test/averager.cc

Code
Comments
Other
Rev Date Author Line
224 28 Jan 05 peter 1 // $Id$
224 28 Jan 05 peter 2
675 10 Oct 06 jari 3 /*
831 27 Mar 07 peter 4   Copyright (C) 2005 Peter Johansson
2121 13 Dec 09 peter 5   Copyright (C) 2006 Jari Häkkinen, Peter Johansson, Markus Ringnér
2121 13 Dec 09 peter 6   Copyright (C) 2007, 2008 Jari Häkkinen, Peter Johansson
4359 23 Aug 23 peter 7   Copyright (C) 2009, 2010, 2012, 2014 Peter Johansson
224 28 Jan 05 peter 8
1437 25 Aug 08 peter 9   This file is part of the yat library, http://dev.thep.lu.se/yat
675 10 Oct 06 jari 10
675 10 Oct 06 jari 11   The yat library is free software; you can redistribute it and/or
675 10 Oct 06 jari 12   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 13   published by the Free Software Foundation; either version 3 of the
675 10 Oct 06 jari 14   License, or (at your option) any later version.
675 10 Oct 06 jari 15
675 10 Oct 06 jari 16   The yat library is distributed in the hope that it will be useful,
675 10 Oct 06 jari 17   but WITHOUT ANY WARRANTY; without even the implied warranty of
675 10 Oct 06 jari 18   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
675 10 Oct 06 jari 19   General Public License for more details.
675 10 Oct 06 jari 20
675 10 Oct 06 jari 21   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 22   along with yat. If not, see <http://www.gnu.org/licenses/>.
675 10 Oct 06 jari 23 */
675 10 Oct 06 jari 24
2881 18 Nov 12 peter 25 #include <config.h>
2881 18 Nov 12 peter 26
1214 07 Mar 08 peter 27 #include "Suite.h"
1214 07 Mar 08 peter 28
675 10 Oct 06 jari 29 #include "yat/statistics/Averager.h"
675 10 Oct 06 jari 30 #include "yat/statistics/AveragerPair.h"
675 10 Oct 06 jari 31 #include "yat/statistics/AveragerPairWeighted.h"
675 10 Oct 06 jari 32 #include "yat/statistics/AveragerWeighted.h"
1120 21 Feb 08 peter 33 #include "yat/utility/Vector.h"
675 10 Oct 06 jari 34
2159 19 Jan 10 peter 35 #include <boost/concept_archetype.hpp>
2159 19 Jan 10 peter 36
848 28 Apr 07 jari 37 #include <cmath>
401 18 Oct 05 peter 38 #include <fstream>
499 18 Jan 06 peter 39 #include <limits>
401 18 Oct 05 peter 40 #include <iostream>
887 24 Sep 07 peter 41 #include <set>
224 28 Jan 05 peter 42
1214 07 Mar 08 peter 43 using namespace theplu::yat;
680 11 Oct 06 jari 44 using namespace theplu::yat::statistics;
484 03 Jan 06 peter 45
1214 07 Mar 08 peter 46 using theplu::yat::test::Suite;
1214 07 Mar 08 peter 47
490 04 Jan 06 peter 48 //Forward declarations
1273 10 Apr 08 jari 49 bool equal(const Averager&, const Averager&, unsigned int, Suite&);
4200 19 Aug 22 peter 50 bool equal(const AveragerWeighted&, const AveragerWeighted&, unsigned int,
1214 07 Mar 08 peter 51            Suite& suite);
1273 10 Apr 08 jari 52 bool equal(const Averager&, const AveragerWeighted&, unsigned int, Suite& suite);
1273 10 Apr 08 jari 53 bool equal(const AveragerPair&, const AveragerPair&, unsigned int, Suite& suite);
4200 19 Aug 22 peter 54 bool equal(const AveragerPair&, const AveragerPairWeighted&, unsigned int,
1214 07 Mar 08 peter 55            Suite& suite);
1273 10 Apr 08 jari 56 bool equal(const AveragerPairWeighted&, const AveragerPairWeighted&,
1273 10 Apr 08 jari 57            unsigned int, Suite& suite);
484 03 Jan 06 peter 58
2798 27 Jul 12 peter 59 void test_plus_assign(test::Suite& suite);
2809 06 Aug 12 peter 60 void test_rescale(test::Suite& suite);
486 04 Jan 06 peter 61
1214 07 Mar 08 peter 62 int main(int argc, char* argv[])
2809 06 Aug 12 peter 63 {
484 03 Jan 06 peter 64
1214 07 Mar 08 peter 65   Suite suite(argc, argv);
399 14 Oct 05 peter 66
486 04 Jan 06 peter 67   // Testing Averager
1214 07 Mar 08 peter 68   suite.out() << "testing Averager" << std::endl;
484 03 Jan 06 peter 69   Averager a;
484 03 Jan 06 peter 70   a.add(1);
484 03 Jan 06 peter 71   a.add(3);
484 03 Jan 06 peter 72   a.add(5);
484 03 Jan 06 peter 73   if (a.n()!=3 || a.mean()!=3 || a.sum_xx()!=35){
1232 15 Mar 08 peter 74     suite.add(false);
1214 07 Mar 08 peter 75     suite.err() << "error: add\n";
484 03 Jan 06 peter 76   }
1295 12 May 08 jari 77   Averager b(a);
1295 12 May 08 jari 78   b.add(5,-1);
1295 12 May 08 jari 79   if (b.n()!=2 || b.mean()!=2 || b.sum_xx()!=10){
1295 12 May 08 jari 80     suite.add(false);
1295 12 May 08 jari 81     suite.err() << "error: add with negative n\n";
1295 12 May 08 jari 82   }
484 03 Jan 06 peter 83
620 04 Sep 06 markus 84   Averager* a1 = new Averager(1.0+3+5,1.0+9+25,3);
1273 10 Apr 08 jari 85   unsigned int tol = 10;
1214 07 Mar 08 peter 86   if (!equal(a,*a1, tol, suite)){
1232 15 Mar 08 peter 87     suite.add(false);
1214 07 Mar 08 peter 88     suite.err().precision(25);
1214 07 Mar 08 peter 89     suite.err() << a.sum_x() << '\t' << a1->sum_x() << std::endl;
1214 07 Mar 08 peter 90     suite.err() << a.sum_xx() << '\t' << a1->sum_xx() << std::endl;
1214 07 Mar 08 peter 91     suite.err() << a.n() << '\t' << a1->n() << std::endl;
1214 07 Mar 08 peter 92     suite.err() << a.variance() << '\t' << a1->variance() << std::endl;
1214 07 Mar 08 peter 93     suite.err() << a.mean() << '\t' << a1->mean() << std::endl;
1214 07 Mar 08 peter 94     suite.err() << a.mean() - a1->mean() << std::endl;
1214 07 Mar 08 peter 95     suite.err() << a.variance() - a1->variance() << std::endl;
1296 13 May 08 jari 96     suite.err() << "error: Averager(double x, double xx, long n)\n";
484 03 Jan 06 peter 97   }
484 03 Jan 06 peter 98   delete a1;
484 03 Jan 06 peter 99
484 03 Jan 06 peter 100   a1 = new Averager(a);
1214 07 Mar 08 peter 101   if (!equal(a,*a1, tol, suite)){
1232 15 Mar 08 peter 102     suite.add(false);
1214 07 Mar 08 peter 103     suite.err() << "error: Copy constructor\n";
484 03 Jan 06 peter 104   }
484 03 Jan 06 peter 105   delete a1;
484 03 Jan 06 peter 106
484 03 Jan 06 peter 107   a.add(3,5);
1704 08 Jan 09 peter 108   if (!suite.equal_sqrt(a.standard_error(), sqrt(a.variance()/a.n()),1)) {
1232 15 Mar 08 peter 109     suite.add(false);
1214 07 Mar 08 peter 110     suite.err() << "error: standard_error\n";
484 03 Jan 06 peter 111   }
2809 06 Aug 12 peter 112
1664 20 Dec 08 peter 113   if ( !suite.equal(a.variance(), a.std()*a.std()) ) {
1232 15 Mar 08 peter 114     suite.add(false);
1214 07 Mar 08 peter 115     suite.err() << "error: std squared should be variance" << std::endl;
1214 07 Mar 08 peter 116     suite.err() << "std2: " << a.std()*a.std() << std::endl;
1214 07 Mar 08 peter 117     suite.err() << "variance: " << a.variance() << std::endl;
2809 06 Aug 12 peter 118     suite.err() << "difference is: " << a.std()*a.std()-a.variance()
2809 06 Aug 12 peter 119                 << std::endl;
484 03 Jan 06 peter 120   }
2809 06 Aug 12 peter 121
484 03 Jan 06 peter 122   if ( a.variance() != a.variance(a.mean()) ){
1232 15 Mar 08 peter 123     suite.add(false);
1214 07 Mar 08 peter 124     suite.err() << "error: variance incorrect\n" << std::endl;
1214 07 Mar 08 peter 125     suite.err() << "variance: " << a.variance() << std::endl;
1214 07 Mar 08 peter 126     suite.err() << "mean: " << a.mean() << std::endl;
1214 07 Mar 08 peter 127     suite.err() << "variance(mean) " << a.variance(a.mean()) << std::endl;
484 03 Jan 06 peter 128   }
1120 21 Feb 08 peter 129   theplu::yat::utility::Vector* tmp_vec = new theplu::yat::utility::Vector(10);
887 24 Sep 07 peter 130   add(a, tmp_vec->begin(), tmp_vec->end());
887 24 Sep 07 peter 131   delete tmp_vec;
887 24 Sep 07 peter 132   std::set<double>* tmp_set = new std::set<double>;
887 24 Sep 07 peter 133   tmp_set->insert(1.0);
887 24 Sep 07 peter 134   tmp_set->insert(2.0);
887 24 Sep 07 peter 135   add(a, tmp_set->begin(), tmp_set->end());
887 24 Sep 07 peter 136   delete tmp_set;
484 03 Jan 06 peter 137
887 24 Sep 07 peter 138
486 04 Jan 06 peter 139   // Testing AveragerWeighted
1214 07 Mar 08 peter 140   suite.err() << "testing AveragerWeighted" << std::endl;
1120 21 Feb 08 peter 141   theplu::yat::utility::Vector x(3,0);
486 04 Jan 06 peter 142   x(0)=0;
486 04 Jan 06 peter 143   x(1)=1;
486 04 Jan 06 peter 144   x(2)=2;
1120 21 Feb 08 peter 145   theplu::yat::utility::Vector w(3,1);
680 11 Oct 06 jari 146   theplu::yat::statistics::AveragerWeighted aw;
3266 03 Jul 14 peter 147   if (!suite.add(aw.n()==0))
3266 03 Jul 14 peter 148     suite.err() << "error: AveragerWeighted::n() returns "
3266 03 Jul 14 peter 149                 << aw.n() << " expected 0\n";
916 30 Sep 07 peter 150   add(aw, x.begin(), x.end(), w.begin());
486 04 Jan 06 peter 151   a.reset();
916 30 Sep 07 peter 152   add(a, x.begin(), x.end());
1214 07 Mar 08 peter 153   if (!equal(a,aw,tol,suite)){
2809 06 Aug 12 peter 154     suite.err() << "error: AveragerWeighted with unitary weights should "
2809 06 Aug 12 peter 155                 << "be equal to Averager" << std::endl;
1232 15 Mar 08 peter 156     suite.add(false);
486 04 Jan 06 peter 157   }
486 04 Jan 06 peter 158
486 04 Jan 06 peter 159   AveragerWeighted* aw2 = new AveragerWeighted(aw);
1214 07 Mar 08 peter 160   if (!equal(aw,*aw2,tol,suite)){
1214 07 Mar 08 peter 161     suite.err() << "error: AveragerWeighted copy constructor " << std::endl;
1232 15 Mar 08 peter 162     suite.add(false);
486 04 Jan 06 peter 163   }
2809 06 Aug 12 peter 164
486 04 Jan 06 peter 165   aw2->add(12,0);
1214 07 Mar 08 peter 166   if (!equal(aw,*aw2,tol,suite)){
2809 06 Aug 12 peter 167     suite.err() << "error: AveragerWeighted adding a data point with weight=0 "
2809 06 Aug 12 peter 168                 << "should make no change " << std::endl;
1232 15 Mar 08 peter 169     suite.add(false);
486 04 Jan 06 peter 170   }
2809 06 Aug 12 peter 171
486 04 Jan 06 peter 172   aw2->reset();
775 01 Mar 07 jari 173   w*=17;
916 30 Sep 07 peter 174   add(*aw2, x.begin(), x.end(), w.begin());
1214 07 Mar 08 peter 175   if (!equal(aw,*aw2,tol,suite)){
2809 06 Aug 12 peter 176     suite.err() << "error: AveragerWeighted rescaling weights "
486 04 Jan 06 peter 177            << "should make no change " << std::endl;
1232 15 Mar 08 peter 178     suite.add(false);
486 04 Jan 06 peter 179   }
486 04 Jan 06 peter 180   delete aw2;
887 24 Sep 07 peter 181   {
1120 21 Feb 08 peter 182     theplu::yat::utility::Vector tmp(10);
887 24 Sep 07 peter 183     add(aw, tmp.begin(), tmp.end());
887 24 Sep 07 peter 184   }
887 24 Sep 07 peter 185   {
887 24 Sep 07 peter 186     std::set<double> tmp;
887 24 Sep 07 peter 187     tmp.insert(1.0);
887 24 Sep 07 peter 188     tmp.insert(2.0);
887 24 Sep 07 peter 189     add(aw, tmp.begin(), tmp.end());
887 24 Sep 07 peter 190   }
486 04 Jan 06 peter 191
2809 06 Aug 12 peter 192
1214 07 Mar 08 peter 193   suite.out() << "testing AveragerPair" << std::endl;
487 04 Jan 06 peter 194   AveragerPair ap;
224 28 Jan 05 peter 195   for (int i=0; i<10; i++)
583 09 May 06 markus 196     ap.add(static_cast<double>(i),i);
1704 08 Jan 09 peter 197   if (!suite.equal(ap.correlation(),1,tol)){
1232 15 Mar 08 peter 198     suite.add(false);
1214 07 Mar 08 peter 199     suite.err() << "correlation: " << ap.correlation() << std::endl;
2809 06 Aug 12 peter 200     suite.err() << "error: correlation between identical vectors should be unity"
487 04 Jan 06 peter 201            << std::endl;
487 04 Jan 06 peter 202   }
1704 08 Jan 09 peter 203   if (!suite.equal(ap.x_averager().variance(),ap.covariance())) {
1232 15 Mar 08 peter 204     suite.add(false);
2809 06 Aug 12 peter 205     suite.err() << "error: covariance of identical vectors should equal to variance"
487 04 Jan 06 peter 206            << std::endl;
487 04 Jan 06 peter 207   }
487 04 Jan 06 peter 208   AveragerPair* ap2 = new AveragerPair(ap);
487 04 Jan 06 peter 209   delete ap2;
224 28 Jan 05 peter 210
1295 12 May 08 jari 211   for (int i=0; i<8; i++)
1295 12 May 08 jari 212     ap.add(static_cast<double>(i),i,-1);
1704 08 Jan 09 peter 213   if (!suite.equal(ap.correlation(),1, tol)) {
1295 12 May 08 jari 214     suite.add(false);
1295 12 May 08 jari 215     suite.err() << "correlation after removal of data: " << ap.correlation()
1295 12 May 08 jari 216                 << std::endl;
2809 06 Aug 12 peter 217     suite.err() << "error: correlation between identical vectors is unity"
1295 12 May 08 jari 218                 << std::endl;
1295 12 May 08 jari 219   }
2202 21 Feb 10 peter 220   // do not run compile test
2202 21 Feb 10 peter 221   if (false) {
2202 21 Feb 10 peter 222     AveragerPair ap3;
2202 21 Feb 10 peter 223     add(ap3, boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 224         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 225         boost::input_iterator_archetype<double>());
2202 21 Feb 10 peter 226     AveragerPairWeighted apw;
2202 21 Feb 10 peter 227     add(apw, boost::input_iterator_archetype_no_proxy<utility::DataWeight>(),
2202 21 Feb 10 peter 228         boost::input_iterator_archetype_no_proxy<utility::DataWeight>(),
2202 21 Feb 10 peter 229         boost::input_iterator_archetype_no_proxy<utility::DataWeight>());
2202 21 Feb 10 peter 230     add(apw, boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 231         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 232         boost::input_iterator_archetype<double>());
2202 21 Feb 10 peter 233     add(apw, boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 234         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 235         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 236         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 237         boost::input_iterator_archetype<double>());
2202 21 Feb 10 peter 238   }
1295 12 May 08 jari 239
1295 12 May 08 jari 240
1214 07 Mar 08 peter 241   suite.out() << "testing AveragerPairWeighted" << std::endl;
490 04 Jan 06 peter 242   AveragerPairWeighted apw;
490 04 Jan 06 peter 243   x(0)=0; x(1)=1; x(2)=2;
1120 21 Feb 08 peter 244   theplu::yat::utility::Vector y(3,0);
1224 11 Mar 08 peter 245   y(0)=0; y(1)=0; y(2)=2;
1043 06 Feb 08 peter 246   add(apw, x.begin(), x.end(), y.begin(), w.begin(), w.begin());
490 04 Jan 06 peter 247   ap.reset();
1043 06 Feb 08 peter 248   add(ap, x.begin(), x.end(), y.begin());
1214 07 Mar 08 peter 249   if (!equal(ap,apw,tol,suite)){
2809 06 Aug 12 peter 250     suite.err() << "error: AveragerPairWeighted with unitary weights should "
490 04 Jan 06 peter 251            << "be equal to AveragerPair" << std::endl;
1232 15 Mar 08 peter 252     suite.add(false);
490 04 Jan 06 peter 253   }
490 04 Jan 06 peter 254
490 04 Jan 06 peter 255   AveragerPairWeighted* apw2 = new AveragerPairWeighted(apw);
1214 07 Mar 08 peter 256   if (!equal(apw,*apw2,tol,suite)){
1214 07 Mar 08 peter 257     suite.err() << "error: AveragerPairWeighted copy constructor " << std::endl;
1232 15 Mar 08 peter 258     suite.add(false);
490 04 Jan 06 peter 259   }
2809 06 Aug 12 peter 260
490 04 Jan 06 peter 261   apw2->add(12,23222.03,32.3,0);
1214 07 Mar 08 peter 262   if (!equal(apw,*apw2,tol,suite)){
4200 19 Aug 22 peter 263     suite.err() << "error: AveragerWeighted adding a data point with weight=0 "
490 04 Jan 06 peter 264            << "should make no change " << std::endl;
1232 15 Mar 08 peter 265     suite.add(false);
490 04 Jan 06 peter 266   }
2809 06 Aug 12 peter 267
490 04 Jan 06 peter 268   apw2->reset();
775 01 Mar 07 jari 269   w*=17;
1043 06 Feb 08 peter 270   add(*apw2, x.begin(), x.end(), y.begin(), w.begin(), w.begin());
1214 07 Mar 08 peter 271   if (!equal(apw,*apw2,tol,suite)){
2809 06 Aug 12 peter 272     suite.err() << "error: AveragerWeighted rescaling weights "
490 04 Jan 06 peter 273            << "should make no change " << std::endl;
1232 15 Mar 08 peter 274     suite.add(false);
490 04 Jan 06 peter 275   }
490 04 Jan 06 peter 276   delete apw2;
490 04 Jan 06 peter 277
2798 27 Jul 12 peter 278   test_plus_assign(suite);
2809 06 Aug 12 peter 279   test_rescale(suite);
2798 27 Jul 12 peter 280
1251 03 Apr 08 peter 281   return suite.return_value();
224 28 Jan 05 peter 282 }
224 28 Jan 05 peter 283
4200 19 Aug 22 peter 284 bool equal(const Averager& a, const Averager& b, unsigned int tol,
1214 07 Mar 08 peter 285            Suite& suite)
490 04 Jan 06 peter 286 {
1234 15 Mar 08 peter 287   return (a.n()==b.n() && suite.equal(a.mean(),b.mean(),tol) &&
1234 15 Mar 08 peter 288           suite.equal(a.variance(),b.variance(),tol));
490 04 Jan 06 peter 289 }
224 28 Jan 05 peter 290
4200 19 Aug 22 peter 291 bool equal(const AveragerWeighted& a, const AveragerWeighted& b,
1273 10 Apr 08 jari 292            const unsigned int tol, Suite& suite)
490 04 Jan 06 peter 293 {
490 04 Jan 06 peter 294   bool equal = true;
1234 15 Mar 08 peter 295   if  ( !suite.equal(a.mean(),b.mean(),tol)){
490 04 Jan 06 peter 296     equal=false;
1214 07 Mar 08 peter 297     suite.err() << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
1214 07 Mar 08 peter 298     suite.err() << "difference:\t" << a.mean()-b.mean() << std::endl;
490 04 Jan 06 peter 299   }
1234 15 Mar 08 peter 300   if ( !suite.equal(a.variance(),b.variance(),tol) ) {
490 04 Jan 06 peter 301     equal=false;
4200 19 Aug 22 peter 302     suite.err() << "error for variance:\t" << a.variance() << " "
1214 07 Mar 08 peter 303                 << b.variance() << std::endl;
4200 19 Aug 22 peter 304   }
1234 15 Mar 08 peter 305   if ( !suite.equal(a.standard_error(),b.standard_error(),tol) ) {
490 04 Jan 06 peter 306     equal =false;
1214 07 Mar 08 peter 307     suite.err() << "error for standard error:\t" << std::endl;
490 04 Jan 06 peter 308   }
490 04 Jan 06 peter 309   return equal;
490 04 Jan 06 peter 310 }
301 30 Apr 05 peter 311
4200 19 Aug 22 peter 312 bool equal(const Averager& a, const AveragerWeighted& b, const unsigned int tol,
1214 07 Mar 08 peter 313            Suite& suite)
490 04 Jan 06 peter 314 {
490 04 Jan 06 peter 315   bool equal = true;
1234 15 Mar 08 peter 316   if  ( !suite.equal(a.mean(),b.mean(),tol)){
490 04 Jan 06 peter 317     equal=false;
1214 07 Mar 08 peter 318     suite.err() << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
490 04 Jan 06 peter 319   }
1234 15 Mar 08 peter 320   if ( !suite.equal(a.variance(),b.variance(),tol) ) {
490 04 Jan 06 peter 321     equal=false;
4200 19 Aug 22 peter 322     suite.err() << "error for variance:\t" << a.variance() << " "
1214 07 Mar 08 peter 323                 << b.variance() << std::endl;
4200 19 Aug 22 peter 324   }
1234 15 Mar 08 peter 325   if ( !suite.equal(a.standard_error(), b.standard_error(),tol) ) {
490 04 Jan 06 peter 326     equal =false;
1214 07 Mar 08 peter 327     suite.err() << "error for standard error:\t" << std::endl;
490 04 Jan 06 peter 328   }
490 04 Jan 06 peter 329   return equal;
490 04 Jan 06 peter 330 }
490 04 Jan 06 peter 331
4200 19 Aug 22 peter 332 bool equal(const AveragerPair& a, const AveragerPair& b,
1273 10 Apr 08 jari 333            const unsigned int tol, Suite& suite)
490 04 Jan 06 peter 334 {
490 04 Jan 06 peter 335   bool ok = true;
1704 08 Jan 09 peter 336   if (!suite.equal(a.covariance(), b.covariance(), tol)) {
1232 15 Mar 08 peter 337     suite.add(false);
4200 19 Aug 22 peter 338     suite.err() << "error covariance: " << a.covariance() << "\t"
490 04 Jan 06 peter 339            << b.covariance() << std::endl;
490 04 Jan 06 peter 340   }
1704 08 Jan 09 peter 341   if ( !suite.equal(a.correlation(),b.correlation(), tol)) {
1232 15 Mar 08 peter 342     suite.add(false);
1214 07 Mar 08 peter 343     suite.err() << "error correlation" << std::endl;
4200 19 Aug 22 peter 344   }
490 04 Jan 06 peter 345   return ok;
490 04 Jan 06 peter 346 }
490 04 Jan 06 peter 347
4200 19 Aug 22 peter 348 bool equal(const AveragerPair& a, const AveragerPairWeighted& b,
1273 10 Apr 08 jari 349            const unsigned int tol, Suite& suite)
490 04 Jan 06 peter 350 {
490 04 Jan 06 peter 351   bool ok = true;
1704 08 Jan 09 peter 352   if (!suite.equal(a.covariance(), b.covariance(), tol) ) {
1232 15 Mar 08 peter 353     suite.add(false);
4200 19 Aug 22 peter 354     suite.err() << "error covariance: " << a.covariance() << "\t"
490 04 Jan 06 peter 355            << b.covariance() << std::endl;
490 04 Jan 06 peter 356   }
1704 08 Jan 09 peter 357   if ( !suite.equal(a.correlation(),b.correlation(), tol) ) {
1232 15 Mar 08 peter 358     suite.add(false);
1214 07 Mar 08 peter 359     suite.err() << "error correlation" << std::endl;
1224 11 Mar 08 peter 360     suite.err() << "unweighted:" << a.correlation() << std::endl;
1224 11 Mar 08 peter 361     suite.err() << "weighted:" << b.correlation() << std::endl;
4200 19 Aug 22 peter 362     suite.err() << "difference:" << a.correlation()-b.correlation()
1224 11 Mar 08 peter 363                 << std::endl;
4200 19 Aug 22 peter 364   }
1214 07 Mar 08 peter 365   if ( !equal(a.x_averager(),b.x_averager(),tol,suite)) {
490 04 Jan 06 peter 366     ok =false;
1214 07 Mar 08 peter 367     suite.err() << "error for x_averager():\t" << std::endl;
490 04 Jan 06 peter 368   }
490 04 Jan 06 peter 369   return ok;
490 04 Jan 06 peter 370 }
4200 19 Aug 22 peter 371 bool equal(const AveragerPairWeighted& a, const AveragerPairWeighted& b,
1273 10 Apr 08 jari 372            const unsigned int tol, Suite& suite)
490 04 Jan 06 peter 373 {
490 04 Jan 06 peter 374   bool ok = true;
1234 15 Mar 08 peter 375   if  ( !suite.equal(a.covariance(),b.covariance(),tol) ){
1232 15 Mar 08 peter 376     suite.add(false);
4200 19 Aug 22 peter 377     suite.err() << "error covariance: " << a.covariance() << "\t"
490 04 Jan 06 peter 378            << b.covariance() << std::endl;
490 04 Jan 06 peter 379   }
1234 15 Mar 08 peter 380   if ( !suite.equal(a.correlation(), b.correlation(), tol) ) {
1232 15 Mar 08 peter 381     suite.add(false);
1214 07 Mar 08 peter 382     suite.err() << "error correlation" << std::endl;
1224 11 Mar 08 peter 383     suite.err() << "a:" << a.correlation() << std::endl;
1224 11 Mar 08 peter 384     suite.err() << "b:" << b.correlation() << std::endl;
4200 19 Aug 22 peter 385     suite.err() << "difference:" << a.correlation()-b.correlation()
1224 11 Mar 08 peter 386                 << std::endl;
1224 11 Mar 08 peter 387     suite.err() << "tol:" << tol << std::endl;
4200 19 Aug 22 peter 388   }
1214 07 Mar 08 peter 389   if ( !equal(a.x_averager(),b.x_averager(),tol,suite)) {
490 04 Jan 06 peter 390     ok =false;
1214 07 Mar 08 peter 391     suite.err() << "error for x_averager():\t" << std::endl;
490 04 Jan 06 peter 392   }
2159 19 Jan 10 peter 393
2159 19 Jan 10 peter 394   // compile tests should not be run
2159 19 Jan 10 peter 395   if (false) {
2159 19 Jan 10 peter 396     Averager averager;
2159 19 Jan 10 peter 397     add(averager, boost::input_iterator_archetype<double>(),
2159 19 Jan 10 peter 398         boost::input_iterator_archetype<double>());
2159 19 Jan 10 peter 399     AveragerWeighted averagerw;
2159 19 Jan 10 peter 400     using utility::DataWeight;
2202 21 Feb 10 peter 401     typedef boost::input_iterator_archetype_no_proxy<DataWeight> Iter;
2159 19 Jan 10 peter 402     Iter iter;
2159 19 Jan 10 peter 403     utility::iterator_traits<Iter> traits;
2202 21 Feb 10 peter 404     traits.weight(iter);
2202 21 Feb 10 peter 405     add(averagerw, iter, iter);
2202 21 Feb 10 peter 406     add(averagerw, boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 407         boost::input_iterator_archetype<double>(),
2202 21 Feb 10 peter 408         boost::input_iterator_archetype<double>());
2159 19 Jan 10 peter 409   }
2159 19 Jan 10 peter 410
490 04 Jan 06 peter 411   return ok;
490 04 Jan 06 peter 412 }
2798 27 Jul 12 peter 413
2798 27 Jul 12 peter 414
2798 27 Jul 12 peter 415 void test_plus_assign(test::Suite& suite)
2798 27 Jul 12 peter 416 {
2798 27 Jul 12 peter 417   suite.out() << "Testing plus assign\n";
2798 27 Jul 12 peter 418   Averager a1;
2798 27 Jul 12 peter 419   a1.add(1);
2798 27 Jul 12 peter 420   a1.add(3);
2798 27 Jul 12 peter 421   Averager a2;
2798 27 Jul 12 peter 422   a2.add(-1);
2798 27 Jul 12 peter 423   a2.add(-3);
2798 27 Jul 12 peter 424   if (!suite.add(suite.equal(a2.mean(), (-1-3)/2)))
2798 27 Jul 12 peter 425     suite.err() << "mean failed\n";
2798 27 Jul 12 peter 426   if (!suite.add(suite.equal(a2.variance(), 1)))
2798 27 Jul 12 peter 427     suite.err() << "variance failed\n";
2798 27 Jul 12 peter 428   a2+=a1;
2798 27 Jul 12 peter 429   if (!suite.add(suite.equal(a2.mean(), (1+3-1-3)/4)))
2798 27 Jul 12 peter 430     suite.err() << "mean failed\n";
2798 27 Jul 12 peter 431   if (!suite.add(suite.equal(a2.variance(), (1+9+1+9)/4)))
2798 27 Jul 12 peter 432     suite.err() << "variance failed\n";
2798 27 Jul 12 peter 433
2798 27 Jul 12 peter 434   suite.out() << "Testing plus assign weighted\n";
2798 27 Jul 12 peter 435   AveragerWeighted aw1;
2798 27 Jul 12 peter 436   aw1.add(1);
2798 27 Jul 12 peter 437   aw1.add(3);
2798 27 Jul 12 peter 438   AveragerWeighted aw2;
2798 27 Jul 12 peter 439   aw2.add(-1);
2798 27 Jul 12 peter 440   aw2.add(-3,0.5);
2798 27 Jul 12 peter 441   aw2.add(-3,0.5);
2798 27 Jul 12 peter 442   aw2.add(42,0.0);
2798 27 Jul 12 peter 443   aw1 += aw2;
2798 27 Jul 12 peter 444   if (!suite.add(suite.equal(aw1.mean(), (1+3-1-3)/4)))
2798 27 Jul 12 peter 445     suite.err() << "mean failed\n";
2798 27 Jul 12 peter 446   if (!suite.add(suite.equal(aw1.variance(), (1+9+1+9)/4)))
2798 27 Jul 12 peter 447     suite.err() << "variance failed\n";
2798 27 Jul 12 peter 448
2798 27 Jul 12 peter 449   suite.out() << "Testing plus assign AveragerPair\n";
2798 27 Jul 12 peter 450   AveragerPair ap1;
2798 27 Jul 12 peter 451   ap1.add(3,3);
2798 27 Jul 12 peter 452   ap1.add(5,5);
2798 27 Jul 12 peter 453   if (!suite.add(suite.equal(ap1.sum_xy_centered(), 2)))
2798 27 Jul 12 peter 454     suite.err() << "sum_xy_centered failed\n";
2798 27 Jul 12 peter 455   AveragerPair ap2;
2798 27 Jul 12 peter 456   ap2.add(-3,-5);
2798 27 Jul 12 peter 457   ap2.add(-5,-3);
2798 27 Jul 12 peter 458   ap1 += ap2;
2798 27 Jul 12 peter 459   if (!suite.add(suite.equal(ap1.sum_xy_centered(), 9+25+15+15)))
2798 27 Jul 12 peter 460     suite.err() << "sum_xy_centered failed\n";
2798 27 Jul 12 peter 461
2798 27 Jul 12 peter 462   suite.out() << "Testing plus assign AveragerPairWeighted\n";
2798 27 Jul 12 peter 463   AveragerPairWeighted apw1;
2798 27 Jul 12 peter 464   apw1.add(3, 3, 1, 1);
2798 27 Jul 12 peter 465   apw1.add(5, 5, 1, 1);
2798 27 Jul 12 peter 466   if (!suite.add(suite.equal(apw1.sum_xy_centered(), 2)))
2798 27 Jul 12 peter 467     suite.err() << "sum_xy_centered failed\n";
2798 27 Jul 12 peter 468   AveragerPairWeighted apw2;
2798 27 Jul 12 peter 469   apw2.add(1, 3, 1, 1);
2798 27 Jul 12 peter 470   apw2.add(3, 1, 0.5, 1);
2798 27 Jul 12 peter 471   apw2.add(3, 1, 1, 0.5);
2798 27 Jul 12 peter 472   apw2.add(42, 99, 0, 0.5);
2798 27 Jul 12 peter 473   apw1 += apw2;
2798 27 Jul 12 peter 474   if (!suite.add(suite.equal(apw1.x_averager().mean(), 3)))
2798 27 Jul 12 peter 475     suite.err() << "x_averager().mean() failed\n";
2798 27 Jul 12 peter 476   if (!suite.add(suite.equal(apw1.y_averager().mean(), 3)))
2798 27 Jul 12 peter 477     suite.err() << "y_averager().mean() failed\n";
2798 27 Jul 12 peter 478   if (!suite.add(suite.equal(apw1.sum_xy(), 9+25+3+3)))
2798 27 Jul 12 peter 479     suite.err() << "sum_xy failed\n";
2798 27 Jul 12 peter 480   if (!suite.add(suite.equal(apw1.sum_xy_centered(), 0+4+0+0)))
2798 27 Jul 12 peter 481     suite.err() << "sum_xy_centered failed\n";
2798 27 Jul 12 peter 482 }
2809 06 Aug 12 peter 483
2809 06 Aug 12 peter 484
2809 06 Aug 12 peter 485 void test_rescale(test::Suite& suite)
2809 06 Aug 12 peter 486 {
2809 06 Aug 12 peter 487   Averager a;
2809 06 Aug 12 peter 488   a.add(1,1);
2809 06 Aug 12 peter 489   a.rescale(3.14);
2809 06 Aug 12 peter 490   suite.add(a.n()==1);
2809 06 Aug 12 peter 491   suite.add(suite.equal(a.mean(),3.14));
2809 06 Aug 12 peter 492 }