fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fRatioMaxCut(0.), fRatioMinCut(0.),
fConeSize(0.), fPtThresholdInCone(0.),fUseJetRefTracks(kTRUE),
fMakeCorrelationInHistoMaker(kFALSE), fSelectIsolated(kTRUE),
- fJetConeSize(0.4),fJetMinPt(5),fJetAreaFraction(0.6),
+ fJetConeSize(0.4),fJetMinPt(5),fJetMinPtBkgSub(-100.),fJetAreaFraction(0.6),
//fNonStandardJetFromReader(kTRUE),
fJetBranchName("jets"),
fBackgroundJetFromReader(kTRUE),
fhDeltaEta(0), /*fhDeltaPhi(0),*/fhDeltaPhiCorrect(0),fhDeltaPhi0PiCorrect(0), fhDeltaPt(0), fhPtRatio(0), fhPt(0),
fhFFz(0),fhFFxi(0),fhFFpt(0),fhNTracksInCone(0),
fhJetFFz(0),fhJetFFxi(0),fhJetFFpt(0),fhJetFFzCor(0),fhJetFFxiCor(0),
- fhBkgFFz(),fhBkgFFxi(),fhBkgFFpt(),fhBkgNTracksInCone(),fhBkgSumPtInCone(),fhBkgSumPtnTracksInCone(),//<<---new
+ fhGamPtPerTrig(0),fhPtGamPtJet(0),
+ fhBkgFFz(),fhBkgFFxi(),fhBkgFFpt(),fhBkgNTracksInCone(),fhBkgSumPtInCone(),fhBkgSumPtnTracksInCone(),
fhNjetsNgammas(0),fhCuts(0),
fhDeltaEtaBefore(0),fhDeltaPhiBefore(0),fhDeltaPtBefore(0),fhPtRatioBefore(0),
fhPtBefore(0),fhDeltaPhi0PiCorrectBefore(0),
- fhJetPtBefore(0),fhJetPt(0),fhJetPtMostEne(0),fhJetPhi(0),fhJetEta(0),fhJetEtaVsPt(0),
+ fhJetPtBefore(0),fhJetPtBeforeCut(0),fhJetPt(0),fhJetPtMostEne(0),fhJetPhi(0),fhJetEta(0),fhJetEtaVsPt(0),
fhJetPhiVsEta(0),fhJetEtaVsNpartInJet(0),fhJetEtaVsNpartInJetBkg(0),fhJetChBkgEnergyVsPt(0),fhJetChAreaVsPt(0),/*fhJetNjet(0),*/
fhTrackPhiVsEta(0),fhTrackAveTrackPt(0),fhJetNjetOverPtCut(),
/*fhJetChBkgEnergyVsPtEtaGt05(0),fhJetChBkgEnergyVsPtEtaLe05(0),fhJetChAreaVsPtEtaGt05(0),fhJetChAreaVsPtEtaLe05(0),*/
fhJetFFxiCor->SetXTitle("p_{T jet}");
outputContainer->Add(fhJetFFxiCor) ;
+ fhGamPtPerTrig = new TH1F("GamPtPerTrig","GamPtPerTrig", nptbins,ptmin,ptmax);
+ fhGamPtPerTrig->SetYTitle("Counts");
+ fhGamPtPerTrig->SetXTitle("p_{T, #gamma}");
+ outputContainer->Add(fhGamPtPerTrig) ;
+
+ fhPtGamPtJet = new TH2F("PtGamPtJet","p_{T #gamma} vs p_{T jet}", nptbins,ptmin,ptmax,150,-50.,100.);
+ fhPtGamPtJet->SetXTitle("p_{T #gamma}");
+ fhPtGamPtJet->SetYTitle("p_{T jet}");
+ outputContainer->Add(fhPtGamPtJet) ;
+
//background FF
fhBkgFFz[0] = new TH2F("BkgFFzRC", "z = p_{T i charged}/p_{T trigger} vs p_{T trigger} Bkg RC" ,nptbins,ptmin,ptmax,200,0.,2);
fhJetPtBefore->SetXTitle("p_{T jet}(GeV/c)");
outputContainer->Add(fhJetPtBefore) ;
+ fhJetPtBeforeCut = new TH1F("JetPtBeforeCut","JetPtBeforeCut",150,-50,100);
+ fhJetPtBeforeCut->SetYTitle("Counts");
+ fhJetPtBeforeCut->SetXTitle("p_{T jet}(GeV/c)");
+ outputContainer->Add(fhJetPtBeforeCut) ;
+
fhJetPt = new TH1F("JetPt","JetPt",150,-50,100);
fhJetPt->SetYTitle("Counts");
fhJetPt->SetXTitle("p_{T jet}(GeV/c)");
fSelectIsolated = kTRUE;
fJetConeSize = 0.4 ;
fJetMinPt = 15. ; //GeV/c
+ fJetMinPtBkgSub = -100. ;//GeV/c
fJetAreaFraction = 0.6 ;
fJetBranchName = "jets";
fBkgJetBranchName = "jets";
Double_t particlePt=particle->Pt();
if(fUseBackgroundSubtractionGamma) {
- Int_t clusterIDtmp = particle->GetCaloLabel(0) ;
- Int_t nCells=0;
- AliVCluster *cluster=0;
- if(!(clusterIDtmp<0) ){
- Int_t iclustmp = -1;
- TObjArray* clusters = GetEMCALClusters();
- cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
- nCells = cluster->GetNCells();
- }
- particlePt-=(fGamRho*nCells);
+ particlePt-=(fGamRho*particle->GetNCells());
+
+// Int_t clusterIDtmp = particle->GetCaloLabel(0) ;
+// Int_t nCells=0;
+// AliVCluster *cluster=0;
+// if(!(clusterIDtmp<0) ){
+// Int_t iclustmp = -1;
+// TObjArray* clusters = GetEMCALClusters();
+// cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
+// nCells = cluster->GetNCells();
+// }
+// particlePt-=(fGamRho*nCells);
}
if(particlePt<=0) {
//printf("Particle with negative or 0 pt\n");
Double_t deltaPhi=-10000.;// in the range (0; 2*pi)
Double_t jetPt=0.;
+ Bool_t photonOnlyOnce=kTRUE;
+
for(Int_t ijet = 0; ijet < njets ; ijet++){
jet = dynamic_cast<AliAODJet*>(aodRecJets->At(ijet));
if(!jet)
jetPt=jet->Pt();
if(jetPt<fJetMinPt) continue;
fhCuts2->Fill(3.,1.);
- if(fBackgroundJetFromReader ){
- jetPt-= (fJetRho * jet->EffectiveAreaCharged() );
- }
- if(jetPt<0.) continue;
//put jet eta requirement here |eta_jet|<0.9-jet_cone_size
- fhCuts2->Fill(4.,1.);
if(TMath::Abs(jet->Eta()) > (0.9 - fJetConeSize) ) continue;
- fhCuts2->Fill(5.,1.);
+ fhCuts2->Fill(4.,1.);
if(jet->EffectiveAreaCharged()<fJetAreaFraction*TMath::Pi()*fJetConeSize*fJetConeSize) continue;
+ fhCuts2->Fill(5.,1.);
+ if(fBackgroundJetFromReader ){
+ jetPt-= (fJetRho * jet->EffectiveAreaCharged() );
+ }
+
+ if(jetPt<fJetMinPtBkgSub) continue;
fhCuts2->Fill(6.,1.);
//printf("jet found\n");
Double_t deltaPhi0pi = TMath::Abs(particle->Phi()-jet->Phi());
if ( deltaPhi0pi > TMath::Pi() ) deltaPhi0pi = 2. * TMath::Pi() - deltaPhi0pi ;
if(deltaPhi<0) deltaPhi +=(TMath::Pi()*2.);
+ //new histogram for Leticia x-check
+ //isolated photon + jet(s)
+ if(OnlyIsolated() && !particle->IsIsolated() &&
+ (deltaPhi > fDeltaPhiMinCut) && (deltaPhi < fDeltaPhiMaxCut) ){
+ //fill 1D photon + 2D photon+jets
+ if(photonOnlyOnce) {
+ fhGamPtPerTrig->Fill(particlePt);
+ photonOnlyOnce=kFALSE;
+ }
+ fhPtGamPtJet->Fill(particlePt,jetPt);
+ }
+
+
fhDeltaPtBefore ->Fill(particlePt, particlePt - jetPt);
fhDeltaPhiBefore->Fill(particlePt, deltaPhi);
fhDeltaEtaBefore->Fill(particlePt, particle->Eta() - jet->Eta());
//calculate average cell energy without most energetic photon
//
Double_t medianPhotonRho=0.;
- TObjArray* clusters = GetEMCALClusters();
- Int_t clusterIDtmp;
- Int_t iclustmp = -1;
- AliVCluster *cluster=0;
+ //TObjArray* clusters = GetEMCALClusters();
+ //Int_t clusterIDtmp;
+ //Int_t iclustmp = -1;
+ //AliVCluster *cluster=0;
if(IsBackgroundSubtractionGamma()){
//
for(Int_t iaod = 0; iaod < ntrig ; iaod++){
particlecorr = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
if(iaod==maxIndex) continue;
- clusterIDtmp = particlecorr->GetCaloLabel(0) ;
- if(clusterIDtmp < 0) continue;
- cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
- photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ cluster->GetNCells();
- numberOfcells+=cluster->GetNCells();
+// clusterIDtmp = particlecorr->GetCaloLabel(0) ;
+// if(clusterIDtmp < 0) continue;
+// cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
+// photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ cluster->GetNCells();
+// numberOfcells+=cluster->GetNCells();
+ photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ particlecorr->GetNCells();
+ numberOfcells+=particlecorr->GetNCells();
+
photonRhoArrayIndex++;
}
if(photonRhoArrayIndex>0) medianPhotonRho=TMath::Median(photonRhoArrayIndex,photonRhoArr);
Int_t indexMostEnePhoton=-1;
AliAODPWG4ParticleCorrelation* particle =0;
Double_t ptCorrect=0.;
- Int_t nCells=0;
+// Int_t nCells=0;
for(Int_t iaod = 0; iaod < ntrig ; iaod++){
particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
- clusterIDtmp = particle->GetCaloLabel(0) ;
- if(!(clusterIDtmp<0)){
- cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
- nCells = cluster->GetNCells();
- }
- ptCorrect = particle->Pt() - medianPhotonRho * nCells;
+// clusterIDtmp = particle->GetCaloLabel(0) ;
+// if(!(clusterIDtmp<0)){
+// cluster = FindCluster(clusters,clusterIDtmp,iclustmp);
+// nCells = cluster->GetNCells();
+// }
+// ptCorrect = particle->Pt() - medianPhotonRho * nCells;
+ ptCorrect = particle->Pt() - medianPhotonRho * particle->GetNCells();
+
if( ptCorrect > mostEnePhotonPt ){
mostEnePhotonPt = ptCorrect;
indexMostEnePhoton = iaod ;
ptMostEne = jetPttmp;
//indexMostEne=ijet;
}
+ if(jettmp->Pt()>=fJetMinPt)
+ fhJetPtBeforeCut->Fill(jetPttmp);
+
fhJetPt->Fill(jetPttmp);
fhJetChBkgEnergyVsPt->Fill(jetPttmp,effectiveChargedBgEnergy);
fhJetChAreaVsPt->Fill(jetPttmp,jettmp->EffectiveAreaCharged());
}
fGamAvEne=0;
- TObjArray* clusters = GetEMCALClusters();
+ //TObjArray* clusters = GetEMCALClusters();
//printf("calculate median bkg energy for photons ");
Double_t medianPhotonRho=0.;
- Int_t clusterID;
- Int_t iclustmp = -1;
+ //Int_t clusterID;
+ //Int_t iclustmp = -1;
Int_t numberOfcells=0;
- AliVCluster *cluster = 0;
+ //AliVCluster *cluster = 0;
if(ntrig>1){
Double_t *photonRhoArr=new Double_t[ntrig-1];
fhPhotonPtMostEne->Fill(maxPt);
if( particlecorr->Pt() > (sumPt-maxPt)/(ntrig-1) ) counterGammaMinus1++;
if(iaod==maxIndex) continue;
- clusterID = particlecorr->GetCaloLabel(0) ;
- if(clusterID < 0) continue;
- cluster = FindCluster(clusters,clusterID,iclustmp);
- photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ cluster->GetNCells();
- numberOfcells+=cluster->GetNCells();
+// clusterID = particlecorr->GetCaloLabel(0) ;
+// if(clusterID < 0) continue;
+// cluster = FindCluster(clusters,clusterID,iclustmp);
+// photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ cluster->GetNCells();
+// numberOfcells+=cluster->GetNCells();
+ photonRhoArr[photonRhoArrayIndex]=particlecorr->Pt()/ particlecorr->GetNCells();
+ numberOfcells+=particlecorr->GetNCells();
photonRhoArrayIndex++;
}
if(photonRhoArrayIndex>0) medianPhotonRho=TMath::Median(photonRhoArrayIndex,photonRhoArr);
fhPhotonBkgRhoVsNcells->Fill(numberOfcells,medianPhotonRho);
- AliVCluster *cluster2 = 0;
+ //AliVCluster *cluster2 = 0;
Double_t photon2Corrected=0;
Double_t sumPtTmp=0.;
Double_t sumPtCorrectTmp=0.;
for(Int_t iaod = 0; iaod < ntrig ; iaod++){
AliAODPWG4ParticleCorrelation* particlecorr = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
- clusterID = particlecorr->GetCaloLabel(0) ;
- if(clusterID < 0) continue;
- cluster = FindCluster(clusters,clusterID,iclustmp);
+// clusterID = particlecorr->GetCaloLabel(0) ;
+// if(clusterID < 0) continue;
+// cluster = FindCluster(clusters,clusterID,iclustmp);
+// Int_t ncells = cluster->GetNCells();
+ Int_t ncells = particlecorr->GetNCells();
fhPhotonPt->Fill(particlecorr->Pt());
- fhPhotonPtCorrected->Fill(particlecorr->Pt() - cluster->GetNCells() * medianPhotonRho);
- fhPhotonPtDiff->Fill(cluster->GetNCells() * medianPhotonRho);
- fhPhotonPtDiffVsCentrality->Fill(GetEventCentrality(),cluster->GetNCells() * medianPhotonRho);
- fhPhotonPtDiffVsNcells->Fill(numberOfcells,cluster->GetNCells() * medianPhotonRho);
- fhPhotonPtDiffVsNtracks->Fill(GetCTSTracks()->GetEntriesFast(),cluster->GetNCells() * medianPhotonRho);
- fhPhotonPtDiffVsNclusters->Fill(ntrig,cluster->GetNCells() * medianPhotonRho);
-
- fhPhotonPtCorrectedZoom->Fill(particlecorr->Pt() - cluster->GetNCells() * medianPhotonRho);
+ fhPhotonPtCorrected->Fill(particlecorr->Pt() - ncells * medianPhotonRho);
+ fhPhotonPtDiff->Fill(ncells * medianPhotonRho);
+ fhPhotonPtDiffVsCentrality->Fill(GetEventCentrality(),ncells * medianPhotonRho);
+ fhPhotonPtDiffVsNcells->Fill(numberOfcells,ncells * medianPhotonRho);
+ fhPhotonPtDiffVsNtracks->Fill(GetCTSTracks()->GetEntriesFast(),ncells * medianPhotonRho);
+ fhPhotonPtDiffVsNclusters->Fill(ntrig,ncells * medianPhotonRho);
+ fhPhotonPtCorrectedZoom->Fill(particlecorr->Pt() - ncells * medianPhotonRho);
//test: sum_pt in the cone 0.3 for each photon
//should be: random fake gamma from MB
for(Int_t iaod2 = 0; iaod2 < ntrig ; iaod2++){
if(iaod==iaod2) continue;
AliAODPWG4ParticleCorrelation* particlecorr2 = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod2));
- clusterID = particlecorr2->GetCaloLabel(0) ;
- if(clusterID < 0) continue;
- cluster2 = FindCluster(clusters,clusterID,iclustmp);
- photon2Corrected = particlecorr2->Pt() - cluster2->GetNCells() * medianPhotonRho;
-
+// clusterID = particlecorr2->GetCaloLabel(0) ;
+// if(clusterID < 0) continue;
+// cluster2 = FindCluster(clusters,clusterID,iclustmp);
+// photon2Corrected = particlecorr2->Pt() - cluster2->GetNCells() * medianPhotonRho;
+ photon2Corrected = particlecorr2->Pt() - particlecorr2->GetNCells() * medianPhotonRho;
+
//if(Pt()<0.5) continue; //<<hardcoded here //FIXME
if( TMath::Sqrt((particlecorr->Eta()-particlecorr2->Eta())*(particlecorr->Eta()-particlecorr2->Eta()) +
(particlecorr->Phi()-particlecorr2->Phi())*(particlecorr->Phi()-particlecorr2->Phi()) )<fGammaConeSize ){//if(/*cone is correct*/){
//
//Fill Correlation Histograms
//
- clusterID = particlecorr->GetCaloLabel(0) ;
- if(!(clusterID<0)){
- cluster = FindCluster(clusters,clusterID,iclustmp);
- //fill tree variables
- fGamNcells = cluster->GetNCells();
- }
+// clusterID = particlecorr->GetCaloLabel(0) ;
+// if(!(clusterID<0)){
+// cluster = FindCluster(clusters,clusterID,iclustmp);
+// //fill tree variables
+// fGamNcells = cluster->GetNCells();
+// }
+
+ fGamNcells = particlecorr->GetNCells();
+
Double_t ptTrig = particlecorr->Pt() - medianPhotonRho * fGamNcells;//<<---changed here
Double_t ptJet = jet->Pt() - rhoEvent * jet->EffectiveAreaCharged();//<<---changed here
Double_t phiTrig = particlecorr->Phi();
fhSelectedNtracks->Fill(GetCTSTracks()->GetEntriesFast());//to be checked
fhSelectedPhotonNLMVsPt->Fill(ptTrig,particlecorr->GetNLM());
-
- if(clusterID < 0 ){
- fhSelectedPhotonLambda0VsPt->Fill(ptTrig,-1);
- //fill tree variables
- fGamLambda0 = -1;
- fGamTime = -1;
- fGamNcells = 0;
- fGamSumPtNeu=0;
- }
- else{
+// if(clusterID < 0 ){
+// fhSelectedPhotonLambda0VsPt->Fill(ptTrig,-1);
+// //fill tree variables
+// fGamLambda0 = -1;
+// fGamTime = -1;
+// fGamNcells = 0;
+// fGamSumPtNeu=0;
+// }
+// else
+// {
//Int_t iclus = -1;
- // TObjArray* clusters = GetEMCALClusters();
+ TObjArray* clusters = GetEMCALClusters();
//cluster = FindCluster(clusters,clusterID,iclustmp);
- Double_t lambda0=cluster->GetM02();
+ Double_t lambda0=particlecorr->GetM02();
fhSelectedPhotonLambda0VsPt->Fill(ptTrig,lambda0);
//fill tree variables
fGamLambda0 = lambda0;
- fGamTime = cluster->GetTOF();
+ fGamTime = particlecorr->GetTime();
//fGamNcells = cluster->GetNCells();
fGamSumPtNeu=0;
//Double_t vectorLength=particlecorr->P();
for(Int_t icalo=0; icalo <clusters->GetEntriesFast(); icalo++){
AliVCluster* calo = (AliVCluster *) clusters->At(icalo);
- if(clusterID==calo->GetID()) continue;//the same cluster as trigger
+ //if(clusterID==calo->GetID()) continue;//the same cluster as trigger
calo->GetMomentum(fMomentum,vertex) ;//Assume that come from vertex in straight line
//printf("min pt %f\n",GetMinPt());
if(fMomentum.Pt()<GetMinPt()) continue; //<<hardcoded here //FIXME 0.5 check if correct
fGamNclusters++;
}
}
- }
+// }
//sum pt of charged tracks in the gamma isolation cone
//starts here
printf("Isolated Trigger? %d\n", fSelectIsolated) ;
printf("Reconstructed jet cone size = %3.2f\n", fJetConeSize) ;
printf("Reconstructed jet minimum pt before background subtraction = %3.2f\n", fJetMinPt) ;
+ printf("Reconstructed jet minimum pt after background subtraction = %3.2f\n", fJetMinPtBkgSub) ;
printf("Reconstructed jet minimum area fraction = %3.2f\n", fJetAreaFraction) ;
//if(!fNonStandardJetFromReader){