ATLAS Offline Software
SctSensorSD.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 //
6 // SCT Sensitive Detector.
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 headers
12 #include "SctSensorSD.h"
13 
14 // athena includes
16 #include "MCTruth/TrackHelper.h"
17 
18 // Geant4 includes
19 #include "G4Step.hh"
20 #include "G4ThreeVector.hh"
21 #include "G4SDManager.hh"
22 #include "G4Geantino.hh"
23 #include "G4ChargedGeantino.hh"
24 #include <G4EventManager.hh>
25 
26 // CLHEP transform
27 #include "CLHEP/Geometry/Transform3D.h"
28 
29 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
30 
31 SctSensorSD::SctSensorSD( const std::string& name, const std::string& hitCollectionName )
32  : G4VSensitiveDetector( name )
33  , m_HitCollName( hitCollectionName )
34 {
35 }
36 
37 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
38 
39 void SctSensorSD::Initialize(G4HCofThisEvent *)
40 {
41  // ISF calls G4SDManager::PrepareNewEvent() before the Geant4 event loop starts...
42  if(auto* eventManger = G4EventManager::GetEventManager()){
43  if(auto* eventInfo = static_cast<AtlasG4EventUserInfo*>(eventManger->GetUserInformation())){
44  m_HitColl = eventInfo->GetHitCollectionMap()->Find<SiHitCollection>(m_HitCollName);
45  }
46  }
47 }
48 
49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50 
51 G4bool SctSensorSD::ProcessHits(G4Step* aStep, G4TouchableHistory* /*ROhist*/)
52 {
53  double edep = aStep->GetTotalEnergyDeposit();
54  if(edep==0.) {
55  if(aStep->GetTrack()->GetDefinition()!=G4Geantino::GeantinoDefinition() &&
56  aStep->GetTrack()->GetDefinition()!=G4ChargedGeantino::ChargedGeantinoDefinition())
57  return false;
58  }
59  edep *= CLHEP::MeV;
60  //
61  // Get the Touchable History:
62  //
63  const G4TouchableHistory *myTouch = dynamic_cast<const G4TouchableHistory*>(aStep->GetPreStepPoint()->GetTouchable());
64  if (not myTouch) {
65  G4cout << "SctSensorSD::ProcessHits bad dynamic_cast" << G4endl;
66  return false;
67  }
68  //
69  // Get the hit coordinates. Start and End Point
70  //
71  G4ThreeVector coord1 = aStep->GetPreStepPoint()->GetPosition();
72  G4ThreeVector coord2 = aStep->GetPostStepPoint()->GetPosition();
73  //
74  // Calculate the local step begin and end position.
75  // From a G4 FAQ:
76  // http://geant4-hn.slac.stanford.edu:5090/HyperNews/public/get/geometry/17/1.html
77  //
78  const G4AffineTransform transformation = myTouch->GetHistory()->GetTopTransform();
79  G4ThreeVector localPosition1 = transformation.TransformPoint(coord1);
80  G4ThreeVector localPosition2 = transformation.TransformPoint(coord2);
81  //
82  // Get it into a vector in local coords and with the right units:
83  //
84  HepGeom::Point3D<double> lP1,lP2;
85 
86  lP1[SiHit::xEta] = localPosition1[2]*CLHEP::mm;
87  lP1[SiHit::xPhi] = localPosition1[1]*CLHEP::mm;
88  lP1[SiHit::xDep] = localPosition1[0]*CLHEP::mm;
89 
90  lP2[SiHit::xEta] = localPosition2[2]*CLHEP::mm;
91  lP2[SiHit::xPhi] = localPosition2[1]*CLHEP::mm;
92  lP2[SiHit::xDep] = localPosition2[0]*CLHEP::mm;
93 
94  // Now Navigate the history to know in what detector the step is:
95  // and finally set the ID of det element in which the hit is.
96  //
97  //G4int History;
98  //
99  // Is it TB, barrel or endcap?
100  //
101  int brlEcap = 0;
102  int layerDisk = 0;
103  int etaMod = 0;
104  int phiMod = 0;
105  int side = 0;
106  this->indexMethod(myTouch, coord1.z(), brlEcap, layerDisk, etaMod, phiMod, side);
107  // get the HepMcParticleLink from the TrackHelper
108  TrackHelper trHelp(aStep->GetTrack());
109  m_HitColl->Emplace(lP1,
110  lP2,
111  edep,
112  aStep->GetPreStepPoint()->GetGlobalTime(),//use the global time. i.e. the time from the beginning of the event
113  trHelp.GenerateParticleLink(),
114  1,brlEcap,layerDisk,etaMod,phiMod,side);
115  return true;
116 }
117 
118 void SctSensorSD::indexMethod(const G4TouchableHistory *myTouch, double coord1z,
119  int &brlEcap, int &layerDisk, int &etaMod, int &phiMod, int &side) {
120 
121 
122  //
123  // In the case of the TB the positioning integers won't be initialized
124  // and the identifying integer will be zero. There is no need to do
125  // anything else
126 
127  const int sensorCopyNo = myTouch->GetVolume()->GetCopyNo();
128 
129  if (verboseLevel>5){
130  for (int i=0;i<myTouch->GetHistoryDepth();i++){
131  std::string detname=myTouch->GetVolume(i)->GetLogicalVolume()->GetName();
132  int copyno=myTouch->GetVolume(i)->GetCopyNo();
133  G4cout << "Volume "<<detname << " Copy Nr. " << copyno << G4endl;
134  }
135  }
136 
137  // Barrel geometry from GeoModel. Flag it with a 1000 in the IdTag
138  //
139  if(sensorCopyNo/1000) {
140  // Barrel
141  if(sensorCopyNo == 1000) {
142  side = myTouch->GetVolume(1)->GetCopyNo();
143  etaMod = myTouch->GetVolume(2)->GetCopyNo();
144  phiMod = myTouch->GetVolume(3)->GetCopyNo();
145  layerDisk = myTouch->GetVolume(5)->GetCopyNo();
146  } else if(sensorCopyNo == 1100) {
147  // SLHC stave layout
148  int etaModTmp = myTouch->GetVolume(1)->GetCopyNo();
149  int sign = (etaModTmp>=0)? +1 : -1;
150  side = std::abs(etaModTmp)%2;
151  etaMod = sign * std::abs(etaModTmp)/2;
152  phiMod = myTouch->GetVolume(2)->GetCopyNo();
153  layerDisk = myTouch->GetVolume(4)->GetCopyNo();
154  } else if(sensorCopyNo/100 == 12) {
155  // Segmented sensor (SLHC only)
156  int iSegment = sensorCopyNo%100;
157  int sensorEnvCopyNo = myTouch->GetVolume(1)->GetCopyNo();
158  if (sensorEnvCopyNo == 1000) {
159  side = myTouch->GetVolume(2)->GetCopyNo();
160  etaMod = myTouch->GetVolume(3)->GetCopyNo() + iSegment;
161  phiMod = myTouch->GetVolume(4)->GetCopyNo();
162  layerDisk = myTouch->GetVolume(6)->GetCopyNo();
163  } else if (sensorEnvCopyNo == 1100) {
164  // Stave layout
165  int etaModTmp = myTouch->GetVolume(2)->GetCopyNo();
166  int sign = (etaModTmp>=0)? +1 : -1;
167  side = std::abs(etaModTmp)%2;
168  etaMod = sign * std::abs(etaModTmp)/2 + iSegment;
169  phiMod = myTouch->GetVolume(3)->GetCopyNo();
170  layerDisk = myTouch->GetVolume(5)->GetCopyNo();
171  } else {
172  G4ExceptionDescription description;
173  description << "ProcessHits: Unrecognized geometry in SCT sensitive detector. Please contact the maintainer of the SCT Detector Description.";
174  G4Exception("SctSensorSD", "UnrecognizedSCTGeometry1", FatalException, description);
175  abort();
176  }
177  } else {
178  G4ExceptionDescription description;
179  description << "ProcessHits: Unrecognized geometry in SCT sensitive detector. Please contact the maintainer of the SCT Detector Description.";
180  G4Exception("SctSensorSD", "UnrecognizedSCTGeometry2", FatalException, description);
181  abort();
182  }
183  if (verboseLevel>5){
184  G4cout << "In the SCT Barrel" << G4endl;
185  G4cout << "----- side : " << side << G4endl;
186  G4cout << "----- eta_module : " << etaMod << G4endl;
187  G4cout << "----- phi_module : " << phiMod << G4endl;
188  G4cout << "----- layer : " << layerDisk << G4endl;
189  G4cout << "----- barrel_ec : " << brlEcap<< G4endl;
190  }
191  } else {
192  // Endcap geometry from GeoModel. Flag it with a 500 or 600 in the IdTag
193  // The level in the hierarchy is different for different geometries
194  // 500 For pre DC3
195  // 600 For DC3 and later (including SLHC)
196  int sensorCopyNoTest = sensorCopyNo/100;
197  if(sensorCopyNoTest == 5 || sensorCopyNoTest == 6) {
198  int signZ = (coord1z > 0) ? 1 : -1;
199  brlEcap = 2 * signZ;
200  side = myTouch->GetVolume(0)->GetCopyNo() % 2;
201  if (sensorCopyNoTest == 5) { // DC2 and Rome
202  etaMod = myTouch->GetVolume(3)->GetCopyNo();
203  layerDisk = myTouch->GetVolume(4)->GetCopyNo();
204  } else { // DC3 and later and SLHC
205  etaMod = myTouch->GetVolume(2)->GetCopyNo();
206  layerDisk = myTouch->GetVolume(3)->GetCopyNo();
207  }
208  int phiNumber = myTouch->GetVolume(1)->GetCopyNo();
209  // Number of modules in ring and module which becomes 0 is encoded in the copy number.
210  // This is in order to recreate the correct id in the negative endcap.
211  int maxModules = (phiNumber & 0x00ff0000) >> 16;
212  int moduleZero = (phiNumber & 0xff000000) >> 24; // When phiMod = moduleZero -> 0 after inverting
213  phiMod = phiNumber & 0x0000ffff;
214  // If -ve endcap then invert numbering
215  // maxModules is non-zero, but check for maxModules != 0 safeguards
216  // against divide by zero (in case we create a geometry without the encoding).
217  if (maxModules && signZ < 0) phiMod = (maxModules + moduleZero - phiMod) % maxModules;
218 
219  // Some printout
220  if (verboseLevel>5){
221  G4cout << "In the SCT EndCap" << G4endl;
222  G4cout << "----- side : " << side << G4endl;
223  G4cout << "----- phi_module : " << phiMod<< G4endl;
224  G4cout << "----- eta_module : " << etaMod << G4endl;
225  G4cout << "----- layerDisk : " << layerDisk << G4endl;
226  G4cout << "----- barrel_ec : " << brlEcap << G4endl;
227  }
228  } else {
229  // Do not expect other numbers. Need to fix SCT_GeoModel or SCT_SLHC_GeoModel if this occurs.
230  G4ExceptionDescription description;
231  description << "ProcessHits: Unrecognized geometry in SCT sensitive detector. Please contact the maintainer of the SCT Detector Description.";
232  G4Exception("SctSensorSD", "UnrecognizedSCTGeometry3", FatalException, description);
233  abort();
234  }
235  }
236  return;
237 }
SctSensorSD::ProcessHits
G4bool ProcessHits(G4Step *, G4TouchableHistory *) override
Definition: SctSensorSD.cxx:51
AtlasG4EventUserInfo
This class is attached to G4Event objects as UserInformation. It holds a pointer to the HepMC::GenEve...
Definition: AtlasG4EventUserInfo.h:23
SctSensorSD::m_HitColl
SiHitCollection * m_HitColl
Definition: SctSensorSD.h:51
SiHit::xDep
@ xDep
Definition: SiHit.h:162
python.SystemOfUnits.mm
float mm
Definition: SystemOfUnits.py:98
SctSensorSD::SctSensorSD
SctSensorSD(const std::string &name, const std::string &hitCollectionName)
Definition: SctSensorSD.cxx:31
TrackHelper.h
SiHit::xPhi
@ xPhi
Definition: SiHit.h:162
AtlasHitsVector< SiHit >
SctSensorSD.h
python.SystemOfUnits.MeV
float MeV
Definition: SystemOfUnits.py:172
TRT::Hit::side
@ side
Definition: HitInfo.h:83
SctSensorSD::m_HitCollName
std::string m_HitCollName
Definition: SctSensorSD.h:50
AtlasHitsVector::Emplace
void Emplace(Args &&... args)
Definition: AtlasHitsVector.h:80
TrackHelper
Definition: TrackHelper.h:14
lumiFormat.i
int i
Definition: lumiFormat.py:85
sign
int sign(int a)
Definition: TRT_StrawNeighbourSvc.h:107
SctSensorSD::Initialize
void Initialize(G4HCofThisEvent *) override final
Definition: SctSensorSD.cxx:39
SctSensorSD::indexMethod
void indexMethod(const G4TouchableHistory *myTouch, double coord1z, int &brlEcap, int &layerDisk, int &etaMod, int &phiMod, int &side)
Definition: SctSensorSD.cxx:118
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
SiHit::xEta
@ xEta
Definition: SiHit.h:162
AtlasG4EventUserInfo.h
TrackHelper::GenerateParticleLink
HepMcParticleLink GenerateParticleLink()
Generates a creates new HepMcParticleLink object on the stack based on GetUniqueID(),...
Definition: TrackHelper.h:35
fillSCTHists.etaMod
etaMod
Definition: fillSCTHists.py:23
description
std::string description
glabal timer - how long have I taken so far?
Definition: hcg.cxx:88