e041be941436f4f9a8e61caf42e309413984db25
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowMaker.cxx
1 //////////////////////////////////////////////////////////////////////
2 //
3 // $Id$
4 //
5 // Author: Emanuele Simili
6 //
7 //////////////////////////////////////////////////////////////////////
8 //_____________________________________________________________
9 //
10 // Description: 
11 //        AliFlowMaker provides the method to create AliFlowEvent(s) 
12 // from AliESD(s). Very basic track cuts are applyed.
13 // The present class can be used in a simple AliRoot macro or in a 
14 // more complex enviroment such as AliSelector or AliTask.
15 //
16 //////////////////////////////////////////////////////////////////////
17
18 #ifndef ALIFLOWMAKER_CXX
19 #define ALIFLOWMAKER_CXX
20
21 // ROOT things
22 #include <TROOT.h>
23 #include <TFile.h>
24 #include <TString.h>
25 #include <TMath.h>
26 #include "TClonesArray.h"
27
28 // AliRoot things
29 #include "AliESD.h"
30 #include "AliESDVertex.h"
31 #include "AliESDtrack.h"
32 #include "AliESDv0.h"
33 //#include "AliKalmanTrack.h"
34
35 // Flow things
36 #include "AliFlowEvent.h"
37 #include "AliFlowTrack.h"
38 #include "AliFlowV0.h"
39 #include "AliFlowConstants.h"
40 #include "AliFlowMaker.h"
41
42 // ANSI things
43 #include <stdlib.h>
44 using namespace std; //required for resolving the 'cout' symbol
45
46 ClassImp(AliFlowMaker) 
47 //-----------------------------------------------------------------------
48 AliFlowMaker::AliFlowMaker():
49   fEventNumber(0), fTrackNumber(0), fV0Number(0), fGoodTracks(0), fGoodV0s(0),
50   fGoodTracksEta(0), fPosiTracks(0), fNegaTracks(0), fUnconstrained(0),
51   fSumAll(0), fCutEvts(0), fCutTrks(0), fCutV0s(0), fNewAli(kFALSE), 
52   fLoopTrks(kTRUE), fLoopV0s(kTRUE), fCounter(0), 
53   fESD(0x0), fTrack(0x0), fV0(0x0), fVertex(0x0),
54   fRunID(0), fNumberOfEvents(0), fNumberOfTracks(0), 
55   fNumberOfV0s(0), fMagField(0), 
56   fFlowEvent(0x0), fFlowTrack(0x0), fFlowV0(0x0),
57   fNHits(0), fElow(0.001), fEup(1000.)
58 {
59  // default constructor 
60  // resets counters , sets defaults
61  
62  for(Int_t bb=0;bb<5;bb++) { fBayesianAll[bb] = 0 ; } ;
63
64  // trak cuts
65  fLabel[0] = 0 ; fLabel[1] = -1 ;
66 }
67 //-----------------------------------------------------------------------
68 AliFlowMaker::~AliFlowMaker()
69 {
70  // default destructor (no actions) 
71 }
72 //-----------------------------------------------------------------------
73
74 //-----------------------------------------------------------------------
75 AliFlowEvent* AliFlowMaker::FillFlowEvent(AliESD* fESD)
76 {
77  // From the AliESD (input) fills the AliFlowEvent (output) . 
78  // It loops on track & v0 and calls the methods to fill the arrays . 
79  // ... . 
80
81  fFlowEvent = new AliFlowEvent() ; if(!fFlowEvent) { return 0 ; }
82  //cout << " -evt- " << fFlowEvent << endl ;
83
84  fRunID = fESD->GetRunNumber() ;
85  fEventNumber = -1 ;
86  // fEventNumber = fESD->GetEventNumber() ;
87  fNumberOfTracks = fESD->GetNumberOfTracks() ;
88  fNumberOfV0s = fESD->GetNumberOfV0s() ;
89  //
90  cout << " *evt n. " << fEventNumber << " (run " << fRunID << ")  -  tracks: " << fNumberOfTracks << " ,   v0s " << fNumberOfV0s << endl ;
91
92  // Event id 
93  fFlowEvent->SetRunID(fRunID) ;                
94  fFlowEvent->SetEventID(fEventNumber) ;         
95  fFlowEvent->SetOrigMult((UInt_t)fNumberOfTracks) ;
96
97  // Run information (fixed - ???)
98  fMagField = fESD->GetMagneticField() ; // cout << " *fMagField " << fMagField << endl ;
99  fFlowEvent->SetMagneticField(fMagField) ;      
100  fFlowEvent->SetCenterOfMassEnergy(AliFlowConstants::fgCenterOfMassEnergy) ;  
101  fFlowEvent->SetBeamMassNumberEast(AliFlowConstants::fgBeamMassNumberEast) ;  
102  fFlowEvent->SetBeamMassNumberWest(AliFlowConstants::fgBeamMassNumberWest) ;  
103
104  // Trigger information (now is: ULon64_t - some trigger mask)
105  fFlowEvent->SetL0TriggerWord((Int_t)fESD->GetTriggerMask()); 
106  
107  // Get primary vertex position
108  fVertex = (AliESDVertex*)fESD->GetVertex() ;
109  Double_t position[3] ; 
110  fVertex->GetXYZ(position) ;
111  fFlowEvent->SetVertexPos((Float_t)position[0],(Float_t)position[1],(Float_t)position[2]) ; 
112
113  // Zero Degree Calorimeter information
114  Int_t zdcp = fESD->GetZDCParticipants() ; 
115  Float_t zdce[3] ; 
116  zdce[0] = fESD->GetZDCN1Energy() + fESD->GetZDCN2Energy(); 
117  zdce[1] = fESD->GetZDCP1Energy() + fESD->GetZDCP2Energy() ; 
118  zdce[2] = fESD->GetZDCEMEnergy() ;
119  fFlowEvent->SetZDCpart(zdcp);                          
120  fFlowEvent->SetZDCenergy(zdce[0],zdce[1],zdce[2]);     
121
122  // Track loop
123  if(fLoopTrks)
124  {
125   Int_t badTrks = 0 ;   
126   for(fTrackNumber=0;fTrackNumber<fNumberOfTracks;fTrackNumber++) 
127   {
128    fTrack = fESD->GetTrack(fTrackNumber) ;
129    if(CheckTrack(fTrack))
130    {
131     FillFlowTrack(fTrack) ;   
132     fGoodTracks++ ;                                   
133    }
134    else { badTrks++ ; continue ; }
135   }
136   fCutTrks += badTrks ;
137  } // cout << " -track number- :  " << fTrackNumber << endl ;
138
139  // V0 loop
140  if(fLoopV0s)
141  {
142   Int_t badV0s = 0 ;
143   for(fV0Number=0;fV0Number<fNumberOfV0s;fV0Number++) 
144   {
145    fV0 = fESD->GetV0(fV0Number) ;                       
146    if(CheckV0(fV0))
147    {
148     FillFlowV0(fV0) ;
149     fGoodV0s++ ;                        
150    }
151    else { badV0s++ ; continue ; }
152   }
153   fCutV0s += badV0s ;
154  } // cout << " -v0 number- :  " << fV0Number << endl ;
155
156  // Evt setting stuff
157  fFlowEvent->SetCentrality();   
158                                 
159  fCounter++ ; 
160  return fFlowEvent ;
161 }
162 //----------------------------------------------------------------------
163 AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
164 {
165  // From the AliESDtrack (input) fills the AliFlowTrack (output) .
166
167  TString name = "" ; name += fTrackNumber ;
168  Int_t idx = fFlowEvent->TrackCollection()->GetEntries() ;
169  fFlowTrack = (AliFlowTrack*)(fFlowEvent->TrackCollection()->New(idx)) ;
170  fFlowTrack->SetName(name.Data()) ;
171  
172  // cout << " -tr- " << name.Data() << "(" << idx << ")"  << endl ;
173
174  // ESD particle label (link: KineTree-ESD)
175  Int_t label = TMath::Abs(fTrack->GetLabel());
176  fFlowTrack->SetLabel(label) ;                  
177
178  // signed DCA from ESDtrack
179  Float_t xy = 0 ; Float_t z = 0 ; 
180  fTrack->GetImpactParameters(xy,z) ; 
181  fFlowTrack->SetDcaSigned(xy,z) ;                   
182
183  // error on the DCA
184  Float_t dcaBis[2] ; Float_t dcaCov[3] ; 
185  for(Int_t dd=0;dd<3;dd++) { dcaCov[dd] = 0. ; }
186  fTrack->GetImpactParameters(dcaBis, dcaCov) ;
187  fFlowTrack->SetDcaError(dcaCov[0],dcaCov[1],dcaCov[2]) ;                   
188
189  // UnConstrained (global) first
190  Double_t gD[3] ;                               
191  fTrack->GetPxPyPz(gD) ;                        
192  // -
193  Float_t phiGl = (Float_t)Phi(gD) ;  
194  if(phiGl<0) { phiGl += 2*TMath::Pi() ; }
195  fFlowTrack->SetPhiGlobal(phiGl) ;              
196  Float_t ptGl = (Float_t)Pt(gD) ;  if(ptGl<=0) { cout << " !!! ptGlobal = " << ptGl << endl ; }
197  fFlowTrack->SetPtGlobal(ptGl) ;                
198  Float_t etaGl = (Float_t)Eta(gD) ; 
199  fFlowTrack->SetEtaGlobal(etaGl) ;              
200
201  // Constrained (NEW)
202  Double_t cD[3] ;
203  Double_t par1 ; Double_t par2 ;  Double_t par3[3] ;
204  if(fTrack->GetConstrainedExternalParameters(par1,par2,par3))
205  {
206   fTrack->GetConstrainedPxPyPz(cD) ;     
207  }
208  else { for(Int_t iii=0;iii<3;iii++) { cD[iii] =0 ; } }
209
210  if(Norm(cD)!=0.)   // ConstrainedPxPyPz != 0 if ConstrainedChi2 < something ...
211  {                                                      
212   Float_t phi = (Float_t)Phi(cD) ; 
213   if(phi<0) { phi += 2*TMath::Pi() ; }
214   fFlowTrack->SetPhi(phi) ;                             
215   Float_t pt = (Float_t)Pt(cD) ;   if(pt<=0) { cout << " !!! pt = " << pt << endl ; }
216   fFlowTrack->SetPt(pt) ;                               
217   Float_t eta = (Float_t)Eta(cD) ; 
218   fFlowTrack->SetEta(eta) ;                             
219  
220   // number of constrainable tracks with |eta| < AliFlowConstants::fgEtaGood (0.9)
221   if(TMath::Abs(eta) < AliFlowConstants::fgEtaGood)  { fGoodTracksEta++ ; }
222  }
223  else  // in case Constriction impossible for track, fill the UnConstrained (global)
224  {
225   fUnconstrained++ ;    
226   fFlowTrack->SetPhi(0.) ;
227   fFlowTrack->SetPt(0.) ;   
228   fFlowTrack->SetEta(0.) ; 
229  }           
230
231  // positive - negative tracks
232  Int_t trkSign = (Int_t)fTrack->GetSign() ; 
233  fFlowTrack->SetCharge(trkSign) ;               
234  if(trkSign>0)  { fPosiTracks++ ; }
235  else if(trkSign<0)     { fNegaTracks++ ; }
236  else                   { return 0 ; }
237
238  // Tracking parameters (fit , TPC , ITS , dE/dx)
239  fFlowTrack->SetChi2(fTrack->GetConstrainedChi2()) ;            
240  fFlowTrack->SetTrackLength(fTrack->GetIntegratedLength()) ;    
241  // -
242  Int_t idXt[180] ; // used for Cluster Map ( see AliESDtrack::GetTPCclusters() )    // old:    Int
243  Int_t idX[6] ;    // used for Cluster Map ( see AliESDtrack::GetITSclusters() )    // old:    UInt
244  Int_t idxr[130] ; // used for Cluster Map ( see AliESDtrack::GetTRDclusters() )    // old:    UInt
245  Int_t nClus = 0 ;      
246  Int_t fNFound = 0 ;                                    // *!* fNFoundable (in AliTPCtrack) ... added by M.Ianov 
247  // -
248  Double_t detPid[5] ;           
249  Float_t  detPid6[AliFlowConstants::kPid] ; 
250
251  Double_t pVecAt[3] ; 
252  for(Int_t gg=0;gg<3;gg++) { pVecAt[gg] = gD[gg] ; }
253  Bool_t boh ; Float_t pAt = 0 ;                         // to get p at each detector
254  // -
255  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgTPCx, fMagField, pVecAt) ; }
256  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgTPCx, fMagField, pVecAt) ; }
257  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
258  nClus = fTrack->GetTPCclusters(idXt) ;
259  fNFound = fTrack->GetTPCNclsF() ;  // was 160
260  if( (fTrack->GetStatus() & AliESDtrack::kTPCpid) != 0 )
261  {
262   fTrack->GetTPCpid(detPid) ; 
263   for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
264   fFlowTrack->SetRespFunTPC(detPid6) ; 
265  }
266  fFlowTrack->SetMaxPtsTPC(fNFound) ;                            
267  fFlowTrack->SetFitPtsTPC(nClus) ;                              
268  fFlowTrack->SetDedxTPC(fTrack->GetTPCsignal()) ;               
269  fFlowTrack->SetChi2TPC((Float_t)(fTrack->GetTPCchi2())) ;      
270  fFlowTrack->SetPatTPC(pAt) ;                                   
271  // -
272  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgITSx, fMagField, pVecAt) ; }
273  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgITSx, fMagField, pVecAt) ; }
274  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
275  nClus = fTrack->GetITSclusters(idX) ;
276  fNFound = 6 ; // ? fixed
277  if( (fTrack->GetStatus() & AliESDtrack::kITSpid) != 0 )
278  {
279   fTrack->GetITSpid(detPid) ; 
280   for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
281   fFlowTrack->SetRespFunITS(detPid6) ; 
282  } 
283  fFlowTrack->SetMaxPtsITS(fNFound) ;                            
284  fFlowTrack->SetFitPtsITS(nClus) ;                              
285  fFlowTrack->SetDedxITS(fTrack->GetITSsignal()) ;               
286  fFlowTrack->SetChi2ITS((Float_t)(fTrack->GetITSchi2())) ;      
287  fFlowTrack->SetPatITS(pAt) ;                                   
288  // -
289  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgTRDx, fMagField, pVecAt) ; } 
290  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgTRDx, fMagField, pVecAt) ; } 
291  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
292  nClus = fTrack->GetTRDclusters(idxr) ;
293  fNFound = fTrack->GetTRDncls() ;  // was 130
294  if( (fTrack->GetStatus() & AliESDtrack::kTRDpid) != 0 )
295  {
296   fTrack->GetTRDpid(detPid) ; 
297   for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
298   fFlowTrack->SetRespFunTRD(detPid6) ; 
299  }
300  fFlowTrack->SetMaxPtsTRD(fNFound) ;                            
301  fFlowTrack->SetNhitsTRD(nClus) ;                               
302  fFlowTrack->SetSigTRD(fTrack->GetTRDsignal()) ;                
303  fFlowTrack->SetChi2TRD((Float_t)fTrack->GetTRDchi2()) ;        
304  fFlowTrack->SetPatTRD(pAt) ;                                   
305  // -
306  if(fNewAli) { boh = fTrack->GetPxPyPzAt(AliFlowConstants::fgTOFx, fMagField, pVecAt) ; }
307  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(AliFlowConstants::fgTOFx, fMagField, pVecAt) ; }
308  pAt = (Float_t)Norm(pVecAt) ; if(!pAt) { pAt = (Float_t)Norm(gD) ; }
309  nClus = fTrack->GetTOFcluster() ;
310  fNFound = 0 ; if(fTrack->GetTOFCalChannel() > 0) { fNFound = 1 ; }
311  if( (fTrack->GetStatus() & AliESDtrack::kTOFpid) != 0 )
312  {
313   fTrack->GetTOFpid(detPid) ; 
314   for(Int_t bb=0;bb<5;bb++) { detPid6[bb] = detPid[bb] ; } detPid6[5] = 0. ; 
315   fFlowTrack->SetRespFunTOF(detPid6) ; 
316  }
317  fFlowTrack->SetMaxPtsTOF(fNFound) ;                            
318  fFlowTrack->SetNhitsTOF(nClus) ;                               
319  fFlowTrack->SetTofTOF(fTrack->GetTOFsignal()) ;                
320  fFlowTrack->SetChi2TOF(fTrack->GetTOFchi2()) ;                 
321  fFlowTrack->SetPatTOF(pAt) ;                                   
322  // -
323  Double_t rIn[3]  ; rIn[0] = 0.  ;  rIn[1] = 0. ;  rIn[2] = 0. ;
324  Double_t rOut[3] ; rOut[0] = 0. ; rOut[1] = 0. ; rOut[2] = 0. ;
325  // -
326  fTrack->GetInnerXYZ(rIn) ;                                     
327  fFlowTrack->SetZFirstPoint(rIn[2]) ; 
328  //fTrack->GetXYZAt(AliFlowConstants::fgTPCx,fMagField,rOut) ;          
329  fTrack->GetOuterXYZ(rOut) ;                                    
330  fFlowTrack->SetZLastPoint(rOut[2]) ; 
331  
332  // ESD-P.Id. = 5-vector of Best detectors probabilities for [e , mu , pi , K , p] 
333  Double_t trkPid[5] ; fTrack->GetESDpid(trkPid) ;               
334  Double_t trkPid6[AliFlowConstants::kPid] ; 
335  for(Int_t bb=0;bb<5;bb++) { trkPid6[bb] = trkPid[bb] ; } 
336  trkPid6[5] = 0. ;                                              // *!* implement P.Id. for Deuterim
337
338  // Bayesian P.Id. method (weighted probabilities for [e , mu , pi , K , p , d])
339  Double_t bsum = 0 ; 
340  Double_t bayePid[AliFlowConstants::kPid] ;    // normalized P.id
341  Double_t storedPid[AliFlowConstants::kPid] ;  // stored P.id
342  for(Int_t nB=0;nB<AliFlowConstants::kPid;nB++)  { bsum += trkPid6[nB]*AliFlowConstants::fgBayesian[nB] ; }
343  if(bsum)
344  {
345   for(Int_t nB=0;nB<AliFlowConstants::kPid;nB++) 
346   { 
347    bayePid[nB] = trkPid6[nB]*AliFlowConstants::fgBayesian[nB] / bsum ; 
348    storedPid[nB] = trkPid6[nB] ; 
349   }
350  }
351  else { cout << " ERROR - Empty Bayesian Vector !!! " << endl ; }
352  
353  fFlowTrack->SetElectronPositronProb(storedPid[0]);               
354  fFlowTrack->SetMuonPlusMinusProb(storedPid[1]);
355  fFlowTrack->SetPionPlusMinusProb(storedPid[2]);
356  fFlowTrack->SetKaonPlusMinusProb(storedPid[3]);
357  fFlowTrack->SetProtonPbarProb(storedPid[4]);
358  fFlowTrack->SetDeuteriumAntiDeuteriumProb(storedPid[5]);       // *!* implement P.Id. for Deuterim
359
360  // P.id. label given via the weighted prob.
361  const Int_t kCode[]   =  {11,13,211,321,2212,10010020} ;
362  Int_t kkk = 2 ;                        // if No id. -> then is a Pi
363  Float_t pidMax = bayePid[2] ;  // (if all equal, Pi probability get's the advantage to be the first)
364  for(Int_t iii=0; iii<5; iii++) 
365  {
366   if(bayePid[iii]>pidMax) { kkk = iii ; pidMax = bayePid[iii] ; }  // !!! Bayesian as well !!!
367  }
368  fBayesianAll[kkk]++ ; fSumAll++ ;      // goes on filling the vector of observed abundance 
369  //-
370  Int_t pdgCode = trkSign*kCode[kkk] ;
371  fFlowTrack->SetMostLikelihoodPID(pdgCode);
372  
373  return fFlowTrack ; 
374 }
375 //-----------------------------------------------------------------------
376 AliFlowV0* AliFlowMaker::FillFlowV0(AliESDv0* fV0)
377 {
378  // From the AliESDv0 (input) fills the AliFlowV0 (output) .
379
380  TString name = "" ; name += fV0Number ;
381  Int_t idx = fFlowEvent->V0Collection()->GetEntries() ;
382  fFlowV0 = (AliFlowV0*)(fFlowEvent->V0Collection()->New(idx)) ;
383  fFlowV0->SetName(name.Data()) ;
384  
385  // cout << " -v0- " << name.Data() << "(" << idx << ")"  << endl ;
386
387  // ESD particle label (link: KineTree-ESD)
388  Int_t label = -1 ; // TMath::Abs(fV0->GetLabel());
389  fFlowV0->SetLabel(label) ;                     
390
391  // reconstructed momentum of the V0
392  Double_t pxyz[3] ;
393  fV0->GetPxPyPz(pxyz[0],pxyz[1],pxyz[2]) ;
394
395  Float_t phi = (Float_t)Phi(pxyz) ; if(phi<0) { phi += 2*TMath::Pi() ; }
396  fFlowV0->SetPhi(phi) ;         
397  Float_t pt = (Float_t)Pt(pxyz) ; 
398  fFlowV0->SetPt(pt) ;           
399  Float_t eta = (Float_t)Eta(pxyz) ; 
400  fFlowV0->SetEta(eta) ;         
401
402  // reconstructed position of the V0 
403  Double_t xyz[3] ;              
404  fV0->GetXYZ(xyz[0],xyz[1],xyz[2]) ;            
405  fFlowV0->SetCrossPoint(xyz[0],xyz[1],xyz[2]) ;
406
407  // V0's impact parameter & error (chi2 , DCA , sigma , pointing angle)  
408  //fFlowV0->SetDca((Float_t)fV0->GetD()) ;    // GetDistNorm 
409  //fFlowV0->SetSigma((Float_t)fV0->GetDistSigma()) ;
410  //fFlowV0->SetCosPointingAngle((Float_t)fV0->GetV0CosineOfPointingAngle()) ;  
411  //fFlowV0->SetDaughtersDca(fV0->GetDcaV0Daughters()) ;
412  //fFlowV0->SetChi2((Float_t)fV0->GetChi2V0()) ;     // AliRoot v4-04-Release (December 2006)  
413  //fFlowV0->SetChi2((Float_t)fV0->GetChi2()) ;       // AliRoot v4-04-Release (old)
414  // ...when they'll stop changing the methods I'll enable the above lines. For now:
415  fFlowV0->SetDca(0.);
416  fFlowV0->SetSigma(0.1); 
417  fFlowV0->SetCosPointingAngle(1.) ;     
418  fFlowV0->SetDaughtersDca(0.) ;
419  fFlowV0->SetChi2(1.) ;
420
421  // P.id. 
422  Int_t pdgCode = fV0->GetPdgCode() ;
423  fFlowV0->SetMostLikelihoodPID(pdgCode);        
424
425  // mass 
426  fFlowV0->SetVmass((Float_t)fV0->GetEffMass()) ; 
427  
428  // daughters 
429  Int_t pN = fV0->GetPindex() ;
430  Int_t nN = fV0->GetNindex() ;
431  AliFlowTrack* pos = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(pN) ; 
432  AliFlowTrack* neg = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(nN) ; 
433  fFlowV0->SetDaughters(pos,neg) ;
434
435  return fFlowV0 ;
436 }
437 //-----------------------------------------------------------------------
438 Bool_t AliFlowMaker::CheckTrack(AliESDtrack* fTrack) const
439 {
440  // applies track cuts (pE , nHits , label)
441
442  Int_t idXt[180] ;  // used for Cluster Map ( see AliESDtrack::GetTPCclusters() )
443  Int_t   nHits = fTrack->GetTPCclusters(idXt) ;
444  Float_t pE    = fTrack->GetP() ;
445  Int_t   label = fTrack->GetLabel() ;
446  
447  if(fNHits && (nHits<=fNHits))                                           { return kFALSE ; }
448  if((fElow < fEup) && ((pE<fElow) || (pE>fEup)))                         { return kFALSE ; }
449  if((fLabel[0] < fLabel[1]) && ((label<fLabel[0]) || (label>fLabel[1]))) { return kFALSE ; }
450  //if(fPrimary && ...)                                                   { return kFALSE ; }
451
452  return kTRUE ;
453 }
454 //-----------------------------------------------------------------------
455 Bool_t AliFlowMaker::CheckV0(AliESDv0* fV0) const
456 {
457  // applies v0 cuts (dummy)
458  
459  if(!fV0) { return kFALSE ; }
460
461  return kTRUE ;
462 }
463 //-----------------------------------------------------------------------
464 Bool_t AliFlowMaker::CheckEvent(AliESD* fESD) const
465 {
466  // applies event cuts (dummy)
467  
468  if(!fESD) { return kFALSE ; }
469
470  return kTRUE ;
471 }
472 //-----------------------------------------------------------------------
473 void AliFlowMaker::PrintCutList()
474
475  // Prints the list of Cuts 
476
477  cout << " * ESD cuts list * " << endl ; 
478  if(fLabel[0]<fLabel[1])
479  { 
480   cout << "  Track Label [ " << fLabel[0] << " , " << fLabel[1] << " ] " << endl ; 
481  }
482  if(fNHits)
483  { 
484   cout << "  TPC clusters > " << fNHits << endl ; 
485  }
486  if(fElow<fEup)
487  { 
488   cout << "  Track Energy (P_total) [ " << fElow << " , " << fEup << " ] " << endl ; 
489  }
490  cout << " *               * " << endl ;
491 }
492 //-----------------------------------------------------------------------
493
494 //-----------------------------------------------------------------------
495 //*** USEFULL METHODS for a 3-array of double (~ TVector3) ***
496 //-----------------------------------------------------------------------
497 Double_t AliFlowMaker::Norm(Double_t nu[3])
498
499  // returns the norm of a double[3] 
500
501  Double_t norm2 = nu[0]*nu[0] + nu[1]*nu[1] + nu[2]*nu[2] ;
502  return TMath::Sqrt(norm2) ; 
503 }
504 //-----------------------------------------------------------------------
505 Double_t AliFlowMaker::Phi(Double_t nu[3])
506 {
507  // returns the azimuthal angle of a double[3] 
508
509  if(nu[0]==0 && nu[1]==0) { return 0. ; }
510  else                     { return TMath::ATan2(nu[1],nu[0]) ; }
511 }
512 //-----------------------------------------------------------------------
513 Double_t AliFlowMaker::Pt(Double_t nu[3])
514 {
515  // returns the transvers momentum of a double[3] 
516
517  Double_t trans = nu[0]*nu[0] + nu[1]*nu[1] ;
518  return TMath::Sqrt(trans) ; 
519 }
520 //-----------------------------------------------------------------------
521 Double_t AliFlowMaker::Eta(Double_t nu[3])
522 {
523  // returns the PseudoRapidity of a double[3] 
524  // if transvers momentum = 0 --> returns +/- 1.000
525
526  Double_t m = Norm(nu) ;
527  if(nu[0]!=0 || nu[1]!=0) { return 0.5*TMath::Log((m+nu[2])/(m-nu[2])) ; }
528  else                     { return TMath::Sign((Double_t)1000.,nu[2]) ; }
529 }
530 //-----------------------------------------------------------------------
531
532 #endif