1#include <vtkCellData.h>
2#include <vtkDataArray.h>
4#include <vtkInformation.h>
6#include <vtkObjectFactory.h>
7#include <vtkPointData.h>
8#include <vtkSignedCharArray.h>
9#include <vtkSmartPointer.h>
10#include <vtkUnstructuredGrid.h>
21 this->SetNumberOfInputPorts(1);
22 this->SetNumberOfOutputPorts(1);
23 ForceInputOffsetScalarField =
false;
28 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(),
"vtkDataSet");
35 if(port == 0 || port == 1) {
36 info->Set(vtkDataObject::DATA_TYPE_NAME(),
"vtkUnstructuredGrid");
42template <
class dataTypeU,
class dataTypeV>
44 const dataTypeV *
const vField,
48 this->execute(jacobiSet_, uField, vField,
49 *
static_cast<TTK_TT *
>(triangulation->
getData()),
55 vtkInformationVector **inputVector,
56 vtkInformationVector *outputVector) {
58 const auto input = vtkDataSet::GetData(inputVector[0]);
59 auto output = vtkUnstructuredGrid::GetData(outputVector);
61 const auto uComponent = this->GetInputArrayToProcess(0, input);
62 const auto vComponent = this->GetInputArrayToProcess(1, input);
64 if(uComponent ==
nullptr || vComponent ==
nullptr)
67 this->
printMsg(
"U-component: `" + std::string{uComponent->GetName()} +
"'");
68 this->
printMsg(
"V-component: `" + std::string{vComponent->GetName()} +
"'");
71 const auto offsetFieldU
72 = this->
GetOrderArray(input, 0, 2, ForceInputOffsetScalarField);
73 const auto offsetFieldV
74 = this->
GetOrderArray(input, 1, 3, ForceInputOffsetScalarField);
82 if(triangulation ==
nullptr)
86#ifndef TTK_ENABLE_DOUBLE_TEMPLATING
87 if(uComponent->GetDataType() != vComponent->GetDataType()) {
89 "Scalar fields should have same input type. Use TTKPointDataConverter or "
90 "TTKArrayEditor to convert array types.");
93 switch(uComponent->GetDataType()) {
100 switch(vtkTemplate2PackMacro(
101 uComponent->GetDataType(), vComponent->GetDataType())) {
109 vtkNew<vtkSignedCharArray> edgeTypes{};
111 edgeTypes->SetNumberOfComponents(1);
112 edgeTypes->SetNumberOfTuples(2 * jacobiSet_.size());
113 edgeTypes->SetName(
"Critical Type");
115 vtkNew<vtkSignedCharArray> isPareto{};
116 isPareto->SetNumberOfComponents(1);
117 isPareto->SetNumberOfTuples(2 * jacobiSet_.size());
118 isPareto->SetName(
"IsPareto");
120 vtkNew<vtkPoints> pointSet{};
121 pointSet->SetNumberOfPoints(2 * jacobiSet_.size());
123 vtkNew<vtkCellArray> cellArray{};
124 vtkNew<vtkIdList> idList{};
125 idList->SetNumberOfIds(2);
127 size_t pointCount = 0;
128 std::array<double, 3> p{};
129 for(
size_t i = 0; i < jacobiSet_.size(); i++) {
131 int edgeId = jacobiSet_[i].first;
133 triangulation->getEdgeVertex(edgeId, 0, vertexId0);
134 triangulation->getEdgeVertex(edgeId, 1, vertexId1);
136 input->GetPoint(vertexId0, p.data());
137 pointSet->SetPoint(pointCount, p.data());
138 edgeTypes->SetTuple1(pointCount, (
float)jacobiSet_[i].second);
139 isPareto->SetTuple1(pointCount, (
float)isPareto_[i]);
141 idList->SetId(0, pointCount);
144 input->GetPoint(vertexId1, p.data());
145 pointSet->SetPoint(pointCount, p.data());
146 edgeTypes->SetTuple1(pointCount, (
float)jacobiSet_[i].second);
147 isPareto->SetTuple1(pointCount, (
float)isPareto_[i]);
148 idList->SetId(1, pointCount);
151 cellArray->InsertNextCell(idList);
153 output->SetPoints(pointSet);
154 output->SetCells(VTK_LINE, cellArray);
155 output->GetPointData()->AddArray(edgeTypes);
156 output->GetPointData()->AddArray(isPareto);
159 vtkNew<ttkSimplexIdTypeArray> edgeIdArray{};
160 edgeIdArray->SetNumberOfComponents(1);
161 edgeIdArray->SetNumberOfTuples(jacobiSet_.size());
162 edgeIdArray->SetName(
"EdgeIds");
165 for(
size_t i = 0; i < jacobiSet_.size(); i++) {
166 edgeIdArray->SetTuple1(pointCount, (
float)jacobiSet_[i].first);
170 output->GetCellData()->AddArray(edgeIdArray);
172 output->GetCellData()->RemoveArray(
"EdgeIds");
177 for(
int i = 0; i < input->GetPointData()->GetNumberOfArrays(); i++) {
179 const auto scalarField = input->GetPointData()->GetArray(i);
182 scalarArray->SetNumberOfComponents(scalarField->GetNumberOfComponents());
183 scalarArray->SetNumberOfTuples(2 * jacobiSet_.size());
184 scalarArray->SetName(scalarField->GetName());
185 std::vector<double> value(scalarField->GetNumberOfComponents());
187 for(
size_t j = 0; j < jacobiSet_.size(); j++) {
188 int edgeId = jacobiSet_[j].first;
190 triangulation->getEdgeVertex(edgeId, 0, vertexId0);
191 triangulation->getEdgeVertex(edgeId, 1, vertexId1);
193 scalarField->GetTuple(vertexId0, value.data());
194 scalarArray->SetTuple(2 * j, value.data());
196 scalarField->GetTuple(vertexId1, value.data());
197 scalarArray->SetTuple(2 * j + 1, value.data());
199 output->GetPointData()->AddArray(scalarArray);
202 for(
int i = 0; i < input->GetPointData()->GetNumberOfArrays(); i++) {
203 output->GetPointData()->RemoveArray(
204 input->GetPointData()->GetArray(i)->GetName());
#define ttkTemplateMacro(triangulationType, call)
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
vtkDataArray * GetOrderArray(vtkDataSet *const inputData, const int scalarArrayIdx, const int orderArrayIdx=0, const bool enforceOrderArrayIdx=false)
ttk::Triangulation * GetTriangulation(vtkDataSet *dataSet)
TTK VTK-filter that computes the Jacobi set of a bivariate volumetric data-set.
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
int FillOutputPortInformation(int port, vtkInformation *info) override
int FillInputPortInformation(int port, vtkInformation *info) override
int dispatch(const dataTypeU *const uField, const dataTypeV *const vField, ttk::Triangulation *const triangulation)
static void * GetVoidPointer(vtkDataArray *array, vtkIdType start=0)
int printMsg(const std::string &msg, const debug::Priority &priority=debug::Priority::INFO, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cout) const
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
void setSosOffsetsV(const SimplexId *const sosOffsets)
void preconditionTriangulation(AbstractTriangulation *const triangulation)
void setSosOffsetsU(const SimplexId *const sosOffsets)
Triangulation is a class that provides time and memory efficient traversal methods on triangulations ...
AbstractTriangulation * getData()
Triangulation::Type getType() const
int SimplexId
Identifier type for simplices of any dimension.
vtkStandardNewMacro(ttkJacobiSet)