]>
Commit | Line | Data |
---|---|---|
7984e5f7 | 1 | // |
2 | // Class to fit the energy distribution. | |
3 | // | |
7c1a1f1d | 4 | #ifndef ALIFMDENERGYFITTER_H |
5 | #define ALIFMDENERGYFITTER_H | |
ffca499d | 6 | /** |
7 | * @file AliFMDEnergyFitter.h | |
8 | * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk> | |
9 | * @date Wed Mar 23 14:02:23 2011 | |
10 | * | |
11 | * @brief | |
12 | * | |
13 | * | |
bd6f5206 | 14 | * @ingroup pwglf_forward_eloss |
ffca499d | 15 | */ |
f8715167 | 16 | #include <TNamed.h> |
17 | #include <TH1D.h> | |
18 | #include <TAxis.h> | |
19 | #include <TList.h> | |
20 | #include <TObjArray.h> | |
0bd4b00f | 21 | #include <TClonesArray.h> |
22 | #include "AliFMDCorrELossFit.h" | |
f8715167 | 23 | #include "AliForwardUtil.h" |
24 | class AliESDFMD; | |
25 | class TFitResult; | |
26 | class TF1; | |
27 | class TArrayD; | |
28 | ||
29 | /** | |
30 | * Class to fit the energy distribution. | |
31 | * | |
32 | * @par Input: | |
33 | * - AliESDFMD object - from reconstruction | |
34 | * | |
35 | * @par Output: | |
36 | * - Lists of histogram - one per ring. Each list has a number of | |
37 | * histograms corresponding to the number of eta bins defined. | |
38 | * | |
39 | * @par Corrections used: | |
40 | * - None | |
41 | * | |
42 | * | |
bd6f5206 | 43 | * @ingroup pwglf_forward_algo |
44 | * @ingroup pwglf_forward_eloss | |
f8715167 | 45 | */ |
46 | class AliFMDEnergyFitter : public TNamed | |
47 | { | |
48 | public: | |
c389303e | 49 | enum { |
50 | kC = AliForwardUtil::ELossFitter::kC, | |
51 | kDelta = AliForwardUtil::ELossFitter::kDelta, | |
52 | kXi = AliForwardUtil::ELossFitter::kXi, | |
53 | kSigma = AliForwardUtil::ELossFitter::kSigma, | |
54 | kSigmaN = AliForwardUtil::ELossFitter::kSigmaN, | |
55 | kN = AliForwardUtil::ELossFitter::kN, | |
56 | kA = AliForwardUtil::ELossFitter::kA | |
57 | }; | |
58 | ||
f8715167 | 59 | /** |
60 | * Destructor | |
61 | */ | |
62 | virtual ~AliFMDEnergyFitter(); | |
63 | /** | |
64 | * Default Constructor - do not use | |
65 | */ | |
66 | AliFMDEnergyFitter(); | |
67 | /** | |
68 | * Constructor | |
69 | * | |
70 | * @param title Title of object - not significant | |
71 | */ | |
72 | AliFMDEnergyFitter(const char* title); | |
73 | /** | |
74 | * Copy constructor | |
75 | * | |
76 | * @param o Object to copy from | |
77 | */ | |
78 | AliFMDEnergyFitter(const AliFMDEnergyFitter& o); | |
79 | /** | |
80 | * Assignment operator | |
81 | * | |
82 | * @param o Object to assign from | |
83 | * | |
84 | * @return Reference to this | |
85 | */ | |
86 | AliFMDEnergyFitter& operator=(const AliFMDEnergyFitter& o); | |
87 | ||
4b9857f3 | 88 | /** |
89 | * Initialise the task | |
90 | * | |
91 | * @param etaAxis The eta axis to use. Note, that if the eta axis | |
92 | * has already been set (using SetEtaAxis), then this parameter will be | |
93 | * ignored | |
94 | */ | |
f8715167 | 95 | void Init(const TAxis& etaAxis); |
4b9857f3 | 96 | /** |
97 | * Set the eta axis to use. This will force the code to use this | |
98 | * eta axis definition - irrespective of whatever axis is passed to | |
99 | * the Init member function. Therefore, this member function can be | |
100 | * used to force another eta axis than one found in the correction | |
101 | * objects. | |
102 | * | |
103 | * @param nBins Number of bins | |
104 | * @param etaMin Minimum of the eta axis | |
105 | * @param etaMax Maximum of the eta axis | |
106 | */ | |
107 | void SetEtaAxis(Int_t nBins, Double_t etaMin, Double_t etaMax); | |
108 | /** | |
109 | * Set the eta axis to use. This will force the code to use this | |
110 | * eta axis definition - irrespective of whatever axis is passed to | |
111 | * the Init member function. Therefore, this member function can be | |
112 | * used to force another eta axis than one found in the correction | |
113 | * objects. | |
114 | * | |
115 | * @param etaAxis Eta axis to use | |
116 | */ | |
117 | void SetEtaAxis(const TAxis& etaAxis); | |
5e4d905e | 118 | /** |
119 | * Set the centrality bins. E.g., | |
120 | * @code | |
121 | * UShort_t n = 12; | |
122 | * Double_t bins[] = { 0., 5., 10., 15., 20., 30., | |
123 | * 40., 50., 60., 70., 80., 100. }; | |
124 | * task->GetFitter().SetCentralityBins(n, bins); | |
125 | * @endcode | |
126 | * | |
127 | * @param nBins Size of @a bins | |
128 | * @param bins Bin limits. | |
129 | */ | |
130 | void SetCentralityAxis(UShort_t nBins, Double_t* bins); | |
f8715167 | 131 | /** |
132 | * Set the low cut used for energy | |
133 | * | |
134 | * @param lowCut Low cut | |
135 | */ | |
136 | void SetLowCut(Double_t lowCut=0.3) { fLowCut = lowCut; } | |
4b9857f3 | 137 | /** |
138 | * Set the number of bins to subtract | |
139 | * | |
140 | * @param n | |
141 | */ | |
c389303e | 142 | void SetFitRangeBinWidth(UShort_t n=4) { fFitRangeBinWidth = n; } |
f8715167 | 143 | /** |
144 | * Whether or not to enable fitting of the final merged result. | |
145 | * Note, fitting takes quite a while and one should be careful not to do | |
146 | * this needlessly | |
147 | * | |
148 | * @param doFit Whether to do the fits or not | |
149 | */ | |
c389303e | 150 | void SetDoFits(Bool_t doFit=kTRUE) { fDoFits = doFit; } |
0bd4b00f | 151 | /** |
152 | * Set whether to make the corrections object on the output. Note, | |
153 | * fits should be enable for this to have any effect. | |
154 | * | |
155 | * @param doMake If true (false is default), do make the corrections object. | |
156 | */ | |
157 | void SetDoMakeObject(Bool_t doMake=kTRUE) { fDoMakeObject = doMake; } | |
f8715167 | 158 | /** |
c389303e | 159 | * Set how many particles we will try to fit at most to the data |
f8715167 | 160 | * |
c389303e | 161 | * @param n Max number of particle to try to fit |
f8715167 | 162 | */ |
c389303e | 163 | void SetNParticles(UShort_t n) { fNParticles = (n < 1 ? 1 : (n > 5 ? 5 : n)); } |
f8715167 | 164 | /** |
c389303e | 165 | * Set the minimum number of entries each histogram must have |
166 | * before we try to fit our response function to it | |
f8715167 | 167 | * |
c389303e | 168 | * @param n Minimum number of entries |
f8715167 | 169 | */ |
170 | void SetMinEntries(UShort_t n) { fMinEntries = (n < 1 ? 1 : n); } | |
4b9857f3 | 171 | /** |
172 | * Set maximum energy loss to consider | |
173 | * | |
174 | * @param x Maximum energy loss to consider | |
175 | */ | |
176 | void SetMaxE(Double_t x) { fMaxE = x; } | |
177 | /** | |
178 | * Set number of energy loss bins | |
179 | * | |
180 | * @param x Number of energy loss bins | |
181 | */ | |
182 | void SetNEbins(Int_t x) { fNEbins = x; } | |
7c1a1f1d | 183 | /** |
184 | * Set the maximum relative error | |
185 | * | |
186 | * @param e Maximum relative error | |
187 | */ | |
0bd4b00f | 188 | void SetMaxRelativeParameterError(Double_t e=0.2) { fMaxRelParError = e; } |
7c1a1f1d | 189 | /** |
190 | * Set the maximum @f$ \chi^2/\nu@f$ | |
191 | * | |
192 | * @param c Maximum @f$ \chi^2/\nu@f$ | |
193 | */ | |
0bd4b00f | 194 | void SetMaxChi2PerNDF(Double_t c=10) { fMaxChi2PerNDF = c; } |
7c1a1f1d | 195 | /** |
196 | * Set the least weight | |
197 | * | |
198 | * @param c Least weight | |
199 | */ | |
0bd4b00f | 200 | void SetMinWeight(Double_t c=1e-7) { fMinWeight = c; } |
4b9857f3 | 201 | /** |
202 | * Set wheter to use increasing bin sizes | |
203 | * | |
204 | * @param x Wheter to use increasing bin sizes | |
205 | */ | |
206 | void SetUseIncreasingBins(Bool_t x) { fUseIncreasingBins = x; } | |
f8715167 | 207 | /** |
208 | * Fitter the input AliESDFMD object | |
209 | * | |
210 | * @param input Input | |
5e4d905e | 211 | * @param cent Event centrality (or < 0 if not valid) |
f8715167 | 212 | * @param empty Whether the event is 'empty' |
213 | * | |
214 | * @return True on success, false otherwise | |
215 | */ | |
216 | Bool_t Accumulate(const AliESDFMD& input, | |
5e4d905e | 217 | Double_t cent, |
f8715167 | 218 | Bool_t empty); |
219 | /** | |
220 | * Scale the histograms to the total number of events | |
221 | * | |
c389303e | 222 | * @param dir Where the histograms are |
f8715167 | 223 | */ |
7c1a1f1d | 224 | void Fit(const TList* dir); |
225 | /** | |
226 | * Generate the corrections object | |
227 | * | |
228 | * @param dir List to analyse | |
229 | */ | |
0bd4b00f | 230 | void MakeCorrectionsObject(TList* dir); |
f8715167 | 231 | |
232 | /** | |
233 | * Define the output histograms. These are put in a sub list of the | |
234 | * passed list. The histograms are merged before the parent task calls | |
235 | * AliAnalysisTaskSE::Terminate | |
236 | * | |
237 | * @param dir Directory to add to | |
238 | */ | |
239 | void DefineOutput(TList* dir); | |
240 | /** | |
241 | * Set the debug level. The higher the value the more output | |
242 | * | |
243 | * @param dbg Debug level | |
244 | */ | |
245 | void SetDebug(Int_t dbg=1); | |
0bd4b00f | 246 | /** |
247 | * Print information | |
248 | * | |
249 | * @param option Not used | |
250 | */ | |
251 | void Print(Option_t* option="") const; | |
f8715167 | 252 | protected: |
253 | /** | |
254 | * Internal data structure to keep track of the histograms | |
255 | */ | |
256 | struct RingHistos : public AliForwardUtil::RingHistos | |
257 | { | |
258 | /** | |
259 | * Default CTOR | |
260 | */ | |
261 | RingHistos(); | |
262 | /** | |
263 | * Constructor | |
264 | * | |
265 | * @param d detector | |
266 | * @param r ring | |
267 | */ | |
268 | RingHistos(UShort_t d, Char_t r); | |
269 | /** | |
270 | * Copy constructor | |
271 | * | |
272 | * @param o Object to copy from | |
273 | */ | |
274 | RingHistos(const RingHistos& o); | |
275 | /** | |
276 | * Assignment operator | |
277 | * | |
278 | * @param o Object to assign from | |
279 | * | |
280 | * @return Reference to this | |
281 | */ | |
282 | RingHistos& operator=(const RingHistos& o); | |
283 | /** | |
284 | * Destructor | |
285 | */ | |
286 | ~RingHistos(); | |
287 | /** | |
288 | * Define outputs | |
289 | * | |
290 | * @param dir | |
291 | */ | |
292 | void Output(TList* dir); | |
293 | /** | |
294 | * Initialise object | |
295 | * | |
c389303e | 296 | * @param eAxis Eta axis |
5e4d905e | 297 | * @param cAxis Centrality axis |
c389303e | 298 | * @param maxDE Max energy loss to consider |
299 | * @param nDEbins Number of bins | |
300 | * @param useIncrBin Whether to use an increasing bin size | |
f8715167 | 301 | */ |
66b34429 | 302 | void Init(const TAxis& eAxis, |
5e4d905e | 303 | const TAxis& cAxis, |
4b9857f3 | 304 | Double_t maxDE=10, |
305 | Int_t nDEbins=300, | |
306 | Bool_t useIncrBin=true); | |
f8715167 | 307 | /** |
308 | * Fill histogram | |
309 | * | |
310 | * @param empty True if event is empty | |
5e4d905e | 311 | * @param ieta Eta bin (0 based) |
312 | * @param icent Centrality bin (1 based) | |
f8715167 | 313 | * @param mult Signal |
314 | */ | |
5e4d905e | 315 | void Fill(Bool_t empty, Int_t ieta, Int_t icent, Double_t mult); |
f8715167 | 316 | /** |
c389303e | 317 | * Fit each histogram to up to @a nParticles particle responses. |
f8715167 | 318 | * |
c389303e | 319 | * @param dir Output list |
320 | * @param eta Eta axis | |
321 | * @param lowCut Lower cut | |
322 | * @param nParticles Max number of convolved landaus to fit | |
323 | * @param minEntries Minimum number of entries | |
324 | * @param minusBins Number of bins from peak to subtract to | |
325 | * get the fit range | |
326 | * @param relErrorCut Cut applied to relative error of parameter. | |
327 | * Note, for multi-particle weights, the cut | |
328 | * is loosend by a factor of 2 | |
329 | * @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ - | |
330 | * the reduced @f$\chi^2@f$ | |
290052e7 | 331 | * |
332 | * @return List of fits | |
f8715167 | 333 | */ |
4b9857f3 | 334 | TObjArray* Fit(TList* dir, |
335 | const TAxis& eta, | |
336 | Double_t lowCut, | |
c389303e | 337 | UShort_t nParticles, |
4b9857f3 | 338 | UShort_t minEntries, |
c389303e | 339 | UShort_t minusBins, |
340 | Double_t relErrorCut, | |
341 | Double_t chi2nuCut) const; | |
f8715167 | 342 | /** |
7c1a1f1d | 343 | * Fit a signal histogram. First, the bin @f$ b_{min}@f$ with |
c389303e | 344 | * maximum bin content in the range @f$ [E_{min},\infty]@f$ is |
345 | * found. Then the fit range is set to the bin range | |
346 | * @f$ [b_{min}-\Delta b,b_{min}+2\Delta b]@f$, and a 1 | |
347 | * particle signal is fitted to that. The parameters of that fit | |
348 | * is then used as seeds for a fit of the @f$ N@f$ particle response | |
349 | * to the data in the range | |
350 | * @f$ [b_{min}-\Delta b,N(\Delta_1+\xi_1\log(N))+2N\xi@f$ | |
f8715167 | 351 | * |
c389303e | 352 | * @param dist Histogram to fit |
353 | * @param lowCut Lower cut @f$ E_{min}@f$ on signal | |
354 | * @param nParticles Max number @f$ N@f$ of convolved landaus to fit | |
355 | * @param minusBins Number of bins @f$ \Delta b@f$ from peak to | |
356 | * subtract to get the fit range | |
357 | * @param relErrorCut Cut applied to relative error of parameter. | |
358 | * Note, for multi-particle weights, the cut | |
359 | * is loosend by a factor of 2 | |
360 | * @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ - | |
361 | * the reduced @f$\chi^2@f$ | |
f8715167 | 362 | * |
363 | * @return The best fit function | |
364 | */ | |
4b9857f3 | 365 | TF1* FitHist(TH1* dist, |
366 | Double_t lowCut, | |
c389303e | 367 | UShort_t nParticles, |
368 | UShort_t minusBins, | |
369 | Double_t relErrorCut, | |
370 | Double_t chi2nuCut) const; | |
0bd4b00f | 371 | /** |
372 | * Find the best fits | |
373 | * | |
374 | * @param d Parent list | |
375 | * @param obj Object to add fits to | |
376 | * @param eta Eta axis | |
377 | * @param relErrorCut Cut applied to relative error of parameter. | |
378 | * Note, for multi-particle weights, the cut | |
379 | * is loosend by a factor of 2 | |
380 | * @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ - | |
381 | * the reduced @f$\chi^2@f$ | |
382 | * @param minWeightCut Least valid @f$ a_i@f$ | |
383 | */ | |
fb3430ac | 384 | void FindBestFits(const TList* d, |
0bd4b00f | 385 | AliFMDCorrELossFit& obj, |
386 | const TAxis& eta, | |
387 | Double_t relErrorCut, | |
388 | Double_t chi2nuCut, | |
389 | Double_t minWeightCut); | |
390 | /** | |
391 | * Find the best fit | |
392 | * | |
393 | * @param dist Histogram | |
394 | * @param relErrorCut Cut applied to relative error of parameter. | |
395 | * Note, for multi-particle weights, the cut | |
396 | * is loosend by a factor of 2 | |
397 | * @param chi2nuCut Cut on @f$ \chi^2/\nu@f$ - | |
398 | * the reduced @f$\chi^2@f$ | |
399 | * @param minWeightCut Least valid @f$ a_i@f$ | |
400 | * | |
401 | * @return Best fit | |
402 | */ | |
fb3430ac | 403 | AliFMDCorrELossFit::ELossFit* FindBestFit(const TH1* dist, |
0bd4b00f | 404 | Double_t relErrorCut, |
405 | Double_t chi2nuCut, | |
406 | Double_t minWeightCut); | |
f8715167 | 407 | /** |
4b9857f3 | 408 | * Check the result of the fit. Returns true if |
c389303e | 409 | * - @f$ \chi^2/\nu < \max{\chi^2/\nu}@f$ |
410 | * - @f$ \Delta p_i/p_i < \delta_e@f$ for all parameters. Note, | |
411 | * for multi-particle fits, this requirement is relaxed by a | |
412 | * factor of 2 | |
413 | * - @f$ a_{n} > 10^{-7}@f$ when fitting to an @f$ n@f$ | |
414 | * particle response | |
f8715167 | 415 | * |
c389303e | 416 | * @param r Result to check |
417 | * @param relErrorCut Cut @f$ \delta_e@f$ applied to relative error | |
418 | * of parameter. | |
419 | * @param chi2nuCut Cut @f$ \max{\chi^2/\nu}@f$ | |
f8715167 | 420 | * |
4b9857f3 | 421 | * @return true if fit is good. |
f8715167 | 422 | */ |
c389303e | 423 | Bool_t CheckResult(TFitResult* r, |
424 | Double_t relErrorCut, | |
425 | Double_t chi2nuCut) const; | |
f8715167 | 426 | /** |
4b9857f3 | 427 | * Make an axis with increasing bins |
f8715167 | 428 | * |
4b9857f3 | 429 | * @param n Number of bins |
430 | * @param min Minimum | |
431 | * @param max Maximum | |
f8715167 | 432 | * |
4b9857f3 | 433 | * @return An axis with quadratically increasing bin size |
f8715167 | 434 | */ |
66b34429 | 435 | TArrayD MakeIncreasingAxis(Int_t n, Double_t min, Double_t max) const; |
f8715167 | 436 | /** |
437 | * Make E/E_mip histogram | |
438 | * | |
7c1a1f1d | 439 | * @param ieta Eta bin |
440 | * @param eMin Least signal | |
441 | * @param eMax Largest signal | |
442 | * @param deMax Maximum energy loss | |
443 | * @param nDeBins Number energy loss bins | |
444 | * @param incr Whether to make bins of increasing size | |
f8715167 | 445 | */ |
66b34429 | 446 | void Make(Int_t ieta, Double_t eMin, Double_t eMax, |
447 | Double_t deMax=12, Int_t nDeBins=300, Bool_t incr=true); | |
f8715167 | 448 | /** |
449 | * Make a parameter histogram | |
450 | * | |
451 | * @param name Name of histogram. | |
452 | * @param title Title of histogram. | |
453 | * @param eta Eta axis | |
454 | * | |
455 | * @return | |
456 | */ | |
457 | TH1D* MakePar(const char* name, const char* title, const TAxis& eta) const; | |
4b9857f3 | 458 | /** |
459 | * Make a histogram that contains the results of the fit over the full ring | |
460 | * | |
461 | * @param name Name | |
462 | * @param title Title | |
463 | * @param eta Eta axis | |
464 | * @param low Least bin | |
465 | * @param high Largest bin | |
466 | * @param val Value of parameter | |
467 | * @param err Error on parameter | |
468 | * | |
469 | * @return The newly allocated histogram | |
470 | */ | |
f8715167 | 471 | TH1D* MakeTotal(const char* name, |
472 | const char* title, | |
473 | const TAxis& eta, | |
474 | Int_t low, | |
475 | Int_t high, | |
476 | Double_t val, | |
477 | Double_t err) const; | |
0bd4b00f | 478 | TH1D* fEDist; // Ring energy distribution |
479 | TH1D* fEmpty; // Ring energy distribution for empty events | |
0526c533 | 480 | TList* fEtaEDists; // Energy distributions per eta bin. |
0bd4b00f | 481 | TList* fList; |
482 | TClonesArray fFits; | |
483 | Int_t fDebug; | |
0526c533 | 484 | ClassDef(RingHistos,2); |
f8715167 | 485 | }; |
486 | /** | |
487 | * Get the ring histogram container | |
488 | * | |
489 | * @param d Detector | |
490 | * @param r Ring | |
491 | * | |
492 | * @return Ring histogram container | |
493 | */ | |
494 | RingHistos* GetRingHistos(UShort_t d, Char_t r) const; | |
495 | ||
496 | TList fRingHistos; // List of histogram containers | |
497 | Double_t fLowCut; // Low cut on energy | |
c389303e | 498 | UShort_t fNParticles; // Number of landaus to try to fit |
f8715167 | 499 | UShort_t fMinEntries; // Minimum number of entries |
c389303e | 500 | UShort_t fFitRangeBinWidth;// Number of bins to subtract from found max |
0bd4b00f | 501 | Bool_t fDoFits; // Whether to actually do the fits |
502 | Bool_t fDoMakeObject; // Whether to make corrections object | |
f8715167 | 503 | TAxis fEtaAxis; // Eta axis |
5e4d905e | 504 | TAxis fCentralityAxis;// Centrality axis |
4b9857f3 | 505 | Double_t fMaxE; // Maximum energy loss to consider |
506 | Int_t fNEbins; // Number of energy loss bins | |
507 | Bool_t fUseIncreasingBins; // Wheter to use increasing bin sizes | |
c389303e | 508 | Double_t fMaxRelParError;// Relative error cut |
509 | Double_t fMaxChi2PerNDF; // chi^2/nu cit | |
0bd4b00f | 510 | Double_t fMinWeight; // Minimum weight value |
f8715167 | 511 | Int_t fDebug; // Debug level |
0bd4b00f | 512 | |
f8715167 | 513 | |
0526c533 | 514 | ClassDef(AliFMDEnergyFitter,3); // |
f8715167 | 515 | }; |
516 | ||
517 | #endif | |
518 | // Local Variables: | |
519 | // mode: C++ | |
520 | // End: |