VTK-m  2.0
DataSetBuilderExplicit.h
Go to the documentation of this file.
1 //============================================================================
2 // Copyright (c) Kitware, Inc.
3 // All rights reserved.
4 // See LICENSE.txt for details.
5 //
6 // This software is distributed WITHOUT ANY WARRANTY; without even
7 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
8 // PURPOSE. See the above copyright notice for more information.
9 //============================================================================
10 #ifndef vtk_m_cont_DataSetBuilderExplicit_h
11 #define vtk_m_cont_DataSetBuilderExplicit_h
12 
18 #include <vtkm/cont/DataSet.h>
19 
20 namespace vtkm
21 {
22 namespace cont
23 {
24 
25 //Coordinates builder??
26 
27 class VTKM_CONT_EXPORT DataSetBuilderExplicit
28 {
29 public:
30  VTKM_CONT
32 
33  //Single cell explicits.
34  //TODO
35 
36  //Zoo explicit cell
37  template <typename T>
38  VTKM_CONT static vtkm::cont::DataSet Create(const std::vector<T>& xVals,
39  const std::vector<vtkm::UInt8>& shapes,
40  const std::vector<vtkm::IdComponent>& numIndices,
41  const std::vector<vtkm::Id>& connectivity,
42  const std::string& coordsNm = "coords")
43  {
44  std::vector<T> yVals(xVals.size(), 0), zVals(xVals.size(), 0);
46  xVals, yVals, zVals, shapes, numIndices, connectivity, coordsNm);
47  }
48 
49  template <typename T>
50  VTKM_CONT static vtkm::cont::DataSet Create(const std::vector<T>& xVals,
51  const std::vector<T>& yVals,
52  const std::vector<vtkm::UInt8>& shapes,
53  const std::vector<vtkm::IdComponent>& numIndices,
54  const std::vector<vtkm::Id>& connectivity,
55  const std::string& coordsNm = "coords")
56  {
57  std::vector<T> zVals(xVals.size(), 0);
59  xVals, yVals, zVals, shapes, numIndices, connectivity, coordsNm);
60  }
61 
62  template <typename T>
63  VTKM_CONT static vtkm::cont::DataSet Create(const std::vector<T>& xVals,
64  const std::vector<T>& yVals,
65  const std::vector<T>& zVals,
66  const std::vector<vtkm::UInt8>& shapes,
67  const std::vector<vtkm::IdComponent>& numIndices,
68  const std::vector<vtkm::Id>& connectivity,
69  const std::string& coordsNm = "coords");
70 
71  template <typename T>
72  VTKM_CONT static vtkm::cont::DataSet Create(const std::vector<vtkm::Vec<T, 3>>& coords,
73  const std::vector<vtkm::UInt8>& shapes,
74  const std::vector<vtkm::IdComponent>& numIndices,
75  const std::vector<vtkm::Id>& connectivity,
76  const std::string& coordsNm = "coords");
77 
78  template <typename T>
83  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
84  const std::string& coordsNm = "coords")
85  {
86  auto offsets = vtkm::cont::ConvertNumComponentsToOffsets(numIndices);
87  return DataSetBuilderExplicit::BuildDataSet(coords, shapes, offsets, connectivity, coordsNm);
88  }
89 
90  template <typename T, typename CellShapeTag>
91  VTKM_CONT static vtkm::cont::DataSet Create(const std::vector<vtkm::Vec<T, 3>>& coords,
92  CellShapeTag tag,
93  vtkm::IdComponent numberOfPointsPerCell,
94  const std::vector<vtkm::Id>& connectivity,
95  const std::string& coordsNm = "coords");
96 
97  template <typename T, typename CellShapeTag>
100  CellShapeTag tag,
101  vtkm::IdComponent numberOfPointsPerCell,
102  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
103  const std::string& coordsNm = "coords")
104  {
106  coords, tag, numberOfPointsPerCell, connectivity, coordsNm);
107  }
108 
109 private:
110  template <typename T>
111  VTKM_CONT static vtkm::cont::DataSet BuildDataSet(
114  const vtkm::cont::ArrayHandle<vtkm::Id>& offsets,
115  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
116  const std::string& coordsNm);
117 
118  template <typename T, typename CellShapeTag>
119  VTKM_CONT static vtkm::cont::DataSet BuildDataSet(
121  CellShapeTag tag,
122  vtkm::IdComponent numberOfPointsPerCell,
123  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
124  const std::string& coordsNm);
125 };
126 
127 template <typename T>
129  const std::vector<T>& xVals,
130  const std::vector<T>& yVals,
131  const std::vector<T>& zVals,
132  const std::vector<vtkm::UInt8>& shapes,
133  const std::vector<vtkm::IdComponent>& numIndices,
134  const std::vector<vtkm::Id>& connectivity,
135  const std::string& coordsNm)
136 {
137  VTKM_ASSERT(xVals.size() == yVals.size() && yVals.size() == zVals.size() && xVals.size() > 0);
138 
140  coordsArray.Allocate(static_cast<vtkm::Id>(xVals.size()));
141  auto coordsPortal = coordsArray.WritePortal();
142  for (std::size_t index = 0; index < xVals.size(); ++index)
143  {
144  coordsPortal.Set(static_cast<vtkm::Id>(index),
145  vtkm::make_Vec(static_cast<vtkm::FloatDefault>(xVals[index]),
146  static_cast<vtkm::FloatDefault>(yVals[index]),
147  static_cast<vtkm::FloatDefault>(zVals[index])));
148  }
149 
150  auto shapesArray = vtkm::cont::make_ArrayHandle(shapes, vtkm::CopyFlag::On);
151  auto connArray = vtkm::cont::make_ArrayHandle(connectivity, vtkm::CopyFlag::On);
152 
153  auto offsetsArray = vtkm::cont::ConvertNumComponentsToOffsets(
155 
157  coordsArray, shapesArray, offsetsArray, connArray, coordsNm);
158 }
159 
160 template <typename T>
162  const std::vector<vtkm::Vec<T, 3>>& coords,
163  const std::vector<vtkm::UInt8>& shapes,
164  const std::vector<vtkm::IdComponent>& numIndices,
165  const std::vector<vtkm::Id>& connectivity,
166  const std::string& coordsNm)
167 {
168  auto coordsArray = vtkm::cont::make_ArrayHandle(coords, vtkm::CopyFlag::On);
169  auto shapesArray = vtkm::cont::make_ArrayHandle(shapes, vtkm::CopyFlag::On);
170  auto connArray = vtkm::cont::make_ArrayHandle(connectivity, vtkm::CopyFlag::On);
171 
172  auto offsetsArray = vtkm::cont::ConvertNumComponentsToOffsets(
174 
176  coordsArray, shapesArray, offsetsArray, connArray, coordsNm);
177 }
178 
179 template <typename T>
183  const vtkm::cont::ArrayHandle<vtkm::Id>& offsets,
184  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
185  const std::string& coordsNm)
186 {
187  vtkm::cont::DataSet dataSet;
188 
189  dataSet.AddCoordinateSystem(vtkm::cont::CoordinateSystem(coordsNm, coords));
190  vtkm::Id nPts = static_cast<vtkm::Id>(coords.GetNumberOfValues());
192 
193  cellSet.Fill(nPts, shapes, connectivity, offsets);
194  dataSet.SetCellSet(cellSet);
195 
196  return dataSet;
197 }
198 
199 template <typename T, typename CellShapeTag>
201  const std::vector<vtkm::Vec<T, 3>>& coords,
202  CellShapeTag tag,
203  vtkm::IdComponent numberOfPointsPerCell,
204  const std::vector<vtkm::Id>& connectivity,
205  const std::string& coordsNm)
206 {
207  auto coordsArray = vtkm::cont::make_ArrayHandle(coords, vtkm::CopyFlag::On);
208  auto connArray = vtkm::cont::make_ArrayHandle(connectivity, vtkm::CopyFlag::On);
209 
211  coordsArray, tag, numberOfPointsPerCell, connArray, coordsNm);
212 }
213 
214 template <typename T, typename CellShapeTag>
217  CellShapeTag tag,
218  vtkm::IdComponent numberOfPointsPerCell,
219  const vtkm::cont::ArrayHandle<vtkm::Id>& connectivity,
220  const std::string& coordsNm)
221 {
222  (void)tag; //C4100 false positive workaround
223  vtkm::cont::DataSet dataSet;
224 
225  dataSet.AddCoordinateSystem(vtkm::cont::CoordinateSystem(coordsNm, coords));
227 
228  cellSet.Fill(coords.GetNumberOfValues(), tag.Id, numberOfPointsPerCell, connectivity);
229  dataSet.SetCellSet(cellSet);
230 
231  return dataSet;
232 }
233 
234 class VTKM_CONT_EXPORT DataSetBuilderExplicitIterative
235 {
236 public:
237  VTKM_CONT
239 
240  VTKM_CONT
241  void Begin(const std::string& coordName = "coords");
242 
243  //Define points.
244  VTKM_CONT
245  vtkm::cont::DataSet Create();
246 
247  VTKM_CONT
248  vtkm::Id AddPoint(const vtkm::Vec3f& pt);
249 
250  VTKM_CONT
251  vtkm::Id AddPoint(const vtkm::FloatDefault& x,
252  const vtkm::FloatDefault& y,
253  const vtkm::FloatDefault& z = 0);
254 
255  template <typename T>
256  VTKM_CONT vtkm::Id AddPoint(const T& x, const T& y, const T& z = 0)
257  {
258  return AddPoint(static_cast<vtkm::FloatDefault>(x),
259  static_cast<vtkm::FloatDefault>(y),
260  static_cast<vtkm::FloatDefault>(z));
261  }
262 
263  template <typename T>
265  {
266  return AddPoint(static_cast<vtkm::Vec3f>(pt));
267  }
268 
269  //Define cells.
270  VTKM_CONT
271  void AddCell(vtkm::UInt8 shape);
272 
273  VTKM_CONT
274  void AddCell(const vtkm::UInt8& shape, const std::vector<vtkm::Id>& conn);
275 
276  VTKM_CONT
277  void AddCell(const vtkm::UInt8& shape, const vtkm::Id* conn, const vtkm::IdComponent& n);
278 
279  VTKM_CONT
280  void AddCellPoint(vtkm::Id pointIndex);
281 
282 private:
283  std::string coordNm;
284 
285  std::vector<vtkm::Vec3f> points;
286  std::vector<vtkm::UInt8> shapes;
287  std::vector<vtkm::IdComponent> numIdx;
288  std::vector<vtkm::Id> connectivity;
289 };
290 }
291 }
292 
293 #endif //vtk_m_cont_DataSetBuilderExplicit_h
vtkm::cont::make_ArrayHandle
VTKM_CONT vtkm::cont::ArrayHandleBasic< T > make_ArrayHandle(const T *array, vtkm::Id numberOfValues, vtkm::CopyFlag copy)
A convenience function for creating an ArrayHandle from a standard C array.
Definition: ArrayHandleBasic.h:217
vtkm::cont::DataSetBuilderExplicitIterative::AddPoint
VTKM_CONT vtkm::Id AddPoint(const T &x, const T &y, const T &z=0)
Definition: DataSetBuilderExplicit.h:256
vtkm::cont::ArrayHandle
Manages an array-worth of data.
Definition: ArrayHandle.h:283
vtkm::cont::DataSetBuilderExplicit::BuildDataSet
static VTKM_CONT vtkm::cont::DataSet BuildDataSet(const vtkm::cont::ArrayHandle< vtkm::Vec< T, 3 >> &coords, const vtkm::cont::ArrayHandle< vtkm::UInt8 > &shapes, const vtkm::cont::ArrayHandle< vtkm::Id > &offsets, const vtkm::cont::ArrayHandle< vtkm::Id > &connectivity, const std::string &coordsNm)
Definition: DataSetBuilderExplicit.h:180
vtkm::cont::DataSetBuilderExplicit::Create
static VTKM_CONT vtkm::cont::DataSet Create(const std::vector< T > &xVals, const std::vector< T > &yVals, const std::vector< vtkm::UInt8 > &shapes, const std::vector< vtkm::IdComponent > &numIndices, const std::vector< vtkm::Id > &connectivity, const std::string &coordsNm="coords")
Definition: DataSetBuilderExplicit.h:50
vtkm
Groups connected points that have the same field value.
Definition: Atomic.h:19
ArrayHandleCast.h
VTKM_ASSERT
#define VTKM_ASSERT(condition)
Definition: Assert.h:43
vtkm::cont::ArrayHandle::Allocate
VTKM_CONT void Allocate(vtkm::Id numberOfValues, vtkm::CopyFlag preserve, vtkm::cont::Token &token) const
Allocates an array large enough to hold the given number of values.
Definition: ArrayHandle.h:465
vtkm::IdComponent
vtkm::Int32 IdComponent
Represents a component ID (index of component in a vector).
Definition: Types.h:168
ArrayPortalToIterators.h
ArrayHandleCompositeVector.h
vtkm::cont::DataSetBuilderExplicit::DataSetBuilderExplicit
VTKM_CONT DataSetBuilderExplicit()
Definition: DataSetBuilderExplicit.h:31
vtkm::cont::CellSetSingleType
Definition: CastAndCall.h:34
vtkm::cont::DataSet
Definition: DataSet.h:34
vtkm::cont::DataSetBuilderExplicitIterative::AddPoint
VTKM_CONT vtkm::Id AddPoint(const vtkm::Vec< T, 3 > &pt)
Definition: DataSetBuilderExplicit.h:264
CoordinateSystem.h
vtkm::Id
vtkm::Int32 Id
Represents an ID (index into arrays).
Definition: Types.h:191
vtkm::cont::CoordinateSystem
Definition: CoordinateSystem.h:25
vtkm::cont::DataSetBuilderExplicitIterative::points
std::vector< vtkm::Vec3f > points
Definition: DataSetBuilderExplicit.h:285
vtkm::cont::DataSetBuilderExplicitIterative::connectivity
std::vector< vtkm::Id > connectivity
Definition: DataSetBuilderExplicit.h:288
vtkm::cont::DataSetBuilderExplicit::Create
static VTKM_CONT vtkm::cont::DataSet Create(const vtkm::cont::ArrayHandle< vtkm::Vec< T, 3 >> &coords, CellShapeTag tag, vtkm::IdComponent numberOfPointsPerCell, const vtkm::cont::ArrayHandle< vtkm::Id > &connectivity, const std::string &coordsNm="coords")
Definition: DataSetBuilderExplicit.h:98
vtkm::cont::DataSetBuilderExplicit::Create
static VTKM_CONT vtkm::cont::DataSet Create(const vtkm::cont::ArrayHandle< vtkm::Vec< T, 3 >> &coords, const vtkm::cont::ArrayHandle< vtkm::UInt8 > &shapes, const vtkm::cont::ArrayHandle< vtkm::IdComponent > &numIndices, const vtkm::cont::ArrayHandle< vtkm::Id > &connectivity, const std::string &coordsNm="coords")
Definition: DataSetBuilderExplicit.h:79
VTKM_CONT
#define VTKM_CONT
Definition: ExportMacros.h:57
vtkm::cont::ArrayHandle::WritePortal
VTKM_CONT WritePortalType WritePortal() const
Get an array portal that can be used in the control environment.
Definition: ArrayHandle.h:435
vtkm::cont::CellSetExplicit::Fill
VTKM_CONT void Fill(vtkm::Id numPoints, const vtkm::cont::ArrayHandle< vtkm::UInt8, ShapesStorageTag > &cellTypes, const vtkm::cont::ArrayHandle< vtkm::Id, ConnectivityStorageTag > &connectivity, const vtkm::cont::ArrayHandle< vtkm::Id, OffsetsStorageTag > &offsets)
Second method to add cells – all at once.
vtkm::UInt8
uint8_t UInt8
Definition: Types.h:157
vtkm::make_Vec
constexpr VTKM_EXEC_CONT vtkm::Vec< T, vtkm::IdComponent(sizeof...(Ts)+1)> make_Vec(T value0, Ts &&... args)
Initializes and returns a Vec containing all the arguments.
Definition: Types.h:1212
vtkm::CopyFlag::On
@ On
vtkm::Vec< T, 3 >
Definition: Types.h:975
vtkm::Vec< vtkm::FloatDefault, 3 >
vtkm::cont::ConvertNumComponentsToOffsets
VTKM_CONT_EXPORT void ConvertNumComponentsToOffsets(const vtkm::cont::UnknownArrayHandle &numComponentsArray, vtkm::cont::ArrayHandle< vtkm::Id > &offsetsArray, vtkm::Id &componentsArraySize, vtkm::cont::DeviceAdapterId device=vtkm::cont::DeviceAdapterTagAny{})
vtkm::CopyFlag::Off
@ Off
vtkm::cont::DataSetBuilderExplicitIterative::numIdx
std::vector< vtkm::IdComponent > numIdx
Definition: DataSetBuilderExplicit.h:287
vtkm::FloatDefault
vtkm::Float32 FloatDefault
The floating point type to use when no other precision is specified.
Definition: Types.h:198
ConvertNumComponentsToOffsets.h
vtkm::cont::CellSetExplicit
Definition: CastAndCall.h:36
vtkm::cont::DataSetBuilderExplicitIterative::shapes
std::vector< vtkm::UInt8 > shapes
Definition: DataSetBuilderExplicit.h:286
vtkm::cont::DataSet::AddCoordinateSystem
VTKM_CONT vtkm::IdComponent AddCoordinateSystem(const vtkm::cont::CoordinateSystem &cs)
Adds the given CoordinateSystem to the DataSet.
vtkm::cont::DataSetBuilderExplicit
Definition: DataSetBuilderExplicit.h:27
vtkm::cont::DataSet::SetCellSet
VTKM_CONT void SetCellSet(const CellSetType &cellSet)
Definition: DataSet.h:355
vtkm::cont::CellSetSingleType::Fill
VTKM_CONT void Fill(vtkm::Id numPoints, vtkm::UInt8 shapeId, vtkm::IdComponent numberOfPointsPerCell, const vtkm::cont::ArrayHandle< vtkm::Id, ConnectivityStorageTag > &connectivity)
Definition: CellSetSingleType.h:186
vtkm::cont::DataSetBuilderExplicitIterative::coordNm
std::string coordNm
Definition: DataSetBuilderExplicit.h:283
vtkm::cont::DataSetBuilderExplicit::Create
static VTKM_CONT vtkm::cont::DataSet Create(const std::vector< T > &xVals, const std::vector< vtkm::UInt8 > &shapes, const std::vector< vtkm::IdComponent > &numIndices, const std::vector< vtkm::Id > &connectivity, const std::string &coordsNm="coords")
Definition: DataSetBuilderExplicit.h:38
DataSet.h
vtkm::cont::DataSetBuilderExplicitIterative
Definition: DataSetBuilderExplicit.h:234