43 ROOT::Math::PxPyPzEVector& KMomentum,
44 const std::vector<Particle*> daughters,
50 ROOT::Math::PxPyPzEVector klDaughters;
51 ROOT::Math::PxPyPzEVector pDaughters;
53 for (
auto daughter : daughters) {
54 if (daughter->getPDGCode() ==
Const::Klong.getPDGCode()) {
56 B2FATAL(
"More than one K_L is detected! This tool accepts only one K_L in the final state.");
58 const KLMCluster* klm_cluster = daughter->getKLMCluster();
60 klDaughters = klm_cluster->getMomentum();
62 const ECLCluster* ecl_cluster = daughter->getECLCluster();
67 B2ERROR(
"K_L candidate must have related KLM or ECL cluster!");
72 pDaughters += daughter->get4Vector();
75 m_j = daughter->getPDGMass();
76 idx = daughter->getArrayIndex() + idx * 100;
81 B2FATAL(
"This tool is meant to reconstruct decays with a K_L0 in the final state. There is no K_L0 in this decay!");
84 if (daughters.size() == 3) {
90 const double s_p = (klDaughters.Vect().
Unit()).Dot(pDaughters.Vect());
91 const double m_sum = (m_b * m_b) - (m_j * m_j) - m_k2;
93 const double s_p2 = s_p * s_p;
94 const double m_sum2 = m_sum * m_sum;
95 const double s_pm = s_p * m_sum;
96 const double e_j2 = pDaughters.E() * pDaughters.E();
98 const double k_mag1 = (s_pm + std::sqrt(s_p2 * m_sum2 - 4 * (e_j2 - s_p2) * (e_j2 * m_k2 - m_sum2 / 4))) / (2 *
100 const double k_mag2 = (s_pm - std::sqrt(s_p2 * m_sum2 - 4 * (e_j2 - s_p2) * (e_j2 * m_k2 - m_sum2 / 4))) / (2 *
104 ROOT::Math::PxPyPzEVector missDaughters;
107 missDaughters = k_mag1 * klDaughters / klDaughters.P();
109 missDaughters = k_mag2 * klDaughters / klDaughters.P();
111 missDaughters.SetE(std::sqrt(m_k2 + missDaughters.P2()));
113 if (std::isnan(missDaughters.P()))
116 ROOT::Math::PxPyPzEVector mom = pDaughters + missDaughters;
117 mom.SetE(std::sqrt(m_b * m_b + mom.P2()));
118 if (std::isnan(mom.P()))
121 KMomentum = missDaughters;