1 // Copyright (C) 2021 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
19 // Author : Anthony Geay
21 #include "vtkExtractThreeD.h"
23 #include <NormalizedGeometricTypes>
24 #include <CellModel.hxx>
25 #include <InterpKernelException.hxx>
26 #include <MEDFileFieldOverView.hxx>
28 #include <vtkAdjacentVertexIterator.h>
29 #include <vtkAlgorithmOutput.h>
30 #include <vtkCellData.h>
31 #include <vtkCharArray.h>
32 #include <vtkDataArraySelection.h>
33 #include <vtkDataArrayTemplate.h>
34 #include <vtkDataObjectTreeIterator.h>
35 #include <vtkDataSet.h>
36 #include <vtkDataSetAttributes.h>
37 #include <vtkDemandDrivenPipeline.h>
38 #include <vtkDoubleArray.h>
39 #include <vtkExecutive.h>
40 #include <vtkInEdgeIterator.h>
41 #include <vtkInformation.h>
42 #include <vtkInformationDataObjectKey.h>
43 #include <vtkInformationStringKey.h>
44 #include <vtkInformationVector.h>
45 #include <vtkIntArray.h>
46 #include <vtkMultiBlockDataSet.h>
47 #include <vtkMutableDirectedGraph.h>
48 #include <vtkObjectFactory.h>
49 #include <vtkPointData.h>
50 #include <vtkStreamingDemandDrivenPipeline.h>
51 #include <vtkStringArray.h>
52 #include <vtkThreshold.h>
53 #include <vtkTimeStamp.h>
54 #include <vtkUnsignedCharArray.h>
55 #include <vtkUnstructuredGrid.h>
56 #include <vtkVariantArray.h>
62 vtkStandardNewMacro(vtkExtractThreeD);
66 class vtkExtractThreeD::vtkExtractThreeDInternal
69 vtkExtractThreeDInternal() {}
70 std::vector<int> getIdsToKeep() const { return _types; }
71 void setIdsToKeep(const std::vector<int> &types) { _types = types; }
74 std::vector<int> _types;
77 vtkExtractThreeD::vtkExtractThreeD() : Internal(new vtkExtractThreeDInternal)
81 vtkExtractThreeD::~vtkExtractThreeD()
83 delete this->Internal;
86 int vtkExtractThreeD::RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
90 //std::cerr << "########################################## vtkExtractThreeD::RequestInformation ##########################################" << std::endl;
91 vtkInformation *outInfo(outputVector->GetInformationObject(0));
92 vtkInformation *inputInfo(inputVector[0]->GetInformationObject(0));
95 vtkDataObject *inp(inputInfo->Get(vtkDataObject::DATA_OBJECT()));
96 if (vtkDataSet::SafeDownCast(inp))
97 input = vtkDataSet::SafeDownCast(inp);
100 vtkMultiBlockDataSet *inputTmp(vtkMultiBlockDataSet::SafeDownCast(inp));
103 if (inputTmp->GetNumberOfBlocks() != 1)
105 vtkDebugMacro("vtkExtractThreeD::RequestInformation : input vtkMultiBlockDataSet must contain exactly 1 block !");
108 vtkDataSet *blk0(vtkDataSet::SafeDownCast(inputTmp->GetBlock(0)));
111 vtkDebugMacro("vtkExtractThreeD::RequestInformation : the single block in input vtkMultiBlockDataSet must be a vtkDataSet instance !");
118 vtkDebugMacro("vtkExtractThreeD::RequestInformation : supported input are vtkDataSet or vtkMultiBlockDataSet !");
124 vtkIdType nbOfCells(input->GetNumberOfCells());
125 std::map<int, INTERP_KERNEL::NormalizedCellType> m;
126 std::set<int> typesToKeep;
127 for (vtkIdType cellId = 0; cellId < nbOfCells; cellId++)
129 int vtkCt(input->GetCellType(cellId));
130 const std::map<int, INTERP_KERNEL::NormalizedCellType>::const_iterator it(m.find(vtkCt));
133 const unsigned char *pos(std::find(MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE, MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE + MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH, vtkCt));
134 if (pos == MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE + MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
136 vtkDebugMacro("vtkExtractThreeD::RequestInformation : cell #" << cellId << " has unrecognized type !");
139 INTERP_KERNEL::NormalizedCellType mcCtype((INTERP_KERNEL::NormalizedCellType)std::distance(MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE, pos));
140 const INTERP_KERNEL::CellModel &cm(INTERP_KERNEL::CellModel::GetCellModel(mcCtype));
141 if (cm.getDimension() == 3)
142 typesToKeep.insert(vtkCt);
145 std::vector<int> typesToKeep2(typesToKeep.begin(), typesToKeep.end());
146 this->Internal->setIdsToKeep(typesToKeep2);
149 catch (INTERP_KERNEL::Exception &e)
151 std::cerr << "Exception has been thrown in vtkExtractThreeD::RequestInformation : " << e.what() << std::endl;
157 vtkDataSet *FilterFamilies(vtkDataSet *input, const std::vector<int> &idsToKeep)
159 const int VTK_DATA_ARRAY_DELETE = vtkDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
160 const char ZE_SELECTION_ARR_NAME[] = "@@ZeSelection@@";
161 vtkDataSet *output(input->NewInstance());
162 output->ShallowCopy(input);
163 vtkSmartPointer<vtkThreshold> thres(vtkSmartPointer<vtkThreshold>::New());
164 thres->SetInputData(output);
165 vtkDataSetAttributes *dscIn(input->GetCellData()), *dscIn2(input->GetPointData());
166 vtkDataSetAttributes *dscOut(output->GetCellData()), *dscOut2(output->GetPointData());
168 double vMin(1.), vMax(2.);
169 thres->ThresholdBetween(vMin, vMax);
171 vtkIdType nbOfCells(input->GetNumberOfCells());
172 vtkCharArray *zeSelection(vtkCharArray::New());
173 zeSelection->SetName(ZE_SELECTION_ARR_NAME);
174 zeSelection->SetNumberOfComponents(1);
175 char *pt(new char[nbOfCells]);
176 zeSelection->SetArray(pt, nbOfCells, 0, VTK_DATA_ARRAY_DELETE);
177 std::fill(pt, pt + nbOfCells, 0);
178 std::vector<bool> pt2(nbOfCells, false);
179 for (std::vector<int>::const_iterator it = idsToKeep.begin(); it != idsToKeep.end(); it++)
181 for (vtkIdType ii = 0; ii < nbOfCells; ii++)
183 if (input->GetCellType(ii) == *it)
187 for (int ii = 0; ii < nbOfCells; ii++)
190 int idx(output->GetCellData()->AddArray(zeSelection));
191 output->GetCellData()->SetActiveAttribute(idx, vtkDataSetAttributes::SCALARS);
192 output->GetCellData()->CopyScalarsOff();
193 zeSelection->Delete();
195 thres->SetInputArrayToProcess(idx, 0, 0, "vtkDataObject::FIELD_ASSOCIATION_CELLS", ZE_SELECTION_ARR_NAME);
197 vtkUnstructuredGrid *zeComputedOutput(thres->GetOutput());
198 zeComputedOutput->GetCellData()->RemoveArray(idx);
200 zeComputedOutput->Register(0);
201 return zeComputedOutput;
204 int vtkExtractThreeD::RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
208 //std::cerr << "########################################## vtkExtractThreeD::RequestData ##########################################" << std::endl;
209 vtkInformation *inputInfo = inputVector[0]->GetInformationObject(0);
210 vtkDataSet *input(vtkDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
211 vtkInformation *info(input->GetInformation());
212 vtkInformation *outInfo(outputVector->GetInformationObject(0));
213 vtkDataSet *output(vtkDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
214 std::vector<int> idsToKeep(this->Internal->getIdsToKeep());
215 vtkDataSet *tryOnCell(FilterFamilies(input, idsToKeep));
216 // first shrink the input
217 output->ShallowCopy(tryOnCell);
220 catch (INTERP_KERNEL::Exception &e)
222 std::cerr << "Exception has been thrown in vtkExtractThreeD::RequestData : " << e.what() << std::endl;
228 void vtkExtractThreeD::PrintSelf(ostream &os, vtkIndent indent)
230 this->Superclass::PrintSelf(os, indent);