1 // Copyright (C) 2010-2015 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 "vtkIdTypeArray.h"
46 #include "vtkDoubleArray.h"
47 #include "vtkIntArray.h"
48 #include "vtkCellArray.h"
49 #include "vtkPointData.h"
50 #include "vtkFieldData.h"
51 #include "vtkCellData.h"
53 #include "vtkMutableDirectedGraph.h"
55 using namespace MEDCoupling;
57 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
59 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
61 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
63 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
65 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
67 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
69 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
71 const char MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME[]="GlobalNodeIds";// WARNING DO NOT CHANGE IT BEFORE HAVING CHECKED IN PV SOURCES !
73 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
75 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
77 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
79 vtkIdTypeArray *ELGACmp::findOrCreate(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds, bool& isNew) const
81 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
90 return createNew(globs,locsReallyUsed,vtkd,ds);
94 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd) const
96 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
97 if(it==_loc_names.end())
99 std::size_t pos(std::distance(_loc_names.begin(),it));
100 vtkIdTypeArray *ret(_elgas[pos]);
101 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
102 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
104 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
106 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
110 vtkIdTypeArray *ELGACmp::createNew(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds) const
112 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
113 std::vector< std::vector<std::string> > locNames(_loc_names);
114 std::vector<vtkIdTypeArray *> elgas(_elgas);
115 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
117 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
118 if(it!=locNames.end())
119 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
120 locNames.push_back(locsReallyUsed);
121 vtkIdTypeArray *elga(vtkIdTypeArray::New());
122 elga->SetNumberOfComponents(1);
123 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
124 std::map<unsigned char,int> m;
125 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
127 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
128 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
129 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
130 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
131 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
132 // 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.
133 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
134 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
135 double *shape(new double[nbPtsPerCell*nbGaussPt]);
136 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
137 calculator.initLocalInfo();
138 const std::vector<double>& wgths(loc.getGaussWeights());
139 for(int i=0;i<nbGaussPt;i++)
141 const double *pt0(calculator.getFunctionValues(i));
142 if(ct!=INTERP_KERNEL::NORM_HEXA27)
143 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
146 for(int j=0;j<27;j++)
147 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
150 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
151 m[vtkType]=nbGaussPt;
152 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
154 key->Set(elga->GetInformation(),def,vtkType);
155 key->Set(vtkd->GetInformation(),def,vtkType);
156 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
159 vtkIdType ncell(ds->GetNumberOfCells());
160 int *pt(new int[ncell]),offset(0);
161 for(vtkIdType cellId=0;cellId<ncell;cellId++)
163 vtkCell *cell(ds->GetCell(cellId));
164 int delta(m[cell->GetCellType()]);
168 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
169 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
170 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
171 std::string ossStr(oss.str());
172 elga->SetName(ossStr.c_str());
173 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
174 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
175 elgas.push_back(elga);
179 _defs.push_back(defs);
183 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
185 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
186 ds->GetCellData()->AddArray(*it);
191 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
193 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
194 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
195 (*it1).first->Delete();
201 class MEDFileVTKTraits
204 typedef void VtkType;
209 class MEDFileVTKTraits<int>
212 typedef vtkIntArray VtkType;
213 typedef MEDCoupling::DataArrayInt MCType;
217 class MEDFileVTKTraits<double>
220 typedef vtkDoubleArray VtkType;
221 typedef MEDCoupling::DataArrayDouble MCType;
225 void AssignDataPointerToVTK(typename MEDFileVTKTraits<T>::VtkType *vtkTab, typename MEDFileVTKTraits<T>::MCType *mcTab, bool noCpyNumNodes)
228 vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),1,vtkDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE);
230 { vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),0,vtkDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE); mcTab->accessToMemArray().setSpecificDeallocator(0); }
233 // here copy is always assumed.
234 template<class VTKT, class MCT>
235 void AssignDataPointerOther(VTKT *vtkTab, MCT *mcTab, int nbElems)
237 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(mcTab->getPointer()),nbElems,0,VTKT::VTK_DATA_ARRAY_FREE);
238 mcTab->accessToMemArray().setSpecificDeallocator(0);
243 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
247 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>& arr):MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
249 std::vector< std::vector<MEDCoupling::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
251 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
252 if(typs[0].size()!=1)
253 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 !");
254 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
255 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
259 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
261 MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>::operator=(other);
264 _ze_name=other._ze_name;
265 _ze_full_name.clear();
269 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
274 int MEDFileFieldRepresentationLeavesArrays::getId() const
279 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
281 return _ze_full_name;
284 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
286 return _ze_full_name.c_str();
289 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
291 vtkIdType refId(sil->AddChild(root,edge));
292 names.push_back(_ze_name);
294 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
296 sil->AddChild(refId,edge);
297 names.push_back(std::string());
301 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
303 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
304 _ze_full_name=oss3.str();
307 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
312 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
314 bool ret(_activated!=status);
319 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshStruct *mst, vtkDataSet *ds) const
321 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
322 tr->setNumberOfTS((operator->())->getNumberOfTS());
324 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
326 MCAuto<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
327 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
328 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
329 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
330 DataArray *crudeArr(0),*postProcessedArr(0);
332 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
334 crudeArr=f1tsPtrInt->getUndergroundDataArray();
336 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
337 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
338 f1ts->loadArraysIfNecessary();
339 MCAuto<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
342 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
344 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
345 vtkFieldData *att(0);
350 att=ds->GetCellData();
355 att=ds->GetPointData();
360 att=ds->GetFieldData();
365 att=ds->GetFieldData();
369 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
373 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
374 vtkDoubleArray *vtkd(vtkDoubleArray::New());
375 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
376 for(int i=0;i<vd->getNumberOfComponents();i++)
377 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
378 AssignDataPointerToVTK<double>(vtkd,vd,postProcessedArr==crudeArr);
379 std::string name(tr->buildName(f1ts->getName()));
380 vtkd->SetName(name.c_str());
383 if(discs[0]==ON_GAUSS_PT)
386 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
388 if(discs[0]==ON_GAUSS_NE)
390 vtkIdTypeArray *elno(vtkIdTypeArray::New());
391 elno->SetNumberOfComponents(1);
392 vtkIdType ncell(ds->GetNumberOfCells());
393 int *pt(new int[ncell]),offset(0);
394 std::set<int> cellTypes;
395 for(vtkIdType cellId=0;cellId<ncell;cellId++)
397 vtkCell *cell(ds->GetCell(cellId));
398 int delta(cell->GetNumberOfPoints());
399 cellTypes.insert(cell->GetCellType());
403 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
404 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
405 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
406 elno->SetName(nameElno.c_str());
407 ds->GetCellData()->AddArray(elno);
408 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
409 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
411 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
412 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
414 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
415 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
417 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
418 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
419 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
420 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
421 double *shape(new double[nbGaussPt*nbGaussPt]);
423 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(ct,dummy));
424 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
425 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
426 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
427 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
428 calculator.initLocalInfo();
429 for(int i=0;i<nbGaussPt;i++)
431 const double *pt0(calculator.getFunctionValues(i));
432 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
434 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
436 key->Set(elno->GetInformation(),def,*it);
437 key->Set(vtkd->GetInformation(),def,*it);
446 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
447 vtkIntArray *vtkd(vtkIntArray::New());
448 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
449 for(int i=0;i<vi->getNumberOfComponents();i++)
450 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
451 AssignDataPointerToVTK<int>(vtkd,vi,postProcessedArr==crudeArr);
452 std::string name(tr->buildName(f1ts->getName()));
453 vtkd->SetName(name.c_str());
458 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
462 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
464 _elga_cmp.appendELGAIfAny(ds);
469 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
473 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS> >& arr,
474 const MEDCoupling::MCAuto<MEDCoupling::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
476 for(std::size_t i=0;i<arr.size();i++)
477 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
480 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
483 _cached_ds->Delete();
486 bool MEDFileFieldRepresentationLeaves::empty() const
488 const MEDFileFastCellSupportComparator *fcscp(_fsp);
489 return fcscp==0 || _arrays.empty();
492 void MEDFileFieldRepresentationLeaves::setId(int& id) const
494 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
498 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
500 return _arrays[0]->getMeshName();
503 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
505 return (int)_arrays.size();
508 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
510 return _arrays[0]->getNumberOfTS();
513 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
515 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
516 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
520 * \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.
522 void MEDFileFieldRepresentationLeaves::feedSIL(const MEDCoupling::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
524 vtkIdType root2(sil->AddChild(root,edge));
525 names.push_back(std::string("Arrs"));
526 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
527 (*it).feedSIL(sil,root2,edge,names);
529 vtkIdType root3(sil->AddChild(root,edge));
530 names.push_back(std::string("InfoOnGeoType"));
531 const MEDCoupling::MEDFileMesh *m(0);
533 m=ms->getMeshWithName(meshName);
534 const MEDCoupling::MEDFileFastCellSupportComparator *fsp(_fsp);
535 if(!fsp || fsp->getNumberOfTS()==0)
537 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
538 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
540 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
541 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
542 sil->AddChild(root3,edge);
543 names.push_back(cmStr);
547 bool MEDFileFieldRepresentationLeaves::containId(int id) const
549 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
550 if((*it).getId()==id)
555 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
557 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
558 if((*it).getZeName()==name)
566 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
568 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
569 status[(*it).getZeName()]=(*it).getStatus();
572 bool MEDFileFieldRepresentationLeaves::isActivated() const
574 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
575 if((*it).getStatus())
580 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
582 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
584 os << " - " << (*it0).getZeName() << " (";
585 if((*it0).getStatus())
589 os << ")" << std::endl;
593 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
595 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
596 (*it).setStatus(true);
599 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
601 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
602 if((*it).getId()==id)
604 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
607 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
610 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
611 std::vector<double> ret;
612 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
613 return tk.getTimeStepsRegardingPolicy(dtits,ret);
616 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
619 return _arrays[0]->getTimeSteps(ts);
623 return std::vector< std::pair<int,int> >();
627 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
629 std::ostringstream oss;
630 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
631 std::vector<double> ret1;
632 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
633 std::size_t sz(ret1.size());
634 for(std::size_t i=0;i<sz;i++)
636 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
639 std::string tmp(oss.str());
640 if(tmp.size()>200 && i!=sz-1)
650 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshes *meshes, vtkDataSet *ds) const
653 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
654 MCAuto<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
655 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
656 if((*it).getStatus())
658 (*it).appendFields(tr,globs,mml,mst,ds);
659 (*it).appendELGAIfAny(ds);
663 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
665 DataArrayDouble *coordsMC(0);
666 DataArrayByte *typesMC(0);
667 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
668 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
669 MCAuto<DataArrayDouble> coordsSafe(coordsMC);
670 MCAuto<DataArrayByte> typesSafe(typesMC);
671 MCAuto<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
673 int nbOfCells(typesSafe->getNbOfElems());
674 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
675 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
676 AssignDataPointerOther<vtkUnsignedCharArray,DataArrayByte>(cellTypes,typesSafe,nbOfCells);
677 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
678 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(cellLocations,cellLocationsSafe,nbOfCells);
679 vtkCellArray *cells(vtkCellArray::New());
680 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
681 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(cells2,cellsSafe,cellsSafe->getNbOfElems());
682 cells->SetCells(nbOfCells,cells2);
684 if(faceLocationsMC!=0 && facesMC!=0)
686 vtkIdTypeArray *faces(vtkIdTypeArray::New());
687 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(faces,facesSafe,facesSafe->getNbOfElems());
688 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
689 AssignDataPointerOther<vtkIdTypeArray,DataArrayInt>(faceLocations,faceLocationsSafe,faceLocationsSafe->getNbOfElems());
690 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
691 faceLocations->Delete();
695 ret->SetCells(cellTypes,cellLocations,cells);
697 cellLocations->Delete();
699 vtkPoints *pts(vtkPoints::New());
700 vtkDoubleArray *pts2(vtkDoubleArray::New());
701 pts2->SetNumberOfComponents(3);
702 AssignDataPointerToVTK<double>(pts2,coordsSafe,statusOfCoords);
711 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(MEDCoupling::MEDCMeshMultiLev *mm) const
714 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
715 vtkDoubleArray *vtkTmp(0);
716 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
717 std::size_t dim(arrs.size());
719 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
720 int sizePerAxe[3]={1,1,1};
721 sizePerAxe[0]=arrs[0]->getNbOfElems();
723 sizePerAxe[1]=arrs[1]->getNbOfElems();
725 sizePerAxe[2]=arrs[2]->getNbOfElems();
726 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
727 vtkTmp=vtkDoubleArray::New();
728 vtkTmp->SetNumberOfComponents(1);
729 AssignDataPointerToVTK<double>(vtkTmp,arrs[0],isInternal);
730 ret->SetXCoordinates(vtkTmp);
735 vtkTmp=vtkDoubleArray::New();
736 vtkTmp->SetNumberOfComponents(1);
737 AssignDataPointerToVTK<double>(vtkTmp,arrs[1],isInternal);
738 ret->SetYCoordinates(vtkTmp);
744 vtkTmp=vtkDoubleArray::New();
745 vtkTmp->SetNumberOfComponents(1);
746 AssignDataPointerToVTK<double>(vtkTmp,arrs[2],isInternal);
747 ret->SetZCoordinates(vtkTmp);
754 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(MEDCoupling::MEDCurveLinearMeshMultiLev *mm) const
756 int meshStr[3]={1,1,1};
757 DataArrayDouble *coords(0);
758 std::vector<int> nodeStrct;
760 mm->buildVTUArrays(coords,nodeStrct,isInternal);
761 std::size_t dim(nodeStrct.size());
763 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
764 meshStr[0]=nodeStrct[0];
766 meshStr[1]=nodeStrct[1];
768 meshStr[2]=nodeStrct[2];
769 vtkStructuredGrid *ret(vtkStructuredGrid::New());
770 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
771 vtkDoubleArray *da(vtkDoubleArray::New());
772 da->SetNumberOfComponents(3);
773 if(coords->getNumberOfComponents()==3)
774 AssignDataPointerToVTK<double>(da,coords,isInternal);//if isIntenal==True VTK has not the ownership of double * because MEDLoader main struct has it !
777 MCAuto<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
778 AssignDataPointerToVTK<double>(da,coords2,false);//let VTK deal with double *
781 vtkPoints *points=vtkPoints::New();
782 ret->SetPoints(points);
789 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDFileMeshes *meshes) const
792 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
793 MCAuto<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
794 MCAuto<MEDMeshMultiLev> mml2(mml->prepare());
795 MEDMeshMultiLev *ptMML2(mml2);
798 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
799 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
800 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
804 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
808 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
812 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
815 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
816 _cached_ds=ret->NewInstance();
817 _cached_ds->ShallowCopy(ret);
821 ret=_cached_ds->NewInstance();
822 ret->ShallowCopy(_cached_ds);
825 appendFields(tr,globs,mml,meshes,ret);
826 // The arrays links to mesh
827 DataArrayInt *famCells(0),*numCells(0);
828 bool noCpyFamCells(false),noCpyNumCells(false);
829 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
832 vtkIntArray *vtkTab(vtkIntArray::New());
833 vtkTab->SetNumberOfComponents(1);
834 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
835 AssignDataPointerToVTK<int>(vtkTab,famCells,noCpyFamCells);
836 ret->GetCellData()->AddArray(vtkTab);
840 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
843 vtkIntArray *vtkTab(vtkIntArray::New());
844 vtkTab->SetNumberOfComponents(1);
845 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
846 AssignDataPointerToVTK<int>(vtkTab,numCells,noCpyNumCells);
847 ret->GetCellData()->AddArray(vtkTab);
851 // The arrays links to mesh
852 DataArrayInt *famNodes(0),*numNodes(0);
853 bool noCpyFamNodes(false),noCpyNumNodes(false);
854 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
857 vtkIntArray *vtkTab(vtkIntArray::New());
858 vtkTab->SetNumberOfComponents(1);
859 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
860 AssignDataPointerToVTK<int>(vtkTab,famNodes,noCpyFamNodes);
861 ret->GetPointData()->AddArray(vtkTab);
865 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
868 vtkIntArray *vtkTab(vtkIntArray::New());
869 vtkTab->SetNumberOfComponents(1);
870 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
871 AssignDataPointerToVTK<int>(vtkTab,numNodes,noCpyNumNodes);
872 ret->GetPointData()->AddArray(vtkTab);
876 // Global Node Ids if any ! (In // mode)
877 DataArrayInt *gni(ptMML2->retrieveGlobalNodeIdsIfAny());
880 vtkIntArray *vtkTab(vtkIntArray::New());
881 vtkTab->SetNumberOfComponents(1);
882 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME);
883 AssignDataPointerToVTK<int>(vtkTab,gni,false);
884 ret->GetPointData()->AddArray(vtkTab);
891 //////////////////////
893 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
897 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
900 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
901 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
902 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
903 ret+=(*it2).getNumberOfArrays();
907 void MEDFileFieldRepresentationTree::assignIds() const
910 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
911 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
912 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
916 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
918 std::size_t it0Cnt(0);
919 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
921 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
922 std::string tsName(oss.str());
923 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
925 std::string meshName((*it1)[0].getMeshName());
926 std::size_t it2Cnt(0);
927 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
929 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
930 std::string comSupStr(oss2.str());
931 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
937 void MEDFileFieldRepresentationTree::activateTheFirst() const
939 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
940 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
941 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
943 (*it2).activateAllArrays();
948 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
950 std::size_t it0Cnt(0);
951 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
953 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
954 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
956 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
957 std::vector<double> ts;
958 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
959 std::size_t nbOfTS(dtits.size());
960 std::ostringstream oss3; oss3 << nbOfTS;
961 names.push_back(oss3.str());
962 for(std::size_t i=0;i<nbOfTS;i++)
964 std::ostringstream oss4; oss4 << dtits[i].first;
965 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
966 names.push_back(oss4.str());
967 std::ostringstream oss5; oss5 << dtits[i].second;
968 vtkIdType ItId(sil->AddChild(DtId,edge));
969 names.push_back(oss5.str());
970 std::ostringstream oss6; oss6 << ts[i];
971 sil->AddChild(ItId,edge);
972 names.push_back(oss6.str());
975 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
976 std::string tsName(oss.str());
977 vtkIdType typeId0(sil->AddChild(root,edge));
978 names.push_back(tsName);
979 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
981 std::string meshName((*it1)[0].getMeshName());
982 vtkIdType typeId1(sil->AddChild(typeId0,edge));
983 names.push_back(meshName);
984 std::size_t it2Cnt(0);
985 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
987 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
988 std::string comSupStr(oss2.str());
989 vtkIdType typeId2(sil->AddChild(typeId1,edge));
990 names.push_back(comSupStr);
991 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
997 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
999 int dummy0(0),dummy1(0),dummy2(0);
1000 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1001 std::string ret(leaf.getMeshName());
1004 for(;i<_ms->getNumberOfMeshes();i++)
1006 m=_ms->getMeshAtPos(i);
1007 if(m->getName()==ret)
1010 if(i==_ms->getNumberOfMeshes())
1011 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1012 vtkIdType typeId0(sil->AddChild(root,edge));
1013 names.push_back(m->getName());
1015 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1016 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1017 std::vector<std::string> grps(m->getGroupsNames());
1018 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1020 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1021 names.push_back(*it0);
1022 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1023 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1025 sil->AddChild(typeId2,edge);
1026 names.push_back((*it1).c_str());
1030 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1031 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1032 std::vector<std::string> fams(m->getFamiliesNames());
1033 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1035 sil->AddChild(typeId11,edge);
1036 int famId(m->getFamilyId((*it00).c_str()));
1037 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1038 names.push_back(oss.str());
1043 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1045 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1046 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1047 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1048 if((*it2).containId(id))
1049 return (*it2).getLeafArr(id);
1050 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1053 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1055 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1056 return elt.getZeName();
1059 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1061 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1062 return elt.getZeNameC();
1065 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1067 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1068 return elt.getStatus();
1071 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1074 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1075 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1076 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1077 if((*it2).containZeName(name,ret))
1079 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1080 throw INTERP_KERNEL::Exception(msg.str().c_str());
1083 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1085 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1086 bool ret(elt.setStatus(status));//to be implemented
1090 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1093 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1094 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1095 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1096 ret=std::max(ret,(*it2).getNumberOfTS());
1103 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1107 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1109 _ms=MEDFileMeshes::New(fileName);
1110 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1114 #ifdef MEDREADER_USE_MPI
1115 _ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1116 int nbMeshes(_ms->getNumberOfMeshes());
1117 for(int i=0;i<nbMeshes;i++)
1119 MEDCoupling::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1120 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1122 MCAuto<DataArrayInt> tmp3(tmp2->zipCoords());
1124 _fields=MEDFileFields::LoadPartOf(fileName,false,_ms);//false is important to not read the values
1126 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1127 throw INTERP_KERNEL::Exception(oss.str().c_str());
1133 MCAuto<MEDCoupling::SauvReader> sr(MEDCoupling::SauvReader::New(fileName));
1134 MCAuto<MEDCoupling::MEDFileData> mfd(sr->loadInMEDFileDS());
1135 _ms=mfd->getMeshes(); _ms->incrRef();
1136 int nbMeshes(_ms->getNumberOfMeshes());
1137 for(int i=0;i<nbMeshes;i++)
1139 MEDCoupling::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1140 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1142 tmp2->forceComputationOfParts();
1144 _fields=mfd->getFields();
1145 if((MEDCoupling::MEDFileFields *)_fields)
1148 if(!((MEDCoupling::MEDFileFields *)_fields))
1150 _fields=BuildFieldFromMeshes(_ms);
1154 AppendFieldFromMeshes(_ms,_fields);
1156 _ms->cartesianizeMe();
1157 _fields->removeFieldsWithoutAnyTimeStep();
1158 std::vector<std::string> meshNames(_ms->getMeshesNames());
1159 std::vector< MCAuto<MEDFileFields> > fields_per_mesh(meshNames.size());
1160 for(std::size_t i=0;i<meshNames.size();i++)
1162 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1164 std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1166 for(std::vector< MCAuto<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1168 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1170 MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1171 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1173 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1175 if(!(*it)->presenceOfMultiDiscPerGeoType())
1176 allFMTSLeavesToDisplaySafe.push_back(*it);
1178 {// The case of some parts of field have more than one discretization per geo type.
1179 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1180 std::size_t it0Cnt(0);
1181 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1183 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1184 (*it0)->setName(oss.str());
1185 allFMTSLeavesToDisplaySafe.push_back(*it0);
1192 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1193 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1195 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1197 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1198 // memory safety part
1199 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1200 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1202 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1203 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1205 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1206 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1209 // end of memory safety part
1210 // 1st : timesteps, 2nd : meshName, 3rd : common support
1211 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1212 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1214 std::vector< std::string > meshNamesLoc;
1215 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1216 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1218 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1219 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1220 if(it==meshNamesLoc.end())
1222 meshNamesLoc.push_back(meshName);
1223 splitByMeshName.resize(splitByMeshName.size()+1);
1224 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1227 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1229 _data_structure[i].resize(meshNamesLoc.size());
1230 for(std::size_t j=0;j<splitByMeshName.size();j++)
1232 std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsp;
1233 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1234 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1235 sbmn[k]=splitByMeshName[j][k];
1236 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1237 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1238 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1239 this->_data_structure[i][j].resize(commonSupSplit.size());
1240 for(std::size_t k=0;k<commonSupSplit.size();k++)
1242 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1243 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1245 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1246 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1249 for(std::size_t k=0;k<commonSupSplit.size();k++)
1250 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1253 this->removeEmptyLeaves();
1255 this->computeFullNameInLeaves();
1258 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1260 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1261 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1263 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1264 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1266 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1267 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1269 newSD1.push_back(*it2);
1271 newSD0.push_back(newSD1);
1274 newSD.push_back(newSD0);
1278 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1280 if(compInfos.size()!=1)
1282 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1285 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1287 return _data_structure[0][0][0].getMeshName();
1290 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1293 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1294 return leaf.getTimeSteps(tk);
1297 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1300 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1301 meshName=leaf.getMeshName();
1302 std::vector<double> ts(leaf.getTimeSteps(tk));
1303 std::size_t zeTimeId(0);
1306 std::vector<double> ts2(ts.size());
1307 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1308 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1309 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1312 if(zeTimeId==(int)ts.size())
1313 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1314 if(zeTimeId==(int)ts.size())
1315 {//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.
1316 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1318 double valAttachedToPos(-std::numeric_limits<double>::max());
1319 for(std::size_t i=0;i<ts.size();i++)
1323 if(ts[i]>valAttachedToPos)
1326 valAttachedToPos=ts[i];
1331 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1332 valAttachedToPos=std::numeric_limits<double>::max();
1333 for(std::size_t i=0;i<ts.size();i++)
1335 if(ts[i]<valAttachedToPos)
1338 valAttachedToPos=ts[i];
1343 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1344 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1345 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1346 std::cerr << oss.str() << std::endl;
1350 tr=new MEDStdTimeReq((int)zeTimeId);
1352 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1353 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1358 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1360 int nbOfActivated(0);
1361 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1362 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1363 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1364 if((*it2).isActivated())
1366 if(nbOfActivated!=1)
1368 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1369 throw INTERP_KERNEL::Exception(oss.str().c_str());
1371 int i0(0),i1(0),i2(0);
1372 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1373 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1374 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1375 if((*it2).isActivated())
1377 lev0=i0; lev1=i1; lev2=i2;
1380 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1383 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1386 os << "#############################################" << std::endl;
1387 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1390 os << "TS" << i << std::endl;
1391 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1394 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1397 os << " " << (*it2).getMeshName() << std::endl;
1398 os << " Comp" << k << std::endl;
1399 (*it2).printMySelf(os);
1403 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1406 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1408 std::map<std::string,bool> ret;
1409 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1410 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1411 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1412 (*it2).dumpState(ret);
1416 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms, MEDCoupling::MEDFileFields *ret)
1419 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1420 for(int i=0;i<ms->getNumberOfMeshes();i++)
1422 MEDFileMesh *mm(ms->getMeshAtPos(i));
1423 std::vector<int> levs(mm->getNonEmptyLevels());
1424 MEDCoupling::MCAuto<MEDCoupling::MEDFileField1TS> f1tsMultiLev(MEDCoupling::MEDFileField1TS::New());
1425 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1428 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1430 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1431 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1433 MEDCoupling::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1434 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1436 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1437 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1440 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1441 f1tsMultiLev->setFieldNoProfileSBT(f);
1446 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1447 if(levsExt.size()==levs.size()+1)
1449 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1450 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_NODES));
1452 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1453 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1454 arr->iota(); f->setArray(arr);
1455 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1456 f1tsMultiLev->setFieldNoProfileSBT(f);
1464 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1465 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1467 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1468 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1471 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1472 f1tsMultiLev->setFieldNoProfileSBT(f);
1475 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> fmtsMultiLev(MEDCoupling::MEDFileFieldMultiTS::New());
1476 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1477 ret->pushField(fmtsMultiLev);
1481 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const MEDCoupling::MEDFileFields *ret)
1483 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1485 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1486 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1487 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1489 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1490 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1491 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1495 MEDCoupling::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms)
1497 MEDCoupling::MCAuto<MEDCoupling::MEDFileFields> ret(MEDCoupling::MEDFileFields::New());
1498 AppendFieldFromMeshes(ms,ret);
1502 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1504 std::vector<std::string> ret;
1506 while(pos!=std::string::npos)
1508 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1509 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1511 pos=fullFieldName.find_first_not_of(sep,curPos);
1517 * Here the non regression tests.
1518 * const char inp0[]="";
1519 * const char exp0[]="";
1520 * const char inp1[]="field";
1521 * const char exp1[]="field";
1522 * const char inp2[]="_________";
1523 * const char exp2[]="_________";
1524 * const char inp3[]="field_p";
1525 * const char exp3[]="field_p";
1526 * const char inp4[]="field__p";
1527 * const char exp4[]="field_p";
1528 * const char inp5[]="field_p__";
1529 * const char exp5[]="field_p";
1530 * const char inp6[]="field_p_";
1531 * const char exp6[]="field_p";
1532 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1533 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1534 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1535 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1536 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1537 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1538 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1539 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1541 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1543 const char SEP('_');
1544 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1546 return fullFieldName;//should never happen
1550 return fullFieldName;
1554 std::string ret(v[0]);
1555 for(std::size_t i=1;i<v.size();i++)
1560 { ret+=SEP; ret+=v[i]; }
1566 return fullFieldName;
1572 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1576 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1581 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1583 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1585 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1591 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1592 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1593 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1595 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1597 std::size_t sz(ts.size());
1598 bool isInHumanRange(true);
1600 for(std::size_t i=0;i<sz;i++)
1603 if(ts[i]<=-1e299 || ts[i]>=1e299)
1604 isInHumanRange=false;
1607 return processedUsingPairOfIds(tsPairs);
1609 return processedUsingPairOfIds(tsPairs);
1610 _postprocessed_time=ts;
1611 return getPostProcessedTime();
1616 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1618 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1620 std::size_t sz(ts.size());
1621 std::vector<double> ts2(sz);
1623 for(std::size_t i=0;i<sz;i++)
1628 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1631 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1633 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1637 void TimeKeeper::printSelf(std::ostream& oss) const
1639 std::size_t sz(_activated_ts.size());
1640 for(std::size_t i=0;i<sz;i++)
1642 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1646 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1648 std::size_t sz(_activated_ts.size());
1649 std::vector<bool> ret(sz);
1650 for(std::size_t i=0;i<sz;i++)
1652 ret[i]=_activated_ts[i].first;
1657 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1659 std::size_t sz(tsPairs.size());
1660 std::set<int> s0,s1;
1661 for(std::size_t i=0;i<sz;i++)
1662 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1665 _postprocessed_time.resize(sz);
1666 for(std::size_t i=0;i<sz;i++)
1667 _postprocessed_time[i]=(double)tsPairs[i].first;
1668 return getPostProcessedTime();
1672 _postprocessed_time.resize(sz);
1673 for(std::size_t i=0;i<sz;i++)
1674 _postprocessed_time[i]=(double)tsPairs[i].second;
1675 return getPostProcessedTime();
1677 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1678 _postprocessed_time.resize(sz);
1679 for(std::size_t i=0;i<sz;i++)
1680 _postprocessed_time[i]=(double)i;
1681 return getPostProcessedTime();
1684 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1686 _activated_ts.resize(maxNumberOfTS);
1687 for(int i=0;i<maxNumberOfTS;i++)
1689 std::ostringstream oss; oss << "000" << i;
1690 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());