e3954ffc23b96677a5bc6667582e12412739bcad
[u/mrichter/AliRoot.git] / FMD / flow / AliFMDFlowResolution.h
1 // -*- mode: C++ -*-
2 /** @file 
3     @brief Declaration of an Resolution class */
4 #ifndef FLOW_RESOLUTION_H
5 #define FLOW_RESOLUTION_H
6 #include <flow/AliFMDFlowStat.h>
7
8 //______________________________________________________
9 /** @class AliFMDFlowResolution flow/AliFMDFlowResolution.h <flow/AliFMDFlowResolution.h>
10     @brief Class to calculate the event plane resolution based on two
11     sub-events 
12     @ingroup a_basic
13     
14     This class calculates the event plane resolution based on the
15     basic formulas given in Phys. Rev. @b C58, 1671.   That is, the
16     resolution is given by 
17     @f[ 
18     R_{k} = \langle\cos(km(\Psi_m-\Psi_R))\rangle 
19     @f] 
20     where @f$ \Psi_R@f$ is the unknown @e true event plane angle, and
21     @f$ n = km@f$
22
23     Using two random sub-events, @f$ A, B@f$ we get that 
24     @f[
25     \langle\cos(n(\Psi^A_m-\Psi^B_m))\rangle = 
26     \langle\cos(n(\Psi^A_m-\Psi_R))\rangle 
27     \langle\cos(n(\Psi^B_m-\Psi_R))\rangle 
28     @f]
29     If the sub-events are of equal size, and randomly chosen, then we
30     get that 
31     @f[ 
32     \langle\cos(n(\Psi^A_m-\Psi_R))\rangle = 
33     \sqrt{\langle\cos(n(\Psi^A_m-\Psi^B_m))}
34     @f]
35     and it follows that 
36     @f{eqnarray*}
37     \langle\cos(km(\Psi_m-\Psi_R))\rangle & = & 
38     \sqrt{2}\langle\cos(n(\Psi^A_m-\Psi_R))\rangle\\
39     \sqrt{2}\sqrt{\langle\cos(n(\Psi^A_m-\Psi^B_m))}
40     @f}
41
42     Hence, the event-plane resolution is simply the square root of the
43     scaled average distance between the two sub-events, multiplied by
44     @f$ \sqrt{s}@f$ 
45
46     The error is therefor @f$ \sqrt{s}@f$ times the variance of the
47     cosine of the distance between the two sub-events. 
48 */
49 class AliFMDFlowResolution : public AliFMDFlowStat
50 {
51 public:
52   /** Constructor
53       @param n Harmonic order */
54   AliFMDFlowResolution(UShort_t n) : fOrder(n) {}
55   /** Destructor */
56   virtual ~AliFMDFlowResolution() {}
57   /** add data point 
58       @param psiA A sub-event plane angle @f$ \Psi_A \in[0,2\pi]@f$
59       @param psiB B sub-event plane angle @f$ \Psi_B \in[0,2\pi]@f$ */
60   virtual void Add(Double_t psiA, Double_t psiB);
61   /** Get the correction for harmonic strength of order @a k 
62       @param k  The harminic strenght order to get the correction for
63       @param e2 The square error on the correction 
64       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
65   virtual Double_t Correction(UShort_t k, Double_t& e2) const;
66   /** Get the correction for harmonic strength of order @a k 
67       @param k The harminic strenght order to get the correction for
68       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
69   virtual Double_t Correction(UShort_t k) const;
70   /** Get the harmnic order */
71   UShort_t Order() const { return fOrder; }
72 protected:
73   /** Order */
74   UShort_t fOrder;
75   /** Define for ROOT I/O */
76   ClassDef(AliFMDFlowResolution,1);
77 };
78
79 //______________________________________________________
80 /** 
81     @ingroup a_basic
82     The event plane angle resolution function is 
83     @f[
84     R_k(\chi) = \frac{\pi} \chi e^{-\chi^2/4}\left( 
85     I_{\frac{k-1}2}(\chi^2/4) + I_{\frac{k+1}2}(\chi^2/4)\right)
86     @f]
87     Where @f$ I_n(x)@f$ is the modified Bessel function of the first
88     kind.   Identifying 
89     @f[
90     y = \chi^2/4\quad C=\frac{\sqrt{\pi}e^{-y}}{2\sqrt{2}}
91     @f]
92     and the short hands @f$ f2(y) = I_{\frac{k-1}2}@f$ and @f$ f3(y) =
93     I_{\frac{k+1}2}(y)@f$ we can write this more compact as 
94     @f[ 
95     R_k(y) = C y (f2(y) + f3(y))
96     @f] 
97     The derivative of the resolution function is 
98     @f[ 
99     R_k'(y) = \frac{C}{2}\left[4\sqrt{y}\left(f2'(y)-f3'(y)\right) 
100     - (4 y - 2)\left(f2(y) + f3(y)\right)\right]
101     @f]
102     Since 
103     @f[
104     I_\nu'(x) = I_{\nu-1}(x) - \frac{\nu}{x} I_\nu(x)\quad,
105     @f]
106     and setting @f$ f1(y) = I_{\frac{k-3}2}(y)@f$, we get 
107     @f[ 
108     R_k'(y) = \frac{C}{2}\left[4y f1(y) + (4-2k) f2(y) - (4y+2k)
109     f3(y)\right]
110     @f]
111
112     In this class, the argument @f$ \chi@f$ is estimated by finding
113     the minima of @f$ R_k(\chi)@f$ near the average of
114     @f$\cos(n(\Psi_A-\Psi_B))@f$.  The error @f$ \delta\chi@f$ is
115     estimated as the largest step size in the minimisation. 
116
117     The total error on the correction is then 
118     @f[ 
119     \delta R_k = R_k'(\chi) \delta\chi
120     @f] 
121 */
122
123 class AliFMDFlowResolutionStar : public AliFMDFlowResolution
124 {
125 public:
126   /** Constructor
127       @param n Harmonic order */
128   AliFMDFlowResolutionStar(UShort_t n) 
129     : AliFMDFlowResolution(n) {}
130   /** Destructor */
131   ~AliFMDFlowResolutionStar() {}
132   /** Get the correction for harmonic strength of order @a k 
133       @param k The harminic strenght order to get the correction for
134       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
135   virtual Double_t Correction(UShort_t k) const;
136
137   /** Get the correction for harmonic strength of order @a k 
138       @param k  The harminic strenght order to get the correction for
139       @param e2 The square error on the correction 
140       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
141   virtual Double_t Correction(UShort_t k, Double_t& e2) const;
142   /** Get @f$ \chi@f$ 
143       @param res  First shot at the resolution. 
144       @param k    Order 
145       @param delta On return, the last step size in @f$ \chi@f$ -
146       which is taken to be @f$ \delta\chi@f$  
147       @return @f$\chi@f$ */
148   virtual Double_t Chi(Double_t res, UShort_t k, Double_t& delta) const;
149 protected:
150   /** Calculate resolution 
151       @param chi @f$ \chi@f$
152       @param k   Order factor 
153       @param dr  On return, the derivative of @f$ R(\chi)@f$
154       @return 
155       @f[ 
156       \frac{\sqrt{\pi/2}}{2}\chi e^{-\chi^2/4}
157       (I_{\frac{(k-1)}{2}}(\chi^2/4)+ I_{\frac{(k+1)}{2}}(\chi^2/4))
158       @f] */
159   Double_t Res(Double_t chi, UShort_t k, Double_t& dr) const;
160   /** Define for ROOT I/O */
161   ClassDef(AliFMDFlowResolutionStar,1);
162 };
163
164 //______________________________________________________
165 /** 
166     @ingroup a_basic
167
168     For more on the event plane angle resolution function, please
169     refer to the class description of ResolutionStar. 
170
171     In this class @f$ \chi@f$ is calculated from the ratio @f$ k/N@f$
172     of events with @f$ |\Psi_A - \Psi_B| > \pi/2@f$ to the total
173     number of events. 
174
175     The pre-print @c nucl-ex/9711003v2 gives the formula 
176     @f[ 
177     \frac{k}{N} = \frac{e^{-\chi^2/2}}{2}
178     @f] 
179     for @f$ \chi@f$.  Note, that this differs from the @f$ \chi@f$
180     used in ResolutionStar by a factor of @f$ 1/sqrt{2}@f$. 
181     
182     We can isolate @f$ \chi = \mp\sqrt{-2\log(2k/N)}@f$ from the
183     equation above.  
184
185     Since @f$ r=k/N@f$ is obviously a efficiency-like ratio, we get
186     that error @f$ \delta r@f$ is given by Binomial errors 
187     @f[ 
188     \delta r =  \sqrt{r\frac{1 - r}{N}}\quad.
189     @f] 
190     The total error @f$ \delta\chi@f$ then becomes 
191     @f[ 
192     \delta^2\chi = \left(\frac{d\chi}{dr}\right)^2 \delta^2r 
193     = \frac{r - 1}{4 k log(2 r)}
194     @f]
195
196     The total error on the correction is, as in
197     ResolutionStar, then given by 
198     @f[ 
199     \delta R_k = R_k'(\chi) \delta\chi
200     @f]
201 */
202 class AliFMDFlowResolutionTDR : public AliFMDFlowResolution
203 {
204 public:
205   /** Constructor
206       @param n Harmonic order */
207   AliFMDFlowResolutionTDR(UShort_t n) 
208     : AliFMDFlowResolution(n), fLarge(0) {}
209   /** DTOR */
210   ~AliFMDFlowResolutionTDR() {}
211   virtual void Clear();
212   /** add a data  point */
213   virtual void Add(Double_t psi_a, Double_t psi_b);
214   /** Get the correction for harmonic strength of order @a k 
215       @param k The harminic strenght order to get the correction for
216       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
217   virtual Double_t Correction(UShort_t k) const;
218   /** Get the correction for harmonic strength of order @a k 
219       @param e2 The square error on the correction 
220       @param k  The harminic strenght order to get the correction for
221       @return @f$ \langle\cos(n(\psi_n - \psi_R))\rangle@f$ */ 
222   virtual Double_t Correction(UShort_t k, Double_t& e2) const;
223   /** Get @f$ \chi^2/2@f$ 
224       @param e2 The square error on the correction 
225       @return @f$ \chi^2/2@f$ */
226   virtual Double_t Chi2Over2(Double_t& e2) const;
227 protected:
228   /** Number of events with large diviation */
229   ULong_t fLarge;
230   /** Define for ROOT I/O */
231   ClassDef(AliFMDFlowResolutionTDR,1);
232 };
233
234
235 #endif
236 //
237 // EOF
238 //