Coding rule violations corrected.
[u/mrichter/AliRoot.git] / ANALYSIS / AliFlowAnalysis.cxx
1 #include "AliFlowAnalysis.h"
2 //________________________________
3 ///////////////////////////////////////////////////////////
4 //
5 // class AliFlowAnalysis
6 //
7 // Flow Analysis
8 //
9 //
10 // S.Radomski@gsi.de
11 // Piotr.Skowronski@cern.ch
12 //
13 ///////////////////////////////////////////////////////////
14 /*********************************************************/
15
16 #include <TString.h>
17 #include <TParticle.h>
18
19 #include <AliStack.h>
20 #include <AliAOD.h>
21 #include <AliVAODParticle.h>
22 #include <AliAODParticleCut.h>
23
24 #include <AliESDtrack.h>
25 #include <AliESD.h>
26
27 ClassImp(AliFlowAnalysis)
28
29 AliFlowAnalysis::AliFlowAnalysis():
30  fPartCut(0x0)
31 {
32  //ctor
33 }
34 /*********************************************************/
35
36 AliFlowAnalysis::~AliFlowAnalysis()
37 {
38  //dtor
39   delete fPartCut;
40 }
41 /*********************************************************/
42
43 Int_t AliFlowAnalysis::Init()
44 {
45   //Initilizes anaysis
46   Info("Init","");
47   return 0; 
48 }
49 /*********************************************************/
50
51 Int_t AliFlowAnalysis::ProcessEvent(AliAOD* aodrec, AliAOD* aodsim)
52 {
53  
54   Info("ProcessEvent","Sim AOD address %#x",aodsim);
55   Double_t psi = 0, v2 = 0;
56   if (aodrec)
57    {
58      GetFlow(aodrec,v2,psi);
59      Info("ProcessEvent","Reconstructed Event: Event plane is %f, V2 is %f",psi,v2);
60    }  
61
62   if (aodsim)
63    {
64      GetFlow(aodsim,v2,psi);
65      Info("ProcessEvent","Simulated Event: Event plane is %f, V2 is %f",psi,v2);
66    }  
67   
68   return 0;
69   
70 }
71 /*********************************************************/
72
73 Int_t AliFlowAnalysis::Finish()
74 {
75   //Finish analysis and writes results
76   Info("Init","Finish");
77   return 0;
78 }
79 /*********************************************************/
80
81 Double_t AliFlowAnalysis::GetEventPlane(AliAOD* aod)
82 {
83   //returns event plane in degrees
84   if (aod == 0x0)
85    {
86      Error("AliFlowAnalysis::GetFlow","Pointer to AOD is NULL");
87      return -1.0;
88    }
89
90   Double_t psi;
91   Int_t mult = aod->GetNumberOfParticles();
92   
93   Double_t ssin = 0, scos = 0;
94
95   for (Int_t i=0; i<mult; i++) 
96    {
97      AliVAODParticle* aodtrack = aod->GetParticle(i);
98      if (aodtrack == 0x0)
99       {
100         Error("AliFlowAnalysis::GetEventPlane","Can not get track %d", i);
101         continue;
102       }
103      
104      if (fPartCut)
105       if (fPartCut->Pass(aodtrack))
106         continue;
107
108      Double_t phi = TMath::Pi()+TMath::ATan2(-aodtrack->Py(),-aodtrack->Px()); 
109      
110      ssin += TMath::Sin( 2.0 * phi );
111      scos += TMath::Cos( 2.0 * phi );
112    }
113
114   psi = atan2 (ssin, scos) / 2.0;
115   psi = psi * 180. / TMath::Pi(); 
116   
117   return psi;
118
119 }
120 /*********************************************************/
121
122 void AliFlowAnalysis::GetFlow(AliAOD* aod,Double_t& v2,Double_t& psi)
123 {
124 //returns flow parameters: v2 and event plane
125   if (aod == 0x0)
126    {
127      Error("AliFlowAnalysis::GetFlow","Pointer to AOD is NULL");
128      return;
129    }
130    
131   psi = GetEventPlane(aod);
132   Int_t mult = aod->GetNumberOfParticles();
133   
134   Double_t ssin = 0, scos = 0;
135
136   for (Int_t i=0; i<mult; i++) 
137    {
138      AliVAODParticle* aodtrack = aod->GetParticle(i);
139      if (aodtrack == 0x0)
140       {
141         Error("AliFlowAnalysis::GetEventPlane","Can not get track %d", i);
142         continue;
143       }
144      if (fPartCut)
145       if (fPartCut->Pass(aodtrack))
146         continue;
147       
148      Double_t phi = TMath::Pi()+TMath::ATan2(-aodtrack->Py(),-aodtrack->Px()); 
149      ssin += TMath::Sin( 2.0 * (phi - psi));
150      scos += TMath::Cos( 2.0 * (phi - psi));
151    }
152    
153   v2 = TMath::Hypot(ssin,scos);
154 }
155
156
157 /*********************************************************/
158
159 Double_t AliFlowAnalysis::GetEventPlane(AliESD* esd)
160 {
161   //returns event plane
162   if (esd == 0x0)
163    {
164      ::Error("AliFlowAnalysis::GetFlow","Pointer to ESD is NULL");
165      return -1.0;
166    }
167
168   Double_t psi;
169   Int_t mult = esd->GetNumberOfTracks();
170   
171   Double_t ssin = 0, scos = 0;
172
173   for (Int_t i=0; i<mult; i++) 
174    {
175      AliESDtrack* esdtrack = esd->GetTrack(i);
176      if (esdtrack == 0x0)
177       {
178         ::Error("AliFlowAnalysis::GetEventPlane","Can not get track %d", i);
179         continue;
180       }
181       
182      Double_t mom[3];//momentum
183      esdtrack->GetPxPyPz(mom); 
184      Double_t phi = TMath::Pi()+TMath::ATan2(-mom[1],-mom[0]); 
185      
186      ssin += TMath::Sin( 2.0 * phi );
187      scos += TMath::Cos( 2.0 * phi );
188    }
189
190   psi = atan2 (ssin, scos) / 2.0;
191   psi = psi * 180. / TMath::Pi(); 
192   
193   return psi;
194
195 }
196 /*********************************************************/
197
198 void AliFlowAnalysis::GetFlow(AliESD* esd,Double_t& v2,Double_t& psi)
199 {
200 //returns flow parameters: v2 and event plane
201   if (esd == 0x0)
202    {
203      ::Error("AliFlowAnalysis::GetFlow","Pointer to ESD is NULL");
204      return;
205    }
206    
207   psi = GetEventPlane(esd);
208   Int_t mult = esd->GetNumberOfTracks();
209   
210   Double_t ssin = 0, scos = 0;
211
212   for (Int_t i=0; i<mult; i++) 
213    {
214      AliESDtrack* esdtrack = esd->GetTrack(i);
215      if (esdtrack == 0x0)
216       {
217         ::Error("AliFlowAnalysis::GetEventPlane","Can not get track %d", i);
218         continue;
219       }
220       
221      Double_t mom[3];//momentum
222      esdtrack->GetPxPyPz(mom); 
223      Double_t phi = TMath::Pi()+TMath::ATan2(-mom[1],-mom[0]); 
224      
225      ssin += TMath::Sin( 2.0 * (phi - psi));
226      scos += TMath::Cos( 2.0 * (phi - psi));
227    }
228    
229   v2 = TMath::Hypot(ssin,scos);
230 }
231 /*********************************************************/