Skip to content

[PWGHF] Add proper lifetime into THnSparse at taskLc #12902

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

Merged
mfaggin merged 8 commits into from
Sep 19, 2025
Merged

[PWGHF] Add proper lifetime into THnSparse at taskLc #12902

Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
1135b0d
add proper lifetime to THnSparse
lubynets Sep 5, 2025
cc0075b
bugfix: do not write KF minv into non-KF tree
lubynets Sep 5, 2025
54a5a51
get rid of magic numbers in Ml classes enumeration
lubynets Sep 8, 2025
b6ded05
Revert "bugfix: do not write KF minv into non-KF tree"
lubynets Sep 8, 2025
8c3b8ad
make thnAxisProperLifetime optional
lubynets Sep 8, 2025
8d6e98c
linter: Double_t -> double
lubynets Sep 8, 2025
bb3a0a5
swap occupancy and lifetime axes for backward compatibility
lubynets Sep 8, 2025
270def9
add proper lifetime also to axesStd
lubynets Sep 8, 2025
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
170 changes: 108 additions & 62 deletions PWGHF/D2H/Tasks/taskLc.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,17 @@ struct HfTaskLc {
Configurable<bool> fillTHn{"fillTHn", false, "fill THn"};
Configurable<bool> storeOccupancy{"storeOccupancy", true, "Flag to store occupancy information"};
Configurable<int> occEstimator{"occEstimator", 2, "Occupancy estimation (None: 0, ITS: 1, FT0C: 2)"};
Configurable<bool> storeProperLifetime{"storeProperLifetime", false, "Flag to store proper lifetime"};

constexpr static float CtToProperLifetimePs = 1.f / o2::constants::physics::LightSpeedCm2PS;
constexpr static float NanoToPico = 1000.f;

enum MlClasses : int {
MlClassBackground = 0,
MlClassPrompt,
MlClassNonPrompt,
NumberOfMlClasses
};

HfHelper hfHelper;
SliceCache cache;
Expand Down Expand Up @@ -100,6 +111,7 @@ struct HfTaskLc {
ConfigurableAxis thnConfigAxisGenPtB{"thnConfigAxisGenPtB", {1000, 0, 100}, "Gen Pt B"};
ConfigurableAxis thnConfigAxisNumPvContr{"thnConfigAxisNumPvContr", {200, -0.5, 199.5}, "Number of PV contributors"};
ConfigurableAxis thnConfigAxisOccupancy{"thnConfigAxisOccupancy", {14, 0, 14000}, "axis for centrality"};
ConfigurableAxis thnConfigAxisProperLifetime{"thnConfigAxisProperLifetime", {200, 0, 2}, "Proper lifetime, ps"};

HistogramRegistry registry{
"registry",
Expand Down Expand Up @@ -317,6 +329,7 @@ struct HfTaskLc {
const AxisSpec thnAxisPtB{thnConfigAxisGenPtB, "#it{p}_{T}^{B} (GeV/#it{c})"};
const AxisSpec thnAxisTracklets{thnConfigAxisNumPvContr, "Number of PV contributors"};
const AxisSpec thnAxisOccupancy{thnConfigAxisOccupancy, "Occupancy"};
const AxisSpec thnAxisProperLifetime{thnConfigAxisProperLifetime, "T_{proper} (ps)"};

bool isDataWithMl = doprocessDataWithMl || doprocessDataWithMlWithFT0C || doprocessDataWithMlWithFT0M;
bool isMcWithMl = doprocessMcWithMl || doprocessMcWithMlWithFT0C || doprocessMcWithMlWithFT0M;
Expand All @@ -342,12 +355,18 @@ struct HfTaskLc {
}

if (storeOccupancy) {
if (!axesWithBdt.empty())
axesWithBdt.push_back(thnAxisOccupancy);
if (!axesStd.empty())
axesStd.push_back(thnAxisOccupancy);
if (!axesGen.empty())
axesGen.push_back(thnAxisOccupancy);
for (const auto& axes : std::array<std::vector<AxisSpec>*, 3>{&axesWithBdt, &axesStd, &axesGen}) {
if (!axes->empty()) {
axes->push_back(thnAxisOccupancy);
}
}
}
if (storeProperLifetime) {
for (const auto& axes : std::array<std::vector<AxisSpec>*, 3>{&axesWithBdt, &axesStd, &axesGen}) {
if (!axes->empty()) {
axes->push_back(thnAxisProperLifetime);
}
}
}

if (isDataWithMl) {
Expand Down Expand Up @@ -542,55 +561,63 @@ struct HfTaskLc {
}
double massLc(-1);
double outputBkg(-1), outputPrompt(-1), outputFD(-1);
const float properLifetime = hfHelper.ctLc(candidate) * CtToProperLifetimePs;
if ((candidate.isSelLcToPKPi() >= selectionFlagLc) && pdgCodeProng0 == kProton) {
massLc = hfHelper.invMassLcToPKPi(candidate);

if constexpr (fillMl) {
if (candidate.mlProbLcToPKPi().size() == 3) {
outputBkg = candidate.mlProbLcToPKPi()[0]; /// bkg score
outputPrompt = candidate.mlProbLcToPKPi()[1]; /// prompt score
outputFD = candidate.mlProbLcToPKPi()[2]; /// non-prompt score
if (candidate.mlProbLcToPKPi().size() == NumberOfMlClasses) {
outputBkg = candidate.mlProbLcToPKPi()[MlClassBackground]; /// bkg score
outputPrompt = candidate.mlProbLcToPKPi()[MlClassPrompt]; /// prompt score
outputFD = candidate.mlProbLcToPKPi()[MlClassNonPrompt]; /// non-prompt score
}
/// Fill the ML outputScores and variables of candidate
std::vector<double> valuesToFill{massLc, pt, cent, outputBkg, outputPrompt, outputFD, static_cast<double>(numPvContributors), ptRecB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, ptRecB, originType, occ);
} else {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, ptRecB, originType);
valuesToFill.push_back(occ);
}

if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
} else {

std::vector<double> valuesToFill{massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, static_cast<double>(numPvContributors), ptRecB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, ptRecB, originType, occ);

} else {

registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, ptRecB, originType);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(valuesToFill.data());
}
}
if ((candidate.isSelLcToPiKP() >= selectionFlagLc) && pdgCodeProng0 == kPiPlus) {
massLc = hfHelper.invMassLcToPiKP(candidate);

if constexpr (fillMl) {
if (candidate.mlProbLcToPiKP().size() == 3) {
outputBkg = candidate.mlProbLcToPiKP()[0]; /// bkg score
outputPrompt = candidate.mlProbLcToPiKP()[1]; /// prompt score
outputFD = candidate.mlProbLcToPiKP()[2]; /// non-prompt score
if (candidate.mlProbLcToPiKP().size() == NumberOfMlClasses) {
outputBkg = candidate.mlProbLcToPiKP()[MlClassBackground]; /// bkg score
outputPrompt = candidate.mlProbLcToPiKP()[MlClassPrompt]; /// prompt score
outputFD = candidate.mlProbLcToPiKP()[MlClassNonPrompt]; /// non-prompt score
}
/// Fill the ML outputScores and variables of candidate (todo: add multiplicity)
std::vector<double> valuesToFill{massLc, pt, cent, outputBkg, outputPrompt, outputFD, static_cast<double>(numPvContributors), ptRecB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, ptRecB, originType, occ);

} else {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, ptRecB, originType);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
} else {
std::vector<double> valuesToFill{massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, static_cast<double>(numPvContributors), ptRecB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, ptRecB, originType, occ);
} else {
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, ptRecB, originType);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(valuesToFill.data());
}
}
}
Expand Down Expand Up @@ -624,6 +651,11 @@ struct HfTaskLc {
occ = o2::hf_occupancy::getOccupancyGenColl(recoCollsPerMcColl, occEstimator);
}

const auto mcDaughter0 = particle.template daughters_as<soa::Join<aod::McParticles, aod::HfCand3ProngMcGen>>().begin();
const float p2m = particle.p() / o2::constants::physics::MassLambdaCPlus;
const float gamma = std::sqrt(1 + p2m * p2m); // mother's particle Lorentz factor
const float properLifetime = mcDaughter0.vt() * NanoToPico / gamma; // from ns to ps * from lab time to proper time

registry.fill(HIST("MC/generated/signal/hPtGen"), ptGen);
registry.fill(HIST("MC/generated/signal/hEtaGen"), particle.eta());
registry.fill(HIST("MC/generated/signal/hYGen"), yGen);
Expand All @@ -634,12 +666,14 @@ struct HfTaskLc {

if (particle.originMcGen() == RecoDecay::OriginType::Prompt) {
if (fillTHn) {
std::vector<double> valuesToFill{ptGen, cent, yGen, static_cast<double>(numPvContributors), ptGenB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(ptGen, cent, yGen, numPvContributors, ptGenB, originType, occ);

} else {
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(ptGen, cent, yGen, numPvContributors, ptGenB, originType);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(valuesToFill.data());
}
registry.fill(HIST("MC/generated/prompt/hPtGenPrompt"), ptGen);
registry.fill(HIST("MC/generated/prompt/hEtaGenPrompt"), particle.eta());
Expand All @@ -652,12 +686,14 @@ struct HfTaskLc {
if (particle.originMcGen() == RecoDecay::OriginType::NonPrompt) {
ptGenB = mcParticles.rawIteratorAt(particle.idxBhadMotherPart()).pt();
if (fillTHn) {
std::vector<double> valuesToFill{ptGen, cent, yGen, static_cast<double>(numPvContributors), ptGenB, static_cast<double>(originType)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(ptGen, cent, yGen, numPvContributors, ptGenB, originType, occ);

} else {
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(ptGen, cent, yGen, numPvContributors, ptGenB, originType);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsGen"))->Fill(valuesToFill.data());
}
registry.fill(HIST("MC/generated/nonprompt/hPtGenNonPrompt"), ptGen);
registry.fill(HIST("MC/generated/nonprompt/hEtaGenNonPrompt"), particle.eta());
Expand Down Expand Up @@ -748,53 +784,63 @@ struct HfTaskLc {
}
double massLc(-1);
double outputBkg(-1), outputPrompt(-1), outputFD(-1);
const float properLifetime = hfHelper.ctLc(candidate) * CtToProperLifetimePs;
if (candidate.isSelLcToPKPi() >= selectionFlagLc) {
massLc = hfHelper.invMassLcToPKPi(candidate);

if constexpr (fillMl) {
if (candidate.mlProbLcToPKPi().size() == 3) {
outputBkg = candidate.mlProbLcToPKPi()[0]; /// bkg score
outputPrompt = candidate.mlProbLcToPKPi()[1]; /// prompt score
outputFD = candidate.mlProbLcToPKPi()[2]; /// non-prompt score
if (candidate.mlProbLcToPKPi().size() == NumberOfMlClasses) {
outputBkg = candidate.mlProbLcToPKPi()[MlClassBackground]; /// bkg score
outputPrompt = candidate.mlProbLcToPKPi()[MlClassPrompt]; /// prompt score
outputFD = candidate.mlProbLcToPKPi()[MlClassNonPrompt]; /// non-prompt score
}
/// Fill the ML outputScores and variables of candidate
std::vector<double> valuesToFill{massLc, pt, cent, outputBkg, outputPrompt, outputFD, static_cast<double>(numPvContributors)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, occ);

} else {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
} else {
std::vector<double> valuesToFill{massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, static_cast<double>(numPvContributors)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {

registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, occ);
} else {

registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(valuesToFill.data());
}
}
if (candidate.isSelLcToPiKP() >= selectionFlagLc) {
massLc = hfHelper.invMassLcToPiKP(candidate);

if constexpr (fillMl) {
if (candidate.mlProbLcToPiKP().size() == 3) {
outputBkg = candidate.mlProbLcToPiKP()[0]; /// bkg score
outputPrompt = candidate.mlProbLcToPiKP()[1]; /// prompt score
outputFD = candidate.mlProbLcToPiKP()[2]; /// non-prompt score
if (candidate.mlProbLcToPiKP().size() == NumberOfMlClasses) {
outputBkg = candidate.mlProbLcToPiKP()[MlClassBackground]; /// bkg score
outputPrompt = candidate.mlProbLcToPiKP()[MlClassPrompt]; /// prompt score
outputFD = candidate.mlProbLcToPiKP()[MlClassNonPrompt]; /// non-prompt score
}
/// Fill the ML outputScores and variables of candidate
std::vector<double> valuesToFill{massLc, pt, cent, outputBkg, outputPrompt, outputFD, static_cast<double>(numPvContributors)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors, occ);
} else {
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(massLc, pt, cent, outputBkg, outputPrompt, outputFD, numPvContributors);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVarsWithBdt"))->Fill(valuesToFill.data());
} else {
std::vector<double> valuesToFill{massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, static_cast<double>(numPvContributors)};
if (storeOccupancy && occEstimator != o2::hf_occupancy::OccupancyEstimator::None) {
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors, occ);
} else {
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(massLc, pt, cent, ptProng0, ptProng1, ptProng2, chi2PCA, decayLength, cpa, numPvContributors);
valuesToFill.push_back(occ);
}
if (storeProperLifetime) {
valuesToFill.push_back(properLifetime);
}
registry.get<THnSparse>(HIST("hnLcVars"))->Fill(valuesToFill.data());
}
}
}
Expand Down
Loading