扇形NaI(Tl)探测器阵列TGS透射图像的重建方法

扇形nai(tl)探测器阵列tgs透射图像的重建方法
技术领域
1.本发明涉及图像重建，尤其是一种扇形nai(tl)探测器阵列tgs透射图像的重建方法。


背景技术：

2.众所周知的：发展核科学与技术是我国社会经济不断发展和人民生活水平不断提高的需 要。在核技术应用和生产过程中，会产生大量放射性固体废物，通常采用标准核废物桶存放。 对于整备在钢桶内的固体核废物，需要对其进行放射性检测和甄别分类，按相应规定进行贮 存处置。层析γ扫描技术(tomographic gamma scanning,tgs)作为一种无损检测技术，能对核 废物桶内放射性核素进行识别和定量分析，并重建出桶内核素活度分布图，是桶装核废物放 射性检测的主要手段。目前，国内外tgs系统主要采用单个hpge探测器进行核废物桶扫描 测量，检测效率低、耗时长。随着核电子学的发展，对探测器的数据获取与处理技术水平逐 步提高，采用探测器阵列进行核废物桶tgs检测成为可能。相较于单探测器tgs系统，探 测器阵列tgs系统的检测效率更高。
3.21世纪以来，我国核科学与技术得到蓬勃发展，促进了国民经济和社会各方面的发展。 在核工业生产过程中会产生大量核废物，不同于其他废物，核废物具有放射性，若不能进行 有效地处理，会严重危害公共安全。核废物主要来源于核工业的前处理厂矿、乏燃料后处理、 核反应堆运行、放射性同位素生产及应用、核设施退役等(罗上庚，2001)。针对固体核废物 的贮存处置，我国制定了一系列标准(gb11928-1989,1990；gb9132-2018,2019)。固体核废物 通常封装和储存在钢桶内(ej1042-2014,2015)，根据《gb11806-2019，放射性物品安全运输规 程》和《gb12711-2018，低、中水平放射性固体废物包安全标准》的要求，核废物桶外表面 标志和标签需有相应的内容物、活度、临界安全指数等相关信息，故在固体核废物贮存处置 中，需对桶装核废物进行甄别分类，获取相关的桶内介质材料和放射性核素信息，继而按对 应贮存处置方式进行处置。
4.目前，桶装核废物检测方法一般采用非破坏分析技术(nondestructive assay,nda)(reillyet al,1991)，nda技术不需对放射性废物进行破坏性采样，就可对待测核废物桶中放射性核 素的种类、分布和含量进行准确分析，相较于需破坏样品物理结构的破坏性技术(destructiveassay,da)，更安全、经济、有效。分段γ扫描技术(segmented gamma scanning,sgs)(martinet al,1977)作为一种nda技术，广泛应用于桶装核废物检测过程中。在废物桶实际测量中， 核废物桶内含有低、中、高密度的固体废物，且其空间分布不均，不同种类的放射性物质也 分布不均。sgs技术将待测核废物桶按轴向分层，通过匀速旋转将桶内分布不均的点源近似 的看作线源，其较适用于核废物桶内放射性核素活度均匀分布、介质密度均匀分布或低密度 的情况。针对桶内活度与介质分布不均的情况，研究人员结合sgs技术和电子计算机断层扫 描技术(computed tomography,ct)，研究出了层析γ扫描技术(tomographic gamma scanning, tgs)(estep et al,1994；estep et al,1997)，在核废物桶内介质密度较高、放射性核素分布不均 时，应用该技术能对桶内的放射性核素进
行定量分析，并重建出桶内核素活度分布图，大大 提高了在此类复杂条件下的核废物桶测量分析准确性。
5.tgs测量技术分为透射测量和发射测量两部分，其在sgs技术的分层基础上，进一步将 每层划分为若干个体素，并假设每个体素内的介质材料和放射性核素呈均匀分布。透射测量 时，使经准直器约束的外置透射源发射的γ射线穿过核废物桶的每个分层，求得每层体素的 线衰减系数，还原出桶内线衰减系数分布图像；发射测量时，核废物桶内γ放射性核素自发 射的特征γ射线穿过桶内介质材料后发生衰减，结合透射测量阶段得到的线衰减系数对衰减 结果进行校正，还原出桶内放射性活度分布图像。
6.因此，可利用合适的图像重建算法提高核废物桶图像重建过程的精度和效率。目前tgs 系统主要采用单个高能量分辨率的hpge探测器测量γ能谱，每层采用平移和旋转的测量方 式，需要测量大量数据才能完成整个废物桶的扫描成像，导致整个系统的检测效率低下。采 用阵列nai(tl)探测器tgs系统，可以提高检测效率、降低仪器成本。tgs系统的准直器设 计是确保重建图像质量的关键因素之一，合适的准直器系统可以提高空间分辨率和准确度， 故对探测器阵列tgs系统的准直器开展优化设计研究，对提高我国桶装核废物检测水平具有 实际意义。
7.针对sgs技术的不足，20世纪80年代末，国外一些研究团队将ct技术应用到桶装核 废物γ无损检测技术中(illerhaus et al,1988；kawasaki et al,1990；estep,1990；martz et al,1990； eifler et al,1995)，形成了层析γ扫描(tgs)技术。1994年，美国lanl(los alamos nationallaboratory)的estep等人成功研制出了tgs装置，并对准直器和图像重建算法等方面进行了 研究(estep et al,1994)。与此同时，llnl(lawrence livermore national laboratory)的roberson 等人研制出基于a&pct(active and passive computed tomography)技术的impact(isotopemeasurements by passive and active computed tomography)实验装置(roberson et al,1994)。 1995年，lanl成功研制了移动式s/tgs系统并开始用于核废物桶无损检测(prettyman et al, 1995；taggart and vigil,1996)。同年，llnl成功研制出车载式具有透射测量和发射测量的装 置wit(waste inspection tomography)(bernardi and martz,1995)，经过测试研究后形成商业化 wit设备(martz et al,1999；camp et al,2002)。2007年，canberra公司研制出商用s/tgs 系统(venkataraman et al,2007)。
8.国内方面，中国原子能研究院最早对tgs技术进行研究。2001年，肖雪夫等人对tgs 装置的准直器、屏蔽层厚度以对测量效率矩阵的获取等问题进行了研究(肖雪夫等,2001)。 2003年，张全虎对tgs的探测效率刻度、透射和发射图像重建技术进行了研究(张全虎,2003)。 上海交通大学的王德忠、成雨、翁文庆、刘诚等人，研制了tgs装置并进行了一系列研究(成 雨,2007；翁文庆,2008；刘诚,2013)。2014年，张金钊利用蒙特卡罗方法提出一种适用于单 hpge探测器tgs检测系统的准直器设计方案(张金钊,2014)。
9.tgs系统主要由放射源、机械和控制系统、数据处理系统、探测器系统等部分构成。根 据tgs系统的组成部分，目前国内外对tgs技术的性能改良研究主要可分为以下两个方面：
10.(1)对探测器系统进行研究。
11.(2)对数据处理系统进行研究。数据处理方面的研究主要集中在体素划分方式、径
迹长 度计算、图像重建算法等方面。成像准确与否与体素划分的精细程度、径迹长度的计算准确 度息息相关。体素划分越精细，重建图像的分辨率越高，重建速度越低；径迹长度计算的越 准确，重建图像越准确。体素划分方式主要分为两种：在直角坐标系中，将体素划分为每个 体素尺寸大小一样的正方形块(roy et al,2015)，或者采用正方形块体素划分时，对于桶边界 的体素进行按其实际不规则的形状进行处理(han et al,2018)；在极坐标系中，将体素划分为 扇形体形状(成雨,2007；翁文庆,2008；刘诚等,2013)，开展自适应动态网格划分方法研究(guet al,2012)。在径迹长度计算方面，主要采用空间解析法(翁文庆,2008)、cyrus-beck裁剪算 法(张全虎,2003)和体素中心法(张金钊,2014)计算径迹长度。重建算法决定了图像重建质量， 对于投影数据较少的情形，应用诸如滤波反投影算法(filtered back-projection,fbp)之类的解 析类方法无法得到质量较好的重建图像，此时一般采用迭代方法进行图像重建。迭代重建算 法通过重复反馈过程，不断逼近所需的结果，其可分为代数迭代重建算法和统计迭代重建算 法。其中代数迭代重建算法包括art、sart(周志波等,2012)、mart、msart，统计迭代 重建算法包括mlem、osem等。代数迭代重建算法的收敛速度较快，其重建性能取决于多 种因素，如：迭代初值的选取、松弛因子等。统计迭代重建算法根据物理成像过程的特性和 实际测量数据的统计特性，将噪声等因素融入到模型建立中，但其收敛速度较慢。
12.近年来针对tgs技术中的迭代重建算法做了大量研究。如应用径向基函数(radial basisfunction,rbf)神经网络算法(张全虎,2003)、mc统计方法(张全虎,2003)、压缩感知方法(liu, 2012；wang,2014)进行tgs图像重建；对解增加先验条件约束，将一个不适定的问题转化为 一个与其近似的适定的新问题进行求解，提高算法的收敛性能和重建效率(thierry,2003；郭志 恒,2012；景越峰等,2016；张斌等,2014；吴朝霞等,2004)；结合层析成像原理对迭代算法相关 参数进行优化，如在迭代算法中选取适当的松弛参数(zheng,2018；孙杰等,2007)或者在迭代 过程中添加加权因子(雷兢等,2007)等；对迭代初值进行预处理(he et al,2018)，可以加快迭代 收敛速度，并提高重建精度。但是现有技术依然存在一些不足之处，比如：透射图像重建质 量较差，图像重建方法复杂等。
13.tgs技术利用窄束γ射线的衰减，结合ct成像技术，建立了桶装核废物tgs图像重建 方程。根据扇形探测器阵列tgs装置稀疏角度数据采样的特点，选取迭代类图像重建算法—— 迭代类算法求解欠定方程组性能更好，可获得质量较好的重建图像，但其计算过程复杂、耗 时长。


技术实现要素：

14.本发明所要解决的技术问题是提供一种能够提高透射图像重建质量的扇形nai(tl)探测 器阵列tgs透射图像的重建方法。
15.本发明解决其技术问题所采用的技术方案是：扇形nai(tl)探测器阵列tgs透射图像的 重建方法，tgs系统测量时将废物桶按轴向划分为数层，每一层又划分多个个体素，并设每 个体素的活度和密度均匀分布；所述扇形探测器阵列tgs，由于扇形探测器阵列的探测范围 涵盖整个废物桶；且废物桶可旋转360
°
；
16.所述扇形探测器阵列tgs系统包括放射源部件，机械和控制系统，探测器系统以及数据 处理系统；
17.所述放射源部件包括透射测量过程所使用的透射源及其用于提供测量所需的窄束特征γ 射线的准直器；
18.所述机械和控制系统，包括透射源升降装置、核废物桶旋转平台、探测器阵列升降装置。 该系统采用扇形探测器阵列扫描方式，在扫描测量时只需旋转废物桶以及升降探测器即可，
19.所述机械和控制系统用于保证透射源、核废物桶以及探测器在垂直方向的平动自由度以 及水平方向上的旋转自由度；
20.所述探测器系统，包括探测器以及与探测器匹配的准直器装置；所述探测器系统用于获 取待测核废物桶tgs的测量信息；
21.所述数据处理系统，包括用于实现对探测器所得数据的快速分析，完成桶内信息重建； 还包括以下步骤：
22.s1、设废物桶的图像值为f；根据mlem统计迭代算法：
[0023][0024]
s2、将上式转化为最小化问题：
[0025][0026]
tgs每层图像为二维图像，用fm,n表示图像第m行n列的像素点的值；
[0027]
s21、结合各向异性全变分模型得到：
[0028][0029]
式中，为降噪图像，为重建的噪声图像，
[0030]
s211、将步骤s11中得到的公式转化为分裂bregman形式：
[0031][0032]
式中，μ为保证项保真参数；
[0033]
s212、对步骤s12中得到的公式采用bregman迭代；
[0034][0035]
式中，λ为降噪因子；
[0036]
公式(1-1)中舍去不含dm和dn的项；对dm和dn采用交替极小化方法求解：
[0037]
[0038][0039]
其中，
[0040]
s22、根据各向同性全变分模型，得到：
[0041][0042]
令和分离l1范数和l2范数有：
[0043][0044]
式中，
[0045]
s221、对dm和dn进行最小化得到以下公式：
[0046][0047]
应用广义阈值公式得到
[0048][0049][0050]
式中，
[0051]
s3、将步骤步骤s212中得到的公式(3-1)中舍去不含的的项；以及步骤s221中得到 的公式(3-2)中舍去不含的的项进行l1和l2范数分开求解；
[0052][0053]
s4、对步骤s3中得到的公式进行优化，其优化条件为：
[0054]
利用傅里叶变换和傅里叶逆变换， 得到：
[0055]fk+1
＝f-1
(f(μi-λδ)-1
·
f(rhs))
[0056]
其中，
[0057][0058]
进一步的，步骤s1中的mlem统计迭代算法包括以下步骤：
[0059]
s11、设置图像值图像值：其中j＝1，2，...，j，j为划分的废物桶当前层体素总数；
[0060]
s12、对第i个探测器，计算投影估计值：
[0061][0062]
s13、计算投影估计值和投影测量值的误差：
[0063][0064]
s14、将比值反投影到图像空间，得到修正值：
[0065][0066]
s15、对xj的值进行修正，
[0067]
s16、k＝k+1，将上一轮迭代结果作为下次迭代初值，重复步骤s2～s5进行新一轮迭代 直到达到收敛条件或满足指定的迭代次数要求。
[0068]
进一步的，设置图像值：和误差界tol1，tol2；
[0069]
如果采用mlem迭代，得到
[0070]
如果
[0071][0072][0073][0074]
[0075][0076]
进一步的，设置图像值：和误差界tol1，tol2；如果
[0077]
采用mlem迭代，得到
[0078]
如果，
[0079][0080][0081][0082][0083][0084]
本发明的有益效果是：本发明所述的扇形nai(tl)探测器阵列tgs透射图像的重建方法 将分裂bregman算法以及tv降噪融入到传统mlem迭代过程中(mlem-sbtv)，对扇形 nai(tl)探测器阵列tgs系统透射扫描测量过程进行模拟，对模拟结果进行线衰减系数分布图 像重建。采用该扇形nai(tl)探测器阵列tgs透射图像的重建方法能够提高图像的重建质量， 并且简化图像重建的步骤，提高效率。
附图说明
[0085]
图1是本发明实施例中扇形nai(tl)探测器阵列tgs系统结构；
[0086]
图2为本发明实施例中扇形探测器阵列tgs系统测量示意图；
[0087]
图3是本发明实施例中模拟透射扫描测量示意图；
[0088]
图4是本发明实施例中线衰减系数透射测量模拟实验示意图；
[0089]
图5是本发明实施例中极坐标体素划分方式示意图；
[0090]
图6是本发明实施例中特征能量662kev下参考图像与透射重建图像；
[0091]
图7是本发明实施例中特征能量1173kev下参考图像与透射重建图像；
[0092]
图8是本发明实施例中特征能量1332kev下参考图像与透射重建图像；
[0093]
图9是本发明实施例中三种方法重建桶内介质线衰减系数分布图像的归一化均方误差示 意图；
[0094]
图中标识：1-透射源升降平台，2-核废物桶旋转平台，3-阵列探测器系统，4-探测器升降 平台，5-透射源和探测器准直器，6-核废物桶。
具体实施方式
[0095]
下面结合附图和实施例对本发明进一步说明。
[0096]
结合附图1至2所示，如图1和图2所示，tgs系统测量时将废物桶按轴向划分为数层， 每一层又划分多个个体素，并设每个体素的活度和密度均匀分布；所述扇形探测器阵列tgs， 由于扇形探测器阵列的探测范围涵盖整个废物桶；且废物桶可旋转360
°
；所述tgs系统采用 hpge探测器；
[0097]
所述扇形探测器阵列tgs系统包括放射源部件，机械和控制系统，探测器系统以及数据 处理系统；
[0098]
所述放射源部件包括透射测量过程所使用的透射源及其用于提供测量所需的窄束特征γ 射线的准直器；
[0099]
所述机械和控制系统，包括透射源升降装置、核废物桶旋转平台、探测器阵列升降装置。 该系统采用扇形探测器阵列扫描方式，在扫描测量时只需旋转废物桶以及升降探测器即可，
[0100]
所述机械和控制系统用于保证透射源、核废物桶以及探测器在垂直方向的平动自由度以 及水平方向上的旋转自由度；
[0101]
所述探测器系统，包括探测器以及与探测器匹配的准直器装置；所述探测器系统用于获 取待测核废物桶tgs的测量信息；
[0102]
所述数据处理系统，包括用于实现对探测器所得数据的快速分析，完成桶内信息重建。
[0103]
本发明所述的扇形nai(tl)探测器阵列tgs透射图像的重建方法，包括以下步骤：
[0104]
s1、设废物桶的图像值为f；根据mlem统计迭代算法：
[0105][0106]
s2、将上式转化为最小化问题：
[0107][0108]
tgs每层图像为二维图像，用fm,n表示图像第m行n列的像素点的值；
[0109]
s21、结合各向异性全变分模型得到：
[0110][0111]
式中，为降噪图像，为重建的噪声图像，
[0112]
s211、将步骤s11中得到的公式转化为分裂bregman形式：
[0113]
[0114]
式中，μ为保证项保真参数；
[0115]
s212、对步骤s12中得到的公式采用bregman迭代；
[0116][0117]
式中，λ为降噪因子；
[0118]
公式(1-1)中舍去不含dm和dn的项；对dm和dn采用交替极小化方法求解：
[0119][0120][0121]
其中，
[0122]
s22、根据各向同性全变分模型，得到：
[0123][0124]
令和分离范数和范数有：
[0125][0126]
式中，
[0127]
s221、对dm和dn进行最小化得到以下公式：
[0128][0129]
应用广义阈值公式得到
[0130][0131][0132]
式中，
[0133]
s3、将步骤步骤s212中得到的公式(3-1)中舍去不含的的项；以及步骤s221中得到 的公式(3-2)中舍去不含的的项进行和范数分开求解；
[0134][0135]
s4、对步骤s3中得到的公式进行优化，其优化条件为：
[0136]
利用傅里叶变换和傅里叶逆变换， 得到：
[0137]fk+1
＝f-1
(f(μi-λδ)-1
·
f(rhs))
[0138]
其中，
[0139][0140]
其中，mlem统计迭代的具体原理如下：
[0141]
极大似然期望最大化算法(maximum likelihood expectation maximization,mlem)作为一 类通过迭代进行极大似然估计的优化算法，该算法具有算法稳定、收敛性较好的优点，且其 抗噪声能力很好(闫镔等,2014)。mlem算法假定放射源发射的γ光子束服从泊松分布，cij 定义为穿过第j个体素到达第i个探测器的γ光子数，投影数据pi为第i个探测器测得的γ光 子之和，即：
[0142][0143]
假设cij服从均值为aijxj的独立泊松分布，则其似然函数为：
[0144][0145]
对数似然函数为：
[0146][0147]
上式中cij为未知变量，对上述对数似然函数求极大值不便。由于投影数据可将cij用其期望代替，求对数似然函数关于pi和当前参数估计值的条件期望，即：
[0148][0149]
式中，const为常数。
[0150]
cij服从均值为aijxj的独立泊松分布，且则：
[0151][0152]
式(1-4)对xj求偏导，令偏导的值为0得到极大值，有：
[0153][0154]
求解式(4-21)，第k次迭代后的图像值为：
[0155][0156]
综上所述步骤s1中的mlem统计迭代算法包括以下步骤：
[0157]
s11、设置图像值图像值：其中j＝1,2,...,j,j为划分的废物桶当前层体素总数；
[0158]
s12、对第i个探测器，计算投影估计值：
[0159][0160]
s13、计算投影估计值和投影测量值的误差：
[0161][0162]
s14、将比值反投影到图像空间，得到修正值：
[0163][0164]
s15、对xj的值进行修正，
[0165]
s16、k＝k+1，将上一轮迭代结果作为下次迭代初值，重复步骤s2～s5进行新一轮迭代 直到达到收敛条件或满足指定的迭代次数要求。
[0166]
其中，步骤2中的bregman迭代原理如下：
[0167]
在稀疏角度情形下，由于投影数据不完全，该问题是不适定的。通常使用正则化方法求 解不适定问题，其基本思路是对解增加先验条件约束，将一个不适定问题转化为一个与之近 似的适定问题进行求解。2005年，osher等人提出了一种新的基于bregman距离的反问题迭 代正则化过程。
[0168]
bregman距离定义：
[0169][0170]
式中，j:x
→
r为凸函数；即p为次微分中的一个次梯度。通常所以bregman距离并非传统意义的空间距离。由于且对于u,v 连线上的所有点w都有成立，bregman距离是用来衡量u,v两点的远近程 度的。bregman迭代主要求解以下问题：
[0171][0172]
式中，a∈rm×n，u∈rn，f∈rm，u0＝p0＝0，
[0173][0174]
令在osher和yin文中(osher et al,2005；yin et al,2008)表明式(1-10)可 以变为：
[0175][0176]
式中，j(u)＝||u||1；u0＝b0＝0。
[0177]
其中步骤3中的分裂bregman迭代的原理如下：
[0178]
对于j(u)＝||u||1这类基追踪问题，其目标函数仅为向量的范数，但有时需要极小化的 目标可能是一个向量值函数的范数，即：
[0179][0180]
式中，||φ(u)||1和h(u)都是凸函数，且假设φ(u)可微。
[0181]
将式(4-5)写成：
[0182][0183]
将上式转化为无约束最优化问题：
[0184][0185]
令e(u,d)＝||d||1+h(u)，a(u,d)＝d-φ(u)，结合bregman迭代有：
[0186][0187][0188][0189]
根据式(4-4)，将上式转化为：
[0190][0191]
利用交替极小化方法求解：
[0192][0193][0194]
参数d之间没有耦合，可以使用收缩算子求解：
[0195][0196]
其中，
[0197][0198]
综上所述分裂bregman迭代算法包括以下步骤：
[0199]
设置图像值：u0＝d0＝f0＝0和误差界tol＞0；如果||u
k-u
k+1
||＞tol
[0200]bk+1
＝bk+(φ(uk)-dk)
[0201]
for n＝1:n
[0202][0203][0204]
所述的扇形nai(tl)探测器阵列tgs透射图像的重建方法，一种优选的方式为 mlem-sbatv迭代，包括以下步骤：
[0205]
设置图像值：和误差界tol1,tol2；
[0206]
如果采用mlem迭代，得到
[0207]
如果
[0208][0209][0210][0211][0212][0213]
其次，根据各向同性全变分(isotropic total variation)模型，有：
[0214]
令和分离范数和 范数有：
[0215]
式中， [0216]
变量dm和dn不能在各向同性下解耦，为在该问题中正确地应用迭代方法，需解决子问题：
[0217][0218]
应用广义阈值公式求解有：
[0219][0220][0221]
式中，
[0222]
因此，所述的扇形nai(tl)探测器阵列tgs透射图像的重建方法另一种优选的方式为 mlem-sbitv迭代，具体的包括以下步骤：
[0223]
设置图像值：和误差界tol1,tol2；如果
[0224]
采用em迭代，得到
[0225]
如果，
[0226][0227][0228][0229][0230][0231]
实施例
[0232]
1、设置系统；
[0233]
扇形nai(tl)探测器阵列tgs装置采用旋转核废物桶的扫描方式获取单层待测核废物桶 相关信息，利用稀疏信号重建原理，以待测核废物桶几何中心为旋转中心逆时针匀速旋转， 每旋转角度30
°
进行一次透射扫描测量，共旋转360
°
测量13组数据。
[0234]
本实施例共采用7个nai(tl)探测器，以透射源为圆心呈等角(角度为9
°
)扇形排布，每 旋转一次测得7个测量值，最终得到13
×
7＝91个投影测量值。图3为模拟透射扫描测量示意 图。
[0235]
本实施例通过比较各体素的线衰减系数重建值与其对应的填充介质的真实线衰减系数值， 判断各体素的误差大小，将该误差作为透射图像重建结果的质量评价标准。利用mcnp程序 对透射扫描模拟实验中的六种材料进行线衰减系数透射测量模拟实验，线衰减系数透射测量 模拟同样利用该测量装置模型，如图4所示。
[0236]
线衰减系数透射测量需确保当前填充介质材料能完全遮挡住目标γ射线窄束，使得探测 器测得的计数均来自经材料衰减后的γ射线。参考本文扇形nai(tl)探测器阵列tgs系统特点， 仅使用4号探测器计数，作为将待测核废物桶替换为一个长方体小块，其厚度为10cm，长 度为35cm，高度为5cm，在长方体小块中分别填充6种介质材料：木材、铝、聚乙烯、硅 酸铝、混凝土和玻璃，结合式(2-1)的γ射线的衰减规律公式分别计算出6种填充介质材料的 线衰减系数值μ，其中i为4号探测器测得γ射线束经介质材料衰减后的计数，i0为4号探测 器测得γ射线束经空气衰减后的计数，x为介质厚度。将该线衰减系数值作为后续各体素重 建值的参考值，特征能量为662kev、1173kev、1332kev时，6种在材料的线衰减系数参考 值如下表所示：
[0237]
表4-1各特征能量下介质材料线衰减系数参考值
[0238][0239]
2、模拟透射图像重建；
[0240]
tgs图像重建过程中需要对核废物桶分层进行体素划分，结合图像重建算法重建出各体 素的值，进而重构出图像。核废物桶体素划分方式可大致分为两类：正方形体素划分方式和 极坐标体素划分方式。采用正方形体素划分方式虽然便于计算，但不能准确重建出桶边界处 体素值(韩苗苗等,2017)。极坐标体素划分是以角度θ将废物桶划分为360
°
/θ个扇区，然后将 桶半径r等距离l划分成r/l个圈层，形成的每个小网格为一个体素。当废物桶的旋转角度 为n
·
θ(n为正整数)时，各网格与旋转前的位置重合，仅体素编号发生变化，可以结合该特性 与探测器阵列的对称性划分体素，进而计算径迹长度和效率矩阵。
[0241]
根据透射测量稀疏采样方式，本节采用极坐标体素划分方式，以角度5
°
废物桶分层划分 为72个扇区，然后将该分层等分成12个同心圆层，每一个圈层含72个体素，共864个体素。 极坐标体素划分方式与编号方法如图5所示。
[0242]
根据表4-1六种介质的线衰减系数参考值以及其在桶内的空间位置分布，将参考值填入 对应的体素位置，得到各能量线衰减系数参考图像。将模拟透射扫描测量得到的91个投影数 据带入式(2-6)，利用mlem、mlem-sbatv和mlem-sbitv这四种重建算法求解该投影 方程，重建出桶内介质在能量为662kev、1173kev和1332kev下线衰减系数分布图像。3 种特征能量下桶内介质材料的线衰减系数分布图如图6至图8所示。
[0243]
从图6、7、8中的参考图像(a)和mlem算法透射重建图像(b)可以看出，相较于3.5.3节 的较大体素划分方式，本节采用的864个细小体素划分方式虽然提高了图像像素——空间分 辨率，但是会导致建立的透射测量方程组中的未知数个数远大于方程个数，而传统mlem算 法并不能准确地求解该欠定方程组，还原出不同桶内介质材料的形状及空间位置分布。在传 统mlem算法迭代过程中加入全变分降噪模型结合分裂bregman迭代方法，可以对该欠定 方程组的解加入具有稀疏性的先验，可以得到该不适定问题稳定的唯一解，有效地求解了该 欠定方程组，如图6、7(c)和(d)所示。相比于传统mlem算法重建结果，mlem-sbatv和 mlem-sbitv算法在3种透射能量下重建出的透射图像更符合桶内介质真实分布情况，图像 中的伪影更少、准确度更高。这两种算法均无法准确地重建出核废物桶边界体素的线衰减系 数值，这可能是采用“点-点模型”计算的径迹矩阵不够精确。“点-点模型”在计算体素径迹长度 时，忽略了透射γ射线束的发散角度，没有综合考虑该发散角度内的各条γ射线穿过体素的 径迹长度，使得桶内介质材料对γ射线的衰减贡献与真实情况不同，导致重建结果不太准确。 从图6、7还可以看出，mlem-sbatv的重建图像在介质边缘处更为锐
利，且其重建出的低 密度介质桶内分布区域小于该介质实际分布区域，而mlem-sbitv重建图像在介质边缘处 更为平滑，对低密度介质的重建效果更好。
[0244]
3、误差分析；
[0245]
透射重建图像的质量评价除其介质空间位置分布的准确性外，还应包括重建线衰减系数 值的准确性，但对透射图像进行主观评价会受限于人眼视觉系统，无法准确评价出重建线衰 减系数值的精确度，为了更准确评价3种方法重建的线衰减系数分布图像的效果，引入图像 质量客观评价参数：归一化均方误差(normalized mean square error,nmse)，其表达式如下：
[0246][0247]
式中，f(m,n)和g(m,n)分别为透射图像中第m行和第n列体素的线衰减系数重建值和真 实值。本章三种方法重建的核废物桶内介质线衰减系数分布图像的归一化均方误差，如表4-2 和图9所示。
[0248]
表4-2三种方法重建桶内介质线衰减系数分布图像的归一化均方误差
[0249][0250]
重建桶内介质线衰减系数分布图像的归一化均方误差越小，其重建图像质量越好。从表 4-2和图4-8可以看出，在不同透射能量下，mlem-sbitv方法重建图像的归一化均方差最 小，其次是mlem-sbatv方法，传统nlem算法最大。这表明mlem-sbatv和 mlem-sbitv均能有效地提高重建图像质量，且随着透射能量的增加，透射重建图像质量也 会有所提升。