PWG/FLOW/Base/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 using std::cout;
  46 using std::endl;
  47 using std::cerr;
  48 ClassImp(AliFlowEventSimple)
  49
  50 //-----------------------------------------------------------------------
  51 AliFlowEventSimple::AliFlowEventSimple():
  52   fTrackCollection(NULL),
  53   fReferenceMultiplicity(0),
  54   fNumberOfTracks(0),
  55   fUseGlauberMCSymmetryPlanes(kFALSE),
  56   fUseExternalSymmetryPlanes(kFALSE),
  57   fPsi1(0.),
  58   fPsi2(0.),
  59   fPsi3(0.),
  60   fPsi4(0.),
  61   fPsi5(0.),
  62   fPsi1Psi3(0x0),
  63   fPsi2Psi4(0x0),
  64   fPsi3Psi5(0x0),
  65   fMCReactionPlaneAngle(0.),
  66   fMCReactionPlaneAngleIsSet(kFALSE),
  67   fAfterBurnerPrecision(0.001),
  68   fUserModified(kFALSE),
  69   fNumberOfTracksWrap(NULL),
  70   fNumberOfRPsWrap(NULL),
  71   fNumberOfPOIsWrap(NULL),
  72   fMCReactionPlaneAngleWrap(NULL),
  73   fShuffledIndexes(NULL),
  74   fShuffleTracks(kFALSE),
  75   fMothersCollection(NULL),
  76   fCentrality(-1.),
  77   fNumberOfPOItypes(2),
  78   fNumberOfPOIs(NULL)
  79 {
  80   cout << "AliFlowEventSimple: Default constructor to be used only by root for io" << endl;
  81 }
  82
  83 //-----------------------------------------------------------------------
  84 AliFlowEventSimple::AliFlowEventSimple( Int_t n,
  85                                         ConstructionMethod method,
  86                                         TF1* ptDist,
  87                                         Double_t phiMin,
  88                                         Double_t phiMax,
  89                                         Double_t etaMin,
  90                                         Double_t etaMax):
  91   fTrackCollection(new TObjArray(n)),
  92   fReferenceMultiplicity(0),
  93   fNumberOfTracks(0),
  94   fUseGlauberMCSymmetryPlanes(kFALSE),
  95   fUseExternalSymmetryPlanes(kFALSE),
  96   fPsi1(0.),
  97   fPsi2(0.),
  98   fPsi3(0.),
  99   fPsi4(0.),
 100   fPsi5(0.),
 101   fPsi1Psi3(0x0),
 102   fPsi2Psi4(0x0),
 103   fPsi3Psi5(0x0),
 104   fMCReactionPlaneAngle(0.),
 105   fMCReactionPlaneAngleIsSet(kFALSE),
 106   fAfterBurnerPrecision(0.001),
 107   fUserModified(kFALSE),
 108   fNumberOfTracksWrap(NULL),
 109   fNumberOfRPsWrap(NULL),
 110   fNumberOfPOIsWrap(NULL),
 111   fMCReactionPlaneAngleWrap(NULL),
 112   fShuffledIndexes(NULL),
 113   fShuffleTracks(kFALSE),
 114   fMothersCollection(new TObjArray()),
 115   fCentrality(-1.),
 116   fNumberOfPOItypes(2),
 117   fNumberOfPOIs(new Int_t[fNumberOfPOItypes])
 118 {
 119   //ctor
 120   // if second argument is set to AliFlowEventSimple::kGenerate
 121   // it generates n random tracks with given Pt distribution
 122   // (a sane default is provided), phi and eta are uniform
 123
 124   if (method==kGenerate)
 125     Generate(n,ptDist,phiMin,phiMax,etaMin,etaMax);
 126 }
 127
 128 //-----------------------------------------------------------------------
 129 AliFlowEventSimple::AliFlowEventSimple(const AliFlowEventSimple& anEvent):
 130   TObject(anEvent),
 131   fTrackCollection((TObjArray*)(anEvent.fTrackCollection)->Clone()),
 132   fReferenceMultiplicity(anEvent.fReferenceMultiplicity),
 133   fNumberOfTracks(anEvent.fNumberOfTracks),
 134   fUseGlauberMCSymmetryPlanes(anEvent.fUseGlauberMCSymmetryPlanes),
 135   fUseExternalSymmetryPlanes(anEvent.fUseExternalSymmetryPlanes),
 136   fPsi1(anEvent.fPsi1),
 137   fPsi2(anEvent.fPsi2),
 138   fPsi3(anEvent.fPsi3),
 139   fPsi4(anEvent.fPsi4),
 140   fPsi5(anEvent.fPsi5),
 141   fPsi1Psi3(anEvent.fPsi1Psi3),
 142   fPsi2Psi4(anEvent.fPsi2Psi4),
 143   fPsi3Psi5(anEvent.fPsi3Psi5),
 144   fMCReactionPlaneAngle(anEvent.fMCReactionPlaneAngle),
 145   fMCReactionPlaneAngleIsSet(anEvent.fMCReactionPlaneAngleIsSet),
 146   fAfterBurnerPrecision(anEvent.fAfterBurnerPrecision),
 147   fUserModified(anEvent.fUserModified),
 148   fNumberOfTracksWrap(anEvent.fNumberOfTracksWrap),
 149   fNumberOfRPsWrap(anEvent.fNumberOfRPsWrap),
 150   fNumberOfPOIsWrap(anEvent.fNumberOfPOIsWrap),
 151   fMCReactionPlaneAngleWrap(anEvent.fMCReactionPlaneAngleWrap),
 152   fShuffledIndexes(NULL),
 153   fShuffleTracks(anEvent.fShuffleTracks),
 154   fMothersCollection(new TObjArray()),
 155   fCentrality(anEvent.fCentrality),
 156   fNumberOfPOItypes(anEvent.fNumberOfPOItypes),
 157   fNumberOfPOIs(new Int_t[fNumberOfPOItypes])
 158 {
 159   //copy constructor
 160   memcpy(fNumberOfPOIs,anEvent.fNumberOfPOIs,fNumberOfPOItypes*sizeof(Int_t));
 161 }
 162
 163 //-----------------------------------------------------------------------
 164 void AliFlowEventSimple::SetNumberOfPOIs( Int_t numberOfPOIs, Int_t poiType)
 165 {
 166   //set the number of poi classes, resize the array if larger is needed
 167   //never shrink the array
 168   //never decrease the stored number
 169   if (poiType>=fNumberOfPOItypes)
 170   {
 171     Int_t n = poiType+1;
 172     Int_t* tmp = new Int_t[n];
 173     for (Int_t j=0; j<n; j++) { tmp[j]=0; }
 174     memcpy(tmp,fNumberOfPOIs,fNumberOfPOItypes*sizeof(Int_t));
 175     delete [] fNumberOfPOIs;
 176     fNumberOfPOIs = tmp;
 177     fNumberOfPOItypes = n;
 178   }
 179
 180   fNumberOfPOIs[poiType] = numberOfPOIs;
 181 }
 182
 183 //-----------------------------------------------------------------------
 184 void AliFlowEventSimple::IncrementNumberOfPOIs(Int_t poiType)
 185 {
 186
 187   if (poiType>=fNumberOfPOItypes) SetNumberOfPOIs(0,poiType);
 188   fNumberOfPOIs[poiType]++;
 189 }
 190
 191 //-----------------------------------------------------------------------
 192 AliFlowEventSimple& AliFlowEventSimple::operator=(const AliFlowEventSimple& anEvent)
 193 {
 194   //assignment operator
 195   if (&anEvent==this) return *this; //check self-assignment
 196   if (fTrackCollection) fTrackCollection->Delete();
 197   delete fTrackCollection;
 198   fTrackCollection = (TObjArray*)(anEvent.fTrackCollection)->Clone(); //deep copy
 199   fReferenceMultiplicity = anEvent.fReferenceMultiplicity;
 200   fNumberOfTracks = anEvent.fNumberOfTracks;
 201   fNumberOfPOItypes = anEvent.fNumberOfPOItypes;
 202   delete [] fNumberOfPOIs;
 203   fNumberOfPOIs=new Int_t[fNumberOfPOItypes];
 204   memcpy(fNumberOfPOIs,anEvent.fNumberOfPOIs,fNumberOfPOItypes*sizeof(Int_t));
 205   fUseGlauberMCSymmetryPlanes = anEvent.fUseGlauberMCSymmetryPlanes;
 206   fUseExternalSymmetryPlanes = anEvent.fUseExternalSymmetryPlanes;
 207   fPsi1 = anEvent.fPsi1;
 208   fPsi2 = anEvent.fPsi2;
 209   fPsi3 = anEvent.fPsi3;
 210   fPsi4 = anEvent.fPsi4;
 211   fPsi5 = anEvent.fPsi5;
 212   fPsi1Psi3 = anEvent.fPsi1Psi3;
 213   fPsi2Psi4 = anEvent.fPsi2Psi4;
 214   fPsi3Psi5 = anEvent.fPsi3Psi5;
 215   fMCReactionPlaneAngle = anEvent.fMCReactionPlaneAngle;
 216   fMCReactionPlaneAngleIsSet = anEvent.fMCReactionPlaneAngleIsSet;
 217   fAfterBurnerPrecision = anEvent.fAfterBurnerPrecision;
 218   fUserModified=anEvent.fUserModified;
 219   fNumberOfTracksWrap = anEvent.fNumberOfTracksWrap;
 220   fNumberOfRPsWrap = anEvent.fNumberOfRPsWrap;
 221   fNumberOfPOIsWrap = anEvent.fNumberOfPOIsWrap;
 222   fMCReactionPlaneAngleWrap=anEvent.fMCReactionPlaneAngleWrap;
 223   fShuffleTracks=anEvent.fShuffleTracks;
 224   delete [] fShuffledIndexes;
 225   return *this;
 226 }
 227
 228 //-----------------------------------------------------------------------
 229 AliFlowEventSimple::~AliFlowEventSimple()
 230 {
 231   //destructor
 232   if (fTrackCollection) fTrackCollection->Delete();
 233   delete fTrackCollection;
 234   delete fNumberOfTracksWrap;
 235   delete fNumberOfRPsWrap;
 236   delete fNumberOfPOIsWrap;
 237   delete fMCReactionPlaneAngleWrap;
 238   delete fShuffledIndexes;
 239   delete fMothersCollection;
 240   delete [] fNumberOfPOIs;
 241 }
 242
 243 //-----------------------------------------------------------------------
 244 void AliFlowEventSimple::SetUseExternalSymmetryPlanes(TF1 *gPsi1Psi3,
 245                                                       TF1 *gPsi2Psi4,
 246                                                       TF1 *gPsi3Psi5) {
 247   //Use symmetry planes, setup correlations between different Psi_n
 248   fUseExternalSymmetryPlanes = kTRUE;
 249
 250   //Correlations between Psi_1 and Psi_3
 251   if(gPsi1Psi3) fPsi1Psi3 = gPsi1Psi3;
 252   else {
 253     fPsi1Psi3 = new TF1("fPsi1Psi3","[0]*x+[1]",0.,2.*TMath::Pi());
 254     fPsi1Psi3->SetParameter(0,1.);
 255     fPsi1Psi3->SetParameter(1,0.);
 256   }
 257
 258   //Correlations between Psi_2 and Psi_4
 259   if(gPsi2Psi4) fPsi2Psi4 = gPsi2Psi4;
 260   else {
 261     fPsi2Psi4 = new TF1("fPsi2Psi4","[0]*x+[1]",0.,2.*TMath::Pi());
 262     fPsi2Psi4->SetParameter(0,1.);
 263     fPsi2Psi4->SetParameter(1,0.);
 264   }
 265
 266   //Correlations between Psi_3 and Psi_5
 267   if(gPsi3Psi5) fPsi3Psi5 = gPsi3Psi5;
 268   else {
 269     fPsi3Psi5 = new TF1("fPsi3Psi5","[0]*x+[1]",0.,2.*TMath::Pi());
 270     fPsi3Psi5->SetParameter(0,1.);
 271     fPsi3Psi5->SetParameter(1,0.);
 272   }
 273 }
 274
 275 //-----------------------------------------------------------------------
 276 void AliFlowEventSimple::Generate(Int_t nParticles,
 277                                   TF1* ptDist,
 278                                   Double_t phiMin,
 279                                   Double_t phiMax,
 280                                   Double_t etaMin,
 281                                   Double_t etaMax)
 282 {
 283   //generate nParticles random tracks uniform in phi and eta
 284   //according to the specified pt distribution
 285   if (!ptDist)
 286   {
 287     static TF1 ptdistribution("ptSpectra","x*TMath::Exp(-pow(0.13957*0.13957+x*x,0.5)/0.4)",0.1,10.);
 288     ptDist=&ptdistribution;
 289   }
 290
 291   for (Int_t i=0; i<nParticles; i++)
 292   {
 293     AliFlowTrackSimple* track = new AliFlowTrackSimple();
 294     track->SetPhi( gRandom->Uniform(phiMin,phiMax) );
 295     track->SetEta( gRandom->Uniform(etaMin,etaMax) );
 296     track->SetPt( ptDist->GetRandom() );
 297     track->SetCharge( (gRandom->Uniform()-0.5<0)?-1:1 );
 298     AddTrack(track);
 299   }
 300   if(fUseExternalSymmetryPlanes) {
 301     Double_t betaParameter = gRandom->Gaus(0.,1.3);
 302     fPsi1Psi3->SetParameter(1,betaParameter);
 303
 304     betaParameter = gRandom->Gaus(0.,0.9);
 305     fPsi2Psi4->SetParameter(1,betaParameter);
 306
 307     betaParameter = gRandom->Gaus(0.,1.5);
 308     fPsi3Psi5->SetParameter(1,betaParameter);
 309
 310     fPsi1 = gRandom->Uniform(2.*TMath::Pi());
 311     fPsi2 = gRandom->Uniform(2.*TMath::Pi());
 312     fPsi3 = fPsi1Psi3->Eval(fPsi1);
 313     if(fPsi3 < 0) fPsi3 += 2.*TMath::Pi();
 314     else if(fPsi3 > 2.*TMath::Pi()) fPsi3 -= 2.*TMath::Pi();
 315     fPsi4 = fPsi2Psi4->Eval(fPsi2);
 316     if(fPsi4 < 0) fPsi4 += 2.*TMath::Pi();
 317     else if(fPsi4 > 2.*TMath::Pi()) fPsi4 -= 2.*TMath::Pi();
 318     fPsi5 = fPsi3Psi5->Eval(fPsi3);
 319     if(fPsi5 < 0) fPsi5 += 2.*TMath::Pi();
 320     else if(fPsi5 > 2.*TMath::Pi()) fPsi5 -= 2.*TMath::Pi();
 321
 322     fMCReactionPlaneAngle=fPsi2;
 323     fMCReactionPlaneAngleIsSet=kTRUE;
 324   }
 325   else {
 326     fMCReactionPlaneAngle=gRandom->Uniform(0.0,TMath::TwoPi());
 327     fMCReactionPlaneAngleIsSet=kTRUE;
 328   }
 329   SetUserModified();
 330 }
 331
 332 //-----------------------------------------------------------------------
 333 AliFlowTrackSimple* AliFlowEventSimple::GetTrack(Int_t i)
 334 {
 335   //get track i from collection
 336   if (i>=fNumberOfTracks) return NULL;
 337   Int_t trackIndex=i;
 338   //if asked use the shuffled index
 339   if (fShuffleTracks)
 340   {
 341     if (!fShuffledIndexes) ShuffleTracks();
 342     trackIndex=fShuffledIndexes[i];
 343   }
 344   AliFlowTrackSimple* pTrack = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(trackIndex)) ;
 345   return pTrack;
 346 }
 347
 348 //-----------------------------------------------------------------------
 349 void AliFlowEventSimple::ShuffleTracks()
 350 {
 351   //shuffle track indexes
 352   if (!fShuffledIndexes)
 353   {
 354     //initialize the table with shuffled indexes
 355     fShuffledIndexes = new Int_t[fNumberOfTracks];
 356     for (Int_t j=0; j<fNumberOfTracks; j++) { fShuffledIndexes[j]=j; }
 357   }
 358   //shuffle
 359   std::random_shuffle(&fShuffledIndexes[0], &fShuffledIndexes[fNumberOfTracks]);
 360   Printf("Tracks shuffled! tracks: %i",fNumberOfTracks);
 361 }
 362
 363 //-----------------------------------------------------------------------
 364 void AliFlowEventSimple::AddTrack( AliFlowTrackSimple* track )
 365 {
 366   //add a track, delete the old one if necessary
 367   if (fNumberOfTracks < fTrackCollection->GetEntriesFast())
 368   {
 369     TObject* o = fTrackCollection->At(fNumberOfTracks);
 370     delete o;
 371   }
 372   fTrackCollection->AddAtAndExpand(track,fNumberOfTracks);
 373   if (track->GetNDaughters()>0)
 374   {
 375     //if there track has daughters cache in the collection of mothers
 376     fMothersCollection->Add(track);
 377   }
 378   TrackAdded();
 379 }
 380
 381 //-----------------------------------------------------------------------
 382 void AliFlowEventSimple::TrackAdded()
 383 {
 384   //book keeping after a new track has been added
 385   fNumberOfTracks++;
 386   if (fShuffledIndexes)
 387   {
 388     delete [] fShuffledIndexes;
 389     fShuffledIndexes=NULL;
 390   }
 391 }
 392
 393 //-----------------------------------------------------------------------
 394 AliFlowTrackSimple* AliFlowEventSimple::MakeNewTrack()
 395 {
 396    AliFlowTrackSimple *t=dynamic_cast<AliFlowTrackSimple *>(fTrackCollection->RemoveAt(fNumberOfTracks));
 397    if( !t ) {  // If there was no track at the end of the list then create a new track
 398       t=new AliFlowTrackSimple();
 399    }
 400
 401    return t;
 402 }
 403
 404 //-----------------------------------------------------------------------
 405 AliFlowVector AliFlowEventSimple::GetQ( Int_t n,
 406                                         TList *weightsList,
 407                                         Bool_t usePhiWeights,
 408                                         Bool_t usePtWeights,
 409                                         Bool_t useEtaWeights )
 410 {
 411   // calculate Q-vector in harmonic n without weights (default harmonic n=2)
 412   Double_t dQX = 0.;
 413   Double_t dQY = 0.;
 414   AliFlowVector vQ;
 415   vQ.Set(0.,0.);
 416
 417   Int_t iOrder = n;
 418   Double_t sumOfWeights = 0.;
 419   Double_t dPhi = 0.;
 420   Double_t dPt = 0.;
 421   Double_t dEta = 0.;
 422   Double_t dWeight = 1.;
 423
 424   AliFlowTrackSimple* pTrack = NULL;
 425
 426   Int_t nBinsPhi = 0;
 427   Double_t dBinWidthPt = 0.;
 428   Double_t dPtMin = 0.;
 429   Double_t dBinWidthEta = 0.;
 430   Double_t dEtaMin = 0.;
 431
 432   Double_t wPhi = 1.; // weight Phi
 433   Double_t wPt = 1.;  // weight Pt
 434   Double_t wEta = 1.; // weight Eta
 435
 436   TH1F *phiWeights = NULL;
 437   TH1D *ptWeights  = NULL;
 438   TH1D *etaWeights = NULL;
 439
 440   if(weightsList)
 441   {
 442     if(usePhiWeights)
 443     {
 444       phiWeights = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights"));
 445       if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
 446     }
 447     if(usePtWeights)
 448     {
 449       ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
 450       if(ptWeights)
 451       {
 452         dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
 453         dPtMin = (ptWeights->GetXaxis())->GetXmin();
 454       }
 455     }
 456     if(useEtaWeights)
 457     {
 458       etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
 459       if(etaWeights)
 460       {
 461         dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
 462         dEtaMin = (etaWeights->GetXaxis())->GetXmin();
 463       }
 464     }
 465   } // end of if(weightsList)
 466
 467   // loop over tracks
 468   for(Int_t i=0; i<fNumberOfTracks; i++)
 469   {
 470     pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
 471     if(pTrack)
 472     {
 473       if(pTrack->InRPSelection())
 474       {
 475         dPhi = pTrack->Phi();
 476         dPt  = pTrack->Pt();
 477         dEta = pTrack->Eta();
 478               dWeight = pTrack->Weight();
 479
 480         // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
 481         if(phiWeights && nBinsPhi)
 482         {
 483           wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
 484         }
 485         // determine v'(pt) weight:
 486         if(ptWeights && dBinWidthPt)
 487         {
 488           wPt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
 489         }
 490         // determine v'(eta) weight:
 491         if(etaWeights && dBinWidthEta)
 492         {
 493           wEta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
 494         }
 495
 496         // building up the weighted Q-vector:
 497         dQX += dWeight*wPhi*wPt*wEta*TMath::Cos(iOrder*dPhi);
 498         dQY += dWeight*wPhi*wPt*wEta*TMath::Sin(iOrder*dPhi);
 499
 500         // weighted multiplicity:
 501         sumOfWeights += dWeight*wPhi*wPt*wEta;
 502
 503       } // end of if (pTrack->InRPSelection())
 504     } // end of if (pTrack)
 505     else
 506     {
 507       cerr << "no particle!!!"<<endl;
 508     }
 509   } // loop over particles
 510
 511   vQ.Set(dQX,dQY);
 512   vQ.SetMult(sumOfWeights);
 513   vQ.SetHarmonic(iOrder);
 514   vQ.SetPOItype(AliFlowTrackSimple::kRP);
 515   vQ.SetSubeventNumber(-1);
 516
 517   return vQ;
 518
 519 }
 520
 521 //-----------------------------------------------------------------------
 522 void AliFlowEventSimple::Get2Qsub( AliFlowVector* Qarray,
 523                                    Int_t n,
 524                                    TList *weightsList,
 525                                    Bool_t usePhiWeights,
 526                                    Bool_t usePtWeights,
 527                                    Bool_t useEtaWeights )
 528 {
 529
 530   // calculate Q-vector in harmonic n without weights (default harmonic n=2)
 531   Double_t dQX = 0.;
 532   Double_t dQY = 0.;
 533
 534   Int_t iOrder = n;
 535   Double_t sumOfWeights = 0.;
 536   Double_t dPhi = 0.;
 537   Double_t dPt  = 0.;
 538   Double_t dEta = 0.;
 539   Double_t dWeight = 1.;
 540
 541   AliFlowTrackSimple* pTrack = NULL;
 542
 543   Int_t    iNbinsPhiSub0 = 0;
 544   Int_t    iNbinsPhiSub1 = 0;
 545   Double_t dBinWidthPt = 0.;
 546   Double_t dPtMin      = 0.;
 547   Double_t dBinWidthEta= 0.;
 548   Double_t dEtaMin     = 0.;
 549
 550   Double_t dWphi = 1.;  // weight Phi
 551   Double_t dWpt  = 1.;  // weight Pt
 552   Double_t dWeta = 1.;  // weight Eta
 553
 554   TH1F* phiWeightsSub0 = NULL;
 555   TH1F* phiWeightsSub1 = NULL;
 556   TH1D* ptWeights  = NULL;
 557   TH1D* etaWeights = NULL;
 558
 559   if(weightsList)
 560   {
 561     if(usePhiWeights)
 562     {
 563       phiWeightsSub0 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub0"));
 564       if(phiWeightsSub0) {
 565         iNbinsPhiSub0 = phiWeightsSub0->GetNbinsX();
 566       }
 567       phiWeightsSub1 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub1"));
 568       if(phiWeightsSub1) {
 569         iNbinsPhiSub1 = phiWeightsSub1->GetNbinsX();
 570       }
 571     }
 572     if(usePtWeights)
 573     {
 574       ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
 575       if(ptWeights)
 576       {
 577         dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
 578         dPtMin = (ptWeights->GetXaxis())->GetXmin();
 579       }
 580     }
 581     if(useEtaWeights)
 582     {
 583       etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
 584       if(etaWeights)
 585       {
 586         dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
 587         dEtaMin = (etaWeights->GetXaxis())->GetXmin();
 588       }
 589     }
 590   } // end of if(weightsList)
 591
 592   //loop over the two subevents
 593   for (Int_t s=0; s<2; s++)
 594   {
 595     // loop over tracks
 596     for(Int_t i=0; i<fNumberOfTracks; i++)
 597     {
 598       pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
 599       if(!pTrack)
 600       {
 601         cerr << "no particle!!!"<<endl;
 602         continue;
 603       }
 604       if(pTrack->InRPSelection() && (pTrack->InSubevent(s)))
 605       {
 606         dPhi    = pTrack->Phi();
 607         dPt     = pTrack->Pt();
 608         dEta    = pTrack->Eta();
 609         dWeight = pTrack->Weight();
 610
 611         // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
 612         //subevent 0
 613         if(s == 0)  {
 614           if(phiWeightsSub0 && iNbinsPhiSub0)  {
 615             Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub0/TMath::TwoPi()));
 616             //use the phi value at the center of the bin
 617             dPhi  = phiWeightsSub0->GetBinCenter(phiBin);
 618             dWphi = phiWeightsSub0->GetBinContent(phiBin);
 619           }
 620         }
 621         //subevent 1
 622         else if (s == 1) {
 623           if(phiWeightsSub1 && iNbinsPhiSub1) {
 624             Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub1/TMath::TwoPi()));
 625             //use the phi value at the center of the bin
 626             dPhi  = phiWeightsSub1->GetBinCenter(phiBin);
 627             dWphi = phiWeightsSub1->GetBinContent(phiBin);
 628           }
 629         }
 630
 631         // determine v'(pt) weight:
 632         if(ptWeights && dBinWidthPt)
 633         {
 634           dWpt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
 635         }
 636
 637         // determine v'(eta) weight:
 638         if(etaWeights && dBinWidthEta)
 639         {
 640           dWeta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
 641         }
 642
 643         // building up the weighted Q-vector:
 644         dQX += dWeight*dWphi*dWpt*dWeta*TMath::Cos(iOrder*dPhi);
 645         dQY += dWeight*dWphi*dWpt*dWeta*TMath::Sin(iOrder*dPhi);
 646
 647         // weighted multiplicity:
 648         sumOfWeights+=dWeight*dWphi*dWpt*dWeta;
 649
 650       } // end of if (pTrack->InRPSelection())
 651     } // loop over particles
 652
 653     Qarray[s].Set(dQX,dQY);
 654     Qarray[s].SetMult(sumOfWeights);
 655     Qarray[s].SetHarmonic(iOrder);
 656     Qarray[s].SetPOItype(AliFlowTrackSimple::kRP);
 657     Qarray[s].SetSubeventNumber(s);
 658
 659     //reset
 660     sumOfWeights = 0.;
 661     dQX = 0.;
 662     dQY = 0.;
 663   }
 664
 665 }
 666
 667
 668 //-----------------------------------------------------------------------
 669 void AliFlowEventSimple::Print(Option_t *option) const
 670 {
 671   //   -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
 672   //             ===============================================
 673   //   printf( "TH1.Print Name  = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
 674   printf( "Class.Print Name = %s, #tracks= %d, Number of RPs= %d, Number of POIs = %d, MC EventPlaneAngle= %f\n",
 675           GetName(),fNumberOfTracks, fNumberOfPOIs[0], fNumberOfPOIs[1], fMCReactionPlaneAngle );
 676
 677   TString optionstr(option);
 678   if (!optionstr.Contains("all")) return;
 679   if (fTrackCollection)
 680   {
 681     fTrackCollection->Print(option);
 682   }
 683   else
 684   {
 685     printf( "Empty track collection \n");
 686   }
 687 }
 688
 689 //-----------------------------------------------------------------------
 690 void AliFlowEventSimple::Browse(TBrowser *b)
 691 {
 692   //browse in TBrowser
 693   if (!b) return;
 694   if (!fNumberOfTracksWrap)
 695   {
 696     fNumberOfTracksWrap = new TParameter<int>("fNumberOfTracks", fNumberOfTracks);
 697     b->Add(fNumberOfTracksWrap);
 698   }
 699   if (!fNumberOfRPsWrap)
 700   {
 701     fNumberOfRPsWrap = new TParameter<int>("fNumberOfRPs", GetNumberOfRPs());
 702     b->Add(fNumberOfRPsWrap);
 703   }
 704   if (!fNumberOfPOIsWrap)
 705   {
 706     fNumberOfPOIsWrap = new TParameter<int>("fNumberOfPOIs", GetNumberOfPOIs());
 707     b->Add(fNumberOfPOIsWrap);
 708   }
 709   if (!fMCReactionPlaneAngleWrap)
 710   {
 711     fMCReactionPlaneAngleWrap = new TParameter<double>(" fMCReactionPlaneAngle",  fMCReactionPlaneAngle);
 712     b->Add( fMCReactionPlaneAngleWrap);
 713   }
 714   if (fTrackCollection) b->Add(fTrackCollection,"AliFlowTracksSimple");
 715 }
 716
 717 //-----------------------------------------------------------------------
 718 AliFlowEventSimple::AliFlowEventSimple( TTree* inputTree,
 719                                         const AliFlowTrackSimpleCuts* rpCuts,
 720                                         const AliFlowTrackSimpleCuts* poiCuts):
 721   fTrackCollection(NULL),
 722   fReferenceMultiplicity(0),
 723   fNumberOfTracks(0),
 724   fUseGlauberMCSymmetryPlanes(kFALSE),
 725   fUseExternalSymmetryPlanes(kFALSE),
 726   fPsi1(0.),
 727   fPsi2(0.),
 728   fPsi3(0.),
 729   fPsi4(0.),
 730   fPsi5(0.),
 731   fPsi1Psi3(0x0),
 732   fPsi2Psi4(0x0),
 733   fPsi3Psi5(0x0),
 734   fMCReactionPlaneAngle(0.),
 735   fMCReactionPlaneAngleIsSet(kFALSE),
 736   fAfterBurnerPrecision(0.001),
 737   fUserModified(kFALSE),
 738   fNumberOfTracksWrap(NULL),
 739   fNumberOfRPsWrap(NULL),
 740   fNumberOfPOIsWrap(NULL),
 741   fMCReactionPlaneAngleWrap(NULL),
 742   fShuffledIndexes(NULL),
 743   fShuffleTracks(kFALSE),
 744   fMothersCollection(new TObjArray()),
 745   fCentrality(-1.),
 746   fNumberOfPOItypes(2),
 747   fNumberOfPOIs(new Int_t[fNumberOfPOItypes])
 748 {
 749   //constructor, fills the event from a TTree of kinematic.root files
 750   //applies RP and POI cuts, tags the tracks
 751
 752   Int_t numberOfInputTracks = inputTree->GetEntries() ;
 753   fTrackCollection = new TObjArray(numberOfInputTracks/2);
 754
 755   TParticle* pParticle = new TParticle();
 756   inputTree->SetBranchAddress("Particles",&pParticle);
 757
 758   for (Int_t i=0; i<numberOfInputTracks; i++)
 759   {
 760     inputTree->GetEntry(i);   //get input particle
 761
 762     if (!pParticle) continue; //no particle
 763     if (!pParticle->IsPrimary()) continue;
 764
 765     Bool_t rpOK = (rpCuts->PassesCuts(pParticle)>0);
 766     Bool_t poiOK = poiCuts->PassesCuts(pParticle);
 767     Int_t poiType = poiCuts->GetPOItype();
 768
 769     if (rpOK || poiOK)
 770     {
 771       AliFlowTrackSimple* pTrack = new AliFlowTrackSimple(pParticle);
 772
 773       //marking the particles used for int. flow:
 774       if(rpOK)
 775       {
 776         pTrack->TagRP(kTRUE);
 777         IncrementNumberOfPOIs(0);
 778         cout<<"numberOfRPs = "<<fNumberOfPOIs[0]<<endl;
 779       }
 780       //marking the particles used for diff. flow:
 781       if(poiOK)
 782       {
 783         pTrack->Tag(poiType);
 784         IncrementNumberOfPOIs(poiType);
 785         printf("fNumberOfPOIs[%i] = %i",poiType,fNumberOfPOIs[poiType]);
 786       }
 787       //adding a particles which were used either for int. or diff. flow to the list
 788       AddTrack(pTrack);
 789     }
 790   }//for i
 791   delete pParticle;
 792 }
 793
 794 //_____________________________________________________________________________
 795 void AliFlowEventSimple::CloneTracks(Int_t n)
 796 {
 797   //clone every track n times to add non-flow
 798   if (n<=0) return; //no use to clone stuff zero or less times
 799   Int_t ntracks = fNumberOfTracks;
 800   fTrackCollection->Expand((n+1)*fNumberOfTracks);
 801   for (Int_t i=0; i<n; i++)
 802   {
 803     for (Int_t itrack=0; itrack<ntracks; itrack++)
 804     {
 805       AliFlowTrackSimple* track = dynamic_cast<AliFlowTrackSimple*>(fTrackCollection->At(itrack));
 806       if (!track) continue;
 807       AddTrack(static_cast<AliFlowTrackSimple*>(track->Clone()));
 808     }
 809   }
 810   SetUserModified();
 811 }
 812
 813 //_____________________________________________________________________________
 814 void AliFlowEventSimple::ResolutionPt(Double_t res)
 815 {
 816   //smear pt of all tracks by gaussian with sigma=res
 817   for (Int_t i=0; i<fNumberOfTracks; i++)
 818   {
 819     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 820     if (track) track->ResolutionPt(res);
 821   }
 822   SetUserModified();
 823 }
 824
 825 //_____________________________________________________________________________
 826 void AliFlowEventSimple::TagSubeventsInEta( Double_t etaMinA,
 827                                             Double_t etaMaxA,
 828                                             Double_t etaMinB,
 829                                             Double_t etaMaxB )
 830 {
 831   //Flag two subevents in given eta ranges
 832   for (Int_t i=0; i<fNumberOfTracks; i++)
 833   {
 834     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 835     if (!track) continue;
 836     track->ResetSubEventTags();
 837     Double_t eta=track->Eta();
 838     if (eta >= etaMinA && eta <= etaMaxA) track->SetForSubevent(0);
 839     if (eta >= etaMinB && eta <= etaMaxB) track->SetForSubevent(1);
 840   }
 841 }
 842
 843 //_____________________________________________________________________________
 844 void AliFlowEventSimple::TagSubeventsByCharge()
 845 {
 846   //Flag two subevents in given eta ranges
 847   for (Int_t i=0; i<fNumberOfTracks; i++)
 848   {
 849     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 850     if (!track) continue;
 851     track->ResetSubEventTags();
 852     Int_t charge=track->Charge();
 853     if (charge<0) track->SetForSubevent(0);
 854     if (charge>0) track->SetForSubevent(1);
 855   }
 856 }
 857
 858 //_____________________________________________________________________________
 859 void AliFlowEventSimple::AddV1( Double_t v1 )
 860 {
 861   //add v2 to all tracks wrt the reaction plane angle
 862   for (Int_t i=0; i<fNumberOfTracks; i++)
 863   {
 864     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 865     if (track) {
 866       if((fUseExternalSymmetryPlanes)||(fUseGlauberMCSymmetryPlanes))
 867         track->AddV1(v1, fPsi1, fAfterBurnerPrecision);
 868       else
 869         track->AddV1(v1, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 870     }
 871   }
 872   SetUserModified();
 873 }
 874
 875 //_____________________________________________________________________________
 876 void AliFlowEventSimple::AddV2( Double_t v2 )
 877 {
 878   //add v2 to all tracks wrt the reaction plane angle
 879   for (Int_t i=0; i<fNumberOfTracks; i++)
 880   {
 881     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 882     if (track) {
 883       if((fUseExternalSymmetryPlanes)||(fUseGlauberMCSymmetryPlanes))
 884         track->AddV2(v2, fPsi2, fAfterBurnerPrecision);
 885       else
 886         track->AddV2(v2, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 887     }
 888   }
 889   SetUserModified();
 890 }
 891
 892 //_____________________________________________________________________________
 893 void AliFlowEventSimple::AddV3( Double_t v3 )
 894 {
 895   //add v3 to all tracks wrt the reaction plane angle
 896   for (Int_t i=0; i<fNumberOfTracks; i++)
 897   {
 898     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 899     if(track) {
 900       if((fUseExternalSymmetryPlanes)||(fUseGlauberMCSymmetryPlanes))
 901         track->AddV3(v3, fPsi3, fAfterBurnerPrecision);
 902       else
 903         track->AddV3(v3, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 904     }
 905   }
 906   SetUserModified();
 907 }
 908
 909 //_____________________________________________________________________________
 910 void AliFlowEventSimple::AddV4( Double_t v4 )
 911 {
 912   //add v4 to all tracks wrt the reaction plane angle
 913   for (Int_t i=0; i<fNumberOfTracks; i++)
 914   {
 915     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 916     if(track) {
 917       if((fUseExternalSymmetryPlanes)||(fUseGlauberMCSymmetryPlanes))
 918         track->AddV4(v4, fPsi4, fAfterBurnerPrecision);
 919       else
 920         track->AddV4(v4, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 921     }
 922   }
 923   SetUserModified();
 924 }
 925
 926 //_____________________________________________________________________________
 927 void AliFlowEventSimple::AddV5( Double_t v5 )
 928 {
 929   //add v4 to all tracks wrt the reaction plane angle
 930   for (Int_t i=0; i<fNumberOfTracks; i++)
 931   {
 932     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 933     if(track) {
 934       if((fUseExternalSymmetryPlanes)||(fUseGlauberMCSymmetryPlanes))
 935         track->AddV5(v5, fPsi5, fAfterBurnerPrecision);
 936       else
 937         track->AddV5(v5, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 938     }
 939   }
 940   SetUserModified();
 941 }
 942
 943 //_____________________________________________________________________________
 944 void AliFlowEventSimple::AddFlow( Double_t v1, Double_t v2, Double_t v3, Double_t v4, Double_t v5,
 945                                   Double_t rp1, Double_t rp2, Double_t rp3, Double_t rp4, Double_t rp5 )
 946 {
 947   //add flow to all tracks wrt the reaction plane angle, for all harmonic separate angle
 948   for (Int_t i=0; i<fNumberOfTracks; i++)
 949   {
 950     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 951     if (track) track->AddFlow(v1,v2,v3,v4,v5,rp1,rp2,rp3,rp4,rp5,fAfterBurnerPrecision);
 952   }
 953   SetUserModified();
 954 }
 955
 956 //_____________________________________________________________________________
 957 void AliFlowEventSimple::AddFlow( Double_t v1, Double_t v2, Double_t v3, Double_t v4, Double_t v5 )
 958 {
 959   //add flow to all tracks wrt the reaction plane angle
 960   for (Int_t i=0; i<fNumberOfTracks; i++)
 961   {
 962     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 963     if (track) track->AddFlow(v1,v2,v3,v4,v5,fMCReactionPlaneAngle, fAfterBurnerPrecision);
 964   }
 965   SetUserModified();
 966 }
 967
 968 //_____________________________________________________________________________
 969 void AliFlowEventSimple::AddV2( TF1* ptDepV2 )
 970 {
 971   //add v2 to all tracks wrt the reaction plane angle
 972   for (Int_t i=0; i<fNumberOfTracks; i++)
 973   {
 974     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 975     if (!track) continue;
 976     Double_t v2 = ptDepV2->Eval(track->Pt());
 977     track->AddV2(v2, fMCReactionPlaneAngle, fAfterBurnerPrecision);
 978   }
 979   SetUserModified();
 980 }
 981
 982 //_____________________________________________________________________________
 983 void AliFlowEventSimple::TagRP( const AliFlowTrackSimpleCuts* cuts )
 984 {
 985   //tag tracks as reference particles (RPs)
 986   for (Int_t i=0; i<fNumberOfTracks; i++)
 987   {
 988     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
 989     if (!track) continue;
 990     Bool_t pass=cuts->PassesCuts(track);
 991     Bool_t rpTrack=track->InRPSelection();
 992     if (pass)
 993     {
 994       if (!rpTrack) fNumberOfPOIs[0]++; //only increase if not already tagged
 995     }
 996     else
 997     {
 998       if (rpTrack) fNumberOfPOIs[0]--; //only decrease if detagging
 999     }
1000     track->SetForRPSelection(pass);
1001   }
1002 }
1003
1004 //_____________________________________________________________________________
1005 void AliFlowEventSimple::TagPOI( const AliFlowTrackSimpleCuts* cuts, Int_t poiType )
1006 {
1007   //tag tracks as particles of interest (POIs)
1008   for (Int_t i=0; i<fNumberOfTracks; i++)
1009   {
1010     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
1011     if (!track) continue;
1012     Bool_t pass=cuts->PassesCuts(track);
1013     Bool_t poiTrack=track->InPOISelection();
1014     if (pass)
1015     {
1016       if (!poiTrack) fNumberOfPOIs[poiType]++; //only increase if not already tagged
1017     }
1018     else
1019     {
1020       if (poiTrack) fNumberOfPOIs[poiType]--; //only decrease if detagging
1021     }
1022     track->Tag(poiType,pass);
1023   }
1024 }
1025
1026 //_____________________________________________________________________________
1027 void AliFlowEventSimple::TagTracks( const AliFlowTrackSimpleCuts* cutsRP, const AliFlowTrackSimpleCuts* cutsPOI)
1028 {
1029     // simple interface to tagging poi's and rp's
1030     TagPOI(cutsRP, 0);
1031     TagPOI(cutsPOI, 1);
1032 }
1033 //_____________________________________________________________________________
1034 void AliFlowEventSimple::DefineDeadZone( Double_t etaMin,
1035                                          Double_t etaMax,
1036                                          Double_t phiMin,
1037                                          Double_t phiMax )
1038 {
1039   //mark tracks in given eta-phi region as dead
1040   //by resetting the flow bits
1041   for (Int_t i=0; i<fNumberOfTracks; i++)
1042   {
1043     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
1044     Double_t eta = track->Eta();
1045     Double_t phi = track->Phi();
1046     if (eta>etaMin && eta<etaMax && phi>phiMin && phi<phiMax)
1047     {
1048       if (track->InRPSelection()) {fNumberOfPOIs[0]--;}
1049       for (Int_t j=1; j<fNumberOfPOItypes; j++)
1050       {
1051         if (track->CheckTag(j)) {fNumberOfPOIs[j]--;}
1052       }
1053       track->ResetPOItype();
1054     }
1055   }
1056 }
1057
1058 //_____________________________________________________________________________
1059 Int_t AliFlowEventSimple::CleanUpDeadTracks()
1060 {
1061   //remove tracks that have no flow tags set and cleanup the container
1062   //returns number of cleaned tracks
1063   Int_t ncleaned=0;
1064   for (Int_t i=0; i<fNumberOfTracks; i++)
1065   {
1066     AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
1067     if (!track) continue;
1068     if (track->IsDead()) {delete track;track=NULL;ncleaned++;}
1069   }
1070   fTrackCollection->Compress(); //clean up empty slots
1071   fNumberOfTracks-=ncleaned; //update number of tracks
1072   delete [] fShuffledIndexes; fShuffledIndexes=NULL;
1073   return ncleaned;
1074 }
1075
1076 //_____________________________________________________________________________
1077 TF1* AliFlowEventSimple::SimplePtDepV2()
1078 {
1079   //return a standard pt dependent v2 formula, user has to clean up!
1080   return new TF1("StandardPtDepV2","((x<1.0)*(0.05/1.0)*x+(x>=1.0)*0.05)");
1081 }
1082
1083 //_____________________________________________________________________________
1084 TF1* AliFlowEventSimple::SimplePtSpectrum()
1085 {
1086   //return a standard pt spectrum, user has to clean up!
1087   return new TF1("StandardPtSpectrum","x*TMath::Exp(-pow(0.13957*0.13957+x*x,0.5)/0.4)",0.1,10.);
1088 }
1089
1090 //_____________________________________________________________________________
1091 void AliFlowEventSimple::ClearFast()
1092 {
1093   //clear the counters without deleting allocated objects so they can be reused
1094   fReferenceMultiplicity = 0;
1095   fNumberOfTracks = 0;
1096   for (Int_t i=0; i<fNumberOfPOItypes; i++)
1097   {
1098     fNumberOfPOIs[i] = 0;
1099   }
1100   fMCReactionPlaneAngle = 0.0;
1101   fMCReactionPlaneAngleIsSet = kFALSE;
1102   fAfterBurnerPrecision = 0.001;
1103   fUserModified = kFALSE;
1104   delete [] fShuffledIndexes; fShuffledIndexes=NULL;
1105 }