Radiation Protection Dosimetry (2015), Vol. 163, No. 2, pp. 217 –221 Advance Access publication 3 May 2014
doi:10.1093/rpd/ncu145
A METHOD BASED ON MONTE CARLO SIMULATION FOR THE DETERMINATION OF THE G(E) FUNCTION Wei Chen1,2,*, Tiancheng Feng2, Jun Liu2, Chuanying Su2 and Yanjie Tian2 1 Institute of Nuclear Science and Technology. SiChuan University, Chengdu 610051, China 2 Northwest Institute of Nuclear Technology, Xi’an 710024, China *Corresponding author: [email protected] Received 4 December 2013; revised 17 March 2014; accepted 1 April 2014
INTRODUCTION (1)
Dose rate is very important for radiation physics , and the assessment of the environmental gamma-ray dose in some areas, such as the area around nuclear power plants, is a matter of great concern(2). NaI(Tl) detector is known to be highly sensitive to gamma radiation and is widely used for the measurement of gamma spectrum in environmental radiation measurement. Apparently, if NaI(Tl) measurement systems have dose measurement function, the gamma-ray dose could be obtained synchronously as measuring the gamma spectrum, that is useful for mapping and assessing the dose rate in the interested area. There are two representative methods for evaluating dose from measured spectrum(3). One is a group categorised in the unfolding analysis using response functions of the detector(4 – 7); the other is the G(E) function method, which can directly obtain the dose value by applying spectrum-dose conversion operator to a pulse-height spectrum(8 – 10). The unfolding techniques mentioned earlier are effectively used in the fields of environmental radiation monitoring, but accurate and large number of response functions are needed(3). If only the total dose rate is interested, the G(E) function method, which could be used to convert gamma spectra into dose rate(8 – 12, 14), is recommended in ICRU Report 53(11). The G(E) function method was proposed first by Moriuchi et al.(8). Hiromi et al. used the G(E) function method in conversion of gamma spectrum of in situ Ge(Li) detector(12) and obtained better results than the results using HASL method(13). The precision of the G(E) function is very dependent on the number of input
response spectra, and it is desirable to obtain as many response gamma-ray spectra as possible at appropriately homogeneous energy intervals(14). However, in practice, there are some restrictions in choosing suitable gamma-ray source owing to the limited number of monoenergetic sources with proper energies; only a few response functions can be obtained. Kimiaki et al. used Monte Carlo method to obtain simulated response gamma-ray spectra(14), and break through those restrictions. In a study by Kimiaki (14), the simulated spectrums were used according to the method presented in the literatures(10, 12) to determine the coefficients of the G(E) function. In this paper, based on a 400 ` 400 ` 1600 NaI(Tl) detector and Monte Carlo N-particle Transport Code (MCNP), a new method to determine the G(E) function is proposed, using MCNP to simulate spectrums of various monoenergetic gamma rays and deposited energy, which was then converted into absorbed dose rate in air, in the corresponding window of full-energy peak in an air ball. The value of the G(E) function at energy E in spectra was obtained by dividing the absorbed dose rate in air by counts of corresponding full-energy peak. Curve-fitting software 1st0pt was used to fit the G(E) values vs. E using the G(E) function; the coefficients were determined. In the authors’ method, only data of full-energy peak are needed. Compared with using all the data of all standard spectrums in conventional method, the number of the histories and consumed time in simulation processing are decreased. Meanwhile, complexity of numerical calculation is reduced, compared with that using least square method to process while data of all standard spectrums.
# The Author 2014. Published by Oxford University Press. All rights reserved. For Permissions, please email: [email protected]
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
The G(E) function method is a spectrometric method for the exposure dose estimation; this paper describes a method based on Monte Carlo method to determine the G(E) function of a 400 `3 400 `3 1600 NaI(Tl) detector. Simulated spectrums of various monoenergetic gamma rays in the region of 40 – 3200 keV and the corresponding deposited energy in an air ball in the energy region of full-energy peak were obtained using Monte Carlo N-particle Transport Code. Absorbed dose rate in air was obtained according to the deposited energy and divided by counts of corresponding full-energy peak to get the G(E) function value at energy E in spectra. Curve-fitting software 1st0pt was used to determine coefficients of the G(E) function. Experimental results show that the calculated dose rates using the G(E) function determined by the authors’ method are accordant well with those values obtained by ionisation chamber, with a maximum deviation of 6.31 %.
W. CHEN ET AL.
PRINCIPLE AND METHODS In the algorithm of the G(E) function method, there exists a conversion operator between the count rates and the absorbed dose rate in air at channel i of the gamma spectra(8), which is a constant and could be described as in the following equation: Di ni
GðEi Þ ¼
ð1Þ
end X
Di ¼
i¼b
end X
GðEi Þni
ð2Þ
i¼b
where D is the absorbed dose rate in air converted from the spectra (nGy h21); b is the virtual staring channel of a spectra, the energy at which is the first positive value; end is the last channel of the spectra and other notations mean the same as in equation (1). From equation (2), it could be found that the key to obtain the absorbed dose rate in air is to determine the conversion operator at channel i, i.e. G(Ei). The relationship function between G(Ei) and the energy at channel i is defined as the G(E) function, which is given by a polynomial equation(10, 12, 14). GðEÞ ¼
K max X
Ak ðlog EÞk1
2
@S ¼ @Ak
J P j¼1
2
0
ððDj /Dj Þ 1Þ
¼0
@Ak
ð5Þ
0
where Dj is the calculated dose rate corresponding to the jth standard spectrum derived from the activity of the jth standard source(12) (nGy h21) and J is the number of the standard spectrums. In the literature(12), the number of terms in equation (3), i.e. Kmax, is increased one by one; when the value of S 22 converges to a minimum, the procedure is stopped. In the literature(14), simulated spectrums were obtained by MCNP, and the value of Kmax was studied, but the method for determining Ak is still as equation (5). In this paper, a simpler method is proposed. MCNP was used to simulate gamma spectra with monoenergetic incident gamma rays, but only the data of full-energy peaks are used. The region of the full-energy peak could be determined according to the simulated spectra, and then the deposited energy Einair in an air ball, whose radius is 9.8 cm (because the radius of ionisation chamber in this study is 9.8 cm), in the same region of the same incident gamma ray is obtained through simulation, converting Einair into absorbed dose rate in air as:
ð3Þ DðEÞ ¼ Einair f
k¼1
This equation is the basic form of the G(E) function, in which G(E) is the absorbed dose rate in air per unit count rates at the channel corresponding to energy E in a spectra (nGy h21 cps21); Kmax is the maximum term number of the G(E) function; Ak is coefficient for the kth term of the polynomial equation. From equations (2) and (3), it could be known that the key of using the G(E) function method to obtain absorbed dose rate in air is to determine Kmax and Ak. In the conventional method(12, 14) for determining the G(E) function, some radioactive sources or MCNP were used to obtain standard spectrums or simulated standard spectrums, and whole data of each spectrum are used, as in the following equation: " # K end max X X k1 Ak ðn ji log ðEi ÞÞ ð4Þ Dj ¼ k¼1
@
ð6Þ
where f is a conversion factor converting megaelectron volt into nanogray per hour, f ¼ 106` 1.6021` 1027 ` 3600. Then, the G(E) value at the ith channel corresponding to energy Ei is derived as in the following equation: GðEi Þ ¼
DðEi Þ Pi
ð7Þ
where Pi is the counts of the full-energy peak of incident gamma ray with energy Ei. Equation (8) is derived by substituting equation (7) into (2):
i¼b
D¼
end X i¼b
218
GðEi Þni ¼
end X DðEi Þn
i
i¼b
Pi
ð8Þ
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
where E(i) is the energy at channel i (MeV); ni is the count rates at channel i, cps (counts per second); G(Ei) is the conversion operator (nGy h21 cps21); and Di is the corresponding absorbed dose rate at channel i (nGy h21). Using these operators at different channel of the spectrum, the absorbed dose rate in air could be derived, as in the following equation: D¼
where Dj is the dose rate corresponding to the jth standard spectrum (nGy h21); nji is the count rates in channel i of the jth standard spectrum (cps) and other notations mean the same as in the previous equations. Then, the least square method was used, as in the following equation:
DETERMINATION OF THE G(E) FUNCTION Table 1. fwhm of 400 `3 400 `3 1600 NaI(Tl) detector at various energy. Energy/MeV fwhm/MeV
0.0595 0.0111
0.6617 0.0476
1.1732 0.0663
1.3325 0.0722
According to least square algorithm, equation (10) is used to determine a, b and c. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi2 @ X4 fwhmðE Þ a b Ei þ cEi 2 i i¼1 @xm
Figure 1. Model of measurement.
¼0
ð10Þ
ðm ¼ 1; 2; 3; x1 ¼ a; x2 ¼ b; x3 ¼ cÞ
Figure 2. Model of 400 ` 400 ` 1600 NaI(Tl) detector.
Simulation of standard spectrums and deposited energy in the air ball The model used in this study is shown in Figure 1, and the model of 4L NaI(Tl) detector is shown in Figure 2. In order to determine the window of a full-energy peak, the Gaussian energy broadening (GEB) treatment (15) of MCNP was used, which could be used to better simulate a physical radiation detector, in which energy peaks exhibit GEB. The parameters of GEB in MCNP specify the full width at half maximum (fwhm) of the observed energy broadening in a physical radiation detector(15): pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi fwhmðEÞ ¼ a þ b E þ cE 2
ð9Þ
where E is the energy of incident gamma ray, and the units of a, b and c are MeV, MeV1/2 and none(15), respectively. A radioactive point source containing 241 Am, 137Cs and 60Co was placed at a random location (here, the distance from point source to detector has not been measured, because fwhm represents the energy resolution of a detector, and the relative position between detector and source does not affect it) in front of the 400 ` 400 ` 1600 surface of the NaI(Tl) detector to get the fwhm at 59.5, 661.7, 1173.2 and 1332.5 keV, which are listed in Table 1.
Determination of the coefficients of the G(E) function Substitute the deposited energy in the air ball of different incident gamma rays into equation (6) to calculate the absorbed dose rate in air, then substitute the counts of the full-energy peak and the corresponding calculated absorbed dose rate in air into equation (7) to get the value of G(Ei) at the ith channel corresponding to energy Ei, the values obtained were listed in Table 2. Curve-fitting software 1sp0pt (16) was used to fit the relationship between G(E) value and energy E. There are two parameters needed to be determined by users, the fitting function and the fitting algorithm. Here, the basic form of the G(E) function, as in equation (3), was used as the fitting function; the Universal Global Optimisation (UGO)(16) was chosen as fitting algorithm. The Kmax is determined by increasing k one by one, and the procedure stopped when the correlation coefficient R is most close to 1. The determined coefficients of the G(E) function are listed in Table 3. The fitted curve is shown in Figure 3. The accuracy of the fitted curve, i.e. the accuracy of the G(E) function determined by the proposed method, was validated indirectly by experiment.
219
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
Substituting fwhm in Table 1 into equation (10), a ¼ 20.00557, b ¼ 0.06675, c ¼ 0.00028 were obtained. In this paper, an ionisation chamber was used in validation step. The radius of the ionisation chamber is 9.8 cm, so an air ball with 9.8 cm radius was simulated, and the centre of this air ball is located in the position of the centre of 400 ` 400 ` 1600 NaI(Tl) detector in Figure 1. The deposited energy in the air ball is recorded by f6 card(15) of MCNP. In simulating process, the histories are 109 and the maximum estimated relative error of interested simulated data in full-energy peak is ,5 %.
W. CHEN ET AL. Table 2. The G(E) value at energy E of the spectra. E/KeV 40 50 60 80 100
G(E)/nGy h21 cps21
E/KeV
G(E)/nGy h21 cps21
E/KeV
G(E)/nGy h21 cps21
1.89E-08 8.82E-08 1.90E-07 4.70E-07 7.22E-07
200 300 500 800 1200
2.21E-06 4.00E-06 8.61E-06 1.70E-05 2.99E-05
1600 2000 2600 2900 3200
4.37E-05 5.86E-05 8.34E-05 9.58E-05 1.07E-04
Table 3. Value of coefficients of the G(E) function. Coefficient
Coefficient
Coefficient
Coefficient
Coefficient
21.95E-02 6.87E-02 21.04E-01
A4 A5 A6
8.97E-02 24.75E-02 1.59E-02
A7 A8 A9
23.29E-03 3.84E-04 21.95E-05
Table 4. Absorbed dose rate in air contributed by natural radionuclides by the G(E) function and ionisation chamber.
Figure 3. Fitted curve of the G(E) function of 400 ` 400 ` 1600 NaI(Tl) (correlation coefficient R ¼ 0.999).
Validation According to the layout shown in Figure 1, a 400 ` 400 00 ` 16 NaI(Tl) detector was used to measure the spectrum of radioisotopes in environment on concrete square, asphalt road and grass, respectively. The calculated absorbed dose rate in air by the G(E) function and measured results by the ionisation chamber are listed in Table 4. It is shown that the calculated dose rate in air by the G(E) function and the measured results by ionisation chamber coincide well, with a maximum deviation of 25.4 %. But it could be noticed that the absorbed dose rates in Table 4 are from environmental background radioactivity and ,160 nGy h21. It is necessary to validate that the determined G(E) function is still applicable; if there are artificial isotopes, the absorbed dose rate is higher and the uniformity of gamma-ray field is not the same as in the environmental background radioactivity. For this purpose,
Location
Value by G(E)/nGy h21
Value by ionisation chamber/nGy h21
Relative error/%
Concrete square Asphalt road Grass
156.5
157.4+2.6
20.57
152.8
145.0+4.1
5.38
150.7
144.2+4.9
4.51
a hexagon mixed area source (241Am þ 137Cs þ 60Co þ 152Eu) with 50 cm length of each side was used to get spectrums (the layout is shown in Figure 4). The results calculated by the determined G(E) function were compared with the results by ionisation chamber, which is listed in Table 5. The results in Table 5 showed that the G(E) function determined by the proposed method is still applicable, and the uniformity of gamma-ray field did not obviously impact on the accuracy of calculation of absorbed dose rate in air by the determined G(E) function. CONCLUSIONS The method for the determination of the coefficients of the G(E) function is based on MCNP simulation. In the simulated spectrums, only the energy region and count rates of full-energy peak were focussed on, without analysis of whole data of spectrums. Compared with the conventional method, it reduced the strong dependence on radioactive sources and the complexity of numerical calculation. The results obtained in the validation experiment show that the
220
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
A1 A2 A3
Value
DETERMINATION OF THE G(E) FUNCTION
Table 5. Absorbed dose rate in air by the G(E) function and ionisation chamber (using source). d/m
Without source 0 1 3 5
Value by G(E)/nGy h21
Value by ionisation chamber/nGy h21
Relative error/%
153.6
158.4+2.8
23.03
973.3 373.9 209.2 181.6
980.4+10.8 351.8+4.8 198.4+5.1 176.1+3.9
20.72 6.31 5.42 3.14
relative error of absorbed dose rate in air by the determined G(E) function is ,6.4 %, compared with the results by ionisation chamber. The uniformity of gamma-ray field did not obviously impact on the accuracy of calculation of absorbed dose rate in air by the determined G(E) function. REFERENCES 1. Clouvas, A., Xanthos, S., Antonopoulos-Domis, M. and Silva, J. Monte Carlo calculation of dose rate conversion factors for external exposure to photon emitters in soil. Health phys. soc. 78, 295–302 (2000). 2. Chul-Young, Yi et al. Calculation of spectrum to dose conversion factor for a HPGe spectrometer. J. Korean Phys. Soc. 30, 186 –193 (1997).
221
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
Figure 4. Layout of the experiment using hexagon area source.
3. Moriuchi, S., Tsutsumi, M. and Saito, K. Development of a dosimetric system using spectrometric technique suitable for operational radiation dose measurements and evaluation. In: Proceeding of The 10th International Congress of IRPA. Hiroshima, May 2000. P-3b-197, 1–8 (2000). 4. Minato, S. and Kawano, M. Evaluation of exposure due to terrestrial gamma-radiation by response matrix method. J. Nucl. Sci. Technol. 7, 401–406 (1970). 5. Okano, M., Izumo, K., Kumagai, H., Katou, T., Nishida, M., Hamada, T. and Kodama, M. Measurement of environmental radiations with a scintillation spectrometer equipped with a spherical NaI(Tl) scintillator. Natural Radiation Environment III, CONF-780422, 2, 896–911 (1980). 6. Urabe, I., Tsujimoto, T. and Katsurayama, K. Estimation of photon energy spectra based on 34 ` 34 response matrix. J. Radiat. Res. 19, 163–170 (1978). 7. Moriuchi, S., Nagaoka, T. and Saito, K. Measurements, theoretical analysis and dose evaluation of low-level external radiation in environment. Report of the Japan-USSR Seminar on Radiation Effects Research. Tokyo, June 1990, National Institute of Radiological Sciences, (1991). 8. Moriuchi, S. and Miyanage, I. A spectrometric method for measurement of low-level gamma exposure dose. Health Phys. 12, 541–551 (1966). 9. Moriuchi, S. and Miyanaga, I. A method of pulse height weighting using the discrimination bias modulation. Health Phys. 12, 1481– 1486 (1966). 10. Moriuchi, S. A New Method of Dose Evaluation by Spectrum-Dose Conversion Operator and Determination of the Operator. JAERI 1209, Japan Atomic Energy Research Institute (1970). 11. International Commission on Radiation Units and Measurements. ICRU REPORT 53. Gamma-ray Spectrometry in the Environment. ICRU Publication 53 (1994). 12. Hiromi, T. et al. Environmental gamma-ray exposure rates measured by in-situ Ge(Li) spectrometer. J. Nucl. Sci. Technol. 17, 281– 290 (1980). 13. Harold, L. B. et al. In situ Ge(Li) and NaI(Tl) gammaray spectrometry. HASL-258, Health and Safety lab., (1972). 14. Kimiaki, S. et al. Calculation of the spectrum-dose conversion operator for a HPGe detector by means of a Monte Carlo simulation. J. Korean Phys. Soc. 25, 504–509 (1992). 15. Los Alamos National Laboratory. MCNP4C Monte Carlo N-Particle Transport Code System User’s Manual. Los Alamos National Laboratory, (2000). 16. 7D-Soft High Technology, Inc. 1st0pt manual. 7D-Soft High Technology, INC., China (2009).
doi:10.1093/rpd/ncu145
A METHOD BASED ON MONTE CARLO SIMULATION FOR THE DETERMINATION OF THE G(E) FUNCTION Wei Chen1,2,*, Tiancheng Feng2, Jun Liu2, Chuanying Su2 and Yanjie Tian2 1 Institute of Nuclear Science and Technology. SiChuan University, Chengdu 610051, China 2 Northwest Institute of Nuclear Technology, Xi’an 710024, China *Corresponding author: [email protected] Received 4 December 2013; revised 17 March 2014; accepted 1 April 2014
INTRODUCTION (1)
Dose rate is very important for radiation physics , and the assessment of the environmental gamma-ray dose in some areas, such as the area around nuclear power plants, is a matter of great concern(2). NaI(Tl) detector is known to be highly sensitive to gamma radiation and is widely used for the measurement of gamma spectrum in environmental radiation measurement. Apparently, if NaI(Tl) measurement systems have dose measurement function, the gamma-ray dose could be obtained synchronously as measuring the gamma spectrum, that is useful for mapping and assessing the dose rate in the interested area. There are two representative methods for evaluating dose from measured spectrum(3). One is a group categorised in the unfolding analysis using response functions of the detector(4 – 7); the other is the G(E) function method, which can directly obtain the dose value by applying spectrum-dose conversion operator to a pulse-height spectrum(8 – 10). The unfolding techniques mentioned earlier are effectively used in the fields of environmental radiation monitoring, but accurate and large number of response functions are needed(3). If only the total dose rate is interested, the G(E) function method, which could be used to convert gamma spectra into dose rate(8 – 12, 14), is recommended in ICRU Report 53(11). The G(E) function method was proposed first by Moriuchi et al.(8). Hiromi et al. used the G(E) function method in conversion of gamma spectrum of in situ Ge(Li) detector(12) and obtained better results than the results using HASL method(13). The precision of the G(E) function is very dependent on the number of input
response spectra, and it is desirable to obtain as many response gamma-ray spectra as possible at appropriately homogeneous energy intervals(14). However, in practice, there are some restrictions in choosing suitable gamma-ray source owing to the limited number of monoenergetic sources with proper energies; only a few response functions can be obtained. Kimiaki et al. used Monte Carlo method to obtain simulated response gamma-ray spectra(14), and break through those restrictions. In a study by Kimiaki (14), the simulated spectrums were used according to the method presented in the literatures(10, 12) to determine the coefficients of the G(E) function. In this paper, based on a 400 ` 400 ` 1600 NaI(Tl) detector and Monte Carlo N-particle Transport Code (MCNP), a new method to determine the G(E) function is proposed, using MCNP to simulate spectrums of various monoenergetic gamma rays and deposited energy, which was then converted into absorbed dose rate in air, in the corresponding window of full-energy peak in an air ball. The value of the G(E) function at energy E in spectra was obtained by dividing the absorbed dose rate in air by counts of corresponding full-energy peak. Curve-fitting software 1st0pt was used to fit the G(E) values vs. E using the G(E) function; the coefficients were determined. In the authors’ method, only data of full-energy peak are needed. Compared with using all the data of all standard spectrums in conventional method, the number of the histories and consumed time in simulation processing are decreased. Meanwhile, complexity of numerical calculation is reduced, compared with that using least square method to process while data of all standard spectrums.
# The Author 2014. Published by Oxford University Press. All rights reserved. For Permissions, please email: [email protected]
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
The G(E) function method is a spectrometric method for the exposure dose estimation; this paper describes a method based on Monte Carlo method to determine the G(E) function of a 400 `3 400 `3 1600 NaI(Tl) detector. Simulated spectrums of various monoenergetic gamma rays in the region of 40 – 3200 keV and the corresponding deposited energy in an air ball in the energy region of full-energy peak were obtained using Monte Carlo N-particle Transport Code. Absorbed dose rate in air was obtained according to the deposited energy and divided by counts of corresponding full-energy peak to get the G(E) function value at energy E in spectra. Curve-fitting software 1st0pt was used to determine coefficients of the G(E) function. Experimental results show that the calculated dose rates using the G(E) function determined by the authors’ method are accordant well with those values obtained by ionisation chamber, with a maximum deviation of 6.31 %.
W. CHEN ET AL.
PRINCIPLE AND METHODS In the algorithm of the G(E) function method, there exists a conversion operator between the count rates and the absorbed dose rate in air at channel i of the gamma spectra(8), which is a constant and could be described as in the following equation: Di ni
GðEi Þ ¼
ð1Þ
end X
Di ¼
i¼b
end X
GðEi Þni
ð2Þ
i¼b
where D is the absorbed dose rate in air converted from the spectra (nGy h21); b is the virtual staring channel of a spectra, the energy at which is the first positive value; end is the last channel of the spectra and other notations mean the same as in equation (1). From equation (2), it could be found that the key to obtain the absorbed dose rate in air is to determine the conversion operator at channel i, i.e. G(Ei). The relationship function between G(Ei) and the energy at channel i is defined as the G(E) function, which is given by a polynomial equation(10, 12, 14). GðEÞ ¼
K max X
Ak ðlog EÞk1
2
@S ¼ @Ak
J P j¼1
2
0
ððDj /Dj Þ 1Þ
¼0
@Ak
ð5Þ
0
where Dj is the calculated dose rate corresponding to the jth standard spectrum derived from the activity of the jth standard source(12) (nGy h21) and J is the number of the standard spectrums. In the literature(12), the number of terms in equation (3), i.e. Kmax, is increased one by one; when the value of S 22 converges to a minimum, the procedure is stopped. In the literature(14), simulated spectrums were obtained by MCNP, and the value of Kmax was studied, but the method for determining Ak is still as equation (5). In this paper, a simpler method is proposed. MCNP was used to simulate gamma spectra with monoenergetic incident gamma rays, but only the data of full-energy peaks are used. The region of the full-energy peak could be determined according to the simulated spectra, and then the deposited energy Einair in an air ball, whose radius is 9.8 cm (because the radius of ionisation chamber in this study is 9.8 cm), in the same region of the same incident gamma ray is obtained through simulation, converting Einair into absorbed dose rate in air as:
ð3Þ DðEÞ ¼ Einair f
k¼1
This equation is the basic form of the G(E) function, in which G(E) is the absorbed dose rate in air per unit count rates at the channel corresponding to energy E in a spectra (nGy h21 cps21); Kmax is the maximum term number of the G(E) function; Ak is coefficient for the kth term of the polynomial equation. From equations (2) and (3), it could be known that the key of using the G(E) function method to obtain absorbed dose rate in air is to determine Kmax and Ak. In the conventional method(12, 14) for determining the G(E) function, some radioactive sources or MCNP were used to obtain standard spectrums or simulated standard spectrums, and whole data of each spectrum are used, as in the following equation: " # K end max X X k1 Ak ðn ji log ðEi ÞÞ ð4Þ Dj ¼ k¼1
@
ð6Þ
where f is a conversion factor converting megaelectron volt into nanogray per hour, f ¼ 106` 1.6021` 1027 ` 3600. Then, the G(E) value at the ith channel corresponding to energy Ei is derived as in the following equation: GðEi Þ ¼
DðEi Þ Pi
ð7Þ
where Pi is the counts of the full-energy peak of incident gamma ray with energy Ei. Equation (8) is derived by substituting equation (7) into (2):
i¼b
D¼
end X i¼b
218
GðEi Þni ¼
end X DðEi Þn
i
i¼b
Pi
ð8Þ
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
where E(i) is the energy at channel i (MeV); ni is the count rates at channel i, cps (counts per second); G(Ei) is the conversion operator (nGy h21 cps21); and Di is the corresponding absorbed dose rate at channel i (nGy h21). Using these operators at different channel of the spectrum, the absorbed dose rate in air could be derived, as in the following equation: D¼
where Dj is the dose rate corresponding to the jth standard spectrum (nGy h21); nji is the count rates in channel i of the jth standard spectrum (cps) and other notations mean the same as in the previous equations. Then, the least square method was used, as in the following equation:
DETERMINATION OF THE G(E) FUNCTION Table 1. fwhm of 400 `3 400 `3 1600 NaI(Tl) detector at various energy. Energy/MeV fwhm/MeV
0.0595 0.0111
0.6617 0.0476
1.1732 0.0663
1.3325 0.0722
According to least square algorithm, equation (10) is used to determine a, b and c. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi2 @ X4 fwhmðE Þ a b Ei þ cEi 2 i i¼1 @xm
Figure 1. Model of measurement.
¼0
ð10Þ
ðm ¼ 1; 2; 3; x1 ¼ a; x2 ¼ b; x3 ¼ cÞ
Figure 2. Model of 400 ` 400 ` 1600 NaI(Tl) detector.
Simulation of standard spectrums and deposited energy in the air ball The model used in this study is shown in Figure 1, and the model of 4L NaI(Tl) detector is shown in Figure 2. In order to determine the window of a full-energy peak, the Gaussian energy broadening (GEB) treatment (15) of MCNP was used, which could be used to better simulate a physical radiation detector, in which energy peaks exhibit GEB. The parameters of GEB in MCNP specify the full width at half maximum (fwhm) of the observed energy broadening in a physical radiation detector(15): pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi fwhmðEÞ ¼ a þ b E þ cE 2
ð9Þ
where E is the energy of incident gamma ray, and the units of a, b and c are MeV, MeV1/2 and none(15), respectively. A radioactive point source containing 241 Am, 137Cs and 60Co was placed at a random location (here, the distance from point source to detector has not been measured, because fwhm represents the energy resolution of a detector, and the relative position between detector and source does not affect it) in front of the 400 ` 400 ` 1600 surface of the NaI(Tl) detector to get the fwhm at 59.5, 661.7, 1173.2 and 1332.5 keV, which are listed in Table 1.
Determination of the coefficients of the G(E) function Substitute the deposited energy in the air ball of different incident gamma rays into equation (6) to calculate the absorbed dose rate in air, then substitute the counts of the full-energy peak and the corresponding calculated absorbed dose rate in air into equation (7) to get the value of G(Ei) at the ith channel corresponding to energy Ei, the values obtained were listed in Table 2. Curve-fitting software 1sp0pt (16) was used to fit the relationship between G(E) value and energy E. There are two parameters needed to be determined by users, the fitting function and the fitting algorithm. Here, the basic form of the G(E) function, as in equation (3), was used as the fitting function; the Universal Global Optimisation (UGO)(16) was chosen as fitting algorithm. The Kmax is determined by increasing k one by one, and the procedure stopped when the correlation coefficient R is most close to 1. The determined coefficients of the G(E) function are listed in Table 3. The fitted curve is shown in Figure 3. The accuracy of the fitted curve, i.e. the accuracy of the G(E) function determined by the proposed method, was validated indirectly by experiment.
219
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
Substituting fwhm in Table 1 into equation (10), a ¼ 20.00557, b ¼ 0.06675, c ¼ 0.00028 were obtained. In this paper, an ionisation chamber was used in validation step. The radius of the ionisation chamber is 9.8 cm, so an air ball with 9.8 cm radius was simulated, and the centre of this air ball is located in the position of the centre of 400 ` 400 ` 1600 NaI(Tl) detector in Figure 1. The deposited energy in the air ball is recorded by f6 card(15) of MCNP. In simulating process, the histories are 109 and the maximum estimated relative error of interested simulated data in full-energy peak is ,5 %.
W. CHEN ET AL. Table 2. The G(E) value at energy E of the spectra. E/KeV 40 50 60 80 100
G(E)/nGy h21 cps21
E/KeV
G(E)/nGy h21 cps21
E/KeV
G(E)/nGy h21 cps21
1.89E-08 8.82E-08 1.90E-07 4.70E-07 7.22E-07
200 300 500 800 1200
2.21E-06 4.00E-06 8.61E-06 1.70E-05 2.99E-05
1600 2000 2600 2900 3200
4.37E-05 5.86E-05 8.34E-05 9.58E-05 1.07E-04
Table 3. Value of coefficients of the G(E) function. Coefficient
Coefficient
Coefficient
Coefficient
Coefficient
21.95E-02 6.87E-02 21.04E-01
A4 A5 A6
8.97E-02 24.75E-02 1.59E-02
A7 A8 A9
23.29E-03 3.84E-04 21.95E-05
Table 4. Absorbed dose rate in air contributed by natural radionuclides by the G(E) function and ionisation chamber.
Figure 3. Fitted curve of the G(E) function of 400 ` 400 ` 1600 NaI(Tl) (correlation coefficient R ¼ 0.999).
Validation According to the layout shown in Figure 1, a 400 ` 400 00 ` 16 NaI(Tl) detector was used to measure the spectrum of radioisotopes in environment on concrete square, asphalt road and grass, respectively. The calculated absorbed dose rate in air by the G(E) function and measured results by the ionisation chamber are listed in Table 4. It is shown that the calculated dose rate in air by the G(E) function and the measured results by ionisation chamber coincide well, with a maximum deviation of 25.4 %. But it could be noticed that the absorbed dose rates in Table 4 are from environmental background radioactivity and ,160 nGy h21. It is necessary to validate that the determined G(E) function is still applicable; if there are artificial isotopes, the absorbed dose rate is higher and the uniformity of gamma-ray field is not the same as in the environmental background radioactivity. For this purpose,
Location
Value by G(E)/nGy h21
Value by ionisation chamber/nGy h21
Relative error/%
Concrete square Asphalt road Grass
156.5
157.4+2.6
20.57
152.8
145.0+4.1
5.38
150.7
144.2+4.9
4.51
a hexagon mixed area source (241Am þ 137Cs þ 60Co þ 152Eu) with 50 cm length of each side was used to get spectrums (the layout is shown in Figure 4). The results calculated by the determined G(E) function were compared with the results by ionisation chamber, which is listed in Table 5. The results in Table 5 showed that the G(E) function determined by the proposed method is still applicable, and the uniformity of gamma-ray field did not obviously impact on the accuracy of calculation of absorbed dose rate in air by the determined G(E) function. CONCLUSIONS The method for the determination of the coefficients of the G(E) function is based on MCNP simulation. In the simulated spectrums, only the energy region and count rates of full-energy peak were focussed on, without analysis of whole data of spectrums. Compared with the conventional method, it reduced the strong dependence on radioactive sources and the complexity of numerical calculation. The results obtained in the validation experiment show that the
220
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
A1 A2 A3
Value
DETERMINATION OF THE G(E) FUNCTION
Table 5. Absorbed dose rate in air by the G(E) function and ionisation chamber (using source). d/m
Without source 0 1 3 5
Value by G(E)/nGy h21
Value by ionisation chamber/nGy h21
Relative error/%
153.6
158.4+2.8
23.03
973.3 373.9 209.2 181.6
980.4+10.8 351.8+4.8 198.4+5.1 176.1+3.9
20.72 6.31 5.42 3.14
relative error of absorbed dose rate in air by the determined G(E) function is ,6.4 %, compared with the results by ionisation chamber. The uniformity of gamma-ray field did not obviously impact on the accuracy of calculation of absorbed dose rate in air by the determined G(E) function. REFERENCES 1. Clouvas, A., Xanthos, S., Antonopoulos-Domis, M. and Silva, J. Monte Carlo calculation of dose rate conversion factors for external exposure to photon emitters in soil. Health phys. soc. 78, 295–302 (2000). 2. Chul-Young, Yi et al. Calculation of spectrum to dose conversion factor for a HPGe spectrometer. J. Korean Phys. Soc. 30, 186 –193 (1997).
221
Downloaded from http://rpd.oxfordjournals.org/ at EPFL Lausanne on March 21, 2015
Figure 4. Layout of the experiment using hexagon area source.
3. Moriuchi, S., Tsutsumi, M. and Saito, K. Development of a dosimetric system using spectrometric technique suitable for operational radiation dose measurements and evaluation. In: Proceeding of The 10th International Congress of IRPA. Hiroshima, May 2000. P-3b-197, 1–8 (2000). 4. Minato, S. and Kawano, M. Evaluation of exposure due to terrestrial gamma-radiation by response matrix method. J. Nucl. Sci. Technol. 7, 401–406 (1970). 5. Okano, M., Izumo, K., Kumagai, H., Katou, T., Nishida, M., Hamada, T. and Kodama, M. Measurement of environmental radiations with a scintillation spectrometer equipped with a spherical NaI(Tl) scintillator. Natural Radiation Environment III, CONF-780422, 2, 896–911 (1980). 6. Urabe, I., Tsujimoto, T. and Katsurayama, K. Estimation of photon energy spectra based on 34 ` 34 response matrix. J. Radiat. Res. 19, 163–170 (1978). 7. Moriuchi, S., Nagaoka, T. and Saito, K. Measurements, theoretical analysis and dose evaluation of low-level external radiation in environment. Report of the Japan-USSR Seminar on Radiation Effects Research. Tokyo, June 1990, National Institute of Radiological Sciences, (1991). 8. Moriuchi, S. and Miyanage, I. A spectrometric method for measurement of low-level gamma exposure dose. Health Phys. 12, 541–551 (1966). 9. Moriuchi, S. and Miyanaga, I. A method of pulse height weighting using the discrimination bias modulation. Health Phys. 12, 1481– 1486 (1966). 10. Moriuchi, S. A New Method of Dose Evaluation by Spectrum-Dose Conversion Operator and Determination of the Operator. JAERI 1209, Japan Atomic Energy Research Institute (1970). 11. International Commission on Radiation Units and Measurements. ICRU REPORT 53. Gamma-ray Spectrometry in the Environment. ICRU Publication 53 (1994). 12. Hiromi, T. et al. Environmental gamma-ray exposure rates measured by in-situ Ge(Li) spectrometer. J. Nucl. Sci. Technol. 17, 281– 290 (1980). 13. Harold, L. B. et al. In situ Ge(Li) and NaI(Tl) gammaray spectrometry. HASL-258, Health and Safety lab., (1972). 14. Kimiaki, S. et al. Calculation of the spectrum-dose conversion operator for a HPGe detector by means of a Monte Carlo simulation. J. Korean Phys. Soc. 25, 504–509 (1992). 15. Los Alamos National Laboratory. MCNP4C Monte Carlo N-Particle Transport Code System User’s Manual. Los Alamos National Laboratory, (2000). 16. 7D-Soft High Technology, Inc. 1st0pt manual. 7D-Soft High Technology, INC., China (2009).
