Skip to content

Commit 83c483b

Browse files
MarcellocostiMarcellocosti
committed
[WIP] D+ matching
1 parent bb84d2f commit 83c483b

File tree

5 files changed

+397
-136
lines changed

5 files changed

+397
-136
lines changed

‎Common/Core/RecoDecay.h‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -789,6 +789,8 @@ struct RecoDecay {
789789
// Check that the number of direct daughters is not larger than the number of expected final daughters.
790790
if constexpr (!acceptIncompleteReco && !checkProcess) {
791791
if (particleMother.daughtersIds().back() - particleMother.daughtersIds().front() + 1 > static_cast<int>(N)) {
792+
// LOG(info) << "MC Rec: Rejected: too many direct daughters: " << particleMother.daughtersIds().back() - particleMother.daughtersIds().front() + 1
793+
// << " (expected " << N << " final)";
792794
// Printf("MC Rec: Rejected: too many direct daughters: %d (expected %ld final)", particleMother.daughtersIds().back() - particleMother.daughtersIds().front() + 1, N);
793795
return -1;
794796
}
@@ -801,7 +803,9 @@ struct RecoDecay {
801803
// }
802804
// printf("\n");
803805
// Check whether the number of actual final daughters is equal to the number of expected final daughters (i.e. the number of provided prongs).
806+
// Printf("MC Rec: Number of final daughters: %ld (expected %ld)", arrAllDaughtersIndex.size(), N);
804807
if (!acceptIncompleteReco && arrAllDaughtersIndex.size() != N) {
808+
// LOG(info) << "MC Rec: Rejected: incorrect number of final daughters: " << arrAllDaughtersIndex.size() << " (expected " << N << ")";
805809
// Printf("MC Rec: Rejected: incorrect number of final daughters: %ld (expected %ld)", arrAllDaughtersIndex.size(), N);
806810
return -1;
807811
}

‎Common/TableProducer/centralityTable.cxx‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -240,6 +240,7 @@ struct CentralityTable {
240240
mRunNumber = 0;
241241
listCalib.setObject(new TList);
242242
if (!produceHistograms.value) {
243+
LOG(info) << "Debug histograms are disabled. No histograms will be produced.";
243244
return;
244245
}
245246

@@ -259,6 +260,7 @@ struct CentralityTable {
259260
histos.addClone("FT0C/", "sel8FT0C/");
260261
histos.addClone("FT0A/", "sel8FT0A/");
261262

263+
LOG(info) << "Debug histograms are enabled. Histograms will be produced.";
262264
histos.print();
263265
}
264266

@@ -569,6 +571,7 @@ struct CentralityTable {
569571
mftInfo.mCalibrationStored = false;
570572
if (callst != nullptr) {
571573
if (produceHistograms) {
574+
LOG(info) << "Adding calibration list to the list of calibration lists";
572575
listCalib->Add(callst->Clone(Form("%i", bc.runNumber())));
573576
}
574577
LOGF(info, "Getting new histograms with %d run number for %d run number", mRunNumber, bc.runNumber());

‎PWGHF/Core/CorrelatedBkgs.h‎

Lines changed: 194 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -13,25 +13,40 @@
1313
/// \brief Definitions of channels for correlated bkgs
1414
/// \author Marcello Di Costanzo <marcello.di.costanzo@cern.ch>, Polytechnic of Turin and INFN
1515

16+
#include <unordered_map>
17+
#include <array>
18+
#include <variant>
19+
1620
#ifndef PWGHF_CORE_CORRELATEDBKGTAGGER_H_
1721
#define PWGHF_CORE_CORRELATEDBKGTAGGER_H_
1822

23+
using namespace o2::constants::physics;
24+
1925
namespace o2::hf_corrbkg
2026
{
2127

28+
enum FinalStatesDMesons {
29+
KKPi = 0,
30+
KPiPi,
31+
KPiPiPi0,
32+
PiPiPi,
33+
PiPiPiPi0,
34+
NFinalStatesDMesons
35+
};
36+
2237
// D± → K± K∓ π±
2338
enum DecayChannelResoDplus {
24-
PhiPi = 0,
25-
K0starK,
39+
PhiPiDplus = 0,
40+
K0starKDplus,
2641
// K01430K,
2742
// RhoPi,
2843
// f2_1270Pi
2944
};
3045

3146
// Ds± → K± K∓ π±
3247
enum DecayChannelResoDs {
33-
PhiPi = DecayChannelResoDplus::PhiPi,
34-
K0starK,
48+
PhiPiDs = 0,
49+
K0starKDs,
3550
// PhiRho,
3651
// f2_1270Pi,
3752
// K0starPi,
@@ -60,15 +75,183 @@ namespace o2::hf_corrbkg
6075
// PhiP,
6176
// };
6277

63-
// Function to associate vectors with enum values
64-
std::vector<int> getFinalDaughters(DecayType finalState) {
65-
switch(finalState) {
66-
case DecayChannelCorrBkg::DstarToPiKPiBkg:
67-
return std::array{+kPiPlus, +kPiPlus, -kKPlus};
68-
default:
69-
return std::array{};
78+
79+
// using FinalStateVariant = std::variant<
80+
// std::array<int, 3>,
81+
// std::array<int, 4>
82+
// >;
83+
84+
// std::unordered_map<int, FinalStateVariant> finalStates =
85+
// {
86+
// {FinalStatesDMesons::KKPi, std::array<int,3>{+kKPlus, +kKPlus, +kPiPlus}},
87+
// {FinalStatesDMesons::KPiPiPi0, std::array<int,4>{+kKPlus, +kKPlus, +kPiPlus, +kPi0}}
88+
// };
89+
90+
std::unordered_map<FinalStatesDMesons, std::pair<std::array<int,3>, bool>> finalStates =
91+
{
92+
{FinalStatesDMesons::KKPi, {std::array<int,3>{-kKPlus, +kKPlus, +kPiPlus}, false}},
93+
{FinalStatesDMesons::KPiPi, {std::array<int,3>{+kPiPlus, -kKPlus, +kPiPlus}, false}},
94+
{FinalStatesDMesons::KPiPiPi0, {std::array<int,3>{+kPiPlus, -kKPlus, +kPiPlus}, true}},
95+
{FinalStatesDMesons::PiPiPi, {std::array<int,3>{+kPiMinus, +kPiPlus, +kPiPlus}, false}},
96+
{FinalStatesDMesons::PiPiPiPi0, {std::array<int,3>{+kPiMinus, +kPiPlus, +kPiPlus}, true}},
97+
};
98+
99+
100+
template <bool matchKinkedDecayTopology = false, bool matchInteractionsWithMaterial = false>
101+
int matchFinalStateDMeson(aod::McParticles const& mcParticles, auto arrayDaughters, int pdg, int8_t* flag, int8_t* sign, int depth, int8_t* nKinkedTracks, int8_t* nInteractionsWithMaterial) {
102+
int indexRec = -1;
103+
LOG(info) << "matchFinalStateDMeson: " << pdg;
104+
for (const std::pair<int, std::pair<std::array<int,3>, bool>>& finalState : finalStates) {
105+
if (finalState.second.second) { // 4 prong decay, partly reconstructed
106+
LOG(info) << "[4prong] Taking channel with final state: " << finalState.second.first[0] << " " << finalState.second.first[1] << " " << finalState.second.first[2] << " " << finalState.second.first[3] << ", has 4 daughters: " << finalState.second.second;
107+
LOG(info) << "[4prong] Trying 4-prong decay ... ";
108+
LOG(info) << "[4prong] Matching channel with final state: " << finalState.second.first[0] << " " << finalState.second.first[1] << " " << finalState.second.first[2];
109+
if constexpr (matchKinkedDecayTopology && matchInteractionsWithMaterial) {
110+
indexRec = RecoDecay::getMatchedMCRec<false, false, true, true, true>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nKinkedTracks, nInteractionsWithMaterial);
111+
} else if constexpr (matchKinkedDecayTopology && !matchInteractionsWithMaterial) {
112+
indexRec = RecoDecay::getMatchedMCRec<false, false, true, true, false>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nKinkedTracks);
113+
} else if constexpr (!matchKinkedDecayTopology && matchInteractionsWithMaterial) {
114+
indexRec = RecoDecay::getMatchedMCRec<false, false, true, false, true>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nullptr, nInteractionsWithMaterial);
115+
} else {
116+
indexRec = RecoDecay::getMatchedMCRec<false, false, true, false, false>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, &sign, depth);
117+
}
118+
} else {
119+
LOG(info) << "[3prong] Trying 3-prong decay ... ";
120+
LOG(info) << "[3prong] Taking channel with final state: " << finalState.second.first[0] << " " << finalState.second.first[1] << " " << finalState.second.first[2] << ", has 4 daughters: " << finalState.second.second;
121+
if constexpr (matchKinkedDecayTopology && matchInteractionsWithMaterial) {
122+
indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, true>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nKinkedTracks, nInteractionsWithMaterial);
123+
} else if constexpr (matchKinkedDecayTopology && !matchInteractionsWithMaterial) {
124+
indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, false>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nKinkedTracks);
125+
} else if constexpr (!matchKinkedDecayTopology && matchInteractionsWithMaterial) {
126+
indexRec = RecoDecay::getMatchedMCRec<false, false, false, false, true>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, sign, depth, nullptr, nInteractionsWithMaterial);
127+
} else {
128+
indexRec = RecoDecay::getMatchedMCRec<false, false, false, false, false>(mcParticles, arrayDaughters, pdg, finalState.second.first, true, &sign, depth);
129+
}
130+
}
131+
if (indexRec > -1) {
132+
LOG(info) << "Found channel with final state: " << finalState.second.first[0] << " " << finalState.second.first[1] << " " << finalState.second.first[2] << ", has 4 prongs: " << finalState.second.second;
133+
*flag = (*sign) * (1 << finalState.first);
134+
break;
135+
}
70136
}
137+
return indexRec;
71138
}
139+
140+
141+
142+
143+
144+
// Function to associate vectors with enum values
145+
// std::vector<int> getFinalDaughters(DecayType finalState) {
146+
// switch(finalState) {
147+
// case DecayChannelCorrBkg::DstarToPiKPiBkg:
148+
// return std::array{+kPiPlus, +kPiPlus, -kKPlus};
149+
// default:
150+
// return std::array{};
151+
// }
152+
// }
153+
154+
// std::vector<Pdg> getFinalStateDMeson(int index) {
155+
// switch (index) {
156+
// case KKPi:
157+
// return std::vector<Pdg>{+kKPlus, -kKPlus, +kPiPlus};
158+
// case KPiPi:
159+
// return std::vector<Pdg>{+kKPlus, +kPiPlus, +kPiPlus};
160+
// case KPiPiPi0:
161+
// return std::vector<Pdg>{+kKPlus, +kPiPlus, +kPiPlus, +kPi0};
162+
// case PiPiPi:
163+
// return std::vector<Pdg>{+kPiPlus, +kPiPlus, +kPiPlus};
164+
// case PiPiPiPi0:
165+
// return std::vector<Pdg>{+kPiPlus, +kPiPlus, +kPiPlus, +kPi0};
166+
// default:
167+
// throw std::out_of_range("Invalid index for final state D meson");
168+
// }
169+
// }
170+
171+
172+
// template <int N>
173+
// std::array<int, N> getFinalStateDMeson(int index) {
174+
// if constexpr (N == 3) {
175+
// switch (index) {
176+
// case FinalStatesDMesons::KKPi:
177+
// return std::array<int, N>{+kKPlus, -kKPlus, +kPiPlus};
178+
// case FinalStatesDMesons::KPiPi:
179+
// return std::array<int, N>{+kKPlus, +kPiPlus, +kPiPlus};
180+
// case FinalStatesDMesons::PiPiPi:
181+
// return std::array<int, N>{+kPiPlus, +kPiPlus, +kPiPlus};
182+
// default:
183+
// throw std::out_of_range("Invalid index for final state D meson");
184+
// }
185+
// } else if constexpr (N == 4) {
186+
// switch (index) {
187+
// case FinalStatesDMesons::KPiPiPi0:
188+
// return std::array<int, N>{+kKPlus, +kPiPlus, +kPiPlus, +kPi0};
189+
// case FinalStatesDMesons::PiPiPiPi0:
190+
// return std::array<int, N>{+kPiPlus, +kPiPlus, +kPiPlus, +kPi0};
191+
// default:
192+
// throw std::out_of_range("Invalid index for final state D meson");
193+
// }
194+
// } else {
195+
// throw std::invalid_argument("Invalid size for final state D meson");
196+
// }
197+
// }
198+
199+
// int getFinalStateParts(int index) {
200+
// switch (index) {
201+
// case FinalStatesDMesons::KKPi:
202+
// return 3;
203+
// case FinalStatesDMesons::KPiPi:
204+
// return 3;
205+
// case FinalStatesDMesons::KPiPiPi0:
206+
// return 4;
207+
// case FinalStatesDMesons::PiPiPi:
208+
// return 3;
209+
// case FinalStatesDMesons::PiPiPiPi0:
210+
// return 4;
211+
// default:
212+
// throw std::out_of_range("Invalid index for final state D meson");
213+
// }
214+
// }
215+
216+
// template <bool matchKinkedDecayTopology = false, bool matchInteractionsWithMaterial = false>
217+
// int matchFinalStateDMeson(int iChannel, aod::McParticles const& mcParticles, auto arrayDaughters, int pdg, int8_t* sign, int depth, int8_t* nKinkedTracks, int8_t* nInteractionsWithMaterial) {
218+
// int indexRec = -1;
219+
// int nFinalStateParts = getFinalStateParts(iChannel);
220+
// if (nFinalStateParts == 3) {
221+
// auto finalState = getFinalStateDMeson<3>(iChannel);
222+
// if constexpr (matchKinkedDecayTopology && matchInteractionsWithMaterial) {
223+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, true>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nKinkedTracks, nInteractionsWithMaterial);
224+
// } else if constexpr (matchKinkedDecayTopology && !matchInteractionsWithMaterial) {
225+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, false>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nKinkedTracks);
226+
// } else if constexpr (!matchKinkedDecayTopology && matchInteractionsWithMaterial) {
227+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, false, true>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nullptr, nInteractionsWithMaterial);
228+
// } else {
229+
// indexRec = RecoDecay::getMatchedMCRec(mcParticles, arrayDaughters, pdg, finalState, true, &sign, depth);
230+
// }
231+
// if (indexRec > -1) {
232+
// LOG(info) << "Matched final state: " << iChannel;
233+
// return indexRec;
234+
// }
235+
// } else if (nFinalStateParts == 4) {
236+
// auto finalState = getFinalStateDMeson<4>(iChannel);
237+
// if constexpr (matchKinkedDecayTopology && matchInteractionsWithMaterial) {
238+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, true>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nKinkedTracks, nInteractionsWithMaterial);
239+
// } else if constexpr (matchKinkedDecayTopology && !matchInteractionsWithMaterial) {
240+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, true, false>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nKinkedTracks);
241+
// } else if constexpr (!matchKinkedDecayTopology && matchInteractionsWithMaterial) {
242+
// indexRec = RecoDecay::getMatchedMCRec<false, false, false, false, true>(mcParticles, arrayDaughters, pdg, finalState, true, sign, depth, nullptr, nInteractionsWithMaterial);
243+
// } else {
244+
// indexRec = RecoDecay::getMatchedMCRec(mcParticles, arrayDaughters, pdg, finalState, true, &sign, depth);
245+
// }
246+
// if (indexRec > -1) {
247+
// LOG(info) << "Matched final state: " << iChannel;
248+
// return indexRec;
249+
// }
250+
// } else {
251+
// throw std::invalid_argument("Invalid size for final state D meson");
252+
// }
253+
// }
254+
72255
} // namespace o2::hf_corrbkg
73256

74257
#endif // PWGHF_CORE_CORRELATEDBKGTAGGER_H_

0 commit comments

Comments
 (0)