VTK  9.1.0
vtkOctreePointLocator.h
Go to the documentation of this file.
1/*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkOctreePointLocator.h
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14=========================================================================*/
15/*----------------------------------------------------------------------------
16 Copyright (c) Sandia Corporation
17 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
18----------------------------------------------------------------------------*/
19
37#ifndef vtkOctreePointLocator_h
38#define vtkOctreePointLocator_h
39
41#include "vtkCommonDataModelModule.h" // For export macro
42
43class vtkCellArray;
44class vtkIdTypeArray;
46class vtkPoints;
47class vtkPolyData;
48
49class VTKCOMMONDATAMODEL_EXPORT vtkOctreePointLocator : public vtkAbstractPointLocator
50{
51public:
53 void PrintSelf(ostream& os, vtkIndent indent) override;
54
56
58
61 vtkSetMacro(MaximumPointsPerRegion, int);
62 vtkGetMacro(MaximumPointsPerRegion, int);
64
66
69 vtkSetMacro(CreateCubicOctants, int);
70 vtkGetMacro(CreateCubicOctants, int);
72
74
80 vtkGetMacro(FudgeFactor, double);
81 vtkSetMacro(FudgeFactor, double);
83
85
89 double* GetBounds() override;
90 void GetBounds(double* bounds) override;
92
94
97 vtkGetMacro(NumberOfLeafNodes, int);
99
103 void GetRegionBounds(int regionID, double bounds[6]);
104
108 void GetRegionDataBounds(int leafNodeID, double bounds[6]);
109
113 int GetRegionContainingPoint(double x, double y, double z);
114
120 void BuildLocator() override;
121
123
127 vtkIdType FindClosestPoint(const double x[3]) override;
128 vtkIdType FindClosestPoint(double x, double y, double z, double& dist2);
130
136 vtkIdType FindClosestPointWithinRadius(double radius, const double x[3], double& dist2) override;
137
139
144 vtkIdType FindClosestPointInRegion(int regionId, double* x, double& dist2);
145 vtkIdType FindClosestPointInRegion(int regionId, double x, double y, double z, double& dist2);
147
152 void FindPointsWithinRadius(double radius, const double x[3], vtkIdList* result) override;
153
162 void FindClosestNPoints(int N, const double x[3], vtkIdList* result) override;
163
168
172 void FreeSearchStructure() override;
173
178 void GenerateRepresentation(int level, vtkPolyData* pd) override;
179
186 void FindPointsInArea(double* area, vtkIdTypeArray* ids, bool clearArray = true);
187
188protected:
191
193 vtkOctreePointLocatorNode** LeafNodeList; // indexed by region/node ID
194
196
198
202 int FindRegion(vtkOctreePointLocatorNode* node, float x, float y, float z);
203 int FindRegion(vtkOctreePointLocatorNode* node, double x, double y, double z);
205
207
209
216 vtkOctreePointLocatorNode* node, double radiusSquared, const double x[3], vtkIdList* ids);
217
218 // Recursive helper for public FindPointsWithinRadius
220
221 // Recursive helper for public FindPointsInArea
223
224 // Recursive helper for public FindPointsInArea
226
227 void DivideRegion(vtkOctreePointLocatorNode* node, int* ordering, int level);
228
229 int DivideTest(int size, int level);
230
232
237 int _FindClosestPointInRegion(int leafNodeId, double x, double y, double z, double& dist2);
238
247 double x, double y, double z, double radius, int skipRegion, double& dist2);
248
250
256
257 double FudgeFactor; // a very small distance, relative to the dataset's size
261
262 float MaxWidth;
263
272
274 void operator=(const vtkOctreePointLocator&) = delete;
275};
276#endif
abstract class to quickly locate points in 3-space
object to represent cell connectivity
Definition: vtkCellArray.h:181
list of point or cell ids
Definition: vtkIdList.h:31
dynamic, self-adjusting array of vtkIdType
a simple class to control print indentation
Definition: vtkIndent.h:34
Octree node that has 8 children each of equal size.
an octree spatial decomposition of a set of points
void AddPolys(vtkOctreePointLocatorNode *node, vtkPoints *pts, vtkCellArray *polys)
void FindPointsWithinRadius(double radius, const double x[3], vtkIdList *result) override
Find all points within a specified radius of position x.
void PrintSelf(ostream &os, vtkIndent indent) override
Standard type and print methods.
int FindClosestPointInSphere(double x, double y, double z, double radius, int skipRegion, double &dist2)
Given a location and a radiues, find the closest point within this radius.
vtkIdType FindClosestPointInRegion(int regionId, double *x, double &dist2)
Find the Id of the point in the given leaf region which is closest to the given point.
static vtkOctreePointLocator * New()
void FindClosestNPoints(int N, const double x[3], vtkIdList *result) override
Find the closest N points to a position.
void FindPointsInArea(vtkOctreePointLocatorNode *node, double *area, vtkIdTypeArray *ids)
void AddAllPointsInRegion(vtkOctreePointLocatorNode *node, vtkIdTypeArray *ids)
vtkOctreePointLocator(const vtkOctreePointLocator &)=delete
int _FindClosestPointInRegion(int leafNodeId, double x, double y, double z, double &dist2)
Given a leaf node id and point, return the local id and the squared distance between the closest poin...
static void DeleteAllDescendants(vtkOctreePointLocatorNode *octant)
int DivideTest(int size, int level)
int NumberOfLeafNodes
The maximum number of points in a region/octant before it is subdivided.
vtkIdType FindClosestPointWithinRadius(double radius, const double x[3], double &dist2) override
Given a position x and a radius r, return the id of the point closest to the point in that radius.
void FindPointsInArea(double *area, vtkIdTypeArray *ids, bool clearArray=true)
Fill ids with points found in area.
double * GetBounds() override
Get the spatial bounds of the entire octree space.
void BuildLeafNodeList(vtkOctreePointLocatorNode *node, int &index)
void AddAllPointsInRegion(vtkOctreePointLocatorNode *node, vtkIdList *ids)
vtkIdType FindClosestPointInRegion(int regionId, double x, double y, double z, double &dist2)
Find the Id of the point in the given leaf region which is closest to the given point.
int FindRegion(vtkOctreePointLocatorNode *node, double x, double y, double z)
Given a point and a node return the leaf node id that contains the point.
void GetRegionDataBounds(int leafNodeID, double bounds[6])
Get the bounds of the data within the leaf node.
int FindRegion(vtkOctreePointLocatorNode *node, float x, float y, float z)
Given a point and a node return the leaf node id that contains the point.
static void SetDataBoundsToSpatialBounds(vtkOctreePointLocatorNode *node)
vtkIdType FindClosestPoint(double x, double y, double z, double &dist2)
Return the Id of the point that is closest to the given point.
vtkOctreePointLocatorNode * Top
~vtkOctreePointLocator() override
vtkIdType FindClosestPoint(const double x[3]) override
Return the Id of the point that is closest to the given point.
void GenerateRepresentation(int level, vtkPolyData *pd) override
Create a polydata representation of the boundaries of the octree regions.
int GetRegionContainingPoint(double x, double y, double z)
Get the id of the leaf region containing the specified location.
int MaximumPointsPerRegion
The maximum number of points in a region/octant before it is subdivided.
void DivideRegion(vtkOctreePointLocatorNode *node, int *ordering, int level)
int CreateCubicOctants
If CreateCubicOctants is non-zero, the bounding box of the points will be expanded such that all octa...
vtkOctreePointLocatorNode ** LeafNodeList
void GetRegionBounds(int regionID, double bounds[6])
Get the spatial bounds of octree region.
void BuildLocator() override
Create the octree decomposition of the cells of the data set or data sets.
void operator=(const vtkOctreePointLocator &)=delete
void GetBounds(double *bounds) override
Get the spatial bounds of the entire octree space.
void FreeSearchStructure() override
Delete the octree data structure.
void FindPointsWithinRadius(vtkOctreePointLocatorNode *node, double radiusSquared, const double x[3], vtkIdList *ids)
Recursive helper for public FindPointsWithinRadius.
vtkIdTypeArray * GetPointsInRegion(int leafNodeId)
Get a list of the original IDs of all points in a leaf node.
represent and manipulate 3D points
Definition: vtkPoints.h:34
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:86
@ level
Definition: vtkX3D.h:401
@ radius
Definition: vtkX3D.h:258
@ size
Definition: vtkX3D.h:259
@ index
Definition: vtkX3D.h:252
int vtkIdType
Definition: vtkType.h:332
#define VTK_NEWINSTANCE