TTK
Loading...
Searching...
No Matches
ttkPersistenceCurve.cpp
Go to the documentation of this file.
1#include <ttkMacros.h>
4#include <ttkUtils.h>
5
6#include <vtkDoubleArray.h>
7#include <vtkUnstructuredGrid.h>
8
10
12 this->SetNumberOfInputPorts(1);
13 this->SetNumberOfOutputPorts(4);
14}
15
17 return this->GetOutput(0);
18}
19
20vtkTable *ttkPersistenceCurve::GetOutput(int port) {
21 return vtkTable::SafeDownCast(this->GetOutputDataObject(port));
22}
23
25 vtkInformation *info) {
26 if(port == 0) {
27 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkUnstructuredGrid");
28 return 1;
29 }
30 return 0;
31}
32
34 vtkInformation *info) {
35 if(port == 0 || port == 1 || port == 2 || port == 3) {
36 info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkTable");
37 return 1;
38 }
39 return 0;
40}
41
42static void getPersistenceCurve(vtkTable *outputCurve,
44 const std::string &pairsType) {
45
46 if(plot.empty()) {
47 return;
48 }
49
50 vtkNew<vtkDoubleArray> persistenceScalars{};
51 const auto persQualifier = "Persistence (" + pairsType + ")";
52 persistenceScalars->SetName(persQualifier.c_str());
53 persistenceScalars->SetNumberOfTuples(plot.size());
54
55 vtkNew<ttkSimplexIdTypeArray> numberOfPairsScalars{};
56 const auto nPairsQualifier = "Number Of Pairs (" + pairsType + ")";
57 numberOfPairsScalars->SetName(nPairsQualifier.c_str());
58 numberOfPairsScalars->SetNumberOfTuples(plot.size());
59
60 for(size_t i = 0; i < plot.size(); ++i) {
61 persistenceScalars->SetTuple1(i, plot[i].first);
62 numberOfPairsScalars->SetTuple1(i, plot[i].second);
63 }
64
65 vtkNew<vtkTable> persistenceCurve{};
66 persistenceCurve->AddColumn(persistenceScalars);
67 persistenceCurve->AddColumn(numberOfPairsScalars);
68
69 outputCurve->ShallowCopy(persistenceCurve);
70}
71
72int ttkPersistenceCurve::RequestData(vtkInformation *ttkNotUsed(request),
73 vtkInformationVector **inputVector,
74 vtkInformationVector *outputVector) {
75
76 ttk::Timer timer;
77
78 auto *input = vtkUnstructuredGrid::GetData(inputVector[0]);
79
80 auto *minSadTable = vtkTable::GetData(outputVector, 0);
81 auto *sadSadTable = vtkTable::GetData(outputVector, 1);
82 auto *sadMaxTable = vtkTable::GetData(outputVector, 2);
83 auto *allPairsTable = vtkTable::GetData(outputVector, 3);
84
85 ttk::DiagramType diagram{};
86 int ret = VTUToDiagram(diagram, input, *this);
87 if(ret != 0) {
88 this->printErr("Could not read input Persistence Diagram");
89 return 0;
90 }
91
92 std::array<PlotType, 4> plots{};
93 ret = this->execute(plots, diagram);
94 if(ret != 0) {
95 this->printErr("Could not execute base code");
96 return -1;
97 }
98
99 getPersistenceCurve(minSadTable, plots[0], "minimum-saddle pairs");
100 getPersistenceCurve(sadSadTable, plots[1], "saddle-saddle pairs");
101 getPersistenceCurve(sadMaxTable, plots[2], "maximum-saddle pairs");
102 getPersistenceCurve(allPairsTable, plots[3], "all pairs");
103
104 this->printMsg("Completed", 1, timer.getElapsedTime(), threadNumber_);
106
107 return 1;
108}
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
Definition BaseClass.h:47
TTK VTK-filter for the computation of persistence curves.
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
int FillOutputPortInformation(int port, vtkInformation *info) override
int FillInputPortInformation(int port, vtkInformation *info) override
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:149
int execute(std::array< PlotType, 4 > &plots, const DiagramType &diagram) const
Compute the Persistence Curve from the input Diagram.
std::vector< std::pair< double, SimplexId > > PlotType
double getElapsedTime()
Definition Timer.h:15
std::vector< PersistencePair > DiagramType
Persistence Diagram type as a vector of Persistence pairs.
vtkStandardNewMacro(ttkPersistenceCurve)
int VTUToDiagram(ttk::DiagramType &diagram, vtkUnstructuredGrid *vtu, const ttk::Debug &dbg)
Converts a Persistence Diagram in the VTK Unstructured Grid format (as generated by the ttkPersistenc...
printMsg(debug::output::BOLD+" | | | | | . \\ | | (__| | / __/| |_| / __/|__ _|"+debug::output::ENDCOLOR, debug::Priority::PERFORMANCE, debug::LineMode::NEW, stream)