Skip to content

Commit a387656

Browse files
jesgumjesgum
authored
[ALICE3] clean up otf decayer (#13945)
1 parent 8e36764 commit a387656

File tree

1 file changed

+57
-51
lines changed

1 file changed

+57
-51
lines changed

‎ALICE3/TableProducer/OTF/onTheFlyTracker.cxx‎

Lines changed: 57 additions & 51 deletions
Original file line numberDiff line numberDiff line change
@@ -575,52 +575,52 @@ struct OnTheFlyTracker {
575575
/// \param xiDecayVertex the address of the xi decay vertex
576576
/// \param laDecayVertex the address of the la decay vertex
577577
template <typename McParticleType>
578-
void decayParticle(McParticleType particle, o2::track::TrackParCov track, std::vector<TLorentzVector>& decayDaughters, std::vector<double>& xiDecayVertex, std::vector<double>& laDecayVertex)
578+
void decayCascade(McParticleType particle, o2::track::TrackParCov track, std::vector<TLorentzVector>& decayDaughters, std::vector<double>& xiDecayVertex, std::vector<double>& laDecayVertex)
579579
{
580-
double u = rand.Uniform(0, 1);
581-
double xi_mass = o2::constants::physics::MassXiMinus;
582-
double la_mass = o2::constants::physics::MassLambda;
583-
double pi_mass = o2::constants::physics::MassPionCharged;
584-
double pr_mass = o2::constants::physics::MassProton;
585-
586-
double xi_gamma = 1 / std::sqrt(1 + (particle.p() * particle.p()) / (xi_mass * xi_mass));
587-
double xi_ctau = 4.91 * xi_gamma;
588-
double xi_rxyz = (-xi_ctau * std::log(1 - u));
580+
const double uXi = rand.Uniform(0, 1);
581+
const double ctauXi = 4.91; // cm
582+
const double betaGammaXi = particle.p() / o2::constants::physics::MassXiMinus;
583+
const double rxyzXi = (-betaGammaXi * ctauXi * std::log(1 - uXi));
584+
589585
float sna, csa;
590-
o2::math_utils::CircleXYf_t xi_circle;
591-
track.getCircleParams(magneticField, xi_circle, sna, csa);
592-
double xi_rxy = xi_rxyz / std::sqrt(1. + track.getTgl() * track.getTgl());
593-
double theta = xi_rxy / xi_circle.rC;
594-
double newX = ((particle.vx() - xi_circle.xC) * std::cos(theta) - (particle.vy() - xi_circle.yC) * std::sin(theta)) + xi_circle.xC;
595-
double newY = ((particle.vy() - xi_circle.yC) * std::cos(theta) + (particle.vx() - xi_circle.xC) * std::sin(theta)) + xi_circle.yC;
596-
double newPx = particle.px() * std::cos(theta) - particle.py() * std::sin(theta);
597-
double newPy = particle.py() * std::cos(theta) + particle.px() * std::sin(theta);
586+
o2::math_utils::CircleXYf_t circleXi;
587+
track.getCircleParams(magneticField, circleXi, sna, csa);
588+
const double rxyXi = rxyzXi / std::sqrt(1. + track.getTgl() * track.getTgl());
589+
const double theta = rxyXi / circleXi.rC;
590+
const double newX = ((particle.vx() - circleXi.xC) * std::cos(theta) - (particle.vy() - circleXi.yC) * std::sin(theta)) + circleXi.xC;
591+
const double newY = ((particle.vy() - circleXi.yC) * std::cos(theta) + (particle.vx() - circleXi.xC) * std::sin(theta)) + circleXi.yC;
592+
const double newPx = particle.px() * std::cos(theta) - particle.py() * std::sin(theta);
593+
const double newPy = particle.py() * std::cos(theta) + particle.px() * std::sin(theta);
594+
const double newE = std::sqrt(o2::constants::physics::MassXiMinus * o2::constants::physics::MassXiMinus + newPx * newPx + newPy * newPy + particle.pz() * particle.pz());
595+
598596
xiDecayVertex.push_back(newX);
599597
xiDecayVertex.push_back(newY);
600-
xiDecayVertex.push_back(particle.vz() + xi_rxyz * (particle.pz() / particle.p()));
598+
xiDecayVertex.push_back(particle.vz() + rxyzXi * (particle.pz() / particle.p()));
601599

602-
std::vector<double> xiDaughters = {la_mass, pi_mass};
603-
TLorentzVector xi(newPx, newPy, particle.pz(), particle.e());
600+
std::vector<double> xiDaughters = {o2::constants::physics::MassLambda, o2::constants::physics::MassPionCharged};
601+
TLorentzVector xi(newPx, newPy, particle.pz(), newE);
604602
TGenPhaseSpace xiDecay;
605603
xiDecay.SetDecay(xi, 2, xiDaughters.data());
606604
xiDecay.Generate();
607605
decayDaughters.push_back(*xiDecay.GetDecay(1));
608606
TLorentzVector la = *xiDecay.GetDecay(0);
609607

610-
double la_gamma = 1 / std::sqrt(1 + (la.P() * la.P()) / (la_mass * la_mass));
611-
double la_ctau = 7.89 * la_gamma;
612-
std::vector<double> laDaughters = {pi_mass, pr_mass};
613-
double la_rxyz = (-la_ctau * std::log(1 - u));
614-
laDecayVertex.push_back(xiDecayVertex[0] + la_rxyz * (xiDecay.GetDecay(0)->Px() / xiDecay.GetDecay(0)->P()));
615-
laDecayVertex.push_back(xiDecayVertex[1] + la_rxyz * (xiDecay.GetDecay(0)->Py() / xiDecay.GetDecay(0)->P()));
616-
laDecayVertex.push_back(xiDecayVertex[2] + la_rxyz * (xiDecay.GetDecay(0)->Pz() / xiDecay.GetDecay(0)->P()));
608+
const double uLa = rand.Uniform(0, 1);
609+
const double ctauLa = 7.845; // cm
610+
const double betaGammaLa = la.P() / o2::constants::physics::MassLambda;
611+
const double rxyzLa = (-betaGammaLa * ctauLa * std::log(1 - uLa));
612+
laDecayVertex.push_back(xiDecayVertex[0] + rxyzLa * (xiDecay.GetDecay(0)->Px() / xiDecay.GetDecay(0)->P()));
613+
laDecayVertex.push_back(xiDecayVertex[1] + rxyzLa * (xiDecay.GetDecay(0)->Py() / xiDecay.GetDecay(0)->P()));
614+
laDecayVertex.push_back(xiDecayVertex[2] + rxyzLa * (xiDecay.GetDecay(0)->Pz() / xiDecay.GetDecay(0)->P()));
617615

616+
std::vector<double> laDaughters = {o2::constants::physics::MassPionCharged, o2::constants::physics::MassProton};
618617
TGenPhaseSpace laDecay;
619618
laDecay.SetDecay(la, 2, laDaughters.data());
620619
laDecay.Generate();
621620
decayDaughters.push_back(*laDecay.GetDecay(0));
622621
decayDaughters.push_back(*laDecay.GetDecay(1));
623622
}
623+
624624
/// Function to decay the V0
625625
/// \param particle the particle to decay
626626
/// \param decayDaughters the address of resulting daughters
@@ -629,36 +629,36 @@ struct OnTheFlyTracker {
629629
void decayV0Particle(McParticleType particle, std::vector<TLorentzVector>& decayDaughters, std::vector<double>& v0DecayVertex, int pdgCode)
630630
{
631631
double u = rand.Uniform(0, 1);
632-
double v0_mass = -1.;
633-
double negDau_mass = -1.;
634-
double posDau_mass = -1.;
632+
double v0Mass = -1.;
633+
double negDauMass = -1.;
634+
double posDauMass = -1.;
635635
double ctau = -1.;
636+
636637
if (std::abs(pdgCode) == kK0Short) {
637-
v0_mass = o2::constants::physics::MassK0Short;
638-
negDau_mass = o2::constants::physics::MassPionCharged;
639-
posDau_mass = o2::constants::physics::MassPionCharged;
638+
v0Mass = o2::constants::physics::MassK0Short;
639+
negDauMass = o2::constants::physics::MassPionCharged;
640+
posDauMass = o2::constants::physics::MassPionCharged;
640641
ctau = 2.68;
641642
} else if (std::abs(pdgCode) == kLambda0) {
642-
v0_mass = o2::constants::physics::MassLambda;
643-
negDau_mass = o2::constants::physics::MassPionCharged;
644-
posDau_mass = o2::constants::physics::MassProton;
643+
v0Mass = o2::constants::physics::MassLambda;
644+
negDauMass = o2::constants::physics::MassPionCharged;
645+
posDauMass = o2::constants::physics::MassProton;
645646
ctau = 7.845;
646647
} else if (std::abs(pdgCode) == kLambda0Bar) {
647-
v0_mass = o2::constants::physics::MassLambda;
648-
negDau_mass = o2::constants::physics::MassProton;
649-
posDau_mass = o2::constants::physics::MassPionCharged;
648+
v0Mass = o2::constants::physics::MassLambda;
649+
negDauMass = o2::constants::physics::MassProton;
650+
posDauMass = o2::constants::physics::MassPionCharged;
650651
ctau = 7.845;
651652
}
652653

653-
double v0_gamma = 1 / std::sqrt(1 + (particle.p() * particle.p()) / (v0_mass * v0_mass));
654-
double v0_ctau = ctau * v0_gamma;
655-
double v0_rxyz = (-v0_ctau * std::log(1 - u));
654+
const double v0BetaGamma = particle.p() / v0Mass;
655+
const double v0rxyz = (-v0BetaGamma * ctau * std::log(1 - u));
656656
TLorentzVector v0(particle.px(), particle.py(), particle.pz(), particle.e());
657657

658-
v0DecayVertex.push_back(particle.vx() + v0_rxyz * (particle.px() / particle.p()));
659-
v0DecayVertex.push_back(particle.vy() + v0_rxyz * (particle.py() / particle.p()));
660-
v0DecayVertex.push_back(particle.vz() + v0_rxyz * (particle.pz() / particle.p()));
661-
std::vector<double> v0Daughters = {negDau_mass, posDau_mass};
658+
v0DecayVertex.push_back(particle.vx() + v0rxyz * (particle.px() / particle.p()));
659+
v0DecayVertex.push_back(particle.vy() + v0rxyz * (particle.py() / particle.p()));
660+
v0DecayVertex.push_back(particle.vz() + v0rxyz * (particle.pz() / particle.p()));
661+
std::vector<double> v0Daughters = {negDauMass, posDauMass};
662662

663663
TGenPhaseSpace v0Decay;
664664
v0Decay.SetDecay(v0, 2, v0Daughters.data());
@@ -744,7 +744,7 @@ struct OnTheFlyTracker {
744744
if (mcParticle.pdgCode() == kXiMinus) {
745745
o2::track::TrackParCov xiTrackParCov;
746746
o2::upgrade::convertMCParticleToO2Track(mcParticle, xiTrackParCov, pdgDB);
747-
decayParticle(mcParticle, xiTrackParCov, decayProducts, xiDecayVertex, laDecayVertex);
747+
decayCascade(mcParticle, xiTrackParCov, decayProducts, xiDecayVertex, laDecayVertex);
748748
xiDecayRadius2D = std::hypot(xiDecayVertex[0], xiDecayVertex[1]);
749749
laDecayRadius2D = std::hypot(laDecayVertex[0], laDecayVertex[1]);
750750
}
@@ -1214,22 +1214,28 @@ struct OnTheFlyTracker {
12141214
std::array{negP[0], negP[1], negP[2]}},
12151215
std::array{o2::constants::physics::MassPionCharged,
12161216
o2::constants::physics::MassPionCharged});
1217-
} else
1217+
} else {
12181218
thisV0.mK0 = -1;
1219+
}
1220+
12191221
if (isLambda) {
12201222
thisV0.mLambda = RecoDecay::m(std::array{std::array{posP[0], posP[1], posP[2]},
12211223
std::array{negP[0], negP[1], negP[2]}},
12221224
std::array{o2::constants::physics::MassPionCharged,
12231225
o2::constants::physics::MassProton});
1224-
} else
1226+
} else {
12251227
thisV0.mLambda = -1;
1228+
}
1229+
12261230
if (isAntiLambda) {
12271231
thisV0.mAntiLambda = RecoDecay::m(std::array{std::array{posP[0], posP[1], posP[2]},
12281232
std::array{negP[0], negP[1], negP[2]}},
12291233
std::array{o2::constants::physics::MassProton,
12301234
o2::constants::physics::MassPionCharged});
1231-
} else
1235+
} else {
12321236
thisV0.mAntiLambda = -1;
1237+
}
1238+
12331239
if (v0DecaySettings.doV0QA) {
12341240
histos.fill(HIST("V0Building/hV0Building"), 4.0f);
12351241

0 commit comments

Comments
 (0)