基于繆子離散能量的材料鑒別實(shí)驗(yàn)研究*

李雨芃 湯秀章 陳欣南 高春宇 陳雁南 范澄軍 呂建友

(中國(guó)原子能科學(xué)研究院,北京 102413)

宇宙射線繆子具有穿透力強(qiáng)、對(duì)重核材料敏感的特點(diǎn),近年來(lái)被廣泛應(yīng)用于核材料檢查等領(lǐng)域.繆子與不同原子序數(shù)材料發(fā)生的多重庫(kù)侖散射效果不同,利用這點(diǎn)可以對(duì)被測(cè)物體進(jìn)行成像以及材料鑒別,而在該過(guò)程中引入繆子的能量信息可以提高材料鑒別的準(zhǔn)確度.本文基于原子能院研制的繆子成像裝置開(kāi)展了5 種樣品的材料鑒別實(shí)驗(yàn),使用離散能量擬合近似連續(xù)能量的繆子散射角分布,進(jìn)而測(cè)量出各材料的輻射長(zhǎng)度并以此作為特征量進(jìn)行材料鑒別.實(shí)驗(yàn)結(jié)果表明,在約1400 個(gè)有效繆子事例下,Al-Fe 和Fe-Pb 可以在95%的置信水平下有效區(qū)分,該方法相比于不引入繆子能量信息對(duì)Pb-W 鑒別的準(zhǔn)確率提高了18.5%.

1 引言

繆子是一種由宇宙射線與大氣層相互作用產(chǎn)生的高能帶電粒子,質(zhì)量是電子的207 倍,平均能量為3—4 GeV[1].繆子的穿透能力強(qiáng),可以有效地穿過(guò)屏蔽層對(duì)物體內(nèi)部結(jié)構(gòu)進(jìn)行成像探測(cè),并且無(wú)放射性危害.目前繆子成像技術(shù)已被應(yīng)用于核材料檢查[2,3]、金字塔探測(cè)[4,5]、陵墓探測(cè)[6]、火山監(jiān)測(cè)[7,8],以及對(duì)反應(yīng)堆[9,10]和核廢料[11]的監(jiān)測(cè)等多個(gè)領(lǐng)域.

2003 年美國(guó)洛斯·阿拉莫斯國(guó)家實(shí)驗(yàn)室(LANL)首次提出利用繆子散射成像方法進(jìn)行核材料檢查[12],該技術(shù)利用繆子穿過(guò)物體時(shí)發(fā)生多重庫(kù)侖散射的角度偏轉(zhuǎn)與材料原子序數(shù)、繆子能量的規(guī)律對(duì)被測(cè)物進(jìn)行成像[13]和識(shí)別[14].在實(shí)際應(yīng)用中,天然繆子的能量是連續(xù)的且實(shí)時(shí)測(cè)量非常困難,因此大多數(shù)研究使用繆子的平均能量代替未知繆子能量,這種近似導(dǎo)致了成像圖像質(zhì)量和材料鑒別準(zhǔn)確度的降低.Vanini 等[15]開(kāi)展了不同能量精度下的高爐成像模擬研究,Bae 和Chatzidakis [16]開(kāi)展了多種特殊核材料的材料鑒別模擬研究,Morris 等[17]提出了繆子能量的多群模型等.這些模擬結(jié)果均表明繆子能量信息對(duì)于成像和材料鑒別有著重要影響.此外,部分學(xué)者也開(kāi)展了繆子能量測(cè)量的實(shí)驗(yàn)研究.例如,加拿大團(tuán)隊(duì)[18]研制出了帶能量測(cè)量結(jié)構(gòu)的CRIPT 探測(cè)器,通過(guò)測(cè)量已知材料阻擋繆子產(chǎn)生的散射角還原繆子動(dòng)量;清華大學(xué)羅志飛[19]利用多氣隙阻性板室(multi-gap resistive plate chamber,MRPC)時(shí)間分辨率高的特性,使用飛行時(shí)間法計(jì)算出繆子的分段能量信息.

本文基于原子能院研制的繆子成像裝置開(kāi)展了5 種不同原子序數(shù)樣品的材料鑒別實(shí)驗(yàn)研究,采用離散能量繆子散射角分布擬合近似連續(xù)能量的繆子散射角分布,以輻射長(zhǎng)度作為特征量進(jìn)行材料鑒別,并分析實(shí)驗(yàn)測(cè)量誤差對(duì)材料鑒別準(zhǔn)確度的影響.

2 基于繆子離散能量信息的材料鑒別方法

當(dāng)繆子穿過(guò)靶材料時(shí),會(huì)受到原子核庫(kù)侖力的作用發(fā)生多重庫(kù)侖散射,其散射角在同一平面內(nèi)的投影近似服從均值為0、標(biāo)準(zhǔn)差為σθ的高斯分布[20].根據(jù)Molière 理論[21,22],散射角分布寬度與繆子能量、材料輻射長(zhǎng)度有如下關(guān)系:

其中,輻射長(zhǎng)度Lrad是材料的特征量,通常原子序數(shù)越大的材料,輻射長(zhǎng)度越小.βc是繆子的速度,約等于光速,p是繆子的動(dòng)量,L為材料的厚度.(1)式表明繆子散射角與繆子能量有關(guān).天然繆子的能譜是連續(xù)的,實(shí)驗(yàn)測(cè)量到的繆子散射角是不同能量的繆子產(chǎn)生的散射角集合.單能的繆子散射角服從高斯分布,對(duì)天然繆子散射角分布的貢獻(xiàn)為其能量對(duì)應(yīng)的高斯概率密度函數(shù),

由此推廣到整個(gè)繆子能量區(qū)間,天然繆子散射角滿足的分布應(yīng)為全能量區(qū)間上的繆子產(chǎn)生的散射角貢獻(xiàn)加權(quán)求和,權(quán)重為各微分能量點(diǎn)上出現(xiàn)的繆子概率.該分布理論上為高斯分布與繆子能譜的卷積.由于繆子能譜較為復(fù)雜,本研究采用繆子的9 個(gè)特征能量點(diǎn)代替繆子能譜.在天然繆子能量區(qū)間上以對(duì)數(shù)線性選取0.25,0.50,1.00,2.00,4.00,8.00,16.00,32.00,64.00 GeV 這9 個(gè)能量點(diǎn),它們的跨度0.25—64.00 GeV 區(qū)間上包含了94.3%的天然繆子,覆蓋了繆子的主要相互作用能量區(qū)間,具有理論意義.用9 個(gè)特征能量點(diǎn)上的繆子散射角分布加權(quán)求和近似連續(xù)能量繆子的散射角分布,對(duì)實(shí)驗(yàn)測(cè)量到的天然繆子散射角統(tǒng)計(jì)點(diǎn)進(jìn)行擬合,可得

擬合函數(shù)(3)式為9 個(gè)能量點(diǎn)對(duì)應(yīng)的高斯概率密度函數(shù)加權(quán)求和形式,常數(shù)項(xiàng)已歸納到待定系數(shù)中.通過(guò)測(cè)量標(biāo)準(zhǔn)材料的繆子散射角標(biāo)定權(quán)重系數(shù)A1—A9,在已知材料種類時(shí)各能量點(diǎn)對(duì)應(yīng)的σi為已知量,可以由(1)式計(jì)算.標(biāo)定得到的權(quán)重A1—A9代表各自對(duì)應(yīng)能量點(diǎn)上的繆子數(shù)量占比.將A1—A9代入(3)式得到一個(gè)近似天然繆子散射角分布的表達(dá)式.此后測(cè)量未知材料的繆子散射角并與標(biāo)定后的(3)式耦合進(jìn)而推算出未知材料的輻射長(zhǎng)度.該方法從統(tǒng)計(jì)角度引入了繆子的離散能量信息.

綜上所述,基于離散能量信息進(jìn)行材料鑒別分為兩個(gè)步驟: 1)利用已知厚度的鉛作為標(biāo)準(zhǔn)材料標(biāo)定各能量點(diǎn)的權(quán)重系數(shù);2) 測(cè)量待鑒別樣品的散射角,通過(guò)標(biāo)定后的繆子散射角分布計(jì)算出各材料的輻射長(zhǎng)度,并以此對(duì)材料進(jìn)行鑒別.

3 實(shí)驗(yàn)配置與蒙特卡羅模擬

基于繆子離散能量的材料鑒別實(shí)驗(yàn)在中國(guó)原子能科學(xué)研究院研制的宇宙射線繆子成像裝置上進(jìn)行.該裝置的主體部分由六層探測(cè)器陣列構(gòu)成,每層由各兩排垂直相交的氣體漂移管緊密排列組成.上、下部分各三層分別用于測(cè)量入射和出射繆子的徑跡坐標(biāo),中間約1 m 間距的探測(cè)區(qū)域擺放實(shí)驗(yàn)樣品.本文采用Geant4 對(duì)實(shí)驗(yàn)過(guò)程進(jìn)行模擬,Geant4 是歐洲核子中心開(kāi)發(fā)的蒙特卡羅軟件[23],用于模擬微觀粒子與宏觀物質(zhì)相互作用的全過(guò)程.

3.1 離散能量點(diǎn)權(quán)重的實(shí)驗(yàn)標(biāo)定

權(quán)重標(biāo)定實(shí)驗(yàn)選用已知厚度的鉛標(biāo)準(zhǔn)材料測(cè)量散射角并計(jì)算各離散能量所占權(quán)重.由于鉛的厚度過(guò)大會(huì)阻擋低能繆子,厚度過(guò)小又使高能繆子產(chǎn)生的散射角過(guò)小,因此本實(shí)驗(yàn)選用5,10,15 cm 三種厚度的鉛塊,取標(biāo)定結(jié)果的平均值作為最終權(quán)重.將三種厚度的鉛塊同時(shí)放入宇宙射線繆子成像裝置中進(jìn)行測(cè)量,實(shí)驗(yàn)設(shè)置和成像結(jié)果如圖1 所示.

圖1 三種不同厚度鉛塊及散射成像結(jié)果Fig.1.Lead cube with different thickness and the image of scattering tomography.

根據(jù)材料在探測(cè)區(qū)域內(nèi)放置位置的先驗(yàn)知識(shí)篩選出穿過(guò)完整厚度鉛塊的繆子.選取其中穿透每種厚度的繆子各90000 個(gè),計(jì)算上述繆子在XZ和YZ兩正交平面上形成的平面散射角.兩平面內(nèi)的散射角獨(dú)立、同分布,可以混合為一個(gè)樣本集.在–200—200 mrad 的范圍內(nèi),以1 mrad 為間隔劃分區(qū)間,統(tǒng)計(jì)測(cè)量散射角的歸一化頻數(shù)后做出散點(diǎn)圖并用(3)式進(jìn)行擬合.在Geant4 中構(gòu)建與上述實(shí)驗(yàn)相同的鉛塊模型,探測(cè)器陣列面積和位置均與實(shí)際裝置一致.繆子從最頂層探測(cè)器表面1 m ×1 m 的正方形區(qū)域內(nèi)均勻抽樣發(fā)射.入射方向和繆子能量均由CRY 庫(kù)抽樣產(chǎn)生[24].模擬計(jì)算繆子穿透鉛塊產(chǎn)生的散射角分布,圖2 為本實(shí)驗(yàn)的Geant4 模型示意圖.

圖2 權(quán)重標(biāo)定實(shí)驗(yàn)的Geant4 模型Fig.2.Geant4 model for weight calibration experiment.

三種不同厚度鉛塊的散射角分布如圖3(a)所示,其中10 cm 組的實(shí)驗(yàn)結(jié)果與模擬結(jié)果對(duì)比如圖3(b)所示,實(shí)驗(yàn)結(jié)果已扣除本底影響.每個(gè)厚度的實(shí)驗(yàn)組可以擬合出一組權(quán)重A1—A9,將三組權(quán)重取平均值,得到的實(shí)驗(yàn)結(jié)果及模擬結(jié)果見(jiàn)表1,權(quán)重系數(shù)A1—A9代表以對(duì)應(yīng)的特征能量點(diǎn)為中心的能量分段中包含的繆子數(shù)量占比.由表1 可知,A3—A5的數(shù)值在實(shí)驗(yàn)和模擬結(jié)果中均較高,表明在天然繆子中1—4 GeV 的繆子數(shù)量較多,對(duì)總體散射角分布貢獻(xiàn)較大.由于實(shí)驗(yàn)測(cè)量存在一定誤差,部分權(quán)重系數(shù)的實(shí)驗(yàn)與模擬結(jié)果存在差異.

表1 離散能量權(quán)重模擬結(jié)果與實(shí)驗(yàn)結(jié)果Table 1.Discrete energy’s weights of experiment and simulation.

圖3 (a) 繆子穿過(guò)不同厚度鉛塊的散射角分布(點(diǎn))與擬合曲線(線);(b) 10 cm 鉛塊的實(shí)驗(yàn)結(jié)果與模擬結(jié)果對(duì)比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 引入離散能量的材料鑒別實(shí)驗(yàn)與分析

材料鑒別實(shí)驗(yàn)選用C,Al,Fe,Pb,W 5 種樣品作為待測(cè)材料,利用標(biāo)定后的繆子散射角分布和各材料散射角測(cè)量值計(jì)算材料輻射長(zhǎng)度.5 種材料的厚度均為10 cm,其中Al,Fe,Pb,W 為底面積10 cm × 10 cm 的立方體,C 為底面直徑10 cm 的圓柱體.實(shí)驗(yàn)設(shè)置和成像結(jié)果如圖4 所示.將5 種材料樣品放入繆子成像裝置中,經(jīng)一段時(shí)間測(cè)量后,穿透每種材料完整厚度的粒子數(shù)各40000 個(gè).在–200—200 mrad 內(nèi),以1 mrad 為間隔劃分區(qū)間,統(tǒng)計(jì)各材料散射角的歸一化頻數(shù)并繪制散點(diǎn)圖,如圖5 所示.

圖4 測(cè)量五種材料的散射角及散射成像結(jié)果Fig.4.Measurement of scattering angles for the different materials and the image of scattering tomography.

根據(jù)(3)式可知,離散能量點(diǎn)近似的天然繆子散射角分布在θ=0 時(shí)取得最大值,將該值與測(cè)量到的天然繆子散射角頻率峰值N0(即圖5 中的峰值)對(duì)應(yīng),

圖5 繆子穿過(guò)不同材料的散射角分布(點(diǎn))與擬合曲線(線)Fig.5.Measured scattering angle distribution of different materials (dot) and fitting curve (line).

從而計(jì)算出材料的輻射長(zhǎng)度,

將穿過(guò)各材料的繆子事例等分為8 組做平行實(shí)驗(yàn),按上述方法分別計(jì)算5 種材料在平行實(shí)驗(yàn)組的輻射長(zhǎng)度值,對(duì)平行實(shí)驗(yàn)的結(jié)果取均值以減小隨機(jī)誤差.在Geant4 中構(gòu)建C,Al,Fe,Pb,W 5 種材料模型,按照與實(shí)驗(yàn)相同步驟進(jìn)行模擬計(jì)算.表2為5 種材料輻射長(zhǎng)度的實(shí)驗(yàn)結(jié)果、模擬結(jié)果以及與勞倫斯伯克利實(shí)驗(yàn)室提供的輻射長(zhǎng)度標(biāo)準(zhǔn)值[25]比較得到的相對(duì)誤差.由表2 可知,重建5 種材料輻射長(zhǎng)度的模擬結(jié)果與輻射長(zhǎng)度標(biāo)準(zhǔn)值也存在一些偏差.由于模擬中不存在繆子徑跡的測(cè)量誤差,模擬結(jié)果的誤差為重建材料輻射長(zhǎng)度算法帶來(lái)的方法誤差.該方法誤差的來(lái)源為兩部分,一是用離散能量繆子散射角分布近似連續(xù)能量的繆子散射角分布產(chǎn)生的誤差,二是對(duì)散射角測(cè)量值擬合的誤差.方法誤差通常在15%以內(nèi)浮動(dòng).在實(shí)驗(yàn)結(jié)果中,Pb 和W 兩種高原子序數(shù)材料的重建輻射長(zhǎng)度誤差較小,分別為4.7%和9.7%;而低原子序數(shù)材料受到實(shí)驗(yàn)測(cè)量誤差的影響,重建輻射長(zhǎng)度的誤差較大,第4 節(jié)將進(jìn)一步分析測(cè)量誤差對(duì)計(jì)算不同材料輻射長(zhǎng)度值的影響.

表2 不同材料輻射長(zhǎng)度的估計(jì)值與相對(duì)誤差Table 2.Estimated radiation lengths of different materials with relative error.

本文引入受試者工作特征曲線(ROC)評(píng)價(jià)兩種材料之間鑒別的準(zhǔn)確度,該曲線可反映各數(shù)據(jù)集整體的區(qū)分能力.其中AUC 值為ROC 曲線下方面積,AUC 值越接近1,代表各材料之間的鑒別準(zhǔn)確度越高.本實(shí)驗(yàn)以材料輻射長(zhǎng)度作為特征量,用二分類法對(duì)5 種材料中密度相鄰的材料進(jìn)行兩兩鑒別.將5 種材料分為C-Al,Al-Fe,Fe-Pb,Pb-W四個(gè)實(shí)驗(yàn)組,計(jì)算每種材料的輻射長(zhǎng)度樣本值30 個(gè).待鑒別的兩種材料的輻射長(zhǎng)度樣本混合后由小到大進(jìn)行排序,從最小值到最大值之間等間隔選取閾值對(duì)混合樣本逐個(gè)判斷鑒別材料種類.對(duì)于單個(gè)閾值,材料分類準(zhǔn)確率為準(zhǔn)確分類的樣本個(gè)數(shù)與混合樣本總數(shù)之比.隨著閾值的改變,分類準(zhǔn)確率也不斷變化.材料鑒別準(zhǔn)確率定義為在所有閾值下分類準(zhǔn)確率的最大值,對(duì)應(yīng)的閾值為最優(yōu)閾值.

繪制各實(shí)驗(yàn)組材料鑒別的ROC 曲線并計(jì)算曲線下面積AUC 值.四個(gè)實(shí)驗(yàn)組的ROC 曲線及材料鑒別準(zhǔn)確率見(jiàn)圖6 和表3.結(jié)果表明,在各樣本值對(duì)應(yīng)1400 個(gè)繆子事例下,四個(gè)實(shí)驗(yàn)組材料鑒別的AUC 值均在0.9 以上.其中Al-Fe 和Fe-Pb 在95%的置信水平下可以區(qū)分;C-Al,Pb-W 在85%的置信水平下可以區(qū)分.

圖6 材料鑒別實(shí)驗(yàn)ROC 曲線Fig.6.ROC curve of material distinguishment experiment.

表3 各材料鑒別準(zhǔn)確率Table 3.Distinguishment accuracy of the different materials.

材料鑒別的準(zhǔn)確度也與兩種材料間輻射長(zhǎng)度的差值大小(材料密度差異)有關(guān).本研究以Pb 材料為參照,模擬了Ag,Cd,Sn,Cu 四種輻射長(zhǎng)度遞增的材料與Pb 進(jìn)行區(qū)分,各材料的鑒別準(zhǔn)確率與輻射長(zhǎng)度的關(guān)系如表4 所列.使用繆子離散能量的材料鑒別方法可在95%置信水平下區(qū)分輻射長(zhǎng)度差值大于0.7 cm 的兩種材料.

表4 材料鑒別準(zhǔn)確率與輻射長(zhǎng)度的關(guān)系Table 4.Materials distinguishment accuracy versus the radiation length.

對(duì)比分析引入繆子離散能量信息與不引入繆子能量信息對(duì)材料鑒別準(zhǔn)確度的提升效果.在不引入繆子能量信息時(shí),需要假定入射繆子能量為常數(shù),材料鑒別方法只能利用繆子散射角大小作為判斷標(biāo)準(zhǔn).計(jì)算本實(shí)驗(yàn)中繆子穿過(guò)Pb 和W 散射角的平方平均值作為特征量,在不引入繆子能量信息下對(duì)Pb-W 進(jìn)行材料鑒別,并與上文中引入繆子離散能量信息的Pb-W 鑒別進(jìn)行對(duì)比,結(jié)果如圖7 所示.引入離散繆子能量信息后,Pb-W 鑒別的AUC值與準(zhǔn)確率分別相比前者提升了16.7%和18.5%,由此表明引入離散繆子能量信息可以提高重核材料鑒別的準(zhǔn)確度.

圖7 兩種方法鑒別鉛-鎢的ROC 曲線對(duì)比Fig.7.ROC curves of the lead-tungsten distinguishment performed by different method.

4 徑跡測(cè)量誤差對(duì)鑒別準(zhǔn)確度的影響

在實(shí)驗(yàn)測(cè)量過(guò)程中,由于探測(cè)器存在一定的位置分辨率,電子學(xué)讀出以及徑跡擬合等過(guò)程也存在一定誤差,多種因素影響到測(cè)量的繆子坐標(biāo)與繆子真實(shí)位置產(chǎn)生了偏差.為分析探測(cè)器徑跡測(cè)量誤差對(duì)材料鑒別準(zhǔn)確度的影響,本節(jié)模擬計(jì)算徑跡測(cè)量誤差分別為0.5,1.0,1.5,2.0 mm 時(shí)5 種材料的輻射長(zhǎng)度值.各材料輻射長(zhǎng)度相對(duì)誤差與徑跡測(cè)量誤差之間的關(guān)系如圖8 所示.

圖8 徑跡測(cè)量誤差與各材料輻射長(zhǎng)度誤差的關(guān)系Fig.8.Relationship between measurement error and the radiation length error of different materials.

圖8 表明各材料的輻射長(zhǎng)度誤差隨徑跡測(cè)量誤差的增加而增加.在相同水平的測(cè)量誤差下,代入繆子離散能量計(jì)算低原子序數(shù)材料輻射長(zhǎng)度的誤差大于高原子序數(shù)材料輻射長(zhǎng)度的誤差.探測(cè)器的徑跡測(cè)量誤差影響材料輻射長(zhǎng)度的計(jì)算精度,進(jìn)而影響材料鑒別的準(zhǔn)確度.隨著測(cè)量誤差的增大,C-Al 實(shí)驗(yàn)組和Pb-W 實(shí)驗(yàn)組的ROC 曲線變化如圖9 所示.圖9 表明,隨著徑跡測(cè)量誤差的增加,兩實(shí)驗(yàn)組的鑒別準(zhǔn)確率均不斷降低.與精準(zhǔn)測(cè)量相比,當(dāng)繆子徑跡測(cè)量存在1 mm 的誤差時(shí),C-Al 鑒別的AUC 值降低了11.1%;Pb-W 鑒別的AUC值降低了9.8%.當(dāng)徑跡測(cè)量存在2 mm 的誤差時(shí),C-Al 鑒別的AUC 值降低了21.6%;Pb-W 鑒別的AUC 值降低了16.6%.由此表明徑跡測(cè)量誤差對(duì)高原子序數(shù)材料間鑒別的影響較小.

圖9 不同徑跡測(cè)量誤差下鑒別ROC 曲線 (a) C-Al;(b) Pb-WFig.9.ROC curve of distinguishment under the different measurement error: (a) C-Al;(b) Pb-W.

5 總結(jié)

本文介紹了基于繆子離散能量進(jìn)行材料鑒別的實(shí)驗(yàn)方法.基于繆子散射現(xiàn)象,通過(guò)已知厚度的鉛塊標(biāo)定出離散能量的權(quán)重系數(shù),測(cè)量樣品的散射角計(jì)算得到各材料的輻射長(zhǎng)度值,并以此為特征量進(jìn)行不同材料的鑒別.Pb 和W 兩種材料輻射長(zhǎng)度的實(shí)驗(yàn)值與標(biāo)準(zhǔn)值相比,分別相差4.7%和9.7%.在約1400 個(gè)有效繆子事例下,分組對(duì)比實(shí)驗(yàn)表明,Al-Fe 和Fe-Pb 在95%的置信水平下可以區(qū)分.進(jìn)一步模擬表明,引入繆子離散能量可在95%置信水平下區(qū)分輻射長(zhǎng)度差值大于0.7 cm 的兩種材料,并且相比于不引入繆子能量信息時(shí),對(duì)Pb-W鑒別的準(zhǔn)確率提高了18.5%.此外,高原子序數(shù)材料的鑒別準(zhǔn)確度受測(cè)量誤差的影響較小.