yat/statistics/Histogram.h

Code
Comments
Other
Rev Date Author Line
680 11 Oct 06 jari 1 #ifndef _theplu_yat_statistics_histogram_
680 11 Oct 06 jari 2 #define _theplu_yat_statistics_histogram_
195 27 Oct 04 jari 3
675 10 Oct 06 jari 4 // $Id$
295 29 Apr 05 peter 5
675 10 Oct 06 jari 6 /*
2119 12 Dec 09 peter 7   Copyright (C) 2004 Jari Häkkinen
2119 12 Dec 09 peter 8   Copyright (C) 2005 Jari Häkkinen, Peter Johansson
2119 12 Dec 09 peter 9   Copyright (C) 2006 Jari Häkkinen
4359 23 Aug 23 peter 10   Copyright (C) 2007 Peter Johansson
4359 23 Aug 23 peter 11   Copyright (C) 2008 Jari Häkkinen, Peter Johansson
4359 23 Aug 23 peter 12   Copyright (C) 2009, 2010, 2016 Peter Johansson
675 10 Oct 06 jari 13
1437 25 Aug 08 peter 14   This file is part of the yat library, http://dev.thep.lu.se/yat
675 10 Oct 06 jari 15
675 10 Oct 06 jari 16   The yat library is free software; you can redistribute it and/or
675 10 Oct 06 jari 17   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 18   published by the Free Software Foundation; either version 3 of the
675 10 Oct 06 jari 19   License, or (at your option) any later version.
675 10 Oct 06 jari 20
675 10 Oct 06 jari 21   The yat library is distributed in the hope that it will be useful,
675 10 Oct 06 jari 22   but WITHOUT ANY WARRANTY; without even the implied warranty of
675 10 Oct 06 jari 23   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
675 10 Oct 06 jari 24   General Public License for more details.
675 10 Oct 06 jari 25
675 10 Oct 06 jari 26   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 27   along with yat. If not, see <http://www.gnu.org/licenses/>.
675 10 Oct 06 jari 28 */
675 10 Oct 06 jari 29
680 11 Oct 06 jari 30 #include "AveragerWeighted.h"
2263 26 May 10 peter 31 #include "yat/utility/concept_check.h"
1145 25 Feb 08 peter 32 #include "yat/utility/iterator_traits.h"
675 10 Oct 06 jari 33
2202 21 Feb 10 peter 34 #include <boost/concept_check.hpp>
3522 10 Oct 16 peter 35 #include <boost/iterator/iterator_concepts.hpp>
2202 21 Feb 10 peter 36
2333 15 Oct 10 peter 37 #include <iosfwd>
195 27 Oct 04 jari 38 #include <string>
195 27 Oct 04 jari 39 #include <vector>
195 27 Oct 04 jari 40
195 27 Oct 04 jari 41 namespace theplu {
680 11 Oct 06 jari 42 namespace yat {
295 29 Apr 05 peter 43 namespace statistics {
195 27 Oct 04 jari 44
195 27 Oct 04 jari 45   ///
767 22 Feb 07 peter 46   /// @brief Histograms provide a convenient way of presenting the
4200 19 Aug 22 peter 47   /// distribution of a set of data.
767 22 Feb 07 peter 48   ///
767 22 Feb 07 peter 49   /// A histogram consists of a set of
195 27 Oct 04 jari 50   /// bins which count the number of events falling into these
195 27 Oct 04 jari 51   /// bins. Currently only one dimensional histograms with uniformly
195 27 Oct 04 jari 52   /// spaced bins are supported.
195 27 Oct 04 jari 53   ///
195 27 Oct 04 jari 54   class Histogram
195 27 Oct 04 jari 55   {
195 27 Oct 04 jari 56   public:
195 27 Oct 04 jari 57
195 27 Oct 04 jari 58     ///
195 27 Oct 04 jari 59     /// The default constructor.
195 27 Oct 04 jari 60     ///
195 27 Oct 04 jari 61     Histogram(void);
195 27 Oct 04 jari 62
2333 15 Oct 10 peter 63     /**
2333 15 Oct 10 peter 64        \brief The istream constructor.
2333 15 Oct 10 peter 65
2333 15 Oct 10 peter 66        This constructor creates a Histogram from output (possibly)
2333 15 Oct 10 peter 67        created with operator<<. Note that the output from operator<<
2333 15 Oct 10 peter 68        does not contain all information required to restore the
2333 15 Oct 10 peter 69        Histogram. The values of all bins are restored, but the
2333 15 Oct 10 peter 70        averagers, averager_all(void) and averager_histogram(void), and
2333 15 Oct 10 peter 71        more specifically the second moment (sum_xx) in the averagers
2333 15 Oct 10 peter 72        cannot be restored due to lack of information in the output.
2333 15 Oct 10 peter 73
2333 15 Oct 10 peter 74        \since New in yat 0.7
2333 15 Oct 10 peter 75     */
2333 15 Oct 10 peter 76     Histogram(std::istream&);
2333 15 Oct 10 peter 77
195 27 Oct 04 jari 78     ///
195 27 Oct 04 jari 79     /// The copy constructor.
195 27 Oct 04 jari 80     ///
195 27 Oct 04 jari 81     Histogram(const Histogram&);
195 27 Oct 04 jari 82
195 27 Oct 04 jari 83     ///
199 27 Oct 04 jari 84     /// Construct a histogram object that covers \f$(xmin,xmax]\f$
195 27 Oct 04 jari 85     /// with the bin spacing \f$(xmax-xmin)/n\f$.
195 27 Oct 04 jari 86     ///
195 27 Oct 04 jari 87     Histogram(const double xmin, const double xmax, const size_t n);
195 27 Oct 04 jari 88
195 27 Oct 04 jari 89     virtual ~Histogram(void);
195 27 Oct 04 jari 90
195 27 Oct 04 jari 91     ///
195 27 Oct 04 jari 92     /// Update the histogram by adding \a weight to the bin whose
195 27 Oct 04 jari 93     /// range contains the observation \a x. No bins are updated when
195 27 Oct 04 jari 94     /// \a x lies outside the range of the histogram but the value is
195 27 Oct 04 jari 95     /// added to the overall integral of the histogram.
195 27 Oct 04 jari 96     ///
1203 05 Mar 08 peter 97     /// @short Add a data point to the histogram.
195 27 Oct 04 jari 98     ///
195 27 Oct 04 jari 99     /// @return 0 if \a x lies within the range of the histogram, -1
195 27 Oct 04 jari 100     /// if \a x is smaller than the lower limit of the histogram, and
195 27 Oct 04 jari 101     /// similarly, 1 is returned if \a x is greater than or equal to
195 27 Oct 04 jari 102     /// the upper limit.
195 27 Oct 04 jari 103     ///
199 27 Oct 04 jari 104     int add(const double x,const double weight=1.0);
195 27 Oct 04 jari 105
195 27 Oct 04 jari 106     ///
295 29 Apr 05 peter 107     /// Gives access to the AveragerWeighted object that keeps track of
195 27 Oct 04 jari 108     /// average of all events presented to the histogram.
195 27 Oct 04 jari 109     ///
195 27 Oct 04 jari 110     /// @short Average of all events presented to the histogram.
195 27 Oct 04 jari 111     ///
295 29 Apr 05 peter 112     /// @return A const reference to an AveragerWeighted object.
195 27 Oct 04 jari 113     ///
1203 05 Mar 08 peter 114     const AveragerWeighted& averager_all(void) const;
195 27 Oct 04 jari 115
195 27 Oct 04 jari 116     ///
295 29 Apr 05 peter 117     /// Gives access to the AveragerWeighted object that keeps track of
195 27 Oct 04 jari 118     /// average of events that fits within the histogram lower and
195 27 Oct 04 jari 119     /// upper limits. This function is equivalent to averager().
195 27 Oct 04 jari 120     ///
195 27 Oct 04 jari 121     /// @short Average of events fitting within histogram.
195 27 Oct 04 jari 122     ///
295 29 Apr 05 peter 123     /// @return A const reference to an AveragerWeighted object.
195 27 Oct 04 jari 124     ///
1203 05 Mar 08 peter 125     const AveragerWeighted& averager_histogram(void) const;
195 27 Oct 04 jari 126
215 06 Nov 04 jari 127     ///
215 06 Nov 04 jari 128     /// @return The number of bins in the histogram
215 06 Nov 04 jari 129     ///
718 26 Dec 06 jari 130     size_t nof_bins(void) const;
195 27 Oct 04 jari 131
215 06 Nov 04 jari 132     ///
372 05 Aug 05 peter 133     ///  There are two ways to normalize the counts.
215 06 Nov 04 jari 134     ///
372 05 Aug 05 peter 135     /// If choice is true: The normalized count is the count in a
372 05 Aug 05 peter 136     /// bin divided by the total number of observations. In this
372 05 Aug 05 peter 137     /// case the relative counts are normalized to sum to unity (
372 05 Aug 05 peter 138     /// minus values outside histogram).  This is the intuitive case
372 05 Aug 05 peter 139     /// where the height of the histogram bar represents the
372 05 Aug 05 peter 140     /// proportion of the data in each class.
363 04 Aug 05 jari 141     ///
372 05 Aug 05 peter 142     /// If choice is false: The normalized count is the count in the
372 05 Aug 05 peter 143     /// class divided by the number of observations times the bin
372 05 Aug 05 peter 144     /// width. For this normalization, the area (or integral) under
372 05 Aug 05 peter 145     /// the histogram is equal to unity (minus the missing area
372 05 Aug 05 peter 146     /// corresponding to counts outside histogram). From a
372 05 Aug 05 peter 147     /// probabilistic point of view, this normalization results in a
372 05 Aug 05 peter 148     /// relative histogram that is most akin to the probability
372 05 Aug 05 peter 149     /// density function If you want to overlay a probability density
372 05 Aug 05 peter 150     /// on top of the histogram, use this normalization. Although this
372 05 Aug 05 peter 151     /// normalization is less intuitive (relative frequencies greater
372 05 Aug 05 peter 152     /// than 1 are quite permissible), it is the appropriate
372 05 Aug 05 peter 153     /// normalization if you are using the histogram to model a
372 05 Aug 05 peter 154     /// probability density function.
372 05 Aug 05 peter 155     ///
372 05 Aug 05 peter 156     /// @short Normalizing the histogram
4200 19 Aug 22 peter 157     ///
372 05 Aug 05 peter 158     void normalize(bool choice = true);
215 06 Nov 04 jari 159
215 06 Nov 04 jari 160     ///
370 05 Aug 05 peter 161     /// @return The value in the middle of bin \a k.
215 06 Nov 04 jari 162     ///
370 05 Aug 05 peter 163     /// @note No check is done that \a k is within the size of the
215 06 Nov 04 jari 164     /// histogram.
215 06 Nov 04 jari 165     ///
718 26 Dec 06 jari 166     double observation_value(const size_t k) const;
215 06 Nov 04 jari 167
2038 25 Aug 09 peter 168     /**
2038 25 Aug 09 peter 169        multiply value of each bin with \a factor
4200 19 Aug 22 peter 170
2038 25 Aug 09 peter 171        \since New in yat 0.6
2038 25 Aug 09 peter 172      */
2038 25 Aug 09 peter 173     void rescale(double factor);
2038 25 Aug 09 peter 174
215 06 Nov 04 jari 175     ///
215 06 Nov 04 jari 176     /// Set everyting to default values, here it means that everything
215 06 Nov 04 jari 177     /// is set to zero except the boundary values that are kept.
215 06 Nov 04 jari 178     ///
195 27 Oct 04 jari 179     void reset(void);
215 06 Nov 04 jari 180
215 06 Nov 04 jari 181     ///
372 05 Aug 05 peter 182     /// @return The width of the bins in the histogram.
215 06 Nov 04 jari 183     ///
718 26 Dec 06 jari 184     double spacing(void) const;
215 06 Nov 04 jari 185
215 06 Nov 04 jari 186     ///
215 06 Nov 04 jari 187     /// @return The histogram upper boundary.
215 06 Nov 04 jari 188     ///
372 05 Aug 05 peter 189     /// @note The upper boundary value is outside the histogram.
215 06 Nov 04 jari 190     ///
718 26 Dec 06 jari 191     double xmax(void) const;
215 06 Nov 04 jari 192
215 06 Nov 04 jari 193     ///
215 06 Nov 04 jari 194     /// @return The histogram lower boundary.
215 06 Nov 04 jari 195     ///
372 05 Aug 05 peter 196     /// @note The lower boundary value is inside the histogram.
215 06 Nov 04 jari 197     ///
718 26 Dec 06 jari 198     double xmin(void) const;
195 27 Oct 04 jari 199
215 06 Nov 04 jari 200     ///
372 05 Aug 05 peter 201     /// @return The count of bin \a k in the histogram.
372 05 Aug 05 peter 202     ///
718 26 Dec 06 jari 203     double operator[](size_t k) const;
199 27 Oct 04 jari 204
215 06 Nov 04 jari 205     ///
215 06 Nov 04 jari 206     /// The assignment operator
372 05 Aug 05 peter 207     ///
195 27 Oct 04 jari 208     const Histogram& operator=(const Histogram&);
195 27 Oct 04 jari 209
2032 19 Aug 09 peter 210     /**
2032 19 Aug 09 peter 211        \brief Addition operator
2032 19 Aug 09 peter 212
2032 19 Aug 09 peter 213        \note The bins in lhs and rhs must match, i.e., min, max, and
2032 19 Aug 09 peter 214        nof_bins() must be equal in lhs and rhs.
2032 19 Aug 09 peter 215
2032 19 Aug 09 peter 216        \since New in yat 0.6
2032 19 Aug 09 peter 217     */
2032 19 Aug 09 peter 218     const Histogram& operator+=(const Histogram& rhs);
2032 19 Aug 09 peter 219
195 27 Oct 04 jari 220   private:
231 21 Feb 05 peter 221     // Returns zero if outside boundaries
718 26 Dec 06 jari 222     size_t bin(double d);
215 06 Nov 04 jari 223
199 27 Oct 04 jari 224     std::vector<double> histogram_;
195 27 Oct 04 jari 225     double xmax_;
195 27 Oct 04 jari 226     double xmin_;
295 29 Apr 05 peter 227     statistics::AveragerWeighted sum_all_;      // average of all data
295 29 Apr 05 peter 228     statistics::AveragerWeighted sum_histogram_;// average of data in histogram
195 27 Oct 04 jari 229   };
195 27 Oct 04 jari 230
1145 25 Feb 08 peter 231   /**
1145 25 Feb 08 peter 232      Add a range [first, last) of values to Histogram.
1886 31 Mar 09 peter 233
3522 10 Oct 16 peter 234      Type Requirements:
3522 10 Oct 16 peter 235      - \c Iterator models a \ref concept_data_iterator
3522 10 Oct 16 peter 236      - \c Iterator models a \single_pass_iterator
2202 21 Feb 10 peter 237
1886 31 Mar 09 peter 238      \relates Histogram
1145 25 Feb 08 peter 239    */
3522 10 Oct 16 peter 240   template<typename Iterator>
3522 10 Oct 16 peter 241   void add(Histogram& h, Iterator first, Iterator last)
1145 25 Feb 08 peter 242   {
3522 10 Oct 16 peter 243     BOOST_CONCEPT_ASSERT((boost_concepts::SinglePassIterator<Iterator>));
3522 10 Oct 16 peter 244     BOOST_CONCEPT_ASSERT((utility::DataIteratorConcept<Iterator>));
3522 10 Oct 16 peter 245     utility::iterator_traits<Iterator> traits;
1145 25 Feb 08 peter 246     while (first!=last) {
2202 21 Feb 10 peter 247       h.add(traits.data(first), traits.weight(first));
1145 25 Feb 08 peter 248       ++first;
1145 25 Feb 08 peter 249     }
1145 25 Feb 08 peter 250   }
1145 25 Feb 08 peter 251
1886 31 Mar 09 peter 252   /**
1886 31 Mar 09 peter 253      The Histogram output operator
3522 10 Oct 16 peter 254
1886 31 Mar 09 peter 255      \relates Histogram
1886 31 Mar 09 peter 256   */
1886 31 Mar 09 peter 257   std::ostream& operator<<(std::ostream& s,const Histogram&);
195 27 Oct 04 jari 258
683 11 Oct 06 jari 259 }}} // of namespace statistics, yat, and theplu
195 27 Oct 04 jari 260
195 27 Oct 04 jari 261 #endif