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 ParaMEDMEM;
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 MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
73 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
75 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
77 vtkIdTypeArray *ELGACmp::findOrCreate(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds, bool& isNew) const
79 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
88 return createNew(globs,locsReallyUsed,vtkd,ds);
92 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd) const
94 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
95 if(it==_loc_names.end())
97 std::size_t pos(std::distance(_loc_names.begin(),it));
98 vtkIdTypeArray *ret(_elgas[pos]);
99 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
100 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
102 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
104 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
108 vtkIdTypeArray *ELGACmp::createNew(const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDoubleArray *vtkd, vtkDataSet *ds) const
110 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
111 std::vector< std::vector<std::string> > locNames(_loc_names);
112 std::vector<vtkIdTypeArray *> elgas(_elgas);
113 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
115 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
116 if(it!=locNames.end())
117 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
118 locNames.push_back(locsReallyUsed);
119 vtkIdTypeArray *elga(vtkIdTypeArray::New());
120 elga->SetNumberOfComponents(1);
121 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
122 std::map<unsigned char,int> m;
123 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
125 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
126 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
127 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
128 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
129 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
130 // 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.
131 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
132 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
133 double *shape(new double[nbPtsPerCell*nbGaussPt]);
134 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
135 calculator.initLocalInfo();
136 const std::vector<double>& wgths(loc.getGaussWeights());
137 for(int i=0;i<nbGaussPt;i++)
139 const double *pt0(calculator.getFunctionValues(i));
140 if(ct!=INTERP_KERNEL::NORM_HEXA27)
141 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
144 for(int j=0;j<27;j++)
145 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
148 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
149 m[vtkType]=nbGaussPt;
150 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
152 key->Set(elga->GetInformation(),def,vtkType);
153 key->Set(vtkd->GetInformation(),def,vtkType);
154 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
157 vtkIdType ncell(ds->GetNumberOfCells());
158 int *pt(new int[ncell]),offset(0);
159 for(vtkIdType cellId=0;cellId<ncell;cellId++)
161 vtkCell *cell(ds->GetCell(cellId));
162 int delta(m[cell->GetCellType()]);
166 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
167 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
168 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
169 std::string ossStr(oss.str());
170 elga->SetName(ossStr.c_str());
171 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
172 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
173 elgas.push_back(elga);
177 _defs.push_back(defs);
181 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
183 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
184 ds->GetCellData()->AddArray(*it);
189 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
191 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
192 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
193 (*it1).first->Delete();
198 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
202 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>& arr):ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
204 std::vector< std::vector<ParaMEDMEM::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
206 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
207 if(typs[0].size()!=1)
208 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 !");
209 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
210 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
214 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
216 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS>::operator=(other);
219 _ze_name=other._ze_name;
220 _ze_full_name.clear();
224 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
229 int MEDFileFieldRepresentationLeavesArrays::getId() const
234 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
236 return _ze_full_name;
239 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
241 return _ze_full_name.c_str();
244 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
246 vtkIdType refId(sil->AddChild(root,edge));
247 names.push_back(_ze_name);
249 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
251 sil->AddChild(refId,edge);
252 names.push_back(std::string());
256 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
258 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
259 _ze_full_name=oss3.str();
262 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
267 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
269 bool ret(_activated!=status);
274 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshStruct *mst, vtkDataSet *ds) const
276 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
277 const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
278 tr->setNumberOfTS((operator->())->getNumberOfTS());
280 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
282 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
283 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
284 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
285 MEDFileIntField1TS *f1tsPtrInt(dynamic_cast<MEDFileIntField1TS *>(f1tsPtr));
286 DataArray *crudeArr(0),*postProcessedArr(0);
288 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
290 crudeArr=f1tsPtrInt->getUndergroundDataArray();
292 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment !");
293 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
294 f1ts->loadArraysIfNecessary();
295 MEDCouplingAutoRefCountObjectPtr<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
298 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
300 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
301 vtkFieldData *att(0);
306 att=ds->GetCellData();
311 att=ds->GetPointData();
316 att=ds->GetFieldData();
321 att=ds->GetFieldData();
325 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
329 DataArray *vPtr(v); DataArrayDouble *vd(static_cast<DataArrayDouble *>(vPtr));
330 vtkDoubleArray *vtkd(vtkDoubleArray::New());
331 vtkd->SetNumberOfComponents(vd->getNumberOfComponents());
332 for(int i=0;i<vd->getNumberOfComponents();i++)
333 vtkd->SetComponentName(i,vd->getInfoOnComponent(i).c_str());
334 if(postProcessedArr!=crudeArr)
336 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vd->accessToMemArray().setSpecificDeallocator(0);
340 vtkd->SetArray(vd->getPointer(),vd->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
342 std::string name(tr->buildName(f1ts->getName()));
343 vtkd->SetName(name.c_str());
346 if(discs[0]==ON_GAUSS_PT)
349 _elga_cmp.findOrCreate(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp);
351 if(discs[0]==ON_GAUSS_NE)
353 vtkIdTypeArray *elno(vtkIdTypeArray::New());
354 elno->SetNumberOfComponents(1);
355 vtkIdType ncell(ds->GetNumberOfCells());
356 int *pt(new int[ncell]),offset(0);
357 std::set<int> cellTypes;
358 for(vtkIdType cellId=0;cellId<ncell;cellId++)
360 vtkCell *cell(ds->GetCell(cellId));
361 int delta(cell->GetNumberOfPoints());
362 cellTypes.insert(cell->GetCellType());
366 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
367 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
368 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
369 elno->SetName(nameElno.c_str());
370 ds->GetCellData()->AddArray(elno);
371 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
372 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
374 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
375 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
377 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
378 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
380 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
381 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
382 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
383 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
384 double *shape(new double[nbGaussPt*nbGaussPt]);
386 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetLocsFromGeometricType(ct,dummy));
387 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
388 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
389 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
390 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
391 calculator.initLocalInfo();
392 for(int i=0;i<nbGaussPt;i++)
394 const double *pt0(calculator.getFunctionValues(i));
395 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
397 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
399 key->Set(elno->GetInformation(),def,*it);
400 key->Set(vtkd->GetInformation(),def,*it);
409 DataArray *vPtr(v); DataArrayInt *vi(static_cast<DataArrayInt *>(vPtr));
410 vtkIntArray *vtkd(vtkIntArray::New());
411 vtkd->SetNumberOfComponents(vi->getNumberOfComponents());
412 for(int i=0;i<vi->getNumberOfComponents();i++)
413 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
414 if(postProcessedArr!=crudeArr)
416 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); vi->accessToMemArray().setSpecificDeallocator(0);
420 vtkd->SetArray(vi->getPointer(),vi->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
422 std::string name(tr->buildName(f1ts->getName()));
423 vtkd->SetName(name.c_str());
428 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
432 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
434 _elga_cmp.appendELGAIfAny(ds);
439 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
443 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileAnyTypeFieldMultiTS> >& arr,
444 const ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
446 for(std::size_t i=0;i<arr.size();i++)
447 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
450 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
453 _cached_ds->Delete();
456 bool MEDFileFieldRepresentationLeaves::empty() const
458 const MEDFileFastCellSupportComparator *fcscp(_fsp);
459 return fcscp==0 || _arrays.empty();
462 void MEDFileFieldRepresentationLeaves::setId(int& id) const
464 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
468 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
470 return _arrays[0]->getMeshName();
473 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
475 return (int)_arrays.size();
478 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
480 return _arrays[0]->getNumberOfTS();
483 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
485 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
486 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
490 * \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.
492 void MEDFileFieldRepresentationLeaves::feedSIL(const ParaMEDMEM::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
494 vtkIdType root2(sil->AddChild(root,edge));
495 names.push_back(std::string("Arrs"));
496 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
497 (*it).feedSIL(sil,root2,edge,names);
499 vtkIdType root3(sil->AddChild(root,edge));
500 names.push_back(std::string("InfoOnGeoType"));
501 const ParaMEDMEM::MEDFileMesh *m(0);
503 m=ms->getMeshWithName(meshName);
504 const ParaMEDMEM::MEDFileFastCellSupportComparator *fsp(_fsp);
505 if(!fsp || fsp->getNumberOfTS()==0)
507 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
508 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
510 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
511 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
512 sil->AddChild(root3,edge);
513 names.push_back(cmStr);
517 bool MEDFileFieldRepresentationLeaves::containId(int id) const
519 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
520 if((*it).getId()==id)
525 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
527 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
528 if((*it).getZeName()==name)
536 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
538 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
539 status[(*it).getZeName()]=(*it).getStatus();
542 bool MEDFileFieldRepresentationLeaves::isActivated() const
544 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
545 if((*it).getStatus())
550 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
552 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
554 os << " - " << (*it0).getZeName() << " (";
555 if((*it0).getStatus())
559 os << ")" << std::endl;
563 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
565 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
566 (*it).setStatus(true);
569 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
571 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
572 if((*it).getId()==id)
574 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
577 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
580 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
581 std::vector<double> ret;
582 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
583 return tk.getTimeStepsRegardingPolicy(dtits,ret);
586 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
589 return _arrays[0]->getTimeSteps(ts);
593 return std::vector< std::pair<int,int> >();
597 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
599 std::ostringstream oss;
600 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
601 std::vector<double> ret1;
602 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
603 std::size_t sz(ret1.size());
604 for(std::size_t i=0;i<sz;i++)
606 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
609 std::string tmp(oss.str());
610 if(tmp.size()>200 && i!=sz-1)
620 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const ParaMEDMEM::MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDMeshMultiLev *mml, const ParaMEDMEM::MEDFileMeshes *meshes, vtkDataSet *ds) const
623 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
624 MEDCouplingAutoRefCountObjectPtr<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
625 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
626 if((*it).getStatus())
628 (*it).appendFields(tr,globs,mml,mst,ds);
629 (*it).appendELGAIfAny(ds);
633 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
635 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
636 DataArrayDouble *coordsMC(0);
637 DataArrayByte *typesMC(0);
638 DataArrayInt *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
639 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
640 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coordsSafe(coordsMC);
641 MEDCouplingAutoRefCountObjectPtr<DataArrayByte> typesSafe(typesMC);
642 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
644 int nbOfCells(typesSafe->getNbOfElems());
645 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
646 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
647 cellTypes->SetArray(reinterpret_cast<unsigned char *>(typesSafe->getPointer()),nbOfCells,0,VTK_DATA_ARRAY_FREE); typesSafe->accessToMemArray().setSpecificDeallocator(0);
648 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
649 cellLocations->SetVoidArray(cellLocationsSafe->getPointer(),nbOfCells,0,VTK_DATA_ARRAY_FREE); cellLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
650 vtkCellArray *cells(vtkCellArray::New());
651 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
652 cells2->SetVoidArray(cellsSafe->getPointer(),cellsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); cellsSafe->accessToMemArray().setSpecificDeallocator(0);
653 cells->SetCells(nbOfCells,cells2);
655 if(faceLocationsMC!=0 && facesMC!=0)
657 vtkIdTypeArray *faces(vtkIdTypeArray::New());
658 faces->SetVoidArray(facesSafe->getPointer(),facesSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); facesSafe->accessToMemArray().setSpecificDeallocator(0);
659 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
660 faceLocations->SetVoidArray(faceLocationsSafe->getPointer(),faceLocationsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); faceLocationsSafe->accessToMemArray().setSpecificDeallocator(0);
661 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
662 faceLocations->Delete();
666 ret->SetCells(cellTypes,cellLocations,cells);
668 cellLocations->Delete();
670 vtkPoints *pts(vtkPoints::New());
671 vtkDoubleArray *pts2(vtkDoubleArray::New());
672 pts2->SetNumberOfComponents(3);
675 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);
676 coordsSafe->accessToMemArray().setSpecificDeallocator(0);
679 pts2->SetArray(coordsSafe->getPointer(),coordsSafe->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
688 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(ParaMEDMEM::MEDCMeshMultiLev *mm) const
690 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
692 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
693 vtkDoubleArray *vtkTmp(0);
694 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
695 std::size_t dim(arrs.size());
697 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
698 int sizePerAxe[3]={1,1,1};
699 sizePerAxe[0]=arrs[0]->getNbOfElems();
701 sizePerAxe[1]=arrs[1]->getNbOfElems();
703 sizePerAxe[2]=arrs[2]->getNbOfElems();
704 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
705 vtkTmp=vtkDoubleArray::New();
706 vtkTmp->SetNumberOfComponents(1);
708 vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
710 { vtkTmp->SetArray(arrs[0]->getPointer(),arrs[0]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[0]->accessToMemArray().setSpecificDeallocator(0); }
711 ret->SetXCoordinates(vtkTmp);
716 vtkTmp=vtkDoubleArray::New();
717 vtkTmp->SetNumberOfComponents(1);
719 vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
721 { vtkTmp->SetArray(arrs[1]->getPointer(),arrs[1]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[1]->accessToMemArray().setSpecificDeallocator(0); }
722 ret->SetYCoordinates(vtkTmp);
728 vtkTmp=vtkDoubleArray::New();
729 vtkTmp->SetNumberOfComponents(1);
731 vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
733 { vtkTmp->SetArray(arrs[2]->getPointer(),arrs[2]->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); arrs[2]->accessToMemArray().setSpecificDeallocator(0); }
734 ret->SetZCoordinates(vtkTmp);
741 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(ParaMEDMEM::MEDCurveLinearMeshMultiLev *mm) const
743 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
744 int meshStr[3]={1,1,1};
745 DataArrayDouble *coords(0);
746 std::vector<int> nodeStrct;
748 mm->buildVTUArrays(coords,nodeStrct,isInternal);
749 std::size_t dim(nodeStrct.size());
751 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
752 meshStr[0]=nodeStrct[0];
754 meshStr[1]=nodeStrct[1];
756 meshStr[2]=nodeStrct[2];
757 vtkStructuredGrid *ret(vtkStructuredGrid::New());
758 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
759 vtkDoubleArray *da(vtkDoubleArray::New());
760 da->SetNumberOfComponents(3);
761 if(coords->getNumberOfComponents()==3)
764 da->SetArray(coords->getPointer(),coords->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);//VTK has not the ownership of double * because MEDLoader main struct has it !
766 { da->SetArray(coords->getPointer(),coords->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); coords->accessToMemArray().setSpecificDeallocator(0); }
770 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
771 da->SetArray(coords2->getPointer(),coords2->getNbOfElems(),0,VTK_DATA_ARRAY_FREE);//let VTK deal with double *
772 coords2->accessToMemArray().setSpecificDeallocator(0);
775 vtkPoints *points=vtkPoints::New();
776 ret->SetPoints(points);
783 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const ParaMEDMEM::MEDFileMeshes *meshes) const
785 const int VTK_DATA_ARRAY_FREE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_FREE;
787 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
788 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
789 MEDCouplingAutoRefCountObjectPtr<MEDMeshMultiLev> mml2(mml->prepare());
790 MEDMeshMultiLev *ptMML2(mml2);
793 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
794 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
795 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
799 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
803 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
807 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
810 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
811 _cached_ds=ret->NewInstance();
812 _cached_ds->ShallowCopy(ret);
816 ret=_cached_ds->NewInstance();
817 ret->ShallowCopy(_cached_ds);
820 appendFields(tr,globs,mml,meshes,ret);
821 // The arrays links to mesh
822 DataArrayInt *famCells(0),*numCells(0);
823 bool noCpyFamCells(false),noCpyNumCells(false);
824 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
827 vtkIntArray *vtkTab(vtkIntArray::New());
828 vtkTab->SetNumberOfComponents(1);
829 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
831 vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
833 { vtkTab->SetArray(famCells->getPointer(),famCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famCells->accessToMemArray().setSpecificDeallocator(0); }
834 ret->GetCellData()->AddArray(vtkTab);
838 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
841 vtkIntArray *vtkTab(vtkIntArray::New());
842 vtkTab->SetNumberOfComponents(1);
843 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
845 vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
847 { vtkTab->SetArray(numCells->getPointer(),numCells->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numCells->accessToMemArray().setSpecificDeallocator(0); }
848 ret->GetCellData()->AddArray(vtkTab);
852 // The arrays links to mesh
853 DataArrayInt *famNodes(0),*numNodes(0);
854 bool noCpyFamNodes(false),noCpyNumNodes(false);
855 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
858 vtkIntArray *vtkTab(vtkIntArray::New());
859 vtkTab->SetNumberOfComponents(1);
860 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
862 vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
864 { vtkTab->SetArray(famNodes->getPointer(),famNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); famNodes->accessToMemArray().setSpecificDeallocator(0); }
865 ret->GetPointData()->AddArray(vtkTab);
869 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
872 vtkIntArray *vtkTab(vtkIntArray::New());
873 vtkTab->SetNumberOfComponents(1);
874 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
876 vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),1,VTK_DATA_ARRAY_FREE);
878 { vtkTab->SetArray(numNodes->getPointer(),numNodes->getNbOfElems(),0,VTK_DATA_ARRAY_FREE); numNodes->accessToMemArray().setSpecificDeallocator(0); }
879 ret->GetPointData()->AddArray(vtkTab);
886 //////////////////////
888 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
892 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
895 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
896 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
897 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
898 ret+=(*it2).getNumberOfArrays();
902 void MEDFileFieldRepresentationTree::assignIds() const
905 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
906 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
907 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
911 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
913 std::size_t it0Cnt(0);
914 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
916 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
917 std::string tsName(oss.str());
918 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
920 std::string meshName((*it1)[0].getMeshName());
921 std::size_t it2Cnt(0);
922 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
924 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
925 std::string comSupStr(oss2.str());
926 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
932 void MEDFileFieldRepresentationTree::activateTheFirst() const
934 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
935 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
936 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
938 (*it2).activateAllArrays();
943 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
945 std::size_t it0Cnt(0);
946 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
948 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
949 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
951 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
952 std::vector<double> ts;
953 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
954 std::size_t nbOfTS(dtits.size());
955 std::ostringstream oss3; oss3 << nbOfTS;
956 names.push_back(oss3.str());
957 for(std::size_t i=0;i<nbOfTS;i++)
959 std::ostringstream oss4; oss4 << dtits[i].first;
960 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
961 names.push_back(oss4.str());
962 std::ostringstream oss5; oss5 << dtits[i].second;
963 vtkIdType ItId(sil->AddChild(DtId,edge));
964 names.push_back(oss5.str());
965 std::ostringstream oss6; oss6 << ts[i];
966 sil->AddChild(ItId,edge);
967 names.push_back(oss6.str());
970 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
971 std::string tsName(oss.str());
972 vtkIdType typeId0(sil->AddChild(root,edge));
973 names.push_back(tsName);
974 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
976 std::string meshName((*it1)[0].getMeshName());
977 vtkIdType typeId1(sil->AddChild(typeId0,edge));
978 names.push_back(meshName);
979 std::size_t it2Cnt(0);
980 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
982 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
983 std::string comSupStr(oss2.str());
984 vtkIdType typeId2(sil->AddChild(typeId1,edge));
985 names.push_back(comSupStr);
986 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
992 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
994 int dummy0(0),dummy1(0),dummy2(0);
995 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
996 std::string ret(leaf.getMeshName());
999 for(;i<_ms->getNumberOfMeshes();i++)
1001 m=_ms->getMeshAtPos(i);
1002 if(m->getName()==ret)
1005 if(i==_ms->getNumberOfMeshes())
1006 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1007 vtkIdType typeId0(sil->AddChild(root,edge));
1008 names.push_back(m->getName());
1010 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1011 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1012 std::vector<std::string> grps(m->getGroupsNames());
1013 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1015 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1016 names.push_back(*it0);
1017 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1018 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1020 sil->AddChild(typeId2,edge);
1021 names.push_back((*it1).c_str());
1025 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1026 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1027 std::vector<std::string> fams(m->getFamiliesNames());
1028 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1030 sil->AddChild(typeId11,edge);
1031 int famId(m->getFamilyId((*it00).c_str()));
1032 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1033 names.push_back(oss.str());
1038 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1040 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1041 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1042 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1043 if((*it2).containId(id))
1044 return (*it2).getLeafArr(id);
1045 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1048 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1050 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1051 return elt.getZeName();
1054 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1056 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1057 return elt.getZeNameC();
1060 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1062 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1063 return elt.getStatus();
1066 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1069 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1070 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1071 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1072 if((*it2).containZeName(name,ret))
1074 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1075 throw INTERP_KERNEL::Exception(msg.str().c_str());
1078 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1080 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1081 bool ret(elt.setStatus(status));//to be implemented
1085 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1088 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1089 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1090 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1091 ret=std::max(ret,(*it2).getNumberOfTS());
1098 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1102 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1104 _ms=MEDFileMeshes::New(fileName);
1105 _fields=MEDFileFields::New(fileName,false);//false is important to not read the values
1109 #ifdef MEDREADER_USE_MPI
1110 _ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1111 int nbMeshes(_ms->getNumberOfMeshes());
1112 for(int i=0;i<nbMeshes;i++)
1114 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1115 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1117 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> tmp3(tmp2->zipCoords());
1119 _fields=MEDFileFields::LoadPartOf(fileName,false,_ms);//false is important to not read the values
1121 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1122 throw INTERP_KERNEL::Exception(oss.str().c_str());
1128 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::SauvReader> sr(ParaMEDMEM::SauvReader::New(fileName));
1129 MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileData> mfd(sr->loadInMEDFileDS());
1130 _ms=mfd->getMeshes(); _ms->incrRef();
1131 int nbMeshes(_ms->getNumberOfMeshes());
1132 for(int i=0;i<nbMeshes;i++)
1134 ParaMEDMEM::MEDFileMesh *tmp(_ms->getMeshAtPos(i));
1135 ParaMEDMEM::MEDFileUMesh *tmp2(dynamic_cast<ParaMEDMEM::MEDFileUMesh *>(tmp));
1137 tmp2->forceComputationOfParts();
1139 _fields=mfd->getFields();
1140 if((ParaMEDMEM::MEDFileFields *)_fields)
1143 if(!((ParaMEDMEM::MEDFileFields *)_fields))
1145 _fields=BuildFieldFromMeshes(_ms);
1149 AppendFieldFromMeshes(_ms,_fields);
1151 _ms->cartesianizeMe();
1152 _fields->removeFieldsWithoutAnyTimeStep();
1153 std::vector<std::string> meshNames(_ms->getMeshesNames());
1154 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> > fields_per_mesh(meshNames.size());
1155 for(std::size_t i=0;i<meshNames.size();i++)
1157 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1159 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1161 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1163 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1165 MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1166 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1168 for(std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1170 if(!(*it)->presenceOfMultiDiscPerGeoType())
1171 allFMTSLeavesToDisplaySafe.push_back(*it);
1173 {// The case of some parts of field have more than one discretization per geo type.
1174 std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1175 std::size_t it0Cnt(0);
1176 for(std::vector< MEDCouplingAutoRefCountObjectPtr< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1178 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1179 (*it0)->setName(oss.str());
1180 allFMTSLeavesToDisplaySafe.push_back(*it0);
1187 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1188 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1190 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1192 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1193 // memory safety part
1194 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1195 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1197 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1198 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1200 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1201 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1204 // end of memory safety part
1205 // 1st : timesteps, 2nd : meshName, 3rd : common support
1206 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1207 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1209 std::vector< std::string > meshNamesLoc;
1210 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1211 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1213 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1214 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1215 if(it==meshNamesLoc.end())
1217 meshNamesLoc.push_back(meshName);
1218 splitByMeshName.resize(splitByMeshName.size()+1);
1219 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1222 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1224 _data_structure[i].resize(meshNamesLoc.size());
1225 for(std::size_t j=0;j<splitByMeshName.size();j++)
1227 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileFastCellSupportComparator> > fsp;
1228 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1229 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1230 sbmn[k]=splitByMeshName[j][k];
1231 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1232 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1233 std::vector< std::vector< MEDCouplingAutoRefCountObjectPtr<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1234 this->_data_structure[i][j].resize(commonSupSplit.size());
1235 for(std::size_t k=0;k<commonSupSplit.size();k++)
1237 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1238 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1240 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1241 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1244 for(std::size_t k=0;k<commonSupSplit.size();k++)
1245 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1248 this->removeEmptyLeaves();
1250 this->computeFullNameInLeaves();
1253 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1255 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1256 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1258 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1259 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1261 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1262 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1264 newSD1.push_back(*it2);
1266 newSD0.push_back(newSD1);
1269 newSD.push_back(newSD0);
1273 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1275 if(compInfos.size()!=1)
1277 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1280 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1282 return _data_structure[0][0][0].getMeshName();
1285 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1288 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1289 return leaf.getTimeSteps(tk);
1292 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk) const
1295 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1296 meshName=leaf.getMeshName();
1297 std::vector<double> ts(leaf.getTimeSteps(tk));
1298 std::size_t zeTimeId(0);
1301 std::vector<double> ts2(ts.size());
1302 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind2nd(std::plus<double>(),-timeReq));
1303 std::transform(ts2.begin(),ts2.end(),ts2.begin(),std::ptr_fun<double,double>(fabs));
1304 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind2nd(std::less<double>(),1e-14)));
1307 if(zeTimeId==(int)ts.size())
1308 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1309 if(zeTimeId==(int)ts.size())
1310 {//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.
1311 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1313 double valAttachedToPos(-std::numeric_limits<double>::max());
1314 for(std::size_t i=0;i<ts.size();i++)
1318 if(ts[i]>valAttachedToPos)
1321 valAttachedToPos=ts[i];
1326 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1327 valAttachedToPos=std::numeric_limits<double>::max();
1328 for(std::size_t i=0;i<ts.size();i++)
1330 if(ts[i]<valAttachedToPos)
1333 valAttachedToPos=ts[i];
1338 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1339 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1340 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1341 std::cerr << oss.str() << std::endl;
1345 tr=new MEDStdTimeReq((int)zeTimeId);
1347 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1348 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms));
1353 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1355 int nbOfActivated(0);
1356 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1357 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1358 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1359 if((*it2).isActivated())
1361 if(nbOfActivated!=1)
1363 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1364 throw INTERP_KERNEL::Exception(oss.str().c_str());
1366 int i0(0),i1(0),i2(0);
1367 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1368 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1369 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1370 if((*it2).isActivated())
1372 lev0=i0; lev1=i1; lev2=i2;
1375 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1378 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1381 os << "#############################################" << std::endl;
1382 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1385 os << "TS" << i << std::endl;
1386 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1389 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1392 os << " " << (*it2).getMeshName() << std::endl;
1393 os << " Comp" << k << std::endl;
1394 (*it2).printMySelf(os);
1398 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1401 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1403 std::map<std::string,bool> ret;
1404 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1405 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1406 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1407 (*it2).dumpState(ret);
1411 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms, ParaMEDMEM::MEDFileFields *ret)
1414 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1415 for(int i=0;i<ms->getNumberOfMeshes();i++)
1417 MEDFileMesh *mm(ms->getMeshAtPos(i));
1418 std::vector<int> levs(mm->getNonEmptyLevels());
1419 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileField1TS> f1tsMultiLev(ParaMEDMEM::MEDFileField1TS::New());
1420 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1423 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1425 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1426 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1428 ParaMEDMEM::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1429 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1431 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1432 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1435 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1436 f1tsMultiLev->setFieldNoProfileSBT(f);
1441 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1442 if(levsExt.size()==levs.size()+1)
1444 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1445 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_NODES));
1447 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1448 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1449 arr->iota(); f->setArray(arr);
1450 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1451 f1tsMultiLev->setFieldNoProfileSBT(f);
1459 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1460 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDCouplingFieldDouble> f(ParaMEDMEM::MEDCouplingFieldDouble::New(ParaMEDMEM::ON_CELLS));
1462 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::DataArrayDouble> arr(ParaMEDMEM::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1463 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1466 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1467 f1tsMultiLev->setFieldNoProfileSBT(f);
1470 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFieldMultiTS> fmtsMultiLev(ParaMEDMEM::MEDFileFieldMultiTS::New());
1471 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1472 ret->pushField(fmtsMultiLev);
1476 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const ParaMEDMEM::MEDFileFields *ret)
1478 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1480 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1481 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1482 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1484 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1485 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1486 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1490 ParaMEDMEM::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const ParaMEDMEM::MEDFileMeshes *ms)
1492 ParaMEDMEM::MEDCouplingAutoRefCountObjectPtr<ParaMEDMEM::MEDFileFields> ret(ParaMEDMEM::MEDFileFields::New());
1493 AppendFieldFromMeshes(ms,ret);
1497 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1499 std::vector<std::string> ret;
1501 while(pos!=std::string::npos)
1503 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1504 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1506 pos=fullFieldName.find_first_not_of(sep,curPos);
1512 * Here the non regression tests.
1513 * const char inp0[]="";
1514 * const char exp0[]="";
1515 * const char inp1[]="field";
1516 * const char exp1[]="field";
1517 * const char inp2[]="_________";
1518 * const char exp2[]="_________";
1519 * const char inp3[]="field_p";
1520 * const char exp3[]="field_p";
1521 * const char inp4[]="field__p";
1522 * const char exp4[]="field_p";
1523 * const char inp5[]="field_p__";
1524 * const char exp5[]="field_p";
1525 * const char inp6[]="field_p_";
1526 * const char exp6[]="field_p";
1527 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1528 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1529 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1530 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1531 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1532 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1533 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1534 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1536 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1538 const char SEP('_');
1539 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1541 return fullFieldName;//should never happen
1545 return fullFieldName;
1549 std::string ret(v[0]);
1550 for(std::size_t i=1;i<v.size();i++)
1555 { ret+=SEP; ret+=v[i]; }
1561 return fullFieldName;
1567 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1571 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1576 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1578 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1580 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1586 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1587 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1588 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1590 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1592 std::size_t sz(ts.size());
1593 bool isInHumanRange(true);
1595 for(std::size_t i=0;i<sz;i++)
1598 if(ts[i]<=-1e299 || ts[i]>=1e299)
1599 isInHumanRange=false;
1602 return processedUsingPairOfIds(tsPairs);
1604 return processedUsingPairOfIds(tsPairs);
1605 _postprocessed_time=ts;
1606 return getPostProcessedTime();
1611 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1613 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1615 std::size_t sz(ts.size());
1616 std::vector<double> ts2(sz);
1618 for(std::size_t i=0;i<sz;i++)
1623 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1626 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1628 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1632 void TimeKeeper::printSelf(std::ostream& oss) const
1634 std::size_t sz(_activated_ts.size());
1635 for(std::size_t i=0;i<sz;i++)
1637 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1641 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1643 std::size_t sz(_activated_ts.size());
1644 std::vector<bool> ret(sz);
1645 for(std::size_t i=0;i<sz;i++)
1647 ret[i]=_activated_ts[i].first;
1652 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1654 std::size_t sz(tsPairs.size());
1655 std::set<int> s0,s1;
1656 for(std::size_t i=0;i<sz;i++)
1657 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1660 _postprocessed_time.resize(sz);
1661 for(std::size_t i=0;i<sz;i++)
1662 _postprocessed_time[i]=(double)tsPairs[i].first;
1663 return getPostProcessedTime();
1667 _postprocessed_time.resize(sz);
1668 for(std::size_t i=0;i<sz;i++)
1669 _postprocessed_time[i]=(double)tsPairs[i].second;
1670 return getPostProcessedTime();
1672 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1673 _postprocessed_time.resize(sz);
1674 for(std::size_t i=0;i<sz;i++)
1675 _postprocessed_time[i]=(double)i;
1676 return getPostProcessedTime();
1679 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1681 _activated_ts.resize(maxNumberOfTS);
1682 for(int i=0;i<maxNumberOfTS;i++)
1684 std::ostringstream oss; oss << "000" << i;
1685 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());