// --- ROOT system ---
#include "TH2F.h"
#include "TClonesArray.h"
+#include "TClass.h"
//---- ANALYSIS system ----
#include "AliNeutralMesonSelection.h"
#include "AliFidutialCut.h"
#include "AliAODTrack.h"
#include "AliAODCaloCluster.h"
+#include "AliMCAnalysisUtils.h"
+#include "TParticle.h"
+#include "AliStack.h"
ClassImp(AliAnaParticleHadronCorrelation)
fhDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
fhDeltaPhiChargedPt = new TH2F
- ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#p^{#pm}i} vs p_{T h^{#pm}}",
+ ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
nptbins,ptmin,ptmax,200,0,6.4);
fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
fhDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
SetAODRefArrayName("Hadrons");
AddToHistogramsName("AnaHadronCorr_");
+ //Correlation with neutrals
+ //SetOutputAODClassName("AliAODPWG4Particle");
+ //SetOutputAODName("Pi0Correlated");
+
SetPtCutRange(2,300);
fDeltaPhiMinCut = 1.5 ;
fDeltaPhiMaxCut = 4.5 ;
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, ABORT \n",GetInputAODName().Data());
abort();
}
+
+ 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();
+ }
+
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());
Int_t naod = GetInputAODBranch()->GetEntriesFast();
for(Int_t iaod = 0; iaod < naod ; iaod++){
AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
-
- //Make correlation with charged hadrons
+
+ //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)
- MakeNeutralCorrelation(particle, GetAODEMCAL(),kFALSE);
+ MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
//Trigger particle in EMCAL, correlation with PHOS
else if(particle->GetDetector()=="EMCAL" && GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
- MakeNeutralCorrelation(particle, GetAODPHOS(),kFALSE);
+ 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)
- MakeNeutralCorrelation(particle, GetAODPHOS(),kFALSE);
+ MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
if(GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
- MakeNeutralCorrelation(particle, GetAODEMCAL(),kFALSE);
+ MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
}
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s >, ABORT \n",GetInputAODName().Data());
abort();
}
-
+
if(GetDebug() > 1){
printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \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));
-
- //check if the particle is isolated or if we want to take the isolation into account
- if(OnlyIsolated() && !particle->IsIsolated()) continue;
-
- //Make correlation with charged hadrons
- TRefArray * reftracks = particle->GetRefArray(GetAODRefArrayName()+"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);
- }
-
- //Make correlation with neutral pions
- TRefArray * refclusters = particle->GetRefArray(GetAODRefArrayName()+"Clusters");
- if(refclusters && refclusters->GetEntriesFast() > 0){
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Cluster Refs entries %d\n",iaod, refclusters->GetEntriesFast());
- if(refclusters->GetEntriesFast() > 0) MakeNeutralCorrelation(particle,refclusters, kTRUE);
- }
-
+ for(Int_t iaod = 0; iaod < naod ; iaod++){
+ AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
+
+ //check if the particle is isolated or if we want to take the isolation into account
+ if(OnlyIsolated() && !particle->IsIsolated()) continue;
+
+ //Make correlation with charged hadrons
+ TRefArray * reftracks = particle->GetRefArray(GetAODRefArrayName()+"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);
+ }
+
+ //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);
+ }
+
}//Aod branch loop
if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
AliAODTrack * track = (AliAODTrack *) (pl->At(ipr)) ;
track->GetPxPyPz(p) ;
TLorentzVector mom(p[0],p[1],p[2],0);
- pt = mom.Pt();
+ pt = mom.Pt();
eta = mom.Eta();
phi = mom.Phi() ;
- if(phi<0) phi+=TMath::TwoPi();
+ if(phi < 0) phi+=TMath::TwoPi();
rat = pt/ptTrig ;
if(IsFidutialCutOn()){
//Select only hadrons in pt range
if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
+ //Selection within angular range
+ Float_t deltaphi = TMath::Abs(phiTrig-phi);
+ if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
+
if(GetDebug() > 2)
- printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Charged hadron: pt %f, phi %f, phi trigger %f. Cuts: delta phi min %2.2f, max%2.2f, pT min %2.2f \n",
- pt,phi,phiTrig,fDeltaPhiMinCut, fDeltaPhiMaxCut, GetMinPt());
-
+ 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,phiTrig-phi);
- fhDeltaPhiChargedPt->Fill(pt,phiTrig-phi);
- //Selection within angular range
- if(((phiTrig-phi)> fDeltaPhiMinCut) && ((phiTrig-phi)<fDeltaPhiMaxCut) ){
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
- fhPtImbalanceCharged->Fill(ptTrig,rat);
- }
+ fhDeltaPhiCharged->Fill(ptTrig, deltaphi);
+ fhDeltaPhiChargedPt->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);
+ fhPtImbalanceCharged->Fill(ptTrig,rat);
}
else{
//Fill AODs
new (reftracks) TRefArray(TProcessID::GetProcessWithUID(track));
first = kFALSE;
}
-
+
reftracks->Add(track);
}//aod particle loop
}// track loop
-
+
//Fill AOD with reference tracks, if not filling histograms
if(!bFillHisto && reftracks->GetEntriesFast() > 0) {
reftracks->SetName(GetAODRefArrayName()+"Tracks");
aodParticle->AddRefArray(reftracks);
}
-
+
}
+
//____________________________________________________________________________
-void AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliAODPWG4ParticleCorrelation * aodParticle,TRefArray* pl, const Bool_t bFillHisto)
+void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD(AliAODPWG4ParticleCorrelation * aodParticle,TRefArray* pl, TString detector)
{
- // Neutral Pion Correlation Analysis
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Make trigger particle - neutral hadron correlation \n");
+ // 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");
+
+ if(!NewOutputAOD()){
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Output aod not created, set AOD class name and branch name in the configuration file, ABORT! \n");
+ abort();
+ }
- 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();
- Bool_t first = kTRUE;
-
+ Int_t tag = -1;
TLorentzVector gammai;
TLorentzVector gammaj;
Double_t vertex[] = {0,0,0};
if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
- TRefArray * refclusters =0x0;
- if(!bFillHisto)
- refclusters = new TRefArray;
-
//Cluster loop, select pairs with good pt, phi and fill AODs or histograms
- for(Int_t iclus = 0;iclus < pl->GetEntries() ; iclus ++ ){
+ //Int_t iEvent= GetReader()->GetEventNumber() ;
+ Int_t nclus = pl->GetEntriesFast();
+ for(Int_t iclus = 0;iclus < nclus ; iclus ++ ){
AliAODCaloCluster * calo = (AliAODCaloCluster *) (pl->At(iclus)) ;
- if(!bFillHisto){
+
+ //Cluster selection, not charged, with photon or pi0 id and in fidutial cut
+ Int_t pdg=0;
+ if(!SelectCluster(calo, vertex, 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){
- //Cluster selection, not charged, with photon or pi0 id and in fidutial cut
- Int_t pdg=0;
- if(!SelectCluster(calo, vertex, gammai, pdg)) continue ;
+ //Select only hadrons in pt range
+ if(gammai.Pt() < GetMinPt() || gammai.Pt() > GetMaxPt()) continue ;
- if(GetDebug() > 2)
- printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Neutral cluster: pt %f, phi %f, phi trigger %f. Cuts: delta phi min %2.2f, max%2.2f, pT min %2.2f \n",
- gammai.Pt(),gammai.Phi(),phiTrig,fDeltaPhiMinCut, fDeltaPhiMaxCut, GetMinPt());
+ //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 ;
- //2 gamma overlapped, found with PID
- if(pdg == AliCaloPID::kPi0){
-
- //Select only hadrons in pt range
- if(gammai.Pt() < GetMinPt() || gammai.Pt() > GetMaxPt()) continue ;
-
- if(first) {
- new (refclusters) TRefArray(TProcessID::GetProcessWithUID(calo));
- first = kFALSE;
- }
+ 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),GetMCStack()));
+ 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)) ;
- refclusters->Add(calo);
- if(GetDebug() > 2) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Correlated with selected pi0 (pid): pt %f, phi %f\n",gammai.Pt(),gammai.Phi());
+ //Cluster selection, not charged with photon or pi0 id and in fidutial cut
+ Int_t pdgj=0;
+ if(!SelectCluster(calo2,vertex, gammaj, pdgj)) continue ;
- }// 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)) ;
+ if(pdgj == AliCaloPID::kPhoton ){
+
+ if((gammai+gammaj).Pt() < GetMinPt() || (gammai+gammaj).Pt() > GetMaxPt()) continue ;
- //Cluster selection, not charged with photon or pi0 id and in fidutial cut
- Int_t pdgj=0;
- if(!SelectCluster(calo2,vertex, gammaj, pdgj)) 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 ;
- if(pdgj == AliCaloPID::kPhoton ){
+ //Select good pair (aperture and invariant mass)
+ if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
- if((gammai+gammaj).Pt() < GetMinPt() || (gammai+gammaj).Pt() > GetMaxPt()) continue ;
+ 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());
- //Select good pair (aperture and invariant mass)
- if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
-
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - 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());
- Int_t labels[]={calo->GetLabel(0),calo2->GetLabel(0)};
- Float_t pid[]={0,0,0,0,0,0,1,0,0,0,0,0,0};//Pi0 weight 1
- Float_t pos[]={(gammai+gammaj).X(), (gammai+gammaj).Y(), (gammai+gammaj).Z()};
-
- AliAODCaloCluster *caloCluster = new AliAODCaloCluster(0,2,labels,(gammai+gammaj).E(), pos, pid,calo->GetType(),0);
- //Put this new object in file, need to think better how to do this!!!!
- GetReader()->GetAODEMCAL()->Add(caloCluster);
+ 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 tag1 = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(0), GetMCStack());
+ Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(calo2->GetLabel(0), GetMCStack());
- if(first) {
- new (refclusters) TRefArray(TProcessID::GetProcessWithUID(caloCluster));
- first = kFALSE;
+ if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
+ if(tag1 == AliMCAnalysisUtils::kMCPi0Decay && tag2 == AliMCAnalysisUtils::kMCPi0Decay){
+
+ //Check if pi0 mother is the same
+ Int_t label1 = calo->GetLabel(0);
+ TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
+ label1 = mother1->GetFirstMother();
+ //mother1 = GetMCStack()->Particle(label1);//pi0
+
+ Int_t label2 = calo2->GetLabel(0);
+ TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
+ label2 = mother2->GetFirstMother();
+ //mother2 = GetMCStack()->Particle(label2);//pi0
+
+ //printf("mother1 %d, mother2 %d\n",label1,label2);
+ if(label1 == label2)
+ tag = AliMCAnalysisUtils::kMCPi0;
}
-
- refclusters->Add(calo);
-
- }//Pair selected
- }//if pair of gammas
- }//2nd loop
- }// if pdg = 22
- }// Fill AODs
- else{ //Fill histograms
-
- calo->GetMomentum(gammai,vertex);//Assume that come from vertex in straight line
- pt = gammai.Pt();
-
- if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
-
- rat = pt/ptTrig ;
- phi = gammai.Phi() ;
- eta = gammai.Eta() ;
-
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Neutral Hadron Correlation: Histograms selected gamma pair: pt %2.2f, phi %2.2f, eta %2.2f\n",pt,phi,eta);
-
- fhEtaNeutral->Fill(ptTrig,eta);
- fhPhiNeutral->Fill(ptTrig,phi);
- fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
- fhDeltaPhiNeutral->Fill(ptTrig,phiTrig-phi);
- fhDeltaPhiNeutralPt->Fill(pt,phiTrig-phi);
- //Selection within angular range
- if(((phiTrig-phi)> fDeltaPhiMinCut) && ((phiTrig-phi)<fDeltaPhiMaxCut) ){
- if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
- fhPtImbalanceNeutral->Fill(ptTrig,rat);
- }
- }//Fill histograms
+ }//Work with stack 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
- //Fill AOD with reference tracks, if not filling histograms
- if(!bFillHisto && refclusters->GetEntriesFast() > 0) {
- refclusters->SetName(GetAODRefArrayName()+"Clusters");
- aodParticle->AddRefArray(refclusters);
+ if(GetDebug() > 1)
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - End, %d pi0's found \n",GetOutputAODBranch()->GetEntriesFast());
+}
+
+//____________________________________________________________________________
+void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms(AliAODPWG4ParticleCorrelation * aodParticle)
+{
+ // Neutral Pion Correlation Analysis
+ if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistogramS() - Make trigger particle - pi0 correlation, %d pi0's \n",GetOutputAODBranch()->GetEntriesFast());
+
+ 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();
+
+ if(!GetOutputAODBranch()){
+ printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
+ abort();
}
+
+ //Loop on stored AOD pi0
+ Int_t naod = GetOutputAODBranch()->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();
+
+ if(pdg != AliCaloPID::kPi0) continue;
+ pt = pi0->Pt();
+
+ if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
+
+ //Selection within angular range
+ phi = pi0->Phi();
+ Float_t deltaphi = TMath::Abs(phiTrig-phi);
+ if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
+
+ 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);
+ fhDeltaPhiNeutralPt->Fill(pt,deltaphi);
+ fhPtImbalanceNeutral->Fill(ptTrig,rat);
+
+ 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
}
+
//____________________________________________________________________________
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
+ //Select cluster depending on its pid and acceptance selections
+
+ //Skip matched clusters with tracks
if(calo->GetNTracksMatched() > 0) {
- return kFALSE;}
-
- //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.
- TString detector = "";
- if(calo->IsPHOSCluster()) detector= "PHOS";
- else if(calo->IsEMCALCluster()) detector= "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() > 1) 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 ;
- }
- }
-
- //Check acceptance selection
- if(IsFidutialCutOn()){
- Bool_t in = kFALSE;
- if(calo->IsPHOSCluster())
- in = GetFidutialCut()->IsInFidutialCut(mom,"PHOS") ;
- else if(calo->IsEMCALCluster())
- in = GetFidutialCut()->IsInFidutialCut(mom,"EMCAL") ;
- if(! in ) { return kFALSE ;}
- }
-
- if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::SelectCluster() - Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
-
- return kTRUE;
-
- }
+ return kFALSE;}
+
+ //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.
+ TString detector = "";
+ if(calo->IsPHOSCluster()) detector= "PHOS";
+ else if(calo->IsEMCALCluster()) detector= "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 ;
+ }
+ }
+
+ //Check acceptance selection
+ if(IsFidutialCutOn()){
+ Bool_t in = kFALSE;
+ if(calo->IsPHOSCluster())
+ in = GetFidutialCut()->IsInFidutialCut(mom,"PHOS") ;
+ else if(calo->IsEMCALCluster())
+ in = GetFidutialCut()->IsInFidutialCut(mom,"EMCAL") ;
+ if(! in ) { return kFALSE ;}
+ }
+
+ if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
+
+ return kTRUE;
+
+}
// --- ROOT system ---
#include "TH2F.h"
#include "TClonesArray.h"
+#include "TClass.h"
+//#include "Riostream.h"
//---- Analysis system ----
#include "AliAODTrack.h"
fhChargedLeadingRatioPt(0),
fhNeutralLeadingPt(0),fhNeutralLeadingPhi(0),fhNeutralLeadingEta(0),
fhNeutralLeadingDeltaPt(0),fhNeutralLeadingDeltaPhi(0),fhNeutralLeadingDeltaEta(0),
- fhNeutralLeadingRatioPt(0),
+ fhNeutralLeadingRatioPt(0),fhChargedLeadingXi(0), fhNeutralLeadingXi(0),
fhJetPt(0),fhJetRatioPt(0),fhJetDeltaPhi(0), fhJetDeltaEta(0),
fhJetLeadingRatioPt(0),fhJetLeadingDeltaPhi(0),fhJetLeadingDeltaEta(0),
fhJetFFz(0),fhJetFFxi(0),fhJetFFpt(0),fhJetNTracksInCone(0),
fhNeutralLeadingEta(jetlc.fhNeutralLeadingEta), fhNeutralLeadingDeltaPt(jetlc.fhNeutralLeadingDeltaPt),
fhNeutralLeadingDeltaPhi(jetlc.fhNeutralLeadingDeltaPhi),fhNeutralLeadingDeltaEta(jetlc.fhNeutralLeadingDeltaEta),
fhNeutralLeadingRatioPt(jetlc.fhNeutralLeadingRatioPt),
+ fhChargedLeadingXi(jetlc.fhChargedLeadingXi), fhNeutralLeadingXi(jetlc.fhNeutralLeadingXi),
fhJetPt(jetlc.fhJetPt),fhJetRatioPt(jetlc.fhJetRatioPt),fhJetDeltaPhi(jetlc.fhJetDeltaPhi),
fhJetDeltaEta(jetlc.fhJetDeltaEta), fhJetLeadingRatioPt(jetlc.fhJetLeadingRatioPt),
fhJetLeadingDeltaPhi(jetlc.fhJetLeadingDeltaPhi),fhJetLeadingDeltaEta(jetlc.fhJetLeadingDeltaEta),
fhNeutralLeadingEta = jetlc.fhNeutralLeadingEta; fhNeutralLeadingDeltaPt = jetlc.fhNeutralLeadingDeltaPt;
fhNeutralLeadingDeltaPhi = jetlc.fhNeutralLeadingDeltaPhi;fhNeutralLeadingDeltaEta = jetlc.fhNeutralLeadingDeltaEta;
fhNeutralLeadingRatioPt = jetlc.fhNeutralLeadingRatioPt;
+ fhChargedLeadingXi = jetlc.fhChargedLeadingXi;
+ fhNeutralLeadingXi = jetlc.fhNeutralLeadingXi;
+
+
fhJetPt = jetlc.fhJetPt;fhJetRatioPt = jetlc.fhJetRatioPt;fhJetDeltaPhi = jetlc.fhJetDeltaPhi;
fhJetDeltaEta = jetlc.fhJetDeltaEta; fhJetLeadingRatioPt = jetlc.fhJetLeadingRatioPt;
fhJetLeadingDeltaPhi = jetlc.fhJetLeadingDeltaPhi;fhJetLeadingDeltaEta = jetlc.fhJetLeadingDeltaEta;
Double_t AliAnaParticleJetLeadingConeCorrelation::CalculateJetRatioLimit(const Double_t ptg, const Double_t *par, const Double_t *x) const {
//Calculate the ratio of the jet and trigger particle limit for the selection
//WARNING: need to check what it does
- //printf("CalculateLimit: x1 %f, x2%f\n",x[0],x[1]);
+ //printf("CalculateLimit: x1 %2.3f, x2%2.3f\n",x[0],x[1]);
Double_t ePP = par[0] + par[1] * ptg ;
Double_t sPP = par[2] + par[3] * ptg ;
Double_t f = x[0] + x[1] * ptg ;
Double_t ePbPb = ePP + par[4] ;
Double_t sPbPb = TMath::Sqrt(sPP*sPP+ par[5]*par[5]) ;
Double_t rat = (ePbPb - sPbPb * f) / ptg ;
- //printf("CalculateLimit: ePP %f, sPP %f, f %f\n", ePP, sPP, f);
- //printf("CalculateLimit: ePbPb %f, sPbPb %f, rat %f\n", ePbPb, sPbPb, rat);
+ //printf("CalculateLimit: ePP %2.3f, sPP %2.3f, f %2.3f\n", ePP, sPP, f);
+ //printf("CalculateLimit: ePbPb %2.3f, sPbPb %2.3f, rat %2.3f\n", ePbPb, sPbPb, rat);
return rat ;
}
Double_t ptLead = leading.Pt();
Double_t phiTrig = particle->Phi();
Double_t phiJet = jet.Phi();
+ if(phiJet < 0) phiJet+=TMath::TwoPi();
Double_t phiLead = leading.Phi();
+ if(phiLead < 0) phiLead+=TMath::TwoPi();
Double_t etaTrig = particle->Eta();
Double_t etaJet = jet.Eta();
Double_t etaLead = leading.Eta();
-
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Pt"+lastname))->
- Fill(ptTrig,ptJet);
-
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"RatioPt"+lastname))->
- Fill(ptTrig,ptJet/ptTrig);
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingRatioPt"+lastname))->
- Fill(ptTrig,ptLead/ptJet);
-// dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Phi"+lastname))->
-// Fill(ptTrig,phiJet);
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaPhi"+lastname))->
- Fill(ptTrig,phiJet-phiTrig);
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaPhi"+lastname))->
- Fill(ptTrig,phiJet-phiLead);
-
- // dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Eta"+lastname))->
- // Fill(ptTrig,etaJet);
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaEta"+lastname))->
- Fill(ptTrig,etaJet-etaTrig);
- dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaEta"+lastname))->
- Fill(ptTrig,etaJet-etaLead);
+
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Pt"+lastname))->
+ Fill(ptTrig,ptJet);
+
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"RatioPt"+lastname))->
+ Fill(ptTrig,ptJet/ptTrig);
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingRatioPt"+lastname))->
+ Fill(ptTrig,ptLead/ptJet);
+ // dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Phi"+lastname))->
+ // Fill(ptTrig,phiJet);
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaPhi"+lastname))->
+ Fill(ptTrig,phiJet-phiTrig);
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaPhi"+lastname))->
+ Fill(ptTrig,phiJet-phiLead);
+
+ // dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"Eta"+lastname))->
+ // Fill(ptTrig,etaJet);
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"DeltaEta"+lastname))->
+ Fill(ptTrig,etaJet-etaTrig);
+ dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"LeadingDeltaEta"+lastname))->
+ Fill(ptTrig,etaJet-etaLead);
//Construct fragmentation function
TRefArray * pl = new TRefArray;
+
if(type == "Jet") pl = particle->GetRefArray(GetAODRefArrayName()+"Tracks");
else if(type == "Bkg") particle->GetRefArray(GetAODRefArrayName()+"TracksBkg");
-
+
+ if(!pl) return ;
+
//Different pt cut for jet particles in different collisions systems
//Only needed when jet is recalculated from AODs
Float_t ptcut = fJetPtThreshold;
TVector3 p3;
Int_t nTracksInCone = 0;
+
for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ ){
AliAODTrack* track = dynamic_cast<AliAODTrack *>(pl->At(ipr)) ;
p3.SetXYZ(track->Px(),track->Py(),track->Pz());
if(!IsParticleInJetCone(p3.Eta(), p3.Phi(), leading.Eta(), leading.Phi()) ) continue ;
nTracksInCone++;
-
+
dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"FFz"+lastname))
->Fill(ptTrig,p3.Pt()/ptTrig);
dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"FFxi"+lastname))
->Fill(ptTrig,p3.Pt());
}//track loop
-
+
if(nTracksInCone > 0) dynamic_cast<TH2F*>(GetOutputContainer()->FindObject(GetAddedHistogramsStringToName()+type+"NTracksInCone"+lastname))
->Fill(ptTrig, nTracksInCone);
{
// Create histograms to be saved in output file and
// store them in fOutCont
-
+
if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::GetCreateOutputObjects() - Init histograms \n");
fOutCont = new TList() ;
Float_t etamin = GetHistoEtaMin();
fhChargedLeadingPt = new TH2F("ChargedLeadingPt","p_{T leading charge} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
- fhChargedLeadingPt->SetYTitle("p_{T leading charge} /p_{T trigger}");
+ fhChargedLeadingPt->SetYTitle("p_{T leading charge}");
fhChargedLeadingPt->SetXTitle("p_{T trigger} (GeV/c)");
fhChargedLeadingPhi = new TH2F("ChargedLeadingPhi","#phi_{h^{#pm}} vs p_{T trigger}", nptbins,ptmin,ptmax,nphibins,phimin,phimax);
fhChargedLeadingEta->SetYTitle("#eta_{h^{#pm}} ");
fhChargedLeadingEta->SetXTitle("p_{T trigger} (GeV/c)");
- fhChargedLeadingDeltaPt = new TH2F("ChargedLeadingDeltaPt","#p_{T trigger} - #p_{T h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
+ fhChargedLeadingDeltaPt = new TH2F("ChargedLeadingDeltaPt","p_{T trigger} - p_{T h^{#pm}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
fhChargedLeadingDeltaPt->SetYTitle("#Delta p_{T} (GeV/c)");
fhChargedLeadingDeltaPt->SetXTitle("p_{T trigger} (GeV/c)");
fhChargedLeadingRatioPt->SetYTitle("p_{T lead charge} /p_{T trigger}");
fhChargedLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
+ fhChargedLeadingXi = new TH2F("ChargedLeadingXi","ln(p_{T trigger} / p_{T leading charge} ) vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,10);
+ fhChargedLeadingXi->SetYTitle("#xi");
+ fhChargedLeadingXi->SetXTitle("p_{T trigger} (GeV/c)");
+
fOutCont->Add(fhChargedLeadingPt) ;
fOutCont->Add(fhChargedLeadingPhi) ;
fOutCont->Add(fhChargedLeadingEta) ;
fOutCont->Add(fhChargedLeadingDeltaPhi) ;
fOutCont->Add(fhChargedLeadingDeltaEta) ;
fOutCont->Add(fhChargedLeadingRatioPt) ;
-
+ fOutCont->Add(fhChargedLeadingXi) ;
+
if(!fJetsOnlyInCTS){
fhNeutralLeadingPt = new TH2F("NeutralLeadingPt","p_{T leading #pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
- fhNeutralLeadingPt->SetYTitle("p_{T leading #pi^{0}} /p_{T trigger}");
+ fhNeutralLeadingPt->SetYTitle("p_{T leading #pi^{0}}");
fhNeutralLeadingPt->SetXTitle("p_{T trigger} (GeV/c)");
fhNeutralLeadingPhi = new TH2F("NeutralLeadingPhi","#phi_{#pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
fhNeutralLeadingEta->SetYTitle("#eta_{#pi^{0}} ");
fhNeutralLeadingEta->SetXTitle("p_{T trigger} (GeV/c)");
- fhNeutralLeadingDeltaPt = new TH2F("NeutralLeadingDeltaPt","#p_{T trigger} - #p_{T #pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
+ fhNeutralLeadingDeltaPt = new TH2F("NeutralLeadingDeltaPt","p_{T trigger} - p_{T #pi^{0}} vs p_{T trigger}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
fhNeutralLeadingDeltaPt->SetYTitle("#Delta p_{T} (GeV/c)");
fhNeutralLeadingDeltaPt->SetXTitle("p_{T trigger} (GeV/c)");
fhNeutralLeadingRatioPt->SetYTitle("p_{T lead #pi^{0}} /p_{T trigger}");
fhNeutralLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
+ fhNeutralLeadingXi = new TH2F("NeutralLeadingXi","ln(p_{T trigger} / p_{T leading #pi^{0}} ) vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,10);
+ fhNeutralLeadingXi->SetYTitle("#xi");
+ fhNeutralLeadingXi->SetXTitle("p_{T trigger} (GeV/c)");
+
fOutCont->Add(fhNeutralLeadingPt) ;
fOutCont->Add(fhNeutralLeadingPhi) ;
fOutCont->Add(fhNeutralLeadingEta) ;
fOutCont->Add(fhNeutralLeadingDeltaPhi) ;
fOutCont->Add(fhNeutralLeadingDeltaEta) ;
fOutCont->Add(fhNeutralLeadingRatioPt) ;
-
+ fOutCont->Add(fhNeutralLeadingXi) ;
}
if(!fSeveralConeAndPtCuts){// not several cones
fhJetLeadingRatioPt->SetYTitle("p_{T leading}/p_{T jet}");
fhJetLeadingRatioPt->SetXTitle("p_{T trigger} (GeV/c)");
- fhJetLeadingDeltaPhi = new TH2F("JetLeadingDeltaPhi","#phi_{jet} - #phi_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,0,TMath::TwoPi());
+ fhJetLeadingDeltaPhi = new TH2F("JetLeadingDeltaPhi","#phi_{jet} - #phi_{leading} vs p_{T trigger}",nptbins,ptmin,ptmax,120,-TMath::Pi(),TMath::Pi());
fhJetLeadingDeltaPhi->SetYTitle("#Delta #phi (rad)");
fhJetLeadingDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
} //icone
}//If we want to study any cone or pt threshold
+ //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++) fOutCont->Add(nmsHistos->At(i)) ;
+ }
+
+
if(GetDebug() > 2){
printf("AliAnaParticleJetLeadingConeCorrelation::GetCreateOutputObjects() - All histograms names : \n");
for(Int_t i = 0 ; i< fOutCont->GetEntries(); i++)
Double_t ptch = pLeadingCh.Pt();
Double_t ptpi = pLeadingPi0.Pt();
-
if (ptch > 0 || ptpi > 0){
if((ptch >= ptpi)){
if(GetDebug() > 1)printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Leading found in CTS \n");
pLeading = pLeadingCh;
- if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %f, phi %f deg, eta %f\n",
+ if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %2.3f, phi %2.3f deg, eta %2.3f\n",
pLeading.Pt(),pLeading.Phi()*TMath::RadToDeg(),pLeading.Eta()) ;
//Put leading in AOD
particle->SetLeading(pLeadingCh);
return kTRUE;
}
else{
- if(GetDebug() > 1)printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Leading found in EMCAL \n");
+ if(GetDebug() > 1)
+ printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Leading found in EMCAL \n");
pLeading = pLeadingPi0;
- if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %f, phi %f, eta %f\n",
+ if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingParticle() - Found Leading: pt %2.3f, phi %2.3f, eta %2.3f\n",
pLeading.Pt(),pLeading.Phi()*TMath::RadToDeg(),pLeading.Eta()) ;
//Put leading in AOD
particle->SetLeading(pLeadingPi0);
//Phi=Phi_trigger-Pi and pT=0.1E_gamma
if(GetAODCTS()){
- Double_t ptTrig = particle->Pt();
+ Double_t ptTrig = particle->Pt();
Double_t phiTrig = particle->Phi();
- Double_t rat = -100 ;
- Double_t ptl = -100 ;
- Double_t phil = -100 ;
- Double_t pt = -100.;
- Double_t phi = -100.;
+ Double_t rat = -100 ;
+ Double_t ptl = -100 ;
+ Double_t phil = -100 ;
+ Double_t pt = -100.;
+ Double_t phi = -100.;
TVector3 p3;
for(Int_t ipr = 0;ipr < GetAODCTS()->GetEntriesFast() ; ipr ++ ){
AliAODTrack* track = (AliAODTrack *)(GetAODCTS()->At(ipr)) ;
p3.SetXYZ(track->Px(),track->Py(),track->Pz());
- pt = p3.Pt();
+ pt = p3.Pt();
phi = p3.Phi() ;
- if(phi<0) phi+=TMath::TwoPi();
+ if(phi < 0) phi+=TMath::TwoPi();
rat = pt/ptTrig ;
-
- //Selection within angular and energy limits
- if(((phiTrig-phi) > fDeltaPhiMinCut) && ((phiTrig-phi)<fDeltaPhiMaxCut) &&
- (rat > fLeadingRatioMinCut) && (rat < fLeadingRatioMaxCut) && (pt > ptl)) {
+ //printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingCharge() - Tracks: pt %2.3f eta %2.3f phi %2.3f pt/ptTrig %2.3f \n",
+ // pt, p3.Eta(), phi,pt/ptTrig) ;
+ Float_t deltaphi = TMath::Abs(phiTrig-phi);
+ if((deltaphi > fDeltaPhiMinCut) && (deltaphi < fDeltaPhiMaxCut) &&
+ (rat > fLeadingRatioMinCut) && (rat < fLeadingRatioMaxCut) && (pt > ptl)) {
phil = phi ;
ptl = pt ;
pLeading.SetVect(p3);
}
}// track loop
- if(GetDebug() > 1&& ptl>0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingCharge() - Leading in CTS: pt %f eta %f phi %f pt/ptTrig %f \n",
- ptl, pLeading.Eta(), phil,ptl/ptTrig) ;
+ if(GetDebug() > 1 && ptl > 0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingCharge() - Leading in CTS: pt %2.3f eta %2.3f phi %2.3f pt/ptTrig %2.3f, |phiTrig-phi| %2.3f \n",
+ ptl, pLeading.Eta(), phil,ptl/ptTrig, TMath::Abs(phiTrig-phil)) ;
}//CTS list exist
}
//Phi=Phi_trigger-Pi and pT=0.1E_gamma
if(GetAODEMCAL()){
- Double_t ptTrig = particle->Pt();
+ Double_t ptTrig = particle->Pt();
Double_t phiTrig = particle->Phi();
- Double_t rat = -100 ;
- Double_t ptl = -100 ;
- Double_t phil = -100 ;
- Double_t pt = -100.;
- Double_t phi = -100.;
+ Double_t rat = -100 ;
+ Double_t ptl = -100 ;
+ Double_t phil = -100 ;
+ Double_t pt = -100.;
+ Double_t phi = -100.;
TLorentzVector gammai;
TLorentzVector gammaj;
Int_t pdgi=0;
if(!SelectCluster(calo,vertex, gammai, pdgi)) continue ;
- if(GetDebug() > 2) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral cluster: pt %f, phi %f \n",
+ if(GetDebug() > 2) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral cluster: pt %2.3f, phi %2.3f \n",
gammai.Pt(),gammai.Phi());
//2 gamma overlapped, found with PID
if(pdgi == AliCaloPID::kPi0){
- pt = gammai.Pt();
+
+ if(GetDebug() > 2) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral cluster ID as Pi0 \n");
+
+ pt = gammai.Pt();
rat = pt/ptTrig;
phi = gammai.Phi();
- if(phi<0) phi+=TMath::TwoPi();
+ if(phi < 0) phi+=TMath::TwoPi();
//Selection within angular and energy limits
- if(ptl > pt && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
- (phiTrig-phil) > fDeltaPhiMinCut && (phiTrig-phil) < fDeltaPhiMaxCut )
+ Float_t deltaphi = TMath::Abs(phiTrig-phi);
+ if(pt > ptl && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
+ deltaphi > fDeltaPhiMinCut && deltaphi < fDeltaPhiMaxCut )
{
- phi = phil ;
- pt = ptl ;
+ phil = phi ;
+ ptl = pt ;
pLeading.SetPxPyPzE(gammai.Px(),gammai.Py(),gammai.Pz(),gammai.E());
}// cuts
}// pdg = AliCaloPID::kPi0
pt = (gammai+gammaj).Pt();
phi = (gammai+gammaj).Phi();
- rat = pt/ptTrig;
-
- //Selection within angular and energy limits
- if(ptl > pt && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
- (phiTrig-phil) > fDeltaPhiMinCut && (phiTrig-phil) < fDeltaPhiMaxCut ){
- //Select good pair (aperture and invariant mass)
- if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
- phi = phil ;
- pt = ptl ;
- pLeading=(gammai+gammaj);
- }//pi0 selection
-
- if(GetDebug() > 3 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral Hadron Correlation: 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());
- }//Pair selected as leading
+ if(phi < 0) phi+=TMath::TwoPi();
+ rat = pt/ptTrig;
+
+ //Selection within angular and energy limits
+ Float_t deltaphi = TMath::Abs(phiTrig-phi);
+ if(GetDebug() > 3 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral Hadron Correlation: gamma pair: pt %2.2f, phi %2.2f, eta %2.2f, |phiTrig-phi| %2.3f, pt/ptTrig %2.3f, M %2.3f\n",
+ pt,phi,(gammai+gammaj).Eta(), deltaphi, rat, (gammai+gammaj).M());
+
+ if(pt > ptl && rat > fLeadingRatioMinCut && rat < fLeadingRatioMaxCut &&
+ deltaphi > fDeltaPhiMinCut && deltaphi < fDeltaPhiMaxCut ){
+ //Select good pair (aperture and invariant mass)
+ if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
+ phil = phi ;
+ ptl = pt ;
+ pLeading=(gammai+gammaj);
+
+ if(GetDebug() > 3 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Neutral Hadron Correlation: Selected gamma pair: pt %2.2f, phi %2.2f, eta %2.2f, M %2.3f\n",
+ ptl,phil,(gammai+gammaj).Eta(), (gammai+gammaj).M());
+ }//pi0 selection
+
+
+ }//Pair selected as leading
}//if pair of gammas
}//2nd loop
}// if pdg = 22
}// 1st Loop
- if(GetDebug()>2 && pLeading.Pt() >0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Leading EMCAL: pt %f eta %f phi %f pt/Eg %f \n", pLeading.Pt(), pLeading.Eta(), pLeading.Phi(), pLeading.Pt()/ptTrig) ;
+ if(GetDebug() > 2 && pLeading.Pt() > 0 ) printf("AliAnaParticleJetLeadingConeCorrelation::GetLeadingPi0() - Leading EMCAL: pt %2.3f eta %2.3f phi %2.3f pt/Eg %2.3f \n",
+ pLeading.Pt(), pLeading.Eta(), phil, pLeading.Pt()/ptTrig) ;
}//EMCAL list exists
-
}
//____________________________________________________________________________
SetInputAODName("PWG4Particle");
SetAODRefArrayName("JetLeadingCone");
AddToHistogramsName("AnaJetLCCorr_");
-
+
fJetsOnlyInCTS = kFALSE ;
fPbPb = kFALSE ;
fReMakeJet = kFALSE ;
fDeltaPhiMaxCut = 3.4 ;
fLeadingRatioMinCut = 0.1;
fLeadingRatioMaxCut = 1.5;
-
+
//Jet selection parameters
//Fixed cut
fJetRatioMaxCut = 1.2 ;
fJetCTSRatioMaxCut = 1.2 ;
fJetCTSRatioMinCut = 0.3 ;
fSelect = 0 ; //0, Accept all jets, 1, selection depends on energy, 2 fixed selection
-
+
fSelectIsolated = kFALSE;
-
+
//Cut depending on gamma energy
fPtTriggerSelectionCut = 10.; //For Low pt jets+BKG, another limits applied
//Reconstructed jet energy dependence parameters
//Given the pt of the jet and the trigger particle, select the jet or not
//3 options, fSelect=0 accepts all, fSelect=1 selects jets depending on a
//function energy dependent and fSelect=2 selects on simple fixed cuts
-
+
if(ptjet == 0) return kFALSE;
Double_t rat = ptTrig / ptjet ;
Double_t min = CalculateJetRatioLimit(ptTrig, par, xmin);
Double_t max = CalculateJetRatioLimit(ptTrig, par, xmax);
- if(GetDebug() > 3)printf("AliAnaParticleJetLeadingConeCorrelation::IsJetSelected() - Jet selection? : Limits min %f, max %f, pt_jet %f, pt_gamma %f, pt_jet / pt_gamma %f\n",min,max,ptjet,ptTrig,rat);
+ if(GetDebug() > 3)printf("AliAnaParticleJetLeadingConeCorrelation::IsJetSelected() - Jet selection? : Limits min %2.3f, max %2.3f, pt_jet %2.3f, pt_gamma %2.3f, pt_jet / pt_gamma %2.3f\n",min,max,ptjet,ptTrig,rat);
if(( min < rat ) && ( max > ptjet/rat))
return kTRUE;
//___________________________________________________________________
Bool_t AliAnaParticleJetLeadingConeCorrelation::IsParticleInJetCone(const Double_t eta, Double_t phi, const Double_t etal, Double_t phil)
- const {
+ const {
//Check if the particle is inside the cone defined by the leading particle
//WARNING: To be rechecked
-
+
if(phi < 0) phi+=TMath::TwoPi();
if(phil < 0) phil+=TMath::TwoPi();
Double_t rad = 10000 + fJetCone;
GetInputAODName().Data());
abort();
}
+
+ if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
+ printf("AliAnaParticleJetLeadingConeCorrelation::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("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - Begin jet leading cone correlation analysis, fill AODs \n");
printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
Int_t naod = GetInputAODBranch()->GetEntriesFast();
for(Int_t iaod = 0; iaod < naod ; iaod++){
AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
-
+
+ // printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillAOD() - Trigger : pt %3.2f, phi %2.2f, eta %2.2f\n",particle->Pt(), particle->Phi(), particle->Eta());
+
//Search leading particles in CTS and EMCAL
if(GetLeadingParticle(particle, pLeading)){
//Construct the jet around the leading, Fill AOD jet particle list, select jet
//and fill AOD with jet and background
MakeAODJet(particle, pLeading);
-
+
}//Leading found
}//AOD trigger particle loop
if(phiL < 0 ) phiL+=TMath::TwoPi();
Double_t etaL = pLeading.Eta();
- if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - Leading found in %s, with pt %3.2f, phi %2.2f, eta %2.2f\n",
- det.Data(), ptL, phiL, etaL);
+ if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAnalysisFillHistograms() - Trigger with pt %3.2f, phi %2.2f, eta %2.2f\n", pt, phi, eta);
+
if(det == "CTS"){
fhChargedLeadingPt->Fill(pt,ptL);
fhChargedLeadingPhi->Fill(pt,phiL);
fhChargedLeadingEta->Fill(pt,etaL);
fhChargedLeadingDeltaPt->Fill(pt,pt-ptL);
- fhChargedLeadingDeltaPhi->Fill(pt,phi-phiL);
+ fhChargedLeadingDeltaPhi->Fill(pt,TMath::Abs(phi-phiL));
fhChargedLeadingDeltaEta->Fill(pt,eta-etaL);
fhChargedLeadingRatioPt->Fill(pt,ptL/pt);
+ fhChargedLeadingXi->Fill(pt,TMath::Log(pt/ptL));
}
else if(det== "EMCAL"){
fhNeutralLeadingPt->Fill(pt,ptL);
fhNeutralLeadingPhi->Fill(pt,phiL);
fhNeutralLeadingEta->Fill(pt,etaL);
fhNeutralLeadingDeltaPt->Fill(pt,pt-ptL);
- fhNeutralLeadingDeltaPhi->Fill(pt,phi-phiL);
+ fhNeutralLeadingDeltaPhi->Fill(pt,TMath::Abs(phi-phiL));
fhNeutralLeadingDeltaEta->Fill(pt,eta-etaL);
fhNeutralLeadingRatioPt->Fill(pt,ptL/pt);
+ fhNeutralLeadingXi->Fill(pt,TMath::Log(pt/ptL));
}
//Fill Jet histograms
//____________________________________________________________________________
void AliAnaParticleJetLeadingConeCorrelation::MakeAODJet(AliAODPWG4ParticleCorrelation *particle, const TLorentzVector pLeading)
-const {
+ const {
//Fill the jet with the particles around the leading particle with
//R=fJetCone and pt_th = fJetPtThres. Calculate the energy of the jet and
//fill aod with found information
//Initialize reference arrays that will contain jet and background tracks
TRefArray * reftracks = new TRefArray;
TRefArray * reftracksbkg = new TRefArray;
-
+
for(Int_t ipr = 0;ipr < (GetAODCTS())->GetEntriesFast() ; ipr ++ ){
AliAODTrack* track = (AliAODTrack *)((GetAODCTS())->At(ipr)) ;
p3.SetXYZ(track->Px(),track->Py(),track->Pz());
}
reftracks->Add(track);
-
+
if(p3.Pt() > ptcut ){
lv.SetVect(p3);
jet+=lv;
}
}
+
//Background around (phi_gamma-pi, eta_leading)
else if(IsParticleInJetCone(p3.Eta(),p3.Phi(),etal, phiTrig)) {
-
+
if(firstbkg) {
new (reftracksbkg) TRefArray(TProcessID::GetProcessWithUID(track));
firstbkg = kFALSE;
//Add referenced tracks to AOD
if(reftracks->GetEntriesFast() > 0 ){
- reftracks->SetName(GetAODRefArrayName()+"Tracks");
- particle->AddRefArray(reftracks);
+ reftracks->SetName(GetAODRefArrayName()+"Tracks");
+ particle->AddRefArray(reftracks);
}
+ else if(GetDebug() > 2 ) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - No tracks in jet cone\n");
if(reftracksbkg->GetEntriesFast() > 0 ){
- reftracksbkg->SetName(GetAODRefArrayName()+"TracksBkg");
- particle->AddRefArray(reftracksbkg);
+ reftracksbkg->SetName(GetAODRefArrayName()+"TracksBkg");
+ particle->AddRefArray(reftracksbkg);
}
-
+ else if(GetDebug() > 2 ) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - No background tracks in jet cone\n");
+
//Add neutral particles to jet
//Initialize reference arrays that will contain jet and background tracks
TRefArray * refclusters = new TRefArray;
Double_t vertex[] = {0,0,0};
if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
-
+
first = kTRUE;
firstbkg = kTRUE;
//Add referenced clusters to AOD
if(refclusters->GetEntriesFast() > 0 ){
- refclusters->SetName(GetAODRefArrayName()+"Clusters");
- particle->AddRefArray(refclusters);
+ refclusters->SetName(GetAODRefArrayName()+"Clusters");
+ particle->AddRefArray(refclusters);
}
+ else if(GetDebug() > 2 ) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - No clusters in jet cone\n");
if(refclustersbkg->GetEntriesFast() > 0 ){
- refclustersbkg->SetName(GetAODRefArrayName()+"ClustersBkg");
- particle->AddRefArray(refclustersbkg);
+ refclustersbkg->SetName(GetAODRefArrayName()+"ClustersBkg");
+ particle->AddRefArray(refclustersbkg);
}
-
+ else if(GetDebug() > 2 ) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - No background clusters in jet cone\n");
+
//If there is any jet found, select after some criteria and
//and fill AOD with corresponding TLorentzVector kinematics
if(IsJetSelected(particle->Pt(), jet.Pt())) {
particle->SetCorrelatedJet(jet);
particle->SetCorrelatedBackground(bkg);
- if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - Found jet: Trigger pt %f, Jet pt %f, Bkg pt %f\n",ptTrig,jet.Pt(),bkg.Pt());
+ if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeAODJet() - Found jet: Trigger pt %2.3f, Jet pt %2.3f, Bkg pt %2.3f\n",ptTrig,jet.Pt(),bkg.Pt());
}
}
//____________________________________________________________________________
void AliAnaParticleJetLeadingConeCorrelation::MakeJetFromAOD(AliAODPWG4ParticleCorrelation *particle, const TLorentzVector pLeading, TLorentzVector & jet, TLorentzVector & bkg)
-const {
+ const {
//Fill the jet with the particles around the leading particle with
//R=fJetCone and pt_th = fJetPtThres. Calculate the energy of the jet and
//fill aod tlorentzvectors with jet and bakcground found
-
+
TLorentzVector lv (0,0,0,0); //Temporal container for jet particles kinematics
Double_t ptTrig = particle->Pt();
Double_t phiTrig = particle->Phi();
Double_t phil = pLeading.Phi();
+ if(phil < 0) phil+=TMath::TwoPi();
Double_t etal = pLeading.Eta();
TRefArray * refclusters = particle->GetRefArray(GetAODRefArrayName()+"Clusters");
TRefArray * reftracks = particle->GetRefArray(GetAODRefArrayName()+"Tracks");
TRefArray * refclustersbkg = particle->GetRefArray(GetAODRefArrayName()+"ClustersBkg");
TRefArray * reftracksbkg = particle->GetRefArray(GetAODRefArrayName()+"TracksBkg");
-
+
//Different pt cut for jet particles in different collisions systems
Float_t ptcut = fJetPtThreshold;
if(fPbPb && !fSeveralConeAndPtCuts && ptTrig > fPtTriggerSelectionCut) ptcut = fJetPtThresPbPb ;
for(Int_t ipr = 0;ipr < reftracks->GetEntriesFast() ; ipr ++ ){
AliAODTrack* track = (AliAODTrack *) reftracks->At(ipr) ;
p3.SetXYZ(track->Px(),track->Py(),track->Pz());
- if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(), p3.Phi(), etal, phil) ){
+ Float_t phi = p3.Phi();
+ if(phi < 0) phi+=TMath::TwoPi();
+ if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(), phi, etal, phil) ){
lv.SetVect(p3);
jet+=lv;
}
}
//Particles in background
if(reftracksbkg){
- for(Int_t ipr = 0;ipr < reftracksbkg->GetEntriesFast() ; ipr ++ ){
- AliAODTrack* track = (AliAODTrack *) reftracksbkg->At(ipr) ;
- p3.SetXYZ(track->Px(),track->Py(),track->Pz());
- if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(),p3.Phi(),etal, phiTrig) ) {
- lv.SetVect(p3);
- bkg+=lv;
- }
- }//background Track loop
+ for(Int_t ipr = 0;ipr < reftracksbkg->GetEntriesFast() ; ipr ++ ){
+ AliAODTrack* track = (AliAODTrack *) reftracksbkg->At(ipr) ;
+ p3.SetXYZ(track->Px(),track->Py(),track->Pz());
+ if(p3.Pt() > ptcut && IsParticleInJetCone(p3.Eta(),p3.Phi(),etal, phiTrig) ) {
+ lv.SetVect(p3);
+ bkg+=lv;
+ }
+ }//background Track loop
}
-
+
//Add neutral particles to jet
if(!fJetsOnlyInCTS && refclusters){
if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex);
//Loop on jet particles
- if(refclusters){
- for(Int_t iclus = 0;iclus < refclusters->GetEntriesFast() ; iclus ++ ){
- AliAODCaloCluster * calo = (AliAODCaloCluster *) refclusters->At(iclus) ;
- calo->GetMomentum(lv,vertex);
- if(lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(), etal, phil)) jet+=lv;
- }//jet cluster loop
+ if(refclusters){
+ for(Int_t iclus = 0;iclus < refclusters->GetEntriesFast() ; iclus ++ ){
+ AliAODCaloCluster * calo = (AliAODCaloCluster *) refclusters->At(iclus) ;
+ calo->GetMomentum(lv,vertex);
+ if(lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(), etal, phil)) jet+=lv;
+ }//jet cluster loop
+ }
+
+ //Loop on background particles
+ if(refclustersbkg){
+ for(Int_t iclus = 0;iclus < refclustersbkg->GetEntriesFast() ; iclus ++ ){
+ AliAODCaloCluster * calo = (AliAODCaloCluster *) refclustersbkg->At(iclus) ;
+ calo->GetMomentum(lv,vertex);
+ if( lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(),etal, phiTrig)) bkg+=lv;
+ }//background cluster loop
}
-
- //Loop on background particles
- if(refclustersbkg){
- for(Int_t iclus = 0;iclus < refclustersbkg->GetEntriesFast() ; iclus ++ ){
- AliAODCaloCluster * calo = (AliAODCaloCluster *) refclustersbkg->At(iclus) ;
- calo->GetMomentum(lv,vertex);
- if( lv.Pt() > ptcut && IsParticleInJetCone(lv.Eta(),lv.Phi(),etal, phiTrig)) bkg+=lv;
- }//background cluster loop
- }
}//clusters in jet
-
+
//If there is any jet found, leave jet and bkg as they are,
//if not set them to 0.
if(!IsJetSelected(particle->Pt(), jet.Pt())) {
bkg.SetPxPyPzE(0.,0.,0.,0.);
}
else
- if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeJetFromAOD() - Found jet: Trigger pt %f, Jet pt %f, Bkg pt %f\n",ptTrig,jet.Pt(),bkg.Pt());
+ if(GetDebug()>1) printf("AliAnaParticleJetLeadingConeCorrelation::MakeJetFromAOD() - Found jet: Trigger pt %2.3f, Jet pt %2.3f, Bkg pt %2.3f\n",ptTrig,jet.Pt(),bkg.Pt());
}
else
pdg = GetCaloPID()->GetPdg("EMCAL",mom,calo);//PID recalculated
- if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
+ // if(GetDebug() > 3) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
//If it does not pass pid, skip
if(pdg != AliCaloPID::kPhoton || pdg != AliCaloPID::kPi0)
return kFALSE ;
if(! in ) return kFALSE ;
}
- if(GetDebug() > 1) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - Cluster selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
+ //if(GetDebug() > 3) printf("AliAnaParticleJetLeadingConeCorrelation::SelectCluster() - Cluster selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
return kTRUE;
//__________________________________________________________________
void AliAnaParticleJetLeadingConeCorrelation::Print(const Option_t * opt) const
{
-
+
//Print some relevant parameters set for the analysis
if(! opt)
return;
printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
AliAnaPartCorrBaseClass::Print(" ");
-
+
if(fJetsOnlyInCTS)printf("Jets reconstructed in CTS \n");
else printf("Jets reconstructed in CTS+EMCAL \n");
-
+
if(fPbPb) printf("PbPb events, pT cut in jet cone energy reconstruction %2.1f \n", fJetPtThreshold);
else printf("pp events, pT cut in jet cone energy reconstruction %2.1f \n", fJetPtThresPbPb);
-
- printf("If pT of trigger < %f, select jets as in pp? \n", fPtTriggerSelectionCut);
-
+
+ printf("If pT of trigger < %2.3f, select jets as in pp? \n", fPtTriggerSelectionCut);
+
printf("Phi gamma-Leading < %3.2f\n", fDeltaPhiMaxCut) ;
printf("Phi gamma-Leading > %3.2f\n", fDeltaPhiMinCut) ;
printf("pT Leading / pT Trigger < %3.2f\n", fLeadingRatioMaxCut) ;
printf("Accept jets depending on trigger energy \n") ;
else
printf("Wrong jet selection option: %d \n", fSelect) ;
-
+
printf("Isolated Trigger? %d\n", fSelectIsolated) ;
}