VTK  9.3.0
vtkPolyhedron.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
132 #ifndef vtkPolyhedron_h
133 #define vtkPolyhedron_h
134 
135 #include "vtkCell3D.h"
136 #include "vtkCommonDataModelModule.h" // For export macro
137 
138 VTK_ABI_NAMESPACE_BEGIN
139 class vtkIdTypeArray;
140 class vtkCellArray;
141 class vtkTriangle;
142 class vtkQuad;
143 class vtkTetra;
144 class vtkPolygon;
145 class vtkLine;
146 class vtkIdToIdVectorMapType;
147 class vtkIdToIdMapType;
148 class vtkEdgeTable;
149 class vtkPolyData;
150 class vtkCellLocator;
151 class vtkGenericCell;
152 class vtkPointLocator;
153 
154 class VTKCOMMONDATAMODEL_EXPORT vtkPolyhedron : public vtkCell3D
155 {
156 public:
157  typedef std::map<vtkIdType, vtkIdType> vtkPointIdMap;
158 
160 
163  static vtkPolyhedron* New();
164  vtkTypeMacro(vtkPolyhedron, vtkCell3D);
165  void PrintSelf(ostream& os, vtkIndent indent) override;
167 
169 
173  void GetEdgePoints(vtkIdType vtkNotUsed(edgeId), const vtkIdType*& vtkNotUsed(pts)) override
174  {
175  vtkWarningMacro(<< "vtkPolyhedron::GetEdgePoints Not Implemented");
176  }
177  vtkIdType GetFacePoints(vtkIdType vtkNotUsed(faceId), const vtkIdType*& vtkNotUsed(pts)) override
178  {
179  vtkWarningMacro(<< "vtkPolyhedron::GetFacePoints Not Implemented");
180  return 0;
181  }
183  vtkIdType vtkNotUsed(edgeId), const vtkIdType*& vtkNotUsed(pts)) override
184  {
185  vtkWarningMacro(<< "vtkPolyhedron::GetEdgeToAdjacentFaces Not Implemented");
186  }
188  vtkIdType vtkNotUsed(faceId), const vtkIdType*& vtkNotUsed(faceIds)) override
189  {
190  vtkWarningMacro(<< "vtkPolyhedron::GetFaceToAdjacentFaces Not Implemented");
191  return 0;
192  }
194  vtkIdType vtkNotUsed(pointId), const vtkIdType*& vtkNotUsed(edgeIds)) override
195  {
196  vtkWarningMacro(<< "vtkPolyhedron::GetPointToIncidentEdges Not Implemented");
197  return 0;
198  }
199  vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType*& faceIds) override;
201  vtkIdType vtkNotUsed(pointId), const vtkIdType*& vtkNotUsed(pts)) override
202  {
203  vtkWarningMacro(<< "vtkPolyhedron::GetPointToOneRingPoints Not Implemented");
204  return 0;
205  }
206  bool GetCentroid(double vtkNotUsed(centroid)[3]) const override
207  {
208  vtkWarningMacro(<< "vtkPolyhedron::GetCentroid Not Implemented");
209  return false;
210  }
212 
216  double* GetParametricCoords() override;
217 
221  int GetCellType() override { return VTK_POLYHEDRON; }
222 
226  int RequiresInitialization() override { return 1; }
227 
233  void Initialize() override;
234 
236 
240  int GetNumberOfEdges() override;
241  vtkCell* GetEdge(int) override;
242  int GetNumberOfFaces() override;
243  vtkCell* GetFace(int faceId) override;
245 
251  void Contour(double value, vtkDataArray* scalars, vtkIncrementalPointLocator* locator,
252  vtkCellArray* verts, vtkCellArray* lines, vtkCellArray* polys, vtkPointData* inPd,
253  vtkPointData* outPd, vtkCellData* inCd, vtkIdType cellId, vtkCellData* outCd) override;
254 
264  void Clip(double value, vtkDataArray* scalars, vtkIncrementalPointLocator* locator,
265  vtkCellArray* connectivity, vtkPointData* inPd, vtkPointData* outPd, vtkCellData* inCd,
266  vtkIdType cellId, vtkCellData* outCd, int insideOut) override;
267 
275  int EvaluatePosition(const double x[3], double closestPoint[3], int& subId, double pcoords[3],
276  double& dist2, double weights[]) override;
277 
282  void EvaluateLocation(int& subId, const double pcoords[3], double x[3], double* weights) override;
283 
290  int IntersectWithLine(const double p1[3], const double p2[3], double tol, double& t, double x[3],
291  double pcoords[3], int& subId) override;
292 
308  int Triangulate(int index, vtkIdList* ptIds, vtkPoints* pts) override;
309 
316  int TriangulateFaces(vtkIdList* newFaces);
317 
326  int subId, const double pcoords[3], const double* values, int dim, double* derivs) override;
327 
332  int CellBoundary(int subId, const double pcoords[3], vtkIdList* pts) override;
333 
338  int GetParametricCenter(double pcoords[3]) override;
339 
343  int IsPrimaryCell() override { return 1; }
344 
346 
351  void InterpolateFunctions(const double x[3], double* sf) override;
352  void InterpolateDerivs(const double x[3], double* derivs) override;
354 
360  int RequiresExplicitFaceRepresentation() override { return 1; }
361 
378  void SetFaces(vtkIdType* faces) override;
379 
396  vtkIdType* GetFaces() override;
397 
404  int IsInside(const double x[3], double tolerance);
405 
412  bool IsConvex();
413 
418 
419 protected:
421  ~vtkPolyhedron() override;
422 
423  // Internal classes for supporting operations on this cell
429 
430  // Filled with the SetFaces method.
431  // These faces are numbered in global id space
432  // (in the legacy vtkCellArray form)
434 
435  // Filled with the SetFaces method.
436  // Used to to point to each face in the GlobalFaces array.
438 
439  // vtkCell has the data members Points (x,y,z coordinates) and PointIds (global cell ids).
440  // These data members are implicitly organized in canonical space, i.e., where the cell
441  // point ids are (0,1,...,npts-1).
442  // The PointIdMap is constructed during the call of the Initialize() method and maps global
443  // point ids to the canonical point ids.
445 
446  // If edges are needed. Note that the edge numbering is in canonical space.
447  int EdgesGenerated; // true/false
448  vtkEdgeTable* EdgeTable; // keep track of all edges
449  vtkIdTypeArray* Edges; // edge pairs kept in this list, in canonical id space
450  vtkIdTypeArray* EdgeFaces; // face pairs that comprise each edge, with the
451  // same ordering as EdgeTable
452  int GenerateEdges(); // method populates the edge table and edge array
453 
454  // Numerous methods needs faces to be numbered in the canonical space.
455  // This method uses PointIdMap to fill the Faces member (faces described
456  // with canonical IDs) from the GlobalFaces member (faces described with
457  // global IDs).
459  vtkIdTypeArray* Faces; // These are numbered in canonical id space
460  int FacesGenerated; // True when Faces have been successfully constructed
461 
462  // Bounds management
465  void ComputeParametricCoordinate(const double x[3], double pc[3]);
466  void ComputePositionFromParametricCoordinate(const double pc[3], double x[3]);
467 
469 
470  // Members for supporting geometric operations
480 
481  // Members used in GetPointToIncidentFaces
484 
485 private:
486  vtkPolyhedron(const vtkPolyhedron&) = delete;
487  void operator=(const vtkPolyhedron&) = delete;
488 
490 };
491 
492 //----------------------------------------------------------------------------
493 inline int vtkPolyhedron::GetParametricCenter(double pcoords[3])
494 {
495  pcoords[0] = pcoords[1] = pcoords[2] = 0.5;
496  return 0;
497 }
498 
499 VTK_ABI_NAMESPACE_END
500 #endif
abstract class to specify 3D cell interface
Definition: vtkCell3D.h:28
object to represent cell connectivity
Definition: vtkCellArray.h:185
represent and manipulate cell attribute data
Definition: vtkCellData.h:40
octree-based spatial search object to quickly locate cells
abstract class to specify cell behavior
Definition: vtkCell.h:59
virtual int GetParametricCenter(double pcoords[3])
Return center of the cell in parametric coordinates.
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:54
keep track of edges (edge is pair of integer id's)
Definition: vtkEdgeTable.h:30
provides thread-safe access to cells
list of point or cell ids
Definition: vtkIdList.h:32
dynamic, self-adjusting array of vtkIdType
Abstract class in support of both point location and point insertion.
a simple class to control print indentation
Definition: vtkIndent.h:38
cell represents a 1D line
Definition: vtkLine.h:32
represent and manipulate point attribute data
Definition: vtkPointData.h:39
quickly locate points in 3-space
represent and manipulate 3D points
Definition: vtkPoints.h:38
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:89
a cell that represents an n-sided polygon
Definition: vtkPolygon.h:41
vtkPolyhedron utilities
A 3D cell defined by a set of polygonal faces.
vtkIdType GetFacePoints(vtkIdType vtkNotUsed(faceId), const vtkIdType *&vtkNotUsed(pts)) override
See vtkCell3D API for description of these methods.
vtkIdType GetPointToIncidentFaces(vtkIdType pointId, const vtkIdType *&faceIds) override
See vtkCell3D API for description of these methods.
vtkPointIdMap * PointIdMap
int RequiresExplicitFaceRepresentation() override
Satisfy the vtkCell API.
vtkIdList * CellIds
vtkPolyData * GetPolyData()
Construct polydata if no one exist, then return this->PolyData.
int Triangulate(int index, vtkIdList *ptIds, vtkPoints *pts) override
Use vtkOrderedTriangulator to tetrahedralize the polyhedron mesh.
int GenerateEdges()
int GetNumberOfFaces() override
A polyhedron is represented internally by a set of polygonal faces.
void SetFaces(vtkIdType *faces) override
Set the faces of the polyhedron.
int EvaluatePosition(const double x[3], double closestPoint[3], int &subId, double pcoords[3], double &dist2, double weights[]) override
Satisfy the vtkCell API.
int IntersectWithLine(const double p1[3], const double p2[3], double tol, double &t, double x[3], double pcoords[3], int &subId) override
Intersect the line (p1,p2) with a given tolerance tol to determine a point of intersection x[3] with ...
~vtkPolyhedron() override
void GenerateFaces()
vtkIdType GetPointToIncidentEdges(vtkIdType vtkNotUsed(pointId), const vtkIdType *&vtkNotUsed(edgeIds)) override
See vtkCell3D API for description of these methods.
vtkQuad * Quad
void ComputeParametricCoordinate(const double x[3], double pc[3])
vtkIdType * ValenceAtPoint
vtkEdgeTable * EdgeTable
int TriangulateFaces(vtkIdList *newFaces)
Triangulate each face of the polyhedron.
void ComputeBounds()
void GetEdgeToAdjacentFaces(vtkIdType vtkNotUsed(edgeId), const vtkIdType *&vtkNotUsed(pts)) override
See vtkCell3D API for description of these methods.
vtkIdType GetPointToOneRingPoints(vtkIdType vtkNotUsed(pointId), const vtkIdType *&vtkNotUsed(pts)) override
See vtkCell3D API for description of these methods.
vtkIdType GetFaceToAdjacentFaces(vtkIdType vtkNotUsed(faceId), const vtkIdType *&vtkNotUsed(faceIds)) override
See vtkCell3D API for description of these methods.
void GetEdgePoints(vtkIdType vtkNotUsed(edgeId), const vtkIdType *&vtkNotUsed(pts)) override
See vtkCell3D API for description of these methods.
vtkPolygon * Polygon
void EvaluateLocation(int &subId, const double pcoords[3], double x[3], double *weights) override
The inverse of EvaluatePosition.
vtkCell * GetFace(int faceId) override
A polyhedron is represented internally by a set of polygonal faces.
static vtkPolyhedron * New()
Standard new methods.
void InterpolateFunctions(const double x[3], double *sf) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives).
vtkCellLocator * CellLocator
void Contour(double value, vtkDataArray *scalars, vtkIncrementalPointLocator *locator, vtkCellArray *verts, vtkCellArray *lines, vtkCellArray *polys, vtkPointData *inPd, vtkPointData *outPd, vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd) override
Satisfy the vtkCell API.
int IsInside(const double x[3], double tolerance)
A method particular to vtkPolyhedron.
vtkIdTypeArray * EdgeFaces
vtkTriangle * Triangle
int GetParametricCenter(double pcoords[3]) override
Return the center of the cell in parametric coordinates.
std::map< vtkIdType, vtkIdType > vtkPointIdMap
vtkIdTypeArray * GlobalFaces
vtkIdTypeArray * Edges
void InterpolateDerivs(const double x[3], double *derivs) override
Compute the interpolation functions/derivatives (aka shape functions/derivatives).
bool IsConvex()
Determine whether or not a polyhedron is convex.
void ConstructPolyData()
vtkIdType ** PointToIncidentFaces
int GetNumberOfEdges() override
A polyhedron is represented internally by a set of polygonal faces.
vtkCell * GetEdge(int) override
A polyhedron is represented internally by a set of polygonal faces.
void PrintSelf(ostream &os, vtkIndent indent) override
Standard new methods.
void GeneratePointToIncidentFacesAndValenceAtPoint()
double * GetParametricCoords() override
See vtkCell3D API for description of this method.
vtkIdTypeArray * Faces
vtkPolyData * PolyData
int IsPrimaryCell() override
A polyhedron is a full-fledged primary cell.
vtkLine * Line
vtkTetra * Tetra
void ConstructLocator()
void Derivatives(int subId, const double pcoords[3], const double *values, int dim, double *derivs) override
Computes derivatives at the point specified by the parameter coordinate.
bool GetCentroid(double vtkNotUsed(centroid)[3]) const override
See vtkCell3D API for description of these methods.
int RequiresInitialization() override
This cell requires that it be initialized prior to access.
vtkGenericCell * Cell
void ComputePositionFromParametricCoordinate(const double pc[3], double x[3])
vtkCellArray * Polys
vtkIdType * GetFaces() override
Get the faces of the polyhedron.
int GetCellType() override
See the vtkCell API for descriptions of these methods.
int CellBoundary(int subId, const double pcoords[3], vtkIdList *pts) override
Find the boundary face closest to the point defined by the pcoords[3] and subId of the cell (subId ca...
vtkIdTypeArray * FaceLocations
void Clip(double value, vtkDataArray *scalars, vtkIncrementalPointLocator *locator, vtkCellArray *connectivity, vtkPointData *inPd, vtkPointData *outPd, vtkCellData *inCd, vtkIdType cellId, vtkCellData *outCd, int insideOut) override
Satisfy the vtkCell API.
void Initialize() override
The Initialize method builds up internal structures of vtkPolyhedron.
a cell that represents a 2D quadrilateral
Definition: vtkQuad.h:37
a 3D cell that represents a tetrahedron
Definition: vtkTetra.h:43
a cell that represents a triangle
Definition: vtkTriangle.h:37
@ value
Definition: vtkX3D.h:220
@ index
Definition: vtkX3D.h:246
@ VTK_POLYHEDRON
Definition: vtkCellType.h:89
int vtkIdType
Definition: vtkType.h:315