Fixes
[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; }
0bd4b00f 135 /**
136 * Set the lower multiplicity cut. This overrides the setting in
137 * the energy loss fits.
138 *
139 * @param cut Cut to use
140 */
141 void SetMultCut(Double_t cut) { fMultCut = cut; }
b6b35c77 142 /**
143 * Set the luming factors used in the Poisson method
144 *
145 * @param eta Must be 1 or larger
146 * @param phi Must be 1 or larger
147 */
148 void SetLumping(Int_t eta, Int_t phi) {
149 fEtaLumping = (eta < 1 ? 1 : eta);
150 fPhiLumping = (phi < 1 ? 1 : phi);
151 }
0bd4b00f 152 /**
153 * Get the multiplicity cut. If the user has set fMultCut (via
154 * SetMultCut) then that value is used. If not, then the lower
155 * value of the fit range for the enery loss fits is returned.
156 *
157 * @return Lower cut on multiplicity
158 */
159 Double_t GetMultCut() const;
160 /**
161 * Print information
162 *
1174780f 163 * @param option Print options
164 * - max Print max weights
0bd4b00f 165 */
166 void Print(Option_t* option="") const;
7e4038b5 167protected:
1174780f 168 /**
169 * Find the max weight to use for FMD<i>dr</i> in eta bin @a iEta
170 *
171 * @param cor Correction
172 * @param d Detector
173 * @param r Ring
174 * @param iEta Eta bin
175 */
fb3430ac 176 Int_t FindMaxWeight(const AliFMDCorrELossFit* cor,
1174780f 177 UShort_t d, Char_t r, Int_t iEta) const;
178
179 /**
180 * Find the max weights and cache them
181 *
182 */
183 void CacheMaxWeights();
184 /**
185 * Find the (cached) maximum weight for FMD<i>dr</i> in
186 * @f$\eta@f$ bin @a iEta
187 *
188 * @param d Detector
189 * @param r Ring
190 * @param iEta Eta bin
191 *
192 * @return max weight or <= 0 in case of problems
193 */
194 Int_t GetMaxWeight(UShort_t d, Char_t r, Int_t iEta) const;
195 /**
196 * Find the (cached) maximum weight for FMD<i>dr</i> iat
197 * @f$\eta@f$
198 *
199 * @param d Detector
200 * @param r Ring
201 * @param eta Eta bin
202 *
203 * @return max weight or <= 0 in case of problems
204 */
205 Int_t GetMaxWeight(UShort_t d, Char_t r, Float_t eta) const;
206
7e4038b5 207 /**
208 * Get the number of particles corresponding to the signal mult
209 *
210 * @param mult Signal
211 * @param d Detector
212 * @param r Ring
213 * @param s Sector
214 * @param t Strip (not used)
215 * @param v Vertex bin
216 * @param eta Pseudo-rapidity
217 * @param lowFlux Low-flux flag
218 *
219 * @return The number of particles
220 */
221 virtual Float_t NParticles(Float_t mult,
222 UShort_t d, Char_t r, UShort_t s, UShort_t t,
0bd4b00f 223 UShort_t v, Float_t eta, Bool_t lowFlux) const;
7e4038b5 224 /**
225 * Get the inverse correction factor. This consist of
226 *
227 * - acceptance correction (corners of sensors)
228 * - double hit correction (for low-flux events)
229 * - dead strip correction
230 *
231 * @param d Detector
232 * @param r Ring
233 * @param s Sector
234 * @param t Strip (not used)
235 * @param v Vertex bin
236 * @param eta Pseudo-rapidity
237 * @param lowFlux Low-flux flag
238 *
239 * @return
240 */
241 virtual Float_t Correction(UShort_t d, Char_t r, UShort_t s, UShort_t t,
0bd4b00f 242 UShort_t v, Float_t eta, Bool_t lowFlux) const;
7e4038b5 243 /**
244 * Get the acceptance correction for strip @a t in an ring of type @a r
245 *
246 * @param r Ring type ('I' or 'O')
247 * @param t Strip number
248 *
249 * @return Inverse acceptance correction
250 */
251 virtual Float_t AcceptanceCorrection(Char_t r, UShort_t t) const;
0bd4b00f 252 /**
253 * Generate the acceptance corrections
254 *
255 * @param r Ring to generate for
256 *
257 * @return Newly allocated histogram of acceptance corrections
258 */
259 virtual TH1D* GenerateAcceptanceCorrection(Char_t r) const;
7e4038b5 260 /**
261 * Internal data structure to keep track of the histograms
262 */
9d99b0dd 263 struct RingHistos : public AliForwardUtil::RingHistos
7e4038b5 264 {
265 /**
266 * Default CTOR
267 */
268 RingHistos();
269 /**
270 * Constructor
271 *
272 * @param d detector
273 * @param r ring
274 */
275 RingHistos(UShort_t d, Char_t r);
276 /**
277 * Copy constructor
278 *
279 * @param o Object to copy from
280 */
281 RingHistos(const RingHistos& o);
282 /**
283 * Assignment operator
284 *
285 * @param o Object to assign from
286 *
287 * @return Reference to this
288 */
289 RingHistos& operator=(const RingHistos& o);
290 /**
291 * Destructor
292 */
293 ~RingHistos();
9d05ffeb 294 /**
295 * Initialize the object
296 *
297 * @param eAxis
298 */
299 void Init(const TAxis& eAxis);
c389303e 300 /**
301 * Make output
302 *
303 * @param dir Where to put it
304 */
7e4038b5 305 void Output(TList* dir);
9d99b0dd 306 /**
307 * Scale the histograms to the total number of events
308 *
c389303e 309 * @param dir Where the output is
9d99b0dd 310 * @param nEvents Number of events
311 */
312 void ScaleHistograms(TList* dir, Int_t nEvents);
e308a636 313 /**
314 * Create Poisson histograms
315 */
b6b35c77 316 void ResetPoissonHistos(const TH2D* h, Int_t etaLumping, Int_t phiLumping);
317 TH2D* fEvsN; // Correlation of Eloss vs uncorrected Nch
318 TH2D* fEvsM; // Correlation of Eloss vs corrected Nch
319 TProfile* fEtaVsN; // Average uncorrected Nch vs eta
320 TProfile* fEtaVsM; // Average corrected Nch vs eta
321 TProfile* fCorr; // Average correction vs eta
322 TH2D* fDensity; // Distribution inclusive Nch
9d05ffeb 323 TH2D* fELossVsPoisson; // Correlation of energy loss vs Poisson N_ch
b6b35c77 324 TH2D* fTotalStrips; // Total number of strips in a region
325 TH2D* fEmptyStrips; // Total number of strips in a region
326 TH2D* fBasicHits ; // Total number basic hits in a region
327 TH2D* fEmptyVsTotal; // # of empty strips vs total number of strips
9fde7142 328
329
e308a636 330 ClassDef(RingHistos,3);
7e4038b5 331 };
332 /**
333 * Get the ring histogram container
334 *
335 * @param d Detector
336 * @param r Ring
337 *
338 * @return Ring histogram container
339 */
340 RingHistos* GetRingHistos(UShort_t d, Char_t r) const;
341
ea3e5d95 342 TList fRingHistos; // List of histogram containers
343 Double_t fMultCut; // Low cut on scaled energy loss
1174780f 344 TH1D* fSumOfWeights; // Histogram
345 TH1D* fWeightedSum; // Histogram
346 TH1D* fCorrections; // Histogram
0bd4b00f 347 UShort_t fMaxParticles; // Maximum particle weight to use
9d05ffeb 348 Bool_t fUsePoisson; // If true, then use poisson statistics
0bd4b00f 349 TH1D* fAccI; // Acceptance correction for inner rings
350 TH1D* fAccO; // Acceptance correction for outer rings
1174780f 351 TArrayI fFMD1iMax; // Array of max weights
352 TArrayI fFMD2iMax; // Array of max weights
353 TArrayI fFMD2oMax; // Array of max weights
354 TArrayI fFMD3iMax; // Array of max weights
355 TArrayI fFMD3oMax; // Array of max weights
b6b35c77 356 Int_t fEtaLumping; // How to lump eta bins for Poisson
357 Int_t fPhiLumping; // How to lump phi bins for Poisson
ea3e5d95 358 Int_t fDebug; // Debug level
7e4038b5 359
b6b35c77 360
9d05ffeb 361 ClassDef(AliFMDDensityCalculator,2); // Calculate Nch density
7e4038b5 362};
363
364#endif
365// Local Variables:
366// mode: C++
367// End:
368