1 /**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //*****************************************************
17 // Class AliEventplane
18 // author: Alberica Toia, Johanna Gramling
19 //*****************************************************
20 /// A container for the event plane stored in AOD in ESD
23 #include "AliEventplane.h"
25 #include "AliVTrack.h"
26 #include "TObjArray.h"
28 #include "AliVEvent.h"
29 #include "AliVVZERO.h"
31 ClassImp(AliEventplane)
33 AliEventplane::AliEventplane() : TNamed("Eventplane", "Eventplane"),
37 fQContributionXsub1(0),
38 fQContributionYsub1(0),
39 fQContributionXsub2(0),
40 fQContributionYsub2(0),
47 fQContributionX = new TArrayF(0);
48 fQContributionY = new TArrayF(0);
49 fQContributionXsub1 = new TArrayF(0);
50 fQContributionYsub1 = new TArrayF(0);
51 fQContributionXsub2 = new TArrayF(0);
52 fQContributionYsub2 = new TArrayF(0);
53 for(Int_t i = 0; i < 11; ++i) {
54 fMeanX2[i]=fMeanY2[i]=0.;
55 fAPlus[i]=fAMinus[i]=1.;
56 fLambdaPlus[i]=fLambdaMinus[i]=0.;
61 AliEventplane::AliEventplane(const AliEventplane& ep) :
66 fQContributionXsub1(0),
67 fQContributionYsub1(0),
68 fQContributionXsub2(0),
69 fQContributionYsub2(0),
76 ((AliEventplane &) ep).CopyEP(*this);
79 AliEventplane& AliEventplane::operator=(const AliEventplane& ep)
81 /// Assignment operator
83 TNamed::operator=(ep);
84 ((AliEventplane &) ep).CopyEP(*this);
89 void AliEventplane::CopyEP(AliEventplane& ep) const
92 AliEventplane& target = (AliEventplane &) ep;
96 if(ep.fQContributionX)
98 *(ep.fQContributionX) = *fQContributionX;
102 ep.fQContributionX = new TArrayF(*fQContributionX);
107 if(ep.fQContributionX) delete ep.fQContributionX;
108 ep.fQContributionX = 0;
113 if(ep.fQContributionY)
115 *(ep.fQContributionY) = *fQContributionY;
119 ep.fQContributionY = new TArrayF(*fQContributionY);
124 if(ep.fQContributionY) delete ep.fQContributionY;
125 ep.fQContributionY = 0;
129 if(fQContributionXsub1)
131 if(ep.fQContributionXsub1)
133 *(ep.fQContributionXsub1) = *fQContributionXsub1;
137 ep.fQContributionXsub1 = new TArrayF(*fQContributionXsub1);
142 if(ep.fQContributionXsub1) delete ep.fQContributionXsub1;
143 ep.fQContributionXsub1 = 0;
146 if(fQContributionYsub1)
148 if(ep.fQContributionYsub1)
150 *(ep.fQContributionYsub1) = *fQContributionYsub1;
154 ep.fQContributionYsub1 = new TArrayF(*fQContributionYsub1);
159 if(ep.fQContributionYsub1) delete ep.fQContributionYsub1;
160 ep.fQContributionYsub1 = 0;
164 if(fQContributionXsub2)
166 if(ep.fQContributionXsub2)
168 *(ep.fQContributionXsub2) = *fQContributionXsub2;
172 ep.fQContributionXsub2 = new TArrayF(*fQContributionXsub2);
177 if(ep.fQContributionXsub2) delete ep.fQContributionXsub2;
178 ep.fQContributionXsub2 = 0;
181 if(fQContributionYsub2)
183 if(ep.fQContributionYsub2)
185 *(ep.fQContributionYsub2) = *fQContributionYsub2;
189 ep.fQContributionYsub2 = new TArrayF(*fQContributionYsub2);
194 if(ep.fQContributionYsub2) delete ep.fQContributionYsub2;
195 ep.fQContributionYsub2 = 0;
199 target.fEventplaneQ = fEventplaneQ;
202 target.fQVector = dynamic_cast<TVector2*> (fQVector->Clone());
204 target.fQsub1 = dynamic_cast<TVector2*> (fQsub1->Clone());
206 target.fQsub2 = dynamic_cast<TVector2*> (fQsub2->Clone());
208 target.fQsubRes = fQsubRes;
210 for(Int_t i = 0; i < 11; ++i) {
211 target.fMeanX2[i]=fMeanX2[i];
212 target.fMeanY2[i]=fMeanY2[i];
213 target.fAPlus[i]=fAPlus[i];
214 target.fAMinus[i]=fAMinus[i];
215 target.fLambdaPlus[i]=fLambdaPlus[i];
216 target.fLambdaMinus[i]=fLambdaMinus[i];
217 target.fCos8Psi[i]=fCos8Psi[i];
221 AliEventplane::~AliEventplane()
224 if (fQContributionX){
225 delete fQContributionX;
228 if (fQContributionY){
229 delete fQContributionY;
232 if (fQContributionXsub1){
233 delete fQContributionXsub1;
234 fQContributionXsub1 = 0;
236 if (fQContributionYsub1){
237 delete fQContributionYsub1;
238 fQContributionYsub1 = 0;
240 if (fQContributionXsub2){
241 delete fQContributionXsub2;
242 fQContributionXsub2 = 0;
244 if (fQContributionYsub2){
245 delete fQContributionYsub2;
246 fQContributionYsub2 = 0;
262 TVector2* AliEventplane::GetQVector()
267 Double_t AliEventplane::GetEventplane(const char *x, const AliVEvent *event, Int_t harmonic) const
270 Double_t qx = 0, qy = 0;
271 if(method.CompareTo("Q")==0) return fEventplaneQ;
272 else if(method.CompareTo("V0A")==0) return CalculateVZEROEventPlane(event, 8, harmonic, qx, qy);
273 else if(method.CompareTo("V0C")==0) return CalculateVZEROEventPlane(event, 9, harmonic, qx, qy);
274 else if(method.CompareTo("V0")==0) return CalculateVZEROEventPlane(event,10, harmonic, qx, qy);
279 Double_t AliEventplane::CalculateVZEROEventPlane(const AliVEvent * event, Int_t firstRing, Int_t lastRing, Int_t harmonic, Double_t &qxTot, Double_t &qyTot) const
281 // Calculate the VZERO event-plane by summing rings
282 // The flatenning is done on every ring separately
284 AliError("No Event received");
287 AliVVZERO *vzeroData = event->GetVZEROData();
289 AliError("Enable to get VZERO Data");
293 AliError("Required harmonic is less or equal to 0");
299 Double_t totMult = 0.;
300 for(Int_t ring = firstRing; ring <=lastRing; ++ring) {
302 Double_t multRing = 0.;
303 for(Int_t iCh = ring*8; iCh < (ring+1)*8; ++iCh) {
304 Double_t phi = TMath::Pi()/8. + TMath::Pi()/4.*(iCh%8);
305 Double_t mult = event->GetVZEROEqMultiplicity(iCh);
306 qx += mult*TMath::Cos(harmonic*phi);
307 qy += mult*TMath::Sin(harmonic*phi);
311 if (multRing < 1e-6) continue;
316 Double_t qxPrime = qx - fMeanX2[ring];
317 Double_t qyPrime = qy - fMeanY2[ring];
319 Double_t trans[2][2];
320 trans[0][0] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
321 trans[0][1] = -fLambdaMinus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
322 trans[1][0] = -fLambdaPlus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
323 trans[1][1] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
324 Double_t qxSeconde = trans[0][0]*qxPrime + trans[0][1]*qyPrime;
325 Double_t qySeconde = trans[1][0]*qxPrime + trans[1][1]*qyPrime;
327 Double_t qxTierce = qxSeconde/fAPlus[ring];
328 Double_t qyTierce = qySeconde/fAMinus[ring];
329 // 4th Fourier momenta in order to flatten the EP within a sector
330 Double_t psi = TMath::ATan2(qyTierce,qxTierce)/harmonic;
331 Double_t deltaPsi = (fCos8Psi[ring]*TMath::Sin(2.*4.*psi))/2.;
332 Double_t qxQuarte = qxTierce*TMath::Cos(2.*deltaPsi) - qyTierce*TMath::Sin(2.*deltaPsi);
333 Double_t qyQuarte = qxTierce*TMath::Sin(2.*deltaPsi) + qyTierce*TMath::Cos(2.*deltaPsi);
343 // In case of no hits return an invalid event-plane angle
344 if (totMult<1e-6) return -999;
346 return (TMath::ATan2(qyTot,qxTot)/harmonic);
349 Double_t AliEventplane::CalculateVZEROEventPlane(const AliVEvent * event, Int_t ring, Int_t harmonic, Double_t &qx, Double_t &qy) const
351 // Calculate the VZERO event-plane from an individual ring
352 // Ring 8 - total V0A (rings 4-7), Ring 9 - total V0C (rings 0-3)
353 // Ring 10 - total V0 (rings 0-7)
355 AliError("No Event received");
358 AliVVZERO *vzeroData = event->GetVZEROData();
360 AliError("Enable to get VZERO Data");
364 AliError("Required harmonic is less or equal to 0");
368 Int_t firstCh=-1,lastCh=-1;
373 else if (ring == 8) {
377 else if (ring == 9) {
381 else if (ring == 10) {
386 Double_t multRing = 0.;
387 for(Int_t iCh = firstCh; iCh < lastCh; ++iCh) {
388 Double_t phi = TMath::Pi()/8. + TMath::Pi()/4.*(iCh%8);
389 Double_t mult = event->GetVZEROEqMultiplicity(iCh);
390 qx += mult*TMath::Cos(harmonic*phi);
391 qy += mult*TMath::Sin(harmonic*phi);
395 // In case of no hits return an invalid event-plane angle
396 if (multRing < 1e-6) return -999;
400 Double_t qxPrime = qx - fMeanX2[ring];
401 Double_t qyPrime = qy - fMeanY2[ring];
403 Double_t trans[2][2];
404 trans[0][0] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
405 trans[0][1] = -fLambdaMinus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
406 trans[1][0] = -fLambdaPlus[ring]/(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
407 trans[1][1] = 1./(1.-fLambdaPlus[ring]*fLambdaMinus[ring]);
408 Double_t qxSeconde = trans[0][0]*qxPrime + trans[0][1]*qyPrime;
409 Double_t qySeconde = trans[1][0]*qxPrime + trans[1][1]*qyPrime;
411 Double_t qxTierce = qxSeconde/fAPlus[ring];
412 Double_t qyTierce = qySeconde/fAMinus[ring];
413 // 4th Fourier momenta in order to flatten the EP within a sector
414 Double_t psi = TMath::ATan2(qyTierce,qxTierce)/harmonic;
415 Double_t deltaPsi = (fCos8Psi[ring]*TMath::Sin(2.*4.*psi))/2.;
416 Double_t qxQuarte = qxTierce*TMath::Cos(2.*deltaPsi) - qyTierce*TMath::Sin(2.*deltaPsi);
417 Double_t qyQuarte = qxTierce*TMath::Sin(2.*deltaPsi) + qyTierce*TMath::Cos(2.*deltaPsi);
422 return (TMath::ATan2(qy,qx)/harmonic);
425 void AliEventplane::SetVZEROEPParams(Int_t ring,
426 Double_t meanX2, Double_t meanY2,
427 Double_t aPlus, Double_t aMinus,
428 Double_t lambdaPlus, Double_t lambdaMinus,
431 // Set the VZERO event-plane
432 // flatenning parameters
433 if (ring < 0 || ring >= 11) AliFatal(Form("Bad ring index (%d)",ring));
435 fMeanX2[ring] = meanX2;
436 fMeanY2[ring] = meanY2;
437 fAPlus[ring] = aPlus;
438 fAMinus[ring] = aMinus;
439 fLambdaPlus[ring] = lambdaPlus;
440 fLambdaMinus[ring] = lambdaMinus;
441 fCos8Psi[ring] = cos8Psi;
444 TVector2* AliEventplane::GetQsub1()
449 TVector2* AliEventplane::GetQsub2()
454 Double_t AliEventplane::GetQsubRes()
459 Bool_t AliEventplane::IsEventInEventplaneClass(Double_t a, Double_t b, const char *x)
462 if ((method.CompareTo("Q")==0) && (fEventplaneQ >=a && fEventplaneQ < b)) return kTRUE;
466 Double_t AliEventplane::GetQContributionX(AliVTrack* track)
468 return fQContributionX->GetAt(track->GetID());
471 Double_t AliEventplane::GetQContributionY(AliVTrack* track)
473 return fQContributionY->GetAt(track->GetID());
476 Double_t AliEventplane::GetQContributionXsub1(AliVTrack* track)
478 return fQContributionXsub1->GetAt(track->GetID());
481 Double_t AliEventplane::GetQContributionYsub1(AliVTrack* track)
483 return fQContributionYsub1->GetAt(track->GetID());
486 Double_t AliEventplane::GetQContributionXsub2(AliVTrack* track)
488 return fQContributionXsub2->GetAt(track->GetID());
491 Double_t AliEventplane::GetQContributionYsub2(AliVTrack* track)
493 return fQContributionYsub2->GetAt(track->GetID());
496 void AliEventplane::Reset()
498 delete fQVector; fQVector=0;
499 fQContributionX->Reset();
500 fQContributionY->Reset();
501 fQContributionXsub1->Reset();
502 fQContributionYsub1->Reset();
503 fQContributionXsub2->Reset();
504 fQContributionYsub2->Reset();
506 delete fQsub1; fQsub1=0;
507 delete fQsub2; fQsub2=0;