From 4ec9cc8a1882890b5c2efe7a4bc2ce4f5674410a Mon Sep 17 00:00:00 2001 From: Alberto Caliva Date: Wed, 3 Dec 2025 13:11:21 +0100 Subject: [PATCH 1/3] [PWGLF] Add histograms for efficiency factorization --- PWGLF/Tasks/Strangeness/strangenessInJets.cxx | 265 ++++++++++++------ 1 file changed, 176 insertions(+), 89 deletions(-) diff --git a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx index c9e34ae9092..774228cf971 100644 --- a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx +++ b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx @@ -302,6 +302,19 @@ struct StrangenessInJets { registryMC.add("Proton_generated_in_jet", "Proton_generated_in_jet", HistType::kTH2F, {multBinning, ptAxisLongLived}); registryMC.add("Proton_generated_in_ue", "Proton_generated_in_ue", HistType::kTH2F, {multBinning, ptAxisLongLived}); } + + // Histograms to calculate probability of hyperons to be found within jets + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + registryMC.add("K0s_generated_fullevent", "K0s_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("Lambda_generated_fullevent", "Lambda_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("AntiLambda_generated_fullevent", "AntiLambda_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + } + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.add("XiPos_generated_fullevent", "XiPos_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("XiNeg_generated_fullevent", "XiNeg_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("OmegaPos_generated_fullevent", "OmegaPos_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("OmegaNeg_generated_fullevent", "OmegaNeg_generated_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + } } // Histograms for mc reconstructed @@ -350,58 +363,21 @@ struct StrangenessInJets { registryMC.add("Proton_reconstructed_in_jet", "Proton_reconstructed_in_jet", HistType::kTH2F, {multBinning, ptAxisLongLived}); registryMC.add("Proton_reconstructed_in_ue", "Proton_reconstructed_in_ue", HistType::kTH2F, {multBinning, ptAxisLongLived}); } - } - } - - /* - // Calculation of perpendicular axes - void getPerpendicularAxis(TVector3 p, TVector3& u, double sign) - { - // initialization - double ux(0), uy(0), uz(0); - - // components of vector p - const double px = p.X(); - const double py = p.Y(); - const double pz = p.Z(); - - // protection 1 - if (px == 0 && py != 0) { - uy = -(pz * pz) / py; - ux = sign * std::sqrt(py * py - (pz * pz * pz * pz) / (py * py)); - uz = pz; - u.SetXYZ(ux, uy, uz); - return; - } - // protection 2 - if (py == 0 && px != 0) { - ux = -(pz * pz) / px; - uy = sign * std::sqrt(px * px - (pz * pz * pz * pz) / (px * px)); - uz = pz; - u.SetXYZ(ux, uy, uz); - return; - } - - // equation parameters - const double a = px * px + py * py; - const double b = 2.0 * px * pz * pz; - const double c = pz * pz * pz * pz - py * py * py * py - px * px * py * py; - const double delta = b * b - 4.0 * a * c; - - // protection agains delta<0 - if (delta < 0) { - return; + // Histograms to calculate probability of hyperons to be found within jets + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + registryMC.add("K0s_reconstructed_fullevent", "K0s_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("Lambda_reconstructed_fullevent", "Lambda_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("AntiLambda_reconstructed_fullevent", "AntiLambda_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + } + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.add("XiPos_reconstructed_fullevent", "XiPos_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("XiNeg_reconstructed_fullevent", "XiNeg_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("OmegaPos_reconstructed_fullevent", "OmegaPos_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + registryMC.add("OmegaNeg_reconstructed_fullevent", "OmegaNeg_reconstructed_fullevent", HistType::kTH2F, {multBinning, ptAxis}); + } } - - // solutions - ux = (-b + sign * std::sqrt(delta)) / (2.0 * a); - uy = (-pz * pz - px * ux) / py; - uz = pz; - u.SetXYZ(ux, uy, uz); - return; } - */ // Delta phi calculation double getDeltaPhi(double a1, double a2) @@ -1246,14 +1222,58 @@ struct StrangenessInJets { strHadronMomentum.emplace_back(particle.px(), particle.py(), particle.pz()); } - // Select physical primary particles or HF decay products - if (!isPhysicalPrimaryOrFromHF(particle, mcParticles)) - continue; - + // Eta and pt selection double minPtParticle = 0.1; if (particle.eta() < etaMin || particle.eta() > etaMax || particle.pt() < minPtParticle) continue; + // Fill generated histograms for full event + if (particle.isPhysicalPrimary()) { + switch (particle.pdgCode()) { + case kK0Short: + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + registryMC.fill(HIST("K0s_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kLambda0: + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + registryMC.fill(HIST("Lambda_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kLambda0Bar: + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + registryMC.fill(HIST("AntiLambda_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kXiMinus: + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.fill(HIST("XiNeg_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kXiPlusBar: + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.fill(HIST("XiPos_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kOmegaMinus: + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.fill(HIST("OmegaNeg_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + case kOmegaPlusBar: + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + registryMC.fill(HIST("OmegaPos_generated_fullevent"), genMultiplicity, particle.pt()); + } + break; + default: + break; + } + } + + // Select physical primary particles or HF decay products + if (!isPhysicalPrimaryOrFromHF(particle, mcParticles)) + continue; + // Build 4-momentum assuming charged pion mass static constexpr float MassPionChargedSquared = o2::constants::physics::MassPionCharged * o2::constants::physics::MassPionCharged; const double energy = std::sqrt(particle.p() * particle.p() + MassPionChargedSquared); @@ -1472,7 +1492,7 @@ struct StrangenessInJets { // Reconstructed MC events void processMCreconstructed(SimCollisions const& collisions, soa::Join const&, DaughterTracksMC const& mcTracks, aod::V0Datas const& fullV0s, - aod::CascDataExt const& Cascades, const aod::McParticles&) + aod::CascDataExt const& Cascades, aod::McParticles const& mcParticles) { // Define per-event containers std::vector fjParticles; @@ -1520,6 +1540,89 @@ struct StrangenessInJets { auto cascPerColl = Cascades.sliceBy(perCollisionCasc, collision.globalIndex()); auto tracksPerColl = mcTracks.sliceBy(perCollisionTrk, collision.globalIndex()); + // V0 particles + if (enabledSignals.value[ParticleOfInterest::kV0Particles]) { + for (const auto& v0 : v0sPerColl) { + const auto& pos = v0.posTrack_as(); + const auto& neg = v0.negTrack_as(); + + // Get MC particles + if (!pos.has_mcParticle() || !neg.has_mcParticle()) + continue; + auto posParticle = pos.mcParticle_as(); + auto negParticle = neg.mcParticle_as(); + if (!posParticle.has_mothers() || !negParticle.has_mothers()) + continue; + + // Get Mothers + auto motherPos = mcParticles.iteratorAt(posParticle.mothersIds()[0]); + auto motherNeg = mcParticles.iteratorAt(negParticle.mothersIds()[0]); + if (motherPos != motherNeg) + continue; + if (!motherPos.isPhysicalPrimary()) + continue; + + // K0s + if (passedK0ShortSelection(v0, pos, neg) && motherPos.pdgCode() == kK0Short) { + registryMC.fill(HIST("K0s_reconstructed_fullevent"), multiplicity, v0.pt()); + } + // Lambda + if (passedLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0) { + registryMC.fill(HIST("Lambda_reconstructed_fullevent"), multiplicity, v0.pt()); + } + // AntiLambda + if (passedAntiLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0Bar) { + registryMC.fill(HIST("AntiLambda_reconstructed_fullevent"), multiplicity, v0.pt()); + } + } + } + + // Cascades + if (enabledSignals.value[ParticleOfInterest::kCascades]) { + for (const auto& casc : cascPerColl) { + auto bach = casc.bachelor_as(); + auto pos = casc.posTrack_as(); + auto neg = casc.negTrack_as(); + + // Get MC particles + if (!bach.has_mcParticle() || !pos.has_mcParticle() || !neg.has_mcParticle()) + continue; + auto posParticle = pos.mcParticle_as(); + auto negParticle = neg.mcParticle_as(); + auto bachParticle = bach.mcParticle_as(); + if (!posParticle.has_mothers() || !negParticle.has_mothers() || !bachParticle.has_mothers()) + continue; + + // Select particles originating from the same parent + auto motherPos = mcParticles.iteratorAt(posParticle.mothersIds()[0]); + auto motherNeg = mcParticles.iteratorAt(negParticle.mothersIds()[0]); + auto motherBach = mcParticles.iteratorAt(bachParticle.mothersIds()[0]); + if (motherPos != motherNeg) + continue; + if (std::abs(motherPos.pdgCode()) != kLambda0) + continue; + if (!motherBach.isPhysicalPrimary()) + continue; + + // Xi+ + if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && motherBach.pdgCode() == kXiPlusBar) { + registryMC.fill(HIST("XiPos_reconstructed_fullevent"), multiplicity, casc.pt()); + } + // Xi- + if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && motherBach.pdgCode() == kXiMinus) { + registryMC.fill(HIST("XiNeg_reconstructed_fullevent"), multiplicity, casc.pt()); + } + // Omega+ + if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && motherBach.pdgCode() == kOmegaPlusBar) { + registryMC.fill(HIST("OmegaPos_reconstructed_fullevent"), multiplicity, casc.pt()); + } + // Omega- + if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && motherBach.pdgCode() == kOmegaMinus) { + registryMC.fill(HIST("OmegaNeg_reconstructed_fullevent"), multiplicity, casc.pt()); + } + } + } + // Loop over reconstructed tracks for (auto const& track : tracksPerColl) { if (!passedTrackSelectionForJetReconstruction(track)) @@ -1596,19 +1699,12 @@ struct StrangenessInJets { continue; // Select particles originating from the same parent - int pdgParent(0); - bool isPhysPrim = false; - for (const auto& particleMotherOfNeg : negParticle.mothers_as()) { - for (const auto& particleMotherOfPos : posParticle.mothers_as()) { - if (particleMotherOfNeg == particleMotherOfPos) { - pdgParent = particleMotherOfNeg.pdgCode(); - isPhysPrim = particleMotherOfNeg.isPhysicalPrimary(); - } - } - } - if (pdgParent == 0) + auto motherPos = mcParticles.iteratorAt(posParticle.mothersIds()[0]); + auto motherNeg = mcParticles.iteratorAt(negParticle.mothersIds()[0]); + if (motherPos != motherNeg) continue; - + bool isPhysPrim = motherPos.isPhysicalPrimary()); + // Compute distance from jet and UE axes double deltaEtaJet = v0dir.Eta() - selectedJet[i].Eta(); double deltaPhiJet = getDeltaPhi(v0dir.Phi(), selectedJet[i].Phi()); @@ -1621,7 +1717,7 @@ struct StrangenessInJets { double deltaRue2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2); // K0s - if (passedK0ShortSelection(v0, pos, neg) && pdgParent == kK0Short && isPhysPrim) { + if (passedK0ShortSelection(v0, pos, neg) && motherPos.pdgCode() == kK0Short && isPhysPrim) { if (deltaRjet < rJet) { registryMC.fill(HIST("K0s_reconstructed_jet"), multiplicity, v0.pt()); } @@ -1630,7 +1726,7 @@ struct StrangenessInJets { } } // Lambda - if (passedLambdaSelection(v0, pos, neg) && pdgParent == kLambda0 && isPhysPrim) { + if (passedLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0 && isPhysPrim) { if (deltaRjet < rJet) { registryMC.fill(HIST("Lambda_reconstructed_jet"), multiplicity, v0.pt()); } @@ -1639,7 +1735,7 @@ struct StrangenessInJets { } } // AntiLambda - if (passedAntiLambdaSelection(v0, pos, neg) && pdgParent == kLambda0Bar && isPhysPrim) { + if (passedAntiLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0Bar && isPhysPrim) { if (deltaRjet < rJet) { registryMC.fill(HIST("AntiLambda_reconstructed_jet"), multiplicity, v0.pt()); } @@ -1696,23 +1792,14 @@ struct StrangenessInJets { continue; // Select particles originating from the same parent - int pdgParent(0); - bool isPhysPrim = false; - for (const auto& particleMotherOfNeg : negParticle.mothers_as()) { - for (const auto& particleMotherOfPos : posParticle.mothers_as()) { - for (const auto& particleMotherOfBach : bachParticle.mothers_as()) { - if (particleMotherOfNeg != particleMotherOfPos) - continue; - if (std::abs(particleMotherOfNeg.pdgCode()) != kLambda0) - continue; - isPhysPrim = particleMotherOfBach.isPhysicalPrimary(); - pdgParent = particleMotherOfBach.pdgCode(); - } - } - } - if (pdgParent == 0) + auto motherPos = mcParticles.iteratorAt(posParticle.mothersIds()[0]); + auto motherNeg = mcParticles.iteratorAt(negParticle.mothersIds()[0]); + auto motherBach = mcParticles.iteratorAt(bachParticle.mothersIds()[0]); + if (motherPos != motherNeg) + continue; + if (std::abs(motherPos.pdgCode()) != kLambda0) continue; - if (!isPhysPrim) + if (!motherBach.isPhysicalPrimary()) continue; // Compute distances from jet and UE axes @@ -1728,7 +1815,7 @@ struct StrangenessInJets { double deltaRue2 = std::sqrt(deltaEtaUe2 * deltaEtaUe2 + deltaPhiUe2 * deltaPhiUe2); // Xi+ - if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && pdgParent == kXiPlusBar) { + if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && motherBach.pdgCode() == kXiPlusBar) { if (deltaRjet < rJet) { registryMC.fill(HIST("XiPos_reconstructed_jet"), multiplicity, casc.pt()); } @@ -1737,7 +1824,7 @@ struct StrangenessInJets { } } // Xi- - if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && pdgParent == kXiMinus) { + if (passedXiSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && motherBach.pdgCode() == kXiMinus) { if (deltaRjet < rJet) { registryMC.fill(HIST("XiNeg_reconstructed_jet"), multiplicity, casc.pt()); } @@ -1746,7 +1833,7 @@ struct StrangenessInJets { } } // Omega+ - if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && pdgParent == kOmegaPlusBar) { + if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() > 0 && motherBach.pdgCode() == kOmegaPlusBar) { if (deltaRjet < rJet) { registryMC.fill(HIST("OmegaPos_reconstructed_jet"), multiplicity, casc.pt()); } @@ -1755,7 +1842,7 @@ struct StrangenessInJets { } } // Omega- - if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && pdgParent == kOmegaMinus) { + if (passedOmegaSelection(casc, pos, neg, bach, collision) && bach.sign() < 0 && motherBach.pdgCode() == kOmegaMinus) { if (deltaRjet < rJet) { registryMC.fill(HIST("OmegaNeg_reconstructed_jet"), multiplicity, casc.pt()); } From 2a4f6cfb12060c2b0fa46f79994b39ff79635543 Mon Sep 17 00:00:00 2001 From: Alberto Caliva Date: Wed, 3 Dec 2025 13:13:40 +0100 Subject: [PATCH 2/3] removed trailing spaces --- PWGLF/Tasks/Strangeness/strangenessInJets.cxx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx index 774228cf971..f984e580cfb 100644 --- a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx +++ b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx @@ -1561,7 +1561,7 @@ struct StrangenessInJets { continue; if (!motherPos.isPhysicalPrimary()) continue; - + // K0s if (passedK0ShortSelection(v0, pos, neg) && motherPos.pdgCode() == kK0Short) { registryMC.fill(HIST("K0s_reconstructed_fullevent"), multiplicity, v0.pt()); @@ -1704,7 +1704,7 @@ struct StrangenessInJets { if (motherPos != motherNeg) continue; bool isPhysPrim = motherPos.isPhysicalPrimary()); - + // Compute distance from jet and UE axes double deltaEtaJet = v0dir.Eta() - selectedJet[i].Eta(); double deltaPhiJet = getDeltaPhi(v0dir.Phi(), selectedJet[i].Phi()); From 90feb3c5302c56a5f82f6d98cc1daef9df39275d Mon Sep 17 00:00:00 2001 From: Alberto Caliva Date: Wed, 3 Dec 2025 15:57:46 +0100 Subject: [PATCH 3/3] fixed typos --- PWGLF/Tasks/Strangeness/strangenessInJets.cxx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx index f984e580cfb..9a399f58540 100644 --- a/PWGLF/Tasks/Strangeness/strangenessInJets.cxx +++ b/PWGLF/Tasks/Strangeness/strangenessInJets.cxx @@ -1703,7 +1703,7 @@ struct StrangenessInJets { auto motherNeg = mcParticles.iteratorAt(negParticle.mothersIds()[0]); if (motherPos != motherNeg) continue; - bool isPhysPrim = motherPos.isPhysicalPrimary()); + bool isPhysPrim = motherPos.isPhysicalPrimary(); // Compute distance from jet and UE axes double deltaEtaJet = v0dir.Eta() - selectedJet[i].Eta(); @@ -1746,7 +1746,7 @@ struct StrangenessInJets { // Fill inclusive spectra // K0s - if (passedK0ShortSelection(v0, pos, neg) && pdgParent == kK0Short) { + if (passedK0ShortSelection(v0, pos, neg) && motherPos.pdgCode() == kK0Short) { if (deltaRjet < rJet) { registryMC.fill(HIST("K0s_reconstructed_jet_incl"), multiplicity, v0.pt()); } @@ -1755,7 +1755,7 @@ struct StrangenessInJets { } } // Lambda - if (passedLambdaSelection(v0, pos, neg) && pdgParent == kLambda0) { + if (passedLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0) { if (deltaRjet < rJet) { registryMC.fill(HIST("Lambda_reconstructed_jet_incl"), multiplicity, v0.pt()); } @@ -1764,7 +1764,7 @@ struct StrangenessInJets { } } // AntiLambda - if (passedAntiLambdaSelection(v0, pos, neg) && pdgParent == kLambda0Bar) { + if (passedAntiLambdaSelection(v0, pos, neg) && motherPos.pdgCode() == kLambda0Bar) { if (deltaRjet < rJet) { registryMC.fill(HIST("AntiLambda_reconstructed_jet_incl"), multiplicity, v0.pt()); }