diff --git a/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx b/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx index 57fba9c23e9..6ef52e5c8bc 100644 --- a/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx +++ b/PWGCF/FemtoUniverse/TableProducer/femtoUniverseProducerTask.cxx @@ -355,6 +355,7 @@ struct FemtoUniverseProducerTask { Configurable confPhiKaonNsigmaTOFfrom055to15{"confPhiKaonNsigmaTOFfrom055to15", 3.0, "Reject if Kaons in 0.55-1.5 are have TOF n sigma above this value."}; Configurable confPhiKaonNsigmaTPCfrom15{"confPhiKaonNsigmaTPCfrom15", 3.0, "Reject if Kaons above 1.5 are have TPC n sigma above this value."}; Configurable confPhiKaonNsigmaTOFfrom15{"confPhiKaonNsigmaTOFfrom15", 3.0, "Reject if Kaons above 1.5 are have TOF n sigma above this value."}; + Configurable confDoLikeSign{"confDoLikeSign", false, "Enable like-sign (fake Phi) background calculation"}; } ConfPhiSelection; // PDG codes for fillMCParticle function @@ -2229,6 +2230,156 @@ struct FemtoUniverseProducerTask { } } + template + void fillLikeSignPhi(CollisionType const& col, TrackType const& tracks) + { + std::vector childIDs = {0, 0}; + std::vector tmpIDtrack; + + // Loop over all track combinations using FullIndexPolicy + for (const auto& [p1, p2] : combinations(soa::CombinationsFullIndexPolicy(tracks, tracks))) { + + // 1. Check minimal track quality cuts for both tracks + if (!trackCuts.isSelectedMinimal(p1) || !trackCuts.isSelectedMinimal(p2)) { + continue; + } + + if (p1.sign() != p2.sign()) { // K+K- pairs are already handled in fillPhi, so we skip them here, only K+K+ and K-K- + continue; + } + + // avoid double counting and self-pairing when using FullIndexPolicy + if (p1.globalIndex() >= p2.globalIndex()) { + continue; + } + + // cuts for daughters + if (ConfPhiSelection.confPhiDoLFPID4Kaons) { + if (!p1.isGlobalTrackWoDCA() || !p2.isGlobalTrackWoDCA() || !p1.isPVContributor() || !p2.isPVContributor()) { + continue; + } + if (!(isKaonNSigmaLF(p1.pt(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon), p1.hasTOF()))) { + continue; + } + if (!(isKaonNSigmaLF(p2.pt(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon), p2.hasTOF()))) { + continue; + } + } else { + if (!(isKaonNSigma(p1.pt(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p1, o2::track::PID::Kaon)))) { + continue; + } + if (!(isKaonNSigma(p2.pt(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Kaon), trackCuts.getNsigmaTOF(p2, o2::track::PID::Kaon)))) { + continue; + } + } + + // Reject protons and pions contamination + if (isKaonRejected(p1.p(), trackCuts.getNsigmaTPC(p1, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p1, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p1, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p1, o2::track::PID::Pion))) { + continue; + } + if (isKaonRejected(p2.p(), trackCuts.getNsigmaTPC(p2, o2::track::PID::Proton), trackCuts.getNsigmaTOF(p2, o2::track::PID::Proton), trackCuts.getNsigmaTPC(p2, o2::track::PID::Pion), trackCuts.getNsigmaTOF(p2, o2::track::PID::Pion))) { + continue; + } + + const auto mMassOne = o2::constants::physics::MassKPlus; + const auto mMassTwo = o2::constants::physics::MassKPlus; + + ROOT::Math::LorentzVector> part1Vec(p1.pt(), p1.eta(), p1.phi(), mMassOne); + ROOT::Math::LorentzVector> part2Vec(p2.pt(), p2.eta(), p2.phi(), mMassTwo); + + ROOT::Math::LorentzVector> sumVec(part1Vec); + sumVec += part2Vec; + + // kinematic cuts for phi candidate + float phiEta = sumVec.Eta(); + if (std::abs(phiEta) > ConfPhiSelection.confPhiEtaHighLimit.value) { + continue; + } + + float phiPt = sumVec.Pt(); + if ((phiPt < ConfPhiSelection.confPhiPtLowLimit.value) || (phiPt > ConfPhiSelection.confPhiPtHighLimit.value)) { + continue; + } + + float phiPhi = RecoDecay::constrainAngle(sumVec.Phi(), 0); + float phiM = sumVec.M(); + + if ((phiM < ConfPhiSelection.confPhiInvMassLowLimit) || (phiM > ConfPhiSelection.confPhiInvMassUpLimit)) { + continue; + } + + phiCuts.fillQA(col, p1, p1, p2, 321, 321); + + int postrackID = p1.globalIndex(); + int rowInPrimaryTrackTablePos = -1; + rowInPrimaryTrackTablePos = getRowDaughters(postrackID, tmpIDtrack); + childIDs[0] = rowInPrimaryTrackTablePos; + childIDs[1] = 0; + + // first daughter track + outputParts(outputCollision.lastIndex(), p1.pt(), + p1.eta(), p1.phi(), + aod::femtouniverseparticle::ParticleType::kPhiChild, + -999, + -999, + p1.dcaXY(), + childIDs, + 0, + p1.sign()); + const int rowOfPosTrack = outputParts.lastIndex(); + if constexpr (isMC) { + fillMCParticle(p1, o2::aod::femtouniverseparticle::ParticleType::kPhiChild); + } + + int negtrackID = p2.globalIndex(); + int rowInPrimaryTrackTableNeg = -1; + rowInPrimaryTrackTableNeg = getRowDaughters(negtrackID, tmpIDtrack); + childIDs[0] = 0; + childIDs[1] = rowInPrimaryTrackTableNeg; + + // second daughter track + outputParts(outputCollision.lastIndex(), + p2.pt(), + p2.eta(), + p2.phi(), + aod::femtouniverseparticle::ParticleType::kPhiChild, + -999, + -999, + p2.dcaXY(), + childIDs, + 0, + p2.sign()); + const int rowOfNegTrack = outputParts.lastIndex(); + if constexpr (isMC) { + fillMCParticle(p2, o2::aod::femtouniverseparticle::ParticleType::kPhiChild); + } + + std::vector indexChildID = {rowOfPosTrack, rowOfNegTrack}; + + // Define like-sign flag: +1 for (K+K+), -1 for (K-K-) + int likeSignFlag = p1.sign(); + + // Save the Fake Phi candidate from K+K+ or K-K- combination + outputParts(outputCollision.lastIndex(), + phiPt, + phiEta, + phiPhi, + aod::femtouniverseparticle::ParticleType::kPhi, + likeSignFlag, // Flag: +1 for K+K+, -1 for K-K-, before it was empty slot for cutcontainer + 0, + phiM, // Invariant mass of the fake Phi candidate + indexChildID, + phiM, + 0); + + if (confIsDebug) { + fillDebugParticle(p1); // QA for first daughter + fillDebugParticle(p2); // QA for second daughter + fillDebugParticle(p1); // QA for fake Phi mother + } + } + } + template void fillParticles(TrackType const& tracks, std::optional>> recoMcIds = std::nullopt) { @@ -2452,6 +2603,9 @@ struct FemtoUniverseProducerTask { } if (ConfGeneral.confIsActivatePhi) { fillPhi(col, tracks); + if (ConfPhiSelection.confDoLikeSign) { + fillLikeSignPhi(col, tracks); + } } } return colcheck; @@ -2628,6 +2782,9 @@ struct FemtoUniverseProducerTask { if (colcheck) { fillTracks(tracks); fillPhi(col, tracks); + if (ConfPhiSelection.confDoLikeSign) { + fillLikeSignPhi(col, tracks); + } } } PROCESS_SWITCH(FemtoUniverseProducerTask, processTrackPhiData, diff --git a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackPhi.cxx b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackPhi.cxx index 03c24e77153..8b5948e5472 100644 --- a/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackPhi.cxx +++ b/PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackPhi.cxx @@ -153,6 +153,7 @@ struct FemtoUniversePairTaskTrackPhi { Configurable confPhiPtHigh{"confPhiPtHigh", 4.0, "Higher limit of the Phi pT."}; Configurable confInvMassLowLimitPhi{"confInvMassLowLimitPhi", 1.011, "Lower limit of the Phi invariant mass"}; // change that to do invariant mass cut Configurable confInvMassUpLimitPhi{"confInvMassUpLimitPhi", 1.027, "Upper limit of the Phi invariant mass"}; + Configurable confFillMinvBackground{"confFillMinvBackground", false, "Fill background histograms for Phi invariant mass."}; /// Partitions for the Phi meson (particle 2) Partition partsPhi = (aod::femtouniverseparticle::partType == uint8_t(aod::femtouniverseparticle::ParticleType::kPhi)) && @@ -551,37 +552,39 @@ struct FemtoUniversePairTaskTrackPhi { sameEventCont.setPair(track, phicandidate, multCol, confUse3D, mWeight); } - // Used for better fitting of invariant mass background. + if (confFillMinvBackground) { + // Used for better fitting of invariant mass background. - TLorentzVector part1Vec; - TLorentzVector part2Vec; + TLorentzVector part1Vec; + TLorentzVector part2Vec; - float mMassOne = o2::constants::physics::MassKPlus; - float mMassTwo = o2::constants::physics::MassKMinus; + float mMassOne = o2::constants::physics::MassKPlus; + float mMassTwo = o2::constants::physics::MassKMinus; - for (auto const& [kaon1, kaon2] : combinations(CombinationsStrictlyUpperIndexPolicy(partsKaons, partsKaons))) { - if ((kaon1.mAntiLambda() == 1) && (kaon2.mAntiLambda() == 1)) { - part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassOne); - part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassOne); - TLorentzVector sumVec(part1Vec); - sumVec += part2Vec; - registryPhiMinvBackground.fill(HIST("InvariantMassKpKp"), sumVec.M()); - } - if ((kaon1.mAntiLambda() == -1) && (kaon2.mAntiLambda() == -1)) { - part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassTwo); - part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassTwo); + for (auto const& [kaon1, kaon2] : combinations(CombinationsStrictlyUpperIndexPolicy(partsKaons, partsKaons))) { + if ((kaon1.mAntiLambda() == 1) && (kaon2.mAntiLambda() == 1)) { + part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassOne); + part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassOne); + TLorentzVector sumVec(part1Vec); + sumVec += part2Vec; + registryPhiMinvBackground.fill(HIST("InvariantMassKpKp"), sumVec.M()); + } + if ((kaon1.mAntiLambda() == -1) && (kaon2.mAntiLambda() == -1)) { + part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassTwo); + part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassTwo); - TLorentzVector sumVec(part1Vec); - sumVec += part2Vec; - registryPhiMinvBackground.fill(HIST("InvariantMassKmKm"), sumVec.M()); - } - if ((kaon1.mAntiLambda() == 1) && (kaon2.mAntiLambda() == -1)) { - part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassOne); - part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassTwo); + TLorentzVector sumVec(part1Vec); + sumVec += part2Vec; + registryPhiMinvBackground.fill(HIST("InvariantMassKmKm"), sumVec.M()); + } + if ((kaon1.mAntiLambda() == 1) && (kaon2.mAntiLambda() == -1)) { + part1Vec.SetPtEtaPhiM(kaon1.pt(), kaon1.eta(), kaon1.phi(), mMassOne); + part2Vec.SetPtEtaPhiM(kaon2.pt(), kaon2.eta(), kaon2.phi(), mMassTwo); - TLorentzVector sumVec(part1Vec); - sumVec += part2Vec; - registryPhiMinvBackground.fill(HIST("InvariantMassKpKm"), sumVec.M()); + TLorentzVector sumVec(part1Vec); + sumVec += part2Vec; + registryPhiMinvBackground.fill(HIST("InvariantMassKpKm"), sumVec.M()); + } } } }