]>
Commit | Line | Data |
---|---|---|
7dacdf84 | 1 | // Object holding the Energy loss fit 'correction' |
7984e5f7 | 2 | // |
3 | // These are generated from Monte-Carlo or real ESDs. | |
7dacdf84 | 4 | // |
0bd4b00f | 5 | #include "AliFMDCorrELossFit.h" |
6 | #include "AliForwardUtil.h" | |
a19faec0 | 7 | #include "AliLandauGaus.h" |
0bd4b00f | 8 | #include <TF1.h> |
7095962e | 9 | #include <TGraph.h> |
0bd4b00f | 10 | #include <TBrowser.h> |
11 | #include <TVirtualPad.h> | |
12 | #include <THStack.h> | |
8449e3e0 | 13 | #include <TLatex.h> |
14 | #include <TLegend.h> | |
a19faec0 | 15 | #include <TLine.h> |
0bd4b00f | 16 | #include <TH1D.h> |
17 | #include <AliLog.h> | |
7dacdf84 | 18 | #include <TMath.h> |
0bd4b00f | 19 | #include <TList.h> |
20 | #include <iostream> | |
21 | #include <iomanip> | |
7095962e CHC |
22 | namespace { |
23 | Int_t fDebug = 1; | |
24 | } | |
0bd4b00f | 25 | |
81775aba | 26 | Double_t AliFMDCorrELossFit::ELossFit::fgMaxRelError = .25; |
27 | Double_t AliFMDCorrELossFit::ELossFit::fgLeastWeight = 1e-7; | |
0bd4b00f | 28 | Double_t AliFMDCorrELossFit::ELossFit::fgMaxChi2nu = 20; |
29 | ||
30 | //____________________________________________________________________ | |
31 | AliFMDCorrELossFit::ELossFit::ELossFit() | |
32 | : fN(0), | |
33 | fNu(0), | |
34 | fChi2(0), | |
35 | fC(0), | |
36 | fDelta(0), | |
37 | fXi(0), | |
38 | fSigma(0), | |
39 | fSigmaN(0), | |
40 | fA(0), | |
41 | fEC(0), | |
42 | fEDelta(0), | |
43 | fEXi(0), | |
44 | fESigma(0), | |
45 | fESigmaN(0), | |
46 | fEA(0), | |
47 | fQuality(0), | |
48 | fDet(0), | |
49 | fRing('\0'), | |
8449e3e0 | 50 | fBin(0), |
51 | fMaxWeight(0) | |
7984e5f7 | 52 | { |
53 | // | |
54 | // Default constructor | |
55 | // | |
56 | // | |
57 | } | |
0bd4b00f | 58 | //____________________________________________________________________ |
59 | AliFMDCorrELossFit::ELossFit::ELossFit(Int_t quality, const TF1& f) | |
a19faec0 | 60 | : fN(f.GetNpar() > AliLandauGaus::kN ? |
61 | Int_t(f.GetParameter(AliLandauGaus::kN)) : | |
0bd4b00f | 62 | 1), |
63 | fNu(f.GetNDF()), | |
64 | fChi2(f.GetChisquare()), | |
a19faec0 | 65 | fC(f.GetParameter(AliLandauGaus::kC)), |
66 | fDelta(f.GetParameter(AliLandauGaus::kDelta)), | |
67 | fXi(f.GetParameter(AliLandauGaus::kXi)), | |
68 | fSigma(f.GetParameter(AliLandauGaus::kSigma)), | |
69 | fSigmaN(f.GetParameter(AliLandauGaus::kSigmaN)), | |
0bd4b00f | 70 | fA(0), |
a19faec0 | 71 | fEC(f.GetParError(AliLandauGaus::kC)), |
72 | fEDelta(f.GetParError(AliLandauGaus::kDelta)), | |
73 | fEXi(f.GetParError(AliLandauGaus::kXi)), | |
74 | fESigma(f.GetParError(AliLandauGaus::kSigma)), | |
75 | fESigmaN(f.GetParError(AliLandauGaus::kSigmaN)), | |
0bd4b00f | 76 | fEA(0), |
77 | fQuality(quality), | |
78 | fDet(0), | |
79 | fRing('\0'), | |
8449e3e0 | 80 | fBin(0), |
81 | fMaxWeight(0) | |
0bd4b00f | 82 | { |
7984e5f7 | 83 | // |
84 | // Construct from a function | |
85 | // | |
86 | // Parameters: | |
87 | // quality Quality flag | |
88 | // f Function | |
89 | // | |
0bd4b00f | 90 | if (fN <= 0) return; |
91 | fA = new Double_t[fN]; | |
92 | fEA = new Double_t[fN]; | |
93 | for (Int_t i = 0; i < fN-1; i++) { | |
a19faec0 | 94 | fA[i] = f.GetParameter(AliLandauGaus::kA+i); |
95 | fEA[i] = f.GetParError(AliLandauGaus::kA+i); | |
0bd4b00f | 96 | } |
97 | fA[fN-1] = -9999; | |
98 | fEA[fN-1] = -9999; | |
99 | } | |
100 | ||
101 | //____________________________________________________________________ | |
102 | AliFMDCorrELossFit::ELossFit::ELossFit(Int_t quality,UShort_t n, | |
103 | Double_t chi2, UShort_t nu, | |
104 | Double_t c, Double_t ec, | |
105 | Double_t delta, Double_t edelta, | |
106 | Double_t xi, Double_t exi, | |
107 | Double_t sigma, Double_t esigma, | |
108 | Double_t sigman, Double_t esigman, | |
fb3430ac | 109 | const Double_t* a,const Double_t* ea) |
0bd4b00f | 110 | : fN(n), |
111 | fNu(nu), | |
112 | fChi2(chi2), | |
113 | fC(c), | |
114 | fDelta(delta), | |
115 | fXi(xi), | |
116 | fSigma(sigma), | |
117 | fSigmaN(sigman), | |
118 | fA(0), | |
119 | fEC(ec), | |
120 | fEDelta(edelta), | |
121 | fEXi(exi), | |
122 | fESigma(esigma), | |
123 | fESigmaN(esigman), | |
124 | fEA(0), | |
125 | fQuality(quality), | |
126 | fDet(0), | |
127 | fRing('\0'), | |
8449e3e0 | 128 | fBin(0), |
129 | fMaxWeight(0) | |
0bd4b00f | 130 | { |
7984e5f7 | 131 | // |
132 | // Constructor with full parameter set | |
133 | // | |
134 | // Parameters: | |
135 | // quality Quality flag | |
136 | // n @f$ N@f$ - Number of fitted peaks | |
137 | // chi2 @f$ \chi^2 @f$ | |
138 | // nu @f$ \nu @f$ - number degrees of freedom | |
139 | // c @f$ C@f$ - scale constant | |
140 | // ec @f$ \delta C@f$ - error on @f$ C@f$ | |
141 | // delta @f$ \Delta@f$ - Most probable value | |
142 | // edelta @f$ \delta\Delta@f$ - error on @f$\Delta@f$ | |
143 | // xi @f$ \xi@f$ - width | |
144 | // exi @f$ \delta\xi@f$ - error on @f$\xi@f$ | |
145 | // sigma @f$ \sigma@f$ - Width of Gaussian | |
146 | // esigma @f$ \delta\sigma@f$ - error on @f$\sigma@f$ | |
147 | // sigman @f$ \sigma_n@f$ - Noise width | |
148 | // esigman @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$ | |
149 | // a Array of @f$ N-1@f$ weights @f$ a_i@f$ for | |
150 | // @f$ i=2,\ldots@f$ | |
151 | // ea Array of @f$ N-1@f$ error on the weights @f$ a_i@f$ for | |
152 | // @f$ i=2,\ldots@f$ | |
153 | // | |
0bd4b00f | 154 | if (fN <= 0) return; |
155 | fA = new Double_t[fN]; | |
156 | fEA = new Double_t[fN]; | |
157 | for (Int_t i = 0; i < fN-1; i++) { | |
158 | fA[i] = a[i]; | |
159 | fEA[i] = ea[i]; | |
160 | } | |
161 | fA[fN-1] = -9999; | |
162 | fEA[fN-1] = -9999; | |
163 | } | |
164 | //____________________________________________________________________ | |
165 | AliFMDCorrELossFit::ELossFit::ELossFit(const ELossFit& o) | |
166 | : TObject(o), | |
167 | fN(o.fN), | |
168 | fNu(o.fNu), | |
169 | fChi2(o.fChi2), | |
170 | fC(o.fC), | |
171 | fDelta(o.fDelta), | |
172 | fXi(o.fXi), | |
173 | fSigma(o.fSigma), | |
174 | fSigmaN(o.fSigmaN), | |
175 | fA(0), | |
176 | fEC(o.fEC), | |
177 | fEDelta(o.fEDelta), | |
178 | fEXi(o.fEXi), | |
179 | fESigma(o.fESigma), | |
180 | fESigmaN(o.fESigmaN), | |
181 | fEA(0), | |
182 | fQuality(o.fQuality), | |
183 | fDet(o.fDet), | |
184 | fRing(o.fRing), | |
8449e3e0 | 185 | fBin(o.fBin), |
186 | fMaxWeight(o.fMaxWeight) | |
0bd4b00f | 187 | { |
7984e5f7 | 188 | // |
189 | // Copy constructor | |
190 | // | |
191 | // Parameters: | |
192 | // o Object to copy from | |
193 | // | |
0bd4b00f | 194 | if (fN <= 0) return; |
195 | fA = new Double_t[fN]; | |
196 | fEA = new Double_t[fN]; | |
197 | for (Int_t i = 0; i < fN-1; i++) { | |
198 | fA[i] = o.fA[i]; | |
199 | fEA[i] = o.fEA[i]; | |
200 | } | |
201 | fA[fN-1] = -9999; | |
202 | fEA[fN-1] = -9999; | |
203 | } | |
204 | ||
205 | //____________________________________________________________________ | |
206 | AliFMDCorrELossFit::ELossFit& | |
207 | AliFMDCorrELossFit::ELossFit::operator=(const ELossFit& o) | |
208 | { | |
7984e5f7 | 209 | // |
210 | // Assignment operator | |
211 | // | |
212 | // Parameters: | |
213 | // o Object to assign from | |
214 | // | |
215 | // Return: | |
216 | // Reference to this object | |
217 | // | |
d015ecfe | 218 | if (&o == this) return *this; |
8449e3e0 | 219 | fN = o.fN; |
220 | fNu = o.fNu; | |
221 | fChi2 = o.fChi2; | |
222 | fC = o.fC; | |
223 | fDelta = o.fDelta; | |
224 | fXi = o.fXi; | |
225 | fSigma = o.fSigma; | |
226 | fSigmaN = o.fSigmaN; | |
227 | fEC = o.fEC; | |
228 | fEDelta = o.fEDelta; | |
229 | fEXi = o.fEXi; | |
230 | fESigma = o.fESigma; | |
231 | fESigmaN = o.fESigmaN; | |
232 | fQuality = o.fQuality; | |
233 | fDet = o.fDet; | |
234 | fRing = o.fRing; | |
235 | fBin = o.fBin; | |
236 | fMaxWeight = o.fMaxWeight; | |
0bd4b00f | 237 | if (fA) delete [] fA; |
238 | if (fEA) delete [] fEA; | |
239 | fA = 0; | |
240 | fEA = 0; | |
241 | ||
242 | if (fN <= 0) return *this; | |
243 | fA = new Double_t[fN]; | |
244 | fEA = new Double_t[fN]; | |
245 | for (Int_t i = 0; i < fN; i++) { | |
246 | fA[i] = o.fA[i]; | |
247 | fEA[i] = o.fEA[i]; | |
248 | } | |
249 | ||
250 | return *this; | |
251 | } | |
252 | ||
253 | //____________________________________________________________________ | |
254 | AliFMDCorrELossFit::ELossFit::~ELossFit() | |
255 | { | |
256 | if (fA) delete[] fA; | |
257 | if (fEA) delete[] fEA; | |
258 | } | |
259 | ||
260 | ||
261 | //____________________________________________________________________ | |
262 | Int_t | |
263 | AliFMDCorrELossFit::ELossFit::FindMaxWeight(Double_t maxRelError, | |
264 | Double_t leastWeight, | |
265 | UShort_t maxN) const | |
266 | { | |
7984e5f7 | 267 | // |
268 | // Find the maximum weight to use. The maximum weight is the | |
269 | // largest i for which | |
270 | // | |
271 | // - @f$ i \leq \max{N}@f$ | |
272 | // - @f$ a_i > \min{a}@f$ | |
273 | // - @f$ \delta a_i/a_i > \delta_{max}@f$ | |
274 | // | |
275 | // Parameters: | |
276 | // maxRelError @f$ \min{a}@f$ | |
277 | // leastWeight @f$ \delta_{max}@f$ | |
278 | // maxN @f$ \max{N}@f$ | |
279 | // | |
280 | // Return: | |
281 | // The largest index @f$ i@f$ for which the above | |
282 | // conditions hold. Will never return less than 1. | |
283 | // | |
8449e3e0 | 284 | if (fMaxWeight > 0) return fMaxWeight; |
0bd4b00f | 285 | Int_t n = TMath::Min(maxN, UShort_t(fN-1)); |
286 | Int_t m = 1; | |
287 | // fN is one larger than we have data | |
288 | for (Int_t i = 0; i < n-1; i++, m++) { | |
289 | if (fA[i] < leastWeight) break; | |
290 | if (fEA[i] / fA[i] > maxRelError) break; | |
291 | } | |
8449e3e0 | 292 | return fMaxWeight = m; |
0bd4b00f | 293 | } |
294 | ||
295 | //____________________________________________________________________ | |
296 | Double_t | |
297 | AliFMDCorrELossFit::ELossFit::Evaluate(Double_t x, | |
298 | UShort_t maxN) const | |
299 | { | |
7984e5f7 | 300 | // |
301 | // Evaluate | |
302 | // @f[ | |
303 | // f_N(x;\Delta,\xi,\sigma') = | |
304 | // \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i') | |
305 | // @f] | |
306 | // | |
a19faec0 | 307 | // (see AliLandauGaus::NLandauGaus) for the maximum @f$ N @f$ |
7984e5f7 | 308 | // that fulfills the requirements |
309 | // | |
310 | // Parameters: | |
311 | // x Where to evaluate | |
312 | // maxN @f$ \max{N}@f$ | |
313 | // | |
314 | // Return: | |
315 | // @f$ f_N(x;\Delta,\xi,\sigma')@f$ | |
316 | // | |
a19faec0 | 317 | return AliLandauGaus::Fn(x, fDelta, fXi, fSigma, fSigmaN, |
318 | TMath::Min(maxN, UShort_t(fN)), fA); | |
0bd4b00f | 319 | } |
320 | ||
321 | //____________________________________________________________________ | |
322 | Double_t | |
323 | AliFMDCorrELossFit::ELossFit::EvaluateWeighted(Double_t x, | |
324 | UShort_t maxN) const | |
325 | { | |
7984e5f7 | 326 | // |
327 | // Evaluate | |
328 | // @f[ | |
329 | // f_W(x;\Delta,\xi,\sigma') = | |
330 | // \frac{\sum_{i=1}^{n} i a_i f_i(x;\Delta,\xi,\sigma')}{ | |
331 | // f_N(x;\Delta,\xi,\sigma')} = | |
332 | // \frac{\sum_{i=1}^{n} i a_i f(x;\Delta_i,\xi_i,\sigma_i')}{ | |
333 | // \sum_{i=1}^{n} a_i f(x;\Delta_i,\xi_i,\sigma_i')} | |
334 | // @f] | |
335 | // where @f$ n@f$ fulfills the requirements (see FindMaxWeight). | |
336 | // | |
337 | // If the denominator is zero, then 1 is returned. | |
338 | // | |
a19faec0 | 339 | // See also AliLandauGaus::Fi and AliLandauGaus::NLandauGaus |
7984e5f7 | 340 | // for more information on the evaluated functions. |
341 | // | |
342 | // Parameters: | |
343 | // x Where to evaluate | |
344 | // maxN @f$ \max{N}@f$ | |
345 | // | |
346 | // Return: | |
347 | // @f$ f_W(x;\Delta,\xi,\sigma')@f$. | |
348 | // | |
0bd4b00f | 349 | UShort_t n = TMath::Min(maxN, UShort_t(fN-1)); |
350 | Double_t num = 0; | |
351 | Double_t den = 0; | |
352 | for (Int_t i = 1; i <= n; i++) { | |
353 | Double_t a = (i == 1 ? 1 : fA[i-1]); | |
354 | if (fA[i-1] < 0) break; | |
a19faec0 | 355 | Double_t f = AliLandauGaus::Fi(x,fDelta,fXi,fSigma,fSigmaN,i); |
0bd4b00f | 356 | num += i * a * f; |
357 | den += a * f; | |
358 | } | |
359 | if (den <= 0) return 1; | |
360 | return num / den; | |
361 | } | |
362 | ||
363 | ||
364 | #define OUTPAR(N,V,E) \ | |
365 | std::setprecision(9) << \ | |
366 | std::setw(12) << N << ": " << \ | |
367 | std::setw(14) << V << " +/- " << \ | |
368 | std::setw(14) << E << " (" << \ | |
369 | std::setprecision(-1) << \ | |
370 | std::setw(5) << 100*(V>0?E/V:1) << "%)\n" | |
371 | ||
372 | ||
373 | //____________________________________________________________________ | |
374 | Int_t | |
375 | AliFMDCorrELossFit::ELossFit::Compare(const TObject* o) const | |
376 | { | |
7984e5f7 | 377 | // |
378 | // Compare to another ELossFit object. | |
379 | // | |
380 | // - +1, if this quality is better (larger) than other objects quality | |
381 | // - -1, if this quality is worse (smaller) than other objects quality | |
382 | // - +1, if this @f$|\chi^2/\nu-1|@f$ is smaller than the same for other | |
383 | // - -1, if this @f$|\chi^2/\nu-1|@f$ is larger than the same for other | |
384 | // - 0 otherwise | |
385 | // | |
386 | // Parameters: | |
387 | // o Other object to compare to | |
388 | // | |
0bd4b00f | 389 | const ELossFit* other = static_cast<const ELossFit*>(o); |
390 | if (this->fQuality < other->fQuality) return -1; | |
391 | if (this->fQuality > other->fQuality) return +1; | |
392 | Double_t chi2nu = (fNu == 0 ? 1000 : fChi2 / fNu); | |
393 | Double_t oChi2nu = (other->fNu == 0 ? 1000 : other->fChi2 / other->fNu); | |
394 | if (TMath::Abs(chi2nu-1) < TMath::Abs(oChi2nu-1)) return +1; | |
395 | if (TMath::Abs(chi2nu-1) > TMath::Abs(oChi2nu-1)) return -1; | |
396 | return 0; | |
397 | } | |
398 | ||
399 | //____________________________________________________________________ | |
400 | void | |
81775aba | 401 | AliFMDCorrELossFit::ELossFit::Print(Option_t* option) const |
0bd4b00f | 402 | { |
7984e5f7 | 403 | // |
404 | // Information to standard output | |
405 | // | |
406 | // Parameters: | |
407 | // option Not used | |
408 | // | |
81775aba | 409 | TString o(option); |
410 | if (o.Contains("S", TString::kIgnoreCase)) { | |
411 | Printf("%15s: q=%2d n=%1d chi2/nu=%6.3f", | |
412 | GetName(), fQuality, fN, (fNu <= 0 ? 999 : fChi2 / fNu)); | |
413 | return; | |
414 | } | |
415 | ||
0bd4b00f | 416 | std::cout << GetName() << ":\n" |
417 | << " chi^2/nu = " << fChi2 << "/" << fNu << " = " | |
418 | << (fNu == 0 ? 999 : fChi2 / fNu) << "\n" | |
419 | << " Quality: " << fQuality << "\n" | |
420 | << " NParticles: " << fN << " (" << FindMaxWeight() << ")\n" | |
a19faec0 | 421 | << OUTPAR("Delta", fDelta, fEDelta) |
422 | << OUTPAR("xi", fXi, fEXi) | |
423 | << OUTPAR("sigma", fSigma, fESigma) | |
0bd4b00f | 424 | << OUTPAR("sigma_n", fSigmaN, fESigmaN); |
425 | for (Int_t i = 0; i < fN-1; i++) | |
426 | std::cout << OUTPAR(Form("a%d", i+2), fA[i], fEA[i]); | |
427 | std::cout << std::flush; | |
428 | } | |
7095962e CHC |
429 | //____________________________________________________________________ |
430 | TF1* | |
431 | AliFMDCorrELossFit::ELossFit::GetF1(Int_t i, Double_t max) const | |
432 | { | |
433 | const Double_t lowX = 0.001; | |
434 | const Double_t upX = (max < 0 ? 10 : max); | |
435 | Int_t maxW = FindMaxWeight(); | |
436 | TF1* ret = 0; | |
437 | ||
438 | if (i <= 0) | |
a19faec0 | 439 | ret = AliLandauGaus::MakeFn(fC * 1, fDelta, fXi, |
440 | fSigma, fSigmaN, maxW/*fN*/, fA, lowX, upX); | |
441 | else if (i == 1) | |
442 | ret = AliLandauGaus::MakeF1(fC, fDelta, fXi, fSigma, fSigmaN, lowX, upX); | |
7095962e | 443 | else if (i <= maxW) |
a19faec0 | 444 | ret = AliLandauGaus::MakeFi(fC*(i == 1 ? 1 : fA[i-2]), |
445 | fDelta, fXi, fSigma, fSigmaN, i, lowX, upX); | |
446 | ||
7095962e CHC |
447 | return ret; |
448 | } | |
449 | //____________________________________________________________________ | |
450 | Double_t | |
451 | AliFMDCorrELossFit::ELossFit::FindProbabilityCut(Double_t low) const | |
452 | { | |
453 | Double_t ret = 1000; | |
454 | TF1* f = 0; | |
455 | TGraph* g = 0; | |
456 | try { | |
7b2fa237 | 457 | if (!(f = GetF1(1,5))) // First landau up to Delta/Delta_{mip}=5 |
7095962e CHC |
458 | throw TString("Didn't TF1 object"); |
459 | if (!(g = new TGraph(f, "i"))) | |
460 | throw TString("Failed to integrate function"); | |
461 | ||
462 | Int_t n = g->GetN(); | |
463 | Double_t total = g->GetY()[n-1]; | |
464 | if (total <= 0) | |
465 | throw TString::Format("Invalid integral: %lf", total); | |
466 | ||
467 | for (Int_t i = 0; i < n; i++) { | |
468 | Double_t prob = g->GetY()[i] / total; | |
469 | if (prob > low) { | |
470 | ret = g->GetX()[i]; | |
471 | break; | |
472 | } | |
473 | } | |
474 | if (ret >= 1000) | |
475 | throw TString::Format("Couldn't find x value for cut %lf", low); | |
476 | } | |
477 | catch (const TString& str) { | |
478 | AliWarningF("%s: %s", GetName(), str.Data()); | |
479 | } | |
480 | if (f) delete f; | |
481 | if (g) delete g; | |
482 | return ret; | |
483 | } | |
484 | ||
0bd4b00f | 485 | |
486 | //____________________________________________________________________ | |
487 | const Char_t* | |
488 | AliFMDCorrELossFit::ELossFit::GetName() const | |
489 | { | |
7984e5f7 | 490 | // |
491 | // Get the name of this object | |
492 | // | |
493 | // | |
494 | // Return: | |
495 | // | |
496 | // | |
0bd4b00f | 497 | return Form("FMD%d%c_etabin%03d", fDet, fRing, fBin); |
498 | } | |
499 | ||
500 | //____________________________________________________________________ | |
501 | void | |
502 | AliFMDCorrELossFit::ELossFit::Browse(TBrowser* b) | |
503 | { | |
7984e5f7 | 504 | // |
505 | // Browse this object | |
506 | // | |
507 | // Parameters: | |
508 | // b Browser | |
509 | // | |
8449e3e0 | 510 | Draw(b ? b->GetDrawOption() : "comp values"); |
0bd4b00f | 511 | gPad->SetLogy(); |
512 | gPad->Update(); | |
513 | } | |
514 | ||
515 | //____________________________________________________________________ | |
516 | void | |
517 | AliFMDCorrELossFit::ELossFit::Draw(Option_t* option) | |
518 | { | |
7984e5f7 | 519 | // |
520 | // Draw this fit | |
521 | // | |
522 | // Parameters: | |
523 | // option Options | |
524 | // - COMP Draw components too | |
525 | // | |
0bd4b00f | 526 | TString opt(option); |
527 | opt.ToUpper(); | |
528 | bool comp = false; | |
8449e3e0 | 529 | bool good = false; |
530 | bool vals = false; | |
531 | bool legd = false; | |
a19faec0 | 532 | bool peak = false; |
0bd4b00f | 533 | if (opt.Contains("COMP")) { |
534 | opt.ReplaceAll("COMP",""); | |
535 | comp = true; | |
536 | } | |
8449e3e0 | 537 | if (opt.Contains("GOOD")) { |
538 | opt.ReplaceAll("GOOD",""); | |
539 | good = true; | |
540 | } | |
541 | if (opt.Contains("VALUES")) { | |
542 | opt.ReplaceAll("VALUES",""); | |
543 | vals = true; | |
544 | } | |
545 | if (opt.Contains("LEGEND")) { | |
546 | opt.ReplaceAll("LEGEND",""); | |
547 | legd = comp; | |
548 | } | |
0bd4b00f | 549 | if (!opt.Contains("SAME")) { |
550 | gPad->Clear(); | |
551 | } | |
a19faec0 | 552 | if (opt.Contains("PEAK")) { |
553 | peak = true; | |
554 | } | |
8449e3e0 | 555 | TLegend* l = 0; |
556 | if (legd) { | |
557 | l = new TLegend(.3, .5, .59, .94); | |
558 | l->SetBorderSize(0); | |
559 | l->SetFillColor(0); | |
560 | l->SetFillStyle(0); | |
561 | } | |
0bd4b00f | 562 | TObjArray cleanup; |
8449e3e0 | 563 | Int_t maxW = FindMaxWeight(); |
a19faec0 | 564 | TF1* tot = AliLandauGaus::MakeFn(fC * 1, fDelta, fXi, fSigma, fSigmaN, |
565 | maxW/*fN*/, fA, 0.01, 10); | |
0bd4b00f | 566 | tot->SetLineColor(kBlack); |
567 | tot->SetLineWidth(2); | |
568 | tot->SetLineStyle(1); | |
569 | tot->SetTitle(GetName()); | |
8449e3e0 | 570 | if (l) l->AddEntry(tot, "Total", "l"); |
0bd4b00f | 571 | Double_t max = tot->GetMaximum(); |
572 | ||
8449e3e0 | 573 | |
0bd4b00f | 574 | if (!opt.Contains("SAME")) { |
575 | TH1* frame = new TH1F(GetName(), | |
e65b8b56 | 576 | Form("FMD%d%c, eta bin %d",fDet,fRing,fBin), |
577 | 100, 0, 10); | |
0bd4b00f | 578 | frame->SetMinimum(max/10000); |
579 | frame->SetMaximum(max*1.4); | |
580 | frame->SetXTitle("#Delta / #Delta_{mip}"); | |
581 | frame->Draw(); | |
582 | opt.Append(" SAME"); | |
583 | } | |
e65b8b56 | 584 | TF1* cpy = tot->DrawCopy(opt.Data()); |
0bd4b00f | 585 | cleanup.Add(tot); |
586 | ||
8449e3e0 | 587 | if (vals) { |
588 | Double_t x1 = .72; | |
589 | Double_t x2 = .73; | |
590 | Double_t y = .90; | |
591 | Double_t dy = .05; | |
592 | TLatex* ltx1 = new TLatex(x1, y, ""); | |
593 | TLatex* ltx2 = new TLatex(x2, y, ""); | |
594 | ltx1->SetNDC(); | |
595 | ltx1->SetTextAlign(33); | |
596 | ltx1->SetTextFont(132); | |
597 | ltx1->SetTextSize(dy-.01); | |
598 | ltx2->SetNDC(); | |
599 | ltx2->SetTextAlign(13); | |
600 | ltx2->SetTextFont(132); | |
601 | ltx2->SetTextSize(dy-.01); | |
602 | ||
603 | ltx1->DrawLatex(x1, y, "Quality"); | |
604 | ltx2->DrawLatex(x2, y, Form("%d", fQuality)); | |
605 | y -= dy; | |
606 | ||
607 | ltx1->DrawLatex(x1, y, "#chi^{2}/#nu"); | |
608 | ltx2->DrawLatex(x2, y, Form("%7.3f", (fNu > 0 ? fChi2 / fNu : -1))); | |
609 | y -= dy; | |
610 | ||
611 | const Char_t* pn[] = { "C", "#Delta", "#xi", "#sigma" }; | |
612 | Double_t pv[] = { fC, fDelta, fXi, fSigma }; | |
613 | Double_t pe[] = { fEC, fEDelta, fEXi, fESigma }; | |
614 | for (Int_t i = 0; i < 4; i++) { | |
615 | ltx1->DrawLatex(x1, y, pn[i]); | |
616 | ltx2->DrawLatex(x2, y, Form("%6.4f#pm%6.4f", pv[i], pe[i])); | |
617 | y -= dy; | |
618 | } | |
619 | for (Int_t i=2; i <= fN; i++) { | |
620 | if (i > maxW) { | |
621 | ltx1->SetTextColor(kRed+3); | |
622 | ltx2->SetTextColor(kRed+3); | |
623 | } | |
624 | ltx1->DrawLatex(x1, y, Form("a_{%d}", i)); | |
625 | ltx2->DrawLatex(x2, y, Form("%6.4f#pm%6.4f", fA[i-2], fEA[i-2])); | |
626 | y -= dy; | |
627 | } | |
628 | ||
629 | } | |
630 | ||
a19faec0 | 631 | if (peak) { |
632 | TLine* pl = new TLine(fDelta, 0.01*max, fDelta, 1.5*max); | |
633 | pl->SetLineStyle(2); | |
634 | pl->SetLineColor(kBlack); | |
635 | pl->Draw(); | |
636 | } | |
0bd4b00f | 637 | if (!comp) { |
638 | gPad->cd(); | |
639 | return; | |
640 | } | |
641 | ||
642 | Double_t min = max; | |
643 | opt.Append(" same"); | |
0bd4b00f | 644 | for (Int_t i=1; i <= fN; i++) { |
8449e3e0 | 645 | if (good && i > maxW) break; |
a19faec0 | 646 | TF1* f = AliLandauGaus::MakeFi(fC*(i == 1 ? 1 : fA[i-2]), |
0bd4b00f | 647 | fDelta, fXi, |
648 | fSigma, fSigmaN, | |
649 | i, 0.01, 10); | |
650 | f->SetLineWidth(2); | |
651 | f->SetLineStyle(i > maxW ? 2 : 1); | |
652 | min = TMath::Min(f->GetMaximum(), min); | |
653 | f->DrawCopy(opt.Data()); | |
8449e3e0 | 654 | if (l) l->AddEntry(f, Form("%d MIP%s", i, (i>1 ? "s" : "")), "l"); |
655 | ||
0bd4b00f | 656 | cleanup.Add(f); |
657 | } | |
658 | min /= 100; | |
e65b8b56 | 659 | if (max <= 0) max = 0.1; |
660 | if (min <= 0) min = 1e-4; | |
661 | cpy->SetMaximum(max); | |
662 | cpy->SetMinimum(min); | |
663 | cpy->GetHistogram()->SetMaximum(max); | |
664 | cpy->GetHistogram()->SetMinimum(min); | |
665 | cpy->GetHistogram()->GetYaxis()->SetRangeUser(min, max); | |
8449e3e0 | 666 | if (l) l->Draw(); |
0bd4b00f | 667 | |
668 | gPad->cd(); | |
669 | } | |
670 | ||
671 | ||
672 | //____________________________________________________________________ | |
673 | #define CHECKPAR(V,E,T) ((V > 0) && (E / V < T)) | |
674 | ||
16ac05a3 | 675 | //____________________________________________________________________ |
676 | Double_t | |
677 | AliFMDCorrELossFit::ELossFit::GetLowerBound(Double_t f) const | |
678 | { | |
679 | // | |
680 | // Return | |
681 | // Delta * f | |
682 | return f * fDelta; | |
683 | } | |
5bb5d1f6 | 684 | //____________________________________________________________________ |
685 | Double_t | |
686 | AliFMDCorrELossFit::ELossFit::GetLowerBound(Double_t f, | |
687 | Bool_t includeSigma) const | |
688 | { | |
689 | // | |
690 | // Return | |
691 | // Delta - f * (xi + sigma) | |
692 | return fDelta - f * (fXi + (includeSigma ? fSigma : 0)); | |
693 | } | |
694 | ||
0bd4b00f | 695 | //____________________________________________________________________ |
696 | void | |
697 | AliFMDCorrELossFit::ELossFit::CalculateQuality(Double_t maxChi2nu, | |
698 | Double_t maxRelError, | |
699 | Double_t leastWeight) | |
700 | { | |
7984e5f7 | 701 | // |
702 | // Calculate the quality | |
703 | // | |
81775aba | 704 | Double_t decline = maxChi2nu; |
0bd4b00f | 705 | Int_t qual = 0; |
81775aba | 706 | if (fNu > 0) { |
707 | Double_t red = fChi2 / fNu; | |
708 | if (red < maxChi2nu) qual += 4; | |
709 | else { | |
77f97e3f | 710 | Int_t q = Int_t((maxChi2nu+decline - red) / decline * 4); |
81775aba | 711 | if (q > 0) qual += q; |
712 | } | |
713 | } | |
0bd4b00f | 714 | if (CHECKPAR(fDelta, fEDelta, maxRelError)) qual++; |
715 | if (CHECKPAR(fXi, fEXi, maxRelError)) qual++; | |
716 | if (CHECKPAR(fSigma, fESigma, maxRelError)) qual++; | |
717 | if (CHECKPAR(fSigmaN, fESigmaN, maxRelError)) qual++; | |
e65b8b56 | 718 | // Large penalty for large sigma - often a bad fit. |
719 | if (fSigma > 10*fXi) qual -= 4; | |
81775aba | 720 | qual += FindMaxWeight(2*maxRelError, leastWeight, fN); |
0bd4b00f | 721 | fQuality = qual; |
722 | } | |
723 | ||
724 | //____________________________________________________________________ | |
725 | AliFMDCorrELossFit::AliFMDCorrELossFit() | |
726 | : TObject(), | |
727 | fRings(), | |
728 | fEtaAxis(0,0,0), | |
8449e3e0 | 729 | fLowCut(0), |
730 | fCache(0) | |
0bd4b00f | 731 | { |
7984e5f7 | 732 | // |
733 | // Default constructor | |
734 | // | |
0bd4b00f | 735 | fRings.SetOwner(kTRUE); |
736 | fEtaAxis.SetTitle("#eta"); | |
737 | fEtaAxis.SetName("etaAxis"); | |
738 | fRings.SetName("rings"); | |
739 | } | |
740 | ||
741 | //____________________________________________________________________ | |
742 | AliFMDCorrELossFit::AliFMDCorrELossFit(const AliFMDCorrELossFit& o) | |
743 | : TObject(o), | |
744 | fRings(o.fRings), | |
745 | fEtaAxis(o.fEtaAxis.GetNbins(),o.fEtaAxis.GetXmin(),o.fEtaAxis.GetXmax()), | |
8449e3e0 | 746 | fLowCut(0), |
747 | fCache(0) | |
0bd4b00f | 748 | { |
7984e5f7 | 749 | // |
750 | // Copy constructor | |
751 | // | |
752 | // Parameters: | |
753 | // o Object to copy from | |
754 | // | |
0bd4b00f | 755 | fEtaAxis.SetTitle("#eta"); |
756 | fEtaAxis.SetName("etaAxis"); | |
757 | } | |
758 | ||
759 | //____________________________________________________________________ | |
760 | AliFMDCorrELossFit::~AliFMDCorrELossFit() | |
761 | { | |
7984e5f7 | 762 | // |
763 | // Destructor | |
764 | // | |
0bd4b00f | 765 | fRings.Clear(); |
766 | } | |
767 | ||
768 | //____________________________________________________________________ | |
769 | AliFMDCorrELossFit& | |
770 | AliFMDCorrELossFit::operator=(const AliFMDCorrELossFit& o) | |
771 | { | |
7984e5f7 | 772 | // |
773 | // Assignment operator | |
774 | // | |
775 | // Parameters: | |
776 | // o Object to assign from | |
777 | // | |
778 | // Return: | |
779 | // Reference to this object | |
780 | // | |
d015ecfe | 781 | if (&o == this) return *this; |
0bd4b00f | 782 | fRings = o.fRings; |
783 | fLowCut = o.fLowCut; | |
8449e3e0 | 784 | fCache = o.fCache; |
0bd4b00f | 785 | SetEtaAxis(o.fEtaAxis.GetNbins(), o.fEtaAxis.GetXmin(), o.fEtaAxis.GetXmax()); |
786 | ||
787 | return *this; | |
788 | } | |
8449e3e0 | 789 | #define CACHE(BIN,IDX,OFFSET) fCache[IDX*OFFSET+BIN-1] |
7b2fa237 | 790 | #define CACHEDR(BIN,D,R,OFFSET) \ |
791 | CACHE(BIN,(D == 1 ? 0 : (D - 2) * 2 + 1 + (R=='I' || R=='i' ? 0 : 1)),OFFSET) | |
8449e3e0 | 792 | |
793 | //____________________________________________________________________ | |
794 | void | |
795 | AliFMDCorrELossFit::CacheBins(UShort_t minQuality) const | |
796 | { | |
a19faec0 | 797 | AliLandauGaus::EnableSigmaShift(TestBit(kHasShift)); |
8449e3e0 | 798 | if (fCache.GetSize() > 0) return; |
799 | ||
800 | Int_t nRings = fRings.GetEntriesFast(); | |
801 | Int_t offset = fEtaAxis.GetNbins(); | |
802 | ||
803 | fCache.Set(nRings*offset); | |
804 | fCache.Reset(-1); | |
805 | ||
806 | for (Int_t i = 0; i < nRings; i++) { | |
807 | TObjArray* ringArray = static_cast<TObjArray*>(fRings.At(i)); | |
808 | ||
809 | // First loop to find where we actually have fits | |
810 | Int_t nFits = 0; | |
811 | Int_t nGood = 0; | |
812 | Int_t minBin = offset+1; | |
813 | Int_t maxBin = -1; | |
814 | Int_t realMinBin = offset+1; | |
815 | Int_t realMaxBin = -1; | |
816 | for (Int_t j = 1; j < ringArray->GetEntriesFast(); j++) { | |
817 | ELossFit* fit = static_cast<ELossFit*>(ringArray->At(j)); | |
818 | if (!fit) continue; | |
819 | nFits++; | |
820 | ||
821 | // Update our range off possible fits | |
822 | realMinBin = TMath::Min(j, realMinBin); | |
823 | realMaxBin = TMath::Max(j, realMaxBin); | |
824 | ||
825 | // Check the quality of the fit | |
81775aba | 826 | fit->CalculateQuality(AliFMDCorrELossFit::ELossFit::fgMaxChi2nu, |
827 | AliFMDCorrELossFit::ELossFit::fgMaxRelError, | |
828 | AliFMDCorrELossFit::ELossFit::fgLeastWeight); | |
8449e3e0 | 829 | if (minQuality > 0 && fit->fQuality < minQuality) continue; |
830 | nGood++; | |
831 | ||
832 | // Check this bin | |
833 | CACHE(j,i,offset) = j; | |
834 | minBin = TMath::Min(j, minBin); | |
835 | maxBin = TMath::Max(j, maxBin); | |
836 | } | |
837 | AliInfoF("Out of %d bins, %d had fits, of which %d are good (%5.1f%%)", | |
997ba0f4 | 838 | offset, nFits, nGood, (nFits > 0 ? 100*float(nGood)/nFits : 0)); |
839 | ||
8449e3e0 | 840 | // Now loop and set neighbors |
841 | realMinBin = TMath::Max(1, realMinBin-1); // Include one more | |
842 | realMaxBin = TMath::Min(offset, realMaxBin+1); // Include one more | |
843 | for (Int_t j = realMinBin; j <= realMaxBin; j++) { | |
844 | if (CACHE(j,i,offset) > 0) continue; | |
845 | ||
846 | Int_t nK = TMath::Max(realMaxBin - j, j - realMinBin); | |
847 | Int_t found = -1; | |
848 | for (Int_t k = 1; k <= nK; k++) { | |
849 | Int_t left = j - k; | |
850 | Int_t right = j + k; | |
851 | if (left > realMinBin && | |
852 | CACHE(left,i,offset) == left) found = left; | |
853 | else if (right < realMaxBin && | |
854 | CACHE(right,i,offset) == right) found = right; | |
855 | if (found > 0) break; | |
856 | } | |
857 | // Now check that we found something | |
858 | if (found) CACHE(j,i,offset) = CACHE(found,i,offset); | |
859 | else AliWarningF("No fit found for etabin=%d in ring=%d", j, i); | |
860 | } | |
861 | } | |
862 | } | |
863 | ||
864 | ||
0bd4b00f | 865 | //____________________________________________________________________ |
866 | Int_t | |
867 | AliFMDCorrELossFit::FindEtaBin(Double_t eta) const | |
868 | { | |
7984e5f7 | 869 | // |
870 | // Find the eta bin corresponding to the given eta | |
871 | // | |
872 | // Parameters: | |
873 | // eta Eta value | |
874 | // | |
875 | // Return: | |
876 | // Bin (in @f$[1,N_{bins}]@f$) corresponding to the given | |
877 | // eta, or 0 if out of range. | |
878 | // | |
fb3430ac | 879 | if (TMath::Abs(fEtaAxis.GetXmin() - fEtaAxis.GetXmax()) < 1e-6 |
880 | || fEtaAxis.GetNbins() == 0) { | |
0bd4b00f | 881 | AliWarning("No eta axis defined"); |
882 | return -1; | |
883 | } | |
884 | Int_t bin = const_cast<TAxis&>(fEtaAxis).FindBin(eta); | |
885 | if (bin <= 0 || bin > fEtaAxis.GetNbins()) return 0; | |
886 | return bin; | |
887 | } | |
888 | ||
889 | //____________________________________________________________________ | |
890 | Bool_t | |
891 | AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, Int_t etaBin, ELossFit* fit) | |
892 | { | |
7984e5f7 | 893 | // |
894 | // Set the fit parameters from a function | |
895 | // | |
896 | // Parameters: | |
897 | // d Detector | |
898 | // r Ring | |
899 | // etaBin Eta (bin number, 1->nBins) | |
900 | // f ELoss fit result - note, the object will take ownership | |
901 | // | |
0bd4b00f | 902 | TObjArray* ringArray = GetOrMakeRingArray(d, r); |
903 | if (!ringArray) { | |
904 | AliError(Form("Failed to make ring array for FMD%d%c", d, r)); | |
905 | return kFALSE; | |
906 | } | |
907 | if (etaBin <= 0 || etaBin >= fEtaAxis.GetNbins()+1) { | |
908 | AliError(Form("bin=%d is out of range [%d,%d]", | |
909 | etaBin, 1, fEtaAxis.GetNbins())); | |
910 | return kFALSE; | |
911 | } | |
912 | // AliInfo(Form("Adding fit %p at %3d", fit, etaBin)); | |
913 | ringArray->AddAtAndExpand(fit, etaBin); | |
914 | return kTRUE; | |
915 | } | |
916 | ||
917 | //____________________________________________________________________ | |
918 | Bool_t | |
919 | AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, Double_t eta, ELossFit* fit) | |
920 | { | |
7984e5f7 | 921 | // |
922 | // Set the fit parameters from a function | |
923 | // | |
924 | // Parameters: | |
925 | // d Detector | |
926 | // r Ring | |
927 | // eta Eta | |
928 | // f ELoss fit result - note, the object will take ownership | |
929 | // | |
0bd4b00f | 930 | Int_t bin = FindEtaBin(eta); |
931 | if (bin <= 0) { | |
932 | AliError(Form("eta=%f is out of range [%f,%f]", | |
933 | eta, fEtaAxis.GetXmin(), fEtaAxis.GetXmax())); | |
934 | return kFALSE; | |
935 | } | |
936 | ||
937 | return SetFit(d, r, bin, fit); | |
938 | } | |
939 | //____________________________________________________________________ | |
940 | Bool_t | |
941 | AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, | |
942 | Double_t eta, | |
943 | Int_t quality,UShort_t n, | |
944 | Double_t chi2, UShort_t nu, | |
945 | Double_t c, Double_t ec, | |
946 | Double_t delta, Double_t edelta, | |
947 | Double_t xi, Double_t exi, | |
948 | Double_t sigma, Double_t esigma, | |
949 | Double_t sigman, Double_t esigman, | |
950 | Double_t* a, Double_t* ea) | |
951 | { | |
7984e5f7 | 952 | // |
953 | // Set the fit parameters from a function | |
954 | // | |
955 | // Parameters: | |
956 | // d Detector number | |
957 | // r Ring identifier | |
958 | // eta Eta value | |
959 | // quality Quality flag | |
960 | // n @f$ N@f$ - Number of fitted peaks | |
961 | // chi2 @f$ \chi^2 @f$ | |
962 | // nu @f$ \nu @f$ - number degrees of freedom | |
963 | // c @f$ C@f$ - scale constant | |
964 | // ec @f$ \delta C@f$ - error on @f$ C@f$ | |
965 | // delta @f$ \Delta@f$ - most probable value | |
966 | // edelta @f$ \delta\Delta@f$ - error on @f$\Delta@f$ | |
967 | // xi @f$ \xi@f$ - Landau width | |
968 | // exi @f$ \delta\xi@f$ - error on @f$\xi@f$ | |
969 | // sigma @f$ \sigma@f$ - Gaussian width | |
970 | // esigma @f$ \delta\sigma@f$ - error on @f$\sigma@f$ | |
971 | // sigman @f$ \sigma_n@f$ - Noise width | |
972 | // esigman @f$ \delta\sigma_n@f$ - error on @f$\sigma_n@f$ | |
973 | // a Array of @f$ N-1@f$ weights @f$ a_i@f$ for | |
974 | // @f$ i=2,\ldots@f$ | |
975 | // ea Array of @f$ N-1@f$ errors on weights @f$ a_i@f$ for | |
976 | // @f$ i=2,\ldots@f$ | |
977 | // | |
0bd4b00f | 978 | ELossFit* e = new ELossFit(quality, n, |
979 | chi2, nu, | |
980 | c, ec, | |
981 | delta, edelta, | |
982 | xi, exi, | |
983 | sigma, esigma, | |
984 | sigman, esigman, | |
985 | a, ea); | |
986 | if (!SetFit(d, r, eta, e)) { | |
987 | delete e; | |
988 | return kFALSE; | |
989 | } | |
990 | return kTRUE; | |
991 | } | |
992 | //____________________________________________________________________ | |
993 | Bool_t | |
994 | AliFMDCorrELossFit::SetFit(UShort_t d, Char_t r, Double_t eta, | |
995 | Int_t quality, const TF1& f) | |
996 | { | |
7984e5f7 | 997 | // |
998 | // Set the fit parameters from a function | |
999 | // | |
1000 | // Parameters: | |
1001 | // d Detector | |
1002 | // r Ring | |
1003 | // eta Eta | |
1004 | // quality Quality flag | |
1005 | // f Function from fit | |
1006 | // | |
0bd4b00f | 1007 | ELossFit* e = new ELossFit(quality, f); |
1008 | if (!SetFit(d, r, eta, e)) { | |
1009 | delete e; | |
1010 | return kFALSE; | |
1011 | } | |
1012 | return kTRUE; | |
1013 | } | |
1174780f | 1014 | //____________________________________________________________________ |
1015 | AliFMDCorrELossFit::ELossFit* | |
1016 | AliFMDCorrELossFit::GetFit(UShort_t d, Char_t r, Int_t etabin) const | |
1017 | { | |
1018 | // | |
1019 | // Get the fit corresponding to the specified parameters | |
1020 | // | |
1021 | // Parameters: | |
1022 | // d Detector | |
1023 | // r Ring | |
1024 | // etabin Eta bin (1 based) | |
1025 | // | |
1026 | // Return: | |
1027 | // Fit parameters or null in case of problems | |
1028 | // | |
1029 | TObjArray* ringArray = GetRingArray(d, r); | |
1030 | if (!ringArray) return 0; | |
1031 | if (etabin <= 0 || etabin >= fEtaAxis.GetNbins()) return 0; | |
1032 | if (etabin > ringArray->GetEntriesFast()) return 0; | |
1033 | else if (etabin >= ringArray->GetEntriesFast()) etabin--; | |
1034 | else if (!ringArray->At(etabin)) etabin++; | |
1035 | return static_cast<ELossFit*>(ringArray->At(etabin)); | |
1036 | } | |
1037 | //____________________________________________________________________ | |
1038 | AliFMDCorrELossFit::ELossFit* | |
1039 | AliFMDCorrELossFit::GetFit(UShort_t d, Char_t r, Double_t eta) const | |
1040 | { | |
1041 | // | |
1042 | // Find the fit corresponding to the specified parameters | |
1043 | // | |
1044 | // Parameters: | |
1045 | // d Detector | |
1046 | // r Ring | |
1047 | // eta Eta value | |
1048 | // | |
1049 | // Return: | |
1050 | // Fit parameters or null in case of problems | |
1051 | // | |
1052 | Int_t etabin = FindEtaBin(eta); | |
1053 | return GetFit(d, r, etabin); | |
1054 | } | |
1055 | ||
0bd4b00f | 1056 | //____________________________________________________________________ |
1057 | AliFMDCorrELossFit::ELossFit* | |
8449e3e0 | 1058 | AliFMDCorrELossFit::FindFit(UShort_t d, Char_t r, Int_t etabin, |
1059 | UShort_t minQ) const | |
0bd4b00f | 1060 | { |
7984e5f7 | 1061 | // |
1062 | // Find the fit corresponding to the specified parameters | |
1063 | // | |
1064 | // Parameters: | |
1065 | // d Detector | |
1066 | // r Ring | |
1067 | // etabin Eta bin (1 based) | |
1068 | // | |
1069 | // Return: | |
1070 | // Fit parameters or null in case of problems | |
1071 | // | |
0bd4b00f | 1072 | if (etabin <= 0 || etabin >= fEtaAxis.GetNbins()) { |
5bb5d1f6 | 1073 | // AliError(Form("Eta bin=%3d out of bounds [%d,%d] for FMD%d%c", |
1074 | // etabin, 1, fEtaAxis.GetNbins(), d, r)); | |
0bd4b00f | 1075 | return 0; |
1076 | } | |
8449e3e0 | 1077 | |
1078 | TObjArray* ringArray = GetRingArray(d, r); | |
1079 | if (!ringArray) { | |
1080 | AliError(Form("Failed to make ring array for FMD%d%c", d, r)); | |
0bd4b00f | 1081 | return 0; |
1082 | } | |
7095962e | 1083 | DMSG(fDebug, 10, "Got ringArray %s for FMD%d%c", ringArray->GetName(), d, r); |
8449e3e0 | 1084 | if (fCache.GetSize() <= 0) CacheBins(minQ); |
7b2fa237 | 1085 | #if 0 |
8449e3e0 | 1086 | Int_t idx = (d == 1 ? 0 : |
1087 | (d - 2) * 2 + 1 + (r=='I' || r=='i' ? 0 : 1)); | |
7b2fa237 | 1088 | #endif |
1089 | Int_t bin = CACHEDR(etabin, d, r, fEtaAxis.GetNbins()); | |
8449e3e0 | 1090 | |
1091 | if (bin < 0 || bin > ringArray->GetEntriesFast()) return 0; | |
1092 | ||
1093 | return static_cast<ELossFit*>(ringArray->At(bin)); | |
0bd4b00f | 1094 | } |
1095 | //____________________________________________________________________ | |
1096 | AliFMDCorrELossFit::ELossFit* | |
8449e3e0 | 1097 | AliFMDCorrELossFit::FindFit(UShort_t d, Char_t r, Double_t eta, |
1098 | UShort_t minQ) const | |
0bd4b00f | 1099 | { |
7984e5f7 | 1100 | // |
1101 | // Find the fit corresponding to the specified parameters | |
1102 | // | |
1103 | // Parameters: | |
1104 | // d Detector | |
1105 | // r Ring | |
1106 | // eta Eta value | |
1107 | // | |
1108 | // Return: | |
1109 | // Fit parameters or null in case of problems | |
1110 | // | |
0bd4b00f | 1111 | Int_t etabin = FindEtaBin(eta); |
8449e3e0 | 1112 | return FindFit(d, r, etabin, minQ); |
0bd4b00f | 1113 | } |
1114 | //____________________________________________________________________ | |
1115 | TObjArray* | |
1116 | AliFMDCorrELossFit::GetRingArray(UShort_t d, Char_t r) const | |
1117 | { | |
7984e5f7 | 1118 | // |
1119 | // Get the ring array corresponding to the specified ring | |
1120 | // | |
1121 | // Parameters: | |
1122 | // d Detector | |
1123 | // r Ring | |
1124 | // | |
1125 | // Return: | |
1126 | // Pointer to ring array, or null in case of problems | |
1127 | // | |
0bd4b00f | 1128 | Int_t idx = -1; |
1129 | switch (d) { | |
1130 | case 1: idx = 0; break; | |
1131 | case 2: idx = (r == 'i' || r == 'I') ? 1 : 2; break; | |
7095962e | 1132 | case 3: idx = (r == 'i' || r == 'I') ? 3 : 4; break; |
0bd4b00f | 1133 | } |
1134 | if (idx < 0 || idx >= fRings.GetEntriesFast()) return 0; | |
1135 | return static_cast<TObjArray*>(fRings.At(idx)); | |
1136 | } | |
1137 | //____________________________________________________________________ | |
1138 | TObjArray* | |
1139 | AliFMDCorrELossFit::GetOrMakeRingArray(UShort_t d, Char_t r) | |
1140 | { | |
7984e5f7 | 1141 | // |
1142 | // Get the ring array corresponding to the specified ring | |
1143 | // | |
1144 | // Parameters: | |
1145 | // d Detector | |
1146 | // r Ring | |
1147 | // | |
1148 | // Return: | |
1149 | // Pointer to ring array, or newly created container | |
1150 | // | |
0bd4b00f | 1151 | Int_t idx = -1; |
1152 | switch (d) { | |
1153 | case 1: idx = 0; break; | |
1154 | case 2: idx = (r == 'i' || r == 'I') ? 1 : 2; break; | |
1155 | case 3: idx = (r == 'o' || r == 'I') ? 3 : 4; break; | |
1156 | } | |
1157 | if (idx < 0) return 0; | |
1158 | if (idx >= fRings.GetEntriesFast() || !fRings.At(idx)) { | |
1159 | TObjArray* a = new TObjArray(0); | |
1160 | // TOrdCollection* a = new TOrdCollection(fEtaAxis.GetNbins()); | |
1161 | a->SetName(Form("FMD%d%c", d, r)); | |
1162 | a->SetOwner(); | |
1163 | fRings.AddAtAndExpand(a, idx); | |
1164 | } | |
1165 | return static_cast<TObjArray*>(fRings.At(idx)); | |
1166 | } | |
1167 | ||
16ac05a3 | 1168 | //____________________________________________________________________ |
1169 | Double_t | |
1170 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Int_t etabin, | |
1171 | Double_t f) const | |
1172 | { | |
8449e3e0 | 1173 | ELossFit* fit = FindFit(d, r, etabin, 20); |
16ac05a3 | 1174 | if (!fit) return -1024; |
1175 | return fit->GetLowerBound(f); | |
1176 | } | |
1177 | //____________________________________________________________________ | |
1178 | Double_t | |
1179 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Double_t eta, | |
1180 | Double_t f) const | |
1181 | { | |
1182 | Int_t bin = FindEtaBin(eta); | |
1183 | if (bin <= 0) return -1024; | |
1184 | return GetLowerBound(d, r, Int_t(bin), f); | |
1185 | } | |
5bb5d1f6 | 1186 | //____________________________________________________________________ |
1187 | Double_t | |
7095962e CHC |
1188 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Int_t etabin, |
1189 | Double_t p, Bool_t) const | |
1190 | { | |
1191 | DGUARD(fDebug, 10, "Get probability cut for FMD%d%c etabin=%d", d, r, etabin); | |
1192 | ELossFit* fit = FindFit(d, r, etabin, 20); | |
1193 | if (!fit) return -1024; | |
1194 | return fit->FindProbabilityCut(p); | |
1195 | } | |
1196 | //____________________________________________________________________ | |
1197 | Double_t | |
1198 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Double_t eta, | |
1199 | Double_t p, Bool_t dummy) const | |
1200 | { | |
1201 | DGUARD(fDebug, 10, "Get probability cut for FMD%d%c eta=%8.4f", d, r, eta); | |
1202 | Int_t bin = FindEtaBin(eta); | |
1203 | DMSG(fDebug, 10, "bin=%4d", bin); | |
1204 | if (bin <= 0) return -1024; | |
1205 | return GetLowerBound(d, r, Int_t(bin), p, dummy); | |
1206 | } | |
1207 | //____________________________________________________________________ | |
1208 | Double_t | |
5bb5d1f6 | 1209 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Int_t etabin, |
1210 | Double_t f, Bool_t showErrors, | |
1211 | Bool_t includeSigma) const | |
1212 | { | |
8449e3e0 | 1213 | ELossFit* fit = FindFit(d, r, etabin, 20); |
5bb5d1f6 | 1214 | if (!fit) { |
1215 | if (showErrors) { | |
1216 | AliWarning(Form("No fit for FMD%d%c @ etabin=%d", d, r, etabin)); | |
1217 | } | |
1218 | return -1024; | |
1219 | } | |
1220 | return fit->GetLowerBound(f, includeSigma); | |
1221 | } | |
1222 | ||
1223 | //____________________________________________________________________ | |
1224 | Double_t | |
1225 | AliFMDCorrELossFit::GetLowerBound(UShort_t d, Char_t r, Double_t eta, | |
1226 | Double_t f, Bool_t showErrors, | |
1227 | Bool_t includeSigma) const | |
1228 | { | |
1229 | Int_t bin = FindEtaBin(eta); | |
1230 | if (bin <= 0) { | |
1231 | if (showErrors) | |
1232 | AliError(Form("eta=%f out of bounds for FMD%d%c", eta, d, r)); | |
1233 | return -1024; | |
1234 | } | |
1235 | return GetLowerBound(d, r, bin, f, showErrors, includeSigma); | |
1236 | } | |
1237 | ||
1238 | //____________________________________________________________________ | |
0bd4b00f | 1239 | namespace { |
1240 | TH1D* MakeHist(const TAxis& axis, const char* name, const char* title, | |
1241 | Int_t color) | |
1242 | { | |
1243 | TH1D* h = new TH1D(Form("%s_%s", name, title), | |
1244 | Form("%s %s", name, title), | |
1245 | axis.GetNbins(), axis.GetXmin(), axis.GetXmax()); | |
1246 | h->SetDirectory(0); | |
1247 | h->SetMarkerStyle(20); | |
1248 | h->SetMarkerColor(color); | |
1249 | h->SetMarkerSize(0.5); | |
1250 | h->SetFillColor(color); | |
1251 | h->SetFillStyle(3001); | |
1252 | h->SetLineColor(color); | |
1253 | return h; | |
1254 | } | |
1255 | } | |
1256 | ||
5bb5d1f6 | 1257 | |
1258 | ||
7dacdf84 | 1259 | #define IDX2RING(I) (i == 0 || i == 1 || i == 3 ? 'I' : 'O') |
1260 | #define IDX2DET(I) (i == 0 ? 1 : (i == 1 || i == 2 ? 2 : 3)) | |
0bd4b00f | 1261 | //____________________________________________________________________ |
f7cfc454 | 1262 | TList* |
8449e3e0 | 1263 | AliFMDCorrELossFit::GetStacks(Bool_t err, |
1264 | Bool_t rel, | |
7b2fa237 | 1265 | Bool_t good, |
8449e3e0 | 1266 | UShort_t maxN) const |
0bd4b00f | 1267 | { |
f7cfc454 | 1268 | // Get a list of THStacks |
1269 | Int_t nRings = fRings.GetEntriesFast(); | |
8449e3e0 | 1270 | // Int_t nPad = 6+maxN-1; // 7 regular params, and maxN-1 weights |
7dacdf84 | 1271 | |
1272 | enum { | |
1273 | kChi2nu = 0, | |
1274 | kC = 1, | |
1275 | kDelta = 2, | |
1276 | kXi = 3, | |
1277 | kSigma = 4, | |
1278 | kN = 5 | |
1279 | }; | |
0bd4b00f | 1280 | |
7dacdf84 | 1281 | THStack* sChi2nu; |
1282 | THStack* sC; | |
1283 | THStack* sDelta; | |
1284 | THStack* sXi; | |
1285 | THStack* sSigma; | |
1286 | // THStack* sigman; | |
0bd4b00f | 1287 | THStack* n; |
f7cfc454 | 1288 | TList* stacks = new TList; |
1289 | stacks->AddAt(sChi2nu= new THStack("chi2", "#chi^{2}/#nu"), kChi2nu); | |
1290 | stacks->AddAt(sC = new THStack("c", "C"), kC); | |
1291 | stacks->AddAt(sDelta = new THStack("delta", "#Delta_{mp}"), kDelta); | |
1292 | stacks->AddAt(sXi = new THStack("xi", "#xi"), kXi); | |
1293 | stacks->AddAt(sSigma = new THStack("sigma", "#sigma"), kSigma); | |
1294 | //stacks->AddAt(sigman= new THStack("sigman", "#sigma_{n}"), 5); | |
8449e3e0 | 1295 | stacks->AddAt(n = new THStack("n", "N"), kN); |
1296 | if (rel) { | |
1297 | sChi2nu->SetName("qual"); | |
1298 | sChi2nu->SetTitle("Quality"); | |
1299 | n->SetName("good"); | |
1300 | n->SetTitle("Bin map"); | |
1301 | } | |
0bd4b00f | 1302 | for (Int_t i = 1; i <= maxN; i++) { |
f7cfc454 | 1303 | stacks->AddAt(new THStack(Form("a_%02d", i+1), Form("a_{%d}", i+1)), kN+i); |
0bd4b00f | 1304 | } |
f7cfc454 | 1305 | |
8449e3e0 | 1306 | // TArrayD min(nPad); |
1307 | // TArrayD max(nPad); | |
1308 | // min.Reset(100000); | |
1309 | // max.Reset(-100000); | |
7dacdf84 | 1310 | |
0bd4b00f | 1311 | for (Int_t i = 0; i < nRings; i++) { |
1312 | if (!fRings.At(i)) continue; | |
1313 | TObjArray* a = static_cast<TObjArray*>(fRings.At(i)); | |
1314 | Int_t nFits = a->GetEntriesFast(); | |
7b2fa237 | 1315 | UShort_t d = IDX2DET(i); |
1316 | Char_t r = IDX2RING(i); | |
1317 | Int_t color = AliForwardUtil::RingColor(d, r); | |
0bd4b00f | 1318 | |
7b2fa237 | 1319 | TH1* hChi = MakeHist(fEtaAxis,a->GetName(), "chi2", color); |
1320 | TH1* hC = MakeHist(fEtaAxis,a->GetName(), "c", color); | |
1321 | TH1* hDelta = MakeHist(fEtaAxis,a->GetName(), "delta", color); | |
1322 | TH1* hXi = MakeHist(fEtaAxis,a->GetName(), "xi", color); | |
1323 | TH1* hSigma = MakeHist(fEtaAxis,a->GetName(), "sigma", color); | |
7dacdf84 | 1324 | // TH1D* hSigmaN = MakeHist(fEtaAxis,a->GetName(), "sigman", color); |
7b2fa237 | 1325 | TH1* hN = MakeHist(fEtaAxis,a->GetName(), "n", color); |
1326 | TH1* hA[maxN]; | |
a26aa8b2 | 1327 | for (Int_t j = 0; j < maxN; j++) hA[j] = 0; |
0bd4b00f | 1328 | const char* ho = (rel || !err ? "hist" : "e"); |
7dacdf84 | 1329 | sChi2nu->Add(hChi, "hist"); // 0 |
1330 | sC ->Add(hC, ho); // 1 | |
1331 | sDelta ->Add(hDelta, ho); // 2 | |
1332 | sXi ->Add(hXi, ho); // 3 | |
1333 | sSigma ->Add(hSigma, ho); // 4 | |
1334 | // sigman->Add(hSigmaN,ho); // 5 | |
1335 | n ->Add(hN, "hist"); // 5 | |
0bd4b00f | 1336 | hChi->SetFillColor(color); |
1337 | hChi->SetFillStyle(3001); | |
1338 | hN->SetFillColor(color); | |
1339 | hN->SetFillStyle(3001); | |
1340 | ||
1341 | for (Int_t k = 1; k <= maxN; k++) { | |
1342 | hA[k-1] = MakeHist(fEtaAxis,a->GetName(), Form("a%02d",k+1), color); | |
f7cfc454 | 1343 | static_cast<THStack*>(stacks->At(kN+k))->Add(hA[k-1], ho); |
0bd4b00f | 1344 | } |
7b2fa237 | 1345 | |
1346 | if (good) { | |
1347 | for (Int_t j = 1; j <= fEtaAxis.GetNbins(); j++) { | |
1348 | ELossFit* f = FindFit(d, r, j, 20); | |
1349 | if (!f) continue; | |
1350 | ||
1351 | UpdateStackHist(f, rel, j, hChi, hN, hC, hDelta, hXi, hSigma, maxN, hA); | |
0bd4b00f | 1352 | } |
7b2fa237 | 1353 | } |
1354 | else { | |
1355 | for (Int_t j = 0; j < nFits; j++) { | |
1356 | ELossFit* f = static_cast<ELossFit*>(a->At(j)); | |
1357 | if (!f) continue; | |
1358 | ||
1359 | UpdateStackHist(f, rel, CACHE(j,i,fEtaAxis.GetNbins()), | |
1360 | hChi, hN, hC, hDelta, hXi, hSigma, maxN, hA); | |
0bd4b00f | 1361 | } |
1362 | } | |
1363 | } | |
f7cfc454 | 1364 | return stacks; |
1365 | } | |
7b2fa237 | 1366 | //____________________________________________________________________ |
1367 | void | |
1368 | AliFMDCorrELossFit::UpdateStackHist(ELossFit* f, Bool_t rel, | |
1369 | Int_t used, | |
1370 | TH1* hChi, TH1* hN, | |
1371 | TH1* hC, TH1* hDelta, | |
1372 | TH1* hXi, TH1* hSigma, | |
1373 | Int_t maxN, TH1** hA) const | |
1374 | { | |
1375 | Int_t b =f->fBin; | |
1376 | Double_t chi2nu =(rel ? f->fQuality : (f->fNu <= 0 ? 0 : f->fChi2 / f->fNu)); | |
1377 | Double_t c =(rel ? (f->fC >0 ?f->fEC /f->fC :0) : f->fC); | |
1378 | Double_t delta =(rel ? (f->fDelta >0 ?f->fEDelta /f->fDelta :0) : f->fDelta); | |
1379 | Double_t xi =(rel ? (f->fXi >0 ?f->fEXi /f->fXi :0) : f->fXi); | |
1380 | Double_t sigma =(rel ? (f->fSigma >0 ?f->fESigma /f->fSigma :0) : f->fSigma); | |
1381 | Int_t w =(rel ? used : f->FindMaxWeight()); | |
1382 | // Double_t sigman = (rel ? (f->fSigmaN>0 ?f->fESigmaN/f->fSigmaN :0) | |
1383 | // : f->SigmaN); | |
1384 | hChi ->SetBinContent(b, chi2nu); | |
1385 | hN ->SetBinContent(b, w); | |
1386 | hC ->SetBinContent(b, c); | |
1387 | hDelta ->SetBinContent(b, delta); | |
1388 | hXi ->SetBinContent(b, xi); | |
1389 | hSigma ->SetBinContent(b, sigma); | |
1390 | ||
1391 | if (!rel) { | |
1392 | hC ->SetBinError(b, f->fEC); | |
1393 | hDelta ->SetBinError(b, f->fEDelta); | |
1394 | hXi ->SetBinError(b, f->fEXi); | |
1395 | hSigma ->SetBinError(b, f->fESigma); | |
1396 | // hSigmaN->SetBinError(b, f->fESigmaN); | |
1397 | } | |
1398 | for (Int_t k = 0; k < f->fN-1 && k < maxN; k++) { | |
1399 | Double_t a = (rel ? (f->fA[k] > 0 ? f->fEA[k] / f->fA[k] : 0) : f->fA[k]); | |
1400 | hA[k]->SetBinContent(b, a); | |
1401 | if (!rel) hA[k]->SetBinError(b, f->fEA[k]); | |
1402 | } | |
1403 | } | |
f7cfc454 | 1404 | |
1405 | //____________________________________________________________________ | |
1406 | void | |
1407 | AliFMDCorrELossFit::Draw(Option_t* option) | |
1408 | { | |
1409 | // | |
1410 | // Draw this object | |
1411 | // | |
1412 | // Parameters: | |
1413 | // option Options. Possible values are | |
1414 | // - err Plot error bars | |
1415 | // | |
1416 | TString opt(Form("nostack %s", option)); | |
1417 | opt.ToLower(); | |
1418 | Bool_t rel = (opt.Contains("relative")); | |
1419 | Bool_t err = (opt.Contains("error")); | |
8449e3e0 | 1420 | Bool_t clr = (opt.Contains("clear")); |
1421 | Bool_t gdd = (opt.Contains("good")); | |
f7cfc454 | 1422 | if (rel) opt.ReplaceAll("relative",""); |
1423 | if (err) opt.ReplaceAll("error",""); | |
8449e3e0 | 1424 | if (clr) opt.ReplaceAll("clear", ""); |
1425 | if (gdd) opt.ReplaceAll("good", ""); | |
f7cfc454 | 1426 | |
1427 | UShort_t maxN = 0; | |
1428 | Int_t nRings = fRings.GetEntriesFast(); | |
1429 | for (Int_t i = 0; i < nRings; i++) { | |
1430 | if (!fRings.At(i)) continue; | |
1431 | TObjArray* a = static_cast<TObjArray*>(fRings.At(i)); | |
1432 | Int_t nFits = a->GetEntriesFast(); | |
1433 | ||
1434 | for (Int_t j = 0; j < nFits; j++) { | |
1435 | ELossFit* fit = static_cast<ELossFit*>(a->At(j)); | |
1436 | if (!fit) continue; | |
1437 | maxN = TMath::Max(maxN, UShort_t(fit->fN)); | |
1438 | } | |
1439 | } | |
1440 | // AliInfo(Form("Maximum N is %d", maxN)); | |
1441 | Int_t nPad = 6+maxN-1; // 7 regular params, and maxN-1 weights | |
1442 | TVirtualPad* pad = gPad; | |
8449e3e0 | 1443 | if (clr) { |
1444 | pad->Clear(); | |
1445 | pad->SetTopMargin(0.02); | |
1446 | pad->SetRightMargin(0.02); | |
1447 | pad->SetBottomMargin(0.15); | |
1448 | pad->SetLeftMargin(0.10); | |
1449 | } | |
f7cfc454 | 1450 | pad->Divide(2, (nPad+1)/2, 0.1, 0, 0); |
1451 | ||
8449e3e0 | 1452 | TList* stacks = GetStacks(err, rel, gdd, maxN); |
f7cfc454 | 1453 | |
0bd4b00f | 1454 | Int_t nPad2 = (nPad+1) / 2; |
1455 | for (Int_t i = 0; i < nPad; i++) { | |
1456 | Int_t iPad = 1 + i/nPad2 + 2 * (i % nPad2); | |
1457 | TVirtualPad* p = pad->cd(iPad); | |
1458 | p->SetLeftMargin(.15); | |
1459 | p->SetFillColor(0); | |
1460 | p->SetFillStyle(0); | |
1461 | p->SetGridx(); | |
1462 | p->SetGridy(); | |
1463 | if (rel && i != 0 && i != 6 && i != 5 && i != 4) p->SetLogy(); | |
1464 | ||
7dacdf84 | 1465 | |
f7cfc454 | 1466 | THStack* stack = static_cast<THStack*>(stacks->At(i)); |
a19faec0 | 1467 | if (!stack->GetHists() || stack->GetHists()->GetEntries() <= 0) { |
1468 | AliWarningF("No histograms in %s", stack->GetName()); | |
1469 | continue; | |
1470 | } | |
f7cfc454 | 1471 | // Double_t powMax = TMath::Log10(max[i]); |
1472 | // Double_t powMin = min[i] <= 0 ? powMax : TMath::Log10(min[i]); | |
1473 | // if (powMax-powMin > 2. && min[i] != 0) p->SetLogy(); | |
1474 | ||
7dacdf84 | 1475 | // stack->SetMinimum(min[i]); |
1476 | // stack->SetMaximum(max[i]); | |
0bd4b00f | 1477 | stack->Draw(opt.Data()); |
1478 | ||
1479 | TString tit(stack->GetTitle()); | |
8449e3e0 | 1480 | if (rel && i != 0 && i != 5) |
0bd4b00f | 1481 | tit = Form("#delta %s/%s", tit.Data(), tit.Data()); |
1482 | TH1* hist = stack->GetHistogram(); | |
1483 | TAxis* yaxis = hist->GetYaxis(); | |
1484 | yaxis->SetTitle(tit.Data()); | |
1485 | yaxis->SetTitleSize(0.15); | |
1486 | yaxis->SetLabelSize(0.08); | |
1487 | yaxis->SetTitleOffset(0.35); | |
7dacdf84 | 1488 | yaxis->SetTitleFont(132); |
1489 | yaxis->SetLabelFont(132); | |
0bd4b00f | 1490 | yaxis->SetNdivisions(5); |
1491 | ||
7dacdf84 | 1492 | |
0bd4b00f | 1493 | TAxis* xaxis = stack->GetHistogram()->GetXaxis(); |
1494 | xaxis->SetTitle("#eta"); | |
1495 | xaxis->SetTitleSize(0.15); | |
1496 | xaxis->SetLabelSize(0.08); | |
1497 | xaxis->SetTitleOffset(0.35); | |
7dacdf84 | 1498 | xaxis->SetTitleFont(132); |
1499 | xaxis->SetLabelFont(132); | |
0bd4b00f | 1500 | xaxis->SetNdivisions(10); |
1501 | ||
7dacdf84 | 1502 | stack->Draw(opt.Data()); |
0bd4b00f | 1503 | } |
1504 | pad->cd(); | |
1505 | } | |
1506 | ||
1507 | //____________________________________________________________________ | |
1508 | void | |
1509 | AliFMDCorrELossFit::Print(Option_t* option) const | |
1510 | { | |
fb3430ac | 1511 | // |
1512 | // Print this object. | |
1513 | // | |
1514 | // Parameters: | |
1515 | // option Options | |
1516 | // - R Print recursive | |
1517 | // | |
1518 | // | |
0bd4b00f | 1519 | TString opt(option); |
1520 | opt.ToUpper(); | |
8449e3e0 | 1521 | Int_t nRings = fRings.GetEntriesFast(); |
1522 | bool recurse = opt.Contains("R"); | |
1523 | bool cache = opt.Contains("C") && fCache.GetSize() > 0; | |
1524 | Int_t nBins = fEtaAxis.GetNbins(); | |
0bd4b00f | 1525 | |
1526 | std::cout << "Low cut in fit range: " << fLowCut << "\n" | |
8449e3e0 | 1527 | << "Eta axis: " << nBins |
0bd4b00f | 1528 | << " bins, range [" << fEtaAxis.GetXmin() << "," |
1529 | << fEtaAxis.GetXmax() << "]" << std::endl; | |
1530 | ||
1531 | for (Int_t i = 0; i < nRings; i++) { | |
1532 | if (!fRings.At(i)) continue; | |
1533 | TObjArray* a = static_cast<TObjArray*>(fRings.At(i)); | |
1534 | Int_t nFits = a->GetEntriesFast(); | |
1535 | ||
1536 | std::cout << a->GetName() << " [" << nFits << " entries]" | |
1537 | << (recurse ? ":\n" : "\t"); | |
1538 | Int_t min = fEtaAxis.GetNbins()+1; | |
1539 | Int_t max = 0; | |
1540 | for (Int_t j = 0; j < nFits; j++) { | |
1541 | if (!a->At(j)) continue; | |
1542 | ||
1543 | min = TMath::Min(j, min); | |
1544 | max = TMath::Max(j, max); | |
1545 | ||
1546 | if (recurse) { | |
1547 | std::cout << "Bin # " << j << "\t"; | |
1548 | ELossFit* fit = static_cast<ELossFit*>(a->At(j)); | |
1549 | fit->Print(option); | |
1550 | } | |
1551 | } | |
1552 | if (!recurse) | |
1553 | std::cout << " bin range: " << std::setw(3) << min | |
1554 | << "-" << std::setw(3) << max << " " << std::setw(3) | |
1555 | << (max-min+1) << " bins" << std::endl; | |
1556 | } | |
8449e3e0 | 1557 | |
1558 | if (!cache) return; | |
1559 | ||
1560 | std::cout << " eta bin | Fit bin \n" | |
1561 | << " # range | FMD1i FMD2i FMD2o FMD3i FMD3o" | |
1562 | // << "----+-----+++------+-----------------------------------" | |
1563 | << std::endl; | |
1564 | size_t oldPrec = std::cout.precision(); | |
1565 | std::cout.precision(3); | |
1566 | for (Int_t i = 1; i <= nBins; i++) { | |
1567 | std::cout << std::setw(4) << i << " " | |
1568 | << std::setw(5) << std::showpos << fEtaAxis.GetBinLowEdge(i) | |
1569 | << " - " << std::setw(5) << fEtaAxis.GetBinUpEdge(i) | |
1570 | << std::noshowpos << " | "; | |
1571 | for (Int_t j = 0; j < 5; j++) { | |
1572 | Int_t bin = CACHE(i,j,nBins); | |
1573 | if (bin <= 0) std::cout << " "; | |
1574 | else std::cout << std::setw(5) << bin | |
1575 | << (bin == i ? ' ' : '*') << ' '; | |
1576 | } | |
1577 | std::cout << std::endl; | |
1578 | } | |
1579 | std::cout.precision(oldPrec); | |
0bd4b00f | 1580 | } |
1581 | ||
1582 | //____________________________________________________________________ | |
1583 | void | |
1584 | AliFMDCorrELossFit::Browse(TBrowser* b) | |
1585 | { | |
7984e5f7 | 1586 | // |
1587 | // Browse this object | |
1588 | // | |
1589 | // Parameters: | |
1590 | // b | |
1591 | // | |
0bd4b00f | 1592 | b->Add(&fRings); |
1593 | b->Add(&fEtaAxis); | |
1594 | } | |
1595 | ||
1596 | ||
1597 | ||
1598 | //____________________________________________________________________ | |
1599 | // | |
1600 | // EOF | |
1601 | // |