Commit f8f36720 authored by Will Schroeder's avatar Will Schroeder

Fast 2D clipper for image data

The filter takes a single execution to process N contour discrete
contour values. The input is an image label map (i.e., the scalar
values are assumed non-continuous). Threading execution based on
using a flying edges approach is used (i.e., processing dyads at
pixel origins). The output is convex polygons plus optional cell
scalar data which corresponds to the labels from the input image.
parent 3b573f34
......@@ -30,6 +30,7 @@ set(Module_SRCS
vtkDicer.cxx
vtkDiscreteFlyingEdges2D.cxx
vtkDiscreteFlyingEdges3D.cxx
vtkDiscreteFlyingEdgesClipper2D.cxx
vtkDiscreteMarchingCubes.cxx
vtkEdgePoints.cxx
vtkExtractSelectedFrustum.cxx
......
......@@ -9,8 +9,10 @@ vtk_add_test_python(
TestDeformPointSet.py
TestDiscreteFlyingEdges2D.py
TestDiscreteFlyingEdges3D.py
TestDiscreteFlyingEdgesClipper2D.py
TestDiscreteMarchingCubes.py
TestDiscreteMarchingCubesAdjacentScalars.py
TestFEDiscreteClipper2D.py
TestGraphLayoutFilter.py
TestMultiBlockFromTimeSeries.py
TestPointConnectivityFilter.py
......
#!/usr/bin/env python
import vtk
from vtk.util.misc import vtkGetDataRoot
VTK_DATA_ROOT = vtkGetDataRoot()
# Create the RenderWindow, Renderer and both Actors
#
ren1 = vtk.vtkRenderer()
ren1.SetBackground(1,1,1)
renWin = vtk.vtkRenderWindow()
renWin.AddRenderer(ren1)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Read the data. Note this creates a 3-component scalar.
red = vtk.vtkPNGReader()
red.SetFileName(VTK_DATA_ROOT + "/Data/RedCircle.png")
red.Update()
# Next filter can only handle RGB *(&&*@
extract = vtk.vtkImageExtractComponents()
extract.SetInputConnection(red.GetOutputPort())
extract.SetComponents(0,1,2)
# Quantize the image into an index
quantize = vtk.vtkImageQuantizeRGBToIndex()
quantize.SetInputConnection(extract.GetOutputPort())
quantize.SetNumberOfColors(3)
# Create the pipeline
discrete = vtk.vtkDiscreteFlyingEdgesClipper2D()
discrete.SetInputConnection(quantize.GetOutputPort())
discrete.SetValue(0,1)
# Polygons are displayed
polyMapper = vtk.vtkPolyDataMapper()
polyMapper.SetInputConnection(discrete.GetOutputPort())
polyActor = vtk.vtkActor()
polyActor.SetMapper(polyMapper)
ren1.AddActor(polyActor)
renWin.Render()
iren.Start()
#!/usr/bin/env python
import vtk
from vtk.util.misc import vtkGetDataRoot
VTK_DATA_ROOT = vtkGetDataRoot()
# Create the RenderWindow, Renderer and both Actors
#
ren1 = vtk.vtkRenderer()
ren1.SetBackground(1,1,1)
renWin = vtk.vtkRenderWindow()
renWin.AddRenderer(ren1)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Create synthetic image data
VTK_SHORT = 4
img = vtk.vtkImageData()
img.SetDimensions(6,5,1)
img.AllocateScalars(VTK_SHORT,1)
scalars = img.GetPointData().GetScalars()
scalars.SetTuple1(0,0)
scalars.SetTuple1(1,0)
scalars.SetTuple1(2,0)
scalars.SetTuple1(3,0)
scalars.SetTuple1(4,0)
scalars.SetTuple1(5,0)
scalars.SetTuple1(6,0)
scalars.SetTuple1(7,0)
scalars.SetTuple1(8,0)
scalars.SetTuple1(9,0)
scalars.SetTuple1(10,0)
scalars.SetTuple1(11,0)
scalars.SetTuple1(12,0)
scalars.SetTuple1(13,0)
scalars.SetTuple1(14,0)
scalars.SetTuple1(15,2)
scalars.SetTuple1(16,4)
scalars.SetTuple1(17,0)
scalars.SetTuple1(18,0)
scalars.SetTuple1(19,0)
scalars.SetTuple1(20,1)
scalars.SetTuple1(21,1)
scalars.SetTuple1(22,3)
scalars.SetTuple1(23,3)
scalars.SetTuple1(24,0)
scalars.SetTuple1(25,0)
scalars.SetTuple1(26,3)
scalars.SetTuple1(27,0)
scalars.SetTuple1(28,0)
scalars.SetTuple1(29,3)
# Create the pipeline, extract some regions
discrete = vtk.vtkDiscreteFlyingEdgesClipper2D()
discrete.SetInputData(img)
discrete.SetValue(0,1)
discrete.SetValue(1,2)
discrete.SetValue(2,3)
discrete.SetValue(3,4)
discrete.Update()
# Clipped polygons are generated
mapper = vtk.vtkPolyDataMapper()
mapper.SetInputConnection(discrete.GetOutputPort())
mapper.SetScalarModeToUseCellData()
mapper.SetScalarRange(1,4);
actor = vtk.vtkActor()
actor.SetMapper(mapper)
# The image gridlines
gridMapper = vtk.vtkDataSetMapper()
gridMapper.SetInputData(img)
gridMapper.ScalarVisibilityOff()
gridActor = vtk.vtkActor()
gridActor.SetMapper(gridMapper)
gridActor.GetProperty().SetRepresentationToWireframe()
gridActor.GetProperty().SetColor(0,0,1)
ren1.AddActor(actor)
ren1.AddActor(gridActor)
renWin.Render()
iren.Start()
This diff is collapsed.
/*=========================================================================
Program: Visualization Toolkit
Module: vtkDiscreteFlyingEdgesClipper2D.h
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notice for more information.
=========================================================================*/
/**
* @class vtkDiscreteFlyingEdgesClipper2D
* @brief generate filled regions from segmented 2D image data
*
* vtkDiscreteFlyingEdgesClipper2D creates filled polygons from a label map
* (e.g., segmented image) using a variation of the flying edges algorithm
* adapted for 2D clipping. The input is a 2D image where each pixel is
* labeled (integer labels are preferred to real values), and the output data
* is polygonal data representing labeled regions. (Note that on output each
* region [corresponding to a different contour value] may share points on a
* shared boundary.)
*
* While this filter is similar to a contouring operation, label maps do not
* provide continuous function values meaning that usual interpolation along
* edges is not possible. Instead, when the edge endpoints are labeled in
* differing regions, the edge is split at its midpoint. In addition, besides
* producing intersection points at the mid-point of edges, the filter may
* also generate points interior to the pixel cells. For example, if the four
* vertices of a pixel cell are labeled with different regions, then an
* interior point is created and four rectangular "regions" are produced.
*
* Note that one nice feature of this filter is that algorithm execution
* occurs only one time no matter the number of contour values. In many
* contouring-like algorithms, each separate contour value requires an
* additional algorithm execution with a new contour value. So in this filter
* large numbers of contour values do not significantly affect overall speed.
*
* @warning This filter is specialized to 2D images.
*
* @warning
* This class has been threaded with vtkSMPTools. Using TBB or other
* non-sequential type (set in the CMake variable
* VTK_SMP_IMPLEMENTATION_TYPE) may improve performance significantly.
*
* @sa
* vtkDiscreteFlyingEdges2D vtkDiscreteMarchingCubes vtkContourLoopExtraction
* vtkFlyingEdges2D vtkFlyingEdges3D
*/
#ifndef vtkDiscreteFlyingEdgesClipper2D_h
#define vtkDiscreteFlyingEdgesClipper2D_h
#include "vtkFiltersGeneralModule.h" // For export macro
#include "vtkPolyDataAlgorithm.h"
#include "vtkContourValues.h" // Needed for direct access to ContourValues
class vtkImageData;
class VTKFILTERSGENERAL_EXPORT vtkDiscreteFlyingEdgesClipper2D : public vtkPolyDataAlgorithm
{
public:
//@{
/**
* Standard methods for instantiation, printing, and type information.
*/
static vtkDiscreteFlyingEdgesClipper2D *New();
vtkTypeMacro(vtkDiscreteFlyingEdgesClipper2D,vtkPolyDataAlgorithm);
void PrintSelf(ostream& os, vtkIndent indent) override;
//@}
/**
* The modified time is a function of the contour values because we delegate to
* vtkContourValues.
*/
vtkMTimeType GetMTime() override;
/**
* Set a particular contour value at contour number i. The index i ranges
* between 0 <= i <NumberOfContours. (Note: while contour values are
* expressed as doubles, the underlying scalar data may be a different
* type. During execution the contour values are static cast to the type of
* the scalar values.)
*/
void SetValue(int i, double value)
{this->ContourValues->SetValue(i,value);}
/**
* Get the ith contour value.
*/
double GetValue(int i)
{return this->ContourValues->GetValue(i);}
/**
* Get a pointer to an array of contour values. There will be
* GetNumberOfContours() values in the list.
*/
double *GetValues()
{return this->ContourValues->GetValues();}
/**
* Fill a supplied list with contour values. There will be
* GetNumberOfContours() values in the list. Make sure you allocate
* enough memory to hold the list.
*/
void GetValues(double *contourValues)
{this->ContourValues->GetValues(contourValues);}
/**
* Set the number of contours to place into the list. You only really
* need to use this method to reduce list size. The method SetValue()
* will automatically increase list size as needed.
*/
void SetNumberOfContours(int number)
{this->ContourValues->SetNumberOfContours(number);}
/**
* Get the number of contours in the list of contour values.
*/
int GetNumberOfContours()
{return this->ContourValues->GetNumberOfContours();}
//@{
/**
* Generate numContours equally spaced contour values between the specified
* range. Contour values will include min/max range values.
*/
void GenerateValues(int numContours, double range[2])
{this->ContourValues->GenerateValues(numContours, range);}
void GenerateValues(int numContours, double rangeStart, double rangeEnd)
{this->ContourValues->GenerateValues(numContours, rangeStart, rangeEnd);}
//@}
//@{
/**
* Option to set the cell scalars of the output. The scalars will be the
* contour values. By default this flag is on.
*/
vtkSetMacro(ComputeScalars,int);
vtkGetMacro(ComputeScalars,int);
vtkBooleanMacro(ComputeScalars,int);
//@}
//@{
/**
* Set/get which component of a multi-component scalar array to contour on;
* defaults to 0.
*/
vtkSetMacro(ArrayComponent, int);
vtkGetMacro(ArrayComponent, int);
//@}
protected:
vtkDiscreteFlyingEdgesClipper2D();
~vtkDiscreteFlyingEdgesClipper2D() override;
int RequestData(vtkInformation *, vtkInformationVector **,
vtkInformationVector *) override;
int FillInputPortInformation(int port, vtkInformation *info) override;
vtkContourValues *ContourValues;
int ComputeScalars;
int ArrayComponent;
private:
vtkDiscreteFlyingEdgesClipper2D(const vtkDiscreteFlyingEdgesClipper2D&) = delete;
void operator=(const vtkDiscreteFlyingEdgesClipper2D&) = delete;
};
#endif
......@@ -386,4 +386,3 @@ void vtkSelectEnclosedPoints::PrintSelf(ostream& os, vtkIndent indent)
os << indent << "Tolerance: " << this->Tolerance << "\n";
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment