ATLAS Offline Software
Public Member Functions | Private Member Functions | Private Attributes | List of all members
PixelSensorGmxSD Class Reference

#include <PixelSensorGmxSD.h>

Inheritance diagram for PixelSensorGmxSD:
Collaboration diagram for PixelSensorGmxSD:

Public Member Functions

 PixelSensorGmxSD (const std::string &name, const std::string &hitCollectionName, GeoModelIO::ReadGeoModel *sqlreader=nullptr)
 
virtual G4bool ProcessHits (G4Step *, G4TouchableHistory *) override final
 
virtual void Initialize (G4HCofThisEvent *) override final
 
template<class... Args>
void AddHit (Args &&... args)
 Templated method to stuff a single hit into the sensitive detector class. More...
 

Private Member Functions

 FRIEND_TEST (PixelSensorGmxSDtest, ProcessHits)
 
 FRIEND_TEST (PixelSensorGmxSDtest, AddHit)
 

Private Attributes

std::string m_HitCollName
 
SiHitCollectionm_HitColl {nullptr}
 
GeoModelIO::ReadGeoModel * m_sqlreader {nullptr}
 

Detailed Description

Definition at line 27 of file PixelSensorGmxSD.h.

Constructor & Destructor Documentation

◆ PixelSensorGmxSD()

PixelSensorGmxSD::PixelSensorGmxSD ( const std::string &  name,
const std::string &  hitCollectionName,
GeoModelIO::ReadGeoModel *  sqlreader = nullptr 
)

Definition at line 36 of file PixelSensorGmxSD.cxx.

37  : G4VSensitiveDetector( name )
38  , m_HitCollName( hitCollectionName )
39 {
40  m_sqlreader = sqlreader;
41 }

Member Function Documentation

◆ AddHit()

template<class... Args>
void PixelSensorGmxSD::AddHit ( Args &&...  args)
inline

Templated method to stuff a single hit into the sensitive detector class.

This could get rather tricky, but the idea is to allow fast simulations to use the very same SD classes as the standard simulation.

Definition at line 44 of file PixelSensorGmxSD.h.

44 { m_HitColl->Emplace( args... ); }

◆ FRIEND_TEST() [1/2]

PixelSensorGmxSD::FRIEND_TEST ( PixelSensorGmxSDtest  ,
AddHit   
)
private

◆ FRIEND_TEST() [2/2]

PixelSensorGmxSD::FRIEND_TEST ( PixelSensorGmxSDtest  ,
ProcessHits   
)
private

◆ Initialize()

void PixelSensorGmxSD::Initialize ( G4HCofThisEvent *  )
finaloverridevirtual

Definition at line 44 of file PixelSensorGmxSD.cxx.

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 }

◆ ProcessHits()

G4bool PixelSensorGmxSD::ProcessHits ( G4Step *  aStep,
G4TouchableHistory *   
)
finaloverridevirtual

Definition at line 56 of file PixelSensorGmxSD.cxx.

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 }

Member Data Documentation

◆ m_HitColl

SiHitCollection* PixelSensorGmxSD::m_HitColl {nullptr}
private

Definition at line 49 of file PixelSensorGmxSD.h.

◆ m_HitCollName

std::string PixelSensorGmxSD::m_HitCollName
private

Definition at line 48 of file PixelSensorGmxSD.h.

◆ m_sqlreader

GeoModelIO::ReadGeoModel* PixelSensorGmxSD::m_sqlreader {nullptr}
private

Definition at line 50 of file PixelSensorGmxSD.h.


The documentation for this class was generated from the following files:
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
python.CaloAddPedShiftConfig.args
args
Definition: CaloAddPedShiftConfig.py:47
SiHit::xPhi
@ xPhi
Definition: SiHit.h:162
AtlasHitsVector< SiHit >
python.SystemOfUnits.MeV
float MeV
Definition: SystemOfUnits.py:172
AtlasHitsVector::Emplace
void Emplace(Args &&... args)
Definition: AtlasHitsVector.h:80
TrackHelper
Definition: TrackHelper.h:14
lumiFormat.i
int i
Definition: lumiFormat.py:85
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
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
PixelSensorGmxSD::m_sqlreader
GeoModelIO::ReadGeoModel * m_sqlreader
Definition: PixelSensorGmxSD.h:50