Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4HadronStoppingProcess.cc
Go to the documentation of this file.
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// $Id$
27//
28//---------------------------------------------------------------------
29//
30// GEANT4 Class
31//
32// File name: G4HadronStoppingProcess
33//
34// Author V.Ivanchenko 21 April 2012
35//
36//
37// Class Description:
38//
39// Base process class for nuclear capture of negatively charged particles
40//
41// Modifications:
42//
43// 20120522 M. Kelsey -- Set enableAtRestDoIt flag for G4ProcessManager
44// 20120914 M. Kelsey -- Pass subType in base ctor, remove enable flags
45// 20121004 K. Genser -- use G4HadronicProcessType in the constructor
46// 20121016 K. Genser -- Reverting to use one argument c'tor
47//
48//------------------------------------------------------------------------
49
53#include "G4EmCaptureCascade.hh"
54#include "G4Nucleus.hh"
55#include "G4HadFinalState.hh"
56#include "G4HadProjectile.hh"
57#include "G4HadSecondary.hh"
58#include "G4Material.hh"
59#include "G4Element.hh"
60
61//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
62
65{
66 // Modify G4VProcess flags to emulate G4VRest instead of G4VDiscrete
67 enableAtRestDoIt = true;
68 enablePostStepDoIt = false;
69
70 fElementSelector = new G4ElementSelector();
71 fEmCascade = new G4EmCaptureCascade(); // Owned by InteractionRegistry
72 fBoundDecay = 0;
74}
75
76//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
77
79{
81 delete fElementSelector;
82 // NOTE: fEmCascade and fEmBoundDecay owned by registry, not locally
83}
84
85//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
86
88{
89 return (p.GetPDGCharge() < 0.0);
90}
91
92//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
93
95{
97}
98
99//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
100
102{
104}
105
106//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
107
110{
112 return 0.0;
113}
114
115//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
116
119{
121 return DBL_MAX;
122}
123
124//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
125
127 const G4Step&)
128{
129 // if primary is not Alive then do nothing
131
133 G4Element* elm = fElementSelector->SelectZandA(track, nucleus);
134
135 G4HadFinalState* result = 0;
136 thePro.Initialise(track);
137 G4double time0 = track.GetGlobalTime();
138 G4bool nuclearCapture = true;
139
140 // Do the electromagnetic cascade in the nuclear field.
141 // EM cascade should keep G4HadFinalState object,
142 // because it will not be deleted at the end of this method
143 //
144 result = fEmCascade->ApplyYourself(thePro, *nucleus);
145 G4double ebound = result->GetLocalEnergyDeposit();
146 G4double edep = 0.0;
147 G4int nSecondaries = result->GetNumberOfSecondaries();
148
149 // Try decay from bound level
150 // For mu- the time of projectile should be changed.
151 // Decay should keep G4HadFinalState object,
152 // because it will not be deleted at the end of this method
153 //
154 thePro.SetBoundEnergy(ebound);
155 if(fBoundDecay) {
156 G4HadFinalState* resultDecay =
157 fBoundDecay->ApplyYourself(thePro, *nucleus);
158 G4int n = resultDecay->GetNumberOfSecondaries();
159 if(0 < n) {
160 nSecondaries += n;
161 result->AddSecondaries(resultDecay);
162 }
163 if(resultDecay->GetStatusChange() == stopAndKill) {
164 nuclearCapture = false;
165 }
166 resultDecay->Clear();
167 }
168
169 if(nuclearCapture) {
170
171 // select model
172 G4HadronicInteraction* model = 0;
173 try {
174 model = ChooseHadronicInteraction(0.0, track.GetMaterial(), elm);
175 }
176 catch(G4HadronicException & aE) {
178 ed << "Target element "<<elm->GetName()<<" Z= "
179 << nucleus->GetZ_asInt() << " A= "
180 << nucleus->GetA_asInt() << G4endl;
181 DumpState(track,"ChooseHadronicInteraction",ed);
182 ed << " No HadronicInteraction found out" << G4endl;
183 G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had005",
184 FatalException, ed);
185 }
186
187 G4HadFinalState* resultNuc = 0;
188 G4int reentryCount = 0;
189 do {
190 // sample final state
191 // nuclear interaction should keep G4HadFinalState object
192 // model should define time of each secondary particle
193 try {
194 resultNuc = model->ApplyYourself(thePro, *nucleus);
195 ++reentryCount;
196 }
197 catch(G4HadronicException aR) {
199 ed << "Call for " << model->GetModelName() << G4endl;
200 ed << "Target element "<<elm->GetName()<<" Z= "
201 << nucleus->GetZ_asInt()
202 << " A= " << nucleus->GetA_asInt() << G4endl;
203 DumpState(track,"ApplyYourself",ed);
204 ed << " ApplyYourself failed" << G4endl;
205 G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006",
206 FatalException, ed);
207 }
208
209 // Check the result for catastrophic energy non-conservation
210 resultNuc = CheckResult(thePro, *nucleus, resultNuc);
211
212 if(reentryCount>100) {
214 ed << "Call for " << model->GetModelName() << G4endl;
215 ed << "Target element "<<elm->GetName()<<" Z= "
216 << nucleus->GetZ_asInt()
217 << " A= " << nucleus->GetA_asInt() << G4endl;
218 DumpState(track,"ApplyYourself",ed);
219 ed << " ApplyYourself does not completed after 100 attempts" << G4endl;
220 G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006",
221 FatalException, ed);
222 }
223 }
224 while(!resultNuc);
225
226 edep = resultNuc->GetLocalEnergyDeposit();
227 nSecondaries += resultNuc->GetNumberOfSecondaries();
228 result->AddSecondaries(resultNuc);
229 resultNuc->Clear();
230 }
231
232 // Fill results
233 //
237 G4double w = track.GetWeight();
239 for(G4int i=0; i<nSecondaries; ++i) {
240 G4HadSecondary* sec = result->GetSecondary(i);
241
242 // add track global time to the reaction time
243 G4double time = sec->GetTime();
244 if(time < 0.0) { time = 0.0; }
245 time += time0;
246
247 // create secondary track
248 G4Track* t = new G4Track(sec->GetParticle(),
249 time,
250 track.GetPosition());
251 t->SetWeight(w*sec->GetWeight());
254 }
255 result->Clear();
256
257 if (epReportLevel != 0) {
258 CheckEnergyMomentumConservation(track, *nucleus);
259 }
260 return theTotalResult;
261}
262
263//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
264
265void G4HadronStoppingProcess::ProcessDescription(std::ostream& outFile) const
266{
267 outFile << "Base process for negatively charged particle capture at rest.\n";
268}
269
270//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
G4double condition(const G4ErrorSymMatrix &m)
@ FatalException
G4ForceCondition
@ NotForced
@ stopAndKill
@ fHadronAtRest
@ fStopAndKill
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
#define G4endl
Definition: G4ios.hh:52
virtual G4Element * SelectZandA(const G4Track &track, G4Nucleus *)
const G4String & GetName() const
Definition: G4Element.hh:127
G4HadFinalStateStatus GetStatusChange() const
G4double GetLocalEnergyDeposit() const
G4int GetNumberOfSecondaries() const
void AddSecondaries(const std::vector< G4HadSecondary > &addSecs)
G4HadSecondary * GetSecondary(size_t i)
void Initialise(const G4Track &aT)
void SetBoundEnergy(G4double e)
G4DynamicParticle * GetParticle()
G4double GetWeight() const
G4double GetTime() const
virtual G4bool IsApplicable(const G4ParticleDefinition &)
virtual void ProcessDescription(std::ostream &outFile) const
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
G4HadronStoppingProcess(const G4String &name="hadronCaptureAtRest")
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition)
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
virtual G4VParticleChange * AtRestDoIt(const G4Track &, const G4Step &)
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
const G4String & GetModelName() const
void DeRegisterExtraProcess(G4VProcess *)
void RegisterExtraProcess(G4VProcess *)
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)
static G4HadronicProcessStore * Instance()
void PrintInfo(const G4ParticleDefinition *)
G4HadProjectile thePro
G4Nucleus * GetTargetNucleusPointer()
G4ParticleChange * theTotalResult
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
G4HadronicInteraction * ChooseHadronicInteraction(G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement)
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
G4HadFinalState * CheckResult(const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result) const
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void AddSecondary(G4Track *aSecondary)
virtual void Initialize(const G4Track &)
G4double GetPDGCharge() const
Definition: G4Step.hh:78
G4double GetWeight() const
void SetWeight(G4double aValue)
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime() const
G4Material * GetMaterial() const
const G4TouchableHandle & GetTouchableHandle() const
void ProposeTrackStatus(G4TrackStatus status)
void ProposeWeight(G4double finalWeight)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetNumberOfSecondaries(G4int totSecondaries)
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
#define DBL_MAX
Definition: templates.hh:83