Skip to content

[Common] adding median time for occupancy estimation #14202

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
altsybee merged 3 commits into from
Dec 15, 2025
Merged

[Common] adding median time for occupancy estimation #14202

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
19bada7
adding median occupancy column
altsybee Dec 12, 2025
3509a91
adding median time for occupancy
altsybee Dec 12, 2025
b388fe3
adding median time for occupancy (dummy for old evSel task)
altsybee Dec 12, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 3 additions & 1 deletion Common/DataModel/EventSelection.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,8 +63,9 @@
DECLARE_SOA_INDEX_COLUMN_FULL(FoundFV0, foundFV0, int, FV0As, "_foundFV0"); //! FV0 entry index in FV0As table (-1 if doesn't exist)
DECLARE_SOA_INDEX_COLUMN_FULL(FoundFDD, foundFDD, int, FDDs, "_foundFDD"); //! FDD entry index in FDDs table (-1 if doesn't exist)
DECLARE_SOA_INDEX_COLUMN_FULL(FoundZDC, foundZDC, int, Zdcs, "_foundZDC"); //! ZDC entry index in ZDCs table (-1 if doesn't exist)
DECLARE_SOA_COLUMN(NumTracksInTimeRange, trackOccupancyInTimeRange, int); //! Occupancy in specified time interval by a number of tracks from nearby collisions // o2-linter: disable=name/o2-column

Check failure on line 66 in Common/DataModel/EventSelection.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column] 

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.
DECLARE_SOA_COLUMN(SumAmpFT0CInTimeRange, ft0cOccupancyInTimeRange, float); //! Occupancy in specified time interval by a sum of FT0C amplitudes from nearby collisions // o2-linter: disable=name/o2-column

Check failure on line 67 in Common/DataModel/EventSelection.h

View workflow job for this annotation

GitHub Actions / O2 linter

[name/o2-column] 

Use UpperCamelCase for names of O2 columns and matching lowerCamelCase names for their getters.
DECLARE_SOA_COLUMN(OccupancyMedianTime, occupancyMedianTime, float); //! Median time for the occupancy in specified time interval (the delta-time wrt to the collision at which occupancy crosses 1/2 of its full value)
} // namespace evsel

// bc-joinable event selection decisions
Expand All @@ -85,7 +86,8 @@
evsel::FoundFDDId,
evsel::FoundZDCId,
evsel::NumTracksInTimeRange,
evsel::SumAmpFT0CInTimeRange);
evsel::SumAmpFT0CInTimeRange,
evsel::OccupancyMedianTime);
using EvSel = EvSels::iterator;
} // namespace o2::aod

Expand Down
6 changes: 3 additions & 3 deletions Common/TableProducer/eventSelection.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -694,7 +694,7 @@
}
}

evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, 0, 0);
evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, 0, 0, 0);
}
PROCESS_SWITCH(EventSelectionTask, processRun2, "Process Run2 event selection", true);

Expand Down Expand Up @@ -750,7 +750,7 @@
int32_t foundFDD = bc.foundFDDId();
int32_t foundZDC = bc.foundZDCId();
uint32_t rct = 0;
evsel(bc.alias_raw(), bc.selection_raw(), rct, kFALSE, kFALSE, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, -1, -1);
evsel(bc.alias_raw(), bc.selection_raw(), rct, kFALSE, kFALSE, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, -1, -1, -1);
}
return;
}
Expand Down Expand Up @@ -1166,7 +1166,7 @@
}

evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC,
vNumTracksITS567inFullTimeWin[colIndex], vSumAmpFT0CinFullTimeWin[colIndex]);
vNumTracksITS567inFullTimeWin[colIndex], vSumAmpFT0CinFullTimeWin[colIndex], 0);
}
}

Expand Down Expand Up @@ -1373,15 +1373,15 @@
double perBcRateTCE = static_cast<double>(mCounterTCE[i + 1] - mCounterTCE[i]) / nOrbits / nCollidingBCs;
double perBcRateZNC = static_cast<double>(mCounterZNC[i + 1] - mCounterZNC[i]) / nOrbits / nCollidingBCs;
double perBcRateZEM = static_cast<double>(mCounterZEM[i + 1] - mCounterZEM[i]) / nOrbits / nCollidingBCs;
double muTVX = (perBcRateTVX < 1 && perBcRateTVX > 1e-10) ? -std::log(1 - perBcRateTVX) : 0;

Check failure on line 1376 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
double muTCE = (perBcRateTCE < 1 && perBcRateTCE > 1e-10) ? -std::log(1 - perBcRateTCE) : 0;

Check failure on line 1377 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
double muZNC = (perBcRateZNC < 1 && perBcRateZNC > 1e-10) ? -std::log(1 - perBcRateZNC) : 0;

Check failure on line 1378 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
double muZEM = (perBcRateZEM < 1 && perBcRateZEM > 1e-10) ? -std::log(1 - perBcRateZEM) : 0;

Check failure on line 1379 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
LOGP(debug, "orbit={} muTVX={} muTCE={} muZNC={} muZEM={}", mOrbits[i], muTVX, muTCE, muZNC, muZEM);
mPileupCorrectionTVX.push_back(muTVX > 1e-10 ? muTVX / (1 - std::exp(-muTVX)) : 1);

Check failure on line 1381 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
mPileupCorrectionTCE.push_back(muTCE > 1e-10 ? muTCE / (1 - std::exp(-muTCE)) : 1);

Check failure on line 1382 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
mPileupCorrectionZNC.push_back(muZNC > 1e-10 ? muZNC / (1 - std::exp(-muZNC)) : 1);

Check failure on line 1383 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
mPileupCorrectionZEM.push_back(muZEM > 1e-10 ? muZEM / (1 - std::exp(-muZEM)) : 1);

Check failure on line 1384 in Common/TableProducer/eventSelection.cxx

View workflow job for this annotation

GitHub Actions / O2 linter

[magic-number] 

Avoid magic numbers in expressions. Assign the value to a clearly named variable or constant.
}
// filling last orbit range using previous orbit range
mPileupCorrectionTVX.push_back(mPileupCorrectionTVX.back());
Expand Down
61 changes: 41 additions & 20 deletions Common/Tools/EventSelectionModule.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -96,7 +96,7 @@ struct bcselConfigurables : o2::framework::ConfigurableGroup {
o2::framework::Configurable<int> confTimeFrameStartBorderMargin{"TimeFrameStartBorderMargin", -1, "Number of bcs to cut at the start of the Time Frame. Take from CCDB if -1"}; // o2-linter: disable=name/configurable (temporary fix)
o2::framework::Configurable<int> confTimeFrameEndBorderMargin{"TimeFrameEndBorderMargin", -1, "Number of bcs to cut at the end of the Time Frame. Take from CCDB if -1"}; // o2-linter: disable=name/configurable (temporary fix)
o2::framework::Configurable<bool> confCheckRunDurationLimits{"checkRunDurationLimits", false, "Check if the BCs are within the run duration limits"}; // o2-linter: disable=name/configurable (temporary fix)
o2::framework::Configurable<std::vector<int>> maxInactiveChipsPerLayer{"maxInactiveChipsPerLayer", {8, 8, 8, 111, 111, 195, 195}, "Maximum allowed number of inactive ITS chips per layer"};
o2::framework::Configurable<std::vector<int>> confMaxInactiveChipsPerLayer{"maxInactiveChipsPerLayer", {8, 8, 8, 111, 111, 195, 195}, "Maximum allowed number of inactive ITS chips per layer"};
o2::framework::Configurable<int> confNumberOfOrbitsPerTF{"NumberOfOrbitsPerTF", -1, "Number of orbits per Time Frame. Take from CCDB if -1"}; // o2-linter: disable=name/configurable (temporary fix)
};

Expand All @@ -121,6 +121,8 @@ struct evselConfigurables : o2::framework::ConfigurableGroup {
o2::framework::Configurable<bool> confUseWeightsForOccupancyVariable{"UseWeightsForOccupancyEstimator", 1, "Use or not the delta-time weights for the occupancy estimator"}; // o2-linter: disable=name/configurable (temporary fix)
o2::framework::Configurable<int> confNumberOfOrbitsPerTF{"NumberOfOrbitsPerTF", -1, "Number of orbits per Time Frame. Take from CCDB if -1"}; // o2-linter: disable=name/configurable (temporary fix)

// o2::framework::Configurable<std::vector<float>> confTimeIntervalsForGranularOccupancy{"timeIntervalsForGranularOccupancy", {-100, -60, -30, -10, 0, 10, 30, 60, 100}, "Delta-time intervals (wrt given collision) for which we store occupancy"};

// configurables for light-ion event selection
o2::framework::Configurable<float> confLightIonsNsigmaOnVzDiff{"VzDiffNsigma", 3.0, "+/- nSigma on vZ difference by FT0 and by tracks"}; // o2-linter: disable=name/configurable (temporary fix)
o2::framework::Configurable<float> confLightIonsMarginVzDiff{"VzDiffMargin", 0.2, "margin for +/- nSigma cut on vZ difference by FT0 and by tracks, cm"}; // o2-linter: disable=name/configurable (temporary fix)
Expand Down Expand Up @@ -561,14 +563,14 @@ class BcSelectionModule
LOGP(debug, "orbit: {}, previous orbit: {}, next orbit: {} ", orbit, prevOrbitForInactiveChips, nextOrbitForInactiveChips);
LOGP(debug, "next inactive chips: {} {} {} {} {} {} {}", vNextInactiveChips[0], vNextInactiveChips[1], vNextInactiveChips[2], vNextInactiveChips[3], vNextInactiveChips[4], vNextInactiveChips[5], vNextInactiveChips[6]);
LOGP(debug, "prev inactive chips: {} {} {} {} {} {} {}", vPrevInactiveChips[0], vPrevInactiveChips[1], vPrevInactiveChips[2], vPrevInactiveChips[3], vPrevInactiveChips[4], vPrevInactiveChips[5], vPrevInactiveChips[6]);
isGoodITSLayer3 = vPrevInactiveChips[3] <= bcselOpts.maxInactiveChipsPerLayer->at(3) && vNextInactiveChips[3] <= bcselOpts.maxInactiveChipsPerLayer->at(3);
isGoodITSLayer3 = vPrevInactiveChips[3] <= bcselOpts.confMaxInactiveChipsPerLayer->at(3) && vNextInactiveChips[3] <= bcselOpts.confMaxInactiveChipsPerLayer->at(3);
isGoodITSLayer0123 = true;
for (int i = 0; i < 4; i++) { // o2-linter: disable=magic-number (counting first 4 ITS layers)
isGoodITSLayer0123 &= vPrevInactiveChips[i] <= bcselOpts.maxInactiveChipsPerLayer->at(i) && vNextInactiveChips[i] <= bcselOpts.maxInactiveChipsPerLayer->at(i);
isGoodITSLayer0123 &= vPrevInactiveChips[i] <= bcselOpts.confMaxInactiveChipsPerLayer->at(i) && vNextInactiveChips[i] <= bcselOpts.confMaxInactiveChipsPerLayer->at(i);
}
isGoodITSLayersAll = true;
for (int i = 0; i < o2::itsmft::ChipMappingITS::NLayers; i++) {
isGoodITSLayersAll &= vPrevInactiveChips[i] <= bcselOpts.maxInactiveChipsPerLayer->at(i) && vNextInactiveChips[i] <= bcselOpts.maxInactiveChipsPerLayer->at(i);
isGoodITSLayersAll &= vPrevInactiveChips[i] <= bcselOpts.confMaxInactiveChipsPerLayer->at(i) && vNextInactiveChips[i] <= bcselOpts.confMaxInactiveChipsPerLayer->at(i);
}
}

Expand Down Expand Up @@ -875,7 +877,7 @@ class EventSelectionModule
}
}

evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, 0, 0);
evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, 0, 0, 0);
}
} // end processRun2

Expand Down Expand Up @@ -924,7 +926,7 @@ class EventSelectionModule
int32_t foundFDD = bcselbuffer[bc.globalIndex()].foundFDDId;
int32_t foundZDC = bcselbuffer[bc.globalIndex()].foundZDCId;
uint32_t rct = 0;
evsel(bcselbuffer[bc.globalIndex()].alias, bcselbuffer[bc.globalIndex()].selection, rct, kFALSE, kFALSE, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, -1, -1);
evsel(bcselbuffer[bc.globalIndex()].alias, bcselbuffer[bc.globalIndex()].selection, rct, kFALSE, kFALSE, foundBC, foundFT0, foundFV0, foundFDD, foundZDC, -1, -1, -1);
}
return;
}
Expand Down Expand Up @@ -1138,9 +1140,10 @@ class EventSelectionModule

// pre-loop for occupancy calculation
std::vector<bool> vIsFullInfoForOccupancy(cols.size(), 0); // info for occupancy in +/- windows is available (i.e. a given coll is not too close to the TF borders)
std::vector<bool> vIsCollAtROFborder(cols.size(), 0); // collision is close to ITS ROF border
std::vector<bool> vIsCollRejectedByTFborderCut(cols.size(), 0); // helper vector with
std::vector<bool> vCanHaveAssocCollsWithinLastDriftTime(cols.size(), 0); // to see if for some collisions in the occupancy calc (that are close to TF border) we will switch to FT0C based occupancy estimation
std::vector<int> vProxyForCollNtracks(cols.size(), 0); // n tracks, or - if a collision is close to ITS ROF border - scaled FT0C mult
std::vector<float> vMedianTimeForOccupancy(cols.size(), 0); // median time for the occupancy in a defined time window

const float timeWinOccupancyCalcMinNS = evselOpts.confTimeIntervalForOccupancyCalculationMin * 1e3; // ns
const float timeWinOccupancyCalcMaxNS = evselOpts.confTimeIntervalForOccupancyCalculationMax * 1e3; // ns
Expand All @@ -1153,10 +1156,9 @@ class EventSelectionModule

int32_t foundBC = vFoundBCindex[colIndex];
auto bcselEntry = bcselbuffer[foundBC];
// check if we are close to ROF or TF borders => N tracks are not reliable, but FT0 can be used for occupancy estimation
if (!bitcheck64(bcselEntry.selection, aod::evsel::kNoITSROFrameBorder)) {
vIsCollAtROFborder[colIndex] = true;
}
// check if we are close to ITS ROF borders => N ITS tracks is not reliable, and FT0C ampl can be used for occupancy estimation
// denominator for vAmpFT0CperColl is the approximate conversion factor b/n FT0C ampl and number of PV tracks after cuts
vProxyForCollNtracks[colIndex] = bitcheck64(bcselEntry.selection, aod::evsel::kNoITSROFrameBorder) ? vTracksITS567perColl[colIndex] : vAmpFT0CperColl[colIndex] / 10.;

if (!bitcheck64(bcselEntry.selection, aod::evsel::kNoTimeFrameBorder)) {
vIsCollRejectedByTFborderCut[colIndex] = true;
Expand Down Expand Up @@ -1239,6 +1241,7 @@ class EventSelectionModule
// ### for occupancy in time windows
std::vector<int> vAssocToThisCol;
std::vector<float> vCollsTimeDeltaWrtGivenColl;
std::vector<float> vProxyNtracksAssocColls; // mult of associated collisions, for the median time calc
// find all collisions in time window before the current one
minColIndex = colIndex - 1;
while (minColIndex >= 0) {
Expand All @@ -1253,6 +1256,7 @@ class EventSelectionModule
break;
vAssocToThisCol.push_back(minColIndex);
vCollsTimeDeltaWrtGivenColl.push_back(dt);
vProxyNtracksAssocColls.push_back(vProxyForCollNtracks[minColIndex]);
minColIndex--;
}
// find all collisions in time window after the current one
Expand All @@ -1267,10 +1271,31 @@ class EventSelectionModule
break;
vAssocToThisCol.push_back(maxColIndex);
vCollsTimeDeltaWrtGivenColl.push_back(dt);
vProxyNtracksAssocColls.push_back(vProxyForCollNtracks[maxColIndex]);
maxColIndex++;
}
vCollsInTimeWin.push_back(vAssocToThisCol);
vTimeDeltaForColls.push_back(vCollsTimeDeltaWrtGivenColl);

// calculation of the median time for the occupancy in a given time window
std::vector<std::pair<float, float>> pairsDeltaTimeMult;
for (size_t iCol = 0; iCol < vCollsTimeDeltaWrtGivenColl.size(); iCol++)
pairsDeltaTimeMult.emplace_back(vCollsTimeDeltaWrtGivenColl[iCol], vProxyNtracksAssocColls[iCol]);
std::sort(pairsDeltaTimeMult.begin(), pairsDeltaTimeMult.end()); // sorts by first element by default
int proxyTotalMultInTimeWin = std::accumulate(vProxyNtracksAssocColls.begin(), vProxyNtracksAssocColls.end(), 0);

float sumMult = 0.0;
for (size_t iCol = 0; iCol < vCollsTimeDeltaWrtGivenColl.size(); iCol++) {
sumMult += pairsDeltaTimeMult[iCol].second;
if (sumMult > proxyTotalMultInTimeWin / 2.0) {
vMedianTimeForOccupancy[colIndex] = pairsDeltaTimeMult[iCol].first / 1e3; // ns -> us
break;
}
}
for (size_t iCol = 0; iCol < vCollsTimeDeltaWrtGivenColl.size(); iCol++) {
LOGP(debug, "dt={} mult={}", pairsDeltaTimeMult[iCol].first, pairsDeltaTimeMult[iCol].second);
}
LOGP(debug, " --> median time = {}", vMedianTimeForOccupancy[colIndex]);
}

// perform the occupancy calculation per ITS ROF and also in the pre-defined time window
Expand Down Expand Up @@ -1330,15 +1355,11 @@ class EventSelectionModule
for (uint32_t iCol = 0; iCol < vAssocToThisCol.size(); iCol++) {
int thisColIndex = vAssocToThisCol[iCol];
float dt = vCollsTimeDeltaWrtGivenColl[iCol] / 1e3; // ns -> us

// check if we are close to ITS ROF borders => N ITS tracks is not reliable, and FT0C ampl can be used for occupancy estimation
// denominator for vAmpFT0CperColl is the approximate conversion factor b/n FT0C ampl and number of PV tracks after cuts
int nItsTracksAssocColl = !vIsCollAtROFborder[thisColIndex] ? vTracksITS567perColl[thisColIndex] : vAmpFT0CperColl[thisColIndex] / 10.;
// counting tracks from other collisions in fixed time windows
if (std::fabs(dt) < evselOpts.confTimeRangeVetoOnCollNarrow)
nITS567tracksForVetoNarrow += nItsTracksAssocColl;
nITS567tracksForVetoNarrow += vProxyForCollNtracks[thisColIndex];
if (std::fabs(dt) < evselOpts.confTimeRangeVetoOnCollStrict)
nITS567tracksForVetoStrict += nItsTracksAssocColl;
nITS567tracksForVetoStrict += vProxyForCollNtracks[thisColIndex];

// veto on high-mult collisions nearby, where artificial structures in the dt-occupancy plots are observed
if (dt > -4.0 && dt < 2.0 && vAmpFT0CperColl[thisColIndex] > evselOpts.confFT0CamplCutVetoOnCollInTimeRange) { // dt in us // o2-linter: disable=magic-number
Expand All @@ -1359,7 +1380,7 @@ class EventSelectionModule
// weighted occupancy
wOccup = calcWeightForOccupancy(dt);
}
nITS567tracksInFullTimeWindow += wOccup * nItsTracksAssocColl;
nITS567tracksInFullTimeWindow += wOccup * vProxyForCollNtracks[thisColIndex];
sumAmpFT0CInFullTimeWindow += wOccup * vAmpFT0CperColl[thisColIndex];
}
}
Expand Down Expand Up @@ -1410,7 +1431,7 @@ class EventSelectionModule
vNoCollInTimeRangeNarrow[colIndex] = (nITS567tracksForVetoNarrow == 0);
vNoCollInTimeRangeStrict[colIndex] = (nITS567tracksForVetoStrict == 0);
vNoHighMultCollInTimeRange[colIndex] = (nCollsWithFT0CAboveVetoStandard == 0 && nITS567tracksForVetoNarrow == 0);
}
} // end of the occupancy calculation

for (const auto& col : cols) {
int32_t colIndex = col.globalIndex();
Expand Down Expand Up @@ -1490,7 +1511,7 @@ class EventSelectionModule
}

evsel(alias, selection, rct, sel7, sel8, foundBC, foundFT0, foundFV0, foundFDD, foundZDC,
vNumTracksITS567inFullTimeWin[colIndex], vSumAmpFT0CinFullTimeWin[colIndex]);
vNumTracksITS567inFullTimeWin[colIndex], vSumAmpFT0CinFullTimeWin[colIndex], vMedianTimeForOccupancy[colIndex]);
}
} // end processRun3
}; // end EventSelectionModule
Expand Down
Loading