1 // This may look like C code, but it is really -*- C++ -*-
3 Copyright (C) 1988 Free Software Foundation
4 written by Dirk Grunwald (grunwald@cs.uiuc.edu)
6 This file is part of the GNU C++ Library. This library is free
7 software; you can redistribute it and/or modify it under the terms of
8 the GNU Library General Public License as published by the Free
9 Software Foundation; either version 2 of the License, or (at your
10 option) any later version. This library is distributed in the hope
11 that it will be useful, but WITHOUT ANY WARRANTY; without even the
12 implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
13 PURPOSE. See the GNU Library General Public License for more details.
14 You should have received a copy of the GNU Library General Public
15 License along with this library; if not, write to the Free Software
16 Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
19 #pragma implementation
29 #define HUGE_VAL DBL_MAX
36 #include <simgear/debug/logstream.hxx>
37 #include "SGSmplstat.hxx"
40 void SampleStatistic::error (const char *msg)
42 SG_LOG(SG_GENERAL, SG_ALERT, msg);
45 // t-distribution: given p-value and degrees of freedom, return t-value
46 // adapted from Peizer & Pratt JASA, vol63, p1416
48 double tval (double p, int df)
51 int positive = p >= 0.5;
52 p = (positive) ? 1.0 - p : p;
53 if (p <= 0.0 || df <= 0)
58 t = 1.0 / tan ((p + p) * 1.57079633);
60 t = sqrt (1.0 / ((p + p) * (1.0 - p)) - 2.0);
64 double a = sqrt (log (1.0 / (p * p)));
66 a = a - ((2.515517 + (0.802853 * a) + (0.010328 * aa)) /
67 (1.0 + (1.432788 * a) + (0.189269 * aa) +
68 (0.001308 * aa * a)));
69 t = ddf - 0.666666667 + 1.0 / (10.0 * ddf);
70 t = sqrt (ddf * (exp (a * a * (ddf - 0.833333333) / (t * t)) - 1.0));
72 return (positive) ? t : -t;
75 void SampleStatistic::reset ()
83 void SampleStatistic::operator += (double value)
87 x2 += (value * value);
94 double SampleStatistic::mean () const
106 double SampleStatistic::var () const
110 return ((x2 - ((x * x) / n)) / (n - 1));
118 double SampleStatistic::stdDev () const
120 if (n <= 0 || this->var () <= 0)
126 return ((double) sqrt (var ()));
130 double SampleStatistic::confidence (int interval) const
135 double t = tval (double (100 + interval) * 0.005, df);
139 return (t * stdDev ()) / sqrt (double (n));
142 double SampleStatistic::confidence (double p_value) const
147 double t = tval ((1.0 + p_value) * 0.5, df);
151 return (t * stdDev ()) / sqrt (double (n));