From 990f4992062489776a0044ae9d7768d09c14d754 Mon Sep 17 00:00:00 2001 From: Christian Blau Date: Wed, 17 Apr 2019 15:50:35 +0200 Subject: [PATCH] mrc density map reading/writing Add a class that combines header and mrc data handling for float values to enable reading and writing of mrc files. refs #2282 Change-Id: Idd84ab2628a52a2fcafb6121ad66f62e88235c74 --- src/gromacs/fileio/mrcdensitymap.cpp | 179 +++++++++++++++++++++++ src/gromacs/fileio/mrcdensitymap.h | 103 +++++++++++++ src/gromacs/fileio/mrcdensitymapheader.cpp | 67 +++++++++ src/gromacs/fileio/mrcdensitymapheader.h | 7 + src/gromacs/fileio/tests/CMakeLists.txt | 2 + src/gromacs/fileio/tests/mrcdensitymap.cpp | 90 ++++++++++++ src/gromacs/fileio/tests/mrcdensitymapheader.cpp | 75 ++++++++++ 7 files changed, 523 insertions(+) create mode 100644 src/gromacs/fileio/mrcdensitymap.cpp create mode 100644 src/gromacs/fileio/mrcdensitymap.h create mode 100644 src/gromacs/fileio/mrcdensitymapheader.cpp create mode 100644 src/gromacs/fileio/tests/mrcdensitymap.cpp create mode 100644 src/gromacs/fileio/tests/mrcdensitymapheader.cpp diff --git a/src/gromacs/fileio/mrcdensitymap.cpp b/src/gromacs/fileio/mrcdensitymap.cpp new file mode 100644 index 0000000000..9c8385c7c2 --- /dev/null +++ b/src/gromacs/fileio/mrcdensitymap.cpp @@ -0,0 +1,179 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2019, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * + * \brief + * Implements mrc/ccp4-file format handling. + * + * \author Christian Blau + * + * \ingroup module_fileio + */ +#include "gmxpre.h" + +#include "mrcdensitymap.h" + +#include +#include + +#include "gromacs/fileio/mrcdensitymapheader.h" +#include "gromacs/utility/exceptions.h" +#include "gromacs/utility/iserializer.h" + +#include "mrcserializer.h" + +namespace gmx +{ + +/******************************************************************** + * MrcDensityMapOfFloatReader::Impl + */ + +/*! \internal \brief + * Private implementation class for MrcDensityMapOfFloatReader. + */ +class MrcDensityMapOfFloatReader::Impl +{ + public: + //! Build the map reader from a serializer. + explicit Impl(ISerializer *serializer); + ~Impl() {} + //! The header of the read mrc file + MrcDensityMapHeader header_; + //! The data of the mrc file + std::vector data_; +}; + +MrcDensityMapOfFloatReader::Impl::Impl(ISerializer *serializer) +{ + if (!serializer->reading()) + { + GMX_THROW(InternalError("Cannot use writing serializer to read.")); + } + + header_ = deserializeMrcDensityMapHeader(serializer); + const auto dataSize = numberOfExpectedDataItems(header_); + data_.resize(dataSize); + for (auto &value : data_) + { + serializer->doFloat(&value); + } +} + +/******************************************************************** + * MrcDensityMapOfFloatReader + */ + +MrcDensityMapOfFloatReader::MrcDensityMapOfFloatReader(ISerializer *serializer) : impl_(new Impl(serializer)) +{ +} + +ArrayRef MrcDensityMapOfFloatReader::data() const +{ + return impl_->data_; +} + +const MrcDensityMapHeader &MrcDensityMapOfFloatReader::header() const +{ + return impl_->header_; +} + +MrcDensityMapOfFloatReader::~MrcDensityMapOfFloatReader() +{ +} + +/******************************************************************** + * MrcDensityMapOfFloatWriter::Impl + */ + +/*! \internal \brief + * Private implementation class for MrcDensityMapOfFloatWriter. + */ +class MrcDensityMapOfFloatWriter::Impl +{ + public: + //! Construct mrc file writer by providing header and data to be written. + Impl(const MrcDensityMapHeader &header, ArrayRef data); + ~Impl() {} + //! Write the header and data from the writer to a given serialier + void write(ISerializer *serializer) const; + //! The mrc density map header data + const MrcDensityMapHeader header_; + //! The density data + const ArrayRef data_; +}; + +MrcDensityMapOfFloatWriter::Impl::Impl(const MrcDensityMapHeader &header, ArrayRef data) : header_(header), data_(data) +{ +} + +void MrcDensityMapOfFloatWriter::Impl::write(ISerializer *serializer) const +{ + if (serializer->reading()) + { + GMX_THROW(InternalError("Cannot use reading serializer to write.")); + } + + serializeMrcDensityMapHeader(serializer, header_); + + if (numberOfExpectedDataItems(header_) != data_.size()) + { + GMX_THROW(InternalError("Mrc data size does not match header information.")); + } + + for (float value : data_) + { + serializer->doFloat(&value); + } +} + +/******************************************************************** + * MrcDensityMapOfFloatWriter + */ + +MrcDensityMapOfFloatWriter::MrcDensityMapOfFloatWriter(const MrcDensityMapHeader &header, ArrayRef data) : impl_(new Impl(header, data)) +{ +} + +void MrcDensityMapOfFloatWriter::write(ISerializer *serializer) const +{ + impl_->write(serializer); +} + +MrcDensityMapOfFloatWriter::~MrcDensityMapOfFloatWriter() +{ +} + +} // namespace gmx diff --git a/src/gromacs/fileio/mrcdensitymap.h b/src/gromacs/fileio/mrcdensitymap.h new file mode 100644 index 0000000000..164859c309 --- /dev/null +++ b/src/gromacs/fileio/mrcdensitymap.h @@ -0,0 +1,103 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2019, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \libinternal \file + * \brief + * Declars mrc/ccp4-file format handling. + * + * \author Christian Blau + * + * \inlibraryapi + * \ingroup module_fileio + */ +#include "gromacs/utility/arrayref.h" +#include "gromacs/utility/classhelpers.h" + +namespace gmx +{ + +struct MrcDensityMapHeader; +class ISerializer; + +/*! \libinternal \brief Read an mrc/ccp4 file that contains float values. + */ +class MrcDensityMapOfFloatReader +{ + public: + /*! \brief Construct from directly de-serializing data into the object. + * \throws InternalError if serializer is not reading + * \throws InternalError if header is inconsistent + * \throws if serializer throws error upon failed reading + * \param[in] serializer Serializer to read the object data from + */ + explicit MrcDensityMapOfFloatReader(ISerializer * serializer); + + ~MrcDensityMapOfFloatReader(); + + //! Return the data vector that holds the data on the density grid + ArrayRef data() const; + //! Return the header + const MrcDensityMapHeader &header() const; + + private: + class Impl; + PrivateImplPointer impl_; + +}; + +/*! \libinternal \brief Write an mrc/ccp4 file that contains float values. + */ +class MrcDensityMapOfFloatWriter +{ + public: + /*! \brief Construct by setting the data and the header. + * + * \throws if the header data description does not match the provided data + * + * \param[in] header mrc density map header + * \param[in] data the density map data + */ + MrcDensityMapOfFloatWriter(const MrcDensityMapHeader &header, ArrayRef data); + + ~MrcDensityMapOfFloatWriter(); + + //! Serialize the mrc density data. + void write(ISerializer *serializer) const; + + private: + class Impl; + PrivateImplPointer impl_; +}; + +} // namespace gmx diff --git a/src/gromacs/fileio/mrcdensitymapheader.cpp b/src/gromacs/fileio/mrcdensitymapheader.cpp new file mode 100644 index 0000000000..bb8a148336 --- /dev/null +++ b/src/gromacs/fileio/mrcdensitymapheader.cpp @@ -0,0 +1,67 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2019, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief + * Implement mrc/ccp4-file metadata. + * + * \author Christian Blau + * + * \ingroup module_fileio + */ +#include "gmxpre.h" + +#include "mrcdensitymapheader.h" + +#include + +#include "gromacs/utility/exceptions.h" + +namespace gmx +{ +size_t numberOfExpectedDataItems(const MrcDensityMapHeader &header) +{ + if (std::any_of( + std::begin(header.numColumnRowSection_), std::end(header.numColumnRowSection_), + [](auto i) { return i < 0; })) + { + GMX_THROW(InternalError("Cannot determine data size, because the mrc " + "density map header is invalid (Negative number " + "describing data extent in at least one dimension).")); + } + + return header.numColumnRowSection_[XX] * header.numColumnRowSection_[YY] * header.numColumnRowSection_[ZZ]; +} + +} // namespace gmx diff --git a/src/gromacs/fileio/mrcdensitymapheader.h b/src/gromacs/fileio/mrcdensitymapheader.h index 2fbd24b663..46adf70818 100644 --- a/src/gromacs/fileio/mrcdensitymapheader.h +++ b/src/gromacs/fileio/mrcdensitymapheader.h @@ -159,5 +159,12 @@ struct MrcDensityMapHeader{ std::vector extendedHeader_ = {}; }; +/*! \brief Return the number of density data items that are expected + * to follow this header. + * \throws InternalError if the number of data items cannot be determined + * \returns the number of voxels + */ +size_t numberOfExpectedDataItems(const MrcDensityMapHeader &header); + } // namespace gmx #endif /* end of include guard: GMX_FILEIO_MRCDENSITYMAPHEADER_H */ diff --git a/src/gromacs/fileio/tests/CMakeLists.txt b/src/gromacs/fileio/tests/CMakeLists.txt index 6c97b18009..788cdaa6da 100644 --- a/src/gromacs/fileio/tests/CMakeLists.txt +++ b/src/gromacs/fileio/tests/CMakeLists.txt @@ -36,6 +36,8 @@ set(test_sources confio.cpp filemd5.cpp mrcserializer.cpp + mrcdensitymap.cpp + mrcdensitymapheader.cpp readinp.cpp ) if (GMX_USE_TNG) diff --git a/src/gromacs/fileio/tests/mrcdensitymap.cpp b/src/gromacs/fileio/tests/mrcdensitymap.cpp new file mode 100644 index 0000000000..d0fae6903a --- /dev/null +++ b/src/gromacs/fileio/tests/mrcdensitymap.cpp @@ -0,0 +1,90 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2019, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief + * Tests for mrc file data structure. + * + * \author Christian Blau + * \ingroup module_fileio + */ +#include "gmxpre.h" + +#include "gromacs/fileio/mrcdensitymap.h" + +#include +#include + +#include "gromacs/fileio/mrcdensitymapheader.h" +#include "gromacs/utility/inmemoryserializer.h" + +#include "testutils/testasserts.h" + +namespace gmx +{ + +namespace +{ + +TEST(MrcDensityMap, RoundTripIsIdempotent) +{ + // write header and data to serializer, store the serialized data + + MrcDensityMapHeader header {}; + header.numColumnRowSection_ = {1, 1, 1}; + + std::vector data(numberOfExpectedDataItems(header)); + + MrcDensityMapOfFloatWriter mrcDensityMapWriter(header, data); + InMemorySerializer serializer; + mrcDensityMapWriter.write(&serializer); + + const auto serializedFile = serializer.finishAndGetBuffer(); + + // read and write again + InMemoryDeserializer deserializer(serializedFile); + MrcDensityMapOfFloatReader mrcDensityMapReader(&deserializer); + + MrcDensityMapOfFloatWriter writerOfDeserializedOutput(mrcDensityMapReader.header(), mrcDensityMapReader.data()); + writerOfDeserializedOutput.write(&serializer); + + const auto roundTripResult = serializer.finishAndGetBuffer(); + + // compare serialized data after reading and writing again to serialized data after one-time serialization + EXPECT_THAT(serializedFile, testing::Pointwise(testing::Eq(), roundTripResult)); +} + +} // namespace + +} // namespace gmx diff --git a/src/gromacs/fileio/tests/mrcdensitymapheader.cpp b/src/gromacs/fileio/tests/mrcdensitymapheader.cpp new file mode 100644 index 0000000000..8359a84c65 --- /dev/null +++ b/src/gromacs/fileio/tests/mrcdensitymapheader.cpp @@ -0,0 +1,75 @@ +/* + * This file is part of the GROMACS molecular simulation package. + * + * Copyright (c) 2019, by the GROMACS development team, led by + * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, + * and including many others, as listed in the AUTHORS file in the + * top-level source directory and at http://www.gromacs.org. + * + * GROMACS is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public License + * as published by the Free Software Foundation; either version 2.1 + * of the License, or (at your option) any later version. + * + * GROMACS is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with GROMACS; if not, see + * http://www.gnu.org/licenses, or write to the Free Software Foundation, + * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + * + * If you want to redistribute modifications to GROMACS, please + * consider that scientific software is very special. Version + * control is crucial - bugs must be traceable. We will be happy to + * consider code for inclusion in the official distribution, but + * derived work must not be called official GROMACS. Details are found + * in the README & COPYING files - if they are missing, get the + * official version at http://www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the research papers on the package. Check out http://www.gromacs.org. + */ +/*! \internal \file + * \brief + * Tests for mrc file header structure. + * + * \author Christian Blau + * \ingroup module_fileio + */ + +#include "gmxpre.h" + +#include "gromacs/fileio/mrcdensitymapheader.h" + +#include +#include + +#include "testutils/testasserts.h" + +namespace gmx +{ + +TEST(MrcDensityMapHeaderTest, DataSizeIsZeroForDefaultHeader) +{ + MrcDensityMapHeader header; + EXPECT_EQ(0, numberOfExpectedDataItems(header)); +} + +TEST(MrcDensityMapHeaderTest, DataSizeIsCorrect) +{ + MrcDensityMapHeader header; + header.numColumnRowSection_ = {1, 2, 3}; + EXPECT_EQ(6, numberOfExpectedDataItems(header)); +} + +TEST(MrcDensityMapHeaderTest, DataSizeThrowsWhenInvalid) +{ + MrcDensityMapHeader header; + header.numColumnRowSection_ = {-1, 2, 3}; + EXPECT_THROW(numberOfExpectedDataItems(header), InternalError); +} + +} // namespace gmx -- 2.11.4.GIT