基于透镜阵列多元聚焦的连续晶体三维位置灵敏探测器设计和仿真

李 婷,庄 凯,李道武,梁秀佐,刘彦韬,张译文,孔令钦,章志明,4,帅 磊,4,*,魏 龙,4

(1.中国科学院 高能物理研究所,北京市射线成像技术与装备工程技术研究中心,北京 100049;2.中国科学院大学 核科学与技术学院,北京 100049; 3.国家原子能机构核技术(核探测与核成像)研发中心,北京 100049;4.济南中科核技术研究院,山东 济南 250131)

射线成像技术广泛应用于生物医学、天文学、辐射环境监测、国土安全以及公共安全等领域[1-5]。位置灵敏探测器作为射线探测的关键设备,是决定射线成像系统性能指标的重要因素。基于闪烁体的位置灵敏探测器阻止本领高,位置分辨、能量分辨和时间分辨优异,还具有操作简易、稳定性好等优点,使得闪烁体探测器在射线成像应用中起着重要的作用。闪烁体探测器利用闪烁晶体将γ射线能量转化为闪烁光,再由光子计数器转化为电信号,电信号内包含γ射线和晶体作用的位置信息和能量信息。目前大多数高分辨的γ探测器采用晶体阵列耦合光子计数器阵列来实现。当以缩小晶体单元尺寸的方式进一步提高晶体阵列探测器的空间分辨时,受到了机械加工的限制,难以继续提升分辨性能,并且晶体间隔的填充降低了探测器的探测效率,这些是晶体阵列探测器存在的问题。

基于连续晶体的探测器可以避免晶体阵列带来的问题,具有设计结构简单、能量分辨率高、探测器效率高以及低造价等优点。随着新的闪烁晶体、硅光电倍增管(SiPM)和数字电子学技术的发展,连续晶体探测器有了新的可行性方案。连续晶体探测器通过所探测到每个γ事例的闪烁光分布来计算γ射线在晶体内的作用位置。因此当采用连续晶体结构的探测器时,高分辨的γ射线位置定位算法是不可缺少的。目前已有多种高分辨的定位算法被研究和应用开发,如神经网络定位算法[6-7]以及K-NN定位算法[8-10],这些算法都依赖于精确获取闪烁光输出分布。

光学成像技术中,在普通成像系统主透镜的一次像面处插入透镜阵列,每个透镜单元及其后对应的传感器区域记录光线在场景中相同部分在不同视角下所成像的集合,因此采用二维透镜阵列能得到同时包含位置和传播方向在内的四维光场数据,解决光场深度估计和光源定位问题。本文引入透镜阵列耦合连续晶体光输出面,提出一种新型的连续晶体探测器结构,实现对射线作用点的三维位置探测,并有望实现多个作用点的定位重建。在闪烁光传输过程中,透镜阵列中的各个透镜单元对其接收到的闪烁光聚焦,在光子计数器阵列平面上形成含有多个聚焦点的图像,每个聚焦点包含有闪烁光源位置和传播方向的信息,由聚焦点图像通过光路反演即可重建闪烁光源在晶体内的位置。这种结构的探测器不需要对连续晶体闪烁光分布进行精确探测,通过探测聚焦位置即可实现作用点定位重建。本文通过建立模型进行蒙特卡罗模拟研究,提出一种连续晶体对射线作用点的三维定位的解决方法,并通过重建算法结果分析来评估这一技术的可行性,并对晶体内两个作用点的情形进行研究。

1 重建算法原理

透镜阵列应用在光场成像领域中,通常单元透镜直径非常小,因而又可称为微透镜阵列。微透镜阵列作为一种基础的阵列光学元器件,对光信息有很好的聚焦、准直、交换、多重成像和综合成像的能力,已经应用于光场成像3D技术、深度估计等研究中[11-13]。在深度估计应用中,通过透镜阵列中的透镜单元在不同视场对同一个物体成像,相当于将三维物体的深度信息转化为二维透镜阵列所成像的角度信息[14-15],每个透镜单元所成的像都包含着深度信息。经过透镜中心的光线传输方向不变,以此简化光线的传播,可将透镜阵列简化为小孔阵列,光场成像深度估计原理如图1所示。不同深度的光源A、B的光线在经过透镜阵列后,在像平面的单元图像中位置不同,反映出A、B光源在物体空间中的位置差异。

图1 透镜阵列用于光场成像深度估计原理示意图Fig.1 Schematic diagram of depth estimation in light field imaging based on lens array

将透镜阵列应用于连续晶体探测器,射线在晶体内作用转化的闪烁光传输到晶体光输出面,再经由透镜单元传输到透镜与空气的界面上发生折射而改变闪烁光传输方向,闪烁光在空气中向靠近光轴的方向继续传输。经过透镜阵列时闪烁光的传播示意图如图2所示,蓝色光线表示晶体沉积射线能量后转化的闪烁光,红色光线表示从透镜阵列出射的闪烁光,出射闪烁光继续传输直到被光电器件探测,透镜对闪烁光传输方向的调制作用使得其接收到的闪烁光到达光子计数器阵列的位置更加集中。光子计数器阵列的光电器件采集图像并输出电信号,被后续电子学转化为数字信号形成数字图像。当调节光子计数器阵列与透镜阵列的间距,使得闪烁光被探测的位置集中形成的光斑小于光电器件像素大小时,采集图像可视为一组聚焦点的像素化图像,一次事例形成的图像中含有多个成像亮点。这些光斑或亮点仍包含有射线作用点的位置和沉积能量信息,成像亮点在图像中的位置与作用点的位置有关,成像亮点的光子数与射线的沉积能量、透镜对射线作用点的立体角有关。

图2 经过透镜阵列的闪烁光传播示意图Fig.2 Schematic diagram of propagation of scintillation light through lens array

采集到多个聚焦亮点图像进行射线作用点位置重建时,主要基于闪烁光路的光线对应关系和光的可逆性原理。射线能量不大时射线作用点可视为一个点光源,因而在闪烁光的传播路径上,闪烁光源、任一透镜中心和该透镜对应的图像亮点在一条光线上。将光路反演,图像亮点和透镜中心的连线反推到晶体内,反推的光线必定经过闪烁光源。算法中通过预设晶体内不同深度的平面,示意图如图3所示,计算反推光线与平面的交点Si,当各个透镜反推光线的交点在晶体的某一深度平面内汇聚为一点时,该点即为射线与晶体的作用点,此时的平面深度即为作用点的z值,再计算当前深度平面中反推光线与平面的交点即为作用点的x、y值。由此,基于透镜阵列的连续晶体探测器通过多点聚焦和光路反演,实现了晶体探测器对射线作用点的位置重建。

图3 多点聚焦的光路反演重建算法示意图Fig.3 Schematic diagram of reconstruction algorithm based on optical path inversion from multiple focus points

2 探测器模拟研究及算法处理流程

2.1 探测器蒙特卡罗模型及模拟参数设置

使用FLUKA蒙特卡罗模拟软件,建立连续晶体及透镜阵列模型以计算光子计数器阵列的采集图像。闪烁体选用尺寸为48 mm×48 mm×45 mm的连续晶体,材料为LYSO,折射率为1.82。晶体表面除出光面耦合透镜阵列端面外,均涂有黑色吸收层以防止多次反射造成多重影像。透镜采用半球型,半径为15 mm,选取3×3阵列,球心均匀分布,材料与晶体相同,在模型中可将透镜阵列与连续晶体视为一个整体,对闪烁光没有介质界面。光子计数器阵列采用SiPM阵列,像素大小为3 mm,量子效率为30%,设置在距离晶体出光面44.4 mm的位置。连续晶体、透镜阵列、SiPM阵列的对称轴在一条直线上,垂直于晶体出光面。采集图像中每个像素的数值为当前像素探测到的闪烁光子计数。

设置γ射线能量为511 keV,在晶体内仅发生一次作用。探测器的作用过程为:当γ射线与闪烁晶体发生作用时,沉积部分或全部能量产生闪烁光。闪烁光在连续晶体和透镜阵列的一体化结构中传输,当传输到透镜与空气的界面时发生折射,并沿折射方向继续传播至光子计数器阵列被探测成像,成像结果如图4所示。为研究晶体内某选定作用点位置的重建结果统计特性,将γ事例转换为质子事例,质子的出射位置即为γ事例的作用点位置。以晶体入射面中心点(0,0,-45)(单位mm,下同)开始设置作用点位置,在z∈[-45,-20]、x∈[0,8]、y∈[x,8](单位mm,下同)晶体空间内,间隔1 mm设置作用点。在x∈[12,24],y∈[x,24]晶体空间内,间隔4 mm设置作用点,z方向间隔仍为1 mm。本文将设置在晶体边沿的作用点x、y坐标取为23 mm,而非24 mm。

a——中心对称轴;b——偏离中心对称轴图4 作用点在晶体中心对称轴和偏离中心对称轴时成像结果示意图Fig.4 Schematic diagram of multiple focus point image when interaction position is on or departure from axis of crystal symmetry

2.2 多元聚焦的光路反演重建算法流程

算法处理流程分为以下3个步骤。

1) 确定亮点像素位置

采集图像以亮点为中心形成9个分立的区域,区域内除亮点外仍有少量的光子。因采用的半球型透镜对“点”物的成像存在像差,部分闪烁光子分布在焦点区域外,亮度低于焦点。首先对采集图像进行区域划分,区域划分采用K-Means聚类算法。划分区域后,在每个区域内采用峰值法确定亮点像素位置,记录为(xi,yi),i=1～9。

2) 计算重建位置

图5 重建时计算最小距离平方和参数结果示意图Fig.5 Schematic diagram of minimum distance square sum parameter in reconstruction

f(z)=∑‖Si(z)-Sj(z)‖2

(1)

3) 判选重建结果

固定一个作用点位置进行多次事例重建了后,对重建位置结果的x、y、z3个方向分别做拟合分析,x、y方向重建结果为高斯分布,超过3σ的重建结果认定为异常值,在后续分析中剔除。在算法程序中,可用式(1)的f(zSmin)作为参数进行判选。在单个作用点时,f(zSmin)的大小可以表征成像亮点位置确定的准确程度,从而可反映重建结果偏离真实值的程度。文中数据分析时设定f(zSmin)值超过40时,重建结果偏离真实值,不进行统计分析。

3 重建结果分析

3.1 单个作用点时探测器定位性能

以2.2节重建流程计算得到2.1节模型中作用点位置设置下1 000次事例的重建结果,重建结果绘制为散点图,如图6所示,以作用点位置为(4,4,-20)时的重建结果示意,画图时将z取绝对值(下同)。分析重建结果的各方向统计特性,如图7所示。以下对选取位置的重建结果分析时发现z方向的结果不是标准的高斯分布,在下文的计算中仍以高斯分布进行计算,对不同选取位置的重建结果进行对比分析。

a——未进行判选;b——经过判选图6 未进行判选和经过判选的重建结果散点图Fig.6 Scatter diagram of reconstruction result before or after data filtering

选取表1中6个作用点位置,分别计算得出重建结果的x、y、z3个方向平均值和标准差来评估当前结构探测器的定位能力。结果表明,透镜阵列的引入,能解决连续晶体的作用点定位问题,定位精度较好。

表1 固定在不同作用点位置的重建结果Table 1 Reconstruction result of interaction at different positions

对中心透镜所对的z∈[-45,-20]、x∈[0,8]、y∈[x,8]晶体空间区域的作用点重建结果的x、y、z3个方向误差进行统计分析,结果列于表2。结果表明,基于连续晶体和透镜阵列结构的探测器在x、y方向的定位准确度较高,在z方向的定位准确度较弱。

表2 中心区域重建结果定位误差情况Table 2 Position error of reconstruction result of interaction in central region

选取y=0 mm平面,以x=0、2、4、6、8 mm时重建结果的x、z方向平均值做分布图,如图8所示。重建结果的位置是不均匀、非线性的,在虚线框选的位置分别出现了z方向的重建位置重叠和xy平面内的重建位置重叠。这种现象的出现是由于采集图像的像素化和重建算法中对图像像素点的定位网格化带来的。

图8 y=0 mm平面固定作用点重建结果Fig.8 Reconstruction result of interaction in y=0 mm plane

3.2 单个作用点时探测器位置分辨性能

表3列出在不同区域的重建结果,分别计算作用点定位在x、y、z3个方向的误差及半高宽,3个区域分别对应中心、边线、边角位置上的透镜。其中,中心透镜对应坐标值范围为x∈[0,8]、y∈[x,8];边线透镜对应坐标值范围为x∈[0,8]、y∈[12,24];边角透镜对应坐标值范围为x∈[12,24]、y∈[x,24]。

表3 不同区域的作用点定位误差及位置分辨Table 3 Position error and resolution of reconstruction result of interactions in different regions

中心透镜对应的区域内误差及分辨性能优于边线透镜,边角透镜对应的区域内作用点重建结果误差及分辨最差。可重建事例数为数据判选后的事例数,数值高可表征发生在该作用点的事例容易被重建,可见晶体内可重建范围集中在晶体入射面、靠近晶体中心轴。以各方向的半高宽最大值计,xy平面位置分辨优于1.54 mm,z方向位置分辨优于3.13 mm。以各方向的半高宽中位值来看,均匀分布在晶体横向截面上的作用点,一半数量的作用点位置分辨优于0.82 mm;均匀分布在z方向上的作用点,一半数量的作用点位置分辨优于2.01 mm。探测器的位置分辨性能良好。

3.3 两个作用点时探测器定位性能

晶体内发生康普顿散射事例或偶然符合事例时,在两个以上作用点位置产生闪烁光,这些闪烁光在一次采集时间窗内不能被区分,因此光子计数器阵列采集到的图像是闪烁光的叠加图像。连续晶体通过光分布重建作用点定位的重建结果仅能给出一个作用点位置,且可能不能定位到任何一个作用点的真实位置上。透镜阵列的引入,将光子计数器阵列对光分布的精确连续探测转换成亮点的离散探测,不仅可以实现单个作用点的位置重建,还有望实现多作用点的位置重建。

本文研究发生偶然符合事例时双作用点的情况,射线能量为单一能量。2.1节模型中固定作用点的图像为一次光电事例的图像,当晶体内发生偶然符合的两次光电事例,成像结果数据可由单个作用点的图像进行加和来获得。选取3组事例分析双作用点重建的可行性,位置组合分别为(0,1,-40)+(0,1,-30)、(0,1,-40)+(0,4,-30)、(-8,-8,-40)+(8,8,-30),重复100次事例。

表4列出3组事例重建结果的统计情况。(0,1,-40)+(0,1,-30)、(-8,-8,-40)+(8,8,-30)两组事例的重建结果平均值及标准差与单个作用点的情况一致,(0,1,-40)+(0,4,-30)这组事例的位置重建结果较差。通过分析xy重建值的分布发现,按照2.2节重建流程的判选条件仍有偏离较大的结果未被剔除,这部分偏离结果是由于亮点位置定位不准确带来的。重复1 000次事例后分析xy重建值分布,超过3σ的重建结果认定为偏离结果,剔除这些偏离结果后这组事例的重建结果平均值及标准差与单个作用点的情况一致,但可重建事例率仅为7%,即当前组合位置的准确重建程度较低。图9为3组事例的重建结果散点图,其中图9b中已剔除偏差较大的结果。

表4 连续晶体内两个作用点重建结果统计情况Table 4 Reconstruction result of two interactions in monolithic crystal

a——(0,1,-40)+(0,1,-30);b——(0,1,-40)+(0,4,-30);c——(-8,-8,-40)+(8,8,-30)图9 连续晶体内两个作用点重建结果散点图Fig.9 Scatter diagram of reconstruction result of two interactions in monolithic crystal

以上重建是基于两个作用点沉积的能量接近,当能量差异较大时重建准确度会下降。考虑应用到康普顿散射时,如设定射线入射方向为垂直晶体,(0,1,-40)+(0,4,-30)这组事例的散射角约为16.7°,(-8,-8,-40)+(8,8,-30)这组事例的散射角约为66.2°。康普顿散射角较大的事例重建准确度较高,且康普顿散射角较大的事例散射点能量和吸收点能量相近,利于双作用点的定位。本文重点研究探测器的定位性能,后续研究如能进一步解决能量探测问题,透镜阵列结构的连续晶体探测器便可为康普顿射线成像技术提供一种新型的探测器构型。

4 结论

透镜阵列应用于连续晶体探测器时,射线沉积能量转换的闪烁光传输到透镜与空气的界面上发生折射,改变其传输方向,在透镜单元接收范围内的闪烁光出射到空气后均向靠近该透镜单元光轴的方向继续传输,直至被光子计数器阵列采集,形成光斑阵列的离散图像。光子计数器阵列采集位置在透镜焦距附近时,采集的光斑图像会形成聚焦点样式的图像。透镜阵列的引入,改变了连续晶体探测器中对出射闪烁光连续分布的精确探测,转换为对离散的闪烁光聚焦点阵列的探测。通过多元聚焦的光路反演重建能解决连续晶体探测器对射线作用点的定位问题,实现射线作用点的三维位置探测。文中模拟采用像素大小为3 mm的光子计数器阵列采集聚焦点图像,xy平面重建位置误差在1.04 mm内,z方向重建位置误差在3.70 mm内,xy平面位置分辨优于1.54 mm,z方向位置分辨优于3.13 mm,探测器性能良好。本文采用数据为连续晶体和透镜阵列一体化结构的较理想情况下的结果,实际探测时透镜阵列和连续晶体之间存在的界面将影响聚焦点图像的亮暗情况,从而影响“亮点”的定位,后续研究中需要改进“亮点”像素位置算法。透镜阵列的引入还能实现对连续晶体内双作用点的同时定位,有望应用于康普顿散射探测,但同样需要改进“亮点”像素位置算法并解决能量探测的问题。