VTK
Threshold.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_Threshold_h
18#define vtkToDax_Threshold_h
19
24
27
28#include <dax/cont/DispatcherGenerateTopology.h>
29#include <dax/cont/DispatcherMapCell.h>
30#include <dax/worklet/Threshold.h>
31
32namespace
33{
34template <typename T> struct ThresholdOuputType
35{
36 typedef T type;
37};
38template <> struct ThresholdOuputType< dax::CellTagVoxel >
39{
40 typedef dax::CellTagHexahedron type;
41};
42}
43
44namespace vtkToDax
45{
46
47 template<int B>
49 {
50 template<class InGridType,
51 class OutGridType,
52 typename ValueType,
53 class Container,
54 class Adapter>
55 int operator()(const InGridType &,
56 vtkDataSet *,
57 OutGridType &,
59 ValueType,
60 ValueType,
61 const dax::cont::ArrayHandle<ValueType,Container,Adapter> &)
62 {
63 vtkGenericWarningMacro(
64 << "Not calling Dax, GridType-CellType combination not supported");
65 return 0;
66 }
67 };
68 template<>
69 struct DoThreshold<1>
70 {
71 template<class InGridType,
72 class OutGridType,
73 typename ValueType,
74 class Container,
75 class Adapter>
77 const InGridType &inDaxGrid,
78 vtkDataSet *inVTKGrid,
79 OutGridType &outDaxGeom,
80 vtkUnstructuredGrid *outVTKGrid,
81 ValueType thresholdMin,
82 ValueType thresholdMax,
83 const dax::cont::ArrayHandle<ValueType,Container,Adapter> &thresholdHandle)
84 {
85 int result=1;
86 try
87 {
88 typedef dax::worklet::ThresholdCount<ValueType> ThresholdCountType;
89
90 typedef dax::cont::DispatcherGenerateTopology<
91 dax::worklet::ThresholdTopology,
92 dax::cont::ArrayHandle< dax::Id >,
93 Adapter > DispatchGT;
94
95 typedef typename DispatchGT::CountHandleType CountHandleType;
96
97 ThresholdCountType countWorklet(thresholdMin,thresholdMax);
98 dax::cont::DispatcherMapCell<ThresholdCountType,Adapter>
99 dispatchCount( countWorklet );
100
101 CountHandleType count;
102 dispatchCount.Invoke(inDaxGrid, thresholdHandle, count);
103
104 DispatchGT resolveTopology(count);
105 resolveTopology.Invoke(inDaxGrid,outDaxGeom);
106
107 // Convert output geometry to VTK.
108 daxToVtk::dataSetConverter(outDaxGeom, outVTKGrid);
109
110 // Copy arrays where possible.
111 vtkToDax::CompactPointField<DispatchGT> compact(resolveTopology,
112 outVTKGrid);
113
114 vtkDispatcher<vtkAbstractArray,int> compactDispatcher;
115 compactDispatcher.Add<vtkFloatArray>(compact);
116 compactDispatcher.Add<vtkDoubleArray>(compact);
117 compactDispatcher.Add<vtkUnsignedCharArray>(compact);
118 compactDispatcher.Add<vtkIntArray>(compact);
119
120 vtkPointData *pd = inVTKGrid->GetPointData();
121 for (int arrayIndex = 0;
122 arrayIndex < pd->GetNumberOfArrays();
123 arrayIndex++)
124 {
125 vtkDataArray *array = pd->GetArray(arrayIndex);
126 if (array == NULL) { continue; }
127
128 compactDispatcher.Go(array);
129 }
130
131 // Pass information about attributes.
132 for (int attributeType = 0;
134 attributeType++)
135 {
136 vtkDataArray *attribute = pd->GetAttribute(attributeType);
137 if (attribute == NULL) { continue; }
138 outVTKGrid->GetPointData()->SetActiveAttribute(attribute->GetName(),
139 attributeType);
140 }
141 }
142 catch(const dax::cont::ErrorControlOutOfMemory& error)
143 {
144 std::cerr << "Ran out of memory trying to use the GPU" << std::endl;
145 std::cerr << error.GetMessage() << std::endl;
146 result = 0;
147 }
148 catch(const dax::cont::ErrorExecution& error)
149 {
150 std::cerr << "Got ErrorExecution from Dax." << std::endl;
151 std::cerr << error.GetMessage() << std::endl;
152 result = 0;
153 }
154 return result;
155 }
156 };
157
158 template<typename FieldType_>
160 {
161 public:
162 typedef FieldType_ FieldType;
163 //we expect FieldType_ to be an dax::cont::ArrayHandle
164 typedef typename FieldType::ValueType T;
165
166 Threshold(const FieldType& f, T min, T max):
167 Result(NULL),
168 Field(f),
169 Min(min),
170 Max(max),
171 Name()
172 {
173 }
174
176 {
177 Result=grid;
178 }
179
180 void setFieldName(const char* name)
181 {
182 Name=std::string(name);
183 }
184
185 template<typename LHS, typename RHS>
186 int operator()(LHS &dataSet, const RHS&) const
187 {
188 typedef CellTypeToType<RHS> VTKCellTypeStruct;
189 typedef DataSetTypeToType<CellTypeToType<RHS>,LHS> DataSetTypeToTypeStruct;
190
191 //get the mapped output type of this operation(threshold)
192 //todo make this a typedef on the threshold
193 typedef typename ThresholdOuputType< typename VTKCellTypeStruct::DaxCellType >::type OutCellType;
194
195 //get the input dataset type
196 typedef typename DataSetTypeToTypeStruct::DaxDataSetType InputDataSetType;
197
198 //construct the output grid type to use the vtk containers
199 //as we know we are going back to vtk. In a more general framework
200 //we would want a tag to say what the destination container tag types
201 //are
203 dax::cont::UnstructuredGrid<OutCellType,
206
207 InputDataSetType inputDaxData = vtkToDax::dataSetConverter(&dataSet,
208 DataSetTypeToTypeStruct());
209
211 int result = threshold(inputDaxData,
212 &dataSet,
213 resultGrid,
214 this->Result,
215 this->Min,
216 this->Max,
217 this->Field);
218
219 return result;
220
221 }
222 private:
223 vtkUnstructuredGrid* Result;
224 FieldType Field;
225 T Min;
226 T Max;
227 std::string Name;
228
229 };
230}
231
232#endif //vtkToDax_Threshold_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
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:46
represent and manipulate point attribute data
Definition: vtkPointData.h:38
dynamic, self-adjusting array of unsigned char
dataset represents arbitrary combinations of all possible cell types
void dataSetConverter(const dax::cont::UniformGrid<> &grid, vtkImageData *output)
Definition: Containers.h:51
VTKDataSetType::DaxDataSetType dataSetConverter(vtkImageData *input, VTKDataSetType)
@ type
Definition: vtkX3D.h:516
@ name
Definition: vtkX3D.h:219
@ string
Definition: vtkX3D.h:490
int operator()(const InGridType &inDaxGrid, vtkDataSet *inVTKGrid, OutGridType &outDaxGeom, vtkUnstructuredGrid *outVTKGrid, ValueType thresholdMin, ValueType thresholdMax, const dax::cont::ArrayHandle< ValueType, Container, Adapter > &thresholdHandle)
Definition: Threshold.h:76
int operator()(const InGridType &, vtkDataSet *, OutGridType &, vtkUnstructuredGrid *, ValueType, ValueType, const dax::cont::ArrayHandle< ValueType, Container, Adapter > &)
Definition: Threshold.h:55
FieldType::ValueType T
Definition: Threshold.h:164
int operator()(LHS &dataSet, const RHS &) const
Definition: Threshold.h:186
void setOutputGrid(vtkUnstructuredGrid *grid)
Definition: Threshold.h:175
Threshold(const FieldType &f, T min, T max)
Definition: Threshold.h:166
void setFieldName(const char *name)
Definition: Threshold.h:180
FieldType_ FieldType
Definition: Threshold.h:162
VTKCellType
Definition: vtkCellType.h:43
#define max(a, b)