From 73782f2b9a476d3c3cb14df4e3f0b4403073b9bb Mon Sep 17 00:00:00 2001 From: Yunseul Date: Thu, 22 Jan 2026 21:45:30 +0900 Subject: [PATCH] [PWGLF] Refine MC handling and expand QA - Add INEL>0 event selection for data and MC - Use ResoMC collision joins and pass MC FT0M centrality into histogram filling - Expand QA: event/track cutflow, track-stage, PID, and pair-level histograms --- PWGLF/Tasks/Resonances/f0980analysis.cxx | 397 +++++++++++++++++------ 1 file changed, 296 insertions(+), 101 deletions(-) diff --git a/PWGLF/Tasks/Resonances/f0980analysis.cxx b/PWGLF/Tasks/Resonances/f0980analysis.cxx index c8d1873eb22..6fa0d8d49f9 100644 --- a/PWGLF/Tasks/Resonances/f0980analysis.cxx +++ b/PWGLF/Tasks/Resonances/f0980analysis.cxx @@ -15,23 +15,25 @@ /// \since 01/07/2024 #include "PWGLF/DataModel/LFResonanceTables.h" +#include "PWGLF/DataModel/mcCentrality.h" +#include "PWGLF/Utils/inelGt.h" #include "Common/DataModel/Centrality.h" #include "Common/DataModel/EventSelection.h" +#include "Common/DataModel/PIDResponse.h" #include "CommonConstants/MathConstants.h" #include "CommonConstants/PhysicsConstants.h" -#include "DataFormatsParameters/GRPObject.h" #include "Framework/ASoAHelpers.h" #include "Framework/AnalysisHelpers.h" #include "Framework/AnalysisTask.h" #include "Framework/runDataProcessing.h" -#include #include "Math/LorentzVector.h" #include "Math/Vector4D.h" #include "TVector2.h" +#include #include using namespace o2; @@ -45,25 +47,21 @@ struct f0980analysis { SliceCache cache; HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; + using ResoMCCols = o2::soa::Join; + using ResoMCColsEP = o2::soa::Join; + // Event selections Configurable cfgMinpT{"cfgMinpT", 0.15, "Minimum transverse momentum for charged track"}; Configurable cfgMaxEta{"cfgMaxEta", 0.8, "Maximum pseudorapidiy for charged track"}; Configurable cfgMaxDCArToPVcut{"cfgMaxDCArToPVcut", 0.5, "Maximum transverse DCA"}; - Configurable cfgMaxDCAzToPVcut{"cfgMaxDCAzToPVcut", 2.0, "Maximum longitudinal DCA"}; + Configurable cfgMaxDCAzToPVcut{"cfgMaxDCAzToPVcut", 0.2, "Maximum longitudinal DCA"}; Configurable cfgMinRap{"cfgMinRap", -0.5, "Minimum rapidity for pair"}; Configurable cfgMaxRap{"cfgMaxRap", 0.5, "Maximum rapidity for pair"}; // Track selections - Configurable cfgPrimaryTrack{"cfgPrimaryTrack", false, "Primary track selection"}; // kGoldenChi2 | kDCAxy | kDCAz - Configurable cfgGlobalTrack{"cfgGlobalTrack", false, "Global track selection"}; // kGoldenChi2 | kDCAxy | kDCAz - Configurable cfgGlobalWoDCATrack{"cfgGlobalWoDCATrack", true, "Global track selection without DCA"}; // kQualityTracks (kTrackType | - // kTPCNCls | kTPCCrossedRows | - // kTPCCrossedRowsOverNCls | - // kTPCChi2NDF | kTPCRefit | - // kITSNCls | kITSChi2NDF | - // kITSRefit | kITSHits) | - // kInAcceptanceTracks (kPtRange | - // kEtaRange) + Configurable cfgPrimaryTrack{"cfgPrimaryTrack", false, "Primary track selection"}; + Configurable cfgGlobalTrack{"cfgGlobalTrack", true, "Global track selection"}; + Configurable cfgGlobalWoDCATrack{"cfgGlobalWoDCATrack", false, "Global track selection without DCA"}; Configurable cfgPVContributor{"cfgPVContributor", true, "PV contributor track selection"}; Configurable cfgUseTPCRefit{"cfgUseTPCRefit", false, "Require TPC Refit"}; Configurable cfgUseITSRefit{"cfgUseITSRefit", false, "Require ITS Refit"}; @@ -75,7 +73,7 @@ struct f0980analysis { Configurable cMaxTPCnSigmaPion{"cMaxTPCnSigmaPion", 5.0, "TPC nSigma cut for Pion"}; Configurable cMaxTPCnSigmaPionWoTOF{"cMaxTPCnSigmaPionWoTOF", 2.0, "TPC nSigma cut without TOF for Pion"}; Configurable nsigmaCutCombinedPion{"nsigmaCutCombinedPion", -999, "Combined nSigma cut for Pion"}; - Configurable selectType{"selectType", 0, "PID selection type"}; + Configurable selectType{"selectType", 3, "PID selection type"}; // Axis ConfigurableAxis massAxis{"massAxis", {400, 0.2, 2.2}, "Invariant mass axis"}; @@ -84,10 +82,19 @@ struct f0980analysis { Configurable cfgFindRT{"cfgFindRT", false, "boolean for RT analysis"}; Configurable cfgFindEP{"cfgFindEP", false, "boolean for Event plane analysis"}; + + Configurable cRecoINELgt0{"cRecoINELgt0", true, "check if INEL>0 for reco events"}; + Configurable cfgQAEPLT{"cfgQAEPLT", false, "Fill QA histograms for Event Plane and Leading Track"}; + Configurable cfgQACutflow{"cfgQACutflow", true, "Fill cutflow QA histograms"}; + Configurable cfgQAPairQA{"cfgQAPairQA", true, "Fill pair-level QA histograms"}; + Configurable cfgQATrackStages{"cfgQATrackStages", true, "Fill track QA at each stage"}; Configurable cfgQASelection{"cfgQASelection", true, "Fill QA histograms for Selection"}; Configurable cfgQAMCTrue{"cfgQAMCTrue", false, "Fill QA histograms for MC True Selection"}; + // MC event selection + Configurable cZvertCutMC{"cZvertCutMC", 10.0, "MC Z-vertex cut"}; + void init(o2::framework::InitContext&) { std::vector lptBinning = {0, 5.0, 13.0, 20.0, 50.0, 1000.0}; @@ -105,8 +112,12 @@ struct f0980analysis { AxisSpec rapqaAxis = {60, -1.5, 1.5, "#it{y}"}; // Rapidity axis AxisSpec dcaxyAxis = {200, -5.0, 5.0, "DCA_{xy} (cm)"}; // DCAxy axis AxisSpec dcazAxis = {200, -5.0, 5.0, "DCA_{z} (cm)"}; // DCAz axis + AxisSpec invMass1DAxis = {400, 0.2, 2.2, "M_{#pi#pi} (GeV/#it{c}^{2})"}; - AxisSpec collCutAxis = {4, -0.5, 3.5, "Collision cut index for MC"}; + AxisSpec evtCutflowAxis = {5, -1.5, 3.5, "Event cutflow step"}; + AxisSpec trkCutflowAxis = {8, -0.5, 7.5, "Track cutflow step"}; + AxisSpec pairRapAxis = {80, -2.0, 2.0, "y_{#pi#pi}"}; + AxisSpec mcSelAxis = {5, -1.5, 3.5, "MC event selection step"}; if (cfgFindRT) { histos.add("hInvMass_f0980_US_RT", "unlike invariant mass", {HistType::kTHnSparseF, {massAxis, pTAxis, centAxis, rtAxis, lptAxis}}); @@ -132,33 +143,56 @@ struct f0980analysis { histos.add("QA/LTpt", "", {HistType::kTH3F, {pTqaAxis, centAxis, phiqaAxis}}); } if (cfgQASelection) { - // General QA - histos.add("QA/TrackPt", "", {HistType::kTH1F, {pTqaAxis}}); - histos.add("QA/TrackEta", "", {HistType::kTH1F, {etaqaAxis}}); - histos.add("QA/TrackPhi", "", {HistType::kTH1F, {phiqaAxis}}); + // Event QA + histos.add("QAevent/hMultiplicityPercentSameE", "Multiplicity percentile of collision", HistType::kTH1F, {centAxis}); + // General track QA + if (cfgQACutflow) { + histos.add("QAtrack/hTrackCutflow", "Track cutflow; step; counts", HistType::kTH1F, {trkCutflowAxis}); + } + if (cfgQATrackStages) { + // stage: 0=input(initializer tracks), 1=base(kin/quality), 2=type(global/primary), 3=DCA, 4=PID(final pion candidates) + histos.add("QAtrack/stage/hPt_stage", "track pT vs stage; pT; stage", {HistType::kTH2F, {pTqaAxis, trkCutflowAxis}}); + histos.add("QAtrack/stage/hEta_stage", "track eta vs stage; eta; stage", {HistType::kTH2F, {etaqaAxis, trkCutflowAxis}}); + histos.add("QAtrack/stage/hPhi_stage", "track phi vs stage; phi; stage", {HistType::kTH2F, {phiqaAxis, trkCutflowAxis}}); + } + histos.add("QAtrack/TrackPt", "", {HistType::kTH1F, {pTqaAxis}}); + histos.add("QAtrack/TrackEta", "", {HistType::kTH1F, {etaqaAxis}}); + histos.add("QAtrack/TrackPhi", "", {HistType::kTH1F, {phiqaAxis}}); // Track selection QA - histos.add("QA/trkDCAxy_BC", "DCA_{xy} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis}); - histos.add("QA/trkDCAz_BC", "DCA_{z} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis}); - histos.add("QA/trkDCAxy", "DCA_{xy} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis}); - histos.add("QA/trkDCAz", "DCA_{z} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis}); + histos.add("QAtrack/trkDCAxy_BC", "DCA_{xy} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis}); + histos.add("QAtrack/trkDCAz_BC", "DCA_{z} for pion tracks (before cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis}); + histos.add("QAtrack/trkDCAxy", "DCA_{xy} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcaxyAxis}); + histos.add("QAtrack/trkDCAz", "DCA_{z} for pion tracks (after cuts)", HistType::kTH2F, {pTqaAxis, dcazAxis}); // PID QA - histos.add("QA/Nsigma_TPC_BC", "TPC n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); - histos.add("QA/Nsigma_TOF_BC", "TOF n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); - histos.add("QA/Nsigma_TPC_TOF_BC", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}}); - histos.add("QA/Nsigma_TPC", "TPC n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); - histos.add("QA/Nsigma_TOF", "TOF n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); - histos.add("QA/Nsigma_TPC_TOF", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TPC_BC", "TPC n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TOF_BC", "TOF n#sigma^{#pi} (before PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TPC_TOF_BC", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TPC", "TPC n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TPC}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TOF", "TOF n#sigma^{#pi} (after PID cuts); p_{T} (GeV/c); n#sigma_{TOF}^{#pi}", {HistType::kTH2F, {pTqaAxis, pidqaAxis}}); + histos.add("QApid/Nsigma_TPC_TOF", "", {HistType::kTH2F, {pidqaAxis, pidqaAxis}}); + if (cfgQAPairQA) { + histos.add("QApair/hYpipi_BC", "y_{#pi#pi} (before y cut); y; counts", HistType::kTH1F, {pairRapAxis}); + histos.add("QApair/hYpipi_AC", "y_{#pi#pi} (after y cut); y; counts", HistType::kTH1F, {pairRapAxis}); + // low-dimensional mass QA (for quick visual sanity) + histos.add("QApair/hMass_US", "M_{#pi#pi} US; M; counts", HistType::kTH1F, {invMass1DAxis}); + histos.add("QApair/hMass_LSpp", "M_{#pi#pi} ++; M; counts", HistType::kTH1F, {invMass1DAxis}); + histos.add("QApair/hMass_LSmm", "M_{#pi#pi} --; M; counts", HistType::kTH1F, {invMass1DAxis}); + } } if (doprocessMCRec) { - histos.add("MCL/hpT_f0980_REC", "Reconstructed f0 signals", HistType::kTH3F, {massAxis, pTqaAxis, centAxis}); + histos.add("MCL/hMass_f0980_REC", "M_{REC} f0", HistType::kTH3F, {massAxis, pTqaAxis, centAxis}); + histos.add("MCL/hpT_f0980_REC", "REC truth-matched f0; p_{T}; Centrality (%)", HistType::kTH2F, {pTqaAxis, centAxis}); + histos.add("MCL/hpT_f0980_REC_truePt", "REC truth-matched f0; p_{T}^{true mother}; Centrality (%)", HistType::kTH2F, {pTqaAxis, centAxis}); } if (doprocessMCTrue) { - // histos.add("MCL/hpT_f0980_GEN", "Generated f0 signals", HistType::kTH2F, {pTqaAxis, centAxis}); - histos.add("MCL/hpT_f0980_GEN", "Generated f0 signals; selIdx; p_{T} (GeV/c); Centrality (%)", HistType::kTH3F, {collCutAxis, pTqaAxis, centAxis}); + // histos.add("MCL/hpT_f0980_GEN", "GEN f0; p_{T} (GeV/c); Centrality (%)", HistType::kTH2F, {pTqaAxis, centAxis}); + histos.add("MCL/hpT_f0980_GEN", "GEN f0; step; p_{T} (GeV/c); centFT0M (%)", HistType::kTH3F, {mcSelAxis, pTqaAxis, centAxis}); if (cfgQAMCTrue) { histos.add("QAMCTrue/f0_pt_y", "Generated f0 ; #it{p}_{T} (GeV/#it{c}) ; #it{y}", HistType::kTH2F, {pTqaAxis, rapqaAxis}); histos.add("QAMCTrue/f0_pt_cent", "Generated f0 ; #it{p}_{T} (GeV/#it{c}); Centrality (%)", HistType::kTH2F, {pTqaAxis, centAxis}); + histos.add("QAMCTrue/f0YieldSteps", "MC true f0 counts per event-selection; step; counts", HistType::kTH1F, {evtCutflowAxis}); + // histos.add("QAMCTrue/trueINELgt0", "MC true INEL>0 flag; flag; counts", HistType::kTH1F, {inelAxis}); } } histos.print(); @@ -185,29 +219,56 @@ struct f0980analysis { } template - bool selTrack(const TrackType track) + bool selTrackBase(const TrackType& track) { if (std::abs(track.pt()) < cfgMinpT) return false; if (std::fabs(track.eta()) > cfgMaxEta) return false; - if (std::abs(track.dcaXY()) > cfgMaxDCArToPVcut) + if (track.tpcNClsFound() < cfgTPCcluster) return false; - if (std::abs(track.dcaZ()) > cfgMaxDCAzToPVcut) + if (cfgPVContributor && !track.isPVContributor()) return false; - if (track.tpcNClsFound() < cfgTPCcluster) + if (cfgUseITSRefit && !track.passedITSRefit()) + return false; + if (cfgUseTPCRefit && !track.passedTPCRefit()) return false; + + return true; + } + + template + bool selTrackType(const TrackType& track) + { if (cfgPrimaryTrack && !track.isPrimaryTrack()) return false; if (cfgGlobalTrack && !track.isGlobalTrack()) return false; if (cfgGlobalWoDCATrack && !track.isGlobalTrackWoDCA()) return false; - if (cfgPVContributor && !track.isPVContributor()) + + return true; + } + + template + bool selTrackDCA(const TrackType& track) + { + if (std::abs(track.dcaXY()) > cfgMaxDCArToPVcut) return false; - if (cfgUseITSRefit && !track.passedITSRefit()) + if (std::abs(track.dcaZ()) > cfgMaxDCAzToPVcut) return false; - if (cfgUseTPCRefit && !track.passedTPCRefit()) + + return true; + } + + template + bool selTrack(const TrackType& track) + { + if (!selTrackBase(track)) + return false; + if (!selTrackType(track)) + return false; + if (!selTrackDCA(track)) return false; if (cfgHasTOF && !track.hasTOF()) return false; @@ -245,8 +306,14 @@ struct f0980analysis { } template - void fillHistograms(const CollisionType& collision, const TracksType& dTracks) + void fillHistograms(const CollisionType& collision, const TracksType& dTracks, + aod::McParticles const* mcParticles = nullptr, + float centOverride = NAN) { + const float cent = std::isfinite(centOverride) ? centOverride : collision.cent(); + if (cfgQASelection) { + histos.fill(HIST("QAevent/hMultiplicityPercentSameE"), cent); + } double lhpT = 0.; double lhphi = 0.; double relphi = 0.; @@ -258,53 +325,115 @@ struct f0980analysis { } } if (cfgQAEPLT) - histos.fill(HIST("QA/LTpt"), lhpT, collision.cent(), lhphi); + histos.fill(HIST("QA/LTpt"), lhpT, cent, lhphi); } else if (cfgFindEP) { if (cfgQAEPLT) { - histos.fill(HIST("QA/EPhist"), collision.cent(), collision.evtPl()); - histos.fill(HIST("QA/hEPResAB"), collision.cent(), collision.evtPlResAB()); - histos.fill(HIST("QA/hEPResBC"), collision.cent(), collision.evtPlResBC()); - histos.fill(HIST("QA/hEPResAC"), collision.cent(), collision.evtPlResAC()); + histos.fill(HIST("QA/EPhist"), cent, collision.evtPl()); + histos.fill(HIST("QA/hEPResAB"), cent, collision.evtPlResAB()); + histos.fill(HIST("QA/hEPResBC"), cent, collision.evtPlResBC()); + histos.fill(HIST("QA/hEPResAC"), cent, collision.evtPlResAC()); } } ROOT::Math::LorentzVector> pion1, pion2, reco; - for (const auto& [trk1, trk2] : combinations(CombinationsStrictlyUpperIndexPolicy(dTracks, dTracks))) { - if (cfgQASelection) { - histos.fill(HIST("QA/trkDCAxy_BC"), trk1.pt(), trk1.dcaXY()); - histos.fill(HIST("QA/trkDCAz_BC"), trk1.pt(), trk1.dcaZ()); + // Single-track QA + if (cfgQASelection) { + for (const auto& trk : dTracks) { + if (cfgQACutflow) { + histos.fill(HIST("QAtrack/hTrackCutflow"), 0); // all tracks (initializer track) + } + if (cfgQATrackStages) { + histos.fill(HIST("QAtrack/stage/hPt_stage"), trk.pt(), 0); + histos.fill(HIST("QAtrack/stage/hEta_stage"), trk.eta(), 0); + histos.fill(HIST("QAtrack/stage/hPhi_stage"), trk.phi(), 0); + } + histos.fill(HIST("QAtrack/trkDCAxy_BC"), trk.pt(), trk.dcaXY()); + histos.fill(HIST("QAtrack/trkDCAz_BC"), trk.pt(), trk.dcaZ()); + + if (!selTrackBase(trk)) { + continue; + } + if (cfgQACutflow) { + histos.fill(HIST("QAtrack/hTrackCutflow"), 1); // after kinematics cut + } + if (cfgQATrackStages) { + histos.fill(HIST("QAtrack/stage/hPt_stage"), trk.pt(), 1); + histos.fill(HIST("QAtrack/stage/hEta_stage"), trk.eta(), 1); + histos.fill(HIST("QAtrack/stage/hPhi_stage"), trk.phi(), 1); + } + histos.fill(HIST("QApid/Nsigma_TPC_BC"), trk.pt(), trk.tpcNSigmaPi()); + if (trk.hasTOF()) { + histos.fill(HIST("QApid/Nsigma_TOF_BC"), trk.pt(), trk.tofNSigmaPi()); + histos.fill(HIST("QApid/Nsigma_TPC_TOF_BC"), trk.tpcNSigmaPi(), trk.tofNSigmaPi()); + } + + if (!selTrackType(trk)) { + continue; + } + if (cfgQACutflow) { + histos.fill(HIST("QAtrack/hTrackCutflow"), 2); // after track type + } + if (cfgQATrackStages) { + histos.fill(HIST("QAtrack/stage/hPt_stage"), trk.pt(), 2); + histos.fill(HIST("QAtrack/stage/hEta_stage"), trk.eta(), 2); + histos.fill(HIST("QAtrack/stage/hPhi_stage"), trk.phi(), 2); + } + + if (!selTrackDCA(trk)) { + continue; + } + histos.fill(HIST("QAtrack/trkDCAxy"), trk.pt(), trk.dcaXY()); + histos.fill(HIST("QAtrack/trkDCAz"), trk.pt(), trk.dcaZ()); + if (cfgQACutflow) { + histos.fill(HIST("QAtrack/hTrackCutflow"), 3); // after DCA + } + if (cfgQATrackStages) { + histos.fill(HIST("QAtrack/stage/hPt_stage"), trk.pt(), 3); + histos.fill(HIST("QAtrack/stage/hEta_stage"), trk.eta(), 3); + histos.fill(HIST("QAtrack/stage/hPhi_stage"), trk.phi(), 3); + } + + if (!selPion(trk)) { + continue; + } + if (cfgQACutflow) { + histos.fill(HIST("QAtrack/hTrackCutflow"), 4); // after PID for pion selection + } + if (cfgQATrackStages) { + histos.fill(HIST("QAtrack/stage/hPt_stage"), trk.pt(), 4); + histos.fill(HIST("QAtrack/stage/hEta_stage"), trk.eta(), 4); + histos.fill(HIST("QAtrack/stage/hPhi_stage"), trk.phi(), 4); + } + histos.fill(HIST("QApid/Nsigma_TPC"), trk.pt(), trk.tpcNSigmaPi()); + if (trk.hasTOF()) { + histos.fill(HIST("QApid/Nsigma_TOF"), trk.pt(), trk.tofNSigmaPi()); + histos.fill(HIST("QApid/Nsigma_TPC_TOF"), trk.tpcNSigmaPi(), trk.tofNSigmaPi()); + } + + histos.fill(HIST("QAtrack/TrackPt"), trk.pt()); + histos.fill(HIST("QAtrack/TrackEta"), trk.eta()); + histos.fill(HIST("QAtrack/TrackPhi"), trk.phi()); } + } + + for (const auto& [trk1, trk2] : combinations(CombinationsStrictlyUpperIndexPolicy(dTracks, dTracks))) { if (!selTrack(trk1) || !selTrack(trk2)) continue; - if (cfgQASelection) { - histos.fill(HIST("QA/trkDCAxy"), trk1.pt(), trk1.dcaXY()); - histos.fill(HIST("QA/trkDCAz"), trk1.pt(), trk1.dcaZ()); - - histos.fill(HIST("QA/Nsigma_TPC_BC"), trk1.pt(), trk1.tpcNSigmaPi()); - if (trk1.hasTOF()) { - histos.fill(HIST("QA/Nsigma_TOF_BC"), trk1.pt(), trk1.tofNSigmaPi()); - histos.fill(HIST("QA/Nsigma_TPC_TOF_BC"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi()); - } - } if (!selPion(trk1) || !selPion(trk2)) continue; - if (cfgQASelection) { - histos.fill(HIST("QA/Nsigma_TPC"), trk1.pt(), trk1.tpcNSigmaPi()); - if (trk1.hasTOF()) { - histos.fill(HIST("QA/Nsigma_TOF"), trk1.pt(), trk1.tofNSigmaPi()); - histos.fill(HIST("QA/Nsigma_TPC_TOF"), trk1.tpcNSigmaPi(), trk1.tofNSigmaPi()); - } - - histos.fill(HIST("QA/TrackPt"), trk1.pt()); - histos.fill(HIST("QA/TrackEta"), trk1.eta()); - histos.fill(HIST("QA/TrackPhi"), trk1.phi()); - } pion1 = ROOT::Math::PxPyPzMVector(trk1.px(), trk1.py(), trk1.pz(), massPi); pion2 = ROOT::Math::PxPyPzMVector(trk2.px(), trk2.py(), trk2.pz(), massPi); reco = pion1 + pion2; + if (cfgQASelection && cfgQAPairQA) { + histos.fill(HIST("QApair/hYpipi_BC"), reco.Rapidity()); + } + if (reco.Rapidity() > cfgMaxRap || reco.Rapidity() < cfgMinRap) continue; + if (cfgQASelection && cfgQAPairQA) { + histos.fill(HIST("QApair/hYpipi_AC"), reco.Rapidity()); + } if (cfgFindEP) { relphi = TVector2::Phi_0_2pi(reco.Phi() - collision.evtPl()); @@ -314,60 +443,119 @@ struct f0980analysis { } if (trk1.sign() * trk2.sign() < 0) { + if (cfgQASelection && cfgQAPairQA) { + histos.fill(HIST("QApair/hMass_US"), reco.M()); + } if (cfgFindRT) { - histos.fill(HIST("hInvMass_f0980_US_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT); + histos.fill(HIST("hInvMass_f0980_US_RT"), reco.M(), reco.Pt(), cent, rtIndex(reco.Phi(), lhphi), lhpT); } else if (cfgFindEP) { - histos.fill(HIST("hInvMass_f0980_US_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi); + histos.fill(HIST("hInvMass_f0980_US_EPA"), reco.M(), reco.Pt(), cent, relphi); } else { - histos.fill(HIST("hInvMass_f0980_US"), reco.M(), reco.Pt(), collision.cent()); + histos.fill(HIST("hInvMass_f0980_US"), reco.M(), reco.Pt(), cent); } if constexpr (IsMC) { + // truth-matched f0 -> pi+pi- (signal-only REC) if (std::abs(trk1.pdgCode()) != kPiPlus || std::abs(trk2.pdgCode()) != kPiPlus) continue; if (trk1.motherId() != trk2.motherId()) continue; if (std::abs(trk1.motherPDG()) != 9010221) continue; - histos.fill(HIST("MCL/hpT_f0980_REC"), reco.M(), reco.Pt(), collision.cent()); + + histos.fill(HIST("MCL/hMass_f0980_REC"), reco.M(), reco.Pt(), cent); + histos.fill(HIST("MCL/hpT_f0980_REC"), reco.Pt(), cent); + if (mcParticles && trk1.motherId() >= 0) { + auto mother = mcParticles->iteratorAt(trk1.motherId()); + histos.fill(HIST("MCL/hpT_f0980_REC_truePt"), mother.pt(), cent); + } } } else if (trk1.sign() > 0 && trk2.sign() > 0) { + if (cfgQASelection && cfgQAPairQA) { + histos.fill(HIST("QApair/hMass_LSpp"), reco.M()); + } if (cfgFindRT) { - histos.fill(HIST("hInvMass_f0980_LSpp_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT); + histos.fill(HIST("hInvMass_f0980_LSpp_RT"), reco.M(), reco.Pt(), cent, rtIndex(reco.Phi(), lhphi), lhpT); } else if (cfgFindEP) { - histos.fill(HIST("hInvMass_f0980_LSpp_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi); + histos.fill(HIST("hInvMass_f0980_LSpp_EPA"), reco.M(), reco.Pt(), cent, relphi); } else { - histos.fill(HIST("hInvMass_f0980_LSpp"), reco.M(), reco.Pt(), collision.cent()); + histos.fill(HIST("hInvMass_f0980_LSpp"), reco.M(), reco.Pt(), cent); } } else if (trk1.sign() < 0 && trk2.sign() < 0) { + if (cfgQASelection && cfgQAPairQA) { + histos.fill(HIST("QApair/hMass_LSmm"), reco.M()); + } if (cfgFindRT) { - histos.fill(HIST("hInvMass_f0980_LSmm_RT"), reco.M(), reco.Pt(), collision.cent(), rtIndex(reco.Phi(), lhphi), lhpT); + histos.fill(HIST("hInvMass_f0980_LSmm_RT"), reco.M(), reco.Pt(), cent, rtIndex(reco.Phi(), lhphi), lhpT); } else if (cfgFindEP) { - histos.fill(HIST("hInvMass_f0980_LSmm_EPA"), reco.M(), reco.Pt(), collision.cent(), relphi); + histos.fill(HIST("hInvMass_f0980_LSmm_EPA"), reco.M(), reco.Pt(), cent, relphi); } else { - histos.fill(HIST("hInvMass_f0980_LSmm"), reco.M(), reco.Pt(), collision.cent()); + histos.fill(HIST("hInvMass_f0980_LSmm"), reco.M(), reco.Pt(), cent); } } } } - void processData(o2::soa::Join::iterator const& collision, + void processData(o2::soa::Join::iterator const& resoCollision, + o2::aod::ResoCollisionColls const& collisionIndex, + o2::soa::Join const& collisions, o2::aod::ResoTracks const& resotracks) { - fillHistograms(collision, resotracks); + if (cRecoINELgt0) { + auto linkRow = collisionIndex.iteratorAt(resoCollision.globalIndex()); + auto collId = linkRow.collisionId(); + if (collId < 0) + return; + auto coll = collisions.iteratorAt(collId); + if (!coll.isInelGt0()) + return; + } + fillHistograms(resoCollision, resotracks); } PROCESS_SWITCH(f0980analysis, processData, "Process Event for data", true); - void processMCRec(o2::soa::Join::iterator const& collision, - o2::soa::Join const& resotracks) + void processMCRec(ResoMCColsEP::iterator const& collision, + o2::aod::ResoCollisionColls const& collisionIndex, + o2::soa::Join const& collisionsMC, + o2::soa::Join const& resotracks, + o2::aod::McParticles const& mcParticles, + o2::soa::Join const&) { - fillHistograms(collision, resotracks); + auto linkRow = collisionIndex.iteratorAt(collision.globalIndex()); + auto collId = linkRow.collisionId(); + if (collId < 0) { + return; + } + auto coll = collisionsMC.iteratorAt(collId); + if (cRecoINELgt0 && !coll.isInelGt0()) + return; + auto mcColl = coll.mcCollision_as>(); + const float centFT0M = mcColl.centFT0M(); + if (!collision.isInAfterAllCuts() || (std::abs(collision.posZ()) > cZvertCutMC)) { + return; + } + + fillHistograms(collision, resotracks, &mcParticles, centFT0M); } - PROCESS_SWITCH(f0980analysis, processMCRec, "Process Event for MC", false); + PROCESS_SWITCH(f0980analysis, processMCRec, "Process Event for MC (REC)", false); - void processMCTrue(o2::soa::Join::iterator const& resoCollision, - o2::aod::ResoMCParents const& resoParents) + void processMCTrue(ResoMCCols::iterator const& resoCollision, + o2::aod::ResoCollisionColls const& collisionIndex, + o2::aod::ResoMCParents const& resoParents, + o2::aod::ResoCollisionCandidatesMC const& collisionsMC, + o2::soa::Join const&) { + auto linkRow = collisionIndex.iteratorAt(resoCollision.globalIndex()); + auto collId = linkRow.collisionId(); + if (collId < 0) { + return; + } + auto coll = collisionsMC.iteratorAt(collId); + if (cRecoINELgt0 && !coll.isInelGt0()) + return; + auto mcColl = coll.mcCollision_as>(); + const float centFT0M = mcColl.centFT0M(); + for (const auto& part : resoParents) { // loop over all pre-filtered MC particles if (std::abs(part.pdgCode()) != 9010221) continue; @@ -376,29 +564,36 @@ struct f0980analysis { if (part.y() < cfgMinRap || part.y() > cfgMaxRap) { continue; } - bool pass = false; - if ((std::abs(part.daughterPDG1()) == kPiPlus && std::abs(part.daughterPDG2()) == kPiPlus)) { - pass = true; - } - if (!pass) // If we have both decay products + const bool decayPions = (std::abs(part.daughterPDG1()) == kPiPlus && std::abs(part.daughterPDG2()) == kPiPlus); + if (!decayPions) continue; - // no event selection - histos.fill(HIST("MCL/hpT_f0980_GEN"), 0, part.pt(), resoCollision.cent()); - // |zvtx|<10 cm + // no event selection: baseline + histos.fill(HIST("MCL/hpT_f0980_GEN"), -1, part.pt(), centFT0M); + histos.fill(HIST("QAMCTrue/f0YieldSteps"), -1); + // |zvtx|<10 if (resoCollision.isVtxIn10()) { - histos.fill(HIST("MCL/hpT_f0980_GEN"), 1, part.pt(), resoCollision.cent()); + histos.fill(HIST("MCL/hpT_f0980_GEN"), 0, part.pt(), centFT0M); + histos.fill(HIST("QAMCTrue/f0YieldSteps"), 0); + } + // |zvtx|<10 && sel8 + if (resoCollision.isVtxIn10() && resoCollision.isInSel8()) { + histos.fill(HIST("MCL/hpT_f0980_GEN"), 1, part.pt(), centFT0M); + histos.fill(HIST("QAMCTrue/f0YieldSteps"), 1); } - // |zvtx|<10 cm & TVX trigger + // |zvtx|<10 && TVX if (resoCollision.isVtxIn10() && resoCollision.isTriggerTVX()) { - histos.fill(HIST("MCL/hpT_f0980_GEN"), 2, part.pt(), resoCollision.cent()); + histos.fill(HIST("MCL/hpT_f0980_GEN"), 2, part.pt(), centFT0M); + histos.fill(HIST("QAMCTrue/f0YieldSteps"), 2); } + // after all event cuts if (resoCollision.isInAfterAllCuts()) { - histos.fill(HIST("MCL/hpT_f0980_GEN"), 3, part.pt(), resoCollision.cent()); + histos.fill(HIST("MCL/hpT_f0980_GEN"), 3, part.pt(), centFT0M); + histos.fill(HIST("QAMCTrue/f0YieldSteps"), 3); } if (cfgQAMCTrue) { histos.fill(HIST("QAMCTrue/f0_pt_y"), part.pt(), part.y()); - histos.fill(HIST("QAMCTrue/f0_pt_cent"), part.pt(), resoCollision.cent()); + histos.fill(HIST("QAMCTrue/f0_pt_cent"), part.pt(), centFT0M); } } };