source: trunk/src/WeightedAverager.h @ 220

// $Id: WeightedAverager.h 220 2004-12-30 22:31:34Z peter $

#ifndef _theplu_cpptools_weighted_averager_
#define _theplu_cpptools_weighted_averager_

//#include <cmath>
#include "Averager.h"

namespace theplu{
namespace statistics{
  ///
  /// Class to calulate simple (first and second moments) averages
  /// with weights.
  ///
  /// @see Averager
  ///
  class WeightedAverager
  {
  public:

    ///
    /// Default constructor
    ///
    WeightedAverager(void);

    ///
    /// Constructor taking the data point, i.e. the value and its
    /// weight (default = 1)
    ///
    WeightedAverager(const double, const double=1);

    ///
    /// Copy constructor
    ///
    WeightedAverager(const WeightedAverager&);

    ///
    /// adding a data point d, with weight w (default is 1)
    ///
    inline void add(const double d,const double w=1)
    {w_.add(w); wx_.add(w*d); wwx_+=w*w*d;}

    ///
    /// Calculating the weighted mean 
    ///
    /// @return \f$ \frac{\sum w_ix_i}{\sum w_i} \f$ 
    ///
    inline double mean(void) const { return sum_w() ? 
                                       sum_wx()/sum_w() : 0; }
  
    ///
    /// rescale object, i.e. each data point is rescaled
    /// \f$ x = a * x \f$
    ///
    inline void rescale(double a) { wx_.rescale(a); wwx_*=a; }

    ///
    /// resets everything to zero
    ///
    inline void reset(void) { wx_.reset(); w_.reset(); wwx_=0;}

    ///
    /// Calculating the squared error according to: \f$ \frac{\sum
    /// w_i^2(x_i-m)^2}{(\sum w_i)^2} \f$ @return squared error
    ///
    inline double squared_error(void) const 
    { return (squared_sum()) / (sum_w()*sum_w()); }

    ///
    /// The standard deviation is calculated as \f$ \frac{\sum w_i^2}{(\sum
    /// w_i)^2-\sum w_i^2}\frac{\sum w_i(x_i-m)^2}{\sum w_i} \f$
    ///
    /// @return standard deviation of mean
    ///
    inline double standard_error(void)  const 
    { return sqrt(squared_sum()/((sum_w()*sum_w())-squared_sum())
              *squared_error()/sum_w()); }

    ///
    /// Calculating the squared sum  
    ///
    /// @return \f$ \sum w_i^2(x_i-m)^2 \f$
    ///
    inline double squared_sum(void) const 
    { return (sum_wwxx()-2*mean()*wwx_+
              mean()*mean()*sum_ww()) ; }

    ///
    /// Calculating the sum of weights: \f$ \sum
    /// w_i \f$ @return sum of weights
    ///
    inline double sum_w(void) const 
    { return w_.sum_x(); }

    ///
    /// \f$ \sum w_ix_i \f$ @return weighted sum of x
    ///
    inline double sum_wx(void)  const 
    { return wx_.sum_x(); }


  private:
    inline double sum_ww(void)  const 
    { return w_.sum_xsqr(); }
    ///
    ///  @return \f$ \sum w_i^2x_i^2 \f$
    ///
    inline double sum_wwxx(void)  const 
    { return wx_.sum_xsqr(); }
    
    ///
    ///  @return \f$ \sum w_i^2x_i \f$
    ///
    inline double sum_wwx(void) const 
    { return wwx_; }

    ///
    /// operator to add a WeightedAverager
    ///
    inline WeightedAverager operator+=(WeightedAverager& a)
    { wx_+=a.wx(); w_+=a.w(); wwx_+=a.sum_wwx(); return *this; }
    
    Averager w_;
    Averager wx_;
    double wwx_;
    
    inline Averager wx(void) const {return wx_;}
    inline Averager w(void) const {return w_;}


  };

///
/// The WeightedAverager output operator
///
std::ostream& operator<<(std::ostream& s,const WeightedAverager&);

}} // of namespace statistics and namespace theplu

#endif