From ba2795e892f546b344f7cc9b5c85adb2a1727ba2 Mon Sep 17 00:00:00 2001 From: Evgeny Kryshen Date: Thu, 29 Oct 2026 18:53:47 +0300 Subject: [PATCH 1/5] New functionality for Pb-Pb analyses --- PWGUD/TableProducer/upcCandProducerBarrel.cxx | 376 +++++++++++++++--- 1 file changed, 311 insertions(+), 65 deletions(-) diff --git a/PWGUD/TableProducer/upcCandProducerBarrel.cxx b/PWGUD/TableProducer/upcCandProducerBarrel.cxx index 4d2a75be18a..70170fe2367 100644 --- a/PWGUD/TableProducer/upcCandProducerBarrel.cxx +++ b/PWGUD/TableProducer/upcCandProducerBarrel.cxx @@ -23,7 +23,9 @@ #include #include +#include #include +#include #include #include #include @@ -36,6 +38,7 @@ #include #include #include +#include #include #include @@ -58,7 +61,6 @@ using namespace o2::aod::rctsel; namespace { constexpr int NBCsPerOrbit = o2::constants::lhc::LHCMaxBunches; -constexpr int MaxStoredDistance = 15; constexpr float MaxFITTime = 30.f; constexpr int NoSignalMCRef = -1; @@ -80,7 +82,7 @@ using SignalRef = std::pair; using DistanceMap = std::vector>; DistanceMap buildMinimumDistanceMap(std::vector const& tfIDs, - int64_t bcSOR, int64_t nBCsPerTF, + int64_t bcSOR, int32_t nBCsPerTF, std::vector& activeBCs, uint32_t maxDistance, bool skipFirstExactMatch = false) { @@ -95,7 +97,7 @@ DistanceMap buildMinimumDistanceMap(std::vector const& tfIDs, size_t nextIndex = 0; for (size_t iTF = 0; iTF < tfIDs.size(); ++iTF) { const int64_t tfStartBC = bcSOR + tfIDs[iTF] * nBCsPerTF; - for (int64_t bcInTF = 0; bcInTF < nBCsPerTF; ++bcInTF) { + for (int64_t bcInTF = 0; bcInTF < nBCsPerTF; ++bcInTF) { // TODO why int64_t const int64_t currentBC = tfStartBC + bcInTF; while (nextIndex < activeBCs.size() && activeBCs[nextIndex] < currentBC) { ++nextIndex; @@ -118,26 +120,13 @@ DistanceMap buildMinimumDistanceMap(std::vector const& tfIDs, return distances; } -void fillVetoHistograms(std::shared_ptr const& hVetoT00, std::shared_ptr const& hVetoTV0, - std::shared_ptr const& hVetoTVD, int bcInOrbit, +void fillVetoHistograms(std::shared_ptr const& hVetoT00, std::shared_ptr const& hVetoTV0, + std::shared_ptr const& hVetoTVD, int bcInOrbit, int tfID, uint32_t distanceFT0, uint32_t distanceFV0, uint32_t distanceFDD) { - hVetoT00->Fill(bcInOrbit, -1); - hVetoTV0->Fill(bcInOrbit, -1); - hVetoTVD->Fill(bcInOrbit, -1); - for (uint32_t threshold = 0; threshold <= MaxStoredDistance; ++threshold) { - if (distanceFT0 <= threshold) { - continue; - } - hVetoT00->Fill(bcInOrbit, threshold); - if (distanceFV0 <= threshold) { - continue; - } - hVetoTV0->Fill(bcInOrbit, threshold); - if (distanceFDD > threshold) { - hVetoTVD->Fill(bcInOrbit, threshold); - } - } + hVetoT00->Fill(bcInOrbit, tfID, distanceFT0); + hVetoTV0->Fill(bcInOrbit, tfID, std::min(distanceFT0, distanceFV0)); + hVetoTVD->Fill(bcInOrbit, tfID, std::min({distanceFT0, distanceFV0, distanceFDD})); } } // namespace @@ -147,6 +136,8 @@ DECLARE_SOA_COLUMN(RunNumber, runNumber, int32_t); DECLARE_SOA_COLUMN(GlobalBC, globalBC, uint64_t); DECLARE_SOA_COLUMN(TfId, tfId, uint32_t); DECLARE_SOA_COLUMN(Timestamp, timestamp, uint64_t); +DECLARE_SOA_COLUMN(Flags, flags, int8_t); +DECLARE_SOA_COLUMN(OccupancyFT0C, occupancyFT0C, float); DECLARE_SOA_COLUMN(PosX, posX, float); DECLARE_SOA_COLUMN(PosY, posY, float); DECLARE_SOA_COLUMN(PosZ, posZ, float); @@ -184,6 +175,8 @@ DECLARE_SOA_TABLE(UPCBarrelCands, "AOD", "UPCBARRELCANDS", upc_cand_prod_bar::GlobalBC, upc_cand_prod_bar::TfId, upc_cand_prod_bar::Timestamp, + upc_cand_prod_bar::Flags, + upc_cand_prod_bar::OccupancyFT0C, upc_cand_prod_bar::PosX, upc_cand_prod_bar::PosY, upc_cand_prod_bar::PosZ, @@ -243,12 +236,19 @@ struct UpcCandProducerBarrel { int cachedRunNumber = -1; int64_t bcSOR = 0; - int64_t nBCsPerTF = 0; + int32_t nBCsPerTF = 0; int configuredOrbitsPerTF = -1; int configuredTFStartBorder = -1; int configuredTFEndBorder = -1; + int configuredITSROFStartBorder = -1; + int configuredITSROFEndBorder = -1; int tfStartBorder = 0; int tfEndBorder = 0; + int rofOffset = 0; + int rofLength = 0; + int nRofPerTF = 0; + int rofStartBorderMargin = 0; + int rofEndBorderMargin = 0; std::bitset collidingBCs; Configurable nTracks{"nTracks", 2, "Required N tracks in selected collisions"}; @@ -258,8 +258,19 @@ struct UpcCandProducerBarrel { Configurable vetoFT0{"vetoFT0", 1, "FT0 veto: 0=off, 1=on"}; Configurable vetoFV0{"vetoFV0", 0, "FV0 veto: 0=off, 1=on"}; Configurable vetoFDD{"vetoFDD", 0, "FDD veto: 0=off, 1=on"}; + Configurable useTCE{"useTCE", 0, "Use TCE for lumi accounting: 0=off, 1=on"}; Configurable signalSourceId{"signalSourceId", 1, "MC source ID stored as embedded signal"}; + Configurable minAmpFV0A{"minAmpFV0A", 150.f, "Minimum amplitude for FV0A veto"}; + Configurable minAmpFT0A{"minAmpFT0A", 150.f, "Minimum amplitude for FT0A veto"}; + Configurable minAmpFT0C{"minAmpFT0C", 50.f, "Minimum amplitude for FT0C veto"}; + Configurable minAmpFDDA{"minAmpFDDA", 1000.f, "Minimum amplitude for FDDA veto"}; + Configurable minAmpFDDC{"minAmpFDDC", 1000.f, "Minimum amplitude for FDDC veto"}; + + Configurable maxStoredDistance{"maxStoredDistance", 32, "Maximum stored veto distance"}; + Configurable nBinsOccupancy{"nBinsOccupancy", 25, "Number of occupancy bins"}; + Configurable nBinsTF{"nBinsTF", 100, "Number of TF bins"}; + using CollisionsWithSels = soa::Join; using BCsWithSels = soa::Join; using TracksWithPID = soa::Join -40 && dt < -5) // o2-linter: disable=magic-number (collisions in the past) + wOccup = 1. / 1225 * (dt + 40) * (dt + 40); // o2-linter: disable=magic-number (parabolic weight for past collisions) + else if (dt >= -5 && dt < 15) // o2-linter: disable=magic-number (collisions near a given one) + wOccup = 1; + else if (dt >= 15 && dt < 40) // o2-linter: disable=magic-number (collisions from the future) + wOccup = -0.4 / 25 * dt + 1.24; // o2-linter: disable=magic-number (linear weight for future collisions) + else if (dt >= 40 && dt < 100) // o2-linter: disable=magic-number (collisions from the distant future) + wOccup = -0.4 / 60 * dt + 0.6 + 0.8 / 3; // o2-linter: disable=magic-number (linear weight for distant future collisions) + return wOccup; + } + + bool noRofBorder(int64_t bc) + { + int bcInRof = (bc + rofLength - rofOffset) % rofLength; + return bcInRof > rofStartBorderMargin && bcInRof < rofLength - rofEndBorderMargin; + } + + int getRofInTF(int64_t bc) + { + return (bc + rofLength - rofOffset) / rofLength % nRofPerTF; + } + void init(InitContext& initContext) { ccdb->setLocalObjectValidityChecking(); @@ -287,15 +323,49 @@ struct UpcCandProducerBarrel { inheritEventSelectionOption("bcselOpts.NumberOfOrbitsPerTF", configuredOrbitsPerTF); inheritEventSelectionOption("bcselOpts.TimeFrameStartBorderMargin", configuredTFStartBorder); inheritEventSelectionOption("bcselOpts.TimeFrameEndBorderMargin", configuredTFEndBorder); - registry.add("hProcessedTFs", "Processed TFs;;TFs", HistType::kTH1D, {{1, 0., 1.}}); - registry.add("hTVX", "TVX counts per run;run;TVX BCs", HistType::kTH1D, {{NBCsPerOrbit, 0., NBCsPerOrbit}}); - registry.add("hTVXRCT", "TVX counts per run;run;TVX BCs", HistType::kTH1D, {{NBCsPerOrbit, 0., NBCsPerOrbit}}); - registry.add("hVetoT00", ";BC in orbit;veto range (BC);colliding BC count, no RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); - registry.add("hVetoTV0", ";BC in orbit;veto range (BC);colliding BC count, no RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); - registry.add("hVetoTVD", ";BC in orbit;veto range (BC);colliding BC count, no RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); - registry.add("hVetoT00RCT", ";BC in orbit;veto range (BC);colliding BC count, RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); - registry.add("hVetoTV0RCT", ";BC in orbit;veto range (BC);colliding BC count, RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); - registry.add("hVetoTVDRCT", ";BC in orbit;veto range (BC);colliding BC count, RCT mask", HistType::kTH2D, {{NBCsPerOrbit, 0., NBCsPerOrbit}, {17, -1.5, 15.5}}); + inheritEventSelectionOption("bcselOpts.ITSROFrameStartBorderMargin", configuredITSROFStartBorder); + inheritEventSelectionOption("bcselOpts.ITSROFrameEndBorderMargin", configuredITSROFEndBorder); + + const AxisSpec axisBCs{NBCsPerOrbit, 0., NBCsPerOrbit, "BC in orbit"}; + const AxisSpec axisVeto{maxStoredDistance + 1, -0.5, maxStoredDistance + 0.5, "veto range in BCs"}; + const AxisSpec axisOccupancy{nBinsOccupancy, 0., 250000., "FT0-based occupancy"}; + const AxisSpec axisTFs{nBinsTF, 0., 4.e7, "TFs"}; // 100 bins -> 400000 TFs per bin (1 bin ~ 1140 seconds for 32 orbits per TF) + + registry.add("hTVX", "TVX counts", HistType::kTH2D, {axisBCs, axisTFs}); + registry.add("hTVXRCT", "TVX counts after RCT selection", HistType::kTH2D, {axisBCs, axisTFs}); + registry.add("hTVXRCTnoPileupROF", "TVX counts after RCT and ROF pileup selection", HistType::kTH2D, {axisBCs, axisTFs}); + registry.add("hTVXupc", "TVX counts in upc ROFs", HistType::kTH2D, {axisBCs, axisTFs}); + registry.add("hTVXRCTupc", "TVX counts after RCT selection in upc ROFs", HistType::kTH2D, {axisBCs, axisTFs}); + registry.add("hTVXRCTnoPileupROFupc", "TVX counts after RCT and ROF pileup selection in upc ROFs", HistType::kTH2D, {axisBCs, axisTFs}); + + registry.add("hVetoT00", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVD", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoT00RCT", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0RCT", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVDRCT", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoT00RCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0RCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVDRCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + + registry.add("hOccupancyT00RCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + registry.add("hOccupancyTV0RCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + registry.add("hOccupancyTVDRCTnoPileupROF", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + + registry.add("hVetoT00upc", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0upc", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVDupc", "no RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoT00RCTupc", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0RCTupc", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVDRCTupc", "after RCT selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoT00RCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTV0RCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + registry.add("hVetoTVDRCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisVeto}); + + registry.add("hOccupancyT00RCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + registry.add("hOccupancyTV0RCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + registry.add("hOccupancyTVDRCTnoPileupROFupc", "after RCT and ROF pileup selection", HistType::kTH3D, {axisBCs, axisTFs, axisOccupancy}); + LOGF(info, "Veto config: window=%d FT0=%d FV0=%d FDD=%d", vetoBCWindow.value, vetoFT0.value, vetoFV0.value, vetoFDD.value); } @@ -307,7 +377,8 @@ struct UpcCandProducerBarrel { const auto runInfo = o2::parameters::AggregatedRunInfo::buildAggregatedRunInfo(ccdb->instance(), runNumber); const int64_t orbitsPerTF = configuredOrbitsPerTF < 0 ? runInfo.orbitsPerTF : configuredOrbitsPerTF; const bool needsEventSelectionParams = configuredTFStartBorder < 0 || configuredTFEndBorder < 0; - const auto* eventSelectionParams = needsEventSelectionParams ? ccdb->getForTimeStamp("EventSelection/EventSelectionParams", runInfo.sor / 2 + runInfo.eor / 2) : nullptr; + const auto ts = (runInfo.sor + runInfo.eor) / 2; + const auto* eventSelectionParams = needsEventSelectionParams ? ccdb->getForTimeStamp("EventSelection/EventSelectionParams", ts) : nullptr; if (orbitsPerTF <= 0 || !runInfo.grpLHC || (needsEventSelectionParams && !eventSelectionParams)) { LOGF(fatal, "Incomplete run information for run %d", runNumber); } @@ -316,6 +387,12 @@ struct UpcCandProducerBarrel { tfStartBorder = configuredTFStartBorder < 0 ? eventSelectionParams->fTimeFrameStartBorderMargin : configuredTFStartBorder; tfEndBorder = configuredTFEndBorder < 0 ? eventSelectionParams->fTimeFrameEndBorderMargin : configuredTFEndBorder; collidingBCs = runInfo.grpLHC->getBunchFilling().getBCPattern(); + + // get ITS ROF parameters + auto alppar = ccdb->getForTimeStamp>("ITS/Config/AlpideParam", ts); + rofOffset = alppar->roFrameBiasInBC; + rofLength = alppar->roFrameLengthInBC; + nRofPerTF = nBCsPerTF / rofLength; cachedRunNumber = runNumber; } @@ -340,7 +417,7 @@ struct UpcCandProducerBarrel { } localTFIndex.push_back(tfIDs.size() - 1); } - registry.fill(HIST("hProcessedTFs"), 0.5, static_cast(tfIDs.size())); + LOGP(info, "Starting new DF. Number of TFs = {}", tfIDs.size()); std::unordered_map signalParticleRefs; if constexpr (isMC) { @@ -356,7 +433,7 @@ struct UpcCandProducerBarrel { for (const auto& mcCollision : *mcCollisions) { bcsWithMCCollision.push_back(static_cast(mcCollision.template bc_as().globalBC())); } - const auto minDistanceMC = buildMinimumDistanceMap(tfIDs, bcSOR, nBCsPerTF, bcsWithMCCollision, MaxStoredDistance, true); + const auto minDistanceMC = buildMinimumDistanceMap(tfIDs, bcSOR, nBCsPerTF, bcsWithMCCollision, maxStoredDistance, true); for (const auto& mcCollision : *mcCollisions) { const auto bc = mcCollision.template bc_as(); @@ -416,78 +493,246 @@ struct UpcCandProducerBarrel { } } - auto hTVX = registry.get(HIST("hTVX")); - auto hTVXRCT = registry.get(HIST("hTVXRCT")); + // collect number of all collisions and upc-mode collisions per ROF + LOGP(debug, "Collect number of all collisions and upc-mode collisions per ROF"); + std::vector> vMultTracks(tfIDs.size(), std::map()); + std::vector> vCollisionsPerROF(tfIDs.size(), std::vector(nRofPerTF, 0)); + std::vector> vUpcModePerROF(tfIDs.size(), std::vector(nRofPerTF, 0)); + for (const auto& collision : collisions) { + auto bc = collision.template bc_as(); + const int64_t globalBC = bc.globalBC(); + const int32_t bcInTF = (globalBC - bcSOR) % nBCsPerTF; + const auto iTF = localTFIndex[bc.globalIndex()]; + const int32_t rofInTF = getRofInTF(bcInTF); + vCollisionsPerROF[iTF][rofInTF]++; + const bool upcMode = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) != 0; + if (upcMode) { + vUpcModePerROF[iTF][rofInTF]++; + } + } + + auto hTVX = registry.get(HIST("hTVX")); + auto hTVXRCT = registry.get(HIST("hTVXRCT")); + auto hTVXRCTnoPileupROF = registry.get(HIST("hTVXRCTnoPileupROF")); + + auto hTVXupc = registry.get(HIST("hTVXupc")); + auto hTVXRCTupc = registry.get(HIST("hTVXRCTupc")); + auto hTVXRCTnoPileupROFupc = registry.get(HIST("hTVXRCTnoPileupROFupc")); std::vector bcsWithFT0, bcsWithFV0, bcsWithFDD; bcsWithFT0.reserve(ft0s.size()); bcsWithFV0.reserve(fv0s.size()); bcsWithFDD.reserve(fdds.size()); + LOGP(debug, "Starting lumi accounting"); + // lumi accounting + std::vector> vMultFT0C(tfIDs.size(), std::map()); for (const auto& ft0 : ft0s) { const auto bc = ft0.bc_as(); - const auto gbc = bc.globalBC(); - const auto bcInOrbit = gbc % NBCsPerOrbit; - if (ft0.timeA() < MaxFITTime || ft0.timeC() < MaxFITTime) { - bcsWithFT0.push_back(gbc); + const auto globalBC = bc.globalBC(); + const auto bcInOrbit = globalBC % NBCsPerOrbit; + const int32_t bcInTF = (globalBC - bcSOR) % nBCsPerTF; + const auto iTF = localTFIndex[bc.globalIndex()]; + const auto tfID = tfIDs[iTF]; + const int32_t rofInTF = getRofInTF(bcInTF); + bool pileupInROF = vCollisionsPerROF[iTF][rofInTF] > 1; + bool upcMode = vCollisionsPerROF[iTF][rofInTF] == 0 || vUpcModePerROF[iTF][rofInTF] > 0; + + // Multiplicity for FT0C-based occupancy estimates + if (ft0.timeC() < MaxFITTime && ft0.sumAmpC() >= minAmpFT0C) { + vMultFT0C[iTF][bcInTF] = ft0.sumAmpC(); + } + + if (!collidingBCs[bcInOrbit] || !TESTBIT(ft0.triggerMask(), o2::fit::Triggers::bitVertex)) { + continue; } - if (!collidingBCs[bcInOrbit] || !(ft0.triggerMask() & (1U << o2::fit::Triggers::bitVertex))) { + + if (useTCE && !TESTBIT(ft0.triggerMask(), o2::fit::Triggers::bitCen)) { continue; } - hTVX->Fill(bcInOrbit); + + hTVX->Fill(bcInOrbit, tfID); + if (upcMode) { + hTVXupc->Fill(bcInOrbit, tfID); + } + if (!bc.selection_bit(aod::evsel::kNoTimeFrameBorder) || !rctChecker(bc)) { continue; } - hTVXRCT->Fill(bcInOrbit); + + hTVXRCT->Fill(bcInOrbit, tfID); + if (upcMode) { + hTVXRCTupc->Fill(bcInOrbit, tfID); + } + + if (pileupInROF) { + continue; + } + + hTVXRCTnoPileupROF->Fill(bcInOrbit, tfID); + if (upcMode) { + hTVXRCTnoPileupROFupc->Fill(bcInOrbit, tfID); + } + } + + LOGP(debug, "Filling FT0 activity arrays"); + for (const auto& ft0 : ft0s) { + const float sumAmpFT0A = ft0.sumAmpA(); + const float sumAmpFT0C = ft0.sumAmpC(); + if ((ft0.timeA() < MaxFITTime && sumAmpFT0A >= minAmpFT0A) || (ft0.timeC() < MaxFITTime && sumAmpFT0C >= minAmpFT0C)) { + bcsWithFT0.push_back(ft0.bc_as().globalBC()); + } } + + LOGP(debug, "Filling FV0 activity arrays"); for (const auto& fv0 : fv0s) { - if (fv0.time() < MaxFITTime) { + const float sumAmpFV0A = std::accumulate(fv0.amplitude().begin(), fv0.amplitude().end(), 0.f); + if (fv0.time() < MaxFITTime && sumAmpFV0A >= minAmpFV0A) { bcsWithFV0.push_back(fv0.bc_as().globalBC()); } } + + LOGP(debug, "Filling FDD activity arrays"); for (const auto& fdd : fdds) { - if (fdd.timeA() < MaxFITTime || fdd.timeC() < MaxFITTime) { + const float sumAmpFDDA = std::accumulate(fdd.chargeA(), fdd.chargeA() + 8, 0.f); + const float sumAmpFDDC = std::accumulate(fdd.chargeC(), fdd.chargeC() + 8, 0.f); + if ((fdd.timeA() < MaxFITTime && sumAmpFDDA >= minAmpFDDA) || (fdd.timeC() < MaxFITTime && sumAmpFDDC >= minAmpFDDC)) { bcsWithFDD.push_back(fdd.bc_as().globalBC()); } } - const uint32_t scanDistance = std::max(MaxStoredDistance, vetoBCWindow.value); + LOGP(debug, "Starting occupancy calculation"); + std::vector> vOccupancyFT0C(tfIDs.size(), std::vector(nBCsPerTF, 0.f)); + int deltaBcPast = 1600; + int deltaBcFuture = 4000; + // TODO: calculate weights only once per task + std::vector weights(deltaBcPast + deltaBcFuture + 1); + for (int d = -deltaBcPast; d <= deltaBcFuture; ++d) { + weights[d + deltaBcPast] = calcWeightForOccupancy(d * o2::constants::lhc::LHCBunchSpacingMUS); + } + + for (size_t iTF = 0; iTF < tfIDs.size(); ++iTF) { + for (const auto& [bc, multFT0C] : vMultFT0C[iTF]) { + LOGP(debug, "iTF={} bc={} mult={}", iTF, bc, multFT0C); + const int bcInTFMin = std::max(bc - deltaBcFuture, 0); + const int bcInTFMax = std::min(bc + deltaBcPast, nBCsPerTF - 1); + for (int bcInTF = bcInTFMin; bcInTF <= bcInTFMax; ++bcInTF) { + vOccupancyFT0C[iTF][bcInTF] += weights[bc - bcInTF + deltaBcPast] * multFT0C; + } + } + } + + LOGP(debug, "Building minimum distance maps"); + const uint32_t scanDistance = std::max(maxStoredDistance.value, vetoBCWindow.value); const DistanceMap minDistanceFT0 = buildMinimumDistanceMap(tfIDs, bcSOR, nBCsPerTF, bcsWithFT0, scanDistance); const DistanceMap minDistanceFV0 = buildMinimumDistanceMap(tfIDs, bcSOR, nBCsPerTF, bcsWithFV0, scanDistance); const DistanceMap minDistanceFDD = buildMinimumDistanceMap(tfIDs, bcSOR, nBCsPerTF, bcsWithFDD, scanDistance); - auto hVetoT00 = registry.get(HIST("hVetoT00")); - auto hVetoTV0 = registry.get(HIST("hVetoTV0")); - auto hVetoTVD = registry.get(HIST("hVetoTVD")); - auto hVetoT00RCT = registry.get(HIST("hVetoT00RCT")); - auto hVetoTV0RCT = registry.get(HIST("hVetoTV0RCT")); - auto hVetoTVDRCT = registry.get(HIST("hVetoTVDRCT")); + auto hVetoT00 = registry.get(HIST("hVetoT00")); + auto hVetoTV0 = registry.get(HIST("hVetoTV0")); + auto hVetoTVD = registry.get(HIST("hVetoTVD")); + auto hVetoT00RCT = registry.get(HIST("hVetoT00RCT")); + auto hVetoTV0RCT = registry.get(HIST("hVetoTV0RCT")); + auto hVetoTVDRCT = registry.get(HIST("hVetoTVDRCT")); + auto hVetoT00RCTnoPileupROF = registry.get(HIST("hVetoT00RCTnoPileupROF")); + auto hVetoTV0RCTnoPileupROF = registry.get(HIST("hVetoTV0RCTnoPileupROF")); + auto hVetoTVDRCTnoPileupROF = registry.get(HIST("hVetoTVDRCTnoPileupROF")); + + auto hVetoT00upc = registry.get(HIST("hVetoT00upc")); + auto hVetoTV0upc = registry.get(HIST("hVetoTV0upc")); + auto hVetoTVDupc = registry.get(HIST("hVetoTVDupc")); + auto hVetoT00RCTupc = registry.get(HIST("hVetoT00RCTupc")); + auto hVetoTV0RCTupc = registry.get(HIST("hVetoTV0RCTupc")); + auto hVetoTVDRCTupc = registry.get(HIST("hVetoTVDRCTupc")); + auto hVetoT00RCTnoPileupROFupc = registry.get(HIST("hVetoT00RCTnoPileupROFupc")); + auto hVetoTV0RCTnoPileupROFupc = registry.get(HIST("hVetoTV0RCTnoPileupROFupc")); + auto hVetoTVDRCTnoPileupROFupc = registry.get(HIST("hVetoTVDRCTnoPileupROFupc")); + + auto hOccupancyT00RCTnoPileupROF = registry.get(HIST("hOccupancyT00RCTnoPileupROF")); + auto hOccupancyTV0RCTnoPileupROF = registry.get(HIST("hOccupancyTV0RCTnoPileupROF")); + auto hOccupancyTVDRCTnoPileupROF = registry.get(HIST("hOccupancyTVDRCTnoPileupROF")); + + auto hOccupancyT00RCTnoPileupROFupc = registry.get(HIST("hOccupancyT00RCTnoPileupROFupc")); + auto hOccupancyTV0RCTnoPileupROFupc = registry.get(HIST("hOccupancyTV0RCTnoPileupROFupc")); + auto hOccupancyTVDRCTnoPileupROFupc = registry.get(HIST("hOccupancyTVDRCTnoPileupROFupc")); + + LOGP(debug, "Filling veto efficiency histograms"); for (size_t iTF = 0; iTF < tfIDs.size(); ++iTF) { const int64_t tfStartBC = bcSOR + static_cast(tfIDs[iTF]) * nBCsPerTF; - for (int64_t bcInTF = tfStartBorder + 1; bcInTF < nBCsPerTF - tfEndBorder; ++bcInTF) { + const auto tfID = tfIDs[iTF]; + for (int64_t bcInTF = tfStartBorder + 1; bcInTF < nBCsPerTF - tfEndBorder; ++bcInTF) { // TODO why int64_t? const int64_t globalBC = tfStartBC + bcInTF; const int bcInOrbit = globalBC % NBCsPerOrbit; - if (!collidingBCs[bcInOrbit]) { - continue; - } + const int32_t rofInTF = getRofInTF(bcInTF); + bool upcMode = vCollisionsPerROF[iTF][rofInTF] == 0 || vUpcModePerROF[iTF][rofInTF] > 0; + bool pileupInROF = vCollisionsPerROF[iTF][rofInTF] > 1; const auto distanceFT0 = minDistanceFT0[iTF][bcInTF]; const auto distanceFV0 = minDistanceFV0[iTF][bcInTF]; const auto distanceFDD = minDistanceFDD[iTF][bcInTF]; - fillVetoHistograms(hVetoT00, hVetoTV0, hVetoTVD, bcInOrbit, distanceFT0, distanceFV0, distanceFDD); - if (tfPassesRCT[iTF] != 0u) { - fillVetoHistograms(hVetoT00RCT, hVetoTV0RCT, hVetoTVDRCT, bcInOrbit, distanceFT0, distanceFV0, distanceFDD); + + fillVetoHistograms(hVetoT00, hVetoTV0, hVetoTVD, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + if (upcMode) { + fillVetoHistograms(hVetoT00upc, hVetoTV0upc, hVetoTVDupc, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + } + + if (tfPassesRCT[iTF] == 0u) { + continue; + } + + fillVetoHistograms(hVetoT00RCT, hVetoTV0RCT, hVetoTVDRCT, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + if (upcMode) { + fillVetoHistograms(hVetoT00RCTupc, hVetoTV0RCTupc, hVetoTVDRCTupc, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + } + + if (pileupInROF) { + continue; + } + + fillVetoHistograms(hVetoT00RCTnoPileupROF, hVetoTV0RCTnoPileupROF, hVetoTVDRCTnoPileupROF, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + if (upcMode) { + fillVetoHistograms(hVetoT00RCTnoPileupROFupc, hVetoTV0RCTnoPileupROFupc, hVetoTVDRCTnoPileupROFupc, bcInOrbit, tfID, distanceFT0, distanceFV0, distanceFDD); + } + + if (distanceFT0 == 0) + continue; + float occupancy = vOccupancyFT0C[iTF][bcInTF]; + hOccupancyT00RCTnoPileupROF->Fill(bcInOrbit, tfID, occupancy); + if (upcMode) { + hOccupancyT00RCTnoPileupROFupc->Fill(bcInOrbit, tfID, occupancy); + } + + if (distanceFV0 == 0) + continue; + hOccupancyTV0RCTnoPileupROF->Fill(bcInOrbit, tfID, occupancy); + if (upcMode) { + hOccupancyTV0RCTnoPileupROFupc->Fill(bcInOrbit, tfID, occupancy); + } + + if (distanceFDD == 0) + continue; + hOccupancyTVDRCTnoPileupROF->Fill(bcInOrbit, tfID, occupancy); + if (upcMode) { + hOccupancyTVDRCTnoPileupROFupc->Fill(bcInOrbit, tfID, occupancy); } } } + LOGP(debug, "Filling upc trees"); for (const auto& collision : collisions) { if (collision.numContrib() != nTracks || !collision.selection_bit(aod::evsel::kNoTimeFrameBorder) || !rctChecker(collision)) { continue; } auto bc = collision.template bc_as(); - auto gbc = bc.globalBC(); + int64_t globalBC = bc.globalBC(); const auto iTF = localTFIndex[bc.globalIndex()]; - const int64_t bcInTF = (static_cast(gbc) - bcSOR) % nBCsPerTF; + const int64_t bcInTF = (globalBC - bcSOR) % nBCsPerTF; + const int32_t rofInTF = getRofInTF(bcInTF); + const bool upcMode = (collision.flags() & dataformats::Vertex>::Flags::UPCMode) != 0; + const bool noPileupInRof = vCollisionsPerROF[iTF][rofInTF] <= 1; + uint8_t flags = 0; + flags |= upcMode << 0; + flags |= noPileupInRof << 1; + const auto occupancyFT0C = vOccupancyFT0C[iTF][bcInTF]; const auto distFT0 = minDistanceFT0[iTF][bcInTF]; const auto distFV0 = minDistanceFV0[iTF][bcInTF]; const auto distFDD = minDistanceFDD[iTF][bcInTF]; @@ -537,12 +782,13 @@ struct UpcCandProducerBarrel { if (px.size() != static_cast(nTracks)) { continue; } - selectedCandidates(bc.runNumber(), gbc, tfIDs[iTF], bc.timestamp(), collision.posX(), collision.posY(), collision.posZ(), + selectedCandidates(bc.runNumber(), globalBC, tfIDs[iTF], bc.timestamp(), flags, occupancyFT0C, + collision.posX(), collision.posY(), collision.posZ(), px, py, pz, tpcSignal, tpcNSigmaEl, tpcNSigmaPi, tpcNSigmaKa, tpcNSigmaPr, tofNSigmaEl, tofNSigmaPi, tofNSigmaKa, tofNSigmaPr, itsClusterMap, nClusters, sign, - std::min(distFT0, MaxStoredDistance + 1), - std::min(distFV0, MaxStoredDistance + 1), - std::min(distFDD, MaxStoredDistance + 1)); + std::min(distFT0, maxStoredDistance + 1), + std::min(distFV0, maxStoredDistance + 1), + std::min(distFDD, maxStoredDistance + 1)); if constexpr (isMC) { selectedCandidateMCLabels(newMCCollisionIds, newParticleIds); } From a609297907efae21de7a8e35f7d264b8034c0b79 Mon Sep 17 00:00:00 2001 From: Evgeny Kryshen Date: Thu, 29 Oct 2026 19:07:58 +0300 Subject: [PATCH 2/5] Adding map header --- PWGUD/TableProducer/upcCandProducerBarrel.cxx | 1 + 1 file changed, 1 insertion(+) diff --git a/PWGUD/TableProducer/upcCandProducerBarrel.cxx b/PWGUD/TableProducer/upcCandProducerBarrel.cxx index 70170fe2367..a72220b7c30 100644 --- a/PWGUD/TableProducer/upcCandProducerBarrel.cxx +++ b/PWGUD/TableProducer/upcCandProducerBarrel.cxx @@ -53,6 +53,7 @@ #include #include #include +#include using namespace o2; using namespace o2::framework; From f84f28848c2462bdf808ee3105f26e3d12a27c63 Mon Sep 17 00:00:00 2001 From: Evgeny Kryshen Date: Thu, 29 Oct 2026 19:10:24 +0300 Subject: [PATCH 3/5] Megalinter fix --- PWGUD/TableProducer/upcCandProducerBarrel.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGUD/TableProducer/upcCandProducerBarrel.cxx b/PWGUD/TableProducer/upcCandProducerBarrel.cxx index a72220b7c30..c6c439a0c14 100644 --- a/PWGUD/TableProducer/upcCandProducerBarrel.cxx +++ b/PWGUD/TableProducer/upcCandProducerBarrel.cxx @@ -48,12 +48,12 @@ #include #include #include +#include #include #include #include #include #include -#include using namespace o2; using namespace o2::framework; From 400950e3d509c0328e7e173e4cd6a2f053b1bf7c Mon Sep 17 00:00:00 2001 From: Evgeny Kryshen Date: Thu, 29 Oct 2026 19:29:50 +0300 Subject: [PATCH 4/5] Trigger CI From 5d9846b6a9c13ed8277112ee64354602d5ce3e6a Mon Sep 17 00:00:00 2001 From: Evgeny Kryshen Date: Thu, 29 Oct 2026 21:14:25 +0300 Subject: [PATCH 5/5] Trigger CI