diff --git a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp index 937f9e1f7168204c8483a3d01d9596eec381666e..b31b62d8332777e396ce835ee96e814476465b65 100644 --- a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp +++ b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp @@ -1,5 +1,6 @@ #include "Edm4hepWriterAnaElemTool.h" +#include "G4Step.hh" #include "G4Event.hh" #include "G4THitsCollection.hh" #include "G4EventManager.hh" @@ -19,6 +20,23 @@ DECLARE_COMPONENT(Edm4hepWriterAnaElemTool) void Edm4hepWriterAnaElemTool::BeginOfRunAction(const G4Run*) { G4cout << "Begin Run of detector simultion..." << G4endl; + + // access geometry service + m_geosvc = service<IGeomSvc>("GeomSvc"); + if (m_geosvc) { + dd4hep::Detector* dd4hep_geo = m_geosvc->lcdd(); + // try to get the constants + R = dd4hep_geo->constant<double>("tracker_region_rmax")/dd4hep::mm*CLHEP::mm; + Z = dd4hep_geo->constant<double>("tracker_region_zmax")/dd4hep::mm*CLHEP::mm; + + info() << "Tracker Region " + << " R: " << R + << " Z: " << Z + << endmsg; + + } else { + error() << "Failed to find GeomSvc." << endmsg; + } } void @@ -498,8 +516,36 @@ Edm4hepWriterAnaElemTool::PostUserTrackingAction(const G4Track* track) { } void -Edm4hepWriterAnaElemTool::UserSteppingAction(const G4Step*) { - +Edm4hepWriterAnaElemTool::UserSteppingAction(const G4Step* aStep) { + auto aTrack = aStep->GetTrack(); + // try to get user track info + auto trackinfo = dynamic_cast<CommonUserTrackInfo*>(aTrack->GetUserInformation()); + + // ======================================================================== + // Note: + // if there is no track info, then do nothing. + // ======================================================================== + if (trackinfo) { + // back scattering is defined as following: + // - pre point is not in tracker + // - post point is in tracker + // For details, look at Mokka's source code: + // https://llrforge.in2p3.fr/trac/Mokka/browser/trunk/source/Kernel/src/SteppingAction.cc + const auto& pre_pos = aStep->GetPreStepPoint()->GetPosition(); + const auto& post_pos = aStep->GetPostStepPoint()->GetPosition(); + + bool is_pre_in_tracker = pre_pos.perp() < R && std::fabs(pre_pos.z()) < Z; + bool is_post_in_tracker = post_pos.perp() < R && std::fabs(post_pos.z()) < Z; + + if ((!is_pre_in_tracker) and is_post_in_tracker) { + // set back scattering + auto idxedm4hep = trackinfo->idxEdm4hep(); + mcCol->at(idxedm4hep).setBackscatter(true); + info() << " set back scatter for MCParticle " + << " (ID: " << idxedm4hep << ")" + << endmsg; + } + } } StatusCode diff --git a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h index aad11e8d1f0963324abfd9c9039467673fdb67e9..9782cf5f5bd0fbbcca105f9d24b2eb25387ac915 100644 --- a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h +++ b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h @@ -9,6 +9,8 @@ #include "DetSimInterface/CommonUserEventInfo.hh" #include "DetSimInterface/CommonUserTrackInfo.hh" +#include "DetInterface/IGeomSvc.h" + #include "edm4hep/MCParticleCollection.h" #include "edm4hep/SimTrackerHitCollection.h" #include "edm4hep/SimCalorimeterHitCollection.h" @@ -143,6 +145,12 @@ private: std::map<int, int> m_track2primary; CommonUserEventInfo* m_userinfo = nullptr; + + // get the limitation of R/Z in tracker + SmartIF<IGeomSvc> m_geosvc; + double R = 0; + double Z = 0; + }; #endif