diff --git a/arbor/morph/place_pwlin.cpp b/arbor/morph/place_pwlin.cpp index 37df48d8c04578e26f60c2965215d49e37cefa11..8bec27ca6c0e7c50e023e5001bc8e3469eafe0c5 100644 --- a/arbor/morph/place_pwlin.cpp +++ b/arbor/morph/place_pwlin.cpp @@ -62,7 +62,16 @@ std::vector<mpoint> place_pwlin::all_at(mlocation loc) const { double pos = is_degenerate(pw_index)? 0: loc.pos; for (auto [bounds, index]: util::make_range(pw_index.equal_range(pos))) { - result.push_back(interpolate_segment(bounds, data_->segments.at(index), pos)); + auto seg = data_->segments.at(index); + + // Add both ends of zero length segment, if they differ. + if (bounds.first==bounds.second && seg.prox!=seg.dist) { + result.push_back(seg.prox); + result.push_back(seg.dist); + } + else { + result.push_back(interpolate_segment(bounds, seg, pos)); + } } return result; } diff --git a/arborio/include/arborio/neuroml.hpp b/arborio/include/arborio/neuroml.hpp index 401841d37093e0e4cda891785fbdd07de7d4981e..eead3b6c054659d547060b2d906c06e798cb4640 100644 --- a/arborio/include/arborio/neuroml.hpp +++ b/arborio/include/arborio/neuroml.hpp @@ -90,7 +90,19 @@ struct nml_morphology_data { struct neuroml_impl; +// TODO: C++20, replace with enum class and deploy using enum as appropriate. +namespace neuroml_options { + enum values { + none = 0, + allow_spherical_root = 1 + }; +} + struct neuroml { + // Correct interpretation of zero-length segments is currently a bit unclear + // in NeuroML 2.0, hence these options. For further options, use flags in powers of two + // so that we can bitwise combine and test them. + neuroml(); explicit neuroml(std::string nml_document); @@ -108,8 +120,8 @@ struct neuroml { // Parse and retrieve top-level morphology or morphology associated with a cell. // Return nullopt if not found. - std::optional<nml_morphology_data> morphology(const std::string& morph_id) const; - std::optional<nml_morphology_data> cell_morphology(const std::string& cell_id) const; + std::optional<nml_morphology_data> morphology(const std::string& morph_id, enum neuroml_options::values = neuroml_options::none) const; + std::optional<nml_morphology_data> cell_morphology(const std::string& cell_id, enum neuroml_options::values = neuroml_options::none) const; ~neuroml(); diff --git a/arborio/neuroml.cpp b/arborio/neuroml.cpp index 63e1bfc0943ce3cb047a2f147133be29a056a8a3..05d594c0132514f4322d79ab61f61f0f14624b87 100644 --- a/arborio/neuroml.cpp +++ b/arborio/neuroml.cpp @@ -116,15 +116,15 @@ std::vector<std::string> neuroml::morphology_ids() const { return result; } -optional<nml_morphology_data> neuroml::morphology(const std::string& morph_id) const { +optional<nml_morphology_data> neuroml::morphology(const std::string& morph_id, enum neuroml_options::values options) const { xml_error_scope err; auto ctx = impl_->make_context(); auto matches = ctx.query("//nml:neuroml/nml:morphology[@id="+xpath_escape(morph_id)+"]"); - return matches.empty()? nullopt: optional(nml_parse_morphology_element(ctx, matches[0])); + return matches.empty()? nullopt: optional(nml_parse_morphology_element(ctx, matches[0], options)); } -optional<nml_morphology_data> neuroml::cell_morphology(const std::string& cell_id) const { +optional<nml_morphology_data> neuroml::cell_morphology(const std::string& cell_id, enum neuroml_options::values options) const { xml_error_scope err; auto ctx = impl_->make_context(); auto matches = ctx.query( @@ -133,7 +133,7 @@ optional<nml_morphology_data> neuroml::cell_morphology(const std::string& cell_i if (matches.empty()) return nullopt; - nml_morphology_data M = nml_parse_morphology_element(ctx, matches[0]); + nml_morphology_data M = nml_parse_morphology_element(ctx, matches[0], options); M.cell_id = cell_id; return M; } diff --git a/arborio/nml_parse_morphology.cpp b/arborio/nml_parse_morphology.cpp index baff719a928b88cfde9e1c1d2bed1fb22897cfa2..24203ea18451fdcdb7ef3004d00bce8139429503 100644 --- a/arborio/nml_parse_morphology.cpp +++ b/arborio/nml_parse_morphology.cpp @@ -161,6 +161,7 @@ struct neuroml_segment { arb::mpoint distal; optional<non_negative> parent_id; double along = 1; + bool spherical = false; // Data for error reporting: unsigned line = 0; @@ -394,12 +395,30 @@ static arb::stitched_morphology construct_morphology(const neuroml_segment_tree& // Construct result from topologically sorted segments. - for (auto& s: segtree) { + for (const auto& s: segtree) { arb::mstitch stitch(nl_to_string(s.id), s.distal); - stitch.prox = s.proximal; + double along = s.along; + + if (s.spherical) { + arb_assert(!s.parent_id); // root segment only! + arb_assert(s.proximal && s.proximal.value()==s.distal); + + arb::mpoint centre = s.distal; + double radius = centre.radius; + + stitch.prox = arb::mpoint{centre.x, centre.y-radius, centre.z, radius}; + stitch.dist = arb::mpoint{centre.x, centre.y+radius, centre.z, radius}; + } + else if (s.parent_id && segtree[s.parent_id.value()].spherical) { + // Check if _parent_ is spherical. If so, we must attach to its centre. + along = 0.5; + } + else { + stitch.prox = s.proximal; + } if (s.parent_id) { - builder.add(stitch, nl_to_string(s.parent_id.value()), s.along); + builder.add(stitch, nl_to_string(s.parent_id.value()), along); } else { builder.add(stitch); @@ -409,7 +428,8 @@ static arb::stitched_morphology construct_morphology(const neuroml_segment_tree& return arb::stitched_morphology(std::move(builder)); } -nml_morphology_data nml_parse_morphology_element(xml_xpathctx ctx, xml_node morph) { +nml_morphology_data nml_parse_morphology_element(xml_xpathctx ctx, xml_node morph, enum neuroml_options::values options) { + using namespace neuroml_options; nml_morphology_data M; M.id = propx<std::string>(morph, "id", ""s); @@ -460,6 +480,10 @@ nml_morphology_data nml_parse_morphology_element(xml_xpathctx ctx, xml_node morp if (diameter<0) throw nml_bad_segment(seg.id, n.line()); seg.distal = arb::mpoint{x, y, z, diameter/2}; + + // Set spherical flag if we have no parent, options has allow_spherical_root flag, + // and proximal == distal. + seg.spherical = (options & allow_spherical_root) && !seg.parent_id && seg.proximal && seg.proximal.value()==seg.distal; } else { throw nml_bad_segment(seg.id, n.line()); diff --git a/arborio/nml_parse_morphology.hpp b/arborio/nml_parse_morphology.hpp index 3d4b5792c120d458d3c183873ff6a2793ba05240..80e88a0594a52331caf6b316c7d577e9162729a4 100644 --- a/arborio/nml_parse_morphology.hpp +++ b/arborio/nml_parse_morphology.hpp @@ -5,6 +5,6 @@ namespace arborio { -nml_morphology_data nml_parse_morphology_element(xmlwrap::xml_xpathctx ctx, xmlwrap::xml_node morph); +nml_morphology_data nml_parse_morphology_element(xmlwrap::xml_xpathctx ctx, xmlwrap::xml_node morph, enum neuroml_options::values); } // namespace arborio diff --git a/doc/cpp/morphology.rst b/doc/cpp/morphology.rst index 80ab858ef394c416c2111db3be834776491a5698..f70352c7330aa49147fa0237ae937c5825f36f77 100644 --- a/doc/cpp/morphology.rst +++ b/doc/cpp/morphology.rst @@ -482,7 +482,6 @@ Unhandleable exceptions from ``libxml2`` are forwarded via an exception NeuroML2 morphology support ^^^^^^^^^^^^^^^^^^^^^^^^^^^ - NeuroML documents are represented by the ``arborio::neuroml`` class, which in turn provides methods for the identification and translation of morphology data. ``neuroml`` objects are moveable and move-assignable, @@ -496,6 +495,9 @@ the underlying libxml2 library reports a problem that cannot be handled by the ` library. Otherwise, exceptions derived from ``aborio::neuroml_exception`` can be thrown when encountering problems interpreting the NeuroML document (see :ref:`cppneuromlexceptions` below). +Special parsing behaviour can be invoked through the use of an enum value in the `neuroml_options` +namespace. + .. cpp:class:: neuroml .. cpp:function:: neuroml(std::string) @@ -510,16 +512,30 @@ when encountering problems interpreting the NeuroML document (see :ref:`cppneuro Return the id of each top-level ``<morphology>`` element defined in the NeuroML document. - .. cpp:function:: std::optional<nml_morphology_data> morphology(const std::string&) const + .. cpp:function:: std::optional<nml_morphology_data> morphology(const std::string&, enum neuroml_options::value = neuroml_options::none) const Return a representation of the top-level morphology with the supplied identifier, or ``std::nullopt`` if no such morphology could be found. - .. cpp:function:: std::optional<nml_morphology_data> cell_morphology(const std::string&) const + .. cpp:function:: std::optional<nml_morphology_data> cell_morphology(const std::string&, enum neuroml_options::value = neuroml_options::none) const Return a representation of the morphology associated with the cell with the supplied identifier, or ``std::nullopt`` if the cell or its morphology could not be found. +.. cpp:enum:: neuroml_options::value + + .. cpp:enumerator:: none + + Perform no special parsing. + + .. cpp:enumerator:: allow_spherical_root + + Replace a zero-length root segment of constant radius with a Y-axis aligned + cylindrical segment of the same radius and with length twice the radius. This + cylinder will have the equivalent surface area to a sphere of the given radius. + + All child segments will connect to the centre of this cylinder, no matter the value of any ``fractionAlong`` attribute. + The morphology representation contains the corresponding Arbor ``arb::morphology`` object, label dictionaries for regions corresponding to its segments and segment groups by name and id, and a map providing the explicit list of segments contained within each defined diff --git a/doc/python/morphology.rst b/doc/python/morphology.rst index a68b4ee0c812ef3cdf1eb516b78bb88208ef1455..c8dfce3744e02093e8f4e8329159e56a22ebd45d 100644 --- a/doc/python/morphology.rst +++ b/doc/python/morphology.rst @@ -613,6 +613,10 @@ constitute part of the CV boundary point set. An implementation limitation restricts valid segment id values to those which can be represented by an unsigned long long value. + The ``allow_spherical_root`` optional parameter below, if set to true, will instruct the parser to + interpret a zero-length constant radius root segment as denoting a spherical segment, and this will + in turn be represented in the resultant morphology by a cylinder of equivalent surface area. + .. py:method:: neuroml(filename) Build a NeuroML document representation from the supplied file contents. @@ -631,20 +635,22 @@ constitute part of the CV boundary point set. :rtype: list[str] - .. py:method:: morphology(morph_id) + .. py:method:: morphology(morph_id, allow_spherical_root=false) Returns a representation of the top-level morphology with the supplied morph_id if it could be found. Parse errors or an inconsistent representation will raise an exception. :param str morph_id: ID of the top-level morphology. + :param bool allow_spherical_root: Treat zero length root segments especially. :rtype: optional(neuroml_morph_data) - .. py:method:: cell_morphology(cell_id) + .. py:method:: cell_morphology(cell_id, allow_spherical_root=false) Returns a representation of the morphology associated with the cell with the supplied cell_id if it could be found. Parse errors or an inconsistent representation will raise an exception. :param str morph_id: ID of the cell. + :param bool allow_spherical_root: Treat zero length root segments especially. :rtype: optional(neuroml_morph_data) .. _pyasc: diff --git a/python/morphology.cpp b/python/morphology.cpp index 8c736369df5e89485c381ae8c8d4aca13eb4bef8..3d67b41ae2dde89c0669125336cfc2857a83ee86 100644 --- a/python/morphology.cpp +++ b/python/morphology.cpp @@ -413,10 +413,10 @@ void register_morphology(py::module& m) { [](const arborio::nml_morphology_data& md) {return label_dict_proxy(md.segments);}, "Label dictionary containing one region expression for each segment id.") .def("named_segments", - [](const arborio::nml_morphology_data& md) {return label_dict_proxy(md.named_segments);}, + [](const arborio::nml_morphology_data& md) {return label_dict_proxy(md.named_segments);}, "Label dictionary containing one region expression for each name applied to one or more segments.") .def("groups", - [](const arborio::nml_morphology_data& md) {return label_dict_proxy(md.groups);}, + [](const arborio::nml_morphology_data& md) {return label_dict_proxy(md.groups);}, "Label dictionary containing one region expression for each segmentGroup id.") .def_readonly("group_segments", &arborio::nml_morphology_data::group_segments, @@ -427,59 +427,64 @@ void register_morphology(py::module& m) { neuroml // constructors .def(py::init( - [](std::string fname){ - std::ifstream fid{fname}; - if (!fid.good()) { - throw pyarb_error(util::pprintf("can't open file '{}'", fname)); - } - try { - std::string string_data((std::istreambuf_iterator<char>(fid)), - std::istreambuf_iterator<char>()); - return arborio::neuroml(string_data); - } - catch (arborio::neuroml_exception& e) { - // Try to produce helpful error messages for NeuroML parsing errors. - throw pyarb_error( - util::pprintf("NeuroML error processing file {}: ", fname, e.what())); - } + [](std::string fname) { + std::ifstream fid{fname}; + if (!fid.good()) { + throw pyarb_error(util::pprintf("can't open file '{}'", fname)); + } + try { + std::string string_data((std::istreambuf_iterator<char>(fid)), + std::istreambuf_iterator<char>()); + return arborio::neuroml(string_data); + } + catch (arborio::neuroml_exception& e) { + // Try to produce helpful error messages for NeuroML parsing errors. + throw pyarb_error(util::pprintf("NeuroML error processing file {}: {}", fname, e.what())); + } })) .def("cell_ids", - [](const arborio::neuroml& nml) { + [](const arborio::neuroml& nml) { try { return nml.cell_ids(); } catch (arborio::neuroml_exception& e) { throw util::pprintf("NeuroML error: {}", e.what()); } - },"Query top-level cells.") + }, + "Query top-level cells.") .def("morphology_ids", - [](const arborio::neuroml& nml) { + [](const arborio::neuroml& nml) { try { return nml.morphology_ids(); } catch (arborio::neuroml_exception& e) { throw util::pprintf("NeuroML error: {}", e.what()); } - },"Query top-level standalone morphologies.") + }, + "Query top-level standalone morphologies.") .def("morphology", - [](const arborio::neuroml& nml, const std::string& morph_id) { + [](const arborio::neuroml& nml, const std::string& morph_id, bool spherical) { try { - return nml.morphology(morph_id); + using namespace arborio::neuroml_options; + return nml.morphology(morph_id, spherical? allow_spherical_root: none); } catch (arborio::neuroml_exception& e) { throw util::pprintf("NeuroML error: {}", e.what()); } - },"morph_id"_a,"Retrieve top-level nml_morph_data associated with morph_id.") + }, "morph_id"_a, "allow_spherical_root"_a=false, + "Retrieve top-level nml_morph_data associated with morph_id.") .def("cell_morphology", - [](const arborio::neuroml& nml, const std::string& cell_id) { - try { - return nml.cell_morphology(cell_id); - } - catch (arborio::neuroml_exception& e) { - throw util::pprintf("NeuroML error: {}", e.what()); - } - },"morph_id"_a,"Retrieve nml_morph_data associated with cell_id."); -#endif + [](const arborio::neuroml& nml, const std::string& cell_id, bool spherical) { + try { + using namespace arborio::neuroml_options; + return nml.cell_morphology(cell_id, spherical? allow_spherical_root: none); + } + catch (arborio::neuroml_exception& e) { + throw util::pprintf("NeuroML error: {}", e.what()); + } + }, "cell_id"_a, "allow_spherical_root"_a=false, + "Retrieve nml_morph_data associated with cell_id."); +#endif // def ARB_NEUROML_ENABLED } } // namespace pyarb diff --git a/test/unit/morph_pred.hpp b/test/unit/morph_pred.hpp index c66e85cf9c5e96d30258b4a07a911244c73594a7..7bb644df5d4fae1ea4e73381027b674868a69173 100644 --- a/test/unit/morph_pred.hpp +++ b/test/unit/morph_pred.hpp @@ -10,6 +10,7 @@ #include <arbor/morph/region.hpp> #include "util/span.hpp" +#include "common.hpp" namespace testing { @@ -85,6 +86,19 @@ inline ::testing::AssertionResult locset_eq(const arb::mprovider& p, arb::locset return mlocationlist_eq(thingify(a, p), thingify(b, p)); } +inline ::testing::AssertionResult point_almost_eq(arb::mpoint p, arb::mpoint q) { + // Compare almost equal in single precision. + auto result = almost_eq<float>(p.x, q.x); + if (!result) return result; + + result = almost_eq<float>(p.y, q.y); + if (!result) return result; + + result = almost_eq<float>(p.z, q.z); + if (!result) return result; + + return almost_eq<float>(p.radius, q.radius); +} } // namespace testing diff --git a/test/unit/test_morph_place.cpp b/test/unit/test_morph_place.cpp index 6d02a26e73cf87ed60ce10d6c3b1218f3194cda2..e7869fa705f5aca08871550b1e545e17ca2404b8 100644 --- a/test/unit/test_morph_place.cpp +++ b/test/unit/test_morph_place.cpp @@ -385,6 +385,35 @@ TEST(place_pwlin, all_at) { EXPECT_TRUE(mpoint_almost_eq(p1d, points_end_b0[1])); EXPECT_TRUE(mpoint_almost_eq(p2d, points_end_b0[2])); } + + // Zero length branch comprising single zero-length segment with differing radius. + // (Please don't do this either.) + { + mpoint p0p{2, 3, 4, 5}; + mpoint p0d{2, 3, 4, 8}; + + segment_tree tree; + (void)tree.append(mnpos, p0p, p0d, 0); + + morphology m(tree); + mprovider p(m, label_dict{}); + place_pwlin place(m); + + auto points_begin_b0 = place.all_at(mlocation{0, 0}); + ASSERT_EQ(2u, points_begin_b0.size()); + EXPECT_TRUE(mpoint_almost_eq(p0p, points_begin_b0[0])); + EXPECT_TRUE(mpoint_almost_eq(p0d, points_begin_b0[1])); + + auto points_mid_b0 = place.all_at(mlocation{0, 0.5}); + ASSERT_EQ(2u, points_begin_b0.size()); + EXPECT_TRUE(mpoint_almost_eq(p0p, points_begin_b0[0])); + EXPECT_TRUE(mpoint_almost_eq(p0d, points_begin_b0[1])); + + auto points_end_b0 = place.all_at(mlocation{0, 1}); + ASSERT_EQ(2u, points_begin_b0.size()); + EXPECT_TRUE(mpoint_almost_eq(p0p, points_begin_b0[0])); + EXPECT_TRUE(mpoint_almost_eq(p0d, points_begin_b0[1])); + } } TEST(place_pwlin, segments) { diff --git a/test/unit/test_nml_morphology.cpp b/test/unit/test_nml_morphology.cpp index ff69fe25549bbb3ae70e4433340dbe96fb9f8179..28e5af58cd55ad3b65b9137cdf7a3fa08c471dc6 100644 --- a/test/unit/test_nml_morphology.cpp +++ b/test/unit/test_nml_morphology.cpp @@ -277,6 +277,165 @@ R"~( } } +TEST(neuroml, spherical_segments) { + using namespace arb; + using namespace arborio::neuroml_options; + + // Spherical root segments can be translated as equivalent-area + // cylinders oriented along the y-axis in the generated morphology. + + // Points used in morphology definitions below. + + mpoint p0{1, -2, 3.5, 4}; + mpoint p1{1, -2, 3.5, 3}; + mpoint p2{3, -3, 4.5, 5}; + + std::string doc = +R"~( +<neuroml xmlns="http://www.neuroml.org/schema/neuroml2"> +<morphology id="m1"> + <!-- Single zero-length segment at p0 --> + <segment id="0"> + <proximal x="1" y="-2" z="3.5" diameter="8"/> + <distal x="1" y="-2" z="3.5" diameter="8"/> + </segment> +</morphology> +<morphology id="m2"> + <!-- Single zero-length segment defined between colocated p0 and p1; + diameters differ, so should not be treated as spherical. --> + <segment id="0"> + <proximal x="1" y="-2" z="3.5" diameter="8"/> + <distal x="1" y="-2" z="3.5" diameter="6"/> + </segment> +</morphology> +<morphology id="m3"> + <!-- Two segments: first is zero-length with ends colocated at p0; + second has distal point p2, and attaches at fractionAlong=0.2, but should inherit proximal point p0. --> + <segment id="0"> + <proximal x="1" y="-2" z="3.5" diameter="8"/> + <distal x="1" y="-2" z="3.5" diameter="8"/> + </segment> + <segment id="1"> + <parent segment="0" fractionAlong="0.2"/> + <distal x="3" y="-3" z="4.5" diameter="10"/> + </segment> +</morphology> +<morphology id="m4"> + <!-- Two segments: first is between p1 and p2; second is from p2 to p2, but is not the root segment + and should not be interpreted as a sphere. --> + <segment id="0"> + <proximal x="1" y="-2" z="3.5" diameter="6"/> + <distal x="3" y="-3" z="4.5" diameter="10"/> + </segment> + <segment id="1"> + <parent segment="0"/> + <proximal x="3" y="-3" z="4.5" diameter="10"/> + <distal x="3" y="-3" z="4.5" diameter="10"/> + </segment> +</morphology> +</neuroml> +)~"; + + arborio::neuroml N(doc); + + { + arborio::nml_morphology_data m1 = N.morphology("m1", allow_spherical_root).value(); + label_dict labels; + labels.import(m1.segments, "seg:"); + mprovider P(m1.morphology, labels); + + EXPECT_TRUE(region_eq(P, reg::branch(0), reg::all())); + EXPECT_TRUE(region_eq(P, reg::named("seg:0"), reg::all())); + + place_pwlin G(P.morphology()); + EXPECT_EQ(p0.radius, G.at(mlocation{0, 0}).radius); + EXPECT_EQ(p0.radius, G.at(mlocation{0, 1}).radius); + + mpoint centre = G.at(mlocation{0, 0.5}); + EXPECT_EQ(p0, centre); + + // Only y-axis points should differ from centre. + mpoint l0 = G.at(mlocation{0, 0}); + mpoint l1 = G.at(mlocation{0, 1}); + EXPECT_EQ(p0.x, l0.x); + EXPECT_NE(p0.y, l0.y); + EXPECT_EQ(p0.z, l0.z); + EXPECT_EQ(p0.x, l1.x); + EXPECT_NE(p0.y, l1.y); + EXPECT_EQ(p0.z, l1.z); + + EXPECT_DOUBLE_EQ(2*p0.radius, P.embedding().branch_length(0)); + } + { + // With spherical root _not_ provided, treat it just as a simple zero-length segment. + arborio::nml_morphology_data m1 = N.morphology("m1", none).value(); + label_dict labels; + labels.import(m1.segments, "seg:"); + mprovider P(m1.morphology, labels); + + EXPECT_TRUE(region_eq(P, reg::branch(0), reg::all())); + EXPECT_TRUE(region_eq(P, reg::named("seg:0"), reg::all())); + + place_pwlin G(P.morphology()); + EXPECT_EQ(p0, G.at(mlocation{0, 0})); + EXPECT_EQ(p0, G.at(mlocation{0, 1})); + } + { + arborio::nml_morphology_data m2 = N.morphology("m2", allow_spherical_root).value(); + label_dict labels; + labels.import(m2.segments, "seg:"); + mprovider P(m2.morphology, labels); + + EXPECT_TRUE(region_eq(P, reg::branch(0), reg::all())); + EXPECT_TRUE(region_eq(P, reg::named("seg:0"), reg::all())); + + // This one shouldn't be interpreted as a sphere. + place_pwlin G(P.morphology()); + auto points = G.all_at(mlocation{0, 0}); + ASSERT_EQ(2u, points.size()); + EXPECT_TRUE((p0==points[0] && p1==points[1]) || + (p0==points[1] && p1==points[0])); + } + { + arborio::nml_morphology_data m3 = N.morphology("m3", allow_spherical_root).value(); + label_dict labels; + labels.import(m3.segments, "seg:"); + mprovider P(m3.morphology, labels); + place_pwlin G(P.morphology()); + + // With segment 1 attached to spherical segment 0, we should have three branches. + EXPECT_EQ(3u, P.morphology().num_branches()); + + mlocation s0centre = thingify(ls::on_components(0.5, reg::named("seg:0")), P).at(0); + EXPECT_EQ(p0, G.at(s0centre)); + + // Compute locations of ends of segment 1. + mlocation s1ploc = thingify(ls::most_proximal(reg::named("seg:1")), P).at(0); + mlocation s1dloc = thingify(ls::most_distal(reg::named("seg:1")), P).at(0); + mpoint s1p = G.at(s1ploc); + mpoint s1d = G.at(s1dloc); + + EXPECT_TRUE(testing::point_almost_eq(p0, s1p)); + EXPECT_EQ(p2, s1d); + } + { + arborio::nml_morphology_data m4 = N.morphology("m4", allow_spherical_root).value(); + label_dict labels; + labels.import(m4.segments, "seg:"); + mprovider P(m4.morphology, labels); + place_pwlin G(P.morphology()); + + // Segment 1 should have colocated, equal radius endpoints, as it is not the root segment. + mlocation s1ploc = thingify(ls::most_proximal(reg::named("seg:1")), P).at(0); + mlocation s1dloc = thingify(ls::most_distal(reg::named("seg:1")), P).at(0); + mpoint s1p = G.at(s1ploc); + mpoint s1d = G.at(s1dloc); + + EXPECT_EQ(p2, s1p); + EXPECT_EQ(p2, s1d); + } +} + TEST(neuroml, segment_errors) { using namespace arb;