VTK  9.0.1
vtk3DLinearGridInternal.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtk3DLinearGridInternal.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 =========================================================================*/
37 #ifndef vtk3DLinearGridInternal_h
38 #define vtk3DLinearGridInternal_h
39 
40 // Include appropriate cell types
41 #include "vtkCellArray.h"
42 #include "vtkCellArrayIterator.h"
43 #include "vtkHexahedron.h"
44 #include "vtkPyramid.h"
45 #include "vtkSmartPointer.h"
46 #include "vtkTetra.h"
47 #include "vtkVoxel.h"
48 #include "vtkWedge.h"
49 #include <cmath>
50 
51 namespace
52 { // anonymous namespace
53 
54 //========================= CELL MACHINARY ====================================
55 
56 // Implementation note: this filter currently handles 3D linear cells. It
57 // could be extended to handle other 3D cell types.
58 
59 // The maximum number of verts per cell (hexahedron)
60 #define MAX_CELL_VERTS 8
61 
62 // Base class to represent cells
63 struct BaseCell
64 {
65  unsigned char CellType;
66  unsigned char NumVerts;
67  unsigned char NumEdges;
68  unsigned short* Cases;
69  static unsigned char Mask[MAX_CELL_VERTS];
70 
71  BaseCell(int cellType)
72  : CellType(cellType)
73  , NumVerts(0)
74  , NumEdges(0)
75  , Cases(nullptr)
76  {
77  }
78  virtual ~BaseCell() {}
79 
80  // Set up the case table. This is done by accessing standard VTK cells and
81  // repackaging the case table for efficiency. The format of the case table
82  // is as follows: a linear array, organized into two parts: 1) offsets into
83  // the second part, and 2) the cases. The first 2^NumVerts entries are the
84  // offsets which refer to the 2^NumVerts cases in the second part. Each
85  // case is represented by the number of edges, followed by pairs of
86  // vertices (v0,v1) for each edge. Note that groups of three contiguous
87  // edges form a triangle.
88  virtual void BuildCases() = 0;
89  void BuildCases(int numCases, const vtkIdType** edges, int** cases, unsigned short* caseArray);
90 };
91 // Used to generate case mask
92 unsigned char BaseCell::Mask[MAX_CELL_VERTS] = { 1, 2, 4, 8, 16, 32, 64, 128 };
93 // Build repackaged case table and place into cases array.
94 void BaseCell::BuildCases(
95  int numCases, const vtkIdType** edges, int** cases, unsigned short* caseArray)
96 {
97  int caseOffset = numCases;
98  for (int caseNum = 0; caseNum < numCases; ++caseNum)
99  {
100  caseArray[caseNum] = caseOffset;
101  int* triCases = cases[caseNum];
102 
103  // Count the number of edges
104  int count;
105  for (count = 0; triCases[count] != (-1); ++count)
106  {
107  }
108  caseArray[caseOffset++] = count;
109 
110  // Now populate the edges
111  const vtkIdType* edge;
112  for (count = 0; triCases[count] != (-1); ++count)
113  {
114  edge = edges[triCases[count]];
115  caseArray[caseOffset++] = edge[0];
116  caseArray[caseOffset++] = edge[1];
117  }
118  } // for all cases
119 }
120 
121 // Contour tetrahedral cell------------------------------------------------------
122 // Repackages case table for more efficient processing.
123 struct TetraCell : public BaseCell
124 {
125  static unsigned short TetraCases[152];
126 
127  TetraCell()
128  : BaseCell(VTK_TETRA)
129  {
130  this->NumVerts = 4;
131  this->NumEdges = 6;
132  this->BuildCases();
133  this->Cases = this->TetraCases;
134  }
135  ~TetraCell() override {}
136  void BuildCases() override;
137 };
138 // Dummy initialization filled in later at initialization. The lengtth of the
139 // array is determined from the equation length=(2*NumCases + 3*2*NumTris).
140 unsigned short TetraCell::TetraCases[152] = { 0 };
141 // Load and transform vtkTetra case table. The case tables are repackaged for
142 // efficiency (e.g., support the GetCase() method).
143 void TetraCell::BuildCases()
144 {
145  const vtkIdType** edges = new const vtkIdType*[this->NumEdges];
146  int numCases = std::pow(2, this->NumVerts);
147  int** cases = new int*[numCases];
148  for (int i = 0; i < this->NumEdges; ++i)
149  {
151  }
152  for (int i = 0; i < numCases; ++i)
153  {
154  cases[i] = vtkTetra::GetTriangleCases(i);
155  }
156 
157  BaseCell::BuildCases(numCases, edges, cases, this->TetraCases);
158 
159  delete[] edges;
160  delete[] cases;
161 }
162 
163 // Contour hexahedral cell------------------------------------------------------
164 struct HexahedronCell : public BaseCell
165 {
166  static unsigned short HexahedronCases[5432];
167 
168  HexahedronCell()
169  : BaseCell(VTK_HEXAHEDRON)
170  {
171  this->NumVerts = 8;
172  this->NumEdges = 12;
173  this->BuildCases();
174  this->Cases = this->HexahedronCases;
175  }
176  ~HexahedronCell() override {}
177  void BuildCases() override;
178 };
179 // Dummy initialization filled in later at instantiation
180 unsigned short HexahedronCell::HexahedronCases[5432] = { 0 };
181 // Load and transform marching cubes case table. The case tables are
182 // repackaged for efficiency (e.g., support the GetCase() method).
183 void HexahedronCell::BuildCases()
184 {
185  const vtkIdType** edges = new const vtkIdType*[this->NumEdges];
186  int numCases = std::pow(2, this->NumVerts);
187  int** cases = new int*[numCases];
188  for (int i = 0; i < this->NumEdges; ++i)
189  {
191  }
192  for (int i = 0; i < numCases; ++i)
193  {
194  cases[i] = vtkHexahedron::GetTriangleCases(i);
195  }
196 
197  BaseCell::BuildCases(numCases, edges, cases, this->HexahedronCases);
198 
199  delete[] edges;
200  delete[] cases;
201 }
202 
203 // Contour wedge cell ------------------------------------------------------
204 struct WedgeCell : public BaseCell
205 {
206  static unsigned short WedgeCases[968];
207 
208  WedgeCell()
209  : BaseCell(VTK_WEDGE)
210  {
211  this->NumVerts = 6;
212  this->NumEdges = 9;
213  this->BuildCases();
214  this->Cases = this->WedgeCases;
215  }
216  ~WedgeCell() override {}
217  void BuildCases() override;
218 };
219 // Dummy initialization filled in later at instantiation
220 unsigned short WedgeCell::WedgeCases[968] = { 0 };
221 // Load and transform marching cubes case table. The case tables are
222 // repackaged for efficiency (e.g., support the GetCase() method).
223 void WedgeCell::BuildCases()
224 {
225  const vtkIdType** edges = new const vtkIdType*[this->NumEdges];
226  int numCases = std::pow(2, this->NumVerts);
227  int** cases = new int*[numCases];
228  for (int i = 0; i < this->NumEdges; ++i)
229  {
231  }
232  for (int i = 0; i < numCases; ++i)
233  {
234  cases[i] = vtkWedge::GetTriangleCases(i);
235  }
236 
237  BaseCell::BuildCases(numCases, edges, cases, this->WedgeCases);
238 
239  delete[] edges;
240  delete[] cases;
241 }
242 
243 // Contour pyramid cell------------------------------------------------------
244 struct PyramidCell : public BaseCell
245 {
246  static unsigned short PyramidCases[448];
247 
248  PyramidCell()
249  : BaseCell(VTK_PYRAMID)
250  {
251  this->NumVerts = 5;
252  this->NumEdges = 8;
253  this->BuildCases();
254  this->Cases = this->PyramidCases;
255  }
256  ~PyramidCell() override {}
257  void BuildCases() override;
258 };
259 // Dummy initialization filled in later at instantiation
260 unsigned short PyramidCell::PyramidCases[448] = { 0 };
261 // Load and transform marching cubes case table. The case tables are
262 // repackaged for efficiency (e.g., support the GetCase() method).
263 void PyramidCell::BuildCases()
264 {
265  const vtkIdType** edges = new const vtkIdType*[this->NumEdges];
266  int numCases = std::pow(2, this->NumVerts);
267  int** cases = new int*[numCases];
268  for (int i = 0; i < this->NumEdges; ++i)
269  {
271  }
272  for (int i = 0; i < numCases; ++i)
273  {
274  cases[i] = vtkPyramid::GetTriangleCases(i);
275  }
276 
277  BaseCell::BuildCases(numCases, edges, cases, this->PyramidCases);
278 
279  delete[] edges;
280  delete[] cases;
281 }
282 
283 // Contour voxel cell------------------------------------------------------
284 struct VoxelCell : public BaseCell
285 {
286  static unsigned short VoxCases[5432];
287 
288  VoxelCell()
289  : BaseCell(VTK_VOXEL)
290  {
291  this->NumVerts = 8;
292  this->NumEdges = 12;
293  this->BuildCases();
294  this->Cases = this->VoxCases;
295  }
296  ~VoxelCell() override {}
297  void BuildCases() override;
298 };
299 // Dummy initialization filled in later at instantiation
300 unsigned short VoxelCell::VoxCases[5432] = { 0 };
301 // Load and transform marching cubes case table. The case tables are
302 // repackaged for efficiency (e.g., support the GetCase() method). Note that
303 // the MC cases (vtkMarchingCubesTriangleCases) are specified for the
304 // hexahedron; voxels require a transformation to produce correct output.
305 void VoxelCell::BuildCases()
306 {
307  // Map the voxel points consistent with the hex edges and cases, Basically
308  // the hex points (2,3,6,7) are ordered (3,2,7,6) on the voxel.
309  const vtkIdType** edges = new const vtkIdType*[this->NumEdges];
310  constexpr vtkIdType voxEdges[12][2] = { { 0, 1 }, { 1, 3 }, { 2, 3 }, { 0, 2 }, { 4, 5 },
311  { 5, 7 }, { 6, 7 }, { 4, 6 }, { 0, 4 }, { 1, 5 }, { 2, 6 }, { 3, 7 } };
312 
313  for (int i = 0; i < this->NumEdges; ++i)
314  {
315  edges[i] = voxEdges[i];
316  }
317 
318  // Build the voxel cases. Have to shuffle them around due to different
319  // vertex ordering.
320  unsigned int numCases = std::pow(2, this->NumVerts);
321  int** cases = new int*[numCases];
322  unsigned int hexCase, voxCase;
323  for (hexCase = 0; hexCase < numCases; ++hexCase)
324  {
325  voxCase = ((hexCase & BaseCell::Mask[0]) ? 1 : 0) << 0;
326  voxCase |= ((hexCase & BaseCell::Mask[1]) ? 1 : 0) << 1;
327  voxCase |= ((hexCase & BaseCell::Mask[2]) ? 1 : 0) << 3;
328  voxCase |= ((hexCase & BaseCell::Mask[3]) ? 1 : 0) << 2;
329  voxCase |= ((hexCase & BaseCell::Mask[4]) ? 1 : 0) << 4;
330  voxCase |= ((hexCase & BaseCell::Mask[5]) ? 1 : 0) << 5;
331  voxCase |= ((hexCase & BaseCell::Mask[6]) ? 1 : 0) << 7;
332  voxCase |= ((hexCase & BaseCell::Mask[7]) ? 1 : 0) << 6;
333  cases[voxCase] = vtkHexahedron::GetTriangleCases(hexCase);
334  }
335 
336  BaseCell::BuildCases(numCases, edges, cases, this->VoxCases);
337 
338  delete[] edges;
339  delete[] cases;
340 }
341 
342 // Contour empty cell. These cells are skipped.---------------------------------
343 struct EmptyCell : public BaseCell
344 {
345  static unsigned short EmptyCases[2];
346 
347  EmptyCell()
348  : BaseCell(VTK_EMPTY_CELL)
349  {
350  this->NumVerts = 0;
351  this->NumEdges = 0;
352  this->Cases = this->EmptyCases;
353  }
354  ~EmptyCell() override {}
355  void BuildCases() override {}
356 };
357 // No triangles generated
358 unsigned short EmptyCell::EmptyCases[2] = { 0, 0 };
359 
360 // This is a general iterator which assumes that the unstructured grid has a
361 // mix of cells. Any cell that is not processed by this contouring algorithm
362 // (i.e., not one of tet, hex, pyr, wedge, voxel) is skipped.
363 struct CellIter
364 {
365  // Current active cell, and whether it is a copy (which controls
366  // the destruction process).
367  bool Copy;
368  BaseCell* Cell;
369 
370  // The iteration state.
371  unsigned char NumVerts;
372  const unsigned short* Cases;
373 
374  // References to unstructured grid for cell traversal.
375  vtkIdType NumCells;
376  const unsigned char* Types;
379 
380  // All possible cell types. The iterator switches between them when
381  // processing. All unsupported cells are of type EmptyCell.
382  TetraCell* Tetra;
383  HexahedronCell* Hexahedron;
384  PyramidCell* Pyramid;
385  WedgeCell* Wedge;
386  VoxelCell* Voxel;
387  EmptyCell* Empty;
388 
389  CellIter()
390  : Copy(true)
391  , Cell(nullptr)
392  , NumVerts(0)
393  , Cases(nullptr)
394  , NumCells(0)
395  , Types(nullptr)
396  , Tetra(nullptr)
397  , Hexahedron(nullptr)
398  , Pyramid(nullptr)
399  , Wedge(nullptr)
400  , Voxel(nullptr)
401  , Empty(nullptr)
402  {
403  }
404 
405  CellIter(vtkIdType numCells, unsigned char* types, vtkCellArray* cellArray)
406  : Copy(false)
407  , Cell(nullptr)
408  , NumVerts(0)
409  , Cases(nullptr)
410  , NumCells(numCells)
411  , Types(types)
412  , CellArray(cellArray)
413  , ConnIter(vtk::TakeSmartPointer(cellArray->NewIterator()))
414  {
415  this->Tetra = new TetraCell;
416  this->Hexahedron = new HexahedronCell;
417  this->Pyramid = new PyramidCell;
418  this->Wedge = new WedgeCell;
419  this->Voxel = new VoxelCell;
420  this->Empty = new EmptyCell;
421  }
422 
423  ~CellIter()
424  {
425  if (!this->Copy)
426  {
427  delete this->Tetra;
428  delete this->Hexahedron;
429  delete this->Pyramid;
430  delete this->Wedge;
431  delete this->Voxel;
432  delete this->Empty;
433  }
434  }
435 
436  CellIter(const CellIter&) = default; // remove compiler warnings
437 
438  // Shallow copy to avoid new/delete.
439  CellIter& operator=(const CellIter& cellIter)
440  {
441  this->Copy = true;
442  this->Cell = nullptr;
443 
444  this->NumVerts = cellIter.NumVerts;
445  this->Cases = cellIter.Cases;
446 
447  this->NumCells = cellIter.NumCells;
448  this->Types = cellIter.Types;
449  this->CellArray = cellIter.CellArray;
450 
451  // This class is passed around by pointer and only copied deliberately
452  // to create thread-local copies. Since we don't want to share state,
453  // create a new iterator here:
454  if (cellIter.ConnIter)
455  {
456  this->ConnIter = vtk::TakeSmartPointer(this->CellArray->NewIterator());
457  this->ConnIter->GoToCell(cellIter.ConnIter->GetCurrentCellId());
458  }
459  else
460  {
461  this->ConnIter = nullptr;
462  }
463 
464  this->Tetra = cellIter.Tetra;
465  this->Hexahedron = cellIter.Hexahedron;
466  this->Pyramid = cellIter.Pyramid;
467  this->Wedge = cellIter.Wedge;
468  this->Voxel = cellIter.Voxel;
469  this->Empty = cellIter.Empty;
470 
471  return *this;
472  }
473 
474  // Decode the case table. (See previous documentation of case table
475  // organization.) Note that bounds/range chacking is not performed
476  // for efficiency.
477  const unsigned short* GetCase(unsigned char caseNum)
478  {
479  return (this->Cases + this->Cases[caseNum]);
480  }
481 
482  // Methods for caching traversal. Initialize() sets up the traversal
483  // process; Next() advances to the next cell. Note that the public data
484  // members representing the iteration state (NumVerts, Cases, ConnIter) are
485  // modified by these methods, and then subsequently read during iteration.
486  const vtkIdType* Initialize(vtkIdType cellId)
487  {
488  this->Cell = this->GetCell(this->Types[cellId]);
489  this->NumVerts = this->Cell->NumVerts;
490  this->Cases = this->Cell->Cases;
491  this->ConnIter->GoToCell(cellId);
492 
493  vtkIdType dummy;
494  const vtkIdType* conn;
495  this->ConnIter->GetCurrentCell(dummy, conn);
496  return conn;
497  }
498 
499  const vtkIdType* Next()
500  {
501  this->ConnIter->GoToNextCell();
502 
503  if (this->ConnIter->IsDoneWithTraversal())
504  {
505  return nullptr;
506  }
507 
508  const vtkIdType currentCellId = this->ConnIter->GetCurrentCellId();
509 
510  // Only update information if the cell type changes. Note however that
511  // empty cells may have to be treated specially.
512  if (this->Cell->CellType == VTK_EMPTY_CELL ||
513  this->Cell->CellType != this->Types[currentCellId])
514  {
515  this->Cell = this->GetCell(this->Types[currentCellId]);
516  this->NumVerts = this->Cell->NumVerts;
517  this->Cases = this->Cell->Cases;
518  }
519 
520  vtkIdType dummy;
521  const vtkIdType* conn;
522  this->ConnIter->GetCurrentCell(dummy, conn);
523  return conn;
524  }
525 
526  // Method for random access of cell, no caching
527  unsigned char GetCellType(vtkIdType cellId) { return this->Types[cellId]; }
528 
529  // Method for random access of cell, no caching
530  const vtkIdType* GetCellIds(vtkIdType cellId)
531  {
532  this->Cell = this->GetCell(this->Types[cellId]);
533  this->NumVerts = this->Cell->NumVerts;
534  this->Cases = this->Cell->Cases;
535  this->ConnIter->GoToCell(cellId);
536 
537  vtkIdType dummy;
538  const vtkIdType* conn;
539  this->ConnIter->GetCurrentCell(dummy, conn);
540  return conn;
541  }
542 
543  // Switch to the appropriate cell type.
544  BaseCell* GetCell(int cellType)
545  {
546  switch (cellType)
547  {
548  case VTK_TETRA:
549  return this->Tetra;
550  case VTK_HEXAHEDRON:
551  return this->Hexahedron;
552  case VTK_WEDGE:
553  return this->Wedge;
554  case VTK_PYRAMID:
555  return this->Pyramid;
556  case VTK_VOXEL:
557  return this->Voxel;
558  default:
559  return this->Empty;
560  }
561  }
562 };
563 
564 } // anonymous namespace
565 
566 #endif // vtk3DLinearGridInternal_h
567 // VTK-HeaderTest-Exclude: vtk3DLinearGridInternal.h
VTK_VOXEL
@ VTK_VOXEL
Definition: vtkCellType.h:57
vtkPyramid::GetEdgeArray
static const vtkIdType * GetEdgeArray(vtkIdType edgeId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
vtkIdType
int vtkIdType
Definition: vtkType.h:338
vtkCellArray.h
vtkSmartPointer< vtkCellArray >
vtkTetra::GetTriangleCases
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
VTK_PYRAMID
@ VTK_PYRAMID
Definition: vtkCellType.h:60
vtkHexahedron::GetTriangleCases
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
vtkHexahedron.h
VTK_EMPTY_CELL
@ VTK_EMPTY_CELL
Definition: vtkCellType.h:46
vtkWedge::GetEdgeArray
static const vtkIdType * GetEdgeArray(vtkIdType edgeId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
vtkCellArrayIterator.h
vtkPyramid.h
vtkCellArray
object to represent cell connectivity
Definition: vtkCellArray.h:179
vtkSmartPointer.h
VTK_HEXAHEDRON
@ VTK_HEXAHEDRON
Definition: vtkCellType.h:58
VTK_TETRA
@ VTK_TETRA
Definition: vtkCellType.h:56
vtkVoxel.h
vtkWedge::GetTriangleCases
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
vtkWedge.h
vtk
Specialization of tuple ranges and iterators for vtkAOSDataArrayTemplate.
Definition: vtkAtomicTypeConcepts.h:21
vtk::TakeSmartPointer
vtkSmartPointer< T > TakeSmartPointer(T *obj)
Construct a vtkSmartPointer<T> containing obj.
Definition: vtkSmartPointer.h:332
vtkTetra::GetEdgeArray
static const vtkIdType * GetEdgeArray(vtkIdType edgeId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
edges
std::pair< boost::graph_traits< vtkGraph * >::edge_iterator, boost::graph_traits< vtkGraph * >::edge_iterator > edges(vtkGraph *g)
Definition: vtkBoostGraphAdapter.h:988
VTK_WEDGE
@ VTK_WEDGE
Definition: vtkCellType.h:59
MAX_CELL_VERTS
#define MAX_CELL_VERTS
Definition: vtk3DLinearGridInternal.h:60
vtkHexahedron::GetEdgeArray
static const vtkIdType * GetEdgeArray(vtkIdType edgeId)
Return the ids of the vertices defining edge/face (edgeId/‘faceId’).
vtkPyramid::GetTriangleCases
static int * GetTriangleCases(int caseId)
Return the case table for table-based isocontouring (aka marching cubes style implementations).
vtkTetra.h