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