Skip to content

Commit 84bf6aa

Browse files
scannitoscannito
committed
Changes to accomodate online mult-integrated efficiency correction
1 parent e93ccf1 commit 84bf6aa

1 file changed

Lines changed: 169 additions & 39 deletions

File tree

‎PWGLF/Tasks/Strangeness/phiStrangeCorrelation.cxx‎

Lines changed: 169 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -84,14 +84,19 @@ enum AnalysisMode {
8484
kDeltaYvsDeltaPhi
8585
};
8686

87-
enum ParticleOfInterest {
87+
enum AssociatedParticleType {
88+
kK0S = 0,
89+
kXi,
90+
kPion,
91+
kAssocPartSize
92+
};
93+
94+
/*enum ParticleOfInterest {
8895
Phi = 0,
8996
K0S,
9097
Pion,
91-
/*PionTPC,
92-
PionTPCTOF*/
9398
ParticleOfInterestSize
94-
};
99+
};*/
95100

96101
/*
97102
#define LIST_OF_PARTICLES_OF_INTEREST \
@@ -121,7 +126,9 @@ static constexpr auto particleOfInterestLabels = std::to_array<std::string_view>
121126
});
122127
*/
123128

124-
struct BoundEfficiencyMap {
129+
using EffMapPtr = std::variant<std::shared_ptr<TH2>, std::shared_ptr<TH3>>;
130+
131+
/*struct BoundEfficiencyMap {
125132
using CoordsTuple = std::tuple<float, float, float>;
126133
127134
const TH3* effMap;
@@ -149,6 +156,58 @@ struct BoundEfficiencyMap {
149156
const auto& [x, y, z] = coords;
150157
return effMap->Interpolate(x, y, z);
151158
}
159+
};*/
160+
161+
struct BoundEfficiencyMap {
162+
using CoordsTuple = std::tuple<float, float, float>;
163+
164+
const EffMapPtr& effMap;
165+
CoordsTuple coords;
166+
167+
BoundEfficiencyMap(const EffMapPtr& effMap, float x, float y, float z) : effMap(effMap), coords(x, y, z) {}
168+
BoundEfficiencyMap(const EffMapPtr& effMap, const CoordsTuple& coords) : effMap(effMap), coords(coords) {}
169+
170+
float getBinEfficiency() const
171+
{
172+
return std::visit(
173+
[this](auto&& mapPtr) -> float {
174+
if (!mapPtr)
175+
return 1.0f;
176+
177+
const auto& [x, y, z] = coords;
178+
179+
// Extract the actual histogram type (TH2 or TH3) held by the smart pointer
180+
using HistoType = typename std::decay_t<decltype(mapPtr)>::element_type;
181+
182+
// Compile-time branching: generates the exact correct function call
183+
if constexpr (std::is_same_v<HistoType, TH2>) {
184+
return mapPtr->GetBinContent(mapPtr->FindFixBin(y, z)); // 2D case only
185+
} else {
186+
return mapPtr->GetBinContent(mapPtr->FindFixBin(x, y, z)); // Full 3D case
187+
}
188+
},
189+
effMap);
190+
}
191+
192+
float interpolateEfficiency() const
193+
{
194+
return std::visit(
195+
[this](auto&& mapPtr) -> float {
196+
if (!mapPtr)
197+
return 1.0f;
198+
199+
const auto& [x, y, z] = coords;
200+
201+
using HistoType = typename std::decay_t<decltype(mapPtr)>::element_type;
202+
203+
if constexpr (std::is_same_v<HistoType, TH2>) {
204+
return mapPtr->Interpolate(y, z); // Native 2D interpolation
205+
} else {
206+
return mapPtr->Interpolate(x, y, z); // Native 3D interpolation
207+
}
208+
},
209+
effMap);
210+
}
152211
};
153212

154213
struct PhiStrangenessCorrelation {
@@ -195,6 +254,12 @@ struct PhiStrangenessCorrelation {
195254
// Configurable on pion pT bins
196255
Configurable<std::vector<double>> binspTPi{"binspTPi", {0.2, 0.3, 0.4, 0.5, 0.6, 0.8, 1.0, 1.2, 1.5, 2.0, 3.0}, "pT bin limits for pions"};
197256

257+
struct : ConfigurableGroup {
258+
Configurable<bool> doK0SCorrelation{"doK0SCorrelation", true, "Enable Phi-K0S correlation"};
259+
Configurable<bool> doXiCorrelation{"doXiCorrelation", false, "Enable Phi-Xi correlation"};
260+
Configurable<bool> doPionCorrelation{"doPionCorrelation", true, "Enable Phi-Pion correlation"};
261+
} activeCorrelationConfigs;
262+
198263
// Configurables for delta y selection
199264
struct : ConfigurableGroup {
200265
Configurable<int> nBinsY{"nBinsY", 20, "Number of bins in y axis"};
@@ -250,18 +315,22 @@ struct PhiStrangenessCorrelation {
250315
// Necessary service to retrieve efficiency maps from CCDB
251316
Service<ccdb::BasicCCDBManager> ccdb;
252317

253-
std::array<std::shared_ptr<TH3>, ParticleOfInterestSize> effMaps{};
318+
// std::shared_ptr<TH3> effMapPhi{};
319+
// std::array<std::shared_ptr<TH3>, kAssocPartSize> effMapsAssoc{};
320+
321+
EffMapPtr effMapPhi{};
322+
std::array<EffMapPtr, kAssocPartSize> effMapsAssoc{};
254323

255324
// Binning policy and axes for mixed event
256-
ConfigurableAxis axisVertexMixing{"axisVertexMixing", {20, -10, 10}, "Z vertex axis binning for mixing"};
257-
ConfigurableAxis axisCentralityMixing{"axisCentralityMixing", {20, 0, 100}, "Multiplicity percentil binning for mixing"};
325+
ConfigurableAxis axisVertexMixing{"axisVertexMixing", {10, -10.0f, 10.0f}, "Z vertex axis binning for mixing"};
326+
ConfigurableAxis axisCentralityMixing{"axisCentralityMixing", {VARIABLE_WIDTH, 0.0f, 1.0f, 5.0f, 10.0f, 15.0f, 20.0f, 30.0f, 40.0f, 50.0f, 70.0f, 100.0f}, "Multiplicity percentage binning for mixing"};
258327

259328
using BinningTypeVertexCent = ColumnBinningPolicy<aod::collision::PosZ, aod::cent::CentFT0M>;
260329
BinningTypeVertexCent binningOnVertexAndCent{{axisVertexMixing, axisCentralityMixing}, true};
261330

262331
static constexpr std::array<std::string_view, 2> phiMassRegionLabels{"Signal", "Sideband"};
263-
static constexpr std::array<std::string_view, ParticleOfInterestSize> particleOfInterestLabels{"Phi", "K0S", "Pion" /*"PionTPC", "PionTPCTOF"*/};
264-
static constexpr std::array<std::string_view, ParticleOfInterestSize - 1> assocParticleLabels{"K0S", "Pi"};
332+
// static constexpr std::array<std::string_view, ParticleOfInterestSize> particleOfInterestLabels{"Phi", "K0S", "Pion" /*"PionTPC", "PionTPCTOF"*/};
333+
static constexpr std::array<std::string_view, kAssocPartSize> assocParticleLabels{"K0S", "Xi", "Pi"};
265334

266335
// Light structures to store only the necessary information for the correlation analysis at MCGen level
267336
struct MiniParticle {
@@ -274,6 +343,7 @@ struct PhiStrangenessCorrelation {
274343
float multiplicity;
275344
std::vector<MiniParticle> phiParticles;
276345
std::vector<MiniParticle> k0sParticles;
346+
std::vector<MiniParticle> xiParticles;
277347
std::vector<MiniParticle> pionParticles;
278348
};
279349

@@ -348,21 +418,56 @@ struct PhiStrangenessCorrelation {
348418
ccdb->setLocalObjectValidityChecking();
349419
ccdb->setFatalWhenNull(false);
350420

351-
for (int i = 0; i < ParticleOfInterestSize; ++i) {
421+
/*for (int i = 0; i < ParticleOfInterestSize; ++i) {
352422
loadEfficiencyMapFromCCDB(static_cast<ParticleOfInterest>(i));
353-
}
423+
}*/
424+
loadEfficiencyMaps();
354425
}
355426

356427
eventBuffer.resize(binsMult->size() - 1);
357428
}
358429

359-
void loadEfficiencyMapFromCCDB(ParticleOfInterest poi)
430+
void fetchSingleEfficiencyMapFromCCDB(EffMapPtr& effMap, std::string_view particleName)
431+
{
432+
std::string path = fmt::format("{}{}", ccdbEfficiencyPath.value, particleName);
433+
434+
if (auto map3D = std::shared_ptr<TH3>(ccdb->get<TH3D>(path))) {
435+
effMap = map3D;
436+
LOG(info) << "Efficiency map (TH3) for " << particleName << " loaded from CCDB";
437+
return;
438+
}
439+
440+
if (auto map2D = std::shared_ptr<TH2>(ccdb->get<TH2D>(path))) {
441+
effMap = map2D;
442+
LOG(info) << "Efficiency map (TH2) for " << particleName << " loaded from CCDB";
443+
return;
444+
}
445+
446+
LOG(fatal) << "Could not load efficiency map (neither TH3 nor TH2) for " << particleName << " from CCDB!";
447+
}
448+
449+
void loadEfficiencyMaps()
450+
{
451+
// Always load the Trigger (Phi) map
452+
fetchSingleEfficiencyMapFromCCDB(effMapPhi, "Phi");
453+
454+
// Map the user configurations for the associated particles
455+
bool doAssocCorrelations[kAssocPartSize] = {activeCorrelationConfigs.doK0SCorrelation, activeCorrelationConfigs.doXiCorrelation, activeCorrelationConfigs.doPionCorrelation};
456+
457+
// Only load the associated maps that are explicitly enabled
458+
for (size_t i = 0; i < kAssocPartSize; ++i) {
459+
if (doAssocCorrelations[i])
460+
fetchSingleEfficiencyMapFromCCDB(effMapsAssoc[i], assocParticleLabels[i]);
461+
}
462+
}
463+
464+
/*void loadEfficiencyMapFromCCDB(ParticleOfInterest poi)
360465
{
361466
effMaps[poi] = std::shared_ptr<TH3>(ccdb->get<TH3D>(fmt::format("{}{}", ccdbEfficiencyPath.value, particleOfInterestLabels[poi])));
362467
if (!effMaps[poi])
363468
LOG(fatal) << "Could not load efficiency map for " << particleOfInterestLabels[poi] << "!";
364469
LOG(info) << "Efficiency map for " << particleOfInterestLabels[poi] << " loaded from CCDB";
365-
}
470+
}*/
366471

367472
// Compute weight based on efficiencies
368473
template <typename... BoundEffMaps>
@@ -397,10 +502,14 @@ struct PhiStrangenessCorrelation {
397502

398503
const std::array<std::pair<float, float>, 2> phiMassRegions = {phiConfigs.rangeMPhiSignal, phiConfigs.rangeMPhiSideband};
399504

505+
bool doAssocCorrelations[kAssocPartSize] = {activeCorrelationConfigs.doK0SCorrelation,
506+
activeCorrelationConfigs.doXiCorrelation,
507+
activeCorrelationConfigs.doPionCorrelation};
508+
400509
const bool applyK0sMassCut = (analysisMode == kDeltaYvsDeltaPhi) && k0sConfigs.selectK0sInSigRegion;
401510
const auto& [minMassK0s, maxMassK0s] = k0sConfigs.rangeMK0sSignal.value;
402511
auto isK0sValid = [&](const auto& k0s) {
403-
return !applyK0sMassCut || k0s.inMassRegion(minMassK0s, maxMassK0s);
512+
return (!applyEfficiency || k0s.pt() <= binspTK0S->back()) && (!applyK0sMassCut || k0s.inMassRegion(minMassK0s, maxMassK0s));
404513
};
405514

406515
const bool applyPionNSigmaCut = (analysisMode == kDeltaYvsDeltaPhi) && pionConfigs.selectPionInSigRegion;
@@ -409,33 +518,40 @@ struct PhiStrangenessCorrelation {
409518
const float& tofPIDThreshold = pionConfigs.tofPIDThreshold;
410519

411520
auto isPionValid = [&](const auto& pion) {
412-
return !applyPionNSigmaCut || pion.inNSigmaRegion(pidTPCMax, tofPIDThreshold, pidTOFMax);
521+
return (!applyEfficiency || pion.pt() <= binspTPi->back()) && (!applyPionNSigmaCut || pion.inNSigmaRegion(pidTPCMax, tofPIDThreshold, pidTOFMax));
413522
};
414523

415524
for (const auto& phiCand : phiCandidates) {
416-
float weightPhi = computeWeight(BoundEfficiencyMap(effMaps[Phi], multiplicity, phiCand.pt(), phiCand.y()));
525+
if (applyEfficiency && phiCand.pt() > binspTPhi->back())
526+
continue;
527+
528+
float weightPhi = computeWeight(BoundEfficiencyMap(effMapPhi, multiplicity, phiCand.pt(), phiCand.y()));
417529

418530
histos.fill(HIST("phi/h3PhiData"), multiplicity, phiCand.pt(), phiCand.m(), weightPhi);
419531

420532
auto processCorrelations = [&](auto fillK0S, auto fillPion) {
421-
// Loop over all reduced K0S candidates
422-
for (const auto& k0s : k0sReduced) {
423-
if (!isK0sValid(k0s))
424-
continue;
533+
if (doAssocCorrelations[kK0S]) {
534+
// Loop over all reduced K0S candidates
535+
for (const auto& k0s : k0sReduced) {
536+
if (!isK0sValid(k0s))
537+
continue;
425538

426-
float weightPhiK0S = computeWeight(BoundEfficiencyMap(effMaps[Phi], multiplicity, phiCand.pt(), phiCand.y()),
427-
BoundEfficiencyMap(effMaps[K0S], multiplicity, k0s.pt(), k0s.y()));
428-
fillK0S(k0s, weightPhiK0S);
539+
float weightPhiK0S = computeWeight(BoundEfficiencyMap(effMapPhi, multiplicity, phiCand.pt(), phiCand.y()),
540+
BoundEfficiencyMap(effMapsAssoc[kK0S], multiplicity, k0s.pt(), k0s.y()));
541+
fillK0S(k0s, weightPhiK0S);
542+
}
429543
}
430544

431-
// Loop over all primary pion candidates
432-
for (const auto& pionTrack : pionTracks) {
433-
if (!isPionValid(pionTrack))
434-
continue;
545+
if (doAssocCorrelations[kPion]) {
546+
// Loop over all primary pion candidates
547+
for (const auto& pionTrack : pionTracks) {
548+
if (!isPionValid(pionTrack))
549+
continue;
435550

436-
float weightPhiPion = computeWeight(BoundEfficiencyMap(effMaps[Phi], multiplicity, phiCand.pt(), phiCand.y()),
437-
BoundEfficiencyMap(effMaps[Pion], multiplicity, pionTrack.pt(), pionTrack.y()));
438-
fillPion(pionTrack, weightPhiPion);
551+
float weightPhiPion = computeWeight(BoundEfficiencyMap(effMapPhi, multiplicity, phiCand.pt(), phiCand.y()),
552+
BoundEfficiencyMap(effMapsAssoc[kPion], multiplicity, pionTrack.pt(), pionTrack.y()));
553+
fillPion(pionTrack, weightPhiPion);
554+
}
439555
}
440556
};
441557

@@ -533,7 +649,7 @@ struct PhiStrangenessCorrelation {
533649
const auto& [minMassK0s, maxMassK0s] = k0sConfigs.rangeMK0sSignal.value;
534650

535651
auto isK0sValid = [&](const auto& k0s) {
536-
return !applyK0sMassCut || k0s.inMassRegion(minMassK0s, maxMassK0s);
652+
return (!applyEfficiency || k0s.pt() <= binspTK0S->back()) && (!applyK0sMassCut || k0s.inMassRegion(minMassK0s, maxMassK0s));
537653
};
538654

539655
auto tuplePhiK0S = std::make_tuple(phiCandidates, k0sReduced);
@@ -544,12 +660,14 @@ struct PhiStrangenessCorrelation {
544660
float multiplicity = c1.centFT0M();
545661

546662
for (const auto& [phiCand, k0s] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(phiCands, k0sRed))) {
663+
if (applyEfficiency && phiCand.pt() > binspTPhi->back())
664+
continue;
547665
if (!isK0sValid(k0s))
548666
continue;
549667

550668
auto processCorrelations = [&](auto fillK0S) {
551-
float weightPhiK0S = computeWeight(BoundEfficiencyMap(effMaps[Phi], multiplicity, phiCand.pt(), phiCand.y()),
552-
BoundEfficiencyMap(effMaps[K0S], multiplicity, k0s.pt(), k0s.y()));
669+
float weightPhiK0S = computeWeight(BoundEfficiencyMap(effMapPhi, multiplicity, phiCand.pt(), phiCand.y()),
670+
BoundEfficiencyMap(effMapsAssoc[kK0S], multiplicity, k0s.pt(), k0s.y()));
553671
fillK0S(k0s, weightPhiK0S);
554672
};
555673

@@ -626,7 +744,7 @@ struct PhiStrangenessCorrelation {
626744
const float& tofPIDThreshold = pionConfigs.tofPIDThreshold;
627745

628746
auto isPionValid = [&](const auto& pion) {
629-
return !applyPionNSigmaCut || pion.inNSigmaRegion(pidTPCMax, tofPIDThreshold, pidTOFMax);
747+
return (!applyEfficiency || pion.pt() <= binspTPi->back()) && (!applyPionNSigmaCut || pion.inNSigmaRegion(pidTPCMax, tofPIDThreshold, pidTOFMax));
630748
};
631749

632750
auto tuplePhiPion = std::make_tuple(phiCandidates, pionTracks);
@@ -637,12 +755,14 @@ struct PhiStrangenessCorrelation {
637755
float multiplicity = c1.centFT0M();
638756

639757
for (const auto& [phiCand, piTrack] : o2::soa::combinations(o2::soa::CombinationsFullIndexPolicy(phiCands, piTracks))) {
758+
if (applyEfficiency && phiCand.pt() > binspTPhi->back())
759+
continue;
640760
if (!isPionValid(piTrack))
641761
continue;
642762

643763
auto processCorrelations = [&](auto fillPion) {
644-
float weightPhiPion = computeWeight(BoundEfficiencyMap(effMaps[Phi], multiplicity, phiCand.pt(), phiCand.y()),
645-
BoundEfficiencyMap(effMaps[Pion], multiplicity, piTrack.pt(), piTrack.y()));
764+
float weightPhiPion = computeWeight(BoundEfficiencyMap(effMapPhi, multiplicity, phiCand.pt(), phiCand.y()),
765+
BoundEfficiencyMap(effMapsAssoc[kPion], multiplicity, piTrack.pt(), piTrack.y()));
646766
fillPion(piTrack, weightPhiPion);
647767
};
648768

@@ -961,6 +1081,7 @@ struct PhiStrangenessCorrelation {
9611081

9621082
std::vector<MiniParticle> phiParticles;
9631083
std::vector<MiniParticle> k0sParticles;
1084+
std::vector<MiniParticle> xiParticles;
9641085
std::vector<MiniParticle> pionParticles;
9651086

9661087
auto inYAcceptance = [&](const auto& mcParticle) {
@@ -997,14 +1118,20 @@ struct PhiStrangenessCorrelation {
9971118
if (multBin < 0)
9981119
return;
9991120

1121+
bool doAssocCorrelations[kAssocPartSize] = {activeCorrelationConfigs.doK0SCorrelation,
1122+
activeCorrelationConfigs.doXiCorrelation,
1123+
activeCorrelationConfigs.doPionCorrelation};
1124+
10001125
// Same Event Correlations
1001-
std::vector<MiniParticle>* currentAssocParticles[] = {&k0sParticles, &pionParticles};
1126+
std::vector<MiniParticle>* currentAssocParticles[] = {&k0sParticles, &xiParticles, &pionParticles};
10021127

10031128
for (const auto& phiParticle : phiParticles) {
10041129
histos.fill(HIST("phi/h3PhiMCGen"), multiplicity, phiParticle.pt, phiParticle.y);
10051130

10061131
static_for<0, assocParticleLabels.size() - 1>([&](auto i_idx) {
10071132
constexpr unsigned int i = i_idx.value;
1133+
if (!doAssocCorrelations[i])
1134+
return;
10081135

10091136
for (const auto& assocParticle : *(currentAssocParticles[i])) {
10101137
histos.fill(HIST("phi") + HIST(assocParticleLabels[i]) + HIST("/h5Phi") + HIST(assocParticleLabels[i]) + HIST("ClosureMCGen"),
@@ -1017,12 +1144,14 @@ struct PhiStrangenessCorrelation {
10171144

10181145
// Mixed Event Correlations
10191146
for (const auto& pastEvent : eventBuffer[multBin]) {
1020-
const std::vector<MiniParticle>* pastAssocParticles[] = {&pastEvent.k0sParticles, &pastEvent.pionParticles};
1147+
const std::vector<MiniParticle>* pastAssocParticles[] = {&pastEvent.k0sParticles, &pastEvent.xiParticles, &pastEvent.pionParticles};
10211148

10221149
// Loop over past events in the same multiplicity bin and fill histograms with all combinations of current phi particles and past associated particles
10231150
for (const auto& phiParticle : phiParticles) {
10241151
static_for<0, assocParticleLabels.size() - 1>([&](auto i_idx) {
10251152
constexpr unsigned int i = i_idx.value;
1153+
if (!doAssocCorrelations[i])
1154+
return;
10261155

10271156
for (const auto& assocParticle : *(pastAssocParticles[i])) {
10281157
histos.fill(HIST("phi") + HIST(assocParticleLabels[i]) + HIST("/h5Phi") + HIST(assocParticleLabels[i]) + HIST("ClosureMCGenME"),
@@ -1038,6 +1167,7 @@ struct PhiStrangenessCorrelation {
10381167
currentEvent.multiplicity = multiplicity;
10391168
currentEvent.phiParticles = std::move(phiParticles);
10401169
currentEvent.k0sParticles = std::move(k0sParticles);
1170+
currentEvent.xiParticles = std::move(xiParticles);
10411171
currentEvent.pionParticles = std::move(pionParticles);
10421172

10431173
eventBuffer[multBin].push_front(std::move(currentEvent));

0 commit comments

Comments
 (0)