6533524408

Committed 16 Oct 2023 12:24PM UTC coverage: 76.051% (-0.02%) from 76.073%

Build # 6533524408

Build Type

push

github

Committed by

mgeplf

Commit Message

fix test_root_node_split

Run Details

1972 of 2593 relevant lines covered (76.05%)

903.23 hits per line

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

96.65

/src/readers/morphologyHDF5.cpp

/* Copyright (c) 2013-2023, EPFL/Blue Brain Project
 *
 * SPDX-License-Identifier: Apache-2.0
 */

#include <cassert>

#include "morphologyHDF5.h"

#include <highfive/H5Utility.hpp>  // HighFive::SilenceHDF5

#include <morphio/errorMessages.h>

namespace {

constexpr size_t SECTION_START_OFFSET = 0;
constexpr size_t SECTION_TYPE = 1;
constexpr size_t SECTION_PARENT_OFFSET = 2;

constexpr int SOMA_ONLY = -1;

//{v1
const std::string _d_structure("structure");
const std::string _d_points("points");

//{ v1.1
const std::string _a_version("version");
const std::string _g_metadata("metadata");
const std::string _a_family("cell_family");
const std::string _d_perimeters("perimeters");
//} v1.1

//{ v1.2
const std::string _g_mitochondria("organelles/mitochondria");

// endoplasmic reticulum
const std::string _g_endoplasmic_reticulum("organelles/endoplasmic_reticulum");
const std::string _d_section_index("section_index");
const std::string _d_volume("volume");
const std::string _d_surface_area("surface_area");
const std::string _d_filament_count("filament_count");
// } v1.2

//{ v1.3

// Dendritic Spine
const std::string _g_postsynaptic_density("organelles/postsynaptic_density");
const std::string _d_dendritic_spine_section_id("section_id");
const std::string _d_dendritic_spine_segment_id("segment_id");
const std::string _d_dendritic_spine_offset("offset");

//}

// } v1

//{ v2
const std::string _g_v2root("neuron1");
//} v2

}  // namespace

namespace morphio {
namespace readers {
namespace h5 {

MorphologyHDF5::MorphologyHDF5(const HighFive::Group& group)
    : _group(group)
    , _uri("HDF5 Group") {}

Property::Properties load(const std::string& uri) {
    try {
        std::lock_guard<std::recursive_mutex> lock(morphio::readers::h5::global_hdf5_mutex());
        HighFive::SilenceHDF5 silence;
        auto file = HighFive::File(uri, HighFive::File::ReadOnly);
        return MorphologyHDF5(file.getGroup("/")).load();

    } catch (const HighFive::FileException& exc) {
        throw RawDataError("Could not open morphology file " + uri + ": " + exc.what());
    }
}

Property::Properties load(const HighFive::Group& group) {
    std::lock_guard<std::recursive_mutex> lock(morphio::readers::h5::global_hdf5_mutex());
    return MorphologyHDF5(group).load();
}

Property::Properties MorphologyHDF5::load() {
    _readMetadata(_uri);

    int firstSectionOffset = _readSections();

    _readPoints(firstSectionOffset);

    if (_properties._cellLevel.minorVersion() >= 1) {
        _readPerimeters(firstSectionOffset);

        if (_properties._cellLevel.minorVersion() >= 2) {
            _readMitochondria();
            _readEndoplasmicReticulum();
        }

        if (_properties._cellLevel.minorVersion() >= 3 &&
            _properties._cellLevel._cellFamily == CellFamily::SPINE) {
            _readDendriticSpinePostSynapticDensity();
        }
    }

    return _properties;
}

void MorphologyHDF5::_readMetadata(const std::string& source) {
    // default to h5v1.0
    uint32_t majorVersion = 1;
    uint32_t minorVersion = 0;
    _properties._cellLevel._cellFamily = CellFamily::NEURON;

    if (!_group.exist(_d_points) || !_group.exist(_d_structure)) {
        // h5v2 is deprecated, but it can be detected, throw a custom error messages if it is
        if (_group.exist(_g_v2root)) {
            throw RawDataError(
                "Error in " + source +
                "\nh5v2 is no longer supported, see: https://github.com/BlueBrain/MorphIO#H5v2");
        }
        throw RawDataError("Missing " + _d_points + " or " + _d_structure +
                           " datasets, cannot load morphology without them");
    }

    // if there is metadata, perhaps it's h5v1, minor version 1, 2
    if (_group.exist(_g_metadata)) {
        const auto metadata = _group.getGroup(_g_metadata);
        if (metadata.hasAttribute(_a_version)) {
            const auto attr = metadata.getAttribute(_a_version);

            std::array<uint32_t, 2> versions = {0, 0};
            attr.read(versions);
            majorVersion = versions[0];
            minorVersion = versions[1];

            if (majorVersion == 1 &&
                (minorVersion == 1 || minorVersion == 2 || minorVersion == 3)) {
                uint32_t family;
                metadata.getAttribute(_a_family).read(family);
                _properties._cellLevel._cellFamily = static_cast<CellFamily>(family);
            } else {
                std::string msg = "Error in " + source +
                                  "\nUnsupported h5 version: " + std::to_string(majorVersion) +
                                  "." + std::to_string(minorVersion);
                throw RawDataError(msg);
            }
        } else {
            throw RawDataError("Missing " + _a_version +
                               " attribute, cannot load morphology without them");
        }
    }

    _properties._cellLevel._version = {"h5", majorVersion, minorVersion};
}

void MorphologyHDF5::_readPoints(int firstSectionOffset) {
    constexpr size_t pointColumns = 4;

    const auto pointsDataSet = _group.getDataSet(_d_points);
    const auto pointsDims = pointsDataSet.getSpace().getDimensions();
    const size_t numberPoints = pointsDims[0];

    if (pointsDims.size() != 2) {
        throw RawDataError("Opening morphology '" + _uri +
                           "': incorrect number of dimensions in 'points'.");
    } else if (pointsDims[1] != pointColumns) {
        throw RawDataError("Opening morphology '" + _uri +
                           "': incorrect number of columns for points");
    }

    std::vector<std::array<floatType, pointColumns>> hdf5Data(numberPoints);

    if (!hdf5Data.empty()) {
        pointsDataSet.read(hdf5Data.front().data());
    }

    const bool hasSoma = firstSectionOffset != 0;
    const bool hasNeurites = static_cast<size_t>(firstSectionOffset) < numberPoints;
    const size_t somaPointCount = hasNeurites ? static_cast<size_t>(firstSectionOffset)
                                              : hdf5Data.size();

    auto& somaPoints = _properties._somaLevel._points;
    auto& somaDiameters = _properties._somaLevel._diameters;

    if (hasSoma) {
        somaPoints.resize(somaPointCount);
        somaDiameters.resize(somaPointCount);

        for (size_t i = 0; i < somaPointCount; ++i) {
            const auto& p = hdf5Data[i];
            somaPoints[i] = {p[0], p[1], p[2]};
            somaDiameters[i] = p[3];
        }
    }

    auto& points = _properties.get_mut<Property::Point>();
    auto& diameters = _properties.get_mut<Property::Diameter>();

    if (hasNeurites) {
        const size_t size = (hdf5Data.size() - somaPointCount);
        points.resize(size);
        diameters.resize(size);
        for (size_t i = somaPointCount; i < hdf5Data.size(); ++i) {
            const auto& p = hdf5Data[i];
            const size_t section = i - somaPointCount;
            points[section] = {p[0], p[1], p[2]};
            diameters[section] = p[3];
        }
    }
}

int MorphologyHDF5::_readSections() {
    // Important: The code used to split the reading of the sections and types
    //            into two separate fine-grained H5 selections. This does not
    //            reduce the number of I/O operations, but increases them by
    //            forcing HDF5 + MPI-IO to read in 4-byte groups. Thus, we now
    //            read the whole dataset at once, and split it in memory.

    constexpr size_t structureV1Columns = 3;

    const auto structure = _group.getDataSet(_d_structure);
    const auto dims = structure.getSpace().getDimensions();

    if (dims.size() != 2 || dims[1] != structureV1Columns) {
        throw(RawDataError("Error reading morphologies " + _uri +
                           " bad number of dimensions in 'structure' dataspace"));
    }

    std::vector<std::array<int, structureV1Columns>> vec(dims[0]);
    if (dims[0] > 0) {
        structure.read(vec.front().data());
    }

    assert(!vec.empty());

    bool hasSoma = true;
    if (static_cast<SectionType>(vec[0][SECTION_TYPE]) != SECTION_SOMA) {
        hasSoma = false;
    } else if (vec.size() == 1) {
        return SOMA_ONLY;
    }

    const size_t firstSection = hasSoma ? 1 : 0;
    const int firstSectionOffset = vec[firstSection][SECTION_START_OFFSET];

    auto& sections = _properties.get_mut<Property::Section>();
    sections.reserve(vec.size() - firstSection);

    auto& types = _properties.get_mut<Property::SectionType>();
    types.reserve(vec.size() - firstSection);

    // The first section is skipped if it corresponds to a soma
    for (size_t i = firstSection; i < vec.size(); ++i) {
        const auto& section = vec[i];
        SectionType type = static_cast<SectionType>(section[SECTION_TYPE]);

        if (section[SECTION_TYPE] >= SECTION_OUT_OF_RANGE_START || section[SECTION_TYPE] <= 0) {
            ErrorMessages err;
            throw RawDataError(err.ERROR_UNSUPPORTED_SECTION_TYPE(0, type));
        } else if (!hasSoma && type == SECTION_SOMA) {
            throw(RawDataError("Error reading morphology " + _uri +
                               ": it has soma section that doesn't come first"));
        } else if (hasSoma && type == SECTION_SOMA) {
            throw(RawDataError("Error reading morphology " + _uri +
                               ": it has multiple soma sections"));
        }

        sections.emplace_back(
            Property::Section::Type{section[SECTION_START_OFFSET] - firstSectionOffset,
                                    section[SECTION_PARENT_OFFSET] - (hasSoma ? 1 : 0)});
        types.emplace_back(type);
    }

    return firstSectionOffset;
}

void MorphologyHDF5::_readPerimeters(int firstSectionOffset) {
    if (!(_properties._cellLevel.majorVersion() == 1 &&
          _properties._cellLevel.minorVersion() > 0)) {
        throw RawDataError("Perimeter information is available starting at v1.1");
    }

    // soma only, won't have perimeters
    if (firstSectionOffset == SOMA_ONLY) {
        return;
    }

    if (!_group.exist(_d_perimeters)) {
        if (_properties._cellLevel._cellFamily == GLIA) {
            throw RawDataError("No empty perimeters allowed for glia morphology");
        }
        return;
    }

    auto& perimeters = _properties.get_mut<Property::Perimeter>();
    _read("", _d_perimeters, 1, perimeters);
    perimeters.erase(perimeters.begin(), perimeters.begin() + firstSectionOffset);
}


template <typename T>
void MorphologyHDF5::_read(const std::string& groupName,
                           const std::string& datasetName,
                           unsigned int expectedDimension,
                           T& data) {
    if (groupName != "" && !_group.exist(groupName)) {
        throw(
            RawDataError("Reading morphology '" + _uri + "': Missing required group " + groupName));
    }
    const auto group = groupName == "" ? _group : _group.getGroup(groupName);

    if (!group.exist(datasetName)) {
        throw(RawDataError("Reading morphology '" + _uri + "': Missing required dataset " +
                           datasetName));
    }
    const HighFive::DataSet dataset = group.getDataSet(datasetName);

    const auto dims = dataset.getSpace().getDimensions();
    if (dims.size() != expectedDimension) {
        throw(RawDataError("Reading morphology '" + _uri + "': bad number of dimensions in " +
                           datasetName));
    }

    data.resize(dims[0]);
    dataset.read(data);
}

void MorphologyHDF5::_readDendriticSpinePostSynapticDensity() {
    std::vector<morphio::Property::DendriticSpine::SectionId_t> sectionIds;
    _read(_g_postsynaptic_density, _d_dendritic_spine_section_id, 1, sectionIds);

    std::vector<morphio::Property::DendriticSpine::SegmentId_t> segmentIds;
    _read(_g_postsynaptic_density, _d_dendritic_spine_segment_id, 1, segmentIds);

    std::vector<morphio::Property::DendriticSpine::Offset_t> offsets;
    _read(_g_postsynaptic_density, _d_dendritic_spine_offset, 1, offsets);

    if (sectionIds.size() != segmentIds.size() || offsets.size() != segmentIds.size()) {
        throw(RawDataError(
            "Dendritic datasets must match in size:"
            " sectionIds: " +
            std::to_string(sectionIds.size()) + " segmentIds: " +
            std::to_string(segmentIds.size()) + " offsets: " + std::to_string(offsets.size())));
    }

    auto& properties = _properties._dendriticSpineLevel._post_synaptic_density;

    properties.reserve(sectionIds.size());
    for (size_t i = 0; i < sectionIds.size(); ++i) {
        properties.push_back({sectionIds[i], segmentIds[i], offsets[i]});
    }
}

void MorphologyHDF5::_readEndoplasmicReticulum() {
    if (!_group.exist(_g_endoplasmic_reticulum)) {
        return;
    }

    _read(_g_endoplasmic_reticulum,
          _d_section_index,
          1,
          _properties._endoplasmicReticulumLevel._sectionIndices);
    _read(_g_endoplasmic_reticulum,
          _d_volume,
          1,
          _properties._endoplasmicReticulumLevel._volumes);
    _read(_g_endoplasmic_reticulum,
          _d_surface_area,
          1,
          _properties._endoplasmicReticulumLevel._surfaceAreas);
    _read(_g_endoplasmic_reticulum,
          _d_filament_count,
          1,
          _properties._endoplasmicReticulumLevel._filamentCounts);
}

void MorphologyHDF5::_readMitochondria() {
    if (!_group.exist(_g_mitochondria)) {
        return;
    }

    std::vector<std::vector<floatType>> points;
    _read(_g_mitochondria, _d_points, 2, points);

    auto& mitoSectionId = _properties.get_mut<Property::MitoNeuriteSectionId>();
    auto& pathlength = _properties.get_mut<Property::MitoPathLength>();
    auto& diameters = _properties.get_mut<Property::MitoDiameter>();
    mitoSectionId.reserve(mitoSectionId.size() + points.size());
    pathlength.reserve(pathlength.size() + points.size());
    diameters.reserve(diameters.size() + points.size());
    for (const auto& p : points) {
        mitoSectionId.push_back(static_cast<Property::MitoNeuriteSectionId::Type>(p[0]));
        pathlength.push_back(p[1]);
        diameters.push_back(p[2]);
    }

    std::vector<std::vector<int32_t>> structure;
    _read(_g_mitochondria, _d_structure, 2, structure);

    auto& mitoSection = _properties.get_mut<Property::MitoSection>();
    mitoSection.reserve(mitoSection.size() + structure.size());
    for (auto& s : structure)
        mitoSection.emplace_back(Property::MitoSection::Type{s[0], s[1]});
}

}  // namespace h5
}  // namespace readers
}  // namespace morphio

1	/* Copyright (c) 2013-2023, EPFL/Blue Brain Project
2	*
3	* SPDX-License-Identifier: Apache-2.0
4	*/
5
6	#include <cassert>
7
8	#include "morphologyHDF5.h"
9
10	#include <highfive/H5Utility.hpp> // HighFive::SilenceHDF5
11
12	#include <morphio/errorMessages.h>
13
14	namespace {
15
16	constexpr size_t SECTION_START_OFFSET = 0;
17	constexpr size_t SECTION_TYPE = 1;
18	constexpr size_t SECTION_PARENT_OFFSET = 2;
19
20	constexpr int SOMA_ONLY = -1;
21
22	//{v1
23	const std::string _d_structure("structure");
24	const std::string _d_points("points");
25
26	//{ v1.1
27	const std::string _a_version("version");
28	const std::string _g_metadata("metadata");
29	const std::string _a_family("cell_family");
30	const std::string _d_perimeters("perimeters");
31	//} v1.1
32
33	//{ v1.2
34	const std::string _g_mitochondria("organelles/mitochondria");
35
36	// endoplasmic reticulum
37	const std::string _g_endoplasmic_reticulum("organelles/endoplasmic_reticulum");
38	const std::string _d_section_index("section_index");
39	const std::string _d_volume("volume");
40	const std::string _d_surface_area("surface_area");
41	const std::string _d_filament_count("filament_count");
42	// } v1.2
43
44	//{ v1.3
45
46	// Dendritic Spine
47	const std::string _g_postsynaptic_density("organelles/postsynaptic_density");
48	const std::string _d_dendritic_spine_section_id("section_id");
49	const std::string _d_dendritic_spine_segment_id("segment_id");
50	const std::string _d_dendritic_spine_offset("offset");
51
52	//}
53
54	// } v1
55
56	//{ v2
57	const std::string _g_v2root("neuron1");
58	//} v2
59
60	} // namespace
61
62	namespace morphio {
63	namespace readers {
64	namespace h5 {
65
66	MorphologyHDF5::MorphologyHDF5(const HighFive::Group& group)	248✔
67	: _group(group)
68	, _uri("HDF5 Group") {}	248✔
69
70	Property::Properties load(const std::string& uri) {	178✔
71	try {
72	std::lock_guard<std::recursive_mutex> lock(morphio::readers::h5::global_hdf5_mutex());	356✔
73	HighFive::SilenceHDF5 silence;	356✔
74	auto file = HighFive::File(uri, HighFive::File::ReadOnly);	206✔
75	return MorphologyHDF5(file.getGroup("/")).load();	584✔
76
77	} catch (const HighFive::FileException& exc) {	4✔
78	throw RawDataError("Could not open morphology file " + uri + ": " + exc.what());	2✔
79	}
80	}
81
82	Property::Properties load(const HighFive::Group& group) {	70✔
83	std::lock_guard<std::recursive_mutex> lock(morphio::readers::h5::global_hdf5_mutex());	70✔
84	return MorphologyHDF5(group).load();	210✔
85	}
86
87	Property::Properties MorphologyHDF5::load() {	246✔
88	_readMetadata(_uri);	246✔
89
90	int firstSectionOffset = _readSections();	238✔
91
92	_readPoints(firstSectionOffset);	230✔
93
94	if (_properties._cellLevel.minorVersion() >= 1) {	226✔
95	_readPerimeters(firstSectionOffset);	218✔
96
97	if (_properties._cellLevel.minorVersion() >= 2) {	210✔
98	_readMitochondria();	152✔
99	_readEndoplasmicReticulum();	152✔
100	}
101
102	if (_properties._cellLevel.minorVersion() >= 3 &&	238✔
103	_properties._cellLevel._cellFamily == CellFamily::SPINE) {	28✔
104	_readDendriticSpinePostSynapticDensity();	26✔
105	}
106	}
107
108	return _properties;	218✔
109	}
110
111	void MorphologyHDF5::_readMetadata(const std::string& source) {	246✔
112	// default to h5v1.0
113	uint32_t majorVersion = 1;	246✔
114	uint32_t minorVersion = 0;	246✔
115	_properties._cellLevel._cellFamily = CellFamily::NEURON;	246✔
116
117	if (!_group.exist(_d_points) \|\| !_group.exist(_d_structure)) {	246✔
118	// h5v2 is deprecated, but it can be detected, throw a custom error messages if it is
119	if (_group.exist(_g_v2root)) {	4✔
120	throw RawDataError(
121	"Error in " + source +	4✔
122	"\nh5v2 is no longer supported, see: https://github.com/BlueBrain/MorphIO#H5v2");	6✔
123	}
124	throw RawDataError("Missing " + _d_points + " or " + _d_structure +	4✔
125	" datasets, cannot load morphology without them");	6✔
126	}
127
128	// if there is metadata, perhaps it's h5v1, minor version 1, 2
129	if (_group.exist(_g_metadata)) {	242✔
130	const auto metadata = _group.getGroup(_g_metadata);	448✔
131	if (metadata.hasAttribute(_a_version)) {	224✔
132	const auto attr = metadata.getAttribute(_a_version);	444✔
133
134	std::array<uint32_t, 2> versions = {0, 0};	222✔
135	attr.read(versions);	222✔
136	majorVersion = versions[0];	222✔
137	minorVersion = versions[1];	222✔
138
139	if (majorVersion == 1 &&	222✔
140	(minorVersion == 1 \|\| minorVersion == 2 \|\| minorVersion == 3)) {	222✔
141	uint32_t family;
142	metadata.getAttribute(_a_family).read(family);	220✔
143	_properties._cellLevel._cellFamily = static_cast<CellFamily>(family);	220✔
144	} else {
145	std::string msg = "Error in " + source +	4✔
146	"\nUnsupported h5 version: " + std::to_string(majorVersion) +	8✔
147	"." + std::to_string(minorVersion);	6✔
148	throw RawDataError(msg);	2✔
149	}
150	} else {
151	throw RawDataError("Missing " + _a_version +	4✔
152	" attribute, cannot load morphology without them");	6✔
153	}
154	}
155
156	_properties._cellLevel._version = {"h5", majorVersion, minorVersion};	238✔
157	}	238✔
158
159	void MorphologyHDF5::_readPoints(int firstSectionOffset) {	230✔
160	constexpr size_t pointColumns = 4;	230✔
161
162	const auto pointsDataSet = _group.getDataSet(_d_points);	460✔
163	const auto pointsDims = pointsDataSet.getSpace().getDimensions();	460✔
164	const size_t numberPoints = pointsDims[0];	230✔
165
166	if (pointsDims.size() != 2) {	230✔
167	throw RawDataError("Opening morphology '" + _uri +	4✔
168	"': incorrect number of dimensions in 'points'.");	6✔
169	} else if (pointsDims[1] != pointColumns) {	228✔
170	throw RawDataError("Opening morphology '" + _uri +	4✔
171	"': incorrect number of columns for points");	6✔
172	}
173
174	std::vector<std::array<floatType, pointColumns>> hdf5Data(numberPoints);	452✔
175
176	if (!hdf5Data.empty()) {	226✔
177	pointsDataSet.read(hdf5Data.front().data());	226✔
178	}
179
180	const bool hasSoma = firstSectionOffset != 0;	226✔
181	const bool hasNeurites = static_cast<size_t>(firstSectionOffset) < numberPoints;	226✔
182	const size_t somaPointCount = hasNeurites ? static_cast<size_t>(firstSectionOffset)	226✔
183	: hdf5Data.size();	4✔
184
185	auto& somaPoints = _properties._somaLevel._points;	226✔
186	auto& somaDiameters = _properties._somaLevel._diameters;	226✔
187
188	if (hasSoma) {	226✔
189	somaPoints.resize(somaPointCount);	196✔
190	somaDiameters.resize(somaPointCount);	196✔
191
192	for (size_t i = 0; i < somaPointCount; ++i) {	88,596✔
193	const auto& p = hdf5Data[i];	88,400✔
194	somaPoints[i] = {p[0], p[1], p[2]};	88,400✔
195	somaDiameters[i] = p[3];	88,400✔
196	}
197	}
198
199	auto& points = _properties.get_mut<Property::Point>();	226✔
200	auto& diameters = _properties.get_mut<Property::Diameter>();	226✔
201
202	if (hasNeurites) {	226✔
203	const size_t size = (hdf5Data.size() - somaPointCount);	222✔
204	points.resize(size);	222✔
205	diameters.resize(size);	222✔
206	for (size_t i = somaPointCount; i < hdf5Data.size(); ++i) {	48,892✔
207	const auto& p = hdf5Data[i];	48,670✔
208	const size_t section = i - somaPointCount;	48,670✔
209	points[section] = {p[0], p[1], p[2]};	48,670✔
210	diameters[section] = p[3];	48,670✔
211	}
212	}
213	}	226✔
214
215	int MorphologyHDF5::_readSections() {	238✔
216	// Important: The code used to split the reading of the sections and types
217	// into two separate fine-grained H5 selections. This does not
218	// reduce the number of I/O operations, but increases them by
219	// forcing HDF5 + MPI-IO to read in 4-byte groups. Thus, we now
220	// read the whole dataset at once, and split it in memory.
221
222	constexpr size_t structureV1Columns = 3;	238✔
223
224	const auto structure = _group.getDataSet(_d_structure);	476✔
225	const auto dims = structure.getSpace().getDimensions();	476✔
226
227	if (dims.size() != 2 \|\| dims[1] != structureV1Columns) {	238✔
228	throw(RawDataError("Error reading morphologies " + _uri +	4✔
229	" bad number of dimensions in 'structure' dataspace"));	6✔
230	}
231
232	std::vector<std::array<int, structureV1Columns>> vec(dims[0]);	472✔
233	if (dims[0] > 0) {	236✔
234	structure.read(vec.front().data());	236✔
235	}
236
237	assert(!vec.empty());	236✔
238
239	bool hasSoma = true;	236✔
240	if (static_cast<SectionType>(vec[0][SECTION_TYPE]) != SECTION_SOMA) {	236✔
241	hasSoma = false;	34✔
242	} else if (vec.size() == 1) {	202✔
243	return SOMA_ONLY;	6✔
244	}
245
246	const size_t firstSection = hasSoma ? 1 : 0;	230✔
247	const int firstSectionOffset = vec[firstSection][SECTION_START_OFFSET];	230✔
248
249	auto& sections = _properties.get_mut<Property::Section>();	230✔
250	sections.reserve(vec.size() - firstSection);	230✔
251
252	auto& types = _properties.get_mut<Property::SectionType>();	230✔
253	types.reserve(vec.size() - firstSection);	230✔
254
255	// The first section is skipped if it corresponds to a soma
256	for (size_t i = firstSection; i < vec.size(); ++i) {	7,024✔
257	const auto& section = vec[i];	6,800✔
258	SectionType type = static_cast<SectionType>(section[SECTION_TYPE]);	6,800✔
259
260	if (section[SECTION_TYPE] >= SECTION_OUT_OF_RANGE_START \|\| section[SECTION_TYPE] <= 0) {	6,800✔
261	ErrorMessages err;	4✔
262	throw RawDataError(err.ERROR_UNSUPPORTED_SECTION_TYPE(0, type));	2✔
263	} else if (!hasSoma && type == SECTION_SOMA) {	6,798✔
264	throw(RawDataError("Error reading morphology " + _uri +	4✔
265	": it has soma section that doesn't come first"));	6✔
266	} else if (hasSoma && type == SECTION_SOMA) {	6,796✔
267	throw(RawDataError("Error reading morphology " + _uri +	4✔
268	": it has multiple soma sections"));	6✔
269	}
270
271	sections.emplace_back(	6,794✔
272	Property::Section::Type{section[SECTION_START_OFFSET] - firstSectionOffset,	6,794✔
273	section[SECTION_PARENT_OFFSET] - (hasSoma ? 1 : 0)});	6,794✔
274	types.emplace_back(type);	6,794✔
275	}
276
277	return firstSectionOffset;	224✔
278	}
279
280	void MorphologyHDF5::_readPerimeters(int firstSectionOffset) {	218✔
281	if (!(_properties._cellLevel.majorVersion() == 1 &&	436✔
282	_properties._cellLevel.minorVersion() > 0)) {	218✔
283	throw RawDataError("Perimeter information is available starting at v1.1");	×
284	}
285
286	// soma only, won't have perimeters
287	if (firstSectionOffset == SOMA_ONLY) {	218✔
288	return;	4✔
289	}
290
291	if (!_group.exist(_d_perimeters)) {	214✔
292	if (_properties._cellLevel._cellFamily == GLIA) {	126✔
293	throw RawDataError("No empty perimeters allowed for glia morphology");	4✔
294	}
295	return;	122✔
296	}
297
298	auto& perimeters = _properties.get_mut<Property::Perimeter>();	88✔
299	_read("", _d_perimeters, 1, perimeters);	96✔
300	perimeters.erase(perimeters.begin(), perimeters.begin() + firstSectionOffset);	84✔
301	}
302
303
304	template <typename T>
305	void MorphologyHDF5::_read(const std::string& groupName,	494✔
306	const std::string& datasetName,
307	unsigned int expectedDimension,
308	T& data) {
309	if (groupName != "" && !_group.exist(groupName)) {	494✔
310	throw(
311	RawDataError("Reading morphology '" + _uri + "': Missing required group " + groupName));	×
312	}
313	const auto group = groupName == "" ? _group : _group.getGroup(groupName);	988✔
314
315	if (!group.exist(datasetName)) {	494✔
316	throw(RawDataError("Reading morphology '" + _uri + "': Missing required dataset " +	×
317	datasetName));	×
318	}
319	const HighFive::DataSet dataset = group.getDataSet(datasetName);	988✔
320
321	const auto dims = dataset.getSpace().getDimensions();	498✔
322	if (dims.size() != expectedDimension) {	494✔
323	throw(RawDataError("Reading morphology '" + _uri + "': bad number of dimensions in " +	4✔
324	datasetName));	4✔
325	}
326
327	data.resize(dims[0]);	490✔
328	dataset.read(data);	490✔
329	}	490✔
330
331	void MorphologyHDF5::_readDendriticSpinePostSynapticDensity() {	26✔
332	std::vector<morphio::Property::DendriticSpine::SectionId_t> sectionIds;	52✔
333	_read(_g_postsynaptic_density, _d_dendritic_spine_section_id, 1, sectionIds);	26✔
334
335	std::vector<morphio::Property::DendriticSpine::SegmentId_t> segmentIds;	52✔
336	_read(_g_postsynaptic_density, _d_dendritic_spine_segment_id, 1, segmentIds);	26✔
337
338	std::vector<morphio::Property::DendriticSpine::Offset_t> offsets;	52✔
339	_read(_g_postsynaptic_density, _d_dendritic_spine_offset, 1, offsets);	26✔
340
341	if (sectionIds.size() != segmentIds.size() \|\| offsets.size() != segmentIds.size()) {	26✔
342	throw(RawDataError(
343	"Dendritic datasets must match in size:"
344	" sectionIds: " +	×
345	std::to_string(sectionIds.size()) + " segmentIds: " +	×
346	std::to_string(segmentIds.size()) + " offsets: " + std::to_string(offsets.size())));	×
347	}
348
349	auto& properties = _properties._dendriticSpineLevel._post_synaptic_density;	26✔
350
351	properties.reserve(sectionIds.size());	26✔
352	for (size_t i = 0; i < sectionIds.size(); ++i) {	78✔
353	properties.push_back({sectionIds[i], segmentIds[i], offsets[i]});	52✔
354	}
355	}	26✔
356
357	void MorphologyHDF5::_readEndoplasmicReticulum() {	152✔
358	if (!_group.exist(_g_endoplasmic_reticulum)) {	152✔
359	return;	112✔
360	}
361
362	_read(_g_endoplasmic_reticulum,	40✔
363	_d_section_index,
364	1,
365	_properties._endoplasmicReticulumLevel._sectionIndices);	40✔
366	_read(_g_endoplasmic_reticulum,	40✔
367	_d_volume,
368	1,
369	_properties._endoplasmicReticulumLevel._volumes);	40✔
370	_read(_g_endoplasmic_reticulum,	40✔
371	_d_surface_area,
372	1,
373	_properties._endoplasmicReticulumLevel._surfaceAreas);	40✔
374	_read(_g_endoplasmic_reticulum,	40✔
375	_d_filament_count,
376	1,
377	_properties._endoplasmicReticulumLevel._filamentCounts);	40✔
378	}
379
380	void MorphologyHDF5::_readMitochondria() {	152✔
381	if (!_group.exist(_g_mitochondria)) {	152✔
382	return;	68✔
383	}
384
385	std::vector<std::vector<floatType>> points;	168✔
386	_read(_g_mitochondria, _d_points, 2, points);	84✔
387
388	auto& mitoSectionId = _properties.get_mut<Property::MitoNeuriteSectionId>();	84✔
389	auto& pathlength = _properties.get_mut<Property::MitoPathLength>();	84✔
390	auto& diameters = _properties.get_mut<Property::MitoDiameter>();	84✔
391	mitoSectionId.reserve(mitoSectionId.size() + points.size());	84✔
392	pathlength.reserve(pathlength.size() + points.size());	84✔
393	diameters.reserve(diameters.size() + points.size());	84✔
394	for (const auto& p : points) {	924✔
395	mitoSectionId.push_back(static_cast<Property::MitoNeuriteSectionId::Type>(p[0]));	840✔
396	pathlength.push_back(p[1]);	840✔
397	diameters.push_back(p[2]);	840✔
398	}
399
400	std::vector<std::vector<int32_t>> structure;	168✔
401	_read(_g_mitochondria, _d_structure, 2, structure);	84✔
402
403	auto& mitoSection = _properties.get_mut<Property::MitoSection>();	84✔
404	mitoSection.reserve(mitoSection.size() + structure.size());	84✔
405	for (auto& s : structure)	336✔
406	mitoSection.emplace_back(Property::MitoSection::Type{s[0], s[1]});	252✔
407	}
408
409	} // namespace h5
410	} // namespace readers
411	} // namespace morphio

BlueBrain / MorphIO / 6533524408

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