5#include <vtkAbstractArray.h>
6#include <vtkCellData.h>
7#include <vtkCellType.h>
8#include <vtkDoubleArray.h>
9#include <vtkInformation.h>
11#include <vtkObjectFactory.h>
12#include <vtkPointData.h>
15#include <vtkUnstructuredGrid.h>
23 this->SetNumberOfInputPorts(1);
24 this->SetNumberOfOutputPorts(1);
29 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(),
"vtkTable");
37 info->Set(vtkDataObject::DATA_TYPE_NAME(),
"vtkUnstructuredGrid");
44 vtkInformationVector **inputVector,
45 vtkInformationVector *outputVector) {
50 auto *input = vtkTable::GetData(inputVector[0]);
51 auto *output = vtkUnstructuredGrid::GetData(outputVector);
53 if(SelectFieldsWithRegexp) {
56 const auto n = input->GetNumberOfColumns();
57 for(
int i = 0; i < n; ++i) {
58 const auto &name = input->GetColumnName(i);
59 if(std::regex_match(name, std::regex(RegexpString))) {
60 ScalarFields.emplace_back(name);
65 const SimplexId numberOfRows = input->GetNumberOfRows();
66 const SimplexId numberOfColumns = ScalarFields.size();
68 if(numberOfRows <= 0 || numberOfColumns <= 0) {
69 this->
printErr(
"Input matrix has invalid dimensions (rows: "
70 + std::to_string(numberOfRows)
71 +
", columns: " + std::to_string(numberOfColumns) +
")");
75 if(numberOfColumns != numberOfRows) {
76 this->
printErr(
"Input distance matrix is not square (rows: "
77 + std::to_string(numberOfRows)
78 +
", columns: " + std::to_string(numberOfColumns) +
")");
82 std::array<vtkAbstractArray *, 3> components{
83 input->GetColumnByName(this->XColumn.data()),
84 input->GetColumnByName(this->YColumn.data()),
85 input->GetColumnByName(this->ZColumn.data()),
89 vtkNew<vtkPoints> points{};
90 const auto nPoints = components[0]->GetNumberOfTuples();
91 points->SetNumberOfPoints(nPoints);
92#ifdef TTK_ENABLE_OPENMP
93#pragma omp parallel for num_threads(this->threadNumber_)
95 for(vtkIdType i = 0; i < nPoints; ++i) {
96 std::array<float, 3> coords{};
97 if(components[0] !=
nullptr) {
98 coords[0] = components[0]->GetVariantValue(i).ToFloat();
100 if(components[1] !=
nullptr && components[1] != components[0]) {
101 coords[1] = components[1]->GetVariantValue(i).ToFloat();
103 if(components[2] !=
nullptr && components[2] != components[0]
104 && components[2] != components[1]) {
105 coords[2] = components[2]->GetVariantValue(i).ToFloat();
107 points->SetPoint(i, coords.data());
110 std::vector<vtkAbstractArray *> arrays{};
111 for(
const auto &s : ScalarFields) {
112 arrays.push_back(input->GetColumnByName(s.data()));
115 std::vector<std::vector<double>> inputMatrix(numberOfRows);
116#ifdef TTK_ENABLE_OPENMP
117#pragma omp parallel for num_threads(this->threadNumber_)
119 for(SimplexId i = 0; i < numberOfRows; ++i) {
120 for(
size_t j = 0; j < arrays.size(); ++j) {
121 inputMatrix[i].emplace_back(arrays[j]->GetVariantValue(i).ToDouble());
125 std::vector<ttk::SimplexId> vec_connectivity{};
126 std::vector<double> diameters{};
128 vtkNew<vtkDoubleArray> diamMin{}, diamMean{}, diamMax{};
129 diamMin->SetName(
"MinDiameter");
130 diamMean->SetName(
"MeanDiameter");
131 diamMax->SetName(
"MaxDiameter");
132 diamMin->SetNumberOfTuples(numberOfRows);
133 diamMean->SetNumberOfTuples(numberOfRows);
134 diamMax->SetNumberOfTuples(numberOfRows);
135 vtkNew<vtkDoubleArray> gaussianDensity{};
136 gaussianDensity->SetName(
"GaussianDensity");
137 gaussianDensity->SetNumberOfTuples(numberOfRows);
140 = this->
execute(vec_connectivity, diameters,
141 std::array<double *const, 3>{
142 ttkUtils::GetPointer<double>(diamMin),
143 ttkUtils::GetPointer<double>(diamMean),
144 ttkUtils::GetPointer<double>(diamMax),
146 inputMatrix, ttkUtils::GetPointer<double>(gaussianDensity));
152 const auto nCells = vec_connectivity.size() / (this->
OutputDimension + 1);
153 vtkNew<ttkSimplexIdTypeArray> offsets{}, connectivity{};
154 offsets->SetNumberOfComponents(1);
155 offsets->SetNumberOfTuples(nCells + 1);
156 connectivity->SetNumberOfComponents(1);
157 connectivity->SetNumberOfTuples(vec_connectivity.size());
159#ifdef TTK_ENABLE_OPENMP
160#pragma omp parallel for num_threads(this->threadNumber_)
162 for(
size_t i = 0; i < nCells + 1; ++i) {
166#ifdef TTK_ENABLE_OPENMP
167#pragma omp parallel for num_threads(this->threadNumber_)
169 for(
size_t i = 0; i < vec_connectivity.size(); ++i) {
170 connectivity->SetTuple1(i, vec_connectivity[i]);
174 vtkNew<vtkCellArray> cells{};
175#ifndef TTK_ENABLE_64BIT_IDS
176 cells->Use32BitStorage();
178 cells->SetData(offsets, connectivity);
180 const auto getCellType = [](
const SimplexId dim) ->
int {
183 }
else if(dim == 2) {
185 }
else if(dim == 1) {
191 output->SetPoints(points);
195 vtkNew<vtkDoubleArray> cellDiameters{};
196 cellDiameters->SetNumberOfTuples(nCells);
197 cellDiameters->SetName(
"Diameter");
199#ifdef TTK_ENABLE_OPENMP
200#pragma omp parallel for num_threads(this->threadNumber_)
202 for(
size_t i = 0; i < nCells; ++i) {
203 cellDiameters->SetTuple1(i, diameters[i]);
206 output->GetPointData()->AddArray(diamMin);
207 output->GetPointData()->AddArray(diamMean);
208 output->GetPointData()->AddArray(diamMax);
210 output->GetPointData()->AddArray(gaussianDensity);
212 output->GetCellData()->AddArray(cellDiameters);
214 if(this->KeepAllDataArrays) {
215 auto pd = output->GetPointData();
217 for(vtkIdType i = 0; i < input->GetNumberOfColumns(); ++i) {
218 pd->AddArray(input->GetColumn(i));
223 this->
printMsg(
"Produced " + std::to_string(nCells) +
" cells of dimension "
225 1.0, tm.getElapsedTime(), this->threadNumber_);
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
TTK VTK-filter that wraps the ttk::RipsComplex processing package.
int FillInputPortInformation(int port, vtkInformation *info) override
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
int FillOutputPortInformation(int port, vtkInformation *info) override
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
bool ComputeGaussianDensity
int execute(std::vector< SimplexId > &connectivity, std::vector< double > &diameters, std::array< double *const, 3 > diamStats, const std::vector< std::vector< double > > &distanceMatrix, double *const density=nullptr) const
Main entry point.
int SimplexId
Identifier type for simplices of any dimension.
vtkStandardNewMacro(ttkRipsComplex)