source: trunk/documents/UserDoc/DocBookUsersGuides/PhysicsReferenceManual/latex/electromagnetic/standard/msc.tex @ 1211

\section[Multiple Scattering]{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.
Geant4 uses a universal interface {\it G4VMultipleScattering} 
is used by all MSC processes:
\begin{itemize}
\item
G4eMultipleScattering;
\item
G4hMultipleScattering;
\item
G4MuMultipleScattering.
\item
G4MultipleScattering;
\end{itemize}
For concrete simulation the {\it G4VMscModel} interface is used, which is an extension
of the base {\it G4VEmModel} interface. The following models are available:
\begin{itemize}
\item
G4UrbanMscModel;
\item
G4UrbanMscModel2;
\item
G4UrbanMscModel90;
\item
G4WentzelVIModel.
\end{itemize}
The last model is not yet in the production mode, so below we will concentrate 
on models developed by L.~Urban \cite{msc.urban}. Today default model for 
electrons and positrons is {\it G4UrbanMscModel2}, for hadrons - {\it G4UrbanMscModel90}.

\subsection{Introduction}

MSC simulation algorithms can be classified as either {\em detailed} or 
{\em 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
enough kinetic energies.  For larger kinetic energies 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
displacement, energy loss, and change of direction of the charged particle. 
These quantities are computed from the multiple scattering theories used in 
the codes.  The accuracy of the condensed simulations is limited by the 
approximations of the multiple scattering theories. \\

\noindent
Most particle physics simulation codes use the multiple scattering theories
of Moli\`ere \cite{msc.moliere}, Goudsmit and Saunderson \cite{msc.goudsmit}
and Lewis \cite{msc.lewis}.  The theories of Moli\`ere 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.
Therefore each of the MSC simulation codes incorporates its own algorithm to 
determine the 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 in the MSC codes.  Therefore the simulation results can depend 
on the value of the step length and generally one has to select the value of 
the step length carefully. \\

\noindent
A new class of MSC simulation, the {\em mixed} simulation algorithms (see 
e.g.\cite{msc.fernandez}), appeared in the literature recently.  The mixed 
algorithm simulates the {\em hard} collisions one by one and uses a MSC theory to 
treat the effects of the {\em 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. \\ 

\noindent
The MSC model used in Geant4 belongs to the class of condensed 
simulations.  It 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 the Lewis theory \cite{msc.lewis}.

\subsection{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 {\em geometrical path length}, $z$. 
In the absence of a magnetic field, this is a straight line.
For non-zero fields, $z$ is the shortest distance 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 physical interactions
like multiple scattering.  
This distance is called the {\em true path length}, $t$.
Constraints on $t$ are imposed by the physical processes acting on the 
particle. \\

\noindent
The properties of the multiple scattering process are completely determined
by the {\em 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

\begin{equation}
 \frac {1}{\lambda_k} = 2 \pi n_a \int_{-1}^{1} \left[1 - P_k(cos\chi) \right]
                   \frac{d\sigma(\chi)}{d\Omega} d(cos\chi)
\label{msc.a1}
\end{equation}

\noindent
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. \\

\noindent 
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
\begin{equation}
 \langle z \rangle =
            \lambda_1 \left[ 1-\exp \left(-\frac{t}{\lambda_1} \right)\right]
\label{msc.a}
\end{equation}

\noindent
Eq.~\ref{msc.a} 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 {\em path length correction}. This formula and other expressions for 
the first moments of the spatial distribution were taken from either
 \cite{msc.fernandez} or \cite{msc.kawrakow}, but were originally
calculated by Goudsmit and Saunderson \cite{msc.goudsmit} and Lewis
\cite{msc.lewis}. \\

\noindent 
At the end of the true step length, $t$, the scattering angle is $\theta$.
The mean value of $cos\theta$ is
\begin{equation}
  \langle cos\theta \rangle = \exp \left[-\frac{t}{\lambda_1} \right]
\label{msc.c}
\end{equation}
The variance of $cos\theta$ can be written as
\begin{equation}
  \sigma^2 = \langle cos^2\theta \rangle - \langle cos\theta \rangle ^2 = 
             \frac{1 + 2 e^{- 2 \kappa \tau}} {3} - e^{-2 \tau}
\label{msc.c1}
\end{equation}

\noindent
where $\tau = t/\lambda_1$ and $\kappa = \lambda_1/\lambda_2$. \\

\noindent
 The mean lateral displacement is given by a more complicated formula
\cite{msc.fernandez}, but this quantity can also be calculated relatively 
easily and accurately.  The square of the {\em mean lateral displacement} is
\begin{equation}
 \langle x^2 + y^2 \rangle = \frac{4 \lambda_1^2}{3} \ \left[\tau 
                  - \frac{\kappa+1}{\kappa}
                  + \frac{\kappa}{\kappa-1} e^{-\tau} -
                  \frac{1}{\kappa (\kappa -1)} e^{-\kappa \tau} \right]
\label{msc.e1}
\end{equation}

\noindent
Here it is assumed that the initial particle direction is parallel to the 
the $z$ axis.\\

\noindent 
 The lateral correlation is determined by the equation
\begin{equation}
  \langle x v_x+y v_y \rangle = \frac{2 \lambda_1}{3} \ \left[1
                      - \frac{\kappa}{\kappa-1} e^{-\tau} +
                        \frac{1}{\kappa-1} e^{-\kappa \tau} \right]
\label{msc.e2}
\end{equation}

\noindent
 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.\\

\noindent
The transport mean free path values have been calculated by Liljequist et al.
\cite{msc.liljequist1}, \cite{msc.liljequist2} for electrons and positrons in
the kinetic energy range \mbox{100 eV - 20 MeV} in 15 materials.  The 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 (\cite{msc.mayol}).
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 \cite{msc.kawrakow}).  Hence, a constant value of 2.5 is used
in the model.\\

\noindent
  Nuclear size effects are negligible for low energy particles and
they are accounted for in the Born approximation in \cite{msc.mayol},
so there is no need for extra corrections of this kind in the model.

\subsection{Path Length Correction}

As mentioned above, the path length correction refers to the transformation
true path length $\longrightarrow$ geometrical path length and its inverse.
The true path length $\longrightarrow$ geometrical path length transformation
is given by Eq.~\ref{msc.a} 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. \ref{msc.c},\ref{msc.a}
should be modified as
\begin{equation}
 \langle cos\theta \rangle = \exp \left[-\int_0^t \frac{du}{\lambda_1 (u)}
                             \right]
   \label{msc.ax}
\end{equation}

\begin{equation}
 \langle z \rangle =  \int_0^t \langle cos\theta \rangle_u \ du  
       \label{msc.bx}
\end{equation}

\noindent
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. \\

\noindent
In order to compute Eqs. \ref{msc.ax},\ref{msc.bx} 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:
   \begin{equation}
          \lambda_1(t) =  \lambda_{10} ( 1- \alpha t)       \label{msc.cx}
   \end{equation}

\noindent
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.~\ref{msc.cx} is 
\emph{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.~\ref{msc.ax} -  \ref{msc.cx} the explicit formula for 
$\langle cos\theta \rangle$ and $\langle z \rangle$ are :
\begin{equation}
   \langle cos\theta \rangle = 
         \left(1 - \alpha t \right)^{\frac{1}{ \alpha \lambda_{10}}}
  \label{msc.ff}
\end{equation}

\begin{equation}
   \langle z \rangle = \frac{1} { \alpha (1 + \frac{1}{\alpha \lambda_{10}})}
        \left[ 1 - (1 - \alpha t)^{1+ \frac{1}{ \alpha \lambda_{10}}} \right] 
  \label{msc.dx}
\end{equation}

\noindent
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
 \begin{equation}
    \alpha = \frac{\lambda_{10} - \lambda_{11}} {t \lambda_{10}} 
   \label{msc.ex}
 \end{equation}

\noindent
 At low energies ( $T_{kin} < M$ ,  M - particle mass) $\alpha $ has a simpler 
 form:
 \begin{equation}
      \alpha = \frac{1} { r_0}   \label{msc.fx}
 \end{equation}
where $r_0$ denotes the range of the particle at the start of the step.

\noindent
 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
 \begin{equation}
   \langle z \rangle = 
     \lambda_{10} \left[ 1-\exp{\left( -\frac{t}{\lambda_{10}}\right)}\right] 
    \label{msc.gx}
 \end{equation}

\noindent
Eq. \ref{msc.dx} or \ref{msc.gx} gives the mean value of the geometrical step
length for a given true step length. \\

\noindent
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]$ :
 \begin{equation}
   f(v) = (k+1)(k+2) v^k (1-v)  
 \label{msc.d2}
 \end{equation}

\noindent
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.~\ref{msc.dx} or 
\ref{msc.gx}. Hence
\begin{equation}
    k = \frac{3 \langle z \rangle - t}{t - \langle z \rangle}
\label{msc.d3}
\end{equation}

\noindent
The value of $z = v t$ is sampled using  $f(v)$ if $k > 0$, otherwise
$z = \langle z \rangle$ is used. \\

\noindent
The geometrical path length $\longrightarrow$ true path length transformation
is performed using the mean values. The transformation can be written as
\begin{equation}
   t(z) = \langle t \rangle = -\lambda_1 \log\left(1-\frac{z}{\lambda_1}\right)
\label{msc.d4}
\end{equation}
 if the geometrical step is small and
 \begin{equation}
   t(z) = \frac{1}{\alpha} \left[ 1 - (1 - \alpha w z)^{\frac{1}{w}} \right]
 \label{msc.hx}
 \end{equation}
 where  $$w = 1 + \frac{1}{\alpha \lambda_{10}}$$ if the step is not small, i.e.
 the energy loss should be taken into account. \\

\noindent
This transformation is needed when the particle arrives at a volume boundary,
causing the step to be geometry-limited.  In this case the true path length
should be computed in order to have the correct energy loss of the particle
after the step.

\subsection{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. \ref{msc.c} and \ref{msc.c1} are satisfied.  The functional form of $g$ is
\begin{equation}
  g(u) = q [p g_1(u) + (1-p) g_2(u)] + (1-q) g_3(u)
\label{msc.d}
\end{equation}

\noindent
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
\begin{equation}
 g_1(u) = C_{1}\hspace{3mm} e^{-a (1-u)} \hspace{2cm} -1 \leq u_0 \leq u \leq 1
\label{msc.d5}
\end{equation}
\begin{equation}
 g_2(u) =  C_{2}\hspace{3mm} \frac{1} { (b-u)^d}
      \hspace{2cm}  -1 \leq u \leq u_0 \leq 1
\label{msc.d6}
\end{equation}
\begin{equation}
 g_3(u) =  C_{3} \hspace{4.8cm}  -1 \leq u \leq 1
\label{msc.d7}
\end{equation}
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 .


\subsection{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:
\begin{equation}
  p\hspace{1mm} g_1(u_0) = (1-p)\hspace{1mm} g_2(u_0)
\label{msc.p1}
\end{equation}
\begin{equation}
  p\hspace{1mm} a\hspace{1mm} g_1(u_0) = (1-p)\hspace{1mm}
  \frac{d}{b-u_0}\hspace{1mm} g_2(u_0)
\label{msc.p2}
  \end{equation}

\noindent
A third constraint comes from Eq. \ref{msc.ax} : $g(u)$ must give the same
mean value for $u$ as the theory. \\

\noindent
 It follows from Eqs. \ref{msc.ff} and \ref{msc.d} that
\begin{equation}
 q \{ p \langle u \rangle_1 + (1-p) \langle u \rangle_2 \} =
 [ 1 - \alpha \ t ]\ ^{\frac{1}{\alpha \lambda_{10}}}
\label{msc.par1}
\end{equation}
where $\langle u \rangle_i$ denotes the mean value of $u$ computed from the
distribution $g_i(u)$. \\



\noindent
The parameter $a$ was chosen according to a modified Highland-Lynch-Dahl formula
for the width of the angular distribution \cite{msc.highland},
\cite{msc.lynch}.
\begin{equation}
   a = \frac {0.5} {1 -cos(\theta_0)}
\end{equation}
where $\theta_0$ is
\begin{equation}
  \theta_0 = \frac {13.6 MeV}{ \beta c p} z_{ch} \sqrt{\frac{t}{X_0}}
             \ \left[ 1 + h_c  \ln \left(\frac{t}{X_0} \right)\ \right]
\end{equation}

\noindent
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\`ere 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). \\
\noindent
The model uses a slightly modified Highland-Lynch-Dahl formula to
compute $\theta_0$. 

For electrons/positrons the modified $\theta_0$ formula is (see G4UrbanMscModel2)
\begin{equation}
\theta_0 = \frac {13.6 MeV}{\beta c p} z_{ch} \sqrt{y }  c
\end{equation}
    where 
   
\begin{equation}
  y =  \ln \left(\frac{t}{X_0}\right)
\end{equation}

 The correction term $c$ is

\begin{equation}
    c = c_1 + c_2 y  
\end{equation}                 
  for $y > -6.5$  and
\begin{equation}
    c = c_1 + c_2 y - 0.011 (6.5 + y)  
\end{equation}
  for $y <= 6.5$. The coefficients $c_i$ are

\begin{equation}
    c_1 = 0.885 + ln(Z) \left[ 0.104-0.0177 ln(Z) \right]
\end{equation}

\begin{equation}
    c_2 = 0.028 + ln(Z) \left[ 0.012-0.00125 ln(Z) \right]
\end{equation}

In the case of heavy charged particles (muons,hadrons)
the modified formula for $\theta_0$ is (see G4UrbanMscModel90)

\begin{equation}
\theta_0 = \frac {13.6 MeV}{\beta c p} z_{ch} \sqrt{\frac{t}{X_0} }
                \left[ 1 + 0.105 \ln \left(\frac{t}{X_0}\right)
                         + 0.0035 \left(\ln \left(\frac{t}{X_0}\right)\right)^2
               \right] ^{\frac{1}{2}} f \left(Z \right)
\end{equation}

\noindent
 where
\begin{equation}
   f \left(Z \right) = 1 - \frac{0.24}{ Z \left(Z + 1 \right) }
\end{equation}
\noindent
 is an empirical correction factor based on high energy proton scattering
 data \cite{msc.shen}.
This formula gives a much smaller step dependence in the angular
distribution than the Highland form.  \\

\noindent
The value of the parameter $u_0$ has been chosen as
\begin{equation}
  u_0 \hspace{4mm} = \hspace{4mm} 1 - \frac{\xi}{a}
\end{equation}
where $\xi$ is a constant ($\xi = 3$). \\

 For electrons/positrons the parameter $d$ is set to
\begin{equation}
  d \hspace{4mm} = \hspace{4mm} d_1                           
\end{equation}
  for $y >= -13.5$ and
\begin{equation}
  d \hspace{4mm} = \hspace{4mm} d_1 + d_2 (y+13.5)^3                          
\end{equation}
 for $y < 13.5$
 where
\begin{equation}
  d_1 \hspace{4mm} = \hspace{2mm} 2.134-\ln(Z) (0.1045-0.00602 ln(Z))                    
\end{equation}
\begin{equation}
  d_2 \hspace{4mm} = \hspace{2mm} 0.001126-\ln(Z) (0.0001089+0.0000247 \ln(Z)) 
\end{equation}

 For heavy particle the  parameter $d$ is set to
\begin{equation}
  d \hspace{4mm} = \hspace{4mm} 2.40 - 0.027 \ Z^ \frac{2}{3}
\end{equation}
 This (empirical) expression is obtained comparing the simulation
 results to the data of the MuScat experiment \cite{msc.attwood}. \\

\noindent
The remaining three parameters can be computed from
Eqs. \ref{msc.p1} - \ref{msc.par1}.
The numerical value of the parameters can be found in the
code. \\


\noindent
It should be noted that in this model there is no step limitation
originating from the multiple scattering process.  Another important feature
of this model is that the sum of the 'true' step lengths of the particle, that
is, the total true path length, does not depend on the length of the steps.
Most algorithms used in simulations do not have these properties. \\

\noindent
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. \\

\noindent
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.
  
\subsection{The MSC Process in Geant4}

 The step length of the 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. \\
 
 \noindent
 Then, this minimum step length
 must be compared with the length determined by the geometry of the detectors
 and one has to select the minimum of the 'physics step length' and the
 'geometrical step length' as the actual step length. \\
 
 \noindent
 This is the point where the MSC model comes into the game. All the
 physics processes use the true path length $t$ to sample the interaction point,
 while the step limitation originated from the geometry is a 
 geometrical path length $z$. The MSC algorithm transforms the 'physics step 
 length' into a 'geometrical step length' before the comparison of the two 
 lengths. This 't'\(\rightarrow\)'z' transformation can be called the 
 inverse of the path length correction. \\
 
 \noindent
 After the actual step length has been determined and the particle relocation
 has been performed  the MSC performs the transformation 'z'\(\rightarrow\)'t',
 because the energy loss and scattering computation need the true step length 
 't'. \\
 
 \noindent
 The scattering angle $\theta$ of the particle after the step of length 't' is 
 sampled according to the model function given in Eq.~\ref{msc.d} .
 The azimuthal angle $\phi$ is generated uniformly in the range $[0, 2 \pi]$. \\
 
 \noindent
 After the simulation of the scattering angle, the lateral displacement is
 computed using Eq.~\ref{msc.e1}. Then the correlation given by Eq.~\ref{msc.e2}
 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.

\subsection{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.\\
 
 \noindent
 A very simple step limitation algorithm has been implemented in the
 MSC code to cure this situation. At the start of a track or after
 entering in a new volume, the algorithm restricts
 the step size to a value
 \begin{equation}
  f_r \cdot max\{r,\lambda_1\}
 \end{equation}
 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
 \begin{equation}
   f_{reff} = f_r \cdot \left[  1 - sc + sc * \frac{\lambda_1}{\lambda_{lim}} \right]
 \end{equation}
 ( The numerical values $sc = 0.25$ and $\lambda_{lim} = 1 \hspace{2mm} 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
\begin{equation}
 tlimitmin = max\left[ \frac{\lambda_1}{nstepmax}, \lambda_{elastic} \right]
 \end{equation}
 with $nstepmax = 25$ and $\lambda_{elastic}$ is the elastic mean free path
 of the particle (see later).\\

\noindent
 It can be easily seen that this kind of step limitation poses a real constraint
 only for low energy particles. \\
 
\noindent 
 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
 \begin{equation}
  \frac {d_{geom}}{f_g}
 \end{equation}
 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
 \begin{equation}
  \frac {2 d_{geom}}{f_g}
 \end{equation}
 
\noindent
 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.
 
\subsection{Boundary Crossing Algorithm}

 A special stepping algorithm has been implemented recently (Autumn 2006)
 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). \\

\noindent
  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). \\
 
\noindent 
 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'. \\
 
\noindent  
  The elastic mean free path can be estimated as 
 \begin{equation}
  \lambda_{elastic} = \lambda_1 \cdot rat \left( T_{kin} \right)
 \end{equation}
  where $rat(T_{kin})$ a simple empirical function computed from the elastic and
  first transport cross section values of Mayol and Salvat \cite{msc.mayol}
  
 \begin{equation}
  rat \left( T_{kin} \right) = \frac{0.001 (MeV)^2} 
                              { T_{kin} \left(T_{kin} + 10 MeV \right)}
 \end{equation}
  $T_{kin}$ is the kinetic energy of the particle. \\
  
\noindent  
  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. \\
\noindent  
  The single scattering is determined by the distribution
 \begin{equation}
   g(u) = C \frac{1} {(2 a + 1 -u)^2}                          
 \end{equation}
  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).


\subsection{Implementation Details}

The concrete implementation of described algorithms of simulation is
provided in {\it G4UrbanMscModel2, G4UrbanMscModel90}.
Because multiple scattering is very similar for different particles the base 
class {\it G4VMultipleScattering} was created to collect and provide significant
features of the calculations which are common to different particle types. \\

\noindent
In the {\tt AlongStepGetPhysicalInteractionLength} method the minimum step 
size due to the physics processes is compared with the step size constraints 
imposed by the transportation process and the geometry.  In order to do this, 
the \mbox{'t' step $\rightarrow$ 'z' step} transformation must be performed.  
Therefore, the method should be invoked after the 
{\tt GetPhysicalInteractionLength} methods of other physics processes, 
but before the same method of the transportation process.
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.\\ 

\noindent
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. \\

\noindent
The {\tt PostStepGetPhysicalInteractionLength} method of the multiple
scattering process simply sets the force flag to 'Forced' in order to 
ensure that {\tt PostStepDoIt} is called at every step.  It also returns a 
large value for the interaction length so that there is no step limitation 
at this level. \\

\noindent
The {\tt AlongStepDoIt} function of the process performs the inverse,\linebreak
\mbox{'z' $\rightarrow$ 't'} transformation.  This function should be 
invoked after the \linebreak {\tt AlongStepDoIt} method of the transportation 
process, that is, after the particle relocation is determined by the 
geometrical step length, but before applying any other physics 
{\tt AlongStepDoIt}. \\

\noindent
The {\tt PostStepDoIt} method of the process samples the scattering angle 
and performs the lateral displacement when the particle is not near a boundary.

Default MSC parameter values  optimized per particle type are shown in Table \ref{msc.t}.
Note, that there is three types of step limitation by multiple scattering process:
\begin{itemize}
\item
{\tt Minimal} - only $f_r$ parameter is used, was used for g4 7.1 release;
\item 
{\tt UseSafety} or $skin=0$ - uses particle range and
geometrical safety;
\item 
{\tt UseDistanceToBoundary} - uses particle range, geometrical  safety and linear distance to
geometrical boundary.
\end{itemize}
\begin{table}[hbt]
\begin{centering}
\begin{tabular}{|c|c|c|c|} 
\hline
particle & $e^+$, $e^-$   & {\it muons, hadrons} & {\it ions}\\ \hline
 {\it StepLimitType} & {\it fUseSafety} & {\it fMinimal} & {\it fMinimal}\\ \hline
{\it skin}            &  0    & 0   & 0 \\ \hline
$f_r$                 &  0.04 & 0.2 & 0.2\\ \hline
$f_g$                 &  2.5  & 0.1 & 0.1\\ \hline
{\it LateralDisplacement} &  {\it true}  & {\it true}& {\it false}\\ \hline
\end{tabular}
\caption{The default values of parameters for different particle type.}
\label{msc.t}
\end{centering}
\end{table}
The parameters of the model can be changed via public functions of the base
class {\it G4VMultipleSacttering}. They can be changed for all 
multiple scattering processes simultaneously via {\it G4EmProcessOptions} class or 
via Geant4 UI commands. The following commands are available:\\
\\
{\it /process/msc/StepLimit UseDistanceToBoundary\\
/process/msc/LateralDisplacement false\\
/process/msc/RangeFactor  0.02\\
/process/msc/GeomFactor   2.5\\
/process/msc/Skin   2}
\\
 

\subsection{Status of this document}
 09.10.98  created by L. Urb\'an. \\
 15.11.01  major revision by L. Urb\'an.\\
 18.04.02  updated by L. Urb\'an. \\
 25.04.02  re-worded by D.H. Wright \\
 07.06.02  major revision by L. Urb\'an. \\
 18.11.02  updated by L. Urb\'an, now it describes the new angle distribution. \\
 05.12.02  grammar check and parts re-written by D.H. Wright \\
 13.11.03  revision by L. Urb\'an. \\
 01.12.03  revision by V. Ivanchenko. \\
 17.05.04  revision by L.Urb\'an. \\
 01.12.04  updated by L.Urb\'an. \\
 18.03.05  sampling z + mistyping corrections (mma) \\
 22.06.05  grammar, spelling check by D.H. Wright \\
 12.12.05  revised by L. Urb\'an, according to Geant4 V8.0 \\
 14.12.05  updated implementation Details (mma) \\
 08.06.06  revised by L. Urb\'an, according to Geant4 V8.1 \\
 25.11.06  revised by L. Urb\'an, according to Geant4 V8.2 \\
 29.03.07  revised by L. Urb\'an, for Geant4 V8.3 \\
 13.06.07  modified introduction (mma) \\
 17.06.07  explain effective $F_R$ (L. Urb\'an) \\
 25.06.07  update description of options by V. Ivanchenko \\
 05.12.07  revised by L. Urb\'an, for Geant4 V9.1 \\
 08.12.08  revised by L. Urb\'an, for Geant4 V9.2 \\
 11.12.08  minor revision by V. Ivanchenko \\
 
\begin{latexonly}

\begin{thebibliography}{99}

\bibitem{msc.urban} L.~Urban, A multiple scattering model,
   {\em CERN-OPEN-2006-077, Dec 2006. 18 pp.}
\bibitem{msc.moliere} G.Z.~Moli\`ere
   {\em Z. Naturforsch. 3a (1948) 78. }
\bibitem{msc.lewis} H.W.~Lewis.
   {\em Phys. Rev. 78 (1950) 526. }
\bibitem{msc.goudsmit}S.~Goudsmit and J.L.~Saunderson.
   {\em Phys. Rev. 57 (1940) 24. }
\bibitem{msc.fernandez}J.M.~Fernandez-Varea et al.
   {\em NIM B73 (1993) 447.}
\bibitem{msc.kawrakow} I.~Kawrakow and A.F.~Bielajew
   {\em NIM B 142 (1998) 253. }
\bibitem{msc.liljequist1} D.~Liljequist and M.~Ismail.
   {\em J.Appl.Phys. 62 (1987) 342. }
\bibitem{msc.liljequist2} D.~Liljequist et al.
   {\em J.Appl.Phys. 68 (1990) 3061. }
\bibitem{msc.mayol} R.~Mayol and F.~Salvat
   {\em At.Data and Nucl.Data Tables 65 (1997) 55.}.
\bibitem{msc.highland} V.L.~Highland 
   {\em NIM 129 (1975) 497. }
\bibitem{msc.lynch} G.R.~Lynch and O.I.~Dahl
   {\em NIM B58 (1991) 6. }
\bibitem{msc.shen} G.~Shen et al.
   {\em Phys. Rev. D 20 (1979) 1584.}
\bibitem{msc.attwood} D.~Attwood et al.
   {\em NIM B 251 (2006) 41.}
\end{thebibliography}

\end{latexonly}

\begin{htmlonly}

\subsection{Bibliography}

\begin{enumerate}

\item L.~Urban, A multiple scattering model,
   {\em CERN-OPEN-2006-077, Dec 2006. 18 pp.}
\item  G. Z. Moli\`ere
   {\em Z. Naturforsch. 3a (1948) 78. }
\item H. W. Lewis
   {\em Phys. Rev. 78 (1950) 526. }
\item S. Goudsmit and J. L. Saunderson.
   {\em Phys. Rev. 57 (1940) 24. }
\item J. M. Fernandez-Varea et al.
   {\em NIM B73 (1993) 447.}
\item I. Kawrakow and Alex F. Bielajew
   {\em NIM B 142 (1998) 253. }
\item D. Liljequist and M. Ismail.
   {\em J.Appl.Phys. 62 (1987) 342. }
\item D. Liljequist et al.
   {\em J.Appl.Phys. 68 (1990) 3061. }
\item R. Mayol and F. Salvat
   {\em At.Data and Nucl.Data Tables 65 (1997) 55.}.
\item V.L. Highland
   {\em NIM 129 (1975) 497. }
\item G.R. Lynch and O.I. Dahl
   {\em NIM B58 (1991) 6. }
\item{msc.shen} G. Shen et al.
   {\em Phys. Rev. D 20 (1979) 1584.}
\item{msc.attwood} D. Attwood et al.
   {\em NIM B 251 (2006) 41.}
\end{enumerate}

\end{htmlonly}