diff --git a/PWGLF/Tasks/Resonances/doublephimeson.cxx b/PWGLF/Tasks/Resonances/doublephimeson.cxx index c9108099e5b..64b3d5a9b54 100644 --- a/PWGLF/Tasks/Resonances/doublephimeson.cxx +++ b/PWGLF/Tasks/Resonances/doublephimeson.cxx @@ -77,6 +77,8 @@ struct doublephimeson { ConfigurableAxis configThnAxisDeltaR{"configThnAxisDeltaR", {200, 0.0, 2.0}, "#it{k}^{*} (GeV/#it{c})"}; ConfigurableAxis configThnAxisCosTheta{"configThnAxisCosTheta", {160, 0.0, 3.2}, "cos #theta{*}"}; ConfigurableAxis configThnAxisNumPhi{"configThnAxisNumPhi", {101, -0.5, 100.5}, "cos #theta{*}"}; + ConfigurableAxis configThnAxisDeltaPt{"configThnAxisDeltaPt", {100, 0.0, 1.0}, "delta pt"}; + Configurable maxDeltaMPhi{"maxDeltaMPhi", 0.01f, "Delta-m cut on the two phi masses: sqrt((m1-mPDG)^2 + (m2-mPDG)^2) < maxDeltaMPhi (GeV/c^2)"}; // Initialize the ananlysis task void init(o2::framework::InitContext&) @@ -95,15 +97,16 @@ struct doublephimeson { histos.add("hDeltaRkaonplus", "hDeltaRkaonplus", kTH1F, {{800, 0.0, 8.0}}); histos.add("hDeltaRkaonminus", "hDeltaRkaonminus", kTH1F, {{800, 0.0, 8.0}}); + const AxisSpec thnAxisdeltapt{configThnAxisDeltaPt, "Delta pt"}; const AxisSpec thnAxisInvMass{configThnAxisInvMass, "#it{M} (GeV/#it{c}^{2})"}; + const AxisSpec thnAxisPt{configThnAxisPt, "#it{p}_{T} (GeV/#it{c})"}; const AxisSpec thnAxisInvMassPhi{configThnAxisInvMassPhi, "#it{M} (GeV/#it{c}^{2})"}; const AxisSpec thnAxisInvMassDeltaPhi{configThnAxisInvMassDeltaPhi, "#it{M} (GeV/#it{c}^{2})"}; - const AxisSpec thnAxisPt{configThnAxisPt, "#it{p}_{T} (GeV/#it{c})"}; const AxisSpec thnAxisDeltaR{configThnAxisDeltaR, "#Delta R)"}; const AxisSpec thnAxisCosTheta{configThnAxisCosTheta, "cos #theta"}; const AxisSpec thnAxisNumPhi{configThnAxisNumPhi, "Number of phi meson"}; - histos.add("SEMassUnlike", "SEMassUnlike", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisDeltaR, thnAxisInvMassDeltaPhi, thnAxisNumPhi}); + histos.add("SEMassUnlike", "SEMassUnlike", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisdeltapt, thnAxisPt, thnAxisDeltaR, thnAxisInvMassDeltaPhi, thnAxisNumPhi}); // histos.add("SEMassLike", "SEMassLike", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisDeltaR, thnAxisInvMassPhi, thnAxisInvMassPhi, thnAxisNumPhi}); histos.add("MEMassUnlike", "MEMassUnlike", HistType::kTHnSparseF, {thnAxisInvMass, thnAxisPt, thnAxisDeltaR, thnAxisInvMassDeltaPhi}); } @@ -179,6 +182,92 @@ struct doublephimeson { bool selectionPID(float nsigmaTPC, float nsigmaTOF, int TOFHit, int PIDStrategy, float ptcand) { + + if (PIDStrategy == 100) { + if (ptcand < 1.2) { + if (TOFHit == 1 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.5) { + return true; + } else if (TOFHit != 1) { + if (ptcand < 0.5 && nsigmaTPC > -2.0 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.5 && ptcand < 0.6 && nsigmaTPC > -1.5 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.6 && ptcand < 0.7 && nsigmaTPC > -1.0 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 0.7 && ptcand < 0.8 && nsigmaTPC > -0.4 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 0.8 && ptcand < 1.0 && nsigmaTPC > 0.0 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 1.0 && ptcand < 1.2 && nsigmaTPC > -2.5 && nsigmaTPC < 0.5) { + return true; + } + } + } else { + if (TOFHit == 1 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.0) { + return true; + } else if (TOFHit != 1 && nsigmaTPC > -2.5 && nsigmaTPC < 1.0) { + return true; + } + } + } + + if (PIDStrategy == 101) { + if (ptcand < 1.0) { + if (TOFHit == 1 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.5) { + return true; + } else if (TOFHit != 1) { + if (ptcand < 0.5 && nsigmaTPC > -2.0 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.5 && ptcand < 0.6 && nsigmaTPC > -1.5 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.6 && ptcand < 0.7 && nsigmaTPC > -1.0 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 0.7 && ptcand < 0.8 && nsigmaTPC > -0.4 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 0.8 && ptcand < 1.0 && nsigmaTPC > 0.0 && nsigmaTPC < 2.0) { + return true; + } + } + } else if (ptcand >= 1.0 && ptcand < 2.0 && TOFHit == 1 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.5) { + return true; + } else if (ptcand > 2.0) { + if (TOFHit == 1 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.0) { + return true; + } else if (TOFHit != 1 && nsigmaTPC > -2.5 && nsigmaTPC < 1.0) { + return true; + } + } + } + + if (PIDStrategy == 102) { + if (TOFHit != 1) { + if (ptcand < 0.5 && nsigmaTPC > -2.0 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.5 && ptcand < 0.6 && nsigmaTPC > -1.5 && nsigmaTPC < 2.5) { + return true; + } + if (ptcand >= 0.6 && ptcand < 0.7 && nsigmaTPC > -1.0 && nsigmaTPC < 2.0) { + return true; + } + if (ptcand >= 2.2 && nsigmaTPC > -2.5 && nsigmaTPC < 1.0) { + return true; + } + } + if (TOFHit == 1 && ptcand > 0.4 && std::sqrt(nsigmaTOF * nsigmaTOF + nsigmaTPC * nsigmaTPC) < 2.0) { + return true; + } + } + // optimized TPC TOF if (PIDStrategy == 0) { if (ptcand < 0.4) { @@ -428,10 +517,10 @@ struct doublephimeson { continue; } if (!isDeep) { - histos.fill(HIST("SEMassUnlike"), exotic.M(), exotic.Pt(), deltaR, deltam, phimult); + histos.fill(HIST("SEMassUnlike"), exotic.M(), std::abs(Phid1.Pt() - Phid2.Pt()) / exotic.Pt(), exotic.Pt(), deltaR, deltam, phimult); } if (isDeep) { - histos.fill(HIST("SEMassUnlike"), exotic.M(), exotic.Pt(), deltaR, deltam, phimult); + histos.fill(HIST("SEMassUnlike"), exotic.M(), std::abs(Phid1.Pt() - Phid2.Pt()) / exotic.Pt(), exotic.Pt(), deltaR, deltam, phimult); } } } @@ -609,13 +698,13 @@ struct doublephimeson { (d4trackid.at(i5) == d3trackid.at(i6) || d4trackid.at(i5) == d4trackid.at(i6))) { if (deltam2 < deltam1 && deltaRkaonplus2 > daughterDeltaR && deltaRkaonminus2 > daughterDeltaR) { - histos.fill(HIST("SEMassUnlike"), exotic2.M(), exotic2.Pt(), deltaR2, deltam2, phimult); + histos.fill(HIST("SEMassUnlike"), exotic2.M(), std::abs(exotic2phi1.Pt() - exotic2phi2.Pt()) / exotic2.Pt(), exotic2.Pt(), deltaR2, deltam2, phimult); // LOGF(info, "Fill exotic Id %d which is pair of Id %d", i6, i5); } else { - histos.fill(HIST("SEMassUnlike"), exotic1.M(), exotic1.Pt(), deltaR1, deltam1, phimult); + histos.fill(HIST("SEMassUnlike"), exotic1.M(), std::abs(exotic2phi1.Pt() - exotic2phi2.Pt()) / exotic1.Pt(), exotic1.Pt(), deltaR1, deltam1, phimult); } } else { - histos.fill(HIST("SEMassUnlike"), exotic1.M(), exotic1.Pt(), deltaR1, deltam1, phimult); + histos.fill(HIST("SEMassUnlike"), exotic1.M(), std::abs(exotic2phi1.Pt() - exotic2phi2.Pt()) / exotic1.Pt(), exotic1.Pt(), deltaR1, deltam1, phimult); } } } @@ -645,219 +734,166 @@ struct doublephimeson { continue; } - histos.fill(HIST("SEMassUnlike"), exotic1.M(), exotic1.Pt(), deltaR1, deltam1, phimult); + histos.fill(HIST("SEMassUnlike"), exotic1.M(), std::abs(exotic1phi1.Pt() - exotic1phi2.Pt()) / exotic1.Pt(), exotic1.Pt(), deltaR1, deltam1, phimult); } } } PROCESS_SWITCH(doublephimeson, processopti, "Process Optimized same event", false); - void processopti2(aod::RedPhiEvents::iterator const& collision, aod::PhiTracks const& phitracks) + void processopti3(aod::RedPhiEvents::iterator const& collision, aod::PhiTracks const& phitracks) { - if (additionalEvsel && (collision.numPos() < 2 || collision.numNeg() < 2)) { + if (additionalEvsel && (collision.numPos() < 2 || collision.numNeg() < 2)) return; - } - // === Helpers === - // PDG phi mass (centralized constant, easy to vary for systematics) + // --- φ multiplicity with PID --- + int phimult = 0; + for (auto const& t : phitracks) { + if (t.phiMass() < minPhiMass1 || t.phiMass() > maxPhiMass1) + continue; + const double kpluspt = std::hypot(t.phid1Px(), t.phid1Py()); + const double kminuspt = std::hypot(t.phid2Px(), t.phid2Py()); + if (kpluspt > maxKaonPt || kminuspt > maxKaonPt) + continue; + if (!selectionPID(t.phid1TPC(), t.phid1TOF(), t.phid1TOFHit(), strategyPID1, kpluspt)) + continue; + if (!selectionPID(t.phid2TPC(), t.phid2TOF(), t.phid2TOFHit(), strategyPID2, kminuspt)) + continue; + ++phimult; + } + if (phimult < 2) + return; + // --- helpers --- constexpr double mPhiPDG = 1.019461; // GeV/c^2 - // ΔR with proper Δφ wrapping (−π, π] + const auto deltaMPhi = [=](double m1, double m2) { + const double d1 = m1 - mPhiPDG, d2 = m2 - mPhiPDG; + return std::sqrt(d1 * d1 + d2 * d2); + }; + const auto deltaR = [](double phi1, double eta1, double phi2, double eta2) { const double dphi = TVector2::Phi_mpi_pi(phi1 - phi2); const double deta = eta1 - eta2; return std::sqrt(dphi * dphi + deta * deta); }; - // Combined mass offset from PDG (tie-breaker metric) - const auto deltam = [=](double m1, double m2) { - return std::sqrt((m1 - mPhiPDG) * (m1 - mPhiPDG) + (m2 - mPhiPDG) * (m2 - mPhiPDG)); - }; - - // Anti-merging/duplication veto: - // Require same-sign kaons from the two φ's to be separated in (η,φ). - // This protects against split/merged tracks and near-duplicate topologies. + // anti-merging: same-sign kaons must be separated const auto passDaughterDR = [&](const ROOT::Math::PtEtaPhiMVector& kplusA, const ROOT::Math::PtEtaPhiMVector& kplusB, const ROOT::Math::PtEtaPhiMVector& kminusA, - const ROOT::Math::PtEtaPhiMVector& kminusB) { + const ROOT::Math::PtEtaPhiMVector& kminusB, + double thr) { const double dRkplus = deltaR(kplusA.Phi(), kplusA.Eta(), kplusB.Phi(), kplusB.Eta()); const double dRkminus = deltaR(kminusA.Phi(), kminusA.Eta(), kminusB.Phi(), kminusB.Eta()); histos.fill(HIST("hDeltaRkaonplus"), dRkplus); histos.fill(HIST("hDeltaRkaonminus"), dRkminus); - return (dRkplus > daughterDeltaR) && (dRkminus > daughterDeltaR); - // If later needed, make pT-aware: - // const double thr = std::max(0.01, daughterDeltaR - 0.002 * std::min(10.0, exoticPt)); - // return (dRkplus > thr) && (dRkminus > thr); + return (dRkplus > thr) && (dRkminus > thr); }; - - // === Phi multiplicity (uses PID1 for d1, PID2 for d2) === - int phimult = 0; - for (auto const& t : phitracks) { - if (t.phiMass() < minPhiMass1 || t.phiMass() > maxPhiMass1) - continue; - const double kpluspt = std::hypot(t.phid1Px(), t.phid1Py()); - const double kminuspt = std::hypot(t.phid2Px(), t.phid2Py()); - if (kpluspt > maxKaonPt) - continue; - if (kminuspt > maxKaonPt) - continue; - if (!selectionPID(t.phid1TPC(), t.phid1TOF(), t.phid1TOFHit(), strategyPID1, kpluspt)) - continue; // d1 → PID1 - if (!selectionPID(t.phid2TPC(), t.phid2TOF(), t.phid2TOFHit(), strategyPID2, kminuspt)) - continue; // d2 → PID2 - ++phimult; - } - - // === Collect candidates first === - std::vector exoticres, phi1v, phi2v, kplus1, kplus2, kminus1, kminus2; - std::vector d1id, d2id, d3id, d4id; - - const auto n = phitracks.size(); - exoticres.reserve(n); - phi1v.reserve(n); - phi2v.reserve(n); - kplus1.reserve(n); - kplus2.reserve(n); - kminus1.reserve(n); - kminus2.reserve(n); - d1id.reserve(n); - d2id.reserve(n); - d3id.reserve(n); - d4id.reserve(n); + // --- collect candidates once --- + std::vector pairV, phi1V, phi2V, kplus1V, kplus2V, kminus1V, kminus2V; for (auto const& t1 : phitracks) { - // Per-φ selection for the first φ const double kplus1pt = std::hypot(t1.phid1Px(), t1.phid1Py()); const double kminus1pt = std::hypot(t1.phid2Px(), t1.phid2Py()); + // PID QA before histos.fill(HIST("hnsigmaTPCTOFKaonBefore"), t1.phid1TPC(), t1.phid1TOF(), kplus1pt); histos.fill(HIST("hnsigmaTPCKaonPlusBefore"), t1.phid1TPC(), kplus1pt); histos.fill(HIST("hnsigmaTPCKaonMinusBefore"), t1.phid2TPC(), kminus1pt); - if (kplus1pt > maxKaonPt) - continue; - if (kminus1pt > maxKaonPt) + if (kplus1pt > maxKaonPt || kminus1pt > maxKaonPt) continue; if (!selectionPID(t1.phid1TPC(), t1.phid1TOF(), t1.phid1TOFHit(), strategyPID1, kplus1pt)) continue; if (!selectionPID(t1.phid2TPC(), t1.phid2TOF(), t1.phid2TOFHit(), strategyPID2, kminus1pt)) continue; - // Set vectors BEFORE QA fills - Phid1.SetXYZM(t1.phiPx(), t1.phiPy(), t1.phiPz(), t1.phiMass()); - Phi1kaonplus.SetXYZM(t1.phid1Px(), t1.phid1Py(), t1.phid1Pz(), 0.493); - Phi1kaonminus.SetXYZM(t1.phid2Px(), t1.phid2Py(), t1.phid2Pz(), 0.493); + TLorentzVector phi1, k1p, k1m; + phi1.SetXYZM(t1.phiPx(), t1.phiPy(), t1.phiPz(), t1.phiMass()); + k1p.SetXYZM(t1.phid1Px(), t1.phid1Py(), t1.phid1Pz(), 0.493); + k1m.SetXYZM(t1.phid2Px(), t1.phid2Py(), t1.phid2Pz(), 0.493); + + // φ mass windows + if (t1.phiMass() < minPhiMass1 || t1.phiMass() > maxPhiMass1) + continue; - // PID QA + // PID QA after histos.fill(HIST("hnsigmaTPCTOFKaon"), t1.phid1TPC(), t1.phid1TOF(), kplus1pt); histos.fill(HIST("hnsigmaTPCKaonPlus"), t1.phid1TPC(), kplus1pt); histos.fill(HIST("hnsigmaTPCKaonMinus"), t1.phid2TPC(), kminus1pt); - histos.fill(HIST("hPhiMass"), Phid1.M(), Phid1.Pt()); + histos.fill(HIST("hPhiMass"), phi1.M(), phi1.Pt()); const auto id1 = t1.index(); for (auto const& t2 : phitracks) { const auto id2 = t2.index(); if (id2 <= id1) - continue; // unique unordered pairs + continue; - // Per-φ selection for the second φ const double kplus2pt = std::hypot(t2.phid1Px(), t2.phid1Py()); const double kminus2pt = std::hypot(t2.phid2Px(), t2.phid2Py()); - if (kplus2pt > maxKaonPt) - continue; - if (kminus2pt > maxKaonPt) + if (kplus2pt > maxKaonPt || kminus2pt > maxKaonPt) continue; if (!selectionPID(t2.phid1TPC(), t2.phid1TOF(), t2.phid1TOFHit(), strategyPID1, kplus2pt)) continue; if (!selectionPID(t2.phid2TPC(), t2.phid2TOF(), t2.phid2TOFHit(), strategyPID2, kminus2pt)) continue; - // Disallow shared same-sign daughters between the two φ's + // block shared same-sign daughters if ((t1.phid1Index() == t2.phid1Index()) || (t1.phid2Index() == t2.phid2Index())) continue; - Phid2.SetXYZM(t2.phiPx(), t2.phiPy(), t2.phiPz(), t2.phiMass()); - Phi2kaonplus.SetXYZM(t2.phid1Px(), t2.phid1Py(), t2.phid1Pz(), 0.493); - Phi2kaonminus.SetXYZM(t2.phid2Px(), t2.phid2Py(), t2.phid2Pz(), 0.493); - - // Mass windows - if (t1.phiMass() < minPhiMass1 || t1.phiMass() > maxPhiMass1) - continue; + TLorentzVector phi2, k2p, k2m; + phi2.SetXYZM(t2.phiPx(), t2.phiPy(), t2.phiPz(), t2.phiMass()); + k2p.SetXYZM(t2.phid1Px(), t2.phid1Py(), t2.phid1Pz(), 0.493); + k2m.SetXYZM(t2.phid2Px(), t2.phid2Py(), t2.phid2Pz(), 0.493); if (t2.phiMass() < minPhiMass2 || t2.phiMass() > maxPhiMass2) continue; - exotic = Phid1 + Phid2; - if (exotic.M() < minExoticMass || exotic.M() > maxExoticMass) + // Δm cut (configurable) + const double dM = deltaMPhi(phi1.M(), phi2.M()); + if (dM > maxDeltaMPhi) continue; - // Store candidate and bookkeeping - exoticres.emplace_back(exotic.Pt(), exotic.Eta(), exotic.Phi(), exotic.M()); - phi1v.emplace_back(Phid1.Pt(), Phid1.Eta(), Phid1.Phi(), Phid1.M()); - phi2v.emplace_back(Phid2.Pt(), Phid2.Eta(), Phid2.Phi(), Phid2.M()); + TLorentzVector pair = phi1 + phi2; + if (pair.M() < minExoticMass || pair.M() > maxExoticMass) + continue; - kplus1.emplace_back(Phi1kaonplus.Pt(), Phi1kaonplus.Eta(), Phi1kaonplus.Phi(), 0.493); - kminus1.emplace_back(Phi1kaonminus.Pt(), Phi1kaonminus.Eta(), Phi1kaonminus.Phi(), 0.493); - kplus2.emplace_back(Phi2kaonplus.Pt(), Phi2kaonplus.Eta(), Phi2kaonplus.Phi(), 0.493); - kminus2.emplace_back(Phi2kaonminus.Pt(), Phi2kaonminus.Eta(), Phi2kaonminus.Phi(), 0.493); + // daughter anti-merging + ROOT::Math::PtEtaPhiMVector k1pV(k1p.Pt(), k1p.Eta(), k1p.Phi(), 0.493); + ROOT::Math::PtEtaPhiMVector k1mV(k1m.Pt(), k1m.Eta(), k1m.Phi(), 0.493); + ROOT::Math::PtEtaPhiMVector k2pV(k2p.Pt(), k2p.Eta(), k2p.Phi(), 0.493); + ROOT::Math::PtEtaPhiMVector k2mV(k2m.Pt(), k2m.Eta(), k2m.Phi(), 0.493); + if (!passDaughterDR(k1pV, k2pV, k1mV, k2mV, daughterDeltaR)) + continue; - d1id.push_back(t1.phid1Index()); - d2id.push_back(t2.phid1Index()); - d3id.push_back(t1.phid2Index()); - d4id.push_back(t2.phid2Index()); + // store for one-pass fill + pairV.emplace_back(pair.Pt(), pair.Eta(), pair.Phi(), pair.M()); + phi1V.emplace_back(phi1.Pt(), phi1.Eta(), phi1.Phi(), phi1.M()); + phi2V.emplace_back(phi2.Pt(), phi2.Eta(), phi2.Phi(), phi2.M()); + kplus1V.emplace_back(k1p.Pt(), k1p.Eta(), k1p.Phi(), 0.493); + kminus1V.emplace_back(k1m.Pt(), k1m.Eta(), k1m.Phi(), 0.493); + kplus2V.emplace_back(k2p.Pt(), k2p.Eta(), k2p.Phi(), 0.493); + kminus2V.emplace_back(k2m.Pt(), k2m.Eta(), k2m.Phi(), 0.493); } } - if (exoticres.empty()) + if (pairV.empty()) return; - // === Special handling if exactly two candidates were built === - if (exoticres.size() == 2) { - const int i0 = 0, i1 = 1; - - const bool sameFour = - ((d1id[i0] == d1id[i1] || d1id[i0] == d2id[i1]) && - (d2id[i0] == d1id[i1] || d2id[i0] == d2id[i1]) && - (d3id[i0] == d3id[i1] || d3id[i0] == d4id[i1]) && - (d4id[i0] == d3id[i1] || d4id[i0] == d4id[i1])); - - const double deltam0 = deltam(phi1v[i0].M(), phi2v[i0].M()); - const double deltam1 = deltam(phi1v[i1].M(), phi2v[i1].M()); - const bool dr0 = passDaughterDR(kplus1[i0], kplus2[i0], kminus1[i0], kminus2[i0]); - const bool dr1 = passDaughterDR(kplus1[i1], kplus2[i1], kminus1[i1], kminus2[i1]); - - if (sameFour) { - // Choose the candidate closer to mφ (if it passes ΔR) - int keep = (deltam1 < deltam0 && dr1) ? i1 : (dr0 ? i0 : -1); - if (keep >= 0) { - const double dR = deltaR(phi1v[keep].Phi(), phi1v[keep].Eta(), phi2v[keep].Phi(), phi2v[keep].Eta()); - const double dm = (keep == i0) ? deltam0 : deltam1; - histos.fill(HIST("SEMassUnlike"), exoticres[keep].M(), exoticres[keep].Pt(), dR, dm, phimult); - } - } else { - // Independent candidates → fill both (respect ΔR) - if (dr0) { - const double dR0 = deltaR(phi1v[i0].Phi(), phi1v[i0].Eta(), phi2v[i0].Phi(), phi2v[i0].Eta()); - histos.fill(HIST("SEMassUnlike"), exoticres[i0].M(), exoticres[i0].Pt(), dR0, deltam0, phimult); - } - if (dr1) { - const double dR1 = deltaR(phi1v[i1].Phi(), phi1v[i1].Eta(), phi2v[i1].Phi(), phi2v[i1].Eta()); - histos.fill(HIST("SEMassUnlike"), exoticres[i1].M(), exoticres[i1].Pt(), dR1, deltam1, phimult); - } - } - return; - } - - // === General case: fill each candidate once (respect ΔR) === - for (size_t i = 0; i < exoticres.size(); ++i) { - if (!passDaughterDR(kplus1[i], kplus2[i], kminus1[i], kminus2[i])) - continue; - const double dR = deltaR(phi1v[i].Phi(), phi1v[i].Eta(), phi2v[i].Phi(), phi2v[i].Eta()); - const double dm = deltam(phi1v[i].M(), phi2v[i].M()); - histos.fill(HIST("SEMassUnlike"), exoticres[i].M(), exoticres[i].Pt(), dR, dm, phimult); + // --- fill the single THnSparse --- + for (size_t i = 0; i < pairV.size(); ++i) { + TLorentzVector p1, p2, pair; + p1.SetPtEtaPhiM(phi1V[i].Pt(), phi1V[i].Eta(), phi1V[i].Phi(), phi1V[i].M()); + p2.SetPtEtaPhiM(phi2V[i].Pt(), phi2V[i].Eta(), phi2V[i].Phi(), phi2V[i].M()); + pair.SetPtEtaPhiM(pairV[i].Pt(), pairV[i].Eta(), pairV[i].Phi(), pairV[i].M()); + const double dM = deltaMPhi(p1.M(), p2.M()); + const double M = pair.M(); + // const double pT = pair.Pt(); + const double dR = deltaR(p1.Phi(), p1.Eta(), p2.Phi(), p2.Eta()); + histos.fill(HIST("SEMassUnlike"), M, std::abs(p1.Pt() - p2.Pt()) / pair.Pt(), pair.Pt(), dR, dM, pairV.size()); } } - PROCESS_SWITCH(doublephimeson, processopti2, "Process Optimized same event", false); + PROCESS_SWITCH(doublephimeson, processopti3, "Process Optimized same event", false); SliceCache cache; using BinningTypeVertexContributor = ColumnBinningPolicy;