diff --git a/arbor/include/arbor/swcio.hpp b/arbor/include/arbor/swcio.hpp index efcd17fa0a97e3591bbf0389cd84e068002b7872..db5b672f2cfdf39f576699ba39f734d4b9db9944 100644 --- a/arbor/include/arbor/swcio.hpp +++ b/arbor/include/arbor/swcio.hpp @@ -33,11 +33,6 @@ struct swc_duplicate_record_id: swc_error { explicit swc_duplicate_record_id(int record_id); }; -// Irregular record ordering. -struct swc_irregular_id: swc_error { - explicit swc_irregular_id(int record_id); -}; - // Smells like a spherical soma. struct swc_spherical_soma: swc_error { explicit swc_spherical_soma(int record_id); @@ -88,25 +83,22 @@ struct swc_data { // in comments (stripping initial '#' and subsequent whitespace). // Stops at EOF or after reading the first line that does not parse as SWC. // -// Note that 'one-point soma' SWC files are explicitly not supported. -// // In `relaxed` mode, it will check that: // * There are no duplicate record ids. // * All record ids are positive. // * There are no records whose parent id is not less than the record id. // * Only one record has parent id -1; all other parent ids correspond to records. -// * There are at least two records. // -// In `strict` mode, it will additionally check: -// * Record ids are numbered contiguously from 1. -// * The data cannot be interpreted as a 'spherical soma' SWC file. -// Specifically, the root record shares its tag with at least one other -// record with has the root as parent. +// In `strict` mode, it will additionally check that the data cannot be interpreted +// as a 'spherical soma' SWC file: +// * The root record must share its tag with at least one other record +// which has the root as parent. This implies that there must be at least +// two SWC records. // // Throws a corresponding exception of type derived from `swc_error` if any of the // conditions above are encountered. // -// SWC records are stored in id order. +// SWC records are returned in id order. enum class swc_mode { relaxed, strict }; @@ -118,6 +110,14 @@ swc_data parse_swc(const std::string& text, swc_mode mode = swc_mode::strict); swc_data parse_swc(std::vector<swc_record>, swc_mode = swc_mode::strict); // Convert a valid, ordered sequence of SWC records to a morphological segment tree. +// +// Note that 'one-point soma' SWC files are explicitly not supported; the swc_data +// is expected to abide by the restrictions of `strict` mode parsing as described +// above. +// +// The generated segment tree will be contiguous. There will be one segment for +// each SWC record after the first: this record defines the tag and distal point +// of the segment, while the proximal point is taken from the parent record. segment_tree as_segment_tree(const std::vector<swc_record>&); diff --git a/arbor/swcio.cpp b/arbor/swcio.cpp index 1703f1a7e902acaaa00ce668a98ee74aacda0c1c..c17b375ecf14e226f81083a5308be2fc5ee75558 100644 --- a/arbor/swcio.cpp +++ b/arbor/swcio.cpp @@ -33,10 +33,6 @@ swc_duplicate_record_id::swc_duplicate_record_id(int record_id): swc_error("duplicate SWC sample id", record_id) {} -swc_irregular_id::swc_irregular_id(int record_id): - swc_error("SWC record id not numbered consecutively from 1", record_id) -{} - swc_spherical_soma::swc_spherical_soma(int record_id): swc_error("SWC with spherical somata are not supported", record_id) {} @@ -94,10 +90,6 @@ static std::vector<swc_record> sort_and_validate_swc(std::vector<swc_record> rec throw swc_record_precedes_parent(r.id); } - if (mode==swc_mode::strict && r.id != (int)i+1) { - throw swc_irregular_id(r.id); - } - if (!seen.insert(r.id).second) { throw swc_duplicate_record_id(r.id); } diff --git a/test/unit/test_swcio.cpp b/test/unit/test_swcio.cpp index ab890246d7474d7080245d34ba9f40f07152f18b..debcb072109ed22912b6bd68b3ed9b0d497bed93 100644 --- a/test/unit/test_swcio.cpp +++ b/test/unit/test_swcio.cpp @@ -4,6 +4,8 @@ #include <sstream> #include <arbor/cable_cell.hpp> +#include <arbor/morph/primitives.hpp> +#include <arbor/morph/segment_tree.hpp> #include <arbor/swcio.hpp> #include "../gtest.h" @@ -110,16 +112,29 @@ TEST(swc_parser, bad_relaxed) { EXPECT_THROW(parse_swc(bad5, swc_mode::relaxed), swc_no_such_parent); } +} +TEST(swc_parser, bad_strict) { { std::string bad6 = - "1 1 0.1 0.2 0.3 0.4 -1\n"; + "1 7 0.1 0.2 0.3 0.4 -1\n"; // just one record EXPECT_THROW(parse_swc(bad6, swc_mode::relaxed), swc_spherical_soma); } + { + std::string bad3 = + "1 4 0.1 0.2 0.3 0.4 -1\n" // solitary tag + "2 6 0.1 0.2 0.3 0.4 1\n" + "3 6 0.1 0.2 0.3 0.4 2\n" + "4 6 0.1 0.2 0.3 0.4 1\n"; + + EXPECT_THROW(parse_swc(bad3, swc_mode::strict), swc_spherical_soma); + EXPECT_NO_THROW(parse_swc(bad3, swc_mode::relaxed)); + } } -TEST(swc_parser, bad_strict) { +TEST(swc_parser, valid_relaxed) { + // Non-contiguous is okay. { std::string bad1 = "1 1 0.1 0.2 0.3 0.4 -1\n" @@ -127,10 +142,10 @@ TEST(swc_parser, bad_strict) { "3 1 0.1 0.2 0.3 0.4 2\n" "5 1 0.1 0.2 0.3 0.4 3\n"; // non-contiguous - EXPECT_THROW(parse_swc(bad1, swc_mode::strict), swc_irregular_id); EXPECT_NO_THROW(parse_swc(bad1, swc_mode::relaxed)); } + // As is out of order. { std::string bad2 = "1 1 0.1 0.2 0.3 0.4 -1\n" @@ -138,23 +153,9 @@ TEST(swc_parser, bad_strict) { "2 1 0.1 0.2 0.3 0.4 1\n" "4 1 0.1 0.2 0.3 0.4 3\n"; - EXPECT_THROW(parse_swc(bad2, swc_mode::strict), swc_irregular_id); EXPECT_NO_THROW(parse_swc(bad2, swc_mode::relaxed)); } - { - std::string bad3 = - "1 1 0.1 0.2 0.3 0.4 -1\n" // solitary tag - "2 0 0.1 0.2 0.3 0.4 1\n" - "3 0 0.1 0.2 0.3 0.4 2\n" - "4 0 0.1 0.2 0.3 0.4 1\n"; - - EXPECT_THROW(parse_swc(bad3, swc_mode::strict), swc_spherical_soma); - EXPECT_NO_THROW(parse_swc(bad3, swc_mode::relaxed)); - } -} - -TEST(swc_parser, valid_relaxed) { } TEST(swc_parser, valid_strict) { @@ -169,12 +170,13 @@ TEST(swc_parser, valid_strict) { } { + // Non-contiguous, out of order records are fine. std::string valid2 = "# Some people put\n" "# <xml /> in here!\n" "1 1 0.1 0.2 0.3 0.4 -1\n" "2 1 0.3 0.4 0.5 0.3 1\n" - "3 2 0.2 0.6 0.8 0.2 2\n" + "5 2 0.2 0.6 0.8 0.2 2\n" "4 0 0.2 0.8 0.6 0.3 2"; swc_data data = parse_swc(valid2, swc_mode::strict); @@ -182,8 +184,8 @@ TEST(swc_parser, valid_strict) { ASSERT_EQ(4u, data.records.size()); EXPECT_EQ(swc_record(1, 1, 0.1, 0.2, 0.3, 0.4, -1), data.records[0]); EXPECT_EQ(swc_record(2, 1, 0.3, 0.4, 0.5, 0.3, 1), data.records[1]); - EXPECT_EQ(swc_record(3, 2, 0.2, 0.6, 0.8, 0.2, 2), data.records[2]); - EXPECT_EQ(swc_record(4, 0, 0.2, 0.8, 0.6, 0.3, 2), data.records[3]); + EXPECT_EQ(swc_record(4, 0, 0.2, 0.8, 0.6, 0.3, 2), data.records[2]); + EXPECT_EQ(swc_record(5, 2, 0.2, 0.6, 0.8, 0.2, 2), data.records[3]); // Trailing garbage is ignored in data records. std::string valid3 = @@ -199,6 +201,63 @@ TEST(swc_parser, valid_strict) { } } +TEST(swc_parser, segment_tree) { + { + // Missing parent record will throw. + std::vector<swc_record> swc{ + {1, 1, 0., 0., 0., 1., -1}, + {5, 3, 1., 1., 1., 1., 2} + }; + EXPECT_THROW(as_segment_tree(swc), bad_swc_data); + } + { + // A single SWC record will throw. + std::vector<swc_record> swc{ + {1, 1, 0., 0., 0., 1., -1} + }; + EXPECT_THROW(as_segment_tree(swc), bad_swc_data); + } + { + // Otherwise, ensure segment ends and tags correspond. + mpoint p0{0.1, 0.2, 0.3, 0.4}; + mpoint p1{0.3, 0.4, 0.5, 0.3}; + mpoint p2{0.2, 0.8, 0.6, 0.3}; + mpoint p3{0.2, 0.6, 0.8, 0.2}; + mpoint p4{0.4, 0.5, 0.5, 0.1}; + + std::vector<swc_record> swc{ + {1, 1, p0.x, p0.y, p0.z, p0.radius, -1}, + {2, 1, p1.x, p1.y, p1.z, p1.radius, 1}, + {4, 3, p2.x, p2.y, p2.z, p2.radius, 2}, + {5, 2, p3.x, p3.y, p3.z, p3.radius, 2}, + {7, 3, p4.x, p4.y, p4.z, p4.radius, 4} + }; + + segment_tree tree = as_segment_tree(swc); + ASSERT_EQ(4u, tree.segments().size()); + + EXPECT_EQ(mnpos, tree.parents()[0]); + EXPECT_EQ(1, tree.segments()[0].tag); + EXPECT_EQ(p0, tree.segments()[0].prox); + EXPECT_EQ(p1, tree.segments()[0].dist); + + EXPECT_EQ(0u, tree.parents()[1]); + EXPECT_EQ(3, tree.segments()[1].tag); + EXPECT_EQ(p1, tree.segments()[1].prox); + EXPECT_EQ(p2, tree.segments()[1].dist); + + EXPECT_EQ(0u, tree.parents()[2]); + EXPECT_EQ(2, tree.segments()[2].tag); + EXPECT_EQ(p1, tree.segments()[2].prox); + EXPECT_EQ(p3, tree.segments()[2].dist); + + EXPECT_EQ(1u, tree.parents()[3]); + EXPECT_EQ(3, tree.segments()[3].tag); + EXPECT_EQ(p2, tree.segments()[3].prox); + EXPECT_EQ(p4, tree.segments()[3].dist); + } +} + // hipcc bug in reading DATADIR #ifndef ARB_HIP TEST(swc_parser, from_neuromorpho)