Back to index

salome-paravis  6.5.0
VTKMEDCouplingUMeshClient.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "VTKMEDCouplingUMeshClient.hxx"
00021 
00022 #include "vtkPoints.h"
00023 #include "vtkCellArray.h"
00024 #include "vtkDoubleArray.h"
00025 #include "vtkSmartPointer.h"
00026 #include "vtkUnstructuredGrid.h"
00027 
00028 #include <set>
00029 #include <vector>
00030 #include <string>
00031 #include <algorithm>
00032 
00033 static const int ParaMEDMEM2VTKTypeTraducer[32]={1,3,21,5,9,7,22,-1,23,-1,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,-1,-1,-1,25,42};
00034 
00035 void ParaMEDMEM2VTK::FillMEDCouplingUMeshInstanceFrom(SALOME_MED::MEDCouplingUMeshCorbaInterface_ptr meshPtr, vtkUnstructuredGrid *ret, bool& isPolyh)
00036 {
00037   meshPtr->Register();
00038   //
00039   SALOME_TYPES::ListOfDouble *tinyD;
00040   SALOME_TYPES::ListOfLong *tinyI;
00041   SALOME_TYPES::ListOfString *tinyS;
00042   meshPtr->getTinyInfo(tinyD,tinyI,tinyS);
00043   //
00044   int spaceDim=(*tinyI)[1];
00045   int nbOfNodes=(*tinyI)[2];
00046   int meshDim=(*tinyI)[5];
00047   int nbOfCells=(*tinyI)[6];
00048   int meshLength=(*tinyI)[7];
00049   std::string name((*tinyS)[0]);
00050   //std::vector<std::string> compoNames(spaceDim);
00051   //for(int i=0;i<spaceDim;i++)
00052   //  compoNames[i]=(*tinyS)[i+4];
00053   delete tinyD;
00054   delete tinyI;
00055   delete tinyS;
00056   //
00057   ret->Initialize();
00058   ret->Allocate(nbOfCells);
00059   vtkPoints *points=vtkPoints::New();
00060   vtkDoubleArray *da=vtkDoubleArray::New();
00061   da->SetNumberOfComponents(3);
00062   da->SetNumberOfTuples(nbOfNodes);
00063   double *pts=da->GetPointer(0);
00064   //
00065   SALOME_TYPES::ListOfLong *a1Corba;
00066   SALOME_TYPES::ListOfDouble *a2Corba;
00067   meshPtr->getSerialisationData(a1Corba,a2Corba);
00068   if(spaceDim==3)
00069     {
00070       int myLgth=a2Corba->length();
00071       for(int i=0;i<myLgth;i++)
00072         *pts++=(*a2Corba)[i];
00073     }
00074   else
00075     {
00076       int offset=3-spaceDim;
00077       for(int i=0;i<nbOfNodes;i++)
00078         {
00079           for(int j=0;j<spaceDim;j++)
00080             *pts++=(*a2Corba)[spaceDim*i+j];
00081           std::fill(pts,pts+offset,0.);
00082           pts+=offset;
00083         }
00084     }
00085   //
00086   vtkIdType *tmp=new vtkIdType[1000];
00087   isPolyh=false;
00088   for(int i=0;i<nbOfCells;i++)
00089     {
00090       int pos=(*a1Corba)[i];
00091       int pos1=(*a1Corba)[i+1];
00092       int nbOfNodeInCurCell=pos1-pos-1;
00093       int typeOfCell=(*a1Corba)[pos+nbOfCells+1];
00094       for(int j=0;j<nbOfNodeInCurCell;j++)
00095         tmp[j]=(*a1Corba)[pos+1+j+nbOfCells+1];
00096       int vtkType=ParaMEDMEM2VTKTypeTraducer[typeOfCell];
00097       if(vtkType!=42)
00098         ret->InsertNextCell(vtkType,nbOfNodeInCurCell,tmp);
00099       else
00100         {//polyhedron
00101           isPolyh=true;
00102           std::set<vtkIdType> s(tmp,tmp+nbOfNodeInCurCell);
00103           vtkSmartPointer<vtkCellArray> faces=vtkSmartPointer<vtkCellArray>::New();
00104           int nbOfFaces=std::count(tmp,tmp+nbOfNodeInCurCell,-1)+1;
00105           vtkIdType *work=tmp;
00106           for(int i=0;i<nbOfFaces;i++)
00107             {
00108               vtkIdType *work2=std::find(work,tmp+nbOfNodeInCurCell,-1);
00109               int nbOfNodesInFace=std::distance(work,work2);
00110               faces->InsertNextCell(nbOfNodesInFace,work);
00111               work=work2+1;
00112             }
00113           s.erase(-1);
00114           std::vector<vtkIdType> v(s.begin(),s.end());
00115           ret->InsertNextCell(VTK_POLYHEDRON,v.size(),&v[0],nbOfFaces,faces->GetPointer());
00116         }
00117     }
00118   delete [] tmp;
00119   //
00120   delete a1Corba;
00121   delete a2Corba;
00122   //
00123   ret->SetPoints(points);
00124   points->Delete();
00125   points->SetData(da);
00126   da->Delete();
00127   //
00128   meshPtr->UnRegister();
00129 }