TTK
Loading...
Searching...
No Matches
ttkPersistenceDiagramUtils.cpp
Go to the documentation of this file.
1#include <MPIUtils.h>
2#include <ttkMacros.h>
3#include <ttkPersistenceDiagramUtils.h>
4#include <ttkUtils.h>
5
6#include <vtkCellData.h>
7#include <vtkDoubleArray.h>
8#include <vtkFloatArray.h>
9#include <vtkPointData.h>
10#include <vtkThreshold.h>
11#include <vtkTransform.h>
12#include <vtkTransformFilter.h>
13#include <vtkUnsignedCharArray.h>
14#include <vtkUnstructuredGrid.h>
15#include <vtkVersionMacros.h> // for VTK_VERSION_CHECK
16
17int VTUToDiagram(ttk::DiagramType &diagram,
18 vtkUnstructuredGrid *vtu,
19 const ttk::Debug &dbg) {
20
21 const auto pd = vtu->GetPointData();
22 const auto cd = vtu->GetCellData();
23 const auto points = vtu->GetPoints();
24
25 if(pd == nullptr) {
26 dbg.printWrn("VTU diagram with NULL Point Data");
27 return -1;
28 }
29 if(cd == nullptr) {
30 dbg.printWrn("VTU diagram with NULL Cell Data");
31 return -2;
32 }
33 if(points == nullptr) {
34 dbg.printWrn("VTU with no points");
35 return -3;
36 }
37
38 // cell data
39 const auto pairId = ttkSimplexIdTypeArray::SafeDownCast(
40 cd->GetArray(ttk::PersistencePairIdentifierName));
41 const auto pairType
42 = vtkIntArray::SafeDownCast(cd->GetArray(ttk::PersistencePairTypeName));
43 const auto pairPers = cd->GetArray(ttk::PersistenceName);
44 const auto birthScalars = cd->GetArray(ttk::PersistenceBirthName);
45 const auto isFinite = cd->GetArray(ttk::PersistenceIsFinite);
46
47 // point data
48 const auto vertexId = ttkSimplexIdTypeArray::SafeDownCast(
49 pd->GetArray(ttk::VertexScalarFieldName));
50 const auto critType
51 = vtkIntArray::SafeDownCast(pd->GetArray(ttk::PersistenceCriticalTypeName));
52 const auto coords = vtkFloatArray::SafeDownCast(
53 pd->GetArray(ttk::PersistenceCoordinatesName));
54
55 const bool embed = coords == nullptr;
56
57 if(pairId == nullptr) {
58 dbg.printWrn("Missing PairIdentifier cell data array");
59 return -5;
60 }
61 if(pairType == nullptr) {
62 dbg.printWrn("Missing PairType cell data array");
63 return -6;
64 }
65 if(pairPers == nullptr) {
66 dbg.printWrn("Missing Persistence cell data array");
67 return -7;
68 }
69 if(vertexId == nullptr) {
70 dbg.printWrn("Missing ttkVertexScalarField point data array");
71 return -8;
72 }
73 if(critType == nullptr) {
74 dbg.printWrn("Missing CriticalType point data array");
75 return -9;
76 }
77 if(birthScalars == nullptr) {
78 dbg.printWrn("Missing Birth cell data array");
79 return -10;
80 }
81 if(isFinite == nullptr) {
82 dbg.printWrn("Missing IsFinite cell data array");
83 return -12;
84 }
85
86 int nPairs = pairId->GetNumberOfTuples();
87
88 // compact pairIds in [0, nPairs - 1] (diagonal excepted)
89 for(int i = 0; i < nPairs; i++) {
90 if(pairId->GetTuple1(i) != -1) {
91 pairId->SetTuple1(i, i);
92 } else {
93 // detect diagram diagonal
94 nPairs -= 1;
95 }
96 }
97
98 if(nPairs < 1) {
99 dbg.printWrn("Diagram has no pairs");
100 return -4;
101 }
102
103 diagram.resize(nPairs);
104
105 // skip diagonal cell (corresponding points already dealt with)
106 for(int i = 0; i < nPairs; ++i) {
107
108 const auto pId = pairId->GetValue(i);
109 // skip diagram diagonal if present
110 if(pId == -1) {
111 continue;
112 }
113
114 const auto i0 = 2 * i + 0;
115 const auto i1 = 2 * i + 1;
116
117 const auto v0 = vertexId->GetValue(i0);
118 const auto v1 = vertexId->GetValue(i1);
119 const auto ct0 = static_cast<ttk::CriticalType>(critType->GetValue(i0));
120 const auto ct1 = static_cast<ttk::CriticalType>(critType->GetValue(i1));
121
122 const auto pType = pairType->GetValue(i);
123 const auto pers = pairPers->GetTuple1(i);
124 const auto birth = birthScalars->GetTuple1(i);
125 const auto isFin = static_cast<bool>(isFinite->GetTuple1(i));
126
127 std::array<float, 3> coordsBirth{}, coordsDeath{};
128 // no vtkPoints::GetPoint() taking a float array, have to do the
129 // conversion by hand...
130 std::array<double, 3> tmp{};
131
132 if(embed) {
133 points->GetPoint(i0, tmp.data());
134 coordsBirth = {static_cast<float>(tmp[0]), static_cast<float>(tmp[1]),
135 static_cast<float>(tmp[2])};
136 points->GetPoint(i1, tmp.data());
137 coordsDeath = {static_cast<float>(tmp[0]), static_cast<float>(tmp[1]),
138 static_cast<float>(tmp[2])};
139 } else {
140 coords->GetTuple(i0, tmp.data());
141 coordsBirth = {static_cast<float>(tmp[0]), static_cast<float>(tmp[1]),
142 static_cast<float>(tmp[2])};
143 coords->GetTuple(i1, tmp.data());
144 coordsDeath = {static_cast<float>(tmp[0]), static_cast<float>(tmp[1]),
145 static_cast<float>(tmp[2])};
146 }
147
148 // put pairs in diagram
149 diagram[i] = ttk::PersistencePair{
150 ttk::CriticalVertex{v0, birth, {}, coordsBirth, ct0},
151 ttk::CriticalVertex{v1, birth + pers, {}, coordsDeath, ct1}, pType,
152 isFin};
153 }
154
155 return 0;
156}
157
158int DiagramToVTU(vtkUnstructuredGrid *vtu,
159 const ttk::DiagramType &diagram,
160 vtkDataArray *const inputScalars,
161 const ttk::Debug &dbg,
162 const int dim,
163 const bool embedInDomain) {
164
165 if(diagram.empty()) {
166 dbg.printWrn("Empty diagram");
167 return -1;
168 }
169
170 const auto pd = vtu->GetPointData();
171 const auto cd = vtu->GetCellData();
172
173 if(pd == nullptr || cd == nullptr) {
174 dbg.printWrn("Grid has no point data or no cell data");
175 return -2;
176 }
177
178 // point data arrays
179
180 vtkNew<ttkSimplexIdTypeArray> vertsId{};
181 vertsId->SetName(ttk::VertexScalarFieldName);
182 vertsId->SetNumberOfTuples(2 * diagram.size());
183 pd->AddArray(vertsId);
184
185 vtkNew<vtkIntArray> critType{};
186 critType->SetName(ttk::PersistenceCriticalTypeName);
187 critType->SetNumberOfTuples(2 * diagram.size());
188 pd->AddArray(critType);
189
190 vtkNew<vtkFloatArray> coordsScalars{};
191
192 if(!embedInDomain) {
193 coordsScalars->SetNumberOfComponents(3);
194 coordsScalars->SetName(ttk::PersistenceCoordinatesName);
195 coordsScalars->SetNumberOfTuples(2 * diagram.size());
196 pd->AddArray(coordsScalars);
197 }
198
199 // cell data arrays
200
201 vtkNew<ttkSimplexIdTypeArray> pairsId{};
202 pairsId->SetName(ttk::PersistencePairIdentifierName);
203 pairsId->SetNumberOfTuples(diagram.size());
204 cd->AddArray(pairsId);
205
206 vtkNew<vtkIntArray> pairsDim{};
207 pairsDim->SetName(ttk::PersistencePairTypeName);
208 pairsDim->SetNumberOfTuples(diagram.size());
209 cd->AddArray(pairsDim);
210
211 vtkSmartPointer<vtkDataArray> const persistence{inputScalars->NewInstance()};
212 persistence->SetName(ttk::PersistenceName);
213 persistence->SetNumberOfTuples(diagram.size());
214 cd->AddArray(persistence);
215
216 vtkSmartPointer<vtkDataArray> const birthScalars{inputScalars->NewInstance()};
217 birthScalars->SetName(ttk::PersistenceBirthName);
218 birthScalars->SetNumberOfTuples(diagram.size());
219 cd->AddArray(birthScalars);
220
221 vtkNew<vtkUnsignedCharArray> isFinite{};
222 isFinite->SetName(ttk::PersistenceIsFinite);
223 isFinite->SetNumberOfTuples(diagram.size());
224 cd->AddArray(isFinite);
225
226 // grid
227
228 vtkNew<vtkPoints> points{};
229 points->SetNumberOfPoints(2 * diagram.size());
230 vtkNew<vtkIdTypeArray> offsets{}, connectivity{};
231 offsets->SetNumberOfComponents(1);
232 offsets->SetNumberOfTuples(diagram.size() + 1);
233 connectivity->SetNumberOfComponents(1);
234 connectivity->SetNumberOfTuples(2 * diagram.size());
235
236#ifdef TTK_ENABLE_OPENMP
237#pragma omp parallel for num_threads(dbg.getThreadNumber())
238#endif // TTK_ENABLE_OPENMP
239 for(size_t i = 0; i < diagram.size(); ++i) {
240 const auto &pair{diagram[i]};
241 const auto i0{2 * i + 0}, i1{2 * i + 1};
242 if(embedInDomain) {
243 points->SetPoint(
244 i0, pair.birth.coords[0], pair.birth.coords[1], pair.birth.coords[2]);
245 points->SetPoint(
246 i1, pair.death.coords[0], pair.death.coords[1], pair.death.coords[2]);
247 } else {
248 points->SetPoint(i0, pair.birth.sfValue, pair.birth.sfValue, 0);
249 points->SetPoint(i1, pair.birth.sfValue, pair.death.sfValue, 0);
250 }
251
252 connectivity->SetTuple1(i0, i0);
253 connectivity->SetTuple1(i1, i1);
254 offsets->SetTuple1(i, 2 * i);
255
256 // point data
257 vertsId->SetTuple1(i0, pair.birth.id);
258 vertsId->SetTuple1(i1, pair.death.id);
259 critType->SetTuple1(i0, static_cast<ttk::SimplexId>(pair.birth.type));
260 critType->SetTuple1(i1, static_cast<ttk::SimplexId>(pair.death.type));
261
262 if(!embedInDomain) {
263 coordsScalars->SetTuple3(
264 i0, pair.birth.coords[0], pair.birth.coords[1], pair.birth.coords[2]);
265 coordsScalars->SetTuple3(
266 i1, pair.death.coords[0], pair.death.coords[1], pair.death.coords[2]);
267 }
268
269 // cell data
270 pairsId->SetTuple1(i, i);
271 persistence->SetTuple1(i, pair.persistence());
272 birthScalars->SetTuple1(i, pair.birth.sfValue);
273 isFinite->SetTuple1(i, pair.isFinite);
274 pairsDim->SetTuple1(
275 i, (pair.dim == 2 && pair.isFinite) ? dim - 1 : pair.dim);
276 }
277 offsets->SetTuple1(diagram.size(), connectivity->GetNumberOfTuples());
278
279 vtkNew<vtkCellArray> cells{};
280 cells->SetData(offsets, connectivity);
281 vtu->SetPoints(points);
282 vtu->SetCells(VTK_LINE, cells);
283
284 if(!embedInDomain) {
285 // highest birth (last pair)
286 const auto lastPair = std::max_element(diagram.begin(), diagram.end());
287 // add diagonal (first point -> last birth/penultimate point)
288 std::array<vtkIdType, 2> diag{
289 0, 2 * std::distance(diagram.begin(), lastPair)};
290 vtu->InsertNextCell(VTK_LINE, 2, diag.data());
291 pairsId->InsertTuple1(diagram.size(), -1);
292 pairsDim->InsertTuple1(diagram.size(), -1);
293 isFinite->InsertTuple1(diagram.size(), false);
294 // persistence of global min-max pair
295 const auto maxPersistence = diagram[0].persistence();
296 persistence->InsertTuple1(diagram.size(), 2 * maxPersistence);
297 // birth == death == 0
298 birthScalars->InsertTuple1(diagram.size(), 0);
299 }
300
301 return 0;
302}
303
304#if defined(TTK_ENABLE_MPI) && defined(TTK_ENABLE_OPENMP)
305int DiagramToDistributedVTU(vtkUnstructuredGrid *vtu,
306 const ttk::DiagramType &diagram,
307 vtkDataArray *const inputScalars,
308 const ttk::Debug &dbg,
309 const int dim,
310 const bool embedInDomain) {
311
312 // Computation of order of first local pair
313 const ttk::SimplexId nPairs = diagram.size();
314 ttk::SimplexId beginning{0};
315 MPI_Datatype MPI_SimplexId = ttk::getMPIType(beginning);
316 MPI_Exscan(&nPairs, &beginning, 1, MPI_SimplexId, MPI_SUM, ttk::MPIcomm_);
317
318 if(diagram.empty()) {
319 dbg.printMsg("Empty diagram on this rank");
320 return 0;
321 }
322
323 const auto pd = vtu->GetPointData();
324 const auto cd = vtu->GetCellData();
325
326 // point data arrays
327
328 vtkNew<ttkSimplexIdTypeArray> vertsId{};
329 vertsId->SetName(ttk::VertexScalarFieldName);
330 vertsId->SetNumberOfTuples(2 * diagram.size());
331 pd->AddArray(vertsId);
332
333 vtkNew<vtkIntArray> critType{};
334 critType->SetName(ttk::PersistenceCriticalTypeName);
335 critType->SetNumberOfTuples(2 * diagram.size());
336 pd->AddArray(critType);
337
338 vtkNew<vtkFloatArray> coordsScalars{};
339
340 if(!embedInDomain) {
341 coordsScalars->SetNumberOfComponents(3);
342 coordsScalars->SetName(ttk::PersistenceCoordinatesName);
343 coordsScalars->SetNumberOfTuples(2 * diagram.size());
344 pd->AddArray(coordsScalars);
345 }
346
347 // cell data arrays
348
349 vtkNew<ttkSimplexIdTypeArray> pairsId{};
350 pairsId->SetName(ttk::PersistencePairIdentifierName);
351 pairsId->SetNumberOfTuples(diagram.size());
352 cd->AddArray(pairsId);
353
354 vtkNew<vtkIntArray> pairsDim{};
355 pairsDim->SetName(ttk::PersistencePairTypeName);
356 pairsDim->SetNumberOfTuples(diagram.size());
357 cd->AddArray(pairsDim);
358
359 vtkSmartPointer<vtkDataArray> const persistence{inputScalars->NewInstance()};
360 persistence->SetName(ttk::PersistenceName);
361 persistence->SetNumberOfTuples(diagram.size());
362 cd->AddArray(persistence);
363
364 vtkSmartPointer<vtkDataArray> const birthScalars{inputScalars->NewInstance()};
365 birthScalars->SetName(ttk::PersistenceBirthName);
366 birthScalars->SetNumberOfTuples(diagram.size());
367 cd->AddArray(birthScalars);
368
369 vtkNew<vtkUnsignedCharArray> isFinite{};
370 isFinite->SetName(ttk::PersistenceIsFinite);
371 isFinite->SetNumberOfTuples(diagram.size());
372 cd->AddArray(isFinite);
373
374 // grid
375
376 vtkNew<vtkPoints> points{};
377 points->SetNumberOfPoints(2 * diagram.size());
378 vtkNew<vtkIdTypeArray> offsets{}, connectivity{};
379 offsets->SetNumberOfComponents(1);
380 offsets->SetNumberOfTuples(diagram.size() + 1);
381 connectivity->SetNumberOfComponents(1);
382 connectivity->SetNumberOfTuples(2 * diagram.size());
383
384#pragma omp parallel for num_threads(dbg.getThreadNumber())
385 for(size_t i = 0; i < diagram.size(); ++i) {
386 const auto &pair{diagram[i]};
387 const auto i0{2 * i + 0}, i1{2 * i + 1};
388 if(embedInDomain) {
389 points->SetPoint(
390 i0, pair.birth.coords[0], pair.birth.coords[1], pair.birth.coords[2]);
391 points->SetPoint(
392 i1, pair.death.coords[0], pair.death.coords[1], pair.death.coords[2]);
393 } else {
394 points->SetPoint(i0, pair.birth.sfValue, pair.birth.sfValue, 0);
395 points->SetPoint(i1, pair.birth.sfValue, pair.death.sfValue, 0);
396 }
397
398 connectivity->SetTuple1(i0, i0);
399 connectivity->SetTuple1(i1, i1);
400 offsets->SetTuple1(i, 2 * i);
401
402 // point data
403 vertsId->SetTuple1(i0, pair.birth.id);
404 vertsId->SetTuple1(i1, pair.death.id);
405 critType->SetTuple1(i0, static_cast<ttk::SimplexId>(pair.birth.type));
406 critType->SetTuple1(i1, static_cast<ttk::SimplexId>(pair.death.type));
407
408 if(!embedInDomain) {
409 coordsScalars->SetTuple3(
410 i0, pair.birth.coords[0], pair.birth.coords[1], pair.birth.coords[2]);
411 coordsScalars->SetTuple3(
412 i1, pair.death.coords[0], pair.death.coords[1], pair.death.coords[2]);
413 }
414
415 // cell data
416 pairsId->SetTuple1(i, beginning + i);
417 persistence->SetTuple1(i, pair.persistence());
418 birthScalars->SetTuple1(i, pair.birth.sfValue);
419 isFinite->SetTuple1(i, pair.isFinite);
420 pairsDim->SetTuple1(
421 i, (pair.dim == 2 && pair.isFinite) ? dim - 1 : pair.dim);
422 }
423 offsets->SetTuple1(diagram.size(), connectivity->GetNumberOfTuples());
424
425 vtkNew<vtkCellArray> cells{};
426 cells->SetData(offsets, connectivity);
427 vtu->SetPoints(points);
428 vtu->SetCells(VTK_LINE, cells);
429
430 if(!embedInDomain) {
431 const auto lastPair = std::max_element(diagram.begin(), diagram.end());
432 // add diagonal (first point -> last birth/penultimate point)
433 std::array<vtkIdType, 2> diag{
434 0, 2 * std::distance(diagram.begin(), lastPair)};
435 vtu->InsertNextCell(VTK_LINE, 2, diag.data());
436 pairsId->InsertTuple1(diagram.size(), -1);
437 pairsDim->InsertTuple1(diagram.size(), -1);
438 isFinite->InsertTuple1(diagram.size(), false);
439 const auto maxPersistence = diagram[0].persistence();
440 persistence->InsertTuple1(diagram.size(), 2 * maxPersistence);
441 birthScalars->InsertTuple1(diagram.size(), 0);
442 }
443 return 0;
444}
445#endif
446
447int ProjectDiagramInsideDomain(vtkUnstructuredGrid *const inputDiagram,
448 vtkUnstructuredGrid *const outputDiagram,
449 const ttk::Debug &dbg) {
450
451 ttk::Timer tm{};
452
453 // use vtkThreshold to remove diagonal (PairIdentifier == -1)
454 vtkNew<vtkThreshold> threshold{};
455 threshold->SetInputDataObject(0, inputDiagram);
456 threshold->SetInputArrayToProcess(0, 0, 0,
457 vtkDataObject::FIELD_ASSOCIATION_CELLS,
458 ttk::PersistencePairIdentifierName);
459#if VTK_VERSION_NUMBER < VTK_VERSION_CHECK(9, 2, 0)
460 threshold->ThresholdByUpper(0);
461#else
462 threshold->SetThresholdFunction(vtkThreshold::THRESHOLD_UPPER);
463 threshold->SetUpperThreshold(0);
464#endif
465
466 threshold->Update();
467
468 auto diagonalLess = threshold->GetOutput();
469 auto diagonalLessData = diagonalLess->GetPointData();
470
471 const auto critCoordinates = vtkFloatArray::SafeDownCast(
472 diagonalLessData->GetAbstractArray(ttk::PersistenceCoordinatesName));
473
474 // set new points from Coordinates array
475 vtkNew<vtkFloatArray> coords{};
476 coords->DeepCopy(critCoordinates);
477 coords->SetName("Points");
478 diagonalLess->GetPoints()->SetData(coords);
479 diagonalLessData->RemoveArray(ttk::PersistenceCoordinatesName);
480
481 outputDiagram->ShallowCopy(diagonalLess);
482
483 // don't forget to forward the Field Data
484 outputDiagram->GetFieldData()->ShallowCopy(inputDiagram->GetFieldData());
485
486 dbg.printMsg("Projected Persistence Diagram inside domain", 1.0,
487 tm.getElapsedTime(), dbg.getThreadNumber());
488
489 return 0;
490}
491
492template <typename dataType>
493void getCoords(vtkPoints *points,
494 const dataType *const births,
495 const dataType *const perss,
496 vtkIdType nPoints,
497 const int nThreads) {
498
499#ifdef TTK_ENABLE_OPENMP
500#pragma omp parallel for num_threads(nThreads)
501#endif // TTK_ENABLE_OPENMP
502 for(int i = 0; i < nPoints / 2; ++i) {
503 std::array<float, 3> pt0{
504 static_cast<float>(births[i]),
505 static_cast<float>(births[i]),
506 0,
507 };
508 std::array<float, 3> pt1{
509 static_cast<float>(births[i]),
510 static_cast<float>(births[i] + perss[i]),
511 0,
512 };
513 points->SetPoint(2 * i + 0, pt0.data());
514 points->SetPoint(2 * i + 1, pt1.data());
515 }
516
517 TTK_FORCE_USE(nThreads);
518}
519
520int ProjectDiagramIn2D(vtkUnstructuredGrid *const inputDiagram,
521 vtkUnstructuredGrid *const outputDiagram,
522 const ttk::Debug &dbg) {
523
524 ttk::Timer tm{};
525
526 outputDiagram->ShallowCopy(inputDiagram);
527
528 auto pointData = outputDiagram->GetPointData();
529
530 auto birth = inputDiagram->GetCellData()->GetArray(ttk::PersistenceBirthName);
531 auto pers = inputDiagram->GetCellData()->GetArray(ttk::PersistenceName);
532
533 if(birth == nullptr || pers == nullptr) {
534 dbg.printWrn("Missing Birth or Persistence arrays");
535 return 1;
536 }
537
538 // generate a new `Coordinates` pointData array
539 vtkNew<vtkFloatArray> coords{};
540 coords->DeepCopy(inputDiagram->GetPoints()->GetData());
541 coords->SetName(ttk::PersistenceCoordinatesName);
542 pointData->AddArray(coords);
543
544 vtkNew<vtkPoints> points{};
545 const auto nPoints = inputDiagram->GetNumberOfPoints();
546 points->SetNumberOfPoints(nPoints);
547
548 if(birth->GetNumberOfTuples() != nPoints / 2
549 || pers->GetNumberOfTuples() != nPoints / 2) {
550 dbg.printWrn("Wrong number of tuples for Birth or Persistence arrays");
551 return 2;
552 }
553
554 switch(birth->GetDataType()) {
555 vtkTemplateMacro(getCoords(points, ttkUtils::GetPointer<VTK_TT>(birth),
556 ttkUtils::GetPointer<VTK_TT>(pers), nPoints,
557 dbg.getThreadNumber()));
558 }
559
560 outputDiagram->SetPoints(points);
561
562 // add diagonal(first point -> last birth/penultimate point)
563 std::array<vtkIdType, 2> diag{0, 2 * (outputDiagram->GetNumberOfCells() - 1)};
564 outputDiagram->InsertNextCell(VTK_LINE, 2, diag.data());
565
566 // add diagonal data
567 auto cellData = outputDiagram->GetCellData();
568 auto pairIdentifierScalars = vtkIntArray::SafeDownCast(
569 cellData->GetArray(ttk::PersistencePairIdentifierName));
570 auto extremumIndexScalars = vtkIntArray::SafeDownCast(
571 cellData->GetArray(ttk::PersistencePairTypeName));
572 auto persistenceScalars = cellData->GetArray(ttk::PersistenceName);
573 auto birthScalars = cellData->GetArray(ttk::PersistenceBirthName);
574 auto isFinite = cellData->GetArray(ttk::PersistenceIsFinite);
575
576 pairIdentifierScalars->InsertNextTuple1(-1);
577 extremumIndexScalars->InsertNextTuple1(-1);
578 isFinite->InsertNextTuple1(0);
579 // 2 * persistence of min-max pair
580 persistenceScalars->InsertNextTuple1(2 * persistenceScalars->GetTuple1(0));
581 // birth == death == 0
582 birthScalars->InsertNextTuple1(0);
583
584 // don't forget to forward the Field Data
585 outputDiagram->GetFieldData()->ShallowCopy(inputDiagram->GetFieldData());
586
587 dbg.printMsg("Projected Persistence Diagram back to 2D", 1.0,
588 tm.getElapsedTime(), dbg.getThreadNumber());
589
590 return 0;
591}
592
593int TranslateDiagram(vtkUnstructuredGrid *const diagram,
594 const std::array<double, 3> &trans) {
595
596 vtkNew<vtkUnstructuredGrid> tmp{};
597 tmp->ShallowCopy(diagram);
598
599 vtkNew<vtkTransform> tr{};
600 tr->Translate(trans.data());
601
602 vtkNew<vtkTransformFilter> trf{};
603 trf->SetTransform(tr);
604 trf->SetInputData(tmp);
605 trf->Update();
606
607 diagram->ShallowCopy(trf->GetOutputDataObject(0));
608
609 return 0;
610}
611
612int ResetDiagramPosition(vtkUnstructuredGrid *const diagram,
613 const ttk::Debug &dbg) {
614
615 const bool embedded
616 = diagram->GetPointData()->GetArray(ttk::PersistenceCoordinatesName)
617 == nullptr;
618
619 if(embedded) {
620 dbg.printWrn("Cannot reset embedded diagram position");
621 return 1;
622 }
623
624 // position of first point in diagram
625 std::array<double, 3> pos{};
626 diagram->GetPoint(diagram->GetCell(0)->GetPointId(0), pos.data());
627
628 // birth value of the first cell in the diagram
629 const auto firstBirth{
630 diagram->GetCellData()->GetArray(ttk::PersistenceBirthName)->GetTuple1(0)};
631
632 if((pos[0] != pos[1] && pos[0] != firstBirth) || pos[2] != 0) {
633 vtkNew<vtkUnstructuredGrid> tmp{};
634 tmp->ShallowCopy(diagram);
635
636 vtkNew<vtkTransform> tr{};
637 tr->Translate(firstBirth - pos[0], firstBirth - pos[1], -pos[2]);
638
639 vtkNew<vtkTransformFilter> trf{};
640 trf->SetTransform(tr);
641 trf->SetInputData(tmp);
642 trf->Update();
643
644 diagram->ShallowCopy(trf->GetOutputDataObject(0));
645
646 dbg.printMsg("Diagram reset to initial position");
647 }
648
649 return 0;
650}
TTK_FORCE_USE
#define TTK_FORCE_USE(x)
Force the compiler to use the function/method parameter.
Definition BaseClass.h:57
ttkUtils::GetPointer
static DT * GetPointer(vtkDataArray *array, vtkIdType start=0)
Definition ttkUtils.h:59
ttk::BaseClass::getThreadNumber
int getThreadNumber() const
Definition BaseClass.h:76
ttk::Debug
Minimalist debugging class.
Definition Debug.h:88
ttk::Debug::printWrn
int printWrn(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition Debug.h:159
ttk::Debug::printMsg
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
Definition Debug.h:118
ttk::SimplexId
int SimplexId
Identifier type for simplices of any dimension.
Definition DataTypes.h:22
ttk::PersistenceName
const char PersistenceName[]
Definition DataTypes.h:80
ttk::CriticalType
CriticalType
default value for critical index
Definition DataTypes.h:88
ttk::PersistenceCoordinatesName
const char PersistenceCoordinatesName[]
Definition DataTypes.h:78
ttk::PersistencePairTypeName
const char PersistencePairTypeName[]
Definition DataTypes.h:81
ttk::PersistenceCriticalTypeName
const char PersistenceCriticalTypeName[]
Definition DataTypes.h:75
ttk::PersistenceIsFinite
const char PersistenceIsFinite[]
Definition DataTypes.h:82
ttk::VertexScalarFieldName
const char VertexScalarFieldName[]
default name for vertex scalar field
Definition DataTypes.h:35
ttk::PersistencePairIdentifierName
const char PersistencePairIdentifierName[]
Definition DataTypes.h:79
ttk::PersistenceBirthName
const char PersistenceBirthName[]
Definition DataTypes.h:76
ttk::DiagramType
std::vector< PersistencePair > DiagramType
Persistence Diagram type as a vector of Persistence pairs.
Definition PersistenceDiagramUtils.h:60
ProjectDiagramInsideDomain
int ProjectDiagramInsideDomain(vtkUnstructuredGrid *const inputDiagram, vtkUnstructuredGrid *const outputDiagram, const ttk::Debug &dbg)
Generate the spatial embedding of a given Persistence Diagram.
Definition ttkPersistenceDiagramUtils.cpp:447
DiagramToVTU
int DiagramToVTU(vtkUnstructuredGrid *vtu, const ttk::DiagramType &diagram, vtkDataArray *const inputScalars, const ttk::Debug &dbg, const int dim, const bool embedInDomain)
Converts a Persistence Diagram in the ttk::DiagramType format to the VTK Unstructured Grid format (as...
Definition ttkPersistenceDiagramUtils.cpp:158
TranslateDiagram
int TranslateDiagram(vtkUnstructuredGrid *const diagram, const std::array< double, 3 > &trans)
Translate a diagram to a new position.
Definition ttkPersistenceDiagramUtils.cpp:593
VTUToDiagram
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...
Definition ttkPersistenceDiagramUtils.cpp:17
ProjectDiagramIn2D
int ProjectDiagramIn2D(vtkUnstructuredGrid *const inputDiagram, vtkUnstructuredGrid *const outputDiagram, const ttk::Debug &dbg)
Generate the 2D embedding of a given Persistence Diagram.
Definition ttkPersistenceDiagramUtils.cpp:520
getCoords
void getCoords(vtkPoints *points, const dataType *const births, const dataType *const perss, vtkIdType nPoints, const int nThreads)
Definition ttkPersistenceDiagramUtils.cpp:493
ResetDiagramPosition
int ResetDiagramPosition(vtkUnstructuredGrid *const diagram, const ttk::Debug &dbg)
Translate back a canonical diagram into its original position.
Definition ttkPersistenceDiagramUtils.cpp:612
Generated on Fri Apr 10 2026 10:02:09 for TTK by doxygen 1.13.2