Skip to content

Commit f756efc

Browse files
nepeivodaRSnepeivodaRS
committed
Fix UPC timing and generated neutron selection
1 parent 83ba53c commit f756efc

1 file changed

Lines changed: 116 additions & 141 deletions

File tree

‎PWGLF/Tasks/Strangeness/derivedupcanalysis.cxx‎

Lines changed: 116 additions & 141 deletions
Original file line numberDiff line numberDiff line change
@@ -176,10 +176,15 @@ struct Derivedupcanalysis {
176176
Configurable<float> ft0c{"ft0c", 50., "FT0C threshold"};
177177
Configurable<float> zdc{"zdc", 1., "ZDC threshold"};
178178
Configurable<float> fddaTimeCut{"fddaTimeCut", -1., "FDDA timing cut (ns); negative: no cut"};
179+
Configurable<bool> requireFDDATiming{"requireFDDATiming", true, "require valid FDDA timing for gap-side selection when fddaTimeCut is enabled"};
179180
Configurable<float> fddcTimeCut{"fddcTimeCut", -1., "FDDC timing cut (ns); negative: no cut"};
181+
Configurable<bool> requireFDDCTiming{"requireFDDCTiming", true, "require valid FDDC timing for gap-side selection when fddcTimeCut is enabled"};
180182
Configurable<float> fv0aTimeCut{"fv0aTimeCut", -1., "FV0A timing cut (ns); negative: no cut"};
183+
Configurable<bool> requireFV0ATiming{"requireFV0ATiming", true, "require valid FV0A timing for gap-side selection when fv0aTimeCut is enabled"};
181184
Configurable<float> ft0aTimeCut{"ft0aTimeCut", -1., "FT0A timing cut (ns); negative: no cut"};
185+
Configurable<bool> requireFT0ATiming{"requireFT0ATiming", true, "require valid FT0A timing for gap-side selection when ft0aTimeCut is enabled"};
182186
Configurable<float> ft0cTimeCut{"ft0cTimeCut", -1., "FT0C timing cut (ns); negative: no cut"};
187+
Configurable<bool> requireFT0CTiming{"requireFT0CTiming", true, "require valid FT0C timing for gap-side selection when ft0cTimeCut is enabled"};
183188
Configurable<float> zdcTimeCut{"zdcTimeCut", 2., "ZDC timing cut (ns)"};
184189
Configurable<bool> requireZDCTiming{"requireZDCTiming", true, "require valid ZDC timing for gap-side selection"};
185190
Configurable<int> genGapSide{"genGapSide", 0, "0 -- A, 1 -- C, 2 -- double"};
@@ -237,6 +242,7 @@ struct Derivedupcanalysis {
237242
Configurable<bool> calculateFeeddownMatrix{"calculateFeeddownMatrix", true, "fill feeddown matrix if MC"};
238243
ConfigurableAxis axisGeneratorIds{"axisGeneratorIds", {256, -0.5f, 255.5f}, "axis for generatorIds"};
239244
Configurable<bool> checkNeutronsInMC{"checkNeutronsInMC", true, "require no neutrons for single-gap in MC"};
245+
Configurable<bool> requireGeneratedZDCNeutrons{"requireGeneratedZDCNeutrons", false, "require generated ZDC neutron topology matching genGapSide in processGenerated"};
240246
Configurable<float> neutronEtaCut{"neutronEtaCut", 8.8, "ZN acceptance"};
241247

242248
// Occupancy cut
@@ -313,8 +319,6 @@ struct Derivedupcanalysis {
313319
ConfigurableAxis axisCtau{"axisCtau", {200, 0.0f, 20.0f}, "c x tau (cm)"};
314320

315321
static constexpr std::string_view kParticlenames[] = {"K0Short", "Lambda", "AntiLambda", "Xi", "AntiXi", "Omega", "AntiOmega"};
316-
static constexpr uint8_t kFT0TriggerBitIsActiveA = 5;
317-
static constexpr uint8_t kFT0TriggerBitIsActiveC = 6;
318322

319323
void setBits(std::bitset<kSelNum>& mask, std::initializer_list<int> selections)
320324
{
@@ -560,7 +564,7 @@ struct Derivedupcanalysis {
560564
histos.fill(HIST(kParticlenames[partID]) + HIST("/h4dNegDetectPropVsCentrality"), ft0ampl, negTrackExtra.detectorMap(), negTrackExtra.itsClusterMap(), pT);
561565
}
562566
if (doDetectPropQA == 2) {
563-
histos.fill(HIST(kParticlenames[partID]) + HIST("/h7dPosDetectPropVsCentrality"), ft0ampl, posDetMap, posITSclusMap, negDetMap, negITSclusMap, pT, invMass);
567+
histos.fill(HIST(kParticlenames[partID]) + HIST("/h7dDetectPropVsCentrality"), ft0ampl, posDetMap, posITSclusMap, negDetMap, negITSclusMap, pT, invMass);
564568
histos.fill(HIST(kParticlenames[partID]) + HIST("/h5dPosDetectPropVsCentrality"), ft0ampl, posTrackExtra.detectorMap(), posTrackExtra.itsClusterMap(), pT, invMass);
565569
histos.fill(HIST(kParticlenames[partID]) + HIST("/h5dNegDetectPropVsCentrality"), ft0ampl, negTrackExtra.detectorMap(), negTrackExtra.itsClusterMap(), pT, invMass);
566570
}
@@ -844,9 +848,9 @@ struct Derivedupcanalysis {
844848
maskAntiLambdaSpecific.set(selTPCPIDNegativeProton);
845849

846850
maskXiSpecific.set(selTPCPIDNegativePion);
847-
maskAntiXiSpecific.set(selTPCPIDPositiveProton);
851+
maskAntiXiSpecific.set(selTPCPIDNegativeProton);
848852
maskOmegaSpecific.set(selTPCPIDNegativePion);
849-
maskAntiOmegaSpecific.set(selTPCPIDPositiveProton);
853+
maskAntiOmegaSpecific.set(selTPCPIDNegativeProton);
850854
}
851855
// TOF PID
852856
if (PIDConfigurations.tofPidNsigmaCutK0Pi < 1e+5) { // safeguard for no cut
@@ -924,8 +928,8 @@ struct Derivedupcanalysis {
924928
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(10, "kNoSameBunchPileup");
925929
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(11, "kNoCollInTimeRangeStd");
926930
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(12, "kNoCollInTimeRangeNarrow");
927-
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(13, "Below min occup.");
928-
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(14, "Above max occup.");
931+
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(13, "Above min occup.");
932+
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(14, "Below max occup.");
929933
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(15, "RCTFlagsChecker");
930934
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(16, "isUPC");
931935
histos.get<TH1>(HIST("eventQA/hEventSelection"))->GetXaxis()->SetBinLabel(17, "has UPC flag");
@@ -1097,14 +1101,10 @@ struct Derivedupcanalysis {
10971101
const float timeZNC = collision.timeZNC();
10981102
const float cut = upcCuts.zdcTimeCut;
10991103

1100-
auto isInvalidTime = [](float time) {
1101-
return !std::isfinite(time) || (std::abs(time) == 999.f);
1102-
};
1103-
1104-
const bool gapA = isInvalidTime(timeZNA) || (std::abs(timeZNA) > cut);
1105-
const bool gapC = isInvalidTime(timeZNC) || (std::abs(timeZNC) > cut);
1106-
const bool neutronA = !isInvalidTime(timeZNA) && (std::abs(timeZNA) < cut);
1107-
const bool neutronC = !isInvalidTime(timeZNC) && (std::abs(timeZNC) < cut);
1104+
const bool gapA = hasNoZdcNeutrons(timeZNA, cut);
1105+
const bool gapC = hasNoZdcNeutrons(timeZNC, cut);
1106+
const bool neutronA = hasZdcNeutrons(timeZNA, cut);
1107+
const bool neutronC = hasZdcNeutrons(timeZNC, cut);
11081108

11091109
if (selGapSide == o2::aod::sgselector::SingleGapA) { // 0nXn
11101110
if (!(gapA && neutronC)) {
@@ -1124,29 +1124,29 @@ struct Derivedupcanalysis {
11241124
return selGapSide;
11251125
}
11261126

1127-
bool isInvalidTime(float time) const
1127+
bool hasNoZdcNeutrons(float time, float cut) const
11281128
{
1129-
return !std::isfinite(time) || (std::abs(time) >= 998.f);
1129+
return (time != -999.f) && (std::abs(time) > cut);
11301130
}
11311131

1132-
bool isTimingCutEnabled(float cut) const
1132+
bool hasZdcNeutrons(float time, float cut) const
11331133
{
1134-
return cut >= 0.f;
1134+
return std::abs(time) < cut;
11351135
}
11361136

1137-
bool isTimingGap(float time, float cut) const
1137+
bool isTimingCutEnabled(float cut, bool requireTiming) const
11381138
{
1139-
return isInvalidTime(time) || (std::abs(time) > cut);
1139+
return requireTiming && cut >= 0.f;
11401140
}
11411141

1142-
bool isTimingActivity(float time, float cut) const
1142+
bool isTimingGap(float time, float cut) const
11431143
{
1144-
return !isInvalidTime(time) && (std::abs(time) < cut);
1144+
return (time != -999.f) && (std::abs(time) > cut);
11451145
}
11461146

1147-
bool hasFT0Activity(uint8_t triggerMask, uint8_t bit) const
1147+
bool isTimingActivity(float time, float cut) const
11481148
{
1149-
return (triggerMask & (static_cast<uint8_t>(1u) << bit)) != 0;
1149+
return std::abs(time) < cut;
11501150
}
11511151

11521152
template <typename TCollision>
@@ -1158,29 +1158,26 @@ struct Derivedupcanalysis {
11581158
return selGapSide;
11591159
}
11601160

1161-
const bool useFDDA = isTimingCutEnabled(upcCuts.fddaTimeCut);
1162-
const bool useFDDC = isTimingCutEnabled(upcCuts.fddcTimeCut);
1163-
const bool useFV0A = isTimingCutEnabled(upcCuts.fv0aTimeCut);
1164-
const bool useFT0A = isTimingCutEnabled(upcCuts.ft0aTimeCut);
1165-
const bool useFT0C = isTimingCutEnabled(upcCuts.ft0cTimeCut);
1161+
const bool useFDDA = isTimingCutEnabled(upcCuts.fddaTimeCut, upcCuts.requireFDDATiming);
1162+
const bool useFDDC = isTimingCutEnabled(upcCuts.fddcTimeCut, upcCuts.requireFDDCTiming);
1163+
const bool useFV0A = isTimingCutEnabled(upcCuts.fv0aTimeCut, upcCuts.requireFV0ATiming);
1164+
const bool useFT0A = isTimingCutEnabled(upcCuts.ft0aTimeCut, upcCuts.requireFT0ATiming);
1165+
const bool useFT0C = isTimingCutEnabled(upcCuts.ft0cTimeCut, upcCuts.requireFT0CTiming);
11661166

11671167
if (!(useFDDA || useFDDC || useFV0A || useFT0A || useFT0C)) {
11681168
return selGapSide;
11691169
}
11701170

1171-
const bool ft0ActiveA = hasFT0Activity(collision.triggerMaskFT0(), kFT0TriggerBitIsActiveA);
1172-
const bool ft0ActiveC = hasFT0Activity(collision.triggerMaskFT0(), kFT0TriggerBitIsActiveC);
1173-
11741171
const bool gapFDDA = !useFDDA || isTimingGap(collision.timeFDDA(), upcCuts.fddaTimeCut);
11751172
const bool actFDDA = !useFDDA || isTimingActivity(collision.timeFDDA(), upcCuts.fddaTimeCut);
11761173
const bool gapFDDC = !useFDDC || isTimingGap(collision.timeFDDC(), upcCuts.fddcTimeCut);
11771174
const bool actFDDC = !useFDDC || isTimingActivity(collision.timeFDDC(), upcCuts.fddcTimeCut);
11781175
const bool gapFV0A = !useFV0A || isTimingGap(collision.timeFV0A(), upcCuts.fv0aTimeCut);
11791176
const bool actFV0A = !useFV0A || isTimingActivity(collision.timeFV0A(), upcCuts.fv0aTimeCut);
1180-
const bool gapFT0A = !useFT0A || !ft0ActiveA || isTimingGap(collision.timeFT0A(), upcCuts.ft0aTimeCut);
1181-
const bool actFT0A = !useFT0A || (ft0ActiveA && isTimingActivity(collision.timeFT0A(), upcCuts.ft0aTimeCut));
1182-
const bool gapFT0C = !useFT0C || !ft0ActiveC || isTimingGap(collision.timeFT0C(), upcCuts.ft0cTimeCut);
1183-
const bool actFT0C = !useFT0C || (ft0ActiveC && isTimingActivity(collision.timeFT0C(), upcCuts.ft0cTimeCut));
1177+
const bool gapFT0A = !useFT0A || isTimingGap(collision.timeFT0A(), upcCuts.ft0aTimeCut);
1178+
const bool actFT0A = !useFT0A || isTimingActivity(collision.timeFT0A(), upcCuts.ft0aTimeCut);
1179+
const bool gapFT0C = !useFT0C || isTimingGap(collision.timeFT0C(), upcCuts.ft0cTimeCut);
1180+
const bool actFT0C = !useFT0C || isTimingActivity(collision.timeFT0C(), upcCuts.ft0cTimeCut);
11841181

11851182
if (selGapSide == o2::aod::sgselector::SingleGapA) {
11861183
if (!(gapFV0A && gapFDDA && gapFT0A && actFDDC && actFT0C)) {
@@ -1199,6 +1196,68 @@ struct Derivedupcanalysis {
11991196
return selGapSide;
12001197
}
12011198

1199+
template <typename TCollision>
1200+
int getGapSideWithoutRecoZDC(TCollision const& collision)
1201+
{
1202+
int selGapSide = o2::aod::sgselector::NoGap;
1203+
selGapSide = sgSelector.trueGap(collision, upcCuts.fv0a, upcCuts.ft0a, upcCuts.ft0c, std::numeric_limits<float>::infinity());
1204+
return applyFITTiming(selGapSide, collision);
1205+
}
1206+
1207+
template <typename TNeutrons>
1208+
bool generatedNeutronsMatchGapSide(int selGapSide, TNeutrons const& neutrons)
1209+
{
1210+
bool neutronA = false;
1211+
bool neutronC = false;
1212+
for (const auto& neutron : neutrons) {
1213+
neutronA = neutronA || (neutron.eta() > neutronEtaCut);
1214+
neutronC = neutronC || (neutron.eta() < -neutronEtaCut);
1215+
}
1216+
1217+
if (selGapSide == o2::aod::sgselector::SingleGapA) { // 0nXn
1218+
return !neutronA && neutronC;
1219+
}
1220+
if (selGapSide == o2::aod::sgselector::SingleGapC) { // Xn0n
1221+
return neutronA && !neutronC;
1222+
}
1223+
if (selGapSide == o2::aod::sgselector::DoubleGap) {
1224+
return !neutronA && !neutronC;
1225+
}
1226+
1227+
return false;
1228+
}
1229+
1230+
template <typename TNeutrons>
1231+
int applyGeneratedNeutronSelection(int selGapSide, TNeutrons const& neutrons)
1232+
{
1233+
if (!checkNeutronsInMC) {
1234+
return selGapSide;
1235+
}
1236+
if (selGapSide != o2::aod::sgselector::SingleGapA &&
1237+
selGapSide != o2::aod::sgselector::SingleGapC &&
1238+
selGapSide != o2::aod::sgselector::DoubleGap) {
1239+
return selGapSide;
1240+
}
1241+
1242+
if (!generatedNeutronsMatchGapSide(selGapSide, neutrons)) {
1243+
selGapSide = o2::aod::sgselector::NoGap;
1244+
}
1245+
1246+
return selGapSide;
1247+
}
1248+
1249+
template <typename TCollision, typename TNeutrons>
1250+
int getGapSideMC(TCollision const& collision, TNeutrons const& neutrons)
1251+
{
1252+
return applyGeneratedNeutronSelection(getGapSideWithoutRecoZDC(collision), neutrons);
1253+
}
1254+
1255+
template <typename TNeutrons>
1256+
bool acceptGeneratedNeutronSelection(TNeutrons const& neutrons)
1257+
{
1258+
return !requireGeneratedZDCNeutrons || generatedNeutronsMatchGapSide(static_cast<int>(upcCuts.genGapSide), neutrons);
1259+
}
1260+
12021261
template <typename TCollision>
12031262
int getGapSide(TCollision const& collision)
12041263
{
@@ -1962,31 +2021,7 @@ struct Derivedupcanalysis {
19622021
continue;
19632022
}
19642023

1965-
int selGapSide = collision.isUPC() ? getGapSide(collision) : -1;
1966-
if (checkNeutronsInMC) {
1967-
for (const auto& neutron : groupedNeutrons) {
1968-
if (selGapSide < -0.5)
1969-
break;
1970-
1971-
const float eta = neutron.eta();
1972-
switch (selGapSide) {
1973-
case 0: // SGA
1974-
if (eta > neutronEtaCut)
1975-
selGapSide = -1;
1976-
break;
1977-
case 1: // SGC
1978-
if (eta < -neutronEtaCut)
1979-
selGapSide = -1;
1980-
break;
1981-
case 2: // DG
1982-
if (eta > neutronEtaCut)
1983-
selGapSide = 1;
1984-
else if (eta < -neutronEtaCut)
1985-
selGapSide = 0;
1986-
break;
1987-
}
1988-
}
1989-
}
2024+
int selGapSide = collision.isUPC() ? getGapSideMC(collision, groupedNeutrons) : -1;
19902025

19912026
if (evSels.studyUPConly && (selGapSide != static_cast<int>(upcCuts.genGapSide)))
19922027
continue;
@@ -2040,31 +2075,7 @@ struct Derivedupcanalysis {
20402075
continue;
20412076
}
20422077

2043-
int selGapSide = collision.isUPC() ? getGapSide(collision) : -1;
2044-
if (checkNeutronsInMC) {
2045-
for (const auto& neutron : groupedNeutrons) {
2046-
if (selGapSide < -0.5)
2047-
break;
2048-
2049-
const float eta = neutron.eta();
2050-
switch (selGapSide) {
2051-
case 0: // SGA
2052-
if (eta > neutronEtaCut)
2053-
selGapSide = -1;
2054-
break;
2055-
case 1: // SGC
2056-
if (eta < -neutronEtaCut)
2057-
selGapSide = -1;
2058-
break;
2059-
case 2: // DG
2060-
if (eta > neutronEtaCut)
2061-
selGapSide = 1;
2062-
else if (eta < -neutronEtaCut)
2063-
selGapSide = 0;
2064-
break;
2065-
}
2066-
}
2067-
}
2078+
int selGapSide = collision.isUPC() ? getGapSideMC(collision, groupedNeutrons) : -1;
20682079

20692080
const bool passStd = !evSels.studyUPConly || (selGapSide == static_cast<int>(upcCuts.genGapSide));
20702081
if (passStd) {
@@ -2208,34 +2219,9 @@ struct Derivedupcanalysis {
22082219
}
22092220
histos.fill(HIST("eventQA/hRawGapSide"), collision.gapSide());
22102221

2211-
int selGapSide = collision.isUPC() ? getGapSide(collision) : -1;
2212-
int selGapSideNoNeutrons = selGapSide;
2213-
22142222
auto groupedNeutrons = neutrons.sliceBy(neutronsPerMcCollision, mcCollision.globalIndex());
2215-
if (checkNeutronsInMC) {
2216-
for (const auto& neutron : groupedNeutrons) {
2217-
if (selGapSide < -0.5)
2218-
break;
2219-
2220-
const float eta = neutron.eta();
2221-
switch (selGapSide) {
2222-
case 0: // SGA
2223-
if (eta > neutronEtaCut)
2224-
selGapSide = -1;
2225-
break;
2226-
case 1: // SGC
2227-
if (eta < -neutronEtaCut)
2228-
selGapSide = -1;
2229-
break;
2230-
case 2: // DG
2231-
if (eta > neutronEtaCut)
2232-
selGapSide = 1;
2233-
else if (eta < -neutronEtaCut)
2234-
selGapSide = 0;
2235-
break;
2236-
}
2237-
}
2238-
}
2223+
int selGapSideNoNeutrons = collision.isUPC() ? getGapSideWithoutRecoZDC(collision) : -1;
2224+
int selGapSide = applyGeneratedNeutronSelection(selGapSideNoNeutrons, groupedNeutrons);
22392225

22402226
if (evSels.studyUPConly && (selGapSide < -0.5))
22412227
continue;
@@ -2351,34 +2337,9 @@ struct Derivedupcanalysis {
23512337
}
23522338
histos.fill(HIST("eventQA/hRawGapSide"), collision.gapSide());
23532339

2354-
int selGapSide = collision.isUPC() ? getGapSide(collision) : -1;
2355-
int selGapSideNoNeutrons = selGapSide;
2356-
23572340
auto groupedNeutrons = neutrons.sliceBy(neutronsPerMcCollision, mcCollision.globalIndex());
2358-
if (checkNeutronsInMC) {
2359-
for (const auto& neutron : groupedNeutrons) {
2360-
if (selGapSide < -0.5)
2361-
break;
2362-
2363-
const float eta = neutron.eta();
2364-
switch (selGapSide) {
2365-
case 0: // SGA
2366-
if (eta > neutronEtaCut)
2367-
selGapSide = -1;
2368-
break;
2369-
case 1: // SGC
2370-
if (eta < -neutronEtaCut)
2371-
selGapSide = -1;
2372-
break;
2373-
case 2: // DG
2374-
if (eta > neutronEtaCut)
2375-
selGapSide = 1;
2376-
else if (eta < -neutronEtaCut)
2377-
selGapSide = 0;
2378-
break;
2379-
}
2380-
}
2381-
}
2341+
int selGapSideNoNeutrons = collision.isUPC() ? getGapSideWithoutRecoZDC(collision) : -1;
2342+
int selGapSide = applyGeneratedNeutronSelection(selGapSideNoNeutrons, groupedNeutrons);
23822343

23832344
if (evSels.studyUPConly && (selGapSide < -0.5))
23842345
continue;
@@ -2447,6 +2408,11 @@ struct Derivedupcanalysis {
24472408
continue;
24482409
}
24492410

2411+
auto groupedNeutrons = neutrons.sliceBy(neutronsPerMcCollision, mcCollision.globalIndex());
2412+
if (!acceptGeneratedNeutronSelection(groupedNeutrons)) {
2413+
continue;
2414+
}
2415+
24502416
// float centrality = -1.f;
24512417
float ft0ampl = -1.f;
24522418
int nTracksGlobal = -1;
@@ -2503,6 +2469,15 @@ struct Derivedupcanalysis {
25032469
if (std::abs(mcCollision.posZ()) > maxZVtxPosition)
25042470
continue;
25052471

2472+
if (std::find(generatorIds->begin(), generatorIds->end(), mcCollision.generatorsID()) == generatorIds->end()) {
2473+
continue;
2474+
}
2475+
2476+
auto groupedNeutrons = neutrons.sliceBy(neutronsPerMcCollision, mcCollision.globalIndex());
2477+
if (!acceptGeneratedNeutronSelection(groupedNeutrons)) {
2478+
continue;
2479+
}
2480+
25062481
// float centrality = -1.f;
25072482
float ft0ampl = -1.f;
25082483
int nTracksGlobal = -1;

0 commit comments

Comments
 (0)