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 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, const std::string& tsName, const std::string& meshName, const std::string& comSupStr, std::vector<std::string>& names) const
239 vtkIdType refId(sil->AddChild(root,edge));
240 names.push_back(_ze_name);
241 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
242 _ze_full_name=oss3.str();
244 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
246 sil->AddChild(refId,edge);
247 names.push_back(std::string());
251 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
256 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
258 bool ret(_activated!=status);
263 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshStruct *mst, vtkDataSet *ds) const
265 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
266 static const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
267 tr->setNumberOfTS((operator->())->getNumberOfTS());
269 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
271 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
272 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
273 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
274 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
275 DataArray *crudeArr(0),*postProcessedArr(0);
277 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
279 crudeArr=f1tsPtrInt->getUndergroundDataArray();
281 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
282 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
283 f1ts->loadArraysIfNecessary();
284 MEDCouplingAutoRefCountObjectPtr<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
287 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
289 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
290 vtkFieldData *att(0);
295 att=ds->GetCellData();
300 att=ds->GetPointData();
305 att=ds->GetFieldData();
310 att=ds->GetFieldData();
314 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
318 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
319 vtkDoubleArray *vtkd(vtkDoubleArray::New());
320 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
321 for(int i=0;i<vd->getNumberOfComponents();i++)
322 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
323 if(postProcessedArr!=crudeArr)
325 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vd->accessToMemArray().setSpecificDeallocator(0);
329 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
331 std::string name(tr->buildName(f1ts->getName()));
332 vtkd->SetName(name.c_str());
335 if(discs[0]==ON_GAUSS_PT)
338 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
340 if(discs[0]==ON_GAUSS_NE)
342 vtkIdTypeArray *elno(vtkIdTypeArray::New());
343 elno->SetNumberOfComponents(1);
344 vtkIdType ncell(ds->GetNumberOfCells());
345 int *pt(new int[ncell]),offset(0);
346 std::set<int> cellTypes;
347 for(vtkIdType cellId=0;cellId<ncell;cellId++)
349 vtkCell *cell(ds->GetCell(cellId));
350 int delta(cell->GetNumberOfPoints());
351 cellTypes.insert(cell->GetCellType());
355 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
356 elno->SetArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
357 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
358 elno->SetName(nameElno.c_str());
359 ds->GetCellData()->AddArray(elno);
360 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
361 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
363 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
364 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
366 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
367 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
369 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
370 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
371 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
372 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
373 double *shape(new double[nbGaussPt*nbGaussPt]);
375 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(ct,dummy));
376 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
377 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
378 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
379 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
380 calculator.initLocalInfo();
381 for(int i=0;i<nbGaussPt;i++)
383 const double *pt0(calculator.getFunctionValues(i));
384 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
386 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
388 key->Set(elno->GetInformation(),def,*it);
389 key->Set(vtkd->GetInformation(),def,*it);
398 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
399 vtkIntArray *vtkd(vtkIntArray::New());
400 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
401 for(int i=0;i<vi->getNumberOfComponents();i++)
402 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
403 if(postProcessedArr!=crudeArr)
405 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vi->accessToMemArray().setSpecificDeallocator(0);
409 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
411 std::string name(tr->buildName(f1ts->getName()));
412 vtkd->SetName(name.c_str());
417 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
421 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
423 _elga_cmp.appendELGAIfAny(ds);
428 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
432 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS> >& arr,
433 const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
435 for(std::size_t i=0;i<arr.size();i++)
436 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
439 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
442 _cached_ds->Delete();
445 bool MEDFileFieldRepresentationLeaves::empty() const
447 const MEDFileFastCellSupportComparator *fcscp(_fsp);
448 return fcscp==0 || _arrays.empty();
451 void MEDFileFieldRepresentationLeaves::setId(int& id) const
453 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
457 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
459 return _arrays[0]->getMeshName();
462 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
464 return (int)_arrays.size();
467 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
469 return _arrays[0]->getNumberOfTS();
473 * \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.
475 void MEDFileFieldRepresentationLeaves::feedSIL(const ParaMEDMEM::MEDFileMeshes *ms, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, const std::string& tsName, const std::string& meshName, const std::string& comSupStr, std::vector<std::string>& names) const
477 vtkIdType root2(sil->AddChild(root,edge));
478 names.push_back(std::string("Arrs"));
479 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
480 (*it).feedSIL(sil,root2,edge,tsName,meshName,comSupStr,names);
482 vtkIdType root3(sil->AddChild(root,edge));
483 names.push_back(std::string("InfoOnGeoType"));
484 const ParaMEDMEM::MEDFileMesh *m(0);
486 m=ms->getMeshWithName(meshName);
487 const ParaMEDMEM::MEDFileFastCellSupportComparator *fsp(_fsp);
488 if(!fsp || fsp->getNumberOfTS()==0)
490 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
491 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
493 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
494 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
495 sil->AddChild(root3,edge);
496 names.push_back(cmStr);
500 bool MEDFileFieldRepresentationLeaves::containId(int id) const
502 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
503 if((*it).getId()==id)
508 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
510 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
511 if((*it).getZeName()==name)
519 bool MEDFileFieldRepresentationLeaves::isActivated() const
521 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
522 if((*it).getStatus())
527 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
529 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
531 os << " - " << (*it0).getZeName() << " (";
532 if((*it0).getStatus())
536 os << ")" << std::endl;
540 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
542 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
543 (*it).setStatus(true);
546 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
548 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
549 if((*it).getId()==id)
551 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
554 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
557 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
558 std::vector<double> ret;
559 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
560 return tk.getTimeStepsRegardingPolicy(dtits,ret);
563 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
566 return _arrays[0]->getTimeSteps(ts);
570 return std::vector< std::pair<int,int> >();
574 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
576 std::ostringstream oss;
577 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
578 std::vector<double> ret1;
579 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
580 std::size_t sz(ret1.size());
581 for(std::size_t i=0;i<sz;i++)
583 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
586 std::string tmp(oss.str());
587 if(tmp.size()>200 && i!=sz-1)
597 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshes *meshes, vtkDataSet *ds) const
600 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
601 MEDCouplingAutoRefCountObjectPtr<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
602 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
603 if((*it).getStatus())
605 (*it).appendFields(tr,globs,mml,mst,ds);
606 (*it).appendELGAIfAny(ds);
610 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
612 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
613 DataArrayDouble *coordsMC(0);
614 DataArrayByte *typesMC(0);
615 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
616 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
617 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coordsSafe(coordsMC);
618 MEDCouplingAutoRefCountObjectPtr<DataArrayByte> typesSafe(typesMC);
619 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
621 int nbOfCells(typesSafe->getNbOfElems());
622 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
623 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
624 cellTypes->SetArray(reinterpret_cast<unsigned char *>(typesSafe->getPointer()),nbOfCells,0,VTK_DATA_ARRAY_FREE); typesSafe->accessToMemArray().setSpecificDeallocator(0);
625 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
626 cellLocations->SetArray(cellLocationsSafe->getPointer(),nbOfCells,0,VTK_DATA_ARRAY_FREE); cellLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
627 vtkCellArray *cells(vtkCellArray::New());
628 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
629 cells2->SetArray(cellsSafe->getPointer(),cellsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); cellsSafe->accessToMemArray().setSpecificDeallocator(0);
630 cells->SetCells(nbOfCells,cells2);
632 if(faceLocationsMC!=0 && facesMC!=0)
634 vtkIdTypeArray *faces(vtkIdTypeArray::New());
635 faces->SetArray(facesSafe->getPointer(),facesSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); facesSafe->accessToMemArray().setSpecificDeallocator(0);
636 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
637 faceLocations->SetArray(faceLocationsSafe->getPointer(),faceLocationsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); faceLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
638 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
639 faceLocations->Delete();
643 ret->SetCells(cellTypes,cellLocations,cells);
645 cellLocations->Delete();
647 vtkPoints *pts(vtkPoints::New());
648 vtkDoubleArray *pts2(vtkDoubleArray::New());
649 pts2->SetNumberOfComponents(3);
652 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
653 coordsSafe->accessToMemArray().setSpecificDeallocator(0);
656 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
665 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(ParaMEDMEM::MEDCMeshMultiLev *mm) const
667 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
669 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
670 vtkDoubleArray *vtkTmp(0);
671 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
672 std::size_t dim(arrs.size());
674 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
675 int sizePerAxe[3]={1,1,1};
676 sizePerAxe[0]=arrs[0]->getNbOfElems();
678 sizePerAxe[1]=arrs[1]->getNbOfElems();
680 sizePerAxe[2]=arrs[2]->getNbOfElems();
681 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
682 vtkTmp=vtkDoubleArray::New();
683 vtkTmp->SetNumberOfComponents(1);
685 vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
687 { vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[0]->accessToMemArray().setSpecificDeallocator(0); }
688 ret->SetXCoordinates(vtkTmp);
693 vtkTmp=vtkDoubleArray::New();
694 vtkTmp->SetNumberOfComponents(1);
696 vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
698 { vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[1]->accessToMemArray().setSpecificDeallocator(0); }
699 ret->SetYCoordinates(vtkTmp);
705 vtkTmp=vtkDoubleArray::New();
706 vtkTmp->SetNumberOfComponents(1);
708 vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
710 { vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[2]->accessToMemArray().setSpecificDeallocator(0); }
711 ret->SetZCoordinates(vtkTmp);
718 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(ParaMEDMEM::MEDCurveLinearMeshMultiLev *mm) const
720 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
721 int meshStr[3]={1,1,1};
722 DataArrayDouble *coords(0);
723 std::vector<int> nodeStrct;
725 mm->buildVTUArrays(coords,nodeStrct,isInternal);
726 std::size_t dim(nodeStrct.size());
728 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
729 meshStr[0]=nodeStrct[0];
731 meshStr[1]=nodeStrct[1];
733 meshStr[2]=nodeStrct[2];
734 vtkStructuredGrid *ret(vtkStructuredGrid::New());
735 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
736 vtkDoubleArray *da(vtkDoubleArray::New());
737 da->SetNumberOfComponents(3);
738 if(coords->getNumberOfComponents()==3)
741 da->SetArray(coords->getPointer(),coords->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);//VTK has not the ownership of double * because MEDLoader main struct has it !
743 { da->SetArray(coords->getPointer(),coords->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); coords->accessToMemArray().setSpecificDeallocator(0); }
747 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
748 da->SetArray(coords2->getPointer(),coords2->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);//let VTK deal with double *
749 coords2->accessToMemArray().setSpecificDeallocator(0);
752 vtkPoints *points=vtkPoints::New();
753 ret->SetPoints(points);
760 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDFileMeshes *meshes) const
762 static const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
764 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
765 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
766 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml2(mml->prepare());
767 MEDMeshMultiLev *ptMML2(mml2);
770 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
771 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
772 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
776 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
780 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
784 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
787 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
788 _cached_ds=ret->NewInstance();
789 _cached_ds->ShallowCopy(ret);
793 ret=_cached_ds->NewInstance();
794 ret->ShallowCopy(_cached_ds);
797 appendFields(tr,globs,mml,meshes,ret);
798 // The arrays links to mesh
799 DataArrayInt *famCells(0),*numCells(0);
800 bool noCpyFamCells(false),noCpyNumCells(false);
801 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
804 vtkIntArray *vtkTab(vtkIntArray::New());
805 vtkTab->SetNumberOfComponents(1);
806 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
808 vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
810 { vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famCells->accessToMemArray().setSpecificDeallocator(0); }
811 ret->GetCellData()->AddArray(vtkTab);
815 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
818 vtkIntArray *vtkTab(vtkIntArray::New());
819 vtkTab->SetNumberOfComponents(1);
820 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
822 vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
824 { vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numCells->accessToMemArray().setSpecificDeallocator(0); }
825 ret->GetCellData()->AddArray(vtkTab);
829 // The arrays links to mesh
830 DataArrayInt *famNodes(0),*numNodes(0);
831 bool noCpyFamNodes(false),noCpyNumNodes(false);
832 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
835 vtkIntArray *vtkTab(vtkIntArray::New());
836 vtkTab->SetNumberOfComponents(1);
837 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
839 vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
841 { vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famNodes->accessToMemArray().setSpecificDeallocator(0); }
842 ret->GetPointData()->AddArray(vtkTab);
846 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
849 vtkIntArray *vtkTab(vtkIntArray::New());
850 vtkTab->SetNumberOfComponents(1);
851 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
853 vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
855 { vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numNodes->accessToMemArray().setSpecificDeallocator(0); }
856 ret->GetPointData()->AddArray(vtkTab);
863 //////////////////////
865 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
869 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
872 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
873 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
874 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
875 ret+=(*it2).getNumberOfArrays();
879 void MEDFileFieldRepresentationTree::assignIds() const
882 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
883 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
884 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
887 void MEDFileFieldRepresentationTree::activateTheFirst() const
889 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
890 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
891 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
893 (*it2).activateAllArrays();
898 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
900 std::size_t it0Cnt(0);
901 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
903 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
904 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
906 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
907 std::vector<double> ts;
908 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
909 std::size_t nbOfTS(dtits.size());
910 std::ostringstream oss3; oss3 << nbOfTS;
911 names.push_back(oss3.str());
912 for(std::size_t i=0;i<nbOfTS;i++)
914 std::ostringstream oss4; oss4 << dtits[i].first;
915 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
916 names.push_back(oss4.str());
917 std::ostringstream oss5; oss5 << dtits[i].second;
918 vtkIdType ItId(sil->AddChild(DtId,edge));
919 names.push_back(oss5.str());
920 std::ostringstream oss6; oss6 << ts[i];
921 sil->AddChild(ItId,edge);
922 names.push_back(oss6.str());
925 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
926 std::string tsName(oss.str());
927 vtkIdType typeId0(sil->AddChild(root,edge));
928 names.push_back(tsName);
929 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
931 std::string meshName((*it1)[0].getMeshName());
932 vtkIdType typeId1(sil->AddChild(typeId0,edge));
933 names.push_back(meshName);
934 std::size_t it2Cnt(0);
935 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
937 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
938 std::string comSupStr(oss2.str());
939 vtkIdType typeId2(sil->AddChild(typeId1,edge));
940 names.push_back(comSupStr);
941 (*it2).feedSIL(_ms,sil,typeId2,edge,tsName,meshName,comSupStr,names);
947 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
949 int dummy0(0),dummy1(0),dummy2(0);
950 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
951 std::string ret(leaf.getMeshName());
954 for(;i<_ms->getNumberOfMeshes();i++)
956 m=_ms->getMeshAtPos(i);
957 if(m->getName()==ret)
960 if(i==_ms->getNumberOfMeshes())
961 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
962 vtkIdType typeId0(sil->AddChild(root,edge));
963 names.push_back(m->getName());
965 vtkIdType typeId1(sil->AddChild(typeId0,edge));
966 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
967 std::vector<std::string> grps(m->getGroupsNames());
968 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
970 vtkIdType typeId2(sil->AddChild(typeId1,edge));
971 names.push_back(*it0);
972 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
973 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
975 sil->AddChild(typeId2,edge);
976 names.push_back((*it1).c_str());
980 vtkIdType typeId11(sil->AddChild(typeId0,edge));
981 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
982 std::vector<std::string> fams(m->getFamiliesNames());
983 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
985 sil->AddChild(typeId11,edge);
986 int famId(m->getFamilyId((*it00).c_str()));
987 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
988 names.push_back(oss.str());
993 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
995 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
996 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
997 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
998 if((*it2).containId(id))
999 return (*it2).getLeafArr(id);
1000 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1003 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1005 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1006 return elt.getZeName();
1009 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1011 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1012 return elt.getStatus();
1015 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1018 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1019 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1020 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1021 if((*it2).containZeName(name,ret))
1023 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1024 throw INTERP_KERNEL::Exception(msg.str().c_str());
1027 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1029 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1030 bool ret(elt.setStatus(status));//to be implemented
1034 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1037 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1038 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1039 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1040 ret=std::max(ret,(*it2).getNumberOfTS());
1047 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv)
1051 _ms=MEDFileMeshes::New(fileName);
1052 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1056 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::SauvReader> sr(ParaMEDMEM::SauvReader::New(fileName));
1057 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileData> mfd(sr->loadInMEDFileDS());
1058 _ms=mfd->getMeshes(); _ms->incrRef();
1059 int nbMeshes(_ms->getNumberOfMeshes());
1060 for(int i=0;i<nbMeshes;i++)
1062 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1063 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1065 tmp2->forceComputationOfParts();
1067 _fields=mfd->getFields();
1068 if((ParaMEDMEM::MEDFileFields *)_fields)
1071 if(!((ParaMEDMEM::MEDFileFields *)_fields))
1073 _fields=BuildFieldFromMeshes(_ms);
1077 AppendFieldFromMeshes(_ms,_fields);
1079 _fields->removeFieldsWithoutAnyTimeStep();
1080 std::vector<std::string> meshNames(_ms->getMeshesNames());
1081 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> > fields_per_mesh(meshNames.size());
1082 for(std::size_t i=0;i<meshNames.size();i++)
1084 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1086 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1088 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1090 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1092 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1093 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1094 allFMTSLeavesToDisplaySafe.insert(allFMTSLeavesToDisplaySafe.end(),tmp.begin(),tmp.end());
1097 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1098 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1100 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1102 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1103 // memory safety part
1104 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1105 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1107 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1108 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1110 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1111 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1114 // end of memory safety part
1115 // 1st : timesteps, 2nd : meshName, 3rd : common support
1116 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1117 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1119 vtksys_stl::vector< vtksys_stl::string > meshNamesLoc;
1120 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1121 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1123 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1124 vtksys_stl::vector< vtksys_stl::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1125 if(it==meshNamesLoc.end())
1127 meshNamesLoc.push_back(meshName);
1128 splitByMeshName.resize(splitByMeshName.size()+1);
1129 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1132 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1134 _data_structure[i].resize(meshNamesLoc.size());
1135 for(std::size_t j=0;j<splitByMeshName.size();j++)
1137 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFastCellSupportComparator> > fsp;
1138 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1139 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1140 sbmn[k]=splitByMeshName[j][k];
1141 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1142 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1143 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1144 this->_data_structure[i][j].resize(commonSupSplit.size());
1145 for(std::size_t k=0;k<commonSupSplit.size();k++)
1147 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1148 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1150 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1151 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1154 for(std::size_t k=0;k<commonSupSplit.size();k++)
1155 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1158 this->removeEmptyLeaves();
1162 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1164 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1165 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1167 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1168 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1170 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1171 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1173 newSD1.push_back(*it2);
1175 newSD0.push_back(newSD1);
1178 newSD.push_back(newSD0);
1182 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1184 if(compInfos.size()!=1)
1186 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1189 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1191 return _data_structure[0][0][0].getMeshName();
1194 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1197 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1198 return leaf.getTimeSteps(tk);
1201 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1204 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1205 meshName=leaf.getMeshName();
1206 std::vector<double> ts(leaf.getTimeSteps(tk));
1207 std::size_t zeTimeId(0);
1210 std::vector<double> ts2(ts.size());
1211 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1212 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1213 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1216 if(zeTimeId==(int)ts.size())
1217 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1218 if(zeTimeId==(int)ts.size())
1219 {//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.
1220 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1222 double valAttachedToPos(-std::numeric_limits<double>::max());
1223 for(std::size_t i=0;i<ts.size();i++)
1227 if(ts[i]>valAttachedToPos)
1230 valAttachedToPos=ts[i];
1235 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1236 valAttachedToPos=std::numeric_limits<double>::max();
1237 for(std::size_t i=0;i<ts.size();i++)
1239 if(ts[i]<valAttachedToPos)
1242 valAttachedToPos=ts[i];
1247 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1248 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1249 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1250 std::cerr << oss.str() << std::endl;
1254 tr=new MEDStdTimeReq((int)zeTimeId);
1256 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1257 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1262 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1264 int nbOfActivated(0);
1265 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1266 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1267 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1268 if((*it2).isActivated())
1270 if(nbOfActivated!=1)
1272 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1273 throw INTERP_KERNEL::Exception(oss.str().c_str());
1275 int i0(0),i1(0),i2(0);
1276 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1277 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1278 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1279 if((*it2).isActivated())
1281 lev0=i0; lev1=i1; lev2=i2;
1284 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1287 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1290 os << "#############################################" << std::endl;
1291 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1294 os << "TS" << i << std::endl;
1295 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1298 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1301 os << " " << (*it2).getMeshName() << std::endl;
1302 os << " Comp" << k << std::endl;
1303 (*it2).printMySelf(os);
1307 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1310 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms, ParaMEDMEM::MEDFileFields *ret)
1312 for(int i=0;i<ms->getNumberOfMeshes();i++)
1314 MEDFileMesh *mm(ms->getMeshAtPos(i));
1315 std::vector<int> levs(mm->getNonEmptyLevels());
1316 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileField1TS> f1tsMultiLev(ParaMEDMEM::MEDFileField1TS::New());
1317 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1320 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1322 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1323 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1325 ParaMEDMEM::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1326 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1328 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1329 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1332 f->setName(mm->getName());
1333 f1tsMultiLev->setFieldNoProfileSBT(f);
1338 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1339 if(levsExt.size()==levs.size()+1)
1341 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getGenMeshAtLevel(1));
1342 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_NODES));
1344 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1345 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1346 arr->iota(); f->setArray(arr);
1347 f->setName(mm->getName());
1348 f1tsMultiLev->setFieldNoProfileSBT(f);
1356 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getGenMeshAtLevel(0));
1357 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1359 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1360 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1363 f->setName(mm->getName());
1364 f1tsMultiLev->setFieldNoProfileSBT(f);
1367 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFieldMultiTS> fmtsMultiLev(ParaMEDMEM::MEDFileFieldMultiTS::New());
1368 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1369 ret->pushField(fmtsMultiLev);
1373 ParaMEDMEM::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms)
1375 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFields> ret(ParaMEDMEM::MEDFileFields::New());
1376 AppendFieldFromMeshes(ms,ret);
1380 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1382 std::vector<std::string> ret;
1384 while(pos!=std::string::npos)
1386 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1387 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1389 pos=fullFieldName.find_first_not_of(sep,curPos);
1395 * Here the non regression tests.
1396 * const char inp0[]="";
1397 * const char exp0[]="";
1398 * const char inp1[]="field";
1399 * const char exp1[]="field";
1400 * const char inp2[]="_________";
1401 * const char exp2[]="_________";
1402 * const char inp3[]="field_p";
1403 * const char exp3[]="field_p";
1404 * const char inp4[]="field__p";
1405 * const char exp4[]="field_p";
1406 * const char inp5[]="field_p__";
1407 * const char exp5[]="field_p";
1408 * const char inp6[]="field_p_";
1409 * const char exp6[]="field_p";
1410 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1411 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1412 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1413 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1414 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1415 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1416 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1417 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1419 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1421 static const char SEP('_');
1422 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1424 return fullFieldName;//should never happen
1428 return fullFieldName;
1432 std::string ret(v[0]);
1433 for(std::size_t i=1;i<v.size();i++)
1438 { ret+=SEP; ret+=v[i]; }
1444 return fullFieldName;
1450 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1454 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1459 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1461 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1463 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1469 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1470 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1471 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1473 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1475 std::size_t sz(ts.size());
1476 bool isInHumanRange(true);
1478 for(std::size_t i=0;i<sz;i++)
1481 if(ts[i]<=-1e299 || ts[i]>=1e299)
1482 isInHumanRange=false;
1485 return processedUsingPairOfIds(tsPairs);
1487 return processedUsingPairOfIds(tsPairs);
1488 _postprocessed_time=ts;
1489 return getPostProcessedTime();
1494 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1496 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1498 std::size_t sz(ts.size());
1499 std::vector<double> ts2(sz);
1501 for(std::size_t i=0;i<sz;i++)
1506 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1509 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1511 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1515 void TimeKeeper::printSelf(std::ostream& oss) const
1517 std::size_t sz(_activated_ts.size());
1518 for(std::size_t i=0;i<sz;i++)
1520 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1524 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1526 std::size_t sz(_activated_ts.size());
1527 std::vector<bool> ret(sz);
1528 for(std::size_t i=0;i<sz;i++)
1530 ret[i]=_activated_ts[i].first;
1535 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1537 std::size_t sz(tsPairs.size());
1538 std::set<int> s0,s1;
1539 for(std::size_t i=0;i<sz;i++)
1540 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1543 _postprocessed_time.resize(sz);
1544 for(std::size_t i=0;i<sz;i++)
1545 _postprocessed_time[i]=(double)tsPairs[i].first;
1546 return getPostProcessedTime();
1550 _postprocessed_time.resize(sz);
1551 for(std::size_t i=0;i<sz;i++)
1552 _postprocessed_time[i]=(double)tsPairs[i].second;
1553 return getPostProcessedTime();
1555 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1556 _postprocessed_time.resize(sz);
1557 for(std::size_t i=0;i<sz;i++)
1558 _postprocessed_time[i]=(double)i;
1559 return getPostProcessedTime();
1562 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1564 _activated_ts.resize(maxNumberOfTS);
1565 for(int i=0;i<maxNumberOfTS;i++)
1567 std::ostringstream oss; oss << "000" << i;
1568 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());