孫世峰

1)(華北電力大學(xué)核科學(xué)與工程學(xué)院,北京102206)

2)(華北電力大學(xué)非能動(dòng)核能安全技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室,北京102206)

(2020年5月7日收到;2020年6月4日收到修改稿)

1 引 言

基于X射線的分析技術(shù)是非常重要的材料表征方法,它們被廣泛用于元素空間分布、化學(xué)狀態(tài)和質(zhì)構(gòu)特性等的測(cè)量[1].其中,X射線熒光(XRF)技術(shù)因其具有對(duì)多元素樣品進(jìn)行無損分析的能力在許多科學(xué)領(lǐng)域和工業(yè)應(yīng)用中得到了迅速發(fā)展[1?4].在很多情況下,獲得高質(zhì)量的X射線圖像非常困難[5].因此,如何提高圖像質(zhì)量一直是相關(guān)學(xué)者研究的熱點(diǎn)[6].類似于可見光成像,X射線圖像的質(zhì)量在很大程度上取決于所采用的光學(xué)系統(tǒng).但是不同于可見光,X射線穿透性強(qiáng)且折射率近似為1,對(duì)其進(jìn)行反射、折射或聚焦都非常困難;用于X射線成像的光學(xué)器件,如反射鏡或透鏡等,必須滿足非常嚴(yán)格的要求,而且造價(jià)非常昂貴[7].

最簡(jiǎn)單的實(shí)現(xiàn)X射線成像的方法是利用小孔成像原理,設(shè)計(jì)基于單孔準(zhǔn)直器的針孔成像系統(tǒng).針孔的直徑越小,成像系統(tǒng)的空間分辨率越好;但直徑越小,系統(tǒng)的光收集效率越低,圖像的信噪比(SNR)也越差.為了解決上述問題,Dicke[8]和Ables[9]分別獨(dú)立提出了采用多開孔的準(zhǔn)直器(編碼孔徑)取代單一開孔的準(zhǔn)直器;編碼孔徑技術(shù)使系統(tǒng)在保留小孔徑對(duì)應(yīng)的高空間分辨率的同時(shí),大幅提高系統(tǒng)的光收集效率(正比于開孔個(gè)數(shù)),從而提高圖像信噪比.編碼孔徑技術(shù)提出之后得到了極大的發(fā)展,并廣泛用于X射線和g射線成像,包括天文觀測(cè)、醫(yī)學(xué)成像和輻射防護(hù)等[10?12].

傳統(tǒng)上,X射線熒光成像采用基于全反射原理制成的多毛細(xì)管作為光學(xué)聚焦器件,但是此類系統(tǒng)能夠收集的X射線的角度范圍過小.為提高系統(tǒng)的數(shù)值孔徑,Haboub等[13]探索將編碼孔徑技術(shù)用于X射線熒光成像,獲得了幾種標(biāo)準(zhǔn)樣品的實(shí)驗(yàn)圖像,初步驗(yàn)證了該技術(shù)的可行性.Kulow等[14]進(jìn)一步研究了近場(chǎng)成像時(shí)的圖像重建算法,并通過實(shí)驗(yàn)證明了采用編碼孔徑技術(shù)可以大幅提高采集到的光子通量,降低成像時(shí)間.另外,相比于小孔或多毛細(xì)管方案,編碼孔徑技術(shù)的一個(gè)很重要的優(yōu)勢(shì)是擴(kuò)展性更強(qiáng),系統(tǒng)可以在不犧牲分辨率的前提下,通過增加編碼準(zhǔn)直器的編碼數(shù)和面積,增大成像視野范圍.

目前,編碼孔徑技術(shù)應(yīng)用于X射線熒光成像已顯示出良好應(yīng)用前景,但仍存在采用傳統(tǒng)的卷積成像模型時(shí),近場(chǎng)成像偽影大、圖像分辨率較低等不足.近年來,得益于計(jì)算機(jī)視覺算法的進(jìn)步,編碼孔徑成像技術(shù)在光學(xué)成像領(lǐng)域取得了一些值得關(guān)注的進(jìn)展.DeWeert[15]和Asif[16]分別提出采用可分離編碼矩陣的編碼孔徑成像模型,可以極大地降低計(jì)算的復(fù)雜度,從而實(shí)現(xiàn)高編碼數(shù)下的高分辨成像.基于上述工作,Adams和Boominathan進(jìn)一步提出了Texas Two-Step(T2S)模型[17],獲得了高幀頻、分辨率達(dá)微米級(jí)的近場(chǎng)可見光圖像.T2S模型已被證明可用于近場(chǎng)高分辨g射線成像[18].本研究基于可分離編碼的T2S模型,進(jìn)一步對(duì)編碼孔徑技術(shù)應(yīng)用于X射線熒光成像進(jìn)行詳細(xì)的理論分析和模擬計(jì)算,期望為將來的實(shí)驗(yàn)研究和實(shí)際應(yīng)用提供參考依據(jù)和分析數(shù)據(jù).

2 編碼孔徑成像原理

編碼孔徑成像系統(tǒng)通常由編碼準(zhǔn)直器和位置靈敏探測(cè)器組成,兩者相互平行.編碼準(zhǔn)直器為多開孔結(jié)構(gòu),其上的開孔是按照某種規(guī)則排列,開孔部分允許射線通過而非開孔部分阻擋射線.如圖1所示,編碼孔徑成像過程可以分為兩部分:編碼投影和圖像重建.編碼過程即物體發(fā)射的射線經(jīng)過準(zhǔn)直器在探測(cè)器平面進(jìn)行投影,對(duì)于準(zhǔn)直器上的每一個(gè)開孔,其均會(huì)在探測(cè)器平面上產(chǎn)生一個(gè)一定比例的物體的倒立投影;編碼孔徑的多開孔,使探測(cè)器平面上記錄的是若干倒立投影的退化疊加像.圖像重建即計(jì)算機(jī)對(duì)探測(cè)器輸出的原始投影圖像進(jìn)行圖像解碼,以獲得物體的重建圖像.

假設(shè)物體X,探測(cè)器測(cè)量的投影圖像Y均為M×N的二維矩陣.對(duì)于任意的編碼準(zhǔn)直器,編碼孔徑成像可以用廣義線性模型表示為

其中,x和y為分別將X和Y的列連接并矢量化后的一維向量(長(zhǎng)度為MN);Φ為傳遞矩陣(大小為MN×MN),表示物體和投影圖像之間的線性關(guān)系;e是測(cè)量噪聲.對(duì)于任意給定的編碼準(zhǔn)直器和探測(cè)器,可以通過模擬射線從物體傳輸?shù)教綔y(cè)器或通過校準(zhǔn)過程來獲得不同幾何條件下的傳遞矩陣.但當(dāng)成像系統(tǒng)分辨率或編碼數(shù)較高時(shí),傳遞矩陣Φ將非常大,很難通過估算或校準(zhǔn)計(jì)算得到Φ,同時(shí),進(jìn)行圖像重建時(shí)的反問題求解也非常困難.比如,當(dāng)物體圖像和探測(cè)器的分辨率均為300×300時(shí),傳遞矩陣Φ將包含(9×104)×(9×104)=8.1×109個(gè)元素.

為減少從物體到測(cè)量投影的線性映射的復(fù)雜性,一般將編碼矩陣設(shè)計(jì)為可分離的,即二維編碼準(zhǔn)直器的圖樣由兩個(gè)一維序列的外積構(gòu)造而來[16,17].此時(shí),可以將式(1)中的成像模型重寫為

其中,P和Q分別為表示物體行和列的一維卷積的矩陣;E是測(cè)量噪聲和模型失配.當(dāng)物體圖像和探測(cè)器的分辨率均為300×300時(shí),P和Q各只有9×104個(gè)元素,而Φ中有8.1×109個(gè)元素.

在T2 S模型中,考慮到近場(chǎng)效應(yīng)的影響,探測(cè)器投影被進(jìn)一步分解為兩個(gè)可分離函數(shù)的疊加[17].對(duì)于一個(gè)二維的物體Xd,設(shè)其距離成像系統(tǒng)的距離為d,探測(cè)器的測(cè)量投影滿足

圖1編碼孔徑成像示意圖Fig.1.Schematics diagram of coded aperture imaging.

其中,Pod和 Qod為無準(zhǔn)直器時(shí)探測(cè)器響應(yīng)的傳遞矩陣;Pcd和 Qcd為準(zhǔn)直器編碼對(duì)應(yīng)的傳遞矩陣(下標(biāo)“o”代表開孔,“c”代表編碼,“d”代表不同距離d).

對(duì)于三維物體XD,可以將其離散為D個(gè)距離成像系統(tǒng)d的二維物體Xd,則測(cè)量投影滿足

為了減小各種噪聲源項(xiàng)對(duì)重建圖像的影響,同時(shí)避免過擬合,可以通過求解正則化的總體最小二乘問題來進(jìn)行圖像重建.采用FISTA(快速迭代閾值收縮算法)[19]求解Lasso問題對(duì)應(yīng)方程,可以得到重建的物體圖像

3 基于Geant4的成像模擬

3.1 成像系統(tǒng)參數(shù)設(shè)計(jì)

本工作基于Geant4軟件包[20], 選用QGSP_BIC物理模型,通過設(shè)計(jì)基于可分離編碼準(zhǔn)直器的X射線成像模擬程序,驗(yàn)證T2S模型用于高分辨X射線成像的有效性.編碼準(zhǔn)直器采用修正均勻冗余陣列(modified uniformly redundant arrays,MURA),選擇高編碼數(shù)(N=463)以使成像系統(tǒng)具有較高的分辨率.如圖2所示,選取原編碼排布中心區(qū)域的90×90個(gè)像素設(shè)計(jì)編碼準(zhǔn)直器,黑色部分為重金屬(鉑),白色部分為開孔,單位像素尺寸為25μm×25μm.編碼準(zhǔn)直器的厚度應(yīng)能夠阻擋大部分垂直入射的X射線,考慮到準(zhǔn)直效應(yīng)可能影響重建圖像的質(zhì)量,準(zhǔn)直器厚度設(shè)為25—200μm.參考目前pnCCD所用探測(cè)器[2],本工作模擬的探測(cè)器材料設(shè)置為硅, 靈敏體積為2 mm ×2 mm×0.5 mm, 單位像素尺寸為25μm ×25μm,總像素?cái)?shù)為80×80.探測(cè)器與編碼準(zhǔn)直器的距離設(shè)置為固定值2 mm.模擬程序記錄入射X射線在探測(cè)器上的沉積能量及對(duì)應(yīng)像素編號(hào)生成不同情況對(duì)應(yīng)的投影圖像.

圖2 模擬選用的可編碼準(zhǔn)直器圖樣Fig.2.The mask pattern used in the simulation.

3.2 模型校準(zhǔn)

對(duì)于給定的成像系統(tǒng),成像模型中的傳遞矩陣可以通過建立射線經(jīng)準(zhǔn)直器到探測(cè)器的模型求解或經(jīng)過校準(zhǔn)來得到.本工作通過模擬校準(zhǔn)過程,計(jì)算所需的所有傳遞矩陣.當(dāng)物體與系統(tǒng)之間距離固定為任一d時(shí),若物體平面僅第i行對(duì)應(yīng)的像素分布著放射源,即物體圖像可表示為Xi=ei1T.此時(shí)探測(cè)器測(cè)量到的投影圖像為

其中,poi和pci分別為Po和Pc的第i列;qo和qc分別為Qo和Qc的列的和.由于poi和pci是正交的,因此可以進(jìn)行截?cái)嗥娈愔捣纸?truncated singular value decomposition,TSVD)投影Yi,僅保留最大的兩個(gè)奇異值,分別賦值給poi和pci.通過逐行掃描物體平面的每一行,可以計(jì)算出Po和Pc的所有元素.類似地,再通過逐列掃描物體平面的所有列,可以計(jì)算出Qo和 Qc的所有元素.通過改變物體與系統(tǒng)之間距離d,重復(fù)上述掃描過程,可以獲得不同d對(duì)應(yīng)的傳遞矩陣.

對(duì)于上述成像系統(tǒng),采用長(zhǎng)度為2 mm,寬度為0.001 mm的線狀X射線源模擬校準(zhǔn)過程.參考相關(guān)研究工作[13,14],本文研究的熒光X射線能量范圍設(shè)置為5—20 keV,因此計(jì)算了線源能量為6、8、10和16 keV時(shí)對(duì)應(yīng)的傳遞矩陣.如果系統(tǒng)的參數(shù)保持不變,只需執(zhí)行一次校準(zhǔn)過程,對(duì)應(yīng)的傳遞矩陣可重復(fù)使用.

3.3 重建圖像評(píng)價(jià)

本工作主要采用2個(gè)二維平面物體測(cè)試系統(tǒng)的成像能力. 一是長(zhǎng)度為0.4 mm, 寬度為0.001 mm的線源;二是形狀為“180”的物體,考察系統(tǒng)對(duì)復(fù)雜物體的成像能力.兩個(gè)物體均設(shè)置在平行于成像系統(tǒng),且距離為10 mm的平面上.

為了比較不同條件下的重建圖像,選用均方根誤差(root mean square error,RMSE)和通用圖像質(zhì)量指標(biāo)(universal quality index,UQI)2個(gè)指標(biāo)定量評(píng)價(jià)圖像質(zhì)量.RMSE的定義為:

其中,P(i)為原始物體圖像的像素值;P’(i)為重建圖像的像素值;N是圖像總的像素?cái)?shù).RMSE值越大,表明重建圖像與原始圖像的差異越大.

UQI的定義為[21]:

其中,cov是協(xié)方差函數(shù);xˉ和yˉ 分別是原始圖像和重建圖像的均值;分別是原始圖像和重建圖像的方差.UQI的動(dòng)態(tài)范圍為[0,1],UQI的值越接近1,表明重建圖像與原始圖像越相似.

4 模擬結(jié)果分析

4.1 準(zhǔn)直器厚度的影響

采用傳統(tǒng)的卷積成像模型進(jìn)行近場(chǎng)成像時(shí),由于編碼準(zhǔn)直器厚度引起的準(zhǔn)直效應(yīng)會(huì)使重建圖像質(zhì)量嚴(yán)重下降[14].因此,本節(jié)研究采用T2S模型時(shí),準(zhǔn)直器厚度對(duì)重建圖像的影響.將X射線能量設(shè)置為10 keV,編碼準(zhǔn)直器的厚度設(shè)置為25、50、100和200μm等4個(gè)厚度.圖3(a)所示為0.4 mm線源的重建圖像,由于系統(tǒng)的空間分辨率有限,重建圖像在橫縱向均有一定程度的模糊.雖然線源沿橫向均勻分布,但其重建圖像沿橫向并不是均勻分布,且與準(zhǔn)直器厚度無關(guān)(圖3b).準(zhǔn)直器的厚度從25μm增加到200μm,系統(tǒng)沿縱向的分辨率并無明顯變化(圖3c).通過進(jìn)一步模擬點(diǎn)源(10 keV)成像, 計(jì)算得到準(zhǔn)直器的厚度為25、50、100和200 μm時(shí),系統(tǒng)的空間分辨率(半高全寬,FWHM)分別為66.5、65.8、64.9和64.2μm.可以看出,隨著準(zhǔn)直器厚度的增加,空間分辨率緩慢減小.作為對(duì)比,Haboub等人基于20μm準(zhǔn)直器研制的成像系統(tǒng)的空間分辨率約為72μm[13].因此,相比于傳統(tǒng)的卷積模型,采用T2 S模型時(shí)分辨率略優(yōu).

對(duì)于形狀復(fù)雜的成像物體,采用T2S模型時(shí),重建圖像比較準(zhǔn)確地還原了物體的原始分布情況(圖4b).而采用傳統(tǒng)的卷積模型近場(chǎng)成像時(shí),重建圖像通常存在較強(qiáng)的偽影[13,14].當(dāng)準(zhǔn)直器厚度從25μm增加到200μm,重建圖像無明顯變化.在后兩節(jié)的模擬中,準(zhǔn)直器的厚度設(shè)置為固定值25μm.

圖4 物體的原始圖像及重建圖像(a)原始圖像;(b)準(zhǔn)直器厚度25微米時(shí)的重建圖像Fig.4.Original image and reconstructed image of the object:(a)Original image;(b)reconstructed image when the mask thickness was 25 micron.

4.2 X射線能量的影響

X射線的能量越高,穿透能力越強(qiáng),散射也越嚴(yán)重;物體發(fā)射的X射線能量不同時(shí),產(chǎn)生的投影圖像不同.同樣,校準(zhǔn)時(shí)采用的X射線能量不同,計(jì)算得到的傳遞矩陣也不同.因此,本節(jié)研究物體發(fā)射的X射線能量和校準(zhǔn)時(shí)采用的X射線能量對(duì)成像結(jié)果的影響.

如圖5所示,當(dāng)校準(zhǔn)能量為10 keV時(shí),物體發(fā)射的X射線能量為6和8 keV對(duì)應(yīng)的重建圖像質(zhì)量較好;而當(dāng)物體能量增大為16 keV時(shí),重建圖像的質(zhì)量明顯變差.如圖6所示,當(dāng)校準(zhǔn)能量為16 keV時(shí),物體能量為16 keV對(duì)應(yīng)的重建圖像的質(zhì)量較好,而物體能量為6和8 keV時(shí),重建圖像均有一定程度的畸變.圖7展示了形狀為“180”的物體的重建圖像質(zhì)量的定量評(píng)價(jià)結(jié)果隨物體能量和校準(zhǔn)能量的變化情況.可以看出,物體能量不同時(shí),不同校準(zhǔn)能量對(duì)應(yīng)的重建圖像的質(zhì)量所能達(dá)到的最佳值相近,且均出現(xiàn)在校準(zhǔn)能量和物體能量相同時(shí);校準(zhǔn)能量和物體能量不同時(shí),重建圖像的質(zhì)量均有不同程度的變差.整體而言,校準(zhǔn)能量為10 keV時(shí),重建圖像質(zhì)量受物體能量變化的影響最小.

圖5校準(zhǔn)能量10 keV情況下物體不同能量時(shí)的重建圖像(a)6 keV;(b)8 keV;(c)16 keVFig.5.Reconstructed images of the object with different energies at a calibration energy of 10 keV:(a)6 keV;(b)8 keV;(c)16 keV.

圖6校準(zhǔn)能量16 keV情況下物體不同能量時(shí)的重建圖像(a)6 keV;(b)8 keV;(c)16 keVFig.6.Reconstructed images of the object with different energies at a calibration energy of 16 keV:(a)6 keV;(b)8 keV;(c)16 keV.

圖7 重建圖像質(zhì)量的定量評(píng)價(jià)結(jié)果隨物體能量和校準(zhǔn)能量的變化情況(a)RMSE;(b)UQIFig.7.Quantitative evaluation results of reconstructed images change with the object energy and calibration energy: (a) RMSE; (b) UQI.

圖8線源的三維重建圖像Fig.8.3D reconstructed image of the line source.

4.3 三維成像結(jié)果

在近場(chǎng)條件下,編碼孔徑成像系統(tǒng)可以在二維成像的基礎(chǔ)上提取三維和重聚焦信息.為了考察系統(tǒng)的三維成像能力,將兩個(gè)成像物體的能量設(shè)為10 keV,并放置于與系統(tǒng)距離為d=10.0 mm的平面上.如圖8所示,雖然模擬的線源是平面物體,但重建的線源出現(xiàn)在不同的軸向切片上,極大值出現(xiàn)在真實(shí)距離(d=10.0 mm)附近;成像系統(tǒng)的軸向分辨率較差(FWHM=1.1 mm).圖9是形狀為“180”的物體的三維重建圖像.類似于線源的重建結(jié)果,重建物體出現(xiàn)在不同的軸向切片上,并且在真實(shí)距離(d=10.0 mm)附近最清晰.因此,相比于傳統(tǒng)基于卷積模型的成像系統(tǒng)僅能實(shí)現(xiàn)二維成像,采用T2S模型的成像系統(tǒng)可以很容易地實(shí)現(xiàn)距離d的估算.

5 總結(jié)與討論

圖9復(fù)雜物體的三維重建圖像Fig.9.3D reconstructed image of the complex object.

本文利用Geant4蒙特卡羅仿真,研究了基于可分離編碼的新型編碼孔徑成像計(jì)算模型用于高分辨X射線熒光成像的可行性.模擬結(jié)果表明,基于該模型的成像系統(tǒng)能夠?qū)崿F(xiàn)對(duì)線源和形狀復(fù)雜物體的高質(zhì)量成像,系統(tǒng)的空間分辨率約為65μm.與采用傳統(tǒng)卷積模型的成像系統(tǒng)不同,該系統(tǒng)的性能在近場(chǎng)成像時(shí)并不受準(zhǔn)直效應(yīng)的影響.重建圖像的質(zhì)量同時(shí)受物體能量和校準(zhǔn)能量的影響,兩者差異越小,重建圖像的質(zhì)量越高.系統(tǒng)能夠從單次測(cè)量的二維投影圖像,正確地重建出物體與系統(tǒng)的距離,軸向空間分辨率約為1.1 mm.

成像系統(tǒng)的空間分辨率受多個(gè)因素的影響,包括編碼準(zhǔn)直器和探測(cè)器的像素尺寸、物體的放大倍數(shù)(取決于物體到編碼準(zhǔn)直器的距離與編碼準(zhǔn)直器到探測(cè)器的距離之比)、圖像重建算法等.提高系統(tǒng)空間分辨率的途徑包括減小編碼準(zhǔn)直器和探測(cè)器的像素尺寸、增大物體的放大倍數(shù)以及采用更復(fù)雜的圖像重建算法等.未來的工作將集中于優(yōu)化系統(tǒng)設(shè)計(jì)并開發(fā)新的重建算法,以進(jìn)一步提高系統(tǒng)的空間分辨率.