1 // Copyright (C) 2010-2022 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"
31 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
32 #include "vtkQuadratureSchemeDefinition.h"
33 #include "vtkUnstructuredGrid.h"
35 #include "MEDUtilities.hxx"
36 #include "InterpKernelAutoPtr.hxx"
38 //vtkCxxRevisionMacro(vtkELNOFilter, "$Revision: 1.2.2.2 $");
39 vtkStandardNewMacro(vtkELNOFilter)
41 vtkELNOFilter::vtkELNOFilter()
43 this->ShrinkFactor = 0.5;
46 vtkELNOFilter::~vtkELNOFilter()
50 int vtkELNOFilter::RequestData(vtkInformation *request, vtkInformationVector **input, vtkInformationVector *output)
52 vtkUnstructuredGrid *usgIn(vtkUnstructuredGrid::SafeDownCast( input[0]->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT())));
53 vtkPolyData *pdOut(vtkPolyData::SafeDownCast(output->GetInformationObject(0)->Get(vtkDataObject::DATA_OBJECT())));
55 //vtkDataArray *array(this->GetInputArrayToProcess(0, input)); // todo: unused
56 vtkIdTypeArray* offsets(vtkIdTypeArray::SafeDownCast(this->GetInputArrayToProcess(0, input)));
58 if(usgIn == NULL || offsets == NULL || pdOut == NULL)
60 vtkDebugMacro("vtkELNOFilter no correctly configured : offsets = " << offsets);
64 vtkInformation *info(offsets->GetInformation());
65 vtkInformationQuadratureSchemeDefinitionVectorKey *key(vtkQuadratureSchemeDefinition::DICTIONARY());
68 vtkDebugMacro("Dictionary is not present in array " << offsets->GetName() << " " << offsets << " Aborting." );
72 int res(this->Superclass::RequestData(request, input, output));
76 int dictSize(key->Size(info));
77 vtkQuadratureSchemeDefinition **dict = new vtkQuadratureSchemeDefinition *[dictSize];
78 key->GetRange(info, dict, 0, 0, dictSize);
80 vtkIdType ncell(usgIn->GetNumberOfCells());
81 vtkPoints *points(pdOut->GetPoints());
83 for(vtkIdType cellId = 0; cellId < ncell; cellId++)
85 //vtkIdType offset(offsets->GetValue(cellId)); // todo: unused
86 int cellType(usgIn->GetCellType(cellId));
87 // a simple check to see if a scheme really exists for this cell type.
88 // should not happen if the cell type has not been modified.
89 if(dict[cellType] == NULL)
91 int np = dict[cellType]->GetNumberOfQuadraturePoints();
92 double center[3] = {0, 0, 0};
93 for(int id = start; id < start + np; id++)
95 double *position = points->GetPoint(id);
96 center[0] += position[0];
97 center[1] += position[1];
98 center[2] += position[2];
103 for(int id = start; id < start + np; id++)
105 double *position = points->GetPoint(id);
107 newpos[0] = position[0] * this->ShrinkFactor + center[0] * (1 - this->ShrinkFactor);
108 newpos[1] = position[1] * this->ShrinkFactor + center[1] * (1 - this->ShrinkFactor);
109 newpos[2] = position[2] * this->ShrinkFactor + center[2] * (1 - this->ShrinkFactor);
110 points->SetPoint(id, newpos);
114 //// bug EDF 8407 PARAVIS - mantis 22610
115 vtkFieldData *fielddata(usgIn->GetFieldData());
116 for(int index=0;index<fielddata->GetNumberOfArrays();index++)
118 vtkDataArray *data(fielddata->GetArray(index));
119 vtkInformation *info(data->GetInformation());
120 const char *arrayOffsetName(info->Get(vtkQuadratureSchemeDefinition::QUADRATURE_OFFSET_ARRAY_NAME()));
121 vtkIdTypeArray *offData(0);
125 vtkCellData *cellData(usgIn->GetCellData());
126 vtkDataArray *offDataTmp(cellData->GetArray(arrayOffsetName));
127 isELNO=offDataTmp->GetInformation()->Get(MEDUtilities::ELNO())==1;
128 offData=dynamic_cast<vtkIdTypeArray *>(offDataTmp);
130 if(isELNO && offData)
132 vtkIdType nbCellsInput(usgIn->GetNumberOfCells());
134 continue ;//no cells -> no fields
135 // First trying to detected if we are in the special case where data can be used directly. To detect that look at offData. If offData.front()==0 && offData->back()+NbOfNodesInLastCell==data->GetNumberOfTuples() OK.
136 vtkCell *cell(usgIn->GetCell(nbCellsInput-1));
137 bool statement0(offData->GetTuple1(0)==0);
138 bool statement1(offData->GetTuple1(nbCellsInput-1)+cell->GetNumberOfPoints()==data->GetNumberOfTuples());
139 if(statement0 && statement1)
140 pdOut->GetPointData()->AddArray(data);//We are lucky ! No extraction needed.
142 {//not lucky ! Extract values from data. A previous threshold has been done... Bug EDF8662
143 vtkDataArray *newArray(data->NewInstance());
144 newArray->SetName(data->GetName());
145 pdOut->GetPointData()->AddArray(newArray);
146 newArray->SetNumberOfComponents(data->GetNumberOfComponents());
147 newArray->SetNumberOfTuples(pdOut->GetNumberOfPoints());
148 newArray->CopyComponentNames(data);
150 vtkIdType *offsetPtr(offData->GetPointer(0));
152 for(vtkIdType cellId=0;cellId<nbCellsInput;cellId++)
154 vtkCell *cell(usgIn->GetCell(cellId));
155 vtkIdType nbPoints(cell->GetNumberOfPoints())/*,offset(offsetPtr[cellId])*/; // todo: offset is unused
156 for(vtkIdType j=0;j<nbPoints;j++,zeId++)
157 newArray->SetTuple(zeId,offsetPtr[cellId]+j,data);
162 AttachCellFieldsOn(usgIn,pdOut->GetCellData(),pdOut->GetNumberOfCells());
166 void vtkELNOFilter::PrintSelf(ostream& os, vtkIndent indent)
168 this->Superclass::PrintSelf(os, indent);
169 os << indent << "ShrinkFactor : " << this->ShrinkFactor << endl;
173 * This method attach fields on cell of \a inGrid and add it as a point data in \a outData.
175 void vtkELNOFilter::AttachCellFieldsOn(vtkUnstructuredGrid *inGrid, vtkCellData *outData, int nbCellsOut)
177 vtkCellData *cd(inGrid->GetCellData());
178 int nbOfArrays(cd->GetNumberOfArrays());
179 vtkIdType nbCells(inGrid->GetNumberOfCells());
182 INTERP_KERNEL::AutoPtr<vtkIdType> tmpPtr(new vtkIdType[nbCells]);
183 for(vtkIdType cellId=0;cellId<nbCells;cellId++)
185 vtkCell *cell(inGrid->GetCell(cellId));
186 tmpPtr[cellId]=cell->GetNumberOfPoints();
188 for(int index=0;index<nbOfArrays;index++)
190 vtkDataArray *data(cd->GetArray(index));
191 vtkDataArray *newArray(data->NewInstance());
192 newArray->SetName(data->GetName());
193 outData->AddArray(newArray);
194 newArray->SetNumberOfComponents(data->GetNumberOfComponents());
195 newArray->SetNumberOfTuples(nbCellsOut);
196 newArray->CopyComponentNames(data);
199 for(vtkIdType cellId=0;cellId<nbCells;cellId++)
201 for(vtkIdType j=0;j<tmpPtr[cellId];j++,offset++)
202 newArray->SetTuple(offset,cellId,data);