]>
Commit | Line | Data |
---|---|---|
2817d3e2 | 1 | // $Id$ |
2 | // Category: event | |
3 | // | |
4 | // See the class description in the header file. | |
5 | ||
6 | #include "TG4StepManager.h" | |
da99be38 | 7 | #include "TG4GeometryServices.h" |
1d28d5b2 | 8 | #include "TG4ParticlesManager.h" |
2817d3e2 | 9 | #include "TG4PhysicsManager.h" |
10 | #include "TG4VSensitiveDetector.h" | |
8658d09c | 11 | #include "TG4Limits.h" |
2817d3e2 | 12 | #include "TG4Globals.h" |
1d28d5b2 | 13 | #include "TG4G3Units.h" |
2817d3e2 | 14 | |
2817d3e2 | 15 | #include <G4SteppingManager.hh> |
16 | #include <G4UserLimits.hh> | |
2817d3e2 | 17 | #include <G4UImanager.hh> |
18 | #include <G4AffineTransform.hh> | |
f2510df1 | 19 | #include <G4TransportationManager.hh> |
20 | #include <G4Navigator.hh> | |
62c4bd4b | 21 | #include <G4VProcess.hh> |
22 | #include <G4ProcessManager.hh> | |
23 | #include <G4ProcessVector.hh> | |
2817d3e2 | 24 | |
25 | #include <Randomize.hh> | |
2817d3e2 | 26 | #include <TLorentzVector.h> |
27 | ||
28 | TG4StepManager* TG4StepManager::fgInstance = 0; | |
29 | ||
30 | TG4StepManager::TG4StepManager() | |
f2510df1 | 31 | : fTrack(0), |
32 | fStep(0), | |
33 | fStepStatus(kNormalStep), | |
2817d3e2 | 34 | fSteppingManager(0) |
35 | { | |
36 | // | |
37 | if (fgInstance) { | |
38 | TG4Globals::Exception( | |
39 | "TG4StepManager: attempt to create two instances of singleton."); | |
40 | } | |
41 | ||
42 | fgInstance = this; | |
43 | } | |
44 | ||
45 | TG4StepManager::TG4StepManager(const TG4StepManager& right) { | |
46 | // | |
47 | TG4Globals::Exception( | |
48 | "Attempt to copy TG4StepManager singleton."); | |
49 | } | |
50 | ||
51 | TG4StepManager::~TG4StepManager() { | |
52 | // | |
53 | } | |
54 | ||
55 | // operators | |
56 | ||
57 | TG4StepManager& TG4StepManager::operator=(const TG4StepManager& right) | |
58 | { | |
59 | // check assignement to self | |
60 | if (this == &right) return *this; | |
61 | ||
62 | TG4Globals::Exception( | |
63 | "Attempt to assign TG4StepManager singleton."); | |
64 | ||
65 | return *this; | |
66 | } | |
67 | ||
68 | // private methods | |
69 | ||
f2510df1 | 70 | void TG4StepManager::CheckTrack() const |
71 | { | |
72 | // Gives exception in case the track is not defined. | |
73 | // --- | |
74 | ||
75 | if (!fTrack) | |
76 | TG4Globals::Exception("TG4StepManager: Track is not defined."); | |
77 | } | |
78 | ||
79 | ||
12139627 | 80 | void TG4StepManager::CheckStep(const G4String& method) const |
f2510df1 | 81 | { |
82 | // Gives exception in case the step is not defined. | |
83 | // --- | |
84 | ||
12139627 | 85 | if (!fStep) { |
86 | G4String text = "TG4StepManager::"; | |
87 | text = text + method + ": Step is not defined."; | |
88 | TG4Globals::Exception(text); | |
89 | } | |
f2510df1 | 90 | } |
91 | ||
92 | ||
93 | void TG4StepManager::CheckSteppingManager() const | |
94 | { | |
95 | // Gives exception in case the step is not defined. | |
96 | // --- | |
97 | ||
98 | if (!fSteppingManager) | |
99 | TG4Globals::Exception("TG4StepManager: Stepping manager is not defined."); | |
100 | } | |
101 | ||
102 | ||
2817d3e2 | 103 | void TG4StepManager::SetTLorentzVector(G4ThreeVector xyz, G4double t, |
104 | TLorentzVector& lv) const | |
105 | { | |
106 | // Fills TLorentzVector with G4ThreeVector and G4double. | |
107 | // --- | |
108 | ||
109 | lv[0] = xyz.x(); | |
110 | lv[1] = xyz.y(); | |
111 | lv[2] = xyz.z(); | |
112 | lv[3] = t; | |
113 | } | |
114 | ||
f2510df1 | 115 | G4VPhysicalVolume* TG4StepManager::GetCurrentOffPhysicalVolume(G4int off) const |
116 | { | |
117 | // Returns the physical volume of the off-th mother's | |
118 | // of the current volume. | |
119 | // --- | |
120 | ||
121 | G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume(); | |
122 | ||
12139627 | 123 | G4VPhysicalVolume* mother = physVolume; |
124 | ||
f2510df1 | 125 | Int_t level = off; |
126 | while (level > 0) { | |
127 | if (mother) mother = mother->GetMother(); | |
128 | level--; | |
129 | } | |
130 | ||
131 | if (!mother) { | |
132 | G4String text = "TG4StepManager::CurrentVolOff: \n"; | |
133 | text = text + " Volume "; | |
134 | text = text + physVolume->GetName(); | |
135 | text = text + " has not defined mother in the required level."; | |
136 | TG4Globals::Warning(text); | |
137 | } | |
138 | ||
139 | return mother; | |
140 | } | |
141 | ||
2817d3e2 | 142 | // public methods |
143 | ||
144 | void TG4StepManager::StopTrack() | |
145 | { | |
146 | // Stops the current track and skips to the next. | |
147 | // ?? do we want to invoke rest processes? | |
148 | // ?? do we want to stop secondaries too? | |
149 | // possible "stop" track status from G4: | |
150 | // fStopButAlive, // Invoke active rest physics processes and | |
151 | // // and kill the current track afterward | |
152 | // fStopAndKill, // Kill the current track | |
153 | // fKillTrackAndSecondaries, | |
154 | // // Kill the current track and also associated | |
155 | // // secondaries. | |
156 | // --- | |
157 | ||
12139627 | 158 | #ifdef TGEANT4_DEBUG |
f2510df1 | 159 | CheckTrack(); |
12139627 | 160 | #endif |
f2510df1 | 161 | |
162 | fTrack->SetTrackStatus(fStopAndKill); | |
163 | // fTrack->SetTrackStatus(fStopButAlive); | |
164 | // fTrack->SetTrackStatus(fKillTrackAndSecondaries); | |
2817d3e2 | 165 | } |
166 | ||
167 | void TG4StepManager::StopEvent() | |
168 | { | |
169 | // Aborts the current event processing. | |
170 | // --- | |
171 | ||
12139627 | 172 | #ifdef TGEANT4_DEBUG |
f2510df1 | 173 | CheckTrack(); |
12139627 | 174 | #endif |
f2510df1 | 175 | |
176 | fTrack->SetTrackStatus(fKillTrackAndSecondaries); | |
177 | //StopTrack(); // cannot be used as it keeps secondaries | |
178 | G4UImanager::GetUIpointer()->ApplyCommand("/event/abort"); | |
179 | G4UImanager::GetUIpointer()->ApplyCommand("/alStacking/clearStack"); | |
2817d3e2 | 180 | } |
181 | ||
2817d3e2 | 182 | void TG4StepManager::SetMaxStep(Float_t step) |
183 | { | |
184 | // Maximum step allowed in the current logical volume. | |
2817d3e2 | 185 | // --- |
186 | ||
f2510df1 | 187 | G4LogicalVolume* curLogVolume |
188 | = GetCurrentPhysicalVolume()->GetLogicalVolume(); | |
189 | G4UserLimits* userLimits | |
190 | = curLogVolume->GetUserLimits(); | |
8658d09c | 191 | |
f2510df1 | 192 | if (userLimits == 0) { |
8658d09c | 193 | // create new limits |
194 | userLimits = new TG4Limits(); | |
195 | ||
196 | // set limits to all logical volumes | |
197 | // corresponding to the current "G3" volume | |
198 | TG4GeometryServices* geometryServices = TG4GeometryServices::Instance(); | |
199 | G4int nofLV = geometryServices->SetUserLimits(userLimits, curLogVolume); | |
2817d3e2 | 200 | } |
8658d09c | 201 | |
202 | // set max step | |
1d28d5b2 | 203 | userLimits->SetMaxAllowedStep(step*TG4G3Units::Length()); |
2817d3e2 | 204 | } |
205 | ||
206 | void TG4StepManager::SetMaxNStep(Int_t maxNofSteps) | |
207 | { | |
208 | // Not yet implemented. | |
209 | // --- | |
210 | ||
211 | TG4Globals::Warning( | |
212 | "TG4StepManager::SetMaxNStep(..) is not yet implemented."); | |
213 | } | |
214 | ||
5025108a | 215 | void TG4StepManager::SetUserDecay(Int_t pdg) |
2817d3e2 | 216 | { |
217 | // Not yet implemented. | |
218 | // --- | |
219 | ||
220 | TG4Globals::Exception( | |
221 | "TG4StepManager::SetUserDecay(..) is not yet implemented."); | |
222 | } | |
223 | ||
f2510df1 | 224 | G4VPhysicalVolume* TG4StepManager::GetCurrentPhysicalVolume() const |
225 | { | |
226 | // Returns the current physical volume. | |
227 | // According to fStepStatus the volume from track vertex, | |
228 | // pre step point or post step point is returned. | |
229 | // --- | |
230 | ||
231 | G4VPhysicalVolume* physVolume; | |
232 | if (fStepStatus == kNormalStep) { | |
12139627 | 233 | |
234 | #ifdef TGEANT4_DEBUG | |
235 | CheckStep("GetCurrentPhysicalVolume"); | |
236 | #endif | |
237 | ||
f2510df1 | 238 | physVolume = fStep->GetPreStepPoint()->GetPhysicalVolume(); |
239 | } | |
240 | else if (fStepStatus == kBoundary) { | |
12139627 | 241 | |
242 | #ifdef TGEANT4_DEBUG | |
243 | CheckStep("GetCurrentPhysicalVolume"); | |
244 | #endif | |
245 | ||
f2510df1 | 246 | physVolume = fStep->GetPostStepPoint()->GetPhysicalVolume(); |
247 | } | |
248 | else { | |
12139627 | 249 | |
250 | #ifdef TGEANT4_DEBUG | |
f2510df1 | 251 | CheckTrack(); |
12139627 | 252 | #endif |
253 | ||
f2510df1 | 254 | G4ThreeVector position = fTrack->GetPosition(); |
255 | G4Navigator* navigator = | |
256 | G4TransportationManager::GetTransportationManager()-> | |
257 | GetNavigatorForTracking(); | |
258 | physVolume | |
259 | = navigator->LocateGlobalPointAndSetup(position); | |
260 | } | |
261 | ||
262 | return physVolume; | |
263 | } | |
264 | ||
2817d3e2 | 265 | Int_t TG4StepManager::CurrentVolID(Int_t& copyNo) const |
266 | { | |
267 | // Returns the current sensitive detector ID | |
268 | // and the copy number of the current physical volume. | |
269 | // --- | |
270 | ||
f2510df1 | 271 | G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume(); |
feb86ea5 | 272 | copyNo = physVolume->GetCopyNo() + 1; |
2817d3e2 | 273 | |
f2510df1 | 274 | // sensitive detector ID |
da99be38 | 275 | TG4GeometryServices* geometryServices = TG4GeometryServices::Instance(); |
276 | return geometryServices->GetVolumeID(physVolume->GetLogicalVolume()); | |
f2510df1 | 277 | } |
2817d3e2 | 278 | |
279 | Int_t TG4StepManager::CurrentVolOffID(Int_t off, Int_t& copyNo) const | |
280 | { | |
281 | // Returns the off-th mother's of the current volume | |
282 | // the sensitive detector ID and the copy number. | |
283 | // --- | |
284 | ||
285 | if (off == 0) return CurrentVolID(copyNo); | |
286 | ||
f2510df1 | 287 | G4VPhysicalVolume* mother = GetCurrentOffPhysicalVolume(off); |
288 | ||
289 | if (mother) { | |
feb86ea5 | 290 | copyNo = mother->GetCopyNo() + 1; |
f2510df1 | 291 | |
292 | // sensitive detector ID | |
da99be38 | 293 | TG4GeometryServices* geometryServices = TG4GeometryServices::Instance(); |
294 | return geometryServices->GetVolumeID(mother->GetLogicalVolume()); | |
2817d3e2 | 295 | } |
f2510df1 | 296 | else { |
297 | copyNo = 0; | |
298 | return 0; | |
299 | } | |
2817d3e2 | 300 | } |
301 | ||
302 | const char* TG4StepManager::CurrentVolName() const | |
303 | { | |
304 | // Returns the current physical volume name. | |
305 | // --- | |
306 | ||
f2510df1 | 307 | return GetCurrentPhysicalVolume()->GetName(); |
2817d3e2 | 308 | } |
309 | ||
310 | const char* TG4StepManager::CurrentVolOffName(Int_t off) const | |
311 | { | |
312 | // Returns the off-th mother's physical volume name. | |
313 | // --- | |
314 | ||
f2510df1 | 315 | if (off == 0) return CurrentVolName(); |
2817d3e2 | 316 | |
f2510df1 | 317 | G4VPhysicalVolume* mother = GetCurrentOffPhysicalVolume(off); |
318 | ||
319 | if (mother) | |
320 | return mother->GetName(); | |
321 | else | |
2817d3e2 | 322 | return 0; |
2817d3e2 | 323 | } |
324 | ||
325 | Int_t TG4StepManager::CurrentMaterial(Float_t &a, Float_t &z, Float_t &dens, | |
326 | Float_t &radl, Float_t &absl) const | |
327 | { | |
328 | // Returns the parameters of the current material during transport | |
329 | // the return value is the number of elements in the mixture. | |
330 | // --- | |
331 | ||
f2510df1 | 332 | G4VPhysicalVolume* physVolume = GetCurrentPhysicalVolume(); |
333 | ||
334 | G4Material* material | |
335 | = physVolume->GetLogicalVolume()->GetMaterial(); | |
336 | ||
12139627 | 337 | if (material) { |
f2510df1 | 338 | G4int nofElements = material->GetNumberOfElements(); |
da99be38 | 339 | TG4GeometryServices* geometryServices = TG4GeometryServices::Instance(); |
340 | a = geometryServices->GetEffA(material); | |
341 | z = geometryServices->GetEffZ(material); | |
2817d3e2 | 342 | |
f2510df1 | 343 | // density |
344 | dens = material->GetDensity(); | |
1d28d5b2 | 345 | dens /= TG4G3Units::MassDensity(); |
2817d3e2 | 346 | |
f2510df1 | 347 | // radiation length |
348 | radl = material->GetRadlen(); | |
1d28d5b2 | 349 | radl /= TG4G3Units::Length(); |
2817d3e2 | 350 | |
f2510df1 | 351 | absl = 0.; // this parameter is not defined in Geant4 |
352 | return nofElements; | |
353 | } | |
354 | else { | |
355 | TG4Globals::Exception( | |
356 | "TG4StepManager::CurrentMaterial(..): material is not defined."); | |
2817d3e2 | 357 | return 0; |
358 | } | |
359 | } | |
360 | ||
361 | void TG4StepManager::Gmtod(Float_t* xm, Float_t* xd, Int_t iflag) | |
362 | { | |
363 | // Transforms a position from the world reference frame | |
364 | // to the current volume reference frame. | |
365 | // | |
366 | // Geant3 desription: | |
367 | // ================== | |
368 | // Computes coordinates XD (in DRS) | |
369 | // from known coordinates XM in MRS | |
370 | // The local reference system can be initialized by | |
371 | // - the tracking routines and GMTOD used in GUSTEP | |
372 | // - a call to GMEDIA(XM,NUMED) | |
373 | // - a call to GLVOLU(NLEVEL,NAMES,NUMBER,IER) | |
374 | // (inverse routine is GDTOM) | |
375 | // | |
376 | // If IFLAG=1 convert coordinates | |
377 | // IFLAG=2 convert direction cosinus | |
378 | // | |
379 | // --- | |
380 | ||
12139627 | 381 | G4AffineTransform affineTransform; |
f2510df1 | 382 | |
12139627 | 383 | if (fStepStatus == kVertex) { |
384 | G4Navigator* navigator = | |
385 | G4TransportationManager::GetTransportationManager()-> | |
386 | GetNavigatorForTracking(); | |
387 | ||
388 | affineTransform = navigator->GetGlobalToLocalTransform(); | |
389 | } | |
390 | else { | |
f2510df1 | 391 | |
12139627 | 392 | #ifdef TGEANT4_DEBUG |
393 | CheckStep("Gmtod"); | |
394 | #endif | |
395 | ||
396 | affineTransform | |
397 | = fStep->GetPreStepPoint()->GetTouchable()->GetHistory() | |
398 | ->GetTopTransform(); | |
399 | } | |
400 | ||
401 | G4ThreeVector theGlobalPoint(xm[0],xm[1],xm[2]); | |
f2510df1 | 402 | G4ThreeVector theLocalPoint; |
403 | if(iflag == 1) | |
404 | theLocalPoint = affineTransform.TransformPoint(theGlobalPoint); | |
405 | else if ( iflag == 2) | |
406 | theLocalPoint = affineTransform.TransformAxis(theGlobalPoint); | |
407 | else | |
408 | TG4Globals::Exception( | |
409 | "TG4StepManager::Gmtod(..,iflag): iflag is not in 1..2"); | |
2817d3e2 | 410 | |
f2510df1 | 411 | xd[0] = theLocalPoint.x(); |
412 | xd[1] = theLocalPoint.y(); | |
413 | xd[2] = theLocalPoint.z(); | |
2817d3e2 | 414 | } |
415 | ||
416 | void TG4StepManager::Gdtom(Float_t* xd, Float_t* xm, Int_t iflag) | |
417 | { | |
418 | // Transforms a position from the current volume reference frame | |
419 | // to the world reference frame. | |
420 | // | |
421 | // Geant3 desription: | |
422 | // ================== | |
423 | // Computes coordinates XM (Master Reference System | |
424 | // knowing the coordinates XD (Detector Ref System) | |
425 | // The local reference system can be initialized by | |
426 | // - the tracking routines and GDTOM used in GUSTEP | |
427 | // - a call to GSCMED(NLEVEL,NAMES,NUMBER) | |
428 | // (inverse routine is GMTOD) | |
429 | // | |
430 | // If IFLAG=1 convert coordinates | |
431 | // IFLAG=2 convert direction cosinus | |
432 | // | |
433 | // --- | |
434 | ||
12139627 | 435 | G4AffineTransform affineTransform; |
f2510df1 | 436 | |
12139627 | 437 | if (fStepStatus == kVertex) { |
438 | G4Navigator* navigator = | |
439 | G4TransportationManager::GetTransportationManager()-> | |
440 | GetNavigatorForTracking(); | |
441 | ||
442 | affineTransform = navigator->GetLocalToGlobalTransform(); | |
443 | } | |
444 | else { | |
2817d3e2 | 445 | |
12139627 | 446 | #ifdef TGEANT4_DEBUG |
447 | CheckStep("Gdtom"); | |
448 | #endif | |
449 | ||
450 | // check this | |
451 | ||
452 | affineTransform | |
453 | = fStep->GetPreStepPoint()->GetTouchable()->GetHistory() | |
2817d3e2 | 454 | ->GetTopTransform().Inverse(); |
12139627 | 455 | } |
2817d3e2 | 456 | |
12139627 | 457 | G4ThreeVector theLocalPoint(xd[0],xd[1],xd[2]); |
f2510df1 | 458 | G4ThreeVector theGlobalPoint; |
459 | if(iflag == 1) | |
460 | theGlobalPoint = affineTransform.TransformPoint(theLocalPoint); | |
461 | else if( iflag == 2) | |
462 | theGlobalPoint = affineTransform.TransformAxis(theLocalPoint); | |
463 | else | |
464 | TG4Globals::Warning( | |
465 | "TG4StepManager::Gdtom(...,iflag): iflag is not in 1..2"); | |
466 | ||
467 | xm[0] = theGlobalPoint.x(); | |
468 | xm[1] = theGlobalPoint.y(); | |
469 | xm[2] = theGlobalPoint.z(); | |
2817d3e2 | 470 | } |
471 | ||
472 | Float_t TG4StepManager::MaxStep() const | |
473 | { | |
474 | // Returns maximum step allowed in the current logical volume | |
475 | // by User Limits. | |
476 | // --- | |
477 | ||
12139627 | 478 | G4LogicalVolume* curLogVolume |
f2510df1 | 479 | = GetCurrentPhysicalVolume()->GetLogicalVolume(); |
12139627 | 480 | |
481 | // check this | |
f2510df1 | 482 | G4UserLimits* userLimits |
483 | = curLogVolume->GetUserLimits(); | |
484 | ||
485 | G4double maxStep; | |
486 | if (userLimits == 0) { | |
487 | G4String text = "User Limits are not defined for log volume "; | |
488 | text = text + curLogVolume->GetName(); | |
489 | TG4Globals::Warning(text); | |
2acb5b6d | 490 | return FLT_MAX; |
2817d3e2 | 491 | } |
f2510df1 | 492 | else { |
493 | const G4Track& trackRef = *(fTrack); | |
494 | maxStep = userLimits->GetMaxAllowedStep(trackRef); | |
1d28d5b2 | 495 | maxStep /= TG4G3Units::Length(); |
f2510df1 | 496 | return maxStep; |
497 | } | |
2817d3e2 | 498 | } |
499 | ||
500 | Int_t TG4StepManager::GetMaxNStep() const | |
501 | { | |
502 | // Not yet implemented. | |
503 | // --- | |
504 | ||
505 | TG4Globals::Warning( | |
506 | "Method GetMaxNStep is not yet implemented in TG4StepManager."); | |
507 | return 0; | |
508 | } | |
509 | ||
510 | void TG4StepManager::TrackPosition(TLorentzVector& position) const | |
511 | { | |
512 | // Current particle position (in the world reference frame) | |
513 | // and the local time since the current track is created | |
514 | // (position of the PostStepPoint). | |
515 | // --- | |
516 | ||
12139627 | 517 | #ifdef TGEANT4_DEBUG |
f2510df1 | 518 | CheckTrack(); |
12139627 | 519 | #endif |
520 | ||
f2510df1 | 521 | // get position |
522 | // check if this is == to PostStepPoint position !! | |
523 | G4ThreeVector positionVector = fTrack->GetPosition(); | |
1d28d5b2 | 524 | positionVector *= 1./(TG4G3Units::Length()); |
2817d3e2 | 525 | |
f2510df1 | 526 | // local time |
527 | G4double time = fTrack->GetLocalTime(); | |
1d28d5b2 | 528 | time /= TG4G3Units::Time(); |
f2510df1 | 529 | |
530 | SetTLorentzVector(positionVector, time, position); | |
2817d3e2 | 531 | } |
532 | ||
533 | Int_t TG4StepManager::GetMedium() const | |
534 | { | |
535 | // Returns the second index of the current material (corresponding to | |
536 | // G3 tracking medium index). | |
537 | // --- | |
538 | ||
f2510df1 | 539 | // current material |
540 | G4Material* curMaterial | |
541 | = GetCurrentPhysicalVolume()->GetLogicalVolume()->GetMaterial(); | |
2817d3e2 | 542 | |
f2510df1 | 543 | // medium index |
da99be38 | 544 | TG4GeometryServices* geometryServices = TG4GeometryServices::Instance(); |
545 | return geometryServices->GetMediumId(curMaterial); | |
2817d3e2 | 546 | } |
547 | ||
548 | void TG4StepManager::TrackMomentum(TLorentzVector& momentum) const | |
549 | { | |
550 | // Current particle "momentum" (px, py, pz, Etot). | |
551 | // --- | |
552 | ||
12139627 | 553 | #ifdef TGEANT4_DEBUG |
f2510df1 | 554 | CheckTrack(); |
12139627 | 555 | #endif |
2817d3e2 | 556 | |
f2510df1 | 557 | G4ThreeVector momentumVector = fTrack->GetMomentum(); |
1d28d5b2 | 558 | momentumVector *= 1./(TG4G3Units::Energy()); |
2817d3e2 | 559 | |
f2510df1 | 560 | G4double energy = fTrack->GetDynamicParticle()->GetTotalEnergy(); |
1d28d5b2 | 561 | energy /= TG4G3Units::Energy(); |
f2510df1 | 562 | |
563 | SetTLorentzVector(momentumVector, energy, momentum); | |
2817d3e2 | 564 | } |
565 | ||
566 | void TG4StepManager::TrackVertexPosition(TLorentzVector& position) const | |
567 | { | |
568 | // The vertex particle position (in the world reference frame) | |
569 | // and the local time since the current track is created. | |
570 | // --- | |
571 | ||
12139627 | 572 | #ifdef TGEANT4_DEBUG |
f2510df1 | 573 | CheckTrack(); |
12139627 | 574 | #endif |
f2510df1 | 575 | |
576 | // position | |
577 | G4ThreeVector positionVector = fTrack->GetVertexPosition(); | |
1d28d5b2 | 578 | positionVector *= 1./(TG4G3Units::Length()); |
2817d3e2 | 579 | |
f2510df1 | 580 | // local time |
581 | // to be checked | |
582 | G4double time = fTrack->GetLocalTime(); | |
1d28d5b2 | 583 | time /= TG4G3Units::Time(); |
2817d3e2 | 584 | |
f2510df1 | 585 | SetTLorentzVector(positionVector, time, position); |
2817d3e2 | 586 | } |
587 | ||
588 | void TG4StepManager::TrackVertexMomentum(TLorentzVector& momentum) const | |
589 | { | |
590 | // The vertex particle "momentum" (px, py, pz, Ekin) | |
591 | // to do: change Ekin -> Etot | |
592 | // --- | |
593 | ||
12139627 | 594 | #ifdef TGEANT4_DEBUG |
f2510df1 | 595 | CheckTrack(); |
12139627 | 596 | #endif |
597 | ||
f2510df1 | 598 | G4ThreeVector momentumVector = fTrack->GetVertexMomentumDirection(); |
1d28d5b2 | 599 | momentumVector *= 1./(TG4G3Units::Energy()); |
2817d3e2 | 600 | |
f2510df1 | 601 | G4double energy = fTrack->GetVertexKineticEnergy(); |
1d28d5b2 | 602 | energy /= TG4G3Units::Energy(); |
2817d3e2 | 603 | |
f2510df1 | 604 | SetTLorentzVector(momentumVector, energy, momentum); |
2817d3e2 | 605 | } |
606 | ||
607 | Float_t TG4StepManager::TrackStep() const | |
608 | { | |
609 | // Returns the current step length. | |
610 | // --- | |
611 | ||
f2510df1 | 612 | G4double length; |
613 | if (fStepStatus == kNormalStep) { | |
12139627 | 614 | |
615 | #ifdef TGEANT4_DEBUG | |
616 | CheckStep("TrackStep"); | |
617 | #endif | |
618 | ||
f2510df1 | 619 | length = fStep->GetStepLength(); |
1d28d5b2 | 620 | length /= TG4G3Units::Length(); |
f2510df1 | 621 | } |
622 | else | |
623 | length = 0; | |
624 | ||
625 | return length; | |
2817d3e2 | 626 | } |
627 | ||
628 | Float_t TG4StepManager::TrackLength() const | |
629 | { | |
630 | // Returns the length of the current track from its origin. | |
631 | // --- | |
632 | ||
12139627 | 633 | #ifdef TGEANT4_DEBUG |
f2510df1 | 634 | CheckTrack(); |
12139627 | 635 | #endif |
f2510df1 | 636 | |
637 | G4double length = fTrack->GetTrackLength(); | |
1d28d5b2 | 638 | length /= TG4G3Units::Length(); |
f2510df1 | 639 | return length; |
2817d3e2 | 640 | } |
641 | ||
642 | Float_t TG4StepManager::TrackTime() const | |
643 | { | |
644 | // Returns the local time since the current track is created. | |
645 | // Comment: | |
646 | // in Geant4: there is also defined proper time as | |
647 | // the proper time of the dynamical particle of the current track. | |
648 | // --- | |
649 | ||
12139627 | 650 | #ifdef TGEANT4_DEBUG |
f2510df1 | 651 | CheckTrack(); |
12139627 | 652 | #endif |
f2510df1 | 653 | |
654 | G4double time = fTrack->GetLocalTime(); | |
1d28d5b2 | 655 | time /= TG4G3Units::Time(); |
f2510df1 | 656 | return time; |
2817d3e2 | 657 | } |
658 | ||
659 | Float_t TG4StepManager::Edep() const | |
660 | { | |
661 | // Returns total energy deposit in this step. | |
662 | // --- | |
663 | ||
f2510df1 | 664 | G4double energyDeposit; |
665 | if (fStepStatus == kNormalStep) { | |
12139627 | 666 | |
667 | #ifdef TGEANT4_DEBUG | |
668 | CheckStep("Edep"); | |
669 | #endif | |
670 | ||
f2510df1 | 671 | energyDeposit = fStep->GetTotalEnergyDeposit(); |
1d28d5b2 | 672 | energyDeposit /= TG4G3Units::Energy(); |
2817d3e2 | 673 | } |
f2510df1 | 674 | else |
675 | energyDeposit = 0; | |
676 | ||
677 | return energyDeposit; | |
2817d3e2 | 678 | } |
679 | ||
680 | Int_t TG4StepManager::TrackPid() const | |
681 | { | |
682 | // Returns the current particle PDG encoding. | |
683 | // --- | |
684 | ||
12139627 | 685 | #ifdef TGEANT4_DEBUG |
f2510df1 | 686 | CheckTrack(); |
12139627 | 687 | #endif |
2817d3e2 | 688 | |
f2510df1 | 689 | G4ParticleDefinition* particle |
690 | = fTrack->GetDynamicParticle()->GetDefinition(); | |
691 | ||
1d28d5b2 | 692 | // ask TG4ParticlesManager to get PDG encoding |
f2510df1 | 693 | // (in order to get PDG from extended TDatabasePDG |
694 | // in case the standard PDG code is not defined) | |
1d28d5b2 | 695 | G4int pdgEncoding |
696 | = TG4ParticlesManager::Instance()->GetPDGEncodingFast(particle); | |
f2510df1 | 697 | |
698 | return pdgEncoding; | |
2817d3e2 | 699 | } |
700 | ||
701 | Float_t TG4StepManager::TrackCharge() const | |
702 | { | |
703 | // Returns the current particle charge. | |
704 | // --- | |
705 | ||
12139627 | 706 | #ifdef TGEANT4_DEBUG |
f2510df1 | 707 | CheckTrack(); |
12139627 | 708 | #endif |
709 | ||
f2510df1 | 710 | G4double charge |
711 | = fTrack->GetDynamicParticle()->GetDefinition() | |
712 | ->GetPDGCharge(); | |
1d28d5b2 | 713 | charge /= TG4G3Units::Charge(); |
f2510df1 | 714 | return charge; |
2817d3e2 | 715 | } |
716 | ||
717 | Float_t TG4StepManager::TrackMass() const | |
718 | { | |
719 | // Returns current particle rest mass. | |
720 | // --- | |
721 | ||
12139627 | 722 | #ifdef TGEANT4_DEBUG |
f2510df1 | 723 | CheckTrack(); |
12139627 | 724 | #endif |
f2510df1 | 725 | |
726 | G4double mass | |
727 | = fTrack->GetDynamicParticle()->GetDefinition() | |
728 | ->GetPDGMass(); | |
1d28d5b2 | 729 | mass /= TG4G3Units::Mass(); |
f2510df1 | 730 | return mass; |
2817d3e2 | 731 | } |
732 | ||
733 | Float_t TG4StepManager::Etot() const | |
734 | { | |
735 | // Returns total energy of the current particle. | |
736 | // --- | |
737 | ||
12139627 | 738 | #ifdef TGEANT4_DEBUG |
f2510df1 | 739 | CheckTrack(); |
12139627 | 740 | #endif |
f2510df1 | 741 | |
742 | G4double energy | |
743 | = fTrack->GetDynamicParticle()->GetTotalEnergy(); | |
1d28d5b2 | 744 | energy /= TG4G3Units::Energy(); |
f2510df1 | 745 | return energy; |
2817d3e2 | 746 | } |
747 | ||
748 | Bool_t TG4StepManager::IsTrackInside() const | |
749 | { | |
750 | // Returns true if particle does not cross geometrical boundary | |
f2510df1 | 751 | // and is not in vertex. |
2817d3e2 | 752 | // --- |
753 | ||
f2510df1 | 754 | if (fStepStatus == kNormalStep && !(IsTrackExiting()) ) { |
755 | // track is always inside during a normal step | |
756 | return true; | |
757 | } | |
758 | ||
759 | return false; | |
2817d3e2 | 760 | } |
761 | ||
762 | Bool_t TG4StepManager::IsTrackEntering() const | |
763 | { | |
764 | // Returns true if particle cross a geometrical boundary | |
f2510df1 | 765 | // or is in vertex. |
2817d3e2 | 766 | // --- |
767 | ||
f2510df1 | 768 | if (fStepStatus != kNormalStep) { |
769 | // track is entering during a vertex or boundary step | |
770 | return true; | |
771 | } | |
772 | ||
773 | return false; | |
2817d3e2 | 774 | } |
775 | ||
776 | Bool_t TG4StepManager::IsTrackExiting() const | |
777 | { | |
f2510df1 | 778 | // Returns true if particle cross a geometrical boundary. |
2817d3e2 | 779 | // --- |
780 | ||
f2510df1 | 781 | if (fStepStatus == kNormalStep) { |
12139627 | 782 | |
783 | #ifdef TGEANT4_DEBUG | |
784 | CheckStep("IsTrackExiting"); | |
785 | #endif | |
786 | ||
f2510df1 | 787 | if (fStep->GetPostStepPoint()->GetStepStatus() == fGeomBoundary) |
788 | return true; | |
2817d3e2 | 789 | } |
f2510df1 | 790 | |
791 | return false; | |
2817d3e2 | 792 | } |
793 | ||
794 | Bool_t TG4StepManager::IsTrackOut() const | |
795 | { | |
68f491b7 | 796 | // Returns true if particle cross the world boundary |
2817d3e2 | 797 | // at post-step point. |
798 | // --- | |
799 | ||
f2510df1 | 800 | if (fStepStatus == kVertex) return false; |
801 | ||
12139627 | 802 | #ifdef TGEANT4_DEBUG |
803 | CheckStep("IsTrackCut"); | |
804 | #endif | |
f2510df1 | 805 | |
806 | // check | |
807 | G4StepStatus status | |
808 | = fStep->GetPostStepPoint()->GetStepStatus(); | |
809 | if (status == fWorldBoundary) | |
810 | return true; | |
811 | else | |
2817d3e2 | 812 | return false; |
2817d3e2 | 813 | } |
814 | ||
815 | Bool_t TG4StepManager::IsTrackStop() const | |
816 | { | |
817 | // Returns true if particle has stopped | |
818 | // or has been killed, suspended or postponed to next event. | |
819 | // | |
820 | // Possible track status from G4: | |
821 | // fAlive, // Continue the tracking | |
822 | // fStopButAlive, // Invoke active rest physics processes and | |
823 | // // and kill the current track afterward | |
824 | // fStopAndKill, // Kill the current track | |
825 | // fKillTrackAndSecondaries, | |
826 | // // Kill the current track and also associated | |
827 | // // secondaries. | |
828 | // fSuspend, // Suspend the current track | |
829 | // fPostponeToNextEvent | |
830 | // // Postpones the tracking of thecurrent track | |
831 | // // to the next event. | |
832 | // --- | |
833 | ||
12139627 | 834 | #ifdef TGEANT4_DEBUG |
f2510df1 | 835 | CheckTrack(); |
12139627 | 836 | #endif |
f2510df1 | 837 | |
838 | // check | |
839 | G4TrackStatus status | |
840 | = fTrack->GetTrackStatus(); | |
841 | if ((status == fStopAndKill) || | |
842 | (status == fKillTrackAndSecondaries) || | |
843 | (status == fSuspend) || | |
844 | (status == fPostponeToNextEvent)) { | |
845 | return true; | |
2817d3e2 | 846 | } |
f2510df1 | 847 | else |
848 | return false; | |
2817d3e2 | 849 | } |
850 | ||
851 | Bool_t TG4StepManager::IsTrackDisappeared() const | |
852 | { | |
853 | // Returns true if particle has disappeared | |
854 | // (due to any physical process) | |
855 | // or has been killed, suspended or postponed to next event. | |
856 | // --- | |
857 | ||
12139627 | 858 | #ifdef TGEANT4_DEBUG |
f2510df1 | 859 | CheckTrack(); |
12139627 | 860 | #endif |
f2510df1 | 861 | |
862 | // check | |
863 | G4TrackStatus status | |
864 | = fTrack->GetTrackStatus(); | |
865 | if ((status == fStopButAlive) || | |
866 | (status == fKillTrackAndSecondaries) || | |
867 | (status == fSuspend) || | |
868 | (status == fPostponeToNextEvent)) { | |
869 | return true; | |
2817d3e2 | 870 | } |
f2510df1 | 871 | else |
872 | return false; | |
2817d3e2 | 873 | } |
874 | ||
875 | Bool_t TG4StepManager::IsTrackAlive() const | |
876 | { | |
877 | // Returns true if particle continues tracking. | |
878 | // --- | |
879 | ||
12139627 | 880 | #ifdef TGEANT4_DEBUG |
f2510df1 | 881 | CheckTrack(); |
12139627 | 882 | #endif |
f2510df1 | 883 | |
884 | G4TrackStatus status | |
885 | = fTrack->GetTrackStatus(); | |
886 | if (status == fAlive) | |
887 | return true; | |
888 | else | |
889 | return false; | |
2817d3e2 | 890 | } |
891 | ||
892 | Bool_t TG4StepManager::IsNewTrack() const | |
893 | { | |
894 | // Returns true when track performs the first step. | |
895 | // --- | |
896 | ||
f2510df1 | 897 | if (fStepStatus == kVertex) |
2817d3e2 | 898 | return true; |
899 | else | |
900 | return false; | |
901 | } | |
902 | ||
903 | Int_t TG4StepManager::NSecondaries() const | |
904 | { | |
905 | // Returns the number of secondary particles generated | |
906 | // in the current step. | |
907 | // --- | |
908 | ||
12139627 | 909 | #ifdef TGEANT4_DEBUG |
f2510df1 | 910 | CheckSteppingManager(); |
12139627 | 911 | #endif |
2817d3e2 | 912 | |
f2510df1 | 913 | G4int nofSecondaries = 0; |
914 | nofSecondaries += fSteppingManager->GetfN2ndariesAtRestDoIt(); | |
915 | nofSecondaries += fSteppingManager->GetfN2ndariesAlongStepDoIt(); | |
916 | nofSecondaries += fSteppingManager->GetfN2ndariesPostStepDoIt(); | |
917 | ||
918 | return nofSecondaries; | |
2817d3e2 | 919 | } |
920 | ||
921 | void TG4StepManager::GetSecondary(Int_t index, Int_t& particleId, | |
922 | TLorentzVector& position, TLorentzVector& momentum) | |
923 | { | |
924 | // Fills the parameters (particle encoding, position, momentum) | |
925 | // of the generated secondary particle which is specified by index. | |
926 | // !! Check if indexing of secondaries is same !! | |
927 | // --- | |
928 | ||
12139627 | 929 | #ifdef TGEANT4_DEBUG |
f2510df1 | 930 | CheckSteppingManager(); |
12139627 | 931 | #endif |
f2510df1 | 932 | |
933 | G4int nofSecondaries = NSecondaries(); | |
934 | G4TrackVector* secondaryTracks = fSteppingManager->GetSecondary(); | |
935 | ||
936 | if (secondaryTracks){ | |
937 | if (index < nofSecondaries) { | |
938 | ||
939 | // the index of the first secondary of this step | |
940 | G4int startIndex | |
941 | = secondaryTracks->entries() - nofSecondaries; | |
942 | // (the secondaryTracks vector contains secondaries | |
943 | // produced by the track at previous steps, too) | |
944 | G4Track* track | |
945 | = (*secondaryTracks)[startIndex + index]; | |
2817d3e2 | 946 | |
f2510df1 | 947 | // particle encoding |
948 | particleId | |
949 | = track->GetDynamicParticle()->GetDefinition()->GetPDGEncoding(); | |
2817d3e2 | 950 | |
f2510df1 | 951 | // position & time |
952 | G4ThreeVector positionVector = track->GetPosition(); | |
1d28d5b2 | 953 | positionVector *= 1./(TG4G3Units::Length()); |
f2510df1 | 954 | G4double time = track->GetLocalTime(); |
1d28d5b2 | 955 | time /= TG4G3Units::Time(); |
f2510df1 | 956 | SetTLorentzVector(positionVector, time, position); |
957 | ||
958 | // momentum & energy | |
959 | G4ThreeVector momentumVector = track->GetMomentum(); | |
960 | G4double energy = track->GetDynamicParticle()->GetTotalEnergy(); | |
1d28d5b2 | 961 | energy /= TG4G3Units::Energy(); |
f2510df1 | 962 | SetTLorentzVector(momentumVector, energy, momentum); |
2817d3e2 | 963 | } |
964 | else { | |
965 | TG4Globals::Exception( | |
f2510df1 | 966 | "TG4StepManager::GetSecondary(): wrong secondary track index."); |
2817d3e2 | 967 | } |
968 | } | |
f2510df1 | 969 | else { |
970 | TG4Globals::Exception( | |
971 | "TG4StepManager::GetSecondary(): secondary tracks vector is empty"); | |
2817d3e2 | 972 | } |
973 | } | |
974 | ||
da99be38 | 975 | AliMCProcess TG4StepManager::ProdProcess(Int_t isec) const |
2817d3e2 | 976 | { |
da99be38 | 977 | // The process that has produced the secondary particles specified |
978 | // with isec index in the current step. | |
2817d3e2 | 979 | // --- |
980 | ||
12139627 | 981 | G4int nofSecondaries = NSecondaries(); |
7bfb4926 | 982 | if (fStepStatus == kVertex || !nofSecondaries) return kPNoProcess; |
983 | ||
12139627 | 984 | #ifdef TGEANT4_DEBUG |
985 | CheckStep("ProdProcess"); | |
986 | #endif | |
2817d3e2 | 987 | |
7bfb4926 | 988 | G4TrackVector* secondaryTracks = fSteppingManager->GetSecondary(); |
989 | ||
990 | // should never happen | |
991 | if (!secondaryTracks) { | |
992 | TG4Globals::Exception( | |
993 | "TG4StepManager::ProdProcess(): secondary tracks vector is empty."); | |
da99be38 | 994 | |
995 | return kPNoProcess; | |
7bfb4926 | 996 | } |
f2510df1 | 997 | |
da99be38 | 998 | if (isec < nofSecondaries) { |
7bfb4926 | 999 | |
da99be38 | 1000 | // the index of the first secondary of this step |
1001 | G4int startIndex | |
1002 | = secondaryTracks->entries() - nofSecondaries; | |
62c4bd4b | 1003 | // the secondaryTracks vector contains secondaries |
1004 | // produced by the track at previous steps, too | |
da99be38 | 1005 | |
1006 | // the secondary track with specified isec index | |
1007 | G4Track* track = (*secondaryTracks)[startIndex + isec]; | |
1008 | ||
62c4bd4b | 1009 | const G4VProcess* kpProcess = track->GetCreatorProcess(); |
da99be38 | 1010 | |
1d28d5b2 | 1011 | AliMCProcess mcProcess |
1012 | = TG4PhysicsManager::Instance()->GetMCProcess(kpProcess); | |
7bfb4926 | 1013 | |
da99be38 | 1014 | // distinguish kPDeltaRay from kPEnergyLoss |
1015 | if (mcProcess == kPEnergyLoss) mcProcess = kPDeltaRay; | |
7bfb4926 | 1016 | |
da99be38 | 1017 | return mcProcess; |
1018 | } | |
1019 | else { | |
1020 | TG4Globals::Exception( | |
1021 | "TG4StepManager::GetSecondary(): wrong secondary track index."); | |
1022 | ||
1023 | return kPNoProcess; | |
1024 | } | |
1025 | } | |
1026 | ||
1027 | ||
1028 | Int_t TG4StepManager::StepProcesses(TArrayI &proc) const | |
1029 | { | |
62c4bd4b | 1030 | // Fills the array of processes that were active in the current step |
1031 | // and returns the number of them. | |
1032 | // TBD: Distinguish between kPDeltaRay and kPEnergyLoss | |
da99be38 | 1033 | // --- |
1034 | ||
62c4bd4b | 1035 | if (fStepStatus == kVertex) { |
1036 | G4cout << "kVertex" << G4endl; | |
1037 | G4int nofProcesses = 1; | |
1038 | proc.Set(nofProcesses); | |
1039 | proc[0] = kPNull; | |
1040 | return nofProcesses; | |
1041 | } | |
1042 | ||
1043 | #ifdef TGEANT4_DEBUG | |
1044 | CheckSteppingManager(); | |
1045 | CheckStep("StepProcesses"); | |
1046 | #endif | |
1047 | ||
1048 | // along step processes | |
1049 | G4ProcessManager* processManager | |
1050 | = fStep->GetTrack()->GetDefinition()->GetProcessManager(); | |
1051 | G4ProcessVector* alongStepProcessVector | |
1052 | = processManager->GetAlongStepProcessVector(); | |
1053 | G4int nofProcesses = alongStepProcessVector->entries(); | |
1054 | ||
1055 | // process defined step | |
1056 | const G4VProcess* kpLastProcess | |
1057 | = fStep->GetPostStepPoint()->GetProcessDefinedStep(); | |
1058 | ||
1059 | // fill the array of processes | |
1060 | proc.Set(nofProcesses); | |
1061 | TG4PhysicsManager* physicsManager = TG4PhysicsManager::Instance(); | |
1062 | G4int i; | |
1063 | for (i=0; i<nofProcesses-1; i++) { | |
1064 | G4VProcess* g4Process = (*alongStepProcessVector)[i]; | |
1065 | // do not fill transportation along step process | |
1066 | if (g4Process->GetProcessName() != "Transportation") { | |
1067 | physicsManager->GetMCProcess(g4Process); | |
1068 | proc[i] = physicsManager->GetMCProcess(g4Process); | |
1069 | } | |
1070 | } | |
1071 | proc[nofProcesses-1] = physicsManager->GetMCProcess(kpLastProcess); | |
1072 | ||
1073 | return nofProcesses; | |
2817d3e2 | 1074 | } |