00d5f8debcd66162218480013c0936b73eaf1401
[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 <TTree.h>
25 #include <TChain.h>
26 #include <TString.h>
27 #include <TObject.h>
28 #include <TObjArray.h>
29 #include <TParticle.h>
30 #include <TParticlePDG.h>
31 #include <TDatabasePDG.h>
32
33 // AliRoot things
34 #include "AliESD.h"
35 #include "AliESDVertex.h"
36 #include "AliESDtrack.h"
37 #include "AliESDv0.h"
38 #include "AliKalmanTrack.h"
39
40 // Flow things
41 #include "AliFlowEvent.h"
42 #include "AliFlowTrack.h"
43 #include "AliFlowV0.h"
44 #include "AliFlowConstants.h"
45 #include "AliFlowMaker.h"
46
47 // ANSI things
48 #include <math.h>
49 #include <stdio.h>
50 #include <stdlib.h>
51 #include <iostream>
52 using namespace std; //required for resolving the 'cout' symbol
53
54 ClassImp(AliFlowMaker) 
55 //-----------------------------------------------------------------------
56 AliFlowMaker::AliFlowMaker()
57 {
58  // default constructor 
59  // resets counters , sets defaults
60  
61  fNewAli = kFALSE ;
62  // flags
63  fLoopTrks = kTRUE ;
64  fLoopV0s = kTRUE ;     
65  fCounter = 0 ; 
66  // loop variable
67  fRunID = 0 ;       
68  fEventNumber = 0 ; 
69  fTrackNumber = 0 ; 
70  fV0Number = 0 ;    
71  fEventNumber = 0 ;
72  fNumberOfTracks = 0 ;
73  fNumberOfV0s = 0 ; 
74  fMagField = 0. ;
75  // counters
76  fGoodTracks = 0  ;  
77  fGoodV0s    = 0  ;     
78  fGoodTracksEta = 0  ;
79  fPosiTracks = 0  ;  
80  fNegaTracks = 0  ;  
81  fUnconstrained = 0  ;
82  fCutEvts    = 0  ;     
83  fCutTrks    = 0  ;     
84  fCutV0s     = 0  ;
85  for(Int_t bb=0;bb<5;bb++) { fBayesianAll[bb] = 0 ; } ;
86  fSumAll = 0  ; 
87
88  // trak cuts
89  fNHits = 0 ;                   
90  fElow = 0.001 ; fEup = 1000. ;   
91  fLabel[0] = 0 ; fLabel[1] = -1 ;
92 }
93 //-----------------------------------------------------------------------
94 AliFlowMaker::~AliFlowMaker()
95 {
96  // default destructor (no actions) 
97 }
98 //-----------------------------------------------------------------------
99
100 //-----------------------------------------------------------------------
101 AliFlowEvent* AliFlowMaker::FillFlowEvent(AliESD* fESD)
102 {
103  // From the AliESD (input) fills the AliFlowEvent (output) . 
104  // It loops on track & v0 and calls the methods to fill the arrays . 
105  // ... . 
106
107  fFlowEvent = new AliFlowEvent() ; if(!fFlowEvent) { return 0 ; }
108  //cout << " -evt- " << fFlowEvent << endl ;
109
110  fRunID = fESD->GetRunNumber() ;
111  fEventNumber = fESD->GetEventNumber() ;
112  fNumberOfTracks = fESD->GetNumberOfTracks() ;
113  fNumberOfV0s = fESD->GetNumberOfV0s() ;
114  //
115  cout << " *evt n. " << fEventNumber << " (run " << fRunID << ")  -  tracks: " << fNumberOfTracks << " ,   v0s " << fNumberOfV0s << endl ;
116
117  // Event id 
118  fFlowEvent->SetRunID(fRunID) ;                
119  fFlowEvent->SetEventID(fEventNumber) ;         
120  fFlowEvent->SetOrigMult((UInt_t)fNumberOfTracks) ;
121
122  // Run information (fixed - ???)
123  fMagField = fESD->GetMagneticField() ; // cout << " *fMagField " << fMagField << endl ;
124  fFlowEvent->SetMagneticField(fMagField) ;      
125  fFlowEvent->SetCenterOfMassEnergy(Flow::fCenterOfMassEnergy) ; 
126  fFlowEvent->SetBeamMassNumberEast(Flow::fBeamMassNumberEast) ; 
127  fFlowEvent->SetBeamMassNumberWest(Flow::fBeamMassNumberWest) ; 
128
129  // Trigger information (now is: ULon64_t - some trigger mask)
130  fFlowEvent->SetL0TriggerWord((Int_t)fESD->GetTriggerMask()); 
131  
132  // Get primary vertex position
133  fVertex = (AliESDVertex*)fESD->GetVertex() ;
134  Double_t position[3] ; 
135  fVertex->GetXYZ(position) ;
136  fFlowEvent->SetVertexPos((Float_t)position[0],(Float_t)position[1],(Float_t)position[2]) ; 
137
138  // Zero Degree Calorimeter information
139  Int_t zdcp = fESD->GetZDCParticipants() ; 
140  Float_t zdce[3] ; 
141  zdce[0] = fESD->GetZDCN1Energy() + fESD->GetZDCN2Energy(); 
142  zdce[1] = fESD->GetZDCP1Energy() + fESD->GetZDCP2Energy() ; 
143  zdce[2] = fESD->GetZDCEMEnergy() ;
144  fFlowEvent->SetZDCpart(zdcp);                          
145  fFlowEvent->SetZDCenergy(zdce[0],zdce[1],zdce[2]);     
146
147  // Track loop
148  if(fLoopTrks)
149  {
150   Int_t badTrks = 0 ;
151   for(fTrackNumber=0;fTrackNumber<fNumberOfTracks;fTrackNumber++) 
152   {
153    fTrack = fESD->GetTrack(fTrackNumber) ;
154    if(CheckTrack(fTrack))
155    {
156     fFlowTrack = FillFlowTrack(fTrack) ;   
157     fFlowEvent->TrackCollection()->Add(fFlowTrack) ;   
158     fGoodTracks++ ;                                   
159    }
160    else { badTrks++ ; continue ; }
161   }
162   fCutTrks += badTrks ;
163  }
164
165  // V0 loop
166  if(fLoopV0s)
167  {
168   Int_t badV0s = 0 ;
169   for(fV0Number=0;fV0Number<fNumberOfV0s;fV0Number++) 
170   {
171    fV0 = fESD->GetV0(fV0Number) ;                       
172    if(CheckV0(fV0))
173    {
174     fFlowV0 =  FillFlowV0(fV0) ;
175     fFlowEvent->V0Collection()->Add(fFlowV0) ;
176     fGoodV0s++ ;                        
177    }
178    else { badV0s++ ; continue ; }
179   }
180   fCutV0s += badV0s ;
181  }
182
183  // Evt setting stuff
184  fFlowEvent->SetCentrality();   
185                                 
186  fCounter++ ; 
187  return fFlowEvent ;
188 }
189 //----------------------------------------------------------------------
190 AliFlowTrack* AliFlowMaker::FillFlowTrack(AliESDtrack* fTrack)
191 {
192  // From the AliESDtrack (input) fills the AliFlowTrack (output) .
193
194  TString name = "" ; name += fTrackNumber ;
195  fFlowTrack = new AliFlowTrack(name.Data()) ;
196  //cout << " -tr- " << name.Data() << endl ;
197  
198  // ESD particle label (link: KineTree-ESD)
199  Int_t label = TMath::Abs(fTrack->GetLabel());
200  fFlowTrack->SetLabel(label) ;                  
201
202  // signed DCA from ESDtrack
203  Float_t xy = 0 ; Float_t z = 0 ; 
204  fTrack->GetImpactParameters(xy,z) ; 
205  fFlowTrack->SetDcaSigned(xy,z) ;                   
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) ; 
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) ; 
234   fFlowTrack->SetPt(pt) ;                               
235   Float_t eta = (Float_t)Eta(cD) ; 
236   fFlowTrack->SetEta(eta) ;                             
237  
238   // number of constrainable tracks with |eta| < Flow::fEtaGood (0.9)
239   if(TMath::Abs(eta) < Flow::fEtaGood)  { 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 trk_sign = (Int_t)fTrack->GetSign() ; 
251  fFlowTrack->SetCharge(trk_sign) ;              
252  if(trk_sign>0)         { fPosiTracks++ ; }
253  else if(trk_sign<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(Flow::fTPCx, fMagField, pVecAt) ; }
270  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(Flow::fTPCx, 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(Flow::fITSx, fMagField, pVecAt) ; }
281  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(Flow::fITSx, 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(Flow::fTRDx, fMagField, pVecAt) ; } 
292  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(Flow::fTRDx, 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(Flow::fTOFx, fMagField, pVecAt) ; }
303  else        { boh = fTrack->GetInnerParam()->GetPxPyPzAt(Flow::fTOFx, 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(Flow::fTPCx,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[Flow::nPid] ; 
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[Flow::nPid] ;    // normalized P.id
331  Double_t storedPid[Flow::nPid] ;  // stored P.id
332  for(Int_t nB=0;nB<Flow::nPid;nB++)  { Bsum += trkPid6[nB]*Flow::fBayesian[nB] ; }
333  if(Bsum)
334  {
335   for(Int_t nB=0;nB<Flow::nPid;nB++) 
336   { 
337    bayePid[nB] = trkPid6[nB]*Flow::fBayesian[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 code[]   =  {11,13,211,321,2212,10010020} ;
352  Int_t kkk = 2 ;                        // if No id. -> then is a Pi
353  Float_t pid_max = 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]>pid_max) { kkk = iii ; pid_max = bayePid[iii] ; }  // !!! Bayesian as well !!!
357  }
358  fBayesianAll[kkk]++ ; fSumAll++ ;      // goes on filling the vector of observed abundance 
359  //-
360  Int_t pdg_code = trk_sign*code[kkk] ;
361  fFlowTrack->SetMostLikelihoodPID(pdg_code);
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  fFlowV0 = new AliFlowV0(name.Data()) ;
372  //cout << " -v0- " << name.Data() << endl ;
373
374  Double_t Pxyz[3] ;             // reconstructed momentum of the V0
375  fV0->GetPxPyPz(Pxyz[0],Pxyz[1],Pxyz[2]) ;
376
377  Float_t phi = (Float_t)Phi(Pxyz) ; if(phi<0) { phi += 2*TMath::Pi() ; }
378  fFlowV0->SetPhi(phi) ;         
379  Float_t pt = (Float_t)Pt(Pxyz) ; 
380  fFlowV0->SetPt(pt) ;           
381  Float_t eta = (Float_t)Eta(Pxyz) ; 
382  fFlowV0->SetEta(eta) ;         
383
384  Double_t xyz[3] ;              // reconstructed position of the V0 
385  fV0->GetXYZ(xyz[0],xyz[1],xyz[2]) ;            
386  fFlowV0->SetCrossPoint(xyz[0],xyz[1],xyz[2]) ;
387
388  // // V0's impact parameter & error (chi2 , DCA , sigma)  
389  // fFlowV0->SetDca((Float_t)fV0->GetDistNorm()) ;    // GetD()  
390  // fFlowV0->SetSigma((Float_t)fV0->GetDistSigma()) ;
391  // fFlowV0->SetChi2((Float_t)fV0->GetChi2V0()) ;     // AliRoot v4-04-Release (December 2006)  
392  // fFlowV0->SetChi2((Float_t)fV0->GetChi2()) ;       // AliRoot v4-04-Release (old)
393  // // ...when they'll stop changing the methods I may enable the above lines. For now:
394  fFlowV0->SetChi2(1.) ;
395  fFlowV0->SetDca(1.);      // GetD()
396  fFlowV0->SetSigma(1.); 
397
398  // P.id. 
399  Int_t pdg_code = fV0->GetPdgCode() ;
400  fFlowV0->SetMostLikelihoodPID(pdg_code);       
401
402  // mass 
403  fFlowV0->SetVmass((Float_t)fV0->GetEffMass()) ; 
404  
405  Int_t pN = fV0->GetPindex() ;
406  Int_t nN = fV0->GetNindex() ;
407  AliFlowTrack* pos = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(pN) ; 
408  AliFlowTrack* neg = (AliFlowTrack*)fFlowEvent->TrackCollection()->At(nN) ; 
409  fFlowV0->SetDaughters(pos,neg) ;
410
411  return fFlowV0 ;
412 }
413 //-----------------------------------------------------------------------
414 Bool_t AliFlowMaker::CheckTrack(AliESDtrack* fTrack)
415 {
416  // applies track cuts (E , nHits , label)
417
418  Int_t idXt[180] ;  // used for Cluster Map ( see AliESDtrack::GetTPCclusters() )
419  Int_t   nHits = fTrack->GetTPCclusters(idXt) ;
420  Float_t E     = fTrack->GetP() ;
421  Int_t   label = fTrack->GetLabel() ;
422  
423  if(fNHits && (nHits<=fNHits))                                           { return kFALSE ; }
424  if((fElow < fEup) && ((E<fElow) || (E>fEup)))                           { return kFALSE ; }
425  if((fLabel[0] < fLabel[1]) && ((label<fLabel[0]) || (label>fLabel[1]))) { return kFALSE ; }
426  //if(fPrimary && ...)                                                   { return kFALSE ; }
427
428  return kTRUE ;
429 }
430 //-----------------------------------------------------------------------
431 Bool_t AliFlowMaker::CheckV0(AliESDv0* fV0)
432 {
433  // applies v0 cuts (dummy)
434  
435  if(!fV0) { return kFALSE ; }
436
437  return kTRUE ;
438 }
439 //-----------------------------------------------------------------------
440 Bool_t AliFlowMaker::CheckEvent(AliESD* fESD)
441 {
442  // applies event cuts (dummy)
443  
444  if(!fESD) { return kFALSE ; }
445
446  return kTRUE ;
447 }
448 //-----------------------------------------------------------------------
449 void AliFlowMaker::PrintCutList()
450
451  // Prints the list of Cuts 
452
453  cout << " * ESD cuts list * " << endl ; 
454  if(fLabel[0]<fLabel[1])
455  { 
456   cout << "  Track Label [ " << fLabel[0] << " , " << fLabel[1] << " ] " << endl ; 
457  }
458  if(fNHits)
459  { 
460   cout << "  TPC clusters > " << fNHits << endl ; 
461  }
462  if(fElow<fEup)
463  { 
464   cout << "  Track Energy (P_total) [ " << fElow << " , " << fEup << " ] " << endl ; 
465  }
466  cout << " *               * " << endl ;
467 }
468 //-----------------------------------------------------------------------
469
470 //-----------------------------------------------------------------------
471 //*** USEFULL METHODS for a 3-array of double (~ TVector3) ***
472 //-----------------------------------------------------------------------
473 Double_t AliFlowMaker::Norm(Double_t nu[3])
474
475  // returns the norm of a double[3] 
476
477  Double_t norm2 = nu[0]*nu[0] + nu[1]*nu[1] + nu[2]*nu[2] ;
478  return TMath::Sqrt(norm2) ; 
479 }
480 //-----------------------------------------------------------------------
481 Double_t AliFlowMaker::Phi(Double_t nu[3])
482 {
483  // returns the azimuthal angle of a double[3] 
484
485  if(nu[0]==0 && nu[1]==0) { return 0. ; }
486  else                     { return TMath::ATan2(nu[1],nu[0]) ; }
487 }
488 //-----------------------------------------------------------------------
489 Double_t AliFlowMaker::Pt(Double_t nu[3])
490 {
491  // returns the transvers momentum of a double[3] 
492
493  Double_t trans = nu[0]*nu[0] + nu[1]*nu[1] ;
494  return TMath::Sqrt(trans) ; 
495 }
496 //-----------------------------------------------------------------------
497 Double_t AliFlowMaker::Eta(Double_t nu[3])
498 {
499  // returns the PseudoRapidity of a double[3] 
500  // if transvers momentum = 0 --> returns +/- 1.000
501
502  Double_t m = Norm(nu) ;
503  if(nu[0]!=0 || nu[1]!=0) { return 0.5*TMath::Log((m+nu[2])/(m-nu[2])) ; }
504  else                     { return TMath::Sign((Double_t)1000.,nu[2]) ; }
505 }
506 //-----------------------------------------------------------------------
507
508 #endif