From 9dae25f4200d20b63a83736a1b01180fb2a182df Mon Sep 17 00:00:00 2001 From: Daniel Nachbaur Date: Tue, 14 Jun 2016 11:28:28 +0200 Subject: [PATCH] Implement morphology version 1.1 specification --- brain/neuron/morphologyImpl.cpp | 9 +- brion/enums.h | 15 +- brion/morphology.cpp | 22 +- brion/morphology.h | 60 ++++- brion/morphologyPlugin.h | 30 ++- brion/plugin/morphologyHDF5.cpp | 550 +++++++++++++++++++++++++++------------- brion/plugin/morphologyHDF5.h | 27 +- brion/plugin/morphologySWC.cpp | 18 +- brion/plugin/morphologySWC.h | 9 +- doc/Changelog.md | 2 + doc/feature/morphology11.md | 22 +- tests/morphology.cpp | 120 ++++++++- 12 files changed, 667 insertions(+), 217 deletions(-) diff --git a/brain/neuron/morphologyImpl.cpp b/brain/neuron/morphologyImpl.cpp index 1dc4fdc..2818d37 100644 --- a/brain/neuron/morphologyImpl.cpp +++ b/brain/neuron/morphologyImpl.cpp @@ -1,5 +1,5 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Juan Hernando * * This file is part of Brion @@ -73,11 +73,8 @@ uint32_ts Morphology::Impl::getSectionIDs( for( size_t i = 0; i != types->size(); ++i ) { const SectionType type = ( *types )[i]; - if( type == SECTION_ALL ) - LBWARN << "Unknown section type " << int(type) << std::endl; - else - if( bits[size_t( type )] ) - result.push_back( i ); + if( bits[size_t( type )] ) + result.push_back( i ); } return result; } diff --git a/brion/enums.h b/brion/enums.h index 1f9c6f7..58de379 100644 --- a/brion/enums.h +++ b/brion/enums.h @@ -1,4 +1,4 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Daniel Nachbaur * * This file is part of Brion @@ -79,14 +79,20 @@ enum MorphologyRepairStage /** The supported versions for morphology files. */ enum MorphologyVersion { - MORPHOLOGY_VERSION_1 = 1, //!< @deprecated use MORPHOLOGY_VERSION_H5_1 - MORPHOLOGY_VERSION_2 = 2, //!< @deprecated use MORPHOLOGY_VERSION_H5_2 MORPHOLOGY_VERSION_H5_1 = 1, MORPHOLOGY_VERSION_H5_2 = 2, + MORPHOLOGY_VERSION_H5_1_1 = 3, MORPHOLOGY_VERSION_SWC_1 = 101, MORPHOLOGY_VERSION_UNDEFINED }; +/** The cell family represented by brion::Morphology. */ +enum CellFamily +{ + FAMILY_NEURON = 0, + FAMILY_GLIA = 1 +}; + /** Output stream formatter for MorphologyVersion */ inline std::ostream& operator << ( std::ostream& os, const MorphologyVersion v ) { @@ -150,7 +156,8 @@ enum SectionType SECTION_AXON = 2, SECTION_DENDRITE = 3, //!< general or basal dendrite (near to soma) SECTION_APICAL_DENDRITE = 4, //!< apical dendrite (far from soma) - SECTION_ALL //!< @internal must be last + SECTION_GLIA_PROCESS = 2, + SECTION_GLIA_ENDFOOT = 3 }; /** The supported attributes of a synapse. */ diff --git a/brion/morphology.cpp b/brion/morphology.cpp index 463e775..a546535 100644 --- a/brion/morphology.cpp +++ b/brion/morphology.cpp @@ -1,4 +1,4 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Daniel Nachbaur * * This file is part of Brion @@ -57,11 +57,21 @@ Morphology::Morphology( const std::string& target, { } +Morphology::Morphology( const std::string& file, const CellFamily family ) + : _impl( new Impl( MorphologyInitData( URI( file ), family ))) +{ +} + Morphology::~Morphology() { delete _impl; } +CellFamily Morphology::getCellFamily() const +{ + return _impl->plugin->getCellFamily(); +} + Vector4fsPtr Morphology::readPoints( const MorphologyRepairStage stage ) const { return _impl->plugin->readPoints( stage ); @@ -83,6 +93,11 @@ Vector2isPtr Morphology::readApicals() const return _impl->plugin->readApicals(); } +floatsPtr Morphology::readPerimeters() const +{ + return _impl->plugin->readPerimeters(); +} + MorphologyVersion Morphology::getVersion() const { return _impl->plugin->getVersion(); @@ -110,6 +125,11 @@ void Morphology::writeApicals( const Vector2is& apicals ) _impl->plugin->writeApicals( apicals ); } +void Morphology::writePerimeters( const floats& perimeters ) +{ + _impl->plugin->writePerimeters( perimeters ); +} + void Morphology::flush() { _impl->plugin->flush(); diff --git a/brion/morphology.h b/brion/morphology.h index 32d43e1..98829c4 100644 --- a/brion/morphology.h +++ b/brion/morphology.h @@ -1,4 +1,4 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Daniel Nachbaur * * This file is part of Brion @@ -49,13 +49,17 @@ class Morphology : public boost::noncopyable */ BRION_API explicit Morphology( const std::string& source ); + /** @return the cell family of that morphology. @version 1.8 */ + BRION_API CellFamily getCellFamily() const; + /** Read points of morphology, representing x,y,z coordinates + diameter. * * @param stage the repair stage to take the points from * @return x,y,z coords + diameter of all points of the morphology * @version 1.0 */ - BRION_API Vector4fsPtr readPoints( MorphologyRepairStage stage ) const; + BRION_API Vector4fsPtr + readPoints( MorphologyRepairStage stage = MORPHOLOGY_UNDEFINED ) const; /** Read sections of morphology, representing section start index and index * of the parent section. @@ -65,7 +69,8 @@ class Morphology : public boost::noncopyable * morphology. * @version 1.0 */ - BRION_API Vector2isPtr readSections( MorphologyRepairStage stage ) const; + BRION_API Vector2isPtr + readSections( MorphologyRepairStage stage = MORPHOLOGY_UNDEFINED ) const; /** Read section types of morphology. * @@ -81,6 +86,14 @@ class Morphology : public boost::noncopyable */ BRION_API Vector2isPtr readApicals() const; + /** + * @return perimeters of the cross sections for each point of the morphology + * in micrometers. + * @throw std::runtime_error if empty for FAMILY_GLIA + * @version 1.8 + */ + BRION_API floatsPtr readPerimeters() const; + /** @internal */ BRION_API MorphologyVersion getVersion() const; //@} @@ -95,35 +108,49 @@ class Morphology : public boost::noncopyable * @param overwrite true to allow overwrite of existing file * @throw std::runtime_error if file could not be created * @version 1.0 + * @deprecated */ BRION_API Morphology( const std::string& target, MorphologyVersion version, bool overwrite = false ); + /** + * Open the given morphology file for write access. + * + * @param file filename of the output morphology file + * @param family the cell family that this morphology represents + * @throw std::runtime_error if file could not be created + * @version 1.8 + */ + BRION_API Morphology( const std::string& file, CellFamily family ); + /** Write points of morphology, representing x,y,z coordinates + diameter. * * @param points x,y,z coords + diameter of all points of the morphology * @param stage the repair stage to save those points to - * @throw std::runtime_error if object not created with write ctor + * @throw std::runtime_error if object not writable + * @throw std::runtime_error points already written + * @throw std::runtime_error number of points does not match number of + * perimeters * @version 1.0 */ BRION_API void writePoints( const Vector4fs& points, - MorphologyRepairStage stage ); + MorphologyRepairStage stage = MORPHOLOGY_UNDEFINED ); /** Write sections of morphology, representing section start index and index * of parent the section. * * @param sections index and parent index of all sections of the morphology * @param stage the repair stage to save those sections to - * @throw std::runtime_error if object not created with write ctor + * @throw std::runtime_error if object not writable * @version 1.0 */ BRION_API void writeSections( const Vector2is& sections, - MorphologyRepairStage stage ); + MorphologyRepairStage stage = MORPHOLOGY_UNDEFINED ); /** Write section types of morphology. * * @param types type of each section of the morphology - * @throw std::runtime_error if object not created with write ctor + * @throw std::runtime_error if object not writable * @version 1.0 */ BRION_API void writeSectionTypes( const SectionTypes& types ); @@ -132,11 +159,26 @@ class Morphology : public boost::noncopyable * * @param apicals section and point index of all apical points * @throw std::runtime_error if object not created with write ctor - * @throw std::runtime_error if called for version 1 files + * @throw std::runtime_error if not supported by implementation * @version 1.0 + * @deprecated */ BRION_API void writeApicals( const Vector2is& apicals ); + /** + * Write perimeters of morphology. + * + * @param perimeters perimeters of the cross sections for each point of + * the morphology in micrometers + * @throw std::runtime_error if object not writable + * @throw std::runtime_error perimeters already written + * @throw std::runtime_error number of points does not match number of + * perimeters + * @throw std::runtime_error if not supported by implementation + * @version 1.8 + */ + BRION_API void writePerimeters( const floats& perimeters ); + /** Flush data to output. @version 1.0 */ BRION_API void flush(); //@} diff --git a/brion/morphologyPlugin.h b/brion/morphologyPlugin.h index 80cc482..803317a 100644 --- a/brion/morphologyPlugin.h +++ b/brion/morphologyPlugin.h @@ -1,5 +1,6 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Juan Hernando + * Daniel Nachbaur * * This file is part of Brion * @@ -40,6 +41,8 @@ class MorphologyInitData : public PluginInitData public: explicit MorphologyInitData( const URI& uri ) : PluginInitData( uri, MODE_READ ) + , _version( MORPHOLOGY_VERSION_H5_1_1 ) + , _family( FAMILY_NEURON ) {} MorphologyInitData( const URI& uri, @@ -47,6 +50,14 @@ class MorphologyInitData : public PluginInitData const unsigned int accessMode ) : PluginInitData( uri, accessMode ) , _version( version ) + , _family( FAMILY_NEURON ) + {} + + MorphologyInitData( const URI& uri, + const CellFamily family ) + : PluginInitData( uri, MODE_WRITE ) + , _version( MORPHOLOGY_VERSION_H5_1_1 ) + , _family( family ) {} MorphologyVersion getVersion() const @@ -54,8 +65,14 @@ class MorphologyInitData : public PluginInitData return _version; } + CellFamily getFamily() const + { + return _family; + } + protected: - MorphologyVersion _version; + const MorphologyVersion _version; + const CellFamily _family; }; /** Base interface for morphology readers plugins. @@ -91,6 +108,9 @@ class MorphologyPlugin : public boost::noncopyable /** @name Read API */ //@{ + /** @copydoc brion::Morphology::getCellFamily */ + virtual CellFamily getCellFamily() const = 0; + /** @copydoc brion::Morphology::readPoints */ virtual Vector4fsPtr readPoints( MorphologyRepairStage stage ) const = 0; @@ -103,6 +123,9 @@ class MorphologyPlugin : public boost::noncopyable /** @copydoc brion::Morphology::readApicals */ virtual Vector2isPtr readApicals() const = 0; + /** @copydoc brion::Morphology::readPerimeters */ + virtual floatsPtr readPerimeters() const = 0; + virtual MorphologyVersion getVersion() const = 0; //@} @@ -122,6 +145,9 @@ class MorphologyPlugin : public boost::noncopyable /** @copydoc brion::Morphology::writeApicals */ virtual void writeApicals( const Vector2is& apicals ) = 0; + /** @copydoc brion::Morphology::writePerimeters */ + virtual void writePerimeters( const floats& perimeters ) = 0; + /** @copydoc brion::Morphology::flush */ virtual void flush() = 0; //@} diff --git a/brion/plugin/morphologyHDF5.cpp b/brion/plugin/morphologyHDF5.cpp index 175eece..1a1fddb 100644 --- a/brion/plugin/morphologyHDF5.cpp +++ b/brion/plugin/morphologyHDF5.cpp @@ -1,4 +1,4 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Daniel Nachbaur * Juan Hernando * @@ -48,14 +48,25 @@ const size_t _pointColumns = 4; const std::string _d_structure( "/structure" ); const size_t _structureV1Columns = 3; +// v1.1 +const std::string _g_metadata( "/metadata" ); +const std::string _e_family( "cell_family_enum" ); +const std::string _a_family( "cell_family" ); +const std::string _d_perimeters( "/perimeters" ); +const std::string _a_software_version( "software_version" ); +const std::string _a_creation_time( "creation_time" ); + +// v1.1 & v2 +const std::string _a_creator( "creator" ); +const std::string _a_version( "version" ); + // v2 const std::string _g_structure( "structure" ); const size_t _structureV2Columns = 2; const std::string _g_root( "neuron1" ); const std::string _d_type( "sectiontype" ); const std::string _a_apical( "apical" ); -const std::string _a_creator( "creator" ); -const std::string _a_version( "version" ); + std::string toString( const brion::enums::MorphologyRepairStage& s ) { @@ -82,8 +93,9 @@ lunchbox::PluginRegisterer< MorphologyHDF5 > registerer; MorphologyHDF5::MorphologyHDF5( const MorphologyInitData& initData ) : _file() - , _version( MORPHOLOGY_VERSION_H5_2 ) + , _version( MORPHOLOGY_VERSION_UNDEFINED ) , _stage( MORPHOLOGY_UNDEFINED ) + , _family( FAMILY_NEURON ) , _write( false ) { lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); @@ -131,13 +143,16 @@ MorphologyHDF5::MorphologyHDF5( const MorphologyInitData& initData ) switch( _version ) { case MORPHOLOGY_VERSION_H5_1: - _writeV1Header(); + // no metadata for version 1 + break; + case MORPHOLOGY_VERSION_H5_1_1: + _writeV11Metadata( initData ); break; case MORPHOLOGY_VERSION_UNDEFINED: _version = MORPHOLOGY_VERSION_H5_2; // no break; case MORPHOLOGY_VERSION_H5_2: - _writeV2Header(); + _writeV2Metadata(); break; default: LBTHROW( std::runtime_error( "Unknown morphology file format for " @@ -169,43 +184,50 @@ bool MorphologyHDF5::handles( const MorphologyInitData& initData ) return path.substr( pos ) == ".h5"; } +CellFamily MorphologyHDF5::getCellFamily() const +{ + return _family; +} + Vector4fsPtr MorphologyHDF5::readPoints( MorphologyRepairStage stage ) const { lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); - H5::DataSet dataset; - try + if( _version == MORPHOLOGY_VERSION_H5_2 ) { - if( _version == MORPHOLOGY_VERSION_H5_1 ) - dataset = _file.openDataSet( "/" + _d_points ); - else + H5::DataSet dataset; + try { if( stage == MORPHOLOGY_UNDEFINED ) stage = _stage; dataset = _file.openDataSet( "/" + _g_root + "/" + toString( stage ) + "/" + _d_points ); } - } - catch( ... ) - { - LBERROR << "Could not open points dataset for morphology file " - << _file.getFileName() << " repair stage " - << toString( stage ) << std::endl; - return Vector4fsPtr( new Vector4fs ); - } + catch( ... ) + { + LBERROR << "Could not open points dataset for morphology file " + << _file.getFileName() << " repair stage " + << toString( stage ) << std::endl; + return Vector4fsPtr( new Vector4fs ); + } - hsize_t dims[2]; - const H5::DataSpace& dspace = dataset.getSpace(); - if( dspace.getSimpleExtentNdims() != 2 || - dspace.getSimpleExtentDims( dims ) < 0 || dims[1] != _pointColumns ) - { - LBTHROW( std::runtime_error( "Reading morphology file '" + - _file.getFileName() + "': bad number of dimensions in" - " 'points' dataspace")); + hsize_t dims[2]; + const H5::DataSpace& dspace = dataset.getSpace(); + if( dspace.getSimpleExtentNdims() != 2 || + dspace.getSimpleExtentDims( dims ) < 0 || dims[1] != _pointColumns ) + { + LBTHROW( std::runtime_error( "Reading morphology file '" + + _file.getFileName() + "': bad number of dimensions in" + " 'points' dataspace")); + } + + Vector4fsPtr data( new Vector4fs( dims[0] )); + dataset.read( data->data(), H5::PredType::NATIVE_FLOAT ); + return data; } - Vector4fsPtr data( new Vector4fs( dims[0] )); - dataset.read( data->data(), H5::PredType::NATIVE_FLOAT ); + Vector4fsPtr data( new Vector4fs( _pointsDims[0] )); + _points.read( data->data(), H5::PredType::NATIVE_FLOAT ); return data; } @@ -213,67 +235,55 @@ Vector2isPtr MorphologyHDF5::readSections( MorphologyRepairStage stage ) const { lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); - if( _version == MORPHOLOGY_VERSION_H5_1 ) + if( _version == MORPHOLOGY_VERSION_H5_2 ) { - const H5::DataSet& dataset = _file.openDataSet( _d_structure ); + if( stage == MORPHOLOGY_UNDEFINED ) + stage = _stage; + + // fixes BBPSDK-295 by restoring old BBPSDK 0.13 implementation + if( stage == MORPHOLOGY_UNRAVELED ) + stage = MORPHOLOGY_RAW; + + H5::DataSet dataset; + try + { + dataset = _file.openDataSet( "/" + _g_root + "/" + _g_structure + + "/" + toString( stage )); + } + catch( ... ) + { + LBERROR << "Could not open sections dataset for morphology file " + << _file.getFileName() << " repair stage " + << toString( stage ) << std::endl; + return Vector2isPtr( new Vector2is ); + } hsize_t dims[2]; const H5::DataSpace& dspace = dataset.getSpace(); if( dspace.getSimpleExtentNdims() != 2 || dspace.getSimpleExtentDims( dims ) < 0 || - dims[1] != _structureV1Columns ) + dims[1] != _structureV2Columns ) { LBTHROW( std::runtime_error( "Reading morphology file '" + _file.getFileName() + "': bad number of dimensions in" " 'structure' dataspace")); } - const hsize_t readCount[2] = { dims[0], 1 }; - const hsize_t readOffset[2] = { 0, 1 }; - dspace.selectHyperslab( H5S_SELECT_XOR, readCount, readOffset ); - Vector2isPtr data( new Vector2is( dims[0] )); - const hsize_t mdim[2] = { dims[0], 2 }; - const H5::DataSpace mspace( 2, mdim ); - dataset.read( data->data(), H5::PredType::NATIVE_INT, mspace, - dspace ); + dataset.read( data->data(), H5::PredType::NATIVE_INT ); return data; } - if( stage == MORPHOLOGY_UNDEFINED ) - stage = _stage; - - // fixes BBPSDK-295 by restoring old BBPSDK 0.13 implementation - if( stage == MORPHOLOGY_UNRAVELED ) - stage = MORPHOLOGY_RAW; - - H5::DataSet dataset; - try - { - dataset = _file.openDataSet( "/" + _g_root + "/" + _g_structure + - "/" + toString( stage )); - } - catch( ... ) - { - LBERROR << "Could not open sections dataset for morphology file " - << _file.getFileName() << " repair stage " - << toString( stage ) << std::endl; - return Vector2isPtr( new Vector2is ); - } - - hsize_t dims[2]; - const H5::DataSpace& dspace = dataset.getSpace(); - if( dspace.getSimpleExtentNdims() != 2 || - dspace.getSimpleExtentDims( dims ) < 0 || - dims[1] != _structureV2Columns ) - { - LBTHROW( std::runtime_error( "Reading morphology file '" + - _file.getFileName() + "': bad number of dimensions in" - " 'structure' dataspace")); - } + const hsize_t readCount[2] = { _sectionsDims[0], 1 }; + const hsize_t readOffset[2] = { 0, 1 }; + H5::DataSpace dspace = _sections.getSpace(); + dspace.selectHyperslab( H5S_SELECT_XOR, readCount, readOffset ); - Vector2isPtr data( new Vector2is( dims[0] )); - dataset.read( data->data(), H5::PredType::NATIVE_INT ); + Vector2isPtr data( new Vector2is( _sectionsDims[0] )); + const hsize_t mdim[2] = { _sectionsDims[0], 2 }; + const H5::DataSpace mspace( 2, mdim ); + _sections.read( data->data(), H5::PredType::NATIVE_INT, mspace, + dspace ); return data; } @@ -281,48 +291,36 @@ SectionTypesPtr MorphologyHDF5::readSectionTypes() const { lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); - if( _version == MORPHOLOGY_VERSION_H5_1 ) + if( _version == MORPHOLOGY_VERSION_H5_2 ) { - const H5::DataSet& dataset = _file.openDataSet( _d_structure ); + const H5::DataSet& dataset = _file.openDataSet( "/" + _g_root + + "/" + _g_structure + "/" + _d_type ); hsize_t dims[2]; const H5::DataSpace& dspace = dataset.getSpace(); if( dspace.getSimpleExtentNdims() != 2 || - dspace.getSimpleExtentDims( dims ) < 0 || - dims[1] != _structureV1Columns ) + dspace.getSimpleExtentDims( dims ) < 0 || dims[1] != 1 ) { LBTHROW( std::runtime_error( "Reading morphology file '" + _file.getFileName() + "': bad number of dimensions in" - " 'structure' dataspace")); + " 'sectiontype' dataspace")); } - const hsize_t readCount[2] = { dims[0], 1 }; - const hsize_t readOffset[2] = { 0, 1 }; - dspace.selectHyperslab( H5S_SELECT_SET, readCount, readOffset ); - SectionTypesPtr data( new SectionTypes( dims[0] )); - const hsize_t mdim = dims[0]; - const H5::DataSpace mspace( 1, &mdim ); - dataset.read( data->data(), H5::PredType::NATIVE_INT, mspace, - dspace ); + dataset.read( data->data(), H5::PredType::NATIVE_INT ); return data; } - const H5::DataSet& dataset = _file.openDataSet( "/" + _g_root + - "/" + _g_structure + "/" + _d_type ); - - hsize_t dims[2]; - const H5::DataSpace& dspace = dataset.getSpace(); - if( dspace.getSimpleExtentNdims() != 2 || - dspace.getSimpleExtentDims( dims ) < 0 || dims[1] != 1 ) - { - LBTHROW( std::runtime_error( "Reading morphology file '" + - _file.getFileName() + "': bad number of dimensions in" - " 'sectiontype' dataspace")); - } + const hsize_t readCount[2] = { _sectionsDims[0], 1 }; + const hsize_t readOffset[2] = { 0, 1 }; + H5::DataSpace dspace = _sections.getSpace(); + dspace.selectHyperslab( H5S_SELECT_SET, readCount, readOffset ); - SectionTypesPtr data( new SectionTypes( dims[0] )); - dataset.read( data->data(), H5::PredType::NATIVE_INT ); + SectionTypesPtr data( new SectionTypes( _sectionsDims[0] )); + const hsize_t mdim = _sectionsDims[0]; + const H5::DataSpace mspace( 1, &mdim ); + _sections.read( data->data(), H5::PredType::NATIVE_INT, mspace, + dspace ); return data; } @@ -352,6 +350,41 @@ Vector2isPtr MorphologyHDF5::readApicals() const return data; } +floatsPtr MorphologyHDF5::readPerimeters() const +{ + if( _version != MORPHOLOGY_VERSION_H5_1_1 ) + return floatsPtr( new floats( )); + + lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); + + try + { + detail::SilenceHDF5 silence; + H5::DataSet dataset = _file.openDataSet( _d_perimeters ); + + hsize_t dims; + const H5::DataSpace& dspace = dataset.getSpace(); + if( dspace.getSimpleExtentNdims() != 1 || + dspace.getSimpleExtentDims( &dims ) < 0 ) + { + LBTHROW( std::runtime_error( "Reading morphology file '" + + _file.getFileName() + "': bad number of dimensions in" + " 'perimeters' dataspace")); + } + + floatsPtr data( new floats( dims )); + dataset.read( data->data(), H5::PredType::NATIVE_FLOAT ); + return data; + } + catch( ... ) + { + if( getCellFamily() == FAMILY_GLIA ) + LBTHROW( std::runtime_error( "No empty perimeters allowed for glia " + "morphology" )); + return floatsPtr( new floats( )); + } +} + MorphologyVersion MorphologyHDF5::getVersion() const { return _version; @@ -364,20 +397,53 @@ void MorphologyHDF5::writePoints( const Vector4fs& points, ASSERT_WRITE; - H5::FloatType pointsDT( H5::PredType::NATIVE_DOUBLE ); - pointsDT.setOrder( H5T_ORDER_LE ); hsize_t dim[2] = { points.size(), 4 }; H5::DataSpace pointsDS( 2, dim ); H5::DataSet dataset; - if( _version == MORPHOLOGY_VERSION_H5_1 ) - dataset = _file.createDataSet( _d_points, pointsDT, pointsDS ); - else + H5::FloatType pointsDT( H5::PredType::NATIVE_DOUBLE ); + pointsDT.setOrder( H5T_ORDER_LE ); + if( _version == MORPHOLOGY_VERSION_H5_2 ) { H5::Group root = _file.openGroup( _g_root ); H5::Group group = root.createGroup( toString( stage )); + try + { + detail::SilenceHDF5 silence; + group.openDataSet( _d_points ); + LBTHROW( std::runtime_error( "Points were already written" )); + } + catch( const H5::Exception& ) {} dataset = group.createDataSet( _d_points, pointsDT, pointsDS ); } + else + { + try + { + detail::SilenceHDF5 silence; + _file.openDataSet( _d_points ); + LBTHROW( std::runtime_error( "Points were already written" )); + } + catch( const H5::Exception& ) {} + + try + { + detail::SilenceHDF5 silence; + H5::DataSet perimeters = _file.openDataSet( _d_perimeters ); + hsize_t dims[2]; + perimeters.getSpace().getSimpleExtentDims( dims ); + if( dims[0] != points.size( )) + { + std::stringstream msg; + msg << "Number of points does not match number of perimeters, " + << dims[0] << " != " << points.size() << std::endl; + LBTHROW( std::runtime_error( msg.str( ))); + } + } + catch( const H5::Exception& ) {} + + dataset = _file.createDataSet( _d_points, pointsDT, pointsDS ); + } dataset.write( points.data(), H5::PredType::NATIVE_FLOAT ); } @@ -389,37 +455,30 @@ void MorphologyHDF5::writeSections( const Vector2is& sections, ASSERT_WRITE; - H5::FloatType structureDT( H5::PredType::NATIVE_INT ); - structureDT.setOrder( H5T_ORDER_LE ); - - if( _version == MORPHOLOGY_VERSION_H5_1 ) + if( _version == MORPHOLOGY_VERSION_H5_2 ) { - hsize_t dim[2] = { sections.size(), 3 }; + hsize_t dim[2] = { sections.size(), 2 }; H5::DataSpace structureDS( 2, dim ); - H5::DataSet dataset = _file.getNumObjs() == 2 ? - _file.openDataSet( _g_structure ) : - _file.createDataSet( _d_structure, structureDT, structureDS ); - - const H5::DataSpace& dspace = dataset.getSpace(); - const hsize_t count[2] = { sections.size(), 1 }; - const hsize_t offset[2] = { 0, 1 }; - dspace.selectHyperslab( H5S_SELECT_XOR, count, offset ); - - const hsize_t mdim[2] = { sections.size(), 2 }; - const H5::DataSpace mspace( 2, mdim ); - dataset.write( sections.data(), H5::PredType::NATIVE_INT, mspace, - dspace ); + H5::FloatType structureDT( H5::PredType::NATIVE_INT ); + structureDT.setOrder( H5T_ORDER_LE ); + H5::Group root = _file.openGroup( _g_root + "/" + _g_structure ); + H5::DataSet dataset = root.createDataSet( toString( stage ), + structureDT, structureDS ); + dataset.write( sections.data(), H5::PredType::NATIVE_INT ); return; } - hsize_t dim[2] = { sections.size(), 2 }; - H5::DataSpace structureDS( 2, dim ); - - H5::Group root = _file.openGroup( _g_root + "/" + _g_structure ); - H5::DataSet dataset = root.createDataSet( toString( stage ), - structureDT, structureDS ); - dataset.write( sections.data(), H5::PredType::NATIVE_INT ); + H5::DataSet dataset = _getStructureDataSet( sections.size( )); + const H5::DataSpace& dspace = dataset.getSpace(); + const hsize_t count[2] = { sections.size(), 1 }; + const hsize_t offset[2] = { 0, 1 }; + dspace.selectHyperslab( H5S_SELECT_XOR, count, offset ); + + const hsize_t mdim[2] = { sections.size(), 2 }; + const H5::DataSpace mspace( 2, mdim ); + dataset.write( sections.data(), H5::PredType::NATIVE_INT, mspace, + dspace ); } void MorphologyHDF5::writeSectionTypes( const SectionTypes& types ) @@ -428,52 +487,44 @@ void MorphologyHDF5::writeSectionTypes( const SectionTypes& types ) ASSERT_WRITE; - H5::FloatType structureDT( H5::PredType::NATIVE_INT ); - structureDT.setOrder( H5T_ORDER_LE ); - - if( _version == MORPHOLOGY_VERSION_H5_1 ) + if( _version == MORPHOLOGY_VERSION_H5_2 ) { - hsize_t dim[2] = { types.size(), 3 }; + hsize_t dim[2] = { types.size(), 1 }; H5::DataSpace structureDS( 2, dim ); - H5::DataSet dataset = _file.getNumObjs() == 2 ? - _file.openDataSet( _g_structure ) : - _file.createDataSet( _d_structure, structureDT, structureDS ); + H5::FloatType structureDT( H5::PredType::NATIVE_INT ); + structureDT.setOrder( H5T_ORDER_LE ); + H5::Group root = _file.openGroup( _g_root + "/" + _g_structure ); + H5::DataSet dataset = root.createDataSet( _d_type, structureDT, + structureDS ); + dataset.write( types.data(), H5::PredType::NATIVE_INT ); - const H5::DataSpace& dspace = dataset.getSpace(); - const hsize_t count[2] = { types.size(), 1 }; - const hsize_t offset[2] = { 0, 1 }; - dspace.selectHyperslab( H5S_SELECT_SET, count, offset ); - - const hsize_t mdim = types.size(); - const H5::DataSpace mspace( 1, &mdim ); - dataset.write( types.data(), H5::PredType::NATIVE_INT, mspace, - dspace ); return; } - hsize_t dim[2] = { types.size(), 1 }; - H5::DataSpace structureDS( 2, dim ); + H5::DataSet dataset = _getStructureDataSet( types.size( )); + const H5::DataSpace& dspace = dataset.getSpace(); + const hsize_t count[2] = { types.size(), 1 }; + const hsize_t offset[2] = { 0, 1 }; + dspace.selectHyperslab( H5S_SELECT_SET, count, offset ); - H5::Group root = _file.openGroup( _g_root + "/" + _g_structure ); - H5::DataSet dataset = root.createDataSet( _d_type, structureDT, - structureDS ); - dataset.write( types.data(), H5::PredType::NATIVE_INT ); + const hsize_t mdim = types.size(); + const H5::DataSpace mspace( 1, &mdim ); + dataset.write( types.data(), H5::PredType::NATIVE_INT, mspace, dspace ); } void MorphologyHDF5::writeApicals( const Vector2is& apicals ) { - lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); - ASSERT_WRITE; - if( _version == MORPHOLOGY_VERSION_H5_1 ) - LBTHROW( std::runtime_error( "No apical sections for version 1 " - "morphology files." )); + if( _version != MORPHOLOGY_VERSION_H5_2 ) + LBTHROW( std::runtime_error( "Need version 2 to write apicals" )); if( apicals.empty( )) return; + lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); + const H5::Group& root = _file.openGroup( _g_root ); hsize_t dims = 2 * apicals.size(); H5::DataSpace apicalDS( 1, &dims ); @@ -482,6 +533,57 @@ void MorphologyHDF5::writeApicals( const Vector2is& apicals ) attr.write( H5::PredType::NATIVE_INT, apicals.data( )); } +void MorphologyHDF5::writePerimeters( const floats& perimeters ) +{ + ASSERT_WRITE; + + if( _version != MORPHOLOGY_VERSION_H5_1_1 ) + LBTHROW( std::runtime_error( "Need version 1.1 to write perimeters" )); + + if( perimeters.empty( )) + { + if( getCellFamily() == FAMILY_GLIA ) + LBTHROW( std::runtime_error( "No empty perimeters allowed for glia " + "morphology" )); + return; + } + + lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); + + try + { + detail::SilenceHDF5 silence; + _file.openDataSet( _d_perimeters ); + LBTHROW( std::runtime_error( "Perimeters were already written" )); + } + catch( const H5::Exception& ) {} + + const hsize_t dim = perimeters.size(); + + try + { + detail::SilenceHDF5 silence; + H5::DataSet points = _file.openDataSet( _d_points ); + hsize_t dims[2]; + points.getSpace().getSimpleExtentDims( dims ); + if( dims[0] != perimeters.size( )) + { + std::stringstream msg; + msg << "Number of perimeters does not match number of points, " + << dims[0] << " != " << perimeters.size() << std::endl; + LBTHROW( std::runtime_error( msg.str( ))); + } + } + catch( const H5::Exception& ) {} + + H5::DataSpace perimeterDS( 1, &dim ); + + H5::DataSet dataset = _file.createDataSet( _d_perimeters, + H5::PredType::NATIVE_FLOAT, + perimeterDS ); + dataset.write( perimeters.data(), H5::PredType::NATIVE_FLOAT ); +} + void MorphologyHDF5::flush() { lunchbox::ScopedWrite mutex( detail::_hdf5Lock ); @@ -492,9 +594,18 @@ void MorphologyHDF5::flush() void MorphologyHDF5::_checkVersion( const std::string& source ) { - if( _isV1( )) + if( _readV11Metadata( )) + return; + + try + { + _resolveV1(); _version = MORPHOLOGY_VERSION_H5_1; - else if( _isV2( )) + return; + } + catch( ... ) {} + + if( _readV2Metadata( )) _version = MORPHOLOGY_VERSION_H5_2; else LBTHROW( std::runtime_error( "Unknown morphology file format for " @@ -503,7 +614,7 @@ void MorphologyHDF5::_checkVersion( const std::string& source ) void MorphologyHDF5::_selectRepairStage() { - if( _version == MORPHOLOGY_VERSION_H5_1 ) + if( _version != MORPHOLOGY_VERSION_H5_2 ) return; MorphologyRepairStage stages[3] = { MORPHOLOGY_REPAIRED, @@ -523,12 +634,98 @@ void MorphologyHDF5::_selectRepairStage() } } -bool MorphologyHDF5:: _isV1() const +void MorphologyHDF5::_resolveV1() +{ + _points = _file.openDataSet( "/" + _d_points ); + const H5::DataSpace pointsSpace = _points.getSpace(); + if( pointsSpace.getSimpleExtentNdims() != 2 || + pointsSpace.getSimpleExtentDims( _pointsDims ) < 0 || + _pointsDims[1] != _pointColumns ) + { + LBTHROW( std::runtime_error( "Opening morphology file '" + + _file.getFileName() + "': bad number of dimensions in" + " 'points' dataspace")); + } + + _sections = _file.openDataSet( _d_structure ); + const H5::DataSpace sectionsSpace = _sections.getSpace(); + if( sectionsSpace.getSimpleExtentNdims() != 2 || + sectionsSpace.getSimpleExtentDims(_sectionsDims ) < 0 || + _sectionsDims[1] != _structureV1Columns ) + { + LBTHROW( std::runtime_error( "Opening morphology file '" + + _file.getFileName() + "': bad number of dimensions in" + " 'structure' dataspace")); + } +} + +void MorphologyHDF5::_writeV11Metadata( const MorphologyInitData& initData ) +{ + ASSERT_WRITE; + + H5::Group metadata = _file.createGroup( _g_metadata ); + + H5::EnumType familyEnum( H5::PredType::STD_U32LE ); + uint32_t enumValue = brion::FAMILY_NEURON; + familyEnum.insert( "NEURON", &enumValue ); + enumValue = brion::FAMILY_GLIA; + familyEnum.insert( "GLIA", &enumValue ); + familyEnum.commit( metadata, _e_family ); + + const hsize_t dim = 1; + H5::DataSpace familyDS( 1, &dim ); + H5::Attribute familyAttr = metadata.createAttribute( _a_family, familyEnum, + familyDS ); + _family = initData.getFamily(); + familyAttr.write( familyEnum, &_family ); + + const std::string creator = "Brion"; + detail::addStringAttribute( metadata, _a_creator, creator ); + + detail::addStringAttribute( metadata, _a_software_version, + Version::getString( )); + + const time_t now = ::time(0); +#ifdef _WIN32 + char* gmtString = ::ctime( &now ); +#else + char gmtString[32]; + ::ctime_r( &now, gmtString ); +#endif + + std::string creation_time = gmtString; + creation_time = creation_time.substr( 0, creation_time.size() -1 ); // ctime_r ends with \n + detail::addStringAttribute( metadata, _a_creation_time, creation_time ); + + hsize_t dims = 2; + H5::DataSpace versionDS( 1, &dims ); + H5::Attribute versionAttr = metadata.createAttribute( _a_version, + H5::PredType::STD_U32LE, versionDS ); + const uint32_t version[2] = { 1, 1 }; + versionAttr.write( H5::PredType::STD_U32LE, &version[0] ); +} + +bool MorphologyHDF5::_readV11Metadata() { try { detail::SilenceHDF5 silence; - _file.openDataSet( _d_structure ); + const H5::Group& metadata = _file.openGroup( _g_metadata ); + const H5::Attribute& attr = metadata.openAttribute( _a_version ); + + uint32_t version[2]; + attr.read( H5::PredType::STD_U32LE, &version[0] ); + if( version[0] != 1 || version[1] != 1 ) + return false; + + _version = MORPHOLOGY_VERSION_H5_1_1; + + const H5::Attribute& familyAttr = metadata.openAttribute( _a_family ); + H5::EnumType familyEnum = metadata.openEnumType( _e_family ); + + familyAttr.read( familyEnum, &_family ); + + _resolveV1(); return true; } catch( const H5::Exception& ) @@ -537,7 +734,7 @@ bool MorphologyHDF5:: _isV1() const } } -bool MorphologyHDF5:: _isV2() const +bool MorphologyHDF5:: _readV2Metadata() const { try { @@ -565,12 +762,7 @@ bool MorphologyHDF5:: _isV2() const } } -void MorphologyHDF5::_writeV1Header() -{ - ASSERT_WRITE; -} - -void MorphologyHDF5::_writeV2Header() +void MorphologyHDF5::_writeV2Metadata() { ASSERT_WRITE; @@ -596,5 +788,23 @@ void MorphologyHDF5::_writeV2Header() versionAttr.write( H5::PredType::NATIVE_INT, &version ); } +H5::DataSet MorphologyHDF5::_getStructureDataSet( size_t nSections ) +{ + H5::DataSet dataset; + try + { + detail::SilenceHDF5 silence; + return _file.openDataSet( _d_structure ); + } + catch( const H5::Exception& ) + { + hsize_t dim[2] = { nSections, 3 }; + H5::DataSpace structureDS( 2, dim ); + H5::FloatType structureDT( H5::PredType::NATIVE_INT ); + structureDT.setOrder( H5T_ORDER_LE ); + return _file.createDataSet( _d_structure, structureDT, structureDS ); + } +} + } } diff --git a/brion/plugin/morphologyHDF5.h b/brion/plugin/morphologyHDF5.h index ffda85a..6bf0706 100644 --- a/brion/plugin/morphologyHDF5.h +++ b/brion/plugin/morphologyHDF5.h @@ -1,4 +1,4 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Daniel Nachbaur * Juan Hernando * @@ -41,6 +41,8 @@ class MorphologyHDF5 : public MorphologyPlugin /** Check if this plugin can handle the given uri. */ static bool handles( const MorphologyInitData& initData ); + CellFamily getCellFamily() const final; + Vector4fsPtr readPoints( MorphologyRepairStage stage ) const final; Vector2isPtr readSections( MorphologyRepairStage stage ) const final; @@ -49,6 +51,8 @@ class MorphologyHDF5 : public MorphologyPlugin Vector2isPtr readApicals() const final; + floatsPtr readPerimeters() const final; + MorphologyVersion getVersion() const final; void writePoints( const Vector4fs& points, @@ -61,25 +65,36 @@ class MorphologyHDF5 : public MorphologyPlugin void writeApicals( const Vector2is& apicals ) final; + void writePerimeters( const floats& perimeters ) final; + void flush() final; private: H5::H5File _file; + + H5::DataSet _points; + hsize_t _pointsDims[2]; + + H5::DataSet _sections; + hsize_t _sectionsDims[2]; + MorphologyVersion _version; MorphologyRepairStage _stage; + CellFamily _family; bool _write; void _checkVersion( const std::string& source ); - void _selectRepairStage(); - bool _isV1() const; + void _resolveV1(); - bool _isV2() const; + void _writeV11Metadata( const MorphologyInitData& initData ); + bool _readV11Metadata(); - void _writeV1Header(); + bool _readV2Metadata() const; + void _writeV2Metadata(); - void _writeV2Header(); + H5::DataSet _getStructureDataSet( size_t nSections ); }; } diff --git a/brion/plugin/morphologySWC.cpp b/brion/plugin/morphologySWC.cpp index e4844dc..14724db 100644 --- a/brion/plugin/morphologySWC.cpp +++ b/brion/plugin/morphologySWC.cpp @@ -1,5 +1,6 @@ /* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Juan Hernando + * Daniel Nachbaur * * This file is part of Brion * @@ -185,6 +186,11 @@ bool MorphologySWC::handles( const MorphologyInitData& initData ) return path.substr( pos ) == ".swc"; } +CellFamily MorphologySWC::getCellFamily() const +{ + return FAMILY_NEURON; +} + Vector4fsPtr MorphologySWC::readPoints( MorphologyRepairStage ) const { return _points; @@ -203,7 +209,12 @@ SectionTypesPtr MorphologySWC::readSectionTypes() const Vector2isPtr MorphologySWC::readApicals() const { // Can these points be inferred? Should it be done at all? - return Vector2isPtr(new Vector2is()); + return Vector2isPtr( new Vector2is( )); +} + +floatsPtr MorphologySWC::readPerimeters() const +{ + return floatsPtr( new floats( )); } MorphologyVersion MorphologySWC::getVersion() const @@ -233,6 +244,11 @@ void MorphologySWC::writeApicals( const Vector2is& ) LBUNIMPLEMENTED } +void MorphologySWC::writePerimeters( const floats& ) +{ + LBUNIMPLEMENTED +} + void MorphologySWC::flush() { LBUNIMPLEMENTED diff --git a/brion/plugin/morphologySWC.h b/brion/plugin/morphologySWC.h index 5198089..3a8a7a3 100644 --- a/brion/plugin/morphologySWC.h +++ b/brion/plugin/morphologySWC.h @@ -1,5 +1,6 @@ -/* Copyright (c) 2013-2015, EPFL/Blue Brain Project +/* Copyright (c) 2013-2016, EPFL/Blue Brain Project * Juan Hernando + * Daniel Nachbaur * * This file is part of Brion * @@ -35,6 +36,8 @@ class MorphologySWC : public MorphologyPlugin /** Check if this plugin can handle the given uri. */ static bool handles( const MorphologyInitData& initData ); + CellFamily getCellFamily() const final; + Vector4fsPtr readPoints( MorphologyRepairStage stage ) const final; Vector2isPtr readSections( MorphologyRepairStage stage ) const final; @@ -43,6 +46,8 @@ class MorphologySWC : public MorphologyPlugin Vector2isPtr readApicals() const final; + floatsPtr readPerimeters() const final; + MorphologyVersion getVersion() const final; void writePoints( const Vector4fs& points, @@ -55,6 +60,8 @@ class MorphologySWC : public MorphologyPlugin void writeApicals( const Vector2is& apicals ) final; + void writePerimeters( const floats& perimeters ) final; + void flush() final; private: diff --git a/doc/Changelog.md b/doc/Changelog.md index 2630764..b06be3e 100644 --- a/doc/Changelog.md +++ b/doc/Changelog.md @@ -3,6 +3,8 @@ Changelog {#Changelog} # git master +* [75](https://github.com/BlueBrain/Brion/pull/75): + Implement morphology version 1.1 specification * [74](https://github.com/BlueBrain/Brion/pull/74): Remove deprecated enums and functions: - `CompartmentReport( const std::string&, const GIDSet& )` and diff --git a/doc/feature/morphology11.md b/doc/feature/morphology11.md index 70f0a0c..1be488e 100644 --- a/doc/feature/morphology11.md +++ b/doc/feature/morphology11.md @@ -64,7 +64,7 @@ None { ... /** - * Open the given morphology file to write a version 1.1 morphology. + * Open the given morphology file for write access. * * @param file filename of the output morphology file * @param family the cell family that this morphology represents. @@ -76,8 +76,8 @@ None /** * @return perimeters of the cross sections for each point of the - * morphology in micrometers. Can be empty for version 1.1 - * morphology of FAMILY_NEURON. + * morphology in micrometers. + * @throw std::runtime_error if empty for FAMILY_GLIA */ floatsPtr readPerimeters() const; @@ -85,17 +85,23 @@ None * Write perimeters of morphology. * * @param perimeters perimeters of the cross sections for each point of - * the morphology in micrometers. + * the morphology in micrometers + * @throw std::runtime_error if object not created with write ctor + * @throw std::runtime_error if not supported by implementation */ void writePerimeters( const floats& perimeters ); ## File format -* new root 'cell_family' with datatype H5::NATIVE_INT * new root dataset 'perimeters' with dimension 1 and datatype H5::NATIVE_FLOAT -* new root attributes 'creator', 'software_version', 'creation_time' with - datatype H5::C_S1 (string) -* new root attribute 'version' with datatype H5::STD_U32LE and dimensionality 2 +* new root group 'metadata' +* new enum 'cell_family_enum' with datatype H5::NATIVE_INT and values + NEURON=0 and GLIA=1 inside 'metadata' +* new attribute 'cell_family' with datatype 'cell_family_enum' inside 'metadata' +* new attributes 'creator', 'software_version', 'creation_time' with datatype + H5::C_S1 (string) inside 'metadata' +* new attribute 'version' with datatype H5::STD_U32LE and dimensionality 2 + inside 'metadata' ## Examples diff --git a/tests/morphology.cpp b/tests/morphology.cpp index ff66e8c..2f35259 100644 --- a/tests/morphology.cpp +++ b/tests/morphology.cpp @@ -64,6 +64,15 @@ const int DENDRITE = brion::SECTION_DENDRITE; const int APICAL_DENDRITE = brion::SECTION_APICAL_DENDRITE; } +template< typename T > +void checkCloseArrays( const std::vector< T >& array1, + const std::vector< T >& array2 ) +{ + BOOST_CHECK_EQUAL( array1.size(), array2.size() ); + for( size_t i = 0; i != std::min( array1.size(), array2.size( )); ++i ) + BOOST_CHECK_CLOSE( array1[i], array2[i], 2e-5f); +} + BOOST_AUTO_TEST_CASE( invalid_open ) { BOOST_CHECK_THROW( brion::Morphology( "/bla" ), std::runtime_error ); @@ -127,6 +136,7 @@ BOOST_AUTO_TEST_CASE( h5_read_v1 ) path /= "local/morphologies/01.07.08/h5/R-C010306G.h5"; const brion::Morphology morphology( path.string( )); + BOOST_CHECK_EQUAL( morphology.getCellFamily(), brion::FAMILY_NEURON ); const brion::Vector4fsPtr points = morphology.readPoints( brion::MORPHOLOGY_RAW ); @@ -148,6 +158,8 @@ BOOST_AUTO_TEST_CASE( h5_read_v1 ) BOOST_CHECK_EQUAL( types->size(), 138 ); BOOST_CHECK_EQUAL( (*types)[0], 1 ); BOOST_CHECK_EQUAL( (*types)[5], 2 ); + + BOOST_CHECK( morphology.readPerimeters()->empty( )); } BOOST_AUTO_TEST_CASE( h5_write_v1 ) @@ -167,6 +179,10 @@ BOOST_AUTO_TEST_CASE( h5_write_v1 ) a.writePoints( *points, brion::MORPHOLOGY_UNDEFINED ); a.writeSections( *sections, brion::MORPHOLOGY_UNDEFINED ); a.writeSectionTypes( *types ); + BOOST_CHECK_THROW( a.writeApicals( brion::Vector2is( )), + std::runtime_error ); + BOOST_CHECK_THROW( a.writePerimeters( brion::floats( )), + std::runtime_error ); } const brion::Morphology source2( "testv1.h5" ); @@ -180,12 +196,101 @@ BOOST_AUTO_TEST_CASE( h5_write_v1 ) BOOST_CHECK( *types == *types2 ); } +BOOST_AUTO_TEST_CASE( h5_write_v11_glia ) +{ + const std::string file( "glia.h5" ); + boost::filesystem::remove( file ); + + boost::filesystem::path path( BBP_TESTDATA ); + path /= "local/morphologies/01.07.08/h5/R-C010306G.h5"; + + const brion::Morphology morphology( path.string( )); + brion::floats perimeters; + brion::Vector4fsPtr points = morphology.readPoints(); + perimeters.reserve( points->size( )); + for( size_t i = 0; i < points->size(); ++i ) + perimeters.push_back( (*points)[i].w() * 4.f ); + + brion::Morphology glia( file, brion::FAMILY_GLIA ); + glia.writePoints( *points ); + glia.writeSections( *morphology.readSections( )); + glia.writeSectionTypes( *morphology.readSectionTypes( )); + glia.writePerimeters( perimeters ); + BOOST_CHECK_THROW( glia.writeApicals( brion::Vector2is( )), + std::runtime_error ); + + const brion::Morphology gliaRead( file ); + BOOST_CHECK_EQUAL( gliaRead.getCellFamily(), brion::FAMILY_GLIA ); + checkCloseArrays( *gliaRead.readPerimeters(), perimeters ); + + boost::filesystem::remove( file ); +} + +BOOST_AUTO_TEST_CASE( h5_write_invalid_glia ) +{ + const std::string file( "glia.h5" ); + boost::filesystem::remove( file ); + brion::Morphology glia( file, brion::FAMILY_GLIA ); + brion::floats perimeters; + BOOST_CHECK_THROW( glia.writePerimeters( perimeters ), std::runtime_error ); + boost::filesystem::remove( file ); +} + +BOOST_AUTO_TEST_CASE( h5_write_v11_neuron ) +{ + const std::string file( "neuron.h5" ); + boost::filesystem::remove( file ); + + boost::filesystem::path path( BBP_TESTDATA ); + path /= "local/morphologies/01.07.08/h5/R-C010306G.h5"; + + const brion::Morphology morphology( path.string( )); + + brion::Morphology neuron( file, brion::FAMILY_NEURON ); + neuron.writePoints( *morphology.readPoints( )); + neuron.writeSections( *morphology.readSections( )); + neuron.writeSectionTypes( *morphology.readSectionTypes( )); + BOOST_CHECK_THROW( neuron.writeApicals( brion::Vector2is( )), + std::runtime_error ); + + const brion::Morphology neuronRead( file ); + BOOST_CHECK_EQUAL( neuronRead.getCellFamily(), brion::FAMILY_NEURON ); + BOOST_CHECK( neuronRead.readPerimeters()->empty( )); + + boost::filesystem::remove( file ); +} + +BOOST_AUTO_TEST_CASE( h5_write_invalid_neuron ) +{ + const std::string file( "neuron.h5" ); + boost::filesystem::remove( file ); + + { + brion::Morphology neuron( file, brion::FAMILY_NEURON ); + neuron.writePoints( brion::Vector4fs( 5 )); + BOOST_CHECK_THROW( neuron.writePoints( brion::Vector4fs( 1 )), std::runtime_error ); + BOOST_CHECK_THROW( neuron.writePerimeters( brion::floats( 4 )), std::runtime_error ); + neuron.writePerimeters( brion::floats( 5 )); + boost::filesystem::remove( file ); + } + + { + brion::Morphology neuron( file, brion::FAMILY_NEURON ); + neuron.writePerimeters( brion::floats( 5 )); + BOOST_CHECK_THROW( neuron.writePerimeters( brion::floats( 4 )), std::runtime_error ); + BOOST_CHECK_THROW( neuron.writePoints( brion::Vector4fs( 3 )), std::runtime_error ); + neuron.writePoints( brion::Vector4fs( 5 )); + boost::filesystem::remove( file ); + } +} + BOOST_AUTO_TEST_CASE( h5_read_v2 ) { boost::filesystem::path path( BBP_TESTDATA ); path /= "local/morphologies/14.07.10_repaired/v2/C010398B-P2.h5"; const brion::Morphology morphology( path.string( )); + BOOST_CHECK_EQUAL( morphology.getCellFamily(), brion::FAMILY_NEURON ); brion::Vector4fsPtr points = morphology.readPoints( brion::MORPHOLOGY_REPAIRED ); @@ -212,6 +317,8 @@ BOOST_AUTO_TEST_CASE( h5_read_v2 ) BOOST_CHECK_EQUAL( apicals->size(), 1 ); BOOST_CHECK_EQUAL( (*apicals)[0].x(), 67 ); BOOST_CHECK_EQUAL( (*apicals)[0].y(), 76 ); + + BOOST_CHECK( morphology.readPerimeters()->empty( )); } BOOST_AUTO_TEST_CASE( h5_write_v2 ) @@ -234,6 +341,8 @@ BOOST_AUTO_TEST_CASE( h5_write_v2 ) a.writeSections( *sections, brion::MORPHOLOGY_REPAIRED ); a.writeSectionTypes( *types ); a.writeApicals( *apicals ); + BOOST_CHECK_THROW( a.writePerimeters( brion::floats( )), + std::runtime_error ); } const brion::Morphology source2( "testv2.h5" ); @@ -337,15 +446,6 @@ void checkCloseArrays( const std::vector< T >& array, va_end( args ); } -template< typename T > -void checkCloseArrays( const std::vector< T >& array1, - const std::vector< T >& array2 ) -{ - BOOST_CHECK_EQUAL( array1.size(), array2.size() ); - for( size_t i = 0; i != std::min( array1.size(), array2.size( )); ++i ) - BOOST_CHECK_CLOSE( array1[i], array2[i], 2e-5f); -} - template< typename T, long unsigned int M > void checkCloseArrays( const std::vector< vmml::vector< M, T > >& array1, const std::vector< vmml::vector< M, T > >& array2 ) @@ -553,6 +653,8 @@ BOOST_AUTO_TEST_CASE( swc_neuron ) brion::Morphology neuron( path.string( )); BOOST_CHECK_EQUAL( neuron.readPoints( brion::MORPHOLOGY_REPAIRED )->size(), 927 ); + BOOST_CHECK_EQUAL( neuron.getCellFamily(), brion::FAMILY_NEURON ); + BOOST_CHECK( neuron.readPerimeters()->empty( )); } namespace