Warning, file /juggler/JugDigi/src/components/SiliconTrackerDigi.cpp was not indexed
or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).
0001
0002
0003
0004 #include <algorithm>
0005 #include <cmath>
0006
0007 #include "Gaudi/Property.h"
0008 #include "Gaudi/Algorithm.h"
0009 #include "GaudiKernel/PhysicalConstants.h"
0010 #include "GaudiKernel/RndmGenerators.h"
0011
0012 #include <k4FWCore/DataHandle.h>
0013
0014
0015
0016 #include "edm4hep/EDM4hepVersion.h"
0017 #include "edm4hep/MCParticleCollection.h"
0018 #include "edm4hep/SimTrackerHitCollection.h"
0019
0020 #include "edm4eic/RawTrackerHitCollection.h"
0021
0022 namespace Jug::Digi {
0023
0024
0025
0026
0027
0028 class SiliconTrackerDigi : public Gaudi::Algorithm {
0029 private:
0030 Gaudi::Property<double> m_timeResolution{this, "timeResolution", 10};
0031 Gaudi::Property<double> m_threshold{this, "threshold", 0. * Gaudi::Units::keV};
0032 Rndm::Numbers m_gaussDist;
0033 mutable DataHandle<edm4hep::SimTrackerHitCollection> m_inputHitCollection{"inputHitCollection", Gaudi::DataHandle::Reader,
0034 this};
0035 mutable DataHandle<edm4eic::RawTrackerHitCollection> m_outputHitCollection{"outputHitCollection", Gaudi::DataHandle::Writer,
0036 this};
0037
0038 public:
0039
0040 SiliconTrackerDigi(const std::string& name, ISvcLocator* svcLoc) : Gaudi::Algorithm(name, svcLoc) {
0041 declareProperty("inputHitCollection", m_inputHitCollection, "");
0042 declareProperty("outputHitCollection", m_outputHitCollection, "");
0043 }
0044 StatusCode initialize() override {
0045 if (Gaudi::Algorithm::initialize().isFailure()) {
0046 return StatusCode::FAILURE;
0047 }
0048 IRndmGenSvc* randSvc = Gaudi::svcLocator()->service<IRndmGenSvc>("RndmGenSvc", true);
0049 StatusCode sc = m_gaussDist.initialize(randSvc, Rndm::Gauss(0.0, m_timeResolution.value()));
0050 if (!sc.isSuccess()) {
0051 return StatusCode::FAILURE;
0052 }
0053 return StatusCode::SUCCESS;
0054 }
0055 StatusCode execute(const EventContext&) const override {
0056
0057 const auto* const simhits = m_inputHitCollection.get();
0058
0059 auto* rawhits = m_outputHitCollection.createAndPut();
0060
0061 std::map<long long, int> cell_hit_map;
0062 for (const auto& ahit : *simhits) {
0063 if (msgLevel(MSG::DEBUG)) {
0064 debug() << "--------------------" << ahit.getCellID() << endmsg;
0065 debug() << "Hit in cellID = " << ahit.getCellID() << endmsg;
0066 debug() << " position = (" << ahit.getPosition().x << "," << ahit.getPosition().y << ","
0067 << ahit.getPosition().z << ")" << endmsg;
0068 debug() << " xy_radius = " << std::hypot(ahit.getPosition().x, ahit.getPosition().y) << endmsg;
0069 debug() << " momentum = (" << ahit.getMomentum().x << "," << ahit.getMomentum().y << ","
0070 << ahit.getMomentum().z << ")" << endmsg;
0071 }
0072 if (ahit.getEDep() * Gaudi::Units::keV < m_threshold) {
0073 if (msgLevel(MSG::DEBUG)) {
0074 debug() << " edep = " << ahit.getEDep() << " (below threshold of " << m_threshold / Gaudi::Units::keV
0075 << " keV)" << endmsg;
0076 }
0077 continue;
0078 } else {
0079 if (msgLevel(MSG::DEBUG)) {
0080 debug() << " edep = " << ahit.getEDep() << endmsg;
0081 }
0082 }
0083 if (cell_hit_map.count(ahit.getCellID()) == 0) {
0084 cell_hit_map[ahit.getCellID()] = rawhits->size();
0085 rawhits->create(
0086 ahit.getCellID(),
0087 #if EDM4HEP_BUILD_VERSION >= EDM4HEP_VERSION(0, 99, 0)
0088 ahit.getParticle().getTime() * 1e6 + m_gaussDist() * 1e3,
0089 #else
0090 ahit.getMCParticle().getTime() * 1e6 + m_gaussDist() * 1e3,
0091 #endif
0092 std::llround(ahit.getEDep() * 1e6)
0093 );
0094 } else {
0095 auto hit = (*rawhits)[cell_hit_map[ahit.getCellID()]];
0096 #if EDM4HEP_BUILD_VERSION >= EDM4HEP_VERSION(0, 99, 0)
0097 hit.setTimeStamp(ahit.getParticle().getTime() * 1e6 + m_gaussDist() * 1e3);
0098 #else
0099 hit.setTimeStamp(ahit.getMCParticle().getTime() * 1e6 + m_gaussDist() * 1e3);
0100 #endif
0101 auto ch = hit.getCharge();
0102 hit.setCharge(ch + std::llround(ahit.getEDep() * 1e6));
0103 }
0104 }
0105 return StatusCode::SUCCESS;
0106 }
0107 };
0108
0109 DECLARE_COMPONENT(SiliconTrackerDigi)
0110
0111 }