Skip to content

Commit 24e8974

Browse files
fgrosafgrosa
authored
[Common, PWGHF, Event Filtering] Implement specific method to tag prompt/nonprompt charm hadrons (#966)
* Add function to check origin of charm hadrons * Use check origin method to tag prompt/nonprompt charm hadrons * Change name of method and move enum from HFSecondaryVertex to RecoDecay * Pass bool by value * Move enum inside RecoDecay class * Fix comment and remove unnecessary include * Resotre include in HFSecondaryVertex.h
1 parent 454d86b commit 24e8974

13 files changed

Lines changed: 101 additions & 42 deletions

‎Common/Core/RecoDecay.h‎

Lines changed: 70 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,11 @@ class RecoDecay
4949
/// Default destructor
5050
~RecoDecay() = default;
5151

52+
// mapping of charm-hadron origin type
53+
enum OriginType { None = 0,
54+
Prompt,
55+
NonPrompt };
56+
5257
// Auxiliary functions
5358

5459
/// Sums numbers.
@@ -913,6 +918,71 @@ class RecoDecay
913918
return true;
914919
}
915920

921+
/// Finds the origin (from charm hadronisation or beauty-hadron decay) of charm hadrons.
922+
/// \param particlesMC table with MC particles
923+
/// \param particle MC particle
924+
/// \param searchUpToQuark if true tag origin based on charm/beauty quark otherwise on b-hadron
925+
/// \return an integer corresponding to the origin (0: none, 1: prompt, 2: nonprompt) as in OriginType
926+
template <typename T>
927+
static int getCharmHadronOrigin(const T& particlesMC,
928+
const typename T::iterator& particle,
929+
const bool searchUpToQuark = false)
930+
{
931+
int stage = 0; // mother tree level (just for debugging)
932+
933+
// vector of vectors with mother indices; each line corresponds to a "stage"
934+
std::vector<std::vector<long int>> arrayIds{};
935+
std::vector<long int> initVec{particle.globalIndex()};
936+
arrayIds.push_back(initVec); // the first vector contains the index of the original particle
937+
938+
while (arrayIds[-stage].size() > 0) {
939+
// vector of mother indices for the current stage
940+
std::vector<long int> arrayIdsStage{};
941+
for (auto& iPart : arrayIds[-stage]) { // check all the particles that were the mothers at the previous stage
942+
auto particleMother = particlesMC.rawIteratorAt(iPart - particlesMC.offset());
943+
if (particleMother.has_mothers()) {
944+
for (auto iMother = particleMother.mothersIds().front(); iMother <= particleMother.mothersIds().back(); ++iMother) { // loop over the mother particles of the analysed particle
945+
if (std::find(arrayIdsStage.begin(), arrayIdsStage.end(), iMother) != arrayIdsStage.end()) { // if a mother is still present in the vector, do not check it again
946+
continue;
947+
}
948+
auto mother = particlesMC.rawIteratorAt(iMother - particlesMC.offset());
949+
// Check mother's PDG code.
950+
auto PDGParticleIMother = std::abs(mother.pdgCode()); // PDG code of the mother
951+
// printf("getMother: ");
952+
// for (int i = stage; i < 0; i++) // Indent to make the tree look nice.
953+
// printf(" ");
954+
// printf("Stage %d: Mother PDG: %d, Index: %d\n", stage, PDGParticleIMother, iMother);
955+
956+
if (searchUpToQuark) {
957+
if (PDGParticleIMother == 5) { // b quark
958+
return OriginType::NonPrompt;
959+
}
960+
if (PDGParticleIMother == 4) { // c quark
961+
return OriginType::Prompt;
962+
}
963+
} else {
964+
if (
965+
(PDGParticleIMother / 100 == 5 || // b mesons
966+
PDGParticleIMother / 1000 == 5) // b baryons
967+
) {
968+
return OriginType::NonPrompt;
969+
}
970+
}
971+
// add mother index in the vector for the current stage
972+
arrayIdsStage.push_back(iMother);
973+
}
974+
}
975+
}
976+
// add vector of mother indices for the current stage
977+
arrayIds.push_back(arrayIdsStage);
978+
stage--;
979+
}
980+
if (!searchUpToQuark) {
981+
return OriginType::Prompt;
982+
}
983+
return OriginType::None;
984+
}
985+
916986
private:
917987
static std::vector<std::tuple<int, double>> mListMass; ///< list of particle masses in form (PDG code, mass)
918988
};

‎EventFiltering/PWGHF/HFFilter.cxx‎

Lines changed: 5 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -39,7 +39,6 @@
3939
using namespace o2;
4040
using namespace o2::framework;
4141
using namespace o2::framework::expressions;
42-
using namespace o2::aod::hf_cand;
4342
using namespace hf_cuts_single_track;
4443
using namespace hf_cuts_bdt_multiclass;
4544

@@ -982,10 +981,10 @@ struct HfFilter { // Main struct for HF triggers
982981
auto indexRec = RecoDecay::getMatchedMCRec(particlesMC, std::array{trackPos, trackNeg}, pdg::Code::kD0, array{+kPiPlus, -kKPlus}, true, &sign);
983982
if (indexRec > -1) {
984983
auto particle = particlesMC.rawIteratorAt(indexRec);
985-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
986-
if (origin == OriginType::NonPrompt) {
984+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
985+
if (origin == RecoDecay::OriginType::NonPrompt) {
987986
flag = kNonPrompt;
988-
} else {
987+
} else if (origin == RecoDecay::OriginType::Prompt) {
989988
flag = kPrompt;
990989
}
991990
} else {
@@ -1041,8 +1040,8 @@ struct HfFilter { // Main struct for HF triggers
10411040

10421041
if (indexRec > -1) {
10431042
auto particle = particlesMC.rawIteratorAt(indexRec);
1044-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
1045-
if (origin == OriginType::NonPrompt) {
1043+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
1044+
if (origin == RecoDecay::OriginType::NonPrompt) {
10461045
flag = kNonPrompt;
10471046
} else {
10481047
flag = kPrompt;

‎PWGHF/DataModel/HFSecondaryVertex.h‎

Lines changed: 0 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -290,10 +290,6 @@ DECLARE_SOA_DYNAMIC_COLUMN(Ct, ct, //!
290290
[](float xVtxP, float yVtxP, float zVtxP, float xVtxS, float yVtxS, float zVtxS, float px, float py, float pz, double m) -> float { return RecoDecay::ct(array{px, py, pz}, RecoDecay::distance(array{xVtxP, yVtxP, zVtxP}, array{xVtxS, yVtxS, zVtxS}), m); });
291291
DECLARE_SOA_DYNAMIC_COLUMN(ImpactParameterXY, impactParameterXY, //!
292292
[](float xVtxP, float yVtxP, float zVtxP, float xVtxS, float yVtxS, float zVtxS, float px, float py, float pz) -> float { return RecoDecay::impParXY(array{xVtxP, yVtxP, zVtxP}, array{xVtxS, yVtxS, zVtxS}, array{px, py, pz}); });
293-
294-
// mapping of origin type
295-
enum OriginType { Prompt = 1,
296-
NonPrompt };
297293
} // namespace hf_cand
298294

299295
// specific 2-prong decay properties

‎PWGHF/TableProducer/HFCandidateCreator2Prong.cxx‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -210,7 +210,7 @@ struct HFCandidateCreator2ProngExpressions {
210210
// Check whether the particle is non-prompt (from a b quark).
211211
if (flag != 0) {
212212
auto particle = particlesMC.rawIteratorAt(indexRec);
213-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
213+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
214214
}
215215

216216
rowMCMatchRec(flag, origin);
@@ -246,7 +246,7 @@ struct HFCandidateCreator2ProngExpressions {
246246

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

252252
rowMCMatchGen(flag, origin);

‎PWGHF/TableProducer/HFCandidateCreator3Prong.cxx‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -251,7 +251,7 @@ struct HFCandidateCreator3ProngExpressions {
251251
// Check whether the particle is non-prompt (from a b quark).
252252
if (flag != 0) {
253253
auto particle = particlesMC.rawIteratorAt(indexRec);
254-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
254+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
255255
}
256256

257257
rowMCMatchRec(flag, origin, swapping, channel);
@@ -305,7 +305,7 @@ struct HFCandidateCreator3ProngExpressions {
305305

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

311311
rowMCMatchGen(flag, origin, channel);

‎PWGHF/TableProducer/HFCandidateCreatorChic.cxx‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -249,7 +249,7 @@ struct HFCandidateCreatorChicMC {
249249
}
250250
if (flag != 0) {
251251
auto particle = particlesMC.rawIteratorAt(indexRec);
252-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? NonPrompt : Prompt);
252+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
253253
}
254254
rowMCMatchRec(flag, origin, channel);
255255
}

‎PWGHF/TableProducer/HFCandidateCreatorX.cxx‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -305,7 +305,7 @@ struct HFCandidateCreatorXMC {
305305
// Check whether the particle is non-prompt (from a b quark).
306306
if (flag != 0) {
307307
auto particle = particlesMC.rawIteratorAt(indexRec);
308-
origin = (RecoDecay::getMother(particlesMC, particle, 5, true) > -1 ? NonPrompt : Prompt);
308+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
309309
}
310310

311311
rowMCMatchRec(flag, origin, channel);
@@ -338,7 +338,7 @@ struct HFCandidateCreatorXMC {
338338

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

344344
rowMCMatchGen(flag, origin, channel);

‎PWGHF/Tasks/HFMCValidation.cxx‎

Lines changed: 7 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -26,7 +26,6 @@ using namespace o2;
2626
using namespace o2::framework;
2727
using namespace o2::framework::expressions;
2828
using namespace o2::aod;
29-
using namespace o2::aod::hf_cand;
3029
using namespace o2::aod::hf_cand_prong2;
3130
using namespace o2::aod::hf_cand_prong3;
3231

@@ -157,10 +156,10 @@ struct ValidationGenLevel {
157156
std::size_t arrayPDGsize = arrPDGFinal[iD].size() - std::count(arrPDGFinal[iD].begin(), arrPDGFinal[iD].end(), 0);
158157
int origin = -1;
159158
if (listDaughters.size() == arrayPDGsize) {
160-
origin = (RecoDecay::getMother(particlesMC, particle, kBottom, true) > -1 ? OriginType::NonPrompt : OriginType::Prompt);
161-
if (origin == OriginType::Prompt) {
159+
origin = RecoDecay::getCharmHadronOrigin(particlesMC, particle);
160+
if (origin == RecoDecay::OriginType::Prompt) {
162161
counterPrompt[iD]++;
163-
} else if (origin == OriginType::NonPrompt) {
162+
} else if (origin == RecoDecay::OriginType::NonPrompt) {
164163
counterNonPrompt[iD]++;
165164
}
166165
}
@@ -190,13 +189,13 @@ struct ValidationGenLevel {
190189
registry.fill(HIST("hPzDiffMotherDaughterGen"), pzDiff);
191190
registry.fill(HIST("hPDiffMotherDaughterGen"), pDiff);
192191
registry.fill(HIST("hPtDiffMotherDaughterGen"), ptDiff);
193-
if (origin == OriginType::Prompt) {
192+
if (origin == RecoDecay::OriginType::Prompt) {
194193
if (std::abs(particle.y()) < 0.5) {
195194
hPromptCharmHadronsPtDistr->Fill(whichHadron, particle.pt());
196195
}
197196
hPromptCharmHadronsYDistr->Fill(whichHadron, particle.y());
198197
hPromptCharmHadronsDecLenDistr->Fill(whichHadron, decayLength * 10000);
199-
} else if (origin == OriginType::NonPrompt) {
198+
} else if (origin == RecoDecay::OriginType::NonPrompt) {
200199
if (std::abs(particle.y()) < 0.5) {
201200
hNonPromptCharmHadronsPtDistr->Fill(whichHadron, particle.pt());
202201
}
@@ -288,7 +287,7 @@ struct ValidationRecLevel {
288287
whichHad = 6;
289288
}
290289
int whichOrigin = -1;
291-
if (cand2Prong.originMCRec() == OriginType::Prompt) {
290+
if (cand2Prong.originMCRec() == RecoDecay::OriginType::Prompt) {
292291
whichOrigin = 0;
293292
} else {
294293
whichOrigin = 1;
@@ -353,7 +352,7 @@ struct ValidationRecLevel {
353352
whichHad = 5;
354353
}
355354
int whichOrigin = -1;
356-
if (cand3Prong.originMCRec() == OriginType::Prompt) {
355+
if (cand3Prong.originMCRec() == RecoDecay::OriginType::Prompt) {
357356
whichOrigin = 0;
358357
} else {
359358
whichOrigin = 1;

‎PWGHF/Tasks/HFSelOptimisation.cxx‎

Lines changed: 4 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,6 @@
2323
using namespace o2;
2424
using namespace o2::framework;
2525
using namespace o2::framework::expressions;
26-
using namespace o2::aod::hf_cand;
2726

2827
namespace
2928
{
@@ -273,7 +272,7 @@ struct HfSelOptimisation {
273272
for (int iDecay{0}; iDecay < n2Prong; ++iDecay) {
274273
if (TESTBIT(cand2Prong.hfflag(), iDecay)) {
275274
if (std::abs(cand2Prong.flagMCMatchRec()) == BIT(iDecay)) {
276-
if (cand2Prong.originMCRec() == OriginType::Prompt) {
275+
if (cand2Prong.originMCRec() == RecoDecay::OriginType::Prompt) {
277276
isPrompt = true;
278277
switch (iDecay) {
279278
case o2::aod::hf_cand_prong2::DecayType::D0ToPiK:
@@ -283,7 +282,7 @@ struct HfSelOptimisation {
283282
testSelections2Prong<o2::aod::hf_cand_prong2::DecayType::JpsiToEE, 0>(cand2Prong, tracks);
284283
break;
285284
}
286-
} else if (cand2Prong.originMCRec() == OriginType::NonPrompt) {
285+
} else if (cand2Prong.originMCRec() == RecoDecay::OriginType::NonPrompt) {
287286
isNonPrompt = true;
288287
switch (iDecay) {
289288
case o2::aod::hf_cand_prong2::DecayType::D0ToPiK:
@@ -329,7 +328,7 @@ struct HfSelOptimisation {
329328
for (int iDecay{0}; iDecay < n3Prong; ++iDecay) {
330329
if (TESTBIT(cand3Prong.hfflag(), iDecay)) {
331330
if (std::abs(cand3Prong.flagMCMatchRec()) == BIT(iDecay)) {
332-
if (cand3Prong.originMCRec() == OriginType::Prompt) {
331+
if (cand3Prong.originMCRec() == RecoDecay::OriginType::Prompt) {
333332
isPrompt = true;
334333
switch (iDecay) {
335334
case o2::aod::hf_cand_prong3::DecayType::DPlusToPiKPi:
@@ -345,7 +344,7 @@ struct HfSelOptimisation {
345344
testSelections3Prong<o2::aod::hf_cand_prong3::DecayType::XicToPKPi, 0>(cand3Prong, tracks);
346345
break;
347346
}
348-
} else if (cand3Prong.originMCRec() == OriginType::NonPrompt) {
347+
} else if (cand3Prong.originMCRec() == RecoDecay::OriginType::NonPrompt) {
349348
isNonPrompt = true;
350349
switch (iDecay) {
351350
case o2::aod::hf_cand_prong3::DecayType::DPlusToPiKPi:

‎PWGHF/Tasks/taskD0.cxx‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -23,7 +23,6 @@
2323
using namespace o2;
2424
using namespace o2::framework;
2525
using namespace o2::framework::expressions;
26-
using namespace o2::aod::hf_cand;
2726
using namespace o2::aod::hf_cand_prong2;
2827
using namespace o2::analysis::hf_cuts_d0_topik;
2928

@@ -201,7 +200,7 @@ struct TaskD0 {
201200
registry.fill(HIST("hPtvsYRecSig_RecoPID"), ptRec, yRec);
202201
}
203202

204-
if (candidate.originMCRec() == OriginType::Prompt) {
203+
if (candidate.originMCRec() == RecoDecay::OriginType::Prompt) {
205204
registry.fill(HIST("hPtRecSigPrompt"), ptRec); // rec. level pT, prompt
206205
if (candidate.isRecoHFFlag() >= d_selectionHFFlag) {
207206
registry.fill(HIST("hPtvsYRecSigPrompt_RecoHFFlag"), ptRec, yRec);
@@ -314,7 +313,7 @@ struct TaskD0 {
314313
auto yGen = RecoDecay::y(array{particle.px(), particle.py(), particle.pz()}, RecoDecay::getMassPDG(particle.pdgCode()));
315314
registry.fill(HIST("hPtGen"), ptGen);
316315
registry.fill(HIST("hPtvsYGen"), ptGen, yGen);
317-
if (particle.originMCGen() == OriginType::Prompt) {
316+
if (particle.originMCGen() == RecoDecay::OriginType::Prompt) {
318317
registry.fill(HIST("hPtGenPrompt"), ptGen);
319318
registry.fill(HIST("hPtvsYGenPrompt"), ptGen, yGen);
320319
} else {

0 commit comments

Comments
 (0)