Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 0 additions & 52 deletions PWGHF/Core/HfHelper.h
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,6 @@
#include <Math/GenVector/Boost.h>
#include <Math/Vector4D.h>
#include <TPDGCode.h>
#include <vector>

#include "CommonConstants/PhysicsConstants.h"

Expand Down Expand Up @@ -932,57 +931,6 @@ class HfHelper
return true;
}

/// Get the pt of the beauty hadron(non-prompt case), return negative value for prompt case
/// \param particlesMC table with MC particles
/// \param particle MC particle
/// \return a double corresponding to beauty mother pt
template <typename T>
double getBeautyMotherPt(const T& particlesMC,
const typename T::iterator& particle)
{
int stage = 0; // mother tree level

// vector of vectors with mother indices; each line corresponds to a "stage"
std::vector<std::vector<int64_t>> arrayIds{};
std::vector<int64_t> initVec{particle.globalIndex()};
arrayIds.push_back(initVec); // the first vector contains the index of the original particle

while (arrayIds[-stage].size() > 0) {
// vector of mother indices for the current stage
std::vector<int64_t> 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

if ((PDGParticleIMother > 500 && PDGParticleIMother < 600) || (PDGParticleIMother > 5000 && PDGParticleIMother < 6000)) {
return mother.pt();
}
if (PDGParticleIMother == 4) {
return -999.;
}
if (PDGParticleIMother == 5) {
return -1.;
}

// 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--;
}
return -999.;
}

private:
};

Expand Down
25 changes: 5 additions & 20 deletions PWGHF/D2H/Tasks/taskD0.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -436,9 +436,6 @@ struct HfTaskD0 {
auto ctCandidate = hfHelper.ctD0(candidate);
auto cpaCandidate = candidate.cpa();
auto cpaxyCandidate = candidate.cpaXY();
float ptB = -1;
auto indexMother = RecoDecay::getMother(mcParticles, candidate.template prong0_as<aod::TracksWMc>().template mcParticle_as<soa::Join<aod::McParticles, aod::HfCand2ProngMcGen>>(), o2::constants::physics::Pdg::kD0, true);
auto particleMother = mcParticles.rawIteratorAt(indexMother);
if (candidate.isSelD0() >= selectionFlagD0) {
registry.fill(HIST("hMassSigBkgD0"), massD0, ptCandidate, rapidityCandidate);
if (candidate.flagMcMatchRec() == (1 << aod::hf_cand_2prong::DecayType::D0ToPiK)) {
Expand Down Expand Up @@ -466,10 +463,7 @@ struct HfTaskD0 {
registry.fill(HIST("hDecLengthxyVsPtSig"), declengthxyCandidate, ptCandidate);
registry.fill(HIST("hMassSigD0"), massD0, ptCandidate, rapidityCandidate);
if constexpr (applyMl) {
if (candidate.originMcRec() == RecoDecay::OriginType::NonPrompt) {
ptB = hfHelper.getBeautyMotherPt(mcParticles, particleMother);
}
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0()[0], candidate.mlProbD0()[1], massD0, ptCandidate, ptB, rapidityCandidate, candidate.originMcRec(), SigD0);
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0()[0], candidate.mlProbD0()[1], massD0, ptCandidate, candidate.ptBhadMotherPart(), rapidityCandidate, candidate.originMcRec(), SigD0);
}
} else {
registry.fill(HIST("hPtProng0Bkg"), ptProng0, rapidityCandidate);
Expand All @@ -489,10 +483,7 @@ struct HfTaskD0 {
if (candidate.flagMcMatchRec() == -(1 << aod::hf_cand_2prong::DecayType::D0ToPiK)) {
registry.fill(HIST("hMassReflBkgD0"), massD0, ptCandidate, rapidityCandidate);
if constexpr (applyMl) {
if (candidate.originMcRec() == RecoDecay::OriginType::NonPrompt) {
ptB = hfHelper.getBeautyMotherPt(mcParticles, particleMother);
}
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0()[0], candidate.mlProbD0()[1], massD0, ptCandidate, ptB, rapidityCandidate, candidate.originMcRec(), ReflectedD0);
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0()[0], candidate.mlProbD0()[1], massD0, ptCandidate, candidate.ptBhadMotherPart(), rapidityCandidate, candidate.originMcRec(), ReflectedD0);
}
}
}
Expand All @@ -502,20 +493,14 @@ struct HfTaskD0 {
if (candidate.flagMcMatchRec() == -(1 << aod::hf_cand_2prong::DecayType::D0ToPiK)) {
registry.fill(HIST("hMassSigD0bar"), massD0bar, ptCandidate, rapidityCandidate);
if constexpr (applyMl) {
if (candidate.originMcRec() == RecoDecay::OriginType::NonPrompt) {
ptB = hfHelper.getBeautyMotherPt(mcParticles, particleMother);
}
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0bar()[0], candidate.mlProbD0bar()[1], massD0bar, ptCandidate, ptB, rapidityCandidate, candidate.originMcRec(), SigD0bar);
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0bar()[0], candidate.mlProbD0bar()[1], massD0bar, ptCandidate, candidate.ptBhadMotherPart(), rapidityCandidate, candidate.originMcRec(), SigD0bar);
}
} else {
registry.fill(HIST("hMassBkgD0bar"), massD0bar, ptCandidate, rapidityCandidate);
if (candidate.flagMcMatchRec() == (1 << aod::hf_cand_2prong::DecayType::D0ToPiK)) {
registry.fill(HIST("hMassReflBkgD0bar"), massD0bar, ptCandidate, rapidityCandidate);
if constexpr (applyMl) {
if (candidate.originMcRec() == RecoDecay::OriginType::NonPrompt) {
ptB = hfHelper.getBeautyMotherPt(mcParticles, particleMother);
}
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0bar()[0], candidate.mlProbD0bar()[1], massD0bar, ptCandidate, ptB, rapidityCandidate, candidate.originMcRec(), ReflectedD0bar);
registry.fill(HIST("hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type"), candidate.mlProbD0bar()[0], candidate.mlProbD0bar()[1], massD0bar, ptCandidate, candidate.ptBhadMotherPart(), rapidityCandidate, candidate.originMcRec(), ReflectedD0bar);
}
}
}
Expand All @@ -539,7 +524,7 @@ struct HfTaskD0 {
registry.fill(HIST("hPtVsYGenPrompt"), ptGen, yGen);
registry.fill(HIST("hSparseAcc"), ptGen, ptGenB, yGen, 1);
} else {
ptGenB = hfHelper.getBeautyMotherPt(mcParticles, particle);
ptGenB = mcParticles.rawIteratorAt(particle.idxBhadMotherPart()).pt();
registry.fill(HIST("hPtGenNonPrompt"), ptGen);
registry.fill(HIST("hYGenNonPrompt"), yGen);
registry.fill(HIST("hPtVsYGenNonPrompt"), ptGen, yGen);
Expand Down
57 changes: 39 additions & 18 deletions PWGHF/DataModel/CandidateReconstructionTables.h
Original file line number Diff line number Diff line change
Expand Up @@ -523,6 +523,10 @@ DECLARE_SOA_DYNAMIC_COLUMN(Ct, ct, //!
DECLARE_SOA_DYNAMIC_COLUMN(ImpactParameterXY, impactParameterXY, //!
[](float xVtxP, float yVtxP, float zVtxP, float xVtxS, float yVtxS, float zVtxS, float px, float py, float pz) -> float { return RecoDecay::impParXY(std::array{xVtxP, yVtxP, zVtxP}, std::array{xVtxS, yVtxS, zVtxS}, std::array{px, py, pz}); });
DECLARE_SOA_COLUMN(KfTopolChi2OverNdf, kfTopolChi2OverNdf, float); //! chi2overndf of the KFParticle topological constraint
// B-hadron mother information
DECLARE_SOA_COLUMN(PtBhadMotherPart, ptBhadMotherPart, float); //! pt of the first B-hadron mother particle (only in case of non-prompt)
DECLARE_SOA_COLUMN(PdgBhadMotherPart, pdgBhadMotherPart, int); //! pdg of the first B-hadron mother particle (only in case of non-prompt)
DECLARE_SOA_COLUMN(IdxBhadMotherPart, idxBhadMotherPart, int); //! index of the first B-hadron mother particle (only in case of non-prompt)

// method of secondary-vertex reconstruction
enum VertexerType { DCAFitter = 0,
Expand Down Expand Up @@ -638,12 +642,15 @@ DECLARE_SOA_TABLE(HfCand2ProngKF, "AOD", "HFCAND2PKF",
// table with results of reconstruction level MC matching
DECLARE_SOA_TABLE(HfCand2ProngMcRec, "AOD", "HFCAND2PMCREC", //!
hf_cand_2prong::FlagMcMatchRec,
hf_cand_2prong::OriginMcRec);
hf_cand_2prong::OriginMcRec,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart);

// table with results of generator level MC matching
DECLARE_SOA_TABLE(HfCand2ProngMcGen, "AOD", "HFCAND2PMCGEN", //!
hf_cand_2prong::FlagMcMatchGen,
hf_cand_2prong::OriginMcGen);
hf_cand_2prong::OriginMcGen,
hf_cand::IdxBhadMotherPart);

// cascade decay candidate table

Expand Down Expand Up @@ -734,12 +741,15 @@ using HfCandCascade = HfCandCascExt;
// table with results of reconstruction level MC matching for Cascade
DECLARE_SOA_TABLE(HfCandCascadeMcRec, "AOD", "HFCANDCASCMCREC", //!
hf_cand_casc::FlagMcMatchRec,
hf_cand_casc::OriginMcRec);
hf_cand_casc::OriginMcRec,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart);

// table with results of generator level MC matching
DECLARE_SOA_TABLE(HfCandCascadeMcGen, "AOD", "HFCANDCASCMCGEN", //!
hf_cand_casc::FlagMcMatchGen,
hf_cand_casc::OriginMcGen);
hf_cand_casc::OriginMcGen,
hf_cand::IdxBhadMotherPart);

// specific BPlus candidate properties
namespace hf_cand_bplus
Expand Down Expand Up @@ -908,13 +918,16 @@ DECLARE_SOA_TABLE(HfCand3ProngMcRec, "AOD", "HFCAND3PMCREC", //!
hf_cand_3prong::FlagMcMatchRec,
hf_cand_3prong::OriginMcRec,
hf_cand_3prong::IsCandidateSwapped,
hf_cand_3prong::FlagMcDecayChanRec);
hf_cand_3prong::FlagMcDecayChanRec,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart);

// table with results of generator level MC matching
DECLARE_SOA_TABLE(HfCand3ProngMcGen, "AOD", "HFCAND3PMCGEN", //!
hf_cand_3prong::FlagMcMatchGen,
hf_cand_3prong::OriginMcGen,
hf_cand_3prong::FlagMcDecayChanGen);
hf_cand_3prong::FlagMcDecayChanGen,
hf_cand::IdxBhadMotherPart);

namespace hf_cand_casc_lf
{
Expand Down Expand Up @@ -1265,39 +1278,47 @@ DECLARE_SOA_TABLE(HfXicToXiPiMCRec, "AOD", "HFXICXIPIMCREC", //!
hf_cand_xic0_omegac0::DebugMcRec,
hf_cand_xic0_omegac0::OriginRec,
hf_cand_xic0_omegac0::CollisionMatched,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart,
o2::soa::Marker<1>);
DECLARE_SOA_TABLE(HfOmegacToXiPiMCRec, "AOD", "HFOMCXIPIMCREC", //!
hf_cand_xic0_omegac0::FlagMcMatchRec,
hf_cand_xic0_omegac0::DebugMcRec,
hf_cand_xic0_omegac0::OriginRec,
hf_cand_xic0_omegac0::CollisionMatched,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart,
o2::soa::Marker<2>);
DECLARE_SOA_TABLE(HfToOmegaPiMCRec, "AOD", "HFTOOMEPIMCREC", //!
hf_cand_xic0_omegac0::FlagMcMatchRec,
hf_cand_xic0_omegac0::DebugMcRec,
hf_cand_xic0_omegac0::OriginRec,
hf_cand_xic0_omegac0::CollisionMatched,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart,
o2::soa::Marker<3>);
DECLARE_SOA_TABLE(HfToOmegaKMCRec, "AOD", "HFTOOMEKMCREC", //!
hf_cand_xic0_omegac0::FlagMcMatchRec,
hf_cand_xic0_omegac0::DebugMcRec,
hf_cand_xic0_omegac0::OriginRec,
hf_cand_xic0_omegac0::CollisionMatched,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart,
o2::soa::Marker<4>);

// table with results of generator level MC matching
DECLARE_SOA_TABLE(HfXicToXiPiMCGen, "AOD", "HFXICXIPIMCGEN", //!
hf_cand_xic0_omegac0::FlagMcMatchGen, hf_cand_xic0_omegac0::DebugGenCharmBar, hf_cand_xic0_omegac0::DebugGenCasc, hf_cand_xic0_omegac0::DebugGenLambda,
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, o2::soa::Marker<1>);
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, hf_cand::IdxBhadMotherPart, o2::soa::Marker<1>);
DECLARE_SOA_TABLE(HfOmegacToXiPiMCGen, "AOD", "HFOMECXIPIMCGEN", //!
hf_cand_xic0_omegac0::FlagMcMatchGen, hf_cand_xic0_omegac0::DebugGenCharmBar, hf_cand_xic0_omegac0::DebugGenCasc, hf_cand_xic0_omegac0::DebugGenLambda,
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, o2::soa::Marker<2>);
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, hf_cand::IdxBhadMotherPart, o2::soa::Marker<2>);
DECLARE_SOA_TABLE(HfToOmegaPiMCGen, "AOD", "HFTOOMEPIMCGEN", //!
hf_cand_xic0_omegac0::FlagMcMatchGen, hf_cand_xic0_omegac0::DebugGenCharmBar, hf_cand_xic0_omegac0::DebugGenCasc, hf_cand_xic0_omegac0::DebugGenLambda,
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, o2::soa::Marker<3>);
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, hf_cand::IdxBhadMotherPart, o2::soa::Marker<3>);
DECLARE_SOA_TABLE(HfToOmegaKMCGen, "AOD", "HFTOOMEKMCGEN", //!
hf_cand_xic0_omegac0::FlagMcMatchGen, hf_cand_xic0_omegac0::DebugGenCharmBar, hf_cand_xic0_omegac0::DebugGenCasc, hf_cand_xic0_omegac0::DebugGenLambda,
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, o2::soa::Marker<4>);
hf_cand_xic0_omegac0::PtCharmBaryonGen, hf_cand_xic0_omegac0::EtaCharmBaryonGen, hf_cand_xic0_omegac0::OriginGen, hf_cand::IdxBhadMotherPart, o2::soa::Marker<4>);

// specific chic candidate properties
namespace hf_cand_chic
Expand Down Expand Up @@ -1641,12 +1662,15 @@ using HfCandSc = HfCandScExt;
// table with results of reconstruction level MC matching
DECLARE_SOA_TABLE(HfCandScMcRec, "AOD", "HFCANDSCMCREC", //!
hf_cand_sigmac::FlagMcMatchRec,
hf_cand_sigmac::OriginMcRec);
hf_cand_sigmac::OriginMcRec,
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart);

// table with results of generation level MC matching
DECLARE_SOA_TABLE(HfCandScMcGen, "AOD", "HFCANDSCMCGEN", //!
hf_cand_sigmac::FlagMcMatchGen,
hf_cand_sigmac::OriginMcGen);
hf_cand_sigmac::OriginMcGen,
hf_cand::IdxBhadMotherPart);

/// D*± → D0(bar) π±
namespace hf_cand_dstar
Expand Down Expand Up @@ -1755,12 +1779,9 @@ DECLARE_SOA_DYNAMIC_COLUMN(PVecSoftPi, pVecSoftPi, [](float px, float py, float
// MC matching result:
DECLARE_SOA_COLUMN(FlagMcMatchRec, flagMcMatchRec, int8_t); //! reconstruction level
DECLARE_SOA_COLUMN(FlagMcMatchGen, flagMcMatchGen, int8_t); //! generator level
DECLARE_SOA_COLUMN(PtBhadMotherPart, ptBhadMotherPart, float); //! pt of the first B-hadron mother particle (only in case of non-prompt)
DECLARE_SOA_COLUMN(PdgBhadMotherPart, pdgBhadMotherPart, int); //! pdg of the first B-hadron mother particle (only in case of non-prompt)

DECLARE_SOA_COLUMN(OriginMcRec, originMcRec, int8_t); //! particle origin, reconstruction level
DECLARE_SOA_COLUMN(OriginMcGen, originMcGen, int8_t); //! particle origin, generator level
DECLARE_SOA_COLUMN(IdxBhadMotherPart, idxBhadMotherPart, int); //! index of the first B-hadron mother particle (only in case of non-prompt)

enum DecayType {
DstarToD0Pi = 0,
Expand Down Expand Up @@ -1879,14 +1900,14 @@ using HfCandDstar = HfCandDstars::iterator;
DECLARE_SOA_TABLE(HfCandDstarMcRec, "AOD", "HFCANDDSTRMCREC",
hf_cand_dstar::FlagMcMatchRec,
hf_cand_dstar::OriginMcRec,
hf_cand_dstar::PtBhadMotherPart,
hf_cand_dstar::PdgBhadMotherPart);
hf_cand::PtBhadMotherPart,
hf_cand::PdgBhadMotherPart);

// table with results of generator level MC matching
DECLARE_SOA_TABLE(HfCandDstarMcGen, "AOD", "HFCANDDSTRMCGEN",
hf_cand_dstar::FlagMcMatchGen,
hf_cand_dstar::OriginMcGen,
hf_cand_dstar::IdxBhadMotherPart);
hf_cand::IdxBhadMotherPart);

#undef HFCAND_COLUMNS

Expand Down
24 changes: 16 additions & 8 deletions PWGHF/TableProducer/candidateCreator2Prong.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -687,10 +687,10 @@ struct HfCandidateCreator2ProngExpressions {
// Match reconstructed candidates.
// Spawned table can be used directly
for (const auto& candidate : *rowCandidateProng2) {

flag = 0;
origin = 0;
auto arrayDaughters = std::array{candidate.prong0_as<aod::TracksWMc>(), candidate.prong1_as<aod::TracksWMc>()};
std::vector<int> idxBhadMothers{};

// D0(bar) → π± K∓
indexRec = RecoDecay::getMatchedMCRec(mcParticles, arrayDaughters, Pdg::kD0, std::array{+kPiPlus, -kKPlus}, true, &sign);
Expand All @@ -717,21 +717,26 @@ struct HfCandidateCreator2ProngExpressions {
// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
auto particle = mcParticles.rawIteratorAt(indexRec);
origin = RecoDecay::getCharmHadronOrigin(mcParticles, particle);
origin = RecoDecay::getCharmHadronOrigin(mcParticles, particle, false, &idxBhadMothers);
}
if (origin == RecoDecay::OriginType::NonPrompt) {
auto bHadMother = mcParticles.rawIteratorAt(idxBhadMothers[0]);
rowMcMatchRec(flag, origin, bHadMother.pt(), bHadMother.pdgCode());
} else {
rowMcMatchRec(flag, origin, -1.f, 0);
}

rowMcMatchRec(flag, origin);
}

// Match generated particles.
for (const auto& particle : mcParticles) {
flag = 0;
origin = 0;
std::vector<int> idxBhadMothers{};

auto mcCollision = particle.mcCollision();
float zPv = mcCollision.posZ();
if (zPv < -zPvPosMax || zPv > zPvPosMax) { // to avoid counting particles in collisions with Zvtx larger than the maximum, we do not match them
rowMcMatchGen(flag, origin);
rowMcMatchGen(flag, origin, -1);
continue;
}

Expand All @@ -756,10 +761,13 @@ struct HfCandidateCreator2ProngExpressions {

// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
origin = RecoDecay::getCharmHadronOrigin(mcParticles, particle);
origin = RecoDecay::getCharmHadronOrigin(mcParticles, particle, false, &idxBhadMothers);
}
if (origin == RecoDecay::OriginType::NonPrompt) {
rowMcMatchGen(flag, origin, idxBhadMothers[0]);
} else {
rowMcMatchGen(flag, origin, -1);
}

rowMcMatchGen(flag, origin);
}
}

Expand Down
Loading