From 68155fd4e8e4498664de2fc0e28db5dfc09f923a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Nicol=C3=B2=20Jacazio?= Date: Sat, 21 Sep 2024 11:47:43 +0200 Subject: [PATCH 1/2] Update --- PWGLF/DataModel/spectraTOF.h | 60 --------- PWGLF/Tasks/Nuspex/spectraTOF.cxx | 209 +++++++++++------------------- 2 files changed, 75 insertions(+), 194 deletions(-) diff --git a/PWGLF/DataModel/spectraTOF.h b/PWGLF/DataModel/spectraTOF.h index 884003b3951..97d5909332e 100644 --- a/PWGLF/DataModel/spectraTOF.h +++ b/PWGLF/DataModel/spectraTOF.h @@ -88,24 +88,12 @@ static constexpr std::string_view hdcaxy[NpCharge] = {"dcaxy/pos/el", "dcaxy/pos "dcaxy/neg/el", "dcaxy/neg/mu", "dcaxy/neg/pi", "dcaxy/neg/ka", "dcaxy/neg/pr", "dcaxy/neg/de", "dcaxy/neg/tr", "dcaxy/neg/he", "dcaxy/neg/al"}; -static constexpr std::string_view hdcaxytot[NpCharge] = {"dcaxytot/pos/el", "dcaxytot/pos/mu", "dcaxytot/pos/pi", - "dcaxytot/pos/ka", "dcaxytot/pos/pr", "dcaxytot/pos/de", - "dcaxytot/pos/tr", "dcaxytot/pos/he", "dcaxytot/pos/al", - "dcaxytot/neg/el", "dcaxytot/neg/mu", "dcaxytot/neg/pi", - "dcaxytot/neg/ka", "dcaxytot/neg/pr", "dcaxytot/neg/de", - "dcaxytot/neg/tr", "dcaxytot/neg/he", "dcaxytot/neg/al"}; static constexpr std::string_view hdcaz[NpCharge] = {"dcaz/pos/el", "dcaz/pos/mu", "dcaz/pos/pi", "dcaz/pos/ka", "dcaz/pos/pr", "dcaz/pos/de", "dcaz/pos/tr", "dcaz/pos/he", "dcaz/pos/al", "dcaz/neg/el", "dcaz/neg/mu", "dcaz/neg/pi", "dcaz/neg/ka", "dcaz/neg/pr", "dcaz/neg/de", "dcaz/neg/tr", "dcaz/neg/he", "dcaz/neg/al"}; -static constexpr std::string_view hdcaztot[NpCharge] = {"dcaztot/pos/el", "dcaztot/pos/mu", "dcaztot/pos/pi", - "dcaztot/pos/ka", "dcaztot/pos/pr", "dcaztot/pos/de", - "dcaztot/pos/tr", "dcaztot/pos/he", "dcaztot/pos/al", - "dcaztot/neg/el", "dcaztot/neg/mu", "dcaztot/neg/pi", - "dcaztot/neg/ka", "dcaztot/neg/pr", "dcaztot/neg/de", - "dcaztot/neg/tr", "dcaztot/neg/he", "dcaztot/neg/al"}; static constexpr std::string_view hdcaxyphi[NpCharge] = {"dcaxyphi/pos/el", "dcaxyphi/pos/mu", "dcaxyphi/pos/pi", "dcaxyphi/pos/ka", "dcaxyphi/pos/pr", "dcaxyphi/pos/de", "dcaxyphi/pos/tr", "dcaxyphi/pos/he", "dcaxyphi/pos/al", @@ -235,18 +223,6 @@ static constexpr std::string_view hdcaxyprm[NpCharge] = {"dcaxyprm/pos/el", "dca "dcaxyprm/neg/el", "dcaxyprm/neg/mu", "dcaxyprm/neg/pi", "dcaxyprm/neg/ka", "dcaxyprm/neg/pr", "dcaxyprm/neg/de", "dcaxyprm/neg/tr", "dcaxyprm/neg/he", "dcaxyprm/neg/al"}; -static constexpr std::string_view hdcaxyprm2[NpCharge] = {"dcaxyprm2/pos/el", "dcaxyprm2/pos/mu", "dcaxyprm2/pos/pi", - "dcaxyprm2/pos/ka", "dcaxyprm2/pos/pr", "dcaxyprm2/pos/de", - "dcaxyprm2/pos/tr", "dcaxyprm2/pos/he", "dcaxyprm2/pos/al", - "dcaxyprm2/neg/el", "dcaxyprm2/neg/mu", "dcaxyprm2/neg/pi", - "dcaxyprm2/neg/ka", "dcaxyprm2/neg/pr", "dcaxyprm2/neg/de", - "dcaxyprm2/neg/tr", "dcaxyprm2/neg/he", "dcaxyprm2/neg/al"}; -static constexpr std::string_view hdcaxyD0[NpCharge] = {"dcaxyD0/pos/el", "dcaxyD0/pos/mu", "dcaxyD0/pos/pi", - "dcaxyD0/pos/ka", "dcaxyD0/pos/pr", "dcaxyD0/pos/de", - "dcaxyD0/pos/tr", "dcaxyD0/pos/he", "dcaxyD0/pos/al", - "dcaxyD0/neg/el", "dcaxyD0/neg/mu", "dcaxyD0/neg/pi", - "dcaxyD0/neg/ka", "dcaxyD0/neg/pr", "dcaxyD0/neg/de", - "dcaxyD0/neg/tr", "dcaxyD0/neg/he", "dcaxyD0/neg/al"}; static constexpr std::string_view hdcaxyprmgoodevs[NpCharge] = {"dcaxyprmgoodevs/pos/el", "dcaxyprmgoodevs/pos/mu", "dcaxyprmgoodevs/pos/pi", "dcaxyprmgoodevs/pos/ka", "dcaxyprmgoodevs/pos/pr", "dcaxyprmgoodevs/pos/de", "dcaxyprmgoodevs/pos/tr", "dcaxyprmgoodevs/pos/he", "dcaxyprmgoodevs/pos/al", @@ -259,18 +235,6 @@ static constexpr std::string_view hdcazprm[NpCharge] = {"dcazprm/pos/el", "dcazp "dcazprm/neg/el", "dcazprm/neg/mu", "dcazprm/neg/pi", "dcazprm/neg/ka", "dcazprm/neg/pr", "dcazprm/neg/de", "dcazprm/neg/tr", "dcazprm/neg/he", "dcazprm/neg/al"}; -static constexpr std::string_view hdcazprm2[NpCharge] = {"dcazprm2/pos/el", "dcazprm2/pos/mu", "dcazprm2/pos/pi", - "dcazprm2/pos/ka", "dcazprm2/pos/pr", "dcazprm2/pos/de", - "dcazprm2/pos/tr", "dcazprm2/pos/he", "dcazprm2/pos/al", - "dcazprm2/neg/el", "dcazprm2/neg/mu", "dcazprm2/neg/pi", - "dcazprm2/neg/ka", "dcazprm2/neg/pr", "dcazprm2/neg/de", - "dcazprm2/neg/tr", "dcazprm2/neg/he", "dcazprm2/neg/al"}; -static constexpr std::string_view hdcazD0[NpCharge] = {"dcazD0/pos/el", "dcazD0/pos/mu", "dcazD0/pos/pi", - "dcazD0/pos/ka", "dcazD0/pos/pr", "dcazD0/pos/de", - "dcazD0/pos/tr", "dcazD0/pos/he", "dcazD0/pos/al", - "dcazD0/neg/el", "dcazD0/neg/mu", "dcazD0/neg/pi", - "dcazD0/neg/ka", "dcazD0/neg/pr", "dcazD0/neg/de", - "dcazD0/neg/tr", "dcazD0/neg/he", "dcazD0/neg/al"}; static constexpr std::string_view hdcazprmgoodevs[NpCharge] = {"dcazprmgoodevs/pos/el", "dcazprmgoodevs/pos/mu", "dcazprmgoodevs/pos/pi", "dcazprmgoodevs/pos/ka", "dcazprmgoodevs/pos/pr", "dcazprmgoodevs/pos/de", "dcazprmgoodevs/pos/tr", "dcazprmgoodevs/pos/he", "dcazprmgoodevs/pos/al", @@ -283,18 +247,6 @@ static constexpr std::string_view hdcaxystr[NpCharge] = {"dcaxystr/pos/el", "dca "dcaxystr/neg/el", "dcaxystr/neg/mu", "dcaxystr/neg/pi", "dcaxystr/neg/ka", "dcaxystr/neg/pr", "dcaxystr/neg/de", "dcaxystr/neg/tr", "dcaxystr/neg/he", "dcaxystr/neg/al"}; -static constexpr std::string_view hdcaxycharm[NpCharge] = {"dcaxycharm/pos/el", "dcaxycharm/pos/mu", "dcaxycharm/pos/pi", - "dcaxycharm/pos/ka", "dcaxycharm/pos/pr", "dcaxycharm/pos/de", - "dcaxycharm/pos/tr", "dcaxycharm/pos/he", "dcaxycharm/pos/al", - "dcaxycharm/neg/el", "dcaxycharm/neg/mu", "dcaxycharm/neg/pi", - "dcaxycharm/neg/ka", "dcaxycharm/neg/pr", "dcaxycharm/neg/de", - "dcaxycharm/neg/tr", "dcaxycharm/neg/he", "dcaxycharm/neg/al"}; -static constexpr std::string_view hdcaxybeauty[NpCharge] = {"dcaxybeauty/pos/el", "dcaxybeauty/pos/mu", "dcaxybeauty/pos/pi", - "dcaxybeauty/pos/ka", "dcaxybeauty/pos/pr", "dcaxybeauty/pos/de", - "dcaxybeauty/pos/tr", "dcaxybeauty/pos/he", "dcaxybeauty/pos/al", - "dcaxybeauty/neg/el", "dcaxybeauty/neg/mu", "dcaxybeauty/neg/pi", - "dcaxybeauty/neg/ka", "dcaxybeauty/neg/pr", "dcaxybeauty/neg/de", - "dcaxybeauty/neg/tr", "dcaxybeauty/neg/he", "dcaxybeauty/neg/al"}; static constexpr std::string_view hdcaxymat[NpCharge] = {"dcaxymat/pos/el", "dcaxymat/pos/mu", "dcaxymat/pos/pi", "dcaxymat/pos/ka", "dcaxymat/pos/pr", "dcaxymat/pos/de", "dcaxymat/pos/tr", "dcaxymat/pos/he", "dcaxymat/pos/al", @@ -307,18 +259,6 @@ static constexpr std::string_view hdcazstr[NpCharge] = {"dcazstr/pos/el", "dcazs "dcazstr/neg/el", "dcazstr/neg/mu", "dcazstr/neg/pi", "dcazstr/neg/ka", "dcazstr/neg/pr", "dcazstr/neg/de", "dcazstr/neg/tr", "dcazstr/neg/he", "dcazstr/neg/al"}; -static constexpr std::string_view hdcazcharm[NpCharge] = {"dcazcharm/pos/el", "dcazcharm/pos/mu", "dcazcharm/pos/pi", - "dcazcharm/pos/ka", "dcazcharm/pos/pr", "dcazcharm/pos/de", - "dcazcharm/pos/tr", "dcazcharm/pos/he", "dcazcharm/pos/al", - "dcazcharm/neg/el", "dcazcharm/neg/mu", "dcazcharm/neg/pi", - "dcazcharm/neg/ka", "dcazcharm/neg/pr", "dcazcharm/neg/de", - "dcazcharm/neg/tr", "dcazcharm/neg/he", "dcazcharm/neg/al"}; -static constexpr std::string_view hdcazbeauty[NpCharge] = {"dcazbeauty/pos/el", "dcazbeauty/pos/mu", "dcazbeauty/pos/pi", - "dcazbeauty/pos/ka", "dcazbeauty/pos/pr", "dcazbeauty/pos/de", - "dcazbeauty/pos/tr", "dcazbeauty/pos/he", "dcazbeauty/pos/al", - "dcazbeauty/neg/el", "dcazbeauty/neg/mu", "dcazbeauty/neg/pi", - "dcazbeauty/neg/ka", "dcazbeauty/neg/pr", "dcazbeauty/neg/de", - "dcazbeauty/neg/tr", "dcazbeauty/neg/he", "dcazbeauty/neg/al"}; static constexpr std::string_view hdcazmat[NpCharge] = {"dcazmat/pos/el", "dcazmat/pos/mu", "dcazmat/pos/pi", "dcazmat/pos/ka", "dcazmat/pos/pr", "dcazmat/pos/de", "dcazmat/pos/tr", "dcazmat/pos/he", "dcazmat/pos/al", diff --git a/PWGLF/Tasks/Nuspex/spectraTOF.cxx b/PWGLF/Tasks/Nuspex/spectraTOF.cxx index 6b4361cf99c..7efc6821ed2 100644 --- a/PWGLF/Tasks/Nuspex/spectraTOF.cxx +++ b/PWGLF/Tasks/Nuspex/spectraTOF.cxx @@ -36,7 +36,7 @@ #include "Framework/O2DatabasePDGPlugin.h" #include "PWGLF/Utils/inelGt.h" #include "PWGLF/DataModel/mcCentrality.h" - +#include "Common/Core/RecoDecay.h" #include "TPDGCode.h" using namespace o2; @@ -44,91 +44,6 @@ using namespace o2::track; using namespace o2::framework; using namespace o2::framework::expressions; -class RecoDecay -{ - public: - /// Default constructor - RecoDecay() = default; - - /// Default destructor - ~RecoDecay() = default; - - // mapping of charm-hadron origin type - enum OriginType { None = 0, - Prompt, - NonPrompt }; - - template - static int getCharmHadronOrigin(const T& particlesMC, - const typename T::iterator& particle, - const bool searchUpToQuark = false) - { - int stage = 0; // mother tree level (just for debugging) - - // vector of vectors with mother indices; each line corresponds to a "stage" - std::vector> arrayIds{}; - std::vector initVec{particle.globalIndex()}; - arrayIds.push_back(initVec); // the first vector contains the index of the original particle - auto PDGParticle = std::abs(particle.pdgCode()); - bool couldBePrompt = false; - if (PDGParticle / 100 == 4 || PDGParticle / 1000 == 4) { - couldBePrompt = true; - } - while (arrayIds[-stage].size() > 0) { - // vector of mother indices for the current stage - std::vector arrayIdsStage{}; - for (auto& iPart : arrayIds[-stage]) { // check all the particles that were the mothers at the previous stage - auto particleMother = particlesMC.rawIteratorAt(iPart - particlesMC.offset()); - if (particleMother.has_mothers()) { - for (auto iMother = particleMother.mothersIds().front(); iMother <= particleMother.mothersIds().back(); ++iMother) { // loop over the mother particles of the analysed particle - if (std::find(arrayIdsStage.begin(), arrayIdsStage.end(), iMother) != arrayIdsStage.end()) { // if a mother is still present in the vector, do not check it again - continue; - } - auto mother = particlesMC.rawIteratorAt(iMother - particlesMC.offset()); - // Check mother's PDG code. - auto PDGParticleIMother = std::abs(mother.pdgCode()); // PDG code of the mother - // printf("getMother: "); - // for (int i = stage; i < 0; i++) // Indent to make the tree look nice. - // printf(" "); - // printf("Stage %d: Mother PDG: %d, Index: %d\n", stage, PDGParticleIMother, iMother); - - if (searchUpToQuark) { - if (PDGParticleIMother == 5) { // b quark - return OriginType::NonPrompt; - } - if (PDGParticleIMother == 4) { // c quark - return OriginType::Prompt; - } - } else { - if ( - (PDGParticleIMother / 100 == 5 || // b mesons - PDGParticleIMother / 1000 == 5) // b baryons - ) { - return OriginType::NonPrompt; - } - if ( - (PDGParticleIMother / 100 == 4 || // c mesons - PDGParticleIMother / 1000 == 4) // c baryons - ) { - couldBePrompt = true; - } - } - // add mother index in the vector for the current stage - arrayIdsStage.push_back(iMother); - } - } - } - // add vector of mother indices for the current stage - arrayIds.push_back(arrayIdsStage); - stage--; - } - if (!searchUpToQuark && couldBePrompt) { // Returns prompt if it's a charm hadron or a charm-hadron daughter. Note: 1) LF decay particles from cases like -> Lc -> p K0S, K0S -> pi pi are marked as prompt. 2) if particles from HF parton showers have to be searched, switch to option "search up to quark" - return OriginType::Prompt; - } - return OriginType::None; - } -}; - // Histograms std::array, NpCharge> hDcaXYZ; std::array, NpCharge> hDcaXYZPrm; @@ -137,6 +52,17 @@ std::array, NpCharge> hDcaXYZMat; std::array, NpCharge> hDcaXYWrongCollisionPrm; std::array, NpCharge> hDcaXYWrongCollisionStr; std::array, NpCharge> hDcaXYWrongCollisionMat; +std::array, NpCharge> hDcaXYMC; // DCA xy in the MC +std::array, NpCharge> hDcaZMC; // DCA z in the MC +std::array, NpCharge> hDcaXYMCPrm; // DCA xy in the MC +std::array, NpCharge> hDcaXYMCD0; // DCA xy in the MC for particles from D0 +std::array, NpCharge> hDcaZMCD0; // DCA z in the MC for particles from D0 +std::array, NpCharge> hDcaXYMCCharm; // DCA xy in the MC for particles from charm +std::array, NpCharge> hdcaZMCCharm; // DCA z in the MC for particles from charm +std::array, NpCharge> hDcaXYMCBeauty; // DCA xy in the MC for particles from beauty +std::array, NpCharge> hDcaZMCBeauty; // DCA z in the MC for particles from beauty +std::array, NpCharge> hDcaXYMCNotHF; // DCA xy in the MC for particles from not a HF +std::array, NpCharge> hDcaZMCNotHF; // DCA z in the MC for particles from not a HF // Spectra task struct tofSpectra { @@ -168,6 +94,7 @@ struct tofSpectra { Configurable enableTPCTOFHistograms{"enableTPCTOFHistograms", true, "Enables TPC TOF histograms"}; Configurable enableDCAxyzHistograms{"enableDCAxyzHistograms", false, "Enables DCAxyz correlation histograms"}; Configurable enableDCAxyphiHistograms{"enableDCAxyphiHistograms", false, "Enables DCAxyphi correlation histograms"}; + Configurable enableDCAvsmotherHistograms{"enableDCAvsmotherHistograms", false, "Enables DCA vs mother histograms"}; struct : ConfigurableGroup { ConfigurableAxis binsPt{"binsPt", {VARIABLE_WIDTH, 0.0, 0.1, 0.12, 0.14, 0.16, 0.18, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8, 5.0}, "Binning of the pT axis"}; @@ -649,9 +576,7 @@ struct tofSpectra { hDcaXYZ[i] = histos.add(Form("dca/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); } else { histos.add(hdcaxy[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - histos.add(hdcaxytot[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); histos.add(hdcaz[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); - histos.add(hdcaztot[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); } if (doprocessMC) { @@ -698,35 +623,41 @@ struct tofSpectra { histos.add(hpt_den_prm_mcgoodev[i].data(), pTCharge[i], kTH2D, {ptAxis, multAxis}); histos.add(hpt_den_prm_mcbadev[i].data(), pTCharge[i], kTH2D, {ptAxis, multAxis}); + const std::string cpName = Form("/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]); if (enableDCAxyzHistograms) { - hDcaXYZPrm[i] = histos.add(Form("dcaprm/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); - hDcaXYZStr[i] = histos.add(Form("dcastr/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); - hDcaXYZMat[i] = histos.add(Form("dcamat/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); + hDcaXYZPrm[i] = histos.add("dcaprm" + cpName, pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); + hDcaXYZStr[i] = histos.add("dcastr" + cpName, pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); + hDcaXYZMat[i] = histos.add("dcamat" + cpName, pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); if (enableDcaGoodEvents) { histos.add(hdcaxyprmgoodevs[i].data(), pTCharge[i], kTH3D, {ptAxis, dcaXyAxis, dcaZAxis}); } } else { - hDcaXYWrongCollisionPrm[i] = histos.add(Form("dcaxywrongcollprm/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - hDcaXYWrongCollisionStr[i] = histos.add(Form("dcaxywrongcollstr/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - hDcaXYWrongCollisionMat[i] = histos.add(Form("dcaxywrongcollmat/%s/%s", (i < Np) ? "pos" : "neg", pN[i % Np]), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaXYWrongCollisionPrm[i] = histos.add("dcaxywrongcollprm" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaXYWrongCollisionStr[i] = histos.add("dcaxywrongcollstr" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaXYWrongCollisionMat[i] = histos.add("dcaxywrongcollmat" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + histos.add(hdcaxyprm[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); histos.add(hdcazprm[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); - histos.add(hdcaxyprm2[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - histos.add(hdcazprm2[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); - histos.add(hdcaxyD0[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - histos.add(hdcazD0[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); histos.add(hdcaxystr[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); histos.add(hdcazstr[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); - histos.add(hdcaxycharm[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - histos.add(hdcazcharm[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); - histos.add(hdcaxybeauty[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); - histos.add(hdcazbeauty[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); histos.add(hdcaxymat[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); histos.add(hdcazmat[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); if (enableDcaGoodEvents) { histos.add(hdcaxyprmgoodevs[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); histos.add(hdcazprmgoodevs[i].data(), pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); } + if (enableDCAvsmotherHistograms) { + hDcaXYMC[i] = histos.add("dcaxymc" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaZMC[i] = histos.add("dcazmc" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); + hDcaXYMCNotHF[i] = histos.add("dcaxynothf" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaZMCNotHF[i] = histos.add("dcaznothf" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); + hDcaXYMCD0[i] = histos.add("dcaxyD0" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaZMCD0[i] = histos.add("dcazD0" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); + hDcaXYMCCharm[i] = histos.add("dcaxycharm" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hdcaZMCCharm[i] = histos.add("dcazcharm" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); + hDcaXYMCBeauty[i] = histos.add("dcaxybeauty" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaXyAxis}); + hDcaZMCBeauty[i] = histos.add("dcazbeauty" + cpName, pTCharge[i], kTH2D, {ptAxis, dcaZAxis}); + } } // Mismatched info @@ -1620,7 +1551,10 @@ struct tofSpectra { using RecoMCCollisions = soa::Join; // RD template - void fillTrackHistograms_MC(TrackType const& track, ParticleType const& mcParticle, RecoMCCollisions::iterator const& collision) + void fillTrackHistograms_MC(TrackType const& track, + ParticleType::iterator const& mcParticle, + RecoMCCollisions::iterator const& collision, + ParticleType const& mcParticles) { if (!isParticleEnabled()) { // Check if the particle is enabled return; @@ -1645,10 +1579,11 @@ struct tofSpectra { return; } - histos.fill(HIST(hdcaxytot[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcaztot[i]), track.pt(), track.dcaZ()); - - if (!mcParticle.isPhysicalPrimary()) { + if (enableDCAvsmotherHistograms) { + hDcaXYMC[i]->Fill(track.pt(), track.dcaXY()); + hDcaZMC[i]->Fill(track.pt(), track.dcaZ()); + } + if (!mcParticle.isPhysicalPrimary()) { // Secondaries (weak decays and material) if (mcParticle.getProcess() == 4) { if (enableDCAxyzHistograms) { hDcaXYZStr[i]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); @@ -1664,55 +1599,61 @@ struct tofSpectra { histos.fill(HIST(hdcazmat[i]), track.pt(), track.dcaZ()); } } - } else { + } else { // Primaries if (enableDCAxyzHistograms) { hDcaXYZPrm[i]->Fill(track.pt(), track.dcaXY(), track.dcaZ()); if (enableDcaGoodEvents.value && collision.has_mcCollision()) { histos.fill(HIST(hdcaxyprmgoodevs[i]), track.pt(), track.dcaXY(), track.dcaZ()); } } else { - // DCAxy for all non primaries - histos.fill(HIST(hdcaxyprm2[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazprm2[i]), track.pt(), track.dcaZ()); + // DCAxy for all primaries + histos.fill(HIST(hdcaxyprm[i]), track.pt(), track.dcaXY()); + histos.fill(HIST(hdcazprm[i]), track.pt(), track.dcaZ()); + } + if (enableDcaGoodEvents.value && collision.has_mcCollision()) { + histos.fill(HIST(hdcaxyprmgoodevs[i]), track.pt(), track.dcaXY()); + histos.fill(HIST(hdcazprmgoodevs[i]), track.pt(), track.dcaZ()); + } + + if (enableDCAvsmotherHistograms) { + bool IsD0Mother = false; + bool IsCharmMother = false; + bool IsBeautyMother = false; - // bool IsD0Mother{false}; - bool IsCharmMother{false}; - bool IsBeautyMother{false}; if (mcParticle.has_mothers()) { + const int charmOrigin = RecoDecay::getCharmHadronOrigin(mcParticles, mcParticle, false); for (const auto& mother : mcParticle.template mothers_as()) { - int motherPdgCode = mother.pdgCode(); + const int motherPdgCode = mother.pdgCode(); if (motherPdgCode == 421) { - std::cout << "Mother PDG Code: " << motherPdgCode << std::endl; - histos.fill(HIST(hdcaxyD0[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazD0[i]), track.pt(), track.dcaZ()); - // IsD0Mother = true; + IsD0Mother = true; } - if (RecoDecay::OriginType::NonPrompt) { + if (charmOrigin == RecoDecay::OriginType::NonPrompt) { if ((motherPdgCode) / 1000 == 5 || (motherPdgCode) / 100 == 5) { - std::cout << "Mother PDG Code starts with 5: " << motherPdgCode << std::endl; - histos.fill(HIST(hdcaxybeauty[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazbeauty[i]), track.pt(), track.dcaZ()); IsBeautyMother = true; } } - if (RecoDecay::OriginType::Prompt) { + if (charmOrigin == RecoDecay::OriginType::Prompt) { if ((motherPdgCode) / 1000 == 4 || (motherPdgCode) / 100 == 4) { - std::cout << "Mother PDG Code starts with 4: " << motherPdgCode << std::endl; - histos.fill(HIST(hdcaxycharm[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazcharm[i]), track.pt(), track.dcaZ()); IsCharmMother = true; } } } } - if (!IsCharmMother && !IsBeautyMother) { - histos.fill(HIST(hdcaxyprm[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazprm[i]), track.pt(), track.dcaZ()); + if (IsD0Mother) { + hDcaXYMCD0[i]->Fill(track.pt(), track.dcaXY()); + hDcaZMCD0[i]->Fill(track.pt(), track.dcaZ()); } - - if (enableDcaGoodEvents.value && collision.has_mcCollision()) { - histos.fill(HIST(hdcaxyprmgoodevs[i]), track.pt(), track.dcaXY()); - histos.fill(HIST(hdcazprmgoodevs[i]), track.pt(), track.dcaZ()); + if (IsCharmMother) { + hDcaXYMCBeauty[i]->Fill(track.pt(), track.dcaXY()); + hDcaZMCBeauty[i]->Fill(track.pt(), track.dcaZ()); + } + if (IsBeautyMother) { + hDcaXYMCCharm[i]->Fill(track.pt(), track.dcaXY()); + hdcaZMCCharm[i]->Fill(track.pt(), track.dcaZ()); + } + if (!IsCharmMother && !IsBeautyMother) { + hDcaXYMCNotHF[i]->Fill(track.pt(), track.dcaXY()); + hDcaZMCNotHF[i]->Fill(track.pt(), track.dcaZ()); } } } From ab3a615f6dbba0988fa62f0ac95f6965090d21fd Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Nicol=C3=B2=20Jacazio?= Date: Sun, 22 Sep 2024 07:56:08 +0200 Subject: [PATCH 2/2] FIx --- PWGLF/Tasks/Nuspex/spectraTOF.cxx | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/PWGLF/Tasks/Nuspex/spectraTOF.cxx b/PWGLF/Tasks/Nuspex/spectraTOF.cxx index 7efc6821ed2..c6449669458 100644 --- a/PWGLF/Tasks/Nuspex/spectraTOF.cxx +++ b/PWGLF/Tasks/Nuspex/spectraTOF.cxx @@ -54,7 +54,6 @@ std::array, NpCharge> hDcaXYWrongCollisionStr; std::array, NpCharge> hDcaXYWrongCollisionMat; std::array, NpCharge> hDcaXYMC; // DCA xy in the MC std::array, NpCharge> hDcaZMC; // DCA z in the MC -std::array, NpCharge> hDcaXYMCPrm; // DCA xy in the MC std::array, NpCharge> hDcaXYMCD0; // DCA xy in the MC for particles from D0 std::array, NpCharge> hDcaZMCD0; // DCA z in the MC for particles from D0 std::array, NpCharge> hDcaXYMCCharm; // DCA xy in the MC for particles from charm @@ -1644,13 +1643,13 @@ struct tofSpectra { hDcaZMCD0[i]->Fill(track.pt(), track.dcaZ()); } if (IsCharmMother) { - hDcaXYMCBeauty[i]->Fill(track.pt(), track.dcaXY()); - hDcaZMCBeauty[i]->Fill(track.pt(), track.dcaZ()); - } - if (IsBeautyMother) { hDcaXYMCCharm[i]->Fill(track.pt(), track.dcaXY()); hdcaZMCCharm[i]->Fill(track.pt(), track.dcaZ()); } + if (IsBeautyMother) { + hDcaXYMCBeauty[i]->Fill(track.pt(), track.dcaXY()); + hDcaZMCBeauty[i]->Fill(track.pt(), track.dcaZ()); + } if (!IsCharmMother && !IsBeautyMother) { hDcaXYMCNotHF[i]->Fill(track.pt(), track.dcaXY()); hDcaZMCNotHF[i]->Fill(track.pt(), track.dcaZ());