src/Plugins/DevelopedSurface/plugin/DevelopedSurfaceModule/VTKToMEDMem.cxx

   1 // Copyright (C) 2017-2022  CEA/DEN, EDF R&D
   2 //
   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.
   7 //
   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.
  12 //
  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
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (EDF R&D)
  20
  21 #include "VTKToMEDMem.h"
  22
  23 #include "vtkAdjacentVertexIterator.h"
  24 #include "vtkIntArray.h"
  25 #include "vtkLongArray.h"
  26 #include "vtkLongLongArray.h"
  27 #include "vtkCellData.h"
  28 #include "vtkPointData.h"
  29 #include "vtkFloatArray.h"
  30 #include "vtkCellArray.h"
  31
  32 #include "vtkStreamingDemandDrivenPipeline.h"
  33 #include "vtkInformationDataObjectMetaDataKey.h"
  34 #include "vtkUnstructuredGrid.h"
  35 #include "vtkMultiBlockDataSet.h"
  36 #include "vtkRectilinearGrid.h"
  37 #include "vtkInformationStringKey.h"
  38 #include "vtkAlgorithmOutput.h"
  39 #include "vtkObjectFactory.h"
  40 #include "vtkMutableDirectedGraph.h"
  41 #include "vtkMultiBlockDataSet.h"
  42 #include "vtkPolyData.h"
  43 #include "vtkDataSet.h"
  44 #include "vtkInformationVector.h"
  45 #include "vtkInformation.h"
  46 #include "vtkDataArraySelection.h"
  47 #include "vtkTimeStamp.h"
  48 #include "vtkInEdgeIterator.h"
  49 #include "vtkInformationDataObjectKey.h"
  50 #include "vtkExecutive.h"
  51 #include "vtkVariantArray.h"
  52 #include "vtkStringArray.h"
  53 #include "vtkDoubleArray.h"
  54 #include "vtkCharArray.h"
  55 #include "vtkUnsignedCharArray.h"
  56 #include "vtkDataSetAttributes.h"
  57 #include "vtkDemandDrivenPipeline.h"
  58 #include "vtkDataObjectTreeIterator.h"
  59 #include "vtkWarpScalar.h"
  60
  61 #include <map>
  62 #include <deque>
  63 #include <sstream>
  64 #include <cstring>
  65
  66 using VTKToMEDMem::Grp;
  67 using VTKToMEDMem::Fam;
  68
  69 using MEDCoupling::MEDFileData;
  70 using MEDCoupling::MEDFileMesh;
  71 using MEDCoupling::MEDFileCMesh;
  72 using MEDCoupling::MEDFileUMesh;
  73 using MEDCoupling::MEDFileFields;
  74 using MEDCoupling::MEDFileMeshes;
  75
  76 using MEDCoupling::MEDFileInt32Field1TS;
  77 using MEDCoupling::MEDFileInt64Field1TS;
  78 using MEDCoupling::MEDFileField1TS;
  79 using MEDCoupling::MEDFileIntFieldMultiTS;
  80 using MEDCoupling::MEDFileFieldMultiTS;
  81 using MEDCoupling::MEDFileAnyTypeFieldMultiTS;
  82 using MEDCoupling::DataArray;
  83 using MEDCoupling::DataArrayInt32;
  84 using MEDCoupling::DataArrayInt64;
  85 using MEDCoupling::DataArrayFloat;
  86 using MEDCoupling::DataArrayDouble;
  87 using MEDCoupling::MEDCouplingMesh;
  88 using MEDCoupling::MEDCouplingUMesh;
  89 using MEDCoupling::MEDCouplingCMesh;
  90 using MEDCoupling::MEDCouplingFieldDouble;
  91 using MEDCoupling::MEDCouplingFieldFloat;
  92 using MEDCoupling::MEDCouplingFieldInt;
  93 using MEDCoupling::MEDCouplingFieldInt64;
  94 using MEDCoupling::MCAuto;
  95 using MEDCoupling::Traits;
  96 using MEDCoupling::MLFieldTraits;
  97
  98 ///////////////////
  99
 100 Fam::Fam(const std::string& name)
 101 {
 102   static const char ZE_SEP[]="@@][@@";
 103   std::size_t pos(name.find(ZE_SEP));
 104   std::string name0(name.substr(0,pos)),name1(name.substr(pos+strlen(ZE_SEP)));
 105   std::istringstream iss(name1);
 106   iss >> _id;
 107   _name=name0;
 108 }
 109
 110 ///////////////////
 111
 112 #include "VTKMEDTraits.hxx"
 113
 114 std::map<int,int> ComputeMapOfType()
 115 {
 116   std::map<int,int> ret;
 117   int nbOfTypesInMC(sizeof(MEDCOUPLING2VTKTYPETRADUCER)/sizeof( decltype(MEDCOUPLING2VTKTYPETRADUCER[0]) ));
 118   for(int i=0;i<nbOfTypesInMC;i++)
 119     {
 120       auto vtkId(MEDCOUPLING2VTKTYPETRADUCER[i]);
 121       if(vtkId!=MEDCOUPLING2VTKTYPETRADUCER_NONE)
 122         ret[vtkId]=i;
 123     }
 124   return ret;
 125 }
 126
 127 std::string GetMeshNameWithContext(const std::vector<int>& context)
 128 {
 129   static const char DFT_MESH_NAME[]="Mesh";
 130   if(context.empty())
 131     return DFT_MESH_NAME;
 132   std::ostringstream oss; oss << DFT_MESH_NAME;
 133   for(std::vector<int>::const_iterator it=context.begin();it!=context.end();it++)
 134     oss << "_" << *it;
 135   return oss.str();
 136 }
 137
 138 DataArrayIdType *ConvertVTKArrayToMCArrayInt(vtkDataArray *data)
 139 {
 140   if(!data)
 141     throw MZCException("ConvertVTKArrayToMCArrayInt : internal error !");
 142   int nbTuples(data->GetNumberOfTuples()),nbComp(data->GetNumberOfComponents());
 143   std::size_t nbElts(nbTuples*nbComp);
 144   MCAuto<DataArrayIdType> ret(DataArrayIdType::New());
 145   ret->alloc(nbTuples,nbComp);
 146   for(int i=0;i<nbComp;i++)
 147     {
 148       const char *comp(data->GetComponentName(i));
 149       if(comp)
 150         ret->setInfoOnComponent(i,comp);
 151     }
 152   mcIdType *ptOut(ret->getPointer());
 153   vtkIntArray *d0(vtkIntArray::SafeDownCast(data));
 154   if(d0)
 155     {
 156       const int *pt(d0->GetPointer(0));
 157       std::copy(pt,pt+nbElts,ptOut);
 158       return ret.retn();
 159     }
 160   vtkLongArray *d1(vtkLongArray::SafeDownCast(data));
 161   if(d1)
 162     {
 163       const long *pt(d1->GetPointer(0));
 164       std::copy(pt,pt+nbElts,ptOut);
 165       return ret.retn();
 166     }
 167   vtkLongLongArray *d1l(vtkLongLongArray::SafeDownCast(data));
 168   if(d1l)
 169     {
 170       const long long *pt(d1l->GetPointer(0));
 171       std::copy(pt,pt+nbElts,ptOut);
 172       return ret.retn();
 173     }
 174   vtkIdTypeArray *d3(vtkIdTypeArray::SafeDownCast(data));
 175   if(d3)
 176     {
 177       const vtkIdType *pt(d3->GetPointer(0));
 178       std::copy(pt,pt+nbElts,ptOut);
 179       return ret.retn();
 180     }
 181   vtkUnsignedCharArray *d2(vtkUnsignedCharArray::SafeDownCast(data));
 182   if(d2)
 183     {
 184       const unsigned char *pt(d2->GetPointer(0));
 185       std::copy(pt,pt+nbElts,ptOut);
 186       return ret.retn();
 187     }
 188   std::ostringstream oss;
 189   oss << "ConvertVTKArrayToMCArrayInt : unrecognized array \"" << typeid(*data).name() << "\" type !";
 190   throw MZCException(oss.str());
 191 }
 192
 193 template<class T>
 194 typename Traits<T>::ArrayType *ConvertVTKArrayToMCArrayDouble(vtkDataArray *data)
 195 {
 196   if(!data)
 197     throw MZCException("ConvertVTKArrayToMCArrayDouble : internal error !");
 198   int nbTuples(data->GetNumberOfTuples()),nbComp(data->GetNumberOfComponents());
 199   std::size_t nbElts(nbTuples*nbComp);
 200   MCAuto< typename Traits<T>::ArrayType > ret(Traits<T>::ArrayType::New());
 201   ret->alloc(nbTuples,nbComp);
 202   for(int i=0;i<nbComp;i++)
 203     {
 204       const char *comp(data->GetComponentName(i));
 205       if(comp)
 206         ret->setInfoOnComponent(i,comp);
 207       else
 208         {
 209           if(nbComp>1 && nbComp<=3)
 210             {
 211               char tmp[2];
 212               tmp[0]=(char)('X'+i); tmp[1]='\0';
 213               ret->setInfoOnComponent(i,tmp);
 214             }
 215         }
 216     }
 217   T *ptOut(ret->getPointer());
 218   typename MEDFileVTKTraits<T>::VtkType *d0(MEDFileVTKTraits<T>::VtkType::SafeDownCast(data));
 219   if(d0)
 220     {
 221       const T *pt(d0->GetPointer(0));
 222       for(std::size_t i=0;i<nbElts;i++)
 223         ptOut[i]=(T)pt[i];
 224       return ret.retn();
 225     }
 226   std::ostringstream oss;
 227   oss << "ConvertVTKArrayToMCArrayDouble : unrecognized array \"" << data->GetClassName() << "\" type !";
 228   throw MZCException(oss.str());
 229 }
 230
 231 DataArrayDouble *ConvertVTKArrayToMCArrayDoubleForced(vtkDataArray *data)
 232 {
 233   if(!data)
 234     throw MZCException("ConvertVTKArrayToMCArrayDoubleForced : internal error 0 !");
 235   vtkFloatArray *d0(vtkFloatArray::SafeDownCast(data));
 236   if(d0)
 237     {
 238       MCAuto<DataArrayFloat> ret(ConvertVTKArrayToMCArrayDouble<float>(data));
 239       MCAuto<DataArrayDouble> ret2(ret->convertToDblArr());
 240       return ret2.retn();
 241     }
 242   vtkDoubleArray *d1(vtkDoubleArray::SafeDownCast(data));
 243   if(d1)
 244     return ConvertVTKArrayToMCArrayDouble<double>(data);
 245   throw MZCException("ConvertVTKArrayToMCArrayDoubleForced : unrecognized type of data for double !");
 246 }
 247
 248 DataArray *ConvertVTKArrayToMCArray(vtkDataArray *data)
 249 {
 250   if(!data)
 251     throw MZCException("ConvertVTKArrayToMCArray : internal error !");
 252   vtkFloatArray *d0(vtkFloatArray::SafeDownCast(data));
 253   if(d0)
 254     return ConvertVTKArrayToMCArrayDouble<float>(data);
 255   vtkDoubleArray *d1(vtkDoubleArray::SafeDownCast(data));
 256   if(d1)
 257     return ConvertVTKArrayToMCArrayDouble<double>(data);
 258   vtkIntArray *d2(vtkIntArray::SafeDownCast(data));
 259   vtkLongArray *d3(vtkLongArray::SafeDownCast(data));
 260   vtkUnsignedCharArray *d4(vtkUnsignedCharArray::SafeDownCast(data));
 261   vtkIdTypeArray *d5(vtkIdTypeArray::SafeDownCast(data));
 262   vtkLongLongArray *d6(vtkLongLongArray::SafeDownCast(data));
 263   if(d2 || d3 || d4 || d5 || d6)
 264     return ConvertVTKArrayToMCArrayInt(data);
 265   std::ostringstream oss;
 266   oss << "ConvertVTKArrayToMCArray : unrecognized array \"" << typeid(*data).name() << "\" type !";
 267   throw MZCException(oss.str());
 268 }
 269
 270 MEDCouplingUMesh *BuildMeshFromCellArray(vtkCellArray *ca, DataArrayDouble *coords, int meshDim, INTERP_KERNEL::NormalizedCellType type)
 271 {
 272   MCAuto<MEDCouplingUMesh> subMesh(MEDCouplingUMesh::New("",meshDim));
 273   subMesh->setCoords(coords); subMesh->allocateCells();
 274   int nbCells(ca->GetNumberOfCells());
 275   if(nbCells==0)
 276     return 0;
 277   //vtkIdType nbEntries(ca->GetNumberOfConnectivityEntries()); // todo: unused
 278   const vtkIdType *conn(ca->GetData()->GetPointer(0));
 279   for(int i=0;i<nbCells;i++)
 280     {
 281       mcIdType sz(ToIdType(*conn++));
 282       std::vector<mcIdType> conn2(sz);
 283       for(int jj=0;jj<sz;jj++)
 284         conn2[jj]=ToIdType(conn[jj]);
 285       subMesh->insertNextCell(type,sz,&conn2[0]);
 286       conn+=sz;
 287     }
 288   return subMesh.retn();
 289 }
 290
 291 MEDCouplingUMesh *BuildMeshFromCellArrayTriangleStrip(vtkCellArray *ca, DataArrayDouble *coords, MCAuto<DataArrayIdType>& ids)
 292 {
 293   MCAuto<MEDCouplingUMesh> subMesh(MEDCouplingUMesh::New("",2));
 294   subMesh->setCoords(coords); subMesh->allocateCells();
 295   int nbCells(ca->GetNumberOfCells());
 296   if(nbCells==0)
 297     return 0;
 298   //vtkIdType nbEntries(ca->GetNumberOfConnectivityEntries()); // todo: unused
 299   const vtkIdType *conn(ca->GetData()->GetPointer(0));
 300   ids=DataArrayIdType::New() ; ids->alloc(0,1);
 301   for(int i=0;i<nbCells;i++)
 302     {
 303       int sz(*conn++);
 304       int nbTri(sz-2);
 305       if(nbTri>0)
 306         {
 307           for(int j=0;j<nbTri;j++,conn++)
 308             {
 309               mcIdType conn2[3]; conn2[0]=conn[0] ; conn2[1]=conn[1] ; conn2[2]=conn[2];
 310               subMesh->insertNextCell(INTERP_KERNEL::NORM_TRI3,3,conn2);
 311               ids->pushBackSilent(i);
 312             }
 313         }
 314       else
 315         {
 316           std::ostringstream oss; oss << "BuildMeshFromCellArrayTriangleStrip : on cell #" << i << " the triangle stip looks bab !";
 317           throw MZCException(oss.str());
 318         }
 319       conn+=sz;
 320     }
 321   return subMesh.retn();
 322 }
 323
 324 class MicroField
 325 {
 326 public:
 327   MicroField(const MCAuto<MEDCouplingUMesh>& m, const std::vector<MCAuto<DataArray> >& cellFs):_m(m),_cellFs(cellFs) { }
 328   MicroField(const std::vector< MicroField >& vs);
 329   void setNodeFields(const std::vector<MCAuto<DataArray> >& nf) { _nodeFs=nf; }
 330   MCAuto<MEDCouplingUMesh> getMesh() const { return _m; }
 331   std::vector<MCAuto<DataArray> > getCellFields() const { return _cellFs; }
 332 private:
 333   MCAuto<MEDCouplingUMesh> _m;
 334   std::vector<MCAuto<DataArray> > _cellFs;
 335   std::vector<MCAuto<DataArray> > _nodeFs;
 336 };
 337
 338 MicroField::MicroField(const std::vector< MicroField >& vs)
 339 {
 340   std::size_t sz(vs.size());
 341   std::vector<const MEDCouplingUMesh *> vs2(sz);
 342   std::vector< std::vector< MCAuto<DataArray> > > arrs2(sz);
 343   int nbElts(-1);
 344   for(std::size_t ii=0;ii<sz;ii++)
 345     {
 346       vs2[ii]=vs[ii].getMesh();
 347       arrs2[ii]=vs[ii].getCellFields();
 348       if(nbElts<0)
 349         nbElts=(int)arrs2[ii].size();
 350       else
 351         if((int)arrs2[ii].size()!=nbElts)
 352           throw MZCException("MicroField cstor : internal error !");
 353     }
 354   _cellFs.resize(nbElts);
 355   for(int ii=0;ii<nbElts;ii++)
 356     {
 357       std::vector<const DataArray *> arrsTmp(sz);
 358       for(std::size_t jj=0;jj<sz;jj++)
 359         {
 360           arrsTmp[jj]=arrs2[jj][ii];
 361         }
 362       _cellFs[ii]=DataArray::Aggregate(arrsTmp);
 363     }
 364   _m=MEDCouplingUMesh::MergeUMeshesOnSameCoords(vs2);
 365 }
 366
 367 template<class T>
 368 void AppendToFields(MEDCoupling::TypeOfField tf, MEDCouplingMesh *mesh, const DataArrayIdType *n2oPtr, typename MEDFileVTKTraits<T>::MCType *dadPtr, MEDFileFields *fs, double timeStep, int tsId)
 369 {
 370   std::string fieldName(dadPtr->getName());
 371   MCAuto< typename Traits<T>::FieldType > f(Traits<T>::FieldType::New(tf));
 372   f->setTime(timeStep,tsId,0);
 373   {
 374     std::string fieldNameForChuckNorris(MEDCoupling::MEDFileAnyTypeField1TSWithoutSDA::FieldNameToMEDFileConvention(fieldName));
 375     f->setName(fieldNameForChuckNorris);
 376   }
 377   if(!n2oPtr)
 378     f->setArray(dadPtr);
 379   else
 380     {
 381       MCAuto< typename Traits<T>::ArrayType > dad2(dadPtr->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
 382       f->setArray(dad2);
 383     }
 384   f->setMesh(mesh);
 385   MCAuto< typename MLFieldTraits<T>::FMTSType > fmts(MLFieldTraits<T>::FMTSType::New());
 386   MCAuto< typename MLFieldTraits<T>::F1TSType > f1ts(MLFieldTraits<T>::F1TSType::New());
 387   f1ts->setFieldNoProfileSBT(f);
 388   fmts->pushBackTimeStep(f1ts);
 389   fs->pushField(fmts);
 390 }
 391
 392 void AppendMCFieldFrom(MEDCoupling::TypeOfField tf, MEDCouplingMesh *mesh, MEDFileData *mfd, MCAuto<DataArray> da, const DataArrayIdType *n2oPtr, double timeStep, int tsId)
 393 {
 394   static const char FAMFIELD_FOR_CELLS[]="FamilyIdCell";
 395   static const char FAMFIELD_FOR_NODES[]="FamilyIdNode";
 396   if(!da || !mesh || !mfd)
 397     throw MZCException("AppendMCFieldFrom : internal error !");
 398   MEDFileFields *fs(mfd->getFields());
 399   MEDFileMeshes *ms(mfd->getMeshes());
 400   if(!fs || !ms)
 401     throw MZCException("AppendMCFieldFrom : internal error 2 !");
 402   MCAuto<DataArrayDouble> dad(MEDCoupling::DynamicCast<DataArray,DataArrayDouble>(da));
 403   if(dad.isNotNull())
 404     {
 405       AppendToFields<double>(tf,mesh,n2oPtr,dad,fs,timeStep,tsId);
 406       return ;
 407     }
 408   MCAuto<DataArrayFloat> daf(MEDCoupling::DynamicCast<DataArray,DataArrayFloat>(da));
 409   if(daf.isNotNull())
 410     {
 411       AppendToFields<float>(tf,mesh,n2oPtr,daf,fs,timeStep,tsId);
 412       return ;
 413     }
 414   MCAuto<DataArrayInt> dai(MEDCoupling::DynamicCast<DataArray,DataArrayInt>(da));
 415   MCAuto<DataArrayIdType> daId(MEDCoupling::DynamicCast<DataArray,DataArrayIdType>(da));
 416   if(dai.isNotNull() || daId.isNotNull())
 417     {
 418       std::string fieldName(da->getName());
 419       if((fieldName!=FAMFIELD_FOR_CELLS || tf!=MEDCoupling::ON_CELLS) && (fieldName!=FAMFIELD_FOR_NODES || tf!=MEDCoupling::ON_NODES))
 420         {
 421           if(!dai)
 422             AppendToFields<mcIdType>(tf,mesh,n2oPtr,daId,fs,timeStep,tsId);
 423           else
 424             AppendToFields<int>(tf,mesh,n2oPtr,dai,fs,timeStep,tsId);
 425           return ;
 426         }
 427       else if(fieldName==FAMFIELD_FOR_CELLS && tf==MEDCoupling::ON_CELLS)
 428         {
 429           MEDFileMesh *mm(ms->getMeshWithName(mesh->getName()));
 430           if(!mm)
 431             throw MZCException("AppendMCFieldFrom : internal error 3 !");
 432           if(!daId)
 433             throw MZCException("AppendMCFieldFrom : internal error 3 (not mcIdType) !");
 434           if(!n2oPtr)
 435             mm->setFamilyFieldArr(mesh->getMeshDimension()-mm->getMeshDimension(),daId);
 436           else
 437             {
 438               MCAuto<DataArrayIdType> dai2(daId->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
 439               mm->setFamilyFieldArr(mesh->getMeshDimension()-mm->getMeshDimension(),dai2);
 440             }
 441         }
 442       else if(fieldName==FAMFIELD_FOR_NODES || tf==MEDCoupling::ON_NODES)
 443         {
 444           MEDFileMesh *mm(ms->getMeshWithName(mesh->getName()));
 445           if(!mm)
 446             throw MZCException("AppendMCFieldFrom : internal error 4 !");
 447           if(!daId)
 448             throw MZCException("AppendMCFieldFrom : internal error 4 (not mcIdType) !");
 449           if(!n2oPtr)
 450             mm->setFamilyFieldArr(1,daId);
 451           else
 452             {
 453               MCAuto<DataArrayIdType> dai2(daId->selectByTupleId(n2oPtr->begin(),n2oPtr->end()));
 454               mm->setFamilyFieldArr(1,dai2);
 455             }
 456         }
 457     }
 458 }
 459
 460 void PutAtLevelDealOrder(MEDFileData *mfd, int meshDimRel, const MicroField& mf, double timeStep, int tsId)
 461 {
 462   if(!mfd)
 463     throw MZCException("PutAtLevelDealOrder : internal error !");
 464   MEDFileMesh *mm(mfd->getMeshes()->getMeshAtPos(0));
 465   MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
 466   if(!mmu)
 467     throw MZCException("PutAtLevelDealOrder : internal error 2 !");
 468   MCAuto<MEDCouplingUMesh> mesh(mf.getMesh());
 469   mesh->setName(mfd->getMeshes()->getMeshAtPos(0)->getName());
 470   MCAuto<DataArrayIdType> o2n(mesh->sortCellsInMEDFileFrmt());
 471   //const DataArrayIdType *o2nPtr(o2n); // todo: unused
 472   MCAuto<DataArrayIdType> n2o;
 473   mmu->setMeshAtLevel(meshDimRel,mesh);
 474   const DataArrayIdType *n2oPtr(0);
 475   if(o2n)
 476     {
 477       n2o=o2n->invertArrayO2N2N2O(mesh->getNumberOfCells());
 478       n2oPtr=n2o;
 479       if(n2oPtr && n2oPtr->isIota(mesh->getNumberOfCells()))
 480         n2oPtr=0;
 481       if(n2oPtr)
 482         mm->setRenumFieldArr(meshDimRel,n2o);
 483     }
 484   //
 485   std::vector<MCAuto<DataArray> > cells(mf.getCellFields());
 486   for(std::vector<MCAuto<DataArray> >::const_iterator it=cells.begin();it!=cells.end();it++)
 487     {
 488       MCAuto<DataArray> da(*it);
 489       AppendMCFieldFrom(MEDCoupling::ON_CELLS,mesh,mfd,da,n2oPtr,timeStep,tsId);
 490     }
 491 }
 492
 493 void AssignSingleGTMeshes(MEDFileData *mfd, const std::vector< MicroField >& ms, double timeStep, int tsId)
 494 {
 495   if(!mfd)
 496     throw MZCException("AssignSingleGTMeshes : internal error !");
 497   MEDFileMesh *mm0(mfd->getMeshes()->getMeshAtPos(0));
 498   MEDFileUMesh *mm(dynamic_cast<MEDFileUMesh *>(mm0));
 499   if(!mm)
 500     throw MZCException("AssignSingleGTMeshes : internal error 2 !");
 501   int meshDim(-std::numeric_limits<int>::max());
 502   std::map<int, std::vector< MicroField > > ms2;
 503   for(std::vector< MicroField >::const_iterator it=ms.begin();it!=ms.end();it++)
 504     {
 505       const MEDCouplingUMesh *elt((*it).getMesh());
 506       if(elt)
 507         {
 508           int myMeshDim(elt->getMeshDimension());
 509           meshDim=std::max(meshDim,myMeshDim);
 510           ms2[myMeshDim].push_back(*it);
 511         }
 512     }
 513   if(ms2.empty())
 514     return ;
 515   for(std::map<int, std::vector< MicroField > >::const_iterator it=ms2.begin();it!=ms2.end();it++)
 516     {
 517       const std::vector< MicroField >& vs((*it).second);
 518       if(vs.size()==1)
 519         {
 520           PutAtLevelDealOrder(mfd,(*it).first-meshDim,vs[0],timeStep,tsId);
 521         }
 522       else
 523         {
 524           MicroField merge(vs);
 525           PutAtLevelDealOrder(mfd,(*it).first-meshDim,merge,timeStep,tsId);
 526         }
 527     }
 528 }
 529
 530 DataArrayDouble *BuildCoordsFrom(vtkPointSet *ds)
 531 {
 532   if(!ds)
 533     throw MZCException("BuildCoordsFrom : internal error !");
 534   vtkPoints *pts(ds->GetPoints());
 535   if(!pts)
 536     throw MZCException("BuildCoordsFrom : internal error 2 !");
 537   vtkDataArray *data(pts->GetData());
 538   if(!data)
 539     throw MZCException("BuildCoordsFrom : internal error 3 !");
 540   return ConvertVTKArrayToMCArrayDoubleForced(data);
 541 }
 542
 543 void AddNodeFields(MEDFileData *mfd, vtkDataSetAttributes *dsa, double timeStep, int tsId)
 544 {
 545   if(!mfd || !dsa)
 546     throw MZCException("AddNodeFields : internal error !");
 547   MEDFileMesh *mm(mfd->getMeshes()->getMeshAtPos(0));
 548   MEDFileUMesh *mmu(dynamic_cast<MEDFileUMesh *>(mm));
 549   if(!mmu)
 550     throw MZCException("AddNodeFields : internal error 2 !");
 551   MCAuto<MEDCouplingUMesh> mesh;
 552   if(!mmu->getNonEmptyLevels().empty())
 553     mesh=mmu->getMeshAtLevel(0);
 554   else
 555     {
 556       mesh=MEDCouplingUMesh::Build0DMeshFromCoords(mmu->getCoords());
 557       mesh->setName(mmu->getName());
 558     }
 559   int nba(dsa->GetNumberOfArrays());
 560   for(int i=0;i<nba;i++)
 561     {
 562       vtkDataArray *arr(dsa->GetArray(i));
 563       const char *name(arr->GetName());
 564       if(!arr)
 565         continue;
 566       MCAuto<DataArray> da(ConvertVTKArrayToMCArray(arr));
 567       da->setName(name);
 568       AppendMCFieldFrom(MEDCoupling::ON_NODES,mesh,mfd,da,NULL,timeStep,tsId);
 569     }
 570 }
 571
 572 std::vector<MCAuto<DataArray> > AddPartFields(const DataArrayIdType *part, vtkDataSetAttributes *dsa)
 573 {
 574   std::vector< MCAuto<DataArray> > ret;
 575   if(!dsa)
 576     return ret;
 577   int nba(dsa->GetNumberOfArrays());
 578   for(int i=0;i<nba;i++)
 579     {
 580       vtkDataArray *arr(dsa->GetArray(i));
 581       if(!arr)
 582         continue;
 583       const char *name(arr->GetName());
 584       //int nbCompo(arr->GetNumberOfComponents()); // todo: unused
 585       //vtkIdType nbTuples(arr->GetNumberOfTuples()); // todo: unused
 586       MCAuto<DataArray> mcarr(ConvertVTKArrayToMCArray(arr));
 587       if(part)
 588         mcarr=mcarr->selectByTupleId(part->begin(),part->end());
 589       mcarr->setName(name);
 590       ret.push_back(mcarr);
 591     }
 592   return ret;
 593 }
 594
 595 std::vector<MCAuto<DataArray> > AddPartFields2(int bg, int end, vtkDataSetAttributes *dsa)
 596 {
 597   std::vector< MCAuto<DataArray> > ret;
 598   if(!dsa)
 599     return ret;
 600   int nba(dsa->GetNumberOfArrays());
 601   for(int i=0;i<nba;i++)
 602     {
 603       vtkDataArray *arr(dsa->GetArray(i));
 604       if(!arr)
 605         continue;
 606       const char *name(arr->GetName());
 607       //int nbCompo(arr->GetNumberOfComponents()); // todo: unused
 608       //vtkIdType nbTuples(arr->GetNumberOfTuples()); // todo: unused
 609       MCAuto<DataArray> mcarr(ConvertVTKArrayToMCArray(arr));
 610       mcarr=mcarr->selectByTupleIdSafeSlice(bg,end,1);
 611       mcarr->setName(name);
 612       ret.push_back(mcarr);
 613     }
 614   return ret;
 615 }
 616
 617 void ConvertFromRectilinearGrid(MEDFileData *ret, vtkRectilinearGrid *ds, const std::vector<int>& context, double timeStep, int tsId)
 618 {
 619   if(!ds || !ret)
 620     throw MZCException("ConvertFromRectilinearGrid : internal error !");
 621   //
 622   MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
 623   ret->setMeshes(meshes);
 624   MCAuto<MEDFileFields> fields(MEDFileFields::New());
 625   ret->setFields(fields);
 626   //
 627   MCAuto<MEDFileCMesh> cmesh(MEDFileCMesh::New());
 628   meshes->pushMesh(cmesh);
 629   MCAuto<MEDCouplingCMesh> cmeshmc(MEDCouplingCMesh::New());
 630   vtkDataArray *cx(ds->GetXCoordinates()),*cy(ds->GetYCoordinates()),*cz(ds->GetZCoordinates());
 631   if(cx)
 632     {
 633       MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cx));
 634       cmeshmc->setCoordsAt(0,arr);
 635     }
 636   if(cy)
 637     {
 638       MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cy));
 639       cmeshmc->setCoordsAt(1,arr);
 640     }
 641   if(cz)
 642     {
 643       MCAuto<DataArrayDouble> arr(ConvertVTKArrayToMCArrayDoubleForced(cz));
 644       cmeshmc->setCoordsAt(2,arr);
 645     }
 646   std::string meshName(GetMeshNameWithContext(context));
 647   cmeshmc->setName(meshName);
 648   cmesh->setMesh(cmeshmc);
 649   std::vector<MCAuto<DataArray> > cellFs(AddPartFields(0,ds->GetCellData()));
 650   for(std::vector<MCAuto<DataArray> >::const_iterator it=cellFs.begin();it!=cellFs.end();it++)
 651     {
 652       MCAuto<DataArray> da(*it);
 653       AppendMCFieldFrom(MEDCoupling::ON_CELLS,cmeshmc,ret,da,NULL,timeStep,tsId);
 654     }
 655   std::vector<MCAuto<DataArray> > nodeFs(AddPartFields(0,ds->GetPointData()));
 656   for(std::vector<MCAuto<DataArray> >::const_iterator it=nodeFs.begin();it!=nodeFs.end();it++)
 657     {
 658       MCAuto<DataArray> da(*it);
 659       AppendMCFieldFrom(MEDCoupling::ON_NODES,cmeshmc,ret,da,NULL,timeStep,tsId);
 660     }
 661 }
 662
 663 void ConvertFromPolyData(MEDFileData *ret, vtkPolyData *ds, const std::vector<int>& context, double timeStep, int tsId)
 664 {
 665   if(!ds || !ret)
 666     throw MZCException("ConvertFromPolyData : internal error !");
 667   //
 668   MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
 669   ret->setMeshes(meshes);
 670   MCAuto<MEDFileFields> fields(MEDFileFields::New());
 671   ret->setFields(fields);
 672   //
 673   MCAuto<MEDFileUMesh> umesh(MEDFileUMesh::New());
 674   meshes->pushMesh(umesh);
 675   MCAuto<DataArrayDouble> coords(BuildCoordsFrom(ds));
 676   umesh->setCoords(coords);
 677   umesh->setName(GetMeshNameWithContext(context));
 678   //
 679   int offset(0);
 680   std::vector< MicroField > ms;
 681   vtkCellArray *cd(ds->GetVerts());
 682   if(cd)
 683     {
 684       MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArray(cd,coords,0,INTERP_KERNEL::NORM_POINT1));
 685       if((const MEDCouplingUMesh *)subMesh)
 686         {
 687           std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
 688           offset+=subMesh->getNumberOfCells();
 689           ms.push_back(MicroField(subMesh,cellFs));
 690         }
 691     }
 692   vtkCellArray *cc(ds->GetLines());
 693   if(cc)
 694     {
 695       MCAuto<MEDCouplingUMesh> subMesh;
 696       try
 697         {
 698           subMesh=BuildMeshFromCellArray(cc,coords,1,INTERP_KERNEL::NORM_SEG2);
 699         }
 700       catch(INTERP_KERNEL::Exception& e)
 701         {
 702           std::ostringstream oss; oss << "MEDWriter does not manage polyline cell type because MED file format does not support it ! Maybe it is the source of the problem ? The cause of this exception was " << e.what() << std::endl;
 703           throw INTERP_KERNEL::Exception(oss.str());
 704         }
 705       if((const MEDCouplingUMesh *)subMesh)
 706         {
 707           std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
 708           offset+=subMesh->getNumberOfCells();
 709           ms.push_back(MicroField(subMesh,cellFs));
 710         }
 711     }
 712   vtkCellArray *cb(ds->GetPolys());
 713   if(cb)
 714     {
 715       MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArray(cb,coords,2,INTERP_KERNEL::NORM_POLYGON));
 716       if((const MEDCouplingUMesh *)subMesh)
 717         {
 718           std::vector<MCAuto<DataArray> > cellFs(AddPartFields2(offset,offset+subMesh->getNumberOfCells(),ds->GetCellData()));
 719           offset+=subMesh->getNumberOfCells();
 720           ms.push_back(MicroField(subMesh,cellFs));
 721         }
 722     }
 723   vtkCellArray *ca(ds->GetStrips());
 724   if(ca)
 725     {
 726       MCAuto<DataArrayIdType> ids;
 727       MCAuto<MEDCouplingUMesh> subMesh(BuildMeshFromCellArrayTriangleStrip(ca,coords,ids));
 728       if((const MEDCouplingUMesh *)subMesh)
 729         {
 730           std::vector<MCAuto<DataArray> > cellFs(AddPartFields(ids,ds->GetCellData()));
 731           offset+=subMesh->getNumberOfCells();
 732           ms.push_back(MicroField(subMesh,cellFs));
 733         }
 734     }
 735   AssignSingleGTMeshes(ret,ms,timeStep,tsId);
 736   AddNodeFields(ret,ds->GetPointData(),timeStep,tsId);
 737 }
 738
 739 void ConvertFromUnstructuredGrid(MEDFileData *ret, vtkUnstructuredGrid *ds, const std::vector<int>& context, double timeStep, int tsId)
 740 {
 741   if(!ds || !ret)
 742     throw MZCException("ConvertFromUnstructuredGrid : internal error !");
 743   //
 744   MCAuto<MEDFileMeshes> meshes(MEDFileMeshes::New());
 745   ret->setMeshes(meshes);
 746   MCAuto<MEDFileFields> fields(MEDFileFields::New());
 747   ret->setFields(fields);
 748   //
 749   MCAuto<MEDFileUMesh> umesh(MEDFileUMesh::New());
 750   meshes->pushMesh(umesh);
 751   MCAuto<DataArrayDouble> coords(BuildCoordsFrom(ds));
 752   umesh->setCoords(coords);
 753   umesh->setName(GetMeshNameWithContext(context));
 754   vtkIdType nbCells(ds->GetNumberOfCells());
 755   vtkCellArray *ca(ds->GetCells());
 756   if(!ca)
 757     return ;
 758   //vtkIdType nbEnt(ca->GetNumberOfConnectivityEntries()); // todo: unused
 759   //vtkIdType *caPtr(ca->GetData()->GetPointer(0)); // todo: unused
 760   vtkUnsignedCharArray *ct(ds->GetCellTypesArray());
 761   if(!ct)
 762     throw MZCException("ConvertFromUnstructuredGrid : internal error");
 763   vtkIdTypeArray *cla(ds->GetCellLocationsArray());
 764   //const vtkIdType *claPtr(cla->GetPointer(0)); // todo: unused
 765   if(!cla)
 766     throw MZCException("ConvertFromUnstructuredGrid : internal error 2");
 767   const unsigned char *ctPtr(ct->GetPointer(0));
 768   std::map<int,int> m(ComputeMapOfType());
 769   MCAuto<DataArrayInt> lev(DataArrayInt::New()) ;  lev->alloc(nbCells,1);
 770   int *levPtr(lev->getPointer());
 771   for(vtkIdType i=0;i<nbCells;i++)
 772     {
 773       std::map<int,int>::iterator it(m.find(ctPtr[i]));
 774       if(it!=m.end())
 775         {
 776           const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel((INTERP_KERNEL::NormalizedCellType)(*it).second));
 777           levPtr[i]=cm.getDimension();
 778         }
 779       else
 780         {
 781           if(ctPtr[i]==VTK_POLY_VERTEX)
 782             {
 783               const INTERP_KERNEL::CellModel& cm(INTERP_KERNEL::CellModel::GetCellModel(INTERP_KERNEL::NORM_POINT1));
 784               levPtr[i]=cm.getDimension();
 785             }
 786           else
 787             {
 788               std::ostringstream oss; oss << "ConvertFromUnstructuredGrid : at pos #" << i << " unrecognized VTK cell with type =" << ctPtr[i];
 789               throw MZCException(oss.str());
 790             }
 791         }
 792     }
 793   MCAuto<DataArrayInt> levs(lev->getDifferentValues());
 794   std::vector< MicroField > ms;
 795   //vtkIdTypeArray *faces(ds->GetFaces()),*faceLoc(ds->GetFaceLocations()); // todo: unused
 796   for(const int *curLev=levs->begin();curLev!=levs->end();curLev++)
 797     {
 798       MCAuto<MEDCouplingUMesh> m0(MEDCouplingUMesh::New("",*curLev));
 799       m0->setCoords(coords); m0->allocateCells();
 800       MCAuto<DataArrayIdType> cellIdsCurLev(lev->findIdsEqual(*curLev));
 801       for(const mcIdType *cellId=cellIdsCurLev->begin();cellId!=cellIdsCurLev->end();cellId++)
 802         {
 803           int vtkType(ds->GetCellType(*cellId));
 804           std::map<int,int>::iterator it(m.find(vtkType));
 805           INTERP_KERNEL::NormalizedCellType ct=it!=m.end()?(INTERP_KERNEL::NormalizedCellType)((*it).second):INTERP_KERNEL::NORM_POINT1;
 806           if(ct!=INTERP_KERNEL::NORM_POLYHED && vtkType!=VTK_POLY_VERTEX)
 807             {
 808               vtkIdType sz(0);
 809               const vtkIdType *pts(nullptr);
 810               ds->GetCellPoints(*cellId, sz, pts);
 811               std::vector<mcIdType> conn2(pts,pts+sz);
 812               m0->insertNextCell(ct,sz,conn2.data());
 813             }
 814           else if(ct==INTERP_KERNEL::NORM_POLYHED)
 815             {
 816               // # de faces du polyèdre
 817               vtkIdType nbOfFaces(0);
 818               // connectivé des faces (numFace0Pts, id1, id2, id3, numFace1Pts,id1, id2, id3, ...)
 819               const vtkIdType *facPtr(nullptr);
 820               ds->GetFaceStream(*cellId, nbOfFaces, facPtr);
 821               std::vector<mcIdType> conn;
 822               for(vtkIdType k=0;k<nbOfFaces;k++)
 823                 {
 824                   vtkIdType nbOfNodesInFace(*facPtr++);
 825                   std::copy(facPtr,facPtr+nbOfNodesInFace,std::back_inserter(conn));
 826                   if(k<nbOfFaces-1)
 827                     conn.push_back(-1);
 828                   facPtr+=nbOfNodesInFace;
 829                 }
 830               m0->insertNextCell(ct,ToIdType(conn.size()),&conn[0]);
 831             }
 832           else
 833             {
 834               vtkIdType sz(0);
 835               const vtkIdType *pts(nullptr);
 836               ds->GetCellPoints(*cellId, sz, pts);
 837               if(sz!=1)
 838                 throw MZCException("ConvertFromUnstructuredGrid : non single poly vertex not managed by MED !");
 839               m0->insertNextCell(ct,1,(const mcIdType*)pts);
 840             }
 841         }
 842       std::vector<MCAuto<DataArray> > cellFs(AddPartFields(cellIdsCurLev,ds->GetCellData()));
 843       ms.push_back(MicroField(m0,cellFs));
 844     }
 845   AssignSingleGTMeshes(ret,ms,timeStep,tsId);
 846   AddNodeFields(ret,ds->GetPointData(),timeStep,tsId);
 847 }
 848
 849 ///////////////////
 850
 851 void WriteMEDFileFromVTKDataSet(MEDFileData *mfd, vtkDataSet *ds, const std::vector<int>& context, double timeStep, int tsId)
 852 {
 853   if(!ds || !mfd)
 854     throw MZCException("Internal error in WriteMEDFileFromVTKDataSet.");
 855   vtkPolyData *ds2(vtkPolyData::SafeDownCast(ds));
 856   vtkUnstructuredGrid *ds3(vtkUnstructuredGrid::SafeDownCast(ds));
 857   vtkRectilinearGrid *ds4(vtkRectilinearGrid::SafeDownCast(ds));
 858   if(ds2)
 859     {
 860       ConvertFromPolyData(mfd,ds2,context,timeStep,tsId);
 861     }
 862   else if(ds3)
 863     {
 864       ConvertFromUnstructuredGrid(mfd,ds3,context,timeStep,tsId);
 865     }
 866   else if(ds4)
 867     {
 868       ConvertFromRectilinearGrid(mfd,ds4,context,timeStep,tsId);
 869     }
 870   else
 871     throw MZCException("Unrecognized vtkDataSet ! Sorry ! Try to convert it to UnstructuredGrid to be able to write it !");
 872 }
 873
 874 void WriteMEDFileFromVTKMultiBlock(MEDFileData *mfd, vtkMultiBlockDataSet *ds, const std::vector<int>& context, double timeStep, int tsId)
 875 {
 876   if(!ds || !mfd)
 877     throw MZCException("Internal error in WriteMEDFileFromVTKMultiBlock.");
 878   int nbBlocks(ds->GetNumberOfBlocks());
 879   if(nbBlocks==1 && context.empty())
 880     {
 881       vtkDataObject *uniqueElt(ds->GetBlock(0));
 882       if(!uniqueElt)
 883         throw MZCException("Unique elt in multiblock is NULL !");
 884       vtkDataSet *uniqueEltc(vtkDataSet::SafeDownCast(uniqueElt));
 885       if(uniqueEltc)
 886         {
 887           WriteMEDFileFromVTKDataSet(mfd,uniqueEltc,context,timeStep,tsId);
 888           return ;
 889         }
 890     }
 891   for(int i=0;i<nbBlocks;i++)
 892     {
 893       vtkDataObject *elt(ds->GetBlock(i));
 894       std::vector<int> context2;
 895       context2.push_back(i);
 896       if(!elt)
 897         {
 898           std::ostringstream oss; oss << "In context ";
 899           std::copy(context.begin(),context.end(),std::ostream_iterator<int>(oss," "));
 900           oss << " at pos #" << i << " elt is NULL !";
 901           throw MZCException(oss.str());
 902         }
 903       vtkDataSet *elt1(vtkDataSet::SafeDownCast(elt));
 904       if(elt1)
 905         {
 906           WriteMEDFileFromVTKDataSet(mfd,elt1,context,timeStep,tsId);
 907           continue;
 908         }
 909       vtkMultiBlockDataSet *elt2(vtkMultiBlockDataSet::SafeDownCast(elt));
 910       if(elt2)
 911         {
 912           WriteMEDFileFromVTKMultiBlock(mfd,elt2,context,timeStep,tsId);
 913           continue;
 914         }
 915       std::ostringstream oss; oss << "In context ";
 916       std::copy(context.begin(),context.end(),std::ostream_iterator<int>(oss," "));
 917       oss << " at pos #" << i << " elt not recognized data type !";
 918       throw MZCException(oss.str());
 919     }
 920 }
 921
 922 void WriteMEDFileFromVTKGDS(MEDFileData *mfd, vtkDataObject *input, double timeStep, int tsId)
 923 {
 924   if(!input || !mfd)
 925     throw MZCException("WriteMEDFileFromVTKGDS : internal error !");
 926   std::vector<int> context;
 927   vtkDataSet *input1(vtkDataSet::SafeDownCast(input));
 928   if(input1)
 929     {
 930       WriteMEDFileFromVTKDataSet(mfd,input1,context,timeStep,tsId);
 931       return ;
 932     }
 933   vtkMultiBlockDataSet *input2(vtkMultiBlockDataSet::SafeDownCast(input));
 934   if(input2)
 935     {
 936       WriteMEDFileFromVTKMultiBlock(mfd,input2,context,timeStep,tsId);
 937       return ;
 938     }
 939   throw MZCException("WriteMEDFileFromVTKGDS : not recognized data type !");
 940 }
 941
 942 void PutFamGrpInfoIfAny(MEDFileData *mfd, const std::string& meshName, const std::vector<Grp>& groups, const std::vector<Fam>& fams)
 943 {
 944   if(!mfd)
 945     return ;
 946   if(meshName.empty())
 947     return ;
 948   MEDFileMeshes *meshes(mfd->getMeshes());
 949   if(!meshes)
 950     return ;
 951   if(meshes->getNumberOfMeshes()!=1)
 952     return ;
 953   MEDFileMesh *mm(meshes->getMeshAtPos(0));
 954   if(!mm)
 955     return ;
 956   mm->setName(meshName);
 957   for(std::vector<Fam>::const_iterator it=fams.begin();it!=fams.end();it++)
 958     mm->setFamilyId((*it).getName(),(*it).getID());
 959   for(std::vector<Grp>::const_iterator it=groups.begin();it!=groups.end();it++)
 960     mm->setFamiliesOnGroup((*it).getName(),(*it).getFamilies());
 961   MEDFileFields *fields(mfd->getFields());
 962   if(!fields)
 963     return ;
 964   for(int i=0;i<fields->getNumberOfFields();i++)
 965     {
 966       MEDFileAnyTypeFieldMultiTS *fmts(fields->getFieldAtPos(i));
 967       if(!fmts)
 968         continue;
 969       fmts->setMeshName(meshName);
 970     }
 971 }