]> SALOME platform Git repositories - modules/visu.git/blob - src/CONVERTOR/VISU_MedConvertor.cxx
Salome HOME
ab1e1026099b0db74674acd7f62b841498dc2e80
[modules/visu.git] / src / CONVERTOR / VISU_MedConvertor.cxx
1 //  VISU OBJECT : interactive object for VISU entities implementation
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 //  File   : VISU_MedConvertor.cxx
24 //  Author : Alexey PETROV
25 //  Module : VISU
26
27
28 #include "VISU_MedConvertor.hxx"
29 #include "VISU_ConvertorUtils.hxx"
30
31 #include <valarray>     
32 #include <vtkCellType.h>
33
34 #define USER_INTERLACE MED_FULL_INTERLACE
35
36 using namespace std;
37
38 #ifdef _DEBUG_
39 static int MYDEBUG = 0;
40 #else
41 static int MYDEBUG = 0;
42 #endif
43
44
45 static med_err ret = 0;
46
47 typedef map<VISU::TEntity,med_entite_maillage> TVisu2MedEntity;
48
49 static TVisu2MedEntity aVisu2MedEntity;
50
51 static int INIT = (
52                    aVisu2MedEntity[VISU::NODE_ENTITY] = MED_NOEUD,
53                    aVisu2MedEntity[VISU::EDGE_ENTITY] = MED_ARETE,
54                    aVisu2MedEntity[VISU::FACE_ENTITY] = MED_FACE,
55                    aVisu2MedEntity[VISU::CELL_ENTITY] = MED_MAILLE,
56                    1);
57
58 static med_geometrie_element NODEGEOM[1] = {MED_POINT1};
59
60 static med_geometrie_element EDGEGEOM[MED_NBR_GEOMETRIE_ARETE] = {
61   MED_SEG2, MED_SEG3
62   };
63
64 static med_geometrie_element FACEGEOM[MED_NBR_GEOMETRIE_FACE] = {
65   MED_TRIA3, MED_QUAD4, MED_TRIA6, MED_QUAD8
66   };
67
68 static med_geometrie_element  CELLGEOM[MED_NBR_GEOMETRIE_MAILLE] = {
69   MED_POINT1, MED_SEG2, MED_SEG3, MED_TRIA3,
70   MED_QUAD4, MED_TRIA6, MED_QUAD8, MED_TETRA4,
71   MED_PYRA5, MED_PENTA6, MED_HEXA8, MED_TETRA10, 
72   MED_PYRA13, MED_PENTA15, MED_HEXA20
73   };
74
75 void GetEntity2Geom(const VISU::TEntity& theEntity, med_geometrie_element*& theVector, int* theEnd) 
76      throw (std::runtime_error&)
77 {
78   switch(theEntity){
79   case VISU::CELL_ENTITY: theVector = CELLGEOM; *theEnd = MED_NBR_GEOMETRIE_MAILLE; break;
80   case VISU::FACE_ENTITY: theVector = FACEGEOM; *theEnd = MED_NBR_GEOMETRIE_FACE; break;
81   case VISU::EDGE_ENTITY: theVector = EDGEGEOM; *theEnd = MED_NBR_GEOMETRIE_ARETE; break;
82   case VISU::NODE_ENTITY: theVector = NODEGEOM; *theEnd = 1; break;
83   default:
84     throw std::runtime_error("GetEntity2Geom >> theEntity is uncorrect");
85   }
86 }
87
88 extern "C"
89 VISU_Convertor* CreateMedConvertor(const string& theFileName) throw(std::runtime_error&){
90   return new VISU_MedConvertor(theFileName);
91 }
92
93 class MedFile{
94   char* myFileName;
95   med_idt myFid;
96   MedFile();
97   MedFile(const MedFile&);
98 public:
99   MedFile(const char* theFileName) throw(std::runtime_error&) :
100     myFileName(strdup(theFileName))
101   {
102     myFid = MEDouvrir(myFileName,MED_LECT);
103     if(myFid < 0){
104       free(myFileName);
105       throw std::runtime_error(string("MedFile::MedFile >> MEDouvrir(...) - ") + theFileName);
106     }
107   }
108   ~MedFile(){
109     free(myFileName);
110     if(myFid >= 0) 
111       MEDfermer(myFid);
112   }
113   const med_idt& GetFid() const { return myFid;};
114 };
115
116 VISU_MedConvertor::VISU_MedConvertor(const string& theFileName) throw (std::runtime_error&) {
117   myFileInfo.setFile(QString(theFileName.c_str()));
118   myName = myFileInfo.baseName().latin1();
119 }
120
121 VISU_Convertor* VISU_MedConvertor::Build() throw (std::runtime_error&) {
122   MedFile aMedFile(myFileInfo.absFilePath());
123   med_idt fid = aMedFile.GetFid();
124   med_int iMeshEnd = MEDnMaa(fid);  //Get number of meshes
125   if (iMeshEnd < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDnMaa(...)"));
126   if(MYDEBUG) MESSAGE("ImportInfo - MEDnMaa = "<<iMeshEnd<<"; myFileInfo = "<<myFileInfo.filePath());
127   for(int iMesh = 1; iMesh <= iMeshEnd; iMesh++){
128     med_int aMeshDim;
129     char aMeshName[MED_TAILLE_NOM+1] = "";
130     MEDmaaInfo(fid,iMesh,aMeshName,&aMeshDim);
131     if(MYDEBUG) MESSAGE("ImportInfo - aMeshName = '"<<aMeshName<<"'; aMeshDim = "<<aMeshDim);
132     VISU::TMesh &aMesh = myMeshMap[aMeshName];
133     aMesh.myDim = aMeshDim;
134     aMesh.myName = aMeshName;
135     typedef map<med_int,VISU::TEntity> TFamily2EntityMap;
136     TFamily2EntityMap aFamily2EntityMap;
137     typedef map<med_int,med_int> TFamilyCounterMap;
138     TFamilyCounterMap aFamilyNbCellsCounterMap, aFamilyCellsSizeCounterMap;
139     TVisu2MedEntity::const_iterator aVisu2MedEntityIter = aVisu2MedEntity.begin();
140     for(; aVisu2MedEntityIter != aVisu2MedEntity.end(); aVisu2MedEntityIter++) {
141       VISU::TEntity anEntity = aVisu2MedEntityIter->first;
142       int iGeomElemEnd;
143       med_geometrie_element* aGeomElemVector;
144       GetEntity2Geom(anEntity,aGeomElemVector,&iGeomElemEnd);
145       med_entite_maillage& aMedEntity = aVisu2MedEntity[anEntity];
146       for (int iGeomElem = 0; iGeomElem < iGeomElemEnd; iGeomElem++) {
147         int medId = getIdMedType(aGeomElemVector[iGeomElem]);
148         int aVtkType = med2vtk[medId].vtkType;
149         med_geometrie_element aMedType = med2vtk[medId].medType;
150         if(aMedEntity == MED_NOEUD){
151           med_geometrie_element typgeo = (med_geometrie_element)0;
152           med_connectivite typco = (med_connectivite)0;
153           med_int iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_COOR,MED_NOEUD,typgeo,typco);
154           if(iNumElemEnd > 0){
155             VISU::TMeshOnEntity& aMeshOnEntity = aMesh.myMeshOnEntityMap[anEntity];
156             aMeshOnEntity.myEntity = anEntity;
157             aMeshOnEntity.myMeshName = aMeshName;
158             aMeshOnEntity.myNbCells = iNumElemEnd;
159             aMeshOnEntity.myCellsSize = 2*iNumElemEnd;
160             aMesh.myNbPoints = iNumElemEnd;
161             if(MYDEBUG) 
162               MESSAGE("ImportInfo -\t anEntity = "<<anEntity<<
163                       "; medName = "<<med2vtk[medId].medName<<
164                       "; myNbCells = "<<aMeshOnEntity.myNbCells<<
165                       "; myCellsSize = "<<aMeshOnEntity.myCellsSize);
166             med_booleen iname_elem, inum_elem;
167             valarray<med_int> num_elem(iNumElemEnd), num_fam_elem(iNumElemEnd);
168             valarray<char> name_elem('\0',iNumElemEnd*MED_TAILLE_PNOM+1);
169             med_repere rep;
170             valarray<char> name_coord('\0',aMesh.myDim*MED_TAILLE_PNOM+1);
171             valarray<char> unit_coord('\0',aMesh.myDim*MED_TAILLE_PNOM+1);
172             valarray<med_float> coord(aMesh.myDim*iNumElemEnd);
173             ret = MEDnoeudsLire(fid,aMeshName,aMesh.myDim,&coord[0],MED_FULL_INTERLACE,&rep,
174                                 &name_coord[0],&unit_coord[0],&name_elem[0],&iname_elem,
175                                 &num_elem[0],&inum_elem,&num_fam_elem[0],iNumElemEnd);
176             if(ret < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDnoeudsLire(...)"));
177             for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++)
178               if(num_fam_elem[iNumElem] != 0){
179                 aFamily2EntityMap[num_fam_elem[iNumElem]] = anEntity;
180                 aFamilyNbCellsCounterMap[num_fam_elem[iNumElem]] += 1;
181                 aFamilyCellsSizeCounterMap[num_fam_elem[iNumElem]] += 2;
182               }
183           }
184         }
185         //Get number of connectivities
186         med_int iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_CONN,aMedEntity,aMedType,MED_NOD); 
187         if (iNumElemEnd > 0) {
188           VISU::TMeshOnEntity& aMeshOnEntity = aMesh.myMeshOnEntityMap[anEntity];
189           aMeshOnEntity.myEntity = anEntity;
190           aMeshOnEntity.myMeshName = aMeshName;
191           med_booleen iname_elem, inum_elem;
192           valarray<med_int> num_elem(iNumElemEnd), num_fam_elem(iNumElemEnd);
193           valarray<char> name_elem('\0',iNumElemEnd*MED_TAILLE_PNOM+1);
194           med_int aNbConnForElem = getNbMedConnect(aMedType,anEntity,aMesh.myDim);
195           aMeshOnEntity.myNbCells += iNumElemEnd;
196           aMeshOnEntity.myCellsSize += iNumElemEnd*(med2vtk[medId].vtkNbNodes + 1);
197           if(MYDEBUG) 
198             MESSAGE("ImportInfo -\t anEntity = "<<anEntity<<
199                     "; medName = "<<med2vtk[medId].medName<<
200                     "; aNbConnForElem = "<<aNbConnForElem<<
201                     "; myNbCells = "<<aMeshOnEntity.myNbCells<<
202                     "; myCellsSize = "<<aMeshOnEntity.myCellsSize);
203           valarray<med_int> conn(0,aNbConnForElem*iNumElemEnd);
204           ret = MEDelementsLire(fid,aMeshName,aMesh.myDim,&conn[0],MED_FULL_INTERLACE,
205                                 &name_elem[0],&iname_elem,&num_elem[0],&inum_elem,
206                                 &num_fam_elem[0],iNumElemEnd,aMedEntity,aMedType,MED_NOD);
207           if (ret < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDelementsLire(...)"));
208           for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++)
209             if(num_fam_elem[iNumElem] != 0){
210               aFamily2EntityMap[num_fam_elem[iNumElem]] = anEntity;
211               aFamilyNbCellsCounterMap[num_fam_elem[iNumElem]] += 1;
212               aFamilyCellsSizeCounterMap[num_fam_elem[iNumElem]] += med2vtk[medId].vtkNbNodes + 1;
213             } 
214         }
215       }
216     }
217     med_int aNbFamily = MEDnFam(fid,aMeshName,0,MED_FAMILLE);
218     if(MYDEBUG) MESSAGE("ImportInfo - aNbFamily = "<<aNbFamily);
219     for(int aFamInd = 1; aFamInd <= aNbFamily; aFamInd++){
220       med_int aNbAttrib = MEDnFam(fid,aMeshName,aFamInd,MED_ATTR);
221       valarray<med_int> anAttId(aNbAttrib), anAttVal(aNbAttrib);
222       valarray<char> anAttDesc('\0',aNbAttrib*MED_TAILLE_DESC+1);
223       med_int aNbGroup = MEDnFam(fid,aMeshName,aFamInd,MED_GROUPE);
224       if(0 && MYDEBUG) 
225         MESSAGE("ImportInfo - aFamInd = "<<aFamInd<<"; aNbAttrib = "<<aNbAttrib<<"; aNbGroup = "<<aNbGroup);
226       valarray<char> aGroupNames('\0',aNbGroup*MED_TAILLE_LNOM+1);
227       char aFamilyName[MED_TAILLE_NOM+1] = "";
228       med_int aFamilyNum;
229       ret = MEDfamInfo(fid,aMeshName,aFamInd,aFamilyName,&aFamilyNum,
230                        &anAttId[0],&anAttVal[0],&anAttDesc[0],&aNbAttrib,
231                        &aGroupNames[0],&aNbGroup);
232       if(ret < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDfamInfo"));
233       if(0 && MYDEBUG) 
234         MESSAGE("ImportInfo - aFamilyNum = "<<aFamilyNum<<"; aNbGroup = "<<aNbGroup);
235       if(aFamily2EntityMap.find(aFamilyNum) == aFamily2EntityMap.end()) {
236         if(MYDEBUG) MESSAGE("ImportInfo - a Family with name '"<<aFamilyName<<"' are empty !!!");
237         continue;
238       }
239       VISU::TEntity anEntity = aFamily2EntityMap[aFamilyNum];
240       VISU::TMeshOnEntity& aMeshOnEntity = aMesh.myMeshOnEntityMap[anEntity];
241       VISU::TFamily& aFamily = aMeshOnEntity.myFamilyMap[aFamilyName];
242       aFamily.myName = aFamilyName;
243       aFamily.myEntity = anEntity;
244       aFamily.myId = aFamilyNum;
245       aFamily.myNbCells = aFamilyNbCellsCounterMap[aFamilyNum];
246       aFamily.myCellsSize = aFamilyCellsSizeCounterMap[aFamilyNum];
247       if(MYDEBUG) 
248         MESSAGE("ImportInfo - aFamily.myEntity = "<<anEntity<<
249                 "; myName = '"<<aFamilyName<<
250                 "'; myId = "<<aFamilyNum<<
251                 "; myNbCells = "<<aFamily.myNbCells<<
252                 "; myCellsSize = "<<aFamily.myCellsSize);
253       VISU::TBindGroups& aBindGroups = aFamily.myGroups;
254       for(int iGroup = 0, iPos = 0; iGroup < aNbGroup; iGroup++, iPos += MED_TAILLE_LNOM){
255         char aGroupName[MED_TAILLE_LNOM+1];
256         strncpy(aGroupName,&aGroupNames[iPos],MED_TAILLE_LNOM);
257         aGroupName[MED_TAILLE_LNOM] = '\0';
258         if(MYDEBUG) MESSAGE("ImportInfo - aGroupName["<<iGroup<<"] = '"<<aGroupName<<"'");
259         aBindGroups.insert(aGroupName);
260       }
261     }
262     //Calculation of TMesh.TGroupMap
263     const VISU::TMeshOnEntityMap& aMeshOnEntityMap = aMesh.myMeshOnEntityMap;
264     if(aMeshOnEntityMap.empty()) continue;
265     VISU::TGroupMap& aGroupMap = aMesh.myGroupMap;
266     VISU::TMeshOnEntityMap::const_iterator aMeshOnEntityMapIter = aMeshOnEntityMap.begin();
267     for(; aMeshOnEntityMapIter != aMeshOnEntityMap.end(); aMeshOnEntityMapIter++){
268       const VISU::TMeshOnEntity& aMeshOnEntity = aMeshOnEntityMapIter->second;
269       const VISU::TFamilyMap& aFamilyMap = aMeshOnEntity.myFamilyMap;
270       if(aFamilyMap.empty()) continue;
271       VISU::TFamilyMap::const_iterator aFamilyMapIter = aFamilyMap.begin();
272       for(; aFamilyMapIter != aFamilyMap.end(); aFamilyMapIter++){
273         const VISU::TFamily& aFamily = aFamilyMapIter->second;
274         const VISU::TBindGroups& aBindGroups = aFamily.myGroups;
275         if(aBindGroups.empty()) continue;
276         VISU::TBindGroups::const_iterator aBindGroupsIter = aBindGroups.begin();
277         for(; aBindGroupsIter != aBindGroups.end(); aBindGroupsIter++){
278           const string& aGroupName = *aBindGroupsIter;
279           VISU::TGroup& aGroup = aGroupMap[aGroupName];
280           aGroup.myName = aGroupName;
281           aGroup.myMeshName = aMesh.myName;
282           VISU::TFamilyAndEntity aFamilyAndEntity(aFamily.myName,aFamily.myEntity);
283           aGroup.myFamilyAndEntitySet.insert(aFamilyAndEntity);
284           aGroup.myNbCells += aFamily.myNbCells;
285           aGroup.myCellsSize += aFamily.myCellsSize;
286         }
287       }
288     }
289     //Displaing information for the TMesh.TGroupMap
290     VISU::TGroupMap::const_iterator aGroupMapIter = aGroupMap.begin();
291     if(MYDEBUG) MESSAGE("ImportInfo - aGroupMap.size() = "<<aGroupMap.size());
292     for(; aGroupMapIter != aGroupMap.end(); aGroupMapIter++){
293       const VISU::TGroup& aCGroup = aGroupMapIter->second;
294       const string& aGroupName = aGroupMapIter->first;
295       if(MYDEBUG) MESSAGE("ImportInfo - aGroupName = '"<<aGroupName<<
296                 "'; myNbCells = "<<aCGroup.myNbCells<<
297                 "; myCellsSize = "<<aCGroup.myCellsSize);
298       const VISU::TGroup& aGroup = aGroupMapIter->second;
299       const VISU::TFamilyAndEntitySet& aFamilyAndEntitySet = aGroup.myFamilyAndEntitySet;
300       VISU::TFamilyAndEntitySet::const_iterator aFamilyAndEntitySetIter = aFamilyAndEntitySet.begin();
301       for(; aFamilyAndEntitySetIter != aFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
302         const string& aFamilyName = aFamilyAndEntitySetIter->first;
303         if(MYDEBUG) MESSAGE("ImportInfo - \t aFamilyName = '"<<aFamilyName<<"'");
304       }
305     }
306   }
307   //Reading information about fields
308   med_int iFieldEnd = MEDnChamp(fid,0);
309   if (iFieldEnd < 0) throw std::runtime_error(EXCEPTION("ImportChamps >> MEDnChamp(fid,0)"));
310   if(MYDEBUG) MESSAGE("ImportInfo - iFieldEnd = "<<iFieldEnd);
311   for(med_int iField = 1; iField <= iFieldEnd; iField++){
312     med_int ncomp = MEDnChamp(fid,iField);
313     if(ncomp < 0) throw std::runtime_error(EXCEPTION("ImportChamps >> MEDnChamp(fid,i)"));
314     valarray<char> aCompNames('\0',ncomp*MED_TAILLE_PNOM + 1);
315     valarray<char> aUnitNames('\0',ncomp*MED_TAILLE_PNOM + 1);
316     char name_field[MED_TAILLE_NOM + 1] = "";
317     med_type_champ type_field;
318     if(MEDchampInfo(fid,iField,name_field,&type_field,&aCompNames[0],&aUnitNames[0],ncomp) < 0)
319       throw std::runtime_error(EXCEPTION("ImportInfo >> MEDchampInfo(...)"));
320     //if(type_field != MED_REEL64) continue; //There is some problem in reading INTXX values
321     TVisu2MedEntity::const_iterator aVisu2MedEntityIter = aVisu2MedEntity.begin();
322     for(; aVisu2MedEntityIter != aVisu2MedEntity.end(); aVisu2MedEntityIter++) {
323       VISU::TEntity anEntity = aVisu2MedEntityIter->first;
324       int iGeomElemEnd;
325       med_geometrie_element* aGeomElemVector;
326       GetEntity2Geom(anEntity,aGeomElemVector,&iGeomElemEnd);
327       med_entite_maillage& aMedEntity = aVisu2MedEntity[anEntity];
328       for (int iGeomElem = 0; iGeomElem < iGeomElemEnd; iGeomElem++) {
329         med_geometrie_element& aGeom = aGeomElemVector[iGeomElem];
330         med_int iTimeStampEnd = MEDnPasdetemps(fid,name_field,aMedEntity,aGeom);
331         if(iTimeStampEnd < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDnPasdetemps(...)"));
332         if(iTimeStampEnd > 0) {
333           for(med_int iTimeStamp = 1; iTimeStamp <= iTimeStampEnd; iTimeStamp++) {
334             char aMeshName[MED_TAILLE_NOM+1] = "";
335             med_int ngauss = 0, numdt = 0, numo = 0;
336             char dt_unit[MED_TAILLE_PNOM+1] = "";
337             med_float dt = 0;
338             ret = MEDpasdetempsInfo(fid,name_field,aMedEntity,aGeom,iTimeStamp,
339                                     aMeshName,&ngauss,&numdt,dt_unit,&dt,&numo);
340             if(ret < 0) throw std::runtime_error(EXCEPTION("ImportInfo >> MEDpasdetempsInfo(...) < 0"));
341             if(myMeshMap.find(aMeshName) == myMeshMap.end())
342               throw std::runtime_error(EXCEPTION("ImportInfo >> MEDpasdetempsInfo(...)"));
343             VISU::TMesh &aMesh = myMeshMap[aMeshName];
344             VISU::TMeshOnEntity& aMeshOnEntity = aMesh.myMeshOnEntityMap[anEntity];
345             VISU::TFieldMap::iterator aFieldMapIter = aMeshOnEntity.myFieldMap.find(name_field);
346             //if(MYDEBUG && aFieldMapIter == aMeshOnEntity.myFieldMap.end()){
347             VISU::TField& aField = aMeshOnEntity.myFieldMap[name_field];
348             if(iTimeStamp == 1){
349               aField.myId = iField;
350               aField.myName = name_field;
351               aField.myEntity = anEntity;
352               aField.myMeshName = aMeshName;
353               aField.myNbComp = ncomp;
354               aField.myNbValField = iTimeStampEnd;
355               aField.myDataSize = aMeshOnEntity.myNbCells * ncomp;
356               aField.myCompNames.resize(ncomp);
357               aField.myUnitNames.resize(ncomp);
358               if(MYDEBUG)
359                 MESSAGE("ImportInfo - aField.myName = '"<<name_field<<
360                         "'; myMeshName = '"<<aMeshName<<
361                         "'; myEntity = "<<anEntity<<
362                         "; myNbComp = "<<ncomp<<
363                         "; myDataSize = "<<aField.myDataSize);
364               for(int iComp = 0, iPos = 0; iComp < ncomp; iComp++, iPos += MED_TAILLE_PNOM){
365                 char aCompName[MED_TAILLE_PNOM+1], aUnitName[MED_TAILLE_PNOM+1];
366                 strncpy(aCompName,&aCompNames[iPos],MED_TAILLE_PNOM);
367                 aCompName[MED_TAILLE_PNOM] = '\0';
368                 aField.myCompNames[iComp] = aCompName;
369                 strncpy(aUnitName,&aUnitNames[iPos],MED_TAILLE_PNOM);
370                 aUnitName[MED_TAILLE_PNOM] = '\0';
371                 aField.myUnitNames[iComp] = aUnitName;
372                 if(MYDEBUG)
373                   MESSAGE("ImportInfo - aCompName = '"<<aCompName<<"'; aUnitName = '"<<aUnitName<<"'");
374               }
375               
376             }
377             VISU::TField::TValForTime& aValForTime = aField.myValField[iTimeStamp];
378             aValForTime.myId = iTimeStamp;
379             aValForTime.myFieldName = aField.myName;
380             aValForTime.myEntity = aField.myEntity;
381             aValForTime.myMeshName = aField.myMeshName;
382             aValForTime.myNbComp = aField.myNbComp;
383             aValForTime.myTime = VISU::TField::TTime(dt,dt_unit);
384             if(MYDEBUG) 
385               MESSAGE("ImportInfo -\t aValForTime.myTime = "<<dt<<", "<<dt_unit);
386           }
387         }
388       }
389     }
390   }
391   return this; 
392 }
393
394 int VISU_MedConvertor::LoadMeshOnEntity(VISU::TMeshOnEntity& theMeshOnEntity, 
395                                         const string& theFamilyName)
396   throw (std::runtime_error&)
397 {
398   //Open the med file (it will be closed by call of destructor)
399   MedFile aMedFile(myFileInfo.absFilePath());
400   med_idt fid = aMedFile.GetFid();
401   //Main part of code
402   const string& aMeshName = theMeshOnEntity.myMeshName;
403   const VISU::TEntity& anEntity = theMeshOnEntity.myEntity;
404   VISU::TMesh& aMesh = myMeshMap[aMeshName];
405   int isPointsUpdated;
406   if(anEntity == VISU::NODE_ENTITY) 
407     isPointsUpdated = LoadPoints(fid,aMesh,theFamilyName);
408   else 
409     isPointsUpdated = LoadPoints(fid,aMesh);
410   int isCellsOnEntityUpdated = LoadCellsOnEntity(fid,theMeshOnEntity,theFamilyName);
411
412   return (isPointsUpdated || isCellsOnEntityUpdated);
413 }
414
415 int VISU_MedConvertor::LoadMeshOnGroup(VISU::TMesh& theMesh, 
416                                        const VISU::TFamilyAndEntitySet& theFamilyAndEntitySet)
417      throw (std::runtime_error&)
418 {
419   //Open the med file (it will be closed by call of destructor)
420   MedFile aMedFile(myFileInfo.absFilePath());
421   med_idt fid = aMedFile.GetFid();
422   //Main part of code
423   int isPointsUpdated = 0;
424   int isCellsOnEntityUpdated = 0;
425   VISU::TFamilyAndEntitySet::const_iterator aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
426   for(; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
427     const string& aFamilyName = aFamilyAndEntitySetIter->first;
428     const VISU::TEntity& anEntity = aFamilyAndEntitySetIter->second;
429     VISU::TMeshOnEntity& aMeshOnEntity = theMesh.myMeshOnEntityMap[anEntity];
430     if(anEntity == VISU::NODE_ENTITY){
431       isPointsUpdated += LoadPoints(fid,theMesh,aFamilyName);
432       isCellsOnEntityUpdated += LoadCellsOnEntity(fid,aMeshOnEntity);
433     }else{
434       isPointsUpdated += LoadPoints(fid,theMesh);
435       isCellsOnEntityUpdated += LoadCellsOnEntity(fid,aMeshOnEntity,aFamilyName);
436     }
437   }
438
439   return (isPointsUpdated || isCellsOnEntityUpdated);
440 }
441
442 int VISU_MedConvertor::LoadFieldOnMesh(VISU::TMesh& theMesh, 
443                                        VISU::TMeshOnEntity& theMeshOnEntity, 
444                                        VISU::TField& theField, 
445                                        VISU::TField::TValForTime& theValForTime)
446   throw (std::runtime_error&)
447 {
448   //Open the med file (it will be closed by call of destructor)
449   MedFile aMedFile(myFileInfo.absFilePath());
450   med_idt fid = aMedFile.GetFid();
451   //Main part of code
452   int isPointsUpdated = LoadPoints(fid,theMesh);
453   int isCellsOnEntityUpdated = LoadCellsOnEntity(fid,theMeshOnEntity);
454   int isFieldUpdated = LoadField(fid,theMeshOnEntity,theField,theValForTime);
455   
456   return (isPointsUpdated || isCellsOnEntityUpdated || isFieldUpdated);
457 }
458
459 int VISU_MedConvertor::LoadPoints(const med_idt& fid, VISU::TMesh& theMesh, const string& theFamilyName) 
460   throw (std::runtime_error&)
461 {
462   try{
463     //Check on existing family
464     VISU::TMeshOnEntity& aMeshOnEntity = theMesh.myMeshOnEntityMap[VISU::NODE_ENTITY];
465     aMeshOnEntity.myEntity = VISU::NODE_ENTITY;
466     aMeshOnEntity.myMeshName = theMesh.myName;
467     VISU::TFamily* pFamily = VISU::GetFamily(aMeshOnEntity,theFamilyName);
468     bool isFamilyPresent = (pFamily != NULL);
469     VISU::TFamily& aFamily = *pFamily;
470     //Check on loading already done
471     bool isPointsLoaded = !theMesh.myPointsCoord.empty();
472     if(isPointsLoaded) 
473       if(!isFamilyPresent) return 0;
474       else if(!aFamily.mySubMesh.empty()) return 0;
475     if(MYDEBUG) 
476       MESSAGE("LoadPoints - isPointsLoaded = "<<isPointsLoaded<<"; theFamilyName = '"<<theFamilyName<<"'");
477     //Main part of code
478     char aMeshName[MED_TAILLE_NOM+1] = "";
479     strcpy(aMeshName,theMesh.myName.c_str());
480     med_geometrie_element typgeo = (med_geometrie_element)0; //MED_POINT1
481     med_connectivite typco = (med_connectivite)0; //MED_NOD
482     med_int iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_COOR,MED_NOEUD,typgeo,typco);
483     if (iNumElemEnd <= 0) throw std::runtime_error(EXCEPTION("LoadPoints >> MEDnEntMaa(...)"));
484     if(MYDEBUG) MESSAGE("LoadPoints - iNumElemEnd = "<<iNumElemEnd);
485     med_repere rep;
486     med_booleen iname_elem, inum_elem;
487     valarray<med_int> num_elem(iNumElemEnd), num_fam_elem(iNumElemEnd);
488     valarray<char> name_elem('\0',iNumElemEnd*MED_TAILLE_PNOM+1);
489     valarray<char> name_coord('\0',theMesh.myDim*MED_TAILLE_PNOM+1);
490     valarray<char> unit_coord('\0',theMesh.myDim*MED_TAILLE_PNOM+1);
491     valarray<med_float> coord(theMesh.myDim*iNumElemEnd);
492     ret = MEDnoeudsLire(fid,aMeshName,theMesh.myDim,&coord[0],MED_FULL_INTERLACE,&rep,
493                         &name_coord[0],&unit_coord[0],&name_elem[0],&iname_elem,
494                         &num_elem[0],&inum_elem,&num_fam_elem[0],iNumElemEnd);
495     if(ret < 0) throw std::runtime_error(EXCEPTION("LoadPoints >> MEDnoeudsLire(...)"));
496     if(!isPointsLoaded){
497       VISU::TMesh::TPointsCoord& aPointsCoord = theMesh.myPointsCoord;
498       aPointsCoord.resize(iNumElemEnd*theMesh.myDim);
499       if(MYDEBUG) MESSAGE("LoadPoints - Filling coordinates of the mesh - inum_elem = "<<inum_elem);
500       inum_elem = MED_FAUX; // It is workaround 
501       if(inum_elem == MED_FAUX)
502         for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++) 
503           for(int iDim = 0, iNumElem2Dim = iNumElem*theMesh.myDim; iDim < theMesh.myDim; iDim++, iNumElem2Dim++)
504             aPointsCoord[iNumElem2Dim] = coord[iNumElem2Dim];
505       else
506         for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++)
507           for(int iDim = 0, iNumElem2Dim = iNumElem*theMesh.myDim; iDim < theMesh.myDim; iDim++, iNumElem2Dim++)
508             aPointsCoord[num_elem[iNumElem2Dim]] = coord[iNumElem2Dim];
509       if(MYDEBUG) MESSAGE("LoadPoints - Filling aMeshOnEntity with type NODE_ENTITY");
510       VISU::TMeshOnEntity::TConnForCellType& aConnForCellType = aMeshOnEntity.myCellsConn[VTK_VERTEX];
511       aConnForCellType.resize(iNumElemEnd);
512       for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++)
513         aConnForCellType[iNumElem] = VISU::TMeshOnEntity::TConnect(1,iNumElem);
514     }
515     if(isFamilyPresent && iNumElemEnd > 0){
516       if(MYDEBUG) MESSAGE("LoadPoints - Filling aFamily SubMesh");
517       VISU::TFamily::TSubMeshOnCellType& aSubMeshOnCellType = aFamily.mySubMesh[VTK_VERTEX];
518       for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++) 
519         if(num_fam_elem[iNumElem] == aFamily.myId)
520           aSubMeshOnCellType.insert(iNumElem);
521     }
522     return 1;
523   }catch(std::runtime_error& exc){
524     theMesh.myPointsCoord.clear();
525     throw std::runtime_error(exc.what());
526   }catch(...){
527     theMesh.myPointsCoord.clear();
528     throw std::runtime_error(EXCEPTION("Unknown exception !!!"));
529   }
530   return 0;
531 }
532
533 int VISU_MedConvertor::LoadCellsOnEntity(const med_idt& fid, VISU::TMeshOnEntity& theMeshOnEntity,
534                                          const string& theFamilyName)
535   throw (std::runtime_error&)
536 {
537   try{
538     //Check on existing family
539     VISU::TFamily* pFamily = VISU::GetFamily(theMeshOnEntity,theFamilyName);
540     bool isFamilyPresent = (pFamily != NULL);
541     VISU::TFamily& aFamily = *pFamily;
542     //Check on loading already done
543     bool isCellsLoaded = !theMeshOnEntity.myCellsConn.empty();
544     if(isCellsLoaded) 
545       if(!isFamilyPresent) return 0;
546       else if(!aFamily.mySubMesh.empty()) return 0;
547     if(MYDEBUG) {
548       MESSAGE("LoadCellsOnEntity - theFamilyName = '"<<theFamilyName<<"'");
549       MESSAGE("LoadCellsOnEntity - isCellsLoaded = "<<isCellsLoaded<<"; isFamilyPresent = "<<isFamilyPresent);
550     }
551     //Main part of code
552     int iGeomElemEnd;
553     med_geometrie_element* aGeomElemVector;
554     const VISU::TEntity& anEntity = theMeshOnEntity.myEntity;
555     GetEntity2Geom(anEntity,aGeomElemVector,&iGeomElemEnd);
556     const med_entite_maillage& aMedEntity = aVisu2MedEntity[anEntity];
557     char aMeshName[MED_TAILLE_NOM+1] = "";
558     strcpy(aMeshName,theMeshOnEntity.myMeshName.c_str());
559     if(MYDEBUG) 
560       MESSAGE("LoadCellsOnEntity - theMeshOnEntity.myEntity = "<<theMeshOnEntity.myEntity<<
561               "; iGeomElemEnd = "<<iGeomElemEnd<<"; theFamilyName = '"<<theFamilyName<<"'");
562     VISU::TMesh &aMesh = myMeshMap[theMeshOnEntity.myMeshName];
563     int aNbPoints = aMesh.myPointsCoord.size()/aMesh.myDim;
564     valarray<med_int>  num_node(aNbPoints);
565     med_booleen inum_node = 
566       med_booleen(MEDnumLire(fid,aMeshName,&num_node[0],aNbPoints,MED_NOEUD,med_geometrie_element(0)) >= 0);
567     if(MYDEBUG) MESSAGE("LoadCellsOnEntity - inum_node = "<<inum_node);
568     map<med_int,med_int> node_map;
569     if(inum_node) 
570       for(int i = 0; i < aNbPoints; i++) 
571         node_map[num_node[i]-1] = i;
572     for (int iGeomElem = 0; iGeomElem < iGeomElemEnd; iGeomElem++) {
573       int medId = getIdMedType(aGeomElemVector[iGeomElem]);
574       int nbMedNodes = med2vtk[medId].medNbNodes;
575       int nbVtkNodes = med2vtk[medId].vtkNbNodes;
576       int aVtkType = med2vtk[medId].vtkType;
577       med_geometrie_element aMedType = med2vtk[medId].medType;
578       med_int iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_CONN,aMedEntity,aMedType,MED_NOD);
579       if (iNumElemEnd > 0) {
580         med_booleen iname_elem, inum_elem;
581         valarray<med_int> num_elem(iNumElemEnd), num_fam_elem(iNumElemEnd);
582         valarray<char> name_elem('\0',iNumElemEnd*MED_TAILLE_PNOM+1);
583         med_int aNbConnForElem = getNbMedConnect(aMedType,anEntity,aMesh.myDim);
584         if(MYDEBUG) MESSAGE("LoadCellsOnEntity - medName = "<<med2vtk[medId].medName<<
585                             "; iNumElemEnd = "<<iNumElemEnd<<"; aNbConnForElem = "<<aNbConnForElem);
586         valarray<med_int> conn(aNbConnForElem*iNumElemEnd);
587         ret = MEDelementsLire(fid,aMeshName,aMesh.myDim,&conn[0],MED_FULL_INTERLACE,
588                               &name_elem[0],&iname_elem,&num_elem[0],&inum_elem,
589                               &num_fam_elem[0],iNumElemEnd,aMedEntity,aMedType,MED_NOD);
590         if (ret < 0) throw std::runtime_error(EXCEPTION("LoadCellsOnEntity >> MEDelementsLire(...)"));
591         if(!isCellsLoaded){
592           VISU::TMeshOnEntity::TConnForCellType& aConnForCellType = theMeshOnEntity.myCellsConn[aVtkType];
593           aConnForCellType.resize(iNumElemEnd);
594           valarray<med_int> aConnect(nbMedNodes);
595           for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++) {
596             VISU::TMeshOnEntity::TConnect& anArray = aConnForCellType[iNumElem];
597             anArray.resize(nbVtkNodes);
598             if(inum_node)
599               for (int k = 0, kj = iNumElem*aNbConnForElem; k < nbMedNodes; k++)
600                 aConnect[k] = node_map[conn[kj+k]-1];
601             else
602               for (int k = 0, kj = iNumElem*aNbConnForElem; k < nbMedNodes; k++)
603                 aConnect[k] = conn[kj+k] - 1;
604             switch(aMedType){
605             case MED_TETRA4 :
606             case MED_TETRA10 :
607               anArray[0] = aConnect[0];
608               anArray[1] = aConnect[1];
609               anArray[2] = aConnect[3];  
610               anArray[3] = aConnect[2];  
611               break;
612             case MED_PYRA5 :
613             case MED_PYRA13 :
614               anArray[0] = aConnect[0];
615               anArray[1] = aConnect[3];  
616               anArray[2] = aConnect[2];
617               anArray[3] = aConnect[1];  
618               anArray[4] = aConnect[4];
619               break;
620             default:
621               for (int k = 0; k < nbVtkNodes; k++) 
622                 anArray[k] = aConnect[k];
623             }
624             for (int k = 0; k < nbVtkNodes; k++) 
625               if(anArray[k] < 0 || aNbPoints <= anArray[k])
626                 throw std::runtime_error(EXCEPT("ImportCells >> aNbPoints(%1) <= anArray[%2][%3](%4) < 0").
627                                          arg(aNbPoints).arg(iNumElem).arg(k).arg(anArray[k]).latin1());
628           }
629         }
630         //Filling aFamily SubMesh
631         if(isFamilyPresent){
632           VISU::TFamily::TSubMeshOnCellType& aSubMeshOnCellType = aFamily.mySubMesh[aVtkType];
633           for (int iNumElem = 0; iNumElem < iNumElemEnd; iNumElem++) 
634             if(num_fam_elem[iNumElem] == aFamily.myId)
635               aSubMeshOnCellType.insert(iNumElem);
636         }
637       }
638     }
639     return 1;
640   }catch(std::runtime_error& exc){
641     theMeshOnEntity.myCellsConn.clear();
642     throw std::runtime_error(exc.what());
643   }catch(...){
644     theMeshOnEntity.myCellsConn.clear();
645     throw std::runtime_error(EXCEPTION("Unknown exception !!!"));
646   }
647   return 0;
648 }
649
650 int VISU_MedConvertor::LoadField(const med_idt& fid, const VISU::TMeshOnEntity& theMeshOnEntity,
651                                  VISU::TField& theField, VISU::TField::TValForTime& theValForTime)
652      throw (std::runtime_error&)
653 {
654   //Check on loading already done
655   if(!theValForTime.myValForCells.empty()) return 0;
656   //Main part of code
657   med_int ncomp = MEDnChamp(fid,theField.myId);
658   if(ncomp < 0) throw std::runtime_error(EXCEPTION("LoadField >> MEDnChamp(fid,i)"));
659   valarray<char> comp('\0',ncomp*MED_TAILLE_PNOM + 1);
660   valarray<char> unit('\0',ncomp*MED_TAILLE_PNOM + 1);
661   char aFieldName[MED_TAILLE_NOM + 1] = "";
662   char aMeshName[MED_TAILLE_NOM+1] = "";
663   strcpy(aMeshName,theValForTime.myMeshName.c_str());
664   med_type_champ type_field;
665   if(MEDchampInfo(fid,theField.myId,aFieldName,&type_field,&comp[0],&unit[0],ncomp) < 0)
666     throw std::runtime_error(EXCEPTION("LoadField >> MEDchampInfo(...)"));
667   int iGeomElemEnd;
668   med_geometrie_element* aGeomElemVector;
669   const VISU::TEntity& anEntity = theField.myEntity;
670   GetEntity2Geom(anEntity,aGeomElemVector,&iGeomElemEnd);
671   med_entite_maillage& aMedEntity = aVisu2MedEntity[anEntity];
672   const VISU::TMeshOnEntity::TCellsConn& aCellsConn = theMeshOnEntity.myCellsConn;
673   if(MYDEBUG){
674     MESSAGE("LoadField - aMeshName = '"<<aMeshName<<"' aFieldName = '"<<aFieldName<<"'; anEntity = "<<anEntity);
675     MESSAGE("LoadField - iGeomElemEnd = "<<iGeomElemEnd<<"; ncomp = "<<ncomp<<"; type_field = "<<type_field);
676   }
677   for (int iGeomElem = 0; iGeomElem < iGeomElemEnd; iGeomElem++) {
678     const med_geometrie_element& aGeom = aGeomElemVector[iGeomElem];
679     med_int iTimeStampEnd = MEDnPasdetemps(fid,aFieldName,aMedEntity,aGeom);
680     //Checking for accordance between the mesh and the field data
681     med_int iNumElemEnd = 0;
682     if(aMedEntity == MED_NOEUD)
683       iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_COOR,MED_NOEUD,med_geometrie_element(0),med_connectivite(0));
684     else
685       iNumElemEnd = MEDnEntMaa(fid,aMeshName,MED_CONN,aMedEntity,aGeom,MED_NOD);
686     int medId = getIdMedType(aGeomElemVector[iGeomElem]);
687     int aVtkType = med2vtk[medId].vtkType;
688     if(iTimeStampEnd <= 0){
689       if(iNumElemEnd > 0){
690         if(!theField.myIsTrimmed){
691           theField.myIsTrimmed = true;
692           theField.myDataSize -= iNumElemEnd*theField.myNbComp;
693         }
694         //throw std::runtime_error(EXCEPT("VISU_MedConvertor::LoadField - There is no the data "
695         //                              "for cells with type %1 of the mesh !!!").
696         //                       arg(med2vtk[medId].medName).latin1());
697       }
698     }else{
699       med_int ngauss = 0, numdt = 0, numo = 0;
700       char dt_unit[MED_TAILLE_PNOM+1] = "";
701       med_float dt = 0;
702       ret = MEDpasdetempsInfo(fid,aFieldName,aMedEntity,aGeom,theValForTime.myId,
703                               aMeshName,&ngauss,&numdt,dt_unit,&dt,&numo);
704       if(ret < 0) throw std::runtime_error(EXCEPTION("LoadField >> MEDpasdetempsInfo(...)"));
705       med_int nval = MEDnVal(fid,aFieldName,aMedEntity,aGeom,numdt,numo);
706       if(nval <= 0) throw std::runtime_error(EXCEPTION("LoadField >> MEDnVal(...) - nval <= 0"));
707       else{
708         //Checking for accordance between the mesh and the field data
709         if(iNumElemEnd <= 0) 
710           throw std::runtime_error(EXCEPTION("LoadField - There is no the geom. elem. on the mesh !!!"));
711         static int aMaxNbGaussPts = 52;  // the workaround for broken files
712         if(ngauss > aMaxNbGaussPts) ngauss = 1;
713         if(iNumElemEnd*ngauss != nval) 
714           throw std::runtime_error(EXCEPT("LoadField - Size of values (%1) and size of mesh (%2) not equal !!!").
715                                    arg(nval).arg(iNumElemEnd*ngauss).latin1());
716         VISU::TField::TValForCellsWithType &anArray = theValForTime.myValForCells[aVtkType];
717         int jEnd = theField.myNbComp*nval;
718         int kEnd = jEnd/ngauss;
719         anArray.resize(kEnd);
720         char pflname[MED_TAILLE_NOM + 1] = "";
721         switch(type_field){
722         case MED_REEL64 : {
723           valarray<med_float> valr(jEnd);
724           ret = MEDchampLire(fid,aMeshName,aFieldName,(unsigned char*)&valr[0],MED_FULL_INTERLACE,MED_ALL,
725                              pflname,aMedEntity,aGeom,numdt,numo);
726           if(ret < 0) throw std::runtime_error(EXCEPTION("LoadField >> MEDchampLire(...)"));
727           for (med_int k = 0, j = 0; k < kEnd; k++, j += ngauss){
728             for (med_int iGauss = 0; iGauss < ngauss; iGauss++)
729               anArray[k] = valr[j+iGauss];
730             anArray[k] /= ngauss;
731           }
732           break;
733         }
734         //case MED_INT64 : //valarray<long long> valr(jEnd);
735         case MED_INT32 : //valarray<long int> valr(jEnd);
736         case MED_INT : {
737           valarray<med_int> valr(jEnd);
738           ret = MEDchampLire(fid,aMeshName,aFieldName,(unsigned char*)&valr[0],MED_FULL_INTERLACE,MED_ALL,
739                              pflname,aMedEntity,aGeom,numdt,numo);
740           if(ret < 0) throw std::runtime_error(EXCEPTION("LoadField >> MEDchampLire(...)"));
741           for (med_int k = 0, j = 0; k < kEnd; k++, j += ngauss){
742             for (med_int iGauss = 0; iGauss < ngauss; iGauss++)
743               anArray[k] = valr[j+iGauss];
744             anArray[k] /= ngauss;
745           }
746           break;
747         }
748         default :
749           throw std::runtime_error(EXCEPTION("LoadField >> Value of med_type_champ for the field is wrong !!!"));
750         }
751         if(MYDEBUG) MESSAGE("LoadField - aGeom = "<<aGeom<<"; nval = "<<nval<<"; ngauss = "<<ngauss
752                             <<"; iTimeStampEnd = "<<iTimeStampEnd<<"; pflname = '"<<pflname<<"'");
753       }
754     }
755   }
756   return 1; 
757 }