diff --git a/main.cpp b/main.cpp
index 1f70eac1836949b320a1c01d7f8baa0eb16126db..ed8c3307699c64dc2ed19860d54de7739cd88e4e 100644
--- a/main.cpp
+++ b/main.cpp
@@ -1,10 +1,12 @@
#include <iostream>
#include <fstream>
#include <numeric>
+#include <vector>
#include "gtest/gtest.h"
#include "cell_tree.hpp"
+#include "swcio.hpp"
#include "json/src/json.hpp"
using json = nlohmann::json;
@@ -250,6 +252,179 @@ TEST(cell_tree, json_load) {
}
}
+// SWC tests
+void expect_cell_equals(const neuron::io::cell_record &expected,
+ const neuron::io::cell_record &actual)
+{
+ EXPECT_EQ(expected.id(), actual.id());
+ EXPECT_EQ(expected.type(), actual.type());
+ EXPECT_FLOAT_EQ(expected.x(), actual.x());
+ EXPECT_FLOAT_EQ(expected.y(), actual.y());
+ EXPECT_FLOAT_EQ(expected.z(), actual.z());
+ EXPECT_FLOAT_EQ(expected.radius(), actual.radius());
+ EXPECT_EQ(expected.parent(), actual.parent());
+}
+
+TEST(cell_record, construction)
+{
+ using namespace neuron::io;
+
+ {
+ // force an invalid type
+ cell_record::kind invalid_type = static_cast<cell_record::kind>(100);
+ EXPECT_THROW(cell_record cell(invalid_type, 7, 1., 1., 1., 1., 5),
+ std::invalid_argument);
+ }
+
+ {
+ // invalid id
+ EXPECT_THROW(cell_record cell(
+ cell_record::custom, -3, 1., 1., 1., 1., 5),
+ std::invalid_argument);
+ }
+
+ {
+ // invalid parent id
+ EXPECT_THROW(cell_record cell(
+ cell_record::custom, 0, 1., 1., 1., 1., -5),
+ std::invalid_argument);
+ }
+
+ {
+ // invalid radius
+ EXPECT_THROW(cell_record cell(
+ cell_record::custom, 0, 1., 1., 1., -1., -1),
+ std::invalid_argument);
+ }
+
+ {
+ // parent_id > id
+ EXPECT_THROW(cell_record cell(
+ cell_record::custom, 0, 1., 1., 1., 1., 2),
+ std::invalid_argument);
+ }
+
+ {
+ // parent_id == id
+ EXPECT_THROW(cell_record cell(
+ cell_record::custom, 0, 1., 1., 1., 1., 0),
+ std::invalid_argument);
+ }
+
+ {
+ // check standard construction by value
+ cell_record cell(cell_record::custom, 0, 1., 1., 1., 1., -1);
+ EXPECT_EQ(cell.id(), 0);
+ EXPECT_EQ(cell.type(), cell_record::custom);
+ EXPECT_EQ(cell.x(), 1.);
+ EXPECT_EQ(cell.y(), 1.);
+ EXPECT_EQ(cell.z(), 1.);
+ EXPECT_EQ(cell.radius(), 1.);
+ EXPECT_EQ(cell.diameter(), 2*1.);
+ EXPECT_EQ(cell.parent(), -1);
+ }
+
+ {
+ // check copy constructor
+ cell_record cell_orig(cell_record::custom, 0, 1., 1., 1., 1., -1);
+ cell_record cell(cell_orig);
+ expect_cell_equals(cell_orig, cell);
+ }
+}
+
+TEST(swc_parser, invalid_input)
+{
+ using namespace neuron::io;
+
+ {
+ // check incomplete lines; missing parent
+ std::istringstream is("1 1 14.566132 34.873772 7.857000 0.717830\n");
+ cell_record cell;
+ EXPECT_THROW(is >> cell, std::logic_error);
+ }
+
+ {
+ // Check non-parsable values
+ std::istringstream is("1a 1 14.566132 34.873772 7.857000 0.717830 -1\n");
+ cell_record cell;
+ EXPECT_THROW(is >> cell, std::logic_error);
+ }
+
+ {
+ // Check invalid cell type
+ std::istringstream is("1 10 14.566132 34.873772 7.857000 0.717830 -1\n");
+ cell_record cell;
+ EXPECT_THROW(is >> cell, std::invalid_argument);
+ }
+}
+
+
+TEST(swc_parser, valid_input)
+{
+ using namespace neuron::io;
+
+ {
+ // check empty file; no cell may be parsed
+ cell_record cell, cell_orig;
+ std::istringstream is("");
+ EXPECT_NO_THROW(is >> cell);
+ expect_cell_equals(cell_orig, cell);
+ }
+
+ {
+ // check comment-only file not ending with a newline;
+ // no cell may be parsed
+ cell_record cell, cell_orig;
+ std::istringstream is("#comment\n#comment");
+ EXPECT_NO_THROW(is >> cell);
+ expect_cell_equals(cell_orig, cell);
+ }
+
+
+ {
+ // check last line case (no newline at the end)
+ std::istringstream is("1 1 14.566132 34.873772 7.857000 0.717830 -1");
+ cell_record cell;
+ EXPECT_NO_THROW(is >> cell);
+ EXPECT_EQ(0, cell.id()); // zero-based indexing
+ EXPECT_EQ(cell_record::soma, cell.type());
+ EXPECT_FLOAT_EQ(14.566132, cell.x());
+ EXPECT_FLOAT_EQ(34.873772, cell.y());
+ EXPECT_FLOAT_EQ( 7.857000, cell.z());
+ EXPECT_FLOAT_EQ( 0.717830, cell.radius());
+ EXPECT_FLOAT_EQ( -1, cell.parent());
+ }
+
+ {
+ // check valid input with a series of records
+ std::vector<cell_record> cells_orig = {
+ cell_record(cell_record::soma, 0,
+ 14.566132, 34.873772, 7.857000, 0.717830, -1),
+ cell_record(cell_record::dendrite, 1,
+ 14.566132+1, 34.873772+1, 7.857000+1, 0.717830+1, -1)
+ };
+
+ std::stringstream swc_input;
+ swc_input << "# this is a comment\n";
+ swc_input << "# this is a comment\n";
+ for (auto c : cells_orig)
+ swc_input << c << "\n";
+
+ swc_input << "# this is a final comment\n";
+ try {
+ std::size_t nr_records = 0;
+ cell_record cell;
+ while ( !(swc_input >> cell).eof()) {
+ ASSERT_LT(nr_records, cells_orig.size());
+ expect_cell_equals(cells_orig[nr_records], cell);
+ ++nr_records;
+ }
+ } catch (std::exception &e) {
+ ADD_FAILURE();
+ }
+ }
+}
+
int main(int argc, char **argv) {
::testing::InitGoogleTest(&argc, argv);
return RUN_ALL_TESTS();
diff --git a/makefile b/makefile
index 3a75e5b6907eb4c1efd25d95d5ed8e521d911831..736a149f83bcb53f082a9a3656d7f1ef5ca5f321 100644
--- a/makefile
+++ b/makefile
@@ -1,5 +1,5 @@
CC=clang++
-FLAGS=-std=c++11 -g
+FLAGS=-std=c++11 -g -pedantic
test.exe : main.cpp *.hpp makefile gtest.o
$(CC) $(FLAGS) main.cpp -o test.exe gtest.o -pthread
diff --git a/swcio.hpp b/swcio.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..d362109dd15b771ac897bb884380682214f113b0
--- /dev/null
+++ b/swcio.hpp
@@ -0,0 +1,263 @@
+#pragma once
+
+#include <exception>
+#include <iostream>
+#include <sstream>
+#include <type_traits>
+
+namespace neuron
+{
+
+namespace io
+{
+
+
+static bool starts_with(const std::string &str, const std::string &prefix)
+{
+ return (str.find(prefix) == 0);
+}
+
+class cell_record
+{
+public:
+
+ // FIXME: enum's are not completely type-safe, since they can accept anything
+ // that can be casted to their underlying type.
+ //
+ // More on SWC files: http://research.mssm.edu/cnic/swc.html
+ enum kind {
+ undefined = 0,
+ soma,
+ axon,
+ dendrite,
+ apical_dendrite,
+ fork_point,
+ end_point,
+ custom
+ };
+
+ // cell records assume zero-based indexing; root's parent remains -1
+ cell_record(kind type, int id,
+ float x, float y, float z, float r,
+ int parent_id)
+ : type_(type)
+ , id_(id)
+ , x_(x)
+ , y_(y)
+ , z_(z)
+ , r_(r)
+ , parent_id_(parent_id)
+ {
+ // Check cell type as well; enum's do not offer complete type safety,
+ // since you can cast anything that fits to its underlying type
+ if (type_ < 0 || type_ > custom)
+ throw std::invalid_argument("unknown cell type");
+
+ if (id_ < 0)
+ throw std::invalid_argument("negative ids not allowed");
+
+ if (parent_id_ < -1)
+ throw std::invalid_argument("parent_id < -1 not allowed");
+
+ if (parent_id_ >= id_)
+ throw std::invalid_argument("parent_id >= id is not allowed");
+
+ if (r_ < 0)
+ throw std::invalid_argument("negative radii are not allowed");
+ }
+
+ cell_record()
+ : type_(cell_record::undefined)
+ , id_(0)
+ , x_(0)
+ , y_(0)
+ , z_(0)
+ , r_(0)
+ , parent_id_(-1)
+ { }
+
+ cell_record(const cell_record &other) = default;
+ cell_record &operator=(const cell_record &other) = default;
+
+ friend bool operator==(const cell_record &lhs,
+ const cell_record &rhs)
+ {
+ return lhs.id_ == rhs.id_;
+ }
+
+ friend bool operator!=(const cell_record &lhs,
+ const cell_record &rhs)
+ {
+ return !(lhs == rhs);
+ }
+
+ friend std::ostream &operator<<(std::ostream &os, const cell_record &cell);
+
+ kind type() const
+ {
+ return type_;
+ }
+
+ int id() const
+ {
+ return id_;
+ }
+
+ int parent() const
+ {
+ return parent_id_;
+ }
+
+ float x() const
+ {
+ return x_;
+ }
+
+ float y() const
+ {
+ return y_;
+ }
+
+ float z() const
+ {
+ return z_;
+ }
+
+ float radius() const
+ {
+ return r_;
+ }
+
+ float diameter() const
+ {
+ return 2*r_;
+ }
+
+private:
+ kind type_; // cell type
+ int id_; // cell id
+ float x_, y_, z_; // cell coordinates
+ float r_; // cell radius
+ int parent_id_; // cell parent's id
+};
+
+class swc_parser
+{
+public:
+ swc_parser(const std::string &delim,
+ std::string comment_prefix)
+ : delim_(delim)
+ , comment_prefix_(comment_prefix)
+ { }
+
+ swc_parser()
+ : delim_(" ")
+ , comment_prefix_("#")
+ { }
+
+ std::istream &parse_record(std::istream &is, cell_record &cell)
+ {
+ while (!is.eof() && !is.bad()) {
+ // consume empty and comment lines first
+ std::getline(is, linebuff_);
+ if (!linebuff_.empty() && !starts_with(linebuff_, comment_prefix_))
+ break;
+ }
+
+ if (is.bad())
+ // let the caller check for such events
+ return is;
+
+ if (is.eof() &&
+ (linebuff_.empty() || starts_with(linebuff_, comment_prefix_)))
+ // last line is either empty or a comment; don't parse anything
+ return is;
+
+ if (is.fail())
+ throw std::runtime_error("too long line detected");
+
+ std::istringstream line(linebuff_);
+ cell = parse_record(line);
+ return is;
+ }
+
+private:
+ void check_parse_status(const std::istream &is)
+ {
+ if (is.fail())
+ // If we try to read past the eof; fail bit will also be set
+ // FIXME: better throw a custom parse_error exception
+ throw std::logic_error("could not parse value");
+ }
+
+ template<typename T>
+ T parse_value_strict(std::istream &is)
+ {
+ T val;
+ check_parse_status(is >> val);
+
+ // everything's fine
+ return val;
+ }
+
+ // Read the record from a string stream; will be treated like a single line
+ cell_record parse_record(std::istringstream &is);
+
+ std::string delim_;
+ std::string comment_prefix_;
+ std::string linebuff_;
+};
+
+
+// specialize parsing for cell types
+template<>
+cell_record::kind swc_parser::parse_value_strict(std::istream &is)
+{
+ int val;
+ check_parse_status(is >> val);
+
+ // Let cell_record's constructor check for the type validity
+ return static_cast<cell_record::kind>(val);
+}
+
+cell_record swc_parser::parse_record(std::istringstream &is)
+{
+ auto id = parse_value_strict<int>(is);
+ auto type = parse_value_strict<cell_record::kind>(is);
+ auto x = parse_value_strict<float>(is);
+ auto y = parse_value_strict<float>(is);
+ auto z = parse_value_strict<float>(is);
+ auto r = parse_value_strict<float>(is);
+ auto parent_id = parse_value_strict<int>(is);
+
+ // Convert to zero-based, leaving parent_id as-is if -1
+ if (parent_id != -1)
+ parent_id--;
+
+ return cell_record(type, id-1, x, y, z, r, parent_id);
+}
+
+
+std::istream &operator>>(std::istream &is, cell_record &cell)
+{
+ swc_parser parser;
+ parser.parse_record(is, cell);
+ return is;
+}
+
+std::ostream &operator<<(std::ostream &os, const cell_record &cell)
+{
+ // output in one-based indexing
+ os << cell.id_+1 << " "
+ << cell.type_ << " "
+ << std::setprecision(7) << cell.x_ << " "
+ << std::setprecision(7) << cell.y_ << " "
+ << std::setprecision(7) << cell.z_ << " "
+ << std::setprecision(7) << cell.r_ << " "
+ << ((cell.parent_id_ == -1) ? cell.parent_id_ : cell.parent_id_+1);
+
+ return os;
+}
+
+} // end of neuron::io
+} // end of neuron