基于壓縮傳感的水下目標方位估計方法*

（杭州應用聲學研究所聲納技術重點實驗室 杭州 310023）

1 引言

聲信號的波達方向（direction-of-arrival，DOA）估計是聲納在目標探測過程中非常重要的一個環節，它可以估計水下目標的空間方位。常規波束形成（conventional beamforming，CBF）［1］是 DOA 估計的傳統方法，其分辨力受到瑞利限的限制。為了突破這一限制，發展起了一些高分辨的空間譜估計方法，如基于子空間分解的多重信號分類（multiple signal classification，MUSIC）方法［2］。除此之外，為了自適應地抑制干擾，基于最大信噪比準則的最小方差無失真響應（Minimum Variance Distortionless Response，MVDR）方法［3］可以自適應地估計DOA。無論是MUSIC方法還是MVDR方法，都需要計算信號的協方差矩陣，這導致了低快拍數、相干源、低信噪比條件下的性能下降。

2006年前后，壓縮傳感（compressive sensing or compressed sensing）又稱為壓縮采樣（compressive sampling）發展起來。E.Candès和 T.Tao等［4～6］證明了，利用信號的稀疏性，能夠用比Nyquist-Shannon采樣定理所要求的少得多的采樣實現稀疏信號的完全重構，因此取名壓縮采樣。同時期，D.Donoho等［7］提出了一種新的信息獲取指導理論，即壓縮傳感。CS用于求解欠定線性逆問題，致力于發展信號的稀疏性或可壓縮性，用比Nyquist采樣率更低的采樣率做采樣而實現信號的重構。稀疏性或可壓縮性是信號的自然性質，聲信號的稀疏性或可壓縮性表現在時域、頻率域、空間域、波數域、空間波束域等域中。傳統DOA估計方法的目標是使重構信號的能量最小，因而得到的是低分辨力、非稀疏的解。將CS應用到聲信號的DOA估計中正是利用了空間波束域的稀疏性，因此得到的是高分辨的稀疏解。李璇等［8］分別利用對角加載最小二乘、λ1正則化和正交匹配追蹤方法求解單快拍情況下的DOA估計問題，得到了高分辨的方位估計結果。房云飛等［9］提出了一種基于波束域的多重測量向量欠定系統正則化聚焦求解DOA的算法，結果表明有較好的估計精度和角度分辨力。A.Xenaki和P.Gerstoft等［10～11］把單快拍和多快拍情況下的 DOA估計分別建模為單一測量向量和多重測量向量的問題，結果表明CS方法的性能優于傳統的波束形成和高分辨力方法。

本文將針對陣元間距大于半波長的線陣，即空間采樣率低于Nyquist采樣率進行研究。首先把DOA估計建模為一個標準的CS問題，并分析傳感矩陣的相干性（coherence），然后利用凸優化方法［12～13］和基于閾值的方法［14］對DOA估計這一欠定線性逆問題進行求解。數值仿真結果和實驗結果均表明了基于CS的波束形成方法的可行性，其分辨力要高于傳統的DOA估計方法，可在單快拍、相干源、低信噪比條件下估計目標的空間方位。

2 壓縮傳感理論

經典的Shannon采樣定理以帶寬為先驗信息，是一個充分非必要條件，對于頻率低于最高頻率的分量，一個周期上采的樣本多于兩個，因此存在冗余。空間采樣亦是如此，陣元間距為最高頻率的半波長時，那么小于最高頻率的波長，一個周期上有兩個以上的樣本，存在冗余。對于傳統的信號處理，在時域按Nyquist-Shannon采樣較容易實現，但由于實際中信號帶寬較大或最高頻率較高使得采樣頻率也較高，因而采集到的數據量較大；在空間域上按Nyquist-Shannon采樣需要在空間內均勻地布放大量傳感器，成本太高。壓縮傳感從根本上顛覆了Nyquist-Shannon采樣定理，在信號滿足稀疏性或可壓縮性的條件下，直接對信號進行壓縮采樣，通過復雜度較高的重構算法來重構原始信號。

壓縮傳感有三個基本要素，分別是稀疏性和可壓縮性、壓縮采樣和重構算法。

1）稀疏性和可壓縮性

一個信號如果在某一域中的分量大多數為零，則稱這個信號是稀疏的。數學上，用λ0范數表示非零元素的個數，如果向量x∈CN最多有s個元素是非零的，就稱 x是s-稀疏的，表示為‖x‖0≤s。稀疏性是一個強約束，實際中的許多信號是可壓縮的。可壓縮性指信號可以很好地由稀疏信號近似表示，這種表示常常是在經過適當的基變換之后，如Fourier變換、小波變換等。假設信號z∈CN可由稀疏向量x∈CN表示為

其中Φ∈CN×N稱為變換基矩陣。

2）壓縮采樣

利用測量矩陣Ψ∈CM×N對信號z進行采樣，得到測量數據y∈CM為

其中 A∈CM×N稱為傳感矩陣，是變換基矩陣和測量矩陣的結合。測量矩陣Ψ可以是隨機的，也可以是確定的，其行數M一定小于列數N，故該采樣過程把一個高維信號變換為低維數據，因而稱為壓縮采樣。

3）重構算法

經過壓縮采樣后得到了一個欠定線性方程

壓縮傳感的目的是利用稀疏性，從少量的測量數據y中重構x或z，本質是求解該欠定線性方程的稀疏解。重構算法能夠在適當的假設下，求得上述欠定線性方程的稀疏解。要想保證重構算法能準確地和有效地重構信號，傳感矩陣A需要滿足一定的性質，如零空間特性（null space property，NSP）、約束等距特性（restricted isometry property，RIP）和相干性等，其中前兩者在計算上都是NP-hard的，因此論文中考慮相干性。重構算法可以粗略地分成三類：優化方法，貪婪方法和基于閾值的方法。

3 基于壓縮傳感的DOA估計方法

3.1 線陣DOA估計的CS模型

假設聲源位于線陣的遠場，即線陣接收到的聲波可以考慮為平面波，考慮聲源輻射的聲信號是窄帶的。線陣示意圖如圖1所示。

根據圖1示意圖，目標聲源的DOA由相對于陣列方向的 θ∈[-90°，90°]決定。第i個聲源對每個傳感器的傳播時延可以由駕駛向量表示：

圖1 線陣示意圖

令向量x∈CN由感興趣方向θ格點上的源信號幅度組成，y表示聲場在M個傳感器處的測量向量。由于我們感興趣的是角度格點上的高分辨力，因此有M＜N。傳感矩陣則由所有可能方向的駕駛向量組成，即

它是測量矩陣Ψ和基變換矩陣Φ的乘積，其中Ψ表示聲場在傳感器位置處的空間采樣，是一個確定的矩陣，Φ表示一個IDFT變換基。考慮存在噪聲n∈CM的情況，測量向量就可以表示為

這是一個標準的欠定線性系統，我們的目標就是利用x的稀疏性，從y中重構x。

3.2 優化方法和基于閾值的方法

CS中求解式（1）的問題被稱為（P0）問題，它屬于非凸優化問題，且一般是NP-hard的，因此通常把（P0）問題凸松弛為（P1）問題，表示為

（P1）問題是一個凸優化問題，常被稱作λ1最小化或基追蹤（basis pursuit，BP）。對于存在噪聲的情況，可以利用Lagrange乘子法求解式（4），被稱為問題，表示為

其中正則化參數λ＞0控制著稀疏解和測量值擬合之間的相對重要性。求解問題的方法被稱作基追蹤去噪（basis pursuit denoising，BPD）。在Matlab上可以利用CVX工具箱求解上述凸優化問題，它利用內點求解器得到優化問題的全局解。

除了優化方法，還有一種基于閾值的方法。由于λ0最小化是非線性操作，因此基于閾值的方法需要定義一個非線性算子Hs，它稱為稀疏性為s的硬閾值算子。對于一個向量x，Hs（x）保持x的s個絕對值最大的元素不變，其他元素置為零。這里介紹一種已知解是s-稀疏情況下的迭代算法，稱為迭代硬閾值（iterative hard thresholding，IHT）。該算法是通過迭代地求解定點方程xn+1=(I-AHA)xn+AHy得到的，其中I為單位陣，這里的上標“H”表示矩陣的共軛轉置。由于已知所求的解是s-稀疏的，所以只需在每次迭代時保留s個絕對值最大的元素，最終可以得到一個s-稀疏的向量。

3.3 傳感矩陣的相干性分析

傳感矩陣A的特性會影響重構算法的性能。如果不考慮噪聲，NSP是穩健的信號重構的充要條件。如果考慮存在噪聲，Candès和Tao引入了RIP［15］，它對于各種算法成功地從有噪測量中重構一個稀疏信號是充分非必要的。無論是RIP還是NSP，它們在計算上都是NP-hard的，因此這里考慮易于計算的矩陣A的相干性［16］。相干性度量了傳感矩陣任意兩列之間的相關性（correlation），首先給出一個矩陣A的相干性μ(A)的定義式為

對于式（3）所表示的傳感矩陣，有

在DOA估計中，為了實現高分辨力，通常需要精細的格點，因此欠定方程中的傳感矩陣A具有線性相關的列，其相關程度反映在G的非對角元素上。下面給出了陣元間距大于半波長的線陣的Gram矩陣示意圖。

從圖2中可以看出當d=4λ5時，不出現柵瓣的范圍大約是[-15°，15°]。但是由于CS重構算法是非線性的，且稀疏解具有唯一性，因而在一定程度上可以準確重構出稀疏信號。

圖2 十六元均勻線陣的Gram矩陣示意圖，d/λ=4/5

4 數值仿真和實驗結果

4.1 數值仿真結果

考慮一個16元均勻線陣，陣元間距d=4λ/5。首先考慮遠場存在一個目標，入射方向為20°，信噪比為10dB，分別采用CBF、MVDR、MUSIC、BPD和IHT算法進行DOA估計。仿真結果如圖3所示，其中CBF、MVDR和MUSIC算法的快拍數為50，BPD和IHT算法采用單快拍。從圖3中可以看出，CBF和MVDR方法會在-65.3°出現柵瓣，雖然IHT算法仍然存在柵瓣，但BPD算法在一定程度上可以準確估計目標方位。因此基于CS的DOA估計可以在單快拍的條件下準確重構出目標方位，且分辨力高于CBF和MVDR。

圖3 單目標DOA估計

然后考慮存在兩個相干源，入射方位分別為0°和5°，信噪比為-5dB，其余條件與之前一致，結果如圖4所示。從圖4中可以看出基于CS的DOA估計方法可以在相干源、低信噪比的條件下基本準確估計目標方位，且分辨力高于經典方法。

圖4 相干源DOA估計

4.2 實驗結果

在消聲水池中進行水下目標方位估計的實驗，目標是一塊長2.35m，寬0.5m的鐵板。利用收發合置形式的主動聲納對目標進行探測，發射中心頻率為12kHz的PCW信號，脈沖寬度為10ms，用16元垂直線陣接收回波信號，陣元間距為0.1m，靜止目標與接收陣之間的水平距離為25m。實驗的目的是從接收到的回波信號中估計目標的空間方位，即目標的俯仰角。分別利用CBF、MVDR、BPD和IHT算法對目標的俯仰角進行估計，結果如圖5所示。其中CBF和MVDR用了100個快拍，兩者估計得到的俯仰角均為-1.3°，BPD和IHT計算了50個快拍，BPD估計的俯仰角均值為-1.3°，IHT估計的俯仰角為-1.5°，這些方法估計得到的結果與實際相符。可以看出，基于CS的DOA估計方法可以在陣元間距大于半波長的情況下有效估計目標的空間方位，尤其是在低快拍數的情況下，且分辨力高于經典的CBF和MVDR方法。

圖5 實驗結果

5 結語

本文將CS應用到線陣的DOA估計上，尤其是陣元間距大于半波長的線陣，把DOA估計轉變為一個欠定線性方程，闡述了基本原理和基本的重構算法BPD和IHT，并分析了傳感矩陣的相干性。數值仿真結果和實驗處理結果均表明了BPD和IHT方法的有效性，尤其是在相干源、低信噪比、低快拍數的情況下，基于CS的DOA估計方法相較傳統的CBF和MVDR等方法具有更高的分辨力，能較好地估計目標空間方位。