Commit | Line | Data |
---|---|---|
15637ed4 RG |
1 | // @(#)Sample.cc 6.2 (Berkeley) 2/25/91 |
2 | ||
3 | // Modified for Berkeley Unix by Donn Seeley, donn@okeeffe.berkeley.edu | |
4 | ||
5 | // This may look like C code, but it is really -*- C++ -*- | |
6 | /* | |
7 | Copyright (C) 1988 Free Software Foundation | |
8 | written by Dirk Grunwald (grunwald@cs.uiuc.edu) | |
9 | ||
10 | This file is part of GNU CC. | |
11 | ||
12 | GNU CC is distributed in the hope that it will be useful, | |
13 | but WITHOUT ANY WARRANTY. No author or distributor | |
14 | accepts responsibility to anyone for the consequences of using it | |
15 | or for whether it serves any particular purpose or works at all, | |
16 | unless he says so in writing. Refer to the GNU CC General Public | |
17 | License for full details. | |
18 | ||
19 | Everyone is granted permission to copy, modify and redistribute | |
20 | GNU CC, but only under the conditions described in the | |
21 | GNU CC General Public License. A copy of this license is | |
22 | supposed to have been given to you along with GNU CC so you | |
23 | can know your rights and responsibilities. It should be in a | |
24 | file named COPYING. Among other things, the copyright notice | |
25 | and this notice must be preserved on all copies. | |
26 | */ | |
27 | #ifdef __GNUG__ | |
28 | #pragma implementation | |
29 | #endif | |
30 | #include <stream.h> | |
31 | #include <SmplStat.h> | |
32 | #include <math.h> | |
33 | ||
34 | // error handling | |
35 | ||
36 | void default_SampleStatistic_error_handler(const char* msg) | |
37 | { | |
38 | cerr << "Fatal SampleStatistic error. " << msg << "\n"; | |
39 | exit(1); | |
40 | } | |
41 | ||
42 | one_arg_error_handler_t SampleStatistic_error_handler = default_SampleStatistic_error_handler; | |
43 | ||
44 | one_arg_error_handler_t set_SampleStatistic_error_handler(one_arg_error_handler_t f) | |
45 | { | |
46 | one_arg_error_handler_t old = SampleStatistic_error_handler; | |
47 | SampleStatistic_error_handler = f; | |
48 | return old; | |
49 | } | |
50 | ||
51 | void SampleStatistic::error(const char* msg) | |
52 | { | |
53 | (*SampleStatistic_error_handler)(msg); | |
54 | } | |
55 | ||
56 | // t-distribution: given p-value and degrees of freedom, return t-value | |
57 | // adapted from Peizer & Pratt JASA, vol63, p1416 | |
58 | ||
59 | double tval(double p, int df) | |
60 | { | |
61 | double t; | |
62 | int positive = p >= 0.5; | |
63 | p = (positive)? 1.0 - p : p; | |
64 | if (p <= 0.0 || df <= 0) | |
65 | t = HUGE; | |
66 | else if (p == 0.5) | |
67 | t = 0.0; | |
68 | else if (df == 1) | |
69 | t = 1.0 / tan((p + p) * 1.57079633); | |
70 | else if (df == 2) | |
71 | t = sqrt(1.0 / ((p + p) * (1.0 - p)) - 2.0); | |
72 | else | |
73 | { | |
74 | double ddf = df; | |
75 | double a = sqrt(log(1.0 / (p * p))); | |
76 | double aa = a * a; | |
77 | a = a - ((2.515517 + (0.802853 * a) + (0.010328 * aa)) / | |
78 | (1.0 + (1.432788 * a) + (0.189269 * aa) + | |
79 | (0.001308 * aa * a))); | |
80 | t = ddf - 0.666666667 + 1.0 / (10.0 * ddf); | |
81 | t = sqrt(ddf * (exp(a * a * (ddf - 0.833333333) / (t * t)) - 1.0)); | |
82 | } | |
83 | return (positive)? t : -t; | |
84 | } | |
85 | ||
86 | void | |
87 | SampleStatistic::reset() | |
88 | { | |
89 | n = 0; x = x2 = 0.0; | |
90 | maxValue = -HUGE; | |
91 | minValue = HUGE; | |
92 | } | |
93 | ||
94 | void | |
95 | SampleStatistic::operator+=(double value) | |
96 | { | |
97 | n += 1; | |
98 | x += value; | |
99 | x2 += (value * value); | |
100 | if ( minValue > value) minValue = value; | |
101 | if ( maxValue < value) maxValue = value; | |
102 | } | |
103 | ||
104 | double | |
105 | SampleStatistic::mean() | |
106 | { | |
107 | if ( n > 0) { | |
108 | return (x / n); | |
109 | } | |
110 | else { | |
111 | return ( 0.0 ); | |
112 | } | |
113 | } | |
114 | ||
115 | double | |
116 | SampleStatistic::var() | |
117 | { | |
118 | if ( n > 1) { | |
119 | return(( x2 - ((x * x) / n)) / ( n - 1)); | |
120 | } | |
121 | else { | |
122 | return ( 0.0 ); | |
123 | } | |
124 | } | |
125 | ||
126 | double | |
127 | SampleStatistic::stdDev() | |
128 | { | |
129 | if ( n <= 0 || this -> var() <= 0) { | |
130 | return(0); | |
131 | } else { | |
132 | return( (double) sqrt( var() ) ); | |
133 | } | |
134 | } | |
135 | ||
136 | double | |
137 | SampleStatistic::confidence(int interval) | |
138 | { | |
139 | int df = n - 1; | |
140 | if (df <= 0) return HUGE; | |
141 | double t = tval(double(100 + interval) * 0.005, df); | |
142 | if (t == HUGE) | |
143 | return t; | |
144 | else | |
145 | return (t * stdDev()) / sqrt(double(n)); | |
146 | } | |
147 | ||
148 | double | |
149 | SampleStatistic::confidence(double p_value) | |
150 | { | |
151 | int df = n - 1; | |
152 | if (df <= 0) return HUGE; | |
153 | double t = tval((1.0 + p_value) * 0.5, df); | |
154 | if (t == HUGE) | |
155 | return t; | |
156 | else | |
157 | return (t * stdDev()) / sqrt(double(n)); | |
158 | } | |
159 | ||
160 | ||
161 | #include <SmplHist.h> | |
162 | ||
163 | const int SampleHistogramMinimum = -2; | |
164 | const int SampleHistogramMaximum = -1; | |
165 | ||
166 | SampleHistogram::SampleHistogram(double low, double high, double width) | |
167 | { | |
168 | if (high < low) { | |
169 | double t = high; | |
170 | high = low; | |
171 | low = t; | |
172 | } | |
173 | ||
174 | if (width == -1) { | |
175 | width = (high - low) / 10; | |
176 | } | |
177 | ||
178 | howManyBuckets = int((high - low) / width) + 2; | |
179 | bucketCount = new int[howManyBuckets]; | |
180 | bucketLimit = new double[howManyBuckets]; | |
181 | double lim = low; | |
182 | for (int i = 0; i < howManyBuckets; i++) { | |
183 | bucketCount[i] = 0; | |
184 | bucketLimit[i] = lim; | |
185 | lim += width; | |
186 | } | |
187 | bucketLimit[howManyBuckets-1] = HUGE; /* from math.h */ | |
188 | } | |
189 | ||
190 | SampleHistogram::~SampleHistogram() | |
191 | { | |
192 | if (howManyBuckets > 0) { | |
193 | delete bucketCount; | |
194 | delete bucketLimit; | |
195 | } | |
196 | } | |
197 | ||
198 | void | |
199 | SampleHistogram::operator+=(double value) | |
200 | { | |
201 | int i; | |
202 | for (i = 0; i < howManyBuckets; i++) { | |
203 | if (value < bucketLimit[i]) break; | |
204 | } | |
205 | bucketCount[i]++; | |
206 | this->SampleStatistic::operator+=(value); | |
207 | } | |
208 | ||
209 | int | |
210 | SampleHistogram::similarSamples(double d) | |
211 | { | |
212 | int i; | |
213 | for (i = 0; i < howManyBuckets; i++) { | |
214 | if (d < bucketLimit[i]) return(bucketCount[i]); | |
215 | } | |
216 | return(0); | |
217 | } | |
218 | ||
219 | void | |
220 | SampleHistogram::printBuckets(ostream& s) | |
221 | { | |
222 | for(int i = 0; i < howManyBuckets; i++) { | |
223 | if (bucketLimit[i] >= HUGE) { | |
224 | s << "< max : " << bucketCount[i] << "\n"; | |
225 | } else { | |
226 | s << "< " << bucketLimit[i] << " : " << bucketCount[i] << "\n"; | |
227 | } | |
228 | } | |
229 | } | |
230 | ||
231 | void | |
232 | SampleHistogram::reset() | |
233 | { | |
234 | this->SampleStatistic::reset(); | |
235 | if (howManyBuckets > 0) { | |
236 | for (register int i = 0; i < howManyBuckets; i++) { | |
237 | bucketCount[i] = 0; | |
238 | } | |
239 | } | |
240 | } | |
241 |