4 // See the class description in the header file.
6 #include "TG4SpecialCuts.h"
7 #include "TG4G3CutVector.h"
10 #include <G4UserLimits.hh>
12 #include <G4EnergyLossTables.hh>
14 TG4SpecialCuts::TG4SpecialCuts(TG4G3ParticleWSP particle,
15 TG4G3CutVector* cutVector,
16 const G4String& processName)
17 : G4UserSpecialCuts(processName),
23 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForGamma;
24 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForGamma;
26 case kElectron: case kEplus:
27 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForElectron;
28 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForElectron;
31 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForHadron;
32 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForHadron;
35 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForNeutralHadron;
36 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForNeutralHadron;
39 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForMuon;
40 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForMuon;
43 fPtrMinEkineInCutVector = &TG4G3CutVector::GetMinEkineForOther;
44 fPtrMinEkineInLimits = &TG4Limits::GetMinEkineForOther;
46 case kNofParticlesWSP:
47 TG4Globals::Exception("TG4SpecialCuts: Wrong particle specification.");
52 TG4SpecialCuts::TG4SpecialCuts() {
56 TG4SpecialCuts::TG4SpecialCuts(const TG4SpecialCuts& right) {
58 TG4Globals::Exception(
59 "TG4SpecialCuts is protected from copying.");
62 TG4SpecialCuts::~TG4SpecialCuts() {
68 TG4SpecialCuts& TG4SpecialCuts::operator=(const TG4SpecialCuts& right)
70 // check assignement to self
71 if (this == &right) return *this;
73 TG4Globals::Exception(
74 "TG4SpecialCuts is protected from assigning.");
81 G4double TG4SpecialCuts::PostStepGetPhysicalInteractionLength(
82 const G4Track& track, G4double previousStepSize,
83 G4ForceCondition* condition)
85 // Returns the Step-size (actual length) which is allowed
90 *condition = NotForced;
92 G4double proposedStep = DBL_MAX;
94 = track.GetVolume()->GetLogicalVolume()->GetUserLimits();
98 = (limits->GetUserMaxTrackLength(track) - track.GetTrackLength());
99 if (proposedStep < 0.) return 0.;
103 = (track.GetDynamicParticle()->GetTotalMomentum()) /
104 (track.GetTotalEnergy());
106 = (limits->GetUserMaxTime(track) - track.GetGlobalTime());
107 G4double temp = beta*c_light*time;
108 if (temp < 0.) return 0.;
109 if (proposedStep > temp) proposedStep = temp;
111 // min remaining range
112 G4ParticleDefinition* particle = track.GetDefinition();
113 G4double kinEnergy = track.GetKineticEnergy();
114 G4Material* material = track.GetMaterial();
116 = G4EnergyLossTables::GetRange(particle, kinEnergy, material);
117 temp = (rangeNow - limits->GetUserMinRange(track));
118 if (temp < 0.) return 0.;
119 if (proposedStep > temp) proposedStep = temp;
121 // min kinetic energy (from limits)
122 // the kin energy cut can be applied only in case
123 // G4EnergyLossTables are defined for the particle
124 if (G4EnergyLossTables::GetDEDXTable(particle)) {
125 TG4Limits* tg4Limits = dynamic_cast<TG4Limits*>(limits);
127 G4String text = "TG4SpecialCuts::PostStepGetPhysicalInteractionLength:\n";
128 text = text + " Unknown limits type.";
129 TG4Globals::Exception(text);
132 = (tg4Limits->*fPtrMinEkineInLimits)(track);
134 = G4EnergyLossTables::GetRange(particle, minEkine, material);
135 temp = rangeNow - minR;
136 if (temp < 0.) return 0.;
137 if (proposedStep > temp) proposedStep = temp;
141 else if (fCutVector) {
142 // min kinetic energy (from cut vector)
143 // the kin energy cut can be applied only in case
144 // G4EnergyLossTables are defined for the particle
145 G4ParticleDefinition* particle = track.GetDefinition();
146 if (G4EnergyLossTables::GetDEDXTable(particle)) {
147 G4double kinEnergy = track.GetKineticEnergy();
148 G4Material* material = track.GetMaterial();
150 = G4EnergyLossTables::GetRange(particle, kinEnergy, material);
152 = (fCutVector->*fPtrMinEkineInCutVector)(track);
154 = G4EnergyLossTables::GetRange(particle, minEkine, material);
155 G4double temp = rangeNow - minR;
156 if (temp < 0.) return 0.;
157 if (proposedStep > temp) proposedStep = temp;
163 G4VParticleChange* TG4SpecialCuts::PostStepDoIt(const G4Track& track,
166 // Kills the current particle, if requested by G4UserLimits.
169 aParticleChange.Initialize(track);
170 aParticleChange.SetEnergyChange(0.) ;
171 aParticleChange.SetLocalEnergyDeposit(track.GetKineticEnergy()) ;
172 aParticleChange.SetStatusChange(fStopAndKill);
173 return &aParticleChange;