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 "MEDCouplingMemArray.txx"
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 "vtkAOSDataArrayTemplate.h"
46 #include "vtkIdTypeArray.h"
47 #include "vtkDoubleArray.h"
48 #include "vtkIntArray.h"
49 #include "vtkLongArray.h"
51 #include "vtkLongLongArray.h"
53 #include "vtkFloatArray.h"
54 #include "vtkCellArray.h"
55 #include "vtkPointData.h"
56 #include "vtkFieldData.h"
57 #include "vtkCellData.h"
59 #include "vtkMutableDirectedGraph.h"
61 using namespace MEDCoupling;
63 const char MEDFileFieldRepresentationLeavesArrays::ZE_SEP[]="@@][@@";
65 const char MEDFileFieldRepresentationLeavesArrays::TS_STR[]="TS";
67 const char MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR[]="ComSup";
69 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME[]="FamilyIdCell";
71 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME[]="NumIdCell";
73 const char MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME[]="FamilyIdNode";
75 const char MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME[]="NumIdNode";
77 const char MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME[]="GlobalNodeIds";// WARNING DO NOT CHANGE IT BEFORE HAVING CHECKED IN PV SOURCES !
79 const char MEDFileFieldRepresentationTree::ROOT_OF_GRPS_IN_TREE[]="zeGrps";
81 const char MEDFileFieldRepresentationTree::ROOT_OF_FAM_IDS_IN_TREE[]="zeFamIds";
83 const char MEDFileFieldRepresentationTree::COMPO_STR_TO_LOCATE_MESH_DA[]="-@?|*_";
86 vtkIdTypeArray *ELGACmp::findOrCreate(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd, vtkDataSet *ds, bool& isNew, ExportedTinyInfo *internalInfo) const
88 vtkIdTypeArray *try0(isExisting(locsReallyUsed,vtkd));
97 return createNew<T>(globs,locsReallyUsed,vtkd,ds,internalInfo);
101 vtkIdTypeArray *ELGACmp::isExisting(const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd) const
103 std::vector< std::vector<std::string> >::iterator it(std::find(_loc_names.begin(),_loc_names.end(),locsReallyUsed));
104 if(it==_loc_names.end())
106 std::size_t pos(std::distance(_loc_names.begin(),it));
107 vtkIdTypeArray *ret(_elgas[pos]);
108 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
109 for(std::vector<std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it=_defs[pos].begin();it!=_defs[pos].end();it++)
111 key->Set(vtkd->GetInformation(),(*it).first,(*it).second);
113 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),ret->GetName());
118 vtkIdTypeArray *ELGACmp::createNew(const MEDCoupling::MEDFileFieldGlobsReal *globs, const std::vector<std::string>& locsReallyUsed, vtkDataArray *vtkd, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
120 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_DELETE;
121 std::vector< std::vector<std::string> > locNames(_loc_names);
122 std::vector<vtkIdTypeArray *> elgas(_elgas);
123 std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > defs;
125 std::vector< std::vector<std::string> >::const_iterator it(std::find(locNames.begin(),locNames.end(),locsReallyUsed));
126 if(it!=locNames.end())
127 throw INTERP_KERNEL::Exception("ELGACmp::createNew : Method is expected to be called after isExisting call ! Entry already exists !");
128 locNames.push_back(locsReallyUsed);
129 vtkIdTypeArray *elga(vtkIdTypeArray::New());
130 elga->SetNumberOfComponents(1);
131 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
132 std::map<unsigned char,int> m;
133 for(std::vector<std::string>::const_iterator it=locsReallyUsed.begin();it!=locsReallyUsed.end();it++)
135 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
136 const MEDFileFieldLoc& loc(globs->getLocalization((*it).c_str()));
137 INTERP_KERNEL::NormalizedCellType ct(loc.getGeoType());
138 unsigned char vtkType(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE[ct]);
139 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
140 int nbGaussPt(loc.getNbOfGaussPtPerCell()),nbPtsPerCell((int)cm.getNumberOfNodes()),dimLoc(loc.getDimension());
141 // 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.
142 std::vector<double> gsCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getGaussCoords()));
143 std::vector<double> refCoods2(INTERP_KERNEL::GaussInfo::NormalizeCoordinatesIfNecessary(ct,dimLoc,loc.getRefCoords()));
145 internalInfo->pushGaussAdditionnalInfo(vtkType,dimLoc,refCoods2,gsCoods2);
146 double *shape(new double[nbPtsPerCell*nbGaussPt]);
147 INTERP_KERNEL::GaussInfo calculator(ct,gsCoods2,nbGaussPt,refCoods2,nbPtsPerCell);
148 calculator.initLocalInfo();
149 const std::vector<double>& wgths(loc.getGaussWeights());
150 for(int i=0;i<nbGaussPt;i++)
152 const double *pt0(calculator.getFunctionValues(i));
153 if(ct!=INTERP_KERNEL::NORM_HEXA27)
154 std::copy(pt0,pt0+nbPtsPerCell,shape+nbPtsPerCell*i);
157 for(int j=0;j<27;j++)
158 shape[nbPtsPerCell*i+j]=pt0[MEDMeshMultiLev::HEXA27_PERM_ARRAY[j]];
161 m[vtkType]=nbGaussPt;
162 def->Initialize(vtkType,nbPtsPerCell,nbGaussPt,shape,const_cast<double *>(&wgths[0]));
164 key->Set(elga->GetInformation(),def,vtkType);
165 key->Set(vtkd->GetInformation(),def,vtkType);
166 defs.push_back(std::pair< vtkQuadratureSchemeDefinition *, unsigned char >(def,vtkType));
169 vtkIdType ncell(ds->GetNumberOfCells());
170 vtkIdType *pt(new vtkIdType[ncell]),offset(0);
171 for(vtkIdType cellId=0;cellId<ncell;cellId++)
173 vtkCell *cell(ds->GetCell(cellId));
174 vtkIdType delta(m[(unsigned char)cell->GetCellType()]);
178 elga->GetInformation()->Set(MEDUtilities::ELGA(),1);
179 elga->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
180 std::ostringstream oss; oss << "ELGA" << "@" << _loc_names.size();
181 std::string ossStr(oss.str());
182 elga->SetName(ossStr.c_str());
183 elga->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
184 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elga->GetName());
185 elgas.push_back(elga);
189 _defs.push_back(defs);
193 void ELGACmp::appendELGAIfAny(vtkDataSet *ds) const
195 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
196 ds->GetCellData()->AddArray(*it);
201 for(std::vector<vtkIdTypeArray *>::const_iterator it=_elgas.begin();it!=_elgas.end();it++)
203 for(std::vector< std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > > >::const_iterator it0=_defs.begin();it0!=_defs.end();it0++)
204 for(std::vector< std::pair< vtkQuadratureSchemeDefinition *, unsigned char > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
205 (*it1).first->Delete();
211 class MEDFileVTKTraits
214 typedef void VtkType;
219 class MEDFileVTKTraits<int>
222 typedef vtkIntArray VtkType;
223 typedef MEDCoupling::DataArrayInt MCType;
228 class MEDFileVTKTraits<long long>
230 class MEDFileVTKTraits<long>
235 typedef vtkLongLongArray VtkType;
237 typedef vtkLongArray VtkType;
239 typedef MEDCoupling::DataArrayInt64 MCType;
243 class MEDFileVTKTraits<float>
246 typedef vtkFloatArray VtkType;
247 typedef MEDCoupling::DataArrayFloat MCType;
251 class MEDFileVTKTraits<double>
254 typedef vtkDoubleArray VtkType;
255 typedef MEDCoupling::DataArrayDouble MCType;
258 typedef typename MEDFileVTKTraits<mcIdType>::VtkType vtkMCIdTypeArray;
262 void AssignDataPointerToVTK(typename MEDFileVTKTraits<T>::VtkType *vtkTab, typename MEDFileVTKTraits<T>::MCType *mcTab, bool noCpyNumNodes)
265 vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),1,vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE);
267 { vtkTab->SetArray(mcTab->getPointer(),mcTab->getNbOfElems(),0,vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_FREE); mcTab->accessToMemArray().setSpecificDeallocator(0); }
270 // here copy is always assumed.
271 template<class VTKT, class MCT>
272 void AssignDataPointerOther(VTKT *vtkTab, MCT *mcTab, vtkIdType nbElems)
274 typedef typename VTKT::ValueType VTKType;
275 if ( sizeof( VTKType ) == sizeof( typename MCT::Type ))
277 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(mcTab->getPointer()),nbElems,0,VTKT::VTK_DATA_ARRAY_FREE);
278 mcTab->accessToMemArray().setSpecificDeallocator(0);
282 VTKType* newArray = new VTKType[ nbElems ];
283 std::copy( mcTab->begin(), mcTab->begin() + nbElems, newArray );
284 vtkTab->SetVoidArray(reinterpret_cast<unsigned char *>(newArray),nbElems,0,VTKT::VTK_DATA_ARRAY_DELETE);
289 void AssignToFieldData(DataArray *vPtr, const MEDTimeReq *tr, vtkFieldData *att, const std::string& crudeName, bool noCpyNumNodes,
290 const std::vector<TypeOfField>& discs, const ELGACmp& elgaCmp, const MEDCoupling::MEDFileFieldGlobsReal *globs,
291 MEDFileAnyTypeField1TS *f1ts, vtkDataSet *ds, ExportedTinyInfo *internalInfo)
293 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<T>::VTK_DATA_ARRAY_DELETE;
294 typename MEDFileVTKTraits<T>::MCType *vi(static_cast<typename MEDFileVTKTraits<T>::MCType *>(vPtr));
295 typename MEDFileVTKTraits<T>::VtkType *vtkd(MEDFileVTKTraits<T>::VtkType::New());
296 vtkd->SetNumberOfComponents((int)vi->getNumberOfComponents());
297 for(unsigned int i=0;i<vi->getNumberOfComponents();i++)
298 vtkd->SetComponentName(i,vi->getVarOnComponent(i).c_str());
299 AssignDataPointerToVTK<T>(vtkd,vi,noCpyNumNodes);
300 std::string name(tr->buildName(crudeName));
301 vtkd->SetName(name.c_str());
304 if(discs[0]==ON_GAUSS_PT)
307 elgaCmp.findOrCreate<T>(globs,f1ts->getLocsReallyUsed(),vtkd,ds,tmp,internalInfo);
309 if(discs[0]==ON_GAUSS_NE)
311 vtkIdTypeArray *elno(vtkIdTypeArray::New());
312 elno->SetNumberOfComponents(1);
313 vtkIdType ncell(ds->GetNumberOfCells());
314 vtkIdType *pt(new vtkIdType[ncell]),offset(0);
315 std::set<int> cellTypes;
316 for(vtkIdType cellId=0;cellId<ncell;cellId++)
318 vtkCell *cell(ds->GetCell(cellId));
319 vtkIdType delta(cell->GetNumberOfPoints());
320 cellTypes.insert(cell->GetCellType());
324 elno->GetInformation()->Set(MEDUtilities::ELNO(),1);
325 elno->SetVoidArray(pt,ncell,0,VTK_DATA_ARRAY_DELETE);
326 std::string nameElno("ELNO"); nameElno+="@"; nameElno+=name;
327 elno->SetName(nameElno.c_str());
328 ds->GetCellData()->AddArray(elno);
329 vtkd->GetInformation()->Set(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME(),elno->GetName());
330 elno->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(),1);
332 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
333 for(std::set<int>::const_iterator it=cellTypes.begin();it!=cellTypes.end();it++)
335 const unsigned char *pos(std::find(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,*it));
336 if(pos==MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
338 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)std::distance(MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos));
339 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(ct));
340 int nbGaussPt(cm.getNumberOfNodes()),dim(cm.getDimension());
341 vtkQuadratureSchemeDefinition *def(vtkQuadratureSchemeDefinition::New());
342 double *shape(new double[nbGaussPt*nbGaussPt]);
344 const double *gsCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));//GetLocsFromGeometricType
345 const double *refCoords(MEDCouplingFieldDiscretizationGaussNE::GetRefCoordsFromGeometricType(ct,dummy));
346 const double *weights(MEDCouplingFieldDiscretizationGaussNE::GetWeightArrayFromGeometricType(ct,dummy));
347 std::vector<double> gsCoords2(gsCoords,gsCoords+nbGaussPt*dim),refCoords2(refCoords,refCoords+nbGaussPt*dim);
348 INTERP_KERNEL::GaussInfo calculator(ct,gsCoords2,nbGaussPt,refCoords2,nbGaussPt);
349 calculator.initLocalInfo();
350 for(int i=0;i<nbGaussPt;i++)
352 const double *pt0(calculator.getFunctionValues(i));
353 std::copy(pt0,pt0+nbGaussPt,shape+nbGaussPt*i);
355 def->Initialize(*it,nbGaussPt,nbGaussPt,shape,const_cast<double *>(weights));
357 key->Set(elno->GetInformation(),def,*it);
358 key->Set(vtkd->GetInformation(),def,*it);
368 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays():_id(-1)
372 MEDFileFieldRepresentationLeavesArrays::MEDFileFieldRepresentationLeavesArrays(const MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>& arr):MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>(arr),_activated(false),_id(-1)
374 std::vector< std::vector<MEDCoupling::TypeOfField> > typs((operator->())->getTypesOfFieldAvailable());
376 throw INTERP_KERNEL::Exception("There is a big internal problem in MEDLoader ! The field time spitting has failed ! A CRASH will occur soon !");
377 if(typs[0].size()!=1)
378 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 !");
379 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDiscretization> fd(MEDCouplingFieldDiscretization::New(typs[0][0]));
380 std::ostringstream oss2; oss2 << (operator->())->getName() << ZE_SEP << fd->getRepr();
384 MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeavesArrays::operator=(const MEDFileFieldRepresentationLeavesArrays& other)
386 MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS>::operator=(other);
389 _ze_name=other._ze_name;
390 _ze_full_name.clear();
394 void MEDFileFieldRepresentationLeavesArrays::setId(int& id) const
399 int MEDFileFieldRepresentationLeavesArrays::getId() const
404 std::string MEDFileFieldRepresentationLeavesArrays::getZeName() const
406 return _ze_full_name;
409 const char *MEDFileFieldRepresentationLeavesArrays::getZeNameC() const
411 return _ze_full_name.c_str();
414 void MEDFileFieldRepresentationLeavesArrays::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
416 vtkIdType refId(sil->AddChild(root,edge));
417 names.push_back(_ze_name);
419 if(MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(((operator->())->getInfo())))
421 sil->AddChild(refId,edge);
422 names.push_back(std::string());
426 void MEDFileFieldRepresentationLeavesArrays::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
428 std::ostringstream oss3; oss3 << tsName << "/" << meshName << "/" << comSupStr << "/" << _ze_name;
429 _ze_full_name=oss3.str();
432 bool MEDFileFieldRepresentationLeavesArrays::getStatus() const
437 bool MEDFileFieldRepresentationLeavesArrays::setStatus(bool status) const
439 bool ret(_activated!=status);
444 void MEDFileFieldRepresentationLeavesArrays::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshStruct *mst, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
446 //const int VTK_DATA_ARRAY_DELETE=vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE; // todo: unused
447 tr->setNumberOfTS((operator->())->getNumberOfTS());
449 for(int timeStepId=0;timeStepId<tr->size();timeStepId++,++(*tr))
451 MCAuto<MEDFileAnyTypeField1TS> f1ts((operator->())->getTimeStepAtPos(tr->getCurrent()));
452 MEDFileAnyTypeField1TS *f1tsPtr(f1ts);
453 MEDFileField1TS *f1tsPtrDbl(dynamic_cast<MEDFileField1TS *>(f1tsPtr));
454 MEDFileInt32Field1TS *f1tsPtrInt(dynamic_cast<MEDFileInt32Field1TS *>(f1tsPtr));
455 MEDFileInt64Field1TS *f1tsPtrInt64(dynamic_cast<MEDFileInt64Field1TS *>(f1tsPtr));
456 MEDFileFloatField1TS *f1tsPtrFloat(dynamic_cast<MEDFileFloatField1TS *>(f1tsPtr));
457 DataArray *crudeArr(0),*postProcessedArr(0);
459 crudeArr=f1tsPtrDbl->getUndergroundDataArray();
461 crudeArr=f1tsPtrInt->getUndergroundDataArray();
462 else if(f1tsPtrInt64)
463 crudeArr=f1tsPtrInt64->getUndergroundDataArray();
464 else if(f1tsPtrFloat)
465 crudeArr=f1tsPtrFloat->getUndergroundDataArray();
467 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64, FLOAT32 and INT32 fields are dealt for the moment !");
468 MEDFileField1TSStructItem fsst(MEDFileField1TSStructItem::BuildItemFrom(f1ts,mst));
469 f1ts->loadArraysIfNecessary();
470 MCAuto<DataArray> v(mml->buildDataArray(fsst,globs,crudeArr));
473 std::vector<TypeOfField> discs(f1ts->getTypesOfFieldAvailable());
475 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : internal error ! Number of spatial discretizations must be equal to one !");
476 vtkFieldData *att(0);
481 att=ds->GetCellData();
486 att=ds->GetPointData();
491 att=ds->GetFieldData();
496 att=ds->GetFieldData();
500 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only CELL and NODE, GAUSS_NE and GAUSS fields are available for the moment !");
504 AssignToFieldData<double>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
508 AssignToFieldData<int>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
510 else if(f1tsPtrFloat)
512 AssignToFieldData<float>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
514 else if(f1tsPtrInt64)
516 AssignToFieldData<Int64>(v,tr,att,f1ts->getName(),postProcessedArr==crudeArr,discs,_elga_cmp,globs,f1ts,ds,internalInfo);
519 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeavesArrays::appendFields : only FLOAT64 and INT32 fields are dealt for the moment ! Internal Error !");
523 void MEDFileFieldRepresentationLeavesArrays::appendELGAIfAny(vtkDataSet *ds) const
525 _elga_cmp.appendELGAIfAny(ds);
530 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves():_cached_ds(0)
534 MEDFileFieldRepresentationLeaves::MEDFileFieldRepresentationLeaves(const std::vector< MEDCoupling::MCAuto<MEDCoupling::MEDFileAnyTypeFieldMultiTS> >& arr,
535 const MEDCoupling::MCAuto<MEDCoupling::MEDFileFastCellSupportComparator>& fsp):_arrays(arr.size()),_fsp(fsp),_cached_ds(0)
537 for(std::size_t i=0;i<arr.size();i++)
538 _arrays[i]=MEDFileFieldRepresentationLeavesArrays(arr[i]);
541 MEDFileFieldRepresentationLeaves::~MEDFileFieldRepresentationLeaves()
544 _cached_ds->Delete();
547 bool MEDFileFieldRepresentationLeaves::empty() const
549 const MEDFileFastCellSupportComparator *fcscp(_fsp);
550 return fcscp==0 || _arrays.empty();
553 void MEDFileFieldRepresentationLeaves::setId(int& id) const
555 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
559 std::string MEDFileFieldRepresentationLeaves::getMeshName() const
561 return _arrays[0]->getMeshName();
564 int MEDFileFieldRepresentationLeaves::getNumberOfArrays() const
566 return (int)_arrays.size();
569 int MEDFileFieldRepresentationLeaves::getNumberOfTS() const
571 return _arrays[0]->getNumberOfTS();
574 void MEDFileFieldRepresentationLeaves::computeFullNameInLeaves(const std::string& tsName, const std::string& meshName, const std::string& comSupStr) const
576 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
577 (*it).computeFullNameInLeaves(tsName,meshName,comSupStr);
581 * \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.
582 * \param [in] meshName
588 void MEDFileFieldRepresentationLeaves::feedSIL(const MEDCoupling::MEDFileMeshes *ms, const std::string& meshName, vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
590 vtkIdType root2(sil->AddChild(root,edge));
591 names.push_back(std::string("Arrs"));
592 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
593 (*it).feedSIL(sil,root2,edge,names);
595 vtkIdType root3(sil->AddChild(root,edge));
596 names.push_back(std::string("InfoOnGeoType"));
597 const MEDCoupling::MEDFileMesh *m(0);
599 m=ms->getMeshWithName(meshName);
600 const MEDCoupling::MEDFileFastCellSupportComparator *fsp(_fsp);
601 if(!fsp || fsp->getNumberOfTS()==0)
603 std::vector< INTERP_KERNEL::NormalizedCellType > gts(fsp->getGeoTypesAt(0,m));
604 for(std::vector< INTERP_KERNEL::NormalizedCellType >::const_iterator it2=gts.begin();it2!=gts.end();it2++)
606 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(*it2));
607 std::string cmStr(cm.getRepr()); cmStr=cmStr.substr(5);//skip "NORM_"
608 sil->AddChild(root3,edge);
609 names.push_back(cmStr);
613 bool MEDFileFieldRepresentationLeaves::containId(int id) const
615 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
616 if((*it).getId()==id)
621 bool MEDFileFieldRepresentationLeaves::containZeName(const char *name, int& id) const
623 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
624 if((*it).getZeName()==name)
632 void MEDFileFieldRepresentationLeaves::dumpState(std::map<std::string,bool>& status) const
634 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
635 status[(*it).getZeName()]=(*it).getStatus();
638 bool MEDFileFieldRepresentationLeaves::isActivated() const
640 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
641 if((*it).getStatus())
646 void MEDFileFieldRepresentationLeaves::printMySelf(std::ostream& os) const
648 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it0=_arrays.begin();it0!=_arrays.end();it0++)
650 os << " - " << (*it0).getZeName() << " (";
651 if((*it0).getStatus())
655 os << ")" << std::endl;
659 void MEDFileFieldRepresentationLeaves::activateAllArrays() const
661 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
662 (*it).setStatus(true);
665 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationLeaves::getLeafArr(int id) const
667 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
668 if((*it).getId()==id)
670 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getLeafArr ! No such id !");
673 std::vector<double> MEDFileFieldRepresentationLeaves::getTimeSteps(const TimeKeeper& tk) const
676 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::getTimeSteps : the array size must be at least of size one !");
677 std::vector<double> ret;
678 std::vector< std::pair<int,int> > dtits(_arrays[0]->getTimeSteps(ret));
679 return tk.getTimeStepsRegardingPolicy(dtits,ret);
682 std::vector< std::pair<int,int> > MEDFileFieldRepresentationLeaves::getTimeStepsInCoarseMEDFileFormat(std::vector<double>& ts) const
685 return _arrays[0]->getTimeSteps(ts);
689 return std::vector< std::pair<int,int> >();
693 std::string MEDFileFieldRepresentationLeaves::getHumanReadableOverviewOfTS() const
695 std::ostringstream oss;
696 oss << _arrays[0]->getNumberOfTS() << " time steps [" << _arrays[0]->getDtUnit() << "]\n(";
697 std::vector<double> ret1;
698 std::vector< std::pair<int,int> > ret2(getTimeStepsInCoarseMEDFileFormat(ret1));
699 std::size_t sz(ret1.size());
700 for(std::size_t i=0;i<sz;i++)
702 oss << ret1[i] << " (" << ret2[i].first << "," << ret2[i].second << ")";
705 std::string tmp(oss.str());
706 if(tmp.size()>200 && i!=sz-1)
716 void MEDFileFieldRepresentationLeaves::appendFields(const MEDTimeReq *tr, const MEDCoupling::MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDMeshMultiLev *mml, const MEDCoupling::MEDFileMeshes *meshes, vtkDataSet *ds, ExportedTinyInfo *internalInfo) const
719 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::appendFields : internal error !");
720 MCAuto<MEDFileMeshStruct> mst(MEDFileMeshStruct::New(meshes->getMeshWithName(_arrays[0]->getMeshName().c_str())));
721 for(std::vector<MEDFileFieldRepresentationLeavesArrays>::const_iterator it=_arrays.begin();it!=_arrays.end();it++)
722 if((*it).getStatus())
724 (*it).appendFields(tr,globs,mml,mst,ds,internalInfo);
725 (*it).appendELGAIfAny(ds);
729 vtkUnstructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationUnstructured(MEDUMeshMultiLev *mm) const
731 DataArrayDouble *coordsMC(0);
732 DataArrayByte *typesMC(0);
733 DataArrayIdType *cellLocationsMC(0),*cellsMC(0),*faceLocationsMC(0),*facesMC(0);
734 bool statusOfCoords(mm->buildVTUArrays(coordsMC,typesMC,cellLocationsMC,cellsMC,faceLocationsMC,facesMC));
735 MCAuto<DataArrayDouble> coordsSafe(coordsMC);
736 MCAuto<DataArrayByte> typesSafe(typesMC);
737 MCAuto<DataArrayIdType> cellLocationsSafe(cellLocationsMC),cellsSafe(cellsMC),faceLocationsSafe(faceLocationsMC),facesSafe(facesMC);
739 vtkIdType nbOfCells(typesSafe->getNbOfElems());
740 vtkUnstructuredGrid *ret(vtkUnstructuredGrid::New());
741 vtkUnsignedCharArray *cellTypes(vtkUnsignedCharArray::New());
742 AssignDataPointerOther<vtkUnsignedCharArray,DataArrayByte>(cellTypes,typesSafe,nbOfCells);
743 vtkIdTypeArray *cellLocations(vtkIdTypeArray::New());
744 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(cellLocations,cellLocationsSafe,nbOfCells);
745 vtkCellArray *cells(vtkCellArray::New());
746 vtkIdTypeArray *cells2(vtkIdTypeArray::New());
747 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(cells2,cellsSafe,cellsSafe->getNbOfElems());
748 cells->SetCells(nbOfCells,cells2);
750 if(faceLocationsMC!=0 && facesMC!=0)
752 vtkIdTypeArray *faces(vtkIdTypeArray::New());
753 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(faces,facesSafe,facesSafe->getNbOfElems());
754 vtkIdTypeArray *faceLocations(vtkIdTypeArray::New());
755 AssignDataPointerOther<vtkIdTypeArray,DataArrayIdType>(faceLocations,faceLocationsSafe,faceLocationsSafe->getNbOfElems());
756 ret->SetCells(cellTypes,cellLocations,cells,faceLocations,faces);
757 faceLocations->Delete();
761 ret->SetCells(cellTypes,cellLocations,cells);
763 cellLocations->Delete();
765 vtkPoints *pts(vtkPoints::New());
766 vtkDoubleArray *pts2(vtkDoubleArray::New());
767 pts2->SetNumberOfComponents(3);
768 AssignDataPointerToVTK<double>(pts2,coordsSafe,statusOfCoords);
777 vtkRectilinearGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCartesian(MEDCoupling::MEDCMeshMultiLev *mm) const
780 std::vector< DataArrayDouble * > arrs(mm->buildVTUArrays(isInternal));
781 vtkDoubleArray *vtkTmp(0);
782 vtkRectilinearGrid *ret(vtkRectilinearGrid::New());
783 std::size_t dim(arrs.size());
785 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCartesian : dimension must be in [1,3] !");
786 int sizePerAxe[3]={1,1,1};
787 sizePerAxe[0]=arrs[0]->getNbOfElems();
789 sizePerAxe[1]=arrs[1]->getNbOfElems();
791 sizePerAxe[2]=arrs[2]->getNbOfElems();
792 ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
793 vtkTmp=vtkDoubleArray::New();
794 vtkTmp->SetNumberOfComponents(1);
795 AssignDataPointerToVTK<double>(vtkTmp,arrs[0],isInternal);
796 ret->SetXCoordinates(vtkTmp);
801 vtkTmp=vtkDoubleArray::New();
802 vtkTmp->SetNumberOfComponents(1);
803 AssignDataPointerToVTK<double>(vtkTmp,arrs[1],isInternal);
804 ret->SetYCoordinates(vtkTmp);
810 vtkTmp=vtkDoubleArray::New();
811 vtkTmp->SetNumberOfComponents(1);
812 AssignDataPointerToVTK<double>(vtkTmp,arrs[2],isInternal);
813 ret->SetZCoordinates(vtkTmp);
820 vtkStructuredGrid *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolationCurveLinear(MEDCoupling::MEDCurveLinearMeshMultiLev *mm) const
822 int meshStr[3]={1,1,1};
823 DataArrayDouble *coords(0);
824 std::vector<mcIdType> nodeStrct;
826 mm->buildVTUArrays(coords,nodeStrct,isInternal);
827 std::size_t dim(nodeStrct.size());
829 throw INTERP_KERNEL::Exception("buildVTKInstanceNoTimeInterpolationCurveLinear : dimension must be in [1,3] !");
830 meshStr[0]=nodeStrct[0];
832 meshStr[1]=nodeStrct[1];
834 meshStr[2]=nodeStrct[2];
835 vtkStructuredGrid *ret(vtkStructuredGrid::New());
836 ret->SetDimensions(meshStr[0],meshStr[1],meshStr[2]);
837 vtkDoubleArray *da(vtkDoubleArray::New());
838 da->SetNumberOfComponents(3);
839 if(coords->getNumberOfComponents()==3)
840 AssignDataPointerToVTK<double>(da,coords,isInternal);//if isIntenal==True VTK has not the ownership of double * because MEDLoader main struct has it !
843 MCAuto<DataArrayDouble> coords2(coords->changeNbOfComponents(3,0.));
844 AssignDataPointerToVTK<double>(da,coords2,false);//let VTK deal with double *
847 vtkPoints *points=vtkPoints::New();
848 ret->SetPoints(points);
855 vtkDataSet *MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation(const MEDTimeReq *tr, const MEDFileFieldGlobsReal *globs, const MEDCoupling::MEDFileMeshes *meshes, ExportedTinyInfo *internalInfo) const
858 //_fsp->isDataSetSupportEqualToThePreviousOne(i,globs);
859 MCAuto<MEDMeshMultiLev> mml(_fsp->buildFromScratchDataSetSupport(0,globs));//0=timestep Id. Make the hypothesis that support does not change
860 MCAuto<MEDMeshMultiLev> mml2(mml->prepare());
861 MEDMeshMultiLev *ptMML2(mml2);
864 MEDUMeshMultiLev *ptUMML2(dynamic_cast<MEDUMeshMultiLev *>(ptMML2));
865 MEDCMeshMultiLev *ptCMML2(dynamic_cast<MEDCMeshMultiLev *>(ptMML2));
866 MEDCurveLinearMeshMultiLev *ptCLMML2(dynamic_cast<MEDCurveLinearMeshMultiLev *>(ptMML2));
870 ret=buildVTKInstanceNoTimeInterpolationUnstructured(ptUMML2);
874 ret=buildVTKInstanceNoTimeInterpolationCartesian(ptCMML2);
878 ret=buildVTKInstanceNoTimeInterpolationCurveLinear(ptCLMML2);
881 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationLeaves::buildVTKInstanceNoTimeInterpolation : unrecognized mesh ! Supported for the moment unstructured, cartesian, curvelinear !");
882 _cached_ds=ret->NewInstance();
883 _cached_ds->ShallowCopy(ret);
887 ret=_cached_ds->NewInstance();
888 ret->ShallowCopy(_cached_ds);
891 appendFields(tr,globs,mml,meshes,ret,internalInfo);
892 // The arrays links to mesh
893 DataArrayIdType *famCells(0),*numCells(0);
894 bool noCpyFamCells(false),noCpyNumCells(false);
895 ptMML2->retrieveFamilyIdsOnCells(famCells,noCpyFamCells);
898 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
899 vtkTab->SetNumberOfComponents(1);
900 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME);
901 AssignDataPointerToVTK<mcIdType>(vtkTab,famCells,noCpyFamCells);
902 ret->GetCellData()->AddArray(vtkTab);
906 ptMML2->retrieveNumberIdsOnCells(numCells,noCpyNumCells);
909 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
910 vtkTab->SetNumberOfComponents(1);
911 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_CELL_NAME);
912 AssignDataPointerToVTK<mcIdType>(vtkTab,numCells,noCpyNumCells);
913 ret->GetCellData()->AddArray(vtkTab);
917 // The arrays links to mesh
918 DataArrayIdType *famNodes(0),*numNodes(0);
919 bool noCpyFamNodes(false),noCpyNumNodes(false);
920 ptMML2->retrieveFamilyIdsOnNodes(famNodes,noCpyFamNodes);
923 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
924 vtkTab->SetNumberOfComponents(1);
925 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME);
926 AssignDataPointerToVTK<mcIdType>(vtkTab,famNodes,noCpyFamNodes);
927 ret->GetPointData()->AddArray(vtkTab);
931 ptMML2->retrieveNumberIdsOnNodes(numNodes,noCpyNumNodes);
934 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
935 vtkTab->SetNumberOfComponents(1);
936 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::NUM_ID_NODE_NAME);
937 AssignDataPointerToVTK<mcIdType>(vtkTab,numNodes,noCpyNumNodes);
938 ret->GetPointData()->AddArray(vtkTab);
942 // Global Node Ids if any ! (In // mode)
943 DataArrayIdType *gni(ptMML2->retrieveGlobalNodeIdsIfAny());
946 vtkMCIdTypeArray *vtkTab(vtkMCIdTypeArray::New());
947 vtkTab->SetNumberOfComponents(1);
948 vtkTab->SetName(MEDFileFieldRepresentationLeavesArrays::GLOBAL_NODE_ID_NAME);
949 AssignDataPointerToVTK<mcIdType>(vtkTab,gni,false);
950 ret->GetPointData()->AddArray(vtkTab);
957 //////////////////////
959 MEDFileFieldRepresentationTree::MEDFileFieldRepresentationTree()
963 int MEDFileFieldRepresentationTree::getNumberOfLeavesArrays() const
966 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
967 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
968 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
969 ret+=(*it2).getNumberOfArrays();
973 void MEDFileFieldRepresentationTree::assignIds() const
976 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
977 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
978 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
982 void MEDFileFieldRepresentationTree::computeFullNameInLeaves() const
984 std::size_t it0Cnt(0);
985 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
987 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
988 std::string tsName(oss.str());
989 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
991 std::string meshName((*it1)[0].getMeshName());
992 std::size_t it2Cnt(0);
993 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
995 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
996 std::string comSupStr(oss2.str());
997 (*it2).computeFullNameInLeaves(tsName,meshName,comSupStr);
1003 void MEDFileFieldRepresentationTree::activateTheFirst() const
1005 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1006 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1007 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1009 (*it2).activateAllArrays();
1014 void MEDFileFieldRepresentationTree::feedSIL(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1016 std::size_t it0Cnt(0);
1017 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,it0Cnt++)
1019 vtkIdType InfoOnTSId(sil->AddChild(root,edge));
1020 names.push_back((*it0)[0][0].getHumanReadableOverviewOfTS());
1022 vtkIdType NbOfTSId(sil->AddChild(InfoOnTSId,edge));
1023 std::vector<double> ts;
1024 std::vector< std::pair<int,int> > dtits((*it0)[0][0].getTimeStepsInCoarseMEDFileFormat(ts));
1025 std::size_t nbOfTS(dtits.size());
1026 std::ostringstream oss3; oss3 << nbOfTS;
1027 names.push_back(oss3.str());
1028 for(std::size_t i=0;i<nbOfTS;i++)
1030 std::ostringstream oss4; oss4 << dtits[i].first;
1031 vtkIdType DtId(sil->AddChild(NbOfTSId,edge));
1032 names.push_back(oss4.str());
1033 std::ostringstream oss5; oss5 << dtits[i].second;
1034 vtkIdType ItId(sil->AddChild(DtId,edge));
1035 names.push_back(oss5.str());
1036 std::ostringstream oss6; oss6 << ts[i];
1037 sil->AddChild(ItId,edge);
1038 names.push_back(oss6.str());
1041 std::ostringstream oss; oss << MEDFileFieldRepresentationLeavesArrays::TS_STR << it0Cnt;
1042 std::string tsName(oss.str());
1043 vtkIdType typeId0(sil->AddChild(root,edge));
1044 names.push_back(tsName);
1045 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1047 std::string meshName((*it1)[0].getMeshName());
1048 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1049 names.push_back(meshName);
1050 std::size_t it2Cnt(0);
1051 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,it2Cnt++)
1053 std::ostringstream oss2; oss2 << MEDFileFieldRepresentationLeavesArrays::COM_SUP_STR << it2Cnt;
1054 std::string comSupStr(oss2.str());
1055 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1056 names.push_back(comSupStr);
1057 (*it2).feedSIL(_ms,meshName,sil,typeId2,edge,names);
1063 std::string MEDFileFieldRepresentationTree::getActiveMeshName() const
1065 int dummy0(0),dummy1(0),dummy2(0);
1066 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1067 return leaf.getMeshName();
1070 std::string MEDFileFieldRepresentationTree::feedSILForFamsAndGrps(vtkMutableDirectedGraph* sil, vtkIdType root, vtkVariantArray *edge, std::vector<std::string>& names) const
1072 int dummy0(0),dummy1(0),dummy2(0);
1073 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(dummy0,dummy1,dummy2));
1074 std::string ret(leaf.getMeshName());
1077 for(;i<_ms->getNumberOfMeshes();i++)
1079 m=_ms->getMeshAtPos(i);
1080 if(m->getName()==ret)
1083 if(i==_ms->getNumberOfMeshes())
1084 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::feedSILForFamsAndGrps : internal error #0 !");
1085 vtkIdType typeId0(sil->AddChild(root,edge));
1086 names.push_back(m->getName());
1088 vtkIdType typeId1(sil->AddChild(typeId0,edge));
1089 names.push_back(std::string(ROOT_OF_GRPS_IN_TREE));
1090 std::vector<std::string> grps(m->getGroupsNames());
1091 for(std::vector<std::string>::const_iterator it0=grps.begin();it0!=grps.end();it0++)
1093 vtkIdType typeId2(sil->AddChild(typeId1,edge));
1094 names.push_back(*it0);
1095 std::vector<std::string> famsOnGrp(m->getFamiliesOnGroup((*it0).c_str()));
1096 for(std::vector<std::string>::const_iterator it1=famsOnGrp.begin();it1!=famsOnGrp.end();it1++)
1098 sil->AddChild(typeId2,edge);
1099 names.push_back((*it1).c_str());
1103 vtkIdType typeId11(sil->AddChild(typeId0,edge));
1104 names.push_back(std::string(ROOT_OF_FAM_IDS_IN_TREE));
1105 std::vector<std::string> fams(m->getFamiliesNames());
1106 for(std::vector<std::string>::const_iterator it00=fams.begin();it00!=fams.end();it00++)
1108 sil->AddChild(typeId11,edge);
1109 int famId(m->getFamilyId((*it00).c_str()));
1110 std::ostringstream oss; oss << (*it00) << MEDFileFieldRepresentationLeavesArrays::ZE_SEP << famId;
1111 names.push_back(oss.str());
1116 const MEDFileFieldRepresentationLeavesArrays& MEDFileFieldRepresentationTree::getLeafArr(int id) const
1118 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1119 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1120 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1121 if((*it2).containId(id))
1122 return (*it2).getLeafArr(id);
1123 throw INTERP_KERNEL::Exception("Internal error in MEDFileFieldRepresentationTree::getLeafArr !");
1126 std::string MEDFileFieldRepresentationTree::getNameOf(int id) const
1128 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1129 return elt.getZeName();
1132 const char *MEDFileFieldRepresentationTree::getNameOfC(int id) const
1134 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1135 return elt.getZeNameC();
1138 bool MEDFileFieldRepresentationTree::getStatusOf(int id) const
1140 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1141 return elt.getStatus();
1144 int MEDFileFieldRepresentationTree::getIdHavingZeName(const char *name) const
1147 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1148 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1149 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1150 if((*it2).containZeName(name,ret))
1152 std::ostringstream msg; msg << "MEDFileFieldRepresentationTree::getIdHavingZeName : No such a name \"" << name << "\" !";
1153 throw INTERP_KERNEL::Exception(msg.str().c_str());
1156 bool MEDFileFieldRepresentationTree::changeStatusOfAndUpdateToHaveCoherentVTKDataSet(int id, bool status) const
1158 const MEDFileFieldRepresentationLeavesArrays& elt(getLeafArr(id));
1159 bool ret(elt.setStatus(status));//to be implemented
1163 int MEDFileFieldRepresentationTree::getMaxNumberOfTimeSteps() const
1166 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1167 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1168 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1169 ret=std::max(ret,(*it2).getNumberOfTS());
1176 void MEDFileFieldRepresentationTree::loadInMemory(MEDCoupling::MEDFileFields *fields, MEDCoupling::MEDFileMeshes *meshes)
1178 _fields=fields; _ms=meshes;
1179 if(_fields.isNotNull())
1184 if(_fields.isNull())
1186 _fields=BuildFieldFromMeshes(_ms);
1190 AppendFieldFromMeshes(_ms,_fields);
1192 _ms->cartesianizeMe();
1193 _fields->removeFieldsWithoutAnyTimeStep();
1194 std::vector<std::string> meshNames(_ms->getMeshesNames());
1195 std::vector< MCAuto<MEDFileFields> > fields_per_mesh(meshNames.size());
1196 for(std::size_t i=0;i<meshNames.size();i++)
1198 fields_per_mesh[i]=_fields->partOfThisLyingOnSpecifiedMeshName(meshNames[i].c_str());
1200 std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS > > allFMTSLeavesToDisplaySafe;
1201 for(std::vector< MCAuto<MEDFileFields> >::const_iterator fields=fields_per_mesh.begin();fields!=fields_per_mesh.end();fields++)
1203 for(int j=0;j<(*fields)->getNumberOfFields();j++)
1205 MCAuto<MEDFileAnyTypeFieldMultiTS> fmts((*fields)->getFieldAtPos((int)j));
1206 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > tmp(fmts->splitDiscretizations());
1208 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::const_iterator it=tmp.begin();it!=tmp.end();it++)
1210 if(!(*it)->presenceOfMultiDiscPerGeoType())
1211 allFMTSLeavesToDisplaySafe.push_back(*it);
1213 {// The case of some parts of field have more than one discretization per geo type.
1214 std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > > subTmp((*it)->splitMultiDiscrPerGeoTypes());
1215 std::size_t it0Cnt(0);
1216 for(std::vector< MCAuto< MEDFileAnyTypeFieldMultiTS > >::iterator it0=subTmp.begin();it0!=subTmp.end();it0++,it0Cnt++)//not const because setName
1218 std::ostringstream oss; oss << (*it0)->getName() << "_" << std::setfill('M') << std::setw(3) << it0Cnt;
1219 (*it0)->setName(oss.str());
1220 allFMTSLeavesToDisplaySafe.push_back(*it0);
1227 std::vector< MEDFileAnyTypeFieldMultiTS *> allFMTSLeavesToDisplay(allFMTSLeavesToDisplaySafe.size());
1228 for(std::size_t i=0;i<allFMTSLeavesToDisplaySafe.size();i++)
1230 allFMTSLeavesToDisplay[i]=allFMTSLeavesToDisplaySafe[i];
1232 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > allFMTSLeavesPerTimeSeries(MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries(allFMTSLeavesToDisplay));
1233 // memory safety part
1234 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > allFMTSLeavesPerTimeSeriesSafe(allFMTSLeavesPerTimeSeries.size());
1235 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeries.size();j++)
1237 allFMTSLeavesPerTimeSeriesSafe[j].resize(allFMTSLeavesPerTimeSeries[j].size());
1238 for(std::size_t k=0;k<allFMTSLeavesPerTimeSeries[j].size();k++)
1240 allFMTSLeavesPerTimeSeries[j][k]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitIntoCommonTimeSeries do not increments the counter
1241 allFMTSLeavesPerTimeSeriesSafe[j][k]=allFMTSLeavesPerTimeSeries[j][k];
1244 // end of memory safety part
1245 // 1st : timesteps, 2nd : meshName, 3rd : common support
1246 this->_data_structure.resize(allFMTSLeavesPerTimeSeriesSafe.size());
1247 for(std::size_t i=0;i<allFMTSLeavesPerTimeSeriesSafe.size();i++)
1249 std::vector< std::string > meshNamesLoc;
1250 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > splitByMeshName;
1251 for(std::size_t j=0;j<allFMTSLeavesPerTimeSeriesSafe[i].size();j++)
1253 std::string meshName(allFMTSLeavesPerTimeSeriesSafe[i][j]->getMeshName());
1254 std::vector< std::string >::iterator it(std::find(meshNamesLoc.begin(),meshNamesLoc.end(),meshName));
1255 if(it==meshNamesLoc.end())
1257 meshNamesLoc.push_back(meshName);
1258 splitByMeshName.resize(splitByMeshName.size()+1);
1259 splitByMeshName.back().push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1262 splitByMeshName[std::distance(meshNamesLoc.begin(),it)].push_back(allFMTSLeavesPerTimeSeriesSafe[i][j]);
1264 _data_structure[i].resize(meshNamesLoc.size());
1265 for(std::size_t j=0;j<splitByMeshName.size();j++)
1267 std::vector< MCAuto<MEDFileFastCellSupportComparator> > fsp;
1268 std::vector< MEDFileAnyTypeFieldMultiTS *> sbmn(splitByMeshName[j].size());
1269 for(std::size_t k=0;k<splitByMeshName[j].size();k++)
1270 sbmn[k]=splitByMeshName[j][k];
1271 //getMeshWithName does not return a newly allocated object ! It is a true get* method !
1272 std::vector< std::vector<MEDFileAnyTypeFieldMultiTS *> > commonSupSplit(MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport(sbmn,_ms->getMeshWithName(meshNamesLoc[j].c_str()),fsp));
1273 std::vector< std::vector< MCAuto<MEDFileAnyTypeFieldMultiTS> > > commonSupSplitSafe(commonSupSplit.size());
1274 this->_data_structure[i][j].resize(commonSupSplit.size());
1275 for(std::size_t k=0;k<commonSupSplit.size();k++)
1277 commonSupSplitSafe[k].resize(commonSupSplit[k].size());
1278 for(std::size_t l=0;l<commonSupSplit[k].size();l++)
1280 commonSupSplit[k][l]->incrRef();//because MEDFileAnyTypeFieldMultiTS::SplitPerCommonSupport does not increment pointers !
1281 commonSupSplitSafe[k][l]=commonSupSplit[k][l];
1284 for(std::size_t k=0;k<commonSupSplit.size();k++)
1285 this->_data_structure[i][j][k]=MEDFileFieldRepresentationLeaves(commonSupSplitSafe[k],fsp[k]);
1288 this->removeEmptyLeaves();
1290 this->computeFullNameInLeaves();
1293 void MEDFileFieldRepresentationTree::loadMainStructureOfFile(const char *fileName, int iPart, int nbOfParts)
1295 MCAuto<MEDFileMeshes> ms;
1296 MCAuto<MEDFileFields> fields;
1298 if((iPart==-1 && nbOfParts==-1) || (iPart==0 && nbOfParts==1))
1300 MCAuto<MEDFileMeshSupports> msups(MEDFileMeshSupports::New(fileName));
1301 MCAuto<MEDFileStructureElements> mse(MEDFileStructureElements::New(fileName,msups));
1302 ms=MEDFileMeshes::New(fileName);
1303 fields=MEDFileFields::NewWithDynGT(fileName,mse,false);//false is important to not read the values
1304 if(ms->presenceOfStructureElements())
1306 fields->loadArrays();
1307 fields->blowUpSE(ms,mse);
1309 int nbMeshes(ms->getNumberOfMeshes());
1310 for(int i=0;i<nbMeshes;i++)
1312 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1313 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1315 tmp2->forceComputationOfParts();
1320 #ifdef MEDREADER_USE_MPI
1321 ms=ParaMEDFileMeshes::New(iPart,nbOfParts,fileName);
1322 int nbMeshes(ms->getNumberOfMeshes());
1323 for(int i=0;i<nbMeshes;i++)
1325 MEDCoupling::MEDFileMesh *tmp(ms->getMeshAtPos(i));
1326 MEDCoupling::MEDFileUMesh *tmp2(dynamic_cast<MEDCoupling::MEDFileUMesh *>(tmp));
1328 MCAuto<DataArrayIdType> tmp3(tmp2->zipCoords());
1330 fields=MEDFileFields::LoadPartOf(fileName,false,ms);//false is important to not read the values
1332 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::loadMainStructureOfFile : request for iPart/nbOfParts=" << iPart << "/" << nbOfParts << " whereas Plugin not compiled with MPI !";
1333 throw INTERP_KERNEL::Exception(oss.str().c_str());
1337 loadInMemory(fields,ms);
1340 void MEDFileFieldRepresentationTree::removeEmptyLeaves()
1342 std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > > newSD;
1343 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1345 std::vector< std::vector< MEDFileFieldRepresentationLeaves > > newSD0;
1346 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1348 std::vector< MEDFileFieldRepresentationLeaves > newSD1;
1349 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1351 newSD1.push_back(*it2);
1353 newSD0.push_back(newSD1);
1356 newSD.push_back(newSD0);
1360 bool MEDFileFieldRepresentationTree::IsFieldMeshRegardingInfo(const std::vector<std::string>& compInfos)
1362 if(compInfos.size()!=1)
1364 return compInfos[0]==COMPO_STR_TO_LOCATE_MESH_DA;
1367 std::string MEDFileFieldRepresentationTree::getDftMeshName() const
1369 return _data_structure[0][0][0].getMeshName();
1372 std::vector<double> MEDFileFieldRepresentationTree::getTimeSteps(int& lev0, const TimeKeeper& tk) const
1375 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1376 return leaf.getTimeSteps(tk);
1379 vtkDataSet *MEDFileFieldRepresentationTree::buildVTKInstance(bool isStdOrMode, double timeReq, std::string& meshName, const TimeKeeper& tk, ExportedTinyInfo *internalInfo) const
1382 const MEDFileFieldRepresentationLeaves& leaf(getTheSingleActivated(lev0,lev1,lev2));
1383 meshName=leaf.getMeshName();
1384 std::vector<double> ts(leaf.getTimeSteps(tk));
1385 std::size_t zeTimeId(0);
1388 std::vector<double> ts2(ts.size());
1389 std::transform(ts.begin(),ts.end(),ts2.begin(),std::bind(std::plus<double>(),std::placeholders::_1,-timeReq));
1390 std::transform(ts2.begin(),ts2.end(),ts2.begin(),[](double c) {return fabs(c);});
1391 zeTimeId=std::distance(ts2.begin(),std::find_if(ts2.begin(),ts2.end(),std::bind(std::less<double>(),std::placeholders::_1,1e-14)));
1394 if(zeTimeId==ts.size())
1395 zeTimeId=std::distance(ts.begin(),std::find(ts.begin(),ts.end(),timeReq));
1396 if(zeTimeId==ts.size())
1397 {//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.
1398 //In this case the default behaviour is taken. Keep the highest time step in this lower than timeReq.
1400 double valAttachedToPos(-std::numeric_limits<double>::max());
1401 for(std::size_t i=0;i<ts.size();i++)
1405 if(ts[i]>valAttachedToPos)
1408 valAttachedToPos=ts[i];
1413 {// timeReq is lower than all time steps (ts). So let's keep the lowest time step greater than timeReq.
1414 valAttachedToPos=std::numeric_limits<double>::max();
1415 for(std::size_t i=0;i<ts.size();i++)
1417 if(ts[i]<valAttachedToPos)
1420 valAttachedToPos=ts[i];
1425 std::ostringstream oss; oss.precision(15); oss << "request for time " << timeReq << " but not in ";
1426 std::copy(ts.begin(),ts.end(),std::ostream_iterator<double>(oss,","));
1427 oss << " ! Keep time " << valAttachedToPos << " at pos #" << zeTimeId;
1428 std::cerr << oss.str() << std::endl;
1432 tr=new MEDStdTimeReq((int)zeTimeId);
1434 tr=new MEDModeTimeReq(tk.getTheVectOfBool(),tk.getPostProcessedTime());
1435 vtkDataSet *ret(leaf.buildVTKInstanceNoTimeInterpolation(tr,_fields,_ms,internalInfo));
1440 const MEDFileFieldRepresentationLeaves& MEDFileFieldRepresentationTree::getTheSingleActivated(int& lev0, int& lev1, int& lev2) const
1442 int nbOfActivated(0);
1443 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1444 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1445 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1446 if((*it2).isActivated())
1448 if(nbOfActivated!=1)
1450 std::ostringstream oss; oss << "MEDFileFieldRepresentationTree::getTheSingleActivated : Only one leaf must be activated ! Having " << nbOfActivated << " !";
1451 throw INTERP_KERNEL::Exception(oss.str().c_str());
1453 int i0(0),i1(0),i2(0);
1454 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i0++)
1455 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,i1++)
1456 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,i2++)
1457 if((*it2).isActivated())
1459 lev0=i0; lev1=i1; lev2=i2;
1462 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::getTheSingleActivated : Internal error !");
1465 void MEDFileFieldRepresentationTree::printMySelf(std::ostream& os) const
1468 os << "#############################################" << std::endl;
1469 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++,i++)
1472 os << "TS" << i << std::endl;
1473 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++,j++)
1476 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++,k++)
1479 os << " " << (*it2).getMeshName() << std::endl;
1480 os << " Comp" << k << std::endl;
1481 (*it2).printMySelf(os);
1485 os << "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$" << std::endl;
1488 std::map<std::string,bool> MEDFileFieldRepresentationTree::dumpState() const
1490 std::map<std::string,bool> ret;
1491 for(std::vector< std::vector< std::vector< MEDFileFieldRepresentationLeaves > > >::const_iterator it0=_data_structure.begin();it0!=_data_structure.end();it0++)
1492 for(std::vector< std::vector< MEDFileFieldRepresentationLeaves > >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
1493 for(std::vector< MEDFileFieldRepresentationLeaves >::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
1494 (*it2).dumpState(ret);
1498 void MEDFileFieldRepresentationTree::AppendFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms, MEDCoupling::MEDFileFields *ret)
1501 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::AppendFieldFromMeshes : internal error ! NULL ret !");
1502 for(int i=0;i<ms->getNumberOfMeshes();i++)
1504 MEDFileMesh *mm(ms->getMeshAtPos(i));
1505 std::vector<int> levs(mm->getNonEmptyLevels());
1506 MEDCoupling::MCAuto<MEDCoupling::MEDFileField1TS> f1tsMultiLev(MEDCoupling::MEDFileField1TS::New());
1507 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
1510 for(std::vector<int>::const_iterator it=levs.begin();it!=levs.end();it++)
1512 std::vector<INTERP_KERNEL::NormalizedCellType> gts(mmu->getGeoTypesAtLevel(*it));
1513 for(std::vector<INTERP_KERNEL::NormalizedCellType>::const_iterator gt=gts.begin();gt!=gts.end();gt++)
1515 MEDCoupling::MEDCouplingMesh *m(mmu->getDirectUndergroundSingleGeoTypeMesh(*gt));
1516 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1518 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1519 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1522 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1523 f1tsMultiLev->setFieldNoProfileSBT(f);
1528 std::vector<int> levsExt(mm->getNonEmptyLevelsExt());
1529 if(levsExt.size()==levs.size()+1)
1531 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(1));
1532 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_NODES));
1534 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(m->getNumberOfNodes());
1535 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1536 arr->iota(); f->setArray(arr);
1537 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1538 f1tsMultiLev->setFieldNoProfileSBT(f);
1546 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingMesh> m(mm->getMeshAtLevel(0));
1547 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingFieldDouble> f(MEDCoupling::MEDCouplingFieldDouble::New(MEDCoupling::ON_CELLS));
1549 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> arr(MEDCoupling::DataArrayDouble::New()); arr->alloc(f->getNumberOfTuplesExpected());
1550 arr->setInfoOnComponent(0,std::string(COMPO_STR_TO_LOCATE_MESH_DA));
1553 f->setName(BuildAUniqueArrayNameForMesh(mm->getName(),ret));
1554 f1tsMultiLev->setFieldNoProfileSBT(f);
1557 MEDCoupling::MCAuto<MEDCoupling::MEDFileFieldMultiTS> fmtsMultiLev(MEDCoupling::MEDFileFieldMultiTS::New());
1558 fmtsMultiLev->pushBackTimeStep(f1tsMultiLev);
1559 ret->pushField(fmtsMultiLev);
1563 std::string MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh(const std::string& meshName, const MEDCoupling::MEDFileFields *ret)
1565 const char KEY_STR_TO_AVOID_COLLIDE[]="MESH@";
1567 throw INTERP_KERNEL::Exception("MEDFileFieldRepresentationTree::BuildAUniqueArrayNameForMesh : internal error ! NULL ret !");
1568 std::vector<std::string> fieldNamesAlreadyExisting(ret->getFieldsNames());
1569 if(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),meshName)==fieldNamesAlreadyExisting.end())
1571 std::string tmpName(KEY_STR_TO_AVOID_COLLIDE); tmpName+=meshName;
1572 while(std::find(fieldNamesAlreadyExisting.begin(),fieldNamesAlreadyExisting.end(),tmpName)!=fieldNamesAlreadyExisting.end())
1573 tmpName=std::string(KEY_STR_TO_AVOID_COLLIDE)+tmpName;
1577 MEDCoupling::MEDFileFields *MEDFileFieldRepresentationTree::BuildFieldFromMeshes(const MEDCoupling::MEDFileMeshes *ms)
1579 MEDCoupling::MCAuto<MEDCoupling::MEDFileFields> ret(MEDCoupling::MEDFileFields::New());
1580 AppendFieldFromMeshes(ms,ret);
1584 std::vector<std::string> MEDFileFieldRepresentationTree::SplitFieldNameIntoParts(const std::string& fullFieldName, char sep)
1586 std::vector<std::string> ret;
1588 while(pos!=std::string::npos)
1590 std::size_t curPos(fullFieldName.find_first_of(sep,pos));
1591 std::string elt(fullFieldName.substr(pos,curPos!=std::string::npos?curPos-pos:std::string::npos));
1593 pos=fullFieldName.find_first_not_of(sep,curPos);
1599 * Here the non regression tests.
1600 * const char inp0[]="";
1601 * const char exp0[]="";
1602 * const char inp1[]="field";
1603 * const char exp1[]="field";
1604 * const char inp2[]="_________";
1605 * const char exp2[]="_________";
1606 * const char inp3[]="field_p";
1607 * const char exp3[]="field_p";
1608 * const char inp4[]="field__p";
1609 * const char exp4[]="field_p";
1610 * const char inp5[]="field_p__";
1611 * const char exp5[]="field_p";
1612 * const char inp6[]="field_p_";
1613 * const char exp6[]="field_p";
1614 * const char inp7[]="field_____EDFGEG//sdkjf_____PP_______________";
1615 * const char exp7[]="field_EDFGEG//sdkjf_PP";
1616 * const char inp8[]="field_____EDFGEG//sdkjf_____PP";
1617 * const char exp8[]="field_EDFGEG//sdkjf_PP";
1618 * const char inp9[]="_field_____EDFGEG//sdkjf_____PP_______________";
1619 * const char exp9[]="field_EDFGEG//sdkjf_PP";
1620 * const char inp10[]="___field_____EDFGEG//sdkjf_____PP_______________";
1621 * const char exp10[]="field_EDFGEG//sdkjf_PP";
1623 std::string MEDFileFieldRepresentationTree::PostProcessFieldName(const std::string& fullFieldName)
1625 const char SEP('_');
1626 std::vector<std::string> v(SplitFieldNameIntoParts(fullFieldName,SEP));
1628 return fullFieldName;//should never happen
1632 return fullFieldName;
1636 std::string ret(v[0]);
1637 for(std::size_t i=1;i<v.size();i++)
1642 { ret+=SEP; ret+=v[i]; }
1648 return fullFieldName;
1654 TimeKeeper::TimeKeeper(int policy):_policy(policy)
1658 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1663 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1665 return getTimeStepsRegardingPolicy0(tsPairs,ts);
1667 throw INTERP_KERNEL::Exception("TimeKeeper::getTimeStepsRegardingPolicy : only policy 0 and 1 supported presently !");
1673 * if all of ts are in -1e299,1e299 and different each other pairs are ignored ts taken directly.
1674 * if all of ts are in -1e299,1e299 but some are not different each other ts are ignored pairs used
1675 * if some of ts are out of -1e299,1e299 ts are ignored pairs used
1677 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy0(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1679 std::size_t sz(ts.size());
1680 bool isInHumanRange(true);
1682 for(std::size_t i=0;i<sz;i++)
1685 if(ts[i]<=-1e299 || ts[i]>=1e299)
1686 isInHumanRange=false;
1689 return processedUsingPairOfIds(tsPairs);
1691 return processedUsingPairOfIds(tsPairs);
1692 _postprocessed_time=ts;
1693 return getPostProcessedTime();
1698 * idem than 0, except that ts is preaccumulated before invoking policy 0.
1700 std::vector<double> TimeKeeper::getTimeStepsRegardingPolicy1(const std::vector< std::pair<int,int> >& tsPairs, const std::vector<double>& ts) const
1702 std::size_t sz(ts.size());
1703 std::vector<double> ts2(sz);
1705 for(std::size_t i=0;i<sz;i++)
1710 return getTimeStepsRegardingPolicy0(tsPairs,ts2);
1713 int TimeKeeper::getTimeStepIdFrom(double timeReq) const
1715 std::size_t pos(std::distance(_postprocessed_time.begin(),std::find(_postprocessed_time.begin(),_postprocessed_time.end(),timeReq)));
1719 void TimeKeeper::printSelf(std::ostream& oss) const
1721 std::size_t sz(_activated_ts.size());
1722 for(std::size_t i=0;i<sz;i++)
1724 oss << "(" << i << "," << _activated_ts[i].first << "), ";
1728 std::vector<bool> TimeKeeper::getTheVectOfBool() const
1730 std::size_t sz(_activated_ts.size());
1731 std::vector<bool> ret(sz);
1732 for(std::size_t i=0;i<sz;i++)
1734 ret[i]=_activated_ts[i].first;
1739 std::vector<double> TimeKeeper::processedUsingPairOfIds(const std::vector< std::pair<int,int> >& tsPairs) const
1741 std::size_t sz(tsPairs.size());
1742 std::set<int> s0,s1;
1743 for(std::size_t i=0;i<sz;i++)
1744 { s0.insert(tsPairs[i].first); s1.insert(tsPairs[i].second); }
1747 _postprocessed_time.resize(sz);
1748 for(std::size_t i=0;i<sz;i++)
1749 _postprocessed_time[i]=(double)tsPairs[i].first;
1750 return getPostProcessedTime();
1754 _postprocessed_time.resize(sz);
1755 for(std::size_t i=0;i<sz;i++)
1756 _postprocessed_time[i]=(double)tsPairs[i].second;
1757 return getPostProcessedTime();
1759 //TimeKeeper::processedUsingPairOfIds : you are not a lucky guy ! All your time steps info in MEDFile are not discriminant taken one by one !
1760 _postprocessed_time.resize(sz);
1761 for(std::size_t i=0;i<sz;i++)
1762 _postprocessed_time[i]=(double)i;
1763 return getPostProcessedTime();
1766 void TimeKeeper::setMaxNumberOfTimeSteps(int maxNumberOfTS)
1768 _activated_ts.resize(maxNumberOfTS);
1769 for(int i=0;i<maxNumberOfTS;i++)
1771 std::ostringstream oss; oss << "000" << i;
1772 _activated_ts[i]=std::pair<bool,std::string>(true,oss.str());