source: branches/0.3.1/test/averager_test.cc @ 840

Last change on this file since 840 was 840, checked in by Jari Häkkinen, 15 years ago

Addresses #222. Comparing doubles with == is poor programming style.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 9.3 KB
Line 
1// $Id: averager_test.cc 840 2007-04-25 07:55:37Z jari $
2
3/*
4  Copyright (C) 2005 Peter Johansson
5  Copyright (C) 2006 Jari Häkkinen, Markus Ringnér, Peter Johansson
6  Copyright (C) 2007 Jari Häkkinen
7
8  This file is part of the yat library, http://lev.thep.lu.se/trac/yat
9
10  The yat library is free software; you can redistribute it and/or
11  modify it under the terms of the GNU General Public License as
12  published by the Free Software Foundation; either version 2 of the
13  License, or (at your option) any later version.
14
15  The yat library is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18  General Public License for more details.
19
20  You should have received a copy of the GNU General Public License
21  along with this program; if not, write to the Free Software
22  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
23  02111-1307, USA.
24*/
25
26#include "yat/statistics/Averager.h"
27#include "yat/statistics/AveragerPair.h"
28#include "yat/statistics/AveragerPairWeighted.h"
29#include "yat/statistics/AveragerWeighted.h"
30#include "yat/utility/vector.h"
31
32#include <cmath>
33#include <fstream>
34#include <limits>
35#include <iostream>
36
37using namespace theplu::yat::statistics;
38
39//Forward declarations
40bool equal(const Averager&, const Averager&);
41bool equal(const AveragerWeighted&, const AveragerWeighted&, const double, 
42           std::ostream* error);
43bool equal(const Averager&, const AveragerWeighted&, const double, 
44           std::ostream* error);
45bool equal(const AveragerPair&, const AveragerPair&, 
46           const double, std::ostream* error);
47bool equal(const AveragerPair&, const AveragerPairWeighted&, 
48           const double, std::ostream* error);
49bool equal(const AveragerPairWeighted&, const AveragerPairWeighted&, 
50           const double, std::ostream* error);
51
52
53int main(const int argc,const char* argv[])
54{ 
55
56  std::ostream* error;
57  if (argc>1 && argv[1]==std::string("-v"))
58    error = &std::cerr;
59  else {
60    error = new std::ofstream("/dev/null");
61    if (argc>1)
62      std::cout << "averager_test -v : for printing extra information\n";
63  }
64  bool ok = true;
65
66  // Testing Averager
67  *error << "testing Averager" << std::endl;
68  Averager a;
69  a.add(1);
70  a.add(3);
71  a.add(5);
72  if (a.n()!=3 || a.mean()!=3 || a.sum_xx()!=35){
73    ok=false;
74    *error << "error: add\n";
75  }
76
77
78  Averager* a1 = new Averager(1.0+3+5,1.0+9+25,3);
79  if (!equal(a,*a1)){
80    ok=false;
81    std::cout.precision(25);
82    std::cout << equal(a,*a1) << std::endl;
83    std::cout << a.sum_x() << '\t' << a1->sum_x() << std::endl;
84    std::cout << a.sum_xx() << '\t' << a1->sum_xx() << std::endl;
85    std::cout << a.n() << '\t' << a1->n() << std::endl;
86    std::cout << a.variance() << '\t' << a1->variance() << std::endl;
87    std::cout << a.mean() << '\t' << a1->mean() << std::endl;
88    std::cout << a.mean() - a1->mean() << std::endl;
89    std::cout << a.variance() - a1->variance() << std::endl;
90    *error << "error: Averager(double x, double xx, u_long n)\n";
91  }
92  delete a1;
93
94  a1 = new Averager(a);
95  if (!equal(a,*a1)){
96    ok=false;
97    *error << "error: Copy constructor\n";
98  }
99  delete a1;
100
101  a.add(3,5);
102  if (! a.standard_error()==sqrt(a.variance()/a.n())){
103    ok=false;
104    *error << "error: standard_error\n";
105  }
106     
107
108  if ( fabs(a.variance() - a.std()*a.std())>
109       std::numeric_limits<double>().round_error() ){
110    ok=false;
111    *error << "error: std squared should be variance" << std::endl;
112    *error << "std2: " << a.std()*a.std() << std::endl;
113    *error << "variance: " << a.variance() << std::endl;
114    *error << "difference is: " << a.std()*a.std()-a.variance() << std::endl;
115  }
116 
117  if ( a.variance() != a.variance(a.mean()) ){
118    ok=false;
119    *error << "error: variance incorrect\n" << std::endl;
120    *error << "variance: " << a.variance() << std::endl;
121    *error << "mean: " << a.mean() << std::endl;
122    *error << "variance(mean) " << a.variance(a.mean()) << std::endl;
123  }
124
125  // Testing AveragerWeighted
126  *error << "testing AveragerWeighted" << std::endl;
127  theplu::yat::utility::vector x(3,0);
128  x(0)=0;
129  x(1)=1;
130  x(2)=2;
131  theplu::yat::utility::vector w(3,1);
132  theplu::yat::statistics::AveragerWeighted aw;
133  aw.add_values(x,w);
134  a.reset();
135  a.add_values(x);
136  const double tol=std::numeric_limits<double>().round_error();
137  if (!equal(a,aw,tol,error)){
138    *error << "error: AveragerWeighted with unitary weights should " 
139           << "be equal to Averager" << std::endl;
140    ok=false;
141  }
142
143  AveragerWeighted* aw2 = new AveragerWeighted(aw);
144  if (!equal(aw,*aw2,tol,error)){
145    *error << "error: AveragerWeighted copy constructor " << std::endl;
146    ok=false;
147  }
148   
149  aw2->add(12,0);
150  if (!equal(aw,*aw2,tol,error)){
151    *error << "error: AveragerWeighted adding a data point with weight=0 " 
152           << "should make no change " << std::endl;
153    ok=false;
154  }
155     
156  aw2->reset();
157  w*=17;
158  aw2->add_values(x,w);
159  if (!equal(aw,*aw2,tol,error)){
160    *error << "error: AveragerWeighted rescaling weights " 
161           << "should make no change " << std::endl;
162    ok=false;
163  }
164  delete aw2;
165 
166
167  *error << "testing AveragerPair" << std::endl;
168  AveragerPair ap;
169  for (int i=0; i<10; i++)
170    ap.add(static_cast<double>(i),i);
171  if (fabs(ap.correlation()-1)>tol){
172    ok=false;
173    *error << "correlation: " << ap.correlation() << std::endl;
174    *error << "error: correlation between identical vectors should be unity" 
175           << std::endl;
176  }
177  if (ap.x_averager().variance()!=ap.covariance()){
178    ok=false;
179    *error << "error: covariance of identical vectors should equal to variance" 
180           << std::endl;
181  }
182  AveragerPair* ap2 = new AveragerPair(ap);
183  delete ap2;
184
185  *error << "testing AveragerPairWeighted" << std::endl;
186  AveragerPairWeighted apw;
187  x(0)=0; x(1)=1; x(2)=2;
188  theplu::yat::utility::vector y(3,0);
189  x(0)=0; x(1)=0; x(2)=2;
190  apw.add_values(x,y,w,w);
191  ap.reset();
192  ap.add_values(x,y);
193  if (!equal(ap,apw,tol,error)){
194    *error << "error: AveragerPairWeighted with unitary weights should " 
195           << "be equal to AveragerPair" << std::endl;
196    ok=false;
197  }
198
199  AveragerPairWeighted* apw2 = new AveragerPairWeighted(apw);
200  if (!equal(apw,*apw2,tol,error)){
201    *error << "error: AveragerPairWeighted copy constructor " << std::endl;
202    ok=false;
203  }
204   
205  apw2->add(12,23222.03,32.3,0);
206  if (!equal(apw,*apw2,tol,error)){
207    *error << "error: AveragerWeighted adding a data point with weight=0 " 
208           << "should make no change " << std::endl;
209    ok=false;
210  }
211     
212  apw2->reset();
213  w*=17;
214  apw2->add_values(x,y,w,w);
215  if (!equal(apw,*apw2,tol,error)){
216    *error << "error: AveragerWeighted rescaling weights " 
217           << "should make no change " << std::endl;
218    ok=false;
219  }
220  delete apw2;
221
222  if (error!=&std::cerr)
223    delete error;
224
225  if (!ok)
226    return -1;
227  return 0;
228}
229
230bool equal(const Averager& a, const Averager& b)
231{
232//  std::cout << (a.n()==b.n()) << std::endl;
233//  std::cout << (a.mean()==b.mean()) << std::endl;
234//  std::cout << (fabs(a.variance()-b.variance()<1e-15)) << std::endl;
235  return (a.n()==b.n() && a.mean()==b.mean() &&
236          fabs(a.variance()-b.variance()<1e-15));
237}
238
239bool equal(const AveragerWeighted& a, const AveragerWeighted& b, 
240           const double tol, std::ostream* error)
241{
242  bool equal = true;
243  if  ( fabs(a.mean()-b.mean())>tol){
244    equal=false;
245    *error << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
246  }
247  if ( fabs(a.variance()-b.variance())>tol ) {
248    equal=false;
249    *error << "error for variance:\t" << a.variance() << " " << b.variance() 
250           << std::endl;
251  } 
252  if ( fabs(a.standard_error()-b.standard_error())>tol ) {
253    equal =false;
254    *error << "error for standard error:\t" << std::endl;
255  }
256  return equal;
257}
258
259bool equal(const Averager& a, const AveragerWeighted& b, const double tol, 
260           std::ostream* error)
261{
262  bool equal = true;
263  if  ( fabs(a.mean()-b.mean())>tol){
264    equal=false;
265    *error << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
266  }
267  if ( fabs(a.variance()-b.variance())>tol ) {
268    equal=false;
269    *error << "error for variance:\t" << a.variance() << " " << b.variance() 
270           << std::endl;
271  } 
272  if ( fabs(a.standard_error()-b.standard_error())>tol ) {
273    equal =false;
274    *error << "error for standard error:\t" << std::endl;
275  }
276  return equal;
277}
278
279bool equal(const AveragerPair& a, const AveragerPair& b, 
280           const double tol, std::ostream* error)
281{
282  bool ok = true;
283  if  ( fabs(a.covariance()-b.covariance())>tol){
284    ok=false;
285    *error << "error covariance: " << a.covariance() << "\t" 
286           << b.covariance() << std::endl;
287  }
288  if ( fabs(a.correlation()-b.correlation())>tol ) {
289    ok=false;
290    *error << "error correlation" << std::endl;
291  } 
292  return ok;
293}
294
295bool equal(const AveragerPair& a, const AveragerPairWeighted& b, 
296           const double tol, std::ostream* error)
297{
298  bool ok = true;
299  if  ( fabs(a.covariance()-b.covariance())>tol){
300    ok=false;
301    *error << "error covariance: " << a.covariance() << "\t" 
302           << b.covariance() << std::endl;
303  }
304  if ( fabs(a.correlation()-b.correlation())>tol ) {
305    ok=false;
306    *error << "error correlation" << std::endl;
307  } 
308  if ( !equal(a.x_averager(),b.x_averager(),tol,error)) {
309    ok =false;
310    *error << "error for x_averager():\t" << std::endl;
311  }
312  return ok;
313}
314bool equal(const AveragerPairWeighted& a, const AveragerPairWeighted& b, 
315           const double tol, std::ostream* error)
316{
317  bool ok = true;
318  if  ( fabs(a.covariance()-b.covariance())>tol){
319    ok=false;
320    *error << "error covariance: " << a.covariance() << "\t" 
321           << b.covariance() << std::endl;
322  }
323  if ( fabs(a.correlation()-b.correlation())>tol ) {
324    ok=false;
325    *error << "error correlation" << std::endl;
326  } 
327  if ( !equal(a.x_averager(),b.x_averager(),tol,error)) {
328    ok =false;
329    *error << "error for x_averager():\t" << std::endl;
330  }
331  return ok;
332}
333
334
335
Note: See TracBrowser for help on using the repository browser.