/* History of cvs commits:
*
* $Log$
+ * Revision 1.14 2006/08/28 10:01:56 kharlov
+ * Effective C++ warnings fixed (Timur Pocheptsov)
+ *
+ * Revision 1.13 2006/04/26 07:32:37 hristov
+ * Coding conventions, clean-up and related changes
+ *
+ * Revision 1.12 2006/03/10 13:23:36 hristov
+ * Using AliESDCaloCluster instead of AliESDtrack
+ *
+ * Revision 1.11 2006/01/31 20:30:52 hristov
+ * Including TFile.h
+ *
+ * Revision 1.10 2006/01/23 18:04:08 hristov
+ * Removing meaningless const
+ *
+ * Revision 1.9 2006/01/12 16:23:26 schutz
+ * ESD is properly read with methods of macros/AliReadESD.C copied in it
+ *
+ * Revision 1.8 2005/12/20 07:08:32 schutz
+ * corrected error in call AliReadESD
+ *
* Revision 1.6 2005/05/28 14:19:04 schutz
* Compilation warnings fixed by T.P.
*
// --- ROOT system ---
-#include "TParticle.h"
-#include "TCanvas.h"
-#include "TPaveLabel.h"
-#include "TPad.h"
+#include <TFile.h>
+#include <TParticle.h>
+#include <TCanvas.h>
+#include <TPaveLabel.h>
+#include <TPad.h>
+#include <TH2.h>
+
#include "AliPHOSGammaJet.h"
#include "AliPHOSGetter.h"
-#include "TH2.h"
#include "AliPHOSGeometry.h"
#include "AliPHOSFastGlobalReconstruction.h"
#include "AliESD.h"
#include "AliESDtrack.h"
+#include "AliESDCaloCluster.h"
#include "Riostream.h"
-#include "macros/AliReadESD.C"
-//#include "../PYTHIA6/AliGenPythia.h"
+
ClassImp(AliPHOSGammaJet)
//____________________________________________________________________________
-AliPHOSGammaJet::AliPHOSGammaJet() : TTask()
+AliPHOSGammaJet::AliPHOSGammaJet() :
+ TTask(), fAnyConeOrPt(0), fOptionGJ(""),
+ fOutputFile(new TFile(gDirectory->GetName())),
+ fOutputFileName(gDirectory->GetName()),
+ fInputFileName(gDirectory->GetName()),
+ fHIJINGFileName(gDirectory->GetName()),
+ fHIJING(0), fESDdata(0), fEtaCut(0.),
+ fOnlyCharged(0), fPhiMaxCut(0.),
+ fPhiMinCut(0.), fPtCut(0.),
+ fNeutralPtCut(0.), fChargedPtCut(0.),
+ fInvMassMaxCut(0.), fInvMassMinCut(0.),
+ fMinDistance(0.), fRatioMaxCut(0.), fRatioMinCut(0.),
+ fTPCCutsLikeEMCAL(0), fDirName(""), fESDTree(""),
+ fPattern(""), fJetTPCRatioMaxCut(0.),
+ fJetTPCRatioMinCut(0.), fJetRatioMaxCut(0.),
+ fJetRatioMinCut(0.), fNEvent(0), fNCone(0),
+ fNPt(0), fCone(0), fPtThreshold(0),
+ fPtJetSelectionCut(0.0),
+ fListHistos(new TObjArray(100)),
+ fFastRec(0), fOptFast(0),
+ fRan(0), fResPara1(0.), fResPara2(0.), fResPara3(0.),
+ fPosParaA(0.), fPosParaB(0.), fAngleMaxParam(), fSelect(0)
{
// ctor
fAngleMaxParam.Set(4) ;
fAngleMaxParam.Reset(0.);
- fAnyConeOrPt = 0 ;
- fEtaCut = 0. ;
- fESDdata = 0 ;
- fFastRec = 0 ;
- fInvMassMaxCut = 0. ;
- fInvMassMinCut = 0. ;
- fJetRatioMaxCut = 0. ;
- fJetRatioMinCut = 0. ;
- fJetTPCRatioMaxCut = 0. ;
- fJetTPCRatioMinCut = 0. ;
- fMinDistance = 0. ;
- fNEvent = 0 ;
- fNCone = 0 ;
- fNPt = 0 ;
- fOnlyCharged = 0 ;
- fOptFast = 0 ;
- fPhiMaxCut = 0. ;
- fPhiMinCut = 0. ;
- fPtCut = 0. ;
- fNeutralPtCut = 0. ;
- fChargedPtCut = 0. ;
- fRatioMaxCut = 0. ;
- fRatioMinCut = 0. ;
- fCone = 0 ;
- fPtThreshold = 0 ;
- fTPCCutsLikeEMCAL = 0 ;
- fSelect = 0 ;
-
- fDirName = "" ;
- fESDTree = "" ;
- fPattern = "" ;
for(Int_t i = 0; i<10; i++){
fCones[i] = 0.0 ;
fNameCones[i] = "" ;
fPtThres[i] = 0.0 ;
- fPtJetSelectionCut = 0.0 ;
fNamePtThres[i] = "" ;
if( i < 6 ){
fJetXMin1[i] = 0.0 ;
}
}
}
-
- fOptionGJ = "" ;
- fOutputFile = new TFile(gDirectory->GetName()) ;
- fInputFileName = gDirectory->GetName() ;
- fOutputFileName = gDirectory->GetName() ;
- fHIJINGFileName = gDirectory->GetName() ;
- fHIJING = 0 ;
- fPosParaA = 0. ;
- fPosParaB = 0. ;
- fRan = 0 ;
- fResPara1 = 0. ;
- fResPara2 = 0. ;
- fResPara3 = 0. ;
- fListHistos = new TObjArray(100) ;
TList * list = gDirectory->GetListOfKeys() ;
TIter next(list) ;
TH2F * h = 0 ;
//____________________________________________________________________________
AliPHOSGammaJet::AliPHOSGammaJet(const TString inputfilename) :
- TTask("GammaJet","Analysis of gamma-jet correlations")
+ TTask("GammaJet","Analysis of gamma-jet correlations"),
+ fAnyConeOrPt(0), fOptionGJ(),
+ fOutputFile(0),
+ fOutputFileName(),
+ fInputFileName(),
+ fHIJINGFileName(),
+ fHIJING(0), fESDdata(0), fEtaCut(0.),
+ fOnlyCharged(0), fPhiMaxCut(0.),
+ fPhiMinCut(0.), fPtCut(0.),
+ fNeutralPtCut(0.), fChargedPtCut(0.),
+ fInvMassMaxCut(0.), fInvMassMinCut(0.),
+ fMinDistance(0.), fRatioMaxCut(0.), fRatioMinCut(0.),
+ fTPCCutsLikeEMCAL(0), fDirName(), fESDTree(),
+ fPattern(), fJetTPCRatioMaxCut(0.),
+ fJetTPCRatioMinCut(0.), fJetRatioMaxCut(0.),
+ fJetRatioMinCut(0.), fNEvent(0), fNCone(0),
+ fNPt(0), fCone(0), fPtThreshold(0),
+ fPtJetSelectionCut(0.0),
+ fListHistos(0),
+ fFastRec(0), fOptFast(0),
+ fRan(0), fResPara1(0.), fResPara2(0.), fResPara3(0.),
+ fPosParaA(0.), fPosParaB(0.), fAngleMaxParam(), fSelect(0)
+
{
// ctor
-// gSystem->Load("$ALICE_ROOT/PHOS/macros/AliReadESD_C.so");
fInputFileName = inputfilename;
fFastRec = new AliPHOSFastGlobalReconstruction(fInputFileName);
AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileName) ;
fNEvent = gime->MaxEvent();
InitParameters();
- fListHistos = 0 ;
}
//____________________________________________________________________________
-AliPHOSGammaJet::AliPHOSGammaJet(const AliPHOSGammaJet & gj) : TTask(gj)
+AliPHOSGammaJet::AliPHOSGammaJet(const AliPHOSGammaJet & gj) :
+ TTask(gj),
+ fAnyConeOrPt(gj.fAnyConeOrPt), fOptionGJ(gj.fOptionGJ),
+ fOutputFile(gj.fOutputFile),
+ fOutputFileName(gj.fOutputFileName),
+ fInputFileName(gj.fInputFileName),
+ fHIJINGFileName(gj.fHIJINGFileName),
+ fHIJING(gj.fHIJING), fESDdata(gj.fESDdata), fEtaCut(gj.fEtaCut),
+ fOnlyCharged(gj.fOnlyCharged), fPhiMaxCut(gj.fPhiMaxCut),
+ fPhiMinCut(gj.fPhiMinCut), fPtCut(gj.fPtCut),
+ fNeutralPtCut(gj.fNeutralPtCut), fChargedPtCut(gj.fChargedPtCut),
+ fInvMassMaxCut(gj.fInvMassMaxCut), fInvMassMinCut(gj.fInvMassMinCut),
+ fMinDistance(gj.fMinDistance), fRatioMaxCut(gj.fRatioMaxCut),
+ fRatioMinCut(gj.fRatioMinCut), fTPCCutsLikeEMCAL(gj.fTPCCutsLikeEMCAL),
+ fDirName(gj.fDirName), fESDTree(gj.fESDTree),
+ fPattern(gj.fPattern), fJetTPCRatioMaxCut(gj.fJetTPCRatioMaxCut),
+ fJetTPCRatioMinCut(gj.fJetTPCRatioMinCut), fJetRatioMaxCut(gj.fJetRatioMaxCut),
+ fJetRatioMinCut(gj.fJetRatioMinCut), fNEvent(gj.fNEvent), fNCone(gj.fNCone),
+ fNPt(gj.fNPt), fCone(gj.fCone), fPtThreshold(gj.fPtThreshold),
+ fPtJetSelectionCut(gj.fPtJetSelectionCut),
+ fListHistos(0),//?????
+ fFastRec(gj.fFastRec), fOptFast(gj.fOptFast),
+ fRan(0), //???
+ fResPara1(gj.fResPara1), fResPara2(gj.fResPara2), fResPara3(gj.fResPara3),
+ fPosParaA(gj.fPosParaA), fPosParaB(gj.fPosParaB),
+ fAngleMaxParam(gj.fAngleMaxParam), fSelect(gj.fSelect)
{
// cpy ctor
- fAngleMaxParam = gj.fAngleMaxParam;
- fAnyConeOrPt = gj.fAnyConeOrPt;
- fESDdata = gj.fESDdata;
- fEtaCut = gj.fEtaCut ;
- fInvMassMaxCut = gj.fInvMassMaxCut ;
- fInvMassMinCut = gj.fInvMassMinCut ;
- fFastRec = gj.fFastRec ;
- fOptionGJ = gj.fOptionGJ ;
- fMinDistance = gj.fMinDistance ;
- fOptFast = gj.fOptFast ;
- fOnlyCharged = gj.fOnlyCharged ;
- fOutputFile = gj.fOutputFile ;
- fInputFileName = gj.fInputFileName ;
- fOutputFileName = gj.fOutputFileName ;
- fHIJINGFileName = gj.fHIJINGFileName ;
- fHIJING = gj.fHIJING ;
- fRatioMaxCut = gj.fRatioMaxCut ;
- fRatioMinCut = gj.fRatioMinCut ;
- fJetRatioMaxCut = gj.fJetRatioMaxCut ;
- fJetRatioMinCut = gj.fJetRatioMinCut ;
- fJetTPCRatioMaxCut = gj.fJetRatioMaxCut ;
- fJetTPCRatioMinCut = gj.fJetRatioMinCut ;
- fNEvent = gj.fNEvent ;
- fNCone = gj.fNCone ;
- fNPt = gj.fNPt ;
- fResPara1 = gj.fResPara1 ;
- fResPara2 = gj.fResPara2 ;
- fResPara3 = gj.fResPara3 ;
- fPtCut = gj.fPtCut ;
- fNeutralPtCut = gj.fNeutralPtCut ;
- fChargedPtCut = gj.fChargedPtCut ;
- fPtJetSelectionCut = gj.fPtJetSelectionCut ;
- fPhiMaxCut = gj.fPhiMaxCut ;
- fPhiMinCut = gj.fPhiMinCut ;
- fPosParaA = gj.fPosParaA ;
- fPosParaB = gj.fPosParaB ;
- fSelect = gj.fSelect ;
- fTPCCutsLikeEMCAL = gj.fTPCCutsLikeEMCAL ;
- fCone = gj.fCone ;
- fPtThreshold = gj.fPtThreshold ;
-
- fDirName = gj.fDirName ;
- fESDTree = gj.fESDTree ;
- fPattern = gj.fPattern ;
-
SetName (gj.GetName()) ;
SetTitle(gj.GetTitle()) ;
AliPHOSGammaJet::~AliPHOSGammaJet()
{
fOutputFile->Close() ;
-
}
//____________________________________________________________________________
AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
+//DP. To be fixed: extract vertex position
+ Double_t vtx[3]={0.,0.,0.} ;
+
if(!fOptFast){
for (iParticle=0 ; iParticle < t->GetEntries() ; iParticle++) {
t->GetEvent(iParticle) ;
}
}
else if((charge == 0) && (particle->Pt() > fNeutralPtCut) ){
- geom->ImpactOnEmc(particle->Theta(),particle->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,particle->Theta(),particle->Phi(), m,z,x);
if(m != 0)
{//Is in PHOS
indexCh++ ;
}
else if(charge == 0){
- geom->ImpactOnEmc(particle->Theta(),particle->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,particle->Theta(),particle->Phi(), m,z,x);
if((particle->GetPdgCode() != 111) && (particle->Pt() > fNeutralPtCut) &&
(TMath::Abs(particle->Eta())<fEtaCut) ){
TLorentzVector part(particle->Px(),particle->Py(),
TParticle * photon1 =
new TParticle(22,1,0,0,0,0,pGamma1.Px(),pGamma1.Py(),
pGamma1.Pz(),pGamma1.E(),0,0,0,0);
- geom->ImpactOnEmc(photon1->Theta(),photon1->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,photon1->Theta(),photon1->Phi(), m,z,x);
if( photon1->Phi()>fPhiEMCALCut[0] && photon1->Phi()<fPhiEMCALCut[1]
&& m == 0){
if(strstr(fOptionGJ,"deb all") || strstr(fOptionGJ,"deb"))
TParticle * photon2 =
new TParticle(22,1,0,0,0,0,pGamma2.Px(), pGamma2.Py(),
pGamma2.Pz(),pGamma2.E(),0,0,0,0);
- geom->ImpactOnEmc(photon2->Theta(),photon2->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,photon2->Theta(),photon2->Phi(), m,z,x);
if(photon2->Phi()>fPhiEMCALCut[0] &&
photon2->Phi()<fPhiEMCALCut[1] && m == 0){
if(strstr(fOptionGJ,"deb all") || strstr(fOptionGJ,"deb"))
const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
gime->Event(iEvent, "X") ;
+ //DP to be fixed: extract true vertex here
+ Double_t vtx[3]={0.,0.,0.} ;
Int_t index = particleList->GetEntries() ;
Int_t indexCh = plCh->GetEntries() ;
new((*particleList)[index++]) TParticle(*particle) ;
}
else if((charge == 0) && (particle->Pt() > fNeutralPtCut)){
- geom->ImpactOnEmc(particle->Theta(),particle->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,particle->Theta(),particle->Phi(), m,z,x);
if(m != 0)
{//Is in PHOS
if(strstr(fOptionGJ,"deb all")|| strstr(fOptionGJ,"deb"))
new((*particleList)[index++]) TParticle(*particle) ;
}
else if(charge == 0) {
- geom->ImpactOnEmc(particle->Theta(),particle->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,particle->Theta(),particle->Phi(), m,z,x);
if((particle->GetPdgCode() != 111) && particle->Pt() > 0 &&
(TMath::Abs(particle->Eta())<fEtaCut))
{
TParticle * photon1 =
new TParticle(22,1,0,0,0,0,pGamma1.Px(),pGamma1.Py(),
pGamma1.Pz(),pGamma1.E(),0,0,0,0);
- geom->ImpactOnEmc(photon1->Theta(),photon1->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,photon1->Theta(),photon1->Phi(), m,z,x);
if( photon1->Phi()>fPhiEMCALCut[0] && photon1->Phi()<fPhiEMCALCut[1]
&& m == 0){
if(strstr(fOptionGJ,"deb all") || strstr(fOptionGJ,"deb"))
TParticle * photon2 =
new TParticle(22,1,0,0,0,0,pGamma2.Px(), pGamma2.Py(),
pGamma2.Pz(),pGamma2.E(),0,0,0,0);
- geom->ImpactOnEmc(photon2->Theta(),photon2->Phi(), m,z,x);
+ geom->ImpactOnEmc(vtx,photon2->Theta(),photon2->Phi(), m,z,x);
if(photon2->Phi()>fPhiEMCALCut[0] &&
photon2->Phi()<fPhiEMCALCut[1] && m == 0){
if(strstr(fOptionGJ,"deb all") || strstr(fOptionGJ,"deb"))
TClonesArray * plNe,
TClonesArray * plNePHOS,
const AliESD * esd){
- //, Int_t iEvent){
+
+ //Create a list of particles from the ESD. These particles have been measured
+ //by the Central Tracking system (TPC), PHOS and EMCAL
+ //(EMCAL not available for the moment).
//Info("CreateParticleListFromESD","Inside");
- //AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileName) ;
- //gime->Event(iEvent, "P") ;
- Double_t v[3] ; //vertex ;
- esd->GetVertex()->GetXYZ(v) ;
Int_t index = pl->GetEntries() ;
+ Int_t npar = 0 ;
+ Float_t *pid = new Float_t[AliPID::kSPECIESN];
+
+ //########### PHOS ##############
+ //Info("CreateParticleListFromESD","Fill ESD PHOS list");
+ Int_t begphos = esd->GetFirstPHOSCluster();
+ Int_t endphos = esd->GetFirstPHOSCluster() +
+ esd->GetNumberOfPHOSClusters() ;
+ Int_t indexNePHOS = plNePHOS->GetEntries() ;
+ if(strstr(fOptionGJ,"deb all"))
+ Info("CreateParticleListFromESD","PHOS: first particle %d, last particle %d",
+ begphos,endphos);
- //PHOS
- Int_t begphos = esd->GetFirstPHOSParticle();
- Int_t endphos = esd->GetFirstPHOSParticle() + esd->GetNumberOfPHOSParticles() ;
+ for (npar = begphos; npar < endphos; npar++) {//////////////PHOS track loop
+ AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd
+
+ //Create a TParticle to fill the particle list
- Int_t indexNePHOS = plNePHOS->GetEntries() ;
- Int_t nphP ;
- //cout<<"Fill PHOS"<<endl;
- //cout<<"begin "<<begphos<<" end "<<endphos<<endl;
- for (nphP = begphos; nphP < endphos; nphP++) {//////////////PHOS track loop
- AliESDtrack * track = esd->GetTrack(nphP) ; // retrieve track from esd
-
- // PID of this track
- Double_t pid[AliPID::kSPECIESN];
- track->GetPHOSpid(pid);
- TParticle * particle = new TParticle() ;
- Double_t en = track->GetPHOSsignal() ;
- Double_t * p = new Double_t();
- track->GetPHOSposition(p) ;
+ Float_t en = clus->GetClusterEnergy() ;
+ Float_t *p = new Float_t();
+ clus->GetGlobalPosition(p) ;
TVector3 pos(p[0],p[1],p[2]) ;
Double_t phi = pos.Phi();
Double_t theta= pos.Theta();
Double_t px = en*TMath::Cos(phi)*TMath::Sin(theta);;
Double_t py = en*TMath::Sin(phi)*TMath::Sin(theta);
Double_t pz = en*TMath::Cos(theta);
+
+ TParticle * particle = new TParticle() ;
particle->SetMomentum(px,py,pz,en) ;
+
+ //Select only photons
+
+ pid=clus->GetPid();
//cout<<"pid "<<pid[AliPID::kPhoton]<<endl ;
if( pid[AliPID::kPhoton] > 0.75)
- new((*plNePHOS)[indexNePHOS++]) TParticle(*particle) ;
+ new((*plNePHOS)[indexNePHOS++]) TParticle(*particle) ;
}
- //TPC
+
+ //########### TPC #####################
+ //Info("CreateParticleListFromESD","Fill ESD TPC list");
Int_t begtpc = 0 ;
Int_t endtpc = esd->GetNumberOfTracks() ;
Int_t indexCh = plCh->GetEntries() ;
- //cout<<"Fill TPC"<<endl;
- //cout<<"begin "<<begtpc<<" end "<<endtpc<<endl;
- for (nphP = begtpc; nphP < endtpc; nphP++) {////////////// track loop
- AliESDtrack * track = esd->GetTrack(nphP) ; // retrieve track from esd
+ if(strstr(fOptionGJ,"deb all"))
+ Info("CreateParticleListFromESD","TPC: first particle %d, last particle %d",
+ begtpc,endtpc);
+
+ for (npar = begtpc; npar < endtpc; npar++) {////////////// track loop
+ AliESDtrack * track = esd->GetTrack(npar) ; // retrieve track from esd
- // PID of this track
- //Double_t pid[AliPID::kSPECIESN];
- //track->GetEMCALpid(pid);
+ Double_t en = track ->GetTPCsignal() ;
+ Double_t mom[3];
+ track->GetPxPyPz(mom) ;
+ Double_t px = mom[0];
+ Double_t py = mom[1];
+ Double_t pz = mom[2]; //Check with TPC people if this is correct.
+
+ //cout<<"TPC signal "<<en<<endl;
+ //cout<<"px "<<px<<"; py "<<py<<"; pz "<<pz<<endl;
TParticle * particle = new TParticle() ;
- Double_t en = track->GetTPCsignal() ;
- //Double_t *p = 0;
- //track->GetEMCALposition(p) ;
- //TVector3 pos(p[0],p[1],p[2]) ;
- //Double_t phi = pos.Phi();
- //Double_t theta= pos.Theta();
-
- TVector3 mom = track->P3() ;
- Double_t px = mom.Px();
- Double_t py = mom.Py();
- Double_t pz = mom.Pz();
-
particle->SetMomentum(px,py,pz,en) ;
- //if( pid[AliPID::kPhoton] > 0.75)
+
new((*plCh)[indexCh++]) TParticle(*particle) ;
new((*pl)[index++]) TParticle(*particle) ;
}
- //EMCAL
- Int_t begem = esd->GetFirstEMCALParticle();
- Int_t endem = esd->GetFirstEMCALParticle() + esd->GetNumberOfEMCALParticles() ;
- Int_t indexNe = plNe->GetEntries() ;
- //cout<<"Fill EMCAL"<<endl;
- //cout<<"begin "<<begem<<" end "<<endem<<endl;
- for (nphP = begem; nphP < endem; nphP++) {//////////////PHOS track loop
- AliESDtrack * track = esd->GetTrack(nphP) ; // retrieve track from esd
-
- // PID of this track
- Double_t pid[AliPID::kSPECIESN];
- track->GetEMCALpid(pid);
- Double_t en = track->GetEMCALsignal() ;
- //cout<<"EMCAL signal "<<en<<endl;
- Double_t *p = new Double_t();
- track->GetEMCALposition(p) ;
- TVector3 pos(p[0],p[1],p[2]) ;
+ //################ EMCAL ##############
+ Double_t v[3] ; //vertex ;
+ esd->GetVertex()->GetXYZ(v) ;
+ //##########Uncomment when ESD for EMCAL works ##########
+ //Info("CreateParticleListFromESD","Fill ESD EMCAL list");
+
+ Int_t begem = esd->GetFirstEMCALCluster();
+ Int_t endem = esd->GetFirstEMCALCluster() +
+ esd->GetNumberOfEMCALClusters() ;
+ Int_t indexNe = plNe->GetEntries() ;
+ if(strstr(fOptionGJ,"deb all"))
+ Info("CreateParticleListFromESD","EMCAL: first particle %d, last particle %d",
+ begem,endem);
+
+ for (npar = begem; npar < endem; npar++) {//////////////EMCAL track loop
+ AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd
+
+ Float_t en = clus->GetClusterEnergy() ;
+ Float_t *p = new Float_t();
+ clus->GetGlobalPosition(p) ;
+ TVector3 pos(p[0],p[1],p[2]) ;
Double_t phi = pos.Phi();
Double_t theta= pos.Theta();
-
Double_t px = en*TMath::Cos(phi)*TMath::Sin(theta);;
Double_t py = en*TMath::Sin(phi)*TMath::Sin(theta);
Double_t pz = en*TMath::Cos(theta);
- Int_t pdg = 0;
+ //cout<<"EMCAL signal "<<en<<endl;
+ //cout<<"px "<<px<<"; py "<<py<<"; pz "<<pz<<endl;
+ //TParticle * particle = new TParticle() ;
//particle->SetMomentum(px,py,pz,en) ;
- if( pid[AliPID::kPhoton] > 0.75) //This has to be fixen.
- pdg = 22;
- //particle->SetPdgCode(22);
- else if( pid[AliPID::kPi0] > 0.75)
- pdg = 111;
- // particle->SetPdgCode(111);
- pdg = 22 ; //Assume all photons, no PID for EMCAL available.
+
+ Int_t pdg = 0;
+ // //Uncomment if PID IS WORKING, photon and pi0 idenitification.
+ // //if( pid[AliPID::kPhoton] > 0.75) //This has to be fixen.
+ // //pdg = 22;
+ // //else if( pid[AliPID::kPi0] > 0.75)
+ // //pdg = 111;
+ pdg = 22; //No PID, assume all photons
TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
px, py, pz, en, v[0], v[1], v[2], 0);
- //cout<<"pid; ph "<<pid[AliPID::kPhoton]<<" pi0 "<<pid[AliPID::kPi0]<<endl;
+
new((*plNe)[indexNe++]) TParticle(*particle) ;
new((*pl)[index++]) TParticle(*particle) ;
}
- //Info("CreateParticleListFromESD","End Inside");
-
+
+ // Info("CreateParticleListFromESD","End Inside");
+
}
if(fESDdata){
// Create chain of esd trees
- const UInt_t kNevent = static_cast<const UInt_t>(GetNEvent()) ;
- t = AliReadESD(kNevent, fDirName, fESDTree, fPattern) ;
+ const UInt_t kNevent = static_cast<UInt_t>(GetNEvent()) ;
+ t = ReadESD(kNevent, fDirName, fESDTree, fPattern) ;
if(!t) {
AliError("Could not create the TChain") ;
//break ;
}
+//__________________________________________________________________________
+TChain * AliPHOSGammaJet::ReadESDfromdisk(const UInt_t eventsToRead,
+ const TString dirName,
+ const TString esdTreeName,
+ const char * pattern)
+{
+ // Reads ESDs from Disk
+ TChain * rv = 0;
+
+ AliInfo( Form("\nReading files in %s \nESD tree name is %s \nReading %d events",
+ dirName.Data(), esdTreeName.Data(), eventsToRead) ) ;
+
+ // create a TChain of all the files
+ TChain * cESDTree = new TChain(esdTreeName) ;
+
+ // read from the directory file until the require number of events are collected
+ void * from = gSystem->OpenDirectory(dirName) ;
+ if (!from) {
+ AliError( Form("Directory %s does not exist") ) ;
+ rv = 0 ;
+ }
+ else{ // reading file names from directory
+ const char * subdir ;
+ // search all subdirectories witch matching pattern
+ while( (subdir = gSystem->GetDirEntry(from)) &&
+ (cESDTree->GetEntries() < eventsToRead)) {
+ if ( strstr(subdir, pattern) != 0 ) {
+ char file[200] ;
+ sprintf(file, "%s%s/AliESDs.root", dirName.Data(), subdir);
+ AliInfo( Form("Adding %s\n", file) );
+ cESDTree->Add(file) ;
+ }
+ } // while file
+
+ AliInfo( Form(" %d events read", cESDTree->GetEntriesFast()) ) ;
+ rv = cESDTree ;
+
+ } // reading file names from directory
+ return rv ;
+}
+
+//__________________________________________________________________________
+TChain * AliPHOSGammaJet::ReadESD(const UInt_t eventsToRead,
+ const TString dirName,
+ const TString esdTreeName,
+ const char * pattern)
+{
+ // Read AliESDs files and return a Chain of events
+
+ if ( dirName == "" ) {
+ AliError("Give the name of the DIR where to find files") ;
+ return 0 ;
+ }
+ if ( esdTreeName == "" )
+ return ReadESDfromdisk(eventsToRead, dirName) ;
+ else if ( strcmp(pattern, "") == 0 )
+ return ReadESDfromdisk(eventsToRead, dirName, esdTreeName) ;
+ else
+ return ReadESDfromdisk(eventsToRead, dirName, esdTreeName, pattern) ;
+}
+
//___________________________________________________________________
void AliPHOSGammaJet::SetJet(TParticle * part, Bool_t & b, Float_t cone,
Double_t eta, Double_t phi)