From 8551d89a75b0f8d270469939f5b3a31dcb88f2b4 Mon Sep 17 00:00:00 2001
From: Chengdong Fu <fucd@ihep.ac.cn>
Date: Sun, 27 Mar 2022 21:46:14 +0800
Subject: [PATCH] add GenericTrackerSD

---
 Detector/DetCRD/scripts/CRD-Sim.py            | 17 ++++++
 Detector/DetCRD/scripts/CRD_o1_v01-SimRec.py  | 17 ++++++
 Detector/DetCRD/scripts/CRD_o1_v02-SimRec.py  | 17 ++++++
 .../DetSimGeom/src/AnExampleDetElemTool.cpp   | 11 +++-
 .../DetSimGeom/src/AnExampleDetElemTool.h     |  1 +
 Simulation/DetSimSD/CMakeLists.txt            |  2 +
 .../src/GenericTrackerSensDetTool.cpp         | 40 +++++++++++++
 .../DetSimSD/src/GenericTrackerSensDetTool.h  | 35 +++++++++++
 .../src/GenericTrackerSensitiveDetector.cpp   | 59 +++++++++++++++++++
 .../src/GenericTrackerSensitiveDetector.h     | 24 ++++++++
 10 files changed, 222 insertions(+), 1 deletion(-)
 create mode 100644 Simulation/DetSimSD/src/GenericTrackerSensDetTool.cpp
 create mode 100644 Simulation/DetSimSD/src/GenericTrackerSensDetTool.h
 create mode 100644 Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.cpp
 create mode 100644 Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.h

diff --git a/Detector/DetCRD/scripts/CRD-Sim.py b/Detector/DetCRD/scripts/CRD-Sim.py
index a4ad74c4..0e734b78 100644
--- a/Detector/DetCRD/scripts/CRD-Sim.py
+++ b/Detector/DetCRD/scripts/CRD-Sim.py
@@ -89,6 +89,23 @@ detsimalg.AnaElems = [
 ]
 detsimalg.RootDetElem = "WorldDetElemTool"
 
+dedxoption = "BetheBlochEquationDedxSimTool"
+from Configurables import DriftChamberSensDetTool
+dc_sensdettool = DriftChamberSensDetTool("DriftChamberSensDetTool")
+dc_sensdettool.DedxSimTool = dedxoption
+
+from Configurables import DummyDedxSimTool
+from Configurables import BetheBlochEquationDedxSimTool
+
+if dedxoption == "DummyDedxSimTool":
+    dedx_simtool = DummyDedxSimTool("DummyDedxSimTool")
+elif dedxoption == "BetheBlochEquationDedxSimTool":
+    dedx_simtool = BetheBlochEquationDedxSimTool("BetheBlochEquationDedxSimTool")
+    dedx_simtool.material_Z = 2
+    dedx_simtool.material_A = 4
+    dedx_simtool.scale = 10
+    dedx_simtool.resolution = 0.0001
+
 # output
 from Configurables import PodioOutput
 out = PodioOutput("outputalg")
diff --git a/Detector/DetCRD/scripts/CRD_o1_v01-SimRec.py b/Detector/DetCRD/scripts/CRD_o1_v01-SimRec.py
index 1810a574..19836013 100644
--- a/Detector/DetCRD/scripts/CRD_o1_v01-SimRec.py
+++ b/Detector/DetCRD/scripts/CRD_o1_v01-SimRec.py
@@ -84,6 +84,23 @@ detsimalg.AnaElems = [
 ]
 detsimalg.RootDetElem = "WorldDetElemTool"
 
+dedxoption = "BetheBlochEquationDedxSimTool"
+from Configurables import DriftChamberSensDetTool
+dc_sensdettool = DriftChamberSensDetTool("DriftChamberSensDetTool")
+dc_sensdettool.DedxSimTool = dedxoption
+
+from Configurables import DummyDedxSimTool
+from Configurables import BetheBlochEquationDedxSimTool
+
+if dedxoption == "DummyDedxSimTool":
+    dedx_simtool = DummyDedxSimTool("DummyDedxSimTool")
+elif dedxoption == "BetheBlochEquationDedxSimTool":
+    dedx_simtool = BetheBlochEquationDedxSimTool("BetheBlochEquationDedxSimTool")
+    dedx_simtool.material_Z = 2
+    dedx_simtool.material_A = 4
+    dedx_simtool.scale = 10
+    dedx_simtool.resolution = 0.0001
+    
 from Configurables import MarlinEvtSeeder
 evtseeder = MarlinEvtSeeder("EventSeeder")
 
diff --git a/Detector/DetCRD/scripts/CRD_o1_v02-SimRec.py b/Detector/DetCRD/scripts/CRD_o1_v02-SimRec.py
index acee999d..253ea21a 100644
--- a/Detector/DetCRD/scripts/CRD_o1_v02-SimRec.py
+++ b/Detector/DetCRD/scripts/CRD_o1_v02-SimRec.py
@@ -84,6 +84,23 @@ detsimalg.AnaElems = [
 ]
 detsimalg.RootDetElem = "WorldDetElemTool"
 
+dedxoption = "BetheBlochEquationDedxSimTool"
+from Configurables import DriftChamberSensDetTool
+dc_sensdettool = DriftChamberSensDetTool("DriftChamberSensDetTool")
+dc_sensdettool.DedxSimTool = dedxoption
+
+from Configurables import DummyDedxSimTool
+from Configurables import BetheBlochEquationDedxSimTool
+
+if dedxoption == "DummyDedxSimTool":
+    dedx_simtool = DummyDedxSimTool("DummyDedxSimTool")
+elif dedxoption == "BetheBlochEquationDedxSimTool":
+    dedx_simtool = BetheBlochEquationDedxSimTool("BetheBlochEquationDedxSimTool")
+    dedx_simtool.material_Z = 2
+    dedx_simtool.material_A = 4
+    dedx_simtool.scale = 10
+    dedx_simtool.resolution = 0.0001
+
 from Configurables import MarlinEvtSeeder
 evtseeder = MarlinEvtSeeder("EventSeeder")
 
diff --git a/Simulation/DetSimGeom/src/AnExampleDetElemTool.cpp b/Simulation/DetSimGeom/src/AnExampleDetElemTool.cpp
index 001675db..ff20af31 100644
--- a/Simulation/DetSimGeom/src/AnExampleDetElemTool.cpp
+++ b/Simulation/DetSimGeom/src/AnExampleDetElemTool.cpp
@@ -130,7 +130,16 @@ AnExampleDetElemTool::ConstructSDandField() {
 		    warning() << "TimeProjectionChamberSensDetTool is not found, and default tracker SD will be used" << endmsg;
 		  }
 		}
-
+		else {
+		  m_tracker_sdtool = ToolHandle<ISensDetTool>("GenericTrackerSensDetTool");
+		  if (m_tracker_sdtool) {
+		    info() << "Find the GenericTrackerSensDetTool" << endmsg;
+		    g4sd = m_tracker_sdtool->createSD(nam);
+		  }
+		  else {
+		    warning() << "GenericTrackerSensDetTool is not found. " << endmsg;
+		  }
+		}
             }
         }
         
diff --git a/Simulation/DetSimGeom/src/AnExampleDetElemTool.h b/Simulation/DetSimGeom/src/AnExampleDetElemTool.h
index 1ad3702f..0fd9c356 100644
--- a/Simulation/DetSimGeom/src/AnExampleDetElemTool.h
+++ b/Simulation/DetSimGeom/src/AnExampleDetElemTool.h
@@ -35,6 +35,7 @@ private:
     ToolHandle<ISensDetTool> m_calo_sdtool;
     ToolHandle<ISensDetTool> m_driftchamber_sdtool;
     ToolHandle<ISensDetTool> m_tpc_sdtool;
+    ToolHandle<ISensDetTool> m_tracker_sdtool;
 };
 
 #endif
diff --git a/Simulation/DetSimSD/CMakeLists.txt b/Simulation/DetSimSD/CMakeLists.txt
index eed28a78..fa8ba9ef 100644
--- a/Simulation/DetSimSD/CMakeLists.txt
+++ b/Simulation/DetSimSD/CMakeLists.txt
@@ -13,6 +13,8 @@ gaudi_add_module(DetSimSD
                          src/TimeProjectionChamberSensDetTool.cpp
                          src/TimeProjectionChamberSensitiveDetector.cpp
 
+                         src/GenericTrackerSensDetTool.cpp
+                         src/GenericTrackerSensitiveDetector.cpp
                  LINK DetSimInterface
                       DetInterface
                       ${DD4hep_COMPONENT_LIBRARIES} 
diff --git a/Simulation/DetSimSD/src/GenericTrackerSensDetTool.cpp b/Simulation/DetSimSD/src/GenericTrackerSensDetTool.cpp
new file mode 100644
index 00000000..480f509a
--- /dev/null
+++ b/Simulation/DetSimSD/src/GenericTrackerSensDetTool.cpp
@@ -0,0 +1,40 @@
+#include "GenericTrackerSensDetTool.h"
+
+#include "G4VSensitiveDetector.hh"
+
+#include "DD4hep/Detector.h"
+
+#include "GenericTrackerSensitiveDetector.h"
+
+#include "CLHEP/Units/SystemOfUnits.h"
+
+DECLARE_COMPONENT(GenericTrackerSensDetTool)
+
+StatusCode GenericTrackerSensDetTool::initialize() {
+  StatusCode sc;
+  debug() << "initialize() " << endmsg;
+  
+  m_geosvc = service<IGeomSvc>("GeomSvc");
+  if (!m_geosvc) {
+    error() << "Failed to find GeomSvc." << endmsg;
+    return StatusCode::FAILURE;
+  }
+  
+  return AlgTool::initialize();
+}
+
+StatusCode GenericTrackerSensDetTool::finalize() {
+  StatusCode sc;
+  
+  return sc;
+}
+
+G4VSensitiveDetector* GenericTrackerSensDetTool::createSD(const std::string& name) {
+  debug() << "createSD for " << name << endmsg;
+  
+  dd4hep::Detector* dd4hep_geo = m_geosvc->lcdd();
+
+  G4VSensitiveDetector* sd = new GenericTrackerSensitiveDetector(name, *dd4hep_geo);
+
+  return sd;
+}
diff --git a/Simulation/DetSimSD/src/GenericTrackerSensDetTool.h b/Simulation/DetSimSD/src/GenericTrackerSensDetTool.h
new file mode 100644
index 00000000..9f46d381
--- /dev/null
+++ b/Simulation/DetSimSD/src/GenericTrackerSensDetTool.h
@@ -0,0 +1,35 @@
+#ifndef GenericTrackerSensDetTool_h
+#define GenericTrackerSensDetTool_h
+
+/*
+ * GenericTrackerSensDetTool is used to create Time Projection Chamber SD.
+ */
+
+#include "GaudiKernel/AlgTool.h"
+#include "GaudiKernel/ToolHandle.h"
+#include "DetSimInterface/ISensDetTool.h"
+#include "DetInterface/IGeomSvc.h"
+
+#include "DD4hep/DD4hepUnits.h"
+
+class GenericTrackerSensDetTool: public extends<AlgTool, ISensDetTool> {
+  
+ public:
+
+  using extends::extends;
+
+  /// Overriding initialize and finalize
+  StatusCode initialize() override;
+  StatusCode finalize() override;
+
+  /// Override ISensDetTool
+  virtual G4VSensitiveDetector* createSD(const std::string& name) override;
+
+private:
+
+  // in order to initialize SD, we need to get the lcdd()
+  SmartIF<IGeomSvc> m_geosvc;
+
+};
+
+#endif
diff --git a/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.cpp b/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.cpp
new file mode 100644
index 00000000..d250079e
--- /dev/null
+++ b/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.cpp
@@ -0,0 +1,59 @@
+#include "GenericTrackerSensitiveDetector.h"
+
+#include "G4Step.hh"
+#include "G4VProcess.hh"
+//#include "G4SteppingManager.hh"
+#include "G4SDManager.hh"
+//#include "UserTrackInformation.hh"
+#include "DD4hep/DD4hepUnits.h"
+
+GenericTrackerSensitiveDetector::GenericTrackerSensitiveDetector(const std::string& name,
+								 dd4hep::Detector& description)
+  : DDG4SensitiveDetector(name, description),
+    m_hc(nullptr){
+  
+  G4String CollName=name+"Collection";
+  collectionName.insert(CollName);
+}
+
+void GenericTrackerSensitiveDetector::Initialize(G4HCofThisEvent* HCE){
+  m_hc = new HitCollection(GetName(), collectionName[0]);
+  int HCID = G4SDManager::GetSDMpointer()->GetCollectionID(m_hc);
+  HCE->AddHitsCollection( HCID, m_hc ); 
+}
+
+G4bool GenericTrackerSensitiveDetector::ProcessHits(G4Step* step, G4TouchableHistory*){
+  
+  G4TouchableHandle touchPost = step->GetPostStepPoint()->GetTouchableHandle(); 
+  G4TouchableHandle touchPre  = step->GetPreStepPoint()->GetTouchableHandle(); 
+  dd4hep::sim::Geant4StepHandler h(step);
+  dd4hep::Position prePos    = h.prePos();
+  dd4hep::Position postPos   = h.postPos();
+  dd4hep::Position direction = postPos - prePos;
+  dd4hep::Position position  = mean_direction(prePos,postPos);
+  double   hit_len   = direction.R();
+  if (hit_len > 0) {
+    double new_len = mean_length(h.preMom(),h.postMom())/hit_len;
+    direction *= new_len/hit_len;
+  }
+  dd4hep::sim::Geant4TrackerHit* hit = nullptr;
+  hit = new dd4hep::sim::Geant4TrackerHit(h.track->GetTrackID(),
+					  h.track->GetDefinition()->GetPDGEncoding(),
+					  step->GetTotalEnergyDeposit(),
+					  h.track->GetGlobalTime());
+
+  if ( hit )  {
+    hit->cellID  = getCellID( step ) ;
+    hit->energyDeposit =  step->GetTotalEnergyDeposit();
+    hit->position = position;
+    hit->momentum = direction;
+    hit->length   = hit_len;
+    m_hc->insert(hit);
+    return true;
+  }
+  throw std::runtime_error("new() failed: Cannot allocate hit object");
+  return false;
+}
+
+void GenericTrackerSensitiveDetector::EndOfEvent(G4HCofThisEvent* HCE){
+}
diff --git a/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.h b/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.h
new file mode 100644
index 00000000..63257295
--- /dev/null
+++ b/Simulation/DetSimSD/src/GenericTrackerSensitiveDetector.h
@@ -0,0 +1,24 @@
+// *********************************************************
+//
+// $Id: GenericTrackerSensitiveDetector.hh,v 1.0 2022/03/27
+
+#ifndef GenericTrackerSensitiveDetector_h
+#define GenericTrackerSensitiveDetector_h
+
+#include "DetSimSD/DDG4SensitiveDetector.h"
+#include "DDG4/Defs.h"
+
+class GenericTrackerSensitiveDetector: public DDG4SensitiveDetector {
+ public:
+  GenericTrackerSensitiveDetector(const std::string& name, dd4hep::Detector& description);
+  
+  void Initialize(G4HCofThisEvent* HCE);
+  G4bool ProcessHits(G4Step* step, G4TouchableHistory* history);
+  void EndOfEvent(G4HCofThisEvent* HCE);
+  
+ protected:
+
+  HitCollection* m_hc = nullptr;
+  
+};
+#endif
-- 
GitLab