1 // Copyright (C) 2021 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "MeshFormatReader.hxx"
22 #include "MEDFileMesh.hxx"
23 #include "MEDFileField.hxx"
24 #include "MEDFileData.hxx"
25 #include "MEDCouplingFieldDouble.hxx"
26 #include "libmesh5.hxx"
27 #include "MEDMESHConverterUtilities.hxx"
31 namespace MEDCoupling {
32 MeshFormatReader::MeshFormatReader():_myMeshName("MESH")
34 MeshFormatReader::MeshFormatReader(const std::string& meshFileName,
35 const std::vector<std::string>& fieldFileName):_myFile(meshFileName),
36 _myFieldFileNames(fieldFileName),
41 MeshFormatReader::~MeshFormatReader()
44 MEDCoupling::MCAuto<MEDCoupling::MEDFileData> MeshFormatReader::loadInMedFileDS()
46 _myStatus = perform();
47 if(_myStatus != MeshFormat::DRS_OK) return 0;
49 if ( !_uMesh->getName().c_str() || strlen( _uMesh->getName().c_str() ) == 0 )
50 _uMesh->setName( _myMeshName );
51 if ( _myFieldFileNames.size() ) performFields();
54 MEDCoupling::MCAuto< MEDCoupling::MEDFileMeshes > meshes = MEDCoupling::MEDFileMeshes::New();
55 _myMed = MEDCoupling::MEDFileData::New();
56 meshes->pushMesh( _uMesh );
57 _myMed->setMeshes( meshes );
59 if ( _fields ) _myMed->setFields( _fields );
64 //================================================================================
66 * \brief Read a GMF file
68 //================================================================================
70 MeshFormat::Status MeshFormatReader::perform()
72 MeshFormat::Localizer loc;
74 MeshFormat::Status status = MeshFormat::DRS_OK;
78 _reader = MeshFormat::MeshFormatParser();
79 _myCurrentOpenFile = _myFile;
80 _myCurrentFileId = _reader.GmfOpenMesh( _myFile.c_str(), GmfRead, &_version, &_dim );
81 if ( !_myCurrentFileId )
83 if ( MeshFormat::isMeshExtensionCorrect( _myFile ))
86 return addMessage( MeshFormat::Comment("Can't open for reading ") << _myFile,
90 return addMessage( MeshFormat::Comment("Not '.mesh' or '.meshb' extension of file ") << _myFile,
98 MEDCoupling::DataArrayDouble* coordArray = MEDCoupling::DataArrayDouble::New();
102 int nbEdges = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfEdges);
103 int nbTria = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfTriangles);
104 int nbQuad = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfQuadrilaterals);
105 int nbTet = _reader.GmfStatKwd( _myCurrentFileId, MeshFormat::GmfTetrahedra );
106 int nbPyr = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfPyramids);
107 int nbHex = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfHexahedra);
108 int nbPrism = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfPrisms);
111 _dim2NbEl = nbTria + nbQuad;
112 _dim3NbEl = nbTet + nbPyr + nbHex + nbPrism;
113 bool okdim1 = (nbEdges > 0), okdim2= (_dim2NbEl > 0), okdim3= (_dim3NbEl > 0);
115 MEDCoupling::MEDCouplingUMesh* dimMesh1;
116 MEDCoupling::MEDCouplingUMesh* dimMesh2;
117 MEDCoupling::MEDCouplingUMesh* dimMesh3;
121 dimMesh1 = MEDCoupling::MEDCouplingUMesh::New();
122 dimMesh1->setCoords( coordArray );
123 dimMesh1->allocateCells( nbEdges);
124 dimMesh1->setMeshDimension( 1 );
128 dimMesh2 = MEDCoupling::MEDCouplingUMesh::New();
129 dimMesh2->setCoords( coordArray );
130 dimMesh2->allocateCells(_dim2NbEl);
131 dimMesh2->setMeshDimension( 2 );
135 dimMesh3 = MEDCoupling::MEDCouplingUMesh::New();
136 dimMesh3->setCoords( coordArray );
137 dimMesh3->allocateCells(_dim3NbEl);
138 dimMesh3->setMeshDimension( 3 );
143 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
149 setEdges(dimMesh1, nbEdges);
156 setTriangles(dimMesh2, nbTria);
159 /* Read quadrangles */
162 setQuadrangles( dimMesh2, nbQuad);
166 dimMesh2->finishInsertingCells();
167 _uMesh->setMeshAtLevel( 2 - _dim, dimMesh2 );
168 dimMesh2->sortCellsInMEDFileFrmt();
174 setTetrahedras( dimMesh3, nbTet);
180 setPyramids( dimMesh3, nbPyr);
186 setHexahedras(dimMesh3, nbHex);
190 //~const int prismIDShift = myMesh->GetMeshInfo().NbElements();
193 setPrisms(dimMesh3, nbPrism);
197 dimMesh3->finishInsertingCells();
198 _uMesh->setMeshAtLevel( 3 - _dim, dimMesh3 );
203 coordArray->decrRef();
204 _reader.GmfCloseMesh(_myCurrentFileId);
205 _myCurrentFileId = -1;
206 _myCurrentOpenFile ="";
210 MeshFormat::Status MeshFormatReader::performFields()
213 MeshFormat::Status status = MeshFormat::DRS_OK;
214 _fields = MEDCoupling::MEDFileFields::New();
217 const MeshFormat::GmfKwdCod meshFormatSol[1] = { MeshFormat::GmfSolAtVertices}; // ONLY VERTEX SOL FOR NOW
220 std::vector<std::string>::const_iterator fieldFileIt = _myFieldFileNames.begin();
222 for (; fieldFileIt !=_myFieldFileNames.end(); ++fieldFileIt)
224 _myCurrentOpenFile = *fieldFileIt;
225 _myCurrentFileId = _reader.GmfOpenMesh( fieldFileIt->c_str(), GmfRead, &version, &dim );
226 if ( !_myCurrentFileId )
228 if ( MeshFormat::isMeshExtensionCorrect( *fieldFileIt ))
231 return addMessage( MeshFormat::Comment("Can't open for reading ") << *fieldFileIt,
235 return addMessage( MeshFormat::Comment("Not '.sol' or '.solb' extension of file ") << _myFile,
239 if (version != _version)
242 addMessage( MeshFormat::Comment("Warning sol file version is different than mesh file version ") << *fieldFileIt,
248 return addMessage( MeshFormat::Comment("Error sol file must have same dimension as mesh file ") << *fieldFileIt,
253 MeshFormat::GmfKwdCod kwd = meshFormatSol[0];
254 int NmbSol, NmbTypes, NmbReals, TypesTab[ GmfMaxTyp ];
255 NmbSol = _reader.GmfStatKwd( _myCurrentFileId, kwd, &NmbTypes, &NmbReals, TypesTab );
259 for(int i=0; i<NmbTypes; i++)
261 _reader.GmfGotoKwd(_myCurrentFileId, kwd);
267 setFields(kwd, NmbSol, 1);
273 setFields(kwd, NmbSol, dim);
279 setFields(kwd, NmbSol, dim*(dim+1)/2 );
285 setFields(kwd, NmbSol, dim*dim);
298 _reader.GmfCloseMesh(_myCurrentFileId);
299 _myCurrentFileId = -1;
300 _myCurrentOpenFile ="";
305 void MeshFormatReader::setFields( MeshFormat::GmfKwdCod kwd, int nmbSol, int nbComp)
308 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> fieldValues = MEDCoupling::DataArrayDouble::New();
309 fieldValues->alloc(nmbSol, nbComp);
310 double* values = fieldValues->getPointer();
313 double *val = new double[nbComp];
315 bool isOnAll = (_uMesh->getNumberOfNodes()== nmbSol);
318 for(int i = 1; i<= nmbSol; i++)
320 callParserGetLin(kwd, val, nbComp, &ref);
322 std::copy(val, val+nbComp, values);
326 values=fieldValues->getPointer(); // to delete
330 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> timeStamp;
331 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> tsField = MEDCoupling::MEDFileFieldMultiTS::New();
332 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingUMesh> dimMesh ;
334 MEDCoupling::TypeOfField typeOfField;
335 setTypeOfFieldAndDimRel(kwd, &typeOfField, &dimRel );
337 timeStamp = MEDCoupling::MEDCouplingFieldDouble::New(typeOfField);
338 dimMesh = _uMesh->getMeshAtLevel(dimRel);
341 timeStamp->setMesh( dimMesh );
342 std::string name = "Field_on_Vertex";
343 timeStamp->setName(name);
344 timeStamp->setArray(fieldValues);
346 const int nbTS = tsField->getNumberOfTS();
348 timeStamp->setOrder( nbTS );
350 // add the time-stamp
351 timeStamp->checkConsistencyLight();
355 tsField->appendFieldNoProfileSBT( timeStamp );
359 _fields->pushField( tsField);
363 void MeshFormatReader::setMeshName(const std::string& theMeshName)
365 _myMeshName = theMeshName;
368 std::string MeshFormatReader::getMeshName() const
374 void MeshFormatReader::setFile(const std::string& theFileName)
376 _myFile = theFileName;
379 void MeshFormatReader::setFieldFileNames(const std::vector<std::string>& theFieldFileNames)
381 _myFieldFileNames = theFieldFileNames;
383 std::vector<std::string> MeshFormatReader::getFieldFileNames() const
385 return _myFieldFileNames;
388 MeshFormat::Status MeshFormatReader::addMessage(const std::string& msg,
389 const bool isFatal/*=false*/)
392 _myErrorMessages.clear(); // warnings are useless if a fatal error encounters
394 _myErrorMessages.push_back( msg );
398 std::cout << msg << std::endl;
400 return ( _myStatus = isFatal ? MeshFormat::DRS_FAIL : MeshFormat::DRS_WARN_SKIP_ELEM );
403 void MeshFormatReader::backward_shift(mcIdType* tab, int size)
405 for(int i = 0; i<size; i++) tab[i] = tab[i]-1;
409 MeshFormat::Status MeshFormatReader::setNodes( MEDCoupling::DataArrayDouble* coordArray)
412 MeshFormat::Status status = MeshFormat::DRS_OK;
413 int nbNodes = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfVertices);
415 return addMessage( "No nodes in the mesh", /*fatal=*/true );
417 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfVertices);
420 _uMesh = MEDCoupling::MEDFileUMesh::New();
421 coordArray->alloc( nbNodes, _dim );
422 double* coordPrt = coordArray->getPointer();
423 std::vector<double> nCoords (_dim+2, 0.0) ;
425 double* coordPointer = &nCoords[0];
427 if ( _version != GmfFloat )
431 for ( int i = 1; i <= nbNodes; ++i )
433 _dim == 2 ? _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfVertices, &x, &y, &ref) : \
434 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfVertices, &x, &y, &z, &ref);
438 std::copy(coordPointer, coordPointer+_dim, coordPrt);
439 MeshFormatElement e(MeshFormat::GmfVertices, i-1);
440 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 1);
450 for ( int i = 1; i <= nbNodes; ++i )
452 _dim == 2 ? _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfVertices, &x, &y, &ref) :
453 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfVertices, &x, &y, &z, &ref);
457 std::copy(coordPointer, coordPointer+_dim, coordPrt);
458 MeshFormatElement e(MeshFormat::GmfVertices, i-1);
459 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 1);
463 _uMesh->setCoords( coordArray );
470 void MeshFormatReader::setEdges( MEDCoupling::MEDCouplingUMesh* dimMesh1, int nbEdges)
473 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
475 // read extra vertices for quadratic edges
476 std::vector<int> quadNodesAtEdges( nbEdges + 1, -1 );
477 if ( int nbQuadEdges = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtEdges))
479 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtEdges);
480 for ( int i = 1; i <= nbQuadEdges; ++i )
482 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtEdges, &iN[0], &iN[1], &iN[2]);
484 quadNodesAtEdges[ iN[0] ] = iN[2];
488 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfEdges);
489 for ( int i = 1; i <= nbEdges; ++i )
491 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfEdges, &iN[0], &iN[1], &ref);
492 const int midN = quadNodesAtEdges[ i ];
496 mcIdType nodalConnPerCell[3] = {iN[0], iN[1], midN};
497 backward_shift(nodalConnPerCell, 3);
498 dimMesh1->insertNextCell(INTERP_KERNEL::NORM_SEG3, 3,nodalConnPerCell);
499 MeshFormatElement e(MeshFormat::GmfEdges, i-1);
500 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 1 -_dim);
505 mcIdType nodalConnPerCell[2] = {iN[0], iN[1]};
506 backward_shift(nodalConnPerCell, 2);
507 dimMesh1->insertNextCell(INTERP_KERNEL::NORM_SEG2, 2,nodalConnPerCell);
508 MeshFormatElement e(MeshFormat::GmfEdges, i-1);
509 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 1 -_dim);
513 dimMesh1->finishInsertingCells();
514 dimMesh1->sortCellsInMEDFileFrmt();
515 _uMesh->setMeshAtLevel( 1 - _dim, dimMesh1 );
518 void MeshFormatReader::setTriangles(MEDCoupling::MEDCouplingUMesh* dimMesh2, int nbTria)
520 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
522 // read extra vertices for quadratic triangles
523 std::vector< std::vector<int> > quadNodesAtTriangles( nbTria + 1 );
524 if ( int nbQuadTria = _reader.GmfStatKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtTriangles ))
526 _reader.GmfGotoKwd( _myCurrentFileId, MeshFormat::GmfExtraVerticesAtTriangles );
527 for ( int i = 1; i <= nbQuadTria; ++i )
529 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtTriangles,
530 &iN[0], &iN[1], &iN[2], &iN[3], &iN[4],
531 &iN[5]); // iN[5] - preview TRIA7
532 if ( iN[0] <= nbTria )
534 std::vector<int>& nodes = quadNodesAtTriangles[ iN[0] ];
535 nodes.insert( nodes.end(), & iN[2], & iN[5+1] );
536 nodes.resize( iN[1] );
543 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfTriangles);
545 for ( int i = 1; i <= nbTria; ++i )
547 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfTriangles, &iN[0], &iN[1], &iN[2], &ref);
548 std::vector<int>& midN = quadNodesAtTriangles[ i ];
549 if ( midN.size() >= 3 )
552 mcIdType nodalConnPerCell[6] = {iN[0], iN[1], iN[2], midN[0], midN[1], midN[2]};
553 backward_shift(nodalConnPerCell, 6);
554 dimMesh2->insertNextCell(INTERP_KERNEL::NORM_TRI6, 6, nodalConnPerCell);
555 MeshFormatElement e(MeshFormat::GmfTriangles, i-1);
556 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 2 -_dim);
560 mcIdType nodalConnPerCell[3] = {iN[0], iN[1], iN[2]};
561 backward_shift(nodalConnPerCell, 3);
562 dimMesh2->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,nodalConnPerCell);
563 MeshFormatElement e(MeshFormat::GmfTriangles, i-1);
564 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 2 -_dim);
566 if ( !midN.empty() ) MeshFormat::FreeVector( midN );
572 void MeshFormatReader::setQuadrangles( MEDCoupling::MEDCouplingUMesh* dimMesh2, int nbQuad)
574 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
576 // read extra vertices for quadratic quadrangles
577 std::vector< std::vector<int> > quadNodesAtQuadrilaterals( nbQuad + 1 );
578 if ( int nbQuadQuad = _reader.GmfStatKwd( _myCurrentFileId, MeshFormat::GmfExtraVerticesAtQuadrilaterals ))
580 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtQuadrilaterals);
581 for ( int i = 1; i <= nbQuadQuad; ++i )
583 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtQuadrilaterals,
584 &iN[0], &iN[1], &iN[2], &iN[3], &iN[4], &iN[5], &iN[6]);
585 if ( iN[0] <= nbQuad )
587 std::vector<int>& nodes = quadNodesAtQuadrilaterals[ iN[0] ];
588 nodes.insert( nodes.end(), & iN[2], & iN[6+1] );
589 nodes.resize( iN[1] );
593 // create quadrangles
594 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfQuadrilaterals);
595 for ( int i = 1; i <= nbQuad; ++i )
597 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfQuadrilaterals, &iN[0], &iN[1], &iN[2], &iN[3], &ref);
598 std::vector<int>& midN = quadNodesAtQuadrilaterals[ i ];
599 if ( midN.size() == 8-4 ) // QUAD8
601 mcIdType nodalConnPerCell[8] = {iN[0], iN[1], iN[2], iN[3],
602 midN[0], midN[1], midN[2], midN[3]
604 backward_shift(nodalConnPerCell, 8);
605 dimMesh2->insertNextCell(INTERP_KERNEL::NORM_QUAD8,8, nodalConnPerCell);
606 MeshFormatElement e(MeshFormat::GmfQuadrilaterals, i-1);
607 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 2 -_dim);
609 else if ( midN.size() > 8-4 ) // QUAD9
611 mcIdType nodalConnPerCell[9] = {iN[0], iN[1], iN[2], iN[3],
612 midN[0], midN[1], midN[2], midN[3], midN[4]
614 backward_shift(nodalConnPerCell, 9);
615 dimMesh2->insertNextCell(INTERP_KERNEL::NORM_QUAD9,9, nodalConnPerCell);
616 MeshFormatElement e(MeshFormat::GmfQuadrilaterals, i-1);
617 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 2 -_dim);
621 mcIdType nodalConnPerCell[4] = {iN[0], iN[1], iN[2], iN[3]};
622 backward_shift(nodalConnPerCell, 4);
623 dimMesh2->insertNextCell(INTERP_KERNEL::NORM_QUAD4, 4, nodalConnPerCell);
624 MeshFormatElement e(MeshFormat::GmfQuadrilaterals, i-1);
625 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 2 -_dim);
627 if ( !midN.empty() ) MeshFormat::FreeVector( midN );
632 void MeshFormatReader::setTetrahedras( MEDCoupling::MEDCouplingUMesh* dimMesh3, int nbTet)
634 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
636 // read extra vertices for quadratic tetrahedra
637 std::vector< std::vector<int> > quadNodesAtTetrahedra( nbTet + 1 );
638 if ( int nbQuadTetra = _reader.GmfStatKwd( _myCurrentFileId, MeshFormat::GmfExtraVerticesAtTetrahedra ))
640 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtTetrahedra);
641 for ( int i = 1; i <= nbQuadTetra; ++i )
643 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtTetrahedra,
644 &iN[0], &iN[1], &iN[2], &iN[3], &iN[4], &iN[5], &iN[6], &iN[7]);
645 if ( iN[0] <= nbTet )
647 std::vector<int>& nodes = quadNodesAtTetrahedra[ iN[0] ];
648 nodes.insert( nodes.end(), & iN[2], & iN[7+1] );
649 nodes.resize( iN[1] );
654 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfTetrahedra);
655 for ( int i = 1; i <= nbTet; ++i )
657 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfTetrahedra, &iN[0], &iN[1], &iN[2], &iN[3], &ref);
658 std::vector<int>& midN = quadNodesAtTetrahedra[ i ];
659 if ( midN.size() >= 10-4 ) // TETRA10
661 mcIdType nodalConnPerCell[10] = {iN[0], iN[2], iN[1], iN[3],
662 midN[2], midN[1], midN[0], midN[3], midN[5], midN[4]
664 backward_shift(nodalConnPerCell, 10);
665 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_TETRA10, 10,nodalConnPerCell);
666 MeshFormatElement e(MeshFormat::GmfTetrahedra, i-1);
667 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
671 mcIdType nodalConnPerCell[4] = {iN[0], iN[2], iN[1], iN[3]};
672 backward_shift(nodalConnPerCell, 4);
673 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_TETRA4,4, nodalConnPerCell);
674 MeshFormatElement e(MeshFormat::GmfTetrahedra, i-1);
675 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
677 if ( !midN.empty() ) MeshFormat::FreeVector( midN );
683 void MeshFormatReader::setPyramids( MEDCoupling::MEDCouplingUMesh* dimMesh3, int nbPyr)
685 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
687 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfPyramids);
688 for ( int i = 1; i <= nbPyr; ++i )
690 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfPyramids, &iN[0], &iN[1], &iN[2], &iN[3], &iN[4], &ref);
691 mcIdType nodalConnPerCell[5] = {iN[3], iN[2], iN[1], iN[0], iN[4]};
692 backward_shift(nodalConnPerCell, 5);
693 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_PYRA5, 5,nodalConnPerCell);
694 MeshFormatElement e(MeshFormat::GmfPyramids, i-1);
695 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
701 void MeshFormatReader::setHexahedras( MEDCoupling::MEDCouplingUMesh* dimMesh3, int nbHex)
703 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
705 // read extra vertices for quadratic hexahedra
706 std::vector< std::vector<int> > quadNodesAtHexahedra( nbHex + 1 );
707 if ( int nbQuadHexa = _reader.GmfStatKwd( _myCurrentFileId, MeshFormat::GmfExtraVerticesAtHexahedra ))
709 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtHexahedra);
710 for ( int i = 1; i <= nbQuadHexa; ++i )
712 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfExtraVerticesAtHexahedra, &iN[0], &iN[1], // Hexa Id, Nb extra vertices
713 &iN[2], &iN[3], &iN[4], &iN[5],
714 &iN[6], &iN[7], &iN[8], &iN[9],
715 &iN[10], &iN[11], &iN[12], &iN[13], // HEXA20
717 &iN[15], &iN[16], &iN[17], &iN[18],
720 if ( iN[0] <= nbHex )
722 std::vector<int>& nodes = quadNodesAtHexahedra[ iN[0] ];
723 nodes.insert( nodes.end(), & iN[2], & iN[20+1] );
724 nodes.resize( iN[1] );
731 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfHexahedra);
732 for ( int i = 1; i <= nbHex; ++i )
734 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfHexahedra, &iN[0], &iN[1], &iN[2], &iN[3],
735 &iN[4], &iN[5], &iN[6], &iN[7], &ref);
736 std::vector<int>& midN = quadNodesAtHexahedra[ i ];
737 if ( midN.size() == 20-8 ) // HEXA20
739 mcIdType nodalConnPerCell[20] = {iN[0], iN[3], iN[2], iN[1],
740 iN[4], iN[7], iN[6], iN[5],
741 midN[3], midN[2], midN[1], midN[0],
742 midN[7], midN[6], midN[5], midN[4],
743 midN[8], midN[11], midN[10], midN[9]
745 backward_shift(nodalConnPerCell, 20);
746 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_HEXA20,20, nodalConnPerCell);
747 MeshFormatElement e(MeshFormat::GmfHexahedra, i-1);
748 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
752 else if ( midN.size() >= 27-8 ) // HEXA27
754 mcIdType nodalConnPerCell[27] = {iN[0], iN[3], iN[2], iN[1],
755 iN[4], iN[7], iN[6], iN[5],
756 midN[3], midN[2], midN[1], midN[0],
757 midN[7], midN[6], midN[5], midN[4],
758 midN[8], midN[11], midN[10], midN[9],
760 midN[16], midN[15], midN[14], midN[13],
764 backward_shift(nodalConnPerCell, 27);
765 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_HEXA27,27, nodalConnPerCell);
766 MeshFormatElement e(MeshFormat::GmfHexahedra, i-1);
767 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
772 mcIdType nodalConnPerCell[8] = {iN[0], iN[3], iN[2], iN[1],
773 iN[4], iN[7], iN[6], iN[5]
775 backward_shift(nodalConnPerCell, 8);
776 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_HEXA8,8, nodalConnPerCell);
777 MeshFormatElement e(MeshFormat::GmfHexahedra, i-1);
778 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
781 if ( !midN.empty() ) MeshFormat::FreeVector( midN );
789 void MeshFormatReader::setPrisms(MEDCoupling::MEDCouplingUMesh* dimMesh3, int nbPrism)
791 int iN[28]; // 28 - nb nodes in HEX27 (+ 1 for safety :)
793 _reader.GmfGotoKwd(_myCurrentFileId, MeshFormat::GmfPrisms);
794 for ( int i = 1; i <= nbPrism; ++i )
796 _reader.GmfGetLin(_myCurrentFileId, MeshFormat::GmfPrisms, &iN[0], &iN[1], &iN[2], &iN[3], &iN[4], &iN[5], &ref);
797 mcIdType nodalConnPerCell[8] = {iN[0], iN[2], iN[1], iN[3], iN[5], iN[4]};
798 backward_shift(nodalConnPerCell, 8);
799 dimMesh3->insertNextCell(INTERP_KERNEL::NORM_PENTA6, 6,nodalConnPerCell);
800 MeshFormatElement e(MeshFormat::GmfPrisms, i-1);
801 _fams.insert(std::pair <int, MeshFormatElement> (ref, e), 3 -_dim);
808 void MeshFormatReader::callParserGetLin( MeshFormat::GmfKwdCod kwd, double* val, int valSize, int* ref)
814 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], ref);
819 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], &val[1], ref);
824 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], &val[1], &val[2], ref);
829 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], &val[1], &val[2], &val[3], ref);
834 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], &val[1], &val[2], &val[3], &val[4], &val[5], ref);
839 _reader.GmfGetLin(_myCurrentFileId, kwd, &val[0], &val[1], &val[2], &val[3], &val[4], &val[5], &val[6], &val[7], &val[8], ref);
845 void MeshFormatReader::setTypeOfFieldAndDimRel(MeshFormat::GmfKwdCod kwd, MEDCoupling::TypeOfField* typeOfField, int* dimRel )
849 case MeshFormat::GmfSolAtVertices :
851 *typeOfField = MEDCoupling::ON_NODES;
855 case MeshFormat::GmfSolAtEdges :
857 *typeOfField = MEDCoupling::ON_CELLS;
861 case MeshFormat::GmfSolAtTriangles :
862 case MeshFormat::GmfSolAtQuadrilaterals :
864 *typeOfField = MEDCoupling::ON_CELLS;
868 case MeshFormat::GmfSolAtTetrahedra :
869 case MeshFormat::GmfSolAtPrisms :
870 case MeshFormat::GmfSolAtHexahedra :
872 *typeOfField = MEDCoupling::ON_CELLS;
880 INTERP_KERNEL::NormalizedCellType MeshFormatReader::toMedType(MeshFormat::GmfKwdCod kwd)
882 INTERP_KERNEL::NormalizedCellType type;
885 //~case MeshFormat::GmfEdges :
888 //~type = INTERP_KERNEL::NORM_SEG2;
891 //~case MeshFormat::GmfTriangles :
893 //~type = INTERP_KERNEL::NORM_TRI3;
896 //~case MeshFormat::GmfQuadrilaterals :
898 //~type = INTERP_KERNEL::NORM_QUAD;
901 //~case MeshFormat::GmfSolAtTetrahedra :
902 //~case MeshFormat::GmfSolAtPrisms :
903 //~case MeshFormat::GmfSolAtHexahedra :
905 //~*typeOfField = MEDCoupling::ON_CELLS;
906 //~*dimRel = 3 - _dim;
914 void MeshFormatReader::buildFamilies()
916 buildNodesFamilies();
917 buildCellsFamilies();
920 void MeshFormatReader::buildCellsFamilies()
922 std::vector<int> levs = _uMesh->getNonEmptyLevels();
923 for (size_t iDim = 0; iDim<levs.size(); iDim++ )
925 int dimRelMax = levs[iDim];
926 std::map <int, std::vector<MeshFormatElement>* > famDim = _fams.getMapAtLevel(dimRelMax);
927 std::map <int, std::vector<MeshFormatElement>* >::const_iterator _meshFormatFamsIt = famDim.begin();
928 std::vector< const MEDCoupling::DataArrayIdType* > fams;
929 MEDCoupling::DataArrayIdType* cellIds = MEDCoupling::DataArrayIdType::New();
930 cellIds->alloc(_uMesh->getSizeAtLevel(dimRelMax), 1);
931 cellIds->fillWithZero();
933 for(; _meshFormatFamsIt!= famDim.end(); ++_meshFormatFamsIt)
935 const int famId = _meshFormatFamsIt->first;
936 std::string famName ="FromMeshGemsFormatAttributFamily_"+std::to_string(famId);
937 std::vector <MeshFormatElement>* cellsInFam = _meshFormatFamsIt->second;
938 if (!famId) continue;
939 std::vector <MeshFormatElement>::iterator cellsInFamIt = cellsInFam->begin();
941 _uMesh->addFamily(famName, famId);
942 for ( ; cellsInFamIt !=cellsInFam->end(); ++cellsInFamIt)
944 cellIds->setIJ(cellsInFamIt->_id, 0, famId);
948 _uMesh->setFamilyFieldArr(dimRelMax, cellIds->deepCopy());
955 void MeshFormatReader::buildNodesFamilies()
957 std::vector<int> levs = _uMesh->getNonEmptyLevels();
959 std::map <int, std::vector<MeshFormatElement>* > famDim = _fams.getMapAtLevel(dimRelMax);
960 std::map <int, std::vector<MeshFormatElement>* >::const_iterator _meshFormatFamsIt = famDim.begin();
961 std::vector< const MEDCoupling::DataArrayIdType* > fams;
962 MEDCoupling::DataArrayIdType* cellIds = MEDCoupling::DataArrayIdType::New();
963 cellIds->alloc(_uMesh->getSizeAtLevel(dimRelMax), 1);
964 cellIds->fillWithZero();
966 for(; _meshFormatFamsIt!= famDim.end(); ++_meshFormatFamsIt)
968 const int famId = _meshFormatFamsIt->first;
969 if (!famId) continue;
970 bool thisIsACellFamily = false;
972 for (size_t iDim = 0; iDim<levs.size(); iDim++ )
974 int dimMesh = levs[iDim];
975 std::map <int, std::vector<MeshFormatElement>* > famDimAtLevel = _fams.getMapAtLevel(dimMesh);
976 std::map <int, std::vector<MeshFormatElement>* >::iterator famDimAtLevelId = famDimAtLevel.find(famId);
977 if (famDimAtLevelId != famDimAtLevel.end())
979 thisIsACellFamily = true;
985 if (thisIsACellFamily) continue;
986 std::string famName ="FromMeshGemsFormatAttributFamily_"+std::to_string(famId);
987 std::vector <MeshFormatElement>* cellsInFam = _meshFormatFamsIt->second;
988 std::vector <MeshFormatElement>::iterator cellsInFamIt = cellsInFam->begin();
990 _uMesh->addFamily(famName, famId);
991 for ( ; cellsInFamIt !=cellsInFam->end(); ++cellsInFamIt)
993 cellIds->setIJ(cellsInFamIt->_id, 0, famId);
997 _uMesh->setFamilyFieldArr(dimRelMax, cellIds->deepCopy());