1 // Copyright (C) 2010-2014 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 #include "vtkXMLUnstructuredGridWriter.h"//
33 #include "vtkUnstructuredGrid.h"
34 #include "vtkRectilinearGrid.h"
35 #include "vtkStructuredGrid.h"
36 #include "vtkUnsignedCharArray.h"
37 #include "vtkQuadratureSchemeDefinition.h"
38 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
39 #include "vtkInformationIntegerKey.h"
40 #include "vtkInformation.h"
41 #include "vtkIdTypeArray.h"
42 #include "vtkDoubleArray.h"
43 #include "vtkIntArray.h"
44 #include "vtkCellArray.h"
45 #include "vtkPointData.h"
46 #include "vtkFieldData.h"
47 #include "vtkCellData.h"
49 #include "vtksys/stl/string"
50 #include "vtksys/ios/fstream"
51 #include "vtksys/stl/algorithm"
52 #include "vtkMutableDirectedGraph.h"
54 using namespace ParaMEDMEM;
56 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
58 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
60 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
62 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
64 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
66 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
68 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
70 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
72 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
74 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
76 vtkIdTypeArray *ELGACmp::findOrCreate(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds, bool& isNew) const
78 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
87 return createNew(globs,locsReallyUsed,vtkd,ds);
91 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd) const
93 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
94 if(it==_loc_names.end())
96 std::size_t pos(std::distance(_loc_names.begin(),it));
97 vtkIdTypeArray *ret(_elgas[pos]);
98 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
99 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
101 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
103 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
107 vtkIdTypeArray *ELGACmp::createNew(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds) const
109 static const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
110 std::vector< std::vector<std::string> > locNames(_loc_names);
111 std::vector<vtkIdTypeArray *> elgas(_elgas);
112 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
114 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
115 if(it!=locNames.end())
116 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
117 locNames.push_back(locsReallyUsed);
118 vtkIdTypeArray *elga(vtkIdTypeArray::New());
119 elga->SetNumberOfComponents(1);
120 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
121 std::map<unsigned char,int> m;
122 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
124 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
125 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
126 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
127 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
128 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
129 // 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.
130 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
131 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
132 double *shape(new double[nbPtsPerCell*nbGaussPt]);
133 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
134 calculator.initLocalInfo();
135 const std::vector<double>& wgths(loc.getGaussWeights());
136 for(int i=0;i<nbGaussPt;i++)
138 const double *pt0(calculator.getFunctionValues(i));
139 if(ct!=INTERP_KERNEL::NORM_HEXA27)
140 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
143 for(int j=0;j<27;j++)
144 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
147 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
148 m[vtkType]=nbGaussPt;
149 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
151 key->Set(elga->GetInformation(),def,vtkType);
152 key->Set(vtkd->GetInformation(),def,vtkType);
153 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
156 vtkIdType ncell(ds->GetNumberOfCells());
157 int *pt(new int[ncell]),offset(0);
158 for(vtkIdType cellId=0;cellId<ncell;cellId++)
160 vtkCell *cell(ds->GetCell(cellId));
161 int delta(m[cell->GetCellType()]);
165 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
166 elga->SetArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
167 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
168 std::string ossStr(oss.str());
169 elga->SetName(ossStr.c_str());
170 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
171 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
172 elgas.push_back(elga);
176 _defs.push_back(defs);
179 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
181 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
182 ds->GetCellData()->AddArray(*it);
187 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
189 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
190 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
191 (*it1).first->Delete();
196 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
200 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>& arr):ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
202 std::vector< std::vector<ParaMEDMEM::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
204 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
205 if(typs[0].size()!=1)
206 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 !");
207 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
208 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
212 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
214 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>::operator=(other);
217 _ze_name=other._ze_name;
218 _ze_full_name.clear();
222 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
227 int MEDFileFieldRepresentationLeavesArrays::getId() const
232 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
234 return _ze_full_name;
237 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
239 return _ze_full_name.c_str();
242 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
244 vtkIdType refId(sil->AddChild(root,edge));
245 names.push_back(_ze_name);
247 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
249 sil->AddChild(refId,edge);
250 names.push_back(std::string());
254 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
256 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
257 _ze_full_name=oss3.str();
260 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
265 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
267 bool ret(_activated!=status);
272 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshStruct *mst, vtkDataSet *ds) const
274 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
275 static const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
276 tr->setNumberOfTS((operator->())->getNumberOfTS());
278 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
280 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
281 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
282 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
283 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
284 DataArray *crudeArr(0),*postProcessedArr(0);
286 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
288 crudeArr=f1tsPtrInt->getUndergroundDataArray();
290 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
291 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
292 f1ts->loadArraysIfNecessary();
293 MEDCouplingAutoRefCountObjectPtr<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
296 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
298 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
299 vtkFieldData *att(0);
304 att=ds->GetCellData();
309 att=ds->GetPointData();
314 att=ds->GetFieldData();
319 att=ds->GetFieldData();
323 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
327 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
328 vtkDoubleArray *vtkd(vtkDoubleArray::New());
329 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
330 for(int i=0;i<vd->getNumberOfComponents();i++)
331 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
332 if(postProcessedArr!=crudeArr)
334 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vd->accessToMemArray().setSpecificDeallocator(0);
338 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
340 std::string name(tr->buildName(f1ts->getName()));
341 vtkd->SetName(name.c_str());
344 if(discs[0]==ON_GAUSS_PT)
347 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
349 if(discs[0]==ON_GAUSS_NE)
351 vtkIdTypeArray *elno(vtkIdTypeArray::New());
352 elno->SetNumberOfComponents(1);
353 vtkIdType ncell(ds->GetNumberOfCells());
354 int *pt(new int[ncell]),offset(0);
355 std::set<int> cellTypes;
356 for(vtkIdType cellId=0;cellId<ncell;cellId++)
358 vtkCell *cell(ds->GetCell(cellId));
359 int delta(cell->GetNumberOfPoints());
360 cellTypes.insert(cell->GetCellType());
364 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
365 elno->SetArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
366 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
367 elno->SetName(nameElno.c_str());
368 ds->GetCellData()->AddArray(elno);
369 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
370 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
372 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
373 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
375 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
376 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
378 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
379 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
380 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
381 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
382 double *shape(new double[nbGaussPt*nbGaussPt]);
384 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(ct,dummy));
385 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
386 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
387 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
388 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
389 calculator.initLocalInfo();
390 for(int i=0;i<nbGaussPt;i++)
392 const double *pt0(calculator.getFunctionValues(i));
393 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
395 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
397 key->Set(elno->GetInformation(),def,*it);
398 key->Set(vtkd->GetInformation(),def,*it);
407 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
408 vtkIntArray *vtkd(vtkIntArray::New());
409 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
410 for(int i=0;i<vi->getNumberOfComponents();i++)
411 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
412 if(postProcessedArr!=crudeArr)
414 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vi->accessToMemArray().setSpecificDeallocator(0);
418 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
420 std::string name(tr->buildName(f1ts->getName()));
421 vtkd->SetName(name.c_str());
426 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
430 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
432 _elga_cmp.appendELGAIfAny(ds);
437 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
441 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS> >& arr,
442 const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
444 for(std::size_t i=0;i<arr.size();i++)
445 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
448 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
451 _cached_ds->Delete();
454 bool MEDFileFieldRepresentationLeaves::empty() const
456 const MEDFileFastCellSupportComparator *fcscp(_fsp);
457 return fcscp==0 || _arrays.empty();
460 void MEDFileFieldRepresentationLeaves::setId(int& id) const
462 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
466 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
468 return _arrays[0]->getMeshName();
471 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
473 return (int)_arrays.size();
476 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
478 return _arrays[0]->getNumberOfTS();
481 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
483 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
484 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
488 * \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.
490 void MEDFileFieldRepresentationLeaves::feedSIL(const ParaMEDMEM::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
492 vtkIdType root2(sil->AddChild(root,edge));
493 names.push_back(std::string("Arrs"));
494 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
495 (*it).feedSIL(sil,root2,edge,names);
497 vtkIdType root3(sil->AddChild(root,edge));
498 names.push_back(std::string("InfoOnGeoType"));
499 const ParaMEDMEM::MEDFileMesh *m(0);
501 m=ms->getMeshWithName(meshName);
502 const ParaMEDMEM::MEDFileFastCellSupportComparator *fsp(_fsp);
503 if(!fsp || fsp->getNumberOfTS()==0)
505 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
506 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
508 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
509 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
510 sil->AddChild(root3,edge);
511 names.push_back(cmStr);
515 bool MEDFileFieldRepresentationLeaves::containId(int id) const
517 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
518 if((*it).getId()==id)
523 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
525 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
526 if((*it).getZeName()==name)
534 bool MEDFileFieldRepresentationLeaves::isActivated() const
536 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
537 if((*it).getStatus())
542 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
544 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
546 os << " - " << (*it0).getZeName() << " (";
547 if((*it0).getStatus())
551 os << ")" << std::endl;
555 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
557 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
558 (*it).setStatus(true);
561 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
563 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
564 if((*it).getId()==id)
566 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
569 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
572 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
573 std::vector<double> ret;
574 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
575 return tk.getTimeStepsRegardingPolicy(dtits,ret);
578 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
581 return _arrays[0]->getTimeSteps(ts);
585 return std::vector< std::pair<int,int> >();
589 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
591 std::ostringstream oss;
592 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
593 std::vector<double> ret1;
594 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
595 std::size_t sz(ret1.size());
596 for(std::size_t i=0;i<sz;i++)
598 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
601 std::string tmp(oss.str());
602 if(tmp.size()>200 && i!=sz-1)
612 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshes *meshes, vtkDataSet *ds) const
615 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
616 MEDCouplingAutoRefCountObjectPtr<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
617 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
618 if((*it).getStatus())
620 (*it).appendFields(tr,globs,mml,mst,ds);
621 (*it).appendELGAIfAny(ds);
625 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
627 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
628 DataArrayDouble *coordsMC(0);
629 DataArrayByte *typesMC(0);
630 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
631 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
632 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coordsSafe(coordsMC);
633 MEDCouplingAutoRefCountObjectPtr<DataArrayByte> typesSafe(typesMC);
634 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
636 int nbOfCells(typesSafe->getNbOfElems());
637 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
638 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
639 cellTypes->SetArray(reinterpret_cast<unsigned char *>(typesSafe->getPointer()),nbOfCells,0,VTK_DATA_ARRAY_FREE); typesSafe->accessToMemArray().setSpecificDeallocator(0);
640 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
641 cellLocations->SetArray(cellLocationsSafe->getPointer(),nbOfCells,0,VTK_DATA_ARRAY_FREE); cellLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
642 vtkCellArray *cells(vtkCellArray::New());
643 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
644 cells2->SetArray(cellsSafe->getPointer(),cellsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); cellsSafe->accessToMemArray().setSpecificDeallocator(0);
645 cells->SetCells(nbOfCells,cells2);
647 if(faceLocationsMC!=0 && facesMC!=0)
649 vtkIdTypeArray *faces(vtkIdTypeArray::New());
650 faces->SetArray(facesSafe->getPointer(),facesSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); facesSafe->accessToMemArray().setSpecificDeallocator(0);
651 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
652 faceLocations->SetArray(faceLocationsSafe->getPointer(),faceLocationsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); faceLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
653 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
654 faceLocations->Delete();
658 ret->SetCells(cellTypes,cellLocations,cells);
660 cellLocations->Delete();
662 vtkPoints *pts(vtkPoints::New());
663 vtkDoubleArray *pts2(vtkDoubleArray::New());
664 pts2->SetNumberOfComponents(3);
667 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
668 coordsSafe->accessToMemArray().setSpecificDeallocator(0);
671 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
680 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(ParaMEDMEM::MEDCMeshMultiLev *mm) const
682 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
684 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
685 vtkDoubleArray *vtkTmp(0);
686 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
687 std::size_t dim(arrs.size());
689 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
690 int sizePerAxe[3]={1,1,1};
691 sizePerAxe[0]=arrs[0]->getNbOfElems();
693 sizePerAxe[1]=arrs[1]->getNbOfElems();
695 sizePerAxe[2]=arrs[2]->getNbOfElems();
696 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
697 vtkTmp=vtkDoubleArray::New();
698 vtkTmp->SetNumberOfComponents(1);
700 vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
702 { vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[0]->accessToMemArray().setSpecificDeallocator(0); }
703 ret->SetXCoordinates(vtkTmp);
708 vtkTmp=vtkDoubleArray::New();
709 vtkTmp->SetNumberOfComponents(1);
711 vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
713 { vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[1]->accessToMemArray().setSpecificDeallocator(0); }
714 ret->SetYCoordinates(vtkTmp);
720 vtkTmp=vtkDoubleArray::New();
721 vtkTmp->SetNumberOfComponents(1);
723 vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
725 { vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[2]->accessToMemArray().setSpecificDeallocator(0); }
726 ret->SetZCoordinates(vtkTmp);
733 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(ParaMEDMEM::MEDCurveLinearMeshMultiLev *mm) const
735 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
736 int meshStr[3]={1,1,1};
737 DataArrayDouble *coords(0);
738 std::vector<int> nodeStrct;
740 mm->buildVTUArrays(coords,nodeStrct,isInternal);
741 std::size_t dim(nodeStrct.size());
743 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
744 meshStr[0]=nodeStrct[0];
746 meshStr[1]=nodeStrct[1];
748 meshStr[2]=nodeStrct[2];
749 vtkStructuredGrid *ret(vtkStructuredGrid::New());
750 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
751 vtkDoubleArray *da(vtkDoubleArray::New());
752 da->SetNumberOfComponents(3);
753 if(coords->getNumberOfComponents()==3)
756 da->SetArray(coords->getPointer(),coords->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);//VTK has not the ownership of double * because MEDLoader main struct has it !
758 { da->SetArray(coords->getPointer(),coords->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); coords->accessToMemArray().setSpecificDeallocator(0); }
762 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
763 da->SetArray(coords2->getPointer(),coords2->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);//let VTK deal with double *
764 coords2->accessToMemArray().setSpecificDeallocator(0);
767 vtkPoints *points=vtkPoints::New();
768 ret->SetPoints(points);
775 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDFileMeshes *meshes) const
777 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
779 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
780 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
781 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml2(mml->prepare());
782 MEDMeshMultiLev *ptMML2(mml2);
785 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
786 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
787 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
791 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
795 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
799 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
802 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
803 _cached_ds=ret->NewInstance();
804 _cached_ds->ShallowCopy(ret);
808 ret=_cached_ds->NewInstance();
809 ret->ShallowCopy(_cached_ds);
812 appendFields(tr,globs,mml,meshes,ret);
813 // The arrays links to mesh
814 DataArrayInt *famCells(0),*numCells(0);
815 bool noCpyFamCells(false),noCpyNumCells(false);
816 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
819 vtkIntArray *vtkTab(vtkIntArray::New());
820 vtkTab->SetNumberOfComponents(1);
821 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
823 vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
825 { vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famCells->accessToMemArray().setSpecificDeallocator(0); }
826 ret->GetCellData()->AddArray(vtkTab);
830 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
833 vtkIntArray *vtkTab(vtkIntArray::New());
834 vtkTab->SetNumberOfComponents(1);
835 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
837 vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
839 { vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numCells->accessToMemArray().setSpecificDeallocator(0); }
840 ret->GetCellData()->AddArray(vtkTab);
844 // The arrays links to mesh
845 DataArrayInt *famNodes(0),*numNodes(0);
846 bool noCpyFamNodes(false),noCpyNumNodes(false);
847 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
850 vtkIntArray *vtkTab(vtkIntArray::New());
851 vtkTab->SetNumberOfComponents(1);
852 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
854 vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
856 { vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famNodes->accessToMemArray().setSpecificDeallocator(0); }
857 ret->GetPointData()->AddArray(vtkTab);
861 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
864 vtkIntArray *vtkTab(vtkIntArray::New());
865 vtkTab->SetNumberOfComponents(1);
866 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
868 vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
870 { vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numNodes->accessToMemArray().setSpecificDeallocator(0); }
871 ret->GetPointData()->AddArray(vtkTab);
878 //////////////////////
880 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
884 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
887 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
888 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
889 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
890 ret+=(*it2).getNumberOfArrays();
894 void MEDFileFieldRepresentationTree::assignIds() const
897 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
898 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
899 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
903 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
905 std::size_t it0Cnt(0);
906 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
908 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
909 std::string tsName(oss.str());
910 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
912 std::string meshName((*it1)[0].getMeshName());
913 std::size_t it2Cnt(0);
914 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
916 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
917 std::string comSupStr(oss2.str());
918 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
924 void MEDFileFieldRepresentationTree::activateTheFirst() const
926 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
927 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
928 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
930 (*it2).activateAllArrays();
935 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
937 std::size_t it0Cnt(0);
938 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
940 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
941 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
943 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
944 std::vector<double> ts;
945 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
946 std::size_t nbOfTS(dtits.size());
947 std::ostringstream oss3; oss3 << nbOfTS;
948 names.push_back(oss3.str());
949 for(std::size_t i=0;i<nbOfTS;i++)
951 std::ostringstream oss4; oss4 << dtits[i].first;
952 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
953 names.push_back(oss4.str());
954 std::ostringstream oss5; oss5 << dtits[i].second;
955 vtkIdType ItId(sil->AddChild(DtId,edge));
956 names.push_back(oss5.str());
957 std::ostringstream oss6; oss6 << ts[i];
958 sil->AddChild(ItId,edge);
959 names.push_back(oss6.str());
962 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
963 std::string tsName(oss.str());
964 vtkIdType typeId0(sil->AddChild(root,edge));
965 names.push_back(tsName);
966 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
968 std::string meshName((*it1)[0].getMeshName());
969 vtkIdType typeId1(sil->AddChild(typeId0,edge));
970 names.push_back(meshName);
971 std::size_t it2Cnt(0);
972 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
974 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
975 std::string comSupStr(oss2.str());
976 vtkIdType typeId2(sil->AddChild(typeId1,edge));
977 names.push_back(comSupStr);
978 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
984 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
986 int dummy0(0),dummy1(0),dummy2(0);
987 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
988 std::string ret(leaf.getMeshName());
991 for(;i<_ms->getNumberOfMeshes();i++)
993 m=_ms->getMeshAtPos(i);
994 if(m->getName()==ret)
997 if(i==_ms->getNumberOfMeshes())
998 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
999 vtkIdType typeId0(sil->AddChild(root,edge));
1000 names.push_back(m->getName());
1002 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1003 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1004 std::vector<std::string> grps(m->getGroupsNames());
1005 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1007 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1008 names.push_back(*it0);
1009 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1010 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1012 sil->AddChild(typeId2,edge);
1013 names.push_back((*it1).c_str());
1017 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1018 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1019 std::vector<std::string> fams(m->getFamiliesNames());
1020 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1022 sil->AddChild(typeId11,edge);
1023 int famId(m->getFamilyId((*it00).c_str()));
1024 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1025 names.push_back(oss.str());
1030 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1032 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1033 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1034 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1035 if((*it2).containId(id))
1036 return (*it2).getLeafArr(id);
1037 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1040 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1042 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1043 return elt.getZeName();
1046 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1048 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1049 return elt.getZeNameC();
1052 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1054 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1055 return elt.getStatus();
1058 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1061 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1062 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1063 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1064 if((*it2).containZeName(name,ret))
1066 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1067 throw INTERP_KERNEL::Exception(msg.str().c_str());
1070 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1072 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1073 bool ret(elt.setStatus(status));//to be implemented
1077 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1080 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1081 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1082 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1083 ret=std::max(ret,(*it2).getNumberOfTS());
1090 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv)
1094 _ms=MEDFileMeshes::New(fileName);
1095 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1099 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::SauvReader> sr(ParaMEDMEM::SauvReader::New(fileName));
1100 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileData> mfd(sr->loadInMEDFileDS());
1101 _ms=mfd->getMeshes(); _ms->incrRef();
1102 int nbMeshes(_ms->getNumberOfMeshes());
1103 for(int i=0;i<nbMeshes;i++)
1105 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1106 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1108 tmp2->forceComputationOfParts();
1110 _fields=mfd->getFields();
1111 if((ParaMEDMEM::MEDFileFields *)_fields)
1114 if(!((ParaMEDMEM::MEDFileFields *)_fields))
1116 _fields=BuildFieldFromMeshes(_ms);
1120 AppendFieldFromMeshes(_ms,_fields);
1122 _fields->removeFieldsWithoutAnyTimeStep();
1123 std::vector<std::string> meshNames(_ms->getMeshesNames());
1124 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> > fields_per_mesh(meshNames.size());
1125 for(std::size_t i=0;i<meshNames.size();i++)
1127 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1129 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1131 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1133 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1135 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1136 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1137 allFMTSLeavesToDisplaySafe.insert(allFMTSLeavesToDisplaySafe.end(),tmp.begin(),tmp.end());
1140 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1141 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1143 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1145 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1146 // memory safety part
1147 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1148 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1150 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1151 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1153 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1154 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1157 // end of memory safety part
1158 // 1st : timesteps, 2nd : meshName, 3rd : common support
1159 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1160 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1162 vtksys_stl::vector< vtksys_stl::string > meshNamesLoc;
1163 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1164 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1166 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1167 vtksys_stl::vector< vtksys_stl::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1168 if(it==meshNamesLoc.end())
1170 meshNamesLoc.push_back(meshName);
1171 splitByMeshName.resize(splitByMeshName.size()+1);
1172 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1175 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1177 _data_structure[i].resize(meshNamesLoc.size());
1178 for(std::size_t j=0;j<splitByMeshName.size();j++)
1180 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFastCellSupportComparator> > fsp;
1181 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1182 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1183 sbmn[k]=splitByMeshName[j][k];
1184 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1185 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1186 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1187 this->_data_structure[i][j].resize(commonSupSplit.size());
1188 for(std::size_t k=0;k<commonSupSplit.size();k++)
1190 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1191 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1193 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1194 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1197 for(std::size_t k=0;k<commonSupSplit.size();k++)
1198 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1201 this->removeEmptyLeaves();
1203 this->computeFullNameInLeaves();
1206 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1208 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1209 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1211 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1212 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1214 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1215 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1217 newSD1.push_back(*it2);
1219 newSD0.push_back(newSD1);
1222 newSD.push_back(newSD0);
1226 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1228 if(compInfos.size()!=1)
1230 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1233 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1235 return _data_structure[0][0][0].getMeshName();
1238 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1241 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1242 return leaf.getTimeSteps(tk);
1245 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1248 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1249 meshName=leaf.getMeshName();
1250 std::vector<double> ts(leaf.getTimeSteps(tk));
1251 std::size_t zeTimeId(0);
1254 std::vector<double> ts2(ts.size());
1255 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1256 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1257 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1260 if(zeTimeId==(int)ts.size())
1261 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1262 if(zeTimeId==(int)ts.size())
1263 {//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.
1264 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1266 double valAttachedToPos(-std::numeric_limits<double>::max());
1267 for(std::size_t i=0;i<ts.size();i++)
1271 if(ts[i]>valAttachedToPos)
1274 valAttachedToPos=ts[i];
1279 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1280 valAttachedToPos=std::numeric_limits<double>::max();
1281 for(std::size_t i=0;i<ts.size();i++)
1283 if(ts[i]<valAttachedToPos)
1286 valAttachedToPos=ts[i];
1291 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1292 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1293 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1294 std::cerr << oss.str() << std::endl;
1298 tr=new MEDStdTimeReq((int)zeTimeId);
1300 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1301 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1306 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1308 int nbOfActivated(0);
1309 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1310 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1311 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1312 if((*it2).isActivated())
1314 if(nbOfActivated!=1)
1316 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1317 throw INTERP_KERNEL::Exception(oss.str().c_str());
1319 int i0(0),i1(0),i2(0);
1320 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1321 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1322 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1323 if((*it2).isActivated())
1325 lev0=i0; lev1=i1; lev2=i2;
1328 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1331 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1334 os << "#############################################" << std::endl;
1335 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1338 os << "TS" << i << std::endl;
1339 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1342 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1345 os << " " << (*it2).getMeshName() << std::endl;
1346 os << " Comp" << k << std::endl;
1347 (*it2).printMySelf(os);
1351 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1354 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms, ParaMEDMEM::MEDFileFields *ret)
1356 for(int i=0;i<ms->getNumberOfMeshes();i++)
1358 MEDFileMesh *mm(ms->getMeshAtPos(i));
1359 std::vector<int> levs(mm->getNonEmptyLevels());
1360 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileField1TS> f1tsMultiLev(ParaMEDMEM::MEDFileField1TS::New());
1361 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1364 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1366 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1367 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1369 ParaMEDMEM::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1370 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1372 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1373 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1376 f->setName(mm->getName());
1377 f1tsMultiLev->setFieldNoProfileSBT(f);
1382 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1383 if(levsExt.size()==levs.size()+1)
1385 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getGenMeshAtLevel(1));
1386 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_NODES));
1388 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1389 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1390 arr->iota(); f->setArray(arr);
1391 f->setName(mm->getName());
1392 f1tsMultiLev->setFieldNoProfileSBT(f);
1400 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getGenMeshAtLevel(0));
1401 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1403 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1404 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1407 f->setName(mm->getName());
1408 f1tsMultiLev->setFieldNoProfileSBT(f);
1411 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFieldMultiTS> fmtsMultiLev(ParaMEDMEM::MEDFileFieldMultiTS::New());
1412 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1413 ret->pushField(fmtsMultiLev);
1417 ParaMEDMEM::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms)
1419 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFields> ret(ParaMEDMEM::MEDFileFields::New());
1420 AppendFieldFromMeshes(ms,ret);
1424 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1426 std::vector<std::string> ret;
1428 while(pos!=std::string::npos)
1430 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1431 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1433 pos=fullFieldName.find_first_not_of(sep,curPos);
1439 * Here the non regression tests.
1440 * const char inp0[]="";
1441 * const char exp0[]="";
1442 * const char inp1[]="field";
1443 * const char exp1[]="field";
1444 * const char inp2[]="_________";
1445 * const char exp2[]="_________";
1446 * const char inp3[]="field_p";
1447 * const char exp3[]="field_p";
1448 * const char inp4[]="field__p";
1449 * const char exp4[]="field_p";
1450 * const char inp5[]="field_p__";
1451 * const char exp5[]="field_p";
1452 * const char inp6[]="field_p_";
1453 * const char exp6[]="field_p";
1454 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1455 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1456 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1457 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1458 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1459 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1460 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1461 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1463 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1465 static const char SEP('_');
1466 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1468 return fullFieldName;//should never happen
1472 return fullFieldName;
1476 std::string ret(v[0]);
1477 for(std::size_t i=1;i<v.size();i++)
1482 { ret+=SEP; ret+=v[i]; }
1488 return fullFieldName;
1494 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1498 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1503 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1505 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1507 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1513 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1514 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1515 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1517 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1519 std::size_t sz(ts.size());
1520 bool isInHumanRange(true);
1522 for(std::size_t i=0;i<sz;i++)
1525 if(ts[i]<=-1e299 || ts[i]>=1e299)
1526 isInHumanRange=false;
1529 return processedUsingPairOfIds(tsPairs);
1531 return processedUsingPairOfIds(tsPairs);
1532 _postprocessed_time=ts;
1533 return getPostProcessedTime();
1538 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1540 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1542 std::size_t sz(ts.size());
1543 std::vector<double> ts2(sz);
1545 for(std::size_t i=0;i<sz;i++)
1550 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1553 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1555 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1559 void TimeKeeper::printSelf(std::ostream& oss) const
1561 std::size_t sz(_activated_ts.size());
1562 for(std::size_t i=0;i<sz;i++)
1564 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1568 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1570 std::size_t sz(_activated_ts.size());
1571 std::vector<bool> ret(sz);
1572 for(std::size_t i=0;i<sz;i++)
1574 ret[i]=_activated_ts[i].first;
1579 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1581 std::size_t sz(tsPairs.size());
1582 std::set<int> s0,s1;
1583 for(std::size_t i=0;i<sz;i++)
1584 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1587 _postprocessed_time.resize(sz);
1588 for(std::size_t i=0;i<sz;i++)
1589 _postprocessed_time[i]=(double)tsPairs[i].first;
1590 return getPostProcessedTime();
1594 _postprocessed_time.resize(sz);
1595 for(std::size_t i=0;i<sz;i++)
1596 _postprocessed_time[i]=(double)tsPairs[i].second;
1597 return getPostProcessedTime();
1599 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1600 _postprocessed_time.resize(sz);
1601 for(std::size_t i=0;i<sz;i++)
1602 _postprocessed_time[i]=(double)i;
1603 return getPostProcessedTime();
1606 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1608 _activated_ts.resize(maxNumberOfTS);
1609 for(int i=0;i<maxNumberOfTS;i++)
1611 std::ostringstream oss; oss << "000" << i;
1612 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());