1 // VISU OBJECT : interactive object for VISU entities implementation
3 // Copyright (C) 2003 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
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.
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.
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
20 // See http://www.opencascade.org/SALOME/ or email : webmaster.salome@opencascade.org
23 // File : VISU_MedConvertor.cxx
24 // Author : Alexey PETROV
28 #include "VISU_MedConvertor.hxx"
29 #include "VISU_Convertor.hxx"
30 #include "VISU_ConvertorUtils.hxx"
32 #include "MED_Factory.hxx"
33 #include "MED_Algorithm.hxx"
34 #include "MED_Utilities.hxx"
36 #include <vtkCellType.h>
38 #define _EDF_NODE_IDS_
45 static int MYDEBUG = 0;
47 static int MYDEBUG = 0;
56 int MEDGeom2NbNodes(MED::EGeometrieElement theMEDGeomType)
58 return theMEDGeomType % 100;
61 int MEDGeomToVTK(MED::EGeometrieElement theMEDGeomType)
63 switch(theMEDGeomType){
64 case ePOINT1: return VTK_VERTEX;
65 case eSEG2: return VTK_LINE;
66 case eSEG3: return VTK_LINE;
67 case eTRIA3: return VTK_TRIANGLE;
68 case eTRIA6: return VTK_TRIANGLE;
69 case eQUAD4: return VTK_QUAD;
70 case eQUAD8: return VTK_QUAD;
71 case eTETRA4: return VTK_TETRA;
72 case eTETRA10: return VTK_TETRA;
73 case eHEXA8: return VTK_HEXAHEDRON;
74 case eHEXA20: return VTK_HEXAHEDRON;
75 case ePENTA6: return VTK_WEDGE;
76 case ePENTA15: return VTK_WEDGE;
77 case ePYRA5: return VTK_PYRAMID;
78 case ePYRA13: return VTK_PYRAMID;
79 case ePOLYGONE: return VTK_POLYGON;
80 // case ePOLYEDRE: return VTK_POLYEDRE;
85 int VTKGeom2NbNodes(int theVTKGeomType)
87 switch(theVTKGeomType){
88 case VTK_VERTEX: return 1;
89 case VTK_LINE: return 2;
90 case VTK_TRIANGLE: return 3;
91 case VTK_QUAD: return 4;
92 case VTK_TETRA: return 4;
93 case VTK_HEXAHEDRON: return 8;
94 case VTK_WEDGE: return 6;
95 case VTK_PYRAMID: return 5;
100 MED::EGeometrieElement VTKGeomToMED(int theVTKGeomType)
102 switch(theVTKGeomType){
103 case VTK_VERTEX: return ePOINT1;
104 case VTK_LINE: return eSEG2;
105 case VTK_TRIANGLE: return eTRIA3;
106 case VTK_QUAD: return eQUAD4;
107 case VTK_TETRA: return eTETRA4;
108 case VTK_HEXAHEDRON: return eHEXA8;
109 case VTK_WEDGE: return ePENTA6;
110 case VTK_PYRAMID: return ePYRA5;
111 case VTK_POLYGON: return ePOLYGONE;
113 return EGeometrieElement(-1);
116 TEntity MEDEntityToVTK(MED::EEntiteMaillage theMEDEntity)
118 switch(theMEDEntity){
119 case eNOEUD: return NODE_ENTITY;
120 case eARETE: return EDGE_ENTITY;
121 case eFACE: return FACE_ENTITY;
122 case eMAILLE: return CELL_ENTITY;
127 MED::EEntiteMaillage VTKEntityToMED(TEntity theVTKEntity)
129 switch(theVTKEntity){
130 case NODE_ENTITY: return eNOEUD;
131 case EDGE_ENTITY: return eARETE;
132 case FACE_ENTITY: return eFACE;
133 case CELL_ENTITY: return eMAILLE;
135 return MED::EEntiteMaillage(-1);
141 VISU_Convertor* CreateConvertor(const string& theFileName)
143 return new VISU_MedConvertor(theFileName);
146 VISU_MedConvertor::VISU_MedConvertor(const string& theFileName) {
147 myFileInfo.setFile(QString(theFileName.c_str()));
148 myName = myFileInfo.baseName().latin1();
151 VISU_Convertor* VISU_MedConvertor::Build() {
152 PWrapper aMed = CrWrapper(myFileInfo.absFilePath().latin1());
153 TInt aNbMeshes = aMed->GetNbMeshes();
154 TMeshMap& aMeshMap = myMeshMap;
156 MSG(MYDEBUG,"VISU_MedConvertor::Build()");
157 INITMSG(MYDEBUG,"GetNbMeshes() = "<<aNbMeshes<<"\n");
159 for(TInt iMesh = 1; iMesh <= aNbMeshes; iMesh++)
161 PMeshInfo aMeshInfo = aMed->GetPMeshInfo(iMesh);
163 PNodeInfo aNodeInfo = aMed->GetPNodeInfo(aMeshInfo);
165 MED::TEntityInfo aEntityInfo = aMed->GetEntityInfo(aMeshInfo);
167 TElemGroup aElemGroup = GetElemsByEntity(aMed,aMeshInfo,aEntityInfo);
169 TFamilyGroup aFamilyGroup = GetFamilies(aMed,aMeshInfo);
171 TFamilyByEntity aFamilyByEntity = GetFamiliesByEntity(aMed,aNodeInfo,aElemGroup,aFamilyGroup);
173 TGroupInfo aGroupInfo = GetFamiliesByGroup(aFamilyGroup);
175 // creating TMesh structure and TMeshOnEntityMap
176 typedef map<TInt,TInt> TFamilyCounterMap;
177 TFamilyCounterMap aFamilyNbCellsCounterMap, aFamilyCellsSizeCounterMap;
178 TFamilyCounterMap aFamilyNbPolygonesCounterMap, aFamilyPolygonesSizeCounterMap;
180 TInt aDim = aMeshInfo->GetDim();
181 const string& aMeshName = aMeshInfo->GetName();
183 PMEDMesh aMesh = aMeshMap[aMeshName](new TMEDMesh());
185 aMesh->myName = aMeshName;
186 aMesh->myNbPoints = aNodeInfo->GetNbElem();
187 aMesh->myMeshInfo = aMeshInfo;
188 aMesh->myEntityInfo = aEntityInfo;
190 INITMSG(MYDEBUG,"aMeshName = '"<<aMeshName<<
191 "'; myNbPoints = "<<aMesh->myNbPoints<<
192 "; aDim = "<<aDim<<"\n");
194 BEGMSG(MYDEBUG,"aEntityInfo.size() = "<<aEntityInfo.size()<<"\n");
196 TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
197 MED::TEntityInfo::iterator anEntityIter = aEntityInfo.begin();
198 for(; anEntityIter != aEntityInfo.end(); anEntityIter++){
199 const EEntiteMaillage& aMEntity = anEntityIter->first;
200 const MED::TGeom& aTGeom = anEntityIter->second;
202 TEntity aVEntity = MEDEntityToVTK(aMEntity);
203 PMEDMeshOnEntity aMeshOnEntity = aMeshOnEntityMap[aVEntity](new TMEDMeshOnEntity());
204 aMeshOnEntity->myEntity = aVEntity;
205 aMeshOnEntity->myMeshName = aMeshName;
206 aMeshOnEntity->myGeom = aTGeom;
208 INITMSG(MYDEBUG,"aMEntity = "<<aMEntity<<"; aVEntity = "<<aVEntity<<"\n");
210 if(aMEntity == eNOEUD){
211 aMeshOnEntity->myNbCells = aMesh->myNbPoints;
212 aMeshOnEntity->myCellsSize = 2*aMesh->myNbPoints;
214 for(TInt iElem = 0; iElem < aMesh->myNbPoints; iElem++){
215 TInt aFamId = aNodeInfo->GetFamNum(iElem);
217 aFamilyNbCellsCounterMap[aFamId] += 1;
218 aFamilyCellsSizeCounterMap[aFamId] += 2;
222 INITMSG(MYDEBUG,"myNbCells = "<<aMeshOnEntity->myNbCells<<
223 "; myCellsSize = "<<aMeshOnEntity->myCellsSize<<"\n");;
226 MED::TGeom::const_iterator anTGeomIter = aTGeom.begin();
227 aMeshOnEntity->myNbCells = 0;
228 aMeshOnEntity->myCellsSize = 0;
229 for(; anTGeomIter != aTGeom.end(); anTGeomIter++){
230 const EGeometrieElement& aGeom = anTGeomIter->first;
235 PPolygoneInfo aPolygoneInfo = aMed->GetPPolygoneInfo(aMeshInfo,aMEntity,aGeom);
236 TInt aNbElem = aPolygoneInfo->GetNbElem();
237 TElemNum aConn = aPolygoneInfo->GetConnectivite();
238 TElemNum aIndex = aPolygoneInfo->GetIndex();
239 TInt aNbIndex = aIndex.size();
240 TInt aNbConn = aConn.size();
242 aMeshOnEntity->myNbCells += aNbElem;
244 for (int ii = 0; ii<aNbElem ; ii++){
245 int aNbConnii = aPolygoneInfo->GetNbConn(ii);
246 aMeshOnEntity->myCellsSize += aNbConnii;
248 INITMSG(MYDEBUG,"aGeom = "<<aGeom<<"; aNbElem = "<<aNbElem<<
249 "; myNbPolygones = "<<aNbElem<<
250 "; nbConn= "<<aNbConn<<"\n");
252 for(TInt iElem = 0; iElem < aNbElem; iElem++){
253 TInt aFamId = aPolygoneInfo->GetFamNum(iElem);
255 aFamilyNbCellsCounterMap[aFamId] += 1;
256 ADDMSG(MYDEBUG,"aFamId = "<<aFamId<<" ");
257 aFamilyCellsSizeCounterMap[aFamId] += aPolygoneInfo->GetNbConn(iElem) + 1;
260 ADDMSG(MYDEBUG,endl);
265 int aVNbNodes = VTKGeom2NbNodes(MEDGeomToVTK(aGeom));
266 PCellInfo aCellInfo = aMed->GetPCellInfo(aMeshInfo,aMEntity,aGeom);
267 TInt aNbElem = aCellInfo->GetNbElem();
268 aMeshOnEntity->myNbCells += aNbElem;
269 aMeshOnEntity->myCellsSize += aNbElem*(aVNbNodes+1);
270 INITMSG(MYDEBUG,"aGeom = "<<aGeom<<"; aNbElem = "<<aNbElem<<
271 "; myNbCells = "<<aMeshOnEntity->myNbCells<<
272 "; myCellsSize = "<<aMeshOnEntity->myCellsSize<<"\n");
274 for(TInt iElem = 0; iElem < aNbElem; iElem++){
275 TInt aFamId = aCellInfo->GetFamNum(iElem);
277 aFamilyNbCellsCounterMap[aFamId] += 1;
278 ADDMSG(MYDEBUG,"aFamId = "<<aFamId<<" ");
279 aFamilyCellsSizeCounterMap[aFamId] += aVNbNodes + 1;
282 ADDMSG(MYDEBUG,endl);
289 TFamilyByEntity::const_iterator aFamilyByEntityIter = aFamilyByEntity.begin();
290 BEGMSG(MYDEBUG,"TFamilyByEntity:\n");
291 for(; aFamilyByEntityIter != aFamilyByEntity.end(); aFamilyByEntityIter++){
292 const EEntiteMaillage& aMEntity = aFamilyByEntityIter->first;
293 const TFamilyGroup& aFamilyGroup = aFamilyByEntityIter->second;
295 TEntity aVEntity = MEDEntityToVTK(aMEntity);
296 VISU::PMEDMeshOnEntity aMeshOnEntity = aMesh->myMeshOnEntityMap[aVEntity];
297 VISU::TFamilyMap& aFamilyMap = aMeshOnEntity->myFamilyMap;
299 if(aFamilyGroup.empty())
302 INITMSG(MYDEBUG,"aMEntity = "<<aMEntity<<"; aVEntity = "<<aVEntity<<"\n");
303 TFamilyGroup::const_iterator aFamilyGroupIter = aFamilyGroup.begin();
304 for(; aFamilyGroupIter != aFamilyGroup.end(); aFamilyGroupIter++){
305 const PFamilyInfo& aFamilyInfo = *aFamilyGroupIter;
306 if (aFamilyInfo->GetId() == 0)
309 const std::string& aFamilyName = aFamilyInfo->GetName();
310 PMEDFamily aFamily = aFamilyMap[aFamilyName](new TMEDFamily());
312 aFamily->myId = aFamilyInfo->GetId();
313 aFamily->myName = aFamilyInfo->GetName();
314 aFamily->myEntity = aVEntity;
315 aFamily->myNbCells = aFamilyNbCellsCounterMap[aFamily->myId];
316 aFamily->myCellsSize = aFamilyCellsSizeCounterMap[aFamily->myId];
318 INITMSG(MYDEBUG,"aFamilyName = '"<<aFamily->myName<<
319 "'; myId = "<<aFamily->myId<<"; "<<
320 "; aNbAttr = "<<aFamilyInfo->GetNbAttr()<<
321 "; aNbGroup = "<<aFamilyInfo->GetNbGroup()<<
322 "; myEntity = "<<aFamily->myEntity<<
323 "; myNbCells = "<<aFamily->myNbCells<<
324 "; myCellsSize = "<<aFamily->myCellsSize<<"\n");
326 VISU::TBindGroups& aBindGroups = aFamily->myGroups;
327 const TInt aNbGroup = aFamilyInfo->GetNbGroup();
328 for(TInt i = 0; i < aNbGroup; i++){
329 const string& aGroupName = aFamilyInfo->GetGroupName(i);
330 aBindGroups.insert(aGroupName);
331 INITMSG(MYDEBUG,"aGroupName = '"<<aGroupName<<"'\n");
336 BEGMSG(MYDEBUG,"VISU::TGroup:\n");
338 VISU::TGroupMap& aGroupMap = aMesh->myGroupMap;
339 TGroupInfo::const_iterator aGroupInfoIter = aGroupInfo.begin();
340 for(;aGroupInfoIter != aGroupInfo.end(); aGroupInfoIter++){
341 const string& aGroupName = aGroupInfoIter->first;
342 const TFamilyGroup& aFamilyGroup = aGroupInfoIter->second;
343 PMEDGroup aGroup(new TMEDGroup());
344 aGroup->myName = aGroupName;
345 aGroup->myMeshName = aMesh->myName;
347 INITMSG(MYDEBUG,"aGroup->myName = '"<<aGroup->myName<<"'\n");
349 TFamilyGroup::const_iterator aFamilyIter = aFamilyGroup.begin();
350 for(; aFamilyIter != aFamilyGroup.end(); aFamilyIter++){
351 const PFamilyInfo& aFamilyInfo = *aFamilyIter;
352 const string& aFamilyName = aFamilyInfo->GetName();
354 TEntity aVEntity = TEntity(-1);
358 const TMeshOnEntityMap& aMeshOnEntityMap = aMesh->myMeshOnEntityMap;
359 TMeshOnEntityMap::const_iterator aMeshOnEntityIter = aMeshOnEntityMap.begin();
360 for(; aMeshOnEntityIter != aMeshOnEntityMap.end(); aMeshOnEntityIter++){
361 const PMeshOnEntity& aMeshOnEntity = aMeshOnEntityIter->second;
362 const TFamilyMap& aFamilyMap = aMeshOnEntity->myFamilyMap;
363 TFamilyMap::const_iterator aFamilyMapIter = aFamilyMap.begin();
364 for (; aFamilyMapIter != aFamilyMap.end(); aFamilyMapIter++){
365 const string& aName = aFamilyMapIter->first;
366 aFamily = aFamilyMapIter->second;
367 if(aName == aFamilyName){
368 aVEntity = aFamily->myEntity;
374 if(aFamily && aVEntity >= 0){
375 aGroup->myFamilyAndEntitySet.insert(TFamilyAndEntity(aFamilyName,aVEntity));
376 INITMSG(MYDEBUG,"aFamilyName = '"<<aFamilyName<<"'; '"<<aFamily->myName<<"'; aVEntity = "<<aVEntity<<"\n");
378 aGroup->myNbCells += aFamily->myNbCells;
379 aGroup->myCellsSize += aFamily->myCellsSize;
382 if(!aGroup->myFamilyAndEntitySet.empty() && aGroup->myNbCells > 0){
383 BEGMSG(MYDEBUG,"myNbCells = "<<aGroup->myNbCells<<
384 "; myCellsSize = "<<aGroup->myCellsSize<<"\n\n");
385 aGroupMap.insert(VISU::TGroupMap::value_type(aGroupName,aGroup));
389 TInt aNbFields = aMed->GetNbFields();
390 BEGMSG(MYDEBUG,"VISU::TField: NbFields="<<aNbFields<<"\n");
391 for(TInt iField = 1; iField <= aNbFields; iField++){
392 PFieldInfo aFieldInfo = aMed->GetPFieldInfo(aMeshInfo,iField);
393 TInt aNbComp = aFieldInfo->GetNbComp();
394 const string& aFieldName = aFieldInfo->GetName();
397 EEntiteMaillage aMEntity;
398 TInt aNbTimeStamps = aMed->GetNbTimeStamps(aFieldInfo,aEntityInfo,aMEntity,aTGeom);
399 TEntity aVEntity = MEDEntityToVTK(aMEntity);
400 VISU::PMeshOnEntity aMeshOnEntity = aMesh->myMeshOnEntityMap[aVEntity];
401 TFieldMap& aFieldMap = aMeshOnEntity->myFieldMap;
402 PMEDField aField = aFieldMap[aFieldName](new TMEDField());
403 aField->myId = iField;
404 aField->myNbComp = aNbComp;
405 aField->myEntity = aVEntity;
406 aField->myName = aFieldName;
407 aField->myMeshName = aMeshName;
408 aField->myDataSize = aMeshOnEntity->myNbCells * aNbComp;
409 aField->myCompNames.resize(aNbComp);
410 aField->myUnitNames.resize(aNbComp);
412 INITMSG(MYDEBUG,"myName = '"<<aField->myName<<
413 "'; myId = "<<aField->myId<<
414 "; myEntity = "<<aField->myEntity<<
415 "; myDataSize = "<<aField->myDataSize<<
416 "; myNbComp = "<<aField->myNbComp<<"\n");
418 for(TInt iComp = 0; iComp < aNbComp; iComp++){
419 aField->myCompNames[iComp] = aFieldInfo->GetCompName(iComp);
420 aField->myUnitNames[iComp] = aFieldInfo->GetUnitName(iComp);
423 for(TInt iTimeStamp = 1; iTimeStamp <= aNbTimeStamps; iTimeStamp++){
424 PTimeStampInfo aTimeStamp = aMed->GetPTimeStampInfo(aFieldInfo,
428 TFloat aDt = aTimeStamp->GetDt();
429 const string& anUnitDt = aTimeStamp->GetUnitDt();
430 PTimeStampVal aTimeStampVal = aMed->GetPTimeStampVal(aTimeStamp);
431 TValField& aValField = aField->myValField;
432 PMEDValForTime aValForTime = aValField[iTimeStamp](new TMEDValForTime());
433 aValForTime->myId = iTimeStamp;
434 aValForTime->myFieldName = aField->myName;
435 aValForTime->myEntity = aField->myEntity;
436 aValForTime->myMeshName = aField->myMeshName;
437 aValForTime->myNbComp = aField->myNbComp;
438 aValForTime->myTime = VISU::TTime(aDt,anUnitDt);
439 INITMSG(MYDEBUG,"aDt = "<<aDt<<", "<<anUnitDt<<"\n");
442 } catch (std::runtime_error& exc){
443 MSG(MYDEBUG,"Follow exception wqs occured in:\n"<<exc.what());
445 EXCEPTION(runtime_error,"Unknown exception !!!");
451 int VISU_MedConvertor::LoadMeshOnEntity(VISU::PMeshOnEntityImpl theMeshOnEntity,
452 const string& theFamilyName)
454 PWrapper aMed = CrWrapper(myFileInfo.absFilePath().latin1());
455 const string& aMeshName = theMeshOnEntity->myMeshName;
456 const VISU::TEntity& anEntity = theMeshOnEntity->myEntity;
457 PMeshImpl aMesh = myMeshMap[aMeshName];
459 if(anEntity == VISU::NODE_ENTITY)
460 isPointsUpdated = LoadPoints(aMed,aMesh,theFamilyName);
462 isPointsUpdated = LoadPoints(aMed,aMesh);
463 int isCellsOnEntityUpdated = LoadCellsOnEntity(aMed,aMesh,theMeshOnEntity,theFamilyName);
465 return (isPointsUpdated || isCellsOnEntityUpdated);
468 int VISU_MedConvertor::LoadMeshOnGroup(VISU::PMeshImpl theMesh,
469 const VISU::TFamilyAndEntitySet& theFamilyAndEntitySet)
471 PWrapper aMed = CrWrapper(myFileInfo.absFilePath().latin1());
472 int isPointsUpdated = 0, isCellsOnEntityUpdated = 0;
473 VISU::TFamilyAndEntitySet::const_iterator aFamilyAndEntitySetIter = theFamilyAndEntitySet.begin();
474 for(; aFamilyAndEntitySetIter != theFamilyAndEntitySet.end(); aFamilyAndEntitySetIter++){
475 const string& aFamilyName = aFamilyAndEntitySetIter->first;
476 const VISU::TEntity& anEntity = aFamilyAndEntitySetIter->second;
477 const VISU::PMEDMeshOnEntity aMeshOnEntity = theMesh->myMeshOnEntityMap[anEntity];
478 if(anEntity == VISU::NODE_ENTITY){
479 isPointsUpdated += LoadPoints(aMed,theMesh,aFamilyName);
480 isCellsOnEntityUpdated += LoadCellsOnEntity(aMed,theMesh,aMeshOnEntity);
482 isPointsUpdated += LoadPoints(aMed,theMesh);
483 isCellsOnEntityUpdated += LoadCellsOnEntity(aMed,theMesh,aMeshOnEntity,aFamilyName);
487 return (isPointsUpdated || isCellsOnEntityUpdated);
490 int VISU_MedConvertor::LoadFieldOnMesh(VISU::PMeshImpl theMesh,
491 VISU::PMeshOnEntityImpl theMeshOnEntity,
492 VISU::PFieldImpl theField,
493 VISU::PValForTimeImpl theValForTime)
495 PWrapper aMed = CrWrapper(myFileInfo.absFilePath().latin1());
496 int isPointsUpdated = LoadPoints(aMed,theMesh);
497 int isCellsOnEntityUpdated = LoadCellsOnEntity(aMed,theMesh,theMeshOnEntity);
498 int isFieldUpdated = LoadField(aMed,theMesh,theMeshOnEntity,theField,theValForTime);
500 return (isPointsUpdated || isCellsOnEntityUpdated || isFieldUpdated);
505 VISU_MedConvertor::LoadPoints(const MED::PWrapper& theMed,
506 VISU::PMEDMesh theMesh,
507 const string& theFamilyName)
510 //Check on existing family
511 VISU::PMEDMeshOnEntity aMeshOnEntity = theMesh->myMeshOnEntityMap[VISU::NODE_ENTITY];
512 aMeshOnEntity->myEntity = VISU::NODE_ENTITY;
513 aMeshOnEntity->myMeshName = theMesh->myName;
514 PFamilyImpl aFamily = GetFamily(aMeshOnEntity,theFamilyName);
515 //Check on loading already done
516 bool isPointsLoaded = !theMesh->myPointsCoord.empty();
520 else if(!aFamily->mySubMesh.empty())
523 INITMSG(MYDEBUG,"LoadPoints - isPointsLoaded = "<<isPointsLoaded<<"; theFamilyName = '"<<theFamilyName<<"'\n");
526 PNodeInfo aNodeInfo = theMed->GetPNodeInfo(theMesh->myMeshInfo);
527 TInt aNbElem = aNodeInfo->GetNbElem();
530 VISU::TMeshImpl::TPointsDim& aPointsDim = theMesh->myPointsDim;
531 aPointsDim.resize(theMesh->myDim);
532 for(int iDim = 0; iDim < theMesh->myDim; iDim++)
533 aPointsDim[iDim] = aNodeInfo->GetCoordName(iDim);
535 VISU::TMeshImpl::TPointsCoord& aPointsCoord = theMesh->myPointsCoord;
536 aPointsCoord.resize(aNbElem*theMesh->myDim);
537 for (int iElem = 0; iElem < aNbElem; iElem++)
538 for(int iDim = 0, iElem2Dim = iElem*theMesh->myDim; iDim < theMesh->myDim; iDim++, iElem2Dim++)
539 aPointsCoord[iElem2Dim] = aNodeInfo->GetNodeCoord(iElem,iDim);
541 VISU::TMeshOnEntityImpl::TConnForCellType& aConnForCellType = aMeshOnEntity->myCellsConn[VTK_VERTEX];
542 aConnForCellType.resize(aNbElem);
543 for (int iElem = 0; iElem < aNbElem; iElem++)
544 aConnForCellType[iElem] = VISU::TMeshOnEntityImpl::TConnect(1,iElem);
546 if(aFamily && aNbElem > 0){
547 VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aFamily->mySubMesh[VTK_VERTEX];
548 for (int iElem = 0; iElem < aNbElem; iElem++)
549 if(aNodeInfo->GetElemNum(iElem) == aFamily->myId)
550 aSubMeshOnCellType.insert(iElem);
553 }catch(std::runtime_error& exc){
554 theMesh->myPointsCoord.clear();
557 theMesh->myPointsCoord.clear();
558 EXCEPTION(runtime_error,"Unknown exception !!!");
565 VISU_MedConvertor::LoadCellsOnEntity(const MED::PWrapper& theMed,
566 VISU::PMEDMesh theMesh,
567 VISU::PMEDMeshOnEntity theMeshOnEntity,
568 const string& theFamilyName)
571 //Check on existing family
572 PFamilyImpl aFamily = GetFamily(theMeshOnEntity,theFamilyName);
573 //Check on loading already done
574 bool isCellsLoaded = !theMeshOnEntity->myCellsConn.empty();
578 else if(!aFamily->mySubMesh.empty())
581 INITMSG(MYDEBUG,"LoadCellsOnEntity - theFamilyName = '"<<theFamilyName<<"'\n");
582 BEGMSG(MYDEBUG,"LoadCellsOnEntity - isCellsLoaded = "<<isCellsLoaded<<"; isFamilyPresent = "<<bool(aFamily)<<endl);
584 const VISU::TEntity& aVEntity = theMeshOnEntity->myEntity;
585 const EEntiteMaillage& aMEntity = VTKEntityToMED(aVEntity);
587 const PMeshInfo& aMeshInfo = theMesh->myMeshInfo;
588 PNodeInfo aNodeInfo = theMed->GetPNodeInfo(aMeshInfo);
589 TInt aNbPoints = aNodeInfo->GetNbElem();
591 std::map<TInt,TInt> aNodeIdMap;
592 #ifdef _EDF_NODE_IDS_
593 EBooleen anIsNodeNum = eFAUX;
595 EBooleen anIsNodeNum = aNodeInfo->IsElemNum();
597 for(TInt i = 0; i < aNbPoints; i++){
598 aNodeIdMap[aNodeInfo->GetElemNum(i)-1] = i;
603 const MED::TGeom& aTGeom = theMeshOnEntity->myGeom;
604 MED::TGeom::const_iterator anTGeomIter = aTGeom.begin();
605 TMeshOnEntityImpl::TCellsConn& aCellsConn = theMeshOnEntity->myCellsConn;
607 for(; anTGeomIter != aTGeom.end(); anTGeomIter++){
608 const EGeometrieElement& aGeom = anTGeomIter->first;
609 int aVTKGeomType = MEDGeomToVTK(aGeom);
610 ADDMSG(MYDEBUG,"LoadCellsOnEntity aGeom="<<aGeom<<"\n");
614 PPolygoneInfo aPolygoneInfo = theMed->GetPPolygoneInfo(aMeshInfo,aMEntity,aGeom);
615 TInt aNbElem = aPolygoneInfo->GetNbElem();
618 VISU::TMeshOnEntityImpl::TConnForCellType& aConnForPolygoneType = aCellsConn[aVTKGeomType];
619 aConnForPolygoneType.resize(aNbElem);
621 int aMNbNodes = aPolygoneInfo->GetConnDim();
623 vector<TInt> aConnect(aMNbNodes);
624 vector<TInt> aIndex = aPolygoneInfo->GetIndex();
626 for (int iElem = 0; iElem < aNbElem; iElem++) {
627 VISU::TMeshOnEntityImpl::TConnect& anArray = aConnForPolygoneType[iElem];
628 int aNbConn = aPolygoneInfo->GetNbConn(iElem);
630 anArray.resize(aNbConn);
632 aConnect = aPolygoneInfo->GetConnectivite();
634 for (int i=0;i<aNbConn;i++){
635 anArray[i] = aConnect[aIndex[iElem]-1+i]-1;
640 VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aFamily->mySubMesh[aVTKGeomType];
641 for(int iElem = 0; iElem < aNbElem; iElem++)
642 if(aPolygoneInfo->GetFamNum(iElem) == aFamily->myId)
643 aSubMeshOnCellType.insert(iElem);
649 int aVNbNodes = VTKGeom2NbNodes(aVTKGeomType);
651 PCellInfo aCellInfo = theMed->GetPCellInfo(aMeshInfo,aMEntity,aGeom);
652 TInt aNbElem = aCellInfo->GetNbElem();
655 VISU::TMeshOnEntityImpl::TConnForCellType& aConnForCellType = aCellsConn[aVTKGeomType];
656 aConnForCellType.resize(aNbElem);
658 int aMNbNodes = MEDGeom2NbNodes(aGeom);
659 vector<TInt> aConnect(aMNbNodes);
661 for (int iElem = 0; iElem < aNbElem; iElem++) {
662 VISU::TMeshOnEntityImpl::TConnect& anArray = aConnForCellType[iElem];
663 anArray.resize(aVNbNodes);
666 for(int i = 0; i < aMNbNodes; i++){
667 aConnect[i] = aNodeIdMap[aCellInfo->GetConn(iElem,i)-1];
670 for(int i = 0; i < aMNbNodes; i++){
671 aConnect[i] = aCellInfo->GetConn(iElem,i)-1;
678 anArray[0] = aConnect[0];
679 anArray[1] = aConnect[1];
680 anArray[2] = aConnect[3];
681 anArray[3] = aConnect[2];
685 anArray[0] = aConnect[0];
686 anArray[1] = aConnect[3];
687 anArray[2] = aConnect[2];
688 anArray[3] = aConnect[1];
689 anArray[4] = aConnect[4];
692 for(int iNode = 0; iNode < aVNbNodes; iNode++)
693 anArray[iNode] = aConnect[iNode];
695 for(int iNode = 0; iNode < aVNbNodes; iNode++)
696 if(anArray[iNode] < 0 || aNbPoints <= anArray[iNode])
697 EXCEPTION(runtime_error,"ImportCells >> aNbPoints("<<aNbPoints<<") "<<
698 "<= anArray["<<iElem<<"]"<<
700 "("<<anArray[iNode]<<") < 0");
703 //Filling aFamily SubMesh
705 VISU::TFamilyImpl::TSubMeshOnCellType& aSubMeshOnCellType = aFamily->mySubMesh[aVTKGeomType];
706 for(int iElem = 0; iElem < aNbElem; iElem++)
707 if(aCellInfo->GetFamNum(iElem) == aFamily->myId)
708 aSubMeshOnCellType.insert(iElem);
714 }catch(std::runtime_error& exc){
715 theMeshOnEntity->myCellsConn.clear();
718 theMeshOnEntity->myCellsConn.clear();
719 EXCEPTION(runtime_error,"Unknown exception !!!");
726 VISU_MedConvertor::LoadField(const MED::PWrapper& theMed,
727 VISU::PMEDMesh theMesh,
728 VISU::PMEDMeshOnEntity theMeshOnEntity,
729 VISU::PMEDField theField,
730 VISU::PMEDValForTime theValForTime)
732 //Check on loading already done
733 if(!theValForTime->myValForCells.empty()) return 0;
736 const std::string& aMeshName = theMeshOnEntity->myMeshName;
737 const PMeshInfo& aMeshInfo = theMesh->myMeshInfo;
738 PFieldInfo aFieldInfo = theMed->GetPFieldInfo(aMeshInfo,theField->myId);
741 EEntiteMaillage aMEntity;
742 theMed->GetNbTimeStamps(aFieldInfo,theMesh->myEntityInfo,aMEntity,aTGeom);
744 PTimeStampInfo aTimeStampInfo = theMed->GetPTimeStampInfo(aFieldInfo,
747 theValForTime->myId);
748 TInt aNbGauss = aTimeStampInfo->GetNbGauss();
749 TInt aNbComp = theField->myNbComp;
751 PTimeStampVal aTimeStampVal = theMed->GetPTimeStampVal(aTimeStampInfo);
753 INITMSG(MYDEBUG,"LoadField - aMeshName = '"<<aMeshName<<
754 "'; aFieldName = '"<<aFieldInfo->GetName()<<
755 "'; aMEntity = "<<aMEntity<<
756 "; anId = "<<theValForTime->myId<<endl);
757 BEGMSG(MYDEBUG,"LoadField - aNbComp = "<<aNbComp<<
758 "; aNbGauss = "<<aNbGauss<<endl);
760 const MED::TGeom& anEntityTGeom = theMeshOnEntity->myGeom;
761 MED::TGeom::const_iterator aTGeomIter = anEntityTGeom.begin();
762 for(; aTGeomIter != anEntityTGeom.end(); aTGeomIter++){
763 const EGeometrieElement& aGeom = aTGeomIter->first;
764 const TInt& aNbElem = aTGeomIter->second;
766 INITMSG(MYDEBUG,"LoadField - aGeom = "<<aGeom<<"; aNbElem = '"<<aNbElem<<endl);
768 if(aTGeom.find(aGeom) == aTGeom.end()){
769 theField->myDataSize -= aNbElem*theField->myNbComp;
770 theField->myIsTrimmed = true;
772 int aVTKGeomType = MEDGeomToVTK(aGeom);
773 VISU::TValForTimeImpl::TValForCellsWithType& anArray = theValForTime->myValForCells[aVTKGeomType];
774 anArray.resize(aNbComp*aNbElem);
775 for(TInt iElem = 0, anId = 0; iElem < aNbElem; iElem++){
776 for(TInt iComp = 0; iComp < aNbComp; iComp++, anId++){
777 for(TInt iGauss = 0; iGauss < aNbGauss; iGauss++){
778 anArray[anId] += aTimeStampVal->GetVal(aGeom,iElem,iComp,iGauss);
780 anArray[anId] /= aNbGauss;