CN109975860B - 一种基于极大似然迭代法的量化核素能谱成分识别方法 - Google Patents

Abstract

本发明公开了一种基于极大似然迭代法的高灵敏度量化核素能谱成分识别方法，其步骤包括：1)生成待识别的核素能谱库，包括N种待识别核素的归一化处理后的能谱以及归一化处理后的环境本底谱；2)对待解析的能谱进行归一化，得到归一化处理后的测试能谱；根据核素能谱库生成一核素含量占比列表；3)将当前核素含量占比列表与核素能谱库做前投影计算，得到计算谱；4)将该测试能谱与当前计算谱相比较，如果两者差别大于设定限定值则继续执行步骤5)，否则输出当前核素含量占比列表；5)将计算谱与测试能谱相比较的结果数组与核素能谱库做反投影过程，得到更新后的核素含量占比列表；6)执行步骤3)～5)设定次数，输出列表。

Description

一种基于极大似然迭代法的量化核素能谱成分识别方法

技术领域

本发明属于核辐射探测及核技术应用领域，特别是涉及一种基于极大似然迭代法的高灵敏度核素能谱成分识别方法。

背景技术

在核辐射探测领域，核素能谱作为放射性核素的特征谱，能够指示出所测得放射性核素的种类，以利于后续对特殊核材料的分析研究，在核辐射监测、核设施退役、核应急、核反恐等领域具有广阔的应用前景。然而实际测得的能谱复杂，本底噪声干扰大，有效计数率低或者包含的核素成分很多，能谱中的核素成分解析和识别的难度很大，因此众多的核素识别算法也应运而生，但是在能谱总计数很低的情况下，即关于高灵敏度量化能谱中所含核素成分的有效占比的相关方法依然缺失。

总体来看，传统的核素识别算法只能给出能谱中所含核素成分的种类，却不能给出每种成分的占比，而且所能分析的能谱总计数很高，核素能谱的解析不能满足高灵敏度的需求。已有的核素识别算法主要包括：1)寻峰法；2)非线性最小二乘拟合法；3)多元线性回归法；4)人工神经网络算法；5)多变量主成分分析法。极大似然迭代法(MLEM)在核素能谱分析中也有应用，已有的文献参考是仅仅用于能谱的去卷积，即滤除核素能谱的探测器响应函数，使得能谱的峰位更加明显，整个过程可以看作是将能量分辨率差，峰位展宽大的能谱滤波为能量分辨好、峰位展宽小的能谱，继而通过寻峰和峰位匹配法给出测试能谱中的放射性核素种类，并没有通过极大似然迭代法(MLEM)直接对测试能谱进行解谱分析。在做核素能谱分析中，本发明涉及到的极大似然迭代法(MLEM)与已有的工作相比较，除了能谱的整个分析流程不同，还有分析方法学上的优势。

发明内容

作为已有核素能谱分析方法学的一种补充和方法学创新，本发明的目的在于，提出一种基于极大似然迭代法的高灵敏度量化核素能谱成分识别方法，为核素能谱分析提供新的技术手段，实现低计数测试能谱中所含核素成分的高灵敏度解析和量化。

本发明的灵感来源于核医学中的基于极大似然迭代法(MLEM)的放射性发射成像方法学(ET)，比如SPECT，如图1所示。核医学发射成像方法学中，放射性药物的分布即热点是要成像所要定位的对象，每个热点发射的射线会在探测器阵列上留下响应，探测器阵列上每个像素的探测事例数的分布也就是当前热点在探测器阵列上响应的概率密度分布。而核素能谱的探测过程与核医学发射成像类似，待识别的核素相当于核医学成像中的热点，每个核素发出的射线谱也都不同，被探测器所探测产生的响应即能谱也就不同，好比每个热点所产生的探测器响应也是不同的。相当于在核医学中是二维的响应，过渡到核素识别中是一维的能谱响应，能谱中每个通道计数的不同代表着响应概率密度的不同。由此可以看出，每种核素产生的探测器响应即一维概率密度分布的能谱，类比于核医学成像中每个热点的探测器响应即二维概率密度分布的投影。基于极大似然迭代法(MLEM)核医学发射成像的方法自1974年至今已经很成熟，是置信度很高的放射性药物热点分布重建方法，本发明采用的基于极大似然迭代法的高灵敏度量化核素能谱成分解析与识别方法便是借鉴于此，每种核素的有无和含量占比的多少好比于热点的有无和强弱。本发明的基于极大似然迭代法的高灵敏度量化核素能谱成分解析与识别方法便是上文中核医学发射成像方法学的衍生。

本发明提出了一种基于极大似然迭代法的高灵敏度量化核素能谱成分识别方法，该方法的流程图如图2所示，具体实施步骤包括：

1.备好核素识别探头，结合数字化能谱采集系统，用于测量放射性核素的能谱。

2.确定好待识别核素的种类，IAEA-2006标准规定核素识别仪能够有效识别的核素包括24种：¹³³Ba,²⁴¹Am,⁵⁷Co,²²Na,¹³⁷Cs,⁶⁰Co,⁶⁷Ga,¹⁹²Ir,¹⁵²Eu,¹³¹I,⁴⁰K,⁷⁵Se,¹⁸F,¹¹¹In,¹³⁸La,⁵¹Cr,¹⁰³Pd,¹³³Xe,²²⁶Ra,²³⁷Np,²³⁷Pu,²³²Th,²³³U和²³⁵U，涉及到特殊核材料、医用核素、工业用核素和天然核素(后续如有新的核素，只需要更新核素库即可，本方法同样适用)。确定好待识别核素的种类后需要获得当前核素识别探头关于上述所提及核素的能谱，并且对每种核素的能谱做归一化处理，加上归一化后的环境本底谱，即生成待识别的核素能谱库。核素能谱库的具体生成过程如下：

2.1用核素识别探头测量获得(实验室)已有的若干核素的能谱，比如以下6种：¹³³Ba,²⁴¹Am,⁵⁷Co,²²Na,¹³⁷Cs和⁶⁰Co；

2.2对测量获得的几种核素能谱的全能峰进行峰位和半高全宽的拟合，拟合得到的参数为当前核素识别探头的能量刻度和能量分辨率刻度系数，即当前探测器能谱响应下能量与道址的对应关系以及能量分辨率与能量的对应关系，关于核素能谱的拟合和刻度是核素分析必不可少的过程且相关技术很成熟，此处不再赘述；

2.3根据当前核素识别探头的能量刻度和能量分辨率刻度系数，通过Geant4模拟给出其他18种核素的能谱，与实验测量获得的6种核素能谱以及加上归一化后的环境本底谱一并合为待识别核素能谱库，其中的每种核素能谱都需要做归一化处理，即能谱中每通道的计数除以能谱总计数，归一化后的能谱中的每一通道对应的是概率密度。

3.在准备好核素能谱库之后，便可以对实际所测得的测试能谱执行本发明所设计的核素识别算法，即基于极大似然迭代法的高灵敏度量化核素能谱成分识别方法，具体的执行步骤如图2所示：

3.1准备好待解析的归一化处理后的测试能谱，开始执行迭代过程；

3.2准备好初始化列表，是一个包含有25个元素值的一维数组，每个值代表着所属核素在测试能谱中的百分占比，25个元素代表着当前测试环境下的本底谱和24种上文提及的核素谱，初始化列表中的每个值都为1/25。

3.3将当前列表与核素能谱库做前投影计算，得到计算谱；

3.4将测试能谱与计算谱相比较，两者的比较过程在此处做一些额外的说明，在主流的迭代算法中，以SART(联合代数迭代重建算法)为例，测试对象与计算对象的比较过程为两者各个元素的相减，以两者间的差值作为迭代反投影过程中对测试对象的修正值；而在MLEM(极大似然期望最大化估计迭代重建算法)中，测试对象与计算对象的比较过程为两者各个元素的相除，以两者间的比值作为迭代反投影过程中对测试对象的修正值，它的好处在于保证先验性的设定待重建对象的非负性。所以，本方法中测试能谱与计算谱相比较的结果为两者对应数组的相除比值结果。如果两者差别大于设定限定值则继续执行下一环节，其中，两个测试能谱与计算谱的差别可以表述为两者每通道的差值与测试能谱对应通道的比值，所有通道平均后的差值比值结果即可以表征两者的近似度，差值比值结果越小表示两者近似度越高；

3.5将计算谱与测试能谱相比较的结果数组与核素能谱库做反投影过程，得到更新后的当前列表，此时完成一次迭代过程；

3.6如果当前列表不作为最终的输出结果，则继续使用当前列表替代3.2中的列表，逐一执行3.3、3.4和3.5，直到3.4中利用当前列表的结果计算出的计算谱与测试谱差别小于限定值(根据经验值，通常限定值设为尽可能的小，比如0.01，表示待达到的相似度为99％)，或者迭代次数达到一定数量级自行停止，本方法鉴于计算机的计算能力，通常选择迭代次数达到为百次或者千次时自行停止，该数量级下迭代结果已经足够收敛，并且运算时间为百毫秒或者几秒量级；

3.7作为最终输出结果的当前列表即为整个方法学的执行结果，即输出列表，列表中的每个元素的值代表着每种核素在测试谱中的含量占比。

与现有技术相比，本发明的积极效果为：

与现有技术方案相比，本发明能够对总计数低的测试能谱进行解析和核素识别，实现了高灵敏度解谱，并且对测试能谱所含核素成分的含量占比实现了量化。本发明的方法能够：1)对总计数低的能谱进行解析，可以对总计数为100-1000量级的能谱进行量化解析，而传统核素识别方法所需要的能谱总计数为10000量级；2)实现测试能谱中所含核素成分的量化，即在能谱总计数中的百分占比；3)不需要通过寻峰法来进行核素匹配识别，有效地利用了整个能谱的形状信息；4)不需要做本底的扣除，本底也作为待识别测试能谱中的有效成分，在解谱的结果即输出列表中，包含有本底的成分，极大的简化了算法执行过程。

本发明的方法为核素能谱分析领域提供了一种补充和方法学创新，以利于后续对特殊核材料的分析研究，在核辐射监测、核设施退役、核应急、核反恐等领域具有广阔的应用前景。

附图说明

图1为方法流程对比图；

(a)为基于极大似然迭代法(MLEM)的放射性发射成像方法学流程图，示意核医学成像中某位置热点的探测器投影响应；

(b)为本发明基于极大似然迭代法(MLEM)的类比流程图，示意核素能谱解析中某种核素的探测器投影响应；

图2为本发明的具体实现方法流程图；

图3为核素能谱库图；

图4为本方法的实验测试能谱解析结果：

(a)为待分析的包含100个计数的⁶⁰Co测试能谱；

(b)为对待分析测试能谱(a)迭代100次后的解析结果。

具体实施方式

本发明提出一种基于极大似然迭代法的高灵敏度核素能谱成分识别方法，实施步骤如下：

1.备好核素识别探头，结合数字化能谱采集系统，用于测量放射性核素的能谱。其中，核素识别探头是由基于闪烁体(比如Ce:GAGG，BGO，YSO，LYSO，CsI等)耦合光子探测器(比如硅光电倍增管Silicon Photomultipliers，SiPM；光电倍增管PhotomultiplierTube，PMT)形成的射线探测单元，数字化能谱采集系统为数字电子学采集板卡或者数字化多道采集系统。

2.确定好待识别核素的种类，IAEA-2006标准规定核素识别仪能够有效识别的核素包括24种：¹³³Ba,²⁴¹Am,⁵⁷Co,²²Na,¹³⁷Cs,⁶⁰Co,⁶⁷Ga,¹⁹²Ir,¹⁵²Eu,¹³¹I,⁴⁰K,⁷⁵Se,¹⁸F,¹¹¹In,¹³⁸La,⁵¹Cr,¹⁰³Pd,¹³³Xe,²²⁶Ra,²³⁷Np,²³⁷Pu,²³²Th,²³³U和²³⁵U，涉及到特殊核材料、医用核素、工业用核素和天然核素，本发明提出的核素能谱分析方法限于识别分析上述提及的核素。确定好待识别核素的种类后需要获得当前核素识别探头关于上述所提及核素的能谱，并且对每种核素的能谱做归一化处理，加上归一化后的环境本底谱，即生成待识别的核素能谱库。核素能谱库的具体生成过程如下：

2.1用核素识别探头测量获得(实验室)已有的几种核素的能谱，包括6种：¹³³Ba,²⁴¹Am,⁵⁷Co,²²Na,¹³⁷Cs和⁶⁰Co；

2.2对测量获得的几种核素能谱的全能峰进行峰位和半高全宽的拟合，拟合的参数即为当前核素识别探头的能量刻度和能量分辨率刻度系数；

2.3根据当前核素识别探头的能量刻度和能量分辨率刻度系数，通过Geant4模拟给出其他18种核素的能谱，与实验测量获得的6种核素能谱一并合为待识别核素能谱库，其中的每种核素能谱都需要做归一化处理，即能谱中每通道的计数除以能谱总计数，归一化后的能谱中的每一通道对应的是概率密度，每个核素能谱的通道数为2048。

2.4最终生成的核素能谱库如图3所示，用数学式表达为：

因此在下文本发明的核素能谱解析算法执行过程中，用T来表示核素能谱库，t_ij表示第i种核素的第j通道的概率密度值，S的数值为25，表示总共有25种核素(包含了本底)，C的数值为2048，表示每个核素能谱的通道数总共有2048道。t_ij的数值范围为(0,1),第i种核素的能谱总和，即

的值为1，此为归一化后的结果；

3.在准备好核素能谱库之后，执行本发明所设计的核素能谱解析方法还需要归一化处理后的测试能谱和核素含量占比列表两个输入值：

3.1将待解析的测试能谱归一化处理为Y＝{y_i|i＝1,2,3……C}，y_i表示当前待分析的测试能谱的第i通道的概率密度值，待解析的测试能谱中各通道计数是已知的，即y_i等于测试能谱中第i通道的计数值除以能谱总计数；

3.2准备好初始化列表，是一个包含有25个元素值的一维数组，每个值代表着所属核素在测试能谱中的百分占比，25个元素代表着本底和24种上文提及的核素，初始化列表中的每个值都为1/25，核素含量占比列表的数学表达式为：P＝{p_i|i＝1,2,3……S}，p_i表示第i种核素在测试能谱中的含量占比；

4.在准备好核素能谱库T、归一化处理后的测试能谱Y、核素含量占比列表P之后，便可以对实际所测得的测试能谱执行本发明所设计的核素能谱解析方法，即基于极大似然迭代法的高灵敏度量化核素能谱成分解析与识别方法，算法执行过程的数学表达式为：

其中P^new表示当前迭代过程执行完一次后更新的核素含量占比列表，P^old表示上一次迭代执行完后得到的核素含量占比列表，Y表示待解析的归一化处理后的测试能谱，T表示核素能谱库，算法执行过程中的每一步说明为：

其中Y’为当前核素含量占比列表P^old与核素能谱库T的前投影卷积计算后的计算谱；

其中D为归一化处理后的测试能谱Y与计算谱Y’的比较结果，在极大似然法迭代中，待测试对象与前投影计算对象的比较通常为两者的相除，此处即为测试能谱对计算谱中每个i通道概率密度值的修正；

其中A为谱比较结果D和核素能谱库的反投影卷积计算结果；

4.4 P^new＝P^old·A＝{p_i ^new|i＝1,2,3......S}，p_i ^new＝p_i ^old·a_i， (5)

其中P^new表示当前迭代过程执行完一次后更新的核素含量占比列表，P^old表示当前迭代执行过程中的核素含量占比列表，P^old与反投影卷积计算结果的点乘结果即为此次迭代更新后的核素含量占比列表P^new；

5.算法执行过程公式(1)每计算完一次，即迭代一次，将会得到更新后的P^new，当迭代次数足够多时，P^new和当前P^old的差别将会小到没有变化，此时停止迭代过程。

6.作为示范事例，图4表示待分析测试能谱的核素为⁶⁰Co,总计数为100，迭代100次后的解析结果显示解析识别出的⁶⁰Co含量占总测试能谱的86.55％，即量化识别精度为86.55％，实现了高灵敏度量化核素能谱成分解析与识别。

综上所述，以上仅为本发明的较佳实施方法，并非用于限定本发明的保护范围。凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种基于极大似然迭代法的量化核素能谱成分识别方法，其步骤包括：

1)生成待识别的核素能谱库，包括N种待识别核素的归一化处理后的能谱以及归一化处理后的环境本底谱；

2)对待解析的能谱进行归一化，得到归一化处理后的测试能谱；根据核素能谱库生成一核素含量占比列表，包括N+1个值，分别代表该N种核素以及环境本底谱在该测试能谱中的百分占比；

3)将当前核素含量占比列表与核素能谱库做前投影计算，得到计算谱；

4)将该测试能谱与当前计算谱相比较，如果两者差别大于设定限定值则继续执行步骤5)，否则执行步骤7)；

5)将计算谱与测试能谱相比较的结果数组与核素能谱库做反投影过程，得到更新后的核素含量占比列表，完成一次迭代过程；

6)执行步骤3)～5)设定次数，输出最后一次迭代得到的核素含量占比列表；

7)将当前核素含量占比列表输出，核素含量占比列表中的每个元素的值代表着对应核素在该测试能谱中的含量占比。

2.如权利要求1所述的方法，其特征在于，生成所述核素能谱库的方法为：

a)用核素识别探头测量获得该待识别核素中若干核素的能谱；

b)对测量获得的每种核素能谱的全能峰进行峰位和半高全宽的拟合，拟合得到的参数为当前核素识别探头的能量刻度和能量分辨率刻度系数；

c)根据当前核素识别探头的能量刻度和能量分辨率刻度系数，模拟给出其他待识别核素的能谱，将模拟能谱与已测量获得的核素能谱进行归一化，并与归一化后的环境本底谱合并，生成该待识别核素能谱库。

3.如权利要求1所述的方法，其特征在于，将该测试能谱与当前计算谱相比较的方法为：取测试能谱第i通道与计算谱第i通道的差值Δ_i；或者取测试能谱第i通道与计算谱第i通道的比值d_i。

4.如权利要求3所述的方法，其特征在于，将所有通道比较结果的平均值与设定限定值进行比较，如果该平均值大于设定限定值则继续执行步骤5)。

5.如权利要求4所述的方法，其特征在于，核素能谱库

其中的元素t_ij表示第i种核素的第j通道的概率密度值，S＝N+1，C表示每个核素能谱的通道数，

的值为1。

6.如权利要求5所述的方法，其特征在于，每次迭代更新核素含量占比列表的方法为：

其中，P^new表示当前迭代过程执行后更新的核素含量占比列表，P^old表示上一次迭代执行完后得到的核素含量占比列表，Y表示待解析的归一化处理后的测试能谱，T表示核素能谱库，

为卷积计算。

7.如权利要求6所述的方法，其特征在于，所述计算谱

8.如权利要求7所述的方法，其特征在于，所述结果数组为

9.如权利要求8所述的方法，其特征在于，结果数组与核素能谱库做反投影计算的结果

10.如权利要求1所述的方法，其特征在于，初始化列表中的每个值均相同。

Images