1 // Copyright (C) 2010-2016 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay
21 #include "MEDTimeReq.hxx"
22 #include "MEDUtilities.hxx"
24 #include "MEDFileFieldRepresentationTree.hxx"
25 #include "MEDCouplingFieldDiscretization.hxx"
26 #include "MEDCouplingFieldDouble.hxx"
27 #include "InterpKernelGaussCoords.hxx"
28 #include "MEDFileData.hxx"
29 #include "SauvReader.hxx"
31 #ifdef MEDREADER_USE_MPI
32 #include "ParaMEDFileMesh.hxx"
35 #include "vtkXMLUnstructuredGridWriter.h"//
37 #include "vtkUnstructuredGrid.h"
38 #include "vtkRectilinearGrid.h"
39 #include "vtkStructuredGrid.h"
40 #include "vtkUnsignedCharArray.h"
41 #include "vtkQuadratureSchemeDefinition.h"
42 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
43 #include "vtkInformationIntegerKey.h"
44 #include "vtkInformation.h"
45 #include "vtkDataArrayTemplate.h"
46 #include "vtkIdTypeArray.h"
47 #include "vtkDoubleArray.h"
48 #include "vtkIntArray.h"
49 #include "vtkCellArray.h"
50 #include "vtkPointData.h"
51 #include "vtkFieldData.h"
52 #include "vtkCellData.h"
54 #include "vtkMutableDirectedGraph.h"
56 using namespace MEDCoupling;
58 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
60 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
62 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
64 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
66 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
68 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
70 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
72 const char MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME[]="GlobalNodeIds";// WARNING DO NOT CHANGE IT BEFORE HAVING CHECKED IN PV SOURCES !
74 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
76 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
78 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
80 vtkIdTypeArray *ELGACmp::findOrCreate(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds, bool& isNew) const
82 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
91 return createNew(globs,locsReallyUsed,vtkd,ds);
95 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd) const
97 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
98 if(it==_loc_names.end())
100 std::size_t pos(std::distance(_loc_names.begin(),it));
101 vtkIdTypeArray *ret(_elgas[pos]);
102 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
103 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
105 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
107 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
111 vtkIdTypeArray *ELGACmp::createNew(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds) const
113 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
114 std::vector< std::vector<std::string> > locNames(_loc_names);
115 std::vector<vtkIdTypeArray *> elgas(_elgas);
116 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
118 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
119 if(it!=locNames.end())
120 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
121 locNames.push_back(locsReallyUsed);
122 vtkIdTypeArray *elga(vtkIdTypeArray::New());
123 elga->SetNumberOfComponents(1);
124 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
125 std::map<unsigned char,int> m;
126 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
128 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
129 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
130 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
131 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
132 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
133 // 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.
134 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
135 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
136 double *shape(new double[nbPtsPerCell*nbGaussPt]);
137 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
138 calculator.initLocalInfo();
139 const std::vector<double>& wgths(loc.getGaussWeights());
140 for(int i=0;i<nbGaussPt;i++)
142 const double *pt0(calculator.getFunctionValues(i));
143 if(ct!=INTERP_KERNEL::NORM_HEXA27)
144 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
147 for(int j=0;j<27;j++)
148 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
151 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
152 m[vtkType]=nbGaussPt;
153 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
155 key->Set(elga->GetInformation(),def,vtkType);
156 key->Set(vtkd->GetInformation(),def,vtkType);
157 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
160 vtkIdType ncell(ds->GetNumberOfCells());
161 int *pt(new int[ncell]),offset(0);
162 for(vtkIdType cellId=0;cellId<ncell;cellId++)
164 vtkCell *cell(ds->GetCell(cellId));
165 int delta(m[cell->GetCellType()]);
169 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
170 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
171 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
172 std::string ossStr(oss.str());
173 elga->SetName(ossStr.c_str());
174 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
175 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
176 elgas.push_back(elga);
180 _defs.push_back(defs);
184 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
186 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
187 ds->GetCellData()->AddArray(*it);
192 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
194 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
195 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
196 (*it1).first->Delete();
202 class MEDFileVTKTraits
205 typedef void VtkType;
210 class MEDFileVTKTraits<int>
213 typedef vtkIntArray VtkType;
214 typedef MEDCoupling::DataArrayInt MCType;
218 class MEDFileVTKTraits<double>
221 typedef vtkDoubleArray VtkType;
222 typedef MEDCoupling::DataArrayDouble MCType;
226 void AssignDataPointerToVTK(typename MEDFileVTKTraits<T>::VtkType *vtkTab, typename MEDFileVTKTraits<T>::MCType *mcTab, bool noCpyNumNodes)
229 vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),1,vtkDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE);
231 { vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),0,vtkDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE); mcTab->accessToMemArray().setSpecificDeallocator(0); }
234 // here copy is always assumed.
235 template<class VTKT, class MCT>
236 void AssignDataPointerOther(VTKT *vtkTab, MCT *mcTab, int nbElems)
238 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(mcTab->getPointer()),nbElems,0,VTKT::VTK_DATA_ARRAY_FREE);
239 mcTab->accessToMemArray().setSpecificDeallocator(0);
244 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
248 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>& arr):MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
250 std::vector< std::vector<MEDCoupling::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
252 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
253 if(typs[0].size()!=1)
254 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 !");
255 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
256 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
260 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
262 MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>::operator=(other);
265 _ze_name=other._ze_name;
266 _ze_full_name.clear();
270 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
275 int MEDFileFieldRepresentationLeavesArrays::getId() const
280 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
282 return _ze_full_name;
285 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
287 return _ze_full_name.c_str();
290 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
292 vtkIdType refId(sil->AddChild(root,edge));
293 names.push_back(_ze_name);
295 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
297 sil->AddChild(refId,edge);
298 names.push_back(std::string());
302 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
304 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
305 _ze_full_name=oss3.str();
308 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
313 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
315 bool ret(_activated!=status);
320 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshStruct *mst, vtkDataSet *ds) const
322 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
323 tr->setNumberOfTS((operator->())->getNumberOfTS());
325 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
327 MCAuto<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
328 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
329 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
330 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
331 DataArray *crudeArr(0),*postProcessedArr(0);
333 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
335 crudeArr=f1tsPtrInt->getUndergroundDataArray();
337 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
338 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
339 f1ts->loadArraysIfNecessary();
340 MCAuto<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
343 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
345 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
346 vtkFieldData *att(0);
351 att=ds->GetCellData();
356 att=ds->GetPointData();
361 att=ds->GetFieldData();
366 att=ds->GetFieldData();
370 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
374 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
375 vtkDoubleArray *vtkd(vtkDoubleArray::New());
376 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
377 for(int i=0;i<vd->getNumberOfComponents();i++)
378 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
379 AssignDataPointerToVTK<double>(vtkd,vd,postProcessedArr==crudeArr);
380 std::string name(tr->buildName(f1ts->getName()));
381 vtkd->SetName(name.c_str());
384 if(discs[0]==ON_GAUSS_PT)
387 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
389 if(discs[0]==ON_GAUSS_NE)
391 vtkIdTypeArray *elno(vtkIdTypeArray::New());
392 elno->SetNumberOfComponents(1);
393 vtkIdType ncell(ds->GetNumberOfCells());
394 int *pt(new int[ncell]),offset(0);
395 std::set<int> cellTypes;
396 for(vtkIdType cellId=0;cellId<ncell;cellId++)
398 vtkCell *cell(ds->GetCell(cellId));
399 int delta(cell->GetNumberOfPoints());
400 cellTypes.insert(cell->GetCellType());
404 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
405 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
406 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
407 elno->SetName(nameElno.c_str());
408 ds->GetCellData()->AddArray(elno);
409 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
410 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
412 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
413 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
415 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
416 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
418 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
419 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
420 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
421 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
422 double *shape(new double[nbGaussPt*nbGaussPt]);
424 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));//GetLocsFromGeometricType
425 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
426 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
427 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
428 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
429 calculator.initLocalInfo();
430 for(int i=0;i<nbGaussPt;i++)
432 const double *pt0(calculator.getFunctionValues(i));
433 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
435 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
437 key->Set(elno->GetInformation(),def,*it);
438 key->Set(vtkd->GetInformation(),def,*it);
447 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
448 vtkIntArray *vtkd(vtkIntArray::New());
449 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
450 for(int i=0;i<vi->getNumberOfComponents();i++)
451 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
452 AssignDataPointerToVTK<int>(vtkd,vi,postProcessedArr==crudeArr);
453 std::string name(tr->buildName(f1ts->getName()));
454 vtkd->SetName(name.c_str());
459 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
463 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
465 _elga_cmp.appendELGAIfAny(ds);
470 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
474 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS> >& arr,
475 const MEDCoupling::MCAuto<MEDCoupling::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
477 for(std::size_t i=0;i<arr.size();i++)
478 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
481 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
484 _cached_ds->Delete();
487 bool MEDFileFieldRepresentationLeaves::empty() const
489 const MEDFileFastCellSupportComparator *fcscp(_fsp);
490 return fcscp==0 || _arrays.empty();
493 void MEDFileFieldRepresentationLeaves::setId(int& id) const
495 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
499 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
501 return _arrays[0]->getMeshName();
504 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
506 return (int)_arrays.size();
509 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
511 return _arrays[0]->getNumberOfTS();
514 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
516 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
517 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
521 * \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.
523 void MEDFileFieldRepresentationLeaves::feedSIL(const MEDCoupling::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
525 vtkIdType root2(sil->AddChild(root,edge));
526 names.push_back(std::string("Arrs"));
527 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
528 (*it).feedSIL(sil,root2,edge,names);
530 vtkIdType root3(sil->AddChild(root,edge));
531 names.push_back(std::string("InfoOnGeoType"));
532 const MEDCoupling::MEDFileMesh *m(0);
534 m=ms->getMeshWithName(meshName);
535 const MEDCoupling::MEDFileFastCellSupportComparator *fsp(_fsp);
536 if(!fsp || fsp->getNumberOfTS()==0)
538 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
539 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
541 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
542 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
543 sil->AddChild(root3,edge);
544 names.push_back(cmStr);
548 bool MEDFileFieldRepresentationLeaves::containId(int id) const
550 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
551 if((*it).getId()==id)
556 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
558 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
559 if((*it).getZeName()==name)
567 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
569 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
570 status[(*it).getZeName()]=(*it).getStatus();
573 bool MEDFileFieldRepresentationLeaves::isActivated() const
575 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
576 if((*it).getStatus())
581 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
583 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
585 os << " - " << (*it0).getZeName() << " (";
586 if((*it0).getStatus())
590 os << ")" << std::endl;
594 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
596 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
597 (*it).setStatus(true);
600 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
602 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
603 if((*it).getId()==id)
605 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
608 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
611 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
612 std::vector<double> ret;
613 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
614 return tk.getTimeStepsRegardingPolicy(dtits,ret);
617 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
620 return _arrays[0]->getTimeSteps(ts);
624 return std::vector< std::pair<int,int> >();
628 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
630 std::ostringstream oss;
631 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
632 std::vector<double> ret1;
633 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
634 std::size_t sz(ret1.size());
635 for(std::size_t i=0;i<sz;i++)
637 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
640 std::string tmp(oss.str());
641 if(tmp.size()>200 && i!=sz-1)
651 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshes *meshes, vtkDataSet *ds) const
654 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
655 MCAuto<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
656 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
657 if((*it).getStatus())
659 (*it).appendFields(tr,globs,mml,mst,ds);
660 (*it).appendELGAIfAny(ds);
664 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
666 DataArrayDouble *coordsMC(0);
667 DataArrayByte *typesMC(0);
668 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
669 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
670 MCAuto<DataArrayDouble> coordsSafe(coordsMC);
671 MCAuto<DataArrayByte> typesSafe(typesMC);
672 MCAuto<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
674 int nbOfCells(typesSafe->getNbOfElems());
675 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
676 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
677 AssignDataPointerOther<vtkUnsignedCharArray,DataArrayByte>(cellTypes,typesSafe,nbOfCells);
678 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
679 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(cellLocations,cellLocationsSafe,nbOfCells);
680 vtkCellArray *cells(vtkCellArray::New());
681 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
682 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(cells2,cellsSafe,cellsSafe->getNbOfElems());
683 cells->SetCells(nbOfCells,cells2);
685 if(faceLocationsMC!=0 && facesMC!=0)
687 vtkIdTypeArray *faces(vtkIdTypeArray::New());
688 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(faces,facesSafe,facesSafe->getNbOfElems());
689 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
690 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(faceLocations,faceLocationsSafe,faceLocationsSafe->getNbOfElems());
691 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
692 faceLocations->Delete();
696 ret->SetCells(cellTypes,cellLocations,cells);
698 cellLocations->Delete();
700 vtkPoints *pts(vtkPoints::New());
701 vtkDoubleArray *pts2(vtkDoubleArray::New());
702 pts2->SetNumberOfComponents(3);
703 AssignDataPointerToVTK<double>(pts2,coordsSafe,statusOfCoords);
712 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(MEDCoupling::MEDCMeshMultiLev *mm) const
715 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
716 vtkDoubleArray *vtkTmp(0);
717 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
718 std::size_t dim(arrs.size());
720 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
721 int sizePerAxe[3]={1,1,1};
722 sizePerAxe[0]=arrs[0]->getNbOfElems();
724 sizePerAxe[1]=arrs[1]->getNbOfElems();
726 sizePerAxe[2]=arrs[2]->getNbOfElems();
727 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
728 vtkTmp=vtkDoubleArray::New();
729 vtkTmp->SetNumberOfComponents(1);
730 AssignDataPointerToVTK<double>(vtkTmp,arrs[0],isInternal);
731 ret->SetXCoordinates(vtkTmp);
736 vtkTmp=vtkDoubleArray::New();
737 vtkTmp->SetNumberOfComponents(1);
738 AssignDataPointerToVTK<double>(vtkTmp,arrs[1],isInternal);
739 ret->SetYCoordinates(vtkTmp);
745 vtkTmp=vtkDoubleArray::New();
746 vtkTmp->SetNumberOfComponents(1);
747 AssignDataPointerToVTK<double>(vtkTmp,arrs[2],isInternal);
748 ret->SetZCoordinates(vtkTmp);
755 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(MEDCoupling::MEDCurveLinearMeshMultiLev *mm) const
757 int meshStr[3]={1,1,1};
758 DataArrayDouble *coords(0);
759 std::vector<int> nodeStrct;
761 mm->buildVTUArrays(coords,nodeStrct,isInternal);
762 std::size_t dim(nodeStrct.size());
764 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
765 meshStr[0]=nodeStrct[0];
767 meshStr[1]=nodeStrct[1];
769 meshStr[2]=nodeStrct[2];
770 vtkStructuredGrid *ret(vtkStructuredGrid::New());
771 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
772 vtkDoubleArray *da(vtkDoubleArray::New());
773 da->SetNumberOfComponents(3);
774 if(coords->getNumberOfComponents()==3)
775 AssignDataPointerToVTK<double>(da,coords,isInternal);//if isIntenal==True VTK has not the ownership of double * because MEDLoader main struct has it !
778 MCAuto<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
779 AssignDataPointerToVTK<double>(da,coords2,false);//let VTK deal with double *
782 vtkPoints *points=vtkPoints::New();
783 ret->SetPoints(points);
790 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDFileMeshes *meshes) const
793 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
794 MCAuto<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
795 MCAuto<MEDMeshMultiLev> mml2(mml->prepare());
796 MEDMeshMultiLev *ptMML2(mml2);
799 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
800 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
801 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
805 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
809 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
813 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
816 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
817 _cached_ds=ret->NewInstance();
818 _cached_ds->ShallowCopy(ret);
822 ret=_cached_ds->NewInstance();
823 ret->ShallowCopy(_cached_ds);
826 appendFields(tr,globs,mml,meshes,ret);
827 // The arrays links to mesh
828 DataArrayInt *famCells(0),*numCells(0);
829 bool noCpyFamCells(false),noCpyNumCells(false);
830 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
833 vtkIntArray *vtkTab(vtkIntArray::New());
834 vtkTab->SetNumberOfComponents(1);
835 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
836 AssignDataPointerToVTK<int>(vtkTab,famCells,noCpyFamCells);
837 ret->GetCellData()->AddArray(vtkTab);
841 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
844 vtkIntArray *vtkTab(vtkIntArray::New());
845 vtkTab->SetNumberOfComponents(1);
846 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
847 AssignDataPointerToVTK<int>(vtkTab,numCells,noCpyNumCells);
848 ret->GetCellData()->AddArray(vtkTab);
852 // The arrays links to mesh
853 DataArrayInt *famNodes(0),*numNodes(0);
854 bool noCpyFamNodes(false),noCpyNumNodes(false);
855 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
858 vtkIntArray *vtkTab(vtkIntArray::New());
859 vtkTab->SetNumberOfComponents(1);
860 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
861 AssignDataPointerToVTK<int>(vtkTab,famNodes,noCpyFamNodes);
862 ret->GetPointData()->AddArray(vtkTab);
866 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
869 vtkIntArray *vtkTab(vtkIntArray::New());
870 vtkTab->SetNumberOfComponents(1);
871 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
872 AssignDataPointerToVTK<int>(vtkTab,numNodes,noCpyNumNodes);
873 ret->GetPointData()->AddArray(vtkTab);
877 // Global Node Ids if any ! (In // mode)
878 DataArrayInt *gni(ptMML2->retrieveGlobalNodeIdsIfAny());
881 vtkIntArray *vtkTab(vtkIntArray::New());
882 vtkTab->SetNumberOfComponents(1);
883 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME);
884 AssignDataPointerToVTK<int>(vtkTab,gni,false);
885 ret->GetPointData()->AddArray(vtkTab);
892 //////////////////////
894 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
898 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
901 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
902 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
903 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
904 ret+=(*it2).getNumberOfArrays();
908 void MEDFileFieldRepresentationTree::assignIds() const
911 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
912 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
913 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
917 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
919 std::size_t it0Cnt(0);
920 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
922 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
923 std::string tsName(oss.str());
924 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
926 std::string meshName((*it1)[0].getMeshName());
927 std::size_t it2Cnt(0);
928 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
930 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
931 std::string comSupStr(oss2.str());
932 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
938 void MEDFileFieldRepresentationTree::activateTheFirst() const
940 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
941 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
942 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
944 (*it2).activateAllArrays();
949 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
951 std::size_t it0Cnt(0);
952 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
954 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
955 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
957 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
958 std::vector<double> ts;
959 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
960 std::size_t nbOfTS(dtits.size());
961 std::ostringstream oss3; oss3 << nbOfTS;
962 names.push_back(oss3.str());
963 for(std::size_t i=0;i<nbOfTS;i++)
965 std::ostringstream oss4; oss4 << dtits[i].first;
966 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
967 names.push_back(oss4.str());
968 std::ostringstream oss5; oss5 << dtits[i].second;
969 vtkIdType ItId(sil->AddChild(DtId,edge));
970 names.push_back(oss5.str());
971 std::ostringstream oss6; oss6 << ts[i];
972 sil->AddChild(ItId,edge);
973 names.push_back(oss6.str());
976 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
977 std::string tsName(oss.str());
978 vtkIdType typeId0(sil->AddChild(root,edge));
979 names.push_back(tsName);
980 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
982 std::string meshName((*it1)[0].getMeshName());
983 vtkIdType typeId1(sil->AddChild(typeId0,edge));
984 names.push_back(meshName);
985 std::size_t it2Cnt(0);
986 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
988 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
989 std::string comSupStr(oss2.str());
990 vtkIdType typeId2(sil->AddChild(typeId1,edge));
991 names.push_back(comSupStr);
992 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
998 std::string MEDFileFieldRepresentationTree::getActiveMeshName() const
1000 int dummy0(0),dummy1(0),dummy2(0);
1001 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1002 return leaf.getMeshName();
1005 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1007 int dummy0(0),dummy1(0),dummy2(0);
1008 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1009 std::string ret(leaf.getMeshName());
1012 for(;i<_ms->getNumberOfMeshes();i++)
1014 m=_ms->getMeshAtPos(i);
1015 if(m->getName()==ret)
1018 if(i==_ms->getNumberOfMeshes())
1019 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1020 vtkIdType typeId0(sil->AddChild(root,edge));
1021 names.push_back(m->getName());
1023 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1024 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1025 std::vector<std::string> grps(m->getGroupsNames());
1026 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1028 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1029 names.push_back(*it0);
1030 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1031 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1033 sil->AddChild(typeId2,edge);
1034 names.push_back((*it1).c_str());
1038 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1039 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1040 std::vector<std::string> fams(m->getFamiliesNames());
1041 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1043 sil->AddChild(typeId11,edge);
1044 int famId(m->getFamilyId((*it00).c_str()));
1045 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1046 names.push_back(oss.str());
1051 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1053 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1054 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1055 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1056 if((*it2).containId(id))
1057 return (*it2).getLeafArr(id);
1058 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1061 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1063 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1064 return elt.getZeName();
1067 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1069 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1070 return elt.getZeNameC();
1073 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1075 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1076 return elt.getStatus();
1079 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1082 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1083 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1084 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1085 if((*it2).containZeName(name,ret))
1087 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1088 throw INTERP_KERNEL::Exception(msg.str().c_str());
1091 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1093 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1094 bool ret(elt.setStatus(status));//to be implemented
1098 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1101 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1102 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1103 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1104 ret=std::max(ret,(*it2).getNumberOfTS());
1111 void MEDFileFieldRepresentationTree::loadInMemory(MEDCoupling::MEDFileFields *fields, MEDCoupling::MEDFileMeshes *meshes)
1113 _fields=fields; _ms=meshes;
1114 if(_fields.isNotNull())
1119 if(_fields.isNull())
1121 _fields=BuildFieldFromMeshes(_ms);
1125 AppendFieldFromMeshes(_ms,_fields);
1127 _ms->cartesianizeMe();
1128 _fields->removeFieldsWithoutAnyTimeStep();
1129 std::vector<std::string> meshNames(_ms->getMeshesNames());
1130 std::vector< MCAuto<MEDFileFields> > fields_per_mesh(meshNames.size());
1131 for(std::size_t i=0;i<meshNames.size();i++)
1133 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1135 std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1137 for(std::vector< MCAuto<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1139 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1141 MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1142 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1144 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1146 if(!(*it)->presenceOfMultiDiscPerGeoType())
1147 allFMTSLeavesToDisplaySafe.push_back(*it);
1149 {// The case of some parts of field have more than one discretization per geo type.
1150 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1151 std::size_t it0Cnt(0);
1152 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1154 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1155 (*it0)->setName(oss.str());
1156 allFMTSLeavesToDisplaySafe.push_back(*it0);
1163 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1164 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1166 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1168 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1169 // memory safety part
1170 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1171 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1173 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1174 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1176 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1177 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1180 // end of memory safety part
1181 // 1st : timesteps, 2nd : meshName, 3rd : common support
1182 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1183 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1185 std::vector< std::string > meshNamesLoc;
1186 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1187 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1189 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1190 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1191 if(it==meshNamesLoc.end())
1193 meshNamesLoc.push_back(meshName);
1194 splitByMeshName.resize(splitByMeshName.size()+1);
1195 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1198 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1200 _data_structure[i].resize(meshNamesLoc.size());
1201 for(std::size_t j=0;j<splitByMeshName.size();j++)
1203 std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsp;
1204 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1205 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1206 sbmn[k]=splitByMeshName[j][k];
1207 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1208 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1209 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1210 this->_data_structure[i][j].resize(commonSupSplit.size());
1211 for(std::size_t k=0;k<commonSupSplit.size();k++)
1213 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1214 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1216 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1217 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1220 for(std::size_t k=0;k<commonSupSplit.size();k++)
1221 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1224 this->removeEmptyLeaves();
1226 this->computeFullNameInLeaves();
1229 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1231 MCAuto<MEDFileMeshes> ms;
1232 MCAuto<MEDFileFields> fields;
1235 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1237 ms=MEDFileMeshes::New(fileName);
1238 fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1242 #ifdef MEDREADER_USE_MPI
1243 ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1244 int nbMeshes(ms->getNumberOfMeshes());
1245 for(int i=0;i<nbMeshes;i++)
1247 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1248 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1250 MCAuto<DataArrayInt> tmp3(tmp2->zipCoords());
1252 fields=MEDFileFields::LoadPartOf(fileName,false,ms);//false is important to not read the values
1254 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1255 throw INTERP_KERNEL::Exception(oss.str().c_str());
1261 MCAuto<MEDCoupling::SauvReader> sr(MEDCoupling::SauvReader::New(fileName));
1262 MCAuto<MEDCoupling::MEDFileData> mfd(sr->loadInMEDFileDS());
1263 ms=mfd->getMeshes(); ms->incrRef();
1264 int nbMeshes(ms->getNumberOfMeshes());
1265 for(int i=0;i<nbMeshes;i++)
1267 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1268 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1270 tmp2->forceComputationOfParts();
1272 fields=mfd->getFields();
1273 if(fields.isNotNull())
1276 loadInMemory(fields,ms);
1279 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1281 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1282 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1284 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1285 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1287 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1288 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1290 newSD1.push_back(*it2);
1292 newSD0.push_back(newSD1);
1295 newSD.push_back(newSD0);
1299 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1301 if(compInfos.size()!=1)
1303 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1306 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1308 return _data_structure[0][0][0].getMeshName();
1311 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1314 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1315 return leaf.getTimeSteps(tk);
1318 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1321 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1322 meshName=leaf.getMeshName();
1323 std::vector<double> ts(leaf.getTimeSteps(tk));
1324 std::size_t zeTimeId(0);
1327 std::vector<double> ts2(ts.size());
1328 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1329 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1330 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1333 if(zeTimeId==(int)ts.size())
1334 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1335 if(zeTimeId==(int)ts.size())
1336 {//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.
1337 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1339 double valAttachedToPos(-std::numeric_limits<double>::max());
1340 for(std::size_t i=0;i<ts.size();i++)
1344 if(ts[i]>valAttachedToPos)
1347 valAttachedToPos=ts[i];
1352 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1353 valAttachedToPos=std::numeric_limits<double>::max();
1354 for(std::size_t i=0;i<ts.size();i++)
1356 if(ts[i]<valAttachedToPos)
1359 valAttachedToPos=ts[i];
1364 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1365 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1366 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1367 std::cerr << oss.str() << std::endl;
1371 tr=new MEDStdTimeReq((int)zeTimeId);
1373 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1374 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1379 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1381 int nbOfActivated(0);
1382 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1383 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1384 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1385 if((*it2).isActivated())
1387 if(nbOfActivated!=1)
1389 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1390 throw INTERP_KERNEL::Exception(oss.str().c_str());
1392 int i0(0),i1(0),i2(0);
1393 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1394 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1395 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1396 if((*it2).isActivated())
1398 lev0=i0; lev1=i1; lev2=i2;
1401 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1404 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1407 os << "#############################################" << std::endl;
1408 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1411 os << "TS" << i << std::endl;
1412 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1415 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1418 os << " " << (*it2).getMeshName() << std::endl;
1419 os << " Comp" << k << std::endl;
1420 (*it2).printMySelf(os);
1424 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1427 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1429 std::map<std::string,bool> ret;
1430 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1431 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1432 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1433 (*it2).dumpState(ret);
1437 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms, MEDCoupling::MEDFileFields *ret)
1440 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1441 for(int i=0;i<ms->getNumberOfMeshes();i++)
1443 MEDFileMesh *mm(ms->getMeshAtPos(i));
1444 std::vector<int> levs(mm->getNonEmptyLevels());
1445 MEDCoupling::MCAuto<MEDCoupling::MEDFileField1TS> f1tsMultiLev(MEDCoupling::MEDFileField1TS::New());
1446 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1449 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1451 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1452 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1454 MEDCoupling::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1455 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1457 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1458 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1461 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1462 f1tsMultiLev->setFieldNoProfileSBT(f);
1467 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1468 if(levsExt.size()==levs.size()+1)
1470 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1471 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_NODES));
1473 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1474 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1475 arr->iota(); f->setArray(arr);
1476 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1477 f1tsMultiLev->setFieldNoProfileSBT(f);
1485 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1486 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1488 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1489 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1492 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1493 f1tsMultiLev->setFieldNoProfileSBT(f);
1496 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> fmtsMultiLev(MEDCoupling::MEDFileFieldMultiTS::New());
1497 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1498 ret->pushField(fmtsMultiLev);
1502 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const MEDCoupling::MEDFileFields *ret)
1504 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1506 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1507 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1508 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1510 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1511 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1512 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1516 MEDCoupling::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms)
1518 MEDCoupling::MCAuto<MEDCoupling::MEDFileFields> ret(MEDCoupling::MEDFileFields::New());
1519 AppendFieldFromMeshes(ms,ret);
1523 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1525 std::vector<std::string> ret;
1527 while(pos!=std::string::npos)
1529 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1530 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1532 pos=fullFieldName.find_first_not_of(sep,curPos);
1538 * Here the non regression tests.
1539 * const char inp0[]="";
1540 * const char exp0[]="";
1541 * const char inp1[]="field";
1542 * const char exp1[]="field";
1543 * const char inp2[]="_________";
1544 * const char exp2[]="_________";
1545 * const char inp3[]="field_p";
1546 * const char exp3[]="field_p";
1547 * const char inp4[]="field__p";
1548 * const char exp4[]="field_p";
1549 * const char inp5[]="field_p__";
1550 * const char exp5[]="field_p";
1551 * const char inp6[]="field_p_";
1552 * const char exp6[]="field_p";
1553 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1554 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1555 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1556 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1557 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1558 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1559 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1560 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1562 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1564 const char SEP('_');
1565 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1567 return fullFieldName;//should never happen
1571 return fullFieldName;
1575 std::string ret(v[0]);
1576 for(std::size_t i=1;i<v.size();i++)
1581 { ret+=SEP; ret+=v[i]; }
1587 return fullFieldName;
1593 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1597 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1602 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1604 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1606 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1612 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1613 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1614 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1616 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1618 std::size_t sz(ts.size());
1619 bool isInHumanRange(true);
1621 for(std::size_t i=0;i<sz;i++)
1624 if(ts[i]<=-1e299 || ts[i]>=1e299)
1625 isInHumanRange=false;
1628 return processedUsingPairOfIds(tsPairs);
1630 return processedUsingPairOfIds(tsPairs);
1631 _postprocessed_time=ts;
1632 return getPostProcessedTime();
1637 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1639 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1641 std::size_t sz(ts.size());
1642 std::vector<double> ts2(sz);
1644 for(std::size_t i=0;i<sz;i++)
1649 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1652 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1654 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1658 void TimeKeeper::printSelf(std::ostream& oss) const
1660 std::size_t sz(_activated_ts.size());
1661 for(std::size_t i=0;i<sz;i++)
1663 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1667 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1669 std::size_t sz(_activated_ts.size());
1670 std::vector<bool> ret(sz);
1671 for(std::size_t i=0;i<sz;i++)
1673 ret[i]=_activated_ts[i].first;
1678 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1680 std::size_t sz(tsPairs.size());
1681 std::set<int> s0,s1;
1682 for(std::size_t i=0;i<sz;i++)
1683 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1686 _postprocessed_time.resize(sz);
1687 for(std::size_t i=0;i<sz;i++)
1688 _postprocessed_time[i]=(double)tsPairs[i].first;
1689 return getPostProcessedTime();
1693 _postprocessed_time.resize(sz);
1694 for(std::size_t i=0;i<sz;i++)
1695 _postprocessed_time[i]=(double)tsPairs[i].second;
1696 return getPostProcessedTime();
1698 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1699 _postprocessed_time.resize(sz);
1700 for(std::size_t i=0;i<sz;i++)
1701 _postprocessed_time[i]=(double)i;
1702 return getPostProcessedTime();
1705 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1707 _activated_ts.resize(maxNumberOfTS);
1708 for(int i=0;i<maxNumberOfTS;i++)
1710 std::ostringstream oss; oss << "000" << i;
1711 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());