1 // Copyright (C) 2010-2011 CEA/DEN, EDF R&D
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 #include "vtkMedUnstructuredGrid.h"
22 #include "vtkObjectFactory.h"
23 #include "vtkSmartPointer.h"
24 #include "vtkDataArray.h"
25 #include "vtkUnstructuredGrid.h"
26 #include "vtkPointData.h"
27 #include "vtkCellData.h"
28 #include "vtkIdList.h"
29 #include "vtkCellType.h"
30 #include "vtkInformation.h"
32 #include "vtkMedUtilities.h"
33 #include "vtkMedEntityArray.h"
34 #include "vtkMedMesh.h"
35 #include "vtkMedFile.h"
36 #include "vtkMedDriver.h"
37 #include "vtkMedFamilyOnEntityOnProfile.h"
38 #include "vtkMedFamilyOnEntity.h"
39 #include "vtkMedProfile.h"
40 #include "vtkMedUtilities.h"
41 #include "vtkMedStructElement.h"
42 #include "vtkMedVariableAttribute.h"
44 #include "vtkMultiProcessController.h"
46 vtkCxxSetObjectMacro(vtkMedUnstructuredGrid,Coordinates,vtkDataArray);
48 // vtkCxxRevisionMacro(vtkMedUnstructuredGrid, "$Revision$")
49 vtkStandardNewMacro(vtkMedUnstructuredGrid)
51 vtkMedUnstructuredGrid::vtkMedUnstructuredGrid()
53 this->Coordinates = NULL;
54 this->NumberOfPoints = 0;
57 vtkMedUnstructuredGrid::~vtkMedUnstructuredGrid()
59 this->SetCoordinates(NULL);
62 int vtkMedUnstructuredGrid::IsCoordinatesLoaded()
64 return this->Coordinates != NULL && this->Coordinates->GetNumberOfTuples()
65 == this->NumberOfPoints;
68 void vtkMedUnstructuredGrid::InitializeCellGlobalIds()
72 for(int id = 0; id < this->EntityArray->size(); id++)
74 vtkMedEntityArray* array = this->EntityArray->at(id);
79 if(array->GetEntity().EntityType == MED_NODE)
82 array->SetInitialGlobalId(ncells + 1);
83 ncells += array->GetNumberOfEntity();
87 void vtkMedUnstructuredGrid::ClearMedSupports()
89 this->Superclass::ClearMedSupports();
90 for(int id = 0; id < this->EntityArray->size(); id++)
92 vtkMedEntityArray* array = this->EntityArray->at(id);
97 void vtkMedUnstructuredGrid::LoadCoordinates()
99 this->GetParentMesh()->GetParentFile()->GetMedDriver()->LoadCoordinates(this);
102 double* vtkMedUnstructuredGrid::GetCoordTuple(med_int index)
104 return this->Coordinates->GetTuple(index);
107 vtkDataSet* vtkMedUnstructuredGrid::CreateVTKDataSet(
108 vtkMedFamilyOnEntityOnProfile* foep)
110 vtkMultiProcessController* controller =
111 vtkMultiProcessController::GetGlobalController();
113 vtkMedFamilyOnEntity* foe = foep->GetFamilyOnEntity();
114 vtkMedEntityArray* array = foe->GetEntityArray();
115 if(foep->GetValid() == 0)
120 vtkUnstructuredGrid* vtkugrid = vtkUnstructuredGrid::New();
122 // now we copy all the flagged nodes in the grid, shallow copy if possible
123 bool shallowCopyPoints=true;
125 if(this->GetParentMesh()->GetSpaceDimension()!=3)
127 shallowCopyPoints=false;
130 shallowCopyPoints = shallowCopyPoints && foep->CanShallowCopyPointField(NULL);
132 foe->GetParentGrid()->LoadCoordinates();
134 vtkIdType numberOfPoints;
135 vtkPoints* points=vtkPoints::New(this->GetCoordinates()->GetDataType());
136 vtkugrid->SetPoints(points);
139 vtkIdTypeArray* pointGlobalIds=vtkIdTypeArray::New();
140 pointGlobalIds->SetName("MED_POINT_ID");
141 pointGlobalIds->SetNumberOfComponents(1);
142 vtkugrid->GetPointData()->SetGlobalIds(pointGlobalIds);
143 pointGlobalIds->Delete();
145 vtkIdTypeArray* cellGlobalIds=vtkIdTypeArray::New();
146 cellGlobalIds->SetName("MED_CELL_ID");
147 cellGlobalIds->SetNumberOfComponents(1);
148 vtkugrid->GetCellData()->SetGlobalIds(cellGlobalIds);
149 cellGlobalIds->Delete();
151 if (shallowCopyPoints)
153 vtkugrid->GetPoints()->SetDataType(this->GetCoordinates()->GetDataType());
154 vtkugrid->GetPoints()->SetData(this->GetCoordinates());
156 numberOfPoints=this->GetNumberOfPoints();
157 pointGlobalIds->SetNumberOfTuples(numberOfPoints);
158 vtkIdType* ptr=pointGlobalIds->GetPointer(0);
159 for(int pid=0; pid<numberOfPoints; pid++)
164 vtkIdType currentIndex=0;
165 for(vtkIdType index=0; index<this->GetNumberOfPoints(); index++)
167 if(!foep->KeepPoint(index))
170 double coord[3]={0.0, 0.0, 0.0};
171 double * tuple=this->GetCoordinates()->GetTuple(index);
172 for(int dim=0; dim<this->GetParentMesh()->GetSpaceDimension()&&dim<3; dim++)
174 coord[dim]=tuple[dim];
176 vtkugrid->GetPoints()->InsertPoint(currentIndex, coord);
177 pointGlobalIds->InsertNextValue(index+1);
180 vtkugrid->GetPoints()->Squeeze();
181 pointGlobalIds->Squeeze();
182 numberOfPoints=currentIndex;
185 vtkMedStructElement* structelem = NULL;
186 int vtkType = VTK_EMPTY_CELL;
187 int nsupportcell = 1;
188 vtkSmartPointer<vtkIdTypeArray> supportIndex = vtkSmartPointer<vtkIdTypeArray>::New();
189 supportIndex->SetName("STRUCT_ELEMENT_INDEX");
190 supportIndex->SetNumberOfComponents(2);
191 supportIndex->GetInformation()->Set(vtkMedUtilities::STRUCT_ELEMENT_INDEX(), 1);
192 supportIndex->GetInformation()->Set(vtkAbstractArray::GUI_HIDE(), 1);
193 supportIndex->SetComponentName(0, "GLOBAL_ID");
194 supportIndex->SetComponentName(1, "LOCAL_ID");
196 if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
198 structelem = array->GetStructElement();
199 if(structelem == NULL)
205 vtkType = vtkMedUtilities::GetVTKCellType(
206 structelem->GetSupportGeometryType());
207 nsupportcell = structelem->GetSupportNumberOfCell();
208 if(structelem->GetSupportEntityType() != MED_CELL
209 || strcmp(structelem->GetName(), MED_PARTICLE_NAME) == 0)
211 // Special case : the support connectivity is implicit
213 // map this to points
214 vtkType = VTK_VERTEX;
215 nsupportcell = structelem->GetSupportNumberOfNode();
218 vtkugrid->GetInformation()->Set(vtkMedUtilities::STRUCT_ELEMENT(), structelem);
219 std::cout << "structelem->GetNumberOfVariableAttribute() = " << structelem->GetNumberOfVariableAttribute() << std::endl;
220 for(int varattid = 0; varattid<structelem->GetNumberOfVariableAttribute(); varattid++)
222 vtkMedVariableAttribute* varatt = structelem->GetVariableAttribute(varattid);
224 vtkAbstractArray* values = array->GetVariableAttributeValue(varatt);
225 vtkugrid->GetFieldData()->AddArray(values);
227 vtkugrid->GetCellData()->AddArray(supportIndex);
232 vtkType = vtkMedUtilities::GetVTKCellType(
233 array->GetEntity().GeometryType);
236 vtkSmartPointer<vtkIdList> pts = vtkSmartPointer<vtkIdList>::New();
237 vtkSmartPointer<vtkIdList> vtkpts = vtkSmartPointer<vtkIdList>::New();
238 vtkIdType intialGlobalId = array->GetInitialGlobalId();
239 vtkMedIntArray* pids = (foep->GetProfile()!=NULL?
240 foep->GetProfile()->GetIds() : NULL);
241 vtkIdType maxId = (pids!=NULL?
242 pids->GetNumberOfTuples():
243 array->GetNumberOfEntity());
245 int valid = foep->GetValid();
246 if (controller != NULL)
247 if (controller->GetNumberOfProcesses() > 1)
250 for (vtkIdType pindex = 0; pindex<maxId && valid; pindex++)
252 vtkIdType realIndex = (pids!=NULL?
253 pids->GetValue(pindex)-1:
256 if (!foep->KeepCell(realIndex))
259 array->GetCellVertices(realIndex, pts);
261 for(int sid = 0; sid < nsupportcell; sid++)
263 cellGlobalIds->InsertNextValue(intialGlobalId+pindex);
266 // The supportIndex array has 2 component, the first to give the index
267 // of the med cell in the variable attributes array, the second to give
268 // the index of the cell in the n cells composing the structural element
269 if(array->GetEntity().EntityType == MED_STRUCT_ELEMENT)
271 for(int sid = 0; sid < nsupportcell; sid++)
273 supportIndex->InsertNextTuple2(pindex, sid);
277 if (array->GetEntity().GeometryType==MED_POLYHEDRON)
279 if(vtkMedUtilities::FormatPolyhedronForVTK(foep, pindex, pts))
281 vtkugrid->InsertNextCell(VTK_POLYHEDRON, pts);
292 vtkpts->Initialize();
293 vtkpts->SetNumberOfIds(pts->GetNumberOfIds());
295 for(vtkIdType node=0; node<pts->GetNumberOfIds(); node++)
297 vtkIdType pid = pts->GetId(node);
298 vtkIdType ptid=foep->GetVTKPointIndex(pid);
299 if(ptid < 0 || ptid >= vtkugrid->GetNumberOfPoints())
301 vtkDebugMacro("Index error, this cell" <<
302 " is not on this profile");
303 #ifndef MedReader_HAVE_PARALLEL_INFRASTRUCTURE
311 vtkpts->SetId(vtkMedUtilities::MedToVTKIndex(vtkType, node), ptid);
314 // for strutural elements, insert nsupportcell instead of only one
317 vtkIdType npts = vtkpts->GetNumberOfIds() / nsupportcell;
318 for(int sid=0; sid < nsupportcell; sid++)
320 vtkIdType* ptids = vtkpts->GetPointer(sid * npts);
321 vtkugrid->InsertNextCell(vtkType, npts, ptids);
326 vtkugrid->InsertNextCell(vtkType, vtkpts);
331 if(vtkugrid->GetNumberOfCells() == vtkugrid->GetNumberOfPoints())
333 for(int fieldId = 0; fieldId < vtkugrid->GetFieldData()->GetNumberOfArrays(); fieldId++)
335 vtkDataArray* fieldData = vtkugrid->GetFieldData()->GetArray(fieldId);
337 if(fieldData->GetNumberOfTuples() == vtkugrid->GetNumberOfCells())
339 vtkugrid->GetPointData()->AddArray(fieldData);
343 vtkDataArray* real_fieldData = fieldData->NewInstance();
344 real_fieldData->SetName(fieldData->GetName());
345 real_fieldData->SetNumberOfComponents(fieldData->GetNumberOfComponents());
347 for(vtkIdType cellId = 0; cellId < vtkugrid->GetNumberOfPoints(); cellId++)
349 vtkIdType supportId = static_cast<int>(supportIndex->GetTuple2(cellId)[0]);
350 real_fieldData->InsertNextTuple(fieldData->GetTuple(supportId));
353 vtkugrid->GetPointData()->AddArray(real_fieldData);
354 real_fieldData->Delete();
362 void vtkMedUnstructuredGrid::PrintSelf(ostream& os, vtkIndent indent)
364 this->Superclass::PrintSelf(os, indent);
365 PRINT_IVAR(os, indent, NumberOfPoints);