eeb78423cfd928f2b63c91ee0f9eb95c1ef59de8
[u/mrichter/AliRoot.git] / FMD / flow / AliFMDFlowBin.cxx
1 /** @file 
2     @brief implementation of a Bin in a Flow histogram */
3 #include "flow/AliFMDFlowBin.h"
4 #include <cmath>
5 #include <iostream>
6 #include <iomanip>
7
8 //====================================================================
9 void 
10 AliFMDFlowBin::Begin() 
11 {
12   // Clear event plane calculators 
13   fPsi.Clear();
14   fPsiA.Clear();
15   fPsiB.Clear();
16 }
17
18 //____________________________________________________________________
19 void 
20 AliFMDFlowBin::AddToEventPlane(Double_t phi, Double_t w, Bool_t a) 
21 {
22   fPsi.Add(phi, w);
23   if (a) fPsiA.Add(phi, w);
24   else   fPsiB.Add(phi, w);
25 }
26
27 //____________________________________________________________________
28 void 
29 AliFMDFlowBin::AddToHarmonic(Double_t phi, Double_t w)
30 {
31   // Disregard the obervation of phi from the event plane angle. 
32   Double_t psi   = fPsi.Psi(phi, w);
33   fHarmonic.Add(phi, psi);
34 }
35
36 //____________________________________________________________________
37 void 
38 AliFMDFlowBin::End()
39 {
40   Double_t psi_A = fPsiA.Psi();
41   Double_t psi_B = fPsiB.Psi();
42
43   // Update the resolutions 
44   fRes.Add(psi_A, psi_B);
45   fResStar.Add(psi_A, psi_B);
46   fResTdr.Add(psi_A, psi_B);
47 }
48
49 //____________________________________________________________________
50 void 
51 AliFMDFlowBin::Event(Double_t* phis, Double_t* ws, UInt_t n) 
52
53   Begin();
54   
55   // Calculate split. 
56   UInt_t split = n / 2;
57   // First sub-event. 
58   for (UInt_t i = 0; i < split; i++) 
59     AddToEventPlane(phis[i], (ws ? ws[i] : 1), kTRUE);
60   // Second sub-event. 
61   for (UInt_t i = split; i < n; i++) 
62     AddToEventPlane(phis[i], (ws ? ws[i] : 1), kFALSE);
63   // Add contributions to the harmonic. 
64   for (UInt_t i = 0; i < n; i++)     
65     AddToHarmonic(phis[i], (ws ? ws[i] : 1));
66
67   End();
68 }
69
70 //____________________________________________________________________
71 Double_t 
72 AliFMDFlowBin::Value(CorType t) const
73
74   Double_t e;
75   return Value(e, t);
76 }
77
78 //____________________________________________________________________
79 Double_t 
80 AliFMDFlowBin::EValue(CorType t) const 
81
82   Double_t e2;
83   Value(e2, t);
84   return sqrt(e2);
85 }
86
87 //____________________________________________________________________
88 Double_t 
89 AliFMDFlowBin::Value(Double_t& e2, CorType t) const
90
91   Double_t r, er2;
92   r = Correction(er2, t);
93   return fHarmonic.Value(r, er2, e2);
94 }
95
96 //____________________________________________________________________
97 Double_t 
98 AliFMDFlowBin::Correction(Double_t& er2, CorType t) const
99 {
100   Double_t r = 1;
101   UShort_t k = fHarmonic.Order()/fRes.Order();
102   switch (t) { 
103   case naive: r = fRes.Correction(k, er2);     break;
104   case star:  r = fResStar.Correction(k, er2); break;
105   case tdr:   r = fResTdr.Correction(k, er2);  break;
106   default:    r = 1; er2 = 0;                  break;
107   }
108   return r;
109 }
110
111 //____________________________________________________________________
112 void 
113 AliFMDFlowBin::Finish() 
114 {}
115
116 //____________________________________________________________________
117 void 
118 AliFMDFlowBin::Print(Option_t*) const
119 {
120   Double_t e2v[4], v[4], r[4], e2r[4];
121   const char* names[] = { "Bare", "Naive", "STAR", "TDR" };
122   v[0] = 100 * Value(e2v[0], AliFMDFlowBin::none);
123   v[1] = 100 * Value(e2v[1], AliFMDFlowBin::naive);
124   v[2] = 100 * Value(e2v[2], AliFMDFlowBin::star);
125   v[3] = 100 * Value(e2v[3], AliFMDFlowBin::tdr);
126   r[0] = 100 * Correction(e2r[0], AliFMDFlowBin::none);
127   r[1] = 100 * Correction(e2r[1], AliFMDFlowBin::naive);
128   r[2] = 100 * Correction(e2r[2], AliFMDFlowBin::star);
129   r[3] = 100 * Correction(e2r[3], AliFMDFlowBin::tdr);
130   
131   std::streamsize         old_prec  = std::cout.precision(3);
132   std::ios_base::fmtflags old_flags = std::cout.setf(std::ios_base::fixed, 
133                                                      std::ios_base::floatfield);
134   std::cout << "  v" << std::setw(1) << fHarmonic.Order() << ": ";
135   for (UInt_t i = 0; i < 4; i++) 
136     std::cout << std::setw(6+(i == 0 ? 0 : 6)) << names[i] << ": " 
137               << std::setw(6) << v[i] << " +/- " 
138               << std::setw(6) << 100*sqrt(e2v[i]) << " ["
139               << std::setw(7) << r[i] << " +/- " 
140               << std::setw(7) << 100*sqrt(e2r[i]) << "]\n";
141   std::cout << std::flush;
142   std::cout.precision(old_prec);
143   std::cout.setf(old_flags, std::ios_base::floatfield);
144 }
145
146
147 //____________________________________________________________________
148 //
149 // EOF
150 //