coverity + some patch for pp
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowEventSimple.cxx
CommitLineData
f1d945a1 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
929098e4 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
f1d945a1 26#include "Riostream.h"
27#include "TObjArray.h"
26c4cbb9 28#include "TFile.h"
03a02aca 29#include "TList.h"
929098e4 30#include "TTree.h"
31#include "TParticle.h"
f1d945a1 32#include "TMath.h"
26c4cbb9 33#include "TH1F.h"
34#include "TH1D.h"
bc231a12 35#include "TF1.h"
26c4cbb9 36#include "TProfile.h"
929098e4 37#include "TParameter.h"
c076fda8 38#include "TBrowser.h"
f1d945a1 39#include "AliFlowVector.h"
40#include "AliFlowTrackSimple.h"
929098e4 41#include "AliFlowTrackSimpleCuts.h"
701f71c1 42#include "AliFlowEventSimple.h"
7382279b 43#include "TRandom.h"
f1d945a1 44
f1d945a1 45ClassImp(AliFlowEventSimple)
46
47//-----------------------------------------------------------------------
46bec39c 48AliFlowEventSimple::AliFlowEventSimple():
49 fTrackCollection(NULL),
85d2ee8d 50 fReferenceMultiplicity(0),
46bec39c 51 fNumberOfTracks(0),
bc231a12 52 fNumberOfRPs(0),
7183fe85 53 fMCReactionPlaneAngle(0.),
929098e4 54 fMCReactionPlaneAngleIsSet(kFALSE),
bc231a12 55 fAfterBurnerPrecision(0.001),
9fa64edc 56 fUserModified(kFALSE),
c076fda8 57 fNumberOfTracksWrap(NULL),
bc231a12 58 fNumberOfRPsWrap(NULL),
a12990bb 59 fMCReactionPlaneAngleWrap(NULL)
46bec39c 60{
61 cout << "AliFlowEventSimple: Default constructor to be used only by root for io" << endl;
62}
63
64//-----------------------------------------------------------------------
9fa64edc 65AliFlowEventSimple::AliFlowEventSimple( Int_t n,
66 ConstructionMethod method,
bc231a12 67 TF1* ptDist,
68 Double_t phiMin,
69 Double_t phiMax,
70 Double_t etaMin,
71 Double_t etaMax):
9fa64edc 72 fTrackCollection(new TObjArray(n)),
a1c43d26 73 fReferenceMultiplicity(0),
bc231a12 74 fNumberOfTracks(0),
75 fNumberOfRPs(0),
76 fMCReactionPlaneAngle(0.),
77 fMCReactionPlaneAngleIsSet(kFALSE),
78 fAfterBurnerPrecision(0.001),
9fa64edc 79 fUserModified(kFALSE),
bc231a12 80 fNumberOfTracksWrap(NULL),
81 fNumberOfRPsWrap(NULL),
82 fMCReactionPlaneAngleWrap(NULL)
83{
9fa64edc 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);
f1d945a1 91}
92
93//-----------------------------------------------------------------------
e35ddff0 94AliFlowEventSimple::AliFlowEventSimple(const AliFlowEventSimple& anEvent):
bc231a12 95 TObject(anEvent),
96 fTrackCollection((TObjArray*)(anEvent.fTrackCollection)->Clone()),
85d2ee8d 97 fReferenceMultiplicity(anEvent.fReferenceMultiplicity),
e35ddff0 98 fNumberOfTracks(anEvent.fNumberOfTracks),
bc231a12 99 fNumberOfRPs(anEvent.fNumberOfRPs),
a12990bb 100 fMCReactionPlaneAngle(anEvent.fMCReactionPlaneAngle),
929098e4 101 fMCReactionPlaneAngleIsSet(anEvent.fMCReactionPlaneAngleIsSet),
bc231a12 102 fAfterBurnerPrecision(anEvent.fAfterBurnerPrecision),
9fa64edc 103 fUserModified(anEvent.fUserModified),
c076fda8 104 fNumberOfTracksWrap(anEvent.fNumberOfTracksWrap),
bc231a12 105 fNumberOfRPsWrap(anEvent.fNumberOfRPsWrap),
a12990bb 106 fMCReactionPlaneAngleWrap(anEvent.fMCReactionPlaneAngleWrap)
f1d945a1 107{
929098e4 108 //copy constructor
f1d945a1 109}
110
111//-----------------------------------------------------------------------
e35ddff0 112AliFlowEventSimple& AliFlowEventSimple::operator=(const AliFlowEventSimple& anEvent)
f1d945a1 113{
124fb262 114 //assignment operator
2047680a 115 if (&anEvent==this) return *this; //check self-assignment
701f71c1 116 if (fTrackCollection) fTrackCollection->Delete();
124fb262 117 delete fTrackCollection;
bc231a12 118 fTrackCollection = (TObjArray*)(anEvent.fTrackCollection)->Clone(); //deep copy
85d2ee8d 119 fReferenceMultiplicity = anEvent.fReferenceMultiplicity;
e35ddff0 120 fNumberOfTracks = anEvent.fNumberOfTracks;
bc231a12 121 fNumberOfRPs = anEvent.fNumberOfRPs;
a12990bb 122 fMCReactionPlaneAngle = anEvent.fMCReactionPlaneAngle;
929098e4 123 fMCReactionPlaneAngleIsSet = anEvent.fMCReactionPlaneAngleIsSet;
bc231a12 124 fAfterBurnerPrecision = anEvent.fAfterBurnerPrecision;
9fa64edc 125 fUserModified=anEvent.fUserModified;
929098e4 126 fNumberOfTracksWrap = anEvent.fNumberOfTracksWrap;
bc231a12 127 fNumberOfRPsWrap = anEvent.fNumberOfRPsWrap;
a12990bb 128 fMCReactionPlaneAngleWrap=anEvent.fMCReactionPlaneAngleWrap;
f1d945a1 129 return *this;
f1d945a1 130}
131
929098e4 132//-----------------------------------------------------------------------
f1d945a1 133AliFlowEventSimple::~AliFlowEventSimple()
134{
135 //destructor
929098e4 136 if (fTrackCollection) fTrackCollection->Delete();
137 delete fTrackCollection;
7d27a354 138 delete fNumberOfTracksWrap;
139 delete fNumberOfRPsWrap;
140 delete fMCReactionPlaneAngleWrap;
f1d945a1 141}
142
929098e4 143//-----------------------------------------------------------------------
bc231a12 144void AliFlowEventSimple::Generate(Int_t nParticles,
6e027e29 145 TF1* ptDist,
146 Double_t phiMin,
147 Double_t phiMax,
148 Double_t etaMin,
149 Double_t etaMax)
bc231a12 150{
151 //generate nParticles random tracks uniform in phi and eta
152 //according to the specified pt distribution
6e027e29 153 if (!ptDist)
154 {
155 static TF1 ptdistribution("ptSpectra","x*TMath::Exp(-pow(0.13957*0.13957+x*x,0.5)/0.4)",0.1,10.);
156 ptDist=&ptdistribution;
157 }
9fa64edc 158
bc231a12 159 for (Int_t i=0; i<nParticles; i++)
160 {
701f71c1 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);
bc231a12 167 }
6e027e29 168 fMCReactionPlaneAngle=gRandom->Uniform(0.0,TMath::TwoPi());
9fa64edc 169 fMCReactionPlaneAngleIsSet=kTRUE;
170 SetUserModified();
bc231a12 171}
172
173//-----------------------------------------------------------------------
f1d945a1 174AliFlowTrackSimple* AliFlowEventSimple::GetTrack(Int_t i)
175{
176 //get track i from collection
7382279b 177 if (i>=fNumberOfTracks) return NULL;
124fb262 178 AliFlowTrackSimple* pTrack = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i)) ;
e35ddff0 179 return pTrack;
f1d945a1 180}
181
929098e4 182//-----------------------------------------------------------------------
7382279b 183void AliFlowEventSimple::AddTrack( AliFlowTrackSimple* track )
184{
7d27a354 185 //add a track, delete the old one if necessary
186 if (fNumberOfTracks < fTrackCollection->GetEntriesFast())
187 {
188 TObject* o = fTrackCollection->At(fNumberOfTracks);
189 if (o) delete o;
190 }
191 fTrackCollection->AddAtAndExpand(track,fNumberOfTracks);
7382279b 192 fNumberOfTracks++;
193}
194
195//-----------------------------------------------------------------------
929098e4 196AliFlowVector AliFlowEventSimple::GetQ(Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
f1d945a1 197{
929098e4 198 // calculate Q-vector in harmonic n without weights (default harmonic n=2)
e35ddff0 199 Double_t dQX = 0.;
200 Double_t dQY = 0.;
201 AliFlowVector vQ;
202 vQ.Set(0.,0.);
929098e4 203
9825d4a9 204 Int_t iOrder = n;
8eca5d19 205 Double_t sumOfWeights = 0.;
44e060e0 206 Double_t dPhi = 0.;
207 Double_t dPt = 0.;
208 Double_t dEta = 0.;
209 Double_t dWeight = 1.;
929098e4 210
26c4cbb9 211 AliFlowTrackSimple* pTrack = NULL;
929098e4 212
44e060e0 213 Int_t nBinsPhi = 0;
214 Double_t dBinWidthPt = 0.;
215 Double_t dPtMin = 0.;
216 Double_t dBinWidthEta = 0.;
217 Double_t dEtaMin = 0.;
929098e4 218
44e060e0 219 Double_t wPhi = 1.; // weight Phi
220 Double_t wPt = 1.; // weight Pt
221 Double_t wEta = 1.; // weight Eta
929098e4 222
03a02aca 223 TH1F *phiWeights = NULL;
ae733b3b 224 TH1D *ptWeights = NULL;
03a02aca 225 TH1D *etaWeights = NULL;
929098e4 226
03a02aca 227 if(weightsList)
26c4cbb9 228 {
929098e4 229 if(usePhiWeights)
230 {
231 phiWeights = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights"));
232 if(phiWeights) nBinsPhi = phiWeights->GetNbinsX();
233 }
234 if(usePtWeights)
03a02aca 235 {
929098e4 236 ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
237 if(ptWeights)
238 {
239 dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
240 dPtMin = (ptWeights->GetXaxis())->GetXmin();
241 }
242 }
243 if(useEtaWeights)
03a02aca 244 {
929098e4 245 etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
246 if(etaWeights)
247 {
248 dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
249 dEtaMin = (etaWeights->GetXaxis())->GetXmin();
250 }
251 }
03a02aca 252 } // end of if(weightsList)
929098e4 253
254 // loop over tracks
255 for(Int_t i=0; i<fNumberOfTracks; i++)
26c4cbb9 256 {
124fb262 257 pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
929098e4 258 if(pTrack)
f1d945a1 259 {
929098e4 260 if(pTrack->InRPSelection())
261 {
262 dPhi = pTrack->Phi();
263 dPt = pTrack->Pt();
264 dEta = pTrack->Eta();
701f71c1 265 dWeight = pTrack->Weight();
929098e4 266
267 // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
268 if(phiWeights && nBinsPhi)
269 {
270 wPhi = phiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*nBinsPhi/TMath::TwoPi())));
271 }
272 // determine v'(pt) weight:
273 if(ptWeights && dBinWidthPt)
274 {
275 wPt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
276 }
277 // determine v'(eta) weight:
278 if(etaWeights && dBinWidthEta)
279 {
280 wEta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
281 }
282
283 // building up the weighted Q-vector:
44e060e0 284 dQX += dWeight*wPhi*wPt*wEta*TMath::Cos(iOrder*dPhi);
285 dQY += dWeight*wPhi*wPt*wEta*TMath::Sin(iOrder*dPhi);
929098e4 286
287 // weighted multiplicity:
8eca5d19 288 sumOfWeights += dWeight*wPhi*wPt*wEta;
929098e4 289
290 } // end of if (pTrack->InRPSelection())
291 } // end of if (pTrack)
292 else
293 {
294 cerr << "no particle!!!"<<endl;
295 }
ae733b3b 296 } // loop over particles
929098e4 297
e35ddff0 298 vQ.Set(dQX,dQY);
8eca5d19 299 vQ.SetMult(sumOfWeights);
929098e4 300
e35ddff0 301 return vQ;
929098e4 302
5fef318d 303}
304
929098e4 305//-----------------------------------------------------------------------
306void AliFlowEventSimple::Get2Qsub(AliFlowVector* Qarray, Int_t n, TList *weightsList, Bool_t usePhiWeights, Bool_t usePtWeights, Bool_t useEtaWeights)
395fadba 307{
929098e4 308
309 // calculate Q-vector in harmonic n without weights (default harmonic n=2)
395fadba 310 Double_t dQX = 0.;
311 Double_t dQY = 0.;
929098e4 312
395fadba 313 Int_t iOrder = n;
8eca5d19 314 Double_t sumOfWeights = 0.;
395fadba 315 Double_t dPhi = 0.;
316 Double_t dPt = 0.;
317 Double_t dEta = 0.;
44e060e0 318 Double_t dWeight = 1.;
929098e4 319
395fadba 320 AliFlowTrackSimple* pTrack = NULL;
929098e4 321
44e060e0 322 Int_t iNbinsPhiSub0 = 0;
323 Int_t iNbinsPhiSub1 = 0;
395fadba 324 Double_t dBinWidthPt = 0.;
325 Double_t dPtMin = 0.;
326 Double_t dBinWidthEta= 0.;
327 Double_t dEtaMin = 0.;
929098e4 328
329 Double_t dWphi = 1.; // weight Phi
330 Double_t dWpt = 1.; // weight Pt
331 Double_t dWeta = 1.; // weight Eta
332
44e060e0 333 TH1F* phiWeightsSub0 = NULL;
334 TH1F* phiWeightsSub1 = NULL;
29195b69 335 TH1D* ptWeights = NULL;
336 TH1D* etaWeights = NULL;
929098e4 337
338 if(weightsList)
339 {
340 if(usePhiWeights)
341 {
44e060e0 342 phiWeightsSub0 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub0"));
44e060e0 343 if(phiWeightsSub0) {
13ff9ccd 344 iNbinsPhiSub0 = phiWeightsSub0->GetNbinsX();
44e060e0 345 }
13ff9ccd 346 phiWeightsSub1 = dynamic_cast<TH1F *>(weightsList->FindObject("phi_weights_sub1"));
44e060e0 347 if(phiWeightsSub1) {
13ff9ccd 348 iNbinsPhiSub1 = phiWeightsSub1->GetNbinsX();
29195b69 349 }
929098e4 350 }
351 if(usePtWeights)
352 {
29195b69 353 ptWeights = dynamic_cast<TH1D *>(weightsList->FindObject("pt_weights"));
929098e4 354 if(ptWeights)
355 {
356 dBinWidthPt = ptWeights->GetBinWidth(1); // assuming that all bins have the same width
357 dPtMin = (ptWeights->GetXaxis())->GetXmin();
358 }
359 }
360 if(useEtaWeights)
361 {
362 etaWeights = dynamic_cast<TH1D *>(weightsList->FindObject("eta_weights"));
363 if(etaWeights)
364 {
365 dBinWidthEta = etaWeights->GetBinWidth(1); // assuming that all bins have the same width
366 dEtaMin = (etaWeights->GetXaxis())->GetXmin();
367 }
368 }
395fadba 369 } // end of if(weightsList)
929098e4 370
b125a454 371 //loop over the two subevents
929098e4 372 for (Int_t s=0; s<2; s++)
373 {
374 // loop over tracks
375 for(Int_t i=0; i<fNumberOfTracks; i++)
376 {
124fb262 377 pTrack = (AliFlowTrackSimple*)fTrackCollection->At(i);
929098e4 378 if(pTrack)
379 {
380 if(pTrack->InRPSelection())
381 {
382 if (pTrack->InSubevent(s))
383 {
44e060e0 384 dPhi = pTrack->Phi();
385 dPt = pTrack->Pt();
386 dEta = pTrack->Eta();
387 dWeight = pTrack->Weight();
929098e4 388
389 // determine Phi weight: (to be improved, I should here only access it + the treatment of gaps in the if statement)
13ff9ccd 390 //subevent 0
391 if(s == 0) {
392 if(phiWeightsSub0 && iNbinsPhiSub0) {
393 Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub0/TMath::TwoPi()));
394 //use the phi value at the center of the bin
395 dPhi = phiWeightsSub0->GetBinCenter(phiBin);
396 dWphi = phiWeightsSub0->GetBinContent(phiBin);
397 }
398 }
399 //subevent 1
400 else if (s == 1) {
401 if(phiWeightsSub1 && iNbinsPhiSub1) {
402 Int_t phiBin = 1+(Int_t)(TMath::Floor(dPhi*iNbinsPhiSub1/TMath::TwoPi()));
403 //use the phi value at the center of the bin
404 dPhi = phiWeightsSub1->GetBinCenter(phiBin);
405 dWphi = phiWeightsSub1->GetBinContent(phiBin);
44e060e0 406 }
407 }
13ff9ccd 408
929098e4 409 // determine v'(pt) weight:
410 if(ptWeights && dBinWidthPt)
411 {
412 dWpt=ptWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-dPtMin)/dBinWidthPt)));
413 }
44e060e0 414
929098e4 415 // determine v'(eta) weight:
416 if(etaWeights && dBinWidthEta)
417 {
418 dWeta=etaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-dEtaMin)/dBinWidthEta)));
419 }
420
421 // building up the weighted Q-vector:
44e060e0 422 dQX += dWeight*dWphi*dWpt*dWeta*TMath::Cos(iOrder*dPhi);
423 dQY += dWeight*dWphi*dWpt*dWeta*TMath::Sin(iOrder*dPhi);
929098e4 424
425 // weighted multiplicity:
8eca5d19 426 sumOfWeights+=dWeight*dWphi*dWpt*dWeta;
929098e4 427
428 } // end of subevent
429 } // end of if (pTrack->InRPSelection())
29195b69 430 } // end of if (pTrack)
929098e4 431 else
432 {
433 cerr << "no particle!!!"<<endl;
434 }
29195b69 435 } // loop over particles
436 Qarray[s].Set(dQX,dQY);
8eca5d19 437 Qarray[s].SetMult(sumOfWeights);
29195b69 438 //reset
8eca5d19 439 sumOfWeights = 0.;
29195b69 440 dQX = 0.;
441 dQY = 0.;
442 }
929098e4 443
395fadba 444}
445
446
929098e4 447//-----------------------------------------------------------------------
c076fda8 448void AliFlowEventSimple::Print(Option_t *option) const
449{
450 // -*-*-*-*-*Print some global quantities for this histogram collection class *-*-*-*-*-*-*-*
451 // ===============================================
452 // printf( "TH1.Print Name = %s, Entries= %d, Total sum= %g\n",GetName(),Int_t(fEntries),GetSumOfWeights());
f873707d 453 printf( "Class.Print Name = %s, #tracks= %d, Number of RPs= %d, MC EventPlaneAngle= %f\n",
bc231a12 454 GetName(),fNumberOfTracks, fNumberOfRPs, fMCReactionPlaneAngle );
c076fda8 455
b4dba88d 456 TString optionstr(option);
457 if (!optionstr.Contains("all")) return;
929098e4 458 if (fTrackCollection)
459 {
c076fda8 460 fTrackCollection->Print(option);
461 }
929098e4 462 else
463 {
c076fda8 464 printf( "Empty track collection \n");
465 }
466}
467
929098e4 468//-----------------------------------------------------------------------
469void AliFlowEventSimple::Browse(TBrowser *b)
c076fda8 470{
cc0afcfc 471 //browse in TBrowser
c076fda8 472 if (!b) return;
929098e4 473 if (!fNumberOfTracksWrap)
474 {
c076fda8 475 fNumberOfTracksWrap = new TParameter<int>("fNumberOfTracks", fNumberOfTracks);
476 b->Add(fNumberOfTracksWrap);
477 }
bc231a12 478 if (!fNumberOfRPsWrap)
929098e4 479 {
bc231a12 480 fNumberOfRPsWrap = new TParameter<int>("fNumberOfRPs", fNumberOfRPs);
481 b->Add(fNumberOfRPsWrap);
c076fda8 482 }
929098e4 483 if (!fMCReactionPlaneAngleWrap)
484 {
7183fe85 485 fMCReactionPlaneAngleWrap = new TParameter<double>(" fMCReactionPlaneAngle", fMCReactionPlaneAngle);
a12990bb 486 b->Add( fMCReactionPlaneAngleWrap);
487 }
c076fda8 488 if (fTrackCollection) b->Add(fTrackCollection,"AliFlowTracksSimple");
489}
490
929098e4 491//-----------------------------------------------------------------------
492AliFlowEventSimple::AliFlowEventSimple( TTree* inputTree,
493 const AliFlowTrackSimpleCuts* rpCuts,
494 const AliFlowTrackSimpleCuts* poiCuts):
495 fTrackCollection(NULL),
85d2ee8d 496 fReferenceMultiplicity(0),
929098e4 497 fNumberOfTracks(0),
bc231a12 498 fNumberOfRPs(0),
929098e4 499 fMCReactionPlaneAngle(0.),
500 fMCReactionPlaneAngleIsSet(kFALSE),
bc231a12 501 fAfterBurnerPrecision(0.001),
9fa64edc 502 fUserModified(kFALSE),
929098e4 503 fNumberOfTracksWrap(NULL),
bc231a12 504 fNumberOfRPsWrap(NULL),
929098e4 505 fMCReactionPlaneAngleWrap(NULL)
506{
507 //constructor, fills the event from a TTree of kinematic.root files
508 //applies RP and POI cuts, tags the tracks
509
510 Int_t numberOfInputTracks = inputTree->GetEntries() ;
511 fTrackCollection = new TObjArray(numberOfInputTracks/2);
512
513 TParticle* pParticle = new TParticle();
514 inputTree->SetBranchAddress("Particles",&pParticle);
515
516 Int_t iSelParticlesPOI = 0;
517
518 for (Int_t i=0; i<numberOfInputTracks; i++)
519 {
520 inputTree->GetEntry(i); //get input particle
521
522 if (!pParticle) continue; //no particle
523 if (!pParticle->IsPrimary()) continue;
524
7382279b 525 Bool_t rpOK = rpCuts->PassesCuts(pParticle);
526 Bool_t poiOK = poiCuts->PassesCuts(pParticle);
929098e4 527
7382279b 528 if (rpOK || poiOK)
929098e4 529 {
530 AliFlowTrackSimple* pTrack = new AliFlowTrackSimple(pParticle);
531
532 //marking the particles used for int. flow:
7382279b 533 if(rpOK)
929098e4 534 {
535 pTrack->SetForRPSelection(kTRUE);
bc231a12 536 fNumberOfRPs++;
929098e4 537 }
538 //marking the particles used for diff. flow:
7382279b 539 if(poiOK)
929098e4 540 {
541 pTrack->SetForPOISelection(kTRUE);
542 iSelParticlesPOI++;
543 }
544 //adding a particles which were used either for int. or diff. flow to the list
124fb262 545 AddTrack(pTrack);
929098e4 546 }
547 }//for i
548 delete pParticle;
549}
550
551//_____________________________________________________________________________
552void AliFlowEventSimple::CloneTracks(Int_t n)
553{
554 //clone every track n times to add non-flow
34b15925 555 if (n<=0) return; //no use to clone stuff zero or less times
bc231a12 556 Int_t ntracks = fNumberOfTracks;
557 fTrackCollection->Expand((n+1)*fNumberOfTracks);
558 for (Int_t i=0; i<n; i++)
929098e4 559 {
bc231a12 560 for (Int_t itrack=0; itrack<ntracks; itrack++)
929098e4 561 {
562 AliFlowTrackSimple* track = dynamic_cast<AliFlowTrackSimple*>(fTrackCollection->At(itrack));
563 if (!track) continue;
bc231a12 564 AddTrack(static_cast<AliFlowTrackSimple*>(track->Clone()));
929098e4 565 }
566 }
9fa64edc 567 SetUserModified();
929098e4 568}
7382279b 569
570//_____________________________________________________________________________
571void AliFlowEventSimple::ResolutionPt(Double_t res)
572{
573 //smear pt of all tracks by gaussian with sigma=res
574 for (Int_t i=0; i<fNumberOfTracks; i++)
575 {
576 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
124fb262 577 if (track) track->ResolutionPt(res);
7382279b 578 }
9fa64edc 579 SetUserModified();
7382279b 580}
124fb262 581
582//_____________________________________________________________________________
34b15925 583void AliFlowEventSimple::TagSubeventsInEta( Double_t etaMinA,
584 Double_t etaMaxA,
585 Double_t etaMinB,
586 Double_t etaMaxB )
124fb262 587{
588 //Flag two subevents in given eta ranges
589 for (Int_t i=0; i<fNumberOfTracks; i++)
590 {
591 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
3ca4688a 592 if (!track) continue;
593 track->ResetSubEventTags();
124fb262 594 Double_t eta=track->Eta();
595 if (eta >= etaMinA && eta <= etaMaxA) track->SetForSubevent(0);
596 if (eta >= etaMinB && eta <= etaMaxB) track->SetForSubevent(1);
597 }
598}
599
600//_____________________________________________________________________________
076df7bf 601void AliFlowEventSimple::TagSubeventsByCharge()
602{
603 //Flag two subevents in given eta ranges
604 for (Int_t i=0; i<fNumberOfTracks; i++)
605 {
606 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
3ca4688a 607 if (!track) continue;
608 track->ResetSubEventTags();
076df7bf 609 Int_t charge=track->Charge();
610 if (charge<0) track->SetForSubevent(0);
611 if (charge>0) track->SetForSubevent(1);
612 }
613}
614
615//_____________________________________________________________________________
34b15925 616void AliFlowEventSimple::AddV1( Double_t v1 )
bc231a12 617{
618 //add v2 to all tracks wrt the reaction plane angle
619 for (Int_t i=0; i<fNumberOfTracks; i++)
620 {
621 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
622 if (track) track->AddV1(v1, fMCReactionPlaneAngle, fAfterBurnerPrecision);
623 }
9fa64edc 624 SetUserModified();
bc231a12 625}
626
627//_____________________________________________________________________________
34b15925 628void AliFlowEventSimple::AddV2( Double_t v2 )
124fb262 629{
244c607a 630 //add v2 to all tracks wrt the reaction plane angle
124fb262 631 for (Int_t i=0; i<fNumberOfTracks; i++)
632 {
633 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
bc231a12 634 if (track) track->AddV2(v2, fMCReactionPlaneAngle, fAfterBurnerPrecision);
635 }
9fa64edc 636 SetUserModified();
bc231a12 637}
638
639//_____________________________________________________________________________
54089829 640void AliFlowEventSimple::AddV3( Double_t v3 )
641{
642 //add v3 to all tracks wrt the reaction plane angle
643 for (Int_t i=0; i<fNumberOfTracks; i++)
644 {
645 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
646 if (track) track->AddV3(v3, fMCReactionPlaneAngle, fAfterBurnerPrecision);
647 }
648 SetUserModified();
649}
650
651//_____________________________________________________________________________
34b15925 652void AliFlowEventSimple::AddV4( Double_t v4 )
bc231a12 653{
654 //add v4 to all tracks wrt the reaction plane angle
655 for (Int_t i=0; i<fNumberOfTracks; i++)
656 {
657 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
658 if (track) track->AddV4(v4, fMCReactionPlaneAngle, fAfterBurnerPrecision);
659 }
9fa64edc 660 SetUserModified();
bc231a12 661}
662
663//_____________________________________________________________________________
1a80f9f6 664void AliFlowEventSimple::AddV5( Double_t v5 )
665{
666 //add v4 to all tracks wrt the reaction plane angle
667 for (Int_t i=0; i<fNumberOfTracks; i++)
668 {
669 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
670 if (track) track->AddV5(v5, fMCReactionPlaneAngle, fAfterBurnerPrecision);
671 }
672 SetUserModified();
673}
674
675//_____________________________________________________________________________
676void AliFlowEventSimple::AddFlow( Double_t v1, Double_t v2, Double_t v3, Double_t v4, Double_t v5,
677 Double_t rp1, Double_t rp2, Double_t rp3, Double_t rp4, Double_t rp5 )
678{
679 //add flow to all tracks wrt the reaction plane angle, for all harmonic separate angle
680 for (Int_t i=0; i<fNumberOfTracks; i++)
681 {
682 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
683 if (track) track->AddFlow(v1,v2,v3,v4,v5,rp1,rp2,rp3,rp4,rp5,fAfterBurnerPrecision);
684 }
685 SetUserModified();
686}
687
688//_____________________________________________________________________________
689void AliFlowEventSimple::AddFlow( Double_t v1, Double_t v2, Double_t v3, Double_t v4, Double_t v5 )
bc231a12 690{
691 //add flow to all tracks wrt the reaction plane angle
692 for (Int_t i=0; i<fNumberOfTracks; i++)
693 {
694 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
1a80f9f6 695 if (track) track->AddFlow(v1,v2,v3,v4,v5,fMCReactionPlaneAngle, fAfterBurnerPrecision);
124fb262 696 }
9fa64edc 697 SetUserModified();
124fb262 698}
699
dd91b595 700//_____________________________________________________________________________
99ff691b 701void AliFlowEventSimple::AddV2( TF1* ptDepV2 )
702{
703 //add v2 to all tracks wrt the reaction plane angle
704 for (Int_t i=0; i<fNumberOfTracks; i++)
705 {
706 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
2a745a5f 707 if (!track) continue;
99ff691b 708 Double_t v2 = ptDepV2->Eval(track->Pt());
2a745a5f 709 track->AddV2(v2, fMCReactionPlaneAngle, fAfterBurnerPrecision);
99ff691b 710 }
711 SetUserModified();
712}
713
714//_____________________________________________________________________________
cc0afcfc 715void AliFlowEventSimple::TagRP( const AliFlowTrackSimpleCuts* cuts )
dd91b595 716{
717 //tag tracks as reference particles (RPs)
718 for (Int_t i=0; i<fNumberOfTracks; i++)
719 {
720 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
721 if (!track) continue;
701f71c1 722 Bool_t pass=cuts->PassesCuts(track);
e0da66a2 723 Bool_t rpTrack=track->InRPSelection();
701f71c1 724 if (pass)
e0da66a2 725 {
726 if (!rpTrack) fNumberOfRPs++; //only increase if not already tagged
727 }
701f71c1 728 else
e0da66a2 729 {
730 if (rpTrack) fNumberOfRPs--; //only decrease if detagging
731 }
732 track->SetForRPSelection(pass);
dd91b595 733 }
734}
735
736//_____________________________________________________________________________
cc0afcfc 737void AliFlowEventSimple::TagPOI( const AliFlowTrackSimpleCuts* cuts )
dd91b595 738{
739 //tag tracks as particles of interest (POIs)
740 for (Int_t i=0; i<fNumberOfTracks; i++)
741 {
742 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
743 if (!track) continue;
701f71c1 744 Bool_t pass=cuts->PassesCuts(track);
745 track->SetForPOISelection(pass);
dd91b595 746 }
747}
748
34b15925 749//_____________________________________________________________________________
750void AliFlowEventSimple::DefineDeadZone( Double_t etaMin,
751 Double_t etaMax,
752 Double_t phiMin,
753 Double_t phiMax )
754{
755 //mark tracks in given eta-phi region as dead
756 //by resetting the flow bits
757 for (Int_t i=0; i<fNumberOfTracks; i++)
758 {
759 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
760 Double_t eta = track->Eta();
761 Double_t phi = track->Phi();
762 if (eta>etaMin && eta<etaMax && phi>phiMin && phi<phiMax)
1a43aca9 763 {
764 if (track->InRPSelection()) fNumberOfRPs--;
34b15925 765 track->ResetFlowTags();
1a43aca9 766 }
34b15925 767 }
768}
769
770//_____________________________________________________________________________
771Int_t AliFlowEventSimple::CleanUpDeadTracks()
772{
773 //remove tracks that have no flow tags set and cleanup the container
65ef7d1a 774 //returns number of cleaned tracks
775 Int_t ncleaned=0;
34b15925 776 for (Int_t i=0; i<fNumberOfTracks; i++)
777 {
778 AliFlowTrackSimple* track = static_cast<AliFlowTrackSimple*>(fTrackCollection->At(i));
cd62a2a8 779 if (!track) continue;
7d27a354 780 if (track->IsDead()) {delete track;track=NULL;ncleaned++;}
34b15925 781 }
782 fTrackCollection->Compress(); //clean up empty slots
7d27a354 783 fNumberOfTracks-=ncleaned; //update number of tracks
65ef7d1a 784 return ncleaned;
34b15925 785}
99ff691b 786
787//_____________________________________________________________________________
788TF1* AliFlowEventSimple::SimplePtDepV2()
789{
790 //return a standard pt dependent v2 formula, user has to clean up!
791 return new TF1("StandardPtDepV2","((x<1.0)*(0.05/1.0)*x+(x>=1.0)*0.05)");
792}
793
794//_____________________________________________________________________________
795TF1* AliFlowEventSimple::SimplePtSpectrum()
796{
797 //return a standard pt spectrum, user has to clean up!
798 return new TF1("StandardPtSpectrum","x*TMath::Exp(-pow(0.13957*0.13957+x*x,0.5)/0.4)",0.1,10.);
799}
7d27a354 800
801//_____________________________________________________________________________
802void AliFlowEventSimple::ClearFast()
803{
804 //clear the counter without deleting allocated objects so they can be reused
805 fReferenceMultiplicity = 0;
806 fNumberOfTracks = 0;
807 fNumberOfRPs = 0;
808 fMCReactionPlaneAngle = 0.0;
809 fMCReactionPlaneAngleIsSet = kFALSE;
810 fAfterBurnerPrecision = 0.001;
811 fUserModified = kFALSE;
812}