中子照相系统整体空间分辨率计算与分析

贺林峰 韩松柏 王 雨 魏国海 武梅梅 王洪立 刘蕴韬 陈东风

（中国原子能科学研究院核物理所 中核集团中子散射应用技术重点实验室 北京 102413）

中子照相作为一种无损测试手段，在航空、航天、军工、核工业、建筑、考古、生物学、材料学、能源存储等领域发挥了重要作用[1−5]。空间分辨率是中子照相设备的核心参数之一，在开展各种实际测试之前，准确分析和评估特定设备的空间分辨率对于选择测量设备、优化设备参数及提高实验效率具有重要意义。

中子照相图像分辨率的测量方法通常是把具有明锐直边并且强吸收中子的薄片（镉片或钆片）置于样品平面上进行成像。在获得具有清晰边缘的图像之后，再利用图像分辨率计算方法进行计算[6]，如刀口分布函数法(Edge Spread Function, ESF)或者调制传递函数法(Modulation Transfer Function,MTF)。但是，这些实验方法是在设备建成后方能进行相应实验。对在建设备的分辨估算和优化设计还需构建相应的理论模型。

本文通过建立描述中子照相系统整体空间分辨率的二维物理模型，将几何不锐度和探测系统固有不锐度等因素引入系统脉冲响应函数，推导出计算空间分辨率的理论公式。利用该公式推算了中国先进研究堆(CARR)中子实时探测系统在不同准直比(L/D)条件下的空间分辨率，讨论了L/D、闪烁分辨率、CCD（或CMOS）相机分辨率等参数对装置整体分辨率的影响规律，并探讨了如何优化调整各影响因素以调高整体分辨，为中国先进研究堆中子照相平台的优化设计提供了重要的理论依据。

1 中子照相系统空间分辨率计算公式推导

中子照相二维数字图像可用线性移不变函数来表示[7]：

式中， jT( x,y)是中子束流在透射样品后投射在闪烁屏上的理想图像； h (x,y)是系统脉冲响应函数（又称点分布函数），它表示照相系统对像点源的弥散分布； I(x,y)可表示为 jT( x,y)和 h(x,y)的卷积。对系统脉冲响应函数 h (x,y)进行建模分析即可评估仪器的设计分辨率和实际分辨率。

图1表示包含可影响理想中子成像系统分辨率的基本物理参数的模型。中子照相系统的数字图像的空间分辨主要影响为：准直比L/D和样品与闪烁屏之间的距离造成的几何不锐度；闪烁屏处光斑弥散造成的固有不锐度；CCD（或CMOS）芯片像素尺寸受限引起的分辨误差（由于样品的随机性，为了便于统一评估，此模型未引入中子在样品处散射造成的分辨误差）。

由于这三方面对系统的影响是相互独立的，将它们引入系统脉冲响应函数 h(x,y)，可得：

式中， hD(x,y)代表准直系统产生的几何不锐度，它可表达为通过光阑的中子束通过样品处点源在闪烁屏处产生的圆斑(图2)，其直径为d/(L/D)，d为样品到闪烁屏距离，D为准直光阑直径，L为光阑与闪烁屏距离。

图1 中子成像系统分辨率模型[7]Fig.1 Model of neutron radiography special resolution[7].

hD(x,y)的数学表达式可写为：

此处定义函数 C YL[r/r0]为：

hs(x,y)代表闪烁屏本身的固有不锐度，即是单个中子被闪烁屏吸收后所产生的弥散发光区域。闪烁屏发光原理是中子被闪烁屏中的 Li6吸收后产生α粒子，α粒子随后将ZnS(Ag)颗粒激发发光。因此α粒子扩散距离及闪烁发光体本身颗粒度大小都会影响单个中子激发发光区域面积。它可表达为像点源在闪烁屏处产生特征值为d的高斯分布圆斑(图3)。

hs(x,y)可写为二维高斯函数：

此处定义函数 G AUS[x]为：

图2 光阑的中子束通过样品处点源在闪烁屏处产生的投影圆斑Fig.2 Projection of object’s point at scintillator.

图3 像点源在闪烁屏产生的扩散圆斑Fig.3 Spread area of image point at scintillator.

hCCD(x,y)代表 CCD芯片像素尺寸受限造成的分辨误差，它可表示为 CCD单位像素在闪烁屏上边长为 D s/MCCD的投影光斑(图4)。其中 Ds为CCD单位像素尺寸， MCCD为CCD相机对闪烁屏物体图像的放大倍数。

图4 探测器像素尺寸在闪烁屏上的投影光斑Fig.4 Projection of CCD’s pixel at scintillator.

它可表示为矩形函数：

此处定义函数 []RECTx为：

将式(3)、(5)、(7)代入式(2)，可得系统脉冲响应函数 (,)hxy的二维数学模型：

为便于实际测量应用，我们通常将系统脉冲响应函数进行一维空间调制，得到的一维调制传递函数(MTF)即可表征系统分辨率：

首先将上述系统脉冲响应函数 (,)hxy式(9)转换为一维模型函数：

然后将 h (x)的三项卷积进行傅里叶变换得到H(u)。这里将 h(x)的三项卷积进行傅里叶变换可认为此三项傅里叶变换后再点乘[8]。

其中 R ECT[x]、G A US[x]的傅里叶变换形式[8]：

这里定义 S INC(x) = sin(π x) /(πx )。

应用傅里叶变换关系式(12)、(13)可得到 ()Hu：

通过式(14)很容易得到 ()MTFu数学表达式：

对于中子照相而言，一般取MTF值为0.1时，对应横坐标值为系统的空间分辨率（u值），其单位为 cyc·mm−1。若已知上述三项函数的特征参数d/(L/D)、δ、Δs/MCCD，代入式(15)即可求出MTF值为0.1的u值。

2 中国先进研究堆中子实时成像探测系统分辨率计算

中国先进研究堆(China Advanced Research Reactor, CARR)中子实时成像探测系统是用于开展动态快速成像测量的装置，结构如图5所示。其采用了AIT公司生产的EJ426HD2型热中子闪烁屏，尺寸为20 cm×20 cm，厚度为500 μm；光学反射镜为浮法石英玻璃镀铝膜；定制光学镜头焦距 f为50mm 光圈 F为 0.9；PCO公司生产的PCO-edge型科学级CMOS相机，2560×2160像素，像素尺寸6.5 μm×6.5 μm，全幅条件测量速率为 0.5−100 fps[9]。

对于该探测系统而言，式(15)的三项特征参数中Δs/MCCD已知，若通过实验得到表征闪烁屏特性的特征参数δ，即可计算得到不同d/(L/D)数值条件下中子照相系统的整体空间分辨率。

在中国先进研究堆利用中子实时成像探测系统在水平孔道进行中子照相实验，将具有明锐直边的三角镉片贴在闪烁屏处进行测量。中子束流强度约为 105n·cm−2·s−1，所用相机的像素为 2160×2560，曝光时间为2 s，连续拍摄 300张，累积曝光时间600 s。利用国际中子成像领域常用开源图像处理软件ImageJ对单帧图像进行降噪处理。降噪方法采用九格中值法对图像进行处理，阈值设为50。样品镉片边缘白光处灰度不超过 120，因此这一方法主要过滤伽马白斑噪点。由于300张照片中在同一像素处有伽马噪点的照片不超过5张，因此此种降噪方法带来的误差可以忽略。将获得减噪后的图像叠加以增强数据统计性，再进行归一化处理。

根据所获得数据图像（图6左），利用ImageJ软件计算出镉片边缘灰度曲线（图6右）并对其进行 MTF分析处理，得出此中子图像分辨约为0.36mm。

在上述实验中，由于镉的热中子吸收截面很大，可忽略样品厚度对分辨的影响，并且镉片紧贴在闪烁屏上，可以近似认为样品到闪烁屏的距离d为0。因此实验测得的图像分辨率只包含中子照相系统空间分辨率计算公式(12)的后两项，即：

当式(14)为 0.1 时，u=1/0.36 cyc·mm−1，其中Δs/MCCD=0.136 mm，对其数值求解，可得代表闪烁屏分项的高斯函数特征值δ=0.29 mm。

图6 去噪后的镉片中子成像图像（左）和镉片刀口分布曲线（右）Fig.6 Rough neutron imaging without noise (left) and edge spread function (right).

此中子实时成像探测系统应用于实际的中子照相装置时，考虑到样品尺寸及三维成像的实验需要，计算时统一设定样品到闪烁屏距离d为10 cm，将上述参数代入式(12)可求得L/D分别为100、250、500时，系统整体空间分辨率为1.1 mm、0.53 mm、0.40 mm。

3 影响空间分辨率主要参数的优化与分析

在中子照相系统空间分辨率计算公式(15)中，三个特征参数d/(L/D)、δ、Δs/MCCD分别对应着系统几何不锐度（准直比L/D和样品与闪烁屏之间的距离d），闪烁屏的固有不锐度，CCD（或CMOS）芯片像素尺寸等影响因素。此外，中子照相系统的中子通量与L/D的平方成反比关系，即L/D越小则照相测试效率越高。综合考虑上述因素的对空间分辨的影响以及测试效率和分辨之间相互制约的关系，细致优化这三个特征参数将最大发挥中子照相装置的应用效果。

图7 L/D为250 (a)和500 (b)时系统及各分项的MTF曲线Fig.7 MTF of system and components at L/D=250 (a) and 500 (b).

以上述中国先进研究堆中子实时成像探测系统为例，由于采用了特定的闪烁屏和CMOS芯片，此时δ、Δs/MCCD等参数为定值。假如一准直比L/D可在 250−500内调整的中子照相装置采用此探测系统，可以分别计算出L/D（d值定为10 cm）、闪烁屏、相机和系统整体的MTF分辨曲线。

当L/D为250时，各参数对应的分辨曲线如图7(a)所示，此时整体空间分辨率主要受L/D的制约，闪烁屏影响次之，CCD影响最小。可见在此种L/D条件下并没有发挥该探测系统的最佳分辨性能。

当L/D=500时，各参数对应的分辨曲线如图7(b)所示，此时系统整体空间分辨率主要受闪烁屏的制约，L/D的影响次之，CCD的影响最小。虽然，在此种L/D条件下已经发挥该探测系统的最佳分辨性能，但是由于L/D值较大，反而过多损失了中子通量，降低了测试效率。

就该探测系统而言，应该以闪烁屏的分辨能力为基点优化照相装置的 L/D值，即让对应 L/D的MTF曲线与闪烁屏 MTF曲线相交于一点，此时u=3.1 cyc·mm−1，可推算出L/D为341。那么此时，照相装置在发挥探测系统最佳分辩的同时，保证了最优化的测试效率。

通过以上的计算分析可以看出，在使用或设计中子照相装置时需要综合考虑各影响因素，通过优化配置几何光路和探测系统参数，以发挥中子照相装置的最佳性能。

4 结语

本文通过建立涵盖几何光路和探测系统参数的中子照相二维脉冲响应函数理论模型，推导出中子照相系统整体空间分辨率计算公式，探讨了公式中各参数与整体空间分辨的定量关系和影响规律。以中国先进研究堆实时中子成像探测系统为例验证了该公式的可行性，并分析讨论了其在优化设计中子照相装置和确定实验测量参数等方面的应用价值。

1 貊大卫, 刘以思, 金光宇, 等. 中子照相[M]. 北京: 原子能出版社, 1996 MO Dawei, LIU Yisi, JIN Guangyu, et al. Neutron radiography[M]. Beijing: Atomic Energy Press, 1996

2 丁大钊, 叶春堂, 赵志祥, 等. 中子物理学(下册)[M].北京: 原子能出版社, 2001: 718–719 DING Dazhao, YE Chuntang, ZHAO Zhixiang, et al.Neutron physics (Vol.2)[M]. Beijing: Atomic Energy Press, 2001: 718–719

3 Johnson M W. The industrial uses of neutron[J].Applied Radiation and Isotopes, 1995, 46(6/7): 673–680

4 Schillinger B, Brunner J, Calzada E. A study of oil lubrication in a rotating engine using stroboscopic neutron imaging[J]. Physics B, 2006, 385–386: 921–923

5 Lehmann E H, Vontobel P, Hermann A. Non-destructive analysis of nuclear fuel by means of thermal and cold neutrons[J]. Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers,Detectors and Associated Equipment, 2003, A515:745–759

6 王雨，韩松柏，贺林峰，等. 中子照相图像分辨率计算方法[J]. 核技术, 2012, 35(4): 275−280 WANG Yu, HAN Songbai, HE Linfeng, et al. Calculation methods for neutron radiography spatial resolution[J].Nuclear Techniques, 2012, 35(4): 275−280

7 Anderson I S, Hurd A J, McGreevy R L, et al. Neutron imaging and applications[M]. Springer-Verlag New York Inc, 2009

8 Gaskill J D. Linear systems, fourier transforms, and optics[M]. John Wiley & Sons, Inc, 1978

9 He L F, Han S B, Wang H L, et al. Design of real-time neutron radiography at China advanced research reactor[J].Physics Procedia, 2013, 43: 48−53