Mmg
Simplicial remeshers (mesh adaptation, isovalue discretization, lagrangian movement)
vtkparser.hpp
Go to the documentation of this file.
1/* =============================================================================
2** This file is part of the mmg software package for the tetrahedral
3** mesh modification.
4** Copyright (c) Bx INP/Inria/UBordeaux/UPMC, 2004- .
5**
6** mmg is free software: you can redistribute it and/or modify it
7** under the terms of the GNU Lesser General Public License as published
8** by the Free Software Foundation, either version 3 of the License, or
9** (at your option) any later version.
10**
11** mmg is distributed in the hope that it will be useful, but WITHOUT
12** ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13** FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14** License for more details.
15**
16** You should have received a copy of the GNU Lesser General Public
17** License and of the GNU General Public License along with mmg (in
18** files COPYING.LESSER and COPYING). If not, see
19** <http://www.gnu.org/licenses/>. Please read their terms carefully and
20** use this copy of the mmg distribution only if you accept them.
21** =============================================================================
22*/
23
24#ifndef VTKPARSER_HPP
25#define VTKPARSER_HPP
26
27#ifdef USE_VTK
28
29#include "mmgcommon_private.h"
30
31#include <vtkSmartPointer.h>
32#include <vtkXMLReader.h>
33#include <vtkXMLWriter.h>
34#include <vtkXMLUnstructuredGridReader.h>
35#include <vtkXMLUnstructuredGridWriter.h>
36#include <vtkXMLPUnstructuredGridWriter.h>
37#include <vtkXMLPolyDataReader.h>
38#include <vtkXMLPolyDataWriter.h>
39#include <vtkXMLPPolyDataWriter.h>
40#include <vtkDataSetReader.h>
41#include <vtkDataSetWriter.h>
42#include <vtkPDataSetWriter.h>
43#include <vtkDataSet.h>
44#include <vtkUnstructuredGrid.h>
45#include <vtkPolyData.h>
46#include <vtkStructuredGrid.h>
47#include <vtkPointData.h>
48#include <vtkCellData.h>
49#include <vtkFieldData.h>
50#include <vtkCellTypes.h>
51#include <vtkDataArray.h>
52#include <vtkFloatArray.h>
53#include <vtkDoubleArray.h>
54
55#include <vtkLine.h>
56#include <vtkTriangle.h>
57#include <vtkQuad.h>
58#include <vtkTetra.h>
59#include <vtkWedge.h>
60#include <vtkCellArray.h>
61#include <typeinfo>
62
63int MMG5_loadVtpMesh_part1(MMG5_pMesh,const char*,vtkDataSet**,int8_t*,int8_t*,int*,int8_t*,int8_t*);
64int MMG5_loadVtuMesh_part1(MMG5_pMesh,const char*,vtkDataSet**,int8_t*,int8_t*,int*,int8_t*,int8_t*);
65int MMG5_loadVtkMesh_part1(MMG5_pMesh,const char*,vtkDataSet**,int8_t*,int8_t*,int*,int8_t*,int8_t*);
66
67int MMG5_loadVtkMesh_part2(MMG5_pMesh,MMG5_pSol*,vtkDataSet**,int8_t,int8_t,int,int8_t,int8_t);
68
76static void MMG5_internal_VTKSetLine(vtkSmartPointer<vtkPolyData> d,vtkSmartPointer<vtkCellArray> *ca) {
77
78 d->SetLines(*ca);
79
80 *ca = NULL;
81 *ca = vtkSmartPointer< vtkCellArray >::New();
82
83}
84
90static void MMG5_internal_VTKSetLine(vtkSmartPointer<vtkUnstructuredGrid> d,vtkSmartPointer<vtkCellArray> *ca) {
91
92}
93
94
101static void MMG5_internal_VTKSetCells(int * t,vtkSmartPointer<vtkPolyData> d,vtkSmartPointer<vtkCellArray> ca) {
102 d->SetPolys(ca);
103}
104
111static void MMG5_internal_VTKSetCells(int * t,vtkSmartPointer<vtkUnstructuredGrid> d,vtkSmartPointer<vtkCellArray> ca) {
112 d->SetCells(t,ca);
113}
114
120static void MMG5_internal_VTKbinary(vtkXMLUnstructuredGridWriter *w, int binary) {
121 if ( binary ) {
122 w->SetDataModeToBinary();
123 }
124 else {
125 w->SetDataModeToAscii();
126 }
127}
128
134static void MMG5_internal_VTKbinary(vtkXMLPolyDataWriter *w, int binary) {
135 if ( binary ) {
136 w->SetDataModeToBinary();
137 }
138 else {
139 w->SetDataModeToAscii();
140 }
141}
142
148static void MMG5_internal_VTKbinary(vtkDataSetWriter *w, int binary) {
149 return;
150}
151
169template <class T, class TWriter, class PWriter>
171 const char *mfilename,int metricData,
172 int binary,int npart,int myid,int master) {
173 int hasPointRef = 0, hasCellRef = 0;
174
175 // Transfer points from Mmg to VTK
176 vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
177
178 MMG5_int np = 0;
179 for ( MMG5_int k=1; k<=mesh->np; ++k ) {
180 MMG5_pPoint ppt = &mesh->point[k];
181 if ( !MG_VOK(ppt) ) continue;
182
183 ppt->tmp = np++;
184
185 if ( mesh->dim == 2 ) ppt->c[2] = 0.;
186 if ( ppt->ref ) hasPointRef = 1;
187
188 points->InsertNextPoint(ppt->c[0], ppt->c[1], ppt->c[2]);
189 }
190
191
192 auto dataset = vtkSmartPointer<T> :: New();
193 dataset->SetPoints(points);
194
195
196 // VTK types
197 vtkSmartPointer<vtkCellArray> cellArray = vtkSmartPointer<vtkCellArray> ::New();
198 vtkSmartPointer<vtkWedge> wedge = vtkSmartPointer<vtkWedge> ::New();
199 vtkSmartPointer<vtkTetra> tetra = vtkSmartPointer<vtkTetra> ::New();
200 vtkSmartPointer<vtkTriangle> tria = vtkSmartPointer<vtkTriangle> ::New();
201 vtkSmartPointer<vtkQuad> quadra = vtkSmartPointer<vtkQuad> ::New();
202 vtkSmartPointer<vtkLine> edge = vtkSmartPointer<vtkLine> ::New();
203
204 MMG5_int nc = mesh->na+mesh->nt+mesh->nquad+mesh->ne+mesh->nprism;
205 MMG5_int ic = 0;
206
207 int* types = NULL;
208 MMG5_SAFE_MALLOC ( types, nc, int, return 0 );
209
210 // Transfer edges from Mmg to VTK
211 for ( MMG5_int k=1; k<=mesh->na; ++k ) {
212 MMG5_pEdge pa = &mesh->edge[k];
213 if ( !pa || !pa->a ) continue;
214
215 if ( pa->ref ) hasCellRef = 1;
216
217 edge->GetPointIds()->SetId(0, mesh->point[pa->a].tmp);
218 edge->GetPointIds()->SetId(1, mesh->point[pa->b].tmp);
219
220 cellArray->InsertNextCell(edge);
221 types[ic++] = VTK_LINE;
222 }
223
224 MMG5_internal_VTKSetLine(dataset, &cellArray);
225
226 // Transfer trias from Mmg to VTK
227 for ( MMG5_int k=1; k<=mesh->nt; ++k ) {
228 MMG5_pTria ptt = &mesh->tria[k];
229 if ( !MG_EOK(ptt) ) continue;
230
231 if ( ptt->ref ) hasCellRef = 1;
232
233 for ( int i=0; i<3; ++i ) {
234 tria->GetPointIds()->SetId(i, mesh->point[ptt->v[i]].tmp );
235 }
236 cellArray->InsertNextCell(tria);
237 types[ic++] = VTK_TRIANGLE;
238 }
239
240 // Transfer quads from Mmg to VTK
241 for ( MMG5_int k=1; k<=mesh->nquad; ++k ) {
242 MMG5_pQuad pq = &mesh->quadra[k];
243 if ( !MG_EOK(pq) ) continue;
244
245 if ( pq->ref ) hasCellRef = 1;
246
247 for ( int i=0; i<4; ++i ) {
248 quadra->GetPointIds()->SetId(i, mesh->point[pq->v[i]].tmp );
249 }
250 cellArray->InsertNextCell(quadra);
251 types[ic++] = VTK_QUAD;
252 }
253
254 // Transfer tets from Mmg to VTK
255 for ( MMG5_int k=1; k<=mesh->ne; ++k ) {
256 MMG5_pTetra pt = &mesh->tetra[k];
257 if ( !MG_EOK(pt) ) continue;
258
259 if ( pt->ref ) hasCellRef = 1;
260
261 for ( int i=0; i<4; ++i ) {
262 tetra->GetPointIds()->SetId(i, mesh->point[pt->v[i]].tmp );
263 }
264 cellArray->InsertNextCell(tetra);
265 types[ic++] = VTK_TETRA;
266 }
267
268 // Transfer prisms from Mmg to VTK
269 for ( MMG5_int k=1; k<=mesh->nprism; ++k ) {
270 MMG5_pPrism ppr = &mesh->prism[k];
271 if ( !MG_EOK(ppr) ) continue;
272
273 if ( ppr->ref ) hasCellRef = 1;
274
275 for ( int i=0; i<6; ++i ) {
276 wedge->GetPointIds()->SetId(i, mesh->point[ppr->v[i]].tmp );
277 }
278 cellArray->InsertNextCell(wedge);
279 types[ic++] = VTK_WEDGE;
280 }
281
282 assert ( ic == nc );
283
284 // Transfer cell array into data set
285 MMG5_internal_VTKSetCells(types,dataset,cellArray);
286
287 // Transfer references if needed (i.e. the mesh contains non 0 refs)
288 if ( hasPointRef ) {
289 auto *ar = vtkIntArray::New();
290
291 ar->SetNumberOfComponents(1);
292 ar->SetNumberOfTuples(mesh->np);
293 ar->SetName("medit:ref");
294
295 for ( MMG5_int k = 0; k < mesh->np; k++ ) {
296 MMG5_pPoint ppt = &mesh->point[k+1];
297 if ( !MG_VOK(ppt) ) continue;
298
299 ar->SetTuple1(ppt->tmp,ppt->ref);
300 }
301
302 dataset->GetPointData()->AddArray(ar);
303 }
304 if ( hasCellRef ) {
305 auto *ar = vtkIntArray::New();
306
307 ar->SetNumberOfComponents(1);
308 ar->SetNumberOfTuples(nc);
309 ar->SetName("medit:ref");
310
311 MMG5_int idx = 0;
312 for ( MMG5_int k = 1; k <= mesh->na; k++ ) {
313 if ( !mesh->edge[k].a ) continue;
314 ar->SetTuple1(idx++,mesh->edge[k].ref);
315 }
316 for ( MMG5_int k = 1; k <= mesh->nt; k++ ) {
317 if ( !MG_EOK(&mesh->tria[k]) ) continue;
318 ar->SetTuple1(idx++,mesh->tria[k].ref);
319 }
320 for ( MMG5_int k = 1; k <= mesh->nquad; k++ ) {
321 if ( !MG_EOK(&mesh->quadra[k]) ) continue;
322 ar->SetTuple1(idx++,mesh->quadra[k].ref);
323 }
324 for ( MMG5_int k = 1; k <= mesh->ne; k++ ) {
325 if ( !MG_EOK(&mesh->tetra[k]) ) continue;
326 ar->SetTuple1(idx++,mesh->tetra[k].ref);
327 }
328 for ( MMG5_int k = 1; k <= mesh->nprism; k++ ) {
329 if ( !MG_EOK(&mesh->prism[k]) ) continue;
330 ar->SetTuple1(idx++,mesh->prism[k].ref);
331 }
332
333 dataset->GetCellData()->AddArray(ar);
334 }
335
336 // Transfer point solutions into data set
337 MMG5_pSol psl = NULL;
338 int fieldData = 0;
339 int nsols = 0;
340
341 // For functions: MMG<Y>_saveVt<X>Mesh;
342 // PMMG_savePvtuMesh and PMMG_savePvtuMesh_and_allData
343 if ( metricData==1 ) {
344 if ( sol && *sol && sol[0]->np ) {
345 nsols = 1;
346 }
347 }
348
349 // If we have other fields to output
350 // For functions: MMG<Y>_saveVt<X>Mesh_and_allData and PMMG_savePvtuMesh_and_allData
351 if (&sol[1][0]) {
352 nsols += mesh->nsols;
353 fieldData = 1;
354 }
355
356 static int mmgWarn = 0;
357 for ( int isol=0; isol<nsols; ++isol) {
358 //- If metric, it is stored in sol[0]
359 if (metricData && isol==0) {
360 psl = sol[0];
361 }
362 //- If metric and other fields, fields are stored in sol[1]
363 else if (metricData){
364 psl = &sol[1][isol-1];
365 }
366 //- If other fields only, fields are stored in sol[1]
367 else {
368 psl = &sol[1][isol];
369 }
370
371 if ( !psl->m ) {
372 if ( !mmgWarn ) {
373 mmgWarn = 1;
374 fprintf(stderr, " ## Warning: %s: missing data for at least 1 solution."
375 " Skipped.\n",__func__);
376 }
377 continue;
378 }
379
380 int ncp;
381 if ( psl->size == 1 ) {
382 ncp = 1;
383 }
384 else if ( psl->size == psl->dim ) {
385 ncp = psl->dim;
386 }
387 else {
388 ncp = psl->dim*psl->dim;
389 }
390
391 auto *ar = vtkDoubleArray::New();
392
393 ar->SetNumberOfComponents(ncp);
394 ar->SetNumberOfTuples(mesh->np);
395
396 if ( psl->namein ) {
397 char *tmp = MMG5_Get_basename(psl->namein);
398 char *data, *data2;
399
400 MMG5_SAFE_CALLOC(data,strlen(tmp)+8,char,
401 MMG5_SAFE_FREE ( types ); return 0);
402 MMG5_SAFE_CALLOC(data2,strlen(tmp)+8,char,
403 MMG5_SAFE_FREE ( types ); return 0);
404
405 strcpy(data,tmp);
406 free(tmp); tmp = 0;
407
408 if ( metricData && strstr(data,":ls")) {
409 memmove(data2,data,strlen(data)-3);
410 strcat ( data2 , ":metric");
411 ar->SetName(data2);
412 }
413 else if (metricData && !strstr(data,":metric") && isol==0) {
414 strcat ( data , ":metric");
415 ar->SetName(data);
416 }
417 else {
418 ar->SetName(data);
419 }
420
421 MMG5_DEL_MEM(mesh,data);
422 }
423 else {
424 ar->SetName("no_name");
425 }
426
427 double* dfmt = NULL;
428 MMG5_SAFE_MALLOC ( dfmt, ncp, double, MMG5_SAFE_FREE ( types );return 0 );
429
430 if ( psl->size!= (psl->dim*(psl->dim+1))/2 ) {
431 /* scalar or vector field: data order and size isn't modified */
432 for ( MMG5_int k=1; k<=mesh->np; k++) {
433 MMG5_pPoint ppt = &mesh->point[k];
434 if ( !MG_VOK(ppt) ) continue;
435
436 MMG5_int iadr = k*psl->size;
437 for ( int i=0; i<psl->size; ++i ) {
438 dfmt[i] = psl->m[iadr+i];
439 }
440
441 if ( psl->dim==2 && ncp==3 ) {
442 dfmt[2] = 0; // z-component for a vector field
443 }
444 ar->SetTuple(ppt->tmp,dfmt);
445 }
446 }
447 else {
448 /* Tensor field: we need to reconstruct the entire symetric matrix */
449 for ( MMG5_int k=1; k<=mesh->np; k++) {
450 MMG5_pPoint ppt = &mesh->point[k];
451 if ( !MG_VOK(ppt) ) continue;
452
453 MMG5_int iadr = k*psl->size;
454 double *d = &psl->m[iadr];
455 if ( psl->dim==2 ) {
456 dfmt[0] = d[0];
457 dfmt[1] = d[1];
458 dfmt[2] = d[1];
459 dfmt[3] = d[2];
460 }
461 else {
462 double dbuf[6] = {0,0,0,0,0,0};
463 if ( metricData ) {
464 assert(!mesh->nsols);
465 MMG5_build3DMetric(mesh,psl,k,dbuf);
466 }
467 else {
468 for ( int i=0; i<psl->size; i++) dbuf[i] = psl->m[psl->size*k+i];
469 }
470 dfmt[0] = dbuf[0];
471 dfmt[1] = dbuf[1];
472 dfmt[2] = dbuf[2];
473 dfmt[3] = dbuf[1];
474 dfmt[4] = dbuf[3];
475 dfmt[5] = dbuf[4];
476 dfmt[6] = dbuf[2];
477 dfmt[7] = dbuf[4];
478 dfmt[8] = dbuf[5];
479 }
480 ar->SetTuple(ppt->tmp,dfmt);
481 }
482 }
483 dataset->GetPointData()->AddArray(ar);
484
485 MMG5_SAFE_FREE(dfmt);
486 }
487
488 if ( npart ) {
489 // distributed IO
490 vtkSmartPointer<PWriter> writer = vtkSmartPointer<PWriter>::New();
491
492#if VTK_MAJOR_VERSION <= 5
493 writer->SetInput(dataset);
494#else
495 writer->SetInputData(dataset);
496#endif
497
498 writer->SetFileName(mfilename);
499
500 writer->SetNumberOfPieces(npart);
501 writer->SetStartPiece(myid);
502 writer->SetEndPiece(myid);
503 writer->Write();
504 }
505 else {
506 // centralized IO
507 vtkSmartPointer<TWriter> writer = vtkSmartPointer<TWriter>::New();
508
509 writer->SetFileName(mfilename);
510
511#if VTK_MAJOR_VERSION <= 5
512 writer->SetInput(dataset);
513#else
514 writer->SetInputData(dataset);
515#endif
516
517 //MMG5_internal_VTKbinary(writer,binary);
518 writer->Write();
519 }
520
521 MMG5_SAFE_FREE(types);
522
523 return 1;
524}
525
538template <class T, class TWriter,class PWriter>
540 int metricData,int binary) {
541
542 return MMG5_saveVtkMesh_i<T,TWriter,PWriter>( mesh,sol,filename,
543 metricData,binary,0,0,0);
544
545}
546
547#endif
548
549#endif
char * MMG5_Get_basename(char *path)
tmp[*strlen0]
MMG5_pMesh MMG5_pSol * sol
MMG5_pMesh * mesh
MMG5_pMesh char * filename
void MMG5_build3DMetric(MMG5_pMesh mesh, MMG5_pSol sol, MMG5_int ip, double dbuf[6])
Definition: inout.c:1524
#define MMG5_SAFE_CALLOC(ptr, size, type, law)
#define MG_EOK(pt)
#define MG_VOK(ppt)
#define MMG5_SAFE_MALLOC(ptr, size, type, law)
#define MMG5_SAFE_FREE(ptr)
#define MMG5_DEL_MEM(mesh, ptr)
Structure to store edges of am MMG mesh.
Definition: libmmgtypes.h:311
MMG5_int b
Definition: libmmgtypes.h:312
MMG5_int ref
Definition: libmmgtypes.h:313
MMG5_int a
Definition: libmmgtypes.h:312
MMG mesh structure.
Definition: libmmgtypes.h:613
MMG5_pQuad quadra
Definition: libmmgtypes.h:656
MMG5_int ne
Definition: libmmgtypes.h:620
MMG5_pPoint point
Definition: libmmgtypes.h:649
MMG5_pPrism prism
Definition: libmmgtypes.h:653
MMG5_int nquad
Definition: libmmgtypes.h:621
MMG5_pTetra tetra
Definition: libmmgtypes.h:651
MMG5_int nt
Definition: libmmgtypes.h:620
MMG5_pTria tria
Definition: libmmgtypes.h:655
MMG5_int np
Definition: libmmgtypes.h:620
MMG5_pEdge edge
Definition: libmmgtypes.h:657
MMG5_int nprism
Definition: libmmgtypes.h:621
MMG5_int na
Definition: libmmgtypes.h:620
Structure to store vertices of an MMG mesh.
Definition: libmmgtypes.h:276
MMG5_int tmp
Definition: libmmgtypes.h:286
double c[3]
Definition: libmmgtypes.h:277
MMG5_int ref
Definition: libmmgtypes.h:284
Structure to store prsim of a MMG mesh.
Definition: libmmgtypes.h:469
MMG5_int v[6]
Definition: libmmgtypes.h:470
MMG5_int ref
Definition: libmmgtypes.h:471
Structure to store quadrangles of an MMG mesh.
Definition: libmmgtypes.h:372
MMG5_int ref
Definition: libmmgtypes.h:374
MMG5_int v[4]
Definition: libmmgtypes.h:373
char * namein
Definition: libmmgtypes.h:682
double * m
Definition: libmmgtypes.h:680
Structure to store tetrahedra of an MMG mesh.
Definition: libmmgtypes.h:407
MMG5_int v[4]
Definition: libmmgtypes.h:409
MMG5_int ref
Definition: libmmgtypes.h:410
Structure to store triangles of a MMG mesh.
Definition: libmmgtypes.h:338
MMG5_int ref
Definition: libmmgtypes.h:341
MMG5_int v[3]
Definition: libmmgtypes.h:340
int MMG5_loadVtuMesh_part1(MMG5_pMesh, const char *, vtkDataSet **, int8_t *, int8_t *, int *, int8_t *, int8_t *)
Definition: vtkparser.cpp:325
static void MMG5_internal_VTKSetLine(vtkSmartPointer< vtkPolyData > d, vtkSmartPointer< vtkCellArray > *ca)
Definition: vtkparser.hpp:76
int MMG5_loadVtpMesh_part1(MMG5_pMesh, const char *, vtkDataSet **, int8_t *, int8_t *, int *, int8_t *, int8_t *)
Definition: vtkparser.cpp:245
int MMG5_saveVtkMesh(MMG5_pMesh mesh, MMG5_pSol *sol, const char *filename, int metricData, int binary)
Definition: vtkparser.hpp:539
static void MMG5_internal_VTKbinary(vtkXMLUnstructuredGridWriter *w, int binary)
Definition: vtkparser.hpp:120
int MMG5_saveVtkMesh_i(MMG5_pMesh mesh, MMG5_pSol *sol, const char *mfilename, int metricData, int binary, int npart, int myid, int master)
Definition: vtkparser.hpp:170
int MMG5_loadVtkMesh_part2(MMG5_pMesh, MMG5_pSol *, vtkDataSet **, int8_t, int8_t, int, int8_t, int8_t)
Definition: vtkparser.cpp:366
static void MMG5_internal_VTKSetCells(int *t, vtkSmartPointer< vtkPolyData > d, vtkSmartPointer< vtkCellArray > ca)
Definition: vtkparser.hpp:101
int MMG5_loadVtkMesh_part1(MMG5_pMesh, const char *, vtkDataSet **, int8_t *, int8_t *, int *, int8_t *, int8_t *)
Definition: vtkparser.cpp:285