]>
Commit | Line | Data |
---|---|---|
b0f5e3fc | 1 | ////////////////////////////////////////////////////////////////////////// |
2 | // Alice ITS class to help keep statistical information // | |
3 | // // | |
4 | // version: 0.0.0 Draft. // | |
5 | // Date: April 18 1999 // | |
6 | // By: Bjorn S. Nilsen // | |
7 | // // | |
8 | ////////////////////////////////////////////////////////////////////////// | |
9 | #include <stdio.h> | |
10 | #include <math.h> | |
11 | #include "TMath.h" | |
12 | #include "AliITSstatistics2.h" | |
13 | ||
14 | ClassImp(AliITSstatistics2) | |
15 | ||
16 | // | |
17 | AliITSstatistics2::AliITSstatistics2() : TObject(){ | |
18 | // | |
19 | // default constructor | |
20 | // | |
07025d6d | 21 | fX = 0; |
22 | fY = 0; | |
23 | fYx = 0; | |
24 | fW = 0; | |
b0f5e3fc | 25 | fN = 0; |
26 | fOrder = 0; | |
27 | return; | |
28 | } | |
29 | ||
30 | ||
31 | AliITSstatistics2::AliITSstatistics2(Int_t order) : TObject(){ | |
32 | // | |
33 | // constructor to maximum moment/order order | |
34 | // | |
b0f5e3fc | 35 | fOrder = order; |
07025d6d | 36 | fX = new Double_t[order]; |
37 | fY = new Double_t[order]; | |
38 | fYx = new Double_t[order]; | |
39 | fW = new Double_t[order]; | |
40 | for(Int_t i=0;i<order;i++) {fX[i] = 0.0;fY[i] = 0.0; | |
41 | fYx[i] = 0.0; fW[i] = 0.0;} | |
b0f5e3fc | 42 | fN = 0; |
43 | return; | |
44 | } | |
45 | ||
46 | AliITSstatistics2::~AliITSstatistics2(){ | |
47 | // | |
48 | // destructor | |
49 | // | |
07025d6d | 50 | if(fX!=0) delete[] fX; |
51 | if(fY!=0) delete[] fY; | |
52 | if(fYx!=0) delete[] fYx; | |
53 | if(fW!=0) delete[] fW; | |
54 | fX = 0; | |
55 | fY = 0; | |
56 | fYx = 0; | |
57 | fW = 0; | |
b0f5e3fc | 58 | fN = 0; |
59 | fOrder = 0; | |
60 | } | |
61 | ||
62 | //_______________________________________________________________ | |
63 | AliITSstatistics2& AliITSstatistics2::operator=(AliITSstatistics2 &source){ | |
64 | // operator = | |
65 | ||
e8189707 | 66 | if(this==&source) return *this; |
67 | if(source.fOrder!=0){ | |
68 | this->fOrder = source.fOrder; | |
69 | this->fN = source.fN; | |
07025d6d | 70 | this->fX = new Double_t[this->fOrder]; |
71 | this->fW = new Double_t[this->fOrder]; | |
e8189707 | 72 | for(Int_t i=0;i<source.fOrder;i++){ |
07025d6d | 73 | this->fX[i] = source.fX[i]; |
74 | this->fW[i] = source.fW[i]; | |
e8189707 | 75 | } // end for i |
76 | }else{ | |
07025d6d | 77 | this->fX = 0; |
78 | this->fW = 0; | |
e8189707 | 79 | this->fN = 0; |
80 | this->fOrder = 0; | |
81 | }// end if source.fOrder!=0 | |
82 | return *this; | |
b0f5e3fc | 83 | } |
84 | //_______________________________________________________________ | |
ac74f489 | 85 | AliITSstatistics2::AliITSstatistics2(AliITSstatistics2 &source): |
86 | TObject(source){ | |
b0f5e3fc | 87 | // Copy constructor |
88 | ||
e8189707 | 89 | if(this==&source) return; |
90 | if(source.fOrder!=0){ | |
91 | this->fOrder = source.fOrder; | |
92 | this->fN = source.fN; | |
07025d6d | 93 | this->fX = new Double_t[this->fOrder]; |
94 | this->fW = new Double_t[this->fOrder]; | |
e8189707 | 95 | for(Int_t i=0;i<source.fOrder;i++){ |
07025d6d | 96 | this->fX[i] = source.fX[i]; |
97 | this->fW[i] = source.fW[i]; | |
e8189707 | 98 | } // end for i |
99 | }else{ | |
07025d6d | 100 | this->fX = 0; |
101 | this->fW = 0; | |
e8189707 | 102 | this->fN = 0; |
103 | this->fOrder = 0; | |
104 | }// end if source.fOrder!=0 | |
b0f5e3fc | 105 | } |
106 | ||
107 | void AliITSstatistics2::Reset(){ | |
108 | // | |
109 | // Reset/zero statistics | |
110 | // | |
07025d6d | 111 | for(Int_t i=0;i<fOrder;i++) {fX[i] = 0.0;fY[i] = 0.0; |
112 | fYx[i] = 0.0; fW[i] = 0.0;} | |
b0f5e3fc | 113 | fN = 0; |
114 | return; | |
115 | } | |
116 | ||
117 | void AliITSstatistics2::AddValue(Double_t y,Double_t x,Double_t w=1.0){ | |
118 | // | |
119 | // add next x,y pair to statistics | |
120 | // | |
121 | Double_t xs=1.0,ys=1.0,yxs=1.0,ws=1.0; | |
122 | Int_t i; | |
123 | ||
e8189707 | 124 | const Double_t kBig=1.0e+38; |
125 | ||
126 | if(y>kBig || x>kBig || w>kBig) return; | |
127 | ||
128 | ||
b0f5e3fc | 129 | fN++; |
130 | for(i=0;i<fOrder;i++){ | |
131 | xs *= x; | |
132 | ys *= y; | |
133 | yxs *= x*y; | |
134 | ws *= w; | |
07025d6d | 135 | fX[i] += xs*w; |
136 | fY[i] += ys*w; | |
137 | fYx[i] += yxs*w; | |
138 | fW[i] += ws; | |
b0f5e3fc | 139 | } // end for i |
140 | } | |
141 | ||
142 | Double_t AliITSstatistics2::GetXNth(Int_t order){ | |
143 | // | |
144 | // This give the unbiased estimator for the RMS. | |
145 | // | |
146 | ||
147 | Double_t s; | |
148 | ||
07025d6d | 149 | if(fW[0]!=0.0&&order<=fOrder) s = fX[order-1]/fW[0]; |
b0f5e3fc | 150 | else { |
151 | s = 0.0; | |
07025d6d | 152 | printf("AliITSstatistics2: error in GetNth: fOrder=%d fN=%d fW[0]=%f\n", |
153 | fOrder,fN,fW[0]); | |
b0f5e3fc | 154 | } // end else |
155 | return s; | |
156 | } | |
157 | Double_t AliITSstatistics2::GetYNth(Int_t order){ | |
158 | // | |
159 | // This give the unbiased estimator for the RMS. | |
160 | // | |
161 | Double_t s; | |
162 | ||
07025d6d | 163 | if(fW[0]!=0.0&&order<=fOrder) s = fY[order-1]/fW[0]; |
b0f5e3fc | 164 | else { |
165 | s = 0.0; | |
07025d6d | 166 | printf("AliITSstatistics2: error in GetNth: fOrder=%d fN=%d fW[0]=%f\n", |
167 | fOrder,fN,fW[0]); | |
b0f5e3fc | 168 | } // end else |
169 | return s; | |
170 | } | |
171 | Double_t AliITSstatistics2::GetYXNth(Int_t order){ | |
172 | // This give the unbiased estimator for the RMS. | |
173 | Double_t s; | |
174 | ||
07025d6d | 175 | if(fW[0]!=0.0&&order<=fOrder) s = fYx[order-1]/fW[0]; |
b0f5e3fc | 176 | else { |
177 | s = 0.0; | |
07025d6d | 178 | printf("AliITSstatistics2: error in GetNth: fOrder=%d fN=%d fW[0]=%f\n", |
179 | fOrder,fN,fW[0]); | |
b0f5e3fc | 180 | } // end else |
181 | return s; | |
182 | } | |
183 | Double_t AliITSstatistics2::GetRMSX(){ | |
184 | // This give the unbiased estimator for the RMS. | |
185 | Double_t x,x2,w,ww,s; | |
186 | ||
187 | x = GetMeanX(); // first order | |
188 | x2 = GetXNth(2); // second order | |
07025d6d | 189 | w = fW[0]; // first order - 1. |
190 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 191 | |
192 | if(w*w==ww) return (-1.0); | |
193 | s = (x2-x*x)*w*w/(w*w-ww); | |
194 | return TMath::Sqrt(s); | |
195 | } | |
196 | Double_t AliITSstatistics2::GetRMSY(){ | |
197 | // This give the unbiased estimator for the RMS. | |
198 | Double_t x,x2,w,ww,s; | |
199 | ||
200 | x = GetMeanY(); // first order | |
201 | x2 = GetYNth(2); // second order | |
07025d6d | 202 | w = fW[0]; // first order - 1. |
203 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 204 | |
205 | if(w*w==ww) return (-1.0); | |
206 | s = (x2-x*x)*w*w/(w*w-ww); | |
207 | return TMath::Sqrt(s); | |
208 | } | |
209 | Double_t AliITSstatistics2::GetRMSYX(){ | |
210 | // This give the unbiased estimator for the RMS. | |
211 | Double_t x,x2,w,ww,s; | |
212 | ||
213 | x = GetMeanYX(); // first order | |
214 | x2 = GetYXNth(2); // second order | |
07025d6d | 215 | w = fW[0]; // first order - 1. |
216 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 217 | |
218 | if(w*w==ww) return (-1.0); | |
219 | s = (x2-x*x)*w*w/(w*w-ww); | |
220 | return TMath::Sqrt(s); | |
221 | } | |
222 | Double_t AliITSstatistics2::GetErrorMeanY(){ | |
223 | //This is the error in the mean or the square root of the variance of the mean. | |
224 | Double_t rms,w,ww,s; | |
225 | ||
226 | rms = GetRMSY(); | |
07025d6d | 227 | w = fW[0]; |
228 | ww = fW[1]; | |
b0f5e3fc | 229 | s = rms*rms*ww/(w*w); |
230 | return TMath::Sqrt(s); | |
231 | } | |
232 | Double_t AliITSstatistics2::GetErrorMeanX(){ | |
233 | //This is the error in the mean or the square root of the variance of the mean. | |
234 | Double_t rms,w,ww,s; | |
235 | ||
236 | rms = GetRMSX(); | |
07025d6d | 237 | w = fW[0]; |
238 | ww = fW[1]; | |
b0f5e3fc | 239 | s = rms*rms*ww/(w*w); |
240 | return TMath::Sqrt(s); | |
241 | } | |
242 | Double_t AliITSstatistics2::GetErrorMeanYX(){ | |
243 | //This is the error in the mean or the square root of the variance of the mean. | |
244 | Double_t rms,w,ww,s; | |
245 | ||
246 | rms = GetRMSYX(); | |
07025d6d | 247 | w = fW[0]; |
248 | ww = fW[1]; | |
b0f5e3fc | 249 | s = rms*rms*ww/(w*w); |
250 | return TMath::Sqrt(s); | |
251 | } | |
252 | ||
253 | ||
254 | Double_t AliITSstatistics2::GetErrorRMSY(){ | |
255 | //This is the error in the mean or the square root of the variance of the mean. | |
256 | // at this moment this routine is only defined for weights=1. | |
257 | Double_t x,x2,x3,x4,w,ww,m2,m4,n,s; | |
258 | ||
07025d6d | 259 | if(fW[0]!=(Double_t)fN||GetN()<4) return (-1.); |
b0f5e3fc | 260 | x = GetMeanY(); // first order |
261 | x2 = GetYNth(2); // second order | |
07025d6d | 262 | w = fW[0]; // first order - 1. |
263 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 264 | if(w*w==ww) return (-1.0); |
265 | s = (x2-x*x)*w*w/(w*w-ww); | |
266 | ||
267 | m2 = s; | |
268 | n = (Double_t) GetN(); | |
269 | x3 = GetYNth(3); | |
270 | x4 = GetYNth(4); | |
271 | // This equation assumes that all of the weights are equal to 1. | |
272 | m4 = (n/(n-1.))*(x4-3.*x*x3+6.*x*x*x2-2.*x*x*x*x); | |
273 | s = (m4-(n-3.)*m2*m2/(n-1.))/n; | |
274 | return TMath::Sqrt(s); | |
275 | } | |
276 | Double_t AliITSstatistics2::GetErrorRMSX(){ | |
277 | //This is the error in the mean or the square root of the variance of the mean. | |
278 | // at this moment this routine is only defined for weights=1. | |
279 | Double_t x,x2,x3,x4,w,ww,m2,m4,n,s; | |
280 | ||
07025d6d | 281 | if(fW[0]!=(Double_t)fN||GetN()<4) return (-1.); |
b0f5e3fc | 282 | x = GetMeanX(); // first order |
283 | x2 = GetXNth(2); // second order | |
07025d6d | 284 | w = fW[0]; // first order - 1. |
285 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 286 | if(w*w==ww) return (-1.0); |
287 | s = (x2-x*x)*w*w/(w*w-ww); | |
288 | ||
289 | m2 = s; | |
290 | n = (Double_t) GetN(); | |
291 | x3 = GetXNth(3); | |
292 | x4 = GetXNth(4); | |
293 | // This equation assumes that all of the weights are equal to 1. | |
294 | m4 = (n/(n-1.))*(x4-3.*x*x3+6.*x*x*x2-2.*x*x*x*x); | |
295 | s = (m4-(n-3.)*m2*m2/(n-1.))/n; | |
296 | return TMath::Sqrt(s); | |
297 | } | |
298 | Double_t AliITSstatistics2::GetErrorRMSYX(){ | |
299 | //This is the error in the mean or the square root of the variance of the mean. | |
300 | // at this moment this routine is only defined for weights=1. | |
301 | Double_t x,x2,x3,x4,w,ww,m2,m4,n,s; | |
302 | ||
07025d6d | 303 | if(fW[0]!=(Double_t)fN||GetN()<4) return (-1.); |
b0f5e3fc | 304 | x = GetMeanYX(); // first order |
305 | x2 = GetYXNth(2); // second order | |
07025d6d | 306 | w = fW[0]; // first order - 1. |
307 | ww = fW[1]; // second order - 1. | |
b0f5e3fc | 308 | if(w*w==ww) return (-1.0); |
309 | s = (x2-x*x)*w*w/(w*w-ww); | |
310 | ||
311 | m2 = s; | |
312 | n = (Double_t) GetN(); | |
313 | x3 = GetYXNth(3); | |
314 | x4 = GetYXNth(4); | |
315 | // This equation assumes that all of the weights are equal to 1. | |
316 | m4 = (n/(n-1.))*(x4-3.*x*x3+6.*x*x*x2-2.*x*x*x*x); | |
317 | s = (m4-(n-3.)*m2*m2/(n-1.))/n; | |
318 | return TMath::Sqrt(s); | |
319 | } | |
320 | //_______________________________________________________________________ | |
321 | Double_t AliITSstatistics2::FitToLine(Double_t &a,Double_t &b){ | |
322 | // fit to y = a+bx returns Chi squared or -1.0 if an error | |
323 | Double_t c,d,e,f,g,h; | |
324 | ||
325 | a = b = 0.0; | |
326 | if(fOrder<2 || fN<3) return -1.0; | |
327 | c = GetWN(1); | |
328 | d = GetYN(1); | |
329 | e = GetXN(1); | |
330 | f = GetYXN(1); | |
331 | g = GetXN(2); | |
332 | h = c*g-e*e; | |
333 | a = d*g-f*e; | |
334 | b = c*f-d*e; | |
335 | if(h!=0.0){ | |
336 | a = a/h; | |
337 | b = b/h; | |
338 | }else{ | |
339 | printf("AliITSstatistics2: Error in FitToLine vertical line\n"); | |
340 | return -1.0; | |
341 | } // end if h | |
342 | h = GetYN(2)+a*a*c+b*b*g-2.0*a*d-2.0*b*f+2.0*a*b*e; | |
343 | h /= (Double_t)fN - 2.0; | |
344 | return h; | |
345 | } |