From fcfa17b3154282c8a70194d0af431cfef56f8d00 Mon Sep 17 00:00:00 2001 From: Anton Riedel Date: Wed, 1 Apr 2026 12:38:41 +0200 Subject: [PATCH] Feat: Add dalitz plot for track-v0 analysis --- PWGCF/Femto/Core/femtoUtils.h | 4 +- PWGCF/Femto/Core/pairBuilder.h | 21 ++- PWGCF/Femto/Core/pairHistManager.h | 228 +++++++++++++++++--------- PWGCF/Femto/Core/pairProcessHelpers.h | 20 +-- PWGCF/Femto/Core/trackHistManager.h | 2 +- PWGCF/Femto/Core/tripletHistManager.h | 6 +- 6 files changed, 181 insertions(+), 100 deletions(-) diff --git a/PWGCF/Femto/Core/femtoUtils.h b/PWGCF/Femto/Core/femtoUtils.h index 000850f3be9..1da7fde8635 100644 --- a/PWGCF/Femto/Core/femtoUtils.h +++ b/PWGCF/Femto/Core/femtoUtils.h @@ -63,7 +63,7 @@ float itsSignal(T const& track) return static_cast(signal); }; -inline double getMass(int pdgCode) +inline double getPdgMass(int pdgCode) { // use this function instead of TDatabasePDG to return masses defined in the PhysicsConstants.h header // this approach saves a lot of memory and important partilces like deuteron are missing in TDatabasePDG anyway @@ -119,7 +119,7 @@ inline double getMass(int pdgCode) mass = o2::constants::physics::MassOmegaMinus; break; default: - LOG(fatal) << "PDG code is not suppored"; + LOG(warn) << "PDG code is not suppored. Return 0..."; } return mass; } diff --git a/PWGCF/Femto/Core/pairBuilder.h b/PWGCF/Femto/Core/pairBuilder.h index 0a44ba639fc..af6e9aa51bb 100644 --- a/PWGCF/Femto/Core/pairBuilder.h +++ b/PWGCF/Femto/Core/pairBuilder.h @@ -34,6 +34,8 @@ #include #include +#include + #include #include #include @@ -595,13 +597,28 @@ class PairTrackV0Builder mTrackCleaner.init(confTrackCleaner); mV0Cleaner.init(confV0Cleaner); + int pdgCodePosDau = 0; + int pdgCodeNegDau = 0; + if (modes::isEqual(v0Type, modes::V0::kK0short)) { + pdgCodeNegDau = kPiPlus; + pdgCodeNegDau = kPiMinus; + } else if (modes::isEqual(v0Type, modes::V0::kLambda) || modes::isEqual(v0Type, modes::V0::kAntiLambda)) { + if (confV0Selection.sign.value > 0) { + pdgCodeNegDau = kProton; + pdgCodeNegDau = kPiMinus; + } else { + pdgCodeNegDau = kProtonBar; + pdgCodeNegDau = kPiPlus; + } + } + mPairHistManagerSe.template init(registry, pairHistSpec, confPairBinning, confPairCuts); - mPairHistManagerSe.setMass(confTrackSelection.pdgCodeAbs.value, confV0Selection.pdgCodeAbs.value); + mPairHistManagerSe.setMass(confTrackSelection.pdgCodeAbs.value, 0, 0, confV0Selection.pdgCodeAbs.value, pdgCodePosDau, pdgCodeNegDau); mPairHistManagerSe.setCharge(confTrackSelection.chargeAbs.value, 1); mCprSe.init(registry, cprHistSpec, confCpr, confTrackSelection.chargeAbs.value); mPairHistManagerMe.template init(registry, pairHistSpec, confPairBinning, confPairCuts); - mPairHistManagerMe.setMass(confTrackSelection.pdgCodeAbs.value, confV0Selection.pdgCodeAbs.value); + mPairHistManagerMe.setMass(confTrackSelection.pdgCodeAbs.value, 0, 0, confV0Selection.pdgCodeAbs.value, pdgCodePosDau, pdgCodeNegDau); mPairHistManagerMe.setCharge(confTrackSelection.chargeAbs.value, 1); mCprMe.init(registry, cprHistSpec, confCpr, confTrackSelection.chargeAbs.value); mPc.template init(confPairCuts); diff --git a/PWGCF/Femto/Core/pairHistManager.h b/PWGCF/Femto/Core/pairHistManager.h index 2e431427f9c..31efc4779a4 100644 --- a/PWGCF/Femto/Core/pairHistManager.h +++ b/PWGCF/Femto/Core/pairHistManager.h @@ -77,6 +77,8 @@ enum PairHist { kKstarVsMtVsMass1VsMass2VsPt1VsPt2, kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, + // dalitz plots + kDalitz, // between a track and pos/neg daughter of another particle // mc kTrueKstarVsKstar, kTrueKtVsKt, @@ -108,6 +110,7 @@ struct ConfMixing : o2::framework::ConfigurableGroup { struct ConfPairBinning : o2::framework::ConfigurableGroup { std::string prefix = std::string("PairBinning"); + o2::framework::Configurable usePdgMass{"usePdgMass", true, "Use PDF masses for 4-vectors. If false, use reconstructed mass (if available)"}; o2::framework::Configurable plot1D{"plot1D", true, "Enable 1D histograms"}; o2::framework::Configurable plot2D{"plot2D", true, "Enable 2D histograms"}; o2::framework::Configurable plotKstarVsMtVsMult{"plotKstarVsMtVsMult", false, "Enable 3D histogram (Kstar Vs Mt Vs Mult)"}; @@ -121,6 +124,7 @@ struct ConfPairBinning : o2::framework::ConfigurableGroup { o2::framework::Configurable plotKstarVsMtVsMass1VsMass2VsPt1VsPt2{"plotKstarVsMtVsMass1VsMass2VsPt1VsPt2", false, "Enable 6D histogram (Kstar Vs Mt Vs Pt1 Vs Pt2 Vs Mass1 Vs Mass2)"}; o2::framework::Configurable plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult{"plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult", false, "Enable 7D histogram (Kstar Vs Mt Vs Pt1 Vs Pt2 Vs Mass1 Vs Mass2 Vs Mult)"}; o2::framework::Configurable plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent{"plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent", false, "Enable 8D histogram (Kstar Vs Mt Vs Pt1 Vs Pt2 Vs Mass1 Vs Mass2 Vs Mult Vs Cent)"}; + o2::framework::Configurable plotDalitz{"plotDalitz", false, "Enable dalitz plot"}; o2::framework::ConfigurableAxis kstar{"kstar", {{600, 0, 6}}, "kstar"}; o2::framework::ConfigurableAxis kt{"kt", {{600, 0, 6}}, "kt"}; o2::framework::ConfigurableAxis mt{"mt", {{500, 0.8, 5.8}}, "mt"}; @@ -130,6 +134,9 @@ struct ConfPairBinning : o2::framework::ConfigurableGroup { o2::framework::ConfigurableAxis pt2{"pt2", {{100, 0, 6}}, "Pt binning for particle 2"}; o2::framework::ConfigurableAxis mass1{"mass1", {{100, 0, 2}}, "Mass binning for particle 1 (if particle has mass getter, otherwise PDG mass)"}; o2::framework::ConfigurableAxis mass2{"mass2", {{100, 0, 2}}, "Mass binning for particle 2 (if particle has mass getter, otherwise PDG mass)"}; + o2::framework::ConfigurableAxis dalitzMtot{"dalitzMtot", {{100, 0, 10}}, "Total invariant mass squared binning in darlitz plot"}; + o2::framework::ConfigurableAxis dalitzM12{"dalitzM12", {{100, 0, 10}}, "Mass12 binning of darlitz plot"}; + o2::framework::ConfigurableAxis dalitzM13{"dalitzM13", {{100, 0, 10}}, "Mass13 binning of darlitz plot"}; o2::framework::Configurable transverseMassType{"transverseMassType", static_cast(modes::TransverseMassType::kAveragePdgMass), "Type of transverse mass (0-> Average Pdg Mass, 1-> Reduced Pdg Mass, 2-> Mt from combined 4 vector)"}; }; @@ -170,20 +177,21 @@ constexpr std::array, kPairHistogramLast> {kKstarVsMass2, o2::framework::kTH2F, "hKstarVsMass2", "k* vs m_{2}; k* (GeV/#it{c}); m_{2} (GeV/#it{c}^{2})"}, {kMass1VsMass2, o2::framework::kTH2F, "hMass1VsMass2", "m_{1} vs m_{2}; m_{1} (GeV/#it{c}^{2}); m_{2} (GeV/#it{c}^{2})"}, // n-D - {kKstarVsMtVsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMult", "k* vs m_{T} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); Multiplicity"}, - {kKstarVsMtVsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMultVsCent", "k* vs m_{T} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); Multiplicity; Centrality (%)"}, + {kKstarVsMtVsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMult", "k* vs m_{T} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); Multiplicity;"}, + {kKstarVsMtVsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMultVsCent", "k* vs m_{T} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); Multiplicity; Centrality (%);"}, // n-D with pt - {kKstarVsMtVsPt1VsPt2, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2", "k* vs m_{T} vs p_{T,1} vs p_{T,2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c})"}, - {kKstarVsMtVsPt1VsPt2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2VsMult", "k* vs m_{T} vs p_{T,1} vs p_{T,2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity"}, - {kKstarVsMtVsPt1VsPt2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2VsMultVsCent", "k* vs m_{T} vs p_{T,1} vs p_{T,2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity; Centrality"}, + {kKstarVsMtVsPt1VsPt2, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2", "k* vs m_{T} vs p_{T,1} vs p_{T,2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c});"}, + {kKstarVsMtVsPt1VsPt2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2VsMult", "k* vs m_{T} vs p_{T,1} vs p_{T,2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity;"}, + {kKstarVsMtVsPt1VsPt2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsPt1VsPt2VsMultVsCent", "k* vs m_{T} vs p_{T,1} vs p_{T,2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity; Centrality;"}, // n-D with mass - {kKstarVsMtVsMass1VsMass2, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2", "k* vs m_{T} vs m_{1} vs m_{2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2})"}, - {kKstarVsMtVsMass1VsMass2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsMult", "k* vs m_{T} vs m_{1} vs m_{2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); Multiplicity"}, - {kKstarVsMtVsMass1VsMass2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsMultVsCent", "k* vs m_{T} vs m_{1} vs m_{2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); Multiplicity; Centrality (%)"}, + {kKstarVsMtVsMass1VsMass2, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2", "k* vs m_{T} vs m_{1} vs m_{2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2});"}, + {kKstarVsMtVsMass1VsMass2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsMult", "k* vs m_{T} vs m_{1} vs m_{2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); Multiplicity;"}, + {kKstarVsMtVsMass1VsMass2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsMultVsCent", "k* vs m_{T} vs m_{1} vs m_{2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); Multiplicity; Centrality (%);"}, // n-D with pt and mass - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c})"}, - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity"}, - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity; Centrality (%)"}, + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2}; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c});"}, + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2} vs multiplicity; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity;"}, + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, o2::framework::kTHnSparseF, "hKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent", "k* vs m_{T} vs m_{1} vs m_{2} vs p_{T,1} vs p_{T,2} vs multiplicity vs centrality; k* (GeV/#it{c}); m_{T} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); m_{1} (GeV/#it{c}^{2}); p_{T,1} (GeV/#it{c}); p_{T,2} (GeV/#it{c}); Multiplicity; Centrality (%);"}, + {kDalitz, o2::framework::kTHnSparseF, "hDalitz", "Dalitz plot; k* (GeV/#it{c}); m^{2}_{123} (GeV/#it{c}^{2})^{2}; m^{2}_{12} (GeV/#it{c}^{2})^{2}; m^{2}_{13} (GeV/#it{c}^{2})^{2};"}, {kTrueKstarVsKstar, o2::framework::kTH2F, "hTrueKstarVsKstar", "k*_{True} vs k*; k*_{True} (GeV/#it{c}); k* (GeV/#it{c})"}, {kTrueKtVsKt, o2::framework::kTH2F, "hTrueKtVsKt", "k_{T,True} vs k_{T}; k_{T,True} (GeV/#it{c}); k_{T} (GeV/#it{c})"}, {kTrueMtVsMt, o2::framework::kTH2F, "hTrueMtVsMt", "m_{T,True} vs m_{T}; m_{T,True} (GeV/#it{c}^{2}); m_{T} (GeV/#it{c}^{2})"}, @@ -191,35 +199,36 @@ constexpr std::array, kPairHistogramLast> {kTrueCentVsCent, o2::framework::kTH2F, "hTrueCentVsCent", "Centrality_{True} vs Centrality; Centrality_{True} (%); Centrality (%)"}, }}; -#define PAIR_HIST_ANALYSIS_MAP(conf) \ - {kKstar, {conf.kstar}}, \ - {kKt, {conf.kt}}, \ - {kMt, {conf.mt}}, \ - {kPt1VsPt2, {conf.pt1, conf.pt2}}, \ - {kPt1VsKstar, {conf.pt1, conf.kstar}}, \ - {kPt2VsKstar, {conf.pt2, conf.kstar}}, \ - {kPt1VsKt, {conf.pt1, conf.kt}}, \ - {kPt2VsKt, {conf.pt2, conf.kt}}, \ - {kPt1VsMt, {conf.pt1, conf.mt}}, \ - {kPt2VsMt, {conf.pt2, conf.mt}}, \ - {kKstarVsKt, {conf.kstar, conf.kt}}, \ - {kKstarVsMt, {conf.kstar, conf.mt}}, \ - {kKstarVsMult, {conf.kstar, conf.multiplicity}}, \ - {kKstarVsCent, {conf.kstar, conf.centrality}}, \ - {kKstarVsMass1, {conf.kstar, conf.mass1}}, \ - {kKstarVsMass2, {conf.kstar, conf.mass2}}, \ - {kMass1VsMass2, {conf.mass1, conf.mass2}}, \ - {kKstarVsMtVsMult, {conf.kstar, conf.mt, conf.multiplicity}}, \ - {kKstarVsMtVsMultVsCent, {conf.kstar, conf.mt, conf.multiplicity, conf.centrality}}, \ - {kKstarVsMtVsPt1VsPt2, {conf.kstar, conf.mt, conf.pt1, conf.pt2}}, \ - {kKstarVsMtVsPt1VsPt2VsMult, {conf.kstar, conf.mt, conf.pt1, conf.pt2, conf.multiplicity}}, \ - {kKstarVsMtVsPt1VsPt2VsMultVsCent, {conf.kstar, conf.mt, conf.pt1, conf.pt2, conf.multiplicity, conf.centrality}}, \ - {kKstarVsMtVsMass1VsMass2, {conf.kstar, conf.mt, conf.mass1, conf.mass2}}, \ - {kKstarVsMtVsMass1VsMass2VsMult, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.multiplicity}}, \ - {kKstarVsMtVsMass1VsMass2VsMultVsCent, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.multiplicity, conf.centrality}}, \ - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2}}, \ - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2, conf.multiplicity}}, \ - {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2, conf.multiplicity, conf.centrality}}, +#define PAIR_HIST_ANALYSIS_MAP(conf) \ + {kKstar, {conf.kstar}}, \ + {kKt, {conf.kt}}, \ + {kMt, {conf.mt}}, \ + {kPt1VsPt2, {conf.pt1, conf.pt2}}, \ + {kPt1VsKstar, {conf.pt1, conf.kstar}}, \ + {kPt2VsKstar, {conf.pt2, conf.kstar}}, \ + {kPt1VsKt, {conf.pt1, conf.kt}}, \ + {kPt2VsKt, {conf.pt2, conf.kt}}, \ + {kPt1VsMt, {conf.pt1, conf.mt}}, \ + {kPt2VsMt, {conf.pt2, conf.mt}}, \ + {kKstarVsKt, {conf.kstar, conf.kt}}, \ + {kKstarVsMt, {conf.kstar, conf.mt}}, \ + {kKstarVsMult, {conf.kstar, conf.multiplicity}}, \ + {kKstarVsCent, {conf.kstar, conf.centrality}}, \ + {kKstarVsMass1, {conf.kstar, conf.mass1}}, \ + {kKstarVsMass2, {conf.kstar, conf.mass2}}, \ + {kMass1VsMass2, {conf.mass1, conf.mass2}}, \ + {kKstarVsMtVsMult, {conf.kstar, conf.mt, conf.multiplicity}}, \ + {kKstarVsMtVsMultVsCent, {conf.kstar, conf.mt, conf.multiplicity, conf.centrality}}, \ + {kKstarVsMtVsPt1VsPt2, {conf.kstar, conf.mt, conf.pt1, conf.pt2}}, \ + {kKstarVsMtVsPt1VsPt2VsMult, {conf.kstar, conf.mt, conf.pt1, conf.pt2, conf.multiplicity}}, \ + {kKstarVsMtVsPt1VsPt2VsMultVsCent, {conf.kstar, conf.mt, conf.pt1, conf.pt2, conf.multiplicity, conf.centrality}}, \ + {kKstarVsMtVsMass1VsMass2, {conf.kstar, conf.mt, conf.mass1, conf.mass2}}, \ + {kKstarVsMtVsMass1VsMass2VsMult, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.multiplicity}}, \ + {kKstarVsMtVsMass1VsMass2VsMultVsCent, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.multiplicity, conf.centrality}}, \ + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2}}, \ + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2, conf.multiplicity}}, \ + {kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, {conf.kstar, conf.mt, conf.mass1, conf.mass2, conf.pt1, conf.pt2, conf.multiplicity, conf.centrality}}, \ + {kDalitz, {conf.kstar, conf.dalitzMtot, conf.dalitzM12, conf.dalitzM13}}, #define PAIR_HIST_MC_MAP(conf) \ {kTrueKstarVsKstar, {conf.kstar, conf.kstar}}, \ @@ -285,6 +294,8 @@ class PairHistManager { mHistogramRegistry = registry; + mUsePdgMass = ConfPairBinning.usePdgMass.value; + // flags for histograms mPlot1d = ConfPairBinning.plot1D.value; mPlot2d = ConfPairBinning.plot2D.value; @@ -303,6 +314,8 @@ class PairHistManager mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult = ConfPairBinning.plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult.value; mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent = ConfPairBinning.plotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent.value; + mPlotDalitz = ConfPairBinning.plotDalitz.value; + // transverse mass type mMtType = static_cast(ConfPairBinning.transverseMassType.value); @@ -325,8 +338,18 @@ class PairHistManager void setMass(int PdgParticle1, int PdgParticle2) { - mPdgMass1 = o2::analysis::femto::utils::getMass(PdgParticle1); - mPdgMass2 = o2::analysis::femto::utils::getMass(PdgParticle2); + mPdgMass1 = utils::getPdgMass(PdgParticle1); + mPdgMass2 = utils::getPdgMass(PdgParticle2); + } + + void setMass(int PdgParticle1, int PdgPosDauParticle1, int PdgNegDauParticle1, int PdgParticle2, int PdgPosDauParticle2, int PdgNegDauParticle2) + { + mPdgMass1 = utils::getPdgMass(PdgParticle1); + mPdgMassPosDau1 = utils::getPdgMass(PdgPosDauParticle1); + mPdgMassNegDau1 = utils::getPdgMass(PdgNegDauParticle1); + mPdgMass2 = utils::getPdgMass(PdgParticle2); + mPdgMassPosDau2 = utils::getPdgMass(PdgPosDauParticle2); + mPdgMassNegDau2 = utils::getPdgMass(PdgNegDauParticle2); } void setCharge(int chargeAbsParticle1, int chargeAbsParticle2) { @@ -335,13 +358,37 @@ class PairHistManager mAbsCharge2 = std::abs(chargeAbsParticle2); } - template - void setPair(const T1& particle1, const T2& particle2) + template + void setPair(T1 const& particle1, T2 const& particle2, T3 const& trackTable) { + // if one of the particles has a mass getter (like lambda), we cache the value for the filling later + // otherwise set it to the pdg mass + if constexpr (utils::HasMass) { + mRecoMass1 = particle1.mass(); + } else { + mRecoMass1 = mPdgMass1; + } + if constexpr (utils::HasMass) { + mRecoMass2 = particle2.mass(); + } else { + mRecoMass2 = mPdgMass2; + } + + // get mass for 4-vectors + double mass1 = 0.f; + double mass2 = 0.f; + if (mUsePdgMass) { + mass1 = mPdgMass1; + mass2 = mPdgMass2; + } else { + mass1 = mRecoMass1; + mass2 = mRecoMass2; + } + // pt in track table is calculated from 1/signedPt from the original track table // in case of He with Z=2, we have to rescale the pt with the absolute charge - mParticle1 = ROOT::Math::PtEtaPhiMVector(mAbsCharge1 * particle1.pt(), particle1.eta(), particle1.phi(), mPdgMass1); - mParticle2 = ROOT::Math::PtEtaPhiMVector(mAbsCharge2 * particle2.pt(), particle2.eta(), particle2.phi(), mPdgMass2); + mParticle1 = ROOT::Math::PtEtaPhiMVector(mAbsCharge1 * particle1.pt(), particle1.eta(), particle1.phi(), mass1); + mParticle2 = ROOT::Math::PtEtaPhiMVector(mAbsCharge2 * particle2.pt(), particle2.eta(), particle2.phi(), mass2); // set kT mKt = getKt(mParticle1, mParticle2); @@ -352,36 +399,37 @@ class PairHistManager // set kstar mKstar = getKstar(mParticle1, mParticle2); - // if one of the particles has a mass getter (like lambda), we cache the value for the filling later - // otherwise we continue to use the pdg mass - mMass1 = mPdgMass1; - if constexpr (utils::HasMass) { - mMass1 = particle1.mass(); - } - mMass2 = mPdgMass2; - if constexpr (utils::HasMass) { - mMass2 = particle2.mass(); + if (mPlotDalitz) { + if constexpr (modes::isEqual(particleType1, modes::Particle::kTrack) && modes::isEqual(particleType2, modes::Particle::kV0)) { + auto posDaughter = trackTable.rawIteratorAt(particle2.posDauId() - trackTable.offset()); + auto negDaughter = trackTable.rawIteratorAt(particle2.negDauId() - trackTable.offset()); + ROOT::Math::PtEtaPhiMVector posDau4v = ROOT::Math::PtEtaPhiMVector(posDaughter.pt(), posDaughter.eta(), posDaughter.phi(), mPdgMassPosDau2); + ROOT::Math::PtEtaPhiMVector negDau4v = ROOT::Math::PtEtaPhiMVector(negDaughter.pt(), negDaughter.eta(), negDaughter.phi(), mPdgMassNegDau2); + mMassTot2 = (mParticle1 + posDau4v + negDau4v).M2(); + mMass12 = (mParticle1 + posDau4v).M2(); + mMass13 = (mParticle1 + negDau4v).M2(); + } } } - template - void setPair(const T1& particle1, const T2& particle2, const T3& col) + template + void setPair(T1 const& particle1, T2 const& particle2, T3 const& trackTable, T4 const& col) { - setPair(particle1, particle2); + setPair(particle1, particle2, trackTable); mMult = col.mult(); mCent = col.cent(); } - template - void setPair(const T1& particle1, const T2& particle2, const T3& col1, const T4& col2) + template + void setPair(T1 const& particle1, T2 const& particle2, T3 const& trackTable, T4 const& col1, T5 const& col2) { - setPair(particle1, particle2); + setPair(particle1, particle2, trackTable); mMult = 0.5f * (col1.mult() + col2.mult()); // if mixing with multiplicity, should be in the same mixing bin mCent = 0.5f * (col1.cent() + col2.cent()); // if mixing with centrality, should be in the same mixing bin } template - void setPairMc(const T1& particle1, const T2& particle2, const T3& /*mcParticles*/) + void setPairMc(T1 const& particle1, T2 const& particle2, const T3& /*mcParticles*/) { if (!particle1.has_fMcParticle() || !particle2.has_fMcParticle()) { mHasMcPair = false; @@ -404,33 +452,33 @@ class PairHistManager mTrueKstar = getKstar(mTrueParticle1, mTrueParticle2); } - template - void setPairMc(const T1& particle1, const T2& particle2, const T3& mcParticles, const T4& col, const T5& /*mcCols*/) + template + void setPairMc(T1 const& particle1, T2 const& particle2, T3 const& trackTable, T4 const& mcParticles, T5 const& col, T6 const& /*mcCols*/) { - setPair(particle1, particle2, col); + setPair(particle1, particle2, trackTable, col); setPairMc(particle1, particle2, mcParticles); if (!col.has_fMcCol()) { mHasMcCol = false; return; } mHasMcCol = true; - auto mcCol = col.template fMcCol_as(); + auto mcCol = col.template fMcCol_as(); mTrueMult = mcCol.mult(); mTrueCent = mcCol.cent(); } - template - void setPairMc(const T1& particle1, const T2& particle2, const T3& mcParticles, const T4& col1, const T5& col2, const T6& /*mcCols*/) + template + void setPairMc(T1 const& particle1, T2 const& particle2, T3 const& trackTable, T4 const& mcParticles, T5 const& col1, T6 const& col2, T7 const& /*mcCols*/) { - setPair(particle1, particle2, col1, col2); + setPair(particle1, particle2, trackTable, col1, col2); setPairMc(particle1, particle2, mcParticles); if (!col1.has_fMcCol() || !col2.has_fMcCol()) { mHasMcCol = false; return; } mHasMcCol = true; - auto mcCol1 = col1.template fMcCol_as(); - auto mcCol2 = col2.template fMcCol_as(); + auto mcCol1 = col1.template fMcCol_as(); + auto mcCol2 = col2.template fMcCol_as(); mTrueMult = 0.5f * (mcCol1.mult() + mcCol2.mult()); mTrueCent = 0.5f * (mcCol1.cent() + mcCol2.cent()); } @@ -521,6 +569,9 @@ class PairHistManager if (mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent) { mHistogramRegistry->add(analysisDir + getHistNameV2(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, HistTable), getHistDesc(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, HistTable), getHistType(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, HistTable), {Specs.at(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent)}); } + if (mPlotDalitz) { + mHistogramRegistry->add(analysisDir + getHistNameV2(kDalitz, HistTable), getHistDesc(kDalitz, HistTable), getHistType(kDalitz, HistTable), {Specs.at(kDalitz)}); + } } void initMc(std::map> const& Specs) @@ -553,9 +604,9 @@ class PairHistManager mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMult, HistTable)), mKstar, mMult); mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsCent, HistTable)), mKstar, mCent); - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMass1, HistTable)), mKstar, mMass1); - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMass2, HistTable)), mKstar, mMass2); - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kMass1VsMass2, HistTable)), mMass1, mMass2); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMass1, HistTable)), mKstar, mRecoMass1); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMass2, HistTable)), mKstar, mRecoMass2); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kMass1VsMass2, HistTable)), mRecoMass1, mRecoMass2); } // n-D histograms are only filled if enabled @@ -577,22 +628,25 @@ class PairHistManager mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsPt1VsPt2VsMultVsCent, HistTable)), mKstar, mMt, mParticle1.Pt(), mParticle2.Pt(), mMult, mCent); } if (mPlotKstarVsMtVsMass1VsMass2) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2, HistTable)), mKstar, mMt, mMass1, mMass2); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2); } if (mPlotKstarVsMtVsMass1VsMass2VsMult) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsMult, HistTable)), mKstar, mMt, mMass1, mMass2, mMult); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsMult, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2, mMult); } if (mPlotKstarVsMtVsMass1VsMass2VsMultVsCent) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsMultVsCent, HistTable)), mKstar, mMt, mMass1, mMass2, mMult, mCent); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsMultVsCent, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2, mMult, mCent); } if (mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2, HistTable)), mKstar, mMt, mMass1, mMass2, mParticle1.Pt(), mParticle2.Pt()); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2, mParticle1.Pt(), mParticle2.Pt()); } if (mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, HistTable)), mKstar, mMt, mMass1, mMass2, mParticle1.Pt(), mParticle2.Pt(), mMult); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2, mParticle1.Pt(), mParticle2.Pt(), mMult); } if (mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent) { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, HistTable)), mKstar, mMt, mMass1, mMass2, mParticle1.Pt(), mParticle2.Pt(), mMult, mCent); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent, HistTable)), mKstar, mMt, mRecoMass1, mRecoMass2, mParticle1.Pt(), mParticle2.Pt(), mMult, mCent); + } + if (mPlotDalitz) { + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kDalitz, HistTable)), mKstar, mMassTot2, mMass12, mMass13); } } @@ -654,8 +708,13 @@ class PairHistManager } o2::framework::HistogramRegistry* mHistogramRegistry = nullptr; + bool mUsePdgMass = true; double mPdgMass1 = 0.; + double mPdgMassPosDau1 = 0; + double mPdgMassNegDau1 = 0; double mPdgMass2 = 0.; + double mPdgMassPosDau2 = 0; + double mPdgMassNegDau2 = 0; modes::TransverseMassType mMtType = modes::TransverseMassType::kAveragePdgMass; @@ -663,13 +722,16 @@ class PairHistManager int mAbsCharge2 = 1; ROOT::Math::PtEtaPhiMVector mParticle1{}; ROOT::Math::PtEtaPhiMVector mParticle2{}; - float mMass1 = 0.f; - float mMass2 = 0.f; + float mRecoMass1 = 0.f; + float mRecoMass2 = 0.f; float mKstar = 0.f; float mKt = 0.f; float mMt = 0.f; float mMult = 0.f; float mCent = 0.f; + double mMass12 = 0.; + double mMass13 = 0.; + double mMassTot2 = 0.; // mc ROOT::Math::PtEtaPhiMVector mTrueParticle1{}; @@ -708,6 +770,8 @@ class PairHistManager bool mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2 = false; bool mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMult = false; bool mPlotKstarVsMtVsMass1VsMass2VsPt1VsPt2VsMultVsCent = false; + + bool mPlotDalitz = false; }; }; // namespace pairhistmanager diff --git a/PWGCF/Femto/Core/pairProcessHelpers.h b/PWGCF/Femto/Core/pairProcessHelpers.h index 50d00196341..9279921b00b 100644 --- a/PWGCF/Femto/Core/pairProcessHelpers.h +++ b/PWGCF/Femto/Core/pairProcessHelpers.h @@ -67,13 +67,13 @@ void processSameEvent(T1 const& SliceParticle, // Randomize pair order if enabled switch (pairOrder) { case kOrder12: - PairHistManager.setPair(p1, p2, Collision); + PairHistManager.setPair(p1, p2, TrackTable, Collision); break; case kOrder21: - PairHistManager.setPair(p2, p1, Collision); + PairHistManager.setPair(p2, p1, TrackTable, Collision); break; default: - PairHistManager.setPair(p1, p2, Collision); + PairHistManager.setPair(p1, p2, TrackTable, Collision); } // fill deta-dphi histograms with kstar cutoff CprManager.fill(PairHistManager.getKstar()); @@ -136,13 +136,13 @@ void processSameEvent(T1 const& SliceParticle, // Randomize pair order if enabled switch (pairOrder) { case kOrder12: - PairHistManager.setPairMc(p1, p2, mcParticles, Collision, mcCollisions); + PairHistManager.setPairMc(p1, p2, TrackTable, mcParticles, Collision, mcCollisions); break; case kOrder21: - PairHistManager.setPairMc(p2, p1, mcParticles, Collision, mcCollisions); + PairHistManager.setPairMc(p2, p1, TrackTable, mcParticles, Collision, mcCollisions); break; default: - PairHistManager.setPairMc(p1, p2, mcParticles, Collision, mcCollisions); + PairHistManager.setPairMc(p1, p2, TrackTable, mcParticles, Collision, mcCollisions); } // fill deta-dphi histograms with kstar cutoff CprManager.fill(PairHistManager.getKstar()); @@ -191,7 +191,7 @@ void processSameEvent(T1 const& SliceParticle1, if (CprManager.isClosePair()) { continue; } - PairHistManager.setPair(p1, p2, Collision); + PairHistManager.setPair(p1, p2, TrackTable, Collision); CprManager.fill(PairHistManager.getKstar()); if (PairHistManager.checkPairCuts()) { PairHistManager.template fill(); @@ -260,7 +260,7 @@ void processSameEvent(T1 const& SliceParticle1, if (CprManager.isClosePair()) { continue; } - PairHistManager.setPairMc(p1, p2, mcParticles, Collision, mcCollisions); + PairHistManager.setPairMc(p1, p2, TrackTable, mcParticles, Collision, mcCollisions); CprManager.fill(PairHistManager.getKstar()); if (PairHistManager.checkPairCuts()) { PairHistManager.template fill(); @@ -311,7 +311,7 @@ void processMixedEvent(T1 const& Collisions, if (CprManager.isClosePair()) { continue; } - PairHistManager.setPair(p1, p2, collision1, collision2); + PairHistManager.setPair(p1, p2, TrackTable, collision1, collision2); CprManager.fill(PairHistManager.getKstar()); if (PairHistManager.checkPairCuts()) { PairHistManager.template fill(); @@ -380,7 +380,7 @@ void processMixedEvent(T1 const& Collisions, if (CprManager.isClosePair()) { continue; } - PairHistManager.setPairMc(p1, p2, mcParticles, collision1, collision2, mcCollisions); + PairHistManager.setPairMc(p1, p2, TrackTable, mcParticles, collision1, collision2, mcCollisions); CprManager.fill(PairHistManager.getKstar()); if (PairHistManager.checkPairCuts()) { PairHistManager.template fill(); diff --git a/PWGCF/Femto/Core/trackHistManager.h b/PWGCF/Femto/Core/trackHistManager.h index 12f69af380b..dbff8bb2dd1 100644 --- a/PWGCF/Femto/Core/trackHistManager.h +++ b/PWGCF/Femto/Core/trackHistManager.h @@ -794,7 +794,7 @@ class TrackHistManager if constexpr (utils::HasMass) { mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kMass, HistTable)), track.mass()); } else { - mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kMass, HistTable)), utils::getMass(mPdgCode)); + mHistogramRegistry->fill(HIST(prefix) + HIST(AnalysisDir) + HIST(getHistName(kMass, HistTable)), utils::getPdgMass(mPdgCode)); } } diff --git a/PWGCF/Femto/Core/tripletHistManager.h b/PWGCF/Femto/Core/tripletHistManager.h index 8eca807c5aa..672e9f68055 100644 --- a/PWGCF/Femto/Core/tripletHistManager.h +++ b/PWGCF/Femto/Core/tripletHistManager.h @@ -252,9 +252,9 @@ class TripletHistManager void setMass(int PdgParticle1, int PdgParticle2, int PdgParticle3) { - mPdgMass1 = o2::analysis::femto::utils::getMass(PdgParticle1); - mPdgMass2 = o2::analysis::femto::utils::getMass(PdgParticle2); - mPdgMass3 = o2::analysis::femto::utils::getMass(PdgParticle3); + mPdgMass1 = utils::getPdgMass(PdgParticle1); + mPdgMass2 = utils::getPdgMass(PdgParticle2); + mPdgMass3 = utils::getPdgMass(PdgParticle3); } void setCharge(int chargeAbsParticle1, int chargeAbsParticle2, int chargeAbsParticle3) {