劉 拓,王洪金, 何赟澤, TOMASZ Chady, 許 健, 李 響

(1.湖南大學(xué) 電氣與信息工程學(xué)院，長沙 410082；2.天津大學(xué) 精密儀器與光電子工程學(xué)院，天津 300000；3.波蘭西波莫瑞工業(yè)大學(xué) 電氣工程系，什切青 70-310)

將太赫茲波用于樹脂玻璃材料的無損檢測研究，所研究的問題是如何利用圖像處理方法提高太赫茲成像分辨率。太赫茲波束廣泛應(yīng)用于非金屬材料以及非極性材料中，具有寬頻帶、低能量、高穿透等特性，在無損檢測領(lǐng)域有著重要的應(yīng)用[1]。太赫茲波束的高穿透特性在同等條件下比渦流攜帶了更多信息；而其光子低能的特性，又使其比X射線更難影響材料內(nèi)部的組織結(jié)構(gòu)[2]。太赫茲時(shí)域光譜(Terahertz Time-domain Spectroscopy, THz-TDS)則是材料與太赫茲脈沖相互作用后，測量得到的太赫茲場強(qiáng)隨著時(shí)間變動的曲線，其信噪比較高、穩(wěn)定性好，可以在室溫下工作[3]。在醫(yī)學(xué)領(lǐng)域，太赫茲圖像可以對微波、紅外、可見光、紫外線或X射線檢測圖像所獲得的信息進(jìn)行互補(bǔ)[4]；在安防領(lǐng)域，適用于安檢的太赫茲設(shè)備已經(jīng)逐步應(yīng)用在車站、機(jī)場等，與太赫茲設(shè)備相匹配的算法也在不斷地出現(xiàn)，不僅可以節(jié)省人力資源，還能降低安檢成本[5]；在軍事領(lǐng)域，太赫茲合成孔徑雷達(dá)在重量、體積、功耗等方面具有較大的優(yōu)勢，更適合于極大帶寬信號的產(chǎn)生與處理，更容易實(shí)現(xiàn)高空間分辨力的成像[6]。但另一方面，太赫茲成像系統(tǒng)的低分辨率一直是其應(yīng)用于無損檢測或其他領(lǐng)域的巨大障礙，大量相關(guān)的研究工作都致力于提高太赫茲系統(tǒng)的成像質(zhì)量[7]。筆者以樹脂玻璃為試驗(yàn)材料，研究了一種基于盲去卷積算法的圖像處理方法，通過圖像灰度范圍調(diào)整、盲去卷積、濾波去噪等步驟處理其太赫茲時(shí)域光譜，以提高其成像的分辨率，進(jìn)一步拓展了太赫茲在無損檢測領(lǐng)域的應(yīng)用范圍。

1 試驗(yàn)系統(tǒng)介紹

使用Picometrix公司的型號為T-ray 4000的脈沖太赫茲時(shí)域頻譜系統(tǒng)進(jìn)行逐點(diǎn)掃描成像。太赫茲時(shí)域頻譜系統(tǒng)工作原理示意如圖1所示。通過飛秒激光發(fā)射的脈沖分為兩支，其中一支通過光纖作用于太赫茲發(fā)射器中的太赫茲發(fā)射材料，促使其發(fā)射頻率約為1.8 THz的太赫茲波束。該系統(tǒng)為反射式，即太赫茲源和太赫茲探測器位于試驗(yàn)材料的同一側(cè)，波束從太赫茲源出發(fā)，因其波長在毫米到亞毫米級別，太赫茲波在大多數(shù)絕緣體中可無損失傳播或存在衰減傳播。當(dāng)太赫茲波遇到具有不同介電常數(shù)的兩個(gè)物體表面時(shí)，會在邊界界面的上下表面分別發(fā)生反射，根據(jù)分界界面與被測物體表面的距離，太赫茲場強(qiáng)峰值會在皮秒級的時(shí)域頻譜中發(fā)生變化，從而得知分界界面的存在。通過對被檢測物體表面進(jìn)行逐點(diǎn)掃描，最終能實(shí)現(xiàn)被檢測物體內(nèi)部的分層成像。太赫茲時(shí)域頻譜系統(tǒng)外觀如圖2所示。

圖2 太赫茲時(shí)域頻譜系統(tǒng)外觀

樹脂玻璃刻槽能夠在平底刻槽的底端形成一個(gè)良好的反射界面，因此文中采用樹脂玻璃刻槽來標(biāo)定基于太赫茲時(shí)域頻譜的成像系統(tǒng)的空間分辨率。試驗(yàn)材料為樹脂玻璃，厚度為5.08 mm。在其表面從左到右共設(shè)置12組人工缺陷，這些表面人工缺陷是利用數(shù)控機(jī)床銑削加工的矩形槽。對于該類型的樣本，目的不是模擬真正的缺陷，而是用于操作校準(zhǔn)、構(gòu)建校準(zhǔn)曲線、確認(rèn)探頭的設(shè)置或驗(yàn)證算法。在試驗(yàn)中，此類缺陷用于評估所用算法對最大分辨率或靈敏度的處理效果。每一組缺陷之間的距離為10 mm，每組兩個(gè)缺陷之間的距離逐漸增大，且所有缺陷的深度都是1 mm。缺陷的寬度為2.05 mm，缺陷的長度為5 mm。每組每兩個(gè)缺陷的距離逐漸增加(見表1)。試驗(yàn)材料結(jié)構(gòu)示意如圖3所示。

表1 每組缺陷之間的距離 mm

圖3 試驗(yàn)材料結(jié)構(gòu)示意

2 盲去卷積算法

實(shí)際的太赫茲成像結(jié)果可以看成是材料的樣品函數(shù)與太赫茲波束的點(diǎn)擴(kuò)散函數(shù)(Point Spread Function, PSF)的卷積，成像過程如式(1)所示[8]

T(x,y)=I(x,y)*F(x,y)+N(x,y)

(1)

式中：T(x,y)為觀測圖像；I(x,y)為真實(shí)圖像函數(shù)；F(x,y)為太赫茲波束的點(diǎn)擴(kuò)散函數(shù)；N(x,y)為圖像的噪聲函數(shù)。

顯然，當(dāng)太赫茲波束的點(diǎn)擴(kuò)散函數(shù)已知時(shí)，可以將其代入式(1)，由于T(x,y)也是已知量，所以可以反導(dǎo)出I(x,y)的數(shù)值。在太赫茲圖像中，I(x,y)是太赫茲波在物體中某一時(shí)刻的傳播方程在垂直于傳播路徑平面上的投影[7]。

然而，太赫茲時(shí)域頻譜在絕緣介質(zhì)中的傳播時(shí)空方程極為復(fù)雜，其不僅與光程有關(guān)，而且與太赫茲的波長有關(guān)。對于寬譜太赫茲時(shí)域頻譜儀得到的數(shù)據(jù)，使用物理擴(kuò)散模型解逆卷積則會非常復(fù)雜，產(chǎn)生不必要的莫爾條紋[7]，因此采用盲去卷積算法進(jìn)行處理。

盲去卷積算法的基本原理是：根據(jù)最大似然圖像模型以及其稀疏特性[9-10]，不斷迭代解卷積核，從而達(dá)到在不需要知道PSF具體數(shù)據(jù)的情況下，通過迭代以得到去卷積后的理想圖像與PSF具體的數(shù)值。該算法基于這樣一種假設(shè)：成像的光子滿足泊松分布，因此理想圖像將會是一個(gè)滿足泊松分布的模型。通過原先輸入圖像的先驗(yàn)概率，以及模糊圖像的先驗(yàn)概率，運(yùn)用大似然估計(jì)獲得原始清晰圖像的后驗(yàn)概率，通過這種方式讓復(fù)原圖像盡可能接近原始圖像[11]。假定理想圖像各像素點(diǎn)獨(dú)立，可以聯(lián)立所有像素點(diǎn)構(gòu)造關(guān)于所求圖像、所求PSF與當(dāng)前圖像的似然函數(shù)。該似然函數(shù)的值是從理想圖像推導(dǎo)出的當(dāng)前圖像的概率，當(dāng)此函數(shù)取得最大值時(shí)，就說明當(dāng)前求出的理想圖像和當(dāng)前求出的PSF能夠最大可能地從這一理想圖像推導(dǎo)出當(dāng)前的圖像。也就是說，能從當(dāng)前圖像最大程度地復(fù)原出理想的圖像。具體的處理步驟如下所述。

(1) 讀取圖像，并對圖像進(jìn)行預(yù)處理。因?yàn)镻SF收斂時(shí)，圖像一般已經(jīng)卷積過度，故需要提前設(shè)定迭代次數(shù)。

(2) 由于PSF未知，保守地估計(jì)其為一個(gè)所有值均相等的矩陣，以此矩陣作為初始PSF，以預(yù)處理后的圖像作為初始圖像。

(3) 使用最大似然估計(jì)的方法，即在頻率域中構(gòu)造卷積核作為PSF，使得當(dāng)前圖像可以反向?qū)ふ页隼硐雸D像的概率最大。

(4) 使用此PSF，反向迭代出新的當(dāng)前圖像。

(5) 使用新的當(dāng)前圖像，反向迭代出新的PSF。

(6) 重復(fù)第(4),(5)步的操作，直到達(dá)到預(yù)先設(shè)置的次數(shù)為止，輸出最終的PSF結(jié)果與圖像，完成整個(gè)盲去卷積的過程。

將這些步驟以流程圖的形式表示，如圖4所示。

圖4 盲去卷積算法流程

3 成像結(jié)果與處理

由于太赫茲波的特殊物理性質(zhì)，其在刻槽分界表面反射時(shí)會產(chǎn)生明顯的瑞利衍射。為了能夠直觀地觀測這些現(xiàn)象對太赫茲檢測結(jié)果的影響,檢測結(jié)果均以B掃描結(jié)果呈現(xiàn)(見圖5)。

圖5 材料的B掃描結(jié)果

圖5中，橫坐標(biāo)為當(dāng)前探頭的水平位置，且只截取了與缺陷處相關(guān)的坐標(biāo)，縱坐標(biāo)為時(shí)間軸，且只截取了與材料有關(guān)的時(shí)間范圍。此外，其為一幅偽彩色圖，顏色坐標(biāo)代表太赫茲的場強(qiáng)，為相對單位。A處的藍(lán)色線段傾斜的原因是材料沒有被水平放置，而B區(qū)域的一對色塊對應(yīng)于材料中的一對缺陷。色塊成對出現(xiàn)，且從左至右，每對色塊之間的距離越來越大，所有的色塊都在一條時(shí)間坐標(biāo)的直線上，這也與材料中缺陷的特征相吻合。C區(qū)域所在直線為材料的下表面，在垂直方向上，這條紅線淺色的部分都與中間的紅色色塊一一對應(yīng)。淺色的部分其實(shí)是該缺陷邊緣周圍的太赫茲衍射造成的，波束在入射到缺陷表面周圍的下表面時(shí)，會對缺陷產(chǎn)生影響，即雖然淺色部分對應(yīng)的區(qū)域?qū)嶋H上是空氣而不是表面，卻依然有場強(qiáng)存在。如果太赫茲系統(tǒng)的定向性足夠優(yōu)越，那么圖5中底部紅線所在區(qū)域只會出現(xiàn)斷續(xù)的紅色與背景色。

分別取橫坐標(biāo)47.5，45.5，153 mm，將3條曲線的谷值對齊，觀察試驗(yàn)材料對應(yīng)的A掃描結(jié)果(見圖6)。太赫茲探測器在每個(gè)點(diǎn)記錄從太赫茲脈沖發(fā)出后320 ps內(nèi)返回的太赫茲波在探測器內(nèi)引起的光電流所代表的電場強(qiáng)度波形圖，如圖6所示。儀器輸出的太赫茲場強(qiáng)度用太赫茲脈沖的場強(qiáng)作為參考，以相對場強(qiáng)的方式進(jìn)行記錄，因而以相對單位進(jìn)行表示。

圖6 試件2在不同x位置時(shí)的A掃描結(jié)果

圖6的橫坐標(biāo)為時(shí)間(進(jìn)行了一定的截取),縱坐標(biāo)為太赫茲場強(qiáng)。緊鄰場強(qiáng)波谷的波峰的物理意義為太赫茲波束從材料的上表面反射，而右起的波峰物理意義為太赫茲波束從材料的下表面反射。

黑色、紅色、黃色曲線分別對應(yīng)試驗(yàn)材料中遠(yuǎn)離缺陷處、缺陷縫隙處與缺陷處。由于黑色曲線遠(yuǎn)離缺陷，太赫茲波衍射的影響較小，位于兩個(gè)波峰之間，可以理解為波束穿過了試驗(yàn)材料；紅色曲線在缺陷縫隙處，與缺陷的距離接近，因此曲線在兩個(gè)表面峰值之間出現(xiàn)了衍射造成的峰值；黃色曲線沒有下表面，只有缺陷表面，同樣是由于衍射，在下表面的位置出現(xiàn)了一個(gè)本來不應(yīng)存在的波峰。

以圖5中藍(lán)色線段為基準(zhǔn)，將圖像對齊后僅裁剪缺陷部分的圖像。又由于圖像處理算法中一般不支持小于零的像素值，故亦需要調(diào)整圖像的灰度，將灰度轉(zhuǎn)化為0到255，調(diào)整后的圖像如圖7所示。

圖7 材料在缺陷附近的灰度圖像

圖7中紅色部位對應(yīng)反射最為強(qiáng)烈的區(qū)域，即分界層表面。從圖7可以看出，在太赫茲B掃描圖中，分界上下出現(xiàn)明顯的瑞利衍射現(xiàn)象；太赫茲的反射強(qiáng)度也并非迅速衰減，從而導(dǎo)致相近刻槽之間的邊界模糊，左起第二組刻槽間的間隙不可分辨，第三組刻槽間的間隙勉強(qiáng)可以區(qū)別。由于只是改變了圖像的灰度范圍，故圖像的形態(tài)沒有改變，此外，很難忽略缺陷周圍“X”狀的深藍(lán)色圖形，而該圖形由太赫茲波在有限尺寸透鏡中的傳播本質(zhì)引起[12]。又由于縱坐標(biāo)是時(shí)間軸，對應(yīng)深度尺寸，無需進(jìn)行去卷積，分別使用3×1、4×1與5×1大小的PSF對圖7進(jìn)行處理，迭代10次，結(jié)果如圖8所示。

圖8 材料使用不同大小PSF去卷積后的圖像

從圖8可以看出，3×1大小的初始PSF處理后的圖像不如4×1大小初始PSF處理后的圖像改善明顯，而5×1大小的初始PSF去卷積已經(jīng)過度，在缺陷的縫隙處出現(xiàn)了錯(cuò)誤的深紅色形狀。因此，最終采用4×1大小初始PSF對圖7進(jìn)行處理。也可以通過預(yù)估，設(shè)置初始PSF為4×1大小，來觀察圖像在不同迭代次數(shù)下去卷積的結(jié)果(見圖9)。

圖9 材料經(jīng)過不同次數(shù)迭代去卷積后的圖像

從圖9可以發(fā)現(xiàn)，去卷積過程中引入了不希望被看到的高頻噪聲。當(dāng)?shù)螖?shù)過少時(shí)，去卷積達(dá)不到應(yīng)有的效果，而在迭代次數(shù)過多時(shí)，缺陷部分的圖像則會發(fā)生混疊，高頻噪聲的影響也過大，難以被去除。因此，選擇初始PSF大小為4×1，迭代次數(shù)為10次。去卷積后，還要對圖像進(jìn)行濾波操作，削弱或消除高頻噪聲的影響，去卷積后的圖像濾波頻譜如圖10所示。

將去卷積后的圖像進(jìn)行快速傅里葉變換(FFT)， 圖像中心的位置是低頻部分，圖像的4個(gè)角是高頻部分。將圖像頻譜4個(gè)角處的值設(shè)置為0，僅僅保留圖像中心處的十字部分，就達(dá)到了低通濾波的效果，將此頻譜進(jìn)行反FFT，即可得到最后的圖像，將圖7與處理后的最終圖像進(jìn)行對比(見圖11)。

空間分辨率的比較是通過橫坐標(biāo)的變化來實(shí)現(xiàn)的，橫坐標(biāo)有變化就意味著分辨率有變化。由圖11可知，在對圖像進(jìn)行處理后，缺陷部分的分辨率明顯提高，每對缺陷的中間縫隙都更容易辨認(rèn)，說明處理效果良好。證實(shí)了盲去卷積的算法可以提高分辨率、改善圖像質(zhì)量。

圖10 去卷積后的圖像濾波頻譜

圖11 材料處理前后的缺陷處圖像

為了對圖像的清晰程度進(jìn)行評價(jià)和估計(jì)，也可以引入灰度平均梯度(GMG)這一指標(biāo)[13]。對于M×N大小的圖像，其定義如式(2)所示。

(2)

一般情況下，GMG與圖像的清晰度成正相關(guān)，可以較好地體現(xiàn)出圖像的邊緣特性。使用式(2)分別求取處理前后圖像的GMG，得到圖像在進(jìn)行處理之前的本指標(biāo)值為3.567 0，而在處理之后的本指標(biāo)值上升到了4.243 6。顯然，經(jīng)過這些處理步驟之后的圖像邊緣特性比處理之前的更好，處理過程有效地降低了圖像的模糊程度，縮小了在太赫茲B掃描圖像中可區(qū)分兩區(qū)域間的間隙。此外，將缺陷所在時(shí)間剖面處的處理前和處理后的曲線進(jìn)行對比(見圖12)。

圖12 材料處理前后的缺陷處時(shí)間剖面圖像

圖12是將圖11上下兩幅圖像缺陷處的時(shí)間坐標(biāo)剖開的結(jié)果，橫坐標(biāo)為探頭水平位置，縱坐標(biāo)為該處的灰度。且黑色曲線是圖像處理之前的缺陷剖面，紅色曲線是圖像處理之后的缺陷剖面。需要說明的是，圖12并非在圖11上的某一固定時(shí)間坐標(biāo)取值，而是依照以下原則：當(dāng)橫坐標(biāo)對應(yīng)的值為缺陷時(shí)，就直接取灰度峰值作為函數(shù)值；當(dāng)橫坐標(biāo)對應(yīng)的值沒有缺陷時(shí)，就取其臨近缺陷處的圖像縱坐標(biāo)作為縱坐標(biāo)值，而取該處橫縱坐標(biāo)的灰度。

通過圖12也可以清晰看到，紅線的每對峰值相較于黑線的每對峰值更容易分辨，處理后的曲線峰值較大，而兩峰之間鞍部處的谷值較小，尤其在左側(cè)，每對缺陷的谷值跟峰值都拉開了更大的差距，這進(jìn)一步證實(shí)了文中圖像處理方法有著良好的效果。

4 結(jié)語

利用以盲去卷積為主的太赫茲時(shí)域光譜處理方法，不需要得到太赫茲波束的點(diǎn)擴(kuò)散函數(shù)，只需要提供初始PSF的規(guī)模與迭代次數(shù)，具備良好的實(shí)用價(jià)值。將太赫茲時(shí)域光譜進(jìn)行了灰度范圍調(diào)整、去卷積處理、濾波去噪等操作，分別對處理前后的圖像進(jìn)行直觀對比、處理前后的圖像灰度平均梯度對比、處理前后的圖像缺陷時(shí)間剖面處對比，證實(shí)了算法能夠減弱太赫茲成像中由于瑞利散射引起的模糊，但是對于太赫茲傳播中的x-t(空間-時(shí)間)耦合分離的作用尚有待提高。