1. INTRODUCTION
For more than half a century, the main tool for hydrodynamic experiment was pulsed X-ray radiography. It allows one to see the inside of an object with a complex structure without disturbing it. However, the dose limitations, position resolution, and backgrounds still limit the utility of X-ray radiography for obtaining the pictures with higher resolution. About 50 years ago, proton radiography technology, which is predicted to achieve higher spatial resolution, was put forward. However, the experimental results of earlier research show that proton radiography has a lower spatial resolution comparing with X-ray radiography (Koehler, Reference Koehler1968; Cookson et al., Reference Cookson, Armitage and Ferguson1972; West & Sherwood, Reference West and Sherwood1972). The main reason is that the interaction between the protons and the sample atoms, especially multiple Coulomb scattering (MCS) lead to the image blur. The problem has not been solved until Mottershead and Zumbro (Reference Mottershead and Zumbro1998) proposed a typical magnetic quadrupole lens system for the correction of image blur in 1998. The proton radiography system has been successfully developed at some laboratories, such as LANSCE (Ziock et al., Reference Ziock, Adams, Alrick, Amann, Boissevain, Crow, Cushing, Eddleman, Espinoza, Fife, Gallegos, Gomez, Gorman, Gray, Hogan, Holmes, Jaramillo, King, Knudson, London, Lopez, McClelland, Merrill, Morley, Morris, Mottershead, Mueller, Neri, Numkena, Pazuchanics, Pillai, Prael, Riedel, Sarracino, Stacy, Takala, Thiessen, Tucker, Walstrom, Yates, Zumbro, Ables, Aufderheide, Barnes, Bionta, Fujino, Hartouni, Park, Soltz, Wright, Balzer, Flores, Thompson, Prigl, Scaduto, Schwaner, Saunders and O'Donnell1998; King et al., Reference King, Ables, Adams, Alrick, Amann, Balzar, Barnes, Crow, Cushing, Eddleman, Fife, Flores, Fujino, Gallegos, Gray, Hartouni, Hogan, Holmes, Jaramillo, Knudsson, London, Lopes, McDonald, McClelland, Merril, Morley, Morris, Naivar, Parker, Park, Pazuchanics, Pillai, Riedel, Sarracino, Shelley, Stacy, Takala, Thompson, Tucker, Yates and Ziock1999; Morris et al., Reference Morris, Ables, Alrick, Aufderheide, Barnes, Buescher, Cagliostro, Clark, Clark, Espinoza, Ferm, Gallegos, Gardner, Gomez, Greene, Hanson, Hartouni, Hogan, King, Kwiatkowski, Liljestrand, Mariam, Merrill, Morgan, Morley, Mottershead, Murray, Pazuchanics, Pearson, Sarracino, Saunders, Scaduto, Schach von Wittenau, Soltz, Sterbenz, Thompson, Vixie, Wilke, Wright and Zumbro2011) in US, ITEP, and IHEP in Russia (Antipov et al., Reference Antipov, Afonin, Vasilevskii, Gusev, Demyanchuk, Zyat'kov, Ignashin, Karshev, Larionov, Maksimov, Matyushin, Minchenko, Minkeev, Mirgorodskii, Peleshko, Rud'ko, Terekhov, Tyurin, Fedotov, Trutnev, Burtsev, Volkov, Ivanin, Kartanov, Kuropatkin, Mikhailov, Mikhailyukov, Oreshkov, Rudnev, Spirov, Syrunin, Tatsenko, Tkachenko and Khramov2010; Kolesnikov et al., Reference Kolesnikov, Golubev, Demidov, Dudin, Kantsyrev, Mintsev, Smirnov, Turtikov, Utkin, Sharkov and Fortov2010) and GSI in Germany (Tahir et al., Reference Tahir, Shutov, Varentsov, Hoffmann, Spiller, Lomonosov, Wieser, Jacoby and Fortov2002). In addition, some new proton radiography systems are being developed at FAIR in Germany (Merrill et al., Reference Merrill, Golubev, Mariam, Turtikov and Varentsov2009) and Institute of Modern Physics (IMP) in China (Zhao, Reference Zhao2011). The Institute of Modern Physics of Chinese Academy of Science (IMP-CAS) has proposed two proton radiography beam line (Sheng et al., Reference Sheng, Zhao, Yang, Wei, Jiang, Zhou, Cheng, Yan, Li, Yang, Yuan, Xia and Xiao2014) to investigate proton radiography technique. This paper is mainly to develop the proton radiography terminal of IMP, including calculation of the beam optic parameters and the simulation of proton radiography by the Monte Carlo method.
2. THE ZUMBRO LENS SYSTEM
The attenuation of a proton beam through a given density of material is described by the equation (Ziock et al., Reference Ziock, Adams, Alrick, Amann, Boissevain, Crow, Cushing, Eddleman, Espinoza, Fife, Gallegos, Gomez, Gorman, Gray, Hogan, Holmes, Jaramillo, King, Knudson, London, Lopez, McClelland, Merrill, Morley, Morris, Mottershead, Mueller, Neri, Numkena, Pazuchanics, Pillai, Prael, Riedel, Sarracino, Stacy, Takala, Thiessen, Tucker, Walstrom, Yates, Zumbro, Ables, Aufderheide, Barnes, Bionta, Fujino, Hartouni, Park, Soltz, Wright, Balzer, Flores, Thompson, Prigl, Scaduto, Schwaner, Saunders and O'Donnell1998):
$$N/N_0 = e^{ - l/\lambda} $$where N 0 and N are respectively the incident number and the transmitted number of proton, and l is the path length of the object, and λ is the mean free path of proton. The transmitted protons undergo MCS with the atom of target, so it will produce an angular distribution after the target. This angular distribution can be approximately described by the Gauss distribution, the half-width of this Gauss distribution is expressed as (Wohl, Reference Wohl1984)
$${\rm \phi} _{\rm c} = \displaystyle{{0.0136\,{\rm GeV}} \over {{\rm \beta} cp}}\sqrt {\displaystyle{z \over {X_0}}} \left[ {1 + 0.038{\rm In}\left( {\displaystyle{z \over {X_0}}} \right)} \right] $$where p is the beam momentum, βc is the velocity of the proton, z is the thickness of the object, X 0 is the proton's radiation length for the material (Tsai, Reference Tsai1974), and its value is related to the material's mass number and atomic number.
The angle distribution produced by MCS will lead to the image blur seriously. So to cancel the effect by MCS, a typical lens was designed by Mottershead and Zumbro, which is called Zumbro lens. The Zumbro lens system (Mottershead & Zumbro, Reference Mottershead and Zumbro1998) is the type of angle-matching lens which is widely used to focus the particles to achieve point-to-point imaging for removing most image blur by MCS. Furthermore, the magnetic lens system can provide angle sorting which allows one to insert an angle cut aperture to aid material identification.
Zumbro lens system consists of two identical quadrupole doublet cells; each cell includes a positive quadrupole and a negative quadrupole, as shown in Figure 1. This quadrupole beam line is reflection symmetric structure because its second cell's parameters are identical with the first one. The reflection symmetric lens focuses the particles at image plane with the magnification of 1 along the x- and y-axes. To achieve above aim, the lens transfer matrix should be 
$ - \overrightarrow I$. If the transfer matrix of one cell is M, then the whole lens system is given by the 2 × 2 matrix identity
$$\overrightarrow R = \overrightarrow M ^2 \, = \, - \overrightarrow I \det (\overrightarrow M ) + {\rm \tau} \overrightarrow M $$
where 
${\rm \tau} \equiv {\rm Tr}(\overrightarrow M ) = {M_{11}} + {M_{22}}$ is the trace of 
$\overrightarrow M$, let τ = 0, the determinant 
$\det \overrightarrow M = 1$ for beam line matrices. So the lens transfer matrix 
$\overrightarrow R = - \overrightarrow I $, which achieves an inverse image with the same size. Considering the particle's momentum deviation Δ = δp/p, the final beam spot size of this particle at the image plane is (to first order in Δ)
$$x_{\rm f} = R_{11} x + R_{12} {\rm \theta} + (R_{11} ^{\prime} + {\rm \omega} R_{12} ^{\prime})x\Delta + R_{12} ^{\prime}{\rm \varphi} \Delta $$where (x, θ) is the particle's initial coordinate, ω = L −1 is the correlation coefficient, L is the distance between the particle source and the entry plane of the lens. A beam would lie along the line θ = ωx in phase space, here if we make
$$R_{11} ^{\prime} + {\rm \omega} R_{12} ^{\prime} = 0, \; {\rm or} \; {\rm \omega} = - R_{11} ^{\prime}/R_{12} ^{\prime} $$all position-dependent chromatic aberrations vanish, so the angle-sorting function can be carried out. With R 12 = 0, R 11 = −1, the final position is
$$x_{\rm f} = - x + R_{12} ^{\prime}{\rm \varphi} {\Delta} $$
Fig. 1. The schematic presentation of the Zumbro lens.
The remaining chromatic aberration depends only on the deviation angle φ.
3. DESIGN OF THE PROTON RADIOGRAPHY BEAM LINE
Before design the proton radiography beam line, the appropriate energy of proton should be decided first. If only considering the penetrating power of proton, the 1 GeV proton is enough for the target. But the spatial resolution must also be considered. The main factor is the MCS effect, which will reduce the spatial resolution, the value of position error can be estimated by Xu (Reference Xu2006)
$$\Delta r = {\rm \theta} _0 l/\sqrt 3 $$where θ0 is the root-mean square (RMS) value of the scattering angle by MCS, l is the thickness of the target. In Eq. (7), Δr∝θ0, θ0 can be derived by
$${\rm \theta} _0 \approx \displaystyle{{14.1} \over {\,p\rm \beta}} \sqrt {\sum\limits_i^n {\displaystyle{{{\rm \rho} _i l_i} \over {R_i}}}} $$where p is the momentum of proton, β = v/c is the relative velocity of proton, ρi and l i are the density and thickness of the target, respectively, R i is the radiation length of proton in the target material. For a determined target, ρi, li, and R i are constants. So θ0∝1/pβ. According to Eqs. (7) and (8), we can derive Δr∝1/pβ, That is to say the choice of higher energy protons can effectively reduce Δr value. The maximum design energy of proton is 2.6 GeV in Heavy Ion Research Facility in Lanzhou-Cooling Storage Ring (HIRFL-CSR), so to reduce Δr we choose 2.6 GeV in this paper to simulate the proton radiography.
According to the beam condition of HIRFL-CSR at IMP, a 2.6 GeV proton radiography terminal is designed, which consists of a matching section, an imaging lens system, and proton detector. Figure 2 shows the schematic of the 2.6 GeV proton radiography beam line scenario. Q1, Q2, and Q3 are the matching lens system, and Q4, Q5, Q6, and Q7 are the imaging lens system.
Fig. 2. The schematic presentation of 2.6 GeV proton radiography beam line scenario (unit/m).
The matching lens section consists of three quadrupole lens to provide the required phase space correlation upstream of the object. Moreover, this section also expands the beam's transverse size to illuminate the field-of-view of the object fully. The Zumbro magnets structure is used here in the imaging lens system, which will eliminate a major part of the chromatic and geometric aberrations, and match the rays with different MCS angles. Moreover, a collimator at the Fourier plane is used to eliminate the larger angle scattered particle to enhance the spatial resolution and the capability of material identification.
According to adjust the parameters from Eqs. (4) to (6), makes the imaging lens system to meet the Zumbro magnets requirements, and the magnifications are both 1 along the x- and y-axes. To meet the above requirements, the corresponding parameters for the 2.6 GeV proton radiography beam line are calculated. The gradients and geometrical parameters of the imaging lens are calculated by WINAGILE code (Bryant, Reference Bryant2000), Figure 3 shows the beam envelop of the beam line by WINAGILE. The values of the parameters are listed in Table 1. With the parameters in Table 1, the particle transport trajectory in the imaging lens is simulated by My-BOC code (Zhang et al., Reference Zhang, Yang and Lv2010); it can be seen that the proton trajectory in image lens achieve the point-to-point imaging in Figure 4.
Fig. 3. Beam envelope of the beam line.
Fig. 4. The particle tracking of proton radiography.
Table 1. Parameters of the imaging lens system of 2.6 GeV proton radiography beam line
4. MONTE CARLO SIMULATION
4.1. Monte Carlo Model
To test the designed proton radiography system (see Fig. 2), the beam parameter and the image were respectively simulated using Monte Carlo codes Geant4 (Agostinelli et al., Reference Agostinelli, Allison, Amako, Apostolakis, Aroujo, Arce, Asai, Axen, Banerjee, Barrand, Behner, Bellagamba, Boudreau, Broglia, Brunengo, Burkhardt, Chauvie, Chuma, Chytracek, Cooperman, Cosmo, Degtyarenko, Dell'Acqua, Depaola, Dietrich, Enami, Feliciello, Ferguson, Fesefeldt, Folger, Foppiano, Forti, Garelli, Giani, Giannitrapani, Gibin, Cadenas, Gonzalez, Abril, Greeniaus, Greiner, Grichine, Grossheim, Guatelli, Gumplinger, Hamatsu, Hashimoto, Hasui, Heikkinen, Howard, Ivanchenko, Johnson, Jones, Kallenbach, Kanaya, Kawabata, Kawabata, Kawaguti, Kelner, Kent, Kimura, Kodama, Kokoulin, Kossov, Kurashige, Lamanna, Lampen, Lara, Lefebure, Lei, Liendl, Lockman, Longo, Magni, Maire, Medernach, Minamimoto, Mora de Freitas, Morita, Murakami, Nagamatu, Nartallo, Nieminen, Nishimura, Ohtsubo, Okamura, O'Neale, Oohata, Paech, Perl, Pfeiffer, Pia, Ranjard, Rybin, Sadilov, Di Salvo, Santin, Sasaki, Savvas, Sawada, Scherer, Sei, Sirotenko, Smith, Starkov, Stoecker, Sulkimo, Takahata, Tanaka, Tcherniaev, Safai Tehrani, Tropeano, Truscott, Uno, Urban, Urban, Verderi, Walkden, Wander, Weber, Wellisch, Wenaus, Williams, Wright, Yamada, Yoshida and Zschiesche2003) and G4beamline (Roberts et al., Reference Roberts, Beard, Ahmed, Huang and Kaplan2011). The parameters of the proton radiography system were set up according the schematic in Figure 2 and the data in Table 1. During the simulation process, the class of QGSP_BIC in Geant4 is invoked to achieve all essential interactions, including Coulomb interaction between the protons and the electrons, strong interaction between the protons and the atomic nuclei. Some self-field effects were ignored, such as the space charge effect and fringe field effect. The ideal detectors are respectively located at the Fourier plane and the image plane to obtain the beam parameters and image results.
4.2. Simulation Results of the Beam Parameters
Based on the above model, the protons angle distribution is firstly simulation using Geant4 code. Here the incident number of proton is 1 × 106. The simulated results are compared with the computed data by Eqs. (1) and (2), they agree with each other. So self-field can be ignored, and the simulation code can be trusted.
The beam parameters are also investigated using G4beamline code. The results are shown in Figure 5, the phase space diagrams at the object plane (Fig. 5a) and the image plane (Fig. 5c) are identical with each other in the horizontal and vertical directions. It means that when proton pass through the imaging lens, the magnification of proton at the image plane is 1. And the horizontal and vertical phase spaces are superposition at the Fourier plane in Figure 5b. The parameters of the beam are totally the same as the calculated results, and meet the design requirement.
Fig. 5. Horizontal and vertical phase space diagrams (the blue and red parts are the phase space in horizontal and vertical directions, respectively): (a) at the object plane, (b) at the Fourier plane, and (c) at the image plane.
To get the high-resolution image, high transportation efficiency of the beam line for 2.6 GeV protons must be ensured from the particle source to image plane without target. So we input the optic parameters into G4beamline to simulate the beam line without target to check the transmission efficiency. In the simulation, the number of incident particles is 106. An ideal detector is located at the image plane to collect protons, and the results show all incident protons are collected by the ideal detector, it means that the transmission efficiency is 100%. This result by Monte Carlo code agrees with the result using WINAGILE and My-BOC codes.
When the protons pass through the target, a lot of secondary particles, such as gamma, neutron, proton, electron, π meson, and k meson, will be produced. These secondary particles may influence the image quality. So the secondary particles need to be investigated. The secondary particles transportation in the proton beam line is simulated using Geant4 code. The results show that only 11 particles arrive at the image detector with a 100 mm diameter under 105 incident protons. A reasonable explanation is that almost all of the charged particles are unable to reach the image plane for these particles, which cannot pass through the focus lens, and a majority of neutral particles are also lost in the beam pipe. Thus, the effect of the secondary particles can be neglected.
4.3. Simulation Results of the Sample Image
To test the image effect of the designed proton radiography beam line, the images of two samples are simulated using Geant4 code. The first sample for testing the spatial resolution is a circular aluminum object (2 cm thickness) with two 10 × 3 mm2 rectangular orifices in the center. Another sample for testing the thick resolution is a step brass plate with three layers to test the thick resolution. Figure 6 shows the structures and geometry sizes of two samples.
Fig. 6. (a) Drawing of orifices sample (unit/mm). (b) Drawing of steps sample (unit/mm).
Figure 7a shows the two-dimensional image of the first sample on the x–y plane, and the two rectangular orifices can be seen clearly. Figure 7b is the proton counts distribution along the x-axis; the two steps in the distribution curve are corresponding to the position of the rectangular orifices. Figure 7c is obtained by a differential operation for the data in Figure 7b. Four peaks marked as 1, 2, 3, and 4 can be clearly seen in Figure 7c. The Sigma of four peaks, which is the standard deviation of the Gauss fitting function, is defined as the position resolution corresponding to four edgy of two rectangular orifices. The data analysis shows that the best resolution is about 36 μm for the edge of rectangular orifices.
Fig. 7. The simulation results with orifices sample by Geant4 code. (a) Image on the x–y plane; (b) transmission in the edge of the x-direction; (c) Gauss fit for the derivative of edge transmission.
According to the same method, the simulation results for the steps sample are obtained, and the results are shown in Figure 8. Figure 8a shows the two-dimensional image of the steps sample on the x–y plane. Every step edge can be seen clearly in the image. The same analysis process is also made for the data of Figure 8b. In Figure 8c, four peaks marked by 1, 2, 3, and 4 are the positions of the step edge. The spatial resolution between the hole at the center and brass is 83 μm. And the best spatial resolution between different thicknesses of brass steps is about 161 μm.
Fig. 8. The simulation results with step sample by Geant4 code. (a) Image on the x–y plane; (b) transmission in the edge of the x-direction; (c) Gauss fit for the derivative of edge transmission.
Moreover, further simulations are operated to estimate the thickness of brass steps by the imaging data in Figure 8b. Brass plates with the thickness of 2, 4, 8, and 12 mm are respectively put into the object plane of the beam line as the target to get the transmission ratio of proton. The simulation results are showed in the second column in Table 2. According to the simulation results, a fitting formula between the transmission ratio of proton and the thickness is established:
$$f\,(x) = p_1 x^2 - p_2 x + p_3 $$where p 1 = 17.39, p 2 = −35.38, and p 3 = 19.49 are fitting parameters, x is the transmission ratio of proton, and f(x) is the areal density of the brass target, which is corresponding to the target thickness. The transmission ratio of steps sample are obtained by calculating the average number of proton each steps in Figure 8b, and are showed in the third column in Table 2. The calculated thicknesses of the steps sample are obtained by a combination of the third data in Table 2 and Eq. (9). The maximum differences are about 7.3% between the real thickness and the calculated thickness.
Table 2. Calculation results of step thickness
5. CONCLUSION
A 2.6 GeV proton beam line for proton radiography, which consists of the matching lens and the Zumbro lens systems, has been designed. The optic parameters of the matching magnetic lens and the Zumbro lens systems are presented by simulating the proton transportation using WINAGILE and My-BOC codes. The transmission efficiency and secondary particles are respectively investigated by Monte Carlo simulation using Geant4 and G4beamline codes. The results show that the transmission efficiency of proton is 100%, the effect of the secondary particles can be neglected. The simulation test of the proton images for the orifices sample and the steps sample have been carried out using Geant4 code. The structure of two samples can be identified clearly in the images and the best spatial resolution is 36 μm. The difference of the thickness is not >10%. Compared with the previous work (Borghesi et al., Reference Borghesi, Audebert, Bulanov, Cowan, Fuchs, Gauthier, Mackinnon, Patel, Pretzler, Romagnani, Schiavi, Toncian and Willi2005, Reference Borghesi, Sarri, Cecchetti, Kourakis, Hoarty, Stevenson, James, Brown, Hobbs, Lockyear, Morton, Willi, Jung and Dieckmann2010; Zhao et al., Reference Zhao, Hu, Cheng, Wang, Peng, Golubev, Zhang, Lu, Zhang, Zhou, Wang, Xu, Ren, Li, Lei, Sun, Zhao, Wang, Wang and Xiao2012), the proton radiography beam line is the first design scheme of high-energy proton radiography which includes the Zumbro lens system in China. According to simulation results, the resolution will be improved comparing with the heavy-ion radiography experiment we have made at IMP (Sheng et al., Reference Sheng, Zhao, Yang, Wei, Jiang, Zhou, Cheng, Yan, Li, Yang, Yuan, Xia and Xiao2014).
Considering the actual experiment, the detect efficiency will influence the final resolution, so the simulation for detector is the next step of work. Furthermore, the collimator in Fourier plane will be studied, including its aperture and material, which will improve the spatial resolution of the terminal. Moreover, further simulation will be made with higher energy protons in High Intensity heavy-ion Accelerator Facility (HIAF) proposed by IMP; higher spatial resolution can be expected due to higher energy of proton. In addition, the experiment will be carried out at IMP in future.
ACKNOWLEDGMENTS
This work is supported by the NSAF Joint Funds of National Natural Science Foundation of China under Grant No. 11176001. The authors would like to thank the accelerator staffs of IMP.
