30#include <boost/variant.hpp>
32#include <unordered_map>
74 using Edges = std::vector<idType>;
75 using Nodes = std::vector<Node>;
83 inline bool operator()(
const size_t &i,
const size_t &j) {
97 std::vector<size_t> &sortedIndices)
const {
101 sortedIndices.resize(nPoints);
102 for(
size_t i = 0; i < nPoints; i++)
103 sortedIndices[i] = i;
105 sort(sortedIndices.begin(), sortedIndices.end(), c);
107 std::stringstream msg;
115 std::vector<Nodes> &timeNodesMap)
const {
119 size_t nT = timeNodesMap.size();
122 for(
size_t t = 1; t < nT; t++) {
123 auto &nodes0 = timeNodesMap[t - 1];
124 auto &nodes1 = timeNodesMap[t];
125 auto &edges = timeEdgesMap[t - 1];
127 size_t nE = edges.size();
129 for(
size_t i = 0; i < nE; i += 4) {
130 auto n0Index = edges[i];
131 auto n1Index = edges[i + 1];
132 auto &n0 = nodes0[n0Index];
133 auto &n1 = nodes1[n1Index];
136 = n0.maxSuccID != -1 ? nodes1[n0.maxSuccID].size : 0;
138 = n1.maxPredID != -1 ? nodes0[n1.maxPredID].size : 0;
139 if(n0MaxSuccSize < n1.size)
140 n0.maxSuccID = n1Index;
141 if(n1MaxPredSize < n0.size)
142 n1.maxPredID = n0Index;
149 for(
size_t t = 0; t < nT; t++)
150 for(
auto &n : timeNodesMap[t])
151 n.branchID = n.maxPredID == -1 ? branchCounter++ : -1;
153 for(
size_t t = 1; t < nT; t++) {
154 auto &nodes0 = timeNodesMap[t - 1];
155 auto &nodes1 = timeNodesMap[t];
157 for(
size_t i = 0; i < nodes1.size(); i++) {
158 auto &n1 = nodes1[i];
159 if(n1.maxPredID != -1
160 && ((
idType)i) != nodes0[n1.maxPredID].maxSuccID)
161 n1.branchID = branchCounter++;
166 for(
size_t t = 1; t < nT; t++) {
167 auto &nodes0 = timeNodesMap[t - 1];
168 auto &nodes1 = timeNodesMap[t];
169 auto &edges = timeEdgesMap[t - 1];
171 size_t nE = edges.size();
173 for(
size_t i = 0; i < nE; i += 4) {
174 auto n0Index = edges[i];
175 auto n1Index = edges[i + 1];
176 auto &n0 = nodes0[n0Index];
177 auto &n1 = nodes1[n1Index];
179 if(n1.branchID == -1 && n0Index == n1.maxPredID)
180 n1.branchID = n0.branchID;
185 for(
size_t t = 1; t < nT; t++) {
186 auto &nodes0 = timeNodesMap[t - 1];
187 auto &nodes1 = timeNodesMap[t];
188 auto &edges = timeEdgesMap[t - 1];
190 size_t nE = edges.size();
192 for(
size_t i = 0; i < nE; i += 4) {
193 auto n0Index = edges[i];
194 auto n1Index = edges[i + 1];
195 auto &n0 = nodes0[n0Index];
196 auto &n1 = nodes1[n1Index];
198 edges[i + 3] = n0.branchID == n1.branchID ? n0.branchID
199 : n0.maxSuccID == n1Index ? n0.branchID
204 std::stringstream msg;
213 template <
typename labelType>
215 const size_t nPoints,
216 std::map<labelType, size_t> &labelIndexMap)
const;
220 template <
typename labelType>
222 const labelType *pointLabels,
223 const size_t nPoints,
227 template <
typename labelType>
229 const float *pointCoordinates1,
230 const labelType *pointLabels0,
231 const labelType *pointLabels1,
232 const size_t nPoints0,
233 const size_t nPoints1,
244template <
typename labelType>
246 const labelType *pointLabels,
247 const size_t nPoints,
248 std::map<labelType, size_t> &labelIndexMap)
const {
249 for(
size_t i = 0; i < nPoints; i++)
250 labelIndexMap[pointLabels[i]] = 0;
252 for(
auto &it : labelIndexMap)
260template <
typename labelType>
262 const labelType *pointLabels,
263 const size_t nPoints,
264 Nodes &nodes)
const {
269 std::map<labelType, size_t> labelIndexMap;
270 this->computeLabelIndexMap(pointLabels, nPoints, labelIndexMap);
272 size_t nNodes = labelIndexMap.size();
274 nodes.resize(nNodes);
275 for(
size_t i = 0, q = 0; i < nPoints; i++) {
276 labelType label = pointLabels[i];
277 Node &n = nodes[labelIndexMap[label]];
280 n.
x += pointCoordinates[q++];
281 n.
y += pointCoordinates[q++];
282 n.
z += pointCoordinates[q++];
285 for(
size_t i = 0; i < nNodes; i++) {
287 float size = (float)n.
size;
295 std::stringstream msg;
296 msg <<
"done (#" << nNodes <<
" in " << t.
getElapsedTime() <<
" s).";
306template <
typename labelType>
308 const float *pointCoordinates1,
309 const labelType *pointLabels0,
310 const labelType *pointLabels1,
311 const size_t nPoints0,
312 const size_t nPoints1,
314 Edges &edges)
const {
318 std::map<labelType, size_t> labelIndexMap0;
319 std::map<labelType, size_t> labelIndexMap1;
320 this->computeLabelIndexMap<labelType>(pointLabels0, nPoints0, labelIndexMap0);
321 this->computeLabelIndexMap<labelType>(pointLabels1, nPoints1, labelIndexMap1);
326 std::vector<size_t> sortedIndices0;
327 std::vector<size_t> sortedIndices1;
328 this->sortCoordinates(pointCoordinates0, nPoints0, sortedIndices0);
329 this->sortCoordinates(pointCoordinates1, nPoints1, sortedIndices1);
342 auto compare = [&](
size_t p0,
size_t p1) {
343 size_t p0CoordIndex = p0 * 3;
344 size_t p1CoordIndex = p1 * 3;
346 float p0_X = pointCoordinates0[p0CoordIndex++];
347 float p0_Y = pointCoordinates0[p0CoordIndex++];
348 float p0_Z = pointCoordinates0[p0CoordIndex];
350 float p1_X = pointCoordinates1[p1CoordIndex++];
351 float p1_Y = pointCoordinates1[p1CoordIndex++];
352 float p1_Z = pointCoordinates1[p1CoordIndex];
354 return p0_X == p1_X ? p0_Y == p1_Y ? p0_Z == p1_Z ? 0
367 std::unordered_map<size_t, std::unordered_map<size_t, size_t>> edgesMap;
369 while(i < nPoints0 && j < nPoints1) {
370 size_t pointIndex0 = sortedIndices0[i];
371 size_t pointIndex1 = sortedIndices1[j];
374 int c = compare(pointIndex0, pointIndex1);
377 labelType label0 = pointLabels0[pointIndex0];
378 labelType label1 = pointLabels1[pointIndex1];
380 size_t &nodeIndex0 = labelIndexMap0[label0];
381 size_t &nodeIndex1 = labelIndexMap1[label1];
384 auto edges0 = edgesMap.find(nodeIndex0);
387 if(edges0 == edgesMap.end()) {
388 edgesMap[nodeIndex0] = std::unordered_map<size_t, size_t>();
389 edges0 = edgesMap.find(nodeIndex0);
393 auto edge = edges0->second.find(nodeIndex1);
396 if(edge == edges0->second.end()) {
397 edges0->second[nodeIndex1] = 0;
398 edge = edges0->second.find(nodeIndex1);
418 edges.resize(nEdges * 4);
420 for(
auto &it0 : edgesMap) {
421 for(
auto &it1 : it0.second) {
422 edges[q++] = it0.first;
423 edges[q++] = it1.first;
424 edges[q++] = it1.second;
432 std::stringstream msg;
433 msg <<
"done (#" << nEdges <<
" in " << t.
getElapsedTime() <<
" s).";
boost::variant< double, float, long long, unsigned long long, long, unsigned long, int, unsigned int, short, unsigned short, char, signed char, unsigned char > labelTypeVariant
unsigned char topologyType
Minimalist debugging class.
void setDebugMsgPrefix(const std::string &prefix)
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
TTK trackingFromOverlap processing package that tracks labeled point sets.
std::vector< idType > Edges
int computeNodes(const float *pointCoordinates, const labelType *pointLabels, const size_t nPoints, Nodes &nodes) const
std::vector< Node > Nodes
int sortCoordinates(const float *pointCoordinates, const size_t nPoints, std::vector< size_t > &sortedIndices) const
int computeBranches(std::vector< Edges > &timeEdgesMap, std::vector< Nodes > &timeNodesMap) const
~TrackingFromOverlap() override=default
int computeOverlap(const float *pointCoordinates0, const float *pointCoordinates1, const labelType *pointLabels0, const labelType *pointLabels1, const size_t nPoints0, const size_t nPoints1, Edges &edges) const
int computeLabelIndexMap(const labelType *pointLabels, const size_t nPoints, std::map< labelType, size_t > &labelIndexMap) const
bool operator()(const size_t &i, const size_t &j)
CoordinateComparator(const float *coords)
const float * coordinates
printMsg(debug::output::GREEN+" "+debug::output::ENDCOLOR+debug::output::GREEN+"▒"+debug::output::ENDCOLOR+debug::output::GREEN+"▒▒▒▒▒▒▒▒▒▒▒▒▒░"+debug::output::ENDCOLOR, debug::Priority::PERFORMANCE, debug::LineMode::NEW, stream)