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);
182 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
184 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
185 ds->GetCellData()->AddArray(*it);
190 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
192 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
193 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
194 (*it1).first->Delete();
199 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
203 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>& arr):MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
205 std::vector< std::vector<MEDCoupling::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
207 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
208 if(typs[0].size()!=1)
209 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field spitting by spatial discretization has failed ! A CRASH will occur soon !");
210 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
211 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
215 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
217 MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>::operator=(other);
220 _ze_name=other._ze_name;
221 _ze_full_name.clear();
225 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
230 int MEDFileFieldRepresentationLeavesArrays::getId() const
235 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
237 return _ze_full_name;
240 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
242 return _ze_full_name.c_str();
245 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
247 vtkIdType refId(sil->AddChild(root,edge));
248 names.push_back(_ze_name);
250 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
252 sil->AddChild(refId,edge);
253 names.push_back(std::string());
257 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
259 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
260 _ze_full_name=oss3.str();
263 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
268 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
270 bool ret(_activated!=status);
275 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshStruct *mst, vtkDataSet *ds) const
277 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
278 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
279 tr->setNumberOfTS((operator->())->getNumberOfTS());
281 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
283 MCAuto<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
284 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
285 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
286 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
287 DataArray *crudeArr(0),*postProcessedArr(0);
289 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
291 crudeArr=f1tsPtrInt->getUndergroundDataArray();
293 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
294 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
295 f1ts->loadArraysIfNecessary();
296 MCAuto<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
299 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
301 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
302 vtkFieldData *att(0);
307 att=ds->GetCellData();
312 att=ds->GetPointData();
317 att=ds->GetFieldData();
322 att=ds->GetFieldData();
326 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
330 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
331 vtkDoubleArray *vtkd(vtkDoubleArray::New());
332 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
333 for(int i=0;i<vd->getNumberOfComponents();i++)
334 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
335 if(postProcessedArr!=crudeArr)
337 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vd->accessToMemArray().setSpecificDeallocator(0);
341 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
343 std::string name(tr->buildName(f1ts->getName()));
344 vtkd->SetName(name.c_str());
347 if(discs[0]==ON_GAUSS_PT)
350 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
352 if(discs[0]==ON_GAUSS_NE)
354 vtkIdTypeArray *elno(vtkIdTypeArray::New());
355 elno->SetNumberOfComponents(1);
356 vtkIdType ncell(ds->GetNumberOfCells());
357 int *pt(new int[ncell]),offset(0);
358 std::set<int> cellTypes;
359 for(vtkIdType cellId=0;cellId<ncell;cellId++)
361 vtkCell *cell(ds->GetCell(cellId));
362 int delta(cell->GetNumberOfPoints());
363 cellTypes.insert(cell->GetCellType());
367 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
368 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
369 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
370 elno->SetName(nameElno.c_str());
371 ds->GetCellData()->AddArray(elno);
372 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
373 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
375 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
376 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
378 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
379 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
381 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
382 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
383 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
384 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
385 double *shape(new double[nbGaussPt*nbGaussPt]);
387 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(ct,dummy));
388 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
389 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
390 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
391 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
392 calculator.initLocalInfo();
393 for(int i=0;i<nbGaussPt;i++)
395 const double *pt0(calculator.getFunctionValues(i));
396 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
398 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
400 key->Set(elno->GetInformation(),def,*it);
401 key->Set(vtkd->GetInformation(),def,*it);
410 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
411 vtkIntArray *vtkd(vtkIntArray::New());
412 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
413 for(int i=0;i<vi->getNumberOfComponents();i++)
414 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
415 if(postProcessedArr!=crudeArr)
417 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vi->accessToMemArray().setSpecificDeallocator(0);
421 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
423 std::string name(tr->buildName(f1ts->getName()));
424 vtkd->SetName(name.c_str());
429 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
433 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
435 _elga_cmp.appendELGAIfAny(ds);
440 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
444 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS> >& arr,
445 const MEDCoupling::MCAuto<MEDCoupling::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
447 for(std::size_t i=0;i<arr.size();i++)
448 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
451 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
454 _cached_ds->Delete();
457 bool MEDFileFieldRepresentationLeaves::empty() const
459 const MEDFileFastCellSupportComparator *fcscp(_fsp);
460 return fcscp==0 || _arrays.empty();
463 void MEDFileFieldRepresentationLeaves::setId(int& id) const
465 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
469 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
471 return _arrays[0]->getMeshName();
474 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
476 return (int)_arrays.size();
479 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
481 return _arrays[0]->getNumberOfTS();
484 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
486 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
487 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
491 * \param [in] ms is the meshes pointer. It can be used only for information of geometric types. No special processing will be requested on ms.
493 void MEDFileFieldRepresentationLeaves::feedSIL(const MEDCoupling::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
495 vtkIdType root2(sil->AddChild(root,edge));
496 names.push_back(std::string("Arrs"));
497 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
498 (*it).feedSIL(sil,root2,edge,names);
500 vtkIdType root3(sil->AddChild(root,edge));
501 names.push_back(std::string("InfoOnGeoType"));
502 const MEDCoupling::MEDFileMesh *m(0);
504 m=ms->getMeshWithName(meshName);
505 const MEDCoupling::MEDFileFastCellSupportComparator *fsp(_fsp);
506 if(!fsp || fsp->getNumberOfTS()==0)
508 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
509 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
511 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
512 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
513 sil->AddChild(root3,edge);
514 names.push_back(cmStr);
518 bool MEDFileFieldRepresentationLeaves::containId(int id) const
520 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
521 if((*it).getId()==id)
526 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
528 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
529 if((*it).getZeName()==name)
537 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
539 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
540 status[(*it).getZeName()]=(*it).getStatus();
543 bool MEDFileFieldRepresentationLeaves::isActivated() const
545 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
546 if((*it).getStatus())
551 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
553 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
555 os << " - " << (*it0).getZeName() << " (";
556 if((*it0).getStatus())
560 os << ")" << std::endl;
564 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
566 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
567 (*it).setStatus(true);
570 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
572 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
573 if((*it).getId()==id)
575 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
578 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
581 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
582 std::vector<double> ret;
583 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
584 return tk.getTimeStepsRegardingPolicy(dtits,ret);
587 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
590 return _arrays[0]->getTimeSteps(ts);
594 return std::vector< std::pair<int,int> >();
598 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
600 std::ostringstream oss;
601 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
602 std::vector<double> ret1;
603 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
604 std::size_t sz(ret1.size());
605 for(std::size_t i=0;i<sz;i++)
607 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
610 std::string tmp(oss.str());
611 if(tmp.size()>200 && i!=sz-1)
621 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshes *meshes, vtkDataSet *ds) const
624 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
625 MCAuto<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
626 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
627 if((*it).getStatus())
629 (*it).appendFields(tr,globs,mml,mst,ds);
630 (*it).appendELGAIfAny(ds);
634 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
636 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
637 DataArrayDouble *coordsMC(0);
638 DataArrayByte *typesMC(0);
639 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
640 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
641 MCAuto<DataArrayDouble> coordsSafe(coordsMC);
642 MCAuto<DataArrayByte> typesSafe(typesMC);
643 MCAuto<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
645 int nbOfCells(typesSafe->getNbOfElems());
646 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
647 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
648 cellTypes->SetArray(reinterpret_cast<unsigned char *>(typesSafe->getPointer()),nbOfCells,0,VTK_DATA_ARRAY_FREE); typesSafe->accessToMemArray().setSpecificDeallocator(0);
649 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
650 cellLocations->SetVoidArray(cellLocationsSafe->getPointer(),nbOfCells,0,VTK_DATA_ARRAY_FREE); cellLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
651 vtkCellArray *cells(vtkCellArray::New());
652 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
653 cells2->SetVoidArray(cellsSafe->getPointer(),cellsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); cellsSafe->accessToMemArray().setSpecificDeallocator(0);
654 cells->SetCells(nbOfCells,cells2);
656 if(faceLocationsMC!=0 && facesMC!=0)
658 vtkIdTypeArray *faces(vtkIdTypeArray::New());
659 faces->SetVoidArray(facesSafe->getPointer(),facesSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); facesSafe->accessToMemArray().setSpecificDeallocator(0);
660 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
661 faceLocations->SetVoidArray(faceLocationsSafe->getPointer(),faceLocationsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); faceLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
662 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
663 faceLocations->Delete();
667 ret->SetCells(cellTypes,cellLocations,cells);
669 cellLocations->Delete();
671 vtkPoints *pts(vtkPoints::New());
672 vtkDoubleArray *pts2(vtkDoubleArray::New());
673 pts2->SetNumberOfComponents(3);
676 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
677 coordsSafe->accessToMemArray().setSpecificDeallocator(0);
680 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
689 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(MEDCoupling::MEDCMeshMultiLev *mm) const
691 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
693 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
694 vtkDoubleArray *vtkTmp(0);
695 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
696 std::size_t dim(arrs.size());
698 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
699 int sizePerAxe[3]={1,1,1};
700 sizePerAxe[0]=arrs[0]->getNbOfElems();
702 sizePerAxe[1]=arrs[1]->getNbOfElems();
704 sizePerAxe[2]=arrs[2]->getNbOfElems();
705 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
706 vtkTmp=vtkDoubleArray::New();
707 vtkTmp->SetNumberOfComponents(1);
709 vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
711 { vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[0]->accessToMemArray().setSpecificDeallocator(0); }
712 ret->SetXCoordinates(vtkTmp);
717 vtkTmp=vtkDoubleArray::New();
718 vtkTmp->SetNumberOfComponents(1);
720 vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
722 { vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[1]->accessToMemArray().setSpecificDeallocator(0); }
723 ret->SetYCoordinates(vtkTmp);
729 vtkTmp=vtkDoubleArray::New();
730 vtkTmp->SetNumberOfComponents(1);
732 vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
734 { vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[2]->accessToMemArray().setSpecificDeallocator(0); }
735 ret->SetZCoordinates(vtkTmp);
742 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(MEDCoupling::MEDCurveLinearMeshMultiLev *mm) const
744 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
745 int meshStr[3]={1,1,1};
746 DataArrayDouble *coords(0);
747 std::vector<int> nodeStrct;
749 mm->buildVTUArrays(coords,nodeStrct,isInternal);
750 std::size_t dim(nodeStrct.size());
752 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
753 meshStr[0]=nodeStrct[0];
755 meshStr[1]=nodeStrct[1];
757 meshStr[2]=nodeStrct[2];
758 vtkStructuredGrid *ret(vtkStructuredGrid::New());
759 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
760 vtkDoubleArray *da(vtkDoubleArray::New());
761 da->SetNumberOfComponents(3);
762 if(coords->getNumberOfComponents()==3)
765 da->SetArray(coords->getPointer(),coords->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);//VTK has not the ownership of double * because MEDLoader main struct has it !
767 { da->SetArray(coords->getPointer(),coords->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); coords->accessToMemArray().setSpecificDeallocator(0); }
771 MCAuto<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
772 da->SetArray(coords2->getPointer(),coords2->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);//let VTK deal with double *
773 coords2->accessToMemArray().setSpecificDeallocator(0);
776 vtkPoints *points=vtkPoints::New();
777 ret->SetPoints(points);
784 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDFileMeshes *meshes) const
786 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
788 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
789 MCAuto<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
790 MCAuto<MEDMeshMultiLev> mml2(mml->prepare());
791 MEDMeshMultiLev *ptMML2(mml2);
794 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
795 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
796 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
800 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
804 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
808 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
811 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
812 _cached_ds=ret->NewInstance();
813 _cached_ds->ShallowCopy(ret);
817 ret=_cached_ds->NewInstance();
818 ret->ShallowCopy(_cached_ds);
821 appendFields(tr,globs,mml,meshes,ret);
822 // The arrays links to mesh
823 DataArrayInt *famCells(0),*numCells(0);
824 bool noCpyFamCells(false),noCpyNumCells(false);
825 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
828 vtkIntArray *vtkTab(vtkIntArray::New());
829 vtkTab->SetNumberOfComponents(1);
830 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
832 vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
834 { vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famCells->accessToMemArray().setSpecificDeallocator(0); }
835 ret->GetCellData()->AddArray(vtkTab);
839 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
842 vtkIntArray *vtkTab(vtkIntArray::New());
843 vtkTab->SetNumberOfComponents(1);
844 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
846 vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
848 { vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numCells->accessToMemArray().setSpecificDeallocator(0); }
849 ret->GetCellData()->AddArray(vtkTab);
853 // The arrays links to mesh
854 DataArrayInt *famNodes(0),*numNodes(0);
855 bool noCpyFamNodes(false),noCpyNumNodes(false);
856 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
859 vtkIntArray *vtkTab(vtkIntArray::New());
860 vtkTab->SetNumberOfComponents(1);
861 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
863 vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
865 { vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famNodes->accessToMemArray().setSpecificDeallocator(0); }
866 ret->GetPointData()->AddArray(vtkTab);
870 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
873 vtkIntArray *vtkTab(vtkIntArray::New());
874 vtkTab->SetNumberOfComponents(1);
875 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
877 vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
879 { vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numNodes->accessToMemArray().setSpecificDeallocator(0); }
880 ret->GetPointData()->AddArray(vtkTab);
884 // Global Node Ids if any ! (In // mode)
885 DataArrayInt *gni(ptMML2->retrieveGlobalNodeIdsIfAny());
888 vtkIntArray *vtkTab(vtkIntArray::New());
889 vtkTab->SetNumberOfComponents(1);
890 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME);
891 vtkTab->SetArray(gni->getPointer(),gni->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
892 gni->accessToMemArray().setSpecificDeallocator(0);
898 //////////////////////
900 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
904 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
907 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
908 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
909 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
910 ret+=(*it2).getNumberOfArrays();
914 void MEDFileFieldRepresentationTree::assignIds() const
917 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
918 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
919 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
923 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
925 std::size_t it0Cnt(0);
926 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
928 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
929 std::string tsName(oss.str());
930 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
932 std::string meshName((*it1)[0].getMeshName());
933 std::size_t it2Cnt(0);
934 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
936 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
937 std::string comSupStr(oss2.str());
938 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
944 void MEDFileFieldRepresentationTree::activateTheFirst() const
946 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
947 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
948 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
950 (*it2).activateAllArrays();
955 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
957 std::size_t it0Cnt(0);
958 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
960 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
961 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
963 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
964 std::vector<double> ts;
965 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
966 std::size_t nbOfTS(dtits.size());
967 std::ostringstream oss3; oss3 << nbOfTS;
968 names.push_back(oss3.str());
969 for(std::size_t i=0;i<nbOfTS;i++)
971 std::ostringstream oss4; oss4 << dtits[i].first;
972 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
973 names.push_back(oss4.str());
974 std::ostringstream oss5; oss5 << dtits[i].second;
975 vtkIdType ItId(sil->AddChild(DtId,edge));
976 names.push_back(oss5.str());
977 std::ostringstream oss6; oss6 << ts[i];
978 sil->AddChild(ItId,edge);
979 names.push_back(oss6.str());
982 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
983 std::string tsName(oss.str());
984 vtkIdType typeId0(sil->AddChild(root,edge));
985 names.push_back(tsName);
986 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
988 std::string meshName((*it1)[0].getMeshName());
989 vtkIdType typeId1(sil->AddChild(typeId0,edge));
990 names.push_back(meshName);
991 std::size_t it2Cnt(0);
992 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
994 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
995 std::string comSupStr(oss2.str());
996 vtkIdType typeId2(sil->AddChild(typeId1,edge));
997 names.push_back(comSupStr);
998 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
1004 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1006 int dummy0(0),dummy1(0),dummy2(0);
1007 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1008 std::string ret(leaf.getMeshName());
1011 for(;i<_ms->getNumberOfMeshes();i++)
1013 m=_ms->getMeshAtPos(i);
1014 if(m->getName()==ret)
1017 if(i==_ms->getNumberOfMeshes())
1018 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1019 vtkIdType typeId0(sil->AddChild(root,edge));
1020 names.push_back(m->getName());
1022 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1023 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1024 std::vector<std::string> grps(m->getGroupsNames());
1025 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1027 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1028 names.push_back(*it0);
1029 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1030 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1032 sil->AddChild(typeId2,edge);
1033 names.push_back((*it1).c_str());
1037 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1038 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1039 std::vector<std::string> fams(m->getFamiliesNames());
1040 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1042 sil->AddChild(typeId11,edge);
1043 int famId(m->getFamilyId((*it00).c_str()));
1044 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1045 names.push_back(oss.str());
1050 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1052 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1053 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1054 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1055 if((*it2).containId(id))
1056 return (*it2).getLeafArr(id);
1057 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1060 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1062 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1063 return elt.getZeName();
1066 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1068 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1069 return elt.getZeNameC();
1072 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1074 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1075 return elt.getStatus();
1078 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1081 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1082 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1083 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1084 if((*it2).containZeName(name,ret))
1086 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1087 throw INTERP_KERNEL::Exception(msg.str().c_str());
1090 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1092 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1093 bool ret(elt.setStatus(status));//to be implemented
1097 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1100 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1101 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1102 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1103 ret=std::max(ret,(*it2).getNumberOfTS());
1110 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1114 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1116 _ms=MEDFileMeshes::New(fileName);
1117 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1121 #ifdef MEDREADER_USE_MPI
1122 _ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1123 int nbMeshes(_ms->getNumberOfMeshes());
1124 for(int i=0;i<nbMeshes;i++)
1126 MEDCoupling::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1127 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1129 MCAuto<DataArrayInt> tmp3(tmp2->zipCoords());
1131 _fields=MEDFileFields::LoadPartOf(fileName,false,_ms);//false is important to not read the values
1133 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1134 throw INTERP_KERNEL::Exception(oss.str().c_str());
1140 MCAuto<MEDCoupling::SauvReader> sr(MEDCoupling::SauvReader::New(fileName));
1141 MCAuto<MEDCoupling::MEDFileData> mfd(sr->loadInMEDFileDS());
1142 _ms=mfd->getMeshes(); _ms->incrRef();
1143 int nbMeshes(_ms->getNumberOfMeshes());
1144 for(int i=0;i<nbMeshes;i++)
1146 MEDCoupling::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1147 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1149 tmp2->forceComputationOfParts();
1151 _fields=mfd->getFields();
1152 if((MEDCoupling::MEDFileFields *)_fields)
1155 if(!((MEDCoupling::MEDFileFields *)_fields))
1157 _fields=BuildFieldFromMeshes(_ms);
1161 AppendFieldFromMeshes(_ms,_fields);
1163 _ms->cartesianizeMe();
1164 _fields->removeFieldsWithoutAnyTimeStep();
1165 std::vector<std::string> meshNames(_ms->getMeshesNames());
1166 std::vector< MCAuto<MEDFileFields> > fields_per_mesh(meshNames.size());
1167 for(std::size_t i=0;i<meshNames.size();i++)
1169 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1171 std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1173 for(std::vector< MCAuto<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1175 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1177 MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1178 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1180 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1182 if(!(*it)->presenceOfMultiDiscPerGeoType())
1183 allFMTSLeavesToDisplaySafe.push_back(*it);
1185 {// The case of some parts of field have more than one discretization per geo type.
1186 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1187 std::size_t it0Cnt(0);
1188 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1190 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1191 (*it0)->setName(oss.str());
1192 allFMTSLeavesToDisplaySafe.push_back(*it0);
1199 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1200 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1202 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1204 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1205 // memory safety part
1206 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1207 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1209 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1210 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1212 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1213 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1216 // end of memory safety part
1217 // 1st : timesteps, 2nd : meshName, 3rd : common support
1218 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1219 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1221 std::vector< std::string > meshNamesLoc;
1222 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1223 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1225 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1226 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1227 if(it==meshNamesLoc.end())
1229 meshNamesLoc.push_back(meshName);
1230 splitByMeshName.resize(splitByMeshName.size()+1);
1231 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1234 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1236 _data_structure[i].resize(meshNamesLoc.size());
1237 for(std::size_t j=0;j<splitByMeshName.size();j++)
1239 std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsp;
1240 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1241 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1242 sbmn[k]=splitByMeshName[j][k];
1243 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1244 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1245 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1246 this->_data_structure[i][j].resize(commonSupSplit.size());
1247 for(std::size_t k=0;k<commonSupSplit.size();k++)
1249 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1250 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1252 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1253 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1256 for(std::size_t k=0;k<commonSupSplit.size();k++)
1257 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1260 this->removeEmptyLeaves();
1262 this->computeFullNameInLeaves();
1265 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1267 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1268 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1270 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1271 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1273 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1274 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1276 newSD1.push_back(*it2);
1278 newSD0.push_back(newSD1);
1281 newSD.push_back(newSD0);
1285 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1287 if(compInfos.size()!=1)
1289 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1292 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1294 return _data_structure[0][0][0].getMeshName();
1297 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1300 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1301 return leaf.getTimeSteps(tk);
1304 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1307 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1308 meshName=leaf.getMeshName();
1309 std::vector<double> ts(leaf.getTimeSteps(tk));
1310 std::size_t zeTimeId(0);
1313 std::vector<double> ts2(ts.size());
1314 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1315 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1316 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1319 if(zeTimeId==(int)ts.size())
1320 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1321 if(zeTimeId==(int)ts.size())
1322 {//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.
1323 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1325 double valAttachedToPos(-std::numeric_limits<double>::max());
1326 for(std::size_t i=0;i<ts.size();i++)
1330 if(ts[i]>valAttachedToPos)
1333 valAttachedToPos=ts[i];
1338 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1339 valAttachedToPos=std::numeric_limits<double>::max();
1340 for(std::size_t i=0;i<ts.size();i++)
1342 if(ts[i]<valAttachedToPos)
1345 valAttachedToPos=ts[i];
1350 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1351 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1352 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1353 std::cerr << oss.str() << std::endl;
1357 tr=new MEDStdTimeReq((int)zeTimeId);
1359 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1360 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1365 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1367 int nbOfActivated(0);
1368 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1369 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1370 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1371 if((*it2).isActivated())
1373 if(nbOfActivated!=1)
1375 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1376 throw INTERP_KERNEL::Exception(oss.str().c_str());
1378 int i0(0),i1(0),i2(0);
1379 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1380 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1381 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1382 if((*it2).isActivated())
1384 lev0=i0; lev1=i1; lev2=i2;
1387 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1390 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1393 os << "#############################################" << std::endl;
1394 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1397 os << "TS" << i << std::endl;
1398 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1401 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1404 os << " " << (*it2).getMeshName() << std::endl;
1405 os << " Comp" << k << std::endl;
1406 (*it2).printMySelf(os);
1410 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1413 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1415 std::map<std::string,bool> ret;
1416 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1417 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1418 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1419 (*it2).dumpState(ret);
1423 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms, MEDCoupling::MEDFileFields *ret)
1426 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1427 for(int i=0;i<ms->getNumberOfMeshes();i++)
1429 MEDFileMesh *mm(ms->getMeshAtPos(i));
1430 std::vector<int> levs(mm->getNonEmptyLevels());
1431 MEDCoupling::MCAuto<MEDCoupling::MEDFileField1TS> f1tsMultiLev(MEDCoupling::MEDFileField1TS::New());
1432 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1435 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1437 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1438 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1440 MEDCoupling::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1441 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1443 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1444 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1447 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1448 f1tsMultiLev->setFieldNoProfileSBT(f);
1453 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1454 if(levsExt.size()==levs.size()+1)
1456 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1457 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_NODES));
1459 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1460 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1461 arr->iota(); f->setArray(arr);
1462 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1463 f1tsMultiLev->setFieldNoProfileSBT(f);
1471 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1472 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1474 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1475 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1478 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1479 f1tsMultiLev->setFieldNoProfileSBT(f);
1482 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> fmtsMultiLev(MEDCoupling::MEDFileFieldMultiTS::New());
1483 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1484 ret->pushField(fmtsMultiLev);
1488 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const MEDCoupling::MEDFileFields *ret)
1490 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1492 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1493 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1494 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1496 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1497 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1498 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1502 MEDCoupling::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms)
1504 MEDCoupling::MCAuto<MEDCoupling::MEDFileFields> ret(MEDCoupling::MEDFileFields::New());
1505 AppendFieldFromMeshes(ms,ret);
1509 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1511 std::vector<std::string> ret;
1513 while(pos!=std::string::npos)
1515 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1516 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1518 pos=fullFieldName.find_first_not_of(sep,curPos);
1524 * Here the non regression tests.
1525 * const char inp0[]="";
1526 * const char exp0[]="";
1527 * const char inp1[]="field";
1528 * const char exp1[]="field";
1529 * const char inp2[]="_________";
1530 * const char exp2[]="_________";
1531 * const char inp3[]="field_p";
1532 * const char exp3[]="field_p";
1533 * const char inp4[]="field__p";
1534 * const char exp4[]="field_p";
1535 * const char inp5[]="field_p__";
1536 * const char exp5[]="field_p";
1537 * const char inp6[]="field_p_";
1538 * const char exp6[]="field_p";
1539 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1540 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1541 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1542 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1543 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1544 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1545 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1546 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1548 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1550 const char SEP('_');
1551 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1553 return fullFieldName;//should never happen
1557 return fullFieldName;
1561 std::string ret(v[0]);
1562 for(std::size_t i=1;i<v.size();i++)
1567 { ret+=SEP; ret+=v[i]; }
1573 return fullFieldName;
1579 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1583 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1588 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1590 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1592 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1598 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1599 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1600 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1602 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1604 std::size_t sz(ts.size());
1605 bool isInHumanRange(true);
1607 for(std::size_t i=0;i<sz;i++)
1610 if(ts[i]<=-1e299 || ts[i]>=1e299)
1611 isInHumanRange=false;
1614 return processedUsingPairOfIds(tsPairs);
1616 return processedUsingPairOfIds(tsPairs);
1617 _postprocessed_time=ts;
1618 return getPostProcessedTime();
1623 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1625 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1627 std::size_t sz(ts.size());
1628 std::vector<double> ts2(sz);
1630 for(std::size_t i=0;i<sz;i++)
1635 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1638 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1640 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1644 void TimeKeeper::printSelf(std::ostream& oss) const
1646 std::size_t sz(_activated_ts.size());
1647 for(std::size_t i=0;i<sz;i++)
1649 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1653 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1655 std::size_t sz(_activated_ts.size());
1656 std::vector<bool> ret(sz);
1657 for(std::size_t i=0;i<sz;i++)
1659 ret[i]=_activated_ts[i].first;
1664 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1666 std::size_t sz(tsPairs.size());
1667 std::set<int> s0,s1;
1668 for(std::size_t i=0;i<sz;i++)
1669 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1672 _postprocessed_time.resize(sz);
1673 for(std::size_t i=0;i<sz;i++)
1674 _postprocessed_time[i]=(double)tsPairs[i].first;
1675 return getPostProcessedTime();
1679 _postprocessed_time.resize(sz);
1680 for(std::size_t i=0;i<sz;i++)
1681 _postprocessed_time[i]=(double)tsPairs[i].second;
1682 return getPostProcessedTime();
1684 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1685 _postprocessed_time.resize(sz);
1686 for(std::size_t i=0;i<sz;i++)
1687 _postprocessed_time[i]=(double)i;
1688 return getPostProcessedTime();
1691 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1693 _activated_ts.resize(maxNumberOfTS);
1694 for(int i=0;i<maxNumberOfTS;i++)
1696 std::ostringstream oss; oss << "000" << i;
1697 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());