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