diff --git a/PWGCF/EbyEFluctuations/Tasks/partNumFluc.cxx b/PWGCF/EbyEFluctuations/Tasks/partNumFluc.cxx index 6aff03adde6..2bcecf4fb12 100644 --- a/PWGCF/EbyEFluctuations/Tasks/partNumFluc.cxx +++ b/PWGCF/EbyEFluctuations/Tasks/partNumFluc.cxx @@ -48,7 +48,6 @@ #include #include #include -#include #include #include @@ -57,6 +56,7 @@ #include #include #include +#include #include #include #include @@ -70,16 +70,16 @@ #include #include -#define C_CS(cs) /*NOLINT(cppcoreguidelines-macro-usage)*/ \ - [](std::index_sequence) { \ - static_assert(std::is_array_v> && \ - std::same_as>, char>); \ - return ConstStr<(cs)[indices]...>{}; \ +#define C_CS(cs) /* NOLINT(cppcoreguidelines-macro-usage) */ \ + [](std::index_sequence) constexpr -> ConstStr<(cs)[indices]...> { \ + static_assert(std::is_array_v> && \ + std::same_as>, char>); \ + return ConstStr<(cs)[indices]...>{}; \ }(std::make_index_sequence{}) -#define C_SV(sv) /*NOLINT(cppcoreguidelines-macro-usage)*/ \ - [](std::index_sequence) { \ - static_assert(std::same_as, std::string_view>); \ - return ConstStr<(sv)[indices]...>{}; \ +#define C_SV(sv) /* NOLINT(cppcoreguidelines-macro-usage) */ \ + [](std::index_sequence) constexpr -> ConstStr<(sv)[indices]...> { \ + static_assert(std::same_as, std::string_view>); \ + return ConstStr<(sv)[indices]...>{}; \ }(std::make_index_sequence<(sv).size()>{}) using namespace o2; @@ -88,7 +88,7 @@ using namespace o2::framework::expressions; namespace o2::aod { -using JoinedCollisions = soa::Join; +using JoinedCollisions = soa::Join; using JoinedTracks = soa::Join; using JoinedCollisionsWithMc = soa::Join; using JoinedTracksWithMc = soa::Join; @@ -139,8 +139,8 @@ namespace fluctuation_calculator_base { inline constexpr std::int8_t MaxOrder{8}; inline constexpr std::int32_t NExponentKeys{MaxOrder * (MaxOrder + 1) / 2}; -inline constexpr std::array, NExponentKeys> ExponentKeys{[] { - std::array, NExponentKeys> result{}; +inline constexpr std::array, NExponentKeys> ExponentKeys{[]() constexpr -> std::array, NExponentKeys> { + std::array, NExponentKeys> result{}; std::int32_t index{}; for (std::int32_t const& iExponent : std::views::iota(1, MaxOrder + 1)) { for (std::int32_t const& jExponent : std::views::iota(1, iExponent + 1)) { @@ -149,25 +149,21 @@ inline constexpr std::array, NExponentKeys> ExponentK } return result; }()}; -inline constexpr std::int32_t NOrderKeys{[] { +inline constexpr std::int32_t NOrderKeys{[]() constexpr -> std::int32_t { std::array counts{1}; - for (std::array const& exponentKey /*o2-linter: disable=const-ref-in-for-loop*/ : ExponentKeys) { - const std::int32_t weight{exponentKey[0]}; + for (std::pair const& exponentKey /* o2-linter: disable=const-ref-in-for-loop */ : ExponentKeys) { + const std::int32_t weight{exponentKey.first}; for (std::int32_t const& sum : std::views::iota(weight, MaxOrder + 1)) { counts[sum] += counts[sum - weight]; } } - std::int32_t result{}; - for (std::int32_t const& count : counts) { - result += count; - } - return result; + return std::accumulate(counts.begin(), counts.end(), 0); }()}; -inline constexpr std::array, NOrderKeys> OrderKeys{[] { +inline constexpr std::array, NOrderKeys> OrderKeys{[]() constexpr -> std::array, NOrderKeys> { std::array, NOrderKeys> result{}; std::array current{}; std::int32_t index{}; - const auto fillOrderKeys{[&](const auto& self, const std::int32_t position, const std::int32_t sum, const std::int32_t target) constexpr -> void { + const auto fillOrderKeys{[¤t, &index, &result](const auto& self, const std::int32_t position, const std::int32_t sum, const std::int32_t target) constexpr -> void { if (position == NExponentKeys) { if (sum == target) { result[index++] = current; @@ -175,9 +171,9 @@ inline constexpr std::array, NOrderKeys> return; } - const std::int32_t weight{ExponentKeys[position][0]}; + const std::int32_t weight{ExponentKeys[position].first}; for (std::int32_t power{}; sum + power * weight <= target; ++power) { - current[position] = power; + current[position] = static_cast(power); self(self, position + 1, sum + power * weight, target); } }}; @@ -192,23 +188,46 @@ class FluctuationCalculatorTrack { public: FluctuationCalculatorTrack() = default; + FluctuationCalculatorTrack(const FluctuationCalculatorTrack&) = default; + FluctuationCalculatorTrack(FluctuationCalculatorTrack&&) noexcept = default; + FluctuationCalculatorTrack& operator=(const FluctuationCalculatorTrack&) = default; + FluctuationCalculatorTrack& operator=(FluctuationCalculatorTrack&&) noexcept = default; virtual ~FluctuationCalculatorTrack() = default; - [[nodiscard]] double getProduct(const std::int32_t orderKeyIndex, const double weight = 1.) const + [[nodiscard]] std::array getProducts(const double weight = 1.) const { - double product{1.}; - for (std::int32_t const& iExponentKey : std::views::iota(0, static_cast(fluctuation_calculator_base::NExponentKeys))) { - if (fluctuation_calculator_base::OrderKeys[orderKeyIndex][iExponentKey] != 0) { - product *= std::pow(mQs[iExponentKey], fluctuation_calculator_base::OrderKeys[orderKeyIndex][iExponentKey]); + std::array, fluctuation_calculator_base::NExponentKeys> powersQ{}; + for (std::int32_t const& iExponentKey : std::views::iota(0, fluctuation_calculator_base::NExponentKeys)) { + powersQ[iExponentKey][1] = mQs[iExponentKey]; + for (std::int32_t const& power : std::views::iota(2, fluctuation_calculator_base::MaxOrder / fluctuation_calculator_base::ExponentKeys[iExponentKey].first + 1)) { + powersQ[iExponentKey][power] = powersQ[iExponentKey][power - 1] * mQs[iExponentKey]; } } - return product * weight; + + std::array products{}; + products.fill(weight); + for (std::int32_t const& iOrderKey : std::views::iota(0, fluctuation_calculator_base::NOrderKeys)) { + for (std::int32_t const& iExponentKey : std::views::iota(0, fluctuation_calculator_base::NExponentKeys)) { + if (const std::int32_t power{fluctuation_calculator_base::OrderKeys[iOrderKey][iExponentKey]}; power > 0) { + products[iOrderKey] *= powersQ[iExponentKey][power]; + } + } + } + return products; } - void init() { mQs.fill(0.); } + void clear() { mQs.fill({}); } void fill(const double charge, const double efficiency, const double weight = 1.) { - for (std::int32_t const& iExponentKey : std::views::iota(0, static_cast(fluctuation_calculator_base::NExponentKeys))) { - mQs[iExponentKey] += std::pow(charge, fluctuation_calculator_base::ExponentKeys[iExponentKey][0]) / std::pow(efficiency, fluctuation_calculator_base::ExponentKeys[iExponentKey][1]) * weight; + const double inverseEfficiency{1. / efficiency}; + std::array powersCharge{1.}; + std::array powersInverseEfficiency{1.}; + for (std::int32_t const& exponent : std::views::iota(1, fluctuation_calculator_base::MaxOrder + 1)) { + powersCharge[exponent] = powersCharge[exponent - 1] * charge; + powersInverseEfficiency[exponent] = powersInverseEfficiency[exponent - 1] * inverseEfficiency; + } + for (std::int32_t const& iExponentKey : std::views::iota(0, fluctuation_calculator_base::NExponentKeys)) { + const auto& [exponentCharge, exponentEfficiency]{fluctuation_calculator_base::ExponentKeys[iExponentKey]}; + mQs[iExponentKey] += weight * powersCharge[exponentCharge] * powersInverseEfficiency[exponentEfficiency]; } } @@ -223,14 +242,14 @@ concept IsValid = IsValidEnum...> && ((EV template requires IsValidEnum -inline constexpr std::int32_t toI(const E e) +constexpr std::int32_t toI(const E e) { return static_cast(e); } template requires IsValidEnum -inline constexpr std::int32_t NEs{toI(E::N)}; +constexpr std::int32_t NEs{toI(E::N)}; template requires IsValidEnum @@ -245,7 +264,7 @@ enum class NameKind { }; template requires IsValidEnum -inline constexpr std::string_view getName(const std::int32_t index) +constexpr std::string_view getName(const std::int32_t index) { if constexpr (NameKindValue == NameKind::Lower) { return EnumInfo>::NamesLower.at(index); @@ -259,15 +278,19 @@ inline constexpr std::string_view getName(const std::int32_t index) } template requires IsValidEnum -inline constexpr std::string_view getName(const E e) +constexpr std::string_view getName(const E e) { return getName(toI(e)); } template requires IsValidEnum -inline std::vector getDisplayNames() +std::vector getDisplayNames() { - return std::vector(EnumInfo>::DisplayNames.begin(), EnumInfo>::DisplayNames.end()); + if constexpr (requires { EnumInfo>::DisplayNames; }) { + return {EnumInfo>::DisplayNames.begin(), EnumInfo>::DisplayNames.end()}; + } else { + return {EnumInfo>::Names.begin(), EnumInfo>::Names.end()}; + } } enum class DataMode { @@ -282,7 +305,7 @@ enum class CentralityDefinition { Ft0m, N }; -enum class DcaKind { +enum class DcaMeasure { Mean = 0, Sigma, N @@ -342,119 +365,164 @@ enum class ChargeNumber { }; template <> -struct EnumInfo { - inline static constexpr std::array> Names{"Mean", "Sigma"}; +struct EnumInfo { + static constexpr std::array> Names{"Mean", "Sigma"}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"Xy", "Z"}; - inline static constexpr std::array> DisplayNames{"Xy", "Z"}; + static constexpr std::array> Names{"Xy", "Z"}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"Tpc", "Tof"}; - inline static constexpr std::array> NamesLower{"tpc", "tof"}; - inline static constexpr std::array> DisplayNames{"TPC", "TOF"}; - inline static constexpr std::array> PidStrategyAllValues{PidStrategyAll::Tpc, PidStrategyAll::Tof}; + static constexpr std::array> Names{"Tpc", "Tof"}; + static constexpr std::array> NamesLower{"tpc", "tof"}; + static constexpr std::array> DisplayNames{"TPC", "TOF"}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"Tpc", "TpcTof"}; - inline static constexpr std::array> NamesLower{"tpc", "tpcTof"}; - inline static constexpr std::array> DisplayNames{"TPC", "TPC+TOF"}; - inline static constexpr std::array> PidStrategyAllValues{PidStrategyAll::Tpc, PidStrategyAll::TpcTofCombined}; + static constexpr std::array> Names{"Tpc", "TpcTof"}; + static constexpr std::array> NamesLower{"tpc", "tpcTof"}; + static constexpr std::array> DisplayNames{"TPC", "TPC+TOF"}; + static constexpr std::array> PidStrategyAllValues{PidStrategyAll::Tpc, PidStrategyAll::TpcTofCombined}; +}; +template <> +struct EnumInfo { + static constexpr std::array> DetectorValues{Detector::Tpc, Detector::Tof, Detector::N, Detector::N}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"Pi", "Ka", "Pr"}; - inline static constexpr std::array> DisplayNames{"Pion", "Kaon", "Proton"}; - inline static constexpr std::array> DisplayNamesLower{"pion", "kaon", "proton"}; - inline static constexpr std::array> ParticleSpeciesAllValues{ParticleSpeciesAll::Pion, ParticleSpeciesAll::Kaon, ParticleSpeciesAll::Proton}; + static constexpr std::array> Names{"Pi", "Ka", "Pr"}; + static constexpr std::array> DisplayNames{"Pion", "Kaon", "Proton"}; + static constexpr std::array> DisplayNamesLower{"pion", "kaon", "proton"}; + static constexpr std::array> ParticleSpeciesAllValues{ParticleSpeciesAll::Pion, ParticleSpeciesAll::Kaon, ParticleSpeciesAll::Proton}; + static constexpr std::array>, NEs> PdgCodes{{{PDG_t::kPiPlus, PDG_t::kPiMinus}, {PDG_t::kKPlus, PDG_t::kKMinus}, {PDG_t::kProton, PDG_t::kProtonBar}}}; + static constexpr std::array> Masses{constants::physics::MassPiPlus, constants::physics::MassKPlus, constants::physics::MassProton}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"", "Pi", "Ka", "Pr"}; - inline static constexpr std::array> DisplayNames{"All", "Pion", "Kaon", "Proton"}; - inline static constexpr std::array> DisplayNamesLower{"all", "pion", "kaon", "proton"}; - inline static constexpr std::array> ParticleSpeciesValues{ParticleSpecies::N, ParticleSpecies::Pion, ParticleSpecies::Kaon, ParticleSpecies::Proton}; - inline static constexpr std::array>, NEs> PdgCodes{{{0, 0}, {PDG_t::kPiPlus, PDG_t::kPiMinus}, {PDG_t::kKPlus, PDG_t::kKMinus}, {PDG_t::kProton, PDG_t::kProtonBar}}}; - inline static constexpr std::array> Masses{0., constants::physics::MassPiPlus, constants::physics::MassKPlus, constants::physics::MassProton}; + static constexpr std::array> Names{"", "Pi", "Ka", "Pr"}; + static constexpr std::array> DisplayNames{"All", "Pion", "Kaon", "Proton"}; + static constexpr std::array> DisplayNamesLower{"all", "pion", "kaon", "proton"}; + static constexpr std::array> Titles{"", "#pi", "K", "p"}; + static constexpr std::array> ParticleSpeciesValues{ParticleSpecies::N, ParticleSpecies::Pion, ParticleSpecies::Kaon, ParticleSpecies::Proton}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"Ch", "Ka", "Pr"}; - inline static constexpr std::array> DisplayNames{"Charge", "Kaon", "Proton"}; - inline static constexpr std::array> DisplayNamesLower{"charge", "kaon", "proton"}; + static constexpr std::array> Names{"Ch", "Ka", "Pr"}; + static constexpr std::array> DisplayNames{"Charge", "Kaon", "Proton"}; + static constexpr std::array> DisplayNamesLower{"charge", "kaon", "proton"}; + static constexpr std::array> Titles{"h", "K", "p"}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"P", "M"}; - inline static constexpr std::array> NamesLower{"p", "m"}; + static constexpr std::array> Names{"P", "M"}; + static constexpr std::array> NamesLower{"p", "m"}; + static constexpr std::array> Titles{"+", "#minus"}; }; template <> struct EnumInfo { - inline static constexpr std::array> Names{"P", "M", "T", "N"}; + static constexpr std::array> Names{"P", "M", "T", "N"}; }; +template + requires IsValidEnum +constexpr std::string_view getTitle(const std::int32_t index) +{ + return EnumInfo>::Titles.at(index); +} template requires IsValidEnum -inline constexpr To getValue(const E e) +constexpr To getValue(const E e) { - if constexpr (std::is_same_v, PidStrategyAll>) { + if constexpr (std::is_same_v, Detector>) { + return EnumInfo::DetectorValues[toI(e)]; + } else if constexpr (std::is_same_v, PidStrategyAll>) { return EnumInfo::PidStrategyAllValues[toI(e)]; } else if constexpr (std::is_same_v, ParticleSpecies>) { return EnumInfo::ParticleSpeciesValues[toI(e)]; } else if constexpr (std::is_same_v, ParticleSpeciesAll>) { return EnumInfo::ParticleSpeciesAllValues[toI(e)]; + } else { + return To::N; } - return To::N; } -inline constexpr std::int32_t getPdgCode(const ParticleSpeciesAll e, const ChargeSpecies chargeSpecies) +constexpr std::int32_t getPdgCode(const ParticleSpecies particleSpecies, const ChargeSpecies chargeSpecies) { - return EnumInfo::PdgCodes[toI(e)][toI(chargeSpecies)]; + return EnumInfo::PdgCodes[toI(particleSpecies)][toI(chargeSpecies)]; } -inline constexpr double getMass(const ParticleSpeciesAll e) +constexpr double getMass(const ParticleSpecies particleSpecies) { - return EnumInfo::Masses[toI(e)]; + return EnumInfo::Masses[toI(particleSpecies)]; } template requires std::is_arithmetic_v -inline constexpr std::int32_t nEnabled(const Configurable>& cfg) +constexpr std::int32_t nEnabled(const Configurable>& cfg) { const LabeledArray& la{cfg.value}; - return std::count_if(la[0], la[0] + la.rows() * la.cols(), [](T x) { return x != T{}; }); + return std::count_if(la[0], la[0] + la.rows() * la.cols(), [](const T x) constexpr -> bool { return x != T{}; }); } template requires std::is_arithmetic_v -inline constexpr bool isEnabled(const Configurable>& cfg) +constexpr bool isEnabled(const Configurable>& cfg) { return nEnabled(cfg) > 0; } + +double interpolate(const TH3* const h, const double x, const double y, const double z) +{ + if (!h) { + return 0.; + } + + static constexpr auto GetBinIndicesWeights{[](const TAxis* const axis, const double position) -> std::array, 2> { + if (!axis) { + return {}; + } + + const std::int32_t n{axis->GetNbins()}; + if (n == 1 || position <= axis->GetBinCenter(1)) { + return {{{1, 1.}, {1, 0.}}}; + } + if (position >= axis->GetBinCenter(n)) { + return {{{n, 1.}, {n, 0.}}}; + } + + const std::int32_t bin{axis->FindFixBin(position)}; + const std::int32_t lower{position < axis->GetBinCenter(bin) ? bin - 1 : bin}; + const std::int32_t upper{lower + 1}; + const double fraction{(position - axis->GetBinCenter(lower)) / (axis->GetBinCenter(upper) - axis->GetBinCenter(lower))}; + + return {{{lower, 1. - fraction}, {upper, fraction}}}; + }}; + + const std::array, 2> xb{GetBinIndicesWeights(h->GetXaxis(), x)}; + const std::array, 2> yb{GetBinIndicesWeights(h->GetYaxis(), y)}; + const std::array, 2> zb{GetBinIndicesWeights(h->GetZaxis(), z)}; + + double result{}; + for (const auto& [ix, wx] : xb) { + for (const auto& [iy, wy] : yb) { + for (const auto& [iz, wz] : zb) { + result += wx * wy * wz * h->GetBinContent(ix, iy, iz); + } + } + } + return result; +} } // namespace struct PartNumFluc { struct HolderCcdb { - [[maybe_unused]] inline static constexpr std::int32_t NDimensionsEfficiency{4}; + [[maybe_unused]] static constexpr std::int32_t NDimensionsEfficiency{4}; std::map> runNumbersIndicesGroupIndices; - std::vector>, NEs>, NEs>> fPtDca; + std::vector, NEs>, NEs>, NEs>> fPtMeasureDca; std::vector>, NEs>, NEs>> hCentralityPtEtaShiftNSigmaPid; std::vector>, NEs>, NEs>> hVzCentralityPtEtaEfficiency; - - void clear() - { - runNumbersIndicesGroupIndices.clear(); - fPtDca.clear(); - hCentralityPtEtaShiftNSigmaPid.clear(); - hVzCentralityPtEtaEfficiency.clear(); - } } holderCcdb{}; struct HolderMcEvent { - std::int32_t runNumber{}; - std::int32_t runIndex{}; - std::int32_t runGroupIndex{}; double vz{}; std::array>, NEs> numbers{}; std::array>, NEs> numbersEff{}; @@ -463,7 +531,8 @@ struct PartNumFluc { } holderMcEvent{}; struct HolderEvent { - inline static constexpr std::pair RangeCentrality{0., 100.}; + static constexpr std::pair RangeCentrality{0., 100.}; + static constexpr bool isValidCentrality(const double value) { return RangeCentrality.first <= value && value <= RangeCentrality.second; } std::int32_t runNumber{}; std::int32_t runIndex{}; @@ -471,7 +540,7 @@ struct PartNumFluc { double vz{}; std::array> nGlobalTracks{}; std::array> nPvContributors{}; - std::array>, NEs>, NEs> dca{}; + std::array>, NEs>, NEs> measureDca{}; std::array> nTofBeta{}; double centralityCalibration{}; double centrality{}; @@ -481,9 +550,9 @@ struct PartNumFluc { void clear() { *this = {}; } [[nodiscard]] std::int32_t getNGlobalTracks() const { return std::accumulate(nGlobalTracks.begin(), nGlobalTracks.end(), 0); } [[nodiscard]] std::int32_t getNPvContributors() const { return std::accumulate(nPvContributors.begin(), nPvContributors.end(), 0); } - template - requires IsValid - [[nodiscard]] double getDca() const + template + requires IsValid + [[nodiscard]] double getMeasureDca() const { const std::int32_t sumNGlobalTracks{getNGlobalTracks()}; if (sumNGlobalTracks == 0) { @@ -491,18 +560,18 @@ struct PartNumFluc { } double sumDca{}; - if constexpr (DcaKindValue == DcaKind::Sigma) { - const double meanDca{getDca()}; + if constexpr (DcaMeasureValue == DcaMeasure::Sigma) { + const double meanDca{getMeasureDca()}; for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - sumDca += (std::pow(dca[toI(DcaKind::Sigma)][toI(DcaAxisValue)][iChargeSpecies], 2.) + std::pow(dca[toI(DcaKind::Mean)][toI(DcaAxisValue)][iChargeSpecies] - meanDca, 2.)) * nGlobalTracks[iChargeSpecies]; + sumDca += (std::pow(measureDca[toI(DcaMeasure::Sigma)][toI(DcaAxisValue)][iChargeSpecies], 2.) + std::pow(measureDca[toI(DcaMeasure::Mean)][toI(DcaAxisValue)][iChargeSpecies] - meanDca, 2.)) * nGlobalTracks[iChargeSpecies]; } return std::sqrt(sumDca / sumNGlobalTracks); + } else { + for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { + sumDca += measureDca[toI(DcaMeasure::Mean)][toI(DcaAxisValue)][iChargeSpecies] * nGlobalTracks[iChargeSpecies]; + } + return sumDca / sumNGlobalTracks; } - - for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - sumDca += dca[toI(DcaKind::Mean)][toI(DcaAxisValue)][iChargeSpecies] * nGlobalTracks[iChargeSpecies]; - } - return sumDca / sumNGlobalTracks; } [[nodiscard]] std::int32_t getNTofBeta() const { return std::accumulate(nTofBeta.begin(), nTofBeta.end(), 0); } } holderEvent{}; @@ -517,19 +586,26 @@ struct PartNumFluc { } holderMcParticle{}; struct HolderTrack { - inline static constexpr double TruncationAbsNSigmaPid{999.}; - static constexpr double truncateNSigmaPid(const double value) { return std::abs(value) < TruncationAbsNSigmaPid ? value : -TruncationAbsNSigmaPid; } + static constexpr double TruncationAbsNSigmaPid{999.}; + static constexpr double truncateNSigmaPid(const double value, const double shift = {}) + { + const double valueShifted{value - shift}; + return std::abs(value) < TruncationAbsNSigmaPid && std::abs(valueShifted) < TruncationAbsNSigmaPid ? valueShifted : -TruncationAbsNSigmaPid; + } + + [[nodiscard]] double getNSigmaPidCombined(const std::int32_t particleSpeciesIndex) const + { + return truncateNSigmaPid(std::copysign(std::hypot(nSigmaPid[toI(Detector::Tpc)][particleSpeciesIndex], nSigmaPid[toI(Detector::Tof)][particleSpeciesIndex]), nSigmaPid[toI(Detector::Tpc)][particleSpeciesIndex] + nSigmaPid[toI(Detector::Tof)][particleSpeciesIndex])); + } std::array> dca{}; std::int32_t sign{}; - double p{}; double pt{}; double eta{}; double phi{}; - double phiIu{}; std::array> hasPid{}; - std::array>, NEs> nSigmaPid{[] { - std::array>, NEs> a{}; + std::array>, NEs> nSigmaPid{[]() constexpr -> std::array>, NEs> { + std::array>, NEs> a{}; std::array> b{}; b.fill(-TruncationAbsNSigmaPid); a.fill(b); @@ -553,8 +629,8 @@ struct PartNumFluc { std::array> nMcParticles{}; std::array> nTracks{}; - std::vector signedEfficienciesMcParticle{[] {std::vector v{}; v.reserve(256); return v; }()}; - std::vector signedEfficienciesTrack{[] {std::vector v{}; v.reserve(256); return v; }()}; + std::vector signedEfficienciesMcParticle{[]() constexpr -> std::vector { std::vector v{}; v.reserve(256); return v; }()}; + std::vector signedEfficienciesTrack{[]() constexpr -> std::vector { std::vector v{}; v.reserve(256); return v; }()}; void clear() { @@ -579,14 +655,14 @@ struct PartNumFluc { Configurable cfgFlagQaCentrality{"cfgFlagQaCentrality", false, "Centrality QA flag"}; Configurable cfgFlagQaTrack{"cfgFlagQaTrack", false, "Track QA flag"}; Configurable cfgFlagQaDca{"cfgFlagQaDca", false, "DCA QA flag"}; - Configurable> cfgFlagsQaAcceptance{"cfgFlagsQaAcceptance", {std::array>{0, 0, 0, 0}.data(), NEs, getDisplayNames()}, "Acceptance QA flags"}; - Configurable> cfgFlagsQaPhi{"cfgFlagsQaPhi", {std::array>{0, 0, 0, 0}.data(), NEs, getDisplayNames()}, "Phi QA flags"}; - Configurable> cfgFlagsQaPid{"cfgFlagsQaPid", {std::array>{0, 0, 0, 0}.data(), NEs, getDisplayNames()}, "PID QA flags"}; + Configurable> cfgFlagsQaAcceptance{"cfgFlagsQaAcceptance", {std::array>{}.data(), NEs, getDisplayNames()}, "Acceptance QA flags"}; + Configurable> cfgFlagsQaPhi{"cfgFlagsQaPhi", {std::array>{}.data(), NEs, getDisplayNames()}, "Phi QA flags"}; + Configurable> cfgFlagsQaPid{"cfgFlagsQaPid", {std::array>{}.data(), NEs, getDisplayNames()}, "PID QA flags"}; Configurable cfgFlagQaMc{"cfgFlagQaMc", false, "MC QA flag"}; - Configurable> cfgFlagsCalculationYield{"cfgFlagsCalculationYield", {std::array>{0, 0, 0}.data(), NEs, getDisplayNames()}, "Yield calculation flags"}; - Configurable> cfgFlagsCalculationPurity{"cfgFlagsCalculationPurity", {std::array>{0, 0, 0}.data(), NEs, getDisplayNames()}, "Purity calculation flags"}; - Configurable> cfgFlagsCalculationFractionPrimary{"cfgFlagsCalculationFractionPrimary", {std::array>{0, 0, 0}.data(), NEs, getDisplayNames()}, "Primary fraction calculation flags"}; - Configurable> cfgFlagsCalculationFluctuation{"cfgFlagsCalculationFluctuation", {std::array>{0, 0, 0}.data(), NEs, getDisplayNames()}, "Fluctuation calculation flags"}; + Configurable> cfgFlagsCalculationYield{"cfgFlagsCalculationYield", {std::array>{}.data(), NEs, getDisplayNames()}, "Yield calculation flags"}; + Configurable> cfgFlagsCalculationPurity{"cfgFlagsCalculationPurity", {std::array>{}.data(), NEs, getDisplayNames()}, "Purity calculation flags"}; + Configurable> cfgFlagsCalculationFractionPrimary{"cfgFlagsCalculationFractionPrimary", {std::array>{}.data(), NEs, getDisplayNames()}, "Primary fraction calculation flags"}; + Configurable> cfgFlagsCalculationFluctuation{"cfgFlagsCalculationFluctuation", {std::array>{}.data(), NEs, getDisplayNames()}, "Fluctuation calculation flags"}; } groupAnalysis{}; struct : ConfigurableGroup { @@ -615,6 +691,7 @@ struct PartNumFluc { Configurable cfgCutMinItsNCls{"cfgCutMinItsNCls", 5, "Minimum number of clusters ITS"}; Configurable cfgCutMaxItsChi2NCls{"cfgCutMaxItsChi2NCls", 25., "Maximum chi2 per cluster ITS"}; Configurable cfgCutMinTpcNCls{"cfgCutMinTpcNCls", 55, "Minimum number of clusters TPC"}; + Configurable cfgCutMinTpcChi2NCls{"cfgCutMinTpcChi2NCls", 0., "Minimum chi2 per cluster TPC"}; Configurable cfgCutMaxTpcChi2NCls{"cfgCutMaxTpcChi2NCls", 3.5, "Maximum chi2 per cluster TPC"}; Configurable cfgCutMaxTpcNClsSharedRatio{"cfgCutMaxTpcNClsSharedRatio", 0.4, "Maximum ratio of shared clusters over clusters TPC"}; Configurable cfgCutMinTpcNCrossedRows{"cfgCutMinTpcNCrossedRows", 80, "Minimum number of crossed rows TPC"}; @@ -625,13 +702,21 @@ struct PartNumFluc { Configurable cfgCutMaxPt{"cfgCutMaxPt", 2., "Maximum pT (GeV/c)"}; Configurable cfgCutMaxAbsEta{"cfgCutMaxAbsEta", 0.8, "Maximum absolute eta"}; Configurable> cfgThresholdsPtTofPid{"cfgThresholdsPtTofPid", {std::array>{0.5, 0.5, 0.8}.data(), NEs, getDisplayNames()}, "pT (GeV/c) thresholds for TOF PID"}; - Configurable> cfgFlagsRecalibrationNSigmaPid{"cfgFlagsRecalibrationNSigmaPid", {std::array>{0, 0, 0}.data(), NEs, getDisplayNames()}, "nSigma PID recalibration flags"}; + Configurable> cfgFlagsRecalibrationNSigmaPid{"cfgFlagsRecalibrationNSigmaPid", {std::array>{}.data(), NEs, getDisplayNames()}, "nSigma PID recalibration flags"}; Configurable cfgFlagRejectionOthers{"cfgFlagRejectionOthers", false, "Other particle species rejection flag"}; Configurable> cfgCutsMaxAbsNSigmaPid{"cfgCutsMaxAbsNSigmaPid", {std::array>{2., 2., 2.}.data(), NEs, getDisplayNames()}, "Maximum absolute nSigma values for PID"}; Configurable cfgFlagMcParticlePhysicalPrimary{"cfgFlagMcParticlePhysicalPrimary", true, "Flag of requiring physical primary MC particle"}; Configurable cfgFlagMcParticleMomentum{"cfgFlagMcParticleMomentum", true, "Flag of using momentum of MC particle"}; } groupTrack{}; + bool doQaAcceptance{}; + bool doQaPhi{}; + bool doQaPid{}; + bool doCalculationYield{}; + bool doCalculationPurity{}; + bool doCalculationFractionPrimary{}; + bool doCalculationFluctuation{}; + HistogramRegistry hrCalculationFluctuation{"hrCalculationFluctuation", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry hrCalculationFractionPrimary{"hrCalculationFractionPrimary", {}, OutputObjHandlingPolicy::AnalysisObject}; HistogramRegistry hrCalculationPurity{"hrCalculationPurity", {}, OutputObjHandlingPolicy::AnalysisObject}; @@ -663,7 +748,7 @@ struct PartNumFluc { Produces miniMcParticle{}; Produces miniTrack{}; - void init(InitContext&) + void init(const InitContext&) { gRandom->SetSeed(0); @@ -676,6 +761,14 @@ struct PartNumFluc { LOG(info) << "Enabling raw data process."; } + doQaAcceptance = isEnabled(groupAnalysis.cfgFlagsQaAcceptance); + doQaPhi = isEnabled(groupAnalysis.cfgFlagsQaPhi); + doQaPid = isEnabled(groupAnalysis.cfgFlagsQaPid); + doCalculationYield = isEnabled(groupAnalysis.cfgFlagsCalculationYield); + doCalculationPurity = isEnabled(groupAnalysis.cfgFlagsCalculationPurity); + doCalculationFractionPrimary = isEnabled(groupAnalysis.cfgFlagsCalculationFractionPrimary); + doCalculationFluctuation = isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation); + rctFlagsChecker.init(groupEvent.cfgLabelFlagsRct.value, static_cast(groupEvent.cfgFlagsRct.value.get("ZDC")), static_cast(groupEvent.cfgFlagsRct.value.get("Acceptance")), static_cast(groupEvent.cfgFlagsRct.value.get("Table"))); ccdb->setURL(groupCcdb.cfgCcdbUrl.value); @@ -689,7 +782,6 @@ struct PartNumFluc { if (!ccdbObject || ccdbObject->IsA() != TList::Class()) { LOG(fatal) << "Invalid CCDB object!"; } - holderCcdb.clear(); std::int32_t nRunsBad{}; std::int32_t nRunGroups{}; @@ -714,10 +806,11 @@ struct PartNumFluc { } for (std::int32_t const& iRun : std::views::iota(0, gRunNumberGroupIndex->GetN())) { if (std::llrint(gRunNumberGroupIndex->GetY()[iRun]) <= 0) { - auto iter = holderCcdb.runNumbersIndicesGroupIndices.find(std::llrint(gRunNumberGroupIndex->GetX()[iRun])); - if (iter != holderCcdb.runNumbersIndicesGroupIndices.end() && iter->second.second > 0) { + if (const auto iter{holderCcdb.runNumbersIndicesGroupIndices.find(std::llrint(gRunNumberGroupIndex->GetX()[iRun]))}; iter != holderCcdb.runNumbersIndicesGroupIndices.end() && iter->second.second > 0) { iter->second.second = -iter->second.second; - ++nRunsBad; + if (groupEvent.cfgFlagRejectionRunBad.value) { + ++nRunsBad; + } } } } @@ -778,9 +871,9 @@ struct PartNumFluc { break; } - const auto readListRunGroup{[&](const std::int32_t runGroupIndex) -> const TList* { - const char* const name{Form("lRunGroup_%d", runGroupIndex)}; - const TList* const lRunGroup{dynamic_cast(ccdbObject->FindObject(name))}; + static constexpr auto ReadListRunGroup{[](const TList* const ccdbObject, const std::int32_t runGroupIndex) -> const TList* { + const std::string name{std::format("lRunGroup_{}", runGroupIndex)}; + const TList* const lRunGroup{dynamic_cast(ccdbObject->FindObject(name.c_str()))}; if (!lRunGroup) { LOG(fatal) << "Invalid " << name << "!"; } @@ -790,18 +883,27 @@ struct PartNumFluc { if (groupTrack.cfgFlagRecalibrationDca.value) { LOG(info) << "Enabling DCA recalibration."; - holderCcdb.fPtDca.resize(nRunGroups); + holderCcdb.fPtMeasureDca.resize(nRunGroups); for (std::int32_t const& iRunGroup : std::views::iota(0, nRunGroups)) { - const TList* const lRunGroup{readListRunGroup(iRunGroup + 1)}; - for (std::int32_t const& iDcaKind : std::views::iota(0, NEs)) { + const TList* const lRunGroup{ReadListRunGroup(ccdbObject, iRunGroup + 1)}; + for (std::int32_t const& iDcaMeasure : std::views::iota(0, NEs)) { for (std::int32_t const& iDcaAxis : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - const char* const name{Form("fPt%sDca%s%s%s_runGroup%d", getName(iDcaKind).data(), getName(iDcaAxis).data(), getName(iChargeSpecies).data(), doProcessMc.value ? "_mc" : "", iRunGroup + 1)}; - holderCcdb.fPtDca[iRunGroup][iDcaKind][iDcaAxis][iChargeSpecies] = dynamic_cast(lRunGroup->FindObject(name)); - if (!holderCcdb.fPtDca[iRunGroup][iDcaKind][iDcaAxis][iChargeSpecies]) { - LOG(fatal) << "Invalid " << name << "!"; + std::pair& calibration{holderCcdb.fPtMeasureDca[iRunGroup][iDcaMeasure][iDcaAxis][iChargeSpecies]}; + const std::string nameFormula{std::format("fPt{}Dca{}{}{}_runGroup{}", getName(iDcaMeasure), getName(iDcaAxis), getName(iChargeSpecies), doProcessMc.value ? "_mc" : "", iRunGroup + 1)}; + calibration.first = dynamic_cast(lRunGroup->FindObject(nameFormula.c_str())); + if (!calibration.first || calibration.first->GetNdim() != 1 || calibration.first->GetNpar() <= 0) { + LOG(fatal) << "Invalid " << nameFormula << "!"; + } + LOG(info) << "Reading from CCDB: " << nameFormula << " \"" << calibration.first->GetExpFormula() << "\""; + const std::int32_t nParameters{calibration.first->GetNpar()}; + + const std::string nameHistogram{std::format("hCentralityEtaParameterPt{}Dca{}{}{}_runGroup{}", getName(iDcaMeasure), getName(iDcaAxis), getName(iChargeSpecies), doProcessMc.value ? "_mc" : "", iRunGroup + 1)}; + calibration.second = dynamic_cast(lRunGroup->FindObject(nameHistogram.c_str())); + if (calibration.second == nullptr || calibration.second->GetNbinsZ() != nParameters || std::ranges::any_of(std::views::iota(0, nParameters), [zAxis = calibration.second->GetZaxis()](const std::int32_t binIndex) -> bool { return zAxis->GetBinCenter(binIndex + 1) != binIndex; })) { + LOG(fatal) << "Invalid " << nameHistogram << "!"; } - LOG(info) << "Reading from CCDB: " << name << " \"" << holderCcdb.fPtDca[iRunGroup][iDcaKind][iDcaAxis][iChargeSpecies]->GetExpFormula("clingp") << "\""; + LOG(info) << "Reading from CCDB: " << nameHistogram; } } } @@ -817,11 +919,11 @@ struct PartNumFluc { holderCcdb.hCentralityPtEtaShiftNSigmaPid.resize(nRunGroups); for (std::int32_t const& iRunGroup : std::views::iota(0, nRunGroups)) { - const TList* const lRunGroup{readListRunGroup(iRunGroup + 1)}; + const TList* const lRunGroup{ReadListRunGroup(ccdbObject, iRunGroup + 1)}; for (std::int32_t const& iDetector : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - const char* const name{Form("hCentralityPtEtaShift%sNSigma%s%s%s_runGroup%d", getName(iDetector).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data(), doProcessMc.value ? "_mc" : "", iRunGroup + 1)}; - holderCcdb.hCentralityPtEtaShiftNSigmaPid[iRunGroup][iDetector][iParticleSpecies][iChargeSpecies] = dynamic_cast(lRunGroup->FindObject(name)); + const std::string name{std::format("hCentralityPtEtaShift{}NSigma{}{}{}_runGroup{}", getName(iDetector), getName(iParticleSpecies), getName(iChargeSpecies), doProcessMc.value ? "_mc" : "", iRunGroup + 1)}; + holderCcdb.hCentralityPtEtaShiftNSigmaPid[iRunGroup][iDetector][iParticleSpecies][iChargeSpecies] = dynamic_cast(lRunGroup->FindObject(name.c_str())); if (!holderCcdb.hCentralityPtEtaShiftNSigmaPid[iRunGroup][iDetector][iParticleSpecies][iChargeSpecies]) { LOG(fatal) << "Invalid " << name << "!"; } @@ -839,22 +941,24 @@ struct PartNumFluc { if (groupAnalysis.cfgFlagQaRun.value) { LOG(info) << "Enabling run QA."; - const HistogramConfigSpec hcsQaRun(HistType::kTProfile, {{static_cast(holderCcdb.runNumbersIndicesGroupIndices.size()), -0.5, holderCcdb.runNumbersIndicesGroupIndices.size() - 0.5, "Run Index"}}); + const HistogramConfigSpec hcsQaRun{HistType::kTProfile, {{static_cast(holderCcdb.runNumbersIndicesGroupIndices.size()), -0.5, holderCcdb.runNumbersIndicesGroupIndices.size() - 0.5, "Run Index"}}}; for (const auto& [name, title, isChargeSpeciesSeparated] : std::to_array>( {{"Vx", "#LT#it{V}_{#it{x}}#GT (cm)", false}, {"Vy", "#LT#it{V}_{#it{y}}#GT (cm)", false}, {"Vz", "#LT#it{V}_{#it{z}}#GT (cm)", false}, + {"MultiplicityFt0a", "FT0A #LTMultiplicity#GT", false}, + {"MultiplicityFt0c", "FT0C #LTMultiplicity#GT", false}, {"CentralityNtpv", "NTPV #LTCentrality#GT", false}, {"CentralityFt0a", "FT0A #LTCentrality#GT", false}, {"CentralityFt0c", "FT0C #LTCentrality#GT", false}, {"CentralityFt0m", "FT0M #LTCentrality#GT", false}, {"NGlobalTracks", "#LTnGlobalTracks#GT", true}, {"NPvContributors", "#LTnPvContributors#GT", true}, - {Form("%sDca%s", getName(DcaKind::Mean).data(), getName(DcaAxis::Xy).data()), "#LT#LTDCA_{#it{xy}}#GT_{event}#GT (cm)", true}, - {Form("%sDca%s", getName(DcaKind::Sigma).data(), getName(DcaAxis::Xy).data()), "#LT#it{#sigma}(DCA_{#it{xy}})_{event}#GT (cm)", true}, - {Form("%sDca%s", getName(DcaKind::Mean).data(), getName(DcaAxis::Z).data()), "#LT#LTDCA_{#it{z}}#GT_{event}#GT (cm)", true}, - {Form("%sDca%s", getName(DcaKind::Sigma).data(), getName(DcaAxis::Z).data()), "#LT#it{#sigma}(DCA_{#it{z}})_{event}#GT (cm)", true}, + {std::format("{}Dca{}", getName(DcaMeasure::Mean), getName(DcaAxis::Xy)), "#LT#LTDCA_{#it{xy}}#GT_{event}#GT (cm)", true}, + {std::format("{}Dca{}", getName(DcaMeasure::Sigma), getName(DcaAxis::Xy)), "#LT#it{#sigma}(DCA_{#it{xy}})_{event}#GT (cm)", true}, + {std::format("{}Dca{}", getName(DcaMeasure::Mean), getName(DcaAxis::Z)), "#LT#LTDCA_{#it{z}}#GT_{event}#GT (cm)", true}, + {std::format("{}Dca{}", getName(DcaMeasure::Sigma), getName(DcaAxis::Z)), "#LT#it{#sigma}(DCA_{#it{z}})_{event}#GT (cm)", true}, {"NTofBeta", "#LTnTofBeta#GT", true}, {"ItsNCls", "ITS #LTnClusters#GT", true}, {"ItsChi2NCls", "ITS #LT#it{#chi}^{2}/nClusters#GT", true}, @@ -867,18 +971,18 @@ struct PartNumFluc { {"Eta", "#LT#it{#eta}#GT", true}, {"Phi", "#LT#it{#varphi}#GT", true}, {"TpcDeDx", "TPC #LTd#it{E}/d#it{x}#GT (a.u.)", true}, - {Form("%sNSigma%s", getName(Detector::Tpc).data(), getName(ParticleSpecies::Pion).data()), Form("%s #LT#it{n}#it{#sigma}_{#pi}#GT", getName(Detector::Tpc).data()), true}, - {Form("%sNSigma%s", getName(Detector::Tpc).data(), getName(ParticleSpecies::Kaon).data()), Form("%s #LT#it{n}#it{#sigma}_{K}#GT", getName(Detector::Tpc).data()), true}, - {Form("%sNSigma%s", getName(Detector::Tpc).data(), getName(ParticleSpecies::Proton).data()), Form("%s #LT#it{n}#it{#sigma}_{p}#GT", getName(Detector::Tpc).data()), true}, + {std::format("{}NSigma{}", getName(Detector::Tpc), getName(ParticleSpecies::Pion)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{#pi}}#GT", getName(Detector::Tpc)), true}, + {std::format("{}NSigma{}", getName(Detector::Tpc), getName(ParticleSpecies::Kaon)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{K}}#GT", getName(Detector::Tpc)), true}, + {std::format("{}NSigma{}", getName(Detector::Tpc), getName(ParticleSpecies::Proton)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{p}}#GT", getName(Detector::Tpc)), true}, {"TofInverseBeta", "TOF #LT1/#it{#beta}#GT", true}, - {Form("%sNSigma%s", getName(Detector::Tof).data(), getName(ParticleSpecies::Pion).data()), Form("%s #LT#it{n}#it{#sigma}_{#pi}#GT", getName(Detector::Tof).data()), true}, - {Form("%sNSigma%s", getName(Detector::Tof).data(), getName(ParticleSpecies::Kaon).data()), Form("%s #LT#it{n}#it{#sigma}_{K}#GT", getName(Detector::Tof).data()), true}, - {Form("%sNSigma%s", getName(Detector::Tof).data(), getName(ParticleSpecies::Proton).data()), Form("%s #LT#it{n}#it{#sigma}_{p}#GT", getName(Detector::Tof).data()), true}})) { + {std::format("{}NSigma{}", getName(Detector::Tof), getName(ParticleSpecies::Pion)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{#pi}}#GT", getName(Detector::Tof)), true}, + {std::format("{}NSigma{}", getName(Detector::Tof), getName(ParticleSpecies::Kaon)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{K}}#GT", getName(Detector::Tof)), true}, + {std::format("{}NSigma{}", getName(Detector::Tof), getName(ParticleSpecies::Proton)), std::format("{} #LT#it{{n}}#it{{#sigma}}_{{p}}#GT", getName(Detector::Tof)), true}})) { if (!isChargeSpeciesSeparated) { - hrQaRun.add(Form("QaRun/pRunIndex%s", name.data()), Form(";;%s", title.data()), hcsQaRun); + hrQaRun.add(std::format("QaRun/pRunIndex{}", name).c_str(), std::format(";;{}", title).c_str(), hcsQaRun); } else { - hrQaRun.add(Form("QaRun/pRunIndex%s_%s", name.data(), getName(ChargeSpecies::Plus).data()), Form(";;%s (#it{q}>0)", title.data()), hcsQaRun); - hrQaRun.add(Form("QaRun/pRunIndex%s_%s", name.data(), getName(ChargeSpecies::Minus).data()), Form(";;%s (#it{q}<0)", title.data()), hcsQaRun); + hrQaRun.add(std::format("QaRun/pRunIndex{}_{}", name, getName(ChargeSpecies::Plus)).c_str(), std::format(";;{} (#it{{q}}>0)", title).c_str(), hcsQaRun); + hrQaRun.add(std::format("QaRun/pRunIndex{}_{}", name, getName(ChargeSpecies::Minus)).c_str(), std::format(";;{} (#it{{q}}<0)", title).c_str(), hcsQaRun); } } } @@ -886,8 +990,8 @@ struct PartNumFluc { if (groupAnalysis.cfgFlagQaEvent.value) { LOG(info) << "Enabling event QA."; - const AxisSpec asNTracks(200, -0.5, 199.5); - const HistogramConfigSpec hcsQaEvent(HistType::kTHnSparseD, {asNTracks, asNTracks}); + const AxisSpec asNTracks{200, -0.5, 199.5}; + const HistogramConfigSpec hcsQaEvent{HistType::kTHnSparseD, {asNTracks, asNTracks}}; hrQaEvent.add("QaEvent/hVxVy", "", {HistType::kTHnSparseD, {{150, -0.15, 0.15, "#it{V}_{#it{x}} (cm)"}, {150, -0.15, 0.15, "#it{V}_{#it{y}} (cm)"}}}); hrQaEvent.add("QaEvent/hVz", "", {HistType::kTH1D, {{300, -15., 15., "#it{V}_{#it{z}} (cm)"}}}); @@ -917,7 +1021,7 @@ struct PartNumFluc { {"TpcChi2NCls", {HistType::kTH1D, {{100, 0., 5., "TPC #it{#chi}^{2}/nClusters"}}}}, {"TpcNClsFindableNCrossedRows", {HistType::kTHnSparseD, {{180, -0.5, 179.5, "TPC nFindableClusters"}, {180, -0.5, 179.5, "TPC nCrossedRows"}}}}})) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaTrack.add(Form("QaTrack/h%s_%s", name.data(), getName(iChargeSpecies).data()), "", configSpec); + hrQaTrack.add(std::format("QaTrack/h{}_{}", name, getName(iChargeSpecies)).c_str(), "", configSpec); } } } @@ -925,8 +1029,8 @@ struct PartNumFluc { if (groupAnalysis.cfgFlagQaDca.value) { LOG(info) << "Enabling DCA QA."; - const AxisSpec asPt(40, 0., 2., "#it{p}_{T} (GeV/#it{c})"); - const HistogramConfigSpec hcsQaDcaProfile(HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, asPt, {24, -1.2, 1.2, "#it{#eta}"}}); + const AxisSpec asPt{40, 0., 2., "#it{p}_{T} (GeV/#it{c})"}; + const HistogramConfigSpec hcsQaDcaProfile{HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, asPt, {24, -1.2, 1.2, "#it{#eta}"}}}; for (const auto& [name, title, configSpec] : std::to_array>( {{"hPtDcaXy", "", {HistType::kTHnSparseD, {asPt, {250, -0.25, 0.25, "DCA_{#it{xy}} (cm)"}}}}, @@ -934,7 +1038,7 @@ struct PartNumFluc { {"hPtDcaZ", "", {HistType::kTHnSparseD, {asPt, {250, -0.5, 0.5, "DCA_{#it{z}} (cm)"}}}}, {"pCentralityPtEtaDcaZ", ";;#LTDCA_{#it{z}}#GT (cm)", hcsQaDcaProfile}})) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaDca.add(Form("QaDca/%s_%s", name.data(), getName(iChargeSpecies).data()), title.data(), configSpec); + hrQaDca.add(std::format("QaDca/{}_{}", name, getName(iChargeSpecies)).c_str(), title.data(), configSpec); } } } @@ -948,7 +1052,7 @@ struct PartNumFluc { for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaAcceptance.add(Form("QaAcceptance/h%sPt_%sEdge%s%s", iParticleSpeciesAll == toI(ParticleSpeciesAll::All) ? "Eta" : "Rapidity", getName(iPidStrategy).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), "", {HistType::kTHnSparseD, {{300, -1.5, 1.5, iParticleSpeciesAll == toI(ParticleSpeciesAll::All) ? "#it{#eta}" : "#it{y}"}, {250, 0., 2.5, "#it{p}_{T} (GeV/#it{c})"}}}); + hrQaAcceptance.add(std::format("QaAcceptance/h{}Pt_{}Edge{}{}", iParticleSpeciesAll == toI(ParticleSpeciesAll::All) ? "Eta" : "Rapidity", getName(iPidStrategy), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), "", {HistType::kTHnSparseD, {{300, -1.5, 1.5, iParticleSpeciesAll == toI(ParticleSpeciesAll::All) ? "#it{#eta}" : "#it{y}"}, {250, 0., 2.5, "#it{p}_{T} (GeV/#it{c})"}}}); } } } @@ -960,12 +1064,11 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleSpeciesAll) << " phi QA."; - const HistogramConfigSpec hcsQaPhi(HistType::kTHnSparseF, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {24, -1.2, 1.2, "#it{#eta}"}, {360, 0., constants::math::TwoPI, "#it{#varphi} (rad)"}}); + const HistogramConfigSpec hcsQaPhi{HistType::kTHnSparseF, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {24, -1.2, 1.2, "#it{#eta}"}, {360, 0., constants::math::TwoPI, "#it{#varphi} (rad)"}}}; for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPhi.add(Form("QaPhi/hCentralityPtEtaPhi_%s%s%s", doProcessMc.value ? Form("mc%s", getName(iPidStrategy).data()) : getName(iPidStrategy).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), "", hcsQaPhi); - hrQaPhi.add(Form("QaPhi/hCentralityPtEtaPhiIu_%s%s%s", doProcessMc.value ? Form("mc%s", getName(iPidStrategy).data()) : getName(iPidStrategy).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), "", hcsQaPhi); + hrQaPhi.add(std::format("QaPhi/hCentralityPtEtaPhi_{}{}{}{}", doProcessMc.value ? "mc" : "", doProcessMc.value ? getName(iPidStrategy) : getName(iPidStrategy), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), "", hcsQaPhi); } } } @@ -977,36 +1080,35 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleSpeciesAll) << " PID QA."; - const AxisSpec asCentrality(groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"); + const AxisSpec asCentrality{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}; if (iParticleSpeciesAll == toI(ParticleSpeciesAll::All)) { - const AxisSpec asPOverQ(350, -3.5, 3.5, "#it{p}/#it{q} (GeV/#it{c})"); - const AxisSpec asEta(48, -1.2, 1.2, "#it{#eta}"); + const AxisSpec asPOverQ{350, -3.5, 3.5, "#it{p}/#it{q} (GeV/#it{c})"}; + const AxisSpec asEta{48, -1.2, 1.2, "#it{#eta}"}; hrQaPid.add("QaPid/hCentralityPOverQEtaTpcLnDeDx", "", {HistType::kTHnSparseF, {asCentrality, asPOverQ, asEta, {240, 3., 9., "TPC ln(d#it{E}/d#it{x} (a.u.))"}}}); hrQaPid.add("QaPid/hCentralityPOverQEtaTofInverseBeta", "", {HistType::kTHnSparseF, {asCentrality, asPOverQ, asEta, {120, 0.5, 3.5, "TOF 1/#it{#beta}"}}}); } else { - const HistogramConfigSpec hcsQaPid(HistType::kTHnSparseF, {asCentrality, {40, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {32, -0.8, 0.8, "#it{#eta}"}, {300, -30., 30.}}); + const HistogramConfigSpec hcsQaPid{HistType::kTHnSparseF, {asCentrality, {40, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {32, -0.8, 0.8, "#it{#eta}"}, {300, -30., 30.}}}; - constexpr std::array> ParticleSpeciesAllTitles{"", "#pi", "K", "p"}; if (doProcessMc.value) { for (std::int32_t const& iDetector : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPid.add(Form("QaPid/hCentralityPtEta%sNSigma%s_mc%s%s", getName(iDetector).data(), getName(iParticleSpeciesAll).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), Form(";;;;%s #it{n}#it{#sigma}_{%s};", getName(iDetector).data(), ParticleSpeciesAllTitles[iParticleSpeciesAll].data()), hcsQaPid); + hrQaPid.add(std::format("QaPid/hCentralityPtEta{}NSigma{}_mc{}{}", getName(iDetector), getName(iParticleSpeciesAll), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), std::format(";;;;{} #it{{n}}#it{{#sigma}}_{{{}}};", getName(iDetector), getTitle(iParticleSpeciesAll)).c_str(), hcsQaPid); } } } else { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPid.add(Form("QaPid/hCentralityPtEta%sNSigma%s_%s", getName(Detector::Tpc).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), Form(";;;;%s #it{n}#it{#sigma}_{%s};", getName(Detector::Tpc).data(), ParticleSpeciesAllTitles[iParticleSpeciesAll].data()), hcsQaPid); + hrQaPid.add(std::format("QaPid/hCentralityPtEta{}NSigma{}_{}", getName(Detector::Tpc), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), std::format(";;;;{} #it{{n}}#it{{#sigma}}_{{{}}};", getName(Detector::Tpc), getTitle(iParticleSpeciesAll)).c_str(), hcsQaPid); } for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPid.add(Form("QaPid/hCentralityPtEta%sNSigma%s_%s%s%s", getName(Detector::Tpc).data(), getName(iParticleSpeciesAll).data(), getName(Detector::Tof).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), Form(";;;;%s #it{n}#it{#sigma}_{%s};", getName(Detector::Tpc).data(), ParticleSpeciesAllTitles[iParticleSpeciesAll].data()), hcsQaPid); + hrQaPid.add(std::format("QaPid/hCentralityPtEta{}NSigma{}_{}{}{}", getName(Detector::Tpc), getName(iParticleSpeciesAll), getName(Detector::Tof), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), std::format(";;;;{} #it{{n}}#it{{#sigma}}_{{{}}};", getName(Detector::Tpc), getTitle(iParticleSpeciesAll)).c_str(), hcsQaPid); } for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPid.add(Form("QaPid/hCentralityPtEta%sNSigma%s_%s%s%s", getName(Detector::Tof).data(), getName(iParticleSpeciesAll).data(), getName(Detector::Tpc).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), Form(";;;;%s #it{n}#it{#sigma}_{%s};", getName(Detector::Tof).data(), ParticleSpeciesAllTitles[iParticleSpeciesAll].data()), hcsQaPid); + hrQaPid.add(std::format("QaPid/hCentralityPtEta{}NSigma{}_{}{}{}", getName(Detector::Tof), getName(iParticleSpeciesAll), getName(Detector::Tpc), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), std::format(";;;;{} #it{{n}}#it{{#sigma}}_{{{}}};", getName(Detector::Tof), getTitle(iParticleSpeciesAll)).c_str(), hcsQaPid); } for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrQaPid.add(Form("QaPid/hCentralityPtEta%sNSigma%s_%s", getName(PidStrategy::TpcTof).data(), getName(iParticleSpeciesAll).data(), getName(iChargeSpecies).data()), Form(";;;;%s #it{n}#it{#sigma}_{%s};", getName(PidStrategy::TpcTof).data(), ParticleSpeciesAllTitles[iParticleSpeciesAll].data()), hcsQaPid); + hrQaPid.add(std::format("QaPid/hCentralityPtEta{}NSigma{}_{}", getName(PidStrategy::TpcTof), getName(iParticleSpeciesAll), getName(iChargeSpecies)).c_str(), std::format(";;;;{} #it{{n}}#it{{#sigma}}_{{{}}};", getName(PidStrategy::TpcTof), getTitle(iParticleSpeciesAll)).c_str(), hcsQaPid); } } } @@ -1017,7 +1119,7 @@ struct PartNumFluc { LOG(info) << "Enabling MC QA."; const double maxAbsVz{std::ceil(groupEvent.cfgFlagMcCollisionVz.value ? groupEvent.cfgCutMaxAbsVzMc.value : groupEvent.cfgCutMaxAbsVz.value)}; - const AxisSpec asCentrality(20, 0., 100., "Centrality (%)"); + const AxisSpec asCentrality{20, 0., 100., "Centrality (%)"}; hrQaMc.add("QaMc/hCentralityVzMcDeltaVz", "", {HistType::kTHnSparseF, {asCentrality, {static_cast(maxAbsVz) * 20, -maxAbsVz, maxAbsVz, "#it{V}_{#it{z}}^{Gen} (cm)"}, {200, -0.2, 0.2, "#it{V}_{#it{z}}^{Rec}#minus#it{V}_{#it{z}}^{Gen} (cm)"}}}); hrQaMc.add("QaMc/hCentralityPtMcEtaMcDeltaPt", "", {HistType::kTHnSparseF, {asCentrality, {200, 0., 2., "#it{p}_{T}^{Gen} (GeV/#it{c})"}, {24, -1.2, 1.2, "#it{#eta}_{Gen}"}, {320, -0.8, 0.8, "#it{p}_{T}^{Rec}#minus#it{p}_{T}^{Gen} (GeV/#it{c})"}}}); @@ -1032,25 +1134,25 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleSpecies) << " yield calculation."; const double maxAbsVz{std::ceil(doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? groupEvent.cfgCutMaxAbsVzMc.value : groupEvent.cfgCutMaxAbsVz.value)}; - const HistogramConfigSpec hcsCalculationYield(HistType::kTHnSparseF, {{static_cast(maxAbsVz) * 2, -maxAbsVz, maxAbsVz, "#it{V}_{#it{z}} (cm)"}, {groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {40, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {32, -0.8, 0.8, "#it{#eta}"}}); + const HistogramConfigSpec hcsCalculationYield{HistType::kTHnSparseF, {{static_cast(maxAbsVz) * 2, -maxAbsVz, maxAbsVz, "#it{V}_{#it{z}} (cm)"}, {groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {40, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {32, -0.8, 0.8, "#it{#eta}"}}}; if (doProcessMc.value) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrCalculationYield.add(Form("CalculationYield/hVzCentralityPtMcEtaMc_mc%s%s", getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationYield); + hrCalculationYield.add(std::format("CalculationYield/hVzCentralityPtMcEtaMc_mc{}{}", getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationYield); } for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { if (groupTrack.cfgFlagMcParticleMomentum.value) { - hrCalculationYield.add(Form("CalculationYield/hVzCentralityPtMcEtaMc_mc%s%s%s", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationYield); + hrCalculationYield.add(std::format("CalculationYield/hVzCentralityPtMcEtaMc_mc{}{}{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationYield); } else { - hrCalculationYield.add(Form("CalculationYield/hVzCentralityPtEta_mc%s%s%s", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationYield); + hrCalculationYield.add(std::format("CalculationYield/hVzCentralityPtEta_mc{}{}{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationYield); } } } } else { for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrCalculationYield.add(Form("CalculationYield/hVzCentralityPtEta_%s%s%s", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationYield); + hrCalculationYield.add(std::format("CalculationYield/hVzCentralityPtEta_{}{}{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationYield); } } } @@ -1064,11 +1166,11 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleSpecies) << " purity calculation."; - const HistogramConfigSpec hcsCalculationPurity(HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {16, -0.8, 0.8, "#it{#eta}"}}); + const HistogramConfigSpec hcsCalculationPurity{HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {16, -0.8, 0.8, "#it{#eta}"}}}; for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrCalculationPurity.add(Form("CalculationPurity/pCentralityPtEtaPurity%s%s%s", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationPurity); + hrCalculationPurity.add(std::format("CalculationPurity/pCentralityPtEtaPurity{}{}{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationPurity); } } } @@ -1082,11 +1184,11 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleSpecies) << " primary fraction calculation."; - const HistogramConfigSpec hcsCalculationFractionPrimary(HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {16, -0.8, 0.8, "#it{#eta}"}}); + const HistogramConfigSpec hcsCalculationFractionPrimary{HistType::kTProfile3D, {{groupEvent.cfgAxisCentralityCalibration, "Centrality (%)"}, {20, 0., 2., "#it{p}_{T} (GeV/#it{c})"}, {16, -0.8, 0.8, "#it{#eta}"}}}; for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - hrCalculationFractionPrimary.add(Form("CalculationFractionPrimary/pCentralityPtEtaFractionPrimary%s%s%s", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data()), "", hcsCalculationFractionPrimary); + hrCalculationFractionPrimary.add(std::format("CalculationFractionPrimary/pCentralityPtEtaFractionPrimary{}{}{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies)).c_str(), "", hcsCalculationFractionPrimary); } } } @@ -1095,7 +1197,7 @@ struct PartNumFluc { if (nEnabled(groupAnalysis.cfgFlagsCalculationFluctuation) > 1) { LOG(fatal) << "Invalid " << groupAnalysis.cfgFlagsCalculationFluctuation.name << "!"; } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation) && groupEvent.cfgNSubgroups.value <= 0) { + if (doCalculationFluctuation && groupEvent.cfgNSubgroups.value <= 0) { LOG(fatal) << "Invalid " << groupEvent.cfgNSubgroups.name << "!"; } for (std::int32_t const& iParticleNumber : std::views::iota(0, NEs)) { @@ -1105,26 +1207,24 @@ struct PartNumFluc { LOG(info) << "Enabling " << getName(iParticleNumber) << " number fluctuation calculation."; - const AxisSpec asCentrality(groupEvent.cfgAxisCentrality, "Centrality (%)"); - const HistogramConfigSpec hcsDistribution(HistType::kTHnSparseD, {asCentrality, {200, -0.5, 199.5}, {200, -0.5, 199.5}}); - const HistogramConfigSpec hcsFluctuationCalculator(HistType::kTH3D, {asCentrality, {groupEvent.cfgNSubgroups.value, -0.5, groupEvent.cfgNSubgroups.value - 0.5, "Subgroup Index"}, {fluctuation_calculator_base::NOrderKeys, -0.5, fluctuation_calculator_base::NOrderKeys - 0.5, "Order Key Index"}}); + const AxisSpec asCentrality{groupEvent.cfgAxisCentrality, "Centrality (%)"}; + const HistogramConfigSpec hcsDistribution{HistType::kTHnSparseD, {asCentrality, {200, -0.5, 199.5}, {200, -0.5, 199.5}}}; + const HistogramConfigSpec hcsFluctuationCalculator{HistType::kTH3D, {asCentrality, {groupEvent.cfgNSubgroups.value, -0.5, groupEvent.cfgNSubgroups.value - 0.5, "Subgroup Index"}, {fluctuation_calculator_base::NOrderKeys, -0.5, fluctuation_calculator_base::NOrderKeys - 0.5, "Order Key Index"}}}; for (std::int32_t const& iChargeNumber : std::views::iota(0, NEs)) { fluctuationCalculatorTrack[iParticleNumber][iChargeNumber] = std::make_unique(); } - constexpr std::array> ParticleNumberTitles{"h", "K", "p"}; - constexpr std::array> ChargeSpeciesTitles{"+", "#minus"}; if (doProcessMc.value) { - hrCalculationFluctuation.add(Form("CalculationFluctuation/hCentralityN%s%sN%s%s_mc", getName(iParticleNumber).data(), getName(ChargeSpecies::Plus).data(), getName(iParticleNumber).data(), getName(ChargeSpecies::Minus).data()), Form(";;#it{N}(%s^{%s});#it{N}(%s^{%s});", ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Plus)].data(), ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Minus)].data()), hcsDistribution); - hrCalculationFluctuation.add(Form("CalculationFluctuation/hCentralityN%s%sN%s%s_mcEff", getName(iParticleNumber).data(), getName(ChargeSpecies::Plus).data(), getName(iParticleNumber).data(), getName(ChargeSpecies::Minus).data()), Form(";;#it{N}(%s^{%s});#it{N}(%s^{%s});", ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Plus)].data(), ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Minus)].data()), hcsDistribution); + hrCalculationFluctuation.add(std::format("CalculationFluctuation/hCentralityN{}{}N{}{}_mc", getName(iParticleNumber), getName(ChargeSpecies::Plus), getName(iParticleNumber), getName(ChargeSpecies::Minus)).c_str(), std::format(";;#it{{N}}({}^{{{}}});#it{{N}}({}^{{{}}});", getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Plus)), getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Minus))).c_str(), hcsDistribution); + hrCalculationFluctuation.add(std::format("CalculationFluctuation/hCentralityN{}{}N{}{}_mcEff", getName(iParticleNumber), getName(ChargeSpecies::Plus), getName(iParticleNumber), getName(ChargeSpecies::Minus)).c_str(), std::format(";;#it{{N}}({}^{{{}}});#it{{N}}({}^{{{}}});", getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Plus)), getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Minus))).c_str(), hcsDistribution); for (std::int32_t const& iChargeNumber : std::views::iota(0, NEs)) { - hrCalculationFluctuation.add(Form("CalculationFluctuation/hFluctuationCalculator%s%s_mc", getName(iParticleNumber).data(), getName(iChargeNumber).data()), "", hcsFluctuationCalculator); + hrCalculationFluctuation.add(std::format("CalculationFluctuation/hFluctuationCalculator{}{}_mc", getName(iParticleNumber), getName(iChargeNumber)).c_str(), "", hcsFluctuationCalculator); } } - hrCalculationFluctuation.add(Form("CalculationFluctuation/hCentralityN%s%sN%s%s", getName(iParticleNumber).data(), getName(ChargeSpecies::Plus).data(), getName(iParticleNumber).data(), getName(ChargeSpecies::Minus).data()), Form(";;#it{N}(%s^{%s});#it{N}(%s^{%s});", ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Plus)].data(), ParticleNumberTitles[iParticleNumber].data(), ChargeSpeciesTitles[toI(ChargeSpecies::Minus)].data()), hcsDistribution); + hrCalculationFluctuation.add(std::format("CalculationFluctuation/hCentralityN{}{}N{}{}", getName(iParticleNumber), getName(ChargeSpecies::Plus), getName(iParticleNumber), getName(ChargeSpecies::Minus)).c_str(), std::format(";;#it{{N}}({}^{{{}}});#it{{N}}({}^{{{}}});", getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Plus)), getTitle(iParticleNumber), getTitle(toI(ChargeSpecies::Minus))).c_str(), hcsDistribution); for (std::int32_t const& iChargeNumber : std::views::iota(0, NEs)) { - hrCalculationFluctuation.add(Form("CalculationFluctuation/hFluctuationCalculator%s%s", getName(iParticleNumber).data(), getName(iChargeNumber).data()), "", hcsFluctuationCalculator); + hrCalculationFluctuation.add(std::format("CalculationFluctuation/hFluctuationCalculator{}{}", getName(iParticleNumber), getName(iChargeNumber)).c_str(), "", hcsFluctuationCalculator); } } @@ -1135,11 +1235,11 @@ struct PartNumFluc { holderCcdb.hVzCentralityPtEtaEfficiency.resize(nRunGroups); for (std::int32_t const& iRunGroup : std::views::iota(0, nRunGroups)) { - const TList* const lRunGroup{readListRunGroup(iRunGroup + 1)}; + const TList* const lRunGroup{ReadListRunGroup(ccdbObject, iRunGroup + 1)}; for (std::int32_t const& iPidStrategy : std::views::iota(0, NEs)) { for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { - const char* const name{Form("hVzCentralityPtEtaEfficiency%s%s%s_runGroup%d", getName(iPidStrategy).data(), getName(iParticleSpecies).data(), getName(iChargeSpecies).data(), iRunGroup + 1)}; - holderCcdb.hVzCentralityPtEtaEfficiency[iRunGroup][iPidStrategy][iParticleSpecies][iChargeSpecies] = dynamic_cast(lRunGroup->FindObject(name)); + const std::string name{std::format("hVzCentralityPtEtaEfficiency{}{}{}_runGroup{}", getName(iPidStrategy), getName(iParticleSpecies), getName(iChargeSpecies), iRunGroup + 1)}; + holderCcdb.hVzCentralityPtEtaEfficiency[iRunGroup][iPidStrategy][iParticleSpecies][iChargeSpecies] = dynamic_cast(lRunGroup->FindObject(name.c_str())); if (!holderCcdb.hVzCentralityPtEtaEfficiency[iRunGroup][iPidStrategy][iParticleSpecies][iChargeSpecies] || holderCcdb.hVzCentralityPtEtaEfficiency[iRunGroup][iPidStrategy][iParticleSpecies][iChargeSpecies]->GetNdimensions() != HolderCcdb::NDimensionsEfficiency) { LOG(fatal) << "Invalid " << name << "!"; } @@ -1152,67 +1252,60 @@ struct PartNumFluc { template requires IsValid - double getEfficiency(const bool doUsingMcParticleMomentum) + double getEfficiency(const bool doUseMcParticleMomentum) const { - const THnBase* const hVzCentralityPtEtaEfficiency{holderCcdb.hVzCentralityPtEtaEfficiency.at(std::abs(holderEvent.runGroupIndex) - 1)[toI(PidStrategyValue)][toI(ParticleSpeciesValue)][toI(ChargeSpeciesValue)]}; - return hVzCentralityPtEtaEfficiency ? hVzCentralityPtEtaEfficiency->GetBinContent(hVzCentralityPtEtaEfficiency->GetBin(std::array{doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? holderMcEvent.vz : holderEvent.vz, holderEvent.centrality, doUsingMcParticleMomentum ? holderMcParticle.pt : holderTrack.pt, doUsingMcParticleMomentum ? holderMcParticle.eta : holderTrack.eta}.data())) : 0.; + const THnBase* const hVzCentralityPtEtaEfficiency{holderCcdb.hVzCentralityPtEtaEfficiency[std::abs(holderEvent.runGroupIndex) - 1][toI(PidStrategyValue)][toI(ParticleSpeciesValue)][toI(ChargeSpeciesValue)]}; + return hVzCentralityPtEtaEfficiency ? hVzCentralityPtEtaEfficiency->GetBinContent(hVzCentralityPtEtaEfficiency->GetBin(std::array{doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? holderMcEvent.vz : holderEvent.vz, holderEvent.centrality, doUseMcParticleMomentum ? holderMcParticle.pt : holderTrack.pt, doUseMcParticleMomentum ? holderMcParticle.eta : holderTrack.eta}.data())) : 0.; } - template + template requires IsValid - double getShiftNSigmaPid() + double getShiftNSigmaPid() const { - if (!groupTrack.cfgFlagsRecalibrationNSigmaPid.value.get(toI(ParticleSpeciesValue)) || holderTrack.sign == 0) { - return 0.; - } - - static const auto interpolate{[](const TH3* h, double x, double y, double z) { - if (!h) { - return 0.; + if constexpr (DoRecalibration) { + if (groupTrack.cfgFlagsRecalibrationNSigmaPid.value.get(toI(ParticleSpeciesValue))) { + return interpolate(holderCcdb.hCentralityPtEtaShiftNSigmaPid[std::abs(holderEvent.runGroupIndex) - 1][toI(DetectorValue)][toI(ParticleSpeciesValue)][holderTrack.sign > 0 ? toI(ChargeSpecies::Plus) : toI(ChargeSpecies::Minus)], holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta); } + } + return 0.; + } - const auto getBinIndicesWeights{[](const TAxis* axis, double position) -> std::array, 2> { - if (!axis) { - return {}; - } - - const std::int32_t n{axis->GetNbins()}; - if (n == 1 || position <= axis->GetBinCenter(1)) { - return {{{1, 1.}, {1, 0.}}}; - } - if (position >= axis->GetBinCenter(n)) { - return {{{n, 1.}, {n, 0.}}}; - } - - const std::int32_t bin{axis->FindFixBin(position)}; - const std::int32_t lower{position < axis->GetBinCenter(bin) ? bin - 1 : bin}; - const std::int32_t upper{lower + 1}; - const double fraction{(position - axis->GetBinCenter(lower)) / (axis->GetBinCenter(upper) - axis->GetBinCenter(lower))}; - - return {{{lower, 1. - fraction}, {upper, fraction}}}; - }}; - - const std::array, 2> xb{getBinIndicesWeights(h->GetXaxis(), x)}; - const std::array, 2> yb{getBinIndicesWeights(h->GetYaxis(), y)}; - const std::array, 2> zb{getBinIndicesWeights(h->GetZaxis(), z)}; + template + void setNSigmaPid(const T& track) + { + if (holderTrack.hasPid[toI(Detector::Tpc)]) { + holderTrack.nSigmaPid[toI(Detector::Tpc)][toI(ParticleSpecies::Pion)] = HolderTrack::truncateNSigmaPid(track.tpcNSigmaPi(), getShiftNSigmaPid()); + holderTrack.nSigmaPid[toI(Detector::Tpc)][toI(ParticleSpecies::Kaon)] = HolderTrack::truncateNSigmaPid(track.tpcNSigmaKa(), getShiftNSigmaPid()); + holderTrack.nSigmaPid[toI(Detector::Tpc)][toI(ParticleSpecies::Proton)] = HolderTrack::truncateNSigmaPid(track.tpcNSigmaPr(), getShiftNSigmaPid()); + } + if (holderTrack.hasPid[toI(Detector::Tof)]) { + holderTrack.nSigmaPid[toI(Detector::Tof)][toI(ParticleSpecies::Pion)] = HolderTrack::truncateNSigmaPid(track.tofNSigmaPi(), getShiftNSigmaPid()); + holderTrack.nSigmaPid[toI(Detector::Tof)][toI(ParticleSpecies::Kaon)] = HolderTrack::truncateNSigmaPid(track.tofNSigmaKa(), getShiftNSigmaPid()); + holderTrack.nSigmaPid[toI(Detector::Tof)][toI(ParticleSpecies::Proton)] = HolderTrack::truncateNSigmaPid(track.tofNSigmaPr(), getShiftNSigmaPid()); + } + } - double result{}; - for (const auto& [ix, wx] : xb) { - for (const auto& [iy, wy] : yb) { - for (const auto& [iz, wz] : zb) { - result += wx * wy * wz * h->GetBinContent(ix, iy, iz); - } - } - } - return result; - }}; + double getMeasureDca(const std::pair& calibration) const + { + static thread_local std::vector parametersPtMeasureDcaScratch{}; - return interpolate(holderCcdb.hCentralityPtEtaShiftNSigmaPid.at(std::abs(holderEvent.runGroupIndex) - 1)[toI(DetectorValue)][toI(ParticleSpeciesValue)][holderTrack.sign > 0 ? toI(ChargeSpecies::Plus) : toI(ChargeSpecies::Minus)], holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta); + const TFormula* const fPtMeasureDca{calibration.first}; + const TH3* const hCentralityEtaParameterPtMeasureDca{calibration.second}; + const std::int32_t nParameters{fPtMeasureDca->GetNpar()}; + if (static_cast(parametersPtMeasureDcaScratch.size()) < nParameters) { + parametersPtMeasureDcaScratch.resize(nParameters); + } + const std::int32_t centralityBinIndex{std::clamp(hCentralityEtaParameterPtMeasureDca->GetXaxis()->FindFixBin(holderEvent.centralityCalibration), 1, hCentralityEtaParameterPtMeasureDca->GetNbinsX())}; + const std::int32_t etaBinIndex{std::clamp(hCentralityEtaParameterPtMeasureDca->GetYaxis()->FindFixBin(holderTrack.eta), 1, hCentralityEtaParameterPtMeasureDca->GetNbinsY())}; + for (std::int32_t const& iParameter : std::views::iota(0, nParameters)) { + parametersPtMeasureDcaScratch[iParameter] = hCentralityEtaParameterPtMeasureDca->GetBinContent(centralityBinIndex, etaBinIndex, iParameter + 1); + } + return fPtMeasureDca->EvalPar(&holderTrack.pt, parametersPtMeasureDcaScratch.data()); } template requires IsValid - bool isPid(const bool doRejectingOthers) + bool isPid(const bool doRejectOthers) const { if constexpr (ParticleSpeciesAllValue == ParticleSpeciesAll::All) { if constexpr (PidStrategyAllValue == PidStrategyAll::Tpc) { @@ -1231,20 +1324,29 @@ struct PartNumFluc { } else { constexpr std::int32_t ParticleSpeciesIndex{toI(getValue(ParticleSpeciesAllValue))}; if constexpr (PidStrategyAllValue == PidStrategyAll::TpcTofSeparated) { - if (!(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { + if (!(std::abs(holderTrack.nSigmaPid[toI(Detector::Tpc)][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { return false; } - if (!(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { + if (!(std::abs(holderTrack.nSigmaPid[toI(Detector::Tof)][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { return false; } - if (doRejectingOthers && !(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][ParticleSpeciesIndex]) < std::min(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][(ParticleSpeciesIndex + 1) % NEs]), std::abs(holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][(ParticleSpeciesIndex + 2) % NEs])))) { + if (doRejectOthers && !(std::abs(holderTrack.nSigmaPid[toI(Detector::Tof)][ParticleSpeciesIndex]) < std::min(std::abs(holderTrack.nSigmaPid[toI(Detector::Tof)][(ParticleSpeciesIndex + 1) % NEs]), std::abs(holderTrack.nSigmaPid[toI(Detector::Tof)][(ParticleSpeciesIndex + 2) % NEs])))) { + return false; + } + } else if constexpr (PidStrategyAllValue == PidStrategyAll::TpcTofCombined) { + const double absNSigmaPidCombined{std::abs(holderTrack.getNSigmaPidCombined(ParticleSpeciesIndex))}; + if (!(absNSigmaPidCombined < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { + return false; + } + if (doRejectOthers && !(absNSigmaPidCombined < std::min(std::abs(holderTrack.getNSigmaPidCombined((ParticleSpeciesIndex + 1) % NEs)), std::abs(holderTrack.getNSigmaPidCombined((ParticleSpeciesIndex + 2) % NEs))))) { return false; } } else { - if (!(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAllValue)][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { + constexpr std::int32_t DetectorIndex{toI(getValue(PidStrategyAllValue))}; + if (!(std::abs(holderTrack.nSigmaPid[DetectorIndex][ParticleSpeciesIndex]) < groupTrack.cfgCutsMaxAbsNSigmaPid.value.get(ParticleSpeciesIndex))) { return false; } - if (doRejectingOthers && !(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAllValue)][ParticleSpeciesIndex]) < std::min(std::abs(holderTrack.nSigmaPid[toI(PidStrategyAllValue)][(ParticleSpeciesIndex + 1) % NEs]), std::abs(holderTrack.nSigmaPid[toI(PidStrategyAllValue)][(ParticleSpeciesIndex + 2) % NEs])))) { + if (doRejectOthers && !(std::abs(holderTrack.nSigmaPid[DetectorIndex][ParticleSpeciesIndex]) < std::min(std::abs(holderTrack.nSigmaPid[DetectorIndex][(ParticleSpeciesIndex + 1) % NEs]), std::abs(holderTrack.nSigmaPid[DetectorIndex][(ParticleSpeciesIndex + 2) % NEs])))) { return false; } } @@ -1254,47 +1356,29 @@ struct PartNumFluc { template requires IsValid - bool isPid() + bool isPid() const { if constexpr (ParticleSpeciesAllValue == ParticleSpeciesAll::All) { - if constexpr (ChargeSpeciesValue == ChargeSpecies::Plus) { - if (holderMcParticle.charge <= 0) { - return false; - } - } else { // ChargeSpeciesValue == ChargeSpecies::Minus - if (holderMcParticle.charge >= 0) { - return false; - } - } + return ChargeSpeciesValue == ChargeSpecies::Plus ? holderMcParticle.charge > 0 : holderMcParticle.charge < 0; } else { - if (holderMcParticle.pdgCode != getPdgCode(ParticleSpeciesAllValue, ChargeSpeciesValue)) { - return false; - } + return holderMcParticle.pdgCode == getPdgCode(getValue(ParticleSpeciesAllValue), ChargeSpeciesValue); } - return true; } - bool isGoodMomentum(const bool doUsingMcParticleMomentum) + bool isGoodMomentum(const bool doUseMcParticleMomentum) const { - if (doUsingMcParticleMomentum) { - if (!(groupTrack.cfgCutMinPt.value < holderMcParticle.pt) || !(holderMcParticle.pt < groupTrack.cfgCutMaxPt.value)) { - return false; - } - if (!(std::abs(holderMcParticle.eta) < groupTrack.cfgCutMaxAbsEta.value)) { - return false; - } - } else { - if (!(groupTrack.cfgCutMinPt.value < holderTrack.pt) || !(holderTrack.pt < groupTrack.cfgCutMaxPt.value)) { - return false; - } - if (!(std::abs(holderTrack.eta) < groupTrack.cfgCutMaxAbsEta.value)) { - return false; - } + const double pt{doUseMcParticleMomentum ? holderMcParticle.pt : holderTrack.pt}; + const double eta{doUseMcParticleMomentum ? holderMcParticle.eta : holderTrack.eta}; + if (!(groupTrack.cfgCutMinPt.value < pt) || !(pt < groupTrack.cfgCutMaxPt.value)) { + return false; + } + if (!(std::abs(eta) < groupTrack.cfgCutMaxAbsEta.value)) { + return false; } return true; } - bool isGoodDca() + bool isGoodDca() const { if (!groupTrack.cfgFlagRecalibrationDca.value) { for (std::int32_t const& iDcaAxis : std::views::iota(0, NEs)) { @@ -1303,13 +1387,12 @@ struct PartNumFluc { } } } else { - if (holderTrack.sign == 0) { - return false; - } const std::int32_t chargeSpeciesIndex{holderTrack.sign > 0 ? toI(ChargeSpecies::Plus) : toI(ChargeSpecies::Minus)}; - const std::array>, NEs>, NEs>& fPtDcaGroup{holderCcdb.fPtDca.at(std::abs(holderEvent.runGroupIndex) - 1)}; + const std::array, NEs>, NEs>, NEs>& fPtMeasureDcaGroup{holderCcdb.fPtMeasureDca[std::abs(holderEvent.runGroupIndex) - 1]}; for (std::int32_t const& iDcaAxis : std::views::iota(0, NEs)) { - if (!fPtDcaGroup[toI(DcaKind::Mean)][iDcaAxis][chargeSpeciesIndex] || !fPtDcaGroup[toI(DcaKind::Sigma)][iDcaAxis][chargeSpeciesIndex] || !(std::abs(holderTrack.dca[iDcaAxis] - fPtDcaGroup[toI(DcaKind::Mean)][iDcaAxis][chargeSpeciesIndex]->Eval(holderTrack.pt)) < groupTrack.cfgCutsMaxAbsNSigmaDca.value.get(iDcaAxis) * fPtDcaGroup[toI(DcaKind::Sigma)][iDcaAxis][chargeSpeciesIndex]->Eval(holderTrack.pt))) { + const double mean{getMeasureDca(fPtMeasureDcaGroup[toI(DcaMeasure::Mean)][iDcaAxis][chargeSpeciesIndex])}; + const double sigma{getMeasureDca(fPtMeasureDcaGroup[toI(DcaMeasure::Sigma)][iDcaAxis][chargeSpeciesIndex])}; + if (!(std::abs(holderTrack.dca[iDcaAxis] - mean) < groupTrack.cfgCutsMaxAbsNSigmaDca.value.get(iDcaAxis) * sigma)) { return false; } } @@ -1318,7 +1401,7 @@ struct PartNumFluc { } template - bool isGoodTrack(const T& track) + bool isGoodTrack(const T& track) const { if (groupTrack.cfgFlagPvContributor.value && !track.isPVContributor()) { return false; @@ -1332,7 +1415,7 @@ struct PartNumFluc { if (!(track.tpcNClsFound() > groupTrack.cfgCutMinTpcNCls.value)) { return false; } - if (!(track.tpcChi2NCl() < groupTrack.cfgCutMaxTpcChi2NCls.value)) { + if (!(groupTrack.cfgCutMinTpcChi2NCls.value < track.tpcChi2NCl()) || !(track.tpcChi2NCl() < groupTrack.cfgCutMaxTpcChi2NCls.value)) { return false; } if (!(track.tpcFractionSharedCls() < groupTrack.cfgCutMaxTpcNClsSharedRatio.value)) { @@ -1348,7 +1431,7 @@ struct PartNumFluc { } template - bool isGoodMcParticle(const MP& mcParticle) + bool isGoodMcParticle(const MP& mcParticle) const { if constexpr (IsMc) { if (!mcParticle.isPhysicalPrimary()) { @@ -1362,14 +1445,14 @@ struct PartNumFluc { requires IsValid void fillQaRunByTrackByChargeSpecies(const T& track) { - const auto fill{[&](const auto& name, const auto value) { + const auto fill{[this](const auto& name, const auto value) -> void { hrQaRun.fill(C_CS("QaRun/pRunIndex") + name + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.runIndex, value); }}; const auto fillNSigmaPidByDetectorParticleSpecies{ - [&] - requires IsValid(DetectorValue), ParticleSpeciesValue> - () { - const double nSigmaPid{holderTrack.nSigmaPid[toI(getValue(DetectorValue))][toI(ParticleSpeciesValue)]}; + [this, &fill] + requires IsValid + () -> void { + const double nSigmaPid{holderTrack.nSigmaPid[toI(DetectorValue)][toI(ParticleSpeciesValue)]}; if (std::abs(nSigmaPid) < HolderTrack::TruncationAbsNSigmaPid) { fill(C_SV(getName(DetectorValue)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesValue)), nSigmaPid); } @@ -1403,17 +1486,17 @@ struct PartNumFluc { requires IsValid void fillQaRunByEventByChargeSpecies() { - const auto fill{[&](const auto& name, const auto value) { + const auto fill{[this](const auto& name, const auto value) -> void { hrQaRun.fill(C_CS("QaRun/pRunIndex") + name + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.runIndex, value); }}; fill(C_CS("NGlobalTracks"), holderEvent.nGlobalTracks[toI(ChargeSpeciesValue)]); fill(C_CS("NPvContributors"), holderEvent.nPvContributors[toI(ChargeSpeciesValue)]); if (holderEvent.nGlobalTracks[toI(ChargeSpeciesValue)] > 0) { - fill(C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.dca[toI(DcaKind::Mean)][toI(DcaAxis::Xy)][toI(ChargeSpeciesValue)]); - fill(C_SV(getName(DcaKind::Sigma)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.dca[toI(DcaKind::Sigma)][toI(DcaAxis::Xy)][toI(ChargeSpeciesValue)]); - fill(C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.dca[toI(DcaKind::Mean)][toI(DcaAxis::Z)][toI(ChargeSpeciesValue)]); - fill(C_SV(getName(DcaKind::Sigma)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.dca[toI(DcaKind::Sigma)][toI(DcaAxis::Z)][toI(ChargeSpeciesValue)]); + fill(C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.measureDca[toI(DcaMeasure::Mean)][toI(DcaAxis::Xy)][toI(ChargeSpeciesValue)]); + fill(C_SV(getName(DcaMeasure::Sigma)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.measureDca[toI(DcaMeasure::Sigma)][toI(DcaAxis::Xy)][toI(ChargeSpeciesValue)]); + fill(C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.measureDca[toI(DcaMeasure::Mean)][toI(DcaAxis::Z)][toI(ChargeSpeciesValue)]); + fill(C_SV(getName(DcaMeasure::Sigma)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.measureDca[toI(DcaMeasure::Sigma)][toI(DcaAxis::Z)][toI(ChargeSpeciesValue)]); } fill(C_CS("NTofBeta"), holderEvent.nTofBeta[toI(ChargeSpeciesValue)]); } @@ -1422,7 +1505,7 @@ struct PartNumFluc { requires IsValid void fillQaTrackByChargeSpecies(const T& track) { - const auto fill{[&](const auto& name, const auto... positionAndWeight) { + const auto fill{[this](const auto& name, const auto... positionAndWeight) -> void { hrQaTrack.fill(C_CS("QaTrack/h") + name + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), positionAndWeight...); }}; @@ -1438,9 +1521,9 @@ struct PartNumFluc { void fillQaDcaByChargeSpecies() { const auto fillByDcaAxis{ - [&] + [this] requires IsValid - () { + () -> void { hrQaDca.fill(C_CS("QaDca/hPtDca") + C_SV(getName(DcaAxisValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderTrack.pt, holderTrack.dca[toI(DcaAxisValue)]); hrQaDca.fill(C_CS("QaDca/pCentralityPtEtaDca") + C_SV(getName(DcaAxisValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.dca[toI(DcaAxisValue)]); }}; @@ -1453,18 +1536,18 @@ struct PartNumFluc { requires IsValid void fillQaAcceptancebyParticleSpeciesAll(const T& track) { - if (!groupAnalysis.cfgFlagsQaAcceptance.value.get(toI(ParticleSpeciesAllValue)) || holderTrack.sign == 0) { + if (!groupAnalysis.cfgFlagsQaAcceptance.value.get(toI(ParticleSpeciesAllValue))) { return; } const auto fillByChargeSpecies{ - [&] + [this] requires IsValid - (const auto& name, const auto value) { + (const auto& name, const double value) -> void { const auto fillByPidStrategy{ - [&] + [this, &name, value] requires IsValid - () { + () -> void { if (isPid(PidStrategyValue), ParticleSpeciesAllValue>(false)) { hrQaAcceptance.fill(C_CS("QaAcceptance/h") + name + C_CS("Pt_") + C_SV(getName(PidStrategyValue)) + C_CS("Edge") + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), value, holderTrack.pt); // NOLINT(clang-analyzer-core.NonNullParamChecker) } @@ -1482,9 +1565,9 @@ struct PartNumFluc { } } else { if (holderTrack.sign > 0) { - fillByChargeSpecies.template operator()(C_CS("Rapidity"), track.rapidity(getMass(ParticleSpeciesAllValue))); + fillByChargeSpecies.template operator()(C_CS("Rapidity"), track.rapidity(getMass(getValue(ParticleSpeciesAllValue)))); } else { - fillByChargeSpecies.template operator()(C_CS("Rapidity"), track.rapidity(getMass(ParticleSpeciesAllValue))); + fillByChargeSpecies.template operator()(C_CS("Rapidity"), track.rapidity(getMass(getValue(ParticleSpeciesAllValue)))); } } } @@ -1497,27 +1580,21 @@ struct PartNumFluc { return; } - if (holderTrack.sign == 0) { // DataModeValue != DataMode::McMcParticle - return; - } - const auto fillByChargeSpecies{ - [&] + [this] requires IsValid - () { + () -> void { const auto fillByPidStrategy{ - [&] + [this] requires IsValid - () { + () -> void { if constexpr (DataModeValue == DataMode::McTrack) { if (isPid() && isPid(PidStrategyValue), ParticleSpeciesAllValue>(false)) { hrQaPhi.fill(C_CS("QaPhi/hCentralityPtEtaPhi_mc") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.phi); - hrQaPhi.fill(C_CS("QaPhi/hCentralityPtEtaPhiIu_mc") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.phiIu); } } else { // DataModeValue == DataMode::RawTrack if (isPid(PidStrategyValue), ParticleSpeciesAllValue>(false)) { hrQaPhi.fill(C_CS("QaPhi/hCentralityPtEtaPhi_") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.phi); - hrQaPhi.fill(C_CS("QaPhi/hCentralityPtEtaPhiIu_") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.phiIu); } } }}; @@ -1537,36 +1614,37 @@ struct PartNumFluc { requires IsValid && (DataModeValue != DataMode::McMcParticle) void fillQaPidByParticleSpeciesAll(const T& track) { - if (!groupAnalysis.cfgFlagsQaPid.value.get(toI(ParticleSpeciesAllValue)) || holderTrack.sign == 0) { + if (!groupAnalysis.cfgFlagsQaPid.value.get(toI(ParticleSpeciesAllValue))) { return; } if constexpr (ParticleSpeciesAllValue == ParticleSpeciesAll::All) { if (isPid(PidStrategy::Tpc), ParticleSpeciesAll::All>(false)) { - hrQaPid.fill(C_CS("QaPid/hCentralityPOverQEtaTpcLnDeDx"), holderEvent.centralityCalibration, holderTrack.p / holderTrack.sign, holderTrack.eta, track.tpcSignal()); + hrQaPid.fill(C_CS("QaPid/hCentralityPOverQEtaTpcLnDeDx"), holderEvent.centralityCalibration, track.p() / holderTrack.sign, holderTrack.eta, std::log(track.tpcSignal())); } if (isPid(PidStrategy::TpcTof), ParticleSpeciesAll::All>(false)) { - hrQaPid.fill(C_CS("QaPid/hCentralityPOverQEtaTofInverseBeta"), holderEvent.centralityCalibration, holderTrack.p / holderTrack.sign, holderTrack.eta, 1. / track.beta()); + hrQaPid.fill(C_CS("QaPid/hCentralityPOverQEtaTofInverseBeta"), holderEvent.centralityCalibration, track.p() / holderTrack.sign, holderTrack.eta, 1. / track.beta()); } } else { + constexpr std::int32_t ParticleSpeciesIndex{toI(getValue(ParticleSpeciesAllValue))}; const auto fillByChargeSpecies{ - [&] + [this] requires IsValid - () { + () -> void { if constexpr (DataModeValue == DataMode::McTrack) { if (isPid()) { - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_mc") + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][toI(ParticleSpeciesAllValue)]); - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_mc") + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][toI(ParticleSpeciesAllValue)]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_mc") + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(Detector::Tpc)][ParticleSpeciesIndex]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_mc") + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(Detector::Tof)][ParticleSpeciesIndex]); } } else { // DataModeValue == DataMode::RawTrack - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][toI(ParticleSpeciesAllValue)]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(Detector::Tpc)][ParticleSpeciesIndex]); if (isPid(false)) { - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(Detector::Tof)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][toI(ParticleSpeciesAllValue)]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tpc)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(Detector::Tof)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(Detector::Tpc)][ParticleSpeciesIndex]); } if (isPid(false)) { - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(Detector::Tpc)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][toI(ParticleSpeciesAllValue)]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(Detector::Tof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(Detector::Tpc)) + C_SV(getName(ParticleSpeciesAllValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(Detector::Tof)][ParticleSpeciesIndex]); } - hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(PidStrategy::TpcTof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.nSigmaPid[toI(PidStrategyAll::TpcTofCombined)][toI(ParticleSpeciesAllValue)]); + hrQaPid.fill(C_CS("QaPid/hCentralityPtEta") + C_SV(getName(PidStrategy::TpcTof)) + C_CS("NSigma") + C_SV(getName(ParticleSpeciesAllValue)) + C_CS("_") + C_SV(getName(ChargeSpeciesValue)), holderEvent.centralityCalibration, holderTrack.pt, holderTrack.eta, holderTrack.getNSigmaPidCombined(ParticleSpeciesIndex)); } }}; @@ -1586,30 +1664,32 @@ struct PartNumFluc { return; } - const std::int32_t chargeSign{[&] { + const std::int32_t chargeSign{[](const std::int32_t chargeMcParticle, const std::int32_t signTrack) constexpr -> std::int32_t { if constexpr (DataModeValue == DataMode::McMcParticle) { - return holderMcParticle.charge; + return chargeMcParticle; } else { - return holderTrack.sign; + return signTrack; + } + }(holderMcParticle.charge, holderTrack.sign)}; + if constexpr (DataModeValue == DataMode::McMcParticle) { + if (chargeSign == 0) { + return; } - }()}; - if (chargeSign == 0) { - return; } const auto fillByChargeSpecies{ - [&] + [this] requires IsValid - () { + () -> void { if constexpr (DataModeValue == DataMode::McMcParticle) { if (isPid(ParticleSpeciesValue), ChargeSpeciesValue>()) { hrCalculationYield.fill(C_CS("CalculationYield/hVzCentralityPtMcEtaMc_mc") + C_SV(getName(ParticleSpeciesValue)) + C_SV(getName(ChargeSpeciesValue)), groupEvent.cfgFlagMcCollisionVz.value ? holderMcEvent.vz : holderEvent.vz, holderEvent.centrality, holderMcParticle.pt, holderMcParticle.eta); } } else { const auto fillByPidStrategy{ - [&] + [this] requires IsValid - () { + () -> void { if constexpr (DataModeValue == DataMode::McTrack) { if (isPid(ParticleSpeciesValue), ChargeSpeciesValue>() && isPid(PidStrategyValue), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value)) { if (groupTrack.cfgFlagMcParticleMomentum.value) { @@ -1641,18 +1721,18 @@ struct PartNumFluc { requires IsValid void fillCalculationPurityByParticleSpecies() { - if (!groupAnalysis.cfgFlagsCalculationPurity.value.get(toI(ParticleSpeciesValue)) || holderTrack.sign == 0) { + if (!groupAnalysis.cfgFlagsCalculationPurity.value.get(toI(ParticleSpeciesValue))) { return; } const auto fillByChargeSpecies{ - [&] + [this] requires IsValid - () { + () -> void { const auto fillByPidStrategy{ - [&] + [this] requires IsValid - () { + () -> void { if (isPid(PidStrategyValue), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value)) { hrCalculationPurity.fill(C_CS("CalculationPurity/pCentralityPtEtaPurity") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centrality, holderTrack.pt, holderTrack.eta, isPid(ParticleSpeciesValue), ChargeSpeciesValue>() ? 1. : 0.); } @@ -1673,18 +1753,18 @@ struct PartNumFluc { requires IsValid void fillCalculationFractionPrimaryByParticleSpecies(const MP& mcParticle) { - if (!groupAnalysis.cfgFlagsCalculationFractionPrimary.value.get(toI(ParticleSpeciesValue)) || holderTrack.sign == 0) { + if (!groupAnalysis.cfgFlagsCalculationFractionPrimary.value.get(toI(ParticleSpeciesValue))) { return; } const auto fillByChargeSpecies{ - [&] + [this, &mcParticle] requires IsValid - () { + () -> void { const auto fillByPidStrategy{ - [&] + [this, &mcParticle] requires IsValid - () { + () -> void { if (isPid(ParticleSpeciesValue), ChargeSpeciesValue>() && isPid(PidStrategyValue), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value)) { hrCalculationFractionPrimary.fill(C_CS("CalculationFractionPrimary/pCentralityPtEtaFractionPrimary") + C_SV(getName(PidStrategyValue)) + C_SV(getName(ParticleSpeciesValue)) + C_SV(getName(ChargeSpeciesValue)), holderEvent.centrality, holderTrack.pt, holderTrack.eta, mcParticle.isPhysicalPrimary() ? 1. : 0.); } @@ -1706,7 +1786,7 @@ struct PartNumFluc { for (std::int32_t const& iParticleNumber : std::views::iota(0, NEs)) { if (static_cast(groupAnalysis.cfgFlagsCalculationFluctuation.value.get(iParticleNumber))) { for (std::int32_t const& iChargeNumber : std::views::iota(0, NEs)) { - fluctuationCalculatorTrack[iParticleNumber][iChargeNumber]->init(); + fluctuationCalculatorTrack[iParticleNumber][iChargeNumber]->clear(); } } } @@ -1720,27 +1800,29 @@ struct PartNumFluc { return; } - const std::int32_t chargeSign{[&] { + const std::int32_t chargeSign{[](const std::int32_t chargeMcParticle, const std::int32_t signTrack) constexpr -> std::int32_t { if constexpr (DataModeValue == DataMode::McMcParticle) { - return holderMcParticle.charge; + return chargeMcParticle; } else { - return holderTrack.sign; + return signTrack; + } + }(holderMcParticle.charge, holderTrack.sign)}; + if constexpr (DataModeValue == DataMode::McMcParticle) { + if (chargeSign == 0) { + return; } - }()}; - if (chargeSign == 0) { - return; } - const bool doUsingMcParticleMomentum{[&] { + const bool doUseMcParticleMomentum{[](const bool flagMcParticleMomentum) constexpr -> bool { if constexpr (DataModeValue == DataMode::McMcParticle) { return true; } else if constexpr (DataModeValue == DataMode::McTrack) { - return groupTrack.cfgFlagMcParticleMomentum.value; - } else { // DataModeValue == DataMode::RawTrack + return flagMcParticleMomentum; + } else { return false; } - }()}; - if (!isGoodMomentum(doUsingMcParticleMomentum)) { + }(groupTrack.cfgFlagMcParticleMomentum.value)}; + if (!isGoodMomentum(doUseMcParticleMomentum)) { return; } @@ -1751,39 +1833,37 @@ struct PartNumFluc { } const auto calculateByParticleSpecies{ - [&] + [this, chargeSign, doUseMcParticleMomentum] requires IsValid && (ParticleNumberValue == ParticleNumber::Charge || (ParticleNumberValue == ParticleNumber::Kaon && ParticleSpeciesValue == ParticleSpecies::Kaon) || (ParticleNumberValue == ParticleNumber::Proton && ParticleSpeciesValue == ParticleSpecies::Proton)) - () { + () -> void { const auto calculateByChargeSpecies{ - [&] + [this, chargeSign, doUseMcParticleMomentum] requires IsValid - () { + () -> void { if constexpr (DataModeValue != DataMode::RawTrack) { if (!isPid(ParticleSpeciesValue), ChargeSpeciesValue>()) { return; } } - const bool doUsingTofPid{[&] { - return [&] { - if constexpr (DataModeValue == DataMode::McMcParticle) { - return holderMcParticle.pt; - } else if constexpr (DataModeValue == DataMode::McTrack) { - return groupTrack.cfgFlagMcParticleMomentum.value ? holderMcParticle.pt : holderTrack.pt; - } else { // DataModeValue == DataMode::RawTrack - return holderTrack.pt; - } - }() >= groupTrack.cfgThresholdsPtTofPid.value.get(toI(ParticleSpeciesValue)); - }()}; + const bool doUseTofPid{[](const bool doUseMcParticleMomentumValue, const double ptMcParticle, const double ptTrack, const double thresholdPtTofPid) constexpr -> bool { + if constexpr (DataModeValue == DataMode::McMcParticle) { + return ptMcParticle >= thresholdPtTofPid; + } else if constexpr (DataModeValue == DataMode::McTrack) { + return (doUseMcParticleMomentumValue ? ptMcParticle : ptTrack) >= thresholdPtTofPid; + } else { + return ptTrack >= thresholdPtTofPid; + } + }(doUseMcParticleMomentum, holderMcParticle.pt, holderTrack.pt, groupTrack.cfgThresholdsPtTofPid.value.get(toI(ParticleSpeciesValue)))}; if constexpr (DataModeValue != DataMode::McMcParticle) { - if (!(doUsingTofPid ? isPid(PidStrategy::TpcTof), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value) : isPid(PidStrategy::Tpc), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value))) { + if (!(doUseTofPid ? isPid(PidStrategy::TpcTof), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value) : isPid(PidStrategy::Tpc), getValue(ParticleSpeciesValue)>(groupTrack.cfgFlagRejectionOthers.value))) { return; } } - const double efficiency{doUsingTofPid ? getEfficiency(doUsingMcParticleMomentum) : getEfficiency(doUsingMcParticleMomentum)}; // NOLINT(clang-analyzer-core.NullDereference) + const double efficiency{doUseTofPid ? getEfficiency(doUseMcParticleMomentum) : getEfficiency(doUseMcParticleMomentum)}; // NOLINT(clang-analyzer-core.NullDereference) const auto fill{ - [&] { + [this, efficiency]() -> void { if constexpr (ChargeSpeciesValue == ChargeSpecies::Plus) { fluctuationCalculatorTrack[toI(ParticleNumberValue)][toI(ChargeNumber::Plus)]->fill(1., efficiency); fluctuationCalculatorTrack[toI(ParticleNumberValue)][toI(ChargeNumber::Net)]->fill(1., efficiency); @@ -1841,14 +1921,15 @@ struct PartNumFluc { } const auto fillByChargeNumber{ - [&] + [this] requires IsValid - () { + () -> void { + const std::array products{fluctuationCalculatorTrack[toI(ParticleNumberValue)][toI(ChargeNumberValue)]->getProducts()}; for (std::int32_t const& iOrderKey : std::views::iota(0, fluctuation_calculator_base::NOrderKeys)) { if constexpr (DataModeValue == DataMode::McMcParticle) { - hrCalculationFluctuation.fill(C_CS("CalculationFluctuation/hFluctuationCalculator") + C_SV(getName(ParticleNumberValue)) + C_SV(getName(ChargeNumberValue)) + C_CS("_mc"), holderEvent.centrality, holderEvent.subgroupIndex, iOrderKey, fluctuationCalculatorTrack[toI(ParticleNumberValue)][toI(ChargeNumberValue)]->getProduct(iOrderKey)); + hrCalculationFluctuation.fill(C_CS("CalculationFluctuation/hFluctuationCalculator") + C_SV(getName(ParticleNumberValue)) + C_SV(getName(ChargeNumberValue)) + C_CS("_mc"), holderEvent.centrality, holderEvent.subgroupIndex, iOrderKey, products[iOrderKey]); } else { - hrCalculationFluctuation.fill(C_CS("CalculationFluctuation/hFluctuationCalculator") + C_SV(getName(ParticleNumberValue)) + C_SV(getName(ChargeNumberValue)), holderEvent.centrality, holderEvent.subgroupIndex, iOrderKey, fluctuationCalculatorTrack[toI(ParticleNumberValue)][toI(ChargeNumberValue)]->getProduct(iOrderKey)); + hrCalculationFluctuation.fill(C_CS("CalculationFluctuation/hFluctuationCalculator") + C_SV(getName(ParticleNumberValue)) + C_SV(getName(ChargeNumberValue)), holderEvent.centrality, holderEvent.subgroupIndex, iOrderKey, products[iOrderKey]); } } }}; @@ -1859,110 +1940,86 @@ struct PartNumFluc { fillByChargeNumber.template operator()(); } - template - bool initTrack(const T& track, const TIs& tracksIu) + template + bool initTrack(const T& track) { + if (std::abs(track.sign()) != 1) { + return false; + } + holderTrack.clear(); holderTrack.dca[toI(DcaAxis::Xy)] = track.dcaXY(); holderTrack.dca[toI(DcaAxis::Z)] = track.dcaZ(); holderTrack.sign = track.sign(); - holderTrack.p = track.p(); holderTrack.pt = track.pt(); holderTrack.eta = track.eta(); holderTrack.phi = track.phi(); - { - const std::int64_t localIndexTrackIu{track.globalIndex() - static_cast(tracksIu.offset())}; - if (0 <= localIndexTrackIu && localIndexTrackIu < tracksIu.size()) { - const auto& trackIu{tracksIu.iteratorAt(localIndexTrackIu)}; - if (track.globalIndex() == trackIu.globalIndex()) { - holderTrack.phiIu = trackIu.phi(); - } else { - LOG(warning) << "Mismatched track " << track.globalIndex() << " and trackIu " << trackIu.globalIndex(); - } - } else { - LOG(warning) << "Invalid trackIu " << track.globalIndex(); - } - } holderTrack.hasPid[toI(Detector::Tpc)] = (track.hasTPC() && track.tpcSignal() > 0.); - if (holderTrack.hasPid[toI(Detector::Tpc)]) { - holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)] = {HolderTrack::truncateNSigmaPid(track.tpcNSigmaPi() - getShiftNSigmaPid()), HolderTrack::truncateNSigmaPid(track.tpcNSigmaKa() - getShiftNSigmaPid()), HolderTrack::truncateNSigmaPid(track.tpcNSigmaPr() - getShiftNSigmaPid())}; - } holderTrack.hasPid[toI(Detector::Tof)] = (track.hasTOF() && track.beta() > 0.); - if (holderTrack.hasPid[toI(Detector::Tof)]) { - holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)] = {HolderTrack::truncateNSigmaPid(track.tofNSigmaPi() - getShiftNSigmaPid()), HolderTrack::truncateNSigmaPid(track.tofNSigmaKa() - getShiftNSigmaPid()), HolderTrack::truncateNSigmaPid(track.tofNSigmaPr() - getShiftNSigmaPid())}; - if (holderTrack.hasPid[toI(Detector::Tpc)]) { - for (std::int32_t const& iParticleSpecies : std::views::iota(0, NEs)) { - holderTrack.nSigmaPid[toI(PidStrategyAll::TpcTofCombined)][iParticleSpecies] = HolderTrack::truncateNSigmaPid(std::copysign(std::hypot(holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][iParticleSpecies], holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][iParticleSpecies]), holderTrack.nSigmaPid[toI(PidStrategyAll::Tpc)][iParticleSpecies] + holderTrack.nSigmaPid[toI(PidStrategyAll::Tof)][iParticleSpecies])); - } - } - } + setNSigmaPid(track); if constexpr (DoInitEvent) { - if (track.isPrimaryTrack() && holderTrack.sign != 0) { + if (track.isPrimaryTrack()) { const std::int32_t chargeSpeciesIndex{toI(holderTrack.sign > 0 ? ChargeSpecies::Plus : ChargeSpecies::Minus)}; ++holderEvent.nGlobalTracks[chargeSpeciesIndex]; if (track.isPVContributor()) { ++holderEvent.nPvContributors[chargeSpeciesIndex]; } for (std::int32_t const& iDcaAxis : std::views::iota(0, NEs)) { - holderEvent.dca[toI(DcaKind::Mean)][iDcaAxis][chargeSpeciesIndex] += holderTrack.dca[iDcaAxis]; - holderEvent.dca[toI(DcaKind::Sigma)][iDcaAxis][chargeSpeciesIndex] += std::pow(holderTrack.dca[iDcaAxis], 2.); + holderEvent.measureDca[toI(DcaMeasure::Mean)][iDcaAxis][chargeSpeciesIndex] += holderTrack.dca[iDcaAxis]; + holderEvent.measureDca[toI(DcaMeasure::Sigma)][iDcaAxis][chargeSpeciesIndex] += std::pow(holderTrack.dca[iDcaAxis], 2.); } if (holderTrack.hasPid[toI(Detector::Tof)]) { ++holderEvent.nTofBeta[chargeSpeciesIndex]; } } - } - if constexpr (DoInitEvent) { if (groupAnalysis.cfgFlagQaRun.value && track.isPrimaryTrack()) { if (holderTrack.sign > 0) { fillQaRunByTrackByChargeSpecies(track); - } else if (holderTrack.sign < 0) { + } else { fillQaRunByTrackByChargeSpecies(track); } } - } - if constexpr (!DoInitEvent) { + return true; + } else { if (groupAnalysis.cfgFlagQaTrack.value && track.isPrimaryTrack()) { if (holderTrack.sign > 0) { fillQaTrackByChargeSpecies(track); - } else if (holderTrack.sign < 0) { + } else { fillQaTrackByChargeSpecies(track); } } - } - if (!isGoodTrack(track)) { - return false; - } + if (!isGoodTrack(track)) { + return false; + } - if constexpr (!DoInitEvent) { if (groupAnalysis.cfgFlagQaDca.value) { if (holderTrack.sign > 0) { fillQaDcaByChargeSpecies(); - } else if (holderTrack.sign < 0) { + } else { fillQaDcaByChargeSpecies(); } } - } - if (!isGoodDca()) { - return false; - } + if (!isGoodDca()) { + return false; + } - if constexpr (!DoInitEvent) { - const double vz{doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? holderMcEvent.vz : holderEvent.vz}; - if (isEnabled(groupAnalysis.cfgFlagsQaAcceptance) && (holderTrack.eta * vz > 0. && std::abs(vz) > (doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? groupEvent.cfgCutMaxAbsVzMc.value : groupEvent.cfgCutMaxAbsVz.value) - 1.)) { - fillQaAcceptancebyParticleSpeciesAll(track); - fillQaAcceptancebyParticleSpeciesAll(track); - fillQaAcceptancebyParticleSpeciesAll(track); - fillQaAcceptancebyParticleSpeciesAll(track); + { + const double vz{doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? holderMcEvent.vz : holderEvent.vz}; + if (doQaAcceptance && (holderTrack.eta * vz > 0. && std::abs(vz) > (doProcessMc.value && groupEvent.cfgFlagMcCollisionVz.value ? groupEvent.cfgCutMaxAbsVzMc.value : groupEvent.cfgCutMaxAbsVz.value) - 1.)) { + fillQaAcceptancebyParticleSpeciesAll(track); + fillQaAcceptancebyParticleSpeciesAll(track); + fillQaAcceptancebyParticleSpeciesAll(track); + fillQaAcceptancebyParticleSpeciesAll(track); + } } - } - return true; + return true; + } } template @@ -1988,8 +2045,8 @@ struct PartNumFluc { return isGoodMcParticle(mcParticle); } - template - bool initEvent(const C& collision, const Ts& tracks, const TIs& tracksIu) + template + bool initEvent(const C& collision, const Ts& tracks) { holderEvent.clear(); holderEvent.vz = collision.posZ(); @@ -2006,65 +2063,56 @@ struct PartNumFluc { } holderEvent.centralityCalibration = (groupEvent.cfgFlagDefinitionCentralitySameQa.value ? holderEvent.centrality : collision.centNTPV()); - hrCounter.fill(C_CS("hNEvents"), 0.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 0.); - } - - if (!collision.has_foundBC()) { - hrCounter.fill(C_CS("hNEvents"), 2.); + const auto fillEventSelection{[this](const auto... selections) -> void { + (hrCounter.fill(C_CS("hNEvents"), selections), ...); if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 2.); + (hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, selections), ...); } + }}; + + fillEventSelection(0.); + + if (!collision.has_foundBC()) { + fillEventSelection(2.); return false; } - const auto& bc{collision.template bc_as()}; - holderEvent.runNumber = bc.runNumber(); + const auto& foundBc{collision.template foundBC_as()}; + holderEvent.runNumber = foundBc.runNumber(); - if (!holderCcdb.runNumbersIndicesGroupIndices.contains(holderEvent.runNumber)) { - hrCounter.fill(C_CS("hNEvents"), 2.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 2.); + { + const auto iter{holderCcdb.runNumbersIndicesGroupIndices.find(holderEvent.runNumber)}; + if (iter == holderCcdb.runNumbersIndicesGroupIndices.end()) { + fillEventSelection(2.); + return false; } - return false; + std::tie(holderEvent.runIndex, holderEvent.runGroupIndex) = iter->second; } - std::tie(holderEvent.runIndex, holderEvent.runGroupIndex) = holderCcdb.runNumbersIndicesGroupIndices.at(holderEvent.runNumber); - if (holderEvent.runGroupIndex == 0 || (groupEvent.cfgFlagRejectionRunBad.value && holderEvent.runGroupIndex < 0)) { - hrCounter.fill(C_CS("hNEvents"), 2.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 2.); - } + fillEventSelection(2.); return false; } - if (!rctFlagsChecker.checkTable(collision)) { - hrCounter.fill(C_CS("hNEvents"), 3.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 3.); - } + if (rctFlagsChecker.any() && !rctFlagsChecker.checkTable(collision)) { + fillEventSelection(3.); + return false; + } + + if (groupEvent.cfgFlagInelEvent.value && !collision.isInelGt0()) { + fillEventSelection(4.); return false; } for (std::int32_t const& iEvSel : std::views::iota(0, aod::evsel::EventSelectionFlags::kNsel)) { if (((groupEvent.cfgBitsSelectionEvent.value >> iEvSel) & 1) && !collision.selection_bit(iEvSel)) { - hrCounter.fill(C_CS("hNEvents"), 4.); - hrCounter.fill(C_CS("hNEvents"), 10. + iEvSel); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 4.); - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 10. + iEvSel); - } + fillEventSelection(5., 10. + iEvSel); return false; } } - if (groupEvent.cfgFlagInelEvent.value && !collision.isInelGt0()) { - hrCounter.fill(C_CS("hNEvents"), 5.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 5.); - } + if (!HolderEvent::isValidCentrality(holderEvent.centrality) || !HolderEvent::isValidCentrality(holderEvent.centralityCalibration)) { + fillEventSelection(6.); return false; } @@ -2074,10 +2122,7 @@ struct PartNumFluc { } if (!(std::abs(holderEvent.vz) < groupEvent.cfgCutMaxAbsVz.value)) { - hrCounter.fill(C_CS("hNEvents"), 6.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 6.); - } + fillEventSelection(7.); return false; } @@ -2085,32 +2130,30 @@ struct PartNumFluc { hrQaRun.fill(C_CS("QaRun/pRunIndexVx"), holderEvent.runIndex, collision.posX()); hrQaRun.fill(C_CS("QaRun/pRunIndexVy"), holderEvent.runIndex, collision.posY()); hrQaRun.fill(C_CS("QaRun/pRunIndexVz"), holderEvent.runIndex, holderEvent.vz); - if (const double centrality = collision.centNTPV(); HolderEvent::RangeCentrality.first <= centrality && centrality <= HolderEvent::RangeCentrality.second) { + hrQaRun.fill(C_CS("QaRun/pRunIndexMultiplicityFt0a"), holderEvent.runIndex, collision.multZeqFT0A()); + hrQaRun.fill(C_CS("QaRun/pRunIndexMultiplicityFt0c"), holderEvent.runIndex, collision.multZeqFT0C()); + if (const double centrality{collision.centNTPV()}; HolderEvent::isValidCentrality(centrality)) { hrQaRun.fill(C_CS("QaRun/pRunIndexCentralityNtpv"), holderEvent.runIndex, centrality); } - if (const double centrality = collision.centFT0A(); HolderEvent::RangeCentrality.first <= centrality && centrality <= HolderEvent::RangeCentrality.second) { + if (const double centrality{collision.centFT0A()}; HolderEvent::isValidCentrality(centrality)) { hrQaRun.fill(C_CS("QaRun/pRunIndexCentralityFt0a"), holderEvent.runIndex, centrality); } - if (const double centrality = collision.centFT0C(); HolderEvent::RangeCentrality.first <= centrality && centrality <= HolderEvent::RangeCentrality.second) { + if (const double centrality{collision.centFT0C()}; HolderEvent::isValidCentrality(centrality)) { hrQaRun.fill(C_CS("QaRun/pRunIndexCentralityFt0c"), holderEvent.runIndex, centrality); } - if (const double centrality = collision.centFT0M(); HolderEvent::RangeCentrality.first <= centrality && centrality <= HolderEvent::RangeCentrality.second) { + if (const double centrality{collision.centFT0M()}; HolderEvent::isValidCentrality(centrality)) { hrQaRun.fill(C_CS("QaRun/pRunIndexCentralityFt0m"), holderEvent.runIndex, centrality); } } for (const auto& track : tracks) { - if (!track.has_collision() || track.collisionId() != collision.globalIndex()) { - continue; - } - - initTrack(track, tracksIu); + initTrack(track); } for (std::int32_t const& iChargeSpecies : std::views::iota(0, NEs)) { if (holderEvent.nGlobalTracks[iChargeSpecies] > 0) { for (std::int32_t const& iDcaAxis : std::views::iota(0, NEs)) { - holderEvent.dca[toI(DcaKind::Mean)][iDcaAxis][iChargeSpecies] /= holderEvent.nGlobalTracks[iChargeSpecies]; - holderEvent.dca[toI(DcaKind::Sigma)][iDcaAxis][iChargeSpecies] = std::sqrt(std::max(0., holderEvent.dca[toI(DcaKind::Sigma)][iDcaAxis][iChargeSpecies] / holderEvent.nGlobalTracks[iChargeSpecies] - std::pow(holderEvent.dca[toI(DcaKind::Mean)][iDcaAxis][iChargeSpecies], 2.))); + holderEvent.measureDca[toI(DcaMeasure::Mean)][iDcaAxis][iChargeSpecies] /= holderEvent.nGlobalTracks[iChargeSpecies]; + holderEvent.measureDca[toI(DcaMeasure::Sigma)][iDcaAxis][iChargeSpecies] = std::sqrt(std::max(0., holderEvent.measureDca[toI(DcaMeasure::Sigma)][iDcaAxis][iChargeSpecies] / holderEvent.nGlobalTracks[iChargeSpecies] - std::pow(holderEvent.measureDca[toI(DcaMeasure::Mean)][iDcaAxis][iChargeSpecies], 2.))); } } } @@ -2123,29 +2166,23 @@ struct PartNumFluc { if (groupAnalysis.cfgFlagQaEvent.value) { hrQaEvent.fill(C_CS("QaEvent/hNPvContributorsNGlobalTracks"), holderEvent.getNPvContributors(), holderEvent.getNGlobalTracks()); if (holderEvent.getNGlobalTracks() > 0) { - hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.getNGlobalTracks(), holderEvent.getDca()); - hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.getNGlobalTracks(), holderEvent.getDca()); + hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)), holderEvent.getNGlobalTracks(), holderEvent.getMeasureDca()); + hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)), holderEvent.getNGlobalTracks(), holderEvent.getMeasureDca()); } hrQaEvent.fill(C_CS("QaEvent/hNTofBetaNGlobalTracks"), holderEvent.getNTofBeta(), holderEvent.getNGlobalTracks()); } if (!(holderEvent.getNPvContributors() - holderEvent.getNGlobalTracks() > groupEvent.cfgCutMinDeviationNPvContributors.value)) { - hrCounter.fill(C_CS("hNEvents"), 7.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 7.); - } + fillEventSelection(8.); return false; } - hrCounter.fill(C_CS("hNEvents"), 1.); - if (groupAnalysis.cfgFlagQaCentrality.value) { - hrQaCentrality.fill(C_CS("QaCentrality/hCentralitySelection"), holderEvent.centrality, 1.); - } + fillEventSelection(1.); if (groupAnalysis.cfgFlagQaEvent.value) { if (holderEvent.getNGlobalTracks() > 0) { - hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)) + C_CS("_nPvContributorsCut"), holderEvent.getNGlobalTracks(), holderEvent.getDca()); - hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaKind::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)) + C_CS("_nPvContributorsCut"), holderEvent.getNGlobalTracks(), holderEvent.getDca()); + hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Xy)) + C_CS("_nPvContributorsCut"), holderEvent.getNGlobalTracks(), holderEvent.getMeasureDca()); + hrQaEvent.fill(C_CS("QaEvent/hNGlobalTracks") + C_SV(getName(DcaMeasure::Mean)) + C_CS("Dca") + C_SV(getName(DcaAxis::Z)) + C_CS("_nPvContributorsCut"), holderEvent.getNGlobalTracks(), holderEvent.getMeasureDca()); } hrQaEvent.fill(C_CS("QaEvent/hNTofBetaNGlobalTracks_nPvContributorsCut"), holderEvent.getNTofBeta(), holderEvent.getNGlobalTracks()); } @@ -2165,38 +2202,18 @@ struct PartNumFluc { hrCounter.fill(C_CS("hNMcEvents"), 0.); - if (!mcCollision.has_bc()) { - hrCounter.fill(C_CS("hNMcEvents"), 2.); - return false; - } - - const auto& bc{mcCollision.template bc_as()}; - holderMcEvent.runNumber = bc.runNumber(); - - if (!holderCcdb.runNumbersIndicesGroupIndices.contains(holderMcEvent.runNumber)) { - hrCounter.fill(C_CS("hNMcEvents"), 2.); - return false; - } - - std::tie(holderMcEvent.runIndex, holderMcEvent.runGroupIndex) = holderCcdb.runNumbersIndicesGroupIndices.at(holderMcEvent.runNumber); - - if (holderMcEvent.runGroupIndex == 0 || (groupEvent.cfgFlagRejectionRunBad.value && holderMcEvent.runGroupIndex < 0)) { - hrCounter.fill(C_CS("hNMcEvents"), 2.); - return false; - } - if (groupEvent.cfgFlagInelEventMc.value && !mcCollision.isInelGt0()) { - hrCounter.fill(C_CS("hNMcEvents"), 3.); + hrCounter.fill(C_CS("hNMcEvents"), 2.); return false; } if (!(std::abs(holderMcEvent.vz) < groupEvent.cfgCutMaxAbsVzMc.value)) { - hrCounter.fill(C_CS("hNMcEvents"), 4.); + hrCounter.fill(C_CS("hNMcEvents"), 3.); return false; } if (groupEvent.cfgFlagSingleCollisionMc.value && mcCollision.numRecoCollision() != 1) { - hrCounter.fill(C_CS("hNMcEvents"), 5.); + hrCounter.fill(C_CS("hNMcEvents"), 4.); return false; } @@ -2205,55 +2222,51 @@ struct PartNumFluc { return true; } - void processRaw(const soa::Filtered::iterator& collision, const soa::Filtered& tracks, const aod::TracksIU& tracksIu, const aod::BCsWithTimestamps&) + void processRaw(const soa::Filtered::iterator& collision, const soa::Filtered& tracks, const aod::BCsWithTimestamps&) { - if (!initEvent(collision, tracks, tracksIu) || (!groupAnalysis.cfgFlagQaTrack.value && !groupAnalysis.cfgFlagQaDca.value && !isEnabled(groupAnalysis.cfgFlagsQaAcceptance) && !isEnabled(groupAnalysis.cfgFlagsQaPhi) && !isEnabled(groupAnalysis.cfgFlagsQaPid) && !isEnabled(groupAnalysis.cfgFlagsCalculationYield) && !isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation))) { + if (!initEvent(collision, tracks) || (!groupAnalysis.cfgFlagQaTrack.value && !groupAnalysis.cfgFlagQaDca.value && !doQaAcceptance && !doQaPhi && !doQaPid && !doCalculationYield && !doCalculationFluctuation)) { return; } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { holderEvent.subgroupIndex = gRandom->Integer(groupEvent.cfgNSubgroups.value); initCalculationFluctuation(); holderDerivedData.clear(); } for (const auto& track : tracks) { - if (!track.has_collision() || track.collisionId() != collision.globalIndex()) { + if (!initTrack(track)) { continue; } - if (!initTrack(track, tracksIu)) { - continue; - } - - if (isEnabled(groupAnalysis.cfgFlagsQaPhi)) { + if (doQaPhi) { fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); } - if (isEnabled(groupAnalysis.cfgFlagsQaPid)) { + if (doQaPid) { fillQaPidByParticleSpeciesAll(track); fillQaPidByParticleSpeciesAll(track); fillQaPidByParticleSpeciesAll(track); fillQaPidByParticleSpeciesAll(track); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationYield)) { + if (doCalculationYield) { fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); } } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber(); @@ -2264,7 +2277,7 @@ struct PartNumFluc { } } - void processMc(const soa::Filtered::iterator& mcCollision, const aod::McParticles& mcParticles, const soa::SmallGroups& collisions, const soa::Filtered& tracksUngrouped, const aod::TracksIU& tracksIuUngrouped, const aod::BCsWithTimestamps&) + void processMc(const soa::Filtered::iterator& mcCollision, const aod::McParticles& mcParticles, const soa::SmallGroups& collisions, const soa::Filtered& tracksUngrouped, const aod::BCsWithTimestamps&) { if (!initMcEvent(mcCollision)) { return; @@ -2276,9 +2289,8 @@ struct PartNumFluc { } const auto& tracks{tracksUngrouped.sliceBy(presliceTracksPerCollision, collision.globalIndex())}; - const auto& tracksIu{tracksIuUngrouped.sliceBy(presliceTracksPerCollision, collision.globalIndex())}; - if (!initEvent(collision, tracks, tracksIu)) { + if (!initEvent(collision, tracks)) { continue; } @@ -2286,49 +2298,45 @@ struct PartNumFluc { hrQaMc.fill(C_CS("QaMc/hCentralityVzMcDeltaVz"), holderEvent.centralityCalibration, holderMcEvent.vz, holderEvent.vz - holderMcEvent.vz); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationYield) || isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationYield || doCalculationFluctuation) { + if (doCalculationFluctuation) { holderEvent.subgroupIndex = gRandom->Integer(groupEvent.cfgNSubgroups.value); initCalculationFluctuation(); holderDerivedData.clear(); } for (const auto& mcParticle : mcParticles) { - if (!mcParticle.has_mcCollision() || mcParticle.mcCollisionId() != mcCollision.globalIndex()) { - continue; - } - if (!initMcParticle(mcParticle)) { continue; } - if (isEnabled(groupAnalysis.cfgFlagsCalculationYield)) { + if (doCalculationYield) { fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); } } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber(); } } - if (groupAnalysis.cfgFlagQaTrack.value || groupAnalysis.cfgFlagQaDca.value || isEnabled(groupAnalysis.cfgFlagsQaAcceptance) || isEnabled(groupAnalysis.cfgFlagsQaPhi) || isEnabled(groupAnalysis.cfgFlagsQaPid) || isEnabled(groupAnalysis.cfgFlagsCalculationYield) || isEnabled(groupAnalysis.cfgFlagsCalculationPurity) || isEnabled(groupAnalysis.cfgFlagsCalculationFractionPrimary) || isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (groupAnalysis.cfgFlagQaTrack.value || groupAnalysis.cfgFlagQaDca.value || doQaAcceptance || doQaPhi || doQaPid || groupAnalysis.cfgFlagQaMc.value || doCalculationYield || doCalculationPurity || doCalculationFractionPrimary || doCalculationFluctuation) { + if (doCalculationFluctuation) { initCalculationFluctuation(); } for (const auto& track : tracks) { - if (!track.has_collision() || track.collisionId() != collision.globalIndex() || !track.has_mcParticle()) { + if (!track.has_mcParticle()) { continue; } @@ -2337,18 +2345,18 @@ struct PartNumFluc { continue; } - if (!initTrack(track, tracksIu) || !initMcParticle(mcParticle)) { + if (!initTrack(track) || !initMcParticle(mcParticle)) { continue; } - if (isEnabled(groupAnalysis.cfgFlagsQaPhi)) { + if (doQaPhi) { fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); fillQaPhiByParticleSpeciesAll(); } - if (isEnabled(groupAnalysis.cfgFlagsQaPid)) { + if (doQaPid) { fillQaPidByParticleSpeciesAll(track); fillQaPidByParticleSpeciesAll(track); fillQaPidByParticleSpeciesAll(track); @@ -2360,32 +2368,32 @@ struct PartNumFluc { hrQaMc.fill(C_CS("QaMc/hCentralityPtMcEtaMcDeltaEta"), holderEvent.centralityCalibration, holderMcParticle.pt, holderMcParticle.eta, holderTrack.eta - holderMcParticle.eta); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationYield) && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { + if (doCalculationYield && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); fillCalculationYieldByParticleSpecies(); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationPurity) && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { + if (doCalculationPurity && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { fillCalculationPurityByParticleSpecies(); fillCalculationPurityByParticleSpecies(); fillCalculationPurityByParticleSpecies(); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFractionPrimary)) { + if (doCalculationFractionPrimary) { fillCalculationFractionPrimaryByParticleSpecies(mcParticle); fillCalculationFractionPrimaryByParticleSpecies(mcParticle); fillCalculationFractionPrimaryByParticleSpecies(mcParticle); } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation) && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { + if (doCalculationFluctuation && (!groupTrack.cfgFlagMcParticlePhysicalPrimary.value || mcParticle.isPhysicalPrimary())) { calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); calculateFluctuationByParticleNumber(); } } - if (isEnabled(groupAnalysis.cfgFlagsCalculationFluctuation)) { + if (doCalculationFluctuation) { fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber(); fillCalculationFluctuationByParticleNumber();