ViennaLS
Loading...
Searching...
No Matches
lsVTKWriter.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <fstream>
4#include <string>
5#include <unordered_map>
6#include <utility>
7
8#include <lsFileFormats.hpp>
9#include <lsMesh.hpp>
10
11#include <vcLogger.hpp>
12#include <vcSmartPointer.hpp>
13
14#ifdef VIENNALS_USE_VTK
15#include <vtkCellArray.h>
16#include <vtkCellData.h>
17#include <vtkFloatArray.h>
18#include <vtkPointData.h>
19#include <vtkPoints.h>
20#include <vtkPolyData.h>
21#include <vtkSmartPointer.h>
22#include <vtkXMLPolyDataWriter.h>
23
24#include <vtkUnstructuredGrid.h>
25#include <vtkXMLUnstructuredGridWriter.h>
26#endif // VIENNALS_USE_VTK
27
28namespace viennals {
29
30using namespace viennacore;
31
33template <class T> class VTKWriter {
34 using MetaDataType = std::unordered_map<std::string, std::vector<double>>;
35
36 SmartPointer<Mesh<T>> mesh = nullptr;
38 std::string fileName;
39 MetaDataType metaData;
40
41#ifdef VIENNALS_USE_VTK
42 template <class In, class Out>
43 void addDataFromMesh(const In &inData, Out outData) const {
44 // now add pointData
45 for (unsigned i = 0; i < inData.getScalarDataSize(); ++i) {
46 vtkSmartPointer<vtkFloatArray> pointData =
47 vtkSmartPointer<vtkFloatArray>::New();
48 pointData->SetNumberOfComponents(1);
49 pointData->SetName(inData.getScalarDataLabel(i).c_str());
50 auto scalars = *(inData.getScalarData(i));
51 for (unsigned j = 0; j < inData.getScalarData(i)->size(); ++j) {
52 pointData->InsertNextValue(scalars[j]);
53 }
54 outData->AddArray(pointData);
55 }
56
57 // now add vector data
58 for (unsigned i = 0; i < inData.getVectorDataSize(); ++i) {
59 vtkSmartPointer<vtkFloatArray> vectorData =
60 vtkSmartPointer<vtkFloatArray>::New();
61 vectorData->SetNumberOfComponents(3);
62 vectorData->SetName(inData.getVectorDataLabel(i).c_str());
63 auto vectors = *(inData.getVectorData(i));
64 for (unsigned j = 0; j < inData.getVectorData(i)->size(); ++j) {
65 vectorData->InsertNextTuple3(vectors[j][0], vectors[j][1],
66 vectors[j][2]);
67 }
68 outData->AddArray(vectorData);
69 }
70 }
71
72 void addMetaDataToVTK(vtkDataSet *data) const {
73 if (metaData.empty()) {
74 return;
75 }
76
77 // add metadata to field data
78 vtkSmartPointer<vtkFieldData> fieldData = data->GetFieldData();
79 for (const auto &meta : metaData) {
80 if (meta.second.empty())
81 continue; // skip empty metadata
82
83 vtkSmartPointer<vtkFloatArray> metaDataArray =
84 vtkSmartPointer<vtkFloatArray>::New();
85 metaDataArray->SetName(meta.first.c_str());
86 metaDataArray->SetNumberOfValues(meta.second.size());
87 for (size_t i = 0; i < meta.second.size(); ++i) {
88 metaDataArray->SetValue(i, meta.second[i]);
89 }
90 fieldData->AddArray(metaDataArray);
91 }
92 }
93#endif // VIENNALS_USE_VTK
94
95public:
96 VTKWriter() = default;
97
98 VTKWriter(SmartPointer<Mesh<T>> passedMesh) : mesh(passedMesh) {}
99
100 VTKWriter(SmartPointer<Mesh<T>> passedMesh, std::string passedFileName)
101 : mesh(passedMesh), fileName(std::move(passedFileName)) {}
102
103 VTKWriter(SmartPointer<Mesh<T>> passedMesh, FileFormatEnum passedFormat,
104 std::string passedFileName)
105 : mesh(passedMesh), fileFormat(passedFormat),
106 fileName(std::move(passedFileName)) {}
107
108 void setMesh(SmartPointer<Mesh<T>> passedMesh) { mesh = passedMesh; }
109
111 void setFileFormat(FileFormatEnum passedFormat) { fileFormat = passedFormat; }
112
114 void setFileName(std::string passedFileName) {
115 fileName = std::move(passedFileName);
116 }
117
118 void setMetaData(const MetaDataType &passedMetaData) {
119 metaData = passedMetaData;
120 }
121
122 void addMetaData(const std::string &key, double value) {
123 metaData[key] = std::vector<double>{value};
124 }
125
126 void addMetaData(const std::string &key, const std::vector<double> &values) {
127 metaData[key] = values;
128 }
129
130 void addMetaData(const MetaDataType &newMetaData) {
131 for (const auto &pair : newMetaData) {
132 metaData[pair.first] = pair.second;
133 }
134 }
135
136 void apply() {
137 // check mesh
138 if (mesh == nullptr) {
139 Logger::getInstance()
140 .addError("No mesh was passed to VTKWriter.")
141 .print();
142 return;
143 }
144 // check filename
145 if (fileName.empty()) {
146 Logger::getInstance()
147 .addError("No file name specified for VTKWriter.")
148 .print();
149 return;
150 }
151 if (mesh->nodes.empty()) {
152 Logger::getInstance().addWarning("Writing empty mesh.").print();
153 return;
154 }
155
156 if (fileFormat == FileFormatEnum::VTK_AUTO) {
157 auto dotPos = fileName.rfind('.');
158 if (dotPos == std::string::npos) {
159 fileFormat = FileFormatEnum::VTP;
160 } else {
161 auto ending = fileName.substr(dotPos);
162 if (ending == ".vtk") {
163 fileFormat = FileFormatEnum::VTK_LEGACY;
164 } else if (ending == ".vtp") {
165 fileFormat = FileFormatEnum::VTP;
166 } else if (ending == ".vtu") {
167 fileFormat = FileFormatEnum::VTU;
168 } else {
169 Logger::getInstance()
170 .addError("No valid file format found based on the file ending "
171 "passed to VTKWriter.")
172 .print();
173 return;
174 }
175 }
176 }
177
178 // check file format
179 switch (fileFormat) {
181 writeVTKLegacy(fileName);
182 break;
183#ifdef VIENNALS_USE_VTK
185 writeVTP(fileName);
186 break;
188 writeVTU(fileName);
189 break;
190#else
193 Logger::getInstance()
194 .addWarning("VTKWriter was built without VTK support. Falling back "
195 "to VTK_LEGACY.")
196 .print();
197 writeVTKLegacy(fileName);
198 break;
199#endif
200 default:
201 Logger::getInstance()
202 .addError("No valid file format set for VTKWriter.")
203 .print();
204 }
205 }
206
207private:
208#ifdef VIENNALS_USE_VTK
209 void writeVTP(std::string filename) const {
210 if (mesh == nullptr) {
211 Logger::getInstance()
212 .addError("No mesh was passed to VTKWriter.")
213 .print();
214 return;
215 }
216
217 if (filename.find(".vtp") != filename.size() - 4)
218 filename += ".vtp";
219 vtkSmartPointer<vtkPolyData> polyData = vtkSmartPointer<vtkPolyData>::New();
220
221 // Points
222 vtkSmartPointer<vtkPoints> polyPoints = vtkSmartPointer<vtkPoints>::New();
223 for (auto it = mesh->getNodes().begin(); it != mesh->getNodes().end();
224 ++it) {
225 polyPoints->InsertNextPoint((*it)[0], (*it)[1], (*it)[2]);
226 }
227 polyData->SetPoints(polyPoints);
228
229 // Vertices
230 if (mesh->vertices.size() > 0) {
231 vtkSmartPointer<vtkCellArray> polyCells =
232 vtkSmartPointer<vtkCellArray>::New();
233 for (auto it = mesh->vertices.begin(); it != mesh->vertices.end(); ++it) {
234 polyCells->InsertNextCell(1);
235 polyCells->InsertCellPoint((*it)[0]);
236 }
237 polyData->SetVerts(polyCells);
238 }
239
240 // Lines
241 if (mesh->lines.size() > 0) {
242 vtkSmartPointer<vtkCellArray> polyCells =
243 vtkSmartPointer<vtkCellArray>::New();
244 for (auto it = mesh->lines.begin(); it != mesh->lines.end(); ++it) {
245 polyCells->InsertNextCell(2);
246 for (unsigned i = 0; i < 2; ++i) {
247 polyCells->InsertCellPoint((*it)[i]);
248 }
249 }
250 polyData->SetLines(polyCells);
251 }
252
253 // Triangles
254 if (mesh->triangles.size() > 0) {
255 vtkSmartPointer<vtkCellArray> polyCells =
256 vtkSmartPointer<vtkCellArray>::New();
257 for (auto it = mesh->triangles.begin(); it != mesh->triangles.end();
258 ++it) {
259 polyCells->InsertNextCell(3);
260 for (unsigned i = 0; i < 3; ++i) {
261 polyCells->InsertCellPoint((*it)[i]);
262 }
263 }
264 polyData->SetPolys(polyCells);
265 }
266
267 addDataFromMesh(mesh->pointData, polyData->GetPointData());
268 addDataFromMesh(mesh->cellData, polyData->GetCellData());
269 addMetaDataToVTK(polyData);
270
271 vtkSmartPointer<vtkXMLPolyDataWriter> pwriter =
272 vtkSmartPointer<vtkXMLPolyDataWriter>::New();
273 pwriter->SetFileName(filename.c_str());
274 pwriter->SetInputData(polyData);
275 pwriter->Write();
276 }
277
278 void writeVTU(std::string filename) const {
279 if (mesh == nullptr) {
280 Logger::getInstance()
281 .addError("No mesh was passed to VTKWriter.")
282 .print();
283 return;
284 }
285
286 if (filename.find(".vtu") != filename.size() - 4)
287 filename += ".vtu";
288
289 vtkSmartPointer<vtkUnstructuredGrid> uGrid =
290 vtkSmartPointer<vtkUnstructuredGrid>::New();
291
292 // Points
293 vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
294 for (auto it = mesh->getNodes().begin(); it != mesh->getNodes().end();
295 ++it) {
296 points->InsertNextPoint((*it)[0], (*it)[1], (*it)[2]);
297 }
298 uGrid->SetPoints(points);
299
300 // Now set all cells
301 vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
302 std::vector<int> cellTypes;
303 cellTypes.reserve(mesh->vertices.size() + mesh->lines.size() +
304 mesh->triangles.size() + mesh->tetras.size() +
305 mesh->hexas.size());
306
307 // Vertices
308 if (mesh->vertices.size() > 0) {
309 for (auto it = mesh->vertices.begin(); it != mesh->vertices.end(); ++it) {
310 cells->InsertNextCell(1);
311 cells->InsertCellPoint((*it)[0]);
312 cellTypes.push_back(1); // vtk Vertex
313 }
314 }
315
316 // Lines
317 if (mesh->lines.size() > 0) {
318 for (auto it = mesh->lines.begin(); it != mesh->lines.end(); ++it) {
319 cells->InsertNextCell(2);
320 for (unsigned i = 0; i < 2; ++i) {
321 cells->InsertCellPoint((*it)[i]);
322 }
323 cellTypes.push_back(3); // vtk Line
324 }
325 }
326
327 // Triangles
328 if (mesh->triangles.size() > 0) {
329 for (auto it = mesh->triangles.begin(); it != mesh->triangles.end();
330 ++it) {
331 cells->InsertNextCell(3);
332 for (unsigned i = 0; i < 3; ++i) {
333 cells->InsertCellPoint((*it)[i]);
334 }
335 cellTypes.push_back(5); // vtk Triangle
336 }
337 }
338
339 // Tetras
340 if (mesh->tetras.size() > 0) {
341 for (auto it = mesh->tetras.begin(); it != mesh->tetras.end(); ++it) {
342 cells->InsertNextCell(4);
343 for (unsigned i = 0; i < 4; ++i) {
344 cells->InsertCellPoint((*it)[i]);
345 }
346 cellTypes.push_back(10); // vtk Tetra
347 }
348 }
349
350 // Hexas
351 if (mesh->hexas.size() > 0) {
352 for (auto it = mesh->hexas.begin(); it != mesh->hexas.end(); ++it) {
353 cells->InsertNextCell(8);
354 for (unsigned i = 0; i < 8; ++i) {
355 cells->InsertCellPoint((*it)[i]);
356 }
357 cellTypes.push_back(12); // vtk Hexahedron
358 }
359 }
360
361 // set cells
362 uGrid->SetCells(&(cellTypes[0]), cells);
363
364 addDataFromMesh(mesh->pointData, uGrid->GetPointData());
365 addDataFromMesh(mesh->cellData, uGrid->GetCellData());
366 addMetaDataToVTK(uGrid);
367
368 // // now add pointData
369 // for (unsigned i = 0; i < mesh->cellData.getScalarDataSize(); ++i) {
370 // vtkSmartPointer<vtkFloatArray> pointData =
371 // vtkSmartPointer<vtkFloatArray>::New();
372 // pointData->SetNumberOfComponents(1);
373 // pointData->SetName(mesh->cellData.getScalarDataLabel(i).c_str());
374 // auto scalars = *(mesh->cellData.getScalarData(i));
375 // for (unsigned j = 0; j < scalars.size(); ++j) {
376 // pointData->InsertNextValue(scalars[j]);
377 // }
378 // uGrid->GetCellData()->AddArray(pointData);
379 // }
380
381 // // now add vector data
382 // for (unsigned i = 0; i < mesh->cellData.getVectorDataSize(); ++i) {
383 // vtkSmartPointer<vtkFloatArray> vectorData =
384 // vtkSmartPointer<vtkFloatArray>::New();
385 // vectorData->SetNumberOfComponents(3);
386 // vectorData->SetName(mesh->cellData.getVectorDataLabel(i).c_str());
387 // auto vectors = *(mesh->cellData.getVectorData(i));
388 // for (unsigned j = 0; j < vectors.size(); ++j) {
389 // vectorData->InsertNextTuple3(vectors[j][0], vectors[j][1],
390 // vectors[j][2]);
391 // }
392 // uGrid->GetCellData()->AddArray(vectorData);
393 // }
394
395 vtkSmartPointer<vtkXMLUnstructuredGridWriter> owriter =
396 vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
397 owriter->SetFileName(filename.c_str());
398 owriter->SetInputData(uGrid);
399 owriter->Write();
400 }
401
402#endif // VIENNALS_USE_VTK
403
404 void writeVTKLegacy(const std::string &filename) {
405 if (mesh == nullptr) {
406 Logger::getInstance()
407 .addError("No mesh was passed to VTKWriter.")
408 .print();
409 return;
410 }
411
412 std::ofstream f(filename.c_str());
413
414 f << "# vtk DataFile Version 2.0" << std::endl;
415 f << ((!mesh->lines.empty()) ? 2 : 3) << "D Surface" << std::endl;
416 f << "ASCII" << std::endl;
417 f << "DATASET UNSTRUCTURED_GRID" << std::endl;
418 f << "POINTS " << mesh->nodes.size() << " float" << std::endl;
419
420 // print nodes
421 for (unsigned int i = 0; i < mesh->nodes.size(); i++) {
422 for (int j = 0; j < 3; j++)
423 f << static_cast<float>(mesh->nodes[i][j]) << " ";
424 f << std::endl;
425 }
426
427 const unsigned numberOfCells = mesh->vertices.size() + mesh->lines.size() +
428 mesh->triangles.size() +
429 mesh->tetras.size() + mesh->hexas.size();
430 const unsigned cellDataSize =
431 2 * mesh->vertices.size() + 3 * mesh->lines.size() +
432 4 * mesh->triangles.size() + 5 * mesh->tetras.size() +
433 9 * mesh->hexas.size();
434
435 f << "CELLS " << numberOfCells << " " << cellDataSize << std::endl;
436
437 // print elements
438 for (unsigned int i = 0; i < mesh->vertices.size(); i++) {
439 f << 1 << " " << mesh->vertices[i][0] << std::endl;
440 }
441 for (unsigned int i = 0; i < mesh->lines.size(); i++) {
442 f << 2 << " ";
443 for (int j = 0; j < 2; j++)
444 f << mesh->lines[i][j] << " ";
445 f << std::endl;
446 }
447
448 for (unsigned int i = 0; i < mesh->triangles.size(); i++) {
449 f << 3 << " ";
450 for (int j = 0; j < 3; j++)
451 f << mesh->triangles[i][j] << " ";
452 f << std::endl;
453 }
454
455 for (unsigned int i = 0; i < mesh->tetras.size(); i++) {
456 f << 4 << " ";
457 for (int j = 0; j < 4; j++)
458 f << mesh->tetras[i][j] << " ";
459 f << std::endl;
460 }
461
462 for (unsigned int i = 0; i < mesh->hexas.size(); i++) {
463 f << 8 << " ";
464 for (int j = 0; j < 8; j++)
465 f << mesh->hexas[i][j] << " ";
466 f << std::endl;
467 }
468
469 f << "CELL_TYPES " << numberOfCells << std::endl;
470 for (unsigned i = 0; i < mesh->vertices.size(); ++i)
471 f << 1 << std::endl;
472
473 for (unsigned i = 0; i < mesh->lines.size(); ++i)
474 f << 3 << std::endl;
475
476 for (unsigned i = 0; i < mesh->triangles.size(); ++i)
477 f << 5 << std::endl;
478
479 for (unsigned i = 0; i < mesh->tetras.size(); ++i)
480 f << 10 << std::endl;
481
482 for (unsigned i = 0; i < mesh->hexas.size(); ++i)
483 f << 12 << std::endl;
484
485 // WRITE POINT DATA
486 if (mesh->pointData.getScalarDataSize() ||
487 mesh->pointData.getVectorDataSize()) {
488 Logger::getInstance()
489 .addWarning("Point data output not supported for legacy VTK output. "
490 "Point data is ignored.")
491 .print();
492 }
493
494 // WRITE SCALAR DATA
495 if (mesh->cellData.getScalarDataSize()) {
496 f << "CELL_DATA " << mesh->cellData.getScalarData(0)->size() << std::endl;
497 for (unsigned i = 0; i < mesh->cellData.getScalarDataSize(); ++i) {
498 auto scalars = *(mesh->cellData.getScalarData(i));
499 f << "SCALARS " << mesh->cellData.getScalarDataLabel(i) << " float"
500 << std::endl;
501 f << "LOOKUP_TABLE default" << std::endl;
502 for (unsigned j = 0; j < scalars.size(); ++j) {
503 f << ((std::abs(scalars[j]) < 1e-6) ? 0.0 : scalars[j]) << std::endl;
504 }
505 }
506 }
507
508 // WRITE VECTOR DATA
509 if (mesh->cellData.getVectorDataSize()) {
510 if (!mesh->cellData.getScalarDataSize())
511 f << "CELL_DATA " << mesh->cellData.getVectorData(0)->size()
512 << std::endl;
513 for (unsigned i = 0; i < mesh->cellData.getVectorDataSize(); ++i) {
514 auto vectors = *(mesh->cellData.getVectorData(i));
515 f << "VECTORS " << mesh->cellData.getVectorDataLabel(i) << " float"
516 << std::endl;
517 for (unsigned j = 0; j < vectors.size(); ++j) {
518 for (unsigned k = 0; k < 3; ++k) {
519 f << vectors[j][k] << " ";
520 }
521 f << std::endl;
522 }
523 }
524 }
525
526 f.close();
527 }
528};
529
530} // namespace viennals
This class holds an explicit mesh, which is always given in 3 dimensions. If it describes a 2D mesh,...
Definition lsMesh.hpp:22
void setFileName(std::string passedFileName)
set file name for file to write
Definition lsVTKWriter.hpp:114
VTKWriter(SmartPointer< Mesh< T > > passedMesh, FileFormatEnum passedFormat, std::string passedFileName)
Definition lsVTKWriter.hpp:103
VTKWriter(SmartPointer< Mesh< T > > passedMesh, std::string passedFileName)
Definition lsVTKWriter.hpp:100
void addMetaData(const std::string &key, double value)
Definition lsVTKWriter.hpp:122
void apply()
Definition lsVTKWriter.hpp:136
void setMetaData(const MetaDataType &passedMetaData)
Definition lsVTKWriter.hpp:118
void setMesh(SmartPointer< Mesh< T > > passedMesh)
Definition lsVTKWriter.hpp:108
void setFileFormat(FileFormatEnum passedFormat)
set file format for file to write. Defaults to VTK_LEGACY.
Definition lsVTKWriter.hpp:111
void addMetaData(const std::string &key, const std::vector< double > &values)
Definition lsVTKWriter.hpp:126
void addMetaData(const MetaDataType &newMetaData)
Definition lsVTKWriter.hpp:130
VTKWriter(SmartPointer< Mesh< T > > passedMesh)
Definition lsVTKWriter.hpp:98
Definition lsAdvect.hpp:36
FileFormatEnum
Definition lsFileFormats.hpp:4
@ VTK_AUTO
Definition lsFileFormats.hpp:8
@ VTK_LEGACY
Definition lsFileFormats.hpp:5
@ VTP
Definition lsFileFormats.hpp:6
@ VTU
Definition lsFileFormats.hpp:7