From 377d819ba6c5c9574c4eef0b0ce968ae29680fb5 Mon Sep 17 00:00:00 2001 From: Jseo Date: Fri, 11 Sep 2026 15:05:56 +0200 Subject: [PATCH 1/4] Fix event mixing bug --- PWGDQ/Core/MixingHandler.cxx | 42 +++- PWGDQ/Core/MixingHandler.h | 58 ++++- PWGDQ/Tasks/tableReader.cxx | 4 + PWGDQ/Tasks/tableReader_withAssoc.h | 228 +++++++++++++------ PWGDQ/Tasks/tableReader_withAssoc_direct.cxx | 14 +- 5 files changed, 268 insertions(+), 78 deletions(-) diff --git a/PWGDQ/Core/MixingHandler.cxx b/PWGDQ/Core/MixingHandler.cxx index 7d3a8e49f62..e1f82cb6287 100644 --- a/PWGDQ/Core/MixingHandler.cxx +++ b/PWGDQ/Core/MixingHandler.cxx @@ -62,6 +62,8 @@ void MixingHandler::AddMixingVariable(int var, const std::vector& binLims { fVariables[var] = fVariableLimits.size(); fVariableLimits.push_back(binLims); + // FillEvent() only fills variables marked as used + VarManager::SetUseVariable(var); } /* @@ -125,6 +127,8 @@ int MixingHandler::FindEventCategory(float* values) // loop over the variables and find out in which bin the value of the variable for the event is located std::vector bin; + // number of bins per variable in the iteration order of fVariables (fVariableLimits is in insertion order) + std::vector nBins; for (auto [var, pos] : fVariables) { // check that the value is within limits, if not return -1 to exclude the event from mixing size_t binValue = std::distance(fVariableLimits[pos].begin(), std::upper_bound(fVariableLimits[pos].begin(), fVariableLimits[pos].end(), values[var])); @@ -132,6 +136,7 @@ int MixingHandler::FindEventCategory(float* values) return -1; // all variables must be inside limits } bin.push_back(binValue - 1); + nBins.push_back(fVariableLimits[pos].size() - 1); } // Hash the bin values to define a unique category @@ -149,7 +154,7 @@ int MixingHandler::FindEventCategory(float* values) if (iv2 == iv1) { tempCategory *= bin[iv2]; } else { - tempCategory *= (fVariableLimits[iv2].size() - 1); + tempCategory *= nBins[iv2]; } } category += tempCategory; @@ -167,15 +172,40 @@ int MixingHandler::GetBinFromCategory(VarManager::Variables var, int category) c return -1; } - // Search for the position of the variable "var" in the internal variable list of the handler - int ivar = fVariables.at(var); + // number of bins and position of var in the iteration order of fVariables, as used by FindEventCategory() + std::vector nBins; + int ivar = -1; + for (auto [v, pos] : fVariables) { + if (v == var) { + ivar = static_cast(nBins.size()); + } + nBins.push_back(fVariableLimits[pos].size() - 1); + } + if (ivar < 0) { + return -1; + } // extract the bin position in variable "var" from the category int norm = 1; - for (int i = fVariables.size() - 1; i > ivar; --i) { - norm *= (fVariableLimits[i].size() - 1); + for (size_t i = nBins.size() - 1; i > static_cast(ivar); --i) { + norm *= nBins[i]; } int truncatedCategory = category - (category % norm); truncatedCategory /= norm; - return truncatedCategory % (fVariableLimits[ivar].size() - 1); + return truncatedCategory % nBins[ivar]; +} + +//_________________________________________________________________________ +void MixingHandler::SetCategoryBinCenters(int category, float* values) const +{ + // + // set the mixing variables to the bin centers of this category (used for the leftover mixing) + // + for (auto [var, pos] : fVariables) { + int bin = GetBinFromCategory(static_cast(var), category); + if (bin < 0) { + continue; + } + values[var] = 0.5 * (fVariableLimits[pos][bin] + fVariableLimits[pos][bin + 1]); + } } diff --git a/PWGDQ/Core/MixingHandler.h b/PWGDQ/Core/MixingHandler.h index 931db50f8a2..3fc6d4bcbc2 100644 --- a/PWGDQ/Core/MixingHandler.h +++ b/PWGDQ/Core/MixingHandler.h @@ -23,6 +23,7 @@ #include +#include #include #include #include @@ -39,11 +40,19 @@ class MixingHandler : public TNamed float eta; float phi; uint32_t filteringFlags; + // globalIndex is unique only within a dataframe, so the dataframe sequence is part of the track identity + uint64_t dataFrameSequence = 0; + uint64_t trackGlobalIndex = 0; + bool IsSamePhysicalTrack(const MixingTrack& other) const + { + return dataFrameSequence == other.dataFrameSequence && trackGlobalIndex == other.trackGlobalIndex; + } // Clear a bit once the track was used in mixing for that bit for the required pool depth. void ClearBit(uint32_t mask) { filteringFlags &= ~mask; } void Print() const { - std::cout << "pt: " << pt << ", eta: " << eta << ", phi: " << phi << ", filteringFlags: " << filteringFlags << std::endl; + std::cout << "pt: " << pt << ", eta: " << eta << ", phi: " << phi << ", filteringFlags: " << filteringFlags + << ", dataframe: " << dataFrameSequence << ", track: " << trackGlobalIndex << std::endl; } }; @@ -70,6 +79,19 @@ class MixingHandler : public TNamed } // Clear bits in the filtering mask. void ClearFilteringMask(uint32_t mask) { filteringMask &= ~mask; } + // clear the cut bits from all tracks and the filtering mask; remove tracks with no active bits left + void ClearBits(uint32_t mask) + { + for (auto& track : tracks1) { + track.ClearBit(mask); + } + tracks1.erase(std::remove_if(tracks1.begin(), tracks1.end(), [](auto const& track) { return track.filteringFlags == 0; }), tracks1.end()); + for (auto& track : tracks2) { + track.ClearBit(mask); + } + tracks2.erase(std::remove_if(tracks2.begin(), tracks2.end(), [](auto const& track) { return track.filteringFlags == 0; }), tracks2.end()); + ClearFilteringMask(mask); + } // 1) increment the counters for a given track cut bit mask and if the counters reached the pool depth, // 2) clear the corresponding bit in the tracks filtering flags to exclude them from further mixing // 3) for each track, if there are no more active bits in the filtering mask, then remove the track from the event @@ -157,6 +179,35 @@ class MixingHandler : public TNamed CleanPool(); events.push_back(event); } + // fixed-block mixing: AddEvent() until GetMixingMask() reports full cuts, mix, then ClearBits() + void AddEvent(const MixingEvent& event) { events.push_back(event); } + // bit mask of the cuts for which at least poolDepth events are in the pool + uint32_t GetMixingMask(int16_t poolDepth) const + { + std::array counts = {0}; + for (auto const& event : events) { + for (int icut = 0; icut < 32; ++icut) { + if (event.filteringMask & (static_cast(1) << icut)) { + counts[icut]++; + } + } + } + uint32_t fullMask = 0; + for (int icut = 0; icut < 32; ++icut) { + if (counts[icut] >= poolDepth) { + fullMask |= static_cast(1) << icut; + } + } + return fullMask; + } + // clear the given cut bits from all events in the pool and remove the events with no tracks left + void ClearBits(uint32_t mask) + { + for (auto& event : events) { + event.ClearBits(mask); + } + CleanPool(); + } // getter for the events in the pool const std::vector& GetEvents() const { return events; } @@ -176,17 +227,22 @@ class MixingHandler : public TNamed // setters void AddMixingVariable(int var, const std::vector& binLims); void SetPoolDepth(int16_t depth) { fPoolDepth = depth; } + // remove all pools (e.g. at a run change) + void ClearPools() { fPools.clear(); } // getters // int GetNMixingVariables() const { return fVariables.size(); } // int GetMixingVariable(VarManager::Variables var); // returns the position in the internal varible list of the handler. Useful for checks, mostly // std::vector GetMixingVariableLimits(VarManager::Variables var); MixingPool& GetPool(int category) { return fPools[category]; } + std::map& GetPools() { return fPools; } int16_t GetPoolDepth() const { return fPoolDepth; } void Init(); int FindEventCategory(float* values); int GetBinFromCategory(VarManager::Variables var, int category) const; + // set the mixing variables to the bin centers of the given category + void SetCategoryBinCenters(int category, float* values) const; private: MixingHandler(const MixingHandler& handler); diff --git a/PWGDQ/Tasks/tableReader.cxx b/PWGDQ/Tasks/tableReader.cxx index 8501e0792e5..ef3e6713b2f 100644 --- a/PWGDQ/Tasks/tableReader.cxx +++ b/PWGDQ/Tasks/tableReader.cxx @@ -273,6 +273,10 @@ struct AnalysisEventSelection { if (fMixHandler != nullptr) { int hh = fMixHandler->FindEventCategory(VarManager::fgValues); + // events outside the mixing limits (-1) get a distinct negative hash so that they are not mixed with each other + if (hh < 0) { + hh = -1 - static_cast(event.globalIndex()); + } hash(hh); } } diff --git a/PWGDQ/Tasks/tableReader_withAssoc.h b/PWGDQ/Tasks/tableReader_withAssoc.h index 541a8f039cb..983dfca7d1b 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.h +++ b/PWGDQ/Tasks/tableReader_withAssoc.h @@ -73,6 +73,7 @@ #include #include #include +#include #include #include #include @@ -81,6 +82,26 @@ #include #include +// OutputObj with a callback run on the first dereference, i.e. right before the end-of-stream snapshot +// (tasks cannot register their own EndOfStream callback, CallbackService::set replaces the framework one) +template +struct FinalizingOutputObj : public o2::framework::OutputObj { + using o2::framework::OutputObj::OutputObj; + std::function finalizeBeforeSnapshot; + bool finalized = false; + + T& operator*() + { + if (!finalized) { + finalized = true; + if (finalizeBeforeSnapshot) { + finalizeBeforeSnapshot(); + } + } + return o2::framework::OutputObj::operator*(); + } +}; + // Some definitions namespace o2::aod { @@ -494,6 +515,10 @@ struct AnalysisEventSelection { // create the mixing hash and publish it into the hash table if (fMixHandler != nullptr) { int hh = fMixHandler->FindEventCategory(dqtablereader_helpers::varValues()); + // events outside the mixing limits (-1) get a distinct negative hash so that they are not mixed with each other + if (hh < 0) { + hh = -1 - static_cast(event.globalIndex()); + } hash(hh); } } @@ -1351,7 +1376,7 @@ struct AnalysisSameEventPairing { TH1D* ResoFlowEP = nullptr; int fCurrentRun = -1; // needed to detect if the run changed and trigger update of calibrations etc. - o2::framework::OutputObj fOutputList{"output"}; + FinalizingOutputObj fOutputList{"output"}; struct : o2::framework::ConfigurableGroup { o2::framework::Configurable track{"cfgTrackCuts", "jpsiO2MCdebugCuts2", "Comma separated list of barrel track cuts"}; @@ -1427,6 +1452,8 @@ struct AnalysisSameEventPairing { HistogramManager* fHistMan = nullptr; MixingHandler fMixingHandler; + // dataframe counter, part of the track identity in the mixing pools + uint64_t fMixingDataFrameSequence = 0; o2::analysis::DQMlResponse fDQMlResponse; std::vector fOutputMlPsi2ee; // TODO: check this is needed or not @@ -1714,6 +1741,12 @@ struct AnalysisSameEventPairing { } if (fConfigRunMixingAcrossTFs) { + if (fConfigMixingDepth.value < 2) { + LOGF(fatal, "cfgMixingDepth must be at least 2"); + } + if (fNCutsBarrel > 32) { + LOGF(fatal, "Across-TF mixing supports at most 32 barrel track-cut bits, got %d", fNCutsBarrel); + } TString mixVarsString = fConfigMixingVariables.value; TString mixVarsJsonString = fConfigMixingVariablesJson.value; std::unique_ptr objArray(mixVarsString.Tokenize(",")); @@ -1816,6 +1849,12 @@ struct AnalysisSameEventPairing { o2::aod::dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill fOutputList.setObject(fHistMan->GetMainHistogramList()); + // mix the events left in the pools at the end of the stream + fOutputList.finalizeBeforeSnapshot = [this]() { + if (fConfigRunMixingAcrossTFs) { + runLeftoverMixing(); + } + }; } } @@ -1893,12 +1932,110 @@ struct AnalysisSameEventPairing { } } + // Mix all the events in the pool for the cuts in mixingMask (event-wise variables are those of the current event) + void runEventMixing(MixingHandler::MixingPool& pool, uint32_t mixingMask) + { + auto const& events = pool.GetEvents(); + // each pair of events is mixed once + for (size_t iev1 = 0; iev1 < events.size(); iev1++) { + for (size_t iev2 = iev1 + 1; iev2 < events.size(); iev2++) { + auto const& mixingEvent = events[iev1]; + auto const& poolEvent = events[iev2]; + if (!(mixingEvent.filteringMask & poolEvent.filteringMask & mixingMask)) { + continue; + } + for (auto const& t1 : mixingEvent.tracks1) { + // run +- pairing + for (auto const& t2 : poolEvent.tracks2) { + // check the two-track filter for the mixed pair + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; + if (!mixedTwoTrackFilter) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < fNCutsBarrel; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + // run ++ pairing + for (auto const& t2 : poolEvent.tracks1) { + // check the two-track filter for the mixed pair and skip the same track associated to both collisions + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; + if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < fNCutsBarrel; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPP_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + } + for (auto const& t1 : mixingEvent.tracks2) { + // run -+ pairing + for (auto const& t2 : poolEvent.tracks1) { + // check the two-track filter for the mixed pair + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; + if (!mixedTwoTrackFilter) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < fNCutsBarrel; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + // run -- pairing + for (auto const& t2 : poolEvent.tracks2) { + // check the two-track filter for the mixed pair and skip the same track associated to both collisions + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; + if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < fNCutsBarrel; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEMM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + } + } + } + } + + // Mix the events left in the pools (end of stream or run change), with the mixing variables set to the category bin centers + void runLeftoverMixing() + { + for (auto& [category, pool] : fMixingHandler.GetPools()) { + const uint32_t mixingMask = pool.GetMixingMask(2); + if (!mixingMask) { + continue; + } + VarManager::ResetValues(0, VarManager::kNEventWiseVariables); + fMixingHandler.SetCategoryBinCenters(category, dqtablereader_helpers::varValues()); + runEventMixing(pool, mixingMask); + pool.ClearBits(mixingMask); + } + } + // Template function to run same event pairing (barrel-barrel, muon-muon, barrel-muon) template void runSameEventPairing(TEvents const& events, o2::framework::Preslice& preslice, TTrackAssocs const& assocs, TTracks const& /*tracks*/) { if (events.size() > 0) { // Additional protection to avoid crashing of events.begin().runNumber() if (fCurrentRun != events.begin().runNumber()) { + if (fConfigRunMixingAcrossTFs) { + // do not mix events from different runs + if (fCurrentRun >= 0) { + runLeftoverMixing(); + } + fMixingHandler.ClearPools(); + } initParamsFromCCDB(events.begin().timestamp(), events.begin().runNumber(), TTwoProngFitter); fCurrentRun = events.begin().runNumber(); } @@ -1981,6 +2118,10 @@ struct AnalysisSameEventPairing { // constexpr bool fillFlowReso = eventHasQvector || eventHasQvectorCentr; bool isSelectedBDT = false; fNPairPerEvent = 0; + uint64_t currentMixingDataFrameSequence = 0; + if (fConfigRunMixingAcrossTFs) { + currentMixingDataFrameSequence = ++fMixingDataFrameSequence; + } for (auto const& event : events) { if (!event.isEventSelected_bit(0)) { @@ -2006,6 +2147,10 @@ struct AnalysisSameEventPairing { } VarManager::FillEventFlowResoFactor(ResoFlowSP, ResoFlowEP); } + int mixingCategory = -1; + if (fConfigRunMixingAcrossTFs) { + mixingCategory = fMixingHandler.FindEventCategory(dqtablereader_helpers::varValues()); + } bool isFirst = true; for (auto const& [a1, a2] : o2::soa::combinations(groupedAssocs, groupedAssocs)) { @@ -2478,6 +2623,9 @@ struct AnalysisSameEventPairing { if (fConfigRunMixingAcrossTFs) { // run event mixing across TFs + if (mixingCategory < 0) { + continue; + } // 1) create a MixingEvent and fill it with the relevant tracks MixingHandler::MixingEvent mixingEvent; uint32_t trackFilterForMixing = 0; @@ -2488,7 +2636,7 @@ struct AnalysisSameEventPairing { continue; } auto t1 = assoc.template reducedtrack_as(); - MixingHandler::MixingTrack mixingTrack(t1.pt(), t1.eta(), t1.phi(), trackFilterForMixing); + MixingHandler::MixingTrack mixingTrack(t1.pt(), t1.eta(), t1.phi(), trackFilterForMixing, currentMixingDataFrameSequence, static_cast(assoc.reducedtrackId())); if (t1.sign() > 0) { mixingEvent.AddTrack1(mixingTrack); } else { @@ -2496,72 +2644,18 @@ struct AnalysisSameEventPairing { } } } - // 2) run the mixing with the events in the pool corresponding to this event - auto& pool = fMixingHandler.GetPool(fMixingHandler.FindEventCategory(dqtablereader_helpers::varValues())); - for (auto const& poolEvent : pool.GetEvents()) { - for (auto const& t1 : mixingEvent.tracks1) { - // run +- pairing - for (auto const& t2 : poolEvent.tracks2) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < ncuts; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - // run ++ pairing - for (auto const& t2 : poolEvent.tracks1) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < ncuts; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPP_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - } - for (auto const& t1 : mixingEvent.tracks2) { - // run -+ pairing - for (auto const& t2 : poolEvent.tracks1) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < ncuts; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - // run -- pairing - for (auto const& t2 : poolEvent.tracks2) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < ncuts; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEMM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - } + if (mixingEvent.tracks1.empty() && mixingEvent.tracks2.empty()) { + continue; + } + // 2) add the event to the pool corresponding to this event + auto& pool = fMixingHandler.GetPool(mixingCategory); + pool.AddEvent(mixingEvent); + // 3) mix all the events in the pool for the cuts which reached the pool depth + uint32_t mixingMask = pool.GetMixingMask(fMixingHandler.GetPoolDepth()); + if (mixingMask) { + runEventMixing(pool, mixingMask); + pool.ClearBits(mixingMask); } - // 3) add the current event to the pool - pool.UpdatePool(mixingEvent, fMixingHandler.GetPoolDepth()); // pool.Print(); } } // end loop over events diff --git a/PWGDQ/Tasks/tableReader_withAssoc_direct.cxx b/PWGDQ/Tasks/tableReader_withAssoc_direct.cxx index 781930985c7..95fceb1e889 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc_direct.cxx +++ b/PWGDQ/Tasks/tableReader_withAssoc_direct.cxx @@ -494,6 +494,16 @@ struct AnalysisEventSelection { VarManager::FillBC(bc); VarManager::FillEvent(event); + // the hash table is joined to the events by row order, so publish one row per event before any event selection + if (fMixHandler != nullptr) { + int hh = fMixHandler->FindEventCategory(VarManager::fgValues); + // events outside the mixing limits (-1) get a distinct negative hash so that they are not mixed with each other + if (hh < 0) { + hh = -1 - static_cast(event.globalIndex()); + } + hash(hh); + } + bool decision = false; if (fConfigQA) { fHistMan->FillHistClass("Event_BeforeCuts", VarManager::fgValues); @@ -533,10 +543,6 @@ struct AnalysisEventSelection { auto& evIndices = fBCCollMap[bc.globalBC()]; evIndices.push_back(event.globalIndex()); } - if (fMixHandler != nullptr) { - int hh = fMixHandler->FindEventCategory(VarManager::fgValues); - hash(hh); - } } } From 287a879e5d3ee917acacb3377422b80bf1e0acf5 Mon Sep 17 00:00:00 2001 From: Jseo Date: Sat, 12 Sep 2026 16:01:29 +0200 Subject: [PATCH 2/4] o2 linter --- PWGDQ/Core/MixingHandler.cxx | 4 ++-- PWGDQ/Core/MixingHandler.h | 11 ++++++++--- PWGDQ/Tasks/tableReader_withAssoc.h | 12 ++++++------ 3 files changed, 16 insertions(+), 11 deletions(-) diff --git a/PWGDQ/Core/MixingHandler.cxx b/PWGDQ/Core/MixingHandler.cxx index e1f82cb6287..06497b1b2cb 100644 --- a/PWGDQ/Core/MixingHandler.cxx +++ b/PWGDQ/Core/MixingHandler.cxx @@ -175,7 +175,7 @@ int MixingHandler::GetBinFromCategory(VarManager::Variables var, int category) c // number of bins and position of var in the iteration order of fVariables, as used by FindEventCategory() std::vector nBins; int ivar = -1; - for (auto [v, pos] : fVariables) { + for (auto const& [v, pos] : fVariables) { if (v == var) { ivar = static_cast(nBins.size()); } @@ -201,7 +201,7 @@ void MixingHandler::SetCategoryBinCenters(int category, float* values) const // // set the mixing variables to the bin centers of this category (used for the leftover mixing) // - for (auto [var, pos] : fVariables) { + for (auto const& [var, pos] : fVariables) { int bin = GetBinFromCategory(static_cast(var), category); if (bin < 0) { continue; diff --git a/PWGDQ/Core/MixingHandler.h b/PWGDQ/Core/MixingHandler.h index 3fc6d4bcbc2..fb81eb42174 100644 --- a/PWGDQ/Core/MixingHandler.h +++ b/PWGDQ/Core/MixingHandler.h @@ -34,6 +34,11 @@ class MixingHandler : public TNamed { public: + // number of track cuts which fit in the 32-bit filtering masks + static constexpr int NMaxCuts = 32; + // smallest pool depth for which a mixed pair can be built + static constexpr int16_t MinPoolDepth = 2; + // Struct to define track properties relevant for mixing and few utility functions struct MixingTrack { float pt; @@ -184,16 +189,16 @@ class MixingHandler : public TNamed // bit mask of the cuts for which at least poolDepth events are in the pool uint32_t GetMixingMask(int16_t poolDepth) const { - std::array counts = {0}; + std::array counts = {0}; for (auto const& event : events) { - for (int icut = 0; icut < 32; ++icut) { + for (int icut = 0; icut < NMaxCuts; ++icut) { if (event.filteringMask & (static_cast(1) << icut)) { counts[icut]++; } } } uint32_t fullMask = 0; - for (int icut = 0; icut < 32; ++icut) { + for (int icut = 0; icut < NMaxCuts; ++icut) { if (counts[icut] >= poolDepth) { fullMask |= static_cast(1) << icut; } diff --git a/PWGDQ/Tasks/tableReader_withAssoc.h b/PWGDQ/Tasks/tableReader_withAssoc.h index 983dfca7d1b..ed81099e0da 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.h +++ b/PWGDQ/Tasks/tableReader_withAssoc.h @@ -1741,11 +1741,11 @@ struct AnalysisSameEventPairing { } if (fConfigRunMixingAcrossTFs) { - if (fConfigMixingDepth.value < 2) { - LOGF(fatal, "cfgMixingDepth must be at least 2"); + if (fConfigMixingDepth.value < MixingHandler::MinPoolDepth) { + LOGF(fatal, "cfgMixingDepth must be at least %d", MixingHandler::MinPoolDepth); } - if (fNCutsBarrel > 32) { - LOGF(fatal, "Across-TF mixing supports at most 32 barrel track-cut bits, got %d", fNCutsBarrel); + if (fNCutsBarrel > MixingHandler::NMaxCuts) { + LOGF(fatal, "Across-TF mixing supports at most %d barrel track-cut bits, got %d", MixingHandler::NMaxCuts, fNCutsBarrel); } TString mixVarsString = fConfigMixingVariables.value; TString mixVarsJsonString = fConfigMixingVariablesJson.value; @@ -2011,8 +2011,8 @@ struct AnalysisSameEventPairing { // Mix the events left in the pools (end of stream or run change), with the mixing variables set to the category bin centers void runLeftoverMixing() { - for (auto& [category, pool] : fMixingHandler.GetPools()) { - const uint32_t mixingMask = pool.GetMixingMask(2); + for (auto& [category, pool] : fMixingHandler.GetPools()) { // o2-linter: disable=const-ref-in-for-loop (the pools are modified) + const uint32_t mixingMask = pool.GetMixingMask(MixingHandler::MinPoolDepth); if (!mixingMask) { continue; } From e92b77d6d89169540e6df1a89e4f10719fbf16d4 Mon Sep 17 00:00:00 2001 From: Jseo Date: Tue, 15 Sep 2026 13:29:32 +0200 Subject: [PATCH 3/4] back to rolling mixing --- PWGDQ/Core/MixingHandler.cxx | 15 --- PWGDQ/Core/MixingHandler.h | 56 ++------ PWGDQ/Tasks/tableReader_withAssoc.h | 199 +++++++++------------------- 3 files changed, 77 insertions(+), 193 deletions(-) diff --git a/PWGDQ/Core/MixingHandler.cxx b/PWGDQ/Core/MixingHandler.cxx index 06497b1b2cb..7bdaa62f072 100644 --- a/PWGDQ/Core/MixingHandler.cxx +++ b/PWGDQ/Core/MixingHandler.cxx @@ -194,18 +194,3 @@ int MixingHandler::GetBinFromCategory(VarManager::Variables var, int category) c truncatedCategory /= norm; return truncatedCategory % nBins[ivar]; } - -//_________________________________________________________________________ -void MixingHandler::SetCategoryBinCenters(int category, float* values) const -{ - // - // set the mixing variables to the bin centers of this category (used for the leftover mixing) - // - for (auto const& [var, pos] : fVariables) { - int bin = GetBinFromCategory(static_cast(var), category); - if (bin < 0) { - continue; - } - values[var] = 0.5 * (fVariableLimits[pos][bin] + fVariableLimits[pos][bin + 1]); - } -} diff --git a/PWGDQ/Core/MixingHandler.h b/PWGDQ/Core/MixingHandler.h index fb81eb42174..53c93be09fd 100644 --- a/PWGDQ/Core/MixingHandler.h +++ b/PWGDQ/Core/MixingHandler.h @@ -23,7 +23,6 @@ #include -#include #include #include #include @@ -36,8 +35,6 @@ class MixingHandler : public TNamed public: // number of track cuts which fit in the 32-bit filtering masks static constexpr int NMaxCuts = 32; - // smallest pool depth for which a mixed pair can be built - static constexpr int16_t MinPoolDepth = 2; // Struct to define track properties relevant for mixing and few utility functions struct MixingTrack { @@ -48,6 +45,8 @@ class MixingHandler : public TNamed // globalIndex is unique only within a dataframe, so the dataframe sequence is part of the track identity uint64_t dataFrameSequence = 0; uint64_t trackGlobalIndex = 0; + // electric charge of the track; 0 means "not set" and disables the charge dependent pair variables + int8_t sign = 0; bool IsSamePhysicalTrack(const MixingTrack& other) const { return dataFrameSequence == other.dataFrameSequence && trackGlobalIndex == other.trackGlobalIndex; @@ -56,7 +55,8 @@ class MixingHandler : public TNamed void ClearBit(uint32_t mask) { filteringFlags &= ~mask; } void Print() const { - std::cout << "pt: " << pt << ", eta: " << eta << ", phi: " << phi << ", filteringFlags: " << filteringFlags + std::cout << "pt: " << pt << ", eta: " << eta << ", phi: " << phi << ", sign: " << static_cast(sign) + << ", filteringFlags: " << filteringFlags << ", dataframe: " << dataFrameSequence << ", track: " << trackGlobalIndex << std::endl; } }; @@ -84,19 +84,6 @@ class MixingHandler : public TNamed } // Clear bits in the filtering mask. void ClearFilteringMask(uint32_t mask) { filteringMask &= ~mask; } - // clear the cut bits from all tracks and the filtering mask; remove tracks with no active bits left - void ClearBits(uint32_t mask) - { - for (auto& track : tracks1) { - track.ClearBit(mask); - } - tracks1.erase(std::remove_if(tracks1.begin(), tracks1.end(), [](auto const& track) { return track.filteringFlags == 0; }), tracks1.end()); - for (auto& track : tracks2) { - track.ClearBit(mask); - } - tracks2.erase(std::remove_if(tracks2.begin(), tracks2.end(), [](auto const& track) { return track.filteringFlags == 0; }), tracks2.end()); - ClearFilteringMask(mask); - } // 1) increment the counters for a given track cut bit mask and if the counters reached the pool depth, // 2) clear the corresponding bit in the tracks filtering flags to exclude them from further mixing // 3) for each track, if there are no more active bits in the filtering mask, then remove the track from the event @@ -184,34 +171,16 @@ class MixingHandler : public TNamed CleanPool(); events.push_back(event); } - // fixed-block mixing: AddEvent() until GetMixingMask() reports full cuts, mix, then ClearBits() - void AddEvent(const MixingEvent& event) { events.push_back(event); } - // bit mask of the cuts for which at least poolDepth events are in the pool - uint32_t GetMixingMask(int16_t poolDepth) const + // Same, but the stored events are aged only for the cuts in agingMask. Passing the filtering mask of the + // incoming event ages an event only for the cuts for which a mixed pair was actually produced, so that the + // pool depth is a number of mixed partners and not a number of arrivals. + void UpdatePool(const MixingEvent& event, int16_t poolDepth, uint32_t agingMask) { - std::array counts = {0}; - for (auto const& event : events) { - for (int icut = 0; icut < NMaxCuts; ++icut) { - if (event.filteringMask & (static_cast(1) << icut)) { - counts[icut]++; - } - } - } - uint32_t fullMask = 0; - for (int icut = 0; icut < NMaxCuts; ++icut) { - if (counts[icut] >= poolDepth) { - fullMask |= static_cast(1) << icut; - } - } - return fullMask; - } - // clear the given cut bits from all events in the pool and remove the events with no tracks left - void ClearBits(uint32_t mask) - { - for (auto& event : events) { - event.ClearBits(mask); + for (auto& poolEvent : events) { // o2-linter: disable=const-ref-in-for-loop (the events are modified) + poolEvent.IncrementCounters(agingMask, poolDepth); } CleanPool(); + events.push_back(event); } // getter for the events in the pool const std::vector& GetEvents() const { return events; } @@ -240,14 +209,11 @@ class MixingHandler : public TNamed // int GetMixingVariable(VarManager::Variables var); // returns the position in the internal varible list of the handler. Useful for checks, mostly // std::vector GetMixingVariableLimits(VarManager::Variables var); MixingPool& GetPool(int category) { return fPools[category]; } - std::map& GetPools() { return fPools; } int16_t GetPoolDepth() const { return fPoolDepth; } void Init(); int FindEventCategory(float* values); int GetBinFromCategory(VarManager::Variables var, int category) const; - // set the mixing variables to the bin centers of the given category - void SetCategoryBinCenters(int category, float* values) const; private: MixingHandler(const MixingHandler& handler); diff --git a/PWGDQ/Tasks/tableReader_withAssoc.h b/PWGDQ/Tasks/tableReader_withAssoc.h index ed81099e0da..b8028f2b773 100644 --- a/PWGDQ/Tasks/tableReader_withAssoc.h +++ b/PWGDQ/Tasks/tableReader_withAssoc.h @@ -73,7 +73,6 @@ #include #include #include -#include #include #include #include @@ -82,26 +81,6 @@ #include #include -// OutputObj with a callback run on the first dereference, i.e. right before the end-of-stream snapshot -// (tasks cannot register their own EndOfStream callback, CallbackService::set replaces the framework one) -template -struct FinalizingOutputObj : public o2::framework::OutputObj { - using o2::framework::OutputObj::OutputObj; - std::function finalizeBeforeSnapshot; - bool finalized = false; - - T& operator*() - { - if (!finalized) { - finalized = true; - if (finalizeBeforeSnapshot) { - finalizeBeforeSnapshot(); - } - } - return o2::framework::OutputObj::operator*(); - } -}; - // Some definitions namespace o2::aod { @@ -1376,7 +1355,7 @@ struct AnalysisSameEventPairing { TH1D* ResoFlowEP = nullptr; int fCurrentRun = -1; // needed to detect if the run changed and trigger update of calibrations etc. - FinalizingOutputObj fOutputList{"output"}; + o2::framework::OutputObj fOutputList{"output"}; struct : o2::framework::ConfigurableGroup { o2::framework::Configurable track{"cfgTrackCuts", "jpsiO2MCdebugCuts2", "Comma separated list of barrel track cuts"}; @@ -1741,9 +1720,6 @@ struct AnalysisSameEventPairing { } if (fConfigRunMixingAcrossTFs) { - if (fConfigMixingDepth.value < MixingHandler::MinPoolDepth) { - LOGF(fatal, "cfgMixingDepth must be at least %d", MixingHandler::MinPoolDepth); - } if (fNCutsBarrel > MixingHandler::NMaxCuts) { LOGF(fatal, "Across-TF mixing supports at most %d barrel track-cut bits, got %d", MixingHandler::NMaxCuts, fNCutsBarrel); } @@ -1849,12 +1825,6 @@ struct AnalysisSameEventPairing { o2::aod::dqhistograms::AddHistogramsFromJSON(fHistMan, fConfigAddJSONHistograms.value.c_str()); // ad-hoc histograms via JSON VarManager::SetUseVars(fHistMan->GetUsedVars()); // provide the list of required variables so that VarManager knows what to fill fOutputList.setObject(fHistMan->GetMainHistogramList()); - // mix the events left in the pools at the end of the stream - fOutputList.finalizeBeforeSnapshot = [this]() { - if (fConfigRunMixingAcrossTFs) { - runLeftoverMixing(); - } - }; } } @@ -1932,97 +1902,6 @@ struct AnalysisSameEventPairing { } } - // Mix all the events in the pool for the cuts in mixingMask (event-wise variables are those of the current event) - void runEventMixing(MixingHandler::MixingPool& pool, uint32_t mixingMask) - { - auto const& events = pool.GetEvents(); - // each pair of events is mixed once - for (size_t iev1 = 0; iev1 < events.size(); iev1++) { - for (size_t iev2 = iev1 + 1; iev2 < events.size(); iev2++) { - auto const& mixingEvent = events[iev1]; - auto const& poolEvent = events[iev2]; - if (!(mixingEvent.filteringMask & poolEvent.filteringMask & mixingMask)) { - continue; - } - for (auto const& t1 : mixingEvent.tracks1) { - // run +- pairing - for (auto const& t2 : poolEvent.tracks2) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < fNCutsBarrel; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - // run ++ pairing - for (auto const& t2 : poolEvent.tracks1) { - // check the two-track filter for the mixed pair and skip the same track associated to both collisions - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; - if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < fNCutsBarrel; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPP_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - } - for (auto const& t1 : mixingEvent.tracks2) { - // run -+ pairing - for (auto const& t2 : poolEvent.tracks1) { - // check the two-track filter for the mixed pair - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; - if (!mixedTwoTrackFilter) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < fNCutsBarrel; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - // run -- pairing - for (auto const& t2 : poolEvent.tracks2) { - // check the two-track filter for the mixed pair and skip the same track associated to both collisions - uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags & mixingMask; - if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { - continue; - } - VarManager::FillPairMEAcrossTFs(t1, t2); - for (int icut = 0; icut < fNCutsBarrel; icut++) { - if (mixedTwoTrackFilter & (static_cast(1) << icut)) { - fHistMan->FillHistClass(Form("PairsBarrelMEMM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); - } - } - } - } - } - } - } - - // Mix the events left in the pools (end of stream or run change), with the mixing variables set to the category bin centers - void runLeftoverMixing() - { - for (auto& [category, pool] : fMixingHandler.GetPools()) { // o2-linter: disable=const-ref-in-for-loop (the pools are modified) - const uint32_t mixingMask = pool.GetMixingMask(MixingHandler::MinPoolDepth); - if (!mixingMask) { - continue; - } - VarManager::ResetValues(0, VarManager::kNEventWiseVariables); - fMixingHandler.SetCategoryBinCenters(category, dqtablereader_helpers::varValues()); - runEventMixing(pool, mixingMask); - pool.ClearBits(mixingMask); - } - } - // Template function to run same event pairing (barrel-barrel, muon-muon, barrel-muon) template void runSameEventPairing(TEvents const& events, o2::framework::Preslice& preslice, TTrackAssocs const& assocs, TTracks const& /*tracks*/) @@ -2031,9 +1910,6 @@ struct AnalysisSameEventPairing { if (fCurrentRun != events.begin().runNumber()) { if (fConfigRunMixingAcrossTFs) { // do not mix events from different runs - if (fCurrentRun >= 0) { - runLeftoverMixing(); - } fMixingHandler.ClearPools(); } initParamsFromCCDB(events.begin().timestamp(), events.begin().runNumber(), TTwoProngFitter); @@ -2636,7 +2512,7 @@ struct AnalysisSameEventPairing { continue; } auto t1 = assoc.template reducedtrack_as(); - MixingHandler::MixingTrack mixingTrack(t1.pt(), t1.eta(), t1.phi(), trackFilterForMixing, currentMixingDataFrameSequence, static_cast(assoc.reducedtrackId())); + MixingHandler::MixingTrack mixingTrack(t1.pt(), t1.eta(), t1.phi(), trackFilterForMixing, currentMixingDataFrameSequence, static_cast(assoc.reducedtrackId()), static_cast(t1.sign())); if (t1.sign() > 0) { mixingEvent.AddTrack1(mixingTrack); } else { @@ -2647,15 +2523,72 @@ struct AnalysisSameEventPairing { if (mixingEvent.tracks1.empty() && mixingEvent.tracks2.empty()) { continue; } - // 2) add the event to the pool corresponding to this event + // 2) run the mixing with the events in the pool corresponding to this event auto& pool = fMixingHandler.GetPool(mixingCategory); - pool.AddEvent(mixingEvent); - // 3) mix all the events in the pool for the cuts which reached the pool depth - uint32_t mixingMask = pool.GetMixingMask(fMixingHandler.GetPoolDepth()); - if (mixingMask) { - runEventMixing(pool, mixingMask); - pool.ClearBits(mixingMask); + for (auto const& poolEvent : pool.GetEvents()) { + for (auto const& t1 : mixingEvent.tracks1) { + // run +- pairing + for (auto const& t2 : poolEvent.tracks2) { + // check the two-track filter for the mixed pair + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; + if (!mixedTwoTrackFilter) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < ncuts; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + // run ++ pairing + for (auto const& t2 : poolEvent.tracks1) { + // check the two-track filter for the mixed pair and skip the same track associated to both collisions + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; + if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < ncuts; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPP_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + } + for (auto const& t1 : mixingEvent.tracks2) { + // run -+ pairing + for (auto const& t2 : poolEvent.tracks1) { + // check the two-track filter for the mixed pair + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; + if (!mixedTwoTrackFilter) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < ncuts; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEPM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + // run -- pairing + for (auto const& t2 : poolEvent.tracks2) { + // check the two-track filter for the mixed pair and skip the same track associated to both collisions + uint32_t mixedTwoTrackFilter = t1.filteringFlags & t2.filteringFlags; + if (!mixedTwoTrackFilter || t1.IsSamePhysicalTrack(t2)) { + continue; + } + VarManager::FillPairMEAcrossTFs(t1, t2); + for (int icut = 0; icut < ncuts; icut++) { + if (mixedTwoTrackFilter & (static_cast(1) << icut)) { + fHistMan->FillHistClass(Form("PairsBarrelMEMM_%s", fTrackCuts[icut].Data()), dqtablereader_helpers::varValues()); + } + } + } + } } + // 3) add the current event to the pool + pool.UpdatePool(mixingEvent, fMixingHandler.GetPoolDepth(), mixingEvent.filteringMask); // pool.Print(); } } // end loop over events From 3804722d041b00abc84189ae63058813941adb5e Mon Sep 17 00:00:00 2001 From: Jseo Date: Tue, 15 Sep 2026 13:39:12 +0200 Subject: [PATCH 4/4] fix --- PWGDQ/Core/MixingHandler.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGDQ/Core/MixingHandler.h b/PWGDQ/Core/MixingHandler.h index 53c93be09fd..839666f36ec 100644 --- a/PWGDQ/Core/MixingHandler.h +++ b/PWGDQ/Core/MixingHandler.h @@ -176,7 +176,7 @@ class MixingHandler : public TNamed // pool depth is a number of mixed partners and not a number of arrivals. void UpdatePool(const MixingEvent& event, int16_t poolDepth, uint32_t agingMask) { - for (auto& poolEvent : events) { // o2-linter: disable=const-ref-in-for-loop (the events are modified) + for (auto& poolEvent : events) { poolEvent.IncrementCounters(agingMask, poolDepth); } CleanPool();