戴光智 韓國強(qiáng) 林偉毅 歐陽顯躍

(1.華南理工大學(xué) 計算機(jī)科學(xué)與工程學(xué)院，廣東 廣州 510006;2.華南理工大學(xué) 自動化科學(xué)與工程學(xué)院，廣東 廣州 510640;3.湖南大學(xué) 嵌入式與網(wǎng)絡(luò)計算湖南省重點實驗室，湖南 長沙 410082)

Donoho 和Candès 等［1-4］提出的采樣與壓縮同步進(jìn)行的壓縮感知(CS)理論與傳統(tǒng)的“先采樣后壓縮”不同，它采用“邊采樣邊壓縮”的方法.CS 理論指出，只要信號是稀疏的，那么信號的采樣頻率可以遠(yuǎn)遠(yuǎn)低于信號頻譜中最高頻率的2 倍而不影響信號的精確重建.將CS 理論應(yīng)用于工業(yè)超聲波成像系統(tǒng)中，可以有效地減少采樣點和采樣頻率.目前，國內(nèi)外將CS 應(yīng)用于超聲領(lǐng)域的研究［5-13］主要集中在醫(yī)用超聲方面，在工業(yè)超聲方面的研究很少.有限新息率(FRI)模型是傳統(tǒng)采樣理論與CS 相結(jié)合的信號采集方法［14］.具有FRI 性質(zhì)的信號可以由各個短脈沖信號的延時和幅度進(jìn)行完備的表示，因此可以有效地減少采樣數(shù)據(jù)和采樣頻率.工業(yè)超聲波成像系統(tǒng)的反射信號可以看成一系列不同延時和幅度的高斯脈沖信號的疊加，因此，工業(yè)超聲波成像系統(tǒng)的反射信號具有FRI 性質(zhì)，可以利用FRI 模型中的采樣方法進(jìn)行信號采樣和重建.但現(xiàn)有的FRI 信號采集方法對高頻信號進(jìn)行采樣時很不穩(wěn)定.通常，工業(yè)超聲波成像系統(tǒng)的反射信號頻率很高，導(dǎo)致現(xiàn)有的FRI 模型并不能很好地應(yīng)用于工業(yè)超聲波成像系統(tǒng)中.因此，文中研究能夠采集高頻信號的FRI 信號采集方法，對CS 理論在工業(yè)超聲波成像系統(tǒng)中的應(yīng)用具有重大的意義.

為有效地減少采樣數(shù)據(jù)和采樣頻率，文中提出了一種適用于工業(yè)超聲信號處理的新型FRI 采集方法——辛格函數(shù)之和(SoS)濾波器法，推導(dǎo)了這種采集方法應(yīng)具備的條件，給出了SoS 濾波器的構(gòu)造方式，并利用SoS 濾波器對一維超聲波成像原始數(shù)據(jù)進(jìn)行處理，最后通過仿真實驗驗證了該方法的可行性.

1 壓縮感知理論和FRI 模型

1.1 壓縮感知理論

CS 理論指出，當(dāng)信號滿足稀疏或可壓縮性時，可以采集少量觀測信號，然后通過求解一個優(yōu)化問題準(zhǔn)確重建原始信號.由信號理論可知，時域信號X(XRn×1)能夠用一組基Ψ(Ψ=(ψ1，ψ2，…，ψn))線性表示，即

式中，α 為n 維系數(shù)向量，α=(α1，α2，…，αn).如果式(1)中系數(shù)αi僅有k 個非0，則稱信號X 是稀疏的，稱Ψ 為信號X 的稀疏基或稀疏字典.如果對X的測量也是在時域上進(jìn)行的，則有

其中Θ(Θ=ΦΨ)為m×n 階矩陣.由于式(3)中的α 是k 稀疏的，且k ＜m ＜n，故α 可以通過式(4)的最優(yōu)化問題獲得，即

可以通過基于l1范數(shù)的最優(yōu)化方法(如貝葉斯方法)求解目標(biāo)信號α，但其前提是測量矩陣必須滿足RIP 性質(zhì).

1.2 FRI 模型

假設(shè)得到如下一個周期為 的無限周期信號:

式中，h(t)為形狀和周期已知的脈沖，tl和al(l=1，2，…，L)分別為未知的延時和幅度，L 為一個周期內(nèi)短脈沖的個數(shù).假設(shè)一個周期內(nèi)得到M 個連續(xù)的傅里葉系數(shù)，則x(t)的傅里葉變換為

選擇M個連續(xù)的整數(shù)，使得矩陣H 是可逆的，定義y=H-1x ，則有

式中，V為Vandermonde 矩陣，因此只要L ≤M，則V 是列滿秩的.

將方程組(8)進(jìn)行重寫，有

由式(7)-(9)可知，一旦得到傅里葉系數(shù)向量x，并保證2L≤M，就可以根據(jù)方程組(8)解出tl和al(l= 1，2，…，L).文中采用零化濾波重建算法［15］求解方程組(8)，當(dāng)然也可以使用其他的方法，如矩陣束算法［16］.

2 FRI 濾波器組及SoS 濾波器

2.1 FRI 濾波器組

傅里葉系數(shù)的采集是稀疏重建成像的關(guān)鍵.在時間域上對FRI 信號x(t)進(jìn)行采樣(采樣函數(shù)為s*())，可以得到M(2L≤M)個傅里葉系數(shù)，能實現(xiàn)此功能的方法的集合稱為FRI 濾波器組.假設(shè)待采樣信號為周期信號，其采樣框架如圖1 所示.

圖1 單通道采樣框架Fig.1 Framework of single channel sampling

由圖1 可知，對于FRI 信號x(t)及采樣周期T，實際采集到的向量為c= (c［1］，c［2］，…，c［N］)，其中

式中:n=1，2，…，N;S*(ω)為S(ω)的伴隨矩陣的元素，S(ω)為s(t)的連續(xù)傅里葉變換，ω=2 k/ .實際上，只要選擇

可見，對濾波后得到的采樣向量c 進(jìn)行離散傅里葉變換DFT{c}，再乘以S-1，就能得到傅里葉系數(shù)x.

在文獻(xiàn)［5］的FRI 濾波器模型中，s*(- t)=Bsin(-cBt)，B=M/ ，2L ≤M ≤N，并且s(t)是一個帶寬為B 的理想低通濾波器:

對于一個有限信號，采用FRI 濾波器的采樣模型定義如下:

假設(shè)h(t)滿足條件

經(jīng)過推算可以驗證，當(dāng)R ≤ 時，有以下公式成立:

其中，s3p(t)=s(t - )+ s(t)+ s(t + ).

2.2 SoS 濾波器

由于理想低通濾波器帶寬是有限的，而脈沖信號在時域上是無限的，因此文獻(xiàn)［5］的FRI 濾波器具有局限性，它不能擴(kuò)展到無限脈沖信號.通常，工業(yè)超聲回波信號的頻率很高，但超聲探頭及轉(zhuǎn)換硬件性能的限制使所獲得的超聲信號帶寬遠(yuǎn)遠(yuǎn)低于實際的信號帶寬，故存在較大的信息丟失.為此，基于FRI 濾波器組，在文獻(xiàn)［6］的基礎(chǔ)之上，文中提出了帶寬不受限制的SoS 濾波器:

式中，bk≠0，bk=，

將它轉(zhuǎn)化到時域上，有

同樣，有限信號采用SoS 濾波器采樣模型，定義為

選擇正確的參數(shù)bk(kκ)，能夠增強(qiáng)濾波器的抗噪性能，文中取長度為M 的漢明窗系數(shù):

可知，bk=，因此濾波后的采樣系數(shù)c［n］必定為實數(shù).

3 實驗及結(jié)果分析

3.1 Matlab 仿真結(jié)果分析

SoS 濾波法的實現(xiàn)步驟如下:

(1)輸入包含L 個高斯脈沖的信號x(t)，高斯脈沖為

其中σ=7 ×10-3，周期= 1 s，時間延時和幅度是隨機(jī)選擇的;

(2)根據(jù)實際的超聲檢測，取脈沖個數(shù)L=5，κ={- L，- L +1，…，L -1，L}，M== 11 ;

(3)采用漢明窗系數(shù)滿足式(23)(計算結(jié)果見圖2)的g(t)對輸入信號進(jìn)行濾波;

圖2 bk 的選取Fig.2 Selection of bk

根據(jù)式(7)-(9)對圖3 所示信號進(jìn)行濾波，SoS濾波器的實部如圖4 所示，采樣結(jié)果如圖5 所示，原始信號的估計延時和幅度如圖6 所示，所得到的重建信號如圖7 所示，由圖可見重建效果很好.

圖3 包含不同周期的原始信號Fig.3 Original signals containing different periods

圖4 SoS 濾波器的實部Fig.4 Real component of SoS filter

圖5 濾波輸出信號中采樣的21 個值Fig.5 Twenty-one sampling values of filter output signal

圖6 原始信號的估計延時和幅度Fig.6 Estimation time delay and amplitude of original signal

圖7 仿真重建信號與原始信號對比Fig.7 Comparison of reconstruction signal with original signal in simulation

3.2 實測結(jié)果分析

(1)利用工業(yè)超聲波成像系統(tǒng)采集一組一維原始信號(如圖8 所示)，其中單一探頭的中心頻率fc=3.4235 MHz，系統(tǒng)的采樣頻率fs= 40 MHz，成像深度Rmax= 0.16 m，在被測物體中，超聲波速c=1540 m/s，可以計算出信號的持續(xù)時間 '=2Rmax/c=2.08×10-4s;

(2)根據(jù)超聲波工業(yè)探頭的特性，將h(t)定為高斯脈沖是比較合適的，選擇高斯參數(shù)σ=3 ×10-7s，因σ ≤ '，故可以使用濾波器(21)對采樣數(shù)據(jù)進(jìn)行濾波，結(jié)果見圖9;

(3)取L=5，為了增強(qiáng)抗噪性能，設(shè)M=2 ×5 ×2+1=21，N=M=21，κ= {- 10，- 9，…，0，…，9，10}，采樣周期T= '/N=2.08 ×10-4/21 s=9.91 μs;

(4)利用式(20)-(22)構(gòu)造濾波器，計算采樣向量c;

(5)利用式(14)計算傅里葉系數(shù)向量x;

(6)利用y=H-1x 求出向量y，從而得到方程組(9);

(7)利用零化濾波重建算法求解方程(9)，得到估計的延時和幅度如圖10 所示，重建信號見圖11.

圖8 實測超聲信號Fig.8 Actually measured ultrasonic signals

圖9 濾波輸出信號的采樣值Fig.9 Sampling values of filter output signal

圖10 原始超聲信號的估計延時和幅度Fig.10 Estimation time delay and amplitude of original ultrasonic signal

圖11 實測原始超聲信號和重建信號比較Fig.11 Comparison of original ultrasonic signal with reconstruction signal in actual measurment

由重建效果可以看出，延時系數(shù)是精確的.在超聲波成像系統(tǒng)中，延遲系數(shù)代表著缺陷的位置，因此，只要延時系數(shù)精確，就能準(zhǔn)確地確定缺陷的位置.

4 結(jié)語

根據(jù)工業(yè)超聲波信號頻率較高的特點，基于FRI 濾波器組，文中提出了一種新的適用于工業(yè)超聲信號處理的FRI 信號采集方法——SoS 濾波器法，Matlab 仿真和實測結(jié)果表明該方法是有效的，從而驗證了CS 理論在工業(yè)超聲波檢測系統(tǒng)中應(yīng)用的可行性.文中就CS 理論在工業(yè)超聲波檢測系統(tǒng)中的應(yīng)用做了有益的探討和實測數(shù)據(jù)仿真，進(jìn)一步豐富了CS 理論.

［1］Donoho D L.Compressed sensing［J］.IEEE Transactions on Information Theory，2006，52(4):1289-1306.

［2］Candès E J，Romberg J，Tao T.Robust uncertainty principles:exact signal reconstruction from highly incomplete frequency information［J］.IEEE Transactions on Information Theory，2006，52(2):489-509.

［3］Candès E J.Compressive sampling［C］∥Proceedings of the International Congress of Mathematics.Madrid:European Mathematical Society，2006:1433-1452.

［4］Candès E J，Wakin M B.An introduction to compressive sampling［J］.IEEE Signal Processing Magazine，2008，25(2):21-30.

［5］Tur R，Eldar Y C，F(xiàn)riedman Z.Innovation rate sampling of pulse streams with application to ultrasound imaging［J］.IEEE Transactions on Signal Processing，2011，59(4):1827-1842.

［6］Wagner N，Eldar Y C，F(xiàn)euer A，et al.Compressed beamforming in ultrasound imaging［J］.IEEE Transactions on Signal Processing，2012，60(9):4643-4657.

［7］Chernyakova T，Eldar Y C.Fourier domain beamforming:the path to compressed ultrasound imaging［J/OL］.Arxiv，2013:11(2013-07-24).http:∥arxiv.org/abs/1307.6345.

［8］Wagner N，Eldar Y C，F(xiàn)euer A，et al.Compressed beamforming with applications to ultrasound imaging ［C］∥Proceedings of 2012 IEEE International Conference on Acoustics，Speech and Signal Processing.Kyoto:IEEE，2012:3641-3644.

［9］Friboulet D，Liebgott H，Prost R.Compressive sensing for raw RF signals reconstruction in ultrasound ［C］∥Proceedings of 2010 IEEE Ultrasonics Symposium.San Diego:IEEE，2010:367-370.

［10］Quinsac C，Basarab A，Kouamé D.Frequency domain compressive sampling for ultrasound imaging［J/OL］.Advances in Acoustics and Vibration，2012，2012:16(2012-04-05).http:∥dx.doi.org/10.1155/2012/231317.

［11］Quinsac C，Basarab A，Kouame D，et al.3D compressed sensing ultrasound imaging［C］∥Proceedings of 2010 IEEE Ultrasonics Symposium.San Diego:IEEE，2010:363-366.

［12］Quinsac C，Basarab A，Girault J，et al.Compressed sensing of ultrasound images:sampling of spatial and frequency domains ［C］∥Proceedings of 2010 IEEE Workshop on Signal Processing Systems.San Diego:IEEE，2010:231-236.

［13］Shi Guangming，Chen Chongyu，Lin Jie，et al.Narrowband ultrasonic detection with high range resolution:separating echoes via compressed sensing and singular value decomposition［J］.IEEE Transactions on Ultrasonics，F(xiàn)erroelectrics and Frequency Control，2012，59(10):2237-2253.

［14］Vetterli M，Marziliano P，Blu T.Sampling signals with finite rate of innovation［J］.IEEE Transactions on Signal Processing，2002，50(6):1417-1428.

［15］Stoica P，Moses R.Introduction to spectral analysis［M］.Englewood Cliffs:Prentice-Hall，2000.

［16］Sarkar T K，Pereira O.Using the matrix pencil method to estimate the parameters of a sum of complex exponentials［J］.IEEE Antennas and Propagation Magazine，1995，37(1):48-55.