Skip to content

[PWGLF] Fix Omega2012 and K1 MC efficiency denominators and truth matching #18431

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
BongHwi wants to merge 3 commits into
base: master
Choose a base branch
from
+288 −49
Open

[PWGLF] Fix Omega2012 and K1 MC efficiency denominators and truth matching #18431

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
4f7fdfb
[PWGLF] Align Omega2012 MC efficiency and exporter event selection
Oct 23, 2026
132b8b8
[PWGLF] Align K1 MC parent truth and efficiency response
Oct 23, 2026
116f23d
[PWGLF] Fix clang-tidy findings in Omega2012 and K1 tables
Oct 23, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
57 changes: 56 additions & 1 deletion PWGLF/Core/K1AnalysisMicroCore.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,7 @@
#include <limits>
#include <set>
#include <type_traits>
#include <unordered_map>
#include <utility>
#include <vector>

Expand All @@ -64,6 +65,16 @@ enum class K1TruthChannel {
KStarPi = 2
};

struct K1ParentTruth {
double pt = std::numeric_limits<double>::quiet_NaN();
double y = std::numeric_limits<double>::quiet_NaN();
int pdgCode = 0;
K1TruthChannel channel = K1TruthChannel::None;
};
using K1ParentTruthMap = std::unordered_map<int64_t, K1ParentTruth>;

inline K1TruthChannel classifyGeneratedK1(int charge, int daughter1, int daughter2);

// The value is the species index of the PID configuration in ResoAnalysisSelectionCore.
enum class Species : int {
Pion = 0,
Expand Down Expand Up @@ -217,6 +228,49 @@ K1TruthChannel classifyK1Truth(const Track& pion1, const Track& pion2, const Kao
return K1TruthChannel::None;
}

template <typename Kaon, typename Pion>
int k1MotherId(K1TruthChannel channel, const Kaon& kaon, const Pion& pion1, const Pion& pion2)
{
if (channel == K1TruthChannel::RhoK) {
return kaon.motherId();
}
if (channel == K1TruthChannel::KStarPi) {
return std::abs(pion1.motherPDG()) == o2::constants::physics::Pdg::kK1_1270Plus ? pion1.motherId() : pion2.motherId();
}
return -1;
}

template <typename ParentsType>
K1ParentTruthMap makeK1ParentTruthMap(const ParentsType& parents)
{
K1ParentTruthMap result;
for (const auto& part : parents) {
if (std::abs(part.pdgCode()) != o2::constants::physics::Pdg::kK1_1270Plus || part.originalMcParticleId() < 0) {
continue;
}
const auto id = static_cast<int64_t>(part.originalMcParticleId());
const int charge = part.pdgCode() > 0 ? 1 : -1;
result.emplace(id, K1ParentTruth{part.pt(), part.y(), part.pdgCode(), classifyGeneratedK1(charge, part.daughterPDG1(), part.daughterPDG2())});
}
return result;
}

inline const K1ParentTruth* findK1ParentTruth(const K1ParentTruthMap& parents, int64_t motherId,
K1TruthChannel channel, int charge)
{
const auto it = parents.find(motherId);
if (it == parents.end()) {
return nullptr;
}
const auto& parent = it->second;
const int expectedPdg = charge * o2::constants::physics::Pdg::kK1_1270Plus;
if (channel == K1TruthChannel::None || parent.channel != channel || parent.pdgCode != expectedPdg ||
!std::isfinite(parent.pt) || parent.pt < 0. || !std::isfinite(parent.y)) {
return nullptr;
}
return &parent;
}

// Immediate decay channel of a generated K1 from the PDG codes of its two daughters.
inline K1TruthChannel classifyGeneratedK1(int charge, int daughter1, int daughter2)
{
Expand Down Expand Up @@ -635,7 +689,8 @@ class K1AnalysisMicroCore
const int charge = part.pdgCode() > 0 ? 1 : -1;
const K1TruthChannel channel = classifyGeneratedK1(charge, part.daughterPDG1(), part.daughterPDG2());
histos.fill(HIST("CutFlow/generated"), 0, static_cast<int>(channel));
if (part.y() < mCandidateCuts.cK1MinRap || part.y() > mCandidateCuts.cK1MaxRap) {
if (!std::isfinite(part.pt()) || part.pt() < 0. || !std::isfinite(part.y()) ||
part.y() < mCandidateCuts.cK1MinRap || part.y() > mCandidateCuts.cK1MaxRap) {
continue;
}
histos.fill(HIST("CutFlow/generated"), 1, static_cast<int>(channel));
Expand Down
22 changes: 17 additions & 5 deletions PWGLF/Core/Omega2012AnalysisCore.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -58,6 +58,17 @@ using o2::analysis::resonance::ResoAnalysisSelectionCore;
using o2::analysis::resonance::TrackCuts;
using o2::analysis::resonance::TrackStage;

// Omega(2012) analysis and ML export share these event-level selection defaults.
// Keep the generic EventCuts defaults unchanged for other resonance tasks.
inline EventCuts makeOmega2012EventCuts()
{
EventCuts cuts;
cuts.cRecoINELgt0.value = true;
cuts.cMCINELgt0.value = true;
cuts.cMCVtxIn10.value = true;
return cuts;
}

enum class DecayMode : uint8_t {
XiK0s = 0, // Omega(2012)- -> Xi- K0S
Xi1530K = 1 // Omega(2012)- -> Xi(1530)0 K- -> Xi- pi+ K-
Expand Down Expand Up @@ -932,19 +943,20 @@ class Omega2012AnalysisCore
}
}

// Generated Omega(2012) parents of a selected reconstructed MC collision (ResoMCParents_001). The optional callback
// receives (parent, immediate channel) for the parents inside the rapidity window. Split reconstructed collisions
// repeat parent sets: this is not an unconditional generated denominator.
// Generated Omega(2012) parents of a selected reconstructed MC collision (ResoMCParents_001). The caller supplies
// the denominator rapidity limit so it can remain independent of a wider candidate-export window. The optional
// callback receives (parent, immediate channel) for the parents inside that rapidity window. Split reconstructed
// collisions repeat parent sets: this is not an unconditional generated denominator.
template <typename ParentsType, typename Callback = std::nullptr_t>
void forEachGeneratedOmega2012(o2::framework::HistogramRegistry& histos, const ParentsType& resoParents, Callback callback = nullptr)
void forEachGeneratedOmega2012(o2::framework::HistogramRegistry& histos, const ParentsType& resoParents, float rapidityMax, Callback callback = nullptr)
{
for (const auto& part : resoParents) {
if (std::abs(part.pdgCode()) != o2::constants::physics::Pdg::kOmega2012Minus) {
continue;
}
const GeneratedChannel channel = classifyGeneratedOmega2012(part.pdgCode(), part.daughterPDG1(), part.daughterPDG2());
histos.fill(HIST("CutFlow/generated"), 0, static_cast<int>(channel));
if (!(std::abs(part.y()) < mCandidateCuts.cfgRapidityCut.value)) {
if (!(std::abs(part.y()) < rapidityMax)) {
continue;
}
histos.fill(HIST("CutFlow/generated"), 1, static_cast<int>(channel));
Expand Down
5 changes: 3 additions & 2 deletions PWGLF/DataModel/LFK1MlTables.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,14 +63,15 @@ DECLARE_SOA_COLUMN(K1Charge, k1Charge, int8_t);
// 16 = candidate cuts.
// Candidates written at the "selected" export stage carry 31.
DECLARE_SOA_COLUMN(K1BaselinePassBits, k1BaselinePassBits, uint16_t);
// NOLINTNEXTLINE(modernize-avoid-c-arrays) -- persistent fixed-size AOD column
DECLARE_SOA_COLUMN(K1MasterFeatures, k1MasterFeatures, float[125]); //! 125 master features of the K1 ML feature contract (PWGLF/Core/K1MlFeatures.h)
DECLARE_SOA_COLUMN(K1FeatureStatus, k1FeatureStatus, uint8_t); //! o2::analysis::k1ml::BuildStatus; only Ok (0) rows are written
DECLARE_SOA_COLUMN(K1TruthStatus, k1TruthStatus, uint8_t); //! 0 data, 1 matched, 2 unmatched
DECLARE_SOA_COLUMN(K1TruthChannel, k1TruthChannel, uint8_t); //! 0 none, 1 rho K, 2 K* pi
DECLARE_SOA_COLUMN(K1MotherPdg, k1MotherPdg, int32_t);
DECLARE_SOA_COLUMN(K1MotherId, k1MotherId, int64_t);
DECLARE_SOA_COLUMN(K1GeneratedPt, k1GeneratedPt, float); //! NaN in K1MlTruth: not filled for reconstructed candidates
DECLARE_SOA_COLUMN(K1GeneratedY, k1GeneratedY, float); //! NaN in K1MlTruth: not filled for reconstructed candidates
DECLARE_SOA_COLUMN(K1GeneratedPt, k1GeneratedPt, float); //! Parent-table pT for a validated reconstructed MC match; NaN otherwise
DECLARE_SOA_COLUMN(K1GeneratedY, k1GeneratedY, float); //! Parent-table y for a validated reconstructed MC match; NaN otherwise
DECLARE_SOA_COLUMN(K1OriginalMcParticleId, k1OriginalMcParticleId, int64_t);
DECLARE_SOA_COLUMN(K1DaughterPdg1, k1DaughterPdg1, int32_t);
DECLARE_SOA_COLUMN(K1DaughterPdg2, k1DaughterPdg2, int32_t);
Expand Down
40 changes: 30 additions & 10 deletions PWGLF/Tasks/Resonances/k1AnalysisMicro.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ struct HistogramOptions : ConfigurableGroup {
Configurable<int> cNbinsDiv{"cNbinsDiv", 1, "Integer to divide the number of bins"};
Configurable<bool> additionalQAplots{"additionalQAplots", true, "Additional QA plots"};
Configurable<int> cfgTruthDebug{"cfgTruthDebug", 0, "Maximum logged matched candidates per truth channel"};
Configurable<bool> cUseTruthRapidity{"cUseTruthRapidity", true, "Require generated K1 rapidity in the configured window for MC efficiency numerators"};
};

enum BinAnti : unsigned int {
Expand Down Expand Up @@ -231,12 +232,15 @@ struct K1AnalysisMicro {
AxisSpec channelAxis = {3, -0.5, 2.5, "0: non-K1, 1: rho K, 2: K* pi"};
histos.add("MCReco/collisions", "Selected reconstructed MC collisions", HistType::kTH1D, {{1, 0, 1}});
histos.add("MCReco/microTracks", "Input micro tracks in selected MC collisions", HistType::kTH1D, {{1, 0, 1}});
histos.add("MCReco/channel", "All selected pi-pi-K combinations by truth channel", HistType::kTH1D, {channelAxis});
// These channel and pT histograms are pre-truth-fiducial QA; MCReco/response is the validated numerator.
histos.add("MCReco/channel", "Pre-truth-fiducial QA: selected pi-pi-K combinations by truth channel", HistType::kTH1D, {channelAxis});
histos.add("MCReco/mass", "Reconstructed mass by truth channel", HistType::kTH2D, {channelAxis, invMassAxisReso});
histos.add("MCReco/pt", "Reconstructed pT by truth channel", HistType::kTH2D, {channelAxis, ptAxis});
histos.add("MCReco/pt", "Pre-truth-fiducial QA: reconstructed pT by truth channel", HistType::kTH2D, {channelAxis, ptAxis});
histos.add("MCReco/piPiMass", "pi-pi mass by truth channel", HistType::kTH2D, {channelAxis, invMassAxisRho});
histos.add("MCReco/pi1KMass", "First pion-kaon mass by truth channel", HistType::kTH2D, {channelAxis, invMassAxisK892});
histos.add("MCReco/pi2KMass", "Second pion-kaon mass by truth channel", HistType::kTH2D, {channelAxis, invMassAxisK892});
histos.add("MCReco/response", "Selected MC candidates;reconstructed pT;generated pT;channel;charge;centrality", HistType::kTHnSparseD,
{ptAxis, ptAxis, channelAxis, AxisSpec{2, -1.5, 1.5, "K1 charge"}, centAxis});
histos.add("k1invmass_MC", "Invariant mass of K1(1270)", HistType::kTH1F, {invMassAxisReso});
histos.add("k1invmass_MC_noK1", "Invariant mass of K1(1270)", HistType::kTH1F, {invMassAxisReso});

Expand Down Expand Up @@ -272,6 +276,8 @@ struct K1AnalysisMicro {
AxisSpec channelAxis = {3, -0.5, 2.5, "0: other/unresolved, 1: rho K, 2: K* pi"};
histos.add("MCGen/chargeChannel", "K1 parents in selected reconstructed events, inside the K1 rapidity window", HistType::kTH2D, {{2, -1.5, 1.5, "K1 charge"}, channelAxis});
histos.add("MCGen/ptChannel", "Generated K1 pT by immediate decay channel", HistType::kTH2D, {channelAxis, ptAxis});
histos.add("MCGen/response", "Generated K1 parents;generated pT;channel;charge;centrality", HistType::kTHnSparseD,
{ptAxis, channelAxis, AxisSpec{2, -1.5, 1.5, "K1 charge"}, centAxis});
}
}

Expand Down Expand Up @@ -380,7 +386,8 @@ struct K1AnalysisMicro {
// (or one mixed pair of collisions). The selection itself is the shared K1 core.
// dTracks1: bachelor kaons, dTracks2: pions.
template <bool IsMC, bool IsMix, bool IsResoMicrotrack, typename CollisionType, typename TracksType>
void fillHistograms(const CollisionType& collision, const TracksType& dTracks1, const TracksType& dTracks2)
void fillHistograms(const CollisionType& collision, const TracksType& dTracks1, const TracksType& dTracks2,
const K1ParentTruthMap* parentTruth = nullptr)
{
const bool fillQA = !IsMix && histogramOptions.additionalQAplots;
const auto multiplicity = collision.cent();
Expand Down Expand Up @@ -474,9 +481,17 @@ struct K1AnalysisMicro {
pion1.motherPDG(), pion2.motherPDG(), kaon.motherPDG(),
pion1.siblingIds()[0], pion1.siblingIds()[1], pion2.siblingIds()[0], pion2.siblingIds()[1], kaon.siblingIds()[0], kaon.siblingIds()[1]);
}
typeK1 = kaon.sign() > 0 ? BinType::kK1P_Rec : BinType::kK1N_Rec;
histos.fill(HIST("hInvmass_K1"), typeNormal, typeK1, multiplicity, c.k1.Pt(), c.k1.M());
histos.fill(HIST("k1invmass_MC"), c.k1.M());
const int motherId = k1MotherId(channel, kaon, pion1, pion2);
const int charge = kaon.pdgCode() > 0 ? 1 : -1;
const auto* parent = parentTruth ? findK1ParentTruth(*parentTruth, motherId, channel, charge) : nullptr;
const bool truthRapidity = c.isUnlikeSign && parent && (!histogramOptions.cUseTruthRapidity || (parent->y >= candidateCuts.cK1MinRap && parent->y <= candidateCuts.cK1MaxRap));
if (truthRapidity) {
typeK1 = kaon.sign() > 0 ? BinType::kK1P_Rec : BinType::kK1N_Rec;
histos.fill(HIST("hInvmass_K1"), typeNormal, typeK1, multiplicity, c.k1.Pt(), c.k1.M());
histos.fill(HIST("k1invmass_MC"), c.k1.M());
const int parentCharge = parent->pdgCode > 0 ? 1 : -1;
histos.fill(HIST("MCReco/response"), c.k1.Pt(), parent->pt, channelBin, parentCharge, multiplicity);
}
histos.fill(HIST("QAMC/K1OA"), c.angle);
histos.fill(HIST("QAMC/K1PairAsym"), c.pairAsym);
histos.fill(HIST("QAMC/hInvmassK892_Rho"), c.mass13, c.secondary.M());
Expand Down Expand Up @@ -515,17 +530,18 @@ struct K1AnalysisMicro {
PROCESS_SWITCH(K1AnalysisMicro, processResoMicroTracks, "Process ResoMicroTracks", true);

void processMC(ResoMCCols::iterator const& collision,
ResoMCTracks const& resotracks)
ResoMCTracks const& resotracks, ResoMCParents const& resoParents)
{
if (!core.passesEventCuts(collision) || !core.passesMCEventCuts(collision)) {
return;
}
histos.fill(HIST("MCReco/collisions"), 0.5);
fillHistograms<true, false, false>(collision, resotracks, resotracks);
const auto parentTruth = makeK1ParentTruthMap(resoParents);
fillHistograms<true, false, false>(collision, resotracks, resotracks, &parentTruth);
}
PROCESS_SWITCH(K1AnalysisMicro, processMC, "Process Event for MC", false);

void processMCMicro(ResoMCCols::iterator const& collision, ResoMCMicroTracks const& tracks)
void processMCMicro(ResoMCCols::iterator const& collision, ResoMCMicroTracks const& tracks, ResoMCParents const& resoParents)
{
// The modular producer already selected these reconstructed collisions.
// Apply precisely the same reconstruction loop as the data baseline.
Expand All @@ -534,7 +550,8 @@ struct K1AnalysisMicro {
}
histos.fill(HIST("MCReco/collisions"), 0.5);
histos.fill(HIST("MCReco/microTracks"), 0.5, tracks.size());
fillHistograms<true, false, true>(collision, tracks, tracks);
const auto parentTruth = makeK1ParentTruthMap(resoParents);
fillHistograms<true, false, true>(collision, tracks, tracks, &parentTruth);
}
PROCESS_SWITCH(K1AnalysisMicro, processMCMicro, "Process reconstructed MC with micro v001 tables", false);

Expand All @@ -548,6 +565,9 @@ struct K1AnalysisMicro {
const int charge = part.pdgCode() > 0 ? 1 : -1;
histos.fill(HIST("MCGen/chargeChannel"), charge, static_cast<int>(channel));
histos.fill(HIST("MCGen/ptChannel"), static_cast<int>(channel), part.pt());
if (channel != K1TruthChannel::None) {
histos.fill(HIST("MCGen/response"), part.pt(), static_cast<int>(channel), charge, collision.cent());
}
});
}
PROCESS_SWITCH(K1AnalysisMicro, processMCTrue, "Process generated K1 in selected events with v001 parents", false);
Expand Down
22 changes: 13 additions & 9 deletions PWGLF/Tasks/Resonances/k1TrainingTable.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -182,7 +182,7 @@ struct K1TrainingTable {

template <bool IsMC, typename Collision, typename Kaon, typename Pion>
void writeK1MlCandidate(Collision const& collision, Kaon const& kaon, Pion const& samePion, Pion const& oppPion,
K1TruthChannel channel, uint16_t passBits)
K1TruthChannel channel, uint16_t passBits, const K1ParentTruthMap* parentTruth = nullptr)
{
using o2::analysis::k1ml::BuildStatus;
const auto canonical = o2::analysis::k1ml::canonicalizeUS(o2::analysis::k1ml::makeTrackSnapshot(kaon), o2::analysis::k1ml::makeTrackSnapshot(samePion), o2::analysis::k1ml::makeTrackSnapshot(oppPion));
Expand Down Expand Up @@ -215,14 +215,17 @@ struct K1TrainingTable {
if constexpr (IsMC) {
const bool matched = channel != K1TruthChannel::None;
int motherId = -1;
const K1ParentTruth* parent = nullptr;
if (matched) {
motherId = channel == K1TruthChannel::RhoK ? kaon.motherId() : std::abs(samePion.motherPDG()) == Pdg::kK1_1270Plus ? samePion.motherId()
: oppPion.motherId();
motherId = k1MotherId(channel, kaon, samePion, oppPion);
const int charge = kaon.pdgCode() > 0 ? 1 : -1;
parent = parentTruth ? findK1ParentTruth(*parentTruth, motherId, channel, charge) : nullptr;
}
// Immediate-mother information alone does not prove a physical UID.
k1MlTruth(row, matched ? uint8_t{1} : uint8_t{2}, static_cast<uint8_t>(channel),
matched ? kaon.sign() * Pdg::kK1_1270Plus : 0, static_cast<int64_t>(motherId),
std::numeric_limits<float>::quiet_NaN(), std::numeric_limits<float>::quiet_NaN());
parent ? parent->pdgCode : (matched ? kaon.sign() * Pdg::kK1_1270Plus : 0), static_cast<int64_t>(motherId),
parent ? static_cast<float>(parent->pt) : std::numeric_limits<float>::quiet_NaN(),
parent ? static_cast<float>(parent->y) : std::numeric_limits<float>::quiet_NaN());
} else {
k1MlTruth(row, uint8_t{0}, uint8_t{0}, 0, int64_t{-1},
std::numeric_limits<float>::quiet_NaN(), std::numeric_limits<float>::quiet_NaN());
Expand All @@ -244,7 +247,7 @@ struct K1TrainingTable {
}
PROCESS_SWITCH(K1TrainingTable, processResoMicroTracks, "Write K1 candidates from data micro v001 tables", true);

void processMCMicro(ResoMCCols::iterator const& collision, ResoMCMicroTracks const& tracks)
void processMCMicro(ResoMCCols::iterator const& collision, ResoMCMicroTracks const& tracks, ResoMCParents const& resoParents)
{
// The modular producer already selected these reconstructed collisions.
if (!core.passesEventCuts(collision) || !core.passesMCEventCuts(collision)) {
Expand All @@ -253,10 +256,11 @@ struct K1TrainingTable {
histos.fill(HIST("MCReco/collisions"), 0.5);
histos.fill(HIST("MCReco/microTracks"), 0.5, tracks.size());
writeK1MlEvent(collision);
const auto parentTruth = makeK1ParentTruthMap(resoParents);
core.forEachCandidate<true, false, true>(histos, collision, tracks, tracks, false, nullptr, nullptr,
[this](auto const& coll, auto const& kaon, auto const& samePion, auto const& oppPion,
K1TruthChannel channel, uint16_t passBits) {
writeK1MlCandidate<true>(coll, kaon, samePion, oppPion, channel, passBits);
[this, &parentTruth](auto const& coll, auto const& kaon, auto const& samePion, auto const& oppPion,
K1TruthChannel channel, uint16_t passBits) {
writeK1MlCandidate<true>(coll, kaon, samePion, oppPion, channel, passBits, &parentTruth);
});
}
PROCESS_SWITCH(K1TrainingTable, processMCMicro, "Write K1 candidates with truth from reconstructed MC micro v001 tables", false);
Expand Down
Loading
Loading