Method

The G4eIonisation class provides the continuous and discrete energy losses of electrons and positrons due to ionisation in a material according to the approach described in Mean Energy Loss. The value of the maximum energy transferable to a free electron \(T_{max}\) is given by the following relation:

(123)\[\begin{split}T_{max} = \left\{ \begin{array}{ll} E-mc^2 & {\mbox{for} \hspace{.2cm} e^+} \\ (E-mc^2)/2 & {\mbox{for} \hspace{.2cm} e^- } \\ \end{array} \right .\end{split}\]

where \(mc^2\) is the electron mass. Above a given threshold energy the energy loss is simulated by the explicit production of delta rays by Möller scattering (\(e^- e^-\)), or Bhabha scattering (\(e^+ e^-\)). Below the threshold the soft electrons ejected are simulated as continuous energy loss by the incident \({e^{\pm}}\).

Continuous Energy Loss

The integration of (29) leads to the Berger-Seltzer formula [MC70]:

(124)\[\left. \frac{dE}{dx} \right]_{T < T_{cut}} = 2 \pi r_e^2 mc^2 n_{el} \frac{1}{\beta^2} \left [\ln \frac{2(\gamma + 1)} {(I/mc^2)^2}+ F^{\pm} (\tau , \tau_{up}) - \delta \right ]\]

with

\[\begin{split}r_e &= \mbox{classical electron radius:} \quad e^2/(4 \pi \epsilon_0 mc^2 ) \\ mc^2 &= \mbox{mass energy of the electron} \\ n_{el} &= \mbox{electron density in the material} \\ I &= \mbox{mean excitation energy in the material}\\ \gamma &= E/mc^2 \\ \beta^2 &= 1-(1/\gamma^2) \\ \tau &= \gamma-1 \\ T_{cut} &= \mbox{minimum energy cut for } \delta \mbox{-ray production} \\ \tau_c &= T_{cut}/mc^2 \\ \tau_{max} &= \mbox{maximum energy transfer:} \tau \mbox{ for } e^+, \tau/2 \mbox{ for } e^- \\ \tau_{up} &= \min(\tau_c,\tau_{max}) \\ \delta &= \mbox{density effect function} .\end{split}\]

In an elemental material the electron density is

\[n_{el} = Z \: n_{at} = Z \: \frac{\mathcal{N}_{av} \rho}{A} .\]

\(\mathcal{N}_{av}\) is Avogadro’s number, \(\rho\) is the material density, and \(A\) is the mass of a mole. In a compound material

\[n_{el} = \sum_i Z_i \: n_{ati} = \sum_i Z_i \: \frac{\mathcal{N}_{av} w_i \rho}{A_i} ,\]

where \(w_i\) is the proportion by mass of the \(i^{th}\) element, with molar mass \(A_i\) .

The mean excitation energies \(I\) for all elements are taken from [eal84].

The functions \(F^{\pm}\) are given by :

\[\begin{split}F^+ (\tau,\tau_{up}) &= \ln(\tau\tau_{up} ) -\frac{\tau_{up}^2}{\tau}\left[\tau + 2 \tau_{up} - \frac{3\tau_{up}^2 y } {2} -\left(\tau_{up} - \frac{\tau_{up}^3 }{3} \right) y^2 - \left (\frac{\tau_{up}^2}{2} - \tau \frac{\tau_{up}^3}{3} + \frac{\tau_{up}^4 } {4} \right) y^3 \right] \\ F^- (\tau,\tau_{up} ) &= -1 -\beta^2 +\ln \left [(\tau - \tau_{up}) \tau_{up} \right ] + \frac{\tau}{\tau -\tau_{up}} + \frac{1}{\gamma^2} \left [ \frac{\tau_{up}^2}{2} + ( 2\tau +1) \ln \left (1- \frac{\tau_{up}}{\tau} \right ) \right ]\end{split}\]

where \(y = 1/(\gamma+1)\).

The density effect correction is calculated according to the formalism of Sternheimer [SP71b]:

\(x\) is a kinetic variable of the particle : \(x = \log_{10}(\gamma \beta) = \ln(\gamma^{2} \beta^{2})/4.606\), and \(\delta(x)\) is defined by

\[\begin{split}\begin{array}{rll} \mbox{for} & x < x_0 : & \delta(x) = 0 \\ \mbox{for} & x \in [x_0,\ x_1] : & \delta(x) = 4.606 x - C + a(x_1 - x)^m \\ \mbox{for} & x > x_1 : & \delta(x) = 4.606 x - C \end{array}\end{split}\]

where the matter-dependent constants are calculated as follows:

\[\begin{split}\begin{array}{lcl} h\nu_p & = & \mbox{ plasma energy of the medium } = \sqrt{4\pi n_{el} r_e^3} mc^2/\alpha = \sqrt{4\pi n_{el} r_e} \hbar c \\ C & = & 1 + 2 \ln (I/h\nu_p) \\ x_a & = & C/4.606 \\ a & = & 4.606(x_a - x_0)/(x_1 - x_0)^m \\ m & = & 3 . \end{array}\end{split}\]

For condensed media

\[\begin{split}\begin{array}{ll} I < 100 \: \mbox{eV} & \left \{ \begin{array}{rll} \mbox{for } C \leq 3.681 & x_0 = 0.2 & x_1 = 2 \\ \mbox{for } C > 3.681 & x_0 = 0.326 C - 1.0 & x_1 = 2 \end{array} \right . \\ I \geq 100 \: \mbox{eV} & \left \{ \begin{array}{rll} \mbox{for } C \leq 5.215 & x_0 = 0.2 & x_1 = 3 \\ \mbox{for } C > 5.215 & x_0 = 0.326 C - 1.5 & x_1 = 3 \end{array} \right . \end{array}\end{split}\]

and for gaseous media

\[\begin{split}\begin{array}{rlll} \mbox{for} & C < 10. & x_0 = 1.6 & x_1 = 4 \\ \mbox{for} & C \in [10.0,\ 10.5[ & x_0 = 1.7 & x_1 = 4 \\ \mbox{for} & C \in [10.5,\ 11.0[ & x_0 = 1.8 & x_1 = 4 \\ \mbox{for} & C \in [11.0,\ 11.5[ & x_0 = 1.9 & x_1 = 4 \\ \mbox{for} & C \in [11.5,\ 12.25[ & x_0 = 2. & x_1 = 4 \\ \mbox{for} & C \in [12.25,\ 13.804[ & x_0 = 2. & x_1 = 5 \\ \mbox{for} & C \geq 13.804 & x_0 = 0.326 C -2.5 & x_1 = 5 . \end{array}\end{split}\]

Total Cross Section per Atom and Mean Free Path

The total cross section per atom for Möller scattering (\(e^- e^-\)) and Bhabha scattering (\(e^+ e^-\)) is obtained by integrating Eq. (30). In Geant4 \(T_{cut}\) is always 1 keV or larger. For delta ray energies much larger than the excitation energy of the material (\(T \gg I\)), the total cross section becomes [MC70] for Möller scattering,

\[\sigma ( Z,E,T_{cut} ) = \frac {2 \pi r_e^2 Z}{\beta^2(\gamma -1)} \left[\frac{(\gamma-1)^2} {\gamma^2}\left(\frac{1}{2}-x\right) +\frac{1}{x}-\frac{1}{1-x}-\frac{2\gamma-1}{\gamma^2}\ln \frac{1-x}{x}\right] ,\]

and for Bhabha scattering (\(e^+ e^-\)),

\[\sigma (Z,E,T_{cut}) = \frac{ 2 \pi r_e^2 Z }{(\gamma -1)} \left [\frac {1 }{\beta^2} \left(\frac{1}{x}-1\right) + B_1 \ln x + B_2 (1-x) - \frac {B_3 } {2} ( 1-x^2 ) +\frac{B_4}{3}(1-x^3)\right] .\]

Here

\[\begin{split}\begin{array}{lcllcl} \gamma & = & E/mc^2 & B_1 & = & 2-y^2 \\ \beta^2 & = & 1-(1/\gamma^2) & B_2 & = & (1-2y)(3+y^2 ) \\ x & = & T_{cut}/(E-mc^2) & B_3 & = & (1-2y)^2+(1-2y)^3 \\ y & = & 1/(\gamma + 1) & B_4 & = & (1-2y)^3 . \end{array}\end{split}\]

The above formulas give the total cross section for scattering above the threshold energies

\[T_{\rm Moller}^{\rm thr} =2T_{cut} \hspace{2cm} \mbox{and} \hspace{2cm} T_{\rm Bhabha}^{\rm thr} = T_{cut} .\]

In a given material the mean free path is then

\[\begin{array}{lll} \lambda = (n_{at} \cdot \sigma)^{-1} & \mbox{or} & \lambda = \left( \sum_i n_{ati} \cdot \sigma_i \right)^{-1} . \end{array}\]

Simulation of Delta-ray Production

Differential Cross Section

For \(T \gg I\) the differential cross section per atom becomes [MC70] for Möller scattering,

(125)\[\frac{d\sigma }{d \epsilon } = \frac{2 \pi r_e^2 Z}{\beta^2 (\gamma -1)} \left[ \frac{(\gamma -1 )^2} {\gamma^2 }+\frac{1}{\epsilon} \left(\frac{1}{\epsilon}-\frac{2 \gamma -1 } {\gamma^2 } \right) + \frac{1}{1- \epsilon}\left(\frac{1} {1- \epsilon} - \frac{2 \gamma - 1} {\gamma^2 }\right) \right]\]

and for Bhabha scattering,

(126)\[\frac{d \sigma}{d \epsilon}=\frac{2 \pi r_e^2 Z}{(\gamma -1)}\left[ \frac{1} {\beta^2 \epsilon^2}-\frac{B_1}{\epsilon}+B_2 - B_3 \epsilon + B_4 \epsilon^2\right] .\]

Here \(\epsilon = T/(E-mc^2)\). The kinematical limits of \(\epsilon\) are

\[\epsilon_0 = \frac{T_{cut}}{E-mc^2} \leq \epsilon \leq \frac{1}{2} \hspace{.2cm} \mbox{ for $e^- e^-$} \hspace{2cm} \epsilon_0 = \frac{T_{cut}}{E-mc^2} \leq \epsilon \leq 1 \hspace{.2cm} \mbox{ for $e^+ e^-$} .\]

Sampling

The delta ray energy is sampled according to methods discussed in Monte Carlo Methods. Apart from normalization, the cross section can be factorized as

\[\frac{d\sigma}{d\epsilon}=f(\epsilon) g(\epsilon) .\]

For \(e^- e^-\) scattering

\[\begin{split}f(\epsilon) &= \frac{1}{\epsilon^2} \frac{\epsilon_0 }{1- 2\epsilon_0} \\ g(\epsilon) &= \frac{4}{9\gamma^2 - 10 \gamma + 5}\left[(\gamma -1)^2 \epsilon^2 - (2 \gamma^2 +2\gamma -1) \frac{\epsilon} {1- \epsilon }+ \frac{\gamma^2}{(1- \epsilon )^2 }\right]\end{split}\]

and for \(e^+ e^-\) scattering

\[\begin{split}f(\epsilon) &= \frac{1}{\epsilon^2} \frac{\epsilon_0}{1- \epsilon_0 } \\ g(\epsilon) &= \frac{B_0 -B_1 \epsilon +B_2 \epsilon^2 -B_3 \epsilon^3 +B_4 \epsilon ^4}{B_ 0-B_1\epsilon_0 +B_2\epsilon_0^2 -B_3 \epsilon_0^3 +B_4 \epsilon_0^4} .\end{split}\]

Here \(B_0=\gamma^2/(\gamma^2-1)\) and all other quantities have been defined above.

To choose \(\epsilon\), and hence the delta ray energy,

  1. \(\epsilon\) is sampled from \(f(\epsilon)\)

  2. the rejection function \(g(\epsilon)\) is calculated using the sampled value of \(\epsilon\)

  3. \(\epsilon\) is accepted with probability \(g(\epsilon)\).

After the successful sampling of \(\epsilon\), the direction of the ejected electron is generated with respect to the direction of the incident particle. The azimuthal angle \(\phi\) is generated isotropically and the polar angle \(\theta\) is calculated from energy-momentum conservation. This information is used to calculate the energy and momentum of both the scattered incident particle and the ejected electron, and to transform them to the global coordinate system.

Bibliography

eal84

M.J. Berger et al. Icru report 37. Journal of the International Commission on Radiation Units and Measurements, os19(2):, dec 1984. URL: https://doi.org/10.1093/jicru/os19.2.Report37, doi:10.1093/jicru/os19.2.report37.

MC70(1,2,3)

H. Messel and D. Crawford. Electron-Photon shower distribution. Pergamon Press, 1970.

SP71

R.M. Sternheimer and R.F. Peierls. General expression for the density effect for the ionization loss of charged particles. Physical Review B, 3(11):3681–3692, jun 1971. URL: https://doi.org/10.1103/PhysRevB.3.3681, doi:10.1103/physrevb.3.3681.