src/gromacs/analysisdata/modules/frameaverager.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2013,2014, 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  * \brief
  37  * Declares gmx::AnalysisDataFrameAverager.
  38  *
  39  * \author Teemu Murtola <teemu.murtola@gmail.com>
  40  * \ingroup module_analysisdata
  41  */
  42 #ifndef GMX_ANALYSISDATA_MODULES_FRAMEAVERAGER_H
  43 #define GMX_ANALYSISDATA_MODULES_FRAMEAVERAGER_H
  44
  45 #include <vector>
  46
  47 #include "gromacs/utility/gmxassert.h"
  48 #include "gromacs/utility/real.h"
  49
  50 namespace gmx
  51 {
  52
  53 class AnalysisDataPointSetRef;
  54
  55 /*! \internal
  56  * \brief
  57  * Helper class for modules that average values over frames.
  58  *
  59  * This class implements common functionality for analysis data modules that
  60  * need to average a set of values over frames.  Currently, it is designed for
  61  * computing averages for each input column independently, but should be
  62  * relatively easy to make more general if required.
  63  *
  64  * This class takes care of accumulating the values and computing their
  65  * variance.  It allows different number of samples for each input column.
  66  * Accumulation is always in double precision and uses a formula that is
  67  * relatively stable numerically.  For now, does nothing fancy,
  68  * but provides ground for other implementation (e.g., related to
  69  * parallelization) that would benefit all such modules.
  70  *
  71  * Methods in this class do not throw unless otherwise indicated.
  72  *
  73  * \ingroup module_analysisdata
  74  */
  75 class AnalysisDataFrameAverager
  76 {
  77     public:
  78         AnalysisDataFrameAverager() : bFinished_(false) {}
  79
  80         /*! \brief
  81          * Returns the number of columns in this averager.
  82          */
  83         int columnCount() const { return values_.size(); }
  84
  85         /*! \brief
  86          * Sets the number of columns in the input data.
  87          *
  88          * \throws std::bad_alloc if out of memory.
  89          *
  90          * Typically called from AnalysisDataModuleInterface::dataStarted().
  91          *
  92          * Must be called exactly once, before setting calling any other method
  93          * in the class.
  94          */
  95         void setColumnCount(int columnCount);
  96         /*! \brief
  97          * Adds a single value to the average for a given column.
  98          *
  99          * \param[in] index  Index of the column to add the value to.
 100          * \param[in] value  Value to add to the sample.
 101          */
 102         void addValue(int index, real value);
 103         /*! \brief
 104          * Accumulates data from a given point set into the average.
 105          *
 106          * Typically called from AnalysisDataModuleInterface::pointsAdded().
 107          *
 108          * Each call accumulates the values for those columns that are present
 109          * in the point set.  Can be called multiple times for a frame, and
 110          * does not need to be called for every frame.
 111          */
 112         void addPoints(const AnalysisDataPointSetRef &points);
 113         /*! \brief
 114          * Finalizes the calculation of the averages and variances.
 115          *
 116          * Does any computation that is not done during the accumulation in
 117          * addPoints().  Currently, does nothing, but provided as a placeholder
 118          * for more complex implementation.
 119          *
 120          * Typically called from AnalysisDataModuleInterface::dataFinished().
 121          */
 122         void finish();
 123
 124         /*! \brief
 125          * Returns the computed average for a given column.
 126          *
 127          * If called before finish(), the results are undefined.
 128          */
 129         real average(int index) const
 130         {
 131             GMX_ASSERT(index >= 0 && index < columnCount(),
 132                        "Invalid column index");
 133             GMX_ASSERT(bFinished_,
 134                        "Values available only after finished() has been called");
 135             return values_[index].average;
 136         }
 137         /*! \brief
 138          * Returns the computed (sample) variance for a given column.
 139          *
 140          * If called before finish(), the results are undefined.
 141          */
 142         real variance(int index) const
 143         {
 144             GMX_ASSERT(index >= 0 && index < columnCount(),
 145                        "Invalid column index");
 146             GMX_ASSERT(bFinished_,
 147                        "Values available only after finished() has been called");
 148             const AverageItem &item = values_[index];
 149             return item.samples > 1 ? item.squaredSum / (item.samples - 1) : 0.0;
 150         }
 151         /*! \brief
 152          * Returns the number of samples for a given column.
 153          *
 154          * If called before finish(), the results are undefined.
 155          */
 156         int sampleCount(int index) const
 157         {
 158             GMX_ASSERT(index >= 0 && index <= static_cast<int>(values_.size()),
 159                        "Invalid column index");
 160             GMX_ASSERT(bFinished_,
 161                        "Values available only after finished() has been called");
 162             return values_[index].samples;
 163         }
 164
 165     private:
 166         struct AverageItem
 167         {
 168             AverageItem() : average(0.0), squaredSum(0.0), samples(0) {}
 169
 170             //! Average of the values so far.
 171             double               average;
 172             //! Sum of squared deviations from the average for values so far.
 173             double               squaredSum;
 174             //! Number of values so far.
 175             int                  samples;
 176         };
 177
 178         std::vector<AverageItem> values_;
 179         bool                     bFinished_;
 180 };
 181
 182 } // namespace gmx
 183
 184 #endif