Analytical computation of prompt gamma ray emission and detection for proton range verification.

Institute of Physics and Engineering in Medicine Phys. Med. Biol. 60 (2015) 4915–4946

Physics in Medicine & Biology doi:10.1088/0031-9155/60/12/4915

Analytical computation of prompt gamma ray emission and detection for proton range verification E Sterpin1, G Janssens2, J Smeets2, François Vander Stappen2, D Prieels2, Marlen Priegnitz3, Irene Perali4,5 and S Vynckier1 1

Université catholique de Louvain, Center of Molecular Imaging, Radiotherapy and Oncology, Institut de recherche expérimentale et clinique, Avenue Hippocrate 54, 1200 Brussels, Belgium 2 Ion Beam Applications SA, Louvain-la-Neuve, Belgium 3 Helmholtz-Zentrum Dresden-Rossendorf, Institute of Radiation Physics, PO Box 510119, 01314 Dresden, Germany 4 Politecnico di Milano, Dipartimento di Elettronica, Informazione e Bioingegneria, Milano, Italy 5 Istituto Nazionale di Fisica Nucleare, Sezione di Milano, Milano, Italy E-mail: [email protected] Received 21 April 2015 Accepted for publication 11 May 2015 Published 9 June 2015 Abstract

A prompt gamma (PG) slit camera prototype recently demonstrated that Bragg Peak position in a clinical proton scanned beam could be measured with 1–2 mm accuracy by comparing an expected PG detection profile to a measured one. The computation of the expected PG detection profile in the context of a clinical framework is challenging but must be solved before clinical implementation. Obviously, Monte Carlo methods (MC) can simulate the expected PG profile but at prohibitively long calculation times. We implemented a much faster method that is based on analytical processing of precomputed MC data that would allow practical evaluation of this range monitoring approach in clinical conditions. Reference PG emission profiles were generated with MC simulations (PENH) in targets consisting of either 12C, 14N, 16O, 31P or 40Ca, with 10% of 1H. In a given geometry, the local PG emission can then be derived by adding the contribution of each element, according to the local energy of the proton obtained by continuous slowing down approximation and the local composition. The actual incident spot size is taken into account using an optical model fitted to measurements and by super sampling the spot with several rays (up to 113). PG transport in the patient/camera geometries and 0031-9155/15/124915+32$33.00 © 2015 Institute of Physics and Engineering in Medicine Printed in the UK

4915

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

the detector response are modelled by convolving the PG production profile with a transfer function. The latter is interpolated from a database of transfer functions fitted to MC data (PENELOPE) generated for a photon source in a cylindrical phantom with various radiuses and a camera placed at various positions. As a benchmark, the analytical model was compared to MC and experiments in homogeneous and heterogeneous phantoms. Comparisons with MC were also performed in a thoracic CT. For all cases, the analytical model reproduced the prediction of the position of the Bragg peak computed with MC within 1 mm for the camera in nominal configuration. When compared to measurements, the shape of the profiles was well reproduced and agreement for the estimation of the position of the Bragg peak was within 2.7 mm on average (1.4 mm standard deviation). On a non-optimized MATLAB code, computation time with the analytical model is between 0.3 to 10 s depending on the number of rays simulated per spot. The analytical model can be further used to determine which spots are the best candidates to evaluate the range in clinical conditions and eventually correct for over- and under-shoots depending on the acquired PG profiles. Keywords: prompt gamma, range monitoring, Monte Carlo (Some figures may appear in colour only in the online journal) I. Introduction Proton range uncertainties may lead to delivered dose distributions that deviate significantly from planned dose distributions with negative clinical outcomes (Paganetti 2012). The treatment planning process can cope to some extent with range uncertainties (robust treatment planning), often at the cost of sub-optimal treatment planning quality. Range uncertainties can also be reduced by measuring indirectly the range of the protons in actual treatment conditions (Knopf and Lomax 2013). Prompt gamma (PG) imaging is a promising technique for on-line monitoring of the proton range in the context of spot scanned beams. However, such technology requires for clinical implementation an integration of efficient modelling tools and hardware. PGs emitted by protons along their track in patient geometry can be detected by appropriate camera designs (Lu 2008, Peterson et al 2010, Testa et al 2010, Bom et al 2012, Min et al 2012, Smeets et al 2012, Golnik et al 2014, Pinto et al 2014, Verburg and Seco 2014). Smeets et al (2012) proposed a slit camera geometry that allows for the detection of shifts of Bragg peak positions in homogeneous targets with an accuracy of 1–2 mm (one standard deviation, considering an actual match of expected and measured Bragg peak positions). Measured PG detection profiles with this camera design show typically a sharp fall-off that is well correlated with the position of the Bragg peak. Measurements and Monte Carlo (MC) simulations were in very good agreement (within 5%) for the PG contribution to the total signal detected (Smeets et al 2012). Therefore, MC computed profiles could be used as the reference to compare with measured profiles and estimate the shifts between expected and actual Bragg peak positions, if any. Shifts between expected and actual ranges can then be derived by acquiring the image of the patient and determining the point of entrance of the beam using for instance an optical surface-imaging system (Schöffel et al 2007). 4916

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 1. Schematic description of the prompt gamma imaging system. Emitted prompt

gammas are collimated by a knife-edge slit collimator and impact LYSO crystals. The nominal magnification factor is 0.8.

However, the requirements are challenging, both for the measurement and the calculation. On the one hand, classical gamma cameras are not suited for the measurement of high-energy gammas at a very high-count rate with an important neutron background and, on the other hand, MC methods are too consuming for calculation of the expected signal in clinical routine. The first issue has been recently cleared with the report of measurements in spot scanning mode at clinical dose rate (Perali et al 2014). The present paper in turn aims at demonstrating the feasibility of the calculation of the expected PG signal acquired by this prototype with sufficient speed for clinical routine. The results will allow practical evaluation of this range monitoring approach in clinical conditions. The detection of a PG by the camera is the result of two main processes: (1) the emission of PGs by protons undergoing nuclear reactions; (2) the transport of the PG through patient, collimator and detector geometries. We propose to model both aspects analytically using extensive MC precomputed data. The applicability of the model was tested on realistic clinical configurations and benchmarked with MC simulations and experiments. II. Material and methods II.A. PG camera prototype

The current camera design was described in Perali et al (2014) and is schematically illustrated in figure 1. The full-size prompt gamma camera prototype consisted of a photodetection system and a knife-edge slit tungsten collimator mounted on a dedicated trolley positioning system. The photodetection systems rely on two rows of 20 LYSO crystal slabs (each slab is 4 mm wide along beam axis, 10 cm high and 3.15 cm deep resulting in 504 cm³ for the 40 slabs) with independent SiPM light readout on one extremity. The field-of-view of the camera can be adjusted by modifying the distance from centre of collimator to the centre of detector crystals (d coll − detectors) and/or the distance from the centre of the collimator to the beam axis (d collimator). In the nominal configuration, the crystals cover a 10 cm field-of-view on the beam axis (magnification factor (M) of 0.8). The energy detection window of the camera is set to 3 to 6 MeV. 4917

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

II.B. MC engine for proton simulation (PENH)

The MC engine used in this study is the recently validated extension of PENELOPE (Salvat et al 2011) to protons called PENH (Salvat 2013, Sterpin et al 2013, 2014). In PENH, nuclear reactions are simulated from ICRU 63 (ICRU 2000) cross sections for 12C, 14N, 16O, 31P and 40Ca, which are the only elements (with 1H) that have a weight contribution of more than 1% in human tissues. Secondary protons are simulated in detail. Neutrons are neglected (the impact of this approximation is discussed in the discussion later and appendix A). Other secondary charged particles are simulated as protons with corrected energy to ensure consistency of their range (see Sterpin et al 2013 for more details). Elastic and inelastic electromagnetic collisions are described by a mixed (class II) algorithm (Salvat et al 2011, Salvat 2013). Spatial displacements are generated by the random hinge method (Salvat et al 2011). The user-input parameters influencing proton transport were C1 = C2 = 0.2, Wcc = 10 keV and EABS = 1 MeV. EABS is the absorption energy cutoff. C1 and C2 control the maximum energy loss and angular deflection per step. Wcc provides the minimum energy that can be transferred to an electron. For electrons, positrons and photons, similar parameters were defined. The definition of these parameters may be found in the PENELOPE reference manual. For electrons and positrons, EABS = 100 keV, C1 = C2 = 0.2, Wcc = 10 keV and Wcr = 10 keV (energy cut-off for bremsstrahlung emission), unless otherwise mentioned. For photons, EABS = 10 keV. II.C. MC modelling of spot scanning delivery, CT and PG camera geometries

Spot scanning delivery is modelled according to the optical model of Grevillot et al (2012). The model is implemented in a dedicated main program that generates the spatial and energy characteristics of the spots. The main program is called 4P (extension of Penelope to Protons for Pencil beam scanning treatment Planning). The optical model used in the present study is based on the beam data library of Gantry 4 in the West German Proton therapy centre Essen (WPE). 4P integrates also PENCT extended to protons. Thus, it is possible to track protons in CT geometries. Finally, variance reduction techniques are used in 4P to speedup the PG generation process. The probability of a nuclear event is increased by a factor of 10. Moreover, at each nuclear event, the number of emitted PGs is multiplied by 50 (energy and direction of the PG are resampled 50 times). The weights of the PGs and other particles involved are corrected accordingly to ensure the consistency of the simulation, according to the recipes provided by Salvat et al (2011). The PG camera is modelled using the PENGEOM package and simulated simultaneously with the CT geometry. For all simulations shown in the document and appendixes, the number of protons simulated equalled 107 with the variance reduction settings as described in section II.C, leading to an effective simulated protons of 5 × 109. 4P was benchmarked with MCNPX for PG emission and detection. The scope of the present paper being the design and benchmark of the analytical model, details on the benchmark of 4P are provided in the appendix A. The choice of 4P, over MCNPX, Geant4 and FluKa, was not only based on past demonstrated performances, but also on very practical reasons (deep knowledge of the source code in order to tailor it in an optimal way for this specific application; easy integration and update of tabulated prompt gamma emission data; already integrated capability to simulate CT geometries; already integrated and benchmarked beam model for proton pencil beam scanning). Overall, 4P is capable to generate PG detection profiles in close agreement with MCNPX that lead to the same 1–2 mm accuracy for detecting shifts of Bragg peak position as described in Smeets et al (2012). The fact that neutrons are not simulated with 4P does not impact significantly the estimation of Bragg peak position 4918

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

for the considered camera design because neutrons are contributing as a flat background to the acquired signal. Quantitative information and additional justifications for using 4P are provided in appendix A. II.D. Analytical model

The main idea is that PG production profiles in depth can be treated almost the same way as dose profiles in a pencil beam dose calculation algorithm. The analytical model is based on the following steps that will be discussed below: • Generation of a proton ray according to the beam model • Curvilinear raytracing in the CT geometry to determine the local composition and the energy of the proton in the continuous slowing down approximation (CSDA). • At the corresponding equivalent depth, look-up in MC precomputed PG emission profiles for various nuclear species • Addition of each individual contribution according to local atomic composition • Convolution of the PG emission profile with the response of the camera to obtain the PG detection profile. The response of the camera is determined from an interpolation within a database of fits to precomputed MC simulations of the response of the camera for various configurations (various patient thicknesses, camera positions, etc). II.D.1. Precomputation of PG emission. Using precomputed PG emission profiles in depth for various nuclear species, one could predict the effective PG emission at a given depth, thus at a given proton energy and for a given local material. As in Sterpin et al (2013), it is assumed that human tissue can be completely described by the nuclear species 1H, 12C, 14N, 16 O, 31P and 40Ca. PG emission profiles in depth were computed with PENH in homogeneous phantoms (60 × 60 × 60 cm3) for incident protons of energies between 40 to 230 MeV with 1 MeV bin. PG emission was scored over the surface perpendicular to the direction of the protons, as the PG camera does. The resolution in depth is variable and depends on the initial energy to keep the number of voxels constant, which helps post processing (around 0.3 mm and 3 mm for 40 MeV and 230 MeV, respectively). Only PGs with energies between 3 and 8 MeV were counted, PG of energies higher than 8 MeV being unlikely for the considered elements (below 4% of the PG emission cross section integrated from 3 MeV for protons energies below 100 MeV). The spectrum can be recovered afterwards using ICRU 63 cross-sections differential in energy for each element. Most tissues have a weight fraction of 1H of about 10% (ICRU 44), excepted cortical bone (3.4% (ICRU 44)). Although 1H does not contribute to the prompt gamma signal, it modifies the proton flux because of the proton-proton elastic interaction. Therefore, all PG production profiles were computed in 5 phantoms with 10% 1H and 90% of either 12C, 14N, 16O, 31P or 40 Ca with a forced mass density of 1 g cm−3 leading to 5 sets of precomputed profiles. II.D.2. Determination of local PG emission in CT geometry. Raytracing is performed to deter-

mine the voxels and their associated material that intersect the proton trajectory. In order to simulate the (Gaussian (Grevillot et al 2012)) energy spread of the incident proton beam, the mean proton energy at nozzle exit, plus two other energies (mean energy minus half full width half maximum and mean energy plus half full width half maximum) are simulated and then summed with the proper weight. Also, the spot size and shape are modelled as in 4P and sampled with several curvilinear rays (up to the arbitrary number 113). In order to mimic the 4919

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 2. Illustration of the geometry used for the MC simulation of photon detection by the PG camera. 5 variables were considered and are shown on the figure (energy, diameter, distance axis—collimator (d collimator), distance collimator—detectors (d coll − detectors) and position along axis). The impulse responses were computed for fixed sets of the 5 variables.

increase of the beam size along depth due to multiple scattering, the distance between each curvilinear ray and the spot axis varies isotropically as a function of the water equivalent depth, following the size of such spot in water. Thus, the materials intersected by the protons may be different depending whether the ray is sampled at the centre of the spot (straight trajectory), or from the border of the spot (‘curved’ trajectory). At a given depth, total equivalent path-lengths are computed for each element (leading to the ‘carbon’ equivalent path-length, the ‘nitrogen’ equivalent path-length, etc). To ensure consistency, the equivalent path-lengths were computed using the same stopping powers as PENH. The curvilinear raytracing assumes that protons follow a monotonously diverging trajectory, which is only an approximation of multiple scattering since the same proton can diverge and converge back to its initial trajectory. The actual PG production in tissue can then be derived by combining the PG production of each element according to its weight fraction and at each equivalent depth. This operation assumes that the proton flux in the considered material is similar to the proton flux in each element. The PG production profiles from the different curvilinear rays are weighted according to their respective distance from spot axis assuming a Gaussian dose profile and then combined. Local information on the spectrum of the produced PGs is obtained using the PG emission cross-sections differential in energy of ICRU 63 (from 3 to 8 MeV). The conversion from depth to local proton energy is computed from CSDA approximation. II.D.3. Transport of PGs from production site to camera detector. We propose to model the PG profile acquired by the detectors by a single convolution of the PG production with a transfer function with several parameters fitted with data generated by the MC code PENELOPE (Salvat et al 2011) with the same settings as defined in section II.A. The response of the camera was simulated with MC considering a photon source on the central axis of a cylindrical water phantom (see illustration in figure 2). Four variables were considered: (1) the energy of the photon source; (2) the distance of the collimator to the cylinder axis (d collimator); (3) the distance between the collimator and the detectors (d coll − detectors); (4) the radius of the cylinder. For these simulations, the source was placed right below the centre of the opening of the collimator. However, the position of the source on phantom axis with respect to the centre of the opening of the camera was also varied in order to validate the use of a single convolution in the computation of the detection profile. 4920

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

We observed that the impulse response from emission to detection by the camera can be approximated by a transfer function F(x) (where x is the position along beam line, with the origin placed at the projection of the centre of the opening of the collimator on the beam line) such that: F (x ) = a (Ro(x ) ⊗ (Gσ1(x ) + Gσ 2(x ))) + b,

where the five parameters are: the amplitude (a), the baseline (b), the opening (o) of the rectangular function Ro and the standard deviations (σ1 and σ2) of the Gaussian functions G. The ⊗ operator represents the convolution. Examples of fits to MC data computed for the five aforementioned parameters are shown in figure 3. These five parameters were computed for different configurations by optimizing the fit between the function F(x) and the corresponding impulse response computed with MC. The fit was optimized using the least-square method. For each parameter, a table covering the range of configurations that were simulated was built. The appropriate value of the parameter can then be interpolated (range of values in the table are shown in brackets) using the actual energy of the PGs (3 to 8 MeV by step of 1 MeV), the actual distance between beam and collimator (5 to 35 cm, by step of 2.5 cm), the actual distance between collimator and detectors (12 to 20 cm by step of 4 cm) and the integral density crossed by the gammas, perpendicular to beam axis towards the centre of the detectors (0 to 30 radiological cm by step of 5 cm, the radiological path being the integral electronic density (relative to water) along the photon ray path). Eventually, the detection profile is computed by convolving the PG emission profile with the impulse response F(x). Because both the PG emission and the impulse response are expressed along the beam axis, the convolution is straightforward. The final result is also expressed onto beam-axis. To obtain the same result but this time in the coordinates of the detectors, an inversion of the profile and a correction by the magnifying factor are necessary. A single convolution is an approximation. Strictly speaking, there is a dependency of the impulse response with the position of the emitted photon, as shown in figure 3(e). Using only the response function for the centralized source leads to an overestimation of the amplitude of the response for outer detectors. Figure 3 displays fits of MC data generated for the transfer function modelling the detection by the PG camera of a photon produced in a water cylinder. 5 parameters were independently varied: (1) the energy of the photon (a); (2) the radius of the cylinder (b); (3) the distance from the centre of the collimator to the axis of the phantom cylinder (c); the distance from the centre of the collimator to the centre of the detectors; and (5) the position of the source along the cylinder axis (e). The origin on the abscissa corresponds to the projection of the centre of the camera to the beam axis. In figure 3(e), the MC computed profiles were shifted to the origin to superimpose and compare the fit and the MC profiles. The fit was computed based on MC data generated for the source aligned with the centre of the collimator and the detectors. II.D.4. Implementation in MATLAB. The analytical model has been programmed in MAT-

LAB within a graphical extension of a MATLAB software called REGGUI and developed in-house. Among other features, REGGUI can import images and contours and process them for various applications. The graphical user interface, called ‘PG simulator’ (see figure 4), can import CT images, phantom geometries and RT plans describing spot scanning treatments. For the latter, a spot can be selected for PG simulation graphically among the imported spots. The configuration of the camera can also be defined. The software processes the data in order to start 4P and analytical simulations. Measured data can also be imported. Analysis of the results is performed 4921

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

1.2 Energies in MEV

6 7

1

5 8

0.8 0.6

Impulse response for various phantom radiuses (5 MeV photon energy, 25 cm from collimator to beam axis, 20 cm from collimator to detector) −4 x 10 1.6 MC Fit 1.4 Yield of detected prompt gammas (counts/emitted prompt gamma)

Yield of detected prompt gammas (counts/emitted prompt gamma)

Impulse response for various photon energies (25 cm from collimator to beam axis, 20 cm from collimator to detector, 10 cm phantom radius) −4 x 10 1.6 MC Fit 1.4

4

0.4 0.2

1

1.2 Radiuses in cm 1

5 10 15

0.8

20

0.6 0.4 0.2

3

0 −60 −50 −40 −30 −20 −10 0 10 20 Beam axis (mm) (a)

30

40

50

0 −60 −50 −40 −30 −20 −10 0 10 20 Beam axis (mm) (b)

60

1 0.8 0.6

20

Distances in cm 25

30 35

0.4 0.2

40

50

60

Impulse response for various distances from collimator to detector (5 MeV photon energy, 25 cm from collimator to beam axis, 10 cm phantom radius) −4 x 10 1.6 MC Fit 1.4 12 Yield of detected prompt gammas (counts/emitted prompt gamma)

1.2

1.2

Distances in cm

1

16

20

0.8 0.6 0.4 0.2

0 −60 −50 −40 −30 −20 −10 0 10 20 Beam axis (mm) (c)



Impulse response for various distances from collimator to beam axis (5 MeV photon energy, 20 cm from collimator to detector, 10 cm phantom radius) −4 x 10 1.6 15 MC Fit 1.4

30

30

40

50

60

0

−50 −40 −30 −20 −10 0 10 20 Beam axis (mm) (d)

30

40

50

Impulse response for various position along beam axis (5 MeV photon energy, 25 cm from collimator to beam axis, 20 cm from collimator to detector, 10 cm phantom radius) −4 x 10 1.6 MC centred source Fit 1.4 MC 1 cm off centre MC 2 cm off centre MC 3 cm off centre 1.2 MC 4 cm off centre MC 5 cm off centre 1 0.8 0.6 0.4 0.2 0 0

10

20 30 Beam axis (mm) (e)

40

50

Figure 3. Fits of MC data generated for the transfer function modelling the detection

by the PG camera of a photon produced in a water cylinder. 5 parameters were independently varied: 1) the energy of the photon (a); 2) the radius of the cylinder (b); 3) the distance from the centre of the collimator to the axis of the phantom cylinder (c); the distance from the centre of the collimator to the centre of the detectors; and 5) the position of the source along the cylinder axis (e). The origin on the abscissa corresponds to the projection of the centre of the camera to the beam axis. In figure (e), the MC computed profiles were shifted to the origin to superimpose and compare the fit and the MC profiles. The fit was computed based on MC data generated for the source aligned with the centre of the collimator and the detectors. 4922

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 4. Snapshot of the graphical user interface ‘PG simulator’ for starting MC (4P) and analytical simulations of prompt gamma emission and detection. 3D distributions of MC simulated prompt gamma emission are also illustrated on top of the CT images.

through visual inspection of the generated data or using more quantitative estimations as discussed in the following sections. II.E. Benchmark of the analytical model II.E.1. General philosophy. As mentioned before, the purpose of the analytical model is to compute a reference signal for comparison with the equivalent measurement in clinical conditions. Discrepancies with measurements may have several sources of uncertainties, among which: (1) the stopping powers (mean excitation energy) (Paganetti 2012); (2) the PG cross sections (Verburg et al 2012); (3) the modelling of the treatment unit; (4) residual uncertainties in the positioning of the camera; (5) the modelling of the camera response; (6) the transport of protons in the patient and photons in CT geometry and camera. 4P and the analytical model evenly share the first 4 sources of uncertainties, that are independent of the proton and photon transport algorithm, which is not the case for the latter two ones. The validation of the analytical model with MC simulations checks that the approximated particle transport in the patient and simplified modelling of the camera do not jeopardize the capability to detect accurately shifts of Bragg peak positions. In case of a successful comparison between 4P and the analytical model, we ensure that any remaining discrepancy between the measurements and the analytical model are caused by inaccuracies in physical data, knowledge of patient anatomy and in the beam model that are the same as 4P. Research is needed to reduce these uncertainties, but this is beyond the scope of the present study. Therefore, the primary aim is to validate the analytical model with MC computations. The experiments described a few sections below help to illustrate the actual performance of the analytical model and to identify the potential areas for future improvement. II.E.2. MC benchmark of the analytical model in phantoms. A phantom with a cylindrical geometry was used and is illustrated for several configurations in figure 5. It was made of 4923

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 5. Configuration of the phantoms used in this study for experiments and simulations. All values are given in centimetres. AP6, SB3, LN300 and G452 correspond to GAMMEX adipose, bone, lung and muscle tissue equivalent materials, respectively. (a) corresponds to the configuration used for homogeneous phantoms. (b) mimics a brain case (succession of adipose, bone, muscle and water-equivalent tissues). (c) mimics a lung case for which range mixing effects are expected because of the lateral heterogeneity muscle–lung. The (lateral) interface muscle–lung is aligned with beam-axis.

two co-axial cylinders. The outer cylinder was made of PMMA, the inner cylinder of slab inserts of various compositions and sizes. Details on the incident energies, the configurations chosen for the camera and the materials used are given in table 1. In all simulations, the slit of the camera was placed around the expected CSDA range of the protons. This CSDA range was computed according to all voxels densities and materials intersecting the proton incident trajectory. Thus, it may differ slightly from tabulated CSDA ranges because there are a few air voxels crossed by the proton before entering the phantom. The voxel resolution was set to 1 × 1 × 1 mm3. It is worth to mention here that tissue-equivalent plastics do not reproduce actual tissue composition. To determine whether the analytical model has similar performance in biological tissue, we performed a simulation with phantom configuration of figure 5(a), but with an ICRU cortical bone insert instead. Results are presented in appendix B. Figure 5(b) mimics a brain tumour with the quick succession of GAMMEX adipose and bone followed by water (‘brain’ configuration). The upstream PMMA slab serves as a range shifter. The small downstream PMMA slab corresponds to the wall of the water container. The two energy chosen (115 and 145 MeV) allow for studying the influence of a succession of heterogeneous materials either close or relatively far from the Bragg peak. Figure 5(c) attempts to reproducing a situation that can likely be met in lung tumours (‘lung’ configuration), i.e. the presence in the beam path of a lateral heterogeneity with a significant change of density (typically from 0.3 to 1 g cm−3). For an incident energy of 110 MeV, 4924

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Table 1. Description of the sets of tests performed in this study. d collimator is the distance

between the beam axis and the centre of the collimator. d coll − detectors is the distance between the centre of the collimator and the centre of the detectors. M is the magnification factor (d coll − detectors/d collimator). The first three sets correspond to a homogeneous inner cylinder. For sets marked witha, measured data is also presented. For phantom configurations marked withb more details are provided in figures 5(b) and (c). Set

Incident energies (MeV)

Inner cylinder materials

d collimator (cm)

d coll − detectors (cm)

1

150

Water

20

2

150

Water

15 25 35 25

3a

100/150/200

4a 5a

115/145 115

Adipose (AP6) Bone (SB3) Water Brain phantomb Lung phantomb

M

25

20 16 12 20

1.33 0.80 0.57 0.80 0.64 0.48 0.80

25 25

20 20

0.80 0.80

Expected CSDA ranges (cm) 15.9

15.9

8.6/17.3/28.5 4.9/9.8/16.0 7.9/15.9/26.1 9.2/14.2 11.4

the protons are assumed to stop around the middle of the muscle heterogeneity. However, the Bragg peak is not uniquely located, because of the so-called ‘range mixing’ effect with the apparition of two Bragg peaks in the dose distribution for the same incident energy. Such phenomenon is difficult to handle for pencil beam convolution algorithms (Paganetti 2012), mainly because of two reasons: (1) insufficient sampling of the incident spot and/or; (2) poor management of multiple scattering. The first issue can be prevented by the use of several rays to sample an incident spot with a known spatial distribution (here, up to 113 rays). For multiple scattering, the analytical model should suffer limitations comparable to pencil beam convolution algorithms. To quantify the accuracy of the analytical model for emission, the symmetric Γ-index was computed (Low et al 1998, Dhakal and Yepes 2014). The emission difference (ED, in % of maximum emission) and distance-to-agreement (DTA) criteria were 3%/0.5 mm, 3%/1 mm and 5%/2 mm. In principle, we are more interested by inaccuracies on the range than on the intensity. However, the emission distribution has a relatively low gradient before the end of the range. Therefore, a composite metric like the Γ-index is well adapted to evaluate the entire distribution. For the detection, the root-mean-square value of the differences between the analytical model and 4P was computed. However, what is the most important here is how accurate the analytical model is to predict a shift of Bragg peak position. Thus, the difference in the estimation of the Bragg peak position was also quantified by estimating the shift between both profiles. The shift was estimated as the one that minimizes the root-squared difference between the first profile and the second shifted profile. The minimization was implemented using an exhaustive search by steps of 0.2 mm. The PG gamma camera design presented here can detect shifts of Bragg peak position with an accuracy of about 1–2 mm (Perali et al 2014). Therefore, an agreement between 4P and the analytical model within 0.5 mm does not contribute significantly to the uncertainty in the estimation of the Bragg peak position. II.E.3. Experimental setup. To illustrate the potential of the analytical model, comparisons with measured profiles were performed for most of the configurations described in table 1 (sets 3, 4 and 5). 4925

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 6. Experimental setup. The target is made of two coaxial cylinders. The outer cylinder was always made of PMMA while the inner cylinder was made of homogeneous or heterogeneous tissue-equivalent plastic materials (see table 1 for more details). The collimator was pegged to the target support and mechanically aligned. The target was positioned on its support using the lasers.

Measurements were performed at the Proton Therapy Centre (PTC) in Prague, Czech Republic, in a gantry treatment room equipped with an IBA spot scanning dedicated nozzle. The spot was always incident along the axis of the cylindrical phantom (15 cm diameter, 20 or 40 cm axis depending on beam energy). Alignment of the beam with the target axis was ensured with radiochromic films. The PG camera collimator aperture was mechanically aligned with the target support and the target itself was then positioned on its support using the lasers. Distance was by default 25 cm between beam axis and the centre of the collimator and 20 cm between the centre of the collimator and the centre of the crystals. The beam time structure was that of the IBA C230 isochronous cyclotron, resulting in proton bunches of the order of 1 ns width every 10 ns. Charge delivery was monitored using the nozzle ionization chambers in clinical mode. When available, measured PG detection profiles were compared to computed profiles. Estimation of shifts of Bragg peak positions were performed using the same methodology as described above (minimization of the root-squared difference). In illustrations, measured values were corrected because these take into account the neutron contribution to the signal. A constant value was subtracted from measured values such that the average of measured values corresponded to the average of values computed by 4P. The experimental setup is illustrated in figure 6. II.E.4. Computation of PG detection profiles in CT geometries. In order to illustrate the capabilities of the analytical model in patient geometries, a CT of a patient with a lung tumour was selected (see figure 7). CT numbers were converted into mass density and materials using an arbitrary CT calibration curve. For the present simulation, 5 materials were used (lung, adipose, muscle and cortical ICRP tissues plus water). Voxel resolution was 2.1 × 2.1 × 2 mm3. The energies were selected such that the CSDA ranges corresponded approximately to the centre of the tumour (100 and 92 MeV for figures 7(a) and (b), respectively). The PG camera was set to its nominal configuration, i.e. M = 0.8, d collimator = 25 cm, d coll − detectors = 20 cm. 4926

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure 7. Beam incidences chosen for the computation of PG detection profiles in CT

images.

III. Results III.A. Validation of the analytical model in phantoms

Quantitative results for all sets of simulations listed in table 1 are regrouped in table 2. III.A.1 Homogeneous water with various camera configurations (sets 1 and 2). Figure 8 and table 2 display the results corresponding to the sets 1 and 2 listed in table 1. Excellent agreement between 4P and the analytical model is achieved for PG emission (figure 8(a)) with 99.4% of the points satisfying the 3%/0.5 mm Γ criteria (100% passing rate for 3%/1 mm). In the nominal configuration (figure 8(b)), excellent agreement is also observed for PG detection (estimations on the Bragg peak position within 0.2 mm). Agreement is still good for d collimator − d coll − detectors equalling 15–20 cm and 25–16 cm (figures 8(c) and (e), respectively). The estimation on the Bragg peak position is within 0.4 mm for both configurations. However, agreement is poorer for the configurations that lead to the largest field-of-views (figures 8(d) and (f)), with a deviation of 1.2 mm. III.A.2. Fixed camera configuration with various homogeneous materials (set 3). For the sake of conciseness, we illustrate only the results for bone (SB3) in figure 9. Figures for water and adipose (AP6) are shown in the appendix B. All results are summarized in table 2 for all materials of the simulation set 3. For PG emission, good agreement is obtained for almost all materials and energies excepted bone (SB3) for the 200 MeV beam with only 91.5% of the points passing the 3%/0.5 mm Γ criteria. For the 3%/1 mm Γ criteria, the passing rate is above 99% for all materials, excepted AP6 for the 200 MeV beam (97.6 % passing rate). For PG detection, good agreement is observed between 4P and the analytical model, with accuracy on the estimation of Bragg peak position within 0.4 mm. Compared to measurements, both 4P and the analytical model show similar accuracy on the estimation of the range, with a maximum deviation of −5.6 mm for AP6. III.A.3. Fixed camera configuration with heterogeneous phantom with various materials (sets 4 and 5). Results for PG detection and emission computed for the ‘brain’ configuration (set 4,

figure 5(b)) are summarized in table 2 and illustrated in figure 10. For both energies (115 and 4927

Energy (MeV)

150 150 150 2 150 150 3 100 150 200 100 150 200 100 150 200 4 115 145 5 110 No spot sampling

1

Set

Water Water Water Water Water AP6 AP6 AP6 SB3 SB3 SB3 Water Water Water ‘Brain’ ‘Brain’ ‘Lung’

Material 25–20 15–20 35–20 25–16 25–12 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20 25–20

d collimator− d coll − detectors (cm) 99.4 99.4 99.4 99.4 99.4 96.5 96.6 97.0 100.0 98.5 91.5 97.9 99.3 99.0 97.7 98.8 97.9 80.7

3/0.5

100.0 100.0 100.0 100.0 100.0 100.0 99.4 97.6 100.0 100.0 99.0 100.0 100.0 99.8 100.0 100.0 99.7 80.7

AM versus 4P Γ ≤ 1 (%)

3/1

100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 100.0 82.0

3/2

Γ-index ED(%)/DTA (mm)

Emission

4.9 5.4 7.5 7.0 11.0 6.9 5.3 9.1 4.8 6.9 5.3 7.0 5.0 9.3 7.9 5.4 6.2 19.0

AM versus 4P (%)

Rms

− 0.2 0.4 1.2 0.4 1.2 − 0.2 − 0.4 0.2 0.0 0.2 0.2 − 0.2 − 0.2 − 0.4 0.4 0.6 0.6 1.8

AM versus 4P (mm)

− 0.4 − 2.2 − 5.6 2.8 3.2 3.6 2.2 2.2 2.4 3.4 2.8 0.4 4.0

AM versus MS (mm)

0 − 1.6 − 5.6 2.6 3.2 4.2 2.2 2.6 2.8 3.0 2.6 − 1.8 − 1.8

4P versus MS (mm)

Accuracy of estimation of Bragg peak position

Detection

emission, the (symmetric) Γ-index was computed taking 4P as the reference and AM as the evaluated distribution, considering only 4P values above 10% of the maximum 4P value. The emission difference (ED) criterion is given in % of maximum emission. DTA stands for distance-toagreement. For detection, the root-mean-square value of the differences between AM and 4P was computed. The accuracy of the estimation of the range was estimated for AM taking the results computed with 4P as the reference. Finally, the accuracy of the estimation of the Bragg peak position of both AM and 4P were both estimated taking the measurements (MS) as the reference.

Table 2. Comparisons of prompt gamma emissions and detections computed with the analytical model (AM), and MC simulations (4P). For

Phys. Med. Biol. 60 (2015) 4915 E Sterpin et al

4928

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

7

-7

16

1 Water target 150 MeV incident protons

6

0.5 Γ value

No. of prompt gamma/(proton mm)

Prompt gamma emission (3-8 MeV)

4P Analytical Γ (3%/1 mm)

5 4 3 0

2


-4

x 10

1 0 0

30

60 90 120 Depth (mm) (a)

-6 Prompt

x 10

1.5

1

0.5 -40

-20 0 20 Detector axis (mm) (c) -6


2.5

x 10

2

8 6 4

1

0.5

-20 0 20 Detector axis (mm) (e)

40

40

Prompt gamma detection (3-6 MeV) 4P Analytical

12

dcollimator = 35 cm dcoll-detectors = 20 cm

10 8 6 4

-20 0 20 Detector axis (mm) (d) -6

3

4P Analytical

1.5

x 10

2 -40

Prompt gamma detection (3-6 MeV)

dcollimator = 25 cm = 16 cm d

-20 0 20 Detector axis (mm) (b) -7

40

coll-detectors

0 -40

10

14 No. of prompt gamma/(proton mm)


2


12

gamma detection (3-6 MeV) 4P Analytical

Prompt gamma detection (3-6 MeV) 4P Analytical

14

2 -40

150



2.5

x 10

x 10

2.5

40

Prompt gamma detection (3-6 MeV) 4P Analytical d = 25 cm collimator dcoll-detectors = 12 cm

2

1.5

1

0.5 -40

-20 0 20 Detector axis (mm) (f)

40

Figure 8. Comparisons of prompt gamma emission and detection profiles computed

with the MC model 4P and the analytical model. The emission was computed/simulated in a water phantom (a). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected CSDA range is 15.9 cm. Detection was computed for several configurations of the camera, with the couples of values distance from beam axis to centre of collimator (d collimator)—distance from centre of collimator to centre of detectors (d coll − detectors) equalling 25–20 (b), 15–20 (c), 35–20 (d), 25–16 (e) and 25–12 cm (f).

4929

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

-6

0

10

-3

20

30 40 Depth (mm)

50


2.5

Bone target 100 MeV incident protons

Γ value


1.1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0

x 10

1.5

1

0.5

-20 0 20 Detector axis (mm)

(a)


-6

3

0.8 0.7 0.6 0.5 0.4 0

0.3

x 10

0.2

40

Prompt gamma detection (3-6 MeV) 4P Analytical MS raw MS corrected

1

0.9


2



x 10

0 -40

60

Γ value





-3

x 10

2.5 2 1.5 1 0.5

0.1 20

-3

0 -40

100



-6

4


0.7 0.6 0.5 0.4 0.3

0

0.2

x 10

3.5 1

0.8

40


3 2.5 2 1.5 1 0.5

0.1 0 0


(b)

Γ value


x 10

40 60 80 Depth (mm)


0 0

40

80 120 Depth (mm)

0 -40

160


40

(c)


with MC simulations (4P) and the analytical model. The inner cylinder consists of homogeneous Gammex compact bone (SB3) and for the incident energies 100, 150 and 200 MeV ((a), (b) and (c) respectively). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges for SB3 are 4.9 (a), 9.8 (b) and 16.0 cm (c). For detection, raw measured data is also shown (MS). In MS corrected, a constant value is subtracted from MS raw such that the average of MS corrected values corresponds to the average of 4P values.

4930

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Prompt gamma emission (3−8 MeV)

Prompt gamma detection (3−6 MeV)

−3

−6

Brain configuration 115 MeV incident protons

2.5

A A SB3 P Water M i P M r 6 A

PMMA

0.8 0.7 0.6 0.5 0.4 0.3

x 10

4P Analytical MS raw MS corrected

1

Γ value



0

0.2


x 10

2

1.5

1

0.5

0.1 0 0

20

40 60 Depth (mm)

80

0 −40

100

Prompt gamma emission (3−8 MeV)

40

−6

PMMA

AASP Water i PBM r 6 3 M A

Γ value

6

3

Brain configuration 1 145 MeV incident protons

5 4 3

0

2




0 20 Detector axis (mm)

Prompt gamma detection (3−6 MeV)

−4

x 10

7

−20

(a)

x 10


2.5

2

1.5

1

0.5

1 0 0

30

60 90 Depth (mm)

120

0 −40

150

−20

0 20 Detector axis (mm)

40

(b)


with the MC model 4P and the analytical model. The inner cylinder consists of slabs of different materials (‘brain’ configuration, see figure 5(b)) and for the incident energies 115 and 145 MeV (a) and (b). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges are 9.2 (a) and 14.2 cm (b). For detection, raw measured data is also shown (MS). In MS corrected, a constant value is subtracted from MS raw such that the average of MS corrected corresponds to the average of 4P values.

145 MeV), there is good agreement between 4P and the analytical model for emission (97.7 and 98.8% of the points satisfying 3%/0.5 mm Γ criteria, respectively; 100% of the points satisfying 3%/1 mm for both energies). For detection, agreement on estimation of Bragg peak position is within 0.4 mm and 0.6 mm, respectively. With respect to measurements, similar agreement is obtained for 4P and the analytical model (from 2.6 to 3.4 mm). For the phantom in ‘lung’ configuration (set 5, figure 5(c)), agreement is much worse if super sampling of the spot is turned off, as shown in table 2 and figure 11. For emission, the passing rate for 3%/0.5 mm Γ criteria goes from 97.9 to 80.7% (figures 11(a) and (b), for super sampling of the spot turned on and off, respectively). For detection, similar observations can 4931

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Prompt gamma emission (3-8 MeV) -4

W a t e r

4

3

LN 300

G452

LN 300

P M1 M A

2 0 1

0 0

30

60 90 120 Depth (mm)



1.4

Lung configuration 110 MeV incident protons LN300 + G452

Γ value


x 10

Lung configuration 11 110 MeV incident protons 10 LN300 + G452 LN 300

LN 300

9

P 8 M M7 A 6

5 4 3 2 1 0

3 2 1 0 0

30

60 90 120 Depth (mm)

1 0.8 0.6 0.4 0.2

1.4


G452

1.2

(a)

Γ value


4

W a t e r


0 -40

150

Prompt gamma emission (3-8 MeV) -4 x 10 4P Analytical Γ (3%/1 mm)

Prompt gamma detection (3-6 MeV) -6 x 10

(b)

40

Prompt gamma detection (3-6 MeV) -6 x 10

1.2


1 0.8 0.6 0.4 0.2 0 -40

150



40


with the MC model 4P and the analytical model. The inner cylinder consists of slabs of different materials (‘lung’ configuration, see figure 5(c)), for the incident energy of 110 MeV, with and without super sampling of the incident spot ((a) and (b), respectively). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected range is 11.2 cm. For detection, raw measured data is also shown (MS). In MS corrected, a constant value is subtracted from MS raw such that the average of MS corrected corresponds to the average of 4P values.

be made, with the analytical model reproducing the estimation of the Bragg peak position of 4P within 0.6 and 1.8 mm for super sampling turned on and off, respectively. III.B. Clinical example

Figure 12 shows PG emission and detection for the clinical case illustrated in figure 7. For the lateral incidence (figure 12(a)), the symmetric Γ passing rate equal 88.8, 90.2 and 97.7 % for 3%/0.5, 3%/1 and 3%/2 mm, respectively. For detection, agreement on the estimation of the Bragg peak position is within 0.8 mm. For the posterior incidence (figure 12(b)): passing rates

4932

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

−

−

−

−

Γ

Γ

−

−

−

−

−

−

Γ

Γ

−

−


with the MC model 4P and the analytical model. The beam incidences, lateral right (a) and posterior (b) are illustrated along with a CT slice in figure 7. Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges are 25.0 (a) and 28.1 cm (b).

corresponding to the same criteria equal 82.6%, 88.5% and 98.0%, while agreement for the estimation of the Bragg peak position is within 0.6 mm. III.C. Simulation speed

Simulation time for 4P is typically around 15 h on a single CPU (Intel Core i7 2.3 GHz) considering 1 × 107 incident protons, which is equivalent to 5 × 109 protons for the variance parameters used. The analytical model is much faster: for a 200 MeV beam (worst case), around 0.3 s without super sampling and around 10 s with super sampling (113 rays). Computing emission is by far the longest, the detection process taking only 0.1 s.

4933

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

IV. Discussion IV.A. Validity of the analytical model with respect to MC simulations

The analytical model performs several approximations, among which: • The proton flux and the proton energy spectrum are assumed identical in the reference PG emission profile and in the considered geometry. Also, protons scattered out of the spot due to nuclear reactions are assumed to stay in the spot and so a wrong material could be assigned to these protons. • Multiple scattering is partially modelled • All the simplifications for PG detection: single convolution of PG emission with a response function, transport of the PG in the phantom/patient described by a single thickness parameter, etc. • No simulations of the neutrons (as PENH does). The effects of these hypotheses are now discussed in detail in the following sub-sections. IV.A.1. Constancy of proton flux and effect of protons scattered out of the spot by nuclear reactions. In a pure carbon material or in a pure oxygen material, proton flux and proton

energy spectrum will be different in depth because of the difference in secondary proton flux produced by nuclear reactions in carbon and oxygen. Thus the proton flux at a given depth in the individual simulations (pure target) may not be representative of the proton flux in mixed materials at the equivalent depth. The simulations in homogeneous materials (set 3) showed that no significant differences could be observed linked to this issue. Excellent agreement for all depths, all energies and all materials were observed. The only exception is for bone (SB3) for 200 MeV (figure 9(c)). However, the discrepancy for bone and 200 MeV is likely due to protons produced by nuclear reactions that escape the inner cylinder. These protons may also produce PGs. The analytical model assumes that they interact with bone, which has a much higher density (1.819 g cm−3), leading to a slight overestimation of the global PG emission. IV.A.2. Influence of multiple scattering. The effect of multiple scattering is modelled by curvilinear raytracing. Also, the ray-tracing algorithm takes into account the exact atomic composition and density laterally in the widening spot. Such strategy cannot model completely the effect of interfaces. For instance in a lateral heterogeneity involving two materials, a proton can travel in the first material, scatter in the second and scatter back in the first. This is not considered in the analytical model. However, it is a second order effect compared to the importance of sampling the spot size, even in the extreme cases shown in figures 11(a) and 12. Significant deviations (up to 40%) do occur locally in the emission profile but these are smeared out by the response function of the camera, leading to sub-millimetre accuracy of estimation of Bragg peak position (0.6, 0.8 and 0.6 mm for figures 11(a), 12(a) and (b), respectively). IV.A.3. Hypotheses for PG detection. In the nominal configuration of the camera, very good

agreement was observed overall between 4P and the analytical model (see figures 8–12 and table 2). However, it can be observed that for the initial millimetres in the detection profile, the analytical model overestimates the signal, especially if the field of view increases (figures 8(d) and (f)). This is expected from the profiles shown in figure 3(e): for a photon emitted 5 cm off camera centre, the response is reduced, violating the approximation of a single convolution. PG attenuation through the patient is modelled using a single thickness parameter (the integral radiological path length seen by the PGs, perpendicularly from beam-axis to the centre 4934

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

of the camera). This approximation should not cause significant deviations for most cases, because a relatively large 5 cm change of the pathlength only leads to a reduction of roughly 10% of the signal collected in the detector (figure 3(b)). For the considered clinical case (figures 7 and 12), this hypothesis seems robust enough (accuracy on the estimation of Bragg peak position within 0.8 and 0.6 mm). However, the validity of these approximations should be confirmed on multiple clinical cases. If significant deviations are observed, a solution would be to compute the radiological path length for every detector and correct the amplitude of the signal accordingly. IV.A.4. Neglecting neutron contribution. As observed in figures 9–11, the raw data do not

match 4P and the analytical model. For 150 MeV protons in Gammex SB3 (figure 9), the measured raw data is roughly twice superior to the results predicted by both 4P and the analytical model, in agreement with previous experiments performed by our group (Perali et al 2014). There is a roughly flat background contribution (see for instance Smeets et al (2012), figure 28 and Perali et al (2014) figure 16) to the total signal that is ignored in the simulations, hence the correction for the baseline (‘MS corrected’) and a lateral matching to predict the shift of the Bragg peak position (i.e. a 2D matching). Ultimately, a 1D matching (only laterally) would be preferable but would demand a very accurate prediction of the baseline as well. Such accurate prediction of the baseline requires an extensive simulation of the contribution of the neutrons (at least). This can be performed advantageously with MCNPX. Our latest prompt gamma experimental results in Perali et al (2014) reached a good agreement (of the order of 10%) with our MCNPX simulations in terms of uncorrelated neutron signal. This agreement contradicts the factor 2–3 difference observed in our earliest experimental results in Smeets et al (2012), but is unfortunately still not sufficient to allow 1D matching. This explains why the expected accuracy on the estimation of the position of the Bragg peak is similar when using 4P or MCNPX, as shown in the appendix A. Because we so far failed to predict with sufficient accuracy the background signal even when simulating neutrons, these can be safely ignored at the present time without jeopardizing the accuracy of the estimation of the Bragg peak position, as shown in the appendix A, figure A4. Such approximation simplifies considerably the design of a fast and accurate analytical model, which is the main aim of the present study. There is still left to understand the issues related to the accurate estimation of the background signal, which is mandatory to enable 1D matching for this particular camera design. In any case, the uncertainties related to nuclear cross section data (more than 5% according to ICRU 63) and to the experiments make unlikely a direct 1D matching of the measured and full MC simulated PG detection profiles. Thus, a 2D matching would be required anyway, which cancels the potential gain in accuracy one can expect by simulating neutrons and further justifies the implementation of an analytical model that neglects neutrons, leading to large gains in speed and simplicity. Obviously, fast models that would integrate neutrons and simulate background radiation accurately in order to enable 1D matching would be a meaningful progress compared to our current solution. However, our current solution was shown to be accurate enough for future clinical testing, which was the main objective. IV.B. Accuracy of Bragg peak position retrieval with respect to measurements

In all figures, the shape of the computed profiles visually matches well corrected measured profiles. However, deviations larger than 2 mm were observed for the estimation of the Bragg peak position (table 2). Many sources of uncertainties may affect this estimation. 4935

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

First, there are uncertainties on the mean excitation energies I used in the simulations. In such case, the error on the range should show a positive or negative trend with the incident energy. This is likely the case for set 3, material AP6, with an error for the analytical model of 0 and −5.6 mm, respectively for 100 MeV and 200 MeV, respectively. For other cases of set 3 and for set 4, there are errors for all energies between 2.2 and 3.6 mm for the analytical model and between 2.2 and 4.2 mm for 4P. Thus either the Bragg peak position is systematically underestimated in the model or the camera is systematically positioned slightly downstream of the expected range. Such deviations may have several causes: inaccurate computing of the actual mean energy of the protons from their nominal energy; incomplete calibration of the camera; uncertainties in the positioning of the phantom and the camera… Measurements were performed in PTC (Prague), while the beam model used in 4P and the analytical model is the one of Gantry 4, WPE (Essen). Although geometrical details of the incident spots have negligible influence for the considered experiments, the actual mean energy of the incident protons may be slightly inaccurate, leading to a systematic range shift. The positioning of the camera is unlikely the cause of the discrepancy. The collimator was pegged to the target support and the camera secured to the collimator. However, the target was positioned on its support using the lasers, which may cause a discrepancy within 1 mm. However, the approximations made in the analytical model are validated with comparisons with MC simulations. Thus, any remaining difference with experimental data will affect both MC and analytical models. The measurements are shown here to illustrate that the shape of expected PG detection profiles reproduce reasonably measured ones, which is the case as shown in figures 9–12. The observed differences can be taken care of by additional calibration of the camera (Perali et al 2014) and accurate characterization of the considered proton unit in both MC and the analytical model. IV.C. Comparison with previous work

To our knowledge, there are no previous attempts for comprehensive analytical simulation of PG emission and detection. However, various approaches are described in literature for proton dose monitoring with PET. For instance, Oelfke et al (1996) and Lopatiuk-Tirpak et al (2011) compute the activity by scaling an analytical expression of the proton flux with the isotope production cross-section and a correction for decay. In their approach, the secondary proton flux due to elastic and inelastic nuclear interactions is not taken into account, while it is implicitly included in the precomputed emission profiles of our model. Parodi and Bortfeld (2006) use analytical fits to MC simulations of isotope production and to dose distributions provided by a clinical dose engine. A single filter per reference material, valid for all incident energies, is derived and used to convolve dose distributions to evaluate the activity. Further results show that the method is able to accurately predict the activity near the distal edge for clinical cases (Attanasi et al 2011). It would be interesting to determine whether such methodology can be transposed to PG imaging. However, computed dose distributions are needed to apply the filter, which is impractical for our specific case where we need to compute PG emission spot by spot. Thus, the dose distributions for potentially each spot need to be provided by the treatment planning system, feature that is usually not available and would be cumbersome due to the number of dose maps needed. IV.D. Practical perspectives of the analytical model

The main purpose of the analytical model is to simulate all the spots in a given treatment plan in order to compare the measurement to the expected profiles and estimate the range shifts. 4936

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Alternatively, the simulation could be used to select some probe spots that would be good candidates to verify the range prior to treatment delivery, in terms of localization, sensitivity to specific range uncertainties and expected profile quality. In either case, simulating more than thousands spots is unrealistic for MC algorithms considering clinical time constraints, hence the implementation of an analytical solution. However, MC simulations can still be used in the probe spot approach. For instance, the analytical model could be used for a first pre-selection of the probe spots and then the profiles recomputed with MC in order to ensure the best accuracy possible for the expected detection profiles. MC codes faster than 4P could be used, like MCsquare (Souris et al 2014) that includes PG computing capability. It is worth mentioning that the speed of the analytical simulation is decent in the present version (of order of 10 s for 113 rays on a single CPU) but can be further improved using faster codes than MATLAB. Simulating all the spots at this speed is still a lengthy process and probably some compromise might be needed for the number of rays simulated (the simulation time being roughly proportional to the number of rays). Although the simulation of the PG detection phase is tailored to the considered prompt gamma geometry, the model for PG emission is universal for all types of treatment units (considering spot scanning). The only thing needed is the precomputed PG emission profiles for the 5 elements. These are valid for any machine provided that the incident energy spectrum is available. Any set of PG cross sections can be used to compute these profiles. Of course, the chosen PG acquisition technique must be robust to the absence of neutrons in the simulated emission. It is important to mention that the 1–2 mm accuracy of the camera only refers to the ability of the camera to estimate the Bragg peak position when the position of the camera (the expected range) matches the actual range. However, a measured shift of the Bragg peak position by the camera may not correspond to an actual shift of the range in the patient. For instance, if the patient is positioned 5 mm upstream the beam, everything else being equal, the camera will detect a 5 mm shift, even though the proton is reaching the correct depth in the patient. This example illustrates that a PG camera solution must be linked to other imaging devices to distinguish between true or false measured shifts (or matches), like cone-beam CT, flat panel detectors or an optical surface-imaging system (Schöffel et al 2007). V. Conclusions We have successfully modelled analytically PG emission and detection for a PG camera prototype. The analytical model is based on the hypothesis that PG emission can be computed as the dose in a pencil beam convolution algorithm. The detection of PGs is modelled using a convolution with a transfer function. Overall, excellent agreement was achieved with MC simulations and promising agreement with measurements was observed. The analytical model can potentially be used to compute the expected PG detection for all spots of a given treatment plan at practical speed. Appendix A: benchmark of 4P with MCNPX for prompt gamma simulation and detection In this appendix, we describe the benchmark of 4P with MCNPX for PG emission and detection. 4P was used instead of Geant4 or MCNPX for several reasons: 1. In Smeets et al paper (2012), MCNPX uses ENDF/B-VI (McLane 2001) data tabulated up to 150 MeV for prompt gamma emission, which is the same data used in ICRU 63 (although extended up to 250 MeV in ICRU 63). 4937

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure A1. Reference setup for prompt gamma detection. The length of the phantom is extended to 40 cm for an incident energy of 230 MeV. The collimator is long enough to cover the entire phantom +1 cm at each edge. More details on the geometry are available in Smeets et al (2012).

2. At the time of the design of this study, the beam optical model of Grevillot et al (2012) was already implemented in 4P. 3. We do not have access to MCNPX code, only executables. 4. 4P uses PENCT that is specific to track particles in CT structures, with trivial preparation of the CT input file. This feature is well mastered by our group from previous work (for instance, Sterpin et al 2009). In MCNPX, the translation of CT geometry to MCNPX input is possible but not trivial (Jabbari et al 2012). 5. Because of our full access to the entire 4P code, we could implement efficient variance reduction techniques as described in the main document.

A.1. Material and methods MCNPX (version 2.5.0) was used as a benchmark. The simulation settings were exactly the same as in Smeets et al (2012), that is, a monoenergetic proton pencil beam (with no lateral spatial distribution) impinging a PMMA phantom with the camera placed right at the position of the Bragg peak. It is important to mention that the geometry for the benchmark of 4P is the one of the first prototype of the slit camera detailed in Smeets et al (2012). In the main article, all experiments and modelling are related to the current prototype (Perali et al 2014). Three sets of tests were performed to benchmark 4P. First, a PMMA cylinder as shown in figure A1 was modelled with an incident proton spot of 150 MeV. The energy spectrum of all photons exiting the PMMA target computed with both MCNPX and 4P were compared. 150 MeV was chosen because it is the maximum energy for which MCNPX uses tabulated data. Second, the reference setup for the prompt gamma camera as described in Smeets et al (2012) was modelled in 4P. The PG detection profiles were simulated with MCNPX and 4P taking into account only the photon component. The profiles were computed and compared for 100, 160 and 230 MeV incident proton energies. Finally, the impact of neglecting neutrons was assessed by comparing the accuracy of the Bragg peak position retrieving method for MCNPX simulated PG profiles with neutron transport turned on, MCNPX PG profiles with 4938

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure A2. Energy spectrum of photons exiting the PMMA phantom geometry for

150 MeV incident protons computed by MCNPX (black line) and 4P (red solid circles). The inset shows the 511 keV lines.

neutron transport turned off and 4P PG profiles. For measurements, data was available for the first version of the prototype for 100 and 160 MeV only. With 230 MeV, data was not usable because of the low current necessary to avoid saturating the camera of the old prototype, which degrades the stability of the beam. A.2. Results and discussion A.2.1. Energy spectrum

Figure A2 shows a comparison between MCNPX and 4P of energy spectra of photons exiting the PMMA phantom for a 150 MeV incident proton beam. It is worth mentioning that the simulation of neutrons is turned off in MCNPX. Very good agreement is overall achieved. For 12 C there are rays in the prompt gamma emission spectrum at 2.00 and 4.44 MeV. For 16O, there are rays at 2.74, 6.13, 6.92 and 7.12 MeV. During a nuclear event, 16O may breakdown in 14N and 12C that can emit also prompt gammas. The rays for 14N are at 1.64, 2.31 and 5.11 MeV. At low proton energies, the most important is the 4.44 MeV line for 12C. For 16O, it is the 6.13 MeV line (Verburg et al 2012). Because 12C is the dominant element in PMMA and the 4.44 MeV is an important spectrum line, it appears clearly in the energy spectrum. The 6.13 and 7.12 MeV lines of 16O can also be distinguished. Other lines cannot be distinguished from the continuous background. The line at 511 keV comes from electron-positron annihilation in the PMMA target. 4939

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure A3. Comparison of prompt gamma detection profiles simulated by 4P and

MCNPX for energies of 100, 160 and 230 MeV ((a), (b) and (c), respectively). In MCNPX, the simulation of neutrons was turned off.

A.2.2. PG detection profiles for 100, 160 and 230 MeV

Figure A3 shows comparisons of PG detection profiles according to the prototype design showed in figure A1. Profiles are within 6%, 10% and 18% on average for 100 MeV, 160 MeV and 230 MeV respectively. Overall, the main difference is a global shift of 4P values towards higher values for higher energies. This could be explained by the different nuclear models used for both codes. Indeed, the same tabulated cross sections (including gamma-production data) are used by MCNPX up to 150 MeV whenever available (la150h library), but the default physics model (Bertini intranuclear cascade model) with default settings is used otherwise. However, the shapes are very similar, which is the most important to detect shifts of Bragg peak positions. Indeed, profiles predicted by 4P lead to similar accuracy for the estimation of the Bragg peak shift than MCNPX (figure A4). The results indicate also that simulating neutrons do not improve the accuracy for the detection of Bragg peak shifts. Neglecting neutrons leads to the necessity of a 2D matching of the PG detection profile (both laterally and vertically (amplitude)) to determine the Bragg peak shift. This is less accurate than a 1D matching (only lateral matching), because of the increased sensitivity to statistical fluctuations from bin to bin. However, 1D matching can be performed only if the amplitude is accurately predicted by both the simulation model and accurately measured. For the simulation, on one hand, there are uncertainties on the cross-sections for nuclear reactions (above 5% for the total cross-section (ICRU 63)). Also, all the potential sources of neutrons 4940

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Figure A4. Comparison of the accuracy on the estimation on the Bragg peak position (here called range) for three incident energies (100, 160 and 230 MeV) using either MCNPX as the reference (pink triangles) or 4P (green triangles). As a reference, the accuracy using MCNPX profiles only is shown as well (blue dots). Simulation of neutrons was turned on in MCNPX. The measurements were performed for the first camera prototype described in Smeets et al (2012). For each curve in the figure, the camera is positioned at a reference, expected beam penetration depth and profiles are acquired with range shifts (with respect to the expected beam penetration depth) and compared to a reference profile with no range shift in order to estimate the value of the shift. In the curve with blue dots, the profiles acquired with range shifts are MC simulations with MCNPX and the reference profile is a MC simulation with MCNPX. In the curve with red triangles, the profiles acquired with range shifts are measurements with the former camera prototype based on the HiCam system and the reference profile is a MC simulation with MCNPX. In the curve with green triangles, the profiles acquired with range shifts are measurements with the former camera prototype based on the HiCam system and the reference profile is a MC simulation with 4P.

that could contaminate the camera signal should be modelled. On the other hand, measurements still need to be corrected for dead-time, which is not the case for both prototypes (see Perali et al (2014) for more details). In other words, it is unlikely that a correspondence in amplitude can be achieved between simulated and measured signals, which forces the use of a 2D matching. The neutrons contributing as a flat background, neglecting the neutrons should not have a significant impact on a 2D matching procedure in this case. Appendix B: results for ICRP cortical bone, water and Gammex adipose (AP6) B.1. ICRP cortical bone The objective is to verify the accuracy of the analytical model for human tissues instead of tissue equivalent plastic, for the same settings as set 3 in table 1 of the main document. For PG emission, similar agreement was achieved than for SB3, (Γ passing rate of 100%, 100% and 99.8% for criteria 3%/1 mm). For PG detection, good agreement was also achieved with an accuracy for the estimation of the Bragg peak position equalling 0,−0.4 and −0.6 mm (100 MeV, 150 MeV and 200 MeV, respectively).

4941

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Prompt gamma emission (3–8 MeV)

Prompt gamma detection (3–6 MeV)

-4

14



1

7 6 Γ value


x 10

5 4 3

0

2 1 0 0

20

40 60 Depth (mm)

x 10

-7

4P Analytical

12 10 8 6 4 2 -40

80


40

(a)

Prompt gamma emission (3– 8 MeV)


-4



6 5

Γ value


x 10

4 3

0

2 1 0 0

30

60 90 120 Depth (mm)

x 10

-7

4P Analytical

12 10 8 6 4 2 -40

150


(b)

40

Prompt gamma emission (3– 8 MeV) -4


13


1

5 Γ value


x 10

4 3 2

0

1 0 0

60

120 180 Depth (mm)

x 10

12

4P Analytical

11 10 9 8 7 6 5 4 -40

240

-7


40

(c)

Figure B1. Comparisons of prompt gamma emission and detection profiles computed

by the MC model 4P and the analytical model. The emission was computed/simulated with the inner cylinder made of homogeneous ICRP compact bone and for the incident energies 100, 150 and 200 MeV ((a), (b) and (c) respectively). For emission, symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges are 4.9 (a), 9.8 (b) and 16.0 cm (c). 4942

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915



-3

-6

2.5 No. of prompt gamma/(proton mm)

1


0.7 0.6

Γ value


x 10

0.5 0.4 0.3

0

0.2 0.1 0 0

20

40 60 Depth (mm)

x 10


2 1.5 1 0.5 0 -40

80


40

(a) Prompt gamma emission (3–8 MeV)


-4

-6

6 5 4 3 0

2 1 40

80 120 Depth (mm)


7

0 0

3 x 10

1


Γ value


x 10

2.5 2 1.5 1 0.5 0 -40

160



40

(b) Prompt gamma emission (3–8 MeV)

Prompt gamma detection (3– 6 MeV)

-4

-6

1

6 5 4 3 2

0

1 0 0

60

120 180 Depth (mm)


3.5


Γ value


x 10

x 10

3 2.5 2 1.5 1 0.5 0 -40

240



40

(c)

Figure B2. Comparisons of prompt gamma emission and detection profiles computed by the MC model 4P and the analytical model. The emission was computed/simulated with the inner cylinder made of homogeneous water and for the incident energies 100, 150 and 200 MeV ((a), (b) and (c) respectively). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges are 7.9, 15.9 and 26.1 cm. For detection, raw measured data is also shown (MS). In MS corrected, values are corrected such that their average corresponds to the average of 4P values. 4943

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915



-4

-6

1.4



1

5 4 Γ value


x 10

3 2

0

1 0 0

20

40 60 Depth (mm)

80

x 10


1.2 1 0.8 0.6 0.4 0.2 0 -40

100


40

(a)



-4

-6

2 x 10 No. of prompt gamma/(proton mm)


1 3

Γ value


x 10

2

0

1 0 0

40

80 120 Depth (mm)


1.5

1

0.5

0 -40

160


40

(b)



-4

-6


2


3 1 Γ value


x 10

2

0

1

0 0

60

120 180 Depth (mm)

240

x 10

2.5 2 1.5 1 0.5 0 -40

300 (c)



40

Figure B3. Comparisons of prompt gamma emission and detection profiles computed

by the MC model 4P and the analytical model. The emission was computed/simulated with the inner cylinder made of homogeneous Gammex adipose tissue (AP6) and for the incident energies 100, 150 and 200 MeV ((a), (b) and (c) respectively). Symmetric Γ-indexes (analytical versus 4P) are shown one point over two for the visibility of the graph. The expected ranges are 8.6, 17.3 and 28.5 cm. For detection, raw measured data is also shown (MS). In MS corrected, values are corrected such that their average corresponds to the average of 4P values.

4944

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

B.2. Water Results for measured and simulated profiles are illustrated in figure B2. Quantitative results may be found in the main document in table 2. B.3. Gammex adipose (AP6) Similarly, figure B3 shows the results for Gammex adipose (AP6). Quantitative results can be found in table 2 of the main document. References Attanasi F, Knopf A, Parodi K, Paganetti H, Bortfeld T, Rosso V and Del Guerra A 2011 Extension and validation of an analytical model for in vivo PET verification of proton therapy—a phantom and clinical study Phys. Med. Biol. 56 5079–98 Bom V, Joulaeizadeh L and Beekman F 2012 Real-time prompt gamma monitoring in spot-scanning proton therapy using imaging through a knife-edge shaped slit Phys. Med. Biol. 57 297–308 Dhakal T R and Yepes P 2014 A symmetric probabilistic γ-index for Monte Carlo dose comparisons Phys. Med. Biol. 59 N153–61 Golnik C G et al 2014 Range assessment in particle therapy based on prompt γ-ray timing measurements Phys. Med. Biol. 59 5399–422 Grevillot L, Bertrand D, Dessy F, Freud N and Sarrut D 2012 GATE as a GEANT4-based Monte Carlo platform for the evaluation of proton pencil beam scanning treatment plans Phys. Med. Biol. 57 4223–44 ICRU 1989 Tissue Substitutes in Radiation Dosimetry and Measurements (ICRU Report 44) (Bethesda, MD: International Commission on Radiation Units and Measurements) ICRU 2000 Nuclear Data for Neutron and Proton Radiotherapy and for Radiation Protection (ICRUReport 63) (Bethesda, MD: International Commission on Radiation Units and Measurements) Jabbari I, Shahriari M, Aghamiri S M R and Monadi S 2012 Advantages of mesh tallying in MCNPX for 3D dose calculations in radiotherapy J. Radioanal. Nucl. Chem. 291 831–7 Knopf A-C and Lomax A 2013 In vivo proton range verification: a review Phys. Med. Biol. 58 R131–60 Lopatiuk-Tirpak O, Su Z, Li Z, Hsi W, Meeks S and Zeidan O 2011 Spatial correlation of proton irradiation-induced activity and dose in polymer gel phantoms for PET/CT delivery verification studies Med. Phys. 38 6483–8 Low D A, Harms W B, Mutic S and Purdy J A 1998 A technique for the quantitative evaluation of dose distributions Med. Phys. 25 656–61 Lu H M 2008 A point dose method for in vivo range verification in proton therapy Phys. Med. Biol. 53 N415–22 McLane V 2001 ENDF-102 Data Formats and Procedures for the Evaluated Nuclear Data File ENDF6 Technical Report BNL-NCS-44945-01/04-Rev (Upton, NY: Brookhaven National Laboratory, National Nuclear Data Center) Min C H, Lee H R, Kim C H and Lee S B 2012 Development of array-type prompt gamma measurement system for in vivo range verification in proton therapy Med. Phys. 39 2100–7 Oelfke U, Lam G K Y and Atkins M S 1996 Proton dose monitoring with PET: quantitative studies in Lucite Proton dose monitoring with PET: quantitative studies in Lucite Phys. Med. Biol. 41 177–96 Paganetti H 2012 Range uncertainties in proton therapy and the role of Monte Carlo simulations Phys. Med. Biol. 57 R99–117 Parodi K and Bortfeld T 2006 A filtering approach based on Gaussian-powerlaw convolutions for local PET verification of proton radiotherapy Phys. Med. Biol. 51 1991–2009 Perali I et al 2014 Prompt gamma imaging of proton pencil beams at clinical dose rate Phys. Med. Biol. 59 5849–71 Perl J et al 2012 TOPAS: an innovative proton Monte Carlo platform for research and clinical applications Med. Phys. 39 6818–37

4945

E Sterpin et al

Phys. Med. Biol. 60 (2015) 4915

Peterson S W, Robertson D and Polf J 2010 Optimizing a three-stage Compton camera for measuring prompt gamma rays emitted during proton radiotherapy Phys. Med. Biol. 55 6841–56 Pinto M, Dauvergne D, Freud N, Krimmer J, Letang J M, Ray C, Roellinghoff F and Testa E 2014 Design optimisation of a TOF-based collimated camera prototype for online hadrontherapy monitoring Phys. Med. Biol. 59 7653–74 Salvat F 2013 A generic algorithm for Monte Carlo simulation of proton transport Nucl. Instru. Methods Phys. Res. B 316 144–59 Salvat F, Fernandez-Varea J M and Sempau J 2011 PENELOPE-2011: a Code System for Monte Carlo Simulation of Electron and Photon Transport (Issy-les-Moulineaux: OECD Nuclear Energy Agency) Schöffel P J, Harms W, Sroka-Perez G, Schlegel W and Karger C P 2007 Accuracy of a commercial optical 3D surface imaging system for realignment of patients for radiotherapy of the thorax Phys. Med. Biol. 52 3949–63 Smeets J et al 2012 Prompt gamma imaging with a slit camera for real-time range control in proton therapy Phys. Med. Biol. 57 3371–405 Souris K, Lee J A and Sterpin E 2014 Intel Xeon Phi implementation of a fast multi-purpose Monte Carlo simulation for proton therapy Med. Phys. 41 535 Sterpin E, Salvat F, Olivera G and Vynckier S 2009 Monte Carlo evaluation of the convolution/ superposition algorithm of Hi-ArtTM tomotherapy in heterogeneous phantoms and clinical cases Med. Phys. 36 1566–75 Sterpin E, Sorriaux J and Vynckier S 2013 Extension of PENELOPE to protons: simulation of nuclear reactions and benchmark with Geant4 Med. Phys. 40 111705 Sterpin E, Sorriaux J, Souris K, Vynckier S and Bouchard H 2014 A Fano cavity test for Monte Carlo proton transport algorithms Med. Phys. 41 011706 Testa M et al 2010 Real-time monitoring of the Bragg-peak position in ion therapy by means of single photon detection Radiat. Environ. Biophys. 49 337–43 Verburg J M, Riley K, Bortfeld T and Seco J 2013 Energy- and time-resolved detection of prompt gamma-rays for proton range verification Phys. Med. Biol. 58 L37–49 Verburg J M and Seco J 2014 Proton range verification through prompt gamma-ray spectroscopy Phys. Med. Biol. 59 7089–106 Verburg J M, Shih H A and Seco J 2012 Simulation of prompt gamma-ray emission during proton radiotherapy Phys. Med. Biol. 57 5459–72

4946

Copyright of Physics in Medicine & Biology is the property of IOP Publishing and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.

Proton range verification through prompt gamma-ray spectroscopy.

Range verification of proton radiotherapy with prompt gamma rays.

Simulation and experimental verification of prompt gamma-ray emissions during proton irradiation.

Range verification of passively scattered proton beams based on prompt gamma time patterns.

Factors influencing the accuracy of beam range estimation in proton therapy using prompt gamma emission.

Compton Camera and Prompt Gamma Ray Timing: Two Methods for In Vivo Range Assessment in Proton Therapy.

Theoretical detection threshold of the proton-acoustic range verification technique.

Detecting prompt gamma emission during proton therapy: the effects of detector size and distance from the patient.

Verifying performance characteristics of quantitative analytical systems: calibration verification, linearity, and analytical measurement range.

Technical Note: Experimental carbon ion range verification in inhomogeneous phantoms using prompt gammas.

Feasibility of Using Distal Endpoints for In-room PET Range Verification of Proton Therapy.

Active neutron and gamma-ray imaging of highly enriched uranium for treaty verification.

Prompt-Gamma Activation Analysis.

Prompt gamma ray evaluation for chlorine analysis in blended cement concrete.

Prompt gamma imaging of proton pencil beams at clinical dose rate.

Experimental Comparison of Knife-Edge and Multi-Parallel Slit Collimators for Prompt Gamma Imaging of Proton Pencil Beams.

Neutrino and cosmic-ray emission from multiple internal shocks in gamma-ray bursts.

Range assessment in particle therapy based on prompt γ-ray timing measurements.

Prompt gamma-ray neutron activation analysis of boron-10 in biological materials.

Monte Carlo study on the sensitivity of prompt gamma imaging to proton range variations due to interfractional changes in prostate cancer patients.

Development of a new deuterium-deuterium (D-D) neutron generator for prompt gamma-ray neutron activation analysis.

Absolute prompt-gamma yield measurements for ion beam therapy monitoring.

Multielement proton-induced X-ray emission analysis of betel leaves, betel nuts, and lime.

Trace-element analysis of human teeth and bone by proton-induced X-ray emission.