From 43405e4c9e8d9ed807c8777142de881d0eb9e5bf Mon Sep 17 00:00:00 2001 From: blacwovie Date: Fri, 22 May 2026 15:17:48 +0800 Subject: [PATCH 1/2] add CPR switch --- .../TableProducer/HadNucleiFemto.cxx | 111 ++++++++++++++++-- 1 file changed, 104 insertions(+), 7 deletions(-) diff --git a/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx b/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx index 88540514c26..409eaf5c75d 100644 --- a/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx +++ b/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx @@ -83,6 +83,7 @@ namespace { constexpr double betheBlochDefault[1][6]{{-136.71, 0.441, 0.2269, 1.347, 0.8035, 0.09}}; static const std::vector betheBlochParNames{"p0", "p1", "p2", "p3", "p4", "resolution"}; +static constexpr std::array tmpRadiiTPC{{85.f, 105.f, 125.f, 145.f, 165.f, 185.f, 205.f, 225.f, 245.f}}; enum Selections { kNoCuts = 0, @@ -189,6 +190,8 @@ struct HadNucleiFemto { Configurable settingCutNsigTOFPrMin{"settingCutNsigTOFPrMin", 3.0f, "Minimum TOF Pr Nsigma cut for rejection"}; Configurable settingCutNsigTOFPiMin{"settingCutNsigTOFPiMin", 3.0f, "Minimum TOF Pi Nsigma cut for rejection"}; Configurable settingCutNsigTOFKaMin{"settingCutNsigTOFKaMin", 3.0f, "Minimum TOF Ka Nsigma cut for rejection"}; + Configurable settingEnablePionProtonRejection{"settingEnablePionProtonRejection", true, "If true, apply proton rejection in the pion PID"}; + Configurable settingEnablePionKaonRejection{"settingEnablePionKaonRejection", true, "If true, apply kaon rejection in the pion PID"}; Configurable settingUsePionReferencePIDCuts{"settingUsePionReferencePIDCuts", false, "If true, use the reference pion track/PID cuts from the pi-p study"}; // Deuteron purity and PID cuts Configurable settingCutPinMinDe{"settingCutPinMinDe", 0.0f, "Minimum Pin for De"}; @@ -199,6 +202,12 @@ struct HadNucleiFemto { Configurable settingCutNsigmaTPCDe{"settingCutNsigmaTPCDe", 2.5f, "Value of the TPC Nsigma cut on De"}; Configurable settingCutNsigmaITSDe{"settingCutNsigmaITSDe", 2.5f, "Value of the ITD Nsigma cut on De"}; Configurable settingCutNsigmaTOFTPCDe{"settingCutNsigmaTOFTPCDe", 2.5f, "Value of the De TOF TPC combNsigma cut"}; + Configurable settingUseProtonMassForKstarMt{"settingUseProtonMassForKstarMt", false, "If true, use proton mass instead of deuteron mass for kstar and mT"}; + Configurable settingEnableClosePairRejection{"settingEnableClosePairRejection", false, "Enable close pair rejection for deuteron-hadron track pairs"}; + Configurable settingClosePairDeltaPhiMax{"settingClosePairDeltaPhiMax", 0.01f, "Maximum delta phi star for close pair rejection"}; + Configurable settingClosePairDeltaEtaMax{"settingClosePairDeltaEtaMax", 0.01f, "Maximum delta eta for close pair rejection"}; + Configurable settingClosePairRadiusMode{"settingClosePairRadiusMode", 1, "Close pair rejection mode: 0 = PV, 1 = average phi star, 2 = specific TPC radius"}; + Configurable settingClosePairSpecificRadius{"settingClosePairSpecificRadius", 85.f, "TPC radius in cm used when close pair rejection mode is 2"}; // Hypertriton-specific cuts Configurable settingCutTPCChi2He{"settingCutTPCChi2He", 0.0f, "Minimum tpcChi2He for Hyper He3"}; Configurable settingCutAverClsSizeHe{"settingCutAverClsSizeHe", 0.0f, "Minimum averClusSizeHe for Hyper He3"}; @@ -278,6 +287,8 @@ struct HadNucleiFemto { {"hHadPin", "P distribution; #it{p} (GeV/#it{c})", {HistType::kTH1F, {{120, -4.0f, 4.0f}}}}, {"hHadEta", "eta distribution; #eta(had)", {HistType::kTH1F, {{200, -1.0f, 1.0f}}}}, {"hHadPhi", "phi distribution; phi(had)", {HistType::kTH1F, {{600, -4.0f, 4.0f}}}}, + {"h2CPRBefore", "Close pair rejection before cut; #Delta#eta; #Delta#phi^{*}", {HistType::kTH2F, {{160, -2.0f, 2.0f}, {160, -3.2f, 3.2f}}}}, + {"h2CPRAfter", "Close pair rejection after cut; #Delta#eta; #Delta#phi^{*}", {HistType::kTH2F, {{160, -2.0f, 2.0f}, {160, -3.2f, 3.2f}}}}, // dE/dx {"h2dEdxNucandidates", "dEdx distribution; #it{p} (GeV/#it{c}); dE/dx (a.u.)", {HistType::kTH2F, {{200, -5.0f, 5.0f}, {100, 0.0f, 2000.0f}}}}, @@ -602,6 +613,85 @@ struct HadNucleiFemto { return true; } + template + float phiAtSpecificRadiiTPC(const Ttrack& track, float radius) const + { + const float absPt = std::abs(track.pt()); + if (absPt <= 0.f) { + return 999.f; + } + const float arg = 0.3f * static_cast(track.sign()) * 0.1f * mDbz * radius * 0.01f / (2.f * absPt); + if (std::fabs(arg) >= 1.f) { + return 999.f; + } + return track.phi() - std::asin(arg); + } + + float wrapDeltaPhi(float dphi) const + { + return std::atan2(std::sin(dphi), std::cos(dphi)); + } + + template + float averagePhiStar(const Ttrack1& track1, const Ttrack2& track2) const + { + constexpr float invalidPhiStar = 999.f; + float dPhiAvg = 0.f; + int meaningfulEntries = 0; + for (const auto& radius : tmpRadiiTPC) { + const float phi1 = phiAtSpecificRadiiTPC(track1, radius); + const float phi2 = phiAtSpecificRadiiTPC(track2, radius); + if (phi1 == invalidPhiStar || phi2 == invalidPhiStar) { + continue; + } + dPhiAvg += wrapDeltaPhi(phi1 - phi2); + meaningfulEntries++; + } + if (meaningfulEntries == 0) { + return invalidPhiStar; + } + return dPhiAvg / static_cast(meaningfulEntries); + } + + template + bool isClosePair(const Ttrack1& track1, const Ttrack2& track2) + { + constexpr int closePairRadiusModePv = 0; + constexpr int closePairRadiusModeSpecificTpc = 2; + constexpr float invalidPhiStar = 999.f; + if (!settingEnableClosePairRejection.value) { + return false; + } + if (track1.sign() != track2.sign()) { + return false; + } + + const float deta = track1.eta() - track2.eta(); + const float dphiAtPV = wrapDeltaPhi(track1.phi() - track2.phi()); + const float dphiAtSpecificRadius = wrapDeltaPhi(phiAtSpecificRadiiTPC(track1, settingClosePairSpecificRadius.value) - phiAtSpecificRadiiTPC(track2, settingClosePairSpecificRadius.value)); + const float dphiAvg = averagePhiStar(track1, track2); + + float dphiToCut = dphiAvg; + if (settingClosePairRadiusMode.value == closePairRadiusModePv) { + dphiToCut = dphiAtPV; + } else if (settingClosePairRadiusMode.value == closePairRadiusModeSpecificTpc) { + dphiToCut = dphiAtSpecificRadius; + } + + if (dphiToCut == invalidPhiStar) { + return false; + } + + mQaRegistry.fill(HIST("h2CPRBefore"), deta, dphiToCut); + const bool isRejected = std::pow(dphiToCut, 2.f) / std::pow(settingClosePairDeltaPhiMax.value, 2.f) + + std::pow(deta, 2.f) / std::pow(settingClosePairDeltaEtaMax.value, 2.f) < + 1.f; + if (!isRejected) { + mQaRegistry.fill(HIST("h2CPRAfter"), deta, dphiToCut); + } + return isRejected; + } + template bool selectionPIDKaon(const Ttrack& candidate) { @@ -700,15 +790,15 @@ struct HadNucleiFemto { if (std::abs(candidate.pt()) < settingCutHadptMin || std::abs(candidate.pt()) > settingCutHadptMax) return false; // reject protons and kaons - if (std::abs(candidate.tpcNSigmaPr()) < settingCutNsigTPCPrMin) + if (settingEnablePionProtonRejection && std::abs(candidate.tpcNSigmaPr()) < settingCutNsigTPCPrMin) return false; - if (std::abs(candidate.tpcNSigmaKa()) < settingCutNsigTPCKaMin) + if (settingEnablePionKaonRejection && std::abs(candidate.tpcNSigmaKa()) < settingCutNsigTPCKaMin) return false; mQaRegistry.fill(HIST("h2NsigmaHadPrTPC"), candidate.tpcNSigmaPr()); mQaRegistry.fill(HIST("h2NsigmaHadKaTPC"), candidate.tpcNSigmaKa()); - if (candidate.hasTOF() && std::abs(candidate.tofNSigmaPr()) < settingCutNsigTOFPrMin) + if (settingEnablePionProtonRejection && candidate.hasTOF() && std::abs(candidate.tofNSigmaPr()) < settingCutNsigTOFPrMin) return false; - if (candidate.hasTOF() && std::abs(candidate.tofNSigmaKa()) < settingCutNsigTOFKaMin) + if (settingEnablePionKaonRejection && candidate.hasTOF() && std::abs(candidate.tofNSigmaKa()) < settingCutNsigTOFKaMin) return false; mQaRegistry.fill(HIST("h2NsigmaHadPrTOF"), candidate.tofNSigmaPr()); if (candidate.hasTOF()) { @@ -794,7 +884,7 @@ struct HadNucleiFemto { const float deDCAzMax = 0.004f + 0.013f / absPt; if (std::abs(candidate.dcaXY()) > deDCAxyMax || std::abs(candidate.dcaZ()) > deDCAzMax) return false; - + if (candidate.hasTOF() && candidate.tpcInnerParam() > settingCutPinMinTOFITSDe) { auto tofNSigmaDe = candidate.tofNSigmaDe(); auto combNsigma = std::sqrt(tofNSigmaDe * tofNSigmaDe + tpcNSigmaDe * tpcNSigmaDe); @@ -967,8 +1057,9 @@ struct HadNucleiFemto { hadNucand.massTOFHad = trackHad.tpcInnerParam() * std::sqrt(1.f / (beta * beta) - 1.f); } - hadNucand.kstar = o2::analysis::femtoWorld::FemtoWorldMath::getkstar(trackHad, MassHad, trackDe, o2::constants::physics::MassDeuteron); - hadNucand.mT = o2::analysis::femtoWorld::FemtoWorldMath::getmT(trackHad, MassHad, trackDe, o2::constants::physics::MassDeuteron); + const float massLightNucleusForKstarMt = settingUseProtonMassForKstarMt ? static_cast(o2::constants::physics::MassProton) : static_cast(o2::constants::physics::MassDeuteron); + hadNucand.kstar = o2::analysis::femtoWorld::FemtoWorldMath::getkstar(trackHad, MassHad, trackDe, massLightNucleusForKstarMt); + hadNucand.mT = o2::analysis::femtoWorld::FemtoWorldMath::getmT(trackHad, MassHad, trackDe, massLightNucleusForKstarMt); return true; } @@ -1075,6 +1166,9 @@ struct HadNucleiFemto { if (!selectTrackHadron(track1) || !selectionPIDHadron(track1)) { continue; } + if (isClosePair(track0, track1)) { + continue; + } SVCand trackPair; trackPair.tr0Idx = track0.globalIndex(); @@ -1130,6 +1224,9 @@ struct HadNucleiFemto { if (!selectTrackHadron(hadCand) || !selectionPIDHadron(hadCand)) { continue; } + if (isClosePair(DeCand, hadCand)) { + continue; + } SVCand trackPair; trackPair.tr0Idx = DeCand.globalIndex(); From 9a0d467bf3a5da9c56362d5e676b03c0709c59d3 Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Fri, 22 May 2026 07:18:43 +0000 Subject: [PATCH 2/2] Please consider the following formatting changes --- PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx b/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx index 409eaf5c75d..c132ac4c3bf 100644 --- a/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx +++ b/PWGCF/Femto/FemtoNuclei/TableProducer/HadNucleiFemto.cxx @@ -884,7 +884,7 @@ struct HadNucleiFemto { const float deDCAzMax = 0.004f + 0.013f / absPt; if (std::abs(candidate.dcaXY()) > deDCAxyMax || std::abs(candidate.dcaZ()) > deDCAzMax) return false; - + if (candidate.hasTOF() && candidate.tpcInnerParam() > settingCutPinMinTOFITSDe) { auto tofNSigmaDe = candidate.tofNSigmaDe(); auto combNsigma = std::sqrt(tofNSigmaDe * tofNSigmaDe + tpcNSigmaDe * tpcNSigmaDe);