diff --git a/EventFiltering/PWGHF/HFFilter.cxx b/EventFiltering/PWGHF/HFFilter.cxx index 6d08f725dcf..2a5e69af2b1 100644 --- a/EventFiltering/PWGHF/HFFilter.cxx +++ b/EventFiltering/PWGHF/HFFilter.cxx @@ -79,6 +79,13 @@ struct HfFilter { // Main struct for HF triggers Configurable requireCharmMassForFemto{"requireCharmMassForFemto", false, "Flags to enable/disable cut on charm-hadron invariant-mass window for femto"}; Configurable ptThresholdForFemtoPid{"ptThresholdForFemtoPid", 8., "pT threshold for changing strategy of proton PID in femto triggers"}; + // parameters for triggers with SigmaC0,++ + Configurable deltaMassLcForSigmaC{"deltaMassLcForSigmaC", 0.04, "Absolute delta cut from PDG Lc mass for Lc candidates to build SigmaC0, ++ candidates"}; + Configurable deltaMassMinSigmaC{"deltaMassMinSigmaC", 0.155, "Minimum delta mass M(pKpipi)-M(pKpi) of SigmaC0,++ candidates"}; + Configurable deltaMassMaxSigmaC{"deltaMassMaxSigmaC", 0.18, "Maximum delta mass M(pKpipi)-M(pKpi) of SigmaC0,++ candidates"}; + Configurable massMinSigmaCKaonPair{"massMinSigmaCKaonPair", 2.92, "Minimum invariant mass for M(SigmaC K) candidates"}; + Configurable massMaxSigmaCKaonPair{"massMaxSigmaCKaonPair", 3.25, "Maximum invariant mass for M(SigmaC K) candidates"}; + // double charm Configurable> enableDoubleCharmChannels{"enableDoubleCharmChannels", {activeDoubleCharmChannels[0], 1, 3, labelsEmpty, labelsColumnsDoubleCharmChannels}, "Flags to enable/disable double charm channels"}; @@ -143,7 +150,7 @@ struct HfFilter { // Main struct for HF triggers std::array, kNCharmParticles> hCharmHighPt{}; std::array, kNCharmParticles> hCharmProtonKstarDistr{}; std::array, kNBeautyParticles> hMassVsPtB{}; - std::array, kNCharmParticles + 9> hMassVsPtC{}; // +6 for resonances (D*+, D*0, Ds*+, Ds1+, Ds2*+, Xic* right sign, Xic* wrong sign) +2 for charm baryons (Xi+Pi, Xi+Ka) + std::array, kNCharmParticles + 13> hMassVsPtC{}; // +6 for resonances (D*+, D*0, Ds*+, Ds1+, Ds2*+, Xic* right sign, Xic* wrong sign) +2 for charm baryons (Xi+Pi, Xi+Ka) +2 for SigmaC (SigmaC++, SigmaC0) +2 for SigmaCK pairs (SigmaC++K-, SigmaC0K0s) std::shared_ptr hProtonTPCPID, hProtonTOFPID; std::array, kNCharmParticles> hBDTScoreBkg{}; std::array, kNCharmParticles> hBDTScorePrompt{}; @@ -191,6 +198,7 @@ struct HfFilter { // Main struct for HF triggers helper.setXiSelections(cutsXiCascades->get(0u, 0u), cutsXiCascades->get(0u, 1u), cutsXiCascades->get(0u, 2u), cutsXiCascades->get(0u, 3u), cutsXiCascades->get(0u, 4u), cutsXiCascades->get(0u, 5u), cutsXiCascades->get(0u, 6u), cutsXiCascades->get(0u, 7u)); helper.setNsigmaPiCutsForCharmBaryonBachelor(nSigmaPidCuts->get(0u, 4u), nSigmaPidCuts->get(1u, 4u)); helper.setTpcPidCalibrationOption(setTPCCalib); + helper.setDeltaMassRangeSigmaC(deltaMassMinSigmaC, deltaMassMaxSigmaC); hProcessedEvents = registry.add("fProcessedEvents", "HF - event filtered;;counts", HistType::kTH1F, {{kNtriggersHF + 2, -0.5, kNtriggersHF + 1.5}}); for (auto iBin = 0; iBin < kNtriggersHF + 2; ++iBin) { @@ -221,6 +229,12 @@ struct HfFilter { // Main struct for HF triggers // charm baryons to LF cascades hMassVsPtC[kNCharmParticles + 7] = registry.add("fMassVsPtCharmBaryonToXiPi", "#it{M} vs. #it{p}_{T} distribution of triggered #Xi+#pi candidates;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 7]}); hMassVsPtC[kNCharmParticles + 8] = registry.add("fMassVsPtCharmBaryonToXiKa", "#it{M} vs. #it{p}_{T} distribution of triggered #Xi+K candidates;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 8]}); + // SigmaC0,++ + hMassVsPtC[kNCharmParticles + 9] = registry.add("fDeltaMassVsPtSigmaCPlusPlus", "#it{M}(pK#pi#pi)-M(pK#pi) vs. #it{p}_{T} distribution of #Sigma_{c}^{++} candidates for triggers;#it{p}_{T}(#Sigma_{c}^{++}) (GeV/#it{c});#it{M}(pK#pi#pi)-M(pK#pi);counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 9]}); + hMassVsPtC[kNCharmParticles + 10] = registry.add("fDeltaMassVsPtSigmaC0", "#it{M}(pK#pi#pi)-M(pK#pi) vs. #it{p}_{T} distribution of #Sigma_{c}^{0} candidates for triggers;#it{p}_{T}(#Sigma_{c}^{0}) (GeV/#it{c});#it{M}(pK#pi#pi)-M(pK#pi);counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 10]}); + // SigmaCKaon pairs + hMassVsPtC[kNCharmParticles + 11] = registry.add("fMassVsPtSigmaCPlusPlusKaMinus", "#it{M}(#Sigma_{c}^{++}K^{-}) vs. #it{p}_{T} distribution of of triggered #Sigma_{c}^{++}K^{-} pairs;#it{p}_{T} (GeV/#it{c});#it{M}(#Sigma_{c}^{++}K^{-});counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 11]}); + hMassVsPtC[kNCharmParticles + 12] = registry.add("fMassVsPtSigmaC0Ka0s", "#it{M}(#Sigma_{c}^{0}K^{0}_{s}) vs. #it{p}_{T} distribution of of triggered #Sigma_{c}^{0}K^{0}_{s} pairs;#it{p}_{T} (GeV/#it{c});#it{M}(#Sigma_{c}^{++}K^{-});counts", HistType::kTH2F, {ptAxis, massAxisC[kNCharmParticles + 12]}); for (int iBeautyPart{0}; iBeautyPart < kNBeautyParticles; ++iBeautyPart) { hMassVsPtB[iBeautyPart] = registry.add(Form("fMassVsPt%s", beautyParticleNames[iBeautyPart].data()), Form("#it{M} vs. #it{p}_{T} distribution of triggered %s candidates;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", beautyParticleNames[iBeautyPart].data()), HistType::kTH2F, {ptAxis, massAxisB[iBeautyPart]}); } @@ -871,56 +885,78 @@ struct HfFilter { // Main struct for HF triggers // SigmaC++ K- trigger if (is3Prong[2] > 0) { // we need a candidate Lc->pKpi - std::array massDausLcPKPi{massProton, massKa, massPi}; - std::array massDausLcPiKP{massPi, massKa, massProton}; - float massLcPKPi = RecoDecay::m(pVec3Prong, massDausLcPKPi); - float massLcPiKP = RecoDecay::m(pVec3Prong, massDausLcPiKP); - if (!helper.isSelectedLcInMassRangeSigmaC(massLcPKPi, massLcPiKP)) { - /// Lc outside the allowed mass range + // select candidate Lc->pKpi in mass range + if (!is3ProngInMass[2]) { continue; } - // we first look to SigmaC++ + // look for SigmaC++ candidates for (const auto& trackSoftPiId : trackIdsThisCollision) { // start loop over tracks (soft pi) // soft pion candidates auto trackSoftPi = trackSoftPiId.track_as(); - // check the candidate SigmaC++ charge - std::array charges = {trackFirst.sign(), trackSecond.sign(), trackThird.sign(), trackSoftPi.sign()}; - if (std::abs(std::accumulate(charges.begin(), charges.end(), 0)) != 2) { - continue; - } - // exclude tracks already used to build the 3-prong candidate auto globalIndexSoftPi = trackSoftPi.globalIndex(); - if (cand3Prong.prong0Id() == trackFirst.globalIndex() || cand3Prong.prong0Id() == trackSecond.globalIndex() || cand3Prong.prong0Id() == trackThird.globalIndex()) { + if (globalIndexSoftPi == trackFirst.globalIndex() || globalIndexSoftPi == trackSecond.globalIndex() || globalIndexSoftPi == trackThird.globalIndex()) { // do not consider as candidate soft pion a track already used to build the current 3-prong candidate continue; } + // check the candidate SigmaC++ charge + std::array chargesSc = {trackFirst.sign(), trackSecond.sign(), trackThird.sign(), trackSoftPi.sign()}; + int chargeSc = std::accumulate(chargesSc.begin(), chargesSc.end(), 0); // SIGNED electric charge of SigmaC candidate + if (std::abs(chargeSc) != 2) { + continue; + } + // select soft pion candidates auto trackParSoftPi = getTrackPar(trackSoftPi); o2::gpu::gpustd::array dcaSoftPi{trackSoftPi.dcaXY(), trackSoftPi.dcaZ()}; std::array pVecSoftPi = {trackSoftPi.px(), trackSoftPi.py(), trackSoftPi.pz()}; if (trackSoftPi.collisionId() != thisCollId) { + // This is a track reassociated to this PV by the track-to-collision-associator + // Let's propagate this track to it, and calculate dcaXY, dcaZ o2::base::Propagator::Instance()->propagateToDCABxByBz({collision.posX(), collision.posY(), collision.posZ()}, trackParSoftPi, 2.f, noMatCorr, &dcaSoftPi); getPxPyPz(trackParSoftPi, pVecSoftPi); } - int isSoftPionSelected = helper.isSelectedTrackForSoftPionOrBeauty(trackSoftPi, trackParSoftPi, dcaSoftPi, kSigmaCppKminus); + int8_t isSoftPionSelected = helper.isSelectedTrackForSoftPionOrBeauty(trackSoftPi, trackParSoftPi, dcaSoftPi, kSigmaCppKminus); if (TESTBIT(isSoftPionSelected, kSoftPionForSigmaC) /*&& (TESTBIT(is3Prong[2], 0) || TESTBIT(is3Prong[2], 1))*/) { // check the mass of the SigmaC++ candidate - std::array pVecDausSigmaC{pVecFirst, pVecSecond, pVecThird, pVecSoftPi}; - std::array massDausSigmaCToLcPKPi{massProton, massKa, massPi, massPi}; - std::array massDausSigmaCToLcPiKP{massPi, massKa, massProton, massPi}; - float massSigmaCToLcPKPi = RecoDecay::m(pVecDausSigmaC, massDausSigmaCToLcPKPi) : float massSigmaCToLcPiKP = RecoDecay::m(pVecDausSigmaC, massDausSigmaCToLcPiKP) : if (!helper.isSelectedSigmaCInDeltaMassRange(massSigmaCToLcPKPi - massLcPKPi, massSigmaCToLcPiKP - massLcPiKP)) + auto pVecSigmaC = RecoDecay::pVec(pVecFirst, pVecSecond, pVecThird, pVecSoftPi); + auto ptSigmaC= RecoDecay::pt(pVecSigmaC); + if (!helper.isSelectedSigmaCInDeltaMassRange(pVecFirst, pVecThird, pVecSecond, pVecSoftPi, ptSigmaC, is3Prong[2], hMassVsPtC[kNCharmParticles + 9], activateQA)) { + /// let's build a candidate SigmaC++K- pair + /// and keep it only if: + /// - it has the correct charge (±1) + /// - it is in the correct mass range + + // check the charge for SigmaC++K- candidates + if( std::abs(chargeSc + track.sign()) != 1 ) { + continue; + } - /// let's build a candidate SigmaC++K- and keep it if it is in the correct mass range - /// TODO (check the charge first) - } + // TODO: any PID selection on kaon? + // In case, this selection can be put even before the scope for this trigger, + // to avoid all the computations for SigmaC + + // check the invariant mass + std::array pVecSigmaC = RecoDecay::pVec(pVecFirst, pVecSecond, pVecThird, pVecSoftPi); + float massSigmaCKaonPair = RecoDecay::m( std::array{pVecSigmaC, pVecFourth}, std::array{massSigmaCPlusPlus, massKa} ); + if(massMinSigmaCKaonPair < massSigmaCKaonPair && massSigmaCKaonPair < massMaxSigmaCKaonPair) { + /// This is a good SigmaC++K- event + /// Let's flag it together with SigmaC0K0s + keepEvent[kV0Charm3P] = true; + /// QA plot + if(activateQA) { + float ptSigmaCKaon = RecoDecay::pt(pVecSigmaC, pVecFourth); + hMassVsPtC[kNCharmParticles + 9]->Fill(ptSigmaCKaon, massSigmaCKaonPair); + } + } + } } // end SigmaC++ candidate } // end loop over tracks (soft pi) } // end candidate Lc->pKpi diff --git a/EventFiltering/PWGHF/HFFilterHelpers.h b/EventFiltering/PWGHF/HFFilterHelpers.h index a1003e485f2..db7854511b0 100644 --- a/EventFiltering/PWGHF/HFFilterHelpers.h +++ b/EventFiltering/PWGHF/HFFilterHelpers.h @@ -156,6 +156,8 @@ constexpr float massB0 = o2::constants::physics::MassB0; constexpr float massBs = o2::constants::physics::MassBS; constexpr float massLb = o2::constants::physics::MassLambdaB0; constexpr float massXib = o2::constants::physics::MassXiB0; +constexpr float massSigmaCPlusPlus = o2::constants::physics::MassSigmaCPlusPlus; +constexpr float massSigmaC0 = o2::constants::physics::MassSigmaC0; static const o2::framework::AxisSpec ptAxis{50, 0.f, 50.f}; static const o2::framework::AxisSpec pAxis{50, 0.f, 10.f}; @@ -166,7 +168,7 @@ static const o2::framework::AxisSpec alphaAxis{100, -1.f, 1.f}; static const o2::framework::AxisSpec qtAxis{100, 0.f, 0.25f}; static const o2::framework::AxisSpec bdtAxis{100, 0.f, 1.f}; static const o2::framework::AxisSpec phiAxis{36, 0., o2::constants::math::TwoPI}; -static const std::array massAxisC = {o2::framework::AxisSpec{100, 1.65f, 2.05f}, o2::framework::AxisSpec{100, 1.65f, 2.05f}, o2::framework::AxisSpec{100, 1.75f, 2.15f}, o2::framework::AxisSpec{100, 2.05f, 2.45f}, o2::framework::AxisSpec{100, 2.25f, 2.65f}, o2::framework::AxisSpec{100, 0.139f, 0.159f}, o2::framework::AxisSpec{100, 0.f, 0.25f}, o2::framework::AxisSpec{100, 0.f, 0.25f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{100, 1.1f, 1.4f}, o2::framework::AxisSpec{100, 1.1f, 1.4f}, o2::framework::AxisSpec{100, 2.3f, 2.9f}, o2::framework::AxisSpec{100, 2.3f, 2.9f}}; +static const std::array massAxisC = {o2::framework::AxisSpec{100, 1.65f, 2.05f}, o2::framework::AxisSpec{100, 1.65f, 2.05f}, o2::framework::AxisSpec{100, 1.75f, 2.15f}, o2::framework::AxisSpec{100, 2.05f, 2.45f}, o2::framework::AxisSpec{100, 2.25f, 2.65f}, o2::framework::AxisSpec{100, 0.139f, 0.159f}, o2::framework::AxisSpec{100, 0.f, 0.25f}, o2::framework::AxisSpec{100, 0.f, 0.25f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{100, 1.1f, 1.4f}, o2::framework::AxisSpec{100, 1.1f, 1.4f}, o2::framework::AxisSpec{100, 2.3f, 2.9f}, o2::framework::AxisSpec{100, 2.3f, 2.9f}, o2::framework::AxisSpec{240, 0.14f, 0.2f}, o2::framework::AxisSpec{240, 0.14f, 0.2f}, o2::framework::AxisSpec{330, 2.92, 3.25}, o2::framework::AxisSpec{330, 2.92, 3.25}}; static const std::array massAxisB = {o2::framework::AxisSpec{240, 4.8f, 6.0f}, o2::framework::AxisSpec{240, 4.8f, 6.0f}, o2::framework::AxisSpec{240, 4.8f, 6.0f}, o2::framework::AxisSpec{240, 4.8f, 6.0f}, o2::framework::AxisSpec{240, 5.0f, 6.2f}, o2::framework::AxisSpec{240, 5.0f, 6.2f}}; // default values for configurables @@ -280,7 +282,6 @@ class HfFilterHelper mNSigmaTofPiKaCutForDzero = nSigmaTof; } void setDeltaMassCharmHadForBeauty(float delta) { mDeltaMassCharmHadForBeauty = delta; } - void setDeltaMassLcForSigmaC(float delta) { mDeltaMassLcForSigmaC = delta; } void setV0Selections(float minGammaCosPa, float minK0sLambdaCosPa, float minK0sLambdaRadius, float nSigmaPrFromLambda, float deltaMassK0s, float deltaMassLambda) { mMinGammaCosinePa = minGammaCosPa; @@ -332,6 +333,8 @@ class HfFilterHelper template int8_t isSelectedLcInMassRange(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const float& ptLc, const int8_t isSelected, const int& activateQA, H2 hMassVsPt); template + int8_t isSelectedSigmaCInDeltaMassRange(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const T& pTrackSoftPi, const float ptSigmaC, const int8_t isSelectedLc, H2 hMassVsPt, const int& activateQA); + template int8_t isSelectedXicInMassRange(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const float& ptXic, const int8_t isSelected, const int& activateQA, H2 hMassVsPt); template int8_t isSelectedV0(const V0& v0, const std::array& dauTracks, const Coll& collision, const int& activateQA, H2 hV0Selected, std::array& hArmPod); @@ -396,9 +399,8 @@ class HfFilterHelper float mNSigmaTpcPiKaCutForDzero{3.}; // maximum Nsigma TPC for pions/kaons in D0 decays float mNSigmaTofPiKaCutForDzero{3.}; // maximum Nsigma TOF for pions/kaons in D0 decays float mDeltaMassCharmHadForBeauty{0.04}; // delta mass cut for charm hadrons in B and charm excited decays - float mDeltaMassLcForSigmaC{0.04}; // delta mass cut for Lc candidates to build SigmaC0, ++ candidates float mDeltaMassMinSigmaC{0.155}; // minimum delta mass M(pKpipi)-M(pKpi) of SigmaC0,++ candidates - float mDeltaMassMaxSigmaC{0.18}; // minimum delta mass M(pKpipi)-M(pKpi) of SigmaC0,++ candidates + float mDeltaMassMaxSigmaC{0.18}; // maximum delta mass M(pKpipi)-M(pKpi) of SigmaC0,++ candidates float mMinGammaCosinePa{0.85}; // minimum cosp for gammas float mMinK0sLambdaCosinePa{0.97}; // minimum cosp for K0S and Lambda in charm excited decays float mMinK0sLambdaRadius{0.5}; // minimum radius for K0S and Lambda in charm excited decays @@ -798,31 +800,41 @@ inline int8_t HfFilterHelper::isSelectedLcInMassRange(const T& pTrackSameChargeF return retValue; } -/// Mass selection of Lc candidates to build SigmaC candidates -/// \param massPKPi is the inv. mass calculated assuming p-K-pi daughters -/// \param massPiKP is the inv. mass calculated assuming pi-K-p daughters -inline int8_t HfFilterHelper::isSelectedLcInMassRangeSigmaC(const float& massPKPi, const float& massPiKP) -{ - int8_t retValue = 0; - if (std::fabs(massPKPi - massLc) < mDeltaMassLcForSigmaC) { - retValue |= BIT(0); - } - if (std::fabs(massPiKP - massLc) < mDeltaMassLcForSigmaC) { - retValue |= BIT(1); - } - return retValue; -} - /// Delta mass selection on SigmaC candidates -inline int8_t HfFilterHelper::isSelectedSigmaCInDeltaMassRange(const float& deltaMassPKPi, const float& deltaMassPiKP) +template +inline int8_t HfFilterHelper::isSelectedSigmaCInDeltaMassRange(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const T& pTrackSoftPi, const float ptSigmaC, const int8_t isSelectedLc, H2 hMassVsPt, const int& activateQA) { int8_t retValue = 0; - if (deltaMassPKPi < mDeltaMassMinSigmaC) { - retValue |= BIT(0); + if (TESTBIT(isSelectedLc, 0)) { + /// Lc->pKpi case + auto invMassLcToPKPi = RecoDecay::m(std::array{pTrackSameChargeFirst, pTrackOppositeCharge, pTrackSameChargeSecond}, std::array{massProton, massKa, massPi}); + std::array massDausSigmaCToLcPKPi{massProton, massKa, massPi, massPi}; + float invMassSigmaCToLcPKPi = RecoDecay::m(std::array{pTrackSameChargeFirst, pTrackOppositeCharge, pTrackSameChargeSecond, pTrackSoftPi}, massDausSigmaCToLcPKPi); + float deltaMassPKPi = invMassSigmaCToLcPKPi - invMassLcToPKPi; + if (mDeltaMassMinSigmaC < deltaMassPKPi && deltaMassPKPi < mDeltaMassMaxSigmaC) { + retValue |= BIT(0); + /// QA plot + if(activateQA) { + hMassVsPt->Fill(ptSigmaC, deltaMassPKPi); + } + } } - if (deltaMassPiKP < mDeltaMassMaxSigmaC) { - retValue |= BIT(1); + if (TESTBIT(isSelectedLc, 1)) { + /// Lc->piKp case + auto invMassLcToPiKP = RecoDecay::m(std::array{pTrackSameChargeFirst, pTrackOppositeCharge, pTrackSameChargeSecond}, std::array{massPi, massKa, massProton}); + std::array massDausSigmaCToLcPiKP{massPi, massKa, massProton, massPi}; + float invMassSigmaCToLcPiKP = RecoDecay::m(std::array{pTrackSameChargeFirst, pTrackOppositeCharge, pTrackSameChargeSecond, pTrackSoftPi}, massDausSigmaCToLcPiKP); + float deltaMassPiKP = invMassSigmaCToLcPiKP - invMassLcToPiKP; + if (mDeltaMassMinSigmaC < deltaMassPiKP && deltaMassPiKP < mDeltaMassMaxSigmaC) { + retValue |= BIT(1); + /// QA plot + if(activateQA) { + hMassVsPt->Fill(ptSigmaC, deltaMassPiKP); + } + } } + /// TODO: add QA plot + return retValue; }