ATLAS Offline Software
PixelSensorGmxSD.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 //
6 // Pixel Sensitive Detector - specialisation for GeoModelXml
7 // The Hits are processed here. For every hit I get the position and
8 // an information on the sensor in which the interaction happened
9 //
10 
11 // Class header
12 #include "PixelSensorGmxSD.h"
13 
14 // Athena headers
16 #include "MCTruth/TrackHelper.h"
17 
18 // Geant4 headers
19 #include "G4ChargedGeantino.hh"
20 #include <G4EventManager.hh>
21 #include "G4Geantino.hh"
22 #include "G4SDManager.hh"
23 #include "G4Step.hh"
24 #include "G4ThreeVector.hh"
25 
26 // CLHEP headers
27 #include "CLHEP/Geometry/Transform3D.h"
28 #include "CLHEP/Units/SystemOfUnits.h"
29 
30 #include <GeoModelKernel/GeoFullPhysVol.h>
31 #include <GeoModelRead/ReadGeoModel.h>
32 
34 
35 
36 PixelSensorGmxSD::PixelSensorGmxSD(const std::string& name, const std::string& hitCollectionName,GeoModelIO::ReadGeoModel * sqlreader)
37  : G4VSensitiveDetector( name )
38  , m_HitCollName( hitCollectionName )
39 {
40  m_sqlreader = sqlreader;
41 }
42 
43 // Initialize from G4 - cache the hit collection for the current event
44 void PixelSensorGmxSD::Initialize(G4HCofThisEvent *)
45 {
46  // ISF calls G4SDManager::PrepareNewEvent() before the Geant4 event loop starts...
47  if(auto* eventManger = G4EventManager::GetEventManager()){
48  if(auto* eventInfo = static_cast<AtlasG4EventUserInfo*>(eventManger->GetUserInformation())){
49  m_HitColl = eventInfo->GetHitCollectionMap()->Find<SiHitCollection>(m_HitCollName);
50  }
51  }
52 
53 }
54 
55 
56 G4bool PixelSensorGmxSD::ProcessHits(G4Step* aStep, G4TouchableHistory* /*ROhist*/)
57 {
58  if (verboseLevel>5) G4cout << "Process Hit" << G4endl;
59 
60  G4double edep = aStep->GetTotalEnergyDeposit();
61  edep *= CLHEP::MeV;
62  if(edep==0.) {
63  if(aStep->GetTrack()->GetDefinition() != G4Geantino::GeantinoDefinition() &&
64  aStep->GetTrack()->GetDefinition() != G4ChargedGeantino::ChargedGeantinoDefinition())
65  return false;
66  }
67 
68  //use the global time. i.e. the time from the beginning of the event
69  //
70  // Get the Touchable History:
71  //
72  const G4TouchableHistory *myTouch = dynamic_cast<const G4TouchableHistory*>(aStep->GetPreStepPoint()->GetTouchable());
73  if (not myTouch) {
74  G4cout << "PixelSensorGmxSD::ProcessHits bad dynamic_cast" << G4endl;
75  return false;
76  }
77  if(verboseLevel>5){
78  for (int i=0;i<myTouch->GetHistoryDepth();i++){
79  std::string detname=myTouch->GetVolume(i)->GetLogicalVolume()->GetName();
80  int copyno=myTouch->GetVolume(i)->GetCopyNo();
81  G4cout << "Volume " <<detname <<" Copy Nr. " << copyno << G4endl;
82  }
83  }
84  //
85  // Get the hit coordinates. Start and End Point
86  //
87  G4ThreeVector coord1 = aStep->GetPreStepPoint()->GetPosition();
88  G4ThreeVector coord2 = aStep->GetPostStepPoint()->GetPosition();
89 
90  // Calculate the local step begin and end position.
91  // From a G4 FAQ:
92  // http://geant4-hn.slac.stanford.edu:5090/HyperNews/public/get/geometry/17/1.html
93  //
94  const G4AffineTransform transformation = myTouch->GetHistory()->GetTopTransform();
95  G4ThreeVector localPosition1 = transformation.TransformPoint(coord1);
96  G4ThreeVector localPosition2 = transformation.TransformPoint(coord2);
97 
98  HepGeom::Point3D<double> lP1,lP2;
99  lP1[SiHit::xEta] = localPosition1[2]*CLHEP::mm;
100  lP1[SiHit::xPhi] = localPosition1[1]*CLHEP::mm;
101  lP1[SiHit::xDep] = localPosition1[0]*CLHEP::mm;
102 
103  lP2[SiHit::xEta] = localPosition2[2]*CLHEP::mm;
104  lP2[SiHit::xPhi] = localPosition2[1]*CLHEP::mm;
105  lP2[SiHit::xDep] = localPosition2[0]*CLHEP::mm;
106 
107  TrackHelper trHelp(aStep->GetTrack());
108 
109  if(m_sqlreader){
110  //if sqlite inputs, Identifier indices come from PhysVol Name
111  std::string physVolName = myTouch->GetVolume(0)->GetName();
112 
113  int hitIdOfWafer = SiHitIdHelper::GetHelper()->buildHitIdFromStringITk(0,physVolName);
114 
115 
116  m_HitColl->Emplace(lP1,
117  lP2,
118  edep,
119  aStep->GetPreStepPoint()->GetGlobalTime(),//use the global time. i.e. the time from the beginning of the event
120  trHelp.GenerateParticleLink(),
121  hitIdOfWafer);
122  return true;
123 
124  }
125  // if not from SQLite, we assume that the Identifier has already been written in as the copy number
126  // (it should have done if GeoModel building ran within Athena)
127 
128  int id = myTouch->GetVolume()->GetCopyNo();
129 
130  m_HitColl->Emplace(lP1,
131  lP2,
132  edep,
133  aStep->GetPreStepPoint()->GetGlobalTime(),
134  trHelp.GenerateParticleLink(),
135  id);
136  return true;
137 }
AtlasG4EventUserInfo
This class is attached to G4Event objects as UserInformation. It holds a pointer to the HepMC::GenEve...
Definition: AtlasG4EventUserInfo.h:23
SiHit::xDep
@ xDep
Definition: SiHit.h:162
python.SystemOfUnits.mm
float mm
Definition: SystemOfUnits.py:98
TrackHelper.h
SiHit::xPhi
@ xPhi
Definition: SiHit.h:162
AtlasHitsVector< SiHit >
PixelSensorGmxSD::PixelSensorGmxSD
PixelSensorGmxSD(const std::string &name, const std::string &hitCollectionName, GeoModelIO::ReadGeoModel *sqlreader=nullptr)
Definition: PixelSensorGmxSD.cxx:36
python.SystemOfUnits.MeV
float MeV
Definition: SystemOfUnits.py:172
SiHitIdHelper.h
AtlasHitsVector::Emplace
void Emplace(Args &&... args)
Definition: AtlasHitsVector.h:80
TrackHelper
Definition: TrackHelper.h:14
lumiFormat.i
int i
Definition: lumiFormat.py:85
PixelSensorGmxSD::Initialize
virtual void Initialize(G4HCofThisEvent *) override final
Definition: PixelSensorGmxSD.cxx:44
PixelSensorGmxSD.h
SiHitIdHelper::buildHitIdFromStringITk
int buildHitIdFromStringITk(int part, const std::string &) const
Definition: SiHitIdHelper.cxx:131
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
SiHit::xEta
@ xEta
Definition: SiHit.h:162
AtlasG4EventUserInfo.h
PixelSensorGmxSD::ProcessHits
virtual G4bool ProcessHits(G4Step *, G4TouchableHistory *) override final
Definition: PixelSensorGmxSD.cxx:56
PixelSensorGmxSD::m_HitColl
SiHitCollection * m_HitColl
Definition: PixelSensorGmxSD.h:49
SiHitIdHelper::GetHelper
static const SiHitIdHelper * GetHelper()
Definition: SiHitIdHelper.cxx:19
PixelSensorGmxSD::m_HitCollName
std::string m_HitCollName
Definition: PixelSensorGmxSD.h:48
TrackHelper::GenerateParticleLink
HepMcParticleLink GenerateParticleLink()
Generates a creates new HepMcParticleLink object on the stack based on GetUniqueID(),...
Definition: TrackHelper.h:35
PixelSensorGmxSD::m_sqlreader
GeoModelIO::ReadGeoModel * m_sqlreader
Definition: PixelSensorGmxSD.h:50