PWG2/FLOW/AliFlowCommon/AliFlowEventSimple.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, 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   AliFlowEventSimple: A simple event
  18   for flow analysis
  19
  20   origin: Naomi van der Kolk (kolk@nikhef.nl)
  21           Ante Bilandzic     (anteb@nikhef.nl)
  22           Raimond Snellings  (Raimond.Snellings@nikhef.nl)
  23   mods:   Mikolaj Krzewicki  (mikolaj.krzewicki@cern.ch)
  24 *****************************************************************/
  25
  26 #include "Riostream.h"
  27 #include "TObjArray.h"
  28 #include "TFile.h"
  29 #include "TList.h"
  30 #include "TTree.h"
  31 #include "TParticle.h"
  32 #include "TMath.h"
  33 #include "TH1F.h"
  34 #include "TH1D.h"
  35 #include "TF1.h"
  36 #include "TProfile.h"
  37 #include "TParameter.h"
  38 #include "TBrowser.h"
  39 #include "AliFlowVector.h"
  40 #include "AliFlowTrackSimple.h"
  41 #include "AliFlowTrackSimpleCuts.h"
  42 #include "AliFlowEventSimple.h"
  43 #include "TRandom.h"
  44
  45 ClassImp(AliFlowEventSimple)
  46
  47 //-----------------------------------------------------------------------
  48 AliFlowEventSimple::AliFlowEventSimple():
  49   fTrackCollection(NULL),
  50   fReferenceMultiplicity(0),
  51   fNumberOfTracks(0),
  52   fNumberOfRPs(0),
  53   fMCReactionPlaneAngle(0.),
  54   fMCReactionPlaneAngleIsSet(kFALSE),
  55   fAfterBurnerPrecision(0.001),
  56   fNumberOfTracksWrap(NULL),
  57   fNumberOfRPsWrap(NULL),
  58   fMCReactionPlaneAngleWrap(NULL)
  59 {
  60   cout << "AliFlowEventSimple: Default constructor to be used only by root for io" << endl;
  61 }
  62
  63 //-----------------------------------------------------------------------
  64 AliFlowEventSimple::AliFlowEventSimple(Int_t aLength):
  65   fTrackCollection(new TObjArray(aLength)),
  66   fReferenceMultiplicity(0),
  67   fNumberOfTracks(0),
  68   fNumberOfRPs(0),
  69   fMCReactionPlaneAngle(0.),
  70   fMCReactionPlaneAngleIsSet(kFALSE),
  71   fAfterBurnerPrecision(0.001),
  72   fNumberOfTracksWrap(NULL),
  73   fNumberOfRPsWrap(NULL),
  74   fMCReactionPlaneAngleWrap(NULL)
  75 {
  76   //constructor
  77 }
  78
  79 //-----------------------------------------------------------------------
  80 AliFlowEventSimple::AliFlowEventSimple( Int_t nParticles,
  81                                         TF1* ptDist,
  82                                         Double_t phiMin,
  83                                         Double_t phiMax,
  84                                         Double_t etaMin,
  85                                         Double_t etaMax):
  86   fTrackCollection(new TObjArray(nParticles)),
  87   fReferenceMultiplicity(nParticles),
  88   fNumberOfTracks(0),
  89   fNumberOfRPs(0),
  90   fMCReactionPlaneAngle(0.),
  91   fMCReactionPlaneAngleIsSet(kFALSE),
  92   fAfterBurnerPrecision(0.001),
  93   fNumberOfTracksWrap(NULL),
  94   fNumberOfRPsWrap(NULL),
  95   fMCReactionPlaneAngleWrap(NULL)
  96 {
  97   //ctor. generates nParticles random tracks with given Pt distribution
  98   //phi and eta are uniform
  99   Generate(nParticles,ptDist,phiMin,phiMax,etaMin,etaMax);
 100 }
 101
 102 //-----------------------------------------------------------------------
 103 AliFlowEventSimple::AliFlowEventSimple(const AliFlowEventSimple& anEvent):
 104   TObject(anEvent),
 105   fTrackCollection((TObjArray*)(anEvent.fTrackCollection)->Clone()),
 106   fReferenceMultiplicity(anEvent.fReferenceMultiplicity),
 107   fNumberOfTracks(anEvent.fNumberOfTracks),
 108   fNumberOfRPs(anEvent.fNumberOfRPs),
 109   fMCReactionPlaneAngle(anEvent.fMCReactionPlaneAngle),
 110   fMCReactionPlaneAngleIsSet(anEvent.fMCReactionPlaneAngleIsSet),
 111   fAfterBurnerPrecision(anEvent.fAfterBurnerPrecision),
 112   fNumberOfTracksWrap(anEvent.fNumberOfTracksWrap),
 113   fNumberOfRPsWrap(anEvent.fNumberOfRPsWrap),
 114   fMCReactionPlaneAngleWrap(anEvent.fMCReactionPlaneAngleWrap)
 115 {
 116   //copy constructor
 117 }
 118
 119 //-----------------------------------------------------------------------
 120 AliFlowEventSimple& AliFlowEventSimple::operator=(const AliFlowEventSimple& anEvent)
 121 {
 122   //assignment operator
 123   if (fTrackCollection) fTrackCollection->Delete();
 124   delete fTrackCollection;
 125   fTrackCollection = (TObjArray*)(anEvent.fTrackCollection)->Clone(); //deep copy
 126   fReferenceMultiplicity = anEvent.fReferenceMultiplicity;
 127   fNumberOfTracks = anEvent.fNumberOfTracks;
 128   fNumberOfRPs = anEvent.fNumberOfRPs;
 129   fMCReactionPlaneAngle = anEvent.fMCReactionPlaneAngle;
 130   fMCReactionPlaneAngleIsSet = anEvent.fMCReactionPlaneAngleIsSet;
 131   fAfterBurnerPrecision = anEvent.fAfterBurnerPrecision;
 132   fNumberOfTracksWrap = anEvent.fNumberOfTracksWrap;
 133   fNumberOfRPsWrap = anEvent.fNumberOfRPsWrap;
 134   fMCReactionPlaneAngleWrap=anEvent.fMCReactionPlaneAngleWrap;
 135   return *this;
 136 }
 137
 138 //-----------------------------------------------------------------------
 139 AliFlowEventSimple::~AliFlowEventSimple()
 140 {
 141   //destructor
 142   if (fTrackCollection) fTrackCollection->Delete();
 143   delete fTrackCollection;
 144   if (fNumberOfTracksWrap) delete fNumberOfTracksWrap;
 145   if (fNumberOfRPsWrap) delete fNumberOfRPsWrap;
 146   if (fMCReactionPlaneAngleWrap) delete fMCReactionPlaneAngleWrap;
 147 }
 148
 149 //-----------------------------------------------------------------------
 150 void AliFlowEventSimple::Generate(Int_t nParticles,
 151                                         TF1* ptDist,
 152                                         Double_t phiMin,
 153                                         Double_t phiMax,
 154                                         Double_t etaMin,
 155                                         Double_t etaMax)
 156 {
 157   //generate nParticles random tracks uniform in phi and eta
 158   //according to the specified pt distribution
 159   for (Int_t i=0; i<nParticles; i++)
 160   {
 161     AliFlowTrackSimple* track = new AliFlowTrackSimple();
 162     track->SetPhi( gRandom->Uniform(phiMin,phiMax) );
 163     track->SetEta( gRandom->Uniform(etaMin,etaMax) );
 164     track->SetPt( ptDist->GetRandom() );
 165     track->SetCharge( (gRandom->Uniform()-0.5<0)?-1:1 );
 166     AddTrack(track);
 167   }
 168 }
 169
 170 //-----------------------------------------------------------------------
 171 AliFlowTrackSimple* AliFlowEventSimple::GetTrack(Int_t i)
 172 {
 173   //get track i from collection
 174   if (i>=fNumberOfTracks) return NULL;
 175   AliFlowTrackSimple* pTrack = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i)) ;
 176   return pTrack;
 177 }
 178
 179 //-----------------------------------------------------------------------
 180 void AliFlowEventSimple::AddTrack( AliFlowTrackSimple* track )
 181 {
 182   //add a track
 183   fTrackCollection->AddLast(track);
 184   fNumberOfTracks++;
 185 }
 186
 187 //-----------------------------------------------------------------------
 188 AliFlowVector AliFlowEventSimple::GetQ(Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
 189 {
 190   // calculate Q-vector in harmonic n without weights (default harmonic n=2)
 191   Double_t dQX = 0.;
 192   Double_t dQY = 0.;
 193   AliFlowVector vQ;
 194   vQ.Set(0.,0.);
 195
 196   Int_t iOrder = n;
 197   Double_t iUsedTracks = 0;
 198   Double_t dPhi = 0.;
 199   Double_t dPt = 0.;
 200   Double_t dEta = 0.;
 201   Double_t dWeight = 1.;
 202
 203   AliFlowTrackSimple* pTrack = NULL;
 204
 205   Int_t nBinsPhi = 0;
 206   Double_t dBinWidthPt = 0.;
 207   Double_t dPtMin = 0.;
 208   Double_t dBinWidthEta = 0.;
 209   Double_t dEtaMin = 0.;
 210
 211   Double_t wPhi = 1.; // weight Phi
 212   Double_t wPt = 1.;  // weight Pt
 213   Double_t wEta = 1.; // weight Eta
 214
 215   TH1F *phiWeights = NULL;
 216   TH1D *ptWeights  = NULL;
 217   TH1D *etaWeights = NULL;
 218
 219   if(weightsList)
 220   {
 221     if(usePhiWeights)
 222     {
 223       phiWeights = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights"));
 224       if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
 225     }
 226     if(usePtWeights)
 227     {
 228       ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
 229       if(ptWeights)
 230       {
 231         dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
 232         dPtMin = (ptWeights->GetXaxis())->GetXmin();
 233       }
 234     }
 235     if(useEtaWeights)
 236     {
 237       etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
 238       if(etaWeights)
 239       {
 240         dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
 241         dEtaMin = (etaWeights->GetXaxis())->GetXmin();
 242       }
 243     }
 244   } // end of if(weightsList)
 245
 246   // loop over tracks
 247   for(Int_t i=0; i<fNumberOfTracks; i++)
 248   {
 249     pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
 250     if(pTrack)
 251     {
 252       if(pTrack->InRPSelection())
 253       {
 254         dPhi = pTrack->Phi();
 255         dPt  = pTrack->Pt();
 256         dEta = pTrack->Eta();
 257               dWeight = pTrack->Weight();
 258
 259         // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
 260         if(phiWeights && nBinsPhi)
 261         {
 262           wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
 263         }
 264         // determine v'(pt) weight:
 265         if(ptWeights && dBinWidthPt)
 266         {
 267           wPt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
 268         }
 269         // determine v'(eta) weight:
 270         if(etaWeights && dBinWidthEta)
 271         {
 272           wEta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
 273         }
 274
 275         // building up the weighted Q-vector:
 276         dQX += dWeight*wPhi*wPt*wEta*TMath::Cos(iOrder*dPhi);
 277         dQY += dWeight*wPhi*wPt*wEta*TMath::Sin(iOrder*dPhi);
 278
 279         // weighted multiplicity:
 280         iUsedTracks += dWeight*wPhi*wPt*wEta;
 281
 282       } // end of if (pTrack->InRPSelection())
 283     } // end of if (pTrack)
 284     else
 285     {
 286       cerr << "no particle!!!"<<endl;
 287     }
 288   } // loop over particles
 289
 290   vQ.Set(dQX,dQY);
 291   vQ.SetMult(iUsedTracks);
 292
 293   return vQ;
 294
 295 }
 296
 297 //-----------------------------------------------------------------------
 298 void AliFlowEventSimple::Get2Qsub(AliFlowVector* Qarray, Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
 299 {
 300
 301   // calculate Q-vector in harmonic n without weights (default harmonic n=2)
 302   Double_t dQX = 0.;
 303   Double_t dQY = 0.;
 304
 305   Int_t iOrder = n;
 306   Double_t iUsedTracks = 0;
 307   Double_t dPhi = 0.;
 308   Double_t dPt  = 0.;
 309   Double_t dEta = 0.;
 310   Double_t dWeight = 1.;
 311
 312   AliFlowTrackSimple* pTrack = NULL;
 313
 314   Int_t    iNbinsPhiSub0 = 0;
 315   Int_t    iNbinsPhiSub1 = 0;
 316   Double_t dBinWidthPt = 0.;
 317   Double_t dPtMin      = 0.;
 318   Double_t dBinWidthEta= 0.;
 319   Double_t dEtaMin     = 0.;
 320
 321   Double_t dWphi = 1.;  // weight Phi
 322   Double_t dWpt  = 1.;  // weight Pt
 323   Double_t dWeta = 1.;  // weight Eta
 324
 325   TH1F* phiWeightsSub0 = NULL;
 326   TH1F* phiWeightsSub1 = NULL;
 327   TH1D* ptWeights  = NULL;
 328   TH1D* etaWeights = NULL;
 329
 330   if(weightsList)
 331   {
 332     if(usePhiWeights)
 333     {
 334       phiWeightsSub0 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub0"));
 335       phiWeightsSub1 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub1"));
 336       if(phiWeightsSub0) {
 337         iNbinsPhiSub0 = phiWeightsSub0->GetNbinsX();
 338       }
 339       if(phiWeightsSub1) {
 340         iNbinsPhiSub1 = phiWeightsSub1->GetNbinsX();
 341       }
 342     }
 343     if(usePtWeights)
 344     {
 345       ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
 346       if(ptWeights)
 347       {
 348         dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
 349         dPtMin = (ptWeights->GetXaxis())->GetXmin();
 350       }
 351     }
 352     if(useEtaWeights)
 353     {
 354       etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
 355       if(etaWeights)
 356       {
 357         dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
 358         dEtaMin = (etaWeights->GetXaxis())->GetXmin();
 359       }
 360     }
 361   } // end of if(weightsList)
 362
 363   //loop over the two subevents
 364   for (Int_t s=0; s<2; s++)
 365   {
 366     // loop over tracks
 367     for(Int_t i=0; i<fNumberOfTracks; i++)
 368     {
 369       pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
 370       if(pTrack)
 371       {
 372         if(pTrack->InRPSelection())
 373         {
 374           if (pTrack->InSubevent(s))
 375           {
 376             dPhi    = pTrack->Phi();
 377             dPt     = pTrack->Pt();
 378             dEta    = pTrack->Eta();
 379             dWeight = pTrack->Weight();
 380
 381             // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
 382             if(s == 0) { //subevent 0
 383               if(phiWeightsSub0 && iNbinsPhiSub0){
 384                 dWphi = phiWeightsSub0->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub0/TMath::TwoPi())));
 385               }
 386             } else if (s == 1) {
 387               if(phiWeightsSub1 && iNbinsPhiSub1){
 388                 dWphi = phiWeightsSub1->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub1/TMath::TwoPi())));
 389               }
 390             }
 391
 392             // determine v'(pt) weight:
 393             if(ptWeights && dBinWidthPt)
 394             {
 395               dWpt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
 396             }
 397
 398             // determine v'(eta) weight:
 399             if(etaWeights && dBinWidthEta)
 400             {
 401               dWeta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
 402             }
 403
 404             // building up the weighted Q-vector:
 405             dQX += dWeight*dWphi*dWpt*dWeta*TMath::Cos(iOrder*dPhi);
 406             dQY += dWeight*dWphi*dWpt*dWeta*TMath::Sin(iOrder*dPhi);
 407
 408             // weighted multiplicity:
 409             iUsedTracks+=dWeight*dWphi*dWpt*dWeta;
 410
 411           } // end of subevent
 412         } // end of if (pTrack->InRPSelection())
 413       } // end of if (pTrack)
 414       else
 415       {
 416         cerr << "no particle!!!"<<endl;
 417       }
 418     } // loop over particles
 419     Qarray[s].Set(dQX,dQY);
 420     Qarray[s].SetMult(iUsedTracks);
 421     //reset
 422     iUsedTracks = 0;
 423     dQX = 0.;
 424     dQY = 0.;
 425   }
 426
 427 }
 428
 429
 430 //-----------------------------------------------------------------------
 431 void AliFlowEventSimple::Print(Option_t *option) const
 432 {
 433   //   -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
 434   //             ===============================================
 435   //   printf( "TH1.Print Name  = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
 436   printf( "Class.Print Name = %s, #tracks= %d, Number of RPs= %d, MC EventPlaneAngle= %f\n",
 437           GetName(),fNumberOfTracks, fNumberOfRPs, fMCReactionPlaneAngle );
 438
 439   TString optionstr(option);
 440   if (!optionstr.Contains("all")) return;
 441   if (fTrackCollection)
 442   {
 443     fTrackCollection->Print(option);
 444   }
 445   else
 446   {
 447     printf( "Empty track collection \n");
 448   }
 449 }
 450
 451 //-----------------------------------------------------------------------
 452 void AliFlowEventSimple::Browse(TBrowser *b)
 453 {
 454   if (!b) return;
 455   if (!fNumberOfTracksWrap)
 456   {
 457     fNumberOfTracksWrap = new TParameter<int>("fNumberOfTracks", fNumberOfTracks);
 458     b->Add(fNumberOfTracksWrap);
 459   }
 460   if (!fNumberOfRPsWrap)
 461   {
 462     fNumberOfRPsWrap = new TParameter<int>("fNumberOfRPs", fNumberOfRPs);
 463     b->Add(fNumberOfRPsWrap);
 464   }
 465   if (!fMCReactionPlaneAngleWrap)
 466   {
 467     fMCReactionPlaneAngleWrap = new TParameter<double>(" fMCReactionPlaneAngle",  fMCReactionPlaneAngle);
 468     b->Add( fMCReactionPlaneAngleWrap);
 469   }
 470   if (fTrackCollection) b->Add(fTrackCollection,"AliFlowTracksSimple");
 471 }
 472
 473 //-----------------------------------------------------------------------
 474 AliFlowEventSimple::AliFlowEventSimple( TTree* inputTree,
 475                                         const AliFlowTrackSimpleCuts* rpCuts,
 476                                         const AliFlowTrackSimpleCuts* poiCuts):
 477   fTrackCollection(NULL),
 478   fReferenceMultiplicity(0),
 479   fNumberOfTracks(0),
 480   fNumberOfRPs(0),
 481   fMCReactionPlaneAngle(0.),
 482   fMCReactionPlaneAngleIsSet(kFALSE),
 483   fAfterBurnerPrecision(0.001),
 484   fNumberOfTracksWrap(NULL),
 485   fNumberOfRPsWrap(NULL),
 486   fMCReactionPlaneAngleWrap(NULL)
 487 {
 488   //constructor, fills the event from a TTree of kinematic.root files
 489   //applies RP and POI cuts, tags the tracks
 490
 491   Int_t numberOfInputTracks = inputTree->GetEntries() ;
 492   fTrackCollection = new TObjArray(numberOfInputTracks/2);
 493
 494   TParticle* pParticle = new TParticle();
 495   inputTree->SetBranchAddress("Particles",&pParticle);
 496
 497   Int_t iSelParticlesPOI = 0;
 498
 499   for (Int_t i=0; i<numberOfInputTracks; i++)
 500   {
 501     inputTree->GetEntry(i);   //get input particle
 502
 503     if (!pParticle) continue; //no particle
 504     if (!pParticle->IsPrimary()) continue;
 505
 506     Bool_t rpOK = rpCuts->PassesCuts(pParticle);
 507     Bool_t poiOK = poiCuts->PassesCuts(pParticle);
 508
 509     if (rpOK || poiOK)
 510     {
 511       AliFlowTrackSimple* pTrack = new AliFlowTrackSimple(pParticle);
 512
 513       //marking the particles used for int. flow:
 514       if(rpOK)
 515       {
 516         pTrack->SetForRPSelection(kTRUE);
 517         fNumberOfRPs++;
 518       }
 519       //marking the particles used for diff. flow:
 520       if(poiOK)
 521       {
 522         pTrack->SetForPOISelection(kTRUE);
 523         iSelParticlesPOI++;
 524       }
 525       //adding a particles which were used either for int. or diff. flow to the list
 526       AddTrack(pTrack);
 527     }
 528   }//for i
 529   delete pParticle;
 530 }
 531
 532 //_____________________________________________________________________________
 533 void AliFlowEventSimple::CloneTracks(Int_t n)
 534 {
 535   //clone every track n times to add non-flow
 536   if (n<=0) return; //no use to clone stuff zero or less times
 537   Int_t ntracks = fNumberOfTracks;
 538   fTrackCollection->Expand((n+1)*fNumberOfTracks);
 539   for (Int_t i=0; i<n; i++)
 540   {
 541     for (Int_t itrack=0; itrack<ntracks; itrack++)
 542     {
 543       AliFlowTrackSimple* track = dynamic_cast<AliFlowTrackSimple*>(fTrackCollection->At(itrack));
 544       if (!track) continue;
 545       AddTrack(static_cast<AliFlowTrackSimple*>(track->Clone()));
 546     }
 547   }
 548 }
 549
 550 //_____________________________________________________________________________
 551 void AliFlowEventSimple::ResolutionPt(Double_t res)
 552 {
 553   //smear pt of all tracks by gaussian with sigma=res
 554   for (Int_t i=0; i<fNumberOfTracks; i++)
 555   {
 556     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 557     if (track) track->ResolutionPt(res);
 558   }
 559 }
 560
 561 //_____________________________________________________________________________
 562 void AliFlowEventSimple::TagSubeventsInEta( Double_t etaMinA,
 563                                             Double_t etaMaxA,
 564                                             Double_t etaMinB,
 565                                             Double_t etaMaxB )
 566 {
 567   //Flag two subevents in given eta ranges
 568   for (Int_t i=0; i<fNumberOfTracks; i++)
 569   {
 570     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 571     Double_t eta=track->Eta();
 572     if (eta >= etaMinA && eta <= etaMaxA) track->SetForSubevent(0);
 573     if (eta >= etaMinB && eta <= etaMaxB) track->SetForSubevent(1);
 574   }
 575 }
 576
 577 //_____________________________________________________________________________
 578 void AliFlowEventSimple::AddV1( Double_t v1 )
 579 {
 580   //add v2 to all tracks wrt the reaction plane angle
 581   for (Int_t i=0; i<fNumberOfTracks; i++)
 582   {
 583     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 584     if (track) track->AddV1(v1, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 585   }
 586 }
 587
 588 //_____________________________________________________________________________
 589 void AliFlowEventSimple::AddV2( Double_t v2 )
 590 {
 591   //add v2 to all tracks wrt the reaction plane angle
 592   for (Int_t i=0; i<fNumberOfTracks; i++)
 593   {
 594     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 595     if (track) track->AddV2(v2, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 596   }
 597 }
 598
 599 //_____________________________________________________________________________
 600 void AliFlowEventSimple::AddV4( Double_t v4 )
 601 {
 602   //add v4 to all tracks wrt the reaction plane angle
 603   for (Int_t i=0; i<fNumberOfTracks; i++)
 604   {
 605     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 606     if (track) track->AddV4(v4, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 607   }
 608 }
 609
 610 //_____________________________________________________________________________
 611 void AliFlowEventSimple::AddFlow( Double_t v1, Double_t v2, Double_t v3, Double_t v4 )
 612 {
 613   //add flow to all tracks wrt the reaction plane angle
 614   for (Int_t i=0; i<fNumberOfTracks; i++)
 615   {
 616     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 617     if (track) track->AddFlow(v1,v2,v3,v4,fMCReactionPlaneAngle, fAfterBurnerPrecision);
 618   }
 619 }
 620
 621 //_____________________________________________________________________________
 622 void AliFlowEventSimple::TagRP( AliFlowTrackSimpleCuts* cuts )
 623 {
 624   //tag tracks as reference particles (RPs)
 625   for (Int_t i=0; i<fNumberOfTracks; i++)
 626   {
 627     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 628     if (!track) continue;
 629     Bool_t pass=cuts->PassesCuts(track);
 630     track->SetForRPSelection(pass);
 631     if (pass)
 632       fNumberOfRPs++;
 633     else
 634       fNumberOfRPs--;
 635   }
 636 }
 637
 638 //_____________________________________________________________________________
 639 void AliFlowEventSimple::TagPOI( AliFlowTrackSimpleCuts* cuts )
 640 {
 641   //tag tracks as particles of interest (POIs)
 642   for (Int_t i=0; i<fNumberOfTracks; i++)
 643   {
 644     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 645     if (!track) continue;
 646     Bool_t pass=cuts->PassesCuts(track);
 647     track->SetForPOISelection(pass);
 648   }
 649 }
 650
 651 //_____________________________________________________________________________
 652 void AliFlowEventSimple::DefineDeadZone( Double_t etaMin,
 653                                          Double_t etaMax,
 654                                          Double_t phiMin,
 655                                          Double_t phiMax )
 656 {
 657   //mark tracks in given eta-phi region as dead
 658   //by resetting the flow bits
 659   for (Int_t i=0; i<fNumberOfTracks; i++)
 660   {
 661     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 662     Double_t eta = track->Eta();
 663     Double_t phi = track->Phi();
 664     if (eta>etaMin && eta<etaMax && phi>phiMin && phi<phiMax)
 665     {
 666       if (track->InRPSelection()) fNumberOfRPs--;
 667       track->ResetFlowTags();
 668     }
 669   }
 670 }
 671
 672 //_____________________________________________________________________________
 673 Int_t AliFlowEventSimple::CleanUpDeadTracks()
 674 {
 675   //remove tracks that have no flow tags set and cleanup the container
 676   //returns number of cleaned tracks
 677   Int_t ncleaned=0;
 678   for (Int_t i=0; i<fNumberOfTracks; i++)
 679   {
 680     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 681     if (track->IsDead()) {delete track;ncleaned++;}
 682   }
 683   fTrackCollection->Compress(); //clean up empty slots
 684   return ncleaned;
 685 }