From c94c85dcaed7bc8132a3df39eb4a4d35ee141f69 Mon Sep 17 00:00:00 2001 From: Christian Blau Date: Fri, 8 Mar 2019 16:06:01 +0100 Subject: [PATCH] densityfitting - density file format (mrc) header. Implements mrc/cc4/map format file headers as decribed in "EMDB Map Distribution Format Description Version 1.01 (c) emdatabank.org 2014" Reading and writing is part of a follow-up patch. refs #2282 and #1877 Change-Id: Ifd43cedb137c0beeb320e1299df8bae9916c8ea8 --- src/gromacs/fileio/mrcdensitymapheader.h | 163 +++++++++++++++++++++++++++++++ 1 file changed, 163 insertions(+) create mode 100644 src/gromacs/fileio/mrcdensitymapheader.h diff --git a/src/gromacs/fileio/mrcdensitymapheader.h b/src/gromacs/fileio/mrcdensitymapheader.h new file mode 100644 index 0000000000..b73281e18b --- /dev/null +++ b/src/gromacs/fileio/mrcdensitymapheader.h @@ -0,0 +1,163 @@ +/* + * 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 + * Implement mrc/ccp4-file metadata. + * + * \author Christian Blau + * + * \inlibraryapi + * \ingroup module_fileio + */ +#ifndef GMX_FILEIO_MRCDENSITYMAPHEADER_H +#define GMX_FILEIO_MRCDENSITYMAPHEADER_H + +#include +#include + +#include "gromacs/math/vectypes.h" + +namespace gmx +{ + +/*! \brief Space group in three dimensions. + * + * Currently only "no symmetry" is supported, the complete enum class would hold + * 230 symbols. + * + * Table 12.3.4.1 Standard space-group symbols, pages 824-831, + * International Tables for Crystallography, Volume A, fifth edition + */ +enum class SpaceGroup : int32_t +{ + P1 = 1, //!< no symmetry +}; + +/*! \brief The type of density data stored in an mrc file. + * As named in + * "EMDB Map Distribution Format Description Version 1.01 (c) emdatabank.org 2014" + * Modes 0-4 are defined by the standard. + * NOTE only mode 2 is currently implemented and used. + */ +enum class MrcDataMode : int32_t +{ + uInt8 = 0, //!< compressed data mode, 8 bits, signed byte (range -128 to 127, ISO/IEC 10967) + int16 = 1, //!< 16 bits, signed integer (range -32768 to 32767, ISO/IEC 10967) + float32 = 2, //!< 32 bits, floating point number (IEEE 754) + complexInt32 = 3, //!< 32 bits, complex signed integers (ISO/IEC 10967) + complexFloat64 = 4, //!< 64 bits, complex floating point numbers (IEEE 754) +}; + +/*! \libinternal + * \brief Statistics about mrc data arrays. + */ +struct MrcDataStatistics +{ + float min_ = 0.; //!< Minimum data value scales values in (currently unsupported) compressed data mode. + float max_ = 0.; //!< Maximum data value scales values in (currently unsupported) compressed data mode. + float mean_ = 0.; //!< mean of the data + float rms_ = 0.; //!< rms of the data +}; + +/*! \libinternal + * \brief Skew matrix and translation. + * As named in + * "EMDB Map Distribution Format Description Version 1.01 (c) emdatabank.org 2014" + */ +struct MrcDensitySkewData +{ + bool valid_ = false; //!< True if skew matrix is stored. + std::array matrix_ = {}; //!< Skew matrix for crystallographic unit cell in Ångström + std::array translation_ = {}; //!< Translation of crystallographic unit cell in Ångström +}; + +/*! \libinternal + * \brief Crystallographic labels for mrc data. + */ +struct CrystallographicLabels +{ + static constexpr int c_labelSize = 80; //!< Length of crystallographic labels is eighty. + //! Number of used crystallographic labels, 0 for imagestacks, 1 for emdb data + int32_t numUsedLabels_ = 0; + + //! Crystallographic labels or "::::EMDataBank.org::::EMD-1234::::" for EMDB entries + std::array, 10> labels_ = {}; +}; + +/*! \libinternal + * \brief A container for the data in mrc density map file formats. + * + * Mrc files are a widely used file format in crystallography and cryo electron + * microscopy to represent volumetric data. Covers ccp4, map and imod. + * + * For a detailed description see + * "EMDB Map Distribution Format Description Version 1.01 (c) emdatabank.org 2014" + */ +struct MrcDensityMapHeader{ + + //! Space group of stored data. + SpaceGroup spaceGroup_ = SpaceGroup::P1; + //! Data mode, currently only mode 2 is supported + MrcDataMode dataMode_ = MrcDataMode::float32; + //! Identifies file format, expected to be "MAP " + std::array formatIdentifier_ = {'M', 'A', 'P', ' '}; + //! 15 unspecified float numbers + std::array userDefinedFloat_ = {}; + + //!Labels for crystallographic data + CrystallographicLabels labels_ = {}; + + //! Length of the crystallographic unit cell in Ångström + std::array cellLength_ = {1., 1., 1.}; + //! crystallographic unit cell angles + std::array cellAngles_ = {90., 90., 90.}; + + //! Data axis order with columns varying the fastest, and sections the slowest. + std::array columnRowSectionToXyz_ = {0, 1, 2}; + + std::array numColumnRowSection_ = {}; //!< Column, row and section count + std::array columnRowSectionStart_ = {}; //!< Start of values in grid + std::array extent_ = {}; //!< The number of grid points in the crystall cell + + //! Statistics about the data stored in the file. + MrcDataStatistics dataStatistics_ = {}; + //! Data to perform crystallographic unit cell skewing + MrcDensitySkewData skewData_ = {}; + //! Extended header with symmetry tables + std::vector extendedHeader_ = {}; +}; + +} // namespace gmx +#endif /* end of include guard: GMX_FILEIO_MRCDENSITYMAPHEADER_H */ -- 2.11.4.GIT