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