一种利用康普顿散射事例统计进行放射性核素识别的方法与流程

本发明属于核辐射探测及核技术应用领域，特别是涉及一种利用康普顿散射事例统计进行放射性核素识别的方法。

背景技术：

在核科学技术大力发展的今天，对放射性核素的定性识别与定量测定是核辐射探测的核心任务之一。具有识别放射性核素功能的伽玛谱仪、伽玛相机等设备，在军控核查、核设施监控、核安检与核应急等领域发挥着重要的作用。放射性核素识别的方法一直以来是伽玛能谱和伽玛成像领域的研究热点，其主要原理是利用不同放射性核素的γ能谱具有不同的特征峰，通过对采集的γ能谱进行处理，识别出对应的放射性核素。目前公知的核素识别方法主要包括：1)特征峰匹配法；2)多元线性回归法；3)多变量主成分分析法；4)人工神经网络算法等。这些方法都没有脱离γ能谱这一先决条件，需要高能量分辨的探测器获取完整的能谱数据来进行放射性核素识别。

伽玛成像技术相比于剂量仪、伽玛谱仪等非成像式探测技术，能够在远距离给出核辐射热点分布情况的二维图像，可大幅降低操作人员所受照射剂量。康普顿成像是目前伽玛成像技术中新兴发展并展现出广泛应用的一种技术，基于入射γ光子的康普顿散射原理，通过测量散射点、吸收点的位置以及所沉积的能量，回溯重建出γ光子来源的方向，具有视野范围宽、能量范围广、体积紧凑等优势。应用康普顿成像技术的设备通常称为康普顿相机。康普顿相机在对辐射热点进行重建定位时需要已知入射射线的能量，即依赖于识别辐射环境中的放射性核素。目前康普顿相机中采用的核素识别方法仍是基于获取的γ能谱，分别进行核素成像与核素识别两个环节，没有有效的利用收集康普顿散射事例的过程，存在核素成像和核素识别不同步、需要增加高能量分辨的探测器等问题。本发明方法正是基于康普顿相机对放射性核素识别的特殊需求而提出的。

技术实现要素：

为了解决康普顿相机核素成像与核素识别不同步的问题，本发明的目的在于，提出一种利用康普顿散射事例统计进行放射性核素识别的方法，实现康普顿相机对未知放射性核素的量化识别。

区别于伽玛谱仪以及其他类型的伽玛相机，康普顿相机的主要特征在于其利用康普顿散射事例进行成像，且需要已知入射γ射线的能量，即内置包含待识别放射性核素特征能量的核素库。考虑到其成像过程依赖于获取核素库中满足各核素条件的康普顿散射事例，因此本发明利用对各核素康普顿散射事例的统计信息，实现对未知放射性核素的识别，并通过极大似然期望最大化算法(mlem)，计算得到识别出的核素量化组分率。mlem算法是一种在核成像领域公知的迭代算法，具有成熟度较高的应用。本发明选用mlem算法作为计算核素量化组分率的方法，是因为本发明方法中获取的康普顿散射事例统计信息与待识别核素一一对应，这与mlem算法中系统响应矩阵(tij)的设定相符。本发明方法有效解决了康普顿相机核素成像与核素识别不同步的问题，利用成像过程中收集的康普顿散射事例的统计信息，经mlem算法计算核素量化组分率，实现对未知放射性核素的识别。

本发明方法首先利用康普顿相机所确定的康普顿散射事例统计数组特征，测量生成该康普顿相机的待识别核素库；应用时，通过测量包含在待识别核素库中的未知核素康普顿散射事例统计数组，依据极大似然期望最大化法得到待测核素的种类及量化组分率。该方法的基本流程图如图1所示，具体步骤包括：

1.依据康普顿相机的能量分辨率、散射能谱范围、以及成像架构等因素，确定康普顿散射事例统计数组的特征，所述康普顿散射事例统计数组中每个数组元素表示一种待识别核素的康普顿散射事例计数率，所述特征包括可识别放射性核素的种类、各核素的判断条件；

2.依据步骤1中该康普顿相机可识别放射性核素的种类，确定其核素库可识别核素的个数s(包含本底，将本底看做一种待识别核素情况)，以j(j＝{0,1,2,…,s-1})表示第j种待识别核素，则康普顿散射事例统计数组包含的数组元素个数为s-1，以i(i＝{0,1,2,…,s-2})表示第i个数组元素。依据步骤1中该康普顿相机各核素的判断条件，测量生成康普顿相机的待识别核素库，具体步骤如下：

2.1测量康普顿散射事例统计数组在本底情况下的各数组元素计数率，记第i元素计数率为本底计数率bi；

2.2分别测量康普顿散射事例统计数组在每种待识别核素单独存在情况下(此时不需再测量本底情况)的各数组元素计数率，记每种核素测量得到的数组中第i元素计数率为ni；

2.3对于每种待识别核素，将步骤2.2中的ni分别除以步骤2.1中的bi，并经归一化运算，即得到康普顿散射事例统计数组中每一元素计数率相对于本底情况的倍数的概率密度；

2.4将步骤2.3的运算结果表示为系统响应矩阵，即得到该康普顿相机的待识别核素库；系统响应矩阵中的元素tij(i＝{0,1,2,…,s-2}，j＝{0,1,2,…,s-1})，其物理意义为在第j种核素单独存在时康普顿散射事例统计数组第i元素的概率密度。

3.测量未知核素(包含在待识别核素库的种类中)的康普顿散射事例统计数组，依据极大似然期望最大法(mlem)得到待测核素的种类及量化组分率，具体步骤如下：

3.1对包含在待识别核素库中的未知核素进行测试，得到待测康普顿散射事例统计数组，记待测康普顿散射事例统计数组中第i元素计数率为mi；

3.2将步骤3.1中的mi分别除以步骤2.1中的本底计数率bi，经归一化运算后，得到待测康普顿散射事例统计数组中每一元素计数率相对于本底情况的倍数的概率密度，记为

3.3依据极大似然期望最大化法(mlem)对核素组分数组进行迭代，得到待测核素的种类及量化组分率，先初始化核素组分数组；然后按照步骤3.3.2至3.3.5的循环执行，每循环一遍即完成一次迭代；每次迭代后，核素组分数组完成更新。mlem算法的执行过程分为以下步骤：

3.3.1根据步骤2中康普顿相机可识别核素的个数s，初始化核素组分数组；初始的核素组分数组中的第j个数组元素λj＝1/s(j＝{0,1,2,…,s-1})，表示第j种核素在核素组分数组中的初始含量占比为1/s；

3.3.2计算核素组分数组与步骤2.4中的待识别核素库的前投影卷积，前投影卷积中的第i个元素(i＝{0,1,2,…,s-2})，表示为此为计算期望过程；

3.3.3将步骤3.2中的待测康普顿散射事例统计数组每一元素概率密度vi分别除以步骤3.3.2中的前投影卷积中的每一元素后，与步骤2.4中的待识别核素库进行反投影卷积运算，反投影卷积中的第j个元素pj(j＝{0,1,2,…,s-1})，表示为

3.3.4计算当前核素组分数组与步骤3.3.3中的反投影卷积的点乘结果，得到此次迭代更新的核素组分数组，迭代更新的核素组分数组中的数组元素λj^new(j＝{0,1,2,…,s-1})，表示为λj^new＝λj^old·pj，其中，λj^old为当前核素组分数组的第j个元素，此为计算期望的最大化过程；

3.3.5循环执行步骤3.3.2至3.3.4，当步骤3.3.2中的前投影卷积与步骤3.2中的以概率密度表示的待测康普顿散射事例统计数组的差别小于限定值(通常为经验值，反映由系统响应矩阵计算的前投影卷积与实际测试的待测康普顿散射事例统计数组之间的相似度)，停止循环迭代，此时的核素组分数组中的λj^new反映了第j种核素占总核素组分的概率，从概率可判断待测核素的种类，并得到识别精度。

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

本发明方法充分发挥了康普顿相机的优势，利用其需要已知待识别核素的能量、和按照待识别核素将康普顿散射事例进行分类存储进而成像的特性，发掘了一种利用康普顿散射事例统计进行放射性核素识别的方法。与康普顿相机中现有的依据γ能谱进行核素识别的方法相比，本发明实现了康普顿相机核素识别与核素成像的同步，不需要依赖高能量分辨的探测器采集γ能谱，直接在康普顿相机获取散射事例数据的过程中实现放射性核素识别的功能，并能给出识别核素的量化组分率，提高了康普顿相机核素识别的效率。本发明方法适用于能够获取康普顿散射事例的设备，主要应用于康普顿相机，是对现有技术的补充和方法学创新。

附图说明

图1为本发明的方法流程图。

图2为本发明实施例所示的双层探测器架构康普顿相机。

图3为本发明实施例所示的康普顿散射事例统计数组来源的说明示意图。

图4为本发明实施例所示的康普顿相机待识别核素库。

图5为一组待测康普顿散射事例统计数组各元素概率密度的分布图。

图6为对图5利用mlem算法计算得到的核素组分数组。

具体实施方式

本发明提出一种利用康普顿散射事例统计进行放射性核素识别的方法，实现康普顿相机对未知放射性核素的量化识别。该方法首先利用康普顿相机所确定的康普顿散射事例统计数组特征，测量生成该康普顿相机的待识别核素库；应用时，通过测量包含在待识别核素库中的未知核素康普顿散射事例统计数组，依据极大似然期望最大化法得到待测核素的种类及量化组分率。

本实施例以双层探测器架构的康普顿相机为例，下面结合附图对本发明作进一步的阐述。

步骤s1：依据康普顿相机的能量分辨率、散射能谱范围、以及成像架构等因素，确定康普顿散射事例统计数组的特征(可识别放射性核素的种类、各核素的判断条件)；

需要说明的是，康普顿散射事例统计数组的特征是根据不同康普顿相机成像能力来确定的，与康普顿相机对其内置核素库中的放射性核素成像所依据的分类条件相一致，即康普顿散射事例统计数组的特征来源于不同类型康普顿相机的核素成像判据。

在本实施例中，所示的双层探测器架构康普顿相机，是康普顿相机的最基本成像架构，包含散射层探测器和吸收层探测器，每层探测器均为位置灵敏探测器，能够探测γ射线作用的位置和能量，如图2所示。因康普顿散射发生在中高能区(几百kev～几mev)，在实施例中，以在300kev～1.5mev能量区间内实验室常见的四种放射性核素¹³³ba(钡133)、²²na(钠22)、¹³⁷cs(铯137)和⁶⁰co(钴60)进行说明。本发明方法所适用的放射性核素，包括但不限于上述4种核素。

在本实施例中，由代表上述四种核素的4个数组元素构成康普顿散射事例统计数组，每个数组元素(n[0],n[1],n[2],n[3])存储一种待识别核素的康普顿散射事例计数率。

所述康普顿散射事例统计数组的特征，包括了康普顿相机可识别放射性核素的种类、以及各放射性核素的判断条件。

在本实施例中，所述可识别放射性核素的种类包含了康普顿散射事例统计数组元素所对应的4种核素(¹³³ba、²²na、¹³⁷cs、⁶⁰co)以及本底，共5种，因此实施例的待识别核素库可识别核素的个数s＝5。所述各放射性核素的判断条件是依据康普顿相机的能量分辨率、散射能谱范围以及成像架构等因素确定的，具体的判断方法为本领域研究人员所熟知，此处不再累述。

表1列出了本实施例所示的康普顿散射事例统计数组及其特征，其中n[2]所存储的为对应于¹³⁷cs核素的康普顿散射事例计数率。在实施例中，n[2]记录的康普顿散射事例为康普顿相机探测器测量数据中满足622kev<e1+e2<702kev且10kev<e1<165kev的事例，这部分事例展示在图3中e1e2能量分布图的实线四边形所圈区域中，其中e1为散射层探测器沉积能量，e2为吸收层探测器沉积能量。

表1为康普顿散射事例统计数组及其特征

步骤s2：依据表1所示康普顿散射事例统计数组的特征，测量生成实施例中康普顿相机的待识别核素库。在实施例中，待识别核素库可识别核素的个数s＝5，以j(j＝{0,1,2,3,4})表示第j种待识别核素。康普顿散射事例统计数组包含的数组元素个数为4，以i(i＝{0,1,2,3})表示第i个数组元素。测量生成待识别核素库的具体步骤如下：

2.1测量康普顿散射事例统计数组在本底情况下的各数组元素计数率，记第i元素计数率为本底计数率bi；

2.2分别测量康普顿散射事例统计数组在每种待识别核素单独存在情况下(此时不需再测量本底情况)的各数组元素计数率，记每种核素测量得到的数组中第i元素计数率为ni；

2.3对于每种待识别核素，将步骤2.2中的ni分别除以步骤2.1中的bi，并经归一化运算，即得到康普顿散射事例统计数组中每一元素计数率相对于本底情况的倍数的概率密度；

2.4将步骤2.3的运算结果表示为系统响应矩阵，即得到该康普顿相机的待识别核素库；系统响应矩阵中的元素tij(i＝{0,1,2,3}，j＝{0,1,2,3,4})，其物理意义为在第j种核素单独存在时康普顿散射事例统计数组第i元素的概率密度，如图4所示；

步骤s3：测量未知核素(包含在待识别核素库的种类中)的康普顿散射事例统计数组，依据极大似然期望最大法(mlem)得到待测核素的种类及量化组分率，具体步骤如下：

3.1对包含在待识别核素库中的未知核素进行测试，得到待测康普顿散射事例统计数组，记数组中第i元素计数率为mi；

3.2将步骤3.1中的mi分别除以步骤2.1中的本底计数率bi，经归一化运算后，得到待测康普顿散射事例统计数组中每一元素计数率相对于本底情况的倍数的概率密度，记为如图5所示；

3.3依据极大似然期望最大化法(mlem)对核素组分数组进行迭代，得到待测核素的种类及量化组分率，先初始化核素组分数组；然后按照步骤3.3.2至3.3.5的循环执行，每循环一遍即完成一次迭代；每次迭代后，核素组分数组完成更新。mlem算法的执行过程分为以下步骤：

3.3.1根据实施例步骤s1中康普顿相机可识别核素的个数s＝5，初始化核素组分数组，初始的核素组分数组中的数组元素λj＝1/s＝0.2(j＝{0,1,2,3,4})，表示5种可识别核素情况中(本底、¹³³ba、²²na、¹³⁷cs、⁶⁰co)第j种核素在核素组分数组中的初始含量占比为0.2；

3.3.2计算核素组分数组与步骤2.4中的待识别核素库的前投影卷积，前投影卷积中的第i个元素(i＝{0,1,2,3})，表示为此为计算期望过程；

3.3.3将步骤3.2中的待测康普顿散射事例统计数组每一元素概率密度vi分别除以步骤3.3.2中的前投影卷积中的每一元素后，与步骤2.4中的待识别核素库进行反投影卷积运算，反投影卷积中的第j个元素pj(j＝{0,1,2,3,4})，表示为

3.3.4计算当前核素组分数组与步骤3.3.3中的反投影卷积的点乘结果，得到此次迭代更新的核素组分数组，迭代更新的核素组分数组中的数组元素λj^new(j＝{0,1,2,3,4})，表示为λj^new＝λj^old·pj，其中，λj^old为当前核素组分数组的第j个元素，此为计算期望的最大化过程；

3.3.5循环执行步骤3.3.2至3.3.4，当步骤3.3.2中的前投影卷积与步骤3.2中的以概率密度表示的待测康普顿散射事例统计数组的差别小于限定值(通常为经验值，反映由系统响应矩阵计算的前投影卷积与实际测试的待测康普顿散射事例统计数组之间的相似度)，停止循环迭代，此时的核素组分数组中的λj^new反映了第j种核素占总核素组分的概率，从概率可判断待测核素的种类，并得到识别精度。

在本实施例中，假定取限定值为0.1，则当估计值(由系统响应矩阵计算的前投影卷积每一元素)与真实值(实际测试的以概率密度表示的待测康普顿散射事例统计数组每一元素vi)之间满足(i＝{0,1,2,3})，表示两者之间的相似度>90％，此时停止循环迭代，输出核素组分数组。

在本实施例中，图5显示了一组对¹³⁷cs核素进行实测后的待测康普顿散射事例统计数组各元素概率密度分布图，其纵轴表示步骤3.2中经归一化后的数组各元素相对于本底情况的倍数的概率密度。图6为利用步骤3.3所述的mlem算法在限定值为0.1时(此实施例中的循环迭代次数为10次)输出的核素组分数组，显示了实施例中的5种待识别核素情况(本底、¹³³ba、²²na、¹³⁷cs、⁶⁰co)的含量占比，其中¹³⁷cs组分含量占总核素组分的比例最高，为95％，表示识别得到的放射性核素为¹³⁷cs，量化的识别精度为95％。

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