]> SALOME platform Git repositories - modules/paravis.git/blob - src/Plugins/ParaMEDCorba/VTKMEDCouplingUMeshClient.cxx
Salome HOME
Add missing MEDFile dependency in DevelopedSurface filter CMakeLists.txt
[modules/paravis.git] / src / Plugins / ParaMEDCorba / VTKMEDCouplingUMeshClient.cxx
1 // Copyright (C) 2010-2016  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
20 #include "VTKMEDCouplingUMeshClient.hxx"
21
22 #include "vtkPoints.h"
23 #include "vtkCellArray.h"
24 #include "vtkDoubleArray.h"
25 #include "vtkSmartPointer.h"
26 #include "vtkUnstructuredGrid.h"
27
28 #include <set>
29 #include <vector>
30 #include <string>
31 #include <algorithm>
32
33 static const int ParaMEDMEM2VTKTypeTraducer[34]={1,3,21,5,9,7,22,34,23,28,-1,-1,-1,-1,10,14,13,-1,12,-1,24,-1,16,27,-1,26,-1,29,-1,-1,25,42,-1,4};
34
35 void ParaMEDMEM2VTK::FillMEDCouplingUMeshInstanceFrom(SALOME_MED::MEDCouplingUMeshCorbaInterface_ptr meshPtr, vtkUnstructuredGrid *ret, bool& isPolyh)
36 {
37   meshPtr->Register();
38   //
39   SALOME_TYPES::ListOfDouble *tinyD;
40   SALOME_TYPES::ListOfLong *tinyI;
41   SALOME_TYPES::ListOfString *tinyS;
42   meshPtr->getTinyInfo(tinyD,tinyI,tinyS);
43   //
44   int spaceDim=(*tinyI)[1];
45   int nbOfNodes=(*tinyI)[2];
46   int meshDim=(*tinyI)[5];
47   int nbOfCells=(*tinyI)[6];
48   int meshLength=(*tinyI)[7];
49   std::string name((*tinyS)[0]);
50   //std::vector<std::string> compoNames(spaceDim);
51   //for(int i=0;i<spaceDim;i++)
52   //  compoNames[i]=(*tinyS)[i+4];
53   delete tinyD;
54   delete tinyI;
55   delete tinyS;
56   //
57   ret->Initialize();
58   ret->Allocate(nbOfCells);
59   vtkPoints *points=vtkPoints::New();
60   vtkDoubleArray *da=vtkDoubleArray::New();
61   da->SetNumberOfComponents(3);
62   da->SetNumberOfTuples(nbOfNodes);
63   double *pts=da->GetPointer(0);
64   //
65   SALOME_TYPES::ListOfLong *a1Corba;
66   SALOME_TYPES::ListOfDouble *a2Corba;
67   meshPtr->getSerialisationData(a1Corba,a2Corba);
68   if(spaceDim==3)
69     {
70       int myLgth=a2Corba->length();
71       for(int i=0;i<myLgth;i++)
72         *pts++=(*a2Corba)[i];
73     }
74   else
75     {
76       int offset=3-spaceDim;
77       for(int i=0;i<nbOfNodes;i++)
78         {
79           for(int j=0;j<spaceDim;j++)
80             *pts++=(*a2Corba)[spaceDim*i+j];
81           std::fill(pts,pts+offset,0.);
82           pts+=offset;
83         }
84     }
85   //
86   vtkIdType *tmp=new vtkIdType[1000];
87   isPolyh=false;
88   for(int i=0;i<nbOfCells;i++)
89     {
90       int pos=(*a1Corba)[i];
91       int pos1=(*a1Corba)[i+1];
92       int nbOfNodeInCurCell=pos1-pos-1;
93       int typeOfCell=(*a1Corba)[pos+nbOfCells+1];
94       for(int j=0;j<nbOfNodeInCurCell;j++)
95         tmp[j]=(*a1Corba)[pos+1+j+nbOfCells+1];
96       int vtkType=ParaMEDMEM2VTKTypeTraducer[typeOfCell];
97       if(vtkType!=42)
98         ret->InsertNextCell(vtkType,nbOfNodeInCurCell,tmp);
99       else
100         {//polyhedron
101           isPolyh=true;
102           std::set<vtkIdType> s(tmp,tmp+nbOfNodeInCurCell);
103           vtkSmartPointer<vtkCellArray> faces=vtkSmartPointer<vtkCellArray>::New();
104           int nbOfFaces=std::count(tmp,tmp+nbOfNodeInCurCell,-1)+1;
105           vtkIdType *work=tmp;
106           for(int i=0;i<nbOfFaces;i++)
107             {
108               vtkIdType *work2=std::find(work,tmp+nbOfNodeInCurCell,-1);
109               int nbOfNodesInFace=std::distance(work,work2);
110               faces->InsertNextCell(nbOfNodesInFace,work);
111               work=work2+1;
112             }
113           s.erase(-1);
114           std::vector<vtkIdType> v(s.begin(),s.end());
115           ret->InsertNextCell(VTK_POLYHEDRON,v.size(),&v[0],nbOfFaces,faces->GetPointer());
116         }
117     }
118   delete [] tmp;
119   //
120   delete a1Corba;
121   delete a2Corba;
122   //
123   ret->SetPoints(points);
124   points->Delete();
125   points->SetData(da);
126   da->Delete();
127   //
128   meshPtr->UnRegister();
129 }