diff --git a/PWGLF/Tasks/Resonances/xi1530Analysisqa.cxx b/PWGLF/Tasks/Resonances/xi1530Analysisqa.cxx index d4169378e9f..3ed15c1328d 100644 --- a/PWGLF/Tasks/Resonances/xi1530Analysisqa.cxx +++ b/PWGLF/Tasks/Resonances/xi1530Analysisqa.cxx @@ -20,7 +20,6 @@ #include #include #include -#include #include #include #include @@ -28,7 +27,7 @@ #include #include #include -#include +#include #include #include #include @@ -36,12 +35,20 @@ #include #include // IWYU pragma: keep (do not replace with Math/Vector4Dfwd.h) #include -#include #include -#include + +#include +#include +#include +#include #include +#include #include +#include +#include +#include +#include #include using namespace o2; @@ -52,34 +59,55 @@ using namespace o2::constants::physics; using LorentzVectorPtEtaPhiMass = ROOT::Math::PtEtaPhiMVector; using LorentzVectorSetXYZM = ROOT::Math::LorentzVector>; -Service pdgDB; - enum { kData = 0, - kLS, - kMixing, - kMCReco, - kMCTrue, - kMCTruePS, - kINEL10, - kINELg010, - kAllType, + kLS = 1, + kMixing = 2, + kMCReco = 3, + kMCTrue = 4, + kMCTruePS = 5, + kINEL10 = 6, + kINELg010 = 7, + kAllType = 8, kXiStar = 3324 }; struct Xi1530Analysisqa { + // Module-initializer tables; full tracks and cascades retain their schema. + using ResoCollisions = aod::ResoCollisions_001; + // Read original IDs only; no track_as() or original AO2D is needed here. + using ResoTracks = soa::Join; + using ResoMicroTracks = aod::ResoMicroTracks_001; + using ResoMCCols = soa::Join; + + // Own just the event rows: zero-copy slices would retain entire DF buffers. + struct MixingCollision { + float x, y, z, centrality, multiplicity; + int64_t index; + bool recINELgt0; + [[nodiscard]] float posX() const { return x; } + [[nodiscard]] float posY() const { return y; } + [[nodiscard]] float posZ() const { return z; } + [[nodiscard]] int64_t globalIndex() const { return index; } + }; + struct MixingEvent { + MixingCollision collision{}; + std::shared_ptr tracks; + }; + std::map> mixingPools; + float mixingBField = 0.f; + // Basic set-up // SliceCache cache; Preslice perRCol = aod::resodaughter::resoCollisionId; Preslice perCollision = aod::track::collisionId; - Preslice perResoCollision = + Preslice perResoCollision = aod::resodaughter::resoCollisionId; Preslice perResoCollisionCasc = aod::resodaughter::resoCollisionId; HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; - using ResoMCCols = soa::Join; // Associated with histograms struct : ConfigurableGroup { ConfigurableAxis binsPt{"binsPt", {VARIABLE_WIDTH, 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3.0, 3.1, 3.2, 3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4.0, 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5.0, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8, 5.9, 6.0, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7, 6.8, 6.9, 7.0, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 8.0, 8.1, 8.2, 8.3, 8.4, 8.5, 8.6, 8.7, 8.8, 8.9, 9.0, 9.1, 9.2, 9.3, 9.4, 9.5, 9.6, 9.7, 9.8, 9.9, 10.0, 10.1, 10.2, 10.3, 10.4, 10.5, 10.6, 10.7, 10.8, 10.9, 11.0, 11.1, 11.2, 11.3, 11.4, 11.5, 11.6, 11.7, 11.8, 11.9, 12.0, 12.1, 12.2, 12.3, 12.4, 12.5, 12.6, 12.7, 12.8, 12.9, 13.0, 13.1, 13.2, 13.3, 13.4, 13.5, 13.6, 13.7, 13.8, 13.9, 14.0, 14.1, 14.2, 14.3, 14.4, 14.5, 14.6, 14.7, 14.8, 14.9, 15.0}, "Binning of the pT axis"}; @@ -124,23 +152,24 @@ struct Xi1530Analysisqa { Configurable cfgTPCRows{"cfgTPCRows", 80, "Minimum Number of TPC Crossed Rows "}; - Configurable cfgHasTOF{"cfgHasTOF", false, "Require TOF"}; - Configurable cfgRapidityMinCut{"cfgRapidityMinCut", -0.5, "Rapidity cut for tracks"}; Configurable cfgRapidityMaxCut{"cfgRapidityMaxCut", 0.5, "Rapidity cut for tracks"}; // Primary track DCAxy to PV ConfigurableAxis cDCAtoPVBins{"cDCAtoPVBins", {1500, 0, 0.3}, "Bins for track DCA to PV"}; - Configurable cDCAxytoPVByPtPiFirstP0{"cDCAxytoPVByPtPiFirstP0", 0.004, "Coeff. Track DCAxy cut to PV by pt for Pion First (p0)"}; - Configurable cDCAxyToPVByPtPiFirstExp{"cDCAxyToPVByPtPiFirstExp", 0.013, "Coeff. Track DCAxy cut to PV by pt for Pion First (exp)"}; + Configurable cDCAxyToPVAsPt{"cDCAxyToPVAsPt", true, "Use the pT-dependent pion DCAxy cut; otherwise use cMaxDCAxyToPVCut"}; + Configurable cMaxDCAxyToPVCut{"cMaxDCAxyToPVCut", 0.5f, "Constant maximum pion |DCAxy| (cm) when its pT-dependent cut is disabled"}; + Configurable cDCAxytoPVByPtPiFirstP0{"cDCAxytoPVByPtPiFirstP0", 0.004, "Constant term of the pion pT-dependent DCAxy/z cut (cm)"}; + Configurable cDCAxyToPVByPtPiFirstExp{"cDCAxyToPVByPtPiFirstExp", 0.013, "Coefficient in pion DCAxy/z cut = P0 + coefficient / pT^power (legacy key)"}; + Configurable cDCAxyToPVByPtPiFirstPower{"cDCAxyToPVByPtPiFirstPower", 1.f, "Power in pion DCAxy/z cut = P0 + coefficient / pT^power"}; Configurable cDCAxyToPVAsPtForCasc{"cDCAxyToPVAsPtForCasc", true, "Set DCAxy to PV selection as pt cut"}; Configurable cDCAxyToPVByPtCascP0{"cDCAxyToPVByPtCascP0", 999., "Coeff. for Track DCAxy cut to PV by pt for Cascade (p0)"}; Configurable cDCAxyToPVByPtCascExp{"cDCAxyToPVByPtCascExp", 1., "Coeff. Track DCAxy cut to PV by pt for Cascade (exp)"}; // Primary track DCAz to PV - Configurable cDCAzToPVAsPt{"cDCAzToPVAsPt", true, "DCAz to PV selection as pt"}; + Configurable cDCAzToPVAsPt{"cDCAzToPVAsPt", true, "Use the shared pion pT-dependent DCA cut for DCAz; otherwise use cMaxDCAzToPVCut"}; Configurable cDCAzToPVAsPtForCasc{"cDCAzToPVAsPtForCasc", true, "Set DCA to PV selection as pt cut"}; Configurable cMaxDCAzToPVCut{"cMaxDCAzToPVCut", 0.5, "Track DCAz cut to PV Maximum"}; Configurable cMaxDCAzToPVCutCasc{"cMaxDCAzToPVCutCasc", 0.5, "Track DCAz cut to PV Maximum for casc"}; @@ -173,7 +202,7 @@ struct Xi1530Analysisqa { // Topological selections for Cascades ConfigurableAxis cDCASecondaryBins{"cDCASecondaryBins", {1000, 0, 0.1}, "Bins for DCA to Cascade secondary"}; ConfigurableAxis cProperLifetimeBins{"cProperLifetimeBins", {300, 0, 30}, "Bins for proper lifetime"}; - ConfigurableAxis cDCABachelorToPVBins{"cDCABachelorToPVcutBins", {1000, 0, 0.1}, "Bins for DCA bachelor to PV cut"}; + ConfigurableAxis cDCABachelorToPVcutBins{"cDCABachelorToPVcutBins", {1000, 0, 0.1}, "Bins for DCA bachelor to PV cut"}; Configurable cDCABachlorToPVcut{"cDCABachlorToPVcut", 0.06, "Bachelor DCA cut to PV"}; Configurable cDCAXiDaugthersCutPtRangeLower{"cDCAXiDaugthersCutPtRangeLower", 1., "Xi- DCA cut to PV as pt range lower"}; Configurable cDCAXiDaugthersCutPtRangeUpper{"cDCAXiDaugthersCutPtRangeUpper", 4., "Xi- DCA cut to PV as pt range upper"}; @@ -207,43 +236,55 @@ struct Xi1530Analysisqa { Configurable cPIDBound{"cPIDBound", 6.349, "configurable for replacing to .has"}; ConfigurableAxis cPIDnSigmaBins{"cPIDnSigmaBins", {131, -6.5, 6.5}, "Bins for nSigma PID"}; - Configurable tofAtHighPt{"tofAtHighPt", false, "Use TOF at high pT"}; - Configurable cMinTOFpt{"cMinTOFpt", 0.5, "Maximum TOF pt cut"}; + Configurable tofAtHighPt{"tofAtHighPt", false, "Apply TOF PID only above cMinTOFpt (cascade pT for its daughters)"}; + Configurable cMinTOFpt{"cMinTOFpt", 0.5, "Lower pT threshold for TOF PID when tofAtHighPt is enabled"}; // PID Selections for Pion First - Configurable cMaxtpcnSigmaPionFirst{"cMaxtpcnSigmaPionFirst", 4.0, "TPC nSigma cut for Pion First"}; - Configurable cMaxtofnSigmaPionFirst{"cMaxtofnSigmaPionFirst", 3.0, "TOF nSigma cut for Pion First"}; + Configurable cMintpcnSigmaPionFirst{"cMintpcnSigmaPionFirst", -4.0, "Strict lower TPC (nSigma - mean) bound for Pion First"}; + Configurable cMaxtpcnSigmaPionFirst{"cMaxtpcnSigmaPionFirst", 4.0, "Strict upper TPC (nSigma - mean) bound for Pion First"}; + Configurable cMeantpcnSigmaPionFirst{"cMeantpcnSigmaPionFirst", 0.0, "TPC nSigma mean subtracted for Pion First PID"}; + Configurable cMintofnSigmaPionFirst{"cMintofnSigmaPionFirst", -3.0, "Strict lower TOF (nSigma - mean) bound for Pion First"}; + Configurable cMaxtofnSigmaPionFirst{"cMaxtofnSigmaPionFirst", 3.0, "Strict upper TOF (nSigma - mean) bound for Pion First"}; + Configurable cMeantofnSigmaPionFirst{"cMeantofnSigmaPionFirst", 0.0, "TOF nSigma mean subtracted for Pion First PID"}; Configurable nsigmaCutCombinedPionFirst{"nsigmaCutCombinedPionFirst", -4.0, "Combined nSigma cut for Pion First"}; - Configurable cUseOnlyTOFTrackPionFirst{"cUseOnlyTOFTrackPionFirst", false, "Use only TOF track for PID selection Pion First"}; Configurable cByPassTOFPionFirst{"cByPassTOFPionFirst", true, "By pass TOF Pion First PID selection"}; // PID Selections for Pion Bachelor - Configurable cMaxtpcnSigmaPionBachelor{"cMaxtpcnSigmaPionBachelor", 4.0, "TPC nSigma cut for Pion Bachelor"}; - Configurable cMaxtofnSigmaPionBachelor{"cMaxtofnSigmaPionBachelor", 3.0, "TOF nSigma cut for Pion Bachelor"}; + Configurable cMintpcnSigmaPionBachelor{"cMintpcnSigmaPionBachelor", -4.0, "Strict lower TPC (nSigma - mean) bound for Pion Bachelor"}; + Configurable cMaxtpcnSigmaPionBachelor{"cMaxtpcnSigmaPionBachelor", 4.0, "Strict upper TPC (nSigma - mean) bound for Pion Bachelor"}; + Configurable cMeantpcnSigmaPionBachelor{"cMeantpcnSigmaPionBachelor", 0.0, "TPC nSigma mean subtracted for Pion Bachelor PID"}; + Configurable cMintofnSigmaPionBachelor{"cMintofnSigmaPionBachelor", -3.0, "Strict lower TOF (nSigma - mean) bound for Pion Bachelor"}; + Configurable cMaxtofnSigmaPionBachelor{"cMaxtofnSigmaPionBachelor", 3.0, "Strict upper TOF (nSigma - mean) bound for Pion Bachelor"}; + Configurable cMeantofnSigmaPionBachelor{"cMeantofnSigmaPionBachelor", 0.0, "TOF nSigma mean subtracted for Pion Bachelor PID"}; Configurable nsigmaCutCombinedPionBachelor{"nsigmaCutCombinedPionBachelor", -4.0, "Combined nSigma cut for Pion Bachelor"}; - Configurable cUseOnlyTOFTrackPionBachelor{"cUseOnlyTOFTrackPionBachelor", false, "Use only TOF track for PID selection Pion Bachelor"}; Configurable cByPassTOFPionBachelor{"cByPassTOFPionBachelor", true, "By pass TOF Pion Bachelor PID selection"}; // PID Selections for Pion - Configurable cMaxtpcnSigmaPion{"cMaxtpcnSigmaPion", 4.0, "TPC nSigma cut for Pion"}; - Configurable cMaxtofnSigmaPion{"cMaxtofnSigmaPion", 3.0, "TOF nSigma cut for Pion"}; + Configurable cMintpcnSigmaPion{"cMintpcnSigmaPion", -4.0, "Strict lower TPC (nSigma - mean) bound for Pion"}; + Configurable cMaxtpcnSigmaPion{"cMaxtpcnSigmaPion", 4.0, "Strict upper TPC (nSigma - mean) bound for Pion"}; + Configurable cMeantpcnSigmaPion{"cMeantpcnSigmaPion", 0.0, "TPC nSigma mean subtracted for Pion PID"}; + Configurable cMintofnSigmaPion{"cMintofnSigmaPion", -3.0, "Strict lower TOF (nSigma - mean) bound for Pion"}; + Configurable cMaxtofnSigmaPion{"cMaxtofnSigmaPion", 3.0, "Strict upper TOF (nSigma - mean) bound for Pion"}; + Configurable cMeantofnSigmaPion{"cMeantofnSigmaPion", 0.0, "TOF nSigma mean subtracted for Pion PID"}; Configurable nsigmaCutCombinedPion{"nsigmaCutCombinedPion", -4.0, "Combined nSigma cut for Pion"}; - Configurable cUseOnlyTOFTrackPion{"cUseOnlyTOFTrackPion", false, "Use only TOF track for PID selection Pion"}; Configurable cByPassTOFPion{"cByPassTOFPion", true, "By pass TOF Pion PID selection"}; // PID Selections for Proton - Configurable cMaxtpcnSigmaProton{"cMaxtpcnSigmaProton", 4.0, "TPC nSigma cut for Proton"}; - Configurable cMaxtofnSigmaProton{"cMaxtofnSigmaProton", 3.0, "TOF nSigma cut for Proton"}; + Configurable cMintpcnSigmaProton{"cMintpcnSigmaProton", -4.0, "Strict lower TPC (nSigma - mean) bound for Proton"}; + Configurable cMaxtpcnSigmaProton{"cMaxtpcnSigmaProton", 4.0, "Strict upper TPC (nSigma - mean) bound for Proton"}; + Configurable cMeantpcnSigmaProton{"cMeantpcnSigmaProton", 0.0, "TPC nSigma mean subtracted for Proton PID"}; + Configurable cMintofnSigmaProton{"cMintofnSigmaProton", -3.0, "Strict lower TOF (nSigma - mean) bound for Proton"}; + Configurable cMaxtofnSigmaProton{"cMaxtofnSigmaProton", 3.0, "Strict upper TOF (nSigma - mean) bound for Proton"}; + Configurable cMeantofnSigmaProton{"cMeantofnSigmaProton", 0.0, "TOF nSigma mean subtracted for Proton PID"}; Configurable nsigmaCutCombinedProton{"nsigmaCutCombinedProton", -4.0, "Combined nSigma cut for Proton"}; - Configurable cUseOnlyTOFTrackProton{"cUseOnlyTOFTrackProton", false, "Use only TOF track for PID selection Proton"}; Configurable cByPassTOFProton{"cByPassTOFProton", true, "By pass TOF Proton PID selection"}; } pidConfig; @@ -262,16 +303,12 @@ struct Xi1530Analysisqa { Configurable studyStableXi{"studyStableXi", false, "Study stable Xi"}; - Configurable cMultNTracksPVFull{"cMultNTracksPVFull", false, "Use full PV track multiplicity"}; - Configurable cMultNTracksPVeta1{"cMultNTracksPVeta1", false, "Use PV track multiplicity within |eta|<1"}; - Configurable cMultNTracksPVetaHalf{"cMultNTracksPVetaHalf", true, "Use PV track multiplicity within |eta|<0.5"}; - Configurable cRecoINELgt0{"cRecoINELgt0", true, "check if INEL>0 for reco events"}; + Configurable cMCINELgt0{"cMCINELgt0", true, "Require generator INEL>0 in reconstructed MC and generated-parent processes"}; + Configurable cMCVtxIn10{"cMCVtxIn10", true, "Require generator |vertex z| < 10 cm using isVtxIn10 in reconstructed MC and generated-parent processes"}; Configurable cUseFixedMassXi{"cUseFixedMassXi", false, "Use fixed mass for Xi-"}; - Configurable cUseTruthRapidity{"cUseTruthRapidity", false, "Use truth rapidity for Xi*"}; - - Configurable cConsiderPairOnly{"cConsiderPairOnly", true, "Consider only existing particle pairs in the event"}; + Configurable cUseTruthRapidity{"cUseTruthRapidity", true, "Use truth rapidity for Xi*"}; Configurable cfgFillRotBkg{"cfgFillRotBkg", true, "Fill rotated background"}; Configurable cfgMinRot{"cfgMinRot", 5.0 * constants::math::PI / 6.0, "Minimum of rotation"}; @@ -280,27 +317,32 @@ struct Xi1530Analysisqa { Configurable cfgNrotBkg{"cfgNrotBkg", 4, "Number of rotated copies (background) per each original candidate"}; } additionalConfig; - TRandom* rn = new TRandom(); - //*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*//*// struct PidSelectionParam { + float cMinTPCnSigma; float cMaxTPCnSigma; + float cMeanTPCnSigma; + float cMinTOFnSigma; float cMaxTOFnSigma; + float cMeanTOFnSigma; bool cByPassTOF; float nsigmaCutCombined; }; void init(o2::framework::InitContext&) { + if (doprocessMEDF && doprocessMEMicro) { + LOGF(fatal, "Enable only one mixing process: processMEDF or processMEMicro"); + } AxisSpec centAxis = {histoConfig.binsCent, "FT0M (%)"}; - AxisSpec dcaxyAxis = {primarytrackConfig.cDCAtoPVBins, "DCA_{#it{xy}} (cm)"}; - AxisSpec dcazAxis = {primarytrackConfig.cDCAtoPVBins, "DCA_{#it{z}} (cm)"}; + AxisSpec dcaxyAxis = {primarytrackConfig.cDCAtoPVBins, "|DCA_{#it{xy}}| (cm)"}; + AxisSpec dcazAxis = {primarytrackConfig.cDCAtoPVBins, "|DCA_{#it{z}}| (cm)"}; AxisSpec dcaSecondaryAxis = {cascadeConfig.cDCASecondaryBins, "DCA_{#it{Secondary}} (cm)"}; - AxisSpec dcaBachAxis = {cascadeConfig.cDCABachelorToPVBins, "DCA_{#it{Bach}} (cm)"}; + AxisSpec dcaBachAxis = {cascadeConfig.cDCABachelorToPVcutBins, "DCA_{#it{Bach}} (cm)"}; AxisSpec dcaDaugAxis = {v0sConfig.cDCADaughtersBins, "DCA_{#it{Daughter}} (cm)"}; AxisSpec cosPAAxis = {cascadeConfig.cCosPABins, "1-cos(PA)"}; - AxisSpec properLifetimeAxis = {cascadeConfig.cProperLifetimeBins, "Proper lifetime (fm/c)"}; + AxisSpec properLifetimeAxis = {cascadeConfig.cProperLifetimeBins, "c#tau (cm)"}; AxisSpec mcLabelAxis = {6, -1.5, 4.5, "MC Label"}; AxisSpec ptAxis = {histoConfig.binsPt, "#it{p}_{T} (GeV/#it{c})"}; AxisSpec ptAxisQA = {histoConfig.binsPtQA, "#it{p}_{T} (GeV/#it{c})"}; @@ -313,7 +355,10 @@ struct Xi1530Analysisqa { if (histoConfig.multQA) { // multiplicity histograms - histos.add("multQA/h2MultCent", "Multiplicity vs Centrality", HistType::kTH2F, {centAxis, histoConfig.multNTracksAxis}); + histos.add("multQA/h2MultCent", "Same-event multiplicity vs centrality (one entry per event)", HistType::kTH2F, {centAxis, histoConfig.multNTracksAxis}); + if (doprocessMEMicro || doprocessMEDF) { + histos.add("multQA/h2MultCentME", "Mixed-event anchor multiplicity vs centrality (one entry per event pair)", HistType::kTH2F, {centAxis, histoConfig.multNTracksAxis}); + } histos.add("multQA/h2MultCentMC", "Multiplicity vs Centrality MC", HistType::kTH2F, {centAxis, histoConfig.multNTracksAxis}); } @@ -403,24 +448,24 @@ struct Xi1530Analysisqa { } // 3d histogram + Flags - histos.add("h3Xi1530invmassDS", "Invariant mass of Xi- differnt sign", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3XiinvmassDS", "Invariant mass of Xi- differnt sign", kTHnSparseF, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); + histos.add("h3Xi1530invmassDS", "Invariant mass of Xi- differnt sign", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3XiinvmassDS", "Invariant mass of Xi- differnt sign", kTHnSparseD, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); - histos.add("h3Xi1530invmassLS", "Invariant mass of Xi(1530)0 same sign", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3Xi1530invmassRotDS", "Invariant mass of Xi(1530)0 rotated DS", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3Xi1530invmassLS", "Invariant mass of Xi(1530)0 same sign", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3Xi1530invmassRotDS", "Invariant mass of Xi(1530)0 rotated DS", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3Xi1530invmassDSAnti", "Invariant mass of Anti-Xi(1530)0 differnt sign", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3XiinvmassDSAnti", "Invariant mass of Anti-Xi- differnt sign", kTHnSparseF, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); + histos.add("h3Xi1530invmassDSAnti", "Invariant mass of Anti-Xi(1530)0 differnt sign", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3XiinvmassDSAnti", "Invariant mass of Anti-Xi- differnt sign", kTHnSparseD, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); - histos.add("h3Xi1530invmassLSAnti", "Invariant mass of Anti-Xi(1530)0 same sign", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3Xi1530invmassRotDSAnti", "Invariant mass of Anti-Xi(1530)0 rotated DS", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3Xi1530invmassLSAnti", "Invariant mass of Anti-Xi(1530)0 same sign", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3Xi1530invmassRotDSAnti", "Invariant mass of Anti-Xi(1530)0 rotated DS", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); - if (doprocessMEMicro) { - histos.add("h3Xi1530invmassME_DS", "Invariant mass of Xi(1530)0 mixed event DS", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); - histos.add("h3Xi1530invmassME_DSAnti", "Invariant mass of Xi(1530)0 mixed event DSAnti", kTHnSparseF, {centAxis, ptAxis, invMassAxis, flagAxis}); + if (doprocessMEMicro || doprocessMEDF) { + histos.add("h3Xi1530invmassME_DS", "Invariant mass of Xi(1530)0 mixed event DS", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); + histos.add("h3Xi1530invmassME_DSAnti", "Invariant mass of Xi(1530)0 mixed event DSAnti", kTHnSparseD, {centAxis, ptAxis, invMassAxis, flagAxis}); } - if (doprocessMC) { + if (doprocessMC || doprocessMCMicro) { // MC QA histos.add("QAMCTrue/trkDCAxy_pi", "DCAxy distribution of pion track candidates", HistType::kTH2F, {ptAxis, dcaxyAxis}); histos.add("QAMCTrue/trkDCAxy_xi", "DCAxy distribution of Xi- track candidates", HistType::kTH2F, {ptAxis, dcaxyAxis}); @@ -449,19 +494,20 @@ struct Xi1530Analysisqa { histos.add("QAMCTrue/TPC_Nsigma_piminus_all", "TPC NSigma for Pion -;#it{p}_{T} (GeV/#it{c});#sigma_{TPC}^{Pion};", {HistType::kTH3F, {centAxis, ptAxisQA, pidQAAxis}}); } - histos.add("h3RecXi1530invmass", "Invariant mass of Reconstructed MC Xi(1530)0", kTHnSparseF, {centAxis, ptAxis, invMassAxis, ptAxis}); - histos.add("h3RecXiinvmass", "Invariant mass of Reconstructed MC Xi-", kTHnSparseF, {centAxis, ptAxis, invMassAxisCasc, ptAxis}); + histos.add("h3RecXi1530invmass", "Invariant mass of Reconstructed MC Xi(1530)0", kTHnSparseD, {centAxis, ptAxis, invMassAxis, ptAxis}); + histos.add("h3RecXiinvmass", "Invariant mass of Reconstructed MC Xi-", kTHnSparseD, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); - histos.add("h3RecXi1530invmassAnti", "Invariant mass of Reconstructed MC Anti-Xi(1530)0", kTHnSparseF, {centAxis, ptAxis, invMassAxis, ptAxis}); - histos.add("h3RecXiinvmassAnti", "Invariant mass of Reconstructed MC Anti-Xi-", kTHnSparseF, {centAxis, ptAxis, invMassAxisCasc, ptAxis}); + histos.add("h3RecXi1530invmassAnti", "Invariant mass of Reconstructed MC Anti-Xi(1530)0", kTHnSparseD, {centAxis, ptAxis, invMassAxis, ptAxis}); + histos.add("h3RecXiinvmassAnti", "Invariant mass of Reconstructed MC Anti-Xi-", kTHnSparseD, {centAxis, ptAxis, invMassAxisCasc, flagAxis}); - histos.add("h3Xi1530Gen", "pT distribution of True MC Xi(1530)0", kTHnSparseF, {ptAxis, centAxis, histoConfig.multNTracksAxis}); - histos.add("h3Xi1530GenAnti", "pT distribution of True MC Anti-Xi(1530)0", kTHnSparseF, {ptAxis, centAxis, histoConfig.multNTracksAxis}); - - histos.add("Xi1530Rec", "pT distribution of Reconstructed MC Xi(1530)0", kTH2F, {ptAxis, centAxis}); - histos.add("Xi1530RecAnti", "pT distribution of Reconstructed MC Anti-Xi(1530)0", kTH2F, {ptAxis, centAxis}); + histos.add("Xi1530Rec", "pT distribution of Reconstructed MC Xi(1530)0", kTH2D, {ptAxis, centAxis}); + histos.add("Xi1530RecAnti", "pT distribution of Reconstructed MC Anti-Xi(1530)0", kTH2D, {ptAxis, centAxis}); histos.add("Xi1530Recinvmass", "Inv mass distribution of Reconstructed MC Xi(1530)0", kTH1F, {invMassAxis}); } + if (doprocessMC || doprocessMCMicro || doprocessMCTrue) { + histos.add("h3Xi1530Gen", "pT distribution of True MC Xi(1530)0", kTHnSparseD, {ptAxis, centAxis, histoConfig.multNTracksAxis}); + histos.add("h3Xi1530GenAnti", "pT distribution of True MC Anti-Xi(1530)0", kTHnSparseD, {ptAxis, centAxis, histoConfig.multNTracksAxis}); + } // QA for topological, kinematical cut for cascades if (histoConfig.additionalQAplots) { histos.add("QAbefore/V0DCATopPV", "V0s DCA to PV distribution as pt", HistType::kTH2F, {ptAxis, dcaxyAxis}); @@ -521,43 +567,42 @@ struct Xi1530Analysisqa { // Primary track selection for the first pion // template - bool primaryTrackCut(const TrackType track) + bool primaryTrackCut(const TrackType& track) { - if (std::abs(track.eta()) >= primarytrackConfig.cMaxetacut) + if (!std::isfinite(track.pt()) || !std::isfinite(track.eta()) || + !std::isfinite(track.dcaXY()) || !std::isfinite(track.dcaZ())) { + return false; + } + if (std::abs(track.eta()) >= primarytrackConfig.cMaxetacut) { return false; - if (std::abs(track.pt()) <= primarytrackConfig.cMinPtcut) + } + if (std::abs(track.pt()) <= primarytrackConfig.cMinPtcut) { return false; - if constexpr (IsResoMicrotrack) { - if (std::abs(o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(track.trackSelectionFlags())) >= (primarytrackConfig.cDCAxytoPVByPtPiFirstP0 + primarytrackConfig.cDCAxyToPVByPtPiFirstExp * std::pow(track.pt(), -1.))) + } + // Keep double precision and Power=1 equivalent to the previous pow(pt, -1.) cut. + // Evaluate the common threshold only when at least one axis needs it. + const double dcaPtCut = (primarytrackConfig.cDCAxyToPVAsPt || primarytrackConfig.cDCAzToPVAsPt) + ? primarytrackConfig.cDCAxytoPVByPtPiFirstP0.value + primarytrackConfig.cDCAxyToPVByPtPiFirstExp.value * std::pow(track.pt(), -static_cast(primarytrackConfig.cDCAxyToPVByPtPiFirstPower.value)) + : 0.; + const double dcaXYCut = primarytrackConfig.cDCAxyToPVAsPt ? dcaPtCut : primarytrackConfig.cMaxDCAxyToPVCut.value; + const double dcaZCut = primarytrackConfig.cDCAzToPVAsPt ? dcaPtCut : primarytrackConfig.cMaxDCAzToPVCut.value; + if (std::abs(track.dcaXY()) >= dcaXYCut || std::abs(track.dcaZ()) >= dcaZCut) { + return false; + } + if constexpr (!IsResoMicrotrack) { + if (track.tpcNClsFound() <= primarytrackConfig.cfgTPCcluster) { return false; - if (primarytrackConfig.cDCAzToPVAsPt) { - if (std::abs(o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(track.trackSelectionFlags())) >= (primarytrackConfig.cDCAxytoPVByPtPiFirstP0 + primarytrackConfig.cDCAxyToPVByPtPiFirstExp * std::pow(track.pt(), -1.))) - return false; - } else { - if (std::abs(o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(track.trackSelectionFlags())) >= primarytrackConfig.cMaxDCAzToPVCut) - return false; } - } else { - if (std::abs(track.dcaXY()) >= (primarytrackConfig.cDCAxytoPVByPtPiFirstP0 + primarytrackConfig.cDCAxyToPVByPtPiFirstExp * std::pow(track.pt(), -1.))) + if (track.tpcNClsCrossedRows() <= primarytrackConfig.cfgTPCRows) { return false; - if (primarytrackConfig.cDCAzToPVAsPt) { - if (std::abs(track.dcaZ()) >= (primarytrackConfig.cDCAxytoPVByPtPiFirstP0 + primarytrackConfig.cDCAxyToPVByPtPiFirstExp * std::pow(track.pt(), -1.))) - return false; - } else { - if (std::abs(track.dcaZ()) >= primarytrackConfig.cMaxDCAzToPVCut) - return false; } - if (track.tpcNClsFound() <= primarytrackConfig.cfgTPCcluster) - return false; - if (track.tpcNClsCrossedRows() <= primarytrackConfig.cfgTPCRows) - return false; } - if (primarytrackConfig.cfgHasTOF && !track.hasTOF()) - return false; - if (primarytrackConfig.cfgPVContributor && !track.isPVContributor()) + if (primarytrackConfig.cfgPVContributor && !track.isPVContributor()) { return false; - if (primarytrackConfig.cfgPrimaryTrack && !track.isPrimaryTrack()) + } + if (primarytrackConfig.cfgPrimaryTrack && !track.isPrimaryTrack()) { return false; + } return true; } @@ -568,25 +613,34 @@ struct Xi1530Analysisqa { // Primary track selection for cascades, Need to more informations for cascades // template - bool cascprimaryTrackCut(const TracksTypeCasc track) + bool cascprimaryTrackCut(const TracksTypeCasc& track) { - if (std::abs(track.eta()) >= primarytrackConfig.cMaxetacut) + if (std::abs(track.eta()) >= primarytrackConfig.cMaxetacut) { return false; - if (std::abs(track.pt()) <= primarytrackConfig.cMinPtcut) + } + if (std::abs(track.pt()) <= primarytrackConfig.cMinPtcut) { return false; - if (track.nCrossedRowsPos() <= cascadeConfig.cMinNCrossedRowsTPCPos) + } + if (track.nCrossedRowsPos() <= cascadeConfig.cMinNCrossedRowsTPCPos) { return false; - if (track.nCrossedRowsNeg() <= cascadeConfig.cMinNCrossedRowsTPCNeg) + } + if (track.nCrossedRowsNeg() <= cascadeConfig.cMinNCrossedRowsTPCNeg) { return false; - if (track.nCrossedRowsBach() <= cascadeConfig.cMinNCrossedRowsTPCBach) + } + if (track.nCrossedRowsBach() <= cascadeConfig.cMinNCrossedRowsTPCBach) { return false; + } if (primarytrackConfig.cDCAxyToPVAsPtForCasc) { - if (std::abs(track.dcaXYCascToPV()) >= (primarytrackConfig.cDCAxyToPVByPtCascP0 + primarytrackConfig.cDCAxyToPVByPtCascExp * track.pt())) + if (std::abs(track.dcaXYCascToPV()) >= (primarytrackConfig.cDCAxyToPVByPtCascP0 + primarytrackConfig.cDCAxyToPVByPtCascExp * track.pt())) { return false; + } } if (primarytrackConfig.cDCAzToPVAsPtForCasc) { - if (std::abs(track.dcaZCascToPV()) >= (primarytrackConfig.cDCAxyToPVByPtCascP0 + primarytrackConfig.cDCAxyToPVByPtCascExp * std::pow(track.pt(), -1.))) + if (std::abs(track.dcaZCascToPV()) >= (primarytrackConfig.cDCAxyToPVByPtCascP0 + primarytrackConfig.cDCAxyToPVByPtCascExp * std::pow(track.pt(), -1.))) { return false; + } + } else if (std::abs(track.dcaZCascToPV()) >= primarytrackConfig.cMaxDCAzToPVCutCasc) { + return false; } return true; @@ -596,127 +650,105 @@ struct Xi1530Analysisqa { // Topological cuts for cascades template - bool casctopCut(const TracksTypeCasc track) + bool casctopCut(const TracksTypeCasc& track) { // Topological cuts for V0s - if (std::abs(track.daughDCA()) >= v0sConfig.cDCALambdaDaugtherscut) + if (std::abs(track.daughDCA()) >= v0sConfig.cDCALambdaDaugtherscut) { return false; - if (std::abs(track.dcav0topv()) <= v0sConfig.cDCALambdaToPVcut) + } + if (std::abs(track.dcav0topv()) <= v0sConfig.cDCALambdaToPVcut) { return false; + } if (track.sign() < 0) { - if (std::abs(track.dcanegtopv()) <= v0sConfig.cDCAPionToPVcut) + if (std::abs(track.dcanegtopv()) <= v0sConfig.cDCAPionToPVcut) { return false; - if (std::abs(track.dcapostopv()) <= v0sConfig.cDCAProtonToPVcut) + } + if (std::abs(track.dcapostopv()) <= v0sConfig.cDCAProtonToPVcut) { return false; + } } else { - if (std::abs(track.dcanegtopv()) <= v0sConfig.cDCAProtonToPVcut) + if (std::abs(track.dcanegtopv()) <= v0sConfig.cDCAProtonToPVcut) { return false; - if (std::abs(track.dcapostopv()) <= v0sConfig.cDCAPionToPVcut) + } + if (std::abs(track.dcapostopv()) <= v0sConfig.cDCAPionToPVcut) { return false; + } } - if (track.v0CosPA() <= std::cos(v0sConfig.cV0CosPACutPtDepP0 - v0sConfig.cV0CosPACutPtDepP1 * track.pt())) + if (track.v0CosPA() <= std::cos(v0sConfig.cV0CosPACutPtDepP0 - v0sConfig.cV0CosPACutPtDepP1 * track.pt())) { return false; - if (track.transRadius() >= v0sConfig.cMaxV0radiuscut || track.transRadius() <= v0sConfig.cMinV0radiuscut) + } + if (track.transRadius() >= v0sConfig.cMaxV0radiuscut || track.transRadius() <= v0sConfig.cMinV0radiuscut) { return false; - if (std::abs(track.mLambda() - MassLambda) >= v0sConfig.cMasswindowV0cut) + } + if (std::abs(track.mLambda() - MassLambda) >= v0sConfig.cMasswindowV0cut) { return false; + } // Topological Cuts for Cascades - if (std::abs(track.dcabachtopv()) <= cascadeConfig.cDCABachlorToPVcut) + if (std::abs(track.dcabachtopv()) <= cascadeConfig.cDCABachlorToPVcut) { return false; + } if (track.pt() <= cascadeConfig.cDCAXiDaugthersCutPtRangeLower) { - if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepLower) + if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepLower) { return false; + } } if (track.pt() >= cascadeConfig.cDCAXiDaugthersCutPtRangeLower && track.pt() <= cascadeConfig.cDCAXiDaugthersCutPtRangeUpper) { - if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepMiddle) + if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepMiddle) { return false; + } } if (track.pt() >= cascadeConfig.cDCAXiDaugthersCutPtRangeUpper) { - if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepUpper) + if (track.cascDaughDCA() >= cascadeConfig.cDCAXiDaugthersCutPtDepUpper) { return false; + } } - if (track.cascCosPA() <= std::cos(cascadeConfig.cCosPACascCutPtDepP0 - cascadeConfig.cCosPACascCutPtDepP1 * track.pt())) + if (track.cascCosPA() <= std::cos(cascadeConfig.cCosPACascCutPtDepP0 - cascadeConfig.cCosPACascCutPtDepP1 * track.pt())) { return false; + } - if (track.cascTransRadius() >= cascadeConfig.cMaxCascradiuscut || track.cascTransRadius() <= cascadeConfig.cMinCascradiuscut) + if (track.cascTransRadius() >= cascadeConfig.cMaxCascradiuscut || track.cascTransRadius() <= cascadeConfig.cMinCascradiuscut) { return false; - if (std::abs(track.mXi() - cascadeConfig.cMassXiminus) >= cascadeConfig.cMasswindowCasccut) + } + if (std::abs(track.mXi() - cascadeConfig.cMassXiminus) >= cascadeConfig.cMasswindowCasccut) { return false; + } return true; } - bool pidSelector(float TPCNsigma, float TOFNsigma, const PidSelectionParam& params, bool tofAtHighPt, float trackPt) + bool pidSelector(float TPCNsigma, float TOFNsigma, const PidSelectionParam& params, bool tofAtHighPt, float trackPt, bool hasTOF) { - bool tpcPIDPassed{false}, tofPIDPassed{false}; - - if (tofAtHighPt && trackPt > pidConfig.cMinTOFpt) { - if (std::abs(TPCNsigma) < params.cMaxTPCnSigma) { - tpcPIDPassed = true; - } - - if (params.cByPassTOF && tpcPIDPassed) { - return true; - } - - if (hasSubsystemInfo(TOFNsigma)) { - if (std::abs(TOFNsigma) < params.cMaxTOFnSigma) { - tofPIDPassed = true; - } - if ((params.nsigmaCutCombined > 0) && - (TPCNsigma * TPCNsigma + TOFNsigma * TOFNsigma < params.nsigmaCutCombined * params.nsigmaCutCombined)) { - tofPIDPassed = true; - } - } else { - tofPIDPassed = true; - } - return tpcPIDPassed && tofPIDPassed; - } else { - - if (std::abs(TPCNsigma) < params.cMaxTPCnSigma) { - tpcPIDPassed = true; - } - - if (params.cByPassTOF && tpcPIDPassed) { - return true; - } - - if (hasSubsystemInfo(TOFNsigma)) { - if (std::abs(TOFNsigma) < params.cMaxTOFnSigma) { - tofPIDPassed = true; - } - if ((params.nsigmaCutCombined > 0) && - (TPCNsigma * TPCNsigma + TOFNsigma * TOFNsigma < params.nsigmaCutCombined * params.nsigmaCutCombined)) { - tofPIDPassed = true; - } - } else { - tofPIDPassed = true; - } - - return tpcPIDPassed && tofPIDPassed; + // Bounds are signed offsets from the detector-specific mean; exclude both edges. + const float centeredTPC = TPCNsigma - params.cMeanTPCnSigma; + const bool passesTPCWindow = params.cMinTPCnSigma < centeredTPC && centeredTPC < params.cMaxTPCnSigma; + if (!passesTPCWindow) { + return false; } + // Missing or bypassed TOF never vetoes a track that passes the TPC window. + if (params.cByPassTOF || (tofAtHighPt && trackPt <= pidConfig.cMinTOFpt) || !hasTOF) { + return true; + } + // Retain TOF-window OR combined acceptance; center both detectors for the latter. + // A present Micro001 TOF value at +/-infinity is overflow, not missing TOF. + const float centeredTOF = TOFNsigma - params.cMeanTOFnSigma; + return (params.cMinTOFnSigma < centeredTOF && centeredTOF < params.cMaxTOFnSigma) || + (params.nsigmaCutCombined > 0 && + centeredTPC * centeredTPC + centeredTOF * centeredTOF < params.nsigmaCutCombined * params.nsigmaCutCombined); } // PID selection for the First Pion // - template + template bool selectionPIDPionFirst(const T& candidate) { - float tpcNsigmaPionFirst, tofNsigmaPionFirst; + const float tpcNsigmaPionFirst = candidate.tpcNSigmaPi(); + const float tofNsigmaPionFirst = candidate.tofNSigmaPi(); float trackPt = candidate.pt(); - if constexpr (IsResoMicrotrack) { - tpcNsigmaPionFirst = o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(candidate.pidNSigmaPiFlag()); - tofNsigmaPionFirst = o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(candidate.pidNSigmaPiFlag()); - } else { - tpcNsigmaPionFirst = candidate.tpcNSigmaPi(); - tofNsigmaPionFirst = candidate.tofNSigmaPi(); - } + PidSelectionParam pionFirstParams = {.cMinTPCnSigma = pidConfig.cMintpcnSigmaPionFirst, .cMaxTPCnSigma = pidConfig.cMaxtpcnSigmaPionFirst, .cMeanTPCnSigma = pidConfig.cMeantpcnSigmaPionFirst, .cMinTOFnSigma = pidConfig.cMintofnSigmaPionFirst, .cMaxTOFnSigma = pidConfig.cMaxtofnSigmaPionFirst, .cMeanTOFnSigma = pidConfig.cMeantofnSigmaPionFirst, .cByPassTOF = pidConfig.cByPassTOFPionFirst, .nsigmaCutCombined = pidConfig.nsigmaCutCombinedPionFirst}; - PidSelectionParam pionFirstParams = {pidConfig.cMaxtpcnSigmaPionFirst, pidConfig.cMaxtofnSigmaPionFirst, pidConfig.cByPassTOFPionFirst, pidConfig.nsigmaCutCombinedPionFirst}; - - return pidSelector(tpcNsigmaPionFirst, tofNsigmaPionFirst, pionFirstParams, pidConfig.tofAtHighPt, trackPt); + return pidSelector(tpcNsigmaPionFirst, tofNsigmaPionFirst, pionFirstParams, pidConfig.tofAtHighPt, trackPt, candidate.hasTOF()); } template @@ -724,9 +756,9 @@ struct Xi1530Analysisqa { { bool lConsistentWithXi{false}, lConsistentWithLambda{false}, lConsistentWithPion{false}, lConsistentWithProton{false}; - float tpcNsigmaBachelor, tofNsigmaBachelor; - float tpcNsigmaPion, tofNsigmaPion; - float tpcNsigmaProton, tofNsigmaProton; + float tpcNsigmaBachelor = 0.f, tofNsigmaBachelor = 0.f; + float tpcNsigmaPion = 0.f, tofNsigmaPion = 0.f; + float tpcNsigmaProton = 0.f, tofNsigmaProton = 0.f; float trackPt = candidate.pt(); if (candidate.sign() < 0) { // Xi- candidates @@ -750,13 +782,13 @@ struct Xi1530Analysisqa { tofNsigmaProton = candidate.daughterTOFNSigmaNegPr(); } - PidSelectionParam bachelorParams = {pidConfig.cMaxtpcnSigmaPionBachelor, pidConfig.cMaxtofnSigmaPionBachelor, pidConfig.cByPassTOFPionBachelor, pidConfig.nsigmaCutCombinedPionBachelor}; - PidSelectionParam pionParams = {pidConfig.cMaxtpcnSigmaPion, pidConfig.cMaxtofnSigmaPion, pidConfig.cByPassTOFPion, pidConfig.nsigmaCutCombinedPion}; - PidSelectionParam protonParams = {pidConfig.cMaxtpcnSigmaProton, pidConfig.cMaxtofnSigmaProton, pidConfig.cByPassTOFProton, pidConfig.nsigmaCutCombinedProton}; + PidSelectionParam bachelorParams = {.cMinTPCnSigma = pidConfig.cMintpcnSigmaPionBachelor, .cMaxTPCnSigma = pidConfig.cMaxtpcnSigmaPionBachelor, .cMeanTPCnSigma = pidConfig.cMeantpcnSigmaPionBachelor, .cMinTOFnSigma = pidConfig.cMintofnSigmaPionBachelor, .cMaxTOFnSigma = pidConfig.cMaxtofnSigmaPionBachelor, .cMeanTOFnSigma = pidConfig.cMeantofnSigmaPionBachelor, .cByPassTOF = pidConfig.cByPassTOFPionBachelor, .nsigmaCutCombined = pidConfig.nsigmaCutCombinedPionBachelor}; + PidSelectionParam pionParams = {.cMinTPCnSigma = pidConfig.cMintpcnSigmaPion, .cMaxTPCnSigma = pidConfig.cMaxtpcnSigmaPion, .cMeanTPCnSigma = pidConfig.cMeantpcnSigmaPion, .cMinTOFnSigma = pidConfig.cMintofnSigmaPion, .cMaxTOFnSigma = pidConfig.cMaxtofnSigmaPion, .cMeanTOFnSigma = pidConfig.cMeantofnSigmaPion, .cByPassTOF = pidConfig.cByPassTOFPion, .nsigmaCutCombined = pidConfig.nsigmaCutCombinedPion}; + PidSelectionParam protonParams = {.cMinTPCnSigma = pidConfig.cMintpcnSigmaProton, .cMaxTPCnSigma = pidConfig.cMaxtpcnSigmaProton, .cMeanTPCnSigma = pidConfig.cMeantpcnSigmaProton, .cMinTOFnSigma = pidConfig.cMintofnSigmaProton, .cMaxTOFnSigma = pidConfig.cMaxtofnSigmaProton, .cMeanTOFnSigma = pidConfig.cMeantofnSigmaProton, .cByPassTOF = pidConfig.cByPassTOFProton, .nsigmaCutCombined = pidConfig.nsigmaCutCombinedProton}; - lConsistentWithXi = pidSelector(tpcNsigmaBachelor, tofNsigmaBachelor, bachelorParams, pidConfig.tofAtHighPt, trackPt); - lConsistentWithPion = pidSelector(tpcNsigmaPion, tofNsigmaPion, pionParams, pidConfig.tofAtHighPt, trackPt); - lConsistentWithProton = pidSelector(tpcNsigmaProton, tofNsigmaProton, protonParams, pidConfig.tofAtHighPt, trackPt); + lConsistentWithXi = pidSelector(tpcNsigmaBachelor, tofNsigmaBachelor, bachelorParams, pidConfig.tofAtHighPt, trackPt, hasSubsystemInfo(tofNsigmaBachelor)); + lConsistentWithPion = pidSelector(tpcNsigmaPion, tofNsigmaPion, pionParams, pidConfig.tofAtHighPt, trackPt, hasSubsystemInfo(tofNsigmaPion)); + lConsistentWithProton = pidSelector(tpcNsigmaProton, tofNsigmaProton, protonParams, pidConfig.tofAtHighPt, trackPt, hasSubsystemInfo(tofNsigmaProton)); lConsistentWithLambda = lConsistentWithProton && lConsistentWithPion; @@ -777,29 +809,27 @@ struct Xi1530Analysisqa { } template - void fillHistograms(const CollisionType& collision, const centType& inCent, const TracksType& dTracks1, const TracksTypeCasc& dTracks2) // Order: ResoColl, ResoTrack, ResoCascTrack + void fillHistograms(const CollisionType& collision, const centType& inCent, const TracksType& dTracks1, const TracksTypeCasc& dTracks2, + const MixingCollision* cascadeCollision = nullptr) // Order: ResoColl, ResoTrack, ResoCascTrack { - auto Cent = inCent; + auto cent = inCent; if (histoConfig.eventQA) { if constexpr (!IsMix) { histos.fill(HIST("QAevent/hVertexZSameE"), collision.posZ()); - histos.fill(HIST("QAevent/hMultiplicityPercentSameE"), Cent); + histos.fill(HIST("QAevent/hMultiplicityPercentSameE"), cent); histos.fill(HIST("QAevent/hCollisionIndexSameE"), collision.globalIndex()); histos.fill(HIST("QAevent/hnTrksSameE"), dTracks1.size()); histos.fill(HIST("QAevent/hnCascsSameE"), dTracks2.size()); } else { histos.fill(HIST("QAevent/hVertexZMixedE"), collision.posZ()); - histos.fill(HIST("QAevent/hMultiplicityPercentMixedE"), Cent); + histos.fill(HIST("QAevent/hMultiplicityPercentMixedE"), cent); histos.fill(HIST("QAevent/hCollisionIndexMixedE"), collision.globalIndex()); histos.fill(HIST("QAevent/hnTrksMixedE"), dTracks1.size()); histos.fill(HIST("QAevent/hnCascsMixedE"), dTracks2.size()); } } - if (additionalConfig.cConsiderPairOnly && (dTracks2.size() < 1 || dTracks1.size() < 1)) - return; - LorentzVectorPtEtaPhiMass lDecayDaughter1, lDecayDaughter2, lResonance, lDaughterRot, lResonanceRot; std::vector pionCandateIndicies = {}; std::vector xiCandateIndicies = {}; @@ -808,55 +838,44 @@ struct Xi1530Analysisqa { for (const auto& trk1 : dTracks1) { auto trk1ptPi = trk1.pt(); - float trk1DCAXY = -1.f; - float trk1DCAZ = -1.f; - float trk1NSigmaPiTPC = -999.f; - float trk1NSigmaPiTOF = -999.f; - if constexpr (IsResoMicrotrack) { - trk1DCAXY = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(trk1.trackSelectionFlags()); - trk1DCAZ = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(trk1.trackSelectionFlags()); - trk1NSigmaPiTPC = o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(trk1.pidNSigmaPiFlag()); - trk1NSigmaPiTOF = o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(trk1.pidNSigmaPiFlag()); - } else { - trk1DCAXY = trk1.dcaXY(); - trk1DCAZ = trk1.dcaZ(); - trk1NSigmaPiTPC = trk1.tpcNSigmaPi(); - trk1NSigmaPiTOF = trk1.tofNSigmaPi(); - } + const float trk1DCAXY = trk1.dcaXY(); + const float trk1DCAZ = trk1.dcaZ(); + const float trk1NSigmaPiTPC = trk1.tpcNSigmaPi(); + const float trk1NSigmaPiTOF = trk1.tofNSigmaPi(); // QA before if constexpr (!IsMix) { if (histoConfig.pidPlots) { - histos.fill(HIST("QAbefore/TPC_Nsigma_pi_first_all"), Cent, trk1ptPi, trk1NSigmaPiTPC); - if (hasSubsystemInfo(trk1NSigmaPiTOF)) { - histos.fill(HIST("QAbefore/TOF_Nsigma_pi_first_all"), Cent, trk1ptPi, trk1NSigmaPiTOF); + histos.fill(HIST("QAbefore/TPC_Nsigma_pi_first_all"), cent, trk1ptPi, trk1NSigmaPiTPC); + if (trk1.hasTOF() && !std::isnan(trk1NSigmaPiTOF)) { + histos.fill(HIST("QAbefore/TOF_Nsigma_pi_first_all"), cent, trk1ptPi, trk1NSigmaPiTOF); histos.fill(HIST("QAbefore/TOF_TPC_Map_pi_first_all"), trk1NSigmaPiTOF, trk1NSigmaPiTPC); } } if (histoConfig.additionalQAplots) { - histos.fill(HIST("QAbefore/trkDCAxy_pi"), trk1ptPi, trk1DCAXY); - histos.fill(HIST("QAbefore/trkDCAz_pi"), trk1ptPi, trk1DCAZ); + histos.fill(HIST("QAbefore/trkDCAxy_pi"), trk1ptPi, std::abs(trk1DCAXY)); + histos.fill(HIST("QAbefore/trkDCAz_pi"), trk1ptPi, std::abs(trk1DCAZ)); } } - if (pidConfig.cUseOnlyTOFTrackPionBachelor && hasSubsystemInfo(trk1NSigmaPiTOF)) + if (!selectionPIDPionFirst(trk1)) { // PID selection for the first pion continue; - if (!selectionPIDPionFirst(trk1)) // PID selection for the first pion - continue; - if (!primaryTrackCut(trk1)) // Primary track selections + } + if (!primaryTrackCut(trk1)) { // Primary track selections continue; + } if constexpr (!IsMix) { if (histoConfig.pidPlots) { - histos.fill(HIST("QAafter/TPC_Nsigma_pi_first_all"), Cent, trk1ptPi, trk1NSigmaPiTPC); + histos.fill(HIST("QAafter/TPC_Nsigma_pi_first_all"), cent, trk1ptPi, trk1NSigmaPiTPC); - if (hasSubsystemInfo(trk1NSigmaPiTOF)) { - histos.fill(HIST("QAafter/TOF_Nsigma_pi_first_all"), Cent, trk1ptPi, trk1NSigmaPiTOF); + if (trk1.hasTOF() && !std::isnan(trk1NSigmaPiTOF)) { + histos.fill(HIST("QAafter/TOF_Nsigma_pi_first_all"), cent, trk1ptPi, trk1NSigmaPiTOF); histos.fill(HIST("QAafter/TOF_TPC_Map_pi_first_all"), trk1NSigmaPiTOF, trk1NSigmaPiTPC); } } if (histoConfig.additionalQAplots) { - histos.fill(HIST("QAafter/trkDCAxy_pi"), trk1ptPi, trk1DCAXY); - histos.fill(HIST("QAafter/trkDCAz_pi"), trk1ptPi, trk1DCAZ); + histos.fill(HIST("QAafter/trkDCAxy_pi"), trk1ptPi, std::abs(trk1DCAXY)); + histos.fill(HIST("QAafter/trkDCAz_pi"), trk1ptPi, std::abs(trk1DCAZ)); } } @@ -868,7 +887,7 @@ struct Xi1530Analysisqa { auto trk2ptXi = trk2.pt(); auto massLambdaCand = trk2.mLambda(); auto massXiCand = trk2.mXi(); - auto trk2ProperLifetime = properLifetime(collision, trk2); + auto trk2ProperLifetime = cascadeCollision ? properLifetime(*cascadeCollision, trk2) : properLifetime(collision, trk2); // Topological variables for cascades auto trk2DCAV0TopPV = trk2.dcav0topv(); @@ -909,27 +928,27 @@ struct Xi1530Analysisqa { //// QA plots before the selection // // --- PID QA if (histoConfig.pidPlots) { - histos.fill(HIST("QAbefore/TPC_Nsigma_pi_bachelor_all"), Cent, 0, trk2NSigmaPiBachelorTPC); // can't take pt information for the cascade secondary + histos.fill(HIST("QAbefore/TPC_Nsigma_pi_bachelor_all"), cent, 0, trk2NSigmaPiBachelorTPC); // can't take pt information for the cascade secondary if (hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) { histos.fill(HIST("QAbefore/TOF_TPC_Map_pi_bachelor_all"), trk2NSigmaPiBachelorTOF, trk2NSigmaPiBachelorTPC); } - histos.fill(HIST("QAbefore/TPC_Nsigma_pr_all"), Cent, 0, trk2NSigmaPrPosTPC); + histos.fill(HIST("QAbefore/TPC_Nsigma_pr_all"), cent, 0, trk2NSigmaPrPosTPC); if (hasSubsystemInfo(trk2NSigmaPrPosTOF)) { histos.fill(HIST("QAbefore/TOF_TPC_Map_pr_all"), trk2NSigmaPrPosTOF, trk2NSigmaPrPosTPC); } - histos.fill(HIST("QAbefore/TPC_Nsigma_antipr_all"), Cent, 0, trk2NSigmaPrNegTPC); + histos.fill(HIST("QAbefore/TPC_Nsigma_antipr_all"), cent, 0, trk2NSigmaPrNegTPC); if (hasSubsystemInfo(trk2NSigmaPrNegTOF)) { histos.fill(HIST("QAbefore/TOF_TPC_Map_antipr_all"), trk2NSigmaPrNegTOF, trk2NSigmaPrNegTPC); } - histos.fill(HIST("QAbefore/TPC_Nsigma_pi_all"), Cent, 0, trk2NSigmaPiPosTPC); + histos.fill(HIST("QAbefore/TPC_Nsigma_pi_all"), cent, 0, trk2NSigmaPiPosTPC); if (hasSubsystemInfo(trk2NSigmaPiPosTOF)) { histos.fill(HIST("QAbefore/TOF_TPC_Map_pi_all"), trk2NSigmaPiPosTOF, trk2NSigmaPiPosTPC); } - histos.fill(HIST("QAbefore/TPC_Nsigma_piminus_all"), Cent, 0, trk2NSigmaPiNegTPC); + histos.fill(HIST("QAbefore/TPC_Nsigma_piminus_all"), cent, 0, trk2NSigmaPiNegTPC); if (hasSubsystemInfo(trk2NSigmaPiNegTOF)) { histos.fill(HIST("QAbefore/TOF_TPC_Map_piminus_all"), trk2NSigmaPiNegTOF, trk2NSigmaPiNegTPC); } @@ -937,8 +956,8 @@ struct Xi1530Analysisqa { if (histoConfig.additionalQAplots) { histos.fill(HIST("QAbefore/V0DCATopPV"), trk2ptXi, trk2DCAV0TopPV); - histos.fill(HIST("QAbefore/trkDCAxy_xi"), trk2ptXi, trk2DCAXY); - histos.fill(HIST("QAbefore/trkDCAz_xi"), trk2ptXi, trk2DCAZ); + histos.fill(HIST("QAbefore/trkDCAxy_xi"), trk2ptXi, std::abs(trk2DCAXY)); + histos.fill(HIST("QAbefore/trkDCAz_xi"), trk2ptXi, std::abs(trk2DCAZ)); histos.fill(HIST("QAbefore/V0DCADoughter"), trk2ptXi, trk2DCAV0sDougthers); histos.fill(HIST("QAbefore/CascDCADoughter"), trk2ptXi, trk2DCACascDougthers); histos.fill(HIST("QAbefore/CascDCABachPV"), trk2ptXi, trk2DCABachPV); @@ -957,22 +976,15 @@ struct Xi1530Analysisqa { } } - if (pidConfig.cUseOnlyTOFTrackPionBachelor && !hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) - continue; - if (pidConfig.cUseOnlyTOFTrackProton && !hasSubsystemInfo(trk2NSigmaPrPosTOF)) - continue; - if (pidConfig.cUseOnlyTOFTrackProton && !hasSubsystemInfo(trk2NSigmaPrNegTOF)) - continue; - if (pidConfig.cUseOnlyTOFTrackPion && !hasSubsystemInfo(trk2NSigmaPiPosTOF)) - continue; - if (pidConfig.cUseOnlyTOFTrackPion && !hasSubsystemInfo(trk2NSigmaPiNegTOF)) + if (!selectionPIDCascades(trk2)) { continue; - if (!selectionPIDCascades(trk2)) - continue; - if (!cascprimaryTrackCut(trk2) || !casctopCut(trk2)) // Primary track selections + } + if (!cascprimaryTrackCut(trk2) || !casctopCut(trk2)) { // Primary track selections continue; - if (trk2ProperLifetime >= cascadeConfig.cMaxProperLifetimeCut) + } + if (trk2ProperLifetime >= cascadeConfig.cMaxProperLifetimeCut) { continue; + } // QA after selections if constexpr (!IsMix) { @@ -981,34 +993,34 @@ struct Xi1530Analysisqa { if (histoConfig.pidPlots) { if (trk2.sign() < 0) { - histos.fill(HIST("QAafter/TPC_Nsigma_pi_bachelor_all"), Cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. + histos.fill(HIST("QAafter/TPC_Nsigma_pi_bachelor_all"), cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. if (hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_pi_bachelor_all"), trk2NSigmaPiBachelorTOF, trk2NSigmaPiBachelorTPC); } - histos.fill(HIST("QAafter/TPC_Nsigma_pr_all"), Cent, 0, trk2NSigmaPrPosTPC); + histos.fill(HIST("QAafter/TPC_Nsigma_pr_all"), cent, 0, trk2NSigmaPrPosTPC); if (hasSubsystemInfo(trk2NSigmaPrPosTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_pr_all"), trk2NSigmaPrPosTOF, trk2NSigmaPrPosTPC); } - histos.fill(HIST("QAafter/TPC_Nsigma_piminus_all"), Cent, 0, trk2NSigmaPiNegTPC); + histos.fill(HIST("QAafter/TPC_Nsigma_piminus_all"), cent, 0, trk2NSigmaPiNegTPC); if (hasSubsystemInfo(trk2NSigmaPiNegTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_piminus_all"), trk2NSigmaPiNegTOF, trk2NSigmaPiNegTPC); } } else { - histos.fill(HIST("QAafter/TPC_Nsigma_pi_bachelor_all"), Cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. + histos.fill(HIST("QAafter/TPC_Nsigma_pi_bachelor_all"), cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. if (hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_pi_bachelor_all"), trk2NSigmaPiBachelorTOF, trk2NSigmaPiBachelorTPC); } - histos.fill(HIST("QAafter/TPC_Nsigma_antipr_all"), Cent, 0, trk2NSigmaPrNegTPC); + histos.fill(HIST("QAafter/TPC_Nsigma_antipr_all"), cent, 0, trk2NSigmaPrNegTPC); if (hasSubsystemInfo(trk2NSigmaPrNegTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_antipr_all"), trk2NSigmaPrNegTOF, trk2NSigmaPrNegTPC); } - histos.fill(HIST("QAafter/TPC_Nsigma_pi_all"), Cent, 0, trk2NSigmaPiPosTPC); + histos.fill(HIST("QAafter/TPC_Nsigma_pi_all"), cent, 0, trk2NSigmaPiPosTPC); if (hasSubsystemInfo(trk2NSigmaPiPosTOF)) { histos.fill(HIST("QAafter/TOF_TPC_Map_pi_all"), trk2NSigmaPiPosTOF, trk2NSigmaPiPosTPC); } @@ -1016,8 +1028,8 @@ struct Xi1530Analysisqa { } if (histoConfig.additionalQAplots) { histos.fill(HIST("QAafter/V0DCATopPV"), trk2ptXi, trk2DCAV0TopPV); - histos.fill(HIST("QAafter/trkDCAxy_xi"), trk2ptXi, trk2DCAXY); - histos.fill(HIST("QAafter/trkDCAz_xi"), trk2ptXi, trk2DCAZ); + histos.fill(HIST("QAafter/trkDCAxy_xi"), trk2ptXi, std::abs(trk2DCAXY)); + histos.fill(HIST("QAafter/trkDCAz_xi"), trk2ptXi, std::abs(trk2DCAZ)); histos.fill(HIST("QAafter/V0DCADoughter"), trk2ptXi, trk2DCAV0sDougthers); histos.fill(HIST("QAafter/CascDCADoughter"), trk2ptXi, trk2DCACascDougthers); histos.fill(HIST("QAafter/CascDCABachPV"), trk2ptXi, trk2DCABachPV); @@ -1037,15 +1049,15 @@ struct Xi1530Analysisqa { if (additionalConfig.studyStableXi) { if (trk2.sign() < 0) { - histos.fill(HIST("h3XiinvmassDS"), Cent, trk2ptXi, massXiCand, kData); + histos.fill(HIST("h3XiinvmassDS"), cent, trk2ptXi, massXiCand, kData); } else if (trk2.sign() > 0) { - histos.fill(HIST("h3XiinvmassDSAnti"), Cent, trk2ptXi, massXiCand, kData); + histos.fill(HIST("h3XiinvmassDSAnti"), cent, trk2ptXi, massXiCand, kData); } if constexpr (IsMC) { - if (trk2.motherPDG() > 0) { - histos.fill(HIST("h3RecXiinvmass"), Cent, trk2ptXi, massXiCand, kMCReco); - } else { - histos.fill(HIST("h3RecXiinvmassAnti"), Cent, trk2ptXi, massXiCand, kMCReco); + if (trk2.pdgCode() == kXiMinus) { + histos.fill(HIST("h3RecXiinvmass"), cent, trk2ptXi, massXiCand, kMCReco); + } else if (trk2.pdgCode() == -kXiMinus) { + histos.fill(HIST("h3RecXiinvmassAnti"), cent, trk2ptXi, massXiCand, kMCReco); } } } @@ -1058,11 +1070,19 @@ struct Xi1530Analysisqa { auto pionCand = dTracks1.iteratorAt(trk1cand); for (const auto& trk2cand : xiCandateIndicies) { auto xiCand = dTracks2.iteratorAt(trk2cand); + if constexpr (!IsMix) { + const auto trackId = pionCand.trackId(); + const auto& daughterIds = xiCand.cascadeIndices(); + if (trackId >= 0 && (trackId == daughterIds[0] || trackId == daughterIds[1] || trackId == daughterIds[2])) { + continue; + } + } auto pionCandPt = pionCand.pt(); auto xiCandPt = xiCand.pt(); float massXiCand = xiCand.mXi(); - if (additionalConfig.cUseFixedMassXi) + if (additionalConfig.cUseFixedMassXi) { massXiCand = cascadeConfig.cMassXiminus; + } lDecayDaughter1 = LorentzVectorPtEtaPhiMass(pionCandPt, pionCand.eta(), pionCand.phi(), massPi); lDecayDaughter2 = LorentzVectorPtEtaPhiMass(xiCandPt, xiCand.eta(), xiCand.phi(), massXiCand); @@ -1071,27 +1091,35 @@ struct Xi1530Analysisqa { auto lResonanceMass = lResonance.M(); auto lResonancePt = lResonance.Pt(); - if ((lResonance.Rapidity() <= primarytrackConfig.cfgRapidityMinCut) || (lResonance.Rapidity() >= primarytrackConfig.cfgRapidityMaxCut)) + if ((lResonance.Rapidity() <= primarytrackConfig.cfgRapidityMinCut) || (lResonance.Rapidity() >= primarytrackConfig.cfgRapidityMaxCut)) { continue; + } if (additionalConfig.cfgCutsOnMother) { - if (lResonancePt >= additionalConfig.cMaxPtMotherCut) // excluding candidates in overflow + if (lResonancePt >= additionalConfig.cMaxPtMotherCut) { // excluding candidates in overflow continue; - if (lResonanceMass >= additionalConfig.cMaxMinvMotherCut) // excluding candidates in overflow + } + if (lResonanceMass >= additionalConfig.cMaxMinvMotherCut) { // excluding candidates in overflow continue; + } } if (pionCand.sign() * xiCand.sign() < 0) { // Signal candidates if constexpr (!IsMix) { if (pionCand.sign() > 0) { - if (histoConfig.invMass1D) + if (histoConfig.invMass1D) { histos.fill(HIST("Xi1530invmassDS"), lResonanceMass); - histos.fill(HIST("h3Xi1530invmassDS"), Cent, lResonancePt, lResonanceMass, kData); + } + histos.fill(HIST("h3Xi1530invmassDS"), cent, lResonancePt, lResonanceMass, kData); if (additionalConfig.cfgFillRotBkg) { for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("QAevent/hRotBkg"), lRotAngle); + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (histoConfig.eventQA) { + histos.fill(HIST("QAevent/hRotBkg"), lRotAngle); + } if (additionalConfig.cfgRotPion) { lDaughterRot = lDecayDaughter1; ROOT::Math::RotationZ rot(lRotAngle); @@ -1105,18 +1133,27 @@ struct Xi1530Analysisqa { lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); lResonanceRot = lDecayDaughter1 + lDaughterRot; } - histos.fill(HIST("h3Xi1530invmassRotDS"), Cent, lResonanceRot.Pt(), lResonanceRot.M(), kData); + if (additionalConfig.cfgCutsOnMother && + (lResonanceRot.Pt() >= additionalConfig.cMaxPtMotherCut || lResonanceRot.M() >= additionalConfig.cMaxMinvMotherCut)) { + continue; + } + histos.fill(HIST("h3Xi1530invmassRotDS"), cent, lResonanceRot.Pt(), lResonanceRot.M(), kData); } } } else if (pionCand.sign() < 0) { - if (histoConfig.invMass1D) + if (histoConfig.invMass1D) { histos.fill(HIST("Xi1530invmassDSAnti"), lResonanceMass); - histos.fill(HIST("h3Xi1530invmassDSAnti"), Cent, lResonancePt, lResonanceMass, kData); + } + histos.fill(HIST("h3Xi1530invmassDSAnti"), cent, lResonancePt, lResonanceMass, kData); if (additionalConfig.cfgFillRotBkg) { for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("QAevent/hRotBkg"), lRotAngle); + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (histoConfig.eventQA) { + histos.fill(HIST("QAevent/hRotBkg"), lRotAngle); + } if (additionalConfig.cfgRotPion) { lDaughterRot = lDecayDaughter1; ROOT::Math::RotationZ rot(lRotAngle); @@ -1130,46 +1167,45 @@ struct Xi1530Analysisqa { lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); lResonanceRot = lDecayDaughter1 + lDaughterRot; } - histos.fill(HIST("h3Xi1530invmassRotDSAnti"), Cent, lResonanceRot.Pt(), lResonanceRot.M(), kData); + if (additionalConfig.cfgCutsOnMother && + (lResonanceRot.Pt() >= additionalConfig.cMaxPtMotherCut || lResonanceRot.M() >= additionalConfig.cMaxMinvMotherCut)) { + continue; + } + histos.fill(HIST("h3Xi1530invmassRotDSAnti"), cent, lResonanceRot.Pt(), lResonanceRot.M(), kData); } } } } else { if (pionCand.sign() > 0) { - histos.fill(HIST("h3Xi1530invmassME_DS"), Cent, lResonancePt, lResonanceMass, kData); + histos.fill(HIST("h3Xi1530invmassME_DS"), cent, lResonancePt, lResonanceMass, kData); } else if (pionCand.sign() < 0) { - histos.fill(HIST("h3Xi1530invmassME_DSAnti"), Cent, lResonancePt, lResonanceMass, kData); + histos.fill(HIST("h3Xi1530invmassME_DSAnti"), cent, lResonancePt, lResonanceMass, kData); } } if constexpr (IsMC) { - if (std::abs(xiCand.motherPDG()) != kXiStar) + if (std::abs(xiCand.motherPDG()) != kXiStar || pionCand.motherPDG() != xiCand.motherPDG()) { continue; - if (std::abs(pionCand.pdgCode()) != kPiPlus || std::abs(xiCand.pdgCode()) != kXiMinus) + } + const int motherSign = xiCand.motherPDG() > 0 ? 1 : -1; + if (pionCand.pdgCode() != motherSign * kPiPlus || xiCand.pdgCode() != motherSign * kXiMinus) { continue; - if (pionCand.motherId() != xiCand.motherId()) + } + if (xiCand.motherId() < 0 || pionCand.motherId() != xiCand.motherId()) { continue; + } auto lResonancePtMC = xiCand.motherPt(); - if (additionalConfig.cUseTruthRapidity) - continue; - if ((xiCand.motherRap() >= primarytrackConfig.cfgRapidityMaxCut) || (xiCand.motherRap() <= primarytrackConfig.cfgRapidityMinCut)) - continue; - if (histoConfig.truthQA) { - float trk1DCAXY = -1.f; - float trk1DCAZ = -1.f; - float trk1NSigmaPiTPC = -999.f; - float trk1NSigmaPiTOF = -999.f; - if constexpr (IsResoMicrotrack) { - trk1DCAXY = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(pionCand.trackSelectionFlags()); - trk1DCAZ = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(pionCand.trackSelectionFlags()); - trk1NSigmaPiTPC = o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(pionCand.pidNSigmaPiFlag()); - trk1NSigmaPiTOF = o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(pionCand.pidNSigmaPiFlag()); - } else { - trk1DCAXY = pionCand.dcaXY(); - trk1DCAZ = pionCand.dcaZ(); - trk1NSigmaPiTPC = pionCand.tpcNSigmaPi(); - trk1NSigmaPiTOF = pionCand.tofNSigmaPi(); + if (additionalConfig.cUseTruthRapidity) { + const bool passesTruthRapidity = xiCand.motherRap() > primarytrackConfig.cfgRapidityMinCut && xiCand.motherRap() < primarytrackConfig.cfgRapidityMaxCut; + if (!passesTruthRapidity) { + continue; } + } + if (histoConfig.truthQA) { + const float trk1DCAXY = pionCand.dcaXY(); + const float trk1DCAZ = pionCand.dcaZ(); + const float trk1NSigmaPiTPC = pionCand.tpcNSigmaPi(); + const float trk1NSigmaPiTOF = pionCand.tofNSigmaPi(); auto trk2DCAXY = xiCand.dcaXYCascToPV(); auto trk2DCAZ = xiCand.dcaZCascToPV(); @@ -1207,80 +1243,88 @@ struct Xi1530Analysisqa { float trk2NSigmaPiPosTOF = xiCand.daughterTOFNSigmaPosPi(); float trk2NSigmaPiNegTOF = xiCand.daughterTOFNSigmaNegPi(); - histos.fill(HIST("QAMCTrue/trkDCAxy_pi"), pionCandPt, trk1DCAXY); - histos.fill(HIST("QAMCTrue/trkDCAz_pi"), pionCandPt, trk1DCAZ); - histos.fill(HIST("QAMCTrue/V0DCATopPV"), xiCandPt, trk2DCAV0TopPV); - histos.fill(HIST("QAMCTrue/trkDCAxy_xi"), xiCandPt, trk2DCAXY); - histos.fill(HIST("QAMCTrue/trkDCAz_xi"), xiCandPt, trk2DCAZ); - - histos.fill(HIST("QAMCTrue/V0DCADoughter"), xiCandPt, trk2DCAV0sDougthers); - histos.fill(HIST("QAMCTrue/CascDCADoughter"), xiCandPt, trk2DCACascDougthers); - histos.fill(HIST("QAMCTrue/CascDCABachPV"), xiCandPt, trk2DCABachPV); - histos.fill(HIST("QAMCTrue/CascDCAPosPV"), xiCandPt, trk2DCAPosPV); - histos.fill(HIST("QAMCTrue/CascDCANegPV"), xiCandPt, trk2DCANegPV); - histos.fill(HIST("QAMCTrue/V0CosPA"), xiCandPt, 1. - trk2V0CosPA); - histos.fill(HIST("QAMCTrue/CascCosPA"), xiCandPt, 1. - trk2CascCosPA); - histos.fill(HIST("QAMCTrue/V0Radius"), xiCandPt, trk2V0Radius); - histos.fill(HIST("QAMCTrue/CascRadius"), xiCandPt, trk2CascRadius); - histos.fill(HIST("QAMCTrue/V0Mass"), xiCandPt, massLambdaCand); - histos.fill(HIST("QAMCTrue/CascMass"), xiCandPt, massXiCand); - histos.fill(HIST("QAMCTrue/ProperLifetime"), xiCandPt, trk2ProperLifetime); - histos.fill(HIST("QAMCTrue/NCrossedRowsPos"), xiCandPt, trks2NCrossedRowsPos); - histos.fill(HIST("QAMCTrue/NCrossedRowsNeg"), xiCandPt, trks2NCrossedRowsNeg); - histos.fill(HIST("QAMCTrue/NCrossedRowsBach"), xiCandPt, trks2NCrossedRowsBach); - - histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_first_all"), Cent, pionCandPt, trk1NSigmaPiTPC); - if (hasSubsystemInfo(trk1NSigmaPiTOF)) { - histos.fill(HIST("QAMCTrue/TOF_Nsigma_pi_first_all"), Cent, pionCandPt, trk1NSigmaPiTOF); - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_first_all"), trk1NSigmaPiTOF, trk1NSigmaPiTPC); + if (histoConfig.additionalQAplots) { + histos.fill(HIST("QAMCTrue/trkDCAxy_pi"), pionCandPt, std::abs(trk1DCAXY)); + histos.fill(HIST("QAMCTrue/trkDCAz_pi"), pionCandPt, std::abs(trk1DCAZ)); + histos.fill(HIST("QAMCTrue/V0DCATopPV"), xiCandPt, trk2DCAV0TopPV); + histos.fill(HIST("QAMCTrue/trkDCAxy_xi"), xiCandPt, std::abs(trk2DCAXY)); + histos.fill(HIST("QAMCTrue/trkDCAz_xi"), xiCandPt, std::abs(trk2DCAZ)); + + histos.fill(HIST("QAMCTrue/V0DCADoughter"), xiCandPt, trk2DCAV0sDougthers); + histos.fill(HIST("QAMCTrue/CascDCADoughter"), xiCandPt, trk2DCACascDougthers); + histos.fill(HIST("QAMCTrue/CascDCABachPV"), xiCandPt, trk2DCABachPV); + histos.fill(HIST("QAMCTrue/CascDCAPosPV"), xiCandPt, trk2DCAPosPV); + histos.fill(HIST("QAMCTrue/CascDCANegPV"), xiCandPt, trk2DCANegPV); + histos.fill(HIST("QAMCTrue/V0CosPA"), xiCandPt, 1. - trk2V0CosPA); + histos.fill(HIST("QAMCTrue/CascCosPA"), xiCandPt, 1. - trk2CascCosPA); + histos.fill(HIST("QAMCTrue/V0Radius"), xiCandPt, trk2V0Radius); + histos.fill(HIST("QAMCTrue/CascRadius"), xiCandPt, trk2CascRadius); + histos.fill(HIST("QAMCTrue/V0Mass"), xiCandPt, massLambdaCand); + histos.fill(HIST("QAMCTrue/CascMass"), xiCandPt, xiCand.mXi()); + histos.fill(HIST("QAMCTrue/ProperLifetime"), xiCandPt, trk2ProperLifetime); + histos.fill(HIST("QAMCTrue/NCrossedRowsPos"), xiCandPt, trks2NCrossedRowsPos); + histos.fill(HIST("QAMCTrue/NCrossedRowsNeg"), xiCandPt, trks2NCrossedRowsNeg); + histos.fill(HIST("QAMCTrue/NCrossedRowsBach"), xiCandPt, trks2NCrossedRowsBach); } - histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_bachelor_all"), Cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. - if (hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) { - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_bachelor_all"), trk2NSigmaPiBachelorTOF, trk2NSigmaPiBachelorTPC); - } + if (histoConfig.pidPlots) { + histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_first_all"), cent, pionCandPt, trk1NSigmaPiTPC); + if (pionCand.hasTOF() && !std::isnan(trk1NSigmaPiTOF)) { + histos.fill(HIST("QAMCTrue/TOF_Nsigma_pi_first_all"), cent, pionCandPt, trk1NSigmaPiTOF); + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_first_all"), trk1NSigmaPiTOF, trk1NSigmaPiTPC); + } - histos.fill(HIST("QAMCTrue/TPC_Nsigma_pr_all"), Cent, 0, trk2NSigmaPrPosTPC); - if (hasSubsystemInfo(trk2NSigmaPrPosTOF)) { - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pr_all"), trk2NSigmaPrPosTOF, trk2NSigmaPrPosTPC); - } - histos.fill(HIST("QAMCTrue/TPC_Nsigma_antipr_all"), Cent, 0, trk2NSigmaPrNegTPC); - if (hasSubsystemInfo(trk2NSigmaPrNegTOF)) { - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_antipr_all"), trk2NSigmaPrNegTOF, trk2NSigmaPrNegTPC); - } + histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_bachelor_all"), cent, 0, trk2NSigmaPiBachelorTPC); // not exist pt information in resocascade yet. + if (hasSubsystemInfo(trk2NSigmaPiBachelorTOF)) { + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_bachelor_all"), trk2NSigmaPiBachelorTOF, trk2NSigmaPiBachelorTPC); + } - histos.fill(HIST("QAMCTrue/TPC_Nsigma_piminus_all"), Cent, 0, trk2NSigmaPiNegTPC); - if (hasSubsystemInfo(trk2NSigmaPiNegTOF)) { - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_piminus_all"), trk2NSigmaPiNegTOF, trk2NSigmaPiNegTPC); - } - histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_all"), Cent, 0, trk2NSigmaPiPosTPC); - if (hasSubsystemInfo(trk2NSigmaPiPosTOF)) { - histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_all"), trk2NSigmaPiPosTOF, trk2NSigmaPiPosTPC); + if (xiCand.sign() < 0) { + histos.fill(HIST("QAMCTrue/TPC_Nsigma_pr_all"), cent, 0, trk2NSigmaPrPosTPC); + if (hasSubsystemInfo(trk2NSigmaPrPosTOF)) { + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pr_all"), trk2NSigmaPrPosTOF, trk2NSigmaPrPosTPC); + } + histos.fill(HIST("QAMCTrue/TPC_Nsigma_piminus_all"), cent, 0, trk2NSigmaPiNegTPC); + if (hasSubsystemInfo(trk2NSigmaPiNegTOF)) { + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_piminus_all"), trk2NSigmaPiNegTOF, trk2NSigmaPiNegTPC); + } + } else { + histos.fill(HIST("QAMCTrue/TPC_Nsigma_antipr_all"), cent, 0, trk2NSigmaPrNegTPC); + if (hasSubsystemInfo(trk2NSigmaPrNegTOF)) { + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_antipr_all"), trk2NSigmaPrNegTOF, trk2NSigmaPrNegTPC); + } + histos.fill(HIST("QAMCTrue/TPC_Nsigma_pi_all"), cent, 0, trk2NSigmaPiPosTPC); + if (hasSubsystemInfo(trk2NSigmaPiPosTOF)) { + histos.fill(HIST("QAMCTrue/TOF_TPC_Map_pi_all"), trk2NSigmaPiPosTOF, trk2NSigmaPiPosTPC); + } + } } } // truthQA // MC histograms if (xiCand.motherPDG() > 0) { - histos.fill(HIST("Xi1530Rec"), lResonancePt, Cent); + histos.fill(HIST("Xi1530Rec"), lResonancePt, cent); histos.fill(HIST("Xi1530Recinvmass"), lResonanceMass); - histos.fill(HIST("h3RecXi1530invmass"), Cent, lResonancePt, lResonanceMass, lResonancePtMC); + histos.fill(HIST("h3RecXi1530invmass"), cent, lResonancePt, lResonanceMass, lResonancePtMC); } else { - histos.fill(HIST("Xi1530RecAnti"), lResonancePt, Cent); + histos.fill(HIST("Xi1530RecAnti"), lResonancePt, cent); histos.fill(HIST("Xi1530Recinvmass"), lResonanceMass); - histos.fill(HIST("h3RecXi1530invmassAnti"), Cent, lResonancePt, lResonanceMass, lResonancePtMC); + histos.fill(HIST("h3RecXi1530invmassAnti"), cent, lResonancePt, lResonanceMass, lResonancePtMC); } } // is MC } else { // Bkg candidates if constexpr (!IsMix) { if (pionCand.sign() < 0) { - if (histoConfig.invMass1D) + if (histoConfig.invMass1D) { histos.fill(HIST("Xi1530invmassLS"), lResonanceMass); - histos.fill(HIST("h3Xi1530invmassLS"), Cent, lResonancePt, lResonanceMass, kLS); + } + histos.fill(HIST("h3Xi1530invmassLS"), cent, lResonancePt, lResonanceMass, kLS); } else if (pionCand.sign() > 0) { - if (histoConfig.invMass1D) + if (histoConfig.invMass1D) { histos.fill(HIST("Xi1530invmassLSAnti"), lResonanceMass); - histos.fill(HIST("h3Xi1530invmassLSAnti"), Cent, lResonancePt, lResonanceMass, kLS); + } + histos.fill(HIST("h3Xi1530invmassLSAnti"), cent, lResonancePt, lResonanceMass, kLS); } } } // -> End if signal or bkg @@ -1288,22 +1332,17 @@ struct Xi1530Analysisqa { } // -> End loop over pion and xi candidates } // -> End fillHistograms - void processData(aod::ResoCollision const& resoCollision, - aod::ResoTracks const& resoTracks, + void processData(ResoCollisions::iterator const& resoCollision, + ResoTracks const& resoTracks, aod::ResoCascades const& cascTracks) { auto inCent = resoCollision.cent(); - auto multiplicity = 0.f; + // Version 001 stores the estimator chosen by cfgMultiplicityEstimator in the producer. + const auto multiplicity = resoCollision.multiplicity(); - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision + if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) { // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision return; - if (additionalConfig.cMultNTracksPVFull) - multiplicity = resoCollision.multNTracksPV(); - else if (additionalConfig.cMultNTracksPVeta1) - multiplicity = resoCollision.multNTracksPVeta1(); - else if (additionalConfig.cMultNTracksPVetaHalf) - multiplicity = resoCollision.multNTracksPVetaHalf(); - + } if (histoConfig.multQA) { histos.fill(HIST("multQA/h2MultCent"), inCent, multiplicity); } @@ -1313,74 +1352,83 @@ struct Xi1530Analysisqa { // Calculate numerator for the Acceptance x Efficiency void processMC(ResoMCCols::iterator const& resoCollision, soa::Join const& cascTracks, - soa::Join const& resoTracks) + soa::Join const& resoTracks) { auto inCent = resoCollision.cent(); - float multiplicity = 0.f; - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision - return; - if (additionalConfig.cMultNTracksPVFull) - multiplicity = resoCollision.multNTracksPV(); - else if (additionalConfig.cMultNTracksPVeta1) - multiplicity = resoCollision.multNTracksPVeta1(); - else if (additionalConfig.cMultNTracksPVetaHalf) - multiplicity = resoCollision.multNTracksPVetaHalf(); - - if (!resoCollision.isInAfterAllCuts()) // MC event selection, all cuts missing vtx cut + const auto multiplicity = resoCollision.multiplicity(); + if ((additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoCollision.isVtxIn10())) { return; + } + // The modular producer writes only selected reconstructed collisions. if (histoConfig.multQA) { histos.fill(HIST("multQA/h2MultCent"), inCent, multiplicity); } fillHistograms(resoCollision, inCent, resoTracks, cascTracks); } + void processMCMicro(ResoMCCols::iterator const& resoCollision, + soa::Join const& cascTracks, + soa::Join const& resoTracks) + { + if ((additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoCollision.isVtxIn10())) { + return; + } + if (histoConfig.multQA) { + histos.fill(HIST("multQA/h2MultCent"), resoCollision.cent(), resoCollision.multiplicity()); + } + fillHistograms(resoCollision, resoCollision.cent(), resoTracks, cascTracks); + } + // Calculate denominator for the Acceptance x Efficiency, actually it is not Trueth info... void processMCTrue(ResoMCCols::iterator const& resoCollision, - aod::ResoMCParents const& resoParents) + aod::ResoMCParents_001 const& resoParents) { auto multiplicity = resoCollision.mcMultiplicity(); auto inCent = resoCollision.cent(); - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision - return; - if (!resoCollision.isInAfterAllCuts()) + if ((additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoCollision.isVtxIn10())) { return; + } if (histoConfig.multQA) { histos.fill(HIST("multQA/h2MultCentMC"), inCent, multiplicity); } for (const auto& part : resoParents) { // loop over all pre-filtered MC particles - if (std::abs(part.pdgCode()) != kXiStar) + if (std::abs(part.pdgCode()) != kXiStar) { continue; - if ((part.y() <= primarytrackConfig.cfgRapidityMinCut) || (part.y() >= primarytrackConfig.cfgRapidityMaxCut)) + } + if ((part.y() <= primarytrackConfig.cfgRapidityMinCut) || (part.y() >= primarytrackConfig.cfgRapidityMaxCut)) { continue; + } bool pass1 = std::abs(part.daughterPDG1()) == kPiPlus || std::abs(part.daughterPDG2()) == kPiPlus; bool pass2 = std::abs(part.daughterPDG1()) == kXiMinus || std::abs(part.daughterPDG2()) == kXiMinus; - if (!pass1 || !pass2) + if (!pass1 || !pass2) { continue; + } - if (part.pdgCode() > 0) + if (part.pdgCode() > 0) { histos.fill(HIST("h3Xi1530Gen"), part.pt(), inCent, multiplicity); - else + } else { histos.fill(HIST("h3Xi1530GenAnti"), part.pt(), inCent, multiplicity); + } } } - void processDataMicro(aod::ResoCollision const& resoCollision, - aod::ResoMicroTracks const& resomicrotracks, + void processDataMicro(ResoCollisions::iterator const& resoCollision, + ResoMicroTracks const& resomicrotracks, aod::ResoCascades const& cascTracks) { - float multiplicity = 0.f; + const auto multiplicity = resoCollision.multiplicity(); auto inCent = resoCollision.cent(); - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision + if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) { // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision return; - if (additionalConfig.cMultNTracksPVFull) - multiplicity = resoCollision.multNTracksPV(); - else if (additionalConfig.cMultNTracksPVeta1) - multiplicity = resoCollision.multNTracksPVeta1(); - else if (additionalConfig.cMultNTracksPVetaHalf) - multiplicity = resoCollision.multNTracksPVetaHalf(); - + } if (histoConfig.multQA) { histos.fill(HIST("multQA/h2MultCent"), inCent, multiplicity); } @@ -1388,64 +1436,98 @@ struct Xi1530Analysisqa { } using BinningTypeVtxZT0M = ColumnBinningPolicy; - void processMEMicro(aod::ResoCollisions const& resoCollisions, - aod::ResoMicroTracks const& resomicrotracks, + void processMEMicro(ResoCollisions const& resoCollisions, + ResoMicroTracks const& resomicrotracks, aod::ResoCascades const& cascTracks) { + if (mixingConfig.nEvtMixing <= 0) { + return; + } auto tracksTuple = std::make_tuple(resomicrotracks, cascTracks); BinningTypeVtxZT0M colBinning{{mixingConfig.cfgVtxBins, mixingConfig.cfgMultBins}, true}; - Pair pairs{colBinning, mixingConfig.nEvtMixing, -1, resoCollisions, tracksTuple, &cache}; + Pair pairs{colBinning, mixingConfig.nEvtMixing, -1, resoCollisions, tracksTuple, &cache}; for (const auto& [collision1, tracks1, collision2, tracks2] : pairs) { - float multiplicity = 0.f; + const auto multiplicity = collision1.multiplicity(); auto inCent = collision1.cent(); - if (additionalConfig.cRecoINELgt0 && !collision1.isRecINELgt0()) // Check reco INELgt0 (at least one PV track in |eta| < 1) about the collision + if (additionalConfig.cRecoINELgt0 && (!collision1.isRecINELgt0() || !collision2.isRecINELgt0())) { continue; - if (additionalConfig.cMultNTracksPVFull) - multiplicity = collision1.multNTracksPV(); - else if (additionalConfig.cMultNTracksPVeta1) - multiplicity = collision1.multNTracksPVeta1(); - else if (additionalConfig.cMultNTracksPVetaHalf) - multiplicity = collision1.multNTracksPVetaHalf(); - + } + if (!std::isfinite(collision1.bMagField()) || collision1.bMagField() != collision2.bMagField()) { + continue; + } if (histoConfig.multQA) { - histos.fill(HIST("multQA/h2MultCent"), inCent, multiplicity); + histos.fill(HIST("multQA/h2MultCentME"), inCent, multiplicity); } - fillHistograms(collision1, inCent, tracks1, tracks2); + const MixingCollision cascadeCollision{.x = collision2.posX(), .y = collision2.posY(), .z = collision2.posZ(), .centrality = collision2.cent(), .multiplicity = collision2.multiplicity(), .index = collision2.globalIndex(), .recINELgt0 = collision2.isRecINELgt0()}; + fillHistograms(collision1, inCent, tracks1, tracks2, &cascadeCollision); } } - // void processMEDF(aod::ResoCollisionDFs const& resoCollisions, aod::ResoTrackDFs const& resotracks, aod::ResoCascadeDFs const& cascTracks) - // { - - /* Will be implemented once the DataFrame for cascade is ready. */ - - // auto tracksTuple = std::make_tuple(resotracks, cascTracks); - - // BinningTypeVtxZT0M colBinning{{mixingConfig.cfgVtxBins, mixingConfig.cfgMultBins}, true}; - // Pair pairs{colBinning, mixingConfig.nEvtMixing, -1, resoCollisions, tracksTuple, &cache}; - - // for (const auto& [collision1, tracks1, collision2, tracks2] : pairs) { + static std::shared_ptr copyMixingTracks(ResoMicroTracks const& tracks) + { + const auto source = tracks.asArrowTableConstrained(); + if (source->num_rows() == 0) { + return arrow::Table::MakeEmpty(source->schema()).ValueOrDie(); + } + std::vector> columns; + columns.reserve(source->num_columns()); + for (const auto& column : source->columns()) { + // Concatenate copies only the slice's rows, releasing the source DF buffers. + columns.push_back(std::make_shared(arrow::Concatenate(column->chunks()).ValueOrDie())); + } + return arrow::Table::Make(source->schema(), columns); + } - // float multiplicity = 0.f; - // auto inCent = collision1.cent(); - // if (histoConfig.multQA) { - // histos.fill(HIST("multQA/h2MultCent"), inCent, multiplicity); - // } - // fillHistograms(collision1, inCent, tracks1, tracks2); - // } - // } + // Cross-DF mixing on module v001 input; no legacy Reso*DF merger is needed. + // History is local to this task instance. Use inputs from one run/condition set. + void processMEDF(ResoCollisions::iterator const& collision, + ResoMicroTracks const& tracks, + aod::ResoCascades const& cascades) + { + if (mixingConfig.nEvtMixing <= 0 || !std::isfinite(collision.bMagField())) { + return; + } + if (mixingBField != collision.bMagField()) { + mixingPools.clear(); + mixingBField = collision.bMagField(); + } + BinningTypeVtxZT0M binning{{mixingConfig.cfgVtxBins, mixingConfig.cfgMultBins}, true}; + const int bin = binning.getBin(std::make_tuple(collision.posZ(), collision.cent())); + if (bin < 0) { + return; + } + const MixingCollision current{.x = collision.posX(), .y = collision.posY(), .z = collision.posZ(), .centrality = collision.cent(), .multiplicity = collision.multiplicity(), .index = collision.globalIndex(), .recINELgt0 = collision.isRecINELgt0()}; + auto& pool = mixingPools[bin]; + for (const auto& previous : pool) { + const auto& anchor = previous.collision; + if (additionalConfig.cRecoINELgt0 && (!anchor.recINELgt0 || !current.recINELgt0)) { + continue; + } + if (histoConfig.multQA) { + histos.fill(HIST("multQA/h2MultCentME"), anchor.centrality, anchor.multiplicity); + } + // Preserve the old pion -> new Cascade direction and the Cascade's own PV. + fillHistograms(anchor, anchor.centrality, ResoMicroTracks{previous.tracks}, cascades, ¤t); + } + // Insert after mixing: DF-local row IDs may repeat, but no event mixes with itself. + pool.push_back({current, copyMixingTracks(tracks)}); + while (pool.size() > static_cast(mixingConfig.nEvtMixing.value)) { + pool.pop_front(); + } + } PROCESS_SWITCH(Xi1530Analysisqa, processData, "Process Event for Data", false); PROCESS_SWITCH(Xi1530Analysisqa, processMC, "Process Event for MC (Reconstructed)", false); + PROCESS_SWITCH(Xi1530Analysisqa, processMCMicro, "Process Event for MC (Reconstructed MicroTrack 001)", false); PROCESS_SWITCH(Xi1530Analysisqa, processMCTrue, "Process Event for MC (Generated)", false); PROCESS_SWITCH(Xi1530Analysisqa, processDataMicro, "Process Event for Data (MicroTrack)", false); PROCESS_SWITCH(Xi1530Analysisqa, processMEMicro, "Process EventMixing (MicroTrack) ", false); - // PROCESS_SWITCH(Xi1530Analysisqa, processMEDF, "Process EventMixing (DataFrame) ", false); + PROCESS_SWITCH(Xi1530Analysisqa, processMEDF, "Process cross-DF EventMixing (MicroTrack 001)", false); }; -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +WorkflowSpec defineDataProcessing(ConfigContext const& context) { - return WorkflowSpec{adaptAnalysisTask(cfgc)}; + return WorkflowSpec{adaptAnalysisTask(context)}; } diff --git a/PWGLF/Tasks/Resonances/xi1820Analysis.cxx b/PWGLF/Tasks/Resonances/xi1820Analysis.cxx index 9ed500e286c..12b89e89024 100644 --- a/PWGLF/Tasks/Resonances/xi1820Analysis.cxx +++ b/PWGLF/Tasks/Resonances/xi1820Analysis.cxx @@ -29,6 +29,7 @@ #include #include #include +#include #include #include #include @@ -38,8 +39,18 @@ #include #include +#include +#include +#include + +#include +#include #include #include +#include +#include +#include +#include #include using namespace o2; @@ -50,12 +61,45 @@ using namespace o2::constants::physics; using LorentzVectorSetXYZM = ROOT::Math::LorentzVector>; struct Xi1820Analysis { + // Module-initializer tables; full tracks and V0s retain their schema. + using ResoCollisions = aod::ResoCollisions_001; + // Original track IDs are used as numbers; no parent AO2D dereference is needed. + using ResoTracks = soa::Join; + using ResoMicroTracks = aod::ResoMicroTracks_001; + using ResoMCCols = soa::Join; + SliceCache cache; Preslice perResoCollisionV0 = aod::resodaughter::resoCollisionId; Preslice perResoCollisionTrack = aod::resodaughter::resoCollisionId; - Preslice perResoCollisionMicroTrack = aod::resodaughter::resoCollisionId; + Preslice perResoCollisionMicroTrack = aod::resodaughter::resoCollisionId; HistogramRegistry histos{"histos", {}, OutputObjHandlingPolicy::AnalysisObject}; - using ResoMCCols = soa::Join; + + struct MixingCollision { + float x, y, z, centrality; + bool recINELgt0; + [[nodiscard]] float posX() const { return x; } + [[nodiscard]] float posY() const { return y; } + [[nodiscard]] float posZ() const { return z; } + [[nodiscard]] float cent() const { return centrality; } + }; + struct MixingEvent { + MixingCollision collision{}; + std::shared_ptr daughters; + }; + std::map> mixingPools; + float mixingBField = 0.f; + std::map> neutralMixingPools; + float neutralMixingBField = 0.f; + + template + struct SelectedLambda { + Iterator candidate; + bool isLambda = false, isAntiLambda = false; + }; + + Configurable cfgFillQA{"cfgFillQA", true, "Fill reconstructed daughter QA"}; + Configurable cfgFillEventQA{"cfgFillEventQA", true, "Fill event and rotation QA"}; + Configurable cfgFillTruthQA{"cfgFillTruthQA", true, "Fill MC truth-matched daughter QA"}; // Constants static constexpr float SmallMomentumDenominator = 1e-10f; // Small value to avoid division by zero @@ -86,13 +130,19 @@ struct Xi1820Analysis { Configurable cKaonITSNClusMin{"cKaonITSNClusMin", 2, "Minimum ITS clusters for kaon"}; // Kaon PID selections + Configurable cKaonTPCNSigmaMin{"cKaonTPCNSigmaMin", -3.5f, "Strict lower TPC (nSigma - mean) bound for kaons"}; + Configurable cKaonTPCNSigmaMean{"cKaonTPCNSigmaMean", 0.f, "TPC nSigma mean subtracted for kaons"}; Configurable cKaonTPCNSigmaMax{"cKaonTPCNSigmaMax", 3.5, "Maximum TPC NSigma for kaon (if not using pT-dependent)"}; + Configurable cKaonTOFNSigmaMin{"cKaonTOFNSigmaMin", -999.f, "Strict lower TOF (nSigma - mean) bound for kaons"}; + Configurable cKaonTOFNSigmaMean{"cKaonTOFNSigmaMean", 0.f, "TOF nSigma mean subtracted for kaons"}; + Configurable cKaonBypassTOF{"cKaonBypassTOF", false, "Bypass TOF PID even when TOF information is present"}; Configurable cKaonTOFNSigmaMax{"cKaonTOFNSigmaMax", 999., "Maximum TOF NSigma for kaon (if not using pT-dependent)"}; Configurable cKaonUsePtDepPID{"cKaonUsePtDepPID", false, "Use pT-dependent PID cuts"}; Configurable> cKaonPIDPtBins{"cKaonPIDPtBins", {0.0f, 0.5f, 0.8f, 2.0f, 999.0f}, "pT bin edges for PID cuts (N+1 values for N bins)"}; Configurable> cKaonTPCNSigmaCuts{"cKaonTPCNSigmaCuts", {3.0f, 3.0f, 2.0f, 2.0f}, "TPC NSigma cuts per pT bin (N values)"}; Configurable> cKaonTOFNSigmaCuts{"cKaonTOFNSigmaCuts", {3.0f, 3.0f, 3.0f, 3.0f}, "TOF NSigma cuts per pT bin (N values)"}; - Configurable> cKaonTOFRequired{"cKaonTOFRequired", {0, 0, 0, 0}, "Require TOF per pT bin (N values, 0=false, 1=true)"}; + Configurable> cKaonTPCNSigmaMinCuts{"cKaonTPCNSigmaMinCuts", {-3.f, -3.f, -2.f, -2.f}, "Strict lower centered TPC bounds per pT bin"}; + Configurable> cKaonTOFNSigmaMinCuts{"cKaonTOFNSigmaMinCuts", {-3.f, -3.f, -3.f, -3.f}, "Strict lower centered TOF bounds per pT bin"}; // V0 (Lambda) selections Configurable cV0MinCosPA{"cV0MinCosPA", 0.995, "V0 minimum pointing angle cosine"}; @@ -103,10 +153,14 @@ struct Xi1820Analysis { Configurable cV0RadiusMax{"cV0RadiusMax", 200.0, "V0 decay radius max"}; Configurable cV0DauPosDCAtoPVMin{"cV0DauPosDCAtoPVMin", 0.05, "V0 positive daughter DCA to PV min"}; Configurable cV0DauNegDCAtoPVMin{"cV0DauNegDCAtoPVMin", 0.05, "V0 negative daughter DCA to PV min"}; - Configurable cV0ProperLifetimeMax{"cV0ProperLifetimeMax", 30.0, "Lambda proper lifetime max (cm/c)"}; + Configurable cV0ProperLifetimeMax{"cV0ProperLifetimeMax", 30.0, "Lambda proper lifetime max (cm)"}; Configurable cV0sCrossMassRejection{"cV0sCrossMassRejection", true, "Enable K0s mass rejection for Lambda"}; Configurable cV0sCrossMassRejectionWindow{"cV0sCrossMassRejectionWindow", 0.005, "K0s mass rejection window for Lambda (GeV/c^2)"}; + Configurable cLambdaDaughterPiTPCNSigmaMin{"cLambdaDaughterPiTPCNSigmaMin", -5.f, "Strict lower centered TPC PID bound"}; + Configurable cLambdaDaughterPiTPCNSigmaMean{"cLambdaDaughterPiTPCNSigmaMean", 0.f, "TPC nSigma mean subtracted for the V0 daughter"}; Configurable cLambdaDaughterPiTPCNSigmaMax{"cLambdaDaughterPiTPCNSigmaMax", 5.0, "Maximum TPC NSigma for Lambda daughter pions"}; + Configurable cLambdaDaughterPrTPCNSigmaMin{"cLambdaDaughterPrTPCNSigmaMin", -5.f, "Strict lower centered TPC PID bound"}; + Configurable cLambdaDaughterPrTPCNSigmaMean{"cLambdaDaughterPrTPCNSigmaMean", 0.f, "TPC nSigma mean subtracted for the V0 daughter"}; Configurable cLambdaDaughterPrTPCNSigmaMax{"cLambdaDaughterPrTPCNSigmaMax", 5.0, "Maximum TPC NSigma for Lambda daughter protons"}; Configurable cLambdaPosDaughterMinCrossedRows{"cLambdaPosDaughterMinCrossedRows", 50, "Minimum TPC crossed rows for Lambda positive daughter"}; Configurable cLambdaNegDaughterMinCrossedRows{"cLambdaNegDaughterMinCrossedRows", 50, "Minimum TPC crossed rows for Lambda negative daughter"}; @@ -115,7 +169,7 @@ struct Xi1820Analysis { Configurable cK0sMinCosPA{"cK0sMinCosPA", 0.98, "K0s minimum pointing angle cosine"}; Configurable cK0sMaxDaughDCA{"cK0sMaxDaughDCA", 0.5, "K0s daughter DCA Maximum"}; Configurable cK0sMassWindow{"cK0sMassWindow", 0.025, "Mass window for K0s selection (GeV/c^2)"}; - Configurable cK0sProperLifetimeMax{"cK0sProperLifetimeMax", 20.0, "K0s proper lifetime max (cm/c)"}; + Configurable cK0sProperLifetimeMax{"cK0sProperLifetimeMax", 20.0, "K0s proper lifetime max (cm)"}; Configurable cK0sArmenterosQtMin{"cK0sArmenterosQtMin", 0.0, "K0s Armenteros qt min"}; Configurable cK0sArmenterosAlphaCoeff{"cK0sArmenterosAlphaCoeff", 0.2, "K0s Armenteros alpha max"}; Configurable cK0sDauPosDCAtoPVMin{"cK0sDauPosDCAtoPVMin", 0.05, "K0s positive daughter DCA to PV min"}; @@ -124,6 +178,8 @@ struct Xi1820Analysis { Configurable cK0sRadiusMax{"cK0sRadiusMax", 200.0, "K0s decay radius max"}; Configurable cK0sCrossMassRejection{"cK0sCrossMassRejection", true, "Enable Lambda mass rejection for K0s"}; Configurable cK0sCrossMassRejectionWindow{"cK0sCrossMassRejectionWindow", 0.01, "Lambda mass rejection window for K0s (GeV/c^2)"}; + Configurable cK0sDaughterPiTPCNSigmaMin{"cK0sDaughterPiTPCNSigmaMin", -5.f, "Strict lower centered TPC PID bound"}; + Configurable cK0sDaughterPiTPCNSigmaMean{"cK0sDaughterPiTPCNSigmaMean", 0.f, "TPC nSigma mean subtracted for the V0 daughter"}; Configurable cK0sDaughterPiTPCNSigmaMax{"cK0sDaughterPiTPCNSigmaMax", 5.0, "Maximum TPC NSigma for K0s daughter pions"}; Configurable cK0sPosDaughterMinCrossedRows{"cK0sPosDaughterMinCrossedRows", 50, "Minimum TPC crossed rows for K0s positive daughter"}; Configurable cK0sNegDaughterMinCrossedRows{"cK0sNegDaughterMinCrossedRows", 50, "Minimum TPC crossed rows for K0s negative daughter"}; @@ -135,14 +191,15 @@ struct Xi1820Analysis { // Additional QA and configurations struct : ConfigurableGroup { - Configurable cRecoINELgt0{"cRecoINELgt0", false, "Apply Reco INEL>0 event selection"}; - Configurable cConsiderPairOnly{"cConsiderPairOnly", true, "Consider only the pair of tracks for the charged K + Lambda analysis"}; - Configurable cConsiderHasV0s{"cConsiderHasV0s", true, "Consider only the pair of tracks for the K0s + Lambda analysis"}; + Configurable cRecoINELgt0{"cRecoINELgt0", true, "Apply Reco INEL>0 event selection"}; + Configurable cMCINELgt0{"cMCINELgt0", true, "Require generator INEL>0 in reconstructed MC and generated-parent processes"}; + Configurable cMCVtxIn10{"cMCVtxIn10", true, "Require generator |vertex z| < 10 cm using isVtxIn10 in reconstructed MC and generated-parent processes"}; Configurable cUseTruthRapidity{"cUseTruthRapidity", false, "Use truth rapidity for MC generated target"}; Configurable cUsePtDepDCAForKaons{"cUsePtDepDCAForKaons", true, "Use pT dependent DCA cuts for kaon tracks"}; Configurable cDCAToPVByPtFirstP0{"cDCAToPVByPtFirstP0", 0.004, "pT dependent DCA cut first parameter (cm)"}; - Configurable cDCAToPVByPtFirstExp{"cDCAToPVByPtFirstExp", 0.013, "pT dependent DCA cut second parameter (exponent)"}; + Configurable cDCAToPVByPtFirstExp{"cDCAToPVByPtFirstExp", 0.013, "Coefficient in kaon DCA cut = P0 + coefficient / pT^power (legacy key)"}; + Configurable cDCAToPVByPtFirstPower{"cDCAToPVByPtFirstPower", 1.f, "Power in kaon DCA cut = P0 + coefficient / pT^power"}; Configurable cMaxDcaToPVV0Lambda{"cMaxDcaToPVV0Lambda", 1.0, "Maximum DCA to PV for Lambda candidates (cm)"}; Configurable cMaxDcaToPVV0K0s{"cMaxDcaToPVV0K0s", 1.0, "Maximum DCA to PV for K0s candidates (cm)"}; @@ -165,6 +222,13 @@ struct Xi1820Analysis { void init(InitContext&) { + const int sameEventModes = static_cast(doprocessDataWithTracks) + static_cast(doprocessDataWithMicroTracks) + + static_cast(doprocessMCWithTracks) + static_cast(doprocessMCWithMicroTracks); + const int mixedEventModes = static_cast(doprocessMixedEventWithTracks) + static_cast(doprocessMixedEventWithMicroTracks) + static_cast(doprocessMEDF); + if (sameEventModes > 1 || mixedEventModes > 1 || (doprocessK0sLambda && doprocessMCK0sLambda) || + (doprocessK0sLambdaMixedEvent && doprocessK0sLambdaMEDF)) { + LOG(fatal) << "Enable at most one same-event mode and one mixing mode per charged/neutral channel"; + } AxisSpec centAxis = {binsCent, "V0M (%)"}; AxisSpec ptAxis = {binsPt, "#it{p}_{T} (GeV/#it{c})"}; AxisSpec ptAxisQA = {binsPtQA, "#it{p}_{T} (GeV/#it{c})"}; @@ -175,7 +239,7 @@ struct Xi1820Analysis { AxisSpec dcazAxis = {400, -0.2, 0.2, "DCA_{z} (cm)"}; AxisSpec cosPAAxis = {1000, 0.95, 1.0, "cos(PA)"}; AxisSpec radiusAxis = {200, 0, 200, "Radius (cm)"}; - AxisSpec lifetimeAxis = {200, 0, 50, "Proper lifetime (cm/c)"}; + AxisSpec lifetimeAxis = {200, 0, 50, "Proper lifetime (cm)"}; AxisSpec nsigmaAxis = {100, -5.0, 5.0, "N#sigma"}; AxisSpec armenterosAlphaAxis = {200, -1.0, 1.0, "Armenteros alpha"}; AxisSpec armenterosQtAxis = {500, 0.0, 0.5, "Armenteros qt (GeV/c)"}; @@ -183,16 +247,39 @@ struct Xi1820Analysis { AxisSpec axisDel = {400, -0.5, 7.0, "#DeltaR"}; // Event QA histograms - histos.add("Event/posZ", "Event vertex Z position", kTH1F, {{200, -20., 20., "V_{z} (cm)"}}); - histos.add("Event/centrality", "Event centrality distribution", kTH1D, {centAxis}); - histos.add("Event/posZvsCent", "Vertex Z vs Centrality", kTH2F, {{200, -20., 20., "V_{z} (cm)"}, centAxis}); - histos.add("Event/nV0s", "Number of V0s per event", kTH1F, {{200, 0., 200., "N_{V0s}"}}); - histos.add("Event/nKaons", "Number of kaons per event", kTH1F, {{200, 0., 200., "N_{kaon}"}}); - histos.add("Event/nLambdasAfterCuts", "Number of Lambdas per event after cuts", kTH1F, {{100, 0., 100., "N_{Lambda}"}}); - histos.add("Event/nKaonsAfterCuts", "Number of kaons (or K0s) per event after cuts", kTH1F, {{100, 0., 100., "N_{Kaon}"}}); - histos.add("Event/hRotBkg", "Rotated angle of rotated background", HistType::kTH1F, {axisRot}); - - if (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMCWithTracks || doprocessK0sLambda || doprocessMCK0sLambda) { + if (cfgFillEventQA) { + histos.add("Event/posZ", "Event vertex Z position", kTH1F, {{200, -20., 20., "V_{z} (cm)"}}); + histos.add("Event/centrality", "Event centrality distribution", kTH1D, {centAxis}); + histos.add("Event/posZvsCent", "Vertex Z vs Centrality", kTH2F, {{200, -20., 20., "V_{z} (cm)"}, centAxis}); + histos.add("Event/nV0s", "Number of V0s per event", kTH1F, {{200, 0., 200., "N_{V0s}"}}); + histos.add("Event/nKaons", "Number of kaons per event", kTH1F, {{200, 0., 200., "N_{kaon}"}}); + histos.add("Event/nLambdasAfterCuts", "Number of Lambdas per event after cuts", kTH1F, {{100, 0., 100., "N_{Lambda}"}}); + histos.add("Event/nKaonsAfterCuts", "Number of charged kaons per event after cuts", kTH1F, {{100, 0., 100., "N_{Kaon}"}}); + if (doprocessK0sLambda || doprocessMCK0sLambda) { + histos.add("Event/nK0s", "Number of K0s candidates per event before cuts", kTH1F, {{200, 0., 200., "N_{K^{0}_{S}}"}}); + histos.add("Event/nK0sAfterCuts", "Number of K0s per event after cuts", kTH1F, {{100, 0., 100., "N_{K^{0}_{S}}"}}); + } + histos.add("Event/hRotBkg", "Rotated angle of rotated background", HistType::kTH1F, {axisRot}); + } + + if (cfgFillEventQA && doprocessK0sLambdaMEDF) { + const int mixingDepth = std::max(0, nEvtMixing.value); + const AxisSpec mixingCountAxis{mixingDepth + 1, -0.5, mixingDepth + 0.5, "N_{events}"}; + const AxisSpec mixingVertexAxis{200, -20., 20., "Current-event vertex z (cm)"}; + const AxisSpec mixingCentAxis{110, 0., 110., "Current-event centrality (%)"}; + const AxisSpec mixingV0Axis{201, -0.5, 200.5, "Current-event N_{V0s} before candidate cuts"}; + const AxisSpec previousVertexAxis{200, -20., 20., "Previous-event vertex z (cm)"}; + const AxisSpec previousCentAxis{110, 0., 110., "Previous-event centrality (%)"}; + const AxisSpec previousV0Axis{201, -0.5, 200.5, "Previous-event N_{V0s} before candidate cuts"}; + histos.add("MixingQA/K0sLambdaMEDF/hPosZvsCent", "Current mixing events;Vertex z (cm);Centrality (%)", kTH2F, {mixingVertexAxis, mixingCentAxis}); + histos.add("MixingQA/K0sLambdaMEDF/hPoolSize", "Buffered events before mixing, including INEL-rejected slots (one entry per accepted current event)", kTH1F, {mixingCountAxis}); + histos.add("MixingQA/K0sLambdaMEDF/hNPartners", "Accepted event partners before V0 cuts (one entry per accepted current event)", kTH1F, {mixingCountAxis}); + histos.add("MixingQA/K0sLambdaMEDF/hPosZPair", "Mixed event vertices", kTH2F, {previousVertexAxis, mixingVertexAxis}); + histos.add("MixingQA/K0sLambdaMEDF/hCentPair", "Mixed event centralities", kTH2F, {previousCentAxis, mixingCentAxis}); + histos.add("MixingQA/K0sLambdaMEDF/hV0CountsPair", "Mixed event V0 counts before candidate cuts", kTH2F, {previousV0Axis, mixingV0Axis}); + } + + if (cfgFillQA && (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMCWithTracks || doprocessMCWithMicroTracks || doprocessK0sLambda || doprocessMCK0sLambda)) { // Lambda QA histograms histos.add("QAbefore/lambdaDCAtoPV", "Lambda DCA to PV before cuts", kTH2F, {ptAxisQA, dcaAxis}); histos.add("QAbefore/lambdaMass", "Lambda mass before cuts", kTH1F, {lambdaMassAxis}); @@ -233,7 +320,7 @@ struct Xi1820Analysis { histos.add("QAafter/lambdaArmenterosPodolanski", "Lambda candidate Armenteros-Podolanski after cuts", kTH3F, {armenterosAlphaAxis, armenterosQtAxis, ptAxisQA}); } - if (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMCWithTracks) { + if (cfgFillQA && (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMCWithTracks || doprocessMCWithMicroTracks)) { // Kaon QA histograms histos.add("QAbefore/kaonPt", "Kaon pT before cuts", kTH1F, {ptAxisQA}); histos.add("QAbefore/kaonEta", "Kaon eta before cuts", kTH1F, {{100, -2.0, 2.0, "#eta"}}); @@ -256,7 +343,7 @@ struct Xi1820Analysis { // Resonance histograms - 4 combinations // K+ Lambda - if (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMixedEventWithTracks || doprocessMixedEventWithMicroTracks || doprocessMCWithTracks) { + if (doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMixedEventWithTracks || doprocessMixedEventWithMicroTracks || doprocessMEDF || doprocessMCWithTracks || doprocessMCWithMicroTracks) { histos.add("xi1820/kplus_lambda/hInvMassKplusLambda", "Invariant mass of K^{+} + #Lambda", kTH1F, {invMassAxis}); histos.add("xi1820/kplus_lambda/hInvMassKplusLambda_Mix", "Mixed event Invariant mass of K^{+} + #Lambda", kTH1F, {invMassAxis}); @@ -288,27 +375,27 @@ struct Xi1820Analysis { histos.add("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda", "K^{-} + #bar{#Lambda} mass vs pT vs cent", kTH3D, {invMassAxis, ptAxis, centAxis}); histos.add("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda_Mix", "Mixed event K^{-} + #bar{#Lambda} mass vs pT vs cent", kTH3D, {invMassAxis, ptAxis, centAxis}); } else { - histos.add("xi1820/kplus_lambda/hMassPtCentDelKplusLambda", "K^{+} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kplus_lambda/hMassPtCentDelKplusLambda_Mix", "Mixed event K^{+} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kplus_lambda/hMassPtCentDelKplusLambda", "K^{+} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kplus_lambda/hMassPtCentDelKplusLambda_Mix", "Mixed event K^{+} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); // K+ Anti-Lambda - histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda", "K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Mix", "Mixed event K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Rot", "Rotated background K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda", "K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Mix", "Mixed event K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Rot", "Rotated background K^{+} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); // K- Lambda - histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda", "K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Mix", "Mixed event K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Rot", "Rotated background K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda", "K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Mix", "Mixed event K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Rot", "Rotated background K^{-} + #Lambda mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); // K- Anti-Lambda - histos.add("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda", "K^{-} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); - histos.add("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda_Mix", "Mixed event K^{-} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseF, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda", "K^{-} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); + histos.add("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda_Mix", "Mixed event K^{-} + #bar{#Lambda} mass vs pT vs cent vs #DeltaR", kTHnSparseD, {invMassAxis, ptAxis, centAxis, axisDel}); } } // MC Reco histograms for charged K + Lambda channel - if (doprocessMCWithTracks) { + if (doprocessMCWithTracks || doprocessMCWithMicroTracks) { histos.add("MC/kplus_antilambda/hMCRecoInvMassKplusAntiLambda", "Invariant mass of Xi(1820) to K^{-} + #Lambda (MC Reco)", kTH1F, {invMassAxis}); histos.add("MC/kplus_antilambda/hMCRecoMassPtCentKplusAntiLambda", "Xi(1820) mass vs pT vs cent (K^{-} + #Lambda) (MC Reco)", kTHnSparseD, {invMassAxis, ptAxis, centAxis, ptAxis}); @@ -317,7 +404,7 @@ struct Xi1820Analysis { } // K0s QA histograms - if (doprocessK0sLambda || doprocessMCK0sLambda) { + if (cfgFillQA && (doprocessK0sLambda || doprocessMCK0sLambda)) { histos.add("QAbefore/k0sDCAtoPV", "K0s DCA to PV before cuts", kTH2F, {ptAxisQA, dcaAxis}); histos.add("QAbefore/k0sMass", "K0s mass before cuts", kTH1F, {{100, 0.4, 0.6, "K^{0}_{S} mass (GeV/#it{c}^{2})"}}); histos.add("QAbefore/k0sPt", "K0s pT before cuts", kTH1F, {ptAxisQA}); @@ -348,7 +435,7 @@ struct Xi1820Analysis { } // K0s + Lambda - if (doprocessK0sLambda || doprocessK0sLambdaMixedEvent || doprocessMCK0sLambda) { + if (doprocessK0sLambda || doprocessK0sLambdaMixedEvent || doprocessK0sLambdaMEDF || doprocessMCK0sLambda) { histos.add("xi1820/k0s_lambda/hInvMassK0sLambda", "Invariant mass of Xi(1820) to K^{0}_{S} + #Lambda", kTH1F, {invMassAxis}); histos.add("xi1820/k0s_lambda/hInvMassK0sLambda_Mix", "Mixed event Invariant mass of Xi(1820) to K^{0}_{S} + #Lambda", kTH1F, {invMassAxis}); @@ -365,14 +452,14 @@ struct Xi1820Analysis { histos.add("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Mix", "Mixed event Xi(1820) mass vs pT vs cent (K^{0}_{S}#bar{#Lambda})", kTH3D, {invMassAxis, ptAxis, centAxis}); histos.add("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Rot", "Rotated background Xi(1820) mass vs pT vs cent (K^{0}_{S}#bar{#Lambda})", kTH3D, {invMassAxis, ptAxis, centAxis}); } else { - histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda", "Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); - histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Mix", "Mixed event Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); - histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Rot", "Rotated background Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda", "Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Mix", "Mixed event Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Rot", "Rotated background Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#Lambda)", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); // K0s + Anti-Lambda - histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda", "Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); - histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Mix", "Mixed event Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); - histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Rot", "Rotated background Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseF, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda", "Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Mix", "Mixed event Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); + histos.add("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Rot", "Rotated background Xi(1820) mass vs pT vs cent vs #DeltaR (K^{0}_{S}#bar{#Lambda})", kTHnSparseD, {{invMassAxis, ptAxis, centAxis, axisDel}}); } } @@ -384,6 +471,22 @@ struct Xi1820Analysis { histos.add("MC/k0s_antilambda/hMCRecoMassPtCentK0sAntiLambda", "Xi(1820) mass vs pT vs cent (K^{0}_{S}#bar{#Lambda}) (MC Reco)", kTHnSparseD, {invMassAxis, ptAxis, centAxis, ptAxis}); } + if (cfgFillTruthQA && (doprocessMCWithTracks || doprocessMCWithMicroTracks || doprocessMCK0sLambda)) { + histos.add("QAMCTrue/lambdaPt", "Truth-matched Lambda pT per signal pair", kTH1F, {ptAxisQA}); + histos.add("QAMCTrue/lambdaDaughterTPCNSigmaPi", "Truth-matched Lambda pion TPC PID", kTH2F, {ptAxisQA, nsigmaAxis}); + histos.add("QAMCTrue/lambdaDaughterTPCNSigmaPr", "Truth-matched Lambda proton TPC PID", kTH2F, {ptAxisQA, nsigmaAxis}); + if (doprocessMCWithTracks || doprocessMCWithMicroTracks) { + histos.add("QAMCTrue/kaonPt", "Truth-matched kaon pT per signal pair", kTH1F, {ptAxisQA}); + histos.add("QAMCTrue/kaonDCAxy", "Truth-matched kaon DCAxy", kTH2F, {ptAxisQA, dcaxyAxis}); + histos.add("QAMCTrue/kaonDCAz", "Truth-matched kaon DCAz", kTH2F, {ptAxisQA, dcazAxis}); + histos.add("QAMCTrue/kaonTPCNSigma", "Truth-matched kaon TPC PID", kTH2F, {ptAxisQA, nsigmaAxis}); + histos.add("QAMCTrue/kaonTOFNSigma", "Truth-matched kaon TOF PID", kTH2F, {ptAxisQA, nsigmaAxis}); + } + if (doprocessMCK0sLambda) { + histos.add("QAMCTrue/k0sPt", "Truth-matched K0s pT per signal pair", kTH1F, {ptAxisQA}); + } + } + if (doprocessMCGen) { histos.add("multQA/h2MultCentMC", "Multiplicity vs Centrality MC", HistType::kTH2D, {centAxis, additionalConfig.multNTracksAxis}); // MC truth invariant mass vs pT (2D) @@ -398,21 +501,21 @@ struct Xi1820Analysis { AxisSpec rapidityAxis = {100, -2.0, 2.0, "y"}; if (doprocessMCTruth) { - histos.add("MC/hMCTruthXi1820Pt", "MC Generated Xi(1820) pT", kTH1F, {ptAxis}); - histos.add("MC/hMCTruthXi1820PtEta", "MC Generated Xi(1820) pT vs eta", kTH2F, {ptAxis, etaAxis}); - histos.add("MC/hMCTruthXi1820Y", "MC Generated Xi(1820) rapidity", kTH1F, {rapidityAxis}); + histos.add("MC/hMCTruthXi1820Pt", "MC Generated Xi(1820) pT", kTH1D, {ptAxis}); + histos.add("MC/hMCTruthXi1820PtEta", "MC Generated Xi(1820) pT vs eta", kTH2D, {ptAxis, etaAxis}); + histos.add("MC/hMCTruthXi1820Y", "MC Generated Xi(1820) rapidity", kTH1D, {rapidityAxis}); // MC truth invariant mass (from MC particles) - histos.add("MC/hMCTruthInvMassKminusLambda", "MC Truth Inv Mass K^{-}#Lambda", kTH1F, {invMassAxis}); - histos.add("MC/hMCTruthInvMassKplusAntiLambda", "MC Truth Inv Mass K^{+}#bar{#Lambda}", kTH1F, {invMassAxis}); - histos.add("MC/hMCTruthInvMassK0sLambda", "MC Truth Inv Mass K^{0}_{S}#Lambda", kTH1F, {invMassAxis}); - histos.add("MC/hMCTruthInvMassK0sAntiLambda", "MC Truth Inv Mass K^{0}_{S}#bar{#Lambda}", kTH1F, {invMassAxis}); + histos.add("MC/hMCTruthInvMassKminusLambda", "MC Truth Inv Mass K^{-}#Lambda", kTH1D, {invMassAxis}); + histos.add("MC/hMCTruthInvMassKplusAntiLambda", "MC Truth Inv Mass K^{+}#bar{#Lambda}", kTH1D, {invMassAxis}); + histos.add("MC/hMCTruthInvMassK0sLambda", "MC Truth Inv Mass K^{0}_{S}#Lambda", kTH1D, {invMassAxis}); + histos.add("MC/hMCTruthInvMassK0sAntiLambda", "MC Truth Inv Mass K^{0}_{S}#bar{#Lambda}", kTH1D, {invMassAxis}); // MC truth invariant mass vs pT (2D) - histos.add("MC/hMCTruthMassPtKminusLambda", "MC Truth Mass vs pT K^{-}#Lambda", kTH2F, {invMassAxis, ptAxis}); - histos.add("MC/hMCTruthMassPtKplusAntiLambda", "MC Truth Mass vs pT K^{+}#bar{#Lambda}", kTH2F, {invMassAxis, ptAxis}); - histos.add("MC/hMCTruthMassPtK0sLambda", "MC Truth Mass vs pT K^{0}_{S}#Lambda", kTH2F, {invMassAxis, ptAxis}); - histos.add("MC/hMCTruthMassPtK0sAntiLambda", "MC Truth Mass vs pT K^{0}_{S}#bar{#Lambda}", kTH2F, {invMassAxis, ptAxis}); + histos.add("MC/hMCTruthMassPtKminusLambda", "MC Truth Mass vs pT K^{-}#Lambda", kTH2D, {invMassAxis, ptAxis}); + histos.add("MC/hMCTruthMassPtKplusAntiLambda", "MC Truth Mass vs pT K^{+}#bar{#Lambda}", kTH2D, {invMassAxis, ptAxis}); + histos.add("MC/hMCTruthMassPtK0sLambda", "MC Truth Mass vs pT K^{0}_{S}#Lambda", kTH2D, {invMassAxis, ptAxis}); + histos.add("MC/hMCTruthMassPtK0sAntiLambda", "MC Truth Mass vs pT K^{0}_{S}#bar{#Lambda}", kTH2D, {invMassAxis, ptAxis}); } } template @@ -443,36 +546,75 @@ struct Xi1820Analysis { return std::sqrt(dPhi * dPhi + dEta * dEta); } + bool passesRapidity(float rapidity) const + { + return -additionalConfig.cfgRapidityCut.value < rapidity && rapidity < additionalConfig.cfgRapidityCut.value; + } + + template + void fillTruthLambdaQA(const V0Type& lambda, bool isLambda) + { + if (!cfgFillTruthQA) { + return; + } + histos.fill(HIST("QAMCTrue/lambdaPt"), lambda.pt()); + histos.fill(HIST("QAMCTrue/lambdaDaughterTPCNSigmaPi"), lambda.pt(), isLambda ? lambda.daughterTPCNSigmaNegPi() : lambda.daughterTPCNSigmaPosPi()); + histos.fill(HIST("QAMCTrue/lambdaDaughterTPCNSigmaPr"), lambda.pt(), isLambda ? lambda.daughterTPCNSigmaPosPr() : lambda.daughterTPCNSigmaNegPr()); + } + + template + void fillTruthKaonQA(const TrackType& track) + { + if (!cfgFillTruthQA) { + return; + } + histos.fill(HIST("QAMCTrue/kaonPt"), track.pt()); + histos.fill(HIST("QAMCTrue/kaonDCAxy"), track.pt(), track.dcaXY()); + histos.fill(HIST("QAMCTrue/kaonDCAz"), track.pt(), track.dcaZ()); + histos.fill(HIST("QAMCTrue/kaonTPCNSigma"), track.pt(), track.tpcNSigmaKa()); + if (track.hasTOF() && !std::isnan(track.tofNSigmaKa())) { + histos.fill(HIST("QAMCTrue/kaonTOFNSigma"), track.pt(), track.tofNSigmaKa()); + } + } + // Lambda/Anti-Lambda selection template bool v0Cut(const CollisionType& collision, const V0Type& v0, bool isLambda) { // Basic kinematic cuts - if (std::abs(v0.eta()) > cMaxV0Etacut) + if (std::abs(v0.eta()) >= cMaxV0Etacut) { return false; - if (v0.pt() < cMinPtcut) + } + if (v0.pt() <= cMinPtcut) { return false; + } // DCA to PV - if (std::abs(v0.dcav0topv()) > additionalConfig.cMaxDcaToPVV0Lambda) + if (std::abs(v0.dcav0topv()) >= additionalConfig.cMaxDcaToPVV0Lambda) { return false; + } // Topological cuts - if (v0.v0CosPA() < cV0MinCosPA) + if (v0.v0CosPA() <= cV0MinCosPA) { return false; - if (v0.daughDCA() > cV0MaxDaughDCA) + } + if (v0.daughDCA() >= cV0MaxDaughDCA) { return false; + } // Daughter DCA to PV cuts - if (std::abs(v0.dcapostopv()) < cV0DauPosDCAtoPVMin) + if (std::abs(v0.dcapostopv()) <= cV0DauPosDCAtoPVMin) { return false; - if (std::abs(v0.dcanegtopv()) < cV0DauNegDCAtoPVMin) + } + if (std::abs(v0.dcanegtopv()) <= cV0DauNegDCAtoPVMin) { return false; + } // Radius cuts float radius = v0.transRadius(); - if (radius < cV0RadiusMin || radius > cV0RadiusMax) + if (radius <= cV0RadiusMin || radius >= cV0RadiusMax) { return false; + } // Proper lifetime cut float dx = v0.decayVtxX() - collision.posX(); @@ -481,41 +623,52 @@ struct Xi1820Analysis { float l = std::sqrt(dx * dx + dy * dy + dz * dz); float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; - if (properLifetime > cV0ProperLifetimeMax) + if (properLifetime >= cV0ProperLifetimeMax) { return false; + } // Mass window if (isLambda) { - if (std::abs(v0.mLambda() - MassLambda) > cV0MassWindow) + if (std::abs(v0.mLambda() - MassLambda) >= cV0MassWindow) { return false; + } // TPC Nsigma cuts for daughter tracks - if (std::abs(v0.daughterTPCNSigmaNegPi()) >= cLambdaDaughterPiTPCNSigmaMax) + if (!passesPIDWindow(v0.daughterTPCNSigmaNegPi(), cLambdaDaughterPiTPCNSigmaMin, cLambdaDaughterPiTPCNSigmaMax, cLambdaDaughterPiTPCNSigmaMean)) { return false; - if (std::abs(v0.daughterTPCNSigmaPosPr()) >= cLambdaDaughterPrTPCNSigmaMax) + } + if (!passesPIDWindow(v0.daughterTPCNSigmaPosPr(), cLambdaDaughterPrTPCNSigmaMin, cLambdaDaughterPrTPCNSigmaMax, cLambdaDaughterPrTPCNSigmaMean)) { return false; + } // Note: TOF beta cut that calibrated by primary vertex is not applied for V0 daughters in this moment } else { - if (std::abs(v0.mAntiLambda() - MassLambda) > cV0MassWindow) + if (std::abs(v0.mAntiLambda() - MassLambda) >= cV0MassWindow) { return false; - if (std::abs(v0.daughterTPCNSigmaPosPi()) >= cLambdaDaughterPiTPCNSigmaMax) + } + if (!passesPIDWindow(v0.daughterTPCNSigmaPosPi(), cLambdaDaughterPiTPCNSigmaMin, cLambdaDaughterPiTPCNSigmaMax, cLambdaDaughterPiTPCNSigmaMean)) { return false; - if (std::abs(v0.daughterTPCNSigmaNegPr()) >= cLambdaDaughterPrTPCNSigmaMax) + } + if (!passesPIDWindow(v0.daughterTPCNSigmaNegPr(), cLambdaDaughterPrTPCNSigmaMin, cLambdaDaughterPrTPCNSigmaMax, cLambdaDaughterPrTPCNSigmaMean)) { return false; + } // Note: TOF beta cut that calibrated by primary vertex is not applied for V0 daughters in this moment } if (cV0sCrossMassRejection) { - if (std::abs(v0.mK0Short() - MassK0Short) < cV0sCrossMassRejectionWindow) + if (std::abs(v0.mK0Short() - MassK0Short) <= cV0sCrossMassRejectionWindow) { return false; + } } - if (v0.nCrossedRowsPos() <= cLambdaPosDaughterMinCrossedRows) + if (v0.nCrossedRowsPos() <= cLambdaPosDaughterMinCrossedRows) { return false; - if (v0.nCrossedRowsNeg() <= cLambdaNegDaughterMinCrossedRows) + } + if (v0.nCrossedRowsNeg() <= cLambdaNegDaughterMinCrossedRows) { return false; + } - if (v0.qtarm() > cK0sArmenterosAlphaCoeff * std::fabs(v0.alpha())) + if (v0.qtarm() >= cK0sArmenterosAlphaCoeff * std::fabs(v0.alpha())) { return false; + } return true; } @@ -524,31 +677,39 @@ struct Xi1820Analysis { bool k0sCut(const CollisionType& collision, const V0Type& v0) { // Basic kinematic cuts - if (std::abs(v0.eta()) > cMaxV0Etacut) + if (std::abs(v0.eta()) >= cMaxV0Etacut) { return false; - if (v0.pt() < cMinPtcut) + } + if (v0.pt() <= cMinPtcut) { return false; + } // Topological cuts - if (v0.v0CosPA() < cK0sMinCosPA) + if (v0.v0CosPA() <= cK0sMinCosPA) { return false; - if (v0.daughDCA() > cK0sMaxDaughDCA) + } + if (v0.daughDCA() >= cK0sMaxDaughDCA) { return false; + } // Daughter DCA to PV cuts - if (std::abs(v0.dcapostopv()) < cK0sDauPosDCAtoPVMin) + if (std::abs(v0.dcapostopv()) <= cK0sDauPosDCAtoPVMin) { return false; - if (std::abs(v0.dcanegtopv()) < cK0sDauNegDCAtoPVMin) + } + if (std::abs(v0.dcanegtopv()) <= cK0sDauNegDCAtoPVMin) { return false; + } // Radius cuts float radius = v0.transRadius(); - if (radius < cK0sRadiusMin || radius > cK0sRadiusMax) + if (radius <= cK0sRadiusMin || radius >= cK0sRadiusMax) { return false; + } // DCA to PV - if (std::abs(v0.dcav0topv()) > additionalConfig.cMaxDcaToPVV0K0s) + if (std::abs(v0.dcav0topv()) >= additionalConfig.cMaxDcaToPVV0K0s) { return false; + } // Proper lifetime cut float dx = v0.decayVtxX() - collision.posX(); @@ -557,225 +718,248 @@ struct Xi1820Analysis { float l = std::sqrt(dx * dx + dy * dy + dz * dz); float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassK0Short; - if (properLifetime > cK0sProperLifetimeMax) + if (properLifetime >= cK0sProperLifetimeMax) { return false; + } // Mass window - if (std::abs(v0.mK0Short() - MassK0Short) > cK0sMassWindow) + if (std::abs(v0.mK0Short() - MassK0Short) >= cK0sMassWindow) { return false; + } // Competing V0 rejection: remove (Anti)Λ if (cK0sCrossMassRejection) { - if (std::abs(v0.mLambda() - MassLambda) < cK0sCrossMassRejectionWindow) + if (std::abs(v0.mLambda() - MassLambda) <= cK0sCrossMassRejectionWindow) { return false; - if (std::abs(v0.mAntiLambda() - MassLambda) < cK0sCrossMassRejectionWindow) + } + if (std::abs(v0.mAntiLambda() - MassLambda) <= cK0sCrossMassRejectionWindow) { return false; + } } - if (std::abs(v0.daughterTPCNSigmaPosPi()) >= cK0sDaughterPiTPCNSigmaMax) + if (!passesPIDWindow(v0.daughterTPCNSigmaPosPi(), cK0sDaughterPiTPCNSigmaMin, cK0sDaughterPiTPCNSigmaMax, cK0sDaughterPiTPCNSigmaMean)) { return false; - if (std::abs(v0.daughterTPCNSigmaNegPi()) >= cK0sDaughterPiTPCNSigmaMax) + } + if (!passesPIDWindow(v0.daughterTPCNSigmaNegPi(), cK0sDaughterPiTPCNSigmaMin, cK0sDaughterPiTPCNSigmaMax, cK0sDaughterPiTPCNSigmaMean)) { return false; + } // Note: TOF beta cut that calibrated by primary vertex is not applied for V0 daughters in this moment - if (v0.nCrossedRowsPos() <= cK0sPosDaughterMinCrossedRows) + if (v0.nCrossedRowsPos() <= cK0sPosDaughterMinCrossedRows) { return false; - if (v0.nCrossedRowsNeg() <= cK0sNegDaughterMinCrossedRows) + } + if (v0.nCrossedRowsNeg() <= cK0sNegDaughterMinCrossedRows) { return false; + } - if (v0.qtarm() < cK0sArmenterosAlphaCoeff * std::fabs(v0.alpha())) + if (v0.qtarm() <= cK0sArmenterosAlphaCoeff * std::fabs(v0.alpha())) { return false; + } return true; } - // Helper function to find pT bin index + // Strict signed window, retaining NaN rejection without negated compound comparisons. + static bool passesPIDWindow(float nSigma, float minimum, float maximum, float mean) + { + const float centered = nSigma - mean; + return minimum < centered && centered < maximum; + } + + // Preserve pT-bin membership [low, high); PID window edges themselves are strict. int getPtBinIndex(float pt) { - auto ptBins = static_cast>(cKaonPIDPtBins); - for (size_t i = 0; i < ptBins.size() - 1; i++) { - if (pt >= ptBins[i] && pt < ptBins[i + 1]) { - return i; + const auto& ptBins = cKaonPIDPtBins.value; + for (std::size_t i = 1; i < ptBins.size(); ++i) { + if (pt >= ptBins[i - 1] && pt < ptBins[i]) { + return static_cast(i - 1); } } - return -1; // should not happen if bins are properly configured + return -1; } - // Kaon PID selection - template + template bool kaonPidCut(const TrackType& track) { - float pt = track.pt(); - - if constexpr (IsResoMicrotrack) { - // For ResoMicroTracks - decode PID from flags - float tpcNSigma = o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(track.pidNSigmaKaFlag()); - float tofNSigma = track.hasTOF() ? o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(track.pidNSigmaKaFlag()) : 999.f; - - if (cKaonUsePtDepPID) { - // pT-dependent PID with binning - int ptBin = getPtBinIndex(pt); - if (ptBin < 0) - return false; // safety check - - auto tpcCuts = static_cast>(cKaonTPCNSigmaCuts); - auto tofCuts = static_cast>(cKaonTOFNSigmaCuts); - auto tofRequired = static_cast>(cKaonTOFRequired); - - // Check array sizes - if (ptBin >= static_cast(tpcCuts.size()) || - ptBin >= static_cast(tofCuts.size()) || - ptBin >= static_cast(tofRequired.size())) { - return false; // safety check - } - - // Apply TPC cut - if (std::abs(tpcNSigma) >= tpcCuts[ptBin]) - return false; - - // Apply TOF requirement and cut - if (tofRequired[ptBin] != 0) { - if (!track.hasTOF()) - return false; - if (std::abs(tofNSigma) >= tofCuts[ptBin]) - return false; - } else { - // TOF optional but apply cut if present - if (track.hasTOF() && std::abs(tofNSigma) >= tofCuts[ptBin]) - return false; - } - - return true; + const bool useTOF = !cKaonBypassTOF && track.hasTOF(); + float tpcMin = cKaonTPCNSigmaMin; + float tpcMax = cKaonTPCNSigmaMax; + float tofMin = cKaonTOFNSigmaMin; + float tofMax = cKaonTOFNSigmaMax; + if (cKaonUsePtDepPID) { + const int ptBin = getPtBinIndex(track.pt()); + if (ptBin < 0) { + return false; } - - // Standard PID - bool tpcPass = std::abs(tpcNSigma) < cKaonTPCNSigmaMax; - bool tofPass = track.hasTOF() ? std::abs(tofNSigma) < cKaonTOFNSigmaMax : true; - return tpcPass && tofPass; - } else { - // For ResoTracks - direct access - float tpcNSigma = track.tpcNSigmaKa(); - float tofNSigma = track.hasTOF() ? track.tofNSigmaKa() : 999.f; - - if (cKaonUsePtDepPID) { - // pT-dependent PID with binning - int ptBin = getPtBinIndex(pt); - if (ptBin < 0) - return false; // safety check - - auto tpcCuts = static_cast>(cKaonTPCNSigmaCuts); - auto tofCuts = static_cast>(cKaonTOFNSigmaCuts); - auto tofRequired = static_cast>(cKaonTOFRequired); - - // Check array sizes - if (ptBin >= static_cast(tpcCuts.size()) || - ptBin >= static_cast(tofCuts.size()) || - ptBin >= static_cast(tofRequired.size())) { - return false; // safety check - } - - // Apply TPC cut - if (std::abs(tpcNSigma) >= tpcCuts[ptBin]) + const auto bin = static_cast(ptBin); + if (bin >= cKaonTPCNSigmaCuts.value.size() || bin >= cKaonTPCNSigmaMinCuts.value.size()) { + return false; + } + tpcMin = cKaonTPCNSigmaMinCuts.value[bin]; + tpcMax = cKaonTPCNSigmaCuts.value[bin]; + if (useTOF) { + if (bin >= cKaonTOFNSigmaCuts.value.size() || bin >= cKaonTOFNSigmaMinCuts.value.size()) { return false; - - // Apply TOF requirement and cut - if (tofRequired[ptBin] != 0) { - if (!track.hasTOF()) - return false; - if (std::abs(tofNSigma) >= tofCuts[ptBin]) - return false; - } else { - // TOF optional but apply cut if present - if (track.hasTOF() && std::abs(tofNSigma) >= tofCuts[ptBin]) - return false; } - - return true; + tofMin = cKaonTOFNSigmaMinCuts.value[bin]; + tofMax = cKaonTOFNSigmaCuts.value[bin]; } - - // Standard PID - bool tpcPass = std::abs(tpcNSigma) < cKaonTPCNSigmaMax; - bool tofPass = track.hasTOF() ? std::abs(tofNSigma) < cKaonTOFNSigmaMax : true; - return tpcPass && tofPass; } + if (!passesPIDWindow(track.tpcNSigmaKa(), tpcMin, tpcMax, cKaonTPCNSigmaMean)) { + return false; + } + // No TOF-presence veto: missing or explicitly bypassed TOF uses TPC alone. + return !useTOF || + passesPIDWindow(track.tofNSigmaKa(), tofMin, tofMax, cKaonTOFNSigmaMean); } // Kaon track selection (for both ResoTracks and ResoMicroTracks) template bool kaonCut(const TrackType& track) { - float candPt = track.pt(); - // Basic kinematic cuts - if (candPt < cKaonPtMin) + const float candPt = track.pt(); + if (!std::isfinite(candPt) || !std::isfinite(track.eta()) || + !std::isfinite(track.dcaXY()) || !std::isfinite(track.dcaZ())) { return false; - if (std::abs(track.eta()) > cKaonEtaMax) + } + if (candPt <= cKaonPtMin || std::abs(track.eta()) >= cKaonEtaMax) { + return false; + } + const double dcaPtCut = additionalConfig.cUsePtDepDCAForKaons + ? additionalConfig.cDCAToPVByPtFirstP0.value + additionalConfig.cDCAToPVByPtFirstExp.value * std::pow(candPt, -static_cast(additionalConfig.cDCAToPVByPtFirstPower.value)) + : 0.; + const double dcaXYCut = additionalConfig.cUsePtDepDCAForKaons ? dcaPtCut : cKaonDCAxyMax.value; + const double dcaZCut = additionalConfig.cUsePtDepDCAForKaons ? dcaPtCut : cKaonDCAzMax.value; + if (std::abs(track.dcaXY()) >= dcaXYCut || std::abs(track.dcaZ()) >= dcaZCut) { return false; - - float dcaXY = -999.f; - float dcaZ = -999.f; - - // DCA cuts - different access for ResoMicroTracks - if constexpr (IsResoMicrotrack) { - dcaXY = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(track.trackSelectionFlags()); - dcaZ = o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(track.trackSelectionFlags()); - - if (additionalConfig.cUsePtDepDCAForKaons) { // Insert pT dependent DCAxy,z cut (tighter than global-track w DCA cut) - if (std::abs(dcaXY) > (additionalConfig.cDCAToPVByPtFirstP0 + additionalConfig.cDCAToPVByPtFirstExp * std::pow(candPt, -1.))) - return false; - if (std::abs(dcaZ) > (additionalConfig.cDCAToPVByPtFirstP0 + additionalConfig.cDCAToPVByPtFirstExp * std::pow(candPt, -1.))) - return false; - } else { - if (std::abs(dcaXY) > cKaonDCAxyMax) - return false; - if (std::abs(dcaZ) > cKaonDCAzMax) - return false; - } - } else { - if (additionalConfig.cUsePtDepDCAForKaons) { // Insert pT dependent DCAxy,z cut (tighter than global-track w DCA cut) - if (std::abs(track.dcaXY()) > (additionalConfig.cDCAToPVByPtFirstP0 + additionalConfig.cDCAToPVByPtFirstExp * std::pow(candPt, -1.))) - return false; - if (std::abs(track.dcaZ()) > (additionalConfig.cDCAToPVByPtFirstP0 + additionalConfig.cDCAToPVByPtFirstExp * std::pow(candPt, -1.))) - return false; - } else { - if (std::abs(track.dcaXY()) > cKaonDCAxyMax) - return false; - if (std::abs(track.dcaZ()) > cKaonDCAzMax) - return false; - } } // Track quality cuts - check if fields are available (only for ResoTracks) if constexpr (!IsResoMicrotrack) { if constexpr (requires { track.tpcNClsFound(); }) { - if (track.tpcNClsFound() < cKaonTPCNClusMin) + if (track.tpcNClsFound() <= cKaonTPCNClusMin) { return false; + } } if constexpr (requires { track.itsNCls(); }) { - if (track.itsNCls() < cKaonITSNClusMin) + if (track.itsNCls() <= cKaonITSNClusMin) { return false; + } } } // PID selection - if (!kaonPidCut(track)) + if (!kaonPidCut(track)) { return false; + } // Flag selections for Primary track selection - if (additionalConfig.cfgPVContributor && !track.isPVContributor()) + if (additionalConfig.cfgPVContributor && !track.isPVContributor()) { return false; - if (additionalConfig.cfgPrimaryTrack && !track.isPrimaryTrack()) + } + if (additionalConfig.cfgPrimaryTrack && !track.isPrimaryTrack()) { return false; + } return true; } - template - void fillChargedKLambda(const CollisionT& collision, const V0sT& v0s, const TracksT& tracks) // Xi(1820) analysis: charged K + Lambda channel + bool hasChargedSameEventProcess() + { + return doprocessDataWithTracks || doprocessDataWithMicroTracks || doprocessMCWithTracks || doprocessMCWithMicroTracks; + } + + // Select each V0 once; QA and event counts must not be weighted by kaon multiplicity. + template + auto selectLambdas(const CollisionT& collision, const V0sT& v0s, bool fillSharedQA = true) + { + std::vector> selected; + selected.reserve(v0s.size()); + for (const auto& v0 : v0s) { + // Lambda QA before cuts + if (!IsMix && cfgFillQA && fillSharedQA) { + histos.fill(HIST("QAbefore/lambdaDCAtoPV"), v0.pt(), v0.dcav0topv()); + histos.fill(HIST("QAbefore/lambdaMass"), v0.mLambda()); + histos.fill(HIST("QAbefore/lambdaMassAnti"), v0.mAntiLambda()); + histos.fill(HIST("QAbefore/lambdaPt"), v0.pt()); + histos.fill(HIST("QAbefore/lambdaEta"), v0.eta()); + histos.fill(HIST("QAbefore/lambdaCosPA"), v0.pt(), v0.v0CosPA()); + histos.fill(HIST("QAbefore/lambdaRadius"), v0.pt(), v0.transRadius()); + histos.fill(HIST("QAbefore/lambdaDauDCA"), v0.pt(), v0.daughDCA()); + histos.fill(HIST("QAbefore/lambdaDauPosDCA"), v0.pt(), std::abs(v0.dcapostopv())); + histos.fill(HIST("QAbefore/lambdaDauNegDCA"), v0.pt(), std::abs(v0.dcanegtopv())); + histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaPosPr"), v0.pt(), v0.daughterTPCNSigmaPosPr()); + histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaNegPi"), v0.pt(), v0.daughterTPCNSigmaNegPi()); + histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaPosPi"), v0.pt(), v0.daughterTPCNSigmaPosPi()); + histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaNegPr"), v0.pt(), v0.daughterTPCNSigmaNegPr()); + histos.fill(HIST("QAbefore/lambdaNCrossedRowsPos"), v0.pt(), v0.nCrossedRowsPos()); + histos.fill(HIST("QAbefore/lambdaNCrossedRowsNeg"), v0.pt(), v0.nCrossedRowsNeg()); + + // Calculate proper lifetime manually + float dx = v0.decayVtxX() - collision.posX(); + float dy = v0.decayVtxY() - collision.posY(); + float dz = v0.decayVtxZ() - collision.posZ(); + float l = std::sqrt(dx * dx + dy * dy + dz * dz); + float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); + auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; + histos.fill(HIST("QAbefore/lambdaProperLifetime"), v0.pt(), properLifetime); + histos.fill(HIST("QAbefore/lambdaArmenterosPodolanski"), v0.alpha(), v0.qtarm(), v0.pt()); + } + + // Try Lambda + bool isLambda = v0Cut(collision, v0, true); + // Try Anti-Lambda + bool isAntiLambda = v0Cut(collision, v0, false); + + if (!isLambda && !isAntiLambda) { + continue; + } + + if (!IsMix && cfgFillQA && fillSharedQA) { + // QA after cuts (fill for whichever passes) + if (isLambda) { + histos.fill(HIST("QAafter/lambdaMass"), v0.mLambda()); + histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaPosPr"), v0.pt(), v0.daughterTPCNSigmaPosPr()); + histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaNegPi"), v0.pt(), v0.daughterTPCNSigmaNegPi()); + } + if (isAntiLambda) { + histos.fill(HIST("QAafter/lambdaMassAnti"), v0.mAntiLambda()); + histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaPosPi"), v0.pt(), v0.daughterTPCNSigmaPosPi()); + histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaNegPr"), v0.pt(), v0.daughterTPCNSigmaNegPr()); + } + histos.fill(HIST("QAafter/lambdaDCAtoPV"), v0.pt(), v0.dcav0topv()); + histos.fill(HIST("QAafter/lambdaPt"), v0.pt()); + histos.fill(HIST("QAafter/lambdaEta"), v0.eta()); + histos.fill(HIST("QAafter/lambdaCosPA"), v0.pt(), v0.v0CosPA()); + histos.fill(HIST("QAafter/lambdaRadius"), v0.pt(), v0.transRadius()); + histos.fill(HIST("QAafter/lambdaDauDCA"), v0.pt(), v0.daughDCA()); + histos.fill(HIST("QAafter/lambdaDauPosDCA"), v0.pt(), std::abs(v0.dcapostopv())); + histos.fill(HIST("QAafter/lambdaDauNegDCA"), v0.pt(), std::abs(v0.dcanegtopv())); + histos.fill(HIST("QAafter/lambdaNCrossedRowsPos"), v0.pt(), v0.nCrossedRowsPos()); + histos.fill(HIST("QAafter/lambdaNCrossedRowsNeg"), v0.pt(), v0.nCrossedRowsNeg()); + + float dx = v0.decayVtxX() - collision.posX(); + float dy = v0.decayVtxY() - collision.posY(); + float dz = v0.decayVtxZ() - collision.posZ(); + float l = std::sqrt(dx * dx + dy * dy + dz * dz); + float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); + auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; + histos.fill(HIST("QAafter/lambdaProperLifetime"), v0.pt(), properLifetime); + histos.fill(HIST("QAafter/lambdaArmenterosPodolanski"), v0.alpha(), v0.qtarm(), v0.pt()); + } + + selected.push_back({v0, isLambda, isAntiLambda}); + } + return selected; + } + + template + void fillChargedKLambda(const CollisionT& collision, const V0sT& v0s, const TracksT& tracks, const LambdaCollisionT& lambdaCollision) // Xi(1820) analysis: charged K + Lambda channel { auto cent = collision.cent(); // Fill event QA histograms (only for same-event) - if constexpr (!IsMix) { + if (!IsMix && cfgFillEventQA) { histos.fill(HIST("Event/posZ"), collision.posZ()); histos.fill(HIST("Event/centrality"), cent); histos.fill(HIST("Event/posZvsCent"), collision.posZ(), cent); @@ -783,12 +967,8 @@ struct Xi1820Analysis { histos.fill(HIST("Event/nKaons"), tracks.size()); } - if (additionalConfig.cConsiderPairOnly && (v0s.size() < 1 || tracks.size() < 1)) - return; // skip events that cannot form pairs if the option is enabled (for increasing processing speed when only pairs are of interest) - - // Count candidates after cuts - int nV0sAfterCuts = 0; - int nKaonsAfterCuts = 0; + std::vector selectedKaons; + selectedKaons.reserve(tracks.size()); // Build 4 combinations ROOT::Math::PxPyPzEVector pKaon, pLambda, pRes, lDaughterRot, lResonanceRot; @@ -796,23 +976,16 @@ struct Xi1820Analysis { // Loop over kaon candidates for (const auto& kaon : tracks) { // QA before cuts - if constexpr (!IsMix) { + if (!IsMix && cfgFillQA) { histos.fill(HIST("QAbefore/kaonPt"), kaon.pt()); histos.fill(HIST("QAbefore/kaonEta"), kaon.eta()); - if constexpr (IsResoMicrotrack) { - histos.fill(HIST("QAbefore/kaonDCAxy"), kaon.pt(), o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(kaon.trackSelectionFlags())); - histos.fill(HIST("QAbefore/kaonDCAz"), kaon.pt(), o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(kaon.trackSelectionFlags())); - histos.fill(HIST("QAbefore/kaonTPCNSigma"), kaon.pt(), o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(kaon.pidNSigmaKaFlag())); - if (kaon.hasTOF()) { - histos.fill(HIST("QAbefore/kaonTOFNSigma"), kaon.pt(), o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(kaon.pidNSigmaKaFlag())); - } - } else { - histos.fill(HIST("QAbefore/kaonDCAxy"), kaon.pt(), kaon.dcaXY()); - histos.fill(HIST("QAbefore/kaonDCAz"), kaon.pt(), kaon.dcaZ()); - histos.fill(HIST("QAbefore/kaonTPCNSigma"), kaon.pt(), kaon.tpcNSigmaKa()); - if (kaon.hasTOF()) { - histos.fill(HIST("QAbefore/kaonTOFNSigma"), kaon.pt(), kaon.tofNSigmaKa()); - } + histos.fill(HIST("QAbefore/kaonDCAxy"), kaon.pt(), kaon.dcaXY()); + histos.fill(HIST("QAbefore/kaonDCAz"), kaon.pt(), kaon.dcaZ()); + histos.fill(HIST("QAbefore/kaonTPCNSigma"), kaon.pt(), kaon.tpcNSigmaKa()); + if (kaon.hasTOF() && !std::isnan(kaon.tofNSigmaKa())) { + histos.fill(HIST("QAbefore/kaonTOFNSigma"), kaon.pt(), kaon.tofNSigmaKa()); + } + if constexpr (!IsResoMicrotrack) { if constexpr (requires { kaon.tpcNClsFound(); }) { histos.fill(HIST("QAbefore/kaonTPCNcls"), kaon.tpcNClsFound()); } @@ -822,28 +995,21 @@ struct Xi1820Analysis { } } - if (!kaonCut(kaon)) + if (!kaonCut(kaon)) { continue; + } - if constexpr (!IsMix) { - nKaonsAfterCuts++; + if (!IsMix && cfgFillQA) { // QA after cuts histos.fill(HIST("QAafter/kaonPt"), kaon.pt()); histos.fill(HIST("QAafter/kaonEta"), kaon.eta()); - if constexpr (IsResoMicrotrack) { - histos.fill(HIST("QAafter/kaonDCAxy"), kaon.pt(), o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAxy(kaon.trackSelectionFlags())); - histos.fill(HIST("QAafter/kaonDCAz"), kaon.pt(), o2::aod::resomicrodaughter::ResoMicroTrackSelFlag::decodeDCAz(kaon.trackSelectionFlags())); - histos.fill(HIST("QAafter/kaonTPCNSigma"), kaon.pt(), o2::aod::resomicrodaughter::PidNSigma::getTPCnSigma(kaon.pidNSigmaKaFlag())); - if (kaon.hasTOF()) { - histos.fill(HIST("QAafter/kaonTOFNSigma"), kaon.pt(), o2::aod::resomicrodaughter::PidNSigma::getTOFnSigma(kaon.pidNSigmaKaFlag())); - } - } else { - histos.fill(HIST("QAafter/kaonDCAxy"), kaon.pt(), kaon.dcaXY()); - histos.fill(HIST("QAafter/kaonDCAz"), kaon.pt(), kaon.dcaZ()); - histos.fill(HIST("QAafter/kaonTPCNSigma"), kaon.pt(), kaon.tpcNSigmaKa()); - if (kaon.hasTOF()) { - histos.fill(HIST("QAafter/kaonTOFNSigma"), kaon.pt(), kaon.tofNSigmaKa()); - } + histos.fill(HIST("QAafter/kaonDCAxy"), kaon.pt(), kaon.dcaXY()); + histos.fill(HIST("QAafter/kaonDCAz"), kaon.pt(), kaon.dcaZ()); + histos.fill(HIST("QAafter/kaonTPCNSigma"), kaon.pt(), kaon.tpcNSigmaKa()); + if (kaon.hasTOF() && !std::isnan(kaon.tofNSigmaKa())) { + histos.fill(HIST("QAafter/kaonTOFNSigma"), kaon.pt(), kaon.tofNSigmaKa()); + } + if constexpr (!IsResoMicrotrack) { if constexpr (requires { kaon.tpcNClsFound(); }) { histos.fill(HIST("QAafter/kaonTPCNcls"), kaon.tpcNClsFound()); } @@ -853,304 +1019,272 @@ struct Xi1820Analysis { } } - int kaonCharge = kaon.sign(); - - // Loop over V0 candidates - for (const auto& v0 : v0s) { - // Lambda QA before cuts - if constexpr (!IsMix) { - histos.fill(HIST("QAbefore/lambdaDCAtoPV"), v0.pt(), v0.dcav0topv()); - histos.fill(HIST("QAbefore/lambdaMass"), v0.mLambda()); - histos.fill(HIST("QAbefore/lambdaMassAnti"), v0.mAntiLambda()); - histos.fill(HIST("QAbefore/lambdaPt"), v0.pt()); - histos.fill(HIST("QAbefore/lambdaEta"), v0.eta()); - histos.fill(HIST("QAbefore/lambdaCosPA"), v0.pt(), v0.v0CosPA()); - histos.fill(HIST("QAbefore/lambdaRadius"), v0.pt(), v0.transRadius()); - histos.fill(HIST("QAbefore/lambdaDauDCA"), v0.pt(), v0.daughDCA()); - histos.fill(HIST("QAbefore/lambdaDauPosDCA"), v0.pt(), std::abs(v0.dcapostopv())); - histos.fill(HIST("QAbefore/lambdaDauNegDCA"), v0.pt(), std::abs(v0.dcanegtopv())); - histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaPosPr"), v0.pt(), v0.daughterTPCNSigmaPosPr()); - histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaNegPi"), v0.pt(), v0.daughterTPCNSigmaNegPi()); - histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaPosPi"), v0.pt(), v0.daughterTPCNSigmaPosPi()); - histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaNegPr"), v0.pt(), v0.daughterTPCNSigmaNegPr()); - histos.fill(HIST("QAbefore/lambdaNCrossedRowsPos"), v0.pt(), v0.nCrossedRowsPos()); - histos.fill(HIST("QAbefore/lambdaNCrossedRowsNeg"), v0.pt(), v0.nCrossedRowsNeg()); - - // Calculate proper lifetime manually - float dx = v0.decayVtxX() - collision.posX(); - float dy = v0.decayVtxY() - collision.posY(); - float dz = v0.decayVtxZ() - collision.posZ(); - float l = std::sqrt(dx * dx + dy * dy + dz * dz); - float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); - auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; - histos.fill(HIST("QAbefore/lambdaProperLifetime"), v0.pt(), properLifetime); - histos.fill(HIST("QAbefore/lambdaArmenterosPodolanski"), v0.alpha(), v0.qtarm(), v0.pt()); - } - - // Try Lambda - bool isLambda = v0Cut(collision, v0, true); - // Try Anti-Lambda - bool isAntiLambda = v0Cut(collision, v0, false); - - if (!isLambda && !isAntiLambda) - continue; + selectedKaons.push_back(kaon); + } + const auto selectedLambdas = selectLambdas(lambdaCollision, v0s); + for (const auto& kaon : selectedKaons) { + const int kaonCharge = kaon.sign(); + for (const auto& selected : selectedLambdas) { + const auto& v0 = selected.candidate; if constexpr (!IsMix) { - nV0sAfterCuts++; - // QA after cuts (fill for whichever passes) - if (isLambda) { - histos.fill(HIST("QAafter/lambdaMass"), v0.mLambda()); - histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaPosPr"), v0.pt(), v0.daughterTPCNSigmaPosPr()); - histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaNegPi"), v0.pt(), v0.daughterTPCNSigmaNegPi()); - } - if (isAntiLambda) { - histos.fill(HIST("QAafter/lambdaMassAnti"), v0.mAntiLambda()); - histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaPosPi"), v0.pt(), v0.daughterTPCNSigmaPosPi()); - histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaNegPr"), v0.pt(), v0.daughterTPCNSigmaNegPr()); + const auto daughterIds = v0DaughterIds(v0); + if (kaon.trackId() == daughterIds[0] || kaon.trackId() == daughterIds[1]) { + continue; } - histos.fill(HIST("QAafter/lambdaDCAtoPV"), v0.pt(), v0.dcav0topv()); - histos.fill(HIST("QAafter/lambdaPt"), v0.pt()); - histos.fill(HIST("QAafter/lambdaEta"), v0.eta()); - histos.fill(HIST("QAafter/lambdaCosPA"), v0.pt(), v0.v0CosPA()); - histos.fill(HIST("QAafter/lambdaRadius"), v0.pt(), v0.transRadius()); - histos.fill(HIST("QAafter/lambdaDauDCA"), v0.pt(), v0.daughDCA()); - histos.fill(HIST("QAafter/lambdaDauPosDCA"), v0.pt(), std::abs(v0.dcapostopv())); - histos.fill(HIST("QAafter/lambdaDauNegDCA"), v0.pt(), std::abs(v0.dcanegtopv())); - histos.fill(HIST("QAafter/lambdaNCrossedRowsPos"), v0.pt(), v0.nCrossedRowsPos()); - histos.fill(HIST("QAafter/lambdaNCrossedRowsNeg"), v0.pt(), v0.nCrossedRowsNeg()); - - float dx = v0.decayVtxX() - collision.posX(); - float dy = v0.decayVtxY() - collision.posY(); - float dz = v0.decayVtxZ() - collision.posZ(); - float l = std::sqrt(dx * dx + dy * dy + dz * dz); - float p = std::sqrt(v0.px() * v0.px() + v0.py() * v0.py() + v0.pz() * v0.pz()); - auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; - histos.fill(HIST("QAafter/lambdaProperLifetime"), v0.pt(), properLifetime); - histos.fill(HIST("QAafter/lambdaArmenterosPodolanski"), v0.alpha(), v0.qtarm(), v0.pt()); } pKaon = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(kaon.pt(), kaon.eta(), kaon.phi(), MassKaonCharged)); - // K+ + Lambda -> Bkg channel for charged Xi(1820) - if (kaonCharge > 0 && isLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mLambda())); - pRes = pKaon + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut + // A failed rapidity or truth match skips only this mass hypothesis. + for (const auto& isLambda : std::array{true, false}) { + if (isLambda ? !selected.isLambda : !selected.isAntiLambda) { continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/kplus_lambda/hInvMassKplusLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentKplusLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentDelKplusLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } + + // K+ + Lambda -> Bkg channel for charged Xi(1820) + if (kaonCharge > 0 && isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mLambda())); + pRes = pKaon + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } - } else { - histos.fill(HIST("xi1820/kplus_lambda/hInvMassKplusLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentKplusLambda_Mix"), pRes.M(), pRes.Pt(), cent); + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/kplus_lambda/hInvMassKplusLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentKplusLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentDelKplusLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentDelKplusLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/kplus_lambda/hInvMassKplusLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentKplusLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kplus_lambda/hMassPtCentDelKplusLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } - } - } // End of Bkg - - // K+ + Anti-Lambda -> Signal channel for Anti-charged Xi(1820) - if (kaonCharge > 0 && isAntiLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mAntiLambda())); - pRes = pKaon + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut - continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/kplus_antilambda/hInvMassKplusAntiLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } // End of Bkg + + // K+ + Anti-Lambda -> Signal channel for Anti-charged Xi(1820) + if (kaonCharge > 0 && !isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mAntiLambda())); + pRes = pKaon + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } - if (additionalConfig.cfgFillRotBkg) { - for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); - if (additionalConfig.cfgRotKaon) { - lDaughterRot = pKaon; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = lDaughterRot + pLambda; - } else { - lDaughterRot = pLambda; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = pKaon + lDaughterRot; - } - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); - } else { - auto pRotDel = additionalConfig.cfgRotKaon ? deltaR(lDaughterRot, pLambda) : deltaR(lDaughterRot, pKaon); - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, pRotDel); + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/kplus_antilambda/hInvMassKplusAntiLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } + if (additionalConfig.cfgFillRotBkg) { + for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (cfgFillEventQA) { + histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); + } + if (additionalConfig.cfgRotKaon) { + lDaughterRot = pKaon; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = lDaughterRot + pLambda; + } else { + lDaughterRot = pLambda; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = pKaon + lDaughterRot; + } + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); + } else { + auto pRotDel = additionalConfig.cfgRotKaon ? deltaR(lDaughterRot, pLambda) : deltaR(lDaughterRot, pKaon); + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, pRotDel); + } } } - } - if constexpr (IsMC) { // Calculate Acceptance x efficiency for "the particle" channel - if (std::abs(v0.motherPDG()) != PdgChargedXi1820) - continue; - if (kaon.pdgCode() != PDG_t::kKPlus || v0.pdgCode() != PDG_t::kLambda0Bar) - continue; - if (kaon.motherId() != v0.motherId()) - continue; - auto pMCPt = v0.motherPt(); // Check particle's pT resolution - if (additionalConfig.cUseTruthRapidity && std::abs(v0.motherRap()) >= additionalConfig.cfgRapidityCut) // skip candidate if True rapidity of mother particle is outside of cut - continue; - histos.fill(HIST("MC/kplus_antilambda/hMCRecoInvMassKplusAntiLambda"), pRes.M()); - histos.fill(HIST("MC/kplus_antilambda/hMCRecoMassPtCentKplusAntiLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); - - // Detail QA histograms for truth particle -> Will be updated - } - } else { - histos.fill(HIST("xi1820/kplus_antilambda/hInvMassKplusAntiLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + if constexpr (IsMC) { // Calculate Acceptance x efficiency for "the particle" channel + if (v0.motherPDG() != -PdgChargedXi1820 || kaon.motherPDG() != v0.motherPDG()) { + continue; + } + if (kaon.pdgCode() != PDG_t::kKPlus || v0.pdgCode() != PDG_t::kLambda0Bar) { + continue; + } + if (v0.motherId() < 0 || kaon.motherId() != v0.motherId()) { + continue; + } + auto pMCPt = v0.motherPt(); // Check particle's pT resolution + if (additionalConfig.cUseTruthRapidity && !passesRapidity(v0.motherRap())) { // skip candidate if True rapidity of mother particle is outside of cut + continue; + } + histos.fill(HIST("MC/kplus_antilambda/hMCRecoInvMassKplusAntiLambda"), pRes.M()); + histos.fill(HIST("MC/kplus_antilambda/hMCRecoMassPtCentKplusAntiLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); + + fillTruthKaonQA(kaon); + fillTruthLambdaQA(v0, false); + } } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/kplus_antilambda/hInvMassKplusAntiLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentKplusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kplus_antilambda/hMassPtCentDelKplusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } - } - } // End of signal - - // K- + Lambda -> Signal channel for Xi(1820)- - if (kaonCharge < 0 && isLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mLambda())); - pRes = pKaon + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut - continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/kminus_lambda/hInvMassKminusLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } // End of signal + + // K- + Lambda -> Signal channel for Xi(1820)- + if (kaonCharge < 0 && isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mLambda())); + pRes = pKaon + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/kminus_lambda/hInvMassKminusLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } - if (additionalConfig.cfgFillRotBkg) { - for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); - if (additionalConfig.cfgRotKaon) { - lDaughterRot = pKaon; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = lDaughterRot + pLambda; - } else { - lDaughterRot = pLambda; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = pKaon + lDaughterRot; + if (additionalConfig.cfgFillRotBkg) { + for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (cfgFillEventQA) { + histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); + } + if (additionalConfig.cfgRotKaon) { + lDaughterRot = pKaon; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = lDaughterRot + pLambda; + } else { + lDaughterRot = pLambda; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = pKaon + lDaughterRot; + } + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); + } else { + auto pCandDel = additionalConfig.cfgRotKaon ? deltaR(lDaughterRot, pLambda) : deltaR(lDaughterRot, pKaon); + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, pCandDel); + } } - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); - } else { - auto pCandDel = additionalConfig.cfgRotKaon ? deltaR(lDaughterRot, pLambda) : deltaR(lDaughterRot, pKaon); - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, pCandDel); + } + if constexpr (IsMC) { // Calculate Acceptance x efficiency for "the particle" channel + if (v0.motherPDG() != PdgChargedXi1820 || kaon.motherPDG() != v0.motherPDG()) { + continue; + } + if (kaon.pdgCode() != PDG_t::kKMinus || v0.pdgCode() != PDG_t::kLambda0) { + continue; + } + if (v0.motherId() < 0 || kaon.motherId() != v0.motherId()) { + continue; } + auto pMCPt = v0.motherPt(); // Check particle's pT resolution + if (additionalConfig.cUseTruthRapidity && !passesRapidity(v0.motherRap())) { // skip candidate if True rapidity of mother particle is outside of cut + continue; + } + histos.fill(HIST("MC/kminus_lambda/hMCRecoInvMassKminusLambda"), pRes.M()); + histos.fill(HIST("MC/kminus_lambda/hMCRecoMassPtCentKminusLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); + + fillTruthKaonQA(kaon); + fillTruthLambdaQA(v0, true); } - } - if constexpr (IsMC) { // Calculate Acceptance x efficiency for "the particle" channel - if (std::abs(v0.motherPDG()) != PdgChargedXi1820) - continue; - if (kaon.pdgCode() != PDG_t::kKMinus || v0.pdgCode() != PDG_t::kLambda0) - continue; - if (kaon.motherId() != v0.motherId()) - continue; - auto pMCPt = v0.motherPt(); // Check particle's pT resolution - if (additionalConfig.cUseTruthRapidity && std::abs(v0.motherRap()) >= additionalConfig.cfgRapidityCut) // skip candidate if True rapidity of mother particle is outside of cut - continue; - histos.fill(HIST("MC/kminus_lambda/hMCRecoInvMassKminusLambda"), pRes.M()); - histos.fill(HIST("MC/kminus_lambda/hMCRecoMassPtCentKminusLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); - - // Detail QA histograms for the truth particle -> Will be updated - } - } else { - histos.fill(HIST("xi1820/kminus_lambda/hInvMassKminusLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda_Mix"), pRes.M(), pRes.Pt(), cent); } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/kminus_lambda/hInvMassKminusLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentKminusLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kminus_lambda/hMassPtCentDelKminusLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } - } - } // End of Signal - - // K- + Anti-Lambda -> Bkg channel for charged Xi(1820) - if (kaonCharge < 0 && isAntiLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mAntiLambda())); - pRes = pKaon + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut - continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/kminus_antilambda/hInvMassKminusAntiLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } // End of Signal + + // K- + Anti-Lambda -> Bkg channel for charged Xi(1820) + if (kaonCharge < 0 && !isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(v0.pt(), v0.eta(), v0.phi(), v0.mAntiLambda())); + pRes = pKaon + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } - } else { - histos.fill(HIST("xi1820/kminus_antilambda/hInvMassKminusAntiLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/kminus_antilambda/hInvMassKminusAntiLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } else { - auto pCandDel = deltaR(pLambda, pKaon); - histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/kminus_antilambda/hInvMassKminusAntiLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentKminusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pKaon); + histos.fill(HIST("xi1820/kminus_antilambda/hMassPtCentDelKminusAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } - } - } // End of Bkg - + } // End of Bkg + } // End of mass-hypothesis loop } // End of V0 loop } // Fill event QA for after-cuts counters (only for same-event) - if constexpr (!IsMix) { - histos.fill(HIST("Event/nLambdasAfterCuts"), nV0sAfterCuts); - histos.fill(HIST("Event/nKaonsAfterCuts"), nKaonsAfterCuts); + if (!IsMix && cfgFillEventQA) { + histos.fill(HIST("Event/nLambdasAfterCuts"), selectedLambdas.size()); + histos.fill(HIST("Event/nKaonsAfterCuts"), selectedKaons.size()); } } - template - void fillK0sLambda(const CollisionT& collision, const V0sT& k0sCands, const V0sT& lambdaCands) // Xi(1820) analysis: K0s + Lambda channel, No need to MicroTrack! + template + void fillK0sLambda(const CollisionT& collision, const V0sT& k0sCands, const V0sT& lambdaCands, const LambdaCollisionT& lambdaCollision) // Xi(1820) analysis: K0s + Lambda channel, No need to MicroTrack! { auto cent = collision.cent(); + // The charged same-event callback owns common event/Lambda QA when enabled. + const bool fillSharedQA = !hasChargedSameEventProcess(); // Fill event QA histograms - if constexpr (!IsMix) { - histos.fill(HIST("Event/posZ"), collision.posZ()); - histos.fill(HIST("Event/centrality"), cent); - histos.fill(HIST("Event/posZvsCent"), collision.posZ(), cent); - histos.fill(HIST("Event/nV0s"), lambdaCands.size()); - histos.fill(HIST("Event/nKaons"), k0sCands.size()); + if (!IsMix && cfgFillEventQA) { + if (fillSharedQA) { + histos.fill(HIST("Event/posZ"), collision.posZ()); + histos.fill(HIST("Event/centrality"), cent); + histos.fill(HIST("Event/posZvsCent"), collision.posZ(), cent); + histos.fill(HIST("Event/nV0s"), lambdaCands.size()); + } + histos.fill(HIST("Event/nK0s"), k0sCands.size()); } - if (additionalConfig.cConsiderHasV0s && (lambdaCands.size() < 1)) - return; // skip events that do not have V0s if the option is enabled - - int nV0sAfterCuts = 0; - int nKaonsAfterCuts = 0; + std::vector selectedK0s; + selectedK0s.reserve(k0sCands.size()); // Loop over V0s for K0s for (const auto& k0s : k0sCands) { // K0s QA before cuts - if constexpr (!IsMix) { + if (!IsMix && cfgFillQA) { histos.fill(HIST("QAbefore/k0sDCAtoPV"), k0s.pt(), k0s.dcav0topv()); histos.fill(HIST("QAbefore/k0sMass"), k0s.mK0Short()); histos.fill(HIST("QAbefore/k0sPt"), k0s.pt()); @@ -1173,12 +1307,11 @@ struct Xi1820Analysis { histos.fill(HIST("QAbefore/k0sArmenterosPodolanski"), k0s.alpha(), k0s.qtarm(), k0s.pt()); } - if (!k0sCut(collision, k0s)) + if (!k0sCut(collision, k0s)) { continue; - auto indexK0s = k0s.index(); + } - if constexpr (!IsMix) { - nKaonsAfterCuts++; + if (!IsMix && cfgFillQA) { // K0s QA after cuts histos.fill(HIST("QAafter/k0sDCAtoPV"), k0s.pt(), k0s.dcav0topv()); histos.fill(HIST("QAafter/k0sMass"), k0s.mK0Short()); @@ -1202,213 +1335,179 @@ struct Xi1820Analysis { histos.fill(HIST("QAafter/k0sArmenterosPodolanski"), k0s.alpha(), k0s.qtarm(), k0s.pt()); } - // Loop over V0s for Lambda - for (const auto& lambda : lambdaCands) { - - auto indexLambda = lambda.index(); + selectedK0s.push_back(k0s); + } + const auto selectedLambdas = selectLambdas(lambdaCollision, lambdaCands, fillSharedQA); + for (const auto& k0s : selectedK0s) { + for (const auto& selected : selectedLambdas) { + const auto& lambda = selected.candidate; if constexpr (!IsMix) { - if (indexLambda == indexK0s) // Avoid self-combination - continue; - if (sharesAnyDaughterId(v0DaughterIds(k0s), v0DaughterIds(lambda))) { + if (k0s.globalIndex() == lambda.globalIndex() || + sharesAnyDaughterId(v0DaughterIds(k0s), v0DaughterIds(lambda))) { continue; } - histos.fill(HIST("QAbefore/lambdaDCAtoPV"), lambda.pt(), lambda.dcav0topv()); - histos.fill(HIST("QAbefore/lambdaMass"), lambda.mLambda()); - histos.fill(HIST("QAbefore/lambdaMassAnti"), lambda.mAntiLambda()); - histos.fill(HIST("QAbefore/lambdaPt"), lambda.pt()); - histos.fill(HIST("QAbefore/lambdaEta"), lambda.eta()); - histos.fill(HIST("QAbefore/lambdaCosPA"), lambda.pt(), lambda.v0CosPA()); - histos.fill(HIST("QAbefore/lambdaRadius"), lambda.pt(), lambda.transRadius()); - histos.fill(HIST("QAbefore/lambdaDauDCA"), lambda.pt(), lambda.daughDCA()); - histos.fill(HIST("QAbefore/lambdaDauPosDCA"), lambda.pt(), std::abs(lambda.dcapostopv())); - histos.fill(HIST("QAbefore/lambdaDauNegDCA"), lambda.pt(), std::abs(lambda.dcanegtopv())); - histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaPosPr"), lambda.pt(), lambda.daughterTPCNSigmaPosPr()); - histos.fill(HIST("QAbefore/lambdaDaughterTPCNSigmaNegPi"), lambda.pt(), lambda.daughterTPCNSigmaNegPi()); - histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaPosPi"), lambda.pt(), lambda.daughterTPCNSigmaPosPi()); - histos.fill(HIST("QAbefore/lambdaAntiDaughterTPCNSigmaNegPr"), lambda.pt(), lambda.daughterTPCNSigmaNegPr()); - histos.fill(HIST("QAbefore/lambdaNCrossedRowsPos"), lambda.pt(), lambda.nCrossedRowsPos()); - histos.fill(HIST("QAbefore/lambdaNCrossedRowsNeg"), lambda.pt(), lambda.nCrossedRowsNeg()); - - // Calculate proper lifetime manually - float dx = lambda.decayVtxX() - collision.posX(); - float dy = lambda.decayVtxY() - collision.posY(); - float dz = lambda.decayVtxZ() - collision.posZ(); - float l = std::sqrt(dx * dx + dy * dy + dz * dz); - float p = std::sqrt(lambda.px() * lambda.px() + lambda.py() * lambda.py() + lambda.pz() * lambda.pz()); - auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; - histos.fill(HIST("QAbefore/lambdaProperLifetime"), lambda.pt(), properLifetime); - histos.fill(HIST("QAbefore/lambdaArmenterosPodolanski"), lambda.alpha(), lambda.qtarm(), lambda.pt()); - } - - // Try Lambda - bool isLambda = v0Cut(collision, lambda, true); - // Try Anti-Lambda - bool isAntiLambda = v0Cut(collision, lambda, false); - - if (!isLambda && !isAntiLambda) - continue; - - if constexpr (!IsMix) { - nV0sAfterCuts++; - // QA after cuts (fill for whichever passes) - if (isLambda) { - histos.fill(HIST("QAafter/lambdaMass"), lambda.mLambda()); - histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaPosPr"), lambda.pt(), lambda.daughterTPCNSigmaPosPr()); - histos.fill(HIST("QAafter/lambdaDaughterTPCNSigmaNegPi"), lambda.pt(), lambda.daughterTPCNSigmaNegPi()); - } - if (isAntiLambda) { - histos.fill(HIST("QAafter/lambdaMassAnti"), lambda.mAntiLambda()); - histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaPosPi"), lambda.pt(), lambda.daughterTPCNSigmaPosPi()); - histos.fill(HIST("QAafter/lambdaAntiDaughterTPCNSigmaNegPr"), lambda.pt(), lambda.daughterTPCNSigmaNegPr()); - } - histos.fill(HIST("QAafter/lambdaDCAtoPV"), lambda.pt(), lambda.dcav0topv()); - histos.fill(HIST("QAafter/lambdaPt"), lambda.pt()); - histos.fill(HIST("QAafter/lambdaEta"), lambda.eta()); - histos.fill(HIST("QAafter/lambdaCosPA"), lambda.pt(), lambda.v0CosPA()); - histos.fill(HIST("QAafter/lambdaRadius"), lambda.pt(), lambda.transRadius()); - histos.fill(HIST("QAafter/lambdaDauDCA"), lambda.pt(), lambda.daughDCA()); - histos.fill(HIST("QAafter/lambdaDauPosDCA"), lambda.pt(), std::abs(lambda.dcapostopv())); - histos.fill(HIST("QAafter/lambdaDauNegDCA"), lambda.pt(), std::abs(lambda.dcanegtopv())); - histos.fill(HIST("QAafter/lambdaNCrossedRowsPos"), lambda.pt(), lambda.nCrossedRowsPos()); - histos.fill(HIST("QAafter/lambdaNCrossedRowsNeg"), lambda.pt(), lambda.nCrossedRowsNeg()); - - float dx = lambda.decayVtxX() - collision.posX(); - float dy = lambda.decayVtxY() - collision.posY(); - float dz = lambda.decayVtxZ() - collision.posZ(); - float l = std::sqrt(dx * dx + dy * dy + dz * dz); - float p = std::sqrt(lambda.px() * lambda.px() + lambda.py() * lambda.py() + lambda.pz() * lambda.pz()); - auto properLifetime = (l / (p + SmallMomentumDenominator)) * MassLambda; - histos.fill(HIST("QAafter/lambdaProperLifetime"), lambda.pt(), properLifetime); - histos.fill(HIST("QAafter/lambdaArmenterosPodolanski"), lambda.alpha(), lambda.qtarm(), lambda.pt()); } // 4-vectors ROOT::Math::PxPyPzEVector pK0s, pLambda, pRes, lDaughterRot, lResonanceRot; pK0s = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(k0s.pt(), k0s.eta(), k0s.phi(), MassK0Short)); - if (isLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(lambda.pt(), lambda.eta(), lambda.phi(), lambda.mLambda())); - pRes = pK0s + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut + // Lambda and anti-Lambda may both pass; evaluate them independently. + for (const auto& isLambda : std::array{true, false}) { + if (isLambda ? !selected.isLambda : !selected.isAntiLambda) { continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/k0s_lambda/hInvMassK0sLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pK0s); - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } + + if (isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(lambda.pt(), lambda.eta(), lambda.phi(), lambda.mLambda())); + pRes = pK0s + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/k0s_lambda/hInvMassK0sLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pK0s); + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } - if (additionalConfig.cfgFillRotBkg) { - for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); - lDaughterRot = pK0s; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = lDaughterRot + pLambda; - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); - } else { - auto lRotDel = deltaR(lDaughterRot, pLambda); - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, lRotDel); + if (additionalConfig.cfgFillRotBkg) { + for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (cfgFillEventQA) { + histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); + } + lDaughterRot = pK0s; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = lDaughterRot + pLambda; + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); + } else { + auto lRotDel = deltaR(lDaughterRot, pLambda); + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, lRotDel); + } } } - } - if constexpr (IsMC) { // Calculate Acceptance x efficiency - if (std::abs(lambda.motherPDG()) != PdgXi1820Zero) - continue; - if (std::abs(k0s.pdgCode()) != PDG_t::kK0Short || lambda.pdgCode() != PDG_t::kLambda0) - continue; - if (k0s.motherId() != lambda.motherId()) - continue; - auto pMCPt = lambda.motherPt(); // Check particle's pT resolution - if (additionalConfig.cUseTruthRapidity && std::abs(lambda.motherRap()) >= additionalConfig.cfgRapidityCut) // skip candidate if True rapidity of mother particle is outside of cut - continue; - histos.fill(HIST("MC/k0s_lambda/hMCRecoInvMassK0sLambda"), pRes.M()); - histos.fill(HIST("MC/k0s_lambda/hMCRecoMassPtCentK0sLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); - // Detail QA histograms for truth particle -> Will be updated - } - } else { - histos.fill(HIST("xi1820/k0s_lambda/hInvMassK0sLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda_Mix"), pRes.M(), pRes.Pt(), cent); + if constexpr (IsMC) { // Calculate Acceptance x efficiency + if (lambda.motherPDG() != PdgXi1820Zero || k0s.motherPDG() != lambda.motherPDG()) { + continue; + } + if (std::abs(k0s.pdgCode()) != PDG_t::kK0Short || lambda.pdgCode() != PDG_t::kLambda0) { + continue; + } + if (lambda.motherId() < 0 || k0s.motherId() != lambda.motherId()) { + continue; + } + auto pMCPt = lambda.motherPt(); // Check particle's pT resolution + if (additionalConfig.cUseTruthRapidity && !passesRapidity(lambda.motherRap())) { // skip candidate if True rapidity of mother particle is outside of cut + continue; + } + histos.fill(HIST("MC/k0s_lambda/hMCRecoInvMassK0sLambda"), pRes.M()); + histos.fill(HIST("MC/k0s_lambda/hMCRecoMassPtCentK0sLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); + fillTruthLambdaQA(lambda, true); + if (cfgFillTruthQA) { + histos.fill(HIST("QAMCTrue/k0sPt"), k0s.pt()); + } + } } else { - auto pCandDel = deltaR(pLambda, pK0s); - histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/k0s_lambda/hInvMassK0sLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentK0sLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pK0s); + histos.fill(HIST("xi1820/k0s_lambda/hMassPtCentDelK0sLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } - } - } // End of Lambda - - if (isAntiLambda) { - pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(lambda.pt(), lambda.eta(), lambda.phi(), lambda.mAntiLambda())); - pRes = pK0s + pLambda; - auto pCandRapidity = pRes.Rapidity(); - if (std::abs(pCandRapidity) >= additionalConfig.cfgRapidityCut) // skip candidate if reconstructed rapidity is outside of cut - continue; - if constexpr (!IsMix) { - histos.fill(HIST("xi1820/k0s_antilambda/hInvMassK0sAntiLambda"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda"), pRes.M(), pRes.Pt(), cent); - } else { - auto pCandDel = deltaR(pLambda, pK0s); - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } // End of Lambda + + if (!isLambda) { + pLambda = ROOT::Math::PxPyPzEVector(ROOT::Math::PtEtaPhiMVector(lambda.pt(), lambda.eta(), lambda.phi(), lambda.mAntiLambda())); + pRes = pK0s + pLambda; + auto pCandRapidity = pRes.Rapidity(); + if (!passesRapidity(pCandRapidity)) { // skip candidate if reconstructed rapidity is outside of cut + continue; } + if constexpr (!IsMix) { + histos.fill(HIST("xi1820/k0s_antilambda/hInvMassK0sAntiLambda"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pK0s); + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda"), pRes.M(), pRes.Pt(), cent, pCandDel); + } - if (additionalConfig.cfgFillRotBkg) { - for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { - auto lRotAngle = additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); - histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); - lDaughterRot = pK0s; - ROOT::Math::RotationZ rot(lRotAngle); - auto p3 = rot * lDaughterRot.Vect(); - lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); - lResonanceRot = lDaughterRot + pLambda; - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); - } else { - auto lRotDel = deltaR(lDaughterRot, pLambda); - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, lRotDel); + if (additionalConfig.cfgFillRotBkg) { + for (int i = 0; i < additionalConfig.cfgNrotBkg; i++) { + const auto lRotAngle = additionalConfig.cfgNrotBkg == 1 + ? 0.5f * (additionalConfig.cfgMinRot + additionalConfig.cfgMaxRot) + : additionalConfig.cfgMinRot + i * ((additionalConfig.cfgMaxRot - additionalConfig.cfgMinRot) / (additionalConfig.cfgNrotBkg - 1)); + if (cfgFillEventQA) { + histos.fill(HIST("Event/hRotBkg"), lRotAngle - o2::constants::math::PI); + } + lDaughterRot = pK0s; + ROOT::Math::RotationZ rot(lRotAngle); + auto p3 = rot * lDaughterRot.Vect(); + lDaughterRot = LorentzVectorSetXYZM(p3.X(), p3.Y(), p3.Z(), lDaughterRot.M()); + lResonanceRot = lDaughterRot + pLambda; + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent); + } else { + auto lRotDel = deltaR(lDaughterRot, pLambda); + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Rot"), lResonanceRot.M(), lResonanceRot.Pt(), cent, lRotDel); + } } } - } - if constexpr (IsMC) { // Calculate Acceptance x efficiency - if (std::abs(lambda.motherPDG()) != PdgXi1820Zero) - continue; - if (std::abs(k0s.pdgCode()) != PDG_t::kK0Short || lambda.pdgCode() != PDG_t::kLambda0Bar) - continue; - if (k0s.motherId() != lambda.motherId()) - continue; - auto pMCPt = lambda.motherPt(); // Check particle's pT resolution - if (additionalConfig.cUseTruthRapidity && std::abs(lambda.motherRap()) >= additionalConfig.cfgRapidityCut) // skip candidate if True rapidity of mother particle is outside of cut - continue; - histos.fill(HIST("MC/k0s_antilambda/hMCRecoInvMassK0sAntiLambda"), pRes.M()); - histos.fill(HIST("MC/k0s_antilambda/hMCRecoMassPtCentK0sAntiLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); - // Detail QA histograms for truth particle -> Will be updated - } - } else { - histos.fill(HIST("xi1820/k0s_antilambda/hInvMassK0sAntiLambda_Mix"), pRes.M()); - if (!(additionalConfig.cfgDelCheck)) { - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + if constexpr (IsMC) { // Calculate Acceptance x efficiency + if (lambda.motherPDG() != -PdgXi1820Zero || k0s.motherPDG() != lambda.motherPDG()) { + continue; + } + if (std::abs(k0s.pdgCode()) != PDG_t::kK0Short || lambda.pdgCode() != PDG_t::kLambda0Bar) { + continue; + } + if (lambda.motherId() < 0 || k0s.motherId() != lambda.motherId()) { + continue; + } + auto pMCPt = lambda.motherPt(); // Check particle's pT resolution + if (additionalConfig.cUseTruthRapidity && !passesRapidity(lambda.motherRap())) { // skip candidate if True rapidity of mother particle is outside of cut + continue; + } + histos.fill(HIST("MC/k0s_antilambda/hMCRecoInvMassK0sAntiLambda"), pRes.M()); + histos.fill(HIST("MC/k0s_antilambda/hMCRecoMassPtCentK0sAntiLambda"), pRes.M(), pRes.Pt(), cent, pMCPt); + fillTruthLambdaQA(lambda, false); + if (cfgFillTruthQA) { + histos.fill(HIST("QAMCTrue/k0sPt"), k0s.pt()); + } + } } else { - auto pCandDel = deltaR(pLambda, pK0s); - histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + histos.fill(HIST("xi1820/k0s_antilambda/hInvMassK0sAntiLambda_Mix"), pRes.M()); + if (!(additionalConfig.cfgDelCheck)) { + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentK0sAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent); + } else { + auto pCandDel = deltaR(pLambda, pK0s); + histos.fill(HIST("xi1820/k0s_antilambda/hMassPtCentDelK0sAntiLambda_Mix"), pRes.M(), pRes.Pt(), cent, pCandDel); + } } } - } + } // End of mass-hypothesis loop } // End of loop over Lambda candidates } // End of loop over K0s candidates // Fill event QA for after-cuts counters (only for same-event) - if constexpr (!IsMix) { - histos.fill(HIST("Event/nLambdasAfterCuts"), nV0sAfterCuts); - histos.fill(HIST("Event/nKaonsAfterCuts"), nKaonsAfterCuts); + if (!IsMix && cfgFillEventQA) { + if (fillSharedQA) { + histos.fill(HIST("Event/nLambdasAfterCuts"), selectedLambdas.size()); + } + histos.fill(HIST("Event/nK0sAfterCuts"), selectedK0s.size()); } } @@ -1418,137 +1517,265 @@ struct Xi1820Analysis { } PROCESS_SWITCH(Xi1820Analysis, processDummy, "Process Dummy", true); - void processDataWithTracks(const aod::ResoCollision& resoCollision, + void processDataWithTracks(ResoCollisions::iterator const& resoCollision, aod::ResoV0s const& resoV0s, - aod::ResoTracks const& resoTracks) + ResoTracks const& resoTracks) { - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) + if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) { return; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillChargedKLambda(resoCollision, resoV0s, resoTracks); + } + fillChargedKLambda(resoCollision, resoV0s, resoTracks, resoCollision); } PROCESS_SWITCH(Xi1820Analysis, processDataWithTracks, "Process Event with ResoTracks", false); - void processDataWithMicroTracks(const aod::ResoCollision& resoCollision, + void processDataWithMicroTracks(ResoCollisions::iterator const& resoCollision, aod::ResoV0s const& resoV0s, - aod::ResoMicroTracks const& resoMicroTracks) + ResoMicroTracks const& resoMicroTracks) { - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) + if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) { return; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillChargedKLambda(resoCollision, resoV0s, resoMicroTracks); + } + fillChargedKLambda(resoCollision, resoV0s, resoMicroTracks, resoCollision); } PROCESS_SWITCH(Xi1820Analysis, processDataWithMicroTracks, "Process Event with ResoMicroTracks", false); - void processMixedEventWithTracks(const aod::ResoCollisions& resoCollisions, - aod::ResoV0s const& resoV0s, - aod::ResoTracks const& resoTracks) + template + bool acceptsMixedCollisions(const FirstCollision& first, const SecondCollision& second) const { + if (additionalConfig.cRecoINELgt0.value && (!first.isRecINELgt0() || !second.isRecINELgt0())) { + return false; + } + return std::isfinite(first.bMagField()) && first.bMagField() == second.bMagField(); + } - auto v0sTracksTuple = std::make_tuple(resoTracks, resoV0s); - BinningTypeVertexContributor colBinning{{cfgVtxBins, cfgMultBins}, true}; - Pair pairs{colBinning, nEvtMixing, -1, resoCollisions, v0sTracksTuple, &cache}; - + void processMixedEventWithTracks(ResoCollisions const& resoCollisions, + aod::ResoV0s const& resoV0s, + ResoTracks const& resoTracks) + { + if (nEvtMixing <= 0) { + return; + } + auto tracksTuple = std::make_tuple(resoTracks, resoV0s); + BinningTypeVertexContributor binning{{cfgVtxBins, cfgMultBins}, true}; + Pair pairs{binning, nEvtMixing, -1, resoCollisions, tracksTuple, &cache}; for (auto& [collision1, tracks1, collision2, v0s2] : pairs) { // o2-linter: disable=const-ref-in-for-loop (structured bindings from Pair iterator cannot be const) - if (additionalConfig.cRecoINELgt0 && !collision1.isRecINELgt0()) - continue; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillChargedKLambda(collision1, v0s2, tracks1); + if (!acceptsMixedCollisions(collision1, collision2)) { + continue; + } + fillChargedKLambda(collision1, v0s2, tracks1, collision2); } } PROCESS_SWITCH(Xi1820Analysis, processMixedEventWithTracks, "Process Mixed Event with ResoTracks", false); - void processMixedEventWithMicroTracks(const aod::ResoCollisions& resoCollisions, + void processMixedEventWithMicroTracks(ResoCollisions const& resoCollisions, aod::ResoV0s const& resoV0s, - aod::ResoMicroTracks const& resoMicroTracks) + ResoMicroTracks const& resoMicroTracks) { + if (nEvtMixing <= 0) { + return; + } + auto tracksTuple = std::make_tuple(resoMicroTracks, resoV0s); + BinningTypeVertexContributor binning{{cfgVtxBins, cfgMultBins}, true}; + Pair pairs{binning, nEvtMixing, -1, resoCollisions, tracksTuple, &cache}; + for (auto& [collision1, tracks1, collision2, v0s2] : pairs) { // o2-linter: disable=const-ref-in-for-loop (structured bindings from Pair iterator cannot be const) + if (!acceptsMixedCollisions(collision1, collision2)) { + continue; + } + fillChargedKLambda(collision1, v0s2, tracks1, collision2); + } + } + PROCESS_SWITCH(Xi1820Analysis, processMixedEventWithMicroTracks, "Process Mixed Event with ResoMicroTracks 001", false); - auto v0sTracksTuple = std::make_tuple(resoV0s, resoMicroTracks); - BinningTypeVertexContributor colBinning{{cfgVtxBins, cfgMultBins}, true}; - Pair pairs{colBinning, nEvtMixing, -1, resoCollisions, v0sTracksTuple, &cache}; + template + static std::shared_ptr copyMixingRows(DaughterTable const& daughters) + { + // Copy only this collision's rows, without retaining the source DF buffers. + const auto source = daughters.asArrowTableConstrained(); + if (source->num_rows() == 0) { + return arrow::Table::MakeEmpty(source->schema()).ValueOrDie(); + } + std::vector> columns; + columns.reserve(source->num_columns()); + for (const auto& column : source->columns()) { + columns.push_back(std::make_shared(arrow::Concatenate(column->chunks()).ValueOrDie())); + } + return arrow::Table::Make(source->schema(), columns); + } - for (auto& [collision1, v0s1, collision2, tracks2] : pairs) { // o2-linter: disable=const-ref-in-for-loop (structured bindings from Pair iterator cannot be const) - if (additionalConfig.cRecoINELgt0 && !collision2.isRecINELgt0()) - continue; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillChargedKLambda(collision2, v0s1, tracks2); + // FIFO history is task-local and retains only the copied MicroTrack rows. + // Use a single run/condition set; a field change clears all stored event pools. + void processMEDF(ResoCollisions::iterator const& collision, + ResoMicroTracks const& tracks, + aod::ResoV0s const& v0s) + { + if (nEvtMixing <= 0 || !std::isfinite(collision.bMagField())) { + return; + } + if (mixingBField != collision.bMagField()) { + mixingPools.clear(); + mixingBField = collision.bMagField(); + } + BinningTypeVertexContributor binning{{cfgVtxBins, cfgMultBins}, true}; + const int bin = binning.getBin(std::make_tuple(collision.posZ(), collision.cent())); + if (bin < 0) { + return; + } + const MixingCollision current{.x = collision.posX(), .y = collision.posY(), .z = collision.posZ(), .centrality = collision.cent(), .recINELgt0 = collision.isRecINELgt0()}; + auto& pool = mixingPools[bin]; + for (const auto& previous : pool) { + if (additionalConfig.cRecoINELgt0 && (!previous.collision.recINELgt0 || !current.recINELgt0)) { + continue; + } + fillChargedKLambda(previous.collision, v0s, ResoMicroTracks{previous.daughters}, current); + } + // Insert after pairing; DF-local IDs may repeat across frames. + pool.push_back({current, copyMixingRows(tracks)}); + while (pool.size() > static_cast(nEvtMixing.value)) { + pool.pop_front(); } } - PROCESS_SWITCH(Xi1820Analysis, processMixedEventWithMicroTracks, "Process Mixed Event with ResoMicroTracks", false); + PROCESS_SWITCH(Xi1820Analysis, processMEDF, "Process cross-DF mixing with MicroTracks 001", false); // K0s + Lambda analysis - void processK0sLambda(const aod::ResoCollision& resoCollision, + void processK0sLambda(ResoCollisions::iterator const& resoCollision, aod::ResoV0s const& resoV0s) { - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) + if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) { return; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillK0sLambda(resoCollision, resoV0s, resoV0s); + } + fillK0sLambda(resoCollision, resoV0s, resoV0s, resoCollision); } PROCESS_SWITCH(Xi1820Analysis, processK0sLambda, "Process K0s + Lambda", false); // K0s + Lambda mixed event analysis - void processK0sLambdaMixedEvent(const aod::ResoCollisions& resoCollisions, + void processK0sLambdaMixedEvent(ResoCollisions const& resoCollisions, aod::ResoV0s const& resoV0s) { + if (nEvtMixing <= 0) { + return; + } auto v0sV0sTuple = std::make_tuple(resoV0s, resoV0s); BinningTypeVertexContributor colBinning{{cfgVtxBins, cfgMultBins}, true}; - Pair pairs{colBinning, nEvtMixing, -1, resoCollisions, v0sV0sTuple, &cache}; + Pair pairs{colBinning, nEvtMixing, -1, resoCollisions, v0sV0sTuple, &cache}; for (auto& [collision1, k0s1, collision2, lambda2] : pairs) { // o2-linter: disable=const-ref-in-for-loop (structured bindings from Pair iterator cannot be const) - if (additionalConfig.cRecoINELgt0 && !collision1.isRecINELgt0()) - continue; // skip event if RecoINEL>0 selection is enabled and event does not pass it - fillK0sLambda(collision1, k0s1, lambda2); + if (!acceptsMixedCollisions(collision1, collision2)) { + continue; + } + fillK0sLambda(collision1, k0s1, lambda2, collision2); } } PROCESS_SWITCH(Xi1820Analysis, processK0sLambdaMixedEvent, "Process K0s + Lambda Mixed Event", false); + // Cross-DF neutral mixing: old-event K0s with current-event Lambda/anti-Lambda. + // Keep a separate bounded V0 pool; use inputs from one run/condition set. + void processK0sLambdaMEDF(ResoCollisions::iterator const& collision, + aod::ResoV0s const& v0s) + { + if (nEvtMixing <= 0 || !std::isfinite(collision.bMagField())) { + return; + } + if (neutralMixingBField != collision.bMagField()) { + neutralMixingPools.clear(); + neutralMixingBField = collision.bMagField(); + } + BinningTypeVertexContributor binning{{cfgVtxBins, cfgMultBins}, true}; + const int bin = binning.getBin(std::make_tuple(collision.posZ(), collision.cent())); + if (bin < 0) { + return; + } + const MixingCollision current{.x = collision.posX(), .y = collision.posY(), .z = collision.posZ(), .centrality = collision.cent(), .recINELgt0 = collision.isRecINELgt0()}; + auto& pool = neutralMixingPools[bin]; + const bool fillCurrentQA = cfgFillEventQA && (!additionalConfig.cRecoINELgt0 || current.recINELgt0); + if (fillCurrentQA) { + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hPosZvsCent"), current.posZ(), current.cent()); + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hPoolSize"), pool.size()); + } + int nPartners = 0; + for (const auto& previous : pool) { + if (additionalConfig.cRecoINELgt0 && (!previous.collision.recINELgt0 || !current.recINELgt0)) { + continue; + } + if (cfgFillEventQA) { + ++nPartners; + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hPosZPair"), previous.collision.posZ(), current.posZ()); + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hCentPair"), previous.collision.cent(), current.cent()); + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hV0CountsPair"), previous.daughters->num_rows(), v0s.size()); + } + // Each V0 uses its own event's PV; original row IDs are not compared across events. + fillK0sLambda(previous.collision, aod::ResoV0s{previous.daughters}, v0s, current); + } + if (fillCurrentQA) { + histos.fill(HIST("MixingQA/K0sLambdaMEDF/hNPartners"), nPartners); + } + // Insert after pairing to avoid self-mixing, even when DF-local IDs repeat. + pool.push_back({current, copyMixingRows(v0s)}); + while (pool.size() > static_cast(nEvtMixing.value)) { + pool.pop_front(); + } + } + PROCESS_SWITCH(Xi1820Analysis, processK0sLambdaMEDF, "Process cross-DF K0s + Lambda mixing", false); + // MC processes for charged K + Lambda analysis void processMCWithTracks(ResoMCCols::iterator const& resoMCcollision, soa::Join const& resoMCV0s, - soa::Join const& resoMCTracks) + soa::Join const& resoMCTracks) { - if (additionalConfig.cRecoINELgt0 && !resoMCcollision.isRecINELgt0()) - return; // skip event if RecoINEL>0 selection is enabled and event does not pass it - if (!resoMCcollision.isInAfterAllCuts()) // MC event selection + if ((additionalConfig.cRecoINELgt0 && !resoMCcollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoMCcollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoMCcollision.isVtxIn10())) { return; - fillChargedKLambda(resoMCcollision, resoMCV0s, resoMCTracks); + } + fillChargedKLambda(resoMCcollision, resoMCV0s, resoMCTracks, resoMCcollision); } PROCESS_SWITCH(Xi1820Analysis, processMCWithTracks, "Process MC for charged K + Lambda", false); void processMCK0sLambda(ResoMCCols::iterator const& resoMCCollision, soa::Join const& resoMCV0s) { - if (additionalConfig.cRecoINELgt0 && !resoMCCollision.isRecINELgt0()) - return; // skip event if RecoINEL>0 selection is enabled and event does not pass it - if (!resoMCCollision.isInAfterAllCuts()) // MC event selection + if ((additionalConfig.cRecoINELgt0 && !resoMCCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoMCCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoMCCollision.isVtxIn10())) { return; - fillK0sLambda(resoMCCollision, resoMCV0s, resoMCV0s); + } + fillK0sLambda(resoMCCollision, resoMCV0s, resoMCV0s, resoMCCollision); } PROCESS_SWITCH(Xi1820Analysis, processMCK0sLambda, "Process MC K0s + Lambda", false); - void processMCWithMicroTracks(const aod::ResoCollision& /*resoCollision*/, - aod::ResoV0s const& /*resoV0s*/, - aod::ResoMicroTracks const& /*resoMicroTracks*/, - aod::McParticles const& /*mcParticles*/) + void processMCWithMicroTracks(ResoMCCols::iterator const& resoCollision, + soa::Join const& resoV0s, + soa::Join const& resoMicroTracks) { - // TODO: Implement MC truth matching for K± + Lambda with MicroTracks - // But is this really necessary? -> Most of the injected MC sizes are already within small-train limit. + if ((additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoCollision.isVtxIn10())) { + return; + } + // Module 001 contains only selected reconstructed collisions. + fillChargedKLambda(resoCollision, resoV0s, resoMicroTracks, resoCollision); } - PROCESS_SWITCH(Xi1820Analysis, processMCWithMicroTracks, "Process MC with ResoMicroTracks (placeholder)", false); + PROCESS_SWITCH(Xi1820Analysis, processMCWithMicroTracks, "Process MC with ResoMicroTracks 001", false); void processMCGen(ResoMCCols::iterator const& resoCollision, // Calculate denominator for the acceptance x efficiency and a part of Event-factor (for selected evennts) - aod::ResoMCParents const& resoParents) + aod::ResoMCParents_001 const& resoParents) { auto multiplicity = resoCollision.mcMultiplicity(); auto inCent = resoCollision.cent(); - if (additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) // Check reco INELgt0 - return; - if (!resoCollision.isInAfterAllCuts()) + if ((additionalConfig.cRecoINELgt0 && !resoCollision.isRecINELgt0()) || + (additionalConfig.cMCINELgt0 && !resoCollision.isINELgt0()) || + (additionalConfig.cMCVtxIn10 && !resoCollision.isVtxIn10())) { return; + } histos.fill(HIST("multQA/h2MultCentMC"), inCent, multiplicity); for (const auto& part : resoParents) { // loop over all pre-filtered Gen particle on selected events auto pdgMother = part.pdgCode(); - if (std::abs(pdgMother) != PdgChargedXi1820 && std::abs(pdgMother) != PdgXi1820Zero) + if (std::abs(pdgMother) != PdgChargedXi1820 && std::abs(pdgMother) != PdgXi1820Zero) { continue; - if (std::abs(part.y()) >= additionalConfig.cfgRapidityCut) + } + if (!passesRapidity(part.y())) { continue; // skip if rapidity of the particle is outside of cut + } auto motherPt = part.pt(); auto daughter1PDG = part.daughterPDG1(); auto daughter2PDG = part.daughterPDG2(); @@ -1579,17 +1806,16 @@ struct Xi1820Analysis { void processMCTruth(aod::McParticles const& mcParticles) // ->Let's keep it and use for injected MC QA...! { // Process MC generated particles (no reconstruction requirement) - // Xi(1820)0 PDG code: 123314 (neutral, decays to K+ Lambda or K0s Lambda) - // Note: PDG doesn't have separate codes for charge states in this case + // Charged and neutral Xi(1820) states have distinct PDG codes. for (const auto& mcParticle : mcParticles) { // Look for Xi(1820) - PDG code can vary, check for resonance mass ~1820 MeV int pdg = mcParticle.pdgCode(); - // Xi(1820)0: PDG 123314 - // Check if it's Xi(1820) or similar resonance - if (std::abs(pdg) != PdgChargedXi1820 && std::abs(pdg) != PdgXi1820Zero) + // Select only the configured charged and neutral Xi(1820) states. + if (std::abs(pdg) != PdgChargedXi1820 && std::abs(pdg) != PdgXi1820Zero) { continue; + } // Fill generated level histograms auto pt = mcParticle.pt(); @@ -1602,8 +1828,9 @@ struct Xi1820Analysis { // Get daughters auto daughters = mcParticle.daughters_as(); - if (daughters.size() != ExpectedDaughters) + if (daughters.size() != ExpectedDaughters) { continue; + } int daughter1PDG = 0, daughter2PDG = 0; ROOT::Math::PxPyPzEVector p1, p2, pMother; @@ -1657,7 +1884,7 @@ struct Xi1820Analysis { PROCESS_SWITCH(Xi1820Analysis, processMCTruth, "Process MC Truth particles", false); }; -WorkflowSpec defineDataProcessing(ConfigContext const& cfgc) +WorkflowSpec defineDataProcessing(ConfigContext const& context) { - return WorkflowSpec{adaptAnalysisTask(cfgc)}; + return WorkflowSpec{adaptAnalysisTask(context)}; }