Workaround to compute conformal track parameters at vertex.
[u/mrichter/AliRoot.git] / JETAN / AliJetParticlesReaderHLT.cxx
1 // $Id$
2
3 //___________________________________________________________________________
4 /////////////////////////////////////////////////////////////////////////////
5 //                                                                         //
6 // File reader for HLT tracks ESD                                          //
7 //                                                                         //
8 // loizides@ikf.uni-frankfurt.de                                           //
9 /////////////////////////////////////////////////////////////////////////////
10
11 #include <Riostream.h>
12 #include <TMath.h>
13 #include <TString.h>
14 #include <TObjString.h>
15 #include <TTree.h>
16 #include <TFile.h>
17 #include <AliESD.h>
18 #include <AliESDtrack.h>
19 #include <AliESDHLTtrack.h>
20 #include <AliL3Track.h>
21 #include <AliKalmanTrack.h>
22 #include <AliJetEventParticles.h>
23 #include "AliJetParticlesReaderHLT.h"
24
25 ClassImp(AliJetParticlesReaderHLT)
26
27 AliJetParticlesReaderHLT::AliJetParticlesReaderHLT(Bool_t bMapper, const Char_t* esdfilename) :
28   AliJetParticlesReaderESD(0,esdfilename),
29   fTrackerType(bMapper),
30   fMinHits(0),
31   fMinWeight(0)
32 {
33   //constructor
34 }
35
36 /********************************************************************/
37   
38 AliJetParticlesReaderHLT::AliJetParticlesReaderHLT(
39                                       Bool_t bMapper,
40                                       TObjArray* dirs,
41                                       const Char_t* esdfilename) :
42   AliJetParticlesReaderESD(0,dirs,esdfilename),
43   fTrackerType(bMapper),
44   fMinHits(0),
45   fMinWeight(0)
46
47 {
48   //constructor
49 }
50
51 /********************************************************************/
52
53 AliJetParticlesReaderHLT::~AliJetParticlesReaderHLT()
54 {
55   //desctructor
56 }
57
58 Int_t AliJetParticlesReaderHLT::ReadESD(AliESD* esd)
59 {
60   //Reads one ESD
61
62   if (esd == 0)
63    {
64      Error("ReadESD","ESD is NULL");
65      return kFALSE;
66    }
67
68   Float_t mf = esd->GetMagneticField(); 
69   if (mf <= 0.0)  
70   {
71      Error("ReadESD","Magnetic Field is 0. Skipping to next event.");
72      return kFALSE;
73   }
74   AliKalmanTrack::SetMagneticField(mf/10.);
75
76   Info("ReadESD","Reading Event %d",fCurrentDir*1000+fCurrentEvent);
77   if((!fOwner) || (fEventParticles==0)) 
78     fEventParticles = new AliJetEventParticles();
79
80   Int_t ntr=0;
81   if(fTrackerType){
82     ntr =esd->GetNumberOfHLTHoughTracks();
83     Info("ReadESD","Found %d conformal tracks.",ntr);
84   } else {
85     ntr=esd->GetNumberOfHLTHoughTracks();
86     Info("ReadESD","Found %d hough tracks.",ntr);
87   }
88   fEventParticles->Reset(ntr);
89
90   TString headdesc="";
91   headdesc+="Run ";
92   headdesc+=esd->GetRunNumber();
93   headdesc+=": Ev ";
94   headdesc+=esd->GetEventNumber();
95   fEventParticles->SetHeader(headdesc);
96
97   Double_t vertexpos[3];//vertex position
98   const AliESDVertex* kvertex = esd->GetVertex();
99   if (kvertex == 0)
100    {
101      Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
102      vertexpos[0] = 0.0;
103      vertexpos[1] = 0.0;
104      vertexpos[2] = 0.0;
105    }
106   else
107    {
108      kvertex->GetXYZ(vertexpos);
109    }
110
111   fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
112
113   AliL3Track l3;
114   for (Int_t i = 0;i<ntr; i++) {
115     AliESDHLTtrack *kesdtrack;
116     if(fTrackerType){
117       kesdtrack=esd->GetHLTConfMapTrack(i);
118     } else {
119       kesdtrack=esd->GetHLTHoughTrack(i);
120     }
121
122     if (kesdtrack == 0)
123       {
124         Error("ReadESD","Can not get track %d", i);
125         continue;
126       }
127
128     if(fTrackerType){
129       //if(!kesdtrack->ComesFromMainVertex()) continue;
130       l3.SetFirstPoint(kesdtrack->GetFirstPointX(),kesdtrack->GetFirstPointY(),kesdtrack->GetFirstPointZ());
131       l3.SetLastPoint(kesdtrack->GetLastPointX(),kesdtrack->GetLastPointY(),kesdtrack->GetLastPointZ());
132       l3.SetCharge(kesdtrack->GetCharge());
133       l3.SetPt(kesdtrack->GetPt());
134       l3.SetTgl(kesdtrack->GetTgl());
135       l3.SetPsi(kesdtrack->GetPsi());
136       l3.CalculateHelix();
137       l3.SetFirstPoint(0,0,0);
138       l3.UpdateToFirstPoint();
139     }
140
141     //const Float_t kpid=kesdtrack->GetPID();
142     const Int_t knhits=kesdtrack->GetNHits();
143     const Int_t kweight=kesdtrack->GetWeight();
144     //cout << i << " " << kweight << " " << knhits << endl;
145     if((fMinHits>0) && (knhits<fMinHits)) continue;    
146     if(kweight>1000) continue; //avoid ghosts 
147     if((fMinWeight>0) && (kweight<fMinWeight)) continue;
148
149     const Float_t kpx=kesdtrack->GetPx();
150     const Float_t kpy=kesdtrack->GetPy();
151     const Float_t kpz=kesdtrack->GetPz();
152     const Float_t kpt=kesdtrack->GetPt();
153     const Float_t kp=TMath::Sqrt(kpz*kpz+kpt*kpt);
154     const Float_t keta=0.5*TMath::Log((kp+kpz+1e-30)/(kp-kpz+1e-30)); 
155     const Float_t kphi=TMath::Pi()+TMath::ATan2(-kpy,-kpx);
156
157     if(IsAcceptedParticle(kpt,kphi,keta))
158       fEventParticles->AddParticle(kpx,kpy,kpz,kp,i,kesdtrack->GetMCid(),knhits,kpt,kphi,keta);
159   } //loop over tracks
160
161   return kTRUE;
162 }