中国辐射卫生  2024, Vol. 33 Issue (3): 260-266  DOI: 10.13491/j.issn.1004-714X.2024.03.006

引用本文 

王玉剑, 王薇, 李兴隆, 汪传高, 庞洪超, 陈凌. 一种超铀核素致伤口污染成像探测装置的模拟和分析[J]. 中国辐射卫生, 2024, 33(3): 260-266. DOI: 10.13491/j.issn.1004-714X.2024.03.006.
WANG Yujian, WANG Wei, LI Xinglong, WANG Chuangao, PANG Hongchao, CHEN Ling. Simulation and analysis of an imaging detection device for wound contamination caused by transuranic nuclides[J]. Chinese Journal of Radiological Health, 2024, 33(3): 260-266. DOI: 10.13491/j.issn.1004-714X.2024.03.006.

通讯作者

王薇,E-mail:janitor0405@126.com

文章历史

收稿日期:2023-12-27
一种超铀核素致伤口污染成像探测装置的模拟和分析
王玉剑 , 王薇 , 李兴隆 , 汪传高 , 庞洪超 , 陈凌     
中国原子能科学研究院, 北京 102413
摘要目的 为解决超铀核素对伤口的放射性污染问题,开展基于编码孔成像技术的伤口辐射成像研究。方法 通过蒙特卡罗模拟多种源项,比较两种图像重建算法在近场条件下成像效果的差异,以及探测器像素和成像平面像素数量对图像分辨率的影响。结果 根据设计的尺寸进行成像系统模拟工作,模拟成像视野为89.4 mm × 89.4 mm,模拟角分辨率为1.98°,根据对比不同条件下重建点源的平均半高宽,得出增加探测器和成像平面中的像素数量可以优化角分辨率,但显著延长蒙特卡罗模拟时间。结论 根据模拟结果,该成像系统的参数可以有效的对放射性污染进行成像,为超铀核素致伤口污染测量提供了方法支持,同时为后续实际伤口污染成像探测系统的研发奠定基础。
关键词编码孔成像    图像重建    MLEM    伽马相机    伤口放射性    
Simulation and analysis of an imaging detection device for wound contamination caused by transuranic nuclides
WANG Yujian , WANG Wei , LI Xinglong , WANG Chuangao , PANG Hongchao , CHEN Ling     
China Institute of Atomic Energy, Beijing 102413 China
Abstract: Objective To address the radioactive contamination of wounds caused by transuranic nuclides, wound radiation imaging based on coded aperture imaging technology was investigated. Methods By simulating multiple source terms using Monte Carlo method, the differences in imaging performance between two image reconstruction algorithms under near-field conditions were compared. The effects of detector pixels and detection plane pixels on image resolution were investigated. Results The imaging system was simulated based on the designed dimensions. The simulated imaging field of view was 89.4 mm × 89.4 mm and the simulated angular resolution was 1.98°. Based on the comparison of the average width at half height of the reconstructed point sources under different conditions, it was found that increasing the number of pixels in the detector and detection plane optimized the angular resolution but significantly prolonged the Monte Carlo simulation time. Conclusion According to the simulation results, the parameters of the imaging system can be used to effectively image radioactive contamination. Our results provide methodological support for the measurement of wound contamination caused by transuranic nuclides, and lay the foundation for the development of wound contamination imaging detection systems in the future.
Key words: Coded aperture imaging    Image reconstruction    MLEM    Gamma camera    Wound radioactivity    

在乏燃料后处理等涉及超铀核素生产操作过程中,极易发生各种创伤(刺伤、割伤、爆炸伤及酸烧伤等)引起超铀核素侵入人体内。截至2001年相关科学文献中报道了超过2100起放射性核素致伤事件,大多数发生在生产、制造、维修或回收处理核燃料组件及已被沾污的设备中,且在大部分事件中90%以上的伤口为发生在手臂和手(主要是手指)处的穿刺伤和化学烧伤[1]。考虑到超铀核素的较强毒性和长期衰变特征,创伤后会绕过人体天然防御皮肤向淋巴结、血液、肝脏、骨等转移。相比于吸入、食入,将直接入血对人体产生内照射,即使相对低的摄入量也可能导致很高的内照射待积有效剂量。此外,某些核素的严重污染滞留还可能引起伤口处局部组织的癌变,因此创伤之后应尽快进行细致准确的伤口污染测量,当初步清洗之后,沾污量仍大到需要手术切除时,精确的伤口放射性污染定位、定量就更重要。由于乏燃料后处理中主要核素为超铀核素中的钚和镅,考虑到本工作中采用的像素型CZT探测器对低能伽马射线探测效率较低,故采用镅作为主要测量对象。

基于目前缺少对超铀核素致伤口污染分布的有效测量方法和手段的现状,本工作创新性的采用成像方式来实现对超铀核素致伤口污染分布特征的测量。从成像原理及算法选择、成像系统设计选型、成像性能指标模拟计算及成像效果影响因素研究等多方面展开研究,研究结论将为实际超铀核素致伤口污染分布的测量提供有益的方法支持,同时为后续适用于后处理等厂址中实际伤口污染成像探测系统的研发奠定技术基础。

1 材料与方法 1.1 成像原理及图像重建算法

辐射成像技术根据准直原理的不同可以分为编码孔成像技术和康普顿成像技术2种。其中编码孔成像技术采用针孔成像的原理,在针孔成像的基础上增加孔的数量以提高探测效率,通过这一方式编码孔成像一方面能够维持较好的成像视野,同时也能提供较好的角分辨率。康普顿成像采用康普顿散射的原理,其视场十分广阔,通常是全景视野,但康普顿相机要求一个粒子在相机内能够与探测器晶体进行多次碰撞,能量较低的射线通常在第一次碰撞时就已经完全沉积,因此对低能射线无法成像。鉴于后处理等场址中致伤源项核素主要为产生低能γ/X射线的超铀核素,因此选择编码孔成像作为超铀核素致伤口污染成像探测的基本技术。

编码孔成像分为成像和图像重建两部分。其中成像过程的物理模型如图1

图 1 编码孔成像的物理模型 Figure 1 The physical model of coded aperture imaging

根据几何光学的原理,探测器平面的成像结果满足公式[2-3]

$ P(\overrightarrow {{r_p}} ) = \iint_{\overrightarrow {{r_o}} } {\left[O(\overrightarrow {{r_o}} )h\left(\frac{{a\overrightarrow {{r_p}} + b\overrightarrow {{r_o}} }}{{\textit{z}}}\right){{\cos }^3}\theta \frac{1}{{4\pi {{\textit{z}}^2}}}\right]}{d^2}\overrightarrow {{r_o}} $ (1)

其中

$ \theta = {\tan ^{ - 1}}\left(\frac{{\overrightarrow {{r_p}} - \overrightarrow {{r_o}} }}{{\textit{z}}}\right) $ (2)
$ {\textit{z}} = a + b $ (3)

式中P表示成像结果,O表示源项,h表示编码板的编码函数,ab分别表示成像平面与编码板的距离和探测器与编码板的距离,θ表示源项位置与源项通过编码孔后投影位置的夹角,当远场成像(a >> ba近似取无穷远)时,该值可忽略不计,但与一般适用于远场成像的γ相机不同,本工作中对伤口污染分布的成像具有近场成像特征(ab尺度可比,不可认为无穷远),此时由于θ项的存在,编码孔成像的近场成像要比远场成像更为复杂。

适用于编码板成像的图像重建算法包括相关解码算法和最大似然期望最大化(MLEM)算法2种。其中相关解码算法从编码的原理上出发,将成像结果P与解码函数G进行相关运算后得到重建图像。由于只进行单次运算,计算速度快是相关解码算法最大的优点,但在近场情况下需要对该方法进行一些修正;MLEM算法则通过迭代计算的方法,假设源项分布并逐渐逼近真实图像,该算法能很好的抑制噪声,提高信噪比,因此是目前主流的图像重建算法,本工作拟采用MLEM算法进行伤口污染图像重建[4-6]

MLEM算法由泊松分布推导,其核心思想为首先对源项进行假设,将假设源项通过卷积等数学计算得到假设源项的投影项,将投影项与实际的成像结果对比,然后将对比的差别返回到假设源项中进行修正得到迭代源项,重复上述过程,使迭代源项逼近实际的源项,最终得到重建图像。MLEM的迭代公式为[7-8]

$ {x_j}^{k + 1} = \frac{{{x_j}^k}}{{\displaystyle\sum\limits_i {{a_{ij}}} }}\displaystyle\sum\limits_i {\frac{{{y_i}{a_{ij}}}}{{\displaystyle\sum\limits_j {{a_{ij}}{x_j}^k} }}} $ (4)

本文将式4)称为系统矩阵法,其中x表示重建图像,y表示实际的投影图像,a表示射线由成像平面第 j 个像素发射,被探测器平面第 i 个像素接受的概率,k表示迭代次数。

公式也可化为卷积形式[6,9]

$ {x^{k + 1}}(m,n) = {x^k}(m,n) \times \left[\frac{{y(m,n)}}{{{x^k}(m,n) * h(m,n)}} \otimes h(m,n)\right] $ (5)

本文将式5)称为直接卷积法,式中 × 表示乘法运算,*表示卷积运算,$\otimes $表示相关运算,h函数为编码板开孔函数。

1.2 成像编码孔准直器设计

本工作中设计的超铀核素致伤口污染成像系统包括基本成像探测器和编码孔准直器两部分,其中基本成像探测器采用16 × 16路,2.56 cm × 2.56 cm × 0.5 cm像素型CZT探测器。编码孔准直器是伤口辐射成像探测系统的核心部件,设计的好坏直接影响重建图像的质量。根据超铀核素致伤口辐射探测的特殊需求,对该系统编码孔准直器展开设计,设计中考虑的主要因素包括:

(1)考虑到成像效果,保证伤口污染成像探测系统的角分辨率设计值 < 5°;

(2)考虑到超铀核素发射的γ、X射线能量较低,伤口源项活度亦较低(约1 × 104 Bq),实际测量时需尽量提高成像效率以减少测量时间来减轻伤员痛苦,因此限制成像平面与探测器平面垂直距离小于10 cm;

(3)为适配CZT探测器尺寸,同时有效提高透光率,选取修正均匀冗余阵列MURA作为本工作编码孔准直器基本编码方式,同时从提高成像系统灵敏度角度考虑,最大限度的利用探测器的有效灵敏区域探测面积,因此限制编码板像素数量为13 × 13,此时对应CZT探测器有效探测面积为20.8 mm × 20.8 mm,有效利用像素数目为13 × 13=169个;

(4)考虑到编码孔准直器开孔加工难度,设计后的准直器编码孔尺寸应大于等于1 mm;

(5)为了便于实际应用,尽量减少测量时间,编码板设计在保证成像系统角分辨率5°的基础上尽量扩大系统的成像视场FCFOV。

1.3 成像系统性能指标理论计算

编码孔成像视野分为全视野场FCFOV和半视野场PCFOV两部分如图2。当编码板的图案被完整的投影到探测器上对应全视野场FCFOV,此时图像重建算法才能根据完整的投影准确重建出放射源的位置[10-11]

图 2 编码孔成像视野 Figure 2 Field of coded aperture imaging

编码孔成像的全视野理论计算公式为:

$ {{FCFOV}} = \frac{a}{b}{D_d} $ (6)

其中${{{D}}_{{d}}}$为探测器的总边长。

编码孔成像的分辨率可以用几何分辨率或角分辨率来表征。假设在空间中存在2个理想点源,相距为λ,分别透过编码孔后在探测器上形成投影图像,其投影图像刚好分布在2个相邻的探测器像素中时,此时λ值即为编码板的几何分辨率[12-13]。如图3

图 3 编码孔成像角分辨率 Figure 3 Angular resolution of coded aperture imaging

编码孔成像的角分辨率理论计算公式为:

$ \Delta \theta = {\tan ^{ - 1}}\frac{{{{\text{d}}_{\text{m}}}}}{{\text{b}}} $ (7)
$ \frac{{{{{d}}_{{m}}}}}{{{a}}} = \frac{{{{{d}}_{{d}}}}}{{{{a}} + {{b}}}} $ (8)

其中${{{d}}_{{m}}}$为编码孔的孔径。

综合考虑1.2节中各设计指标,利用上式计算并综合评估各设计参数下的成像性能后,最终确定系统设计参数及对应成像性能指标如表1

表 1 伤口成像探测系统初步设计参数及基本成像性能指标 Table 1 Preliminary design parameters and basic imaging performance indicators of wound imaging detection system
2 结 果 2.1 模拟条件

采用Geant4软件对表1所示伤口成像探测系统初步设计参数建模,模拟过程采用241Am源,点源和线源的模拟粒子数为1千万个,面源为1亿个。为减少粒子未经过编码直接射入探测器的影响,用钨在编码板周围建立4 mm厚的屏蔽体,为提高计算效率,设置粒子发射的角度为一锥形区域。建立的Geant4模型如图4。Geant4软件采用11.0.2版本,考虑光电效应、康普顿效应和电子对效应。根据模拟系统响应矩阵的时间计算,每一百万个粒子需要模拟 1 ~ 2 min,每 1000 组数据模拟需要 1 d左右。

图 4 编码孔成像的Geant4模型 Figure 4 Geant4 model of coded aperture imaging

本工作对不同数量源项粒子数的模拟成像效果进行了研究,研究发现粒子数为107个已经能够达到成像要求,因此认为蒙卡模拟造成的统计误差影响可以忽略。

将视野范围内的成像平面划分为52 × 52个像素,探测器平面划分为13 × 13个像素,在上述模拟条件下使用系统矩阵法迭代50次对单个、多个241Am点源、线源和体源的图像重建效果如图5

图 5 不同源项的重建图像、重建的三维图像 Figure 5 Reconstructed images and reconstructed 3D images with different source terms

模拟时点源由成像平面中心开始,延轴线间隔5 mm向成像平面边缘移动,每次移动过程中记录重建图像中峰值位置,并将其作为模拟位置,同时计算重建图像峰半高宽并计算其对应角分辨率,重建图像的计算结果如表2,其中模拟位置以重建矩阵的像素为单位。

表 2 重建图像的实际位置,重建位置,角分辨率 Table 2 Actual position, reconstructed position, and angular resolution of the reconstructed image

根据表2的蒙卡结果与理论设计对比,给出蒙卡和理论的成像视野及角分辨率。模拟成像系统的平均角分辨率为1.98°,视野为89.4 mm × 89.4 mm;理论计算成像系统的角分辨率为3.90°,视野为88.67 mm × 88.67 mm。验证了该成像条件能够满足设计需求。

表2可以看出,成像系统的平均角分辨率为1.98°。该值要优于理论计算的角分辨率3.90°,其主要原因是当利用针孔成像原理对系统分辨率进行理论计算时,认为点源的投影图像完全落在2个探测器像素时才能对源的位置进行区分,而实际上投影图像完全落在一个像素内和分布在2个探测器像素中时,也可以对源的位置进行区分,即利用针孔成像原理计算的分辨率理论值将偏大。

表2给出的由实际位置和模拟位置可计算该系统的成像视野。此时计算可得单位像素对应的实际距离平均值约为1.72 mm,结合成像视野52 × 52个像素,因此推算成像视野约为89.4 mm × 89.4 mm。与理论的88.67 mm × 88.67 mm相近。

2.2 模拟结论及分析 2.2.1 不同算法在近场条件下的系统分辨率对比

将视野范围内的成像平面划分为52 × 52个像素,探测器平面划分为13 × 13个像素,分别用系统矩阵法和直接卷积法模拟研究在近场条件下系统的成像分辨率。模拟时点源由成像平面中心开始,延轴线间隔1 mm向成像平面边缘移动,共移动43次,采用高斯拟合成像半高宽作为位置分辨率,重建时2种算法均迭代50次。计算结果见图6

图 6 模拟点源分辨率与位置的关系 Figure 6 The relationship between resolution and position of simulated point sources

图6可知,系统矩阵法重建后的系统分辨率具有周期性,但分辨率整体趋势不随点源位置变化,而直接卷积法重建后的系统分辨率从成像中心向边沿呈增大趋势,分辨率恶化,且这一问题难以修正。因此考虑到系统在整个成像平面内分辨率的均匀性,相比于直接卷积法,系统矩阵法更适合于近场成像。

2.2.2 系统矩阵法中不同探测器像素数量对系统分辨率的影响

将成像平面划分为52 × 52个像素不变,保证探测器平面大小为20.8 mm × 20.8 mm × 5 mm,像素数量分别取13 × 13、26 × 26、39 × 39、52 × 52,使用系统矩阵法研究在近场条件下不同探测器像素数量对系统的成像分辨率的影响。模拟时点源由成像平面中心开始,延轴线间隔1 mm向成像平面边缘移动,共移动43次,采用高斯拟合成像半高宽作为位置分辨率,重建时迭代50次。计算结果见图7

图 7 不同探测器像素数下系统矩阵法的分辨率与点源位置的关系 Figure 7 The relationship between resolution and point source position measured by the system matrix method with different pixel numbers of detectors

图7中可以看出,探测器平面取不同的像素数量时,将造成成像平面角分辨率的周期性分布,且该周期性与视野大小成正比,与探测器平面像素数目成反比。当探测器像素为13 × 13时,成像平面角分辨率周期约为6.85 mm;探测器像素为26 × 26时,成像平面角分辨率周期约为3.42 mm;探测器像素为39 × 39和52 × 52时,成像平面角分辨率周期分别约为2.28 mm和1.71 mm,由于模拟是以1 mm的间隔进行的,所以在这2种条件下周期性不明显。此外可以看出不同探测器的像素数下成像平面角分辨率的最小值基本不变,但会改变其分布的周期性,从而影响平均的角分辨率。

2.2.3 系统矩阵法中不同成像平面像素数量对系统分辨率的影响

将探测器平面划分为13 × 13个像素不变,保证成像平面总大小不变,像素数量分别取52 × 52、39 × 39,使用系统矩阵法研究在近场条件下不同成像平面像素数量对系统的成像分辨率的影响。模拟时,点源由成像平面中心开始,延轴线间隔1 mm向成像平面边缘移动,共移动43次,采用高斯拟合成像半高宽作为位置分辨率,重建时迭代50次。计算结果见图8

图 8 不同成像平面像素下系统矩阵法的分辨率与点源位置的关系 Figure 8 The relationship between resolution and point source position measured by the system matrix method with different imaging plane pixels

图8计算成像平面像素数为52 × 52时平均分辨率为3.45 mm,像素数为39 × 39时平均分辨率为3.99 mm,成像平面的像素数不会改变分辨率在成像平面分布周期性,但会影响平均角分辨率,成像平面划分像素数量越多,平均角分辨率越好,但相应的图像重建计算时间将延长,实际成像平面像素数目的选择应权衡成像效果和计算时间后最终确定。

3 讨 论

本工作根据理论设计编码板32.5 mm × 32.5 mm × 4 mm,编码板与探测器间距1.9 cm,编码板与放射源间距8.1 cm,利用Geant4模拟241Am源的成像并进行图像重建,模拟结果表明该成像装置设计参数可用于241Am的成像,模拟成像视野为89.4 mm × 89.4 mm,模拟角分辨率为1.98°。研究工作表明采用针孔成像计算的角分辨率结果偏大;在近场条件下相比于直接卷积法,系统矩阵法稳定性较好;探测器平面和成像平面像素数越多,图像的分辨率越好,但也会显著的增加模拟计算的时间。

目前编码板成像的产品及研究主要集中在场所的远场成像,但对于适用于伤口辐射成像的近场成像研究较少,未见成熟产品[14]。本研究从超铀核素伤口测量实际需求出发,提出并实现了以成像编码孔准直器为核心的伤口成像探测装置基本设计思路,结合图像重建算法及模拟仿真计算对不同成像条件展开了成像效果影响研究,为类似成像测量装置的设计积累了经验,同时也为下一步实验装置的搭建奠定了基础。

参考文献
[1]
吕晓凡, 卢丙慧, 冉新泽, 等. 核应急医学救援中放射性核素污染伤口的医学处理思考[J]. 中国辐射卫生, 2023, 32(4): 402-407,412.
Lyu XF, Lu BH, Ran XZ, et al. A reflection on medical treatment of radionuclide-contaminated wounds during medical response to nuclear emergencies[J]. Chin J Radiol Health, 2023, 32(4): 402-407,412. DOI:10.13491/j.issn.1004-714X.2023.04.007
[2]
Accorsi R. Design of a near-field coded aperture cameras for high-resolution medical and industrial gamma-ray imaging[D]. Cambridge: Massachusetts Institute of Technology, 2001.
[3]
李汉平, 王锋, 艾宪芸. 编码板成像系统MLEM算法优化[J]. 核技术, 2017, 40(2): 020404.
Li HP, Wang F, Ai XY. Algorithm optimization of MLEM in coded aperture imaging system[J]. Nucl Tech, 2017, 40(2): 020404. DOI:10.11889/j.0253-3219.2017.hjs.40.020404
[4]
Altomare C, Di Venere L, Fanchini E, et al. A gamma-ray imaging camera for ambient radioactivity detection[J]. Nucl Instrum Methods Phys Res Sect A: Accel, Spectrom, Detect Assoc Equip, 2020, 981: 164492. DOI:10.1016/j.nima.2020.164492
[5]
Hu YF, Fan P, Lyu ZL, et al. Design and performance evaluation of a 4π-view gamma camera with mosaic-patterned 3D position-sensitive scintillators[J]. Nucl Instrum Methods Phys Res Sect A: Accel, Spectrom, Detect Assoc Equip, 2022, 1023: 165971. DOI:10.1016/j.nima.2021.165971
[6]
Lee HS, Kim JH, Lee J, et al. Design and performance prediction of large-area hybrid gamma imaging system (LAHGIS) for localization of low-level radioactive material[J]. Nucl Eng Technol, 2021, 53(4): 1259-1265. DOI:10.1016/j.net.2020.09.008
[7]
Lange K, Carson R. EM reconstruction algorithms for emission and transmission tomography[J]. J Comput Assist Tomogr, 1984, 8(2): 306-316.
[8]
Tain JL, Cano-Ott D. Algorithms for the analysis of β-decay total absorption spectra[J]. Nucl Instrum Methods Phys Res Sect A: Accel, Spectrom, Detect Assoc Equip, 2007, 571(3): 728-738. DOI:10.1016/j.nima.2006.10.098
[9]
Mu ZP, Liu YH. Aperture collimation correction and maximum-likelihood image reconstruction for near-field coded aperture imaging of single photon emission computerized tomography[J]. IEEE Trans Med Imaging, 2006, 25(6): 701-711. DOI:10.1109/TMI.2006.873298
[10]
Sun SF, Zhang ZM, Shuai L, et al. Development of a panorama coded-aperture gamma camera for radiation detection[J]. Radiat Meas, 2015, 77: 34-40. DOI:10.1016/j.radmeas.2015.04.014
[11]
Byard K, Ramsden D. Coded aperture imaging using imperfect detector systems[J]. Nucl Instrum Methods Phys Res Sect A: Accel, Spectrom, Detect Assoc Equip, 1994, 342(2/3): 600-608. DOI:10.1016/0168-9002(94)90292-5
[12]
Shen XL, Gong P, Tang XB, et al. Encoding methods matching the 16 × 16 pixel CZT detector of a coded aperture gamma camera[J]. Nucl Sci Tech, 2020, 31(9): 92. DOI:10.1007/s41365-020-00796-5
[13]
孙世杰. 编码孔成像算法研究[D]. 哈尔滨: 哈尔滨工程大学, 2019. DOI: 10.27060/d.cnki.ghbcu.2019.000699.
Sun SJ. Research on the coded aperture imaging algorithm[D]. Harbin: Harbin Engineering University, 2019. DOI: 10.27060/d.cnki.ghbcu.2019.000699.
[14]
刘崎. 编码孔径成像放射源定位关键技术研究[D]. 成都: 成都理工大学, 2021. DOI: 10.26986/d.cnki.gcdlc.2021.000040.
Liu Q. Research on key techniques of radioactive sources localization in coded-aperture imaging[D]. Chengdu: Chengdu University of Technology, 2021. DOI: 10.26986/d.cnki.gcdlc.2021.000040.