1 // Copyright (C) 2017 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 (EDF R&D)
21 #include "VTKToMEDMem.hxx"
23 #include "vtkAdjacentVertexIterator.h"
24 #include "vtkIntArray.h"
25 #include "vtkLongArray.h"
26 #include "vtkCellData.h"
27 #include "vtkPointData.h"
28 #include "vtkFloatArray.h"
29 #include "vtkCellArray.h"
31 #include "vtkStreamingDemandDrivenPipeline.h"
32 #include "vtkInformationDataObjectMetaDataKey.h"
33 #include "vtkUnstructuredGrid.h"
34 #include "vtkMultiBlockDataSet.h"
35 #include "vtkRectilinearGrid.h"
36 #include "vtkInformationStringKey.h"
37 #include "vtkAlgorithmOutput.h"
38 #include "vtkObjectFactory.h"
39 #include "vtkMutableDirectedGraph.h"
40 #include "vtkMultiBlockDataSet.h"
41 #include "vtkPolyData.h"
42 #include "vtkDataSet.h"
43 #include "vtkInformationVector.h"
44 #include "vtkInformation.h"
45 #include "vtkDataArraySelection.h"
46 #include "vtkTimeStamp.h"
47 #include "vtkInEdgeIterator.h"
48 #include "vtkInformationDataObjectKey.h"
49 #include "vtkExecutive.h"
50 #include "vtkVariantArray.h"
51 #include "vtkStringArray.h"
52 #include "vtkDoubleArray.h"
53 #include "vtkCharArray.h"
54 #include "vtkUnsignedCharArray.h"
55 #include "vtkDataSetAttributes.h"
56 #include "vtkDemandDrivenPipeline.h"
57 #include "vtkDataObjectTreeIterator.h"
58 #include "vtkWarpScalar.h"
65 using VTKToMEDMem::Grp;
66 using VTKToMEDMem::Fam;
68 using MEDCoupling::MEDFileData;
69 using MEDCoupling::MEDFileMesh;
70 using MEDCoupling::MEDFileCMesh;
71 using MEDCoupling::MEDFileUMesh;
72 using MEDCoupling::MEDFileFields;
73 using MEDCoupling::MEDFileMeshes;
75 using MEDCoupling::MEDFileIntField1TS;
76 using MEDCoupling::MEDFileField1TS;
77 using MEDCoupling::MEDFileIntFieldMultiTS;
78 using MEDCoupling::MEDFileFieldMultiTS;
79 using MEDCoupling::MEDFileAnyTypeFieldMultiTS;
80 using MEDCoupling::DataArray;
81 using MEDCoupling::DataArrayInt;
82 using MEDCoupling::DataArrayFloat;
83 using MEDCoupling::DataArrayDouble;
84 using MEDCoupling::MEDCouplingMesh;
85 using MEDCoupling::MEDCouplingUMesh;
86 using MEDCoupling::MEDCouplingCMesh;
87 using MEDCoupling::MEDCouplingFieldDouble;
88 using MEDCoupling::MEDCouplingFieldFloat;
89 using MEDCoupling::MEDCouplingFieldInt;
90 using MEDCoupling::MCAuto;
91 using MEDCoupling::Traits;
92 using MEDCoupling::MLFieldTraits;
96 Fam::Fam(const std::string& name)
98 static const char ZE_SEP[]="@@][@@";
99 std::size_t pos(name.find(ZE_SEP));
100 std::string name0(name.substr(0,pos)),name1(name.substr(pos+strlen(ZE_SEP)));
101 std::istringstream iss(name1);
108 #include "VTKMEDTraits.hxx"
110 std::map<int,int> ComputeMapOfType()
112 std::map<int,int> ret;
113 int nbOfTypesInMC(sizeof(MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER)/sizeof(int));
114 for(int i=0;i<nbOfTypesInMC;i++)
116 int vtkId(MEDCouplingUMesh::MEDCOUPLING2VTKTYPETRADUCER[i]);
123 std::string GetMeshNameWithContext(const std::vector<int>& context)
125 static const char DFT_MESH_NAME[]="Mesh";
127 return DFT_MESH_NAME;
128 std::ostringstream oss; oss << DFT_MESH_NAME;
129 for(std::vector<int>::const_iterator it=context.begin();it!=context.end();it++)
134 DataArrayInt *ConvertVTKArrayToMCArrayInt(vtkDataArray *data)
137 throw MZCException("ConvertVTKArrayToMCArrayInt : internal error !");
138 int nbTuples(data->GetNumberOfTuples()),nbComp(data->GetNumberOfComponents());
139 std::size_t nbElts(nbTuples*nbComp);
140 MCAuto<DataArrayInt> ret(DataArrayInt::New());
141 ret->alloc(nbTuples,nbComp);
142 for(int i=0;i<nbComp;i++)
144 const char *comp(data->GetComponentName(i));
146 ret->setInfoOnComponent(i,comp);
148 int *ptOut(ret->getPointer());
149 vtkIntArray *d0(vtkIntArray::SafeDownCast(data));
152 const int *pt(d0->GetPointer(0));
153 std::copy(pt,pt+nbElts,ptOut);
156 vtkLongArray *d1(vtkLongArray::SafeDownCast(data));
159 const long *pt(d1->GetPointer(0));
160 std::copy(pt,pt+nbElts,ptOut);
163 vtkUnsignedCharArray *d2(vtkUnsignedCharArray::SafeDownCast(data));
166 const unsigned char *pt(d2->GetPointer(0));
167 std::copy(pt,pt+nbElts,ptOut);
170 std::ostringstream oss;
171 oss << "ConvertVTKArrayToMCArrayInt : unrecognized array \"" << typeid(*data).name() << "\" type !";
172 throw MZCException(oss.str());
176 typename Traits<T>::ArrayType *ConvertVTKArrayToMCArrayDouble(vtkDataArray *data)
179 throw MZCException("ConvertVTKArrayToMCArrayDouble : internal error !");
180 int nbTuples(data->GetNumberOfTuples()),nbComp(data->GetNumberOfComponents());
181 std::size_t nbElts(nbTuples*nbComp);
182 MCAuto< typename Traits<T>::ArrayType > ret(Traits<T>::ArrayType::New());
183 ret->alloc(nbTuples,nbComp);
184 for(int i=0;i<nbComp;i++)
186 const char *comp(data->GetComponentName(i));
188 ret->setInfoOnComponent(i,comp);
191 if(nbComp>1 && nbComp<=3)
194 tmp[0]='X'+i; tmp[1]='\0';
195 ret->setInfoOnComponent(i,tmp);
199 T *ptOut(ret->getPointer());
200 typename MEDFileVTKTraits<T>::VtkType *d0(MEDFileVTKTraits<T>::VtkType::SafeDownCast(data));
203 const T *pt(d0->GetPointer(0));
204 for(std::size_t i=0;i<nbElts;i++)
208 std::ostringstream oss;
209 oss << "ConvertVTKArrayToMCArrayDouble : unrecognized array \"" << data->GetClassName() << "\" type !";
210 throw MZCException(oss.str());
213 DataArrayDouble *ConvertVTKArrayToMCArrayDoubleForced(vtkDataArray *data)
216 throw MZCException("ConvertVTKArrayToMCArrayDoubleForced : internal error 0 !");
217 vtkFloatArray *d0(vtkFloatArray::SafeDownCast(data));
220 MCAuto<DataArrayFloat> ret(ConvertVTKArrayToMCArrayDouble<float>(data));
221 MCAuto<DataArrayDouble> ret2(ret->convertToDblArr());
224 vtkDoubleArray *d1(vtkDoubleArray::SafeDownCast(data));
226 return ConvertVTKArrayToMCArrayDouble<double>(data);
227 throw MZCException("ConvertVTKArrayToMCArrayDoubleForced : unrecognized type of data for double !");
230 DataArray *ConvertVTKArrayToMCArray(vtkDataArray *data)
233 throw MZCException("ConvertVTKArrayToMCArray : internal error !");
234 vtkFloatArray *d0(vtkFloatArray::SafeDownCast(data));
236 return ConvertVTKArrayToMCArrayDouble<float>(data);
237 vtkDoubleArray *d1(vtkDoubleArray::SafeDownCast(data));
239 return ConvertVTKArrayToMCArrayDouble<double>(data);
240 vtkIntArray *d2(vtkIntArray::SafeDownCast(data));
241 vtkLongArray *d3(vtkLongArray::SafeDownCast(data));
242 vtkUnsignedCharArray *d4(vtkUnsignedCharArray::SafeDownCast(data));
244 return ConvertVTKArrayToMCArrayInt(data);
245 std::ostringstream oss;
246 oss << "ConvertVTKArrayToMCArray : unrecognized array \"" << typeid(*data).name() << "\" type !";
247 throw MZCException(oss.str());
250 MEDCouplingUMesh *BuildMeshFromCellArray(vtkCellArray *ca, DataArrayDouble *coords, int meshDim, INTERP_KERNEL::NormalizedCellType type)
252 MCAuto<MEDCouplingUMesh> subMesh(MEDCouplingUMesh::New("",meshDim));
253 subMesh->setCoords(coords); subMesh->allocateCells();
254 int nbCells(ca->GetNumberOfCells());
257 vtkIdType nbEntries(ca->GetNumberOfConnectivityEntries());
258 const vtkIdType *conn(ca->GetPointer());
259 for(int i=0;i<nbCells;i++)
262 std::vector<int> conn2(sz);
263 for(int jj=0;jj<sz;jj++)
265 subMesh->insertNextCell(type,sz,&conn2[0]);
268 return subMesh.retn();
271 MEDCouplingUMesh *BuildMeshFromCellArrayTriangleStrip(vtkCellArray *ca, DataArrayDouble *coords, MCAuto<DataArrayInt>& ids)
273 MCAuto<MEDCouplingUMesh> subMesh(MEDCouplingUMesh::New("",2));
274 subMesh->setCoords(coords); subMesh->allocateCells();
275 int nbCells(ca->GetNumberOfCells());
278 vtkIdType nbEntries(ca->GetNumberOfConnectivityEntries());
279 const vtkIdType *conn(ca->GetPointer());
280 ids=DataArrayInt::New() ; ids->alloc(0,1);
281 for(int i=0;i<nbCells;i++)
287 for(int j=0;j<nbTri;j++,conn++)
289 int conn2[3]; conn2[0]=conn[0] ; conn2[1]=conn[1] ; conn2[2]=conn[2];
290 subMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn2);
291 ids->pushBackSilent(i);
296 std::ostringstream oss; oss << "BuildMeshFromCellArrayTriangleStrip : on cell #" << i << " the triangle stip looks bab !";
297 throw MZCException(oss.str());
301 return subMesh.retn();
307 MicroField(const MCAuto<MEDCouplingUMesh>& m, const std::vector<MCAuto<DataArray> >& cellFs):_m(m),_cellFs(cellFs) { }
308 MicroField(const std::vector< MicroField >& vs);
309 void setNodeFields(const std::vector<MCAuto<DataArray> >& nf) { _nodeFs=nf; }
310 MCAuto<MEDCouplingUMesh> getMesh() const { return _m; }
311 std::vector<MCAuto<DataArray> > getCellFields() const { return _cellFs; }
313 MCAuto<MEDCouplingUMesh> _m;
314 std::vector<MCAuto<DataArray> > _cellFs;
315 std::vector<MCAuto<DataArray> > _nodeFs;
318 MicroField::MicroField(const std::vector< MicroField >& vs)
320 std::size_t sz(vs.size());
321 std::vector<const MEDCouplingUMesh *> vs2(sz);
322 std::vector< std::vector< MCAuto<DataArray> > > arrs2(sz);
324 for(std::size_t ii=0;ii<sz;ii++)
326 vs2[ii]=vs[ii].getMesh();
327 arrs2[ii]=vs[ii].getCellFields();
329 nbElts=arrs2[ii].size();
331 if(arrs2[ii].size()!=nbElts)
332 throw MZCException("MicroField cstor : internal error !");
334 _cellFs.resize(nbElts);
335 for(int ii=0;ii<nbElts;ii++)
337 std::vector<const DataArray *> arrsTmp(sz);
338 for(std::size_t jj=0;jj<sz;jj++)
340 arrsTmp[jj]=arrs2[jj][ii];
342 _cellFs[ii]=DataArray::Aggregate(arrsTmp);
344 _m=MEDCouplingUMesh::MergeUMeshesOnSameCoords(vs2);
348 void AppendToFields(MEDCoupling::TypeOfField tf, MEDCouplingMesh *mesh, const DataArrayInt *n2oPtr, typename MEDFileVTKTraits<T>::MCType *dadPtr, MEDFileFields *fs, double timeStep, int tsId)
350 std::string fieldName(dadPtr->getName());
351 MCAuto< typename Traits<T>::FieldType > f(Traits<T>::FieldType::New(tf));
352 f->setTime(timeStep,tsId,0);
354 std::string fieldNameForChuckNorris(MEDCoupling::MEDFileAnyTypeField1TSWithoutSDA::FieldNameToMEDFileConvention(fieldName));
355 f->setName(fieldNameForChuckNorris);
361 MCAuto< typename Traits<T>::ArrayType > dad2(dadPtr->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
365 MCAuto< typename MLFieldTraits<T>::FMTSType > fmts(MLFieldTraits<T>::FMTSType::New());
366 MCAuto< typename MLFieldTraits<T>::F1TSType > f1ts(MLFieldTraits<T>::F1TSType::New());
367 f1ts->setFieldNoProfileSBT(f);
368 fmts->pushBackTimeStep(f1ts);
372 void AppendMCFieldFrom(MEDCoupling::TypeOfField tf, MEDCouplingMesh *mesh, MEDFileData *mfd, MCAuto<DataArray> da, const DataArrayInt *n2oPtr, double timeStep, int tsId)
374 static const char FAMFIELD_FOR_CELLS[]="FamilyIdCell";
375 static const char FAMFIELD_FOR_NODES[]="FamilyIdNode";
376 if(!da || !mesh || !mfd)
377 throw MZCException("AppendMCFieldFrom : internal error !");
378 MEDFileFields *fs(mfd->getFields());
379 MEDFileMeshes *ms(mfd->getMeshes());
381 throw MZCException("AppendMCFieldFrom : internal error 2 !");
382 MCAuto<DataArrayDouble> dad(MEDCoupling::DynamicCast<DataArray,DataArrayDouble>(da));
385 AppendToFields<double>(tf,mesh,n2oPtr,dad,fs,timeStep,tsId);
388 MCAuto<DataArrayFloat> daf(MEDCoupling::DynamicCast<DataArray,DataArrayFloat>(da));
391 AppendToFields<float>(tf,mesh,n2oPtr,daf,fs,timeStep,tsId);
394 MCAuto<DataArrayInt> dai(MEDCoupling::DynamicCast<DataArray,DataArrayInt>(da));
397 std::string fieldName(dai->getName());
398 if((fieldName!=FAMFIELD_FOR_CELLS || tf!=MEDCoupling::ON_CELLS) && (fieldName!=FAMFIELD_FOR_NODES || tf!=MEDCoupling::ON_NODES))
400 AppendToFields<int>(tf,mesh,n2oPtr,dai,fs,timeStep,tsId);
403 else if(fieldName==FAMFIELD_FOR_CELLS && tf==MEDCoupling::ON_CELLS)
405 MEDFileMesh *mm(ms->getMeshWithName(mesh->getName()));
407 throw MZCException("AppendMCFieldFrom : internal error 3 !");
409 mm->setFamilyFieldArr(mesh->getMeshDimension()-mm->getMeshDimension(),dai);
412 MCAuto<DataArrayInt> dai2(dai->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
413 mm->setFamilyFieldArr(mesh->getMeshDimension()-mm->getMeshDimension(),dai2);
416 else if(fieldName==FAMFIELD_FOR_NODES || tf==MEDCoupling::ON_NODES)
418 MEDFileMesh *mm(ms->getMeshWithName(mesh->getName()));
420 throw MZCException("AppendMCFieldFrom : internal error 4 !");
422 mm->setFamilyFieldArr(1,dai);
425 MCAuto<DataArrayInt> dai2(dai->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
426 mm->setFamilyFieldArr(1,dai2);
432 void PutAtLevelDealOrder(MEDFileData *mfd, int meshDimRel, const MicroField& mf, double timeStep, int tsId)
435 throw MZCException("PutAtLevelDealOrder : internal error !");
436 MEDFileMesh *mm(mfd->getMeshes()->getMeshAtPos(0));
437 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
439 throw MZCException("PutAtLevelDealOrder : internal error 2 !");
440 MCAuto<MEDCouplingUMesh> mesh(mf.getMesh());
441 mesh->setName(mfd->getMeshes()->getMeshAtPos(0)->getName());
442 MCAuto<DataArrayInt> o2n(mesh->sortCellsInMEDFileFrmt());
443 const DataArrayInt *o2nPtr(o2n);
444 MCAuto<DataArrayInt> n2o;
445 mmu->setMeshAtLevel(meshDimRel,mesh);
446 const DataArrayInt *n2oPtr(0);
449 n2o=o2n->invertArrayO2N2N2O(mesh->getNumberOfCells());
451 if(n2oPtr && n2oPtr->isIota(mesh->getNumberOfCells()))
454 mm->setRenumFieldArr(meshDimRel,n2o);
457 std::vector<MCAuto<DataArray> > cells(mf.getCellFields());
458 for(std::vector<MCAuto<DataArray> >::const_iterator it=cells.begin();it!=cells.end();it++)
460 MCAuto<DataArray> da(*it);
461 AppendMCFieldFrom(MEDCoupling::ON_CELLS,mesh,mfd,da,n2oPtr,timeStep,tsId);
465 void AssignSingleGTMeshes(MEDFileData *mfd, const std::vector< MicroField >& ms, double timeStep, int tsId)
468 throw MZCException("AssignSingleGTMeshes : internal error !");
469 MEDFileMesh *mm0(mfd->getMeshes()->getMeshAtPos(0));
470 MEDFileUMesh *mm(dynamic_cast<MEDFileUMesh *>(mm0));
472 throw MZCException("AssignSingleGTMeshes : internal error 2 !");
473 int meshDim(-std::numeric_limits<int>::max());
474 std::map<int, std::vector< MicroField > > ms2;
475 for(std::vector< MicroField >::const_iterator it=ms.begin();it!=ms.end();it++)
477 const MEDCouplingUMesh *elt((*it).getMesh());
480 int myMeshDim(elt->getMeshDimension());
481 meshDim=std::max(meshDim,myMeshDim);
482 ms2[myMeshDim].push_back(*it);
487 for(std::map<int, std::vector< MicroField > >::const_iterator it=ms2.begin();it!=ms2.end();it++)
489 const std::vector< MicroField >& vs((*it).second);
492 PutAtLevelDealOrder(mfd,(*it).first-meshDim,vs[0],timeStep,tsId);
496 MicroField merge(vs);
497 PutAtLevelDealOrder(mfd,(*it).first-meshDim,merge,timeStep,tsId);
502 DataArrayDouble *BuildCoordsFrom(vtkPointSet *ds)
505 throw MZCException("BuildCoordsFrom : internal error !");
506 vtkPoints *pts(ds->GetPoints());
508 throw MZCException("BuildCoordsFrom : internal error 2 !");
509 vtkDataArray *data(pts->GetData());
511 throw MZCException("BuildCoordsFrom : internal error 3 !");
512 return ConvertVTKArrayToMCArrayDoubleForced(data);
515 void AddNodeFields(MEDFileData *mfd, vtkDataSetAttributes *dsa, double timeStep, int tsId)
518 throw MZCException("AddNodeFields : internal error !");
519 MEDFileMesh *mm(mfd->getMeshes()->getMeshAtPos(0));
520 MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
522 throw MZCException("AddNodeFields : internal error 2 !");
523 MCAuto<MEDCouplingUMesh> mesh;
524 if(!mmu->getNonEmptyLevels().empty())
525 mesh=mmu->getMeshAtLevel(0);
528 mesh=MEDCouplingUMesh::Build0DMeshFromCoords(mmu->getCoords());
529 mesh->setName(mmu->getName());
531 int nba(dsa->GetNumberOfArrays());
532 for(int i=0;i<nba;i++)
534 vtkDataArray *arr(dsa->GetArray(i));
535 const char *name(arr->GetName());
538 MCAuto<DataArray> da(ConvertVTKArrayToMCArray(arr));
540 AppendMCFieldFrom(MEDCoupling::ON_NODES,mesh,mfd,da,NULL,timeStep,tsId);
544 std::vector<MCAuto<DataArray> > AddPartFields(const DataArrayInt *part, vtkDataSetAttributes *dsa)
546 std::vector< MCAuto<DataArray> > ret;
549 int nba(dsa->GetNumberOfArrays());
550 for(int i=0;i<nba;i++)
552 vtkDataArray *arr(dsa->GetArray(i));
555 const char *name(arr->GetName());
556 int nbCompo(arr->GetNumberOfComponents());
557 vtkIdType nbTuples(arr->GetNumberOfTuples());
558 MCAuto<DataArray> mcarr(ConvertVTKArrayToMCArray(arr));
560 mcarr=mcarr->selectByTupleId(part->begin(),part->end());
561 mcarr->setName(name);
562 ret.push_back(mcarr);
567 std::vector<MCAuto<DataArray> > AddPartFields2(int bg, int end, vtkDataSetAttributes *dsa)
569 std::vector< MCAuto<DataArray> > ret;
572 int nba(dsa->GetNumberOfArrays());
573 for(int i=0;i<nba;i++)
575 vtkDataArray *arr(dsa->GetArray(i));
578 const char *name(arr->GetName());
579 int nbCompo(arr->GetNumberOfComponents());
580 vtkIdType nbTuples(arr->GetNumberOfTuples());
581 MCAuto<DataArray> mcarr(ConvertVTKArrayToMCArray(arr));
582 mcarr=mcarr->selectByTupleIdSafeSlice(bg,end,1);
583 mcarr->setName(name);
584 ret.push_back(mcarr);
589 void ConvertFromRectilinearGrid(MEDFileData *ret, vtkRectilinearGrid *ds, const std::vector<int>& context, double timeStep, int tsId)
592 throw MZCException("ConvertFromRectilinearGrid : internal error !");
594 MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
595 ret->setMeshes(meshes);
596 MCAuto<MEDFileFields> fields(MEDFileFields::New());
597 ret->setFields(fields);
599 MCAuto<MEDFileCMesh> cmesh(MEDFileCMesh::New());
600 meshes->pushMesh(cmesh);
601 MCAuto<MEDCouplingCMesh> cmeshmc(MEDCouplingCMesh::New());
602 vtkDataArray *cx(ds->GetXCoordinates()),*cy(ds->GetYCoordinates()),*cz(ds->GetZCoordinates());
605 MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cx));
606 cmeshmc->setCoordsAt(0,arr);
610 MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cy));
611 cmeshmc->setCoordsAt(1,arr);
615 MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cz));
616 cmeshmc->setCoordsAt(2,arr);
618 std::string meshName(GetMeshNameWithContext(context));
619 cmeshmc->setName(meshName);
620 cmesh->setMesh(cmeshmc);
621 std::vector<MCAuto<DataArray> > cellFs(AddPartFields(0,ds->GetCellData()));
622 for(std::vector<MCAuto<DataArray> >::const_iterator it=cellFs.begin();it!=cellFs.end();it++)
624 MCAuto<DataArray> da(*it);
625 AppendMCFieldFrom(MEDCoupling::ON_CELLS,cmeshmc,ret,da,NULL,timeStep,tsId);
627 std::vector<MCAuto<DataArray> > nodeFs(AddPartFields(0,ds->GetPointData()));
628 for(std::vector<MCAuto<DataArray> >::const_iterator it=nodeFs.begin();it!=nodeFs.end();it++)
630 MCAuto<DataArray> da(*it);
631 AppendMCFieldFrom(MEDCoupling::ON_NODES,cmeshmc,ret,da,NULL,timeStep,tsId);
635 void ConvertFromPolyData(MEDFileData *ret, vtkPolyData *ds, const std::vector<int>& context, double timeStep, int tsId)
638 throw MZCException("ConvertFromPolyData : internal error !");
640 MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
641 ret->setMeshes(meshes);
642 MCAuto<MEDFileFields> fields(MEDFileFields::New());
643 ret->setFields(fields);
645 MCAuto<MEDFileUMesh> umesh(MEDFileUMesh::New());
646 meshes->pushMesh(umesh);
647 MCAuto<DataArrayDouble> coords(BuildCoordsFrom(ds));
648 umesh->setCoords(coords);
649 umesh->setName(GetMeshNameWithContext(context));
652 std::vector< MicroField > ms;
653 vtkCellArray *cd(ds->GetVerts());
656 MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArray(cd,coords,0,INTERP_KERNEL::NORM_POINT1));
657 if((const MEDCouplingUMesh *)subMesh)
659 std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
660 offset+=subMesh->getNumberOfCells();
661 ms.push_back(MicroField(subMesh,cellFs));
664 vtkCellArray *cc(ds->GetLines());
667 MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArray(cc,coords,1,INTERP_KERNEL::NORM_SEG2));
668 if((const MEDCouplingUMesh *)subMesh)
670 std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
671 offset+=subMesh->getNumberOfCells();
672 ms.push_back(MicroField(subMesh,cellFs));
675 vtkCellArray *cb(ds->GetPolys());
678 MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArray(cb,coords,2,INTERP_KERNEL::NORM_POLYGON));
679 if((const MEDCouplingUMesh *)subMesh)
681 std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
682 offset+=subMesh->getNumberOfCells();
683 ms.push_back(MicroField(subMesh,cellFs));
686 vtkCellArray *ca(ds->GetStrips());
689 MCAuto<DataArrayInt> ids;
690 MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArrayTriangleStrip(ca,coords,ids));
691 if((const MEDCouplingUMesh *)subMesh)
693 std::vector<MCAuto<DataArray> > cellFs(AddPartFields(ids,ds->GetCellData()));
694 offset+=subMesh->getNumberOfCells();
695 ms.push_back(MicroField(subMesh,cellFs));
698 AssignSingleGTMeshes(ret,ms,timeStep,tsId);
699 AddNodeFields(ret,ds->GetPointData(),timeStep,tsId);
702 void ConvertFromUnstructuredGrid(MEDFileData *ret, vtkUnstructuredGrid *ds, const std::vector<int>& context, double timeStep, int tsId)
705 throw MZCException("ConvertFromUnstructuredGrid : internal error !");
707 MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
708 ret->setMeshes(meshes);
709 MCAuto<MEDFileFields> fields(MEDFileFields::New());
710 ret->setFields(fields);
712 MCAuto<MEDFileUMesh> umesh(MEDFileUMesh::New());
713 meshes->pushMesh(umesh);
714 MCAuto<DataArrayDouble> coords(BuildCoordsFrom(ds));
715 umesh->setCoords(coords);
716 umesh->setName(GetMeshNameWithContext(context));
717 vtkIdType nbCells(ds->GetNumberOfCells());
718 vtkCellArray *ca(ds->GetCells());
721 vtkIdType nbEnt(ca->GetNumberOfConnectivityEntries());
722 vtkIdType *caPtr(ca->GetPointer());
723 vtkUnsignedCharArray *ct(ds->GetCellTypesArray());
725 throw MZCException("ConvertFromUnstructuredGrid : internal error");
726 vtkIdTypeArray *cla(ds->GetCellLocationsArray());
727 const vtkIdType *claPtr(cla->GetPointer(0));
729 throw MZCException("ConvertFromUnstructuredGrid : internal error 2");
730 const unsigned char *ctPtr(ct->GetPointer(0));
731 std::map<int,int> m(ComputeMapOfType());
732 MCAuto<DataArrayInt> lev(DataArrayInt::New()) ; lev->alloc(nbCells,1);
733 int *levPtr(lev->getPointer());
734 for(vtkIdType i=0;i<nbCells;i++)
736 std::map<int,int>::iterator it(m.find(ctPtr[i]));
739 const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)(*it).second));
740 levPtr[i]=cm.getDimension();
744 std::ostringstream oss; oss << "ConvertFromUnstructuredGrid : at pos #" << i << " unrecognized VTK cell with type =" << ctPtr[i];
745 throw MZCException(oss.str());
749 MCAuto<DataArrayInt> levs(lev->getDifferentValues());
750 std::vector< MicroField > ms;
751 vtkIdTypeArray *faces(ds->GetFaces()),*faceLoc(ds->GetFaceLocations());
752 for(const int *curLev=levs->begin();curLev!=levs->end();curLev++)
754 MCAuto<MEDCouplingUMesh> m0(MEDCouplingUMesh::New("",*curLev));
755 m0->setCoords(coords); m0->allocateCells();
756 MCAuto<DataArrayInt> cellIdsCurLev(lev->findIdsEqual(*curLev));
757 for(const int *cellId=cellIdsCurLev->begin();cellId!=cellIdsCurLev->end();cellId++)
759 std::map<int,int>::iterator it(m.find(ctPtr[*cellId]));
760 vtkIdType offset(claPtr[*cellId]);
761 vtkIdType sz(caPtr[offset]);
762 INTERP_KERNEL::NormalizedCellType ct((INTERP_KERNEL::NormalizedCellType)(*it).second);
763 if(ct!=INTERP_KERNEL::NORM_POLYHED)
765 std::vector<int> conn2(sz);
766 for(int kk=0;kk<sz;kk++)
767 conn2[kk]=caPtr[offset+1+kk];
768 m0->insertNextCell(ct,sz,&conn2[0]);
772 if(!faces || !faceLoc)
773 throw MZCException("ConvertFromUnstructuredGrid : faces are expected when there are polyhedra !");
774 const vtkIdType *facPtr(faces->GetPointer(0)),*facLocPtr(faceLoc->GetPointer(0));
775 std::vector<int> conn;
776 int off0(facLocPtr[*cellId]);
777 int nbOfFaces(facPtr[off0++]);
778 for(int k=0;k<nbOfFaces;k++)
780 int nbOfNodesInFace(facPtr[off0++]);
781 std::copy(facPtr+off0,facPtr+off0+nbOfNodesInFace,std::back_inserter(conn));
782 off0+=nbOfNodesInFace;
786 m0->insertNextCell(ct,conn.size(),&conn[0]);
789 std::vector<MCAuto<DataArray> > cellFs(AddPartFields(cellIdsCurLev,ds->GetCellData()));
790 ms.push_back(MicroField(m0,cellFs));
792 AssignSingleGTMeshes(ret,ms,timeStep,tsId);
793 AddNodeFields(ret,ds->GetPointData(),timeStep,tsId);
798 void WriteMEDFileFromVTKDataSet(MEDFileData *mfd, vtkDataSet *ds, const std::vector<int>& context, double timeStep, int tsId)
801 throw MZCException("Internal error in WriteMEDFileFromVTKDataSet.");
802 vtkPolyData *ds2(vtkPolyData::SafeDownCast(ds));
803 vtkUnstructuredGrid *ds3(vtkUnstructuredGrid::SafeDownCast(ds));
804 vtkRectilinearGrid *ds4(vtkRectilinearGrid::SafeDownCast(ds));
807 ConvertFromPolyData(mfd,ds2,context,timeStep,tsId);
811 ConvertFromUnstructuredGrid(mfd,ds3,context,timeStep,tsId);
815 ConvertFromRectilinearGrid(mfd,ds4,context,timeStep,tsId);
818 throw MZCException("Unrecognized vtkDataSet ! Sorry ! Try to convert it to UnstructuredGrid to be able to write it !");
821 void WriteMEDFileFromVTKMultiBlock(MEDFileData *mfd, vtkMultiBlockDataSet *ds, const std::vector<int>& context, double timeStep, int tsId)
824 throw MZCException("Internal error in WriteMEDFileFromVTKMultiBlock.");
825 int nbBlocks(ds->GetNumberOfBlocks());
826 if(nbBlocks==1 && context.empty())
828 vtkDataObject *uniqueElt(ds->GetBlock(0));
830 throw MZCException("Unique elt in multiblock is NULL !");
831 vtkDataSet *uniqueEltc(vtkDataSet::SafeDownCast(uniqueElt));
834 WriteMEDFileFromVTKDataSet(mfd,uniqueEltc,context,timeStep,tsId);
838 for(int i=0;i<nbBlocks;i++)
840 vtkDataObject *elt(ds->GetBlock(i));
841 std::vector<int> context2;
842 context2.push_back(i);
845 std::ostringstream oss; oss << "In context ";
846 std::copy(context.begin(),context.end(),std::ostream_iterator<int>(oss," "));
847 oss << " at pos #" << i << " elt is NULL !";
848 throw MZCException(oss.str());
850 vtkDataSet *elt1(vtkDataSet::SafeDownCast(elt));
853 WriteMEDFileFromVTKDataSet(mfd,elt1,context,timeStep,tsId);
856 vtkMultiBlockDataSet *elt2(vtkMultiBlockDataSet::SafeDownCast(elt));
859 WriteMEDFileFromVTKMultiBlock(mfd,elt2,context,timeStep,tsId);
862 std::ostringstream oss; oss << "In context ";
863 std::copy(context.begin(),context.end(),std::ostream_iterator<int>(oss," "));
864 oss << " at pos #" << i << " elt not recognized data type !";
865 throw MZCException(oss.str());
869 void WriteMEDFileFromVTKGDS(MEDFileData *mfd, vtkDataObject *input, double timeStep, int tsId)
872 throw MZCException("WriteMEDFileFromVTKGDS : internal error !");
873 std::vector<int> context;
874 vtkDataSet *input1(vtkDataSet::SafeDownCast(input));
877 WriteMEDFileFromVTKDataSet(mfd,input1,context,timeStep,tsId);
880 vtkMultiBlockDataSet *input2(vtkMultiBlockDataSet::SafeDownCast(input));
883 WriteMEDFileFromVTKMultiBlock(mfd,input2,context,timeStep,tsId);
886 throw MZCException("WriteMEDFileFromVTKGDS : not recognized data type !");
889 void PutFamGrpInfoIfAny(MEDFileData *mfd, const std::string& meshName, const std::vector<Grp>& groups, const std::vector<Fam>& fams)
895 MEDFileMeshes *meshes(mfd->getMeshes());
898 if(meshes->getNumberOfMeshes()!=1)
900 MEDFileMesh *mm(meshes->getMeshAtPos(0));
903 mm->setName(meshName);
904 for(std::vector<Fam>::const_iterator it=fams.begin();it!=fams.end();it++)
905 mm->setFamilyId((*it).getName(),(*it).getID());
906 for(std::vector<Grp>::const_iterator it=groups.begin();it!=groups.end();it++)
907 mm->setFamiliesOnGroup((*it).getName(),(*it).getFamilies());
908 MEDFileFields *fields(mfd->getFields());
911 for(int i=0;i<fields->getNumberOfFields();i++)
913 MEDFileAnyTypeFieldMultiTS *fmts(fields->getFieldAtPos(i));
916 fmts->setMeshName(meshName);