yat/statistics/Histogram.cc

Code
Comments
Other
Rev Date Author Line
195 27 Oct 04 jari 1 // $Id$
195 27 Oct 04 jari 2
675 10 Oct 06 jari 3 /*
2119 12 Dec 09 peter 4   Copyright (C) 2004 Jari Häkkinen
831 27 Mar 07 peter 5   Copyright (C) 2005 Peter Johansson
2119 12 Dec 09 peter 6   Copyright (C) 2006 Jari Häkkinen
4359 23 Aug 23 peter 7   Copyright (C) 2007 Peter Johansson
4359 23 Aug 23 peter 8   Copyright (C) 2008 Jari Häkkinen, Peter Johansson
4359 23 Aug 23 peter 9   Copyright (C) 2009, 2010, 2011, 2012 Peter Johansson
195 27 Oct 04 jari 10
1437 25 Aug 08 peter 11   This file is part of the yat library, http://dev.thep.lu.se/yat
675 10 Oct 06 jari 12
675 10 Oct 06 jari 13   The yat library is free software; you can redistribute it and/or
675 10 Oct 06 jari 14   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 15   published by the Free Software Foundation; either version 3 of the
675 10 Oct 06 jari 16   License, or (at your option) any later version.
675 10 Oct 06 jari 17
675 10 Oct 06 jari 18   The yat library is distributed in the hope that it will be useful,
675 10 Oct 06 jari 19   but WITHOUT ANY WARRANTY; without even the implied warranty of
675 10 Oct 06 jari 20   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
675 10 Oct 06 jari 21   General Public License for more details.
675 10 Oct 06 jari 22
675 10 Oct 06 jari 23   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 24   along with yat. If not, see <http://www.gnu.org/licenses/>.
675 10 Oct 06 jari 25 */
675 10 Oct 06 jari 26
2881 18 Nov 12 peter 27 #include <config.h>
2881 18 Nov 12 peter 28
680 11 Oct 06 jari 29 #include "Histogram.h"
675 10 Oct 06 jari 30
2333 15 Oct 10 peter 31 #include <yat/utility/utility.h>
2333 15 Oct 10 peter 32
2032 19 Aug 09 peter 33 #include <algorithm>
2032 19 Aug 09 peter 34 #include <cassert>
199 27 Oct 04 jari 35 #include <cmath>
2032 19 Aug 09 peter 36 #include <functional>
2333 15 Oct 10 peter 37 #include <istream>
2032 19 Aug 09 peter 38 #include <ostream>
195 27 Oct 04 jari 39
195 27 Oct 04 jari 40 namespace theplu {
680 11 Oct 06 jari 41 namespace yat {
298 29 Apr 05 peter 42 namespace statistics {
195 27 Oct 04 jari 43
195 27 Oct 04 jari 44
718 26 Dec 06 jari 45   Histogram::Histogram(void)
718 26 Dec 06 jari 46     : xmax_(0), xmin_(0), sum_all_(), sum_histogram_()
718 26 Dec 06 jari 47   {
718 26 Dec 06 jari 48   }
195 27 Oct 04 jari 49
195 27 Oct 04 jari 50
2333 15 Oct 10 peter 51   Histogram::Histogram(std::istream& is)
2333 15 Oct 10 peter 52   {
2333 15 Oct 10 peter 53     std::string line;
2333 15 Oct 10 peter 54     getline(is, line, ':');
2333 15 Oct 10 peter 55     getline(is, line, '\n');
2333 15 Oct 10 peter 56     xmin_ = utility::convert<double>(line);
2333 15 Oct 10 peter 57     getline(is, line, ':');
2333 15 Oct 10 peter 58     getline(is, line, '\n');
2333 15 Oct 10 peter 59     xmax_ = utility::convert<double>(line);
2333 15 Oct 10 peter 60     getline(is, line, ':');
2333 15 Oct 10 peter 61     getline(is, line, '\n');
2333 15 Oct 10 peter 62     size_t n = utility::convert<size_t>(line);
2333 15 Oct 10 peter 63     histogram_.resize(n);
2333 15 Oct 10 peter 64     getline(is, line, '\n');
2333 15 Oct 10 peter 65     getline(is, line, ':');
2333 15 Oct 10 peter 66     getline(is, line, '\n');
2333 15 Oct 10 peter 67     double total = utility::convert<double>(line);
2333 15 Oct 10 peter 68     getline(is, line, '\n');
2333 15 Oct 10 peter 69     getline(is, line, '\n');
2333 15 Oct 10 peter 70     std::vector<std::vector<double> > data;
2411 13 Jan 11 peter 71     utility::load(is, data, '\0', '\n');
2333 15 Oct 10 peter 72     for (size_t i=0; i<histogram_.size(); ++i) {
2333 15 Oct 10 peter 73       assert(data[i].size()==2);
2333 15 Oct 10 peter 74       add(data[i][0], data[i][1]);
2333 15 Oct 10 peter 75     }
2333 15 Oct 10 peter 76     add(xmax_, total - averager_all().sum_w());
2333 15 Oct 10 peter 77   }
2333 15 Oct 10 peter 78
2333 15 Oct 10 peter 79
718 26 Dec 06 jari 80   Histogram::Histogram(const Histogram& b)
718 26 Dec 06 jari 81   {
718 26 Dec 06 jari 82     *this=b;
718 26 Dec 06 jari 83   }
195 27 Oct 04 jari 84
195 27 Oct 04 jari 85
718 26 Dec 06 jari 86   Histogram::Histogram(const double min, const double max, const size_t n)
718 26 Dec 06 jari 87     : histogram_(std::vector<double>(n,0.0)),
718 26 Dec 06 jari 88       xmax_(max), xmin_(min),
718 26 Dec 06 jari 89       sum_all_(), sum_histogram_()
718 26 Dec 06 jari 90   {
718 26 Dec 06 jari 91   }
195 27 Oct 04 jari 92
195 27 Oct 04 jari 93
718 26 Dec 06 jari 94   Histogram::~Histogram(void)
718 26 Dec 06 jari 95   {
718 26 Dec 06 jari 96   }
195 27 Oct 04 jari 97
195 27 Oct 04 jari 98
718 26 Dec 06 jari 99   int Histogram::add(const double x, const double w)
718 26 Dec 06 jari 100   {
718 26 Dec 06 jari 101     sum_all_.add(x,w);
718 26 Dec 06 jari 102     if (x<xmin_)
718 26 Dec 06 jari 103       return -1;
718 26 Dec 06 jari 104     else if (x>=xmax_)
718 26 Dec 06 jari 105       return 1;
195 27 Oct 04 jari 106
718 26 Dec 06 jari 107     sum_histogram_.add(x,w);
718 26 Dec 06 jari 108     histogram_[bin(x)] += w;
718 26 Dec 06 jari 109     return 0;
718 26 Dec 06 jari 110   }
195 27 Oct 04 jari 111
195 27 Oct 04 jari 112
1203 05 Mar 08 peter 113   const AveragerWeighted& Histogram::averager_all(void) const
718 26 Dec 06 jari 114   {
718 26 Dec 06 jari 115     return sum_all_;
718 26 Dec 06 jari 116   }
195 27 Oct 04 jari 117
195 27 Oct 04 jari 118
1203 05 Mar 08 peter 119   const AveragerWeighted& Histogram::averager_histogram(void) const
718 26 Dec 06 jari 120   {
718 26 Dec 06 jari 121     return sum_histogram_;
718 26 Dec 06 jari 122   }
195 27 Oct 04 jari 123
195 27 Oct 04 jari 124
718 26 Dec 06 jari 125   size_t Histogram::bin(double d)
718 26 Dec 06 jari 126   {
718 26 Dec 06 jari 127     return (((d<xmin_) || (d>xmax_)) ? 0 :
718 26 Dec 06 jari 128             static_cast<size_t>(floor((d-xmin_)/spacing() )));
718 26 Dec 06 jari 129   }
372 05 Aug 05 peter 130
195 27 Oct 04 jari 131
718 26 Dec 06 jari 132   size_t Histogram::nof_bins(void) const
718 26 Dec 06 jari 133   {
718 26 Dec 06 jari 134     return histogram_.size();
718 26 Dec 06 jari 135   }
195 27 Oct 04 jari 136
195 27 Oct 04 jari 137
718 26 Dec 06 jari 138   void Histogram::normalize(bool choice)
718 26 Dec 06 jari 139   {
4200 19 Aug 22 peter 140     if (choice)
2038 25 Aug 09 peter 141       rescale(1.0/sum_all_.sum_w());
718 26 Dec 06 jari 142     else
2038 25 Aug 09 peter 143       rescale(1.0/sum_all_.sum_w()/spacing());
718 26 Dec 06 jari 144   }
195 27 Oct 04 jari 145
195 27 Oct 04 jari 146
718 26 Dec 06 jari 147   double Histogram::observation_value(const size_t k) const
718 26 Dec 06 jari 148   {
718 26 Dec 06 jari 149     return xmin_+spacing()*(k+0.5);
718 26 Dec 06 jari 150   }
195 27 Oct 04 jari 151
195 27 Oct 04 jari 152
2038 25 Aug 09 peter 153   void Histogram::rescale(double factor)
2038 25 Aug 09 peter 154   {
2038 25 Aug 09 peter 155     for (size_t i=0; i<histogram_.size(); i++)
2038 25 Aug 09 peter 156       histogram_[i]*=factor;
2038 25 Aug 09 peter 157   }
2038 25 Aug 09 peter 158
2038 25 Aug 09 peter 159
718 26 Dec 06 jari 160   void Histogram::reset(void)
718 26 Dec 06 jari 161   {
1271 09 Apr 08 peter 162     for (size_t i=0; i<histogram_.size(); i++)
718 26 Dec 06 jari 163       histogram_[i]=0;
718 26 Dec 06 jari 164     sum_all_.reset();
718 26 Dec 06 jari 165     sum_histogram_.reset();
718 26 Dec 06 jari 166   }
195 27 Oct 04 jari 167
195 27 Oct 04 jari 168
718 26 Dec 06 jari 169   double Histogram::spacing(void) const
718 26 Dec 06 jari 170   {
718 26 Dec 06 jari 171     return (xmax_-xmin_)/nof_bins();
718 26 Dec 06 jari 172   }
195 27 Oct 04 jari 173
195 27 Oct 04 jari 174
718 26 Dec 06 jari 175   double Histogram::xmax(void) const
718 26 Dec 06 jari 176   {
718 26 Dec 06 jari 177     return xmax_;
718 26 Dec 06 jari 178   }
718 26 Dec 06 jari 179
718 26 Dec 06 jari 180
718 26 Dec 06 jari 181   double Histogram::xmin(void) const
718 26 Dec 06 jari 182   {
718 26 Dec 06 jari 183     return xmin_;
718 26 Dec 06 jari 184   }
718 26 Dec 06 jari 185
718 26 Dec 06 jari 186
718 26 Dec 06 jari 187   double Histogram::operator[](size_t k) const
718 26 Dec 06 jari 188   {
718 26 Dec 06 jari 189     return histogram_[k];
718 26 Dec 06 jari 190   }
718 26 Dec 06 jari 191
718 26 Dec 06 jari 192
718 26 Dec 06 jari 193   const Histogram& Histogram::operator=(const Histogram& b)
718 26 Dec 06 jari 194   {
718 26 Dec 06 jari 195     if (this==&b)
718 26 Dec 06 jari 196       return *this;
718 26 Dec 06 jari 197     histogram_=b.histogram_;
718 26 Dec 06 jari 198     xmax_=b.xmax_;
718 26 Dec 06 jari 199     xmin_=b.xmin_;
718 26 Dec 06 jari 200     sum_all_=b.sum_all_;
718 26 Dec 06 jari 201     sum_histogram_=b.sum_histogram_;
718 26 Dec 06 jari 202     return *this;
718 26 Dec 06 jari 203   }
718 26 Dec 06 jari 204
718 26 Dec 06 jari 205
2032 19 Aug 09 peter 206   const Histogram& Histogram::operator+=(const Histogram& rhs)
2032 19 Aug 09 peter 207   {
2032 19 Aug 09 peter 208     assert(histogram_.size()==rhs.histogram_.size());
2032 19 Aug 09 peter 209     std::transform(histogram_.begin(), histogram_.end(),
2032 19 Aug 09 peter 210                    rhs.histogram_.begin(), histogram_.begin(),
2032 19 Aug 09 peter 211                    std::plus<double>());
2034 20 Aug 09 peter 212     assert(xmax_==rhs.xmax_);
2034 20 Aug 09 peter 213     assert(xmin_==rhs.xmin_);
2032 19 Aug 09 peter 214     sum_all_ += rhs.sum_all_;
2032 19 Aug 09 peter 215     sum_histogram_ += rhs.sum_histogram_;
2032 19 Aug 09 peter 216     return *this;
2032 19 Aug 09 peter 217   }
2032 19 Aug 09 peter 218
2032 19 Aug 09 peter 219
718 26 Dec 06 jari 220   std::ostream& operator<<(std::ostream& s,const Histogram& histogram)
718 26 Dec 06 jari 221   {
718 26 Dec 06 jari 222     s << "# histogram min : " << histogram.xmin() << '\n';
718 26 Dec 06 jari 223     s << "# histogram max : " << histogram.xmax() << '\n';
718 26 Dec 06 jari 224     s << "# number of bins: " << histogram.nof_bins() << '\n';
4200 19 Aug 22 peter 225     s << "# nof points in histogram : "
718 26 Dec 06 jari 226       << histogram.averager_histogram().sum_w() << '\n';
4200 19 Aug 22 peter 227     s << "# nof points in total:      "
718 26 Dec 06 jari 228       << histogram.averager_all().sum_w() << '\n';
718 26 Dec 06 jari 229     s << "# column 1: center of observation bin\n"
718 26 Dec 06 jari 230       << "# column 2: frequency\n";
718 26 Dec 06 jari 231
1271 09 Apr 08 peter 232     for (size_t i=0; i<histogram.nof_bins(); i++) {
718 26 Dec 06 jari 233       s.width(12);
718 26 Dec 06 jari 234       s << histogram.observation_value(i);
718 26 Dec 06 jari 235       s.width(12);
718 26 Dec 06 jari 236       s << histogram[i] << '\n';
718 26 Dec 06 jari 237     }
718 26 Dec 06 jari 238
718 26 Dec 06 jari 239     return s;
718 26 Dec 06 jari 240   }
718 26 Dec 06 jari 241
680 11 Oct 06 jari 242 }}} // of namespace statistics, yat, and theplu