1 // Copyright (C) 2010-2014 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, or (at your option) any later version.
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 "vtkELNOFilter.h"
21 #include "vtkInformation.h"
22 #include "vtkInformationVector.h"
23 #include "vtkObjectFactory.h"
24 #include "vtkPolyDataAlgorithm.h"
25 #include "vtkPolyData.h"
26 #include "vtkIdTypeArray.h"
27 #include "vtkFieldData.h"
28 #include "vtkCellData.h"
29 #include "vtkPointData.h"
30 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
31 #include "vtkQuadratureSchemeDefinition.h"
32 #include "vtkUnstructuredGrid.h"
34 #include "MEDUtilities.hxx"
36 //vtkCxxRevisionMacro(vtkELNOFilter, "$Revision: 1.2.2.2 $");
37 vtkStandardNewMacro(vtkELNOFilter);
39 vtkELNOFilter::vtkELNOFilter()
41 this->ShrinkFactor = 0.5;
44 vtkELNOFilter::~vtkELNOFilter()
48 int vtkELNOFilter::RequestData(vtkInformation *request,
49 vtkInformationVector **input, vtkInformationVector *output)
51 vtkUnstructuredGrid *usgIn = vtkUnstructuredGrid::SafeDownCast(
52 input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
54 vtkPolyData *pdOut = vtkPolyData::SafeDownCast(
55 output->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT()));
57 vtkDataArray* array = this->GetInputArrayToProcess(0, input);
58 vtkIdTypeArray* offsets = vtkIdTypeArray::SafeDownCast(
59 this->GetInputArrayToProcess(0, input));
61 if(usgIn == NULL || offsets == NULL || pdOut == NULL)
63 vtkDebugMacro("vtkELNOFilter no correctly configured : offsets = " << offsets);
67 vtkInformation *info = offsets->GetInformation();
68 vtkInformationQuadratureSchemeDefinitionVectorKey *key =
69 vtkQuadratureSchemeDefinition::DICTIONARY();
72 vtkDebugMacro("Dictionary is not present in array " << offsets->GetName() << " " << offsets << " Aborting." );
76 int res = this->Superclass::RequestData(request, input, output);
82 int dictSize = key->Size(info);
83 vtkQuadratureSchemeDefinition **dict =
84 new vtkQuadratureSchemeDefinition *[dictSize];
85 key->GetRange(info, dict, 0, 0, dictSize);
87 vtkIdType ncell = usgIn->GetNumberOfCells();
88 vtkPoints *points = pdOut->GetPoints();
90 for(vtkIdType cellId = 0; cellId < ncell; cellId++)
92 vtkIdType offset = offsets->GetValue(cellId);
93 int cellType = usgIn->GetCellType(cellId);
94 // a simple check to see if a scheme really exists for this cell type.
95 // should not happen if the cell type has not been modified.
96 if(dict[cellType] == NULL)
98 int np = dict[cellType]->GetNumberOfQuadraturePoints();
99 double center[3] = {0, 0, 0};
100 for(int id = start; id < start + np; id++)
102 double *position = points->GetPoint(id);
103 center[0] += position[0];
104 center[1] += position[1];
105 center[2] += position[2];
110 for(int id = start; id < start + np; id++)
112 double *position = points->GetPoint(id);
114 newpos[0] = position[0] * this->ShrinkFactor + center[0] * (1
115 - this->ShrinkFactor);
116 newpos[1] = position[1] * this->ShrinkFactor + center[1] * (1
117 - this->ShrinkFactor);
118 newpos[2] = position[2] * this->ShrinkFactor + center[2] * (1
119 - this->ShrinkFactor);
120 points->SetPoint(id, newpos);
124 //// bug EDF 8407 PARAVIS - mantis 22610
125 vtkFieldData *fielddata(usgIn->GetFieldData());
126 for(int index = 0; index < fielddata->GetNumberOfArrays(); index++)
128 vtkDataArray *data(fielddata->GetArray(index));
129 vtkInformation *info(data->GetInformation());
130 const char *arrayOffsetName(info->Get(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME()));
134 vtkCellData *cellData(usgIn->GetCellData());
135 vtkDataArray *offData(cellData->GetArray(arrayOffsetName));
136 isELNO=offData->GetInformation()->Get(MEDUtilities::ELNO())==1;
140 pdOut->GetPointData()->AddArray(data);
146 void vtkELNOFilter::PrintSelf(ostream& os, vtkIndent indent)
148 this->Superclass::PrintSelf(os, indent);
150 os << indent << "ShrinkFactor : " << this->ShrinkFactor << endl;