3 //____________________________________________________________________
4 //////////////////////////////////////////////////////////////////////
6 // class AliHBTReaderESD //
8 // Reader for ALICE Event Summary Data (ESD). //
10 // Piotr.Skowronski@cern.ch //
12 //////////////////////////////////////////////////////////////////////
17 #include <TObjString.h>
23 #include <AliESDtrack.h>
24 #include <AliKalmanTrack.h>
25 #include <AliJetEventParticles.h>
26 #include "AliJetParticlesReaderESD.h"
28 ClassImp(AliJetParticlesReaderESD)
30 AliJetParticlesReaderESD::AliJetParticlesReaderESD(const Char_t* esdfilename) :
31 AliJetParticlesReader(),
32 fESDFileName(esdfilename),
37 fPassFlag(AliESDtrack::kTPCrefit)
42 /********************************************************************/
44 AliJetParticlesReaderESD::AliJetParticlesReaderESD(
46 const Char_t* esdfilename) :
47 AliJetParticlesReader(dirs),
48 fESDFileName(esdfilename),
53 fPassFlag(AliESDtrack::kTPCrefit)
58 /********************************************************************/
60 AliJetParticlesReaderESD::~AliJetParticlesReaderESD()
64 if(fTree) delete fTree;
65 if(fKeyIterator) delete fKeyIterator;
66 if(fFile) delete fFile;
70 Int_t AliJetParticlesReaderESD::ReadNext()
72 //reads next event from fFile
74 do // is OK even if 0 dirs specified,
75 { // in that case we try to read from "./"
78 fFile = OpenFile(fCurrentDir);
81 Error("ReadNext","Can't get fFile for dir no. %d",fCurrentDir);
87 //fFile->GetListOfKeys()->Print();
89 if(fTree) delete fTree;
90 fTree = dynamic_cast<TTree*>(fFile->Get("esdTree"));
92 fTree->SetBranchAddress("ESD",&fESD);
94 fKeyIterator = new TIter(fFile->GetListOfKeys());
99 if(AliKalmanTrack::GetConvConst()<=0.)
100 AliKalmanTrack::SetMagneticField(0.5);
101 if(fCurrentEvent>=fTree->GetEntries())
110 fTree->GetEvent(fCurrentEvent);
113 { // "old" way via ESD objects stored in root file
114 TKey* key = (TKey*)fKeyIterator->Next();
120 delete fFile; //we have to assume there are no more ESD objects in the fFile
124 TString esdname = "ESD";
125 esdname+=fCurrentEvent;
126 if(fESD) delete fESD;
127 fESD = dynamic_cast<AliESD*>(fFile->Get(esdname));
130 Info("ReadNext","Can't find AliESD object named %s",esdname.Data());
134 delete fFile;//we have to assume there is no more ESD objects in the fFile
142 return kTRUE;//success -> read one event
143 } while(fCurrentDir < GetNumberOfDirs());
144 //end of loop over directories specified in fDirs Obj Array
145 return kFALSE; //no more directories to read
148 /**********************************************************/
150 Int_t AliJetParticlesReaderESD::ReadESD(AliESD* esd)
156 Error("ReadESD","ESD is NULL");
160 //TDatabasePDG* pdgdb = TDatabasePDG::Instance();
163 // Error("ReadESD","Can not get PDG Database Instance.");
167 Float_t mf = esd->GetMagneticField();
170 Error("ReadESD","Magnetic Field is 0. Skipping to next event.");
173 AliKalmanTrack::SetMagneticField(mf/10.);
175 Info("ReadESD","Reading Event %d",fCurrentEvent);
176 if((!fOwner) || (fEventParticles==0))
177 fEventParticles = new AliJetEventParticles();
179 const Int_t kntr = esd->GetNumberOfTracks();
180 Info("ReadESD","Found %d tracks.",kntr);
181 fEventParticles->Reset(kntr);
185 headdesc+=esd->GetRunNumber();
187 headdesc+=esd->GetEventNumber();
188 fEventParticles->SetHeader(headdesc);
190 Double_t vertexpos[3];//vertex position
191 const AliESDVertex* kvertex = esd->GetVertex();
194 Info("ReadESD","ESD returned NULL pointer to vertex - assuming (0.0,0.0,0.0)");
201 kvertex->GetXYZ(vertexpos);
203 fEventParticles->SetVertex(vertexpos[0],vertexpos[1],vertexpos[2]);
206 for (Int_t i = 0;i<kntr; i++)
208 const AliESDtrack *kesdtrack = esd->GetTrack(i);
211 Error("ReadESD","Can't get track %d", i);
215 if ((kesdtrack->GetStatus() & fPassFlag) == kFALSE)
217 Info("ReadNext","Particle skipped: %ud.",kesdtrack->GetStatus());
221 Double_t mom[3]; //momentum
222 kesdtrack->GetPxPyPz(mom);
223 //kesdtrack->GetConstrainedPxPyPz(mom);
224 //Double_t pos[3];//position
225 //kesdtrack->GetXYZ(pos);
226 //kesdtrack->GetConstrainedXYZ(pos);
227 const Float_t kmass=kesdtrack->GetMass();
228 const Float_t kp2=mom[0]*mom[0]+mom[1]*mom[1]+mom[2]*mom[2];
229 const Float_t ketot=TMath::Sqrt(kmass*kmass+kp2);
230 const Float_t kpt=TMath::Sqrt(mom[0]*mom[0]+mom[1]*mom[1]);
231 const Float_t kp=TMath::Sqrt(kp2);
232 const Float_t keta=0.5*TMath::Log((kp+mom[2]+1e-30)/(kp-mom[2]+1e-30));
233 const Float_t kphi=TMath::Pi()+TMath::ATan2(-mom[1],-mom[0]);
234 if(IsAcceptedParticle(kpt,kphi,keta))
235 fEventParticles->AddParticle(mom[0],mom[1],mom[2],ketot,i,kesdtrack->GetLabel(),kpt,kphi,keta);
237 } // loop over tracks
242 /**********************************************************/
244 void AliJetParticlesReaderESD::Rewind()
248 if(fFile) delete fFile;
249 if(fKeyIterator) delete fKeyIterator;
256 /**********************************************************/
258 TFile* AliJetParticlesReaderESD::OpenFile(Int_t n)
260 //opens fFile with kine tree
262 const TString& dirname = GetDirName(n);
265 Error("OpenFiles","Can't get directory name");
268 TString filename = dirname +"/"+ fESDFileName;
269 TFile *ret = TFile::Open(filename.Data());
272 Error("OpenFiles","Can't open fFile %s",filename.Data());
277 Error("OpenFiles","Can't open fFile %s",filename.Data());