--- /dev/null
+// File : VTK2MED_Convertor.cxx
+// Created : Fri May 4 15:10:49 2007
+// Author : Eugeny NIKOLAEV (enk)
+// Copyright : Open CASCADE
+
+#include "VTK2MED_Convertor.hxx"
+
+// QT includes
+#include <qdir.h>
+#include <qfileinfo.h>
+#include <qstringlist.h>
+
+// VTK includes
+#include <vtkSmartPointer.h>
+#include <vtkDataReader.h>
+
+#include <vtkStructuredPointsReader.h>
+#include <vtkStructuredGridReader.h>
+#include <vtkRectilinearGridReader.h>
+#include <vtkUnstructuredGridReader.h>
+#include <vtkPolyDataReader.h>
+#include <vtkDataSetReader.h>
+
+#include <vtkPointSet.h>
+#include <vtkDataSet.h>
+#include <vtkPolyData.h>
+#include <vtkDataObject.h>
+#include <vtkCellTypes.h>
+#include <vtkCellType.h>
+#include <vtkCell.h>
+#include <vtkCellData.h>
+#include <vtkPointData.h>
+
+#include <vtkUnstructuredGrid.h>
+#include <vtkFloatArray.h>
+
+// MED Warpper includes
+#include "MED_Factory.hxx"
+
+// STL includes
+#include <string>
+#include <vector>
+#include <iostream>
+#include <sstream>
+#include <set>
+#include <map>
+//#include <pair.h>
+
+#ifdef _DEBUG_
+static int MYDEBUG = 0;
+static int MYDEBUG_VALUES = 0;
+#else
+static int MYDEBUG = 0;
+static int MYDEBUG_VALUES = 0;
+#endif
+
+using namespace std;
+
+/*
+# = dynamic_cast</*
+*>( );define VTK_EMPTY_CELL 0
+#define VTK_VERTEX 1
+#define VTK_POLY_VERTEX 2
+#define VTK_LINE 3
+#define VTK_POLY_LINE 4
+#define VTK_TRIANGLE 5
+#define VTK_TRIANGLE_STRIP 6
+#define VTK_POLYGON 7
+#define VTK_PIXEL 8
+#define VTK_QUAD 9
+#define VTK_TETRA 10
+#define VTK_VOXEL 11
+#define VTK_HEXAHEDRON 12
+#define VTK_WEDGE 13
+#define VTK_PYRAMID 14
+#define VTK_QUADRATIC_EDGE 21
+#define VTK_QUADRATIC_TRiTSIANGLE 22
+#define VTK_QUADRATIC_QUAD 23
+#define VTK_QUADRATIC_TETRA 24
+#define VTK_QUADRATIC_HEXAHEDRON 25
+#define VTK_CONVEX_POINT_SET 41
+#define VTK_PARAMETRIC_CURVE 51
+#define VTK_PARAMETRIC_SURFACE 52
+#define VTK_PARAMETRIC_TRI_SURFACE 53
+#define VTK_PARAMETRIC_QUAD_SURFACE 54
+#define VTK_PARAMETRIC_TETRA_REGION 55
+#define VTK_PARAMETRIC_HEX_REGION 56
+*/
+
+static MED::EGeometrieElement VTK2MED( const int theGeom )
+{
+ // Ignoring vtk types:
+ // VTK_PIXEL,
+ // VTK_VERTEX,
+ // VTK_POLY_VERTEX
+ // VTK_VOXEL
+ // VTK_POLY_LINE
+ // VTK_TRIANGLE_STRIP
+ // VTK_PARAMETRIC_CURVE
+ // VTK_PARAMETRIC_SURFACE
+ // VTK_PARAMETRIC_TRI_SURFACE
+ // VTK_PARAMETRIC_QUAD_SURFACE
+ // VTK_PARAMETRIC_TETRA_REGION
+ // VTK_PARAMETRIC_HEX_REGION
+
+ MED::EGeometrieElement aEmptyGeom = MED::EGeometrieElement(-1);
+ switch(theGeom){
+ case VTK_LINE: return MED::eSEG2;
+ case VTK_TRIANGLE: return MED::eTRIA3;
+ case VTK_POLYGON: return MED::ePOLYGONE;
+ case VTK_QUAD: return MED::eQUAD4;
+ case VTK_TETRA: return MED::eTETRA4;
+ case VTK_HEXAHEDRON: return MED::eHEXA8;
+ case VTK_WEDGE: return MED::ePENTA6;
+ case VTK_PYRAMID: return MED::ePYRA5;
+ // QUADRATIC elements
+ case VTK_QUADRATIC_EDGE: return MED::eSEG3;
+ case VTK_QUADRATIC_TRIANGLE: return MED::eTRIA6;
+ case VTK_QUADRATIC_QUAD: return MED::eQUAD8;
+ case VTK_QUADRATIC_TETRA: return MED::eTETRA10;
+ case VTK_QUADRATIC_HEXAHEDRON:return MED::eHEXA20;
+ case VTK_CONVEX_POINT_SET: return MED::ePOLYEDRE;
+ }
+
+ return aEmptyGeom;
+}
+
+/*!
+ \class VTK2MED_Convertor
+ \brief The general main of the class VTK2MED_Convertor is converting from
+ one or several VTK files to the one MED file...
+
+ The VTK2MED_Convertor interface allows us to create the MED file according
+ to VTK files in next ways:
+ - Extract geometry and fields from one VTK file.
+ - Extract geometry and fields from first VTK file and fields from others
+ VTK files (geometry ignoring). Also the fields which have same names join into
+ corresponding fields with different time stamp
+*/
+
+/*!
+ \brief Constructor
+ - Sets default output mesh name
+ - Sets default version of output MED file.
+ - Sets default ignoring fields list
+ - Sets default points and cells ids mapping field names.
+
+*/
+VTK2MED_Convertor
+::VTK2MED_Convertor():
+ myVersion(MED::eV2_2),
+ myMeshName("vtk2med")
+{
+ myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
+ myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
+ myIgnoringFieldList.insert("VISU_FIELD");
+ setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
+ setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
+}
+
+/*!
+ \brief Constructor
+ - Sets default output mesh name
+ - Sets default version of output MED file.
+ - Sets default ignoring fields list
+ - Sets default points and cells ids mapping field names.
+
+ \param theMEDFileName output med file name
+ \param theFirstVTKFileName first vtk file name
+*/
+VTK2MED_Convertor
+::VTK2MED_Convertor(const string theMEDFileName,
+ const string theFirstVTKFileName):
+ myVersion(MED::eV2_2),
+ myMeshName("vtk2med")
+{
+ myMEDFileName = theMEDFileName;
+ myFirstVTKFileName = theFirstVTKFileName;
+ myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
+ myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
+ myIgnoringFieldList.insert("VISU_FIELD");
+ setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
+ setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
+}
+
+/*!
+ \brief Constructor
+ - Sets default output mesh name
+ - Sets default version of output MED file.
+ - Sets default ignoring fields list
+ - Sets default points and cells ids mapping field names.
+
+ \param theMEDFileName output med file name
+ \param theFirstVTKFileName first vtk file name
+ \param theDataVTKFileNames of vtk file names, which will be using as values on points and cells
+*/
+VTK2MED_Convertor
+::VTK2MED_Convertor(const string theMEDFileName,
+ const string theFirstVTKFileName,
+ const TVectorString theDataVTKFileNames):
+ myVersion(MED::eV2_2),
+ myMeshName("vtk2med")
+{
+ myMEDFileName = theMEDFileName;
+ myFirstVTKFileName = theFirstVTKFileName;
+ myDataVTKFileNames = theDataVTKFileNames;
+ myIgnoringFieldList.insert("VISU_POINTS_MAPPER");
+ myIgnoringFieldList.insert("VISU_CELLS_MAPPER");
+ myIgnoringFieldList.insert("VISU_FIELD");
+ setMeshName("vtk2med");
+ setCellDataFieldNameIDS("VISU_CELLS_MAPPER");
+ setPointDataFieldNameIDS("VISU_POINTS_MAPPER");
+}
+
+/*!
+ \brief Adds field names, which used as specific fields with ids or elements
+ (or something else). (Default: \93VISU_CELLS_MAPPER\94,\94VISU_POINTS_MAPPER\94,\94VISU_FILED\94)
+ \param theFieldName field name
+ \sa eraseFromIgnoringFieldList()
+*/
+void
+VTK2MED_Convertor
+::addToIgnoringFieldList( const string& theFieldName )
+{
+ myIgnoringFieldList.insert(theFieldName);
+}
+
+/*!
+ \brief Sets the output MED file name
+ \param theFileName file name
+ \sa getMEDFileName()
+*/
+void
+VTK2MED_Convertor
+::setMEDFileName( const string theFileName )
+{
+ myMEDFileName=theFileName;
+};
+
+/*!
+ \brief Gets the output MED file name
+ \return output MED file name
+ \sa setMEDFileName()
+*/
+string
+VTK2MED_Convertor
+::getMEDFileName() const
+{
+ return myMEDFileName;
+}
+
+/*!
+ \brief Sets the first input vtk file name
+ \param theFileName file name
+ \sa getFirstVTKFileName()
+*/
+void
+VTK2MED_Convertor
+::setFirstVTKFileName( const string theFileName )
+{
+ myFirstVTKFileName=theFileName;
+};
+
+/*!
+ \brief Fets the first input vtk file name
+ \return first input vtk file name
+ \sa setFirstVTKFileName()
+*/
+string
+VTK2MED_Convertor
+::getFirstVTKFileName() const
+{
+ return myFirstVTKFileName;
+}
+
+/*!
+ \brief Sets list of vtk file names, which will be using as values on points and cells
+ \param theFileNames list of vtk file names
+ \sa getDataVTKFileNames()
+*/
+void
+VTK2MED_Convertor
+::setDataVTKFileNames( const TVectorString theFileNames )
+{
+ myDataVTKFileNames = theFileNames;
+};
+
+/*!
+ \brief Gets list of vtk file names, which will be using as values on points and cells
+ \param theFileNames out list of vtk file names
+ \sa setDataVTKFileNames()
+*/
+void
+VTK2MED_Convertor
+::getDataVTKFileNames( TVectorString& theDataVTKFileNames ) const
+{
+ theDataVTKFileNames = myDataVTKFileNames;
+};
+
+/*!
+ \brief Sets version of the output MED file MED::V2_2(is default) or MED::V2_1
+ \param theVersion version of the output MED file
+*/
+void
+VTK2MED_Convertor
+::setVersion( const MED::EVersion theVersion )
+{
+ myVersion = theVersion;
+}
+
+/*!
+ \brief Gets version of the output MED file MED::V2_2(is default) or MED::V2_1
+ \return version of the output MED file
+*/
+MED::EVersion
+VTK2MED_Convertor
+::getVersion() const
+{
+ return myVersion;
+}
+
+/*!
+ \brief Sets output mesh name. (\93vtk2med\94 - default)
+ \param theMeshName mesh name
+ \sa getMeshName()
+*/
+void
+VTK2MED_Convertor
+::setMeshName( const string theMeshName )
+{
+ myMeshName = theMeshName;
+}
+
+/*!
+ \brief Gets output mesh name. (\93vtk2med\94 - default)
+ \return mesh name
+ \sa setMeshName()
+*/
+string
+VTK2MED_Convertor
+::getMeshName() const
+{
+ return myMeshName;
+}
+
+/*!
+ \brief Sets field name with cell ids (Default - VISU_CELLS_MAPPER)
+ \param theFieldName field name
+ \sa getCellDataFieldNameIDS
+*/
+void
+VTK2MED_Convertor
+::setCellDataFieldNameIDS( const string& theFieldName )
+{
+ myCellDataFieldNameIDS = theFieldName;
+}
+
+/*!
+ \brief Gets field name with cell ids (Default - VISU_CELLS_MAPPER)
+ \return field name
+ \sa setCellDataFieldNameIDS()
+*/
+const string&
+VTK2MED_Convertor
+::getCellDataFieldNameIDS() const
+{
+ return myCellDataFieldNameIDS;
+}
+
+/*!
+ \brief Erases field names which used as specific fields with ids or elements
+ (or something else)
+ \param theFieldName field name
+ \sa addToIgnoringFieldList()
+*/
+void
+VTK2MED_Convertor
+::eraseFromIgnoringFieldList(const string& theFieldName)
+{
+ myIgnoringFieldList.erase(theFieldName);
+}
+
+/*!
+ \brief Gets list of field names which used as specific fields with ids or elements
+ \return list of field names
+*/
+const std::set<std::string>&
+VTK2MED_Convertor
+::getIgnoringFieldList() const
+{
+ return myIgnoringFieldList;
+}
+
+/*!
+ \brief Sets field name with point ids
+ \param theFieldName field name
+ \sa getPointDataFieldNameIDS()
+*/
+void
+VTK2MED_Convertor
+::setPointDataFieldNameIDS( const string& theFieldName )
+{
+ myPointDataFieldNameIDS = theFieldName;
+}
+
+/*!
+ \brief Gets field name with point ids
+ \return field name
+ \sa setPointDataFieldNameIDS()
+*/
+const string&
+VTK2MED_Convertor
+::getPointDataFieldNameIDS() const
+{
+ return myPointDataFieldNameIDS;
+}
+
+/*!
+ \brief Sets values of time stamps If this array is not specified values of time
+ stamps are generated automatically ( 0, 1, 2 ... )
+ \param theTStamps vector of time stamps
+ \sa getTimeStamps()
+*/
+void
+VTK2MED_Convertor
+::setTimeStamps( const TVectorDouble& theTStamps )
+{
+ myTStamps = theTStamps;
+}
+
+/*!
+ \brief Gets values of time stamps If this array is not specified values of time
+ stamps are generated automatically ( 0, 1, 2 ... )
+ \param theTStamps out vector of time stamps
+ \sa setTimeStamps()
+*/
+void
+VTK2MED_Convertor
+::getTimeStamps( TVectorDouble& theTStamps ) const
+{
+ theTStamps = myTStamps;
+}
+
+/*!
+ \brief Retrieves identifiers of cells from input data set corresponding to given type
+ \param theInput input data set
+ \param type type
+ \param array out array
+*/
+void
+VTK2MED_Convertor
+::GetIdsOfCellsOfType(vtkDataSet* theInput,
+ const int type,
+ vtkIntArray *array)
+{
+ for (int cellId = 0; cellId < theInput->GetNumberOfCells(); cellId++)
+ if (theInput->GetCellType(cellId) == type)
+ array->InsertNextValue(cellId);
+}
+
+/*!
+ \brief Creates elements (private auxiliary method)
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::CreateElements(vtkDataSet* theInput,
+ MED::PMeshInfo theMeshInfo,
+ MED::PWrapper theMed,
+ vtkIntArray* theCellsMapper,
+ MED::EEntiteMaillage theEntity,
+ int theVTKGeom,
+ int nbPointsInGeom,
+ std::vector<int>& theNumberingConvertor,
+ TGeom2CellIds& theGeom2CellIdMap)
+{
+ bool aIdsConvert = (theNumberingConvertor.size() > 0);
+ vtkIntArray* aCellIds = vtkIntArray::New();
+ const MED::EConnectivite aConnType = MED::eNOD;
+
+ MED::TIntVector aConn;
+ MED::TIntVector aFamilyNums;// -1
+ MED::TIntVector aElemNums;
+
+ this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
+ int nbElems = aCellIds->GetNumberOfTuples();
+ if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
+ aConn.reserve(nbElems*nbPointsInGeom);
+ if(nbElems>0){
+ TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
+ int* aPointer = aCellIds->GetPointer(0);
+ for(int i=0;i<nbElems;i++,aPointer++){
+ int aCellId = *aPointer;
+ aCellIdsMapper.push_back(aCellId);
+ vtkCell* aCell = theInput->GetCell(aCellId);
+ int nbPointsInCell = aCell->GetNumberOfPoints();
+ if(nbPointsInCell!=nbPointsInGeom){
+ cout << "Error in file=|" << __FILE__<<"| line:[" << __LINE__ << "]" << endl;
+ cout << "Must be "<<nbPointsInGeom<<" nodes in VTK Geometry:"<<theVTKGeom<<" element" << endl;
+ aCellIds->Delete();
+ return 1; // exit
+ }
+ aFamilyNums.push_back(-1);
+ for(int j=0;j<nbPointsInCell;j++){
+ if (aIdsConvert)
+ aConn.push_back(aCell->GetPointId(theNumberingConvertor[j])+1);
+ else
+ aConn.push_back(aCell->GetPointId(j)+1);
+ }
+ if(theCellsMapper){
+ if(theCellsMapper->GetNumberOfComponents()==2)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
+ else if(theCellsMapper->GetNumberOfComponents()==1)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
+ }
+ }
+
+
+ MED::PCellInfo aCellInfo = theMed->CrCellInfo(theMeshInfo,
+ theEntity,
+ VTK2MED(theVTKGeom),
+ aConn,
+ aConnType,
+ aFamilyNums,
+ aElemNums);
+ theMed->SetCellInfo(aCellInfo);
+ }
+
+ aCellIds->Delete();
+
+ return 0;
+}
+
+/*!
+ \brief Creates polygons (private auxiliary method)
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::CreatePolygons(vtkDataSet* theInput,
+ MED::PMeshInfo theMeshInfo,
+ MED::PWrapper theMed,
+ vtkIntArray* theCellsMapper,
+ MED::EEntiteMaillage theEntity,
+ TGeom2CellIds& theGeom2CellIdMap)
+{
+ int theVTKGeom = VTK_POLYGON;
+ vtkIntArray* aCellIds = vtkIntArray::New();
+ const MED::EConnectivite aConnType = MED::eNOD;
+
+ MED::TIntVector aConn;
+ MED::TIntVector aFamilyNums;// -1
+ MED::TIntVector aElemNums;
+ MED::TIntVector aPolygoneInds;
+ aPolygoneInds.push_back(1); // reference on the first element in the connectivities
+
+ this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
+ int nbElems = aCellIds->GetNumberOfTuples();
+ if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
+ if(nbElems>0){
+ TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
+ int* aPointer = aCellIds->GetPointer(0);
+ for(int i=0;i<nbElems;i++,aPointer++){
+ int aCellId = *aPointer;
+ aCellIdsMapper.push_back(aCellId);
+ vtkCell* aCell = theInput->GetCell(aCellId);
+ int nbPointsInCell = aCell->GetNumberOfPoints();
+ aFamilyNums.push_back(-1);
+ int aPrevPos = aPolygoneInds.back();
+ aPolygoneInds.push_back(aPrevPos+nbPointsInCell);
+ for(int j=0;j<nbPointsInCell;j++)
+ aConn.push_back(aCell->GetPointId(j)+1);
+ if(theCellsMapper){
+ if(theCellsMapper->GetNumberOfComponents()==2)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
+ else if(theCellsMapper->GetNumberOfComponents()==1)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
+ }
+ }
+
+
+ MED::PPolygoneInfo aCellInfo = theMed->CrPolygoneInfo(theMeshInfo,
+ theEntity,
+ VTK2MED(theVTKGeom),
+ aPolygoneInds,
+ aConn,
+ aConnType,
+ aFamilyNums,
+ aElemNums);
+ theMed->SetPolygoneInfo(aCellInfo);
+ }
+
+ aCellIds->Delete();
+
+ return 0;
+}
+
+/*!
+ \brief Creates polyedres (private auxiliary method)
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::CreatePolyedres(vtkDataSet* theInput,
+ MED::PMeshInfo theMeshInfo,
+ MED::PWrapper theMed,
+ vtkIntArray* theCellsMapper,
+ MED::EEntiteMaillage theEntity,
+ TGeom2CellIds& theGeom2CellIdMap)
+{
+ int theVTKGeom = VTK_CONVEX_POINT_SET;
+ vtkIntArray* aCellIds = vtkIntArray::New();
+ const MED::EConnectivite aConnType = MED::eNOD;
+
+ MED::TIntVector aConn;
+ MED::TIntVector aFamilyNums;// -1
+ MED::TIntVector aElemNums;
+ MED::TIntVector aPolyedreInds;
+ MED::TIntVector aPolyedreFaces;
+
+ aPolyedreInds.push_back(1); // reference on the first element in the connectivities
+ aPolyedreFaces.push_back(1);
+
+ this->GetIdsOfCellsOfType(theInput,theVTKGeom,aCellIds);
+ int nbElems = aCellIds->GetNumberOfTuples();
+ if(MYDEBUG) cout << "\tnbElems in geom:"<<VTK2MED(theVTKGeom)<<" ="<<nbElems<<endl;
+ if(nbElems>0){
+ TCellIds& aCellIdsMapper = theGeom2CellIdMap[VTK2MED(theVTKGeom)];
+ int* aPointer = aCellIds->GetPointer(0);
+ for(int i=0;i<nbElems;i++,aPointer++){
+ int aCellId = *aPointer;
+ aCellIdsMapper.push_back(aCellId);
+ vtkCell* aCell = theInput->GetCell(aCellId);
+ int nbPointsInCell = aCell->GetNumberOfPoints();
+ for(int j=0;j<nbPointsInCell;j++)
+ aConn.push_back(aCell->GetPointId(j)+1);
+ int aPrevPos = aPolyedreFaces.back();
+ aPolyedreFaces.push_back(aPrevPos + nbPointsInCell);
+ aPrevPos = aPolyedreInds.back();
+ aPolyedreInds.push_back(aPrevPos + 1/*aNbFaces*/);
+ aFamilyNums.push_back(-1);
+
+ if(theCellsMapper){
+ if(theCellsMapper->GetNumberOfComponents()==2)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId*2));
+ else if(theCellsMapper->GetNumberOfComponents()==1)
+ aElemNums.push_back(*theCellsMapper->GetPointer(aCellId));
+ }
+ }
+
+
+ MED::PPolyedreInfo aCellInfo = theMed->CrPolyedreInfo(theMeshInfo,
+ theEntity,
+ VTK2MED(theVTKGeom),
+ aPolyedreInds,
+ aPolyedreFaces,
+ aConn,
+ aConnType,
+ aFamilyNums,
+ aElemNums);
+ theMed->SetPolyedreInfo(aCellInfo);
+ }
+
+ aCellIds->Delete();
+
+ return 0;
+}
+
+/*!
+ \brief Converts geometry to med (private auxiliary method)
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::Geometry2MED(vtkDataSet* aInput,
+ MED::PWrapper myMed,
+ MED::PMeshInfo aMeshInfo,
+ TGeom2CellIds& outGeom2CellIdMap)
+{
+ int aNbNodes = aInput->GetNumberOfPoints();
+ int aMeshDimension = aMeshInfo->GetDim();
+ // ----------------------- NODES -------------------------
+ vtkIntArray* aPointsMapper;
+ if(aInput->GetPointData())
+ aPointsMapper = dynamic_cast<vtkIntArray*>(aInput->GetPointData()->GetArray(myPointDataFieldNameIDS.c_str()));
+
+ MED::TFloatVector aCoordinates(aNbNodes*aMeshDimension);
+ MED::TIntVector anElemNumsNodes; // takes from VISU_POINTS_MAPPER array
+ MED::TIntVector aFamilyNumsNodes;
+ MED::TStringVector aCoordNamesNodes;
+ MED::TStringVector aCoordUnitsNodes;
+
+ vtkFloatingPointType aPntCoord[3];
+ if(aPointsMapper){
+ int nbComp = aPointsMapper->GetNumberOfComponents();
+ int *aPointsMapperPtr = aPointsMapper->GetPointer(0);
+ for(int i=0;i<aNbNodes;i++){
+ aInput->GetPoint(i,aPntCoord);
+ aCoordinates[i*3] = aPntCoord[0];
+ aCoordinates[i*3+1] = aPntCoord[1];
+ aCoordinates[i*3+2] = aPntCoord[2];
+ if(nbComp == 2){
+ anElemNumsNodes.push_back(*aPointsMapperPtr);
+ aPointsMapperPtr++;aPointsMapperPtr++;
+ }
+ else if (nbComp == 1){
+ anElemNumsNodes.push_back(*aPointsMapperPtr);
+ aPointsMapperPtr++;
+ }
+ else{
+ cout << "Error in file=|" << __FILE__<<"| line:[" << __LINE__ << "]" << endl;
+ cout << "Code must be adapted for more than 2 components array |VISU_POINTS_MAPPER|" << endl;
+ return 1;
+ }
+
+ }
+ } else {
+ for(int i=0;i<aNbNodes;i++){
+ aInput->GetPoint(i,aPntCoord);
+ aCoordinates[i*3] = aPntCoord[0];
+ aCoordinates[i*3+1] = aPntCoord[1];
+ aCoordinates[i*3+2] = aPntCoord[2];
+ }
+ }
+
+
+ MED::PNodeInfo aNodeInfo = myMed->CrNodeInfo(aMeshInfo,
+ aCoordinates,
+ MED::eFULL_INTERLACE,
+ MED::eCART,
+ aCoordNamesNodes,
+ aCoordUnitsNodes,
+ aFamilyNumsNodes,
+ anElemNumsNodes);
+ myMed->SetNodeInfo(aNodeInfo);
+
+ vtkIntArray* aCellsMapper;
+ if(vtkCellData* aCD = aInput->GetCellData())
+ aCellsMapper = dynamic_cast<vtkIntArray*>(aCD->GetArray(myCellDataFieldNameIDS.c_str()));
+
+ if(MYDEBUG)
+ {
+ // debug info
+ // print all cell types in the input
+ vtkCellTypes* aCellTypes = vtkCellTypes::New();
+ aInput->GetCellTypes(aCellTypes);
+ cout << "Cell types in the input data:"<<endl;
+ for(int aNbCellType = 0;aNbCellType<aCellTypes->GetNumberOfTypes();aNbCellType++)
+ cout << (int)(aCellTypes->GetCellType(aNbCellType)) << endl;
+ aCellTypes->Delete();
+ }
+
+ //----------------------
+ // Entity EDGES (eARETE)
+ //----------------------
+ vector<int> aNumberingConvertor;
+ {
+ // aVTKGeom->eSEG2
+ MED::EEntiteMaillage aEntity = MED::eMAILLE;//eARETE;
+ int aVTKGeom = 0;
+ int nbPointsInGeom = 0;
+
+ // aNumberingConvertor NULL - OK
+ aVTKGeom = VTK_LINE;
+ nbPointsInGeom = 2;
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ // aNumberingConvertor NULL - OK
+ // debug info: checked - OK
+ aVTKGeom = VTK_QUADRATIC_EDGE;
+ nbPointsInGeom = 3;
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ }
+ //----------------------------
+ // Entity FACES (eFACE)
+ // eTRIA3,eQUAD4,eTRIA6,eQUAD8
+ // ePOLYGONE
+ //----------------------------
+ {
+
+ MED::EEntiteMaillage aEntity = MED::eMAILLE;//MED::eFACE;
+ int aVTKGeom = 0;
+ int nbPointsInGeom = 0;
+
+ // debug info: checked - OK
+ aVTKGeom = VTK_TRIANGLE;
+ aNumberingConvertor.clear();
+ nbPointsInGeom = 3;
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ // debug info: checked - OK
+ aVTKGeom = VTK_QUAD;
+ nbPointsInGeom = 4;
+ aNumberingConvertor.clear();
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+
+ // debug info: checked - OK
+ aVTKGeom = VTK_QUADRATIC_TRIANGLE;
+ nbPointsInGeom = 6;
+ aNumberingConvertor.clear();
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ aVTKGeom = VTK_QUADRATIC_QUAD;
+ nbPointsInGeom = 8;
+ aNumberingConvertor.clear();
+ // 0,3,2,1,7,6,5,4
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ // debug info: checked - OK
+ CreatePolygons(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ outGeom2CellIdMap);
+
+ }
+ //----------------------------
+ // Entity CELLS (eMAILLE)
+ // eTETRA4,eHEXA8,ePENTA6,ePYRA5,
+ // eTETRA10,eHEXA20,ePOLYEDRE
+ //----------------------------
+ {
+
+ MED::EEntiteMaillage aEntity = MED::eMAILLE;
+ int aVTKGeom = 0;
+ int nbPointsInGeom = 0;
+
+ aVTKGeom = VTK_TETRA;
+ nbPointsInGeom = 4;
+ aNumberingConvertor.clear();
+ aNumberingConvertor.push_back(0);
+ aNumberingConvertor.push_back(2);
+ aNumberingConvertor.push_back(1);
+ aNumberingConvertor.push_back(3);
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ aVTKGeom = VTK_HEXAHEDRON;
+ nbPointsInGeom = 8;
+ aNumberingConvertor.clear();
+ // 0,3,2,1,4,7,6,5
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ aVTKGeom = VTK_WEDGE;
+ nbPointsInGeom = 6;
+ aNumberingConvertor.clear();
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ aVTKGeom = VTK_PYRAMID;
+ nbPointsInGeom = 5;
+ aNumberingConvertor.clear();
+ aNumberingConvertor.push_back(0);
+ aNumberingConvertor.push_back(3);
+ aNumberingConvertor.push_back(2);
+ aNumberingConvertor.push_back(1);
+ aNumberingConvertor.push_back(4);
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ aVTKGeom = VTK_QUADRATIC_TETRA;
+ nbPointsInGeom = 10;
+ aNumberingConvertor.clear();
+ // 0,2,1,3,6,5,4,7,9,8
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+
+ aVTKGeom = VTK_QUADRATIC_HEXAHEDRON;
+ nbPointsInGeom = 20;
+ aNumberingConvertor.clear();
+ // 0,3,2,1,4,7,6,5,11,10,9,8,15,14,13,12,16,19,18,17
+ CreateElements(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ aVTKGeom,
+ nbPointsInGeom,
+ aNumberingConvertor,
+ outGeom2CellIdMap);
+
+ // debug info: checked OK
+ CreatePolyedres(aInput,
+ aMeshInfo,
+ myMed,
+ aCellsMapper,
+ aEntity,
+ outGeom2CellIdMap);
+ }
+
+
+ return 0;
+}
+
+/*!
+ \brief Converts data to med (private auxiliary method)
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::Data2MED(std::vector<vtkDataSet*> theListForAdd,
+ MED::PWrapper myMed,
+ MED::PMeshInfo theMeshInfo,
+ TGeom2CellIds& theGeom2CellIdMap)
+{
+ typedef std::vector<vtkPointData*> TPDVec;
+ typedef std::vector<vtkCellData*> TCDVec;
+ typedef std::map<std::string,TPDVec> TNameToPointData;
+ typedef std::map<std::string,TCDVec> TNameToCellData;
+
+ TNameToPointData aName2PointData;
+ TNameToCellData aName2CellData;
+
+ MED::TErr* theErrCode = new MED::TErr();
+ MED::EGeometrieElement geomType;
+ MED::EEntiteMaillage entity;
+
+ // prepare data to create time stamps
+ const MED::TInt aNbGauss = 1;
+ MED::TProfileInfo::TInfo aTInfo("",0);
+ MED::PProfileInfo aPProfileInfo = myMed->CrProfileInfo( aTInfo );
+
+ MED::TGeom2Size aGeom2Size;
+ MED::TGeom2NbGauss aGeom2NbGauss;
+ MED::TGeom2Profile aTGeom2Profile;
+
+ int nbPointsInFirstData = 0;
+ if(theListForAdd.size()>0)
+ nbPointsInFirstData = (theListForAdd[0])->GetNumberOfPoints();
+
+ int nbClasses = theListForAdd.size();
+
+ for(int iClass=0;iClass<nbClasses;iClass++){
+ vtkDataSet* aClassPtr = theListForAdd[iClass];
+ if(aClassPtr->GetNumberOfPoints() != nbPointsInFirstData){
+ cout << "Warning in PointData: Some vtk file consist of number of points( " <<aClassPtr->GetNumberOfPoints()
+ << ") not equal number of points in first file("<<nbPointsInFirstData<<")"<<endl;
+ cout << "This data will be lost." << endl;
+ continue;
+ }
+
+ if(vtkPointData* aPD = aClassPtr->GetPointData()){
+ int nbArrays = aPD->GetNumberOfArrays();
+
+ for(int aArrNum=0;aArrNum<nbArrays;aArrNum++)
+ {
+ vtkDataArray* aArr = aPD->GetArray(aArrNum);
+ std::string aName = aArr->GetName();
+ std::set<std::string>::const_iterator aIgnoreIter = myIgnoringFieldList.find(aName);
+ if(aIgnoreIter!=myIgnoringFieldList.end())
+ continue;
+ (aName2PointData[aName]).push_back(aPD);
+ }
+ }
+
+ if(vtkCellData* aCD = aClassPtr->GetCellData()){
+ int nbArrays = aCD->GetNumberOfArrays();
+
+ for(int aArrNum=0;aArrNum<nbArrays;aArrNum++)
+ {
+ vtkDataArray* aArr = aCD->GetArray(aArrNum);
+ std::string aName = aArr->GetName();
+ std::set<std::string>::const_iterator aIgnoreIter = myIgnoringFieldList.find(aName);
+ if(aIgnoreIter!=myIgnoringFieldList.end())
+ continue;
+ (aName2CellData[aName]).push_back(aCD);
+ }
+ }
+ }
+
+
+ // PointData
+ int aLastField = 0;
+ {
+ geomType = MED::ePOINT1;
+ entity = MED::eNOEUD;
+ aGeom2Size[geomType] = nbPointsInFirstData;
+ aGeom2NbGauss[geomType] = aNbGauss;
+ aTGeom2Profile[geomType]= aPProfileInfo;
+
+ TNameToPointData::const_iterator aIter = aName2PointData.begin();
+ for(int iField=1;aIter!=aName2PointData.end();aIter++,iField++){
+ std::string aFieldName = aIter->first;
+ TPDVec aPD2Vec = (aIter->second);
+ int nbComp = 0;
+ if(aPD2Vec.size() >0 ){
+ if(vtkPointData* aPD0 = aPD2Vec[0]){
+ if(vtkDataArray* aArr0 = aPD0->GetArray(aFieldName.c_str()))
+ nbComp = aArr0->GetNumberOfComponents();
+ }
+ }
+
+ MED::PFieldInfo aFieldInfo = myMed->CrFieldInfo(theMeshInfo,
+ nbComp);
+
+ string aFieldName_PD = "Point " + aFieldName;
+ aFieldInfo->SetName(aFieldName_PD.c_str());
+
+ myMed->SetFieldInfo(aFieldInfo);
+
+ TPDVec::const_iterator aPDIter = aPD2Vec.begin();
+ for(int iTStamp=0;aPDIter!=aPD2Vec.end();aPDIter++,iTStamp++){
+ vtkPointData* aPD = *aPDIter;
+ vtkDataArray* aArr = aPD->GetArray(aFieldName.c_str());
+ MED::TFloat aTS = iTStamp < (int)myTStamps.size() ? myTStamps[ iTStamp ] : iTStamp;
+ MED::PTimeStampInfo aTempTimeStampInfo = myMed->CrTimeStampInfo (aFieldInfo,
+ entity,
+ aGeom2Size,
+ aGeom2NbGauss,
+ iTStamp,
+ iTStamp,
+ aTS);
+
+ MED::PTimeStampVal aTempTimeStampVal = myMed->CrTimeStampVal (aTempTimeStampInfo,
+ aTGeom2Profile);
+
+ MED::TMeshValue& aTMeshValue = aTempTimeStampVal->GetMeshValue(geomType);
+
+ MED::TValue& aValue = aTMeshValue.myValue; // float
+ int nbValues = aValue.size();
+ for(int i=0;i<nbValues;i++){
+ aValue[i] = *(float*)(aArr->GetVoidPointer(i));
+ }
+
+ myMed->SetTimeStamp( aTempTimeStampVal, theErrCode);
+ if(*theErrCode==0){
+ cout << "Error in "<<__FILE__<<"["<<__LINE__<<"] in method SetTimeStamp(...)"<<endl;
+ return 1;
+ }
+ }
+ aLastField = iField;
+ }
+ }
+ // CellData
+ {
+ MED::TEntityInfo aEntityInfo = myMed->GetEntityInfo(theMeshInfo);
+ aGeom2Size.clear();
+ aGeom2NbGauss.clear();
+ aTGeom2Profile.clear();
+ aGeom2Size = aEntityInfo[ MED::eMAILLE ];
+ entity = MED::eMAILLE;
+ MED::TGeom2Size::iterator geom_nb;
+ for ( geom_nb = aGeom2Size.begin(); geom_nb != aGeom2Size.end(); ++geom_nb ) { // loop on geometric types of cell
+ aGeom2NbGauss[ geom_nb->first ] = aNbGauss;
+ aTGeom2Profile[ geom_nb->first ] = aPProfileInfo;
+ }
+
+
+ TNameToCellData::const_iterator aIter = aName2CellData.begin();
+ for(int iField=1;aIter!=aName2CellData.end();aIter++,iField++){
+ std::string aFieldName = aIter->first;
+ TCDVec aCD2Vec = (aIter->second);
+ int nbComp = 0;
+ if(aCD2Vec.size() >0 ){
+ if(vtkCellData* aCD0 = aCD2Vec[0]){
+ if(vtkDataArray* aArr0 = aCD0->GetArray(aFieldName.c_str()))
+ nbComp = aArr0->GetNumberOfComponents();
+ }
+ }
+
+ MED::PFieldInfo aFieldInfo = myMed->CrFieldInfo(theMeshInfo,
+ nbComp);
+
+ string aFieldName_CD = "Cell " + aFieldName;
+ aFieldInfo->SetName(aFieldName_CD.c_str());
+
+ myMed->SetFieldInfo(aFieldInfo);
+
+ TCDVec::const_iterator aCDIter = aCD2Vec.begin();
+ for(int iTStamp=0;aCDIter!=aCD2Vec.end();aCDIter++,iTStamp++){
+ vtkCellData* aCD = *aCDIter;
+ vtkDataArray* aArr = aCD->GetArray(aFieldName.c_str());
+ MED::TFloat aTS = iTStamp < (int)myTStamps.size() ? myTStamps[ iTStamp ] : iTStamp;
+ MED::PTimeStampInfo aTempTimeStampInfo = myMed->CrTimeStampInfo (aFieldInfo,
+ entity,
+ aGeom2Size,
+ aGeom2NbGauss,
+ iTStamp,
+ iTStamp,
+ aTS);
+
+ MED::PTimeStampVal aTempTimeStampVal = myMed->CrTimeStampVal (aTempTimeStampInfo,
+ aTGeom2Profile);
+
+ for ( geom_nb = aGeom2Size.begin(); geom_nb != aGeom2Size.end(); ++geom_nb ) { // loop on geometric types of cell
+ geomType = geom_nb->first;
+ TCellIds& aCellIds = theGeom2CellIdMap[geomType];
+
+ MED::TMeshValue& aTMeshValue = aTempTimeStampVal->GetMeshValue(geomType);
+
+ MED::TValue& aValue = aTMeshValue.myValue; // float
+ int nbValues = aValue.size();
+ int nbCellIds = aCellIds.size();
+ if(nbValues!=nbCellIds*nbComp){
+ cout << "Warning in "<<__FILE__<<"["<<__LINE__<<"] the data for geometry:"<<geomType<<" will be ignored"<<endl;
+ continue;
+ }
+ if(MYDEBUG_VALUES) cout << "Geom["<<geomType<<"]"<<endl;
+ for(int i=0;i<nbCellIds;i++){
+ if(MYDEBUG_VALUES) cout << "\t|";
+ for(int iComp=0;iComp<nbComp;iComp++){
+ aValue[nbComp*i+iComp] = *(float*)(aArr->GetVoidPointer(nbComp*aCellIds[i]+iComp));
+ if(MYDEBUG_VALUES) cout << aValue[nbComp*i+iComp] << " ";
+ }
+ if(MYDEBUG_VALUES) cout << "|" << endl;
+ }
+
+ myMed->SetTimeStamp( aTempTimeStampVal, theErrCode);
+ if(*theErrCode==0){
+ cout << "Error in "<<__FILE__<<"["<<__LINE__<<"] in method SetTimeStamp(...)"<<endl;
+ return 1;
+ }
+ }
+ }
+ }
+ }
+
+
+ return 0;
+}
+
+/*!
+ \brief Writes data to MED file
+ \return 0 if operation has been completed successfully, 1 otherwise
+*/
+int
+VTK2MED_Convertor
+::Execute()
+{
+ int aStatus = 1;
+ if (myFirstVTKFileName.size() == 0 ||
+ myMEDFileName.size() == 0)
+ {
+ cout << "Error! Bad input file names output med file name or vtk file name."<<endl;
+ cout << "Exit."<<endl;
+ return 1;
+ }
+
+ MED::PWrapper myMed;
+ MED::PMeshInfo aMeshInfo;
+ int aMeshDimension = 3;
+ myMed = CrWrapper(myMEDFileName.c_str(),myVersion);
+ aMeshInfo = myMed->CrMeshInfo(aMeshDimension,myMeshName.c_str());
+ myMed->SetMeshInfo(aMeshInfo);
+
+ {
+ typedef vtkDataSetReader TReader;
+ TReader* aReader = TReader::New();
+ aReader->SetFileName(myFirstVTKFileName.c_str());
+ aReader->Update();
+ TGeom2CellIds myGeom2CellIds;
+
+ typedef std::vector<vtkDataSet*> TListUG;
+ TListUG aList;
+
+ if(aReader->IsFilePolyData())
+ {
+ if(MYDEBUG) cout << "PolyData" << endl;
+ typedef vtkPolyData TCommonType;
+ TCommonType* aInput = aReader->GetPolyDataOutput();
+
+ aStatus = Geometry2MED(aInput,
+ myMed,
+ aMeshInfo,
+ myGeom2CellIds);
+
+
+ TCommonType* aUG1 = TCommonType::New();
+ aUG1->ShallowCopy(aInput);
+ vtkDataSet* aTmp1 = dynamic_cast<vtkDataSet*>(aUG1);
+ aList.push_back(aTmp1);
+
+ TVectorString::iterator aFilesIter = myDataVTKFileNames.begin();
+ for(;aFilesIter!=myDataVTKFileNames.end();aFilesIter++){
+ aReader->SetFileName((*aFilesIter).c_str());
+ aReader->Update();
+ TCommonType* aUG2 = TCommonType::New();
+ aUG2->ShallowCopy(aReader->GetPolyDataOutput());
+ vtkDataSet* aTmp2 = dynamic_cast<vtkDataSet*>(aUG2);
+ aList.push_back(aTmp2);
+ }
+ } else if (aReader->IsFileUnstructuredGrid()){
+ if (MYDEBUG) cout << "UnstructuredGrid" << endl;
+ typedef vtkUnstructuredGrid TCommonType;
+ TCommonType* aInput = aReader->GetUnstructuredGridOutput();
+
+ aStatus = Geometry2MED(aInput,
+ myMed,
+ aMeshInfo,
+ myGeom2CellIds);
+
+
+ TCommonType* aUG1 = TCommonType::New();
+ aUG1->ShallowCopy(aInput);
+ vtkDataSet* aTmp1 = dynamic_cast<vtkDataSet*>(aUG1);
+ aList.push_back(aTmp1);
+
+ TVectorString::iterator aFilesIter = myDataVTKFileNames.begin();
+ for(;aFilesIter!=myDataVTKFileNames.end();aFilesIter++){
+ aReader->SetFileName((*aFilesIter).c_str());
+ aReader->Update();
+ TCommonType* aUG2 = TCommonType::New();
+ aUG2->ShallowCopy(aReader->GetUnstructuredGridOutput());
+ vtkDataSet* aTmp2 = dynamic_cast<vtkDataSet*>(aUG2);
+ aList.push_back(aTmp2);
+ }
+ }
+
+ Data2MED(aList,
+ myMed,
+ aMeshInfo,
+ myGeom2CellIds);
+
+
+ // clear aList by removing of unstructured grids
+ TListUG::iterator aIter = aList.begin();
+ for(;aIter!=aList.end();aIter++)
+ (*aIter)->Delete();
+
+ aReader->Delete();
+ }
+
+ return aStatus;
+}
