UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 CONSTRUCTIING THE RESPONSE FUNCTION FOR A BGO DETECTOR USING MCNP5 AND DEVELOPING THE DECONVOLUTION ALGORITHM IN THE LOW GAMMA ENERGY By JANGYONG HUH A THESIS PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE UNIVERSITY OF FLORIDA 200 9 PAGE 2 2 200 9 J angyong H uh PAGE 3 3 To my mother who is in heaven and to my father who has sacrificed everything for me PAGE 4 4 ACKNOWLEDGMENTS The Nuclea r and Radiological Engineering D epartment of University of Florida was very unfamiliar with me when I cam e here for the first time. N ow, I feel it is my home. I think life is a long journey. I have been on a hard but f unn y travel in America. During travel, you come to know the truth: you cannot take one step without other s help no matter how great you are. Theref ore, the first thing you learn is how to appreciate the people who help you keep traveling Fortunate ly, my journey has been wonderful so far even though it was very hard to me sometimes That means I am indeb ted to many people : Professor Allieza Haghighat would not be better for a lifetime teacher as well as an academic supervisor I appreciate my co adv iser, Dr. Baciak for reviewing my paper. I was happy that I knew other faculties and staffs in the department. And also I have made too many good friends to mention in one page: J orge, Vishal, Romel, Benoir, and others Especially I greatly appreciate my roommate Melissa for proofreading my thesis. Along with them, my strongest supporters are my family and friends in Korea. Without their sacrifice, I would not have been here. My mother would be very proud of me even in the heaven. M any events have been m ixed with bad, sad and good things: a bad thing is that my goal in America has yet to be accomplished a sad thing is that my mother passed away two years ago, and a good thing is that I am in good health. I don t need to be disappointed with the past nor too optimistic for the future because my journey is still going on. I will keep my voyage in America with much more efforts until my life ends. In conclusion I might not be that smart because I have already started forgetting what I studied in a class, bu t I will never lose my gratitude of all of them. PAGE 5 5 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ ............... 4 LIST OF TABLES ................................ ................................ ................................ ........................... 7 LIST OF FIGURES ................................ ................................ ................................ ......................... 8 ABSTRACT ................................ ................................ ................................ ................................ ... 11 CHAPTER 1 I NTRODUCTION ................................ ................................ ................................ .................. 13 Charac teristics of the Gamma Spectrum ................................ ................................ ................ 13 General Deconvolution Methods for Gamma Spectrum ................................ ........................ 14 Response Function ................................ ................................ ................................ .................. 17 Bismuth Germinate Oxide (BGO) Detector for the Gamma Spectroscopy ........................... 18 Objectives ................................ ................................ ................................ ............................... 18 2 C ONSTRUCTION O F THE RESPONSE FUNCTION ................................ ........................ 23 Overview ................................ ................................ ................................ ................................ 23 Modeling of the BGO Detector Using the MCNP5 Code ................................ ...................... 26 Characteristics of the MCNP5 Code for Gamma Transport ................................ ........... 27 Crystal Size ................................ ................................ ................................ ...................... 28 Source Position ................................ ................................ ................................ ................ 28 Window Thickness ................................ ................................ ................................ .......... 29 Angular Dependence ................................ ................................ ................................ ....... 30 Comparison of the Experimental and MCNP5 Simulation R e sults ................................ ........ 31 Experimental Setup ................................ ................................ ................................ ......... 31 Energy Broadening ................................ ................................ ................................ .......... 32 Comparison of the Exp erimental Results to the MCNP5 Predictions ............................. 34 Generation of the Response Matrix ................................ ................................ ........................ 36 3 D ECONVOLUTION METHODS ................................ ................................ .......................... 55 Maximum Likelihood Estimation (MLE) ................................ ................................ ............... 57 Expectation Maximization (EM) ................................ ................................ ............................ 59 Maximum Likelihood Expect ation Maximization (MLEM) ................................ .................. 60 Direct Inverse Method & Gaussian Mixture Method (DIM GMM) ................................ ...... 62 Finite Mixture Model (FMM) ................................ ................................ ......................... 62 Gaussian Mixture Model and Parameter Estimation through the EM ............................ 64 Concept of the New DIM GMM ................................ ................................ ............................ 67 PAGE 6 6 4 R ESULTS ................................ ................................ ................................ ............................... 73 Deconvolution of the Sample Spectrum Obtained from the MCNP Simulation .................... 74 Comparison of the Deconvolution E ffec t of the MLEM and the DIM (S ources 1 and 2) ................................ ................................ ................................ ........................... 74 Comparison of the Deconvolution Effect of the MLEM and the DIM GMM (S ource 3) ................................ ................................ ................................ ..................... 75 Deco nvolution of the Sample Spectrum Obtained from E xperiments ................................ ... 76 Comparison of the Deconvolution E ffect of the MLEM and the DIM (S ources 1 and 2) ................................ ................................ ................................ ........................... 76 Comparison of the Deconvolution Effect of the MLEM and the DIM GMM (S ource 3) ................................ ................................ ................................ ..................... 78 5 D ISCUSSION ................................ ................................ ................................ ......................... 88 Response Function ................................ ................................ ................................ .................. 89 Deconvolution Methods ................................ ................................ ................................ .......... 91 6 C ONCLUSION ................................ ................................ ................................ ....................... 94 APPENDIX A MCNP CODE FOR THE BGO D ETECTOR MODELING ................................ .................. 98 B C LANGUAGE CODE FOR ENERGY BROADENING EFFECT ................................ .... 103 C L INEAR INTERPOLATION OF THE RESPONSE FUNCTION ................................ ..... 105 D M ATLAB CODE FOR INTERPOLATION OF THE RESPONSE FUNCTION MATRIX ................................ ................................ ................................ ............................... 108 E M ATLAB CODE FOR THE SPECTRAL DECONVOLUTION USING THE MLEM ALGOR ITHM ................................ ................................ ................................ ...................... 112 F C LANGUAGE CODE FOR THE SPECTRAL DECONVOLUTION USING THE GMM ................................ ................................ ................................ ................................ .... 114 LIST OF REFERENCES ................................ ................................ ................................ ............. 119 BIOGRAPHICAL SKETCH ................................ ................................ ................................ ....... 123 PAGE 7 7 LIST OF TABLES Table page 2 1 Tally of two full energy peaks and their normalized tally by tally of 0 as a function of the angles for the n on shielded case ( 60 Co). ................................ ................................ .. 38 2 2 Tally of two full energy peaks and their normalized tally by tally of 0 as a function of the angles for the lead shielded case ( 60 Co). ................................ ................................ .. 3 8 2 3 List of experiments performed for examining the accuracy of Monte Carlo modeling .... 3 9 4 1 Comparison of deconvolution results of the MLEM and the DIM over the sample spectra which are obtained from MCNP5 simulations for two dif ferent sources, source 1 ( 54 Mn ) and source 2 ( 22 Na ) respectively ................................ .............................. 80 4 2 Comparison of deconvolution results of the MLEM and the DIM & GMM over the sample spectra which MCNP5 simulates a source mixing with 54 Mn 22 Na and 137 Cs (Sou rce 3) ................................ ................................ ................................ ........................... 80 4 3 Comparison of deconvolution results of the MLEM and the DIM over the sample spectra which are obtained from experiments for two different sources, source 1 ( 54 Mn ) and source 2 ( 22 Na ) respectively. ................................ ................................ ........... 81 4 4 Comparison of deconvolution results of the MLEM and the DIM GMM over the sample spectra which are obtained from experiments for a source mixing with 54 Mn 22 Na and 137 Cs (Source 3) ................................ ................................ ................................ 8 1 PAGE 8 8 LIST OF FIGURES Figure page 1 1 T ypical gamma spectrum analyzed according to the spectrum attribute .......................... 2 1 1 2 T ypical forward method for deconvolution ................................ ................................ ...... 22 1 3 T ypical inverse method for deconvo lution ................................ ................................ ....... 22 2 1 Detection mechanisms for an actual scintillator detector and a MCNP5 simulation ....... 40 2 2 Tally and P/T as a function of crystal size ( 137 Cs source) ................................ ................. 41 2 3 Tally and P /T as a function of source positions ( 137 Cs source) ................................ ......... 41 2 4 Change of the solid angle as the source to detector distance increases ............................ 42 2 5 Tally and P/T as a function of detector window thickness ( 137 Cs sou rce ). ....................... 42 2 6 G eometry of a source and a detector is illustrated for the non shield ed case: A 60 Co source is rotated by 20 degree for every simulation ................................ ......................... 43 2 7 G eometry of a source and a detector is illustra ted for the lead shield ed case: A 60 Co source is rotated by 20 degree for every simulation ................................ ......................... 43 2 8 S pectrum of 60 Co as a function of the angles for the non shield ed case ........................... 4 4 2 9 S pectrum of 60 Co as a functi on of the angles for the lead shield ed case .......................... 4 4 2 10 The BGO detector and gamma source are supported with the low density epoxy board ................................ ................................ ................................ ................................ 45 2 11 Experimental setup for measuring the gamma rays from the disk s ource ........................ 45 2 1 2 Three different response function of a 1x1 detector are depicted; MCNP5 without the Gaussian energy broadening, MCNP5 with the Gaussian energy broadening and the experiment ................................ ................................ ................................ ......................... 46 2 1 3 Compariso n of the experimental results to calculated spectra with and without source distances ......... 47 2 1 4 Comparison of the experimental results to calculated spectra with and w ithout source distances ......... 48 2 1 5 Comparison of the experimental results to calculated spectra with and without detector source distances ......... 49 PAGE 9 9 2 1 6 Comparison of the experimental results to calculated spectra with and without source distances ......... 50 2 1 7 C omparison of the experimental results to calculated spectra with and without source distances ......... 51 2 1 8 Comparison of the experimental results to calculated spectra wi th and without source distances ......... 52 2 1 9 Examples of response functions which are constructed using the numerical interpolation with and without the energy broadening ................................ ..................... 5 3 2 20 Comparison of 662 keV spectrum obtained from interpolation and 662 keV spectrum reduced from convolution of the response function and a source ................................ ... 5 4 3 1 H ow statistic methods define source data and measured data and how MLE and EM estimate parameters given a source data ................................ ................................ ............ 70 3 2 S pectr a before and after it is deconvolved. ................................ ................................ ....... 71 3 3 Discrip tions of how the DIM GMM works. ................................ ................................ ..... 72 4 1 C omparison of unfolding ability of two algorithms (the MLEM and the DIM) for the MCNP5 simulated sample spectrum on source set 1, 54 Mn ................................ .............. 82 4 2 Comparison of unfolding ability of two algorithms (the MLEM and the DIM) for the MCNP5 sim ulated sample spectrum on source set 2, 22 Na ................................ ............... 82 4 3 T he MCNP5 simulated sample spectrum for a mixture of 54 Mn, 22 Na, 137 Cs on source set 3 ................................ ................................ ................................ ................................ .... 83 4 4 U nfolding spectrum of source set 3 for the MCNP5 simulated s ample spectrum by the MLEM algorithm ................................ ................................ ................................ ......... 83 4 5 Unfolding spectrum of the MCNP5 simulated sample by the DIM on source set 3 ......... 84 4 6 The GMM separates peaks of the spectrum which was obtained by the DIM at Figur e 4 5 ................................ ................................ ................................ ........................... 84 4 7 Comparison of unfolding ability of two algorithms (the MLEM and the DIM) for the experimentally generated spectrum on source set 1, 54 Mn ................................ ................ 85 4 8 Comparison of unfolding ability of two algorithms (the MLEM and the DIM) for the experimentally generated spectrum on source set 2, 22 Na ................................ ................. 85 4 9 Comparison of the MCNP5 and measured sample spectra for 22 Na and 54 Mn sources ..... 86 4 10 E xperimentally generated sample spectrum for a mixture of 54 Mn, 22 Na, 137 Cs on source set 3 ................................ ................................ ................................ ......................... 8 7 PAGE 10 10 4 11 U nfolding spectrum of source set 3 for the experimentally generated sample spectrum by the MLEM algorithm ................................ ................................ .................... 8 7 4 12 Unfolding spectru m of the experimentally generated sample spectrum by the DIM on source set 3 ................................ ................................ ................................ ......................... 8 7 C 1 Linear expansion of the response function of the low energy ................................ ......... 10 7 C 2 Linear compression of the response function of the l ow energy ................................ ..... 10 7 PAGE 11 11 Abstract of Thesis Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for th e Degree of Master of Science CONSTRUCTIING THE RESPONSE FU NCTION FOR A BGO DETECTOR USING MCNP5 AND DEVELOPING THE DECONVOLUTION ALGORITHM IN THE LOW GAMMA ENERGY By Jangyong Huh May 200 9 Chair: Alireza Haghighat Cochair: James Edward Baciak, Jr. Major: Nuclear Engineering S cience s In our study, all the gen eral procedures necessary for extracting useful information from the measured gamma spectrum of a BGO detector and point radioisotope sources have been discussed with regard to the MCNP modeling and the deconvolution methods T he sensitivity of the BGO det ector was studied by using the MCNP5. From simulation results the MCNP5 code was verified for its feasibility of the spectral analysis. The detector response function is subject to the detector system condition: such as crystal material, crystal size, pho t on energy, and source position. Therefore, various detector responses were modeled on the conditions of the detector system. A response function matrix was constructed on the basis of the MCNP5 simulation results. All the results were compared with the ex perimental ones. In terms of the inverse problem, the deconvolution of the measured gamma spectrum is identical to development of the algorithm for deriving the stable solution from the intrinsically ill posed deconvolution problem. Despite its slow conver gence, d econvolution algorithms based on statistical iteration have recently drawn much attention as the powerful alternative for resolving the above problem. Thus, we have tested two well known statistical deconvolution methods, the MLEM (Maximum Likeliho od Expectation Method) and the GMM (Gaussian PAGE 12 12 Mixture Method) for the deconvolution of the measured gamma spectrum. The new DIM GMM DIM ( D irect I nverse M ethod) and the GMM (Gaussian Mixture M odel) Those methods have shown good performance. PAGE 13 13 CHAPTER 1 INTRODUCTION Characteristics of the G amma S pectrum Gamma rays are high energy photons, which are radiated from de excitation of the nuclei or from sub atomic interactions. Each nuclide has its own characteristic energies corresponding to the energy level of its nuclear state. Due to this property, the gamma energy has been used for a variety of spectral applications. One of the most favorable methods to detect gamma rays is to use scintillation detectors. As with other detector s t he detection mechanism in scintillation detectors relies on indirect measurement. This involve s very complicated processes which consequently cause undesired signals. T he basi c principle can be descri bed as follows : an energetic radiation photon ionizes the detector crystal proportional to its incident energy. Then, ionized electrons produce visual photons through so called fluorescence. T h e se photons are guided into a PMT (Photo Multiplier Tube) and the PMT converts them into an electric signal which contain s information on the radiation source In addition to perturbation from th e detector, other external effect s such as backscatter and cosmic rays contribute to the formation of the gamma spectrum. Th us, the actually measured spectrum makes a more complicated shape which stems from convolution of unwanted peaks and backgrounds as well as photo peaks of interest as shown in Figure 1 1 Accordingly, proper analysis of the gamma spectrum is closely conn characteristics. If the measured spectrum is analyzed in terms of spectr al attributes [1] it can be classified into the following: P hoto peaks : the full energy peak of interest which is generated by the photoelectr ic effect. PAGE 14 14 Compton backgrounds : they are produced by Compton scattering, single and double escape peaks generated by pair production, and x ray escapes created by the photoelectric effect Compton scatters and background radiations : the general gamm a spectrum is inevitably deconvolved with Compton scatters coming from an interaction between the original radiation source and materials surrounding the detector crystal and backgrounds radiations which are due to external factors such as cosmic rays and terrestrial sources Detector resolution : Due to the statistical fluctuations which are primarily related to a collection of visible photons, the full photo peaks broaden to form the Gaussian distribution. Geometry of a source and a detector : Location of a sample source and a detector affects the shape of the spectrum. Because of the above factors, the measured spectrum will not completely coincide wit h the original radiation source T he spectral shape will vary even for the same radi ation source according to the crystal material s and size and the sample detector geometry used. Therefore, many additional processes are required before and after measurements in order to extract only the necessary information from such a complicated spectrum. Gener al Deconvolution Methods for Gamma S pectrum Various unfolding models for the effective data analysis have been introduced depending on objectives of spectral applications [ 2 ]. However, making a correct estimation of the source information from the observed data is not an easy task because obtained data loses its original information in the real detection system due to effects discussed in the previous section. PAGE 15 15 From the classical viewpoint, the methodology to approach the gamma spectrum analysis can be divi ded into two types, forward method and inverse method [ 3 ] [ 4 ]. Mathema tically, it can be expressed as = + ( 1 1 ) where M indicates the measured data, S represents the source data, and R is the detector response function. All measurement errors or uncertainties are symbolized with G iven a model function well characterizing the physical system, th e forward method estimates the true spectrum (source spectrum) by best fitting the model to the measured data. However, the inverse method predicts the true spectrum (source information) given results of the measured spectrum. In the traditional forward m ethod, the procedure f or deconvolution is as follows : first develop a model function which best characterizes the physical system, and estimate the data using the model function, then compare the estimated data with the measured spectrum, and repeat that process by an iteration techni que until it converges. Figure 1 2 depicts the typical mechanism of the forward method. Advantages of the forward method are that it provides stabilit y of a solution and time saving s while its disadvantage is that the forward method has difficulty resolving the overlapped peaks. Total peak area method (TPA) [5], peak search method [6 ] and stripping method [7 ] can be categorized as examples of the forward method. The way to approach deconvolution in the frame of the inverse pro blem is based on correlation of cause (source data) and effect (measured spectrum) in Figure 1 3 Mathematically, it appears to be very simple to solve this inverse problem. H owever, in reality it presents the mathematical challenge because the real world is often ill posed [ 8 ]. The concept of the well posed or ill posed problem was first introduced by Hadamard in his 1902 paper. The wel l posed problem is defined as follows : Given a measured spectrum (M), PAGE 16 16 i) T here exists a solution (S) ii) T he solution is unique iii) T he relation of the solution and the measured data is continuous. If the invers e problem does not satisfy above conditions, it is defined as ill posed. Usually, most in verse problems are ill posed. This means either a solution does not exist or it is not unique unlike the forward problem which provides a unique solution. Moreover, the solution could be unstable even if the inverse problem is well posed. Specially, t his phenomenon is termed as [ 9 ]. In this unstable system, even a small change of the system can cause a large amplification of the original data. Therefore, general algorithms for the inverse problem seek an approximate solution in a numerical way rather than an exact solution. One of those examples is the least squares method (LSM) which was introduced by Carl Friedrich Gauss in 1805 Its basic idea is to find the minimization of a sum of squares of difference between observ ed value and value given by the model: = arg min 2 ( 1 2 ) where M includes measurement error The LSM is the very classical method for deconvolution of the gamma spectrum [ 10, 11 ]. In many cases, it encounters an oscillatory be havior of the deconvolved spectrum due to the ill conditioned system so it often produce s an unphysical result [ 12 ]. Many mathematical algorithms have been developed in order to reduce this oscillatory effect. Regularization is the fundamental method which has improved instability of a solution and accomplished the meaningful approximate solution by introducing additional information [ 13 ]. Eq.1 3 shows the T ikhonov regularization which has been the most commonly used regularization method for ill condition ed problems [1 4 ]. PAGE 17 17 = arg min { 2 + 2 } ( 1 3 ) side of Eq. 1 3 is artificially introduced for enforcing smooth constraints on the system and, as a consequence, obtaini ng the unique stable solution. However this constraint parameter inevitably causes biased estimates. Therefore, the optimized constraint parameter should be determined by taking into consideration the balance between artificial effects and noise probl ems [ 12 ]. A statistical method is an example for the inverse technique [14]. The main idea of an inverse technique is to consider that the measured data is a random variable which follows a Gaussian or Poisson distribution and the technique allows for dete rmination of the necessary parameters associated with these distributions. Thus, they assume that the obtained sample is independent and identical ly distributed Then the best model which characterizes that distribution is designed with parameters. Using the maximum likelihood methods, it estimates the parameters which make that distribution most likely. The Expectation maximization (EM) [ 1 5 ], and the maximum likelihood expectation method (MLEM) [16 ] are typical examples of the statistical methods. Respon se F unction The r esponse of scintillation detectors has been studied for a long time bec ause of their importance to gamma spectr al analysis as well as investigation of detector characteristics [ 17, 1 8, 19 ] An o bse rved spectrum does not illustrate the orig inal data information but rather the distorted data that is perturbed by the detector system. In other words, the measured spectrum could not be used for spectral applications without the deconvolution process. The response function includes information o n the data which is lost through the detection process. Thus the result of the spectral deconvolution greatly depend s on the response function used by the PAGE 18 18 deconvolution algorithm. For that reason, it is important to generate an accurate response function in the spectral deconvolution process In principle, the detector response function is defined as a probability distribution R (E, E ), indicating that a photon source emitted with energy E is measured as a pulse height with energy E. Strictly speaking, the detector response function is involved only in photo n interaction s with in the detector such as photoelectric interaction, Compton scattering, and pair production that are caused by incidence photons of interest. However, the actual response function i ncludes other external factors such as Compton backscatter and cosmic rays [ 20 ] Bismuth Germinate Oxide ( BGO ) D etector for the G am ma S pectroscopy Despite having some drawbacks such as low detection efficiency, and high susceptibility to radiation damag e, a semiconductor detector such as a HPGe (High Purity Germanium) detector has a huge advantage over scintillator detectors in gamma spectroscopy analysis because of its excellent resolution (~2 keV at 1 MeV) [21 ]. However, the high detection efficiency o f a BGO detector and its ease of use outweigh its low energy resolution in a certain physical environment [22 23 ] compared to that of the HPGe detector. Furthermore a certain effective unfolding algorithm even allows the BGO detector to achieve as high a resolution as the HPGe detector [ 24 ]. For that reason a BGO detector of high efficiency has been often selected as the detection s ystem for many applications despite its poor detector resolution [2 5 ]. Objectives The objectives of this thesis are to constr uct a response function for a BGO detector using the MCNP5 code, then find the optimized interpolation technique for the response function over the different detector system conditions based on the response function obtained, and finally develop a new comp utational algorithm for deconvolution of the measured spectrum in the BGO detector for a radioisotope source. A ll the proce s s es ha ve been discussed in this study. PAGE 19 19 The spectrum deconvolution requires a complicated and multistage procedure for the proper da ta anal ysis since the observed spectra arise from convolution of various photon interactions in the detector system. The detector res ponse function contains all th e information concerning transformation of the original source. Therefore, the accurate const ruc tion of the response function requires the most work prior to beginning development of the deconvolution algorithm. In this study, we constructed the resp onse function library using MCNP5 simulation s while t aking into consideration the experimental fac tors that have critical effects on the sensitivity of a detector : the BGO crystal size, the detector radiation source position, radiation energy of the radioactive source. The simulated results were then compared with experimental results, and the feasibil ity of the MCNP5 simulation s were discussed. Based on the MCNP5 simulation results, the response function matrix was built for deconvolution of the measured gamma spectrum. Selection of the deconvolution method depends on the characteristics of the s pectral application s being used. T here has not been such a n algorithm which works well for all spectral applications : for example, l ibrary least square s show s a good performance for identification of the radioisotopes only if an accurate library is provide d [ 2 6 ] The total peak method is a good alternative unless the full energy peaks are densely overlap ped in the measured spectrum [2 7 ]. As long as the accurate response function is given, the peak stripping method is a very powerful method [ 2 8 2 9 ]. D eco nvolution algorithms based on statistical iteration have recently drawn much attention especially in the field of medical imaging [ 30 3 1 ]. Some researchers have obtained a good result PAGE 20 20 when applying it to the energy spectrum analysis Especially Meng et a l. have demonstrated the MLEM algorithm is superior to other deconvolution methods[ 24 ]. Our study introduced a newly developed deconvolution method and compared it to the MLEM algorithm. The new deconvolutio n method consists of two steps In the first st ep, the ill posed response function with energy broadening is decomposed into a response function matrix without energy broadening and an energy broadening mat rix and measured spectrum is reproduced through the inverse calcul ations of the response matrix without energy broadening. In the second step, this preprocessed spectrum is deconvolved using the GMM (Gaussian Mixture Method) [32] PAGE 21 21 Figure 1 1. T ypical gamma spectrum analyzed according to the spectru m attribute: 1) indicates the ideal spectrum of Cs 137, 2) illustrates the real gamma spectrum, and 3) 5) show the spectrum shape when the real spectrum is decomposed corresponding to physical effects. PAGE 22 22 Figure 1 2. T ypical f orward method fo r deconvolution Figure 1 3. T ypical i nverse method for deconvolution PAGE 23 23 CHAPTER 2 CONSTRUCTION OF THE RESPONSE FUNCTION Overview One of the goals of the current research is the develop ment of a response function f or a BGO detector using the MC NP5 code. The detector response function is subject to various physical factors and surrounding environments such as crystal materials, crystal size, detect or source distance and so on. This implies that a different detector response function matrix needs to be constructed whenever the detection system is changed. It is highly irrational and inefficient because it takes enormous time and effort to build the response function matrix for every case. If sensitivity of the detector response has property of line arity over change of the physical parameter, we can build a correlation of the sensitivity and the parameter change and thus estimate the response function of other parameters based on that correlation. For example, suppose we have constructed response fun ctions of two different crystal sizes, 1x1 inch and 3x3 inch. We could interpolate the response function of a crystal size of 2x2 inch without additional sim ulation or experiment by finding a correlation between the detector sensitivity and the crystal siz e. For the first time, we simulated performance of the BGO detector using the MCNP5 code under various physical conditions in the detection system in order to examine feasibility of the MCNP5 code in model ing the BGO detector [APPENDIX A] After that a r esponse function matrix was constructed of 1500 different energies between 1 keV to 1500 keV on the basis of the MCNP5 simulation results Note that t he response f unctions were produced for a few energie s and the response functions for the rest of energies were obtained using an interpolation technique. PAGE 24 24 In general, three methods are known for generating the response function: the experimental method [ 33 34 ], the semi empirical method [ 35 36 ] and the Monte Carlo simulation [ 37 38 ] method. Typically, the experimental method is not practical for development of the response function because the number of single mono energetic radioactive sources required is limited. In addition, most of the radioactive sources emit more than one single gamma energy, which m akes it more difficult to estimate the detector response for a single energy From the 19 70s through the 19 90s the semi empirical approach had been widely used where the analytical function is established regarding all the possible physical characteristics of a mono energetic photon interaction. I ts parameters are calculated from the least square fits to several response function spectra obtained from experiments, and then this semi empirical function is used to generate the response function of other gamm a energies. W ith the advent of remarkable development s in the computer technology and the particle transport algorithm technique, Monte Carlo simulation codes such as MCNP5 [ 39 ] and GEANT4 [ 40 ] have been favored for construction of the response function si nce it saves time and effort in production of the response function. In the frame of the mathematical approach, the observed spectrum can be expressed in the integral form such as: = ` ` ` + ( ) 0 ( 2 1 ) where M ( E ) denotes the me asured spectrum at energy E S ( E` ) indicates the radioisotope source at energy E` and R ( E,E` ) represents the response function operator which transfers a photon from spectrum S ( E` ) to spectrum M ( E ). ( E ) is the measurement error. In a radiation detecti on system energetic photons are recorded through a series of electronic processing. The electric signal whose amplitude is proportional to the photon energy deposited in the detector crystal is generated from the interactio n of an incident gamma ray and PAGE 25 25 t he crystal material. Each signal is digitized by the ADC (Analog Digital Convertor) and then counted in registers (channels) corresponding to the digitized pulse height [4 1]. Consequently, the actual measured spectrum is represented in a discrete form The degree of discretization depends on the number of the ADC channels used. Considering this characteristic of the radiation detection system, the integral form of Eq. 2 1 can be rewritten as the summation of the discrete energies : = ` ` + ` = 0 ( ) ( 2 2 ) Eq. 2 2 can be expressed i n matrix notation as follows : = + ( 2 3 ) w here M is a vector ( ) S is a vector ( ) is a vector ( ) and R is a matrix ( ) The dimension of vector M is not ne cessarily the same as that of vector S but we treat only the same dimension for simplicity of a problem. Thus, the dimension of t he response function matrix R is assigned a n x n matrix. From a viewpoint of the imaging sy stem, the detector response matrix would be the point spread function [ 42 ] which convolves with a point source and spreads it out. The complicated shape of the observed spectrum is subject to the feature s of the response function matrix, i.e., the detecto r system. If the response function matrix can be decomposed according to causes of the convolution, then there is an obvious merit in handling a deconvolution problem. Unlike the experiment based response function matrix the simulation based response func tion matrix c an be decomposed into two parts: the response function matrix withou t energy broadening a nd the energy broadening matrix. PAGE 26 26 1 2 = 1 1 1 2 1 2 1 2 1 0 2 1 0 0 0 0 0 0 1 0 0 0 0 0 1 2 ( 2 4 ) = 1 1 1 2 1 0 2 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 1 1 2 0 0 0 0 2 1 2 1 0 0 0 0 2 1 0 0 0 0 0 0 0 1 0 0 0 0 0 1 2 ( 2 5 ) = + = + ( 2 6 ) where is the detector response matrix without the energy broadening and B is the energy broadening matrix. The reason for splitting the r esponse function into two matrices is associa ted with the methodology for the spectrum deconvolution. The detector response is generally ill posed at Eq. 2 4, but this ill posed response ma trix can be transformed into a well posed response matrix by the decomposition process. As shown in Eq. 2 5, th e response matrix without energy broadening is a triangul ar matrix when it is decomposed; implying that it is well posed. In the present stud y, we produced two types of response function matrix es : with energy broadening and without energy broadening. One i s used for the ML EM method and the other for a new combination of the DIM and the GMM. Modeling of the BGO D etector Using the MCNP5 C ode Primary features of the gamma s pectrum are closely related to the detector properties and geometr y of the detection sy stem. Unfortunately many features are not simulated by the MCNP5 PAGE 27 27 code. F or example, flat continuum and asymmetry of the tail of a photo peak which arise from leak age of electrons and detector imperfection within the crystal [20, 43 ]. Neverthel ess, the MC NP5 code has shown good performance for many detector related simulations. This is because in many case s such detail ed description of the detector physics is not that critical for model ing a detector. There fore, we are interested in critical factors that affect the detector sensitivity such as crystal size, source position detector window thickness, and source detector angle Such a study will contribute to the efficient development of a computational algorithm for construction of the response function ma trix using the MCNP5 code and, at the same time, provide in depth understanding of detection mechanis m of the scintillation detector. In this section, characteristics of the gamma spectrum in the BGO detector are modeled c onsidering aforementioned physical parameter s. The detector performance is a ss essed with the peak to total ratio (P/T). The P/T is a physical quantity referred to as the ratio of the full energy peak count to the total count recorded by a detector, which is generally used to measure the d etector intrinsic efficiency [ 44 ]. Characte ristics of the MCNP5 C ode for Gamma T ransport As shown in Figure 2 1, a gamma ray produces visible light through so called scintillation process within the BGO detector. Those visible lights are low energy photon s whose number is proportional to energy deposited. The energy of a scintillation photon is around 480 nm or 2.58 eV, and the number of scintillation yields is around 9000 photons/MeV for the BGO crystal. The reflector guides them into the PMT where they a re converted into electrons whi ch contribute to the electronic signal. On the other hand, the MCNP5 code is not designed to transport the visible photons. Once the gamma ray is determined to interact with a crystal by the Monte Carlo algorithm, the energy deposited in the crystal is regarded as absorbed at that spot without further transporting the visible photons that are induced from the interaction of the gamma rays and the PAGE 28 28 detector crystal as illustrated in Figure 2 1 b Thus, the MCNP5 code does not acc ount for uncertainties which are involved in scintillation yields, leak of photons, and quantum electrons while they contribute to the energy resolution in the real detector. For this reason, the MCNP5 simulation shows a salient difference from a spectrum obtained from an experiment. The energy broadening effect should be treated after a simulation in order to compare to the experimental output. Crystal S ize Three gamma spectra were generated for three different crystal sizes for a 137 Cs point source loc ated at 4 cm away from the center of the detector window. Their crystal diameters and heights are 1 x 1 2 x 2 and 3 x 3 The P/T (peak to total ratio) is affected much more by the crystal volume to surface ratio rather than the crystal size [ 45 ]. The most influential factor on the P/T is the escape of the Compton scattered photons which lead to enhancing the Compton background. From the geomet rical viewpoint, a large volume to surface ratio suppresses the Compton background and increases the f ull energy peak. Therefore, if the volume to surface increases the P/T increases. Figure 2 2 shows the P/T obtained from the MCNP5 simulation increases (from 0.606 to 0.764) as the crystal volume to surfa ce ratio enlarges. Source P osition For this analys is we are considering a 137 Cs is otropic point source located at 2 cm, 10 cm and 20 cm from the center of the detector window. Figure 2 3 shows the spectra obtained in these three configurations. The shape of a spectrum is highly dependent on the change of the source to detector distance since the source to detector distance has a large effect on the absolute efficiency. On the other hand, t he P/T does not depend significantly on the source to detector distance while it is closely connected to the incident energy, the volume to surface ratio, and the detector dimension PAGE 29 29 [ 45 ]. Such a feature is observe d in Figur e 2 3. As a source is located fa rther away from the detector, intensity of the full energy peak notably decreases (from 0.1404 to 0.0006) due to reduct change significantly, i.e. 0.76 to 0.797 A slight increment of the P/T can be explained by the solid angle difference in F igure 2 4 T he gamma rays which are more parall el to the centerline of the crystal at the lower solid angle result in less leakage out of the crystal and consequently, less Compton background and higher P/T. Window T hickness Three detector window s of thicknesses (0.02 0.04 and 0.12 ) were consi dered for a 137 Cs point source located at 0.1 cm from th e center of the detector window. For t he detector window only aluminum was considered because the reflector of visible photons is transparent to high energy photons. The most prominent feature due t o increment of the window thickness is increment of the spectrum intensity in th e multi scattering region. A considerable number of photons lose their energy in the detector window before being absorbed into the detector crystal. As a result, the photon s l osing their energy contribute to the spectrum intensity in the multi scattering region Figure 2 5 shows the expected window effect based on the MCNP5. As the detector window is thickened, the multi scattering intensity increases, but the photo peak area d ecreases by 35.4 % from 0.16689 to 0.10775. In comparison to the change in photopeak intensity, the change in the window thickness did not make as a large impact on the P/T as on the photopeak intensity for the range of thickness PAGE 30 30 tested P/T decreases f rom 0.662 to 0.611 when the thickness of the detector window is increased. Angular D ependence In order to simulate the angular effect in the photon detection, five angles (0 20 40 60 and 80 ) were selected Unlike in the previous models where pho tons were incident on the front of a BGO detector, in the angular dependence simulation photons can be project ed on to the side of the detector as well as its front face Thus, two cases were considered in Figure s 2 6 and 2 7 : one is to treat photons incide nt on both the front and the side face s and the other is to treat only photons incident on the front face For implementing the la tter model, a couple of lead bricks were placed on the side face of the detector window to block photons from being proj ected on the side face in Figure 2 7 60 Co was considered as a point source located at 25 cm from the center of the detector window and was rotated by 20 degree. On the non shield ed model, the intensity of the full energy peak gradually increase s corresponding to increment of the angle as shown in Figure 2 8 This result arises from the fact a large angle allows for more incident photons. Table 2 1 shows that the peak intensity of two peaks increases by 26 % for an 80 degree ang le as compared to 0 degree angle case A slight difference of the intensity between two peaks (1173.3 keV and 1332.5 keV) can be attributed to the detector efficiency The lead shield ed model shows the apparent discrepancy between the non shield ed spectrum and the shield ed spectrum. Firs t, t he x ray peak arising from the lead shield is observed which is not found in the non shield ed spectrum as shown in Figure 2 9 Secondly, the full energy peak intensit y diminish es as the incident angle is enlarged on the contrary to the non shield ed mod el. This complete different phenomenon stems from the fact that t he decrease of the PAGE 31 31 solid angle gives rise to reduction in the incident photons. Compared to the peak intensity of an an gle of 0 degree, it decreases by 8 5 % at an 80 degree angle. The peak in tensity hardly changes by an angle of 20 degree, but a rate of change in the peak intensity is much higher at the large angle than that in the non shield ed model as shown in Table 2 2 P/T was not treated in this simulation because of the mathematical diff iculty in calculation of the number of the incident photons. Comparison of the e xperimental and MCNP5 simulation results E xperimental S etup Several experiments were done to test the accuracy of the MCNP5 modeling results. Also, various experimental setups were presented to study a ny correlation between the detector sensitivity and considered physical parameters: source type ( 137 Cs, 60 Co), detector size distance (2 cm to 20 cm) between the radioactive source and the BGO detector and angles between the rad ioactive source and the BGO detector centerline (0 to 80 ) Table 2 3 summarizes four experimental sets which are presented for the variety of the physical conditions corresponding to the MCN P5 modeling. Experimental sets 1, 2, and 3 are involved in the study of detector performance, where two isotopes, three BGO crystals of the different sizes, and fourteen different geometries are utilized. E xperimental set 4 is designed for two purposes: v erifying the response function obtained from the MCNP5 modeling and provid ing a test spectrum for spectral deconvolution. The modeling set 4 will be discussed in Chapter 3. Three different BGO scintillators were selected which have the crystal sizes of 1 x 1 2 x 2 and 3 x 3 The same detector window thi ckness of 0.05 cm was considered for all detectors. Four radioisotopes were employed: 137 Cs, 60 Co, 22 Na, and 54 Mn, whose activities are 1 20% Ci [46]. Figure 2 10 depicts the experimental PAGE 32 32 setup in which a BGO detector and a radioisotope source are supporte d with the low density epoxy board. The data collection and processing follow the typical data acquisition system as illustrated in Figure 2 11 The BGO detector was connected to the standard NIMs (Nuclear Instrumentation Modules) [41] and the data was co llected using a multi channel analyzer (MCA) and ORTEC Maestro software. The energy calibration was performed using a 32 keV x ray and a 661.6 keV photopeak s emi t t ed from 137 Cs and 1173.2 keV and 1332.5 keV photopeak s emitt ed from 60 Co in off line analysis T o compare the experimental results to the MCNP prediction the total count of each channel w as normalized by the total number of activation counts of the radioactive source during the given acquisition time. The BGO detector was isolated from its surrou ndings in order to inhibit external backgrounds as much as possible. However, a measured spectrum intrinsically includes following three types of background counts: external backgrounds which originate from any other gamma source except for the radioisotop e source, backscatter due to surrounding materials, and the electronic noise. In this experimental analysis, only the external background was taken into account. External background measurements were performed before and after the main experiment. Their av erage value was subtracted from the measured data. Energy B roadening The measured data is a result of convolution of many physical phenomena and full y absorbed energies of the radioactive sources: Compton scattering, backscatter, cosmic ray, electric noi se, and energy broadening. Since the energy broadening which is involved in statistical fluctuation makes a significant impact on the convoluted spectrum, the energy resolution effect should be accurately treated for modeling the BGO detector. However, MCN P5 PAGE 33 33 does not provide the transport of visible photons which are the primary attribute of the energy broadening. In the current study, the energy broadening was artificially estimated after the MCNP5 simulation on the assumption that the full energy peak fo rms a Gaussian distribution considering the following formula [ 47 ]: = + + 2 ( 2 7 ) where FWHM denotes full width half maximum Note that coefficients of a, b, and c may change depending on the detector type and experimental setup For an example, parameters for a detector of 2x2 BGO crystal are determi ned as a = 4.87125, b = 2.87358, and c = 0.00026. The Gaussian energy broadening effect can be applied to each tally scored in the pulse height tally (F8) by sampling from the following Gaussian function: ( ) = ( ) 2 ( 2 8 ) wher e E is the energy broadening, E o indicates the energy without broadening effect and C represents the normalized constant. The variance ( 2 ) of the energy broadening is obtained from the following equation: = 2 2 ( 2 9 ) Actually, the F8 tally presents the Gaussian energy broadening option (GEB) in MCNP5 [ 39 ]. We, however, did not use the GEB option. Instead, we made the Gaussian energy broadening code which includes the above FWHM formula because of two objectives [APPENDIX B] F irst, the reso lution correction after the MCNP simulation provides a better efficiency for building a large number of the response functions. Secondly, for deconvolution the response function PAGE 34 34 matrix needs to be decomposed into the response matrices without and with the energy broadening. Figure 2 1 2 describes the result of the Gaussian energy broadening correction of the MCNP5 simulation, where three different detector responses are compared to each other: i) MCNP5 without Gaussian energy broadening, ii) MCNP5 with Gaus sian energy broadening, and iii) the experiment. Originally, the MCNP5 code did not give a complete description on the exponential tail in the low energy side of the Gaussian energy peak. Such a deformation of the Gaussian peak can be attributed to the fol lowing; incomplete charge collection due to crystal imperfection, and the leak out of Auger electrons and x ray escape [43, 48 ]. Only the x ray escape is taken into account in the MCNP5 simulations. However, Figure 2 1 2 exhibits that the energy broaden MC NP5 and experimental results are in a good agreement except for backscatter, which is not treated in the MCNP5 modeling. Comparison of the E xperiment al R esults to the MCNP5 P redictions The spectra obtained from the MCNP5 simulation were compared with the experimental results for different BGO crys tal sizes, for different source detector positions for different detector source angles, and for different radioactive source s Figure s 2 1 3 (a f) and Figure s 2 1 4 (a f) show the results of m odeling set 1 which correspond to the BGO crystal sizes of 1x1 and 3x3 for the 137 Cs, and for detector source distance between 2 cm and 20 cm. A small difference is observed around the Compton edge where the spectrum of the experiment is a little bit higher than that of the M CNP5 simulation result. This can be attributed to the Compton backscattering, which i s not modeled in the MCNP5 simulation. This effect magnifies at the crystal size of 3 x 3 because increment s of the d etector size diminish the multi scattering. PAGE 35 35 Figu re s 2 1 5 (a f) and Figure s 2 1 6 (a f) s how the results of modeling set 2 which correspond to the BGO crystal sizes of 1 x 1 and 3 x 3 for the 60 Co source and for detector source distance between 2 cm and 20 cm. In the modeling set 2 for the crysta l size of 3 x 3 the full energy photopeak of the experiment is found low er than that of the MCNP5 prediction in Figure s 2 1 6 (a f). This might be associated with the sum peak It has been reported in literature [49 ] that the loss of the peak efficienc y due to the coincidence summing effect is around 3% for 60 Co. As the number of radiation counts increases the number of radiation count s in the modeling set 2 (two gamma lines, 1173.2 keV and 1332.5 keV) is twice in the modeling set 1 (one gamma line, 6 61.6 keV) for the same elapsed time a large amount of signal processing loads are imposed on the detection system. It, consequently, leads to increment of a data loss due to the sum peak of 1173.2 keV and 1332.5 keV. Such reasoning is supporte d by the fa ct that difference in the full energy peaks between the two results maximizes at the shortest detector source distance (2 cm) where more radiation photons are incident on the detector. The results of modeling set 3 are given in Figures 2 1 7 (a e) and F igur es 2 1 8 (a e), which correspond to the BGO crystal size of 2 x 2 for the 60 Co source, and six detector source angles. In modeling set 3, t he full energy peak of the experiment do e s not exhibit as large difference from that of the MCNP5 prediction for m odeling of set 2. Since the distance between a source and a detector (25 cm) for the modeling set 3 is longer than the distance used for the modeling set 2 (2 cm to 19 cm), the data loss due to the sum peak is not considered large enough to cause the inc onsistency between the experimental and the MCNP5 result. On the whole, t he MCNP5 simulations accurately reproduce the angular dependence characteristics of the BGO detector for two cases, the non shield ed model and the lead shield ed model PAGE 36 36 For the mode l set s 1 2, and 3 all the results show that the Monte Carlo predictions are in good agreement with experiment s The visible light self absorption within the scintillator crystal becomes a dominant factor as the crystal size increases. This behavior canno t be represented by MCNP5 since visible photons are not transported. Nevertheless, the reason that MCNP5 and the experimental results show a good agreement is due to the fact that lo ss of the visible photons during particle transport does not contribute to counts but to the detector resolution. I n other words, the energy broadening of the MCNP5 result is estimated by using the value of the detector resolution obtained from the experimental values. Generation of the R esponse M atrix We constructed the respon se function matrix for the BGO detector with a crystal size of 1 x 1 where the mono energetic source is located 15 cm away from the center of the detector window. Considering the energy range of the sample radioactive sources (Maximum gamma energy is 1332.5 keV emitting from 60 Co), a matrix dimension of the response was develope d as 1500 x 1500 whose elements have 1 keV bin. The response function relies on the geometry of the model as well as on the intrinsic features of the detector. It means that a new response function matrix is required whenever the detector or th e geometry changes. It is a very time consuming work to build the response functions of 1500 different energies. Thus, the numerical interpolation was adopted to cope with that problem [ 50 51 ]. 14 different energies from 200 keV up to 1500 keV at interval s of 100 keV were selected and the response function corresponding to each energy was simulated using MCNP5. Below 2 00 keV, two intervals were used: 10 keV and 20 keV. T hen, the response f unctions of energies between two selections were linearly interpolated [APPENDIX C APPENDIX D ] The response function matrix obtained from the MCNP5 modeling and the PAGE 37 37 linear interpolation does not include energy broadening. Therefore, the Gaussian code was used for building the response matrix with the energy broadening Figure 2 1 9 shows se veral response functions with and with out the energy broadening, which are linearly interpolated. To verify the response matrix built, the monoenergetic energy of 662 keV, which is the equivalent of gamma lines resulting from 137 Cs, was convolved with the response matrix and its result was compared with the spectrum directly produced by the MCNP5 simulation. They are con sistent with each other in Figure 2 20 PAGE 38 38 Table 2 1 Tally of two full energy peaks and their normalized tally by tally of 0 as a function of the angles for the n on shield ed case ( 60 Co). Non shielded 2x2 crystal Source position Peak area Normalized by value of 0 (degree) 1173.2 keV 1332.5 keV 1173.2 keV 1332.5 keV 0 0.00065 0.00060 1.00 1.00 20 0.00067 0.00062 1.03 1.03 40 0.00071 0.00065 1.09 1.09 60 0.00077 0.00070 1.18 1.17 80 0.00083 0.00075 1.26 1.26 Table 2 2 Tally of two full energy peaks and their normalized tally by tally of 0 as a function of the angles for the lead shield ed case ( 60 Co). Lead shielded 2x2 crystal Source position Peak area Normalized by value of 0 (degree) 1173.2 keV 1332.5 keV 1173.2 keV 1332.5 keV 0 0.00066 0.00060 1.00 1.00 20 0.00066 0.00060 1.00 1.00 40 0.00058 0.00054 0.89 0.90 60 0.00039 0.00037 0.60 0.61 80 0.00010 0.00010 0.15 0.16 PAGE 39 39 Table 2 3 List of experime nts performed for examining the accuracy of Monte Carlo modeling Set Source BGO Detector, crystal sizes Distance between a detector and a source 1 Cs 137 2 cm, 5 cm, 10 cm. 15 cm, 19 cm, and 20 cm 2 Co 60 2 cm, 4 cm, 5 cm, 9 cm, 10 cm, 14 cm, 15 cm, 19 cm, and 20 cm 3 Co 60 25 cm Angles between the detector centerline and a source (0 20 40 60 and 80 ) 4 Na 22, Mn 54, Cs 137 15 cm PAGE 40 40 Figure 2 1. Detection mechanisms for an actual scintillator detector and a MCNP5 simulation. PAGE 41 41 Figure 2 2. Tally and P/T as a function of crystal size ( 137 Cs sourc e). Figure 2 3. Tally and P/T as a function of source positions ( 137 Cs source) PAGE 42 42 Figure 2 4 Change of the solid angle as the source to detector distance increases: the smaller a solid angle, the more parallel the incident photon. Figure 2 5. Tally and P/T as a function of detector window thickness ( 137 Cs source) PAGE 43 43 Figure 2 6 G eometry of a source and a detector is illustrated for the non shield ed case: A 60 Co source is rotated by 20 degree for e very simulation. Figure 2 7 G eometry of a source and a detector is illustrated for the lead shield ed case: A 60 Co source is rotated by 20 degree for every simulation. PAGE 44 44 Figure 2 8 S pectrum of 60 Co as a function of the angle s for the non shield ed case. Figure 2 9 S pectrum of 60 Co as a function of the angles for the lead shield ed case. PAGE 45 45 Figure 2 10 A BGO detector and a gamma source are supported with the low density epoxy board. Figure 2 11 Experimental setup for measuring the gamma rays from the disk source. PAGE 46 46 Figure 2 1 2 Three different response function of a 1 x 1 detector are depicted: MCNP5 without the Gaussian energy broadening, MCNP5 with the Gaussian e nerg y broadening and the experiment. PAGE 47 47 (a) Detector source distance of 2 cm (b) detector source distance of 5 cm ( c) Detector source distance of 10 cm (d) Detector source distance of 15 cm (e) Detector source distance of 19 cm (f) Detector source distance of 20 cm Figure 2 1 3 Comparison of the experimental results to calculated spectra with and without broadening for BGO crystal of size 1 x 1 for different d etector source distances PAGE 48 48 (a) Detect or source distance of 2 cm (b) detector source distance of 5 cm ( c) Detector source distance of 10 cm (d) Detector source distance of 15 cm (e) Detector source distance of 19 cm (f) Detector source distance of 20 cm Figure 2 1 4 Comparison of the experimental results to calculated spectra with and without broadening for BGO crystal of size 3 x 3 for diff erent detector source distances PAGE 49 49 (a) Detector source distance of 2 cm (b) Detector source distance of 5 cm (c) Detector source d istance of 10 cm (d) Detector source distance of 15 cm (e) Detector source distance of 19 cm (f) Detector source distance of 20 cm Figure 2 1 5 Comparison of the experim ental results to calculated spectra with and without broadening for BGO crystal of size 1 x 1 for diff erent detector source distances PAGE 50 50 (a) Detector source distance of 2 cm (b) Detect or source distance of 5 cm (c) Detector source distance of 10 cm (d) Detector source distance of 15 cm (e) Detector source distance of 19 cm (f) Detector source distance of 20 cm Figure 2 1 6 Comparison of the experimental results to calculated spectra with and without broadening for BGO crystal of size 3 x 3 for diff erent detector source distances PAGE 51 51 (a) At the angle of 0 (b) At the angle of 20 (c) At the angle of 40 (d) At the angle of 60 (e) At the angle of 80 Figure 2 1 7 Comparison of the experimental results to calculated spectra with and without broadening for BGO crystal of size 2 x 2 for differ ent detector source angles without a shield PAGE 52 52 (a) At the angle of 0 (b) At the angle of 20 (c) At the angle of 40 (d) At the angle of 60 (e) At the angle of 80 Figure 2 1 8 Comparison of the experimental results to calculated spectra with and without broadening for BGO crystal of size 2 x 2 for different detector s ource angles with a lead shield PAGE 53 53 Figure 2 1 9 E xample s of response functions which are constructed usi ng numerical int erpolation with and without ener gy broadening PAGE 54 54 Figure 2 20 Comparison of 662 keV spectrum obtained from interpolation and 662 keV spectrum reduced from convolution of the response function and a source PAGE 55 55 CHAPTER 3 DECONVOLUTION METHOD S The purpose of deconvolution in gamma spectroscopy is to extract the original source information from the measured data which is perturbed by the detector system. From the viewpoint of the inverse problem, the largest difficulty in deconvolution arises from the non exist ence of a unique solution or instability of a solution due to an ill posed or an ill condition ed problem [ 52 ]. Th e se troubles hinder a deconvolution method from deriving the exact solution, and instead lead it to seek the physically meaningful approximate solution within permissible uncertainties. There are two types of approaches to handle with this approximate problem: deterministic deconvolution and statistical deconvolution [ 53 ]. A basic idea of the deterministic method is to construct a model and fit the observed data to estimate the unknown source data. It looks for the m ost probable parameter s by data matching of the measured data and calculated data. As discussed in the previous chapter, LSM or regularization method can be regarded as the determini stic deconvolution method [ 12 ]. In contrast, the statistical model assumes that the observed data is distributed from a random sampling. A model is constructed which characterizes the d istribution and its associated parameters. This model is called the pr obability distribution function. Once the probability distribution function is built, it seeks estimators, i.e., parameters which make the measured spectrum most likelihood [5 4 ]. Both methods do not ex hibit a large difference from the point of model para meterization. The two models are actually related with each other. The distinctive difference comes from interpretation of the model over the observed data. The former seek s the optimum parameter of PAGE 56 56 t he model to best describe the observed spectrum, while t he la tter find s the parameter of the highest probability distribution model to generate the observed spectrum. The biggest advantage of the deterministic approach is its simplicity to treat a pro blem but it has not provided a satisfactory outcome for hi gh dimensionality of the data. This drawback is mostly attributed to the lack of information. The regularization method makes considerable improvement by introducing additional information. However, determination of a functional and a regularization parame ter considerably limit s its application to the spectral deconvolution [ 8 ]. By contrast statistical deconvolution resolves aforementioned issues more efficiently than the deterministic model. Also, it is an appropriate methodology to treat nuisance pro blems such as errors in the measured data and uncertainties of a model since its principle is the inference estimation based on statistics concepts, i.e., variance and means square error. Many methodologies have been developed for the deconvolution of a measured spectrum in terms of statistics: maximum entropy m ethod (MEP) [55], genetic algorithm (GA) [56 ], maximum likelihood exp ectation maximization (MLEM) [16 ] and e xpectation maximization (EM) [57 ]. Statistic al methodologies for deconvolution are based on statistical theory such as maximum l ikelihood estimation (MLE) [ 58 ] and Bayesian model In the present study, we have utilized the MLEM for deconvolution since it has already demonstrated good performance on spectral deconvolution. Gaussian mixture mo del (GMM) ha s been introduced for development of a new model. The GMM has been coupled with the direct inverse method (DIM) for the deconvolution o f the gamma spectrum. We term this new methodology od and Gaussian mixture m ethod ( DIM In t his C hapter we w ill discuss the statistic theories which are applied to the MLEM and the GMM. PAGE 57 57 Maximum Likelihood Estimation (MLE) The objective of spectr al analysis is to infer quantitative or qualitative information from the measured data. I n gamma spectroscopy such information cannot be obtained directly from the measured spectrum. The maximum li kelihood estimation (MLE) is a fundamental data analysis t ool which is utilized for extrac ting information of interest from the observed data. This method is based on a statistical approach which is used to estimate parame ter s of a model. It was originally proposed by Sir Roland A. Fisher in 1912 [ 5 9 ] and has been widely used for various applications including the analysis of gamma spectra [ 60, 61 6 2 ] The basic proced ure of the MLE is illustrated in Figure 3 1. When the source data is collected through measurement, two different types of data are generated: measured data and missing data. The measured data is referred to as incomplete data in terms of statistics, which implies it loses some information when the sou rce data is transformed into measured data. This is why the source data is termed as complete data which can be recovered via only the indirect method. The missing data which is unobserva ble [22] and is often called the hidden or latent data is used by the expectation maximization algorithm which will be discussed later The MLE deals with the measured data or incomplete data. The MLE assumes that the measured spectrum is obtained from a sampling of random variables which are independent and i dentically distributed. Based on this assumption, a statistical parameterization model called a probability density function is constructed The density function characterizes the distribution of the random variables: 1 2 ( 3 1 ) Because the observed samples are independent the probability density function ( pdf ) can be expressed as a product of the probabilit ies of each data point; i.e, PAGE 58 58 1 2 = 1 1 = = 1 ( 3 2 ) The pdf specifies characteristics of the probability that the sample distribution is Y = 1 2 given parameter The Poisson distribution or Gaussian di stribution has been often us ed as the pd f [ 63 ]: = = 1 ( 3 3 ) where denotes the mean of the Poisson distribution and = 1 2 = 1 ( ) 2 2 2 ( 3 4 ) where and indicate the mean and the standard deviation of the Gaussian distribution. Unlike the deterministic method which looks for parameters by fitting a model to the measured data, the MLE method finds the parameters taking advantage of a statistical quantity, i.e., the lik elihood function [ 5 5 ]. Mathematically, the probability density function and the likelihood function are treated as the same: 1 2 1 2 ( 3 5 ) However, the likelihood function searches for the estimator which makes the probability density the highest for an observed data. The estimator is obtained by maximizing the likelihood fu nction i.e., ` = max ( 3 6 ) L og likelihood is often used for calculation of estimators instead of the likelihood for ease of calculation: PAGE 59 59 ` = max log or log = log = 1 = 0 ( 3 7 ) The value ` calculated from Eq. 3 7 is called the maximum likelihood estimator (MLE). Expectation Maximization (EM) The methodology of the MLE presents some advantages in terms of the fact that it uses advanced statistical inference techn ique s Thus, it has been applied to various models and different types of data. Nevertheless, some problems still arise in parameter estimations. Usually, the likelihood function is mathematically intractable to maximize the likelihood or find this likelih ood estimator. This drawback is mainly due to high dimensionality of a data and non linearity of a model. In addition, information of interest is the unobserv able complete data or missing data rather than unknown paramete rs of the incomplete data in gamm a spectroscopy. It means that the incomplete data log likelihood in the MLE does not provide the source data information, and also the complete data log likelihood in the MLE cannot be used because the source data information is not given. The MLE algorith m is not an advisable choice for proble m s for which the measured data is incomplete Therefore, i n order to treat these problem s, expectation maximization (EM) algorithm is used. The EM associates the observed incomplete data with the complete data in a w ay that estimation of the likelihood function is computational ly tractable The general formulation of the EM algorithm was introduced by Dempster, Laird, and Rubin in 1977 [ 15 ]. Its most salient feature distinguishing from the MLE is that it maximizes t he complete data log likelihood via a numerical iteration technique. The general procedure of the EM algorithm is performed by two successive steps, E (Expectation) step and M ( Maximization ) step. The E step is the process to PAGE 60 60 relate the incomplete data and the complete data. G iven the observable incomplete data (Y) and initial parameters ( 0 ) t he conditional expectation value of the complete data log likelihood is defined as + 1 = log  + 1 ( 3 8 ) In the M step, one it erates on the value of parameter until a maximum value of Q is achieved, i.e., + 1 = max ( 3 9 ) Maximum Likelihood Expectation Maximization (MLEM) The maximum likelihood expectation maximization (MLEM) is the special cas e of the EM algorithm. The MLEM algorithm was proposed by Shepp and Vardi [ 16 ] for reconstruction of medical image s in emission tomography. I t has been used in the field of medical imaging [ 31 6 4 ] and has exhibited good performa nce in spectral deconvolut ion [ 24 6 5 ]. The main principle is explained by the following procedure: suppose that the measured spectrum is subject to a Poisson distribution t he pdf characterizing the measured data is constructed on the basis of a Poisson distribution. In the MLEM method the likelihood of the pdf is maximize d using the MLE algorithm. A numerical iterative equation is derived from this MLE process. This equation is identical to the expectation step of the EM method. The final step is to maximize the equation as the EM method does. In the MLEM, it is important to note how to relate the measured data (incomplete data) with unknown source data (complete data) in order to obtain the best estimator. The MLEM defines the source data as the estimator, and utilizes the g iven response function which is associated with the source data by following equation : PAGE 61 61 s d = s b R ( b d ) = 1 ( 3 10 ) w here s(b) indicates the number of photons of the energy bin b emitted from the radioactive source, and the response function mat rix R(b,d) represents the probability that a photon of the energy bin b emitted from the source is detected at energy bin d in the detector. s*(d) is the mean of the photons detected at the energy bin d in the detector The mathematical theory is dev eloped via the following procedure: the pdf is expressed as m s = ( ) ( ) m ( d ) m ( d ) = 1 ( 3 11 ) where m*(d) denotes the measured spectrum at the energy bin d in the detector. The log likelihood associated with the pdf is given as log = log m s = ( ) ( ) m ( d ) m ( d ) = 1 ( 3 12 ) To estimate s corresponding to maximum of the log likelihood, we d ifferentiate Eq. 3 12 with respect to s(b) : log ( ) s ( b ) = ( ) = 1 + ( ) ` ( ` ) ` = 1 = 1 = 0 ( 3 13 ) M ultiplying both sides of Eq.3 13 by s(b), then ( ) log ( ) s ( b ) = ( ) = 1 + ( ) ` ( ` ) ` = 1 = 1 = 0 ( 3 14 ) Eq.3 14 can be rewritten for an iteration scheme that converges to a maximum likelihood: PAGE 62 62 + 1 = ( ) = 1 ( ) ` ( ` ) ` = 1 = 1 ( 3 15 ) In the transition from Eq. 3 14 to Eq. 3 15, the term s s g = s g b ` R ( b ` d ) ` = 1 is referred to as the expectation step, and the process to estimate + 1 is referred to as the maximization step. Direct Inverse Method & Gaussian Mixture Method ( DIM GMM) Finite Mixture M odel (FMM) Finite mixture model (F MM) is a methodology to analyze dat a which is characterized by a random distribution, using statistical modeling [ 32 ]. It has been around one century since the biometrician Karl Person used mixture models for analysis of a data which is mixed of two differ ent normal probability density functions. The enormous computational load s in fitting the mixture model have made it difficult to be applied to a problem with a large number of variables (high degree of freedom) even though the FMMs have huge advantages of modeling the complex mixture distribution. However, with the advent of fast er computer technology, FMMs have been applied to numerous fields, which are involved in the computational pattern recognition, biol ogy, engineering, marketing and so on [ 66, 67 6 8 ] The FMMs are used to identify unobserved heterogeneity in the data, or to cluster underlying groups in the data. To this end, finite mixture models assume that observed data are a mixture of random distributions, and each distribution is characteri zed by the probability density function with its own parameters. A model of the observed data (probability density function) is mathematically defined as = k = 1 ( 3 1 6 ) PAGE 63 63 w here = 1 denotes a random dat a set of size N, and indicates the probability density function of k th component characterized with parameters of a vector The parameter is the k th mixing coefficient or the k th mixing weight whose physical meaning is probability t hat observation comes from the k th component. It has the following features: 0 1 ( = 1 ) = 1 = 1 ( 3 17 ) In the specified vector space the quantity containing all the unkno wn parameters can be defined as = a 1 a k 1 1 k ( 3 18 ) In short, Eq. 3 1 6 implies that the observed data can be expressed as a mixture of the probability density distributions with different parameters on the assumption that the observed data is distributed at rand om. Once the probability density model is constructed for the observed data, the next step is to find unobservable parameters in order to infer which cluster the observed data originated from One advantage of FMMs is that they can model the appreciabl y complicated distributions using the simple equation Eq. 3 1 6 However, it is still an open issue to estimate parameters in the fitting of the FMM to the ob served data. Several methods have been introduced for parameter estimation: minimum distance method s [6 9 ], maximum likelihood [ 70 7 1 ], Bayesian approaches [ 7 2 ] and so on. In fact, the explicit methodology for parameter estimation in the FMM does not exist, but the EM method has been known as an effective and straightforward algorithm for pa rameter est imat ion in the FMM [32 ]. PAGE 64 64 Gaussian Mixture Model and Parameter E stimation through the EM The Gaussian mixture model is a specified case of f inite mixture model (FMM) where the observed data distribution is a mixture of Gaussian distribution s Figure 3 2 a depicts the observed spectrum which is represented as the probability density function, in the mixture model. This mixture distribution can be deconvolved into three Gaussian peaks of different characteristics, variance, mean, and peak amplitu de by taking advantage of the GMM if the observed spectrum originates from a convolution of th e three Gaussian distribution s Figure 3 2 b The GMM starts with the assumption that the observed spectru m is independent and identically distributed i.e., Y = y 1 y 2 y N where N is the sample size and its components are characterized by the Gaussian distribution. Then, each Ga ussian distribution is given as = = 1 2 2 1 2 2 ( ) 2 = 1 ( 3 19 ) where and 2 are mean and variance of the k th component respectively. and the EM is applied to the GMM. The log likelihood should be determined by the corresponding probability density function: log = log = = 1 log a k = 1 = 1 = log a k = 1 = 1 ( 3 20 ) Eq. 3 20 denotes the incomplete data log likelihood. The estimation of the parameter is mathematically a challenge for this incomplete log likelihood function b ecause the equation includes a sum of the log function. Therefore, it would be a better idea to take advantage of the PAGE 65 65 complete log likelihood for the parameter estimation which is mathematically tractable The complete data set for the GMM can be given as Z = Y H = y 1 h 1 y 2 h 2 y N h N ( 3 21 ) where Z is the complete data set consisting of the observed data set Y= y 1 y 2 y N and the unobserved data set H = h 1 h 2 h N wh ich is called hidden data or missing data. In the GMM the hidden data (h i ) is the component identity of the observable data (y i ): that is, it tells which component is linked with the observed data. The complete da ta log likelihood is written as log = log Y H = log = = 1 log = 1 ( 3 22 ) where = a = ( a ). With introduction of the additional information, the hidden data H, the estimator calculation will be mathematically more tractable. Then, the E step of the EM method is given as g = log  g = g log = = 1 = 1 ( 3 23 ) where g = a k g k g g = a k g k g a k g k g = 1 ( 3 24 ) Note that E q. 3 24 If Eq. 3 23 is developed into the simpler equation, then g = log  g ( 3 25 ) PAGE 66 66 = log a k g + = 1 = = 1 log k g = 1 = = 1 Because each parameter is independent from each other, we can maximize Eq. 3 25 by differentiating it with respect to an individual parameter. = max ( 3 26 ) a k a k = 0 ( 3 27 ) k k = 0 ( 3 28 ) k k = 0 ( 3 29 ) From t he above maximization process the new es timators are obtained as follows : + 1 = 1 = 1 ( 3 30 ) + 1 = = 1 = 1 ( 3 31 ) + 1 = = 1 ( + 1 ) 2 = 1 ( 3 32 ) where the hidden value is replaced by k. The new estimator calculation begins with guess of the initial parameters and is repeated until it converges to the expe cted tolerance. PAGE 67 67 Concept of the N ew DIM GMM The DIM GMM algorithm consists of two steps, the direct inverse of the detector response function, and the peak separation of the Gaussian mixture method. The most straightforward and the easiest way to solve th e inverse problem is to calculate the inverse of the response function matrix. However, the detector response matrix is ill posed as previously mentioned several times Some algorithms reach the approximate solution by grouping the ill posed response matr ix into several matri c es which provide stable solution and then by estimating the solution using an iterative scheme. Van Cittert Algorithm [ 7 3 ], and gold deconvolution [ 7 4 ] are those examples. This method, however, has intrinsic limit due to the fact that the decomposition of the response function matrix is performed only in terms of mathematics without taking into account the physical situation. Thus, it may cause artificial effect. By contrast the direct inverse method decomposes the response function m atrix into two matri c es corresponding to physical characteristics: the detector response function matrix ( ) without the energy broadening and the response function matrix ( ) with energy broadening as shown at Eq. 2 6 Physically, t he matrix is i nvolved in any effect which produces distortion of the source data and the matrix B is the detector energy resolution resulting from statistical uncertainty of the collection of visible photons in the detector and photoelectrons generated in the photo mul tiplier tube (PMT). From the mathematical point of view, t he response function matrix is well posed even though the energy broadening matrix B is still ill posed. If the inverse of is performed on both side of Eq. 2 6 then 1 = + ( 3 33 ) PAGE 68 68 1 = = ( 3 34 ) is still subject to the energy broadening effect, but all other effects are removed except for the photo peaks of interest. As a result, Compton backgrounds and x ray escape peaks are all eliminated, and the original data intensity of the source is recovered in Figur e 3 3. This i s due to properties of the response function which includes all the physical information regarding the distortion of the original source information. The reconstructed source data can be treated just as a mixture of the Gaussian peaks corresponding to the energy broadening. This assumption is reasonable as long as the measured data is large enough to fo rm the Gaussian distribution. In the second step, the reproduced spectrum is deconvolved again using the GMM (F igure 3 3 ), in which Eq. 3 3 0 to Eq. 3 32 are used to estimate three different parameter s and In this spectral application, we can benefit from the fact that the energy resolution is given by t he experimental results. In other words, the number of unknowns can be reduced. The reduction in the number of unknowns diminishes the number of parameter estimation s. This means it save s time taken for the parameter estimation and reduces error propagati on of the parameter estimation. As a result, it leads to a large decrease in conversion time and enhancement of the parameter estimation accuracy. In order to apply this advantage to the parameter estimation, estimation equations Eq. 3 31 and Eq. 3 32 sh ould be redefined. The full wi d th half maximum (FWHM) is obtained from the experiment over all energies. = + + 2 ( 3 35 ) where a, b, and c are the user defined coefficients which are determined from the experiment. The detector resolution is related to the FWHM and the standard deviation as follows : PAGE 69 69 = = 2 35 (3 36) From Eq. 3 3 6 the detector resolution is decided over all energies. Then, dependency between t he mean energy and the sta ndard deviation is rewritten as = 2 35 ( 3 37 ) Therefore, t he Gaussian distribution is redefined as = = 1 2 2 1 2 2 ( ) 2 = 1 = 1 2 2 35 1 2 2 35 2 ( ) 2 = 1 ( 3 38 ) If we r ewrite Eq. 3 25 using the above Gaussian distribution and m aximize it with respect to then we obtain + 1 = = 1 + = 1 2 + 4 2 = 1 2 = 1 2 = 1 ( 3 39 ) where = PAGE 70 70 Figure 3 1 H ow statistical methods define source data and measured data and how MLE and EM estimate parameters given a source data. PAGE 71 71 Figure 3 2. S pectr a before and after it is deconvolved : a) S pectrum which is mixed with the Gaussian distribution b) C omponent s of the s pectrum after it is deconvolved PAGE 72 72 Figure 3 3. Descriptions of how the DIM GMM works. The complex sample spectrum is deconvolved by the DIM and becomes the mixture of the Gaussian peak in the upper figure. In the lower figure, the preproce ssed spectrum is separated by the GMM algorithm. PAGE 73 73 CHAPTER 4 RESULTS Two types of sample spectra were generated by MCNP5 simulations and experiments. They were employed for demonstrating performance of deconvolution algorithms, the MLEM algorithm [A PPENDIX E ] the DIM, and the DIM GMM [APPENDIX F ]. The detection sy stem for demonstration was composed of a BGO det ector of 2x2 inch crystal size and a source The detector and the source were positioned 15 cm away from each other Three sources 54 Mn 22 Na and 137 Cs with activity of 1.0 Ci were considered. Th e sources used in the experiment w ere encapsulated in the plastic sea led disk of 2 cm diameter Th e encapsulated source contains a small spot in epoxy in the center of the hole The actual area of the spot is not well defined but can be 2 4 mm in diameter depending on how the liquid wicks during pipetting. However, the source wa s treated as a point source in the MCNP5 since the source area is small In the experiment the BGO detector was connecte d with the standard NIMs (Nuclear Instrumentation Modules) [41] and the data was c ollected for each sample for 180 minutes which is regar ded as to be enough time to achieve good s tatistics. The collected data wa s normalized by the total number of photons emitted from a radioactive source for a given time in order to compare the deconvolution effect of the MCNP5 prediction and the experimentally obtained spectrum Three different sources were considered: source 1 is 54 Mn which includes one gamma line (834 .8 keV), source 2 is 22 Na which consists of two gamma lines (51 1 keV and 1274.5 keV) and source 3 is mixed with three sources 54 Mn, 22 Na and 137 Cs which emit one x ray (32.2 keV) and four gamma lines (511 keV, 661.6 keV, 834.8 keV and 1274.5 keV). S pectr a 1, 2 and 3 were generated orderly from each source using the MCNP5 simulation and spectr a 4, 5, and 6 wer e produced from the experiment over the same source s as shown in Figures 4 1 to 4 3 and Fi gures 4 7, 4 8, and 4 10 respectively Through decon volution of these PAGE 74 74 spectra three most important performances of deconvolution method s were examined separately: i) recover y of the original source intensity ii) peak searching (source identification) and iii) peak separation For deconvolutions of spect r a 1 and 2 (or spectr a 4 and 5 in the experiment) which have one single photo peak and two separated photo peaks respectively, the MLEM and the DIM were employed in order to verify their ability of the source intensity recovery and the peak searching Note that the accuracy of the DIM GMM strongly depends on the recovery ability of the DIM for the source intensity. T he DIM GMM w as not performed since peaks of these spectra were not closely packed so as to require for the peak separation. Accordingly, peak separation ability of the deconvolution methods used was not discussed on spectr a 1 and 2 S pectrum 3 (or spectrum 6 in the experiment) is of a rather complicated shape which includes the overlapped peak. T he DIM i s inappro priate for deconvolution of tha t spectrum due to its limit of peak separation. Therefore, in order to use the high ability of the DIM in the source intensity restoration for deconvolution of spectrum 3, one more step the GMM ( Gaussian mixture model ), was introduced which is able to sep arate the overlapped peak s The deconvolution effects were assessed by comparison of two met hods, the MLEM and the DIM GMM. Deconvolution of the Sample Spectrum Obtained from the MCNP S imulation C omparison of the Deconvolution E ffect of the MLEM and the DIM (source s 1 and 2) As depicted in Figure s 4 1 and 4 2, the Compton backgrounds were remarkably removed by both algorithms, the DIM and the MLEM. This means that the two algorithms improve the signal to noise ratio. In terms of ability to recover the sou rce intensity, both methods recover ed the radiation strength of the MCNP5 simulated sample spectra within 1 % of the true value for spectrum 1 Unfolding of spectrum 2 shows the similar results as that of spectrum 1 even if a PAGE 75 75 slight increase of the relativ e error was observed. Table 4 1 compares true and estimated values by deconvo lution Regarding the detector resolution, the MLEM me thod show s even a better energy resolution than the actual BGO detector. For example, the resolution of the BGO detector use d is 11.6 % at 834.8 keV while t he MLEM leads to a significant reduction in resolution up to 2.3 % at 834.8 keV of spectrum 1 T he DIM does not exhibit any discrepancy of the resolution between before and after the simulation This inefficiency of the DIM in improvement of the energy resolution was anticipated because of its intrinsic problem that the spectrum contains the energy broadening after deconvolution. However, such a drawback is not an issue for the convolved spectrum consisting of separated peak s since the DIM successfully removes background effects such as Compton scattering, which distorted the original source spectrum, and the source intensity can be calculated by the simple integration of a peak area in F igure s 4 1 and 4 2. Comparison of th e Deconvolution E ffect of the MLEM and the DIM GMM (source 3) Over the complex spectrum 3 resulting from source 3 of more gamma lines, the deconvolution ability of the MLEM method declines compared to the case of spectra 1 and 2. A couple of small noise peaks which were not found in deconvolving spectra 1 and 2 are observed around the low energy region (~100 keV) of the reconstructed spectrum in Figure 4 3 Table 4 2 indicates the relative errors of the reconstructed source intensity considerably increas e for the all peaks. The best source restoration of spectrum 3 (7.7 % relative error) occurs at the peak of 1274.5 keV which is positioned away from the overlapped peaks. It implies that the deconvolution of the MLEM method works be tter on the separable ph otopeak s than the close ly packed photopeaks. T here is no large variation in energy broadening for the deconvolved spectrum 3. The energy resolution however, becomes a little bit lower than in the deconvolution of spectr a 1 and 2 ( 2.3 % at 834.8 keV for th e deconvolved spectrum 1 but 2.5 % at 834 keV for PAGE 76 76 the deconvolved spectrum 3 ) The MLEM method shows a good performance on the peak searching. It identifies all the peaks within 1~2 keV. On the other hand, the DIM GMM exhibits a better performance for unf olding of the spectrum 3 than the MLEM as shown at Table 4 2. The reconstructed source intensity is in a good agreement with the true value, within 1.3 % of the relative error. Those values are similar to what were obtained in deconvolution of spectr a 1 an d 2 which are the relatively simple r spectra. It is important to note that, as shown in Figures 4 5 and 4 6 the DIM GMM did not produce noise peaks which were found in the low energy region s in the MLEM deconvolution approach Each f ull energy peak posit ion coincides with the characteristic energy of the source s with uncertainties of ~1 keV. Deconvolution of the sample spectrum obtained from experiments Comparison of the Deconvolution E ffect of the MLEM and the DIM ( source s 1 and 2) Unlike the case of th e MCNP5 simulated sample spectrum, the experimentally generated sample spectr a show un desired peaks in t he low energy region in F igure s 4 7 a and 4 8a This discrepancy arises from the backscattering effect which was not considered in the response function generation. A backscattering peak remains in both folded spectra after the deconvolving process and several small spurious peaks are observed in the spect ra deconvolved by the MLEM ( Figures 4 7b and 4 8b). The recovery of the source intensity of the meas ured sample spectrum significantly deteriorates compared to the deconvolution of the MCNP5 simulated sample spectrum. Especially, the spectrum 5 of two photopeaks exhibits the worse performance than spectrum 4 of one single photopeak. As shown in Table 4 3 t his undesirable e ffect is much higher on the photopeak closer to the energy region of backscatter regardless of the deconvolution algorithm PAGE 77 77 used. For example, the relative error of the recovered source intensity for spectrum 5 increases to ~30 % at the gamma line 551 keV. The magnitude of relative error is almost twice as high as that of the gamma line 1274.5 keV. The two main r easons for a large relative error is inconsistency of the measured sample data and the MCNP5 simulated sample data, and the bac kscatter effect. As shown in Figure 4 9 two peaks (551 keV and 1274.5 keV) of the measured spectrum for the 22 Na source are slightly larger than those of the MCNP5 simulated spectrum while the measured spectrum for the 54 Mn source is in a good agreement w ith the MCNP5 simulation spectrum. One of t he reason s that the experimenta l result is higher than expectation, might be related to uncertainty of the activity of the source (1 20% Ci ). Consequently this discrepancy wa s reflected on the deconvol v ed spectr um during the deconvolution process The other reason is the backscatter effect. Even though the backscatter effect does not directly influence the deconvolution of the measured photopeak it is apparently attributed to the poor recovery of photopeaks of interest For the 834.8 keV peak of Mn 54 source and the 1274.5 keV of the 22 Na source backscatter does not contribute to the intensity of those peaks However backscatter has a significant effect on the intensity of 551 keV of the 22 Na source In other words, the measured 551 keV peak of the 22 Na source contains contribution resulting from backscatter before deconvolution. Therefore, th at increment due to backscatter is transferred to the intensity of the recovered 551 keV through the deconvolution proce ss. Table 4 3 supports this theory. The relative error of the deconvolved 551 keV peak positioned around the backscatter energy region is two times higher than that of the decon v ol v ed 1274.5 keV peak whose measured sample spectrum does not include the back scatter effect even though two peaks originate from the same source, 22 Na PAGE 78 78 Because of this reason, it is important to treat the backscatter effect in generation of the response function. The deconvolution ability of t he MLEM and the DIM is not directly as sociated with the backscatter effect, but they cannot recover the true intensity of the measured spectrum when the measured spectrum includes multi peaks : that is, when the full energy peaks of the measure spectrum are located around the backscatter energy region. Therefore, i n order to prevent worsening the recovery performance of the measured spectrum due to the backscatter, the response function is required to be constructed by considering the backscatter effect. Consequently, backscatter was not removed by both of two deconvolution methods, and the MLEM created the unpredicted spurious peaks in the low energy as illustrated in Figures 4 8b and 4 9b. The deconvolved photopeaks of 22 Na produced a large relative error in the peak intensity since the mea sured photopeak w as intrinsically contaminated by backscatter On the other hands the recovery of 54 Mn in which the intensity of the photopeak is not involved in backscatter shows a good agreement between the deconvolved peak intensity and the expected tru e value. Comparison of the Deconvolution E ff ect of the MLEM and the DIM GMM ( source 3 ) In deconvolving spectrum 6, the MLEM produces several unexpected peaks around the backscatte r energy region, which were not found in the spectrum deconvolved from the MCNP5 simulated sampl e spectrum. Those spurious peaks increases both in size, and in number compared to deconvolution of spectra 1 and 2 in Figure 4 11 It implies that the unfolded spectrum may develop into the more unpredictable form if the original sou rce consists of more complex gamma lines. It is worth nothing that the MLEM is not able to resolve the backscatter effect, and makes some undesired spurious peaks, but nevertheless its source intensity recovery over the PAGE 79 79 experimentally generated sample spe ctrum does not show a large discrepancy compared to the MCNP5 simulated sam ple spectrum ( Tables 4 2 and 4 4). As the first process of the DIM GMM, the DIM was performed for recovery of the source intensity of spectrum 6. The recovered spectrum was obtaine d in Figure 4 12. However, the GMM failed to separate the recovered spectrum by the DIM. This is due to the fundamental principle of the GMM which is based on the assumption that the spectrum should be a mixture of the Gaussian peak. The deconvolved spect rum obtained from the first process is convoluted with five Gaussian shape peaks except for backscatter whose shape cannot be defined. Therefore, the GMM was expected to be able to separate at least the Gaussian shape peaks within the relative error simila r to that of the MLEM. The outcome, however, was not good. It indicates the GMM is very sensitive to the peak shape. PAGE 80 80 Table 4 1. Comparison of deconvolution results of the MLEM and the DIM over the sample spectra which are obtained from MCNP5 si mulations for two different sources, source 1 ( 54 Mn ) and source 2 ( 22 Na ) respectively Peak area (peak strength) of deconvoluted spectrum MCNP5 Source 1 ( 54 Mn) Source 2 ( 22 Na) Simulated 834.8 (keV) 511 (keV) 1274.5 (keV) True value 0.999 0.643 0.357 DIM 0.997 (0.2 %) 0.642 (0.2 %) 0.355 (0.6 %) MLEM 0.994 (0.5 %) 0.641 (0.3 %) 0.353 (1.1 %) Figures in the parenthesis indicate the relative error of the reconstructed source peak. Table 4 2. Comparison of deconvolution results of t he MLEM and the DIM GMM over the sample spectra which are obtained from MCNP5 simulations for a source mixing with 54 Mn 22 Na and 137 Cs (Source 3) Peak area (peak strength) of deconvoluted spectrum ( Source 3) MCNP5 Peak area Peak position Simulated 511 661.6 834.8 1274.5 511 661.6 834.8 1274.5 True value 0.396 0.274 0.093 0.22 511 661.6 834.8 1274.5 MLEM 0.451 0.304 0.104 0.237 509 660 834 1274 (13.9 %)* (11.1%) (11.3%) (7.7%) DIM & GMM 0.396 0.278 0.094 0.22 511 660 632 1275 (0.1%) (1.3%) (0.8%) (0.1%) Figures in the parenthesis indicate the relative error of the reconstructed source peak. PAGE 81 81 Table 4 3. Comparison of deconvolution results of the M LEM and the DIM over the sample spectra which are obtained from experiments for two different sources, source 1 ( 54 Mn ) and source 2 ( 22 Na ) respectively Peak area (peak strength) of deconvoluted spectrum Experimentally Source 1 ( 54 Mn) Source 2 ( 22 Na) obt ained 834.8 (keV) 511 (keV) 1274.5 (keV) True value 0.999 0.643 0.357 DIM 1.015 (1.6 %) 0.841 (30.8 %) 0.411 (15.1 %) MLEM 1.009 (1.0 %) 0.833 (29.5 %) 0.408 (14.3 %) Figures in the parenthesis indicate the relative error of the reco nstructed source peak Table 4 4. Comparison of deconvolution results of the MLEM and the DIM GMM over the sample spectra which are obtained from experiments for a source mixing with 54 Mn 22 Na and 137 Cs (Source 3) Peak area (peak strength) of deconvol uted spectrum ( Source 3) Experimentally Peak area Peak position obtained 511 661.6 834.8 1274.5 511 661.6 834.8 1274.5 True value 0.396 0.274 0.093 0.22 511 661.6 834.8 1274.5 MLEM 0.442 0.252 0.092 0.231 512 662 836 1276 (11.7 % ) (8.0 % ) (1.4 % ) (4.6 % ) DIM & GMM N/A N/A Figures in the parenthesis indicate the relative error of the reconstructed source peak PAGE 82 82 Figure 4 1. Comparison of unfold ing ability of two algorithms (the MLEM and the DIM) for the MCNP5 simulated sample spectrum on source set 1, 54 Mn Figure 4 2. Comparison of unfolding ability of two algorithms (the MLEM and the DIM) for the MCNP5 simulated sam ple spectrum on source set 2, 22 Na PAGE 83 83 Figure 4 3. The MCNP5 simulated sample spectrum for a mixture of 54 Mn, 22 Na, 137 Cs on source set 3 Figure 4 4. Unfolding spectrum of source set 3 for the MCNP5 simula ted sample spectrum by the MLEM algorithm PAGE 84 84 Figure 4 5. Unfolding spectrum of the MCNP5 simulated sample by the DI M on source set 3 Figure 4 6. The GMM separates the unfolding spectrum ob tained by the DIM at Figure 4 5 PAGE 85 85 Figure 4 7. Comparison of unfolding ability of two algorithms (the MLEM and the DIM) for the experimentally generated spectrum on source set 1, 54 Mn Figure 4 8. Comparison of unfoldi ng ability of two algorithms (the MLEM and the DIM) for the experimentally generated spectrum on source set 2, 22 Na PAGE 86 86 Figure 4 9 Comparison of the MCNP5 and measured sample spectra for 22 Na and 54 Mn sources: Two spectra are in a good agreement for 54 Mn but the measured sample spectrum is a little bit larger than the MCNP5 simulated sample spectra for 22 Na. Especially, this phenomenon is larger for 551 keV peak due to backscatter. PAGE 87 87 Figure 4 10. E xp erimentally generated sample spectrum for a mixture of 54 Mn, 22 Na, 137 Cs on source set 3 Figure 4 11. Unfolding spectrum of source set 3 for the experimentally generated sample spectrum by the MLEM algorithm Figure 4 12. Unfolding spectrum of the experimentally generated sample spec trum by the DIM on source set 3 PAGE 88 88 CHAPTER 5 DISCUSSION The in depth study for deconvolution of the low energy gamma spectrum (< 5 MeV) has been presented on the basis of the following three steps: A profound understanding of the detector sensitivity: the measured spectrum reflects physical features of a detector as well as information of the incident gamma rays. The simulation code should implement such a detailed physical phe nomenon. It is a very fundamental task for the spectral deconvolution. The MCNP5 code wa s verified for its feasibility of the spectral analysis. Construction of the detector response function for modeling: all the deconvolution algorithms depend on the response function. Therefore, generation of the correct detector response function is the crucial factor to regulate performance of the deconvolution methods. Regarding the response function, another important factor is to accumulate various detector resp onse functions because t he detector response function varies according to detection system conditions such as the incident gamma energy and geometry as well as the detector crystal property H owever, it takes a lot of time and effort to construct the respo nse function corresponding to a ll the detection environments. To overcome this difficulty, the most effective technique is to gen erate the response functions for several selected detection conditions, and then obtain other response s via interpolation D evelopment of the deconvolution tool: due to the intrinsic ill posedness of deconvolution problems, there is no deconvolution method which is able to extract the exact solution. Therefore, the deconvolution technique inevitably looks for the approximate so lution. The main issue is how to generate a stable and unique solution. The statistical method for PAGE 89 89 deconvolution treats this problem efficiently in spite of its slow convergence. The MLEM method and Response F unction The response function plays a key role in deconvolution of the gamma spectrum. Therefore, many researches have been performed on construction of the response function over years The e xperimental methods and the semi empirical method s which were most widely used, have experimen tal difficulties and inefficiency of use but t hey have advantage in that the response function c an be generated without a deep understanding of the complicated p hysics in the detector since the se method s employ the i nterpo l a tion and the analytic fun ction for generation of the response function on the basis of the mono energ etic photon spectrum obtained from experiments By contrast the Monte Carlo simulated method such as the M CNP5 obviate s those problems that two methods encounter In order to reproduce the true response, the MCNP5 should be able to simulate all physical features of the photo n interaction in the detector. MCNP5 does not give a detailed description on some physics such as the exponential tail an d flat continuum due to imperfection of the crystal or electron leakages [43 ], asymmetry of the full energy peaks due to the non linear ity of the detector response [ 7 5 ], and other effects (dead time, temperat ure sensitivity, etc) [ 7 6 ]. Moreover, energy re solution, one of the biggest characteristics on the detector, is not treated in the MCNP5. Nevertheless, the MCNP5 prediction wa s in a good agreement with the experimental results as shown in Figures 2 1 3 to 2 1 8 Such a result can be explained by two fact s: first, m inor features of the detector response mentioned above have a relatively small effect on peak intensity and peak position of interest. Secondly, t he energy br oadening is constructed using the PAGE 90 90 experimental data. In other words t he ene rgy broaden ing effect is determined by the experiment r ather than the MCNP5 simulation and t he response features, which are not simulated, are offset by the energy broadening effect This is why the simulation results have show n a good agreement with the experimenta l results even if the MCNP5 d id not model some physical features of a detector. It is worth noting that only one issue, i.e. loss of a count due to the dead time is not still accounted for. As shown in Figure 2 1 6 the peak amplitude is slightly higher in the MCNP5 simulation of 60 Co source than the experimental result unlike other modeling of 137 Cs source where the experimental result is in a good agreement with the MCNP5 simulation. This discrepancy is regarded as coming from a loss of a signal due to the sum peak of 1173.2 keV and 1332.5 keV which is not treated in MCNP5 simulation The MCNP5 based simulation exhibits a good performance for generation of the response function. However, it does not guarantee that the desirable solution can be obtained when the response function is applied to the spectral decovolution This is because t he MCNP5 generated response function is an approximate value even though it is very similar to the true response function. Note that even a small discrepancy can lead to u nphysical results in the deconvolution process because the res ponse function matrix is ill posed. In the current study, two deconvolution algorithm s MLEM and DIM GMM reproduce d the source intensity within a low relative error of the true value in Tables 4 1 to 4 4 F ailure of the GMM to separate the overlapped peaks in the experimental ly obtained source spectrum is not due to the response function, but to the mechanism of the GMM. T he unexpected spurious peaks were also found in the deconvolved spectra ( Figures 4 7b, 4 8b and 4 11), which was deconvolved by the MLEM It is unclear whether those peaks are associated with inaccuracy of PAGE 91 91 the MCNP5 generated response function or the deconvolution methodology applied. Regarding that the DIM d id not produce the spurious peaks in the folded spectrum ( Figure 4 12) the spurious peaks seem to be related to the MLEM rather than the inaccuracy of the MCNP5 generated response function. There are some factors that the MCNP5 simulation did not treat in producing t he response function: the minor physical effects in the detector background radiation, and backscatter. The former hardly makes an impact on the spectral deconvolution. However, t he effect of backscatter becomes fairly obvious in the low energy region. Ba ckscatter will not have a direct effect on the spectral deconvolution, but could distort the intensity of the full energy peaks around th e backscattering energy region in the measurement As a result, it causes the poor recovery of the measured spectrum. T herefore, w ithout proper treatment of backscatter in the generation of the response function matrix it could be a very serious problem in deconvolving the measured s pectrum. T his adverse effect was already witnessed in the MELM and the GMM. Therefore, it is more significant in two methods to develop the response function algorithm which can deal with Compton backscatt er rather than to try to construct an accurate response function. Deconvolution M ethods The MLEM exhibits a very good decon volution ability f or the MCNP5 simulated spectra 1 and 2 which do not include backscatter and overlapped peaks. When the MLEM wa s performed on the MCNP5 simulated spectrum 3 which include the overlapped peaks, it was very effective in peak searching and peak separation, but it presented a relatively poor result for the recovery of the source intensity And also, s everal spurious peaks were observed around in the low energy region. The DI M GMM is combined of advantage s of two algorithms the DIM (D irect I nverse M ethod ) and t he GMM ( Gaussian M ixture M ethod ): t he DIM provides a simple and effective PAGE 92 92 way to recover the source intensity, and the GMM a resolving power to separate the overlapped peaks. The DIM GMM decomposes the response function matrix into two other matri c es one of which is well pos ed. The well posed matrix can be directly solved by the simple inverse method. Therefore, t he biggest advantage of the DIM is that its solution is not artificially distorted but physically understandable unlike a solution of the MLEM wh ich may cause the spurious peaks. The DIM eliminates all phenomena occurring in a detector except for the photoelectric absorption such as Compton scattering, x ray escape, and single and double escape peaks which attribute to the spec trum backgrounds. An d then, i t recover s information about photon s including backscatter before enter ing the detector. As described at section 4 3 2 ( Figure 4 11) the DIM can remove the backgrounds mentioned above from even a complicated spectrum which is convolved of backsca tter and overlapped spectra. Its drawback is that it is not able to resolve the energy broadened spectrum. In the DIM GMM, this problem is solved using the GMM which is a technique to separate the overlapped spectrum o n the assumption that the unfolded sp ectrum is convolved with Gaussian peaks This method is very effective to resolve a spectrum of a mix of Gaussian peaks. T he most significant issue of the GMM is how to improve the slow convergence due to the computational complexity, and how to determine th e initial guessing with regard to the local maximum. Especially, determination of the initial value may be a critical problem when there are many peaks to be separated For that reason, the initial value guessing should be carefully treated in deconvolut ion of a spectrum of a mixture of unknown sources. As discussed in Chapter 4, the DIM GMM gave a n excellent accuracy of the spectral deconvolution (the relative error of the true value is less than 1 %) when the spectrum which was processed by the DIM wa s a mix of Gaussian peaks, i.e., it does not contain backscatter PAGE 93 93 However, it had a difficulty dealing with the backscatter effect. Previously mentioned, the backscatter effect may contribute to the intensity of the photopeak in the measured spectrum. Such a phenomenon inevitably brings about the inaccuracy of the recovered peaks in the deconvol ution process Accordingly, for a successful deconvolution, the response function should be generated considering backscatter. PAGE 94 94 CHAPTER 6 CONCLUSION All the general procedures necessary for extracting information from the measured gamma spectrum of a BGO detector and point line radioisotope sources have been discussed with regard to the MCNP modeling and different deconvolution methods. The MCNP5 code was employed to study the detector sensitivity. The code embodies the typical characteristics of the BGO detector except for backscatter which w as not treated in this simulation. Even though we partially treated the exponential tail on the low energy side of the full energy pea k which is a very interesting topic in many literatures, the simulation results are in good agreement with the experimental result s The modeling of the BGO detector exhibits a reasonable spectrum variation corresponding to four phys ical changes of the detector i.e., crystal size, detector to source distance, window thickness and detector source angle The energy resolution was semi empirically calculated through the analytical function since the MCNP5 does not present the transport of visible photons which are the main cause of the energy broadening The artificially added energy broadening effect shows a good correspondence with the experimental result. The energy broadening function of Eq. 2 7 was used to implement the energy bro adening effect in the MCNP5 simulated spectrum. Coefficients of the energy broadening function was determined using three gamma energies, 81 keV, 661.2 keV, and 1274.5 keV coming from 81 keV 133 Ba, 137 Cs, and 22 Na respectively. However, f or the exact est imation of the coefficients, more gamma ray sources should be selected in the gamma spectrum analysis which spans the wide energy range because there is no guarantee that the coefficients will fit the entire region [7 7 ] Actually, in the energy region belo w the 50 keV a slight discrepancy was observed between the simulation result and the experimental result. PAGE 95 95 The detector response matrix was constructed which has a 1500 x 1500 dimension. Each matrix element corresponds to 1 keV. In fact, 1500 response funct ions are required for generating the 1500 x 1500 response function matrix. It would be an unimaginably tedious and laborious job. Thus, most of the response functions were produced using a linear numerical interpolation from 21 selected response functions This is a very effective method to build many response functions. The MLEM method exhibits a very good unfolding ability for the MCNP5 simulated spectra The energy resolution is remarkably improved from ~14 % to 2 %. In the experimentally measured spect r a, some artificial peaks were observed in the Compton backscattering region. The measured spectrum which wa s convolved with a variety of backgrounds reduces the detector resolution compared to the MCNP5 simulated spectrum. However, the overall performance of the MLEM is satisfactory. The new deconvolution algorithm, so called the DIM G MM, was developed, which combines advantage s of the DIM and the GMM The DIM GMM was carried out with a series of two methods. First, the DIM filtered out the Compton backgr ound effectively and recovered the ful l energy peak of the source. It did not, however, exhibit good performance in unfolding of the measured spectrum for closely spaced full energy peaks due to the low resolution of the BGO detector. If the DIM is applied to a detector of high resolution such as the HPGe detector, it could yield a better result Therefore, one more process is required for resolving the spectra obtained from the DIM The overlapping peaks were unfolded using the GMM after deconvolution of t he DIM. The deconvolving power of the DIM GMM is superior to the MELM method for the MCNP5 simulated sample spectrum. However, it failed to separate peaks from the PAGE 96 96 experimentally measured sample spectrum because backscatter w as not treated in constructing the response function. The reason that the DIM & GMM brought about the worse outcome for the experimentally measured sample spectrum than the MELM is that the GMM is based on the assumption that the sample spectrum consist of many Gaussian peaks : that is, backscatter is not a Gaussian form of a peak. The present study has demonstrated that the MCNP5 code is a very good tool to model characteristics of the BGO detector and to build the detector response function even though it did not model some detai led features of the BGO detector The small discrepancy of the MCNP5 calculated response function from the true one did not make a noticeable impact on the spectral deconvolution. On the contrary, t wo deconvolution methods show that treatment of backscatte r in the response function is more significant rather than accuracy of the calculated response function The need for treatment of backscatter in the response function does not arise from a problem of the deconvolution method, but from distortion on the me asured photopeak intensity due to backscatter. The deconvolution r esult has proved that the MLEM and the DIM & GMM are ve ry useful in deconvolution of the MCNP5 simulated spectrum which does not contain backscatter The DIM GMM gives a better performance Such a result is due to the po werful recovery ability of the source information of the DIM, and the effective peak separation ability of the GMM. Of course, the backscatter effect should be well treated for the desirable deconvolution of the experimental d ata because it seriously degrades performance of the spectral deconvolution. Then the DIM GMM will provide a strong tool for unfolding of the spectrum. In order to enhance the deconvolution ability of the DIM GMM in various applications, more in depth stu dies are required in the future with regards to its drawbacks already shown in PAGE 97 97 this research. First development of an algorithm should be carried out, which is able to treat backscatter in construction of the detector response function. Secondly, the peak separating ability o f the IDM GMM should be studied for more complicated situations in which the measured spectrum includes more gamma peaks or overlapped peaks Thirdly, an application method of the DIM GMM to multi sources or bulk source should be inves tigated because the gamma source is not a point source in real spectral applications such as NDA (Non Destructive Analysis). Lastly, it is worth applying the DIM to the measured gamma spectrum obtained by the HPGe detector. The current study has shown that the DIM successfully recovers information of the original source except for energy broadening which is created by the BGO detector. Therefore, if the DIM is applied to the HPGe detector measured spectrum whose energy resolution is very low, it could provi de a simple and powerful tool for spectral deconvolution without using complicated algorithms. PAGE 98 98 APPENDIX A MCNP CODE FOR THE BGO DETECTOR MODELING In this Appendix, we provide sample MCNP5 input files for the 137 Co placed at 5 cm away from a BGO detecto *********************** 60 Co ********************************************* MODEL: 3.0"(diameter) X 3.0"(length) Case 3 c c MODELING OF THE BGO DETECTOR c c CELL CARDS 1 1 7.13 4 5 11 imp:e=1 imp:p=1 $ CRYSTAL 2 0 (3 4 12):(4 5 11 12) imp:e=1 imp:p=1 $ Between CRYSTAL AND WINDOW 3 2 2.7 (2 3 13):(3 5 12 13) imp:e=1 imp:p=1 $ DETECTOR WINDOW 4 0 1 2 13 imp:e=1 imp:p=1 $ SRC CELL (ISOTROPIC DISK SRC) 5 0 1: 5: 13 imp:e=0 imp:p=0 $ REST OF THE UNIVERSE c SURFACE CARDS 1 px 15.1 $ PLANE DELIMITING THE SOURCE REGION 2 px 0.0 $ PLANE DELIMITING THE FRONT OF THE DETECTOR 3 px 0.05 $ PLANE DELIMITING THE INSIDE OF THE DETECTOR WINDOW 4 px 0.25 $ PLANE DELIMITING THE FRONT OF THE CRYSTAL 5 px 7.87 $ PLANE DELIMITING THE BACK OF THE CRYSTAL 11 cx 3.81 $ R ADIUS OF THE CRYSRAL (3.0" DIAMETER) 12 cx 4.01 $ RADIUS OF THE INSIDE OF THE DETECTOR WINDOW 13 cx 4.06 $ RADIUS OF THE DETECTOR c DATA CARDS c SIMPLE GEOMETRY MCNP5 CALCULATES THE VOLUM ES mode p e $ PHOTON/ELECTRON MODE c c SOURCE CARDS sdef erg=d1 pos 2.0 0 0 axs= 1 0 0 rad=d2 par 2 $ ISOTROPIC DISK SRC Co 60 si1 l 1.1732 1.3325 sp1 1.0 1.0 c si2 0 0.15 $ Co 60 source radius sp2 21 1 c c TALLY CARDS f4:p 1 fc4 Crystal total flux f8:p 1 fc8 Pulse height in crystal volume (cell 1) PAGE 99 99 e8 0.0 1.0E 05 1.0E 04 0.001 0 .002 0.003 0.004 0.005 0.006 0.007 0.008 0.009 0.01 0.011 0.012 0.013 0.014 0.015 0.016 0.017 0.018 0.019 0.02 0.021 0.022 0.023 0.024 0.025 0.026 0.027 0.028 0.029 0.03 0.031 0.032 0.033 0.034 0.035 0.036 0.037 0.038 0.039 0.04 0.041 0 .042 0.043 0.044 0.045 0.046 0.047 0.048 0.049 0.05 0.051 0.052 0.053 0.054 0.055 0.056 0.057 0.058 0.059 0.06 0.061 0.062 0.063 0.064 0.065 0.066 0.067 0.068 0.069 0.07 0.071 0.072 0.073 0.074 0.075 0.076 0.077 0.078 0.079 0.08 0.081 0 .082 0.083 0.084 0.085 0.086 0.087 0.088 0.089 0.09 0.091 0.092 0.093 0.094 0.095 0.096 0.097 0.098 0.099 0.1 c 0.101 0.102 0.103 0.104 0.105 0.106 0.107 0.108 0.109 0.11 0.111 0.112 0.113 0.114 0.115 0.116 0.117 0.118 0.119 0.12 0.121 0.122 0.123 0.124 0.125 0.126 0.127 0.128 0.129 0.13 0.131 0.132 0.133 0.134 0.135 0.136 0.137 0.138 0.139 0.14 0.141 0.142 0.143 0.144 0.145 0.146 0.147 0.148 0.149 0.15 0.151 0.152 0.153 0.154 0.155 0.156 0.157 0.158 0.159 0.16 0.161 0.162 0.163 0.164 0.165 0.166 0.167 0.168 0.169 0.17 0.171 0.172 0.173 0.174 0.175 0.176 0.177 0.178 0.179 0.18 0.181 0.182 0.183 0.184 0.185 0.186 0.187 0.188 0.189 0.19 0.191 0.192 0.193 0.194 0.195 0.196 0.197 0.198 0.199 0.2 c 0.201 0.202 0.203 0.204 0.205 0.206 0.207 0.208 0.209 0.21 0.211 0.212 0.213 0.214 0.215 0.216 0.217 0.218 0.219 0.22 0.221 0.222 0.223 0.224 0.225 0.226 0.227 0.228 0.229 0.23 0.231 0.232 0.233 0.234 0.235 0.236 0.237 0.238 0.239 0.24 0.241 0.242 0.243 0.244 0.245 0.246 0.247 0.248 0.249 0.25 0.251 0.252 0.253 0.254 0.255 0.256 0.257 0.258 0.259 0.26 0.261 0.262 0.263 0.264 0.265 0.266 0.267 0.268 0.269 0.27 0.271 0.272 0.273 0.274 0.275 0.276 0.277 0.278 0.279 0.28 0.281 0.282 0.283 0.284 0.285 0.286 0.287 0.288 0.289 0.29 0.291 0.292 0.293 0.294 0.295 0.296 0.297 0.298 0.299 0.3 c 0.301 0.302 0.303 0.304 0.305 0.306 0.307 0.308 0.309 0.31 0.311 0.312 0.313 0.314 0.315 0.316 0.317 0.318 0.319 0.32 0.32 1 0.322 0.323 0.324 0.325 0.326 0.327 0.328 0.329 0.33 0.331 0.332 0.333 0.334 0.335 0.336 0.337 0.338 0.339 0.34 0.341 0.342 0.343 0.344 0.345 0.346 0.347 0.348 0.349 0.35 0.351 0.352 0.353 0.354 0.355 0.356 0.357 0.358 0.359 0.36 0.36 1 0.362 0.363 0.364 0.365 0.366 0.367 0.368 0.369 0.37 0.371 0.372 0.373 0.374 0.375 0.376 0.377 0.378 0.379 0.38 0.381 0.382 0.383 0.384 0.385 0.386 0.387 0.388 0.389 0.39 0.391 0.392 0.393 0.394 0.395 0.396 0.397 0.398 0.399 0.4 c 0.4 01 0.402 0.403 0.404 0.405 0.406 0.407 0.408 0.409 0.41 PAGE 100 100 0.411 0.412 0.413 0.414 0.415 0.416 0.417 0.418 0.419 0.42 0.421 0.422 0.423 0.424 0.425 0.426 0.427 0.428 0.429 0.43 0.431 0.432 0.433 0.434 0.435 0.436 0.437 0.438 0.439 0.44 0.4 41 0.442 0.443 0.444 0.445 0.446 0.447 0.448 0.449 0.45 0.451 0.452 0.453 0.454 0.455 0.456 0.457 0.458 0.459 0.46 0.461 0.462 0.463 0.464 0.465 0.466 0.467 0.468 0.469 0.47 0.471 0.472 0.473 0.474 0.475 0.476 0.477 0.478 0.479 0.48 0.4 81 0.482 0.483 0.484 0.485 0.486 0.487 0.488 0.489 0.49 0.491 0.492 0.493 0.494 0.495 0.496 0.497 0.498 0.499 0.5 c 0.501 0.502 0.503 0.504 0.505 0.506 0.507 0.508 0.509 0.51 0.511 0.512 0.513 0.514 0.515 0.516 0.517 0.518 0.519 0.52 0. 521 0.522 0.523 0.524 0.525 0.526 0.527 0.528 0.529 0.53 0.531 0.532 0.533 0.534 0.535 0.536 0.537 0.538 0.539 0.54 0.541 0.542 0.543 0.544 0.545 0.546 0.547 0.548 0.549 0.55 0.551 0.552 0.553 0.554 0.555 0.556 0.557 0.558 0.559 0.56 0. 561 0.562 0.563 0.564 0.565 0.566 0.567 0.568 0.569 0.57 0.571 0.572 0.573 0.574 0.575 0.576 0.577 0.578 0.579 0.58 0.581 0.582 0.583 0.584 0.585 0.586 0.587 0.588 0.589 0.59 0.591 0.592 0.593 0.594 0.595 0.596 0.597 0.598 0.599 0.6 c 0 .601 0.602 0.603 0.604 0.605 0.606 0.607 0.608 0.609 0.61 0.611 0.612 0.613 0.614 0.615 0.616 0.617 0.618 0.619 0.62 0.621 0.622 0.623 0.624 0.625 0.626 0.627 0.628 0.629 0.63 0.631 0.632 0.633 0.634 0.635 0.636 0.637 0.638 0.639 0.64 0 .641 0.642 0.643 0.644 0.645 0.646 0.647 0.648 0.649 0.65 0.651 0.652 0.653 0.654 0.655 0.656 0.657 0.658 0.659 0.66 0.661 0.662 0.663 0.664 0.665 0.666 0.667 0.668 0.669 0.67 0.671 0.672 0.673 0.674 0.675 0.676 0.677 0.678 0.679 0.68 0 .681 0.682 0.683 0.684 0.685 0.686 0.687 0.688 0.689 0.69 0.691 0.692 0.693 0.694 0.695 0.696 0.697 0.698 0.699 0.7 c 0.701 0.702 0.703 0.704 0.705 0.706 0.707 0.708 0.709 0.71 0.711 0.712 0.713 0.714 0.715 0.716 0.717 0.718 0.719 0.72 0.721 0.722 0.723 0.724 0.725 0.726 0.727 0.728 0.729 0.73 0.731 0.732 0.733 0.734 0.735 0.736 0.737 0.738 0.739 0.74 0.741 0.742 0.743 0.744 0.745 0.746 0.747 0.748 0.749 0.75 0.751 0.752 0.753 0.754 0.755 0.756 0.757 0.758 0.759 0.76 0.761 0.762 0.763 0.764 0.765 0.766 0.767 0.768 0.769 0.77 0.771 0.772 0.773 0.774 0.775 0.776 0.777 0.778 0.779 0.78 0.781 0.782 0.783 0.784 0.785 0.786 0.787 0.788 0.789 0.79 0.791 0.792 0.793 0.794 0.795 0.796 0.797 0.798 0.799 0.8 c 0.801 0.802 0.803 0.804 0.805 0.806 0.807 0.808 0.809 0.81 0.811 0.812 0.813 0.814 0.815 0.816 0.817 0.818 0.819 0.82 0.821 0.822 0.823 0.824 0.825 0.826 0.827 0.828 0.829 0.83 PAGE 101 101 0.831 0.832 0.833 0.834 0.835 0.836 0.837 0.838 0.839 0.84 0.841 0.842 0.843 0.844 0.845 0.846 0.847 0.848 0.849 0.85 0.851 0.852 0.853 0.854 0.855 0.856 0.857 0.858 0.859 0.86 0.861 0.862 0.863 0.864 0.865 0.866 0.867 0.868 0.869 0.87 0.871 0.872 0.873 0.874 0.875 0.876 0.877 0.878 0.879 0.88 0.881 0.882 0.883 0.884 0.885 0.886 0.887 0.888 0.889 0.89 0.891 0.892 0.893 0.894 0.895 0.896 0.897 0.898 0.899 0.9 c 0.901 0.902 0.903 0.904 0.905 0.906 0.907 0.908 0.909 0.91 0.911 0.912 0.913 0.914 0.915 0.916 0.917 0.918 0.919 0.92 0.921 0.922 0.923 0.924 0.925 0.926 0.927 0.928 0.929 0.93 0.931 0.932 0.933 0.934 0.935 0.936 0.937 0.938 0.939 0.94 0.941 0.942 0.943 0.944 0.945 0.946 0.947 0.948 0.949 0.95 0.951 0.952 0.953 0.954 0.955 0.956 0.957 0.958 0.959 0.96 0.961 0.962 0.963 0.964 0.965 0.966 0.967 0.968 0.969 0.97 0.971 0.972 0.973 0.974 0.975 0.976 0.977 0.978 0.979 0.98 0.981 0.982 0.983 0.984 0.985 0.986 0.987 0.988 0.989 0.99 0.991 0.992 0.993 0.994 0.995 0.996 0.997 0.998 0.999 1.0 c 1.001 1.002 1.003 1.004 1.005 1.006 1.007 1.008 1.009 1.01 1.011 1.012 1.013 1.014 1.015 1.016 1.017 1.018 1.019 1.02 1.021 1.022 1.023 1.024 1.025 1.026 1.027 1.028 1.029 1.03 1.031 1.032 1.033 1.034 1.035 1.036 1.037 1.038 1.039 1.04 1.041 1.042 1.043 1.044 1.045 1.046 1.047 1.048 1.049 1.05 1.051 1.052 1.053 1.054 1.055 1.056 1.057 1.058 1.059 1.06 1.061 1.062 1.063 1.064 1.065 1.066 1.067 1.068 1.069 1.07 1.071 1.072 1.073 1.074 1.075 1.076 1.077 1.078 1.079 1.08 1.081 1.082 1.083 1.084 1.085 1.086 1.087 1.088 1.089 1.09 1.091 1.092 1.093 1.094 1.095 1.096 1.097 1.098 1.099 1.1 c 1.101 1.102 1.103 1.104 1.105 1.106 1.107 1.108 1.109 1.11 1.111 1.112 1.113 1.114 1.115 1.116 1.117 1.118 1.119 1.12 1.121 1.122 1.123 1.124 1.125 1.126 1.127 1.128 1.129 1.13 1.131 1.132 1.133 1.134 1.135 1.136 1.137 1.138 1.139 1.14 1.141 1.142 1.143 1.144 1.145 1.146 1.147 1.148 1.149 1.15 1.151 1.152 1.153 1.154 1.155 1.156 1.157 1.158 1.159 1.16 1.161 1.162 1.163 1.164 1.165 1.166 1.167 1.168 1.169 1.17 1.171 1.172 1.173 1.174 1.175 1.176 1.177 1.178 1.179 1.18 1.181 1.182 1.183 1.184 1.185 1.186 1.187 1.188 1.189 1.19 1.191 1.192 1.193 1.194 1.195 1.196 1.197 1.198 1.199 1.2 c 1.201 1.202 1.203 1.204 1.205 1.206 1.207 1.208 1.209 1.21 1.211 1.212 1.213 1.214 1.215 1.216 1.217 1.218 1.219 1.22 1.221 1.222 1.223 1.224 1.225 1.226 1.227 1.228 1.229 1.23 1.231 1.232 1.233 1.234 1.235 1.236 1.237 1.238 1.239 1.24 1.241 1.242 1.243 1.244 1.245 1.246 1.247 1.248 1.249 1.25 PAGE 102 102 1.251 1.252 1.253 1.254 1.255 1.256 1.257 1.258 1.259 1.26 1.261 1.262 1.263 1.264 1.265 1.266 1.267 1.268 1.269 1.27 1.271 1.272 1.273 1.274 1.275 1.276 1.277 1.278 1.279 1.28 1.281 1.282 1.283 1.284 1.285 1.286 1.287 1.288 1.289 1.29 1.291 1.292 1.293 1.294 1.295 1.296 1.297 1.298 1.299 1.3 c 1.301 1.302 1.303 1.304 1.305 1.306 1.307 1.308 1.309 1.31 1.311 1.312 1.313 1.314 1.315 1.316 1.317 1.318 1.319 1.32 1.321 1.322 1.323 1.324 1.325 1.326 1.327 1.328 1.329 1.33 1.331 1.332 1.333 1.334 1.335 1.336 1.337 1.338 1.339 1.34 1.341 1.342 1.343 1.344 1.345 1.346 1.347 1.348 1.349 1.35 1.351 1.352 1.353 1.354 1.355 1.356 1.357 1.358 1.359 1.36 1.361 1.362 1.363 1.364 1.365 1.366 1.367 1.368 1.369 1.37 1.371 1.372 1.373 1.374 1.375 1.376 1.377 1.378 1.379 1.38 1.381 1.382 1.383 1.384 1.385 1.386 1.387 1.388 1.389 1.39 1.391 1.392 1.393 1.394 1.395 1.396 1.397 1.398 1.399 1.4 c 1.401 1.402 1.403 1.404 1.405 1.406 1.407 1.408 1.409 1.41 1.411 1.412 1.413 1.414 1.415 1.416 1.417 1.418 1.419 1.42 1.421 1.422 1.423 1.424 1.425 1.426 1.427 1.428 1.429 1.43 1.431 1.432 1.433 1.434 1.435 1.436 1.437 1.438 1.439 1.4 4 1.441 1.442 1.443 1.444 1.445 1.446 1.447 1.448 1.449 1.45 1.451 1.452 1.453 1.454 1.455 1.456 1.457 1.458 1.459 1.46 1.461 1.462 1.463 1.464 1.465 1.466 1.467 1.468 1.469 1.47 1.471 1.472 1.473 1.474 1.475 1.476 1.477 1.478 1.479 1.4 8 1.481 1.482 1.483 1.484 1.485 1.486 1.487 1.488 1.489 1.49 1.491 1.492 1.493 1.494 1.495 1.496 1.497 1.498 1.499 1.5 c c MATERIAL CARDS c MOST RECENT PHOTOATOMIC LIBRARY MCLIB04 IS USED c M1 = B GO CRYSTAL: Bi4 Ge3 O12 c M2 = ALUMINUM HOUSING (DETECTOR WINDOW): Al m1 83000.04p 4 32000.04p 3 8000.04p 12 m2 13000.04p 1 c ENERGY CARDS phys:p $ DEFAULT VALUES phys:e $ D EFAULT VALUES c CUTOFF CARDS cut:p $ DEFAULT VALUE cut:e $ DEFAULT VALUE nps 100000000 c PERIPHERAL CARDS print prdmp 2j 1 1 5000000 PAGE 103 103 APPENDIX B C LANGUAGE CODE FOR ENERGY BROA DENING EFFECT In this Appendix, we provide a C code for the energy broadening. #include PAGE 104 104 { peak[energy]=1.0; } } for (j=1 ; j <= (SourceE[i]) ; j++) { FWHM = a + b*sqrt(Eo[j]+c*pow(Eo[j],2)); s = FWHM/2.35; // standard deviation E = Eo[j]; gaus = 1; while (gaus > 0.0001) { gaus = exp( pow(E Eo[j],2)/(2*pow(s,2))); E = E + 1; Eho = E; Eh = E; } x=0; sum = 0.0; gaus_n = 0.0; Amp = 0.0; if (E > 1500) Eh = 1500; for ( x = (2*Eo[j] E) ; x <= E; x=x+1) { Amp = peak[j]/sum; s_peak =0.0; for ( x = (2*Eo[j] Eho); x <= Eh; x=x+1) { if ( x >= 0) { intx = x; P[intx] = P[intx] + Amp*exp( pow(x Eo[j],2)/(2*pow(s,2))); s_peak = s_peak + Amp*exp( pow(x Eo[j],2)/(2*pow(s,2))); } } } } fclose(out); } PAGE 105 105 APPENDIX C Linear interpolation of the response function The first procedure for generating the response function determines energy intervals whose response functions are constructed using the MCNP5 simulation. In the present study, 100 keV was used as the energy interval in the range of 200 keV to 1500 keV. Then, a linear interpolation approach is used to obtain the response function for the energies in the selected energies. Suppose that the response function of the 662 keV is obtained from the linear interpolation. Two response functions of 600 keV and 700 keV are selected as shown in Figures C 1 a and C 2 a The response function of the low energy (600 keV) is scaled up and the response functi on of the high energy (700 keV) is scaled down. The reason for interpolating from two response functions is to estimate the average response function so as to minimize the error. For more exact interpolation, the energy scale is divided into two regions co rresponding to physical ch aracteristics of the response: R egion 1 (C ompton scattering region), and R egion 2 (full photopeak region) as shown in Figures C 1 and C 2 Let us start the interpolation with expanding the response function of 600 keV to that of the 662 keV. The scaling of the response f unction is done separately for R egions 1 and 2. The relation of the interpolated response function and the response function of the low energy is given by 1 = 1 1 1 0 < < 1 C 1 ) where 1 indicates the divis ion point of the low energy in R egion 1, and 1 denotes the division point of the interpolated energy in this region. 1 is the interpolated response function o btained from the low energy in R egion 1, and 1 1 1 is the response function of the low energy (600 keV) in R egion 1 ( Figure C 1 ). PAGE 106 106 The scali ng of the response function in R egion 2 is carried out followi ng the same procedure as in R egion 1. 2 = 2 1 + 1 2 1 2 1 1 < < 2 C 2 ) Then, we compress the response funct ion of the high energy (700 keV). The procedure of the scale down of the response function is the same as that of the scale up. Therefore, the response function of scale up can be expressed as 1 = 1 1 1 0 < < 1 C 3 ) 2 = 2 1 + 1 2 1 2 1 1 < < 2 C 4 ) Two response functions obtained from the interpolation is averaged over each reg ion considering weighting factor. 1 = 100 1 + 1 100 C 5 ) 2 = 100 2 + 2 100 C 6 ) where weighting factor is given as k =662 600 = 62. PAGE 107 107 Figure C 1. Linear expansion of the response function of the low energy Figure C 2. Linear compression of the resp onse function of the low energy PAGE 108 108 APPENDIX D MATLAB CODE FOR INTERPOLATION OF THE RESPONSE F UNCTION MATRIX In this Appendix, we provide a C code for the interpolation of the detector response matrix. #include PAGE 109 109 printf (" Type the 1st Source energy:"); scanf("%d", &peakL0); printf (" Type the 2nd Source energy:"); scanf("%d", &peakH0); if ((in1 = fopen(Open1, "rt")) == NULL) { fprintf(stderr, \ nCannot open input file. \ n"); exit(0); } if ((in2 = fopen(Open2, "rt")) == NULL) { fprintf(stderr, \ nCannot open output file. \ n"); exit(0); } for (i=0; i < 101; i++) { if ((out[i] = fopen(Oput[i], "wt")) == NULL) { fprintf(stderr, \ nCa nnot open output file. \ n"); exit(0); } } i = 0; while(!feof(in1)) // (!feof(in)) { fscanf(in1," %f %g %f \ n", &DL0, &DL1, &DL2); if (i >= 3 ) { gL[i 2] = DL1; } i = i+1; } //while i=0; while(!feof(in2)) { fscanf(in2," %f %g %f \ n", &DH0, &DH1, &DH2); if (i>=3) { PAGE 110 110 gH[i 2] = DH1; } i = i+1; } zeroH = 1 + peakH0; peakL1 = pea kL0 75; peakL2 = peakL0 77; peakL3 = peakL0 87; peakL4 = peakL0 90; peakL5 = (int) (2.0*peakL0*peakL0/(moC2 + 2*peakL0) + 0.51); peakL6 = peakL0 511; peakL7 = peakL0 1022; peakH1 = peakH0 75; peakH2 = peakH0 77; peakH3 = peakH0 87; peakH4 = peakH0 90; peakH5 = (int) (2.0*peakH0*peakH0/(moC2 + 2*peakH0) + 0.51); peakH6 = peakH0 511; peakH7 = peakH0 1022; AL = peakL0 92; BL = peakL0 ; AH = peakH0 92; BH = peakH0; /********** ******************************************************************/ for(k=0;k<101;k++) { Ipeak0 = peakL0 + k; Ipeak5 = (int) (2.0*Ipeak0*Ipeak0/(moC2 + 2*Ipeak0) + 0.51); AI = Ipeak0 92; BI = Ipeak0; for(chn=1; chn <= AI; chn++) { EL[k][chn] = gL[(int)(chn*AL/AI +0.5)]; } for(chn= AI; chn <= BI; chn++) { EL[k][chn] = gL[AL + (int)((chn AI)*(peakL0 AL)/(Ipeak0 AI)+ 0.5)]; } PAGE 111 111 for(chn=1; chn <= AI; chn++) { EH[k][chn] = gH[(int)(chn*AH/AI + 0.5)]; } for(chn=AI; chn <= BI; chn++) { EH[k][chn] = gH[AH + (int)((chn AI)*(peakH0 AH)/(Ipeak0 AI) +0.5) ]; } for(chn=1; chn <= BI; chn++) { E[k ][chn] = ((100 k)*EL[k][chn] + k*EH[k][chn])/100; } } for (i=0; i<101; i++) { for(chn=1; chn <= BI; chn++) { fprintf( out[i], "%d %g \ n", chn, E[i][chn]); } } fclose(in1); fclose(in2); for (i=0; i<101; i++) { fclose(out[i]); } } PAGE 112 112 APPENDIX E MATLAB CODE FOR THE SPECTRAL DECONVO LUTION USING THE MLE M ALGORITHM In this Appendix, we provide a Matlab code for the deconvolution of the measured spectrum using t he MLEM algorithm. clear all; in1 = fopen('Rmatrix R00','r'); in2 = fopen('RMn.txt','r'); in3 = fopen('Initial10 7.txt','r'); out1 = fopen('ML MP YB Mn result.txt','w'); e=1.0; A1=fscanf(in1,'%g %g g %g \ n', [1500 inf]); A2=fscanf(in2,'%d %g \ n', [ 2 inf]); A3=fscanf(in3,'%d %g \ n', [2 inf]); MS =zeros(1500,2); initS =zeros(1500,1); convergeS =zeros(1,1500); NewS = zeros(1500,1); P = zeros(1500,1); iNo = zeros(1500,1); e = zeros(1500,1); correction =zeros(1500,1); for (i=1:1500) MS(i) = A2(2,i) ; initS(i) = A3(2,i); end for b=1:1500 for d=1:1500 P(b) = P(b) + A1(d,b); end end ee=1.0; eiNo=50; while (ee > 0.0001 & eiNo < 500) correction = zeros(1500,1); PAGE 113 113 for (d=10:1500) dmean =0.0; for(b=10:1500) dmean = dmean + initS(b)*A1(d,b); end if (dmean ~= 0.0) for (k=10:1500) correction(k) = correction(k) + MS(d)*A1(d,k )/dmean; end end end for (k=10:1500) NewS(k) =(initS(k)/P(k))*correction(k); if k==1170 dd = NewS(k,1); end if initS(k) > 0.0 diff = NewS(k) initS(k); if diff > 0 e(k) = diff/initS(k); else e(k) = 0.0; else e(k) = 0.0; end i nitS(k) = NewS(k); kk = initS(k); iNo(k)= iNo(k)+1; end eiNo = max(iNo(1:1500)); ee = max(e(1:1500)); end for (i=1:1500) convergeS(1,i)= initS(i); end fprintf(out1,'%g \ n', convergeS); fclose(in1); fclo se(in2); fclose(in3); fclose(out1); PAGE 114 114 APPENDIX F C LANGUAGE CODE FOR THE SPECTRAL DECONVO LUTION USING THE GMM In this Appendix, we provide a C code for the peak separation using the GMM. #include PAGE 115 115 printf(" \ n \ n Output file name2 :"); scanf("%s",oput1); printf(" \ n \ n Put the number of peaks :"); scanf("%d", &NoP); K=NoP; if ((out = fopen(Oput, "w t")) == NULL) { fprintf(stderr, \ nCannot open output file. \ n"); exit(0); } if ((out1 = fopen(oput1, "wt")) == NULL) { fprintf(stderr, \ nCannot open output file. \ n"); exit(0); } while(!feof(in)){ i=i+1; fscanf(in,"%f \ n", &measured_x); Sample_X[i]= measured_x; } h=0.0; for (i=1; i <= N; i++){ if (Sample_X[i] < 1e 7 ){ h = h + 1.0; Sample_X[i] = .0; } else { Sample_I nt = (int) 1000000*Sample_X[i]; Sample_X[i] = (float) Sample_Int; TSample_No = TSample_No + Sample_X[i]; } } for (m=1; m<=1500; m=m+1){ FWHM = c1 + c2*sqrt(m + c3*pow(m,2)); R[m] = (1/2.35)*FWHM/m; } PAGE 116 116 while (flag ==1){ for (i=1; i<=K; i=i+1){ a[i] = (float) 1/K; mu[1] = 20; mu[2]= 120; mu[3]= 450; mu[4]=550; mu[5]=900; mu[6]=1350; FWHM = c1 + c2*sqrt(mu[i]+ c3*pow(mu[i],2)); Sig[i] = pow(FWHM/2.35,2); } m=0; Tolerance1 = 1.0; while (m<200){ //&& Tolerance1 > 0.000001){ m = m+1; if ( m != 1){ for (j=1; j<= K; j=j+1){ a[j] = a_new[j]; } } for (i=1; i<=K; i=i+1){ a_new[i] = 0.0; for (x=1; x<=N; x=x+1){ g[x] = Gaussian( (float) x, mu[i], Sig[i], flag0); sumup = 0.0; for (j=1; j<=K; j=j+1){ sumup = sumup + a[j]*Gaussian((float) x, mu[j], Sig[j], flag0); } if (sumup != 0.0){ a_ik = a[i]*g[x]/sumup; if (Sample_X[x] != 0.0) a_new[i] = a_new[i] + Sample_X[x]*a_ ik; else a_new[i] = a_new[i] + a_ik; } else { a_new[i] = a_new[i] + 1; } } a_new[i]= (1/(TSample_No))*a_new[i]; e1[i] = abs(a_new[i] a[i])/a[i]; } PAGE 117 117 if ( m != 1){ for (j=1; j<=K; j=j+1){ mu[j] = mu_new[j]; } } for (i=1; i<=K; i=i+1){ mu_p1=0.0; mu_p2 =0.0; mu_p3 =0.0; for (x=1; x<=N; x=x+1){ g[x] = Gaussian( (float) x, mu[i], Sig[i], flag0); sumup =0.0; for (j=1; j<=K; j=j+1){ sumup = sumup + a[i]*Gaussian( (float) x,mu[j], Sig [j], flag0); } if (sumup != 0.0) Px_mu = a[i]*g[x]/sumup; else Px_mu = 0.0; mu_p1 = mu_p1 + Sample_X[x]*Px_mu; mu_p2 = mu_p2 + Sample_X[x]*x*Px_mu; mu_p3 = mu_p3 + Sample_X[x]*pow(x,2)*Px_mu; } FWHM = c1 + c2*sqrt(mu[i] + c3*pow(mu[i],2)); R0 = (1/2.35)*FWHM/mu[i]; mu_p5 = pow(R0,2); mu_p4 = pow(mu_p2,2) + 4*mu_p1*mu_p3*mu_p5; mu_new[i] = ( mu_p2 + sqrt(mu_p4))/(2*mu_p1*mu_p5); e2[i] = abs(mu_new[i] mu[i])/mu[i]; } for (i=1; i<=K; i=i+1){ FWHM = c1 + c2*sqrt(mu_new[i] + c3*pow(mu_new[i],2)); sig_p1 = FWHM/2.35; Sig_new[i] = pow(sig_p1,2); e3[i] = abs(Sig_new[i] Sig[i])/Sig[i]; Sig[i] = Sig_new[i]; } PAGE 118 118 Tolerance1 = e1[1]; for (i=1; i<=K; i=i+1){ t = Tolerance1 e1[i]; if (t < 0.0) Tolerance1 = e1[i]; t = Tolerance1 e2[i]; if (t < 0.0) Tolerance1 = e2[i]; t = Tolerance1 e3[i]; if (t < 0.0) Tolerance1 = e3[i]; } } mixing = a_new[1]+a_new[2]+a_new[3]+a_new[4]+a_new[5]; for (x=1; x<=N; x=x+1){ for (i=1; i<=K; i=i+1){ Calculated_X[x]= Calculated_X[x] + a[i]*Gaussian( (float) x, mu_new[ i], Sig_new[i], flag0); } peak[1][x] = a[1]*Gaussian((float) x, mu_new[1], Sig_new[1], flag0); peak[2][x] = a[2]*Gaussian((float) x, mu_new[2], Sig_ new[2], flag0); peak[3][x] = a[3]*Gaussian((float) x, mu_new[3], Sig_new[3], flag0); peak[4][x] = a[4]*Gaussian((float) x, mu_new[4], Sig_new[4], flag0); peak[5][x] = a[5]*Gaussian((float) x, mu_new[5], Sig_new[5], flag0); peak[6][x] = a[6]*Gaussian((float) x, mu_new[6], Sig_new[6], flag0); Error[x] = Sample_X[x] Calculated_X[x]; } flag = 0; } for (j=1; j<=1500; j=j+1){ fprintf(out,"%d %f \ n", j, Calculated_X[j]); fprintf(out1,"%d %f %f %f %f %f %f \ n", j, peak[1][j], peak[2][j], peak[3][j], peak[4][j], peak[5][j], peak[6][j]); } fclose(in); fclose(out); fclose(out1); } PAGE 119 119 LIST OF REFERENCES [ 1 ] G. F. Knoll Radiation detection and measurement, John Wiley & Sons, New York 2000 [ 2 ] L. Bouchet Astron. Astrophys. Suppl. Ser. 113 (1995) 167. [ 3 ] V. B. Anikeev, Nucl. Instr. and Meth. A 502 (2003) 792. [ 4 ] W. M. Smith, R. C. Barr, Journal of cardiovascular electrophysiology, 12 (2001) 253. [ 5 ] D. F. Covell, Anal. Chem. 31 (1959) 1785. [ 6 ] A. Likar, T. Vidmar, J. Phys. D, Appl. Phys. 36 (2003) 1903. [ 7 ] L. Loska, Appl. Radiat. Isot. 39 (1988) 475. [ 8 ] A. Mohammad Dj afari, Proceeding of the Maximum Entropy Conference (1996) 77. [ 9 ] A. Tarantola, Inverse problem theory and methods for model parameter estimation, SIAM, Paris, 2005. [ 10 ] Cs. Sukosd, W. Galster, I. Licot, M. P. Simonart, Nucl. Instr. and Meth. A 355 (1995) 552. [ 11 ] C M. Shyu, R. P. Gardner, K. Verghese, Nucl. Geophys. 7 (1993) 241. [ 12 ] O. Helene, C. Takiya, V. R. Vanin, Braz. J. Phys. 31 (2001) 8. [ 13 ] J. Idier, A. Mohammd Djafari, G. Demoment, 2 nd Intern. Conf. on Inverse Problems in Engineering, 1996. [ 14 ] M. Bertero, P. Boccac ci, Introduction to inverse problems in imaging IOP, London, 1998. [ 15 ] A. Dempster, N. Laird, D. Rubin, Journal of the Royal Statistics, Series B, 39 (1997) 1. [ 16 ] L. A. Shepp, Y. Vardi, IEEE Trans. Med. Imag. MI 1 (1982) 113. [ 17 ] H. Shi, B. Chen, T. Li, D. Yun, App l. Radiat. Isot. 57 (2002) 517. [ 18 ] P. H. G. M. Hendriks, M. Maucec, R. J. de Meijer, Appl. Radiat. Isot. 57 (2002) 449. [ 19 ] H. Vincke, E. Gschwendtner, C. W. Fabjan, T. otto, Nucl. Instr. and Meth. A 484 (2002) 102. [ 20 ] A. Sood, Ph. D. Thesis, NCSU, 2000. [ 21 ] P. Sangsing keow, K. D. Berry, E. J. Dumas, T. W. Raudorf, T. A. Underwood, Nucl. Instr. and Meth. A 505 (2003) 183. PAGE 120 120 [ 22 ] N. Inadama, H. Murayama, T. Yamaya, K. Kitamura, T. Yamashita, H. Kawai, T. Tsuda, M. Sato, Y. Ono, M. Hamamoto, IEEE Trans. Nucl. Sci. NS 53 (2006) 30 [ 23 ] coolant pump performance monitoring in IRIS using nitrogen th International Topical Meeting on Nuclear Plant Instrumentation, Control and Human Mac hine Interface Technology, Columbus, Ohio, 2004. [ 24 ] L. J. Meng, D. Ramsden, IEEE Trans. Nucl. Sci. NS 47 (4) (2000) 1011. [ 25 ] M. Moszynski, Nucl. Instr. and Meth. A 505 (2003) 101. [ 26 ] W. Guo, Ph. D. Thesis, NCSU, 2003. [ 27 ] M. Blaauw, V. O. Fernandez, W. Westmeier, Nucl. Instr. and Meth. A 387 (1997) 410. [ 28 ] W. Guo, R. P. Gardner, A. C. Todd, Nucl. Instr. and Meth. A 516 (2004) 586. [ 29 ] P. Quittner, Gamma ray spectroscopy, Adam Hilger, London, 1972. [ 30 ] G. H. Olivera, D. M. Shepard, P. J. Reckwerdt, K. Ruchala, J. Zachman, E. E. Fit chard, T. R. Mackie., Phys. Med. Biol. 43 (1998) 3277. [ 31 ] K. Chuang, M. Jan, J. Way, J. Lu, S. Chen, C. Hsu, Y. Fu, Comput Med Imaging and Graph, (2005) 571. [ 32 ] G. J. McLachlan, D. Peel, Finite mixture models, John Wiley & Sons, New York, 2000. [ 33 ] R. R. Heath, Sci ntillation spectrometry, Idaho National Engineering & Environmental Laboratory, 1964. [ 34 ] R. L. Heath, R.G. Helmer, L.A. Schmittroth, G.A. Cazier, NIM 47 (1967) 281. [ 35 ] M. J. Berger, S. M. Seltzer, NIM 104 (1972) 317. [ 36 ] M. C. Lee, K. Verghese, R.P. Gardner, Nucl. I nstr. and Meth. A 262 (1987) 430. [ 37 ] J. C. Vitorelli, A.X. Silva, V.R. Crispim, E.S. da Fonseca, W.W. Pereira, Appl. Radiat. Isot. 62 (2005) 619. [ 38 ] I. Orion, L. Wielopolski, In Vivo Body Composition Symposium, Ann. New York Acad. Sci. 904 (2000) 271. [ 39 ] X 5 MONTE A General Monte Carlo N Particle Transport LA CP 03 0245, 2003 [ 40 ] S. Agostinelli, et al., Nucl. Instr. Meth., A 506 (2003) 250. PAGE 121 121 [ 41 ] W. R. Leo, Techniques for nucl ear and particle physics experiments, Springer Verlog, 1994 [ 42 ] P. M. Marinkovi, V. Lj. Ljubenov, Jpn. J. Appl. Phys. 41 (2002) 5757. [ 43 ] M. Moszynski, M. Balcerzyk, W. Czarnacki, M. Kapusta, W. Klamra, A. Syntfeld, and M. Szawlowski, 2004, IEEE Trans. Nucl. Sci. NS 51 (3) (2004) 1074. [ 44 ] R. Venkataraman, S. Croft, W. R. Russ, J. Radioanal. Nucl. Chem. 264 (2005) 183. [ 45 ] A. Cesana, M. Terrani, Nucl. Instr. and Meth. A 281 (1989) 172. [ 46 ] J. Lassater, ORTEC Technical Support, Private communication, 2006. [ 47 ] E. P. Nassens, X. G Xu, Health Physics 77 (1) (1999) 76. [ 48 ] A. Sood, P. G. Robin Nucl. Instr. and Meth. B 213 (2004) 100. [ 49 ] Z. Wang, B. Kahn, J. D. Valentine, IEEE Trans. Nucl. Sci. NS 49 (9) (2002) 1925. [ 50 ] R. R. Kiziah, J. R. Lowell, NIM 305 (1991) 129. [ 51 ] M. Jandel, M. Morhac, J. Kliman, L. Krupa, V. Matousek, J. H. Hamilton, A. V. Ramayya, Nucl. Instr. and Meth. A 516 (2004) 172. [ 52 ] M. Bertero, T. A. Poggio, V. Torre, Proceeding of the IEEE 76 (8) (1988) 869. [ 53 ] A. Mohammad Djafari, Proc. SPIE 3459 (1998) 2. [ 54 ] N. S. Evans and P. B. Stark, 2002, Inverse Problems 18 (4) (2002) R55. [ 55 ] J. Aldrich, Statistical Science 12 (3) (1997) 162. [ 56 ] M. Garcia Talavera, B. Ulicny, Nucl. Instr. and Meth. A 286 (2003) 585. [ 57 ] G. J. McLachlan and T. Krishnan, The EM algorithm and extensions, John Wiley & Sons, New Jersey, 2008. [ 58 ] I. J. Myung, J. Mathematical Psychology 47 (2003) 90. [ 59 ] R. A. Fisher, Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character 222 (1922) 309. [ 60 ] E. Pardo Iguzquiza, P. A. Dow d, Computer. & Geosciences 23 (7) (1997) 793. [ 61 ] J.R. Moon, Target tracking 2004: Algorithms and Applications, IEE, (2004) 93. [ 62 ] L.R. Bahl, F. Jelinek, R.L. Mercer, IEEE Transactions on Pattern Analysis and Machine Intelligence, PAMI 5 (2) (1983) 179. PAGE 122 122 [ 63 ] Y. H. Chu ng, T. Y. Song, J. H. Jung, G. Cho, Y. S. Choe, K. Lee, S. E. Kim, B. Kim, 2004, IEEE Trans. Nucl. Sci. NS 51 (1) (2004) 101. [ 64 ] J. E. Bowsher, C. E. Floryd, J. Nucl. Med. 32 (6) (1991) 1285. [ 65 ] J. C. Engdahl and K. Bharwani, IEEE Nuclear Science Symposium Confe rence Record, (2005) 1346. [ 66 ] F. Leisch, Journal of Statistical Software 11 (8) (2004) 1 [ 67 ] M. Aitkin, Biometrics 55 (1999) 117. [ 68 ] M. Wedel, W. A. Kamakura, Market segmentation: conceptual and methodological foundations, Kluwer Academic Publishers, 2000. [ 69 ] J. Wolfowi tz, Ann. Math. Statist. 28 (1) (1957) 75. [ 70 ] S. Dasgupta, L. Schulman, J. of Machine Learning Research 8 (2007) 203. [ 71 ] P. Arcidiacono, J. B. Jones, Econometrica 71 (3) (2003) 933. [ 72 ] S. J. Roberts, D. Husmeier, I. Rezek, and W. Penny, IEEE Trans. Pattern Analysis and Machine Intelligence 20 (11) (1998) 1133. [ 73 ] P. Bandzuch, M. Morhac, J. Kristiak, Nucl. Instr. and Meth. A 384 (1997) 506. [ 74 ] M. Morhac, J. Kliman, V. Matousek, M. Veselsky, I. Turzo, Nucl. Instr. and Meth. A 401 (1997) 385. [ 75 ] P. Dorenbos, J. T. M. de Haas, C. W. van Eijk, IEEE Trans. Nucl. Sci. NS 42 (6) (1995) 2190. [ 76 ] C. W. Van Eijk, P. Dorenbos, E. V. D. van Loef, K. Kramer, H. U. Gudel, Radiation Measurement 33 (2001) 521. [ 77 ] H. M. Hakimabad, H. Panjeh, A. Vejdani Noghreiyan, Asian J. Exp. Sci. 21 (2) (2007) 23 3. PAGE 123 12 3 BIOGRAPHICAL SKETCH Jangyong Huh was born in 1967 in Seoul, South Korea to Jongwook Huh and Chundo L ee and is th eir only son They did not interfere with his life unless he broke the fundamental rules and principles which are required for the human society So, he grew up in a free atmosphere and as a result had a very open mind to the world. He began his academic career at Chung Ang University in Seoul in 1986 where he studied physics. After finishing with a Bachelor of Physics degree in 1992, he c ontinued his master degree for nuclear physics in the same university in 1994. Finally, he completed his doctorate in 2002. Even though he was offered a decent job after his doctorate, he decided to further his education at the University of Florida in pu rsuit of his long cherished desire to study and work in America. Such a hardly understandable decision might be due to his liberal spirit which has been built up through his life. For now, he has completed a Master of Nuclear Science degree in University o f Florida in 2008 and then has been waiting for his another journey in America. 