Skip to content

Commit 72664a5

Browse files
SteffimroSteffimro
authored
[PWGEM] add adjustable QA settings + qinv gate for mixed pairs + fix headers (#17228)
1 parent a3556ab commit 72664a5

1 file changed

Lines changed: 57 additions & 14 deletions

File tree

‎PWGEM/PhotonMeson/Tasks/photonhbt.cxx‎

Lines changed: 57 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -29,6 +29,7 @@
2929
#include <CCDB/BasicCCDBManager.h>
3030
#include <CommonConstants/MathConstants.h>
3131
#include <DataFormatsParameters/GRPMagField.h>
32+
#include <Framework/ASoA.h>
3233
#include <Framework/ASoAHelpers.h>
3334
#include <Framework/AnalysisDataModel.h>
3435
#include <Framework/AnalysisHelpers.h>
@@ -47,6 +48,7 @@
4748
#include <Math/Vector3Dfwd.h>
4849
#include <Math/Vector4D.h> // IWYU pragma: keep
4950
#include <Math/Vector4Dfwd.h>
51+
#include <TH1.h>
5052
#include <TPDGCode.h>
5153
#include <TString.h>
5254

@@ -222,6 +224,9 @@ struct Photonhbt {
222224
Configurable<float> cfgMaxQinvForQA{"cfgMaxQinvForQA", 0.1f, "fill per-step pair QA histograms only when q_inv < this value"};
223225
Configurable<float> cfgMaxQinvForFullRange{"cfgMaxQinvForFullRange", 0.3f, "fill full-range histograms only when q_inv < this value"};
224226
Configurable<float> cfgMaxQinvForMCQA{"cfgMaxQinvForMCQA", 0.3f, "fill MC truth 1D histograms only when q_inv < this value"};
227+
Configurable<int> cfgQaLevel{"cfgQaLevel", 2, "QA: 0 no QA, 1 standard, 2 full diagnostics with leg information etc."};
228+
Configurable<bool> cfgFillDRDZSparse{"cfgFillDRDZSparse", true, "book/fill the FullRange |R1-R2|-Deltaz-qinv sparse (large; needs cfgQaLevel>=2)"};
229+
Configurable<float> cfgMaxQinvForProcessing{"cfgMaxQinvForProcessing", 0.5, "skip mixed pairs with q_inv above this before building observables"};
225230
} qaflags;
226231

227232
// ─── HBT analysis mode ───────────────────────────────────────────────────────────
@@ -507,9 +512,20 @@ struct Photonhbt {
507512
// INITS
508513
/*************************************************/
509514

515+
bool mDoPairQa{true}, mDoSinglePhotonQa{true}, mDoLegPairQA{true};
516+
bool mDoPairSepQA{true}, mFillDRDZSparse{true};
517+
510518
void init(InitContext& context)
511519
{
512520
isMC = context.mOptions.get<bool>("processMC");
521+
const int qaLevel = qaflags.cfgQaLevel.value;
522+
mDoPairQa = qaflags.doPairQa.value && qaLevel >= 1; // o2-linter: disable=magic-number (QA set-up)
523+
mDoSinglePhotonQa = qaflags.doSinglePhotonQa.value && qaLevel >= 1; // o2-linter: disable=magic-number (QA set-up)
524+
mDoLegPairQA = qaflags.doLegPairQA.value && qaLevel >= 2; // o2-linter: disable=magic-number (QA set-up)
525+
mDoPairSepQA = pairsep.cfgDoPairSepQA.value && qaLevel >= 2; // o2-linter: disable=magic-number (QA set-up)
526+
mFillDRDZSparse = qaflags.cfgFillDRDZSparse.value && qaLevel >= 2; // o2-linter: disable=magic-number (QA set-up)
527+
LOGF(info, "photonhbt QA level %d -> pairQA %d, singlePhotonQA %d, legPairQA %d, pairSep %d, dRdZ sparse %d",
528+
qaLevel, mDoPairQa, mDoSinglePhotonQa, mDoLegPairQA, mDoPairSepQA, mFillDRDZSparse);
513529
mRunNumber = 0;
514530
parseBins(mixing.confVtxBins, ztxBinEdges);
515531
parseBins(mixing.confCentBins, centBinEdges);
@@ -676,7 +692,7 @@ struct Photonhbt {
676692

677693
if (isMC) {
678694
addPairMCHistograms();
679-
if (qaflags.doLegPairQA) {
695+
if (mDoLegPairQA) {
680696
addLegPairMCHistograms();
681697
}
682698
addTruthMCHistograms();
@@ -763,8 +779,10 @@ struct Photonhbt {
763779
}
764780

765781
fRegistryPairQA.addClone("Pair/same/QA/", "Pair/mix/QA/");
766-
addLegPairQAForStep("Pair/same/QA/Before/");
767-
addLegPairQAForStep("Pair/same/QA/AfterPairCuts/");
782+
if (mDoLegPairQA) {
783+
addLegPairQAForStep("Pair/same/QA/Before/");
784+
addLegPairQAForStep("Pair/same/QA/AfterPairCuts/");
785+
}
768786
}
769787

770788
void addPairMCHistograms()
@@ -1005,7 +1023,9 @@ struct Photonhbt {
10051023
fRegistryCF.add((path + "hDeltaR3DVsQinv").c_str(), "#Delta r_{3D} vs q_{inv};q_{inv} (GeV/c);#Delta r_{3D} (cm)", kTH2D, {axisQinv, axisDeltaR3D}, true);
10061024
fRegistryCF.add((path + "hQinvVsCent").c_str(), "q_{inv} vs centrality;centrality (%);q_{inv} (GeV/c)", kTH2D, {axisCentQA, axisQinv}, true);
10071025
fRegistryCF.add((path + "hQinvVsOccupancy").c_str(), "q_{inv} vs occupancy;occupancy;q_{inv} (GeV/c)", kTH2D, {axisOccupancy, axisQinv}, true);
1008-
fRegistryCF.add((path + "hSparseDeltaRDeltaZQinv").c_str(), "|R_{1}-R_{2}|,#Delta z,q_{inv}", kTHnSparseD, {axisDeltaR, axisDeltaZ, axisQinv}, true);
1026+
if (mFillDRDZSparse) {
1027+
fRegistryCF.add((path + "hSparseDeltaRDeltaZQinv").c_str(), "|R_{1}-R_{2}|,#Delta z,q_{inv}", kTHnSparseD, {axisDeltaR, axisDeltaZ, axisQinv}, true);
1028+
}
10091029
fRegistryCF.add((path + "hDeltaRCosOAVsQinv").c_str(), "#Delta r/cos(#theta_{op}/2) vs q_{inv};q_{inv} (GeV/c);#Delta r/cos(#theta_{op}/2) (cm)", kTH2D, {axisQinv, {100, 0, 100}}, true);
10101030
}
10111031

@@ -1266,6 +1286,19 @@ struct Photonhbt {
12661286
return s;
12671287
}
12681288

1289+
template <typename TG1, typename TG2>
1290+
[[nodiscard]] inline bool passFastQinvGate(TG1 const& g1, TG2 const& g2) const
1291+
{
1292+
const float qMax = qaflags.cfgMaxQinvForProcessing.value;
1293+
if (qMax > 1e9f) { // o2-linter: disable=magic-number (skip if non-sensical value is chosen)
1294+
return true;
1295+
}
1296+
const float dEta = g1.eta() - g2.eta();
1297+
const float dPhi = RecoDecay::constrainAngle(g1.phi() - g2.phi(), -o2::constants::math::PI);
1298+
const float q2 = 2.f * g1.pt() * g2.pt() * (std::cosh(dEta) - std::cos(dPhi));
1299+
return q2 <= qMax * qMax;
1300+
}
1301+
12691302
[[nodiscard]] inline bool passLegSepCut(PairSep const& s) const
12701303
{
12711304
if (!pairsep.cfgDoLegSepCut.value) {
@@ -1386,13 +1419,15 @@ struct Photonhbt {
13861419
fRegistryCF.fill(HIST(base) + HIST("hDeltaR3DVsQinv"), obs.qinv, obs.deltaR3D);
13871420
fRegistryCF.fill(HIST(base) + HIST("hQinvVsCent"), cent, obs.qinv);
13881421
fRegistryCF.fill(HIST(base) + HIST("hQinvVsOccupancy"), occupancy, obs.qinv);
1389-
fRegistryCF.fill(HIST(base) + HIST("hSparseDeltaRDeltaZQinv"), obs.deltaR, obs.deltaZ, obs.qinv);
1422+
if (mFillDRDZSparse) {
1423+
fRegistryCF.fill(HIST(base) + HIST("hSparseDeltaRDeltaZQinv"), obs.deltaR, obs.deltaZ, obs.qinv);
1424+
}
13901425
}
13911426

13921427
template <int ev_id, bool after_cut = false>
13931428
inline void fillPairSep(PairSep const& s, PairQAObservables const& obs)
13941429
{
1395-
if (!pairsep.cfgDoPairSepQA.value) {
1430+
if (!mDoPairSepQA) {
13961431
return;
13971432
}
13981433
const float limit = pairsep.cfgPairSepMaxQinv.value;
@@ -1422,7 +1457,7 @@ struct Photonhbt {
14221457
template <int step_id, typename TPhoton>
14231458
inline void fillSinglePhotonQAStep(TPhoton const& g)
14241459
{
1425-
if (!qaflags.doSinglePhotonQa) {
1460+
if (!mDoSinglePhotonQa) {
14261461
return;
14271462
}
14281463
constexpr auto base = singlePhotonQAPrefix<step_id>();
@@ -1519,7 +1554,7 @@ struct Photonhbt {
15191554
template <int ev_id, int step_id>
15201555
inline void fillPairQAStep(PairQAObservables const& o, float /*cent*/, float /*occupancy*/)
15211556
{
1522-
if (!qaflags.doPairQa) {
1557+
if (!mDoPairQa) {
15231558
return;
15241559
}
15251560
constexpr auto base = qaPrefix<ev_id, step_id>();
@@ -1559,7 +1594,7 @@ struct Photonhbt {
15591594
template <int step_id>
15601595
inline void fillLegPairQAStep(LegPairObservables const& lo, float kt)
15611596
{
1562-
if (!qaflags.doPairQa) {
1597+
if (!mDoLegPairQA) {
15631598
return;
15641599
}
15651600
constexpr auto base = qaPrefix<0, step_id>();
@@ -1948,7 +1983,7 @@ struct Photonhbt {
19481983
auto keyDFCollision = std::make_pair(ndf, collision.globalIndex());
19491984
auto photons1Coll = photons1.sliceBy(perCollision1, collision.globalIndex());
19501985
auto photons2Coll = photons2.sliceBy(perCollision2, collision.globalIndex());
1951-
if (qaflags.doSinglePhotonQa) {
1986+
if (mDoSinglePhotonQa) {
19521987
for (const auto& g : photons1Coll) {
19531988
if (cut1.template IsSelected<decltype(g), TSubInfos1>(g)) {
19541989
fillSinglePhotonQAStep<0>(g);
@@ -2039,7 +2074,7 @@ struct Photonhbt {
20392074
addToPool(g1, pwl1);
20402075
addToPool(g2, pwl2);
20412076
}
2042-
if (qaflags.doSinglePhotonQa) {
2077+
if (mDoSinglePhotonQa) {
20432078
for (const auto& g : photons1Coll) {
20442079
if (cut1.template IsSelected<decltype(g), TSubInfos1>(g)) {
20452080
if (idsAfterPairCuts.contains(g.globalIndex())) {
@@ -2070,6 +2105,9 @@ struct Photonhbt {
20702105
if (!passAsymmetryCut(g1.pt(), g2.pt())) {
20712106
continue;
20722107
}
2108+
if (!passFastQinvGate(g1, g2)) {
2109+
continue;
2110+
}
20732111
auto obs = buildPairQAObservables(g1, g2);
20742112
if (!obs.valid) {
20752113
continue;
@@ -2170,7 +2208,7 @@ struct Photonhbt {
21702208
auto keyBin = std::make_tuple(zbin, centbin, epbin, occbin);
21712209
auto keyDFCollision = std::make_pair(ndf, collision.globalIndex());
21722210
auto photonsColl = photons.sliceBy(perCollision, collision.globalIndex());
2173-
if (qaflags.doSinglePhotonQa) {
2211+
if (mDoSinglePhotonQa) {
21742212
for (const auto& g : photonsColl) {
21752213
if (cut.template IsSelected<decltype(g), TLegs>(g)) {
21762214
fillSinglePhotonQAStep<0>(g);
@@ -2265,7 +2303,9 @@ struct Photonhbt {
22652303
} else {
22662304
const bool doMCQA = passQinvMCQAGate(obs.qinv);
22672305
fillMCPairQA<false>(truthType, obs, doQA, doMCQA);
2268-
fillLegPairMC(truthType, legObs, obs.kt);
2306+
if (mDoLegPairQA) {
2307+
fillLegPairMC(truthType, legObs, obs.kt);
2308+
}
22692309
if (doFR) {
22702310
fillMCPairQAFullRange<false>(truthType, obs);
22712311
}
@@ -2318,7 +2358,7 @@ struct Photonhbt {
23182358
addToPool(g1, pwl1);
23192359
addToPool(g2, pwl2);
23202360
}
2321-
if (qaflags.doSinglePhotonQa) {
2361+
if (mDoSinglePhotonQa) {
23222362
for (const auto& g : photonsColl) {
23232363
if (cut.template IsSelected<decltype(g), TLegs>(g)) {
23242364
if (idsAfterPairCuts.contains(g.globalIndex())) {
@@ -2349,6 +2389,9 @@ struct Photonhbt {
23492389
if (!passAsymmetryCut(g1.pt(), g2.pt())) {
23502390
continue;
23512391
}
2392+
if (!passFastQinvGate(g1, g2)) {
2393+
continue;
2394+
}
23522395
auto obs = buildPairQAObservables(g1, g2);
23532396
if (!obs.valid) {
23542397
continue;

0 commit comments

Comments
 (0)