algorithms/reco/ElectronReconstruction.cc

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2023 Daniel Brandenburg
0003 #include "ElectronReconstruction.h"
0004
0005 #include <edm4eic/ClusterCollection.h>
0006 #include <edm4eic/ReconstructedParticleCollection.h>
0007 #include <edm4hep/utils/vector_utils.h>
0008 #include <fmt/core.h>
0009 #include <podio/RelationRange.h>
0010
0011 #include "algorithms/reco/ElectronReconstructionConfig.h"
0012
0013 namespace eicrecon {
0014
0015 void ElectronReconstruction::init(std::shared_ptr<spdlog::logger> logger) { m_log = logger; }
0016
0017 std::unique_ptr<edm4eic::ReconstructedParticleCollection>
0018 ElectronReconstruction::execute(const edm4eic::ReconstructedParticleCollection* rcparts) {
0019
0020   // Some obvious improvements:
0021   // - E/p cut from real study optimized for electron finding and hadron rejection
0022   // - use of any HCAL info?
0023   // - check for duplicates?
0024
0025   // output container
0026   auto out_electrons = std::make_unique<edm4eic::ReconstructedParticleCollection>();
0027   out_electrons
0028       ->setSubsetCollection(); // out_electrons is a subset of the ReconstructedParticles collection
0029
0030   for (const auto particle : *rcparts) {
0031     // if we found a reco particle then test for electron compatibility
0032     if (particle.getClusters().empty()) {
0033       continue;
0034     }
0035     if (particle.getCharge() == 0) { // Skip over photons/other particles without a track
0036       continue;
0037     }
0038     double E      = particle.getClusters()[0].getEnergy();
0039     double p      = edm4hep::utils::magnitude(particle.getMomentum());
0040     double EOverP = E / p;
0041
0042     m_log->trace(
0043         "ReconstructedElectron: Energy={} GeV, p={} GeV, E/p = {} for PDG (from truth): {}", E, p,
0044         EOverP, particle.getPDG());
0045     // Apply the E/p cut here to select electons
0046     if (EOverP >= m_cfg.min_energy_over_momentum && EOverP <= m_cfg.max_energy_over_momentum) {
0047       out_electrons->push_back(particle);
0048     }
0049   }
0050   m_log->debug("Found {} electron candidates", out_electrons->size());
0051   return out_electrons;
0052 }
0053 } // namespace eicrecon