Skip to content

Commit ddc75d2

Browse files
MarcellocostiMarcellocosti
committed
Factorize function to change Pi0 sign for antiparticles
1 parent fc47508 commit ddc75d2

File tree

4 files changed

+88
-82
lines changed

4 files changed

+88
-82
lines changed

‎PWGHF/TableProducer/candidateCreator2Prong.cxx‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -798,6 +798,7 @@ struct HfCandidateCreator2ProngExpressions {
798798
if (indexRec > -1) {
799799
auto motherParticle = mcParticles.rawIteratorAt(indexRec);
800800
std::array<int, 3> finalStateParts2ProngAll = std::array{finalState[0], finalState[1], finalState[2]};
801+
convertPi0ToAntiPi0(motherParticle.pdgCode(), finalStateParts2ProngAll);
801802
if (!RecoDecay::isMatchedMCGen(mcParticles, motherParticle, Pdg::kD0, finalStateParts2ProngAll, true, &sign, FinalStateDepth)) {
802803
indexRec = -1; // Reset indexRec if the generated decay does not match the reconstructed one does not match the reconstructed one
803804
}
@@ -820,7 +821,6 @@ struct HfCandidateCreator2ProngExpressions {
820821
flag = sign * (1 << chn);
821822

822823
// Flag the resonant decay channel
823-
int resoMaxDepth = 1;
824824
std::vector<int> arrResoDaughIndex = {};
825825
RecoDecay::getDaughters(mcParticles.rawIteratorAt(indexRec), &arrResoDaughIndex, std::array{0}, ResoDepth);
826826
std::array<int, NDaughtersResonant> arrPDGDaugh = {};
@@ -829,7 +829,7 @@ struct HfCandidateCreator2ProngExpressions {
829829
auto daughI = mcParticles.rawIteratorAt(arrResoDaughIndex[iProng]);
830830
arrPDGDaugh[iProng] = daughI.pdgCode();
831831
}
832-
flagResonantDecay(Pdg::kD0, &channel, arrPDGDaugh);
832+
channel = flagResonantDecay(Pdg::kD0, arrPDGDaugh);
833833
}
834834
break;
835835
}

‎PWGHF/TableProducer/candidateCreator3Prong.cxx‎

Lines changed: 16 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -957,7 +957,7 @@ struct HfCandidateCreator3ProngExpressions {
957957
if (pdg == Pdg::kDStar) {
958958
depth = MaxDepth+1; // D0 resonant decays are active
959959
}
960-
auto finalStates = getDecayChannel3Prong(pdg);
960+
auto finalStates = getDecayChannelMain(pdg);
961961
for (const auto& [chn, finalState] : finalStates) {
962962
std::array<int, 3> finalStateParts3Prong = std::array{finalState[0], finalState[1], finalState[2]};
963963
if (finalState.size() > 3) { // Partly Reco 4-prong decays
@@ -975,25 +975,13 @@ struct HfCandidateCreator3ProngExpressions {
975975
auto motherParticle = mcParticles.rawIteratorAt(indexRec);
976976
if (finalState.size() == 4) { // Check if the final state has 4 particles
977977
std::array<int, 4> finalStateParts3ProngAll = std::array{finalState[0], finalState[1], finalState[2], finalState[3]};
978-
if (sign < 0) {
979-
for (auto& part : finalStateParts3ProngAll) {
980-
if (part == kPi0) {
981-
part = -part; // The Pi0 pdg code does not change between particle and antiparticle
982-
}
983-
}
984-
}
978+
convertPi0ToAntiPi0(motherParticle.pdgCode(), finalStateParts3ProngAll);
985979
if (!RecoDecay::isMatchedMCGen(mcParticles, motherParticle, pdg, finalStateParts3ProngAll, true, &sign, depth)) {
986980
indexRec = -1; // Reset indexRec if the generated decay does not match the reconstructed one is not matched
987981
}
988982
} else if (finalState.size() == 5) { // Check if the final state has 5 particles
989983
std::array<int, 5> finalStateParts3ProngAll = std::array{finalState[0], finalState[1], finalState[2], finalState[3], finalState[4]};
990-
if (sign < 0) {
991-
for (auto& part : finalStateParts3ProngAll) {
992-
if (part == kPi0) {
993-
part = -part; // The Pi0 pdg code does not change between particle and antiparticle
994-
}
995-
}
996-
}
984+
convertPi0ToAntiPi0(motherParticle.pdgCode(), finalStateParts3ProngAll);
997985
if (!RecoDecay::isMatchedMCGen(mcParticles, motherParticle, pdg, finalStateParts3ProngAll, true, &sign, depth)) {
998986
indexRec = -1; // Reset indexRec if the generated decay does not match the reconstructed one is not matched
999987
}
@@ -1017,19 +1005,26 @@ struct HfCandidateCreator3ProngExpressions {
10171005
flag = sign * chn;
10181006

10191007
// Flag the resonant decay channel
1020-
int resoMaxDepth = 1;
1021-
if (std::abs(mcParticles.rawIteratorAt(indexRec).pdgCode()) == Pdg::kDStar) {
1022-
resoMaxDepth = 2; // Flag D0 resonances
1023-
}
10241008
std::vector<int> arrResoDaughIndex = {};
1025-
RecoDecay::getDaughters(mcParticles.rawIteratorAt(indexRec), &arrResoDaughIndex, std::array{0}, resoMaxDepth);
1009+
if (pdg == Pdg::kDStar) {
1010+
std::vector<int> arrResoDaughIndexDStar = {};
1011+
RecoDecay::getDaughters(mcParticles.rawIteratorAt(indexRec), &arrResoDaughIndexDStar, std::array{0}, ResoMaxDepth);
1012+
for (size_t iDaug = 0; iDaug < arrResoDaughIndexDStar.size(); iDaug++) {
1013+
auto daughDstar = mcParticles.rawIteratorAt(arrResoDaughIndexDStar[iDaug]);
1014+
if (std::abs(daughDstar.pdgCode()) == Pdg::kD0 || std::abs(daughDstar.pdgCode()) == Pdg::kDPlus) {
1015+
break;
1016+
}
1017+
}
1018+
} else {
1019+
RecoDecay::getDaughters(mcParticles.rawIteratorAt(indexRec), &arrResoDaughIndex, std::array{0}, ResoMaxDepth);
1020+
}
10261021
std::array<int, NDaughtersResonant> arrPDGDaugh = {};
10271022
if (arrResoDaughIndex.size() == NDaughtersResonant) {
10281023
for (auto iProng = 0u; iProng < NDaughtersResonant; ++iProng) {
10291024
auto daughI = mcParticles.rawIteratorAt(arrResoDaughIndex[iProng]);
10301025
arrPDGDaugh[iProng] = daughI.pdgCode();
10311026
}
1032-
flagResonantDecay<true>(pdg, &channel, arrPDGDaugh);
1027+
channel = flagResonantDecay<true>(pdg, arrPDGDaugh);
10331028
}
10341029
break; // Exit loop if a match is found
10351030
}

‎PWGHF/Utils/utilsMcGen.h‎

Lines changed: 7 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -42,7 +42,6 @@ void fillMcMatchGen2Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
4242
constexpr std::size_t NDaughtersResonant{2u};
4343

4444
// Match generated particles.
45-
int maxDepth = 2; // Default depth for matching
4645
for (const auto& particle : mcParticlesPerMcColl) {
4746
int8_t flag = 0;
4847
int8_t origin = 0;
@@ -62,6 +61,7 @@ void fillMcMatchGen2Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
6261
for (const auto& [chn, finalState] : o2::hf_decay::hf_cand_2prong::daughtersD0Main) {
6362
if (finalState.size() == 3) { // Partly Reco 3-prong decays
6463
std::array<int, 3> finalStateParts = std::array{finalState[0], finalState[1], finalState[2]};
64+
o2::hf_decay::convertPi0ToAntiPi0(particle.pdgCode(), finalStateParts);
6565
matched = RecoDecay::isMatchedMCGen(mcParticles, particle, Pdg::kD0, finalStateParts, true, &sign, MaxDepth);
6666
} else if (finalState.size() == 2) { // Fully Reco 2-prong decays
6767
std::array<int, 2> finalStateParts = std::array{finalState[0], finalState[1]};
@@ -74,7 +74,6 @@ void fillMcMatchGen2Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
7474
flag = sign * (1 << chn);
7575

7676
// Flag the resonant decay channel
77-
int resoMaxDepth = 1;
7877
std::vector<int> arrResoDaughIndex = {};
7978
RecoDecay::getDaughters(particle, &arrResoDaughIndex, std::array{0}, ResoMaxDepth);
8079
std::array<int, NDaughtersResonant> arrPDGDaugh = {};
@@ -83,7 +82,7 @@ void fillMcMatchGen2Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
8382
auto daughI = mcParticles.rawIteratorAt(arrResoDaughIndex[iProng]);
8483
arrPDGDaugh[iProng] = daughI.pdgCode();
8584
}
86-
o2::hf_decay::flagResonantDecay(Pdg::kD0, &channel, arrPDGDaugh);
85+
channel = o2::hf_decay::flagResonantDecay(Pdg::kD0, arrPDGDaugh);
8786
}
8887
break;
8988
}
@@ -164,26 +163,14 @@ void fillMcMatchGen3Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
164163

165164
std::vector<int> arrAllDaughtersIndex;
166165
for (const auto& [chn, finalState] : finalStates) {
167-
if (finalState.size() == 5) { // Partly Reco 3-prong decays
166+
if (finalState.size() == 5) { // Partly Reco 3-prong decays from 5-prong decays
168167
std::array<int, 5> finalStateParts = std::array{finalState[0], finalState[1], finalState[2], finalState[3], finalState[4]};
169-
if (particle.pdgCode() < 0) {
170-
for (auto& part : finalStateParts) {
171-
if (part == kPi0) {
172-
part = -part; // The Pi0 pdg code does not change between particle and antiparticle
173-
}
174-
}
175-
}
168+
o2::hf_decay::convertPi0ToAntiPi0(particle.pdgCode(), finalStateParts);
176169
RecoDecay::getDaughters<false>(particle, &arrAllDaughtersIndex, finalStateParts, maxDepth);
177170
matched = RecoDecay::isMatchedMCGen(mcParticles, particle, motherPdgCode, finalStateParts, true, &sign, -1);
178-
} else if (finalState.size() == 4) { // Partly Reco 3-prong decays
171+
} else if (finalState.size() == 4) { // Partly Reco 3-prong decays from 4-prong decays
179172
std::array<int, 4> finalStateParts = std::array{finalState[0], finalState[1], finalState[2], finalState[3]};
180-
if (particle.pdgCode() < 0) {
181-
for (auto& part : finalStateParts) {
182-
if (part == kPi0) {
183-
part = -part; // The Pi0 pdg code does not change between particle and antiparticle
184-
}
185-
}
186-
}
173+
o2::hf_decay::convertPi0ToAntiPi0(particle.pdgCode(), finalStateParts);
187174
RecoDecay::getDaughters<false>(particle, &arrAllDaughtersIndex, finalStateParts, maxDepth);
188175
matched = RecoDecay::isMatchedMCGen(mcParticles, particle, motherPdgCode, finalStateParts, true, &sign, -1);
189176
} else if (finalState.size() == 3) { // Fully Reco 3-prong decays
@@ -197,7 +184,6 @@ void fillMcMatchGen3Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
197184
if (matched) {
198185
flag = sign * chn;
199186
// Flag the resonant decay channel
200-
int resoMaxDepth = 1;
201187
std::vector<int> arrResoDaughIndex = {};
202188
if (std::abs(motherPdgCode) == Pdg::kDStar) {
203189
std::vector<int> arrResoDaughIndexDStar = {};
@@ -220,7 +206,7 @@ void fillMcMatchGen3Prong(T const& mcParticles, U const& mcParticlesPerMcColl, V
220206
auto daughI = mcParticles.rawIteratorAt(arrResoDaughIndex[iProng]);
221207
arrPDGDaugh[iProng] = daughI.pdgCode();
222208
}
223-
o2::hf_decay::flagResonantDecay<true>(motherPdgCode, &channel, arrPDGDaugh);
209+
channel = o2::hf_decay::flagResonantDecay<true>(motherPdgCode, arrPDGDaugh);
224210
}
225211
break; // Exit loop if a match is found
226212
}

‎PWGHF/Utils/utilsMcMatching.h‎

Lines changed: 63 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -164,34 +164,13 @@ inline std::unordered_map<DecayChannelMain, const std::vector<int>> getDecayChan
164164
}
165165
} // namespace hf_cand_3prong
166166

167-
168-
{
169-
switch (pdgMother) {
170-
case o2::constants::physics::Pdg::kDPlus:
171-
return DaughtersDPlusResonant;
172-
case o2::constants::physics::Pdg::kDS:
173-
return DaughtersDsResonant;
174-
case o2::constants::physics::Pdg::kDStar:
175-
return DaughtersDstarResonant;
176-
case o2::constants::physics::Pdg::kLambdaCPlus:
177-
return DaughtersLcResonant;
178-
case o2::constants::physics::Pdg::kXiCPlus:
179-
return DaughtersXiCResonant;
180-
default:
181-
LOG(error) << "Unknown PDG code for 3-prong final states: " << pdgMother;
182-
return {};
183-
}
184-
}
185-
186-
} // namespace hf_cand_3prong
187-
188167
/// Perform the matching for a single resonant channel
189168
/// \tparam N size of the array of daughter PDG codes
190169
/// \param arrPdgResoChn array of daughter indices
191170
/// \param arrPdgDaugs array of PDG codes for the resonant decay
192171
/// \return true if the resonant channel is matched, false otherwise
193172
template <std::size_t N>
194-
inline bool checkResonantDecay(std::array<int, N> const& arrPdgResoChn, std::array<int, N> arrPdgDaugs)
173+
inline bool checkDecayChannel(std::array<int, N> const& arrPdgResoChn, std::array<int, N> arrPdgDaugs)
195174
{
196175
for (std::size_t i = 0; i < N; i++) {
197176
bool findDaug = false;
@@ -217,24 +196,70 @@ inline bool checkResonantDecay(std::array<int, N> const& arrPdgResoChn, std::arr
217196
/// \param channel decay channel flag to be set
218197
/// \param arrDaughPdgs array of daughter PDG codes
219198
template <bool is3Prong = false, std::size_t N>
220-
inline void flagResonantDecay(int motherPdg, int8_t* channel, std::array<int, N> const& arrDaughPdgs)
199+
inline int8_t flagResonantDecay(int motherPdg, std::array<int, N> const& arrDaughPdgs)
221200
{
222-
if constexpr (is3Prong) {
223-
std::unordered_map<o2::hf_decay::hf_cand_3prong::DecayChannelResonant, const std::array<int, 2>> resoStates = o2::hf_decay::hf_cand_3prong::getResoChannels3Prong(motherPdg);
224-
for (const auto& [flag, pdgCodes] : resoStates) {
225-
if (o2::hf_decay::checkResonantDecay(arrDaughPdgs, pdgCodes)) {
226-
*channel = flag;
227-
break;
201+
switch (motherPdg) {
202+
case o2::constants::physics::Pdg::kD0:
203+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_2prong::daughtersD0Resonant) {
204+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
205+
return flag;
206+
}
228207
}
229-
}
230-
} else {
231-
if (motherPdg != o2::constants::physics::Pdg::kD0) {
232-
return;
233-
}
234-
for (const auto& [flag, pdgCodes] : hf_cand_2prong::DaughtersD0Resonant) {
235-
if (o2::hf_decay::checkResonantDecay(arrDaughPdgs, pdgCodes)) {
236-
*channel = flag;
237-
break;
208+
break;
209+
case o2::constants::physics::Pdg::kDPlus:
210+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_3prong::daughtersDplusResonant) {
211+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
212+
return flag;
213+
}
214+
}
215+
break;
216+
case o2::constants::physics::Pdg::kDS:
217+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_3prong::daughtersDsResonant) {
218+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
219+
return flag;
220+
}
221+
}
222+
break;
223+
case o2::constants::physics::Pdg::kDStar:
224+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_3prong::daughtersDstarResonant) {
225+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
226+
return flag;
227+
}
228+
}
229+
break;
230+
case o2::constants::physics::Pdg::kLambdaCPlus:
231+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_3prong::daughtersLcResonant) {
232+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
233+
return flag;
234+
}
235+
}
236+
break;
237+
case o2::constants::physics::Pdg::kXiCPlus:
238+
for (const auto& [flag, pdgCodes] : o2::hf_decay::hf_cand_3prong::daughtersXicResonant) {
239+
if (o2::hf_decay::checkDecayChannel(arrDaughPdgs, pdgCodes)) {
240+
return flag;
241+
}
242+
}
243+
break;
244+
default:
245+
LOG(fatal) << "Unknown PDG code for 3-prong final states: " << motherPdg;
246+
return {};
247+
}
248+
return 0;
249+
}
250+
251+
/// Convert Pi0 to AntiPi0 in the array of PDG codes for antiparticles since
252+
/// the Pi0s stemming from antiparticle decays keep the pdg code 111.
253+
/// \tparam N size of the array of PDG codes
254+
/// \param partPdgCode PDG code of the particle
255+
/// \param arrPdgIndexes array of PDG codes to be modified
256+
template <std::size_t N>
257+
inline void convertPi0ToAntiPi0(const int partPdgCode, std::array<int, N>& arrPdgIndexes)
258+
{
259+
if (partPdgCode < 0) {
260+
for (auto& part : arrPdgIndexes) { // o2-linter: disable=const-ref-in-for-loop (int elements)
261+
if (part == kPi0) {
262+
part = -part; // The Pi0 pdg code does not change between particle and antiparticle
238263
}
239264
}
240265
}

0 commit comments

Comments
 (0)