- Email: [email protected]
Principal Component Analysis for pulse-shape discrimination of scintillation radiation detectors
Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Q1 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
Contents lists available at ScienceDirect
Nuclear Instruments and Methods in Physics Research A journal homepage: www.elsevier.com/locate/nima
Principal Component Analysis for pulse-shape discrimination of scintillation radiation detectors T. Alharbi Department of Physics, Faculty of Science in Zulfi, Almajmaah University, P.O. Box 1712, 11932, Saudi Arabia
art ic l e i nf o
a b s t r a c t
Article history: Received 21 July 2015 Received in revised form 11 October 2015 Accepted 11 October 2015
In this paper, we report on the application of Principal Component analysis (PCA) for pulse-shape discrimination (PSD) of scintillation radiation detectors. The details of the method are described and the performance of the method is experimentally examined by discriminating between neutrons and gamma-rays with a liquid scintillation detector in a mixed radiation field. The performance of the method is also compared against that of the conventional charge-comparison method, demonstrating the superior performance of the method particularly at low light output range. PCA analysis has the important advantage of automatic extraction of the pulse-shape characteristics which makes the PSD method directly applicable to various scintillation detectors without the need for the adjustment of a PSD parameter. & 2015 Published by Elsevier B.V.
Keywords: Principal component analysis Pulse-shape discrimination Scintillation detectors
1. Introduction In several scintillation detectors sensing differences in the shape of output pulses provide very useful information on the type of radiation interacting with the detector. This is commonly referred to as pulse-shape discrimination (PSD). The most common PSD application is in the discrimination between neutrons and gamma-rays with organic scintillation detectors used as fast neutron detector [1]. The PSD techniques are also widely used in applications involving phoswich scintillation detectors [2]. The traditional PSD methods are the Charge-Comparison (CC) and the Zero-crossing (ZC) methods which have been implemented on analog circuitry for several decades [3,4]. The CC method is based on the comparison of integrals of the pulses for two different time intervals while in the ZC method a zero-crossing time bears the pulse-shape information. In recent years, there has been tremendous interest in using digital pulse-processing methods for PSD of scintillation detectors and a number of techniques such as digital versions of CC and ZC methods have been used for this purpose [5,6]. Moreover, pure digital techniques based on curve fitting [7], frequency domain analysis [8] and neural networks [9] have been used for PSD applications. Furthermore, digital PSD systems with the capability of real-time operation have been developed [10]. However, the available PSD methods are either computationally time intensive or require a PSD parameter to be precisely adjusted in accordance to the pulse-shape characteristics of the detector concerned. In this paper, a new digital PSD method based on the Principal Component Analysis (PCA) of sampled
detector pulses is proposed. PCA is a widely used technique for extracting strong patterns in a dataset and has found applications in various disciplines such as image processing, neuroscience, etc. [11]. Our motivation for employing PCA for PSD applications is that PCA analysis is inherently adaptive to the intrinsic statistical characteristics of the data, and therefore, a PCA-based PSD system would be directly applicable to various scintillation detectors without the need for adjusting a PSD parameter. We describe the details of the method and the performance of the method is experimentally examined by discriminating between neutrons and gamma-rays with a liquid scintillation detector in a mixed radiation field. The performance of the PCA pulse classification algorithm is also compared against that of the conventional CC method.
2. The PCA method In most of the scintillation detectors, a PSD technique makes use of a simple basic fact that the decay-times of a scintillation light pulse are dependent on the type of radiation interacting with the detector. The shape of the pulses is, however, subject to variations due to statistical fluctuations in the number of photoelectrons, time spread in the photomultiplier tube (PMT) and electronic noise [12]. Our PSD method is based on the employment of a PCA analysis for extracting the patterns in the decay-times of sampled PMT pulses, thereby different types of radiation can be identified. A comprehensive review of PCA can be found in many
http://dx.doi.org/10.1016/j.nima.2015.10.030 0168-9002/& 2015 Published by Elsevier B.V.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
2
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
references (e.g. Ref. [11]). Here, we briefly describe the PCA procedure in regard to our PSD application. In computational terms, the PSD method is based on the extraction of the eigenvectors and eigenvalues of a data covariance matrix formed by a set of pulses acquired during a pre-acquisition run. Consider during a preacquisition run a set of n pulses are acquired. Each pulse is composed of m samples covering the whole duration of the pulse. The pulses are then normalized to their amplitude and a data matrix X is formed in such a way that each pulse represents a row of the matrix. The matrix Xn m is described as: 2 3 x1;1 x1;2 : : x1;m 6 : : : : : 7 6 7 6 7 : : : : : 7 ð1Þ X nm ¼ 6 6 7 6 : : : : : 7 4 5 : : : xn;m xn;1
used. The detector container is made of aluminium and the PMT is directly mounted on the back circular surface of the detector. The detector has a 0.5 mm thick aluminium window. The PMT was operated with a negative voltage of 1600 V and the pulses from the anode of the PMT were directly digitized by using a fast waveform digitizer at 4 GS/s sampling rate and with 10-bit resolution (model DC252HF from Agilient Technologies Inc). Tests were carried out using an americium–beryllium (Am–Be) neutron source. A lead brick was placed in front of the detector to increase the ratio of the number of recorded neutrons to gamma-rays. Some tests with standard laboratory sources were also performed to characterise the light output performance of the system. The sampled PMT pulses were analyzed with algorithms written in MATLAB programming environment.
The j-th column of this matrix represents samples of random variable xj which are the normalized amplitudes of the pulses at a specific time during their lifetime. The covariance matrix of the random variables xj will be a matrix Cm m as: 2 3 C 1;1 C 1;2 : : C 1;m 6 : : : : : 7 6 7 6 7 6 : : : : 7: ð2Þ C mm ¼ 6 : 7 6 : : : : : 7 4 5 : : : C m;m C m;1
4. Experimental results
The elements of the covariance matrix Cm m are computed by the following: C j;k ¼
n X ðxi;j x j Þðxi;k x k Þ
n1
i¼1
;
ð3Þ
where x j and x k are the average values of column j and k, respectively. By having the covariance matrix C, the eigenvalues and eigenvectors of the matrix can be easily computed. A nonzero vector q of dimension m is an eigenvector of the square matrix Cm m, if it satisfies the Eigen equation: Cq ¼ λq;
ð4Þ
where λ is a scaler called eigenvalue. Each eigenvalue corresponds to an eigenvector qm 1. In practise, only the eigenvectors associated with the largest eigenvalues are used for PCA analysis. If we use the k eigenvectors corresponding to k largest eigenvalues λ1,… λk, then the PCA transformation reduces the dimension of an input pulse vector from m to k. By using the k eigenvectors, we can create the matrix Qm k as: Q ¼ ½ q1
q2
:::
qk :
In the first step, we calibrate the light output of the system. The calibration was performed by integrating the sampled pulses and shaping the pulses with a CR–RC filter before measuring the amplitude of the pulses [13]. Prior to the integration process, the baseline of each pulse was corrected by subtracting the average of the samples before the trigger from all pulse samples. The shaping time constant of the filter was 0.5 μs. The Compton edges in the pulse-height spectrum of data taken by 137Cs and 22Na were used for the calibration [14], indicating that our system is able of accepting pulses of light outputs from 85 keVee to about 1000 keVee (electron equivalent energy). The eigenvectors were determined by using a data matrix containing 3000 pulses acquired during the pre-acquisition run. Each pulse in the matrix contains 700 samples and the pulses are aligned in time by the trigger level of the digitizer. It is worth mentioning that the number of pulses is sufficiently large to yield a good statistical accuracy. A large number of pulses complicates the calculation of eigenvectors by increasing the size of the initial data matrix while it has little effect of the statistical accuracy of the calculations. Fig. 1 shows a typical input pulse and the first eigenvector extracted from the data matrix. In our study, the first four eigenvectors were examined for the PSD purpose. The examination of the PSD performance of the eigenvectors was performed by using a separate set of pulses. Fig. 2 shows the scatter plots of the first and second components of pulses against the corresponding light outputs. It is seen that, by using the first principal component, the neutron and gamma-rays are distributed in two different regions
ð5Þ
The matrix Q is then used to extract the k principal components of an arbitrary input pulse. To do so, the input pulse is normalized to its amplitude and then the principal components of the normalized pulse y1 m are calculated as: ð6Þ PC ¼ yQ ¼ pc1 pc2 ::: pck The projected vector PC is now of dimension k and its elements are termed the first principal component (pc1), the second principal component (pc2) and so on. In fact, the principal components are a set of uncorrelated variables that reflect most of the statistical variations in the shape of the input pulses. In our study, the principal components of input pulses are examined as a PSD index which determines the type of radiation. 3. Experimental setup In our experiments, a cylindrical NE213 liquid scintillation detector (200 200 ) coupled to a PMT of type R329 Hamamatsu was
Fig. 1. (Top) A sample input pulse and (Bottom) the first eigenvector calculated by using pules acquired during the pre-acquisition run.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
while, by using the second principal component, the two group of pulses are indistinguishable. The same results as that with the second principal component were obtained by using the third and fourth principal components indicating the fact that the pulseshape characteristics of the pulses are only reflected in the first principal component of the pulses and thus the first principal component was used as PSD index. In order to label each group of pulses in Fig. 2(A) as neutron or gamma-rays, a set of pulses from a pure gamma-ray source (137Cs) was also analyzed, indicating that pulses with smaller first principal component are gamma-rays. It has to be mentioned that, by using the first principal component, the PSD method is performed by only a simple multiplication of two vectors, including the sampled pulse vector and the first eigenvector. This means that the method is computationally very economical, comparable or slightly higher than that of the standard CC method. To evaluate the separation of the neutron and gamma-ray pulses in a quantitative way, a figure-of-merit (FOM) for each
Fig. 2. The scatter plot of the first principal component (A) and second principal component (B) against the corresponding light output. The separation between neutrons and gamma-rays can be achieved only by using the first principal component. The lower plume corresponds to gamma-rays and the upper plume corresponds to neutrons.
3
method has been calculated as [15]: FOM ¼
ΔP
FWHMn þ FWHMγ
;
ð7Þ
where ΔP is the separation between the centroids of the neutron and gamma-ray peaks in the distribution of the first principal component and FWHMn and FWHMγ are the full-width at halfmaximum (FWHM) of the neutron and gamma-ray peaks. This definition requires for both peaks being Gaussian in shape. In the case of the presented results in Fig. 2(A), due to the presence of nonlinear regions, particularly in the low light output range, both peaks are asymmetric. Therefore, the FOM values were calculated for narrow intervals of scintillation light output in which the distributions of the PSD index show fairly good Gaussian shapes for both neutron and gamma-ray events. The parameters of Eq. (7) were determined by fitting a double Gaussian function to the distribution of the first principal components for events lying in narrow light output intervals (here 100 keVee wide). Fig. 3 shows the distributions of the first principal component for two different intervals of light output. We also compared the performance of the PSD method against the conventional charge-comparison (CC) method. In general, the CC method is implemented by integrating PMT pulses for two different time intervals. Each pulse was first integrated from the beginning of the pulse to the end of the tail (here 250 ns). This integral is referred to as the total integral. The second integral was taken from an optimized starting position after the pulse maximum to the same ending point with the total integral. This integral is referred to as the tail integral. Then the ratio of the tail integral to the total integral was used as PSD index. In general, the optimization of the total integral has only a minor effect on the overall PSD performance of the CC methods, and therefore, the optimization of the CC method was done by calculating the FOM values for various length of the tail integral. FOM values of the CC method were calculated by fitting a double Gaussian function to the distribution of the PSD index. The results of FOM calculations for various starting points of the tail integral are shown in Fig. 4. The results show that the best performance is achieved by using a tail integral starting 40 ns after the peak of each pulse. After optimizing the CC method, the FOM values calculated for different light output intervals with the CC method were compared against the FOM values calculated for our PSD method. The results of the
Fig. 3. The distributions of the first principal components for events in two different light output ranges of 85–185 keVee (right) and 200–300 keVee (left).
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
4
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
5. Summary and conclusion
Fig. 4. Optimization of the CC method. The highest FOM is achieved by using a tail integral starting about 40 ns after a pulse's peak.
In this work, we applied the PCA technique for PSD application of scintillation radiation detectors. The PCA analysis was experimentally performed by discriminating between neutron and gamma-ray events with a liquid scintillation detector in a mixed radiation field. Our results show that the pulse-shape characteristics of the scintillation pulses are reflected in the first principal component of the pulses and, by using it as PSD index, a clear discrimination of neutron and gamma-rays was achieved. The discrimination performance of the method in the light output range of 60–1000 keVee was quantified by calculating the average of FOM values as FOM ¼1.7423. The performance of the method was also compared with that of the conventional CC method indicating the superior performance of the method particularly at light output range of 85–600 keVee. The PSD method based on the PCA analysis also offers the advantage that is self-calibrating allowing optimal discrimination for any type of scintillation detector. The method is also computationally cheap which makes it a powerful method for pulse-shape discrimination of various scintillation detectors.
Acknowledgment The author acknowledges the fund and support of Majmaah University, Saudi Arabia.
References
Fig. 5. A comparison of the FOM values of the proposed PSD method with the conventional CC method. The better performance of our method is apparent at light output range 85–600 keVee. The average of FOM value for the CC method is 1.5393 and for the PCA method is 1.7423.
comparison are shown in Fig. 5. In this figure the FOM values are given as a function of the beginning of 50 keVee wide light output intervals. The superior performance of the PCA method is apparent, particularly at low light output range of 85–600 keVee. The averages of FOM values over the whole examined light output range were calculated to be 1.7423 and 1.5393, respectively, for the PCA and the CC methods.
[1] F.D. Brooks, Nuclear Instruments and Methods 162 (1979) 477. [2] G.F. Knoll, Radiation Detection and Measurements, John Wiley & Sons, New York, 2000 third edition. [3] F.D. Brooks, Nuclear Instruments and Methods 4 (1959) 151. [4] T.K. Alexander, F.S. Goudling, Nuclear Instruments and Methods A 13 (1961) 244. [5] M. Nakhostin, P.M. Walker, Nuclear Instruments and Methods A 621 (2010) 498. [6] Y. Kaschuck, B. Esposito, Nuclear Instruments and Methods A 551 (2005) 420. [7] S. Marrone, et al., Nuclear Instruments and Methods A 490 (2002) 299. [8] S. Yousefi, et al., Nuclear Instruments and Methods A 598 (2002) 551. [9] G. Liu, et al., Nuclear Instruments and Methods A 607 (2009) 620. [10] M.J. Joyce, et al., IEEE Transactions on Nuclear Science NS59 (2012) 1245. [11] I.T. Jolliffe, Principal Component Analysis, Springer-Verlag, 2002, second edition. [12] Z. Cao, L.F. Miller, Nuclear Instruments and Methods A 416 (1998) 32. [13] M. Nakhostin, IEEE Transactions on Nuclear Science NS58 (2011) 2378. [14] A.A. Naqvi, et al., Nuclear Instruments and Methods A 306 (1991) 267. [15] R.A. Winyard, et al., Nuclear Instruments and Methods 95 (1971) 141.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 Q2 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 Q3108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Q1 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
Contents lists available at ScienceDirect
Nuclear Instruments and Methods in Physics Research A journal homepage: www.elsevier.com/locate/nima
Principal Component Analysis for pulse-shape discrimination of scintillation radiation detectors T. Alharbi Department of Physics, Faculty of Science in Zulfi, Almajmaah University, P.O. Box 1712, 11932, Saudi Arabia
art ic l e i nf o
a b s t r a c t
Article history: Received 21 July 2015 Received in revised form 11 October 2015 Accepted 11 October 2015
In this paper, we report on the application of Principal Component analysis (PCA) for pulse-shape discrimination (PSD) of scintillation radiation detectors. The details of the method are described and the performance of the method is experimentally examined by discriminating between neutrons and gamma-rays with a liquid scintillation detector in a mixed radiation field. The performance of the method is also compared against that of the conventional charge-comparison method, demonstrating the superior performance of the method particularly at low light output range. PCA analysis has the important advantage of automatic extraction of the pulse-shape characteristics which makes the PSD method directly applicable to various scintillation detectors without the need for the adjustment of a PSD parameter. & 2015 Published by Elsevier B.V.
Keywords: Principal component analysis Pulse-shape discrimination Scintillation detectors
1. Introduction In several scintillation detectors sensing differences in the shape of output pulses provide very useful information on the type of radiation interacting with the detector. This is commonly referred to as pulse-shape discrimination (PSD). The most common PSD application is in the discrimination between neutrons and gamma-rays with organic scintillation detectors used as fast neutron detector [1]. The PSD techniques are also widely used in applications involving phoswich scintillation detectors [2]. The traditional PSD methods are the Charge-Comparison (CC) and the Zero-crossing (ZC) methods which have been implemented on analog circuitry for several decades [3,4]. The CC method is based on the comparison of integrals of the pulses for two different time intervals while in the ZC method a zero-crossing time bears the pulse-shape information. In recent years, there has been tremendous interest in using digital pulse-processing methods for PSD of scintillation detectors and a number of techniques such as digital versions of CC and ZC methods have been used for this purpose [5,6]. Moreover, pure digital techniques based on curve fitting [7], frequency domain analysis [8] and neural networks [9] have been used for PSD applications. Furthermore, digital PSD systems with the capability of real-time operation have been developed [10]. However, the available PSD methods are either computationally time intensive or require a PSD parameter to be precisely adjusted in accordance to the pulse-shape characteristics of the detector concerned. In this paper, a new digital PSD method based on the Principal Component Analysis (PCA) of sampled
detector pulses is proposed. PCA is a widely used technique for extracting strong patterns in a dataset and has found applications in various disciplines such as image processing, neuroscience, etc. [11]. Our motivation for employing PCA for PSD applications is that PCA analysis is inherently adaptive to the intrinsic statistical characteristics of the data, and therefore, a PCA-based PSD system would be directly applicable to various scintillation detectors without the need for adjusting a PSD parameter. We describe the details of the method and the performance of the method is experimentally examined by discriminating between neutrons and gamma-rays with a liquid scintillation detector in a mixed radiation field. The performance of the PCA pulse classification algorithm is also compared against that of the conventional CC method.
2. The PCA method In most of the scintillation detectors, a PSD technique makes use of a simple basic fact that the decay-times of a scintillation light pulse are dependent on the type of radiation interacting with the detector. The shape of the pulses is, however, subject to variations due to statistical fluctuations in the number of photoelectrons, time spread in the photomultiplier tube (PMT) and electronic noise [12]. Our PSD method is based on the employment of a PCA analysis for extracting the patterns in the decay-times of sampled PMT pulses, thereby different types of radiation can be identified. A comprehensive review of PCA can be found in many
http://dx.doi.org/10.1016/j.nima.2015.10.030 0168-9002/& 2015 Published by Elsevier B.V.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
2
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
references (e.g. Ref. [11]). Here, we briefly describe the PCA procedure in regard to our PSD application. In computational terms, the PSD method is based on the extraction of the eigenvectors and eigenvalues of a data covariance matrix formed by a set of pulses acquired during a pre-acquisition run. Consider during a preacquisition run a set of n pulses are acquired. Each pulse is composed of m samples covering the whole duration of the pulse. The pulses are then normalized to their amplitude and a data matrix X is formed in such a way that each pulse represents a row of the matrix. The matrix Xn m is described as: 2 3 x1;1 x1;2 : : x1;m 6 : : : : : 7 6 7 6 7 : : : : : 7 ð1Þ X nm ¼ 6 6 7 6 : : : : : 7 4 5 : : : xn;m xn;1
used. The detector container is made of aluminium and the PMT is directly mounted on the back circular surface of the detector. The detector has a 0.5 mm thick aluminium window. The PMT was operated with a negative voltage of 1600 V and the pulses from the anode of the PMT were directly digitized by using a fast waveform digitizer at 4 GS/s sampling rate and with 10-bit resolution (model DC252HF from Agilient Technologies Inc). Tests were carried out using an americium–beryllium (Am–Be) neutron source. A lead brick was placed in front of the detector to increase the ratio of the number of recorded neutrons to gamma-rays. Some tests with standard laboratory sources were also performed to characterise the light output performance of the system. The sampled PMT pulses were analyzed with algorithms written in MATLAB programming environment.
The j-th column of this matrix represents samples of random variable xj which are the normalized amplitudes of the pulses at a specific time during their lifetime. The covariance matrix of the random variables xj will be a matrix Cm m as: 2 3 C 1;1 C 1;2 : : C 1;m 6 : : : : : 7 6 7 6 7 6 : : : : 7: ð2Þ C mm ¼ 6 : 7 6 : : : : : 7 4 5 : : : C m;m C m;1
4. Experimental results
The elements of the covariance matrix Cm m are computed by the following: C j;k ¼
n X ðxi;j x j Þðxi;k x k Þ
n1
i¼1
;
ð3Þ
where x j and x k are the average values of column j and k, respectively. By having the covariance matrix C, the eigenvalues and eigenvectors of the matrix can be easily computed. A nonzero vector q of dimension m is an eigenvector of the square matrix Cm m, if it satisfies the Eigen equation: Cq ¼ λq;
ð4Þ
where λ is a scaler called eigenvalue. Each eigenvalue corresponds to an eigenvector qm 1. In practise, only the eigenvectors associated with the largest eigenvalues are used for PCA analysis. If we use the k eigenvectors corresponding to k largest eigenvalues λ1,… λk, then the PCA transformation reduces the dimension of an input pulse vector from m to k. By using the k eigenvectors, we can create the matrix Qm k as: Q ¼ ½ q1
q2
:::
qk :
In the first step, we calibrate the light output of the system. The calibration was performed by integrating the sampled pulses and shaping the pulses with a CR–RC filter before measuring the amplitude of the pulses [13]. Prior to the integration process, the baseline of each pulse was corrected by subtracting the average of the samples before the trigger from all pulse samples. The shaping time constant of the filter was 0.5 μs. The Compton edges in the pulse-height spectrum of data taken by 137Cs and 22Na were used for the calibration [14], indicating that our system is able of accepting pulses of light outputs from 85 keVee to about 1000 keVee (electron equivalent energy). The eigenvectors were determined by using a data matrix containing 3000 pulses acquired during the pre-acquisition run. Each pulse in the matrix contains 700 samples and the pulses are aligned in time by the trigger level of the digitizer. It is worth mentioning that the number of pulses is sufficiently large to yield a good statistical accuracy. A large number of pulses complicates the calculation of eigenvectors by increasing the size of the initial data matrix while it has little effect of the statistical accuracy of the calculations. Fig. 1 shows a typical input pulse and the first eigenvector extracted from the data matrix. In our study, the first four eigenvectors were examined for the PSD purpose. The examination of the PSD performance of the eigenvectors was performed by using a separate set of pulses. Fig. 2 shows the scatter plots of the first and second components of pulses against the corresponding light outputs. It is seen that, by using the first principal component, the neutron and gamma-rays are distributed in two different regions
ð5Þ
The matrix Q is then used to extract the k principal components of an arbitrary input pulse. To do so, the input pulse is normalized to its amplitude and then the principal components of the normalized pulse y1 m are calculated as: ð6Þ PC ¼ yQ ¼ pc1 pc2 ::: pck The projected vector PC is now of dimension k and its elements are termed the first principal component (pc1), the second principal component (pc2) and so on. In fact, the principal components are a set of uncorrelated variables that reflect most of the statistical variations in the shape of the input pulses. In our study, the principal components of input pulses are examined as a PSD index which determines the type of radiation. 3. Experimental setup In our experiments, a cylindrical NE213 liquid scintillation detector (200 200 ) coupled to a PMT of type R329 Hamamatsu was
Fig. 1. (Top) A sample input pulse and (Bottom) the first eigenvector calculated by using pules acquired during the pre-acquisition run.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
while, by using the second principal component, the two group of pulses are indistinguishable. The same results as that with the second principal component were obtained by using the third and fourth principal components indicating the fact that the pulseshape characteristics of the pulses are only reflected in the first principal component of the pulses and thus the first principal component was used as PSD index. In order to label each group of pulses in Fig. 2(A) as neutron or gamma-rays, a set of pulses from a pure gamma-ray source (137Cs) was also analyzed, indicating that pulses with smaller first principal component are gamma-rays. It has to be mentioned that, by using the first principal component, the PSD method is performed by only a simple multiplication of two vectors, including the sampled pulse vector and the first eigenvector. This means that the method is computationally very economical, comparable or slightly higher than that of the standard CC method. To evaluate the separation of the neutron and gamma-ray pulses in a quantitative way, a figure-of-merit (FOM) for each
Fig. 2. The scatter plot of the first principal component (A) and second principal component (B) against the corresponding light output. The separation between neutrons and gamma-rays can be achieved only by using the first principal component. The lower plume corresponds to gamma-rays and the upper plume corresponds to neutrons.
3
method has been calculated as [15]: FOM ¼
ΔP
FWHMn þ FWHMγ
;
ð7Þ
where ΔP is the separation between the centroids of the neutron and gamma-ray peaks in the distribution of the first principal component and FWHMn and FWHMγ are the full-width at halfmaximum (FWHM) of the neutron and gamma-ray peaks. This definition requires for both peaks being Gaussian in shape. In the case of the presented results in Fig. 2(A), due to the presence of nonlinear regions, particularly in the low light output range, both peaks are asymmetric. Therefore, the FOM values were calculated for narrow intervals of scintillation light output in which the distributions of the PSD index show fairly good Gaussian shapes for both neutron and gamma-ray events. The parameters of Eq. (7) were determined by fitting a double Gaussian function to the distribution of the first principal components for events lying in narrow light output intervals (here 100 keVee wide). Fig. 3 shows the distributions of the first principal component for two different intervals of light output. We also compared the performance of the PSD method against the conventional charge-comparison (CC) method. In general, the CC method is implemented by integrating PMT pulses for two different time intervals. Each pulse was first integrated from the beginning of the pulse to the end of the tail (here 250 ns). This integral is referred to as the total integral. The second integral was taken from an optimized starting position after the pulse maximum to the same ending point with the total integral. This integral is referred to as the tail integral. Then the ratio of the tail integral to the total integral was used as PSD index. In general, the optimization of the total integral has only a minor effect on the overall PSD performance of the CC methods, and therefore, the optimization of the CC method was done by calculating the FOM values for various length of the tail integral. FOM values of the CC method were calculated by fitting a double Gaussian function to the distribution of the PSD index. The results of FOM calculations for various starting points of the tail integral are shown in Fig. 4. The results show that the best performance is achieved by using a tail integral starting 40 ns after the peak of each pulse. After optimizing the CC method, the FOM values calculated for different light output intervals with the CC method were compared against the FOM values calculated for our PSD method. The results of the
Fig. 3. The distributions of the first principal components for events in two different light output ranges of 85–185 keVee (right) and 200–300 keVee (left).
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
4
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
T. Alharbi / Nuclear Instruments and Methods in Physics Research A ∎ (∎∎∎∎) ∎∎∎–∎∎∎
5. Summary and conclusion
Fig. 4. Optimization of the CC method. The highest FOM is achieved by using a tail integral starting about 40 ns after a pulse's peak.
In this work, we applied the PCA technique for PSD application of scintillation radiation detectors. The PCA analysis was experimentally performed by discriminating between neutron and gamma-ray events with a liquid scintillation detector in a mixed radiation field. Our results show that the pulse-shape characteristics of the scintillation pulses are reflected in the first principal component of the pulses and, by using it as PSD index, a clear discrimination of neutron and gamma-rays was achieved. The discrimination performance of the method in the light output range of 60–1000 keVee was quantified by calculating the average of FOM values as FOM ¼1.7423. The performance of the method was also compared with that of the conventional CC method indicating the superior performance of the method particularly at light output range of 85–600 keVee. The PSD method based on the PCA analysis also offers the advantage that is self-calibrating allowing optimal discrimination for any type of scintillation detector. The method is also computationally cheap which makes it a powerful method for pulse-shape discrimination of various scintillation detectors.
Acknowledgment The author acknowledges the fund and support of Majmaah University, Saudi Arabia.
References
Fig. 5. A comparison of the FOM values of the proposed PSD method with the conventional CC method. The better performance of our method is apparent at light output range 85–600 keVee. The average of FOM value for the CC method is 1.5393 and for the PCA method is 1.7423.
comparison are shown in Fig. 5. In this figure the FOM values are given as a function of the beginning of 50 keVee wide light output intervals. The superior performance of the PCA method is apparent, particularly at low light output range of 85–600 keVee. The averages of FOM values over the whole examined light output range were calculated to be 1.7423 and 1.5393, respectively, for the PCA and the CC methods.
[1] F.D. Brooks, Nuclear Instruments and Methods 162 (1979) 477. [2] G.F. Knoll, Radiation Detection and Measurements, John Wiley & Sons, New York, 2000 third edition. [3] F.D. Brooks, Nuclear Instruments and Methods 4 (1959) 151. [4] T.K. Alexander, F.S. Goudling, Nuclear Instruments and Methods A 13 (1961) 244. [5] M. Nakhostin, P.M. Walker, Nuclear Instruments and Methods A 621 (2010) 498. [6] Y. Kaschuck, B. Esposito, Nuclear Instruments and Methods A 551 (2005) 420. [7] S. Marrone, et al., Nuclear Instruments and Methods A 490 (2002) 299. [8] S. Yousefi, et al., Nuclear Instruments and Methods A 598 (2002) 551. [9] G. Liu, et al., Nuclear Instruments and Methods A 607 (2009) 620. [10] M.J. Joyce, et al., IEEE Transactions on Nuclear Science NS59 (2012) 1245. [11] I.T. Jolliffe, Principal Component Analysis, Springer-Verlag, 2002, second edition. [12] Z. Cao, L.F. Miller, Nuclear Instruments and Methods A 416 (1998) 32. [13] M. Nakhostin, IEEE Transactions on Nuclear Science NS58 (2011) 2378. [14] A.A. Naqvi, et al., Nuclear Instruments and Methods A 306 (1991) 267. [15] R.A. Winyard, et al., Nuclear Instruments and Methods 95 (1971) 141.
Please cite this article as: T. Alharbi, Nuclear Instruments & Methods in Physics Research A (2015), http://dx.doi.org/10.1016/j. nima.2015.10.030i
67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 Q2 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 Q3108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132