Skip to content

Commit 6900742

Browse files
jesgumjesgum
authored
[ALICE3] Add procces function to onTheFlyTracker.cxx for mcParticle with daughters (#14508)
1 parent 55af663 commit 6900742

File tree

5 files changed

+364
-37
lines changed

5 files changed

+364
-37
lines changed

‎ALICE3/Core/Decayer.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -69,7 +69,7 @@ class Decayer
6969
vy = pos[1] + rxyz * (mom[1] / track.getP());
7070
vz = pos[2] + rxyz * (mom[2] / track.getP());
7171
px = mom[0];
72-
py = mom[2];
72+
py = mom[1];
7373
} else {
7474
float sna, csa;
7575
o2::math_utils::CircleXYf_t circle;

‎ALICE3/Core/TrackUtilities.h‎

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,8 +33,10 @@ struct OTFParticle {
3333
float mE;
3434
float mVx, mVy, mVz;
3535
float mPx, mPy, mPz;
36+
bool mIsAlive;
3637

3738
// Setters
39+
void setIsAlive(bool isAlive) { mIsAlive = isAlive; }
3840
void setPDG(int pdg) { mPdgCode = pdg; }
3941
void setVxVyVz(float vx, float vy, float vz)
4042
{
@@ -52,6 +54,7 @@ struct OTFParticle {
5254

5355
// Getters
5456
int pdgCode() const { return mPdgCode; }
57+
int isAlive() const { return mIsAlive; }
5558
float vx() const { return mVx; }
5659
float vy() const { return mVy; }
5760
float vz() const { return mVz; }
@@ -97,7 +100,9 @@ void convertTLorentzVectorToO2Track(int pdgCode,
97100
/// \param o2track the address of the resulting TrackParCov
98101
/// \param pdg the pdg service
99102
template <typename PdgService>
100-
void convertOTFParticleToO2Track(const OTFParticle& particle, o2::track::TrackParCov& o2track, const PdgService& pdg)
103+
void convertOTFParticleToO2Track(const OTFParticle& particle,
104+
o2::track::TrackParCov& o2track,
105+
const PdgService& pdg)
101106
{
102107
static TLorentzVector tlv;
103108
tlv.SetPxPyPzE(particle.px(), particle.py(), particle.pz(), particle.e());

‎ALICE3/DataModel/OTFMCParticle.h‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,6 +32,9 @@ DECLARE_SOA_COLUMN(Eta, eta, float);
3232
DECLARE_SOA_COLUMN(Pt, pt, float);
3333
DECLARE_SOA_COLUMN(P, p, float);
3434
DECLARE_SOA_COLUMN(Y, y, float);
35+
DECLARE_SOA_COLUMN(IsAlive, isAlive, bool);
36+
DECLARE_SOA_COLUMN(IsPrimary, isPrimary, bool);
37+
3538
} // namespace otfmcparticle
3639

3740
DECLARE_SOA_TABLE_FULL(McParticlesWithDau, "McParticlesWithDau", "AOD", "MCPARTICLEWITHDAU",
@@ -56,6 +59,8 @@ DECLARE_SOA_TABLE_FULL(McParticlesWithDau, "McParticlesWithDau", "AOD", "MCPARTI
5659
otfmcparticle::Pt,
5760
otfmcparticle::P,
5861
otfmcparticle::Y,
62+
otfmcparticle::IsAlive,
63+
otfmcparticle::IsPrimary,
5964
mcparticle::PVector<mcparticle::Px, mcparticle::Py, mcparticle::Pz>,
6065
mcparticle::ProducedByGenerator<mcparticle::Flags>,
6166
mcparticle::FromBackgroundEvent<mcparticle::Flags>,

‎ALICE3/TableProducer/OTF/onTheFlyDecayer.cxx‎

Lines changed: 74 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -51,7 +51,7 @@ using namespace o2::framework;
5151
static constexpr int NumDecays = 7;
5252
static constexpr int NumParameters = 1;
5353
static constexpr int DefaultParameters[NumDecays][NumParameters]{{1}, {1}, {1}, {1}, {1}, {1}, {1}};
54-
static constexpr float VerySmall = 1e-7f;
54+
static constexpr float Tolerance = 1e-7f;
5555
static const std::vector<std::string> parameterNames{"enable"};
5656
static const std::vector<std::string> particleNames{"K0s",
5757
"Lambda",
@@ -102,6 +102,12 @@ struct OnTheFlyDecayer {
102102
histos.add("AntiXi/hGenAntiXi", "hGenAntiXi", kTH1D, {axisPt});
103103
histos.add("Omega/hGenOmega", "hGenOmega", kTH1D, {axisPt});
104104
histos.add("AntiOmega/hGenAntiOmega", "hGenAntiOmega", kTH1D, {axisPt});
105+
106+
histos.add("GeneratedElectron/hGenEl", "hGenEl", kTH1D, {axisPt});
107+
histos.add("GeneratedMuon/hGenMu", "hGenMu", kTH1D, {axisPt});
108+
histos.add("GeneratedPion/hGenPi", "hGenPi", kTH1D, {axisPt});
109+
histos.add("GeneratedKaon/hGenKa", "hGenKa", kTH1D, {axisPt});
110+
histos.add("GeneratedProton/hGenPr", "hGenPr", kTH1D, {axisPt});
105111
}
106112

107113
bool canDecay(const int pdgCode)
@@ -123,40 +129,60 @@ struct OnTheFlyDecayer {
123129
o2::upgrade::convertMCParticleToO2Track(particle, o2track, pdgDB);
124130
decayDaughters = decayer.decayParticle(pdgDB, o2track, particle.pdgCode());
125131
for (size_t idau{0}; idau < decayDaughters.size(); ++idau) {
126-
const auto& dau = decayDaughters[idau];
132+
o2::upgrade::OTFParticle& dau = decayDaughters[idau];
127133
o2::track::TrackParCov dauTrack;
128134
o2::upgrade::convertOTFParticleToO2Track(dau, dauTrack, pdgDB);
129135
if (canDecay(dau.pdgCode())) {
136+
dau.setIsAlive(false);
130137
std::vector<o2::upgrade::OTFParticle> cascadingDaughers = decayer.decayParticle(pdgDB, dauTrack, dau.pdgCode());
131138
for (const auto& daudau : cascadingDaughers) {
132139
decayDaughters.push_back(daudau);
133140
if (MaxNestedDecays < ++nDecays) {
134141
LOG(error) << "Seemingly stuck trying to perpetually decay products from pdg: " << particle.pdgCode();
135142
}
136143
}
144+
} else {
145+
dau.setIsAlive(true);
137146
}
138147
}
139148

140149
if (decayDaughters.empty()) {
141150
LOG(error) << "Attempted to decay " << particle.pdgCode() << " but resulting vector of daugthers were empty";
142151
continue;
143152
}
144-
}
145153

146-
if (particle.pdgCode() == kK0Short) {
147-
histos.fill(HIST("K0S/hGenK0S"), particle.pt());
148-
} else if (particle.pdgCode() == kLambda0) {
149-
histos.fill(HIST("Lambda/hGenLambda"), particle.pt());
150-
} else if (particle.pdgCode() == kLambda0Bar) {
151-
histos.fill(HIST("AntiLambda/hGenAntiLambda"), particle.pt());
152-
} else if (particle.pdgCode() == kXiMinus) {
153-
histos.fill(HIST("Xi/hGenXi"), particle.pt());
154-
} else if (particle.pdgCode() == kXiPlusBar) {
155-
histos.fill(HIST("AntiXi/hGenAntiXi"), particle.pt());
156-
} else if (particle.pdgCode() == kOmegaMinus) {
157-
histos.fill(HIST("Omega/hGenOmega"), particle.pt());
158-
} else if (particle.pdgCode() == kOmegaPlusBar) {
159-
histos.fill(HIST("AntiOmega/hGenAntiOmega"), particle.pt());
154+
switch (particle.pdgCode()) {
155+
case kK0Short:
156+
histos.fill(HIST("K0S/hGenK0S"), particle.pt());
157+
break;
158+
159+
case kLambda0:
160+
histos.fill(HIST("Lambda/hGenLambda"), particle.pt());
161+
break;
162+
163+
case kLambda0Bar:
164+
histos.fill(HIST("AntiLambda/hGenAntiLambda"), particle.pt());
165+
break;
166+
167+
case kXiMinus:
168+
histos.fill(HIST("Xi/hGenXi"), particle.pt());
169+
break;
170+
171+
case kXiPlusBar:
172+
histos.fill(HIST("AntiXi/hGenAntiXi"), particle.pt());
173+
break;
174+
175+
case kOmegaMinus:
176+
histos.fill(HIST("Omega/hGenOmega"), particle.pt());
177+
break;
178+
179+
case kOmegaPlusBar:
180+
histos.fill(HIST("AntiOmega/hGenAntiOmega"), particle.pt());
181+
break;
182+
183+
default:
184+
break;
185+
}
160186
}
161187

162188
int daughtersIdSlice[2];
@@ -178,13 +204,13 @@ struct OnTheFlyDecayer {
178204
float p = std::sqrt(particle.px() * particle.px() + particle.py() * particle.py() + particle.pz() * particle.pz());
179205
float y; // As https://github.com/AliceO2Group/AliceO2/blob/dev/Framework/Core/include/Framework/AnalysisDataModel.h#L1943
180206

181-
if ((p - particle.pz()) < VerySmall) {
207+
if ((p - particle.pz()) < Tolerance) {
182208
eta = (particle.pz() < 0.0f) ? -100.0f : 100.0f;
183209
} else {
184210
eta = 0.5f * std::log((p + particle.pz()) / (p - particle.pz()));
185211
}
186212

187-
if ((particle.e() - particle.pz()) < VerySmall) {
213+
if ((particle.e() - particle.pz()) < Tolerance) {
188214
y = (particle.pz() < 0.0f) ? -100.0f : 100.0f;
189215
} else {
190216
y = 0.5f * std::log((particle.e() + particle.pz()) / (particle.e() - particle.pz()));
@@ -196,7 +222,7 @@ struct OnTheFlyDecayer {
196222
particle.flags(), motherSpan, daughtersIdSlice, particle.weight(),
197223
particle.px(), particle.py(), particle.pz(), particle.e(),
198224
particle.vx(), particle.vy(), particle.vz(), particle.vt(),
199-
phi, eta, pt, p, y);
225+
phi, eta, pt, p, y, !canDecay(particle.pdgCode()), true);
200226
}
201227

202228
int daughtersIdSlice[2] = {-1, -1};
@@ -217,26 +243,51 @@ struct OnTheFlyDecayer {
217243
float p = std::sqrt(dau.px() * dau.px() + dau.py() * dau.py() + dau.pz() * dau.pz());
218244
float y; // Conditional as https://github.com/AliceO2Group/AliceO2/blob/dev/Framework/Core/include/Framework/AnalysisDataModel.h#L1943
219245

220-
if ((p - dau.pz()) < VerySmall) {
246+
if ((p - dau.pz()) < Tolerance) {
221247
eta = (dau.pz() < 0.0f) ? -100.0f : 100.0f;
222248
} else {
223249
eta = 0.5f * std::log((p + dau.pz()) / (p - dau.pz()));
224250
}
225251

226-
if ((dau.e() - dau.pz()) < VerySmall) {
252+
if ((dau.e() - dau.pz()) < Tolerance) {
227253
y = (dau.pz() < 0.0f) ? -100.0f : 100.0f;
228254
} else {
229255
y = 0.5f * std::log((dau.e() + dau.pz()) / (dau.e() - dau.pz()));
230256
}
231257

258+
switch (dau.pdgCode()) {
259+
case kElectron:
260+
histos.fill(HIST("GeneratedElectron/hGenEl"), pt);
261+
break;
262+
263+
case kMuonMinus:
264+
histos.fill(HIST("GeneratedMuon/hGenMu"), pt);
265+
break;
266+
267+
case kPiPlus:
268+
histos.fill(HIST("GeneratedPion/hGenPi"), pt);
269+
break;
270+
271+
case kKPlus:
272+
histos.fill(HIST("GeneratedKaon/hGenKa"), pt);
273+
break;
274+
275+
case kProton:
276+
histos.fill(HIST("GeneratedProton/hGenPr"), pt);
277+
break;
278+
279+
default:
280+
break;
281+
}
282+
232283
// TODO: Particle status code
233284
// TODO: Expression columns
234285
// TODO: vt
235286
tableMcParticlesWithDau(mother.mcCollisionId(), dau.pdgCode(), 1,
236287
mother.flags(), motherSpan, daughtersIdSlice, mother.weight(),
237288
dau.px(), dau.py(), dau.pz(), dau.e(),
238289
dau.vx(), dau.vy(), dau.vz(), mother.vt(),
239-
phi, eta, pt, p, y);
290+
phi, eta, pt, p, y, dau.isAlive(), false);
240291
}
241292
}
242293
}

0 commit comments

Comments
 (0)