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