OGS
MeshLib::VtkMeshConverter Class Reference

Detailed Description

Converter for VtkUnstructured Grids to OGS meshes.

Definition at line 39 of file VtkMeshConverter.h.

#include <VtkMeshConverter.h>

Static Public Member Functions

static MeshLib::MeshconvertUnstructuredGrid (vtkUnstructuredGrid *grid, std::string const &mesh_name="vtkUnstructuredGrid")
 Converts a vtkUnstructuredGrid object to a Mesh. More...
 

Static Private Member Functions

static void convertScalarArrays (vtkUnstructuredGrid &grid, MeshLib::Mesh &mesh)
 
static std::vector< MeshLib::Node * > createNodeVector (std::vector< double > const &elevation, std::vector< int > &node_idx_map, GeoLib::RasterHeader const &header, bool use_elevation)
 
static std::vector< MeshLib::Element * > createElementVector (std::vector< double > const &pix_val, std::vector< bool > const &pix_vis, std::vector< MeshLib::Node * > const &nodes, std::vector< int > const &node_idx_map, std::size_t const imgHeight, std::size_t const imgWidth, MeshElemType elem_type)
 Creates a mesh element vector based on image data. More...
 
template<typename T >
static void fillPropertyVector (MeshLib::PropertyVector< T > &prop_vec, std::vector< double > const &pix_val, std::vector< bool > const &pix_vis, const std::size_t &imgHeight, const std::size_t &imgWidth, MeshElemType elem_type)
 Creates a scalar array/mesh property based on pixel values. More...
 
static void convertArray (vtkDataArray &array, MeshLib::Properties &properties, MeshLib::MeshItemType type)
 
template<typename T >
static void convertTypedArray (vtkDataArray &array, MeshLib::Properties &properties, MeshLib::MeshItemType type)
 

Member Function Documentation

◆ convertArray()

void MeshLib::VtkMeshConverter::convertArray ( vtkDataArray &  array,
MeshLib::Properties properties,
MeshLib::MeshItemType  type 
)
staticprivate

Definition at line 273 of file VtkMeshConverter.cpp.

276{
277 if (vtkDoubleArray::SafeDownCast(&array))
278 {
279 VtkMeshConverter::convertTypedArray<double>(array, properties, type);
280 return;
281 }
282
283 if (vtkFloatArray::SafeDownCast(&array))
284 {
285 VtkMeshConverter::convertTypedArray<float>(array, properties, type);
286 return;
287 }
288
289 if (vtkBitArray::SafeDownCast(&array))
290 {
291 VtkMeshConverter::convertTypedArray<bool>(array, properties, type);
292 return;
293 }
294
295 // This also covers the vtkTypeInt8Array type, which is derived from the
296 // vtkCharArray type.
297 if (vtkCharArray::SafeDownCast(&array) ||
298 vtkSignedCharArray::SafeDownCast(&array))
299 {
300 if constexpr (sizeof(std::int8_t) != sizeof(char))
301 {
302 OGS_FATAL(
303 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
304 "not convertible to char.",
305 array.GetName(), array.GetDataTypeAsString());
306 }
307 VtkMeshConverter::convertTypedArray<char>(array, properties, type);
308 return;
309 }
310
311 // This also covers the vtkTypeInt16Array type, which is derived from the
312 // vtkShortArray type.
313 if (vtkShortArray::SafeDownCast(&array))
314 {
315 if constexpr (sizeof(std::int16_t) != sizeof(short))
316 {
317 OGS_FATAL(
318 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
319 "not convertible to short.",
320 array.GetName(), array.GetDataTypeAsString());
321 }
322 VtkMeshConverter::convertTypedArray<short>(array, properties, type);
323 return;
324 }
325
326 // This also covers the vtkTypeInt32Array type, which is derived from the
327 // vtkIntArray type.
328 if (vtkIntArray::SafeDownCast(&array))
329 {
330 if constexpr (sizeof(std::int32_t) != sizeof(int))
331 {
332 OGS_FATAL(
333 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
334 "not convertible to int.",
335 array.GetName(), array.GetDataTypeAsString());
336 }
337 VtkMeshConverter::convertTypedArray<int>(array, properties, type);
338 return;
339 }
340
341 // This is not a unique conversion depending on the sizes of int, long, and
342 // long long. Converting to the apparently smaller type.
343 if (vtkLongArray::SafeDownCast(&array))
344 {
345 if constexpr (sizeof(long) == sizeof(int))
346 {
347 VtkMeshConverter::convertTypedArray<int>(array, properties, type);
348 return;
349 }
350 if constexpr (sizeof(long long) == sizeof(long))
351 {
352 VtkMeshConverter::convertTypedArray<long>(array, properties, type);
353 return;
354 }
355
356 OGS_FATAL(
357 "Array '{:s}' in VTU file uses unsupported data type '{:s}' ({:d}) "
358 "not convertible to int ({:d}), or long ({:d}) types.",
359 array.GetName(), array.GetDataTypeAsString(),
360 array.GetDataTypeSize(), sizeof(int), sizeof(long));
361 }
362
363 // Ensure, vtkTypeInt64Array and vtkLongLongArray are using the same type.
364 static_assert(sizeof(std::int64_t) == sizeof(long long));
365
366 // This also covers the vtkTypeInt64Array type, which is derived from the
367 // vtkLongLongArray type.
368 // Converting to the apparently smaller type.
369 if (vtkLongLongArray::SafeDownCast(&array))
370 {
371 if constexpr (sizeof(long long) == sizeof(long))
372 {
373 VtkMeshConverter::convertTypedArray<long>(array, properties, type);
374 return;
375 }
376 else // ll > l. Other cases are not possible because ll >= l in c++.
377 {
378 VtkMeshConverter::convertTypedArray<long long>(array, properties,
379 type);
380 return;
381 }
382 }
383
384 // This also covers the vtkTypeUInt8Array type, which is derived from the
385 // vtkUnsignedCharArray type.
386 if (vtkUnsignedCharArray::SafeDownCast(&array))
387 {
388 if constexpr (sizeof(std::uint8_t) != sizeof(unsigned char))
389 {
390 OGS_FATAL(
391 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
392 "not convertible to unsigned char.",
393 array.GetName(), array.GetDataTypeAsString());
394 }
395 VtkMeshConverter::convertTypedArray<unsigned char>(array, properties,
396 type);
397 return;
398 }
399
400 // This also covers the vtkTypeUInt16Array type, which is derived from the
401 // vtkUnsignedShortArray type.
402 if (vtkUnsignedShortArray::SafeDownCast(&array))
403 {
404 if constexpr (sizeof(std::uint16_t) != sizeof(unsigned short))
405 {
406 OGS_FATAL(
407 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
408 "not convertible to unsigned short.",
409 array.GetName(), array.GetDataTypeAsString());
410 }
411 VtkMeshConverter::convertTypedArray<unsigned short>(array, properties,
412 type);
413 return;
414 }
415
416 // This also covers the vtkTypeUInt32Array type, which is derived from the
417 // vtkUnsignedIntArray type.
418 if (vtkUnsignedIntArray::SafeDownCast(&array))
419 {
420 if constexpr (sizeof(std::uint32_t) != sizeof(unsigned))
421 {
422 OGS_FATAL(
423 "Array '{:s}' in VTU file uses unsupported data type '{:s}' "
424 "not convertible to unsigned.",
425 array.GetName(), array.GetDataTypeAsString());
426 }
427
428 // MaterialIDs are assumed to be integers
429 if (std::strncmp(array.GetName(), "MaterialIDs", 11) == 0)
430 {
431 VtkMeshConverter::convertTypedArray<int>(array, properties, type);
432 }
433 else
434 {
435 VtkMeshConverter::convertTypedArray<unsigned>(array, properties,
436 type);
437 }
438
439 return;
440 }
441
442 // This is not a unique conversion depending on the sizes of unsigned,
443 // unsigned long, and unsigned long long. Converting to the apparently
444 // smaller type.
445 if (vtkUnsignedLongArray::SafeDownCast(&array))
446 {
447 if constexpr (sizeof(unsigned long) == sizeof(unsigned))
448 {
449 VtkMeshConverter::convertTypedArray<unsigned>(array, properties,
450 type);
451 return;
452 }
453 if constexpr (sizeof(unsigned long long) == sizeof(unsigned long))
454 {
455 VtkMeshConverter::convertTypedArray<unsigned long>(
456 array, properties, type);
457 return;
458 }
459
460 OGS_FATAL(
461 "Array '{:s}' in VTU file uses unsupported data type '{:s}' ({:d}) "
462 "not convertible to unsigned ({:d}), or unsigned long ({:d}) "
463 "types.",
464 array.GetName(), array.GetDataTypeAsString(),
465 array.GetDataTypeSize(), sizeof(unsigned), sizeof(unsigned long));
466 }
467
468 // Ensure, vtkTypeUInt64Array and vtkUnsignedLongLongArray are using the
469 // same type.
470 static_assert(sizeof(std::uint64_t) == sizeof(unsigned long long));
471
472 // This also covers the vtkTypeUInt64Array type, which is derived from the
473 // vtkUnsignedLongLongArray type.
474 // Converting to the apparently smaller type.
475 if (vtkUnsignedLongLongArray::SafeDownCast(&array))
476 {
477 if constexpr (sizeof(unsigned long long) == sizeof(unsigned long))
478 {
479 VtkMeshConverter::convertTypedArray<unsigned long>(
480 array, properties, type);
481 return;
482 }
483 else // ull > ul. Other cases are not possible because ull >= ul in
484 // c++.
485 {
486 VtkMeshConverter::convertTypedArray<unsigned long long>(
487 array, properties, type);
488 return;
489 }
490 }
491
492 WARN(
493 "Array '{:s}' in VTU file uses unsupported data type '{:s}' of size "
494 "{:d}. The data array will not be available.",
495 array.GetName(), array.GetDataTypeAsString(), array.GetDataTypeSize());
496}
#define OGS_FATAL(...)
Definition: Error.h:26
void WARN(char const *fmt, Args const &... args)
Definition: Logging.h:39

References OGS_FATAL, and WARN().

Referenced by convertScalarArrays().

◆ convertScalarArrays()

void MeshLib::VtkMeshConverter::convertScalarArrays ( vtkUnstructuredGrid &  grid,
MeshLib::Mesh mesh 
)
staticprivate

Definition at line 238 of file VtkMeshConverter.cpp.

240{
241 vtkPointData* point_data = grid.GetPointData();
242 auto const n_point_arrays =
243 static_cast<unsigned>(point_data->GetNumberOfArrays());
244 for (unsigned i = 0; i < n_point_arrays; ++i)
245 {
246 convertArray(*point_data->GetArray(i),
247 mesh.getProperties(),
249 }
250
251 vtkCellData* cell_data = grid.GetCellData();
252 auto const n_cell_arrays =
253 static_cast<unsigned>(cell_data->GetNumberOfArrays());
254 for (unsigned i = 0; i < n_cell_arrays; ++i)
255 {
256 convertArray(*cell_data->GetArray(i),
257 mesh.getProperties(),
259 }
260
261 vtkFieldData* field_data = grid.GetFieldData();
262 auto const n_field_arrays =
263 static_cast<unsigned>(field_data->GetNumberOfArrays());
264 for (unsigned i = 0; i < n_field_arrays; ++i)
265 {
267 *vtkDataArray::SafeDownCast(field_data->GetAbstractArray(i)),
268 mesh.getProperties(),
270 }
271}
Properties & getProperties()
Definition: Mesh.h:128
static void convertArray(vtkDataArray &array, MeshLib::Properties &properties, MeshLib::MeshItemType type)

References MeshLib::Cell, convertArray(), MeshLib::Mesh::getProperties(), MeshLib::IntegrationPoint, and MeshLib::Node.

Referenced by convertUnstructuredGrid().

◆ convertTypedArray()

template<typename T >
static void MeshLib::VtkMeshConverter::convertTypedArray ( vtkDataArray &  array,
MeshLib::Properties properties,
MeshLib::MeshItemType  type 
)
inlinestaticprivate

Definition at line 103 of file VtkMeshConverter.h.

106 {
107 // Keep for debugging purposes, since the vtk array type and the end
108 // type T are not always the same.
109 // DBUG(
110 // "Converting a vtk array '{:s}' with data type '{:s}' of "
111 // "size {:d} to a type '{:s}' of size {:d}.",
112 // array.GetName(), array.GetDataTypeAsString(),
113 // array.GetDataTypeSize(), typeid(T).name(), sizeof(T));
114
115 if (sizeof(T) != array.GetDataTypeSize())
116 {
117 OGS_FATAL(
118 "Trying to convert a vtk array '{:s}' with data type '{:s}' of "
119 "size {:d} to a different sized type '{:s}' of size {:d}.",
120 array.GetName(), array.GetDataTypeAsString(),
121 array.GetDataTypeSize(), typeid(T).name(), sizeof(T));
122 }
123
124 vtkIdType const nTuples(array.GetNumberOfTuples());
125 int const nComponents(array.GetNumberOfComponents());
126 char const* const array_name(array.GetName());
127
128 auto* const vec = properties.createNewPropertyVector<T>(
129 array_name, type, nComponents);
130 if (!vec)
131 {
132 WARN("Array {:s} could not be converted to PropertyVector.",
133 array_name);
134 return;
135 }
136 vec->reserve(nTuples * nComponents);
137 auto* data_array = static_cast<T*>(array.GetVoidPointer(0));
138 std::copy(&data_array[0],
139 &data_array[nTuples * nComponents],
140 std::back_inserter(*vec));
141 return;
142 }
PropertyVector< T > * createNewPropertyVector(std::string const &name, MeshItemType mesh_item_type, std::size_t n_components=1)
void copy(PETScVector const &x, PETScVector &y)
Definition: LinAlg.cpp:37

References MathLib::LinAlg::copy(), MeshLib::Properties::createNewPropertyVector(), OGS_FATAL, and WARN().

◆ convertUnstructuredGrid()

MeshLib::Mesh * MeshLib::VtkMeshConverter::convertUnstructuredGrid ( vtkUnstructuredGrid *  grid,
std::string const &  mesh_name = "vtkUnstructuredGrid" 
)
static

Converts a vtkUnstructuredGrid object to a Mesh.

Definition at line 65 of file VtkMeshConverter.cpp.

67{
68 if (!grid)
69 {
70 return nullptr;
71 }
72
73 // set mesh nodes
74 const std::size_t nNodes = grid->GetPoints()->GetNumberOfPoints();
75 std::vector<MeshLib::Node*> nodes(nNodes);
76 double* coords = nullptr;
77 for (std::size_t i = 0; i < nNodes; i++)
78 {
79 coords = grid->GetPoints()->GetPoint(i);
80 nodes[i] = new MeshLib::Node(coords[0], coords[1], coords[2], i);
81 }
82
83 // set mesh elements
84 const std::size_t nElems = grid->GetNumberOfCells();
85 std::vector<MeshLib::Element*> elements(nElems);
86 auto node_ids = vtkSmartPointer<vtkIdList>::New();
87 for (std::size_t i = 0; i < nElems; i++)
88 {
89 MeshLib::Element* elem;
90 grid->GetCellPoints(i, node_ids);
91
92 int cell_type = grid->GetCellType(i);
93 switch (cell_type)
94 {
95 case VTK_VERTEX:
96 {
97 elem = detail::createElementWithSameNodeOrder<MeshLib::Point>(
98 nodes, node_ids, i);
99 break;
100 }
101 case VTK_LINE:
102 {
103 elem = detail::createElementWithSameNodeOrder<MeshLib::Line>(
104 nodes, node_ids, i);
105 break;
106 }
107 case VTK_TRIANGLE:
108 {
109 elem = detail::createElementWithSameNodeOrder<MeshLib::Tri>(
110 nodes, node_ids, i);
111 break;
112 }
113 case VTK_QUAD:
114 {
115 elem = detail::createElementWithSameNodeOrder<MeshLib::Quad>(
116 nodes, node_ids, i);
117 break;
118 }
119 case VTK_PIXEL:
120 {
121 auto** quad_nodes = new MeshLib::Node*[4];
122 quad_nodes[0] = nodes[node_ids->GetId(0)];
123 quad_nodes[1] = nodes[node_ids->GetId(1)];
124 quad_nodes[2] = nodes[node_ids->GetId(3)];
125 quad_nodes[3] = nodes[node_ids->GetId(2)];
126 elem = new MeshLib::Quad(quad_nodes, i);
127 break;
128 }
129 case VTK_TETRA:
130 {
131 elem = detail::createElementWithSameNodeOrder<MeshLib::Tet>(
132 nodes, node_ids, i);
133 break;
134 }
135 case VTK_HEXAHEDRON:
136 {
137 elem = detail::createElementWithSameNodeOrder<MeshLib::Hex>(
138 nodes, node_ids, i);
139 break;
140 }
141 case VTK_VOXEL:
142 {
143 auto** voxel_nodes = new MeshLib::Node*[8];
144 voxel_nodes[0] = nodes[node_ids->GetId(0)];
145 voxel_nodes[1] = nodes[node_ids->GetId(1)];
146 voxel_nodes[2] = nodes[node_ids->GetId(3)];
147 voxel_nodes[3] = nodes[node_ids->GetId(2)];
148 voxel_nodes[4] = nodes[node_ids->GetId(4)];
149 voxel_nodes[5] = nodes[node_ids->GetId(5)];
150 voxel_nodes[6] = nodes[node_ids->GetId(7)];
151 voxel_nodes[7] = nodes[node_ids->GetId(6)];
152 elem = new MeshLib::Hex(voxel_nodes, i);
153 break;
154 }
155 case VTK_PYRAMID:
156 {
157 elem = detail::createElementWithSameNodeOrder<MeshLib::Pyramid>(
158 nodes, node_ids, i);
159 break;
160 }
161 case VTK_WEDGE:
162 {
163 auto** prism_nodes = new MeshLib::Node*[6];
164 for (unsigned j = 0; j < 3; ++j)
165 {
166 prism_nodes[j] = nodes[node_ids->GetId(j + 3)];
167 prism_nodes[j + 3] = nodes[node_ids->GetId(j)];
168 }
169 elem = new MeshLib::Prism(prism_nodes, i);
170 break;
171 }
172 case VTK_QUADRATIC_EDGE:
173 {
174 elem = detail::createElementWithSameNodeOrder<MeshLib::Line3>(
175 nodes, node_ids, i);
176 break;
177 }
178 case VTK_QUADRATIC_TRIANGLE:
179 {
180 elem = detail::createElementWithSameNodeOrder<MeshLib::Tri6>(
181 nodes, node_ids, i);
182 break;
183 }
184 case VTK_QUADRATIC_QUAD:
185 {
186 elem = detail::createElementWithSameNodeOrder<MeshLib::Quad8>(
187 nodes, node_ids, i);
188 break;
189 }
190 case VTK_BIQUADRATIC_QUAD:
191 {
192 elem = detail::createElementWithSameNodeOrder<MeshLib::Quad9>(
193 nodes, node_ids, i);
194 break;
195 }
196 case VTK_QUADRATIC_TETRA:
197 {
198 elem = detail::createElementWithSameNodeOrder<MeshLib::Tet10>(
199 nodes, node_ids, i);
200 break;
201 }
202 case VTK_QUADRATIC_HEXAHEDRON:
203 {
204 elem = detail::createElementWithSameNodeOrder<MeshLib::Hex20>(
205 nodes, node_ids, i);
206 break;
207 }
208 case VTK_QUADRATIC_PYRAMID:
209 {
210 elem =
211 detail::createElementWithSameNodeOrder<MeshLib::Pyramid13>(
212 nodes, node_ids, i);
213 break;
214 }
215 case VTK_QUADRATIC_WEDGE:
216 {
217 elem = detail::createElementWithSameNodeOrder<MeshLib::Prism15>(
218 nodes, node_ids, i);
219
220 break;
221 }
222 default:
223 ERR("VtkMeshConverter::convertUnstructuredGrid(): Unknown mesh "
224 "element type '{:d}'.",
225 cell_type);
226 return nullptr;
227 }
228
229 elements[i] = elem;
230 }
231
232 MeshLib::Mesh* mesh = new MeshLib::Mesh(mesh_name, nodes, elements);
233 convertScalarArrays(*grid, *mesh);
234
235 return mesh;
236}
void ERR(char const *fmt, Args const &... args)
Definition: Logging.h:44
static void convertScalarArrays(vtkUnstructuredGrid &grid, MeshLib::Mesh &mesh)
TemplateElement< MeshLib::PrismRule6 > Prism
Definition: Prism.h:25
TemplateElement< MeshLib::HexRule8 > Hex
Definition: Hex.h:25
TemplateElement< MeshLib::QuadRule4 > Quad
Definition: Quad.h:28

References convertScalarArrays(), and ERR().

Referenced by VtkVisPipelineView::convertVTKToOGSMesh(), MeshLib::IO::VtuInterface::readVTKFile(), and MeshLib::IO::VtuInterface::readVTUFile().

◆ createElementVector()

static std::vector< MeshLib::Element * > MeshLib::VtkMeshConverter::createElementVector ( std::vector< double > const &  pix_val,
std::vector< bool > const &  pix_vis,
std::vector< MeshLib::Node * > const &  nodes,
std::vector< int > const &  node_idx_map,
std::size_t const  imgHeight,
std::size_t const  imgWidth,
MeshElemType  elem_type 
)
staticprivate

Creates a mesh element vector based on image data.

◆ createNodeVector()

static std::vector< MeshLib::Node * > MeshLib::VtkMeshConverter::createNodeVector ( std::vector< double > const &  elevation,
std::vector< int > &  node_idx_map,
GeoLib::RasterHeader const &  header,
bool  use_elevation 
)
staticprivate

◆ fillPropertyVector()

template<typename T >
static void MeshLib::VtkMeshConverter::fillPropertyVector ( MeshLib::PropertyVector< T > &  prop_vec,
std::vector< double > const &  pix_val,
std::vector< bool > const &  pix_vis,
const std::size_t &  imgHeight,
const std::size_t &  imgWidth,
MeshElemType  elem_type 
)
inlinestaticprivate

Creates a scalar array/mesh property based on pixel values.

Definition at line 69 of file VtkMeshConverter.h.

75 {
76 for (std::size_t i = 0; i < imgHeight; i++)
77 {
78 for (std::size_t j = 0; j < imgWidth; j++)
79 {
80 if (!pix_vis[i * imgWidth + j])
81 {
82 continue;
83 }
84 T val(static_cast<T>(pix_val[i * (imgWidth + 1) + j]));
85 if (elem_type == MeshElemType::TRIANGLE)
86 {
87 prop_vec.push_back(val);
88 prop_vec.push_back(val);
89 }
90 else if (elem_type == MeshElemType::QUAD)
91 {
92 prop_vec.push_back(val);
93 }
94 }
95 }
96 }

References MeshLib::QUAD, and MeshLib::TRIANGLE.


The documentation for this class was generated from the following files: