ATLAS Offline Software
Loading...
Searching...
No Matches
RPCSensitiveDetector Class Reference

#include <RPCSensitiveDetector.h>

Inheritance diagram for RPCSensitiveDetector:
Collaboration diagram for RPCSensitiveDetector:

Public Member Functions

 RPCSensitiveDetector (const std::string &name, const std::string &hitCollectionName, unsigned int nGasGaps)
 construction/destruction
void Initialize (G4HCofThisEvent *) override final
 member functions
G4bool ProcessHits (G4Step *, G4TouchableHistory *) override final

Private Member Functions

 FRIEND_TEST (RPCSensitiveDetectortest, Initialize)
 FRIEND_TEST (RPCSensitiveDetectortest, ProcessHits)

Private Attributes

std::string m_hitCollectionName
 member data
RPCSimHitCollectionm_myRPCHitColl {nullptr}
AtlasG4EventUserInfom_g4UserEventInfo {nullptr}
const RpcHitIdHelperm_muonHelper {nullptr}

Detailed Description

Author
Andre.nosp@m.a.Di.nosp@m.Simon.nosp@m.e@ce.nosp@m.rn.ch

Class methods and properties

The method RPCSensitiveDetector::ProcessHits is executed by the G4 kernel each time a particle crosses one of the RPC gas gaps. Navigating with the touchableHistory method GetHistoryDepth() through the hierarchy of volumes crossed by the particle, the Sensitive Detector determines the correct set of Simulation Identifiers to associate to each hit. The RPC SimIDs are 32-bit unsigned integers, built using the MuonSimEvent/RpcHitIdHelper class which inherits from the MuonHitIdHelper base class.

We describe here how each field of the identifier is determined.

1) stationName, stationEta, stationPhi: when a volume is found in the hierarchy whose name contains the substring "station", the stationName is extracted from the volume's name. stationEta and stationPhi values are calculated starting from the volume's copy number.

2) doubletR: when a volume is found in the hierarchy whose name contains the substring "rpccomponent", its copy number is used to calulate the doubletR of the hit.

3) gasGap: the substring "layer" is searched through the names of the volumes in the hierarchy. When a the volume is found, its copy number is used, together with the sign of the copy number of the "rpccomponent" to decide what is the correct gasGap value.

4) doubletPhi: when the volume with the substring "gas volume" in its name is found, its copy number is enough to determine the doubletPhi of the hit if the hit is registered in a standard chamber. For special chambers (i.e. chambers with only one gas gap per layer readout by two strip panels per direction, or chambers in the ribs) some a special attribution is done using also the stationName.

5) doubletZ: "gazGap" is the required substring in the volume's name. doubletZ attribution is fairly complicated. For standard chambers it just uses the sign of the copy number of the "rpccomponent", but severa arrangements (based on the stationName and the technology name) are needed for special chambers.

6) the doubletPhi calculated as described above needs one further correction for some chambers to take into account chamber rotation before placement in the spectrometer. This is done just before creating the hit.

notes:

1) presently, chamber efficiency is assumed to be 100% at the level of the Sensitive Detector, i.e. two hits (one in eta and one in phi) are created each time the sensitive detector is created.

2) the hits created as described above contain information about their local position in the gas gap and the simulation identifier of the strip panel

3) the strip number is not assigned by the sensitive detector, and must be calculated by the digitization algorithm.

4) points 2 and 3 are due to the necessity of not introducing any dependency on the geometry description in the Sensitive Detector.

5) the present version of the Sensitive Detector produces hits which are different depending on whether hand coded G4 geometry or GeoModel is used for the geometrical description of the setup. When hand coded geometry is used, both GeoModel and the old MuonDetDescr can be used for digitization (using a proper tag of RPC_Digitization), while if GeoModel is used in the simulation, it must be used also in the digitization.

6) for each hit, the time of flight (the G4 globalTime), is recorded and associated to the hit.

7) the RPCHit object contains: the SimID, the globalTime, the hit local position and the track barcode.

Class methods and properties

The method RPCSensitiveDetector::ProcessHits is executed by the G4 kernel each time a particle crosses one of the RPC gas gaps. Navigating with the touchableHistory method GetHistoryDepth() through the hierarchy of volumes crossed by the particle, the Sensitive Detector determines the correct set of Simulation Identifiers to associate to each hit. The RPC SimIDs are 32-bit unsigned integers, built using the MuonSimEvent/RpcHitIdHelper class which inherits from the MuonHitIdHelper base class.

We describe here how each field of the identifier is determined.

1) stationName, stationEta, stationPhi: when a volume is found in the hierarchy whose name contains the substring "station", the stationName is extracted from the volume's name. stationEta and stationPhi values are calculated starting from the volume's copy number.

2) doubletR: when a volume is found in the hierarchy whose name contains the substring "rpccomponent", its copy number is used to calulate the doubletR of the hit.

3) gasGap: the substring "layer" is searched through the names of the volumes in the hierarchy. When a the volume is found, its copy number is used, together with the sign of the copy number of the "rpccomponent" to decide what is the correct gasGap value.

4) doubletPhi: when the volume with the substring "gas volume" in its name is found, its copy number is enough to determine the doubletPhi of the hit if the hit is registered in a standard chamber. For special chambers (i.e. chambers with only one gas gap per layer readout by two strip panels per direction, or chambers in the ribs) some a special attribution is done using also the stationName.

5) doubletZ: "gazGap" is the required substring in the volume's name. doubletZ attribution is fairly complicated. For standard chambers it just uses the sign of the copy number of the "rpccomponent", but severa arrangements (based on the stationName and the technology name) are needed for special chambers.

6) the doubletPhi calculated as described above needs one further correction for some chambers to take into account chamber rotation before placement in the spectrometer. This is done just before creating the hit.

notes:

1) presently, chamber efficiency is assumed to be 100% at the level of the Sensitive Detector, i.e. two hits (one in eta and one in phi) are created each time the sensitive detector is created.

2) the hits created as described above contain information about their local position in the gas gap and the simulation identifier of the strip panel

3) the strip number is not assigned by the sensitive detector, and must be calculated by the digitization algorithm.

4) points 2 and 3 are due to the necessity of not introducing any dependency on the geometry description in the Sensitive Detector.

5) the present version of the Sensitive Detector produces hits which are different depending on whether hand coded G4 geometry or GeoModel is used for the geometrical description of the setup. When hand coded geometry is used, both GeoModel and the old MuonDetDescr can be used for digitization (using a proper tag of RPC_Digitization), while if GeoModel is used in the simulation, it must be used also in the digitization.

6) for each hit, the time of flight (the G4 globalTime), is recorded and associated to the hit.

7) the RPCHit object contains: the SimID, the globalTime, the hit local position and the track barcode.

Definition at line 98 of file RPCSensitiveDetector.h.

Constructor & Destructor Documentation

◆ RPCSensitiveDetector()

RPCSensitiveDetector::RPCSensitiveDetector ( const std::string & name,
const std::string & hitCollectionName,
unsigned int nGasGaps )

construction/destruction

Definition at line 25 of file RPCSensitiveDetector.cxx.

26 : G4VSensitiveDetector( name )
27 , m_hitCollectionName( hitCollectionName )
28{
30}
std::string m_hitCollectionName
member data
const RpcHitIdHelper * m_muonHelper
static const RpcHitIdHelper * GetHelper(unsigned int nGasGaps=2)

Member Function Documentation

◆ FRIEND_TEST() [1/2]

RPCSensitiveDetector::FRIEND_TEST ( RPCSensitiveDetectortest ,
Initialize  )
private

◆ FRIEND_TEST() [2/2]

RPCSensitiveDetector::FRIEND_TEST ( RPCSensitiveDetectortest ,
ProcessHits  )
private

◆ Initialize()

void RPCSensitiveDetector::Initialize ( G4HCofThisEvent * )
finaloverride

member functions

Definition at line 32 of file RPCSensitiveDetector.cxx.

33{
34 m_myRPCHitColl = nullptr;
35 if (auto* eventInfo = AtlasG4EventUserInfo::GetEventUserInfo()) {
36 m_myRPCHitColl = eventInfo->GetHitCollectionMap()->Find<RPCSimHitCollection>(m_hitCollectionName);
37 m_g4UserEventInfo = eventInfo;
38 }
39 //FIXME probably only need to call this bit at start of the event
40 //loop rather than the start of each G4Event.
41 if (verboseLevel>5) G4cout << "Initializing SD" << G4endl;
42
43}
AtlasHitsVector< RPCSimHit > RPCSimHitCollection
static AtlasG4EventUserInfo * GetEventUserInfo()
AtlasG4EventUserInfo * m_g4UserEventInfo
RPCSimHitCollection * m_myRPCHitColl

◆ ProcessHits()

G4bool RPCSensitiveDetector::ProcessHits ( G4Step * aStep,
G4TouchableHistory *  )
finaloverride

RPCs sensitive to charged particle only

Reject hits that are parallel through the gas gap

Definition at line 45 of file RPCSensitiveDetector.cxx.

45 {
46
47 if (!m_myRPCHitColl) {
48 G4Exception("RPCSensitiveDetector::ProcessHits", "RPCHitCollectionMissing", FatalException,
49 "Hit collection not initialized; did SetupEvent run?");
50 return false;
51 }
52
53 G4Track* track = aStep->GetTrack();
54
56 if (track->GetDefinition()->GetPDGCharge() == 0.0) {
57 if (track->GetDefinition()!=G4Geantino::GeantinoDefinition()) {
58 return true;
59 }
60 }
61
62 const G4TouchableHistory* touchHist = static_cast<const G4TouchableHistory*>(aStep->GetPreStepPoint()->GetTouchable());
63 G4ThreeVector position = aStep->GetPreStepPoint()->GetPosition();
64 G4ThreeVector postPosition = aStep->GetPostStepPoint()->GetPosition();
65 const G4AffineTransform trans = track->GetTouchable()->GetHistory()->GetTopTransform(); // from global to local
66
67 // necessary to assign correct identifiers
68 int rpcIsRotated = 0;
69 // int layer = 0;
70
71 std::string tech;
72 // fields for the RPC identifier construction
73 std::string stationName;
74 int stationEta= 0;
75 int stationPhi= 0;
76 int doubletR= 0;
77 // int doubletZ= 0;
78 int doubletPhi= 0;
79 int gasGap= 0;
80 // int eta_inversed=0;
81 // int phi_inversed=0;
82 int station_rotated=0; // tells us if the station was rotated before being positioned.
83 // if =1 we have to correct some IDs
84 int zNeg_original=0; // tells if the station at z<0 was obtained duplicating a station at z>0
85
86 // RPC hit information
87 double globalTime = aStep->GetPreStepPoint()->GetGlobalTime();
88 Amg::Vector3D localPosition = Amg::Hep3VectorToEigen( trans.TransformPoint(position) );
89 Amg::Vector3D localPostPosition = Amg::Hep3VectorToEigen( trans.TransformPoint(postPosition) );
90 {
91 const Amg::Vector3D stepVector = localPostPosition - localPosition;
93 if (stepVector.mag()>std::numeric_limits<float>::epsilon() &&
94 std::abs(std::abs(Amg::angle(stepVector, Amg::Vector3D::UnitX()))
95 - 90.*Gaudi::Units::deg) < 0.0001* Gaudi::Units::deg) {
96 return true;
97 }
98 }
99 int mydbZ=0;
100 int mydbPMod=0;
101 int mydbP=0;
102 bool isAssembly = false;
103 // scan geometry tree to identify the hit channel
104 for (int i=touchHist->GetHistoryDepth();i>=0;i--) {
105
106 std::string::size_type npos;
107 std::string volName = touchHist->GetVolume(i)->GetName();
108 std::string num=volName.substr(3,2);
109 if(num[0]==' ') num[0]=0;
110
111 // check if this station is an assembly
112 if ((npos = volName.find("av_")) != std::string::npos &&
113 (npos = volName.find("impr_")) != std::string::npos) isAssembly = true;
114
115 // stationName, stationEta, stationPhi
116 if ((npos = volName.find("station")) != std::string::npos && (!isAssembly)) {
117
118 stationName = volName.substr(0,npos-1);
119
120 int volCopyNo = touchHist->GetVolume(i)->GetCopyNo();
121
122 if(abs(volCopyNo/1000)==1){
123 zNeg_original=1;
124 volCopyNo=volCopyNo%1000;
125 }
126
127 stationEta = volCopyNo/100;
128 stationPhi = abs(volCopyNo%100);
129
130 if(stationEta<0&&!zNeg_original) station_rotated=1;
131
132 // stationName, stationEta, stationPhi
133 } else if ((npos = volName.find("RPC")) != std::string::npos && isAssembly) {
134 // vol name for Assembly components are
135 // av_WWW_impr_XXX_Muon::BMSxMDTxx_pv_ZZZ_NAME
136 // where WWW is ass. istance nr.
137 // XXX is comp. imprint nr.
138 // BMSxMDTxx is the name of the comp. log.Vol.
139 // x station sub-type; xx technology subtype
140 // ZZZ is the comp. number of order
141 // NAME is the comp. tag (geoIdentifierTag)
142 // for RPCs it is SwapdbPdbRdbZ[aaa]
143 // with aaa = doubletZ + doubletR*100 + dbphi*1000;
144 // aaa = -aaa if iswap == -1
145 // dbphi = 1 always but for S1 or S2 at doubletPhi = 2
146 // copy numbers for Ass.components are =
147 // CopyNoBase(= geoIdentifierTag of the assembly) + geoIdentifierTag of the component
148 // geoIdentifierTag of the component = Job
149 // geoIdentifierTag of the assembly = (sideC*10000 +
150 // mirsign*1000 + abs(zi)*100 + fi+1)*100000;
151 // mirsign = 1 if zi<0 and !mirrored; otherwise 0
152 std::string::size_type loc1,loc2;
153 if ((loc1 = volName.find("Muon::")) != std::string::npos) {
154 stationName = volName.substr(loc1+6,4); //type and subtype
155 }
156
157 int volCopyNo = touchHist->GetVolume(i)->GetCopyNo();
158 int copyNrBase = int(volCopyNo/100000);
159 int sideC = int(copyNrBase/10000);
160 int zi = int((copyNrBase%1000)/100);
161 int mirfl = int((copyNrBase%10000)/1000);
162 int fi = int(copyNrBase%100);
163 if (sideC == 1) zi = -zi;
164 stationEta = zi;
165 stationPhi = fi;
166 zNeg_original = mirfl;
167
168 if(stationEta<0&&!zNeg_original) station_rotated=1;
169
170 // doubletZ = 1;
171 tech=volName.substr(npos,5);
172
173 // now get the geoIdentifierTag of the rpc components
174 int gmID = 0;
175 if ((loc1 = volName.find('[')) != std::string::npos) {
176 if ((loc2 = volName.find(']', loc1+1)) != std::string::npos) {
177 //G4cout << "first [ is at "<<loc1<<" first ] at "<<loc2 << G4endl;
178 std::istringstream istrvar(volName.substr(loc1+1,loc2-loc1-1));
179 istrvar>>gmID;
180 }
181 }
182 int kk=gmID;
183
184 if (kk < 0) rpcIsRotated=1;
185
186 doubletR =(abs(kk)%1000)/100;
187 mydbZ = abs(int(kk%10));
188 mydbPMod = abs(int(kk/1000));
189 // doubleR, doubletZ
190 } else if ((npos = volName.find("rpccomponent")) != std::string::npos && (!isAssembly)) {
191
192 std::string::size_type loc1,loc2;
193 tech=volName.substr(npos-5,5);
194 int gmID = 0;
195 if ((loc1 = volName.find('[')) != std::string::npos) {
196 if ((loc2 = volName.find(']', loc1+1)) != std::string::npos) {
197 std::istringstream istrvar(volName.substr(loc1+1,loc2-loc1-1));
198 istrvar>>gmID;
199 }
200 }
201 mydbZ = abs(int(gmID%10));
202 mydbPMod = abs(int(gmID/1000));
203
204 int kk=touchHist->GetVolume(i)->GetCopyNo();
205
206 if (kk < 0) rpcIsRotated=1;
207
208 doubletR=(abs(kk)%1000)/100;
209
210 } else if ((npos = volName.find("layer")) != std::string::npos) {
211
212 int copyNo = touchHist->GetVolume(i)->GetCopyNo();
213
214 if (copyNo == 1) {
215 rpcIsRotated ? gasGap = 2 : gasGap = 1;
216 } else if (copyNo ==2) {
217 rpcIsRotated ? gasGap = 1 : gasGap = 2;
218 } else if (copyNo ==3) {
219 gasGap = 3;
220 }
221 } else if((npos = volName.find("gas volume")) != std::string::npos) {
222
223 int copyNo = touchHist->GetVolume(i)->GetCopyNo();
224 if (copyNo == 0) {
225 doubletPhi = 1;
226 } else if (copyNo == 10) {
227 doubletPhi = 2;
228 }
229 mydbP = doubletPhi;
230 int ngap_in_s=0;
231 int nstrippanel_in_s=0;
232 std::string::size_type loc1;
233 if ((loc1 = volName.find("gg_in_s")) != std::string::npos) {
234 std::istringstream istrvar(volName.substr(loc1-1,1));
235 istrvar>>ngap_in_s;
236 }
237 if ((loc1 = volName.find("sp_in_s")) != std::string::npos) {
238 std::istringstream istrvar(volName.substr(loc1-1,1));
239 istrvar>>nstrippanel_in_s;
240 }
241 if (ngap_in_s == 1 && nstrippanel_in_s == 2) {
242 if(localPosition.y()>0) mydbP=2;
243 else mydbP=1;
244 } else if (ngap_in_s == 1 && nstrippanel_in_s == 1) {
245 mydbP = mydbPMod;
246 }
247 }
248 }
249
251 // correct wrong IDs due to station rotation. only doubletZ and doubletPhi are affected at this point
253
254 if(station_rotated){
255 // doubletPhi correction
256 // note that this is correct also for ribs chambers
257 mydbP++;
258 if (mydbP>2) mydbP=1;
259 // G4cout << "Rcd PC - swap dbPhi because ch is MIRRORED "<<doubletPhi << G4endl;
260
261 // strip numbering: if the station is rotated both eta and phi directions get inversed.
262 // commented out for geomodel!
263 }
264
266
267 // now we have the position in the gas gap, with correct axis orientation, and the offlineIDs *of the strip panel*
268 // nstrip is supposed to be calculated by RPC_Digitizer.
269
270 // construct two (one in eta and one in phi) new RPC hits and store them in hit collection
271
272 //construct the hit identifiers
273 if (verboseLevel>5) {
274 G4cout << "hit in station "<<stationName<< " on technology "<<tech << G4endl;
275 G4cout << "constructing ids (stName, stEta, stPhi, dr, dZ, dPhi)= "<<stationName<< " "<< stationEta<<" " << stationPhi<< " "<<doubletR<< " "<< mydbZ<< " "<<mydbP << G4endl;
276 }
277
278 HitID RPCid_eta = m_muonHelper->BuildRpcHitId(stationName, stationPhi, stationEta,
279 mydbZ, doubletR, gasGap, mydbP,0);
280
281 HitID RPCid_phi = m_muonHelper->BuildRpcHitId(stationName, stationPhi, stationEta,
282 mydbZ, doubletR, gasGap, mydbP,1);
283
284 // retrieve track barcode
285 TrackHelper trHelp(aStep->GetTrack());
286
287 //construct new rpc hit
288 m_myRPCHitColl->Emplace(RPCid_eta, globalTime,
289 localPosition,
290 trHelp.GenerateParticleLink(m_g4UserEventInfo ? m_g4UserEventInfo->GetEventStore() : nullptr),
291 localPostPosition,
292 aStep->GetTotalEnergyDeposit(),
293 aStep->GetStepLength(),
294 track->GetDefinition()->GetPDGEncoding(),
295 aStep->GetPreStepPoint()->GetKineticEnergy());
296 m_myRPCHitColl->Emplace(RPCid_phi, globalTime,
297 localPosition,
298 trHelp.GenerateParticleLink(m_g4UserEventInfo ? m_g4UserEventInfo->GetEventStore() : nullptr),
299 localPostPosition,
300 aStep->GetTotalEnergyDeposit(),
301 aStep->GetStepLength(),
302 track->GetDefinition()->GetPDGEncoding(),
303 aStep->GetPreStepPoint()->GetKineticEnergy());
304
305 return true;
306}
int HitID
double angle(const Amg::Vector3D &v1, const Amg::Vector3D &v2)
calculates the opening angle between two vectors
Amg::Vector3D Hep3VectorToEigen(const CLHEP::Hep3Vector &CLHEPvector)
Converts a CLHEP-based CLHEP::Hep3Vector into an Eigen-based Amg::Vector3D.
Eigen::Matrix< double, 3, 1 > Vector3D
constexpr uint8_t stationPhi
station Phi 1 to 8
unsigned int constexpr sideC
Definition RPDUtils.h:15
@ loc2
generic first and second local coordinate
Definition ParamDefs.h:35
@ loc1
Definition ParamDefs.h:34

Member Data Documentation

◆ m_g4UserEventInfo

AtlasG4EventUserInfo* RPCSensitiveDetector::m_g4UserEventInfo {nullptr}
private

Definition at line 114 of file RPCSensitiveDetector.h.

114{nullptr};

◆ m_hitCollectionName

std::string RPCSensitiveDetector::m_hitCollectionName
private

member data

Definition at line 112 of file RPCSensitiveDetector.h.

◆ m_muonHelper

const RpcHitIdHelper* RPCSensitiveDetector::m_muonHelper {nullptr}
private

Definition at line 115 of file RPCSensitiveDetector.h.

115{nullptr};

◆ m_myRPCHitColl

RPCSimHitCollection* RPCSensitiveDetector::m_myRPCHitColl {nullptr}
private

Definition at line 113 of file RPCSensitiveDetector.h.

113{nullptr};

The documentation for this class was generated from the following files: