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