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
118#ifndef vtkDelaunay2D_h
119#define vtkDelaunay2D_h
120
121#include "vtkAbstractTransform.h" // For point transformation
122#include "vtkFiltersCoreModule.h" // For export macro
123#include "vtkPolyDataAlgorithm.h"
124
125VTK_ABI_NAMESPACE_BEGIN
126class vtkCellArray;
127class vtkIdList;
128class vtkPointSet;
129
130#define VTK_DELAUNAY_XY_PLANE 0
131#define VTK_SET_TRANSFORM_PLANE 1
132#define VTK_BEST_FITTING_PLANE 2
133
134class VTKFILTERSCORE_EXPORT vtkDelaunay2D : public vtkPolyDataAlgorithm
135{
136public:
138 void PrintSelf(ostream& os, vtkIndent indent) override;
139
145
156
166
171
173
179 vtkSetClampMacro(Alpha, double, 0.0, VTK_DOUBLE_MAX);
180 vtkGetMacro(Alpha, double);
182
184
189 vtkSetClampMacro(Tolerance, double, 0.0, 1.0);
190 vtkGetMacro(Tolerance, double);
192
194
198 vtkSetClampMacro(Offset, double, 0.75, VTK_DOUBLE_MAX);
199 vtkGetMacro(Offset, double);
201
203
209 vtkSetMacro(BoundingTriangulation, vtkTypeBool);
210 vtkGetMacro(BoundingTriangulation, vtkTypeBool);
211 vtkBooleanMacro(BoundingTriangulation, vtkTypeBool);
213
215
228
230
238 vtkSetClampMacro(ProjectionPlaneMode, int, VTK_DELAUNAY_XY_PLANE, VTK_BEST_FITTING_PLANE);
239 vtkGetMacro(ProjectionPlaneMode, int);
241
249
251
256 vtkSetMacro(RandomPointInsertion, vtkTypeBool);
257 vtkGetMacro(RandomPointInsertion, vtkTypeBool);
258 vtkBooleanMacro(RandomPointInsertion, vtkTypeBool);
260
261protected:
263
265
266 double Alpha;
267 double Tolerance;
269 double Offset;
271
272 // Transform input points (if necessary)
274
275 int ProjectionPlaneMode; // selects the plane in 3D where the Delaunay triangulation will be
276 // computed.
277
278private:
279 vtkSmartPointer<vtkPolyData> Mesh; // the created mesh
280
281 // the raw points in double precision, and methods to access them
282 double* Points;
283 void SetPoint(vtkIdType id, double* x)
284 {
285 vtkIdType idx = 3 * id;
286 this->Points[idx] = x[0];
287 this->Points[idx + 1] = x[1];
288 this->Points[idx + 2] = x[2];
289 }
290 void GetPoint(vtkIdType id, double x[3])
291 {
292 double* ptr = this->Points + 3 * id;
293 x[0] = *ptr++;
294 x[1] = *ptr++;
295 x[2] = *ptr;
296 }
297
298 // Keep track of the bounding radius of all points (including the eight bounding points).
299 // This is used occasionally for numerical sanity checks to determine whether a point is
300 // within a circumcircle.
301 double BoundingRadius2;
302
303 int NumberOfDuplicatePoints;
304 int NumberOfDegeneracies;
305
306 // Various methods to support the Delaunay algorithm
307 int* RecoverBoundary(vtkPolyData* source);
308 int RecoverEdge(vtkPolyData* source, vtkIdType p1, vtkIdType p2);
309 void FillPolygons(vtkCellArray* polys, int* triUse);
310
311 int InCircle(double x[3], double x1[3], double x2[3], double x3[3]);
312 vtkIdType FindTriangle(double x[3], vtkIdType ptIds[3], vtkIdType tri, double tol,
313 vtkIdType nei[3], vtkIdList* neighbors);
314
315 // CheckEdge() is a recursive function to determine if triangles satisfy the Delaunay
316 // criterion. To prevent segfaults due to excessive recursion, recursion depth is limited.
317 bool CheckEdge(vtkIdType ptId, double x[3], vtkIdType p1, vtkIdType p2, vtkIdType tri,
318 bool recursive, unsigned int depth);
319
320 int FillInputPortInformation(int, vtkInformation*) override;
321
322 vtkDelaunay2D(const vtkDelaunay2D&) = delete;
323 void operator=(const vtkDelaunay2D&) = delete;
324};
325
326VTK_ABI_NAMESPACE_END
327#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
create 2D Delaunay triangulation of input points
vtkPolyData * GetSource()
Get a pointer to the source object.
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.
vtkSmartPointer< vtkAbstractTransform > Transform
static vtkAbstractTransform * ComputeBestFittingPlane(vtkPointSet *input)
This method computes the best fit plane to a set of points represented by a vtkPointSet.
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:23
a simple class to control print indentation
Definition vtkIndent.h:29
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
concrete class for storing a set of points
Definition vtkPointSet.h:59
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:80
Hold a reference to a vtkObjectBase instance.
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