Salome HOME
Added support for polygones and polyhedres
[modules/smesh.git] / src / DriverMED / DriverMED_R_SMESHDS_Mesh.cxx
1 //  SMESH DriverMED : driver to read and write 'med' files
2 //
3 //  Copyright (C) 2003  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.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org 
21 //
22 //
23 //
24 //  File   : DriverMED_R_SMESHDS_Mesh.cxx
25 //  Module : SMESH
26
27 #include "DriverMED_R_SMESHDS_Mesh.h"
28 #include "DriverMED_R_SMDS_Mesh.h"
29 #include "SMESHDS_Mesh.hxx"
30 #include "utilities.h"
31
32 #include "DriverMED_Family.h"
33
34 #include "SMESHDS_Group.hxx"
35
36 #include "MED_Factory.hxx"
37 #include "MED_Utilities.hxx"
38
39 #include <stdlib.h>
40
41 #ifdef _DEBUG_
42 static int MYDEBUG = 0;
43 #else
44 static int MYDEBUG = 0;
45 #endif
46
47 #define _EDF_NODE_IDS_
48
49 using namespace MED;
50
51 void DriverMED_R_SMESHDS_Mesh::SetMeshName(string theMeshName)
52 {
53   myMeshName = theMeshName;
54 }
55
56 static const SMDS_MeshNode* 
57 FindNode(const SMDS_Mesh* theMesh, TInt theId){
58   const SMDS_MeshNode* aNode = theMesh->FindNode(theId);
59   if(aNode) return aNode;
60   EXCEPTION(runtime_error,"SMDS_Mesh::FindNode - cannot find a SMDS_MeshNode for ID = "<<theId);
61 }
62
63
64 enum ECoordName{eX, eY, eZ, eNone};
65 typedef TFloat (*TGetCoord)(MED::PNodeInfo&, TInt);
66
67 template<ECoordName TheCoordId>
68 TFloat GetCoord(MED::PNodeInfo& thePNodeInfo, TInt theElemId){
69   return thePNodeInfo->GetNodeCoord(theElemId,TheCoordId);
70 }
71
72 template<>
73 TFloat GetCoord<eNone>(MED::PNodeInfo& thePNodeInfo, TInt theElemId){
74   return 0.0;
75 }
76
77
78 static TGetCoord aXYZGetCoord[3] = {
79   &GetCoord<eX>, 
80   &GetCoord<eY>, 
81   &GetCoord<eZ>
82 };
83
84
85 static TGetCoord aXYGetCoord[3] = {
86   &GetCoord<eX>, 
87   &GetCoord<eY>, 
88   &GetCoord<eNone>
89 };
90
91 static TGetCoord aYZGetCoord[3] = {
92   &GetCoord<eNone>,
93   &GetCoord<eX>, 
94   &GetCoord<eY>
95 };
96
97 static TGetCoord aXZGetCoord[3] = {
98   &GetCoord<eX>, 
99   &GetCoord<eNone>,
100   &GetCoord<eY>
101 };
102
103
104 static TGetCoord aXGetCoord[3] = {
105   &GetCoord<eX>, 
106   &GetCoord<eNone>,
107   &GetCoord<eNone>
108 };
109
110 static TGetCoord aYGetCoord[3] = {
111   &GetCoord<eNone>,
112   &GetCoord<eX>, 
113   &GetCoord<eNone>
114 };
115
116 static TGetCoord aZGetCoord[3] = {
117   &GetCoord<eNone>,
118   &GetCoord<eNone>,
119   &GetCoord<eX>
120 };
121
122
123 class TCoordHelper{
124   MED::PNodeInfo myPNodeInfo;
125   TGetCoord* myGetCoord;
126 public:
127   TCoordHelper(const MED::PNodeInfo& thePNodeInfo,
128                TGetCoord* theGetCoord):
129     myPNodeInfo(thePNodeInfo),
130     myGetCoord(theGetCoord)
131   {}
132   virtual ~TCoordHelper(){}
133   TFloat GetCoord(TInt theElemId, TInt theCoodId){
134     return (*myGetCoord[theCoodId])(myPNodeInfo,theElemId);
135   }
136 };
137 typedef boost::shared_ptr<TCoordHelper> TCoordHelperPtr;
138
139
140 Driver_Mesh::Status DriverMED_R_SMESHDS_Mesh::Perform()
141 {
142   Status aResult = DRS_FAIL;
143   try{
144     myFamilies.clear();
145     if(MYDEBUG) MESSAGE("Perform - myFile : "<<myFile);
146     PWrapper aMed = CrWrapper(myFile);
147
148     aResult = DRS_EMPTY;
149     if(TInt aNbMeshes = aMed->GetNbMeshes()){
150       for(int iMesh = 0; iMesh < aNbMeshes; iMesh++){
151         // Reading the MED mesh
152         //---------------------
153         PMeshInfo aMeshInfo = aMed->GetPMeshInfo(iMesh+1);
154         string aMeshName;
155         if (myMeshId != -1) {
156           ostringstream aMeshNameStr;
157           aMeshNameStr<<myMeshId;
158           aMeshName = aMeshNameStr.str();
159         } else {
160           aMeshName = myMeshName;
161         }
162         if(MYDEBUG) MESSAGE("Perform - aMeshName : "<<aMeshName<<"; "<<aMeshInfo->GetName());
163         if(aMeshName != aMeshInfo->GetName()) continue;
164         aResult = DRS_OK;
165         //TInt aMeshDim = aMeshInfo->GetDim();
166         
167         // Reading MED families to the temporary structure
168         //------------------------------------------------
169         TErr anErr;
170         TInt aNbFams = aMed->GetNbFamilies(aMeshInfo);
171         if(MYDEBUG) MESSAGE("Read " << aNbFams << " families");
172         for (TInt iFam = 0; iFam < aNbFams; iFam++) {
173           PFamilyInfo aFamilyInfo = aMed->GetPFamilyInfo(aMeshInfo,iFam+1,&anErr);
174           if(anErr >= 0){
175             TInt aFamId = aFamilyInfo->GetId();
176             if(MYDEBUG) MESSAGE("Family " << aFamId << " :");
177             
178             DriverMED_FamilyPtr aFamily (new DriverMED_Family);
179             
180             TInt aNbGrp = aFamilyInfo->GetNbGroup();
181             if(MYDEBUG) MESSAGE("belong to " << aNbGrp << " groups");
182             for (TInt iGr = 0; iGr < aNbGrp; iGr++) {
183               string aGroupName = aFamilyInfo->GetGroupName(iGr);
184               if(MYDEBUG) MESSAGE(aGroupName);
185               aFamily->AddGroupName(aGroupName);
186             }
187             myFamilies[aFamId] = aFamily;
188           }
189         }
190
191         // Reading MED nodes to the corresponding SMDS structure
192         //------------------------------------------------------
193         PNodeInfo aNodeInfo = aMed->GetPNodeInfo(aMeshInfo);
194
195         TCoordHelperPtr aCoordHelperPtr;
196         {
197           TInt aMeshDimension = aMeshInfo->GetDim();
198           bool anIsDimPresent[3] = {false, false, false};
199           for(TInt iDim = 0; iDim < aMeshDimension; iDim++){
200             string aDimName = aNodeInfo->GetCoordName(iDim);
201             if(aDimName == "x" || aDimName == "X")
202               anIsDimPresent[eX] = true;
203             else if(aDimName == "y" || aDimName == "Y")
204               anIsDimPresent[eY] = true;
205             else if(aDimName == "z" || aDimName == "Z")
206               anIsDimPresent[eZ] = true;
207           }
208           switch(aMeshDimension){
209           case 3:
210             aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aXYZGetCoord));
211             break;
212           case 2:
213             if(anIsDimPresent[eY] && anIsDimPresent[eZ])
214               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aYZGetCoord));
215             else if(anIsDimPresent[eX] && anIsDimPresent[eZ])
216               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aXZGetCoord));
217             else
218               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aXYGetCoord));
219             break;
220           case 1:
221             if(anIsDimPresent[eY])
222               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aYGetCoord));
223             else if(anIsDimPresent[eZ])
224               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aZGetCoord));
225             else
226               aCoordHelperPtr.reset(new TCoordHelper(aNodeInfo,aXGetCoord));
227             break;
228           }
229         }
230
231         EBooleen anIsNodeNum = aNodeInfo->IsElemNum();
232         TInt aNbElems = aNodeInfo->GetNbElem();
233         if(MYDEBUG) MESSAGE("Perform - aNodeInfo->GetNbElem() = "<<aNbElems<<"; anIsNodeNum = "<<anIsNodeNum);
234         for(TInt iElem = 0; iElem < aNbElems; iElem++){
235           double aCoords[3] = {0.0, 0.0, 0.0};
236           for(TInt iDim = 0; iDim < 3; iDim++)
237             aCoords[iDim] = aCoordHelperPtr->GetCoord(iElem,iDim);
238           const SMDS_MeshNode* aNode;
239           if(anIsNodeNum) {
240             aNode = myMesh->AddNodeWithID
241               (aCoords[0],aCoords[1],aCoords[2],aNodeInfo->GetElemNum(iElem));
242           } else {
243             aNode = myMesh->AddNode
244               (aCoords[0],aCoords[1],aCoords[2]);
245           }
246           //cout<<aNode->GetID()<<": "<<aNode->X()<<", "<<aNode->Y()<<", "<<aNode->Z()<<endl;
247
248           // Save reference to this node from its family
249           TInt aFamNum = aNodeInfo->GetFamNum(iElem);
250           if (myFamilies.find(aFamNum) != myFamilies.end())
251           {
252             myFamilies[aFamNum]->AddElement(aNode);
253             myFamilies[aFamNum]->SetType(SMDSAbs_Node);
254           }
255         }
256
257         // Reading pre information about all MED cells
258         //--------------------------------------------
259         bool takeNumbers = true;  // initially we trust the numbers from file
260         MED::TEntityInfo aEntityInfo = aMed->GetEntityInfo(aMeshInfo);
261         MED::TEntityInfo::iterator anEntityIter = aEntityInfo.begin();
262         for(; anEntityIter != aEntityInfo.end(); anEntityIter++){
263           const EEntiteMaillage& anEntity = anEntityIter->first;
264           if(anEntity == eNOEUD) continue;
265           // Reading MED cells to the corresponding SMDS structure
266           //------------------------------------------------------
267           const MED::TGeom& aTGeom = anEntityIter->second;
268           MED::TGeom::const_iterator anTGeomIter = aTGeom.begin();
269           for(; anTGeomIter != aTGeom.end(); anTGeomIter++){
270             const EGeometrieElement& aGeom = anTGeomIter->first;
271
272             if (aGeom == ePOINT1) {
273               continue;
274
275             } else if (aGeom == ePOLYGONE) {
276               PPolygoneInfo aPolygoneInfo = aMed->GetPPolygoneInfo(aMeshInfo,anEntity,aGeom);
277               EBooleen anIsElemNum = takeNumbers ? aPolygoneInfo->IsElemNum() : eFAUX;
278
279               TElemNum aConn  = aPolygoneInfo->GetConnectivite();
280               TElemNum aIndex = aPolygoneInfo->GetIndex();
281
282               TInt nbPolygons = aPolygoneInfo->GetNbElem();
283               for (TInt iPG = 0; iPG < nbPolygons; iPG++) {
284                 // get nodes
285                 TInt aCurrPG_FirstNodeIndex = aIndex[iPG] - 1;
286                 int aNbNodes = aPolygoneInfo->GetNbConn(iPG);
287                 std::vector<int> nodes_ids (aNbNodes);
288                 for (TInt inode = 0; inode < aNbNodes; inode++) {
289                   nodes_ids[inode] = aConn[aCurrPG_FirstNodeIndex + inode];
290                 }
291
292                 bool isRenum = false;
293                 SMDS_MeshElement* anElement = NULL;
294                 TInt aFamNum = aPolygoneInfo->GetFamNum(iPG);
295
296                 try {
297                   if (anIsElemNum) {
298                     anElement = myMesh->AddPolygonalFaceWithID
299                       (nodes_ids, aPolygoneInfo->GetElemNum(iPG));
300                   }
301                   if (!anElement) {
302                     std::vector<const SMDS_MeshNode*> nodes (aNbNodes);
303                     for (int inode = 0; inode < aNbNodes; inode++) {
304                       nodes[inode] = FindNode(myMesh, nodes_ids[inode]);
305                     }
306                     anElement = myMesh->AddPolygonalFace(nodes);
307                     isRenum = anIsElemNum;
308                   }
309                   
310                 } catch (const std::exception& exc) {
311                   aResult = DRS_FAIL;
312                   INFOS("Follow exception was cought:\n\t"<<exc.what());
313                 } catch (...) {
314                   aResult = DRS_FAIL;
315                   INFOS("Follow unknown exception was cought!");
316                 }
317
318                 if (!anElement) {
319                   aResult = DRS_WARN_SKIP_ELEM;
320                 } else {
321                   if (isRenum) {
322                     anIsElemNum = eFAUX;
323                     takeNumbers = false;
324                     if (aResult < DRS_WARN_RENUMBER)
325                       aResult = DRS_WARN_RENUMBER;
326                   }
327                   if (myFamilies.find(aFamNum) != myFamilies.end()) {
328                     // Save reference to this element from its family
329                     if (MYDEBUG){ 
330                       cout<<"myFamilies["<<aFamNum
331                           <<"] IsPoly()="<<anElement->IsPoly()
332                           <<"; GetType="<<anElement->GetType()
333                           <<"; aNbNodes="<<aNbNodes<<endl;                    
334                     }
335                     myFamilies[aFamNum]->AddElement(anElement);
336                     myFamilies[aFamNum]->SetType(anElement->GetType());
337                   }
338                 }
339               }
340               continue;
341
342             } else if (aGeom == ePOLYEDRE) {
343               PPolyedreInfo aPolyedreInfo = aMed->GetPPolyedreInfo(aMeshInfo,anEntity,aGeom);
344               EBooleen anIsElemNum = takeNumbers ? aPolyedreInfo->IsElemNum() : eFAUX;
345
346               TElemNum aConn       = aPolyedreInfo->GetConnectivite();
347               TElemNum aFacesIndex = aPolyedreInfo->GetFacesIndex();
348               TElemNum aIndex      = aPolyedreInfo->GetIndex();
349
350               TInt nbPolyedres = aPolyedreInfo->GetNbElem();
351               
352               for (int iPE = 0; iPE < nbPolyedres; iPE++) {
353                 // get faces
354                 int aCurrPE_FirstFaceIndex = aIndex[iPE] - 1;
355                 int aNextPE_FirstFaceIndex = aIndex[iPE + 1] - 1;
356                 int nbFaces = aNextPE_FirstFaceIndex - aCurrPE_FirstFaceIndex;
357                 std::vector<int> quantities (nbFaces);
358                 for (int iFa = 0; iFa < nbFaces; iFa++) {
359                   int aCurrFace_FirstNodeIndex = aFacesIndex[aCurrPE_FirstFaceIndex + iFa] - 1;
360                   int aNextFace_FirstNodeIndex = aFacesIndex[aCurrPE_FirstFaceIndex + iFa + 1] - 1;
361
362                   int aNbNodes = aNextFace_FirstNodeIndex - aCurrFace_FirstNodeIndex;
363                   quantities[iFa] = aNbNodes;
364                 }
365
366                 // get nodes
367                 int aCurrPE_FirstNodeIndex = aFacesIndex[aCurrPE_FirstFaceIndex] - 1;
368                 int nbPENodes = aPolyedreInfo->GetNbConn(iPE);
369                 std::vector<int> nodes_ids (nbPENodes);
370                 for (int inode = 0; inode < nbPENodes; inode++) {
371                   nodes_ids[inode] = aConn[aCurrPE_FirstNodeIndex + inode];
372                 }
373
374                 bool isRenum = false;
375                 SMDS_MeshElement* anElement = NULL;
376                 TInt aFamNum = aPolyedreInfo->GetFamNum(iPE);
377
378                 try {
379                   if (anIsElemNum) {
380                     anElement = myMesh->AddPolyhedralVolumeWithID
381                       (nodes_ids, quantities, aPolyedreInfo->GetElemNum(iPE));
382                   }
383                   if (!anElement) {
384                     std::vector<const SMDS_MeshNode*> nodes (nbPENodes);
385                     for (int inode = 0; inode < nbPENodes; inode++) {
386                       nodes[inode] = FindNode(myMesh, nodes_ids[inode]);
387                     }
388                     anElement = myMesh->AddPolyhedralVolume(nodes, quantities);
389                     isRenum = anIsElemNum;
390                   }
391                 } catch (const std::exception& exc) {
392                   aResult = DRS_FAIL;
393                   INFOS("Follow exception was cought:\n\t"<<exc.what());
394                 } catch (...) {
395                   aResult = DRS_FAIL;
396                   INFOS("Follow unknown exception was cought!");
397                 }
398
399                 if (!anElement) {
400                   aResult = DRS_WARN_SKIP_ELEM;
401                 } else {
402                   if (isRenum) {
403                     anIsElemNum = eFAUX;
404                     takeNumbers = false;
405                     if (aResult < DRS_WARN_RENUMBER)
406                       aResult = DRS_WARN_RENUMBER;
407                   }
408                   if (myFamilies.find(aFamNum) != myFamilies.end()) {
409                     // Save reference to this element from its family
410                     if (MYDEBUG){
411                       cout<<"myFamilies["<<aFamNum
412                           <<"] IsPoly()="<<anElement->IsPoly()
413                           <<"; GetType="<<anElement->GetType()
414                           <<"; aNbNodes="<<nbPENodes<<endl;
415                     }
416                     myFamilies[aFamNum]->AddElement(anElement);
417                     myFamilies[aFamNum]->SetType(anElement->GetType());
418                   }
419                 }
420               }
421               continue;
422
423             }
424
425             PCellInfo aCellInfo = aMed->GetPCellInfo(aMeshInfo,anEntity,aGeom);
426             EBooleen anIsElemNum = takeNumbers ? aCellInfo->IsElemNum() : eFAUX;
427             TInt aNbElems = aCellInfo->GetNbElem();
428             if(MYDEBUG) MESSAGE("Perform - anEntity = "<<anEntity<<"; anIsElemNum = "<<anIsElemNum);
429             if(MYDEBUG) MESSAGE("Perform - aGeom = "<<aGeom<<"; aNbElems = "<<aNbElems);
430
431             for(int iElem = 0; iElem < aNbElems; iElem++){
432               TInt aNbNodes = -1;
433               switch(aGeom){
434               case eSEG2:
435               case eSEG3:
436                 aNbNodes = 2;
437                 break;
438               case eTRIA3:
439               case eTRIA6:
440                 aNbNodes = 3;
441                 break;
442                 break;
443               case eQUAD4:
444               case eQUAD8:
445                 aNbNodes = 4;
446                 break;
447               case eTETRA4:
448               case eTETRA10:
449                 aNbNodes = 4;
450                 break;
451               case ePYRA5:
452               case ePYRA13:
453                 aNbNodes = 5;
454                 break;
455               case ePENTA6:
456               case ePENTA15:
457                 aNbNodes = 6;
458                 break;
459               case eHEXA8:
460               case eHEXA20:
461                 aNbNodes = 8;
462                 break;
463               }
464               vector<TInt> aNodeIds(aNbNodes);
465               bool anIsValidConnect = false;
466
467               try{
468 #ifdef _EDF_NODE_IDS_
469                 if(anIsNodeNum) {
470                   for(int i = 0; i < aNbNodes; i++){
471                     aNodeIds[i] = aNodeInfo->GetElemNum(aCellInfo->GetConn(iElem,i)-1);
472                   }
473                 }else{
474                   for(int i = 0; i < aNbNodes; i++){
475                     aNodeIds[i] = aCellInfo->GetConn(iElem,i);
476                   }
477                 }
478 #else
479                 for(int i = 0; i < aNbNodes; i++){
480                   aNodeIds[i] = aCellInfo->GetConn(iElem,i);
481                 }
482 #endif
483                 anIsValidConnect = true;
484               }catch(const std::exception& exc){
485                 //INFOS("Follow exception was cought:\n\t"<<exc.what());
486                 aResult = DRS_FAIL;
487               }catch(...){
488                 //INFOS("Unknown exception was cought !!!");
489                 aResult = DRS_FAIL;
490               }
491               
492               if(!anIsValidConnect)
493                 continue;
494
495               bool isRenum = false;
496               SMDS_MeshElement* anElement = NULL;
497               TInt aFamNum = aCellInfo->GetFamNum(iElem);
498               try{
499                 switch(aGeom){
500                 case eSEG2:
501                 case eSEG3:
502                   if(anIsElemNum)
503                     anElement = myMesh->AddEdgeWithID(aNodeIds[0],
504                                                       aNodeIds[1],
505                                                       aCellInfo->GetElemNum(iElem));
506                   if (!anElement) {
507                     anElement = myMesh->AddEdge(FindNode(myMesh,aNodeIds[0]),
508                                                 FindNode(myMesh,aNodeIds[1]));
509                     isRenum = anIsElemNum;
510                   }
511                   break;
512                 case eTRIA3:
513                 case eTRIA6:
514                   aNbNodes = 3;
515                   if(anIsElemNum)
516                     anElement = myMesh->AddFaceWithID(aNodeIds[0],
517                                                       aNodeIds[1],
518                                                       aNodeIds[2],
519                                                       aCellInfo->GetElemNum(iElem));
520                   if (!anElement) {
521                     anElement = myMesh->AddFace(FindNode(myMesh,aNodeIds[0]),
522                                                 FindNode(myMesh,aNodeIds[1]),
523                                                 FindNode(myMesh,aNodeIds[2]));
524                     isRenum = anIsElemNum;
525                   }
526                   break;
527                 case eQUAD4:
528                 case eQUAD8:
529                   aNbNodes = 4;
530                   // There is some differnce between SMDS and MED
531                   if(anIsElemNum)
532                     anElement = myMesh->AddFaceWithID(aNodeIds[0],
533                                                       aNodeIds[1],
534                                                       aNodeIds[2],
535                                                       aNodeIds[3],
536                                                       aCellInfo->GetElemNum(iElem));
537                   if (!anElement) {
538                     anElement = myMesh->AddFace(FindNode(myMesh,aNodeIds[0]),
539                                                 FindNode(myMesh,aNodeIds[1]),
540                                                 FindNode(myMesh,aNodeIds[2]),
541                                                 FindNode(myMesh,aNodeIds[3]));
542                     isRenum = anIsElemNum;
543                   }
544                   break;
545                 case eTETRA4:
546                 case eTETRA10:
547                   aNbNodes = 4;
548                   if(anIsElemNum)
549                     anElement = myMesh->AddVolumeWithID(aNodeIds[0],
550                                                         aNodeIds[1],
551                                                         aNodeIds[2],
552                                                         aNodeIds[3],
553                                                         aCellInfo->GetElemNum(iElem));
554                   if (!anElement) {
555                     anElement = myMesh->AddVolume(FindNode(myMesh,aNodeIds[0]),
556                                                   FindNode(myMesh,aNodeIds[1]),
557                                                   FindNode(myMesh,aNodeIds[2]),
558                                                   FindNode(myMesh,aNodeIds[3]));
559                     isRenum = anIsElemNum;
560                   }
561                   break;
562                 case ePYRA5:
563                 case ePYRA13:
564                   aNbNodes = 5;
565                   // There is some differnce between SMDS and MED
566                   if(anIsElemNum)
567                     anElement = myMesh->AddVolumeWithID(aNodeIds[0],
568                                                         aNodeIds[1],
569                                                         aNodeIds[2],
570                                                         aNodeIds[3],
571                                                         aNodeIds[4],
572                                                         aCellInfo->GetElemNum(iElem));
573                   if (!anElement) {
574                     anElement = myMesh->AddVolume(FindNode(myMesh,aNodeIds[0]),
575                                                   FindNode(myMesh,aNodeIds[1]),
576                                                   FindNode(myMesh,aNodeIds[2]),
577                                                   FindNode(myMesh,aNodeIds[3]),
578                                                   FindNode(myMesh,aNodeIds[4]));
579                     isRenum = anIsElemNum;
580                   }
581                   break;
582                 case ePENTA6:
583                 case ePENTA15:
584                   aNbNodes = 6;
585                   if(anIsElemNum)
586                     anElement = myMesh->AddVolumeWithID(aNodeIds[0],
587                                                         aNodeIds[1],
588                                                         aNodeIds[2],
589                                                         aNodeIds[3],
590                                                         aNodeIds[4],
591                                                         aNodeIds[5],
592                                                         aCellInfo->GetElemNum(iElem));
593                   if (!anElement) {
594                     anElement = myMesh->AddVolume(FindNode(myMesh,aNodeIds[0]),
595                                                   FindNode(myMesh,aNodeIds[1]),
596                                                   FindNode(myMesh,aNodeIds[2]),
597                                                   FindNode(myMesh,aNodeIds[3]),
598                                                   FindNode(myMesh,aNodeIds[4]),
599                                                   FindNode(myMesh,aNodeIds[5]));
600                     isRenum = anIsElemNum;
601                   }
602                   break;
603                 case eHEXA8:
604                 case eHEXA20:
605                   aNbNodes = 8;
606                   if(anIsElemNum)
607                     anElement = myMesh->AddVolumeWithID(aNodeIds[0],
608                                                         aNodeIds[1],
609                                                         aNodeIds[2],
610                                                         aNodeIds[3],
611                                                         aNodeIds[4],
612                                                         aNodeIds[5],
613                                                         aNodeIds[6],
614                                                         aNodeIds[7],
615                                                         aCellInfo->GetElemNum(iElem));
616                   if (!anElement) {
617                     anElement = myMesh->AddVolume(FindNode(myMesh,aNodeIds[0]),
618                                                   FindNode(myMesh,aNodeIds[1]),
619                                                   FindNode(myMesh,aNodeIds[2]),
620                                                   FindNode(myMesh,aNodeIds[3]),
621                                                   FindNode(myMesh,aNodeIds[4]),
622                                                   FindNode(myMesh,aNodeIds[5]),
623                                                   FindNode(myMesh,aNodeIds[6]),
624                                                   FindNode(myMesh,aNodeIds[7]));
625                     isRenum = anIsElemNum;
626                   }
627                   break;
628                 }
629               }catch(const std::exception& exc){
630                 //INFOS("Follow exception was cought:\n\t"<<exc.what());
631                 aResult = DRS_FAIL;
632               }catch(...){
633                 //INFOS("Unknown exception was cought !!!");
634                 aResult = DRS_FAIL;
635               }
636                 
637               if (!anElement) {
638                 aResult = DRS_WARN_SKIP_ELEM;
639               }
640               else {
641                 if (isRenum) {
642                   anIsElemNum = eFAUX;
643                   takeNumbers = false;
644                   if (aResult < DRS_WARN_RENUMBER)
645                     aResult = DRS_WARN_RENUMBER;
646                 }
647                 if (myFamilies.find(aFamNum) != myFamilies.end()) {
648                   // Save reference to this element from its family
649                   myFamilies[aFamNum]->AddElement(anElement);
650                   myFamilies[aFamNum]->SetType(anElement->GetType());
651                 }
652               }
653             }
654           }
655         }
656         break;
657       }
658     }
659   }catch(const std::exception& exc){
660     INFOS("Follow exception was cought:\n\t"<<exc.what());
661     aResult = DRS_FAIL;
662   }catch(...){
663     INFOS("Unknown exception was cought !!!");
664     aResult = DRS_FAIL;
665   }
666   if(MYDEBUG) MESSAGE("Perform - aResult status = "<<aResult);
667   return aResult;
668 }
669
670 list<string> DriverMED_R_SMESHDS_Mesh::GetMeshNames(Status& theStatus)
671 {
672   list<string> aMeshNames;
673
674   try {
675     if(MYDEBUG) MESSAGE("GetMeshNames - myFile : " << myFile);
676     theStatus = DRS_OK;
677     PWrapper aMed = CrWrapper(myFile);
678
679     if (TInt aNbMeshes = aMed->GetNbMeshes()) {
680       for (int iMesh = 0; iMesh < aNbMeshes; iMesh++) {
681         // Reading the MED mesh
682         //---------------------
683         PMeshInfo aMeshInfo = aMed->GetPMeshInfo(iMesh+1);
684         aMeshNames.push_back(aMeshInfo->GetName());
685       }
686     }
687   }catch(const std::exception& exc){
688     INFOS("Follow exception was cought:\n\t"<<exc.what());
689     theStatus = DRS_FAIL;
690   }catch(...){
691     INFOS("Unknown exception was cought !!!");
692     theStatus = DRS_FAIL;
693   }
694
695   return aMeshNames;
696 }
697
698 list<string> DriverMED_R_SMESHDS_Mesh::GetGroupNames()
699 {
700   list<string> aResult;
701   set<string> aResGroupNames;
702
703   map<int, DriverMED_FamilyPtr>::iterator aFamsIter = myFamilies.begin();
704   for (; aFamsIter != myFamilies.end(); aFamsIter++)
705   {
706     DriverMED_FamilyPtr aFamily = (*aFamsIter).second;
707     const MED::TStringSet& aGroupNames = aFamily->GetGroupNames();
708     set<string>::iterator aGrNamesIter = aGroupNames.begin();
709     for (; aGrNamesIter != aGroupNames.end(); aGrNamesIter++)
710     {
711       string aName = *aGrNamesIter;
712       // Check, if this is a Group or SubMesh name
713 //if (aName.substr(0, 5) == string("Group")) {
714         if (aResGroupNames.find(aName) == aResGroupNames.end()) {
715           aResGroupNames.insert(aName);
716           aResult.push_back(aName);
717         }
718 //    }
719     }
720   }
721
722   return aResult;
723 }
724
725 void DriverMED_R_SMESHDS_Mesh::GetGroup(SMESHDS_Group* theGroup)
726 {
727   string aGroupName (theGroup->GetStoreName());
728   if(MYDEBUG) MESSAGE("Get Group " << aGroupName);
729
730   map<int, DriverMED_FamilyPtr>::iterator aFamsIter = myFamilies.begin();
731   for (; aFamsIter != myFamilies.end(); aFamsIter++)
732   {
733     DriverMED_FamilyPtr aFamily = (*aFamsIter).second;
734     if (aFamily->MemberOf(aGroupName))
735     {
736       const set<const SMDS_MeshElement *>& anElements = aFamily->GetElements();
737       set<const SMDS_MeshElement *>::iterator anElemsIter = anElements.begin();
738       const SMDS_MeshElement * element = 0;
739       for (; anElemsIter != anElements.end(); anElemsIter++)
740       {
741         element = *anElemsIter;
742         theGroup->SMDSGroup().Add(element);
743       }
744       if ( element )
745         theGroup->SetType( element->GetType() );
746     }
747   }
748 }
749
750 void DriverMED_R_SMESHDS_Mesh::GetSubMesh (SMESHDS_SubMesh* theSubMesh,
751                                            const int theId)
752 {
753   char submeshGrpName[ 30 ];
754   sprintf( submeshGrpName, "SubMesh %d", theId );
755   string aName (submeshGrpName);
756   map<int, DriverMED_FamilyPtr>::iterator aFamsIter = myFamilies.begin();
757   for (; aFamsIter != myFamilies.end(); aFamsIter++)
758   {
759     DriverMED_FamilyPtr aFamily = (*aFamsIter).second;
760     if (aFamily->MemberOf(aName))
761     {
762       const set<const SMDS_MeshElement *>& anElements = aFamily->GetElements();
763       set<const SMDS_MeshElement *>::iterator anElemsIter = anElements.begin();
764       if (aFamily->GetType() == SMDSAbs_Node)
765       {
766         for (; anElemsIter != anElements.end(); anElemsIter++)
767         {
768           const SMDS_MeshNode* node = static_cast<const SMDS_MeshNode*>(*anElemsIter);
769           theSubMesh->AddNode(node);
770         }
771       }
772       else
773       {
774         for (; anElemsIter != anElements.end(); anElemsIter++)
775         {
776           theSubMesh->AddElement(*anElemsIter);
777         }
778       }
779     }
780   }
781 }
782
783 void DriverMED_R_SMESHDS_Mesh::CreateAllSubMeshes ()
784 {
785   map<int, DriverMED_FamilyPtr>::iterator aFamsIter = myFamilies.begin();
786   for (; aFamsIter != myFamilies.end(); aFamsIter++)
787   {
788     DriverMED_FamilyPtr aFamily = (*aFamsIter).second;
789     MED::TStringSet aGroupNames = aFamily->GetGroupNames();
790     set<string>::iterator aGrNamesIter = aGroupNames.begin();
791     for (; aGrNamesIter != aGroupNames.end(); aGrNamesIter++)
792     {
793       string aName = *aGrNamesIter;
794       // Check, if this is a Group or SubMesh name
795       if (aName.substr(0, 7) == string("SubMesh"))
796       {
797         int Id = atoi(string(aName).substr(7).c_str());
798         set<const SMDS_MeshElement *> anElements = aFamily->GetElements();
799         set<const SMDS_MeshElement *>::iterator anElemsIter = anElements.begin();
800         if (aFamily->GetType() == SMDSAbs_Node)
801         {
802           for (; anElemsIter != anElements.end(); anElemsIter++)
803           {
804             SMDS_MeshNode* node = const_cast<SMDS_MeshNode*>
805               ( static_cast<const SMDS_MeshNode*>( *anElemsIter ));
806             // find out a shape type
807             TopoDS_Shape aShape = myMesh->IndexToShape( Id );
808             int aShapeType = ( aShape.IsNull() ? -1 : aShape.ShapeType() );
809             switch ( aShapeType ) {
810             case TopAbs_FACE:
811               myMesh->SetNodeOnFace(node, Id); break;
812             case TopAbs_EDGE:
813               myMesh->SetNodeOnEdge(node, Id); break;
814             case TopAbs_VERTEX:
815               myMesh->SetNodeOnVertex(node, Id); break;
816             default:
817               myMesh->SetNodeInVolume(node, Id);
818             }
819           }
820         }
821         else
822         {
823           for (; anElemsIter != anElements.end(); anElemsIter++)
824           {
825             myMesh->SetMeshElementOnShape(*anElemsIter, Id);
826           }
827         }
828       }
829     }
830   }
831 }