4 // Author: I. Hrivnacova
6 // Class TG4SpecialCuts
7 // --------------------
8 // See the class description in the header file.
10 #include "TG4SpecialCuts.h"
11 #include "TG4G3CutVector.h"
12 #include "TG4Limits.h"
14 #include <G4UserLimits.hh>
16 #include <G4EnergyLossTables.hh>
18 //_____________________________________________________________________________
19 TG4SpecialCuts::TG4SpecialCuts(TG4G3ParticleWSP particle,
20 TG4G3CutVector* cutVector,
21 const G4String& processName)
22 : G4UserSpecialCuts(processName),
28 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForGamma;
29 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForGamma;
31 case kElectron: case kEplus:
32 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForElectron;
33 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForElectron;
36 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForChargedHadron;
37 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForChargedHadron;
40 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForNeutralHadron;
41 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForNeutralHadron;
44 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForMuon;
45 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForMuon;
48 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForOther;
49 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForOther;
51 case kNofParticlesWSP:
52 TG4Globals::Exception("TG4SpecialCuts: Wrong particle specification.");
57 //_____________________________________________________________________________
58 TG4SpecialCuts::TG4SpecialCuts() {
62 //_____________________________________________________________________________
63 TG4SpecialCuts::TG4SpecialCuts(const TG4SpecialCuts& right) {
65 TG4Globals::Exception(
66 "TG4SpecialCuts is protected from copying.");
69 //_____________________________________________________________________________
70 TG4SpecialCuts::~TG4SpecialCuts() {
76 //_____________________________________________________________________________
77 TG4SpecialCuts& TG4SpecialCuts::operator=(const TG4SpecialCuts& right)
79 // check assignement to self
80 if (this == &right) return *this;
82 TG4Globals::Exception(
83 "TG4SpecialCuts is protected from assigning.");
90 //_____________________________________________________________________________
91 G4double TG4SpecialCuts::PostStepGetPhysicalInteractionLength(
92 const G4Track& track, G4double previousStepSize,
93 G4ForceCondition* condition)
95 // Returns the Step-size (actual length) which is allowed
100 *condition = NotForced;
102 G4double proposedStep = DBL_MAX;
104 = track.GetVolume()->GetLogicalVolume()->GetUserLimits();
108 = (limits->GetUserMaxTrackLength(track) - track.GetTrackLength());
109 if (proposedStep < 0.) return 0.;
113 = (track.GetDynamicParticle()->GetTotalMomentum()) /
114 (track.GetTotalEnergy());
116 = (limits->GetUserMaxTime(track) - track.GetGlobalTime());
117 G4double temp = beta*c_light*time;
118 if (temp < 0.) return 0.;
119 if (proposedStep > temp) proposedStep = temp;
121 // min remaining range
122 G4ParticleDefinition* particle = track.GetDefinition();
123 if (particle->GetPDGCharge() != 0.) {
124 G4double kinEnergy = track.GetKineticEnergy();
125 G4Material* material = track.GetMaterial();
127 = G4EnergyLossTables::GetRange(particle, kinEnergy, material);
128 temp = (rangeNow - limits->GetUserMinRange(track));
129 if (temp < 0.) return 0.;
130 if (proposedStep > temp) proposedStep = temp;
132 // min kinetic energy (from limits)
133 // the kin energy cut can be applied only in case
134 // G4EnergyLossTables are defined for the particle
135 if (G4EnergyLossTables::GetDEDXTable(particle)) {
136 TG4Limits* tg4Limits = dynamic_cast<TG4Limits*>(limits);
138 G4String text = "TG4SpecialCuts::PostStepGetPhysicalInteractionLength:\n";
139 text = text + " Unknown limits type.";
140 TG4Globals::Exception(text);
143 = (tg4Limits->*fPtrMinEkineInLimits)(track);
144 G4EnergyLossTables::GetTables(particle);
146 = G4EnergyLossTables::GetRange(particle, minEkine, material);
147 temp = rangeNow - minR;
148 if (temp < 0.) return 0.;
149 if (proposedStep > temp) proposedStep = temp;
154 else if (fCutVector) {
155 // min kinetic energy (from cut vector)
156 // the kin energy cut can be applied only in case
157 // G4EnergyLossTables are defined for the particle
158 G4ParticleDefinition* particle = track.GetDefinition();
159 if (G4EnergyLossTables::GetDEDXTable(particle)) {
160 G4double kinEnergy = track.GetKineticEnergy();
161 G4Material* material = track.GetMaterial();
163 = G4EnergyLossTables::GetRange(particle, kinEnergy, material);
165 = (fCutVector->*fPtrMinEkineInCutVector)(track);
167 = G4EnergyLossTables::GetRange(particle, minEkine, material);
168 G4double temp = rangeNow - minR;
169 if (temp < 0.) return 0.;
170 if (proposedStep > temp) proposedStep = temp;
176 //_____________________________________________________________________________
177 G4VParticleChange* TG4SpecialCuts::PostStepDoIt(const G4Track& track,
180 // Kills the current particle, if requested by G4UserLimits.
183 aParticleChange.Initialize(track);
184 aParticleChange.SetEnergyChange(0.) ;
185 aParticleChange.SetLocalEnergyDeposit(track.GetKineticEnergy()) ;
186 aParticleChange.SetStatusChange(fStopAndKill);
187 return &aParticleChange;