Skip to content

[PWGCF] Add lite colshape table #17410

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 4 commits into from
Aug 11, 2026
Merged

[PWGCF] Add lite colshape table #17410

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
4dd5374
Feat: add FLiteColShapes table
ariedel-cern Aug 11, 2026
9b2f986
Feat: cleanup derived to derived builder/producer
ariedel-cern Aug 11, 2026
0484b3d
Feat: use common daughterRow function
ariedel-cern Aug 11, 2026
ce4dc92
Feat: add FLiteColShapes in collisionBuilder
ariedel-cern Aug 11, 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
56 changes: 31 additions & 25 deletions PWGCF/Femto/Core/cascadeBuilder.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -793,15 +793,15 @@ struct CascadeBuilderDerivedToDerivedProducts : o2::framework::ProducesGroup {
class CascadeBuilderDerivedToDerived
{
public:
CascadeBuilderDerivedToDerived() = default;
~CascadeBuilderDerivedToDerived() = default;

template <typename T>
void init(T& config)
{
mLimitXi = config.limitXi.value;
mLimitOmega = config.limitOmega.value;

if (mLimitXi < 0 || mLimitOmega < 0) {
LOG(fatal) << "Cascade limits must be non-negative (got " << mLimitXi << " and " << mLimitOmega << "). Breaking...";
}
if (mLimitXi == 0 && mLimitOmega == 0) {
LOG(fatal) << "Both xi limit and omega limit are 0. Breaking...";
}
Expand All @@ -810,34 +810,40 @@ class CascadeBuilderDerivedToDerived
template <typename T1, typename T2, typename T3, typename T4>
bool collisionHasTooFewXis(T1 const& col, T2 const& /*xiTable*/, T3& partitionXi, T4& cache)
{
if (mLimitXi == 0) { // xis disabled, cannot reject on them
return false;
}
auto xiSlice = partitionXi->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);
return xiSlice.size() < mLimitXi;
return xiSlice.size() < static_cast<int64_t>(mLimitXi);
}

template <typename T1, typename T2, typename T3, typename T4>
bool collisionHasTooFewOmegas(T1 const& col, T2 const& /*omegaTable*/, T3& partitionOmega, T4& cache)
{
if (mLimitOmega == 0) { // omegas disabled, cannot reject on them
return false;
}
auto omegaSlice = partitionOmega->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);
return omegaSlice.size() < mLimitOmega;
return omegaSlice.size() < static_cast<int64_t>(mLimitOmega);
}

template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
void processXis(T1 const& col, T2 const& /*xiTable*/, T3 const& oldTrackTable, T4& partitionXi, T5& trackBuilder, T6& cache, T7& newXiTable, T8& newTrackTable, T9& newCollisionTable)
{
if (mLimitXi == 0) { // xis disabled
return;
}

auto xiSlice = partitionXi->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);

for (auto const& xi : xiSlice) {
auto bachelor = oldTrackTable.rawIteratorAt(utils::daughterRow(xi.bachelorId(), oldTrackTable));
auto posDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(xi.posDauId(), oldTrackTable));
auto negDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(xi.negDauId(), oldTrackTable));

// auto bachelor = xi.template bachelor_as<T3>();
// auto posDaughter = xi.template posDau_as<T3>();
// auto negDaughter = xi.template negDau_as<T3>();
auto bachelor = oldTrackTable.rawIteratorAt(xi.bachelorId() - oldTrackTable.offset());
auto posDaughter = oldTrackTable.rawIteratorAt(xi.posDauId() - oldTrackTable.offset());
auto negDaughter = oldTrackTable.rawIteratorAt(xi.negDauId() - oldTrackTable.offset());

int bachelorIndex = trackBuilder.getDaughterIndex(bachelor, newTrackTable, newCollisionTable);
int posDaughterIndex = trackBuilder.getDaughterIndex(posDaughter, newTrackTable, newCollisionTable);
int negDaughterIndex = trackBuilder.getDaughterIndex(negDaughter, newTrackTable, newCollisionTable);
int64_t bachelorIndex = trackBuilder.getDaughterIndex(bachelor, newTrackTable, newCollisionTable);
int64_t posDaughterIndex = trackBuilder.getDaughterIndex(posDaughter, newTrackTable, newCollisionTable);
int64_t negDaughterIndex = trackBuilder.getDaughterIndex(negDaughter, newTrackTable, newCollisionTable);

newXiTable.producedXis(newCollisionTable.producedCollision.lastIndex(),
xi.signedPt(),
Expand All @@ -854,20 +860,20 @@ class CascadeBuilderDerivedToDerived
template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
void processOmegas(T1 const& col, T2 const& /*omegaTable*/, T3 const& oldTrackTable, T4& partitionOmega, T5& trackBuilder, T6& cache, T7& newOmegaTable, T8& newTrackTable, T9& newCollisionTable)
{
if (mLimitOmega == 0) { // omegas disabled
return;
}

auto omegaSlice = partitionOmega->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);

for (auto const& omega : omegaSlice) {
auto bachelor = oldTrackTable.rawIteratorAt(utils::daughterRow(omega.bachelorId(), oldTrackTable));
auto posDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(omega.posDauId(), oldTrackTable));
auto negDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(omega.negDauId(), oldTrackTable));

// auto bachelor = omega.template bachelor_as<T3>();
// auto posDaughter = omega.template posDau_as<T3>();
// auto negDaughter = omega.template negDau_as<T3>();
auto bachelor = oldTrackTable.rawIteratorAt(omega.bachelorId() - oldTrackTable.offset());
auto posDaughter = oldTrackTable.rawIteratorAt(omega.posDauId() - oldTrackTable.offset());
auto negDaughter = oldTrackTable.rawIteratorAt(omega.negDauId() - oldTrackTable.offset());

int bachelorIndex = trackBuilder.getDaughterIndex(bachelor, newTrackTable, newCollisionTable);
int posDaughterIndex = trackBuilder.getDaughterIndex(posDaughter, newTrackTable, newCollisionTable);
int negDaughterIndex = trackBuilder.getDaughterIndex(negDaughter, newTrackTable, newCollisionTable);
int64_t bachelorIndex = trackBuilder.getDaughterIndex(bachelor, newTrackTable, newCollisionTable);
int64_t posDaughterIndex = trackBuilder.getDaughterIndex(posDaughter, newTrackTable, newCollisionTable);
int64_t negDaughterIndex = trackBuilder.getDaughterIndex(negDaughter, newTrackTable, newCollisionTable);

newOmegaTable.producedOmegas(newCollisionTable.producedCollision.lastIndex(),
omega.signedPt(),
Expand Down
38 changes: 26 additions & 12 deletions PWGCF/Femto/Core/collisionBuilder.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -291,7 +291,7 @@ class CollisionSelection : public baseselection::BaseSelection<float, o2::analys
{collisionFilterNames.at(kFilterSphericityMax), mSphericityMax},
{collisionFilterNames.at(kFilterRctFlags), mUseRctFlags ? 1.f : 0.f},
});
};
}

/// \brief Initialize the Zorro trigger machinery for a new run. No-op if no triggers configured.
template <typename T1, typename T2>
Expand Down Expand Up @@ -320,7 +320,7 @@ class CollisionSelection : public baseselection::BaseSelection<float, o2::analys
mMagField = MagField;
}

float getMagneticField()
[[nodiscard]] int getMagneticField()
{
return mMagField;
}
Expand Down Expand Up @@ -491,14 +491,14 @@ class CollisionSelection : public baseselection::BaseSelection<float, o2::analys
}

this->assembleBitmask<SelectionHistName>();
};
}

protected:
template <typename T>
float computeSphericity(T const& tracks)
{
int minNumberTracks = 2;
double maxSphericity = 2.f;
const int64_t minNumberTracks = 2;
const double maxSphericity = 2.f;
if (tracks.size() <= minNumberTracks) {
return maxSphericity;
}
Expand Down Expand Up @@ -529,7 +529,7 @@ class CollisionSelection : public baseselection::BaseSelection<float, o2::analys

// filter cuts
float mVtxZMin = -12.f;
float mVtxZMax = -12.f;
float mVtxZMax = 12.f;
float mSphericityMin = 0.f;
float mSphericityMax = 1.f;
float mMagFieldMin = -5.f;
Expand All @@ -539,7 +539,7 @@ class CollisionSelection : public baseselection::BaseSelection<float, o2::analys
float mCentMin = 0.f;
float mCentMax = 100.f;

int mMagField = 0.f;
int mMagField = 0;
float mSphericity = 0.f;
float mCentrality = 0.f;
float mMultiplicity = 0.f;
Expand Down Expand Up @@ -568,6 +568,7 @@ struct CollisionBuilderProducts : o2::framework::ProducesGroup {
o2::framework::Produces<o2::aod::FColMults> producedMultiplicityEstimators;
o2::framework::Produces<o2::aod::FColCents> producedCentralityEstimators;
o2::framework::Produces<o2::aod::FColShapes> producedShapes;
o2::framework::Produces<o2::aod::FLiteColShapes> producedLiteShapes;
};

struct ConfCollisionTables : o2::framework::ConfigurableGroup {
Expand All @@ -580,6 +581,7 @@ struct ConfCollisionTables : o2::framework::ConfigurableGroup {
o2::framework::Configurable<int> produceMults{"produceMults", -1, "Produce Multiplicities (-1: auto; 0 off; 1 on)"};
o2::framework::Configurable<int> produceCents{"produceCents", -1, "Produce Centralities (-1: auto; 0 off; 1 on)"};
o2::framework::Configurable<int> produceShapes{"produceShapes", -1, "Produce Event shape variables (-1: auto; 0 off; 1 on)"};
o2::framework::Configurable<int> produceLiteShapes{"produceLiteShapes", -1, "Produce Lite Event shape variables (-1: auto; 0 off; 1 on)"};
};

template <auto& SelectionHistName, auto& FilterHistName>
Expand All @@ -606,16 +608,22 @@ class CollisionBuilder
mProducedMultiplicities = utils::enableTable("FColMults_001", confTable.produceMults.value, initContext);
mProducedCentralities = utils::enableTable("FColCents_001", confTable.produceCents.value, initContext);
mProducedShapes = utils::enableTable("FColShapes_001", confTable.produceShapes.value, initContext);
mProducedLiteShapes = utils::enableTable("FLiteColShapes_001", confTable.produceLiteShapes.value, initContext);

if (mProducedCollisions && mProducedLiteCollisions) {
LOG(fatal) << "FCols and FLiteCols are mutually exclusive -- enable only one. "
<< "FLiteCols is meant to only replace FCols at the producer stage (for better compression in derived data); "
<< "use the dedicated converter task to reconstruct FCols from FLiteCols downstream.";
}
if (mProducedShapes && mProducedLiteShapes) {
LOG(fatal) << "FColShapes and FLiteColShapes are mutually exclusive -- enable only one. "
<< "FLiteColShapes is meant to only replace FColShapes at the producer stage (for better compression in derived data); "
<< "use the dedicated converter task to reconstruct FColShapes from FLiteColShapes downstream.";
}

if (mProducedCollisions || mProducedLiteCollisions || mProducedCollisionMasks ||
mProducedPositions || mProducedSphericities || mProducedMultiplicities ||
mProducedCentralities) {
mProducedCentralities || mProducedShapes || mProducedLiteShapes) {
mFillAnyTable = true;
} else {
LOG(info) << "No tables configured, Selection object will not be configured...";
Expand All @@ -634,15 +642,14 @@ class CollisionBuilder
if (mRunNumber != bc.runNumber()) {
mRunNumber = bc.runNumber();
if (mMagFieldForced == 0) {
static o2::parameters::GRPMagField* grpo = nullptr;
LOG(info) << "Get magentic field with Path: " << mGrpPath << "; Run number: " << mRunNumber;
grpo = ccdb->template getForRun<o2::parameters::GRPMagField>(mGrpPath, mRunNumber);
o2::parameters::GRPMagField* grpo = ccdb->template getForRun<o2::parameters::GRPMagField>(mGrpPath, mRunNumber);
LOG(info) << "Get magnetic field with Path: " << mGrpPath << "; Run number: " << mRunNumber;
if (grpo == nullptr) {
LOG(fatal) << "GRP object not found for Run " << mRunNumber;
}
mMagField = static_cast<int>(grpo->getNominalL3Field()); // get magnetic field in kG
} else {
LOG(info) << "Force magentic field to " << mMagFieldForced << "kG";
LOG(info) << "Force magnetic field to " << mMagFieldForced << "kG";
mMagField = mMagFieldForced;
}

Expand Down Expand Up @@ -752,6 +759,12 @@ class CollisionBuilder
mCollisionSelection.getEventPlane());
}

if (mProducedLiteShapes) {
collisionProducts.producedLiteShapes(
o2::aod::femtocollisions::lite::binQvec(mCollisionSelection.getQvector()),
o2::aod::femtocollisions::lite::binEventPlaneAngle(mCollisionSelection.getEventPlane()));
}

mCollisionAlreadyFilled = true;
}

Expand Down Expand Up @@ -795,6 +808,7 @@ class CollisionBuilder
bool mProducedMultiplicities = false;
bool mProducedCentralities = false;
bool mProducedShapes = false;
bool mProducedLiteShapes = false;
};

struct CollisionBuilderDerivedToDerivedProducts : o2::framework::ProducesGroup {
Expand Down
68 changes: 41 additions & 27 deletions PWGCF/Femto/Core/femtoUtils.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,20 @@ inline std::optional<T2> getIndex(const T1& index, const std::unordered_map<T1,
return std::nullopt;
}

/// Translate a global daughter index into a row of the current track table frame.
/// Aborts if the index does not fall inside the frame, which would otherwise
/// silently produce an out-of-range iterator.
template <typename T>
int64_t daughterRow(int64_t daughterId, T const& trackTable)
{
const int64_t row = daughterId - trackTable.offset();
if (daughterId < 0 || row < 0 || row >= static_cast<int64_t>(trackTable.size())) {
LOG(fatal) << "Daughter index " << daughterId << " out of range for track table (offset "
<< trackTable.offset() << ", size " << trackTable.size() << "). Breaking...";
}
return row;
}

template <typename T>
float itsSignal(T const& track)
{
Expand Down Expand Up @@ -243,7 +257,7 @@ inline int signum(T x)
}

template <typename T>
inline T binLinear(float value, float lo, float hi, float step)
T binLinear(float value, float lo, float hi, float step)
{
float v = std::clamp(value, lo, hi);
auto idx = static_cast<int64_t>(std::round((v - lo) / step));
Expand All @@ -253,46 +267,46 @@ inline T binLinear(float value, float lo, float hi, float step)
}

template <typename T>
inline float unBinLinear(T binned, float lo, float step)
float unBinLinear(T binned, float lo, float step)
{
auto idx = static_cast<int64_t>(binned) - static_cast<int64_t>(std::numeric_limits<T>::min());
return lo + static_cast<float>(idx) * step;
}

template <typename T>
inline T binLogSigned(float signedValue, float magMin, float magMax)
T binLogSigned(float signedValue, float magMin, float magMax)
{
static_assert(std::is_unsigned_v<T>, "binLogSigned requires an unsigned storage type");
constexpr uint32_t TotalBits = sizeof(T) * 8;
constexpr uint32_t HalfLevels = 1u << (TotalBits - 1);
uint32_t sign = (signedValue < 0.f) ? 1u : 0u;
float mag = std::clamp(std::fabs(signedValue), magMin, magMax);
float logLo = std::log(magMin);
float logHi = std::log(magMax);
float step = (logHi - logLo) / static_cast<float>(HalfLevels - 1);
auto idx = static_cast<uint32_t>(std::round((std::log(mag) - logLo) / step));
idx = std::clamp(idx, 0u, HalfLevels - 1);
return static_cast<T>((sign << (TotalBits - 1)) | idx);
constexpr uint64_t TotalBits = sizeof(T) * 8;
constexpr uint64_t HalfLevels = uint64_t{1} << (TotalBits - 1);
const uint64_t sign = (signedValue < 0.f) ? uint64_t{1} : uint64_t{0};
const float mag = std::clamp(std::fabs(signedValue), magMin, magMax);
const float logLo = std::log(magMin);
const float logHi = std::log(magMax);
const float step = (logHi - logLo) / static_cast<float>(HalfLevels - 1);
auto idx = static_cast<int64_t>(std::round((std::log(mag) - logLo) / step));
idx = std::clamp(idx, int64_t{0}, static_cast<int64_t>(HalfLevels - 1));
return static_cast<T>((sign << (TotalBits - 1)) | static_cast<uint64_t>(idx));
}

template <typename T>
inline float unBinLogSigned(T binned, float magMin, float magMax)
float unBinLogSigned(T binned, float magMin, float magMax)
{
constexpr uint32_t TotalBits = sizeof(T) * 8;
constexpr uint32_t HalfLevels = 1u << (TotalBits - 1);
constexpr T SignMask = static_cast<T>(1u << (TotalBits - 1));
static_assert(std::is_unsigned_v<T>, "unBinLogSigned requires an unsigned storage type");
constexpr uint64_t TotalBits = sizeof(T) * 8;
constexpr uint64_t HalfLevels = uint64_t{1} << (TotalBits - 1);
constexpr T SignMask = static_cast<T>(uint64_t{1} << (TotalBits - 1));
constexpr T MagMask = static_cast<T>(SignMask - 1);
float sign = (binned & SignMask) ? -1.f : 1.f;
uint32_t idx = binned & MagMask;
float logLo = std::log(magMin);
float logHi = std::log(magMax);
float step = (logHi - logLo) / static_cast<float>(HalfLevels - 1);
float mag = std::exp(logLo + static_cast<float>(idx) * step);
return sign * mag;
const float sign = (binned & SignMask) ? -1.f : 1.f;
const auto idx = static_cast<uint64_t>(binned & MagMask);
const float logLo = std::log(magMin);
const float logHi = std::log(magMax);
const float step = (logHi - logLo) / static_cast<float>(HalfLevels - 1);
return sign * std::exp(logLo + static_cast<float>(idx) * step);
}

template <typename T>
inline int unBinSign(T binned)
int unBinSign(T binned)
{
static_assert(std::is_unsigned_v<T>, "unBinSign requires an unsigned storage type");
constexpr uint64_t TotalBits = sizeof(T) * 8;
Expand All @@ -301,7 +315,7 @@ inline int unBinSign(T binned)
}

template <typename T>
inline T binLogUnsigned(float value, float magMin, float magMax)
T binLogUnsigned(float value, float magMin, float magMax)
{
static_assert(std::is_unsigned_v<T>, "binLogUnsigned requires an unsigned storage type");
constexpr uint64_t TotalBits = sizeof(T) * 8;
Expand All @@ -316,7 +330,7 @@ inline T binLogUnsigned(float value, float magMin, float magMax)
}

template <typename T>
inline float unBinLogUnsigned(T binned, float magMin, float magMax)
float unBinLogUnsigned(T binned, float magMin, float magMax)
{
constexpr uint64_t TotalBits = sizeof(T) * 8;
constexpr uint64_t Levels = uint64_t{1} << TotalBits;
Expand Down
Loading
Loading