diff --git a/cc3d.hpp b/cc3d.hpp index 5fb4d74..d813d3d 100644 --- a/cc3d.hpp +++ b/cc3d.hpp @@ -918,36 +918,180 @@ OUT* connected_components2d_4( // Raster Scan 1: Set temporary labels and // record equivalences in a disjoint set. +#define NEW_LABEL() \ + next_label++; \ + out_labels[loc] = next_label; \ + equivalences.add(out_labels[loc]); + +#define INIT_TREE() \ + if (++x >= xend - 1) { goto FINAL_TREE; } \ + loc = x + sx * y; \ + cur = in_labels[loc]; \ + next = in_labels[loc + 1]; + T cur = 0; + T next = 0; for (int64_t y = 0; y < sy; y++, row++) { const int64_t xstart = runs[row << 1]; const int64_t xend = runs[(row << 1) + 1]; - for (int64_t x = xstart; x < xend; x++) { - loc = x + sx * y; - cur = in_labels[loc]; + int64_t x = xstart - 1; + goto NO_B; - if (cur == 0) { - continue; + FULLTREE: + INIT_TREE() + if (cur == 0) { + goto NO_B; + } + else if (x > 0 && cur == in_labels[loc + B]) { + out_labels[loc + A] = out_labels[loc + B]; + if (y > 0 && cur == in_labels[loc + C]) { + if (cur != in_labels[loc + D]) { + equivalences.unify(out_labels[loc + A], out_labels[loc + C]); + } + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; + } + } + else if (cur == next) { + goto NO_D; + } + else { + goto NO_B; + } + } + else if (y > 0 && cur == in_labels[loc + C]) { + out_labels[loc + A] = out_labels[loc + C]; + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; } + } + else { + NEW_LABEL() + goto NO_D; + } - if (x > 0 && cur == in_labels[loc + B]) { - out_labels[loc + A] = out_labels[loc + B]; - if (y > 0 && cur != in_labels[loc + D] && cur == in_labels[loc + C]) { - equivalences.unify(out_labels[loc + A], out_labels[loc + C]); + NO_B: + INIT_TREE() + if (cur == 0) { + goto NO_B; + } + else if (y > 0 && cur == in_labels[loc + C]) { + out_labels[loc + A] = out_labels[loc + C]; + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; + } + } + else { + NEW_LABEL() + goto NO_D; + } + + ASSUME_BD: + INIT_TREE() + if (cur == 0) { + goto NO_B; + } + else { + out_labels[loc + A] = out_labels[loc + B]; + if (cur == next) { + goto ASSUME_B; + } + goto NO_B; + } + + ASSUME_B: + INIT_TREE() + out_labels[loc + A] = out_labels[loc + B]; + if (y > 0 && cur == in_labels[loc + C]) { + if (cur != in_labels[loc + D]) { + equivalences.unify(out_labels[loc + A], out_labels[loc + C]); + } + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; + } + } + else if (cur == next) { + goto NO_D; + } + else { + goto NO_B; + } + + NO_D: + INIT_TREE() + if (cur == 0) { + goto NO_B; + } + else if (x > 0 && cur == in_labels[loc + B]) { + out_labels[loc + A] = out_labels[loc + B]; + if (y > 0 && cur == in_labels[loc + C]) { + equivalences.unify(out_labels[loc + A], out_labels[loc + C]); + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; } } - else if (y > 0 && cur == in_labels[loc + C]) { - out_labels[loc + A] = out_labels[loc + C]; + else if (cur == next) { + goto NO_D; } else { - next_label++; - out_labels[loc + A] = next_label; - equivalences.add(out_labels[loc + A]); + goto NO_B; + } + } + else if (y > 0 && cur == in_labels[loc + C]) { + out_labels[loc + A] = out_labels[loc + C]; + if (cur == next) { + goto ASSUME_BD; + } + else { + goto NO_B; + } + } + else { + NEW_LABEL() + goto NO_D; + } + + FINAL_TREE: + if (x >= xend) { continue; } + loc = x + sx * y; + cur = in_labels[loc]; + + if (cur == 0) { + continue; + } + else if (x > 0 && cur == in_labels[loc + B]) { + out_labels[loc + A] = out_labels[loc + B]; + if (y > 0 && cur != in_labels[loc + D] && cur == in_labels[loc + C]) { + equivalences.unify(out_labels[loc + A], out_labels[loc + C]); } } + else if (y > 0 && cur == in_labels[loc + C]) { + out_labels[loc + A] = out_labels[loc + C]; + } + else { + NEW_LABEL() + } } +#undef NEW_LABEL +#undef INIT_TREE + out_labels = relabel(out_labels, sx, sy, /*sz=*/1, next_label, equivalences, N, runs); delete[] runs; return out_labels; diff --git a/test.cpp b/test.cpp index 677b386..b00e59e 100644 --- a/test.cpp +++ b/test.cpp @@ -39,10 +39,7 @@ uint32_t* read_subvol() { return input; } -uint32_t* randomvol() { - size_t sx = 512; - size_t sy = 512; - size_t sz = 512; +uint32_t* randomvol(size_t sx = 512, size_t sy = 512, size_t sz = 512) { size_t voxels = sx * sy * sz; uint32_t *big = new uint32_t[sx*sy*sz](); @@ -52,23 +49,18 @@ uint32_t* randomvol() { return big; } -uint32_t* black() { - size_t sx = 512; - size_t sy = 512; - size_t sz = 512; - size_t voxels = sx * sy * sz; - +uint32_t* black(size_t sx = 512, size_t sy = 512, size_t sz = 512) { return new uint32_t[sx*sy*sz](); } int main () { - size_t sx = 512; - size_t sy = 512; - size_t sz = 512; + size_t sx = 10000; + size_t sy = 10000; + size_t sz = 1; size_t voxels = sx * sy * sz; - uint32_t* subvol = black(); + uint32_t* subvol = randomvol(sx,sy,sz); typedef std::chrono::high_resolution_clock Time; typedef std::chrono::milliseconds ms; @@ -76,7 +68,7 @@ int main () { auto t0 = Time::now(); int N = 1; for (int i = 0; i < N; i++) { - cc3d::connected_components3d(subvol, sx,sy,sz, 26); + cc3d::connected_components3d(subvol, sx,sy,sz, 4); } auto t1 = Time::now();