TTK
Loading...
Searching...
No Matches
ttkTrackingFromOverlap.cpp
Go to the documentation of this file.
2
3#include <vtkPointSet.h>
4#include <vtkStreamingDemandDrivenPipeline.h>
5#include <vtkUnstructuredGrid.h>
6
7#include <vtkCellData.h>
8#include <vtkCharArray.h>
9#include <vtkDoubleArray.h>
10#include <vtkFloatArray.h>
11#include <vtkIdTypeArray.h>
12#include <vtkInformation.h>
13#include <vtkInformationVector.h>
14#include <vtkLongLongArray.h>
15#include <vtkPointData.h>
16
17#include <ttkMacros.h>
18#include <ttkUtils.h>
19
20using namespace std;
21using namespace ttk;
22
24
25// Function to fetch data of a specificd block
26static void getData(vtkMultiBlockDataSet *mb,
27 size_t time,
28 size_t level,
29 const string &labelFieldName,
30 vtkPointSet *&pointSet,
31 vtkDataArray *&labels) {
32 auto timesteps = vtkMultiBlockDataSet::SafeDownCast(mb->GetBlock(level));
33 pointSet = vtkPointSet::SafeDownCast(timesteps->GetBlock(time));
34 if(pointSet == nullptr) {
35 return;
36 }
37 const auto pd = pointSet->GetPointData();
38 labels = pd->GetArray(labelFieldName.data());
39};
40
41// Function to compute number of levels and timesteps contained in a
42// vtkMultiBlockDataSet
43static void getNumberOfLevelsAndTimesteps(vtkMultiBlockDataSet *mb,
44 size_t &nL,
45 size_t &nT) {
46 nL = mb->GetNumberOfBlocks();
47 auto timesteps = vtkMultiBlockDataSet::SafeDownCast(mb->GetBlock(0));
48 nT = timesteps->GetNumberOfBlocks();
49};
50
51// =============================================================================
52// Finalize
53// =============================================================================
54template <typename labelType>
55int finalize(vector<vector<TrackingFromOverlap::Nodes>> &levelTimeNodesMap,
56 vector<vector<TrackingFromOverlap::Edges>> &levelTimeEdgesTMap,
57 vector<vector<TrackingFromOverlap::Edges>> &timeLevelEdgesNMap,
58 int labelTypeId,
59 const string &labelFieldName,
60
61 vtkDataObject *trackingGraphObject) {
62 auto trackingGraph = vtkUnstructuredGrid::SafeDownCast(trackingGraphObject);
63
64 size_t nL = levelTimeNodesMap.size();
65 size_t nT = levelTimeNodesMap[0].size();
66
67 auto prepArray = [](vtkAbstractArray *array, const string &name,
68 size_t nComponents, size_t nValues) {
69 array->SetName(name.data());
70 array->SetNumberOfComponents(nComponents);
71 array->SetNumberOfTuples(nValues);
72 };
73
74 // Add Points
75 {
76 size_t nNodes = 0;
77 for(size_t t = 0; t < nT; t++)
78 for(size_t l = 0; l < nL; l++)
79 nNodes += levelTimeNodesMap[l][t].size();
80
81 auto points = vtkSmartPointer<vtkPoints>::New();
82 points->SetNumberOfPoints(nNodes);
83 auto pointCoords = (float *)ttkUtils::GetVoidPointer(points);
84
86 prepArray(sequence, "SequenceIndex", 1, nNodes);
87 auto sequenceData = (long long *)ttkUtils::GetVoidPointer(sequence);
88
90 prepArray(level, "LevelIndex", 1, nNodes);
91 auto levelData = (long long *)ttkUtils::GetVoidPointer(level);
92
94 prepArray(size, "Size", 1, nNodes);
95 auto sizeData = (float *)ttkUtils::GetVoidPointer(size);
96
98 prepArray(branch, "BranchId", 1, nNodes);
99 auto branchData = (long long *)ttkUtils::GetVoidPointer(branch);
100
102 vtkDataArray::CreateDataArray(labelTypeId));
103 prepArray(label, labelFieldName, 1, nNodes);
104 auto labelData = (labelType *)ttkUtils::GetVoidPointer(label);
105
106 size_t q1 = 0, q2 = 0;
107 for(size_t t = 0; t < nT; t++) {
108 for(size_t l = 0; l < nL; l++) {
109 for(auto &node : levelTimeNodesMap[l][t]) {
110 pointCoords[q1++] = node.x;
111 pointCoords[q1++] = node.y;
112 pointCoords[q1++] = node.z;
113
114 sequenceData[q2] = t;
115 levelData[q2] = l;
116 sizeData[q2] = node.size;
117 branchData[q2] = node.branchID;
118 labelData[q2] = boost::get<labelType>(node.label);
119
120 q2++;
121 }
122 }
123 }
124
125 trackingGraph->SetPoints(points);
126
127 auto pointData = trackingGraph->GetPointData();
128 pointData->AddArray(sequence);
129 pointData->AddArray(level);
130 pointData->AddArray(size);
131 pointData->AddArray(label);
132 pointData->AddArray(branch);
133 }
134
135 // Add Cells
136 {
137 if(nT * nL + 1 <= 0) {
138 return 0;
139 }
140
141 // Build node index offset vector
142 vector<size_t> timeLevelOffsetMap(nT * nL + 1);
143 {
144 timeLevelOffsetMap[0] = 0;
145 size_t q = 1;
146 for(size_t t = 0; t < nT; t++)
147 for(size_t l = 0; l < nL; l++) {
148 timeLevelOffsetMap[q]
149 = timeLevelOffsetMap[q - 1] + levelTimeNodesMap[l][t].size();
150 q++;
151 }
152 }
153
154 size_t nEdgesT = 0;
155 if(nT > 1)
156 for(size_t t = 0; t < nT - 1; t++)
157 for(size_t l = 0; l < nL; l++)
158 nEdgesT += levelTimeEdgesTMap[l][t].size() / 4;
159
160 size_t nEdgesN = 0;
161 if(nL > 1)
162 for(size_t l = 0; l < nL - 1; l++)
163 for(size_t t = 0; t < nT; t++)
164 nEdgesN += timeLevelEdgesNMap[t][l].size() / 4;
165
167 cells->SetNumberOfValues(3 * nEdgesT + 3 * nEdgesN);
168 auto cellIds = (vtkIdType *)ttkUtils::GetVoidPointer(cells);
169
171 prepArray(overlap, "Overlap", 1, nEdgesT + nEdgesN);
172 auto overlapData = (float *)ttkUtils::GetVoidPointer(overlap);
173
175 prepArray(branch, "BranchId", 1, nEdgesT + nEdgesN);
176 auto branchData = (long long *)ttkUtils::GetVoidPointer(branch);
177
179 prepArray(type, "Type", 1, nEdgesT + nEdgesN);
180 auto typeData = (char *)ttkUtils::GetVoidPointer(type);
181
182 size_t q0 = 0, q1 = 0;
183
184 // Tracking graphs
185 if(nT > 1)
186 for(size_t t = 1; t < nT; t++) {
187 for(size_t l = 0; l < nL; l++) {
188 auto &edges = levelTimeEdgesTMap[l][t - 1];
189 for(size_t i = 0, j = edges.size(); i < j;) {
190 cellIds[q0++] = 2;
191 cellIds[q0++]
192 = (vtkIdType)(timeLevelOffsetMap[(t - 1) * nL + l] + edges[i++]);
193 cellIds[q0++]
194 = (vtkIdType)(timeLevelOffsetMap[(t)*nL + l] + edges[i++]);
195 typeData[q1] = 0;
196 overlapData[q1] = edges[i++];
197 branchData[q1] = edges[i++];
198 q1++;
199 }
200 }
201 }
202
203 // Nesting trees
204 if(nL > 1)
205 for(size_t l = 1; l < nL; l++) {
206 for(size_t t = 0; t < nT; t++) {
207 auto &edges = timeLevelEdgesNMap[t][l - 1];
208 size_t temp = t * nL;
209 for(size_t i = 0, j = edges.size(); i < j;) {
210 cellIds[q0++] = 2;
211 cellIds[q0++]
212 = (vtkIdType)(timeLevelOffsetMap[temp + (l - 1)] + edges[i++]);
213 cellIds[q0++]
214 = (vtkIdType)(timeLevelOffsetMap[temp + (l)] + edges[i++]);
215 typeData[q1] = 1;
216 overlapData[q1] = edges[i++];
217 branchData[q1] = edges[i++];
218 q1++;
219 }
220 }
221 }
222
223 auto cellArray = vtkSmartPointer<vtkCellArray>::New();
224 cellArray->SetCells(nEdgesT + nEdgesN, cells);
225 trackingGraph->SetCells(VTK_LINE, cellArray);
226
227 auto cellData = trackingGraph->GetCellData();
228 cellData->AddArray(type);
229 cellData->AddArray(overlap);
230 cellData->AddArray(branch);
231 }
232
233 return 1;
234}
235
236// =============================================================================
237// Mesh Nested Tracking Graph
238// =============================================================================
240 vtkDataObject *trackingGraph) {
241 Timer t;
242
243 printMsg("=======================================================",
244 debug::Priority::INFO);
245 printMsg("Meshing nested tracking graph", debug::Priority::INFO);
246
247 switch(this->LabelDataType) {
248 vtkTemplateMacro(
249 finalize<VTK_TT>(this->levelTimeNodesMap, this->levelTimeEdgesTMap,
250 this->timeLevelEdgesNMap, this->LabelDataType,
251 this->GetLabelFieldName(), trackingGraph));
252 }
253
254 printMsg("-------------------------------------------------------");
255 stringstream msg;
256 msg << "Nested tracking graph meshed in " << t.getElapsedTime() << " s. ("
257 << threadNumber_ << " thread(s)).";
258 printMsg(msg.str(), debug::Priority::PERFORMANCE);
259
260 return 1;
261}
262
263// =============================================================================
264// Reset
265// =============================================================================
267 this->LabelDataType = -1;
268
269 this->levelTimeNodesMap.clear();
270 this->levelTimeEdgesTMap.clear();
271 this->timeLevelEdgesNMap.clear();
272
273 this->previousIterationData = nullptr;
274
275 return 1;
276}
277
278// =============================================================================
279// Pack Data
280// =============================================================================
282 vtkDataObject *inputDataObject, vtkMultiBlockDataSet *packedInput) const {
283 /* Enforce following vtkMultiBlockDataSet structure:
284 {
285 level_0: {
286 time_0: vtkPointSet,
287 ...
288 time_nT: vtkPointSet
289 },
290 ...
291 level_nL: {
292 time_0: vtkPointSet,
293 ...
294 time_nT: vtkPointSet
295 }
296 }
297 */
298
299 // Check inputDataObject depth: 2->(level->time), 1->(-,time), 0->(-,-)
300 auto inputAsMB = vtkMultiBlockDataSet::SafeDownCast(inputDataObject);
301 auto inputAsPS = vtkPointSet::SafeDownCast(inputDataObject);
302 bool error = false;
303 if(inputAsMB) {
304 size_t n = inputAsMB->GetNumberOfBlocks();
305
306 // Check if blocks are vtkPointSets or vtkMultiBlockDataSets ...
307 size_t psCounter = 0;
308 size_t mbCounter = 0;
309 for(size_t i = 0; i < n; i++) {
310 auto block = inputAsMB->GetBlock(i);
311 auto blockAsMB = vtkMultiBlockDataSet::SafeDownCast(block);
312 auto blockAsPS = vtkPointSet::SafeDownCast(block);
313 if(blockAsMB)
314 mbCounter++;
315 if(blockAsPS)
316 psCounter++;
317 }
318
319 if(mbCounter
320 == n) { // input already contains timesteps per level -> nothing to do
321 packedInput->ShallowCopy(inputAsMB);
322 } else if(psCounter
323 == n) { // if input is a single list of vtkPointSets over time
325 for(size_t i = 0; i < n; i++) {
326 level->SetBlock(i, vtkPointSet::SafeDownCast(inputAsMB->GetBlock(i)));
327 }
328 packedInput->SetBlock(0, level);
329 } else { // Unexpected input structure
330 error = true;
331 }
332 } else if(inputAsPS) {
334 level->SetBlock(0, inputAsPS);
335 packedInput->SetBlock(0, level);
336 } else { // Unexpected input structure
337 error = true;
338 }
339
340 if(error) {
341 printErr("Unable to convert input into "
342 "'vtkPointSet' collection.");
343 return 0;
344 }
345
346 return 1;
347}
348
349// =============================================================================
350// Check Data
351// =============================================================================
352int ttkTrackingFromOverlap::checkData(vtkMultiBlockDataSet *data) {
353 size_t nL = data->GetNumberOfBlocks();
354 size_t nT = 0;
355
356 if(nL < 1) {
357 printErr("Input must have at least one "
358 "vtkPointSet.");
359 return 0;
360 }
361
362 for(size_t l = 0; l < nL; l++) {
363 auto timesteps = vtkMultiBlockDataSet::SafeDownCast(data->GetBlock(l));
364 size_t n = timesteps->GetNumberOfBlocks();
365 if(n < 1) {
366 printErr("Input must have at least one "
367 "vtkPointSet.");
368 return 0;
369 }
370 if(nT == 0)
371 nT = n;
372 if(nT != n) {
373 printErr("Timeseries have unequal length.");
374 return 0;
375 }
376
377 for(size_t t = 0; t < nT; t++) {
378 auto pointSet = vtkPointSet::SafeDownCast(timesteps->GetBlock(t));
379 if(pointSet == nullptr) {
380 return 0;
381 }
382 size_t nPoints = pointSet->GetNumberOfPoints();
383 auto labels = pointSet->GetPointData()->GetAbstractArray(
384 this->GetLabelFieldName().data());
385
386 if(nPoints > 0 && labels == nullptr) {
387 printErr("Point labels '" + this->GetLabelFieldName() + "' not found.");
388 return 0;
389 }
390 if(labels == nullptr)
391 continue;
392
393 int labelDataType = labels->GetDataType();
394 if(this->LabelDataType < 0)
395 this->LabelDataType = labelDataType;
396 if(this->LabelDataType != labelDataType) {
397 printErr("Point labels do not have same "
398 "type across point sets.");
399 return 0;
400 }
401 }
402 }
403
404 return 1;
405}
406
407// =============================================================================
408// Pack Streamed Data
409// =============================================================================
410int ttkTrackingFromOverlap::packStreamedData(vtkMultiBlockDataSet *streamedData,
411 vtkMultiBlockDataSet *data) const {
412 size_t nL_PI = this->previousIterationData->GetNumberOfBlocks();
413 size_t nL_CI = streamedData->GetNumberOfBlocks();
414 if(nL_PI != nL_CI) {
415 printErr("Number of levels differ over time.");
416 return 0;
417 }
418 for(size_t l = 0; l < nL_PI; l++) {
419 auto timestepsPI = vtkMultiBlockDataSet::SafeDownCast(
420 this->previousIterationData->GetBlock(l));
421 auto timestepsCI
422 = vtkMultiBlockDataSet::SafeDownCast(streamedData->GetBlock(l));
423 size_t nT_CI = timestepsCI->GetNumberOfBlocks();
424
426
427 // First timestep is previous iteration
428 timesteps->SetBlock(0, timestepsPI->GetBlock(0));
429
430 // Remaining timesteps are from current iteration
431 for(size_t t = 0; t < nT_CI; t++)
432 timesteps->SetBlock(t + 1, timestepsCI->GetBlock(t));
433
434 data->SetBlock(l, timesteps);
435 }
436
437 return 1;
438}
439
440// =============================================================================
441// Store Streamed Data
442// =============================================================================
444 vtkMultiBlockDataSet *streamedData) {
445
447 size_t nL = streamedData->GetNumberOfBlocks();
448 for(size_t l = 0; l < nL; l++) {
449 auto timestepsSD
450 = vtkMultiBlockDataSet::SafeDownCast(streamedData->GetBlock(l));
451 size_t nT = timestepsSD->GetNumberOfBlocks();
452
455 lastTimestep->SetBlock(0, timestepsSD->GetBlock(nT - 1));
456
457 temp->SetBlock(l, lastTimestep);
458 }
459
460 this->previousIterationData = vtkSmartPointer<vtkMultiBlockDataSet>::New();
461 this->previousIterationData->DeepCopy(temp);
462
463 return 1;
464}
465
466// =============================================================================
467// Compute Nodes
468// =============================================================================
469int ttkTrackingFromOverlap::computeNodes(vtkMultiBlockDataSet *data) {
470
471 Timer timer;
472
473 size_t nL, nT;
474 getNumberOfLevelsAndTimesteps(data, nL, nT);
475
476 // Reusable variables
477 vtkPointSet *pointSet = nullptr;
478 vtkDataArray *labels = nullptr;
479
480 // Compute Nodes
481 {
482 printMsg("=======================================================",
483 debug::Priority::INFO);
484 printMsg("Computing nodes", debug::Priority::INFO);
485
486 if(this->levelTimeNodesMap.size() != nL)
487 this->levelTimeNodesMap.resize(nL);
488
489 for(size_t l = 0; l < nL; l++) {
490 {
491 printMsg("-------------------------------------------------------");
492 stringstream msg;
493 msg << "Level Index: " << l;
494 printMsg(msg.str(), debug::Priority::INFO);
495 }
496
497 vector<Nodes> &timeNodesMap = this->levelTimeNodesMap[l];
498 size_t timeOffset = timeNodesMap.size();
499 timeNodesMap.resize(timeOffset + nT);
500
501 for(size_t t = 0; t < nT; t++) {
502 getData(data, t, l, this->GetLabelFieldName(), pointSet, labels);
503
504 size_t nPoints = pointSet->GetNumberOfPoints();
505 if(nPoints < 1)
506 continue;
507
508 switch(this->LabelDataType) {
509 vtkTemplateMacro(this->ttk::TrackingFromOverlap::computeNodes<VTK_TT>(
510 (float *)ttkUtils::GetVoidPointer(pointSet->GetPoints()),
511 (VTK_TT *)ttkUtils::GetVoidPointer(labels), nPoints,
512 timeNodesMap[timeOffset + t]));
513 }
514 }
515 }
516
517 {
518 printMsg("-------------------------------------------------------");
519 stringstream msg;
520 msg << "Nodes computed in " << timer.getElapsedTime() << " s. ("
521 << threadNumber_ << " thread(s)).";
522 printMsg(msg.str(), debug::Priority::PERFORMANCE);
523 }
524 }
525
526 return 1;
527}
528
529// =============================================================================
530// Compute Tracking Graphs
531// =============================================================================
532int ttkTrackingFromOverlap::computeTrackingGraphs(vtkMultiBlockDataSet *data) {
533
534 Timer timer;
535
536 size_t nL, nT;
537 getNumberOfLevelsAndTimesteps(data, nL, nT);
538
539 if(nT < 2)
540 return 1;
541
542 // Reusable variables
543 vtkPointSet *pointSet0 = nullptr;
544 vtkPointSet *pointSet1 = nullptr;
545 vtkDataArray *labels0 = nullptr;
546 vtkDataArray *labels1 = nullptr;
547
548 printMsg("=======================================================",
549 debug::Priority::INFO);
550 printMsg("Computing tracking graphs", debug::Priority::INFO);
551
552 if(this->levelTimeEdgesTMap.size() != nL)
553 this->levelTimeEdgesTMap.resize(nL);
554
555 for(size_t l = 0; l < nL; l++) {
556 {
557 printMsg("-------------------------------------------------------");
558 stringstream msg;
559 msg << "Level Index: " << l;
560 printMsg(msg.str(), debug::Priority::INFO);
561 }
562
563 vector<Edges> &timeEdgesTMap = this->levelTimeEdgesTMap[l];
564 size_t timeOffset = timeEdgesTMap.size();
565 timeEdgesTMap.resize(timeOffset + nT - 1);
566
567 for(size_t t = 1; t < nT; t++) {
568 getData(data, t - 1, l, this->GetLabelFieldName(), pointSet0, labels0);
569 getData(data, t, l, this->GetLabelFieldName(), pointSet1, labels1);
570
571 size_t nPoints0 = pointSet0->GetNumberOfPoints();
572 size_t nPoints1 = pointSet1->GetNumberOfPoints();
573 if(nPoints0 < 1 || nPoints1 < 1)
574 continue;
575
576 switch(this->LabelDataType) {
577 vtkTemplateMacro(this->computeOverlap<VTK_TT>(
578 (float *)ttkUtils::GetVoidPointer(pointSet0->GetPoints()),
579 (float *)ttkUtils::GetVoidPointer(pointSet1->GetPoints()),
580 (VTK_TT *)ttkUtils::GetVoidPointer(labels0),
581 (VTK_TT *)ttkUtils::GetVoidPointer(labels1), nPoints0, nPoints1,
582 timeEdgesTMap[timeOffset + t - 1]));
583 }
584 }
585 }
586
587 {
588 printMsg("-------------------------------------------------------");
589 stringstream msg;
590 msg << "Tracking graphs computed in " << timer.getElapsedTime() << " s. ("
591 << threadNumber_ << " thread(s)).";
592 printMsg(msg.str(), debug::Priority::PERFORMANCE);
593 }
594 return 1;
595}
596
597// =============================================================================
598// Compute Nesting Trees
599// =============================================================================
600int ttkTrackingFromOverlap::computeNestingTrees(vtkMultiBlockDataSet *data) {
601
602 Timer timer;
603
604 size_t nL, nT;
605 getNumberOfLevelsAndTimesteps(data, nL, nT);
606
607 if(nL < 2)
608 return 1;
609
610 // Reusable variables
611 vtkPointSet *pointSet0 = nullptr;
612 vtkPointSet *pointSet1 = nullptr;
613 vtkDataArray *labels0 = nullptr;
614 vtkDataArray *labels1 = nullptr;
615
616 printMsg("=======================================================",
617 debug::Priority::INFO);
618 printMsg("Computing nesting trees", debug::Priority::INFO);
619
620 size_t timeOffset = this->timeLevelEdgesNMap.size();
621 this->timeLevelEdgesNMap.resize(timeOffset + nT);
622
623 for(size_t t = 0; t < nT; t++) {
624 {
625 printMsg("-------------------------------------------------------");
626 stringstream msg;
627 msg << "Time Index: " << t;
628 printMsg(msg.str(), debug::Priority::INFO);
629 }
630
631 vector<Edges> &levelEdgesNMap = this->timeLevelEdgesNMap[timeOffset + t];
632 levelEdgesNMap.resize(nL - 1);
633
634 for(size_t l = 1; l < nL; l++) {
635 getData(data, t, l - 1, this->GetLabelFieldName(), pointSet0, labels0);
636 getData(data, t, l, this->GetLabelFieldName(), pointSet1, labels1);
637
638 size_t nPoints0 = pointSet0->GetNumberOfPoints();
639 size_t nPoints1 = pointSet1->GetNumberOfPoints();
640 if(nPoints0 < 1 || nPoints1 < 1)
641 continue;
642
643 switch(this->LabelDataType) {
644 vtkTemplateMacro(this->computeOverlap<VTK_TT>(
645 (float *)ttkUtils::GetVoidPointer(pointSet0->GetPoints()),
646 (float *)ttkUtils::GetVoidPointer(pointSet1->GetPoints()),
647 (VTK_TT *)ttkUtils::GetVoidPointer(labels0),
648 (VTK_TT *)ttkUtils::GetVoidPointer(labels1), nPoints0, nPoints1,
649 levelEdgesNMap[l - 1]));
650 }
651 }
652 }
653
654 {
655 printMsg("-------------------------------------------------------");
656 stringstream msg;
657 msg << "Nesting trees computed in " << timer.getElapsedTime() << " s. ("
658 << threadNumber_ << " thread(s)).";
659 printMsg(msg.str(), debug::Priority::PERFORMANCE);
660 }
661
662 return 1;
663}
664
665// =============================================================================
666// Compute Branches
667// =============================================================================
669
670 size_t nL = this->levelTimeEdgesTMap.size();
671
672 for(size_t l = 0; l < nL; l++)
674 this->levelTimeEdgesTMap[l], this->levelTimeNodesMap[l]);
675
676 return 1;
677}
678
679// =============================================================================
680// Request Data
681// =============================================================================
683 vtkInformationVector **inputVector,
684 vtkInformationVector *outputVector) {
685 Timer timer;
686
687 // Print Status
688 {
689 printMsg(
690 "===================================================================");
691 printMsg("RequestData", debug::Priority::INFO);
692 }
693
694 // Get Input Object
695 vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
696 auto inputObject = inInfo->Get(vtkDataObject::DATA_OBJECT());
697
698 // Get iteration information
699 auto iterationInformation = vtkDoubleArray::SafeDownCast(
700 inputObject->GetFieldData()->GetAbstractArray("_ttk_IterationInfo"));
701
702 bool useStreamingOverTime = iterationInformation != nullptr;
703
704 double iteration = 0;
705 double nIterations = 0;
706 if(useStreamingOverTime) {
707 iteration = iterationInformation->GetValue(0);
708 nIterations = iterationInformation->GetValue(1);
709 }
710
711 // On first iteration reset
712 if(!useStreamingOverTime || iteration == 0)
713 this->reset();
714
715 // -------------------------------------------------------------------------
716 // Prepare Input
717 // -------------------------------------------------------------------------
719
720 if(!this->packInputData(inputObject, packedInput))
721 return 0;
722
723 // Check input integrity
724 if(!this->checkData(packedInput))
725 return 0;
726
727 // During streaming append data to previous timestep; Otherwise pass through
729 if(useStreamingOverTime && this->previousIterationData != nullptr) {
730 if(!this->packStreamedData(packedInput, data))
731 return 0;
732 } else
733 data->ShallowCopy(packedInput);
734
735 // -------------------------------------------------------------------------
736 // Process Input
737 // -------------------------------------------------------------------------
738 // Compute nodes only for current input (previous iterations have already been
739 // processed)
740 if(!this->computeNodes(packedInput))
741 return 0;
742
743 // Compute tracking graphs
744 if(!this->computeTrackingGraphs(data))
745 return 0;
746
747 // Compute nesting trees
748 if(!this->computeNestingTrees(packedInput))
749 return 0;
750
751 // Store last timestep for next iteration
752 if(useStreamingOverTime && !this->storeStreamedData(packedInput))
753 return 0;
754
755 // -------------------------------------------------------------------------
756 // Generate Output
757 // -------------------------------------------------------------------------
758 if(!useStreamingOverTime || iteration == nIterations - 1) {
759 // Get Output
760 vtkInformation *outInfo = outputVector->GetInformationObject(0);
761 auto trackingGraph = outInfo->Get(vtkDataObject::DATA_OBJECT());
762
763 // Compute Branches
764 if(!this->computeBranches())
765 return 0;
766
767 // Mesh Graph
768 if(!this->meshNestedTrackingGraph(trackingGraph))
769 return 0;
770 }
771
772 // -------------------------------------------------------------------------
773 // Print total performance
774 // -------------------------------------------------------------------------
775 if(!useStreamingOverTime) {
776 printMsg("=======================================================");
777 stringstream msg;
778 msg << "Nested tracking graph generated in " << timer.getElapsedTime()
779 << " s. (" << threadNumber_ << " thread(s)).";
780 printMsg(msg.str(), debug::Priority::PERFORMANCE);
781 }
782
783 return 1;
784}
#define ttkNotUsed(x)
Mark function/method parameters that are not used in the function body at all.
Definition: BaseClass.h:47
TTK VTK-filter that computes the overlap between labeled vtkPointSets.
int checkData(vtkMultiBlockDataSet *data)
int computeTrackingGraphs(vtkMultiBlockDataSet *data)
int storeStreamedData(vtkMultiBlockDataSet *data)
int computeNestingTrees(vtkMultiBlockDataSet *data)
int meshNestedTrackingGraph(vtkDataObject *trackingGraph)
int computeNodes(vtkMultiBlockDataSet *data)
int packInputData(vtkDataObject *inputDataObject, vtkMultiBlockDataSet *packedData) const
int packStreamedData(vtkMultiBlockDataSet *streamedData, vtkMultiBlockDataSet *packedData) const
virtual std::string GetLabelFieldName()
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
static void * GetVoidPointer(vtkDataArray *array, vtkIdType start=0)
Definition: ttkUtils.cpp:225
int threadNumber_
Definition: BaseClass.h:95
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
int printErr(const std::string &msg, const debug::LineMode &lineMode=debug::LineMode::NEW, std::ostream &stream=std::cerr) const
Definition: Debug.h:149
double getElapsedTime()
Definition: Timer.h:15
int computeBranches(std::vector< Edges > &timeEdgesMap, std::vector< Nodes > &timeNodesMap) const
The Topology ToolKit.
int finalize(vector< vector< TrackingFromOverlap::Nodes > > &levelTimeNodesMap, vector< vector< TrackingFromOverlap::Edges > > &levelTimeEdgesTMap, vector< vector< TrackingFromOverlap::Edges > > &timeLevelEdgesNMap, int labelTypeId, const string &labelFieldName, vtkDataObject *trackingGraphObject)
vtkStandardNewMacro(ttkTrackingFromOverlap)