基于能谱重建技术的全身计数器γ能谱分析

赵 日，刘立业，曹勤剑

(中国辐射防护研究院，山西 太原 030006)

全身计数器是专用于测量人体体内放射性的仪器。它从体外直接探测体内发出的γ射线，通过分析γ射线能谱实现对放射性核素种类和含量的定量测量。全身计数器在辐射防护、核安全等领域均有重要应用[1]。

目前，全身计数器的γ能谱分析最常用的方法是全能峰法[2-5]。该方法先用滤波算法寻找能谱中由光电效应形成的全能峰，然后对每个峰进行曲线拟合以求得其净计数(峰面积)，最后，结合能量和效率刻度函数确定测量对象发出的γ射线的能量和数量。然而，由于全身计数器通常采用NaI(Tl)探测器以提高探测效率，而该类型探测器能量分辨率差，且存在温漂等现象，再加上人体对γ射线的散射以及天然40K干扰等因素，全身计数器获取的γ能谱常呈现全能峰辨识度低、重峰现象严重、峰下侧基底高、有假峰等不利性态[6]。对于这种γ能谱，全能峰法得到的结果通常准确性很低：其在寻峰时极易误识别和漏识别，峰拟合时因重峰和高基底的原因而误差较大。如何提高全身计数器γ能谱分析精度一直是长期受关注的课题。

近期，一种基于能谱重建算法的分析技术在低分辨率γ能谱中得到越来越多的应用。该技术首先建立核素出射谱到探测器测量谱直接的测量矩阵，然后使用特殊迭代算法从测量谱反解出出射谱。由于出射谱在全谱范围内只在若干能量点处有幅值，这些能量点与γ射线能量一一对应，而幅值等于测量时间内该γ射线总发射数，因此，理论上可直接从反解结果中读取γ射线能量和数量，以此实现对放射性核素的定量分析。该技术在NaI探测器解谱领域的应用始于20世纪末：1997年Bandzuch等[7]对正电子湮灭谱使用了能谱重建技术；2000年Meng等[8]在对7.62 cm×7.62 cm NaI(Tl)探测器的γ能谱进行分析时比较了3种不同重建算法的效果。2004年Jandel等[9]对152Eu和56Co的测量谱进行了能谱重建，结果显示，全能峰被压缩入孤立的道，而基底部分被全部移除；2009年Rahman等[10]则将该技术应用于核电厂流出物监测系统(使用NaI(Tl)探测器)所获取的γ能谱，所得到的重建能谱中峰康比大幅提高；2010年Hanka[11]提出了基于动态指标集和梯度法的γ能谱重建算法；Morhac等[12]综述了目前常用的若干种γ能谱重建算法。本文为提高全身计数器γ能谱分析精度，开展基于能谱重建技术的γ能谱分析方法研究。

1 原理和实验

1.1 能谱重建技术原理

假设核素出射的γ能谱为x，探测器测得的γ能谱为y，则理论上有：

(1)

其中，r(E′→E)表示出射谱中能量为E′的单位数量γ射线在测量谱E处的相应计数。

当能谱计数率较高时，会出现符合相加效应和死时间现象，式(1)并不严格成立；不过对于全身计数器测量，人体体内放射性一般非常弱，加上本底辐射是天然环境辐射，此时γ能谱计数率较低，符合相加效应和死时间现象几乎可忽略，因此可认为式(1)是γ射线出射谱和探测谱之间映射关系的合理描述。

将式(1)写为离散形式，并考虑到各道计数存在统计涨落，则有下式：

y=Rx+ε

(2)

其中：R为测量矩阵，其元素Rij表示出射谱中单位数量的第i道γ射线在测量谱第j道处的相应计数；ε为由统计涨落带来的误差。

图1 直接矩阵求逆所得的解Fig.1 Solution from direct matrix invert method

由于测量矩阵R高度病态，且式中存在误差项，式(2)若用一般的矩阵求逆算法求解，只能得到剧烈振荡的无意义解，如图1所示。为从式(2)出发，稳定求解出x且保证x的所有分量均为非负(能谱幅值不可能是负值)，提出了若干精细算法，其中最常用、最有效的主要有两种，即Gold算法[12-13]和EM算法[14-16]。两者均属于迭代算法，且均具有全局收敛性和解非负性的优点，但Gold算法较EM算法收敛更快。本文将采用基于Gold算法的能谱重建技术。

1.2 Gold算法原理

Gold算法的目标是求解式(2)的最小二乘法解，即求解式(3)：

RTy=RTRx

(3)

令y*=RTy，M=RTR，则式(3)可改写为：

Mx=y*

(4)

基于式(4)，Gold算法的迭代流程可表示如下。

1) 取初值x(0)=[1,1,…,1]T，确定最大迭代次数K和停机准则δ。

在上述Gold算法流程基础上，Jandel等[9]又进一步提出了boosted-Gold算法，其原理是将Gold算法作为内循环，将若干步内循环得到的x进行激励操作，即令xi=|xi|p(i=1,…,n)，其中p(p>1)为用户选择的激励因子，然后将新的x作为下一次内循环的初始值，如此往复，直到外循环达到预设次数。

1.3 实验过程

本文选择了中国辐射防护研究院剂量学实验室自主研制的一台立式NaI型全身计数器作为实验对象。首先，利用点源对该全身计数器进行了能量和全能峰峰形刻度，获得了该设备的能量刻度函数和全能峰展宽函数。然后，使用该全身计数器对两套完全相同但灌装不同核素的BOMAB人体模型(下文以BOMAB-A和BOMAB-B代表)进行测量。BOMAB模型由中国辐射防护研究院研制，每套由10个不同形状大小的空心柱形聚乙烯容器组成以代表人体，聚乙烯容器内充满放射性水溶液。BOMAB-A灌装了4种活度已知的57Co、134Cs、137Cs、60Co放射性核素水溶液，而BOMAB-B仅灌装了水。测量两套BOMAB模型是为了更准确扣减γ能谱中天然辐射本底的贡献：将BOMAB-A测量谱减去BOMAB-B测量谱(按测量时间的比例折算)则可得到只由BOMAB-A中4种核素贡献的净谱，而将天然本底的计数贡献消除。本文分别对两套模型进行2 000 s的测量，经过扣减后，最终所得BOMAB-A净γ能谱如图2所示。能谱共2 048道，能量范围为0～2 MeV，图2还标注了能谱中各全能峰对应的核素名称和γ射线能量。

图2 实测BOMAB-A净γ能谱Fig.2 Measured net γ spectrum for BOMAB-A

a——探头接PMT端的透视图；b——探头远离PMT端的透视图图3 NaI探头的X光透视成像Fig.3 X perspective imaging of NaI detector

采用蒙特卡罗模拟方法计算式(2)中的测量矩阵R。计算方法如下：对该全身计数器和被测BOMAB人体模型进行完整的数字建模，其中全身计数器的支撑结构等对γ射线探测过程影响较小的部分未进行建模；为保证建模精度，对全身计数器中两个大体积NaI探头进行X光透视成像，如图3所示，通过影像获得细节结构的准确尺寸。图4示出了该测量场景的数字模型，其中，正对体模的为探测模块，由两个纵向排列的探测器和屏蔽结构组成。基于数字模型，使用MCNPX 2.5.1进行测量矩阵的计算。计算时，源粒子在BOMAB模型中均匀分布，且能量在50～1 500 keV等间隔变化，共计算800个能量点，每次模拟的源粒子数为107个，能谱的展宽函数使用了实测函数，形式如下：

(5)

其中：FWHM为全能峰半高宽，keV；E为全能峰峰中心的能量，keV。

图4 探测系统完整数字模型Fig.4 Complete digital model of detecting system

图5 蒙特卡罗模拟得到的响应能谱Fig.5 Measure matrix constructed by Monte Carlo simulation

另外，F8计数卡的能量范围为20～1 600 keV，能量节点数为1 500个，节点能量值则与实测能谱在20～1 600 keV区间内的各道能量对应。最终计算所得的测量矩阵R(即响应能谱)如图5所示，其维数为1 500×800。

2 实验结果

利用所得的BOMAB-A净γ能谱和测量矩阵R，使用Gold算法求解核素出射谱。图6比较了实测的BOMAB-A净γ能谱与重建能谱，其中，实测能谱的幅值放大了500倍。可见，重建能谱中高幅值区域与实测能谱中的全能峰位置吻合，其中，57Co核素的122 keV的全能峰收缩至单独1道，其余核素的全能峰也均有不同程度的收缩，同时，全能峰区域外则几乎不再有计数。这表明了重建能谱在全能峰收缩和基底消除方面体现了一定优势。然而，这一重建结果仍不够理想，因为核素出射谱的形态应为：所有全能峰均只占据1道而其余道计数均为0，从而能从中直接读取γ射线能量和数量，但基于目前的重建结果却无法实现直接读取定量结果。

图6 Gold算法的重建结果Fig.6 Reconstruction result given by Gold algorithm

为改善重建效果，本文进一步使用boosted-Gold算法进行计算。

当参数p=3时，boosted-Gold算法在第15次外迭代时到达停机准则。图7对比显示了实测BOMAB-A净γ能谱与boosted-Gold算法的重建能谱。由图可见，此时重建能谱已达理想状态，即全谱只有全能峰位置有计数，且所有全能峰均收缩至1道。

另外，实验发现，当1

图7 boosted-Gold算法的重建结果Fig.7 Reconstruction result given by boosted-Gold algorithm

由boosted-Gold算法的重建能谱，可直接读取各γ射线的能量和计数，而核素的活度可通过式(6)得到：

(6)

其中：A为核素活度；N为出射能谱中γ射线的计数；T为测量时长；η为该γ射线的分支比。

当核素有多个γ射线时，其活度可由最小二乘法确定。如134Cs有604.7、797.0 keV两条射线，60Co有1 173.2、1 332.5 keV两条射线，则它们的活度可由式(7)求得：

(7)

其中，下标1、2分别代表不同射线。

为比较能谱重建法与全能峰法的分析结果，使用Genie 2000软件对实测BOMAB-A净γ能谱进行分析。Genie 2000由Canberra公司开发并在核探测领域被广泛使用，其算法采用全能峰法[17]。

最终的计算结果列于表1、2。

表1 两种方法的γ射线识别能力比较Table 1 Comparison of γ identification ability between two methods

表2 两种方法的核素活度计算能力比较Table 2 Comparison of nuclide activity calculation ability between two methods

由表1可见，在γ射线能量确定方面，能谱重建法与全能峰分析法精度相当，给出的γ射线能量与实际能量误差在±3 keV以内。但值得注意的是，Genie 2000软件给出的结果中还包括一实际并不存在的84.2 keV的峰值。而由表2可见，在核素活度确定方面，能谱重建算法则显著地优于全能峰分析法，具体来说，对于BOMAB-A中的4种核素，能谱重建算法的活度计算值与实际活度的相对误差小于10%，最小仅为-2.2%，而Genie 2000软件给出的活度结果的误差除60Co外均大于10%，最大达到17%。另外，两种方法在57Co核素的活度计算结果上差异最大。

3 讨论

表1显示Genie 2000软件误识别出84.2 keV能量的全能峰，这印证了全能峰法在峰识别环节的弊端。全能峰法使用滤波法来识别全能峰，具体对于Genie 2000软件，其使用高斯函数二阶导函数作为滤波器[17]，滤波器展宽由全能峰展宽函数(由峰形刻度实验给出)确定。当能谱局部的形态与该处高斯函数相近时，算法判定存在全能峰。参照图2可见，在84.2 keV能量附近，实测能谱呈现近似全能峰的形态，这也是Genie 2000误识别的原因。然而，这一假峰结构事实上是由大体积源(人体模型)的低能散射峰效应和57Co核素122 keV能量γ射线的康普顿散射效应等综合因素造成的。由于原理所限，全能峰法无法辨别其真伪。

表2中，Genie 2000软件计算核素的活度时误差较大，这也是由全能峰法的原理所致。该方法计算全能峰净计数时需扣减峰下侧的基底，当扣减的基底与实际情况不符时，就会导致峰净计数偏大或偏小，从而造成对核素活度估算的不准确。对于实测的BOMAB-A净γ能谱，由图2可知，57Co核素122 keV全能峰处于康普顿散射计数堆积区域，下侧基底幅值较高且形态复杂，难以准确扣减，而一旦扣减不准确就会导致峰净计数误差很大，这也就是其对57Co核素的活度计算相对误差达到17.0%的原因；同时，能谱中134Cs的604.7 keV全能峰与137Cs的661.7 keV全能峰相互重叠，形成重峰，其下侧基底同样较为复杂，因此，其给出的这两个核素活度的相对误差也都超过了10%；而60Co核素的两个能量全能峰相互间重叠较少，且处在能谱高能端，下侧基底幅值低、形态简单，所以对其活度估算的误差最小。

能谱重建法不仅能准确识别所有主要γ射线，且没有任何误识别情况，同时所有核素的活度计算相对误差均小于10%。总而言之，在全身计数器γ能谱分析中，能谱重建法明显优于全能峰法。这是易于理解的，因为全能峰法只着眼于能谱的局部形态和计数，且掌握的其他信息极其有限，仅为能量、效率、峰形刻度函数，而能谱重建法不局限于部分能谱，而是将全谱作为分析对象，并通过测量矩阵将探测场景所有已知的物理、几何信息(如探测器、被测对象的结构、尺寸、材料等及射线与物质相互作用过程等)均包含其中。基于更丰富的全系统信息，其分析结果必然更准确。

当然，能谱重建法也有其明显的缺点：它的实施依赖于高精度的测量矩阵，而这必须通过细致的探测器表征及建模工作来实现，但该工作极为耗时耗力，使得该方法的应用成本相对传统方法要大；测量矩阵包含了被测对象的物理、几何信息以及其与探测器的相对位置信息，这意味着上述条件稍加改变，测量矩阵就会发生变化，故该方法只适合于被测对象单一且其与探测器相对位置固定的测量场景；复杂的迭代计算过程导致该方法无法从理论上对其计算结果给出定量的误差估计，而仅能通过与真值的比较而给出后验分析，这无疑也影响了其实施的可靠性。

4 总结

本文开展了将能谱重建技术应用于全身计数器γ能谱分析的实验研究。实验中，首先使用一台立式NaI型全身计数器测量两套BOMAB人体模型(一套为含源模型，一套为本底模型)，通过扣减本底谱得到了只含模型中放射性核素计数贡献的γ能谱；然后对该全身计数器和人体模型进行精细数字建模，并使用蒙特卡罗模拟方法计算了该系统的测量矩阵；基于实测能谱和测量矩阵，使用Gold算法和boosted-Gold算法分别求解出射能谱；最后，将求解结果与Genie 2000软件分析结果进行了比较。实验结果显示，单纯的Gold算法效果不佳，而boosted-Gold算法在参数p取1～10时能较准确重建出理论的出射能谱。根据boosted-Gold算法所求结果，能实现对所有核素的准确识别，且活度计算相对误差均小于10%。这一结果显著优于Genie 2000软件。

本文还深入分析了能谱重建法在核素定量分析方面优于全能峰法的机制。这些结论显示了基于能谱重建技术的γ能谱分析方法在全身计数器测量以及其他核探测领域进行实际应用的潜力。