4 // Author: I. Hrivnacova
6 // Class TG4VSpecialCuts
7 // ---------------------
8 // See the class description in the header file.
10 #include "TG4VSpecialCuts.h"
11 #include "TG4G3CutVector.h"
12 #include "TG4GeometryServices.h"
13 #include "TG4Limits.h"
15 #include <G4UserLimits.hh>
16 #include <G4EnergyLossTables.hh>
18 //_____________________________________________________________________________
19 TG4VSpecialCuts::TG4VSpecialCuts(const G4String& processName)
20 : G4VProcess(processName) {
24 //_____________________________________________________________________________
25 TG4VSpecialCuts::TG4VSpecialCuts() {
29 //_____________________________________________________________________________
30 TG4VSpecialCuts::TG4VSpecialCuts(const TG4VSpecialCuts& right) {
32 TG4Globals::Exception(
33 "TG4VSpecialCuts is protected from copying.");
36 //_____________________________________________________________________________
37 TG4VSpecialCuts::~TG4VSpecialCuts() {
43 //_____________________________________________________________________________
44 TG4VSpecialCuts& TG4VSpecialCuts::operator=(const TG4VSpecialCuts& right)
46 // check assignement to self
47 if (this == &right) return *this;
49 TG4Globals::Exception(
50 "TG4VSpecialCuts is protected from assigning.");
57 //_____________________________________________________________________________
58 G4double TG4VSpecialCuts::PostStepGetPhysicalInteractionLength(
59 const G4Track& track, G4double previousStepSize,
60 G4ForceCondition* condition)
62 // Returns the Step-size (actual length) which is allowed
67 *condition = NotForced;
69 G4double proposedStep = DBL_MAX;
70 G4double minStep = (1.0e-9)*m;
75 = TG4GeometryServices::Instance()
76 ->GetLimits(track.GetVolume()->GetLogicalVolume()->GetUserLimits());
79 G4String text = "TG4VSpecialCuts::PostStepGetPhysicalInteractionLength:\n";
80 text = text + " " + track.GetVolume()->GetLogicalVolume()->GetName();
81 text = text + " has not limits.";
82 TG4Globals::Exception(text);
86 = (TG4Limits*) track.GetVolume()->GetLogicalVolume()->GetUserLimits();
91 = (limits->GetUserMaxTrackLength(track) - track.GetTrackLength());
92 if (proposedStep < 0.) return minStep;
96 = (track.GetDynamicParticle()->GetTotalMomentum()) /
97 (track.GetTotalEnergy());
99 = (limits->GetUserMaxTime(track) - track.GetGlobalTime());
100 G4double temp = beta*c_light*time;
101 if (temp < 0.) return minStep;
102 if (proposedStep > temp) proposedStep = temp;
104 G4ParticleDefinition* particle = track.GetDefinition();
105 if (particle->GetPDGCharge() != 0.) {
106 // min remaining range
107 G4double kinEnergy = track.GetKineticEnergy();
108 G4Material* material = track.GetMaterial();
110 = G4EnergyLossTables::GetRange(particle, kinEnergy, material);
112 temp = (rangeNow - limits->GetUserMinRange(track));
113 if (temp < 0.) return minStep;
114 if (proposedStep > temp) proposedStep = temp;
116 // min kinetic energy (from limits)
117 // the kin energy cut can be applied only in case
118 // G4EnergyLossTables are defined for the particle
119 if (G4EnergyLossTables::GetDEDXTable(particle)) {
120 G4double minEkine = GetMinEkine(*limits, track);
121 G4EnergyLossTables::GetTables(particle);
123 = G4EnergyLossTables::GetRange(particle, minEkine, material);
124 temp = rangeNow - minR;
125 if (temp < 0.) return minStep;
126 if (proposedStep > temp) proposedStep = temp;
130 // min kinetic energy (from limits)
131 // for neutral particles
132 G4double minEkine = GetMinEkine(*limits, track);
133 if (track.GetKineticEnergy() <= minEkine) return minStep;
139 //_____________________________________________________________________________
140 G4VParticleChange* TG4VSpecialCuts::PostStepDoIt(const G4Track& track,
143 // Kills the current particle, if requested by G4UserLimits.
146 aParticleChange.Initialize(track);
147 aParticleChange.SetEnergyChange(0.) ;
148 aParticleChange.SetLocalEnergyDeposit(track.GetKineticEnergy()) ;
149 aParticleChange.SetStatusChange(fStopAndKill);
150 return &aParticleChange;