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 #include <TClass.h>
  17
  18 ClassImp(AliReducedTrack)
  19 ClassImp(AliReducedPair)
  20 ClassImp(AliReducedFMD)
  21 ClassImp(AliReducedEvent)
  22 ClassImp(AliReducedEventFriend)
  23 ClassImp(AliReducedCaloCluster)
  24
  25 TClonesArray* AliReducedEvent::fgTracks = 0;
  26 TClonesArray* AliReducedEvent::fgCandidates = 0;
  27 TClonesArray* AliReducedEvent::fgCaloClusters = 0;
  28 TClonesArray* AliReducedEvent::fgFMD1 = 0;
  29 TClonesArray* AliReducedEvent::fgFMD2I = 0;
  30 TClonesArray* AliReducedEvent::fgFMD2O = 0;
  31 TClonesArray* AliReducedEvent::fgFMD3I = 0;
  32 TClonesArray* AliReducedEvent::fgFMD3O = 0;
  33
  34 //_______________________________________________________________________________
  35 AliReducedTrack::AliReducedTrack() :
  36   fTrackId(-1),
  37   fStatus(0),
  38   fGlobalPhi(0.0),
  39   fGlobalPt(0.0),
  40   fGlobalEta(0.0),
  41   fTPCPhi(0.0),
  42   fTPCPt(0.0),
  43   fTPCEta(0.0),
  44   fMomentumInner(0.0),
  45   fDCA(),
  46   fTrackLength(0.0),
  47   fITSclusterMap(0),
  48   fITSsignal(0.0),
  49   fITSnSig(),
  50   fITSchi2(0.0),
  51   fTPCNcls(0),
  52   fTPCCrossedRows(0),
  53   fTPCNclsF(0),
  54   fTPCNclsIter1(0),
  55   fTPCClusterMap(0),
  56   fTPCsignal(0),
  57   fTPCsignalN(0),
  58   fTPCnSig(),
  59   fTPCchi2(0.0),
  60   fTOFbeta(0.0),
  61   fTOFnSig(),
  62   fTOFdeltaBC(0),
  63   fTRDntracklets(),
  64   fTRDpid(),
  65   fTRDpidLQ2D(),
  66   fCaloClusterId(-999),
  67   fBayesPID(),
  68   fFlags(0),
  69   fMoreFlags(0)
  70 {
  71   //
  72   // Constructor
  73   //
  74   fDCA[0] = 0.0; fDCA[1]=0.0;
  75   for(Int_t i=0; i<4; ++i) {fTPCnSig[i]=-999.; fTOFnSig[i]=-999.; fITSnSig[i]=-999.;}
  76   for(Int_t i=0; i<3; ++i) {fBayesPID[i]=-999.;}
  77   fTRDpid[0]=-999.; fTRDpid[1]=-999.;
  78   fTRDpidLQ2D[0] = -999.; fTRDpidLQ2D[1] = -999.;
  79 }
  80
  81
  82 //_______________________________________________________________________________
  83 AliReducedTrack::~AliReducedTrack()
  84 {
  85   //
  86   // De-Constructor
  87   //
  88 }
  89
  90
  91 //_______________________________________________________________________________
  92 AliReducedPair::AliReducedPair() :
  93   fCandidateId(-1),
  94   fPairType(-1),
  95   fLegIds(),
  96   fMass(),
  97   fPhi(0.0),
  98   fPt(0.0),
  99   fEta(0.0),
 100   fLxy(0.0),
 101   fLxyErr(0.0),
 102   fPointingAngle(0.0),
 103   fChisquare(0.0),
 104   fMCid(0)
 105 {
 106   //
 107   // Constructor
 108   //
 109   fLegIds[0] = -1; fLegIds[1] = -1;
 110   fMass[0]=-999.; fMass[1]=-999.; fMass[2]=-999.; fMass[3]=-999.;
 111 }
 112
 113
 114 //_______________________________________________________________________________
 115 AliReducedPair::AliReducedPair(const AliReducedPair &c) :
 116   TObject(c),
 117   fCandidateId(c.CandidateId()),
 118   fPairType(c.PairType()),
 119   fLegIds(),
 120   fMass(),
 121   fPhi(c.Phi()),
 122   fPt(c.Pt()),
 123   fEta(c.Eta()),
 124   fLxy(c.Lxy()),
 125   fLxyErr(c.LxyErr()),
 126   fPointingAngle(c.PointingAngle()),
 127   fChisquare(c.Chi2()),
 128   fMCid(c.MCid())
 129 {
 130   //
 131   // copy constructor
 132   //
 133   fLegIds[0] = c.LegId(0);
 134   fLegIds[1] = c.LegId(1);
 135   fMass[0] = c.Mass(0); fMass[1] = c.Mass(1); fMass[2] = c.Mass(2); fMass[3] = c.Mass(3);
 136 }
 137
 138
 139 //_______________________________________________________________________________
 140 AliReducedPair::~AliReducedPair() {
 141   //
 142   // destructor
 143   //
 144 }
 145
 146
 147 //_______________________________________________________________________________
 148 AliReducedFMD::AliReducedFMD() :
 149   fMultiplicity(0),
 150   //fEta(0.0),
 151   fId(0)
 152 {
 153   AliReducedFMD::Class()->IgnoreTObjectStreamer();
 154   //
 155   // Constructor
 156   //
 157 }
 158
 159
 160 //_______________________________________________________________________________
 161 AliReducedFMD::~AliReducedFMD()
 162 {
 163   //
 164   // De-Constructor
 165   //
 166 }
 167
 168
 169 //____________________________________________________________________________
 170 AliReducedEvent::AliReducedEvent() :
 171   TObject(),
 172   fEventTag(0),
 173   fEventNumberInFile(0),
 174   fL0TriggerInputs(0),
 175   fL1TriggerInputs(0),
 176   fL2TriggerInputs(0),
 177   fRunNo(0),
 178   fBC(0),
 179   fTimeStamp(0),
 180   fEventType(0),
 181   fTriggerMask(0),
 182   fIsPhysicsSelection(kTRUE),
 183   fIsSPDPileup(kFALSE),
 184   fIsSPDPileupMultBins(kFALSE),
 185   fIRIntClosestIntMap(),
 186   fIsFMDReduced(kTRUE),
 187   fVtx(),
 188   fNVtxContributors(0),
 189   fVtxTPC(),
 190   fNVtxTPCContributors(0),
 191   fNpileupSPD(0),
 192   fNpileupTracks(0),
 193   fNPMDtracks(0),
 194   fNTRDtracks(0),
 195   fNTRDtracklets(0),
 196   fCentrality(),
 197   fCentQuality(0),
 198   fNV0candidates(),
 199   fNDielectronCandidates(0),
 200   fNtracks(),
 201   fSPDntracklets(0),
 202   fVZEROMult(),
 203   fZDCnEnergy(),
 204   fZDCpEnergy(),
 205   fT0amplitude(),
 206   fT0TOF(),
 207   fT0TOFbest(),
 208   fT0zVertex(0),
 209   fT0start(0),
 210   fT0pileup(kFALSE),
 211   fT0sattelite(kFALSE),
 212   fTracks(0x0),
 213   fCandidates(0x0),
 214   fNFMDchannels(),
 215   fFMDtotalMult(),
 216   fFMD1(0x0),
 217   fFMD2I(0x0),
 218   fFMD2O(0x0),
 219   fFMD3I(0x0),
 220   fFMD3O(0x0),
 221   fNCaloClusters(0),
 222   fCaloClusters(0x0)
 223 {
 224   //
 225   // Constructor
 226   //
 227   for(Int_t i=0; i<2; ++i) fIRIntClosestIntMap[i] = 0;
 228   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.;}
 229   for(Int_t i=0; i<4; ++i) fCentrality[i]=-1.;
 230   fNV0candidates[0]=0; fNV0candidates[1]=0;
 231   fNtracks[0]=0; fNtracks[1]=0;
 232   for(Int_t i=0; i<32; ++i) fSPDntrackletsEta[i]=0;
 233   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i]=0;
 234   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 235   for(Int_t i=0; i<10; ++i) fZDCnEnergy[i]=0.0;
 236   for(Int_t i=0; i<10; ++i) fZDCpEnergy[i]=0.0;
 237   for(Int_t i=0; i<26; ++i) fT0amplitude[i]=0.0;
 238   for(Int_t i=0; i<3; ++i)  fT0TOF[i]=0.0;
 239   for(Int_t i=0; i<3; ++i)  fT0TOFbest[i]=0.0;
 240   for(Int_t i=0; i<5; ++i)  fNFMDchannels[i]=0.0;
 241   for(Int_t i=0; i<5; ++i)  fFMDtotalMult[i]=0.0;
 242 }
 243
 244
 245 //____________________________________________________________________________
 246 AliReducedEvent::AliReducedEvent(const Char_t* /*name*/) :
 247   TObject(),
 248   fEventTag(0),
 249   fEventNumberInFile(0),
 250   fL0TriggerInputs(0),
 251   fL1TriggerInputs(0),
 252   fL2TriggerInputs(0),
 253   fRunNo(0),
 254   fBC(0),
 255   fTimeStamp(0),
 256   fEventType(0),
 257   fTriggerMask(0),
 258   fIsPhysicsSelection(kTRUE),
 259   fIsSPDPileup(kFALSE),
 260   fIsSPDPileupMultBins(kFALSE),
 261   fIRIntClosestIntMap(),
 262   fIsFMDReduced(kTRUE),
 263   fVtx(),
 264   fNVtxContributors(0),
 265   fVtxTPC(),
 266   fNVtxTPCContributors(0),
 267   fNpileupSPD(0),
 268   fNpileupTracks(0),
 269   fNPMDtracks(0),
 270   fNTRDtracks(0),
 271   fNTRDtracklets(0),
 272   fCentrality(),
 273   fCentQuality(0),
 274   fNV0candidates(),
 275   fNDielectronCandidates(0),
 276   fNtracks(),
 277   fSPDntracklets(0),
 278   fVZEROMult(),
 279   fZDCnEnergy(),
 280   fZDCpEnergy(),
 281   fT0amplitude(),
 282   fT0TOF(),
 283   fT0TOFbest(),
 284   fT0zVertex(0),
 285   fT0start(0),
 286   fT0pileup(kFALSE),
 287   fT0sattelite(kFALSE),
 288   fTracks(0x0),
 289   fCandidates(0x0),
 290   fNFMDchannels(),
 291   fFMDtotalMult(),
 292   fFMD1(0x0),
 293   fFMD2I(0x0),
 294   fFMD2O(0x0),
 295   fFMD3I(0x0),
 296   fFMD3O(0x0),
 297   fNCaloClusters(0),
 298   fCaloClusters(0x0)
 299 //fFMDMult()
 300 {
 301   //
 302   // Constructor
 303   //
 304   for(Int_t i=0; i<2; ++i) fIRIntClosestIntMap[i] = 0;
 305   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.;}
 306   for(Int_t i=0; i<4; ++i) fCentrality[i]=-1.;
 307   fNV0candidates[0]=0; fNV0candidates[1]=0;
 308   fNtracks[0]=0; fNtracks[1]=0;
 309   for(Int_t i=0; i<32; ++i) fSPDntrackletsEta[i]=0;
 310   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i]=0;
 311   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 312   for(Int_t i=0; i<10; ++i) fZDCnEnergy[i]=0.0;
 313   for(Int_t i=0; i<10; ++i) fZDCpEnergy[i]=0.0;
 314   for(Int_t i=0; i<26; ++i) fT0amplitude[i]=0.0;
 315   for(Int_t i=0; i<3; ++i)  fT0TOF[i]=0.0;
 316   for(Int_t i=0; i<3; ++i)  fT0TOFbest[i]=0.0;
 317   for(Int_t i=0; i<5; ++i)  fNFMDchannels[i]=0.0;
 318   for(Int_t i=0; i<5; ++i)  fFMDtotalMult[i]=0.0;
 319
 320   if(!fgTracks) fgTracks = new TClonesArray("AliReducedTrack", 100000);
 321   fTracks = fgTracks;
 322   if(!fgCandidates) fgCandidates = new TClonesArray("AliReducedPair", 100000);
 323   fCandidates = fgCandidates;
 324   if(!fgCaloClusters) fgCaloClusters = new TClonesArray("AliReducedCaloCluster", 50000);
 325   fCaloClusters = fgCaloClusters;
 326   if(!fgFMD1) fgFMD1 = new TClonesArray("AliReducedFMD", 10240);
 327   fFMD1 = fgFMD1;
 328   if(!fgFMD2I) fgFMD2I = new TClonesArray("AliReducedFMD", 10240);
 329   fFMD2I = fgFMD2I;
 330   if(!fgFMD2O) fgFMD2O = new TClonesArray("AliReducedFMD", 10240);
 331   fFMD2O = fgFMD2O;
 332   if(!fgFMD3I) fgFMD3I = new TClonesArray("AliReducedFMD", 10240);
 333   fFMD3I = fgFMD3I;
 334   if(!fgFMD3O) fgFMD3O = new TClonesArray("AliReducedFMD", 10240);
 335   fFMD3O = fgFMD3O;
 336 }
 337
 338
 339 //____________________________________________________________________________
 340 AliReducedEvent::~AliReducedEvent()
 341 {
 342   //
 343   // De-Constructor
 344   //
 345   //ClearEvent();
 346 }
 347
 348
 349 //____________________________________________________________________________
 350 Float_t AliReducedEvent::MultVZEROA() const
 351 {
 352   //
 353   // Total VZERO multiplicity in A side
 354   //
 355   Float_t mult=0.0;
 356   for(Int_t i=32;i<64;++i)
 357     mult += fVZEROMult[i];
 358   return mult;
 359 }
 360
 361
 362 //____________________________________________________________________________
 363 Float_t AliReducedEvent::MultVZEROC() const
 364 {
 365   //
 366   // Total VZERO multiplicity in C side
 367   //
 368   Float_t mult=0.0;
 369   for(Int_t i=0;i<32;++i)
 370     mult += fVZEROMult[i];
 371   return mult;
 372 }
 373
 374
 375 //____________________________________________________________________________
 376 Float_t AliReducedEvent::MultVZERO() const
 377 {
 378   //
 379   // Total VZERO multiplicity
 380   //
 381   return MultVZEROA()+MultVZEROC();
 382 }
 383
 384
 385 //____________________________________________________________________________
 386 Float_t AliReducedEvent::MultRingVZEROA(Int_t ring) const
 387 {
 388   //
 389   //  VZERO multiplicity in a given ring on A side
 390   //
 391   if(ring<0 || ring>3) return -1.0;
 392
 393   Float_t mult=0.0;
 394   for(Int_t i=32+8*ring; i<32+8*(ring+1); ++i)
 395     mult += fVZEROMult[i];
 396   return mult;
 397 }
 398
 399
 400 //____________________________________________________________________________
 401 Float_t AliReducedEvent::MultRingVZEROC(Int_t ring) const
 402 {
 403   //
 404   //  VZERO multiplicity in a given ring on C side
 405   //
 406   if(ring<0 || ring>3) return -1.0;
 407
 408   Float_t mult=0.0;
 409   for(Int_t i=8*ring; i<8*(ring+1); ++i)
 410     mult += fVZEROMult[i];
 411   return mult;
 412 }
 413
 414 //____________________________________________________________________________
 415 Float_t AliReducedEvent::AmplitudeTZEROA() const
 416 {
 417   //
 418   // Total TZERO multiplicity in A side
 419   //
 420   Float_t mult=0.0;
 421   for(Int_t i=12;i<24;++i)
 422     mult += fT0amplitude[i];
 423   return mult;
 424 }
 425
 426
 427 //____________________________________________________________________________
 428 Float_t AliReducedEvent::AmplitudeTZEROC() const
 429 {
 430   //
 431   // Total TZERO multiplicity in C side
 432   //
 433   Float_t mult=0.0;
 434   for(Int_t i=0;i<12;++i)
 435     mult += fT0amplitude[i];
 436   return mult;
 437 }
 438
 439
 440 //____________________________________________________________________________
 441 Float_t AliReducedEvent::AmplitudeTZERO() const
 442 {
 443   //
 444   // Total TZERO multiplicity
 445   //
 446   return AmplitudeTZEROA()+AmplitudeTZEROC();
 447 }
 448
 449
 450
 451 //____________________________________________________________________________
 452 Float_t AliReducedEvent::EnergyZDCA() const
 453 {
 454   //
 455   // Total ZDC energy in A side
 456   //
 457   Float_t energy=0.0;
 458   for(Int_t i=6;i<10;++i){
 459       if(fZDCnEnergy[i]>0.) {
 460   Float_t zdcNenergyAlpha = TMath::Power(fZDCnEnergy[i], gZdcNalpha);
 461     energy += zdcNenergyAlpha;}
 462   }
 463   return energy;
 464 }
 465
 466
 467 //____________________________________________________________________________
 468 Float_t AliReducedEvent::EnergyZDCC() const
 469 {
 470   //
 471   // Total ZDC energy in C side
 472   //
 473   Float_t energy=0.0;
 474   for(Int_t i=1;i<5;++i){
 475       if(fZDCnEnergy[i]>0.) {
 476     Float_t zdcNenergyAlpha = TMath::Power(fZDCnEnergy[i], gZdcNalpha);
 477     energy += zdcNenergyAlpha;}
 478   }
 479   return energy;
 480 }
 481
 482
 483 //____________________________________________________________________________
 484 Float_t AliReducedEvent::EnergyZDC() const
 485 {
 486   //
 487   // Total ZDC energy
 488   //
 489   return EnergyZDCA()+EnergyZDCC();
 490 }
 491
 492
 493 //____________________________________________________________________________
 494 Float_t AliReducedEvent::EnergyZDCn(Int_t ch) const
 495 {
 496   //
 497   // Total ZDC energy in channel
 498   //
 499   Float_t energy=0;
 500   if(ch<0 || ch>9) return -999.;
 501       if(fZDCnEnergy[ch]>0.) {
 502                         energy = TMath::Power(fZDCnEnergy[ch], gZdcNalpha);
 503                 }
 504   return energy;
 505
 506 }
 507
 508
 509 //_____________________________________________________________________________
 510 void AliReducedEvent::ClearEvent() {
 511   //
 512   // clear the event
 513   //
 514   if(fTracks) fTracks->Clear("C");
 515   if(fCandidates) fCandidates->Clear("C");
 516   if(fCaloClusters) fCaloClusters->Clear("C");
 517   if(fFMD1) fFMD1->Clear("C");
 518   if(fFMD2I) fFMD2I->Clear("C");
 519   if(fFMD2O) fFMD2O->Clear("C");
 520   if(fFMD3I) fFMD3I->Clear("C");
 521   if(fFMD3O) fFMD3O->Clear("C");
 522   fEventTag = 0;
 523   fEventNumberInFile = -999;
 524   fL0TriggerInputs=0;
 525   fL1TriggerInputs=0;
 526   fL2TriggerInputs=0;
 527   fIRIntClosestIntMap[0] = 0; fIRIntClosestIntMap[1] = 0;
 528   fRunNo = 0;
 529   fBC = 0;
 530   fTimeStamp = 0;
 531   fEventType = 0;
 532   fTriggerMask = 0;
 533   fIsPhysicsSelection = kTRUE;
 534   fIsSPDPileup = kFALSE;
 535   fIsSPDPileupMultBins = kFALSE;
 536   fIsFMDReduced = kTRUE;
 537   fNVtxContributors = 0;
 538   fNVtxTPCContributors = 0;
 539   fCentQuality = 0;
 540   for(Int_t i=0;i<4;++i) fCentrality[i] = -1.0;
 541   fNV0candidates[0] = 0; fNV0candidates[1] = 0;
 542   fNpileupSPD=0;
 543   fNpileupTracks=0;
 544   fNPMDtracks=0;
 545   fNTRDtracks=0;
 546   fNTRDtracklets=0;
 547   fNDielectronCandidates = 0;
 548   fNtracks[0] = 0; fNtracks[1] = 0;
 549   for(Int_t i=0; i<32; ++i) fSPDntrackletsEta[i] = 0;
 550   for(Int_t i=0; i<32; ++i) fNtracksPerTrackingFlag[i] = 0;
 551   for(Int_t i=0; i<3; ++i) {fVtx[i]=-999.; fVtxTPC[i]=-999.;}
 552   for(Int_t i=0; i<64; ++i) fVZEROMult[i] = 0.0;
 553   for(Int_t i=0; i<10; ++i) fZDCnEnergy[i]=0.0;
 554   for(Int_t i=0; i<10; ++i) fZDCpEnergy[i]=0.0;
 555   for(Int_t i=0; i<26; ++i) fT0amplitude[i]=0.0;
 556   for(Int_t i=0; i<3; ++i)  fT0TOF[i]=0.0;
 557   for(Int_t i=0; i<3; ++i)  fT0TOFbest[i]=0.0;
 558   fT0pileup = kFALSE;
 559   fT0zVertex = -999.;
 560   fT0start = -999.;
 561   fT0sattelite = kFALSE;
 562   for(Int_t i=0; i<5; ++i)  fNFMDchannels[i]=0.0;
 563   for(Int_t i=0; i<5; ++i)  fFMDtotalMult[i]=0.0;
 564 }
 565
 566
 567 //_______________________________________________________________________________
 568 AliReducedEventFriend::AliReducedEventFriend() :
 569  fQvector(),
 570  fEventPlaneStatus()
 571 {
 572   //
 573   // default constructor
 574   //
 575   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 576     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 577       fEventPlaneStatus[idet][ih] = 0;
 578       for(Int_t ic=0; ic<2; ++ic)
 579         fQvector[idet][ih][ic] = 0.0;
 580     }
 581   }
 582 }
 583
 584
 585 //____________________________________________________________________________
 586 AliReducedEventFriend::~AliReducedEventFriend()
 587 {
 588   //
 589   // De-Constructor
 590   //
 591   ClearEvent();
 592 }
 593
 594
 595 //_____________________________________________________________________________
 596 void AliReducedEventFriend::ClearEvent() {
 597   //
 598   // clear the event
 599   //
 600   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 601     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 602       fEventPlaneStatus[idet][ih] = 0;
 603       for(Int_t ic=0; ic<2; ++ic)
 604         fQvector[idet][ih][ic] = 0.0;
 605     }
 606   }
 607 }
 608
 609
 610 //____________________________________________________________________________
 611 void AliReducedEventFriend::CopyEvent(const AliReducedEventFriend* event) {
 612   //
 613   // copy information from another event to this one
 614   //
 615   for(Int_t idet=0; idet<kNdetectors; ++idet) {
 616     for(Int_t ih=0; ih<fgkNMaxHarmonics; ++ih) {
 617       fQvector[idet][ih][0] = event->Qx(idet, ih+1);
 618       fQvector[idet][ih][1] = event->Qy(idet, ih+1);
 619       fEventPlaneStatus[idet][ih] = event->GetEventPlaneStatus(idet, ih+1);
 620     }
 621   }
 622 }
 623
 624
 625 //_____________________________________________________________________________
 626 AliReducedCaloCluster::AliReducedCaloCluster() :
 627  fType(kUndefined),
 628  fEnergy(-999.),
 629  fTrackDx(-999.),
 630  fTrackDz(-999.),
 631  fM20(-999.),
 632  fM02(-999.),
 633  fDispersion(-999.)
 634 {
 635   //
 636   // default constructor
 637   //
 638 }
 639
 640
 641 //_____________________________________________________________________________
 642 AliReducedCaloCluster::~AliReducedCaloCluster()
 643 {
 644   //
 645   // destructor
 646   //
 647 }
 648
 649
 650 //_______________________________________________________________________________
 651 void AliReducedEvent::GetQvector(Double_t Qvec[][2], Int_t det,
 652                                  Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 653                                  Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 654   //
 655   // Get the event plane for a specified detector
 656   //
 657   if(det==AliReducedEventFriend::kTPC ||
 658      det==AliReducedEventFriend::kTPCptWeights ||
 659      det==AliReducedEventFriend::kTPCpos ||
 660      det==AliReducedEventFriend::kTPCneg) {
 661     GetTPCQvector(Qvec, det, etaMin, etaMax, IsTrackSelected);
 662     return;
 663   }
 664   if(det==AliReducedEventFriend::kVZEROA ||
 665      det==AliReducedEventFriend::kVZEROC) {
 666     GetVZEROQvector(Qvec, det);
 667     return;
 668   }
 669   if(det==AliReducedEventFriend::kZDCA ||
 670      det==AliReducedEventFriend::kZDCC) {
 671     GetZDCQvector(Qvec, det);
 672     return;
 673   }
 674   if(det==AliReducedEventFriend::kFMD) {
 675     //TODO implementation
 676     return;
 677   }
 678   return;
 679 }
 680
 681
 682 //_________________________________________________________________
 683 Int_t AliReducedEvent::GetTPCQvector(Double_t Qvec[][2], Int_t det,
 684                                      Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 685                                      Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 686   //
 687   // Construct the event plane using tracks in the barrel
 688   //
 689   if(!(det==AliReducedEventFriend::kTPC ||
 690        det==AliReducedEventFriend::kTPCpos ||
 691        det==AliReducedEventFriend::kTPCneg))
 692     return 0;
 693   Int_t nUsedTracks = 0;
 694   Short_t charge = 0;
 695   AliReducedTrack* track=0x0;
 696   Double_t weight=0.0; Double_t absWeight = 0.0; Double_t x=0.0; Double_t y=0.0;
 697   TIter nextTrack(fTracks);
 698   while((track=static_cast<AliReducedTrack*>(nextTrack()))) {
 699     if(track->Eta()<etaMin) continue;
 700     if(track->Eta()>etaMax) continue;
 701     charge = track->Charge();
 702     if(det==AliReducedEventFriend::kTPCpos && charge<0) continue;
 703     if(det==AliReducedEventFriend::kTPCneg && charge>0) continue;
 704
 705     if(IsTrackSelected && !IsTrackSelected(track)) continue;
 706     absWeight = 1.0;
 707     if(det==AliReducedEventFriend::kTPCptWeights) {
 708       absWeight = track->Pt();
 709       if(absWeight>2.0) absWeight = 2.0;    // pt is the weight used for the event plane
 710     }
 711     weight = absWeight;
 712
 713     ++nUsedTracks;
 714     x = TMath::Cos(track->Phi());
 715     y = TMath::Sin(track->Phi());
 716     //  1st harmonic
 717     Qvec[0][0] += weight*x;
 718     Qvec[0][1] += weight*y;
 719     //  2nd harmonic
 720     Qvec[1][0] += absWeight*(2.0*TMath::Power(x,2.0)-1);
 721     Qvec[1][1] += absWeight*(2.0*x*y);
 722     //  3rd harmonic
 723     Qvec[2][0] += weight*(4.0*TMath::Power(x,3.0)-3.0*x);
 724     Qvec[2][1] += weight*(3.0*y-4.0*TMath::Power(y,3.0));
 725     //  4th harmonic
 726     Qvec[3][0] += absWeight*(1.0-8.0*TMath::Power(x*y,2.0));
 727     Qvec[3][1] += absWeight*(4.0*x*y-8.0*x*TMath::Power(y,3.0));
 728     //  5th harmonic
 729     Qvec[4][0] += weight*(16.0*TMath::Power(x,5.0)-20.0*TMath::Power(x, 3.0)+5.0*x);
 730     Qvec[4][1] += weight*(16.0*TMath::Power(y,5.0)-20.0*TMath::Power(y, 3.0)+5.0*y);
 731     //  6th harmonic
 732     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);
 733     Qvec[5][1] += absWeight*(x*y*(32.0*TMath::Power(y,4.0)-32.0*TMath::Power(y, 2.0)+6.0));
 734   }
 735   return nUsedTracks;
 736 }
 737
 738
 739 //____________________________________________________________________________
 740 void AliReducedEvent::SubtractParticleFromQvector(
 741         AliReducedTrack* particle, Double_t Qvec[][2], Int_t det,
 742         Float_t etaMin/*=-0.8*/, Float_t etaMax/*=+0.8*/,
 743         Bool_t (*IsTrackSelected)(AliReducedTrack*)/*=NULL*/) {
 744   //
 745   // subtract a particle from the event Q-vector
 746   //
 747   Float_t eta = particle->Eta();
 748   if(eta<etaMin) return;
 749   if(eta>etaMax) return;
 750
 751   Float_t charge = particle->Charge();
 752   if(det==AliReducedEventFriend::kTPCpos && charge<0) return;
 753   if(det==AliReducedEventFriend::kTPCneg && charge>0) return;
 754
 755   if(IsTrackSelected && !IsTrackSelected(particle)) return;
 756
 757   Double_t weight=0.0; Double_t absWeight = 0.0;
 758   if(det==AliReducedEventFriend::kTPCptWeights) {
 759     absWeight = particle->Pt();
 760     if(absWeight>2.0) absWeight = 2.0;
 761   }
 762   weight = absWeight;
 763   //  if(eta<0.0) weight *= -1.0;
 764
 765   Float_t x = TMath::Cos(particle->Phi());
 766   Float_t y = TMath::Sin(particle->Phi());
 767
 768   //  1st harmonic
 769   Qvec[0][0] -= weight*x;
 770   Qvec[0][1] -= weight*y;
 771   //  2nd harmonic
 772   Qvec[1][0] -= absWeight*(2.0*TMath::Power(x,2.0)-1);
 773   Qvec[1][1] -= absWeight*(2.0*x*y);
 774   //  3rd harmonic
 775   Qvec[2][0] -= weight*(4.0*TMath::Power(x,3.0)-3.0*x);
 776   Qvec[2][1] -= weight*(3.0*y-4.0*TMath::Power(y,3.0));
 777   //  4th harmonic
 778   Qvec[3][0] -= absWeight*(1.0-8.0*TMath::Power(x*y,2.0));
 779   Qvec[3][1] -= absWeight*(4.0*x*y-8.0*x*TMath::Power(y,3.0));
 780   //  5th harmonic
 781   Qvec[4][0] -= weight*(16.0*TMath::Power(x,5.0)-20.0*TMath::Power(x, 3.0)+5.0*x);
 782   Qvec[4][1] -= weight*(16.0*TMath::Power(y,5.0)-20.0*TMath::Power(y, 3.0)+5.0*y);
 783   //  6th harmonic
 784   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);
 785   Qvec[5][1] -= absWeight*(x*y*(32.0*TMath::Power(y,4.0)-32.0*TMath::Power(y, 2.0)+6.0));
 786
 787   return;
 788 }
 789
 790
 791 //_________________________________________________________________
 792 void AliReducedEvent::GetVZEROQvector(Double_t Qvec[][2], Int_t det) {
 793   //
 794   // Get the reaction plane from the VZERO detector for a given harmonic
 795   //
 796   GetVZEROQvector(Qvec, det, fVZEROMult);
 797 }
 798
 799
 800 //_________________________________________________________________
 801 void AliReducedEvent::GetVZEROQvector(Double_t Qvec[][2], Int_t det, Float_t* vzeroMult) {
 802   //
 803   // Get the reaction plane from the VZERO detector for a given harmonic
 804   //
 805   //  Q{x,y} = SUM_i mult(i) * {cos(n*phi_i), sin(n*phi_i)}
 806   //  phi_i - phi angle of the VZERO sector i
 807   //          Each sector covers 45 degrees(8 sectors per ring). Middle of sector 0 is at 45/2
 808   //        channel 0: 22.5
 809   //                1: 22.5+45
 810   //                2: 22.5+45*2
 811   //               ...
 812   //        at the next ring continues the same
 813   //        channel 8: 22.5
 814   //        channel 9: 22.5 + 45
 815   //
 816   if(!(det==AliReducedEventFriend::kVZEROA ||
 817        det==AliReducedEventFriend::kVZEROC))
 818     return;
 819
 820   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)
 821   const Double_t kY[8] = {0.38268, 0.92388, 0.92388, 0.38268, -0.38268, -0.92388, -0.92388, -0.38268};    // sines     -- " --
 822   Int_t phi;
 823
 824   for(Int_t iChannel=0; iChannel<64; ++iChannel) {
 825     if(iChannel<32 && det==AliReducedEventFriend::kVZEROA) continue;
 826     if(iChannel>=32 && det==AliReducedEventFriend::kVZEROC) continue;
 827     phi=iChannel%8;
 828     //  1st harmonic
 829     Qvec[0][0] += vzeroMult[iChannel]*kX[phi];
 830     Qvec[0][1] += vzeroMult[iChannel]*kY[phi];
 831     //  2nd harmonic
 832     Qvec[1][0] += vzeroMult[iChannel]*(2.0*TMath::Power(kX[phi],2.0)-1);
 833     Qvec[1][1] += vzeroMult[iChannel]*(2.0*kX[phi]*kY[phi]);
 834     //  3rd harmonic
 835     Qvec[2][0] += vzeroMult[iChannel]*(4.0*TMath::Power(kX[phi],3.0)-3.0*kX[phi]);
 836     Qvec[2][1] += vzeroMult[iChannel]*(3.0*kY[phi]-4.0*TMath::Power(kY[phi],3.0));
 837     //  4th harmonic
 838     Qvec[3][0] += vzeroMult[iChannel]*(1.0-8.0*TMath::Power(kX[phi]*kY[phi],2.0));
 839     Qvec[3][1] += vzeroMult[iChannel]*(4.0*kX[phi]*kY[phi]-8.0*kX[phi]*TMath::Power(kY[phi],3.0));
 840     //  5th harmonic
 841     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]);
 842     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]);
 843     //  6th harmonic
 844     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);
 845     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));
 846   }    // end loop over channels
 847 }
 848
 849
 850 //_________________________________________________________________
 851 void AliReducedEvent::GetZDCQvector(Double_t Qvec[][2], Int_t det) const {
 852   //
 853   // Get the reaction plane from the ZDC detector for a given harmonic
 854   //
 855   GetZDCQvector(Qvec, det, fZDCnEnergy);
 856 }
 857
 858
 859 //_________________________________________________________________
 860 void AliReducedEvent::GetZDCQvector(Double_t Qvec[][2], Int_t det, const Float_t* zdcEnergy) const {
 861   //
 862   // Construct the event plane using the ZDC
 863   // ZDC has 2 side (A and C) with 4 calorimeters on each side
 864   // The XY position of each calorimeter is specified by the
 865   // zdcNTowerCenters_x and zdcNTowerCenters_y arrays
 866   if(!(det==AliReducedEventFriend::kZDCA ||
 867        det==AliReducedEventFriend::kZDCC ))
 868     return;
 869
 870
 871   //Int_t minTower,maxTower;
 872   //Float_t totalEnergy = 0.0;
 873
 874   const Float_t zdcTowerCenter = 1.75;
 875
 876
 877   Float_t zdcNTowerCentersX[4] = {0.0};
 878   Float_t zdcNTowerCentersY[4] = {0.0};
 879
 880
 881   if(det==AliReducedEventFriend::kZDCA){
 882     zdcNTowerCentersX[0] =  zdcTowerCenter; zdcNTowerCentersX[1] =-zdcTowerCenter;
 883     zdcNTowerCentersX[2] =  zdcTowerCenter; zdcNTowerCentersX[3] =-zdcTowerCenter;
 884   }
 885
 886   if(det==AliReducedEventFriend::kZDCC){
 887     zdcNTowerCentersX[0] = -zdcTowerCenter; zdcNTowerCentersX[1] = zdcTowerCenter;
 888     zdcNTowerCentersX[2] = -zdcTowerCenter; zdcNTowerCentersX[3] = zdcTowerCenter;
 889   }
 890
 891   zdcNTowerCentersY[0] = -zdcTowerCenter; zdcNTowerCentersY[1] =-zdcTowerCenter;
 892   zdcNTowerCentersY[2] =  zdcTowerCenter; zdcNTowerCentersY[3] = zdcTowerCenter;
 893
 894
 895   for(Int_t ih=0; ih<6; ih++) {Qvec[ih][0] = 0.0; Qvec[ih][1] = 0.0;}   // harmonic Q-vector
 896   Float_t zdcNCentroidSum = 0;
 897
 898   for(Int_t iTow=(det==AliReducedEventFriend::kZDCA ? 6 : 1); iTow<(det==AliReducedEventFriend::kZDCA ? 10 : 5); ++iTow) {
 899     //if(fZDCnEnergy[i+(det==AliReducedEventFriend::kZDCA ? 4 : 0)]>0.0) {
 900     //  1st harmonic
 901     Qvec[0][0] += zdcEnergy[iTow]*zdcNTowerCentersX[(iTow-1)%5];
 902     Qvec[0][1] += zdcEnergy[iTow]*zdcNTowerCentersY[(iTow-1)%5];
 903     //  2nd harmonic
 904     Qvec[1][0] += zdcEnergy[iTow]*(2.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5],2.0)-1);
 905     Qvec[1][1] += zdcEnergy[iTow]*(2.0*zdcNTowerCentersX[(iTow-1)%5]*zdcNTowerCentersY[(iTow-1)%5]);
 906     //  3rd harmonic
 907     Qvec[2][0] += zdcEnergy[iTow]*(4.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5],3.0)-3.0*zdcNTowerCentersX[(iTow-1)%5]);
 908     Qvec[2][1] += zdcEnergy[iTow]*(3.0*zdcNTowerCentersY[(iTow-1)%5]-4.0*TMath::Power(zdcNTowerCentersY[(iTow-1)%5],3.0));
 909     //  4th harmonic
 910     Qvec[3][0] += zdcEnergy[iTow]*(1.0-8.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5]*zdcNTowerCentersY[(iTow-1)%5],2.0));
 911     Qvec[3][1] += zdcEnergy[iTow]*(4.0*zdcNTowerCentersX[(iTow-1)%5]*zdcNTowerCentersY[(iTow-1)%5]-8.0*zdcNTowerCentersX[(iTow-1)%5]*TMath::Power(zdcNTowerCentersY[(iTow-1)%5],3.0));
 912     //  5th harmonic
 913     Qvec[4][0] += zdcEnergy[iTow]*(16.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5],5.0)-20.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5], 3.0)+5.0*zdcNTowerCentersX[(iTow-1)%5]);
 914     Qvec[4][1] += zdcEnergy[iTow]*(16.0*TMath::Power(zdcNTowerCentersY[(iTow-1)%5],5.0)-20.0*TMath::Power(zdcNTowerCentersY[(iTow-1)%5], 3.0)+5.0*zdcNTowerCentersY[(iTow-1)%5]);
 915     //  6th harmonic
 916     Qvec[5][0] += zdcEnergy[iTow]*(32.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5],6.0)-48.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5], 4.0)+18.0*TMath::Power(zdcNTowerCentersX[(iTow-1)%5], 2.0)-1.0);
 917     Qvec[5][1] += zdcEnergy[iTow]*(zdcNTowerCentersX[(iTow-1)%5]*zdcNTowerCentersY[(iTow-1)%5]*(32.0*TMath::Power(zdcNTowerCentersY[(iTow-1)%5],4.0)-32.0*TMath::Power(zdcNTowerCentersY[(iTow-1)%5], 2.0)+6.0));
 918
 919     zdcNCentroidSum += zdcEnergy[iTow];
 920
 921     if(zdcNCentroidSum>0.0) {
 922       Qvec[0][0] /= zdcNCentroidSum;
 923       Qvec[0][1] /= zdcNCentroidSum;
 924       Qvec[1][0] /= zdcNCentroidSum;
 925       Qvec[1][1] /= zdcNCentroidSum;
 926       Qvec[2][0] /= zdcNCentroidSum;
 927       Qvec[2][1] /= zdcNCentroidSum;
 928       Qvec[3][0] /= zdcNCentroidSum;
 929       Qvec[3][1] /= zdcNCentroidSum;
 930       Qvec[4][0] /= zdcNCentroidSum;
 931       Qvec[4][1] /= zdcNCentroidSum;
 932       Qvec[5][0] /= zdcNCentroidSum;
 933       Qvec[5][1] /= zdcNCentroidSum;
 934     }
 935   }   // end loop over channels
 936   return;
 937 }