Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4UrbanMscModel.hh
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//
27// -------------------------------------------------------------------
28//
29//
30// GEANT4 Class header file
31//
32//
33// File name: G4UrbanMscModel
34//
35// Author: Laszlo Urban
36//
37// Creation date: 19.02.2013
38//
39// Created from G4UrbanMscModel96
40//
41// New parametrization for theta0
42// Correction for very small step length
43//
44// Class Description:
45//
46// Implementation of the model of multiple scattering based on
47// H.W.Lewis Phys Rev 78 (1950) 526 and L.Urban model
48
49// -------------------------------------------------------------------
50//
51
52#ifndef G4UrbanMscModel_h
53#define G4UrbanMscModel_h 1
54
55//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56
58
59#include "G4VMscModel.hh"
60#include "G4MscStepLimitType.hh"
61#include "G4Log.hh"
62#include "G4Exp.hh"
63
65class G4SafetyHelper;
67
68//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69
71{
72
73public:
74
75 explicit G4UrbanMscModel(const G4String& nam = "UrbanMsc");
76
77 ~G4UrbanMscModel() override;
78
80 const G4DataVector&) override;
81
82 void StartTracking(G4Track*) override;
83
86 G4double KineticEnergy,
87 G4double AtomicNumber,
88 G4double AtomicWeight=0.,
89 G4double cut =0.,
90 G4double emax=DBL_MAX) override;
91
93 G4double safety) override;
94
96 G4double& currentMinimalStep) override;
97
98 G4double ComputeGeomPathLength(G4double truePathLength) override;
99
100 G4double ComputeTrueStepLength(G4double geomStepLength) override;
101
102 G4double ComputeTheta0(G4double truePathLength, G4double KineticEnergy);
103
104 // hide assignment operator
105 G4UrbanMscModel & operator=(const G4UrbanMscModel &right) = delete;
107
108private:
109
110 G4double SampleCosineTheta(G4double trueStepLength, G4double KineticEnergy);
111
112 void SampleDisplacement(G4double sinTheta, G4double phi);
113
114 void SampleDisplacementNew(G4double sinTheta, G4double phi);
115
116 void InitialiseModelCache();
117
118 inline void SetParticle(const G4ParticleDefinition*);
119
120 inline G4double Randomizetlimit();
121
122 inline G4double SimpleScattering(G4double xmeanth, G4double x2meanth);
123
124 inline G4double ComputeStepmin();
125
126 inline G4double ComputeTlimitmin();
127
128 CLHEP::HepRandomEngine* rndmEngineMod;
129
130 const G4ParticleDefinition* particle;
131 const G4ParticleDefinition* positron;
132 G4ParticleChangeForMSC* fParticleChange;
133
134 const G4MaterialCutsCouple* couple;
135 G4LossTableManager* theManager;
136
137 G4double mass;
138 G4double charge,ChargeSquare;
139 G4double masslimite,fr;
140
141 G4double taubig;
142 G4double tausmall;
143 G4double taulim;
144 G4double currentTau;
145 G4double tlimit;
146 G4double tlimitmin;
147 G4double tlimitminfix,tlimitminfix2;
148 G4double tgeom;
149
150 G4double geombig;
151 G4double geommin;
152 G4double geomlimit;
153 G4double skindepth;
154 G4double smallstep;
155
156 G4double presafety;
157
158 G4double lambda0;
159 G4double lambdaeff;
160 G4double tPathLength;
161 G4double zPathLength;
162 G4double par1,par2,par3;
163
164 G4double stepmin;
165
166 G4double currentKinEnergy;
167 G4double currentLogKinEnergy;
168 G4double currentRange;
169 G4double rangeinit;
170 G4double currentRadLength;
171
172 G4double rangecut;
173 G4double drr,finalr;
174
175 G4double tlow;
176 G4double invmev;
177
178 struct mscData {
179 G4double ecut, Zeff, Z23, sqrtZ;
180 G4double coeffth1, coeffth2;
181 G4double coeffc1, coeffc2, coeffc3, coeffc4;
182 G4double stepmina, stepminb;
183 G4double doverra, doverrb;
184 };
185 static std::vector<mscData*> msc;
186
187 // index of G4MaterialCutsCouple
188 G4int idx;
189
190 G4bool firstStep;
191 G4bool insideskin;
192
193 G4bool latDisplasmentbackup;
194 G4bool dispAlg96;
195
196};
197
198//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
199//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
200
201inline
202void G4UrbanMscModel::SetParticle(const G4ParticleDefinition* p)
203{
204 if (p != particle) {
205 particle = p;
206 mass = p->GetPDGMass();
207 charge = p->GetPDGCharge()/CLHEP::eplus;
208 ChargeSquare = charge*charge;
209 }
210}
211
212//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
213
214inline G4double G4UrbanMscModel::Randomizetlimit()
215{
216 G4double res = tlimitmin;
217 if(tlimit > tlimitmin)
218 {
219 res = G4RandGauss::shoot(rndmEngineMod,tlimit,0.1*(tlimit-tlimitmin));
220 res = std::max(res, tlimitmin);
221 }
222 return res;
223}
224
225//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
226
227inline
228G4double G4UrbanMscModel::SimpleScattering(G4double xmeanth, G4double x2meanth)
229{
230 // 'large angle scattering'
231 // 2 model functions with correct xmean and x2mean
232 G4double a = (2.*xmeanth+9.*x2meanth-3.)/(2.*xmeanth-3.*x2meanth+1.);
233 G4double prob = (a+2.)*xmeanth/a;
234
235 // sampling
236 G4double rdm = rndmEngineMod->flat();
237 G4double cth = (rndmEngineMod->flat() < prob)
238 ? -1.+2.*G4Exp(G4Log(rdm)/(a+1.)) : -1.+2.*rdm;
239 return cth;
240}
241
242inline G4double G4UrbanMscModel::ComputeStepmin()
243{
244 // define stepmin using estimation of the ratio
245 // of lambda_elastic/lambda_transport
246 G4double rat = currentKinEnergy*invmev;
247 return lambda0*1.e-3/(2.e-3+rat*(msc[idx]->stepmina+msc[idx]->stepminb*rat));
248}
249
250inline G4double G4UrbanMscModel::ComputeTlimitmin()
251{
252 G4double x = (particle == positron) ?
253 0.7*msc[idx]->sqrtZ*stepmin : 0.87*msc[idx]->Z23*stepmin;
254 if(currentKinEnergy < tlow) { x *= 0.5*(1.+currentKinEnergy/tlow); }
255 return std::max(x, tlimitminfix);
256}
257
258//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
259
260
261#endif
262
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179
G4double G4Log(G4double x)
Definition: G4Log.hh:226
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
virtual double flat()=0
G4double GetPDGCharge() const
G4ThreeVector & SampleScattering(const G4ThreeVector &, G4double safety) override
G4double ComputeTrueStepLength(G4double geomStepLength) override
G4double ComputeTheta0(G4double truePathLength, G4double KineticEnergy)
void StartTracking(G4Track *) override
~G4UrbanMscModel() override
void Initialise(const G4ParticleDefinition *, const G4DataVector &) override
G4double ComputeGeomPathLength(G4double truePathLength) override
G4UrbanMscModel & operator=(const G4UrbanMscModel &right)=delete
G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *particle, G4double KineticEnergy, G4double AtomicNumber, G4double AtomicWeight=0., G4double cut=0., G4double emax=DBL_MAX) override
G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &currentMinimalStep) override
G4UrbanMscModel(const G4UrbanMscModel &)=delete
#define DBL_MAX
Definition: templates.hh:62