基于壓縮感知的低復(fù)雜度超分辨角度估計(jì)方法

摘要： 在空域目標(biāo)角度估計(jì)中，分辨率受陣列孔徑的限制，依靠增加陣元數(shù)量提高分辨率會(huì)增加系統(tǒng)成本。為了在受限陣列尺寸中減少陣元數(shù)量，基于壓縮感知（compressive sensing， CS）理論，提出了一種超分辨角度估計(jì)算法。首先建立陣列接收信號(hào)模型并構(gòu)造冗余字典，然后利用目標(biāo)空間域稀疏先驗(yàn)信息將目標(biāo)角度估計(jì)問(wèn)題建模為優(yōu)化問(wèn)題，最后設(shè)計(jì)低復(fù)雜度算法對(duì)優(yōu)化問(wèn)題求解。所提算法通過(guò)貝葉斯CS理論推導(dǎo)了正則化系數(shù)，保證了算法的噪聲抑制性能，通過(guò)共軛梯度運(yùn)算及Hadamard乘積提高了算法效率。所提算法可利用較少快拍在信號(hào)數(shù)目未知的條件下，實(shí)現(xiàn)高精度角度估計(jì)。仿真結(jié)果和實(shí)測(cè)數(shù)據(jù)驗(yàn)證了所提算法的有效性。

關(guān)鍵詞： 陣列信號(hào)處理; 壓縮感知; 波達(dá)方向估計(jì); 波束形成

中圖分類號(hào)： TN 911.7

文獻(xiàn)標(biāo)志碼： A

DOI：10.12305/j.issn.1001-506X.2024.06.01

Low complexity super-resolution angle estimation method based on compressive sensing

WU Min^1，2， LI Zihao^1，2， HAO Chengpeng^{1，2，*}， HU Qiao³

（1. Institute of Acoustics， Chinese Academy of Sciences， Beijing 100190， China; 2. University of Chinese Academy of Sciences， Beijing 100049， China; 3. School of Mechanical Engineering， Xi’an Jiaotong University， Xi’an 710049， China）

Abstract： In aerial target angle estimation， the resolution is constrained by the aperture length. Increasing the number of array elements to improve the resolution will increase the system lost. To reduce the number of elements in the limited array size， a novel algorithm of super-resolution angle estimation is addressed based on compressive sensing （CS） theory. The array received signal model is established and the redundant dictionary is formed. By exploiting the sparse prior information of the observation area， the target angle estimation problem is converted into the optimization problem. The direction of arrival of targets can be estimated with accuracy via an iterative optimization algorithm. In the proposed algorithm， the regularization coefficient is derived by Bayesian CS theory to ensure the noise robustness of the algorithm. Besides， the efficiency of the proposed algorithm is improved by using the conjugate gradient algorithm and Hadamard product. The effectiveness of the proposed algorithm is verified by simulation and measured data.

Keywords： array signal processing; compressive sensing （CS）; direction of arrival （DOA） estimation; beam forming

0 引 言

隨著聲納技術(shù)的高速發(fā)展，陣列信號(hào)處理可實(shí)現(xiàn)對(duì)水下目標(biāo)的高分辨探測(cè)和定位，其基本原理是利用陣元傳感器接收各個(gè)方向的信號(hào)，估計(jì)信號(hào)的到達(dá)方向（direction of arrival， DOA）角是陣列信號(hào)處理的重要任務(wù)，以獲取信號(hào)源的信源數(shù)、角度等信息為目的^［1］。常見的DOA估計(jì)方法有時(shí)延求和（delay-and-sum， DAS）和最小方差無(wú)畸變（minimum variance distortion-less response， MVDR）算法。DAS通過(guò)補(bǔ)償各通道時(shí)延差后求和得到波束輸出，其計(jì)算復(fù)雜度低，穩(wěn)健性高，但受孔徑大小的限制，空間分辨率較低。Capon^［2］提出的MVDR波束形成器，采用最小化噪聲方差準(zhǔn)則，在感興趣角度無(wú)失真輸出的前提下最小化輸出功率，形成空間譜，采用的線性模型對(duì)于噪聲的利用率較低，導(dǎo)致算法具有一定的局限性。20世紀(jì)70年代末涌現(xiàn)的一些子空間類方法對(duì)DOA估計(jì)研究具有劃時(shí)代意義，其中最具代表性的是多重信號(hào)分類（multiple signal classification， MUSIC）算法^［3］。該算法通過(guò)對(duì)接收數(shù)據(jù)進(jìn)行空間分解，將空間劃分為信號(hào)子空間和噪聲子空間，利用子空間的正交特性構(gòu)造出空間譜圖并進(jìn)行譜峰搜索，從理論上克服了瑞利極限，可獲得目標(biāo)角度的超分辨估計(jì)，但自相關(guān)矩陣特征值分解以及譜峰搜索導(dǎo)致算法計(jì)算復(fù)雜度高?；谛D(zhuǎn)不變的信號(hào)參數(shù)估計(jì)（estimation of signal parameters by rotational invariance techniques， ESPRIT）的子空間類方法^［4］，此算法無(wú)需進(jìn)行譜峰搜索，降低了算法運(yùn)算量，增加了DOA估計(jì)方法的實(shí)用性。近年來(lái)，學(xué)者們?cè)贛USIC算法和ESPRIT算法基礎(chǔ)上進(jìn)行改進(jìn)，提出了許多改進(jìn)的子空間類方法^［5］。然而，子空間類算法通常需要滿足信號(hào)源數(shù)目已知、多快拍以及高信噪比（signal to noise ratio， SNR）的條件，在實(shí)際水下測(cè)向環(huán)境中應(yīng)用受限。同時(shí)，為避免角度估計(jì)模糊，要求陣列中相鄰兩個(gè)陣元間隔不大于信號(hào)半倍波長(zhǎng)。

現(xiàn)代聲納工作頻率范圍較寬，工作頻率決定了陣元間隔以及陣元物理尺寸，當(dāng)信號(hào)頻率過(guò)高時(shí)，較小的陣元間距使得半波長(zhǎng)物理布陣難以實(shí)現(xiàn)，且會(huì)產(chǎn)生較大的陣元互耦。隨著海洋探索的不斷深入，對(duì)于多目標(biāo)、集群目標(biāo)的水下探測(cè)需求日益增加，對(duì)角度分辨率要求越來(lái)越高，高分辨率意味著較大的陣列孔徑，一般需要大量的陣元來(lái)擴(kuò)大陣列的有效孔徑，大量的物理陣元會(huì)增加系統(tǒng)成本和復(fù)雜度。針對(duì)上述關(guān)于陣列尺寸、系統(tǒng)成本以及DOA估計(jì)精度之間存在的問(wèn)題，稀疏陣列因其能夠克服陣元間隔的限制而得到巨大發(fā)展及廣泛應(yīng)用。稀疏陣列構(gòu)型方法打破了傳統(tǒng)均勻陣列構(gòu)型規(guī)則，按照某種規(guī)則進(jìn)行非均勻且稀疏放置，其陣元間隔不受半波長(zhǎng)的限制，因而在同樣陣元條件下，稀疏陣列可獲得比傳統(tǒng)滿陣更大的陣列孔徑，提高空間角度分辨率和測(cè)角精度;換言之，在陣列尺寸一定時(shí)，采用稀疏陣列可以顯著減少陣元數(shù)量，從而降低系統(tǒng)成本。

當(dāng)陣元數(shù)目有限時(shí)，陣元對(duì)空域信號(hào)的采樣將不滿足奈奎斯特采樣定理，常規(guī)DOA估計(jì)算法很難保證角度估計(jì)的性能。因此，如何提高有限元陣列或稀疏陣列的空間增益和角度分辨率成為DOA估計(jì)領(lǐng)域的熱點(diǎn)問(wèn)題。針對(duì)互質(zhì)陣中的DOA估計(jì)，研究人員在傳統(tǒng)MVDR算法的基礎(chǔ)上，對(duì)干擾噪聲協(xié)方差矩陣進(jìn)行重構(gòu)，取得了更優(yōu)的角度估計(jì)性能，但在低SNR下，算法穩(wěn)健性較差^［6］。He等^［7］提出了基于MUSIC算法的空間任意四元陣測(cè)向技術(shù)，但此類算法依然存在運(yùn)算量和快拍數(shù)過(guò)大的問(wèn)題。隨著對(duì)壓縮感知（compressive sensing， CS）理論^［8］和稀疏信號(hào)表征研究的不斷深入，在信號(hào)具有稀疏先驗(yàn)的情況下，可以極大概率在欠采樣條件下恢復(fù)目標(biāo)信號(hào)^［9］，學(xué)者們將稀疏重構(gòu)理論應(yīng)用于超分辨DOA估計(jì)領(lǐng)域。當(dāng)目標(biāo)在空間分布滿足稀疏先驗(yàn)條件、測(cè)量矩陣滿足約束等容（restricted isometry property， RIP）條件時(shí)，可實(shí)現(xiàn)超分辨DOA估計(jì)^［10-12］。文獻(xiàn)［13］提出了一種單快拍自適應(yīng)角度估計(jì)算法，可在單快拍情況下，得到較為理想的估計(jì)結(jié)果，但是單快拍DOA估計(jì)在低SNR條件下角度估計(jì)穩(wěn)定性差。為了將稀疏重構(gòu)思想應(yīng)用于多快拍數(shù)據(jù)場(chǎng)景中，Cotter等^［14］提出基于多測(cè)量矢量（multiple measurement vector， MMV）模型的角度估計(jì)算法，由于陣列獲取的是多快拍數(shù)據(jù)，算法在低SNR的情況下能以更高概率得到高精度DOA估計(jì)，但其運(yùn)算復(fù)雜度隨著快拍數(shù)的增加而增大。Malioutov等^［15］提出了基于奇異值分解（singular value decomposition， SVD）的l₁-SVD的DOA估計(jì)算法。該算法首先將矩陣進(jìn)行SVD，再將其轉(zhuǎn)換為二階錐規(guī)劃問(wèn)題，最后采用內(nèi)點(diǎn)法等相關(guān)凸優(yōu)化方法求解得到角度估計(jì)值。算法可解決多快拍情況下復(fù)雜度過(guò)大的問(wèn)題。然而，由于SVD技術(shù)需要正確的信源個(gè)數(shù)給予支撐，在信源個(gè)數(shù)未知的情況下，該算法的性能將嚴(yán)重下降甚至失效。

本文提出了一種基于CS的超分辨DOA估計(jì)算法，利用場(chǎng)景稀疏性，構(gòu)建基于CS的DOA估計(jì)模型，構(gòu)造測(cè)量矩陣和l₁范數(shù)的優(yōu)化函數(shù)，設(shè)計(jì)低復(fù)雜度優(yōu)化算法對(duì)空間譜進(jìn)行重構(gòu)，算法可在信源個(gè)數(shù)未知和有限陣元條件下得到超分辨角度估計(jì)值。在構(gòu)造優(yōu)化問(wèn)題約束條件時(shí)，利用了噪聲背景下的多快拍數(shù)據(jù)，算法可有效抑制噪聲;利用貝葉斯CS推導(dǎo)了正則化系數(shù)，算法能自適應(yīng)獲得空間譜，算法穩(wěn)健;利用修正牛頓法求解優(yōu)化問(wèn)題，算法運(yùn)算效率大大提高。

1 陣列接收信號(hào)模型

假設(shè)接收陣列為M陣元的隨機(jī)線性陣列，采用直角坐標(biāo)表示。假設(shè)空間中存在N個(gè)方向的遠(yuǎn)場(chǎng)信號(hào)，所有傳感器均為各向同性。第m個(gè)陣元的輸出可表示為

y_m（t）=∑N－1n=0s（t－τ_m，n）+ε_m（t）（1）

式中：τ_m，n為遠(yuǎn)場(chǎng)信號(hào)抵達(dá)第m個(gè)陣元相較于抵達(dá)參考點(diǎn)的時(shí)延，τ_m，n=R_m，n/c，c為信號(hào)傳播速度;ε_m（t）為噪聲，這里假定噪聲與信號(hào)是統(tǒng)計(jì)獨(dú)立的。對(duì)于窄帶情形，時(shí)延等價(jià)于相移，即s（t－τ_m，n）=s（t）e^－j2πf_c^τ_m，n，f_c為窄帶信號(hào)的載頻，則第m個(gè)陣元輸出為

y_m（t）=∑N－1n=0s（t）e^－j2πf_c^τ_m，n+ε_m（t）（2）

M個(gè)陣元輸出信號(hào)的矢量形式為

y（t）=∑N－1n=0a（θ_n）s（t）+ε（t）=As（t）+ε（t）（3）

式中：y（t）=［y₁（t），y₂（t），…，y_M（t）］^T為陣列輸出;s（t）=［s₀（t），s₁（t），…，s_N-1（t）］^T代表空間譜;ε（t）=［ε₁（t），ε₂（t），…，ε_M（t）］^T為噪聲;A=［a（θ₀），a（θ₁），…，a（θ_N－1）］為導(dǎo)向矩陣，a（θ_n）=［e^－j2πf_c^τ_1，n，e^－j2πf_c^τ_2，n，…，e^－j2πf_c^τ_M，n］^T代表導(dǎo)向矢量，表征了目標(biāo)方向和各陣元位置之間的相位差異。

2 基于CS的超分辨角度估計(jì)

2.1 陣列角度估計(jì)模型

觀測(cè)角度范圍被分成H個(gè)角度單元，對(duì)應(yīng)的導(dǎo)向矩陣變?yōu)榇笮镸×H的矩陣A=［a（θ₀），a（θ₁），…，a（θ_H－1）］，導(dǎo)向矢量為a（θ_h）=［e^－j2πf_c^τ_1，h，e^－j2πf_c^τ_2，h，…，e^－j2πf_c^τ_M，h］^T，當(dāng)接收數(shù)據(jù)有K個(gè)快拍時(shí)，接收信號(hào)表示為

Y=AS+E（4）

式中：E=［ε（1），ε（2），…，ε（K）］代表大小為M×K的矩陣;Y=［y（1），y（2），…，y（K）］是大小為M×K的矩陣;S代表需要恢復(fù)的二維空間譜，第n行第k列元素為s_n（k）。通過(guò)向量化二維信號(hào)模型轉(zhuǎn)化為一維形式：

y=A-s+ε（5）

式中：A-=diag{A，A，…，A}是分塊矩陣，大小為MK×HK;y=vec（Y）;s=vec（S）;ε=vec（E），vec（·）指矩陣向量化，表示將矩陣的各列向量依次堆疊形成一個(gè)列向量。假設(shè)在多次快拍觀測(cè)中，目標(biāo)角度不發(fā)生變化，即空時(shí)不耦合。因此，對(duì)空時(shí)譜在時(shí)間上做累積可作為空間譜估計(jì)s-（h）=∑_k|s_h（k）|，通過(guò)對(duì)譜估計(jì)結(jié)果s-（h）做檢測(cè)，可得到目標(biāo)的角度估計(jì)值和個(gè)數(shù)。目標(biāo)場(chǎng)景中只在少數(shù)觀測(cè)角度中有目標(biāo)，即Nlt;lt;H，空間譜具有強(qiáng)稀疏特性，利用這一特性可將空間譜估計(jì)過(guò)程轉(zhuǎn)化為稀疏重構(gòu)過(guò)程。在實(shí)際應(yīng)用中，陣元數(shù)通常小于目標(biāo)區(qū)域角度劃分網(wǎng)格數(shù)，A-是一個(gè)高維到低維的測(cè)量矩陣，從有限回波數(shù)據(jù)中恢復(fù)空間譜是一個(gè)病態(tài)過(guò)程。CS理論表明，在s為稀疏先驗(yàn)條件時(shí)，如果測(cè)量矩陣A-滿足RIP特性，則可以從接收信號(hào)y中恢復(fù)未知空時(shí)譜矩陣s。滿足RIP特性，即存在限制等距常數(shù)δ∈（0，1），使條件（1－δ）s22≤A -s22≤（1+δ）s22成立，RIP特性與構(gòu)成測(cè)量矩陣的基向量的相關(guān)性有關(guān)，相關(guān)系數(shù)越小，基的正交性越好。A～={a～（θ_h）}^H－1_h=0為歸一化的基向量集合，相互之間的相關(guān)系數(shù)為

R（A～）=maxp≠q0≤p，q≤H－1|lt;a～（θ_p），a～（θ_q）gt;|（6）

當(dāng)空間譜矩陣滿足以下條件時(shí)，滿足RIP特性。

在滿足稀疏先驗(yàn)以及RIP特性條件下，可將上述病態(tài)求解問(wèn)題轉(zhuǎn)化為l₁范數(shù)的優(yōu)化問(wèn)題，從而得到目標(biāo)角度信息?？紤]噪聲影響，構(gòu)建的l₁范數(shù)優(yōu)化問(wèn)題為

s=arg min{y-A-s22+ρs₁}（8）

式中：ρ為正則化系數(shù);·2和·₁分別表示l₂范數(shù)和l₁范數(shù);min（·）是最小化運(yùn)算。

2.2 稀疏優(yōu)化求解算法實(shí)現(xiàn)

利用優(yōu)化求解算法，如正交匹配追蹤（orthogonal matching pursuit， OMP）算法^［16］、BFGS（Broyden，F(xiàn)letcher，Gold-farb，Shanno）算法等^［17］可求解l₁范數(shù)優(yōu)化問(wèn)題，通常采用CVX、SEDUMI等內(nèi)點(diǎn)法求解器，不利于工程應(yīng)用^［18-20］，本文推導(dǎo)基于擬牛頓算法的低復(fù)雜度方法求解式（8）中的優(yōu)化問(wèn)題。為提高算法的噪聲穩(wěn)健性和環(huán)境適應(yīng)能力，首先推導(dǎo)了優(yōu)化問(wèn)題中的正則化系數(shù)ρ，假設(shè)噪聲矩陣E服從復(fù)高斯分布，均值為零，方差為σ²，ε的概率密度函數(shù)（probability density function， PDF）可表示為

因此，當(dāng)已知s時(shí)，y的概率密度函數(shù)可表示為

基于貝葉斯CS原理，在空間譜具有強(qiáng)稀疏性的前提下，空間譜服從Laplace分布，s的概率密度函數(shù)可表示為

式中：V是Laplace函數(shù)尺度系數(shù)。在已知陣列接收矩陣y時(shí)，可利用最大后驗(yàn)概率準(zhǔn)則估計(jì)空時(shí)譜，s的最大后驗(yàn)概率估計(jì)為

將式（10）和式（11）代入式（12）推導(dǎo)可得

通過(guò)以上的貝葉斯CS的推導(dǎo)，正則化系數(shù)為ρ=2σ²V，下面介紹估計(jì)σ²和V的方法。首先，利用常規(guī)恒虛警檢測(cè)進(jìn)行粗檢測(cè)，超過(guò)檢測(cè)門限的區(qū)域判定為目標(biāo)區(qū)域，低于檢測(cè)門限的區(qū)域判定為噪聲區(qū)域。由于假設(shè)噪聲滿足零均值的高斯分布，其方差σ²可通過(guò)對(duì)噪聲區(qū)域進(jìn)行最大似然估計(jì)確定^［21］。由于s服從Laplace分布，因此V的最大似然可通過(guò)求max ln［PDF（s|V）］確定，也就等價(jià)于求max{ ln［PDF（s|V）］}，對(duì)ln［PDF（s|V）］關(guān)于V求導(dǎo)可得

?ln［PDF（s|V）］?σ=HKσ-∑HK-1i=0|s（i）|（14）

則V的最大似然估計(jì)為

V=HK∑HK－1i=0|s（i）|（15）

通過(guò)以上的推導(dǎo)，估計(jì)了式（13）中正則化系數(shù)ρ，由于在推導(dǎo)中考慮了目標(biāo)和噪聲的分布特性，因此算法具有較強(qiáng)的噪聲抑制性能。

下面對(duì)式（13）的優(yōu)化函數(shù)進(jìn)行求解，對(duì)代價(jià)函數(shù)J（s）=y－A-s22+ρs₁關(guān)于s求共軛梯度為

ΔJ_s*（s）=H（s）s-2A-^Hy（16）

式中：H（s）矩陣是s的系數(shù)矩陣，可表示為

H（s）=2A-^HA-+ρΛ（s）（17）

式中：Λ（s）=diag［1/|s（hk）|²+τ］。傳統(tǒng)的牛頓法和擬牛頓法不適合處理包含非二次型的優(yōu)化問(wèn)題^［22］，特別是代價(jià)函數(shù)中包含了s的模，更增加了求解的復(fù)雜性。本文中用修正的牛頓迭代算法來(lái)求解該優(yōu)化問(wèn)題^［23］。由牛頓法可知，s的第l+1次迭代式為

s^l+1=s^l－β［H（s^l）］^－1ΔJ_s*（s^l）（18）

式中：β為迭代步長(zhǎng)，［·］^－1表示矩陣的逆。如迭代步長(zhǎng)設(shè)為β=1，把共軛梯度式（16）代入式（18）推導(dǎo)得

s^l+1=2［H（s^l）］^－1A-Hy（19）

當(dāng)?shù)鷿M足以下條件時(shí)，算法終止

s^l+1－s^l22s^l22lt;δ（20）

將二維空時(shí)譜在時(shí)間域上做累積即可得到目標(biāo)角度空間譜

s-（h）=∑k|s_h（k）|（21）

對(duì)s-進(jìn)行常規(guī)的恒虛警檢測(cè)能獲得目標(biāo)的角度估計(jì)結(jié)果，其位置和個(gè)數(shù)分別表示目標(biāo)的角度估計(jì)值和目標(biāo)個(gè)數(shù)。

2.3 快速算法及運(yùn)算復(fù)雜度分析

以上優(yōu)化求解算法計(jì)算量主要集中于Hessian矩陣H（s^l）的求逆，由上面的推導(dǎo)可知，H（s^l）是大小為HK×HK的矩陣，用常用的Cholesky分解，每次迭代中均需要運(yùn)行（HK）³/3+2（HK）²次乘法運(yùn)算，算法運(yùn)算量過(guò)大。

為解決Hessian矩陣求逆所引起的運(yùn)算量大的問(wèn)題，本文算法在第l+1次迭代運(yùn)算中求解s^l+1時(shí)，采用共軛梯度法（conjugate gradient algorithm， CGA）求解。

H（s^l）s^l+1=2A-^Hy（22）

Hessian矩陣在每次迭代中均需要進(jìn)行更新H（s^l）=2A-^HA-+ρΛ（s^l），在第l+1次迭代中，CGA需要完成下列方程的多次重復(fù)計(jì)算：

H（s^l）s^=2A-^HA-s^+ρΛ（s^l）s^（23）

式中：s^為迭代更新值。利用分塊對(duì)角矩陣的性質(zhì)，式（23）第1項(xiàng)可表示為

2A-^HA-s^=2diag{A^H …A^H}diag{A …A}s^=

2diag{A^HA …A^HA}s^=2s^（24）

式（23）第2項(xiàng)中Λ（s^l）可表示為對(duì)角矩陣Λ（s^l）=diag［1/|s（hk）|²+τ］，因此Λ（s^l）s^是對(duì)角線上的HK個(gè)元素組成的向量與向量s^的Hadamard乘積，通過(guò)HK次乘法即能獲得矩陣相乘結(jié)果。

由于優(yōu)化算法中主要運(yùn)算量為乘法運(yùn)算，利用所提出的快速求解算法，算法運(yùn)算復(fù)雜度為O（HK），相比于直接矩陣運(yùn)算，運(yùn)算復(fù)雜度降低了兩個(gè)數(shù)量級(jí)。

BFGS和OMP也可以利用上述的快速矩陣乘法來(lái)提高運(yùn)算效率。對(duì)于BFGS算法來(lái)說(shuō)，算法主要運(yùn)算量為每次迭代中的內(nèi)積運(yùn)算以及ΔJ_s*（s^l）的更新，由以上推導(dǎo)可知，ΔJ_s*（s）=H（s）s-2A-^Hy，H（s）s可利用上述的快速矩陣運(yùn)算提高運(yùn)算效率，BFGS算法的運(yùn)算量為O（HK²）。OMP算法的主要運(yùn)算量集中于最大相關(guān)運(yùn)算和最小二乘估計(jì)，在每次迭代中，最大相關(guān)法運(yùn)算量為O（HK），在第l次迭代運(yùn)算中，最小二乘法運(yùn)算量為O（HKl），因此OMP的運(yùn)算量為O（HK）。

3 實(shí)驗(yàn)仿真與結(jié)果分析

仿真實(shí)驗(yàn)中聲速為1 500 m/s，信號(hào)載頻為1 kHz，半波長(zhǎng)為0.75 m。假定陣列為長(zhǎng)度為24 m的隨機(jī)布放直線陣列，陣元個(gè)數(shù)M=8。各陣元接收數(shù)據(jù)中加入均值為0的高斯白噪聲，SNR為0～10 dB。

首先進(jìn)行單目標(biāo)仿真，水聲目標(biāo)的真實(shí)角度為0°，快拍次數(shù)K=16，SNR為10 dB，圖1給出了DAS、MVDR、MUSIC以及所提CS算法的空間譜估計(jì)結(jié)果。本實(shí)驗(yàn)的CS算法以DAS算法結(jié)果為初始值，通過(guò)第2.2節(jié)介紹的正則化系數(shù)估計(jì)方法可得ρ=1.5，設(shè)τ=0.01，［-90°， 90°］角度范圍內(nèi)劃分網(wǎng)格數(shù)為360。從圖1中可以看出，對(duì)于稀疏布放陣列，基于CS的DOA估計(jì)算法能準(zhǔn)確估計(jì)出目標(biāo)角度，峰值旁瓣電平約為-70 dB，遠(yuǎn)低于DAS、MVDR和MUSIC算法。表1給出了各算法進(jìn)行200次蒙特卡羅實(shí)驗(yàn)的平均運(yùn)算時(shí)間，結(jié)果表明，DAS算法運(yùn)算復(fù)雜度最低，MVDR算法由于矩陣求逆操作，運(yùn)算復(fù)雜度高于DAS算法，MUSIC算法由于譜峰搜索，所需計(jì)算時(shí)間遠(yuǎn)高于其他算法，而所提CS算法利用了CGA和Hadamard乘積，可獲得較高的運(yùn)算效率。

圖2給出了角度均方誤差隨SNR的變化，實(shí)驗(yàn)次數(shù)為200次，SNR在0 dB到10 dB之間，間隔為2 dB。由于陣列是稀疏的，在低SNR下，DAS算法失效，無(wú)法準(zhǔn)確估計(jì)目標(biāo)角度。MUSIC算法由于快拍數(shù)較少以及目標(biāo)數(shù)量未知，無(wú)法獲得有效的目標(biāo)角度，MVDR算法和所提CS算法的角度估計(jì)均方誤差較低，在低SNR下依然能獲得高精度角度估計(jì)結(jié)果。

接下來(lái)，為驗(yàn)證算法的稀疏恢復(fù)性能，在不同陣元數(shù)下對(duì)各算法進(jìn)行實(shí)驗(yàn)分析。圖3給出了角度估計(jì)均方誤差與陣元數(shù)之間的關(guān)系，選取陣元數(shù)6到16進(jìn)行實(shí)驗(yàn)，陣元隨機(jī)布放，SNR為10 dB，快拍數(shù)為16。從結(jié)果可以看出，在陣元數(shù)大于10時(shí)，DAS、MVDR以及CS算法角度估計(jì)均方誤差小，MUSIC在低快拍條件下欠穩(wěn)健，而本文所提CS算法在陣元數(shù)為6時(shí)仍能以較高精度獲得目標(biāo)角度。結(jié)果表明，所提算法在保證高精度角度估計(jì)的前提下可降低陣元數(shù)量。

最后，通過(guò)多目標(biāo)角度估計(jì)驗(yàn)證算法的超分辨性能，圖4（a）中的兩個(gè)目標(biāo)角度分別為-5°和5°，圖4（b）中兩個(gè)目標(biāo)角度分別為-3°和3°。從圖4（a）中的結(jié)果可以看出，CS算法估計(jì)結(jié)果主瓣最尖銳，角度分辨率最高，MUSIC由于目標(biāo)數(shù)量未知以及快拍數(shù)較少，無(wú)法準(zhǔn)確恢復(fù)目標(biāo)角度，雖然DAS算法和MVDR算法可恢復(fù)出兩個(gè)目標(biāo)，但是由于旁瓣較高，在兩個(gè)目標(biāo)間產(chǎn)生了一個(gè)虛假目標(biāo)，當(dāng)兩個(gè)目標(biāo)進(jìn)一步靠近，角度為-3°和3°時(shí)，DAS和MVDR算法均已無(wú)法分辨兩個(gè)目標(biāo)，而本文所提算法依然能獲得良好的角度估計(jì)結(jié)果。結(jié)果表明，對(duì)于多目標(biāo)應(yīng)用場(chǎng)景，所提算法可實(shí)現(xiàn)臨近目標(biāo)的超分辨估計(jì)。圖5給出了角度估計(jì)的均方誤差隨兩個(gè)角度間隔（角度分辨率）的變化曲線，實(shí)驗(yàn)次數(shù)為200次，角度分辨率在6°到18°之間，間隔為4°。隨著角度間隔的降低，算法的估計(jì)精度越來(lái)越低，所提CS算法利用稀疏優(yōu)化重構(gòu)，在兩個(gè)目標(biāo)相鄰較近時(shí)，依然能獲得高精度角度估計(jì)結(jié)果，算法超分辨效果顯著。

4 實(shí)驗(yàn)數(shù)據(jù)處理分析

為驗(yàn)證提出的基于CS的超分辨角度估計(jì)算法的性能，采用千島湖的某次湖試的數(shù)據(jù)進(jìn)行處理和分析。聲速剖面如圖6所示，實(shí)驗(yàn)場(chǎng)景中分別有單個(gè)和兩個(gè)實(shí)測(cè)目標(biāo)，采用的換能器組組成了一個(gè)6元均勻直線陣列，陣列水平排布，陣元間隔為0.2 m。試驗(yàn)當(dāng)天水面浪高小于0.3 m且陣風(fēng)，無(wú)過(guò)往船只，SNR約為-8 dB。系統(tǒng)發(fā)射信號(hào)為連續(xù)波，利用16個(gè)快拍數(shù)據(jù)對(duì)目標(biāo)進(jìn)行角度估計(jì)。正則化系數(shù)估計(jì)值為ρ=2.3。

圖7給出的是實(shí)驗(yàn)數(shù)據(jù)的角度估計(jì)的結(jié)果。由圖7可知，當(dāng)場(chǎng)景中只有單個(gè)目標(biāo)時(shí)，采用傳統(tǒng)的DAS估計(jì)技術(shù)時(shí)分辨率低，水下強(qiáng)噪聲使算法很難產(chǎn)生足夠明顯的峰值，從而嚴(yán)重影響了MUSIC算法的角度估計(jì)結(jié)果。相比于MUSIC算法，MVDR算法可得到高峰值空間譜，峰值旁瓣電平約為-50 dB，通過(guò)利用稀疏重構(gòu)方法，基于CS的算法即使在強(qiáng)噪聲水下環(huán)境中，依然能獲得目標(biāo)的高峰值空間譜，峰值旁瓣電平約為-70 dB。

圖8（a）給出了角度為-5°和5°的多目標(biāo)角度估計(jì)結(jié)果，圖8（b）給出了角度為-3°和3°的多目標(biāo)角度估計(jì)結(jié)果，通過(guò)多目標(biāo)角度估計(jì)結(jié)果驗(yàn)證算法的超分辨性能。當(dāng)兩個(gè)目標(biāo)角度為-5°和5°時(shí)，MUSIC算法無(wú)法確定目標(biāo)點(diǎn)位置峰值，DAS算法和MVDR算法雖然在目標(biāo)位置有峰值，但是在其他角度也同樣產(chǎn)生了高峰值空間譜，這會(huì)造成后續(xù)檢測(cè)中出現(xiàn)虛警。本文所提算法可獲得兩個(gè)目標(biāo)角度的高分辨估計(jì)結(jié)果，峰值旁瓣電平約為-60 dB。當(dāng)兩個(gè)目標(biāo)角度繼續(xù)靠近，角度為-3°和3°時(shí)，圖6（b）中DAS、MVDR以及MUSIC算法在非目標(biāo)位置產(chǎn)生了高峰值空間譜，大量的虛假峰值造成檢測(cè)算法無(wú)法獲得正確的目標(biāo)角度信息，所提CS算法依然可以分辨出兩個(gè)目標(biāo)，峰值旁瓣電平約為-50 dB。

5 結(jié)束語(yǔ)

本文提出了基于CS的超分辨DOA估計(jì)算法，能有效提高目標(biāo)角度估計(jì)精度，具有復(fù)雜度低、分辨率高、魯棒性強(qiáng)的優(yōu)點(diǎn)。該算法利用目標(biāo)空間域稀疏先驗(yàn)將目標(biāo)角度估計(jì)問(wèn)題建模為優(yōu)化問(wèn)題，并通過(guò)貝葉斯CS理論推導(dǎo)了正則化系數(shù)，保證了算法的噪聲抑制性能，設(shè)計(jì)了高效的迭代算法，運(yùn)算復(fù)雜度低。算法可利用較少陣元數(shù)和較少快拍實(shí)現(xiàn)高精度的目標(biāo)角度估計(jì)。該算法性能優(yōu)于傳統(tǒng)算法，特別是在低SNR、有限陣元以及較少快拍條件下優(yōu)于MVDR算法和MUSIC算法，有更尖銳的空間譜峰值。仿真和實(shí)測(cè)實(shí)驗(yàn)驗(yàn)證了該算法的有效性。未來(lái)將會(huì)根據(jù)水下探測(cè)的實(shí)際需求，將算法推廣到多種陣型，并且聯(lián)合距離維信息，實(shí)現(xiàn)對(duì)多目標(biāo)、群目標(biāo)的高精度位置估計(jì)。

參考文獻(xiàn)

［1］陳世進(jìn)， 閆晟， 郝程鵬， 等. 一種適用于多輸入多輸出聲吶的穩(wěn)健空時(shí)自適應(yīng)檢測(cè)方法［J］. 聲學(xué)學(xué)報(bào)， 2022， 47（6）： 777-788.

CHEN S J， YAN S， HAO C P， et al. A robust space-time adaptive detection method for multiple-input multiple-output sonar［J］. Chinese Journal of Acoustics， 2022， 47（6）： 777-788.

［2］CAPON J. High-resolution frequency-wavenumber spectrum analysis［J］. Proceedings of the IEEE，1969， 57（8）： 1408-1418.

［3］SCHMIDT R. Multiple emitter location and signal parameter estimation［J］. IEEE Trans.on Antennas and Propagation， 1986， 34（3）： 276-280.

［4］VANDERVEEN M C， VEEN A J， PAULRAJ A. Estimation of multipath parameters in wireless communications［J］. IEEE Trans.on Signal Processing， 1998， 46（3）： 682-690.

［5］RAO B D， HARI K V S. Performance analysis of root-MUSIC［J］. IEEE Trans.on Acoustics Speech and Signal Processing， 1989， 37（12）： 1939-1949.

［6］ZHOU C， SHI Z， GU Y. Coprime array adaptive beamforming with enhanced degrees of freedom capability［C］∥Proc.of the IEEE Radar Conference， 2017.

［7］HE K， XIANG J， HAN P， et al. Performance analysis of MUSIC for arbitary array［J］. Computer Measurement and Control， 2010， 18（3）： 688-690.

［8］DONOHO D. Compressed sensing［J］. IEEE Trans.on Information Theory， 2006， 52（4）： 5406-5425.

［9］CANDES E， ROMBERG J， TAO T. Stable signal recovery from incomplete and inaccurate measurements［J］. Communications on Pure and Applied Mathematics， 2006， 59（8）： 1207-1223.

［10］ZHENG Z， YANG T C， GERSTOFT P. Joint towed array shape and directions of arrival estimation using sparse Bayesian learning during maneuvering［J］. The Journal of the Acoustical Society of America， 2020， 147（3）： 1738-1751.

［11］PING G， FERNANDEZ-GRANDE E， GERSTOFT P. Three-dimensional source localization using sparse Bayesian learning on a spherical microphone array［J］. The Journal of the Acoustical Society of America， 2020， 147（6）： 3895-3904.

［12］虞飛， 余赟， 周利輝， 等. 一種不依賴超參數(shù)的稀疏信號(hào)單快拍DOA估計(jì)方法［J］. 系統(tǒng)工程與電子技術(shù)， 2021， 43（4）： 894-900.

YU F， YU Y， ZHOU L H， et al. Hyperparameter-free sparse signal direction-of-arrival estimation method with single-snapshot［J］. Systems Engineering and Electronics， 2021， 43（4）： 894-900.

［13］李洪濤， 賀亞鵬， 顧陳， 等， 基于壓縮感知的單快拍自適應(yīng)波束形成算法［J］. 南京理工大學(xué)學(xué)報(bào)， 2012， 36（1）： 91-95.

LI H T， HE Y P， GU C， et al. Compressive sensing based single-snapshot adaptive beamforming algorithm［J］. Journal of Nanjing University of Science and Technology， 2012， 36（1）： 91-95.

［14］COTTER S F， RAO B D， ENGAN K， et al. Sparse solutions to linear inverse problems with multiple measurement vectors［J］. IEEE Trans.on Signal Processing， 2005， 53（7）： 2477-2488.

［15］MALIOUTOV D， ?ETIN M， WILLSKY A S. A sparse signal reconstruction perspective for source localization with sensor arrays［J］. IEEE Trans.on Signal Processing， 2005， 53（8）： 3010-3022.

［16］TROPP J A， GILBERT A C. Signal recovery from random measurements via orthogonal matching pursuit［J］. IEEE Trans.on Information Theory， 2007， 53（12）： 4655-4666.

［17］NOCEDAL J， WRIGHT S J. Numerical optimization［M］. 2nd ed. Springer： New York， 2000.

［18］陳力恒， 馬曉川， 李璇， 等， 結(jié)合壓縮感知模型的稀疏陣列波束形成方法［J］. 信號(hào)處理， 2020， 36（4）： 475-485.

CHEN L H， MA X C， LI X， et al. Sparse array beamforming method combined with compressed sensing model［J］. Journal of Signal Processing， 2020， 36（4）： 475-485.

［19］MALLAT S， ZHANG Z. Matching pursuits with time frequency dictionaries［J］. IEEE Trans.on Signal Processing， 1993， 41（12）： 3397-3415.

［20］任明健， 胡國(guó)平， 周豪， 等. 基于耦合張量分解的稀疏陣列二維DOA估計(jì)算法［J］. 系統(tǒng)工程與電子技術(shù)， 2023， 45（4）： 958-964.

REN M J， HU G P， ZHOU H， et al. DOA estimation algorithm with two-dimensional sparse array based on coupling tensor decomposition［J］. Systems Engineering and Electronics， 2023， 45（4）： 958-964.

［21］JI S H， XUE Y， CARIN L. Bayesian compressive sensing［J］. IEEE Trans.on Signal Processing， 2008， 56（6）： 2346-2356.

［22］ZHOU H， ALEXANDER D， LANGE K. A quasi-Newton acceleration for high-dimensional optimization algorithms［J］. Statistics and Computing， 2011， 21（2）： 261-273.

［23］WU M， ZHANG L， XIA X G， et al. Phase adjustment for polarimetric ISAR with compressive sensing［J］. IEEE Trans.on Aerospace and Electronic Systems， 2016， 52（4）： 1592-1606.

作者簡(jiǎn)介

吳 敏（1988—），女，副研究員，博士，主要研究方向?yàn)樗曅盘?hào)處理、目標(biāo)探測(cè)與成像、新體制聲納信號(hào)處理。

黎子浩（1996—），男，博士研究生，主要研究方向?yàn)樗曅盘?hào)處理。

郝程鵬（1975—），男，研究員，博士，主要研究方向?yàn)樗曅盘?hào)處理、陣列信號(hào)處理、信號(hào)檢測(cè)。

胡 橋（1977—），男，研究員，教授，博士研究生導(dǎo)師，博士，主要研究方向?yàn)樗憙蓷律鷻C(jī)器人、水下機(jī)器人集群信號(hào)處理。