TTK
Loading...
Searching...
No Matches
ttkQuadrangulationSubdivision.cpp
Go to the documentation of this file.
1#include <ttkMacros.h>
3#include <ttkUtils.h>
4
5#include <vtkCellData.h>
6#include <vtkFloatArray.h>
7#include <vtkIdTypeArray.h>
8#include <vtkInformation.h>
9#include <vtkPointData.h>
10#include <vtkPolyData.h>
11
13
15 // MSC quadrangulation + initial 2D mesh
16 SetNumberOfInputPorts(2);
17 // quad mesh (containing ttkVertexIdentifiers of critical points)
18 SetNumberOfOutputPorts(1);
19}
20
22 int port, vtkInformation *info) {
23 if(port == 0) { // input quadrangulation
24 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkPolyData");
25 return 1;
26 } else if(port == 1) { // triangulated domain
27 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
28 return 1;
29 }
30 return 0;
31}
32
34 int port, vtkInformation *info) {
35 if(port == 0) {
36 info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkPolyData");
37 return 1;
38 }
39 return 0;
40}
41
43 vtkInformation *ttkNotUsed(request),
44 vtkInformationVector **inputVector,
45 vtkInformationVector *outputVector) {
46
47 auto quads = vtkPolyData::GetData(inputVector[0]);
48 auto mesh = vtkDataSet::GetData(inputVector[1]);
49 auto output = vtkPolyData::GetData(outputVector);
50
51 auto triangulation = ttkAlgorithm::GetTriangulation(mesh);
52 if(triangulation == nullptr) {
53 return 0;
54 }
55 this->preconditionTriangulation(triangulation);
56
57 auto inputCells = quads->GetPolys();
58 if(inputCells == nullptr || inputCells->GetData() == nullptr) {
59 this->printErr("Invalid input quadrangle cells");
60 return 0;
61 }
62
63 auto inputPoints = quads->GetPoints();
64 auto pointData = quads->GetPointData();
65 if(inputPoints == nullptr || inputPoints->GetData() == nullptr
66 || pointData == nullptr) {
67 this->printErr("Invalid input quadrangle points");
68 return 0;
69 }
70
71 auto identifiers = pointData->GetArray(
72 static_cast<const char *>(ttk::VertexScalarFieldName));
73 if(identifiers == nullptr) {
74 this->printErr("Missing point data array named "
75 + std::string(ttk::VertexScalarFieldName));
76 return 0;
77 }
78
79 this->setInputQuads(
80 // get quads from PolyData's connectivity array
81 ttkUtils::GetVoidPointer(inputCells->GetConnectivityArray()),
82 inputCells->GetNumberOfCells());
83 this->setInputVertices(ttkUtils::GetVoidPointer(inputPoints->GetData()),
84 inputPoints->GetNumberOfPoints());
86 ttkUtils::GetVoidPointer(identifiers), identifiers->GetNumberOfTuples());
87
88 int res{-1};
90 triangulation->getType(),
91 res = this->execute(*static_cast<TTK_TT *>(triangulation->getData())));
92
93 if(res != 0) {
94 this->printWrn("Please increase the number of relaxation iterations, of "
95 "subdivision levels or consider another function (higher "
96 "eigenfunctions).");
97 if(!ShowResError) {
98 return 0;
99 }
100 }
101
102 vtkNew<vtkCellArray> cells{};
103
104 for(size_t i = 0; i < outputQuads_.size(); i++) {
105 cells->InsertNextCell(4, this->outputQuads_[i].data());
106 }
107
108 // update output: get quadrangle values
109 output->SetPolys(cells);
110
111 vtkNew<vtkPoints> points{};
112 for(size_t i = 0; i < outputPoints_.size(); ++i) {
113 points->InsertNextPoint(outputPoints_[i].data());
114 }
115
116 // update output: get quadrangle vertices
117 output->SetPoints(points);
118
119 // add data array of points valences
120 vtkNew<ttkSimplexIdTypeArray> valences{};
121 valences->SetName("Valence");
123 valences, outputValences_.data(), outputValences_.size(), 1);
124 output->GetPointData()->AddArray(valences);
125
126 vtkNew<vtkFloatArray> density{};
127 density->SetName("Density");
129 density, outputDensity_.data(), outputDensity_.size(), 1);
130 output->GetPointData()->AddArray(density);
131
132 vtkNew<vtkFloatArray> deformity{};
133 deformity->SetName("Deformity");
135 deformity, outputDifformity_.data(), outputDifformity_.size(), 1);
136 output->GetPointData()->AddArray(deformity);
137
138 // add data array of points infos
139 vtkNew<ttkSimplexIdTypeArray> infos{};
140 infos->SetName("Type");
142 infos, outputVertType_.data(), outputVertType_.size(), 1);
143 output->GetPointData()->AddArray(infos);
144
145 vtkNew<ttkSimplexIdTypeArray> subd{};
146 subd->SetName("Subdivision");
148 subd, outputSubdivision_.data(), outputSubdivision_.size(), 1);
149 output->GetPointData()->AddArray(subd);
150
151 vtkNew<ttkSimplexIdTypeArray> nearestVert{};
152 nearestVert->SetName(ttk::VertexScalarFieldName);
154 nearestVertexIdentifier_.size(), 1);
155 output->GetPointData()->AddArray(nearestVert);
156
157 if(QuadStatistics) {
158 vtkNew<vtkFloatArray> quadArea{};
159 quadArea->SetName("Quad Area");
160 ttkUtils::SetVoidArray(quadArea, quadArea_.data(), quadArea_.size(), 1);
161 output->GetCellData()->AddArray(quadArea);
162
163 vtkNew<vtkFloatArray> diagsRatio{};
164 diagsRatio->SetName("Diagonals Ratio");
166 diagsRatio, quadDiagsRatio_.data(), quadDiagsRatio_.size(), 1);
167 output->GetCellData()->AddArray(diagsRatio);
168
169 vtkNew<vtkFloatArray> edgesRatio{};
170 edgesRatio->SetName("Edges Ratio");
172 edgesRatio, quadEdgesRatio_.data(), quadEdgesRatio_.size(), 1);
173 output->GetCellData()->AddArray(edgesRatio);
174
175 vtkNew<vtkFloatArray> anglesRatio{};
176 anglesRatio->SetName("Angles Ratio");
178 anglesRatio, quadAnglesRatio_.data(), quadAnglesRatio_.size(), 1);
179 output->GetCellData()->AddArray(anglesRatio);
180
181 vtkNew<vtkFloatArray> hausDist{};
182 hausDist->SetName("Hausdorff");
183 ttkUtils::SetVoidArray(hausDist, hausdorff_.data(), hausdorff_.size(), 1);
184 output->GetPointData()->AddArray(hausDist);
185 }
186
187 // shallow copy input field data
188 output->GetFieldData()->ShallowCopy(mesh->GetFieldData());
189
190 return 1;
191}
#define ttkTemplateMacro(triangulationType, call)
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
Definition BaseClass.h:47
ttk::Triangulation * GetTriangulation(vtkDataSet *dataSet)
TTK VTK-filter for surface quadrangulation.
int FillInputPortInformation(int port, vtkInformation *info) override
int FillOutputPortInformation(int port, vtkInformation *info) override
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
static void * GetVoidPointer(vtkDataArray *array, vtkIdType start=0)
Definition ttkUtils.cpp:226
static void SetVoidArray(vtkDataArray *array, void *data, vtkIdType size, int save)
Definition ttkUtils.cpp:280
int printWrn(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:159
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:149
void setInputVertices(void *const address, unsigned int size)
std::vector< SimplexId > nearestVertexIdentifier_
void preconditionTriangulation(AbstractTriangulation *const triangl)
void setInputVertexIdentifiers(void *const address, unsigned int size)
void setInputQuads(void *const address, unsigned int size)
const char VertexScalarFieldName[]
default name for vertex scalar field
Definition DataTypes.h:35
vtkStandardNewMacro(ttkQuadrangulationSubdivision)