Added some more scripts
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis2 / AliForwardUtil.h
CommitLineData
7e4038b5 1#ifndef ALIROOT_PWG2_FORWARD_ALIFORWARDUTIL_H
2#define ALIROOT_PWG2_FORWARD_ALIFORWARDUTIL_H
3#include <TObject.h>
9d99b0dd 4#include <TString.h>
7f759bb7 5#include <TObjArray.h>
7e4038b5 6class TH2D;
9d99b0dd 7class TH1I;
8class TH1;
7f759bb7 9class TF1;
7e4038b5 10class TAxis;
9d99b0dd 11class AliESDEvent;
7e4038b5 12
13/**
14 * Utilities used in the forward multiplcity analysis
15 *
16 * @ingroup pwg2_forward_analysis
17 */
18class AliForwardUtil : public TObject
19{
9d99b0dd 20public:
7f759bb7 21 //__________________________________________________________________
22 /**
23 * Number of steps to do in the Landau, Gaussiam convolution
24 */
25 static Int_t fgConvolutionSteps;
26 //------------------------------------------------------------------
27 /**
28 * How many sigma's of the Gaussian in the Landau, Gaussian
29 * convolution to integrate over
30 */
31 static Double_t fgConvolutionNSigma;
32 //------------------------------------------------------------------
33 /**
34 * Calculate the shifted Landau
35 * @f[
36 * f'_{L}(x;\Delta,\xi) = f_L(x;\Delta+0.22278298\xi)
37 * @f]
38 *
39 * where @f$ f_{L}@f$ is the ROOT implementation of the Landau
40 * distribution (known to have @f$ \Delta_{p}=-0.22278298@f$ for
41 * @f$\Delta=0,\xi=1@f$.
42 *
43 * @param x Where to evaluate @f$ f'_{L}@f$
44 * @param delta Most probable value
45 * @param xi The 'width' of the distribution
46 *
c389303e 47 * @return @f$ f'_{L}(x;\Delta,\xi) @f$
7f759bb7 48 */
49 static Double_t Landau(Double_t x, Double_t delta, Double_t xi);
50
51 //------------------------------------------------------------------
9d99b0dd 52 /**
7f759bb7 53 * Calculate the value of a Landau convolved with a Gaussian
9d99b0dd 54 *
7f759bb7 55 * @f[
c389303e 56 * f(x;\Delta,\xi,\sigma') = \frac{1}{\sigma' \sqrt{2 \pi}}
7f759bb7 57 * \int_{-\infty}^{+\infty} d\Delta' f'_{L}(x;\Delta',\xi)
c389303e 58 * \exp{-\frac{(\Delta-\Delta')^2}{2\sigma'^2}}
7f759bb7 59 * @f]
9d99b0dd 60 *
c389303e 61 * where @f$ f'_{L}@f$ is the Landau distribution, @f$ \Delta@f$ the
62 * energy loss, @f$ \xi@f$ the width of the Landau, and
63 * @f$ \sigma'^2=\sigma^2-\sigma_n^2 @f$. Here, @f$\sigma@f$ is the
7f759bb7 64 * variance of the Gaussian, and @f$\sigma_n@f$ is a parameter modelling
65 * noise in the detector.
66 *
67 * Note that this function uses the constants fgConvolutionSteps and
68 * fgConvolutionNSigma
69 *
70 * References:
71 * - <a href="http://dx.doi.org/10.1016/0168-583X(84)90472-5">Nucl.Instrum.Meth.B1:16</a>
72 * - <a href="http://dx.doi.org/10.1103/PhysRevA.28.615">Phys.Rev.A28:615</a>
73 * - <a href="http://root.cern.ch/root/htmldoc/tutorials/fit/langaus.C.html">ROOT implementation</a>
74 *
75 * @param x where to evaluate @f$ f@f$
76 * @param delta @f$ \Delta@f$ of @f$ f(x;\Delta,\xi,\sigma')@f$
77 * @param xi @f$ \xi@f$ of @f$ f(x;\Delta,\xi,\sigma')@f$
c389303e 78 * @param sigma @f$ \sigma@f$ of @f$\sigma'^2=\sigma^2-\sigma_n^2 @f$
79 * @param sigma_n @f$ \sigma_n@f$ of @f$\sigma'^2=\sigma^2-\sigma_n^2 @f$
7f759bb7 80 *
81 * @return @f$ f@f$ evaluated at @f$ x@f$.
9d99b0dd 82 */
7f759bb7 83 static Double_t LandauGaus(Double_t x, Double_t delta, Double_t xi,
84 Double_t sigma, Double_t sigma_n);
85
86 //------------------------------------------------------------------
9d99b0dd 87 /**
7f759bb7 88 * Evaluate
c389303e 89 * @f[
90 f_N(x;\Delta,\xi,\sigma') = \sum_{i=1}^N a_i f(x;\Delta_i,\xi_i,\sigma'_i)
91 @f]
9d99b0dd 92 *
7f759bb7 93 * where @f$ f(x;\Delta,\xi,\sigma')@f$ is the convolution of a
94 * Landau with a Gaussian (see LandauGaus). Note that
c389303e 95 * @f$ a_1 = 1@f$, @f$\Delta_i = i(\Delta_1 + \xi\log(i))@f$,
96 * @f$\xi_i=i\xi_1@f$, and @f$\sigma_i'^2 = \sigma_n^2 + i\sigma_1^2@f$.
7f759bb7 97 *
98 * References:
99 * - <a href="http://dx.doi.org/10.1016/0168-583X(84)90472-5">Nucl.Instrum.Meth.B1:16</a>
100 * - <a href="http://dx.doi.org/10.1103/PhysRevA.28.615">Phys.Rev.A28:615</a>
101 * - <a href="http://root.cern.ch/root/htmldoc/tutorials/fit/langaus.C.html">ROOT implementation</a>
9d99b0dd 102 *
7f759bb7 103 * @param x Where to evaluate @f$ f_N@f$
104 * @param delta @f$ \Delta_1@f$
105 * @param xi @f$ \xi_1@f$
106 * @param sigma @f$ \sigma_1@f$
107 * @param sigma_n @f$ \sigma_n@f$
108 * @param n @f$ N@f$ in the sum above.
109 * @param a Array of size @f$ N-1@f$ of the weights @f$ a_i@f$ for
110 * @f$ i > 1@f$
111 *
112 * @return @f$ f_N(x;\Delta,\xi,\sigma')@f$
9d99b0dd 113 */
7f759bb7 114 static Double_t NLandauGaus(Double_t x, Double_t delta, Double_t xi,
115 Double_t sigma, Double_t sigma_n, Int_t n,
116 Double_t* a);
117 //__________________________________________________________________
118 /**
119 * Structure to do fits to the energy loss spectrum
120 *
121 */
122 struct ELossFitter
123 {
c389303e 124 enum {
125 kC = 0,
126 kDelta,
127 kXi,
128 kSigma,
129 kSigmaN,
130 kN,
131 kA
132 };
7f759bb7 133 /**
134 * Constructor
135 *
136 * @param lowCut Lower cut of spectrum - data below this cuts is ignored
137 * @param maxRange Maximum range to fit to
138 * @param minusBins The number of bins below maximum to use
139 */
140 ELossFitter(Double_t lowCut, Double_t maxRange, UShort_t minusBins);
141 virtual ~ELossFitter();
142 /**
143 * Clear internal arrays
144 *
145 */
146 void Clear();
147 /**
148 * Fit a 1-particle signal to the passed energy loss distribution
149 *
150 * Note that this function clears the internal arrays first
151 *
152 * @param dist Data to fit the function to
153 * @param sigman If larger than zero, the initial guess of the
154 * detector induced noise. If zero or less, then this
155 * parameter is ignored in the fit (fixed at 0)
156 *
157 * @return The function fitted to the data
158 */
159 TF1* Fit1Particle(TH1* dist, Double_t sigman=-1);
160 /**
161 * Fit a N-particle signal to the passed energy loss distribution
162 *
163 * If there's no 1-particle fit present, it does that first
164 *
165 * @param dist Data to fit the function to
c389303e 166 * @param n Number of particle signals to fit
7f759bb7 167 * @param sigman If larger than zero, the initial guess of the
168 * detector induced noise. If zero or less, then this
169 * parameter is ignored in the fit (fixed at 0)
170 *
171 * @return The function fitted to the data
172 */
173 TF1* FitNParticle(TH1* dist, UShort_t n, Double_t sigman=-1);
174
175
176 const Double_t fLowCut; // Lower cut on data
177 const Double_t fMaxRange; // Maximum range to fit
178 const UShort_t fMinusBins; // Number of bins from maximum to fit 1st peak
179 TObjArray fFitResults; // Array of fit results
180 TObjArray fFunctions; // Array of functions
181 };
182
183
9d99b0dd 184 //__________________________________________________________________
7e4038b5 185 /**
186 * Structure to hold histograms
187 *
188 * @ingroup pwg2_forward_analysis
189 */
190 struct Histos : public TObject
191 {
192 /**
193 * Constructor
194 *
195 *
196 */
197 Histos() : fFMD1i(0), fFMD2i(0), fFMD2o(0), fFMD3i(0), fFMD3o(0) {}
198 /**
199 * Copy constructor
200 *
201 * @param o Object to copy from
202 */
203 Histos(const Histos& o)
204 : TObject(o),
205 fFMD1i(o.fFMD1i),
206 fFMD2i(o.fFMD2i),
207 fFMD2o(o.fFMD2o),
208 fFMD3i(o.fFMD3i),
209 fFMD3o(o.fFMD3o)
210 {}
211 /**
212 * Assignement operator
213 *
214 * @return Reference to this
215 */
216 Histos& operator=(const Histos&) { return *this;}
217 /**
218 * Destructor
219 */
220 ~Histos();
221 /**
222 * Initialize the object
223 *
224 * @param etaAxis Eta axis to use
225 */
226 void Init(const TAxis& etaAxis);
227 /**
228 * Make a histogram
229 *
230 * @param d Detector
231 * @param r Ring
232 * @param etaAxis Eta axis to use
233 *
234 * @return Newly allocated histogram
235 */
236 TH2D* Make(UShort_t d, Char_t r, const TAxis& etaAxis) const;
237 /**
238 * Clear data
239 *
240 * @param option Not used
241 */
242 void Clear(Option_t* option="");
243 // const TH2D* Get(UShort_t d, Char_t r) const;
244 /**
245 * Get the histogram for a particular detector,ring
246 *
247 * @param d Detector
248 * @param r Ring
249 *
250 * @return Histogram for detector,ring or nul
251 */
252 TH2D* Get(UShort_t d, Char_t r) const;
253 TH2D* fFMD1i; // Histogram for FMD1i
254 TH2D* fFMD2i; // Histogram for FMD2i
255 TH2D* fFMD2o; // Histogram for FMD2o
256 TH2D* fFMD3i; // Histogram for FMD3i
257 TH2D* fFMD3o; // Histogram for FMD3o
9d99b0dd 258
259 ClassDef(Histos,1)
7e4038b5 260 };
261
9d99b0dd 262 //__________________________________________________________________
263 struct RingHistos : public TObject
264 {
265 RingHistos() : fDet(0), fRing('\0'), fName("") {}
266 RingHistos(UShort_t d, Char_t r)
267 : fDet(d), fRing(r), fName(TString::Format("FMD%d%c", d, r))
268 {}
269 RingHistos(const RingHistos& o)
270 : TObject(o), fDet(o.fDet), fRing(o.fRing), fName(o.fName)
271 {}
272 virtual ~RingHistos() {}
273 RingHistos& operator=(const RingHistos& o)
274 {
275 TObject::operator=(o);
276 fDet = o.fDet;
277 fRing = o.fRing;
278 fName = o.fName;
279 return *this;
280 }
281 TList* DefineOutputList(TList* d) const;
282 TList* GetOutputList(TList* d) const;
283 TH1* GetOutputHist(TList* d, const char* name) const;
7f759bb7 284 Color_t Color() const
285 {
286 return ((fDet == 1 ? kRed : (fDet == 2 ? kGreen : kBlue))
287 + ((fRing == 'I' || fRing == 'i') ? 2 : -2));
288 }
289
9d99b0dd 290 UShort_t fDet;
291 Char_t fRing;
292 TString fName;
293
294 ClassDef(RingHistos,1)
295 };
296
7e4038b5 297};
298
299#endif
300// Local Variables:
301// mode: C++
302// End:
303