基于時(shí)間反轉(zhuǎn)聲學(xué)理論的超聲相控陣超分辨率成像方法*

樊程廣，余孫全，趙 勇，楊 磊

(國(guó)防科技大學(xué) 空天科學(xué)學(xué)院， 湖南 長(zhǎng)沙 410073)

超聲相控陣是一種新型的超聲傳感器，包含多個(gè)陣元，每一個(gè)陣元可以發(fā)射和接收聲波信號(hào)。通過(guò)計(jì)算機(jī)技術(shù)控制陣列中各個(gè)陣元激勵(lì)脈沖的時(shí)間延遲，可以靈活控制合成波陣面的偏轉(zhuǎn)和聚焦，達(dá)到掃描成像的目的。因此，將超聲相控陣置于被測(cè)對(duì)象表面的固定檢測(cè)位置，可以對(duì)被測(cè)對(duì)象內(nèi)部的某一個(gè)區(qū)域進(jìn)行掃描成像，其已被廣泛應(yīng)用于工業(yè)無(wú)損檢測(cè)領(lǐng)域[1-3]。利用超聲相控陣對(duì)被測(cè)對(duì)象進(jìn)行檢測(cè)時(shí)，首先，采用全矩陣采集(Full Matrix Capture, FMC)方法在線獲取超聲陣列數(shù)據(jù)。FMC方法利用陣列的每一對(duì)陣元組合獲取超聲回波信號(hào)，因此數(shù)據(jù)中包含了與被測(cè)對(duì)象相關(guān)的最大量信息，同時(shí)數(shù)據(jù)可以被不同的成像方法處理，避免了重復(fù)采集，降低了實(shí)驗(yàn)成本[4-5]；其次，利用不同成像方法對(duì)超聲陣列數(shù)據(jù)進(jìn)行離線處理，得到超聲圖像。成像方法是時(shí)域信號(hào)轉(zhuǎn)換為圖像的“橋梁”，相同的數(shù)據(jù)被不同的成像方法處理，得到不同質(zhì)量的超聲圖像；最后，提取超聲圖像特征，評(píng)估被測(cè)對(duì)象內(nèi)部狀況。分辨率是重要的評(píng)價(jià)指標(biāo)，表征成像方法區(qū)分被測(cè)對(duì)象包含相鄰點(diǎn)目標(biāo)的能力。

瑞利準(zhǔn)則定義了常規(guī)成像方法的分辨率極限，若被測(cè)對(duì)象中相鄰點(diǎn)目標(biāo)的距離小于分辨率極限時(shí)，根據(jù)超聲圖像無(wú)法區(qū)分相鄰點(diǎn)目標(biāo)。其中，分辨率極限與超聲波工作波長(zhǎng)相關(guān)，為了提高成像分辨率，需要使用更高頻率的超聲波；但是，超聲波的工作頻率越高，其在介質(zhì)中傳播時(shí)的衰減越厲害，影響探測(cè)深度[6]。如何在保持超聲波工作頻率不變且不影響系統(tǒng)探測(cè)深度的前提下，提高成像分辨率，即為本文研究的超分辨率成像。

Fink等基于聲學(xué)波動(dòng)方程的時(shí)間反轉(zhuǎn)不變性，提出了時(shí)間反轉(zhuǎn)聲學(xué)(Time Reversal Acoustics, TRA)理論，并將其應(yīng)用于無(wú)損檢測(cè)領(lǐng)域[7]。Lev-Ari等在時(shí)間反轉(zhuǎn)聲學(xué)理論的基礎(chǔ)上，結(jié)合多信號(hào)分類(MUltiple SIgnal Classification, MUSIC)法，提出了單陣元脈沖回波工作模式下的時(shí)間反轉(zhuǎn)多信號(hào)分類(Time Reversal-MUltiple SIgnal Classification, TR-MUSIC)法[8]，并隨后擴(kuò)展至雙陣元傳輸工作模式[9]。利用TR-MUSIC法對(duì)介質(zhì)中相鄰點(diǎn)目標(biāo)進(jìn)行成像，在理論模型中考慮聲波在不同目標(biāo)之間的多次散射(Multiple Scattering, MS)，仿真結(jié)果驗(yàn)證了TR-MUSIC對(duì)相鄰點(diǎn)目標(biāo)的超分辨率特性[10-11]。Simonetti等通過(guò)理論分析指出，聲波在不同目標(biāo)之間的MS包含了目標(biāo)的次波長(zhǎng)信息，因此，MS是實(shí)現(xiàn)超分辨率成像的關(guān)鍵因素，但是，de Rosny等對(duì)該結(jié)論提出了異議，故而MS在超分辨率成像中的具體作用有待進(jìn)一步研究[12-13]。由于TR-MUSIC具備超分辨率特性，其已經(jīng)在醫(yī)學(xué)成像領(lǐng)域得到了研究，并被應(yīng)用于乳腺癌細(xì)胞的成像和定位[14-17]。此外，由于醫(yī)學(xué)成像領(lǐng)域和無(wú)損檢測(cè)領(lǐng)域存在差異[18]，需要進(jìn)一步研究TR-MUSIC的超分辨率特性及其在無(wú)損檢測(cè)領(lǐng)域的應(yīng)用。

1 超聲相控陣檢測(cè)原理

利用超聲相控陣對(duì)被測(cè)對(duì)象進(jìn)行檢測(cè)，其原理如圖1所示。在笛卡爾坐標(biāo)系中，一維線性陣列的中心位置為坐標(biāo)原點(diǎn)o，其所包含的N個(gè)陣元沿x軸方向均勻分布，陣元寬度為a。激勵(lì)陣元發(fā)出的超聲波在x-z平面?zhèn)鞑ィ山邮贞囋邮眨瑇-z平面中包含理想點(diǎn)散射體。

圖1 超聲相控陣檢測(cè)原理圖Fig.1 Detection schematic diagram of ultrasonic phased array

對(duì)于陣列的任意一組陣元組合Rtx-Rrx，超聲回波信號(hào)的頻域表達(dá)如式(1)所示[19]。

Htx,rx(ω)=F(ω)Dtx(θtx,ω)G(Rtx,rs,ω)

fsG(rs,Rrx,ω)Drx(θrx,ω)

(1)

式中，F(xiàn)(ω)是陣元輸出信號(hào)頻譜，fs為理想點(diǎn)散射體的散射系數(shù)。

Dtx(θtx,ω)和Drx(θrx,ω)是陣元在固體介質(zhì)中的遠(yuǎn)場(chǎng)指向性函數(shù)，其定義分別如式(2)～(3)所示。

(2)

(3)

其中，λ(ω)為超聲波波長(zhǎng)，cL和cS分別為超聲縱波和超聲橫波在介質(zhì)中的傳播速度，θtx和θrx的定義如圖1所示。F0函數(shù)的定義如式(4)所示。

F0(ξ)=[2ξ2-(cL/cS)2]2-

4ξ2(ξ2-1)1/2[ξ2-(cL/cS)2]1/2

(4)

G(Rtx,rs,ω)和G(rs,Rrx,ω)為介質(zhì)的格林函數(shù)，定義如式(5)所示。

(5)

利用FMC方法獲取超聲陣列數(shù)據(jù)，如圖2所示。由于FMC方法利用陣列的任意一對(duì)陣元組合獲取超聲回波信號(hào)，因此，由其獲取的超聲陣列數(shù)據(jù)包含N2個(gè)信號(hào)。

圖2 超聲陣列數(shù)據(jù)和陣列響應(yīng)矩陣Fig.2 Ultrasonic array data and array response matrix

2 時(shí)間反轉(zhuǎn)多信號(hào)分類法數(shù)學(xué)定義

2.1 陣列響應(yīng)矩陣

陣列響應(yīng)矩陣描述了超聲相控陣檢測(cè)系統(tǒng)中不同陣元之間的脈沖響應(yīng)。對(duì)于一個(gè)工作在回波模式的N元陣列，其任意陣元均可用于激勵(lì)和接收聲波信號(hào)，依次激勵(lì)陣列中的陣元，由其引起的聲波信號(hào)同時(shí)被陣列的N個(gè)陣元接收，即可得到陣列響應(yīng)矩陣K[20]。在實(shí)際應(yīng)用過(guò)程中，通常通過(guò)超聲陣列數(shù)據(jù)構(gòu)建陣列響應(yīng)矩陣K，其計(jì)算過(guò)程如圖2所示。給定工作頻率，通常選擇陣列中心頻率ωc，提取其對(duì)應(yīng)的H(ωc)值，作為矩陣K(ωc)的元素，構(gòu)建陣列響應(yīng)矩陣，其為N×N復(fù)矩陣。

2.2 奇異值分解和子空間劃分

TR-MUSIC法通過(guò)對(duì)陣列響應(yīng)矩陣進(jìn)行奇異值分解(Singular Value Decomposition, SVD)實(shí)現(xiàn)成像目的[11, 21]。

K(ω)的奇異值分解方程如式(6)所示。

K(ω)=U(ω)Σ(ω)VH(ω)

(6)

式中，U(ω)和V(ω)由奇異向量組成，Σ(ω)由奇異值組成，上標(biāo)H表示復(fù)矩陣的共軛轉(zhuǎn)置。

根據(jù)奇異值分布特征，將奇異向量(μi和νi)劃分為信號(hào)子空間(US(ω),VS(ω))和噪聲子空間(UN(ω),VN(ω))，如式(7)～(8)所示。

U(ω)=[US(ω)|UN(ω)]

=[μ1(ω),…,μm(ω),|μm+1(ω),…,μN(yùn)(ω)]

(7)

V(ω)=[VS(ω)|VN(ω)]

=[ν1(ω),…,νm(ω),|νm+1(ω),…,νN(ω)]

(8)

其中，m為信號(hào)子空間的維度。

2.3 TR-MUSIC成像函數(shù)

給定成像區(qū)域，對(duì)于成像區(qū)域中的任意像素點(diǎn)r，定義方向向量g(r,ω)，如式(9)所示。

g(r,ω)=[G(R1,r,ω),…,G(RN,r,ω)]T

(9)

式中，Rl(l=1～N)表示陣元位置，上標(biāo)T表示復(fù)矩陣的轉(zhuǎn)置。

在陣列中心頻率點(diǎn)ωc處，根據(jù)方向向量和信號(hào)子空間，定義TR-MUSIC成像函數(shù)如式(10)所示。

(10)

式中，‖·‖表示范數(shù)。

2.4 點(diǎn)擴(kuò)散函數(shù)

點(diǎn)擴(kuò)散函數(shù)(Point Spread Function, PSF)描述超聲成像系統(tǒng)對(duì)理想點(diǎn)散射體的響應(yīng)[19]。給定陣列和理想點(diǎn)散射體，利用FMC方法獲取超聲陣列數(shù)據(jù)，通過(guò)成像方法得到理想點(diǎn)散射體的超聲圖像，提取圖像特征用于評(píng)價(jià)不同成像方法。

數(shù)值仿真參數(shù)設(shè)置為：一維線性陣列包含64個(gè)均勻分布的陣元，所有陣元可以激勵(lì)和接收超聲波，超聲波在介質(zhì)中的傳播速度為6 300 m/s，成像區(qū)域?yàn)?0 mm×20 mm，在成像區(qū)域中有一個(gè)理想點(diǎn)散射體，其位置坐標(biāo)為(x=0,z=20 mm)，其散射系數(shù)定義為1。利用FMC方法獲取超聲陣列數(shù)據(jù)，由于線性陣列包含64個(gè)陣元，因此，超聲陣列數(shù)據(jù)共包含4 096個(gè)時(shí)域信號(hào)。利用平面B掃描方法和TR-MUSIC方法處理超聲陣列數(shù)據(jù)，得到理想點(diǎn)散射體的二維和三維超聲圖像，如圖3所示。平面B掃描是常用的基于延遲和疊加原理(Delay And Sum, DAS)的波束形成技術(shù)，由于其操作簡(jiǎn)單，對(duì)噪聲具備魯棒性，已經(jīng)被廣泛用于工業(yè)無(wú)損檢測(cè)領(lǐng)域[4]。對(duì)于TR-MUSIC方法，利用其對(duì)介質(zhì)中的理想點(diǎn)散射體成像時(shí)，在不考慮噪聲的前提下，其信號(hào)子空間的維度與理想點(diǎn)散射體的數(shù)目一致[22]，因此，在仿真過(guò)程中，信號(hào)子空間的維度選擇為1。

圖3 點(diǎn)擴(kuò)散函數(shù)Fig.3 Point spread function

由圖3可知，基于平面B掃描和TR-MUSIC圖像均可定位理想點(diǎn)散射體的位置。相較于平面B掃描方法，TR-MUSIC方法的成像結(jié)果中的理想點(diǎn)散射體的圖像區(qū)域明顯變小，說(shuō)明其成像分辨率得到提高。

3 實(shí)驗(yàn)驗(yàn)證

3.1 實(shí)驗(yàn)系統(tǒng)

實(shí)驗(yàn)系統(tǒng)原理圖如圖4所示，主要包括陣列控制器、線性陣列傳感器以及測(cè)試對(duì)象。陣列控制器的一端與PC機(jī)相連，一端與線性陣列相連，陣列通過(guò)耦合劑放置在測(cè)試對(duì)象表面。PC機(jī)發(fā)出指令，通過(guò)陣列控制器控制陣列傳感器中各個(gè)陣元的激勵(lì)和接收。

圖4 實(shí)驗(yàn)系統(tǒng)原理Fig.4 Principle of experimental system

圖5 實(shí)驗(yàn)設(shè)備Fig.5 Experimental equipment

實(shí)驗(yàn)設(shè)備實(shí)物圖如圖5所示，其中陣列控制器選擇英國(guó)Peak NDT公司的MicroPulse 5PA，其參數(shù)見(jiàn)表1，它包括128個(gè)獨(dú)立通道，用于進(jìn)行數(shù)據(jù)傳輸，其中每一個(gè)通道可以進(jìn)行16位的數(shù)據(jù)采集。線性陣列傳感器選擇法國(guó)Imasonic公司的5 MHz、64元線性陣列，其參數(shù)見(jiàn)表2。

表1 陣列控制器參數(shù)

表2 線性陣列參數(shù)

測(cè)試對(duì)象選擇不銹鋼試塊，在其內(nèi)部加工6個(gè)直徑為1 mm的貫通孔。超聲波在不銹鋼試塊中的傳播速度為6 150 m/s，在陣列傳感器中心頻率5 MHz處，其對(duì)應(yīng)的波長(zhǎng)為1.23 mm，大于貫通孔直徑1 mm，因此，這6個(gè)貫通孔可視為點(diǎn)散射體。

3.2 實(shí)驗(yàn)結(jié)果分析

首先，利用FMC方法獲取超聲陣列數(shù)據(jù)。由于超聲波在固體介質(zhì)中傳播時(shí)，存在超聲縱波、橫波等，因此，對(duì)于超聲陣列數(shù)據(jù)中的每一個(gè)時(shí)域信號(hào)，需要從中提取與缺陷相關(guān)的散射信號(hào)，選擇在時(shí)域加窗函數(shù)的方法進(jìn)行提取，部分原始信號(hào)和提取的散射信號(hào)如圖6所示。

其次，對(duì)預(yù)處理后的時(shí)域信號(hào)進(jìn)行時(shí)域-頻域轉(zhuǎn)換，在陣列傳感器的中心頻率點(diǎn)處，構(gòu)建陣列響應(yīng)矩陣。由于陣列包含64個(gè)陣元，因此，陣列響應(yīng)矩陣為64×64的復(fù)矩陣。

再次，對(duì)陣列響應(yīng)矩陣進(jìn)行奇異值分解，在陣列中心頻率點(diǎn)處，得到64個(gè)奇異值，其分布曲線如圖7所示。在64個(gè)奇異值中，前11個(gè)奇異值大于剩下的奇異值，因此，信號(hào)子空間的維度設(shè)定為11。信號(hào)子空間維度會(huì)影響超聲成像結(jié)果質(zhì)量。信號(hào)子空間維度過(guò)小，基于超聲圖像無(wú)法定位所有目標(biāo)；信號(hào)子空間維度過(guò)大，超聲成像結(jié)果中存在偽目標(biāo)，無(wú)法準(zhǔn)確評(píng)估被測(cè)對(duì)象內(nèi)部狀況。

圖7 奇異值分布曲線Fig.7 Distribution of singular values

最后，設(shè)定成像區(qū)域20 mm×20 mm，對(duì)于其中的任意像素點(diǎn)，定義方向向量，基于信號(hào)子空間得到TR-MUSIC圖像，如圖8所示。為了進(jìn)行比較，圖8給出了基于相同超聲陣列數(shù)據(jù)得到的常規(guī)平面B掃描圖像。

圖8 實(shí)驗(yàn)結(jié)果Fig.8 Experimental results

由圖8可知，基于平面B掃描圖像，無(wú)法區(qū)分不銹鋼試塊包含的6個(gè)相鄰貫通孔，基于TR-MUSIC圖像，則可以區(qū)分這6個(gè)相鄰貫通孔，并且對(duì)其進(jìn)行準(zhǔn)確定位。實(shí)驗(yàn)結(jié)果驗(yàn)證了TR-MUSIC方法的超分辨率特性。

4 結(jié)論

本文研究基于時(shí)間反轉(zhuǎn)聲學(xué)理論的超分辨率成像方法，用于解決常規(guī)超聲成像分辨率受到超聲波波長(zhǎng)限制的問(wèn)題。詳細(xì)闡述了超聲相控陣檢測(cè)原理；給出了TR-MUSIC法的數(shù)學(xué)定義，通過(guò)數(shù)值仿真分析了其點(diǎn)擴(kuò)散函數(shù)；搭建了實(shí)驗(yàn)系統(tǒng)，通過(guò)實(shí)驗(yàn)驗(yàn)證了其超分辨率特性。實(shí)驗(yàn)在不銹鋼試塊上完成，在試塊的內(nèi)部加工6個(gè)可視為點(diǎn)散射體、直徑為1 mm的貫通孔。實(shí)驗(yàn)結(jié)果表明，相比較常規(guī)的平面B掃描方法，TR-MUSIC方法可以克服瑞利準(zhǔn)則，實(shí)現(xiàn)超分辨率成像，提高超聲圖像質(zhì)量。