基于缪子离散能量的材料鉴别实验研究*

李雨芃 汤秀章 陈欣南 高春宇 陈雁南 范澄军 吕建友

(中国原子能科学研究院,北京 102413)

宇宙射线缪子具有穿透力强、对重核材料敏感的特点,近年来被广泛应用于核材料检查等领域.缪子与不同原子序数材料发生的多重库仑散射效果不同,利用这点可以对被测物体进行成像以及材料鉴别,而在该过程中引入缪子的能量信息可以提高材料鉴别的准确度.本文基于原子能院研制的缪子成像装置开展了5 种样品的材料鉴别实验,使用离散能量拟合近似连续能量的缪子散射角分布,进而测量出各材料的辐射长度并以此作为特征量进行材料鉴别.实验结果表明,在约1400 个有效缪子事例下,Al-Fe 和Fe-Pb 可以在95%的置信水平下有效区分,该方法相比于不引入缪子能量信息对Pb-W 鉴别的准确率提高了18.5%.

1 引言

缪子是一种由宇宙射线与大气层相互作用产生的高能带电粒子,质量是电子的207 倍,平均能量为3—4 GeV[1].缪子的穿透能力强,可以有效地穿过屏蔽层对物体内部结构进行成像探测,并且无放射性危害.目前缪子成像技术已被应用于核材料检查[2,3]、金字塔探测[4,5]、陵墓探测[6]、火山监测[7,8],以及对反应堆[9,10]和核废料[11]的监测等多个领域.

2003 年美国洛斯·阿拉莫斯国家实验室(LANL)首次提出利用缪子散射成像方法进行核材料检查[12],该技术利用缪子穿过物体时发生多重库仑散射的角度偏转与材料原子序数、缪子能量的规律对被测物进行成像[13]和识别[14].在实际应用中,天然缪子的能量是连续的且实时测量非常困难,因此大多数研究使用缪子的平均能量代替未知缪子能量,这种近似导致了成像图像质量和材料鉴别准确度的降低.Vanini 等[15]开展了不同能量精度下的高炉成像模拟研究,Bae 和Chatzidakis [16]开展了多种特殊核材料的材料鉴别模拟研究,Morris 等[17]提出了缪子能量的多群模型等.这些模拟结果均表明缪子能量信息对于成像和材料鉴别有着重要影响.此外,部分学者也开展了缪子能量测量的实验研究.例如,加拿大团队[18]研制出了带能量测量结构的CRIPT 探测器,通过测量已知材料阻挡缪子产生的散射角还原缪子动量;清华大学罗志飞[19]利用多气隙阻性板室(multi-gap resistive plate chamber,MRPC)时间分辨率高的特性,使用飞行时间法计算出缪子的分段能量信息.

本文基于原子能院研制的缪子成像装置开展了5 种不同原子序数样品的材料鉴别实验研究,采用离散能量缪子散射角分布拟合近似连续能量的缪子散射角分布,以辐射长度作为特征量进行材料鉴别,并分析实验测量误差对材料鉴别准确度的影响.

2 基于缪子离散能量信息的材料鉴别方法

当缪子穿过靶材料时,会受到原子核库仑力的作用发生多重库仑散射,其散射角在同一平面内的投影近似服从均值为0、标准差为σθ的高斯分布[20].根据Molière 理论[21,22],散射角分布宽度与缪子能量、材料辐射长度有如下关系:

其中,辐射长度Lrad是材料的特征量,通常原子序数越大的材料,辐射长度越小.βc是缪子的速度,约等于光速,p是缪子的动量,L为材料的厚度.(1)式表明缪子散射角与缪子能量有关.天然缪子的能谱是连续的,实验测量到的缪子散射角是不同能量的缪子产生的散射角集合.单能的缪子散射角服从高斯分布,对天然缪子散射角分布的贡献为其能量对应的高斯概率密度函数,

由此推广到整个缪子能量区间,天然缪子散射角满足的分布应为全能量区间上的缪子产生的散射角贡献加权求和,权重为各微分能量点上出现的缪子概率.该分布理论上为高斯分布与缪子能谱的卷积.由于缪子能谱较为复杂,本研究采用缪子的9 个特征能量点代替缪子能谱.在天然缪子能量区间上以对数线性选取0.25,0.50,1.00,2.00,4.00,8.00,16.00,32.00,64.00 GeV 这9 个能量点,它们的跨度0.25—64.00 GeV 区间上包含了94.3%的天然缪子,覆盖了缪子的主要相互作用能量区间,具有理论意义.用9 个特征能量点上的缪子散射角分布加权求和近似连续能量缪子的散射角分布,对实验测量到的天然缪子散射角统计点进行拟合,可得

拟合函数(3)式为9 个能量点对应的高斯概率密度函数加权求和形式,常数项已归纳到待定系数中.通过测量标准材料的缪子散射角标定权重系数A1—A9,在已知材料种类时各能量点对应的σi为已知量,可以由(1)式计算.标定得到的权重A1—A9代表各自对应能量点上的缪子数量占比.将A1—A9代入(3)式得到一个近似天然缪子散射角分布的表达式.此后测量未知材料的缪子散射角并与标定后的(3)式耦合进而推算出未知材料的辐射长度.该方法从统计角度引入了缪子的离散能量信息.

综上所述,基于离散能量信息进行材料鉴别分为两个步骤: 1)利用已知厚度的铅作为标准材料标定各能量点的权重系数;2) 测量待鉴别样品的散射角,通过标定后的缪子散射角分布计算出各材料的辐射长度,并以此对材料进行鉴别.

3 实验配置与蒙特卡罗模拟

基于缪子离散能量的材料鉴别实验在中国原子能科学研究院研制的宇宙射线缪子成像装置上进行.该装置的主体部分由六层探测器阵列构成,每层由各两排垂直相交的气体漂移管紧密排列组成.上、下部分各三层分别用于测量入射和出射缪子的径迹坐标,中间约1 m 间距的探测区域摆放实验样品.本文采用Geant4 对实验过程进行模拟,Geant4 是欧洲核子中心开发的蒙特卡罗软件[23],用于模拟微观粒子与宏观物质相互作用的全过程.

3.1 离散能量点权重的实验标定

权重标定实验选用已知厚度的铅标准材料测量散射角并计算各离散能量所占权重.由于铅的厚度过大会阻挡低能缪子,厚度过小又使高能缪子产生的散射角过小,因此本实验选用5,10,15 cm 三种厚度的铅块,取标定结果的平均值作为最终权重.将三种厚度的铅块同时放入宇宙射线缪子成像装置中进行测量,实验设置和成像结果如图1 所示.

图1 三种不同厚度铅块及散射成像结果Fig.1.Lead cube with different thickness and the image of scattering tomography.

根据材料在探测区域内放置位置的先验知识筛选出穿过完整厚度铅块的缪子.选取其中穿透每种厚度的缪子各90000 个,计算上述缪子在XZ和YZ两正交平面上形成的平面散射角.两平面内的散射角独立、同分布,可以混合为一个样本集.在–200—200 mrad 的范围内,以1 mrad 为间隔划分区间,统计测量散射角的归一化频数后做出散点图并用(3)式进行拟合.在Geant4 中构建与上述实验相同的铅块模型,探测器阵列面积和位置均与实际装置一致.缪子从最顶层探测器表面1 m ×1 m 的正方形区域内均匀抽样发射.入射方向和缪子能量均由CRY 库抽样产生[24].模拟计算缪子穿透铅块产生的散射角分布,图2 为本实验的Geant4 模型示意图.

图2 权重标定实验的Geant4 模型Fig.2.Geant4 model for weight calibration experiment.

三种不同厚度铅块的散射角分布如图3(a)所示,其中10 cm 组的实验结果与模拟结果对比如图3(b)所示,实验结果已扣除本底影响.每个厚度的实验组可以拟合出一组权重A1—A9,将三组权重取平均值,得到的实验结果及模拟结果见表1,权重系数A1—A9代表以对应的特征能量点为中心的能量分段中包含的缪子数量占比.由表1 可知,A3—A5的数值在实验和模拟结果中均较高,表明在天然缪子中1—4 GeV 的缪子数量较多,对总体散射角分布贡献较大.由于实验测量存在一定误差,部分权重系数的实验与模拟结果存在差异.

表1 离散能量权重模拟结果与实验结果Table 1.Discrete energy’s weights of experiment and simulation.

图3 (a) 缪子穿过不同厚度铅块的散射角分布(点)与拟合曲线(线);(b) 10 cm 铅块的实验结果与模拟结果对比Fig.3.(a) Measured scattering angle distributions of lead cubes under different thicknesses (dot) and fitting curve (line);(b) the comparison between experiments and simulations of 10 cm Pb.

3.2 引入离散能量的材料鉴别实验与分析

材料鉴别实验选用C,Al,Fe,Pb,W 5 种样品作为待测材料,利用标定后的缪子散射角分布和各材料散射角测量值计算材料辐射长度.5 种材料的厚度均为10 cm,其中Al,Fe,Pb,W 为底面积10 cm × 10 cm 的立方体,C 为底面直径10 cm 的圆柱体.实验设置和成像结果如图4 所示.将5 种材料样品放入缪子成像装置中,经一段时间测量后,穿透每种材料完整厚度的粒子数各40000 个.在–200—200 mrad 内,以1 mrad 为间隔划分区间,统计各材料散射角的归一化频数并绘制散点图,如图5 所示.

图4 测量五种材料的散射角及散射成像结果Fig.4.Measurement of scattering angles for the different materials and the image of scattering tomography.

根据(3)式可知,离散能量点近似的天然缪子散射角分布在θ=0 时取得最大值,将该值与测量到的天然缪子散射角频率峰值N0(即图5 中的峰值)对应,

图5 缪子穿过不同材料的散射角分布(点)与拟合曲线(线)Fig.5.Measured scattering angle distribution of different materials (dot) and fitting curve (line).

从而计算出材料的辐射长度,

将穿过各材料的缪子事例等分为8 组做平行实验,按上述方法分别计算5 种材料在平行实验组的辐射长度值,对平行实验的结果取均值以减小随机误差.在Geant4 中构建C,Al,Fe,Pb,W 5 种材料模型,按照与实验相同步骤进行模拟计算.表2为5 种材料辐射长度的实验结果、模拟结果以及与劳伦斯伯克利实验室提供的辐射长度标准值[25]比较得到的相对误差.由表2 可知,重建5 种材料辐射长度的模拟结果与辐射长度标准值也存在一些偏差.由于模拟中不存在缪子径迹的测量误差,模拟结果的误差为重建材料辐射长度算法带来的方法误差.该方法误差的来源为两部分,一是用离散能量缪子散射角分布近似连续能量的缪子散射角分布产生的误差,二是对散射角测量值拟合的误差.方法误差通常在15%以内浮动.在实验结果中,Pb 和W 两种高原子序数材料的重建辐射长度误差较小,分别为4.7%和9.7%;而低原子序数材料受到实验测量误差的影响,重建辐射长度的误差较大,第4 节将进一步分析测量误差对计算不同材料辐射长度值的影响.

表2 不同材料辐射长度的估计值与相对误差Table 2.Estimated radiation lengths of different materials with relative error.

本文引入受试者工作特征曲线(ROC)评价两种材料之间鉴别的准确度,该曲线可反映各数据集整体的区分能力.其中AUC 值为ROC 曲线下方面积,AUC 值越接近1,代表各材料之间的鉴别准确度越高.本实验以材料辐射长度作为特征量,用二分类法对5 种材料中密度相邻的材料进行两两鉴别.将5 种材料分为C-Al,Al-Fe,Fe-Pb,Pb-W四个实验组,计算每种材料的辐射长度样本值30 个.待鉴别的两种材料的辐射长度样本混合后由小到大进行排序,从最小值到最大值之间等间隔选取阈值对混合样本逐个判断鉴别材料种类.对于单个阈值,材料分类准确率为准确分类的样本个数与混合样本总数之比.随着阈值的改变,分类准确率也不断变化.材料鉴别准确率定义为在所有阈值下分类准确率的最大值,对应的阈值为最优阈值.

绘制各实验组材料鉴别的ROC 曲线并计算曲线下面积AUC 值.四个实验组的ROC 曲线及材料鉴别准确率见图6 和表3.结果表明,在各样本值对应1400 个缪子事例下,四个实验组材料鉴别的AUC 值均在0.9 以上.其中Al-Fe 和Fe-Pb 在95%的置信水平下可以区分;C-Al,Pb-W 在85%的置信水平下可以区分.

图6 材料鉴别实验ROC 曲线Fig.6.ROC curve of material distinguishment experiment.

表3 各材料鉴别准确率Table 3.Distinguishment accuracy of the different materials.

材料鉴别的准确度也与两种材料间辐射长度的差值大小(材料密度差异)有关.本研究以Pb 材料为参照,模拟了Ag,Cd,Sn,Cu 四种辐射长度递增的材料与Pb 进行区分,各材料的鉴别准确率与辐射长度的关系如表4 所列.使用缪子离散能量的材料鉴别方法可在95%置信水平下区分辐射长度差值大于0.7 cm 的两种材料.

表4 材料鉴别准确率与辐射长度的关系Table 4.Materials distinguishment accuracy versus the radiation length.

对比分析引入缪子离散能量信息与不引入缪子能量信息对材料鉴别准确度的提升效果.在不引入缪子能量信息时,需要假定入射缪子能量为常数,材料鉴别方法只能利用缪子散射角大小作为判断标准.计算本实验中缪子穿过Pb 和W 散射角的平方平均值作为特征量,在不引入缪子能量信息下对Pb-W 进行材料鉴别,并与上文中引入缪子离散能量信息的Pb-W 鉴别进行对比,结果如图7 所示.引入离散缪子能量信息后,Pb-W 鉴别的AUC值与准确率分别相比前者提升了16.7%和18.5%,由此表明引入离散缪子能量信息可以提高重核材料鉴别的准确度.

图7 两种方法鉴别铅-钨的ROC 曲线对比Fig.7.ROC curves of the lead-tungsten distinguishment performed by different method.

4 径迹测量误差对鉴别准确度的影响

在实验测量过程中,由于探测器存在一定的位置分辨率,电子学读出以及径迹拟合等过程也存在一定误差,多种因素影响到测量的缪子坐标与缪子真实位置产生了偏差.为分析探测器径迹测量误差对材料鉴别准确度的影响,本节模拟计算径迹测量误差分别为0.5,1.0,1.5,2.0 mm 时5 种材料的辐射长度值.各材料辐射长度相对误差与径迹测量误差之间的关系如图8 所示.

图8 径迹测量误差与各材料辐射长度误差的关系Fig.8.Relationship between measurement error and the radiation length error of different materials.

图8 表明各材料的辐射长度误差随径迹测量误差的增加而增加.在相同水平的测量误差下,代入缪子离散能量计算低原子序数材料辐射长度的误差大于高原子序数材料辐射长度的误差.探测器的径迹测量误差影响材料辐射长度的计算精度,进而影响材料鉴别的准确度.随着测量误差的增大,C-Al 实验组和Pb-W 实验组的ROC 曲线变化如图9 所示.图9 表明,随着径迹测量误差的增加,两实验组的鉴别准确率均不断降低.与精准测量相比,当缪子径迹测量存在1 mm 的误差时,C-Al 鉴别的AUC 值降低了11.1%;Pb-W 鉴别的AUC值降低了9.8%.当径迹测量存在2 mm 的误差时,C-Al 鉴别的AUC 值降低了21.6%;Pb-W 鉴别的AUC 值降低了16.6%.由此表明径迹测量误差对高原子序数材料间鉴别的影响较小.

图9 不同径迹测量误差下鉴别ROC 曲线 (a) C-Al;(b) Pb-WFig.9.ROC curve of distinguishment under the different measurement error: (a) C-Al;(b) Pb-W.

5 总结

本文介绍了基于缪子离散能量进行材料鉴别的实验方法.基于缪子散射现象,通过已知厚度的铅块标定出离散能量的权重系数,测量样品的散射角计算得到各材料的辐射长度值,并以此为特征量进行不同材料的鉴别.Pb 和W 两种材料辐射长度的实验值与标准值相比,分别相差4.7%和9.7%.在约1400 个有效缪子事例下,分组对比实验表明,Al-Fe 和Fe-Pb 在95%的置信水平下可以区分.进一步模拟表明,引入缪子离散能量可在95%置信水平下区分辐射长度差值大于0.7 cm 的两种材料,并且相比于不引入缪子能量信息时,对Pb-W鉴别的准确率提高了18.5%.此外,高原子序数材料的鉴别准确度受测量误差的影响较小.