復(fù)域壓縮感知近場(chǎng)聲圖測(cè)量方法

陳新華， 鄭恩明， 李嶷， 楊鶴,2， 周權(quán)斌,3

(1.中國(guó)科學(xué)院 聲學(xué)研究所, 北京 100190; 2.哈爾濱工程大學(xué) 水聲工程學(xué)院, 黑龍江 哈爾濱 150001;3.中國(guó)科學(xué)院大學(xué), 北京 100049)

0 引言

聲圖測(cè)量是一種通過(guò)聚焦波束形成方法實(shí)現(xiàn)近場(chǎng)精確定位的技術(shù)，被廣泛應(yīng)用于空氣聲學(xué)、水聲學(xué)和醫(yī)學(xué)超聲領(lǐng)域的聲學(xué)成像技術(shù)[1-3]。為得到性能較優(yōu)的聲圖測(cè)量方法，一些學(xué)者從常規(guī)波束形成線性模型到子空間分類非線性模型進(jìn)行了深入研究，并提出了一些相應(yīng)改善方法[4-7]。雖然子空間分類方法所得聲圖具有較高的空間分辨率，但需要事先預(yù)估目標(biāo)信源個(gè)數(shù)，信噪比較低時(shí)，存在目標(biāo)信源個(gè)數(shù)估計(jì)不準(zhǔn)問(wèn)題，致使子空間分類方法存在漏報(bào)或虛報(bào)目標(biāo)，限制了子空間分類方法在實(shí)際工程中的應(yīng)用[8-11]。

壓縮感知作為一種改變“奈奎斯特”采樣理論的新興理論，已被廣泛應(yīng)用到水下目標(biāo)方位估計(jì)和定位中。例如：陳力恒等[12]將稀疏陣列接收數(shù)據(jù)模型轉(zhuǎn)化為更高自由度下的單快拍接收數(shù)據(jù)模型，并結(jié)合壓縮感知模型實(shí)現(xiàn)了稀疏陣列波束形成；李賀等[13]基于壓縮感知和約束隨機(jī)線陣實(shí)現(xiàn)了目標(biāo)方位估計(jì)；沈小正等[14]基于聲矢量傳感器陣空間稀疏模型和壓縮感知快速實(shí)現(xiàn)了小樣本條件下的方位估計(jì)；郭雙樂(lè)等[15]在水下目標(biāo)定位的稀疏數(shù)學(xué)模型基礎(chǔ)上，結(jié)合壓縮感知理論，實(shí)現(xiàn)了高效的水下目標(biāo)定位；時(shí)潔等[16]以壓縮感知空間目標(biāo)空域稀疏性作基礎(chǔ)，通過(guò)構(gòu)造相應(yīng)的感知矩陣和觀測(cè)序列實(shí)現(xiàn)了聲圖測(cè)量。但受復(fù)雜海洋環(huán)境影響，采用壓縮感知進(jìn)行水下目標(biāo)方位估計(jì)或定位中，會(huì)存在性能下降問(wèn)題，對(duì)此：康春玉等[17-18]基于壓縮感知采用單快拍數(shù)據(jù)實(shí)現(xiàn)了目標(biāo)方位估計(jì)和信號(hào)恢復(fù)，并采用重構(gòu)頻域陣列信號(hào)提升壓縮感知在目標(biāo)方位估計(jì)中的性能；周明陽(yáng)等[19]基于改進(jìn)的高斯先驗(yàn)?zāi)Ｐ吞嵘素惾~斯壓縮感知在目標(biāo)方位估計(jì)中的性能。上述基于壓縮感知進(jìn)行的目標(biāo)方位估計(jì)或定位均是在頻域?qū)崿F(xiàn)的，構(gòu)造觀測(cè)序列和感知矩陣方法相似，低信噪比條件情況下，均存在一定的性能退化問(wèn)題。因此，必須建立一種其他變換域壓縮感知近場(chǎng)聲圖測(cè)量模型方法，通過(guò)改變觀測(cè)序列和感知矩陣構(gòu)建方式，提升壓縮感知近場(chǎng)聲圖測(cè)量對(duì)信噪比寬容性。

1 頻域壓縮感知方法

1.1 信號(hào)模型

線列陣近場(chǎng)信號(hào)模型如圖1所示，線列陣布放在空間x軸上，掃描平面上有K個(gè)聲源，K∈[1,M]，M為線列陣傳感器個(gè)數(shù)，圖1中展示了第k個(gè)聲源與線列陣之間的空間幾何關(guān)系，k∈[1,K]。

圖1 均勻線列陣布放示意圖Fig.1 Schematic diagram of liner array layout

圖1中：(0,0)為線列陣中心位置；1#為線列陣第1個(gè)傳感器編號(hào)；M#為線列陣第M個(gè)傳感器編號(hào)；R為聲源相對(duì)線列陣中心位置距離標(biāo)識(shí)；θ為聲源相對(duì)線列陣整橫方向夾角標(biāo)識(shí)；N為空間掃描點(diǎn)數(shù)，N?K.

如圖1所示，令K個(gè)近場(chǎng)聲源相對(duì)該線列陣中心位置(參考點(diǎn))所在空間位置為(RK,ΘK)=[(R1,θ1),(R2,θ2),…，(RK,θK)]，則線列陣拾取數(shù)據(jù)X(f)(后續(xù)簡(jiǎn)稱陣列數(shù)據(jù))可表示為

X(f)=A(RK,ΘK)S(f)+V(f)，

(1)

式中：f為頻率變量；S(f)為K個(gè)聲源輻射信號(hào)，S(f)=[S1(f),S2(f),…，SK(f)]T；V(f)為線列陣拾取背景噪聲數(shù)據(jù)，V(f)=[V1(f),V2(f),…，VM(f)]T；A(RK,ΘK)為線列陣陣列流形矩陣，A(RK,ΘK)=[a(R1,θ1),a(R2,θ2),…，a(RK,θK)]；(·)T為矩陣轉(zhuǎn)置運(yùn)算符。

在近場(chǎng)聲傳播中，a(Rk,θk)的具體形式為：

a(Rk,θk)=[e-j2πfτk,1,e-j2πfτk,2,…，e-j2πfτk,M]T，

(2)

1.2 聲圖測(cè)量原理

a(fj,(n,n))=
[e-j(2πfj/c)τn,1,e-j(2πfj/c)τn,2,…,e-j(2πfj/c)τn,M]T.

(3)

陣列數(shù)據(jù)的頻域窄帶模型可表示為

(fj)=A(fj,N,N)S(fj)+V(fj).

(4)

很顯然，A(RK,ΘK)∈A(fj,N,N)；在(RK,ΘK)位置上，S(fj)值較大，而其他位置上，S(fj)較小，即S(fj)是信號(hào)空間域的一種稀疏表示。此時(shí)，如果將(fj)作為觀測(cè)序列，A(fj,N,N)作為感知矩陣W(fj,N,N)，S(fj)為待求解稀疏系數(shù)分量，則可按壓縮感知模型通過(guò)求解如下凸優(yōu)化問(wèn)題實(shí)現(xiàn)空間信號(hào)稀疏系數(shù)S(fj)求解。

(5)

通過(guò)S(fj)計(jì)算出每個(gè)掃描點(diǎn)上對(duì)應(yīng)的信號(hào)能量，如(6)式所示，即第j子帶的聲圖測(cè)量，且在存在目標(biāo)位置(RK,ΘK)上有最大值，而無(wú)目標(biāo)位置則為一個(gè)足夠小的值。

P(fj,N,N)=|S(fj)|2.

(6)

所有子帶重復(fù)上述過(guò)程，并對(duì)P(fj,N,N)累加求和，可得寬帶聲圖測(cè)量結(jié)果P(N,N)，搜索P(N,N)峰值位置即可完成對(duì)(RK,ΘK)估計(jì)。

(7)

(5)式是基于陣列數(shù)據(jù)頻域形式實(shí)現(xiàn)近場(chǎng)聲圖測(cè)量[17-18]，本文將該方法稱為頻域壓縮感知(FCS)方法。

由凸優(yōu)化問(wèn)題來(lái)求解過(guò)程可知，F(xiàn)CS方法凸優(yōu)化求解過(guò)程核心步驟是根據(jù)觀測(cè)序列(fj)和感知矩陣W(fj,N,N)按S(fj)=(WT(fj,N,N)·W(fj,N,N))-1WT(fj,N,N)(fj)求解空間信號(hào)稀疏系數(shù)，該結(jié)果只是對(duì)觀測(cè)序列采用線性累加方式近似求解得到空間信號(hào)稀疏系數(shù)，未對(duì)觀測(cè)序列做相應(yīng)變換處理以提高觀測(cè)序列數(shù)據(jù)所含信噪比，提升凸優(yōu)化求解穩(wěn)定性。

2 復(fù)域壓縮感知方法

2.1 問(wèn)題解析

由FCS近場(chǎng)聲圖測(cè)量模型可知，利用頻域數(shù)據(jù)構(gòu)建觀測(cè)序列時(shí)，只是利用各傳感器數(shù)據(jù)自身信息，并未利用各傳感器數(shù)據(jù)中信號(hào)與噪聲相關(guān)性差異特性實(shí)現(xiàn)處理數(shù)據(jù)所含信噪比的增強(qiáng)。對(duì)此，本文提出一種復(fù)域處理方法，在復(fù)域通過(guò)對(duì)各傳感器數(shù)據(jù)進(jìn)行相關(guān)、累加處理，提升壓縮感知重構(gòu)聲圖的穩(wěn)定性，稱之為復(fù)域壓縮感知(CCS)方法。

2.2 原理分析

可得經(jīng)時(shí)延補(bǔ)償后數(shù)據(jù)為

(8)

(9)

式中：I=[1,1,…,1].

為了實(shí)現(xiàn)CCS近場(chǎng)聲圖測(cè)量模型構(gòu)建，對(duì)(9)式進(jìn)行變換處理，將I作為觀測(cè)序列，P(N,N)作為待求解稀疏系數(shù)分量，則感知矩陣W(N,N)=[W(1,1)W(2,2) …W(N,N)]，W(n,n)可設(shè)計(jì)為

W(n,n)=

(10)

采用(10)式處理后，可使陣列數(shù)據(jù)相關(guān)所得結(jié)果變?yōu)榕c掃描位置信息有關(guān)的列數(shù)據(jù)，賦值到感知矩陣相應(yīng)列中，完成感知矩陣構(gòu)建。

最后，可通過(guò)求解如下凸優(yōu)化問(wèn)題來(lái)求解空間信號(hào)稀疏系數(shù)S(t)：

(11)

根據(jù)空間信號(hào)稀疏系數(shù)S(t)，此時(shí)可得近場(chǎng)聲圖測(cè)量結(jié)果P(N,N)，搜索P(N,N)峰值位置即可完成對(duì)(RK,ΘK)估計(jì)。

P(RN,N)=|S(t)|2.

(12)

2.3 能力分析

(13)

(fj)=[X1(fj),X2(fj),…，XM(fj)]T=
[S(fj)+V1(fj),S(fj)+V2(fj),…，S(fj)+VM(fj)]T，

(14)

式中：S(fj)為代表空間目標(biāo)信號(hào)；Vm(fj)為第m個(gè)傳感器噪聲數(shù)據(jù)。觀測(cè)序列每個(gè)位置數(shù)據(jù)所含信噪比為

(15)

同樣，由(13)式和(10)式可知，事先將陣列信號(hào)處理引入到感知矩陣構(gòu)造中，充分利用了信號(hào)噪聲之間相關(guān)性差異特性，通過(guò)對(duì)陣列數(shù)據(jù)進(jìn)行相關(guān)、累積處理，可使感知矩陣各位置數(shù)據(jù)具有一定陣增益，因而CCS方法中感知矩陣每個(gè)位置元數(shù)所含信噪比為

(16)

對(duì)比(15)式和(16)式可知，通過(guò)將FCS方法中觀測(cè)序列數(shù)據(jù)移植到感知矩陣構(gòu)建中，CCS方法提升了感知矩陣每個(gè)位置數(shù)據(jù)所含信噪比，改善了FCS方法對(duì)聲圖重構(gòu)效果，對(duì)最低信噪比的適應(yīng)性得到近10lgMdB提升。

由此可知，CCS方法與FCS方法本質(zhì)區(qū)別是壓縮感知數(shù)學(xué)模型構(gòu)建中對(duì)陣列數(shù)據(jù)處理方面的不同。FCS方法是直接將陣列數(shù)據(jù)作為觀測(cè)序列；然后再采用完備陣列流形或陣列某一傳感器數(shù)據(jù)形成感知矩陣；在壓縮感知模型構(gòu)建中未對(duì)陣列數(shù)據(jù)之間進(jìn)行相應(yīng)處理。而CCS方法首先對(duì)FCS中觀測(cè)序列數(shù)據(jù)進(jìn)行相關(guān)、累加等處理，即將FCS方法中觀測(cè)序列數(shù)據(jù)事先進(jìn)行非線性處理，將其移植到感知矩陣構(gòu)建中，即利用觀測(cè)序列各位置數(shù)據(jù)所含信號(hào)、噪聲相關(guān)性差異特性提升感知矩陣各位置數(shù)據(jù)所含信噪比；之后采用凸優(yōu)化求解方法對(duì)空間信號(hào)稀疏系數(shù)求解過(guò)程進(jìn)一步對(duì)陣列數(shù)據(jù)作了二次非線性處理，提升了求解空間信號(hào)稀疏系數(shù)的穩(wěn)定性。

2.4 CCS方法實(shí)現(xiàn)

由2.3節(jié)理論分析可知，CCS方法實(shí)現(xiàn)流程如圖2所示，具體過(guò)程可分為如下步驟實(shí)現(xiàn)：

圖2 復(fù)域壓縮感知方法實(shí)現(xiàn)流程圖Fig.2 Flow chart of compressed sensing in complex domain

輸入：陣列數(shù)據(jù)x(t)，處理頻帶[wd,wu]，wd為處理頻帶下限，wu為處理頻帶上限，L為分幀數(shù)。

輸出：聲圖P(N,N)，目標(biāo)位置估計(jì)值(K,K)。

5) 更新處理幀數(shù)據(jù)，l=l+1，重復(fù)執(zhí)行步驟2至步驟4，直到l=L；

6) 對(duì)L個(gè)聲圖測(cè)量值Pl(N,N)進(jìn)行累加處理，得到本次處理最終聲圖測(cè)量值：

(17)

7) 對(duì)P(N,N)進(jìn)行峰值篩選，可得到目標(biāo)位置估計(jì)值(K,K)。

3 數(shù)據(jù)處理分析

由于本文只討論CCS近場(chǎng)聲圖測(cè)量模型搭建思路，與凸優(yōu)化算法無(wú)關(guān)，所以后續(xù)數(shù)據(jù)處理分析中，無(wú)論是FCS方法還是CCS方法，均采用正交匹配跟蹤(OMP)算法實(shí)現(xiàn)對(duì)聲圖重構(gòu)。

3.1 數(shù)值仿真分析

3.1.1 單目標(biāo)情況

在該仿真中首先按表1所示參數(shù)進(jìn)行設(shè)置，然后采用不同方法進(jìn)行處理分析。

表1 水平陣仿真參數(shù)Tab.1 Numerical simulation parameters of horizontal array

掃描平面為水平距離[20 m,100 m]、方位角度[-90°,90°]形成的區(qū)域，將該區(qū)域按掃描網(wǎng)格劃分，垂直方向網(wǎng)格間距為2 m，水平角度網(wǎng)格間距為1°.在最小方差無(wú)畸變響應(yīng)(MVDR)方法、多重信號(hào)分類(MUSIC)方法、FCS方法和CCS方法實(shí)現(xiàn)中首先對(duì)單次處理數(shù)據(jù)分64幀處理，相鄰幀按半幀采樣數(shù)據(jù)重疊處理。圖3～圖5為目標(biāo)信噪比(SNR)為-25～0 dB情況下，由MVDR方法、FCS方法和CCS方法通過(guò)200次獨(dú)立統(tǒng)計(jì)所得目標(biāo)檢測(cè)概率和估計(jì)均方根誤差(RMSE).

圖3 4種方法的目標(biāo)檢測(cè)概率Fig.3 Target detection probabilities of four methods

圖4 4種方法的距離估計(jì)RMSEFig.4 RMSEs of space estimation of four methods

圖5 4種方法的方位估計(jì)RMSEFig.5 RMSEs of bearing estimation of four methods

由圖3～圖5的仿真結(jié)果可知，相比FCS方法，CCS方法對(duì)最低信噪比要求上降低了8 dB，提升了壓縮感知重構(gòu)聲圖穩(wěn)定性，低信噪比情況下實(shí)現(xiàn)了對(duì)目標(biāo)位置有效估計(jì)。

3.1.2 多目標(biāo)情況

陣列參數(shù)、數(shù)據(jù)處理參數(shù)以及4種方法實(shí)現(xiàn)與單目標(biāo)情況一致，目標(biāo)信號(hào)為兩個(gè)等強(qiáng)度信號(hào)，頻帶為200～500 Hz，目標(biāo)分別位于相對(duì)線列陣陣中心(40 m,-2°)和(40 m,2°)位置處。 圖6～圖8為不同信噪比情況下，由4種方法所得聲圖。

圖6 4種方法所得40 m處不同方位聲圖(SNR=5 dB)Fig.6 Output acoustic images at 40 m obtained by four methods (SNR=5 dB)

圖7 4種方法所得40 m處不同方位聲圖(SNR=-5 dB)Fig.7 Output acoustic images at 40 m obtained by four methods (SNR=-5 dB)

圖8 4種方法所得40 m處不同方位聲圖(SNR=-15 dB)Fig.8 Output acoustic images at 40 m obtained by four methods (SNR=-15 dB)

由圖6～圖8仿真結(jié)果可得：相比MVDR方法和MUSIC方法，F(xiàn)CS方法和CCS方法具有更窄的主瓣寬度，可對(duì)相鄰目標(biāo)實(shí)現(xiàn)高分辨檢測(cè)和分辨；但隨著信噪比的降低，MVDR方法和MUSIC方法分辨能力下降較為厲害，在SNR為-5 dB時(shí)，MVDR方法和MUSIC方法已不能對(duì)兩目標(biāo)的實(shí)現(xiàn)有效分辨，而FCS方法和CCS方法同樣保持著高信噪比情況下的目標(biāo)分辨能力；但在SNR為-15 dB時(shí)，F(xiàn)CS方法所得聲圖中兩目標(biāo)之間的背景級(jí)受噪聲污染比較嚴(yán)重，已無(wú)法對(duì)兩目標(biāo)實(shí)現(xiàn)有效檢測(cè)。

3.1.3 運(yùn)算量分析

由于CCS方法可事先求取所需復(fù)解析小波，所以影響CCS方法運(yùn)算量主要因素為復(fù)解析小波函數(shù)與陣列數(shù)據(jù)卷積運(yùn)算、協(xié)方差矩陣求取以及OMP求解，F(xiàn)CS方法的運(yùn)算量主要在于傅里葉變換、協(xié)方差矩陣及其求逆運(yùn)算以及OMP求解，由于FCS方法和CCS方法都有OMP求解過(guò)程，后續(xù)兩種方法計(jì)算復(fù)雜度比較過(guò)程不再將OMP求解過(guò)程考慮進(jìn)去。

3.2 海試數(shù)據(jù)分析

本次處理數(shù)據(jù)由進(jìn)行近場(chǎng)目標(biāo)定位試驗(yàn)所得，相關(guān)參數(shù)如表2所示。

表2 水平陣數(shù)據(jù)處理參數(shù)Tab.2 Parameters of processing data

水平陣布放示意圖如圖9所示，在小于100 m、-90°～90°范圍內(nèi)，存在兩個(gè)近場(chǎng)目標(biāo)，兩個(gè)目標(biāo)輻射信號(hào)頻帶主要為200～500 Hz.掃描平面為水平距離[20 m,100 m]、方位角[-90°,90°]形成的區(qū)域，將該區(qū)域按掃描網(wǎng)格劃分，垂直方向網(wǎng)格間距為2 m，水平角度網(wǎng)格間距為1°.在該處理數(shù)據(jù)時(shí)段，目標(biāo)1相對(duì)水平陣中心位置(52 m，-16°)附近處，目標(biāo)2相對(duì)水平陣中心位置(68 m，-9°)附近處。

圖9 陣列布放示意圖Fig.9 Schematic diagram of array layout

4種方法處理過(guò)程與3.1.1節(jié)數(shù)值仿真一致。圖10～圖13為4種方法所得某一時(shí)聲圖。

圖10 MVDR方法所得聲圖Fig.10 Output acoustic image obtained by MVDR method

圖11 MUSIC方法所得聲圖Fig.11 Output acoustic image obtained by MUSIC method

圖12 FCS方法所得聲圖Fig.12 Output acoustic image obtained by FCS method

圖13 CCS方法所得聲圖Fig.13 Output acoustic image obtained by CCS method

由圖10～圖13所示結(jié)果可知：MVDR方法和MUSIC方法背景級(jí)較高，顯示效果較差，影響最后目標(biāo)位置估計(jì)；FCS方法受信噪比影響比較嚴(yán)重，主要原因在于采用FCS模型構(gòu)建觀測(cè)序列時(shí)，其穩(wěn)定性差于經(jīng)過(guò)多個(gè)時(shí)域采樣點(diǎn)數(shù)據(jù)相關(guān)處理后數(shù)據(jù)，在低信噪比下，再采用凸優(yōu)化進(jìn)行空間信號(hào)稀疏系數(shù)求解所得結(jié)果穩(wěn)定性不如高信噪比下穩(wěn)定，容易致使其他位置處出現(xiàn)虛假目標(biāo)；而CCS方法所得聲圖能清晰顯示兩目標(biāo)位置，且目標(biāo)方位明晰可辨，背景級(jí)遠(yuǎn)低于另外3種方法，該現(xiàn)象進(jìn)一步說(shuō)明了通過(guò)對(duì)復(fù)域數(shù)據(jù)進(jìn)行相關(guān)、累加處理，提升了基于壓縮感知模型的近場(chǎng)聲圖測(cè)量性能。

實(shí)測(cè)數(shù)據(jù)處理結(jié)果再次驗(yàn)證：CCS方法通過(guò)對(duì)陣列數(shù)據(jù)進(jìn)行復(fù)解析變換和數(shù)據(jù)相關(guān)處理，在低信噪比下對(duì)目標(biāo)位置實(shí)現(xiàn)了有效估計(jì)。

4 結(jié)論

針對(duì)FCS近場(chǎng)下降的問(wèn)題，本文從觀測(cè)序列和感知矩陣構(gòu)建方面入手，通過(guò)對(duì)壓縮感知觀測(cè)序列進(jìn)行相關(guān)、累加預(yù)處理，將觀測(cè)序列數(shù)據(jù)移植到感知矩陣構(gòu)建中，并根據(jù)目標(biāo)在空間域的稀疏性，構(gòu)建了復(fù)域壓縮感知凸優(yōu)化模型，實(shí)現(xiàn)了CCS近場(chǎng)聲圖測(cè)量。

本文方法的關(guān)鍵是采用相關(guān)、累加處理對(duì)陣列數(shù)據(jù)進(jìn)行預(yù)處理，構(gòu)建感知矩陣和凸優(yōu)化問(wèn)題，得到了具有更高信噪比的感知矩陣，保證了后續(xù)聲圖測(cè)量方法性能提高。數(shù)值仿真和實(shí)測(cè)數(shù)據(jù)處理結(jié)果均證明了CCS方法在聲圖測(cè)量方面的有效性，待別是在弱目標(biāo)檢測(cè)方面明顯優(yōu)于FCS方法。

另外，可以預(yù)見(jiàn)，如果進(jìn)一步采用改善FCS的方法來(lái)改善CCS方法，則CCS方法的聲圖測(cè)量性能也將會(huì)進(jìn)一步得到改善，這也是后續(xù)需要進(jìn)一步研究的問(wèn)題。