Skip to content

[PWGCF/FemtoUniverse] Efficiency correction in 3D #11270

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
zchochul merged 9 commits into from
May 21, 2025
Merged

[PWGCF/FemtoUniverse] Efficiency correction in 3D #11270

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
1be3775
fix: correct classification of material
davkk May 12, 2025
cbd8db0
feat: fill 3d histograms for efficiency correction
davkk May 15, 2025
7a90125
feat: macro with projection and working efficiency correction calcula…
davkk May 18, 2025
0c9edb8
fix: project the other way around
davkk May 19, 2025
8b0d824
feat: add 3d hist to loading from ccdb
davkk May 19, 2025
7067772
feat: use configurables for axis specifications
davkk May 19, 2025
5760877
feat: add flags to macro, add axes titles
davkk May 20, 2025
7489af0
feat: change variables configurable to string
davkk May 20, 2025
3ba9c08
fix: megalinter
davkk May 20, 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
188 changes: 156 additions & 32 deletions PWGCF/FemtoUniverse/Core/FemtoUniverseEfficiencyCorrection.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,15 +16,22 @@
#ifndef PWGCF_FEMTOUNIVERSE_CORE_FEMTOUNIVERSEEFFICIENCYCORRECTION_H_
#define PWGCF_FEMTOUNIVERSE_CORE_FEMTOUNIVERSEEFFICIENCYCORRECTION_H_

#include <CCDB/BasicCCDBManager.h>
#include <Framework/Configurable.h>
#include <Framework/ConfigurableKinds.h>
#include <Framework/HistogramRegistry.h>
#include <Framework/HistogramSpec.h>

#include <Framework/O2DatabasePDGPlugin.h>
#include <TH1.h>
#include <TH2.h>

#include <vector>
#include <string>
#include <algorithm>
#include <ranges>

#include "Framework/Configurable.h"
#include "CCDB/BasicCCDBManager.h"
#include "TH1.h"
#include "TH2.h"
#include "TH3.h"
#include "PWGCF/FemtoUniverse/DataModel/FemtoDerived.h"

namespace o2::analysis::femto_universe::efficiency_correction
{
Expand All @@ -34,39 +41,41 @@ enum ParticleNo : size_t {
};

template <size_t T>
concept isOneOrTwo = T == ParticleNo::ONE || T == ParticleNo::TWO;

template <typename T>
consteval auto getHistDim() -> int
{
if (std::is_same_v<T, TH1>)
return 1;
else if (std::is_same_v<T, TH2>)
return 2;
else if (std::is_same_v<T, TH3>)
return 3;
else
return -1;
}
concept IsOneOrTwo = T == ParticleNo::ONE || T == ParticleNo::TWO;

struct EffCorConfigurableGroup : framework::ConfigurableGroup {
framework::Configurable<bool> confEffCorApply{"confEffCorApply", false, "[Efficiency Correction] Should apply efficiency corrections"};
framework::Configurable<bool> confEffCorFillHist{"confEffCorFillHist", false, "[Efficiency Correction] Should fill histograms for efficiency corrections"};
framework::Configurable<std::string> confEffCorCCDBUrl{"confEffCorCCDBUrl", "http://alice-ccdb.cern.ch", "[Efficiency Correction] CCDB URL to use"};
framework::Configurable<std::string> confEffCorCCDBPath{"confEffCorCCDBPath", "", "[Efficiency Correction] CCDB path to histograms"};
framework::Configurable<std::vector<std::string>> confEffCorCCDBTimestamps{"confEffCorCCDBTimestamps", {}, "[Efficiency Correction] Timestamps of histograms in CCDB (0 can be used as a placeholder, e.g. when running subwagons)"};
framework::Configurable<std::string> confEffCorVariables{"confEffCorVariables", "pt", "[Efficiency Correction] Variables for efficiency correction histogram dimensions (available: 'pt'; 'pt,eta'; 'pt,mult'; 'pt,eta,mult')"};
};

template <typename HistType>
requires std::is_base_of_v<TH1, HistType>
class EfficiencyCorrection
{
public:
explicit EfficiencyCorrection(EffCorConfigurableGroup* config) : config(config) // o2-linter: disable=name/function-variable
{
}

auto init() -> void
auto init(framework::HistogramRegistry* registry, std::vector<framework::AxisSpec> axisSpecs) -> void
{
shouldFillHistograms = config->confEffCorFillHist;

histRegistry = registry;
if (shouldFillHistograms) {
for (const auto& suffix : histSuffix) {
auto path = std::format("{}/{}", histDirectory, suffix);
registry->add((path + "/hMCTruth").c_str(), "MCTruth; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
registry->add((path + "/hPrimary").c_str(), "Primary; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
registry->add((path + "/hSecondary").c_str(), "Secondary; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
registry->add((path + "/hMaterial").c_str(), "Material; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
registry->add((path + "/hFake").c_str(), "Fake; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
registry->add((path + "/hOther").c_str(), "Other; #it{p}_{T} (GeV/#it{c}); #it{#eta}; Mult", framework::kTH3F, axisSpecs);
}
}

ccdb.setURL(config->confEffCorCCDBUrl);
ccdb.setLocalObjectValidityChecking();
ccdb.setFatalWhenNull(false);
Expand All @@ -86,20 +95,123 @@ class EfficiencyCorrection
continue;
}

hLoaded[idx] = timestamp > 0 ? loadHistFromCCDB(timestamp) : nullptr;
if (timestamp > 0) {
switch (getDimensionFromVariables()) {
case 1:
hLoaded[idx] = loadHistFromCCDB<TH1>(timestamp);
break;
case 2:
hLoaded[idx] = loadHistFromCCDB<TH2>(timestamp);
break;
case 3:
hLoaded[idx] = loadHistFromCCDB<TH3>(timestamp);
break;
default:
LOGF(fatal, notify("Unknown configuration for efficiency variables"));
break;
}
} else {
hLoaded[idx] = nullptr;
}
}
}
}

template <typename... BinVars>
requires(sizeof...(BinVars) == getHistDim<HistType>())
auto getWeight(ParticleNo partNo, const BinVars&... binVars) const -> float
template <uint8_t N>
requires IsOneOrTwo<N>
void fillTruthHist(auto particle)
{
if (!shouldFillHistograms) {
return;
}

histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hMCTruth"),
particle.pt(),
particle.eta(),
particle.fdCollision().multV0M());
}

template <uint8_t N>
requires IsOneOrTwo<N>
void fillRecoHist(auto particle, int particlePDG)
{
if (!shouldFillHistograms) {
return;
}

if (!particle.has_fdMCParticle()) {
return;
}

auto mcParticle = particle.fdMCParticle();

if (mcParticle.pdgMCTruth() == particlePDG) {
switch (mcParticle.partOriginMCTruth()) {
case (o2::aod::femtouniverse_mc_particle::kPrimary):
histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hPrimary"),
mcParticle.pt(),
mcParticle.eta(),
particle.fdCollision().multV0M());
break;

case (o2::aod::femtouniverse_mc_particle::kDaughter):
case (o2::aod::femtouniverse_mc_particle::kDaughterLambda):
case (o2::aod::femtouniverse_mc_particle::kDaughterSigmaplus):
histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hSecondary"),
mcParticle.pt(),
mcParticle.eta(),
particle.fdCollision().multV0M());
break;

case (o2::aod::femtouniverse_mc_particle::kMaterial):
histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hMaterial"),
mcParticle.pt(),
mcParticle.eta(),
particle.fdCollision().multV0M());
break;

case (o2::aod::femtouniverse_mc_particle::kFake):
histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hFake"),
mcParticle.pt(),
mcParticle.eta(),
particle.fdCollision().multV0M());
break;

default:
histRegistry->fill(HIST(histDirectory) + HIST("/") + HIST(histSuffix[N - 1]) + HIST("/hOther"),
mcParticle.pt(),
mcParticle.eta(),
particle.fdCollision().multV0M());
break;
}
}
}

auto getWeight(ParticleNo partNo, auto particle) -> float
{
auto weight = 1.0f;
auto hWeights = hLoaded[partNo - 1];

if (shouldApplyCorrection && hWeights) {
auto bin = hWeights->FindBin(binVars...);
auto dim = static_cast<size_t>(hWeights->GetDimension());
if (dim != getDimensionFromVariables()) {
LOGF(fatal, notify("Histogram \"%s\" has wrong dimension %d != %d"), config->confEffCorCCDBPath.value, dim, config->confEffCorVariables.value.size());
return weight;
}

auto bin = -1;
if (config->confEffCorVariables.value == "pt") {
bin = hLoaded[partNo - 1]->FindBin(particle.pt());
} else if (config->confEffCorVariables.value == "pt,eta") {
bin = hLoaded[partNo - 1]->FindBin(particle.pt(), particle.eta());
} else if (config->confEffCorVariables.value == "pt,mult") {
bin = hLoaded[partNo - 1]->FindBin(particle.pt(), particle.fdCollision().multV0M());
} else if (config->confEffCorVariables.value == "pt,eta,mult") {
bin = hLoaded[partNo - 1]->FindBin(particle.pt(), particle.eta(), particle.fdCollision().multV0M());
} else {
LOGF(fatal, notify("Unknown configuration for efficiency variables"));
return weight;
}
weight = hWeights->GetBinContent(bin);
}

Expand All @@ -112,7 +224,7 @@ class EfficiencyCorrection
return fmt::format("[EFFICIENCY CORRECTION] {}", msg);
}

static auto isHistEmpty(HistType* hist) -> bool
static auto isHistEmpty(TH1* hist) -> bool
{
if (!hist) {
return true;
Expand All @@ -125,9 +237,10 @@ class EfficiencyCorrection
return true;
}

auto loadHistFromCCDB(const int64_t timestamp) const -> HistType*
template <typename H>
auto loadHistFromCCDB(const int64_t timestamp) const -> H*
{
auto hWeights = ccdb.getForTimeStamp<HistType>(config->confEffCorCCDBPath, timestamp);
auto hWeights = ccdb.getForTimeStamp<H>(config->confEffCorCCDBPath, timestamp);
if (!hWeights || hWeights->IsZombie()) {
LOGF(error, notify("Could not load histogram \"%s/%ld\""), config->confEffCorCCDBPath.value, timestamp);
return nullptr;
Expand All @@ -141,12 +254,23 @@ class EfficiencyCorrection
return hWeights;
}

auto getDimensionFromVariables() -> size_t
{
auto parts = std::views::split(config->confEffCorVariables.value, ',');
return std::ranges::distance(parts);
}

EffCorConfigurableGroup* config{};

bool shouldApplyCorrection = false;
bool shouldApplyCorrection{false};
bool shouldFillHistograms{false};

o2::ccdb::BasicCCDBManager& ccdb{o2::ccdb::BasicCCDBManager::instance()};
std::array<HistType*, 2> hLoaded{};
std::array<TH1*, 2> hLoaded{};

framework::HistogramRegistry* histRegistry{};
static constexpr std::string_view histDirectory{"EfficiencyCorrection"};
static constexpr std::string_view histSuffix[2]{"one", "two"};
};

} // namespace o2::analysis::femto_universe::efficiency_correction
Expand Down
Loading