Skip to content

[PWGCF] Femto: Use mom. at PV for Sigma kstar calculation #18661

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
ariedel-cern merged 3 commits into from
Nov 9, 2026
Merged

[PWGCF] Femto: Use mom. at PV for Sigma kstar calculation #18661

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
c5d18d8
Fix wrong phi calculation for Sigmas
Nov 9, 2026
5b0cace
Include kinkBuilder.h changes
Nov 9, 2026
6eaed6a
Fix: Use epsilon and reco decay + add PhiAtPv in KinkFilters
Nov 9, 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
14 changes: 14 additions & 0 deletions PWGCF/Femto/Core/femtoUtils.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,8 +16,10 @@
#ifndef PWGCF_FEMTO_CORE_FEMTOUTILS_H_
#define PWGCF_FEMTO_CORE_FEMTOUTILS_H_

#include "Common/Core/RecoDecay.h"
#include "Common/Core/TableHelper.h"

#include <CommonConstants/MathConstants.h>
#include <CommonConstants/PhysicsConstants.h>
#include <Framework/InitContext.h>
#include <Framework/Logger.h>
Expand Down Expand Up @@ -173,6 +175,18 @@ concept HasEventShapeRow = requires(T row) {
template <typename T>
concept HasEventShape = HasEventShapeRow<std::decay_t<T>> || (requires { typename std::decay_t<T>::iterator; } && HasEventShapeRow<typename std::decay_t<T>::iterator>);

/// Azimuth of the momentum of a charged particle at the primary vertex, from its azimuth at the decay vertex.
/// The decay vertex (xDecay, yDecay) has to be given relative to the primary vertex. In a solenoidal field the transverse
/// trajectory is a circle through both points, and the tangents at the two ends of a chord make the same angle with it:
/// phiPV = 2*phiChord - phiDecay (independent of polarity).
inline float phiAtPrimaryVertex(float phiDecay, float xDecay, float yDecay)
{
if (std::abs(xDecay) < o2::constants::math::Epsilon && std::abs(yDecay) < o2::constants::math::Epsilon) {
return phiDecay;
}
return RecoDecay::constrainAngle(2.f * std::atan2(yDecay, xDecay) - phiDecay);
}

/// Recalculate the pT of a kink mother from its direction and the charged daughter momentum
inline float calcPtnew(float pxMother, float pyMother, float pzMother, float pxDaughter, float pyDaughter, float pzDaughter,
float massMother, float massChargedDaughter, float massNeutralDaughter)
Expand Down
28 changes: 24 additions & 4 deletions PWGCF/Femto/Core/kinkBuilder.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -53,6 +53,7 @@ struct ConfKinkFilters : o2::framework::ConfigurableGroup {
o2::framework::Configurable<float> etaMax{"etaMax", 10.f, "Maximum eta"};
o2::framework::Configurable<float> phiMin{"phiMin", 0.f, "Minimum phi"};
o2::framework::Configurable<float> phiMax{"phiMax", 1.f * o2::constants::math::TwoPI, "Maximum phi"};
o2::framework::Configurable<bool> kinkPhiAtPv{"kinkPhiAtPv", true, "Kinks and PCM SigmaPlus: store the azimuth of the Sigma momentum at the primary vertex"};
Comment thread
ariedel-cern marked this conversation as resolved.
Resolved
o2::framework::Configurable<float> massMinSigma{"massMinSigma", 1.1f, "Minimum mass for Sigma hypothesis"};
o2::framework::Configurable<float> massMaxSigma{"massMaxSigma", 1.3f, "Maximum mass for Sigma hypothesis"};
o2::framework::Configurable<float> massMinSigmaPlus{"massMinSigmaPlus", 1.1f, "Minimum mass for SigmaPlus hypothesis"};
Expand Down Expand Up @@ -234,6 +235,7 @@ enum KinkFilters {
kPhiMax,
kMassMin,
kMassMax,
kPhiAtPv,
kMlScoreMin,
kKinkFiltersMax
};
Expand All @@ -250,6 +252,7 @@ const std::unordered_map<KinkFilters, std::string> kinkFilterNames = {
{kPhiMax, "phiMax"},
{kMassMin, "massMin"},
{kMassMax, "massMax"},
{kPhiAtPv, "kinkPhiAtPv"},
{kMlScoreMin, "mlScoreMin"}};

/// \class KinkCuts
Expand All @@ -272,6 +275,7 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f
mEtaMax = filter.etaMax.value;
mPhiMin = filter.phiMin.value;
mPhiMax = filter.phiMax.value;
mKinkPhiAtPv = filter.kinkPhiAtPv.value;

if constexpr (modes::isEqual(kinkType, modes::Kink::kSigma)) {
mMassSigmaLowerLimit = filter.massMinSigma.value;
Expand Down Expand Up @@ -317,6 +321,7 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f
{kinkFilterNames.at(kPhiMax), mPhiMax},
{kinkFilterNames.at(kMassMin), mMassSigmaLowerLimit},
{kinkFilterNames.at(kMassMax), mMassSigmaUpperLimit},
{kinkFilterNames.at(kPhiAtPv), static_cast<float>(mKinkPhiAtPv)},
});
}
if constexpr (modes::isEqual(kinkType, modes::Kink::kSigmaPlus)) {
Expand All @@ -331,6 +336,7 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f
{kinkFilterNames.at(kPhiMax), mPhiMax},
{kinkFilterNames.at(kMassMin), mMassSigmaPlusLowerLimit},
{kinkFilterNames.at(kMassMax), mMassSigmaPlusUpperLimit},
{kinkFilterNames.at(kPhiAtPv), static_cast<float>(mKinkPhiAtPv)},
});
}
};
Expand Down Expand Up @@ -421,6 +427,9 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f

mKinkMotherEta = RecoDecay::eta(momMother);
mKinkMotherPhi = RecoDecay::phi(momMother);
if (mKinkPhiAtPv) {
mKinkMotherPhi = utils::phiAtPrimaryVertex(mKinkMotherPhi, kinkCand.xDecVtx(), kinkCand.yDecVtx());
}

// Recalculate pT using kinematic constraints of the decay channel
float ptRecalc = -999.f;
Expand Down Expand Up @@ -468,6 +477,8 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f
this->template fillFilter<FilterHistName>(kPhiMax, p);
pass &= p;

this->template fillFilter<FilterHistName>(kPhiAtPv, mKinkPhiAtPv);

if constexpr (modes::isEqual(kinkType, modes::Kink::kSigma)) {
float sigmaMass = kinkCand.mSigmaMinus();

Expand Down Expand Up @@ -517,6 +528,7 @@ class KinkSelection : public baseselection::BaseSelection<float, o2::analysis::f
float mEtaMax = 1.f;
float mPhiMin = 0.f;
float mPhiMax = o2::constants::math::TwoPI;
bool mKinkPhiAtPv = true;

// mother kinematic
float mKinkMotherPt = 0.f; // recalculated pT (falls back to the original pT if the recalculation fails)
Expand Down Expand Up @@ -827,6 +839,7 @@ class SigmaPlusPcmSelection : public baseselection::BaseSelection<float, o2::ana
mEtaMax = filter.etaMax.value;
mPhiMin = filter.phiMin.value;
mPhiMax = filter.phiMax.value;
mPhiAtPv = filter.kinkPhiAtPv.value;
mMassLowerLimit = filter.massMinSigmaPlusPcm.value;
mMassUpperLimit = filter.massMaxSigmaPlusPcm.value;
mMlScoreMin = filter.mlScoreMinSigmaPlusPcm.value;
Expand All @@ -852,17 +865,21 @@ class SigmaPlusPcmSelection : public baseselection::BaseSelection<float, o2::ana
{kinkFilterNames.at(kPhiMax), mPhiMax},
{kinkFilterNames.at(kMassMin), mMassLowerLimit},
{kinkFilterNames.at(kMassMax), mMassUpperLimit},
{kinkFilterNames.at(kPhiAtPv), static_cast<float>(mPhiAtPv)},
{kinkFilterNames.at(kMlScoreMin), mMlScoreMin},
});
};

template <typename T>
void computeKinematics(T const& sigmaPlusPcmCand)
template <typename T1, typename T2>
void computeKinematics(T1 const& sigmaPlusPcmCand, T2 const& col)
{
std::array<float, 3> momMother = {sigmaPlusPcmCand.pxSigmaPlus(), sigmaPlusPcmCand.pySigmaPlus(), sigmaPlusPcmCand.pzSigmaPlus()};
mPt = RecoDecay::pt(momMother);
mEta = RecoDecay::eta(momMother);
mPhi = RecoDecay::phi(momMother);
if (mPhiAtPv) {
mPhi = utils::phiAtPrimaryVertex(mPhi, sigmaPlusPcmCand.xDecVtx() - col.posX(), sigmaPlusPcmCand.yDecVtx() - col.posY());
}
mMass = sigmaPlusPcmCand.massSigmaPlus();
mMlScore = sigmaPlusPcmCand.mlScore();
}
Expand Down Expand Up @@ -897,6 +914,8 @@ class SigmaPlusPcmSelection : public baseselection::BaseSelection<float, o2::ana
this->template fillFilter<FilterHistName>(kPhiMax, p);
pass &= p;

this->template fillFilter<FilterHistName>(kPhiAtPv, mPhiAtPv);

p = mMass > mMassLowerLimit;
this->template fillFilter<FilterHistName>(kMassMin, p);
pass &= p;
Expand Down Expand Up @@ -951,6 +970,7 @@ class SigmaPlusPcmSelection : public baseselection::BaseSelection<float, o2::ana
float mEtaMin = -1.f;
float mEtaMax = 1.f;
float mPhiMin = 0.f;
bool mPhiAtPv = true;
float mPhiMax = o2::constants::math::TwoPI;

// mother kinematic
Expand Down Expand Up @@ -1015,7 +1035,7 @@ class SigmaPlusPcmBuilder
int64_t daughterIndex = 0;

for (const auto& sigmaPlusPcm : sigmaPlusPcms) {
mSigmaPlusPcmSelection.computeKinematics(sigmaPlusPcm);
mSigmaPlusPcmSelection.computeKinematics(sigmaPlusPcm, col);
if (!mSigmaPlusPcmSelection.checkFilters(sigmaPlusPcm)) {
continue;
}
Expand All @@ -1041,7 +1061,7 @@ class SigmaPlusPcmBuilder
int64_t daughterIndex = 0;

for (const auto& sigmaPlusPcm : sigmaPlusPcms) {
mSigmaPlusPcmSelection.computeKinematics(sigmaPlusPcm);
mSigmaPlusPcmSelection.computeKinematics(sigmaPlusPcm, col);
if (!mSigmaPlusPcmSelection.checkFilters(sigmaPlusPcm)) {
continue;
}
Expand Down
64 changes: 61 additions & 3 deletions PWGCF/Femto/TableProducer/femtoProducer.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -270,11 +270,13 @@ struct FemtoProducer {
}
if ((builders.sigmaBuilder.fillAnyTable() || builders.sigmaPlusBuilder.fillAnyTable()) &&
(!doprocessTracksKinksRun3pp && !doprocessTracksV0sCascadesKinksRun3pp && !doprocessTracksV0sCascadesKinksPhotonsRun3pp &&
!doprocessTracksKinksRun3ppMc && !doprocessTracksV0sKinksRun3ppMc)) {
!doprocessTracksKinksRun3ppMc && !doprocessTracksV0sKinksRun3ppMc &&
!doprocessTracksKinksSigmaPlusPcmsRun3pp && !doprocessTracksKinksSigmaPlusPcmsRun3ppMc)) {
LOG(fatal) << "At least one kink table is enabled, but wrong process function is enabled. Breaking...";
}
if (builders.sigmaPlusPcmBuilder.fillAnyTable() &&
(!doprocessTracksSigmaPlusPcmsRun3pp && !doprocessTracksSigmaPlusPcmsRun3ppMc)) {
(!doprocessTracksSigmaPlusPcmsRun3pp && !doprocessTracksSigmaPlusPcmsRun3ppMc &&
!doprocessTracksKinksSigmaPlusPcmsRun3pp && !doprocessTracksKinksSigmaPlusPcmsRun3ppMc)) {
LOG(fatal) << "At least one sigma plus (PCM) table is enabled, but wrong process function is enabled. Breaking...";
}
if ((builders.d0Builder.fillAnyTable() || builders.d0barBuilder.fillAnyTable()) &&
Expand All @@ -292,7 +294,7 @@ struct FemtoProducer {
!doprocessTracksV0sRun3ppMc && !doprocessTracksV0sRun3PbPbMc &&
!doprocessTracksV0sCascadesRun3ppMc && !doprocessTracksV0sCascadesRun3PbPbMc &&
!doprocessTracksKinksRun3ppMc && !doprocessTracksV0sKinksRun3ppMc &&
!doprocessTracksSigmaPlusPcmsRun3ppMc &&
!doprocessTracksSigmaPlusPcmsRun3ppMc && !doprocessTracksKinksSigmaPlusPcmsRun3ppMc &&
!doprocessTracksD0sRun3ppMc && !doprocessTracksD0sRun3PbPbMc &&
!doprocessTracksLcsRun3ppMc && !doprocessTracksLcsRun3PbPbMc &&
!doprocessMcOnly)) {
Expand Down Expand Up @@ -379,6 +381,7 @@ struct FemtoProducer {
static_cast<int>(doprocessTracksV0sCascadesKinksRun3pp) +
static_cast<int>(doprocessTracksV0sCascadesKinksPhotonsRun3pp) +
static_cast<int>(doprocessTracksSigmaPlusPcmsRun3pp) +
static_cast<int>(doprocessTracksKinksSigmaPlusPcmsRun3pp) +
static_cast<int>(doprocessTracksD0sRun3pp) +
static_cast<int>(doprocessTracksD0sRun3PbPb) +
static_cast<int>(doprocessTracksRun3ppMc) +
Expand All @@ -390,6 +393,7 @@ struct FemtoProducer {
static_cast<int>(doprocessTracksKinksRun3ppMc) +
static_cast<int>(doprocessTracksV0sKinksRun3ppMc) +
static_cast<int>(doprocessTracksSigmaPlusPcmsRun3ppMc) +
static_cast<int>(doprocessTracksKinksSigmaPlusPcmsRun3ppMc) +
static_cast<int>(doprocessTracksD0sRun3ppMc) +
static_cast<int>(doprocessTracksD0sRun3PbPbMc) +
static_cast<int>(doprocessTracksLcsRun3pp) +
Expand Down Expand Up @@ -833,6 +837,29 @@ struct FemtoProducer {
}
PROCESS_SWITCH(FemtoProducer, processTracksSigmaPlusPcmsRun3pp, "Provide tracks and sigma plus reconstructed with a PCM photon", false);

void processTracksKinksSigmaPlusPcmsRun3pp(rawinputs::Run3PpCollisions const& cols,
o2::aod::BCsWithTimestamps const& bcs,
rawinputs::Run3FullPidTracks const& tracks,
rawinputs::Run3Kinks const& kinks,
rawinputs::Run3SigmaPlusPcms const& sigmaPlusPcms)
{
builders.trackBuilder.reset(tracks);
for (auto const& col : cols) {
auto tracksThisCol = tracks.sliceBy(preslices.perColTracks, col.globalIndex());
if (!processCollisions<modes::System::kPP_Run3>(col, bcs, tracksThisCol)) {
continue;
}
auto tracksWithItsPid = o2::soa::Attach<rawinputs::Run3FullPidTracks, o2::aod::pidits::ITSNSigmaEl, o2::aod::pidits::ITSNSigmaPi, o2::aod::pidits::ITSNSigmaKa,
o2::aod::pidits::ITSNSigmaPr, o2::aod::pidits::ITSNSigmaDe, o2::aod::pidits::ITSNSigmaTr, o2::aod::pidits::ITSNSigmaHe>(tracksThisCol);
processTracks<modes::System::kPP_Run3>(col, tracksWithItsPid);
auto kinksThisCol = kinks.sliceBy(preslices.perColRecoKinks, col.globalIndex());
processKinks<modes::System::kPP_Run3>(col, tracksThisCol, kinksThisCol);
auto sigmaPlusPcmsThisCol = sigmaPlusPcms.sliceBy(preslices.perColSigmaPlusPcms, col.globalIndex());
processSigmaPlusPcms<modes::System::kPP_Run3>(col, tracksThisCol, sigmaPlusPcmsThisCol);
}
}
PROCESS_SWITCH(FemtoProducer, processTracksKinksSigmaPlusPcmsRun3pp, "Provide tracks, kinks and sigma plus reconstructed with a PCM photon", false);

void processTracksV0sCascadesKinksRun3pp(rawinputs::Run3PpCollisions const& cols,
o2::aod::BCsWithTimestamps const& bcs,
rawinputs::Run3FullPidTracks const& tracks,
Expand Down Expand Up @@ -1168,6 +1195,37 @@ struct FemtoProducer {
}
PROCESS_SWITCH(FemtoProducer, processTracksSigmaPlusPcmsRun3ppMc, "Provide reconstructed and generated tracks and sigma plus reconstructed with a PCM photon", false);

void processTracksKinksSigmaPlusPcmsRun3ppMc(rawinputs::Run3PpMcGenCollisions const& mcCols,
rawinputs::Run3PpMcRecoCollisions const& cols,
o2::aod::BCsWithTimestamps const& bcs,
rawinputs::Run3McRecoTracks const& tracks,
rawinputs::Run3Kinks const& kinks,
rawinputs::Run3RecoSigmaPlusPcms const& sigmaPlusPcms,
rawinputs::Run3McGenParticles const& mcParticles)
{
builders.mcBuilder.reset(mcCols, mcParticles);
builders.trackBuilder.reset(tracks);
if (builders.mcBuilder.isPassThrough()) {
builders.mcBuilder.fillMcPassThrough<modes::System::kPP_Run3_MC>(mcCols, mcParticles, preslices.perMcCollision, mcProducts);
}

for (const auto& col : cols) {
auto tracksThisCol = tracks.sliceBy(preslices.perColRecoTracks, col.globalIndex());
if (!processMcCollisions<modes::System::kPP_Run3_MC>(col, mcCols, bcs, tracksThisCol)) {
continue;
}
auto tracksWithItsPid = o2::soa::Attach<rawinputs::Run3McRecoTracks, o2::aod::pidits::ITSNSigmaEl, o2::aod::pidits::ITSNSigmaPi, o2::aod::pidits::ITSNSigmaKa,
o2::aod::pidits::ITSNSigmaPr, o2::aod::pidits::ITSNSigmaDe, o2::aod::pidits::ITSNSigmaTr, o2::aod::pidits::ITSNSigmaHe>(tracksThisCol);
processMcTracks<modes::System::kPP_Run3_MC>(col, mcCols, tracksThisCol, tracksWithItsPid, mcParticles);

auto kinksThisCol = kinks.sliceBy(preslices.perColRecoKinks, col.globalIndex());
processMcKinks<modes::System::kPP_Run3_MC>(col, mcCols, tracks, kinksThisCol, mcParticles);
auto sigmaPlusPcmsThisCol = sigmaPlusPcms.sliceBy(preslices.perColRecoSigmaPlusPcms, col.globalIndex());
processMcSigmaPlusPcms<modes::System::kPP_Run3_MC>(col, mcCols, tracks, sigmaPlusPcmsThisCol, mcParticles);
}
}
PROCESS_SWITCH(FemtoProducer, processTracksKinksSigmaPlusPcmsRun3ppMc, "Provide reconstructed and generated tracks, kinks and sigma plus reconstructed with a PCM photon", false);

void processTracksV0sKinksRun3ppMc(rawinputs::Run3PpMcGenCollisions const& mcCols,
rawinputs::Run3PpMcRecoCollisions const& cols,
o2::aod::BCsWithTimestamps const& bcs,
Expand Down
Loading