Added ESD reader for HLT tracks.
[u/mrichter/AliRoot.git] / JETAN / AliJetParticlesReaderESD.cxx
CommitLineData
d7c6ab14 1// $Id$
2
3//____________________________________________________________________
4//////////////////////////////////////////////////////////////////////
5// //
6// class AliHBTReaderESD //
7// //
8// Reader for ALICE Event Summary Data (ESD). //
9// //
10// Piotr.Skowronski@cern.ch //
11// //
12//////////////////////////////////////////////////////////////////////
13
14#include <TMath.h>
15#include <TPDGCode.h>
16#include <TString.h>
17#include <TObjString.h>
18#include <TTree.h>
19#include <TFile.h>
20#include <TKey.h>
21#include <TH1.h>
22#include <AliESD.h>
23#include <AliESDtrack.h>
24#include <AliJetEventParticles.h>
25#include "AliJetParticlesReaderESD.h"
26
27ClassImp(AliJetParticlesReaderESD)
28
29AliJetParticlesReaderESD::AliJetParticlesReaderESD(const Char_t* esdfilename) :
30 AliJetParticlesReader(),
31 fESDFileName(esdfilename),
5d60c8f9 32 fESD(0),
d7c6ab14 33 fFile(0),
5d60c8f9 34 fTree(0),
d7c6ab14 35 fKeyIterator(0),
36 fPassFlag(AliESDtrack::kTPCrefit)
37{
38 //constructor
39}
40
41/********************************************************************/
42
43AliJetParticlesReaderESD::AliJetParticlesReaderESD(
44 TObjArray* dirs,
45 const Char_t* esdfilename) :
46 AliJetParticlesReader(dirs),
47 fESDFileName(esdfilename),
5d60c8f9 48 fESD(0),
d7c6ab14 49 fFile(0),
5d60c8f9 50 fTree(0),
d7c6ab14 51 fKeyIterator(0),
52 fPassFlag(AliESDtrack::kTPCrefit)
53{
54 //constructor
55}
56
57/********************************************************************/
58
59AliJetParticlesReaderESD::~AliJetParticlesReaderESD()
60{
61 //desctructor
62 if(fFile) delete fFile;
5d60c8f9 63 if(fESD) delete fESD;
64 if(fTree) delete fTree;
d7c6ab14 65 if(fKeyIterator) delete fKeyIterator;
66}
67
68
69Int_t AliJetParticlesReaderESD::ReadNext()
70{
71 //reads next event from fFile
72
73 do // is OK even if 0 dirs specified,
74 { // in that case we try to read from "./"
75 if (fFile == 0)
76 {
77 fFile = OpenFile(fCurrentDir);
78 if (fFile == 0)
79 {
80 Error("ReadNext","Cannot get fFile for dir no. %d",fCurrentDir);
81 fCurrentDir++;
82 continue;
83 }
84
d7c6ab14 85 fCurrentEvent = 0;
d7c6ab14 86 //fFile->GetListOfKeys()->Print();
5d60c8f9 87
88 if(fTree) delete fTree;
89 fTree = dynamic_cast<TTree*>(fFile->Get("esdTree"));
90 if(fTree)
91 fTree->SetBranchAddress("ESD",&fESD);
92 else
93 fKeyIterator = new TIter(fFile->GetListOfKeys());
d7c6ab14 94 }
95
5d60c8f9 96 if(fTree)
d7c6ab14 97 {
5d60c8f9 98 if(fCurrentEvent>=fTree->GetEntries())
99 {
100 fCurrentDir++;
101 delete fTree;
102 fTree = 0;
103 delete fFile;
104 fFile = 0;
105 continue;
106 }
107 fTree->GetEvent(fCurrentEvent);
108 }
109 else
110 { // "old" way via ESD objects stored in root file
111 TKey* key = (TKey*)fKeyIterator->Next();
112 if (key == 0)
113 {
114 fCurrentDir++;
115 delete fKeyIterator;
116 fKeyIterator = 0;
117 delete fFile; //we have to assume there are no more ESD objects in the fFile
118 fFile = 0;
119 continue;
120 }
121 TString esdname = "ESD";
122 esdname+=fCurrentEvent;
123 if(fESD) delete fESD;
124 fESD = dynamic_cast<AliESD*>(fFile->Get(esdname));
125 if (fESD == 0)
126 {
127 Info("ReadNext","Can not find AliESD object named %s",esdname.Data());
128 fCurrentDir++;
129 delete fKeyIterator;
130 fKeyIterator = 0;
131 delete fFile;//we have to assume there is no more ESD objects in the fFile
132 fFile = 0;
133 continue;
134 }
d7c6ab14 135 }
5d60c8f9 136 ReadESD(fESD);
d7c6ab14 137 fCurrentEvent++;
138 fNEventsRead++;
d7c6ab14 139 return kTRUE;//success -> read one event
140 } while(fCurrentDir < GetNumberOfDirs());
5d60c8f9 141 //end of loop over directories specified in fDirs Obj Array
d7c6ab14 142 return kFALSE; //no more directories to read
143}
144
145/**********************************************************/
146
147Int_t AliJetParticlesReaderESD::ReadESD(AliESD* esd)
148{
149 //Reads one ESD
150
151 if (esd == 0)
152 {
153 Error("ReadESD","ESD is NULL");
154 return kFALSE;
155 }
156
157 //TDatabasePDG* pdgdb = TDatabasePDG::Instance();
158 //if (pdgdb == 0)
159 //{
160 // Error("ReadESD","Can not get PDG Database Instance.");
161 // return kFALSE;
162 // }
163
164 //Float_t mf = esd->GetMagneticField();
165 //if ( (mf == 0.0) && (fNTrackPoints > 0) )
166 //{
167 // Error("ReadESD","Magnetic Field is 0 and Track Points demanded. Skipping to next event.");
168 // return kFALSE;
169 //}
170
171 Info("ReadESD","Reading Event %d",fCurrentEvent);
172 if((!fOwner) || (fEventParticles==0))
173 fEventParticles = new AliJetEventParticles();
174
175 Double_t vertexpos[3];//vertex position
176 const AliESDVertex* kvertex = esd->GetVertex();
177 if (kvertex == 0)
178 {
179 Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
180 vertexpos[0] = 0.0;
181 vertexpos[1] = 0.0;
182 vertexpos[2] = 0.0;
183 }
184 else
185 {
186 kvertex->GetXYZ(vertexpos);
187 }
188 fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
189
190 const Int_t kntr = esd->GetNumberOfTracks();
191 Info("ReadESD","Found %d tracks.",kntr);
192 for (Int_t i = 0;i<kntr; i++)
193 {
194 const AliESDtrack *kesdtrack = esd->GetTrack(i);
195 if (kesdtrack == 0)
196 {
197 Error("ReadESD","Can not get track %d", i);
198 continue;
199 }
200
201 if ((kesdtrack->GetStatus() & fPassFlag) == kFALSE)
202 {
5d60c8f9 203 Info("ReadNext","Particle skipped: %ud.",kesdtrack->GetStatus());
d7c6ab14 204 continue;
205 }
206
207 Double_t mom[3]; //momentum
208 kesdtrack->GetPxPyPz(mom);
209 //kesdtrack->GetConstrainedPxPyPz(mom);
210 //Double_t pos[3];//position
211 //kesdtrack->GetXYZ(pos);
212 //kesdtrack->GetConstrainedXYZ(pos);
213 const Float_t kmass=kesdtrack->GetMass();
214 const Float_t kp2=mom[0]*mom[0]+mom[1]*mom[1]+mom[2]*mom[2];
215 const Float_t ketot=TMath::Sqrt(kmass*kmass+kp2);
216 const Float_t kpt=TMath::Sqrt(mom[0]*mom[0]+mom[1]*mom[1]);
217 const Float_t kp=TMath::Sqrt(kp2);
218 const Float_t keta=0.5*TMath::Log((kp+mom[2]+1e-30)/(kp-mom[2]+1e-30));
219 const Float_t kphi=TMath::Pi()+TMath::ATan2(-mom[1],-mom[0]);
220 if(IsAcceptedParticle(kpt,kphi,keta))
221 fEventParticles->AddParticle(mom[0],mom[1],mom[2],ketot,i,kesdtrack->GetLabel(),kpt,kphi,keta);
222
223 } // loop over tracks
224
225 return kTRUE;
226}
227
228/**********************************************************/
229
230void AliJetParticlesReaderESD::Rewind()
231{
232 //rewinds reading
233
234 if(fFile) delete fFile;
235 if(fKeyIterator) delete fKeyIterator;
236 fKeyIterator = 0;
237 fFile = 0;
238 fCurrentDir = 0;
239 fNEventsRead = 0;
240}
241
242/**********************************************************/
243
244TFile* AliJetParticlesReaderESD::OpenFile(Int_t n)
245{
246 //opens fFile with kine tree
247
248 const TString& dirname = GetDirName(n);
249 if (dirname == "")
250 {
251 Error("OpenFiles","Can not get directory name");
252 return 0;
253 }
254 TString filename = dirname +"/"+ fESDFileName;
255 TFile *ret = TFile::Open(filename.Data());
256 if (ret == 0)
257 {
258 Error("OpenFiles","Can't open fFile %s",filename.Data());
259 return 0;
260 }
261 if (!ret->IsOpen())
262 {
263 Error("OpenFiles","Can't open fFile %s",filename.Data());
264 return 0;
265 }
266
267 return ret;
268}