Skip to content

Commit 178c15c

Browse files
arvindkhuntiaArvindiiti
arvindkhuntia
and
Arvindiiti
authored
[PWGJE] Corrected posz selection, histograms for gen spectra, and centrality info added (#16792)
Co-authored-by: Arvind Khuntia <arvind.khuntia@cern.ch>
1 parent e06f9bf commit 178c15c

1 file changed

Lines changed: 92 additions & 52 deletions

File tree

‎PWGJE/Tasks/nucleiInJets.cxx‎

Lines changed: 92 additions & 52 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@
1717
#include "PWGJE/DataModel/JetSubtraction.h"
1818
//
1919
#include "PWGLF/DataModel/LFParticleIdentification.h"
20+
#include "PWGLF/DataModel/mcCentrality.h"
2021
#include "PWGLF/Utils/inelGt.h"
2122

2223
#include "Common/Core/RecoDecay.h"
@@ -190,6 +191,7 @@ struct nucleiInJets {
190191
// using EventTable = soa::Join<aod::JetCollisions, aod::EvSels, aod::CentFT0Ms, aod::CentFV0As, aod::CentFT0Cs>;
191192
using EventTable = aod::JetCollisions;
192193
using EventTableMC = soa::Join<EventTable, aod::JMcCollisionLbs>;
194+
using JetMcCollisionWithCent = soa::Join<aod::JetMcCollisions, aod::McCentFT0Cs, aod::McCentFV0As>::iterator;
193195
using JetCollWithLabel = o2::soa::Join<aod::JetCollisions, aod::JMcCollisionLbs, aod::BkgChargedRhos>::iterator;
194196
using TrackCandidates = soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullPi, aod::pidTPCFullKa,
195197
aod::pidTPCFullPr, aod::pidTPCFullDe, aod::pidTPCFullHe, aod::pidTPCFullTr, aod::pidTOFFullPi, aod::pidTOFFullKa,
@@ -600,13 +602,16 @@ struct nucleiInJets {
600602
// Event and signal loss analysis histograms (inclusive)
601603
jetHist.add("eventLoss/hEventStatistics", "Event Statistics for Loss Analysis", kTH1F, {{10, 0.f, 10.f}});
602604
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(1, "All Generated");
603-
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(2, "Gen |Vz|<10");
605+
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(2, "Gen |Vz| cut");
604606
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(3, "Gen True INEL>0");
605607
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(4, "Has Reco Coll");
606608
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(5, "Pass Sel8");
607-
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(6, "Pass |Vz|<10");
608-
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(7, "Pass rec INEL>0");
609-
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(8, "EvSelPassedRecINELgt0");
609+
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(6, "Pass reco |Vz| cut");
610+
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(7, "Pass extra event sel");
611+
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(8, "Reco selected");
612+
jetHist.get<TH1>(HIST("eventLoss/hEventStatistics"))->GetXaxis()->SetBinLabel(9, "Reco selected + true INEL>0");
613+
jetHist.add<TH2>("eventLoss/hRecoCollPerMCCollVsCent", "Reco collisions per MC collision;N_{reco coll};Centrality", HistType::kTH2F, {{10, -0.5f, 9.5f}, {100, 0.f, 100.f}});
614+
jetHist.add<TH2>("eventLoss/hSelectedRecoCollPerMCCollVsCent", "Selected reco collisions per MC collision;N_{selected reco coll};Centrality", HistType::kTH2F, {{10, -0.5f, 9.5f}, {100, 0.f, 100.f}});
610615

611616
// Signal loss histograms (only the ones that are actually used)
612617
jetHist.add<TH3>("eventLoss/signalLoss/h3GenParticlesPtVsEtaVsCent_INELgt0", "Generated Particles p_{T} vs #eta vs Centrality", HistType::kTH3F, {{PtAxis}, {100, -1.5f, 1.5f}, {100, 0, 100}});
@@ -733,6 +738,10 @@ struct nucleiInJets {
733738
jetHist.add<TH2>("eff/recmatched/perpCone/pt/PtParticleTypeTPCTOFVeto", "Pt (p) vs particletype", HistType::kTH2D, {{PtAxis}, {14, -7, 7}});
734739
jetHist.add<TH3>("eff/recmatched/gen/pt/PtParticleType", "Pt (p) vs jetflag vs particletype", HistType::kTH3D, {{PtAxis}, {2, 0, 2}, {14, -7, 7}});
735740
jetHist.add<TH2>("eff/recmatched/gen/perpCone/pt/PtParticleType", "Pt (p) vs particletype", HistType::kTH2D, {{PtAxis}, {14, -7, 7}});
741+
jetHist.add<TH3>("eff/recmatched/mcC/gen/pt/PtParticleType", "Pt (gen, mcC) vs jetflag vs particletype", HistType::kTH3D, {{PtAxis}, {2, 0, 2}, {14, -7, 7}});
742+
jetHist.add<TH3>("eff/recmatched/mcCSpectra/gen/pt/PtParticleType", "Pt (gen, mcCSpectra) vs jetflag vs particletype", HistType::kTH3D, {{PtAxis}, {2, 0, 2}, {14, -7, 7}});
743+
jetHist.add<TH2>("eff/recmatched/mcC/gen/perpCone/pt/PtParticleType", "Pt (gen, mcC, perp cone) vs particletype", HistType::kTH2D, {{PtAxis}, {14, -7, 7}});
744+
jetHist.add<TH2>("eff/recmatched/mcCSpectra/gen/perpCone/pt/PtParticleType", "Pt (gen, mcCSpectra, perp cone) vs particletype", HistType::kTH2D, {{PtAxis}, {14, -7, 7}});
736745
// gen matched
737746
jetHist.add<TH2>("genmatched/hRecMatchedJetPt", "matched jet pT (Rec level);#it{p}_{T,jet part} (GeV/#it{c}); #it{p}_{T,jet part} - #it{p}_{T,jet det}", HistType::kTH2F, {{100, 0., 100.}, {400, -20., 20.}});
738747
jetHist.add<TH2>("genmatched/hRecMatchedVsGenJetPt", "matched jet pT (Rec level);#it{p}_{T,jet det}; #it{p}_{T,jet part} (GeV/#it{c})", HistType::kTH2F, {{100, 0., 100.}, {100, 0., 100.}});
@@ -2414,8 +2423,22 @@ struct nucleiInJets {
24142423

24152424
if (mapPDGToValue(mcParticle.pdgCode()) != 0) {
24162425
jetHist.fill(HIST("eff/recmatched/gen/pt/PtParticleType"), mcParticle.pt(), jetFlagMC, mapPDGToValue(mcParticle.pdgCode()));
2426+
if (useMcC) {
2427+
if (useDataLikeHist) {
2428+
jetHist.fill(HIST("eff/recmatched/mcC/gen/pt/PtParticleType"), mcParticle.pt(), jetFlagMC, mapPDGToValue(mcParticle.pdgCode()));
2429+
} else {
2430+
jetHist.fill(HIST("eff/recmatched/mcCSpectra/gen/pt/PtParticleType"), mcParticle.pt(), jetFlagMC, mapPDGToValue(mcParticle.pdgCode()));
2431+
}
2432+
}
24172433
if (jetFlagPerpConeMC) {
24182434
jetHist.fill(HIST("eff/recmatched/gen/perpCone/pt/PtParticleType"), mcParticle.pt(), mapPDGToValue(mcParticle.pdgCode()));
2435+
if (useMcC) {
2436+
if (useDataLikeHist) {
2437+
jetHist.fill(HIST("eff/recmatched/mcC/gen/perpCone/pt/PtParticleType"), mcParticle.pt(), mapPDGToValue(mcParticle.pdgCode()));
2438+
} else {
2439+
jetHist.fill(HIST("eff/recmatched/mcCSpectra/gen/perpCone/pt/PtParticleType"), mcParticle.pt(), mapPDGToValue(mcParticle.pdgCode()));
2440+
}
2441+
}
24192442
}
24202443
}
24212444
} // mcParticle
@@ -2714,96 +2737,113 @@ struct nucleiInJets {
27142737
}
27152738

27162739
// Process function for event and signal loss analysis (inclusive)
2717-
void processEventSignalLoss(aod::JetMcCollision const& mcCollision,
2740+
void processEventSignalLoss(JetMcCollisionWithCent const& mcCollision,
27182741
soa::SmallGroups<soa::Join<aod::JMcCollisionLbs, aod::JetCollisions>> const& recoColls,
27192742
aod::JetParticles const& mcParticles,
27202743
TrackCandidates const&)
27212744
{
2722-
2723-
// Fill generated event statistics
27242745
jetHist.fill(HIST("eventLoss/hEventStatistics"), 0.5); // All Generated
27252746

2726-
// Check if we have a reconstructed collision
2727-
bool hasRecoColl = false;
2728-
bool passSel8 = false;
2729-
bool passVz = false;
2730-
bool passINELgt0 = false;
2731-
bool isSel8 = false;
2732-
2733-
float centrality = -999;
2747+
auto mcParticles_perColl = mcParticles.sliceBy(perMCCol, mcCollision.globalIndex());
2748+
float centrality = -999.f;
27342749
switch (centralityType) {
2735-
case 0: // FT0M
2750+
case 0:
27362751
centrality = mcCollision.centFT0M();
27372752
break;
2738-
case 1: // FT0C
2739-
centrality = mcCollision.multFT0C();
2753+
case 1:
2754+
centrality = mcCollision.centFT0C();
27402755
break;
2741-
case 2: // V0A
2742-
centrality = mcCollision.multFV0A();
2756+
case 2:
2757+
centrality = mcCollision.centFV0A();
27432758
break;
27442759
default:
2745-
centrality = -999;
2760+
centrality = mcCollision.centFT0M();
2761+
break;
27462762
}
2763+
const bool passGenVz = std::abs(mcCollision.posZ()) < cfgMaxZVertex;
2764+
const bool mcINELgt0 = o2::pwglf::isINELgt0mc(mcParticles_perColl, pdgDB);
27472765

2748-
// Check INEL>0 at MC level using PWGLF functionality
2749-
bool mcINELgt0 = o2::pwglf::isINELgt0mc(mcParticles, pdgDB);
2750-
if (mcCollision.posZ() < 10) {
2751-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 1.5);
2752-
if (mcINELgt0) {
2753-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 2.5);
2754-
}
2766+
if (!passGenVz) {
2767+
return;
2768+
}
2769+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 1.5); // Gen |Vz| cut
2770+
2771+
if (!mcINELgt0) {
2772+
return;
27552773
}
2774+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 2.5); // Gen True INEL>0
27562775

2776+
int nRecoColls = 0;
2777+
int nSelectedRecoColls = 0;
2778+
bool hasRecoColl = false;
2779+
bool passSel8 = false;
2780+
bool passRecoVz = false;
2781+
bool passExtraEventSel = false;
2782+
bool hasSelectedRecoColl = false;
27572783
for (const auto& recoColl : recoColls) {
2784+
++nRecoColls;
27582785
hasRecoColl = true;
2759-
if (jetderiveddatautilities::selectCollision(recoColl, jetderiveddatautilities::initialiseEventSelectionBits("sel8")))
2760-
isSel8 = true;
2761-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 3.5); // Has Reco Coll
2786+
const bool isSel8 = jetderiveddatautilities::selectCollision(recoColl, jetderiveddatautilities::initialiseEventSelectionBits("sel8"));
2787+
const bool isRecoVz = std::abs(recoColl.posZ()) < cfgMaxZVertex;
2788+
const bool isNoSameBunchPileup = !selNoSameBunchPileup || jetderiveddatautilities::selectCollision(recoColl, jetderiveddatautilities::initialiseEventSelectionBits("NoSameBunchPileup"));
2789+
const bool isGoodZvtxFT0vsPV = !selIsGoodZvtxFT0vsPV || jetderiveddatautilities::selectCollision(recoColl, jetderiveddatautilities::initialiseEventSelectionBits("IsGoodZvtxFT0vsPV"));
2790+
const bool isOccupancy = !useOccupancy || isOccupancyAccepted(recoColl);
2791+
const bool isExtraEventSel = isNoSameBunchPileup && isGoodZvtxFT0vsPV && isOccupancy;
27622792

27632793
if (isSel8) {
27642794
passSel8 = true;
2765-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 4.5); // Pass Sel8
27662795
}
2767-
2768-
if (std::abs(recoColl.posZ()) < 10.0) {
2769-
passVz = true;
2770-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 5.5); // Pass |Vz|<10
2796+
if (isRecoVz) {
2797+
passRecoVz = true;
27712798
}
2772-
2773-
if (mcINELgt0) {
2774-
passINELgt0 = true;
2775-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 6.5); // Pass rec INEL>0
2799+
if (isExtraEventSel) {
2800+
passExtraEventSel = true;
27762801
}
27772802

2778-
break; // Only first reco collision
2803+
if (isSel8 && isRecoVz && isExtraEventSel) {
2804+
++nSelectedRecoColls;
2805+
hasSelectedRecoColl = true;
2806+
}
27792807
}
27802808

2781-
// Final selection (all cuts passed)
2782-
if (hasRecoColl && passSel8 && passVz && passINELgt0) {
2783-
jetHist.fill(HIST("eventLoss/hEventStatistics"), 7.5); // Final Selection
2809+
jetHist.fill(HIST("eventLoss/hRecoCollPerMCCollVsCent"), nRecoColls, centrality);
2810+
jetHist.fill(HIST("eventLoss/hSelectedRecoCollPerMCCollVsCent"), nSelectedRecoColls, centrality);
2811+
2812+
if (hasRecoColl) {
2813+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 3.5); // Has Reco Coll
2814+
}
2815+
if (passSel8) {
2816+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 4.5); // Pass Sel8
2817+
}
2818+
if (passRecoVz) {
2819+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 5.5); // Pass reco |Vz| cut
2820+
}
2821+
if (passExtraEventSel) {
2822+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 6.5); // Pass extra event sel
2823+
}
2824+
if (hasSelectedRecoColl) {
2825+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 7.5); // Reco selected
2826+
jetHist.fill(HIST("eventLoss/hEventStatistics"), 8.5); // Reco selected + true INEL>0
27842827
}
27852828

2786-
auto mcParticles_perColl = mcParticles.sliceBy(perMCCol, mcCollision.globalIndex());
27872829
for (const auto& mcParticle : mcParticles_perColl) {
27882830
if (!mcParticle.isPhysicalPrimary())
27892831
continue;
27902832

27912833
// Apply kinematic cuts similar to track selection
27922834
if (std::fabs(mcParticle.eta()) > cfgtrkMaxEta)
27932835
continue;
2836+
if (std::fabs(mcParticle.y()) > cfgtrkMaxRap)
2837+
continue;
27942838

27952839
int particleType = mapPDGToValue(mcParticle.pdgCode());
27962840
if (particleType == 0)
27972841
continue; // Only interested particles
27982842

2799-
// Fill INEL>0 specific histograms
2800-
if (mcINELgt0) {
2801-
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticlesPtVsEtaVsCent_TrueINELgt0"), mcParticle.pt(), mcParticle.eta(), centrality);
2802-
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticleTypeVsPtVsCent_TrueINELgt0"), mcParticle.pt(), particleType, centrality);
2803-
}
2843+
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticlesPtVsEtaVsCent_TrueINELgt0"), mcParticle.pt(), mcParticle.eta(), centrality);
2844+
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticleTypeVsPtVsCent_TrueINELgt0"), mcParticle.pt(), particleType, centrality);
28042845

2805-
// Fill generated particle histograms (rec events)
2806-
if (hasRecoColl && passSel8 && passVz && passINELgt0) {
2846+
if (hasSelectedRecoColl) {
28072847
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticlesPtVsEtaVsCent_INELgt0"), mcParticle.pt(), mcParticle.eta(), centrality);
28082848
jetHist.fill(HIST("eventLoss/signalLoss/h3GenParticleTypeVsPtVsCent_INELgt0"), mcParticle.pt(), particleType, centrality);
28092849
}

0 commit comments

Comments
 (0)