// Base class for analysis algorithms
//-- Author: Gustavo Conesa (LNF-INFN)
//_________________________________________________________________________
+// --Yaxian Mao: Add the possibality for event selection analysis based on vertex and multiplicity bins (10/10/2010)
// --- ROOT system ---
AliAnaPartCorrBaseClass::AliAnaPartCorrBaseClass() :
TObject(), fDataMC(0), fDebug(0), fCheckFidCut(0),
fCheckCaloPID(0), fRecalculateCaloPID(0), fMinPt(0), fMaxPt(0),
+ fMultiBin(0),fNZvertBin(0),fNrpBin(0),fZvtxCut(0.), fMaxMulti(0),fMinMulti(0),
+ fUseSelectEvent(kFALSE),
fReader(0x0), fInputAODBranch(0x0), fInputAODName(""),
fOutputAODBranch(0x0), fNewAOD(kFALSE),
fOutputAODName(""), fOutputAODClassName(""),
if(fIC) delete fIC ;
if(fMCUtils) delete fMCUtils ;
if(fNMS) delete fNMS ;
-
-// printf("--- analysis deleted \n");
+
+ // printf("--- analysis deleted \n");
}
//____________________________________________________________________________
fRecalculateCaloPID = kFALSE ;
fMinPt = 0.1 ; //Min pt in particle analysis
fMaxPt = 300. ; //Max pt in particle analysis
-
+ fMultiBin = 1;
+ fNZvertBin = 1;
+ fNrpBin = 1;
+ fZvtxCut = 40;
+ fMaxMulti = 1000;
+ fMinMulti = 0;
+ fUseSelectEvent = kFALSE ;
+
//fReader = new AliCaloTrackReader(); //Initialized in maker
//fCaloUtils = new AliCalorimeterUtils();//Initialized in maker
//_________________________________________________________________________
// Base class for analysis algorithms
//-- Author: Gustavo Conesa (INFN-LNF)
-
+//-Add the possibality for event selection analysis based on vertex and multiplicity bins (Yaxian Mao, 10/10/2010)
#include <cstdlib>
//ROOT
virtual void SetMinPt(Float_t pt) {fMinPt = pt ; }
virtual void SetPtCutRange(Double_t ptmin, Double_t ptmax)
{ fMaxPt=ptmax; fMinPt=ptmin;}
+ //Setters for parameters of event buffers
+ virtual void SetMultiBin(Int_t n=1) {fMultiBin=n ;} //number of bins in Multiplicity
+ virtual void SetNZvertBin(Int_t n=1) {fNZvertBin=n ;} //number of bins for vertex position
+ virtual void SetNRPBin(Int_t n=1) {fNrpBin=n ;} //number of bins in reaction plain
+ virtual void SetZvertexCut(Float_t zcut=40.){fZvtxCut=zcut ;} //cut on vertex position
+ virtual void SetMultiplicity(Int_t multimin, Int_t multimax) {fMinMulti = multimin ; fMaxMulti = multimax ; }
+ virtual void SwitchOnEventSelection() {fUseSelectEvent = kTRUE ; }
+ virtual void SwitchOffEventSelection() {fUseSelectEvent = kFALSE ; }
+ //Getters for event selection
+ virtual Int_t GetMultiBin() const {return fMultiBin ;} //number of bins in Multiplicity
+ virtual Int_t GetNZvertBin() const {return fNZvertBin ;} //number of bins in vertex
+ virtual Int_t GetNRPBin() const {return fNrpBin ;} //number of bins in reaction plain
+ //Getters for event selection
+ virtual Float_t GetZvertexCut() const {return fZvtxCut ;} //cut on vertex position
+ virtual Int_t GetMaxMulti() const {return fMaxMulti ; }
+ virtual Int_t GetMinMulti() const {return fMinMulti ; }
+
+ // Do correlation analysis with different event buffers
+ virtual Bool_t DoEventSelect() const {return fUseSelectEvent ; }
//Histogrammes setters and getters
//Pt, Energy
Bool_t fRecalculateCaloPID ; // Recalculate PID or use PID weights in calorimeters
Float_t fMinPt ; // Maximum pt of (trigger) particles in the analysis
Float_t fMaxPt ; // Minimum pt of (trigger) particles in the analysis
+ Int_t fMultiBin ; // Number of bins in event container for multiplicity
+ Int_t fNZvertBin ; // Number of bins in event container for vertex position
+ Int_t fNrpBin ; // Number of bins in event container for reaction plain
+ Float_t fZvtxCut ; // Cut on vertex position
+ Int_t fMaxMulti ; // Maximum multiplicity of particles in the analysis
+ Int_t fMinMulti ; // Maximum multiplicity of particles in the analysis
+ Bool_t fUseSelectEvent ; // Select events based on multiplicity and vertex cuts
+
AliCaloTrackReader * fReader; // Acces to ESD/AOD/MC data
AliMCAnalysisUtils * fMCUtils; //! MonteCarlo Analysis utils
AliNeutralMesonSelection * fNMS; //! Neutral Meson Selection
AliCalorimeterUtils * fCaloUtils ; // Pointer to CalorimeterUtils
-
//Histograms binning and range
Int_t fHistoPtBins ; // Number of bins in pt axis
Float_t fHistoAsymMax ; // Maximum value of asymmetry histogram range
Float_t fHistoAsymMin ; // Minimum value of asymmetry histogram range
- ClassDef(AliAnaPartCorrBaseClass,9)
+ ClassDef(AliAnaPartCorrBaseClass,10)
} ;
//
//
//*-- Author: Gustavo Conesa (LNF-INFN)
+
+//-Yaxian Mao (add the possibility for different IC method with different pt range, 01/10/2010)
+//-Yaxian Mao (check the candidate particle is the leading particle or not at the same hemishere)
+
//////////////////////////////////////////////////////////////////////////////
//____________________________________________________________________________
AliIsolationCut::AliIsolationCut() :
TObject(),
- fConeSize(0.),fPtThreshold(0.), fPtFraction(0.), fICMethod(0),fPartInCone(0)
+ fConeSize(0.),fPtThreshold(0.), fSumPtThreshold(0.), fPtFraction(0.), fICMethod(0),fPartInCone(0)
{
//default ctor
fConeSize = 0.4 ;
fPtThreshold = 1. ;
+ fSumPtThreshold = 0.5 ;
fPtFraction = 0.1 ;
fPartInCone = kNeutralAndCharged;
fICMethod = kPtThresIC; // 0 pt threshol method, 1 cone pt sum method
n = 0 ;
coneptsum = 0.;
isolated = kFALSE;
-
+
//Initialize the array with refrences
TObjArray * refclusters = 0x0;
TObjArray * reftracks = 0x0;
Int_t ntrackrefs = 0;
Int_t nclusterrefs = 0;
-
//Check charged particles in cone.
if(plCTS && (fPartInCone==kOnlyCharged || fPartInCone==kNeutralAndCharged)){
TVector3 p3;
phi = p3.Phi() ;
if(phi<0) phi+=TMath::TwoPi();
+ //only loop the particle at the same side of candidate
+ if(TMath::Abs(phi-phiC)>TMath::PiOver2()) continue ;
+ //if at the same side has particle larger than candidate, then candidate can not be the leading, skip such events
+ if(pt > ptC){
+ n = -1;
+ nfrac = -1;
+ coneptsum = -1;
+ isolated = kFALSE;
+ if(fillAOD && reftracks) reftracks->Clear();
+ return ;
+ }
//Check if there is any particle inside cone with pt larger than fPtThreshold
rad = TMath::Sqrt((eta-etaC)*(eta-etaC)+ (phi-phiC)*(phi-phiC));
eta = mom.Eta();
phi = mom.Phi() ;
if(phi<0) phi+=TMath::TwoPi();
+ //only loop the particle at the same side of candidate
+
+ if(TMath::Abs(phi-phiC)>TMath::PiOver2()) continue ;
+ //if at the same side has particle larger than candidate, then candidate can not be the leading, skip such events
+ if(pt > ptC){
+ n = -1;
+ nfrac = -1;
+ coneptsum = -1;
+ isolated = kFALSE;
+ if(fillAOD){
+ if(reftracks) reftracks->Clear();
+ if(refclusters)refclusters->Clear();
+ }
+ return ;
+ }
//Check if there is any particle inside cone with pt larger than fPtThreshold
rad = TMath::Sqrt((eta-etaC)*(eta-etaC)+ (phi-phiC)*(phi-phiC));
//printf("neutral in isolation cone pt %f, phi %f, eta %f, R %f \n",pt,phi,eta,rad);
coneptsum+=pt;
if(pt > fPtThreshold ) n++;
- if(pt > fPtFraction*ptC ) nfrac++;
+ //if fPtFraction*ptC<fPtThreshold then consider the fPtThreshold directly
+ if(pt > fPtFraction*ptC && pt>fPtThreshold) nfrac++;
}//in cone
}// neutral particle loop
}//neutrals
-
+
//printf("Isolation Cut: in cone with: pT>pTthres %d, pT > pTfrac*pTcandidate %d \n",n,nfrac);
//Add reference arrays to AOD when filling AODs only
if(refclusters) pCandidate->AddObjArray(refclusters);
if(reftracks) pCandidate->AddObjArray(reftracks);
}
-
//Check isolation, depending on method.
if( fICMethod == kPtThresIC){
if(n==0) isolated = kTRUE ;
}
else if( fICMethod == kSumPtIC){
- if(coneptsum < fPtThreshold)
+ if(coneptsum < fSumPtThreshold)
isolated = kTRUE ;
}
else if( fICMethod == kPtFracIC){
if(nfrac==0) isolated = kTRUE ;
}
else if( fICMethod == kSumPtFracIC){
- if(coneptsum < fPtFraction*ptC)
- isolated = kTRUE ;
+ //when the fPtFraction*ptC < fSumPtThreshold then consider the later case
+ if(coneptsum < fPtFraction*ptC && coneptsum < fSumPtThreshold) isolated = kTRUE ;
}
//if(refclusters) delete refclusters;
//
// -- Author: Gustavo Conesa (INFN-LNF)
+//-Yaxian Mao (add the possibility for different IC method with different pt range, 01/10/2010)
+
// --- ROOT system ---
#include <TObject.h>
class TObjArray ;
Float_t GetConeSize() const {return fConeSize ;}
Float_t GetPtThreshold() const {return fPtThreshold;}
+ Float_t GetSumPtThreshold() const {return fSumPtThreshold;}
Float_t GetPtFraction() const {return fPtFraction ;}
Int_t GetICMethod() const {return fICMethod ;}
Int_t GetParticleTypeInCone() const {return fPartInCone ;}
void SetConeSize(Float_t r) {fConeSize = r ; }
void SetPtThreshold(Float_t pt) {fPtThreshold = pt; }
+ void SetSumPtThreshold(Float_t sumpt) {fSumPtThreshold = sumpt; }
void SetPtFraction(Float_t pt) {fPtFraction = pt; }
void SetICMethod(Int_t i ) {fICMethod = i ; }
void SetParticleTypeInCone(Int_t i){fPartInCone = i;}
private:
Float_t fConeSize ; //Size of the isolation cone
- Float_t fPtThreshold ; //Mimium pt of the particles in the cone or sum in cone
+ Float_t fPtThreshold ; //Mimium pt of the particles in the cone or sum in cone (UE pt mean in the forward region cone)
+ Float_t fSumPtThreshold ; //Minium of sum pt of the particles in the cone (UE sum in the forward region cone)
Float_t fPtFraction ; //Fraction of the momentum of particles in cone or sum in cone
Int_t fICMethod ; //Isolation cut method to be used
// kPtIC: Pt threshold method
- ClassDef(AliIsolationCut,2)
+ ClassDef(AliIsolationCut,3)
} ;
Int_t trackIndex = 0;
Int_t nModule = -1;
+ //Get vertex for photon momentum calculation and event selection
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+ if (TMath::Abs(v[2]) > GetZvertexCut()) return ;
+
//Play with the MC stack if available
//Get the MC arrays and do some checks
if(IsDataMC()){
if(GetReader()->ReadStack()){
if(!GetMCStack()) {
- printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
- abort();
+ printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
+ abort();
}
//Fill some pure MC histograms, only primaries.
for(Int_t i=0 ; i<GetMCStack()->GetNprimary(); i++){//Only primary particles, for all MC transport put GetNtrack()
- TParticle *primary = GetMCStack()->Particle(i) ;
+ TParticle *primary = GetMCStack()->Particle(i) ;
//printf("i %d, %s: status = %d, primary? %d\n",i, primary->GetName(), primary->GetStatusCode(), primary->IsPrimary());
- if (primary->GetStatusCode() > 11) continue; //Working for PYTHIA and simple generators, check for HERWIG
- primary->Momentum(mom);
- MCHistograms(mom,TMath::Abs(primary->GetPdgCode()));
+ if (primary->GetStatusCode() > 11) continue; //Working for PYTHIA and simple generators, check for HERWIG
+ primary->Momentum(mom);
+ MCHistograms(mom,TMath::Abs(primary->GetPdgCode()));
} //primary loop
}
else if(GetReader()->ReadAODMCParticles()){
if(!GetReader()->GetAODMCParticles(0)) {
- printf("AliAnaPhoton::MakeAnalysisFillHistograms() - AODMCParticles not available!\n");
- abort();
+ printf("AliAnaPhoton::MakeAnalysisFillHistograms() - AODMCParticles not available!\n");
+ abort();
}
//Fill some pure MC histograms, only primaries.
for(Int_t i=0 ; i < (GetReader()->GetAODMCParticles(0))->GetEntriesFast(); i++){
- AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ;
+ AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ;
//printf("i %d, %s: primary? %d physical primary? %d, flag %d\n",
// i,(TDatabasePDG::Instance()->GetParticle(aodprimary->GetPdgCode()))->GetName(),
// aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), aodprimary->GetFlag());
- if (!aodprimary->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
+ if (!aodprimary->IsPrimary()) continue; //accept all which is not MC transport generated. Don't know how to avoid partons
//aodprimary->Momentum(mom);
- mom.SetPxPyPzE(aodprimary->Px(),aodprimary->Py(),aodprimary->Pz(),aodprimary->E());
- MCHistograms(mom,TMath::Abs(aodprimary->GetPdgCode()));
+ mom.SetPxPyPzE(aodprimary->Px(),aodprimary->Py(),aodprimary->Pz(),aodprimary->E());
+ MCHistograms(mom,TMath::Abs(aodprimary->GetPdgCode()));
} //primary loop
}
if(GetDebug() > 0)
printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In %s there are %d clusters \n", fCalorimeter.Data(), nCaloClusters);
- //Get vertex for photon momentum calculation
- Double_t v[3] = {0,0,0}; //vertex ;
- GetReader()->GetVertex(v);
AliVTrack * track = 0x0;
-
Float_t pos[3] ;
Float_t showerShape[3] ;
Double_t tof = 0;
Bool_t in = kTRUE;
if(IsFiducialCutOn()) in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
if(!in) continue;
-
+
//Get module of cluster
nCaloClustersAccepted++;
nModule = GetModuleNumber(clus);
//Cells per cluster
nCaloCellsPerCluster = clus->GetNCells();
//if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster);
-
+
//matched cluster with tracks
nTracksMatched = clus->GetNTracksMatched();
trackIndex = clus->GetTrackMatchedIndex();
if(GetDebug()>1) printf("Invariant mass \n");
//do not do for bad vertex
- Float_t fZvtxCut = 40. ;
- if(v[2]<-fZvtxCut || v[2]> fZvtxCut) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
+ // Float_t fZvtxCut = 40. ;
+ if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
Int_t nModule2 = -1;
Int_t nCaloCellsPerCluster2=0;
fhYE ->Fill(pos[1],e);
fhZE ->Fill(pos[2],e);
fhXYZ ->Fill(pos[0], pos[1],pos[2]);
-
+
fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
if(GetReader()->ReadStack() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){ //it MC stack and known tag
if( label >= GetMCStack()->GetNtrack()) {
- if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** large label ***: label %d, n tracks %d \n", label, GetMCStack()->GetNtrack());
- return ;
+ if(GetDebug() >= 0) printf("AliAnaCalorimeterQA::ClusterHistograms() *** large label ***: label %d, n tracks %d \n", label, GetMCStack()->GetNtrack());
+ return ;
}
primary = GetMCStack()->Particle(label);
iParent = primary->GetFirstMother();
if(GetDebug() > 1 ) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
- printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg0, primary->GetName(),status, iParent);
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
+ printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg0, primary->GetName(),status, iParent);
}
//Get final particle, no conversion products
if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
//Get the parent
- primary = GetMCStack()->Particle(iParent);
- pdg = TMath::Abs(primary->GetPdgCode());
- if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
- while((pdg == 22 || pdg == 11) && status != 1){
- iMother = iParent;
- primary = GetMCStack()->Particle(iMother);
- status = primary->GetStatusCode();
- iParent = primary->GetFirstMother();
- pdg = TMath::Abs(primary->GetPdgCode());
- if(GetDebug() > 1 )printf("\t pdg %d, index %d, %s, status %d \n",pdg, iMother, primary->GetName(),status);
- }
+ primary = GetMCStack()->Particle(iParent);
+ pdg = TMath::Abs(primary->GetPdgCode());
+ if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
+ while((pdg == 22 || pdg == 11) && status != 1){
+ iMother = iParent;
+ primary = GetMCStack()->Particle(iMother);
+ status = primary->GetStatusCode();
+ iParent = primary->GetFirstMother();
+ pdg = TMath::Abs(primary->GetPdgCode());
+ if(GetDebug() > 1 )printf("\t pdg %d, index %d, %s, status %d \n",pdg, iMother, primary->GetName(),status);
+ }
- if(GetDebug() > 1 ) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
- printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg, primary->GetName(), status, iParent);
- }
+ if(GetDebug() > 1 ) {
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
+ printf("\t Mother label %d, pdg %d, %s, status %d, parent %d \n",iMother, pdg, primary->GetName(), status, iParent);
+ }
}
//Overlapped pi0 (or eta, there will be very few), get the meson
if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
- if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: \n");
- while(pdg != 111 && pdg != 221){
- iMother = iParent;
- primary = GetMCStack()->Particle(iMother);
- status = primary->GetStatusCode();
- iParent = primary->GetFirstMother();
- pdg = TMath::Abs(primary->GetPdgCode());
- if(GetDebug() > 1 ) printf("\t pdg %d, %s, index %d\n",pdg, primary->GetName(),iMother);
- if(iMother==-1) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
+ if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: \n");
+ while(pdg != 111 && pdg != 221){
+ iMother = iParent;
+ primary = GetMCStack()->Particle(iMother);
+ status = primary->GetStatusCode();
+ iParent = primary->GetFirstMother();
+ pdg = TMath::Abs(primary->GetPdgCode());
+ if(GetDebug() > 1 ) printf("\t pdg %d, %s, index %d\n",pdg, primary->GetName(),iMother);
+ if(iMother==-1) {
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
//break;
- }
- }
+ }
+ }
- if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
+ if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
primary->GetName(),iMother);
}
else if(GetReader()->ReadAODMCParticles() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){//it MC AOD and known tag
//Get the list of MC particles
if(!GetReader()->GetAODMCParticles(0)) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - MCParticles not available!\n");
- abort();
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - MCParticles not available!\n");
+ abort();
}
aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(label);
iParent = aodprimary->GetMother();
if(GetDebug() > 1 ) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
- printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d \n",
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Cluster most contributing mother: \n");
+ printf("\t Mother label %d, pdg %d, Primary? %d, Physical Primary? %d, parent %d \n",
iMother, pdg0, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary(), iParent);
}
//Get final particle, no conversion products
if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCConversion)){
- if(GetDebug() > 1 )
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
+ if(GetDebug() > 1 )
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted cluster!. Find before conversion: \n");
//Get the parent
- aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iParent);
- pdg = TMath::Abs(aodprimary->GetPdgCode());
- while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary()) {
- iMother = iParent;
- aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
- status = aodprimary->IsPrimary();
- iParent = aodprimary->GetMother();
- pdg = TMath::Abs(aodprimary->GetPdgCode());
- if(GetDebug() > 1 )
- printf("\t pdg %d, index %d, Primary? %d, Physical Primary? %d \n",
+ aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iParent);
+ pdg = TMath::Abs(aodprimary->GetPdgCode());
+ while ((pdg == 22 || pdg == 11) && !aodprimary->IsPhysicalPrimary()) {
+ iMother = iParent;
+ aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
+ status = aodprimary->IsPrimary();
+ iParent = aodprimary->GetMother();
+ pdg = TMath::Abs(aodprimary->GetPdgCode());
+ if(GetDebug() > 1 )
+ printf("\t pdg %d, index %d, Primary? %d, Physical Primary? %d \n",
pdg, iMother, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
- }
-
- if(GetDebug() > 1 ) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
- printf("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d \n",
+ }
+
+ if(GetDebug() > 1 ) {
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Converted Cluster mother before conversion: \n");
+ printf("\t Mother label %d, pdg %d, parent %d, Primary? %d, Physical Primary? %d \n",
iMother, pdg, iParent, aodprimary->IsPrimary(), aodprimary->IsPhysicalPrimary());
- }
-
+ }
+
}
//Overlapped pi0 (or eta, there will be very few), get the meson
if(GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCPi0) ||
GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCEta)){
- if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: PDG %d, mom %d \n",pdg, iMother);
- while(pdg != 111 && pdg != 221){
-
- iMother = iParent;
- aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
- status = aodprimary->IsPrimary();
- iParent = aodprimary->GetMother();
- pdg = TMath::Abs(aodprimary->GetPdgCode());
+ if(GetDebug() > 1 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped Meson decay!, Find it: PDG %d, mom %d \n",pdg, iMother);
+ while(pdg != 111 && pdg != 221){
- if(GetDebug() > 1 ) printf("\t pdg %d, index %d\n",pdg, iMother);
-
- if(iMother==-1) {
- printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
+ iMother = iParent;
+ aodprimary = (AliAODMCParticle*)(GetReader()->GetAODMCParticles(0))->At(iMother);
+ status = aodprimary->IsPrimary();
+ iParent = aodprimary->GetMother();
+ pdg = TMath::Abs(aodprimary->GetPdgCode());
+
+ if(GetDebug() > 1 ) printf("\t pdg %d, index %d\n",pdg, iMother);
+
+ if(iMother==-1) {
+ printf("AliAnaCalorimeterQA::ClusterHistograms() - Tagged as Overlapped photon but meson not found, why?\n");
//break;
- }
- }
-
- if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
+ }
+ }
+
+ if(GetDebug() > 2 ) printf("AliAnaCalorimeterQA::ClusterHistograms() - Overlapped %s decay, label %d \n",
aodprimary->GetName(),iMother);
}
fhEMVxyz ->Fill(vxMC,vyMC);//,vz);
fhEMR ->Fill(e,rVMC);
if( nTracksMatched > 0){
- fhEleECharged ->Fill(e,eMC);
- fhElePtCharged ->Fill(pt,ptMC);
- fhElePhiCharged ->Fill(phi,phiMC);
- fhEleEtaCharged ->Fill(eta,etaMC);
+ fhEleECharged ->Fill(e,eMC);
+ fhElePtCharged ->Fill(pt,ptMC);
+ fhElePhiCharged ->Fill(phi,phiMC);
+ fhEleEtaCharged ->Fill(eta,etaMC);
}
}
else if(charge == 0){
fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
fhHaR ->Fill(e,rVMC);
if( nTracksMatched > 0){
- fhNeHadECharged ->Fill(e,eMC);
- fhNeHadPtCharged ->Fill(pt,ptMC);
- fhNeHadPhiCharged ->Fill(phi,phiMC);
- fhNeHadEtaCharged ->Fill(eta,etaMC);
+ fhNeHadECharged ->Fill(e,eMC);
+ fhNeHadPtCharged ->Fill(pt,ptMC);
+ fhNeHadPhiCharged ->Fill(phi,phiMC);
+ fhNeHadEtaCharged ->Fill(eta,etaMC);
}
}
else if(charge!=0){
fhHaVxyz ->Fill(vxMC,vyMC);//,vz);
fhHaR ->Fill(e,rVMC);
if( nTracksMatched > 0){
- fhChHadECharged ->Fill(e,eMC);
- fhChHadPtCharged ->Fill(pt,ptMC);
- fhChHadPhiCharged ->Fill(phi,phiMC);
- fhChHadEtaCharged ->Fill(eta,etaMC);
+ fhChHadECharged ->Fill(e,eMC);
+ fhChHadPtCharged ->Fill(pt,ptMC);
+ fhChHadPhiCharged ->Fill(phi,phiMC);
+ fhChHadEtaCharged ->Fill(eta,etaMC);
}
}
}//Work with MC
}
pLegendPhiCl.Draw();
- //ETA
+ //ETA
cetaphic->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhEtaPhiE->ProjectionX("heta_cluster_nocut",0,-1,0,-1);
htmp ->SetLineColor(1);
rbEta = GetNewRebinForRePlotting(htmp,etamin, etamax,netabins) ;
pLegendPhiCell.SetTextSize(0.03);
pLegendPhiCell.SetFillColor(10);
pLegendPhiCell.SetBorderSize(1);
-
+
+ delete htmp;
htmp = fhEtaPhiAmp->ProjectionY("hphi_cell_nocut",0,-1,0,-1);
if(htmp){
htmp->SetMinimum(1);
cetaphicell->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
-
+
+ delete htmp;
htmp = fhEtaPhiAmp->ProjectionX("heta_cell_nocut",0,-1,0,-1);
if(htmp){
htmp ->SetLineColor(1);
snprintf(name,buffersize,"QA_%s_ClusterX_Y_Z.eps",fCalorimeter.Data());
cx->Print(name); printf("Create plot %s\n",name);
}
- //CELLS
+ //CELLS
if(fFillAllPosHisto)
{
snprintf(cname,buffersize,"%s_QA_CellXY",fCalorimeter.Data());
pLegendXCl.SetFillColor(10);
pLegendXCl.SetBorderSize(1);
+ delete htmp;
htmp = fhRE->ProjectionX("hre_cluster_nocut",0,-1);
Int_t rbR=1;
if(htmp){
cxe->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhXE->ProjectionX("hxe_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
cxe->cd(3) ;
gPad->SetLogy();
gPad->SetGridy();
+
+ delete htmp;
htmp = fhYE->ProjectionX("hye_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhZE->ProjectionX("hze_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
pLegendXClN.SetFillColor(10);
pLegendXClN.SetBorderSize(1);
+ delete htmp;
htmp = fhRNCells->ProjectionX("hrn_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
cxn->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
+
+ delete htmp;
htmp = fhXNCells->ProjectionX("hxn_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
cxn->cd(3) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhYNCells->ProjectionX("hyn_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhZNCells->ProjectionX("hzn_cluster_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
pLegendXCell.SetFillColor(10);
pLegendXCell.SetBorderSize(1);
+ delete htmp;
htmp = fhRCellE->ProjectionX("hre_cell_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhXCellE->ProjectionX("hxe_cells_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
cxecell->cd(3) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhYCellE->ProjectionX("hye_cells_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhYCellE->ProjectionX(Form("hye_cells_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbY);
cxecell->cd(4) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhZCellE->ProjectionX("hze_cells_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhZCellE->ProjectionX(Form("hze_cells_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbZ);
pLegendXClD.SetFillColor(10);
pLegendXClD.SetBorderSize(1);
+ delete htmp;
htmp = fhDeltaCellClusterRE->ProjectionX("hrde_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhDeltaCellClusterRE->ProjectionX(Form("hrde_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbDR);
cxde->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhDeltaCellClusterXE->ProjectionX("hxde_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhDeltaCellClusterXE->ProjectionX(Form("hxde_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbDX);
cxde->cd(3) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhDeltaCellClusterYE->ProjectionX("hyde_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhDeltaCellClusterYE->ProjectionX(Form("hyde_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbDY);
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhDeltaCellClusterZE->ProjectionX("hzde_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
for (Int_t i = 0; i < ncuts; i++) {
binmin = hE->FindBin(ecut[i]);
//printf(" bins %d for e %f\n",binmin[i],ecut[i]);
+ delete htmp;
htmp = fhDeltaCellClusterZE->ProjectionX(Form("hzde_cut%d",i),binmin,-1);
htmp->SetLineColor(ecutcolor[i]);
htmp->Rebin(rbZ);
pLegendXClDN.SetTextSize(0.03);
pLegendXClDN.SetFillColor(10);
pLegendXClDN.SetBorderSize(1);
-
+ delete htmp;
htmp = fhDeltaCellClusterRNCells->ProjectionX("hrdn_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
pLegendXClDN.AddEntry(htmp,"No cut","L");
for (Int_t i = 0; i < ncellcuts; i++) {
+ delete htmp;
if(i < ncellcuts-1) htmp = fhDeltaCellClusterRNCells->ProjectionX(Form("hrdn_cut%d",i),ncellcut[i],ncellcut[i]);
else htmp = fhDeltaCellClusterRNCells->ProjectionX(Form("hrdn_cut%d",i),ncellcut[i],-1);
htmp->SetLineColor(ecutcolor[i]);
cxdn->cd(2) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhDeltaCellClusterXNCells->ProjectionX("hxdn_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
htmp->Draw("HE");
for (Int_t i = 0; i < ncellcuts; i++) {
+ delete htmp;
if(i < ncellcuts-1)htmp = fhDeltaCellClusterXNCells->ProjectionX(Form("hxdn_cut%d",i),ncellcut[i],ncellcut[i]);
else htmp = fhDeltaCellClusterXNCells->ProjectionX(Form("hxdn_cut%d",i),ncellcut[i],-1);
htmp->SetLineColor(ecutcolor[i]);
cxdn->cd(3) ;
gPad->SetLogy();
gPad->SetGridy();
+ delete htmp;
htmp = fhDeltaCellClusterYNCells->ProjectionX("hydn_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
htmp->Draw("HE");
for (Int_t i = 0; i < ncellcuts; i++) {
+ delete htmp;
if(i < ncellcuts-1) htmp = fhDeltaCellClusterYNCells->ProjectionX(Form("hydn_cut%d",i),ncellcut[i],ncellcut[i]);
else htmp = fhDeltaCellClusterYNCells->ProjectionX(Form("hydn_cut%d",i),ncellcut[i],-1);
htmp->SetLineColor(ecutcolor[i]);
cxdn->cd(4) ;
gPad->SetLogy();
gPad->SetGridy();
-
+ delete htmp;
htmp = fhDeltaCellClusterZNCells->ProjectionX("hzdn_nocut",0,-1);
if(htmp){
htmp->SetMinimum(1);
htmp->Draw("HE");
for (Int_t i = 0; i < ncellcuts; i++) {
+ delete htmp;
if(i < ncellcuts-1)htmp = fhDeltaCellClusterZNCells->ProjectionX(Form("hzdn_cut%d",i),ncellcut[i],ncellcut[i]);
else htmp = fhDeltaCellClusterZNCells->ProjectionX(Form("hzdn_cut%d",i),ncellcut[i],-1);
htmp->SetLineColor(ecutcolor[i]);
cN->cd(2) ;
gPad->SetLogx();
for(Int_t imod = 1; imod < fNModules; imod++){
+ delete htmp;
htmp = (TH1D*)fhNClustersMod[imod]->Clone(Form("hNClustersRat%d",imod));
htmp->Divide(fhNClustersMod[0]);
htmp->SetLineColor(modColorIndex[imod]);
cN->cd(4) ;
gPad->SetLogx();
for(Int_t imod = 1; imod < fNModules; imod++){
+ delete htmp;
htmp = (TH1D*)fhNCellsMod[imod]->Clone(Form("hNCellsRat%d",imod));
htmp->Divide(fhNCellsMod[0]);
htmp->SetLineColor(modColorIndex[imod]);
cN->cd(6) ;
gPad->SetLogx();
for(Int_t imod = 1; imod < fNModules; imod++){
+ delete htmp;
htmp = (TH1D*)hNCellsCluster1D[imod]->Clone(Form("hNClustersCells1DRat%d",imod));
htmp->Divide(hNCellsCluster1D[0]);
htmp->SetLineColor(modColorIndex[imod]);
// --- ROOT system ---
#include "TParticle.h"
#include "TH2F.h"
-
+#include "TH3D.h"
//---- AliRoot system ----
#include "AliAnaChargedParticles.h"
#include "AliCaloTrackReader.h"
//____________________________________________________________________________
AliAnaChargedParticles::AliAnaChargedParticles() :
- AliAnaPartCorrBaseClass(),fPdg(0), fhPt(0),fhPhi(0),fhEta(0),
+ AliAnaPartCorrBaseClass(),fPdg(0), fhNtracks(0),fhVertex(0), fhPt(0),fhPhi(0),fhEta(0),
+ fhPtEtaPhiPos(0), fhPtEtaPhiNeg(0),
fhPtPion(0),fhPhiPion(0),fhEtaPion(0),
fhPtProton(0),fhPhiProton(0),fhEtaProton(0),
fhPtElectron(0),fhPhiElectron(0),fhEtaElectron(0),
Float_t ptmin = GetHistoPtMin();
Float_t phimin = GetHistoPhiMin();
Float_t etamin = GetHistoEtaMin();
+
+ fhNtracks = new TH1F ("hNtracks","# of tracks", 1000,0,1000);
+ fhNtracks->SetXTitle("# of tracks");
+ outputContainer->Add(fhNtracks);
+
+ fhVertex = new TH3D ("Vertex","vertex position", 100,-50.,50., 100,-50.,50., 100,-50.,50.);
+ fhVertex->SetXTitle("X");
+ fhVertex->SetYTitle("Y");
+ fhVertex->SetZTitle("Z");
+ outputContainer->Add(fhVertex);
fhPt = new TH1F ("hPtCharged","p_T distribution", nptbins,ptmin,ptmax);
fhPt->SetXTitle("p_{T} (GeV/c)");
fhEta->SetXTitle("#eta ");
outputContainer->Add(fhEta);
+ fhPtEtaPhiPos = new TH3D ("hPtEtaPhiPos","pt/eta/phi of positive charge",nptbins,ptmin,ptmax, netabins,etamin,etamax, nphibins,phimin,phimax);
+ fhPtEtaPhiPos->SetXTitle("p_{T}^{h^{+}} (GeV/c)");
+ fhPtEtaPhiPos->SetYTitle("#eta ");
+ fhPtEtaPhiPos->SetZTitle("#phi (rad)");
+ outputContainer->Add(fhPtEtaPhiPos);
+
+ fhPtEtaPhiNeg = new TH3D ("hPtEtaPhiNeg","pt/eta/phi of negative charge",nptbins,ptmin,ptmax, netabins,etamin,etamax, nphibins,phimin,phimax);
+ fhPtEtaPhiNeg->SetXTitle("p_{T}^{h^{-}} (GeV/c)");
+ fhPtEtaPhiNeg->SetYTitle("#eta ");
+ fhPtEtaPhiNeg->SetZTitle("#phi (rad)");
+ outputContainer->Add(fhPtEtaPhiNeg);
if(IsDataMC()){
//Fill AODParticle with CTS aods
TVector3 p3;
+ Int_t evtIndex = 0;
+ Double_t vert[3]={0,0,0};
for(Int_t i = 0; i < ntracks; i++){
AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(i));
//Keep only particles identified with fPdg
//Selection not done for the moment
//Should be done here.
+ if (GetMixedEvent()){
+ evtIndex = GetMixedEvent()->EventIndex(track->GetID()) ;
+ }
+ GetVertex(vert,evtIndex);
+ if(TMath::Abs(vert[2])> GetZvertexCut()) continue;
AliAODPWG4Particle tr = AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
tr.SetDetector("CTS");
tr.SetLabel(track->GetLabel());
tr.SetTrackLabel(track->GetID(),-1);
+ // tr.SetChargedBit(track->Charge());
//Input from second AOD?
//if(GetReader()->GetAODCTSNormalInputEntries() <= i) tr.SetInputFileIndex(1);
//Loop on stored AODParticles
Int_t naod = GetOutputAODBranch()->GetEntriesFast();
+ if(naod!=0)fhNtracks->Fill(GetAODCTS()->GetEntriesFast()) ;
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+ fhVertex->Fill(v[0],v[1],v[2]);
if(GetDebug() > 0) printf("AliAnaChargedParticles::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
for(Int_t iaod = 0; iaod < naod ; iaod++){
AliAODPWG4Particle* tr = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
fhPt->Fill(tr->Pt());
fhPhi->Fill(tr->Pt(), tr->Phi());
fhEta->Fill(tr->Pt(), tr->Eta());
+ //for charge information
+ AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(iaod));
+ if(track->Charge()>0)fhPtEtaPhiPos->Fill(tr->Pt(), tr->Eta(),tr->Phi());
+ if(track->Charge()<0)fhPtEtaPhiNeg->Fill(tr->Pt(), tr->Eta(),tr->Phi());
if(IsDataMC()){
//Play with the MC stack if available
- Int_t mompdg = -1;
- Int_t label = tr->GetLabel();
- if(GetReader()->ReadStack()){
- TParticle * mom = GetMCStack()->Particle(label);
- mompdg =TMath::Abs(mom->GetPdgCode());
+ Int_t mompdg = -1;
+ Int_t label = tr->GetLabel();
+ if(GetReader()->ReadStack()){
+ TParticle * mom = GetMCStack()->Particle(label);
+ mompdg =TMath::Abs(mom->GetPdgCode());
}
- else if(GetReader()->ReadAODMCParticles()){
- AliAODMCParticle * aodmom = 0;
- //Get the list of MC particles
- aodmom = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(tr->GetInputFileIndex()))->At(label);
- mompdg =TMath::Abs(aodmom->GetPdgCode());
- }
-
- if(mompdg==211){
- fhPtPion->Fill(tr->Pt());
- fhPhiPion->Fill(tr->Pt(), tr->Phi());
- fhEtaPion->Fill(tr->Pt(), tr->Eta());
+ else if(GetReader()->ReadAODMCParticles()){
+ AliAODMCParticle * aodmom = 0;
+ //Get the list of MC particles
+ aodmom = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(tr->GetInputFileIndex()))->At(label);
+ mompdg =TMath::Abs(aodmom->GetPdgCode());
+ }
+
+ if(mompdg==211){
+ fhPtPion->Fill(tr->Pt());
+ fhPhiPion->Fill(tr->Pt(), tr->Phi());
+ fhEtaPion->Fill(tr->Pt(), tr->Eta());
}
else if(mompdg==2212){
- fhPtProton->Fill(tr->Pt());
- fhPhiProton->Fill(tr->Pt(), tr->Phi());
- fhEtaProton->Fill(tr->Pt(), tr->Eta());
+ fhPtProton->Fill(tr->Pt());
+ fhPhiProton->Fill(tr->Pt(), tr->Phi());
+ fhEtaProton->Fill(tr->Pt(), tr->Eta());
}
else if(mompdg==321){
- fhPtKaon->Fill(tr->Pt());
- fhPhiKaon->Fill(tr->Pt(), tr->Phi());
- fhEtaKaon->Fill(tr->Pt(), tr->Eta());
+ fhPtKaon->Fill(tr->Pt());
+ fhPhiKaon->Fill(tr->Pt(), tr->Phi());
+ fhEtaKaon->Fill(tr->Pt(), tr->Eta());
}
else if(mompdg==11){
- fhPtElectron->Fill(tr->Pt());
- fhPhiElectron->Fill(tr->Pt(), tr->Phi());
- fhEtaElectron->Fill(tr->Pt(), tr->Eta());
+ fhPtElectron->Fill(tr->Pt());
+ fhPhiElectron->Fill(tr->Pt(), tr->Phi());
+ fhEtaElectron->Fill(tr->Pt(), tr->Eta());
}
else {
- //printf("unknown pdg %d\n",mompdg);
- fhPtUnknown->Fill(tr->Pt());
- fhPhiUnknown->Fill(tr->Pt(), tr->Phi());
- fhEtaUnknown->Fill(tr->Pt(), tr->Eta());
+ //printf("unknown pdg %d\n",mompdg);
+ fhPtUnknown->Fill(tr->Pt());
+ fhPhiUnknown->Fill(tr->Pt(), tr->Phi());
+ fhEtaUnknown->Fill(tr->Pt(), tr->Eta());
}
}//Work with stack also
}// aod branch loop
//-- Author: Gustavo Conesa (INFN-LNF)
// Root system
class TH2F;
-
+class TH3D;
// Analysis system
#include "AliAnaPartCorrBaseClass.h"
Int_t fPdg ; //identified particle id
//Histograms
+ TH1F * fhNtracks; //! track multiplicity distribution
+ TH3D * fhVertex; //! vertex distribution
TH1F * fhPt; //! pT distribution
TH2F * fhPhi; //! phi distribution vs pT
TH2F * fhEta; //! eta distribution vs pT
-
+ TH3D * fhPtEtaPhiPos; //! pT and phi distribution of positive charge
+ TH3D * fhPtEtaPhiNeg; //! pT and phi distribution of positive charge
+
//MC
TH1F * fhPtPion; //! pT distribution
TH2F * fhPhiPion; //! phi distribution vs pT
TH2F * fhPhiUnknown; //! phi distribution vs pT
TH2F * fhEtaUnknown; //! eta distribution vs pT
- ClassDef(AliAnaChargedParticles,1)
+ ClassDef(AliAnaChargedParticles,2)
} ;
// 2. change the correlation variable
// 3. Only use leading particle(cluster/track) as trigger for correlation (2010/07/02)
// 4. Make decay photon-hadron correlations where decay contribute pi0 mass (2010/09/09)
+// 5. fill the pout to extract kt at the end, also to study charge asymmetry(2010/10/06)
//////////////////////////////////////////////////////////////////////////////
AliAnaPartCorrBaseClass(),
fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fSelectIsolated(0),
fMakeSeveralUE(0), fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(0.),
- fPi0AODBranchName(""),fPi0Trigger(0),
+ fPi0AODBranchName(""),fNeutralCorr(0), fPi0Trigger(0),
+ // fMultiBin(0),fNZvertBin(0),fNrpBin(0),fZvtxCut(0.),
+ // fUseSelectEvent(kFALSE),
+ fhNclustersNtracks(0), fhVertex(0),
fhPtLeading(0),fhPhiLeading(0),fhEtaLeading(0),
- fhDeltaPhiDeltaEtaCharged(0),fhDeltaPhiDeltaEtaNeutral(0),
- fhPhiCharged(0), fhPhiNeutral(0), fhEtaCharged(0), fhEtaNeutral(0),
- fhDeltaPhiCharged(0), fhDeltaPhiNeutral(0),
- fhDeltaEtaCharged(0), fhDeltaEtaNeutral(0),
- fhDeltaPhiChargedPt(0), fhDeltaPhiNeutralPt(0),
- fhDeltaPhiUeChargedPt(0), fhDeltaPhiUeNeutralPt(0),
- fhPtImbalanceNeutral(0), fhPtImbalanceCharged(0),
- fhPtImbalanceUeCharged(0),fhPtImbalanceUeNeutral(0),
+ fhDeltaPhiDeltaEtaCharged(0),
+ fhPhiCharged(0), fhEtaCharged(0),
+ fhDeltaPhiCharged(0),
+ fhDeltaEtaCharged(0),
+ fhDeltaPhiChargedPt(0),
+ fhDeltaPhiUeChargedPt(0),
+ fhPtImbalanceCharged(0),
+ fhPtImbalanceUeCharged(0),
+ fhPtImbalancePosCharged(0),fhPtImbalanceNegCharged(0),
fhPtHbpCharged(0), fhPtHbpUeCharged(0),
- fhPtHbpNeutral(0), fhPtHbpUeNeutral(0),
fhDeltaPhiUeLeftCharged(0),fhDeltaPhiUeRightCharged(0),
- fhDeltaPhiUeLeftNeutral(0),fhDeltaPhiUeRightNeutral(0),
fhPtImbalanceUeLeftCharged(0),fhPtImbalanceUeRightCharged(0),
+ fhPtHbpUeLeftCharged(0),fhPtHbpUeRightCharged(0),
+ fhPoutTrig(0), fhPtTrigCharged(0),
+ fhTrigDeltaPhiDeltaEtaCharged(0x0), fhTrigCorr(0x0),fhTrigUeCorr(0x0),
+ fhDeltaPhiDeltaEtaNeutral(0),
+ fhPhiNeutral(0), fhEtaNeutral(0),
+ fhDeltaPhiNeutral(0), fhDeltaEtaNeutral(0),
+ fhDeltaPhiNeutralPt(0),fhDeltaPhiUeNeutralPt(0),
+ fhPtImbalanceNeutral(0),fhPtImbalanceUeNeutral(0),
+ fhPtHbpNeutral(0), fhPtHbpUeNeutral(0),
+ fhDeltaPhiUeLeftNeutral(0),fhDeltaPhiUeRightNeutral(0),
fhPtImbalanceUeLeftNeutral(0),fhPtImbalanceUeRightNeutral(0),
- fhPtHbpUeLeftCharged(0),fhPtHbpUeRightCharged(0),
fhPtHbpUeLeftNeutral(0),fhPtHbpUeRightNeutral(0),
fhPtPi0DecayRatio(0),
- fhDeltaPhiDecay1Charged(0),fhDeltaPhiDecay2Charged(0),
- fhPtImbalanceDecay1Charged(0), fhPtImbalanceDecay2Charged(0),
- fhDeltaPhiDecay1Neutral(0),fhDeltaPhiDecay2Neutral(0),
- fhPtImbalanceDecay1Neutral(0), fhPtImbalanceDecay2Neutral(0)
+ fhDeltaPhiDecayCharged(0),
+ fhPtImbalanceDecayCharged(0),
+ fhDeltaPhiDecayNeutral(0),
+ fhPtImbalanceDecayNeutral(0)
{
//Default Ctor
Float_t ptmin = GetHistoPtMin();
Float_t phimin = GetHistoPhiMin();
Float_t etamin = GetHistoEtaMin();
+
+
+ fhNclustersNtracks = new TH2F ("Multiplicity","Neutral cluster and charged track multiplicity",1000, 0., 1000.,1000, 0., 1000.);
+ fhNclustersNtracks->SetYTitle("# of tracks");
+ fhNclustersNtracks->SetXTitle("# of clusters");
+
+ fhVertex = new TH3D ("Vertex","vertex position", 100,-50.,50., 100,-50.,50., 100,-50.,50.);
+ fhVertex->SetXTitle("X");
+ fhVertex->SetYTitle("Y");
+ fhVertex->SetZTitle("Z");
fhPtLeading = new TH1F ("hPtLeading","p_T distribution of leading particles", nptbins,ptmin,ptmax);
fhPtLeading->SetXTitle("p_{T}^{trig} (GeV/c)");
fhEtaLeading = new TH2F ("hEtaLeading","#eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax);
fhEtaLeading->SetYTitle("#eta ");
+ outputContainer->Add(fhNclustersNtracks);
+ outputContainer->Add(fhVertex);
outputContainer->Add(fhPtLeading);
outputContainer->Add(fhPhiLeading);
outputContainer->Add(fhEtaLeading);
fhPtImbalanceUeCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
fhPtImbalanceUeCharged->SetXTitle("p_{T trigger}");
+ fhPtImbalancePosCharged =
+ new TH2F("CorrelationPositiveCharged","z_{trigger h^{+}} = p_{T h^{+}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
+ fhPtImbalancePosCharged->SetYTitle("z_{trigger h^{+}}");
+ fhPtImbalancePosCharged->SetXTitle("p_{T trigger}");
+
+ fhPtImbalanceNegCharged =
+ new TH2F("CorrelationNegativeCharged","z_{trigger h^{-}} = p_{T h^{-}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
+ fhPtImbalanceNegCharged->SetYTitle("z_{trigger h^{-}}");
+ fhPtImbalanceNegCharged->SetXTitle("p_{T trigger}");
+
fhPtHbpCharged =
new TH2F("HbpCharged","#xi = ln(1/x_{E}) with charged hadrons",
nptbins,ptmin,ptmax,200,0.,10.);
nptbins,ptmin,ptmax,200,0.,10.);
fhPtHbpUeCharged->SetYTitle("ln(1/x_{E})");
fhPtHbpUeCharged->SetXTitle("p_{T trigger}");
+
+ fhPoutTrig =
+ new TH2F("PoutTrigPt","Pout with charged hadrons",
+ nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
+ fhPoutTrig->SetYTitle("p_{out} (GeV/c)");
+ fhPoutTrig->SetXTitle("p_{T trigger} (GeV/c)");
+
+ fhPtTrigCharged =
+ new TH2F("PtTrigCharged","trgger and charged tracks pt distribution",
+ nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
+ fhPtTrigCharged->SetYTitle("p_{T asso} (GeV/c)");
+ fhPtTrigCharged->SetXTitle("p_{T trigger} (GeV/c)");
outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
outputContainer->Add(fhPhiCharged) ;
outputContainer->Add(fhDeltaPhiChargedPt) ;
outputContainer->Add(fhDeltaPhiUeChargedPt) ;
outputContainer->Add(fhPtImbalanceCharged) ;
+ outputContainer->Add(fhPtImbalancePosCharged) ;
+ outputContainer->Add(fhPtImbalanceNegCharged) ;
outputContainer->Add(fhPtImbalanceUeCharged) ;
outputContainer->Add(fhPtHbpCharged) ;
outputContainer->Add(fhPtHbpUeCharged) ;
+ outputContainer->Add(fhPoutTrig) ;
+ outputContainer->Add(fhPtTrigCharged) ;
+ if(DoEventSelect()){
+ Int_t nMultiBins = GetMultiBin();
+ fhTrigDeltaPhiDeltaEtaCharged = new TH3D*[nMultiBins] ;
+ fhTrigCorr = new TH2F*[nMultiBins];
+ fhTrigUeCorr = new TH2F*[nMultiBins];
+ for(Int_t im=0; im<nMultiBins; im++){
+ fhTrigDeltaPhiDeltaEtaCharged[im] = new TH3D
+ (Form("fhTrigDeltaPhiDeltaEtaCharged_%d",im),Form("fhTrigDeltaPhiDeltaEtaCharged_%d",im), nptbins,ptmin,ptmax, 140,-2.,5., 200,-2,2);
+ fhTrigDeltaPhiDeltaEtaCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
+ fhTrigDeltaPhiDeltaEtaCharged[im]->SetYTitle("#Delta #phi");
+ fhTrigDeltaPhiDeltaEtaCharged[im]->SetZTitle("#Delta #eta");
+ fhTrigCorr[im] = new TH2F
+ (Form("fhTrigPtCorr_%d",im),Form("fhTrigPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
+ fhTrigCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
+ fhTrigCorr[im]->SetXTitle("p_{T trigger}");
+ fhTrigUeCorr[im] = new TH2F
+ (Form("fhTrigPtUeCorr_%d",im),Form("fhTrigPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
+ fhTrigUeCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
+ fhTrigUeCorr[im]->SetXTitle("p_{T trigger}");
+
+ outputContainer->Add(fhTrigDeltaPhiDeltaEtaCharged[im]) ;
+ outputContainer->Add(fhTrigCorr[im]);
+ outputContainer->Add(fhTrigUeCorr[im]);
+
+ }
+ }
+
if(fPi0Trigger){
fhPtPi0DecayRatio = new TH3D
("hPtPi0DecayRatio","p_T of #pi^{0} and the ratio of pt for two decay",
fhPtPi0DecayRatio->SetYTitle("p_{T}^{Decay1}/p_{T}^{#pi^{0}}");
fhPtPi0DecayRatio->SetZTitle("p_{T}^{Decay2}/p_{T}^{#pi^{0}}");
- fhDeltaPhiDecay1Charged = new TH2F
- ("DeltaPhiDecay1Charged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
- nptbins,ptmin,ptmax,140,-2.,5.);
- fhDeltaPhiDecay1Charged->SetYTitle("#Delta #phi");
- fhDeltaPhiDecay1Charged->SetXTitle("p_{T Decay} (GeV/c)");
-
- fhDeltaPhiDecay2Charged = new TH2F
- ("DeltaPhiDecay2Charged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
+ fhDeltaPhiDecayCharged = new TH2F
+ ("DeltaPhiDecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
nptbins,ptmin,ptmax,140,-2.,5.);
- fhDeltaPhiDecay2Charged->SetYTitle("#Delta #phi");
- fhDeltaPhiDecay2Charged->SetXTitle("p_{T Decay} (GeV/c)");
+ fhDeltaPhiDecayCharged->SetYTitle("#Delta #phi");
+ fhDeltaPhiDecayCharged->SetXTitle("p_{T Decay} (GeV/c)");
- fhPtImbalanceDecay1Charged =
- new TH2F("CorrelationDecay1Charged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
+ fhPtImbalanceDecayCharged =
+ new TH2F("CorrelationDecayCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
nptbins,ptmin,ptmax,200,0.,2.);
- fhPtImbalanceDecay1Charged->SetYTitle("z_{decay h^{#pm}}");
- fhPtImbalanceDecay1Charged->SetXTitle("p_{T decay}");
+ fhPtImbalanceDecayCharged->SetYTitle("z_{decay h^{#pm}}");
+ fhPtImbalanceDecayCharged->SetXTitle("p_{T decay}");
- fhPtImbalanceDecay2Charged =
- new TH2F("CorrelationDecay2Charged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
- nptbins,ptmin,ptmax,200,0.,2.);
- fhPtImbalanceDecay2Charged->SetYTitle("z_{decay h^{#pm}}");
- fhPtImbalanceDecay2Charged->SetXTitle("p_{T decay}");
outputContainer->Add(fhPtPi0DecayRatio) ;
- outputContainer->Add(fhDeltaPhiDecay1Charged) ;
- outputContainer->Add(fhDeltaPhiDecay2Charged) ;
- outputContainer->Add(fhPtImbalanceDecay1Charged) ;
- outputContainer->Add(fhPtImbalanceDecay2Charged) ;
+ outputContainer->Add(fhDeltaPhiDecayCharged) ;
+ outputContainer->Add(fhPtImbalanceDecayCharged) ;
}
} //Correlation with charged hadrons
//Correlation with neutral hadrons
- if(GetReader()->IsEMCALSwitchedOn() || GetReader()->IsPHOSSwitchedOn()){
+ if(fNeutralCorr){
fhDeltaPhiDeltaEtaNeutral = new TH2F
("DeltaPhiDeltaEtaNeutral","#phi_{trigger} - #phi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
outputContainer->Add(fhPtHbpUeNeutral) ;
if(fPi0Trigger){
- fhDeltaPhiDecay1Neutral = new TH2F
- ("DeltaPhiDecay1Neutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
- nptbins,ptmin,ptmax,140,-2.,5.);
- fhDeltaPhiDecay1Neutral->SetYTitle("#Delta #phi");
- fhDeltaPhiDecay1Neutral->SetXTitle("p_{T Decay} (GeV/c)");
-
- fhDeltaPhiDecay2Neutral = new TH2F
- ("DeltaPhiDecay2Neutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
+ fhDeltaPhiDecayNeutral = new TH2F
+ ("DeltaPhiDecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
nptbins,ptmin,ptmax,140,-2.,5.);
- fhDeltaPhiDecay2Neutral->SetYTitle("#Delta #phi");
- fhDeltaPhiDecay2Neutral->SetXTitle("p_{T Decay} (GeV/c)");
-
- fhPtImbalanceDecay1Neutral =
- new TH2F("CorrelationDecay1Neutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
- nptbins,ptmin,ptmax,200,0.,2.);
- fhPtImbalanceDecay1Neutral->SetYTitle("z_{decay h^{0}}");
- fhPtImbalanceDecay1Neutral->SetXTitle("p_{T decay}");
+ fhDeltaPhiDecayNeutral->SetYTitle("#Delta #phi");
+ fhDeltaPhiDecayNeutral->SetXTitle("p_{T Decay} (GeV/c)");
- fhPtImbalanceDecay2Neutral =
- new TH2F("CorrelationDecay2Neutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
+ fhPtImbalanceDecayNeutral =
+ new TH2F("CorrelationDecayNeutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
nptbins,ptmin,ptmax,200,0.,2.);
- fhPtImbalanceDecay2Neutral->SetYTitle("z_{decay h^{0}}");
- fhPtImbalanceDecay2Neutral->SetXTitle("p_{T decay}");
+ fhPtImbalanceDecayNeutral->SetYTitle("z_{decay h^{0}}");
+ fhPtImbalanceDecayNeutral->SetXTitle("p_{T decay}");
- outputContainer->Add(fhDeltaPhiDecay1Neutral) ;
- outputContainer->Add(fhDeltaPhiDecay2Neutral) ;
- outputContainer->Add(fhPtImbalanceDecay1Neutral) ;
- outputContainer->Add(fhPtImbalanceDecay2Neutral) ;
+ outputContainer->Add(fhDeltaPhiDecayNeutral) ;
+ outputContainer->Add(fhPtImbalanceDecayNeutral) ;
}
fMakeSeveralUE = kFALSE;
fUeDeltaPhiMinCut = 1. ;
fUeDeltaPhiMaxCut = 1.5 ;
+ fNeutralCorr = kFALSE ;
fPi0Trigger = kFALSE ;
}
MakeChargedCorrelation(particle, GetAODCTS(),kFALSE);
TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
- if(pi0list && pi0list->GetEntriesFast() > 0)
+ if(fNeutralCorr && pi0list && pi0list->GetEntriesFast() > 0)
MakeNeutralCorrelation(particle, pi0list,kFALSE);
}//Correlate leading
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
}
+ //Get the vertex and check it is not too large in z
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+ if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;
+
//Loop on stored AOD particles, find leading
Int_t naod = GetInputAODBranch()->GetEntriesFast();
+ if(naod!=0)fhNclustersNtracks->Fill(naod, GetAODCTS()->GetEntriesFast());
+ if(naod!=0)fhVertex->Fill(v[0],v[1],v[2]);
Double_t ptTrig = 0.;
Int_t trigIndex = -1 ;
for(Int_t iaod = 0; iaod < naod ; iaod++){ //loop on input trigger AOD file
AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
-
+ //vertex cut in case of mixing
+ if (GetMixedEvent()) {
+ Int_t evt=-1;
+ Int_t id =-1;
+ if (particle->GetCaloLabel(0) >= 0 ){
+ id=particle->GetCaloLabel(0);
+ if(id >= 0 )evt= GetMixedEvent()-> EventIndexForCaloCluster(id) ;
+ }
+ else if(particle->GetTrackLabel(0) >= 0 ){
+ id=particle->GetTrackLabel(0);
+ if(id >= 0 )evt= GetMixedEvent()->EventIndex(id) ;
+ }
+ else continue;
+
+ if (TMath::Abs(GetVertex(evt)[2]) > GetZvertexCut())
+ continue ;
+ }
+
//check if the particle is isolated or if we want to take the isolation into account
if(OnlyIsolated() && !particle->IsIsolated()) continue;
- //find the leading particles with highest momentum
+ //check if inside the vertex cut
+ //find the leading particles with highest momentum
if (particle->Pt()>ptTrig) {
ptTrig = particle->Pt() ;
trigIndex = iaod ;
if(trigIndex!=-1){ //using trigger partilce to do correlations
AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(trigIndex));
+ if (GetMixedEvent()) {
+ Int_t evt=-1;
+ Int_t id = 0;
+ if (particle->GetCaloLabel(0) >= 0 ){
+ id=particle->GetCaloLabel(0);
+ if(id >= 0 )evt= GetMixedEvent()-> EventIndexForCaloCluster(id) ;
+ }
+ else if(particle->GetTrackLabel(0) >= 0 ){
+ id=particle->GetTrackLabel(0);
+ if(id >= 0 )evt= GetMixedEvent()->EventIndex(id) ;
+ }
+ else return;
+
+ if (TMath::Abs(GetVertex(evt)[2]) > GetZvertexCut())
+ return ;
+ }
+
//Fill leading particle histogram
fhPtLeading->Fill(particle->Pt());
Float_t phi = particle->Phi();
//Make correlation with neutral pions
TObjArray * refpi0 = particle->GetObjArray(GetAODObjArrayName()+"Pi0s");
- if(refpi0){
+ if(refpi0 && fNeutralCorr){
if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Leading %d, In Pi0 Refs entries %d\n",trigIndex, refpi0->GetEntriesFast());
MakeNeutralCorrelation(particle, refpi0, kTRUE);
}
Int_t evtIndex12 = 0 ;
Int_t indexPhoton1 = -1 ;
Int_t indexPhoton2 = -1 ;
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+ if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;
+
+ Int_t nTracks = GetAODCTS()->GetEntriesFast() ;
if (GetMixedEvent()) {
evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
Double_t pyTrig = aodParticle->Py();
Double_t phiTrig = aodParticle->Phi();
+ Double_t etaTrig = aodParticle->Eta();
Double_t pt = -100.;
Double_t px = -100.;
//Int_t currentIndex = -1 ;
for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ ){
AliAODTrack * track = (AliAODTrack *) (pl->At(ipr)) ;
+ //check if inside the vertex cut
+ //printf("charge = %d\n", track->Charge());
Int_t evtIndex2 = 0 ;
if (GetMixedEvent()) {
evtIndex2 = GetMixedEvent()->EventIndex(track->GetID()) ;
if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2) // photon and track from different events
continue ;
+ //vertex cut
+ if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut())
+ continue ;
// if(currentIndex == evtIndex2) // tracks from different event
// continue ;
// currentIndex = evtIndex2 ;
if(track->GetID()==aodParticle->GetTrackLabel(0))
continue ;
- //if(track->Pt()>ptTrig)
- //continue ;
-
+ if(track->Pt()==ptTrig && track->Phi()==phiTrig && track->Eta()==etaTrig)
+ continue ;
+
Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
p3.SetXYZ(mom[0],mom[1],mom[2]);
pt = p3.Pt();
fhDeltaPhiDeltaEtaCharged->Fill(deltaphi,aodParticle->Eta()-eta);
if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
-
+ //fill different multiplicity histogram
+ if(DoEventSelect()){
+ for(Int_t im=0; im<GetMultiBin(); im++){
+ if(nTracks < GetMaxMulti()/GetMultiBin()*(im+1))fhTrigDeltaPhiDeltaEtaCharged[im]->Fill(ptTrig,deltaphi,aodParticle->Eta()-eta);
+ }
+ }
//delta phi cut for correlation
if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
fhDeltaPhiChargedPt->Fill(pt,deltaphi);
- fhPtImbalanceCharged->Fill(ptTrig,rat);
+ fhPtImbalanceCharged->Fill(ptTrig,xE);
fhPtHbpCharged->Fill(ptTrig,cosi);
- }
+ Double_t pout = pt*TMath::Sin(deltaphi) ;
+ fhPoutTrig->Fill(ptTrig, pout) ;
+ fhPtTrigCharged->Fill(ptTrig, pt) ;
+ if(track->Charge()>0) fhPtImbalancePosCharged->Fill(ptTrig,xE) ;
+ else fhPtImbalanceNegCharged->Fill(ptTrig,xE) ;
+ //fill different multiplicity histogram
+ if(DoEventSelect()){
+ for(Int_t im=0; im<GetMultiBin(); im++){
+ if(nTracks < ( GetMaxMulti() - 0 )/GetMultiBin()*(im+1))
+ fhTrigCorr[im]->Fill(ptTrig,xE);
+ }
+ } //multiplicity events selection
+ } //delta phi cut for correlation
else {
fhDeltaPhiUeChargedPt->Fill(pt,deltaphi);
- fhPtImbalanceUeCharged->Fill(ptTrig,rat);
+ fhPtImbalanceUeCharged->Fill(ptTrig,xE);
fhPtHbpUeCharged->Fill(ptTrig,cosi);
+ if(DoEventSelect()){
+ for(Int_t im=0; im<GetMultiBin(); im++){
+ if(nTracks < ( GetMaxMulti() - 0 )/GetMultiBin()*(im+1))
+ fhTrigUeCorr[im]->Fill(ptTrig,xE);
+ }
+ } //multiplicity events selection
+
}
if(fPi0Trigger){
if(indexPhoton1!=-1 && indexPhoton2!=-1){
- fhDeltaPhiDecay1Charged->Fill(ptDecay1, deltaphiDecay1);
- fhDeltaPhiDecay2Charged->Fill(ptDecay2, deltaphiDecay2);
+ fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaphiDecay1);
+ fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaphiDecay2);
if(GetDebug() > 1)printf("deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaphiDecay1, deltaphiDecay2);
if( (deltaphiDecay1 > fDeltaPhiMinCut) && ( deltaphiDecay1 < fDeltaPhiMaxCut) )
- fhPtImbalanceDecay1Charged->Fill(ptDecay1,ratDecay1);
+ fhPtImbalanceDecayCharged->Fill(ptDecay1,ratDecay1);
if( (deltaphiDecay2 > fDeltaPhiMinCut) && ( deltaphiDecay2 < fDeltaPhiMaxCut) )
- fhPtImbalanceDecay2Charged->Fill(ptDecay2,ratDecay2);
+ fhPtImbalanceDecayCharged->Fill(ptDecay2,ratDecay2);
if(GetDebug() > 1)printf("ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
} //index of decay photons found
} //make decay-hadron correlation
if(deltaphiDecay1>3*TMath::PiOver2()) deltaphiDecay1-=TMath::TwoPi();
if(deltaphiDecay2< -TMath::PiOver2()) deltaphiDecay2+=TMath::TwoPi();
if(deltaphiDecay2>3*TMath::PiOver2()) deltaphiDecay2-=TMath::TwoPi();
- fhDeltaPhiDecay1Neutral->Fill(ptDecay1, deltaphiDecay1);
- fhDeltaPhiDecay2Neutral->Fill(ptDecay2, deltaphiDecay2);
+ fhDeltaPhiDecayNeutral->Fill(ptDecay1, deltaphiDecay1);
+ fhDeltaPhiDecayNeutral->Fill(ptDecay2, deltaphiDecay2);
if(GetDebug() > 1)printf("deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaphiDecay1, deltaphiDecay2);
if( (deltaphiDecay1 > fDeltaPhiMinCut) && ( deltaphiDecay1 < fDeltaPhiMaxCut) )
- fhPtImbalanceDecay1Neutral->Fill(ptDecay1,ratDecay1);
+ fhPtImbalanceDecayNeutral->Fill(ptDecay1,ratDecay1);
if( (deltaphiDecay2 > fDeltaPhiMinCut) && ( deltaphiDecay2 < fDeltaPhiMaxCut) )
- fhPtImbalanceDecay2Neutral->Fill(ptDecay2,ratDecay2);
+ fhPtImbalanceDecayNeutral->Fill(ptDecay2,ratDecay2);
if(GetDebug() > 1)printf("ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
}
} //do decay-hadron correlation
// 2. change the correlation variable
// 3. Only use leading particle(cluster/track) as trigger for correlation (2010/07/02)
// 4. Make decay photon-hadron correlations where decay contribute pi0 mass (2010/09/09)
-
+// 5. fill the pout to extract kt at the end, also to study charge asymmetry(2010/10/06)
+// 6. Add the possibality for event selection analysis based on vertex and multiplicity bins (10/10/2010)
// --- ROOT system ---
class TH3D;
void SwitchOnSeveralUECalculation() { fMakeSeveralUE = kTRUE;}
void SwitchOffSeveralUECalculation() { fMakeSeveralUE = kFALSE;}
+ // Do trigger-neutral correlation
+ Bool_t DoNeutralCorr() const {return fNeutralCorr ; }
+ void SwitchOnNeutralCorr() { fNeutralCorr = kTRUE;}
+ void SwitchOffNeutralCorr() { fNeutralCorr = kFALSE;}
+
// Do decay-hadron correlation if it is pi0 trigger
Bool_t IsPi0Trigger() const {return fPi0Trigger ; }
void SwitchOnDecayCorr() { fPi0Trigger = kTRUE;}
Bool_t OnlyIsolated() const {return fSelectIsolated ; }
void SelectIsolated(Bool_t select) {fSelectIsolated = select ; }
+
+// //Setters for parameters of event buffers
+// void SetMultiBin(Int_t n=1) {fMultiBin=n ;} //number of bins in Multiplicity
+// void SetNRPBin(Int_t n=1) {fNrpBin=n ;} //number of bins in reaction plain
+// //Setters for event selection
+// void SetZvertexCut(Float_t zcut=40.){fZvtxCut=zcut ;} //cut on vertex position
+// Int_t GetMultiBin() const {return fMultiBin ;} //number of bins in Multiplicity
+// Int_t GetNRPBin() const {return fNrpBin=n ;} //number of bins in reaction plain
+// //Getters for event selection
+// Float_t GetZvertexCut() const {return fZvtxCut ;} //cut on vertex position
+// void SwitchOnEventSelection() {fUseSelectEvent = kTRUE ; }
+// void SwitchOffEventSelection() {fUseSelectEvent = kFALSE ; } s
+// // Do correlation analysis with different event buffers
+// Bool_t IsEventSelect() const {return fUseSelectEvent ; }
void InitParameters();
Double_t fUeDeltaPhiMaxCut ; // Minimum Delta Phi Gamma-Underlying Hadron
Double_t fUeDeltaPhiMinCut ; // Maximum Delta Phi Gamma-Underlying Hadron
TString fPi0AODBranchName; // Name of AOD branch with pi0, not trigger
- Bool_t fPi0Trigger ; // Do analysis with decay photon from pi0 trigger
-
+ Bool_t fNeutralCorr ; // switch the analysis with neutral particles
+ Bool_t fPi0Trigger ; // switch the analysis with decay photon from pi0 trigger
+// Int_t fMultiBin ; // Number of bins in event container for multiplicity
+// Int_t fNZvertBin ; // Number of bins in event container for vertex position
+// Int_t fNrpBin ; // Number of bins in event container for reaction plain
+// Float_t fZvtxCut ; // Cut on vertex position
+// Bool_t fUseSelectEvent ; // Select events based on multiplicity and vertex cuts
+
+
//Histograms
+ TH2F * fhNclustersNtracks; //charge and cluster multiplicity distribution
+ TH3D * fhVertex; //vertex position
//leading particles
TH1F * fhPtLeading; //! pT distribution of leading particles
TH2F * fhPhiLeading; //! phi distribution vs pT of leading particles
TH2F * fhEtaLeading; //! eta distribution vs pT of leading particles
+
+ //trigger-charged histograms
TH2F * fhDeltaPhiDeltaEtaCharged ; //! differences of eta and phi between trigger and charged hadrons
- TH2F * fhDeltaPhiDeltaEtaNeutral ; //! differences of eta and phi between trigger and neutral hadrons (pi0)
-
TH2F * fhPhiCharged ; //! Phi distribution of charged particles
- TH2F * fhPhiNeutral ; //! Phi distribution of neutral particles
TH2F * fhEtaCharged ; //! Eta distribution of charged particles
- TH2F * fhEtaNeutral ; //! Eta distribution of neutral particles
TH2F * fhDeltaPhiCharged ; //! Difference of charged particle phi and trigger particle phi as function of trigger particle pT
- TH2F * fhDeltaPhiNeutral ; //! Difference of neutral particle phi and trigger particle phi as function of trigger particle pT
TH2F * fhDeltaEtaCharged ; //! Difference of charged particle eta and trigger particle eta as function of trigger particle pT
- TH2F * fhDeltaEtaNeutral ; //! Difference of neutral particle eta and trigger particle eta as function of trigger particle pT
TH2F * fhDeltaPhiChargedPt ; //! Difference of charged particle phi and trigger particle phi as function of charged particle pT
- TH2F * fhDeltaPhiNeutralPt ; //! Difference of neutral particle phi and trigger particle phi as function of neutral particle particle pT
TH2F * fhDeltaPhiUeChargedPt ; //! Difference of charged particle from underlying events phi and trigger particle phi as function of charged particle pT
- TH2F * fhDeltaPhiUeNeutralPt ; //! Difference of neutral particle phi and trigger particle phi as function of neutral particle particle pT
-
- TH2F * fhPtImbalanceNeutral ; //! Trigger particle - neutral hadron momentum imbalance histogram
TH2F * fhPtImbalanceCharged ; //! Trigger particle -charged hadron momentim imbalance histogram
TH2F * fhPtImbalanceUeCharged ; //! Trigger particle -underlying charged hadron momentim imbalance histogram
- TH2F * fhPtImbalanceUeNeutral ; //! Trigger particle - neutral hadron momentum imbalance histogram
-
+ TH2F * fhPtImbalancePosCharged ; //! Trigger particle -positive charged hadron momentim imbalance histogram
+ TH2F * fhPtImbalanceNegCharged ; //! Trigger particle -negative charged hadron momentim imbalance histogram
//with different imblance varible defination HBP distribution
TH2F * fhPtHbpCharged ; //! Trigger particle -charged hadron momentim HBP histogram
TH2F * fhPtHbpUeCharged ; //! Trigger particle -underlying charged hadron momentim HBP histogram
+
+ //if several UE calculation is on, most useful for jet-jet events contribution
+ TH2F * fhDeltaPhiUeLeftCharged ; //! Difference of charged particle from underlying events phi and trigger particle phi as function of charged particle pT
+ TH2F * fhDeltaPhiUeRightCharged ; //! Difference of charged particle from underlying events phi and trigger particle phi
+ TH2F * fhPtImbalanceUeLeftCharged ; //! Trigger particle -underlying charged hadron momentim imbalance histogram
+ TH2F * fhPtImbalanceUeRightCharged ; //! Trigger particle -underlying charged hadron momentim imbalance histogram
+ TH2F * fhPtHbpUeLeftCharged ; //! Trigger particle -underlying charged hadron momentim HBP histogram
+ TH2F * fhPtHbpUeRightCharged ; //! Trigger particle -underlying charged hadron momentim HBP histogram
+
+ //for pout and kt extraction
+ TH2F * fhPoutTrig ; // Pout =associated pt*sin(delta phi) distribution vs trigger pt
+ TH2F * fhPtTrigCharged ; //trigger and correlated particl pt, to be used for mean value for kt
+
+ //if different multiplicity analysis asked
+ TH3D ** fhTrigDeltaPhiDeltaEtaCharged ; //! differences of eta and phi between trigger and charged hadrons
+ TH2F ** fhTrigCorr ; //! Trigger particle -charged hadron momentim imbalance histogram
+ TH2F ** fhTrigUeCorr ; //! Trigger particle -UE charged hadron momentim imbalance histogram
+
+ //trigger-neutral histograms
+ TH2F * fhDeltaPhiDeltaEtaNeutral ; //! differences of eta and phi between trigger and neutral hadrons (pi0)
+ TH2F * fhPhiNeutral ; //! Phi distribution of neutral particles
+ TH2F * fhEtaNeutral ; //! Eta distribution of neutral particles
+ TH2F * fhDeltaPhiNeutral ; //! Difference of neutral particle phi and trigger particle phi as function of trigger particle pT
+ TH2F * fhDeltaEtaNeutral ; //! Difference of neutral particle eta and trigger particle eta as function of trigger particle pT
+ TH2F * fhDeltaPhiNeutralPt ; //! Difference of neutral particle phi and trigger particle phi as function of neutral particle particle pT
+ TH2F * fhDeltaPhiUeNeutralPt ; //! Difference of neutral particle phi and trigger particle phi as function of neutral particle particle pT
+ TH2F * fhPtImbalanceNeutral ; //! Trigger particle - neutral hadron momentum imbalance histogram
+ TH2F * fhPtImbalanceUeNeutral ; //! Trigger particle - neutral hadron momentum imbalance histogram
+ //with different imblance varible defination HBP distribution
TH2F * fhPtHbpNeutral ; //! Trigger particle -neutral particle momentim HBP histogram
TH2F * fhPtHbpUeNeutral ; //! Trigger particle -underlying neutral hadron momentim HBP histogram
//if several UE calculation is on, most useful for jet-jet events contribution
- TH2F * fhDeltaPhiUeLeftCharged ; //! Difference of charged particle from underlying events phi and trigger particle phi as function of charged particle pT
- TH2F * fhDeltaPhiUeRightCharged ; //! Difference of charged particle from underlying events phi and trigger particle phi
TH2F * fhDeltaPhiUeLeftNeutral ; //! Difference of charged particle from underlying events phi and trigger particle phi as function of neutral particle pT
TH2F * fhDeltaPhiUeRightNeutral ; //! Difference of charged particle from underlying events phi and trigger particle phi
- TH2F * fhPtImbalanceUeLeftCharged ; //! Trigger particle -underlying charged hadron momentim imbalance histogram
- TH2F * fhPtImbalanceUeRightCharged ; //! Trigger particle -underlying charged hadron momentim imbalance histogram
TH2F * fhPtImbalanceUeLeftNeutral ; //! Trigger particle -underlying neutral hadron momentim imbalance histogram
TH2F * fhPtImbalanceUeRightNeutral ; //! Trigger particle -underlying neutral hadron momentim imbalance histogram
- TH2F * fhPtHbpUeLeftCharged ; //! Trigger particle -underlying charged hadron momentim HBP histogram
- TH2F * fhPtHbpUeRightCharged ; //! Trigger particle -underlying charged hadron momentim HBP histogram
TH2F * fhPtHbpUeLeftNeutral ; //! Trigger particle -underlying neutral hadron momentim HBP histogram
- TH2F * fhPtHbpUeRightNeutral ; //! Trigger particle -underlying neutral hadron momentim HBP histogram
-
+ TH2F * fhPtHbpUeRightNeutral ; //! Trigger particle -underlying neutral hadron momentim HBP histogram
+
//for decay photon trigger correlation
TH3D * fhPtPi0DecayRatio ; //! for pi0 pt and ratio of decay photon pt
- TH2F * fhDeltaPhiDecay1Charged ; //! Difference of charged particle phi and decay trigger
- TH2F * fhDeltaPhiDecay2Charged ; //! Difference of charged particle phi and decay trigger
- TH2F * fhPtImbalanceDecay1Charged ; //! Trigger particle (decay 1 from pi0)-charged hadron momentim imbalance histogram
- TH2F * fhPtImbalanceDecay2Charged ; //! Trigger particle (decay 2 from pi0) -charged hadron momentim imbalance histogram
-
- TH2F * fhDeltaPhiDecay1Neutral ; //! Difference of neutral particle phi and decay trigger
- TH2F * fhDeltaPhiDecay2Neutral ; //! Difference of neutral particle phi and decay trigger
- TH2F * fhPtImbalanceDecay1Neutral ; //! Trigger particle (decay 1 from pi0)-neutral hadron momentim imbalance histogram
- TH2F * fhPtImbalanceDecay2Neutral ; //! Trigger particle (decay 2 from pi0) -neutral hadron momentim imbalance histogram
+ TH2F * fhDeltaPhiDecayCharged ; //! Difference of charged particle phi and decay trigger
+ TH2F * fhPtImbalanceDecayCharged ; //! Trigger particle (decay from pi0)-charged hadron momentim imbalance histogram
+ TH2F * fhDeltaPhiDecayNeutral ; //! Difference of neutral particle phi and decay trigger
+ TH2F * fhPtImbalanceDecayNeutral ; //! Trigger particle (decay from pi0)-neutral hadron momentim imbalance histogram
+
- ClassDef(AliAnaParticleHadronCorrelation,4)
+ ClassDef(AliAnaParticleHadronCorrelation,5)
} ;
// (see AliRoot versions previous Release 4-09)
//
// -- Author: Gustavo Conesa (LNF-INFN)
+
+//-Yaxian Mao (add the possibility for different IC method with different pt range, 01/10/2010)
//////////////////////////////////////////////////////////////////////////////
AliAnaPartCorrBaseClass(), fCalorimeter(""),
fReMakeIC(0), fMakeSeveralIC(0), fMakeInvMass(0),
fhPtIso(0),fhPhiIso(0),fhEtaIso(0), fhConeSumPt(0), fhPtInCone(0),
+ fhFRConeSumPt(0), fhPtInFRCone(0),
//Several IC
fNCones(0),fNPtThresFrac(0), fConeSizes(), fPtThresholds(), fPtFractions(),
//MC
fhPtInCone->SetXTitle("p_{T} (GeV/c)");
outputContainer->Add(fhPtInCone) ;
+ fhFRConeSumPt = new TH2F
+ ("hFRConePtSum","#Sigma p_{T} in the froward region isolation cone ",nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
+ fhFRConeSumPt->SetYTitle("#Sigma p_{T}");
+ fhFRConeSumPt->SetXTitle("p_{T} (GeV/c)");
+ outputContainer->Add(fhFRConeSumPt) ;
+
+ fhPtInFRCone = new TH2F
+ ("hPtInFRCone","p_{T} in froward region isolation cone ",nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
+ fhPtInFRCone->SetYTitle("p_{T in cone} (GeV/c)");
+ fhPtInFRCone->SetXTitle("p_{T} (GeV/c)");
+ outputContainer->Add(fhPtInFRCone) ;
+
fhPtIso = new TH1F("hPt","Isolated Number of particles",nptbins,ptmin,ptmax);
fhPtIso->SetYTitle("N");
fhPtIso->SetXTitle("p_{T}(GeV/c)");
//If too small or too large pt, skip
if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
+ //check if it is low pt trigger particle, then adjust the isolation method
+ if(aodinput->Pt() < GetIsolationCut()->GetPtThreshold() || aodinput->Pt() < GetIsolationCut()->GetSumPtThreshold())
+ continue ; //trigger should not from underlying event
+ // if(GetIsolationCut()->GetICMethod()!=AliIsolationCut::kPtThresIC && aodinput->Pt()!=0.){
+// //low pt trigger use pt threshold IC instead of other method
+// if (aodinput->Pt()*(GetIsolationCut()->GetPtFraction())<GetIsolationCut()->GetPtThreshold() || aodinput->Pt()*(GetIsolationCut()->GetPtFraction())<GetIsolationCut()->GetSumPtThreshold()) {
+// // printf("change the IC method to PtThresIC\n") ;
+// GetIsolationCut()->SetICMethod(AliIsolationCut::kPtThresIC);
+// }
+// }
+
//Check invariant mass, if pi0, skip.
Bool_t decay = kFALSE ;
if(fMakeInvMass) decay = CheckInvMass(iaod,aodinput);
Float_t ptcluster = aod->Pt();
Float_t phicluster = aod->Phi();
Float_t etacluster = aod->Eta();
+ // Float_t phiForward = aod->Phi()+TMath::PiOver2() ;
+ Float_t conesize = GetIsolationCut()->GetConeSize();
+
//Recover reference arrays with clusters and tracks
TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
-
//If too small or too large pt, skip
if(aod->Pt() < GetMinPt() || aod->Pt() > GetMaxPt() ) continue ;
-
+
+
if(fMakeSeveralIC) {
//Analysis of multiple IC at same time
MakeSeveralICAnalysis(aod);
//Fill pt distribution of particles in cone
//Tracks
coneptsum=0;
+ Double_t sumptFR = 0. ;
+ TObjArray * trackList = GetAODCTS() ;
+ for(Int_t itrack=0; itrack < trackList->GetEntriesFast(); itrack++){
+ AliAODTrack* track = (AliAODTrack *) trackList->At(itrack);
+ //fill the histograms at forward range
+ if(TMath::Sqrt((phicluster+TMath::PiOver2()-track->Phi())*(phicluster+TMath::PiOver2()-track->Phi())+(etacluster-track->Eta())*(etacluster-track->Eta())) < conesize){
+ fhPtInFRCone->Fill(ptcluster,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
+ sumptFR+=track->Pt();
+ }
+ if(TMath::Sqrt((phicluster-TMath::PiOver2()-track->Phi())*(phicluster-TMath::PiOver2()-track->Phi())+(etacluster-track->Eta())*(etacluster-track->Eta())) < conesize){
+ fhPtInFRCone->Fill(ptcluster,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
+ sumptFR+=track->Pt();
+ }
+ }
+ fhFRConeSumPt->Fill(ptcluster,sumptFR);
if(reftracks){
for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++){
AliAODTrack* track = (AliAODTrack *) reftracks->At(itrack);
coneptsum+=track->Pt();
}
}
-
+
//CaloClusters
if(refclusters){
TLorentzVector mom ;
TH2F * fhEtaIso ; //! eta of isolated particles
TH2F * fhConeSumPt ; //! Sum Pt in the cone
TH2F * fhPtInCone ; //! Particle Pt in the cone
+ TH2F * fhFRConeSumPt ; //! Sum Pt in the forward region cone (phi +90)
+ TH2F * fhPtInFRCone ; //! Particle Pt in the forward region cone (phi +90 )
//Prompt photon analysis data members for multiple cones and pt thresholds
Int_t fNCones ; //! Number of cone sizes to test
//Do photon analysis and fill aods
// Double_t vertex2[] = {0,0,0} ; //vertex from second input aod
+ //Get the vertex and check it is not too large in z, cut for SE
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+ if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;
//Select the Calorimeter of the photon
TObjArray * pl = 0x0;
Int_t evtIndex = 0 ;
if (GetMixedEvent()) {
evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
+ //Get the vertex and check it is not too large in z
+ if(TMath::Abs(GetVertex(evtIndex)[2])> GetZvertexCut()) continue ;
}
//Cluster selection, not charged, with photon id and in fiducial cut
//________________________________________________________________________________________________________________________________________________
AliAnaPi0::AliAnaPi0() : AliAnaPartCorrBaseClass(),
-fDoOwnMix(kFALSE),fNCentrBin(0),fNZvertBin(0),fNrpBin(0),
-fNPID(0),fNmaxMixEv(0), fZvtxCut(0.),fCalorimeter(""),
+fDoOwnMix(kFALSE),fNCentrBin(0),//fNZvertBin(0),fNrpBin(0),
+fNPID(0),fNmaxMixEv(0), fCalorimeter(""),
fNModules(12), fUseAngleCut(kFALSE), fEventsList(0x0), fMultiCutAna(kFALSE),
fNPtCuts(0),fPtCuts(0x0),fNAsymCuts(0),fAsymCuts(0x0),
fNCellNCuts(0),fCellNCuts(0x0),fNPIDBits(0),fPIDBits(0x0),fhReMod(0x0),
if(fDoOwnMix && fEventsList){
for(Int_t ic=0; ic<fNCentrBin; ic++){
- for(Int_t iz=0; iz<fNZvertBin; iz++){
- for(Int_t irp=0; irp<fNrpBin; irp++){
- fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp]->Delete() ;
- delete fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp] ;
+ for(Int_t iz=0; iz<GetNZvertBin(); iz++){
+ for(Int_t irp=0; irp<GetNRPBin(); irp++){
+ fEventsList[ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp]->Delete() ;
+ delete fEventsList[ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp] ;
}
}
}
AddToHistogramsName("AnaPi0_");
fNModules = 12; // set maximum to maximum number of EMCAL modules
fNCentrBin = 1;
- fNZvertBin = 1;
- fNrpBin = 1;
+// fNZvertBin = 1;
+// fNrpBin = 1;
fNPID = 9;
fNmaxMixEv = 10;
- fZvtxCut = 40;
+
fCalorimeter = "PHOS";
fUseAngleCut = kFALSE;
parList+=onePar ;
snprintf(onePar,buffersize,"Number of bins in Centrality: %d \n",fNCentrBin) ;
parList+=onePar ;
- snprintf(onePar,buffersize,"Number of bins in Z vert. pos: %d \n",fNZvertBin) ;
+ snprintf(onePar,buffersize,"Number of bins in Z vert. pos: %d \n",GetNZvertBin()) ;
parList+=onePar ;
- snprintf(onePar,buffersize,"Number of bins in Reac. Plain: %d \n",fNrpBin) ;
+ snprintf(onePar,buffersize,"Number of bins in Reac. Plain: %d \n",GetNRPBin()) ;
parList+=onePar ;
snprintf(onePar,buffersize,"Depth of event buffer: %d \n",fNmaxMixEv) ;
parList+=onePar ;
parList+=onePar ;
snprintf(onePar,buffersize,"Cuts: \n") ;
parList+=onePar ;
- snprintf(onePar,buffersize,"Z vertex position: -%f < z < %f \n",fZvtxCut,fZvtxCut) ;
+ snprintf(onePar,buffersize,"Z vertex position: -%f < z < %f \n",GetZvertexCut(),GetZvertexCut()) ;
parList+=onePar ;
snprintf(onePar,buffersize,"Calorimeter: %s \n",fCalorimeter.Data()) ;
parList+=onePar ;
// store them in fOutputContainer
//create event containers
- fEventsList = new TList*[fNCentrBin*fNZvertBin*fNrpBin] ;
+ fEventsList = new TList*[fNCentrBin*GetNZvertBin()*GetNRPBin()] ;
for(Int_t ic=0; ic<fNCentrBin; ic++){
- for(Int_t iz=0; iz<fNZvertBin; iz++){
- for(Int_t irp=0; irp<fNrpBin; irp++){
- fEventsList[ic*fNZvertBin*fNrpBin+iz*fNrpBin+irp] = new TList() ;
+ for(Int_t iz=0; iz<GetNZvertBin(); iz++){
+ for(Int_t irp=0; irp<GetNRPBin(); irp++){
+ fEventsList[ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp] = new TList() ;
}
}
}
}// multi cuts analysis
fhEvents=new TH3D("hEvents","Number of events",fNCentrBin,0.,1.*fNCentrBin,
- fNZvertBin,0.,1.*fNZvertBin,fNrpBin,0.,1.*fNrpBin) ;
+ GetNZvertBin(),0.,1.*GetNZvertBin(),GetNRPBin(),0.,1.*GetNRPBin()) ;
outputContainer->Add(fhEvents) ;
fhRealOpeningAngle = new TH2D
AliAnaPartCorrBaseClass::Print(" ");
printf("Number of bins in Centrality: %d \n",fNCentrBin) ;
- printf("Number of bins in Z vert. pos: %d \n",fNZvertBin) ;
- printf("Number of bins in Reac. Plain: %d \n",fNrpBin) ;
+ printf("Number of bins in Z vert. pos: %d \n",GetNZvertBin()) ;
+ printf("Number of bins in Reac. Plain: %d \n",GetNRPBin()) ;
printf("Depth of event buffer: %d \n",fNmaxMixEv) ;
printf("Number of different PID used: %d \n",fNPID) ;
printf("Cuts: \n") ;
- printf("Z vertex position: -%2.3f < z < %2.3f \n",fZvtxCut,fZvtxCut) ;
+ printf("Z vertex position: -%2.3f < z < %2.3f \n",GetZvertexCut(),GetZvertexCut()) ;
printf("Number of modules: %d \n",fNModules) ;
printf("Select pairs with their angle: %d \n",fUseAngleCut) ;
if(fMultiCutAna){
// get the event index in the mixed buffer where the photon comes from
// in case of mixing with analysis frame, not own mixing
evtIndex1 = GetEventIndex(p1, vert) ;
+ //printf("charge = %d\n", track->Charge());
if ( evtIndex1 == -1 )
return ;
if ( evtIndex1 == -2 )
//does not exist in ESD yet????
curCentrBin = 0 ;
curRPBin = 0 ;
- curZvertBin = (Int_t)(0.5*fNZvertBin*(vert[2]+fZvtxCut)/fZvtxCut) ;
+ curZvertBin = (Int_t)(0.5*GetNZvertBin()*(vert[2]+GetZvertexCut())/GetZvertexCut()) ;
fhEvents->Fill(curCentrBin+0.5,curZvertBin+0.5,curRPBin+0.5) ;
currentEvtIndex = evtIndex1 ;
}
if(fDoOwnMix){
//Fill mixed
- TList * evMixList=fEventsList[curCentrBin*fNZvertBin*fNrpBin+curZvertBin*fNrpBin+curRPBin] ;
+ TList * evMixList=fEventsList[curCentrBin*GetNZvertBin()*GetNRPBin()+curZvertBin*GetNRPBin()+curRPBin] ;
Int_t nMixed = evMixList->GetSize() ;
for(Int_t ii=0; ii<nMixed; ii++){
TClonesArray* ev2= (TClonesArray*) (evMixList->At(ii));
evtIndex = GetMixedEvent()->EventIndexForCaloCluster(part->GetCaloLabel(0)) ;
GetVertex(vert,evtIndex);
- if(vert[2] < -fZvtxCut || vert[2] > fZvtxCut)
+ if(TMath::Abs(vert[2])> GetZvertexCut())
evtIndex = -2 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
} else {// Single event
GetVertex(vert);
- if(vert[2] < -fZvtxCut || vert[2] > fZvtxCut)
+ if(TMath::Abs(vert[2])> GetZvertexCut())
evtIndex = -1 ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
else
evtIndex = 0 ;
class TObjString;
//Analysis
+#include "AliAnaPartCorrBaseClass.h"
class AliAODEvent ;
class AliESDEvent ;
class AliAODPWG4Particle ;
-#include "AliAnaPartCorrBaseClass.h"
class AliAnaPi0 : public AliAnaPartCorrBaseClass {
//Setters for parameters of event buffers
void SetNCentrBin(Int_t n=5) {fNCentrBin=n ;} //number of bins in centrality
- void SetNZvertBin(Int_t n=5) {fNZvertBin=n ;} //number of bins for vertex position
- void SetNRPBin(Int_t n=6) {fNrpBin=n ;} //number of bins in reaction plain
+// void SetNZvertBin(Int_t n=5) {fNZvertBin=n ;} //number of bins for vertex position
+// void SetNRPBin(Int_t n=6) {fNrpBin=n ;} //number of bins in reaction plain
void SetNMaxEvMix(Int_t n=20){fNmaxMixEv=n ;} //Maximal number of events for mixing
//Setters for event selection
- void SetZvertexCut(Float_t zcut=40.){fZvtxCut=zcut ;} //cut on vertex position
+ // void SetZvertexCut(Float_t zcut=40.){fZvtxCut=zcut ;} //cut on vertex position
TString GetCalorimeter() const {return fCalorimeter ; }
void SetCalorimeter(TString det) {fCalorimeter = det ; }
private:
Bool_t fDoOwnMix; // Do combinatorial background not the one provided by the frame
Int_t fNCentrBin ; // Number of bins in event container for centrality
- Int_t fNZvertBin ; // Number of bins in event container for vertex position
- Int_t fNrpBin ; // Number of bins in event container for reaction plain
+ // Int_t fNZvertBin ; // Number of bins in event container for vertex position
+// Int_t fNrpBin ; // Number of bins in event container for reaction plain
Int_t fNPID ; // Number of possible PID combinations
Int_t fNmaxMixEv ; // Maximal number of events stored in buffer for mixing
- Float_t fZvtxCut ; // Cut on vertex position
+// Float_t fZvtxCut ; // Cut on vertex position
TString fCalorimeter ; // Select Calorimeter for IM
Int_t fNModules ; // Number of EMCAL/PHOS modules, set as many histogras as modules
Bool_t fUseAngleCut ; // Select pairs depending on their opening angle
Int_t evtIndex1 = 0 ;
if(GetMixedEvent())
evtIndex1 = GetMixedEvent()->EventIndexForCaloCluster(photon1->GetCaloLabel(0)) ;
-
+ if(TMath::Abs(GetVertex(evtIndex1)[2]) > GetZvertexCut()) continue ; //vertex cut
mom1 = *(photon1->Momentum());
for(Int_t jphoton = iphoton+1; jphoton < GetInputAODBranch()->GetEntriesFast()-1; jphoton++){
evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
if(GetMixedEvent() && (evtIndex1 == evtIndex2))
continue ;
-
+ if(TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut()) continue ; //vertex cut
mom2 = *(photon2->Momentum());
//Int_t input = -1; //if -1 photons come from different files, not a pi0
//if(photon1->GetInputFileIndex() == photon2->GetInputFileIndex())
Int_t tag1 = 0;
Int_t tag2 = 0;
Int_t tag = 0;
-
+ Int_t evtIndex = 0;
if(!GetInputAODBranch()){
printf("AliAnaPi0EbE::MakeInvMassInCalorimeterAndCTS() - No input calo photons in AOD branch with name < %s > , STOP\n",GetInputAODName().Data());
abort();
}
for(Int_t jphoton = iphoton+1; jphoton < fInputAODGammaConv->GetEntriesFast()-1; jphoton++){
AliAODPWG4Particle * photon2 = (AliAODPWG4Particle*) (fInputAODGammaConv->At(jphoton));
+ if(GetMixedEvent())
+ evtIndex = GetMixedEvent()->EventIndexForCaloCluster(photon2->GetCaloLabel(0)) ;
+ if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
+
mom2 = *(photon2->Momentum());
//Int_t input = -1; //if -1 photons come from different files, not a pi0
if (GetMixedEvent()) {
evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ;
}
-
+ if(TMath::Abs(GetVertex(evtIndex)[2]) > GetZvertexCut()) continue ; //vertex cut
+
//Cluster selection, not charged, with pi0 id and in fiducial cut
//Input from second AOD?
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff