Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
G4UrbanMscModel93.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//
27// $Id$
28//
29// -------------------------------------------------------------------
30//
31// GEANT4 Class file
32//
33//
34// File name: G4UrbanMscModel93
35//
36// Author: Laszlo Urban
37//
38// Creation date: 06.03.2008
39//
40// Modifications:
41//
42// 06-03-2008 starting point : G4UrbanMscModel2 = G4UrbanMscModel 9.1 ref 02
43//
44// 13-03-08 Bug in SampleScattering (which caused lateral asymmetry) fixed
45// (L.Urban)
46//
47// 14-03-08 Simplification of step limitation in ComputeTruePathLengthLimit,
48// + tlimitmin is the same for UseDistancetoBoundary and
49// UseSafety (L.Urban)
50//
51// 16-03-08 Reorganization of SampleCosineTheta + new method SimpleScattering
52// SimpleScattering is used if the relative energy loss is too big
53// or theta0 is too big (see data members rellossmax, theta0max)
54// (L.Urban)
55//
56// 17-03-08 tuning of the correction factor in ComputeTheta0 (L.Urban)
57//
58// 19-03-08 exponent c of the 'tail' model function is not equal to 2 any more,
59// value of c has been extracted from some e- scattering data (L.Urban)
60//
61// 24-03-08 Step limitation in ComputeTruePathLengthLimit has been
62// simplified further + some data members have been removed (L.Urban)
63//
64// 24-07-08 central part of scattering angle (theta0) has been tuned
65// tail of the scattering angle distribution has been tuned
66// using some e- and proton scattering data
67//
68// 05-08-08 bugfix in ComputeTruePathLengthLimit (L.Urban)
69//
70// 09-10-08 theta0 and tail have been retuned using some e-,mu,proton
71// scattering data (L.Urban)
72// + single scattering without path length correction for
73// small steps (t < tlimitmin, for UseDistanceToBoundary only)
74//
75// 15-10-08 Moliere-Bethe screening in the single scattering part(L.Urban)
76//
77// 17-10-08 stepping similar to that in model (9.1) for UseSafety case
78// for e+/e- in order to speed up the code for calorimeters
79//
80// 23-10-08 bugfix in the screeningparameter of the single scattering part,
81// some technical change in order to speed up the code (UpdateCache)
82//
83// 27-10-08 bugfix in ComputeTruePathLengthLimit (affects UseDistanceToBoundary
84// stepping type only) (L.Urban)
85//
86// 28-10-09 V.Ivanchenko moved G4UrbanMscModel2 to G4UrbanMscModel93,
87// now it is a frozen version of the Urban model corresponding
88// to g4 9.3
89
90// Class Description:
91//
92// Implementation of the model of multiple scattering based on
93// H.W.Lewis Phys Rev 78 (1950) 526 and others
94
95// -------------------------------------------------------------------
96// In its present form the model can be used for simulation
97// of the e-/e+, muon and charged hadron multiple scattering
98//
99
100
101//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
102//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
103
104#include "G4UrbanMscModel93.hh"
105#include "G4PhysicalConstants.hh"
106#include "G4SystemOfUnits.hh"
107#include "Randomize.hh"
108#include "G4Electron.hh"
109#include "G4LossTableManager.hh"
110#include "G4ParticleChangeForMSC.hh"
111
112#include "G4Poisson.hh"
113#include "globals.hh"
114
115//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
116
117using namespace std;
118
119G4UrbanMscModel93::G4UrbanMscModel93(const G4String& nam)
120 : G4VMscModel(nam)
121{
122 masslimite = 0.6*MeV;
123 lambdalimit = 1.*mm;
124 fr = 0.02;
125 taubig = 8.0;
126 tausmall = 1.e-16;
127 taulim = 1.e-6;
128 currentTau = taulim;
129 tlimitminfix = 1.e-6*mm;
130 stepmin = tlimitminfix;
131 smallstep = 1.e10;
132 currentRange = 0. ;
133 rangeinit = 0.;
134 tlimit = 1.e10*mm;
135 tlimitmin = 10.*tlimitminfix;
136 tgeom = 1.e50*mm;
137 geombig = 1.e50*mm;
138 geommin = 1.e-3*mm;
139 geomlimit = geombig;
140 presafety = 0.*mm;
141
142 y = 0.;
143
144 Zold = 0.;
145 Zeff = 1.;
146 Z2 = 1.;
147 Z23 = 1.;
148 lnZ = 0.;
149 coeffth1 = 0.;
150 coeffth2 = 0.;
151 coeffc1 = 0.;
152 coeffc2 = 0.;
153 scr1ini = fine_structure_const*fine_structure_const*
154 electron_mass_c2*electron_mass_c2/(0.885*0.885*4.*pi);
155 scr2ini = 3.76*fine_structure_const*fine_structure_const;
156 scr1 = 0.;
157 scr2 = 0.;
158
159 theta0max = pi/6.;
160 rellossmax = 0.50;
161 third = 1./3.;
162 particle = 0;
163 theManager = G4LossTableManager::Instance();
164 firstStep = true;
165 inside = false;
166 insideskin = false;
167
168 skindepth = skin*stepmin;
169
170 mass = proton_mass_c2;
171 charge = ChargeSquare = 1.0;
172 currentKinEnergy = currentRadLength = lambda0 = lambdaeff = tPathLength
173 = zPathLength = par1 = par2 = par3 = 0;
174
175 currentMaterialIndex = -1;
176 fParticleChange = 0;
177 couple = 0;
178}
179
180//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
181
182G4UrbanMscModel93::~G4UrbanMscModel93()
183{}
184
185//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
186
187void G4UrbanMscModel93::Initialise(const G4ParticleDefinition* p,
188 const G4DataVector&)
189{
190 skindepth = skin*stepmin;
191
192 // set values of some data members
193 SetParticle(p);
194
195 if(p->GetPDGMass() > MeV) {
196 G4cout << "### WARNING: G4UrbanMscModel93 model is used for "
197 << p->GetParticleName() << " !!! " << G4endl;
198 G4cout << "### This model should be used only for e+-"
199 << G4endl;
200 }
201
202 fParticleChange = GetParticleChangeForMSC(p);
203}
204
205//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
206
207G4double G4UrbanMscModel93::ComputeCrossSectionPerAtom(
208 const G4ParticleDefinition* part,
209 G4double KineticEnergy,
210 G4double AtomicNumber,G4double,
211 G4double, G4double)
212{
213 const G4double sigmafactor = twopi*classic_electr_radius*classic_electr_radius;
214 const G4double epsfactor = 2.*electron_mass_c2*electron_mass_c2*
215 Bohr_radius*Bohr_radius/(hbarc*hbarc);
216 const G4double epsmin = 1.e-4 , epsmax = 1.e10;
217
218 const G4double Zdat[15] = { 4., 6., 13., 20., 26., 29., 32., 38., 47.,
219 50., 56., 64., 74., 79., 82. };
220
221 const G4double Tdat[22] = { 100*eV, 200*eV, 400*eV, 700*eV,
222 1*keV, 2*keV, 4*keV, 7*keV,
223 10*keV, 20*keV, 40*keV, 70*keV,
224 100*keV, 200*keV, 400*keV, 700*keV,
225 1*MeV, 2*MeV, 4*MeV, 7*MeV,
226 10*MeV, 20*MeV};
227
228 // corr. factors for e-/e+ lambda for T <= Tlim
229 G4double celectron[15][22] =
230 {{1.125,1.072,1.051,1.047,1.047,1.050,1.052,1.054,
231 1.054,1.057,1.062,1.069,1.075,1.090,1.105,1.111,
232 1.112,1.108,1.100,1.093,1.089,1.087 },
233 {1.408,1.246,1.143,1.096,1.077,1.059,1.053,1.051,
234 1.052,1.053,1.058,1.065,1.072,1.087,1.101,1.108,
235 1.109,1.105,1.097,1.090,1.086,1.082 },
236 {2.833,2.268,1.861,1.612,1.486,1.309,1.204,1.156,
237 1.136,1.114,1.106,1.106,1.109,1.119,1.129,1.132,
238 1.131,1.124,1.113,1.104,1.099,1.098 },
239 {3.879,3.016,2.380,2.007,1.818,1.535,1.340,1.236,
240 1.190,1.133,1.107,1.099,1.098,1.103,1.110,1.113,
241 1.112,1.105,1.096,1.089,1.085,1.098 },
242 {6.937,4.330,2.886,2.256,1.987,1.628,1.395,1.265,
243 1.203,1.122,1.080,1.065,1.061,1.063,1.070,1.073,
244 1.073,1.070,1.064,1.059,1.056,1.056 },
245 {9.616,5.708,3.424,2.551,2.204,1.762,1.485,1.330,
246 1.256,1.155,1.099,1.077,1.070,1.068,1.072,1.074,
247 1.074,1.070,1.063,1.059,1.056,1.052 },
248 {11.72,6.364,3.811,2.806,2.401,1.884,1.564,1.386,
249 1.300,1.180,1.112,1.082,1.073,1.066,1.068,1.069,
250 1.068,1.064,1.059,1.054,1.051,1.050 },
251 {18.08,8.601,4.569,3.183,2.662,2.025,1.646,1.439,
252 1.339,1.195,1.108,1.068,1.053,1.040,1.039,1.039,
253 1.039,1.037,1.034,1.031,1.030,1.036 },
254 {18.22,10.48,5.333,3.713,3.115,2.367,1.898,1.631,
255 1.498,1.301,1.171,1.105,1.077,1.048,1.036,1.033,
256 1.031,1.028,1.024,1.022,1.021,1.024 },
257 {14.14,10.65,5.710,3.929,3.266,2.453,1.951,1.669,
258 1.528,1.319,1.178,1.106,1.075,1.040,1.027,1.022,
259 1.020,1.017,1.015,1.013,1.013,1.020 },
260 {14.11,11.73,6.312,4.240,3.478,2.566,2.022,1.720,
261 1.569,1.342,1.186,1.102,1.065,1.022,1.003,0.997,
262 0.995,0.993,0.993,0.993,0.993,1.011 },
263 {22.76,20.01,8.835,5.287,4.144,2.901,2.219,1.855,
264 1.677,1.410,1.224,1.121,1.073,1.014,0.986,0.976,
265 0.974,0.972,0.973,0.974,0.975,0.987 },
266 {50.77,40.85,14.13,7.184,5.284,3.435,2.520,2.059,
267 1.837,1.512,1.283,1.153,1.091,1.010,0.969,0.954,
268 0.950,0.947,0.949,0.952,0.954,0.963 },
269 {65.87,59.06,15.87,7.570,5.567,3.650,2.682,2.182,
270 1.939,1.579,1.325,1.178,1.108,1.014,0.965,0.947,
271 0.941,0.938,0.940,0.944,0.946,0.954 },
272 {55.60,47.34,15.92,7.810,5.755,3.767,2.760,2.239,
273 1.985,1.609,1.343,1.188,1.113,1.013,0.960,0.939,
274 0.933,0.930,0.933,0.936,0.939,0.949 }};
275
276 G4double cpositron[15][22] = {
277 {2.589,2.044,1.658,1.446,1.347,1.217,1.144,1.110,
278 1.097,1.083,1.080,1.086,1.092,1.108,1.123,1.131,
279 1.131,1.126,1.117,1.108,1.103,1.100 },
280 {3.904,2.794,2.079,1.710,1.543,1.325,1.202,1.145,
281 1.122,1.096,1.089,1.092,1.098,1.114,1.130,1.137,
282 1.138,1.132,1.122,1.113,1.108,1.102 },
283 {7.970,6.080,4.442,3.398,2.872,2.127,1.672,1.451,
284 1.357,1.246,1.194,1.179,1.178,1.188,1.201,1.205,
285 1.203,1.190,1.173,1.159,1.151,1.145 },
286 {9.714,7.607,5.747,4.493,3.815,2.777,2.079,1.715,
287 1.553,1.353,1.253,1.219,1.211,1.214,1.225,1.228,
288 1.225,1.210,1.191,1.175,1.166,1.174 },
289 {17.97,12.95,8.628,6.065,4.849,3.222,2.275,1.820,
290 1.624,1.382,1.259,1.214,1.202,1.202,1.214,1.219,
291 1.217,1.203,1.184,1.169,1.160,1.151 },
292 {24.83,17.06,10.84,7.355,5.767,3.707,2.546,1.996,
293 1.759,1.465,1.311,1.252,1.234,1.228,1.238,1.241,
294 1.237,1.222,1.201,1.184,1.174,1.159 },
295 {23.26,17.15,11.52,8.049,6.375,4.114,2.792,2.155,
296 1.880,1.535,1.353,1.281,1.258,1.247,1.254,1.256,
297 1.252,1.234,1.212,1.194,1.183,1.170 },
298 {22.33,18.01,12.86,9.212,7.336,4.702,3.117,2.348,
299 2.015,1.602,1.385,1.297,1.268,1.251,1.256,1.258,
300 1.254,1.237,1.214,1.195,1.185,1.179 },
301 {33.91,24.13,15.71,10.80,8.507,5.467,3.692,2.808,
302 2.407,1.873,1.564,1.425,1.374,1.330,1.324,1.320,
303 1.312,1.288,1.258,1.235,1.221,1.205 },
304 {32.14,24.11,16.30,11.40,9.015,5.782,3.868,2.917,
305 2.490,1.925,1.596,1.447,1.391,1.342,1.332,1.327,
306 1.320,1.294,1.264,1.240,1.226,1.214 },
307 {29.51,24.07,17.19,12.28,9.766,6.238,4.112,3.066,
308 2.602,1.995,1.641,1.477,1.414,1.356,1.342,1.336,
309 1.328,1.302,1.270,1.245,1.231,1.233 },
310 {38.19,30.85,21.76,15.35,12.07,7.521,4.812,3.498,
311 2.926,2.188,1.763,1.563,1.484,1.405,1.382,1.371,
312 1.361,1.330,1.294,1.267,1.251,1.239 },
313 {49.71,39.80,27.96,19.63,15.36,9.407,5.863,4.155,
314 3.417,2.478,1.944,1.692,1.589,1.480,1.441,1.423,
315 1.409,1.372,1.330,1.298,1.280,1.258 },
316 {59.25,45.08,30.36,20.83,16.15,9.834,6.166,4.407,
317 3.641,2.648,2.064,1.779,1.661,1.531,1.482,1.459,
318 1.442,1.400,1.354,1.319,1.299,1.272 },
319 {56.38,44.29,30.50,21.18,16.51,10.11,6.354,4.542,
320 3.752,2.724,2.116,1.817,1.692,1.554,1.499,1.474,
321 1.456,1.412,1.364,1.328,1.307,1.282 }};
322
323 //data/corrections for T > Tlim
324 G4double Tlim = 10.*MeV;
325 G4double beta2lim = Tlim*(Tlim+2.*electron_mass_c2)/
326 ((Tlim+electron_mass_c2)*(Tlim+electron_mass_c2));
327 G4double bg2lim = Tlim*(Tlim+2.*electron_mass_c2)/
328 (electron_mass_c2*electron_mass_c2);
329
330 G4double sig0[15] = {0.2672*barn, 0.5922*barn, 2.653*barn, 6.235*barn,
331 11.69*barn , 13.24*barn , 16.12*barn, 23.00*barn ,
332 35.13*barn , 39.95*barn , 50.85*barn, 67.19*barn ,
333 91.15*barn , 104.4*barn , 113.1*barn};
334
335 G4double hecorr[15] = {120.70, 117.50, 105.00, 92.92, 79.23, 74.510, 68.29,
336 57.39, 41.97, 36.14, 24.53, 10.21, -7.855, -16.84,
337 -22.30};
338
339 G4double sigma;
340 SetParticle(part);
341
342 Z23 = pow(AtomicNumber,2./3.);
343
344 // correction if particle .ne. e-/e+
345 // compute equivalent kinetic energy
346 // lambda depends on p*beta ....
347
348 G4double eKineticEnergy = KineticEnergy;
349
350 if(mass > electron_mass_c2)
351 {
352 G4double TAU = KineticEnergy/mass ;
353 G4double c = mass*TAU*(TAU+2.)/(electron_mass_c2*(TAU+1.)) ;
354 G4double w = c-2. ;
355 G4double tau = 0.5*(w+sqrt(w*w+4.*c)) ;
356 eKineticEnergy = electron_mass_c2*tau ;
357 }
358
359 G4double eTotalEnergy = eKineticEnergy + electron_mass_c2 ;
360 G4double beta2 = eKineticEnergy*(eTotalEnergy+electron_mass_c2)
361 /(eTotalEnergy*eTotalEnergy);
362 G4double bg2 = eKineticEnergy*(eTotalEnergy+electron_mass_c2)
363 /(electron_mass_c2*electron_mass_c2);
364
365 G4double eps = epsfactor*bg2/Z23;
366
367 if (eps<epsmin) sigma = 2.*eps*eps;
368 else if(eps<epsmax) sigma = log(1.+2.*eps)-2.*eps/(1.+2.*eps);
369 else sigma = log(2.*eps)-1.+1./eps;
370
371 sigma *= ChargeSquare*AtomicNumber*AtomicNumber/(beta2*bg2);
372
373 // interpolate in AtomicNumber and beta2
374 G4double c1,c2,cc1,cc2,corr;
375
376 // get bin number in Z
377 G4int iZ = 14;
378 while ((iZ>=0)&&(Zdat[iZ]>=AtomicNumber)) iZ -= 1;
379 if (iZ==14) iZ = 13;
380 if (iZ==-1) iZ = 0 ;
381
382 G4double ZZ1 = Zdat[iZ];
383 G4double ZZ2 = Zdat[iZ+1];
384 G4double ratZ = (AtomicNumber-ZZ1)*(AtomicNumber+ZZ1)/
385 ((ZZ2-ZZ1)*(ZZ2+ZZ1));
386
387 if(eKineticEnergy <= Tlim)
388 {
389 // get bin number in T (beta2)
390 G4int iT = 21;
391 while ((iT>=0)&&(Tdat[iT]>=eKineticEnergy)) iT -= 1;
392 if(iT==21) iT = 20;
393 if(iT==-1) iT = 0 ;
394
395 // calculate betasquare values
396 G4double T = Tdat[iT], E = T + electron_mass_c2;
397 G4double b2small = T*(E+electron_mass_c2)/(E*E);
398
399 T = Tdat[iT+1]; E = T + electron_mass_c2;
400 G4double b2big = T*(E+electron_mass_c2)/(E*E);
401 G4double ratb2 = (beta2-b2small)/(b2big-b2small);
402
403 if (charge < 0.)
404 {
405 c1 = celectron[iZ][iT];
406 c2 = celectron[iZ+1][iT];
407 cc1 = c1+ratZ*(c2-c1);
408
409 c1 = celectron[iZ][iT+1];
410 c2 = celectron[iZ+1][iT+1];
411 cc2 = c1+ratZ*(c2-c1);
412
413 corr = cc1+ratb2*(cc2-cc1);
414
415 sigma *= sigmafactor/corr;
416 }
417 else
418 {
419 c1 = cpositron[iZ][iT];
420 c2 = cpositron[iZ+1][iT];
421 cc1 = c1+ratZ*(c2-c1);
422
423 c1 = cpositron[iZ][iT+1];
424 c2 = cpositron[iZ+1][iT+1];
425 cc2 = c1+ratZ*(c2-c1);
426
427 corr = cc1+ratb2*(cc2-cc1);
428
429 sigma *= sigmafactor/corr;
430 }
431 }
432 else
433 {
434 c1 = bg2lim*sig0[iZ]*(1.+hecorr[iZ]*(beta2-beta2lim))/bg2;
435 c2 = bg2lim*sig0[iZ+1]*(1.+hecorr[iZ+1]*(beta2-beta2lim))/bg2;
436 if((AtomicNumber >= ZZ1) && (AtomicNumber <= ZZ2))
437 sigma = c1+ratZ*(c2-c1) ;
438 else if(AtomicNumber < ZZ1)
439 sigma = AtomicNumber*AtomicNumber*c1/(ZZ1*ZZ1);
440 else if(AtomicNumber > ZZ2)
441 sigma = AtomicNumber*AtomicNumber*c2/(ZZ2*ZZ2);
442 }
443 return sigma;
444
445}
446
447//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
448
449void G4UrbanMscModel93::StartTracking(G4Track* track)
450{
451 SetParticle(track->GetDynamicParticle()->GetDefinition());
452 firstStep = true;
453 inside = false;
454 insideskin = false;
455 tlimit = geombig;
456}
457
458//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
459
460G4double G4UrbanMscModel93::ComputeTruePathLengthLimit(
461 const G4Track& track,
462 G4double& currentMinimalStep)
463{
464 tPathLength = currentMinimalStep;
465 const G4DynamicParticle* dp = track.GetDynamicParticle();
466 G4StepPoint* sp = track.GetStep()->GetPreStepPoint();
467 G4StepStatus stepStatus = sp->GetStepStatus();
468 couple = track.GetMaterialCutsCouple();
469 SetCurrentCouple(couple);
470 currentMaterialIndex = couple->GetIndex();
471 currentKinEnergy = dp->GetKineticEnergy();
472 currentRange = GetRange(particle,currentKinEnergy,couple);
473 lambda0 = GetTransportMeanFreePath(particle,currentKinEnergy);
474
475 // stop here if small range particle
476 if(inside || tPathLength < tlimitminfix) {
477 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
478 }
479
480 if(tPathLength > currentRange) { tPathLength = currentRange; }
481
482 presafety = sp->GetSafety();
483
484 // G4cout << "G4Urban2::StepLimit tPathLength= "
485 // <<tPathLength<<" safety= " << presafety
486 // << " range= " <<currentRange<< " lambda= "<<lambda0
487 // << " Alg: " << steppingAlgorithm <<G4endl;
488
489 // far from geometry boundary
490 if(currentRange < presafety)
491 {
492 inside = true;
493 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
494 }
495
496 // standard version
497 //
498 if (steppingAlgorithm == fUseDistanceToBoundary)
499 {
500 //compute geomlimit and presafety
501 geomlimit = ComputeGeomLimit(track, presafety, currentRange);
502
503 // is it far from boundary ?
504 if(currentRange < presafety)
505 {
506 inside = true;
507 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
508 }
509
510 smallstep += 1.;
511 insideskin = false;
512
513 if(firstStep || stepStatus == fGeomBoundary)
514 {
515 rangeinit = currentRange;
516 if(firstStep) smallstep = 1.e10;
517 else smallstep = 1.;
518
519 //define stepmin here (it depends on lambda!)
520 //rough estimation of lambda_elastic/lambda_transport
521 G4double rat = currentKinEnergy/MeV ;
522 rat = 1.e-3/(rat*(10.+rat)) ;
523 //stepmin ~ lambda_elastic
524 stepmin = rat*lambda0;
525 skindepth = skin*stepmin;
526 //define tlimitmin
527 tlimitmin = 10.*stepmin;
528 if(tlimitmin < tlimitminfix) tlimitmin = tlimitminfix;
529 //G4cout << "rangeinit= " << rangeinit << " stepmin= " << stepmin
530 // << " tlimitmin= " << tlimitmin << " geomlimit= " << geomlimit <<G4endl;
531 // constraint from the geometry
532 if((geomlimit < geombig) && (geomlimit > geommin))
533 {
534 // geomlimit is a geometrical step length
535 // transform it to true path length (estimation)
536 if((1.-geomlimit/lambda0) > 0.)
537 geomlimit = -lambda0*log(1.-geomlimit/lambda0)+tlimitmin ;
538
539 if(stepStatus == fGeomBoundary)
540 tgeom = geomlimit/facgeom;
541 else
542 tgeom = 2.*geomlimit/facgeom;
543 }
544 else
545 tgeom = geombig;
546 }
547
548
549 //step limit
550 tlimit = facrange*rangeinit;
551 if(tlimit < facsafety*presafety)
552 tlimit = facsafety*presafety;
553
554 //lower limit for tlimit
555 if(tlimit < tlimitmin) tlimit = tlimitmin;
556
557 if(tlimit > tgeom) tlimit = tgeom;
558
559 //G4cout << "tgeom= " << tgeom << " geomlimit= " << geomlimit
560 // << " tlimit= " << tlimit << " presafety= " << presafety << G4endl;
561
562 // shortcut
563 if((tPathLength < tlimit) && (tPathLength < presafety) &&
564 (smallstep >= skin) && (tPathLength < geomlimit-0.999*skindepth))
565 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
566
567 // step reduction near to boundary
568 if(smallstep < skin)
569 {
570 tlimit = stepmin;
571 insideskin = true;
572 }
573 else if(geomlimit < geombig)
574 {
575 if(geomlimit > skindepth)
576 {
577 if(tlimit > geomlimit-0.999*skindepth)
578 tlimit = geomlimit-0.999*skindepth;
579 }
580 else
581 {
582 insideskin = true;
583 if(tlimit > stepmin) tlimit = stepmin;
584 }
585 }
586
587 if(tlimit < stepmin) tlimit = stepmin;
588
589 // randomize 1st step or 1st 'normal' step in volume
590 if(firstStep || ((smallstep == skin) && !insideskin))
591 {
592 G4double temptlimit = tlimit;
593 if(temptlimit > tlimitmin)
594 {
595 do {
596 temptlimit = G4RandGauss::shoot(tlimit,0.3*tlimit);
597 } while ((temptlimit < tlimitmin) ||
598 (temptlimit > 2.*tlimit-tlimitmin));
599 }
600 else
601 temptlimit = tlimitmin;
602 if(tPathLength > temptlimit) tPathLength = temptlimit;
603 }
604 else
605 {
606 if(tPathLength > tlimit) tPathLength = tlimit ;
607 }
608
609 }
610 // for 'normal' simulation with or without magnetic field
611 // there no small step/single scattering at boundaries
612 else if(steppingAlgorithm == fUseSafety)
613 {
614 // compute presafety again if presafety <= 0 and no boundary
615 // i.e. when it is needed for optimization purposes
616 if((stepStatus != fGeomBoundary) && (presafety < tlimitminfix))
617 presafety = ComputeSafety(sp->GetPosition(),tPathLength);
618
619 // is far from boundary
620 if(currentRange < presafety)
621 {
622 inside = true;
623 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
624 }
625
626 if(firstStep || stepStatus == fGeomBoundary)
627 {
628 rangeinit = currentRange;
629 fr = facrange;
630 // 9.1 like stepping for e+/e- only (not for muons,hadrons)
631 if(mass < masslimite)
632 {
633 if(lambda0 > currentRange)
634 rangeinit = lambda0;
635 if(lambda0 > lambdalimit)
636 fr *= 0.75+0.25*lambda0/lambdalimit;
637 }
638
639 //lower limit for tlimit
640 G4double rat = currentKinEnergy/MeV ;
641 rat = 1.e-3/(rat*(10.+rat)) ;
642 tlimitmin = 10.*lambda0*rat;
643 if(tlimitmin < tlimitminfix) tlimitmin = tlimitminfix;
644 }
645 //step limit
646 tlimit = fr*rangeinit;
647
648 if(tlimit < facsafety*presafety)
649 tlimit = facsafety*presafety;
650
651 //lower limit for tlimit
652 if(tlimit < tlimitmin) tlimit = tlimitmin;
653
654 if(tPathLength > tlimit) tPathLength = tlimit;
655
656 }
657
658 // version similar to 7.1 (needed for some experiments)
659 else
660 {
661 if (stepStatus == fGeomBoundary)
662 {
663 if (currentRange > lambda0) tlimit = facrange*currentRange;
664 else tlimit = facrange*lambda0;
665
666 if(tlimit < tlimitmin) tlimit = tlimitmin;
667 if(tPathLength > tlimit) tPathLength = tlimit;
668 }
669 }
670 //G4cout << "tPathLength= " << tPathLength
671 // << " currentMinimalStep= " << currentMinimalStep << G4endl;
672 return ConvertTrueToGeom(tPathLength, currentMinimalStep);
673}
674
675//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
676
677G4double G4UrbanMscModel93::ComputeGeomPathLength(G4double)
678{
679 firstStep = false;
680 lambdaeff = lambda0;
681 par1 = -1. ;
682 par2 = par3 = 0. ;
683
684 // do the true -> geom transformation
685 zPathLength = tPathLength;
686
687 // z = t for very small tPathLength
688 if(tPathLength < tlimitminfix) return zPathLength;
689
690 // this correction needed to run MSC with eIoni and eBrem inactivated
691 // and makes no harm for a normal run
692 if(tPathLength > currentRange)
693 tPathLength = currentRange ;
694
695 G4double tau = tPathLength/lambda0 ;
696
697 if ((tau <= tausmall) || insideskin) {
698 zPathLength = tPathLength;
699 if(zPathLength > lambda0) zPathLength = lambda0;
700 return zPathLength;
701 }
702
703 G4double zmean = tPathLength;
704 if (tPathLength < currentRange*dtrl) {
705 if(tau < taulim) zmean = tPathLength*(1.-0.5*tau) ;
706 else zmean = lambda0*(1.-exp(-tau));
707 } else if(currentKinEnergy < mass || tPathLength == currentRange) {
708 par1 = 1./currentRange ;
709 par2 = 1./(par1*lambda0) ;
710 par3 = 1.+par2 ;
711 if(tPathLength < currentRange)
712 zmean = (1.-exp(par3*log(1.-tPathLength/currentRange)))/(par1*par3) ;
713 else
714 zmean = 1./(par1*par3) ;
715 } else {
716 G4double T1 = GetEnergy(particle,currentRange-tPathLength,couple);
717 G4double lambda1 = GetTransportMeanFreePath(particle,T1);
718
719 par1 = (lambda0-lambda1)/(lambda0*tPathLength) ;
720 par2 = 1./(par1*lambda0) ;
721 par3 = 1.+par2 ;
722 zmean = (1.-exp(par3*log(lambda1/lambda0)))/(par1*par3) ;
723 }
724
725 zPathLength = zmean ;
726
727 // sample z
728 if(samplez)
729 {
730 const G4double ztmax = 0.99 ;
731 G4double zt = zmean/tPathLength ;
732
733 if (tPathLength > stepmin && zt < ztmax)
734 {
735 G4double u,cz1;
736 if(zt >= third)
737 {
738 G4double cz = 0.5*(3.*zt-1.)/(1.-zt) ;
739 cz1 = 1.+cz ;
740 G4double u0 = cz/cz1 ;
741 G4double grej ;
742 do {
743 u = exp(log(G4UniformRand())/cz1) ;
744 grej = exp(cz*log(u/u0))*(1.-u)/(1.-u0) ;
745 } while (grej < G4UniformRand()) ;
746 }
747 else
748 {
749 cz1 = 1./zt-1.;
750 u = 1.-exp(log(G4UniformRand())/cz1) ;
751 }
752 zPathLength = tPathLength*u ;
753 }
754 }
755
756 if(zPathLength > lambda0) zPathLength = lambda0;
757 //G4cout << "zPathLength= " << zPathLength << " lambda1= " << lambda0 << G4endl;
758 return zPathLength;
759}
760
761//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
762
763G4double G4UrbanMscModel93::ComputeTrueStepLength(G4double geomStepLength)
764{
765 // step defined other than transportation
766 if(geomStepLength == zPathLength && tPathLength <= currentRange)
767 return tPathLength;
768
769 // t = z for very small step
770 zPathLength = geomStepLength;
771 tPathLength = geomStepLength;
772 if(geomStepLength < tlimitminfix) return tPathLength;
773
774 // recalculation
775 if((geomStepLength > lambda0*tausmall) && !insideskin)
776 {
777 if(par1 < 0.)
778 tPathLength = -lambda0*log(1.-geomStepLength/lambda0) ;
779 else
780 {
781 if(par1*par3*geomStepLength < 1.)
782 tPathLength = (1.-exp(log(1.-par1*par3*geomStepLength)/par3))/par1 ;
783 else
784 tPathLength = currentRange;
785 }
786 }
787 if(tPathLength < geomStepLength) tPathLength = geomStepLength;
788
789 //G4cout << "tPathLength= " << tPathLength << " step= " << geomStepLength << G4endl;
790
791 return tPathLength;
792}
793
794//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
795
796G4double G4UrbanMscModel93::ComputeTheta0(G4double trueStepLength,
797 G4double KineticEnergy)
798{
799 // for all particles take the width of the central part
800 // from a parametrization similar to the Highland formula
801 // ( Highland formula: Particle Physics Booklet, July 2002, eq. 26.10)
802 const G4double c_highland = 13.6*MeV ;
803 G4double betacp = sqrt(currentKinEnergy*(currentKinEnergy+2.*mass)*
804 KineticEnergy*(KineticEnergy+2.*mass)/
805 ((currentKinEnergy+mass)*(KineticEnergy+mass)));
806 y = trueStepLength/currentRadLength;
807 G4double theta0 = c_highland*std::abs(charge)*sqrt(y)/betacp;
808 y = log(y);
809 // correction factor from e- scattering data
810 G4double corr = coeffth1+coeffth2*y;
811
812 theta0 *= corr ;
813
814 return theta0;
815}
816
817//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
818
819G4ThreeVector&
820G4UrbanMscModel93::SampleScattering(const G4DynamicParticle* dynParticle,
821 G4double safety)
822{
823 fDisplacement.set(0.0,0.0,0.0);
824 G4double kineticEnergy = currentKinEnergy;
825 if (tPathLength > currentRange*dtrl) {
826 kineticEnergy = GetEnergy(particle,currentRange-tPathLength,couple);
827 } else {
828 kineticEnergy -= tPathLength*GetDEDX(particle,currentKinEnergy,couple);
829 }
830 if((kineticEnergy <= eV) || (tPathLength <= tlimitminfix) ||
831 (tPathLength/tausmall < lambda0)) { return fDisplacement; }
832
833 G4double cth = SampleCosineTheta(tPathLength,kineticEnergy);
834
835 // protection against 'bad' cth values
836 if(std::fabs(cth) > 1.) { return fDisplacement; }
837
838 // extra protection agaist high energy particles backscattered
839 if(cth < 1.0 - 1000*tPathLength/lambda0 && kineticEnergy > 20*MeV) {
840 //G4cout << "Warning: large scattering E(MeV)= " << kineticEnergy
841 // << " s(mm)= " << tPathLength/mm
842 // << " 1-cosTheta= " << 1.0 - cth << G4endl;
843 // do Gaussian central scattering
844 if(kineticEnergy > GeV && cth < 0.0) {
845 G4ExceptionDescription ed;
846 ed << dynParticle->GetDefinition()->GetParticleName()
847 << " E(MeV)= " << kineticEnergy/MeV
848 << " Step(mm)= " << tPathLength/mm
849 << " in " << CurrentCouple()->GetMaterial()->GetName()
850 << " CosTheta= " << cth
851 << " is too big - the angle is resampled" << G4endl;
852 G4Exception("G4UrbanMscModel93::SampleScattering","em0004",
853 JustWarning, ed,"");
854 }
855 do {
856 cth = 1.0 + 2*log(G4UniformRand())*tPathLength/lambda0;
857 } while(cth < -1.0);
858 }
859
860 G4double sth = sqrt((1.0 - cth)*(1.0 + cth));
861 G4double phi = twopi*G4UniformRand();
862 G4double dirx = sth*cos(phi);
863 G4double diry = sth*sin(phi);
864
865 G4ThreeVector oldDirection = dynParticle->GetMomentumDirection();
866 G4ThreeVector newDirection(dirx,diry,cth);
867 newDirection.rotateUz(oldDirection);
868 fParticleChange->ProposeMomentumDirection(newDirection);
869
870 if (latDisplasment && safety > tlimitminfix) {
871
872 G4double r = SampleDisplacement();
873 /*
874 G4cout << "G4UrbanMscModel93::SampleSecondaries: e(MeV)= " << kineticEnergy
875 << " sinTheta= " << sth << " r(mm)= " << r
876 << " trueStep(mm)= " << tPathLength
877 << " geomStep(mm)= " << zPathLength
878 << G4endl;
879 */
880 if(r > 0.)
881 {
882 G4double latcorr = LatCorrelation();
883 if(latcorr > r) latcorr = r;
884
885 // sample direction of lateral displacement
886 // compute it from the lateral correlation
887 G4double Phi = 0.;
888 if(std::abs(r*sth) < latcorr)
889 Phi = twopi*G4UniformRand();
890 else
891 {
892 G4double psi = std::acos(latcorr/(r*sth));
893 if(G4UniformRand() < 0.5)
894 Phi = phi+psi;
895 else
896 Phi = phi-psi;
897 }
898
899 dirx = r*std::cos(Phi);
900 diry = r*std::sin(Phi);
901
902 fDisplacement.set(dirx,diry,0.0);
903 fDisplacement.rotateUz(oldDirection);
904 }
905 }
906 return fDisplacement;
907}
908
909//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
910
911G4double G4UrbanMscModel93::SampleCosineTheta(G4double trueStepLength,
912 G4double KineticEnergy)
913{
914 G4double cth = 1. ;
915 G4double tau = trueStepLength/lambda0 ;
916
917 Zeff = couple->GetMaterial()->GetTotNbOfElectPerVolume()/
918 couple->GetMaterial()->GetTotNbOfAtomsPerVolume() ;
919
920 if(Zold != Zeff)
921 UpdateCache();
922
923 if(insideskin)
924 {
925 //no scattering, single or plural scattering
926 G4double mean = trueStepLength/stepmin ;
927
928 G4int n = G4Poisson(mean);
929 if(n > 0)
930 {
931 //screening (Moliere-Bethe)
932 G4double mom2 = KineticEnergy*(2.*mass+KineticEnergy);
933 G4double beta2 = mom2/((KineticEnergy+mass)*(KineticEnergy+mass));
934 G4double ascr = scr1/mom2;
935 ascr *= 1.13+scr2/beta2;
936 G4double ascr1 = 1.+2.*ascr;
937 G4double bp1=ascr1+1.;
938 G4double bm1=ascr1-1.;
939
940 // single scattering from screened Rutherford x-section
941 G4double ct,st,phi;
942 G4double sx=0.,sy=0.,sz=0.;
943 for(G4int i=1; i<=n; i++)
944 {
945 ct = ascr1-bp1*bm1/(2.*G4UniformRand()+bm1);
946 if(ct < -1.) ct = -1.;
947 if(ct > 1.) ct = 1.;
948 st = sqrt(1.-ct*ct);
949 phi = twopi*G4UniformRand();
950 sx += st*cos(phi);
951 sy += st*sin(phi);
952 sz += ct;
953 }
954 cth = sz/sqrt(sx*sx+sy*sy+sz*sz);
955 }
956 }
957 else
958 {
959 if(trueStepLength >= currentRange*dtrl)
960 {
961 if(par1*trueStepLength < 1.)
962 tau = -par2*log(1.-par1*trueStepLength) ;
963 // for the case if ioni/brems are inactivated
964 // see the corresponding condition in ComputeGeomPathLength
965 else if(1.-KineticEnergy/currentKinEnergy > taulim)
966 tau = taubig ;
967 }
968 currentTau = tau ;
969 lambdaeff = trueStepLength/currentTau;
970 currentRadLength = couple->GetMaterial()->GetRadlen();
971
972 if (tau >= taubig) cth = -1.+2.*G4UniformRand();
973 else if (tau >= tausmall)
974 {
975 G4double xsi = 3.;
976 G4double x0 = 1.;
977 G4double a = 1., ea = 0., eaa = 1.;
978 G4double b=2.,b1=3.,bx=1.,eb1=3.,ebx=1.;
979 G4double prob = 1. , qprob = 1. ;
980 G4double xmean1 = 1., xmean2 = 0.;
981 G4double xmeanth = exp(-tau);
982 G4double x2meanth = (1.+2.*exp(-2.5*tau))/3.;
983
984 G4double relloss = 1.-KineticEnergy/currentKinEnergy;
985 if(relloss > rellossmax)
986 return SimpleScattering(xmeanth,x2meanth);
987
988 G4double theta0 = ComputeTheta0(trueStepLength,KineticEnergy);
989
990 // protection for very small angles
991 if(theta0*theta0 < tausmall) return cth;
992
993 if(theta0 > theta0max)
994 return SimpleScattering(xmeanth,x2meanth);
995 G4double sth = sin(0.5*theta0);
996 a = 0.25/(sth*sth);
997
998 ea = exp(-xsi);
999 eaa = 1.-ea ;
1000 xmean1 = 1.-(1.-(1.+xsi)*ea)/(a*eaa);
1001 x0 = 1.-xsi/a;
1002
1003 if(xmean1 <= 0.999*xmeanth)
1004 return SimpleScattering(xmeanth,x2meanth);
1005
1006 // from e- and muon scattering data
1007 G4double c = coeffc1+coeffc2*y; ;
1008
1009 if(abs(c-3.) < 0.001) c = 3.001;
1010 if(abs(c-2.) < 0.001) c = 2.001;
1011 if(abs(c-1.) < 0.001) c = 1.001;
1012
1013 G4double c1 = c-1.;
1014
1015 //from continuity of derivatives
1016 b = 1.+(c-xsi)/a;
1017
1018 b1 = b+1.;
1019 bx = c/a;
1020 eb1 = exp(c1*log(b1));
1021 ebx = exp(c1*log(bx));
1022
1023 xmean2 = (x0*eb1+ebx-(eb1*bx-b1*ebx)/(c-2.))/(eb1-ebx);
1024
1025 G4double f1x0 = a*ea/eaa;
1026 G4double f2x0 = c1*eb1/(bx*(eb1-ebx));
1027 prob = f2x0/(f1x0+f2x0);
1028
1029 qprob = xmeanth/(prob*xmean1+(1.-prob)*xmean2);
1030
1031 // sampling of costheta
1032 if(G4UniformRand() < qprob)
1033 {
1034 if(G4UniformRand() < prob)
1035 cth = 1.+log(ea+G4UniformRand()*eaa)/a ;
1036 else
1037 cth = b-b1*bx/exp(log(ebx+(eb1-ebx)*G4UniformRand())/c1) ;
1038 }
1039 else
1040 cth = -1.+2.*G4UniformRand();
1041 }
1042 }
1043 return cth ;
1044}
1045
1046//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1047
1048G4double G4UrbanMscModel93::SimpleScattering(G4double xmeanth,G4double x2meanth)
1049{
1050 // 'large angle scattering'
1051 // 2 model functions with correct xmean and x2mean
1052 G4double a = (2.*xmeanth+9.*x2meanth-3.)/(2.*xmeanth-3.*x2meanth+1.);
1053 G4double prob = (a+2.)*xmeanth/a;
1054
1055 // sampling
1056 G4double cth = 1.;
1057 if(G4UniformRand() < prob)
1058 cth = -1.+2.*exp(log(G4UniformRand())/(a+1.));
1059 else
1060 cth = -1.+2.*G4UniformRand();
1061 return cth;
1062}
1063
1064//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1065
1066G4double G4UrbanMscModel93::SampleDisplacement()
1067{
1068 // compute rmean = sqrt(<r**2>) from theory
1069 const G4double kappa = 2.5;
1070 const G4double kappapl1 = kappa+1.;
1071 const G4double kappami1 = kappa-1.;
1072 // Compute rmean = sqrt(<r**2>) from theory
1073 G4double rmean = 0.0;
1074 if ((currentTau >= tausmall) && !insideskin) {
1075 if (currentTau < taulim) {
1076 rmean = kappa*currentTau*currentTau*currentTau*
1077 (1.-kappapl1*currentTau*0.25)/6. ;
1078
1079 } else {
1080 G4double etau = 0.0;
1081 if (currentTau<taubig) etau = exp(-currentTau);
1082 rmean = -kappa*currentTau;
1083 rmean = -exp(rmean)/(kappa*kappami1);
1084 rmean += currentTau-kappapl1/kappa+kappa*etau/kappami1;
1085 }
1086 if (rmean>0.) rmean = 2.*lambdaeff*sqrt(rmean/3.0);
1087 else rmean = 0.;
1088 }
1089
1090 if(rmean == 0.) return rmean;
1091
1092 // protection against z > t ...........................
1093 G4double rmax = (tPathLength-zPathLength)*(tPathLength+zPathLength);
1094 if(rmax <= 0.)
1095 rmax = 0.;
1096 else
1097 rmax = sqrt(rmax);
1098
1099 if(rmean >= rmax) return rmax;
1100
1101 return rmean;
1102 // VI comment out for the time being
1103 /*
1104 //sample r (Gaussian distribution with a mean of rmean )
1105 G4double r = 0.;
1106 G4double sigma = min(rmean,rmax-rmean);
1107 sigma /= 3.;
1108 G4double rlow = rmean-3.*sigma;
1109 G4double rhigh = rmean+3.*sigma;
1110 do {
1111 r = G4RandGauss::shoot(rmean,sigma);
1112 } while ((r < rlow) || (r > rhigh));
1113
1114 return r;
1115 */
1116}
1117
1118//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1119
1120G4double G4UrbanMscModel93::LatCorrelation()
1121{
1122 const G4double kappa = 2.5;
1123 const G4double kappami1 = kappa-1.;
1124
1125 G4double latcorr = 0.;
1126 if((currentTau >= tausmall) && !insideskin)
1127 {
1128 if(currentTau < taulim)
1129 latcorr = lambdaeff*kappa*currentTau*currentTau*
1130 (1.-(kappa+1.)*currentTau/3.)/3.;
1131 else
1132 {
1133 G4double etau = 0.;
1134 if(currentTau < taubig) etau = exp(-currentTau);
1135 latcorr = -kappa*currentTau;
1136 latcorr = exp(latcorr)/kappami1;
1137 latcorr += 1.-kappa*etau/kappami1 ;
1138 latcorr *= 2.*lambdaeff/3. ;
1139 }
1140 }
1141
1142 return latcorr;
1143}
1144
1145//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
JustWarning
@ JustWarning
Definition: G4ExceptionSeverity.hh:64
fUseSafety
@ fUseSafety
Definition: G4MscStepLimitType.hh:50
fUseDistanceToBoundary
@ fUseDistanceToBoundary
Definition: G4MscStepLimitType.hh:51
G4Poisson
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:50
G4StepStatus
G4StepStatus
Definition: G4StepStatus.hh:51
fGeomBoundary
@ fGeomBoundary
Definition: G4StepStatus.hh:54
G4double
double G4double
Definition: G4Types.hh:64
G4int
int G4int
Definition: G4Types.hh:66
G4endl
#define G4endl
Definition: G4ios.hh:52
G4cout
G4DLLIMPORT std::ostream G4cout
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:53
CLHEP::Hep3Vector::set
void set(double x, double y, double z)
CLHEP::Hep3Vector::rotateUz
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4LossTableManager::Instance
static G4LossTableManager * Instance()
Definition: G4LossTableManager.cc:107
G4MaterialCutsCouple::GetMaterial
const G4Material * GetMaterial() const
Definition: G4MaterialCutsCouple.hh:150
G4Material::GetTotNbOfAtomsPerVolume
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:208
G4Material::GetTotNbOfElectPerVolume
G4double GetTotNbOfElectPerVolume() const
Definition: G4Material.hh:211
G4Material::GetRadlen
G4double GetRadlen() const
Definition: G4Material.hh:219
G4Material::GetName
const G4String & GetName() const
Definition: G4Material.hh:177
G4ParticleChangeForMSC::ProposeMomentumDirection
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:115
G4Step::GetPreStepPoint
G4StepPoint * GetPreStepPoint() const
G4Track::GetDynamicParticle
const G4DynamicParticle * GetDynamicParticle() const
G4Track::GetMaterialCutsCouple
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4Track::GetStep
const G4Step * GetStep() const
G4UrbanMscModel93::Initialise
void Initialise(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4UrbanMscModel93.cc:187
G4UrbanMscModel93::SampleScattering
G4ThreeVector & SampleScattering(const G4DynamicParticle *, G4double safety)
Definition: G4UrbanMscModel93.cc:820
G4UrbanMscModel93::StartTracking
void StartTracking(G4Track *)
Definition: G4UrbanMscModel93.cc:449
G4UrbanMscModel93::ComputeGeomPathLength
G4double ComputeGeomPathLength(G4double truePathLength)
Definition: G4UrbanMscModel93.cc:677
G4UrbanMscModel93::G4UrbanMscModel93
G4UrbanMscModel93(const G4String &nam="UrbanMsc93")
Definition: G4UrbanMscModel93.cc:119
G4UrbanMscModel93::ComputeCrossSectionPerAtom
G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *particle, G4double KineticEnergy, G4double AtomicNumber, G4double AtomicWeight=0., G4double cut=0., G4double emax=DBL_MAX)
Definition: G4UrbanMscModel93.cc:207
G4UrbanMscModel93::ComputeTruePathLengthLimit
G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &currentMinimalStep)
Definition: G4UrbanMscModel93.cc:460
G4UrbanMscModel93::ComputeTheta0
G4double ComputeTheta0(G4double truePathLength, G4double KineticEnergy)
Definition: G4UrbanMscModel93.cc:796
G4UrbanMscModel93::ComputeTrueStepLength
G4double ComputeTrueStepLength(G4double geomStepLength)
Definition: G4UrbanMscModel93.cc:763
G4VEmModel::SetCurrentCouple
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:370
G4VEmModel::CurrentCouple
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:377
G4VMscModel::dtrl
G4double dtrl
Definition: G4VMscModel.hh:180
G4VMscModel::GetDEDX
G4double GetDEDX(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:273
G4VMscModel::facrange
G4double facrange
Definition: G4VMscModel.hh:176
G4VMscModel::ComputeGeomLimit
G4double ComputeGeomLimit(const G4Track &, G4double &presafety, G4double limit)
Definition: G4VMscModel.hh:256
G4VMscModel::samplez
G4bool samplez
Definition: G4VMscModel.hh:188
G4VMscModel::skin
G4double skin
Definition: G4VMscModel.hh:179
G4VMscModel::GetTransportMeanFreePath
G4double GetTransportMeanFreePath(const G4ParticleDefinition *part, G4double kinEnergy)
Definition: G4VMscModel.hh:332
G4VMscModel::GetEnergy
G4double GetEnergy(const G4ParticleDefinition *part, G4double range, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:304
G4VMscModel::GetRange
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:288
G4VMscModel::ComputeSafety
G4double ComputeSafety(const G4ThreeVector &position, G4double limit)
Definition: G4VMscModel.hh:238
G4VMscModel::steppingAlgorithm
G4MscStepLimitType steppingAlgorithm
Definition: G4VMscModel.hh:186
G4VMscModel::GetParticleChangeForMSC
G4ParticleChangeForMSC * GetParticleChangeForMSC(const G4ParticleDefinition *p=0)
Definition: G4VMscModel.cc:89
G4VMscModel::ConvertTrueToGeom
G4double ConvertTrueToGeom(G4double &tLength, G4double &gLength)
Definition: G4VMscModel.hh:246
G4VMscModel::latDisplasment
G4bool latDisplasment
Definition: G4VMscModel.hh:189
G4VMscModel::facsafety
G4double facsafety
Definition: G4VMscModel.hh:178
G4VMscModel::fDisplacement
G4ThreeVector fDisplacement
Definition: G4VMscModel.hh:185
G4VMscModel::facgeom
G4double facgeom
Definition: G4VMscModel.hh:177
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4ExceptionDescription
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76