diff --git a/PWGLF/DataModel/LFResonanceTables.h b/PWGLF/DataModel/LFResonanceTables.h index 11fd45842ef..673d2ae8393 100644 --- a/PWGLF/DataModel/LFResonanceTables.h +++ b/PWGLF/DataModel/LFResonanceTables.h @@ -38,6 +38,7 @@ #include #include #include +#include namespace o2::aod { @@ -66,6 +67,7 @@ enum { }; DECLARE_SOA_INDEX_COLUMN_FULL(Collision, collision, int, Collisions, "_Col"); //! DECLARE_SOA_COLUMN(Cent, cent, float); //! Centrality (Multiplicity) percentile (Default: FT0M) +DECLARE_SOA_COLUMN(Multiplicity, multiplicity, float); //! Configurable reconstructed multiplicity estimator DECLARE_SOA_COLUMN(Spherocity, spherocity, float); //! Spherocity of the event DECLARE_SOA_COLUMN(EvtPl, evtPl, float); //! Second harmonic event plane DECLARE_SOA_COLUMN(EvtPlResAB, evtPlResAB, float); //! Second harmonic event plane resolution of A-B sub events @@ -84,16 +86,33 @@ DECLARE_SOA_COLUMN(MCMultiplicity, mcMultiplicity, float); //! MC Multiplicit } // namespace resocollision -// Keep the established ResoCollisionColls schema above unchanged. Automatic -// GroupSlicer association to aod::Collisions requires the canonical physical -// column name fIndexCollisions, so the modular initializer writes this small -// companion table for the hybrid daughter process. +// Legacy version-0 companion schema. Its canonical fIndexCollisions column +// enables GroupSlicer association to aod::Collisions, but also makes the table +// dependent on that source parent. The modular initializer therefore writes +// the scalar-only version 001 declared below instead. namespace resocollisiongroup { DECLARE_SOA_INDEX_COLUMN_FULL_CUSTOM(OriginalCollision, originalCollision, int, Collisions, "Collisions", ""); //! } // namespace resocollisiongroup -DECLARE_SOA_TABLE(ResoCollisions, "AOD", "RESOCOLLISION", +// Scalar-only original-collision mapping for standalone modular output. Unlike +// resocollisiongroup::OriginalCollisionId above, this column carries no index +// target metadata and therefore does not require the source Collisions table +// to be present while merging the derived AO2D. +namespace resocollisiongroup001 +{ +DECLARE_SOA_COLUMN(OriginalCollisionId, originalCollisionId, int); //! Original aod::Collision row number +} // namespace resocollisiongroup001 + +// Optional soft link from a positional ResoMCCollisions_001 row to its source +// generator collision. It is filled only when the source AO2D is retained as +// a linked-derived-data parent. +namespace resomccollision +{ +DECLARE_SOA_INDEX_COLUMN_FULL_CUSTOM(OriginalMcCollision, originalMcCollision, int, McCollisions, "McCollisions", "_MC"); //! +} // namespace resomccollision + +DECLARE_SOA_TABLE(ResoCollisions_000, "AOD", "RESOCOLLISION", o2::soa::Index<>, o2::aod::mult::MultNTracksPV, o2::aod::mult::MultNTracksPVeta1, @@ -104,6 +123,23 @@ DECLARE_SOA_TABLE(ResoCollisions, "AOD", "RESOCOLLISION", resocollision::Cent, resocollision::BMagField, resocollision::IsRecINELgt0); + +// Version 001 stores one configurable multiplicity estimator instead of three +// fixed PV-track multiplicities. The producer configuration determines which +// reconstructed estimator is persisted in the common float payload. +DECLARE_SOA_TABLE_VERSIONED(ResoCollisions_001, "AOD", "RESOCOLLISION", 1, + o2::soa::Index<>, + resocollision::Multiplicity, + collision::PosX, + collision::PosY, + collision::PosZ, + resocollision::Cent, + resocollision::BMagField, + resocollision::IsRecINELgt0); + +// Keep the established generic collision-table alias on version 000. The +// modular producer and its consumers request version 001 explicitly. +using ResoCollisions = ResoCollisions_000; using ResoCollision = ResoCollisions::iterator; DECLARE_SOA_TABLE(ResoCollisionColls, "AOD", "RESOCOLLISIONCOL", @@ -114,6 +150,12 @@ DECLARE_SOA_TABLE(ResoCollisionGroups, "AOD", "RESOCOLLGROUP", resocollisiongroup::OriginalCollisionId); using ResoCollisionGroup = ResoCollisionGroups::iterator; +// Version 001 replaces the hard relation to the source Collisions table with +// a scalar row number. The legacy version-0 schema remains unchanged. +DECLARE_SOA_TABLE_VERSIONED(ResoCollisionGroups_001, "AOD", "RESOCOLLGROUP", 1, + resocollisiongroup001::OriginalCollisionId); +using ResoCollisionGroup_001 = ResoCollisionGroups_001::iterator; + DECLARE_SOA_TABLE(ResoMCCollisions, "AOD", "RESOMCCOLLISION", o2::soa::Index<>, resocollision::IsVtxIn10, @@ -125,6 +167,20 @@ DECLARE_SOA_TABLE(ResoMCCollisions, "AOD", "RESOMCCOLLISION", resocollision::MCMultiplicity); using ResoMCCollision = ResoMCCollisions::iterator; +// Version 001 deliberately persists generator-collision properties only. +// Reconstructed event-selection decisions remain available only in legacy +// version-0 output and must not be mixed into this generator-level payload. +DECLARE_SOA_TABLE_VERSIONED(ResoMCCollisions_001, "AOD", "RESOMCCOLLISION", 1, + o2::soa::Index<>, + resocollision::IsVtxIn10, + resocollision::IsINELgt0, + resocollision::ImpactParameter, + resocollision::MCMultiplicity); + +DECLARE_SOA_TABLE(ResoMCCollisionIds, "AOD", "RESOMCCOLLID", + resomccollision::OriginalMcCollisionId); +using ResoMCCollisionId = ResoMCCollisionIds::iterator; + DECLARE_SOA_TABLE(ResoSpheroCollisions, "AOD", "RESOSPHEROCOLLISION", o2::soa::Index<>, resocollision::Spherocity); @@ -204,6 +260,12 @@ struct ResoTrackFlags { #define DECLARE_DYN_TRKSEL_COLUMN(_Name_, _Getter_, _Mask_) \ DECLARE_SOA_DYNAMIC_COLUMN(_Name_, _Getter_, [](ResoTrackFlags::flagtype flags) -> bool { return ResoTrackFlags::checkFlag(flags, _Mask_); }); +// Keep the default foreign-key target on the legacy ResoCollisions_000 alias. +// Modular v001 consumers use resoCollisionId() as the scalar row reference and +// explicitly call resoCollision_as() only when the parent row must be +// dereferenced, with T equal to the exact bound v001 table type (plain +// aod::ResoCollisions_001 or its analysis Join). The version-equivalence +// declaration below permits binding the same physical index column to v001. DECLARE_SOA_INDEX_COLUMN(ResoCollision, resoCollision); DECLARE_SOA_INDEX_COLUMN(ResoCollisionDF, resoCollisionDF); DECLARE_SOA_INDEX_COLUMN_FULL(Track, track, int, Tracks, "_Trk"); //! Soft link to the original track @@ -246,7 +308,7 @@ DECLARE_SOA_COLUMN(DecayVtxY, decayVtxY, float); DECLARE_SOA_COLUMN(DecayVtxZ, decayVtxZ, float); //! Z position of the decay vertex DECLARE_SOA_COLUMN(Alpha, alpha, float); //! Alpha of the decay vertex DECLARE_SOA_COLUMN(QtArm, qtarm, float); //! Armenteros Qt of the decay vertex, o2-linter: disable=name/o2-column (pre-existing public column name kept for schema and API compatibility) -DECLARE_SOA_COLUMN(TpcSignal10, tpcSignal10, int16_t); //! TPC signal of the track x10 +DECLARE_SOA_COLUMN(TpcSignal10, tpcSignal10, int16_t); //! TPC signal of the track x100 (public column name retained for compatibility) DECLARE_SOA_COLUMN(DaughterTPCNSigmaPosPi10, daughterTPCNSigmaPosPi10, int8_t); //! TPC PID x10 of the positive daughter as Pion DECLARE_SOA_COLUMN(DaughterTPCNSigmaPosKa10, daughterTPCNSigmaPosKa10, int8_t); //! TPC PID x10 of the positive daughter as Kaon DECLARE_SOA_COLUMN(DaughterTPCNSigmaPosPr10, daughterTPCNSigmaPosPr10, int8_t); //! TPC PID x10 of the positive daughter as Proton @@ -340,7 +402,7 @@ DECLARE_SOA_DYNAMIC_COLUMN(DaughterTOFNSigmaBachKa, daughterTOFNSigmaBachKa, [](int8_t daughterTOFNSigmaBachKa10) { return (float)daughterTOFNSigmaBachKa10 / 10.f; }); DECLARE_SOA_DYNAMIC_COLUMN(DaughterTOFNSigmaBachPr, daughterTOFNSigmaBachPr, [](int8_t daughterTOFNSigmaBachPr10) { return (float)daughterTOFNSigmaBachPr10 / 10.f; }); -// TPC signal x10 +// TPC signal x100 DECLARE_SOA_DYNAMIC_COLUMN(TpcSignal, tpcSignal, [](int16_t tpcSignal10) { return (float)tpcSignal10 / 100.f; }); // pT, Eta, Phi @@ -554,6 +616,192 @@ struct ResoMicroTrackSelFlag { DECLARE_SOA_DYNAMIC_COLUMN(Pt, pt, [](float px, float py) -> float { return RecoDecay::sqrtSumOfSquares(px, py); }); } // namespace resomicrodaughter +// Version 001 uses signed, lower-inclusive PID bins and keeps the DCAxy/DCAz +// pT-dependent selection results independently from the quantised DCA values. +namespace resomicrodaughter001 +{ +// Keep the original row number as a scalar for standalone pair comparisons. +// The separate typed relation remains available to legacy version-0 users. +DECLARE_SOA_COLUMN(TrackId, trackId, int); //! Original track row number for pair-level comparisons + +/// @brief Compact signed TPC/TOF n-sigma values into two four-bit fields. +/// The upper bit of each field stores the sign. Magnitudes 1..6 represent +/// lower-inclusive 0.25-sigma-wide bins [2.0, 2.25), ..., [3.25, 3.5), while +/// magnitude 0 is |n-sigma| < 2 and magnitude 7 is |n-sigma| >= 3.5. Overflow +/// decodes to signed infinity, while code 8 is reserved for invalid values and +/// decodes to NaN. Missing TOF information is also carried independently by +/// resodaughter::TrackFlags::kHasTOF. +struct PidNSigma { + static constexpr float MinFineNSigma = 2.f; + static constexpr float Step = 0.25f; + static constexpr float MaxNSigma = 3.5f; + static constexpr uint8_t SignMask = 0x08; + static constexpr uint8_t MagnitudeMask = 0x07; + static constexpr uint8_t MaxRegularCode = 6; + static constexpr uint8_t AboveRangeCode = 7; + static constexpr uint8_t InvalidCode = SignMask; + + uint8_t flag; + + PidNSigma(float tpcNSigma, float tofNSigma, bool hasTOF) + { + const uint8_t tpcEncoded = encodeNSigma(tpcNSigma); + const uint8_t tofEncoded = hasTOF ? encodeNSigma(tofNSigma) : InvalidCode; + flag = (tpcEncoded << 4) | tofEncoded; + } + + static uint8_t encodeNSigma(float nSigma) + { + if (!std::isfinite(nSigma)) { + return InvalidCode; + } + const float value = std::abs(nSigma); + if (value < MinFineNSigma) { + return 0; + } + uint8_t magnitude = AboveRangeCode; + if (value < MaxNSigma) { + const int encoded = 1 + static_cast(std::floor((value - MinFineNSigma) / Step)); + magnitude = static_cast(std::clamp(encoded, 1, static_cast(MaxRegularCode))); + } + const uint8_t sign = std::signbit(nSigma) ? SignMask : 0; + return static_cast(sign | magnitude); + } + + static float decodeNSigma(uint8_t encoded) + { + const uint8_t code = encoded & 0x0F; + if (code == InvalidCode) { + return NAN; + } + const uint8_t magnitude = code & MagnitudeMask; + if (magnitude == 0) { + return 0.f; + } + const float value = magnitude == AboveRangeCode + ? std::numeric_limits::infinity() + : MinFineNSigma + static_cast(magnitude - 1) * Step; + return (code & SignMask) != 0 ? -value : value; + } + + static float getTPCNSigma(uint8_t encoded) + { + return decodeNSigma((encoded >> 4) & 0x0F); + } + + static float getTOFNSigma(uint8_t encoded, bool hasTOF) + { + return hasTOF ? decodeNSigma(encoded & 0x0F) : NAN; + } + + operator uint8_t() const { return flag; } +}; + +/// @brief Store two pT-dependent DCA pass bits and two three-bit DCA values. +/// Bits 7 and 6 store the DCAxy and DCAz pass results, respectively. Bits 5..3 +/// and 2..0 store |DCAxy| and |DCAz| in lower-inclusive 0.025 cm bins from 0 +/// to 0.15 cm. Code 6 represents overflow and code 7 an invalid value. +struct DCAEncoding { + static constexpr float MaxDCA = 0.15f; + static constexpr float Step = 0.025f; + static constexpr float InverseStep = 40.f; + static constexpr uint8_t MaxRegularCode = 5; + static constexpr uint8_t AboveRangeCode = 6; + static constexpr uint8_t InvalidCode = 7; + static constexpr uint8_t PassedPtDependentDCAxyMask = 0x80; + static constexpr uint8_t PassedPtDependentDCAzMask = 0x40; + static constexpr uint8_t DCAxyMask = 0x38; + static constexpr uint8_t DCAzMask = 0x07; + static constexpr uint8_t DCAxyShift = 3; + + uint8_t flag = 0; + + DCAEncoding() = default; + DCAEncoding(float dcaXY, float dcaZ, bool passedPtDependentDCAxy, bool passedPtDependentDCAz) + : flag(static_cast((passedPtDependentDCAxy ? PassedPtDependentDCAxyMask : 0) | + (passedPtDependentDCAz ? PassedPtDependentDCAzMask : 0) | + (encodeDCA(dcaXY) << DCAxyShift) | + encodeDCA(dcaZ))) + { + } + + static uint8_t encodeDCA(float dca) + { + const float value = std::abs(dca); + if (!std::isfinite(value)) { + return InvalidCode; + } + if (value >= MaxDCA) { + return AboveRangeCode; + } + const int encoded = static_cast(std::floor(value * InverseStep)); + return static_cast(std::clamp(encoded, 0, static_cast(MaxRegularCode))); + } + + static float decodeDCA(uint8_t encoded) + { + const uint8_t code = encoded & DCAzMask; + if (code == InvalidCode) { + return NAN; + } + return code == AboveRangeCode ? MaxDCA : static_cast(code) * Step; + } + + static float decodeDCAxy(uint8_t encoded) + { + return decodeDCA((encoded & DCAxyMask) >> DCAxyShift); + } + + static float decodeDCAz(uint8_t encoded) + { + return decodeDCA(encoded & DCAzMask); + } + + static bool testPtDependentDCAxy(uint8_t encoded) + { + return (encoded & PassedPtDependentDCAxyMask) != 0; + } + + static bool testPtDependentDCAz(uint8_t encoded) + { + return (encoded & PassedPtDependentDCAzMask) != 0; + } + + operator uint8_t() const { return flag; } +}; + +DECLARE_SOA_COLUMN(TrackSelectionFlags, trackSelectionFlags, uint8_t); //! Packed DCA selection and absolute DCAxy/DCAz values +DECLARE_SOA_DYNAMIC_COLUMN(TpcNSigmaPi, tpcNSigmaPi, + [](uint8_t pidNSigmaPiFlag) { return PidNSigma::getTPCNSigma(pidNSigmaPiFlag); }); +DECLARE_SOA_DYNAMIC_COLUMN(TpcNSigmaKa, tpcNSigmaKa, + [](uint8_t pidNSigmaKaFlag) { return PidNSigma::getTPCNSigma(pidNSigmaKaFlag); }); +DECLARE_SOA_DYNAMIC_COLUMN(TpcNSigmaPr, tpcNSigmaPr, + [](uint8_t pidNSigmaPrFlag) { return PidNSigma::getTPCNSigma(pidNSigmaPrFlag); }); +DECLARE_SOA_DYNAMIC_COLUMN(TofNSigmaPi, tofNSigmaPi, + [](uint8_t pidNSigmaPiFlag, uint8_t trackFlags) -> float { + const bool hasTOF = resodaughter::ResoTrackFlags::checkFlag(trackFlags, resodaughter::ResoTrackFlags::kHasTOF); + return PidNSigma::getTOFNSigma(pidNSigmaPiFlag, hasTOF); + }); +DECLARE_SOA_DYNAMIC_COLUMN(TofNSigmaKa, tofNSigmaKa, + [](uint8_t pidNSigmaKaFlag, uint8_t trackFlags) -> float { + const bool hasTOF = resodaughter::ResoTrackFlags::checkFlag(trackFlags, resodaughter::ResoTrackFlags::kHasTOF); + return PidNSigma::getTOFNSigma(pidNSigmaKaFlag, hasTOF); + }); +DECLARE_SOA_DYNAMIC_COLUMN(TofNSigmaPr, tofNSigmaPr, + [](uint8_t pidNSigmaPrFlag, uint8_t trackFlags) -> float { + const bool hasTOF = resodaughter::ResoTrackFlags::checkFlag(trackFlags, resodaughter::ResoTrackFlags::kHasTOF); + return PidNSigma::getTOFNSigma(pidNSigmaPrFlag, hasTOF); + }); +DECLARE_SOA_DYNAMIC_COLUMN(DcaXY, dcaXY, + [](uint8_t trackSelectionFlags) { return DCAEncoding::decodeDCAxy(trackSelectionFlags); }); +DECLARE_SOA_DYNAMIC_COLUMN(DcaZ, dcaZ, + [](uint8_t trackSelectionFlags) { return DCAEncoding::decodeDCAz(trackSelectionFlags); }); +DECLARE_SOA_DYNAMIC_COLUMN(PassedPtDependentDCAxy, passedPtDependentDCAxy, + [](uint8_t trackSelectionFlags) { return DCAEncoding::testPtDependentDCAxy(trackSelectionFlags); }); +DECLARE_SOA_DYNAMIC_COLUMN(PassedPtDependentDCAz, passedPtDependentDCAz, + [](uint8_t trackSelectionFlags) { return DCAEncoding::testPtDependentDCAz(trackSelectionFlags); }); +} // namespace resomicrodaughter001 + // Ultra-micro track representation. The momentum components are quantised // to 1 MeV/c and only one (pion/kaon/proton) PID flag is retained. namespace resoultramicrodaughter @@ -786,6 +1034,45 @@ DECLARE_SOA_TABLE(ResoMicroTracks, "AOD", "RESOMICROTRACK", resodaughter::Sign); using ResoMicroTrack = ResoMicroTracks::iterator; +// Keep ResoMicroTracks as the version-0 API for existing producers and +// consumers. Version-1 users must request ResoMicroTracks_001 explicitly. +DECLARE_SOA_TABLE_VERSIONED(ResoMicroTracks_001, "AOD", "RESOMICROTRACK", 1, + o2::soa::Index<>, + resodaughter::ResoCollisionId, + resomicrodaughter001::TrackId, + resodaughter::Px, + resodaughter::Py, + resodaughter::Pz, + resomicrodaughter::PidNSigmaPiFlag, + resomicrodaughter::PidNSigmaKaFlag, + resomicrodaughter::PidNSigmaPrFlag, + resomicrodaughter001::TrackSelectionFlags, + resodaughter::TrackFlags, + // Dynamic columns + resomicrodaughter::Pt, + resodaughter::Eta, + resodaughter::Phi, + resomicrodaughter001::TpcNSigmaPi, + resomicrodaughter001::TpcNSigmaKa, + resomicrodaughter001::TpcNSigmaPr, + resomicrodaughter001::TofNSigmaPi, + resomicrodaughter001::TofNSigmaKa, + resomicrodaughter001::TofNSigmaPr, + resomicrodaughter001::DcaXY, + resomicrodaughter001::DcaZ, + resomicrodaughter001::PassedPtDependentDCAxy, + resomicrodaughter001::PassedPtDependentDCAz, + resodaughter::PassedITSRefit, + resodaughter::PassedTPCRefit, + resodaughter::IsGlobalTrackWoDCA, + resodaughter::IsGlobalTrack, + resodaughter::IsPrimaryTrack, + resodaughter::IsPVContributor, + resodaughter::HasTOF, + resodaughter::Sign); +// Positional soft-link side table retained for ResoMicroTracks version 000. +// Version 001 stores the same row number as a scalar and should be consumed +// without joining this side table, since both columns expose trackId(). DECLARE_SOA_TABLE(ResoMicroTrackTracks, "AOD", "RESOMICROTRACKTRACK", resodaughter::TrackId); using ResoMicroTrackTrack = ResoMicroTrackTracks::iterator; @@ -1074,6 +1361,17 @@ DECLARE_SOA_TABLE(ResoMCTracks, "AOD", "RESOMCTRACK", resodaughter::ProducedByGenerator); using ResoMCTrack = ResoMCTracks::iterator; +// Positional MC extension for ResoMicroTracks_001. One row must be written for +// every ResoMicroTracks_001 row, including tracks without an MC particle label. +DECLARE_SOA_TABLE_VERSIONED(ResoMCMicroTracks_001, "AOD", "RESOMCMICROTRACK", 1, + mcparticle::PdgCode, + resodaughter::MotherId, + resodaughter::MotherPDG, + resodaughter::SiblingIds, + resodaughter::IsPhysicalPrimary, + resodaughter::ProducedByGenerator); +using ResoMCMicroTrack = ResoMCMicroTracks_001::iterator; + DECLARE_SOA_TABLE(ResoMCV0s, "AOD", "RESOMCV0", mcparticle::PdgCode, resodaughter::MotherId, @@ -1123,6 +1421,34 @@ DECLARE_SOA_TABLE(ResoMCParents, "AOD", "RESOMCPARENT", resodaughter::Phi); using ResoMCParent = ResoMCParents::iterator; +// Module-specific parent rows keep the source MC-particle row number as a +// scalar. This avoids a hard relation to McParticles in standalone derived +// AO2D while preserving the legacy ResoMCParents schema and typed API. +namespace resomcparent001 +{ +DECLARE_SOA_COLUMN(OriginalMcParticleId, originalMcParticleId, int); //! Original aod::McParticle row number +} // namespace resomcparent001 + +DECLARE_SOA_TABLE_VERSIONED(ResoMCParents_001, "AOD", "RESOMCPARENT", 1, + o2::soa::Index<>, + resodaughter::ResoCollisionId, + resomcparent001::OriginalMcParticleId, + mcparticle::PdgCode, + resodaughter::DaughterPDG1, + resodaughter::DaughterPDG2, + resodaughter::IsPhysicalPrimary, + resodaughter::ProducedByGenerator, + resodaughter::Pt, + resodaughter::Px, + resodaughter::Py, + resodaughter::Pz, + mcparticle::Y, + mcparticle::E, + mcparticle::StatusCode, + resodaughter::Eta, + resodaughter::Phi); +using ResoMCParent_001 = ResoMCParents_001::iterator; + using Reso2TracksExt = soa::Join; // without Extra using Reso2TracksMC = soa::Join; using Reso2TracksPID = soa::Join; @@ -1141,4 +1467,11 @@ using ResoCascadesCandidatesMC = soa::Join; } // namespace o2::aod + +namespace o2::soa +{ +// Preserve the legacy v000 default target while allowing an explicitly typed +// resoCollision_as() binding to the exact v001 parent in modular analyses. +DECLARE_EQUIVALENT_FOR_INDEX(aod::ResoCollisions_000, aod::ResoCollisions_001); +} // namespace o2::soa #endif // PWGLF_DATAMODEL_LFRESONANCETABLES_H_ diff --git a/PWGLF/TableProducer/Resonances/resonanceModuleInitializer.cxx b/PWGLF/TableProducer/Resonances/resonanceModuleInitializer.cxx index 1d748682f6c..21e534db7fc 100644 --- a/PWGLF/TableProducer/Resonances/resonanceModuleInitializer.cxx +++ b/PWGLF/TableProducer/Resonances/resonanceModuleInitializer.cxx @@ -13,15 +13,15 @@ /// \brief Initializes variables for the resonance candidate producers /// /// \author Bong-Hwi Lim , Minjae Kim -/// \since Aug.18 2026 +/// \since Aug.31 2026 #include "PWGLF/DataModel/LFResonanceTables.h" #include "PWGLF/DataModel/LFStrangenessTables.h" #include "PWGLF/DataModel/mcCentrality.h" #include "PWGLF/Utils/collisionCuts.h" -#include "Common/CCDB/EventSelectionParams.h" #include "Common/CCDB/RCTSelectionFlags.h" +#include "Common/DataModel/EventSelection.h" #include "Common/DataModel/Multiplicity.h" #include "Common/DataModel/PIDResponseTOF.h" #include "Common/DataModel/PIDResponseTPC.h" @@ -40,10 +40,14 @@ #include #include #include -#include #include #include +#include +#include +#include + +#include #include #include #include @@ -73,23 +77,60 @@ using namespace o2::aod::rctsel; struct ResonanceModuleInitializer { static constexpr double BzOverrideThreshold = -990.; static constexpr double MinimumNonzeroBz = 1.e-5; - static constexpr double MinimumChargedParticleCharge = 3.; // ROOT particle charge is stored in units of e/3 - static constexpr int MCCentralityRecoEstimator = 0; - static constexpr int MCCentralityGeneratorEstimator = 1; - static constexpr int MCCentralityImpactParameterEstimator = 2; + static constexpr int CentralityFT0M = 0; + static constexpr int CentralityFT0C = 1; + static constexpr int CentralityFT0A = 2; + static constexpr int CentralityFV0A = 3; + static constexpr int MultiplicityNTracksPV = 0; + static constexpr int MultiplicityNTracksPVeta1 = 1; + static constexpr int MultiplicityNTracksPVetaHalf = 2; + static constexpr int MultiplicityFT0M = 3; + static constexpr int MultiplicityFT0A = 4; + static constexpr int MultiplicityFT0C = 5; + static constexpr int MultiplicityFV0A = 6; + static constexpr int DetailedQARCTStage = o2::analysis::CollisonCuts::kAllpassed + 1; + static constexpr int DetailedQAStages = DetailedQARCTStage + 1; static constexpr float MCVertexZMax = 10.f; - - int mRunNumber = 0; ///< Run number for the current data - int multEstimator = 0; ///< Multiplicity estimator type - float dBz = 0.f; ///< Magnetic field value - float centrality = 0.f; ///< Centrality value for the event - Service ccdb; ///< CCDB manager service - Service pdg; ///< PDG database service - - Produces resoCollisions; ///< Output table for resonance collisions - Produces resoCollisionColls; ///< Output table for collision references - Produces resoCollisionGroups; ///< Canonical original-collision grouping references - Produces resoMCCollisions; ///< Output table for MC resonance collisions + // PDG codes used by the persistent resonance-parent selection. Named O2/ROOT + // values are preferred where available; the remaining resonances are kept as + // local named constants because neither PDG_t nor PhysicsConstants defines them. + static constexpr int PdgKStar0 = o2::constants::physics::Pdg::kK0Star892; + static constexpr int PdgKStarCharged = o2::constants::physics::Pdg::kKPlusStar892; + static constexpr int PdgPhi = o2::constants::physics::Pdg::kPhi; + static constexpr int F0Code980 = 9010221; + static constexpr int F0Code1370 = 10221; + static constexpr int F0Code1500 = 9030221; + static constexpr int F0Code1710 = 10331; + static constexpr int F1Code1285 = 20223; + static constexpr int F1Code1420 = 20333; + static constexpr int F2PrimeCode1525 = 335; + static constexpr int PdgRho0 = PDG_t::kRho770_0; + static constexpr int PdgRhoCharged = PDG_t::kRho770Plus; + static constexpr int SigmaStarPlusCode = 3224; + static constexpr int PdgLambda1520 = o2::constants::physics::Pdg::kLambda1520_Py; + static constexpr int Xi1530Code = 3324; + static constexpr int PdgK1Plus1270 = o2::constants::physics::Pdg::kK1_1270Plus; + static constexpr int Xi1820NeutralCode = 123314; + static constexpr int Xi1820MinusCode = 123324; + static constexpr int Omega2012MinusCode = 123334; + static constexpr int PdgProton = PDG_t::kProton; + static constexpr int PdgLambda0 = PDG_t::kLambda0; + static constexpr int PdgXiMinus = PDG_t::kXiMinus; + static constexpr int PdgXi0 = o2::constants::physics::Pdg::kXi0; + static constexpr int PdgOmegaMinus = PDG_t::kOmegaMinus; + + int mRunNumber = 0; ///< Run number for the current data + int multEstimator = CentralityFT0M; ///< Centrality estimator type + float dBz = 0.f; ///< Magnetic field value + float centrality = 0.f; ///< Centrality value for the event + Service ccdb; ///< CCDB manager service + + Produces resoCollisions; ///< Output table for resonance collisions + Produces resoCollisionColls; ///< Optional source collision soft links + Produces resoCollisionGroups001; ///< Scalar original-collision grouping keys + Produces resoMCCollisions001; ///< Generator-only MC collision extension + Produces resoMCCollisionIds; ///< Optional source generator-collision soft links + Produces reso2mcparents; ///< Generated parents with scalar source-particle IDs // CCDB options struct : ConfigurableGroup { @@ -99,18 +140,21 @@ struct ResonanceModuleInitializer { Configurable lutPath{"lutPath", "GLO/Param/MatLUT", "Path of the Lut parametrization"}; Configurable geoPath{"geoPath", "GLO/Config/GeometryAligned", "Path of the geometry file"}; Configurable cfgFatalWhenNull{"cfgFatalWhenNull", true, "Fatal when null on ccdb access"}; - Configurable cfgBypassCollIndexFill{"cfgBypassCollIndexFill", false, "Unsupported in the modular workflow; must remain false"}; + Configurable cfgBypassCollIndexFill{"cfgBypassCollIndexFill", false, "Deprecated compatibility option; collision mapping tables are always written"}; } CCDB; // General event options struct : ConfigurableGroup { Configurable dBzInput{"dBzInput", -999, "bz field, -999 is automatic"}; Configurable cfgFillQA{"cfgFillQA", true, "Fill QA histograms"}; + Configurable cfgFillDetailedQA{"cfgFillDetailedQA", true, "Fill the Run 3 event-selection stage vs vertex-z vs centrality vs multiplicity THnSparse"}; Configurable cfgBypassCCDB{"cfgBypassCCDB", true, "Bypass loading CCDB part to save CPU time and memory"}; // will be affected to b_z value. - Configurable cfgMultName{"cfgMultName", "FT0M", "The name of multiplicity estimator"}; - Configurable cfgCentralityMC{"cfgCentralityMC", 0, "Centrality estimator for MC (0: Reco, 1: MC, 2: impact parameter)"}; - ConfigurableAxis binsCent{"binsCent", {VARIABLE_WIDTH, 0., 0.01, 0.1, 1.0, 5.0, 10., 15., 20., 30., 40., 50., 70., 100.0, 105.}, "Binning of the centrality axis"}; - ConfigurableAxis cfgVtxBins{"cfgVtxBins", {VARIABLE_WIDTH, -20, -15, -10, -7, -5, -3, -2, -1, 0, 1, 2, 3, 5, 7, 10, 15, 20}, "Mixing bins - z-vertex"}; + Configurable cfgMultName{"cfgMultName", "FT0M", "Centrality estimator: FT0M, FT0C, FT0A, or FV0A"}; + Configurable cfgMultiplicityEstimator{"cfgMultiplicityEstimator", 3, + "Stored multiplicity (NOT percentile): 0 -> NTracksPV, 1 -> NTracksPVeta1, 2 -> NTracksPVetaHalf, 3 -> FT0M, 4 -> FT0A, 5 -> FT0C, 6 -> FV0A"}; + ConfigurableAxis binsCent{"binsCent", {VARIABLE_WIDTH, 0., 0.01, 0.1, 1., 5., 10., 15., 20., 30., 40., 50., 60., 70., 80., 90., 100., 105.}, "Binning of the centrality axis"}; + ConfigurableAxis binsMultiplicity{"binsMultiplicity", {500, 0.f, 5000.f}, "Binning of the reconstructed multiplicity axis for detailed collision QA"}; + ConfigurableAxis cfgVtxBins{"cfgVtxBins", {400, -20.f, 20.f}, "Binning of the collision vertex-z axis for detailed QA"}; } EventConfig; /// Event cuts @@ -119,26 +163,86 @@ struct ResonanceModuleInitializer { Configurable cfgEvtZvtx{"cfgEvtZvtx", 10.f, "Evt sel: Max. z-Vertex (cm)"}; Configurable cfgEvtOccupancyInTimeRange{"cfgEvtOccupancyInTimeRange", -1, "Evt sel: maximum track occupancy"}; Configurable cfgEvtTriggerCheck{"cfgEvtTriggerCheck", false, "Evt sel: check for trigger"}; - Configurable cfgEvtOfflineCheck{"cfgEvtOfflineCheck", true, "Evt sel: check for offline selection"}; - Configurable cfgEvtTriggerTVXSel{"cfgEvtTriggerTVXSel", false, "Evt sel: triggerTVX selection (MB)"}; - Configurable cfgEvtTFBorderCut{"cfgEvtTFBorderCut", false, "Evt sel: apply TF border cut"}; + Configurable cfgEvtOfflineCheck{"cfgEvtOfflineCheck", false, "Evt sel: check for offline selection (sel8)"}; + Configurable cfgEvtTriggerTVXSel{"cfgEvtTriggerTVXSel", true, "Evt sel: triggerTVX selection (MB)"}; + Configurable cfgEvtTFBorderCut{"cfgEvtTFBorderCut", true, "Evt sel: apply TF border cut"}; Configurable cfgEvtUseITSTPCvertex{"cfgEvtUseITSTPCvertex", false, "Evt sel: use at lease on ITS-TPC track for vertexing"}; Configurable cfgEvtCollInTimeRangeNarrow{"cfgEvtCollInTimeRangeNarrow", false, "Evt sel: apply NoCollInTimeRangeNarrow"}; Configurable cfgEvtZvertexTimedifference{"cfgEvtZvertexTimedifference", false, "Evt sel: apply Z-vertex time difference"}; Configurable cfgEvtPileupRejection{"cfgEvtPileupRejection", false, "Evt sel: apply pileup rejection"}; Configurable cfgEvtNoITSROBorderCut{"cfgEvtNoITSROBorderCut", false, "Evt sel: apply NoITSRO border cut"}; - Configurable cfgEvtRun2AliEventCuts{"cfgEvtRun2AliEventCuts", true, "Evt sel: apply Run2 AliEventCuts"}; + Configurable cfgEvtRun2AliEventCuts{"cfgEvtRun2AliEventCuts", false, "Evt sel: apply Run2 AliEventCuts"}; Configurable cfgEvtRun2INELgtZERO{"cfgEvtRun2INELgtZERO", false, "Evt sel: apply Run2 INELgtZERO"}; - Configurable cfgEvtUseRCTFlagChecker{"cfgEvtUseRCTFlagChecker", false, "Evt sel: use RCT flag checker"}; + Configurable cfgEvtUseRCTFlagChecker{"cfgEvtUseRCTFlagChecker", true, "Evt sel: use RCT flag checker"}; + Configurable cfgEvtBCRCT{"cfgEvtBCRCT", false, "Evt sel: check RCT on the nominal associated BCSEL instead of the collision EVSEL"}; Configurable cfgEvtRCTFlagCheckerLabel{"cfgEvtRCTFlagCheckerLabel", "CBT_hadronPID", "Evt sel: RCT flag checker label"}; Configurable cfgEvtRCTFlagCheckerZDCCheck{"cfgEvtRCTFlagCheckerZDCCheck", false, "Evt sel: RCT flag checker ZDC check"}; Configurable cfgEvtRCTFlagCheckerLimitAcceptAsBad{"cfgEvtRCTFlagCheckerLimitAcceptAsBad", false, "Evt sel: RCT flag checker treat Limited Acceptance As Bad"}; + Configurable cfgEvtRCTCheckTableValidity{"cfgEvtRCTCheckTableValidity", false, "Evt sel: reject collisions when the RCT CCDB payload is unavailable"}; } EventCuts; - RCTFlagsChecker rctChecker; + RCTFlagsChecker recoRCTChecker; - HistogramRegistry qaRegistry{"QAHistos", {}, OutputObjHandlingPolicy::AnalysisObject}; + // Generator-level event and resonance QA + struct : ConfigurableGroup { + Configurable cfgGenBCRCT{"cfgGenBCRCT", false, "GenEvent: apply the RCT flag checker to the associated BC"}; + Configurable cfgGenRCTCheckTableValidity{"cfgGenRCTCheckTableValidity", false, "GenEvent: reject MC collisions when the RCT CCDB payload is unavailable"}; + Configurable cfgGenMult05{"cfgGenMult05", true, "GenEvent: multiplicity in |eta| < 0.5"}; + Configurable cfgGenMult10{"cfgGenMult10", false, "GenEvent: multiplicity in |eta| < 1.0"}; + Configurable cfgGenMultFT0M{"cfgGenMultFT0M", false, "GenEvent: generated charged-particle multiplicity in the FT0A + FT0C acceptance"}; + Configurable cfgGenMultFT0C{"cfgGenMultFT0C", false, "GenEvent: generated charged-particle multiplicity in the FT0C acceptance"}; + Configurable cfgGenMultFV0A{"cfgGenMultFV0A", false, "GenEvent: generated charged-particle multiplicity in the FV0A acceptance"}; + Configurable cfgGenMultPercentile{"cfgGenMultPercentile", true, "Use the configured FT0M, FT0C, or FV0A percentile from the MC centrality wagon"}; + Configurable cfgFillMCCollisionSoftLink{"cfgFillMCCollisionSoftLink", false, + "Write the source MC-collision soft link; enable only with the original AO2D as linked parent"}; + Configurable isZvtxcutGen{"isZvtxcutGen", true, "Apply the generator-collision z-vertex cut"}; + Configurable cutzvertexGen{"cutzvertexGen", 10.f, "Maximum absolute generator-collision z vertex (cm)"}; + Configurable checkIsTrueINELgt0{"checkIsTrueINELgt0", true, "Classify true INEL>0 generator collisions"}; + ConfigurableAxis binsCentGen{"binsCentGen", + {VARIABLE_WIDTH, 0., 0.01, 0.1, 1., 5., 10., 15., 20., 30., 40., 50., 60., 70., 80., 90., 100., 105.}, + "Generator centrality axis"}; + ConfigurableAxis ptAxisGen{"ptAxisGen", {400, 0.f, 20.f}, "#it{p}_{T} (GeV/#it{c})"}; + ConfigurableAxis multNTracksAxis{"multNTracksAxis", {500, 0.f, 5000.f}, "Number of charged particles"}; + ConfigurableAxis impactParameterAxis{"impactParameterAxis", {500, 0.f, 50.f}, "Impact parameter (fm)"}; + Configurable isDaughterCheck{"isDaughterCheck", true, "Require the configured two-body decay"}; + Configurable cfgRapidityCutMinGen{"cfgRapidityCutMinGen", -0.5f, "Minimum generated-particle rapidity"}; + Configurable cfgRapidityCutMaxGen{"cfgRapidityCutMaxGen", 0.5f, "Maximum generated-particle rapidity"}; + Configurable pdgTruthMother{"pdgTruthMother", static_cast(Xi1530Code), "Absolute PDG code of the generated mother"}; + Configurable pdgTruthDaughter1{"pdgTruthDaughter1", static_cast(PdgXiMinus), "Absolute PDG code of the first daughter"}; + Configurable pdgTruthDaughter2{"pdgTruthDaughter2", PDG_t::kPiPlus, "Absolute PDG code of the second daughter"}; + Configurable cfgDoSignalLoss{"cfgDoSignalLoss", false, "Save reference particles for mT-scaling signal-loss studies"}; + } GenCuts; + RCTFlagsChecker genRCTChecker; + + // Keep the established ResoMCParents content compatible with the legacy + // initializer. The additional stable-particle species are written only for + // signal-loss studies and are filtered in fillMCParents. + Partition selectedMCParticles = (nabs(aod::mcparticle::pdgCode) == PdgKStar0) // K*(892)0 + || (nabs(aod::mcparticle::pdgCode) == PdgKStarCharged) // K*(892)+ + || (nabs(aod::mcparticle::pdgCode) == PdgPhi) // phi(1020) + || (nabs(aod::mcparticle::pdgCode) == F0Code980) // f0(980) + || (nabs(aod::mcparticle::pdgCode) == F0Code1370) // f0(1370) + || (nabs(aod::mcparticle::pdgCode) == F0Code1500) // f0(1500) + || (nabs(aod::mcparticle::pdgCode) == F0Code1710) // f0(1710) + || (nabs(aod::mcparticle::pdgCode) == F1Code1285) // f1(1285) + || (nabs(aod::mcparticle::pdgCode) == F1Code1420) // f1(1420) + || (nabs(aod::mcparticle::pdgCode) == F2PrimeCode1525) // f2'(1525) + || (nabs(aod::mcparticle::pdgCode) == PdgRho0) // rho(770)0 + || (nabs(aod::mcparticle::pdgCode) == PdgRhoCharged) // rho(770)+ + || (nabs(aod::mcparticle::pdgCode) == SigmaStarPlusCode) // Sigma(1385)+ + || (nabs(aod::mcparticle::pdgCode) == PdgLambda1520) // Lambda(1520) + || (nabs(aod::mcparticle::pdgCode) == Xi1530Code) // Xi(1530)0 + || (nabs(aod::mcparticle::pdgCode) == PdgK1Plus1270) // K1(1270)+ + || (nabs(aod::mcparticle::pdgCode) == Xi1820NeutralCode) // Xi(1820)0 + || (nabs(aod::mcparticle::pdgCode) == Xi1820MinusCode) // Xi(1820)- + || (nabs(aod::mcparticle::pdgCode) == Omega2012MinusCode) // Omega(2012)- + || (nabs(aod::mcparticle::pdgCode) == PdgProton) // proton + || (nabs(aod::mcparticle::pdgCode) == PdgLambda0) // Lambda0 + || (nabs(aod::mcparticle::pdgCode) == PdgXiMinus) // Xi- + || (nabs(aod::mcparticle::pdgCode) == PdgXi0) // Xi0 + || (nabs(aod::mcparticle::pdgCode) == PdgOmegaMinus); // Omega- + Preslice mcParticlesPerMcCollision = aod::mcparticle::mcCollisionId; - Filter collisionFilter = nabs(aod::collision::posZ) < EventCuts.cfgEvtZvtx; + HistogramRegistry qaRegistry{"QAHistos", {}, OutputObjHandlingPolicy::AnalysisObject}; /** * @brief Initializes the task @@ -150,32 +254,65 @@ struct ResonanceModuleInitializer { mRunNumber = 0; dBz = 0; centrality = 0; - // Determine the multiplicity estimator based on the configuration - multEstimator = 0; + // Determine the centrality estimator based on the configuration. if (EventConfig.cfgMultName.value == "FT0M") { - multEstimator = 0; + multEstimator = CentralityFT0M; } else if (EventConfig.cfgMultName.value == "FT0C") { - multEstimator = 1; + multEstimator = CentralityFT0C; } else if (EventConfig.cfgMultName.value == "FT0A") { - multEstimator = 2; + multEstimator = CentralityFT0A; + } else if (EventConfig.cfgMultName.value == "FV0A") { + multEstimator = CentralityFV0A; + } else { + LOGF(fatal, "Unsupported cfgMultName '%s'; choose FT0M, FT0C, FT0A, or FV0A", EventConfig.cfgMultName.value.c_str()); + } + LOGF(info, "Centrality estimator: %d, %s", multEstimator, EventConfig.cfgMultName.value.c_str()); + if (EventConfig.cfgMultiplicityEstimator.value < MultiplicityNTracksPV || + EventConfig.cfgMultiplicityEstimator.value > MultiplicityFV0A) { + LOG(fatal) << "cfgMultiplicityEstimator must be in the range [0, 6]"; } - LOGF(info, "Mult estimator: %d, %s", multEstimator, EventConfig.cfgMultName.value.c_str()); + LOGF(info, "Stored collision multiplicity estimator: %d", EventConfig.cfgMultiplicityEstimator.value); - // Ensure that only one process type is active at a time - if (doprocessRun3 && doprocessRun2) { - LOG(fatal) << "You cannot run both Run2 and Run3 processes at the same time"; + // Run 2 and Run 3 callbacks require different event-selection semantics. + const bool anyRun2Process = doprocessRun2 || doprocessRun2MC; + const bool anyRun3Process = doprocessRun3 || doprocessRun3MC || doprocessMCgen; + if (anyRun2Process && anyRun3Process) { + LOG(fatal) << "Run 2 and Run 3 processes cannot be enabled in the same ResonanceModuleInitializer"; } - if (doprocessRun2MC && doprocessRun3MC) { - LOG(fatal) << "You cannot run both Run2 and Run3 MC processes at the same time"; + if (doprocessRun2 && doprocessRun2MC) { + LOG(fatal) << "processRun2MC writes both ResoCollisions and ResoMCCollisions_001; do not enable processRun2 with it"; } - if (CCDB.cfgBypassCollIndexFill) { - LOG(fatal) << "cfgBypassCollIndexFill is incompatible with ResonanceDaughterInitializer"; + if (doprocessRun3 && doprocessRun3MC) { + LOG(fatal) << "processRun3MC writes both ResoCollisions and ResoMCCollisions_001; do not enable processRun3 with it"; + } + const int enabledGenMultiplicityEstimators = static_cast(GenCuts.cfgGenMult05.value) + + static_cast(GenCuts.cfgGenMult10.value) + + static_cast(GenCuts.cfgGenMultFT0M.value) + + static_cast(GenCuts.cfgGenMultFT0C.value) + + static_cast(GenCuts.cfgGenMultFV0A.value); + if ((doprocessMCgen || doprocessRun2MC || doprocessRun3MC) && enabledGenMultiplicityEstimators > 1) { + LOG(fatal) << "Only one generator multiplicity estimator can be enabled: cfgGenMult05, cfgGenMult10, cfgGenMultFT0M, cfgGenMultFT0C, or cfgGenMultFV0A"; + } + if (doprocessMCgen) { + if (GenCuts.cfgGenMultPercentile && multEstimator != CentralityFT0M && + multEstimator != CentralityFT0C && multEstimator != CentralityFV0A) { + LOGF(fatal, "cfgGenMultPercentile supports cfgMultName=FT0M, FT0C, or FV0A"); + } + if (GenCuts.isZvtxcutGen && + (!std::isfinite(GenCuts.cutzvertexGen.value) || GenCuts.cutzvertexGen.value <= 0.f)) { + LOG(fatal) << "cutzvertexGen must be finite and positive when the generator vertex cut is enabled"; + } + if (!std::isfinite(GenCuts.cfgRapidityCutMinGen.value) || + !std::isfinite(GenCuts.cfgRapidityCutMaxGen.value) || + GenCuts.cfgRapidityCutMinGen.value >= GenCuts.cfgRapidityCutMaxGen.value) { + LOG(fatal) << "Generator rapidity limits must be finite and satisfy cfgRapidityCutMinGen < cfgRapidityCutMaxGen"; + } } // Initialize event selection cuts based on the process type - if (doprocessRun2) { + if (anyRun2Process) { colCuts.setCuts(EventCuts.cfgEvtZvtx, EventCuts.cfgEvtTriggerCheck, EventCuts.cfgEvtOfflineCheck, false); - } else if (doprocessRun3) { + } else if (anyRun3Process) { colCuts.setCuts(EventCuts.cfgEvtZvtx, EventCuts.cfgEvtTriggerCheck, EventCuts.cfgEvtOfflineCheck, true, false, EventCuts.cfgEvtOccupancyInTimeRange); } colCuts.init(&qaRegistry); @@ -189,7 +326,31 @@ struct ResonanceModuleInitializer { colCuts.setApplyRun2AliEventCuts(EventCuts.cfgEvtRun2AliEventCuts); colCuts.setApplyRun2INELgtZERO(EventCuts.cfgEvtRun2INELgtZERO); - rctChecker.init(EventCuts.cfgEvtRCTFlagCheckerLabel, EventCuts.cfgEvtRCTFlagCheckerZDCCheck, EventCuts.cfgEvtRCTFlagCheckerLimitAcceptAsBad); + if (EventConfig.cfgFillDetailedQA && (doprocessRun3 || doprocessRun3MC)) { + AxisSpec selectionStageAxis{DetailedQAStages, -0.5f, static_cast(DetailedQAStages) - 0.5f, "Passed event-selection stage"}; + AxisSpec vertexAxis{EventConfig.cfgVtxBins, "Collision vertex z (cm)"}; + AxisSpec centralityAxis{EventConfig.binsCent, "Centrality (%)"}; + AxisSpec multiplicityAxis{EventConfig.binsMultiplicity, "Multiplicity"}; + qaRegistry.add("Event/h4EventSelectionDetail", "Event-selection cut flow", kTHnSparseD, + {selectionStageAxis, vertexAxis, centralityAxis, multiplicityAxis}); + + auto detailedQA = qaRegistry.get(HIST("Event/h4EventSelectionDetail")); + auto cutCounts = qaRegistry.get(HIST("CollCutCounts")); + for (int stage = o2::analysis::CollisonCuts::kAllEvent; + stage <= o2::analysis::CollisonCuts::kAllpassed; ++stage) { + detailedQA->GetAxis(0)->SetBinLabel(stage + 1, cutCounts->GetXaxis()->GetBinLabel(colCuts.binLabel(stage))); + } + detailedQA->GetAxis(0)->SetBinLabel(DetailedQARCTStage + 1, "RCT"); + } + + recoRCTChecker.init(EventCuts.cfgEvtRCTFlagCheckerLabel, + EventCuts.cfgEvtRCTFlagCheckerZDCCheck, + EventCuts.cfgEvtRCTFlagCheckerLimitAcceptAsBad, + EventCuts.cfgEvtRCTCheckTableValidity.value); + genRCTChecker.init(EventCuts.cfgEvtRCTFlagCheckerLabel, + EventCuts.cfgEvtRCTFlagCheckerZDCCheck, + EventCuts.cfgEvtRCTFlagCheckerLimitAcceptAsBad, + GenCuts.cfgGenRCTCheckTableValidity.value); // Configure CCDB access if not bypassed if (!EventConfig.cfgBypassCCDB) { @@ -201,11 +362,26 @@ struct ResonanceModuleInitializer { ccdb->setCreatedNotAfter(now); // TODO must become global parameter from the train creation time } - // Initialize QA histograms if required - if (EventConfig.cfgFillQA && (doprocessRun3MC || doprocessRun2MC)) { - AxisSpec centAxis = {EventConfig.binsCent, "Centrality (%)"}; - AxisSpec idxMCAxis = {26, -0.5, 25.5, "Index"}; - qaRegistry.add("Event/hMCEventIndices", "hMCEventIndices", kTH2D, {centAxis, idxMCAxis}); + if (doprocessMCgen) { + constexpr std::array MCEventLabels{"All", "z vertex", "BC RCT", "INEL", "INEL>0"}; + AxisSpec centAxisGen = {GenCuts.binsCentGen, "Centrality (%)"}; + AxisSpec eventTypeAxis = {2, 0.f, 2.f, "Event type"}; + qaRegistry.add("EventGen/hNEventsMC", "Generator event selection", kTH1D, {{5, 0.f, 5.f}}); + auto eventCounter = qaRegistry.get(HIST("EventGen/hNEventsMC")); + for (std::size_t i = 0; i < MCEventLabels.size(); ++i) { + eventCounter->GetXaxis()->SetBinLabel(static_cast(i + 1), MCEventLabels[i]); + } + qaRegistry.add("EventGen/h5ResonanceTruth", "Generated resonance", kTHnSparseD, + {eventTypeAxis, GenCuts.ptAxisGen, centAxisGen, GenCuts.multNTracksAxis, GenCuts.impactParameterAxis}); + qaRegistry.add("EventGen/h5ResonanceTruthAnti", "Generated anti-resonance", kTHnSparseD, + {eventTypeAxis, GenCuts.ptAxisGen, centAxisGen, GenCuts.multNTracksAxis, GenCuts.impactParameterAxis}); + qaRegistry.add("EventGen/hZCollisionGen", "Generator collision z vertex", kTH1D, {{100, -20.f, 20.f}}); + qaRegistry.add("EventGen/h4MultCent_genMC", "Generator-event multiplicity and centrality", kTHnSparseD, + {eventTypeAxis, centAxisGen, GenCuts.multNTracksAxis, GenCuts.impactParameterAxis}); + qaRegistry.add("EventGen/h4MultCent_recMC", "Reconstructed-event multiplicity and centrality", kTHnSparseD, + {eventTypeAxis, centAxisGen, GenCuts.multNTracksAxis, GenCuts.impactParameterAxis}); + qaRegistry.add("EventGen/h2CentralityVsMultMC", "Representative reconstructed centrality vs generator multiplicity", kTH2D, + {centAxisGen, GenCuts.multNTracksAxis}); } } @@ -214,7 +390,8 @@ struct ResonanceModuleInitializer { * * @param bc BC iterator */ - void initCCDB(aod::BCsWithTimestamps::iterator const& bc) // Simple copy from LambdaKzeroFinder.cxx + template + void initCCDB(BCType const& bc) // Simple copy from LambdaKzeroFinder.cxx { if (EventConfig.cfgBypassCCDB) { return; @@ -258,203 +435,352 @@ struct ResonanceModuleInitializer { LOGF(info, "Bz set to %f for run: ", dBz, mRunNumber); } - /** - * @brief Checks if the collision is INEL>0 - * - * @tparam MCPart Type of MC particles - * @param mcparts MC particles - * @return true if INEL>0, false otherwise - */ - template - bool isTrueINEL0(MCPart const& mcparts) - { - for (auto const& mcparticle : mcparts) { - if (!mcparticle.isPhysicalPrimary()) { - continue; - } - auto p = pdg->GetParticle(mcparticle.pdgCode()); - if (p != nullptr) { - if (std::abs(p->Charge()) >= MinimumChargedParticleCharge) { - if (std::abs(mcparticle.eta()) < 1) { - return true; - } - } - } - } - return false; - } - /** * @brief Centrality estimator selection * * @tparam ResoColl Type of resonance collision - * @tparam isMC Boolean indicating if it's MC * @param resoEvents Resonance events * @return Centrality value */ - template + template float centEst(ResoColl const& resoEvents) { - float returnValue = -999.0; switch (multEstimator) { - case 0: - returnValue = resoEvents.centFT0M(); - break; - case 1: - if constexpr (isMC) { - LOG(fatal) << "CentFT0C is not available for MC"; - return returnValue; - } else { - returnValue = resoEvents.centFT0C(); - break; - } - case 2: - if constexpr (isMC) { - LOG(fatal) << "CentFT0A is not available for MC"; - return returnValue; - } else { - returnValue = resoEvents.centFT0A(); - break; - } + case CentralityFT0M: + return resoEvents.centFT0M(); + case CentralityFT0C: + return resoEvents.centFT0C(); + case CentralityFT0A: + return resoEvents.centFT0A(); + case CentralityFV0A: + return resoEvents.centFV0A(); default: - returnValue = resoEvents.centFT0M(); - break; + return -999.f; + } + } + + /** + * @brief Returns the configured reconstructed multiplicity estimator + */ + template + float collisionMultiplicity(CollisionType const& collision) + { + switch (EventConfig.cfgMultiplicityEstimator.value) { + case MultiplicityNTracksPV: + return collision.multNTracksPV(); + case MultiplicityNTracksPVeta1: + return collision.multNTracksPVeta1(); + case MultiplicityNTracksPVetaHalf: + return collision.multNTracksPVetaHalf(); + case MultiplicityFT0M: + return collision.multFT0M(); + case MultiplicityFT0A: + return collision.multFT0A(); + case MultiplicityFT0C: + return collision.multFT0C(); + case MultiplicityFV0A: + return collision.multFV0A(); + default: + return -1.f; + } + } + + /// Fill the detailed Run 3 collision QA at one passed selection stage. + template + void fillDetailedCollisionQA(CollisionType const& collision, int selectionStage) + { + qaRegistry.fill(HIST("Event/h4EventSelectionDetail"), + selectionStage, + collision.posZ(), + centEst(collision), + collisionMultiplicity(collision)); + } + + /// Reproduce the configured Run 3 collision cut flow for detailed QA only. + /// The authoritative event decision remains CollisonCuts::isSelected(). + template + void fillDetailedRun3SelectionQA(CollisionType const& collision) + { + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::kAllEvent); + if (std::abs(collision.posZ()) > EventCuts.cfgEvtZvtx.value) { + return; + } + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::kFlagZvertex); + +// Keep this QA-only mapping synchronized with CollisonCuts without exposing +// its internal selection registry. +// NOLINTNEXTLINE(cppcoreguidelines-macro-usage) -- X-macro interface of EventSelectionFlagsMapping.def +#define EVSEL_FLAG(enumVal, member, defaultVal, evtSelEnum, setter, getter, label, desc) \ + if (colCuts.getSelection(o2::analysis::CollisonCuts::evtSelEnum)) { \ + if (!collision.selection_bit(o2::aod::evsel::enumVal)) { \ + return; \ + } \ + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::evtSelEnum); \ + } +#include "PWGLF/Utils/EventSelectionFlagsMapping.def" // NOLINT(build/include) +#undef EVSEL_FLAG + + if (EventCuts.cfgEvtOfflineCheck.value && !collision.sel8()) { + return; + } + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::kFlagSel8); + + if (EventCuts.cfgEvtOccupancyInTimeRange.value > 0 && + collision.trackOccupancyInTimeRange() > EventCuts.cfgEvtOccupancyInTimeRange.value) { + return; } - return returnValue; + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::kFlagOccupancy); + fillDetailedCollisionQA(collision, o2::analysis::CollisonCuts::kAllpassed); } - using GenMCCollisions = soa::Join; - float centEstMC(const GenMCCollisions::iterator& collision) { return centEst(collision); } + + using GenMCCollisions = soa::Join; + using Run3MCCollisions = soa::Join; + using Run2MCCollisions = soa::Join; + using GenRecoCollisions = soa::Join; + using BCsWithRCT = soa::Join; /** - * @brief Fills MC particles + * @brief Applies the configured RCT selection to a reconstructed collision * - * @tparam CollisionType Type of collision - * @tparam SelectedMCPartType Type of selected MC particles - * @tparam TotalMCParts Type of total MC particles - * @param collision Collision data - * @param mcParts Selected MC particles - * @param mcParticles Total MC particles + * The collision-level mode reads the RCT value stored in EVSEL. The BC-level + * mode reads the value directly from the nominal BC referenced by collision. */ - template - void fillMCParticles(CollisionType collision, SelectedMCPartType const& mcParts, TotalMCParts const& mcParticles) + template + bool isRecoRCTSelected(CollisionType const& collision) { - for (auto const& mcPart : mcParts) { - std::vector daughterPDGs; - if (mcPart.has_daughters()) { - auto daughter01 = mcParticles.rawIteratorAt(mcPart.daughtersIds()[0] - mcParticles.offset()); - auto daughter02 = mcParticles.rawIteratorAt(mcPart.daughtersIds()[1] - mcParticles.offset()); + if (!EventCuts.cfgEvtUseRCTFlagChecker.value) { + return true; + } + if (!EventCuts.cfgEvtBCRCT.value) { + return recoRCTChecker(collision); + } + const auto bc = collision.template bc_as(); + return recoRCTChecker(bc); + } + + /// Apply the reconstructed Run 3 event and RCT selections. Detailed cut-flow + /// QA is evaluated independently and never controls the event decision. + template + bool isRun3CollisionSelected(CollisionType const& collision) + { + if (EventConfig.cfgFillDetailedQA) { + fillDetailedRun3SelectionQA(collision); + } + if (!colCuts.isSelected(collision, EventConfig.cfgFillQA) || !isRecoRCTSelected(collision)) { + return false; + } + if (EventConfig.cfgFillDetailedQA) { + fillDetailedCollisionQA(collision, DetailedQARCTStage); + } + return true; + } + + /** + * @brief Fills generator-level resonance QA + * + * @tparam MCParticlesType Type of MC-particle group + * @param mcParticles MC particles grouped by generator collision + * @param generatorCentrality Generator or representative reconstructed centrality + * @param multiplicity Generator-level charged-particle multiplicity + * @param impactParameter Generator collision impact parameter + * @param eventType INEL/INEL>0 category + */ + template + void fillMCGenParticles(MCParticlesType const& mcParticles, + float generatorCentrality, + float multiplicity, + float impactParameter, + int eventType) + { + for (auto const& mcPart : mcParticles) { + if (std::abs(mcPart.pdgCode()) != std::abs(GenCuts.pdgTruthMother.value)) { + continue; + } + if (mcPart.y() <= GenCuts.cfgRapidityCutMinGen || mcPart.y() >= GenCuts.cfgRapidityCutMaxGen) { + continue; + } + + std::array daughterPDGs{-1, -1}; + const bool hasDaughters = mcPart.has_daughters(); + const auto daughterIds = mcPart.daughtersIds(); + if (hasDaughters) { + auto daughter01 = mcParticles.rawIteratorAt(daughterIds[0] - mcParticles.offset()); + auto daughter02 = mcParticles.rawIteratorAt(daughterIds[1] - mcParticles.offset()); daughterPDGs = {daughter01.pdgCode(), daughter02.pdgCode()}; + } + if (GenCuts.isDaughterCheck) { + const int daughter1 = std::abs(GenCuts.pdgTruthDaughter1.value); + const int daughter2 = std::abs(GenCuts.pdgTruthDaughter2.value); + const int firstDaughterPDG = std::abs(daughterPDGs[0]); + const int secondDaughterPDG = std::abs(daughterPDGs[1]); + // Match the configured decay to two distinct daughter slots. This is + // important when both configured absolute PDGs are equal (e.g. K+K-): + // one matching daughter must not satisfy both requirements. + const bool matchesConfiguredDecay = + (firstDaughterPDG == daughter1 && secondDaughterPDG == daughter2) || + (firstDaughterPDG == daughter2 && secondDaughterPDG == daughter1); + if (!matchesConfiguredDecay) { + continue; + } + } + + if (mcPart.pdgCode() > 0) { + qaRegistry.fill(HIST("EventGen/h5ResonanceTruth"), eventType, mcPart.pt(), generatorCentrality, multiplicity, impactParameter); } else { - daughterPDGs = {-1, -1}; + qaRegistry.fill(HIST("EventGen/h5ResonanceTruthAnti"), eventType, mcPart.pt(), generatorCentrality, multiplicity, impactParameter); + } + } + } + + /** + * @brief Fills generated resonance parents for one reduced collision + * + * @tparam SelectedMCParticlesType Type of the selected MC-particle slice + * @tparam MCParticlesType Type of the complete MC-particle table + * @param reducedCollisionId Reduced collision referenced by the output rows + * @param selectedParents Selected parent particles in the associated MC collision + * @param mcParticles Complete MC-particle table used to resolve daughter indices + * + * The source MC-particle global index is persisted as a scalar in + * ResoMCParents_001; it is not a relation requiring McParticles at merge time. + */ + template + void fillMCParents(int64_t reducedCollisionId, + SelectedMCParticlesType const& selectedParents, + MCParticlesType const& mcParticles) + { + for (auto const& mcPart : selectedParents) { + if (!GenCuts.cfgDoSignalLoss) { + const int absPdg = std::abs(mcPart.pdgCode()); + if (absPdg == PdgProton || absPdg == PdgLambda0 || absPdg == PdgXiMinus || absPdg == PdgXi0 || absPdg == PdgOmegaMinus) { + continue; + } + } + + std::array daughterPDGs{-1, -1}; + if (mcPart.has_daughters()) { + const auto daughter1 = mcParticles.rawIteratorAt(mcPart.daughtersIds()[0] - mcParticles.offset()); + const auto daughter2 = mcParticles.rawIteratorAt(mcPart.daughtersIds()[1] - mcParticles.offset()); + daughterPDGs = {daughter1.pdgCode(), daughter2.pdgCode()}; } - reso2mcparents(collision.globalIndex(), + reso2mcparents(reducedCollisionId, mcPart.globalIndex(), mcPart.pdgCode(), - daughterPDGs[0], daughterPDGs[1], + daughterPDGs[0], + daughterPDGs[1], mcPart.isPhysicalPrimary(), mcPart.producedByGenerator(), mcPart.pt(), mcPart.px(), mcPart.py(), mcPart.pz(), - mcPart.eta(), - mcPart.phi(), - mcPart.y()); - daughterPDGs.clear(); + mcPart.y(), + mcPart.e(), + mcPart.statusCode()); } } /** - * @brief Fills MC collision data - * - * @tparam isRun2 Boolean indicating if it's Run2 - * @tparam MCCol Type of MC collision - * @tparam MCPart Type of MC particles - * @param mccol MC collision data - * @param mcparts MC particles + * @brief Returns the configured generator-level multiplicity estimator */ - template - void fillMCCollision(MCCol const& mccol, MCPart const& mcparts) - { - const auto& mcColg = mccol.template mcCollision_as(); - float mcCent = 999.0; - if constexpr (isRun2) { - if (EventConfig.cfgCentralityMC == MCCentralityRecoEstimator) { - mcCent = mccol.centRun2V0M(); - } else { - mcCent = mcColg.impactParameter(); - } - } else { - if (EventConfig.cfgCentralityMC == MCCentralityRecoEstimator) { - mcCent = centEst(mccol); - } else if (EventConfig.cfgCentralityMC == MCCentralityGeneratorEstimator) { - mcCent = centEstMC(mcColg); - } else if (EventConfig.cfgCentralityMC == MCCentralityImpactParameterEstimator) { - mcCent = mcColg.impactParameter(); - } + template + float getMCMultiplicity(MCCollision const& mcCollision) + { + if (GenCuts.cfgGenMult05) { + return mcCollision.multMCNParticlesEta05(); + } + if (GenCuts.cfgGenMult10) { + return mcCollision.multMCNParticlesEta10(); + } + if (GenCuts.cfgGenMultFT0M) { + return mcCollision.multMCFT0A() + mcCollision.multMCFT0C(); + } + if (GenCuts.cfgGenMultFT0C) { + return mcCollision.multMCFT0C(); } - const bool inVtx10 = !(std::abs(mcColg.posZ()) > MCVertexZMax); - bool isTrueINELgt0 = isTrueINEL0(mcparts); - bool isTriggerTVX = mccol.selection_bit(aod::evsel::kIsTriggerTVX); - bool isSel8 = mccol.sel8(); - bool isSelected = colCuts.isSelected(mccol, EventConfig.cfgFillQA); - resoMCCollisions(inVtx10, isTrueINELgt0, isTriggerTVX, isSel8, isSelected, mcCent, -1.0f); + if (GenCuts.cfgGenMultFV0A) { + return mcCollision.multMCFV0A(); + } + return -1.f; + } - if (EventConfig.cfgFillQA) { - // QA for trigger efficiency - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kINEL); - if (inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kINEL10); - } - if (isTrueINELgt0) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kINELg0); - } - if (inVtx10 && isTrueINELgt0) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kINELg010); - } + /** + * @brief Returns the configured generator-level centrality percentile + */ + template + float getMCCentrality(MCCollision const& mcCollision) const + { + if (!GenCuts.cfgGenMultPercentile.value) { + return 100.5f; + } + switch (multEstimator) { + case CentralityFT0M: + return mcCollision.centFT0M(); + case CentralityFT0C: + return mcCollision.centFT0C(); + case CentralityFV0A: + return mcCollision.centFV0A(); + default: + return 100.5f; + } + } - // TVX MB trigger - if (isTriggerTVX) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kTrig); - } - if (isTriggerTVX && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kTrig10); - } - if (isTriggerTVX && isTrueINELgt0) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kTrigINELg0); - } - if (isTriggerTVX && isTrueINELgt0 && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kTrigINELg010); - } + /** + * @brief Fills the generator-only MC collision extension and optional source link + * + * Accepting only the generator collision here prevents reconstructed event + * selection state from accidentally entering ResoMCCollisions_001. + */ + template + void fillMCCollision001(MCCollision const& mcCollision) + { + const float mcMultiplicity = getMCMultiplicity(mcCollision); + const bool inTrueVtx10 = std::abs(mcCollision.posZ()) < MCVertexZMax; + // Keep the generator-level definition identical to MultsExtraMC: + // at least one physical-primary charged particle within |eta| < 1. + const bool isTrueINELgt0 = mcCollision.isInelGt0(); + resoMCCollisions001(inTrueVtx10, + isTrueINELgt0, + mcCollision.impactParameter(), + mcMultiplicity); + if (GenCuts.cfgFillMCCollisionSoftLink) { + resoMCCollisionIds(mcCollision.globalIndex()); + } + } - // Sel8 event selection - if (isSel8) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kSel8); - } - if (isSel8 && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kSel810); - } - if (isSel8 && isTrueINELgt0) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kSel8INELg0); - } - if (isSel8 && isTrueINELgt0 && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kSel8INELg010); - } + /// Write the Run 3 base collision, optional soft link, and scalar grouping key. + template + void fillRun3Collision(CollisionType const& collision) + { + const bool isRecINELgt0 = collision.isInelGt0(); + centrality = centEst(collision); + resoCollisions(collisionMultiplicity(collision), + collision.posX(), + collision.posY(), + collision.posZ(), + centrality, + dBz, + isRecINELgt0); + resoCollisionColls(collision.globalIndex()); + resoCollisionGroups001(collision.globalIndex()); + } - // CollisionCuts selection - if (isSelected) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kAllCuts); - } - if (isSelected && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kAllCuts10); - } - if (isSelected && isTrueINELgt0) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kAllCutsINELg0); - } - if (isSelected && isTrueINELgt0 && inVtx10) { - qaRegistry.fill(HIST("Event/hMCEventIndices"), mcCent, aod::resocollision::kAllCutsINELg010); - } - } + /// Write the Run 2 base collision, optional soft link, and scalar grouping key. + template + void fillRun2Collision(CollisionType const& collision) + { + centrality = collision.centRun2V0M(); + // The configured Run 3 estimators are not available in the Run 2 input + // schema. Preserve the previous zero-filled Run 2 behaviour. + resoCollisions(0.f, + collision.posX(), + collision.posY(), + collision.posZ(), + centrality, + dBz, + 0); + resoCollisionColls(collision.globalIndex()); + resoCollisionGroups001(collision.globalIndex()); } /** @@ -473,27 +799,19 @@ struct ResonanceModuleInitializer { * @param collision Collision data * @param bc BC data */ - void processRun3(soa::Filtered::iterator const& collision, - aod::BCsWithTimestamps const&) + void processRun3(aod::ResoCollisionCandidates::iterator const& collision, + BCsWithRCT const&) { - auto bc = collision.bc_as(); + auto bc = collision.bc_as(); initCCDB(bc); // Default event selection - if (!colCuts.isSelected(collision, EventConfig.cfgFillQA)) { - return; - } - if (EventCuts.cfgEvtUseRCTFlagChecker && !rctChecker(collision)) { + if (!isRun3CollisionSelected(collision)) { return; } if (EventConfig.cfgFillQA) { colCuts.fillQA(collision); } - const bool isRecINELgt0 = collision.isInelGt0(); - centrality = centEst(collision); - - resoCollisions(collision.multNTracksPV(), collision.multNTracksPVeta1(), collision.multNTracksPVetaHalf(), collision.posX(), collision.posY(), collision.posZ(), centEst(collision), dBz, isRecINELgt0); - resoCollisionColls(collision.globalIndex()); - resoCollisionGroups(collision.globalIndex()); + fillRun3Collision(collision); } PROCESS_SWITCH(ResonanceModuleInitializer, processRun3, "Default process for RUN3", false); @@ -503,7 +821,7 @@ struct ResonanceModuleInitializer { * @param collision Collision data * @param bc BC data */ - void processRun2(soa::Filtered::iterator const& collision, + void processRun2(aod::ResoRun2CollisionCandidates::iterator const& collision, aod::BCsWithRun2Info const&) { // auto bc = collision.bc_as(); @@ -514,41 +832,154 @@ struct ResonanceModuleInitializer { if (EventConfig.cfgFillQA) { colCuts.fillQARun2(collision); } - centrality = collision.centRun2V0M(); - - resoCollisions(0, 0, 0, collision.posX(), collision.posY(), collision.posZ(), centrality, dBz, 0); - resoCollisionColls(collision.globalIndex()); - resoCollisionGroups(collision.globalIndex()); + fillRun2Collision(collision); } PROCESS_SWITCH(ResonanceModuleInitializer, processRun2, "process for RUN2", false); + /** + * @brief Processes generator-level MC event and resonance QA + * + * The original MC collision is the grouping key. MC particles are therefore + * grouped automatically, while reconstructed collisions arrive as a 0..N + * SmallGroup through their McCollisionLabels relation. This auxiliary + * callback intentionally fills generator-level QA only and does not write + * reduced AOD tables. RCT quality is evaluated through the generator + * collision's associated BC because it is a run-condition property. + */ + void processMCgen(GenMCCollisions::iterator const& mcCollision, + aod::McParticles const& mcParticles, + soa::SmallGroups const& collisions, + BCsWithRCT const&) + { + auto bc = mcCollision.bc_as(); + initCCDB(bc); + + const auto getReconstructedCentrality = [&](auto const& collision) { + return centEst(collision); + }; + + const float generatorCentrality = getMCCentrality(mcCollision); + const float impactParameter = mcCollision.impactParameter(); + const float multiplicity = getMCMultiplicity(mcCollision); + + qaRegistry.fill(HIST("EventGen/hNEventsMC"), 0.5); + if (GenCuts.isZvtxcutGen && std::abs(mcCollision.posZ()) > GenCuts.cutzvertexGen) { + return; + } + qaRegistry.fill(HIST("EventGen/hNEventsMC"), 1.5); + if (GenCuts.cfgGenBCRCT && !genRCTChecker(bc)) { + return; + } + qaRegistry.fill(HIST("EventGen/hNEventsMC"), 2.5); + qaRegistry.fill(HIST("EventGen/hZCollisionGen"), mcCollision.posZ()); + + int eventType = 0; + qaRegistry.fill(HIST("EventGen/hNEventsMC"), 3.5); + if (GenCuts.checkIsTrueINELgt0 && mcCollision.isInelGt0()) { + eventType = 1; + qaRegistry.fill(HIST("EventGen/hNEventsMC"), 4.5); + } + + bool hasSelectedRecoCollision = false; + int largestNContributors = -1; + float reconstructedCentrality = 100.5f; + for (auto const& collision : collisions) { + if (!isRecoRCTSelected(collision)) { + continue; + } + if (!colCuts.isSelected(collision, false)) { + continue; + } + const int nContributors = static_cast(collision.multPVTotalContributors()); + if (nContributors > largestNContributors) { + largestNContributors = nContributors; + reconstructedCentrality = getReconstructedCentrality(collision); + } + hasSelectedRecoCollision = true; + } + + if (GenCuts.cfgGenMultPercentile) { + fillMCGenParticles(mcParticles, generatorCentrality, multiplicity, impactParameter, eventType); + qaRegistry.fill(HIST("EventGen/h4MultCent_genMC"), eventType, generatorCentrality, multiplicity, impactParameter); + } else { + fillMCGenParticles(mcParticles, reconstructedCentrality, multiplicity, impactParameter, eventType); + qaRegistry.fill(HIST("EventGen/h4MultCent_genMC"), eventType, reconstructedCentrality, multiplicity, impactParameter); + qaRegistry.fill(HIST("EventGen/h2CentralityVsMultMC"), reconstructedCentrality, multiplicity); + } + if (hasSelectedRecoCollision) { + qaRegistry.fill(HIST("EventGen/h4MultCent_recMC"), eventType, reconstructedCentrality, multiplicity, impactParameter); + } + } + PROCESS_SWITCH(ResonanceModuleInitializer, processMCgen, "Process generator-level MC QA", false); + /** * @brief Processes Run3 MC data * + * Reconstructed collisions pass the same event selection and QA sequence as + * Run 3 data before the reduced collision and its MC extension are written. + * * @param collision Collision data - * @param mcParticles MC particles * @param mcCollisions MC collisions + * @param mcParticles MC particles used to fill reduced resonance parents */ - void processRun3MC(soa::Filtered::iterator const& collision, - aod::McParticles const& mcParticles, GenMCCollisions const&) + void processRun3MC(aod::ResoCollisionCandidatesMC::iterator const& collision, + Run3MCCollisions const&, + aod::McParticles const& mcParticles, + BCsWithRCT const&) { - if (EventCuts.cfgEvtUseRCTFlagChecker && !rctChecker(collision)) { + auto bc = collision.bc_as(); + initCCDB(bc); + if (!isRun3CollisionSelected(collision)) { return; } - fillMCCollision(collision, mcParticles); + if (EventConfig.cfgFillQA) { + colCuts.fillQA(collision); + } + if (!collision.has_mcCollision()) { + return; + } + const auto& mcCollision = collision.mcCollision_as(); + + // ResoMCCollisions_001 is a positional extension of ResoCollisions. When + // enabled, ResoMCCollisionIds is written in the same callback and is 1:1. + fillRun3Collision(collision); + const int64_t reducedCollisionId = resoCollisions.lastIndex(); + fillMCCollision001(mcCollision); + + // A generator collision can be associated with more than one reconstructed + // collision. Write one parent set for each reduced collision, as in the + // legacy collision-wise producer, so ResoCollisionId remains unambiguous. + auto selectedParents = selectedMCParticles->sliceBy(mcParticlesPerMcCollision, mcCollision.globalIndex()); + fillMCParents(reducedCollisionId, selectedParents, mcParticles); } PROCESS_SWITCH(ResonanceModuleInitializer, processRun3MC, "process MC for RUN3", false); /** * @brief Processes Run2 MC data * + * Reconstructed collisions pass the same event selection and QA sequence as + * Run 2 data before the reduced collision and its MC extension are written. + * * @param collision Collision data - * @param mcParticles MC particles */ - void processRun2MC(soa::Filtered::iterator const& collision, - aod::McParticles const& mcParticles) + void processRun2MC(aod::ResoRun2CollisionCandidatesMC::iterator const& collision, + Run2MCCollisions const&) { - fillMCCollision(collision, mcParticles); + if (!colCuts.isSelected(collision, EventConfig.cfgFillQA)) { + return; + } + if (EventConfig.cfgFillQA) { + colCuts.fillQARun2(collision); + } + if (!collision.has_mcCollision()) { + return; + } + const auto& mcCollision = collision.mcCollision_as(); + + // Keep the base and MC extension one-to-one; the optional source link is + // written by fillMCCollision001 in the same order. + fillRun2Collision(collision); + fillMCCollision001(mcCollision); } PROCESS_SWITCH(ResonanceModuleInitializer, processRun2MC, "process MC for RUN2", false); }; @@ -572,8 +1003,107 @@ struct ResonanceDaughterInitializer { static constexpr int TrackSelectionGlobalWoDCA = 3; static constexpr int TrackSelectionQuality = 4; static constexpr int TrackSelectionInAcceptance = 5; + static constexpr int PairGateModeConfigured = 0; + static constexpr int PairGateModeEither = 1; static constexpr float MomentumQuantizationScale = 1000.f; static constexpr std::size_t StoredMCRelationCount = 2; + static constexpr std::size_t MaxCandidateDaughters = 3; + + /// Selected-candidate state and the optional global daughter-ID veto set. + /// By default only candidate existence is recorded and daughter reuse is + /// rejected later for each concrete pair through the stored trackId. + struct SelectedCandidateDaughters { + std::vector allDaughterIds; + bool hasSelectedCandidate = false; + bool useGlobalDaughterVeto = false; + + explicit SelectedCandidateDaughters(bool globalDaughterVeto = false) + : useGlobalDaughterVeto(globalDaughterVeto) + { + } + + void addCandidate() + { + hasSelectedCandidate = true; + } + + template + void addCandidate(std::array const& daughterIds) + { + static_assert(DaughterCount <= MaxCandidateDaughters); + allDaughterIds.insert(allDaughterIds.end(), daughterIds.begin(), daughterIds.end()); + hasSelectedCandidate = true; + } + + void finalize() + { + if (!useGlobalDaughterVeto) { + return; + } + std::sort(allDaughterIds.begin(), allDaughterIds.end()); + allDaughterIds.erase(std::unique(allDaughterIds.begin(), allDaughterIds.end()), allDaughterIds.end()); + } + + [[nodiscard]] bool accepts(int64_t trackId) const + { + if (!hasSelectedCandidate) { + return false; + } + if (useGlobalDaughterVeto) { + return !std::binary_search(allDaughterIds.begin(), allDaughterIds.end(), trackId); + } + return true; + } + }; + + /// Producer-side track retention for the candidate types enabled by a callback. + struct PairTrackSelection { + SelectedCandidateDaughters v0Candidates; + SelectedCandidateDaughters cascadeCandidates; + bool useV0Candidates = false; + bool useCascadeCandidates = false; + bool useGlobalDaughterVeto = false; + + explicit PairTrackSelection(bool globalDaughterVeto = false) + : v0Candidates(globalDaughterVeto), + cascadeCandidates(globalDaughterVeto), + useGlobalDaughterVeto(globalDaughterVeto) + { + } + + template + bool operator()(TrackType const& track) const + { + if (!useGlobalDaughterVeto) { + return (useV0Candidates && v0Candidates.hasSelectedCandidate) || + (useCascadeCandidates && cascadeCandidates.hasSelectedCandidate); + } + + const auto trackId = static_cast(track.globalIndex()); + bool hasCandidate = false; + if (useV0Candidates && v0Candidates.hasSelectedCandidate) { + hasCandidate = true; + if (!v0Candidates.accepts(trackId)) { + return false; + } + } + if (useCascadeCandidates && cascadeCandidates.hasSelectedCandidate) { + hasCandidate = true; + if (!cascadeCandidates.accepts(trackId)) { + return false; + } + } + return hasCandidate; + } + }; + + struct KeepAllTracks { + template + bool operator()(TrackType const&) const + { + return true; + } + }; UltraMicroPidSpecies ultraMicroPidSpecies = UltraMicroPidSpecies::Pion; bool warnedUltraMicroMomentumRange = false; @@ -586,8 +1116,9 @@ struct ResonanceDaughterInitializer { Preslice cascadesMCPerCollision = aod::cascdata::collisionId; Produces reso2trks; ///< Output table for resonance tracks Produces resoTrackTracks; ///< Output table for original track row IDs - Produces reso2microtrks; ///< Output table for resonance microtracks - Produces resoMicroTrackTracks; ///< Output table for original microtrack row IDs + Produces reso2microtrks; ///< Output table for resonance microtracks + Produces resoMicroTrackTracks; ///< Positional original-track soft links for microtracks + Produces reso2mcmicrotrks; ///< Positional MC extension for resonance microtracks Produces reso2ultramicrotrks; ///< Output table for resonance ultra-microtracks Produces resoUltraMicroTrackTracks; ///< Output table for original ultra-microtrack row IDs Produces reso2mctracks; ///< Output table for MC resonance tracks @@ -606,86 +1137,96 @@ struct ResonanceDaughterInitializer { struct : ConfigurableGroup { Configurable cfgCutEta{"cfgCutEta", 0.8f, "Eta range for tracks"}; Configurable cfgCutMinPt{"cfgCutMinPt", 0.1f, "Minimum pT for tracks (GeV/c)"}; - Configurable cfgCutMaxPt{"cfgCutMaxPt", 32.767f, "Maximum pT for tracks (GeV/c)"}; - Configurable pidnSigmaPreSelectionCut{"pidnSigmaPreSelectionCut", 5.0f, "TPC PID cut (loose, improve performance)"}; - Configurable mincrossedrows{"mincrossedrows", 70, "Minimum crossed rows for V0 daughter tracks"}; + Configurable cfgCutMaxPt{"cfgCutMaxPt", 999.0f, "Maximum pT for tracks (GeV/c)"}; + Configurable pidnSigmaPreSelectionCut{"pidnSigmaPreSelectionCut", 5.0f, "TPC PID half-width around the configured species mean (loose preselection)"}; + Configurable pidnSigmaPreSelectionCutTOF{"pidnSigmaPreSelectionCutTOF", 5.0f, "TOF PID half-width around the configured species mean (loose preselection)"}; + Configurable pidnSigmaPreSelectionMeanPion{"pidnSigmaPreSelectionMeanPion", 0.000f, "Offset for TPC PID mean for pions"}; + Configurable pidnSigmaPreSelectionMeanKaon{"pidnSigmaPreSelectionMeanKaon", 0.000f, "Offset for TPC PID mean for kaons"}; + Configurable pidnSigmaPreSelectionMeanProton{"pidnSigmaPreSelectionMeanProton", 0.000f, "Offset for TPC PID mean for protons"}; + Configurable pidnSigmaPreSelectionMeanTOFPion{"pidnSigmaPreSelectionMeanTOFPion", 0.000f, "Offset for TOF PID mean for pions"}; + Configurable pidnSigmaPreSelectionMeanTOFKaon{"pidnSigmaPreSelectionMeanTOFKaon", 0.000f, "Offset for TOF PID mean for kaons"}; + Configurable pidnSigmaPreSelectionMeanTOFProton{"pidnSigmaPreSelectionMeanTOFProton", 0.000f, "Offset for TOF PID mean for protons"}; + Configurable cfgUseTOFPIDPreSelection{"cfgUseTOFPIDPreSelection", false, + "Apply the TOF PID cut to tracks with TOF; false ignores TOF PID, and tracks without TOF use TPC only"}; Configurable trackSelection{"trackSelection", 3, "Track selection: 0 -> No Cut, 1 -> kGlobalTrack, 2 -> kGlobalTrackWoPtEta, 3 -> kGlobalTrackWoDCA, 4 -> kQualityTracks, 5 -> kInAcceptanceTracks"}; - Configurable cMaxDCArToPVcut{"cMaxDCArToPVcut", 2.0, "Track DCAr cut to PV Maximum"}; - Configurable cMaxDCAzToPVcut{"cMaxDCAzToPVcut", 2.0, "Track DCAz cut to PV Maximum"}; - Configurable cMinDCAzToPVcut{"cMinDCAzToPVcut", 0.0, "Track DCAz cut to PV Minimum"}; + Configurable cMaxDCArToPVcut{"cMaxDCArToPVcut", 0.5f, "Track DCAr cut to PV Maximum"}; + Configurable cMaxDCAzToPVcut{"cMaxDCAzToPVcut", 1.0f, "Track DCAz cut to PV Maximum"}; + Configurable cMinDCAzToPVcut{"cMinDCAzToPVcut", 0.0f, "Track DCAz cut to PV Minimum"}; Configurable cfgApplyTightDCAPtDepSelection{"cfgApplyTightDCAPtDepSelection", true, "Apply the pT-dependent tight DCA selection"}; Configurable cfgTightDCAOffset{"cfgTightDCAOffset", 0.004f, "Constant term of the tight DCA threshold (cm)"}; Configurable cfgTightDCAPtCoefficient{"cfgTightDCAPtCoefficient", 0.013f, "Coefficient of the pT-dependent tight DCA threshold"}; Configurable cfgTightDCAPtPower{"cfgTightDCAPtPower", 1.f, "Power in tight DCA = offset + coefficient / pT^power"}; } TrackCuts; - // V0 and V0-daughter cuts + // V0 and V0-daughter cuts : based on loose cuts in V0s production analysis in pp 13.6 TeV struct : ConfigurableGroup { + Configurable mincrossedrowsV0s{"mincrossedrowsV0s", 70, "Minimum crossed rows for V0 daughter tracks"}; Configurable cMinV0PosDCArToPVcut{"cMinV0PosDCArToPVcut", 0.05f, "V0 Positive Track DCAr cut to PV Minimum"}; Configurable cMinV0NegDCArToPVcut{"cMinV0NegDCArToPVcut", 0.05f, "V0 Negative Track DCAr cut to PV Minimum"}; - Configurable cMinV0Radius{"cMinV0Radius", 0.0, "Minimum V0 radius from PV"}; - Configurable cMaxV0Radius{"cMaxV0Radius", 200.0, "Maximum V0 radius from PV"}; - Configurable cMinV0CosPA{"cMinV0CosPA", 0.995, "Minimum V0 CosPA to PV"}; + Configurable cMinV0Radius{"cMinV0Radius", 0.9f, "Minimum V0 radius from PV"}; + Configurable cMaxV0Radius{"cMaxV0Radius", 200.0f, "Maximum V0 radius from PV"}; + Configurable cMinV0CosPA{"cMinV0CosPA", 0.95f, "Minimum V0 CosPA to PV"}; } V0Cuts; - // Cascade and cascade-daughter cuts + // Additional K0s and Lambda0 selections not covered by V0Cuts: based on loose cuts in V0s production analysis in pp 13.6 TeV struct : ConfigurableGroup { - Configurable cfgMinCrossedRowsCascBach{"cfgMinCrossedRowsCascBach", 70, "min crossed rows for bachelor track from cascade"}; - Configurable cMinCascBachDCArToPVcut{"cMinCascBachDCArToPVcut", 0.05f, "Cascade Bachelor Track DCAr cut to PV Minimum"}; - Configurable cMaxCascBachDCArToPVcut{"cMaxCascBachDCArToPVcut", 999.0f, "Cascade Bachelor Track DCAr cut to PV Maximum"}; - Configurable cMaxCascDCAV0Daughters{"cMaxCascDCAV0Daughters", 1.6, "Cascade DCA between V0 daughters Maximum"}; - Configurable cMaxCascDCACascDaughters{"cMaxCascDCACascDaughters", 1.6, "Cascade DCA between Casc daughters Maximum"}; - Configurable cMinCascV0CosPA{"cMinCascV0CosPA", 0.97, "Minimum Cascade V0 CosPA to PV"}; - Configurable cMaxCascV0Radius{"cMaxCascV0Radius", 200.0, "Maximum Cascade V0 radius from PV"}; - Configurable cMinCascV0Radius{"cMinCascV0Radius", 0.0, "Minimum Cascade V0 radius from PV"}; - Configurable cMinCascRadius{"cMinCascRadius", 0.0, "Minimum Cascade radius from PV"}; - Configurable cMaxCascRadius{"cMaxCascRadius", 200.0, "Maximum Cascade radius from PV"}; - Configurable cMinCascCosPA{"cMinCascCosPA", 0.97, "Minimum Cascade CosPA to PV"}; - Configurable cCascMassResol{"cCascMassResol", 999, "Cascade mass resolution"}; + Configurable cfgSecondaryRequire{"cfgSecondaryRequire", true, "Secondary cuts on/off"}; + Configurable cfgSecondaryArmenterosCut{"cfgSecondaryArmenterosCut", false, "cut on Armenteros-Podolanski graph"}; + Configurable cfgSecondaryCrossMassHypothesisCut{"cfgSecondaryCrossMassHypothesisCut", false, "Apply cut based on the lambda mass hypothesis"}; + Configurable cfgByPassDauPIDSelection{"cfgByPassDauPIDSelection", false, "Bypass TPC PID preselection for V0 daughters"}; + Configurable cfgSecondaryDauDCAMax{"cfgSecondaryDauDCAMax", 1.0f, "Maximum DCA between V0 daughters"}; + Configurable cfgSecondaryPtMin{"cfgSecondaryPtMin", 0.0f, "Minimum transverse momentum of Secondary"}; + Configurable cfgSecondaryRapidityMax{"cfgSecondaryRapidityMax", 0.5f, "Maximum rapidity of Secondary"}; + Configurable cfgSecondaryDCAtoPVMax{"cfgSecondaryDCAtoPVMax", 0.4f, "Maximum DCA Secondary to PV"}; + Configurable cfgSecondaryProperLifetimeMax{"cfgSecondaryProperLifetimeMax", 40.f, "Maximum Secondary Lifetime"}; + Configurable cfgSecondaryparamArmenterosCut{"cfgSecondaryparamArmenterosCut", 0.2f, "parameter for Armenteros Cut"}; + Configurable cfgSecondaryMassWindow{"cfgSecondaryMassWindow", 0.03f, "Secondary inv mass selection window (GeV/c^2)"}; + Configurable cfgSecondaryCrossMassCutWindow{"cfgSecondaryCrossMassCutWindow", 0.02f, "Secondary inv mass selection window with (anti)lambda hypothesis (GeV/c^2)"}; + } SecondaryCuts; + + // Cascade and cascade-daughter cuts: based on loose cuts in cascades production analysis in pp 13.6 TeV + struct : ConfigurableGroup { + Configurable cfgMinCrossedRowsCascBach{"cfgMinCrossedRowsCascBach", 50, "min crossed rows for bachelor track from cascade"}; + Configurable cMinCascBachDCArToPVcut{"cMinCascBachDCArToPVcut", 0.05f, "Cascade Bachelor Track DCAr cut to PV Minimum"}; + Configurable cMaxCascBachDCArToPVcut{"cMaxCascBachDCArToPVcut", 200.0f, "Cascade Bachelor Track DCAr cut to PV Maximum"}; + Configurable cMaxCascDCAV0Daughters{"cMaxCascDCAV0Daughters", 0.5f, "Cascade DCA between V0 daughters Maximum"}; + Configurable cMaxCascDCACascDaughters{"cMaxCascDCACascDaughters", 1.2f, "Cascade DCA between Casc daughters Maximum"}; + Configurable cMinCascV0CosPA{"cMinCascV0CosPA", 0.98f, "Minimum Cascade V0 CosPA to PV"}; + Configurable cMaxCascV0Radius{"cMaxCascV0Radius", 200.0f, "Maximum Cascade V0 radius from PV"}; + Configurable cMinCascV0Radius{"cMinCascV0Radius", 0.4f, "Minimum Cascade V0 radius from PV"}; + Configurable cMinCascRadius{"cMinCascRadius", 0.4f, "Minimum Cascade radius from PV"}; + Configurable cMaxCascRadius{"cMaxCascRadius", 200.0f, "Maximum Cascade radius from PV"}; + Configurable cMinCascCosPA{"cMinCascCosPA", 0.99, "Minimum Cascade CosPA to PV"}; + Configurable cMaxXiMassWindow{"cMaxXiMassWindow", 0.02f, "Xi mass Window (GeV/c^2)"}; } CascadeCuts; // Derived dataset selections struct : ConfigurableGroup { - Configurable cfgFillPionTracks{"cfgFillPionTracks", false, "Fill pion tracks"}; - Configurable cfgFillKaonTracks{"cfgFillKaonTracks", false, "Fill kaon tracks"}; - Configurable cfgFillProtonTracks{"cfgFillProtonTracks", false, "Fill proton tracks"}; - Configurable cfgFillPionMicroTracks{"cfgFillPionMicroTracks", false, "Fill pion micro tracks"}; - Configurable cfgFillKaonMicroTracks{"cfgFillKaonMicroTracks", false, "Fill kaon micro tracks"}; - Configurable cfgFillProtonMicroTracks{"cfgFillProtonMicroTracks", false, "Fill proton micro tracks"}; - Configurable cfgFillPionUltraMicroTracks{"cfgFillPionUltraMicroTracks", true, "Fill pion ultra micro tracks"}; - Configurable cfgFillKaonUltraMicroTracks{"cfgFillKaonUltraMicroTracks", false, "Fill kaon ultra micro tracks"}; - Configurable cfgFillProtonUltraMicroTracks{"cfgFillProtonUltraMicroTracks", false, "Fill proton ultra micro tracks"}; - Configurable cfgFillK0s{"cfgFillK0s", false, "Fill K0s"}; - Configurable cfgFillLambda0{"cfgFillLambda0", false, "Fill Lambda0"}; - Configurable cfgBypassNoPairV0s{"cfgBypassNoPairV0s", false, "In a *WithPairGate process, bypass track fill if no V0 passes the configured selections"}; - Configurable cfgBypassNoPairCascades{"cfgBypassNoPairCascades", true, "In a *WithPairGate process, bypass track fill if no cascade passes the configured selections"}; + Configurable cfgFillPionTracks{"cfgFillPionTracks", false, "Apply the pion PID filter to every enabled track table"}; + Configurable cfgFillKaonTracks{"cfgFillKaonTracks", false, "Apply the kaon PID filter to every enabled track table"}; + Configurable cfgFillProtonTracks{"cfgFillProtonTracks", false, "Apply the proton PID filter to every enabled track table"}; + Configurable cfgFillK0s{"cfgFillK0s", true, "Fill K0s"}; + Configurable cfgFillLambda0{"cfgFillLambda0", true, "Fill Lambda0"}; + Configurable cfgPairGateMode{"cfgPairGateMode", 0, + "Combined *WithPairGate mode: 0 -> require every enabled candidate type, " + "1 -> require a selected V0 or cascade"}; + Configurable cfgBypassNoPairV0s{"cfgBypassNoPairV0s", false, + "Require a selected V0 and at least one selected collision track in pair-gate mode 0"}; + Configurable cfgBypassNoPairCascades{"cfgBypassNoPairCascades", true, + "Require a selected cascade and at least one selected collision track in pair-gate mode 0"}; + Configurable cfgGlobalDaughterVeto{"cfgGlobalDaughterVeto", false, + "Pair-gate callbacks only: reject a track if its original ID is a daughter " + "of any selected V0 or cascade; false stores all selected collision tracks " + "for mandatory candidate-local rejection through trackId"}; Configurable cfgFillMicroTracks{"cfgFillMicroTracks", false, "Fill micro tracks"}; Configurable cfgFillUltraMicroTracks{"cfgFillUltraMicroTracks", false, "Fill ultra micro tracks"}; Configurable cfgBypassTrackFill{"cfgBypassTrackFill", false, "Bypass the full ResoTracks table fill"}; - Configurable cfgBypassTrackIndexFill{"cfgBypassTrackIndexFill", false, "Bypass original daughter ID table fill"}; + Configurable cfgBypassTrackIndexFill{"cfgBypassTrackIndexFill", false, + "Bypass optional source-object soft-link side tables; " + "ResoMicroTracks_001 keeps its inline scalar trackId, but Full/Ultra " + "outputs then cannot perform candidate-local daughter rejection"}; } FilterForDerivedTables; - // Secondary selections for K0s and Lambda0 - struct : ConfigurableGroup { - Configurable cfgSecondaryRequire{"cfgSecondaryRequire", false, "Secondary cuts on/off"}; - Configurable cfgSecondaryArmenterosCut{"cfgSecondaryArmenterosCut", false, "cut on Armenteros-Podolanski graph"}; - Configurable cfgSecondaryCrossMassHypothesisCut{"cfgSecondaryCrossMassHypothesisCut", false, "Apply cut based on the lambda mass hypothesis"}; - Configurable cfgByPassDauPIDSelection{"cfgByPassDauPIDSelection", true, "Bypass TPC PID preselection for V0 daughters"}; - Configurable cfgSecondaryDauDCAMax{"cfgSecondaryDauDCAMax", 0.2, "Maximum DCA Secondary daughters to PV"}; - Configurable cfgSecondaryDauPosDCAtoPVMin{"cfgSecondaryDauPosDCAtoPVMin", 0.0, "Minimum DCA Secondary positive daughters to PV"}; - Configurable cfgSecondaryDauNegDCAtoPVMin{"cfgSecondaryDauNegDCAtoPVMin", 0.0, "Minimum DCA Secondary negative daughters to PV"}; - Configurable cfgSecondaryPtMin{"cfgSecondaryPtMin", 0.f, "Minimum transverse momentum of Secondary"}; - Configurable cfgSecondaryRapidityMax{"cfgSecondaryRapidityMax", 0.5, "Maximum rapidity of Secondary"}; - Configurable cfgSecondaryRadiusMin{"cfgSecondaryRadiusMin", 0.0, "Minimum transverse radius of Secondary"}; - Configurable cfgSecondaryRadiusMax{"cfgSecondaryRadiusMax", 999.9, "Maximum transverse radius of Secondary"}; - Configurable cfgSecondaryCosPAMin{"cfgSecondaryCosPAMin", 0.998, "Mininum cosine pointing angle of Secondary"}; - Configurable cfgSecondaryDCAtoPVMax{"cfgSecondaryDCAtoPVMax", 0.4, "Maximum DCA Secondary to PV"}; - Configurable cfgSecondaryProperLifetimeMax{"cfgSecondaryProperLifetimeMax", 20., "Maximum Secondary Lifetime"}; - Configurable cfgSecondaryparamArmenterosCut{"cfgSecondaryparamArmenterosCut", 0.2, "parameter for Armenteros Cut"}; - Configurable cfgSecondaryMassWindow{"cfgSecondaryMassWindow", 0.03, "Secondary inv mass selection window"}; - Configurable cfgSecondaryCrossMassCutWindow{"cfgSecondaryCrossMassCutWindow", 0.05, "Secondary inv mass selection window with (anti)lambda hypothesis"}; - } SecondaryCuts; - // Track selection filter based on configuration Filter trackFilter = (TrackCuts.trackSelection.node() == TrackSelectionNone) || ((TrackCuts.trackSelection.node() == TrackSelectionGlobal) && requireGlobalTrackInFilter()) || // kGlobalTrack = kQualityTracks | kPrimaryTracks | kInAcceptanceTracks @@ -698,11 +1239,17 @@ struct ResonanceDaughterInitializer { aod::track::pt <= TrackCuts.cfgCutMaxPt; HistogramRegistry qaRegistry{"QAHistos", {}, OutputObjHandlingPolicy::AnalysisObject}; - // The daughter task needs this row-wise mapping back to the original collision. - // Keep ResonanceModuleInitializer::cfgBypassCollIndexFill disabled and enable - // the matching Run 2/Run 3 base event process for MC workflows. - using ResoCollisionWithIndex = soa::Join; - using SelectedResoCollisions = soa::Join; + // The daughter task treats fIndexResoCollisions as a scalar v001 row ID. The + // data-model relation keeps its legacy v000 default target; v001 consumers + // must use resoCollision_as() with the exact bound v001 table type for + // dereferencing (aod::ResoCollisions_001 when the parent is not joined). + // Original-collision grouping is provided by the scalar-only version-1 + // ResoCollisionGroups_001 table, avoiding a hard source-AO2D relation. + // Collision mappings are always written; cfgBypassCollIndexFill is retained + // only so existing configuration files remain accepted. + using ResoCollisionWithIndex = soa::Join; + using SelectedResoCollisions = soa::Join; + PresliceUnsorted reducedCollisionsPerOriginalCollision = aod::resocollisiongroup001::originalCollisionId; /** * @brief Initializes the task @@ -711,36 +1258,107 @@ struct ResonanceDaughterInitializer { */ void init(InitContext&) { - const bool processTrackDataEnabled = doprocessData || doprocessDataHybrid || doprocessDataWithPairGate; - const bool processTrackMCEnabled = doprocessMC || doprocessMCWithPairGate; - const bool processV0DataEnabled = doprocessV0Data || doprocessV0DataHybrid; - const bool processCascDataEnabled = doprocessCascData || doprocessCascDataHybrid; + if (FilterForDerivedTables.cfgPairGateMode.value < PairGateModeConfigured || + FilterForDerivedTables.cfgPairGateMode.value > PairGateModeEither) { + LOGF(fatal, "cfgPairGateMode must be 0 (configured gates) or 1 (V0 or cascade)"); + } + const bool useEitherPairGate = FilterForDerivedTables.cfgPairGateMode.value == PairGateModeEither; + const bool processTrackDataEnabled = doprocessData || doprocessDataHybrid || doprocessDataWithPairGate || + doprocessDataWithV0PairGate || doprocessDataWithCascPairGate; + const bool processTrackMCEnabled = doprocessMC || doprocessMCWithPairGate || + doprocessMCWithV0PairGate || doprocessMCWithCascPairGate; + const bool pairTrackProcessEnabled = doprocessDataWithPairGate || doprocessDataWithV0PairGate || + doprocessDataWithCascPairGate || doprocessMCWithPairGate || + doprocessMCWithV0PairGate || doprocessMCWithCascPairGate; + const bool processV0DataEnabled = doprocessV0Data; + const bool processCascDataEnabled = doprocessCascData; + const bool anyDataProcessEnabled = processTrackDataEnabled || processV0DataEnabled || processCascDataEnabled; + const bool anyMCProcessEnabled = processTrackMCEnabled || doprocessV0MC || doprocessCascMC; const int enabledTrackProcesses = static_cast(doprocessData) + static_cast(doprocessDataHybrid) + static_cast(doprocessDataWithPairGate) + + static_cast(doprocessDataWithV0PairGate) + + static_cast(doprocessDataWithCascPairGate) + static_cast(doprocessMC) + - static_cast(doprocessMCWithPairGate); + static_cast(doprocessMCWithPairGate) + + static_cast(doprocessMCWithV0PairGate) + + static_cast(doprocessMCWithCascPairGate); if (enabledTrackProcesses > 1) { LOGF(fatal, "Only one track process can be enabled in ResonanceDaughterInitializer"); } - if (static_cast(doprocessV0Data) + static_cast(doprocessV0DataHybrid) + static_cast(doprocessV0MC) > 1) { + if (pairTrackProcessEnabled && + FilterForDerivedTables.cfgBypassTrackFill && + !FilterForDerivedTables.cfgFillMicroTracks && + !FilterForDerivedTables.cfgFillUltraMicroTracks) { + LOGF(fatal, "A pair-gate process requires at least one enabled Full, Micro, or UltraMicro track output"); + } + if (pairTrackProcessEnabled && !FilterForDerivedTables.cfgGlobalDaughterVeto && + FilterForDerivedTables.cfgBypassTrackIndexFill && + (!FilterForDerivedTables.cfgBypassTrackFill || FilterForDerivedTables.cfgFillUltraMicroTracks)) { + LOGF(warn, + "Default pair-gate mode defers daughter reuse rejection to analysis trackId comparisons, " + "but cfgBypassTrackIndexFill removes that ID from Full/Ultra track outputs"); + } + if (useEitherPairGate && + (doprocessDataWithV0PairGate || doprocessDataWithCascPairGate || + doprocessMCWithV0PairGate || doprocessMCWithCascPairGate)) { + LOGF(fatal, "cfgPairGateMode 1 requires the combined processDataWithPairGate or processMCWithPairGate callback"); + } + if (static_cast(doprocessV0Data) + static_cast(doprocessV0MC) > 1) { LOGF(fatal, "Only one V0 process can be enabled in ResonanceDaughterInitializer"); } - if (static_cast(doprocessCascData) + static_cast(doprocessCascDataHybrid) + static_cast(doprocessCascMC) > 1) { + if (static_cast(doprocessCascData) + static_cast(doprocessCascMC) > 1) { LOGF(fatal, "Only one cascade process can be enabled in ResonanceDaughterInitializer"); } if ((doprocessData || doprocessDataHybrid || doprocessMC) && - (FilterForDerivedTables.cfgBypassNoPairV0s || FilterForDerivedTables.cfgBypassNoPairCascades)) { + (useEitherPairGate || FilterForDerivedTables.cfgBypassNoPairV0s || FilterForDerivedTables.cfgBypassNoPairCascades || + FilterForDerivedTables.cfgGlobalDaughterVeto)) { LOGF(warn, "Pair-gate options are ignored by processData/processDataHybrid/processMC; enable the matching *WithPairGate process to apply them"); } - if (doprocessDataWithPairGate && FilterForDerivedTables.cfgBypassNoPairV0s && !processV0DataEnabled) { - LOGF(fatal, "cfgBypassNoPairV0s requires processV0Data or processV0DataHybrid so an accepted V0 is written for every retained collision"); - } - if (doprocessDataWithPairGate && FilterForDerivedTables.cfgBypassNoPairCascades && !processCascDataEnabled) { - LOGF(fatal, "cfgBypassNoPairCascades requires processCascData or processCascDataHybrid so an accepted cascade is written for every retained collision"); + const auto validatePairGateOutputs = [&](bool pairProcessEnabled, + bool v0OutputEnabled, + bool cascadeOutputEnabled, + char const* processName) { + if (!pairProcessEnabled) { + return; + } + if (useEitherPairGate) { + if (!v0OutputEnabled || !cascadeOutputEnabled) { + LOGF(fatal, "%s with pair-gate mode 1 requires both V0 and cascade output processes", processName); + } + } else { + if (!FilterForDerivedTables.cfgBypassNoPairV0s && !FilterForDerivedTables.cfgBypassNoPairCascades) { + LOGF(fatal, "%s with pair-gate mode 0 requires at least one enabled V0/cascade gate", processName); + } + if (FilterForDerivedTables.cfgBypassNoPairV0s && !v0OutputEnabled) { + LOGF(fatal, "%s requires a V0 output process when cfgBypassNoPairV0s is enabled", processName); + } + if (FilterForDerivedTables.cfgBypassNoPairCascades && !cascadeOutputEnabled) { + LOGF(fatal, "%s requires a cascade output process when cfgBypassNoPairCascades is enabled", processName); + } + } + }; + validatePairGateOutputs(doprocessDataWithPairGate, + processV0DataEnabled, + processCascDataEnabled, + "processDataWithPairGate"); + validatePairGateOutputs(doprocessMCWithPairGate, + doprocessV0MC, + doprocessCascMC, + "processMCWithPairGate"); + if (doprocessDataWithV0PairGate && !processV0DataEnabled) { + LOGF(fatal, "processDataWithV0PairGate requires processV0Data"); + } + if (doprocessDataWithCascPairGate && !processCascDataEnabled) { + LOGF(fatal, "processDataWithCascPairGate requires processCascData"); + } + if (doprocessMCWithV0PairGate && !doprocessV0MC) { + LOGF(fatal, "processMCWithV0PairGate requires processV0MC"); + } + if (doprocessMCWithCascPairGate && !doprocessCascMC) { + LOGF(fatal, "processMCWithCascPairGate requires processCascMC"); } - if (!std::isfinite(TrackCuts.cfgCutMinPt.value) || !std::isfinite(TrackCuts.cfgCutMaxPt.value) || TrackCuts.cfgCutMinPt.value < 0.f || @@ -765,7 +1383,28 @@ struct ResonanceDaughterInitializer { TrackCuts.pidnSigmaPreSelectionCut.value < 0.f) { LOGF(fatal, "pidnSigmaPreSelectionCut must be finite and non-negative"); } - + const std::array tpcPidMeans{TrackCuts.pidnSigmaPreSelectionMeanPion.value, + TrackCuts.pidnSigmaPreSelectionMeanKaon.value, + TrackCuts.pidnSigmaPreSelectionMeanProton.value}; + if (!std::all_of(tpcPidMeans.begin(), tpcPidMeans.end(), [](float mean) { + return std::isfinite(mean); + })) { + LOGF(fatal, "All TPC PID preselection means must be finite"); + } + if (TrackCuts.cfgUseTOFPIDPreSelection.value) { + if (!std::isfinite(TrackCuts.pidnSigmaPreSelectionCutTOF.value) || + TrackCuts.pidnSigmaPreSelectionCutTOF.value < 0.f) { + LOGF(fatal, "pidnSigmaPreSelectionCutTOF must be finite and non-negative when TOF PID selection is enabled"); + } + const std::array tofPidMeans{TrackCuts.pidnSigmaPreSelectionMeanTOFPion.value, + TrackCuts.pidnSigmaPreSelectionMeanTOFKaon.value, + TrackCuts.pidnSigmaPreSelectionMeanTOFProton.value}; + if (!std::all_of(tofPidMeans.begin(), tofPidMeans.end(), [](float mean) { + return std::isfinite(mean); + })) { + LOGF(fatal, "All TOF PID preselection means must be finite when TOF PID selection is enabled"); + } + } if (TrackCuts.cfgApplyTightDCAPtDepSelection.value && (!std::isfinite(TrackCuts.cfgTightDCAOffset.value) || !std::isfinite(TrackCuts.cfgTightDCAPtCoefficient.value) || @@ -775,28 +1414,38 @@ struct ResonanceDaughterInitializer { TrackCuts.cfgTightDCAPtPower.value < 0.f)) { LOGF(fatal, "Tight-DCA offset, pT coefficient, and power must be finite and non-negative"); } - + if (!std::isfinite(V0Cuts.cMinV0CosPA.value) || + V0Cuts.cMinV0CosPA.value < -1. || V0Cuts.cMinV0CosPA.value >= 1.) { + LOGF(fatal, "cMinV0CosPA must be finite and satisfy -1 <= cMinV0CosPA < 1"); + } + if (!std::isfinite(CascadeCuts.cMinCascCosPA.value) || + CascadeCuts.cMinCascCosPA.value < -1. || CascadeCuts.cMinCascCosPA.value >= 1.) { + LOGF(fatal, "cMinCascCosPA must be finite and satisfy -1 <= cMinCascCosPA < 1"); + } + if (!std::isfinite(CascadeCuts.cMinCascRadius.value) || + !std::isfinite(CascadeCuts.cMaxCascRadius.value) || + CascadeCuts.cMinCascRadius.value < 0. || + CascadeCuts.cMaxCascRadius.value <= CascadeCuts.cMinCascRadius.value) { + LOGF(fatal, "Cascade radius limits must be finite and satisfy 0 <= cMinCascRadius < cMaxCascRadius"); + } if (FilterForDerivedTables.cfgFillUltraMicroTracks) { - constexpr float MaxQuantizedMomentum = static_cast(std::numeric_limits::max()) / MomentumQuantizationScale; - const float maxLongitudinalMomentum = TrackCuts.cfgCutMaxPt.value * std::sinh(TrackCuts.cfgCutEta.value); - if (TrackCuts.cfgCutMaxPt.value > MaxQuantizedMomentum || - !std::isfinite(maxLongitudinalMomentum) || - maxLongitudinalMomentum > MaxQuantizedMomentum) { - LOGF(fatal, "Ultra-micro cfgCutMaxPt/cfgCutEta allow momentum components beyond the int16_t quantization range"); - } int enabledUltraMicroSpecies = 0; - enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillPionUltraMicroTracks ? 1 : 0; - enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillKaonUltraMicroTracks ? 1 : 0; - enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillProtonUltraMicroTracks ? 1 : 0; + enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillPionTracks ? 1 : 0; + enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillKaonTracks ? 1 : 0; + enabledUltraMicroSpecies += FilterForDerivedTables.cfgFillProtonTracks ? 1 : 0; if (enabledUltraMicroSpecies != 1) { - LOGF(fatal, "Exactly one pion/kaon/proton PID species must be enabled when filling ultra-micro tracks"); + LOGF(fatal, "Exactly one of cfgFillPionTracks, cfgFillKaonTracks, or cfgFillProtonTracks must be enabled when filling ultra-micro tracks"); } if (TrackCuts.pidnSigmaPreSelectionCut.value > o2::aod::resoultramicrodaughter::PidNSigma::MaxNSigma) { - LOGF(fatal, "Ultra-micro PID encoding requires pidnSigmaPreSelectionCut <= 5"); + LOGF(fatal, "Ultra-micro TPC PID encoding requires pidnSigmaPreSelectionCut <= 5"); } - if (FilterForDerivedTables.cfgFillKaonUltraMicroTracks) { + if (TrackCuts.cfgUseTOFPIDPreSelection.value && + TrackCuts.pidnSigmaPreSelectionCutTOF.value > o2::aod::resoultramicrodaughter::PidNSigma::MaxNSigma) { + LOGF(fatal, "Ultra-micro TOF PID encoding requires pidnSigmaPreSelectionCutTOF <= 5"); + } + if (FilterForDerivedTables.cfgFillKaonTracks) { ultraMicroPidSpecies = UltraMicroPidSpecies::Kaon; - } else if (FilterForDerivedTables.cfgFillProtonUltraMicroTracks) { + } else if (FilterForDerivedTables.cfgFillProtonTracks) { ultraMicroPidSpecies = UltraMicroPidSpecies::Proton; } else { ultraMicroPidSpecies = UltraMicroPidSpecies::Pion; @@ -804,29 +1453,14 @@ struct ResonanceDaughterInitializer { } if (cfgFillQA) { - AxisSpec idxAxis = {8, 0.0, 8.0, "Index"}; + AxisSpec idxAxis = {8, 0.0, 8.0, "Cumulative selection stage"}; AxisSpec ptAxis = {100, 0.0f, 10.0f, "#it{p}_{T} (GeV/#it{c})"}; - // The DCA maps cover the full configured pT range in 0.1 GeV/c steps. - constexpr float DcaPtBinWidth = 0.1f; - const int maxDcaPtBin = static_cast(std::ceil(TrackCuts.cfgCutMaxPt.value / DcaPtBinWidth)); - const int nDcaPtBins = maxDcaPtBin + 1; - const float dcaPtAxisHalfBin = 0.5f * DcaPtBinWidth; - AxisSpec dcaPtAxis = {nDcaPtBins, - -dcaPtAxisHalfBin, - maxDcaPtBin * DcaPtBinWidth + dcaPtAxisHalfBin, - "#it{p}_{T} (GeV/#it{c})"}; + AxisSpec dcaPtAxis = {300, 0.f, 30.f, "#it{p}_{T} (GeV/#it{c})"}; AxisSpec etaAxis = {100, -1.0f, 1.0f, "#eta"}; - AxisSpec phiAxis = {100, 0.0f, TwoPI, "#phi"}; - // Keep the configured signed DCA limits at bin centres so tracks exactly - // on an accepted selection boundary are not moved to overflow. - constexpr int NDcaBins = 201; - constexpr float DcaAxisPaddingFraction = 1.f / 200.f; - const auto configuredDcaXYMax = static_cast(TrackCuts.cMaxDCArToPVcut.value); - const auto configuredDcaZMax = static_cast(TrackCuts.cMaxDCAzToPVcut.value); - const float dcaXYAxisHalfRange = configuredDcaXYMax > 0.f ? configuredDcaXYMax * (1.f + DcaAxisPaddingFraction) : 1.e-4f; - const float dcaZAxisHalfRange = configuredDcaZMax > 0.f ? configuredDcaZMax * (1.f + DcaAxisPaddingFraction) : 1.e-4f; - AxisSpec dcaXYAxis = {NDcaBins, -dcaXYAxisHalfRange, dcaXYAxisHalfRange, "DCA_{xy} (cm)"}; - AxisSpec dcaZAxis = {NDcaBins, -dcaZAxisHalfRange, dcaZAxisHalfRange, "DCA_{z} (cm)"}; + // Keep the azimuthal QA bins aligned with the 18 TPC sectors (10 bins per sector). + AxisSpec phiAxis = {180, 0.0f, TwoPI, "#phi"}; + AxisSpec dcaXYAxis = {1000, -0.5, 0.5, "DCA_{xy} (cm)"}; + AxisSpec dcaZAxis = {1000, -0.5, 0.5, "DCA_{z} (cm)"}; // Tiny PID has 255 discrete values from -6.35 to 6.35 in 0.05 steps. // Extend the histogram limits by half a step so every decoded value is // located at a bin centre instead of on a floating-point bin boundary. @@ -845,28 +1479,34 @@ struct ResonanceDaughterInitializer { if (processTrackDataEnabled || processTrackMCEnabled) { qaRegistry.add("QA/hGoodTrackIndices", "hGoodTrackIndices", kTH1D, {idxAxis}); + auto trackSelection = qaRegistry.get(HIST("QA/hGoodTrackIndices")); + trackSelection->GetXaxis()->SetBinLabel(1, "Before DCA cuts"); + trackSelection->GetXaxis()->SetBinLabel(2, "Finite DCA, |DCA_{xy}| #leq max"); + trackSelection->GetXaxis()->SetBinLabel(3, "min #leq |DCA_{z}| #leq max"); + trackSelection->GetXaxis()->SetBinLabel(8, "Pass DCA selection"); if (processTrackMCEnabled) { qaRegistry.add("QA/hGoodMCTrackIndices", "hGoodMCTrackIndices", kTH1D, {idxAxis}); + qaRegistry.get(HIST("QA/hGoodMCTrackIndices"))->GetXaxis()->SetBinLabel(1, "MC path: before DCA cuts"); } if (cfgDetailTrackQA) { if (!FilterForDerivedTables.cfgBypassTrackFill) { - qaRegistry.add("QA/h4TrackPtEtaPhi", "ResoTracks pT, eta, phi", kTH3F, {ptAxis, etaAxis, phiAxis}); - qaRegistry.add("QA/h2TrackDCAxyVsPt", "ResoTracks DCAxy vs pT", kTH2F, {dcaPtAxis, dcaXYAxis}); - qaRegistry.add("QA/h2TrackDCAzVsPt", "ResoTracks DCAz vs pT", kTH2F, {dcaPtAxis, dcaZAxis}); + qaRegistry.add("QA/h4TrackPtEtaPhi", "ResoTracks pT, eta, phi", kTHnSparseD, {ptAxis, etaAxis, phiAxis}); + qaRegistry.add("QA/h2TrackDCAxyVsPt", "ResoTracks DCAxy vs pT", kTH2D, {dcaPtAxis, dcaXYAxis}); + qaRegistry.add("QA/h2TrackDCAzVsPt", "ResoTracks DCAz vs pT", kTH2D, {dcaPtAxis, dcaZAxis}); qaRegistry.add("QA/h4TrackTPCnSigma", "ResoTracks TPC nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTPCAxis, nSigmaTPCAxis, nSigmaTPCAxis}); qaRegistry.add("QA/h4TrackTOFnSigma", "ResoTracks TOF nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTOFAxis, nSigmaTOFAxis, nSigmaTOFAxis}); } if (FilterForDerivedTables.cfgFillMicroTracks) { - qaRegistry.add("QA/h4MicroTrackPtEtaPhi", "ResoMicroTracks pT, eta, phi", kTH3F, {ptAxis, etaAxis, phiAxis}); - qaRegistry.add("QA/h2MicroTrackDCAxyVsPt", "ResoMicroTracks DCAxy vs pT", kTH2F, {dcaPtAxis, dcaXYAxis}); - qaRegistry.add("QA/h2MicroTrackDCAzVsPt", "ResoMicroTracks DCAz vs pT", kTH2F, {dcaPtAxis, dcaZAxis}); + qaRegistry.add("QA/h4MicroTrackPtEtaPhi", "ResoMicroTracks pT, eta, phi", kTHnSparseD, {ptAxis, etaAxis, phiAxis}); + qaRegistry.add("QA/h2MicroTrackDCAxyVsPt", "ResoMicroTracks DCAxy vs pT", kTH2D, {dcaPtAxis, dcaXYAxis}); + qaRegistry.add("QA/h2MicroTrackDCAzVsPt", "ResoMicroTracks DCAz vs pT", kTH2D, {dcaPtAxis, dcaZAxis}); qaRegistry.add("QA/h4MicroTrackTPCnSigma", "ResoMicroTracks TPC nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTPCAxis, nSigmaTPCAxis, nSigmaTPCAxis}); qaRegistry.add("QA/h4MicroTrackTOFnSigma", "ResoMicroTracks TOF nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTOFAxis, nSigmaTOFAxis, nSigmaTOFAxis}); } if (FilterForDerivedTables.cfgFillUltraMicroTracks) { - qaRegistry.add("QA/h4UltraMicroTrackPtEtaPhi", "ResoUltraMicroTracks pT, eta, phi", kTH3F, {ptAxis, etaAxis, phiAxis}); - qaRegistry.add("QA/h2UltraMicroTrackDCAxyVsPt", "ResoUltraMicroTracks DCAxy vs pT", kTH2F, {dcaPtAxis, dcaXYAxis}); - qaRegistry.add("QA/h2UltraMicroTrackDCAzVsPt", "ResoUltraMicroTracks DCAz vs pT", kTH2F, {dcaPtAxis, dcaZAxis}); + qaRegistry.add("QA/h4UltraMicroTrackPtEtaPhi", "ResoUltraMicroTracks pT, eta, phi", kTHnSparseD, {ptAxis, etaAxis, phiAxis}); + qaRegistry.add("QA/h2UltraMicroTrackDCAxyVsPt", "ResoUltraMicroTracks DCAxy vs pT", kTH2D, {dcaPtAxis, dcaXYAxis}); + qaRegistry.add("QA/h2UltraMicroTrackDCAzVsPt", "ResoUltraMicroTracks DCAz vs pT", kTH2D, {dcaPtAxis, dcaZAxis}); qaRegistry.add("QA/h4UltraMicroTrackTPCnSigma", "ResoUltraMicroTracks TPC nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTPCAxis, nSigmaTPCAxis, nSigmaTPCAxis}); qaRegistry.add("QA/h4UltraMicroTrackTOFnSigma", "ResoUltraMicroTracks TOF nSigma Pi, Ka, Pr as pT", kTHnSparseD, {ptAxis, nSigmaTOFAxis, nSigmaTOFAxis, nSigmaTOFAxis}); } @@ -875,24 +1515,33 @@ struct ResonanceDaughterInitializer { if (processV0DataEnabled || doprocessV0MC) { qaRegistry.add("QA/hGoodV0Indices", "hGoodV0Indices", kTH1D, {idxAxis}); + auto v0Selection = qaRegistry.get(HIST("QA/hGoodV0Indices")); + v0Selection->GetXaxis()->SetBinLabel(1, "Before V0 cuts"); + v0Selection->GetXaxis()->SetBinLabel(2, "Daughter TPC rows"); + v0Selection->GetXaxis()->SetBinLabel(3, "Daughter |DCA_{xy}| to PV"); + v0Selection->GetXaxis()->SetBinLabel(4, "V0 radius window"); + v0Selection->GetXaxis()->SetBinLabel(5, "V0 cosPA cut passed"); if (doprocessV0MC) { qaRegistry.add("QA/hGoodMCV0Indices", "hGoodMCV0Indices", kTH1D, {idxAxis}); + qaRegistry.get(HIST("QA/hGoodMCV0Indices"))->GetXaxis()->SetBinLabel(1, "Reco V0 cuts passed (MC)"); } - AxisSpec radiusAxis = {100, 0.0, 200.0, "V0 Radius"}; - AxisSpec cosPAAxis = {100, 0.995, 1.0, "V0 CosPA"}; - qaRegistry.add("QA/hV0Radius", "V0 Radius", kTH1F, {radiusAxis}); - qaRegistry.add("QA/hV0CosPA", "V0 CosPA", kTH1F, {cosPAAxis}); } if (processCascDataEnabled || doprocessCascMC) { - AxisSpec radiusAxis = {100, 0.0, 200.0, "Cascade Radius"}; - AxisSpec cosPAAxis = {100, 0.97, 1.0, "Cascade CosPA"}; qaRegistry.add("QA/hGoodCascIndices", "hGoodCascIndices", kTH1D, {idxAxis}); + auto cascadeSelection = qaRegistry.get(HIST("QA/hGoodCascIndices")); + cascadeSelection->GetXaxis()->SetBinLabel(1, "Before cascade cuts"); + cascadeSelection->GetXaxis()->SetBinLabel(2, "Bachelor TPC rows"); + cascadeSelection->GetXaxis()->SetBinLabel(3, "Bachelor DCA_{xy} window"); + cascadeSelection->GetXaxis()->SetBinLabel(4, "V0/cascade daughter DCA"); + cascadeSelection->GetXaxis()->SetBinLabel(5, "Cascade/V0 cosPA"); + cascadeSelection->GetXaxis()->SetBinLabel(6, "V0 radius window"); + cascadeSelection->GetXaxis()->SetBinLabel(7, "Cascade radius window"); + cascadeSelection->GetXaxis()->SetBinLabel(8, "#Xi mass window"); if (doprocessCascMC) { qaRegistry.add("QA/hGoodMCCascIndices", "hGoodMCCascIndices", kTH1D, {idxAxis}); + qaRegistry.get(HIST("QA/hGoodMCCascIndices"))->GetXaxis()->SetBinLabel(1, "Reco cascade cuts passed (MC)"); } - qaRegistry.add("QA/hCascRadius", "Cascade Radius", kTH1F, {radiusAxis}); - qaRegistry.add("QA/hCascCosPA", "Cascade CosPA", kTH1F, {cosPAAxis}); } } if (processTrackDataEnabled || processTrackMCEnabled) { @@ -905,52 +1554,29 @@ struct ResonanceDaughterInitializer { LOGF(info, "ResonanceDaughterInitializer initialized with cascades"); } - // Check if the module is initialized with both data and MC - if ((processTrackDataEnabled && processTrackMCEnabled) || (processV0DataEnabled && doprocessV0MC) || (processCascDataEnabled && doprocessCascMC)) { - LOGF(fatal, "ResonanceDaughterInitializer initialized with both data and MC"); - } // Check if none of the processes are enabled - if (!doprocessDummy && !processTrackDataEnabled && !processTrackMCEnabled && !processV0DataEnabled && !doprocessV0MC && !processCascDataEnabled && !doprocessCascMC) { + if (!doprocessDummy && !anyDataProcessEnabled && !anyMCProcessEnabled) { LOGF(fatal, "ResonanceDaughterInitializer not initialized, enable at least one process"); } } - template - bool filterMicroTrack(T const& track) + static bool passesCenteredPID(float nSigma, float mean, float cut) { - // if no selection is requested, return true - if (!FilterForDerivedTables.cfgFillPionMicroTracks && !FilterForDerivedTables.cfgFillKaonMicroTracks && !FilterForDerivedTables.cfgFillProtonMicroTracks) { - return true; - } - if (FilterForDerivedTables.cfgFillPionMicroTracks) { - if (std::abs(track.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } - } - if (FilterForDerivedTables.cfgFillKaonMicroTracks) { - if (std::abs(track.tpcNSigmaKa()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } - } - if (FilterForDerivedTables.cfgFillProtonMicroTracks) { - if (std::abs(track.tpcNSigmaPr()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } - } - return false; + return std::isfinite(nSigma) && std::abs(nSigma - mean) < cut; } - template - bool filterUltraMicroTrack(T const& track) + bool passesPIDPreSelection(float tpcNSigma, + float tofNSigma, + bool hasTOF, + float tpcMean, + float tofMean, + float tpcCut, + float tofCut) const { - switch (ultraMicroPidSpecies) { - case UltraMicroPidSpecies::Pion: - return std::abs(track.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut; - case UltraMicroPidSpecies::Kaon: - return std::abs(track.tpcNSigmaKa()) < TrackCuts.pidnSigmaPreSelectionCut; - case UltraMicroPidSpecies::Proton: - return std::abs(track.tpcNSigmaPr()) < TrackCuts.pidnSigmaPreSelectionCut; + if (!passesCenteredPID(tpcNSigma, tpcMean, tpcCut)) { + return false; } - return false; + return !TrackCuts.cfgUseTOFPIDPreSelection.value || + !hasTOF || passesCenteredPID(tofNSigma, tofMean, tofCut); } template @@ -960,20 +1586,25 @@ struct ResonanceDaughterInitializer { if (!FilterForDerivedTables.cfgFillPionTracks && !FilterForDerivedTables.cfgFillKaonTracks && !FilterForDerivedTables.cfgFillProtonTracks) { return true; } - if (FilterForDerivedTables.cfgFillPionTracks) { - if (std::abs(track.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } + const float tpcCut = TrackCuts.pidnSigmaPreSelectionCut.value; + const float tofCut = TrackCuts.pidnSigmaPreSelectionCutTOF.value; + if (FilterForDerivedTables.cfgFillPionTracks && + passesPIDPreSelection(track.tpcNSigmaPi(), track.tofNSigmaPi(), track.hasTOF(), + TrackCuts.pidnSigmaPreSelectionMeanPion.value, + TrackCuts.pidnSigmaPreSelectionMeanTOFPion.value, tpcCut, tofCut)) { + return true; } - if (FilterForDerivedTables.cfgFillKaonTracks) { - if (std::abs(track.tpcNSigmaKa()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } + if (FilterForDerivedTables.cfgFillKaonTracks && + passesPIDPreSelection(track.tpcNSigmaKa(), track.tofNSigmaKa(), track.hasTOF(), + TrackCuts.pidnSigmaPreSelectionMeanKaon.value, + TrackCuts.pidnSigmaPreSelectionMeanTOFKaon.value, tpcCut, tofCut)) { + return true; } - if (FilterForDerivedTables.cfgFillProtonTracks) { - if (std::abs(track.tpcNSigmaPr()) < TrackCuts.pidnSigmaPreSelectionCut) { - return true; - } + if (FilterForDerivedTables.cfgFillProtonTracks && + passesPIDPreSelection(track.tpcNSigmaPr(), track.tofNSigmaPr(), track.hasTOF(), + TrackCuts.pidnSigmaPreSelectionMeanProton.value, + TrackCuts.pidnSigmaPreSelectionMeanTOFProton.value, tpcCut, tofCut)) { + return true; } return false; } @@ -988,13 +1619,8 @@ struct ResonanceDaughterInitializer { return true; } if (v0.dcaV0daughters() > SecondaryCuts.cfgSecondaryDauDCAMax || - std::abs(v0.dcapostopv()) < SecondaryCuts.cfgSecondaryDauPosDCAtoPVMin || - std::abs(v0.dcanegtopv()) < SecondaryCuts.cfgSecondaryDauNegDCAtoPVMin || v0.pt() < SecondaryCuts.cfgSecondaryPtMin || - v0.v0radius() < SecondaryCuts.cfgSecondaryRadiusMin || - v0.v0radius() > SecondaryCuts.cfgSecondaryRadiusMax || - v0.dcav0topv() > SecondaryCuts.cfgSecondaryDCAtoPVMax || - v0.v0cosPA() < SecondaryCuts.cfgSecondaryCosPAMin) { + v0.dcav0topv() > SecondaryCuts.cfgSecondaryDCAtoPVMax) { return false; } if (SecondaryCuts.cfgSecondaryArmenterosCut && @@ -1006,11 +1632,12 @@ struct ResonanceDaughterInitializer { const auto posTrack = v0.template posTrack_as(); const auto negTrack = v0.template negTrack_as(); const bool bypassDaughterPID = SecondaryCuts.cfgByPassDauPIDSelection; + const float daughterPIDCut = TrackCuts.pidnSigmaPreSelectionCut.value; bool selected = false; if (FilterForDerivedTables.cfgFillK0s) { const bool passesK0DaughterPID = bypassDaughterPID || - (std::abs(posTrack.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut && - std::abs(negTrack.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut); + (passesCenteredPID(posTrack.tpcNSigmaPi(), TrackCuts.pidnSigmaPreSelectionMeanPion.value, daughterPIDCut) && + passesCenteredPID(negTrack.tpcNSigmaPi(), TrackCuts.pidnSigmaPreSelectionMeanPion.value, daughterPIDCut)); const bool passesK0 = std::fabs(v0.yK0Short()) <= SecondaryCuts.cfgSecondaryRapidityMax && decayLengthOverMomentum * MassK0Short <= SecondaryCuts.cfgSecondaryProperLifetimeMax && std::fabs(v0.mK0Short() - MassK0Short) <= SecondaryCuts.cfgSecondaryMassWindow && @@ -1022,11 +1649,11 @@ struct ResonanceDaughterInitializer { } if (FilterForDerivedTables.cfgFillLambda0) { const bool passesLambdaPID = bypassDaughterPID || - (std::abs(posTrack.tpcNSigmaPr()) < TrackCuts.pidnSigmaPreSelectionCut && - std::abs(negTrack.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut); + (passesCenteredPID(posTrack.tpcNSigmaPr(), TrackCuts.pidnSigmaPreSelectionMeanProton.value, daughterPIDCut) && + passesCenteredPID(negTrack.tpcNSigmaPi(), TrackCuts.pidnSigmaPreSelectionMeanPion.value, daughterPIDCut)); const bool passesAntiLambdaPID = bypassDaughterPID || - (std::abs(posTrack.tpcNSigmaPi()) < TrackCuts.pidnSigmaPreSelectionCut && - std::abs(negTrack.tpcNSigmaPr()) < TrackCuts.pidnSigmaPreSelectionCut); + (passesCenteredPID(posTrack.tpcNSigmaPi(), TrackCuts.pidnSigmaPreSelectionMeanPion.value, daughterPIDCut) && + passesCenteredPID(negTrack.tpcNSigmaPr(), TrackCuts.pidnSigmaPreSelectionMeanProton.value, daughterPIDCut)); const bool passesLambdaMassAndPID = (std::fabs(v0.mLambda() - MassLambda0) <= SecondaryCuts.cfgSecondaryMassWindow && passesLambdaPID) || (std::fabs(v0.mAntiLambda() - MassLambda0Bar) <= SecondaryCuts.cfgSecondaryMassWindow && passesAntiLambdaPID); @@ -1057,9 +1684,9 @@ struct ResonanceDaughterInitializer { } template - bool isTrackSelected(CollisionType const&, TrackType const& track) + bool isTrackSelected(CollisionType const&, TrackType const& track, bool fillSelectionQA = true) { - if (cfgFillQA) { + if (cfgFillQA && fillSelectionQA) { qaRegistry.fill(HIST("QA/hGoodTrackIndices"), 0.5); if constexpr (isMC) { qaRegistry.fill(HIST("QA/hGoodMCTrackIndices"), 0.5); @@ -1071,13 +1698,13 @@ struct ResonanceDaughterInitializer { if (std::fabs(track.dcaXY()) > TrackCuts.cMaxDCArToPVcut) { return false; } - if (cfgFillQA) { + if (cfgFillQA && fillSelectionQA) { qaRegistry.fill(HIST("QA/hGoodTrackIndices"), 1.5); } if (std::fabs(track.dcaZ()) > TrackCuts.cMaxDCAzToPVcut || std::fabs(track.dcaZ()) < TrackCuts.cMinDCAzToPVcut) { return false; } - if (cfgFillQA) { + if (cfgFillQA && fillSelectionQA) { qaRegistry.fill(HIST("QA/hGoodTrackIndices"), 2.5); qaRegistry.fill(HIST("QA/hGoodTrackIndices"), 7.5); } @@ -1093,7 +1720,7 @@ struct ResonanceDaughterInitializer { auto posTrack = v0.template posTrack_as(); auto negTrack = v0.template negTrack_as(); - if (posTrack.tpcNClsCrossedRows() < TrackCuts.mincrossedrows || negTrack.tpcNClsCrossedRows() < TrackCuts.mincrossedrows) { + if (posTrack.tpcNClsCrossedRows() < V0Cuts.mincrossedrowsV0s || negTrack.tpcNClsCrossedRows() < V0Cuts.mincrossedrowsV0s) { return false; } if (cfgFillQA && fillSelectionQA) { @@ -1173,7 +1800,7 @@ struct ResonanceDaughterInitializer { if (cfgFillQA && fillSelectionQA) { qaRegistry.fill(HIST("QA/hGoodCascIndices"), 6.5); } - if (std::abs(casc.mXi() - MassXiMinus) > CascadeCuts.cCascMassResol) { + if (std::abs(casc.mXi() - MassXiMinus) > CascadeCuts.cMaxXiMassWindow) { return false; } if (cfgFillQA && fillSelectionQA) { @@ -1185,28 +1812,122 @@ struct ResonanceDaughterInitializer { return true; } - /// @brief Check whether a collision has at least one V0 that would be written + /// @brief Find a selected V0 and collect daughter IDs only for the optional global veto template - bool hasSelectedV0(CollisionType const& collision, V0Type const& v0s, TrackType const& tracks) + SelectedCandidateDaughters collectSelectedV0Daughters(CollisionType const& collision, + V0Type const& v0s, + TrackType const& tracks, + bool useGlobalDaughterVeto) { + SelectedCandidateDaughters selectedCandidates{useGlobalDaughterVeto}; for (auto const& v0 : v0s) { - if (isV0Selected(collision, v0, tracks, false) && filterV0(collision, v0, tracks)) { - return true; + if (!isV0Selected(collision, v0, tracks, false) || + !filterV0(collision, v0, tracks)) { + continue; + } + if (!useGlobalDaughterVeto) { + selectedCandidates.addCandidate(); + break; } + selectedCandidates.addCandidate( + std::array{static_cast(v0.posTrackId()), + static_cast(v0.negTrackId())}); } - return false; + selectedCandidates.finalize(); + return selectedCandidates; } - /// @brief Check whether a collision has at least one cascade that would be written + /// @brief Find a selected cascade and collect daughter IDs only for the optional global veto template - bool hasSelectedCascade(CollisionType const& collision, CascType const& cascades, TrackType const& tracks) + SelectedCandidateDaughters collectSelectedCascadeDaughters(CollisionType const& collision, + CascType const& cascades, + TrackType const& tracks, + bool useGlobalDaughterVeto) { + SelectedCandidateDaughters selectedCandidates{useGlobalDaughterVeto}; for (auto const& casc : cascades) { - if (isCascSelected(collision, casc, tracks, false)) { - return true; + if (!isCascSelected(collision, casc, tracks, false)) { + continue; + } + if (!useGlobalDaughterVeto) { + selectedCandidates.addCandidate(); + break; } + selectedCandidates.addCandidate( + std::array{static_cast(casc.posTrackId()), + static_cast(casc.negTrackId()), + static_cast(casc.bachelorId())}); } - return false; + selectedCandidates.finalize(); + return selectedCandidates; + } + + /// @brief Build the per-track selection while preserving the configured collision gate + template + bool preparePairTrackSelection(CollisionType const& collision, + TrackType const& tracks, + V0Type const& v0s, + CascType const& cascades, + V0PresliceType const& v0Preslice, + CascPresliceType const& cascPreslice, + PairTrackSelection& selection) + { + const bool useEitherPairGate = FilterForDerivedTables.cfgPairGateMode.value == PairGateModeEither; + selection.useV0Candidates = useEitherPairGate || FilterForDerivedTables.cfgBypassNoPairV0s; + selection.useCascadeCandidates = useEitherPairGate || FilterForDerivedTables.cfgBypassNoPairCascades; + + if (selection.useV0Candidates) { + auto v0sThisCollision = v0s.sliceBy(v0Preslice, collision.collisionId()); + selection.v0Candidates = + collectSelectedV0Daughters(collision, v0sThisCollision, tracks, selection.useGlobalDaughterVeto); + } + if (selection.useCascadeCandidates && + (!useEitherPairGate || selection.useGlobalDaughterVeto || !selection.v0Candidates.hasSelectedCandidate)) { + auto cascadesThisCollision = cascades.sliceBy(cascPreslice, collision.collisionId()); + selection.cascadeCandidates = + collectSelectedCascadeDaughters(collision, cascadesThisCollision, tracks, selection.useGlobalDaughterVeto); + } + + if (useEitherPairGate) { + return selection.v0Candidates.hasSelectedCandidate || + selection.cascadeCandidates.hasSelectedCandidate; + } + + if (FilterForDerivedTables.cfgBypassNoPairV0s && !selection.v0Candidates.hasSelectedCandidate) { + return false; + } + if (FilterForDerivedTables.cfgBypassNoPairCascades && !selection.cascadeCandidates.hasSelectedCandidate) { + return false; + } + return true; + } + + /// @brief Round and saturate a floating-point value into a persistent integer column + template + static Integer quantizeSaturated(float value, double scale) + { + constexpr Integer LowerLimit = std::numeric_limits::lowest(); + constexpr Integer UpperLimit = std::numeric_limits::max(); + if (std::isnan(value)) { + return UpperLimit; + } + if (!std::isfinite(value)) { + return std::signbit(value) ? LowerLimit : UpperLimit; + } + const double rounded = std::round(static_cast(value) * scale); + if (rounded <= static_cast(LowerLimit)) { + return LowerLimit; + } + if (rounded >= static_cast(UpperLimit)) { + return UpperLimit; + } + return static_cast(rounded); } static bool quantizeP(float p, int16_t& quantized) @@ -1237,28 +1958,137 @@ struct ResonanceDaughterInitializer { return std::isfinite(threshold) ? threshold : -1.f; } - template - void fillUltraMicroTracks(CollisionType const& collision, TrackType const& tracks) + template + bool evaluatePtDependentDCA(TrackType const& track, + bool& passedPtDependentDCAxy, + bool& passedPtDependentDCAz) const { - // Loop over tracks + passedPtDependentDCAxy = false; + passedPtDependentDCAz = false; + if (!TrackCuts.cfgApplyTightDCAPtDepSelection.value) { + return true; + } + const float dcaThreshold = tightDCAThreshold(track.pt()); + if (dcaThreshold < 0.f) { + return false; + } + passedPtDependentDCAxy = std::isfinite(track.dcaXY()) && std::abs(track.dcaXY()) < dcaThreshold; + passedPtDependentDCAz = std::isfinite(track.dcaZ()) && std::abs(track.dcaZ()) < dcaThreshold; + return passedPtDependentDCAxy && passedPtDependentDCAz; + } + + template + bool isFullTrackOutputSelected(CollisionType const& collision, + TrackType const& track, + bool fillSelectionQA = true) + { + return isTrackSelected(collision, track, fillSelectionQA) && filterTrack(track); + } + + template + bool isMicroTrackOutputSelected(CollisionType const& collision, + TrackType const& track, + bool& passedPtDependentDCAxy, + bool& passedPtDependentDCAz) + { + return isMicroTrackSelected(collision, track) && + filterTrack(track) && + evaluatePtDependentDCA(track, passedPtDependentDCAxy, passedPtDependentDCAz); + } + + template + bool isUltraMicroTrackOutputSelected(CollisionType const& collision, TrackType const& track) + { + return isMicroTrackSelected(collision, track) && + filterTrack(track) && + o2::aod::resoultramicrodaughter::DCAEncoding::isValid(track.dcaXY()) && + o2::aod::resoultramicrodaughter::DCAEncoding::isValid(track.dcaZ()); + } + + template + static bool quantizeUltraMicroMomentum(TrackType const& track, + int16_t& px1000, + int16_t& py1000, + int16_t& pz1000) + { + return quantizeP(track.px(), px1000) && + quantizeP(track.py(), py1000) && + quantizeP(track.pz(), pz1000); + } + + /// Check that every enabled track output will receive at least one row. + template + bool hasTracksForEnabledOutputs(CollisionType const& collision, + TrackTableType const& tracks, + TrackPredicate const& keepTrack) + { + constexpr uint8_t FullTrackOutput = 1u << 0; + constexpr uint8_t MicroTrackOutput = 1u << 1; + constexpr uint8_t UltraMicroTrackOutput = 1u << 2; + uint8_t requiredOutputs = 0; + if (!FilterForDerivedTables.cfgBypassTrackFill.value) { + requiredOutputs |= FullTrackOutput; + } + if (FilterForDerivedTables.cfgFillMicroTracks.value) { + requiredOutputs |= MicroTrackOutput; + } + if (FilterForDerivedTables.cfgFillUltraMicroTracks.value) { + requiredOutputs |= UltraMicroTrackOutput; + } + + uint8_t availableOutputs = 0; for (auto const& track : tracks) { - if (!isMicroTrackSelected(collision, track)) { + if (!keepTrack(track)) { continue; } - if (!filterUltraMicroTrack(track)) { + if ((requiredOutputs & FullTrackOutput) != 0 && + (availableOutputs & FullTrackOutput) == 0 && + isFullTrackOutputSelected(collision, track, false)) { + availableOutputs |= FullTrackOutput; + } + if ((requiredOutputs & MicroTrackOutput) != 0 && + (availableOutputs & MicroTrackOutput) == 0) { + bool passedPtDependentDCAxy = false; + bool passedPtDependentDCAz = false; + if (isMicroTrackOutputSelected(collision, track, + passedPtDependentDCAxy, + passedPtDependentDCAz)) { + availableOutputs |= MicroTrackOutput; + } + } + if ((requiredOutputs & UltraMicroTrackOutput) != 0 && + (availableOutputs & UltraMicroTrackOutput) == 0 && + isUltraMicroTrackOutputSelected(collision, track)) { + int16_t px1000 = 0; + int16_t py1000 = 0; + int16_t pz1000 = 0; + if (quantizeUltraMicroMomentum(track, px1000, py1000, pz1000)) { + availableOutputs |= UltraMicroTrackOutput; + } + } + if (availableOutputs == requiredOutputs) { + return true; + } + } + return false; + } + + template + void fillUltraMicroTracks(CollisionType const& collision, TrackType const& tracks, TrackPredicate const& keepTrack) + { + // Loop over tracks + for (auto const& track : tracks) { + if (!keepTrack(track)) { continue; } - if (!o2::aod::resoultramicrodaughter::DCAEncoding::isValid(track.dcaXY()) || - !o2::aod::resoultramicrodaughter::DCAEncoding::isValid(track.dcaZ())) { + if (!isUltraMicroTrackOutputSelected(collision, track)) { continue; } o2::aod::resoultramicrodaughter::DCAEncoding dcaEncoding(track.dcaXY(), track.dcaZ()); int16_t px1000 = 0; int16_t py1000 = 0; int16_t pz1000 = 0; - if (!quantizeP(track.px(), px1000) || - !quantizeP(track.py(), py1000) || - !quantizeP(track.pz(), pz1000)) { + if (!quantizeUltraMicroMomentum(track, px1000, py1000, pz1000)) { if (!warnedUltraMicroMomentumRange) { LOGF(warn, "Skipping ultra-micro tracks with non-finite or out-of-range momentum components"); warnedUltraMicroMomentumRange = true; @@ -1318,28 +2148,36 @@ struct ResonanceDaughterInitializer { * @tparam CollisionType Type of collision * @param collision Collision data * @param tracks Track data + * @note ResoMicroTracks_001 is intended for fixed producer-side selections. + * Exact downstream cuts on its decoded values require zero-centred + * strict PID windows |nSigma| < C with + * C in {2.0, 2.25, 2.5, 2.75, 3.0, 3.25, 3.5}, and strict DCA windows + * |DCA| < C with C = N * 0.025 cm (N = 1, ..., 6). The listed PID + * guarantee assumes a zero mean; shifted PID windows require both + * interval edges to align with encoding boundaries and separate + * validation. PID cuts below 2 sigma and off-grid DCA/PID cuts cannot + * be reconstructed exactly. In particular, pT-dependent DCA cuts must + * be applied to the unquantised values here via + * cfgApplyTightDCAPtDepSelection; consumers should retain that producer + * decision rather than retune the decoded lower-edge DCA values. */ - template - void fillMicroTracks(CollisionType const& collision, TrackType const& tracks) + template + void fillMicroTracks(CollisionType const& collision, TrackType const& tracks, TrackPredicate const& keepTrack) { // Loop over tracks for (auto const& track : tracks) { - if (!isMicroTrackSelected(collision, track)) { + if (!keepTrack(track)) { continue; } - if (!filterMicroTrack(track)) { + bool passedPtDependentDCAxy = false; + bool passedPtDependentDCAz = false; + if (!isMicroTrackOutputSelected(collision, track, + passedPtDependentDCAxy, + passedPtDependentDCAz)) { continue; } - o2::aod::resomicrodaughter::ResoMicroTrackSelFlag trackSelFlag(track.dcaXY(), track.dcaZ()); - if (TrackCuts.cfgApplyTightDCAPtDepSelection) { - const float dcaThreshold = tightDCAThreshold(track.pt()); - if (dcaThreshold >= 0.f && std::abs(track.dcaXY()) < dcaThreshold) { - trackSelFlag.setDCAxy0(); - } - if (dcaThreshold >= 0.f && std::abs(track.dcaZ()) < dcaThreshold) { - trackSelFlag.setDCAz0(); - } - } + const o2::aod::resomicrodaughter001::DCAEncoding trackSelFlag( + track.dcaXY(), track.dcaZ(), passedPtDependentDCAxy, passedPtDependentDCAz); uint8_t trackFlags = (track.passedITSRefit() << 0) | (track.passedTPCRefit() << 1) | (track.isGlobalTrackWoDCA() << 2) | @@ -1358,17 +2196,21 @@ struct ResonanceDaughterInitializer { } } reso2microtrks(collision.globalIndex(), + track.globalIndex(), track.px(), track.py(), track.pz(), - static_cast(o2::aod::resomicrodaughter::PidNSigma(std::abs(track.tpcNSigmaPi()), std::abs(track.tofNSigmaPi()), track.hasTOF())), - static_cast(o2::aod::resomicrodaughter::PidNSigma(std::abs(track.tpcNSigmaKa()), std::abs(track.tofNSigmaKa()), track.hasTOF())), - static_cast(o2::aod::resomicrodaughter::PidNSigma(std::abs(track.tpcNSigmaPr()), std::abs(track.tofNSigmaPr()), track.hasTOF())), + static_cast(o2::aod::resomicrodaughter001::PidNSigma(track.tpcNSigmaPi(), track.tofNSigmaPi(), track.hasTOF())), + static_cast(o2::aod::resomicrodaughter001::PidNSigma(track.tpcNSigmaKa(), track.tofNSigmaKa(), track.hasTOF())), + static_cast(o2::aod::resomicrodaughter001::PidNSigma(track.tpcNSigmaPr(), track.tofNSigmaPr(), track.hasTOF())), static_cast(trackSelFlag), trackFlags); if (!FilterForDerivedTables.cfgBypassTrackIndexFill) { resoMicroTrackTracks(track.globalIndex()); } + if constexpr (isMC) { + fillMCTrack(track, reso2mcmicrotrks); + } } } @@ -1381,18 +2223,18 @@ struct ResonanceDaughterInitializer { * @param collision Collision data * @param tracks Track data */ - template - void fillTracks(CollisionType const& collision, TrackType const& tracks) + template + void fillTracks(CollisionType const& collision, TrackType const& tracks, TrackPredicate const& keepTrack) { if (FilterForDerivedTables.cfgBypassTrackFill) { return; } // Loop over tracks for (auto const& track : tracks) { - if (!isTrackSelected(collision, track)) { + if (!keepTrack(track)) { continue; } - if (!filterTrack(track)) { + if (!isFullTrackOutputSelected(collision, track)) { continue; } uint8_t trackFlags = (track.passedITSRefit() << 0) | @@ -1419,21 +2261,21 @@ struct ResonanceDaughterInitializer { track.pz(), static_cast(track.tpcNClsCrossedRows()), static_cast(track.tpcNClsFound()), - static_cast(std::round(track.dcaXY() * 10000)), - static_cast(std::round(track.dcaZ() * 10000)), - static_cast(std::round(track.tpcNSigmaPi() * 10)), - static_cast(std::round(track.tpcNSigmaKa() * 10)), - static_cast(std::round(track.tpcNSigmaPr() * 10)), - static_cast(std::round(track.tofNSigmaPi() * 10)), - static_cast(std::round(track.tofNSigmaKa() * 10)), - static_cast(std::round(track.tofNSigmaPr() * 10)), - static_cast(std::round(track.tpcSignal() * 100)), + quantizeSaturated(track.dcaXY(), 10000.), + quantizeSaturated(track.dcaZ(), 10000.), + quantizeSaturated(track.tpcNSigmaPi(), 10.), + quantizeSaturated(track.tpcNSigmaKa(), 10.), + quantizeSaturated(track.tpcNSigmaPr(), 10.), + quantizeSaturated(track.tofNSigmaPi(), 10.), + quantizeSaturated(track.tofNSigmaKa(), 10.), + quantizeSaturated(track.tofNSigmaPr(), 10.), + quantizeSaturated(track.tpcSignal(), 100.), trackFlags); if (!FilterForDerivedTables.cfgBypassTrackIndexFill) { resoTrackTracks(track.globalIndex()); } if constexpr (isMC) { - fillMCTrack(track); + fillMCTrack(track, reso2mctracks); } } } @@ -1442,10 +2284,12 @@ struct ResonanceDaughterInitializer { * @brief Fills MC track data * * @tparam TrackType Type of track + * @tparam MCOutputTable Type of positional MC extension to fill * @param track Track data + * @param output Positional MC extension writer */ - template - void fillMCTrack(TrackType const& track) + template + void fillMCTrack(TrackType const& track, Produces& output) { // ------ Temporal lambda function to prevent error in build auto getMothersIndeces = [&](auto const& theMcParticle) { @@ -1498,20 +2342,20 @@ struct ResonanceDaughterInitializer { if (siblingsTemp.size() > 1) { siblings[1] = siblingsTemp[1]; } - reso2mctracks(particle.pdgCode(), - mothers[0], - motherPDGs[0], - siblings.data(), - particle.isPhysicalPrimary(), - particle.producedByGenerator()); + output(particle.pdgCode(), + mothers[0], + motherPDGs[0], + siblings.data(), + particle.isPhysicalPrimary(), + particle.producedByGenerator()); } else { // No MC particle associated - reso2mctracks(0, - mothers[0], - motherPDGs[0], - siblings.data(), - 0, - 0); + output(0, + mothers[0], + motherPDGs[0], + siblings.data(), + 0, + 0); } } @@ -1536,10 +2380,8 @@ struct ResonanceDaughterInitializer { if (!filterV0(collision, v0, tracks)) { continue; } - if (cfgFillQA) { - qaRegistry.fill(HIST("QA/hV0Radius"), v0.v0radius()); - qaRegistry.fill(HIST("QA/hV0CosPA"), v0.v0cosPA()); - } + const auto posTrack = v0.template posTrack_as(); + const auto negTrack = v0.template negTrack_as(); const std::array childIDs{v0.posTrackId(), v0.negTrackId()}; // Original track IDs for downstream pair-level shared-daughter rejection reso2v0s(collision.globalIndex(), v0.pt(), @@ -1547,25 +2389,25 @@ struct ResonanceDaughterInitializer { v0.py(), v0.pz(), childIDs.data(), - (int8_t)(v0.template posTrack_as().tpcNSigmaPi() * 10), - (int8_t)(v0.template posTrack_as().tpcNSigmaKa() * 10), - (int8_t)(v0.template posTrack_as().tpcNSigmaPr() * 10), - (int8_t)(v0.template negTrack_as().tpcNSigmaPi() * 10), - (int8_t)(v0.template negTrack_as().tpcNSigmaKa() * 10), - (int8_t)(v0.template negTrack_as().tpcNSigmaPr() * 10), - (int8_t)(v0.template posTrack_as().tofNSigmaPi() * 10), - (int8_t)(v0.template posTrack_as().tofNSigmaKa() * 10), - (int8_t)(v0.template posTrack_as().tofNSigmaPr() * 10), - (int8_t)(v0.template negTrack_as().tofNSigmaPi() * 10), - (int8_t)(v0.template negTrack_as().tofNSigmaKa() * 10), - (int8_t)(v0.template negTrack_as().tofNSigmaPr() * 10), + quantizeSaturated(posTrack.tpcNSigmaPi(), 10.), + quantizeSaturated(posTrack.tpcNSigmaKa(), 10.), + quantizeSaturated(posTrack.tpcNSigmaPr(), 10.), + quantizeSaturated(negTrack.tpcNSigmaPi(), 10.), + quantizeSaturated(negTrack.tpcNSigmaKa(), 10.), + quantizeSaturated(negTrack.tpcNSigmaPr(), 10.), + quantizeSaturated(posTrack.tofNSigmaPi(), 10.), + quantizeSaturated(posTrack.tofNSigmaKa(), 10.), + quantizeSaturated(posTrack.tofNSigmaPr(), 10.), + quantizeSaturated(negTrack.tofNSigmaPi(), 10.), + quantizeSaturated(negTrack.tofNSigmaKa(), 10.), + quantizeSaturated(negTrack.tofNSigmaPr(), 10.), v0.v0cosPA(), v0.dcaV0daughters(), v0.dcapostopv(), v0.dcanegtopv(), v0.dcav0topv(), - static_cast(v0.template posTrack_as().tpcNClsCrossedRows()), - static_cast(v0.template negTrack_as().tpcNClsCrossedRows()), + static_cast(posTrack.tpcNClsCrossedRows()), + static_cast(negTrack.tpcNClsCrossedRows()), v0.mLambda(), v0.mAntiLambda(), v0.mK0Short(), @@ -1661,9 +2503,9 @@ struct ResonanceDaughterInitializer { daughters = getDaughtersIndeces(v0mc); daughterPDGs = getDaughtersPDGCodes(v0mc); } - if (daughters.size() > StoredMCRelationCount) { - LOGF(info, "daughters.size() is larger than 2"); - } + // if (daughters.size() > StoredMCRelationCount) { + // LOGF(info, "daughters.size() is larger than 2"); + // } daughters.resize(StoredMCRelationCount, -1); daughterPDGs.resize(StoredMCRelationCount, -1); reso2mcv0s(v0mc.pdgCode(), @@ -1710,10 +2552,9 @@ struct ResonanceDaughterInitializer { if (!isCascSelected(collision, casc, tracks)) { continue; } - if (cfgFillQA) { - qaRegistry.fill(HIST("QA/hCascRadius"), casc.cascradius()); - qaRegistry.fill(HIST("QA/hCascCosPA"), casc.casccosPA(collision.posX(), collision.posY(), collision.posZ())); - } + const auto posTrack = casc.template posTrack_as(); + const auto negTrack = casc.template negTrack_as(); + const auto bachelor = casc.template bachelor_as(); const std::array childIDs{casc.posTrackId(), casc.negTrackId(), casc.bachelorId()}; // Original track IDs for downstream pair-level shared-daughter rejection reso2cascades(collision.globalIndex(), casc.pt(), @@ -1721,24 +2562,24 @@ struct ResonanceDaughterInitializer { casc.py(), casc.pz(), childIDs.data(), - (int8_t)(casc.template posTrack_as().tpcNSigmaPi() * 10), - (int8_t)(casc.template posTrack_as().tpcNSigmaKa() * 10), - (int8_t)(casc.template posTrack_as().tpcNSigmaPr() * 10), - (int8_t)(casc.template negTrack_as().tpcNSigmaPi() * 10), - (int8_t)(casc.template negTrack_as().tpcNSigmaKa() * 10), - (int8_t)(casc.template negTrack_as().tpcNSigmaPr() * 10), - (int8_t)(casc.template bachelor_as().tpcNSigmaPi() * 10), - (int8_t)(casc.template bachelor_as().tpcNSigmaKa() * 10), - (int8_t)(casc.template bachelor_as().tpcNSigmaPr() * 10), - (int8_t)(casc.template posTrack_as().tofNSigmaPi() * 10), - (int8_t)(casc.template posTrack_as().tofNSigmaKa() * 10), - (int8_t)(casc.template posTrack_as().tofNSigmaPr() * 10), - (int8_t)(casc.template negTrack_as().tofNSigmaPi() * 10), - (int8_t)(casc.template negTrack_as().tofNSigmaKa() * 10), - (int8_t)(casc.template negTrack_as().tofNSigmaPr() * 10), - (int8_t)(casc.template bachelor_as().tofNSigmaPi() * 10), - (int8_t)(casc.template bachelor_as().tofNSigmaKa() * 10), - (int8_t)(casc.template bachelor_as().tofNSigmaPr() * 10), + quantizeSaturated(posTrack.tpcNSigmaPi(), 10.), + quantizeSaturated(posTrack.tpcNSigmaKa(), 10.), + quantizeSaturated(posTrack.tpcNSigmaPr(), 10.), + quantizeSaturated(negTrack.tpcNSigmaPi(), 10.), + quantizeSaturated(negTrack.tpcNSigmaKa(), 10.), + quantizeSaturated(negTrack.tpcNSigmaPr(), 10.), + quantizeSaturated(bachelor.tpcNSigmaPi(), 10.), + quantizeSaturated(bachelor.tpcNSigmaKa(), 10.), + quantizeSaturated(bachelor.tpcNSigmaPr(), 10.), + quantizeSaturated(posTrack.tofNSigmaPi(), 10.), + quantizeSaturated(posTrack.tofNSigmaKa(), 10.), + quantizeSaturated(posTrack.tofNSigmaPr(), 10.), + quantizeSaturated(negTrack.tofNSigmaPi(), 10.), + quantizeSaturated(negTrack.tofNSigmaKa(), 10.), + quantizeSaturated(negTrack.tofNSigmaPr(), 10.), + quantizeSaturated(bachelor.tofNSigmaPi(), 10.), + quantizeSaturated(bachelor.tofNSigmaKa(), 10.), + quantizeSaturated(bachelor.tofNSigmaPr(), 10.), casc.v0cosPA(collision.posX(), collision.posY(), collision.posZ()), casc.casccosPA(collision.posX(), collision.posY(), collision.posZ()), casc.dcaV0daughters(), @@ -1750,9 +2591,9 @@ struct ResonanceDaughterInitializer { casc.dcaXYCascToPV(), casc.dcaZCascToPV(), casc.sign(), - static_cast(casc.template posTrack_as().tpcNClsCrossedRows()), - static_cast(casc.template negTrack_as().tpcNClsCrossedRows()), - static_cast(casc.template bachelor_as().tpcNClsCrossedRows()), + static_cast(posTrack.tpcNClsCrossedRows()), + static_cast(negTrack.tpcNClsCrossedRows()), + static_cast(bachelor.tpcNClsCrossedRows()), casc.mLambda(), casc.mXi(), casc.v0radius(), casc.cascradius(), casc.x(), casc.y(), casc.z()); @@ -1846,9 +2687,9 @@ struct ResonanceDaughterInitializer { daughters = getDaughtersIndeces(cascmc); daughterPDGs = getDaughtersPDGCodes(cascmc); } - if (daughters.size() > StoredMCRelationCount) { - LOGF(info, "daughters.size() is larger than 2"); - } + // if (daughters.size() > StoredMCRelationCount) { + // LOGF(info, "daughters.size() is larger than 2"); + // } daughters.resize(StoredMCRelationCount, -1); daughterPDGs.resize(StoredMCRelationCount, -1); reso2mccascades(cascmc.pdgCode(), @@ -1882,7 +2723,7 @@ struct ResonanceDaughterInitializer { * * @param collision Collision data */ - void processDummy(aod::ResoCollision const&) + void processDummy(aod::ResoCollisions_001::iterator const&) { } PROCESS_SWITCH(ResonanceDaughterInitializer, processDummy, "Process dummy", true); @@ -1894,18 +2735,24 @@ struct ResonanceDaughterInitializer { * @param collision Reduced collision used as the output foreign key * @param tracks Tracks belonging to the corresponding original collision */ - template - void fillTrackTables(CollisionType const& collision, TrackTableType const& tracks) + template + void fillTrackTables(CollisionType const& collision, TrackTableType const& tracks, TrackPredicate const& keepTrack) { - fillTracks(collision, tracks); + fillTracks(collision, tracks, keepTrack); if (FilterForDerivedTables.cfgFillMicroTracks) { - fillMicroTracks(collision, tracks); + fillMicroTracks(collision, tracks, keepTrack); } if (FilterForDerivedTables.cfgFillUltraMicroTracks) { - fillUltraMicroTracks(collision, tracks); + fillUltraMicroTracks(collision, tracks, keepTrack); } } + template + void fillTrackTables(CollisionType const& collision, TrackTableType const& tracks) + { + fillTrackTables(collision, tracks, KeepAllTracks{}); + } + /** * @brief Fills track tables for one original collision * @@ -1939,49 +2786,94 @@ struct ResonanceDaughterInitializer { /** * @brief Processes data tracks using the two-stage hybrid grouping * - * The canonical fIndexCollisions column in ResoCollisionGroups lets - * GroupSlicer associate both reduced collisions and tracks to the same - * original aod::Collision. The tracks argument is therefore already the - * selected slice for this collision and must not be sliced again. + * GroupSlicer associates tracks automatically to the original + * aod::Collision. Reduced collisions retain a scalar original-collision row + * number and are explicitly sliced from the much smaller mapping table. The + * tracks argument is already the selected slice and must not be sliced again. */ - void processDataHybrid(aod::Collision const&, - soa::SmallGroups const& reducedCollisions, + void processDataHybrid(aod::Collision const& originalCollision, + SelectedResoCollisions const& reducedCollisions, soa::Filtered const& tracks) { - if (reducedCollisions.size() == 0) { + auto reducedCollisionsThisCollision = reducedCollisions.sliceBy(reducedCollisionsPerOriginalCollision, originalCollision.globalIndex()); + if (reducedCollisionsThisCollision.size() == 0) { return; } - if (reducedCollisions.size() != 1) { - LOGF(fatal, "Expected exactly one reduced collision for an original collision, found %zu", reducedCollisions.size()); + if (reducedCollisionsThisCollision.size() > 1) { + LOGF(error, "Found %zu reduced collisions for one original collision; skipping the ambiguous association", reducedCollisionsThisCollision.size()); + return; } - auto reducedCollision = reducedCollisions.begin(); + auto reducedCollision = reducedCollisionsThisCollision.begin(); fillTrackTables(reducedCollision, tracks); } PROCESS_SWITCH(ResonanceDaughterInitializer, processDataHybrid, "Process data tracks with the two-stage hybrid grouping", false); /** - * @brief Processes data tracks with configurable selected-V0 and selected-cascade gates + * @brief Processes data tracks with configurable selected V0 and cascade gates */ void processDataWithPairGate(ResoCollisionWithIndex::iterator const& collision, soa::Filtered const& tracks, aod::ResoV0Candidates const& v0s, aod::ResoCascadesCandidates const& cascades) { - if (FilterForDerivedTables.cfgBypassNoPairV0s) { - auto v0sThisCollision = v0s.sliceBy(v0sPerCollision, collision.collisionId()); - if (!hasSelectedV0(collision, v0sThisCollision, tracks)) { - return; - } + auto tracksThisCollision = tracks.sliceBy(tracksPerCollision, collision.collisionId()); + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + if (!preparePairTrackSelection(collision, tracks, v0s, cascades, + v0sPerCollision, cascadesPerCollision, pairSelection)) { + return; } - if (FilterForDerivedTables.cfgBypassNoPairCascades) { - auto cascadesThisCollision = cascades.sliceBy(cascadesPerCollision, collision.collisionId()); - if (!hasSelectedCascade(collision, cascadesThisCollision, tracks)) { - return; - } + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { + return; } - fillTrackTablesForCollision(collision, tracks, tracksPerCollision); + fillTrackTables(collision, tracksThisCollision, pairSelection); + } + PROCESS_SWITCH(ResonanceDaughterInitializer, processDataWithPairGate, "Process data tracks with the configured pair-gate mode", false); + + /** + * @brief Processes data tracks when a selectex0 and collision track exist + * + * This dedicated callback deliberately has no cascade input, so enabling + * the V0 pair gate does not activate the cascade upstream dependency. + */ + void processDataWithV0PairGate(ResoCollisionWithIndex::iterator const& collision, + soa::Filtered const& tracks, + aod::ResoV0Candidates const& v0s) + { + auto tracksThisCollision = tracks.sliceBy(tracksPerCollision, collision.collisionId()); + auto v0sThisCollision = v0s.sliceBy(v0sPerCollision, collision.collisionId()); + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + pairSelection.useV0Candidates = true; + pairSelection.v0Candidates = + collectSelectedV0Daughters(collision, v0sThisCollision, tracks, pairSelection.useGlobalDaughterVeto); + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { + return; + } + fillTrackTables(collision, tracksThisCollision, pairSelection); } - PROCESS_SWITCH(ResonanceDaughterInitializer, processDataWithPairGate, "Process data tracks with configurable pair gates", false); + PROCESS_SWITCH(ResonanceDaughterInitializer, processDataWithV0PairGate, "Process data tracks requiring a selected V0", false); + + /** + * @brief Processes data tracks when a selected cascade and collision track exist + * + * This dedicated callback deliberately has no V0 input, so enabling the + * cascade pair gate does not activate the V0 upstream dependency. + */ + void processDataWithCascPairGate(ResoCollisionWithIndex::iterator const& collision, + soa::Filtered const& tracks, + aod::ResoCascadesCandidates const& cascades) + { + auto tracksThisCollision = tracks.sliceBy(tracksPerCollision, collision.collisionId()); + auto cascadesThisCollision = cascades.sliceBy(cascadesPerCollision, collision.collisionId()); + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + pairSelection.useCascadeCandidates = true; + pairSelection.cascadeCandidates = + collectSelectedCascadeDaughters(collision, cascadesThisCollision, tracks, pairSelection.useGlobalDaughterVeto); + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { + return; + } + fillTrackTables(collision, tracksThisCollision, pairSelection); + } + PROCESS_SWITCH(ResonanceDaughterInitializer, processDataWithCascPairGate, "Process data tracks requiring a selected cascade", false); /** * @brief Processes MC tracks @@ -1999,7 +2891,7 @@ struct ResonanceDaughterInitializer { PROCESS_SWITCH(ResonanceDaughterInitializer, processMC, "Process tracks for MC", false); /** - * @brief Processes MC tracks with configurable V0 and cascade candidate gates + * @brief Processes MC tracks with configurable selected V0 and cascade gates */ void processMCWithPairGate(ResoCollisionWithIndex::iterator const& collision, soa::Filtered const& tracks, @@ -2007,17 +2899,60 @@ struct ResonanceDaughterInitializer { aod::ResoCascadesCandidatesMC const& cascades, aod::McParticles const&) { + auto tracksThisCollision = tracks.sliceBy(tracksMCPerCollision, collision.collisionId()); + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + if (!preparePairTrackSelection(collision, tracks, v0s, cascades, + v0sMCPerCollision, cascadesMCPerCollision, pairSelection)) { + return; + } + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { + return; + } + fillTrackTables(collision, tracksThisCollision, pairSelection); + } + PROCESS_SWITCH(ResonanceDaughterInitializer, processMCWithPairGate, "Process MC tracks with the configured pair-gate mode", false); + + /** + * @brief Processes MC tracks when a selected V0 and collision track exist + */ + void processMCWithV0PairGate(ResoCollisionWithIndex::iterator const& collision, + soa::Filtered const& tracks, + aod::ResoV0CandidatesMC const& v0s, + aod::McParticles const&) + { + auto tracksThisCollision = tracks.sliceBy(tracksMCPerCollision, collision.collisionId()); auto v0sThisCollision = v0s.sliceBy(v0sMCPerCollision, collision.collisionId()); - if (FilterForDerivedTables.cfgBypassNoPairV0s && v0sThisCollision.size() < 1) { + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + pairSelection.useV0Candidates = true; + pairSelection.v0Candidates = + collectSelectedV0Daughters(collision, v0sThisCollision, tracks, pairSelection.useGlobalDaughterVeto); + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { return; } + fillTrackTables(collision, tracksThisCollision, pairSelection); + } + PROCESS_SWITCH(ResonanceDaughterInitializer, processMCWithV0PairGate, "Process MC tracks requiring a selected V0", false); + + /** + * @brief Processes MC tracks when a selected cascade and collision track exist + */ + void processMCWithCascPairGate(ResoCollisionWithIndex::iterator const& collision, + soa::Filtered const& tracks, + aod::ResoCascadesCandidatesMC const& cascades, + aod::McParticles const&) + { + auto tracksThisCollision = tracks.sliceBy(tracksMCPerCollision, collision.collisionId()); auto cascadesThisCollision = cascades.sliceBy(cascadesMCPerCollision, collision.collisionId()); - if (FilterForDerivedTables.cfgBypassNoPairCascades && cascadesThisCollision.size() < 1) { + PairTrackSelection pairSelection{FilterForDerivedTables.cfgGlobalDaughterVeto.value}; + pairSelection.useCascadeCandidates = true; + pairSelection.cascadeCandidates = + collectSelectedCascadeDaughters(collision, cascadesThisCollision, tracks, pairSelection.useGlobalDaughterVeto); + if (!hasTracksForEnabledOutputs(collision, tracksThisCollision, pairSelection)) { return; } - fillTrackTablesForCollision(collision, tracks, tracksMCPerCollision); + fillTrackTables(collision, tracksThisCollision, pairSelection); } - PROCESS_SWITCH(ResonanceDaughterInitializer, processMCWithPairGate, "Process MC tracks with configurable pair gates", false); + PROCESS_SWITCH(ResonanceDaughterInitializer, processMCWithCascPairGate, "Process MC tracks requiring a selected cascade", false); /** * @brief Processes V0 data @@ -2033,29 +2968,6 @@ struct ResonanceDaughterInitializer { } PROCESS_SWITCH(ResonanceDaughterInitializer, processV0Data, "Process V0s for data", false); - /** - * @brief Processes data V0s grouped automatically by their original collision - * - * Both V0s and tracks are already restricted to the current original - * aod::Collision by GroupSlicer. The unfiltered track table is required for - * resolving the positive and negative daughter indices. - */ - void processV0DataHybrid(aod::Collision const&, - soa::SmallGroups const& reducedCollisions, - aod::ResoV0Candidates const& v0s, - aod::ResoTrackCandidates const& tracks) - { - if (reducedCollisions.size() == 0) { - return; - } - if (reducedCollisions.size() != 1) { - LOGF(fatal, "Expected exactly one reduced collision for an original collision, found %zu", reducedCollisions.size()); - } - auto reducedCollision = reducedCollisions.begin(); - fillV0s(reducedCollision, v0s, tracks); - } - PROCESS_SWITCH(ResonanceDaughterInitializer, processV0DataHybrid, "Process data V0s with the two-stage hybrid grouping", false); - /** * @brief Processes MC V0 data * @@ -2084,29 +2996,6 @@ struct ResonanceDaughterInitializer { } PROCESS_SWITCH(ResonanceDaughterInitializer, processCascData, "Process Cascades for data", false); - /** - * @brief Processes data cascades grouped automatically by their original collision - * - * Cascades and tracks arrive as original-collision groups. Keeping this as a - * separate callback from V0 processing avoids enabling either upstream input - * dependency unless its process switch is selected. - */ - void processCascDataHybrid(aod::Collision const&, - soa::SmallGroups const& reducedCollisions, - aod::ResoCascadesCandidates const& cascades, - aod::ResoTrackCandidates const& tracks) - { - if (reducedCollisions.size() == 0) { - return; - } - if (reducedCollisions.size() != 1) { - LOGF(fatal, "Expected exactly one reduced collision for an original collision, found %zu", reducedCollisions.size()); - } - auto reducedCollision = reducedCollisions.begin(); - fillCascades(reducedCollision, cascades, tracks); - } - PROCESS_SWITCH(ResonanceDaughterInitializer, processCascDataHybrid, "Process data cascades with the two-stage hybrid grouping", false); - /** * @brief Processes MC cascade data *