辐射场中NaI（TI）闪烁体探测器响应全过程模拟

俞赛云，仇怀利，李佳，朱涛，宋逢泉

（1.合肥工业大学 物理学院，安徽 合肥 230601；2.中国科学技术大学 核科学技术学院，安徽 合肥 230031；3.福建省辐射环境监督站，福建 福州 350013）

引言

闪烁体探测器具有探测效率高、响应时间快等特点[1]，是目前应用较为广泛的核辐射探测仪器[2-4]。典型的闪烁体探测器系统主要由闪烁晶体探头、光电转换器及信号处理电路等组成，如图1所示。

图1 NaI(TI)晶体闪烁体探测器结构简图Fig.1 Brief diagram of NaI (TI) crystal scintillator detector structure

闪烁体探测器的基本工作原理为：射线粒子进入闪烁体时把部分或全部能量转移给电子，这些次级电子使闪烁材料分子电离或激发，接着被激发的分子退激，此过程以光辐射方式向外发出荧光光子，荧光光子在出射闪烁晶体后被光电倍增管光阴极收集，通过光电效应转化为光电子；光电子在电场中逐级倍增，在阳极聚集产生电脉冲信号，最后通过信号处理电路对输出电脉冲信号进行处理分析。

对上述过程的模拟是一项重要的技术手段，不仅有助于深入理解闪烁体辐射探测器的性能特征，对于发展新的辐射探测方法也都具有重要的指导意义。针对射线粒子在闪烁体中的输运，目前模拟大多是通过蒙特卡罗程序，如FLUKA、MCNP(Monte Carlo N particle transport code)、GEANT4 等。例如DEMIR N[5]、MOUHTI I[6]分别采用FLUKA、MCNPX 模拟光子在闪烁体中的沉积能量即脉冲高度光谱，与NaI(TI)探测器测量结果取得良好的一致性；HIRANO Y 采用Geant4 蒙特卡罗程序，模拟了射线与材料的相互作用和光子在闪烁体中传播的光学过程，并将仿真获得的空间和能量分辨率与探测器探测的实验值进行了比较[7]。对于粒子在闪烁体中沉积能量与光产额之间的关系，BIRKS J B 于1965 年提出著名的Birks 公式[8-9]。除了基础的模拟研究外，PAYNE S E、BIZARRI G、GAO F 等人还进行了进一步研究，针对不同类型的闪烁体开展了荧光产额的非线性研究[10-12]；WIRTH S 等人在对闪烁体阵列进行模拟时得到荧光逃逸效应对闪烁像素探测器的噪声性能起主要作用的结论，并发现光的传输过程可能会对信噪比产生影响[13]。关于可见光在介质中传播的模拟人们也常常用到蒙特卡罗方法，MENG F 等人在研究浑浊介质中极化光的传播时，使用蒙特卡罗方法进行数值计算，得到散射光的空间强度分布和Stocks 矢量模式，并结合实验测量结果得出散射光的空间因子，可用于标定散射系数、吸收系数等组织光学参数[14]；王建岗等人也利用蒙特卡罗程序模拟了光子在高散射介质中的迁移过程，包括光子每一步与介质的相互作用，即光子被散射和吸收，以及根据菲涅尔定律和斯涅尔定律分析光子在折射率不匹配边界上的反射和折射[15-16]；郑兴荣等人运用Matlab 软件对光传播特性进行了数据计算和仿真模拟，得到了光在介质面上传播的各种特性曲线，以此研究电矢量在介质表面的光传播特性[17]。

上述模拟大多数是针对单个物理过程，且依赖自研程序，通用性难以保证。因此，需要建立闪烁体探测器辐射场响应的通用方法模型，针对辐射场中晶体产生的输出信号进行模拟，为相关核探测技术提供基础理论模型。本文在此方面开展了初步的工作，以常用的NaI(TI)晶体为例，模拟分析了辐射场下闪烁体探测器输出信号，首先采用MCNP 程序[18-20]，获得粒子在NaI(TI)中沉积能量分布，然后结合Birks 公式和晶体光学参数，通过射线追迹程序，将晶体中沉积能量转换为晶体可见光信号输出，结合光电倍增管参数以及电路分析，将可见光信号转换为电信号，完成NaI(TI)探测器输出信号的模拟。本文通过实验验证了模拟方法的结果，分析了实验结果和模拟结果差异的原因，为下一步模型改进积累相关经验。

1 分析方法与模型

针对辐射场中NaI(TI)探测器的全过程响应进行模拟，分析流程如图2 所示。主要分为3 个过程进行模拟：射线粒子在闪烁体中输运过程、沉积能量转换为晶体输出光信号过程、晶体输出光信号转换为探测器输出电信号过程。

图2 全过程模拟方法的描述Fig.2 Description of whole-process simulation method

各个过程的模拟方法情况如下。

1）射线粒子在闪烁体中输运过程模拟，目的在于获得入射粒子在晶体中能量沉积的空间和时间分布。根据闪烁体探测器的几何与材料信息，尤其是晶体探头的参数，采用由美国洛斯阿拉莫斯国家实验室开发的MCNP 程序完成模拟[18]。MCNP 程序是通过模拟单粒子在指定系统内的径迹、核反应等行为，以获得相关变量的期望值，具有几何包容性强，减方差技巧丰富以及材料数据截面齐全等优点，是一种常用的解决中子、光子、电子及其联合输运问题的蒙特卡罗程序[19]。

2）沉积能量转换为晶体输出光信号过程模拟，既涉及可见光的产生，又涉及可见光的输运。闪烁体的荧光产额与沉积能量间存在非线性关系，沉积能量只有一小部分转化为可见光能量，剩余能量主要通过内部淬灭效应，激发能向邻近分子辐射和非辐射迁移以及自吸收进行消散。因此，BIRKS J B 提出以下关系[8-9]：

式中：dL/dx为单位路径的荧光产额；dE/dx为单位路径上带电粒子的能量损失；S为荧光发光效率，仅与闪烁体材料有关。根据BIRKS J B 的理论，沿着粒子电离路径的闪烁体分子分为损坏和未损坏两部分，B表征损伤分子和未受损分子的比例，K表征受损伤分子的淬灭概率，但在具体实验中KB只能作为一个参数进行测量[8]。

γ 射线在闪烁体中发生康普顿效应等，在产生荧光光子的同时，闪烁体中产生的次级电子伴随着致轫辐射、散射和正电子湮灭等作用沉积在闪烁体中，生成的荧光光谱范围一般包括可见光到紫外光。这里仅对荧光光子进行可见光传输模拟，采用射线追迹程序完成。使用射线追迹程序模拟γ 辐射粒子所激发的荧光可见光在晶体中输运过程时，计算中将荧光可见光作为光射线处理，通过求解一组光射线传输位置和波矢微分方程，追踪射线在模拟区域内的运行轨迹，获取晶体输出的荧光可见光功率参数。可见光传播模型如图3(a)所示。

图3 可见光传播示意图Fig.3 Schematic diagram of visible light propagation model

当光线入射到界面时与界面法线方向夹角不同时，光线在不同介质接触面上会存在折射、反射及全反射现象，如图3(b)所示。由于该模型中求解了射线功率，则“材料不连续性”特征将计算反射和折射光线这些变量的更新值。重新初始化使用的菲涅耳方程[21]为

式中：rp和rs是p 光和s 光的反射系数；tp和ts是p 光和s 光的反射系数；n1和n2分别是入射介质和折射介质的折射率；θ是入射角；θ′是反射角,θ=θ′；β是折射角。

射线追迹程序得出的数据仅是离散时间点的光输出，而可见光的发光是服从指数衰减的湮灭过程，所以要通过离散时间得到连续衰减，需要对已知数据进行展宽处理，将展宽后的数据整合在一起，得到光响应数据。

3）晶体输出光信号转换为探测器输出电信号过程模拟，涉及光电转换、电路仿真。光电转换过程可根据光电倍增管参数，如量子效率、电子渡越时间和时间分散、增益等，对光电转换后电脉冲信号进行评估，然后根据探测器电路原理图，对电信号处理过程进行仿真。由于本文使用的探测器的运放等器件来自ADI 公司，因此采用该公司开发的电路仿真软件LTspice 完成电路模拟仿真[22]。

2 模型测试与验证

2.1 模拟

以最为常见的NaI(TI)闪烁体探测器作为模拟对象，探头核心为Φ50 mm×50 mm NaI(TI)晶体，其侧面和前端分别覆盖2.8 mm 和1 mm 厚的MgO反射层，探头后端为1.8 mm 厚硅胶，再后边为4 mm厚的光学玻璃，最外侧为承担保护和降低噪声作用的铝壳，其结构为外径增加内径固定的空心圆柱，探测器的外部环境设置为空气，如图4 所示。放射源的类型为指向铝壳前端轴心发射的137Cs 点源，距离铝壳前端轴心距离为10 cm，与探头对称轴所成的角度为0°，这就是MCNP 模拟中模型设置。MCNP 模拟出的是射线照射到闪烁体能量沉积数据，通过Birks 公式可将沉积能量转化为荧光产额。

图4 NaI(TI)晶体探头的MCNP 计算模型Fig.4 MCNP calculation model of NaI (TI) crystal probe

本文利用MCNP 中Fmesh 计数功能计算了不同入射角度下光子在NaI(TI)晶体内的能量沉积情况，并通过公式(1)将光子产生的能量沉积转化为光产额。实验常用的137Cs源，能量为0.662 MeV，与NaI(TI)晶体主要发生康普顿散射和光电效应，产生的次级电子沿自身径迹损失能量，在此过程中NaI(TI)分子吸收能量并发出荧光。因此，射线主要沿其入射方向沉积能量，生成荧光光子，并且沿着入射方向不同位置的发光时间有先后。光子产生后，发光是一个指数衰减的湮灭过程[23]，满足以下公式：

式中：N为该单元的荧光产额；f(t)与N对应，t单位为ns 且t＞ 0，约1 500 ns 后可近似认为衰减完毕。

为了直观分析荧光产生分布与源入射角度的关系，对MCNP 的能量沉积输出结果进行可视化，图5 给出了射线0 °、30 °、60 °入射时荧光产生分布的三维切片图。

图5 不同角度入射能量沉积分布Fig.5 Distribution of incident energy deposition at different angles

当源入射角度为0 °时，光子沿探头的对称轴入射，切片图采用直角坐标系，是为了便于观察荧光产额进行对数转换，颜色越鲜亮，表示荧光产额越高。从图5 中可以看处，荧光产额自中心向外迅速衰减，并且随着入射深度增加，荧光产额也逐渐衰减。同理，当源入射角度为30 °和60 °时，光子也沿着入射方向的法向迅速衰减，并随着入射深度缓慢衰减。结合模拟数据与理论分析可知，射线入射到晶体后有一定的穿透深度，所以能量沉积的最大值并非晶体表面，而是与表面有一段距离。

接下来利用射线追迹程序模拟可见光在闪烁体中的传播。如图3(a)所示，直径为50 mm、高为50 mm 的圆柱形NaI(TI)晶体(图3（a）中2)，可见光出射处连接着光学玻璃(3-K9)，射线入射端面的反射层厚度为1.8 mm(1-MgO)，晶体和玻璃表面是厚度为3 mm 的反射层(4-MgO)。可见光在闪烁体和反射层的界面遵循菲涅尔定律，MgO 与空气接触面反射率设为0.9。可见光从光学玻璃的另一端出射，可将光学玻璃出射面设置为冻结，同时在此节点下添加一个沉积射线功率节点，以计算可见光的输出。图6 为某一时刻可见光在闪烁体中传播轨迹（图示仅选取部分射线）。

图6 部分可见光在闪烁体中的传播轨迹Fig.6 Propagation trajectories of partial visible light in scintillator

前文中提到荧光可见光产生时间与位置有关，且发光时间服从指数衰减，而射线追迹程序得出的数据仅是离散时间点的光输出，所以根据发光服从的指数衰减公式，将每个角度的荧光产额数据划分为离散时间的数据，导入射线追迹程序中，为方便之后Matlab 对数据插值展宽处理，这里取离散时间间隔为10 ns。对同一角度多次模拟得到的离散时间的光输出数据，利用Matlab 的插值函数spline 对这些数据进行一维数据插值再展宽拟合后，可得到连续的光输出，图7 为拟合的不同角度的曲线。

图7 不同角度入射时晶体输出可见光功率参数Fig.7 Visible light power parameters for crystal output at different incidence angles

根据光电倍增管的量子效率和增益可获得光电转换后的电信号。根据光电倍增管的渡越时间离散，对获得的电流脉冲信号，基于图8 电路原理图，采用LTspice 软件得到NaI(TI)探测器的电信号输出。图8 电路中，输入点为I1，即为光电倍增管输出的脉冲电流信号；输出点为R10 与R9 之间的电势差，输出端负载如图8 所示。

图8 NaI(TI)探测器中电路原理图Fig.8 Circuit schematic diagram of NaI (TI) detector

2.2 实验

基于137Cs 标准源进行验证实验，采用Φ50 mm×50 mm NaI(TI)晶体耦合光电倍增管（海南展创XP3240），结合数字化分析电路进行实验。使用的光电倍增管的量子效率为75 μA/lm，电子渡越时间和渡越分散时间分别为48 ns 和12 ns，光电倍增管的增益为7×105。以射线0°角入射为例，根据光电倍增管的渡越时间离散12 ns，通过Matlab 对0°数据做高斯展宽，高斯半高宽为12 ns，获得的脉冲信号如图9 所示。

图9 脉冲信号曲线Fig.9 Pulse signal curve

图10 为NaI(TI)探测器单个输出脉冲信号的实验结果和模拟结果比较。用单个脉冲的上升/下降时间的比值反映脉冲形状，实验结果和模拟结果分别为0.39 和0.36，相差7.69%，总体上脉冲信号的实验结果与模拟结果符合良好，为进一步模型开发提供了分析基础。

另外，从模拟曲线与实验曲线分布形态的细节上看，两者并不是完全吻合，存有一定差异，对于单个脉冲的上升时间，实验结果和模拟结果分别为215.2 ns 和85 ns；对于单个脉冲的下降时间，实验结果和模拟结果分别为550.1 ns 和236.3 ns。分析其原因，主要有以下3 点：1）在射线光学模拟计算中，采用了比较简化的模型，没有充分考虑光射线在晶体闪烁体中的散射过程和吸收过程；2）射线光学程序本身不能直接设定可见光输出功率随时间连续变化的函数，需通过对离散输出数据进行展宽和插值，数值拟合的时间连续的光输出函数曲线，在拟合过程产生了一定的误差；3）在晶体输出荧光可见光信号转换为探测器输出脉冲电压信号的模拟中，对电子在光电倍增管电路中渡越过程的计算，没有构建详细的计算模型。

3 结果与分析

本文对γ 辐射场中NaI(Tl)晶体闪烁体探测器响应信号输出的全过程开展了模拟计算分析，使用MCNP 程序、Birks 公式和射线追踪程序等模拟分析工具，模拟计算出γ 辐射场中典型NaI(Tl)晶体闪烁体探测器的光响应曲线，并通过137Cs 放射源开展验证实验，所得测试数据与模拟计算结果大体一致，表明本文所采用的模拟分析方法基本正确。将模拟结果对比实验测量的NaI(TI)闪烁体探测器的光响应曲线，可以看到整体趋势一致，但依然存在一定的误差，这是因为：首先在各阶段模拟过程中，为了模拟方便将模型简化，并且光传播模拟时没有考虑光在闪烁晶体中的散射、吸收；其次，射线追踪自身局限性也导致了模拟结果与实验测量结果存在差异；最后，晶体输出光信号转换为探测器输出电信号过程中没有建立光电倍增管的详细模型，对电子渡越过程的模拟不够详细。

本文的研究工作，对于深入理解NaI(Tl)晶体闪烁体中可见光射线的传输规律，开展闪烁体辐射探测器系统的优化设计具有一定参考意义。考虑以后将MCNP 与射线追踪结合使用，可以利用Geant4 模拟验证闪烁体光输出信号，同时对提高闪烁体探测器光提取效率做进一步的研究。