]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FORWARD/analysis2/AliFMDDensityCalculator.h
Mega commit.
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliFMDDensityCalculator.h
CommitLineData
7984e5f7 1// This class calculates the inclusive charged particle density
2// in each for the 5 FMD rings.
3//
4#ifndef ALIFMDDENSITYCALCULATOR_H
5#define ALIFMDDENSITYCALCULATOR_H
ffca499d 6/**
7 * @file AliFMDDensityCalculator.h
8 * @author Christian Holm Christensen <cholm@dalsgaard.hehi.nbi.dk>
9 * @date Wed Mar 23 14:02:09 2011
10 *
11 * @brief
12 *
13 *
14 * @ingroup pwg2_forward_aod
15 */
7e4038b5 16#include <TNamed.h>
17#include <TList.h>
1174780f 18#include <TArrayI.h>
7e4038b5 19#include "AliForwardUtil.h"
20class AliESDFMD;
21class TH2D;
dd497217 22class TH1D;
0bd4b00f 23class TProfile;
1174780f 24class AliFMDCorrELossFit;
7e4038b5 25
26/**
27 * This class calculates the inclusive charged particle density
28 * in each for the 5 FMD rings.
29 *
30 * @par Input:
31 * - AliESDFMD object possibly corrected for sharing
32 *
33 * @par Output:
34 * - 5 RingHistos objects - each with a number of vertex dependent
35 * 2D histograms of the inclusive charge particle density
36 *
37 * @par Corrections used:
38 * - AliFMDAnaCalibEnergyDistribution
39 * - AliFMDDoubleHitCorrection
40 * - AliFMDDeadCorrection
41 *
8eb443e1 42 * @ingroup pwg2_forward_algo
ffca499d 43 * @ingroup pwg2_forward_aod
7e4038b5 44 */
45class AliFMDDensityCalculator : public TNamed
46{
47public:
48 /**
49 * Constructor
50 */
51 AliFMDDensityCalculator();
52 /**
53 * Constructor
54 *
55 * @param name Name of object
56 */
57 AliFMDDensityCalculator(const char* name);
58 /**
59 * Copy constructor
60 *
61 * @param o Object to copy from
62 */
63 AliFMDDensityCalculator(const AliFMDDensityCalculator& o);
64 /**
65 * Destructor
66 */
67 virtual ~AliFMDDensityCalculator();
68 /**
69 * Assignement operator
70 *
71 * @param o Object to assign from
72 *
73 * @return Reference to this object
74 */
7c1a1f1d 75 AliFMDDensityCalculator& operator=(const AliFMDDensityCalculator& o);
1174780f 76 /**
77 * Initialize this sub-algorithm
78 *
79 * @param etaAxis Not used
80 */
81 virtual void Init(const TAxis& etaAxis);
7e4038b5 82 /**
83 * Do the calculations
84 *
85 * @param fmd AliESDFMD object (possibly) corrected for sharing
86 * @param hists Histogram cache
87 * @param vtxBin Vertex bin
88 * @param lowFlux Low flux flag.
89 *
90 * @return true on successs
91 */
92 virtual Bool_t Calculate(const AliESDFMD& fmd,
93 AliForwardUtil::Histos& hists,
0bd4b00f 94 UShort_t vtxBin, Bool_t lowFlux);
7e4038b5 95 /**
96 * Scale the histograms to the total number of events
97 *
c389303e 98 * @param dir where to put the output
7e4038b5 99 * @param nEvents Number of events
100 */
fb3430ac 101 virtual void ScaleHistograms(const TList* dir, Int_t nEvents);
7e4038b5 102 /**
103 * Output diagnostic histograms to directory
104 *
105 * @param dir List to write in
106 */
8eb443e1 107 virtual void DefineOutput(TList* dir);
ea3e5d95 108 /**
109 * Set the debug level. The higher the value the more output
110 *
111 * @param dbg Debug level
112 */
113 void SetDebug(Int_t dbg=1) { fDebug = dbg; }
0bd4b00f 114 /**
115 * Maximum particle weight to use
116 *
117 * @param m
118 */
119 void SetMaxParticles(UShort_t m) { fMaxParticles = m; }
9d05ffeb 120 /**
121 * Set whether to use poisson statistics to estimate the
122 * number of particles that has hit within a region. If this is true,
123 * then the average charge particle density is given by
124 * @f[
ffca499d 125 * \lambda = -\log\left(\frac{N_e}{N_t}\right)
9d05ffeb 126 * @f]
127 * where $N_e$ is the number of strips within the region that has no
128 * hits over threshold, and $N_t$ is the total number of strips in the
129 * region/
130 *
131 * @param u Whether to use poisson statistics to estimate the
132 * number of particles that has hit within a region.
133 */
134 void SetUsePoisson(Bool_t u) { fUsePoisson = u; }
5bb5d1f6 135 /**
136 * Set whether to use the phi acceptance correction.
137 *
138 * @param u If true, use the phi acceptance (default is false)
139 */
140 void SetUsePhiAcceptance(Bool_t u) { fUsePhiAcceptance = u; }
0bd4b00f 141 /**
142 * Set the lower multiplicity cut. This overrides the setting in
143 * the energy loss fits.
144 *
145 * @param cut Cut to use
146 */
5bb5d1f6 147 void SetMultCut(Double_t cut) { SetMultCuts(cut,cut,cut,cut,cut); }
148 /**
149 * Set the lower multiplicity cuts
150 *
151 * @param fmd1i Lower mulitplicyt cut for FMD1i
152 * @param fmd2i Lower mulitplicyt cut for FMD2i
153 * @param fmd2o Lower mulitplicyt cut for FMD2o
154 * @param fmd3i Lower mulitplicyt cut for FMD3i
155 * @param fmd3o Lower mulitplicyt cut for FMD3o
156 */
157 void SetMultCuts(Double_t fmd1i,
158 Double_t fmd2i,
159 Double_t fmd2o,
160 Double_t fmd3i,
161 Double_t fmd3o);
b6b35c77 162 /**
163 * Set the luming factors used in the Poisson method
164 *
165 * @param eta Must be 1 or larger
166 * @param phi Must be 1 or larger
167 */
168 void SetLumping(Int_t eta, Int_t phi) {
169 fEtaLumping = (eta < 1 ? 1 : eta);
170 fPhiLumping = (phi < 1 ? 1 : phi);
171 }
5bb5d1f6 172 /**
173 * Set the number of landau width to subtract from the most probably
174 * value to get the low cut.
175 *
176 * @param n Number of @f$ \xi@f$
177 */
178 void SetNXi(Double_t nXi) { fNXi = nXi; }
179 /**
180 * Whether to include sigma in the number subtracted from the MPV to
181 * get the low cut
182 *
183 * @param u If true, then low cut is @f$ \Delta_{mp} - n(\xi+\sigma)@f$
184 */
185 void SetIncludeSigma(Bool_t u) { fIncludeSigma = u; }
186 /**
187 * Get the multiplicity cut. If the user has set fMultCut (via
188 * SetMultCut) then that value is used. If not, then the lower
189 * value of the fit range for the enery loss fits is returned.
190 *
191 * @return Lower cut on multiplicity
192 */
193 Double_t GetMultCut(UShort_t d, Char_t r, Double_t eta,
194 Bool_t errors=true) const;
0bd4b00f 195 /**
196 * Get the multiplicity cut. If the user has set fMultCut (via
197 * SetMultCut) then that value is used. If not, then the lower
198 * value of the fit range for the enery loss fits is returned.
199 *
200 * @return Lower cut on multiplicity
201 */
5bb5d1f6 202 Double_t GetMultCut(UShort_t d, Char_t r, Int_t ieta,
203 Bool_t errors=true) const;
0bd4b00f 204 /**
205 * Print information
206 *
1174780f 207 * @param option Print options
208 * - max Print max weights
0bd4b00f 209 */
210 void Print(Option_t* option="") const;
7e4038b5 211protected:
1174780f 212 /**
213 * Find the max weight to use for FMD<i>dr</i> in eta bin @a iEta
214 *
215 * @param cor Correction
216 * @param d Detector
217 * @param r Ring
218 * @param iEta Eta bin
219 */
fb3430ac 220 Int_t FindMaxWeight(const AliFMDCorrELossFit* cor,
1174780f 221 UShort_t d, Char_t r, Int_t iEta) const;
222
223 /**
224 * Find the max weights and cache them
225 *
226 */
227 void CacheMaxWeights();
228 /**
229 * Find the (cached) maximum weight for FMD<i>dr</i> in
230 * @f$\eta@f$ bin @a iEta
231 *
232 * @param d Detector
233 * @param r Ring
234 * @param iEta Eta bin
235 *
236 * @return max weight or <= 0 in case of problems
237 */
238 Int_t GetMaxWeight(UShort_t d, Char_t r, Int_t iEta) const;
239 /**
240 * Find the (cached) maximum weight for FMD<i>dr</i> iat
241 * @f$\eta@f$
242 *
243 * @param d Detector
244 * @param r Ring
245 * @param eta Eta bin
246 *
247 * @return max weight or <= 0 in case of problems
248 */
249 Int_t GetMaxWeight(UShort_t d, Char_t r, Float_t eta) const;
250
7e4038b5 251 /**
252 * Get the number of particles corresponding to the signal mult
253 *
254 * @param mult Signal
255 * @param d Detector
256 * @param r Ring
257 * @param s Sector
258 * @param t Strip (not used)
259 * @param v Vertex bin
260 * @param eta Pseudo-rapidity
261 * @param lowFlux Low-flux flag
262 *
263 * @return The number of particles
264 */
265 virtual Float_t NParticles(Float_t mult,
266 UShort_t d, Char_t r, UShort_t s, UShort_t t,
0bd4b00f 267 UShort_t v, Float_t eta, Bool_t lowFlux) const;
7e4038b5 268 /**
269 * Get the inverse correction factor. This consist of
270 *
271 * - acceptance correction (corners of sensors)
272 * - double hit correction (for low-flux events)
273 * - dead strip correction
274 *
275 * @param d Detector
276 * @param r Ring
277 * @param s Sector
278 * @param t Strip (not used)
279 * @param v Vertex bin
280 * @param eta Pseudo-rapidity
281 * @param lowFlux Low-flux flag
282 *
283 * @return
284 */
285 virtual Float_t Correction(UShort_t d, Char_t r, UShort_t s, UShort_t t,
0bd4b00f 286 UShort_t v, Float_t eta, Bool_t lowFlux) const;
7e4038b5 287 /**
288 * Get the acceptance correction for strip @a t in an ring of type @a r
289 *
290 * @param r Ring type ('I' or 'O')
291 * @param t Strip number
292 *
293 * @return Inverse acceptance correction
294 */
295 virtual Float_t AcceptanceCorrection(Char_t r, UShort_t t) const;
0bd4b00f 296 /**
297 * Generate the acceptance corrections
298 *
299 * @param r Ring to generate for
300 *
301 * @return Newly allocated histogram of acceptance corrections
302 */
303 virtual TH1D* GenerateAcceptanceCorrection(Char_t r) const;
7e4038b5 304 /**
305 * Internal data structure to keep track of the histograms
306 */
9d99b0dd 307 struct RingHistos : public AliForwardUtil::RingHistos
7e4038b5 308 {
309 /**
310 * Default CTOR
311 */
312 RingHistos();
313 /**
314 * Constructor
315 *
316 * @param d detector
317 * @param r ring
318 */
319 RingHistos(UShort_t d, Char_t r);
320 /**
321 * Copy constructor
322 *
323 * @param o Object to copy from
324 */
325 RingHistos(const RingHistos& o);
326 /**
327 * Assignment operator
328 *
329 * @param o Object to assign from
330 *
331 * @return Reference to this
332 */
333 RingHistos& operator=(const RingHistos& o);
334 /**
335 * Destructor
336 */
337 ~RingHistos();
9d05ffeb 338 /**
339 * Initialize the object
340 *
341 * @param eAxis
342 */
343 void Init(const TAxis& eAxis);
c389303e 344 /**
345 * Make output
346 *
347 * @param dir Where to put it
348 */
7e4038b5 349 void Output(TList* dir);
9d99b0dd 350 /**
351 * Scale the histograms to the total number of events
352 *
c389303e 353 * @param dir Where the output is
9d99b0dd 354 * @param nEvents Number of events
355 */
356 void ScaleHistograms(TList* dir, Int_t nEvents);
e308a636 357 /**
358 * Create Poisson histograms
359 */
b6b35c77 360 void ResetPoissonHistos(const TH2D* h, Int_t etaLumping, Int_t phiLumping);
361 TH2D* fEvsN; // Correlation of Eloss vs uncorrected Nch
362 TH2D* fEvsM; // Correlation of Eloss vs corrected Nch
363 TProfile* fEtaVsN; // Average uncorrected Nch vs eta
364 TProfile* fEtaVsM; // Average corrected Nch vs eta
365 TProfile* fCorr; // Average correction vs eta
366 TH2D* fDensity; // Distribution inclusive Nch
9d05ffeb 367 TH2D* fELossVsPoisson; // Correlation of energy loss vs Poisson N_ch
b6b35c77 368 TH2D* fTotalStrips; // Total number of strips in a region
369 TH2D* fEmptyStrips; // Total number of strips in a region
370 TH2D* fBasicHits ; // Total number basic hits in a region
371 TH2D* fEmptyVsTotal; // # of empty strips vs total number of strips
9fde7142 372
373
1f480471 374 ClassDef(RingHistos,4);
7e4038b5 375 };
376 /**
377 * Get the ring histogram container
378 *
379 * @param d Detector
380 * @param r Ring
381 *
382 * @return Ring histogram container
383 */
384 RingHistos* GetRingHistos(UShort_t d, Char_t r) const;
385
ea3e5d95 386 TList fRingHistos; // List of histogram containers
387 Double_t fMultCut; // Low cut on scaled energy loss
5bb5d1f6 388 Double_t fNXi; // Delta-n(xi+sigma) factor
389 Bool_t fIncludeSigma; // Whether or not to include sigma in cut
1174780f 390 TH1D* fSumOfWeights; // Histogram
391 TH1D* fWeightedSum; // Histogram
392 TH1D* fCorrections; // Histogram
0bd4b00f 393 UShort_t fMaxParticles; // Maximum particle weight to use
9d05ffeb 394 Bool_t fUsePoisson; // If true, then use poisson statistics
5bb5d1f6 395 Bool_t fUsePhiAcceptance; // Whether to correct for corners
0bd4b00f 396 TH1D* fAccI; // Acceptance correction for inner rings
397 TH1D* fAccO; // Acceptance correction for outer rings
1174780f 398 TArrayI fFMD1iMax; // Array of max weights
399 TArrayI fFMD2iMax; // Array of max weights
400 TArrayI fFMD2oMax; // Array of max weights
401 TArrayI fFMD3iMax; // Array of max weights
402 TArrayI fFMD3oMax; // Array of max weights
5bb5d1f6 403 TH2D* fMaxWeights; // Histogram of max weights
404 TH2D* fLowCuts; // Histogram of low cuts
b6b35c77 405 Int_t fEtaLumping; // How to lump eta bins for Poisson
406 Int_t fPhiLumping; // How to lump phi bins for Poisson
ea3e5d95 407 Int_t fDebug; // Debug level
5bb5d1f6 408 Double_t fMultCuts[5]; // Per-ring multiplicity cuts
7e4038b5 409
5bb5d1f6 410 ClassDef(AliFMDDensityCalculator,3); // Calculate Nch density
7e4038b5 411};
412
413#endif
414// Local Variables:
415// mode: C++
416// End:
417