src/DriverMED/DriverMED_W_SMESHDS_Mesh.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 //  SMESH DriverMED : driver to read and write 'med' files
  24 //  File   : DriverMED_W_SMESHDS_Mesh.cxx
  25 //  Module : SMESH
  26 //
  27 #include <sstream>
  28
  29 #include "DriverMED_W_SMESHDS_Mesh.h"
  30 #include "DriverMED_Family.h"
  31
  32 #include "SMESHDS_Mesh.hxx"
  33 #include "SMDS_MeshElement.hxx"
  34 #include "SMDS_MeshNode.hxx"
  35 #include "SMDS_PolyhedralVolumeOfNodes.hxx"
  36
  37 #include "utilities.h"
  38
  39 #include "MED_Utilities.hxx"
  40
  41 #define _EDF_NODE_IDS_
  42 //#define _ELEMENTS_BY_DIM_
  43
  44 using namespace std;
  45 using namespace MED;
  46
  47
  48 DriverMED_W_SMESHDS_Mesh::DriverMED_W_SMESHDS_Mesh():
  49   myAllSubMeshes (false),
  50   myDoGroupOfNodes (false),
  51   myDoGroupOfEdges (false),
  52   myDoGroupOfFaces (false),
  53   myDoGroupOfVolumes (false),
  54   myDoGroupOf0DElems(false),
  55   myDoGroupOfBalls(false)
  56 {}
  57
  58 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName,
  59                                        MED::EVersion theId)
  60 {
  61   myMed = CrWrapper(theFileName,theId);
  62   Driver_SMESHDS_Mesh::SetFile(theFileName);
  63 }
  64
  65 void DriverMED_W_SMESHDS_Mesh::SetFile(const std::string& theFileName)
  66 {
  67   return SetFile(theFileName,MED::eV2_2);
  68 }
  69
  70 string DriverMED_W_SMESHDS_Mesh::GetVersionString(const MED::EVersion theVersion, int theNbDigits)
  71 {
  72   TInt majeur, mineur, release;
  73   majeur =  mineur = release = 0;
  74 //   if ( theVersion == eV2_1 )
  75 //     MED::GetVersionRelease<eV2_1>(majeur, mineur, release);
  76 //   else
  77     MED::GetVersionRelease<eV2_2>(majeur, mineur, release);
  78   ostringstream name;
  79   if ( theNbDigits > 0 )
  80     name << majeur;
  81   if ( theNbDigits > 1 )
  82     name << "." << mineur;
  83   if ( theNbDigits > 2 )
  84     name << "." << release;
  85   return name.str();
  86 }
  87
  88 void DriverMED_W_SMESHDS_Mesh::AddGroup(SMESHDS_GroupBase* theGroup)
  89 {
  90   myGroups.push_back(theGroup);
  91 }
  92
  93 void DriverMED_W_SMESHDS_Mesh::AddAllSubMeshes()
  94 {
  95   myAllSubMeshes = true;
  96 }
  97
  98 void DriverMED_W_SMESHDS_Mesh::AddSubMesh(SMESHDS_SubMesh* theSubMesh, int theID)
  99 {
 100   mySubMeshes[theID] = theSubMesh;
 101 }
 102
 103 void DriverMED_W_SMESHDS_Mesh::AddGroupOfNodes()
 104 {
 105   myDoGroupOfNodes = true;
 106 }
 107
 108 void DriverMED_W_SMESHDS_Mesh::AddGroupOfEdges()
 109 {
 110   myDoGroupOfEdges = true;
 111 }
 112
 113 void DriverMED_W_SMESHDS_Mesh::AddGroupOfFaces()
 114 {
 115   myDoGroupOfFaces = true;
 116 }
 117
 118 void DriverMED_W_SMESHDS_Mesh::AddGroupOfVolumes()
 119 {
 120   myDoGroupOfVolumes = true;
 121 }
 122
 123 namespace{
 124   typedef double (SMDS_MeshNode::* TGetCoord)() const;
 125   typedef const char* TName;
 126   typedef const char* TUnit;
 127
 128   // name length in a mesh must be equal to 16 :
 129   //         1234567890123456
 130   TName M = "m               ";
 131   TName X = "x               ";
 132   TName Y = "y               ";
 133   TName Z = "z               ";
 134
 135   TUnit aUnit[3] = {M,M,M};
 136
 137   // 3 dim
 138   TGetCoord aXYZGetCoord[3] = {
 139     &SMDS_MeshNode::X,
 140     &SMDS_MeshNode::Y,
 141     &SMDS_MeshNode::Z
 142   };
 143   TName aXYZName[3] = {X,Y,Z};
 144
 145   // 2 dim
 146   TGetCoord aXYGetCoord[2] = {
 147     &SMDS_MeshNode::X,
 148     &SMDS_MeshNode::Y
 149   };
 150   TName aXYName[2] = {X,Y};
 151
 152   TGetCoord aYZGetCoord[2] = {
 153     &SMDS_MeshNode::Y,
 154     &SMDS_MeshNode::Z
 155   };
 156   TName aYZName[2] = {Y,Z};
 157
 158   TGetCoord aXZGetCoord[2] = {
 159     &SMDS_MeshNode::X,
 160     &SMDS_MeshNode::Z
 161   };
 162   TName aXZName[2] = {X,Z};
 163
 164   // 1 dim
 165   TGetCoord aXGetCoord[1] = {
 166     &SMDS_MeshNode::X
 167   };
 168   TName aXName[1] = {X};
 169
 170   TGetCoord aYGetCoord[1] = {
 171     &SMDS_MeshNode::Y
 172   };
 173   TName aYName[1] = {Y};
 174
 175   TGetCoord aZGetCoord[1] = {
 176     &SMDS_MeshNode::Z
 177   };
 178   TName aZName[1] = {Z};
 179
 180
 181   class TCoordHelper{
 182     SMDS_NodeIteratorPtr myNodeIter;
 183     const SMDS_MeshNode* myCurrentNode;
 184     TGetCoord* myGetCoord;
 185     TName* myName;
 186     TUnit* myUnit;
 187   public:
 188     TCoordHelper(const SMDS_NodeIteratorPtr& theNodeIter,
 189                  TGetCoord* theGetCoord,
 190                  TName* theName,
 191                  TUnit* theUnit = aUnit):
 192       myNodeIter(theNodeIter),
 193       myGetCoord(theGetCoord),
 194       myName(theName),
 195       myUnit(theUnit)
 196     {}
 197     virtual ~TCoordHelper(){}
 198     bool Next(){
 199       return myNodeIter->more() &&
 200         (myCurrentNode = myNodeIter->next());
 201     }
 202     const SMDS_MeshNode* GetNode(){
 203       return myCurrentNode;
 204     }
 205     MED::TIntVector::value_type GetID(){
 206       return myCurrentNode->GetID();
 207     }
 208     MED::TFloatVector::value_type GetCoord(TInt theCoodId){
 209       return (myCurrentNode->*myGetCoord[theCoodId])();
 210     }
 211     MED::TStringVector::value_type GetName(TInt theDimId){
 212       return myName[theDimId];
 213     }
 214     MED::TStringVector::value_type GetUnit(TInt theDimId){
 215       return myUnit[theDimId];
 216     }
 217   };
 218   typedef boost::shared_ptr<TCoordHelper> TCoordHelperPtr;
 219
 220   //-------------------------------------------------------
 221   /*!
 222    * \brief Structure describing element type
 223    */
 224   //-------------------------------------------------------
 225   struct TElemTypeData
 226   {
 227     EEntiteMaillage     _entity;
 228     EGeometrieElement   _geomType;
 229     TInt                _nbElems;
 230     SMDSAbs_ElementType _smdsType;
 231
 232     TElemTypeData (EEntiteMaillage entity, EGeometrieElement geom, TInt nb, SMDSAbs_ElementType type)
 233       : _entity(entity), _geomType(geom), _nbElems( nb ), _smdsType( type ) {}
 234   };
 235
 236
 237   typedef NCollection_DataMap< Standard_Address, int > TElemFamilyMap;
 238
 239   //================================================================================
 240   /*!
 241    * \brief Fills element to famaly ID map for element type.
 242    * Removes all families of anElemType
 243    */
 244   //================================================================================
 245
 246   void fillElemFamilyMap( TElemFamilyMap &            anElemFamMap,
 247                           list<DriverMED_FamilyPtr> & aFamilies,
 248                           const SMDSAbs_ElementType   anElemType)
 249   {
 250     anElemFamMap.Clear();
 251     //anElemFamMap.clear();
 252     list<DriverMED_FamilyPtr>::iterator aFamsIter = aFamilies.begin();
 253     while ( aFamsIter != aFamilies.end() )
 254     {
 255       if ((*aFamsIter)->GetType() != anElemType) {
 256         aFamsIter++;
 257       }
 258       else {
 259         int aFamId = (*aFamsIter)->GetId();
 260         const set<const SMDS_MeshElement *>& anElems = (*aFamsIter)->GetElements();
 261         set<const SMDS_MeshElement *>::const_iterator anElemsIter = anElems.begin();
 262         for (; anElemsIter != anElems.end(); anElemsIter++)
 263         {
 264           anElemFamMap.Bind( (Standard_Address)*anElemsIter, aFamId );
 265           //anElemFamMap[*anElemsIter] = aFamId;
 266         }
 267         // remove a family from the list
 268         aFamilies.erase( aFamsIter++ );
 269       }
 270     }
 271   }
 272
 273   //================================================================================
 274   /*!
 275    * \brief For an element, return family ID found in the map or a default one
 276    */
 277   //================================================================================
 278
 279   int getFamilyId( const TElemFamilyMap &  anElemFamMap,
 280                    const SMDS_MeshElement* anElement,
 281                    const int               aDefaultFamilyId)
 282   {
 283     if ( anElemFamMap.IsBound( (Standard_Address) anElement ))
 284       return anElemFamMap( (Standard_Address) anElement );
 285 //     TElemFamilyMap::iterator elem_famNum = anElemFamMap.find( anElement );
 286 //     if ( elem_famNum != anElemFamMap.end() )
 287 //       return elem_famNum->second;
 288     return aDefaultFamilyId;
 289   }
 290 }
 291
 292 Driver_Mesh::Status DriverMED_W_SMESHDS_Mesh::Perform()
 293 {
 294   Status aResult = DRS_OK;
 295   if (myMesh->hasConstructionEdges() || myMesh->hasConstructionFaces()) {
 296     INFOS("SMDS_MESH with hasConstructionEdges() or hasConstructionFaces() do not supports!!!");
 297     return DRS_FAIL;
 298   }
 299   try {
 300     MESSAGE("Perform - myFile : "<<myFile);
 301
 302     // Creating the MED mesh for corresponding SMDS structure
 303     //-------------------------------------------------------
 304     string aMeshName;
 305     if (myMeshId != -1) {
 306       ostringstream aMeshNameStr;
 307       aMeshNameStr<<myMeshId;
 308       aMeshName = aMeshNameStr.str();
 309     } else {
 310       aMeshName = myMeshName;
 311     }
 312
 313     // Mesh dimension definition
 314     TInt aSpaceDimension;
 315     TCoordHelperPtr aCoordHelperPtr;
 316     {
 317       bool anIsXDimension = false;
 318       bool anIsYDimension = false;
 319       bool anIsZDimension = false;
 320       {
 321         SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator();
 322         double aBounds[6];
 323         if(aNodesIter->more()){
 324           const SMDS_MeshNode* aNode = aNodesIter->next();
 325           aBounds[0] = aBounds[1] = aNode->X();
 326           aBounds[2] = aBounds[3] = aNode->Y();
 327           aBounds[4] = aBounds[5] = aNode->Z();
 328         }
 329         while(aNodesIter->more()){
 330           const SMDS_MeshNode* aNode = aNodesIter->next();
 331           aBounds[0] = min(aBounds[0],aNode->X());
 332           aBounds[1] = max(aBounds[1],aNode->X());
 333
 334           aBounds[2] = min(aBounds[2],aNode->Y());
 335           aBounds[3] = max(aBounds[3],aNode->Y());
 336
 337           aBounds[4] = min(aBounds[4],aNode->Z());
 338           aBounds[5] = max(aBounds[5],aNode->Z());
 339         }
 340
 341         double EPS = 1.0E-7;
 342         anIsXDimension = (aBounds[1] - aBounds[0]) + abs(aBounds[1]) + abs(aBounds[0]) > EPS;
 343         anIsYDimension = (aBounds[3] - aBounds[2]) + abs(aBounds[3]) + abs(aBounds[2]) > EPS;
 344         anIsZDimension = (aBounds[5] - aBounds[4]) + abs(aBounds[5]) + abs(aBounds[4]) > EPS;
 345         aSpaceDimension = anIsXDimension + anIsYDimension + anIsZDimension;
 346         if(!aSpaceDimension)
 347           aSpaceDimension = 3;
 348         // PAL16857(SMESH not conform to the MED convention):
 349         if ( aSpaceDimension == 2 && anIsZDimension ) // 2D only if mesh is in XOY plane
 350           aSpaceDimension = 3;
 351         // PAL18941(a saved study with a mesh belong Z is opened and the mesh is belong X)
 352         if ( aSpaceDimension == 1 && !anIsXDimension ) {// 1D only if mesh is along OX
 353           if ( anIsYDimension ) {
 354             aSpaceDimension = 2;
 355             anIsXDimension = true;
 356           } else {
 357             aSpaceDimension = 3;
 358           }
 359         }
 360       }
 361
 362       SMDS_NodeIteratorPtr aNodesIter = myMesh->nodesIterator(/*idInceasingOrder=*/true);
 363       switch(aSpaceDimension){
 364       case 3:
 365         aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYZGetCoord,aXYZName));
 366         break;
 367       case 2:
 368         if(anIsXDimension && anIsYDimension)
 369           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXYGetCoord,aXYName));
 370         if(anIsYDimension && anIsZDimension)
 371           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYZGetCoord,aYZName));
 372         if(anIsXDimension && anIsZDimension)
 373           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXZGetCoord,aXZName));
 374         break;
 375       case 1:
 376         if(anIsXDimension)
 377           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aXGetCoord,aXName));
 378         if(anIsYDimension)
 379           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aYGetCoord,aYName));
 380         if(anIsZDimension)
 381           aCoordHelperPtr.reset(new TCoordHelper(aNodesIter,aZGetCoord,aZName));
 382         break;
 383       }
 384     }
 385     TInt aMeshDimension = 0;
 386     if ( myMesh->NbEdges() > 0 )
 387       aMeshDimension = 1;
 388     if ( myMesh->NbFaces() > 0 )
 389       aMeshDimension = 2;
 390     if ( myMesh->NbVolumes() > 0 )
 391       aMeshDimension = 3;
 392
 393     PMeshInfo aMeshInfo = myMed->CrMeshInfo(aMeshDimension,aSpaceDimension,aMeshName);
 394     MESSAGE("Add - aMeshName : "<<aMeshName<<"; "<<aMeshInfo->GetName());
 395     myMed->SetMeshInfo(aMeshInfo);
 396
 397     // Storing SMDS groups and sub-meshes as med families
 398     //----------------------------------------------------
 399     int myNodesDefaultFamilyId      = 0;
 400     int my0DElementsDefaultFamilyId = 0;
 401     int myBallsDefaultFamilyId      = 0;
 402     int myEdgesDefaultFamilyId      = 0;
 403     int myFacesDefaultFamilyId      = 0;
 404     int myVolumesDefaultFamilyId    = 0;
 405     int nbNodes      = myMesh->NbNodes();
 406     int nb0DElements = myMesh->Nb0DElements();
 407     int nbBalls      = myMesh->NbBalls();
 408     int nbEdges      = myMesh->NbEdges();
 409     int nbFaces      = myMesh->NbFaces();
 410     int nbVolumes    = myMesh->NbVolumes();
 411     if (myDoGroupOfNodes && nbNodes) myNodesDefaultFamilyId = REST_NODES_FAMILY;
 412     if (myDoGroupOfEdges && nbEdges) myEdgesDefaultFamilyId = REST_EDGES_FAMILY;
 413     if (myDoGroupOfFaces && nbFaces) myFacesDefaultFamilyId = REST_FACES_FAMILY;
 414     if (myDoGroupOfVolumes && nbVolumes) myVolumesDefaultFamilyId = REST_VOLUMES_FAMILY;
 415     if (myDoGroupOf0DElems && nb0DElements) my0DElementsDefaultFamilyId = REST_0DELEM_FAMILY;
 416     if (myDoGroupOfBalls && nbBalls) myBallsDefaultFamilyId = REST_BALL_FAMILY;
 417
 418     MESSAGE("Perform - aFamilyInfo");
 419     //cout << " DriverMED_Family::MakeFamilies() " << endl;
 420     list<DriverMED_FamilyPtr> aFamilies;
 421     if (myAllSubMeshes) {
 422       aFamilies = DriverMED_Family::MakeFamilies
 423         (myMesh->SubMeshes(), myGroups,
 424          myDoGroupOfNodes   && nbNodes,
 425          myDoGroupOfEdges   && nbEdges,
 426          myDoGroupOfFaces   && nbFaces,
 427          myDoGroupOfVolumes && nbVolumes,
 428          myDoGroupOf0DElems && nb0DElements,
 429          myDoGroupOfBalls   && nbBalls);
 430     } else {
 431       aFamilies = DriverMED_Family::MakeFamilies
 432         (mySubMeshes, myGroups,
 433          myDoGroupOfNodes   && nbNodes,
 434          myDoGroupOfEdges   && nbEdges,
 435          myDoGroupOfFaces   && nbFaces,
 436          myDoGroupOfVolumes && nbVolumes,
 437          myDoGroupOf0DElems && nb0DElements,
 438          myDoGroupOfBalls   && nbBalls);
 439     }
 440     //cout << " myMed->SetFamilyInfo() " << endl;
 441     list<DriverMED_FamilyPtr>::iterator aFamsIter;
 442     for (aFamsIter = aFamilies.begin(); aFamsIter != aFamilies.end(); aFamsIter++)
 443     {
 444       PFamilyInfo aFamilyInfo = (*aFamsIter)->GetFamilyInfo(myMed,aMeshInfo);
 445       myMed->SetFamilyInfo(aFamilyInfo);
 446     }
 447
 448     // Storing SMDS nodes to the MED file for the MED mesh
 449     //----------------------------------------------------
 450 #ifdef _EDF_NODE_IDS_
 451     typedef map<TInt,TInt> TNodeIdMap;
 452     TNodeIdMap aNodeIdMap;
 453 #endif
 454     const EModeSwitch   theMode        = eFULL_INTERLACE;
 455     const ERepere       theSystem      = eCART;
 456     const EBooleen      theIsElemNum   = eVRAI;
 457     const EBooleen      theIsElemNames = eFAUX;
 458     const EConnectivite theConnMode    = eNOD;
 459
 460     TInt aNbNodes = myMesh->NbNodes();
 461     //cout << " myMed->CrNodeInfo() aNbNodes = " << aNbNodes << endl;
 462     PNodeInfo aNodeInfo = myMed->CrNodeInfo(aMeshInfo, aNbNodes,
 463                                             theMode, theSystem, theIsElemNum, theIsElemNames);
 464
 465     //cout << " fillElemFamilyMap( SMDSAbs_Node )" << endl;
 466     // find family numbers for nodes
 467     TElemFamilyMap anElemFamMap;
 468     fillElemFamilyMap( anElemFamMap, aFamilies, SMDSAbs_Node );
 469
 470     for (TInt iNode = 0; aCoordHelperPtr->Next(); iNode++)
 471     {
 472       // coordinates
 473       TCoordSlice aTCoordSlice = aNodeInfo->GetCoordSlice( iNode );
 474       for(TInt iCoord = 0; iCoord < aSpaceDimension; iCoord++){
 475         aTCoordSlice[iCoord] = aCoordHelperPtr->GetCoord(iCoord);
 476       }
 477       // node number
 478       int aNodeID = aCoordHelperPtr->GetID();
 479       aNodeInfo->SetElemNum( iNode, aNodeID );
 480 #ifdef _EDF_NODE_IDS_
 481       aNodeIdMap[aNodeID] = iNode+1;
 482 #endif
 483       // family number
 484       const SMDS_MeshNode* aNode = aCoordHelperPtr->GetNode();
 485       int famNum = getFamilyId( anElemFamMap, aNode, myNodesDefaultFamilyId );
 486       aNodeInfo->SetFamNum( iNode, famNum );
 487     }
 488     anElemFamMap.Clear();
 489
 490     // coordinate names and units
 491     for (TInt iCoord = 0; iCoord < aSpaceDimension; iCoord++) {
 492       aNodeInfo->SetCoordName( iCoord, aCoordHelperPtr->GetName(iCoord));
 493       aNodeInfo->SetCoordUnit( iCoord, aCoordHelperPtr->GetUnit(iCoord));
 494     }
 495
 496     //cout << " SetNodeInfo(aNodeInfo)" << endl;
 497     MESSAGE("Perform - aNodeInfo->GetNbElem() = "<<aNbNodes);
 498     myMed->SetNodeInfo(aNodeInfo);
 499     aNodeInfo.reset(); // free memory used for arrays
 500
 501
 502     // Storing SMDS elements to the MED file for the MED mesh
 503     //-------------------------------------------------------
 504     // Write one element type at once in order to minimize memory usage (PAL19276)
 505
 506     const SMDS_MeshInfo& nbElemInfo = myMesh->GetMeshInfo();
 507
 508     // poly elements are not supported by med-2.1
 509     bool polyTypesSupported = myMed->CrPolygoneInfo(aMeshInfo,eMAILLE,ePOLYGONE,0,0);
 510     TInt nbPolygonNodes = 0, nbPolyhedronNodes = 0, nbPolyhedronFaces = 0;
 511
 512     // collect info on all geom types
 513
 514     list< TElemTypeData > aTElemTypeDatas;
 515
 516     EEntiteMaillage anEntity = eMAILLE;
 517 #ifdef _ELEMENTS_BY_DIM_
 518     anEntity = eNOEUD_ELEMENT;
 519 #endif
 520     aTElemTypeDatas.push_back(TElemTypeData(anEntity,
 521                                             ePOINT1,
 522                                             nbElemInfo.Nb0DElements(),
 523                                             SMDSAbs_0DElement));
 524 #ifdef _ELEMENTS_BY_DIM_
 525     anEntity = eSTRUCT_ELEMENT;
 526 #endif
 527     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 528                                              eBALL,
 529                                              nbElemInfo.NbBalls(),
 530                                              SMDSAbs_Ball));
 531 #ifdef _ELEMENTS_BY_DIM_
 532     anEntity = eARETE;
 533 #endif
 534     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 535                                              eSEG2,
 536                                              nbElemInfo.NbEdges( ORDER_LINEAR ),
 537                                              SMDSAbs_Edge));
 538     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 539                                              eSEG3,
 540                                              nbElemInfo.NbEdges( ORDER_QUADRATIC ),
 541                                              SMDSAbs_Edge));
 542 #ifdef _ELEMENTS_BY_DIM_
 543     anEntity = eFACE;
 544 #endif
 545     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 546                                              eTRIA3,
 547                                              nbElemInfo.NbTriangles( ORDER_LINEAR ),
 548                                              SMDSAbs_Face));
 549     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 550                                              eTRIA6,
 551                                              nbElemInfo.NbTriangles( ORDER_QUADRATIC ) -
 552                                              nbElemInfo.NbBiQuadTriangles(),
 553                                              SMDSAbs_Face));
 554     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 555                                              eTRIA7,
 556                                              nbElemInfo.NbBiQuadTriangles(),
 557                                              SMDSAbs_Face));
 558     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 559                                              eQUAD4,
 560                                              nbElemInfo.NbQuadrangles( ORDER_LINEAR ),
 561                                              SMDSAbs_Face));
 562     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 563                                              eQUAD8,
 564                                              nbElemInfo.NbQuadrangles( ORDER_QUADRATIC ) -
 565                                              nbElemInfo.NbBiQuadQuadrangles(),
 566                                              SMDSAbs_Face));
 567     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 568                                              eQUAD9,
 569                                              nbElemInfo.NbBiQuadQuadrangles(),
 570                                              SMDSAbs_Face));
 571     if ( polyTypesSupported ) {
 572       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 573                                                ePOLYGONE,
 574                                                nbElemInfo.NbPolygons(),
 575                                                SMDSAbs_Face));
 576       // we need one more loop on poly elements to count nb of their nodes
 577       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 578                                                ePOLYGONE,
 579                                                nbElemInfo.NbPolygons(),
 580                                                SMDSAbs_Face));
 581     }
 582 #ifdef _ELEMENTS_BY_DIM_
 583     anEntity = eMAILLE;
 584 #endif
 585     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 586                                              eTETRA4,
 587                                              nbElemInfo.NbTetras( ORDER_LINEAR ),
 588                                              SMDSAbs_Volume));
 589     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 590                                              eTETRA10,
 591                                              nbElemInfo.NbTetras( ORDER_QUADRATIC ),
 592                                              SMDSAbs_Volume));
 593     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 594                                              ePYRA5,
 595                                              nbElemInfo.NbPyramids( ORDER_LINEAR ),
 596                                              SMDSAbs_Volume));
 597     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 598                                              ePYRA13,
 599                                              nbElemInfo.NbPyramids( ORDER_QUADRATIC ),
 600                                              SMDSAbs_Volume));
 601     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 602                                              ePENTA6,
 603                                              nbElemInfo.NbPrisms( ORDER_LINEAR ),
 604                                              SMDSAbs_Volume));
 605     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 606                                              ePENTA15,
 607                                              nbElemInfo.NbPrisms( ORDER_QUADRATIC ),
 608                                              SMDSAbs_Volume));
 609     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 610                                              eHEXA8,
 611                                              nbElemInfo.NbHexas( ORDER_LINEAR ),
 612                                              SMDSAbs_Volume));
 613     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 614                                              eHEXA20,
 615                                              nbElemInfo.NbHexas( ORDER_QUADRATIC )-
 616                                              nbElemInfo.NbTriQuadHexas(),
 617                                              SMDSAbs_Volume));
 618     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 619                                              eHEXA27,
 620                                              nbElemInfo.NbTriQuadHexas(),
 621                                              SMDSAbs_Volume));
 622     aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 623                                              eOCTA12,
 624                                              nbElemInfo.NbHexPrisms(),
 625                                              SMDSAbs_Volume));
 626     if ( polyTypesSupported ) {
 627       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 628                                                ePOLYEDRE,
 629                                                nbElemInfo.NbPolyhedrons(),
 630                                                SMDSAbs_Volume));
 631       // we need one more loop on poly elements to count nb of their nodes
 632       aTElemTypeDatas.push_back( TElemTypeData(anEntity,
 633                                                ePOLYEDRE,
 634                                                nbElemInfo.NbPolyhedrons(),
 635                                                SMDSAbs_Volume));
 636     }
 637
 638     vector< bool > isElemFamMapBuilt( SMDSAbs_NbElementTypes, false );
 639
 640     // loop on all geom types of elements
 641
 642     list< TElemTypeData >::iterator aElemTypeData = aTElemTypeDatas.begin();
 643     for ( ; aElemTypeData != aTElemTypeDatas.end(); ++aElemTypeData )
 644     {
 645       if ( aElemTypeData->_nbElems == 0 )
 646         continue;
 647
 648       int defaultFamilyId = 0;
 649       switch ( aElemTypeData->_smdsType ) {
 650       case SMDSAbs_0DElement:
 651         defaultFamilyId = my0DElementsDefaultFamilyId;
 652         break;
 653       case SMDSAbs_Ball:
 654         defaultFamilyId = myBallsDefaultFamilyId;
 655         break;
 656       case SMDSAbs_Edge:
 657         defaultFamilyId = myEdgesDefaultFamilyId;
 658         break;
 659       case SMDSAbs_Face:
 660         defaultFamilyId = myFacesDefaultFamilyId;
 661         break;
 662       case SMDSAbs_Volume:
 663         defaultFamilyId = myVolumesDefaultFamilyId;
 664         break;
 665       default:
 666         continue;
 667       }
 668
 669       // iterator on elements of a current type
 670       SMDS_ElemIteratorPtr elemIterator;
 671       int iElem = 0;
 672
 673       // Treat POLYGONs
 674       // ---------------
 675       if ( aElemTypeData->_geomType == ePOLYGONE )
 676       {
 677         elemIterator = myMesh->elementGeomIterator( SMDSGeom_POLYGON );
 678         if ( nbPolygonNodes == 0 ) {
 679           // Count nb of nodes
 680           while ( elemIterator->more() ) {
 681             const SMDS_MeshElement* anElem = elemIterator->next();
 682             nbPolygonNodes += anElem->NbNodes();
 683             if ( ++iElem == aElemTypeData->_nbElems )
 684               break;
 685           }
 686         }
 687         else {
 688           // Store in med file
 689           PPolygoneInfo aPolygoneInfo = myMed->CrPolygoneInfo(aMeshInfo,
 690                                                               aElemTypeData->_entity,
 691                                                               aElemTypeData->_geomType,
 692                                                               aElemTypeData->_nbElems,
 693                                                               nbPolygonNodes,
 694                                                               theConnMode, theIsElemNum,
 695                                                               theIsElemNames);
 696           TElemNum & index = *(aPolygoneInfo->myIndex.get());
 697           index[0] = 1;
 698
 699           while ( elemIterator->more() )
 700           {
 701             const SMDS_MeshElement* anElem = elemIterator->next();
 702             // index
 703             TInt aNbNodes = anElem->NbNodes();
 704             index[ iElem+1 ] = index[ iElem ] + aNbNodes;
 705
 706             // connectivity
 707             TConnSlice aTConnSlice = aPolygoneInfo->GetConnSlice( iElem );
 708             for(TInt iNode = 0; iNode < aNbNodes; iNode++) {
 709               const SMDS_MeshElement* aNode = anElem->GetNode( iNode );
 710 #ifdef _EDF_NODE_IDS_
 711               aTConnSlice[ iNode ] = aNodeIdMap[aNode->GetID()];
 712 #else
 713               aTConnSlice[ iNode ] = aNode->GetID();
 714 #endif
 715             }
 716             // element number
 717             aPolygoneInfo->SetElemNum( iElem, anElem->GetID() );
 718
 719             // family number
 720             int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 721             aPolygoneInfo->SetFamNum( iElem, famNum );
 722
 723             if ( ++iElem == aPolygoneInfo->GetNbElem() )
 724               break;
 725           }
 726           myMed->SetPolygoneInfo(aPolygoneInfo);
 727         }
 728
 729       }
 730
 731       // Treat POLYEDREs
 732       // ----------------
 733       else if (aElemTypeData->_geomType == ePOLYEDRE )
 734       {
 735         elemIterator = myMesh->elementGeomIterator( SMDSGeom_POLYHEDRA );
 736
 737         if ( nbPolyhedronNodes == 0 ) {
 738           // Count nb of nodes
 739           while ( elemIterator->more() ) {
 740             const SMDS_MeshElement* anElem = elemIterator->next();
 741             const SMDS_VtkVolume *aPolyedre = dynamic_cast<const SMDS_VtkVolume*>(anElem);
 742             if ( !aPolyedre ) continue;
 743             nbPolyhedronNodes += aPolyedre->NbNodes();
 744             nbPolyhedronFaces += aPolyedre->NbFaces();
 745             if ( ++iElem == aElemTypeData->_nbElems )
 746               break;
 747           }
 748         }
 749         else {
 750           // Store in med file
 751           PPolyedreInfo aPolyhInfo = myMed->CrPolyedreInfo(aMeshInfo,
 752                                                            aElemTypeData->_entity,
 753                                                            aElemTypeData->_geomType,
 754                                                            aElemTypeData->_nbElems,
 755                                                            nbPolyhedronFaces+1,
 756                                                            nbPolyhedronNodes,
 757                                                            theConnMode,
 758                                                            theIsElemNum,
 759                                                            theIsElemNames);
 760           TElemNum & index = *(aPolyhInfo->myIndex.get());
 761           TElemNum & faces = *(aPolyhInfo->myFaces.get());
 762           TElemNum & conn  = *(aPolyhInfo->myConn.get());
 763           index[0] = 1;
 764           faces[0] = 1;
 765
 766           TInt iFace = 0, iNode = 0;
 767           while ( elemIterator->more() )
 768           {
 769             const SMDS_MeshElement* anElem = elemIterator->next();
 770             const SMDS_VtkVolume *aPolyedre = dynamic_cast<const SMDS_VtkVolume*>(anElem);
 771             if ( !aPolyedre ) continue;
 772             // index
 773             TInt aNbFaces = aPolyedre->NbFaces();
 774             index[ iElem+1 ] = index[ iElem ] + aNbFaces;
 775
 776             // face index
 777             for (TInt f = 1; f <= aNbFaces; ++f, ++iFace ) {
 778               int aNbFaceNodes = aPolyedre->NbFaceNodes( f );
 779               faces[ iFace+1 ] = faces[ iFace ] + aNbFaceNodes;
 780             }
 781             // connectivity
 782             SMDS_ElemIteratorPtr nodeIt = anElem->nodesIterator();
 783             while ( nodeIt->more() ) {
 784               const SMDS_MeshElement* aNode = nodeIt->next();
 785 #ifdef _EDF_NODE_IDS_
 786               conn[ iNode ] = aNodeIdMap[aNode->GetID()];
 787 #else
 788               conn[ iNode ] = aNode->GetID();
 789 #endif
 790               ++iNode;
 791             }
 792             // element number
 793             aPolyhInfo->SetElemNum( iElem, anElem->GetID() );
 794
 795             // family number
 796             int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 797             aPolyhInfo->SetFamNum( iElem, famNum );
 798
 799             if ( ++iElem == aPolyhInfo->GetNbElem() )
 800               break;
 801           }
 802           myMed->SetPolyedreInfo(aPolyhInfo);
 803         }
 804       } // if (aElemTypeData->_geomType == ePOLYEDRE )
 805
 806       // Treat BALLs
 807       // ----------------
 808       else if (aElemTypeData->_geomType == eBALL )
 809       {
 810         // allocate data arrays
 811         PBallInfo aBallInfo = myMed->CrBallInfo( aMeshInfo,
 812                                                  aElemTypeData->_nbElems );
 813
 814         // build map of family numbers for this type
 815         if ( !isElemFamMapBuilt[ aElemTypeData->_smdsType ])
 816         {
 817           fillElemFamilyMap( anElemFamMap, aFamilies, aElemTypeData->_smdsType );
 818           isElemFamMapBuilt[ aElemTypeData->_smdsType ] = true;
 819         }
 820
 821         elemIterator = myMesh->elementsIterator( SMDSAbs_Ball );
 822         while ( elemIterator->more() )
 823         {
 824           const SMDS_MeshElement* anElem = elemIterator->next();
 825           // connectivity
 826           const SMDS_MeshElement* aNode = anElem->GetNode( 0 );
 827 #ifdef _EDF_NODE_IDS_
 828           (*aBallInfo->myConn)[ iElem ] = aNodeIdMap[aNode->GetID()];
 829 #else
 830           (*aBallInfo->myConn)[ iElem ] = aNode->GetID();
 831 #endif
 832           // element number
 833           aBallInfo->SetElemNum( iElem, anElem->GetID() );
 834
 835           // diameter
 836           aBallInfo->myDiameters[ iElem ] =
 837             static_cast<const SMDS_BallElement*>( anElem )->GetDiameter();
 838
 839           // family number
 840           int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 841           aBallInfo->SetFamNum( iElem, famNum );
 842           ++iElem;
 843         }
 844         // store data in a file
 845         myMed->SetBallInfo(aBallInfo);
 846       }
 847
 848       else
 849       {
 850         // Treat standard types
 851         // ---------------------
 852
 853         // allocate data arrays
 854         PCellInfo aCellInfo = myMed->CrCellInfo( aMeshInfo,
 855                                                  aElemTypeData->_entity,
 856                                                  aElemTypeData->_geomType,
 857                                                  aElemTypeData->_nbElems,
 858                                                  theConnMode,
 859                                                  theIsElemNum,
 860                                                  theIsElemNames);
 861         // build map of family numbers for this type
 862         if ( !isElemFamMapBuilt[ aElemTypeData->_smdsType ])
 863         {
 864           //cout << " fillElemFamilyMap()" << endl;
 865           fillElemFamilyMap( anElemFamMap, aFamilies, aElemTypeData->_smdsType );
 866           isElemFamMapBuilt[ aElemTypeData->_smdsType ] = true;
 867         }
 868
 869         TInt aNbNodes = MED::GetNbNodes(aElemTypeData->_geomType);
 870         elemIterator = myMesh->elementsIterator( aElemTypeData->_smdsType );
 871         while ( elemIterator->more() )
 872         {
 873           const SMDS_MeshElement* anElem = elemIterator->next();
 874           if ( anElem->NbNodes() != aNbNodes || anElem->IsPoly() )
 875             continue; // other geometry
 876
 877           // connectivity
 878           TConnSlice aTConnSlice = aCellInfo->GetConnSlice( iElem );
 879           for (TInt iNode = 0; iNode < aNbNodes; iNode++) {
 880             const SMDS_MeshElement* aNode = anElem->GetNode( iNode );
 881 #ifdef _EDF_NODE_IDS_
 882             aTConnSlice[ iNode ] = aNodeIdMap[aNode->GetID()];
 883 #else
 884             aTConnSlice[ iNode ] = aNode->GetID();
 885 #endif
 886           }
 887           // element number
 888           aCellInfo->SetElemNum( iElem, anElem->GetID() );
 889
 890           // family number
 891           int famNum = getFamilyId( anElemFamMap, anElem, defaultFamilyId );
 892           aCellInfo->SetFamNum( iElem, famNum );
 893
 894           if ( ++iElem == aCellInfo->GetNbElem() )
 895             break;
 896         }
 897         // store data in a file
 898         myMed->SetCellInfo(aCellInfo);
 899       }
 900
 901     } // loop on geom types
 902
 903
 904   }
 905   catch(const std::exception& exc) {
 906     INFOS("The following exception was caught:\n\t"<<exc.what());
 907     throw;
 908   }
 909   catch(...) {
 910     INFOS("Unknown exception was caught !!!");
 911     throw;
 912   }
 913
 914   myMeshId = -1;
 915   myGroups.clear();
 916   mySubMeshes.clear();
 917   return aResult;
 918 }