Multiple Scattering
Elastic scattering of electrons and other charged particles is an important component of any transport code. Elastic cross section is huge when particle energy decreases, so multiple scattering (MSC) approach should be introduced in order to have acceptable CPU performance of the simulation. A universal interface G4VMultipleScattering is used by all Geant4 MSC processes [eal09, VNIU10]:
G4eMultipleScattering;
G4hMultipleScattering;
G4MuMultipleScattering.
For concrete simulation the G4VMscModel interface is used, which is an extension of the base G4VEmModel interface. The following models are available:
G4UrbanMscModel - it is applicable to all types of particles and is the default model for electrons and positrons below 100 MeV [eal17];
G4GoudsmitSaundersonModel - for electrons and positrons [OKT09], this model was recently rewritten [eal17] and provides the best accuracy of electron transport below 100 MeV [SIN18];
G4WentzelVIModel - is the default model for all changed particles [eal17, eal16, VNIU10] including electrons and positrons above 100 MeV, it is included in Physics List together with G4CoulombScattering process, which is responsible for large angle scattering;
G4LowEWentzelVIModel - for all particles with low-energy limit 10 eV - is a less accurate variant of the above mode.
The discussion on Geant4 MSC models is available in Ref. [VNIU10]. Below we will describe models developed by L. Urban [Urb06], because these models are used in many Geant4 applications and have general components reused by other models.
Introduction
MSC simulation algorithms can be classified as either detailed or condensed. In the detailed algorithms, all the collisions/interactions experienced by the particle are simulated. This simulation can be considered as exact, it gives the same results as the solution of the transport equation. However, it can be used only if the number of collisions is not too large, a condition fulfilled only for special geometries (such as thin foils, or low density gas). In solid or liquid media the average number of collisions is very large and the detailed simulation becomes very inefficient. High energy simulation codes use condensed simulation algorithms, in which the global effects of the collisions are simulated at the end of a track segment. The global effects generally computed in these codes are the net energy loss, displacement, and change of direction of the charged particle. The last two quantities are computed from MSC theories used in the codes and the accuracy of the condensed simulations is limited by accuracy of MSC approximation.
Most particle physics simulation codes use the multiple scattering theories of Molière [Moliere48], Goudsmit and Saunderson [GS40] and Lewis [Lew50]. The theories of Molière and Goudsmit-Saunderson give only the angular distribution after a step, while the Lewis theory computes the moments of the spatial distribution as well. None of these MSC theories gives the probability distribution of the spatial displacement. Each of the MSC simulation codes incorporates its own algorithm to determine the angular deflection, true path length correction, and spatial displacement of the charged particle after a given step. These algorithms are not exact, of course, and are responsible for most of the uncertainties of the transport codes. Also due to inaccuracy of MSC the simulation results can depend on the value of the step length and generally user has to select the value of the step length carefully.
A new class of MSC simulation, the mixed simulation algorithms (see e.g.[JMFernandezVS93]), appeared in the literature recently. The mixed algorithm simulates the hard collisions one by one and uses a MSC theory to treat the effects of the soft collisions at the end of a given step. Such algorithms can prevent the number of steps from becoming too large and also reduce the dependence on the step length. Geant4 original implementation of a similar approach is realized in G4WentzelVIModel [VNIU10].
The Urban MSC models used in Geant4 belongs to the class of condensed simulations. Urban uses model functions to determine the angular and spatial distributions after a step. The functions have been chosen in such a way as to give the same moments of the (angular and spatial) distributions as are given by the Lewis theory [Lew50].
Definition of Terms
In simulation, a particle is transported by steps through the detector geometry. The shortest distance between the endpoints of a step is called the geometrical path length, \(z\). In the absence of a magnetic field, this is a straight line. For non-zero fields, \(z\) is the length along a curved trajectory. Constraints on \(z\) are imposed when particle tracks cross volume boundaries. The path length of an actual particle, however, is usually longer than the geometrical path length, due to multiple scattering. This distance is called the true path length, \(t\). Constraints on \(t\) are imposed by the physical processes acting on the particle.
The properties of the MSC process are determined by the transport mean free paths, \(\lambda_k\), which are functions of the energy in a given material. The \(k\)-th transport mean free path is defined as
where \(d\sigma(\chi)/d\Omega\) is the differential cross section of the scattering, \(P_k(\cos\chi)\) is the \(k\)-th Legendre polynomial, and \(n_a\) is the number of atoms per volume.
Most of the mean properties of MSC computed in the simulation codes depend only on the first and second transport mean free paths. The mean value of the geometrical path length (first moment) corresponding to a given true path length \(t\) is given by
Eq.(47) is an exact result for the mean value of \(z\) if the differential cross section has axial symmetry and the energy loss can be neglected. The transformation between true and geometrical path lengths is called the path length correction. This formula and other expressions for the first moments of the spatial distribution were taken from either [JMFernandezVS93] or [KB98], but were originally calculated by Goudsmit and Saunderson [GS40] and Lewis [Lew50].
At the end of the true step length, \(t\), the scattering angle is \(\theta\). The mean value of \(\cos\theta\) is
The variance of \(\cos\theta\) can be written as
where \(\tau = t/\lambda_1\) and \(\kappa = \lambda_1/\lambda_2\). The mean lateral displacement is given by a more complicated formula [JMFernandezVS93], but this quantity can also be calculated relatively easily and accurately. The square of the mean lateral displacement is
Here it is assumed that the initial particle direction is parallel to the the \(z\) axis. The lateral correlation is determined by the equation
where \(v_x\) and \(v_y\) are the x and y components of the direction unit vector. This equation gives the correlation strength between the final lateral position and final direction.
The transport mean free path values have been calculated in Refs. [LI87], [eal90] for electrons and positrons in the kinetic energy range in 15 materials. The Urban MSC model in Geant4 uses these values for kinetic energies below 10 MeV. For high energy particles (above 10 MeV) the transport mean free path values have been taken from a paper of R. Mayol and F. Salvat [MS97]. When necessary, the model linearly interpolates or extrapolates the transport cross section, \(\sigma_1 = 1 / \lambda_1\), in atomic number \(Z\) and in the square of the particle velocity, \(\beta^2\). The ratio \(\kappa\) is a very slowly varying function of the energy: \(\kappa > 2\) for \(T >\) a few keV, and \(\kappa \rightarrow 3\) for very high energies (see [KB98]). Hence, a constant value of 2.5 is used in the model.
Nuclear size effects are negligible for low energy particles and they are accounted for in the Born approximation in [MS97], so there is no need for extra corrections of this kind in the Urban model.
Path Length Correction
As mentioned above, the path length correction refers to the transformation \(t \longrightarrow g\) and its inverse. The \(t \longrightarrow g\) transformation is given by Eq.(47) if the step is small and the energy loss can be neglected. If the step is not small the energy dependence makes the transformation more complicated. For this case Eqs.(48),(47) should be modified as
where \(\theta\) is the scattering angle, \(t\) and \(z\) are the true and geometrical path lengths, and \(\lambda_1\) is the transport mean free path.
In order to compute Eqs.(52),(53) the \(t\) dependence of the transport mean free path must be known. \(\lambda_1\) depends on the kinetic energy of the particle which decreases along the step. All computations in the model use a linear approximation for this \(t\) dependence:
Here \(\lambda_{10}\) denotes the value of \(\lambda_1\) at the start of the step, and \(\alpha\) is a constant. It is worth noting that Eq.(54) is not a crude approximation. It is rather good at low (\(<\) 1 MeV) energy. At higher energies the step is generally much smaller than the range of the particle, so the change in energy is small and so is the change in \(\lambda_1\). Using Eqs.(52) - (54) the explicit formula for \(\langle \cos\theta \rangle\) and \(\langle z \rangle\) are:
The value of the constant \(\alpha\) can be expressed using \(\lambda_{10}\) and \(\lambda_{11}\) where \(\lambda_{11}\) is the value of the transport mean free path at the end of the step
At low energies ( \(T_{kin} < M\) , where \(M\) is the particle mass) \(\alpha\) has a simpler form:
where \(r_0\) denotes the range of the particle at the start of the step. It can easily be seen that for a small step (i.e. for a step with small relative energy loss) the formula of \(\langle z \rangle\) is
Eq. (56) or (57) gives the mean value of the geometrical step length for a given true step length. The actual geometrical path length is sampled in the model according to the simple probability density function defined for \(v = z/t \in [0 , 1]\) :
The value of the exponent \(k\) is computed from the requirement that \(f(v)\) must give the same mean value for \(z = v t\) as Eq. (56) or (57). Hence
The value of \(z = v t\) is sampled using \(f(v)\) if \(k > 0\), otherwise \(z = \langle z \rangle\) is used. The \(g \longrightarrow t\) transformation is performed using the mean values. The transformation can be written as
if the geometrical step is small and
where
if the step is not small, i.e.the energy loss should be taken into account.
Angular Distribution
The quantity \(u = \cos\theta\) is sampled according to a model function \(g(u)\). The shape of this function has been chosen such that Eqs. (48) and (49) are satisfied. The functional form of \(g\) is
where \(0 \leq p,q \leq 1\), and the \(g_i\) are simple functions of \(u = \cos\theta\), normalized over the range \(u \in [-1,\ 1]\). The functions \(g_i\) have been chosen as
where \(a > 0\), \(b > 0\), \(d > 0\) and \(u_0\) are model parameters, and the \(C_{i}\) are normalization constants. It is worth noting that for small scattering angles, \(\theta\), \(g_1(u)\) is nearly Gaussian (\(\exp(-\theta^2/2 \theta_0^2)\)) if \(\theta_0^2 \approx 1 / a\), while \(g_2(u)\) has a Rutherford-like tail for large \(\theta\), if \(b \approx 1\) and \(d\) is not far from 2 .
Determination of the Model Parameters
The parameters \(a\), \(b\), \(d\), \(u_0\) and \(p\), \(q\) are not independent. The requirement that the angular distribution function \(g(u)\) and its first derivative be continuous at \(u = u_0\) imposes two constraints on the parameters:
A third constraint comes from Eq. (52) : \(g(u)\) must give the same mean value for \(u\) as the theory. It follows from Eqs. (55) and (58) that
where \(\langle u \rangle_i\) denotes the mean value of \(u\) computed from the distribution \(g_i(u)\). The parameter \(a\) was chosen according to a modified Highland-Lynch-Dahl formula for the width of the angular distribution [Hig75], [LD91].
where \(\theta_0\) is
when the original Highland-Lynch-Dahl formula is used. Here \(\theta_0 = \theta^{rms}_{plane}\) is the width of the approximate Gaussian projected angle distribution, \(p\), \(\beta c\) and \(z_{ch}\) are the momentum, velocity and charge number of the incident particle, and \(t/X_0\) is the true path length in radiation length unit. The correction term \(h_c\) = 0.038 in the formula. This value of \(\theta_0\) is from a fit to the Molière distribution for singly charged particles with \(\beta = 1\) for all Z, and is accurate to 11% or better for \(10^{-3} \leq t/X_0 \leq 100\) (see e.g. Rev. of Particle Properties, section 23.3).
The model uses a slightly modified Highland-Lynch-Dahl formula to compute \(\theta_0\). For electrons/positrons the modified \(\theta_0\) formula is
where
The correction term \(c\) and coefficients \(c_i\) are
This formula gives a much smaller step dependence in the angular distribution than the Highland form. The value of the parameter \(u_0\) has been chosen as
where
with
The parameters \(d_i\)-s have the form
The numerical values of the \(d_{ij}\) constants can be found in the code.
The tail parameter \(d\) is the same as the parameter \(\xi\) .
This (empirical) expression is obtained comparing the simulation results to the data of the MuScat experiment [eal06]. The remaining three parameters can be computed from Eqs. (59) - (60). The numerical value of the parameters can be found in the code.
In the case of heavy charged particles (\(\mu\), \(\pi\), \(p\), etc.) the mean transport free path is calculated from the electron or positron \(\lambda_1\) values with a ’scaling’ applied. This is possible because the transport mean free path \(\lambda_1\) depends only on the variable \(P \beta c\), where \(P\) is the momentum, and \(\beta c\) is the velocity of the particle.
In its present form the model samples the path length correction and angular distribution from model functions, while for the lateral displacement and the lateral correlation only the mean values are used and all the other correlations are neglected. However, the model is general enough to incorporate other random quantities and correlations in the future.
Step Limitation Algorithm
In Geant4 the boundary crossing is treated by the transportation process. The transportation ensures that the particle does not penetrate in a new volume without stopping at the boundary, it restricts the step size when the particle leaves a volume. However, this step restriction can be rather weak in big volumes and this fact can result a not very good angular distribution after the volume. At the same time, there is no similar step limitation when a particle enters a volume and this fact does not allow a good backscattering simulation for low energy particles. Low energy particles penetrate too deeply into the volume in the first step and then, because of energy loss, they are not able to reach again the boundary in backward direction.
MSC step limitation algorithm has been developed [Urb06] in order to achieve optimal balance between simulation precision and CPU performance of simulation for different applications. At the start of a track or after entering in a new volume, the algorithm restricts the step size to a value
where \(r\) is the range of the particle, \(f_r\) is a parameter \(\in [0, 1]\), taking the max of \(r\) and \(\lambda_1\) is an empirical choice. The value of \(f_r\) is constant for low energy particles while for particles with \(\lambda_1 > \lambda_{lim}\) an effective value is used given by the scaling equation
(The numerical values \(sc = 0.25\) and \(\lambda_{lim} = 1\) mm are used in the equation.) In order not to use very small - unphysical - step sizes a lower limit is given for the step size as
with \(nstepmax = 25\) and \(\lambda_{elastic}\) is the elastic mean free path of the particle (see later). It can be easily seen that this kind of step limitation poses a real constraint only for low energy particles. In order to prevent a particle from crossing a volume in just one step, an additional limitation is imposed: after entering a volume the step size cannot be bigger than
where \(d_{geom}\) is the distance to the next boundary (in the direction of the particle) and \(f_g\) is a constant parameter. A similar restriction at the start of a track is
At this point the program also checks whether the particle has entered a new volume. If it has, the particle steps cannot be bigger than \(t_{lim} = f_r\hspace{2mm} \max( r, \lambda )\). This step limitation is governed by the physics, because \(t_{lim}\) depends on the particle energy and the material.
The choice of the parameters \(f_r\) and \(f_g\) is also related to performance. By default \(f_r = 0.02\) and \(f_g = 2.5\) are used, but these may be set to any other value in a simple way. One can get an approximate simulation of the backscattering with the default value, while if a better backscattering simulation is needed it is possible to get it using a smaller value for \(f_r\). However, this model is very simple and it can only approximately reproduce the backscattering data.
Boundary Crossing Algorithm
A special stepping algorithm has been implemented in order to improve the simulation around interfaces. This algorithm does not allow ‘big’ last steps in a volume and ‘big’ first steps in the next volume. The step length of these steps around a boundary crossing can not be bigger than the mean free path of the elastic scattering of the particle in the given volume (material). After these small steps the particle scattered according to a single scattering law (i.e.there is no multiple scattering very close to the boundary or at the boundary).
The key parameter of the algorithm is the variable called \(skin\). The algorithm is not active for \(skin \leq 0\), while for \(skin > 0\) it is active in layers of thickness \(skin \cdot \lambda_{elastic}\) before boundary crossing and of thickness \((skin-1) \cdot \lambda_{elastic}\) after boundary crossing (for \(skin=1\) there is only one small step just before the boundary). In this active area the particle performs steps of length \(\lambda_{elastic}\) (or smaller if the particle reaches the boundary traversing a smaller distance than this value).
The scattering at the end of a small step is single or plural and for these small steps there are no path length correction and lateral displacement computation. In other words the program works in this thin layer in ‘microscopic mode’. The elastic mean free path can be estimated as
where \(rat(T_{kin})\) a simple empirical function computed from the elastic and first transport cross section values of Mayol and Salvat [MS97]
\(T_{kin}\) is the kinetic energy of the particle.
At the end of a small step the number of scatterings is sampled according to the Poisson’s distribution with a mean value \(t/\lambda_{elastic}\) and in the case of plural scattering the final scattering angle is computed by summing the contributions of the individual scatterings. The single scattering is determined by the distribution
where \(u = \cos(\theta)\) , \(a\) is the screening parameter, \(C\) is a normalization constant. The form of the screening parameter is the same as in the single scattering (see there).
Implementation Details
The step length of a particles is determined by the physics processes or the geometry of the detectors. The tracking/stepping algorithm checks all the step lengths demanded by the (continuous or discrete) physics processes and determines the minimum of these step lengths (see True Step Length). The MSC model should be called to compute step limit after all processes except the transportation process. The following sequence of computations are performed to make the step:
the minimum of all processes true step length limit \(t\) including one of the MSC process is selected;
The conversion \(t \longrightarrow g\) (geometrical step limit) is performed;
the minimum of obtained value \(g\) and the transportation step limit is selected;
The final conversion \(g \longrightarrow t\) is performed.
The reason for this ordering is that the physics processes ‘feel’ the true path length \(t\) traveled by the particle, while the transportation process (geometry) uses the \(z\) step length.
A new optional mechanism was recently introduced allowing sample displacement in the vicinity of geometry boundary. If it is enabled and transportation limits the step due to a geometry boundary, then after initial sampling of the displacement an additional ‘push’ of the track is applied forcing the end point be at the boundary. Corresponding correction to the true step length is applied according to the value of the ‘push’.
After the actual step of the particle is done, the MSC model is responsible for sampling of scattering angle and relocation of the end-point of the step. The scattering angle \(\theta\) of the particle after the step of length t is sampled according to the model function given in Eq. (58) . The azimuthal angle \(\phi\) is generated uniformly in the range \([0, 2 \pi]\).
After the simulation of the scattering angle, the lateral displacement is computed using Eq. (50). Then the correlation given by Eq. (51) is used to determine the direction of the lateral displacement. Before ’moving’ the particle according to the displacement a check is performed to ensure that the relocation of the particle with the lateral displacement does not take the particle beyond the volume boundary.
Default MSC parameter values optimized per particle type are shown in Table 15. Note, that there are four types of step limitation by multiple scattering process:
Minimal
- only \(f_r\) parameter and range are used;UseSafety
- \(f_r\) parameter, range and geometrical safety are used;UseSafetyPlus
- \(f_r\) parameter, range and geometrical safety are used;UseDistanceToBoundary
- uses particle range, geometrical safety and linear distance to geometrical boundary.
particle |
\(e^+\), \(e^-\) |
muons, hadrons |
ions |
---|---|---|---|
StepLimitType |
fUseSafety |
fMinimal |
fMinimal |
skin |
0 |
0 |
0 |
\(f_r\) |
0.04 |
0.2 |
0.2 |
\(f_g\) |
2.5 |
0.1 |
0.1 |
LateralDisplacement |
true |
false |
false |
The parameters of the model can be changed via public functions of the base class G4VMultipleSacttering. They can be changed for all multiple scattering processes simultaneously via G4EmParameters class, G4EmProcessOptions class, or via Geant4 UI commands. The following commands are available:
/process/msc/StepLimit UseDistanceToBoundary
/process/msc/LateralDisplacement false
/process/msc/MuHadLateralDisplacement false
/process/msc/DisplacementBeyondSafety true
/process/msc/RangeFactor 0.02
/process/msc/GeomFactor 2.5
/process/msc/Skin 2
Bibliography
- eal17(1,2,3)
A. Bagulya et al. Recent progress of geant4 electromagnetic physics for lhc and other applications. Journal of Physics: Conference Series, 898():042032, 2017. URL: https://doi.org/10.1088/1742-6596/898/4/042032, doi:10.1088/1742-6596/898/4/042032.
- eal06
D. Attwood et al. The scattering of muons in low-z materials. Nucl. Instr. and Meth. in Phys. Research B, 251(1):41–55, sep 2006. URL: https://doi.org/10.1016/j.nimb.2006.05.006, doi:10.1016/j.nimb.2006.05.006.
- eal90
D. Liljequist et al. Transport mean free path tabulated for the multiple elastic scattering of electrons and positrons at energies ≤20 mev. Journal of Applied Physics, 68(7):3061–3065, oct 1990. URL: https://doi.org/10.1063/1.346399, doi:10.1063/1.346399.
- eal16
J. Allison et al. Recent developments in geant4. Nucl. Instr. and Meth. in Phys. Research Section A, 835:186–225, nov 2016. URL: https://doi.org/10.1016/j.nima.2016.06.125, doi:10.1016/j.nima.2016.06.125.
- eal09
J. Apostolakis et al. Geometry and physics of the geant4 toolkit for high and medium energy applications. Radiation Physics and Chemistry, 78(10):859–873, oct 2009. URL: https://doi.org/10.1016/j.radphyschem.2009.04.026, doi:10.1016/j.radphyschem.2009.04.026.
- GS40(1,2)
S. Goudsmit and J. L. Saunderson. Multiple scattering of electrons. Physical Review, 57(1):24–29, jan 1940. URL: https://doi.org/10.1103/PhysRev.57.24, doi:10.1103/physrev.57.24.
- Hig75
V.L. Highland. Some practical remarks on multiple scattering. Nuclear Instruments and Methods, 129:497–499, November 1975. doi:10.1016/0029-554X(75)90743-0.
- JMFernandezVS93(1,2,3)
J. Baró J.M. Fernández-Varea, R. Mayol and F. Salvat. On the theory and simulation of multiple elastic scattering of electrons. Nucl. Instrum. and Meth. in Phys. Research B, 73:447–473, apr 1993. doi:10.1016/0168-583X(93)95827-R.
- KB98(1,2)
I. Kawrakow and A.F. Bielajew. On the condensed history technique for electron transport. Nucl. Instr. and Meth. in Phys. Research Section B, 142(3):253–280, jul 1998. URL: https://doi.org/10.1016/S0168-583X(98)00274-2, doi:10.1016/s0168-583x(98)00274-2.
- Lew50(1,2,3)
H. W. Lewis. Multiple scattering in an infinite medium. Physical Review, 78(5):526–529, jun 1950. URL: https://doi.org/10.1103/PhysRev.78.526, doi:10.1103/physrev.78.526.
- LI87
D. Liljequist and M. Ismail. Transport mean free path related to trajectory patterns: comparison of nonrelativistic and highly relativistic electron penetration through matter. Journal of Applied Physics, 62(2):342–350, jul 1987. URL: http://doi.org/10.1063/1.339802, doi:10.1063/1.339802.
- LD91
G.R. Lynch and O.I. Dahl. Approximations to multiple Coulomb scattering. Nucl. Instr. and Meth. in Phys. Research B, 58:6–10, May 1991. doi:10.1016/0168-583X(91)95671-Y.
- MS97(1,2,3)
Ricard Mayol and Francesc Salvat. Total and Transport Cross Sections for Elastic Scattering of Electrons by Atoms. Atomic Data and Nuclear Data Tables, 65(1):55–154, jan 1997. URL: https://doi.org/10.1006/adnd.1997.0734, doi:10.1006/adnd.1997.0734.
- Moliere48
G. Molière. Theorie der Streuung schneller geladener Teilchen II. Mehrfach- und Vielfachstreuung. Zeitschrift Naturforschung Teil A, 3:78–97, February 1948. doi:10.1515/zna-1948-0203.
- OKT09
F. Gharbi O. Kadri, V. Ivanchenko and A. Trabelsi. Incorporation of the goudsmit–saunderson electron transport theory in the geant4 monte carlo code. Nucl. Instr. and Meth. in Phys. Research Section B, 267(23-24):3624–3632, dec 2009. URL: https://doi.org/10.1016/j.nimb.2009.09.015, doi:10.1016/j.nimb.2009.09.015.
- SIN18
V Ivanchenko S Incerti and M Novak. Recent progress of geant4 electromagnetic physics for calorimeter simulation. Journal of Instrumentation, 13():C02054, feb 2018. URL: https://doi.org/10.1088/1748-0221/13/02/C02054, doi:10.1088/1748-0221/13/02/C02054.
- Urb06(1,2)
László Urbán. A model for multiple scattering in GEANT4. Technical Report, CERN, Geneva, Dec 2006. URL: http://cds.cern.ch/record/1004190.
- VNIU10(1,2,3,4)
M. Maire V.N. Ivanchenko, O. Kadri and L. Urban. Geant4 models for simulation of multiple scattering. Journal of Physics: Conference Series, 219(3):032045, apr 2010. URL: https://doi.org/10.1088/1742-6596/219/3/032045, doi:10.1088/1742-6596/219/3/032045.