VTK
Contour.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 vtkToDax_Contour_h
18#define vtkToDax_Contour_h
19
20#include "vtkDispatcher.h"
21#include "vtkPolyData.h"
22
27
30
31#include <dax/cont/DispatcherGenerateInterpolatedCells.h>
32#include <dax/cont/DispatcherMapCell.h>
33#include <dax/worklet/MarchingCubes.h>
34#include <dax/worklet/MarchingTetrahedra.h>
35
36namespace
37{
38template <typename T> struct MarchingCubesOuputType
39{
40 typedef dax::CellTagTriangle type;
41};
42
43}
44
45namespace vtkToDax
46{
47
48template<int B>
50{
51 template<class InGridType,
52 class OutGridType,
53 typename ValueType,
54 class Container1,
55 class Adapter>
56 int operator()(const InGridType &,
57 vtkDataSet *,
58 OutGridType &,
60 ValueType,
61 const dax::cont::ArrayHandle<ValueType,Container1,Adapter> &,
62 bool)
63 {
64 vtkGenericWarningMacro(
65 << "Not calling Dax, GridType-CellType combination not supported");
66 return 0;
67 }
68};
69template<>
70struct DoContour<1>
71{
72 template<class InGridType,
73 class OutGridType,
74 typename ValueType,
75 class Container1,
76 class Adapter>
78 const InGridType &inDaxGrid,
79 vtkDataSet *inVTKGrid,
80 OutGridType &outDaxGeom,
81 vtkPolyData *outVTKGrid,
82 ValueType isoValue,
83 const dax::cont::ArrayHandle<ValueType,Container1,Adapter> &mcHandle,
84 bool computeScalars)
85 {
86 dax::Scalar isoValueT(isoValue);
87
88 dax::worklet::MarchingCubesCount countWorklet(isoValueT);
89 dax::worklet::MarchingCubesGenerate generateWorklet(isoValueT);
90
91 return this->DispatchWork(inDaxGrid,
92 inVTKGrid,
93 outDaxGeom,
94 outVTKGrid,
95 countWorklet,
96 generateWorklet,
97 mcHandle,
98 computeScalars);
99 }
100
101 template<class GridCellContainer,
102 class GridPointContainer,
103 class OutGridType,
104 typename ValueType,
105 class Container1,
106 class Adapter>
108 const dax::cont::UnstructuredGrid<
109 dax::CellTagTetrahedron,GridCellContainer,GridPointContainer,Adapter>
110 &inDaxGrid,
111 vtkDataSet *inVTKGrid,
112 OutGridType &outDaxGeom,
113 vtkPolyData *outVTKGrid,
114 ValueType isoValue,
115 const dax::cont::ArrayHandle<ValueType,Container1,Adapter> &mcHandle,
116 bool computeScalars)
117 {
118 dax::Scalar isoValueT(isoValue);
119
120 dax::worklet::MarchingTetrahedraCount countWorklet(isoValueT);
121 dax::worklet::MarchingTetrahedraGenerate generateWorklet(isoValueT);
122
123 return this->DispatchWork(inDaxGrid,
124 inVTKGrid,
125 outDaxGeom,
126 outVTKGrid,
127 countWorklet,
128 generateWorklet,
129 mcHandle,
130 computeScalars);
131 }
132
133 template<class InGridType,
134 class OutGridType,
135 typename ValueType,
136 class Container1,
137 class Adapter,
138 class CountWorkletType,
139 class GenerateWorkletType>
141 const InGridType &inDaxGrid,
142 vtkDataSet *inVTKGrid,
143 OutGridType &outDaxGeom,
144 vtkPolyData *outVTKGrid,
145 CountWorkletType &countWorklet,
146 GenerateWorkletType &generateWorklet,
147 const dax::cont::ArrayHandle<ValueType,Container1,Adapter> &mcHandle,
148 bool computeScalars)
149 {
150 int result=1;
151
152 try
153 {
154
155 typedef dax::cont::DispatcherGenerateInterpolatedCells<
156 GenerateWorkletType,
157 dax::cont::ArrayHandle< dax::Id >,
158 Adapter > DispatchIC;
159
160 typedef typename DispatchIC::CountHandleType CountHandleType;
161
162 dax::cont::DispatcherMapCell<CountWorkletType,Adapter>
163 dispatchCount( countWorklet );
164
165 CountHandleType count;
166 dispatchCount.Invoke(inDaxGrid, mcHandle, count);
167
168
169 DispatchIC generateSurface(count, generateWorklet);
170 generateSurface.SetRemoveDuplicatePoints(true);
171 generateSurface.Invoke(inDaxGrid,outDaxGeom,mcHandle);
172
173 // Convert output geometry to VTK.
174 daxToVtk::dataSetConverter(outDaxGeom, outVTKGrid);
175
176 // Interpolate arrays where possible.
177 if (computeScalars)
178 {
179 vtkToDax::CompactPointField<DispatchIC> compact(generateSurface,
180 outVTKGrid);
181 vtkDispatcher<vtkAbstractArray,int> compactDispatcher;
182 compactDispatcher.Add<vtkFloatArray>(compact);
183 compactDispatcher.Add<vtkDoubleArray>(compact);
184
185 vtkPointData *pd = inVTKGrid->GetPointData();
186 for (int arrayIndex = 0;
187 arrayIndex < pd->GetNumberOfArrays();
188 arrayIndex++)
189 {
190 vtkDataArray *array = pd->GetArray(arrayIndex);
191 if (array == NULL) { continue; }
192
193 compactDispatcher.Go(array);
194 }
195
196 // Pass information about attributes.
197 for (int attributeType = 0;
199 attributeType++)
200 {
201 vtkDataArray *attribute = pd->GetAttribute(attributeType);
202 if (attribute == NULL) { continue; }
203 outVTKGrid->GetPointData()->SetActiveAttribute(attribute->GetName(),
204 attributeType);
205 }
206 } //computeScalars
207 }
208 catch(const dax::cont::ErrorControlOutOfMemory& error)
209 {
210 std::cerr << "Ran out of memory trying to use the GPU" << std::endl;
211 std::cerr << error.GetMessage() << std::endl;
212 result = 0;
213 }
214 catch(const dax::cont::ErrorExecution& error)
215 {
216 std::cerr << "Got ErrorExecution from Dax." << std::endl;
217 std::cerr << error.GetMessage() << std::endl;
218 result = 0;
219 }
220 return result;
221 }
222};
223
224 template<typename FieldType_>
225 struct Contour
226 {
227 public:
228 typedef FieldType_ FieldType;
229 //we expect FieldType_ to be an dax::cont::ArrayHandle
230 typedef typename FieldType::ValueType T;
231
232 Contour(const FieldType& f, T value, bool computeScalars):
233 Result(NULL),
234 Field(f),
235 Value(value),
236 ComputeScalars(computeScalars),
237 Name()
238 {
239 }
240
242 {
243 this->Result=grid;
244 }
245
246 void setFieldName(const char* name)
247 {
248 this->Name=std::string(name);
249 }
250
251 template<typename LHS, typename RHS>
252 int operator()(LHS &dataSet, const RHS&) const
253 {
254 typedef CellTypeToType<RHS> VTKCellTypeStruct;
255 typedef DataSetTypeToType<CellTypeToType<RHS>,LHS> DataSetTypeToTypeStruct;
256
257 //get the mapped output type of this operation(MarchingCubes)
258 //todo make this a typedef on the MarchingCubes
259 typedef typename MarchingCubesOuputType< typename VTKCellTypeStruct::DaxCellType >::type OutCellType;
260
261 //get the input dataset type
262 typedef typename DataSetTypeToTypeStruct::DaxDataSetType InputDataSetType;
263
264 //construct the output grid type to use the vtk containers
265 //as we know we are going back to vtk. In a more general framework
266 //we would want a tag to say what the destination container tag types
267 //are. We don't need the points container be
269 dax::cont::UnstructuredGrid<OutCellType,
272
273 InputDataSetType inputDaxData = vtkToDax::dataSetConverter(&dataSet,
274 DataSetTypeToTypeStruct());
275
277 int result = mc(inputDaxData,
278 &dataSet,
279 resultGrid,
280 this->Result,
281 this->Value,
282 this->Field,
283 this->ComputeScalars);
284
285 return result;
286 }
287 private:
288 vtkPolyData* Result;
289 FieldType Field;
290 T Value;
291 bool ComputeScalars;
292 std::string Name;
293
294 };
295}
296
297#endif //vtkToDax_Contour_h
virtual char * GetName()
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:55
int SetActiveAttribute(const char *name, int attributeType)
Make the array with the given name the active attribute.
abstract class to specify dataset behavior
Definition: vtkDataSet.h:63
vtkPointData * GetPointData()
Return a pointer to this dataset's point data.
Definition: vtkDataSet.h:256
Dispatch to functor based on a pointer type.
Definition: vtkDispatcher.h:92
ReturnType Go(BaseLhs *lhs)
Given a pointer to an object that derives from the BaseLhs we find the matching functor that was adde...
void Add(Functor fun)
Add in a functor that is mapped to the template SomeLhs parameter.
dynamic, self-adjusting array of double
dynamic, self-adjusting array of float
Definition: vtkFloatArray.h:42
represent and manipulate point attribute data
Definition: vtkPointData.h:38
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:86
void dataSetConverter(const dax::cont::UniformGrid<> &grid, vtkImageData *output)
VTKDataSetType::DaxDataSetType dataSetConverter(vtkImageData *input, VTKDataSetType)
@ value
Definition: vtkX3D.h:220
@ type
Definition: vtkX3D.h:516
@ name
Definition: vtkX3D.h:219
@ string
Definition: vtkX3D.h:490
FieldType_ FieldType
Definition: Contour.h:228
void setOutputGrid(vtkPolyData *grid)
Definition: Contour.h:241
FieldType::ValueType T
Definition: Contour.h:230
Contour(const FieldType &f, T value, bool computeScalars)
Definition: Contour.h:232
void setFieldName(const char *name)
Definition: Contour.h:246
int operator()(LHS &dataSet, const RHS &) const
Definition: Contour.h:252
int operator()(const dax::cont::UnstructuredGrid< dax::CellTagTetrahedron, GridCellContainer, GridPointContainer, Adapter > &inDaxGrid, vtkDataSet *inVTKGrid, OutGridType &outDaxGeom, vtkPolyData *outVTKGrid, ValueType isoValue, const dax::cont::ArrayHandle< ValueType, Container1, Adapter > &mcHandle, bool computeScalars)
Definition: Contour.h:107
int DispatchWork(const InGridType &inDaxGrid, vtkDataSet *inVTKGrid, OutGridType &outDaxGeom, vtkPolyData *outVTKGrid, CountWorkletType &countWorklet, GenerateWorkletType &generateWorklet, const dax::cont::ArrayHandle< ValueType, Container1, Adapter > &mcHandle, bool computeScalars)
Definition: Contour.h:140
int operator()(const InGridType &inDaxGrid, vtkDataSet *inVTKGrid, OutGridType &outDaxGeom, vtkPolyData *outVTKGrid, ValueType isoValue, const dax::cont::ArrayHandle< ValueType, Container1, Adapter > &mcHandle, bool computeScalars)
Definition: Contour.h:77
int operator()(const InGridType &, vtkDataSet *, OutGridType &, vtkPolyData *, ValueType, const dax::cont::ArrayHandle< ValueType, Container1, Adapter > &, bool)
Definition: Contour.h:56
VTKCellType
Definition: vtkCellType.h:43