Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[RecoDecay, Event Filtering, PWGHF] Revert to MC particle slices #695

Merged
merged 1 commit into from
Apr 28, 2022
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
26 changes: 9 additions & 17 deletions Common/Core/RecoDecay.h
Original file line number Diff line number Diff line change
Expand Up @@ -562,16 +562,14 @@ class RecoDecay
}

/// Finds the mother of an MC particle by looking for the expected PDG code in the mother chain.
/// \param particlesMC table with MC particles
/// \param particle MC particle
/// \param PDGMother expected mother PDG code
/// \param acceptAntiParticles switch to accept the antiparticle of the expected mother
/// \param sign antiparticle indicator of the found mother w.r.t. PDGMother; 1 if particle, -1 if antiparticle, 0 if mother not found
/// \param depthMax maximum decay tree level to check; Mothers up to this level will be considered. If -1, all levels are considered.
/// \return index of the mother particle if found, -1 otherwise
template <typename T>
static int getMother(const T& particlesMC,
const typename T::iterator& particle,
static int getMother(const T& particle,
int PDGMother,
bool acceptAntiParticles = false,
int8_t* sign = nullptr,
Expand All @@ -588,8 +586,8 @@ class RecoDecay
if (depthMax > -1 && -stage >= depthMax) { // Maximum depth has been reached.
return -1;
}
auto indexMotherTmp = particleMother.mothersIds().front();
particleMother = particlesMC.rawIteratorAt(indexMotherTmp - particlesMC.offset());
particleMother = particleMother.template mothers_first_as<typename std::decay_t<T>::parent_t>(); // get mother 0
auto indexMotherTmp = particleMother.globalIndex();
// Check mother's PDG code.
auto PDGParticleIMother = particleMother.pdgCode(); // PDG code of the mother
//printf("getMother: ");
Expand All @@ -614,7 +612,6 @@ class RecoDecay
}

/// Gets the complete list of indices of final-state daughters of an MC particle.
/// \param particlesMC table with MC particles
/// \param particle MC particle
/// \param list vector where the indices of final-state daughters will be added
/// \param arrPDGFinal array of PDG codes of particles to be considered final if found
Expand All @@ -623,8 +620,7 @@ class RecoDecay
/// \note Final state is defined as particles from arrPDGFinal plus final daughters of any other decay branch.
/// \note Antiparticles of particles in arrPDGFinal are accepted as well.
template <std::size_t N, typename T>
static void getDaughters(const T& particlesMC,
const typename T::iterator& particle,
static void getDaughters(const T& particle,
std::vector<int>* list,
const array<int, N>& arrPDGFinal,
int8_t depthMax = -1,
Expand Down Expand Up @@ -675,12 +671,8 @@ class RecoDecay
//printf("Stage %d: %d (PDG %d) -> %d-%d\n", stage, index, PDGParticle, indexDaughterFirst, indexDaughterLast);
// Call itself to get daughters of daughters recursively.
stage++;
if (particle.daughtersIds().front() != particle.daughtersIds().back()) {
for (auto& idxDau : particle.daughtersIds()) {
getDaughters(particlesMC, particlesMC.rawIteratorAt(idxDau - particlesMC.offset()), list, arrPDGFinal, depthMax, stage);
}
} else {
getDaughters(particlesMC, particlesMC.rawIteratorAt(particle.daughtersIds().front() - particlesMC.offset()), list, arrPDGFinal, depthMax, stage);
for (auto& dau : particle.template daughters_as<typename std::decay_t<T>::parent_t>()) {
getDaughters(dau, list, arrPDGFinal, depthMax, stage);
}
}

Expand Down Expand Up @@ -721,7 +713,7 @@ class RecoDecay
if (iProng == 0) {
// Get the mother index and its sign.
// PDG code of the first daughter's mother determines whether the expected mother is a particle or antiparticle.
indexMother = getMother(particlesMC, particleI, PDGMother, acceptAntiParticles, &sgn, depthMax);
indexMother = getMother(particleI, PDGMother, acceptAntiParticles, &sgn, depthMax);
// Check whether mother was found.
if (indexMother <= -1) {
//Printf("MC Rec: Rejected: bad mother index or PDG");
Expand All @@ -740,7 +732,7 @@ class RecoDecay
return -1;
}
// Get the list of actual final daughters.
getDaughters(particlesMC, particleMother, &arrAllDaughtersIndex, arrPDGDaughters, depthMax);
getDaughters(particleMother, &arrAllDaughtersIndex, arrPDGDaughters, depthMax);
//printf("MC Rec: Mother %d has %d final daughters:", indexMother, arrAllDaughtersIndex.size());
//for (auto i : arrAllDaughtersIndex) {
// printf(" %d", i);
Expand Down Expand Up @@ -858,7 +850,7 @@ class RecoDecay
return false;
}
// Get the list of actual final daughters.
getDaughters(particlesMC, candidate, &arrAllDaughtersIndex, arrPDGDaughters, depthMax);
getDaughters(candidate, &arrAllDaughtersIndex, arrPDGDaughters, depthMax);
//printf("MC Gen: Mother %ld has %ld final daughters:", candidate.globalIndex(), arrAllDaughtersIndex.size());
//for (auto i : arrAllDaughtersIndex) {
// printf(" %d", i);
Expand Down
4 changes: 2 additions & 2 deletions EventFiltering/PWGHF/HFFilter.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -978,7 +978,7 @@ struct HfFilter { // Main struct for HF triggers
auto indexRec = RecoDecay::getMatchedMCRec(particlesMC, std::array{trackPos, trackNeg}, pdg::Code::kD0, array{+kPiPlus, -kKPlus}, true, &sign);
if (indexRec > -1) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
if (origin == OriginType::NonPrompt) {
flag = kNonPrompt;
} else {
Expand Down Expand Up @@ -1037,7 +1037,7 @@ struct HfFilter { // Main struct for HF triggers

if (indexRec > -1) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
if (origin == OriginType::NonPrompt) {
flag = kNonPrompt;
} else {
Expand Down
4 changes: 2 additions & 2 deletions PWGHF/TableProducer/HFCandidateCreator2Prong.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -185,7 +185,7 @@ struct HFCandidateCreator2ProngExpressions {
// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
}

rowMCMatchRec(flag, origin);
Expand Down Expand Up @@ -221,7 +221,7 @@ struct HFCandidateCreator2ProngExpressions {

// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
}

rowMCMatchGen(flag, origin);
Expand Down
8 changes: 4 additions & 4 deletions PWGHF/TableProducer/HFCandidateCreator3Prong.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -195,7 +195,7 @@ struct HFCandidateCreator3ProngExpressions {
if (arrayDaughters[0].has_mcParticle()) {
swapping = int8_t(std::abs(arrayDaughters[0].mcParticle().pdgCode()) == kPiPlus);
}
RecoDecay::getDaughters(particlesMC, particlesMC.rawIteratorAt(indexRec), &arrDaughIndex, array{0}, 1);
RecoDecay::getDaughters(particlesMC.rawIteratorAt(indexRec), &arrDaughIndex, array{0}, 1);
if (arrDaughIndex.size() == 2) {
for (auto iProng = 0u; iProng < arrDaughIndex.size(); ++iProng) {
auto daughI = particlesMC.rawIteratorAt(arrDaughIndex[iProng]);
Expand Down Expand Up @@ -224,7 +224,7 @@ struct HFCandidateCreator3ProngExpressions {
// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
}

rowMCMatchRec(flag, origin, swapping, channel);
Expand All @@ -251,7 +251,7 @@ struct HFCandidateCreator3ProngExpressions {
flag = sign * (1 << DecayType::LcToPKPi);

//Printf("Flagging the different Λc± → p± K∓ π± decay channels");
RecoDecay::getDaughters(particlesMC, particle, &arrDaughIndex, array{0}, 1);
RecoDecay::getDaughters(particle, &arrDaughIndex, array{0}, 1);
if (arrDaughIndex.size() == 2) {
for (auto jProng = 0u; jProng < arrDaughIndex.size(); ++jProng) {
auto daughJ = particlesMC.rawIteratorAt(arrDaughIndex[jProng]);
Expand All @@ -278,7 +278,7 @@ struct HFCandidateCreator3ProngExpressions {

// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
}

rowMCMatchGen(flag, origin, channel);
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/TableProducer/HFCandidateCreatorCascade.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -275,7 +275,7 @@ struct HFCandidateCreatorCascadeMC {
MY_DEBUG_MSG(sign, LOG(info) << "Lc in K0S p");
arrDaughLcIndex.clear();
// checking that the final daughters (decay depth = 3) are p, pi+, pi-
RecoDecay::getDaughters(particlesMC, particle, &arrDaughLcIndex, arrDaughLcPDGRef, 3); // best would be to check the K0S daughters
RecoDecay::getDaughters(particle, &arrDaughLcIndex, arrDaughLcPDGRef, 3); // best would be to check the K0S daughters
if (arrDaughLcIndex.size() == 3) {
for (std::size_t iProng = 0; iProng < arrDaughLcIndex.size(); ++iProng) {
auto daughI = particlesMC.rawIteratorAt(arrDaughLcIndex[iProng]);
Expand Down
10 changes: 5 additions & 5 deletions PWGHF/TableProducer/HFCandidateCreatorChic.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -231,15 +231,15 @@ struct HFCandidateCreatorChicMC {
if (indexRec > -1) {
hMassJpsiToMuMuMatched->Fill(InvMassJpsiToMuMu(candidate.index0()));

int indexMother = RecoDecay::getMother(particlesMC, particlesMC.rawIteratorAt(indexRec), pdg::Code::kChic1);
int indexMotherGamma = RecoDecay::getMother(particlesMC, particlesMC.rawIteratorAt(candidate.index1().mcparticleId()), pdg::Code::kChic1);
int indexMother = RecoDecay::getMother(particlesMC.rawIteratorAt(indexRec), pdg::Code::kChic1);
int indexMotherGamma = RecoDecay::getMother(particlesMC.rawIteratorAt(candidate.index1().mcparticleId()), pdg::Code::kChic1);
if (indexMother > -1 && indexMotherGamma == indexMother && candidate.index1().mcparticle().pdgCode() == kGamma) {
auto particleMother = particlesMC.rawIteratorAt(indexMother);
hEphotonMatched->Fill(candidate.index1().e());
hMassEMatched->Fill(sqrt(candidate.index1().px() * candidate.index1().px() + candidate.index1().py() * candidate.index1().py() + candidate.index1().pz() * candidate.index1().pz()));
if (particleMother.has_daughters()) {
std::vector<int> arrAllDaughtersIndex;
RecoDecay::getDaughters(particlesMC, particleMother, &arrAllDaughtersIndex, array{(int)(kGamma), (int)(pdg::Code::kJpsi)}, 1);
RecoDecay::getDaughters(particleMother, &arrAllDaughtersIndex, array{(int)(kGamma), (int)(pdg::Code::kJpsi)}, 1);
if (arrAllDaughtersIndex.size() == 2) {
flag = 1 << hf_cand_chic::DecayType::ChicToJpsiToMuMuGamma;
hMassChicToJpsiToMuMuGammaMatched->Fill(InvMassChicToJpsiGamma(candidate));
Expand All @@ -249,7 +249,7 @@ struct HFCandidateCreatorChicMC {
}
if (flag != 0) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? NonPrompt : Prompt);
origin = (RecoDecay::getMother(particle, kBottom, true) > -1 ? NonPrompt : Prompt);
}
rowMCMatchRec(flag, origin, channel);
}
Expand All @@ -265,7 +265,7 @@ struct HFCandidateCreatorChicMC {
if (RecoDecay::isMatchedMCGen(particlesMC, particle, pdg::Code::kChic1, array{(int)(pdg::Code::kJpsi), (int)(kGamma)}, true)) {
// Match J/psi --> e+e-
std::vector<int> arrDaughter;
RecoDecay::getDaughters(particlesMC, particle, &arrDaughter, array{(int)(pdg::Code::kJpsi)}, 1);
RecoDecay::getDaughters(particle, &arrDaughter, array{(int)(pdg::Code::kJpsi)}, 1);
auto jpsiCandMC = particlesMC.rawIteratorAt(arrDaughter[0]);
if (RecoDecay::isMatchedMCGen(particlesMC, jpsiCandMC, pdg::Code::kJpsi, array{+kElectron, -kElectron}, true)) {
flag = 1 << hf_cand_chic::DecayType::ChicToJpsiToEEGamma;
Expand Down
6 changes: 3 additions & 3 deletions PWGHF/TableProducer/HFCandidateCreatorX.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -305,7 +305,7 @@ struct HFCandidateCreatorXMC {
// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
auto particle = particlesMC.rawIteratorAt(indexRec);
origin = (RecoDecay::getMother(particlesMC, particle, 5, true) > -1 ? NonPrompt : Prompt);
origin = (RecoDecay::getMother(particle, 5, true) > -1 ? NonPrompt : Prompt);
}

rowMCMatchRec(flag, origin, channel);
Expand All @@ -323,7 +323,7 @@ struct HFCandidateCreatorXMC {
if (RecoDecay::isMatchedMCGen(particlesMC, particle, pdgCodeX, array{pdgCodeJpsi, +kPiPlus, -kPiPlus}, true)) {
// Match J/psi --> e+e-
std::vector<int> arrDaughter;
RecoDecay::getDaughters(particlesMC, particle, &arrDaughter, array{pdgCodeJpsi}, 1);
RecoDecay::getDaughters(particle, &arrDaughter, array{pdgCodeJpsi}, 1);
auto jpsiCandMC = particlesMC.rawIteratorAt(arrDaughter[0]);
if (RecoDecay::isMatchedMCGen(particlesMC, jpsiCandMC, pdgCodeJpsi, array{+kElectron, -kElectron}, true)) {
flag = 1 << hf_cand_x::DecayType::XToJpsiToEEPiPi;
Expand All @@ -338,7 +338,7 @@ struct HFCandidateCreatorXMC {

// Check whether the particle is non-prompt (from a b quark).
if (flag != 0) {
origin = (RecoDecay::getMother(particlesMC, particle, 5, true) > -1 ? NonPrompt : Prompt);
origin = (RecoDecay::getMother(particle, 5, true) > -1 ? NonPrompt : Prompt);
}

rowMCMatchGen(flag, origin, channel);
Expand Down
6 changes: 3 additions & 3 deletions PWGHF/Tasks/HFMCValidation.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -133,7 +133,7 @@ struct ValidationGenLevel {
int whichHadron = -1;
if (std::abs(particlePdgCode) == PDGArrayParticle[iD]) {
whichHadron = iD;
RecoDecay::getDaughters(particlesMC, particle, &listDaughters, arrPDGFinal[iD], -1);
RecoDecay::getDaughters(particle, &listDaughters, arrPDGFinal[iD], -1);
std::size_t arrayPDGsize = arrPDGFinal[iD].size() - std::count(arrPDGFinal[iD].begin(), arrPDGFinal[iD].end(), 0);
if (listDaughters.size() == arrayPDGsize) {
counter[iD]++;
Expand Down Expand Up @@ -232,7 +232,7 @@ struct ValidationRecLevel {
if (whichHad >= 0) {
int indexParticle = 0;
if (cand2Prong.index0_as<aod::BigTracksMC>().has_mcParticle()) {
indexParticle = RecoDecay::getMother(particlesMC, cand2Prong.index0_as<aod::BigTracksMC>().mcParticle(), PDGArrayParticle[whichHad], true);
indexParticle = RecoDecay::getMother(cand2Prong.index0_as<aod::BigTracksMC>().mcParticle(), PDGArrayParticle[whichHad], true);
}
auto mother = particlesMC.rawIteratorAt(indexParticle);
histPt[whichHad]->Fill(cand2Prong.pt() - mother.pt());
Expand Down Expand Up @@ -278,7 +278,7 @@ struct ValidationRecLevel {
if (whichHad >= 0) {
int indexParticle = 0;
if (cand3Prong.index0_as<aod::BigTracksMC>().has_mcParticle()) {
indexParticle = RecoDecay::getMother(particlesMC, cand3Prong.index0_as<aod::BigTracksMC>().mcParticle(), PDGArrayParticle[whichHad], true);
indexParticle = RecoDecay::getMother(cand3Prong.index0_as<aod::BigTracksMC>().mcParticle(), PDGArrayParticle[whichHad], true);
}
auto mother = particlesMC.rawIteratorAt(indexParticle);
histPt[whichHad]->Fill(cand3Prong.pt() - mother.pt());
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/Tasks/taskBPlus.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -182,7 +182,7 @@ struct HfTaskBplusMc {
}
if (std::abs(candidate.flagMCMatchRec()) == 1 << hf_cand_bplus::DecayType::BPlusToD0Pi) {

auto indexMother = RecoDecay::getMother(particlesMC, candidate.index1_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandBPMCGen>>(), pdg::Code::kBPlus, true);
auto indexMother = RecoDecay::getMother(candidate.index1_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandBPMCGen>>(), pdg::Code::kBPlus, true);
auto particleMother = particlesMC.rawIteratorAt(indexMother);
registry.fill(HIST("hPtGenSig"), particleMother.pt());
registry.fill(HIST("hPtRecSig"), candidate.pt());
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/Tasks/taskChic.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -165,7 +165,7 @@ struct TaskChicMC {
}
if (candidate.flagMCMatchRec() == 1 << decayMode) {
//FIXME the access to the MC particle gen not yet functional
//int indexMother = RecoDecay::getMother(particlesMC, particlesMC.rawIteratorAt(candidate.index1().mcParticle_as<aod::McParticles_000>().globalIndex()), 20443);
//int indexMother = RecoDecay::getMother(particlesMC.rawIteratorAt(candidate.index1().mcParticle_as<aod::McParticles_000>().globalIndex()), 20443);
//auto particleMother = particlesMC.rawIteratorAt(indexMother);
//registry.fill(HIST("hPtGenSig"), particleMother.pt());
registry.fill(HIST("hPtRecSig"), candidate.pt());
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/Tasks/taskD0.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -182,7 +182,7 @@ struct TaskD0 {
}
if (std::abs(candidate.flagMCMatchRec()) == 1 << DecayType::D0ToPiK) {
// Get the corresponding MC particle.
auto indexMother = RecoDecay::getMother(particlesMC, candidate.index0_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandProng2MCGen>>(), pdg::Code::kD0, true);
auto indexMother = RecoDecay::getMother(candidate.index0_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandProng2MCGen>>(), pdg::Code::kD0, true);
auto particleMother = particlesMC.rawIteratorAt(indexMother);
registry.fill(HIST("hPtGenSig"), particleMother.pt()); // gen. level pT
auto ptRec = candidate.pt();
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/Tasks/taskDPlus.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -144,7 +144,7 @@ struct TaskDPlus {
}
if (std::abs(candidate.flagMCMatchRec()) == 1 << DecayType::DPlusToPiKPi) {
// Get the corresponding MC particle.
auto indexMother = RecoDecay::getMother(particlesMC, candidate.index0_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandProng3MCGen>>(), pdg::Code::kDPlus, true);
auto indexMother = RecoDecay::getMother(candidate.index0_as<aod::BigTracksMC>().mcParticle_as<soa::Join<aod::McParticles, aod::HfCandProng3MCGen>>(), pdg::Code::kDPlus, true);
auto particleMother = particlesMC.rawIteratorAt(indexMother);
registry.fill(HIST("hPtGenSig"), particleMother.pt()); // gen. level pT
auto ptRec = candidate.pt();
Expand Down
2 changes: 1 addition & 1 deletion PWGHF/Tasks/taskJpsi.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -231,7 +231,7 @@ struct TaskJpsiMC {
}
if (candidate.flagMCMatchRec() == 1 << decayMode) {
//Get the corresponding MC particle.
auto indexMother = RecoDecay::getMother(particlesMC, candidate.index0_as<aod::BigTracksMC>().mcParticle_as<McParticlesHf>(), pdg::Code::kJpsi, true);
auto indexMother = RecoDecay::getMother(candidate.index0_as<aod::BigTracksMC>().mcParticle_as<McParticlesHf>(), pdg::Code::kJpsi, true);
auto particleMother = particlesMC.rawIteratorAt(indexMother);
registry.fill(HIST("hPtGenSig"), particleMother.pt()); // gen. level pT
registry.fill(HIST("hPtRecSig"), candidate.pt()); // rec. level pT
Expand Down
Loading