Skip to content
Merged
Changes from all commits
Commits
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
130 changes: 128 additions & 2 deletions PWGEM/PhotonMeson/Tasks/photonhbt.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -149,6 +149,11 @@ struct Photonhbt {
float u{999.f};
};

struct CrossObs {
std::array<float, 2> mee{999.f, 999.f}; // crossed m(e+e-) per combo
std::array<float, 2> dist{999.f, 999.f}; // crossed line-line distance (cm)
};

struct TruthGamma {
int id = -1, posId = -1, negId = -1;
float eta = 0.f, phi = 0.f, pt = 0.f;
Expand Down Expand Up @@ -336,6 +341,15 @@ struct Photonhbt {
Configurable<float> cfgUCut{"cfgUCut", 0.f, "min u (sigma units); pairs with u < cut rejected in SE AND ME"};
} pairsep;

struct : ConfigurableGroup {
std::string prefix = "crosspair_group";
Configurable<bool> cfgDoCrossPairQA{"cfgDoCrossPairQA", false, "fill crossed-hypothesis QA (mee_cross, dist_cross vs qinv) for pairs passing all other pair cuts; SE and ME"};
Configurable<bool> cfgDoCrossPairCut{"cfgDoCrossPairCut", false, "reject pairs with a photon-like crossed combination; applied to SE AND ME"};
Configurable<float> cfgCrossMaxMee{"cfgCrossMaxMee", 0.04f, "veto: crossed m_ee below this (GeV/c^2)"};
Configurable<float> cfgCrossMaxDist{"cfgCrossMaxDist", 1.0f, "veto: crossed line-line distance below this (cm)"};
Configurable<float> cfgCrossMaxQinvQA{"cfgCrossMaxQinvQA", 0.3f, "fill crossed QA sparses only below this qinv"};
} crosspair;

struct : ConfigurableGroup {
std::string prefix = "eventcut_group";
Configurable<float> cfgZvtxMin{"cfgZvtxMin", -10.f, "min. Zvtx"};
Expand Down Expand Up @@ -656,6 +670,7 @@ struct Photonhbt {
addEventHistograms();
addPairCFHistograms();
addPairSepHistograms();
addCrossPairHistograms();
addSinglePhotonQAHistograms();
addPairQAHistograms();

Expand Down Expand Up @@ -704,8 +719,6 @@ struct Photonhbt {

fRegistryCF.add("Pair/mix/hDiffBC", "diff. global BC in mixed event;|BC_{current}-BC_{mixed}|", kTH1D, {{10001, -0.5, 10000.5}}, true);
}

// ─── CF: PairSep ────────────────────────
// ─── CF: PairSep ────────────────────────
void addPairSepHistograms()
{
Expand All @@ -717,6 +730,21 @@ struct Photonhbt {
}
}

void addCrossPairHistograms()
{
if (!(crosspair.cfgDoCrossPairQA.value || crosspair.cfgDoCrossPairCut.value)) {
return;
}
for (const auto& sm : {std::string("Pair/same/CrossPair/"), std::string("Pair/mix/CrossPair/")}) {
fRegistryCF.add((sm + "hSparse_Mee_Dist_Qinv").c_str(),
"crossed-hypothesis QA;m_{ee}^{cross} (GeV/c^{2});crossed line distance (cm);q_{inv} (GeV/c)",
kTHnSparseF, {{100, 0.f, 0.1f}, {100, 0.f, 10.f}, {60, 0.f, 0.3f}}, true);
auto h = fRegistryCF.add<TH1>((sm + "hVetoCounter").c_str(), "crossed-pair veto;;pairs", kTH1D, {{2, -0.5f, 1.5f}}, true);
h->GetXaxis()->SetBinLabel(1, "evaluated");
h->GetXaxis()->SetBinLabel(2, "vetoed");
}
}

// ─── pairQA: single-photon QA per cut step ─────────────────────────────────
void addSinglePhotonQAHistograms()
{
Expand Down Expand Up @@ -1246,6 +1274,72 @@ struct Photonhbt {
return s.u > pairsep.cfgUCut.value;
}

[[nodiscard]] static float lineLineDistance(std::array<float, 3> const& o1, std::array<float, 3> const& d1,
std::array<float, 3> const& o2, std::array<float, 3> const& d2)
{
const std::array<float, 3> w = {o2[0] - o1[0], o2[1] - o1[1], o2[2] - o1[2]};
const std::array<float, 3> n = {d1[1] * d2[2] - d1[2] * d2[1],
d1[2] * d2[0] - d1[0] * d2[2],
d1[0] * d2[1] - d1[1] * d2[0]};
const float nMag = std::sqrt(n[0] * n[0] + n[1] * n[1] + n[2] * n[2]);
const float d1Mag2 = d1[0] * d1[0] + d1[1] * d1[1] + d1[2] * d1[2];
if (nMag < 1e-6f * d1Mag2) { // o2-linter: disable=magic-number (near-parallel: point-to-line distance)
const float wDotD = (w[0] * d1[0] + w[1] * d1[1] + w[2] * d1[2]) / d1Mag2;
const std::array<float, 3> perp = {w[0] - wDotD * d1[0], w[1] - wDotD * d1[1], w[2] - wDotD * d1[2]};
return std::sqrt(perp[0] * perp[0] + perp[1] * perp[1] + perp[2] * perp[2]);
}
return std::fabs(w[0] * n[0] + w[1] * n[1] + w[2] * n[2]) / nMag;
}

[[nodiscard]] CrossObs computeCrossObs(PhotonWithLegs const& a, PhotonWithLegs const& b) const
{
constexpr float kMe = 0.000510999f; // electron mass, GeV/c^2
CrossObs c;
auto legVec = [](PhotonWithLegs const& p, int i) {
return ROOT::Math::PtEtaPhiMVector(p.fLegPt[i], p.fLegEta[i], p.fLegPhi[i], kMe);
};
auto legDir = [](PhotonWithLegs const& p, int i) -> std::array<float, 3> {
return {p.fLegPt[i] * std::cos(p.fLegPhi[i]), p.fLegPt[i] * std::sin(p.fLegPhi[i]),
p.fLegPt[i] * std::sinh(p.fLegEta[i])};
};
// Leg index 0 = e+, 1 = e- (makePhotonWithLegs filling order).
// Combo 0: a's e+ with b's e-. Combo 1: b's e+ with a's e-.
for (int combo = 0; combo < 2; ++combo) { // o2-linter: disable=magic-number (combinations of fake photons)
PhotonWithLegs const& pPos = (combo == 0) ? a : b;
PhotonWithLegs const& pEle = (combo == 0) ? b : a;
c.mee[combo] = static_cast<float>((legVec(pPos, 0) + legVec(pEle, 1)).M());
c.dist[combo] = lineLineDistance({pPos.fVx, pPos.fVy, pPos.fVz}, legDir(pPos, 0),
{pEle.fVx, pEle.fVy, pEle.fVz}, legDir(pEle, 1));
}
return c;
}

[[nodiscard]] inline bool passCrossPairVeto(CrossObs const& c) const
{
if (!crosspair.cfgDoCrossPairCut.value) {
return true;
}
for (int i = 0; i < 2; ++i) { // o2-linter: disable=magic-number (combinations of fake photons)
if (c.mee[i] < crosspair.cfgCrossMaxMee.value && c.dist[i] < crosspair.cfgCrossMaxDist.value) {
return false;
}
}
return true;
}

template <int ev_id>
inline void fillCrossPair(CrossObs const& c, float qinv, bool vetoed)
{
constexpr auto dir = (ev_id == 0) ? "Pair/same/CrossPair/" : "Pair/mix/CrossPair/";
fRegistryCF.fill(HIST(dir) + HIST("hVetoCounter"), vetoed ? 1.0 : 0.0);
if (qinv > crosspair.cfgCrossMaxQinvQA.value) {
return;
}
for (int i = 0; i < 2; ++i) { // o2-linter: disable=magic-number (combinations of fake photons)
fRegistryCF.fill(HIST(dir) + HIST("hSparse_Mee_Dist_Qinv"), c.mee[i], c.dist[i], qinv);
}
}

[[nodiscard]] float uTrue(TruthGamma const& g1, TruthGamma const& g2) const
{
const float sE = pairsep.cfgUSigEta.value, sP = pairsep.cfgUSigPhi.value;
Expand Down Expand Up @@ -1913,6 +2007,14 @@ struct Photonhbt {
if (isInsideEllipse(obs.deta, obs.dphi)) {
continue;
}
if (crosspair.cfgDoCrossPairQA.value || crosspair.cfgDoCrossPairCut.value) {
const auto cross = computeCrossObs(pwl1, pwl2);
const bool vetoed = !passCrossPairVeto(cross);
fillCrossPair<0>(cross, obs.qinv, vetoed);
if (vetoed) {
continue;
}
}

// ───after pair cuts ──────────────────────────────────────
fillPairSep<0, true>(sep, obs);
Expand Down Expand Up @@ -2000,6 +2102,14 @@ struct Photonhbt {
if (isInsideEllipse(obs.deta, obs.dphi)) {
continue;
}
if (crosspair.cfgDoCrossPairQA.value || crosspair.cfgDoCrossPairCut.value) {
const auto cross = computeCrossObs(g1, g2);
const bool vetoed = !passCrossPairVeto(cross);
fillCrossPair<1>(cross, obs.qinv, vetoed);
if (vetoed) {
continue;
}
}

fillPairSep<1, true>(sep, obs);
if (doQA) {
Expand Down Expand Up @@ -2122,6 +2232,14 @@ struct Photonhbt {
if (isInsideEllipse(obs.deta, obs.dphi)) {
continue;
}
if (crosspair.cfgDoCrossPairQA.value || crosspair.cfgDoCrossPairCut.value) {
const auto cross = computeCrossObs(pwl1, pwl2);
const bool vetoed = !passCrossPairVeto(cross);
fillCrossPair<0>(cross, obs.qinv, vetoed);
if (vetoed) {
continue;
}
}

// ─── after pair cuts ──────────────────────────────────────
fillPairSep<0, true>(sep, obs);
Expand Down Expand Up @@ -2260,6 +2378,14 @@ struct Photonhbt {
if (isInsideEllipse(obs.deta, obs.dphi)) {
continue;
}
if (crosspair.cfgDoCrossPairQA.value || crosspair.cfgDoCrossPairCut.value) {
const auto cross = computeCrossObs(g1, g2);
const bool vetoed = !passCrossPairVeto(cross);
fillCrossPair<1>(cross, obs.qinv, vetoed);
if (vetoed) {
continue;
}
}

// ─after pair cuts ──────────────────────────────────
fillPairSep<1, true>(sep, obs);
Expand Down
Loading