Cite this article as: Ulmer W. A new calculation formula of the nuclear cross-section of therapeutic protons. Int J Cancer Ther Oncol 2014; 2(2):020211.
DOI: http://dx.doi.org/10.14319/ijcto.0202.11

Introduction

The knowledge of the total nuclear cross-section Q^tot of protons is an important impact with regard to sophisticated features of therapy planning, since Q^tot provides essential information of the following aspects: Decrease of the fluence of primary protons Φ_pp and release of secondary particles and their transport (secondary protons, neutrons, clusters like H₂¹, H₃¹, He₃², He₄³, heavy recoil nuclei, which usually undergo either a β⁺ or β^- decay with additional emission of a γ-quant). With regard to secondary protons we have to differ between two kinds, namely protons resulting from nuclear reactions with production of heavy recoils and those protons, which are, in reality, primary protons and have undergone elastic and inelastic scatter by strong interactions according to the Breit-Wigner formula. Elastic scatter by nuclear forces is only a deflection of the projectile protons with additional energy-momentum transfer to the whole target nucleus, whereas inelastic scatter is connected to excitations of nuclear vibrations, rotations and transitions to excited states without releasing other nuclear particles, i.e. some quantum number will be changed. The resonances due to Breit-Wigner formula represent the main part of the resonance domain E_res according to Figure 1. Nuclear reaction types cannot be regarded as simple resonances, they mainly occur for proton energies E > E_res.

In previous publications ^{2, 3}, we have developed a calculation method based on a nonlinear and nonlocal Schrödinger equation with a Gaussian kernel and on an interacting many-body system containing strong interactions, spin-orbit coupling and Coulomb interaction with inclusion of various excited configurations. The results of these calculation methods, which provide excited nuclear states, virtual compounds and nuclear reactions via S-matrix, transition probabilities and finally, total nuclear cross-sections Q^tot(E), can be translated to the collective nuclear excitation model. This model only uses Z and A_N as parameters and suitable analytic functions to describe all properties of Q^tot(E). The complete contents of this figure have been previously discussed. For protons, a threshold energy E_Th exists to surmount the potential barrier of the oxygen nucleus. Thus for proton energies less than E_Th nuclear reactions cannot occur. At E = E_res = 20.12 MeV Q^tot(E) shows a maximum value, and thereafter, it decreases exponentially to reach the asymptotic behavior at about E = 100 MeV. According to an integration procedure ^{2, 6} the analytic version of Figure 1 (equation (1)) provides the decrease of primary protons. The following integration procedure has to be carried out:

FIG. 1: Total nuclear cross-section of the proton – oxygen interaction.^{4, 5}

With regard to oxygen we have to put E_Th = 7 MeV (O) and Mc² = 938.27 MeV. Equation (2) can be summarized by Figure 2.

FIG. 2: Decrease of the fluence of primary protons in water.

FIG. 3: Nuclear potential energy of oxygen.

A further result is that in the environment of a nucleus the effective potential can be calculated by a linear combination of two Gaussians. This is not true for longer ranges of Coulomb forces. However, this is not a simple task, since the shielding of the nuclear repulsion by the shell electrons has to be accounted for. The nuclear potential according to Figure 3 and the related parameters can be best calculated by equation (3), which assumes the shape:

The units of this formula are stated in terms of units of Rstrong (range of strong interaction) according to Figure 3.

Methods and Materials

In previous publications ^{2, 3}, we have performed the task to calculate all nuclear properties to determine finally with the help of the S-matrix and transition probability of all possible configuration states the total nuclear cross-section Q^tot(E) in terms of Z and A_N. The results are closely related to the Bethe-Weizsäcker formula for the nuclear binding energy. The old formula for the calculation of E_Th and Q^tot(E) according to a previous paper will be stated in section “Summary of previous investigations”. In order to determine all these properties we have to know the threshold energy E^Th, i.e. that energy necessary to surmount the nuclear potential mount according to Figure 3. There exists also the possibility to circumvent this threshold restriction by the quantum mechanical tunneling process with energy E < E_Th, but this is only a rather small effect, which leads to a little roundness of Q^tot(E) at E < E_Th. This aspect will be considered in section “Calculation procedure of the energy levels of excited states and quantum mechanical tunneling effect for E < E_th”.

Summary of previous investigations

According to the results presented in ^{2, 3} we need for the calculation of Q^tot(E), at first, the threshold energy E_Th as a function of Z and A_N. This function can be obtained via the formulas (4 - 6). The second step provides the determination of the total nuclear cross-section Q^tot(E), which can be accounted for by preceding studies ^{2, 6} with the help of equations (7 - 9). Thus the formulas (4 - 6) are necessary to determine the threshold energy ETh by a balance equation of Coulomb repulsion and oscillator model for strong interactions.

In the domain of the resonance energy E_res we obtain:

TABLE 1: Parameters for the evaluation of equations (7 – 9).

In the domain of the transition to the asymptotic behavior of the total nuclear cross-section the following formula is applicable:

In order to give a qualitative motivation of the new calculation method we consider again Figure 1. At first, we have to repeat that Q_^tot(E) is not only restricted to proper nuclear reactions of protons with release of secondary particles: Thus for proton energies E > E_Th (the threshold energy E_Th amounts to 7 MeV) we can verify a rapid increase of Q^tot up to a resonance maximum E_res = 20.12 MeV. This behavior up to the environment of the resonance maximum can be described by a proper Gaussian distribution. What happens in this domain? Elastic scatter of proton at the nuclear potential is dominated by strong interaction and mediated mainly by mesons, if the quantum state of the nucleus is not changed. This implies that in order to satisfy energy-momentum relation only the whole nuclear adopts energy and momentum (kinetic energy), the impinging proton is slightly deflected. If the nucleus is also excited by vibrations, rotations or transitions to an excited configuration, then the whole process is inelastic and by emitting γ-quanta it is damped to finally return to the ground state. These effects are mainly described by the Breit-Wigner formula and its generalization.^{6 - 8} It has to be pointed out that the secondary proton under these conditions is still the primary proton, which is deflected by a slightly higher scatter angle compared to Molière multiple scatter theory.⁹

New calculation formulas for Q^tot(E) and E_Th

Figure 1 presents the total nuclear cross-section of oxygen; we can verify that, besides the Gaussian distribution of the environment of the maximal value E_max, after a slower decrease of Q^tot, the asymptotic behavior is reached (this is certainly valid for proton energies E < 300 MeV). The asymptotic branch can be represented by a suitable error function erf(E), which has to satisfy some further boundary conditions. This error function erf(E) results by an integration over Gaussian resonance distributions of the energy within finite boundaries. The transition from the domain Emax to the asymptotic domain can either be represented by a sum of exponential functions or by a Gaussian distribution with an additional energy shift. In every case, this consideration indicates that there exists an alternative representation of the total nuclear cross-section Q^tot(E) besides the previously studied one according to equations (7 – 9). Thus equation (10) provides the new formula for Q^tot(E), which appears to lead to a better access to quantum mechanical resonance mechanisms expressed by harmonic oscillators. Equation (11) provides all terms necessary for the determination of boundary conditions of the whole problem, and equation (12) an alternative calculation procedure for E_Th. A_boundary has the purpose to ensure that Q^tot(E) assumes zero at E = E_Th, since the Gaussians do not vanish at this position. A modification will be accounted for, when we shall include the tunneling effect.

Table 2 presents the corresponding parameter values to perform calculations of equation (10) with the help of equation (11). All necessary parameters of equations (10, 11) can be determined via Table 2 and the formula:

TABLE 2: Parameters of the new cross-section formula (10) and equations (11) using equation (10a).

The total nuclear cross-section requires also the threshold energy E_Th, which can be calculated by the following formula (12), which is easier to handle than formulas (4 – 6):

The parameters of formula (12) are: C = 6.565304 MeV, p = -0.10368 and q = 0.00481, D = -1.2889, p1 = -0.6597, q1 = -0.6601.

TABLE 3: Calculation of the threshold energy E_Th according to formula 4 – 6 (B) and present formula 12 (A).

Calculation procedure of the energy levels of excited states and quantum mechanical tunneling effect for E < E_Th

The investigations referring to the tunneling effect are as old as the quantum mechanics itself, since this theory has been used to explain the α-decay of heavy nuclei by Gamow. In many textbooks of quantum mechanics the tunnel effect is treated in detail of a particle passing with energy E through a potential box V with E < V. In the meantime, some other important tunnel effects have been studied, e.g. Josephson junctions and, recently, the H bonds between the complementary DNA strands through a double minimum potential.¹⁴ [Also, see reference list in Ulmer et al. ¹⁴]. After this digression, we intend to return to Figure 3, which refers to the nuclear potential of O₈ ¹⁶, but all formulas developed here can be applied to other nuclei, which exhibit similar properties of their nuclear potential. In following, we need the deflection point ξ of this figure. If we reduce formula (3) to one Gaussian with V₁ = 0, this point is easy to determine, since vanishing of the second derivation provides ξ² = σ₀²/2. The deflection point of the nuclear potential according to formula (3) can only be obtained by an iteration procedure:

In the next step, we have to solve a quadratic equation, too, since we insert the solution with ξ₀ into the Gaussians of equation (13) and determine ξ₁² according this equation by putting ξ₁² = x², i.e. the following step is given by:

In principle, this iterative procedure can be repeated, but stopping after step 2 is very sufficient. We now determine the energy levels of the potential type like that of formula (3). It offers to approximate the potential by a 3D harmonic oscillator. However, this is not sufficient, since the energy levels of the excited states are not equidistant.² Therefore we use a nonlinear and nonlocal Schrödinger equation, which incorporates nuclear interactions as a self-interacting field. With the help of equation 3 this generalized Schrödinger equations assumes the form:

The coupling constant λ has to be chosen such that the dimension agrees on both sides, but it can be put λ = 1; the magnetic interaction, i.e. spin-orbit-coupling, can be added using principles previously given.^{2, 3} With the help of the generating functions of Hermite polynomials we are able to write equation (15) in the form:

The equation above represents a highly inharmonic oscillator equation of a self-interacting field. Since the square of the wave-function is always positive definite, all terms with odd numbers of n1, n2, and n3 vanish due to the anti-symmetric properties of those Hermite polynomials. For r_c ≤ ξ (domain with positive curvature), the whole equation is reduced to a harmonic oscillator with self-interaction; the higher-order terms are small perturbations:

The solutions of this equation are those of a 3D harmonic oscillator; the classification of the states by SU₃ and all previously developed statements with regard to the angular momentum are still valid. The only difference is that the energy levels are not equidistant; this property can easily be verified in one dimension. The usual ground state energy is ћω₀/2. This energy level is lowered by the term ~Φ_0,0,0, depending on the ground-state wave-function. The energy difference between the ground and the first excited state amounts to ћω0; this is not true in the case above, since the energy levels of all excited states depend on the corresponding eigen-functions themselves (these are still the oscillator eigen-functions!). Next, we will include the terms of the next order, which are of the form ~ λ∙(Φ_0,2,2, Φ_2,2,0, Φ_2,0,2):

The additional term T represents tensor forces. The whole problem is still exact soluble. In further extensions of the nonlinear/nonlocal Schrödinger equation, we are able to account for spin, isotopic spin, and spin-orbit coupling. The spin-orbit coupling, as an effect of an internal field with nonlocal self-interaction, is plausible, since the extended nucleonic particle has internal structure; consequently, we have to add H_so to the nonlinear term:

Ψ is now (at least) a Pauli spinor (i.e. a two-component wave-function), and due to the term H_so the SU₃ symmetry is broken. We should like to point out that the operation grad φ acts on the Gaussian kernel K:

The expression in the bracket of the previous equation represents a vector, and p (p → -iћgrad) acts on the wave-function. Since the neutron is not a charged particle, the spin-orbit coupling of a neutron can only involve the angular momentum of a proton. In nuclear physics, these nonlinear fields are adequate for the analysis of clusters (deuteron, He, etc.). The complete wave-function Ψ_c is now given by the product of a function in configuration space Ψ multiplied with the total spin and isotopic-spin functions.

A further access to solve the above nonlinear/nonlocal equation is obtained by the Ritz’ variation principle applied to energy minimum. For this purpose, we start with the following set of functions:

With the help of the function set (21) all integrals can be evaluated analytically. By that, the task remains to determine the coefficients of the function set by an iterative procedure, since equation (23) yields a third order equation of all coefficients forming the wave-function. However, this so-called self-consistent field method is very familiar in many-particle problems, and, therefore, a detailed description is superfluous. The performance of the above task provides, besides the ground state energy, excited states for E < E_Th and virtually excited states, if E > E_Th. However, we are interested in the role of excited states of a nucleus, when a proton can penetrate the potential wall for E < E_Th and virtually excited states are only important with regard to nuclear reaction. Therefore we now turn our interest to the tunneling effect of external protons in connection with potential types like Figure 3.

The quantum mechanical tunneling effect must be treated in three dimensions. We denote by μ the reduced mass ‘proton – nucleus’, then the Schrödinger equation of the problem reads:

The boundary condition for the proton energy is 0 < E < E_Th. Since a closed analytical solution is unknown, we shall make use here of the Dyson series; this is similar to the Feynman propagator method, and easy to manipulate for plain waves and Gaussians (kernels, potentials, etc.). In absence of external perturbation, the time evolution operator reads:

Using the interaction picture of quantum mechanics and the definition of U (t, t₀) according to equation (25), the interaction term can be transformed to the expression:

Within the restriction of our task, namely that H₀ is the free-particle operator and H₁ a linear combination of 2 Gaussian functions (see equation (25)), the expansion (28) is easy to solve: The unitary operator U related to H₀ applied to the terms of H₁ again provides Gaussian terms. Since the operator U also acts on plain waves, too, the similar behavior is valid. Thus the following expressions have to subject to iterations (the factor α has to associate with the proper dimension):

With respect to the Dyson expansion the terms according to equations (29, 29a) are very convenient, since iterations always yield terms of the same structure. Please note that the kinetic energy term of the proton assumes complex values, if the impinging energy is lower than the potential wall, i.e. the position probability within the positive, repulsive potential suffers damping. Some consequences will be given in the result section.

Results

The following Figures 4 – 7 deal with the nuclei of carbon, oxygen, calcium and copper. They show that both calculation formulas provide equivalent results; we use the following abbreviations: FM (former method) and PM (present method). However, the chosen type of functions used in this communication offers new ways to analyze nuclear cross-sections due to the advantages of Gaussian functions and kernels.

FIG. 4: Nuclear cross-section of the interaction proton – carbon nucleus; solid line: PM, dots: FM.

FIG.5: Nuclear cross-section of the interaction proton – oxygen nucleus; solid line: PM, dots: FM.

FIG. 6: Nuclear cross-section of the interaction proton – calcium nucleus; solid line: PM, dots: FM.

FIG. 7: Nuclear cross-section of the interaction proton – copper nucleus; solid line: PM, dots: FM.

The following Figures 8 – 11 show with regard to the most important nuclei in radiotherapy with protons roundness in the environment of the threshold energy E_Th and the previously assumed condition Q^tot(E) for E < E_Th does not hold due to the quantum mechanical tunneling effect. This roundness can be amplified by the energy spectrum of the impinging protons resulting from the beam-line and a jump of the fluence decrease of primary protons Φ_pp at E < E_Th is prevented (see also Figure 2 with regard to the passage of protons through water). The passage of protons through other media, e.g. calcium, leads to a similar roundness.³ By that, the quantum mechanical tunneling effect can be best studied by really mono-energetic protons with initial energy E < E_Th.

FIG. 8: Proton – nucleus (oxygen) interaction by proton energies less than 7 MeV via tunnel ling through potential wall according to Figure 4.

FIG. 9: Comparison of the tunneling probability for oxygen (solid line), carbon (dots) and calcium (dashes).

Both Figures 8 and 9 also present a way to check ranges of protons in biologically significant media, since at the end of the particle track very specific nuclear interactions occur with the help of the quantum mechanical tunneling effect. The intermediary existence of the isotopes N₇¹³ and S_c21⁴¹ provides a tool to study proton interactions for very low proton energies in that media. The quantum mechanical tunnel effect represents the necessary basis for these studies.

FIG. 10: Total nuclear cross-section of oxygen with in inclusion of the tunneling effect.

FIG. 11: See Figure 10 in the low energy domain with proton energy < 14 MeV.

Discussion

Let us first consider the low energy domain, i.e. up to the resonance maximum E_res. Thus for E < E_Th only due to the tunneling effect the projectile proton can enter the interior of the potential wall according to Figures 9 - 11. This provides in the case of oxygen the isotope F₉¹⁷ as intermediary excited state (or N₇¹³ and S_c21⁴¹ with respect to carbon or calcium). Since E < E_Th the lifetime of the excited states depends on the transition probability to a lower energy state to finally produce a γ-quant and a neutron by an exchange interaction (Pauli principle) within the nuclei via a π⁺-meson to yield a neutron, which can easily escape the potential wall. The resulting processes are also described by equations (30 – 32).

These equations are also valid for comparably low proton energies, but with E > E_Th. However, we are now located in the energy domain E > E_Th up to the resonance domain E_res. This domain, which is referred to as Breit-Wigner resonances, is characterized by various excitations, e.g. rotations of the whole nuclei, vibrations via proper deformations of the nuclei and excitations to excited configurations. All these processes are damped by emission of photons with proper energies. It should also be pointed out that these resonances belong to the category of inelastic scatter, but the resulting secondary proton is identical with the primary proton. Only the lateral scatter and energy losses are different to Molière scatter processes. If these scattering processes occur without nuclear excitations, rotations and vibrations, then they are referred to as elastic nuclear scatter of the protons, which are also contained in the Breit-Wigner resonance formula and in its generalization ⁷, given by Flügge. (The original Breit-Wigner formula ⁶ is restricted to scatter processes induced by so-called ‘s-states’ of the S-matrix; the inclusion of ‘p-states’ and states of higher order have later been accounted for 7, 8). The very often used way of notation ‘inelastic nuclear scatter of protons’ is not quite correct, since elastic scatter processes of protons at nuclei even occur in the asymptotic domain of Q^tot(E).  A very important nuclear process occurring in the same energy domain must be mentioned again, namely the exchange interaction of the proton with the mesons of the nucleus, i.e. we consider the reaction:

The recoil nucleus F₉¹⁶ undergoes ß⁺-decay of electron capture (EC) to become finally again O₈¹⁶. With regard to nuclear reactions of therapeutic protons the cross-section of oxygen is certainly most important. However, the interaction with calcium (bone) and proteins (carbon) seems also to be worthy of interest. Below we present lists of the most important reactions with carbon and calcium. For brevity, we do not state the decay reactions of the isotopes (detailed listings can be found in web), but for most of them electron capture (EC) and ß⁺-decay with emission of γ-quanta are preferred reaction channels. For these reasons we have state the symbols X in equations (31 – 31g) and equations (32 – 32h).

As already verified ^{2 – 3} in the case of oxygen, the nuclear reactions presented in the listings 4.1.1 and 4.1.2 depend on the available or residual proton energy. For this purpose we consider the domain E > E_res up to the beginning of the asymptotic behavior, i.e. the proton energy amounts to ca. 100 MeV. We refer to this energy branch of Q^tot(E) as E_temperate. (In the case of copper, the asymptotic behavior is reached at a somewhat higher energy (ca. 150 MeV), but this element is mainly of importance with regard to the determination of beam-line properties). Thus in the energy domain E_temperate the nuclear reactions according to relations (31 – 31c) and (32 – 32d) preferably occur, whereas for E > E_temperate or E >> E_temperate the release of clusters according to relations (31d – 31g) and (32e – 32h) are more probable. A typical case of a cluster release is the α-particle as a secondary particle. According to ³ this threshold energy amounts to 100 MeV. In the case of carbon it is increased to be 101 MeV, for calcium it amounts to 98 MeV and for copper to 97 MeV. At about E = 190 MeV the probability of this release of α-particles as outcome of nuclear reactions vanishes for all cases considered here. The main difference, however, exists in the yield of this secondary particle. We normalize total yield of α-particles to ‘1’ for the reaction of proton with oxygen. Then in the case of carbon we obtain a total yield of 0.97. For calcium it amounts to 1.68 and for copper to 1.93. Therefore it is generally correct that mainly the yield of clusters is increasing, if the nuclear charge Z and mass number A_N correspondingly increase. A previous check of GEANT4 ¹⁰ in the papers ^{1–3, 13, 15} revealed that this Monte-Carlo code provides default reaction channels, which underrate cluster formations. However, since this Monte-Carlo code is an open system, it can be improved by the user. A further lack of this Monte-Carlo code is that the release of neutrons is not sufficiently accounted for, as we could verify at a glance at the papers.^{11, 12} According to the present results the release of low energy neutrons should be much more in focus in dose calculations due to the rather high RBE.

The calculation method for the total nuclear cross-section Qtot(E) provides a significant advantage compared to the previously published method ^{2, 3} : Each Gaussian distribution and the error-function distribution with a specific energy shift can be associated to probability distributions for the occurrence of some nuclear reaction types. Thus the energy shifts now refer to threshold energies with regard to the corresponding maxima and the half-breadths represent a measure for the yield of some reaction types. Furthermore, it is possible to analyze measured curves of Q^tot(E) by suitable deconvolution procedures as worked out in.¹⁶ These convolution methods are developed with regard to linear combinations of Gaussian kernels and shifts. Since there exists a connection between the number of kernels and the underlying statistics, which is strictly a non-relativistic Boltzmann distribution in the case of one Gaussian kernel and a linear combination of different Gaussian kernels with further shifts in the case of the Fermi-Dirac statistics, we have also obtained a way to interpret nuclear reactions by the operator formulation of Fermi-Dirac statistics.^{2, 16} Such statistical models can help to calculate nuclear interaction processes in a simpler way compared to many-body-problems of relativistic quantum mechanics. The evaluation of measurement data in clinical practice of proton therapy, e.g. the contributions of secondary particles, requires relatively simple but reliable methods.

Listings of possible nuclear reactions of some nuclei of particular interest

Nuclear reaction types of protons with carbon and calcium (with oxygen already reported in ^{1 - 3}).

Interactions proton – carbon:

Interaction proton – calcium:

Please note that X stands for further nuclear decay reactions, which preferably incorporate β⁺, β^- decay and electron capture (EC). In particular, the emission of positrons is associated for further radiation processes by pair annihilations.

Conflict of Interest

The authors declare that they have no conflicts of interest. The authors alone are responsible for the content and writing of the paper.

References

1. Ulmer W. Theoretical aspects of energy range relations, stopping power and energy straggling of protons. Radiat Phys Chem 2007; 76:1089. [CrossRef]
2. Ulmer W, Matsinos E. Theoretical methods for the calculation of Bragg curves and 3D distribution of proton beams. Eur Phys J Special Topics 2010; 190: 1-81. [CrossRef]
3. Ulmer W, Matsinos E. A Calculation Method of Nuclear Cross-Sections of Proton Beams by the Collective Model and the Extended Nuclear-Shell Theory with Applications to Radiotherapy and Technical Problems. Journal of Nuclear and Particle Physics 2012; 2: 42-56. [CrossRef]
4. Chadwick MB, YoungPG. Los Alamos National Laboratory Report LAUR-96-1649, 1996.
5. Berger MJ, Coursey JS, Zucker MA, ESTAR, PSTAR and ASTAR. Computer Programs for calculating Stopping-Power and Range Tables for Electrons, Protons and a-particles (version 1.2.2) National Institute of Standards and Technology, Gaithersburg, MD 2000.
6. Segrè E. Nuclei and Particles. (Benjamin, New York) 1964.
7. Flügge S. Erweiterung der Breit-Wigner Formel auf höhere Quantenzustände. Z Naturforschung 3a, 1948; 97–110.
8. Ramond T. Semi-classical study of quantum scattering on the line. Communications in Mathematical Physics 1996;177: 221-54. [CrossRef]
9. Molière GZ. Theorie der Streuung schneller geladener Teilchen III. Z Naturforschung 1955; 177 – 200.
10. GEANT4 Documents 2005; Available from:
    http: //geant4.web.cern/geant4/G4UsersDocuments /Overview/ html/
11. Paganetti H. Nuclear interactions in proton therapy: dose and relative biological effect distributions originating from primary and secondary particles. Phys Med Biol 2002; 47: 747 - 64. [CrossRef]
12. Boon SN. Dosimetry and quality control of scanning proton beams. PhD Thesis, Rijks University Groningen, 1998. 
13. Gottschalk B, Platais R, Paganetti H. Nuclear interactions of 160 MeV protons stopping in copper: a test of Monte Carlo nuclear models. Med Phys 1999; 26: 2597– 601. [CrossRef]
14. Ulmer W, Cornelissen G. Coupled Electromagnetic Circuits and Their Connection to Quantum Mechanical Resonance Interactions and Biorhythms. Open Journal of Biophysics 2013; 3: 253-74. [CrossRef]
15. Ulmer W, Schaffner B. Foundation of an analytical proton beamlet model for inclusion in a general proton dose calculation system. Radiation Physics and Chemistry 2011; 80: 378–402. [CrossRef]
16. Ulmer W. Inverse problem of a linear combination of Gaussian convolution kernels (deconvolution) and some applications to proton/photon dosimetry and image processing. Inverse Problems 26, 2010; 085002. [CrossRef]