MCEventsAlg.cxx

#include "MCEventsAlg.h"
#include "TrackerIdentifier/FaserSCT_ID.h"


MCEventsAlg::MCEventsAlg(const std::string& name, ISvcLocator* pSvcLocator)
    : AthReentrantAlgorithm(name, pSvcLocator) {}


StatusCode MCEventsAlg::initialize() {
  ATH_CHECK(detStore()->retrieve(m_idHelper,"FaserSCT_ID"));
  ATH_CHECK(m_mcEventKey.initialize());
  ATH_CHECK(m_faserSiHitKey.initialize());

  return StatusCode::SUCCESS;
}


StatusCode MCEventsAlg::execute(const EventContext& ctx) const {
  SG::ReadHandle<FaserSiHitCollection> siHitCollection (m_faserSiHitKey, ctx);
  std::map<int, std::vector<const FaserSiHit*>> particle_map {};
  ATH_CHECK(siHitCollection.isValid());
  for (const FaserSiHit& hit : *siHitCollection) {
    int barcode = hit.particleLink().barcode();
    particle_map[barcode].push_back(&hit);
  }

  SG::ReadHandle<McEventCollection> mcEvents (m_mcEventKey, ctx);
  ATH_CHECK(mcEvents.isValid());
  for (const HepMC::GenEvent* evt : *mcEvents) {
    for (const HepMC::GenParticle* p : evt->particle_range()) {
      const HepMC::FourVector vertex = p->production_vertex()->position();
      const HepMC::FourVector momentum = p->momentum();
      ATH_MSG_DEBUG("pdg: " << p->pdg_id() << ", barcode: " << p->barcode());
      ATH_MSG_DEBUG("vertex: " << vertex.x() << ", " << vertex.y() << ", " << vertex.z());
      ATH_MSG_DEBUG("momentum: " << momentum.px() << ", " << momentum.py() << ", " << momentum.pz());
      auto hits = particle_map[p->barcode()];
      for (const FaserSiHit* hit : hits) {
        int station = hit->getStation();
        int layer = hit->getPlane();
        int phi = hit->getRow();
        int eta = hit->getModule();
        int side = hit->getSensor();
        ATH_MSG_DEBUG(station << "/" << layer << "/" << phi << "/" << eta << "/" << side);
      }
    }
  }
  return StatusCode::SUCCESS;
}


StatusCode MCEventsAlg::finalize() {
  return StatusCode::SUCCESS;
}