1 // Copyright (C) 2010-2016 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay
21 #include "MEDTimeReq.hxx"
22 #include "MEDUtilities.hxx"
24 #include "MEDFileFieldRepresentationTree.hxx"
25 #include "MEDCouplingFieldDiscretization.hxx"
26 #include "MEDCouplingFieldDouble.hxx"
27 #include "InterpKernelGaussCoords.hxx"
28 #include "MEDFileData.hxx"
29 #include "SauvReader.hxx"
31 #ifdef MEDREADER_USE_MPI
32 #include "ParaMEDFileMesh.hxx"
35 #include "vtkXMLUnstructuredGridWriter.h"//
37 #include "vtkUnstructuredGrid.h"
38 #include "vtkRectilinearGrid.h"
39 #include "vtkStructuredGrid.h"
40 #include "vtkUnsignedCharArray.h"
41 #include "vtkQuadratureSchemeDefinition.h"
42 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
43 #include "vtkInformationIntegerKey.h"
44 #include "vtkInformation.h"
45 #include "vtkDataArrayTemplate.h"
46 #include "vtkIdTypeArray.h"
47 #include "vtkDoubleArray.h"
48 #include "vtkIntArray.h"
49 #include "vtkCellArray.h"
50 #include "vtkPointData.h"
51 #include "vtkFieldData.h"
52 #include "vtkCellData.h"
54 #include "vtkMutableDirectedGraph.h"
56 using namespace ParaMEDMEM;
58 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
60 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
62 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
64 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
66 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
68 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
70 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
72 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
74 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
76 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
78 vtkIdTypeArray *ELGACmp::findOrCreate(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds, bool& isNew) const
80 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
89 return createNew(globs,locsReallyUsed,vtkd,ds);
93 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd) const
95 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
96 if(it==_loc_names.end())
98 std::size_t pos(std::distance(_loc_names.begin(),it));
99 vtkIdTypeArray *ret(_elgas[pos]);
100 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
101 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
103 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
105 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
109 vtkIdTypeArray *ELGACmp::createNew(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds) const
111 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
112 std::vector< std::vector<std::string> > locNames(_loc_names);
113 std::vector<vtkIdTypeArray *> elgas(_elgas);
114 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
116 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
117 if(it!=locNames.end())
118 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
119 locNames.push_back(locsReallyUsed);
120 vtkIdTypeArray *elga(vtkIdTypeArray::New());
121 elga->SetNumberOfComponents(1);
122 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
123 std::map<unsigned char,int> m;
124 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
126 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
127 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
128 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
129 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
130 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
131 // WARNING : these 2 lines are a workaround, due to users that write a ref element with dimension not equal to dimension of the geometric element.
132 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
133 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
134 double *shape(new double[nbPtsPerCell*nbGaussPt]);
135 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
136 calculator.initLocalInfo();
137 const std::vector<double>& wgths(loc.getGaussWeights());
138 for(int i=0;i<nbGaussPt;i++)
140 const double *pt0(calculator.getFunctionValues(i));
141 if(ct!=INTERP_KERNEL::NORM_HEXA27)
142 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
145 for(int j=0;j<27;j++)
146 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
149 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
150 m[vtkType]=nbGaussPt;
151 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
153 key->Set(elga->GetInformation(),def,vtkType);
154 key->Set(vtkd->GetInformation(),def,vtkType);
155 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
158 vtkIdType ncell(ds->GetNumberOfCells());
159 int *pt(new int[ncell]),offset(0);
160 for(vtkIdType cellId=0;cellId<ncell;cellId++)
162 vtkCell *cell(ds->GetCell(cellId));
163 int delta(m[cell->GetCellType()]);
167 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
168 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
169 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
170 std::string ossStr(oss.str());
171 elga->SetName(ossStr.c_str());
172 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
173 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
174 elgas.push_back(elga);
178 _defs.push_back(defs);
182 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
184 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
185 ds->GetCellData()->AddArray(*it);
190 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
192 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
193 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
194 (*it1).first->Delete();
199 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
203 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>& arr):ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
205 std::vector< std::vector<ParaMEDMEM::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
207 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
208 if(typs[0].size()!=1)
209 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field spitting by spatial discretization has failed ! A CRASH will occur soon !");
210 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
211 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
215 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
217 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>::operator=(other);
220 _ze_name=other._ze_name;
221 _ze_full_name.clear();
225 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
230 int MEDFileFieldRepresentationLeavesArrays::getId() const
235 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
237 return _ze_full_name;
240 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
242 return _ze_full_name.c_str();
245 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
247 vtkIdType refId(sil->AddChild(root,edge));
248 names.push_back(_ze_name);
250 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
252 sil->AddChild(refId,edge);
253 names.push_back(std::string());
257 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
259 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
260 _ze_full_name=oss3.str();
263 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
268 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
270 bool ret(_activated!=status);
275 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshStruct *mst, vtkDataSet *ds) const
277 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
278 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
279 tr->setNumberOfTS((operator->())->getNumberOfTS());
281 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
283 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
284 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
285 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
286 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
287 DataArray *crudeArr(0),*postProcessedArr(0);
289 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
291 crudeArr=f1tsPtrInt->getUndergroundDataArray();
293 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
294 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
295 f1ts->loadArraysIfNecessary();
296 MEDCouplingAutoRefCountObjectPtr<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
299 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
301 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
302 vtkFieldData *att(0);
307 att=ds->GetCellData();
312 att=ds->GetPointData();
317 att=ds->GetFieldData();
322 att=ds->GetFieldData();
326 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
330 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
331 vtkDoubleArray *vtkd(vtkDoubleArray::New());
332 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
333 for(int i=0;i<vd->getNumberOfComponents();i++)
334 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
335 if(postProcessedArr!=crudeArr)
337 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vd->accessToMemArray().setSpecificDeallocator(0);
341 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
343 std::string name(tr->buildName(f1ts->getName()));
344 vtkd->SetName(name.c_str());
347 if(discs[0]==ON_GAUSS_PT)
350 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
352 if(discs[0]==ON_GAUSS_NE)
354 vtkIdTypeArray *elno(vtkIdTypeArray::New());
355 elno->SetNumberOfComponents(1);
356 vtkIdType ncell(ds->GetNumberOfCells());
357 int *pt(new int[ncell]),offset(0);
358 std::set<int> cellTypes;
359 for(vtkIdType cellId=0;cellId<ncell;cellId++)
361 vtkCell *cell(ds->GetCell(cellId));
362 int delta(cell->GetNumberOfPoints());
363 cellTypes.insert(cell->GetCellType());
367 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
368 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
369 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
370 elno->SetName(nameElno.c_str());
371 ds->GetCellData()->AddArray(elno);
372 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
373 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
375 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
376 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
378 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
379 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
381 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
382 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
383 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
384 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
385 double *shape(new double[nbGaussPt*nbGaussPt]);
387 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));//GetLocsFromGeometricType
388 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
389 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
390 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
391 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
392 calculator.initLocalInfo();
393 for(int i=0;i<nbGaussPt;i++)
395 const double *pt0(calculator.getFunctionValues(i));
396 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
398 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
400 key->Set(elno->GetInformation(),def,*it);
401 key->Set(vtkd->GetInformation(),def,*it);
410 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
411 vtkIntArray *vtkd(vtkIntArray::New());
412 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
413 for(int i=0;i<vi->getNumberOfComponents();i++)
414 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
415 if(postProcessedArr!=crudeArr)
417 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vi->accessToMemArray().setSpecificDeallocator(0);
421 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
423 std::string name(tr->buildName(f1ts->getName()));
424 vtkd->SetName(name.c_str());
429 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
433 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
435 _elga_cmp.appendELGAIfAny(ds);
440 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
444 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS> >& arr,
445 const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
447 for(std::size_t i=0;i<arr.size();i++)
448 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
451 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
454 _cached_ds->Delete();
457 bool MEDFileFieldRepresentationLeaves::empty() const
459 const MEDFileFastCellSupportComparator *fcscp(_fsp);
460 return fcscp==0 || _arrays.empty();
463 void MEDFileFieldRepresentationLeaves::setId(int& id) const
465 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
469 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
471 return _arrays[0]->getMeshName();
474 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
476 return (int)_arrays.size();
479 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
481 return _arrays[0]->getNumberOfTS();
484 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
486 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
487 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
491 * \param [in] ms is the meshes pointer. It can be used only for information of geometric types. No special processing will be requested on ms.
493 void MEDFileFieldRepresentationLeaves::feedSIL(const ParaMEDMEM::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
495 vtkIdType root2(sil->AddChild(root,edge));
496 names.push_back(std::string("Arrs"));
497 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
498 (*it).feedSIL(sil,root2,edge,names);
500 vtkIdType root3(sil->AddChild(root,edge));
501 names.push_back(std::string("InfoOnGeoType"));
502 const ParaMEDMEM::MEDFileMesh *m(0);
504 m=ms->getMeshWithName(meshName);
505 const ParaMEDMEM::MEDFileFastCellSupportComparator *fsp(_fsp);
506 if(!fsp || fsp->getNumberOfTS()==0)
508 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
509 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
511 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
512 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
513 sil->AddChild(root3,edge);
514 names.push_back(cmStr);
518 bool MEDFileFieldRepresentationLeaves::containId(int id) const
520 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
521 if((*it).getId()==id)
526 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
528 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
529 if((*it).getZeName()==name)
537 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
539 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
540 status[(*it).getZeName()]=(*it).getStatus();
543 bool MEDFileFieldRepresentationLeaves::isActivated() const
545 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
546 if((*it).getStatus())
551 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
553 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
555 os << " - " << (*it0).getZeName() << " (";
556 if((*it0).getStatus())
560 os << ")" << std::endl;
564 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
566 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
567 (*it).setStatus(true);
570 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
572 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
573 if((*it).getId()==id)
575 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
578 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
581 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
582 std::vector<double> ret;
583 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
584 return tk.getTimeStepsRegardingPolicy(dtits,ret);
587 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
590 return _arrays[0]->getTimeSteps(ts);
594 return std::vector< std::pair<int,int> >();
598 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
600 std::ostringstream oss;
601 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
602 std::vector<double> ret1;
603 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
604 std::size_t sz(ret1.size());
605 for(std::size_t i=0;i<sz;i++)
607 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
610 std::string tmp(oss.str());
611 if(tmp.size()>200 && i!=sz-1)
621 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshes *meshes, vtkDataSet *ds) const
624 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
625 MEDCouplingAutoRefCountObjectPtr<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
626 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
627 if((*it).getStatus())
629 (*it).appendFields(tr,globs,mml,mst,ds);
630 (*it).appendELGAIfAny(ds);
634 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
636 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
637 DataArrayDouble *coordsMC(0);
638 DataArrayByte *typesMC(0);
639 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
640 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
641 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coordsSafe(coordsMC);
642 MEDCouplingAutoRefCountObjectPtr<DataArrayByte> typesSafe(typesMC);
643 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
645 int nbOfCells(typesSafe->getNbOfElems());
646 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
647 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
648 cellTypes->SetArray(reinterpret_cast<unsigned char *>(typesSafe->getPointer()),nbOfCells,0,VTK_DATA_ARRAY_FREE); typesSafe->accessToMemArray().setSpecificDeallocator(0);
649 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
650 cellLocations->SetVoidArray(cellLocationsSafe->getPointer(),nbOfCells,0,VTK_DATA_ARRAY_FREE); cellLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
651 vtkCellArray *cells(vtkCellArray::New());
652 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
653 cells2->SetVoidArray(cellsSafe->getPointer(),cellsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); cellsSafe->accessToMemArray().setSpecificDeallocator(0);
654 cells->SetCells(nbOfCells,cells2);
656 if(faceLocationsMC!=0 && facesMC!=0)
658 vtkIdTypeArray *faces(vtkIdTypeArray::New());
659 faces->SetVoidArray(facesSafe->getPointer(),facesSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); facesSafe->accessToMemArray().setSpecificDeallocator(0);
660 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
661 faceLocations->SetVoidArray(faceLocationsSafe->getPointer(),faceLocationsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); faceLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
662 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
663 faceLocations->Delete();
667 ret->SetCells(cellTypes,cellLocations,cells);
669 cellLocations->Delete();
671 vtkPoints *pts(vtkPoints::New());
672 vtkDoubleArray *pts2(vtkDoubleArray::New());
673 pts2->SetNumberOfComponents(3);
676 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
677 coordsSafe->accessToMemArray().setSpecificDeallocator(0);
680 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
689 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(ParaMEDMEM::MEDCMeshMultiLev *mm) const
691 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
693 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
694 vtkDoubleArray *vtkTmp(0);
695 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
696 std::size_t dim(arrs.size());
698 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
699 int sizePerAxe[3]={1,1,1};
700 sizePerAxe[0]=arrs[0]->getNbOfElems();
702 sizePerAxe[1]=arrs[1]->getNbOfElems();
704 sizePerAxe[2]=arrs[2]->getNbOfElems();
705 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
706 vtkTmp=vtkDoubleArray::New();
707 vtkTmp->SetNumberOfComponents(1);
709 vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
711 { vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[0]->accessToMemArray().setSpecificDeallocator(0); }
712 ret->SetXCoordinates(vtkTmp);
717 vtkTmp=vtkDoubleArray::New();
718 vtkTmp->SetNumberOfComponents(1);
720 vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
722 { vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[1]->accessToMemArray().setSpecificDeallocator(0); }
723 ret->SetYCoordinates(vtkTmp);
729 vtkTmp=vtkDoubleArray::New();
730 vtkTmp->SetNumberOfComponents(1);
732 vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
734 { vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[2]->accessToMemArray().setSpecificDeallocator(0); }
735 ret->SetZCoordinates(vtkTmp);
742 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(ParaMEDMEM::MEDCurveLinearMeshMultiLev *mm) const
744 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
745 int meshStr[3]={1,1,1};
746 DataArrayDouble *coords(0);
747 std::vector<int> nodeStrct;
749 mm->buildVTUArrays(coords,nodeStrct,isInternal);
750 std::size_t dim(nodeStrct.size());
752 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
753 meshStr[0]=nodeStrct[0];
755 meshStr[1]=nodeStrct[1];
757 meshStr[2]=nodeStrct[2];
758 vtkStructuredGrid *ret(vtkStructuredGrid::New());
759 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
760 vtkDoubleArray *da(vtkDoubleArray::New());
761 da->SetNumberOfComponents(3);
762 if(coords->getNumberOfComponents()==3)
765 da->SetArray(coords->getPointer(),coords->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);//VTK has not the ownership of double * because MEDLoader main struct has it !
767 { da->SetArray(coords->getPointer(),coords->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); coords->accessToMemArray().setSpecificDeallocator(0); }
771 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
772 da->SetArray(coords2->getPointer(),coords2->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);//let VTK deal with double *
773 coords2->accessToMemArray().setSpecificDeallocator(0);
776 vtkPoints *points=vtkPoints::New();
777 ret->SetPoints(points);
784 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDFileMeshes *meshes) const
786 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
788 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
789 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
790 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml2(mml->prepare());
791 MEDMeshMultiLev *ptMML2(mml2);
794 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
795 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
796 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
800 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
804 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
808 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
811 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
812 _cached_ds=ret->NewInstance();
813 _cached_ds->ShallowCopy(ret);
817 ret=_cached_ds->NewInstance();
818 ret->ShallowCopy(_cached_ds);
821 appendFields(tr,globs,mml,meshes,ret);
822 // The arrays links to mesh
823 DataArrayInt *famCells(0),*numCells(0);
824 bool noCpyFamCells(false),noCpyNumCells(false);
825 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
828 vtkIntArray *vtkTab(vtkIntArray::New());
829 vtkTab->SetNumberOfComponents(1);
830 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
832 vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
834 { vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famCells->accessToMemArray().setSpecificDeallocator(0); }
835 ret->GetCellData()->AddArray(vtkTab);
839 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
842 vtkIntArray *vtkTab(vtkIntArray::New());
843 vtkTab->SetNumberOfComponents(1);
844 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
846 vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
848 { vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numCells->accessToMemArray().setSpecificDeallocator(0); }
849 ret->GetCellData()->AddArray(vtkTab);
853 // The arrays links to mesh
854 DataArrayInt *famNodes(0),*numNodes(0);
855 bool noCpyFamNodes(false),noCpyNumNodes(false);
856 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
859 vtkIntArray *vtkTab(vtkIntArray::New());
860 vtkTab->SetNumberOfComponents(1);
861 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
863 vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
865 { vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famNodes->accessToMemArray().setSpecificDeallocator(0); }
866 ret->GetPointData()->AddArray(vtkTab);
870 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
873 vtkIntArray *vtkTab(vtkIntArray::New());
874 vtkTab->SetNumberOfComponents(1);
875 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
877 vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
879 { vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numNodes->accessToMemArray().setSpecificDeallocator(0); }
880 ret->GetPointData()->AddArray(vtkTab);
887 //////////////////////
889 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
893 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
896 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
897 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
898 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
899 ret+=(*it2).getNumberOfArrays();
903 void MEDFileFieldRepresentationTree::assignIds() const
906 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
907 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
908 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
912 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
914 std::size_t it0Cnt(0);
915 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
917 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
918 std::string tsName(oss.str());
919 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
921 std::string meshName((*it1)[0].getMeshName());
922 std::size_t it2Cnt(0);
923 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
925 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
926 std::string comSupStr(oss2.str());
927 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
933 void MEDFileFieldRepresentationTree::activateTheFirst() const
935 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
936 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
937 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
939 (*it2).activateAllArrays();
944 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
946 std::size_t it0Cnt(0);
947 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
949 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
950 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
952 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
953 std::vector<double> ts;
954 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
955 std::size_t nbOfTS(dtits.size());
956 std::ostringstream oss3; oss3 << nbOfTS;
957 names.push_back(oss3.str());
958 for(std::size_t i=0;i<nbOfTS;i++)
960 std::ostringstream oss4; oss4 << dtits[i].first;
961 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
962 names.push_back(oss4.str());
963 std::ostringstream oss5; oss5 << dtits[i].second;
964 vtkIdType ItId(sil->AddChild(DtId,edge));
965 names.push_back(oss5.str());
966 std::ostringstream oss6; oss6 << ts[i];
967 sil->AddChild(ItId,edge);
968 names.push_back(oss6.str());
971 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
972 std::string tsName(oss.str());
973 vtkIdType typeId0(sil->AddChild(root,edge));
974 names.push_back(tsName);
975 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
977 std::string meshName((*it1)[0].getMeshName());
978 vtkIdType typeId1(sil->AddChild(typeId0,edge));
979 names.push_back(meshName);
980 std::size_t it2Cnt(0);
981 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
983 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
984 std::string comSupStr(oss2.str());
985 vtkIdType typeId2(sil->AddChild(typeId1,edge));
986 names.push_back(comSupStr);
987 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
993 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
995 int dummy0(0),dummy1(0),dummy2(0);
996 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
997 std::string ret(leaf.getMeshName());
1000 for(;i<_ms->getNumberOfMeshes();i++)
1002 m=_ms->getMeshAtPos(i);
1003 if(m->getName()==ret)
1006 if(i==_ms->getNumberOfMeshes())
1007 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1008 vtkIdType typeId0(sil->AddChild(root,edge));
1009 names.push_back(m->getName());
1011 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1012 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1013 std::vector<std::string> grps(m->getGroupsNames());
1014 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1016 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1017 names.push_back(*it0);
1018 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1019 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1021 sil->AddChild(typeId2,edge);
1022 names.push_back((*it1).c_str());
1026 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1027 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1028 std::vector<std::string> fams(m->getFamiliesNames());
1029 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1031 sil->AddChild(typeId11,edge);
1032 int famId(m->getFamilyId((*it00).c_str()));
1033 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1034 names.push_back(oss.str());
1039 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1041 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1042 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1043 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1044 if((*it2).containId(id))
1045 return (*it2).getLeafArr(id);
1046 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1049 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1051 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1052 return elt.getZeName();
1055 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1057 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1058 return elt.getZeNameC();
1061 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1063 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1064 return elt.getStatus();
1067 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1070 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1071 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1072 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1073 if((*it2).containZeName(name,ret))
1075 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1076 throw INTERP_KERNEL::Exception(msg.str().c_str());
1079 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1081 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1082 bool ret(elt.setStatus(status));//to be implemented
1086 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1089 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1090 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1091 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1092 ret=std::max(ret,(*it2).getNumberOfTS());
1099 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1103 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1105 _ms=MEDFileMeshes::New(fileName);
1106 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1110 #ifdef MEDREADER_USE_MPI
1111 _ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1112 int nbMeshes(_ms->getNumberOfMeshes());
1113 for(int i=0;i<nbMeshes;i++)
1115 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1116 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1118 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp3(tmp2->zipCoords());
1120 _fields=MEDFileFields::LoadPartOf(fileName,false,_ms);//false is important to not read the values
1122 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1123 throw INTERP_KERNEL::Exception(oss.str().c_str());
1129 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::SauvReader> sr(ParaMEDMEM::SauvReader::New(fileName));
1130 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileData> mfd(sr->loadInMEDFileDS());
1131 _ms=mfd->getMeshes(); _ms->incrRef();
1132 int nbMeshes(_ms->getNumberOfMeshes());
1133 for(int i=0;i<nbMeshes;i++)
1135 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1136 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1138 tmp2->forceComputationOfParts();
1140 _fields=mfd->getFields();
1141 if((ParaMEDMEM::MEDFileFields *)_fields)
1144 if(!((ParaMEDMEM::MEDFileFields *)_fields))
1146 _fields=BuildFieldFromMeshes(_ms);
1150 AppendFieldFromMeshes(_ms,_fields);
1152 _ms->cartesianizeMe();
1153 _fields->removeFieldsWithoutAnyTimeStep();
1154 std::vector<std::string> meshNames(_ms->getMeshesNames());
1155 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> > fields_per_mesh(meshNames.size());
1156 for(std::size_t i=0;i<meshNames.size();i++)
1158 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1160 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1162 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1164 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1166 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1167 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1169 for(std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1171 if(!(*it)->presenceOfMultiDiscPerGeoType())
1172 allFMTSLeavesToDisplaySafe.push_back(*it);
1174 {// The case of some parts of field have more than one discretization per geo type.
1175 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1176 std::size_t it0Cnt(0);
1177 for(std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1179 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1180 (*it0)->setName(oss.str());
1181 allFMTSLeavesToDisplaySafe.push_back(*it0);
1188 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1189 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1191 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1193 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1194 // memory safety part
1195 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1196 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1198 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1199 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1201 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1202 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1205 // end of memory safety part
1206 // 1st : timesteps, 2nd : meshName, 3rd : common support
1207 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1208 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1210 std::vector< std::string > meshNamesLoc;
1211 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1212 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1214 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1215 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1216 if(it==meshNamesLoc.end())
1218 meshNamesLoc.push_back(meshName);
1219 splitByMeshName.resize(splitByMeshName.size()+1);
1220 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1223 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1225 _data_structure[i].resize(meshNamesLoc.size());
1226 for(std::size_t j=0;j<splitByMeshName.size();j++)
1228 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFastCellSupportComparator> > fsp;
1229 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1230 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1231 sbmn[k]=splitByMeshName[j][k];
1232 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1233 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1234 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1235 this->_data_structure[i][j].resize(commonSupSplit.size());
1236 for(std::size_t k=0;k<commonSupSplit.size();k++)
1238 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1239 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1241 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1242 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1245 for(std::size_t k=0;k<commonSupSplit.size();k++)
1246 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1249 this->removeEmptyLeaves();
1251 this->computeFullNameInLeaves();
1254 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1256 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1257 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1259 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1260 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1262 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1263 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1265 newSD1.push_back(*it2);
1267 newSD0.push_back(newSD1);
1270 newSD.push_back(newSD0);
1274 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1276 if(compInfos.size()!=1)
1278 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1281 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1283 return _data_structure[0][0][0].getMeshName();
1286 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1289 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1290 return leaf.getTimeSteps(tk);
1293 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1296 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1297 meshName=leaf.getMeshName();
1298 std::vector<double> ts(leaf.getTimeSteps(tk));
1299 std::size_t zeTimeId(0);
1302 std::vector<double> ts2(ts.size());
1303 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1304 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1305 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1308 if(zeTimeId==(int)ts.size())
1309 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1310 if(zeTimeId==(int)ts.size())
1311 {//OK the time requested does not fit time series given to ParaView. It is typically the case if more than one MEDReader instance are created or TimeInspector in real time mode.
1312 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1314 double valAttachedToPos(-std::numeric_limits<double>::max());
1315 for(std::size_t i=0;i<ts.size();i++)
1319 if(ts[i]>valAttachedToPos)
1322 valAttachedToPos=ts[i];
1327 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1328 valAttachedToPos=std::numeric_limits<double>::max();
1329 for(std::size_t i=0;i<ts.size();i++)
1331 if(ts[i]<valAttachedToPos)
1334 valAttachedToPos=ts[i];
1339 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1340 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1341 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1342 std::cerr << oss.str() << std::endl;
1346 tr=new MEDStdTimeReq((int)zeTimeId);
1348 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1349 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1354 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1356 int nbOfActivated(0);
1357 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1358 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1359 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1360 if((*it2).isActivated())
1362 if(nbOfActivated!=1)
1364 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1365 throw INTERP_KERNEL::Exception(oss.str().c_str());
1367 int i0(0),i1(0),i2(0);
1368 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1369 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1370 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1371 if((*it2).isActivated())
1373 lev0=i0; lev1=i1; lev2=i2;
1376 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1379 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1382 os << "#############################################" << std::endl;
1383 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1386 os << "TS" << i << std::endl;
1387 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1390 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1393 os << " " << (*it2).getMeshName() << std::endl;
1394 os << " Comp" << k << std::endl;
1395 (*it2).printMySelf(os);
1399 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1402 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1404 std::map<std::string,bool> ret;
1405 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1406 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1407 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1408 (*it2).dumpState(ret);
1412 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms, ParaMEDMEM::MEDFileFields *ret)
1415 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1416 for(int i=0;i<ms->getNumberOfMeshes();i++)
1418 MEDFileMesh *mm(ms->getMeshAtPos(i));
1419 std::vector<int> levs(mm->getNonEmptyLevels());
1420 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileField1TS> f1tsMultiLev(ParaMEDMEM::MEDFileField1TS::New());
1421 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1424 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1426 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1427 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1429 ParaMEDMEM::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1430 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1432 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1433 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1436 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1437 f1tsMultiLev->setFieldNoProfileSBT(f);
1442 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1443 if(levsExt.size()==levs.size()+1)
1445 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1446 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_NODES));
1448 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1449 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1450 arr->iota(); f->setArray(arr);
1451 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1452 f1tsMultiLev->setFieldNoProfileSBT(f);
1460 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1461 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1463 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1464 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1467 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1468 f1tsMultiLev->setFieldNoProfileSBT(f);
1471 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFieldMultiTS> fmtsMultiLev(ParaMEDMEM::MEDFileFieldMultiTS::New());
1472 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1473 ret->pushField(fmtsMultiLev);
1477 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const ParaMEDMEM::MEDFileFields *ret)
1479 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1481 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1482 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1483 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1485 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1486 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1487 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1491 ParaMEDMEM::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms)
1493 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFields> ret(ParaMEDMEM::MEDFileFields::New());
1494 AppendFieldFromMeshes(ms,ret);
1498 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1500 std::vector<std::string> ret;
1502 while(pos!=std::string::npos)
1504 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1505 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1507 pos=fullFieldName.find_first_not_of(sep,curPos);
1513 * Here the non regression tests.
1514 * const char inp0[]="";
1515 * const char exp0[]="";
1516 * const char inp1[]="field";
1517 * const char exp1[]="field";
1518 * const char inp2[]="_________";
1519 * const char exp2[]="_________";
1520 * const char inp3[]="field_p";
1521 * const char exp3[]="field_p";
1522 * const char inp4[]="field__p";
1523 * const char exp4[]="field_p";
1524 * const char inp5[]="field_p__";
1525 * const char exp5[]="field_p";
1526 * const char inp6[]="field_p_";
1527 * const char exp6[]="field_p";
1528 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1529 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1530 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1531 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1532 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1533 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1534 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1535 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1537 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1539 const char SEP('_');
1540 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1542 return fullFieldName;//should never happen
1546 return fullFieldName;
1550 std::string ret(v[0]);
1551 for(std::size_t i=1;i<v.size();i++)
1556 { ret+=SEP; ret+=v[i]; }
1562 return fullFieldName;
1568 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1572 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1577 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1579 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1581 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1587 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1588 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1589 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1591 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1593 std::size_t sz(ts.size());
1594 bool isInHumanRange(true);
1596 for(std::size_t i=0;i<sz;i++)
1599 if(ts[i]<=-1e299 || ts[i]>=1e299)
1600 isInHumanRange=false;
1603 return processedUsingPairOfIds(tsPairs);
1605 return processedUsingPairOfIds(tsPairs);
1606 _postprocessed_time=ts;
1607 return getPostProcessedTime();
1612 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1614 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1616 std::size_t sz(ts.size());
1617 std::vector<double> ts2(sz);
1619 for(std::size_t i=0;i<sz;i++)
1624 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1627 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1629 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1633 void TimeKeeper::printSelf(std::ostream& oss) const
1635 std::size_t sz(_activated_ts.size());
1636 for(std::size_t i=0;i<sz;i++)
1638 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1642 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1644 std::size_t sz(_activated_ts.size());
1645 std::vector<bool> ret(sz);
1646 for(std::size_t i=0;i<sz;i++)
1648 ret[i]=_activated_ts[i].first;
1653 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1655 std::size_t sz(tsPairs.size());
1656 std::set<int> s0,s1;
1657 for(std::size_t i=0;i<sz;i++)
1658 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1661 _postprocessed_time.resize(sz);
1662 for(std::size_t i=0;i<sz;i++)
1663 _postprocessed_time[i]=(double)tsPairs[i].first;
1664 return getPostProcessedTime();
1668 _postprocessed_time.resize(sz);
1669 for(std::size_t i=0;i<sz;i++)
1670 _postprocessed_time[i]=(double)tsPairs[i].second;
1671 return getPostProcessedTime();
1673 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1674 _postprocessed_time.resize(sz);
1675 for(std::size_t i=0;i<sz;i++)
1676 _postprocessed_time[i]=(double)i;
1677 return getPostProcessedTime();
1680 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1682 _activated_ts.resize(maxNumberOfTS);
1683 for(int i=0;i<maxNumberOfTS;i++)
1685 std::ostringstream oss; oss << "000" << i;
1686 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());