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