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