VTK
vtkDaxMarchingCubesImpl.h
Go to the documentation of this file.
1 //=============================================================================
2 //
3 // Copyright (c) Kitware, Inc.
4 // All rights reserved.
5 // See LICENSE.txt for details.
6 //
7 // This software is distributed WITHOUT ANY WARRANTY; without even
8 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
9 // PURPOSE. See the above copyright notice for more information.
10 //
11 // Copyright 2012 Sandia Corporation.
12 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
13 // the U.S. Government retains certain rights in this software.
14 //
15 //=============================================================================
16 
17 #ifndef __vtkDaxMarchingCubesImpl_h
18 #define __vtkDaxMarchingCubesImpl_h
19 
20 // Common code
21 #include "vtkDaxConfig.h"
22 #include "vtkDaxDetailCommon.h"
23 
24 #include "vtkDispatcher.h"
25 #include "vtkDoubleDispatcher.h"
26 #include "vtkNew.h"
27 
28 //fields we support
29 #include "vtkFloatArray.h"
30 
31 //cell types we support
32 #include "vtkCellTypes.h"
33 #include "vtkGenericCell.h"
34 #include "vtkTriangle.h"
35 #include "vtkVoxel.h"
36 
37 //datasets we support
38 #include "vtkDataObjectTypes.h"
39 #include "vtkImageData.h"
40 #include "vtkUniformGrid.h"
41 #include "vtkPolyData.h"
42 
43 //helpers that convert vtk to dax
44 #include "vtkToDax/Portals.h"
45 #include "vtkToDax/Containers.h"
46 #include "vtkToDax/CellTypeToType.h"
47 #include "vtkToDax/DataSetTypeToType.h"
48 #include "vtkToDax/FieldTypeToType.h"
49 #include "vtkToDax/MarchingCubes.h"
50 
51 namespace vtkDax {
52 namespace detail {
54  {
55  typedef int ReturnType;
58  double IsoValue;
59 
61 
63  vtkCell* cell, double isoValue):
64  Input(in),Cell(cell),IsoValue(isoValue),Result(out){}
65 
66  template<typename LHS>
67  int operator()(LHS &arrayField) const
68  {
69  //we can derive the type of the field at compile time, but not the
70  //length
71  if (arrayField.GetNumberOfComponents() == 1)
72  {
73  //first we extract the field type of the array
74  //second we extract the number of components
75  typedef typename vtkToDax::FieldTypeToType<LHS,1>::FieldType FT1;
76  return dispatchOnFieldType<LHS,FT1>(arrayField);
77  }
78  return 0;
79  }
80 
81  template<typename VTKArrayType, typename DaxFieldType>
82  int dispatchOnFieldType(VTKArrayType& vtkField) const
83  {
84  typedef DaxFieldType FieldType;
85  typedef vtkToDax::vtkArrayContainerTag<VTKArrayType> FieldTag;
86  typedef dax::cont::ArrayHandle<FieldType,FieldTag> FieldHandle;
87 
88  typedef typename dax::cont::ArrayHandle
89  <FieldType, FieldTag>::PortalConstControl PortalType;
90 
91  FieldHandle field = FieldHandle( PortalType(&vtkField,
92  vtkField.GetNumberOfTuples() ) );
93  vtkToDax::MarchingCubes<FieldHandle> marching(field,
94  FieldType(IsoValue));
95  marching.setFieldName(vtkField.GetName());
96  marching.setOutputGrid(this->Result);
97 
98  // see if we have a valid data set type if so will perform the
99  // marchingcubes if possible
101  dataDispatcher.Add<vtkImageData,vtkVoxel>(marching);
102  dataDispatcher.Add<vtkUniformGrid,vtkVoxel>(marching);
103 
104  int validMC = dataDispatcher.Go(this->Input,this->Cell);
105  return validMC;
106  }
107  private:
108  void operator=(const ValidMarchingCubesInput&);
109  };
110 } //namespace detail
111 
112 
113 //------------------------------------------------------------------------------
115  vtkDataArray* field, float isoValue)
116 {
117  //we are doing a point threshold now verify we have suitable cells
118  //Dax currently supports: hexs,lines,quads,tets,triangles,vertex,voxel,wedge
119  //if something a cell that doesn't match that list we punt to the
120  //VTK implementation.
122 
123  //construct the object that holds all the state needed to do the MC
124  vtkDax::detail::ValidMarchingCubesInput validInput(input,output,cType.Cell,
125  isoValue);
126 
127 
128  //setup the dispatch to only allow float and int array to go to the next step
129  vtkDispatcher<vtkAbstractArray,int> fieldDispatcher;
130  fieldDispatcher.Add<vtkFloatArray>(validInput);
131  return fieldDispatcher.Go(field);
132 }
133 
134 } //end vtkDax namespace
135 // VTK-HeaderTest-Exclude: vtkDaxMarchingCubesImpl.h
136 #endif
ReturnType Go(BaseLhs *lhs, BaseRhs *rhs)
void Add(Functor fun)
GLenum GLenum GLenum input
Definition: vtkgl.h:15941
abstract class to specify dataset behavior
Definition: vtkDataSet.h:60
dynamic, self-adjusting array of float
Definition: vtkFloatArray.h:45
GLuint in
Definition: vtkgl.h:16905
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:83
void Add(Functor fun)
ValidMarchingCubesInput(vtkDataSet *in, vtkPolyData *out, vtkCell *cell, double isoValue)
abstract class to specify cell behavior
Definition: vtkCell.h:58
Dispatch to functor based on two pointer types.
a cell that represents a 3D orthogonal parallelepiped
Definition: vtkVoxel.h:43
topologically and geometrically regular array of data
Definition: vtkImageData.h:44
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:53
int MarchingCubes(vtkDataSet *input, vtkPolyData *output, vtkDataArray *field, float isoValue)
Dispatch to functor based on a pointer type.
Definition: vtkDispatcher.h:90
CellTypeInDataSet cellType(vtkDataSet *input)
image data with blanking
int dispatchOnFieldType(VTKArrayType &vtkField) const
ReturnType Go(BaseLhs *lhs)