Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
10 changes: 5 additions & 5 deletions PWGCF/Femto/Core/baseSelection.h
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ class BaseSelection
{
// check index
if (static_cast<std::size_t>(observableIndex) >= NumObservables) {
LOG(fatal) << "Observable is not valid. Observable (index) has to be smaller than " << NumObservables;
LOG(fatal) << "Selection '" << selectionName << "': observable is not valid. Observable (index) has to be smaller than " << NumObservables;
}
// init selection container for selection at given index
mSelectionContainers.at(observableIndex) = selectioncontainer::SelectionContainer<T, BitmaskType>(selectionName,
Expand Down Expand Up @@ -112,7 +112,7 @@ class BaseSelection
bool isOptionalCut)
{
if (static_cast<std::size_t>(observableIndex) >= NumObservables) {
LOG(fatal) << "Observable is not valid. Observable (index) has to be smaller than " << NumObservables;
LOG(fatal) << "Selection '" << selectionName << "': observable is not valid. Observable (index) has to be smaller than " << NumObservables;
}
mSelectionContainers.at(observableIndex) = selectioncontainer::SelectionContainer<T, BitmaskType>(selectionName,
lowerLimit,
Expand Down Expand Up @@ -146,7 +146,7 @@ class BaseSelection
{
// check index
if (static_cast<std::size_t>(observableIndex) >= NumObservables) {
LOG(fatal) << "Observable is not valid. Observable (index) has to be smaller than " << NumObservables;
LOG(fatal) << "Selection '" << selectionName << "': observable is not valid. Observable (index) has to be smaller than " << NumObservables;
}
// init selection container for selection at given index
mSelectionContainers.at(observableIndex) = selectioncontainer::SelectionContainer<T, BitmaskType>(selectionName,
Expand Down Expand Up @@ -191,7 +191,7 @@ class BaseSelection
mSelectionContainers.at(observableIndex) = selectioncontainer::SelectionContainer<T, BitmaskType>(selectionName, std::vector<T>{1}, limits::LimitType::kEqual, false, false, false);
break;
default:
LOG(fatal) << "Invalid switch for boolean selection";
LOG(fatal) << "Selection '" << selectionName << "': invalid switch for boolean selection";
}
addSelection(observableIndex);
}
Expand Down Expand Up @@ -577,7 +577,7 @@ class BaseSelection
}

if (mNSelectionBits > sizeof(BitmaskType) * CHAR_BIT) {
LOG(fatal) << "Too many selections. At most " << sizeof(BitmaskType) * CHAR_BIT << " number of bits are supported";
LOG(fatal) << "Selection '" << mSelectionContainers.at(observableIndex).getSelectionName() << "': too many selections. At most " << sizeof(BitmaskType) * CHAR_BIT << " number of bits are supported";
}
}

Expand Down
58 changes: 42 additions & 16 deletions PWGCF/Femto/Core/selectionContainer.h
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@
#include <climits>
#include <cmath>
#include <cstddef>
#include <exception>
#include <limits>
#include <sstream>
#include <string>
Expand Down Expand Up @@ -122,10 +123,10 @@ class SelectionContainer
mIsOptionalCut(isOptionalCut)
{
if (mSelectionValues.size() > sizeof(BitmaskType) * CHAR_BIT) {
LOG(fatal) << "Too many selections for single a observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
LOG(fatal) << "Selection '" << mSelectionName << "': too many selections for a single observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
}
if (isMinimalCut && isOptionalCut) {
LOG(fatal) << "A selection cannot be both minimal and optional";
LOG(fatal) << "Selection '" << mSelectionName << "': a selection cannot be both minimal and optional";
}
// values for selection are not necessarily ordered correctly
sortSelections();
Expand Down Expand Up @@ -157,13 +158,13 @@ class SelectionContainer
mIsOptionalCut(isOptionalCut)
{
if (mLimitType != limits::kInsideRange && mLimitType != limits::kOutsideRange) {
LOG(fatal) << "Range-string selections only support kInsideRange or kOutsideRange";
LOG(fatal) << "Selection '" << mSelectionName << "': range-string selections only support kInsideRange or kOutsideRange";
}
if (rangeStrings.size() > sizeof(BitmaskType) * CHAR_BIT) {
LOG(fatal) << "Too many selections for a single observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
LOG(fatal) << "Selection '" << mSelectionName << "': too many selections for a single observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
}
if (isMinimalCut && isOptionalCut) {
LOG(fatal) << "A selection cannot be both minimal and optional";
LOG(fatal) << "Selection '" << mSelectionName << "': a selection cannot be both minimal and optional";
}
for (auto const& rangeStr : rangeStrings) {
T lower, upper;
Expand Down Expand Up @@ -196,7 +197,7 @@ class SelectionContainer
: (mSelectionRanges[i].first <= mSelectionRanges[i + 1].first &&
mSelectionRanges[i].second >= mSelectionRanges[i + 1].second);
if (!nested) {
LOG(fatal) << "Ranges for minimal cut " << selectionName
LOG(fatal) << "Selection '" << selectionName << "': ranges for minimal cut"
<< " are not nested. Range [" << mSelectionRanges[i].first << ";" << mSelectionRanges[i].second
<< "] is not " << (mLimitType == limits::kOutsideRange ? "contained in" : "containing")
<< " [" << mSelectionRanges[i + 1].first << ";" << mSelectionRanges[i + 1].second << "]";
Expand Down Expand Up @@ -229,10 +230,10 @@ class SelectionContainer
mIsOptionalCut(isOptionalCut)
{
if (functions.size() > sizeof(BitmaskType) * CHAR_BIT) {
LOG(fatal) << "Too many selections for single a observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
LOG(fatal) << "Selection '" << mSelectionName << "': too many selections for a single observable. Limit is " << sizeof(BitmaskType) * CHAR_BIT;
}
if (isMinimalCut && isOptionalCut) {
LOG(fatal) << "A selection cannot be both minimal and optional";
LOG(fatal) << "Selection '" << mSelectionName << "': a selection cannot be both minimal and optional";
}
for (std::size_t i = 0; i < functions.size(); i++) {
mSelectionFunctions.emplace_back((mSelectionName + std::to_string(i)).c_str(), functions.at(i).c_str(), lowerLimit, upperLimit);
Expand All @@ -255,7 +256,7 @@ class SelectionContainer
// note: threshold values may be reordered by sortSelections() or sortFunctions(),
// so comments must be provided in the already-sorted order
if (comments.size() != getNSelections()) {
LOG(fatal) << "Number of comments and number of selections are inconsistent";
LOG(fatal) << "Selection '" << mSelectionName << "': number of comments and number of selections are inconsistent";
}
mComments = comments;
}
Expand Down Expand Up @@ -382,7 +383,7 @@ class SelectionContainer
void evaluate(std::vector<T> const& values)
{
if (values.size() != mSelectionValues.size()) {
LOG(fatal) << "Wrong number of values have been passed";
LOG(fatal) << "Selection '" << mSelectionName << "': wrong number of values have been passed";
}
// reset the bitmask before evaluating
mBitmask.reset();
Expand Down Expand Up @@ -462,7 +463,7 @@ class SelectionContainer
T getLoosestSelection() const
{
if (mSelectionValues.empty()) {
LOG(fatal) << "No selections configured";
LOG(fatal) << "Selection '" << mSelectionName << "': no selections configured";
}
return mSelectionValues.at(0);
}
Expand Down Expand Up @@ -605,7 +606,7 @@ class SelectionContainer
[[nodiscard]] int getBitPosition(int selectionIndex) const
{
if (selectionIndex == 0 && mSkipMostPermissiveBit) {
LOG(fatal) << "Trying to accessed the bit position of a skipped selection. Breaking...";
LOG(fatal) << "Selection '" << mSelectionName << "': trying to access the bit position of a skipped selection. Breaking...";
return -1;
}
if (mSkipMostPermissiveBit) {
Expand Down Expand Up @@ -706,7 +707,10 @@ class SelectionContainer
{
auto pos = rangeStr.find(';');
if (pos == std::string::npos) {
LOG(fatal) << "Range string '" << rangeStr << "' is missing ';' separator. Expected format: 'lower;upper'";
LOG(fatal) << "Selection '" << mSelectionName << "': range string '" << rangeStr << "' is missing ';' separator. Expected format: 'lower;upper'";
}
if (rangeStr.find(';', pos + 1) != std::string::npos) {
LOG(fatal) << "Selection '" << mSelectionName << "': range string '" << rangeStr << "' has more than one ';' separator. Expected format: 'lower;upper'";
}
std::string lowerStr = rangeStr.substr(0, pos);
std::string upperStr = rangeStr.substr(pos + 1);
Expand All @@ -723,14 +727,36 @@ class SelectionContainer
return;
}

lower = lowerStr.empty() ? -std::numeric_limits<T>::max() : static_cast<T>(std::stod(lowerStr));
upper = upperStr.empty() ? std::numeric_limits<T>::max() : static_cast<T>(std::stod(upperStr));
lower = lowerStr.empty() ? -std::numeric_limits<T>::max() : parseBound(rangeStr, lowerStr);
upper = upperStr.empty() ? std::numeric_limits<T>::max() : parseBound(rangeStr, upperStr);

if (lower >= upper) {
LOG(fatal) << "Range string '" << rangeStr << "' has lower bound >= upper bound";
LOG(fatal) << "Selection '" << mSelectionName << "': range string '" << rangeStr << "' has lower bound >= upper bound";
}
}

/// \brief Parse a single range bound and make sure the whole string was consumed, so that
/// malformed input such as "3;5" (extra ';') or "-3,5" (comma used as decimal
/// separator) is rejected instead of silently truncated by std::stod (which happily
/// parses just the numeric prefix "3" or "-3" and ignores the rest).
/// \param rangeStr Full "lower;upper" string, only used for the error message.
/// \param boundStr The individual bound substring (lower or upper) to parse.
T parseBound(std::string const& rangeStr, std::string const& boundStr) const
{
std::size_t idx = 0;
double value = 0.;
try {
value = std::stod(boundStr, &idx);
} catch (std::exception const& e) {
LOG(fatal) << "Selection '" << mSelectionName << "': range string '" << rangeStr << "' contains bound '" << boundStr << "' that is not a valid number (" << e.what() << ")";
}
if (idx != boundStr.size()) {
LOG(fatal) << "Selection '" << mSelectionName << "': range string '" << rangeStr << "' contains bound '" << boundStr
<< "' with trailing characters ('" << boundStr.substr(idx) << "') that could not be parsed";
}
return static_cast<T>(value);
}

std::string mSelectionName;
std::vector<T> mSelectionValues = {}; ///< Threshold values, sorted from most permissive to most restrictive
std::vector<TF1> mSelectionFunctions; ///< TF1 threshold functions (empty for static selections)
Expand Down
116 changes: 116 additions & 0 deletions PWGCF/Femto/Core/twoTrackResonanceBuilder.h
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,7 @@
#include <Math/Vector4D.h> // IWYU pragma: keep (do not replace with Math/Vector4Dfwd.h)
#include <Math/Vector4Dfwd.h>

#include <cstdint>
#include <string>

namespace o2::analysis::femto::twotrackresonancebuilder
Expand Down Expand Up @@ -341,5 +342,120 @@ class TwoTrackResonanceBuilder

}; // namespace twotrackresonancebuilder

struct ConfTwoTrackResonanceTablesDerivedToDerived : o2::framework::ConfigurableGroup {
std::string prefix = std::string("TwoTrackResonanceTables");
o2::framework::Configurable<int> limitRho0{"limitRho0", 1, "Require at least this many Rho0s in the collision. Set to 0 to skip Rho0s entirely (not written to output)."};
o2::framework::Configurable<int> limitPhi{"limitPhi", 0, "Require at least this many Phis in the collision. Set to 0 to skip Phis entirely (not written to output)."};
};

struct TwoTrackResonanceBuilderDerivedToDerivedProducts : o2::framework::ProducesGroup {
o2::framework::Produces<o2::aod::StoredFRho0s> producedRho0s;
o2::framework::Produces<o2::aod::StoredFRho0Masks> producedRho0Masks;
o2::framework::Produces<o2::aod::StoredFPhis> producedPhis;
o2::framework::Produces<o2::aod::StoredFPhiMasks> producedPhiMasks;
};

class TwoTrackResonanceBuilderDerivedToDerived
{
public:
template <typename T>
void init(T& config)
{
mLimitRho0 = config.limitRho0.value;
mLimitPhi = config.limitPhi.value;

if (mLimitRho0 < 0 || mLimitPhi < 0) {
LOG(fatal) << "Two-track resonance limits must be non-negative (got " << mLimitRho0 << " and " << mLimitPhi << "). Breaking...";
}
if (mLimitRho0 == 0 && mLimitPhi == 0) {
LOG(fatal) << "Both Rho0 limit and Phi limit are 0. Breaking...";
}
}

template <typename T1, typename T2, typename T3, typename T4>
bool collisionHasTooFewRho0s(T1 const& col, T2 const& /*rho0Table*/, T3& partitionRho0, T4& cache)
{
if (mLimitRho0 == 0) { // Rho0s disabled, cannot reject on them
return false;
}
auto rho0Slice = partitionRho0->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);
return rho0Slice.size() < static_cast<int64_t>(mLimitRho0);
}

template <typename T1, typename T2, typename T3, typename T4>
bool collisionHasTooFewPhis(T1 const& col, T2 const& /*phiTable*/, T3& partitionPhi, T4& cache)
{
if (mLimitPhi == 0) { // Phis disabled, cannot reject on them
return false;
}
auto phiSlice = partitionPhi->sliceByCached(o2::aod::femtobase::stored::fColId, col.globalIndex(), cache);
return phiSlice.size() < static_cast<int64_t>(mLimitPhi);
}

template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
void processRho0s(T1 const& col, T2 const& /*rho0Table*/, T3 const& oldTrackTable, T4& partitionRho0, T5& trackBuilder, T6& cache, T7& newRho0Table, T8& newTrackTable, T9& newCollisionTable)
{
if (mLimitRho0 == 0) { // Rho0s disabled
return;
}

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

for (auto const& rho0 : rho0Slice) {
auto posDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(rho0.posDauId(), oldTrackTable));
auto negDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(rho0.negDauId(), oldTrackTable));

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

newRho0Table.producedRho0s(newCollisionTable.producedCollision.lastIndex(),
rho0.pt(),
rho0.eta(),
rho0.phi(),
rho0.mass(),
posDaughterIndex,
negDaughterIndex);
newRho0Table.producedRho0Masks(rho0.maskPosDau(),
rho0.posDauHasHighMomentum(),
rho0.maskNegDau(),
rho0.negDauHasHighMomentum());
}
}

template <typename T1, typename T2, typename T3, typename T4, typename T5, typename T6, typename T7, typename T8, typename T9>
void processPhis(T1 const& col, T2 const& /*phiTable*/, T3 const& oldTrackTable, T4& partitionPhi, T5& trackBuilder, T6& cache, T7& newPhiTable, T8& newTrackTable, T9& newCollisionTable)
{
if (mLimitPhi == 0) { // Phis disabled
return;
}

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

for (auto const& phi : phiSlice) {
auto posDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(phi.posDauId(), oldTrackTable));
auto negDaughter = oldTrackTable.rawIteratorAt(utils::daughterRow(phi.negDauId(), oldTrackTable));

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

newPhiTable.producedPhis(newCollisionTable.producedCollision.lastIndex(),
phi.pt(),
phi.eta(),
phi.phi(),
phi.mass(),
posDaughterIndex,
negDaughterIndex);
newPhiTable.producedPhiMasks(phi.maskPosDau(),
phi.posDauHasHighMomentum(),
phi.maskNegDau(),
phi.negDauHasHighMomentum());
}
}

private:
int mLimitRho0 = 0;
int mLimitPhi = 0;
};

} // namespace o2::analysis::femto::twotrackresonancebuilder
#endif // PWGCF_FEMTO_CORE_TWOTRACKRESONANCEBUILDER_H_
4 changes: 4 additions & 0 deletions PWGCF/Femto/DataModel/FemtoTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -552,12 +552,14 @@ DECLARE_SOA_TABLE_STAGED_VERSIONED(FPhis_001, "FPHI", 1, //! femto phis
femtobase::dynamic::Pz<femtobase::stored::Pt, femtobase::stored::Eta>,
femtobase::dynamic::Theta<femtobase::stored::Eta>);
using FPhis = FPhis_001;
using StoredFPhis = StoredFPhis_001;
DECLARE_SOA_TABLE_STAGED_VERSIONED(FPhiMasks_001, "FPHIMASK", 1, //! mask for phis
femtotwotrackresonances::MaskPosDau,
femtotwotrackresonances::PosDauHasHighMomentum,
femtotwotrackresonances::MaskNegDau,
femtotwotrackresonances::NegDauHasHighMomentum);
using FPhiMasks = FPhiMasks_001;
using StoredFPhiMasks = StoredFPhiMasks_001;

// table for kstars
DECLARE_SOA_TABLE_STAGED_VERSIONED(FKstar0s_001, "FKSTAR0", 1, //! femto k0star
Expand Down Expand Up @@ -599,12 +601,14 @@ DECLARE_SOA_TABLE_STAGED_VERSIONED(FRho0s_001, "FRHO0", 1, //! femto rho0s
femtobase::dynamic::Pz<femtobase::stored::Pt, femtobase::stored::Eta>,
femtobase::dynamic::Theta<femtobase::stored::Eta>);
using FRho0s = FRho0s_001;
using StoredFRho0s = StoredFRho0s_001;
DECLARE_SOA_TABLE_STAGED_VERSIONED(FRho0Masks_001, "FRHO0MASK", 1, //! rho0s masks
femtotwotrackresonances::MaskPosDau,
femtotwotrackresonances::PosDauHasHighMomentum,
femtotwotrackresonances::MaskNegDau,
femtotwotrackresonances::NegDauHasHighMomentum);
using FRho0Masks = FRho0Masks_001;
using StoredFRho0Masks = StoredFRho0Masks_001;

namespace femtov0s
{
Expand Down
Loading
Loading