File indexing completed on 2024-09-27 07:02:31
0001
0002
0003
0004 #include <algorithms/truth/MC2SmearedParticle.h>
0005
0006 #include <cmath>
0007 #include <edm4hep/utils/vector_utils.h>
0008
0009 namespace algorithms::truth {
0010
0011 void MC2SmearedParticle::init() {
0012 ;
0013 }
0014
0015 void MC2SmearedParticle::process(const MC2SmearedParticle::Input& input,
0016 const MC2SmearedParticle::Output& output) const {
0017 const auto [parts] = input;
0018 auto [out_parts] = output;
0019
0020 for (const auto& p : *parts) {
0021 if (p.getGeneratorStatus() > 1) {
0022 if (aboveDebugThreshold()) {
0023 debug() << "ignoring particle with generatorStatus = " << p.getGeneratorStatus() << endmsg;
0024 }
0025 continue;
0026 }
0027
0028
0029
0030
0031 const auto pvec = p.getMomentum();
0032 const auto pgen = std::hypot(pvec.x, pvec.y, pvec.z);
0033 const auto momentum = pgen * m_rng.gaussian<double>(0., m_smearing);
0034
0035 using MomType = decltype(edm4eic::ReconstructedParticle().getMomentum().x);
0036 const MomType energy =
0037 std::sqrt(p.getEnergy() * p.getEnergy() - pgen * pgen + momentum * momentum);
0038 const MomType px = p.getMomentum().x * momentum / pgen;
0039 const MomType py = p.getMomentum().y * momentum / pgen;
0040 const MomType pz = p.getMomentum().z * momentum / pgen;
0041
0042 const MomType dpx = m_smearing * px;
0043 const MomType dpy = m_smearing * py;
0044 const MomType dpz = m_smearing * pz;
0045 const MomType dE = m_smearing * energy;
0046
0047
0048 const MomType vx = p.getVertex().x;
0049 const MomType vy = p.getVertex().y;
0050 const MomType vz = p.getVertex().z;
0051
0052 auto rec_part = out_parts->create();
0053 rec_part.setType(-1);
0054 rec_part.setEnergy(energy);
0055 rec_part.setMomentum({px, py, pz});
0056 rec_part.setReferencePoint({vx, vy, vz});
0057 rec_part.setCharge(p.getCharge());
0058 rec_part.setMass(p.getMass());
0059 rec_part.setGoodnessOfPID(1);
0060 rec_part.setCovMatrix({dpx, dpy, dpz, dE});
0061 rec_part.setPDG(p.getPDG());
0062 }
0063 }
0064
0065 }
0066