source: trunk/test/averager_test.cc @ 831

Last change on this file since 831 was 831, checked in by Peter, 16 years ago

Refs #185.

  • 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 831 2007-03-27 13:22:09Z peter $
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 <fstream>
33#include <limits>
34#include <iostream>
35
36using namespace theplu::yat::statistics;
37
38//Forward declarations
39bool equal(const Averager&, const Averager&);
40bool equal(const AveragerWeighted&, const AveragerWeighted&, const double, 
41           std::ostream* error);
42bool equal(const Averager&, const AveragerWeighted&, const double, 
43           std::ostream* error);
44bool equal(const AveragerPair&, const AveragerPair&, 
45           const double, std::ostream* error);
46bool equal(const AveragerPair&, const AveragerPairWeighted&, 
47           const double, std::ostream* error);
48bool equal(const AveragerPairWeighted&, const AveragerPairWeighted&, 
49           const double, std::ostream* error);
50
51
52int main(const int argc,const char* argv[])
53{ 
54
55  std::ostream* error;
56  if (argc>1 && argv[1]==std::string("-v"))
57    error = &std::cerr;
58  else {
59    error = new std::ofstream("/dev/null");
60    if (argc>1)
61      std::cout << "averager_test -v : for printing extra information\n";
62  }
63  bool ok = true;
64
65  // Testing Averager
66  *error << "testing Averager" << std::endl;
67  Averager a;
68  a.add(1);
69  a.add(3);
70  a.add(5);
71  if (a.n()!=3 || a.mean()!=3 || a.sum_xx()!=35){
72    ok=false;
73    *error << "error: add\n";
74  }
75
76
77  Averager* a1 = new Averager(1.0+3+5,1.0+9+25,3);
78  if (!equal(a,*a1)){
79    ok=false;
80    std::cout.precision(25);
81    std::cout << equal(a,*a1) << std::endl;
82    std::cout << a.sum_x() << '\t' << a1->sum_x() << std::endl;
83    std::cout << a.sum_xx() << '\t' << a1->sum_xx() << std::endl;
84    std::cout << a.n() << '\t' << a1->n() << std::endl;
85    std::cout << a.variance() << '\t' << a1->variance() << std::endl;
86    std::cout << a.mean() << '\t' << a1->mean() << std::endl;
87    std::cout << a.mean() - a1->mean() << std::endl;
88    std::cout << a.variance() - a1->variance() << std::endl;
89    *error << "error: Averager(double x, double xx, u_long n)\n";
90  }
91  delete a1;
92
93  a1 = new Averager(a);
94  if (!equal(a,*a1)){
95    ok=false;
96    *error << "error: Copy constructor\n";
97  }
98  delete a1;
99
100  a.add(3,5);
101  if (! a.standard_error()==sqrt(a.variance()/a.n())){
102    ok=false;
103    *error << "error: standard_error\n";
104  }
105     
106
107  if ( fabs(a.variance() - a.std()*a.std())>
108       std::numeric_limits<double>().round_error() ){
109    ok=false;
110    *error << "error: std squared should be variance" << std::endl;
111    *error << "std2: " << a.std()*a.std() << std::endl;
112    *error << "variance: " << a.variance() << std::endl;
113    *error << "difference is: " << a.std()*a.std()-a.variance() << std::endl;
114  }
115 
116  if ( a.variance() != a.variance(a.mean()) ){
117    ok=false;
118    *error << "error: variance incorrect\n" << std::endl;
119    *error << "variance: " << a.variance() << std::endl;
120    *error << "mean: " << a.mean() << std::endl;
121    *error << "variance(mean) " << a.variance(a.mean()) << std::endl;
122  }
123
124  // Testing AveragerWeighted
125  *error << "testing AveragerWeighted" << std::endl;
126  theplu::yat::utility::vector x(3,0);
127  x(0)=0;
128  x(1)=1;
129  x(2)=2;
130  theplu::yat::utility::vector w(3,1);
131  theplu::yat::statistics::AveragerWeighted aw;
132  aw.add_values(x,w);
133  a.reset();
134  a.add_values(x);
135  const double tol=std::numeric_limits<double>().round_error();
136  if (!equal(a,aw,tol,error)){
137    *error << "error: AveragerWeighted with unitary weights should " 
138           << "be equal to Averager" << std::endl;
139    ok=false;
140  }
141
142  AveragerWeighted* aw2 = new AveragerWeighted(aw);
143  if (!equal(aw,*aw2,tol,error)){
144    *error << "error: AveragerWeighted copy constructor " << std::endl;
145    ok=false;
146  }
147   
148  aw2->add(12,0);
149  if (!equal(aw,*aw2,tol,error)){
150    *error << "error: AveragerWeighted adding a data point with weight=0 " 
151           << "should make no change " << std::endl;
152    ok=false;
153  }
154     
155  aw2->reset();
156  w*=17;
157  aw2->add_values(x,w);
158  if (!equal(aw,*aw2,tol,error)){
159    *error << "error: AveragerWeighted rescaling weights " 
160           << "should make no change " << std::endl;
161    ok=false;
162  }
163  delete aw2;
164 
165
166  *error << "testing AveragerPair" << std::endl;
167  AveragerPair ap;
168  for (int i=0; i<10; i++)
169    ap.add(static_cast<double>(i),i);
170  if (fabs(ap.correlation()-1)>tol){
171    ok=false;
172    *error << "correlation: " << ap.correlation() << std::endl;
173    *error << "error: correlation between identical vectors should be unity" 
174           << std::endl;
175  }
176  if (ap.x_averager().variance()!=ap.covariance()){
177    ok=false;
178    *error << "error: covariance of identical vectors should equal to variance" 
179           << std::endl;
180  }
181  AveragerPair* ap2 = new AveragerPair(ap);
182  delete ap2;
183
184  *error << "testing AveragerPairWeighted" << std::endl;
185  AveragerPairWeighted apw;
186  x(0)=0; x(1)=1; x(2)=2;
187  theplu::yat::utility::vector y(3,0);
188  x(0)=0; x(1)=0; x(2)=2;
189  apw.add_values(x,y,w,w);
190  ap.reset();
191  ap.add_values(x,y);
192  if (!equal(ap,apw,tol,error)){
193    *error << "error: AveragerPairWeighted with unitary weights should " 
194           << "be equal to AveragerPair" << std::endl;
195    ok=false;
196  }
197
198  AveragerPairWeighted* apw2 = new AveragerPairWeighted(apw);
199  if (!equal(apw,*apw2,tol,error)){
200    *error << "error: AveragerPairWeighted copy constructor " << std::endl;
201    ok=false;
202  }
203   
204  apw2->add(12,23222.03,32.3,0);
205  if (!equal(apw,*apw2,tol,error)){
206    *error << "error: AveragerWeighted adding a data point with weight=0 " 
207           << "should make no change " << std::endl;
208    ok=false;
209  }
210     
211  apw2->reset();
212  w*=17;
213  apw2->add_values(x,y,w,w);
214  if (!equal(apw,*apw2,tol,error)){
215    *error << "error: AveragerWeighted rescaling weights " 
216           << "should make no change " << std::endl;
217    ok=false;
218  }
219  delete apw2;
220
221  if (error!=&std::cerr)
222    delete error;
223
224  if (!ok)
225    return -1;
226  return 0;
227}
228
229bool equal(const Averager& a, const Averager& b)
230{
231//  std::cout << (a.n()==b.n()) << std::endl;
232//  std::cout << (a.mean()==b.mean()) << std::endl;
233//  std::cout << (a.variance()==b.variance()) << std::endl;
234  return (a.n()==b.n() && a.mean()==b.mean() && a.variance()==b.variance());
235}
236
237bool equal(const AveragerWeighted& a, const AveragerWeighted& b, 
238           const double tol, std::ostream* error)
239{
240  bool equal = true;
241  if  ( fabs(a.mean()-b.mean())>tol){
242    equal=false;
243    *error << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
244  }
245  if ( fabs(a.variance()-b.variance())>tol ) {
246    equal=false;
247    *error << "error for variance:\t" << a.variance() << " " << b.variance() 
248           << std::endl;
249  } 
250  if ( fabs(a.standard_error()-b.standard_error())>tol ) {
251    equal =false;
252    *error << "error for standard error:\t" << std::endl;
253  }
254  return equal;
255}
256
257bool equal(const Averager& a, const AveragerWeighted& b, const double tol, 
258           std::ostream* error)
259{
260  bool equal = true;
261  if  ( fabs(a.mean()-b.mean())>tol){
262    equal=false;
263    *error << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
264  }
265  if ( fabs(a.variance()-b.variance())>tol ) {
266    equal=false;
267    *error << "error for variance:\t" << a.variance() << " " << b.variance() 
268           << std::endl;
269  } 
270  if ( fabs(a.standard_error()-b.standard_error())>tol ) {
271    equal =false;
272    *error << "error for standard error:\t" << std::endl;
273  }
274  return equal;
275}
276
277bool equal(const AveragerPair& a, const AveragerPair& b, 
278           const double tol, std::ostream* error)
279{
280  bool ok = true;
281  if  ( fabs(a.covariance()-b.covariance())>tol){
282    ok=false;
283    *error << "error covariance: " << a.covariance() << "\t" 
284           << b.covariance() << std::endl;
285  }
286  if ( fabs(a.correlation()-b.correlation())>tol ) {
287    ok=false;
288    *error << "error correlation" << std::endl;
289  } 
290  return ok;
291}
292
293bool equal(const AveragerPair& a, const AveragerPairWeighted& b, 
294           const double tol, std::ostream* error)
295{
296  bool ok = true;
297  if  ( fabs(a.covariance()-b.covariance())>tol){
298    ok=false;
299    *error << "error covariance: " << a.covariance() << "\t" 
300           << b.covariance() << std::endl;
301  }
302  if ( fabs(a.correlation()-b.correlation())>tol ) {
303    ok=false;
304    *error << "error correlation" << std::endl;
305  } 
306  if ( !equal(a.x_averager(),b.x_averager(),tol,error)) {
307    ok =false;
308    *error << "error for x_averager():\t" << std::endl;
309  }
310  return ok;
311}
312bool equal(const AveragerPairWeighted& a, const AveragerPairWeighted& b, 
313           const double tol, std::ostream* error)
314{
315  bool ok = true;
316  if  ( fabs(a.covariance()-b.covariance())>tol){
317    ok=false;
318    *error << "error covariance: " << a.covariance() << "\t" 
319           << b.covariance() << std::endl;
320  }
321  if ( fabs(a.correlation()-b.correlation())>tol ) {
322    ok=false;
323    *error << "error correlation" << std::endl;
324  } 
325  if ( !equal(a.x_averager(),b.x_averager(),tol,error)) {
326    ok =false;
327    *error << "error for x_averager():\t" << std::endl;
328  }
329  return ok;
330}
331
332
333
Note: See TracBrowser for help on using the repository browser.