diff --git a/Analysis/DumpEvent/src/DumpSimHitAlg.cpp b/Analysis/DumpEvent/src/DumpSimHitAlg.cpp
index 5c9a3d1f22f99b75f7d4e987e4befa9e5066e30e..9ea53b5fb0e6f0596a3215a4433f3f623001b11d 100644
--- a/Analysis/DumpEvent/src/DumpSimHitAlg.cpp
+++ b/Analysis/DumpEvent/src/DumpSimHitAlg.cpp
@@ -51,6 +51,13 @@ StatusCode DumpSimHitAlg::initialize() {
StatusCode DumpSimHitAlg::execute() {
auto mcCol = m_mcParCol.get();
+ for (auto particle: *mcCol) {
+ info() << "mc particle -> "
+ << " (ID: " << particle.getObjectID().index << ") "
+ << " (simulator status: " << particle.getSimulatorStatus() << ") "
+ << endmsg;
+ }
+
auto vxdCol = m_VXDCol.get();
for (auto hit: *vxdCol) {
diff --git a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.cpp
index ca3bf7b63705b2b7ecf841a30d9e2b089c31cfff..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
@@ -472,6 +490,15 @@ Edm4hepWriterAnaElemTool::PostUserTrackingAction(const G4Track* track) {
mcp.addToParents(primary_particle);
primary_particle.addToDaughters(mcp);
+
+ // store the edm4hep obj idx in track info.
+ // using this idx, the MCParticle object could be modified later.
+ auto trackinfo = new CommonUserTrackInfo();
+ trackinfo->setIdxEdm4hep(mcp.getObjectID().index);
+ sectrk->SetUserInformation(trackinfo);
+ info() << " Appending MCParticle: (id: "
+ << mcp.getObjectID().index << ")"
+ << endmsg;
}
}
}
@@ -489,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 b290c40bd8f10545ce16d1c8588eb12f5aa54998..9782cf5f5bd0fbbcca105f9d24b2eb25387ac915 100644
--- a/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h
+++ b/Simulation/DetSimAna/src/Edm4hepWriterAnaElemTool.h
@@ -7,6 +7,9 @@
#include "k4FWCore/DataHandle.h"
#include "DetSimInterface/IAnaElemTool.h"
#include "DetSimInterface/CommonUserEventInfo.hh"
+#include "DetSimInterface/CommonUserTrackInfo.hh"
+
+#include "DetInterface/IGeomSvc.h"
#include "edm4hep/MCParticleCollection.h"
#include "edm4hep/SimTrackerHitCollection.h"
@@ -142,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
diff --git a/Simulation/DetSimInterface/CMakeLists.txt b/Simulation/DetSimInterface/CMakeLists.txt
index b5fd14616a4bf34b98f3d4545d6c398954423bda..21963df1a8e0fd3289849181995d25908f69641e 100644
--- a/Simulation/DetSimInterface/CMakeLists.txt
+++ b/Simulation/DetSimInterface/CMakeLists.txt
@@ -5,6 +5,7 @@
gaudi_add_library(DetSimInterface
SOURCES src/IDetSimSvc.cpp
src/CommonUserEventInfo.cc
+ src/CommonUserTrackInfo.cc
LINK Gaudi::GaudiKernel
${Geant4_LIBRARIES}
)
diff --git a/Simulation/DetSimInterface/include/DetSimInterface/CommonUserTrackInfo.hh b/Simulation/DetSimInterface/include/DetSimInterface/CommonUserTrackInfo.hh
new file mode 100644
index 0000000000000000000000000000000000000000..ce63e602d26267dd87465406d454c7d558eb83ad
--- /dev/null
+++ b/Simulation/DetSimInterface/include/DetSimInterface/CommonUserTrackInfo.hh
@@ -0,0 +1,35 @@
+#ifndef CommonUserTrackInfo_hh
+#define CommonUserTrackInfo_hh
+
+/* Description:
+ * This class is a part of simulation framework to extend the G4Track.
+ *
+ * Some secondaries are created due to decay. However, their G4 Track IDs are
+ * not valid until the tracks are tracking by geant4. In order to associate
+ * these tracks and their edm4hep MC particle, we use the track information
+ * to record the extra track information.
+ *
+ * Author:
+ * Tao Lin <lintao AT ihep.ac.cn>
+ */
+
+#include "G4VUserTrackInformation.hh"
+
+class CommonUserTrackInfo: public G4VUserTrackInformation {
+public:
+ CommonUserTrackInfo();
+ ~CommonUserTrackInfo();
+
+public:
+
+ virtual void Print() const;
+
+ // get the idx in the EDM4hep MC particle collection
+ bool setIdxEdm4hep(int idxEdm4hep);
+ int idxEdm4hep() const;
+
+private:
+ int m_idxEdm4hep = -1;
+};
+
+#endif
diff --git a/Simulation/DetSimInterface/src/CommonUserTrackInfo.cc b/Simulation/DetSimInterface/src/CommonUserTrackInfo.cc
new file mode 100644
index 0000000000000000000000000000000000000000..9a923a35a22e7fafe496c827e066e591864f3c50
--- /dev/null
+++ b/Simulation/DetSimInterface/src/CommonUserTrackInfo.cc
@@ -0,0 +1,22 @@
+#include "DetSimInterface/CommonUserTrackInfo.hh"
+#include <iostream>
+
+CommonUserTrackInfo::CommonUserTrackInfo() {
+
+}
+
+CommonUserTrackInfo::~CommonUserTrackInfo() {
+
+}
+
+void CommonUserTrackInfo::Print() const {
+
+}
+
+bool CommonUserTrackInfo::setIdxEdm4hep(int idxEdm4hep) {
+ m_idxEdm4hep = idxEdm4hep;
+}
+
+int CommonUserTrackInfo::idxEdm4hep() const {
+ return m_idxEdm4hep;
+}