// Class for the analysis of particle - hadron correlations
// Particle (for example direct gamma) must be found in a previous analysis
//-- Author: Gustavo Conesa (LNF-INFN)
+
+// Modified by Yaxian Mao:
+// 1. add the UE subtraction for corrlation study
+// 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 possibility for event selection analysis based on vertex and multiplicity bins (10/10/2010)
+// 7. change the way of delta phi cut for UE study due to memory issue (reduce histograms)
+// 8. Add the possibility to request the absolute leading particle at the near side or not, set trigger bins, general clean-up (08/2011)
//////////////////////////////////////////////////////////////////////////////
// --- ROOT system ---
-#include "TH2F.h"
-#include "TClonesArray.h"
+//#include "TClonesArray.h"
#include "TClass.h"
+#include "TMath.h"
+#include "TH3D.h"
+#include "TDatabasePDG.h"
//---- ANALYSIS system ----
#include "AliNeutralMesonSelection.h"
#include "AliCaloTrackReader.h"
#include "AliCaloPID.h"
#include "AliAODPWG4ParticleCorrelation.h"
-#include "AliFidutialCut.h"
-#include "AliAODTrack.h"
-#include "AliAODCaloCluster.h"
+#include "AliFiducialCut.h"
+#include "AliVTrack.h"
+#include "AliVCluster.h"
#include "AliMCAnalysisUtils.h"
#include "TParticle.h"
#include "AliStack.h"
#include "AliAODMCParticle.h"
+#include "AliMixedEvent.h"
ClassImp(AliAnaParticleHadronCorrelation)
AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation():
AliAnaPartCorrBaseClass(),
fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fSelectIsolated(0),
- fMakeSeveralUE(0), fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(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),
+ fMakeSeveralUE(0),
+ fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(0.),
+ fPi0AODBranchName(""),fNeutralCorr(0), fPi0Trigger(0),
+ fMakeAbsoluteLeading(0), fLeadingTriggerIndex(-1),
+ fNTriggerPtBins(0),
+ fhPtLeading(0),fhPhiLeading(0),fhEtaLeading(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),
fhDeltaPhiUeLeftCharged(0),fhDeltaPhiUeRightCharged(0),
- fhDeltaPhiUeLeftNeutral(0),fhDeltaPhiUeRightNeutral(0),
fhPtImbalanceUeLeftCharged(0),fhPtImbalanceUeRightCharged(0),
- fhPtImbalanceUeLeftNeutral(0),fhPtImbalanceUeRightNeutral(0)
+ fhPtHbpUeLeftCharged(0),fhPtHbpUeRightCharged(0),
+ fhPtTrigPout(0), fhPtAssocDeltaPhi(0),
+ fhPtTrigCharged(0),
+ fhTrigDeltaPhiCharged(0x0), fhTrigDeltaEtaCharged(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),
+ fhPtHbpUeLeftNeutral(0),fhPtHbpUeRightNeutral(0),
+ fhPtPi0DecayRatio(0),
+ fhDeltaPhiDecayCharged(0),
+ fhPtImbalanceDecayCharged(0),
+ fhDeltaPhiDecayNeutral(0),
+ fhPtImbalanceDecayNeutral(0),
+ fh2phiLeadingParticle(0x0),
+ fhMCLeadingCount(0),
+ fhMCEtaCharged(0),
+ fhMCPhiCharged(0),
+ fhMCDeltaEtaCharged(0),
+ fhMCDeltaPhiCharged(0x0),
+ fhMCDeltaPhiDeltaEtaCharged(0),
+ fhMCDeltaPhiChargedPt(0),
+ fhMCPtImbalanceCharged(0),
+ fhMCPtHbpCharged(0),
+ fhMCPtTrigPout(0),
+ fhMCPtAssocDeltaPhi(0)
{
//Default Ctor
InitParameters();
}
-//____________________________________________________________________________
-AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation(const AliAnaParticleHadronCorrelation & g):
- AliAnaPartCorrBaseClass(g),
- fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut),
- fSelectIsolated(g.fSelectIsolated),
- fMakeSeveralUE(g.fMakeSeveralUE), fUeDeltaPhiMaxCut(g.fUeDeltaPhiMaxCut),
- fUeDeltaPhiMinCut(g.fUeDeltaPhiMinCut),
- fhPhiCharged(g.fhPhiCharged), fhPhiNeutral(g.fhPhiNeutral),
- fhEtaCharged(g.fhEtaCharged), fhEtaNeutral(g.fhEtaNeutral),
- fhDeltaPhiCharged(g.fhDeltaPhiCharged),
- fhDeltaPhiNeutral(g.fhDeltaPhiNeutral),
- fhDeltaEtaCharged(g.fhDeltaEtaCharged),
- fhDeltaEtaNeutral(g.fhDeltaEtaNeutral),
- fhDeltaPhiChargedPt(g.fhDeltaPhiChargedPt),
- fhDeltaPhiNeutralPt(g.fhDeltaPhiNeutralPt),
- fhDeltaPhiUeChargedPt(g.fhDeltaPhiUeChargedPt),
- fhDeltaPhiUeNeutralPt(g.fhDeltaPhiUeNeutralPt),
- fhPtImbalanceNeutral(g.fhPtImbalanceNeutral),
- fhPtImbalanceCharged(g.fhPtImbalanceCharged),
- fhPtImbalanceUeCharged(g.fhPtImbalanceUeCharged),
- fhPtImbalanceUeNeutral(g.fhPtImbalanceUeNeutral),
- fhDeltaPhiUeLeftCharged(g.fhDeltaPhiUeLeftCharged),
- fhDeltaPhiUeRightCharged(g.fhDeltaPhiUeRightCharged),
- fhDeltaPhiUeLeftNeutral(g.fhDeltaPhiUeLeftNeutral),
- fhDeltaPhiUeRightNeutral(g.fhDeltaPhiUeRightNeutral),
- fhPtImbalanceUeLeftCharged(g.fhPtImbalanceUeLeftCharged),
- fhPtImbalanceUeRightCharged(g.fhPtImbalanceUeRightCharged),
- fhPtImbalanceUeLeftNeutral(g.fhPtImbalanceUeLeftNeutral),
- fhPtImbalanceUeRightNeutral(g.fhPtImbalanceUeRightNeutral)
-{
- // cpy ctor
-
-}
-
-//_________________________________________________________________________
-AliAnaParticleHadronCorrelation & AliAnaParticleHadronCorrelation::operator = (const AliAnaParticleHadronCorrelation & source)
-{
- // assignment operator
-
- if(this == &source)return *this;
- ((AliAnaPartCorrBaseClass *)this)->operator=(source);
-
- fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ;
- fDeltaPhiMinCut = source.fDeltaPhiMinCut ;
- fSelectIsolated = source.fSelectIsolated ;
- fMakeSeveralUE = source.fMakeSeveralUE ;
- fUeDeltaPhiMaxCut = source.fUeDeltaPhiMaxCut ;
- fUeDeltaPhiMinCut = source.fUeDeltaPhiMinCut ;
- fhPhiCharged = source.fhPhiCharged ; fhPhiNeutral = source.fhPhiNeutral ;
- fhEtaCharged = source.fhEtaCharged ; fhEtaNeutral = source.fhEtaNeutral ;
- fhDeltaPhiCharged = source.fhDeltaPhiCharged ;
- fhDeltaPhiNeutral = source.fhDeltaPhiNeutral ;
- fhDeltaEtaCharged = source.fhDeltaEtaCharged ;
- fhDeltaEtaNeutral = source.fhDeltaEtaNeutral ;
- fhDeltaPhiChargedPt = source.fhDeltaPhiChargedPt ;
- fhDeltaPhiNeutralPt = source.fhDeltaPhiNeutralPt ;
- fhDeltaPhiUeChargedPt = source.fhDeltaPhiUeChargedPt ;
- fhDeltaPhiUeNeutralPt = source.fhDeltaPhiUeNeutralPt ;
- fhPtImbalanceNeutral = source.fhPtImbalanceNeutral ;
- fhPtImbalanceCharged = source.fhPtImbalanceCharged ;
- fhPtImbalanceUeCharged = source.fhPtImbalanceUeCharged ;
- fhPtImbalanceUeNeutral = source.fhPtImbalanceUeNeutral ;
- fhDeltaPhiUeLeftCharged = source.fhDeltaPhiUeLeftCharged ;
- fhDeltaPhiUeRightCharged = source.fhDeltaPhiUeRightCharged ;
- fhDeltaPhiUeLeftNeutral = source.fhDeltaPhiUeLeftNeutral ;
- fhDeltaPhiUeRightNeutral = source.fhDeltaPhiUeRightNeutral ;
- fhPtImbalanceUeLeftCharged = source.fhPtImbalanceUeLeftCharged ;
- fhPtImbalanceUeRightCharged = source.fhPtImbalanceUeRightCharged ;
- fhPtImbalanceUeLeftNeutral = source.fhPtImbalanceUeLeftNeutral ;
- fhPtImbalanceUeRightNeutral = source.fhPtImbalanceUeRightNeutral ;
- return *this;
-
-}
-
//________________________________________________________________________
TList * AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
{
TList * outputContainer = new TList() ;
outputContainer->SetName("CorrelationHistos") ;
- Int_t nptbins = GetHistoNPtBins();
- Int_t nphibins = GetHistoNPhiBins();
- Int_t netabins = GetHistoNEtaBins();
+ Int_t nptbins = GetHistoPtBins();
+ Int_t nphibins = GetHistoPhiBins();
+ Int_t netabins = GetHistoEtaBins();
Float_t ptmax = GetHistoPtMax();
Float_t phimax = GetHistoPhiMax();
Float_t etamax = GetHistoEtaMax();
Float_t phimin = GetHistoPhiMin();
Float_t etamin = GetHistoEtaMin();
+ fhPtLeading = new TH1F ("hPtLeading","p_T distribution of leading particles", nptbins,ptmin,ptmax);
+ fhPtLeading->SetXTitle("p_{T}^{trig} (GeV/c)");
+
+ fhPhiLeading = new TH2F ("hPhiLeading","#phi distribution of leading Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
+ fhPhiLeading->SetYTitle("#phi (rad)");
+
+ fhEtaLeading = new TH2F ("hEtaLeading","#eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax);
+ fhEtaLeading->SetYTitle("#eta ");
+
+ outputContainer->Add(fhPtLeading);
+ outputContainer->Add(fhPhiLeading);
+ outputContainer->Add(fhEtaLeading);
+
//Correlation with charged hadrons
if(GetReader()->IsCTSSwitchedOn()) {
+ fhDeltaPhiDeltaEtaCharged = new TH2F
+ ("DeltaPhiDeltaEtaCharged","#phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
+ 140,-2.,5.,200,-2,2);
+ fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
+ fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
+
fhPhiCharged = new TH2F
- ("PhiCharged","#phi_{h^{#pm}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,nphibins,phimin,phimax);
+ ("PhiCharged","#phi_{h^{#pm}} vs p_{T #pm}",
+ nptbins,ptmin,ptmax,nphibins,phimin,phimax);
fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
- fhPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
+ fhPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
fhEtaCharged = new TH2F
- ("EtaCharged","#eta_{h^{#pm}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,netabins,etamin,etamax);
+ ("EtaCharged","#eta_{h^{#pm}} vs p_{T #pm}",
+ nptbins,ptmin,ptmax,netabins,etamin,etamax);
fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
- fhEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
+ fhEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
fhDeltaPhiCharged = new TH2F
- ("DeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiCharged->SetYTitle("#Delta #phi");
fhDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
fhDeltaPhiChargedPt = new TH2F
- ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
fhDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
-
+
fhDeltaPhiUeChargedPt = new TH2F
- ("DeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}}",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi");
fhDeltaPhiUeChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
fhDeltaEtaCharged = new TH2F
- ("DeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,200,-2,2);
+ ("DeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",
+ nptbins,ptmin,ptmax,200,-2,2);
fhDeltaEtaCharged->SetYTitle("#Delta #eta");
fhDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
fhPtImbalanceCharged =
- new TH2F("CorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceCharged->SetYTitle("z_{trigger h^{#pm}}");
fhPtImbalanceCharged->SetXTitle("p_{T trigger}");
fhPtImbalanceUeCharged =
- new TH2F("CorrelationUeCharged","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger}",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeCharged","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
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.);
+ fhPtHbpCharged->SetYTitle("ln(1/x_{E})");
+ fhPtHbpCharged->SetXTitle("p_{T trigger}");
+
+ fhPtHbpUeCharged =
+ new TH2F("HbpUeCharged","#xi = ln(1/x_{E}) with charged hadrons",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeCharged->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeCharged->SetXTitle("p_{T trigger}");
+
+ fhPtTrigPout =
+ new TH2F("PtTrigPout","Pout with triggers",
+ nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
+ fhPtTrigPout->SetYTitle("p_{out} (GeV/c)");
+ fhPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)");
+
+ fhPtAssocDeltaPhi =
+ new TH2F("fhPtAssocDeltaPhi"," charged hadrons vs. delta phi",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhPtAssocDeltaPhi->SetXTitle("p_{T h^{#pm}} (GeV/c)");
+ fhPtAssocDeltaPhi->SetYTitle("#Delta #phi (GeV/c)");
+
+ fhPtTrigCharged =
+ new TH2F("PtTrigCharged","trgger and charged tracks pt distribution",
+ nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
+ fhPtTrigCharged->SetYTitle("p_{T h^{#pm}} (GeV/c)");
+ fhPtTrigCharged->SetXTitle("p_{T trigger} (GeV/c)");
+
+ outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
outputContainer->Add(fhPhiCharged) ;
outputContainer->Add(fhEtaCharged) ;
outputContainer->Add(fhDeltaPhiCharged) ;
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(fhPtTrigPout) ;
+ outputContainer->Add(fhPtAssocDeltaPhi) ;
+ outputContainer->Add(fhPtTrigCharged) ;
+
+ if(DoEventSelect()){
+ Int_t nMultiBins = GetMultiBin();
+ fhTrigDeltaPhiCharged = new TH2F*[nMultiBins] ;
+ fhTrigDeltaEtaCharged = new TH2F*[nMultiBins] ;
+ fhTrigCorr = new TH2F*[nMultiBins];
+ fhTrigUeCorr = new TH2F*[nMultiBins];
+ for(Int_t im=0; im<nMultiBins; im++){
+ fhTrigDeltaPhiCharged[im] = new TH2F
+ (Form("fhTrigDeltaPhiCharged_%d",im),Form("fhTrigDeltaPhiCharged_%d",im), nptbins,ptmin,ptmax, 140,-2.,5.);
+ fhTrigDeltaPhiCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
+ fhTrigDeltaPhiCharged[im]->SetYTitle("#Delta #phi");
+ fhTrigDeltaEtaCharged[im] = new TH2F
+ (Form("fhTrigDeltaEtaCharged_%d",im),Form("fhTrigDeltaEtaCharged_%d",im), nptbins,ptmin,ptmax, 200,-2,2);
+ fhTrigDeltaEtaCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
+ fhTrigDeltaEtaCharged[im]->SetYTitle("#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(fhTrigDeltaPhiCharged[im]) ;
+ outputContainer->Add(fhTrigDeltaEtaCharged[im]) ;
+ outputContainer->Add(fhTrigCorr[im]);
+ outputContainer->Add(fhTrigUeCorr[im]);
+
+ }
+ }
+
+ for(Int_t i = 0 ; i < fNTriggerPtBins ; i++){
+ fhTrigPt[i] = new TH1F(Form("fhTrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ Form("fhTrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ 20, fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]);
+ fhDphiTrigPtAssocPt[i] = new TH1F(Form("fhDphiTrigPtAssocPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ Form("fhDphiTrigPtAssocPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ 288, -1.0/3.0, 5.0/3.0);
+ fhAssocPtBkg[i] = new TH1F(Form("fhAssocPtBkg%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ Form("fhAssocPtBkg%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ 600, 0., 30.0);
+ fhXE[i] = new TH1F(Form("fhXETrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ Form("fhXETrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
+ 50, 0.0, 2.0);
+ outputContainer->Add(fhTrigPt[i]) ;
+ outputContainer->Add(fhDphiTrigPtAssocPt[i]) ;
+ outputContainer->Add(fhAssocPtBkg[i]);
+ outputContainer->Add(fhXE[i]);
+ }
+
+ if(fPi0Trigger){
+ fhPtPi0DecayRatio = new TH2F
+ ("hPtPi0DecayRatio","p_T of #pi^{0} and the ratio of pt for two decay",
+ nptbins,ptmin,ptmax, 100,0.,2.);
+ fhPtPi0DecayRatio->SetXTitle("p_{T}^{#pi^{0}} (GeV/c)");
+ fhPtPi0DecayRatio->SetYTitle("p_{T}^{Decay}/p_{T}^{#pi^{0}}");
+
+ fhDeltaPhiDecayCharged = new TH2F
+ ("DeltaPhiDecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhDeltaPhiDecayCharged->SetYTitle("#Delta #phi");
+ fhDeltaPhiDecayCharged->SetXTitle("p_{T Decay} (GeV/c)");
+
+ fhPtImbalanceDecayCharged =
+ new TH2F("CorrelationDecayCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
+ nptbins,ptmin,ptmax,200,0.,2.);
+ fhPtImbalanceDecayCharged->SetYTitle("z_{decay h^{#pm}}");
+ fhPtImbalanceDecayCharged->SetXTitle("p_{T decay}");
+
+ outputContainer->Add(fhPtPi0DecayRatio) ;
+ outputContainer->Add(fhDeltaPhiDecayCharged) ;
+ outputContainer->Add(fhPtImbalanceDecayCharged) ;
+ }
+
+
if(fMakeSeveralUE){
fhDeltaPhiUeLeftCharged = new TH2F
- ("DeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE left side range of trigger particles",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE left side range of trigger particles",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi");
fhDeltaPhiUeLeftCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
-
+
fhDeltaPhiUeRightCharged = new TH2F
- ("DeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE right side range of trigger particles",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE right side range of trigger particles",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeRightCharged->SetYTitle("#Delta #phi");
fhDeltaPhiUeRightCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
outputContainer->Add(fhDeltaPhiUeRightCharged) ;
fhPtImbalanceUeLeftCharged =
- new TH2F("CorrelationUeChargedLeft","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE left side of trigger",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeChargedLeft","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE left side of trigger",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceUeLeftCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
fhPtImbalanceUeLeftCharged->SetXTitle("p_{T trigger}");
outputContainer->Add(fhPtImbalanceUeLeftCharged) ;
fhPtImbalanceUeRightCharged =
- new TH2F("CorrelationUeChargedRight","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE right side of trigger",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeChargedRight","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE right side of trigger",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceUeRightCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
fhPtImbalanceUeRightCharged->SetXTitle("p_{T trigger}");
outputContainer->Add(fhPtImbalanceUeRightCharged) ;
+
+ fhPtHbpUeLeftCharged =
+ new TH2F("HbpUeChargedLeft","#xi = ln(1/x_{E}) with charged UE left side of trigger",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeLeftCharged->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeLeftCharged->SetXTitle("p_{T trigger}");
+ outputContainer->Add(fhPtHbpUeLeftCharged) ;
+
+ fhPtHbpUeRightCharged =
+ new TH2F("HbpUeChargedRight","#xi = ln(1/x_{E}) with charged UE right side of trigger",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeRightCharged->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeRightCharged->SetXTitle("p_{T trigger}");
+ outputContainer->Add(fhPtHbpUeRightCharged) ;
+
}
} //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}}",
+ 140,-2.,5.,200,-2,2);
+ fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #phi");
+ fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");
+
fhPhiNeutral = new TH2F
- ("PhiNeutral","#phi_{#pi^{0}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,nphibins,phimin,phimax);
+ ("PhiNeutral","#phi_{#pi^{0}} vs p_{T #pi^{0}}",
+ nptbins,ptmin,ptmax,nphibins,phimin,phimax);
fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
- fhPhiNeutral->SetXTitle("p_{T trigger} (GeV/c)");
+ fhPhiNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
fhEtaNeutral = new TH2F
- ("EtaNeutral","#eta_{#pi^{0}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,netabins,etamin,etamax);
+ ("EtaNeutral","#eta_{#pi^{0}} vs p_{T #pi^{0}}",
+ nptbins,ptmin,ptmax,netabins,etamin,etamax);
fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
- fhEtaNeutral->SetXTitle("p_{T trigger} (GeV/c)");
+ fhEtaNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
fhDeltaPhiNeutral = new TH2F
- ("DeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,nphibins,phimin,phimax);
+ ("DeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",
+ nptbins,ptmin,ptmax,nphibins,phimin,phimax);
fhDeltaPhiNeutral->SetYTitle("#Delta #phi");
fhDeltaPhiNeutral->SetXTitle("p_{T trigger} (GeV/c)");
fhDeltaPhiNeutralPt = new TH2F
- ("DeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi");
- fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0} (GeV/c)");
-
+ fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
+
fhDeltaPhiUeNeutralPt = new TH2F
- ("DeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi");
- fhDeltaPhiUeNeutralPt->SetXTitle("p_{T h^{0} (GeV/c)");
+ fhDeltaPhiUeNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
fhDeltaEtaNeutral = new TH2F
- ("DeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",
- nptbins,ptmin,ptmax,200,-2,2);
+ ("DeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",
+ nptbins,ptmin,ptmax,200,-2,2);
fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
fhDeltaEtaNeutral->SetXTitle("p_{T trigger} (GeV/c)");
fhPtImbalanceNeutral =
- new TH2F("CorrelationNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceNeutral->SetYTitle("z_{trigger #pi^{0}}");
fhPtImbalanceNeutral->SetXTitle("p_{T trigger}");
-
+
fhPtImbalanceUeNeutral =
- new TH2F("CorrelationUeNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceUeNeutral->SetYTitle("z_{trigger #pi^{0}}");
fhPtImbalanceUeNeutral->SetXTitle("p_{T trigger}");
-
+
+ fhPtHbpNeutral =
+ new TH2F("HbpNeutral","#xi = ln(1/x_{E}) with neutral particles",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpNeutral->SetYTitle("ln(1/x_{E})");
+ fhPtHbpNeutral->SetXTitle("p_{T trigger}");
+
+ fhPtHbpUeNeutral =
+ new TH2F("HbpUeNeutral","#xi = ln(1/x_{E}) with neutral particles",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeNeutral->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeNeutral->SetXTitle("p_{T trigger}");
+
+
+ outputContainer->Add(fhDeltaPhiDeltaEtaNeutral);
outputContainer->Add(fhPhiNeutral) ;
outputContainer->Add(fhEtaNeutral) ;
outputContainer->Add(fhDeltaPhiNeutral) ;
outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
outputContainer->Add(fhDeltaEtaNeutral) ;
outputContainer->Add(fhPtImbalanceNeutral) ;
- outputContainer->Add(fhPtImbalanceUeNeutral) ;
-
+ outputContainer->Add(fhPtImbalanceUeNeutral) ;
+ outputContainer->Add(fhPtHbpNeutral) ;
+ outputContainer->Add(fhPtHbpUeNeutral) ;
+
+ if(fPi0Trigger){
+ fhDeltaPhiDecayNeutral = new TH2F
+ ("DeltaPhiDecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhDeltaPhiDecayNeutral->SetYTitle("#Delta #phi");
+ fhDeltaPhiDecayNeutral->SetXTitle("p_{T Decay} (GeV/c)");
+
+ fhPtImbalanceDecayNeutral =
+ new TH2F("CorrelationDecayNeutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
+ nptbins,ptmin,ptmax,200,0.,2.);
+ fhPtImbalanceDecayNeutral->SetYTitle("z_{decay h^{0}}");
+ fhPtImbalanceDecayNeutral->SetXTitle("p_{T decay}");
+
+ outputContainer->Add(fhDeltaPhiDecayNeutral) ;
+ outputContainer->Add(fhPtImbalanceDecayNeutral) ;
+ }
+
+
if(fMakeSeveralUE){
fhDeltaPhiUeLeftNeutral = new TH2F
- ("DeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T h^{0}} with UE left side range of trigger particles",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T h^{0}} with neutral UE left side range of trigger particles",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi");
fhDeltaPhiUeLeftNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
-
+
fhDeltaPhiUeRightNeutral = new TH2F
- ("DeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T Ueh^{0}} with UE right side range of trigger particles",
- nptbins,ptmin,ptmax,700,-2.,5.);
+ ("DeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T Ueh^{0}} with neutral UE right side range of trigger particles",
+ nptbins,ptmin,ptmax,140,-2.,5.);
fhDeltaPhiUeRightNeutral->SetYTitle("#Delta #phi");
fhDeltaPhiUeRightNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
outputContainer->Add(fhDeltaPhiUeRightNeutral) ;
fhPtImbalanceUeLeftNeutral =
- new TH2F("CorrelationUeNeutralLeft","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with UE left side of trigger",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeNeutralLeft","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE left side of trigger",
+ nptbins,ptmin,ptmax,140,0.,2.);
fhPtImbalanceUeLeftNeutral->SetYTitle("z_{trigger Ueh^{0}}");
fhPtImbalanceUeLeftNeutral->SetXTitle("p_{T trigger}");
outputContainer->Add(fhPtImbalanceUeLeftNeutral) ;
fhPtImbalanceUeRightNeutral =
- new TH2F("CorrelationUeNeutralRight","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with UE right side of trigger",
- nptbins,ptmin,ptmax,1000,0.,2.);
+ new TH2F("CorrelationUeNeutralRight","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE right side of trigger",
+ nptbins,ptmin,ptmax,200,0.,2.);
fhPtImbalanceUeRightNeutral->SetYTitle("z_{trigger Ueh^{0}}");
fhPtImbalanceUeRightNeutral->SetXTitle("p_{T trigger}");
outputContainer->Add(fhPtImbalanceUeRightNeutral) ;
+
+ fhPtHbpUeLeftNeutral =
+ new TH2F("HbpUeNeutralLeft","#xi = ln(1/x_{E}) with neutral UE left side of trigger",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeLeftNeutral->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeLeftNeutral->SetXTitle("p_{T trigger}");
+ outputContainer->Add(fhPtHbpUeLeftNeutral) ;
+
+ fhPtHbpUeRightNeutral =
+ new TH2F("HbpUeNeutralRight","#xi = ln(1/x_{E}) with neutral UE right side of trigger",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhPtHbpUeRightNeutral->SetYTitle("ln(1/x_{E})");
+ fhPtHbpUeRightNeutral->SetXTitle("p_{T trigger}");
+ outputContainer->Add(fhPtHbpUeRightNeutral) ;
+
}
-
+
+ //if data is MC, fill more histograms
+ if(IsDataMC()){
+ fh2phiLeadingParticle=new TH2F("fh2phiLeadingParticle","#phi resolustion for trigger particles",nptbins,ptmin,ptmax,100,-1,1);
+ fh2phiLeadingParticle->GetXaxis()->SetTitle("p_{T gen Leading} (GeV/c)");
+ fh2phiLeadingParticle->GetYaxis()->SetTitle("(#phi_{rec}-#phi_{gen})/#phi_{gen}");
+
+ fhMCLeadingCount=new TH1F("MCLeadingTriggerCount","MCLeadingTriggerCount",nptbins,ptmin,ptmax);
+ fhMCLeadingCount->SetXTitle("p_{T trig}");
+
+ fhMCEtaCharged = new TH2F
+ ("MCEtaCharged","MC #eta_{h^{#pm}} vs p_{T #pm}",
+ nptbins,ptmin,ptmax,netabins,etamin,etamax);
+ fhMCEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
+ fhMCEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
+
+ fhMCPhiCharged = new TH2F
+ ("MCPhiCharged","#MC phi_{h^{#pm}} vs p_{T #pm}",
+ 200,ptmin,ptmax,nphibins,phimin,phimax);
+ fhMCPhiCharged->SetYTitle("MC #phi_{h^{#pm}} (rad)");
+ fhMCPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
+
+ fhMCDeltaPhiDeltaEtaCharged = new TH2F
+ ("MCDeltaPhiDeltaEtaCharged","#MC phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
+ 140,-2.,5.,200,-2,2);
+ fhMCDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
+ fhMCDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
+
+ fhMCDeltaEtaCharged = new TH2F
+ ("MCDeltaEtaCharged","MC #eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger} and p_{T assoc}",
+ nptbins,ptmin,ptmax,200,-2,2);
+ fhMCDeltaEtaCharged->SetYTitle("#Delta #eta");
+ fhMCDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
+
+ fhMCDeltaPhiCharged = new TH2F
+ ("MCDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhMCDeltaPhiCharged->SetYTitle("#Delta #phi");
+ fhMCDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
+
+ fhMCDeltaPhiChargedPt = new TH2F
+ ("MCDeltaPhiChargedPt","MC #phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhMCDeltaPhiChargedPt->SetYTitle("#Delta #phi");
+ fhMCDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
+
+ fhMCPtImbalanceCharged =
+ new TH2F("MCCorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
+ nptbins,ptmin,ptmax,200,0.,2.);
+ fhMCPtImbalanceCharged->SetYTitle("z_{trigger h^{#pm}}");
+ fhMCPtImbalanceCharged->SetXTitle("p_{T trigger}");
+
+ fhMCPtHbpCharged =
+ new TH2F("MCHbpCharged","MC #xi = ln(1/x_{E}) with charged hadrons",
+ nptbins,ptmin,ptmax,200,0.,10.);
+ fhMCPtHbpCharged->SetYTitle("ln(1/x_{E})");
+ fhMCPtHbpCharged->SetXTitle("p_{T trigger}");
+
+ fhMCPtTrigPout =
+ new TH2F("MCPtTrigPout","AOD MC Pout with triggers",
+ nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
+ fhMCPtTrigPout->SetYTitle("p_{out} (GeV/c)");
+ fhMCPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)");
+
+ fhMCPtAssocDeltaPhi =
+ new TH2F("fhMCPtAssocDeltaPhi","AOD MC delta phi with associated charged hadrons",
+ nptbins,ptmin,ptmax,140,-2.,5.);
+ fhMCPtAssocDeltaPhi->SetYTitle("#Delta #phi");
+ fhMCPtAssocDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
+
+ outputContainer->Add(fh2phiLeadingParticle);
+ outputContainer->Add(fhMCLeadingCount);
+ outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged);
+ outputContainer->Add(fhMCPhiCharged) ;
+ outputContainer->Add(fhMCEtaCharged) ;
+ outputContainer->Add(fhMCDeltaEtaCharged) ;
+ outputContainer->Add(fhMCDeltaPhiCharged) ;
+
+ outputContainer->Add(fhMCDeltaPhiChargedPt) ;
+ outputContainer->Add(fhMCPtImbalanceCharged) ;
+ outputContainer->Add(fhMCPtHbpCharged) ;
+ outputContainer->Add(fhMCPtTrigPout) ;
+ outputContainer->Add(fhMCPtAssocDeltaPhi) ;
+ } //for MC histogram
+
//Keep neutral meson selection histograms if requiered
//Setting done in AliNeutralMesonSelection
if(GetNeutralMesonSelection()){
TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
- for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
+ for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
+ delete nmsHistos;
}
-
+
}//Correlation with neutral hadrons
return outputContainer;
-
+
}
//____________________________________________________________________________
SetInputAODName("PWG4Particle");
SetAODObjArrayName("Hadrons");
AddToHistogramsName("AnaHadronCorr_");
-
- //Correlation with neutrals
- //SetOutputAODClassName("AliAODPWG4Particle");
- //SetOutputAODName("Pi0Correlated");
-
+
SetPtCutRange(0.,300);
- fDeltaPhiMinCut = 1.5 ;
- fDeltaPhiMaxCut = 4.5 ;
- fSelectIsolated = kFALSE;
- fMakeSeveralUE = kFALSE;
- fUeDeltaPhiMinCut = 1. ;
- fUeDeltaPhiMaxCut = 1.5 ;
+ fDeltaPhiMinCut = 1.5 ;
+ fDeltaPhiMaxCut = 4.5 ;
+ fSelectIsolated = kFALSE;
+ fMakeSeveralUE = kFALSE;
+ fUeDeltaPhiMinCut = 1. ;
+ fUeDeltaPhiMaxCut = 1.5 ;
+ fNeutralCorr = kFALSE ;
+ fPi0Trigger = kFALSE ;
+ fMakeAbsoluteLeading = kTRUE;
+
+ fNTriggerPtBins = 6 ;
+ fTriggerPtBinLimit[0] = 4.;
+ fTriggerPtBinLimit[1] = 5.;
+ fTriggerPtBinLimit[2] = 6.;
+ fTriggerPtBinLimit[3] = 7.;
+ fTriggerPtBinLimit[4] = 8.;
+ fTriggerPtBinLimit[5] =10.;
+ fTriggerPtBinLimit[6] =12.;
+ fTriggerPtBinLimit[7] =15.;
+ fTriggerPtBinLimit[8] =20.;
+ fTriggerPtBinLimit[9] =25.;
+ fTriggerPtBinLimit[6] =100.;
+
}
//__________________________________________________________________
printf("Phi trigger particle-UeHadron < %3.2f\n", fUeDeltaPhiMaxCut) ;
printf("Phi trigger particle-UeHadron > %3.2f\n", fUeDeltaPhiMinCut) ;
printf("Several UE? %d\n", fMakeSeveralUE) ;
+ printf("Name of AOD Pi0 Branch %s \n",fPi0AODBranchName.Data());
+ printf("Do Decay-hadron correlation ? %d\n", fPi0Trigger) ;
+ printf("Select absolute leading for cluster triggers ? %d\n", fMakeAbsoluteLeading) ;
+ printf("Trigger pt bins %d\n", fNTriggerPtBins) ;
+ for (Int_t ibin = 0; ibin<fNTriggerPtBins; ibin++) {
+ printf("\t bin %d = [%2.1f,%2.1f]\n", ibin, fTriggerPtBinLimit[ibin], fTriggerPtBinLimit[ibin+1]) ;
+ }
+
+}
+
+//__________________________________________________________________
+TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
+{
+ //Save parameters used for analysis
+ TString parList ; //this will be list of parameters used for this analysis.
+ const Int_t buffersize = 560;
+ char onePar[buffersize] ;
+
+ snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---\n") ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Phi trigger particle-Hadron < %3.2f ", fDeltaPhiMaxCut) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Phi trigger particle-Hadron > %3.2f ", fDeltaPhiMinCut) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Isolated Trigger? %d\n", fSelectIsolated) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Phi trigger particle-UeHadron < %3.2f ", fUeDeltaPhiMaxCut) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Phi trigger particle-UeHadron > %3.2f ", fUeDeltaPhiMinCut) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Several UE? %d\n", fMakeSeveralUE) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s ",fPi0AODBranchName.Data());
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Do Decay-hadron correlation ? %d", fPi0Trigger) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d\n", fMakeAbsoluteLeading) ;
+ parList+=onePar ;
+ snprintf(onePar,buffersize,"Trigger pt bins %d: ", fNTriggerPtBins) ;
+ parList+=onePar ;
+ for (Int_t ibin = 0; ibin<fNTriggerPtBins; ibin++) {
+ snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fTriggerPtBinLimit[ibin], fTriggerPtBinLimit[ibin+1]) ;
+ }
+ parList+=onePar ;
+
+ //Get parameters set in base class.
+ parList += GetBaseParametersList() ;
+
+ //Get parameters set in PID class.
+ //parList += GetCaloPID()->GetPIDParametersList() ;
+
+ //Get parameters set in FiducialCut class (not available yet)
+ //parlist += GetFidCut()->GetFidCutParametersList()
+
+ return new TObjString(parList) ;
+
}
+
//____________________________________________________________________________
void AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()
{
}
if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
-
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
- abort();
+ printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
+ abort();
}
if(GetDebug() > 1){
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
- printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetAODCTS()->GetEntriesFast());
- printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetAODEMCAL()->GetEntriesFast());
- printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetAODPHOS()->GetEntriesFast());
+ printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetCTSTracks()->GetEntriesFast());
+ printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetEMCALClusters()->GetEntriesFast());
+ printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetPHOSClusters()->GetEntriesFast());
}
- //Loop on stored AOD particles, trigger
- Int_t naod = 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 trigger
+ Double_t ptTrig = GetMinPt() ;
+ fLeadingTriggerIndex = -1 ;
+ Int_t naod = GetInputAODBranch()->GetEntriesFast() ;
for(Int_t iaod = 0; iaod < naod ; iaod++){
AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
-
- //Make correlation with charged hadrons
- if(GetReader()->IsCTSSwitchedOn() )
- MakeChargedCorrelation(particle, GetAODCTS(),kFALSE);
-
- //Make correlation with neutral pions
- //Trigger particle in PHOS, correlation with EMCAL
- if(particle->GetDetector()=="PHOS" && GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
- MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
- //Trigger particle in EMCAL, correlation with PHOS
- else if(particle->GetDetector()=="EMCAL" && GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
- MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
- //Trigger particle in CTS, correlation with PHOS, EMCAL and CTS
- else if(particle->GetDetector()=="CTS" ){
- if(GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
- MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
- if(GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
- MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
+
+ //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())
+ return ;
}
-
+
+ // find the leading particles with highest momentum
+ if (particle->Pt() > ptTrig) {
+ ptTrig = particle->Pt() ;
+ fLeadingTriggerIndex = iaod ;
+ }
+ }// finish search of leading trigger particle
- }//Aod branch loop
+ //Do correlation with leading particle
+ if(fLeadingTriggerIndex >= 0){
+
+ AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
+
+ //check if the particle is isolated or if we want to take the isolation into account
+ if(OnlyIsolated() && !particle->IsIsolated()) return;
+
+ //Make correlation with charged hadrons
+ Bool_t okcharged = kTRUE;
+ Bool_t okneutral = kTRUE;
+ if(GetReader()->IsCTSSwitchedOn() )
+ okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kFALSE);
+
+ TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
+ if(fNeutralCorr && pi0list && pi0list->GetEntriesFast() > 0)
+ okneutral = MakeNeutralCorrelation(particle, pi0list,kFALSE);
+
+ }//Correlate leading
if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
}
- //Loop on stored AOD particles
- Int_t naod = GetInputAODBranch()->GetEntriesFast();
- for(Int_t iaod = 0; iaod < naod ; iaod++){
- AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
-
+ //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
+ Double_t ptTrig = GetMinPt() ;
+ if(fLeadingTriggerIndex < 0){//Search leading if not done before
+ Int_t naod = GetInputAODBranch()->GetEntriesFast() ;
+ 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())
+ return ;
+ }
+ //check if the particle is isolated or if we want to take the isolation into account
+ if(OnlyIsolated() && !particle->IsIsolated()) continue;
+
+ //check if inside the vertex cut
+ //find the leading particles with highest momentum
+ if (particle->Pt() > ptTrig) {
+ ptTrig = particle->Pt() ;
+ fLeadingTriggerIndex = iaod ;
+ }
+ }//finish search of leading trigger particle
+ }//Search leading if not done before
+
+ if(fLeadingTriggerIndex >= 0 ){ //using trigger particle to do correlations
+ AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
+
//check if the particle is isolated or if we want to take the isolation into account
- if(OnlyIsolated() && !particle->IsIsolated()) continue;
+ if(OnlyIsolated() && !particle->IsIsolated()) return;
//Make correlation with charged hadrons
- TObjArray * reftracks = particle->GetObjArray(GetAODObjArrayName()+"Tracks");
- if(reftracks){
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Track Refs entries %d\n", iaod, reftracks->GetEntriesFast());
- if(reftracks->GetEntriesFast() > 0) MakeChargedCorrelation(particle, reftracks,kTRUE);
- }
+ Bool_t okcharged = kTRUE;
+ Bool_t okneutral = kTRUE;
+ if(GetReader()->IsCTSSwitchedOn() ){
+ okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kTRUE);
+ if(IsDataMC()){
+ MakeMCChargedCorrelation(particle);
+ }
+ }
- //Make correlation with neutral pions
- if(GetOutputAODBranch() && GetOutputAODBranch()->GetEntriesFast() > 0){
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Cluster Refs entries %d\n",iaod, GetOutputAODBranch()->GetEntriesFast());
- MakeNeutralCorrelationFillHistograms(particle);
+ TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
+ if(fNeutralCorr && pi0list){
+ if(pi0list->GetEntriesFast() > 0)
+ okneutral = MakeNeutralCorrelation(particle, pi0list,kTRUE);
}
+ // Fill leading particle histogram if correlation went well and
+ // no problem was found, like not absolute leading, or bad vertex in mixing.
+ if(okcharged && okneutral){
+ fhPtLeading->Fill(particle->Pt());
+ Float_t phi = particle->Phi();
+ if(phi<0)phi+=TMath::TwoPi();
+ fhPhiLeading->Fill(particle->Pt(), phi);
+ fhEtaLeading->Fill(particle->Pt(), particle->Eta());
+ }//ok charged && neutral
}//Aod branch loop
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
+ //Reinit for next event
+ fLeadingTriggerIndex = -1;
+ if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
}
//____________________________________________________________________________
-void AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle, TObjArray* const pl, const Bool_t bFillHisto)
+Bool_t AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle,
+ TObjArray* pl, const Bool_t bFillHisto)
{
- // Charged Hadron Correlation Analysis
- if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
-
- Double_t ptTrig = aodParticle->Pt();
- Double_t phiTrig = aodParticle->Phi();
- Double_t pt = -100.;
- Double_t rat = -100.;
- Double_t phi = -100. ;
- Double_t eta = -100. ;
- Double_t p[3];
-
- TObjArray * reftracks =0x0;
- if(!bFillHisto)
- reftracks = new TObjArray;
-
- //Track loop, select tracks with good pt, phi and fill AODs or histograms
- for(Int_t ipr = 0;ipr < pl->GetEntries() ; ipr ++ ){
- AliAODTrack * track = (AliAODTrack *) (pl->At(ipr)) ;
- track->GetPxPyPz(p) ;
- TLorentzVector mom(p[0],p[1],p[2],0);
- pt = mom.Pt();
- eta = mom.Eta();
- phi = mom.Phi() ;
- if(phi < 0) phi+=TMath::TwoPi();
- rat = pt/ptTrig ;
-
- if(IsFidutialCutOn()){
- Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,"CTS") ;
- if(! in ) continue ;
- }
-
- //Select only hadrons in pt range
- if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
-
- //Selection within angular range
- Float_t deltaphi = phiTrig-phi;
- if(deltaphi<-TMath::PiOver2()) deltaphi+=TMath::TwoPi();
- if(deltaphi>3*TMath::PiOver2()) deltaphi-=TMath::TwoPi();
-
- if(GetDebug() > 2)
- printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Charged hadron: pt %f, phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
- pt,phi, phiTrig,fDeltaPhiMinCut, deltaphi, fDeltaPhiMaxCut, GetMinPt());
-
- if(bFillHisto){
- // Fill Histograms
- fhEtaCharged->Fill(ptTrig,eta);
- fhPhiCharged->Fill(ptTrig,phi);
- fhDeltaEtaCharged->Fill(ptTrig,aodParticle->Eta()-eta);
- fhDeltaPhiCharged->Fill(ptTrig, deltaphi);
-
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
-
- //delta phi cut for correlation
- if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
- fhDeltaPhiChargedPt->Fill(pt,deltaphi);
- fhPtImbalanceCharged->Fill(ptTrig,rat);
- }
- else {
- fhDeltaPhiUeChargedPt->Fill(pt,deltaphi);
- fhPtImbalanceUeCharged->Fill(ptTrig,rat);
- }
- //several UE calculation
- if(fMakeSeveralUE){
- if((deltaphi<-fUeDeltaPhiMinCut) && (deltaphi >-fUeDeltaPhiMaxCut)){
- fhDeltaPhiUeLeftCharged->Fill(pt,deltaphi);
- fhPtImbalanceUeLeftCharged->Fill(ptTrig,rat);
- }
- if((deltaphi>fUeDeltaPhiMinCut) && (deltaphi <fUeDeltaPhiMaxCut)){
- fhDeltaPhiUeRightCharged->Fill(pt,deltaphi);
- fhPtImbalanceUeRightCharged->Fill(ptTrig,rat);
- }
- } //several UE calculation
-
- }
- else{
- reftracks->Add(track);
- }//aod particle loop
- }// track loop
-
- //Fill AOD with reference tracks, if not filling histograms
- if(!bFillHisto && reftracks->GetEntriesFast() > 0) {
- reftracks->SetName(GetAODObjArrayName()+"Tracks");
- aodParticle->AddObjArray(reftracks);
- }
+ // Charged Hadron Correlation Analysis
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
-}
-
-//____________________________________________________________________________
-void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD(AliAODPWG4ParticleCorrelation* const aodParticle,TObjArray* const pl, TString detector)
-{
- // Neutral Pion Correlation Analysis, find pi0, put them in new output aod, if correlation cuts passed
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Make trigger particle - neutral hadron correlation \n");
+ Int_t evtIndex11 = -1 ; //cluster trigger or pi0 trigger
+ Int_t evtIndex12 = -1 ; // pi0 trigger
+ Int_t evtIndex13 = -1 ; // charged trigger
+ Int_t indexPhoton1 = -1 ;
+ Int_t indexPhoton2 = -1 ;
- if(!NewOutputAOD()){
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Output aod not created, set AOD class name and branch name in the configuration file, STOP! \n");
- abort();
+ Double_t v[3] = {0,0,0}; //vertex ;
+ GetReader()->GetVertex(v);
+
+ if (GetMixedEvent()) {
+ evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
+ evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
+ evtIndex13 = GetMixedEvent()->EventIndex(aodParticle->GetTrackLabel(0)) ;
}
Double_t phiTrig = aodParticle->Phi();
- Int_t tag = 0;
- TLorentzVector gammai;
- TLorentzVector gammaj;
-
- //Get vertex for photon momentum calculation
- Double_t vertex[] = {0,0,0} ; //vertex
- Double_t vertex2[] = {0,0,0} ; //vertex of second input aod
- if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC)
- {
- GetReader()->GetVertex(vertex);
- if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex2);
+ Double_t etaTrig = aodParticle->Eta();
+ Double_t ptTrig = aodParticle->Pt();
+ Int_t trBin = -1;
+ for(Int_t i = 0 ; i < fNTriggerPtBins ; i++){
+ if(ptTrig > fTriggerPtBinLimit[i] && ptTrig < fTriggerPtBinLimit[i+1]) trBin= i;
}
-
- //Cluster loop, select pairs with good pt, phi and fill AODs or histograms
- //Int_t iEvent= GetReader()->GetEventNumber() ;
- Int_t nclus = pl->GetEntriesFast();
- for(Int_t iclus = 0;iclus < nclus ; iclus ++ ){
- AliAODCaloCluster * calo = (AliAODCaloCluster *) (pl->At(iclus)) ;
-
- //Input from second AOD?
- Int_t inputi = 0;
- if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= iclus) inputi = 1 ;
- else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= iclus) inputi = 1;
-
- //Cluster selection, not charged, with photon or pi0 id and in fidutial cut
- Int_t pdg=0;
- if (inputi == 0 && !SelectCluster(calo, vertex, gammai, pdg)) continue ;
- else if(inputi == 1 && !SelectCluster(calo, vertex2, gammai, pdg)) continue ;
-
- if(GetDebug() > 2)
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Neutral cluster in %s: pt %f, phi %f, phi trigger %f. Cuts: delta phi min %2.2f, max %2.2f, pT min %2.2f \n",
- detector.Data(), gammai.Pt(),gammai.Phi(),phiTrig,fDeltaPhiMinCut, fDeltaPhiMaxCut, GetMinPt());
-
- //2 gamma overlapped, found with PID
- if(pdg == AliCaloPID::kPi0){
-
- //Select only hadrons in pt range
- if(gammai.Pt() < GetMinPt() || gammai.Pt() > GetMaxPt()) continue ;
-
- //Selection within angular range
- Float_t phi = gammai.Phi();
- if(phi < 0) phi+=TMath::TwoPi();
- //Float_t deltaphi = TMath::Abs(phiTrig-phi);
- //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
-
- AliAODPWG4Particle pi0 = AliAODPWG4Particle(gammai);
- //pi0.SetLabel(calo->GetLabel(0));
- pi0.SetPdg(AliCaloPID::kPi0);
- pi0.SetDetector(detector);
-
- if(IsDataMC()){
- pi0.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(0),GetReader(),inputi));
- if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of candidate %d\n",pi0.GetTag());
- }//Work with stack also
- //Set the indeces of the original caloclusters
- pi0.SetCaloLabel(calo->GetID(),-1);
- AddAODParticle(pi0);
-
- if(GetDebug() > 2)
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Correlated with selected pi0 (pid): pt %f, phi %f\n",pi0.Pt(),pi0.Phi());
-
- }// pdg = 111
-
- //Make invariant mass analysis
- else if(pdg == AliCaloPID::kPhoton){
- //Search the photon companion in case it comes from a Pi0 decay
- //Apply several cuts to select the good pair;
- for(Int_t jclus = iclus+1; jclus < pl->GetEntries() ; jclus ++ ){
- AliAODCaloCluster * calo2 = (AliAODCaloCluster *) (pl->At(jclus)) ;
-
- //Input from second AOD?
- Int_t inputj = 0;
- if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= jclus) inputj = 1;
- else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= jclus) inputj = 1;
-
- //Cluster selection, not charged with photon or pi0 id and in fidutial cut
- Int_t pdgj=0;
- if (inputj == 0 && !SelectCluster(calo2, vertex, gammaj, pdgj)) continue ;
- else if(inputj == 1 && !SelectCluster(calo2, vertex2, gammaj, pdgj)) continue ;
-
- if(!SelectCluster(calo2,vertex, gammaj, pdgj)) continue ;
-
- if(pdgj == AliCaloPID::kPhoton ){
-
- if((gammai+gammaj).Pt() < GetMinPt() || (gammai+gammaj).Pt() > GetMaxPt()) continue ;
-
- //Selection within angular range
- Float_t phi = (gammai+gammaj).Phi();
- if(phi < 0) phi+=TMath::TwoPi();
- //Float_t deltaphi = TMath::Abs(phiTrig-phi);
- //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
-
- //Select good pair (aperture and invariant mass)
- if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
-
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Neutral Hadron Correlation: AOD Selected gamma pair: pt %2.2f, phi %2.2f, eta %2.2f, M %2.3f\n",
- (gammai+gammaj).Pt(),(gammai+gammaj).Phi(),(gammai+gammaj).Eta(), (gammai+gammaj).M());
-
- TLorentzVector pi0mom = gammai+gammaj;
- AliAODPWG4Particle pi0 = AliAODPWG4Particle(pi0mom);
- //pi0.SetLabel(calo->GetLabel(0));
- pi0.SetPdg(AliCaloPID::kPi0);
- pi0.SetDetector(detector);
- if(IsDataMC()){
- //Check origin of the candidates
-
- Int_t label1 = calo->GetLabel(0);
- Int_t label2 = calo2->GetLabel(0);
- Int_t tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader(), inputi);
- Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader(), inputj);
-
- if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
- if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) && GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay)){
-
- //Check if pi0 mother is the same
- if(GetReader()->ReadStack()){
- TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
- label1 = mother1->GetFirstMother();
- //mother1 = GetMCStack()->Particle(label1);//pi0
-
- TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
- label2 = mother2->GetFirstMother();
- //mother2 = GetMCStack()->Particle(label2);//pi0
- }
- else if(GetReader()->ReadAODMCParticles()){
- AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputi))->At(label1);//photon in kine tree
- label1 = mother1->GetMother();
- //mother1 = GetMCStack()->Particle(label1);//pi0
- AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputj))->At(label2);//photon in kine tree
- label2 = mother2->GetMother();
- //mother2 = GetMCStack()->Particle(label2);//pi0
- }
-
- //printf("mother1 %d, mother2 %d\n",label1,label2);
- if(label1 == label2)
- GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
- }
- }//Work with mc information also
- pi0.SetTag(tag);
- //Set the indeces of the original caloclusters
- pi0.SetCaloLabel(calo->GetID(), calo2->GetID());
- AddAODParticle(pi0);
-
-
- }//Pair selected
- }//if pair of gammas
- }//2nd loop
- }// if pdg = 22
- }//1st loop
-
- if(GetDebug() > 1)
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - End, %d pi0's found \n",GetOutputAODBranch()->GetEntriesFast());
-}
+ if(trBin==-1) return kFALSE;
+ if(OnlyIsolated() && aodParticle->IsIsolated()){
+ fhTrigPt[trBin]->Fill(ptTrig);
+ }
+
+ Double_t pt = -100. ;
+ Double_t px = -100. ;
+ Double_t py = -100. ;
+ Double_t rat = -100. ;
+ Double_t xE = -100. ;
+ Double_t cosi = -100. ;
+ Double_t phi = -100. ;
+ Double_t eta = -100. ;
+ Double_t pout = -100. ;
+
+ Double_t ptDecay1 = 0. ;
+ Double_t pxDecay1 = 0. ;
+ Double_t pyDecay1 = 0. ;
+ Double_t phiDecay1 = 0. ;
+ Double_t ptDecay2 = 0. ;
+ Double_t pxDecay2 = 0. ;
+ Double_t pyDecay2 = 0. ;
+ Double_t phiDecay2 = 0. ;
+
+ Double_t ratDecay1 = -100. ;
+ Double_t ratDecay2 = -100. ;
+ Double_t deltaPhi = -100. ;
+ Double_t deltaPhiOrg = -100. ;
+ Double_t deltaPhiDecay1 = -100. ;
+ Double_t deltaPhiDecay2 = -100. ;
+
+ TVector3 p3;
+ TLorentzVector photonMom ;
+ TObjArray * clusters = 0x0 ;
+ TObjArray * reftracks = 0x0;
+ Int_t nrefs = 0;
+ Int_t nTracks = GetCTSTracks()->GetEntriesFast() ;
+
+ if(fPi0Trigger){
+ indexPhoton1 = aodParticle->GetCaloLabel (0);
+ indexPhoton2 = aodParticle->GetCaloLabel (1);
+ if(GetDebug() > 1)printf("indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
+
+ if(indexPhoton1!=-1 && indexPhoton2!=-1){
+ if(aodParticle->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
+ else clusters = GetPHOSClusters() ;
+ for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++){
+ AliVCluster * photon = (AliVCluster*) (clusters->At(iclus));
+ photon->GetMomentum(photonMom,GetVertex(0)) ;
+ if(photon->GetID()==indexPhoton1) {
+ ptDecay1 = photonMom.Pt();
+ pxDecay1 = photonMom.Px();
+ pyDecay1 = photonMom.Py();
+ phiDecay1 = photonMom.Phi();
+ if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay1/ptTrig);
+ }
+ if(photon->GetID()==indexPhoton2) {
+ ptDecay2 = photonMom.Pt();
+ pxDecay2 = photonMom.Px();
+ pyDecay2 = photonMom.Py();
+ phiDecay2 = photonMom.Phi();
+ if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay2/ptTrig);
+ }
+ if(GetDebug() > 1)printf("Photon1 = %f, Photon2 = %f \n", ptDecay1, ptDecay2);
+ } //cluster loop
+ } //index of decay photons found
+ } //make decay-hadron correlation
+
+ //Track loop, select tracks with good pt, phi and fill AODs or histograms
+ for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ ){
+ AliVTrack * track = (AliVTrack *) (pl->At(ipr)) ;
+
+ Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
+ p3.SetXYZ(mom[0],mom[1],mom[2]);
+ pt = p3.Pt();
+ px = p3.Px();
+ py = p3.Py();
+ eta = p3.Eta();
+ phi = p3.Phi() ;
+ if(phi < 0) phi+=TMath::TwoPi();
+
+ //Select only hadrons in pt range
+ if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
+ //remove trigger itself for correlation when use charged triggers
+ if( track->GetID() == aodParticle->GetTrackLabel(0) || track->GetID() == aodParticle->GetTrackLabel(1) ||
+ track->GetID() == aodParticle->GetTrackLabel(2) || track->GetID() == aodParticle->GetTrackLabel(3) )
+ continue ;
+
+ if(IsFiducialCutOn()){
+ Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
+ if(! in ) continue ;
+ }
+
+ //jump out this event if near side associated particle pt larger than trigger
+ if (fMakeAbsoluteLeading){
+ if(pt > ptTrig && TMath::Abs(phi-phiTrig)<TMath::PiOver2()) return kFALSE;
+ }
+
+ //Only for mixed event
+ Int_t evtIndex2 = 0 ;
+ if (GetMixedEvent()) {
+ evtIndex2 = GetMixedEvent()->EventIndex(track->GetID()) ;
+ if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || evtIndex13 == evtIndex2 ) // photon and track from different events
+ continue ;
+ //vertex cut
+ if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut())
+ return kFALSE;
+ }
+
+ if(fPi0Trigger){
+ if(indexPhoton1!=-1 && indexPhoton2!=-1){
+ if(ptDecay1) ratDecay1 = pt/ptDecay1 ;
+ if(ptDecay2) ratDecay2 = pt/ptDecay2 ;
+ deltaPhiDecay1 = phiDecay1-phi;
+ deltaPhiDecay2 = phiDecay2-phi;
+ if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
+ if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
+ if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
+ if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
+ }
+ } //do decay-hadron correlation
+
+ //Selection within angular range
+ deltaPhi = phiTrig-phi;
+ deltaPhiOrg = deltaPhi;
+ if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
+ if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
+
+ pout = pt*TMath::Sin(deltaPhi) ;
+ rat = pt/ptTrig ;
+ xE =-pt/ptTrig*TMath::Cos(deltaPhi);
+ cosi = TMath::Log(1/xE);
+
+ if(GetDebug() > 2)
+ printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Charged hadron: pt %f, phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
+ pt,phi, phiTrig,fDeltaPhiMinCut, deltaPhi, fDeltaPhiMaxCut, GetMinPt());
+
+ // Fill Histograms
+ if(bFillHisto){
+
+ if(OnlyIsolated() && aodParticle->IsIsolated()){ //FILIP, if only wantisolated histograms
+ if(TMath::Cos(deltaPhi)<0){//away side
+ fhXE[trBin]->Fill(xE) ;
+ }
+
+ //Hardcoded values, BAD, FIXME
+ Double_t dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
+ if(TMath::Abs(dphiBrad)>0.325 && TMath::Abs(dphiBrad)<0.475){
+ fhAssocPtBkg[trBin]->Fill(pt);
+ }
+ if(dphiBrad<-1./3) dphiBrad += 2;
+ fhDphiTrigPtAssocPt[trBin]->Fill(dphiBrad);
+ }
+
+ fhEtaCharged->Fill(pt,eta);
+ fhPhiCharged->Fill(pt,phi);
+ fhDeltaEtaCharged->Fill(ptTrig,aodParticle->Eta()-eta);
+ fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
+ fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi,aodParticle->Eta()-eta);
+ fhPtAssocDeltaPhi->Fill(pt, deltaPhi);
+
+ 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() - GetMinMulti() )/GetMultiBin()*(im+1)){
+ fhTrigDeltaPhiCharged[im]->Fill(ptTrig,deltaPhi);
+ fhTrigDeltaEtaCharged[im]->Fill(ptTrig,aodParticle->Eta()-eta);
+ }
+ }
+ }
+ //delta phi cut for correlation
+ if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) ) {
+ fhDeltaPhiChargedPt->Fill(pt,deltaPhi);
+ fhPtImbalanceCharged->Fill(ptTrig,xE);
+ fhPtHbpCharged->Fill(ptTrig,cosi);
+ fhPtTrigPout->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() - GetMinMulti() )/GetMultiBin()*(im+1))
+ fhTrigCorr[im]->Fill(ptTrig,xE);
+ }
+ } //multiplicity events selection
+ } //delta phi cut for correlation
+ else if ((deltaPhi > fUeDeltaPhiMinCut) && ( deltaPhi < fUeDeltaPhiMaxCut)) { //UE study
+ fhDeltaPhiUeChargedPt->Fill(pt,deltaPhi);
+ Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
+ Double_t uexE = -(pt/ptTrig)*TMath::Cos(randomphi);
+ if(uexE < 0.) uexE = -uexE;
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - xe = %f, uexE = %f \n", xE, uexE);
+ fhPtImbalanceUeCharged->Fill(ptTrig,uexE);
+ fhPtHbpUeCharged->Fill(ptTrig,TMath::Log(1/uexE));
+ if(DoEventSelect()){
+ for(Int_t im=0; im<GetMultiBin(); im++){
+ if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
+ fhTrigUeCorr[im]->Fill(ptTrig,xE);
+ }
+ } //multiplicity events selection
+
+ } //UE study
+
+ if(fPi0Trigger){
+ if(indexPhoton1!=-1 && indexPhoton2!=-1){
+ fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
+ fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
+ if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
+ fhPtImbalanceDecayCharged->Fill(ptDecay1,ratDecay1);
+ if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
+ fhPtImbalanceDecayCharged->Fill(ptDecay2,ratDecay2);
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
+ } //index of decay photons found
+ } //make decay-hadron correlation
+
+ //several UE calculation
+ if(fMakeSeveralUE){
+ if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut)){
+ fhDeltaPhiUeLeftCharged->Fill(pt,deltaPhi);
+ fhPtImbalanceUeLeftCharged->Fill(ptTrig,rat);
+ fhPtHbpUeLeftCharged->Fill(ptTrig,cosi);
+ }
+ if((deltaPhi>fUeDeltaPhiMinCut) && (deltaPhi <fUeDeltaPhiMaxCut)){
+ fhDeltaPhiUeRightCharged->Fill(pt,deltaPhi);
+ fhPtImbalanceUeRightCharged->Fill(ptTrig,rat);
+ fhPtHbpUeRightCharged->Fill(ptTrig,cosi);
+
+ }
+ } //several UE calculation
+
+ //Fill leading particle histogram
+ fhPtLeading->Fill(ptTrig);
+ if(phiTrig<0)phiTrig+=TMath::TwoPi();
+ fhPhiLeading->Fill(ptTrig, phiTrig);
+ fhEtaLeading->Fill(ptTrig, etaTrig);
+
+ } //Fill histogram
+ else{
+ nrefs++;
+ if(nrefs==1){
+ reftracks = new TObjArray(0);
+ TString trackname = Form("%s+Tracks", GetAODObjArrayName().Data());
+ reftracks->SetName(trackname.Data());
+ reftracks->SetOwner(kFALSE);
+ }
+ reftracks->Add(track);
+ }//aod particle loop
+ }// track loop
+
+ //Fill AOD with reference tracks, if not filling histograms
+ if(!bFillHisto && reftracks) {
+ aodParticle->AddObjArray(reftracks);
+ }
+
+ return kTRUE;
+
+}
//____________________________________________________________________________
-void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms(AliAODPWG4ParticleCorrelation * const aodParticle)
+Bool_t AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliAODPWG4ParticleCorrelation * const aodParticle,
+ TObjArray* pi0list, const Bool_t bFillHisto)
{
// Neutral Pion Correlation Analysis
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistogramS() - Make trigger particle - pi0 correlation, %d pi0's \n",GetOutputAODBranch()->GetEntriesFast());
+ if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Make trigger particle - pi0 correlation, %d pi0's \n",pi0list->GetEntriesFast());
+
+ Int_t evtIndex11 = 0 ;
+ Int_t evtIndex12 = 0 ;
+ if (GetMixedEvent()) {
+ evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
+ evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
+ }
+
+ Double_t pt = -100. ;
+ Double_t px = -100. ;
+ Double_t py = -100. ;
+ Double_t rat = -100. ;
+ Double_t phi = -100. ;
+ Double_t eta = -100. ;
+ Double_t xE = -100. ;
+ Double_t cosi = -100. ;
- Double_t pt = -100.;
- Double_t rat = -100.;
- Double_t phi = -100.;
- Double_t eta = -100.;
Double_t ptTrig = aodParticle->Pt();
Double_t phiTrig = aodParticle->Phi();
Double_t etaTrig = aodParticle->Eta();
+ Double_t pxTrig = aodParticle->Px();
+ Double_t pyTrig = aodParticle->Py();
- if(!GetOutputAODBranch()){
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
- abort();
- }
+ Int_t indexPhoton1 =-1 ;
+ Int_t indexPhoton2 =-1 ;
+ Double_t ptDecay1 = 0. ;
+ Double_t pxDecay1 = 0. ;
+ Double_t pyDecay1 = 0. ;
+ Double_t phiDecay1 = 0. ;
+ Double_t ptDecay2 = 0. ;
+ Double_t pxDecay2 = 0. ;
+ Double_t pyDecay2 = 0. ;
+ Double_t phiDecay2 = 0. ;
+
+ Double_t ratDecay1 = -100. ;
+ Double_t ratDecay2 = -100. ;
+ Double_t deltaPhi = -100. ;
+ Double_t deltaPhiDecay1 = -100. ;
+ Double_t deltaPhiDecay2 = -100. ;
+
+ TObjArray * clusters = 0x0 ;
+ TLorentzVector photonMom ;
+
+ if(fPi0Trigger){
+ indexPhoton1 = aodParticle->GetCaloLabel (0);
+ indexPhoton2 = aodParticle->GetCaloLabel (1);
+ if(GetDebug() > 1)
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
+
+ if(indexPhoton1!=-1 && indexPhoton2!=-1){
+ if(aodParticle->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
+ else clusters = GetPHOSClusters() ;
+ for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++){
+ AliVCluster * photon = (AliVCluster*) (clusters->At(iclus));
+ photon->GetMomentum(photonMom,GetVertex(0)) ;
+ if(photon->GetID()==indexPhoton1) {
+ ptDecay1 = photonMom.Pt();
+ pxDecay1 = photonMom.Px();
+ pyDecay1 = photonMom.Py();
+ phiDecay1 = photonMom.Phi();
+ }
+ if(photon->GetID()==indexPhoton2) {
+ ptDecay2 = photonMom.Pt();
+ pxDecay2 = photonMom.Px();
+ pyDecay2 = photonMom.Py();
+ phiDecay2 = photonMom.Phi();
+ }
+ if(GetDebug() > 1)
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - Photon1 = %f, Photon2 = %f \n", ptDecay1, ptDecay2);
+ } //photonAOD loop
+ } //index of decay photons found
+ if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay1/ptTrig, ptDecay2/ptTrig);
+ } //make decay-hadron correlation
+
+ TObjArray * refpi0 =0x0;
+ Int_t nrefs = 0;
//Loop on stored AOD pi0
- Int_t naod = GetOutputAODBranch()->GetEntriesFast();
- if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - aod branch entries %d\n", naod);
+ Int_t naod = pi0list->GetEntriesFast();
+ if(GetDebug() > 0)
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - aod branch entries %d\n", naod);
for(Int_t iaod = 0; iaod < naod ; iaod++){
- AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
- Int_t pdg = pi0->GetPdg();
+ AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (pi0list->At(iaod));
- if(pdg != AliCaloPID::kPi0) continue;
- pt = pi0->Pt();
+ Int_t evtIndex2 = 0 ;
+ Int_t evtIndex3 = 0 ;
+ if (GetMixedEvent()) {
+ evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(0)) ;
+ evtIndex3 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(1)) ;
+
+ if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 ||
+ evtIndex11 == evtIndex3 || evtIndex12 == evtIndex3) // trigger and pi0 are not from different events
+ continue ;
+ }
+ //Int_t pdg = pi0->GetPdg();
+ //if(pdg != AliCaloPID::kPi0) continue;
+
+ pt = pi0->Pt();
+ px = pi0->Px();
+ py = pi0->Py();
if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
+ //jumped out this event if near side associated partile pt larger than trigger
+ if(pt > ptTrig && TMath::Abs(phi-phiTrig)<TMath::PiOver2()) break ;
//Selection within angular range
phi = pi0->Phi();
- //Float_t deltaphi = TMath::Abs(phiTrig-phi);
- //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
- Float_t deltaphi = phiTrig-phi;
- if(deltaphi<-TMath::PiOver2()) deltaphi+=TMath::TwoPi();
- if(deltaphi>3*TMath::PiOver2()) deltaphi-=TMath::TwoPi();
-
- rat = pt/ptTrig ;
- phi = pi0->Phi() ;
- eta = pi0->Eta() ;
-
- fhEtaNeutral->Fill(ptTrig,eta);
- fhPhiNeutral->Fill(ptTrig,phi);
- fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
- fhDeltaPhiNeutral->Fill(ptTrig,deltaphi);
-
- //delta phi cut for correlation
- if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
- fhDeltaPhiNeutralPt->Fill(pt,deltaphi);
- fhPtImbalanceNeutral->Fill(ptTrig,rat);
- }
+ //Float_t deltaPhi = TMath::Abs(phiTrig-phi);
+ //if( (deltaPhi < fDeltaPhiMinCut) || ( deltaPhi > fDeltaPhiMaxCut) ) continue ;
+
+ if(bFillHisto){
+
+ deltaPhi = phiTrig-phi;
+ if(deltaPhi<-TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
+ if(deltaPhi>3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
+
+ rat = pt/ptTrig ;
+ phi = pi0->Phi() ;
+ eta = pi0->Eta() ;
+ xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
+ if(xE <0.)xE =-xE;
+ cosi = TMath::Log(1/xE);
+
+ if(fPi0Trigger){
+ if(indexPhoton1!=-1 && indexPhoton2!=-1){
+ if(ptDecay1) ratDecay1 = pt/ptDecay1 ;
+ if(ptDecay2) ratDecay2 = pt/ptDecay2 ;
+ deltaPhiDecay1 = phiDecay1-phi;
+ deltaPhiDecay2 = phiDecay2-phi;
+ if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
+ if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
+ if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
+ if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
+ fhDeltaPhiDecayNeutral->Fill(ptDecay1, deltaPhiDecay1);
+ fhDeltaPhiDecayNeutral->Fill(ptDecay2, deltaPhiDecay2);
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
+ if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
+ fhPtImbalanceDecayNeutral->Fill(ptDecay1,ratDecay1);
+ if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
+ fhPtImbalanceDecayNeutral->Fill(ptDecay2,ratDecay2);
+ if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
+ }
+ } //do decay-hadron correlation
+
+ fhEtaNeutral->Fill(pt,eta);
+ fhPhiNeutral->Fill(pt,phi);
+ fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
+ fhDeltaPhiNeutral->Fill(ptTrig,deltaPhi);
+ fhDeltaPhiDeltaEtaNeutral->Fill(deltaPhi,etaTrig-eta);
+
+ //delta phi cut for correlation
+ if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) ) {
+ fhDeltaPhiNeutralPt->Fill(pt,deltaPhi);
+ fhPtImbalanceNeutral->Fill(ptTrig,rat);
+ fhPtHbpNeutral->Fill(ptTrig,cosi);
+ }
else {
- fhDeltaPhiUeNeutralPt->Fill(pt,deltaphi);
- fhPtImbalanceUeNeutral->Fill(ptTrig,rat);
+ fhDeltaPhiUeNeutralPt->Fill(pt,deltaPhi);
+ fhPtImbalanceUeNeutral->Fill(ptTrig,rat);
+ fhPtHbpUeNeutral->Fill(ptTrig,cosi);
}
- //several UE calculation
+ //several UE calculation
if(fMakeSeveralUE){
- if((deltaphi<-fUeDeltaPhiMinCut) && (deltaphi >-fUeDeltaPhiMaxCut)){
- fhDeltaPhiUeLeftNeutral->Fill(pt,deltaphi);
- fhPtImbalanceUeLeftNeutral->Fill(ptTrig,rat);
- }
- if((deltaphi>fUeDeltaPhiMinCut) && (deltaphi <fUeDeltaPhiMaxCut)){
- fhDeltaPhiUeRightNeutral->Fill(pt,deltaphi);
- fhPtImbalanceUeRightNeutral->Fill(ptTrig,rat);
- }
+ if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut)){
+ fhDeltaPhiUeLeftNeutral->Fill(pt,deltaPhi);
+ fhPtImbalanceUeLeftNeutral->Fill(ptTrig,rat);
+ fhPtHbpUeLeftNeutral->Fill(ptTrig,cosi);
+ }
+ if((deltaPhi>fUeDeltaPhiMinCut) && (deltaPhi <fUeDeltaPhiMaxCut)){
+ fhDeltaPhiUeRightNeutral->Fill(pt,deltaPhi);
+ fhPtImbalanceUeRightNeutral->Fill(ptTrig,rat);
+ fhPtHbpUeRightNeutral->Fill(ptTrig,cosi);
+ }
} //several UE calculation
-
-
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
+ }
+ else{
+ nrefs++;
+ if(nrefs==1){
+ refpi0 = new TObjArray(0);
+ refpi0->SetName(GetAODObjArrayName()+"Pi0s");
+ refpi0->SetOwner(kFALSE);
+ }
+ refpi0->Add(pi0);
+ }//put references in trigger AOD
+
+ if(GetDebug() > 2 )
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
}//loop
+
+ return kTRUE;
}
-
-
+
//____________________________________________________________________________
-Bool_t AliAnaParticleHadronCorrelation::SelectCluster(AliAODCaloCluster * calo, Double_t *vertex, TLorentzVector & mom, Int_t & pdg) const {
- //Select cluster depending on its pid and acceptance selections
-
- //Skip matched clusters with tracks
- if(calo->GetNTracksMatched() > 0) return kFALSE;
-
- TString detector = "";
- if(calo->IsPHOSCluster()) detector= "PHOS";
- else if(calo->IsEMCALCluster()) detector= "EMCAL";
-
- //Check PID
- calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line
- pdg = AliCaloPID::kPhoton;
- if(IsCaloPIDOn()){
- //Get most probable PID, 2 options check PID weights (in MC this option is mandatory)
- //or redo PID, recommended option for EMCal.
-
- if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC )
- pdg = GetCaloPID()->GetPdg(detector,calo->PID(),mom.E());//PID with weights
- else
- pdg = GetCaloPID()->GetPdg(detector,mom,calo);//PID recalculated
-
- if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
-
- //If it does not pass pid, skip
- if(pdg != AliCaloPID::kPhoton && pdg != AliCaloPID::kPi0) {
- return kFALSE ;
- }
- }//PID on
+void AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
+{
+ // Charged Hadron Correlation Analysis with MC information
+ if(GetDebug()>1)
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Make trigger particle - charged hadron correlation in AOD MC level\n");
+
+ AliStack * stack = 0x0 ;
+ TParticle * primary = 0x0 ;
+ TClonesArray * mcparticles0 = 0x0 ;
+ TClonesArray * mcparticles = 0x0 ;
+ AliAODMCParticle * aodprimary = 0x0 ;
+
+ Double_t eprim = 0 ;
+ Double_t ptprim = 0 ;
+ Double_t phiprim = 0 ;
+ Double_t etaprim = 0 ;
+ Double_t pxprim = 0 ;
+ Double_t pyprim = 0 ;
+ Double_t pzprim = 0 ;
+ Int_t nTracks = 0 ;
+ Int_t iParticle = 0 ;
+ Double_t charge = 0.;
+
+ Double_t mcrat =-100 ;
+ Double_t mcxE =-100 ;
+ Double_t mccosi =-100 ;
+
+ //Track loop, select tracks with good pt, phi and fill AODs or histograms
+ //Int_t currentIndex = -1 ;
+ Double_t mcTrackPt = 0 ;
+ Double_t mcTrackPhi = 0 ;
+ Double_t mcTrackEta = 0 ;
+ Double_t mcTrackPx = 0 ;
+ Double_t mcTrackPy = 0 ;
+ Double_t mcTrackPz = 0 ;
- //Check acceptance selection
- if(IsFidutialCutOn()){
- Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,detector) ;
- if(! in ) return kFALSE ;
+ if(GetReader()->ReadStack()){
+ nTracks = GetMCStack()->GetNtrack() ;
}
+ else nTracks = GetReader()->GetAODMCParticles()->GetEntriesFast() ;
+ //Int_t trackIndex[nTracks];
- if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
+ Int_t label= aodParticle->GetLabel();
+ if(label<0){
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** bad label ***: label %d \n", label);
+ return;
+ }
- return kTRUE;
+ if(GetReader()->ReadStack()){
+ stack = GetMCStack() ;
+ if(!stack) {
+ printf(" AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation- Stack not available, is the MC handler called? STOP\n");
+ abort();
+ }
+
+ nTracks=stack->GetNprimary();
+ if(label >= stack->GetNtrack()) {
+ if(GetDebug() > 2) printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
+ return ;
+ }
+ primary = stack->Particle(label);
+ if(!primary){
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no primary ***: label %d \n", label);
+ return;
+ }
+
+ eprim = primary->Energy();
+ ptprim = primary->Pt();
+ phiprim = primary->Phi();
+ etaprim = primary->Eta();
+ pxprim = primary->Px();
+ pyprim = primary->Py();
+ pzprim = primary->Pz();
+
+ if(primary){
+
+ for (iParticle = 0 ; iParticle < nTracks ; iParticle++) {
+ TParticle * particle = stack->Particle(iParticle);
+ TLorentzVector momentum;
+ //keep only final state particles
+ if(particle->GetStatusCode()!=1) continue ;
+ Int_t pdg = particle->GetPdgCode();
+ charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
+ particle->Momentum(momentum);
+
+ //---------- Charged particles ----------------------
+ if(charge != 0){
+ //Particles in CTS acceptance
+ Bool_t inCTS = GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
+ if(TMath::Abs(pdg) == 11 && stack->Particle(particle->GetFirstMother())->GetPdgCode()==22) continue ;
+ if(inCTS&&momentum.Pt() >GetMinPt())
+ {
+ mcTrackPt = particle->Pt();
+ mcTrackPhi = particle->Phi();
+ mcTrackEta = particle->Eta();
+ mcTrackPx = particle->Px();
+ mcTrackPy = particle->Py();
+ mcTrackPz = particle->Pz();
+ if(mcTrackPhi < 0) mcTrackPhi+=TMath::TwoPi();
+ //Select only hadrons in pt range
+ if(mcTrackPt < GetMinPt() || mcTrackPt > GetMaxPt()) continue ;
+ //remove trigger itself for correlation when use charged triggers
+ if(label==iParticle && mcTrackPt==ptprim && mcTrackPhi==phiprim && mcTrackEta==etaprim)
+ continue ;
+ //jumped out this event if near side associated partile pt larger than trigger
+ if( mcTrackPt> ptprim && TMath::Abs(mcTrackPhi-phiprim)<TMath::PiOver2())
+ return ;
+
+ mcrat = mcTrackPt/ptprim ;
+ mcxE = -(mcTrackPx*pxprim+mcTrackPy*pyprim)/(ptprim*ptprim);
+ if(mcxE <0.) mcxE =-mcxE;
+ mccosi = TMath::Log(1/mcxE);
+ // printf("rat = %f, xE = %f, cosi =%f \n", rat, xE, cosi);
+ // printf("phi = %f \n", phi);
+
+ //Selection within angular range
+ Double_t mcdeltaPhi = phiprim-mcTrackPhi;
+ if( mcdeltaPhi< -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
+ if( mcdeltaPhi>3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
+ Double_t mcpout = mcTrackPt*TMath::Sin(mcdeltaPhi) ;
+ if(GetDebug()>0 )
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
+ mcTrackPt,mcTrackPhi, phiprim,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt());
+ // Fill Histograms
+ fhMCEtaCharged->Fill(mcTrackPt,mcTrackEta);
+ fhMCPhiCharged->Fill(mcTrackPt,mcTrackPhi);
+ fhMCDeltaEtaCharged->Fill(ptprim,etaprim-mcTrackEta);
+ fhMCDeltaPhiCharged->Fill(ptprim,mcdeltaPhi);
+ fhMCPtAssocDeltaPhi->Fill(mcTrackPt, mcdeltaPhi);
+ // fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
+ fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,etaprim-mcTrackEta);
+
+ //delta phi cut for correlation
+ if( (mcdeltaPhi > fDeltaPhiMinCut) && ( mcdeltaPhi < fDeltaPhiMaxCut) ) {
+ fhMCDeltaPhiChargedPt->Fill(mcTrackPt,mcdeltaPhi);
+ fhMCPtImbalanceCharged->Fill(ptprim,mcxE);
+ fhMCPtHbpCharged->Fill(ptprim,mccosi);
+ fhMCPtTrigPout->Fill(ptprim, mcpout) ;
+ }//delta phi cut for correlation
+ } //tracks after cuts
+ }//Charged
+ } //track loop
+ } //when the leading particles could trace back to MC
+ } //ESD MC
+ else if(GetReader()->ReadAODMCParticles()){
+ //Get the list of MC particles
+ mcparticles0 = GetReader()->GetAODMCParticles(0);
+ if(!mcparticles0) return;
+ if(label >=mcparticles0->GetEntriesFast()){
+ if(GetDebug() > 2)
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***: label %d, n tracks %d \n", label,mcparticles0->GetEntriesFast());
+ return;
+ }
+ //Get the particle
+ aodprimary = (AliAODMCParticle*) mcparticles0->At(label);
+ if(!aodprimary) {
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no AOD primary ***: label %d \n", label);
+ return;
+ }
+
+ ptprim = aodprimary->Pt();
+ phiprim = aodprimary->Phi();
+ etaprim =aodprimary->Eta();
+ pxprim =aodprimary->Px();
+ pyprim =aodprimary->Py();
+ pzprim =aodprimary->Pz();
+ if(aodprimary){
+ mcparticles= GetReader()->GetAODMCParticles();
+ for (Int_t i=0; i<nTracks;i++) {
+ AliAODMCParticle *part = (AliAODMCParticle*) mcparticles->At(i);
+ if (!part->IsPhysicalPrimary()) continue;
+ Int_t pdg = part->GetPdgCode();
+ charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
+ TLorentzVector momentum(part->Px(),part->Py(),part->Pz(),part->E());
+ if(charge != 0){
+ if(part->Pt()> GetReader()->GetCTSPtMin()){
+ //Particles in CTS acceptance
+ Bool_t inCTS = GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
+ Int_t indexmother=part->GetMother();
+ if(indexmother>-1)
+ {
+ Int_t mPdg = ((AliAODMCParticle*) mcparticles->At(indexmother)) ->GetPdgCode();
+ if(TMath::Abs(pdg) == 11 && mPdg == 22) continue;
+ }
+ if(inCTS&&momentum.Pt() >GetMinPt())
+ {
+ mcTrackPt=part->Pt();
+ mcTrackPhi=part->Phi();
+ mcTrackEta=part->Eta();
+ mcTrackPx=part->Px();
+ mcTrackPy=part->Py();
+ mcTrackPz=part->Pz();
+ if(mcTrackPhi < 0) mcTrackPhi+=TMath::TwoPi();
+ //Select only hadrons in pt range
+ if(mcTrackPt < GetMinPt() || mcTrackPt > GetMaxPt()) continue ;
+ //remove trigger itself for correlation when use charged triggers
+ if(label==i && mcTrackPt==ptprim && mcTrackPhi==phiprim && mcTrackEta==etaprim)
+ continue ;
+ //jumped out this event if near side associated partile pt larger than trigger
+ if( mcTrackPt> ptprim && TMath::Abs(mcTrackPhi-phiprim)<TMath::PiOver2())
+ return ;
+
+ mcrat = mcTrackPt/ptprim ;
+ mcxE = -(mcTrackPx*pxprim+mcTrackPy*pyprim)/(ptprim*ptprim);
+ if(mcxE <0.)mcxE =-mcxE;
+ mccosi = TMath::Log(1/mcxE);
+
+ //Selection within angular range
+ Double_t mcdeltaPhi = phiprim-mcTrackPhi;
+ if( mcdeltaPhi< -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
+ if( mcdeltaPhi>3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
+ Double_t mcpout = mcTrackPt*TMath::Sin(mcdeltaPhi) ;
+ if(GetDebug()>0)
+ printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
+ mcTrackPt,mcTrackPhi, phiprim,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt());
+ // Fill Histograms
+ fhMCEtaCharged->Fill(mcTrackPt,mcTrackEta);
+ fhMCPhiCharged->Fill(mcTrackPt,mcTrackPhi);
+ fhMCDeltaEtaCharged->Fill(ptprim,etaprim-mcTrackEta);
+ fhMCDeltaPhiCharged->Fill(ptprim,mcdeltaPhi);
+ fhMCPtAssocDeltaPhi->Fill(mcTrackPt, mcdeltaPhi);
+ // fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
+ fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,etaprim-mcTrackEta);
+
+ //delta phi cut for correlation
+ if( (mcdeltaPhi > fDeltaPhiMinCut) && ( mcdeltaPhi < fDeltaPhiMaxCut) ) {
+ fhMCDeltaPhiChargedPt->Fill(mcTrackPt,mcdeltaPhi);
+ fhMCPtImbalanceCharged->Fill(ptprim,mcxE);
+ fhMCPtHbpCharged->Fill(ptprim,mccosi);
+ fhMCPtTrigPout->Fill(ptprim, mcpout) ;
+ }//delta phi cut for correlation
+
+ } //tracks after cuts
+
+ } //with minimum pt cut
+ } //only charged particles
+ } //MC particle loop
+ } //when the leading particles could trace back to MC
+ }// AOD MC
+}
+
+//____________________________________________________________________________
+void AliAnaParticleHadronCorrelation::SetNTriggerPtBins(Int_t n)
+{
+ // Set number of bins
+ if(n < 10 && n > 0)
+ {
+ fNTriggerPtBins = n ;
+ }
+ else
+ {
+ printf("n = larger than 9 or too small, set to 9 \n");
+ fNTriggerPtBins = 9;
+ }
}
+
+//____________________________________________________________________________
+void AliAnaParticleHadronCorrelation::SetTriggerPtBinLimit(Int_t ibin, Float_t pt)
+{
+ // Set the list of limits for the trigger pt bins
+
+ if(ibin <= fNTriggerPtBins || ibin >= 0)
+ {
+ fTriggerPtBinLimit[ibin] = pt ;
+ }
+ else {
+ printf("AliAnaParticleHadronCorrelation::SetTriggerPtBinLimit() - bin number too large %d > %d or small, nothing done\n", ibin, fNTriggerPtBins) ;
+
+ }
+}
+
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff