New method and system for HPGe gamma spectrum analysis

Abstract

The application provides a novel HPGe gamma spectrum analysis method and system, relates to the technical field of spectrum analysis, and comprises the following steps: S1, smoothing the spectrum through a low-pass filter; S2, peak searching of the smoothed spectrum through a CWT method; S3, after the peak searching, calculating the peak width by using an improved secondary convolution method; S4, background deduction through four-order filter adaptive SNIP; S5, calibrating the peak shape filter through an experimental method; and S6, spectrum analysis through regularized deconvolution solving. The application can realize the process of overlapping peak decomposition without the aid of a nuclide library, improves the measurement accuracy and precision, and provides key technical support for the safe disposal of LILW.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of energy spectrum analysis, in particular to a novel HPGe gamma spectrum analysis method and system. BACKGROUND

[0002] A million kilowatt nuclear power unit will produce about 13,000 cubic meters of low and medium level radioactive waste (LILW) in 60 years, which is solidified or compressed into 200L, 400L waste barrels, with a total volume of about 3000 cubic meters. IAEA pointed out in the technical report that in order to monitor safety, LILW needs to be characterized throughout its life cycle. At the same time, according to IAEA safety standards, LILW of different activity ranges is disposed of at different depths, so it is necessary to quantitatively analyze LILW. Since destructive analysis needs to reopen the container and it takes time to drill holes in the container, non-destructive analysis (NDA) is currently the best choice. Some important LILW nuclides are alpha or beta nuclides, which cannot be subjected to NDA, so these nuclides are called difficult-to-measure nuclides (DTM). ISO stipulates in the standard that by carrying out non-destructive gamma scanning on key nuclides such as Cs-137 and Co-60, and by establishing a proportionality factor between the key nuclides and the DTM, the activity of the DTM can be calculated. Therefore, for the characterization of LILW, gamma scanning is first required.

[0003] ASTM stipulates that NDA includes segmented gamma scanning (SGS) and tomographic gamma scanning (TGS), and TGS is suitable for LILW with non-uniform distribution of density and nuclides. Regardless of which method is used, first of all, the gamma spectrum of LILW needs to be qualitatively and quantitatively analyzed by using a gamma detector. Since LILW has many types of nuclides, the types of nuclides from different waste streams are also different, and the nuclides are difficult to identify, so ASTM requires the use of a high-resolution detector. ISO points out in the international standard that high-resolution detectors include high-purity germanium (HPGe) detectors and cadmium zinc telluride (CdZnTe or CZT) detectors. However, the size of CZT is difficult to make large, so the detection efficiency of CZT is only 1% of that of HPGe, and it is more difficult for high-energy gamma rays to deposit in CZT, so HPGe detectors are selected internationally for LILW measurement.

[0004] The spectrum measured by HPGe detector needs to be analyzed to obtain the counts of different energy rays. Spectrum analysis methods mainly include: smoothing, background subtraction, peak searching and peak area calculation, among which nuclide identification relies on the peak searching process. The most mature peak searching method is the derivative method proposed by Mariscotti, which is adopted by Gamma Vision (GV), Genie 2000 and the like. However, due to the limited peak searching capability of the derivative method, for an overlapping peak composed of two full-energy peaks with a very close distance, a nuclide library needs to be set to assist in identification. This method assumes that all γ rays in the nuclide library exist in the spectrum, and therefore attempts to match each energy peak listed in the library. Therefore, for the same peak, different nuclide libraries often give different results. For example, the common LILW nuclides have Cs-137 emitting 662 keV and Ag-110m emitting 658 keV, and before LILW measurement, we cannot accurately know the nuclide composition. For the same 662 keV full-energy peak, if the nuclide library only sets 662 keV, the quantitative analysis result is accurate; but if the nuclide library gives both 658 keV and 662 keV, it will certainly give a 658 keV full-energy peak count, causing nuclide misidentification, and the calculated 662 keV peak area is reduced, causing inaccuracy in activity calculation, which poses a safety hazard for LILW disposal. The IAEA report also points out that GV needs a good nuclide library. For other peak searching methods such as IF function, symmetric zero area transformation, Gaussian integral function, covariance, continuous wavelet transform (CWT), neural network, compared with the derivative method, there is certain improvement, but ultimately still need the assistance of nuclide library for identification.

[0005] Firstly, several basic steps of traditional spectrum analysis method are introduced, and the existing problems are expounded. For the smoothing process of the spectrum, the five-point smoothing formula is used:

[0006] S i =(O i-2 +4O i-1 +6O i +4O i+1 + i+2 ) / 16

[0007] In the formula, i is the channel address, O i is the count of each channel in the original spectrum, and S i is the count of each channel after smoothing

[0008] The problem here is that this algorithm is smoothed for each channel of the spectrum, so each time the processing will result in the peak value of the full energy peak to be reduced, and the counts on both sides of the peak bottom to be increased, so that the peak area calculation is inaccurate. At the same time, for the overlapping peak composed of two full energy peaks with similar energies, it is more difficult to identify the nuclide after several smoothings using this method.

[0009] The background subtraction object includes single peak and multi-peak. For single peak, linear background is used; for multi-peak, according to the different trends of background and peak height, it is divided into step background and parabolic background, as shown in Figure 1 The left is linear background, the middle is step background, and the right is parabolic background.

[0010] The problem here is that linear background and step background are only rough approximations of the background. At the same time, studies have shown that parabolic background has poor robustness. Moreover, regardless of which method is used, for the full energy peak falling on the Compton edge position, accurate background fitting cannot be performed, resulting in increased peak area calculation error.

[0011] The peak searching method is the second derivative method. The problem here is that in addition to the limited ability of overlapping peak separation, this peak searching method is very susceptible to noise. Therefore, peak searching often needs to smooth the spectrum several times, resulting in the problems of peak value reduction, peak area calculation error increase, and overlapping peak difficulty in effective separation mentioned above.

[0012] The last step is to quantitatively analyze the overlapping peak. If the derivative method's peak searching result matches the nuclide library, the peak can be fitted and calculated by the function. If the derivative method can only identify a single peak, it can only rely on the auxiliary identification of the nuclide library, resulting in the instability of the analysis result caused by the uncertainty of the nuclide library. Therefore, each analysis step in the traditional spectrum analysis method has deficiencies, and the forced matching of the nuclide library for the overlapping peak that cannot be decomposed is the biggest problem.

[0013] In summary, the problems of the traditional smoothing method are:

[0014] 1) For each channel of the spectrum, the data is smoothed, so each time the processing will result in the peak value of the full energy peak to be reduced, and the counts on both sides of the peak bottom to be increased, so that the peak area calculation is inaccurate. At the same time, for the overlapping peak composed of two full energy peaks with similar energies, it is more difficult to identify the nuclide after several smoothings using this method.

[0015] 2) The problem of the traditional background subtraction method: poor robustness; for the full energy peak falling on the Compton edge position, accurate background fitting cannot be performed, resulting in increased peak area calculation error.

[0016] 3) Traditional peak searching method has the following problems: limited separation ability of overlapping peaks; vulnerable to noise interference; peak searching often needs multiple smoothing of the energy spectrum, which leads to the problems of peak value drop, peak area calculation error increase and difficulty in effectively separating overlapping peaks.

[0017] 4) Overlapping peak separation: depends on the accuracy of the nuclide library.

[0018] Term explanation:

[0019] High purity germanium (HPGe) detector: HPGe refers to a single crystal of germanium with stable net concentration of electrical active impurities at room temperature, and the typical value of the net concentration of impurities is less than 3*1010 cm-3. Under appropriate bias, a detector made of a high purity germanium single crystal of a conventional size can be fully depleted. The semiconductor detector made of HPGe is called HPGe detector, which is the best energy resolution detector on the market.

[0020] Energy resolution: the ratio of the full width at half maximum (FWHM) of the full energy peak to the peak energy, which represents the ability of the detector to distinguish different energy rays, and is therefore the most important performance indicator of the spectrometer detector. The actual measured energy resolution is related to the processes of signal generation, transmission, conversion, amplification and collection of the detector. The stronger the useful signal and the weaker the interference factors, the better the energy resolution.

[0021] Gamma spectrum: the energy spectrum refers to the distribution curve of the count rate with particle energy after the pulse amplitude is energy scaled. The gamma spectrum refers to the information of the reaction in the energy spectrum is gamma rays.

[0022] Radioactive waste: radioactive waste is a substance containing radioactive nuclides or contaminated by radioactive nuclides, whose concentration or specific activity is greater than the clean control level specified by the national regulatory department, and which is not expected to be used again. According to the form, it can be divided into gaseous waste, liquid waste and solid waste.

[0023] Low and medium level solid waste: radioactive waste is divided into five categories: extremely short-lived radioactive waste, extremely low level radioactive waste, low level radioactive waste, medium level radioactive waste and high level radioactive waste. Radioactive waste solidification is to convert gaseous, liquid or solid waste into a monolithic solidified body with performance indicators meeting the disposal requirements. The purpose is to form an object suitable for loading, transportation and temporary storage, and the performance meets the disposal requirements. Low and medium level solid waste refers to the waste formed by using solidification technology to solidify and prepare low level and medium level radioactive waste. Summary of the application

[0024] In view of the defects in the prior art, the application provides a novel HPGe gamma spectrum analysis method and system.

[0025] According to the application, a new HPGe gamma spectrum analysis method and system are provided, and the method comprises the following steps:

[0026] In a first aspect, a new HPGe gamma spectrum analysis method is provided, and the method comprises the following steps:

[0027] Step S1: smoothing the original spectrum through a low-pass filter to obtain a smoothed spectrum;

[0028] Step S2: performing peak searching on the smoothed spectrum through a CWT method;

[0029] Step S3: calculating the peak width based on the peak position obtained after the peak searching, and using an improved double convolution method;

[0030] Step S4: combining the peak position and the peak width, and performing background deduction through a fourth-order filter adaptive SNIP to obtain a net spectrum;

[0031] Step S5: calibrating a peak shape filter through an experimental method;

[0032] Step S5: using the peak shape filter to perform deconvolution on the net spectrum, and solving through a regularization iteration method to obtain a result.

[0033] Preferably, the step S1 comprises: performing smoothing step by step for each channel i in the spectrum;

[0034] For a certain channel x, the corresponding count is y(i), and an interval is taken on the left and right sides of x, and the counts in the interval are added to obtain total counts L(i) and R(i);

[0035] The interval is selected as 1.25 FWHM, and the total count in the peak area is S, so S=y(i)+L(i)+R(i), and a judgment is made for S as follows:

[0036] If S≥N, the conventional five-point smoothing is performed; if S

[0037] If 1 / <(+1) / (+1)<r, smoothing is performed; wherein r is a constant, and different detector parameters are different.

[0038] Preferably, the step S2 comprises:

[0039] Calculate the scale image:

[0040]

[0041] In the formula, T(E,s) is the scale image, E is the ray energy, s is the scale of the wavelet, and ψ E,s() is the mother wavelet; scale s represents stretching and shrinking of the wavelet, as shown in the following formula:

[0042]

[0043] The selection of the wavelet is a key to the result of the spectrum analysis, and the gamma spectrum total peak is a Gaussian function, and the Mexh is similar to the Gaussian peak shape, so the Mexh is selected, and the definition of the Mexh is as follows:

[0044]

[0045] In the formula, ψ(t) represents the Mexh wavelet function; t represents the position of the wavelet; σ represents a standard deviation, and e is a constant;

[0046] The position information of the peak is obtained by analyzing the scale image.

[0047] Preferably, the step S3 comprises:

[0048] The expression of the total peak is:

[0049]

[0050] In the formula, x is a channel address, A is a peak value count, a is a peak position, and σ is a standard deviation of the total peak;

[0051] The calculation result of the second convolution of the Gaussian first derivative of the above formula is:

[0052]

[0053] In the formula, δ represents the width of the filter;

[0054] For the formula, when dC() / ≠0, δ=σ is easily obtained;

[0055] The absolute value of C(x) is taken, and by traversing the possible values of δ, a matrix M(x, δ) can be obtained, wherein x is a channel;

[0056] The position of the peak is calculated by using the CWT, each column curve corresponding to the peak position in the matrix M(x, δ) is analyzed, there is a maximum value, the δ value corresponding to the maximum value is the σ value of the peak position, and [a-3σ, a+3σ] is the width of the peak region.

[0057] Preferably, the step S4 comprises:

[0058] The dynamic range of each channel count in the compressed spectrum is compressed, so as to accelerate iteration, and the LLS operator is applied to the spectrum data:

[0059]

[0060] where i is a channel index, y(i) is the count of each channel, v i is a working vector;

[0061] For the count v(i) of each channel, a parameter m is selected, then w=2m+1;

[0062] The second-order filter shown in the following formula is used to iteratively calculate v1(i), 2(i), and so on until v m (i):

[0063]

[0064] Finally, the LLS inverse operator is operated on v m (i) to obtain the background spectrum;

[0065] A fourth-order clipping filter is applied to the background deduction in the peak region, and the definition of the fourth-order clipping filter is as follows:

[0066]

[0067] Wherein, the range of p is 1~m.

[0068] Preferably, the step S5 comprises:

[0069] An experimental method is used to calibrate FWHM; since the FWHM at different energies are different, the following formula is used for parameter fitting:

[0070]

[0071] Wherein, Energy is energy.

[0072] Based on the FWHM information, a peak shape filter is calculated by combining a Gaussian function.

[0073] Preferably, the step S6 comprises:

[0074] Through the response function of the detector, the relationship between the ray incident spectrum and the energy deposition spectrum is established, as shown in the following formula:

[0075]

[0076] Wherein, x(τ) is the incident spectrum, y(t) is the energy deposition spectrum, h(t) is the response function, t is the channel index, and n(t) is the system noise; the peak-shaped response function h(t) is used as a filter, and the incident spectrum is obtained through the deconvolution process to realize the separation of overlapping peaks;

[0077] For a discrete system, h(t) is expressed as a matrix:

[0078]

[0079] In the formula, each column of the matrix is the probability density function of the full-energy peak response, and the subscript of h represents the ray emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a-n represent the peak widths of different full-energy peaks; and finally, the Boosted-MLEM algorithm is used to solve the formula y(t).

[0080] In a second aspect, a novel HPGe gamma spectrum analysis system is provided, and the system comprises:

[0081] Module M1: smoothing the original spectrum through a low-pass filter to obtain a smoothed spectrum;

[0082] Module M2: peak searching on the smoothed spectrum through a CWT method;

[0083] Module M3: based on the peak positions obtained after peak searching, using an improved two-convolution method to calculate the peak width;

[0084] Module M4: combining the peak positions and the peak width, performing background subtraction through a fourth-order filter adaptive SNIP to obtain a net spectrum;

[0085] Module M5: calibrating the peak shape filter through an experimental method;

[0086] Module M5: using the peak shape filter to perform deconvolution on the net spectrum, and solving through a regularization iteration method to obtain a result.

[0087] Preferably, the module M1 comprises: step-by-step smoothing for each channel i in the spectrum;

[0088] For a certain channel x, the corresponding count is y(i), and an interval is taken on the left and right sides of x, and the counts in the interval are added to obtain total counts L(i) and R(i);

[0089] The interval is selected as 1.25 FWHM; let the total count in the peak region be S, then S=y(i)+L(i)+R(i), and a judgment is made for S as follows:

[0090] If S≥N, perform regular five-point smoothing; if S<N, do not perform smoothing; wherein A is a constant and is different for different detector parameters;

[0091] If 1 / <(+1) / (+1)<r, perform smoothing; wherein r is a constant and is different for different detector parameters;

[0092] The module M2 comprises:

[0093] Calculate the scale image:

[0094]

[0095] In the formula, T(E, s) is a scale image; E is a ray energy, s is a scale of a wavelet, ψ E,s () is a mother wavelet; the scale s represents stretching and shrinking of the wavelet, as shown in the following formula:

[0096]

[0097] Selection of the wavelet is a key of a spectrum analysis result, a gamma spectrum total peak is a Gaussian function, and the Mexh is similar to the Gaussian peak shape, so the Mexh is selected, and the Mexh is defined as follows:

[0098]

[0099] In the formula, ψ(t) represents a Mexh wavelet function; t represents a position of the wavelet; and sigma represents a standard deviation;

[0100] Position information of the peak is obtained by analyzing the scale image;

[0101] The module M3 comprises:

[0102] The total peak expression is:

[0103]

[0104] In the formula, x is a channel, A is a peak value count, a is a peak position, sigma is a standard deviation of the total peak, and e is a constant;

[0105] Position information of the peak is obtained by analyzing the scale image;

[0106] The Gaussian first derivative second convolution calculation result of the above formula is:

[0107]

[0108] In the formula, delta represents a width of a filter;

[0109] For the formula, when dC() / ≠0, delta is equal to sigma;

[0110] The absolute value of C(x) is taken, and by traversing possible values of delta, a matrix M(x, delta) can be obtained, wherein x is a channel;

[0111] The peak position is calculated by using the CWT, each column curve corresponding to the peak position in the matrix M(x, delta) is analyzed, there is a maximum value, the delta value corresponding to the maximum value is the sigma value of the peak position, and [a-3sigma, a+3sigma] is a width of a peak region;

[0112] The module M4 comprises:

[0113] Compress the dynamic range of each channel count in the spectrum to accelerate iteration, and apply the LLS operator to the spectrum data:

[0114]

[0115] In the formula, i is a channel address, y(i) is the count of each channel, v i is a working vector;

[0116] For the count v(i) of each channel, a parameter m is selected, then w=2m+1;

[0117] Iteratively calculate v1(i), 2(i), and so on to v m (i) using the second-order filter shown in the following formula:

[0118]

[0119] Finally, perform LLS inverse operator operation on v m (i) to obtain the background spectrum;

[0120] In the peak region, a fourth-order clipping filter is applied for background subtraction, and the fourth-order clipping filter is defined as:

[0121]

[0122] Wherein, the range of p is 1~m;

[0123] The module M5 comprises:

[0124] An experimental method is used to calibrate FWHM; since FWHM at different energies are different, a parameter fitting is performed using the following formula:

[0125]

[0126] Wherein, Energy is energy;

[0127] Based on FWHM information, a peak shape filter is calculated in combination with a Gaussian function;

[0128] The module M6 comprises:

[0129] Through the response function of the detector, the relationship between the ray incident spectrum and the energy deposition spectrum is established, as shown in the following formula:

[0130]

[0131] In the formula, x (tau) is incident spectrum, y (t) is energy deposition spectrum, h (t) is response function, t is channel address, and n (t) is system noise; the peak-shaped response function h (t) is taken as a filter, the incident spectrum is obtained through a deconvolution process, and overlapping peaks are separated;

[0132] For a discrete system, h (t) is expressed as a matrix:

[0133]

[0134] In the formula, each column of the matrix is a probability density function of a full-energy peak response, and the subscript of h indicates a ray emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a to n represent the peak widths of different full-energy peaks; finally, the Boosted-MLEM algorithm is used to solve formula y (t).

[0135] In a third aspect, a computer readable storage medium storing a computer program is provided, and the computer program is executed by a processor to implement the steps of the novel HPGe gamma spectrum analysis method.

[0136] Compared with the prior art, the present application has the following beneficial effects:

[0137] 1. The present application and the GV with an accurate nuclide library have close spectrum analysis capabilities. Compared with the results of the GV with an accurate nuclide library, through the analysis of Cs-137 spectra without attenuation and after water and sand attenuation, the errors of the new method and the GV are both less than 3 %; through the analysis of Eu-152 different energy rays, the average errors of the two methods are both less than 4 %;

[0138] 2. The present application does not need a nuclide library, and the analysis result is more stable. For the same Cs-137 spectrum, the GV uses different nuclide libraries to give different analysis results, and the errors are 2.7 %, 32.3 % and cannot be identified, respectively, and the analysis result is very unstable; and the new method maintains a calculation error of 0.4 %;

[0139] 3. The present application can quantitatively analyze all nuclides of LILW without a nuclide library. The overlapping peaks composed of three full-energy peaks (ray energy is 601 keV, 603 keV and 605 keV) with the closest energy in the key nuclides of LILW are analyzed. For the overlapping peaks with peak area ratios of 5:1:10, 10:1:20 and 20:1:30, the nuclides are accurately identified, and the error of a strong peak is less than 0.5 %, and the maximum error of a weak peak is 5.74 %. Even in the limit condition of a peak area ratio of 50:1:100, the weak peak can still be accurately identified;

[0140] 4、The application improves the accuracy and precision of LILW measurement, provides key technical support for reasonable disposal of LILW, and guarantees the safety of human and environment.

[0141] Other beneficial effects of the application will be described in the specific embodiments by introducing specific technical features and technical solutions, and those skilled in the art should understand the beneficial technical effects brought by the technical features and technical solutions through the introduction. BRIEF DESCRIPTION OF DRAWINGS

[0142] Other features, objects and advantages of the application will become more apparent through reading the following detailed description of non-restrictive embodiments with reference to the attached drawings:

[0143] Figure 1 Schematic diagram of background deduction method for GammaVision;

[0144] Figure 2 Comparison of the method of the application and the traditional method;

[0145] Figure 3 Mexh Wavele waveform diagram;

[0146] Figure 4 Comparison of Cs-137 original energy spectrum and normalized energy spectrum of different media attenuation;

[0147] Figure 5 Analysis results of Cs-137 (air) energy spectrum by different methods;

[0148] Figure 6 Analysis results of Cs-137 (water) energy spectrum by different methods;

[0149] Figure 7 Analysis results of Cs-137 (sand) energy spectrum by different methods;

[0150] Figure 8 Eu-152 energy spectrum;

[0151] Figure 9 Analysis results of Eu-152 energy spectrum by the new method;

[0152] Figure 10 Comparison of Cs-137 original energy spectrum, GV results based on nuclide library "Cs-137", GV results based on nuclide library "Cs-137+Ag-110m", and GV results based on nuclide library "Ag-110m";

[0153] Figure 11 Analysis results of overlapping peaks with peak area ratio of 5:1:10;

[0154] Figure 12 Analysis result of overlapping peaks with peak area ratio of 10:1:20;

[0155] Figure 13 Analysis result of overlapping peaks with peak area ratio of 20:1:30;

[0156] Figure 14 Analysis result of overlapping peaks with peak area ratio of 50:1:100. DETAILED DESCRIPTION

[0157] The present application will be described in detail below with specific examples. The following examples will help those skilled in the art to further understand the present application, but do not limit the present application in any form. It should be noted that those skilled in the art can make several changes and improvements without departing from the concept of the present application. These are within the scope of protection of the present application.

[0158] The embodiment of the present application provides a new HPGe gamma spectrum analysis method. In view of the problem that the uncertainty of the nuclide library causes instability of the analysis result and reduces the activity reconstruction accuracy in the analysis process of overlapping peaks in the LILW spectrum by using the existing method, the method specifically includes the following contents. Figure 2

[0159] Step S1: Smooth the original spectrum by a low-pass filter to obtain a smoothed spectrum.

[0160] Basic idea: Only the low count (high statistical fluctuation) part of the spectrum is smoothed, and the peak region is not smoothed.

[0161] Smooth each channel i in the spectrum. For a certain channel x, the corresponding count is y(i), and an interval is taken on the left and right of x, and the counts in the interval are added to obtain the total count L(i) and R(i).

[0162] The interval is selected as 1.25FWHM; let the total count in the peak region be S, then S=y(i)+L(i)+R(i), and the judgment is made as follows:

[0163] If S≥N, the conventional five-point smoothing is performed; if S<N, no smoothing is performed; wherein A is a constant, different for different detector parameters;

[0164] If 1 / <(+1) / (+1)<r, smoothing is performed; wherein r is a constant, different for different detector parameters.

[0165] Step S2: Peak position recognition is performed on the smoothed spectrum by the CWT method.

[0166] ​The scale image is calculated:

[0167]

[0168] In the formula, T(E, s) is the scale image; E is the ray energy, s is the scale of the wavelet, ψ E,s is the mother wavelet; the scale s represents stretching and shrinking of the wavelet, as shown in the following formula:

[0169]

[0170] The selection of the wavelet is the key to the result of the spectral analysis. Considering that the full-energy peak of the gamma spectrum is a Gaussian function, because the Mexican hat wavelet (Mexh) is similar to the Gaussian peak shape, the Mexh is selected. The definition of the Mexh is shown in the following formula, and the waveform is shown in Figure 3 .

[0171]

[0172] In the formula, ψ(t) represents the Mexh wavelet function; t represents the position of the wavelet; σ represents the standard deviation, and e is a constant;

[0173] The position information of the peak is obtained by analyzing the scale image.

[0174] Step S3: Based on the peak position obtained in step S2, the peak width is calculated using the improved quadratic convolution method.

[0175] The expression of the full-energy peak is:

[0176]

[0177] In the formula, x is the channel address, A is the peak count, a is the peak position, and σ is the standard deviation of the full-energy peak;

[0178] The calculation result of the Gaussian first derivative quadratic convolution on the above formula is:

[0179]

[0180] In the formula: δ represents the width of the filter;

[0181] For the formula, when dC() / ≠0, it is easy to obtain δ=σ. By taking the absolute value of C(x) and traversing the possible values of δ, a matrix M(x, δ) can be obtained, where x is the channel.

[0182] The position of the peak is calculated using the CWT, and each column curve corresponding to the peak position in the matrix M(x, δ) is analyzed. There is a maximum value, and the δ value corresponding to the maximum value is the σ value under the peak position. Therefore, [a-3σ, a+3σ] is the width of the peak region.

[0183] Step S4: combining the peak position and width information of steps S2 and S3, background subtraction is performed by fourth-order filter adaptive SNIP to obtain a net spectrum.

[0184] The dynamic range of each channel count in the compressed spectrum is compressed to accelerate iteration, and the LLS operator is applied to the energy spectrum data:

[0185]

[0186] In the formula, i is a channel address, y(i) is the count of each channel, v i is a working vector;

[0187] For the count v(i) of each channel, a parameter m is selected, then w=2m+1; the second-order filter shown in the following formula is used to iteratively calculate v1(i), 2(i), and so on until v m (i):

[0188]

[0189] Finally, the LLS inverse operator is operated on v m (i) to obtain a background spectrum; the second-order filter shown in the formula can effectively remove the linear background. However, in order to solve the problem of being unable to subtract the background under the Compton edge of the full-energy peak, the present application proposes a new method.

[0190] Since the Compton edge can be approximated by a cubic function to some extent, and the cubic function can be removed by a fourth-order clipping filter, the present application proposes to apply a fourth-order clipping filter in the peak region to perform background subtraction, and the formula derivation parameters come from Pascal’s triangle formula. The definition of the fourth-order clipping filter is:

[0191]

[0192] Wherein, the range of p is 1~m;

[0193] Step S5: the peak shape filter is calibrated by an experimental method.

[0194] The FWHM is calibrated by an experimental method; since the FWHM under different energies is different, the following formula is used for parameter fitting:

[0195]

[0196] Wherein, Energy is energy;

[0197] Based on the FWHM information, the peak shape filter is calculated by combining the Gaussian function.

[0198] Step S6: deconvolution of the net spectrum of step S4 using the peak shape filter of step S5, and the result is obtained by solving the regularized iterative method.

[0199] The relationship between the incident spectrum and the energy deposition spectrum is established by the response function of the detector, as shown in the following formula:

[0200]

[0201] In the formula, x (τ) is the incident spectrum, y (t) is the energy deposition spectrum, h (t) is the response function, t is the channel address, and n (t) is the system noise; therefore, the essence of the method is to use the peak-shaped response function h (t) as a filter to obtain the incident spectrum through the deconvolution process, and realize the separation of overlapping peaks.

[0202] For a discrete system, h (t) is expressed as a matrix:

[0203]

[0204] In the formula, each column of the matrix is the probability density function of the full-energy peak response, and the subscript of h indicates the rays emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a-n represent the peak widths of different full-energy peaks; finally, the Boosted-MLEM algorithm is used to solve the formula y (t).

[0205] The application also provides a novel HPGe gamma spectrum analysis system, which can be realized by executing the process steps of the novel HPGe gamma spectrum analysis method, that is, the novel HPGe gamma spectrum analysis method can be understood by those skilled in the art as the preferred embodiment of the novel HPGe gamma spectrum analysis system.

[0206] Module M1: smoothing the spectrum by a low-pass filter.

[0207] Basic idea: only the low-count (high statistical fluctuation) part of the spectrum is smoothed, and the peak region is not smoothed.

[0208] For each channel i in the spectrum, smoothing is performed step by step. For a certain channel x, the corresponding count is y (i), and an interval is taken on the left and right of x, and the counts in the interval are added to obtain the total counts L (i) and R (i).

[0209] The selection of the interval is 1.25 FWHM; let the total count in the peak region be S, then S=y (i)+L (i)+R (i), and the judgment is made for S as follows:

[0210] If S≥N, the five-point smoothing is performed; if S<N, no smoothing is performed; where A is a constant, different for different detector parameters.

[0211] If 1 / <(+1) / (+1)<r, smoothing is performed; where r is a constant, different for different detector parameters.

[0212] Module M2: Peak searching is performed on the smoothed energy spectrum by the CWT method.

[0213] The scale image is calculated:

[0214]

[0215] In the formula, T(E,s) is the scale image; E is the ray energy, s is the scale of the wavelet, ψ E,s is the mother wavelet; the scale s represents stretching and shrinking of the wavelet, as shown in the following formula:

[0216]

[0217] The selection of the wavelet is the key to the spectral analysis result. Considering that the full-energy peak of the gamma spectrum is a Gaussian function, because the Mexican hat wavelet (Mexh) is similar to the Gaussian peak shape, the Mexh is selected. The definition of the Mexh is shown in the following formula, and the waveform is shown in Figure 3 .

[0218]

[0219] In the formula, ψ(t) represents the Mexh wavelet function; t represents the position of the wavelet; σ represents the standard deviation;

[0220] The position information of the peak is obtained by analyzing the scale image.

[0221] Module M3: After the peak searching, the peak width is calculated by using the improved quadratic convolution method.

[0222] The full-energy peak expression is:

[0223]

[0224] In the formula, x is the channel address, A is the peak count, a is the peak position, and σ is the standard deviation of the full-energy peak;

[0225] The calculation result of the first-order derivative quadratic convolution of the above formula is:

[0226]

[0227] In the formula: δ represents the width of the filter;

[0228] For this formula, when dC() / =0, it is easy to get delta = sigma. Taking the absolute value of C(x), by traversing the possible values of delta, a matrix M(x, delta) can be obtained, where x is the channel.

[0229] The position of the peak is calculated by CWT, and each column curve corresponding to the peak position in the matrix M(x, delta) is analyzed, and there is a maximum value, then the delta value corresponding to the maximum value is the sigma value under the peak position, and [a-3sigma, a+3sigma] is the width of the peak region.

[0230] Module M4: background subtraction by fourth-order filter adaptive SNIP.

[0231] Compress the dynamic range of each channel count in the spectrum to speed up iteration, and apply LLS operator to the spectrum data:

[0232]

[0233] In the formula, i is the channel address, y(i) is the count of each channel, v i is the working vector;

[0234] For the count v(i) of each channel, select a parameter m, then w = 2m+1; in turn, use the second-order filter shown in the following formula to iteratively calculate v1(i), 2(i), and so on until v m (i):

[0235]

[0236] Finally, the LLS inverse operator is operated on v m (i) to obtain the background spectrum; the formula shows a second-order filter, which can effectively remove the linear background. However, in order to solve the problem of not being able to subtract the background of the Compton edge under the full energy peak, the present application proposes a new method.

[0237] Since the Compton edge can be approximated by a cubic function to some extent, and the cubic function can be removed by a fourth-order clipping filter, the present application proposes to apply a fourth-order clipping filter in the peak region to subtract the background, and the formula derivation parameters come from Pascal's triangle formula. The definition of the fourth-order clipping filter is:

[0238]

[0239] Wherein, the range of p is 1~m.

[0240] Module M5: calibrate the peak shape filter by experimental method.

[0241] The FWHM is scaled by experiment; since the FWHM under different energy is different, the following formula is used for parameter fitting:

[0242]

[0243] Wherein, Energy is energy;

[0244] Based on the FWHM information, the peak shape filter is calculated by combining the Gaussian function.

[0245] Module M6: energy spectrum analysis is carried out by regularized deconvolution solution.

[0246] Through the response function of the detector, the relationship between the incident spectrum and the energy deposition spectrum is established, as shown in the following formula:

[0247]

[0248] In the formula, x(t) is the incident spectrum, y(t) is the energy deposition spectrum, h(t) is the response function, t is the channel address, and n(t) is the system noise; therefore, the essence of the method is to take the peak-shaped response function h(t) as a filter, and obtain the incident spectrum by deconvolution process to realize the separation of overlapping peaks.

[0249] For a discrete system, h(t) is expressed as a matrix:

[0250]

[0251] In the formula, each column of the matrix is the probability density function of the full-energy peak response, and the subscript of h indicates the ray emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a-n represent the peak widths of different full-energy peaks; finally, the Boosted-MLEM algorithm is used to solve the formula y(t).

[0252] Next, the present application is more specifically described.

[0253] The verification experiment of the energy spectrum analysis method includes preliminary verification experiment of single nuclide, nuclide library instability experiment, and overlapping peak simulation energy spectrum experiment.

[0254] 1. Preliminary verification experiment of single nuclide

[0255] a) The original energy spectrum of the emitted single-energy gamma ray nuclide attenuated by different media:

[0256] The original energy spectrum of the Cs-137 nuclide without attenuation, and the energy spectrum after attenuation by water and sand respectively is shown in the left of Figure 4 The normalized energy spectrum is shown in Figure 4The Compton plateau of the decayed spectrum is much higher than that of the un-decayed spectrum. The three spectra were analyzed by the method proposed in this study. For comparison, an accurate nuclide library was set for GV, and the nuclide analysis results of GV were given. The standard peak areas were obtained by manually selecting the region of interest (ROI) using GV software, and the spectrum analysis results are shown in Figures 5 to 7 , and the quantitative results are shown in Table 1. Among them, Figures 5 to 7 are the spectrum analysis results of Cs-137 (air), Cs-137 (water), and Cs-137 (sand) by different methods, left: standard count of ROI; middle: GV; right: new method.

[0257] Table 1. Spectrum analysis results of Cs-137 by different methods after decayed by different media

[0258]

[0259] From Figures 5 to 7 , it can be found that the method proposed in this study converts the full-energy peak channel number to 1 or 2 channels by deconvolution, which is close to the δ function, so it is very easy to identify the nuclide. As shown in the results of Table 1, since the single-peak nuclide identification task is simple, the results obtained are similar to those of GV, and the errors of the two methods are less than 3%.

[0260] b) Spectrum of multi-energy γ-ray emitting nuclide:

[0261] Since the LILW nuclides generally emit multi-energy γ-rays, the Eu-152 spectrum was analyzed and verified. The original spectrum of Eu-152 is shown in Figure 8 , and 8 rays with the largest peak area were selected for result analysis. Since the energy resolution of HPGe is very good, the overall spectrum is very clear. However, if we zoom in on the details, such as the 778.90 keV full-energy peak in Figure 8 , we can find that the peak width is about 20 channels. The analysis results of the new method for this spectrum are shown in Figure 9 , in which only one channel is used for the rest of the energy rays except for 1112.07 keV, which occupies 2 channels, thus achieving accurate nuclide identification. The comparison of the analysis results of the new method and GV is shown in Table 2. It can be found that for the calculation error of different energy full-energy peaks, the new method is more stable, with a maximum error of 5.9%, while the maximum error of GV is 10.0%. Overall, the average errors of the two methods are less than 4%.

[0262] 2. Nuclide library instability experiment:

[0263] The GV analysis results above are based on an accurate nuclide library. However, in practical applications, it is generally difficult to make a very accurate judgment on the setting of the nuclide library. Therefore, for the same Cs-137 (662 keV) energy spectrum after sand attenuation, different nuclide libraries were used to analyze the GV results. Similarly, results from a new method that does not require setting a nuclide library are presented. Since Ag-110m (658 keV) is the most common corrosion activation product in LILW and its energy is similar to that of Cs-137, the nuclide library settings are as follows: Figure 10 The analysis results are shown in Table 3, where Figure 10 In the middle, from left to right, are: the original energy spectrum of Cs-137; the GV result based on the nuclide library "Cs-137"; the GV result based on the nuclide library "Cs-137+Ag-110m"; and the GV result based on the nuclide library "Ag-110m".

[0264] Table 3 Comparison of energy spectrum analysis results with different nuclide libraries

[0265]

[0266] from Figure 10 It can be observed that for Cs-137, if only Cs-137 is provided in the nuclide library, the analysis results are accurate, with a peak area error of only 2.7%. However, when Ag-110m is added to the nuclide library, according to the energy dispersive spectroscopy (EDS) mechanism, a peak area for Ag-110m will definitely be provided. This means that the 662 keV full-energy peak of Cs-137 is fitted with Gaussian functions located at 658 keV and 662 keV, thus obtaining the area of ​​the two full-energy peaks. Therefore, the peak area error of Cs-137 increases to 32.3%. Furthermore, when Cs-137 is not present in the nuclide library, the report does not provide a count of Cs-137, meaning different nuclide libraries give completely different analysis results. The new method, because it does not require setting a nuclide library, maintains a stable analysis result error of 0.4%.

[0267] 3. Energy spectrum analysis of overlapping peaks of typical LILW nuclides:

[0268] Based on the key nuclide types for LILW provided by Qinshan Nuclear Power Plant, and through comparative analysis of the nuclides and their characteristic energies, three nuclides with the closest characteristic energies were selected: Sb-125 (601 keV, 17.65%), Sb-124 (603 keV, 98.30%), and Cs-134 (605 keV, 97.62%). If the new method can accurately analyze the overlapping peaks composed of these three energies, it can be proven that the new method is fully applicable to LILW measurements.

[0269] Different peak area ratios are needed for the three full-energy peaks. Because it is more difficult to separate the overlapping peaks, it is easier to identify the nuclide as a full-energy peak in the middle channel of the overlapping peak, so if the peak area ratio on both sides of the middle energy is low, a weak peak is formed, and it is easier to verify the reliability of the energy spectrum analysis method. Therefore, for 601 keV, 603 keV, and 605 keV, the peak area ratios of 5:1:10, 10:1:20, 20:1:30, and 50:1:100 are set respectively. At the same time, in order to illustrate the difficulty of overlapping peak recognition, the strongest performance CWT in the traditional peak search method is used for comparison. The energy spectrum analysis results are shown in Figures 11 to 14 , and the quantitative results are shown in Table 4. Specifically, Figures 11 to 14 The peak area ratios are 5:1:10, 10:1:20, 20:1:30, and 50:1:100 respectively. In the four figures, the left side is the original energy spectrum, energy1, energy2, and energy3 represent the energies 601 keV, 603 keV, and 605 keV respectively; the middle is the energy spectrum analysis result of the new method; and the right side is the scale image of CWT peak search.

[0270] Table 4 Energy spectrum analysis results of different overlapping peaks

[0271]

[0272]

[0273] In Figure 11 , it can be found that the new method accurately decomposes the overlapping peak into three incident spectra shaped as δ functions through the precise deconvolution process, and accurately identifies the positions of the three full-energy peaks without energy deviation; and the maximum error of the strong peak is 0.24%, and the error of the weak peak is 1.64%. Using CWT to analyze the overlapping peak, it is found that because the peak area difference is not very large, in the CWT scale image, three full-energy peaks can be found, and the middle one is a weak peak. In Figure 12 , similarly, the positions of the three full-energy peaks are accurately located without energy deviation, and the error of the strong peak is less than 0.5%, and the error of the weak peak is 3.48%.

[0274] In Figure 13 and Figure 14In the middle, through the detail enlargement of the original energy spectrum, it can be seen that because the area ratio of the three full energy peaks is very different, the 603keV full energy peak in the middle is completely submerged by the high count peaks on both sides. In the scale image of CWT, the weak peak cannot be identified, and only the 601keV and 605keV strong peaks can be identified. For the new method, accurate nuclide identification can still be achieved for such extreme overlapping peaks. When the peak area ratio is 20:1:30, the strong peak error is less than 0.5%, and the weak peak error is 5.74%. When the peak area ratio is 50:1:100, the strong peak error is less than 0.5%, and the weak peak error is up to 18.18%. However, because the peak area difference is two orders of magnitude, the activity calculated by the weak peak has little contribution to the total activity of LILW, and does not affect the final disposal of LILW.

[0275] Finally, it needs to be pointed out that, by Figures 11 to 14 The results show that for some extreme overlapping peak conditions, CWT cannot achieve accurate nuclide analysis. However, in the new method, CWT is used because we only use it to find the position of the overlapping peak and perform subsequent peak width calculation, and we do not care what the nuclide is that constitutes the overlapping peak.

[0276] The embodiment of the present application provides a new type of HPGe gamma spectrum analysis method and system, which improves each step of spectrum analysis, realizes the assistance of nuclide library in the overlapping peak decomposition process, improves the measurement accuracy and precision, and provides key technical support for the safe disposal of LILW.

[0277] Those skilled in the art know that, in addition to implementing the system and each device, module and unit thereof provided by the present application in the form of pure computer readable program code, the same function can also be achieved by logically programming the method steps to make the system and each device, module and unit thereof provided by the present application in the form of logic gates, switches, application specific integrated circuits, programmable logic controllers and embedded microcontrollers. Therefore, the system and each device, module and unit thereof provided by the present application can be considered as a hardware component, and the devices, modules and units included therein for achieving various functions can also be considered as structures within the hardware component. The devices, modules and units for achieving various functions can also be considered as both software modules for implementing methods and structures within hardware components.

[0278] The specific embodiments of the present application are described above. It should be understood that the present application is not limited to the above specific embodiments, and those skilled in the art can make various changes or modifications within the scope of the claims, which does not affect the essential content of the present application. The embodiments of the present application and the features in the embodiments can be arbitrarily combined with each other without conflict.

Claims

1. A novel method of HPGe gamma spectrum analysis, characterized by, The method comprises the following steps: Step S1: smoothing the original energy spectrum by a low-pass filter to obtain a smoothed energy spectrum; Step S2: peak searching on the smoothed energy spectrum by a CWT method; Step S3: calculating the peak width by using an improved double convolution method based on the peak position obtained after the peak searching; Step S4: background deduction by a fourth-order filter adaptive SNIP in combination with the peak position and the peak width to obtain a net spectrum; Step S5: calibrating a peak shape filter by an experimental method; Step S5: deconvolution of the net spectrum by using the peak shape filter and solving by a regularization iteration method to obtain a result; The step S3 comprises: The full-energy peak expression is as follows: In the formula, x is a channel address, A is a peak value count, a is a peak position, and sigma is a standard deviation of the full-energy peak; The result of the Gaussian first derivative double convolution calculation on the above formula is as follows: In the formula, delta represents the width of the filter; For the formula, when dC(x) / d delta = 0, delta = sigma is easily obtained; Taking the absolute value of C(x), a matrix M(x, delta) can be obtained by traversing the possible values of delta, wherein x is a channel; The peak position is calculated by using the CWT, and each column curve corresponding to the peak position in the matrix M(x, delta) is analyzed, and there is a maximum value, so the delta value corresponding to the maximum value is the sigma value of the peak position, and [a-3sigma, a+3sigma] is the width of the peak region; The step S5 comprises: The FWHM is calibrated by an experimental method, and since the FWHM at different energies is different, the following formula is used for parameter fitting: In the formula, Energy is energy; The peak shape filter is calculated based on the FWHM information and in combination with a Gaussian function; The method further comprises a step S6: The relationship between the incident spectrum and the energy deposition spectrum is established by using the response function of the detector, as shown in the following formula: In the formula, x(t) is an incident spectrum, y(t) is an energy deposition spectrum, h(t) is a response function, t is a channel address, and n(t) is system noise; the peak-shaped response function h(t) is used as a filter, the incident spectrum is obtained by a deconvolution process, and overlapping peaks are separated; For a discrete system, h(t) is expressed as a matrix: In the formula, each column of the matrix is a probability density function of a full-energy peak response, and the subscript of h indicates a ray emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a-n represent the peak widths of different full-energy peaks; finally, the Boosted-MLEM algorithm is used to solve the formula y(t).

2. The novel HPGe gamma spectrum analysis method according to claim 1, characterized in that, The step S1 comprises: gradually smoothing each channel i in the energy spectrum; For a certain channel x, the corresponding count is y(i), an interval is taken on the left and right sides of x, the counts in the interval are added, and the total counts are L(i) and R(i); The interval is selected as 1.25 FWHM; the total count in the peak region is S, S=y(i)+L(i)+R(i), and the following judgment is made for S: If S > N, then a regular five-point smoothing is performed; if S < N, then no smoothing is performed; where A is a constant, different for different detector parameters. If 1 / r<(R+1) / (L+1)<r, smoothing is performed; wherein r is a constant, and different detectors have different parameters.

3. The novel method of HPGe gamma spectrum analysis according to claim 1, characterized by, The step S2 comprises: Calculating a scale image: where T(E,s) is the scaled image; E is the energy of the ray, s is the scale of the wavelet, ψ E,s (t) is the mother wavelet; the scale s indicates stretching and shrinking of the wavelet as shown in the following equation: The selection of wavelet is the key of the result of spectrum analysis. The full energy peak of gamma spectrum is Gaussian function, and the shape of Mexh is similar to Gaussian peak, so Mexh is selected. The definition of Mexh is as follows: In the formula, ψ(t) represents the Mexh wavelet function; t represents the position of wavelet; σ represents standard deviation, and e is a constant; The position information of the peak is obtained by analyzing the scale image.

4. The novel HPGe gamma spectrum analysis method according to claim 1, characterized in that, The step S4 comprises: The dynamic range of each channel count in the compressed energy spectrum is compressed so as to accelerate iteration, and the LLS operator is applied to the energy spectrum data: where i is the track address, y(i) is the count for each track, v i is the work vector; For the count v(i) of each channel, a parameter m is selected, and then w=2m+1; v1(i), v2(i), and so on are iteratively calculated using second-order filters represented by the following equations in this order m (i): Finally, the v m (i) performing LLS inverse operator operation to obtain the background spectrum; A fourth-order clipping filter is applied in the peak region for background deduction, and the definition of the fourth-order clipping filter is as follows: In the formula, the range of p is 1 to m.

5. A novel HPGe gamma spectrum analysis system characterized by, Comprise: Module M1: smoothing the original energy spectrum through a low-pass filter to obtain a smoothed energy spectrum; Module M2: searching for a peak through the CWT method for the smoothed energy spectrum; Module M3: calculating a peak width by using an improved double convolution method based on the peak position obtained after the peak searching; Module M4: combining the peak position and the peak width, and performing background deduction on the smoothed energy spectrum through a fourth-order filter adaptive SNIP to obtain a net spectrum; Module M5: calibrating a peak shape filter through an experimental method; Module M5: performing deconvolution on the net spectrum by using the peak shape filter, and solving the result through a regularization iteration method; The module M3 comprises: The expression of the full energy peak is as follows: In the formula, x is a channel address, A is a peak value count, a is a peak position, σ is a standard deviation of the full energy peak, and e is a constant; The position information of the peak is obtained by analyzing the scale image; The calculation result of the first-order derivative double convolution of Gaussian is as follows: In the formula, δ represents the width of the filter; For the formula, δ=σ is easily obtained when dC(x) / dδ=0; Taking the absolute value of c(x), a matrix M(x, δ) can be obtained by traversing the possible values of δ, wherein x is a channel; The position of the peak is calculated by using the CWT, and each column curve corresponding to the peak position in the matrix M(x, δ) is analyzed, and there is a maximum value, and the δ value corresponding to the maximum value is the σ value of the peak position, and [a-3σ, a+3σ] is the width of the peak region; The module M5 comprises: The FWHM is calibrated through an experimental method, and since the FWHM at different energies is different, the following formula is used for parameter fitting: In the formula, Energy is energy; The peak shape filter is calculated based on the FWHM information and the Gaussian function; Further comprise module M6: The relationship between the incident spectrum and the energy deposition spectrum is established through the response function of the detector, and is as follows: In the formula, x(τ) is an incident spectrum, y(t) is an energy deposition spectrum, h(t) is a response function, t is a channel address, and n(t) is system noise; the response function h(t) of the peak shape is used as a filter, the incident spectrum is obtained through a deconvolution process, and the separation of overlapping peaks is realized; For a discrete system, h(t) is represented as a matrix: In the formula, each column of the matrix is the probability density function of the full-energy peak response, and the subscript of h represents the ray emitted from the channel; since the FWHM of different full-energy peaks is different, the values of a-n represent the peak widths of different full-energy peaks; finally, the Boosted-MLEM algorithm is used to solve the formula y(t).

6. The novel HPGe gamma spectrum analysis system of claim 5, wherein, The module M1 comprises: step-by-step smoothing for each channel i in the energy spectrum; For a certain channel x, the corresponding count is y(i), and an interval is taken on the left and right sides of x, and the counts in the interval are added, and the total count is L(i) and R(i); The selection of the interval is 1.25 FWHM; let the total count in the peak region be S, then S=y(i)+L(i)+R(i), and the judgment is made for S as follows: If S > N, the conventional five-point smoothing is performed; if S < N, no smoothing is performed; where A is a constant, different for different detector parameters. If 1 / r<(R+1) / (L+1)<r, then smoothing is performed; wherein r is a constant and is different for different detector parameters; The module M2 comprises: Calculation of the scale image: where T(E,s) is the scaled image; E is the energy of the ray, s is the scale of the wavelet, and E,s (t) is the mother wavelet; the scale s indicates stretching and shrinking of the wavelet as shown in the following equation: The selection of the wavelet is the key of the spectral analysis result, the full-energy peak of the gamma spectrum is a Gaussian function, and the Mexh is similar to the Gaussian peak shape, so the Mexh is selected, and the definition of the Mexh is as follows: In the formula, Ψ(t) represents the Mexh wavelet function; t represents the position of the wavelet; and sigma represents the standard deviation; Through the analysis of the scale image, the position information of the peak is obtained; The module M4 comprises: Compressing the dynamic range of the count of each channel in the energy spectrum to accelerate iteration, and applying the LLS operator to the energy spectrum data: where i is the track address, y(i) is the count for each track, v i is the work vector; For the count v(i) of each channel, a parameter m is selected, and then w=2m+1; v1(i), v2(i), and so on are iteratively calculated using second-order filters represented by the following equations in this order m (i): Finally, the v m (i) performing LLS inverse operator operation to obtain the background spectrum; In the peak region, a fourth-order clipping filter is applied for background subtraction, and the definition of the fourth-order clipping filter is as follows: Wherein, the range of p is 1-m.

7. A computer readable storage medium storing a computer program, characterized in that, The computer program is executed by the processor to realize the steps of the novel HPGe gamma spectrum analysis method in any one of claims 1 to 4.

Images