source: trunk/documents/UserDoc/DocBookUsersGuides/PhysicsReferenceManual/latex/electromagnetic/standard/hion.tex @ 1222

\section[Ionization]{Hadron and Ion Ionization} \label{hion}

\subsection{Method}
 
The class {\it G4hIonisation} provides the continuous energy loss due to 
ionization and simulates the 'discrete' part of the ionization, that is, 
delta rays produced by charged hadrons.  The class {\it G4ionIonisation} is
intended for the simulation of energy loss by positive ions with change greater than unit.
Inside these classes the following models are used:
\begin{itemize}
\item
{\it G4BetherBlochModel} (valid for protons with $T > 2\; MeV$)
\item 
{\it G4BraggModel} (valid for protons with $T < 2\; MeV$)
\item
{\it G4BraggIonModel} (valid for protons with $T < 2\; MeV$)
\end{itemize}
The scaling relation (\ref{enloss.sc}) is a basic conception for the description 
of ionization of heavy charged particles. It is used both in energy loss
calculation and in determination of the validity range of models. Namely the
$T_p = 2 MeV$ limit for protons is scaled for a particle with mass $M_i$ 
by the ratio of the particle mass to the proton mass $T_i = T_p M_p/M_i$.

For all ionization models the value of the maximum energy 
transferable to a free electron $T_{max}$ is given by the following relation \cite{hion.pdg}:
\begin{equation}
\label{hion.c}
T_{max} =\frac{2m_ec^2(\gamma^2 -1)}{1+2\gamma (m_e/M)+(m_e/M)^2 },
\end{equation}
where $m_e$ is the electron mass and $M$ is the mass of the incident particle.
The method of calculation of the continuous energy loss and the total 
cross-section are explained below. 

\subsection{Continuous Energy Loss}

The integration of \ref{comion.a} leads to the Bethe-Bloch restricted energy 
loss ($T < T_{cut}$ formula \cite{hion.pdg}, which is modified taken into 
account various corrections \cite{hion.ahlen}:
\begin{equation}
\label{hion.d}
\frac{dE}{dx} =
       2 \pi r_e^2 mc^2 n_{el} \frac{z^2}{\beta^2}
       \left [\ln \left(\frac{2mc^2 \beta^2 \gamma^2 T_{up}} {I^2} \right)
       - \beta^2 \left( 1 + \frac{T_{up}}{T_{max}} \right)
       - \delta - \frac{2C_e}{Z} + F\right ]
\end{equation}
 where
\[
\begin{array}{ll}
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{electrons density in the material}     \\
I            & \mbox{mean excitation energy in the material}\\
Z            & \mbox{atomic number of the material}         \\
$z$          & \mbox{charge of the hadron in units of the electron change} \\ 
\gamma       & \mbox{$E/mc^2$}                              \\
\beta^2      & 1-(1/\gamma^2)                               \\
T_{up}       & \min(T_{cut},T_{max})                        \\
\delta       & \mbox{density effect function}               \\
C_e          & \mbox{shell correction function}             \\
F            & \mbox{high order corrections}
\end{array}
\]
In a single element the electron density is
$$ n_{el} = Z \: n_{at} = Z \: \frac{\mathcal{N}_{av} \rho}{A} $$
($\mathcal{N}_{av}$: Avogadro number, $\rho$: density of the material, 
 $A$: 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} . 
$$
$w_i$ is the proportion by mass of the $i^{th}$ element, with molar mass $A_i$.


The mean excitation energy $I$ for all elements is tabulated according to 
the ICRU recommended values \cite{hion.ICRU37}.

\subsubsection{Shell Correction} 

$2C_e/Z$ is the so-called {\it shell correction term} which accounts for the 
fact of interaction of atomic electrons with atomic nucleus. This term more visible
at low energies and for heavy atoms. The classical 
 expression for the term \cite{hion.ICRU49} is used 
\begin{equation}
\label{hion.dh}
C = \sum{C_{\nu}(\theta_{\nu},\eta_{\nu})}, \;\; \nu=K,L,M,..., 
\; \theta=\frac{J_{\nu}}{\epsilon_{\nu}}, \;\; \eta_{\nu}=\frac{\beta^2}{\alpha^2 Z^2_{\nu}},
\end{equation}
where $\alpha$ is the fine structure constant, $\beta$ is the hadron velocity,
$J_{\nu}$ is the ionisation energy of the shell $\nu$, $\epsilon_{\nu}$ is 
Bohr ionisation energy of the shell $\nu$, $Z_{\nu}$ is the effective charge of the
shell $\nu$. 
First terms $C_K$ and $C_L$ can be analytically computed in using an assumption
non-relativistic hydrogenic wave functions \cite{hion.37,hion.38}.
The results \cite{hion.39} of tabulation of these computations in the interval of parameters
$\eta_{\nu} = 0.005\div 10$ and $\theta_{\nu}=0.25 \div 0.95$ are used directly.
For higher values of $\eta_{\nu}$ the parameterization \cite{hion.39}
is applied:
\begin{equation}
C_{\nu} = \frac{K_1}{\eta} + \frac{K_2}{\eta^2} + \frac{K_3}{\eta^3},
\end{equation}
where coefficients $K_i$ provide smooth shape of the function. 
The effective nuclear charge for the $L$-shell can be reproduced as $Z_L = Z - d$, $d$ 
is a parameter shown in Table \ref{hion.t}.
\begin{table}[hbt]
\begin{centering}
\begin{tabular}{|c|c|c|c|c|c|c|c|c|} 
\hline
$Z$ & 3    & 4    & 5    & 6    & 7    & 8    & 9    & $>$9\\ \hline
$d$ & 1.72 & 2.09 & 2.48 & 2.82 & 3.16 & 3.53 & 3.84 & 4.15\\ \hline
\end{tabular}
\caption{Effective nuclear charge for the $L$-shell \cite{hion.ICRU49}.}
\label{hion.th}
\end{centering}
\end{table}
For outer
shells the calculations are not available, so $L$-shell parameterization is used and the
following scaling relation 
\cite{hion.ICRU49,hion.40} is applied:
\begin{equation}
\label{hion.dd}
C_{\nu} = V_{\nu}C_L(\theta_L,H_{\nu}\eta_L), \;\; V_{\nu}=\frac{n_{\nu}}{n_L}, \;\; H_{\nu}=\frac{J_{\nu}}{J_L},
\end{equation}
where $V_{\nu}$ is a vertical scaling factor proportional to number of
electrons at the shell $n_{\nu}$.
The contribution of the shell correction term is about 10\%  for protons
at $T = 2 MeV$.   

\subsubsection{Density Correction} 

$\delta$ is a correction term which takes into account the reduction in energy 
loss due to the so-called {\it density effect}.  This becomes important at 
high energies because media have a tendency to become polarized as the 
incident particle velocity increases.  As a consequence, the atoms in a 
medium can no longer be considered as isolated.  To correct for this effect 
the formulation of Sternheimer~\cite{hion.sternheimer} is used:
\input{electromagnetic/utils/densityeffect}

\subsubsection{High Order Corrections} 

High order corrections term to Bethe-Bloch formula (\ref{hion.d}) can 
be expressed as
\begin{equation}
\label{hion.cor}
F = G - S + 2(z L_1 + z^2 L_2),
\end{equation}
where G is the Mott correction term, S is the finite size correction term, 
$L_1$ is the Barkas correction, $L_2$ is the Bloch correction. The Mott term 
\cite{hion.ahlen} describes the close-collision corrections tend to become
more important at large velocities and higher charge of projectile.
The Fermi result is used:
\begin{equation}
G = \pi\alpha z\beta.
\end{equation}
The Barkas correction term describes distant collisions. The parameterization
of Ref. is expressed in the form:
\begin{equation}
L_1 = \frac{1.29 F_A(b/x^{1/2})}{Z^{1/2}x^{3/2}}, \;\; x = \frac{\beta^2}{Z\alpha^2},
\end{equation}
where $F_A$ is tabulated function \cite{hion.Ashley}, b is scaled minimum impact parameter 
shown in Table \ref{hion.t1}. This and other corrections depending on atomic 
properties are assumed to be additive for mixtures and compounds.
\begin{table}[hbt]
\begin{centering}
\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|} 
\hline
$Z$ & 1 ($H_2$ gas) & 1   & 2   & 3 - 10 & 11 - 17 & 18  & 19 - 25 & 26 - 50 & $>$ 50\\ \hline
$d$ & 0.6         & 1.8 & 0.6 & 1.8    & 1.4     & 1.8 & 1.4     & 1.35    & 1.3\\ \hline
\end{tabular}
\caption{Scaled minimum impact parameter b \cite{hion.ICRU49}.}
\label{hion.t1}
\end{centering}
\end{table}
For the Bloch correction term  the classical expression \cite{hion.ICRU49} is following:
\begin{equation}
z^2L_2 = -y^2 \sum^{\infty}_{n=1} \frac{1}{n(n^2 + y^2)}, \;\; y = \frac{z\alpha}{\beta}.
\end{equation}
The finite size correction term takes into account the space 
distribution of charge of 
the projectile particle. For muon it is zero, for hadrons this term become 
visible at energies above few hundred GeV and the following parameterization 
\cite{hion.ahlen} is used:
\begin{equation}
S = ln(1 + q), \;\; q = \frac{2 m_e T_{max}}{ \varepsilon^2},
\end{equation}
where $T_{max}$ is given in relation (\ref{hion.c}), $\varepsilon$ is proportional to
the inverse effective radius of the projectile (Table \ref{hion.t2}).
\begin{table}[hbt]
\begin{centering}
\begin{tabular}{|c|c|} 
\hline
mesons, spin = 0 ($\pi^{\pm}$, $K^{\pm}$) & $0.736\;GeV$\\ \hline
baryons, spin = 1/2 & $0.843\;GeV$\\ \hline
ions & $0.843\; A^{1/3}\;GeV$\\ \hline
\end{tabular}
\caption{The values of the $\varepsilon$ parameter for different particle types.}
\label{hion.t2}
\end{centering}
\end{table}
All these terms 
break scaling relation (\ref{enloss.sc}) if the projectile particle charge differs from
$\pm$1. To take this circumstance into account  in {\it G4ionIonisation} process 
at initialisation time 
the term $F$ is ignored for the computation of the $dE/dx$ table. 
At run time this term is taken into account by adding to the mean energy loss a value
\begin{equation}
\Delta T' = 2 \pi r_e^2 mc^2 n_{el} \frac{z^2}{\beta^2} F\Delta s,
\end{equation}
where $\Delta s$ is the {\it true step length} and $F$ is the high order correction
term (\ref{hion.cor}).

\subsubsection{Parameterizations at Low Energies} 

For scaled  energies below $T_{lim} = 2\;MeV$ shell correction becomes very 
large and precision of the 
Bethe-Bloch formula degrades, so parameterisation of evaluated data
for stopping powers at low energies is required. 
These parameterisations for  all atoms
is available from ICRU'49 report \cite{hion.ICRU49}. The proton parametrisation is used 
in {\it G4BraggModel}, which is included by default in the process {\it G4hIonisation}.
The alpha particle parameterisation is used in the
{\it G4BraggIonModel}, which is included by default in the process {\it G4ionIonisation}. 
To provide a smooth transition between low-energy and high-energy models the modified
energy loss expression is used for high energy
\begin{equation}
S(T) = S_H (T) + (S_L(T_{lim}) - S_H(T_{lim}))\frac{T_{lim}}{T}, \;\; T > T_{lim},
\end{equation}
where $S$ is smoothed stopping power, $S_H$ is stopping power from formula (\ref{hion.d})
 and $S_L$ is the low-energy parameterisation. 

The precision of Bethe-Bloch formula for $T>10 MeV$ is within 2\%, below the precision 
degrades and at $1 keV$ only 20\% may be garanteed. In the energy interval $1 - 10 MeV$
the quality of description of the stopping power varied from atom to atom. To 
provide more stable and precise parameterisation the data from 
the NIST databases are included inside the standard package. These data are provided for
74 materials of the NIST material database \cite{hion.nist}. The data from the PSTAR database 
are included into {\it G4BraggModel}. The data from the ASTAR database   
are included into {\it G4BraggIonModel}. So, if Geant4 material is defined as a NIST 
material, than NIST data are used for low-energy parameterisation of stopping power.
If material is not from the NIST database, then the ICRU'49 parameterisation is used.


\subsubsection{Nuclear Stopping} 

Nuclear stopping due to elastic ion-ion scattering since Geant4 v9.3
can be simulated with the continuous process
{\it G4NuclearStopping}. By default 
this correction is active and the ICRU'49 parameterisation \cite{hion.ICRU49} is used, 
which is implemented in the model class {\it G4ICRU49NuclearStoppingModel}.  


\subsection{Total Cross Section per Atom and Mean Free Path}

For $T \gg I $ the differential cross section can be written as
\begin{equation}
\label{hion.i}
\frac{d\sigma }{dT} = 2\pi r_e^2 mc^2 Z \frac{z_p^2}{\beta^2} \frac{1}{T^2}
     \left[ 1 - \beta^2 \frac{T}{T_{max}} + \frac{T^2}{2E^2} \right] 
\end{equation}
\cite{hion.pdg}.  In {\sc Geant4} $T_{cut} \geq 1$ keV.  Integrating from
$T_{cut}$ to $T_{max}$ gives the total cross section per atom :
\begin{eqnarray}
\label{hion.j}
\sigma (Z,E,T_{cut}) & = & \frac {2\pi r_e^2 Z z_p^2}{\beta^2} mc^2 \times 
  \\ & &     \left[ \left( \frac{1}{T_{cut}} - \frac{1}{T_{max}} \right)
                   - \frac{\beta^2}{T_{max}} \ln \frac{T_{max}}{T_{cut}}
                   + \frac{T_{max} - T_{cut}}{2E^2} 
             \right]  \nonumber   
\end{eqnarray}
The last term is for spin $1/2$ only.  In a given material the mean free path 
is:
\begin{equation}
\begin{array}{lll} 
\lambda = (n_{at} \cdot \sigma)^{-1} & or &
\lambda = \left( \sum_i n_{ati} \cdot \sigma_i \right)^{-1}
\end{array}
\end{equation}
The mean free path is tabulated during initialization as a function of the
material and of the energy for all kinds of charged particles.

\subsection{Simulating Delta-ray Production}

A short overview of the sampling method is given in Chapter \ref{secmessel}.
Apart from the normalization, the cross section \ref{hion.i} can be 
factorized :
\begin{eqnarray}
\frac{d\sigma}{dT}=f(T) g(T) &with& T \in [T_{cut}, \ T_{max}]
\end{eqnarray}
where
\begin{eqnarray}
f(T) &=& \left(\frac{1}{T_{cut}} - \frac{1}{T_{max}} \right) \frac{1}{T^2} \\
g(T) &=& 1 - \beta^2 \frac{T}{T_{max}} + \frac{T^2}{2E^2} .
\end{eqnarray}
The last term in $g(T)$ is for spin $1/2$ only.  The energy $T$ is chosen by
\begin{enumerate}
\item sampling $T$ from $f(T)$
\item calculating the rejection function $g(T)$ and accepting the sampled 
$T$ with a probability of $g(T)$.
\end{enumerate}
After the successful sampling of the energy, the direction
of the scattered electron is generated with respect to the direction of the
incident particle. The azimuthal angle $\phi$ is generated isotropically. 
The polar angle $\theta$ is calculated from energy-momentum conservation.
This information is used to calculate the energy and momentum of both 
scattered particles and to transform them into the {\em global} coordinate 
system.

\subsubsection{Ion Effective Charge}

As ions penetrate matter they exchange electrons with the medium. In the 
implementation of {\it G4ionIonisation} the effective charge approach is 
used \cite{hion.Ziegler85}. 
A state of equilibrium between the ion and the medium is assumed, so that 
the ion's effective charge can be calculated as a function of its kinetic 
energy in a given material.  Before and after each step the dynamic 
charge of the ion is recalculated and saved in $G4DynamicParticle$, where 
it can be used not only for energy loss calculations but also for the
sampling of transportation in an electromagnetic field.

The ion effective charge is expressed via 
the ion charge $z_i$ and the
fractional effective charge of ion $\gamma_i$:
\begin{equation}
z_{eff} = \gamma_i z_i.
\label{hlei.p}
\end{equation}  
For helium ions 
fractional effective charge
is parameterized for all elements
\begin{eqnarray}
(\gamma_{He})^2 & = &\left (1-\exp\left [-\sum_{j=0}^5{C_jQ^j}\right ]\right)
\left ( 1  + \frac{ 7 + 0.05  Z }{1000} \exp( -(7.6-Q)^2 ) \right )^2, 
\nonumber \\
 Q & = & \max ( 0, \ln T),  
\label{hion.q} 
\end{eqnarray}  
where the coefficients $C_j$ are the same for all elements, and the
helium ion kinetic energy $T$ is in $keV/amu$.


The following expression is used for heavy ions \cite{hion.BK}:  
\begin{equation}
\gamma_i = \left ( q + \frac{1-q}{2} \left (\frac{v_0}{v_F} \right )^2 
\ln {\left ( 1 + \Lambda^2 \right )} \right )
\left ( 1 + \frac{(0.18+0.0015Z)\exp(-(7.6-Q)^2)}{Z_i^2} \right ),
\label{hion.s}
\end{equation}  
where $q$ is 
the fractional average charge of the ion,
$v_0$ is the Bohr velocity,
$v_F$ is the Fermi velocity of  
the electrons in the target medium, and $\Lambda$ is 
the term taking into account the screening effect:
\begin{equation}
\Lambda = 10 \frac{v_F}{v_0} \frac{(1-q)^{2/3}}{Z_i^{1/3}(6+q)}.
\label{hion.t}
\end{equation}  
The Fermi velocity of the medium is of the same order as the Bohr velocity, and
its exact value depends on the detailed electronic structure of the medium.
The expression for the fractional average charge of the ion is the following:
\begin{equation}
q = [1 -\exp(0.803y^{0.3}-1.3167y^{0.6}-0.38157y-0.008983y^2)],
\label{hion.u}
\end{equation}  
where $y$ is a parameter that depends on the ion velocity $v_i$
\begin{equation}
y = \frac{v_i}{v_0Z^{2/3}} \left ( 1 +\frac {v_F^2}{5v_i^2} \right ).
\label{hion.v}
\end{equation}  
The parametrisation of the effective charge of the ion applied 
if the kinetic energy is below limit value
\begin{equation}
T < 10 z_i \frac{M_i}{M_p}\;MeV,
\label{hion.x}
\end{equation}  
where $M_i$ is the ion mass and $M_p$ is the proton mass.


\subsection{Status of this document}
  09.10.98 created by L. Urb\'an. \\
  14.12.01 revised by M.Maire \\
  29.11.02 re-worded by D.H. Wright \\
  01.12.03 revised by V. Ivanchenko     \\
  21.06.07 revised by V. Ivanchenko     \\

\begin{latexonly}

\begin{thebibliography}{99}

\bibitem{hion.pdg}
W.-M.~Yao et al.,  Jour. of Phys. G33 (2006) 1.
\bibitem{hion.ahlen} 
S.P. Ahlen, Rev. Mod. Phys. 52 (1980) 121.
\bibitem{hion.ICRU37} 
ICRU (A.~Allisy et al),
Stopping Powers for Electrons and Positrons,
ICRU Report 37, 1984.
\bibitem{hion.ICRU49}ICRU (A.~Allisy et al),
Stopping Powers and Ranges for Protons and Alpha
Particles,
ICRU Report 49, 1993.
\bibitem{hion.37} 
M.C.~Walske, Phys. Rev. 88 (1952) 1283.
\bibitem{hion.38} 
M.C.~Walske, Phys. Rev. 181 (1956) 940.
\bibitem{hion.39} 
G.S.~Khandelwal, Nucl. Phys. A116 (1968) 97.
\bibitem{hion.40} 
H.~Bichsel, Phys. Rev. A46 (1992) 5761.
\bibitem{hion.sternheimer}
  R.M.~Sternheimer. Phys.Rev. B3 (1971) 3681.
\bibitem{hion.Ashley}
J.C.~Ashley, R.H.~Ritchie and W.~Brandt, Phys. Rev. A8 (1973) 2402.
\bibitem{hion.nist}
http://physics.nist.gov/PhysRevData/contents-radi.html
\bibitem{hion.Ziegler85}
J.F.~Ziegler, J.P.~Biersack, U.~Littmark, The Stopping
and Ranges of Ions in Solids. Vol.1, Pergamon Press, 1985.
\bibitem{hion.BK}
W.~Brandt and M.~Kitagawa, Phys. Rev. B25 (1982) 5631.
\end{thebibliography}

\end{latexonly}

\begin{htmlonly}

\subsection{Bibliography}

\begin{enumerate}
\item{hion.pdg}
W.-M.~Yao et al.,  Jour. of Phys. G33 (2006) 1.
\item{hion.ICRU37} 
S.P. Ahlen, Rev. Mod. Phys. 52 (1980) 121.
\item{hion.ICRU37} 
ICRU (A.~Allisy et al),
Stopping Powers for Electrons and Positrons,
ICRU Report 37, 1984.
\item{hion.ICRU49}ICRU (A.~Allisy et al),
Stopping Powers and Ranges for Protons and Alpha
Particles,
ICRU Report 49, 1993.
\item{hion.37} 
M.C.~Walske, Phys. Rev. 88 (1952) 1283.
\item{hion.38} 
M.C.~Walske, Phys. Rev. 181 (1956) 940.
\item{hion.39} 
G.S.~Khandelwal, Nucl. Phys. A116 (1968) 97.
\item{hion.40} 
H.~Bichsel, Phys. Rev. A46 (1992) 5761.
\item{hion.sternheimer}
  R.M.~Sternheimer. Phys.Rev. B3 (1971) 3681.
\item{hion.Ashley}
J.C.~Ashley, R.H.~Ritchie and W.~Brandt, Phys. Rev. A8 (1973) 2402.
\item{hion.nist}
http://physics.nist.gov/PhysRevData/contents-radi.html
\item{hion.Ziegler85}
J.F.~Ziegler, J.P.~Biersack, U.~Littmark, The Stopping
and Ranges of Ions in Solids. Vol.1, Pergamon Press, 1985.
\item{hion.BK}
W.~Brandt and M.~Kitagawa, Phys. Rev. B25 (1982) 5631.
\end{enumerate}

\end{htmlonly}