diff --git a/EventFiltering/PWGHF/HFFilter.cxx b/EventFiltering/PWGHF/HFFilter.cxx index b5cd0db679f..8bbcbeafcec 100644 --- a/EventFiltering/PWGHF/HFFilter.cxx +++ b/EventFiltering/PWGHF/HFFilter.cxx @@ -147,6 +147,21 @@ struct HfFilter { // Main struct for HF triggers Configurable keepAlsoWrongDmesProtonPairs{"keepAlsoWrongDmesProtonPairs", true, "flat go keep also wrong sign D0p pairs"}; Configurable keepAlsoWrongDstarMesProtonPairs{"keepAlsoWrongDstarMesProtonPairs", true, "flat go keep also wrong sign D*0p pairs"}; + // parameters for Sigma_C + struct : o2::framework::ConfigurableGroup { + Configurable> trkPIDspecies{"trkPIDspecies", std::vector{o2::track::PID::Proton, o2::track::PID::Pion, o2::track::PID::Kaon}, "Trk sel: Particles species for PID, proton, pion, kaon"}; + Configurable> pidTPCMax{"pidTPCMax", std::vector{3., 0., 0.}, "maximum nSigma TPC"}; + Configurable> pidTOFMax{"pidTOFMax", std::vector{3., 0., 0.}, "maximum nSigma TOF"}; + Configurable forceTOF{"forceTOF", false, "fill PID information for associated tracks"}; + Configurable tofPIDThreshold{"tofPIDThreshold", 1.0, "minimum pT after which TOF PID is applicable"}; + Configurable minPtProton{"minPtProton", 0.39, "minimum pT for associated proton"}; + Configurable maxPtProton{"maxPtProton", 4.51, "maximum pT for associated proton"}; + Configurable minMassSigmaCCorr{"minMassSigmaCCorr", 0.15, "minimum mass of SigmaC for correlation with proton"}; + Configurable maxMassSigmaCCorr{"maxMassSigmaCCorr", 0.19, "maximum mass of SigmaC for correlation with proton"}; + Configurable minPtSigmaC{"minPtSigmaC", 4.99, "minimum pT of SigmaC for correlation with proton"}; + Configurable maxPtSigmaC{"maxPtSigmaC", 12.0, "maximum pT of SigmaC for correlation with proton"}; + } configSigmaC; + // parameters for charm baryons to Xi bachelor Configurable> cutsXiCascades{"cutsXiCascades", {cutsCascades[0], 1, 8, labelsEmpty, labelsColumnsCascades}, "Selections for cascades (Xi) for Xi+bachelor triggers"}; Configurable> cutsXiBachelor{"cutsXiBachelor", {cutsCharmBaryons[0], 1, 11, labelsEmpty, labelsColumnsCharmBarCuts}, "Selections for charm baryons (Xi+Pi, Xi+Ka, Xi+Pi+Pi)"}; @@ -161,6 +176,7 @@ struct HfFilter { // Main struct for HF triggers Configurable> thresholdBDTScoreDSToPiKK{"thresholdBDTScoreDSToPiKK", {hf_cuts_bdt_multiclass::Cuts[0], hf_cuts_bdt_multiclass::NBinsPt, hf_cuts_bdt_multiclass::NCutBdtScores, hf_cuts_bdt_multiclass::labelsPt, hf_cuts_bdt_multiclass::labelsCutBdt}, "Threshold values for BDT output scores of Ds+ candidates"}; Configurable> thresholdBDTScoreLcToPiKP{"thresholdBDTScoreLcToPiKP", {hf_cuts_bdt_multiclass::Cuts[0], hf_cuts_bdt_multiclass::NBinsPt, hf_cuts_bdt_multiclass::NCutBdtScores, hf_cuts_bdt_multiclass::labelsPt, hf_cuts_bdt_multiclass::labelsCutBdt}, "Threshold values for BDT output scores of Lc+ candidates"}; Configurable> thresholdBDTScoreXicToPiKP{"thresholdBDTScoreXicToPiKP", {hf_cuts_bdt_multiclass::Cuts[0], hf_cuts_bdt_multiclass::NBinsPt, hf_cuts_bdt_multiclass::NCutBdtScores, hf_cuts_bdt_multiclass::labelsPt, hf_cuts_bdt_multiclass::labelsCutBdt}, "Threshold values for BDT output scores of Xic+ candidates"}; + Configurable> thresholdBDTScoreScLcToPiKP{"thresholdBDTScoreScLcToPiKP", {hf_cuts_bdt_multiclass::Cuts[0], hf_cuts_bdt_multiclass::NBinsPt, hf_cuts_bdt_multiclass::NCutBdtScores, hf_cuts_bdt_multiclass::labelsPt, hf_cuts_bdt_multiclass::labelsCutBdt}, "Threshold values for BDT output scores of Lc<--Sc candidates"}; Configurable acceptBdtBkgOnly{"acceptBdtBkgOnly", true, "Enable / disable selection based on BDT bkg score only"}; @@ -206,7 +222,7 @@ struct HfFilter { // Main struct for HF triggers int currentRun{0}; // needed to detect if the run changed and trigger update of calibrations etc. // array of BDT thresholds - std::array, kNCharmParticles> thresholdBDTScores; + std::array, kNCharmParticles + 1> thresholdBDTScores; o2::vertexing::DCAFitterN<2> df2; // fitter for Charm Hadron vertex (2-prong vertex fitter) o2::vertexing::DCAFitterN<3> df3; // fitter for Charm/Beauty Hadron vertex (3-prong vertex fitter) @@ -224,7 +240,7 @@ struct HfFilter { // Main struct for HF triggers std::array, kNCharmParticles> hCharmProtonKstarDistr{}; std::array, kNCharmParticles> hCharmDeuteronKstarDistr{}; std::array, nTotBeautyParts> hMassVsPtB{}; - std::array, kNCharmParticles + 23> hMassVsPtC{}; // +9 for resonances (D*+, D*0, Ds*+, Ds1+, Ds2*+, Xic+* right sign, Xic+* wrong sign, Xic0* right sign, Xic0* wrong sign) +2 for SigmaC (SigmaC++, SigmaC0) +2 for SigmaCK pairs (SigmaC++K-, SigmaC0K0s) +3 for charm baryons (Xi+Pi, Xi+Ka, Xi+Pi+Pi) + JPsi + 4 for charm baryons (D0+p, D0+pWrongSign, D*0p, D*0+pWrongSign) + std::array, kNCharmParticles + 24> hMassVsPtC{}; // +9 for resonances (D*+, D*0, Ds*+, Ds1+, Ds2*+, Xic+* right sign, Xic+* wrong sign, Xic0* right sign, Xic0* wrong sign) +2 for SigmaC (SigmaC++, SigmaC0) +2 for SigmaCK pairs (SigmaC++K-, SigmaC0K0s) +3 for charm baryons (Xi+Pi, Xi+Ka, Xi+Pi+Pi) + JPsi + 4 for charm baryons (D0+p, D0+pWrongSign, D*0p, D*0+pWrongSign) std::array, 4> hPrDePID; // proton TPC, proton TOF, deuteron TPC, deuteron TOF std::array, kNCharmParticles> hBDTScoreBkg{}; std::array, kNCharmParticles> hBDTScorePrompt{}; @@ -369,6 +385,8 @@ struct HfFilter { // Main struct for HF triggers // ThetaC hMassVsPtC[kNCharmParticles + 21] = registry.add("fMassVsPtCharmBaryonToDstarP", "#it{M} vs. #it{p}_{T} distribution of triggered D^{*0}#p candidates;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", HistType::kTH2D, {ptAxis, massAxisC[kNCharmParticles + 21]}); hMassVsPtC[kNCharmParticles + 22] = registry.add("fMassVsPtCharmBaryonToDstarPWrongSign", "#it{M} vs. #it{p}_{T} distribution of triggered D^{*0}#p wrong sign candidates;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", HistType::kTH2D, {ptAxis, massAxisC[kNCharmParticles + 22]}); + // SigmaC-p + hMassVsPtC[kNCharmParticles + 23] = registry.add("fMassVsPtSigmaCPr", "#it{M} vs. #it{p}_{T} distribution of #Sigma_{c} for SigmaCPr trigger;#it{p}_{T} (GeV/#it{c});#it{M} (GeV/#it{c}^{2});counts", HistType::kTH2D, {ptAxis, massAxisC[kNCharmParticles + 23]}); 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::kTH2D, {ptAxis, massAxisB[iBeautyPart]}); @@ -431,7 +449,7 @@ struct HfFilter { // Main struct for HF triggers ccdbApi.init(url); lut = o2::base::MatLayerCylSet::rectifyPtrFromFile(ccdb->get("GLO/Param/MatLUT")); - thresholdBDTScores = {thresholdBDTScoreD0ToKPi, thresholdBDTScoreDPlusToPiKPi, thresholdBDTScoreDSToPiKK, thresholdBDTScoreLcToPiKP, thresholdBDTScoreXicToPiKP}; + thresholdBDTScores = {thresholdBDTScoreD0ToKPi, thresholdBDTScoreDPlusToPiKPi, thresholdBDTScoreDSToPiKK, thresholdBDTScoreLcToPiKP, thresholdBDTScoreXicToPiKP, thresholdBDTScoreScLcToPiKP}; } void process(CollsWithEvSel const& collisions, @@ -455,7 +473,7 @@ struct HfFilter { // Main struct for HF triggers bool isSelectedITSROFBorder = evSel.applyITSROFBorderCut ? collision.selection_bit(o2::aod::evsel::kNoITSROFrameBorder) : true; bool isSelectedPvZ = (std::fabs(collision.posZ()) < evSel.maxPvZ); if (!isSelectedTVX || !isSelectedTFBorder || !isSelectedITSROFBorder || !isSelectedPvZ) { - tags(keepEvent[kHighPt2P], keepEvent[kHighPt3P], keepEvent[kBeauty3P], keepEvent[kBeauty4P], keepEvent[kFemto2P], keepEvent[kFemto3P], keepEvent[kDoubleCharm2P], keepEvent[kDoubleCharm3P], keepEvent[kDoubleCharmMix], keepEvent[kV0Charm2P], keepEvent[kV0Charm3P], keepEvent[kCharmBarToXiBach], keepEvent[kSigmaCPPK], keepEvent[kSigmaC0K0], keepEvent[kPhotonCharm2P], keepEvent[kPhotonCharm3P], keepEvent[kSingleCharm2P], keepEvent[kSingleCharm3P], keepEvent[kSingleNonPromptCharm2P], keepEvent[kSingleNonPromptCharm3P], keepEvent[kCharmBarToXi2Bach], keepEvent[kPrCharm2P], keepEvent[kBtoJPsiKa], keepEvent[kBtoJPsiKstar], keepEvent[kBtoJPsiPhi], keepEvent[kBtoJPsiPrKa], keepEvent[kBtoJPsiPi]); + tags(keepEvent[kHighPt2P], keepEvent[kHighPt3P], keepEvent[kBeauty3P], keepEvent[kBeauty4P], keepEvent[kFemto2P], keepEvent[kFemto3P], keepEvent[kDoubleCharm2P], keepEvent[kDoubleCharm3P], keepEvent[kDoubleCharmMix], keepEvent[kV0Charm2P], keepEvent[kV0Charm3P], keepEvent[kCharmBarToXiBach], keepEvent[kSigmaCPPK], keepEvent[kSigmaC0K0], keepEvent[kSigmaCPr], keepEvent[kPhotonCharm2P], keepEvent[kPhotonCharm3P], keepEvent[kSingleCharm2P], keepEvent[kSingleCharm3P], keepEvent[kSingleNonPromptCharm2P], keepEvent[kSingleNonPromptCharm3P], keepEvent[kCharmBarToXi2Bach], keepEvent[kPrCharm2P], keepEvent[kBtoJPsiKa], keepEvent[kBtoJPsiKstar], keepEvent[kBtoJPsiPhi], keepEvent[kBtoJPsiPrKa], keepEvent[kBtoJPsiPi]); continue; } @@ -1507,8 +1525,12 @@ struct HfFilter { // Main struct for HF triggers } } // end femto selection - // SigmaC++ K- trigger - if (!keepEvent[kSigmaCPPK] && is3Prong[2] > 0 && is3ProngInMass[2] > 0 && isSignalTagged[2] > 0 && helper.isSelectedKaonFromXicResoToSigmaC(track)) { + // SigmaC++ K- and SigmaC++,0 - p trigger + + bool isTrackKaon = helper.isSelectedKaonFromXicResoToSigmaC(track); + bool isTrackProton = helper.isSelectedTrack4Corr(track, configSigmaC.trkPIDspecies, configSigmaC.pidTPCMax, configSigmaC.pidTOFMax, configSigmaC.minPtProton, configSigmaC.maxPtProton, configSigmaC.tofPIDThreshold, configSigmaC.forceTOF); + + if ((!keepEvent[kSigmaCPPK] || !keepEvent[kSigmaCPr]) && is3Prong[2] > 0 && is3ProngInMass[2] > 0 && isSignalTagged[2] > 0 && (isTrackKaon || isTrackProton)) { // we need a candidate Lc->pKpi and a candidate soft kaon // look for SigmaC++ candidates @@ -1523,8 +1545,7 @@ struct HfFilter { // Main struct for HF triggers // do not consider as candidate soft pion a track already used to build the current 3-prong candidate continue; } - - // exclude already the current track if it corresponds to the K- candidate + // exclude already the current track if it corresponds to the K- or proton candidate if (globalIndexSoftPi == track.globalIndex()) { continue; } @@ -1532,9 +1553,6 @@ struct HfFilter { // Main struct for HF triggers // 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); @@ -1548,8 +1566,7 @@ struct HfFilter { // Main struct for HF triggers } int16_t isSoftPionSelected = helper.isSelectedTrackForSoftPionOrBeauty(trackSoftPi, trackParSoftPi, dcaSoftPi); if (TESTBIT(isSoftPionSelected, kSoftPionForSigmaC) /*&& (TESTBIT(is3Prong[2], 0) || TESTBIT(is3Prong[2], 1))*/) { - - // check the mass of the SigmaC++ candidate + // check the mass of the SigmaC++,0 candidate auto pVecSigmaC = RecoDecay::pVec(pVecFirst, pVecSecond, pVecThird, pVecSoftPi); auto ptSigmaC = RecoDecay::pt(pVecSigmaC); int8_t whichSigmaC = helper.isSelectedSigmaCInDeltaMassRange<2>(pVecFirst, pVecThird, pVecSecond, pVecSoftPi, ptSigmaC, is3Prong[2], hMassVsPtC[kNCharmParticles + 9], activateQA); @@ -1560,51 +1577,59 @@ struct HfFilter { // Main struct for HF triggers /// - it is in the correct mass range // check the charge for SigmaC++K- candidates - if (std::abs(chargeSc + track.sign()) != 1) { - continue; - } - // check the invariant mass - float massSigmaCPKPi{-999.}, massSigmaCPiKP{-999.}, deltaMassXicResoPKPi{-999.}, deltaMassXicResoPiKP{-999.}; - float ptSigmaCKaon = RecoDecay::pt(pVecSigmaC, pVecFourth); + if (!keepEvent[kSigmaCPPK] && (std::abs(chargeSc + track.sign()) == 1 && std::abs(chargeSc) == 2)) { + // check the invariant mass + float massSigmaCPKPi{-999.}, massSigmaCPiKP{-999.}, deltaMassXicResoPKPi{-999.}, deltaMassXicResoPiKP{-999.}; + float ptSigmaCKaon = RecoDecay::pt(pVecSigmaC, pVecFourth); - if (ptSigmaCKaon > cutsPtDeltaMassCharmReso->get(2u, 10u)) { - if (TESTBIT(whichSigmaC, 0)) { - massSigmaCPKPi = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi}, std::array{massProton, massKa, massPi, massPi}); - deltaMassXicResoPKPi = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi, pVecFourth}, std::array{massProton, massKa, massPi, massPi, massKa}) - massSigmaCPKPi; - } - if (TESTBIT(whichSigmaC, 1)) { - massSigmaCPiKP = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi}, std::array{massPi, massKa, massProton, massPi}); - deltaMassXicResoPiKP = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi, pVecFourth}, std::array{massPi, massKa, massProton, massPi, massKa}) - massSigmaCPiKP; - } - bool isPKPiOk = (cutsPtDeltaMassCharmReso->get(0u, 10u) < deltaMassXicResoPKPi && deltaMassXicResoPKPi < cutsPtDeltaMassCharmReso->get(1u, 10u)); - bool isPiKPOk = (cutsPtDeltaMassCharmReso->get(0u, 10u) < deltaMassXicResoPiKP && deltaMassXicResoPiKP < cutsPtDeltaMassCharmReso->get(1u, 10u)); - if (isPKPiOk || isPiKPOk) { - /// This is a good SigmaC++K- event - keepEvent[kSigmaCPPK] = true; + if (ptSigmaCKaon > cutsPtDeltaMassCharmReso->get(2u, 10u)) { + if (TESTBIT(whichSigmaC, 0)) { + massSigmaCPKPi = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi}, std::array{massProton, massKa, massPi, massPi}); + deltaMassXicResoPKPi = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi, pVecFourth}, std::array{massProton, massKa, massPi, massPi, massKa}) - massSigmaCPKPi; + } + if (TESTBIT(whichSigmaC, 1)) { + massSigmaCPiKP = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi}, std::array{massPi, massKa, massProton, massPi}); + deltaMassXicResoPiKP = RecoDecay::m(std::array{pVecFirst, pVecSecond, pVecThird, pVecSoftPi, pVecFourth}, std::array{massPi, massKa, massProton, massPi, massKa}) - massSigmaCPiKP; + } + bool isPKPiOk = (cutsPtDeltaMassCharmReso->get(0u, 10u) < deltaMassXicResoPKPi && deltaMassXicResoPKPi < cutsPtDeltaMassCharmReso->get(1u, 10u)); + bool isPiKPOk = (cutsPtDeltaMassCharmReso->get(0u, 10u) < deltaMassXicResoPiKP && deltaMassXicResoPiKP < cutsPtDeltaMassCharmReso->get(1u, 10u)); + if ((isPKPiOk || isPiKPOk) && isTrackKaon) { + /// This is a good SigmaC++K- event + keepEvent[kSigmaCPPK] = true; - /// QA plot - if (activateQA) { - if (isPKPiOk) { - if (TESTBIT(whichSigmaC, 2)) { - hMassVsPtC[kNCharmParticles + 11]->Fill(ptSigmaCKaon, deltaMassXicResoPKPi); - } - if (TESTBIT(whichSigmaC, 3)) { - hMassVsPtC[kNCharmParticles + 12]->Fill(ptSigmaCKaon, deltaMassXicResoPKPi); - } - } - if (isPiKPOk) { - if (TESTBIT(whichSigmaC, 2)) { - hMassVsPtC[kNCharmParticles + 11]->Fill(ptSigmaCKaon, deltaMassXicResoPiKP); + /// QA plot + if (activateQA) { + if (isPKPiOk) { + if (TESTBIT(whichSigmaC, 2)) { + hMassVsPtC[kNCharmParticles + 11]->Fill(ptSigmaCKaon, deltaMassXicResoPKPi); + } + if (TESTBIT(whichSigmaC, 3)) { + hMassVsPtC[kNCharmParticles + 12]->Fill(ptSigmaCKaon, deltaMassXicResoPKPi); + } } - if (TESTBIT(whichSigmaC, 3)) { - hMassVsPtC[kNCharmParticles + 12]->Fill(ptSigmaCKaon, deltaMassXicResoPiKP); + if (isPiKPOk) { + if (TESTBIT(whichSigmaC, 2)) { + hMassVsPtC[kNCharmParticles + 11]->Fill(ptSigmaCKaon, deltaMassXicResoPiKP); + } + if (TESTBIT(whichSigmaC, 3)) { + hMassVsPtC[kNCharmParticles + 12]->Fill(ptSigmaCKaon, deltaMassXicResoPiKP); + } } } } } } } + float deltaEta = std::abs(RecoDecay::eta(pVecSigmaC) - track.eta()); + if (!keepEvent[kSigmaCPr] && (isTrackProton && deltaEta < 1.0 && pt3Prong > 3.0)) { + + auto tagBDT = helper.isBDTSelected(scores[2], thresholdBDTScores[5]); + + if (helper.selectionSigmaCForScPCorr(pVecFirst, pVecThird, pVecSecond, pVecSoftPi, ptSigmaC, is3Prong[2], hMassVsPtC[kNCharmParticles + 23], activateQA, configSigmaC.minMassSigmaCCorr, configSigmaC.maxMassSigmaCCorr, configSigmaC.minPtSigmaC, configSigmaC.maxPtSigmaC) && TESTBIT(tagBDT, RecoDecay::OriginType::Prompt)) { + keepEvent[kSigmaCPr] = true; + } + } } // end SigmaC++ candidate } // end loop over tracks (soft pi) } // end candidate Lc->pKpi @@ -2001,7 +2026,7 @@ struct HfFilter { // Main struct for HF triggers } } - tags(keepEvent[kHighPt2P], keepEvent[kHighPt3P], keepEvent[kBeauty3P], keepEvent[kBeauty4P], keepEvent[kFemto2P], keepEvent[kFemto3P], keepEvent[kDoubleCharm2P], keepEvent[kDoubleCharm3P], keepEvent[kDoubleCharmMix], keepEvent[kV0Charm2P], keepEvent[kV0Charm3P], keepEvent[kCharmBarToXiBach], keepEvent[kSigmaCPPK], keepEvent[kSigmaC0K0], keepEvent[kPhotonCharm2P], keepEvent[kPhotonCharm3P], keepEvent[kSingleCharm2P], keepEvent[kSingleCharm3P], keepEvent[kSingleNonPromptCharm2P], keepEvent[kSingleNonPromptCharm3P], keepEvent[kCharmBarToXi2Bach], keepEvent[kPrCharm2P], keepEvent[kBtoJPsiKa], keepEvent[kBtoJPsiKstar], keepEvent[kBtoJPsiPhi], keepEvent[kBtoJPsiPrKa], keepEvent[kBtoJPsiPi]); + tags(keepEvent[kHighPt2P], keepEvent[kHighPt3P], keepEvent[kBeauty3P], keepEvent[kBeauty4P], keepEvent[kFemto2P], keepEvent[kFemto3P], keepEvent[kDoubleCharm2P], keepEvent[kDoubleCharm3P], keepEvent[kDoubleCharmMix], keepEvent[kV0Charm2P], keepEvent[kV0Charm3P], keepEvent[kCharmBarToXiBach], keepEvent[kSigmaCPPK], keepEvent[kSigmaC0K0], keepEvent[kSigmaCPr], keepEvent[kPhotonCharm2P], keepEvent[kPhotonCharm3P], keepEvent[kSingleCharm2P], keepEvent[kSingleCharm3P], keepEvent[kSingleNonPromptCharm2P], keepEvent[kSingleNonPromptCharm3P], keepEvent[kCharmBarToXi2Bach], keepEvent[kPrCharm2P], keepEvent[kBtoJPsiKa], keepEvent[kBtoJPsiKstar], keepEvent[kBtoJPsiPhi], keepEvent[kBtoJPsiPrKa], keepEvent[kBtoJPsiPi]); if (!std::accumulate(keepEvent, keepEvent + kNtriggersHF, 0)) { hProcessedEvents->Fill(1); diff --git a/EventFiltering/PWGHF/HFFilterHelpers.h b/EventFiltering/PWGHF/HFFilterHelpers.h index 04297f1ff09..a80d2252d6c 100644 --- a/EventFiltering/PWGHF/HFFilterHelpers.h +++ b/EventFiltering/PWGHF/HFFilterHelpers.h @@ -28,13 +28,15 @@ // #include "Common/Core/RecoDecay.h" #include "Common/Core/trackUtilities.h" +#include "Common/DataModel/PIDResponseTOF.h" +#include "Common/DataModel/PIDResponseTPC.h" +#include "ReconstructionDataFormats/PID.h" #include #include #include #include #include -#include #include #include #include @@ -43,6 +45,7 @@ #include #include #include +#include #include #include // IWYU pragma: keep (do not replace with Math/Vector4Dfwd.h) @@ -97,6 +100,7 @@ enum HfTriggers { kBtoJPsiPhi, kBtoJPsiPrKa, kBtoJPsiPi, + kSigmaCPr, kNtriggersHF }; @@ -245,7 +249,7 @@ static const int nTotBeautyParts = static_cast(kNBeautyParticles) + static_ static const std::array beautyParticleNames{"Bplus", "B0toDStar", "Bc", "B0", "Bs", "Lb", "Xib", "BplusToJPsi", "B0ToJPsi", "BsToJPsi", "LbToJPsi", "BcToJPsi"}; static const std::array pdgCodesCharm{421, 411, 431, 4122, 4232}; static const std::array eventTitles = {"all", "rejected"}; -static const std::vector hfTriggerNames{filtering::HfHighPt2P::columnLabel(), filtering::HfHighPt3P::columnLabel(), filtering::HfBeauty3P::columnLabel(), filtering::HfBeauty4P::columnLabel(), filtering::HfFemto2P::columnLabel(), filtering::HfFemto3P::columnLabel(), filtering::HfDoubleCharm2P::columnLabel(), filtering::HfDoubleCharm3P::columnLabel(), filtering::HfDoubleCharmMix::columnLabel(), filtering::HfV0Charm2P::columnLabel(), filtering::HfV0Charm3P::columnLabel(), filtering::HfCharmBarToXiBach::columnLabel(), filtering::HfSigmaCPPK::columnLabel(), filtering::HfSigmaC0K0::columnLabel(), filtering::HfPhotonCharm2P::columnLabel(), filtering::HfPhotonCharm3P::columnLabel(), filtering::HfSingleCharm2P::columnLabel(), filtering::HfSingleCharm3P::columnLabel(), filtering::HfSingleNonPromptCharm2P::columnLabel(), filtering::HfSingleNonPromptCharm3P::columnLabel(), filtering::HfCharmBarToXi2Bach::columnLabel(), filtering::HfPrCharm2P::columnLabel(), filtering::HfBtoJPsiKa::columnLabel(), filtering::HfBtoJPsiKstar::columnLabel(), filtering::HfBtoJPsiPhi::columnLabel(), filtering::HfBtoJPsiPrKa::columnLabel(), filtering::HfBtoJPsiPi::columnLabel()}; +static const std::vector hfTriggerNames{filtering::HfHighPt2P::columnLabel(), filtering::HfHighPt3P::columnLabel(), filtering::HfBeauty3P::columnLabel(), filtering::HfBeauty4P::columnLabel(), filtering::HfFemto2P::columnLabel(), filtering::HfFemto3P::columnLabel(), filtering::HfDoubleCharm2P::columnLabel(), filtering::HfDoubleCharm3P::columnLabel(), filtering::HfDoubleCharmMix::columnLabel(), filtering::HfV0Charm2P::columnLabel(), filtering::HfV0Charm3P::columnLabel(), filtering::HfCharmBarToXiBach::columnLabel(), filtering::HfSigmaCPPK::columnLabel(), filtering::HfSigmaC0K0::columnLabel(), filtering::HfPhotonCharm2P::columnLabel(), filtering::HfPhotonCharm3P::columnLabel(), filtering::HfSingleCharm2P::columnLabel(), filtering::HfSingleCharm3P::columnLabel(), filtering::HfSingleNonPromptCharm2P::columnLabel(), filtering::HfSingleNonPromptCharm3P::columnLabel(), filtering::HfCharmBarToXi2Bach::columnLabel(), filtering::HfPrCharm2P::columnLabel(), filtering::HfBtoJPsiKa::columnLabel(), filtering::HfBtoJPsiKstar::columnLabel(), filtering::HfBtoJPsiPhi::columnLabel(), filtering::HfBtoJPsiPrKa::columnLabel(), filtering::HfBtoJPsiPi::columnLabel(), filtering::HfSigmaCPr::columnLabel()}; static const std::array v0Labels{"#gamma", "K_{S}^{0}", "#Lambda", "#bar{#Lambda}"}; static const std::array v0Names{"Photon", "K0S", "Lambda", "AntiLambda"}; @@ -293,7 +297,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{250, 1.65f, 2.15f}, o2::framework::AxisSpec{250, 1.65f, 2.15f}, o2::framework::AxisSpec{250, 1.75f, 2.25f}, o2::framework::AxisSpec{250, 2.05f, 2.55f}, o2::framework::AxisSpec{250, 2.25f, 2.75f}, o2::framework::AxisSpec{200, 0.139f, 0.159f}, o2::framework::AxisSpec{250, 0.f, 0.25f}, o2::framework::AxisSpec{250, 0.f, 0.25f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{170, 0.13f, 0.3f}, o2::framework::AxisSpec{170, 0.13f, 0.3f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{240, 2.4f, 3.6f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}}; +static const std::array massAxisC = {o2::framework::AxisSpec{250, 1.65f, 2.15f}, o2::framework::AxisSpec{250, 1.65f, 2.15f}, o2::framework::AxisSpec{250, 1.75f, 2.25f}, o2::framework::AxisSpec{250, 2.05f, 2.55f}, o2::framework::AxisSpec{250, 2.25f, 2.75f}, o2::framework::AxisSpec{200, 0.139f, 0.159f}, o2::framework::AxisSpec{250, 0.f, 0.25f}, o2::framework::AxisSpec{250, 0.f, 0.25f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 0.48f, 0.88f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{200, 1.1f, 1.4f}, o2::framework::AxisSpec{170, 0.13f, 0.3f}, o2::framework::AxisSpec{170, 0.13f, 0.3f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{200, 0.4f, 0.8f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{350, 2.3f, 3.0f}, o2::framework::AxisSpec{240, 2.4f, 3.6f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.7f, 1.3f}, o2::framework::AxisSpec{300, 0.14f, 0.26f}}; static const std::array massAxisB = {o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 5.4f, 7.4f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 4.4f, 6.4f}, o2::framework::AxisSpec{400, 5.0f, 6.6f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{500, 4.2f, 6.2f}, o2::framework::AxisSpec{400, 5.0f, 6.6f}, o2::framework::AxisSpec{240, 5.8f, 7.0f}}; // default values for configurables @@ -625,6 +629,8 @@ class HfFilterHelper int16_t isSelectedTrackForSoftPionOrBeauty(const T& track, const T1& trackPar, const T2& dca); template bool isSelectedTrack4Femto(const T1& track, const T2& trackPar, const int& activateQA, H2 hTPCPID, H2 hTOFPID, const int& trackSpecies); + template + bool isSelectedTrack4Corr(Atrack const& track, SpeciesContainer const mPIDspecies, T1 const maxTPC, T2 const maxTOF, float minPt = 0.39, float maxPt = 4.6, float ptThreshold = 1.0, bool tofForced = false); template int8_t isDzeroPreselected(const T& trackPos, const T& trackNeg); template @@ -644,6 +650,8 @@ class HfFilterHelper 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 + bool selectionSigmaCForScPCorr(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const T& pTrackSoftPi, const float ptSigmaC, const int8_t isSelectedLc, H2 hMassVsPt, const int& activateQA, float mDeltaMassMinSigmaC = 0.155, float mDeltaMassMaxSigmaC = 0.2, float mPtMinSigmaC = 4.99, float mPtMaxSigmaC = 12.0); + 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 int& activateQA, H2 hV0Selected, std::array& hArmPod); @@ -867,7 +875,7 @@ inline int16_t HfFilterHelper::isSelectedTrackForSoftPionOrBeauty(const T& track return kRejected; } - if constexpr (whichTrigger == kSigmaCPPK || whichTrigger == kSigmaC0K0) { + if constexpr (whichTrigger == kSigmaCPPK || whichTrigger == kSigmaC0K0 || whichTrigger == kSigmaCPr) { // SigmaC0,++ soft pion pt cut if (pT < mPtMinSoftPionForSigmaC || pT > mPtMaxSoftPionForSigmaC) { @@ -922,6 +930,110 @@ inline int16_t HfFilterHelper::isSelectedTrackForSoftPionOrBeauty(const T& track return retValue; } +/// Basic selection of proton or deuteron candidates +/// \param track is a track +/// \param mPIDspecies is a vector of different particle species +/// \param activateQA flag to activate the filling of QA histos +/// \param maxTPC is a vector of max TPCnSigma for different particle species +/// \param maxTOF is a vector of max TOFnSigma for different particle species +/// \param hProtonTPCPID histo with NsigmaTPC vs. p +/// \param hProtonTOFPID histo with NsigmaTOF vs. p +/// \return true if track passes all cuts +template +inline bool HfFilterHelper::isSelectedTrack4Corr(Atrack const& track, SpeciesContainer const mPIDspecies, + T1 const maxTPC, T2 const maxTOF, float minPt, float maxPt, float ptThreshold, bool tofForced) +{ + // Ensure size consistency + if (mPIDspecies.value.size() != maxTPC.value.size() || mPIDspecies.value.size() != maxTOF.value.size()) { + LOGF(error, "Size of particle species and corresponding nSigma selection arrays should be the same"); + return false; // Early exit on error + } + + if (!track.isGlobalTrackWoDCA()) { + return false; + } + + if (track.pt() < minPt || track.pt() > maxPt) { + return false; + } + + for (size_t speciesIndex = 0; speciesIndex < mPIDspecies.value.size(); ++speciesIndex) { + float nSigmaTPC; + auto const& pid = mPIDspecies->at(speciesIndex); + + nSigmaTPC = o2::aod::pidutils::tpcNSigma(pid, track); + + if (track.pt() > ptThreshold && tofForced && !track.hasTOF()) + return false; + + int parSpecies = -1; + float tpcNCls = track.tpcNClsFound(); + float tpcPin = track.tpcInnerParam(); + float eta = track.eta(); + + if (pid == o2::track::PID::Proton) { + parSpecies = kPr; + } else if (pid == o2::track::PID::Kaon) { + parSpecies = kKa; + } else if (pid == o2::track::PID::Pion) { + parSpecies = kPi; + } else { + LOGF(fatal, "particle species is not defined in isSelectedTrack4Corr"); + } + + // 2. Apply nSigmaTPC using selected calibration + if (mTpcPidCalibrationOption == 1) { + // Option 1: post-calibration → no sign dependence + nSigmaTPC = getTPCPostCalib(tpcPin, tpcNCls, eta, nSigmaTPC, parSpecies); + } + if (mTpcPidCalibrationOption == 2) { + float dEdx = track.tpcSignal(); + // Option 2: spline calibration → charge-dependent + if (track.sign() > 0) { + // Positive track + nSigmaTPC = getTPCSplineCalib(tpcPin, dEdx, parSpecies); + } else { + // Negative track + if (pid == o2::track::PID::Proton) { + parSpecies = kAntiPr; + } else if (pid == o2::track::PID::Kaon) { + parSpecies = kAntiKa; + } else if (pid == o2::track::PID::Pion) { + parSpecies = kAntiPi; + } else { + LOGF(fatal, "particle species is not defined in isSelectedTrack4Corr"); + } + nSigmaTPC = getTPCSplineCalib(tpcPin, dEdx, parSpecies); + } + } + + if (speciesIndex == 0) { // First species logic + + if (std::abs(nSigmaTPC) > maxTPC->at(speciesIndex)) { + return false; // TPC check failed + } + if (track.hasTOF()) { + auto nSigmaTOF = o2::aod::pidutils::tofNSigma(pid, track); + if (std::abs(nSigmaTOF) > maxTOF->at(speciesIndex)) { + return false; // TOF check failed + } + } + } else { // Other species logic + if (std::abs(nSigmaTPC) < maxTPC->at(speciesIndex)) { // Check TPC nSigma first + if (track.hasTOF()) { + auto nSigmaTOF = o2::aod::pidutils::tofNSigma(pid, track); + if (std::abs(nSigmaTOF) < maxTOF->at(speciesIndex)) { + return false; // Reject if both TPC and TOF are within thresholds + } + } else { + return false; // Reject if only TPC is within threshold and TOF is unavailable + } + } + } + } + return true; // Passed all checks +} + /// Basic selection of proton or deuteron candidates /// \param track is a track /// \param trackPar is a track parameter @@ -1312,6 +1424,48 @@ inline int8_t HfFilterHelper::isSelectedLcInMassRange(const T& pTrackSameChargeF return retValue; } +/// Delta mass selection on SigmaC candidates for correlation +template +inline bool HfFilterHelper::selectionSigmaCForScPCorr(const T& pTrackSameChargeFirst, const T& pTrackSameChargeSecond, const T& pTrackOppositeCharge, const T& pTrackSoftPi, const float ptSigmaC, const int8_t isSelectedLc, H2 hMassVsPt, const int& activateQA, float mDeltaMassMinSigmaC, float mDeltaMassMaxSigmaC, float mPtMinSigmaC, float mPtMaxSigmaC) +{ + if (ptSigmaC < mPtMinSigmaC || ptSigmaC > mPtMaxSigmaC) { + return false; + } + bool isSigmaCSelected{false}; + 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; + isSigmaCSelected = (mDeltaMassMinSigmaC < deltaMassPKPi && deltaMassPKPi < mDeltaMassMaxSigmaC); + if (isSigmaCSelected) { + if (activateQA) { + hMassVsPt->Fill(ptSigmaC, deltaMassPKPi); + } + return true; + } + } + 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; + + isSigmaCSelected = (mDeltaMassMinSigmaC < deltaMassPiKP && deltaMassPiKP < mDeltaMassMaxSigmaC); + if (isSigmaCSelected) { + if (activateQA) { + hMassVsPt->Fill(ptSigmaC, deltaMassPiKP); + } + return true; + } + } + + return isSigmaCSelected; + /// TODO: add QA plot +} + /// Delta mass selection on SigmaC candidates 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) diff --git a/EventFiltering/filterTables.h b/EventFiltering/filterTables.h index 736e397a647..e0908aeb18a 100644 --- a/EventFiltering/filterTables.h +++ b/EventFiltering/filterTables.h @@ -92,6 +92,7 @@ DECLARE_SOA_COLUMN(HfCharmBarToXi2Bach, hasHfCharmBarToXi2Bach, bool); DECLARE_SOA_COLUMN(HfPrCharm2P, hasHfPrCharm2P, bool); //! Charm baryon to 2-prong + bachelors DECLARE_SOA_COLUMN(HfSigmaCPPK, hasHfSigmaCPPK, bool); //! SigmaC(2455)++K- and SigmaC(2520)++K- + c.c. DECLARE_SOA_COLUMN(HfSigmaC0K0, hasHfSigmaC0K0, bool); //! SigmaC(2455)0KS0 and SigmaC(2520)0KS0 +DECLARE_SOA_COLUMN(HfSigmaCPr, hasHfSigmaCPr, bool); //! SigmaC(2455)Proton pairs DECLARE_SOA_COLUMN(HfPhotonCharm2P, hasHfPhotonCharm2P, bool); //! photon with 2-prong charm hadron DECLARE_SOA_COLUMN(HfPhotonCharm3P, hasHfPhotonCharm3P, bool); //! photon with 3-prong charm hadron DECLARE_SOA_COLUMN(HfSingleCharm2P, hasHfSingleCharm2P, bool); //! 2-prong charm hadron (for efficiency studies) @@ -280,6 +281,7 @@ DECLARE_SOA_TABLE(HfFilters, "AOD", "HfFilters", //! filtering::HfCharmBarToXiBach, filtering::HfSigmaCPPK, filtering::HfSigmaC0K0, + filtering::HfSigmaCPr, filtering::HfPhotonCharm2P, filtering::HfPhotonCharm3P, filtering::HfSingleCharm2P,