1 // Copyright (C) 2010-2021 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"
30 #include "MEDCouplingMemArray.txx"
32 #ifdef MEDREADER_USE_MPI
33 #include "ParaMEDFileMesh.hxx"
36 #include "vtkXMLUnstructuredGridWriter.h"//
38 #include "vtkUnstructuredGrid.h"
39 #include "vtkRectilinearGrid.h"
40 #include "vtkStructuredGrid.h"
41 #include "vtkUnsignedCharArray.h"
42 #include "vtkQuadratureSchemeDefinition.h"
43 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
44 #include "vtkInformationIntegerKey.h"
45 #include "vtkInformation.h"
46 #include "vtkAOSDataArrayTemplate.h"
47 #include "vtkIdTypeArray.h"
48 #include "vtkDoubleArray.h"
49 #include "vtkIntArray.h"
50 #include "vtkLongArray.h"
52 #include "vtkLongLongArray.h"
54 #include "vtkFloatArray.h"
55 #include "vtkCellArray.h"
56 #include "vtkPointData.h"
57 #include "vtkFieldData.h"
58 #include "vtkCellData.h"
60 #include "vtkMutableDirectedGraph.h"
62 using namespace MEDCoupling;
64 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
66 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
68 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
70 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
72 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
74 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
76 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
78 const char MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME[]="GlobalNodeIds";// WARNING DO NOT CHANGE IT BEFORE HAVING CHECKED IN PV SOURCES !
80 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
82 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
84 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
87 vtkIdTypeArray *ELGACmp::findOrCreate(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd, vtkDataSet *ds, bool& isNew, ExportedTinyInfo *internalInfo) const
89 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
98 return createNew<T>(globs,locsReallyUsed,vtkd,ds,internalInfo);
102 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd) const
104 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
105 if(it==_loc_names.end())
107 std::size_t pos(std::distance(_loc_names.begin(),it));
108 vtkIdTypeArray *ret(_elgas[pos]);
109 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
110 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
112 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
114 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
119 vtkIdTypeArray *ELGACmp::createNew(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
121 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_DELETE;
122 std::vector< std::vector<std::string> > locNames(_loc_names);
123 std::vector<vtkIdTypeArray *> elgas(_elgas);
124 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
126 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
127 if(it!=locNames.end())
128 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
129 locNames.push_back(locsReallyUsed);
130 vtkIdTypeArray *elga(vtkIdTypeArray::New());
131 elga->SetNumberOfComponents(1);
132 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
133 std::map<unsigned char,int> m;
134 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
136 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
137 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
138 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
139 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
140 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
141 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
142 // 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.
143 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
144 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
146 internalInfo->pushGaussAdditionnalInfo(vtkType,dimLoc,refCoods2,gsCoods2);
147 double *shape(new double[nbPtsPerCell*nbGaussPt]);
148 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
149 calculator.initLocalInfo();
150 const std::vector<double>& wgths(loc.getGaussWeights());
151 for(int i=0;i<nbGaussPt;i++)
153 const double *pt0(calculator.getFunctionValues(i));
154 if(ct!=INTERP_KERNEL::NORM_HEXA27)
155 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
158 for(int j=0;j<27;j++)
159 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
162 m[vtkType]=nbGaussPt;
163 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
165 key->Set(elga->GetInformation(),def,vtkType);
166 key->Set(vtkd->GetInformation(),def,vtkType);
167 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
170 vtkIdType ncell(ds->GetNumberOfCells());
171 vtkIdType *pt(new vtkIdType[ncell]),offset(0);
172 for(vtkIdType cellId=0;cellId<ncell;cellId++)
174 vtkCell *cell(ds->GetCell(cellId));
175 vtkIdType delta(m[(unsigned char)cell->GetCellType()]);
179 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
180 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
181 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
182 std::string ossStr(oss.str());
183 elga->SetName(ossStr.c_str());
184 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
185 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
186 elgas.push_back(elga);
190 _defs.push_back(defs);
194 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
196 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
197 ds->GetCellData()->AddArray(*it);
202 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
204 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
205 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
206 (*it1).first->Delete();
212 class MEDFileVTKTraits
215 typedef void VtkType;
220 class MEDFileVTKTraits<int>
223 typedef vtkIntArray VtkType;
224 typedef MEDCoupling::DataArrayInt MCType;
229 class MEDFileVTKTraits<long long>
231 class MEDFileVTKTraits<long>
236 typedef vtkLongLongArray VtkType;
238 typedef vtkLongArray VtkType;
240 typedef MEDCoupling::DataArrayInt64 MCType;
244 class MEDFileVTKTraits<float>
247 typedef vtkFloatArray VtkType;
248 typedef MEDCoupling::DataArrayFloat MCType;
252 class MEDFileVTKTraits<double>
255 typedef vtkDoubleArray VtkType;
256 typedef MEDCoupling::DataArrayDouble MCType;
259 typedef typename MEDFileVTKTraits<mcIdType>::VtkType vtkMCIdTypeArray;
263 void AssignDataPointerToVTK(typename MEDFileVTKTraits<T>::VtkType *vtkTab, typename MEDFileVTKTraits<T>::MCType *mcTab, bool noCpyNumNodes)
266 vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),1,vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE);
268 { vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),0,vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE); mcTab->accessToMemArray().setSpecificDeallocator(0); }
271 // here copy is always assumed.
272 template<class VTKT, class MCT>
273 void AssignDataPointerOther(VTKT *vtkTab, MCT *mcTab, vtkIdType nbElems)
275 typedef typename VTKT::ValueType VTKType;
276 if ( sizeof( VTKType ) == sizeof( typename MCT::Type ))
278 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(mcTab->getPointer()),nbElems,0,VTKT::VTK_DATA_ARRAY_FREE);
279 mcTab->accessToMemArray().setSpecificDeallocator(0);
283 VTKType* newArray = new VTKType[ nbElems ];
284 std::copy( mcTab->begin(), mcTab->begin() + nbElems, newArray );
285 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(newArray),nbElems,0,VTKT::VTK_DATA_ARRAY_DELETE);
290 void AssignToFieldData(DataArray *vPtr, const MEDTimeReq *tr, vtkFieldData *att, const std::string& crudeName, bool noCpyNumNodes,
291 const std::vector<TypeOfField>& discs, const ELGACmp& elgaCmp, const MEDCoupling::MEDFileFieldGlobsReal *globs,
292 MEDFileAnyTypeField1TS *f1ts, vtkDataSet *ds, ExportedTinyInfo *internalInfo)
294 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_DELETE;
295 typename MEDFileVTKTraits<T>::MCType *vi(static_cast<typename MEDFileVTKTraits<T>::MCType *>(vPtr));
296 typename MEDFileVTKTraits<T>::VtkType *vtkd(MEDFileVTKTraits<T>::VtkType::New());
297 vtkd->SetNumberOfComponents((int)vi->getNumberOfComponents());
298 for(unsigned int i=0;i<vi->getNumberOfComponents();i++)
299 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
300 AssignDataPointerToVTK<T>(vtkd,vi,noCpyNumNodes);
301 std::string name(tr->buildName(crudeName));
302 vtkd->SetName(name.c_str());
305 if(discs[0]==ON_GAUSS_PT)
308 elgaCmp.findOrCreate<T>(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp,internalInfo);
310 if(discs[0]==ON_GAUSS_NE)
312 vtkIdTypeArray *elno(vtkIdTypeArray::New());
313 elno->SetNumberOfComponents(1);
314 vtkIdType ncell(ds->GetNumberOfCells());
315 vtkIdType *pt(new vtkIdType[ncell]),offset(0);
316 std::set<int> cellTypes;
317 for(vtkIdType cellId=0;cellId<ncell;cellId++)
319 vtkCell *cell(ds->GetCell(cellId));
320 vtkIdType delta(cell->GetNumberOfPoints());
321 cellTypes.insert(cell->GetCellType());
325 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
326 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
327 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
328 elno->SetName(nameElno.c_str());
329 ds->GetCellData()->AddArray(elno);
330 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
331 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
333 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
334 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
336 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
337 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
339 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
340 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
341 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
342 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
343 double *shape(new double[nbGaussPt*nbGaussPt]);
345 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));//GetLocsFromGeometricType
346 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
347 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
348 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
349 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
350 calculator.initLocalInfo();
351 for(int i=0;i<nbGaussPt;i++)
353 const double *pt0(calculator.getFunctionValues(i));
354 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
356 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
358 key->Set(elno->GetInformation(),def,*it);
359 key->Set(vtkd->GetInformation(),def,*it);
369 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
373 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>& arr):MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
375 std::vector< std::vector<MEDCoupling::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
377 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
378 if(typs[0].size()!=1)
379 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 !");
380 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
381 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
385 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
387 MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>::operator=(other);
390 _ze_name=other._ze_name;
391 _ze_full_name.clear();
395 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
400 int MEDFileFieldRepresentationLeavesArrays::getId() const
405 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
407 return _ze_full_name;
410 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
412 return _ze_full_name.c_str();
415 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
417 vtkIdType refId(sil->AddChild(root,edge));
418 names.push_back(_ze_name);
420 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
422 sil->AddChild(refId,edge);
423 names.push_back(std::string());
427 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
429 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
430 _ze_full_name=oss3.str();
433 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
438 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
440 bool ret(_activated!=status);
445 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshStruct *mst, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
447 //const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE; // todo: unused
448 tr->setNumberOfTS((operator->())->getNumberOfTS());
450 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
452 MCAuto<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
453 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
454 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
455 MEDFileInt32Field1TS *f1tsPtrInt(dynamic_cast<MEDFileInt32Field1TS *>(f1tsPtr));
456 MEDFileInt64Field1TS *f1tsPtrInt64(dynamic_cast<MEDFileInt64Field1TS *>(f1tsPtr));
457 MEDFileFloatField1TS *f1tsPtrFloat(dynamic_cast<MEDFileFloatField1TS *>(f1tsPtr));
458 DataArray *crudeArr(0),*postProcessedArr(0);
460 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
462 crudeArr=f1tsPtrInt->getUndergroundDataArray();
463 else if(f1tsPtrInt64)
464 crudeArr=f1tsPtrInt64->getUndergroundDataArray();
465 else if(f1tsPtrFloat)
466 crudeArr=f1tsPtrFloat->getUndergroundDataArray();
468 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64, FLOAT32 and INT32 fields are dealt for the moment !");
469 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
470 f1ts->loadArraysIfNecessary();
471 MCAuto<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
474 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
476 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
477 vtkFieldData *att(0);
482 att=ds->GetCellData();
487 att=ds->GetPointData();
492 att=ds->GetFieldData();
497 att=ds->GetFieldData();
501 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
505 AssignToFieldData<double>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
509 AssignToFieldData<int>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
511 else if(f1tsPtrFloat)
513 AssignToFieldData<float>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
515 else if(f1tsPtrInt64)
517 AssignToFieldData<Int64>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
520 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
524 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
526 _elga_cmp.appendELGAIfAny(ds);
531 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
535 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS> >& arr,
536 const MEDCoupling::MCAuto<MEDCoupling::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
538 for(std::size_t i=0;i<arr.size();i++)
539 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
542 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
545 _cached_ds->Delete();
548 bool MEDFileFieldRepresentationLeaves::empty() const
550 const MEDFileFastCellSupportComparator *fcscp(_fsp);
551 return fcscp==0 || _arrays.empty();
554 void MEDFileFieldRepresentationLeaves::setId(int& id) const
556 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
560 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
562 return _arrays[0]->getMeshName();
565 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
567 return (int)_arrays.size();
570 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
572 return _arrays[0]->getNumberOfTS();
575 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
577 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
578 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
582 * \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.
583 * \param [in] meshName
589 void MEDFileFieldRepresentationLeaves::feedSIL(const MEDCoupling::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
591 vtkIdType root2(sil->AddChild(root,edge));
592 names.push_back(std::string("Arrs"));
593 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
594 (*it).feedSIL(sil,root2,edge,names);
596 vtkIdType root3(sil->AddChild(root,edge));
597 names.push_back(std::string("InfoOnGeoType"));
598 const MEDCoupling::MEDFileMesh *m(0);
600 m=ms->getMeshWithName(meshName);
601 const MEDCoupling::MEDFileFastCellSupportComparator *fsp(_fsp);
602 if(!fsp || fsp->getNumberOfTS()==0)
604 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
605 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
607 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
608 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
609 sil->AddChild(root3,edge);
610 names.push_back(cmStr);
614 bool MEDFileFieldRepresentationLeaves::containId(int id) const
616 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
617 if((*it).getId()==id)
622 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
624 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
625 if((*it).getZeName()==name)
633 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
635 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
636 status[(*it).getZeName()]=(*it).getStatus();
639 bool MEDFileFieldRepresentationLeaves::isActivated() const
641 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
642 if((*it).getStatus())
647 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
649 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
651 os << " - " << (*it0).getZeName() << " (";
652 if((*it0).getStatus())
656 os << ")" << std::endl;
660 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
662 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
663 (*it).setStatus(true);
666 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
668 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
669 if((*it).getId()==id)
671 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
674 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
677 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
678 std::vector<double> ret;
679 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
680 return tk.getTimeStepsRegardingPolicy(dtits,ret);
683 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
686 return _arrays[0]->getTimeSteps(ts);
690 return std::vector< std::pair<int,int> >();
694 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
696 std::ostringstream oss;
697 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
698 std::vector<double> ret1;
699 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
700 std::size_t sz(ret1.size());
701 for(std::size_t i=0;i<sz;i++)
703 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
706 std::string tmp(oss.str());
707 if(tmp.size()>200 && i!=sz-1)
717 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshes *meshes, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
720 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
721 MCAuto<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
722 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
723 if((*it).getStatus())
725 (*it).appendFields(tr,globs,mml,mst,ds,internalInfo);
726 (*it).appendELGAIfAny(ds);
730 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
732 DataArrayDouble *coordsMC(0);
733 DataArrayByte *typesMC(0);
734 DataArrayIdType *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
735 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
736 MCAuto<DataArrayDouble> coordsSafe(coordsMC);
737 MCAuto<DataArrayByte> typesSafe(typesMC);
738 MCAuto<DataArrayIdType> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
740 vtkIdType nbOfCells(typesSafe->getNbOfElems());
741 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
742 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
743 AssignDataPointerOther<vtkUnsignedCharArray,DataArrayByte>(cellTypes,typesSafe,nbOfCells);
744 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
745 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(cellLocations,cellLocationsSafe,nbOfCells);
746 vtkCellArray *cells(vtkCellArray::New());
747 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
748 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(cells2,cellsSafe,cellsSafe->getNbOfElems());
749 cells->SetCells(nbOfCells,cells2);
751 if(faceLocationsMC!=0 && facesMC!=0)
753 vtkIdTypeArray *faces(vtkIdTypeArray::New());
754 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(faces,facesSafe,facesSafe->getNbOfElems());
755 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
756 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(faceLocations,faceLocationsSafe,faceLocationsSafe->getNbOfElems());
757 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
758 faceLocations->Delete();
762 ret->SetCells(cellTypes,cellLocations,cells);
764 cellLocations->Delete();
766 vtkPoints *pts(vtkPoints::New());
767 vtkDoubleArray *pts2(vtkDoubleArray::New());
768 pts2->SetNumberOfComponents(3);
769 AssignDataPointerToVTK<double>(pts2,coordsSafe,statusOfCoords);
778 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(MEDCoupling::MEDCMeshMultiLev *mm) const
781 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
782 vtkDoubleArray *vtkTmp(0);
783 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
784 std::size_t dim(arrs.size());
786 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
787 int sizePerAxe[3]={1,1,1};
788 sizePerAxe[0]=arrs[0]->getNbOfElems();
790 sizePerAxe[1]=arrs[1]->getNbOfElems();
792 sizePerAxe[2]=arrs[2]->getNbOfElems();
793 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
794 vtkTmp=vtkDoubleArray::New();
795 vtkTmp->SetNumberOfComponents(1);
796 AssignDataPointerToVTK<double>(vtkTmp,arrs[0],isInternal);
797 ret->SetXCoordinates(vtkTmp);
802 vtkTmp=vtkDoubleArray::New();
803 vtkTmp->SetNumberOfComponents(1);
804 AssignDataPointerToVTK<double>(vtkTmp,arrs[1],isInternal);
805 ret->SetYCoordinates(vtkTmp);
811 vtkTmp=vtkDoubleArray::New();
812 vtkTmp->SetNumberOfComponents(1);
813 AssignDataPointerToVTK<double>(vtkTmp,arrs[2],isInternal);
814 ret->SetZCoordinates(vtkTmp);
821 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(MEDCoupling::MEDCurveLinearMeshMultiLev *mm) const
823 int meshStr[3]={1,1,1};
824 DataArrayDouble *coords(0);
825 std::vector<mcIdType> nodeStrct;
827 mm->buildVTUArrays(coords,nodeStrct,isInternal);
828 std::size_t dim(nodeStrct.size());
830 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
831 meshStr[0]=nodeStrct[0];
833 meshStr[1]=nodeStrct[1];
835 meshStr[2]=nodeStrct[2];
836 vtkStructuredGrid *ret(vtkStructuredGrid::New());
837 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
838 vtkDoubleArray *da(vtkDoubleArray::New());
839 da->SetNumberOfComponents(3);
840 if(coords->getNumberOfComponents()==3)
841 AssignDataPointerToVTK<double>(da,coords,isInternal);//if isIntenal==True VTK has not the ownership of double * because MEDLoader main struct has it !
844 MCAuto<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
845 AssignDataPointerToVTK<double>(da,coords2,false);//let VTK deal with double *
848 vtkPoints *points=vtkPoints::New();
849 ret->SetPoints(points);
856 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDFileMeshes *meshes, ExportedTinyInfo *internalInfo) const
859 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
860 MCAuto<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
861 MCAuto<MEDMeshMultiLev> mml2(mml->prepare());
862 MEDMeshMultiLev *ptMML2(mml2);
865 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
866 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
867 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
871 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
875 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
879 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
882 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
883 _cached_ds=ret->NewInstance();
884 _cached_ds->ShallowCopy(ret);
888 ret=_cached_ds->NewInstance();
889 ret->ShallowCopy(_cached_ds);
892 appendFields(tr,globs,mml,meshes,ret,internalInfo);
893 // The arrays links to mesh
894 DataArrayIdType *famCells(0),*numCells(0);
895 bool noCpyFamCells(false),noCpyNumCells(false);
896 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
899 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
900 vtkTab->SetNumberOfComponents(1);
901 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
902 AssignDataPointerToVTK<mcIdType>(vtkTab,famCells,noCpyFamCells);
903 ret->GetCellData()->AddArray(vtkTab);
907 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
910 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
911 vtkTab->SetNumberOfComponents(1);
912 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
913 AssignDataPointerToVTK<mcIdType>(vtkTab,numCells,noCpyNumCells);
914 ret->GetCellData()->AddArray(vtkTab);
918 // The arrays links to mesh
919 DataArrayIdType *famNodes(0),*numNodes(0);
920 bool noCpyFamNodes(false),noCpyNumNodes(false);
921 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
924 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
925 vtkTab->SetNumberOfComponents(1);
926 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
927 AssignDataPointerToVTK<mcIdType>(vtkTab,famNodes,noCpyFamNodes);
928 ret->GetPointData()->AddArray(vtkTab);
932 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
935 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
936 vtkTab->SetNumberOfComponents(1);
937 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
938 AssignDataPointerToVTK<mcIdType>(vtkTab,numNodes,noCpyNumNodes);
939 ret->GetPointData()->AddArray(vtkTab);
943 // Global Node Ids if any ! (In // mode)
944 DataArrayIdType *gni(ptMML2->retrieveGlobalNodeIdsIfAny());
947 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
948 vtkTab->SetNumberOfComponents(1);
949 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME);
950 AssignDataPointerToVTK<mcIdType>(vtkTab,gni,false);
951 ret->GetPointData()->AddArray(vtkTab);
958 //////////////////////
960 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
964 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
967 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
968 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
969 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
970 ret+=(*it2).getNumberOfArrays();
974 void MEDFileFieldRepresentationTree::assignIds() const
977 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
978 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
979 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
983 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
985 std::size_t it0Cnt(0);
986 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
988 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
989 std::string tsName(oss.str());
990 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
992 std::string meshName((*it1)[0].getMeshName());
993 std::size_t it2Cnt(0);
994 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
996 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
997 std::string comSupStr(oss2.str());
998 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
1004 void MEDFileFieldRepresentationTree::activateTheFirst() const
1006 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1007 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1008 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1010 (*it2).activateAllArrays();
1015 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1017 std::size_t it0Cnt(0);
1018 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
1020 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
1021 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
1023 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
1024 std::vector<double> ts;
1025 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
1026 std::size_t nbOfTS(dtits.size());
1027 std::ostringstream oss3; oss3 << nbOfTS;
1028 names.push_back(oss3.str());
1029 for(std::size_t i=0;i<nbOfTS;i++)
1031 std::ostringstream oss4; oss4 << dtits[i].first;
1032 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
1033 names.push_back(oss4.str());
1034 std::ostringstream oss5; oss5 << dtits[i].second;
1035 vtkIdType ItId(sil->AddChild(DtId,edge));
1036 names.push_back(oss5.str());
1037 std::ostringstream oss6; oss6 << ts[i];
1038 sil->AddChild(ItId,edge);
1039 names.push_back(oss6.str());
1042 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
1043 std::string tsName(oss.str());
1044 vtkIdType typeId0(sil->AddChild(root,edge));
1045 names.push_back(tsName);
1046 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1048 std::string meshName((*it1)[0].getMeshName());
1049 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1050 names.push_back(meshName);
1051 std::size_t it2Cnt(0);
1052 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
1054 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
1055 std::string comSupStr(oss2.str());
1056 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1057 names.push_back(comSupStr);
1058 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
1064 std::string MEDFileFieldRepresentationTree::getActiveMeshName() const
1066 int dummy0(0),dummy1(0),dummy2(0);
1067 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1068 return leaf.getMeshName();
1071 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1073 int dummy0(0),dummy1(0),dummy2(0);
1074 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1075 std::string ret(leaf.getMeshName());
1078 for(;i<_ms->getNumberOfMeshes();i++)
1080 m=_ms->getMeshAtPos(i);
1081 if(m->getName()==ret)
1084 if(i==_ms->getNumberOfMeshes())
1085 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1086 vtkIdType typeId0(sil->AddChild(root,edge));
1087 names.push_back(m->getName());
1089 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1090 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1091 std::vector<std::string> grps(m->getGroupsNames());
1092 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1094 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1095 names.push_back(*it0);
1096 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1097 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1099 sil->AddChild(typeId2,edge);
1100 names.push_back((*it1).c_str());
1104 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1105 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1106 std::vector<std::string> fams(m->getFamiliesNames());
1107 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1109 sil->AddChild(typeId11,edge);
1110 int famId(m->getFamilyId((*it00).c_str()));
1111 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1112 names.push_back(oss.str());
1117 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1119 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1120 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1121 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1122 if((*it2).containId(id))
1123 return (*it2).getLeafArr(id);
1124 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1127 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1129 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1130 return elt.getZeName();
1133 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1135 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1136 return elt.getZeNameC();
1139 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1141 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1142 return elt.getStatus();
1145 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1148 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1149 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1150 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1151 if((*it2).containZeName(name,ret))
1153 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1154 throw INTERP_KERNEL::Exception(msg.str().c_str());
1157 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1159 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1160 bool ret(elt.setStatus(status));//to be implemented
1164 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1167 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1168 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1169 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1170 ret=std::max(ret,(*it2).getNumberOfTS());
1177 void MEDFileFieldRepresentationTree::loadInMemory(MEDCoupling::MEDFileFields *fields, MEDCoupling::MEDFileMeshes *meshes)
1179 _fields=fields; _ms=meshes;
1180 if(_fields.isNotNull())
1185 if(_fields.isNull())
1187 _fields=BuildFieldFromMeshes(_ms);
1191 AppendFieldFromMeshes(_ms,_fields);
1193 _ms->cartesianizeMe();
1194 _fields->removeFieldsWithoutAnyTimeStep();
1195 std::vector<std::string> meshNames(_ms->getMeshesNames());
1196 std::vector< MCAuto<MEDFileFields> > fields_per_mesh(meshNames.size());
1197 for(std::size_t i=0;i<meshNames.size();i++)
1199 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1201 std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1202 for(std::vector< MCAuto<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1204 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1206 MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1207 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1209 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1211 if(!(*it)->presenceOfMultiDiscPerGeoType())
1212 allFMTSLeavesToDisplaySafe.push_back(*it);
1214 {// The case of some parts of field have more than one discretization per geo type.
1215 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1216 std::size_t it0Cnt(0);
1217 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1219 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1220 (*it0)->setName(oss.str());
1221 allFMTSLeavesToDisplaySafe.push_back(*it0);
1228 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1229 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1231 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1233 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1234 // memory safety part
1235 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1236 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1238 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1239 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1241 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1242 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1245 // end of memory safety part
1246 // 1st : timesteps, 2nd : meshName, 3rd : common support
1247 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1248 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1250 std::vector< std::string > meshNamesLoc;
1251 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1252 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1254 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1255 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1256 if(it==meshNamesLoc.end())
1258 meshNamesLoc.push_back(meshName);
1259 splitByMeshName.resize(splitByMeshName.size()+1);
1260 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1263 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1265 _data_structure[i].resize(meshNamesLoc.size());
1266 for(std::size_t j=0;j<splitByMeshName.size();j++)
1268 std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsp;
1269 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1270 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1271 sbmn[k]=splitByMeshName[j][k];
1272 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1273 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1274 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1275 this->_data_structure[i][j].resize(commonSupSplit.size());
1276 for(std::size_t k=0;k<commonSupSplit.size();k++)
1278 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1279 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1281 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1282 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1285 for(std::size_t k=0;k<commonSupSplit.size();k++)
1286 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1289 this->removeEmptyLeaves();
1291 this->computeFullNameInLeaves();
1294 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, bool isMEDOrSauv, int iPart, int nbOfParts)
1296 MCAuto<MEDFileMeshes> ms;
1297 MCAuto<MEDFileFields> fields;
1300 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1302 MCAuto<MEDFileMeshSupports> msups(MEDFileMeshSupports::New(fileName));
1303 MCAuto<MEDFileStructureElements> mse(MEDFileStructureElements::New(fileName,msups));
1304 ms=MEDFileMeshes::New(fileName);
1305 fields=MEDFileFields::NewWithDynGT(fileName,mse,false);//false is important to not read the values
1306 if(ms->presenceOfStructureElements())
1308 fields->loadArrays();
1309 fields->blowUpSE(ms,mse);
1311 int nbMeshes(ms->getNumberOfMeshes());
1312 for(int i=0;i<nbMeshes;i++)
1314 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1315 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1317 tmp2->forceComputationOfParts();
1322 #ifdef MEDREADER_USE_MPI
1323 ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1324 int nbMeshes(ms->getNumberOfMeshes());
1325 for(int i=0;i<nbMeshes;i++)
1327 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1328 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1330 MCAuto<DataArrayIdType> tmp3(tmp2->zipCoords());
1332 fields=MEDFileFields::LoadPartOf(fileName,false,ms);//false is important to not read the values
1334 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1335 throw INTERP_KERNEL::Exception(oss.str().c_str());
1341 MCAuto<MEDCoupling::SauvReader> sr(MEDCoupling::SauvReader::New(fileName));
1342 MCAuto<MEDCoupling::MEDFileData> mfd(sr->loadInMEDFileDS());
1343 ms=mfd->getMeshes(); ms->incrRef();
1344 int nbMeshes(ms->getNumberOfMeshes());
1345 for(int i=0;i<nbMeshes;i++)
1347 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1348 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1350 tmp2->forceComputationOfParts();
1352 fields=mfd->getFields();
1353 if(fields.isNotNull())
1356 loadInMemory(fields,ms);
1359 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1361 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1362 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1364 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1365 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1367 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1368 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1370 newSD1.push_back(*it2);
1372 newSD0.push_back(newSD1);
1375 newSD.push_back(newSD0);
1379 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1381 if(compInfos.size()!=1)
1383 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1386 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1388 return _data_structure[0][0][0].getMeshName();
1391 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1394 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1395 return leaf.getTimeSteps(tk);
1398 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk, ExportedTinyInfo *internalInfo) const
1401 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1402 meshName=leaf.getMeshName();
1403 std::vector<double> ts(leaf.getTimeSteps(tk));
1404 std::size_t zeTimeId(0);
1407 std::vector<double> ts2(ts.size());
1408 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind(std::plus<double>(),std::placeholders::_1,-timeReq));
1409 std::transform(ts2.begin(),ts2.end(),ts2.begin(),[](double c) {return fabs(c);});
1410 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind(std::less<double>(),std::placeholders::_1,1e-14)));
1413 if(zeTimeId==ts.size())
1414 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1415 if(zeTimeId==ts.size())
1416 {//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.
1417 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1419 double valAttachedToPos(-std::numeric_limits<double>::max());
1420 for(std::size_t i=0;i<ts.size();i++)
1424 if(ts[i]>valAttachedToPos)
1427 valAttachedToPos=ts[i];
1432 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1433 valAttachedToPos=std::numeric_limits<double>::max();
1434 for(std::size_t i=0;i<ts.size();i++)
1436 if(ts[i]<valAttachedToPos)
1439 valAttachedToPos=ts[i];
1444 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1445 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1446 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1447 std::cerr << oss.str() << std::endl;
1451 tr=new MEDStdTimeReq((int)zeTimeId);
1453 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1454 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms,internalInfo));
1459 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1461 int nbOfActivated(0);
1462 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1463 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1464 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1465 if((*it2).isActivated())
1467 if(nbOfActivated!=1)
1469 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1470 throw INTERP_KERNEL::Exception(oss.str().c_str());
1472 int i0(0),i1(0),i2(0);
1473 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1474 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1475 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1476 if((*it2).isActivated())
1478 lev0=i0; lev1=i1; lev2=i2;
1481 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1484 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1487 os << "#############################################" << std::endl;
1488 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1491 os << "TS" << i << std::endl;
1492 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1495 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1498 os << " " << (*it2).getMeshName() << std::endl;
1499 os << " Comp" << k << std::endl;
1500 (*it2).printMySelf(os);
1504 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1507 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1509 std::map<std::string,bool> ret;
1510 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1511 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1512 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1513 (*it2).dumpState(ret);
1517 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms, MEDCoupling::MEDFileFields *ret)
1520 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1521 for(int i=0;i<ms->getNumberOfMeshes();i++)
1523 MEDFileMesh *mm(ms->getMeshAtPos(i));
1524 std::vector<int> levs(mm->getNonEmptyLevels());
1525 MEDCoupling::MCAuto<MEDCoupling::MEDFileField1TS> f1tsMultiLev(MEDCoupling::MEDFileField1TS::New());
1526 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1529 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1531 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1532 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1534 MEDCoupling::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1535 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1537 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1538 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1541 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1542 f1tsMultiLev->setFieldNoProfileSBT(f);
1547 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1548 if(levsExt.size()==levs.size()+1)
1550 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1551 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_NODES));
1553 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1554 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1555 arr->iota(); f->setArray(arr);
1556 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1557 f1tsMultiLev->setFieldNoProfileSBT(f);
1565 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1566 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1568 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1569 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1572 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1573 f1tsMultiLev->setFieldNoProfileSBT(f);
1576 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> fmtsMultiLev(MEDCoupling::MEDFileFieldMultiTS::New());
1577 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1578 ret->pushField(fmtsMultiLev);
1582 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const MEDCoupling::MEDFileFields *ret)
1584 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1586 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1587 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1588 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1590 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1591 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1592 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1596 MEDCoupling::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms)
1598 MEDCoupling::MCAuto<MEDCoupling::MEDFileFields> ret(MEDCoupling::MEDFileFields::New());
1599 AppendFieldFromMeshes(ms,ret);
1603 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1605 std::vector<std::string> ret;
1607 while(pos!=std::string::npos)
1609 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1610 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1612 pos=fullFieldName.find_first_not_of(sep,curPos);
1618 * Here the non regression tests.
1619 * const char inp0[]="";
1620 * const char exp0[]="";
1621 * const char inp1[]="field";
1622 * const char exp1[]="field";
1623 * const char inp2[]="_________";
1624 * const char exp2[]="_________";
1625 * const char inp3[]="field_p";
1626 * const char exp3[]="field_p";
1627 * const char inp4[]="field__p";
1628 * const char exp4[]="field_p";
1629 * const char inp5[]="field_p__";
1630 * const char exp5[]="field_p";
1631 * const char inp6[]="field_p_";
1632 * const char exp6[]="field_p";
1633 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1634 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1635 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1636 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1637 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1638 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1639 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1640 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1642 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1644 const char SEP('_');
1645 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1647 return fullFieldName;//should never happen
1651 return fullFieldName;
1655 std::string ret(v[0]);
1656 for(std::size_t i=1;i<v.size();i++)
1661 { ret+=SEP; ret+=v[i]; }
1667 return fullFieldName;
1673 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1677 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1682 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1684 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1686 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1692 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1693 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1694 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1696 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1698 std::size_t sz(ts.size());
1699 bool isInHumanRange(true);
1701 for(std::size_t i=0;i<sz;i++)
1704 if(ts[i]<=-1e299 || ts[i]>=1e299)
1705 isInHumanRange=false;
1708 return processedUsingPairOfIds(tsPairs);
1710 return processedUsingPairOfIds(tsPairs);
1711 _postprocessed_time=ts;
1712 return getPostProcessedTime();
1717 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1719 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1721 std::size_t sz(ts.size());
1722 std::vector<double> ts2(sz);
1724 for(std::size_t i=0;i<sz;i++)
1729 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1732 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1734 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1738 void TimeKeeper::printSelf(std::ostream& oss) const
1740 std::size_t sz(_activated_ts.size());
1741 for(std::size_t i=0;i<sz;i++)
1743 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1747 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1749 std::size_t sz(_activated_ts.size());
1750 std::vector<bool> ret(sz);
1751 for(std::size_t i=0;i<sz;i++)
1753 ret[i]=_activated_ts[i].first;
1758 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1760 std::size_t sz(tsPairs.size());
1761 std::set<int> s0,s1;
1762 for(std::size_t i=0;i<sz;i++)
1763 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1766 _postprocessed_time.resize(sz);
1767 for(std::size_t i=0;i<sz;i++)
1768 _postprocessed_time[i]=(double)tsPairs[i].first;
1769 return getPostProcessedTime();
1773 _postprocessed_time.resize(sz);
1774 for(std::size_t i=0;i<sz;i++)
1775 _postprocessed_time[i]=(double)tsPairs[i].second;
1776 return getPostProcessedTime();
1778 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1779 _postprocessed_time.resize(sz);
1780 for(std::size_t i=0;i<sz;i++)
1781 _postprocessed_time[i]=(double)i;
1782 return getPostProcessedTime();
1785 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1787 _activated_ts.resize(maxNumberOfTS);
1788 for(int i=0;i<maxNumberOfTS;i++)
1790 std::ostringstream oss; oss << "000" << i;
1791 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());