PWGDQ/dielectron/AliReducedEvent.cxx

   1 /*
   2 ***********************************************************
   3   Implementation of reduced ESD information classes for
   4   quick analysis.
   5   Contact: i.c.arsene@gsi.de, i.c.arsene@cern.ch
   6   2012/06/21
   7   *********************************************************
   8 */
   9
  10 #ifndef ALIREDUCEDEVENT_H
  11 #include "AliReducedEvent.h"
  12 #endif
  13
  14 #include <TMath.h>
  15 #include <TClonesArray.h>
  16
  17 ClassImp(AliReducedTrack)
  18 ClassImp(AliReducedPair)
  19 ClassImp(AliReducedEvent)
  20 ClassImp(AliReducedEventFriend)
  21 ClassImp(AliReducedCaloCluster)
  22
  23 TClonesArray* AliReducedEvent::fgTracks = 0;
  24 TClonesArray* AliReducedEvent::fgCandidates = 0;
  25 TClonesArray* AliReducedEvent::fgCaloClusters = 0;
  26
  27 //_______________________________________________________________________________
  28 AliReducedTrack::AliReducedTrack() :
  29   fTrackId(-1),
  30   fStatus(0),
  31   fGlobalPhi(0.0),
  32   fGlobalPt(0.0),
  33   fGlobalEta(0.0),
  34   fTPCPhi(0.0),
  35   fTPCPt(0.0),
  36   fTPCEta(0.0),
  37   fMomentumInner(0.0),
  38   fDCA(),
  39   fITSclusterMap(0),
  40   fITSsignal(0.0),
  41   fITSnSig(),
  42   fTPCNcls(0),
  43   fTPCCrossedRows(0),
  44   fTPCNclsF(0),
  45   fTPCNclsIter1(0),
  46   fTPCClusterMap(0),
  47   fTPCsignal(0),
  48   fTPCnSig(),
  49   fTOFbeta(0.0),
  50   fTOFnSig(),
  51   fTRDntracklets(),
  52   fTRDpid(),
  53   fCaloClusterId(-999),
  54   fBayesPID(),
  55   fFlags(0),
  56   fMoreFlags(0)
  57 {
  58   //
  59   // Constructor
  60   //
  61   fDCA[0] = 0.0; fDCA[1]=0.0;
  62   for(Int_t i=0; i<4; ++i) {fTPCnSig[i]=-999.; fTOFnSig[i]=-999.; fITSnSig[i]=-999.;}
  63   for(Int_t i=0; i<3; ++i) {fBayesPID[i]=-999.;}
  64   fTRDpid[0]=-999.; fTRDpid[1]=-999.;
  65 }
  66
  67
  68 //_______________________________________________________________________________
  69 AliReducedTrack::~AliReducedTrack()
  70 {
  71   //
  72   // De-Constructor
  73   //
  74 }
  75
  76
  77 //_______________________________________________________________________________
  78 AliReducedPair::AliReducedPair() :
  79   fCandidateId(-1),
  80   fPairType(-1),
  81   fLegIds(),
  82   fMass(),
  83   fPhi(0.0),
  84   fPt(0.0),
  85   fEta(0.0),
  86   fLxy(0.0),
  87   fLxyErr(0.0),
  88   fPointingAngle(0.0),
  89   fChisquare(0.0),
  90   fMCid(0)
  91 {
  92   //
  93   // Constructor
  94   //
  95   fLegIds[0] = -1; fLegIds[1] = -1;
  96   fMass[0]=-999.; fMass[1]=-999.; fMass[2]=-999.; fMass[3]=-999.;
  97 }
  98
  99
 100 //_______________________________________________________________________________
 101 AliReducedPair::AliReducedPair(const AliReducedPair &c) :
 102   TObject(c),
 103   fCandidateId(c.CandidateId()),
 104   fPairType(c.PairType()),
 105   fLegIds(),
 106   fMass(),
 107   fPhi(c.Phi()),
 108   fPt(c.Pt()),
 109   fEta(c.Eta()),
 110   fLxy(c.Lxy()),
 111   fLxyErr(c.LxyErr()),
 112   fPointingAngle(c.PointingAngle()),
 113   fChisquare(c.Chi2()),
 114   fMCid(c.MCid())
 115 {
 116   //
 117   // copy constructor
 118   //
 119   fLegIds[0] = c.LegId(0);
 120   fLegIds[1] = c.LegId(1);
 121   fMass[0] = c.Mass(0); fMass[1] = c.Mass(1); fMass[2] = c.Mass(2); fMass[3] = c.Mass(3);
 122 }
 123
 124
 125 //_______________________________________________________________________________
 126 AliReducedPair::~AliReducedPair() {
 127   //
 128   // destructor
 129   //
 130 }
 131
 132 //____________________________________________________________________________
 133 AliReducedEvent::AliReducedEvent() :
 134   TObject(),
 135   fRunNo(0),
 136   fBC(0),
 137   fTimeStamp(0),
 138   fEventType(0),
 139   fTriggerMask(0),
 140   fIsPhysicsSelection(kTRUE),
 141   fIsSPDPileup(kTRUE),
 142   fVtx(),
 143   fNVtxContributors(0),
 144   fVtxTPC(),
 145   fNVtxTPCContributors(0),
 146   fNpileupSPD(0),
 147   fNpileupTracks(0),
 148   fNPMDtracks(0),
 149   fNTRDtracks(0),
 150   fNTRDtracklets(0),
 151   fCentrality(),
 152   fCentQuality(0),
 153   fNV0candidates(),
 154   fNDielectronCandidates(0),
 155   fNtracks(),
 156   fSPDntracklets(0),
 157   fVZEROMult(),
 158   fZDCnEnergy(),
 159   fTracks(0x0),
 160   fCandidates(0x0),
 161   fNCaloClusters(0),
 162   fCaloClusters(0x0)
 163 //fFMDMult()
 164 {
 165   //
 166   // Constructor
 167   //
 168   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.;}
 169   for(Int_t i=0; i<4; ++i) fCentrality[i]=-1.;
 170   fNV0candidates[0]=0; fNV0candidates[1]=0;
 171   fNtracks[0]=0; fNtracks[1]=0;
 172   for(Int_t i=0; i<16; ++i) fSPDntrackletsEta[i]=0;
 173   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i]=0;
 174   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 175   for(Int_t i=0; i<8; ++i) fZDCnEnergy[i]=0.0;
 176 }
 177
 178
 179 //____________________________________________________________________________
 180 AliReducedEvent::AliReducedEvent(const Char_t* /*name*/) :
 181   TObject(),
 182   fRunNo(0),
 183   fBC(0),
 184   fTimeStamp(0),
 185   fEventType(0),
 186   fTriggerMask(0),
 187   fIsPhysicsSelection(kTRUE),
 188   fIsSPDPileup(kTRUE),
 189   fVtx(),
 190   fNVtxContributors(0),
 191   fVtxTPC(),
 192   fNVtxTPCContributors(0),
 193   fNpileupSPD(0),
 194   fNpileupTracks(0),
 195   fNPMDtracks(0),
 196   fNTRDtracks(0),
 197   fNTRDtracklets(0),
 198   fCentrality(),
 199   fCentQuality(0),
 200   fNV0candidates(),
 201   fNDielectronCandidates(0),
 202   fNtracks(),
 203   fSPDntracklets(0),
 204   fVZEROMult(),
 205   fZDCnEnergy(),
 206   fTracks(0x0),
 207   fCandidates(0x0),
 208   fNCaloClusters(0),
 209   fCaloClusters(0x0)
 210 //fFMDMult()
 211 {
 212   //
 213   // Constructor
 214   //
 215   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.;}
 216   for(Int_t i=0; i<4; ++i) fCentrality[i]=-1.;
 217   fNV0candidates[0]=0; fNV0candidates[1]=0;
 218   fNtracks[0]=0; fNtracks[1]=0;
 219   for(Int_t i=0; i<16; ++i) fSPDntrackletsEta[i]=0;
 220   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i]=0;
 221   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 222   for(Int_t i=0; i<8; ++i) fZDCnEnergy[i]=0.0;
 223
 224   if(!fgTracks) fgTracks = new TClonesArray("AliReducedTrack", 100000);
 225   fTracks = fgTracks;
 226   if(!fgCandidates) fgCandidates = new TClonesArray("AliReducedPair", 100000);
 227   fCandidates = fgCandidates;
 228   if(!fgCaloClusters) fgCaloClusters = new TClonesArray("AliReducedCaloCluster", 50000);
 229   fCaloClusters = fgCaloClusters;
 230 }
 231
 232
 233 //____________________________________________________________________________
 234 AliReducedEvent::~AliReducedEvent()
 235 {
 236   //
 237   // De-Constructor
 238   //
 239   //ClearEvent();
 240 }
 241
 242
 243 //____________________________________________________________________________
 244 Float_t AliReducedEvent::MultVZEROA() const
 245 {
 246   //
 247   // Total VZERO multiplicity in A side
 248   //
 249   Float_t mult=0.0;
 250   for(Int_t i=32;i<64;++i)
 251     mult += fVZEROMult[i];
 252   return mult;
 253 }
 254
 255
 256 //____________________________________________________________________________
 257 Float_t AliReducedEvent::MultVZEROC() const
 258 {
 259   //
 260   // Total VZERO multiplicity in C side
 261   //
 262   Float_t mult=0.0;
 263   for(Int_t i=0;i<32;++i)
 264     mult += fVZEROMult[i];
 265   return mult;
 266 }
 267
 268
 269 //____________________________________________________________________________
 270 Float_t AliReducedEvent::MultVZERO() const
 271 {
 272   //
 273   // Total VZERO multiplicity
 274   //
 275   return MultVZEROA()+MultVZEROC();
 276 }
 277
 278
 279 //____________________________________________________________________________
 280 Float_t AliReducedEvent::MultRingVZEROA(Int_t ring) const
 281 {
 282   //
 283   //  VZERO multiplicity in a given ring on A side
 284   //
 285   if(ring<0 || ring>3) return -1.0;
 286
 287   Float_t mult=0.0;
 288   for(Int_t i=32+8*ring; i<32+8*(ring+1); ++i)
 289     mult += fVZEROMult[i];
 290   return mult;
 291 }
 292
 293
 294 //____________________________________________________________________________
 295 Float_t AliReducedEvent::MultRingVZEROC(Int_t ring) const
 296 {
 297   //
 298   //  VZERO multiplicity in a given ring on C side
 299   //
 300   if(ring<0 || ring>3) return -1.0;
 301
 302   Float_t mult=0.0;
 303   for(Int_t i=8*ring; i<8*(ring+1); ++i)
 304     mult += fVZEROMult[i];
 305   return mult;
 306 }
 307
 308 //_____________________________________________________________________________
 309 void AliReducedEvent::ClearEvent() {
 310   //
 311   // clear the event
 312   //
 313   if(fTracks) fTracks->Clear("C");
 314   if(fCandidates) fCandidates->Clear("C");
 315   if(fCaloClusters) fCaloClusters->Clear("C");
 316   fRunNo = 0;
 317   fBC = 0;
 318   fTimeStamp = 0;
 319   fEventType = 0;
 320   fTriggerMask = 0;
 321   fIsPhysicsSelection = kTRUE;
 322   fIsSPDPileup = kFALSE;
 323   fCentQuality = 0;
 324   fNV0candidates[0] = 0; fNV0candidates[1] = 0;
 325   fNpileupSPD=0;
 326   fNpileupTracks=0;
 327   fNPMDtracks=0;
 328   fNTRDtracks=0;
 329   fNTRDtracklets=0;
 330   fNDielectronCandidates = 0;
 331   fNtracks[0] = 0; fNtracks[1] = 0;
 332   for(Int_t i=0; i<16; ++i) fSPDntrackletsEta[i] = 0;
 333   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i] = 0;
 334   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.; fCentrality[i]=-1.;}
 335   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 336   for(Int_t i=0; i<8; ++i) fZDCnEnergy[i]=0.0;
 337 }
 338
 339
 340 //_______________________________________________________________________________
 341 AliReducedEventFriend::AliReducedEventFriend() :
 342  fQvector(),
 343  fEventPlaneStatus()
 344 {
 345   //
 346   // default constructor
 347   //
 348   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 349     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 350       fEventPlaneStatus[idet][ih] = 0;
 351       for(Int_t ic=0; ic<2; ++ic)
 352         fQvector[idet][ih][ic] = 0.0;
 353     }
 354   }
 355 }
 356
 357
 358 //____________________________________________________________________________
 359 AliReducedEventFriend::~AliReducedEventFriend()
 360 {
 361   //
 362   // De-Constructor
 363   //
 364   ClearEvent();
 365 }
 366
 367
 368 //_____________________________________________________________________________
 369 void AliReducedEventFriend::ClearEvent() {
 370   //
 371   // clear the event
 372   //
 373   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 374     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 375       fEventPlaneStatus[idet][ih] = 0;
 376       for(Int_t ic=0; ic<2; ++ic)
 377         fQvector[idet][ih][ic] = 0.0;
 378     }
 379   }
 380 }
 381
 382
 383 //____________________________________________________________________________
 384 void AliReducedEventFriend::CopyEvent(const AliReducedEventFriend* event) {
 385   //
 386   // copy information from another event to this one
 387   //
 388   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 389     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 390       fQvector[idet][ih][0] = event->Qx(idet, ih+1);
 391       fQvector[idet][ih][1] = event->Qy(idet, ih+1);
 392       fEventPlaneStatus[idet][ih] = event->GetEventPlaneStatus(idet, ih+1);
 393     }
 394   }
 395 }
 396
 397
 398 //_____________________________________________________________________________
 399 AliReducedCaloCluster::AliReducedCaloCluster() :
 400  fType(kUndefined),
 401  fEnergy(-999.),
 402  fTrackDx(-999.),
 403  fTrackDz(-999.),
 404  fM20(-999.),
 405  fM02(-999.),
 406  fDispersion(-999.)
 407 {
 408   //
 409   // default constructor
 410   //
 411 }
 412
 413
 414 //_____________________________________________________________________________
 415 AliReducedCaloCluster::~AliReducedCaloCluster()
 416 {
 417   //
 418   // destructor
 419   //
 420 }
 421
 422
 423 //_______________________________________________________________________________
 424 void AliReducedEvent::GetQvector(Double_t Qvec[][2], Int_t det,
 425                                  Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 426                                  Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 427   //
 428   // Get the event plane for a specified detector
 429   //
 430   if(det==AliReducedEventFriend::kTPC ||
 431      det==AliReducedEventFriend::kTPCptWeights ||
 432      det==AliReducedEventFriend::kTPCpos ||
 433      det==AliReducedEventFriend::kTPCneg) {
 434     GetTPCQvector(Qvec, det, etaMin, etaMax, IsTrackSelected);
 435     return;
 436   }
 437   if(det==AliReducedEventFriend::kVZEROA ||
 438      det==AliReducedEventFriend::kVZEROC) {
 439     GetVZEROQvector(Qvec, det);
 440     return;
 441   }
 442   if(det==AliReducedEventFriend::kZDCA ||
 443      det==AliReducedEventFriend::kZDCC) {
 444     GetZDCQvector(Qvec, det);
 445     return;
 446   }
 447   if(det==AliReducedEventFriend::kFMD) {
 448     //TODO implementation
 449     return;
 450   }
 451   return;
 452 }
 453
 454
 455 //_________________________________________________________________
 456 Int_t AliReducedEvent::GetTPCQvector(Double_t Qvec[][2], Int_t det,
 457                                      Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 458                                      Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 459   //
 460   // Construct the event plane using tracks in the barrel
 461   //
 462   if(!(det==AliReducedEventFriend::kTPC ||
 463        det==AliReducedEventFriend::kTPCpos ||
 464        det==AliReducedEventFriend::kTPCneg))
 465     return 0;
 466   Int_t nUsedTracks = 0;
 467   Short_t charge = 0;
 468   AliReducedTrack* track=0x0;
 469   Double_t weight=0.0; Double_t absWeight = 0.0; Double_t x=0.0; Double_t y=0.0;
 470   TIter nextTrack(fTracks);
 471   while((track=static_cast<AliReducedTrack*>(nextTrack()))) {
 472     if(track->Eta()<etaMin) continue;
 473     if(track->Eta()>etaMax) continue;
 474     charge = track->Charge();
 475     if(det==AliReducedEventFriend::kTPCpos && charge<0) continue;
 476     if(det==AliReducedEventFriend::kTPCneg && charge>0) continue;
 477
 478     if(IsTrackSelected && !IsTrackSelected(track)) continue;
 479     absWeight = 1.0;
 480     if(det==AliReducedEventFriend::kTPCptWeights) {
 481       absWeight = track->Pt();
 482       if(absWeight>2.0) absWeight = 2.0;    // pt is the weight used for the event plane
 483     }
 484     weight = absWeight;
 485
 486     ++nUsedTracks;
 487     x = TMath::Cos(track->Phi());
 488     y = TMath::Sin(track->Phi());
 489     //  1st harmonic
 490     Qvec[0][0] += weight*x;
 491     Qvec[0][1] += weight*y;
 492     //  2nd harmonic
 493     Qvec[1][0] += absWeight*(2.0*TMath::Power(x,2.0)-1);
 494     Qvec[1][1] += absWeight*(2.0*x*y);
 495     //  3rd harmonic
 496     Qvec[2][0] += weight*(4.0*TMath::Power(x,3.0)-3.0*x);
 497     Qvec[2][1] += weight*(3.0*y-4.0*TMath::Power(y,3.0));
 498     //  4th harmonic
 499     Qvec[3][0] += absWeight*(1.0-8.0*TMath::Power(x*y,2.0));
 500     Qvec[3][1] += absWeight*(4.0*x*y-8.0*x*TMath::Power(y,3.0));
 501     //  5th harmonic
 502     Qvec[4][0] += weight*(16.0*TMath::Power(x,5.0)-20.0*TMath::Power(x, 3.0)+5.0*x);
 503     Qvec[4][1] += weight*(16.0*TMath::Power(y,5.0)-20.0*TMath::Power(y, 3.0)+5.0*y);
 504     //  6th harmonic
 505     Qvec[5][0] += absWeight*(32.0*TMath::Power(x,6.0)-48.0*TMath::Power(x, 4.0)+18.0*TMath::Power(x, 2.0)-1.0);
 506     Qvec[5][1] += absWeight*(x*y*(32.0*TMath::Power(y,4.0)-32.0*TMath::Power(y, 2.0)+6.0));
 507   }
 508   return nUsedTracks;
 509 }
 510
 511
 512 //____________________________________________________________________________
 513 void AliReducedEvent::SubtractParticleFromQvector(
 514         AliReducedTrack* particle, Double_t Qvec[][2], Int_t det,
 515         Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 516         Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 517   //
 518   // subtract a particle from the event Q-vector
 519   //
 520   Float_t eta = particle->Eta();
 521   if(eta<etaMin) return;
 522   if(eta>etaMax) return;
 523
 524   Float_t charge = particle->Charge();
 525   if(det==AliReducedEventFriend::kTPCpos && charge<0) return;
 526   if(det==AliReducedEventFriend::kTPCneg && charge>0) return;
 527
 528   if(IsTrackSelected && !IsTrackSelected(particle)) return;
 529
 530   Double_t weight=0.0; Double_t absWeight = 0.0;
 531   if(det==AliReducedEventFriend::kTPCptWeights) {
 532     absWeight = particle->Pt();
 533     if(absWeight>2.0) absWeight = 2.0;
 534   }
 535   weight = absWeight;
 536   //  if(eta<0.0) weight *= -1.0;
 537
 538   Float_t x = TMath::Cos(particle->Phi());
 539   Float_t y = TMath::Sin(particle->Phi());
 540
 541   //  1st harmonic
 542   Qvec[0][0] -= weight*x;
 543   Qvec[0][1] -= weight*y;
 544   //  2nd harmonic
 545   Qvec[1][0] -= absWeight*(2.0*TMath::Power(x,2.0)-1);
 546   Qvec[1][1] -= absWeight*(2.0*x*y);
 547   //  3rd harmonic
 548   Qvec[2][0] -= weight*(4.0*TMath::Power(x,3.0)-3.0*x);
 549   Qvec[2][1] -= weight*(3.0*y-4.0*TMath::Power(y,3.0));
 550   //  4th harmonic
 551   Qvec[3][0] -= absWeight*(1.0-8.0*TMath::Power(x*y,2.0));
 552   Qvec[3][1] -= absWeight*(4.0*x*y-8.0*x*TMath::Power(y,3.0));
 553   //  5th harmonic
 554   Qvec[4][0] -= weight*(16.0*TMath::Power(x,5.0)-20.0*TMath::Power(x, 3.0)+5.0*x);
 555   Qvec[4][1] -= weight*(16.0*TMath::Power(y,5.0)-20.0*TMath::Power(y, 3.0)+5.0*y);
 556   //  6th harmonic
 557   Qvec[5][0] -= absWeight*(32.0*TMath::Power(x,6.0)-48.0*TMath::Power(x, 4.0)+18.0*TMath::Power(x, 2.0)-1.0);
 558   Qvec[5][1] -= absWeight*(x*y*(32.0*TMath::Power(y,4.0)-32.0*TMath::Power(y, 2.0)+6.0));
 559
 560   return;
 561 }
 562
 563
 564 //_________________________________________________________________
 565 void AliReducedEvent::GetVZEROQvector(Double_t Qvec[][2], Int_t det) {
 566   //
 567   // Get the reaction plane from the VZERO detector for a given harmonic
 568   //
 569   GetVZEROQvector(Qvec, det, fVZEROMult);
 570 }
 571
 572
 573 //_________________________________________________________________
 574 void AliReducedEvent::GetVZEROQvector(Double_t Qvec[][2], Int_t det, Float_t* vzeroMult) {
 575   //
 576   // Get the reaction plane from the VZERO detector for a given harmonic
 577   //
 578   //  Q{x,y} = SUM_i mult(i) * {cos(n*phi_i), sin(n*phi_i)}
 579   //  phi_i - phi angle of the VZERO sector i
 580   //          Each sector covers 45 degrees(8 sectors per ring). Middle of sector 0 is at 45/2
 581   //        channel 0: 22.5
 582   //                1: 22.5+45
 583   //                2: 22.5+45*2
 584   //               ...
 585   //        at the next ring continues the same
 586   //        channel 8: 22.5
 587   //        channel 9: 22.5 + 45
 588   //
 589   if(!(det==AliReducedEventFriend::kVZEROA ||
 590        det==AliReducedEventFriend::kVZEROC))
 591     return;
 592
 593   const Double_t kX[8] = {0.92388, 0.38268, -0.38268, -0.92388, -0.92388, -0.38268, 0.38268, 0.92388};    // cosines of the angles of the VZERO sectors (8 per ring)
 594   const Double_t kY[8] = {0.38268, 0.92388, 0.92388, 0.38268, -0.38268, -0.92388, -0.92388, -0.38268};    // sines     -- " --
 595   Int_t phi;
 596
 597   for(Int_t iChannel=0; iChannel<64; ++iChannel) {
 598     if(iChannel<32 && det==AliReducedEventFriend::kVZEROA) continue;
 599     if(iChannel>=32 && det==AliReducedEventFriend::kVZEROC) continue;
 600     phi=iChannel%8;
 601     //  1st harmonic
 602     Qvec[0][0] += vzeroMult[iChannel]*kX[phi];
 603     Qvec[0][1] += vzeroMult[iChannel]*kY[phi];
 604     //  2nd harmonic
 605     Qvec[1][0] += vzeroMult[iChannel]*(2.0*TMath::Power(kX[phi],2.0)-1);
 606     Qvec[1][1] += vzeroMult[iChannel]*(2.0*kX[phi]*kY[phi]);
 607     //  3rd harmonic
 608     Qvec[2][0] += vzeroMult[iChannel]*(4.0*TMath::Power(kX[phi],3.0)-3.0*kX[phi]);
 609     Qvec[2][1] += vzeroMult[iChannel]*(3.0*kY[phi]-4.0*TMath::Power(kY[phi],3.0));
 610     //  4th harmonic
 611     Qvec[3][0] += vzeroMult[iChannel]*(1.0-8.0*TMath::Power(kX[phi]*kY[phi],2.0));
 612     Qvec[3][1] += vzeroMult[iChannel]*(4.0*kX[phi]*kY[phi]-8.0*kX[phi]*TMath::Power(kY[phi],3.0));
 613     //  5th harmonic
 614     Qvec[4][0] += vzeroMult[iChannel]*(16.0*TMath::Power(kX[phi],5.0)-20.0*TMath::Power(kX[phi], 3.0)+5.0*kX[phi]);
 615     Qvec[4][1] += vzeroMult[iChannel]*(16.0*TMath::Power(kY[phi],5.0)-20.0*TMath::Power(kY[phi], 3.0)+5.0*kY[phi]);
 616     //  6th harmonic
 617     Qvec[5][0] += vzeroMult[iChannel]*(32.0*TMath::Power(kX[phi],6.0)-48.0*TMath::Power(kX[phi], 4.0)+18.0*TMath::Power(kX[phi], 2.0)-1.0);
 618     Qvec[5][1] += vzeroMult[iChannel]*(kX[phi]*kY[phi]*(32.0*TMath::Power(kY[phi],4.0)-32.0*TMath::Power(kY[phi], 2.0)+6.0));
 619   }    // end loop over channels
 620 }
 621
 622
 623 //_________________________________________________________________
 624 void AliReducedEvent::GetZDCQvector(Double_t Qvec[][2], Int_t det) const {
 625   //
 626   // Construct the event plane using the ZDC
 627   // ZDC has 2 side (A and C) with 4 calorimeters on each side
 628   // The XY position of each calorimeter is specified by the
 629   // zdcNTowerCenters_x and zdcNTowerCenters_y arrays
 630   if(!(det==AliReducedEventFriend::kZDCA ||
 631        det==AliReducedEventFriend::kZDCC)) return;       // bad detector option
 632   const Float_t zdcTowerCenter = 1.75;
 633   const Float_t zdcNTowerCentersX[4] = {-zdcTowerCenter,  zdcTowerCenter, -zdcTowerCenter, zdcTowerCenter};
 634   const Float_t zdcNTowerCentersY[4] = {-zdcTowerCenter, -zdcTowerCenter,  zdcTowerCenter, zdcTowerCenter};
 635
 636   Qvec[0][0] = 0.0; Qvec[0][1] = 0.0;   // first harmonic Q-vector
 637   Float_t zdcNCentroidSum = 0;
 638   Float_t zdcNalpha = 0.5;
 639
 640   for(Int_t i=0; i<4; ++i) {
 641     if(fZDCnEnergy[i+(det==AliReducedEventFriend::kZDCA ? 4 : 0)]>0.0) {
 642       Float_t zdcNenergyAlpha = TMath::Power(fZDCnEnergy[i+(det==AliReducedEventFriend::kZDCA ? 4 : 0)], zdcNalpha);
 643       Qvec[0][0] += zdcNTowerCentersX[i-1]*zdcNenergyAlpha;
 644       Qvec[0][1] += zdcNTowerCentersY[i-1]*zdcNenergyAlpha;
 645       zdcNCentroidSum += zdcNenergyAlpha;
 646     }
 647   }   // end loop over channels
 648
 649   if(zdcNCentroidSum>0.0) {
 650     Qvec[0][0] /= zdcNCentroidSum;
 651     Qvec[0][1] /= zdcNCentroidSum;
 652   }
 653 }