src/CONVERTOR/VISU_Vtk2MedConvertor.cxx

   1 // Copyright (C) 2007-2011  CEA/DEN, EDF R&D, OPEN CASCADE
   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.
   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
  20 //  File   : VISU_Vtk2MedConvertor.cxx
  21 //  Author : Eugeny NIKOLAEV, Open CASCADE SAS
  22 //
  23 #include "VISU_Vtk2MedConvertor.hxx"
  24
  25 // QT includes
  26 #include <qdir.h>
  27 #include <qfileinfo.h>
  28 #include <qstringlist.h>
  29
  30 // VTK includes
  31 #include <vtkSmartPointer.h>
  32 #include <vtkDataReader.h>
  33
  34 #include <vtkStructuredPointsReader.h>
  35 #include <vtkStructuredGridReader.h>
  36 #include <vtkRectilinearGridReader.h>
  37 #include <vtkUnstructuredGridReader.h>
  38 #include <vtkPolyDataReader.h>
  39 #include <vtkDataSetReader.h>
  40
  41 #include <vtkPointSet.h>
  42 #include <vtkDataSet.h>
  43 #include <vtkPolyData.h>
  44 #include <vtkDataObject.h>
  45 #include <vtkCellTypes.h>
  46 #include <vtkCellType.h>
  47 #include <vtkCell.h>
  48 #include <vtkCellData.h>
  49 #include <vtkPointData.h>
  50
  51 #include <vtkUnstructuredGrid.h>
  52 #include <vtkFloatArray.h>
  53
  54 // MED Warpper includes
  55 #include "MED_Factory.hxx"
  56
  57 // STL includes
  58 #include <string>
  59 #include <vector>
  60 #include <iostream>
  61 #include <sstream>
  62 #include <set>
  63 #include <map>
  64 //#include <pair.h>
  65
  66 #ifdef _DEBUG_
  67 static int MYDEBUG = 0;
  68 static int MYDEBUG_VALUES = 0;
  69 #else
  70 static int MYDEBUG = 0;
  71 static int MYDEBUG_VALUES = 0;
  72 #endif
  73
  74 using namespace std;
  75
  76 /*
  77 # = dynamic_cast</*
  78 *>(  );define   VTK_EMPTY_CELL   0
  79 #define         VTK_VERTEX   1
  80 #define         VTK_POLY_VERTEX   2
  81 #define         VTK_LINE   3
  82 #define         VTK_POLY_LINE   4
  83 #define         VTK_TRIANGLE   5
  84 #define         VTK_TRIANGLE_STRIP   6
  85 #define         VTK_POLYGON   7
  86 #define         VTK_PIXEL   8
  87 #define         VTK_QUAD   9
  88 #define         VTK_TETRA   10
  89 #define         VTK_VOXEL   11
  90 #define         VTK_HEXAHEDRON   12
  91 #define         VTK_WEDGE   13
  92 #define         VTK_PYRAMID   14
  93 #define         VTK_QUADRATIC_EDGE   21
  94 #define         VTK_QUADRATIC_TRiTSIANGLE   22
  95 #define         VTK_QUADRATIC_QUAD   23
  96 #define         VTK_QUADRATIC_TETRA   24
  97 #define         VTK_QUADRATIC_HEXAHEDRON   25
  98 #define         VTK_CONVEX_POINT_SET   41
  99 #define         VTK_PARAMETRIC_CURVE   51
 100 #define         VTK_PARAMETRIC_SURFACE   52
 101 #define         VTK_PARAMETRIC_TRI_SURFACE   53
 102 #define         VTK_PARAMETRIC_QUAD_SURFACE   54
 103 #define         VTK_PARAMETRIC_TETRA_REGION   55
 104 #define         VTK_PARAMETRIC_HEX_REGION   56
 105 */
 106
 107 static MED::EGeometrieElement VTK2MED( const int theGeom )
 108 {
 109   // Ignoring vtk types:
 110   // VTK_PIXEL,
 111   // VTK_VERTEX,
 112   // VTK_POLY_VERTEX
 113   // VTK_VOXEL
 114   // VTK_POLY_LINE
 115   // VTK_TRIANGLE_STRIP
 116   // VTK_PARAMETRIC_CURVE
 117   // VTK_PARAMETRIC_SURFACE
 118   // VTK_PARAMETRIC_TRI_SURFACE
 119   // VTK_PARAMETRIC_QUAD_SURFACE
 120   // VTK_PARAMETRIC_TETRA_REGION
 121   // VTK_PARAMETRIC_HEX_REGION
 122
 123   MED::EGeometrieElement aEmptyGeom = MED::EGeometrieElement(-1);
 124   switch(theGeom){
 125   case VTK_LINE:             return MED::eSEG2;
 126   case VTK_TRIANGLE:         return MED::eTRIA3;
 127   case VTK_POLYGON:          return MED::ePOLYGONE;
 128   case VTK_QUAD:             return MED::eQUAD4;
 129   case VTK_TETRA:            return MED::eTETRA4;
 130   case VTK_HEXAHEDRON:       return MED::eHEXA8;
 131   case VTK_WEDGE:            return MED::ePENTA6;
 132   case VTK_PYRAMID:          return MED::ePYRA5;
 133   case VTK_HEXAGONAL_PRISM:  return MED::eOCTA12;
 134   // QUADRATIC elements
 135   case VTK_QUADRATIC_EDGE:      return MED::eSEG3;
 136   case VTK_QUADRATIC_TRIANGLE:  return MED::eTRIA6;
 137   case VTK_QUADRATIC_QUAD:      return MED::eQUAD8;
 138   case VTK_BIQUADRATIC_QUAD:      return MED::eQUAD9;
 139   case VTK_QUADRATIC_TETRA:     return MED::eTETRA10;
 140   case VTK_QUADRATIC_HEXAHEDRON:return MED::eHEXA20;
 141   case VTK_TRIQUADRATIC_HEXAHEDRON :return MED::eHEXA27;
 142   case VTK_CONVEX_POINT_SET:    return MED::ePOLYEDRE;
 143   }
 144
 145   return aEmptyGeom;
 146 }
 147
 148 /*!
 149   \class VISU_Vtk2MedConvertor
 150   \brief The general main of the class VISU_Vtk2MedConvertor is converting from
 151   one or several VTK files to the one MED file...
 152
 153   The VISU_Vtk2MedConvertor interface allows us to create the MED file according
 154   to VTK files in next ways:
 155   - Extract geometry and fields from one VTK file.
 156   - Extract geometry and fields from first VTK file and fields from others
 157   VTK files (geometry ignoring). Also the fields which have same names join into
 158   corresponding fields with different time stamp
 159 */
 160
 161 /*!
 162   \brief Constructor
 163   - Sets default output mesh name
 164   - Sets default version of output MED file.
 165   - Sets default ignoring fields list
 166   - Sets default points and cells ids mapping field names.
 167
 168 */
 169 VISU_Vtk2MedConvertor
 170 ::VISU_Vtk2MedConvertor():
 171   myVersion(MED::eV2_2),
 172   myMeshName("vtk2med")
 173 {
 174   myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
 175   myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
 176   myIgnoringFieldList.insert("VISU_FIELD");
 177   setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
 178   setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
 179 }
 180
 181 /*!
 182   \brief Constructor
 183   - Sets default output mesh name
 184   - Sets default version of output MED file.
 185   - Sets default ignoring fields list
 186   - Sets default points and cells ids mapping field names.
 187
 188   \param theMEDFileName output med file name
 189   \param theFirstVTKFileName first vtk file name
 190 */
 191 VISU_Vtk2MedConvertor
 192 ::VISU_Vtk2MedConvertor(const string theMEDFileName,
 193                     const string theFirstVTKFileName):
 194   myVersion(MED::eV2_2),
 195   myMeshName("vtk2med")
 196 {
 197   myMEDFileName   = theMEDFileName;
 198   myFirstVTKFileName = theFirstVTKFileName;
 199   myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
 200   myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
 201   myIgnoringFieldList.insert("VISU_FIELD");
 202   setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
 203   setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
 204 }
 205
 206 /*!
 207   \brief Constructor
 208   - Sets default output mesh name
 209   - Sets default version of output MED file.
 210   - Sets default ignoring fields list
 211   - Sets default points and cells ids mapping field names.
 212
 213   \param theMEDFileName output med file name
 214   \param theFirstVTKFileName first vtk file name
 215   \param theDataVTKFileNames of vtk file names, which will be using as values on points and cells
 216 */
 217 VISU_Vtk2MedConvertor
 218 ::VISU_Vtk2MedConvertor(const string theMEDFileName,
 219                     const string theFirstVTKFileName,
 220                     const TVectorString theDataVTKFileNames):
 221   myVersion(MED::eV2_2),
 222   myMeshName("vtk2med")
 223 {
 224   myMEDFileName      = theMEDFileName;
 225   myFirstVTKFileName = theFirstVTKFileName;
 226   myDataVTKFileNames = theDataVTKFileNames;
 227   myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
 228   myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
 229   myIgnoringFieldList.insert("VISU_FIELD");
 230   setMeshName("vtk2med");
 231   setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
 232   setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
 233 }
 234
 235 /*!
 236   \brief Adds field names, which used as specific fields with ids or elements
 237   (or something else). (Default: "VISU_CELLS_MAPPER","VISU_POINTS_MAPPER","VISU_FILED")
 238   \param theFieldName field name
 239   \sa eraseFromIgnoringFieldList()
 240 */
 241 void
 242 VISU_Vtk2MedConvertor
 243 ::addToIgnoringFieldList( const string& theFieldName )
 244 {
 245   myIgnoringFieldList.insert(theFieldName);
 246 }
 247
 248 /*!
 249   \brief Sets the output MED file name
 250   \param theFileName file name
 251   \sa getMEDFileName()
 252 */
 253 void
 254 VISU_Vtk2MedConvertor
 255 ::setMEDFileName( const string theFileName )
 256 {
 257   myMEDFileName=theFileName;
 258 };
 259
 260 /*!
 261   \brief Gets the output MED file name
 262   \return output MED file name
 263   \sa setMEDFileName()
 264 */
 265 string
 266 VISU_Vtk2MedConvertor
 267 ::getMEDFileName() const
 268 {
 269   return myMEDFileName;
 270 }
 271
 272 /*!
 273   \brief Sets the first input vtk file name
 274   \param theFileName file name
 275   \sa getFirstVTKFileName()
 276 */
 277 void
 278 VISU_Vtk2MedConvertor
 279 ::setFirstVTKFileName( const string theFileName )
 280 {
 281   myFirstVTKFileName=theFileName;
 282 };
 283
 284 /*!
 285   \brief Fets the first input vtk file name
 286   \return first input vtk file name
 287   \sa setFirstVTKFileName()
 288 */
 289 string
 290 VISU_Vtk2MedConvertor
 291 ::getFirstVTKFileName() const
 292 {
 293   return myFirstVTKFileName;
 294 }
 295
 296 /*!
 297   \brief Sets list of vtk file names, which will be using as values on points and cells
 298   \param theFileNames list of vtk file names
 299   \sa getDataVTKFileNames()
 300 */
 301 void
 302 VISU_Vtk2MedConvertor
 303 ::setDataVTKFileNames( const TVectorString theFileNames )
 304 {
 305   myDataVTKFileNames = theFileNames;
 306 };
 307
 308 /*!
 309   \brief Gets list of vtk file names, which will be using as values on points and cells
 310   \param theFileNames out list of vtk file names
 311   \sa setDataVTKFileNames()
 312 */
 313 void
 314 VISU_Vtk2MedConvertor
 315 ::getDataVTKFileNames( TVectorString& theDataVTKFileNames ) const
 316 {
 317   theDataVTKFileNames = myDataVTKFileNames;
 318 };
 319
 320 /*!
 321   \brief Sets version of the output MED file MED::V2_2(is default) or MED::V2_1
 322   \param theVersion version of the output MED file
 323 */
 324 void
 325 VISU_Vtk2MedConvertor
 326 ::setVersion( const MED::EVersion theVersion )
 327 {
 328   myVersion = theVersion;
 329 }
 330
 331 /*!
 332   \brief Gets version of the output MED file MED::V2_2(is default) or MED::V2_1
 333   \return version of the output MED file
 334 */
 335 MED::EVersion
 336 VISU_Vtk2MedConvertor
 337 ::getVersion() const
 338 {
 339   return myVersion;
 340 }
 341
 342 /*!
 343   \brief Sets output mesh name. ("vtk2med" - default)
 344   \param theMeshName mesh name
 345   \sa getMeshName()
 346 */
 347 void
 348 VISU_Vtk2MedConvertor
 349 ::setMeshName( const string theMeshName )
 350 {
 351   myMeshName = theMeshName;
 352 }
 353
 354 /*!
 355   \brief Gets output mesh name. ("vtk2med" - default)
 356   \return mesh name
 357   \sa setMeshName()
 358 */
 359 string
 360 VISU_Vtk2MedConvertor
 361 ::getMeshName() const
 362 {
 363   return myMeshName;
 364 }
 365
 366 /*!
 367   \brief Sets field name with cell ids (Default - VISU_CELLS_MAPPER)
 368   \param theFieldName field name
 369   \sa getCellDataFieldNameIDS
 370 */
 371 void
 372 VISU_Vtk2MedConvertor
 373 ::setCellDataFieldNameIDS( const string& theFieldName )
 374 {
 375   myCellDataFieldNameIDS = theFieldName;
 376 }
 377
 378 /*!
 379   \brief Gets field name with cell ids (Default - VISU_CELLS_MAPPER)
 380   \return field name
 381   \sa setCellDataFieldNameIDS()
 382 */
 383 const string&
 384 VISU_Vtk2MedConvertor
 385 ::getCellDataFieldNameIDS() const
 386 {
 387   return myCellDataFieldNameIDS;
 388 }
 389
 390 /*!
 391   \brief Erases field names which used as specific fields with ids or elements
 392   (or something else)
 393   \param theFieldName field name
 394   \sa addToIgnoringFieldList()
 395 */
 396 void
 397 VISU_Vtk2MedConvertor
 398 ::eraseFromIgnoringFieldList(const string& theFieldName)
 399 {
 400   myIgnoringFieldList.erase(theFieldName);
 401 }
 402
 403 /*!
 404   \brief Gets list of field names which used as specific fields with ids or elements
 405   \return list of field names
 406 */
 407 const std::set<std::string>&
 408 VISU_Vtk2MedConvertor
 409 ::getIgnoringFieldList() const
 410 {
 411   return myIgnoringFieldList;
 412 }
 413
 414 /*!
 415   \brief Sets field name with point ids
 416   \param theFieldName field name
 417   \sa getPointDataFieldNameIDS()
 418 */
 419 void
 420 VISU_Vtk2MedConvertor
 421 ::setPointDataFieldNameIDS( const string& theFieldName )
 422 {
 423   myPointDataFieldNameIDS = theFieldName;
 424 }
 425
 426 /*!
 427   \brief Gets field name with point ids
 428   \return field name
 429   \sa setPointDataFieldNameIDS()
 430 */
 431 const string&
 432 VISU_Vtk2MedConvertor
 433 ::getPointDataFieldNameIDS() const
 434 {
 435   return myPointDataFieldNameIDS;
 436 }
 437
 438 /*!
 439   \brief Sets values of time stamps If this array is not specified values of time
 440   stamps are generated automatically ( 0, 1, 2 ... )
 441   \param theTStamps vector of time stamps
 442   \sa getTimeStamps()
 443 */
 444 void
 445 VISU_Vtk2MedConvertor
 446 ::setTimeStamps( const TVectorDouble& theTStamps )
 447 {
 448   myTStamps = theTStamps;
 449 }
 450
 451 /*!
 452   \brief Gets values of time stamps If this array is not specified values of time
 453   stamps are generated automatically ( 0, 1, 2 ... )
 454   \param theTStamps out vector of time stamps
 455   \sa setTimeStamps()
 456 */
 457 void
 458 VISU_Vtk2MedConvertor
 459 ::getTimeStamps( TVectorDouble& theTStamps ) const
 460 {
 461   theTStamps = myTStamps;
 462 }
 463
 464 /*!
 465   \brief Retrieves identifiers of cells from input data set corresponding to given type
 466   \param theInput input data set
 467   \param type type
 468   \param array out array
 469 */
 470 void
 471 VISU_Vtk2MedConvertor
 472 ::GetIdsOfCellsOfType(vtkDataSet* theInput,
 473                       const int type,
 474                       vtkIntArray *array)
 475 {
 476   for (int cellId = 0; cellId < theInput->GetNumberOfCells(); cellId++)
 477     if (theInput->GetCellType(cellId) == type)
 478       array->InsertNextValue(cellId);
 479 }
 480
 481 /*!
 482   \brief Creates elements (private auxiliary method)
 483   \return 0 if operation has been completed successfully, 1 otherwise
 484 */
 485 int
 486 VISU_Vtk2MedConvertor
 487 ::CreateElements(vtkDataSet* theInput,
 488                  MED::PMeshInfo theMeshInfo,
 489                  MED::PWrapper  theMed,
 490                  vtkIntArray* theCellsMapper,
 491                  MED::EEntiteMaillage theEntity,
 492                  int theVTKGeom,
 493                  int nbPointsInGeom,
 494                  std::vector<int>& theNumberingConvertor,
 495                  TGeom2CellIds& theGeom2CellIdMap)
 496 {
 497   bool aIdsConvert = (theNumberingConvertor.size() > 0);
 498   vtkIntArray* aCellIds = vtkIntArray::New();
 499   const MED::EConnectivite aConnType = MED::eNOD;
 500
 501   MED::TIntVector aConn;
 502   MED::TIntVector aFamilyNums;// -1
 503   MED::TIntVector aElemNums;
 504
 505   this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
 506   int nbElems = aCellIds->GetNumberOfTuples();
 507   if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
 508   aConn.reserve(nbElems*nbPointsInGeom);
 509   if(nbElems>0){
 510     TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
 511     int* aPointer = aCellIds->GetPointer(0);
 512     for(int i=0;i<nbElems;i++,aPointer++){
 513       int aCellId = *aPointer;
 514       aCellIdsMapper.push_back(aCellId);
 515       vtkCell* aCell = theInput->GetCell(aCellId);
 516       int nbPointsInCell = aCell->GetNumberOfPoints();
 517       if(nbPointsInCell!=nbPointsInGeom){
 518         cout << "Error in file=|" << __FILE__<<"| line:[" << __LINE__ << "]" << endl;
 519         cout << "Must be "<<nbPointsInGeom<<" nodes in VTK Geometry:"<<theVTKGeom<<" element" << endl;
 520         aCellIds->Delete();
 521         return 1; // exit
 522       }
 523       aFamilyNums.push_back(-1);
 524       for(int j=0;j<nbPointsInCell;j++){
 525         if (aIdsConvert)
 526           aConn.push_back(aCell->GetPointId(theNumberingConvertor[j])+1);
 527         else
 528           aConn.push_back(aCell->GetPointId(j)+1);
 529       }
 530       if(theCellsMapper){
 531         if(theCellsMapper->GetNumberOfComponents()==2)
 532           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
 533         else if(theCellsMapper->GetNumberOfComponents()==1)
 534           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
 535       }
 536     }
 537
 538
 539     MED::PCellInfo aCellInfo = theMed->CrCellInfo(theMeshInfo,
 540                                                   theEntity,
 541                                                   VTK2MED(theVTKGeom),
 542                                                   aConn,
 543                                                   aConnType,
 544                                                   aFamilyNums,
 545                                                   aElemNums);
 546     theMed->SetCellInfo(aCellInfo);
 547   }
 548
 549   aCellIds->Delete();
 550
 551   return 0;
 552 }
 553
 554 /*!
 555   \brief Creates polygons (private auxiliary method)
 556   \return 0 if operation has been completed successfully, 1 otherwise
 557 */
 558 int
 559 VISU_Vtk2MedConvertor
 560 ::CreatePolygons(vtkDataSet* theInput,
 561                  MED::PMeshInfo theMeshInfo,
 562                  MED::PWrapper  theMed,
 563                  vtkIntArray* theCellsMapper,
 564                  MED::EEntiteMaillage theEntity,
 565                  TGeom2CellIds& theGeom2CellIdMap)
 566 {
 567   int theVTKGeom = VTK_POLYGON;
 568   vtkIntArray* aCellIds = vtkIntArray::New();
 569   const MED::EConnectivite aConnType = MED::eNOD;
 570
 571   MED::TIntVector aConn;
 572   MED::TIntVector aFamilyNums;// -1
 573   MED::TIntVector aElemNums;
 574   MED::TIntVector aPolygoneInds;
 575   aPolygoneInds.push_back(1); // reference on the first element in the connectivities
 576
 577   this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
 578   int nbElems = aCellIds->GetNumberOfTuples();
 579   if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
 580   if(nbElems>0){
 581     TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
 582     int* aPointer = aCellIds->GetPointer(0);
 583     for(int i=0;i<nbElems;i++,aPointer++){
 584       int aCellId = *aPointer;
 585       aCellIdsMapper.push_back(aCellId);
 586       vtkCell* aCell = theInput->GetCell(aCellId);
 587       int nbPointsInCell = aCell->GetNumberOfPoints();
 588       aFamilyNums.push_back(-1);
 589       int aPrevPos = aPolygoneInds.back();
 590       aPolygoneInds.push_back(aPrevPos+nbPointsInCell);
 591       for(int j=0;j<nbPointsInCell;j++)
 592         aConn.push_back(aCell->GetPointId(j)+1);
 593       if(theCellsMapper){
 594         if(theCellsMapper->GetNumberOfComponents()==2)
 595           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
 596         else if(theCellsMapper->GetNumberOfComponents()==1)
 597           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
 598       }
 599     }
 600
 601
 602     MED::PPolygoneInfo aCellInfo = theMed->CrPolygoneInfo(theMeshInfo,
 603                                                           theEntity,
 604                                                           VTK2MED(theVTKGeom),
 605                                                           aPolygoneInds,
 606                                                           aConn,
 607                                                           aConnType,
 608                                                           aFamilyNums,
 609                                                           aElemNums);
 610     theMed->SetPolygoneInfo(aCellInfo);
 611   }
 612
 613   aCellIds->Delete();
 614
 615   return 0;
 616 }
 617
 618 /*!
 619   \brief Creates polyedres (private auxiliary method)
 620   \return 0 if operation has been completed successfully, 1 otherwise
 621 */
 622 int
 623 VISU_Vtk2MedConvertor
 624 ::CreatePolyedres(vtkDataSet* theInput,
 625                   MED::PMeshInfo theMeshInfo,
 626                   MED::PWrapper  theMed,
 627                   vtkIntArray* theCellsMapper,
 628                   MED::EEntiteMaillage theEntity,
 629                   TGeom2CellIds& theGeom2CellIdMap)
 630 {
 631   int theVTKGeom = VTK_CONVEX_POINT_SET;
 632   vtkIntArray* aCellIds = vtkIntArray::New();
 633   const MED::EConnectivite aConnType = MED::eNOD;
 634
 635   MED::TIntVector aConn;
 636   MED::TIntVector aFamilyNums;// -1
 637   MED::TIntVector aElemNums;
 638   MED::TIntVector aPolyedreInds;
 639   MED::TIntVector aPolyedreFaces;
 640
 641   aPolyedreInds.push_back(1); // reference on the first element in the connectivities
 642   aPolyedreFaces.push_back(1);
 643
 644   this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
 645   int nbElems = aCellIds->GetNumberOfTuples();
 646   if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
 647   if(nbElems>0){
 648     TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
 649     int* aPointer = aCellIds->GetPointer(0);
 650     for(int i=0;i<nbElems;i++,aPointer++){
 651       int aCellId = *aPointer;
 652       aCellIdsMapper.push_back(aCellId);
 653       vtkCell* aCell = theInput->GetCell(aCellId);
 654       int nbPointsInCell = aCell->GetNumberOfPoints();
 655       for(int j=0;j<nbPointsInCell;j++)
 656         aConn.push_back(aCell->GetPointId(j)+1);
 657       int aPrevPos = aPolyedreFaces.back();
 658       aPolyedreFaces.push_back(aPrevPos + nbPointsInCell);
 659       aPrevPos = aPolyedreInds.back();
 660       aPolyedreInds.push_back(aPrevPos + 1/*aNbFaces*/);
 661       aFamilyNums.push_back(-1);
 662
 663       if(theCellsMapper){
 664         if(theCellsMapper->GetNumberOfComponents()==2)
 665           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
 666         else if(theCellsMapper->GetNumberOfComponents()==1)
 667           aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
 668       }
 669     }
 670
 671
 672     MED::PPolyedreInfo aCellInfo = theMed->CrPolyedreInfo(theMeshInfo,
 673                                                           theEntity,
 674                                                           VTK2MED(theVTKGeom),
 675                                                           aPolyedreInds,
 676                                                           aPolyedreFaces,
 677                                                           aConn,
 678                                                           aConnType,
 679                                                           aFamilyNums,
 680                                                           aElemNums);
 681     theMed->SetPolyedreInfo(aCellInfo);
 682   }
 683
 684   aCellIds->Delete();
 685
 686   return 0;
 687 }
 688
 689 /*!
 690   \brief Converts geometry to med (private auxiliary method)
 691   \return 0 if operation has been completed successfully, 1 otherwise
 692 */
 693 int
 694 VISU_Vtk2MedConvertor
 695 ::Geometry2MED(vtkDataSet* aInput,
 696                MED::PWrapper myMed,
 697                MED::PMeshInfo aMeshInfo,
 698                TGeom2CellIds& outGeom2CellIdMap)
 699 {
 700   int aNbNodes = aInput->GetNumberOfPoints();
 701   int aMeshDimension = aMeshInfo->GetDim();
 702   // ----------------------- NODES -------------------------
 703   vtkIntArray* aPointsMapper;
 704   if(aInput->GetPointData())
 705     aPointsMapper = dynamic_cast<vtkIntArray*>(aInput->GetPointData()->GetArray(myPointDataFieldNameIDS.c_str()));
 706
 707   MED::TFloatVector aCoordinates(aNbNodes*aMeshDimension);
 708   MED::TIntVector   anElemNumsNodes; // takes from VISU_POINTS_MAPPER array
 709   MED::TIntVector  aFamilyNumsNodes;
 710   MED::TStringVector aCoordNamesNodes;
 711   MED::TStringVector aCoordUnitsNodes;
 712
 713   vtkFloatingPointType aPntCoord[3];
 714   if(aPointsMapper){
 715     int nbComp = aPointsMapper->GetNumberOfComponents();
 716     int *aPointsMapperPtr = aPointsMapper->GetPointer(0);
 717     for(int i=0;i<aNbNodes;i++){
 718       aInput->GetPoint(i,aPntCoord);
 719       aCoordinates[i*3]   = aPntCoord[0];
 720       aCoordinates[i*3+1] = aPntCoord[1];
 721       aCoordinates[i*3+2] = aPntCoord[2];
 722       if(nbComp == 2){
 723         anElemNumsNodes.push_back(*aPointsMapperPtr);
 724         aPointsMapperPtr++;aPointsMapperPtr++;
 725       }
 726       else if (nbComp == 1){
 727         anElemNumsNodes.push_back(*aPointsMapperPtr);
 728         aPointsMapperPtr++;
 729       }
 730       else{
 731         cout << "Error in file=|" << __FILE__<<"| line:[" << __LINE__ << "]" << endl;
 732         cout << "Code must be adapted for more than 2 components array |VISU_POINTS_MAPPER|" << endl;
 733         return 1;
 734       }
 735
 736     }
 737   } else {
 738     for(int i=0;i<aNbNodes;i++){
 739       aInput->GetPoint(i,aPntCoord);
 740       aCoordinates[i*3]   = aPntCoord[0];
 741       aCoordinates[i*3+1] = aPntCoord[1];
 742       aCoordinates[i*3+2] = aPntCoord[2];
 743     }
 744   }
 745
 746
 747   MED::PNodeInfo aNodeInfo = myMed->CrNodeInfo(aMeshInfo,
 748                                                aCoordinates,
 749                                                MED::eFULL_INTERLACE,
 750                                                MED::eCART,
 751                                                aCoordNamesNodes,
 752                                                aCoordUnitsNodes,
 753                                                aFamilyNumsNodes,
 754                                                anElemNumsNodes);
 755   myMed->SetNodeInfo(aNodeInfo);
 756
 757   vtkIntArray* aCellsMapper;
 758   if(vtkCellData* aCD = aInput->GetCellData())
 759     aCellsMapper = dynamic_cast<vtkIntArray*>(aCD->GetArray(myCellDataFieldNameIDS.c_str()));
 760
 761   if(MYDEBUG)
 762   {
 763     // debug info
 764     // print all cell types in the input
 765     vtkCellTypes* aCellTypes = vtkCellTypes::New();
 766     aInput->GetCellTypes(aCellTypes);
 767     cout << "Cell types in the input data:"<<endl;
 768     for(int aNbCellType = 0;aNbCellType<aCellTypes->GetNumberOfTypes();aNbCellType++)
 769       cout << (int)(aCellTypes->GetCellType(aNbCellType)) << endl;
 770     aCellTypes->Delete();
 771   }
 772
 773   //----------------------
 774   // Entity EDGES (eARETE)
 775   //----------------------
 776   vector<int> aNumberingConvertor;
 777   {
 778     // aVTKGeom->eSEG2
 779     MED::EEntiteMaillage aEntity = MED::eMAILLE;//eARETE;
 780     int aVTKGeom = 0;
 781     int nbPointsInGeom = 0;
 782
 783     // aNumberingConvertor NULL - OK
 784     aVTKGeom       = VTK_LINE;
 785     nbPointsInGeom = 2;
 786     CreateElements(aInput,
 787                    aMeshInfo,
 788                    myMed,
 789                    aCellsMapper,
 790                    aEntity,
 791                    aVTKGeom,
 792                    nbPointsInGeom,
 793                    aNumberingConvertor,
 794                    outGeom2CellIdMap);
 795
 796     // aNumberingConvertor NULL - OK
 797     // debug info: checked - OK
 798     aVTKGeom       = VTK_QUADRATIC_EDGE;
 799     nbPointsInGeom = 3;
 800     CreateElements(aInput,
 801                    aMeshInfo,
 802                    myMed,
 803                    aCellsMapper,
 804                    aEntity,
 805                    aVTKGeom,
 806                    nbPointsInGeom,
 807                    aNumberingConvertor,
 808                    outGeom2CellIdMap);
 809
 810   }
 811   //----------------------------
 812   // Entity FACES (eFACE)
 813   // eTRIA3,eQUAD4,eTRIA6,eQUAD8
 814   // ePOLYGONE
 815   //----------------------------
 816   {
 817
 818     MED::EEntiteMaillage aEntity = MED::eMAILLE;//MED::eFACE;
 819     int aVTKGeom = 0;
 820     int nbPointsInGeom = 0;
 821
 822     // debug info: checked - OK
 823     aVTKGeom       = VTK_TRIANGLE;
 824     aNumberingConvertor.clear();
 825     nbPointsInGeom = 3;
 826     CreateElements(aInput,
 827                    aMeshInfo,
 828                    myMed,
 829                    aCellsMapper,
 830                    aEntity,
 831                    aVTKGeom,
 832                    nbPointsInGeom,
 833                    aNumberingConvertor,
 834                    outGeom2CellIdMap);
 835
 836     // debug info: checked - OK
 837     aVTKGeom       = VTK_QUAD;
 838     nbPointsInGeom = 4;
 839     aNumberingConvertor.clear();
 840     CreateElements(aInput,
 841                    aMeshInfo,
 842                    myMed,
 843                    aCellsMapper,
 844                    aEntity,
 845                    aVTKGeom,
 846                    nbPointsInGeom,
 847                    aNumberingConvertor,
 848                    outGeom2CellIdMap);
 849
 850
 851     // debug info: checked - OK
 852     aVTKGeom       = VTK_QUADRATIC_TRIANGLE;
 853     nbPointsInGeom = 6;
 854     aNumberingConvertor.clear();
 855     CreateElements(aInput,
 856                    aMeshInfo,
 857                    myMed,
 858                    aCellsMapper,
 859                    aEntity,
 860                    aVTKGeom,
 861                    nbPointsInGeom,
 862                    aNumberingConvertor,
 863                    outGeom2CellIdMap);
 864
 865     aVTKGeom       = VTK_QUADRATIC_QUAD;
 866     nbPointsInGeom = 8;
 867     aNumberingConvertor.clear();
 868     // 0,3,2,1,7,6,5,4
 869     CreateElements(aInput,
 870                    aMeshInfo,
 871                    myMed,
 872                    aCellsMapper,
 873                    aEntity,
 874                    aVTKGeom,
 875                    nbPointsInGeom,
 876                    aNumberingConvertor,
 877                    outGeom2CellIdMap);
 878
 879     // debug info: checked - OK
 880     CreatePolygons(aInput,
 881                    aMeshInfo,
 882                    myMed,
 883                    aCellsMapper,
 884                    aEntity,
 885                    outGeom2CellIdMap);
 886
 887   }
 888   //----------------------------
 889   // Entity CELLS (eMAILLE)
 890   // eTETRA4,eHEXA8,ePENTA6,ePYRA5,
 891   // eTETRA10,eHEXA20,ePOLYEDRE
 892   //----------------------------
 893   {
 894
 895     MED::EEntiteMaillage aEntity = MED::eMAILLE;
 896     int aVTKGeom = 0;
 897     int nbPointsInGeom = 0;
 898
 899     aVTKGeom       = VTK_TETRA;
 900     nbPointsInGeom = 4;
 901     aNumberingConvertor.clear();
 902     aNumberingConvertor.push_back(0);
 903     aNumberingConvertor.push_back(2);
 904     aNumberingConvertor.push_back(1);
 905     aNumberingConvertor.push_back(3);
 906     CreateElements(aInput,
 907                    aMeshInfo,
 908                    myMed,
 909                    aCellsMapper,
 910                    aEntity,
 911                    aVTKGeom,
 912                    nbPointsInGeom,
 913                    aNumberingConvertor,
 914                    outGeom2CellIdMap);
 915
 916     aVTKGeom       = VTK_HEXAHEDRON;
 917     nbPointsInGeom = 8;
 918     aNumberingConvertor.clear();
 919     // 0,3,2,1,4,7,6,5
 920     CreateElements(aInput,
 921                    aMeshInfo,
 922                    myMed,
 923                    aCellsMapper,
 924                    aEntity,
 925                    aVTKGeom,
 926                    nbPointsInGeom,
 927                    aNumberingConvertor,
 928                    outGeom2CellIdMap);
 929
 930     aVTKGeom       = VTK_WEDGE;
 931     nbPointsInGeom = 6;
 932     aNumberingConvertor.clear();
 933     CreateElements(aInput,
 934                    aMeshInfo,
 935                    myMed,
 936                    aCellsMapper,
 937                    aEntity,
 938                    aVTKGeom,
 939                    nbPointsInGeom,
 940                    aNumberingConvertor,
 941                    outGeom2CellIdMap);
 942
 943     aVTKGeom       = VTK_PYRAMID;
 944     nbPointsInGeom = 5;
 945     aNumberingConvertor.clear();
 946     aNumberingConvertor.push_back(0);
 947     aNumberingConvertor.push_back(3);
 948     aNumberingConvertor.push_back(2);
 949     aNumberingConvertor.push_back(1);
 950     aNumberingConvertor.push_back(4);
 951     CreateElements(aInput,
 952                    aMeshInfo,
 953                    myMed,
 954                    aCellsMapper,
 955                    aEntity,
 956                    aVTKGeom,
 957                    nbPointsInGeom,
 958                    aNumberingConvertor,
 959                    outGeom2CellIdMap);
 960
 961     aVTKGeom       = VTK_QUADRATIC_TETRA;
 962     nbPointsInGeom = 10;
 963     aNumberingConvertor.clear();
 964     // 0,2,1,3,6,5,4,7,9,8
 965     CreateElements(aInput,
 966                    aMeshInfo,
 967                    myMed,
 968                    aCellsMapper,
 969                    aEntity,
 970                    aVTKGeom,
 971                    nbPointsInGeom,
 972                    aNumberingConvertor,
 973                    outGeom2CellIdMap);
 974
 975
 976     aVTKGeom       = VTK_QUADRATIC_HEXAHEDRON;
 977     nbPointsInGeom = 20;
 978     aNumberingConvertor.clear();
 979     // 0,3,2,1,4,7,6,5,11,10,9,8,15,14,13,12,16,19,18,17
 980     CreateElements(aInput,
 981                    aMeshInfo,
 982                    myMed,
 983                    aCellsMapper,
 984                    aEntity,
 985                    aVTKGeom,
 986                    nbPointsInGeom,
 987                    aNumberingConvertor,
 988                    outGeom2CellIdMap);
 989
 990     // debug info: checked OK
 991     CreatePolyedres(aInput,
 992                     aMeshInfo,
 993                     myMed,
 994                     aCellsMapper,
 995                     aEntity,
 996                     outGeom2CellIdMap);
 997   }
 998
 999
1000   return 0;
1001 }
1002
1003 /*!
1004   \brief Converts data to med (private auxiliary method)
1005   \return 0 if operation has been completed successfully, 1 otherwise
1006 */
1007 int
1008 VISU_Vtk2MedConvertor
1009 ::Data2MED(std::vector<vtkDataSet*> theListForAdd,
1010            MED::PWrapper myMed,
1011            MED::PMeshInfo theMeshInfo,
1012            TGeom2CellIds& theGeom2CellIdMap)
1013 {
1014   typedef std::vector<vtkPointData*> TPDVec;
1015   typedef std::vector<vtkCellData*> TCDVec;
1016   typedef std::map<std::string,TPDVec> TNameToPointData;
1017   typedef std::map<std::string,TCDVec> TNameToCellData;
1018
1019   TNameToPointData aName2PointData;
1020   TNameToCellData  aName2CellData;
1021
1022   MED::TErr* theErrCode = new MED::TErr();
1023   MED::EGeometrieElement geomType;
1024   MED::EEntiteMaillage   entity;
1025
1026   // prepare data to create time stamps
1027   const MED::TInt aNbGauss = 1;
1028   MED::TProfileInfo::TInfo aTInfo("",0);
1029   MED::PProfileInfo aPProfileInfo = myMed->CrProfileInfo( aTInfo );
1030
1031   MED::TGeom2Size    aGeom2Size;
1032   MED::TGeom2NbGauss aGeom2NbGauss;
1033   MED::TGeom2Profile aTGeom2Profile;
1034
1035   int nbPointsInFirstData = 0;
1036   if(theListForAdd.size()>0)
1037     nbPointsInFirstData = (theListForAdd[0])->GetNumberOfPoints();
1038
1039   int nbClasses = theListForAdd.size();
1040
1041   for(int iClass=0;iClass<nbClasses;iClass++){
1042     vtkDataSet* aClassPtr = theListForAdd[iClass];
1043     if(aClassPtr->GetNumberOfPoints() != nbPointsInFirstData){
1044       cout << "Warning in PointData: Some vtk file consist of number of points( " <<aClassPtr->GetNumberOfPoints()
1045            << ") not equal number of points in first file("<<nbPointsInFirstData<<")"<<endl;
1046       cout << "This data will be lost." << endl;
1047       continue;
1048     }
1049
1050     if(vtkPointData* aPD = aClassPtr->GetPointData()){
1051       int nbArrays = aPD->GetNumberOfArrays();
1052
1053       for(int aArrNum=0;aArrNum<nbArrays;aArrNum++)
1054       {
1055         vtkDataArray* aArr = aPD->GetArray(aArrNum);
1056         std::string aName = aArr->GetName();
1057         std::set<std::string>::const_iterator aIgnoreIter = myIgnoringFieldList.find(aName);
1058         if(aIgnoreIter!=myIgnoringFieldList.end())
1059           continue;
1060         (aName2PointData[aName]).push_back(aPD);
1061       }
1062     }
1063
1064     if(vtkCellData* aCD = aClassPtr->GetCellData()){
1065       int nbArrays = aCD->GetNumberOfArrays();
1066
1067       for(int aArrNum=0;aArrNum<nbArrays;aArrNum++)
1068       {
1069         vtkDataArray* aArr = aCD->GetArray(aArrNum);
1070         std::string aName = aArr->GetName();
1071         std::set<std::string>::const_iterator aIgnoreIter = myIgnoringFieldList.find(aName);
1072         if(aIgnoreIter!=myIgnoringFieldList.end())
1073           continue;
1074         (aName2CellData[aName]).push_back(aCD);
1075       }
1076     }
1077   }
1078
1079
1080   // PointData
1081   int aLastField = 0;
1082   {
1083     geomType = MED::ePOINT1;
1084     entity   = MED::eNOEUD;
1085     aGeom2Size[geomType]    = nbPointsInFirstData;
1086     aGeom2NbGauss[geomType] = aNbGauss;
1087     aTGeom2Profile[geomType]= aPProfileInfo;
1088
1089     TNameToPointData::const_iterator aIter = aName2PointData.begin();
1090     for(int iField=1;aIter!=aName2PointData.end();aIter++,iField++){
1091       std::string aFieldName = aIter->first;
1092       TPDVec aPD2Vec = (aIter->second);
1093       int nbComp = 0;
1094       if(aPD2Vec.size() >0 ){
1095         if(vtkPointData* aPD0  = aPD2Vec[0]){
1096           if(vtkDataArray* aArr0 = aPD0->GetArray(aFieldName.c_str()))
1097             nbComp = aArr0->GetNumberOfComponents();
1098         }
1099       }
1100
1101       MED::PFieldInfo aFieldInfo = myMed->CrFieldInfo(theMeshInfo,
1102                                                       nbComp);
1103
1104       string aFieldName_PD = "Point " + aFieldName;
1105       aFieldInfo->SetName(aFieldName_PD.c_str());
1106
1107       myMed->SetFieldInfo(aFieldInfo);
1108
1109       TPDVec::const_iterator aPDIter = aPD2Vec.begin();
1110       for(int iTStamp=0;aPDIter!=aPD2Vec.end();aPDIter++,iTStamp++){
1111         vtkPointData* aPD = *aPDIter;
1112         vtkDataArray* aArr = aPD->GetArray(aFieldName.c_str());
1113         MED::TFloat aTS = iTStamp < (int)myTStamps.size() ? myTStamps[ iTStamp ] : iTStamp;
1114         MED::PTimeStampInfo aTempTimeStampInfo = myMed->CrTimeStampInfo (aFieldInfo,
1115                                                                          entity,
1116                                                                          aGeom2Size,
1117                                                                          aGeom2NbGauss,
1118                                                                          iTStamp,
1119                                                                          iTStamp,
1120                                                                          aTS);
1121
1122         MED::PTimeStampVal aTempTimeStampVal = myMed->CrTimeStampVal (aTempTimeStampInfo,
1123                                                                       aTGeom2Profile);
1124
1125         MED::TMeshValue& aTMeshValue = aTempTimeStampVal->GetMeshValue(geomType);
1126
1127         MED::TValue& aValue = aTMeshValue.myValue; // float
1128         int nbValues = aValue.size();
1129         for(int i=0;i<nbValues;i++){
1130           aValue[i] = *(float*)(aArr->GetVoidPointer(i));
1131         }
1132
1133         myMed->SetTimeStamp( aTempTimeStampVal, theErrCode);
1134         if(*theErrCode==0){
1135           cout << "Error in "<<__FILE__<<"["<<__LINE__<<"] in method SetTimeStamp(...)"<<endl;
1136           return 1;
1137         }
1138       }
1139       aLastField = iField;
1140     }
1141   }
1142   // CellData
1143   {
1144     MED::TEntityInfo aEntityInfo = myMed->GetEntityInfo(theMeshInfo);
1145     aGeom2Size.clear();
1146     aGeom2NbGauss.clear();
1147     aTGeom2Profile.clear();
1148     aGeom2Size = aEntityInfo[ MED::eMAILLE ];
1149     entity   = MED::eMAILLE;
1150     MED::TGeom2Size::iterator geom_nb;
1151     for ( geom_nb = aGeom2Size.begin(); geom_nb != aGeom2Size.end(); ++geom_nb ) { // loop on geometric types of cell
1152       aGeom2NbGauss[ geom_nb->first ] = aNbGauss;
1153       aTGeom2Profile[ geom_nb->first ] = aPProfileInfo;
1154     }
1155
1156
1157     TNameToCellData::const_iterator aIter = aName2CellData.begin();
1158     for(int iField=1;aIter!=aName2CellData.end();aIter++,iField++){
1159       std::string aFieldName = aIter->first;
1160       TCDVec aCD2Vec = (aIter->second);
1161       int nbComp = 0;
1162       if(aCD2Vec.size() >0 ){
1163         if(vtkCellData* aCD0  = aCD2Vec[0]){
1164           if(vtkDataArray* aArr0 = aCD0->GetArray(aFieldName.c_str()))
1165             nbComp = aArr0->GetNumberOfComponents();
1166         }
1167       }
1168
1169       MED::PFieldInfo aFieldInfo = myMed->CrFieldInfo(theMeshInfo,
1170                                                       nbComp);
1171
1172       string aFieldName_CD = "Cell " + aFieldName;
1173       aFieldInfo->SetName(aFieldName_CD.c_str());
1174
1175       myMed->SetFieldInfo(aFieldInfo);
1176
1177       TCDVec::const_iterator aCDIter = aCD2Vec.begin();
1178       for(int iTStamp=0;aCDIter!=aCD2Vec.end();aCDIter++,iTStamp++){
1179         vtkCellData* aCD = *aCDIter;
1180         vtkDataArray* aArr = aCD->GetArray(aFieldName.c_str());
1181         MED::TFloat aTS = iTStamp < (int)myTStamps.size() ? myTStamps[ iTStamp ] : iTStamp;
1182         MED::PTimeStampInfo aTempTimeStampInfo = myMed->CrTimeStampInfo (aFieldInfo,
1183                                                                          entity,
1184                                                                          aGeom2Size,
1185                                                                          aGeom2NbGauss,
1186                                                                          iTStamp,
1187                                                                          iTStamp,
1188                                                                          aTS);
1189
1190         MED::PTimeStampVal aTempTimeStampVal = myMed->CrTimeStampVal (aTempTimeStampInfo,
1191                                                                       aTGeom2Profile);
1192
1193         for ( geom_nb = aGeom2Size.begin(); geom_nb != aGeom2Size.end(); ++geom_nb ) { // loop on geometric types of cell
1194           geomType = geom_nb->first;
1195           TCellIds& aCellIds = theGeom2CellIdMap[geomType];
1196
1197           MED::TMeshValue& aTMeshValue = aTempTimeStampVal->GetMeshValue(geomType);
1198
1199           MED::TValue& aValue = aTMeshValue.myValue; // float
1200           int nbValues = aValue.size();
1201           int nbCellIds = aCellIds.size();
1202           if(nbValues!=nbCellIds*nbComp){
1203             cout << "Warning in "<<__FILE__<<"["<<__LINE__<<"] the data for geometry:"<<geomType<<" will be ignored"<<endl;
1204             continue;
1205           }
1206           if(MYDEBUG_VALUES) cout << "Geom["<<geomType<<"]"<<endl;
1207           for(int i=0;i<nbCellIds;i++){
1208             if(MYDEBUG_VALUES) cout << "\t|";
1209             for(int iComp=0;iComp<nbComp;iComp++){
1210               aValue[nbComp*i+iComp] = *(float*)(aArr->GetVoidPointer(nbComp*aCellIds[i]+iComp));
1211               if(MYDEBUG_VALUES) cout << aValue[nbComp*i+iComp] << "   ";
1212             }
1213             if(MYDEBUG_VALUES) cout << "|" << endl;
1214           }
1215
1216           myMed->SetTimeStamp( aTempTimeStampVal, theErrCode);
1217           if(*theErrCode==0){
1218             cout << "Error in "<<__FILE__<<"["<<__LINE__<<"] in method SetTimeStamp(...)"<<endl;
1219             return 1;
1220           }
1221         }
1222       }
1223     }
1224   }
1225
1226
1227   return 0;
1228 }
1229
1230 /*!
1231   \brief Writes data to MED file
1232   \return 0 if operation has been completed successfully, 1 otherwise
1233 */
1234 int
1235 VISU_Vtk2MedConvertor
1236 ::Execute()
1237 {
1238   int aStatus = 1;
1239   if (myFirstVTKFileName.size() == 0 ||
1240       myMEDFileName.size() == 0)
1241   {
1242     cout << "Error! Bad input file names output med file name or vtk file name."<<endl;
1243     cout << "Exit."<<endl;
1244     return 1;
1245   }
1246
1247   MED::PWrapper myMed;
1248   MED::PMeshInfo aMeshInfo;
1249   int aMeshDimension = 3;
1250   int aSpaceDimension = 3;
1251   myMed = CrWrapper(myMEDFileName.c_str(),myVersion);
1252   aMeshInfo = myMed->CrMeshInfo(aMeshDimension,aSpaceDimension,myMeshName.c_str());
1253   myMed->SetMeshInfo(aMeshInfo);
1254
1255   {
1256     typedef vtkDataSetReader TReader;
1257     TReader* aReader = TReader::New();
1258     aReader->SetFileName(myFirstVTKFileName.c_str());
1259     aReader->Update();
1260     TGeom2CellIds myGeom2CellIds;
1261
1262     typedef std::vector<vtkDataSet*> TListUG;
1263     TListUG aList;
1264
1265     if(aReader->IsFilePolyData())
1266     {
1267       if(MYDEBUG) cout << "PolyData" << endl;
1268       typedef vtkPolyData TCommonType;
1269       TCommonType* aInput = aReader->GetPolyDataOutput();
1270
1271       aStatus = Geometry2MED(aInput,
1272                              myMed,
1273                              aMeshInfo,
1274                              myGeom2CellIds);
1275
1276
1277       TCommonType* aUG1 = TCommonType::New();
1278       aUG1->ShallowCopy(aInput);
1279       vtkDataSet* aTmp1 = dynamic_cast<vtkDataSet*>(aUG1);
1280       aList.push_back(aTmp1);
1281
1282       TVectorString::iterator aFilesIter = myDataVTKFileNames.begin();
1283       for(;aFilesIter!=myDataVTKFileNames.end();aFilesIter++){
1284         aReader->SetFileName((*aFilesIter).c_str());
1285         aReader->Update();
1286         TCommonType* aUG2 = TCommonType::New();
1287         aUG2->ShallowCopy(aReader->GetPolyDataOutput());
1288         vtkDataSet* aTmp2 = dynamic_cast<vtkDataSet*>(aUG2);
1289         aList.push_back(aTmp2);
1290       }
1291     } else if (aReader->IsFileUnstructuredGrid()){
1292       if (MYDEBUG) cout << "UnstructuredGrid" << endl;
1293       typedef vtkUnstructuredGrid TCommonType;
1294       TCommonType* aInput = aReader->GetUnstructuredGridOutput();
1295
1296       aStatus = Geometry2MED(aInput,
1297                              myMed,
1298                              aMeshInfo,
1299                              myGeom2CellIds);
1300
1301
1302       TCommonType* aUG1 = TCommonType::New();
1303       aUG1->ShallowCopy(aInput);
1304       vtkDataSet* aTmp1 = dynamic_cast<vtkDataSet*>(aUG1);
1305       aList.push_back(aTmp1);
1306
1307       TVectorString::iterator aFilesIter = myDataVTKFileNames.begin();
1308       for(;aFilesIter!=myDataVTKFileNames.end();aFilesIter++){
1309         aReader->SetFileName((*aFilesIter).c_str());
1310         aReader->Update();
1311         TCommonType* aUG2 = TCommonType::New();
1312         aUG2->ShallowCopy(aReader->GetUnstructuredGridOutput());
1313         vtkDataSet* aTmp2 = dynamic_cast<vtkDataSet*>(aUG2);
1314         aList.push_back(aTmp2);
1315       }
1316     }
1317
1318     Data2MED(aList,
1319              myMed,
1320              aMeshInfo,
1321              myGeom2CellIds);
1322
1323
1324     // clear aList by removing of unstructured grids
1325     TListUG::iterator aIter = aList.begin();
1326     for(;aIter!=aList.end();aIter++)
1327       (*aIter)->Delete();
1328
1329     aReader->Delete();
1330   }
1331
1332   return aStatus;
1333 }