VTK  9.3.0
vtkDelaunay2D.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-License-Identifier: BSD-3-Clause
127 #ifndef vtkDelaunay2D_h
128 #define vtkDelaunay2D_h
129 
130 #include "vtkAbstractTransform.h" // For point transformation
131 #include "vtkFiltersCoreModule.h" // For export macro
132 #include "vtkPolyDataAlgorithm.h"
133 
134 VTK_ABI_NAMESPACE_BEGIN
135 class vtkCellArray;
136 class vtkIdList;
137 class vtkPointSet;
138 
139 #define VTK_DELAUNAY_XY_PLANE 0
140 #define VTK_SET_TRANSFORM_PLANE 1
141 #define VTK_BEST_FITTING_PLANE 2
142 
143 class VTKFILTERSCORE_EXPORT vtkDelaunay2D : public vtkPolyDataAlgorithm
144 {
145 public:
147  void PrintSelf(ostream& os, vtkIndent indent) override;
148 
153  static vtkDelaunay2D* New();
154 
165 
175 
180 
182 
188  vtkSetClampMacro(Alpha, double, 0.0, VTK_DOUBLE_MAX);
189  vtkGetMacro(Alpha, double);
191 
193 
198  vtkSetClampMacro(Tolerance, double, 0.0, 1.0);
199  vtkGetMacro(Tolerance, double);
201 
203 
207  vtkSetClampMacro(Offset, double, 0.75, VTK_DOUBLE_MAX);
208  vtkGetMacro(Offset, double);
210 
212 
218  vtkSetMacro(BoundingTriangulation, vtkTypeBool);
219  vtkGetMacro(BoundingTriangulation, vtkTypeBool);
220  vtkBooleanMacro(BoundingTriangulation, vtkTypeBool);
222 
224 
237 
239 
247  vtkSetClampMacro(ProjectionPlaneMode, int, VTK_DELAUNAY_XY_PLANE, VTK_BEST_FITTING_PLANE);
248  vtkGetMacro(ProjectionPlaneMode, int);
250 
258 
260 
265  vtkSetMacro(RandomPointInsertion, vtkTypeBool);
266  vtkGetMacro(RandomPointInsertion, vtkTypeBool);
267  vtkBooleanMacro(RandomPointInsertion, vtkTypeBool);
269 
270 protected:
272 
274 
275  double Alpha;
276  double Tolerance;
278  double Offset;
280 
281  // Transform input points (if necessary)
283 
284  int ProjectionPlaneMode; // selects the plane in 3D where the Delaunay triangulation will be
285  // computed.
286 
287 private:
288  vtkSmartPointer<vtkPolyData> Mesh; // the created mesh
289 
290  // the raw points in double precision, and methods to access them
291  double* Points;
292  void SetPoint(vtkIdType id, double* x)
293  {
294  vtkIdType idx = 3 * id;
295  this->Points[idx] = x[0];
296  this->Points[idx + 1] = x[1];
297  this->Points[idx + 2] = x[2];
298  }
299  void GetPoint(vtkIdType id, double x[3])
300  {
301  double* ptr = this->Points + 3 * id;
302  x[0] = *ptr++;
303  x[1] = *ptr++;
304  x[2] = *ptr;
305  }
306 
307  // Keep track of the bounding radius of all points (including the eight bounding points).
308  // This is used occasionally for numerical sanity checks to determine whether a point is
309  // within a circumcircle.
310  double BoundingRadius2;
311 
312  int NumberOfDuplicatePoints;
313  int NumberOfDegeneracies;
314 
315  // Various methods to support the Delaunay algorithm
316  int* RecoverBoundary(vtkPolyData* source);
317  int RecoverEdge(vtkPolyData* source, vtkIdType p1, vtkIdType p2);
318  void FillPolygons(vtkCellArray* polys, int* triUse);
319 
320  int InCircle(double x[3], double x1[3], double x2[3], double x3[3]);
321  vtkIdType FindTriangle(double x[3], vtkIdType ptIds[3], vtkIdType tri, double tol,
322  vtkIdType nei[3], vtkIdList* neighbors);
323 
324  // CheckEdge() is a recursive function to determine if triangles satisfy the Delaunay
325  // criterion. To prevent segfaults due to excessive recursion, recursion depth is limited.
326  bool CheckEdge(vtkIdType ptId, double x[3], vtkIdType p1, vtkIdType p2, vtkIdType tri,
327  bool recursive, unsigned int depth);
328 
329  int FillInputPortInformation(int, vtkInformation*) override;
330 
331  vtkDelaunay2D(const vtkDelaunay2D&) = delete;
332  void operator=(const vtkDelaunay2D&) = delete;
333 };
334 
335 VTK_ABI_NAMESPACE_END
336 #endif
void GetPoint(int i, int j, int k, double pnt[3])
superclass for all geometric transformations
Proxy object to connect input/output ports.
object to represent cell connectivity
Definition: vtkCellArray.h:185
create 2D Delaunay triangulation of input points
static vtkAbstractTransform * ComputeBestFittingPlane(vtkPointSet *input)
This method computes the best fit plane to a set of points represented by a vtkPointSet.
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
This is called by the superclass.
static vtkDelaunay2D * New()
Construct object with Alpha = 0.0; Tolerance = 0.001; Offset = 1.25; BoundingTriangulation turned off...
void SetSourceConnection(vtkAlgorithmOutput *algOutput)
Specify the source object used to specify constrained edges and loops.
vtkPolyData * GetSource()
Get a pointer to the source object.
vtkSmartPointer< vtkAbstractTransform > Transform
vtkTypeBool BoundingTriangulation
vtkTypeBool RandomPointInsertion
vtkSetSmartPointerMacro(Transform, vtkAbstractTransform)
Set / get the transform which is applied to points to generate a 2D problem.
void SetSourceData(vtkPolyData *)
Specify the source object used to specify constrained edges and loops.
vtkGetSmartPointerMacro(Transform, vtkAbstractTransform)
Set / get the transform which is applied to points to generate a 2D problem.
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
list of point or cell ids
Definition: vtkIdList.h:32
a simple class to control print indentation
Definition: vtkIndent.h:38
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
concrete class for storing a set of points
Definition: vtkPointSet.h:68
Superclass for algorithms that produce only polydata as output.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:89
@ Transform
Definition: vtkX3D.h:41
int vtkTypeBool
Definition: vtkABI.h:64
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
#define VTK_BEST_FITTING_PLANE
#define VTK_DELAUNAY_XY_PLANE
int vtkIdType
Definition: vtkType.h:315
#define VTK_DOUBLE_MAX
Definition: vtkType.h:154