source: trunk/test/averager_test.cc @ 1297

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

Fixes #361. Removed tests negative count n in averager.

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.3 KB
Line 
1// $Id: averager_test.cc 1297 2008-05-13 14:20:23Z jari $
2
3/*
4  Copyright (C) 2005 Peter Johansson
5  Copyright (C) 2006 Jari Häkkinen, Peter Johansson, Markus Ringnér
6  Copyright (C) 2007, 2008 Jari Häkkinen, Peter Johansson
7
8  This file is part of the yat library, http://trac.thep.lu.se/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 "Suite.h"
27
28#include "yat/statistics/Averager.h"
29#include "yat/statistics/AveragerPair.h"
30#include "yat/statistics/AveragerPairWeighted.h"
31#include "yat/statistics/AveragerWeighted.h"
32#include "yat/utility/Vector.h"
33
34#include <cmath>
35#include <fstream>
36#include <limits>
37#include <iostream>
38#include <set>
39
40using namespace theplu::yat;
41using namespace theplu::yat::statistics;
42
43using theplu::yat::test::Suite;
44
45//Forward declarations
46bool equal(const Averager&, const Averager&, unsigned int, Suite&);
47bool equal(const AveragerWeighted&, const AveragerWeighted&, unsigned int, 
48           Suite& suite);
49bool equal(const Averager&, const AveragerWeighted&, unsigned int, Suite& suite);
50bool equal(const AveragerPair&, const AveragerPair&, unsigned int, Suite& suite);
51bool equal(const AveragerPair&, const AveragerPairWeighted&, unsigned int, 
52           Suite& suite);
53bool equal(const AveragerPairWeighted&, const AveragerPairWeighted&,
54           unsigned int, Suite& suite);
55
56
57int main(int argc, char* argv[])
58{ 
59
60  Suite suite(argc, argv);
61
62  // Testing Averager
63  suite.out() << "testing Averager" << std::endl;
64  Averager a;
65  a.add(1);
66  a.add(3);
67  a.add(5);
68  if (a.n()!=3 || a.mean()!=3 || a.sum_xx()!=35){
69    suite.add(false);
70    suite.err() << "error: add\n";
71  }
72  Averager b(a);
73  b.add(5,-1);
74  if (b.n()!=2 || b.mean()!=2 || b.sum_xx()!=10){
75    suite.add(false);
76    suite.err() << "error: add with negative n\n";
77  }
78
79  Averager* a1 = new Averager(1.0+3+5,1.0+9+25,3);
80  unsigned int tol = 10;
81  if (!equal(a,*a1, tol, suite)){
82    suite.add(false);
83    suite.err().precision(25);
84    suite.err() << a.sum_x() << '\t' << a1->sum_x() << std::endl;
85    suite.err() << a.sum_xx() << '\t' << a1->sum_xx() << std::endl;
86    suite.err() << a.n() << '\t' << a1->n() << std::endl;
87    suite.err() << a.variance() << '\t' << a1->variance() << std::endl;
88    suite.err() << a.mean() << '\t' << a1->mean() << std::endl;
89    suite.err() << a.mean() - a1->mean() << std::endl;
90    suite.err() << a.variance() - a1->variance() << std::endl;
91    suite.err() << "error: Averager(double x, double xx, long n)\n";
92  }
93  delete a1;
94
95  a1 = new Averager(a);
96  if (!equal(a,*a1, tol, suite)){
97    suite.add(false);
98    suite.err() << "error: Copy constructor\n";
99  }
100  delete a1;
101
102  a.add(3,5);
103  if (std::abs(a.standard_error()-sqrt(a.variance()/a.n()))>
104      std::numeric_limits<double>().round_error() ){
105    suite.add(false);
106    suite.err() << "error: standard_error\n";
107  }
108     
109  if ( std::abs(a.variance() - a.std()*a.std())>
110       std::numeric_limits<double>().round_error() ){
111    suite.add(false);
112    suite.err() << "error: std squared should be variance" << std::endl;
113    suite.err() << "std2: " << a.std()*a.std() << std::endl;
114    suite.err() << "variance: " << a.variance() << std::endl;
115    suite.err() << "difference is: " << a.std()*a.std()-a.variance() << std::endl;
116  }
117 
118  if ( a.variance() != a.variance(a.mean()) ){
119    suite.add(false);
120    suite.err() << "error: variance incorrect\n" << std::endl;
121    suite.err() << "variance: " << a.variance() << std::endl;
122    suite.err() << "mean: " << a.mean() << std::endl;
123    suite.err() << "variance(mean) " << a.variance(a.mean()) << std::endl;
124  }
125  theplu::yat::utility::Vector* tmp_vec = new theplu::yat::utility::Vector(10);
126  add(a, tmp_vec->begin(), tmp_vec->end());
127  delete tmp_vec;
128  std::set<double>* tmp_set = new std::set<double>;
129  tmp_set->insert(1.0);
130  tmp_set->insert(2.0);
131  add(a, tmp_set->begin(), tmp_set->end());
132  delete tmp_set;
133
134
135  // Testing AveragerWeighted
136  suite.err() << "testing AveragerWeighted" << std::endl;
137  theplu::yat::utility::Vector x(3,0);
138  x(0)=0;
139  x(1)=1;
140  x(2)=2;
141  theplu::yat::utility::Vector w(3,1);
142  theplu::yat::statistics::AveragerWeighted aw;
143  add(aw, x.begin(), x.end(), w.begin());
144  a.reset();
145  add(a, x.begin(), x.end());
146  if (!equal(a,aw,tol,suite)){
147    suite.err() << "error: AveragerWeighted with unitary weights should " 
148           << "be equal to Averager" << std::endl;
149    suite.add(false);
150  }
151
152  AveragerWeighted* aw2 = new AveragerWeighted(aw);
153  if (!equal(aw,*aw2,tol,suite)){
154    suite.err() << "error: AveragerWeighted copy constructor " << std::endl;
155    suite.add(false);
156  }
157   
158  aw2->add(12,0);
159  if (!equal(aw,*aw2,tol,suite)){
160    suite.err() << "error: AveragerWeighted adding a data point with weight=0 " 
161           << "should make no change " << std::endl;
162    suite.add(false);
163  }
164     
165  aw2->reset();
166  w*=17;
167  add(*aw2, x.begin(), x.end(), w.begin());
168  if (!equal(aw,*aw2,tol,suite)){
169    suite.err() << "error: AveragerWeighted rescaling weights " 
170           << "should make no change " << std::endl;
171    suite.add(false);
172  }
173  delete aw2;
174  {
175    theplu::yat::utility::Vector tmp(10);
176    add(aw, tmp.begin(), tmp.end());
177  }
178  {
179    std::set<double> tmp;
180    tmp.insert(1.0);
181    tmp.insert(2.0);
182    add(aw, tmp.begin(), tmp.end());
183  }
184 
185
186  suite.out() << "testing AveragerPair" << std::endl;
187  AveragerPair ap;
188  for (int i=0; i<10; i++)
189    ap.add(static_cast<double>(i),i);
190  if (std::abs(ap.correlation()-1)>tol){
191    suite.add(false);
192    suite.err() << "correlation: " << ap.correlation() << std::endl;
193    suite.err() << "error: correlation between identical vectors should be unity" 
194           << std::endl;
195  }
196  if (ap.x_averager().variance()!=ap.covariance()){
197    suite.add(false);
198    suite.err() << "error: covariance of identical vectors should equal to variance" 
199           << std::endl;
200  }
201  AveragerPair* ap2 = new AveragerPair(ap);
202  delete ap2;
203
204  for (int i=0; i<8; i++)
205    ap.add(static_cast<double>(i),i,-1);
206  if (std::abs(ap.correlation()-1)>tol) {
207    suite.add(false);
208    suite.err() << "correlation after removal of data: " << ap.correlation()
209                << std::endl;
210    suite.err() << "error: correlation between identical vectors is unity" 
211                << std::endl;
212  }
213
214
215  suite.out() << "testing AveragerPairWeighted" << std::endl;
216  AveragerPairWeighted apw;
217  x(0)=0; x(1)=1; x(2)=2;
218  theplu::yat::utility::Vector y(3,0);
219  y(0)=0; y(1)=0; y(2)=2;
220  add(apw, x.begin(), x.end(), y.begin(), w.begin(), w.begin());
221  ap.reset();
222  add(ap, x.begin(), x.end(), y.begin());
223  if (!equal(ap,apw,tol,suite)){
224    suite.err() << "error: AveragerPairWeighted with unitary weights should " 
225           << "be equal to AveragerPair" << std::endl;
226    suite.add(false);
227  }
228
229  AveragerPairWeighted* apw2 = new AveragerPairWeighted(apw);
230  if (!equal(apw,*apw2,tol,suite)){
231    suite.err() << "error: AveragerPairWeighted copy constructor " << std::endl;
232    suite.add(false);
233  }
234   
235  apw2->add(12,23222.03,32.3,0);
236  if (!equal(apw,*apw2,tol,suite)){
237    suite.err() << "error: AveragerWeighted adding a data point with weight=0 " 
238           << "should make no change " << std::endl;
239    suite.add(false);
240  }
241     
242  apw2->reset();
243  w*=17;
244  add(*apw2, x.begin(), x.end(), y.begin(), w.begin(), w.begin());
245  if (!equal(apw,*apw2,tol,suite)){
246    suite.err() << "error: AveragerWeighted rescaling weights " 
247           << "should make no change " << std::endl;
248    suite.add(false);
249  }
250  delete apw2;
251
252  return suite.return_value();
253}
254
255bool equal(const Averager& a, const Averager& b, unsigned int tol, 
256           Suite& suite)
257{
258//  std::cout << (a.n()==b.n()) << std::endl;
259//  std::cout << (a.mean()==b.mean()) << std::endl;
260//  std::cout << (std::abs(a.variance()-b.variance()<1e-15)) << std::endl;
261  return (a.n()==b.n() && suite.equal(a.mean(),b.mean(),tol) &&
262          suite.equal(a.variance(),b.variance(),tol));
263}
264
265bool equal(const AveragerWeighted& a, const AveragerWeighted& b, 
266           const unsigned int tol, Suite& suite)
267{
268  bool equal = true;
269  if  ( !suite.equal(a.mean(),b.mean(),tol)){
270    equal=false;
271    suite.err() << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
272    suite.err() << "difference:\t" << a.mean()-b.mean() << std::endl;
273  }
274  if ( !suite.equal(a.variance(),b.variance(),tol) ) {
275    equal=false;
276    suite.err() << "error for variance:\t" << a.variance() << " " 
277                << b.variance() << std::endl;
278  } 
279  if ( !suite.equal(a.standard_error(),b.standard_error(),tol) ) {
280    equal =false;
281    suite.err() << "error for standard error:\t" << std::endl;
282  }
283  return equal;
284}
285
286bool equal(const Averager& a, const AveragerWeighted& b, const unsigned int tol, 
287           Suite& suite)
288{
289  bool equal = true;
290  if  ( !suite.equal(a.mean(),b.mean(),tol)){
291    equal=false;
292    suite.err() << "mean:\t" << a.mean() << "\t" << b.mean() << std::endl;
293  }
294  if ( !suite.equal(a.variance(),b.variance(),tol) ) {
295    equal=false;
296    suite.err() << "error for variance:\t" << a.variance() << " " 
297                << b.variance() << std::endl;
298  } 
299  if ( !suite.equal(a.standard_error(), b.standard_error(),tol) ) {
300    equal =false;
301    suite.err() << "error for standard error:\t" << std::endl;
302  }
303  return equal;
304}
305
306bool equal(const AveragerPair& a, const AveragerPair& b, 
307           const unsigned int tol, Suite& suite)
308{
309  bool ok = true;
310  if  ( std::abs(a.covariance()-b.covariance())>tol){
311    suite.add(false);
312    suite.err() << "error covariance: " << a.covariance() << "\t" 
313           << b.covariance() << std::endl;
314  }
315  if ( std::abs(a.correlation()-b.correlation())>tol ) {
316    suite.add(false);
317    suite.err() << "error correlation" << std::endl;
318  } 
319  return ok;
320}
321
322bool equal(const AveragerPair& a, const AveragerPairWeighted& b, 
323           const unsigned int tol, Suite& suite)
324{
325  bool ok = true;
326  if  ( std::abs(a.covariance()-b.covariance())>tol){
327    suite.add(false);
328    suite.err() << "error covariance: " << a.covariance() << "\t" 
329           << b.covariance() << std::endl;
330  }
331  if ( std::abs(a.correlation()-b.correlation())>tol ) {
332    suite.add(false);
333    suite.err() << "error correlation" << std::endl;
334    suite.err() << "unweighted:" << a.correlation() << std::endl;
335    suite.err() << "weighted:" << b.correlation() << std::endl;
336    suite.err() << "difference:" << a.correlation()-b.correlation() 
337                << std::endl;
338  } 
339  if ( !equal(a.x_averager(),b.x_averager(),tol,suite)) {
340    ok =false;
341    suite.err() << "error for x_averager():\t" << std::endl;
342  }
343  return ok;
344}
345bool equal(const AveragerPairWeighted& a, const AveragerPairWeighted& b, 
346           const unsigned int tol, Suite& suite)
347{
348  bool ok = true;
349  if  ( !suite.equal(a.covariance(),b.covariance(),tol) ){
350    suite.add(false);
351    suite.err() << "error covariance: " << a.covariance() << "\t" 
352           << b.covariance() << std::endl;
353  }
354  if ( !suite.equal(a.correlation(), b.correlation(), tol) ) {
355    suite.add(false);
356    suite.err() << "error correlation" << std::endl;
357    suite.err() << "a:" << a.correlation() << std::endl;
358    suite.err() << "b:" << b.correlation() << std::endl;
359    suite.err() << "difference:" << a.correlation()-b.correlation() 
360                << std::endl;
361    suite.err() << "tol:" << tol << std::endl;
362  } 
363  if ( !equal(a.x_averager(),b.x_averager(),tol,suite)) {
364    ok =false;
365    suite.err() << "error for x_averager():\t" << std::endl;
366  }
367  return ok;
368}
369
370
371
Note: See TracBrowser for help on using the repository browser.