Introduction of VZERO event-plane selection task that can be used in order to flatten...
[u/mrichter/AliRoot.git] / STEER / STEERBase / AliEventplane.cxx
CommitLineData
ce7adfe9 1/**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//*****************************************************
17// Class AliEventplane
18// author: Alberica Toia, Johanna Gramling
19//*****************************************************
20/// A container for the event plane stored in AOD in ESD
21
fd5c821a 22#include "AliLog.h"
ce7adfe9 23#include "AliEventplane.h"
24#include "TVector2.h"
aa260dd1 25#include "AliVTrack.h"
ce7adfe9 26#include "TObjArray.h"
27#include "TArrayF.h"
fd5c821a 28#include "AliVEvent.h"
29#include "AliVVZERO.h"
ce7adfe9 30
31ClassImp(AliEventplane)
32
33AliEventplane::AliEventplane() : TNamed("Eventplane", "Eventplane"),
34 fQVector(0),
35 fQContributionX(0),
36 fQContributionY(0),
e51055a0 37 fQContributionXsub1(0),
38 fQContributionYsub1(0),
39 fQContributionXsub2(0),
40 fQContributionYsub2(0),
3e64f60c 41 fEventplaneQ(-1),
ce7adfe9 42 fQsub1(0),
43 fQsub2(0),
44 fQsubRes(0)
45{
46 /// constructor
3e64f60c 47 fQContributionX = new TArrayF(0);
48 fQContributionY = new TArrayF(0);
e51055a0 49 fQContributionXsub1 = new TArrayF(0);
50 fQContributionYsub1 = new TArrayF(0);
51 fQContributionXsub2 = new TArrayF(0);
52 fQContributionYsub2 = new TArrayF(0);
0643f5a7 53 for(Int_t i = 0; i < 8; ++i) {
54 fMeanX2[i]=fMeanY2[i]=0.;
55 fAPlus[i]=fAMinus[i]=1.;
56 fLambdaPlus[i]=fLambdaMinus[i]=0.;
57 fCos8Psi[i]=0.;
58 }
ce7adfe9 59}
60
61AliEventplane::AliEventplane(const AliEventplane& ep) :
62 TNamed(),
63 fQVector(0),
64 fQContributionX(0),
65 fQContributionY(0),
e51055a0 66 fQContributionXsub1(0),
67 fQContributionYsub1(0),
68 fQContributionXsub2(0),
69 fQContributionYsub2(0),
ce7adfe9 70 fEventplaneQ(0),
71 fQsub1(0),
72 fQsub2(0),
73 fQsubRes(0)
74{
75 /// Copy constructor
76 ((AliEventplane &) ep).CopyEP(*this);
77}
78
79AliEventplane& AliEventplane::operator=(const AliEventplane& ep)
80{
81 /// Assignment operator
82 if (this!=&ep)
83 ((AliEventplane &) ep).CopyEP(*this);
84
85 return *this;
86}
87
88void AliEventplane::CopyEP(AliEventplane& ep) const
89{ // copy function
90
91 AliEventplane& target = (AliEventplane &) ep;
92 if (fQContributionX)
93 target.fQContributionX = fQContributionX;
94 if (fQContributionY)
95 target.fQContributionY = fQContributionY;
e51055a0 96 if (fQContributionXsub1)
97 target.fQContributionXsub1 = fQContributionXsub1;
98 if (fQContributionYsub1)
99 target.fQContributionYsub1 = fQContributionYsub1;
100 if (fQContributionXsub2)
101 target.fQContributionXsub2 = fQContributionXsub2;
102 if (fQContributionYsub2)
103 target.fQContributionYsub2 = fQContributionYsub2;
ce7adfe9 104 if (fEventplaneQ)
105 target.fEventplaneQ = fEventplaneQ;
106 if (fQVector)
107 target.fQVector = dynamic_cast<TVector2*> (fQVector->Clone());
108 if (fQsub1)
109 target.fQsub1 = dynamic_cast<TVector2*> (fQsub1->Clone());
110 if (fQsub2)
111 target.fQsub2 = dynamic_cast<TVector2*> (fQsub2->Clone());
112 if (fQsubRes)
113 target.fQsubRes = fQsubRes;
0643f5a7 114
115 for(Int_t i = 0; i < 8; ++i) {
116 target.fMeanX2[i]=fMeanX2[i];
117 target.fMeanY2[i]=fMeanY2[i];
118 target.fAPlus[i]=fAPlus[i];
119 target.fAMinus[i]=fAMinus[i];
120 target.fLambdaPlus[i]=fLambdaPlus[i];
121 target.fLambdaMinus[i]=fLambdaMinus[i];
122 target.fCos8Psi[i]=fCos8Psi[i];
123 }
ce7adfe9 124}
125
126AliEventplane::~AliEventplane()
127{
128 /// destructor
129 if (fQContributionX){
130 delete fQContributionX;
131 fQContributionX = 0;
132 }
133 if (fQContributionY){
134 delete fQContributionY;
135 fQContributionY = 0;
136 }
e51055a0 137 if (fQContributionXsub1){
138 delete fQContributionXsub1;
139 fQContributionXsub1 = 0;
140 }
141 if (fQContributionYsub1){
142 delete fQContributionYsub1;
143 fQContributionYsub1 = 0;
144 }
145 if (fQContributionXsub2){
146 delete fQContributionXsub2;
147 fQContributionXsub2 = 0;
148 }
149 if (fQContributionYsub2){
150 delete fQContributionYsub2;
151 fQContributionYsub2 = 0;
152 }
ce7adfe9 153 if (fQVector){
154 delete fQVector;
155 fQVector = 0;
156 }
157 if (fQsub1){
158 delete fQsub1;
159 fQsub1 = 0;
160 }
161 if (fQsub2){
162 delete fQsub2;
163 fQsub2 = 0;
164 }
165}
166
167TVector2* AliEventplane::GetQVector()
168{
169 return fQVector;
170}
171
fd5c821a 172Double_t AliEventplane::GetEventplane(const char *x, const AliVEvent *event, Int_t harmonic) const
ce7adfe9 173{
174 TString method = x;
175 if(method.CompareTo("Q")==0) return fEventplaneQ;
fd5c821a 176 else if(method.CompareTo("V0A")==0) return CalculateVZEROEventPlane(event, 4, 7, harmonic);
177 else if(method.CompareTo("V0C")==0) return CalculateVZEROEventPlane(event, 0, 3, harmonic);
178 else if(method.CompareTo("V0")==0) return CalculateVZEROEventPlane(event, 0, 7, harmonic);
179
180 return -1000.;
181}
182
183Double_t AliEventplane::CalculateVZEROEventPlane(const AliVEvent * event, Int_t firstRing, Int_t lastRing, Int_t harmonic) const
184{
185 if(!event) {
186 AliError("No Event received");
187 return -1000.;
188 }
189 AliVVZERO *vzeroData = event->GetVZEROData();
190 if(!vzeroData) {
191 AliError("Enable to get VZERO Data");
192 return -1000.;
193 }
194 if(harmonic <= 0) {
195 AliError("Required harmonic is less or equal to 0");
196 return -1000.;
197 }
198
0643f5a7 199 Double_t qxTot=0., qyTot=0.;
200 Double_t qx, qy;
201 Double_t totMult = 0.;
202 for(Int_t ring = firstRing; ring <=lastRing; ++ring) {
203 qx=qy=0.;
204 Double_t trans[2][2];
205 trans[0][0] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
206 trans[0][1] = -fLambdaMinus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
207 trans[1][0] = -fLambdaPlus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
208 trans[1][1] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
209 Double_t multRing = 0.;
210 for(Int_t iCh = ring*8; iCh < (ring+1)*8; ++iCh) {
211 Double_t phi = TMath::Pi()/8. + TMath::Pi()/4.*(iCh%8);
212 Double_t mult = event->GetVZEROEqMultiplicity(iCh);
213 qx += mult*TMath::Cos(harmonic*phi);
214 qy += mult*TMath::Sin(harmonic*phi);
215 multRing += mult;
216 }
217
218 if (multRing < 1e-6) continue;
219 totMult += multRing;
220
221 if (harmonic == 2) {
222 // Recentering
223 Double_t qxPrime = qx - fMeanX2[ring];
224 Double_t qyPrime = qy - fMeanY2[ring];
225 // Twist
226 Double_t qxSeconde = trans[0][0]*qxPrime + trans[0][1]*qyPrime;
227 Double_t qySeconde = trans[1][0]*qxPrime + trans[1][1]*qyPrime;
228 // Rescaling
229 Double_t qxTierce = qxSeconde/fAPlus[ring];
230 Double_t qyTierce = qySeconde/fAMinus[ring];
231 qxTot += qxTierce;
232 qyTot += qyTierce;
233 }
234 else {
235 qxTot += qx;
236 qyTot += qy;
237 }
238 }
fd5c821a 239
0643f5a7 240 // In case of no hits return an invalid event-plane angle
241 if (totMult<1e-6) return -999;
fd5c821a 242
0643f5a7 243 // 4th Fourier momenta in order to flatten the EP within a sector
244 Double_t psi = TMath::ATan2(qyTot,qxTot)/harmonic;
245 if ((harmonic == 2) && (firstRing == lastRing)) {
246 psi += (fCos8Psi[firstRing]*TMath::Sin(2.*4.*psi))/2.;
247 if (psi > TMath::Pi()/2) psi -= TMath::Pi();
248 if (psi <-TMath::Pi()/2) psi += TMath::Pi();
fd5c821a 249 }
0643f5a7 250
251 return psi;
ce7adfe9 252}
253
0643f5a7 254void AliEventplane::SetVZEROEPParams(Int_t ring,
255 Double_t meanX2, Double_t meanY2,
256 Double_t aPlus, Double_t aMinus,
257 Double_t lambdaPlus, Double_t lambdaMinus,
258 Double_t cos8Psi)
259{
260 // Set the VZERO event-plane
261 // flatenning parameters
262 if (ring < 0 || ring >=8) AliFatal(Form("Bad ring index (%d)",ring));
263
264 fMeanX2[ring] = meanX2;
265 fMeanY2[ring] = meanY2;
266 fAPlus[ring] = aPlus;
267 fAMinus[ring] = aMinus;
268 fLambdaPlus[ring] = lambdaPlus;
269 fLambdaMinus[ring] = lambdaMinus;
270 fCos8Psi[ring] = cos8Psi;
271}
fd5c821a 272
ce7adfe9 273TVector2* AliEventplane::GetQsub1()
274{
275 return fQsub1;
276}
277
278TVector2* AliEventplane::GetQsub2()
279{
280 return fQsub2;
281}
282
283Double_t AliEventplane::GetQsubRes()
284{
285 return fQsubRes;
286}
287
288Bool_t AliEventplane::IsEventInEventplaneClass(Double_t a, Double_t b, const char *x)
289{
290 TString method = x;
291 if ((method.CompareTo("Q")==0) && (fEventplaneQ >=a && fEventplaneQ < b)) return kTRUE;
292 else return kFALSE;
293}
294
aa260dd1 295Double_t AliEventplane::GetQContributionX(AliVTrack* track)
ce7adfe9 296{
297 return fQContributionX->GetAt(track->GetID());
298}
299
aa260dd1 300Double_t AliEventplane::GetQContributionY(AliVTrack* track)
ce7adfe9 301{
302 return fQContributionY->GetAt(track->GetID());
303}
3e64f60c 304
e51055a0 305Double_t AliEventplane::GetQContributionXsub1(AliVTrack* track)
306{
307 return fQContributionXsub1->GetAt(track->GetID());
308}
309
310Double_t AliEventplane::GetQContributionYsub1(AliVTrack* track)
311{
312 return fQContributionYsub1->GetAt(track->GetID());
313}
314
315Double_t AliEventplane::GetQContributionXsub2(AliVTrack* track)
316{
317 return fQContributionXsub2->GetAt(track->GetID());
318}
319
320Double_t AliEventplane::GetQContributionYsub2(AliVTrack* track)
321{
322 return fQContributionYsub2->GetAt(track->GetID());
323}
324
3e64f60c 325void AliEventplane::Reset()
326{
327 delete fQVector; fQVector=0;
328 fQContributionX->Reset();
329 fQContributionY->Reset();
e51055a0 330 fQContributionXsub1->Reset();
331 fQContributionYsub1->Reset();
332 fQContributionXsub2->Reset();
333 fQContributionYsub2->Reset();
3e64f60c 334 fEventplaneQ = -1;
335 delete fQsub1; fQsub1=0;
336 delete fQsub2; fQsub2=0;
337 fQsubRes = 0;
338}