diff --git a/arbor/include/arbor/morph/morphology.hpp b/arbor/include/arbor/morph/morphology.hpp index 4285e3aec034b67f22f6876b8fc535a0d6dc4627..aa9ad61b88e0969acea52976327e7982c97a1548 100644 --- a/arbor/include/arbor/morph/morphology.hpp +++ b/arbor/include/arbor/morph/morphology.hpp @@ -128,5 +128,8 @@ private: mcable_list cables_; }; +// Reduced representation of an extent, excluding zero-length cables +// that are covered by more proximal or non-zero-length cables. +mcable_list canonical(const morphology& m, const mextent& a); } // namespace arb diff --git a/arbor/morph/morphology.cpp b/arbor/morph/morphology.cpp index ead854d28155bde26c0422cc34d4e392a0aaedd8..ad93e472873d9bd3efe6fa6c9e60aab56a67151a 100644 --- a/arbor/morph/morphology.cpp +++ b/arbor/morph/morphology.cpp @@ -481,5 +481,26 @@ mextent join(const mextent& a, const mextent& b) { return m; } +mcable_list canonical(const morphology& m, const mextent& a) { + // For zero-length cables representing isolated points, keep + // only the most proximal. All other zero-length cables should be + // elided. + + mcable_list result; + std::unordered_set<msize_t> remove_set; + + for (auto& c: a.cables()) { + if (c.prox_pos==0 && c.dist_pos>0) { + remove_set.insert(m.branch_parent(c.branch)); + } + } + for (auto& c: a.cables()) { + if (c.prox_pos==1 && remove_set.count(c.branch)) continue; + if (c.dist_pos==0 && m.branch_parent(c.branch)!=mnpos) continue; + result.push_back(c); + } + return result; +} + } // namespace arb diff --git a/test/unit/test_morphology.cpp b/test/unit/test_morphology.cpp index 5bf310274a867acbb928d1d3ae32135cfb0392f8..fd3bece9d4e26dbf4701f53bbe51cea0660116fd 100644 --- a/test/unit/test_morphology.cpp +++ b/test/unit/test_morphology.cpp @@ -734,3 +734,56 @@ TEST(mextent, intersects) { EXPECT_TRUE(x3.intersects(mlocation{3, 0.})); EXPECT_FALSE(x3.intersects(mlocation{3, 1.})); } + +TEST(mextent, canonical) { + using namespace arb; + using testing::cablelist_eq; + + using pvec = std::vector<msize_t>; + using svec = std::vector<msample>; + using cl = mcable_list; + + // Eight samples + // no-sphere + // sample branch + // 0 0 + // 1 3 0 1 + // 2 4 0 1 + // 5 6 2 3 + // 7 3 + pvec parents = {mnpos, 0, 1, 0, 3, 4, 4, 6}; + svec samples = { + {{ 0, 0, 0, 2}, 3}, + {{ 10, 0, 0, 2}, 3}, + {{100, 0, 0, 2}, 3}, + {{ 0, 10, 0, 2}, 3}, + {{ 0,100, 0, 2}, 3}, + {{100,100, 0, 2}, 3}, + {{ 0,200, 0, 2}, 3}, + {{ 0,300, 0, 2}, 3}, + }; + + morphology m(sample_tree(samples, parents)); + + mextent x1(m, cl{{0, 0.3, 0.7}, {1, 0.5, 0.7}, {3, 0.4, 1}}); + EXPECT_TRUE(cablelist_eq(canonical(m, x1), x1.cables())); + + mcable_list cl2{{1, 0.5, 1}, {2, 0, 0.5}}; // a reduced representation + mextent x2(m, cl2); + EXPECT_FALSE(cablelist_eq(cl2, x2.cables())); + EXPECT_TRUE(cablelist_eq(cl2, canonical(m, x2))); + + mcable_list cl3{{3, 0, 0}}; // extent is cover of fork point + mextent x3(m, cl3); + EXPECT_TRUE(cablelist_eq(cl{{1, 1, 1}}, canonical(m, x3))); + + mcable_list cl4{{1, 1, 1}, {3, 0, 0.5}}; // canonical repn excludes zero-length cable + mextent x4(m, cl4); + EXPECT_TRUE(cablelist_eq(cl{{3, 0, 0.5}}, canonical(m, x4))); + + // Heads of top-level branches should be preserved. + mcable_list cl5{{1, 0, 0}}; + mextent x5(m, cl5); + EXPECT_TRUE(cablelist_eq(cl{{0, 0, 0}, {1, 0, 0}}, x5.cables())); + EXPECT_TRUE(cablelist_eq(x5.cables(), canonical(m, x5))); +}