5679884058

pending completion

Build # 5679884058

Build Type

push

github

Committed by

mgeplf

Commit Message

cleanup

Run Details

8 of 8 new or added lines in 1 file covered. (100.0%)

1968 of 2587 relevant lines covered (76.07%)

905.27 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

89.07

/src/readers/morphologySWC.cpp

#include "morphologySWC.h"
#include "utils.h"

#include <cctype>         // isdigit
#include <cstdint>        // uint32_t
#include <memory>         // std::shared_ptr
#include <string>         // std::string
#include <unordered_map>  // std::unordered_map
#include <vector>         // std::vector

#include <morphio/errorMessages.h>
#include <morphio/mut/morphology.h>
#include <morphio/mut/section.h>
#include <morphio/mut/soma.h>
#include <morphio/properties.h>

namespace details {
using namespace morphio;

// It's not clear if -1 is the only way of identifying a root section.
const int SWC_UNDEFINED_PARENT = -1;
const unsigned int SWC_ROOT = 0xFFFFFFFD;

/* simple stream parser for SWC file format which is a line oriented format
 *
 * This parser advances across comments and blank lines, and allows the caller
 * to get integers and floats
 */
class SWCTokenizer
{
public:
  explicit SWCTokenizer(std::string contents, const morphio::readers::ErrorMessages& err)
      : contents_(std::move(contents))
      , err_(err) {
      // ensure null termination
      (void) contents_.c_str();
  }

  bool done() const noexcept {
      return pos_ >= contents_.size();
  }

  size_t lineNumber() const noexcept {
      return line_;
  }

  void skip_to(char value) {
      std::size_t pos = contents_.find_first_of(value, pos_);
      if (pos == std::string::npos) {
          pos_ = contents_.size();
      }
      pos_ = pos;
  }

  void advance_to_non_whitespace() {
      std::size_t pos = contents_.find_first_not_of(" \t", pos_);
      if (pos == std::string::npos) {
          pos_ = contents_.size();
      }
      pos_ = pos;
  }

  void advance_to_number() {
      while (consume_line_and_trailing_comments()) {
      }

      if (done()) {
          throw morphio::RawDataError(err_.EARLY_END_OF_FILE(line_));
      }

      auto c = contents_.at(pos_);
      if (std::isdigit(c) != 0 || c == '-' || c == '+' || c == '.') {
          return;
      }

      throw morphio::RawDataError(err_.ERROR_LINE_NON_PARSABLE(line_));
  }

  int64_t read_int() {
      advance_to_number();
      auto parsed = stn_.toInt(contents_, pos_);
      pos_ = std::get<1>(parsed);
      return std::get<0>(parsed);
  }

  morphio::floatType read_float() {
      advance_to_number();
      auto parsed = stn_.toFloat(contents_, pos_);
      pos_ = std::get<1>(parsed);
      return std::get<0>(parsed);
  }

  bool consume_line_and_trailing_comments() {
      bool found_newline = false;

      advance_to_non_whitespace();
      while (!done() && (contents_.at(pos_) == '#' || contents_.at(pos_) == '\n')) {
          switch (contents_.at(pos_)) {
          case '#':
              skip_to('\n');
              break;
          case '\n':
              ++line_;
              ++pos_;
              found_newline = true;
              break;
          }
          advance_to_non_whitespace();
      }
      return found_newline || done();
  }

private:
  size_t pos_ = 0;
  size_t line_ = 1;
  std::string contents_;
  morphio::StringToNumber stn_{};
  morphio::readers::ErrorMessages err_;
};

std::vector<morphio::readers::Sample> readSamples(const std::string& contents,
                                                  const morphio::readers::ErrorMessages& err) {
    std::vector<morphio::readers::Sample> samples;
    morphio::readers::Sample sample;

    SWCTokenizer tokenizer{contents, err};
    tokenizer.consume_line_and_trailing_comments();

    while (!tokenizer.done()) {
        sample.lineNumber = static_cast<unsigned int>(tokenizer.lineNumber());

        int64_t id = tokenizer.read_int();
        if (id < 0) {
            throw morphio::RawDataError(err.ERROR_NEGATIVE_ID(sample.lineNumber));
        }

        sample.id = static_cast<unsigned int>(id);

        sample.type = static_cast<morphio::SectionType>(tokenizer.read_int());

        for (auto& point : sample.point) {
            point = tokenizer.read_float();
        }

        sample.diameter = 2 * tokenizer.read_float();

        int64_t parentId = tokenizer.read_int();
        if (parentId < -1) {
            throw morphio::RawDataError(err.ERROR_NEGATIVE_ID(sample.lineNumber));
        } else if (parentId == SWC_UNDEFINED_PARENT) {
            sample.parentId = SWC_ROOT;
        } else {
            sample.parentId = static_cast<unsigned int>(parentId);
        }

        if (!tokenizer.consume_line_and_trailing_comments()) {
            throw morphio::RawDataError(err.ERROR_LINE_NON_PARSABLE(sample.lineNumber));
        }
        samples.push_back(sample);
    }
    return samples;
}

/**
  Parsing SWC according to this specification:
http://www.neuronland.org/NLMorphologyConverter/MorphologyFormats/SWC/Spec.html
 **/
using morphio::readers::ErrorMessages;
using morphio::readers::Sample;

enum class DeclaredID : unsigned int {};
unsigned int operator+(DeclaredID type) {
    return static_cast<unsigned int>(type);
}

class SWCBuilder
{
    using MorphioID = uint32_t;
    using PhysicalID = size_t;
    using Samples = std::vector<Sample>;

  public:
    explicit SWCBuilder(const std::string& path)
        : err_(path) {}

    Property::Properties buildProperties(const std::string& contents, unsigned int options) {
        const Samples samples = readSamples(contents, err_);
        buildSWC(samples);
        morph_.applyModifiers(options);
        return morph_.buildReadOnly();
    }

  private:
    void _checkNeuromorph3PointSoma(const Samples& soma_samples){
        // First point is the 'center'; has 2 children
        const Sample& center = soma_samples[0];
        const Sample& child1 = soma_samples[1];
        const Sample& child2 = soma_samples[2];

        floatType x = center.point[0];
        floatType z = center.point[2];
        floatType d = center.diameter;
        floatType r = center.diameter / 2;
        floatType y = center.point[1];

        // whether the soma should be checked for the special case of 3 point soma
        // for details see https://github.com/BlueBrain/MorphIO/issues/273
        // If the 2nd and the 3rd point have the same x, z, d values then the only valid soma
        // is: 1 1 x   y   z r -1 2 1 x (y-r) z r  1 3 1 x (y+r) z r  1
        if ((child1.point[0] != x || child2.point[0] != x || child1.point[1] != y - r ||
             child2.point[1] != y + r || child1.point[2] != z || child2.point[2] != z ||
             child1.diameter != d || child2.diameter != d) &&
            std::fabs(child1.diameter - d) < morphio::epsilon &&
            std::fabs(child2.diameter - d) < morphio::epsilon &&
            std::fabs(child1.point[0] - x) < morphio::epsilon &&
            std::fabs(child2.point[0] - x) < morphio::epsilon &&
            std::fabs(child1.point[2] - z) < morphio::epsilon &&
            std::fabs(child2.point[2] - z) < morphio::epsilon
            ) {
            printError(Warning::SOMA_NON_CONFORM,
                       err_.WARNING_NEUROMORPHO_SOMA_NON_CONFORM(center, child1, child2));
        }
    }

    void build_soma(const Samples& soma_samples) {
        auto& soma = morph_.soma();

        if (soma_samples.empty()) {
            soma->type() = SOMA_UNDEFINED;
            return;
        } else if (soma_samples.size() == 1) {
            Sample sample = soma_samples[0];

            if (sample.parentId != SWC_ROOT && samples_.at(sample.parentId).type != SECTION_SOMA) {
                throw morphio::SomaError(err_.ERROR_SOMA_WITH_NEURITE_PARENT(sample));
            }

            soma->type() = SOMA_SINGLE_POINT;
            soma->points() = {sample.point};
            soma->diameters() = {sample.diameter};
            return;
        } else if (soma_samples.size() == 3) {
            const Sample& center = soma_samples[0];
            const Sample& child1 = soma_samples[1];
            const Sample& child2 = soma_samples[2];
            // All soma that bifurcate with the first parent having two children are considered
            // SOMA_NEUROMORPHO_THREE_POINT_CYLINDERS
            if(center.id == child1.parentId && center.id == child2.parentId){
                soma->type() = SOMA_NEUROMORPHO_THREE_POINT_CYLINDERS;
                soma->points() = {center.point, child1.point, child2.point};
                soma->diameters() = {center.diameter, child1.diameter, child2.diameter};
                _checkNeuromorph3PointSoma(soma_samples);
               return;
           }
        }
        // might also have 3 points at this point, as well

        // a "normal" SWC soma
        soma->type() = SOMA_CYLINDERS;
        auto& points = soma->points();
        auto& diameters = soma->diameters();
        points.reserve(soma_samples.size());
        diameters.reserve(soma_samples.size());

        size_t parent_count = 0;
        for (const auto& s : soma_samples) {
            if (s.parentId == SWC_ROOT) {
                parent_count++;
            } else if (samples_.count(s.parentId) == 0) {
                throw morphio::MissingParentError(err_.ERROR_MISSING_PARENT(s));
            } else if (samples_.at(s.parentId).type != SECTION_SOMA) {
                throw morphio::SomaError(err_.ERROR_SOMA_WITH_NEURITE_PARENT(s));
            }

            if (children_.count(s.id) > 0 && children_.at(s.id).size() > 1) {
                std::vector<Sample> soma_bifurcations;
                for (auto id : children_.at(s.id)) {
                    if (samples_[id].type == SECTION_SOMA) {
                        soma_bifurcations.push_back(samples_[id]);
                    }
                }
                if (soma_bifurcations.size() > 1) {
                    throw morphio::SomaError(err_.ERROR_SOMA_BIFURCATION(s, soma_bifurcations));
                }
            }
            points.push_back(s.point);
            diameters.push_back(s.diameter);
        }

        if (parent_count > 1) {
            throw morphio::SomaError(err_.ERROR_MULTIPLE_SOMATA(soma_samples));
        }
    }

    void buildSWC(const Samples& samples) {
        Samples soma_samples;
        Samples root_samples;

        for (const auto& sample: samples) {
            // { checks
            if (sample.diameter < morphio::epsilon) {
                printError(Warning::ZERO_DIAMETER, err_.WARNING_ZERO_DIAMETER(sample));
            }

            if (sample.parentId == sample.id) {
                throw morphio::RawDataError(err_.ERROR_SELF_PARENT(sample));
            }

            if (sample.type >= morphio::SECTION_OUT_OF_RANGE_START || sample.type <= 0) {
                throw morphio::RawDataError(
                    err_.ERROR_UNSUPPORTED_SECTION_TYPE(sample.lineNumber, sample.type));
            }
            if (sample.parentId == SWC_ROOT && sample.type != SECTION_SOMA) {
                printError(Warning::DISCONNECTED_NEURITE,
                           err_.WARNING_DISCONNECTED_NEURITE(sample));
            }
            // } checks

            if (sample.type == SECTION_SOMA) {
                soma_samples.push_back(sample);
            }

            if (sample.parentId == SWC_ROOT || sample.type == SECTION_SOMA) {
                root_samples.push_back(sample);
            }

            if (!samples_.insert({sample.id, sample}).second) {
                throw RawDataError(err_.ERROR_REPEATED_ID(samples[sample.id], sample));
            }

            children_[sample.parentId].push_back(sample.id);
        }

        // can only check for missing parents once all samples are loaded
        // since it's possible there may be forward references
        for (const auto& sample: samples) {
            if(sample.parentId != SWC_ROOT && samples_.count(sample.parentId) == 0){
                throw morphio::MissingParentError(err_.ERROR_MISSING_PARENT(sample));
            }
        }

        build_soma(soma_samples);

        std::unordered_map<DeclaredID, std::shared_ptr<morphio::mut::Section>> declared_to_swc;
        declared_to_swc.reserve(samples.size());

        for (const Sample& root_sample : root_samples) {
            if (children_.count(root_sample.id) == 0) {
                continue;
            }
            for (unsigned int child_id : children_.at(root_sample.id)) {
                if (samples_.at(child_id).type == SECTION_SOMA) {
                    continue;
                }
                if (root_sample.type == SECTION_SOMA) {
                    assembleSections(child_id,
                                     DeclaredID(root_sample.id),
                                     declared_to_swc,
                                     morph_.soma()->points()[0],
                                     morph_.soma()->diameters()[0],
                                     true);
                } else {
                    // this is neurite as the start
                    assembleSections(root_sample.id,
                                     DeclaredID(SWC_ROOT),
                                     declared_to_swc,
                                     root_sample.point,
                                     root_sample.diameter,
                                     true);
                    break;
                }
            }
        }
    }

    void assembleSections(
        unsigned int id,
        DeclaredID parent_id,
        std::unordered_map<DeclaredID, std::shared_ptr<morphio::mut::Section>>& declared_to_swc,
        const Point& start_point,
        floatType start_diameter,
        bool is_root) {

        morphio::Property::PointLevel properties;
        auto& points = properties._points;
        auto& diameters = properties._diameters;

        auto appendSection = [&](DeclaredID section_id_, DeclaredID parent_id_, SectionType starting_section_type) {
            std::shared_ptr<morphio::mut::Section> new_section;
            if (is_root) {
                new_section = morph_.appendRootSection(properties, starting_section_type);
            } else {
                new_section = declared_to_swc.at(parent_id_)
                                  ->appendSection(properties, starting_section_type);
            }
            declared_to_swc[section_id_] = new_section;
        };

        auto get_child_count = [&](unsigned int child_id) {
            return children_.count(child_id) == 0 ? 0 : children_.at(child_id).size();
        };

        const Sample* sample = &samples_.at(id);

        // create duplicate point if needed
        if (!is_root && sample->point != start_point /*|| sample->diameter != start_diameter */) {
        //if (!(is_root && sample->type == SECTION_SOMA) && sample->point != start_point /*|| sample->diameter != start_diameter */) {
            points.push_back(start_point);
            diameters.push_back(start_diameter);
        }

        // try and combine as many single samples into a single section as possible
        size_t children_count = get_child_count(id);
        while (children_count == 1) {
            sample = &samples_.at(id);
            if(sample->type != samples_.at(children_.at(id)[0]).type){
                break;
            }
            points.push_back(sample->point);
            diameters.push_back(sample->diameter);
            id = children_.at(id)[0];
            children_count = get_child_count(id);
        }
        sample = &samples_.at(id);
        points.push_back(sample->point);
        diameters.push_back(sample->diameter);
        appendSection(DeclaredID(id), parent_id, sample->type);

        if (children_count == 0) {
            // section was already appended above, nothing to do
        } else if (children_count == 1) {
            // section_type changed
            size_t offset = properties._points.size() - 1;
            const Point& new_start_point = properties._points[offset];
            floatType new_start_diameter = properties._diameters[offset];
            assembleSections(children_.at(id)[0], DeclaredID(id), declared_to_swc, new_start_point, new_start_diameter, false);
        } else {
            size_t offset = properties._points.size() - 1;
            const Point& new_start_point = properties._points[offset];
            floatType new_start_diameter = properties._diameters[offset];
            for (unsigned int child_id : children_.at(id)) {
                assembleSections(child_id,
                                 DeclaredID(id),
                                 declared_to_swc,
                                 new_start_point,
                                 new_start_diameter,
                                 false);
            }
        }
    }

    std::unordered_map<unsigned int, std::vector<unsigned int>> children_;
    std::unordered_map<unsigned int, Sample> samples_;
    mut::Morphology morph_;
    ErrorMessages err_;
};


}  // namespace details

namespace morphio {
namespace readers {
namespace swc {
Property::Properties load(const std::string& path,
                          const std::string& contents,
                          unsigned int options) {
    auto properties = details::SWCBuilder(path).buildProperties(contents, options);

    properties._cellLevel._cellFamily = NEURON;
    properties._cellLevel._version = {"swc", 1, 0};
    return properties;
}

}  // namespace swc
}  // namespace readers
}  // namespace morphio

1	#include "morphologySWC.h"
2	#include "utils.h"
3
4	#include <cctype> // isdigit
5	#include <cstdint> // uint32_t
6	#include <memory> // std::shared_ptr
7	#include <string> // std::string
8	#include <unordered_map> // std::unordered_map
9	#include <vector> // std::vector
10
11	#include <morphio/errorMessages.h>
12	#include <morphio/mut/morphology.h>
13	#include <morphio/mut/section.h>
14	#include <morphio/mut/soma.h>
15	#include <morphio/properties.h>
16
17	namespace details {
18	using namespace morphio;
19
20	// It's not clear if -1 is the only way of identifying a root section.
21	const int SWC_UNDEFINED_PARENT = -1;
22	const unsigned int SWC_ROOT = 0xFFFFFFFD;
23
24	/* simple stream parser for SWC file format which is a line oriented format
25	*
26	* This parser advances across comments and blank lines, and allows the caller
27	* to get integers and floats
28	*/
29	class SWCTokenizer
30	{
31	public:
32	explicit SWCTokenizer(std::string contents, const morphio::readers::ErrorMessages& err)	66✔
33	: contents_(std::move(contents))	66✔
34	, err_(err) {	66✔
35	// ensure null termination
36	(void) contents_.c_str();	66✔
37	}	66✔
38
39	bool done() const noexcept {	13,728✔
40	return pos_ >= contents_.size();	13,728✔
41	}
42
43	size_t lineNumber() const noexcept {	506✔
44	return line_;	506✔
45	}
46
47	void skip_to(char value) {	684✔
48	std::size_t pos = contents_.find_first_of(value, pos_);	684✔
49	if (pos == std::string::npos) {	684✔
50	pos_ = contents_.size();	×
51	}
52	pos_ = pos;	684✔
53	}	684✔
54
55	void advance_to_non_whitespace() {	6,072✔
56	std::size_t pos = contents_.find_first_not_of(" \t", pos_);	6,072✔
57	if (pos == std::string::npos) {	6,072✔
58	pos_ = contents_.size();	64✔
59	}
60	pos_ = pos;	6,072✔
61	}	6,072✔
62
63	void advance_to_number() {	3,542✔
64	while (consume_line_and_trailing_comments()) {	3,542✔
65	}
66
67	if (done()) {	3,542✔
68	throw morphio::RawDataError(err_.EARLY_END_OF_FILE(line_));	×
69	}
70
71	auto c = contents_.at(pos_);	3,542✔
72	if (std::isdigit(c) != 0 \|\| c == '-' \|\| c == '+' \|\| c == '.') {	3,542✔
73	return;	3,540✔
74	}
75
76	throw morphio::RawDataError(err_.ERROR_LINE_NON_PARSABLE(line_));	2✔
77	}
78
79	int64_t read_int() {	1,518✔
80	advance_to_number();	1,518✔
81	auto parsed = stn_.toInt(contents_, pos_);	1,516✔
82	pos_ = std::get<1>(parsed);	1,516✔
83	return std::get<0>(parsed);	1,516✔
84	}
85
86	morphio::floatType read_float() {	2,024✔
87	advance_to_number();	2,024✔
88	auto parsed = stn_.toFloat(contents_, pos_);	2,024✔
89	pos_ = std::get<1>(parsed);	2,024✔
90	return std::get<0>(parsed);	2,024✔
91	}
92
93	bool consume_line_and_trailing_comments() {	4,112✔
94	bool found_newline = false;	4,112✔
95
96	advance_to_non_whitespace();	4,112✔
97	while (!done() && (contents_.at(pos_) == '#' \|\| contents_.at(pos_) == '\n')) {	6,072✔
98	switch (contents_.at(pos_)) {	1,960✔
99	case '#':	684✔
100	skip_to('\n');	684✔
101	break;	684✔
102	case '\n':	1,276✔
103	++line_;	1,276✔
104	++pos_;	1,276✔
105	found_newline = true;	1,276✔
106	break;	1,276✔
107	}
108	advance_to_non_whitespace();	1,960✔
109	}
110	return found_newline \|\| done();	4,112✔
111	}
112
113	private:
114	size_t pos_ = 0;
115	size_t line_ = 1;
116	std::string contents_;
117	morphio::StringToNumber stn_{};
118	morphio::readers::ErrorMessages err_;
119	};
120
121	std::vector<morphio::readers::Sample> readSamples(const std::string& contents,	66✔
122	const morphio::readers::ErrorMessages& err) {
123	std::vector<morphio::readers::Sample> samples;	66✔
124	morphio::readers::Sample sample;	66✔
125
126	SWCTokenizer tokenizer{contents, err};	132✔
127	tokenizer.consume_line_and_trailing_comments();	66✔
128
129	while (!tokenizer.done()) {	570✔
130	sample.lineNumber = static_cast<unsigned int>(tokenizer.lineNumber());	506✔
131
132	int64_t id = tokenizer.read_int();	506✔
133	if (id < 0) {	506✔
134	throw morphio::RawDataError(err.ERROR_NEGATIVE_ID(sample.lineNumber));	×
135	}
136
137	sample.id = static_cast<unsigned int>(id);	506✔
138
139	sample.type = static_cast<morphio::SectionType>(tokenizer.read_int());	506✔
140
141	for (auto& point : sample.point) {	2,024✔
142	point = tokenizer.read_float();	1,518✔
143	}
144
145	sample.diameter = 2 * tokenizer.read_float();	506✔
146
147	int64_t parentId = tokenizer.read_int();	506✔
148	if (parentId < -1) {	504✔
149	throw morphio::RawDataError(err.ERROR_NEGATIVE_ID(sample.lineNumber));	×
150	} else if (parentId == SWC_UNDEFINED_PARENT) {	504✔
151	sample.parentId = SWC_ROOT;	66✔
152	} else {
153	sample.parentId = static_cast<unsigned int>(parentId);	438✔
154	}
155
156	if (!tokenizer.consume_line_and_trailing_comments()) {	504✔
157	throw morphio::RawDataError(err.ERROR_LINE_NON_PARSABLE(sample.lineNumber));	×
158	}
159	samples.push_back(sample);	504✔
160	}
161	return samples;	128✔
162	}
163
164	/**
165	Parsing SWC according to this specification:
166	http://www.neuronland.org/NLMorphologyConverter/MorphologyFormats/SWC/Spec.html
167	**/
168	using morphio::readers::ErrorMessages;
169	using morphio::readers::Sample;
170
171	enum class DeclaredID : unsigned int {};
172	unsigned int operator+(DeclaredID type) {	×
173	return static_cast<unsigned int>(type);	×
174	}
175
176	class SWCBuilder
177	{
178	using MorphioID = uint32_t;
179	using PhysicalID = size_t;
180	using Samples = std::vector<Sample>;
181
182	public:
183	explicit SWCBuilder(const std::string& path)	66✔
184	: err_(path) {}	66✔
185
186	Property::Properties buildProperties(const std::string& contents, unsigned int options) {	66✔
187	const Samples samples = readSamples(contents, err_);	130✔
188	buildSWC(samples);	64✔
189	morph_.applyModifiers(options);	50✔
190	return morph_.buildReadOnly();	100✔
191	}
192
193	private:
194	void _checkNeuromorph3PointSoma(const Samples& soma_samples){	4✔
195	// First point is the 'center'; has 2 children
196	const Sample& center = soma_samples[0];	4✔
197	const Sample& child1 = soma_samples[1];	4✔
198	const Sample& child2 = soma_samples[2];	4✔
199
200	floatType x = center.point[0];	4✔
201	floatType z = center.point[2];	4✔
202	floatType d = center.diameter;	4✔
203	floatType r = center.diameter / 2;	4✔
204	floatType y = center.point[1];	4✔
205
206	// whether the soma should be checked for the special case of 3 point soma
207	// for details see https://github.com/BlueBrain/MorphIO/issues/273
208	// If the 2nd and the 3rd point have the same x, z, d values then the only valid soma
209	// is: 1 1 x y z r -1 2 1 x (y-r) z r 1 3 1 x (y+r) z r 1
210	if ((child1.point[0] != x \|\| child2.point[0] != x \|\| child1.point[1] != y - r \|\|	8✔
211	child2.point[1] != y + r \|\| child1.point[2] != z \|\| child2.point[2] != z \|\|	4✔
212	child1.diameter != d \|\| child2.diameter != d) &&	4✔
213	std::fabs(child1.diameter - d) < morphio::epsilon &&	×
214	std::fabs(child2.diameter - d) < morphio::epsilon &&	×
215	std::fabs(child1.point[0] - x) < morphio::epsilon &&	×
216	std::fabs(child2.point[0] - x) < morphio::epsilon &&	×
217	std::fabs(child1.point[2] - z) < morphio::epsilon &&	8✔
218	std::fabs(child2.point[2] - z) < morphio::epsilon	×
219	) {
220	printError(Warning::SOMA_NON_CONFORM,	×
221	err_.WARNING_NEUROMORPHO_SOMA_NON_CONFORM(center, child1, child2));	×
222	}
223	}	4✔
224
225	void build_soma(const Samples& soma_samples) {	56✔
226	auto& soma = morph_.soma();	56✔
227
228	if (soma_samples.empty()) {	56✔
229	soma->type() = SOMA_UNDEFINED;	×
230	return;	×
231	} else if (soma_samples.size() == 1) {	56✔
232	Sample sample = soma_samples[0];	38✔
233
234	if (sample.parentId != SWC_ROOT && samples_.at(sample.parentId).type != SECTION_SOMA) {	38✔
235	throw morphio::SomaError(err_.ERROR_SOMA_WITH_NEURITE_PARENT(sample));	2✔
236	}
237
238	soma->type() = SOMA_SINGLE_POINT;	36✔
239	soma->points() = {sample.point};	36✔
240	soma->diameters() = {sample.diameter};	36✔
241	return;	36✔
242	} else if (soma_samples.size() == 3) {	18✔
243	const Sample& center = soma_samples[0];	4✔
244	const Sample& child1 = soma_samples[1];	4✔
245	const Sample& child2 = soma_samples[2];	4✔
246	// All soma that bifurcate with the first parent having two children are considered
247	// SOMA_NEUROMORPHO_THREE_POINT_CYLINDERS
248	if(center.id == child1.parentId && center.id == child2.parentId){	4✔
249	soma->type() = SOMA_NEUROMORPHO_THREE_POINT_CYLINDERS;	4✔
250	soma->points() = {center.point, child1.point, child2.point};	4✔
251	soma->diameters() = {center.diameter, child1.diameter, child2.diameter};	4✔
252	_checkNeuromorph3PointSoma(soma_samples);	4✔
253	return;	4✔
254	}
255	}
256	// might also have 3 points at this point, as well
257
258	// a "normal" SWC soma
259	soma->type() = SOMA_CYLINDERS;	14✔
260	auto& points = soma->points();	14✔
261	auto& diameters = soma->diameters();	14✔
262	points.reserve(soma_samples.size());	14✔
263	diameters.reserve(soma_samples.size());	14✔
264
265	size_t parent_count = 0;	14✔
266	for (const auto& s : soma_samples) {	56✔
267	if (s.parentId == SWC_ROOT) {	44✔
268	parent_count++;	16✔
269	} else if (samples_.count(s.parentId) == 0) {	28✔
270	throw morphio::MissingParentError(err_.ERROR_MISSING_PARENT(s));	×
271	} else if (samples_.at(s.parentId).type != SECTION_SOMA) {	28✔
272	throw morphio::SomaError(err_.ERROR_SOMA_WITH_NEURITE_PARENT(s));	×
273	}
274
275	if (children_.count(s.id) > 0 && children_.at(s.id).size() > 1) {	44✔
276	std::vector<Sample> soma_bifurcations;	20✔
277	for (auto id : children_.at(s.id)) {	30✔
278	if (samples_[id].type == SECTION_SOMA) {	20✔
279	soma_bifurcations.push_back(samples_[id]);	12✔
280	}
281	}
282	if (soma_bifurcations.size() > 1) {	10✔
283	throw morphio::SomaError(err_.ERROR_SOMA_BIFURCATION(s, soma_bifurcations));	2✔
284	}
285	}
286	points.push_back(s.point);	42✔
287	diameters.push_back(s.diameter);	42✔
288	}
289
290	if (parent_count > 1) {	12✔
291	throw morphio::SomaError(err_.ERROR_MULTIPLE_SOMATA(soma_samples));	2✔
292	}
293	}
294
295	void buildSWC(const Samples& samples) {	64✔
296	Samples soma_samples;	128✔
297	Samples root_samples;	128✔
298
299	for (const auto& sample: samples) {	562✔
300	// { checks
301	if (sample.diameter < morphio::epsilon) {	504✔
302	printError(Warning::ZERO_DIAMETER, err_.WARNING_ZERO_DIAMETER(sample));	10✔
303	}
304
305	if (sample.parentId == sample.id) {	504✔
306	throw morphio::RawDataError(err_.ERROR_SELF_PARENT(sample));	2✔
307	}
308
309	if (sample.type >= morphio::SECTION_OUT_OF_RANGE_START \|\| sample.type <= 0) {	502✔
310	throw morphio::RawDataError(
311	err_.ERROR_UNSUPPORTED_SECTION_TYPE(sample.lineNumber, sample.type));	2✔
312	}
313	if (sample.parentId == SWC_ROOT && sample.type != SECTION_SOMA) {	500✔
314	printError(Warning::DISCONNECTED_NEURITE,	4✔
315	err_.WARNING_DISCONNECTED_NEURITE(sample));	8✔
316	}
317	// } checks
318
319	if (sample.type == SECTION_SOMA) {	500✔
320	soma_samples.push_back(sample);	104✔
321	}
322
323	if (sample.parentId == SWC_ROOT \|\| sample.type == SECTION_SOMA) {	500✔
324	root_samples.push_back(sample);	108✔
325	}
326
327	if (!samples_.insert({sample.id, sample}).second) {	500✔
328	throw RawDataError(err_.ERROR_REPEATED_ID(samples[sample.id], sample));	2✔
329	}
330
331	children_[sample.parentId].push_back(sample.id);	498✔
332	}
333
334	// can only check for missing parents once all samples are loaded
335	// since it's possible there may be forward references
336	for (const auto& sample: samples) {	548✔
337	if(sample.parentId != SWC_ROOT && samples_.count(sample.parentId) == 0){	492✔
338	throw morphio::MissingParentError(err_.ERROR_MISSING_PARENT(sample));	2✔
339	}
340	}
341
342	build_soma(soma_samples);	56✔
343
344	std::unordered_map<DeclaredID, std::shared_ptr<morphio::mut::Section>> declared_to_swc;	100✔
345	declared_to_swc.reserve(samples.size());	50✔
346
347	for (const Sample& root_sample : root_samples) {	134✔
348	if (children_.count(root_sample.id) == 0) {	84✔
349	continue;	18✔
350	}
351	for (unsigned int child_id : children_.at(root_sample.id)) {	218✔
352	if (samples_.at(child_id).type == SECTION_SOMA) {	152✔
353	continue;	34✔
354	}
355	if (root_sample.type == SECTION_SOMA) {	118✔
356	assembleSections(child_id,	236✔
357	DeclaredID(root_sample.id),	118✔
358	declared_to_swc,
359	morph_.soma()->points()[0],	118✔
360	morph_.soma()->diameters()[0],	118✔
361	true);
362	} else {
363	// this is neurite as the start
364	assembleSections(root_sample.id,	×
365	DeclaredID(SWC_ROOT),
366	declared_to_swc,
367	root_sample.point,	×
368	root_sample.diameter,	×
369	true);
370	break;	×
371	}
372	}
373	}
374	}	50✔
375
376	void assembleSections(	238✔
377	unsigned int id,
378	DeclaredID parent_id,
379	std::unordered_map<DeclaredID, std::shared_ptr<morphio::mut::Section>>& declared_to_swc,
380	const Point& start_point,
381	floatType start_diameter,
382	bool is_root) {
383
384	morphio::Property::PointLevel properties;	476✔
385	auto& points = properties._points;	238✔
386	auto& diameters = properties._diameters;	238✔
387
388	auto appendSection = [&](DeclaredID section_id_, DeclaredID parent_id_, SectionType starting_section_type) {	238✔
389	std::shared_ptr<morphio::mut::Section> new_section;	238✔
390	if (is_root) {	238✔
391	new_section = morph_.appendRootSection(properties, starting_section_type);	118✔
392	} else {
393	new_section = declared_to_swc.at(parent_id_)	120✔
394	->appendSection(properties, starting_section_type);	120✔
395	}
396	declared_to_swc[section_id_] = new_section;	238✔
397	};	238✔
398
399	auto get_child_count = [&](unsigned int child_id) {	388✔
400	return children_.count(child_id) == 0 ? 0 : children_.at(child_id).size();	388✔
401	};	238✔
402
403	const Sample* sample = &samples_.at(id);	238✔
404
405	// create duplicate point if needed
406	if (!is_root && sample->point != start_point /\|\| sample->diameter != start_diameter /) {	238✔
407	//if (!(is_root && sample->type == SECTION_SOMA) && sample->point != start_point /\|\| sample->diameter != start_diameter /) {
408	points.push_back(start_point);	112✔
409	diameters.push_back(start_diameter);	112✔
410	}
411
412	// try and combine as many single samples into a single section as possible
413	size_t children_count = get_child_count(id);	238✔
414	while (children_count == 1) {	388✔
415	sample = &samples_.at(id);	150✔
416	if(sample->type != samples_.at(children_.at(id)[0]).type){	150✔
417	break;	×
418	}
419	points.push_back(sample->point);	150✔
420	diameters.push_back(sample->diameter);	150✔
421	id = children_.at(id)[0];	150✔
422	children_count = get_child_count(id);	150✔
423	}
424	sample = &samples_.at(id);	238✔
425	points.push_back(sample->point);	238✔
426	diameters.push_back(sample->diameter);	238✔
427	appendSection(DeclaredID(id), parent_id, sample->type);	238✔
428
429	if (children_count == 0) {	238✔
430	// section was already appended above, nothing to do
431	} else if (children_count == 1) {	60✔
432	// section_type changed
433	size_t offset = properties._points.size() - 1;	×
434	const Point& new_start_point = properties._points[offset];	×
435	floatType new_start_diameter = properties._diameters[offset];	×
436	assembleSections(children_.at(id)[0], DeclaredID(id), declared_to_swc, new_start_point, new_start_diameter, false);	×
437	} else {
438	size_t offset = properties._points.size() - 1;	60✔
439	const Point& new_start_point = properties._points[offset];	60✔
440	floatType new_start_diameter = properties._diameters[offset];	60✔
441	for (unsigned int child_id : children_.at(id)) {	180✔
442	assembleSections(child_id,	120✔
443	DeclaredID(id),
444	declared_to_swc,
445	new_start_point,
446	new_start_diameter,
447	false);
448	}
449	}
450	}	238✔
451
452	std::unordered_map<unsigned int, std::vector<unsigned int>> children_;
453	std::unordered_map<unsigned int, Sample> samples_;
454	mut::Morphology morph_;
455	ErrorMessages err_;
456	};
457
458
459	} // namespace details
460
461	namespace morphio {
462	namespace readers {
463	namespace swc {
464	Property::Properties load(const std::string& path,	66✔
465	const std::string& contents,
466	unsigned int options) {
467	auto properties = details::SWCBuilder(path).buildProperties(contents, options);	82✔
468
469	properties._cellLevel._cellFamily = NEURON;	50✔
470	properties._cellLevel._version = {"swc", 1, 0};	50✔
471	return properties;	50✔
472	}
473
474	} // namespace swc
475	} // namespace readers
476	} // namespace morphio

BlueBrain / MorphIO / 5679884058

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous