src/gromacs/fileio/mrcdensitymap.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2019, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*! \internal \file
  36  *
  37  * \brief
  38  * Implements mrc/ccp4-file format handling.
  39  *
  40  * \author Christian Blau <blau@kth.se>
  41  *
  42  * \ingroup module_fileio
  43  */
  44 #include "gmxpre.h"
  45
  46 #include "mrcdensitymap.h"
  47
  48 #include <algorithm>
  49 #include <vector>
  50
  51 #include "gromacs/fileio/gmxfio.h"
  52 #include "gromacs/fileio/gmxfio_xdr.h"
  53 #include "gromacs/fileio/mrcdensitymapheader.h"
  54 #include "gromacs/utility/exceptions.h"
  55 #include "gromacs/utility/inmemoryserializer.h"
  56 #include "gromacs/utility/iserializer.h"
  57
  58 #include "mrcserializer.h"
  59
  60 namespace gmx
  61 {
  62
  63 namespace
  64 {
  65
  66 /*! \brief Read file into memory as vector of chars
  67  *
  68  * \param[in] filename of the file to be read
  69  * \returns the file contents as a vector
  70  * \throws FileIOError if file not found
  71  * \throws FileIOError if reading was not successful
  72  */
  73 std::vector<char> readCharBufferFromFile(const std::string& filename)
  74 {
  75     if (!gmx_fexist(filename))
  76     {
  77         GMX_THROW(FileIOError("Error while reading '" + filename + "' - file not found."));
  78     }
  79     t_fileio* mrcFile = gmx_fio_open(filename.c_str(), "r");
  80
  81     // Determine file size
  82     gmx_fseek(gmx_fio_getfp(mrcFile), 0, SEEK_END);
  83     gmx_off_t fileSize = gmx_fio_ftell(mrcFile);
  84     gmx_fseek(gmx_fio_getfp(mrcFile), 0, SEEK_SET);
  85     // Read whole file into buffer the size of the file
  86     std::vector<char> fileContentBuffer(fileSize);
  87     size_t readSize = fread(fileContentBuffer.data(), sizeof(char), fileContentBuffer.size(),
  88                             gmx_fio_getfp(mrcFile));
  89     gmx_fio_close(mrcFile);
  90
  91     if (fileContentBuffer.size() != readSize)
  92     {
  93         GMX_THROW(FileIOError("Error while reading '" + filename
  94                               + "' - file size and read buffer size do not match."));
  95     }
  96
  97     return fileContentBuffer;
  98 }
  99
 100 } // namespace
 101
 102 /********************************************************************
 103  * MrcDensityMapOfFloatReader::Impl
 104  */
 105
 106 /*! \internal \brief
 107  * Private implementation class for MrcDensityMapOfFloatReader.
 108  */
 109 class MrcDensityMapOfFloatReader::Impl
 110 {
 111 public:
 112     //! Build the map reader from a serializer.
 113     explicit Impl(ISerializer* serializer);
 114     ~Impl() {}
 115     //! The header of the read mrc file
 116     MrcDensityMapHeader header_;
 117     //! The data of the mrc file
 118     std::vector<float> data_;
 119 };
 120
 121 MrcDensityMapOfFloatReader::Impl::Impl(ISerializer* serializer)
 122 {
 123     if (!serializer->reading())
 124     {
 125         GMX_THROW(InternalError("Cannot use writing serializer to read."));
 126     }
 127
 128     header_             = deserializeMrcDensityMapHeader(serializer);
 129     const auto dataSize = numberOfExpectedDataItems(header_);
 130     data_.resize(dataSize);
 131     for (auto& value : data_)
 132     {
 133         serializer->doFloat(&value);
 134     }
 135 }
 136
 137 /********************************************************************
 138  * MrcDensityMapOfFloatReader
 139  */
 140
 141 MrcDensityMapOfFloatReader::MrcDensityMapOfFloatReader(ISerializer* serializer) :
 142     impl_(new Impl(serializer))
 143 {
 144 }
 145
 146 ArrayRef<const float> MrcDensityMapOfFloatReader::constView() const
 147 {
 148     return impl_->data_;
 149 }
 150
 151 const MrcDensityMapHeader& MrcDensityMapOfFloatReader::header() const
 152 {
 153     return impl_->header_;
 154 }
 155
 156 MrcDensityMapOfFloatReader::~MrcDensityMapOfFloatReader() {}
 157
 158 /********************************************************************
 159  * MrcDensityMapOfFloatFromFileReader::Impl
 160  */
 161
 162
 163 class MrcDensityMapOfFloatFromFileReader::Impl
 164 {
 165 public:
 166     explicit Impl(const std::string& fileName);
 167     ~Impl() = default;
 168     const MrcDensityMapOfFloatReader& reader() const;
 169
 170 private:
 171     const std::vector<char>                     buffer_;
 172     std::unique_ptr<InMemoryDeserializer>       serializer_;
 173     std::unique_ptr<MrcDensityMapOfFloatReader> reader_;
 174 };
 175
 176 MrcDensityMapOfFloatFromFileReader::Impl::Impl(const std::string& filename) :
 177     buffer_(readCharBufferFromFile(filename)),
 178     serializer_(std::make_unique<InMemoryDeserializer>(buffer_, false)),
 179     reader_(std::make_unique<MrcDensityMapOfFloatReader>(serializer_.get()))
 180 {
 181     if (!mrcHeaderIsSane(reader_->header()))
 182     {
 183         serializer_ = std::make_unique<InMemoryDeserializer>(buffer_, false, EndianSwapBehavior::Swap);
 184         reader_     = std::make_unique<MrcDensityMapOfFloatReader>(serializer_.get());
 185         if (!mrcHeaderIsSane(reader_->header()))
 186         {
 187             GMX_THROW(FileIOError(
 188                     "Header of '" + filename
 189                     + "' fails sanity check for little- as well as big-endian reading."));
 190         }
 191     }
 192
 193     layout_right::mapping<dynamicExtents3D> map(getDynamicExtents3D(reader_->header()));
 194     if (map.required_span_size() != reader_->constView().ssize())
 195     {
 196         GMX_THROW(FileIOError("File header density extent information of " + filename
 197                               + "' does not match density data size"));
 198     }
 199 }
 200
 201 const MrcDensityMapOfFloatReader& MrcDensityMapOfFloatFromFileReader::Impl::reader() const
 202 {
 203     return *reader_;
 204 }
 205
 206 /********************************************************************
 207  * MrcDensityMapOfFloatFromFileReader
 208  */
 209
 210 MrcDensityMapOfFloatFromFileReader::MrcDensityMapOfFloatFromFileReader(const std::string& filename) :
 211     impl_(new Impl(filename))
 212 {
 213 }
 214
 215 MrcDensityMapOfFloatFromFileReader::~MrcDensityMapOfFloatFromFileReader() = default;
 216
 217 TranslateAndScale MrcDensityMapOfFloatFromFileReader::transformationToDensityLattice() const
 218 {
 219     return getCoordinateTransformationToLattice(impl_->reader().header());
 220 }
 221
 222 MultiDimArray<std::vector<float>, dynamicExtents3D> MrcDensityMapOfFloatFromFileReader::densityDataCopy() const
 223 {
 224     MultiDimArray<std::vector<float>, dynamicExtents3D> result(
 225             getDynamicExtents3D(impl_->reader().header()));
 226     std::copy(std::begin(impl_->reader().constView()), std::end(impl_->reader().constView()),
 227               begin(result.asView()));
 228     return result;
 229 }
 230
 231 /********************************************************************
 232  * MrcDensityMapOfFloatWriter::Impl
 233  */
 234
 235 /*! \internal \brief
 236  * Private implementation class for MrcDensityMapOfFloatWriter.
 237  */
 238 class MrcDensityMapOfFloatWriter::Impl
 239 {
 240 public:
 241     //! Construct mrc file writer by providing header and data to be written.
 242     Impl(const MrcDensityMapHeader& header, ArrayRef<const float> data);
 243     ~Impl() {}
 244     //! Write the header and data from the writer to a given serialier
 245     void write(ISerializer* serializer) const;
 246     //! The mrc density map header data
 247     const MrcDensityMapHeader header_;
 248     //! The density data
 249     const ArrayRef<const float> data_;
 250 };
 251
 252 MrcDensityMapOfFloatWriter::Impl::Impl(const MrcDensityMapHeader& header, ArrayRef<const float> data) :
 253     header_(header),
 254     data_(data)
 255 {
 256 }
 257
 258 void MrcDensityMapOfFloatWriter::Impl::write(ISerializer* serializer) const
 259 {
 260     if (serializer->reading())
 261     {
 262         GMX_THROW(InternalError("Cannot use reading serializer to write."));
 263     }
 264
 265     serializeMrcDensityMapHeader(serializer, header_);
 266
 267     if (numberOfExpectedDataItems(header_) != data_.size())
 268     {
 269         GMX_THROW(InternalError("Mrc data size does not match header information."));
 270     }
 271
 272     for (float value : data_)
 273     {
 274         serializer->doFloat(&value);
 275     }
 276 }
 277
 278 /********************************************************************
 279  * MrcDensityMapOfFloatWriter
 280  */
 281
 282 MrcDensityMapOfFloatWriter::MrcDensityMapOfFloatWriter(const MrcDensityMapHeader& header,
 283                                                        ArrayRef<const float>      data) :
 284     impl_(new Impl(header, data))
 285 {
 286 }
 287
 288 void MrcDensityMapOfFloatWriter::write(ISerializer* serializer) const
 289 {
 290     impl_->write(serializer);
 291 }
 292
 293 MrcDensityMapOfFloatWriter::~MrcDensityMapOfFloatWriter() {}
 294
 295 } // namespace gmx