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