陣元失效下稀疏陣列的二維DOA估計算法

司偉建 馬萬禹 姚璐 曲明超 梁義魯

摘 要：本文針對二維稀疏陣列在陣元失效條件下， 因數(shù)據(jù)缺失導致虛擬陣列連續(xù)性被破壞及自由度下降的問題， 提出了一種二維DOA估計算法。 首先基于二維差分共陣構建虛擬陣列， 然后利用解耦原子范數(shù)最小化理論， 以矩陣填充的形式恢復協(xié)方差矩陣數(shù)據(jù)， 實現(xiàn)對虛擬陣列中丟失虛擬陣元的內插， 最后采用SS-MUSIC算法進行多信源的二維DOA估計。 所提方法彌補了物理陣元失效所造成的影響， 恢復了原始虛擬陣列的完整孔徑特性， 保持了虛擬陣列的自由度， 從而確保了較高精度的二維DOA估計性能。 仿真實驗結果表明， 在相同陣元數(shù)量及陣元失效情況下， 本文提出的算法相比已有方法能有效地估計更多信源， 并在小快拍數(shù)和低信噪比條件下表現(xiàn)出更高的穩(wěn)健性， 最大限度地保留并利用了稀疏陣列在二維DOA估計中的自由度優(yōu)勢。

關鍵詞：? 二維DOA估計； 稀疏陣列； 差分共陣； 陣元失效； 解耦原子范數(shù)最小化； 矩陣填充

中圖分類號： TJ765； TN911.7

文獻標識碼： A

文章編號：? 1673-5048（2024）02-0114-09

DOI： 10.12132/ISSN.1673-5048.2024.0042

0 引? 言

波達方向（Direction of Arrival， DOA）估計是陣列信號處理領域內的一個極其重要的分支， 長期以來受到學者們的廣泛關注。 在過去的幾十年間， DOA估計理論不斷發(fā)展， 從雷達探測、 無源定位、 無線通信、 衛(wèi)星導航， 乃至物聯(lián)網環(huán)境中的智能感知， DOA估計技術在眾多領域均展現(xiàn)出廣泛而深遠的應用價值［1-3］。 以多重信號分類（Multiple Signal Classification， MUSIC）算法［4］為代表的子空間分解類算法憑借高精度和超分辨的DOA估計能力為人所熟知， 并獲得了廣泛的應用與實踐驗證。 然而這類算法的信號接收模型中通常依賴于均勻陣列， 其理論上可估計的最大信源數(shù)目受限于陣元的個數(shù)。

因此， 在日益復雜的空間電磁環(huán)境下， 針對信號密集場景的精確測向問題， 基于稀疏陣列的欠定DOA估計成為新的研究熱點［5］， 稀疏陣列突破了陣元間距必須為入射信號半波長整數(shù)倍的限制， 通過設計陣列構型以構建大孔徑虛擬陣列， 從而顯著提升了陣列自由度（Degree of Freedom， DOF）， 進而能夠實現(xiàn)超過實際物理陣元數(shù)量的DOA估計［6］。 Pal等學者開創(chuàng)性地提出嵌套陣列的概念， 并基于此發(fā)明出一種欠定DOA估計算法［7］。 該算法首先通過對原始接收數(shù)據(jù)進行協(xié)方差矩陣的向量化處理， 形成與虛擬陣列相對應的單快拍數(shù)據(jù)， 隨后執(zhí)行空間平滑（Spatial Smoothing， SS）操作獲得滿秩的協(xié)方差矩陣， 使用MUSIC算法完成角度估計。 隨后Pal等人進一步拓展， 又提出了適用于二維空間的平面嵌套陣列設計及相應DOA估計算法［8-9］。 盡管SS-MUSIC算法通?？梢詰糜诖蠖鄶?shù)類型的稀疏陣列DOA估計中， 其有效性卻高度依賴于虛擬陣列結構的連續(xù)性［10］， 一旦虛擬陣列呈現(xiàn)出不連續(xù)性， 即其間存在孔洞， 那么陣列的連續(xù)自由度擴展將受到顯著限制。

通過使用協(xié)方差矩陣重構技術， 能夠對虛擬孔洞位置的空缺數(shù)據(jù)信息進行近似化填充， 從而擴大虛擬陣列孔徑。 一系列基于互質線陣內插的DOA估計算法被提出， 文獻［11］通過核范數(shù)最小化構建凸優(yōu)化模型， 填補虛擬陣列孔洞并提升DOA估計精度， 但核范數(shù)方法可能引入一定近似誤差； 文獻［12］則進一步將核范數(shù)最小化問題轉化為跡最小化以求解， 實現(xiàn)離網格DOA估計； 文獻［13］對虛擬陣列內插算法進行了完整分析， 優(yōu)化模型限制了此過程中的誤差積累。 近年來， 原子范數(shù)理論因其對信號空域稀疏特性的有效利用， 在DOA估計領域取得了顯著進展， 原子范數(shù)最小化（Atomic Norm Minimization， ANM）作為一種無網格方法， 避免了一般壓縮感知方法中的基失配問題。 文獻［14］將ANM算法應用于互質線陣的虛擬陣列內插， 實現(xiàn)了高精度的無網格DOA估計。 文獻［15-16］則針對二維DOA估計提出解耦原子范數(shù)最小化（Decoupled ANM， DANM）算法， 基于二維陣列接收數(shù)據(jù)模型， 將二維DOA估計問題轉化為兩個獨立的一維DOA估計； 文獻［17］中對DANM算法做出改進， 使之適用于多快拍情形下的二維DOA估計。 但這兩種方法都局限于均勻矩形陣列， 且最大可估計信源數(shù)受限于較短一邊的陣元數(shù)。 文獻［18］基于互質平面陣列， 采用DANM算法進行二維虛擬陣列內插， 可實現(xiàn)欠定的二維DOA估計。

航空兵器 2024年第31卷第2期

司偉建， 等： 陣元失效下稀疏陣列的二維DOA估計算法

在實際應用環(huán)境中， 由于硬件老化、 環(huán)境干擾等多種難以預測的因素導致性能退化乃至完全失效。 一旦陣元隨機發(fā)生故障， 虛擬陣列的連續(xù)性很可能會因此受到顯著破壞， 形成位置和大小不定的孔洞， 進而大幅降低連續(xù)自由度， 并對DOA估計精度產生嚴重影響［19］。 面對陣元失效情況下的DOA估計難題， 類似地， 可通過內插方法填補陣元失效后虛擬陣列中的孔洞， 近似恢復丟失虛擬陣元的數(shù)據(jù)信息， 以維持虛擬陣列的連續(xù)性， 從而最大限度地利用陣列原始的完整孔徑優(yōu)勢。 文獻［20］采用Toeplitz矩陣重構技術恢復協(xié)方差矩陣； 而文獻［21］則基于重加權l(xiāng)2， 1范數(shù)最小化算法恢復了完整的陣列數(shù)據(jù)矩陣， 實現(xiàn)了一維DOA的有效估計。 然而這兩種方法在將算法擴展到二維場景時仍面臨較大困難。

鑒于當前研究現(xiàn)狀， 陣元失效問題的研究重心主要集中在一維陣列上， 而在二維稀疏陣列DOA估計領域， 關于陣元失效下的虛擬陣列孔洞內插問題， 現(xiàn)有的研究尚不充分。 因此， 本文提出一種陣元失效下的二維稀疏陣列DOA估計算法， 該算法利用二維差分共陣形成虛擬陣列， 通過建立基于DANM理論的矩陣填充問題， 恢復協(xié)方差矩陣數(shù)據(jù)， 對陣元失效導致的孔洞內插， 以恢復原虛擬陣列的完整自由度， 并結合SS-MUSIC算法進行二維DOA估計。 對比傳統(tǒng)算法， 本文方法能夠避免舍棄部分虛擬陣元而造成的自由度損失， 有效維持了稀疏陣列二維DOA估計的優(yōu)勢性能。

1 信號模型

考慮空間中K個窄帶遠場獨立信號入射到XOY平面中一個由M個標量天線構成的稀疏平面陣列上， 如圖 1所示。 在平面直角坐標系中， 陣列中任一陣元的位置都可表示為（xid， yid）的形式， 其中： xi， yi∈瘙綄， i=1， 2， …， M； d=λ/2， 即半波長。 令ni=（xi， yi）∈瘙綄2， 代表單個陣元位置， 構成該稀疏平面陣列的所有物理陣元的二維位置整數(shù)集合記為P。

假設第k個信號的入射角度為（θk， φk）， 其中： θk為仰角， 即入射方向與Z軸的夾角， θk∈［0， π/2］； φk為方位角， 即入射方向在XOY平面上的投影與X軸的夾角，? φk=［0， 2π）。 第k個信號的導向矢量可表示為

aP（θk， φk）=［ax1， y1（θk， φk）， …， axM， yM（θk， φk）］T（1）

式中： axi， yi（θk， φk） = ejπ（xicosφksinθk+yisinφksinθk）， i = 1， 2， …， M， j是虛數(shù)單位， k=1， 2， …， K。

假設噪聲為加性高斯白噪聲， 且噪聲與信號相互統(tǒng)計獨立， 則在t時刻時， 稀疏平面陣列的信號接收數(shù)據(jù)模型為

XP（t）=APS（t）+N（t）=∑Kk=1aP（θk，φk）sk（t）+N（t） （2）

式中：? S（t）=［s1（t）， s2（t）， …， sK（t）］T為信號向量；N（t）=［n1（t）， n2（t）， …， nM（t）］T為噪聲向量； 陣列流形矩陣AP=［aP（θ1， φ1）， …， aP（θk， φk）］∈瘙綇M×K。

當考慮全部T個快拍數(shù)據(jù)時， 式（2）轉化為

XP=APS+N （3）

式中： XP∈瘙綇M×T， S∈瘙綇K×T， N∈瘙綇M×T。

計算陣列接收數(shù)據(jù)的協(xié)方差矩陣RP=E［XPXHP］， 其中E［·］表示數(shù)學期望， 將式（3）代入可得：

RP=APRSAHP+σ2NI=∑Kk=1σ2kaP（θk，φk）aHP（θk，φk）+σ2NI（4）

式中： RS為信號的自協(xié)方差矩陣， 其對角線元素表示信號功率； σ2N為輸入噪聲N的方差， I為單位陣。

由于采樣數(shù)量有限， 通常采用最大似然估計方法對式（4）進行近似計算， 可表示為

R^P=1T∑Tt=1XPXHP（5）

對R^P進行向量化操作， 可得

z=vec（R^P）=（A*P⊙AP）p+σ2Nι（6）

式中： p =［p1， …， pK］T，? ι =vec（ I ）； ⊙表示Khatri-Rao積。 可將 z 視作一個單快拍的虛擬陣列接收信號矢量，? A*P⊙ AP則是虛擬陣列的流形矩陣， 大小為M2×K， 其第k列為 a*P（θk， φk） aP（θk， φk）， 表示Kronecker積， 將式（1）代入后展開可得

a*P（θk， φk）aP（θk， φk）=a*P（～x， k， ～y， k）

aP（～x， k， ～y， k）=［ejπ［（x1－x1）～x， k+（y1－y1）～y， k］， …，

ejπ［（xM-x1）～x， k+（yM-y1）～y， k］， ejπ［（x1-x2）～x， k+（y1-y2）～y， k］， …，

ejπ［（xM-x2）～x， k+（yM-y2）～y， k］， …，? ejπ［（x1-xM）～x， k+（y1-yM）～y， k］， …，

ejπ［（xM－xM）～x， k+（yM－yM）～y， k］］T（7）

式中： 為便于表示， 令～x， k=cosφksinθk， ～y， k=sinφksinθk。 通過觀察式（7）中波程差的形式， 可以將二維虛擬陣列的陣元位置坐標表示為

D={（xi－xj， yi－yj）|（xi， yi）， （xj， yj）∈P，

i， j=1， 2， …， M}（8）

即D為P中任意兩個坐標之間差的集合（包括i=j的情況）， 但集合中元素有唯一性， 同一個虛擬陣元位置可能是由多對物理陣元的位置差而得， 所以D中元素個數(shù)一定小于M2， 但虛擬陣列的陣元數(shù)量較原物理陣列的大幅提升將使二維欠定DOA估計得以實現(xiàn)。

對應地， 式（7）向量中的元素也會有重復， 于是根據(jù)與D中各陣元對應關系， 對式（6）中向量z進行去冗余操作， 即將其中對應同一虛擬陣元的元素取平均值再按對應坐標大小重新排序［10］， 則對應D的單快拍虛擬接收數(shù)據(jù)為

zD=ADp+σ2NιD（9）

式中： AD為虛擬陣列D接收數(shù)據(jù)的流形矩陣， 其大小為|D|×K（|·|表示集合的勢， 即集合中元素的個數(shù)）， ιD表示式（6）中向量化單位陣vec（I）對應D的重排序。

2 陣元失效下二維DOA估計

2.1 二維差分共陣

根據(jù)式（8）對二維虛擬陣列的表示， 可以給出以下概念。 對于由集合P={ni=（xi， yi）|i=1， 2， …， M}指定的二維陣列， 其差分共陣D定義為P陣元位置坐標兩兩之間差的集合［22-23］：

D={m|m=ni－nj， ni， nj∈P}（10）

式中： m=（mx， my）∈瘙綄2表示虛擬陣元的平面位置坐標。

如此由M個物理陣元形成的差分共陣的陣元數(shù)最大可達|D|=M（M－1）+1［7］。 由式（10）易知， 差分共陣的陣元位置是關于原點（0， 0）中心對稱的， 所以， 當分別取得虛擬陣元坐標中mx和my的最大和最小值時， 有mxmax=－mxmin， mymax=－mymin。 于是可以定義一個包含D中所有陣元的最小均勻矩形陣列V：

V={（x， y）|－mxmax≤x≤mxmax，

－mymax≤y≤mymax， （x， y）∈瘙綄2}（11）

但大多數(shù)二維稀疏陣列的D并不是一個完整的均勻矩形陣列（Uniform Rectangular Array， URA）， 所以將屬于V但不屬于D的位置坐標集合稱為“孔洞”， 記為H=V＼D。 而當H=時， 有D=V， 則稱D為無孔的虛擬陣列， 在這種情形下， 無需在虛擬陣列中央劃取一個連續(xù)URA而舍棄掉外圍的虛擬陣元， 能夠最大化利用D的自由度。

當陣列中個別陣元出現(xiàn)故障以至于失效時， 陣元位置集合P中對應元素缺失， 其差分共陣D有可能發(fā)生改變， 破壞虛擬陣列中的連續(xù)結構， 使其連續(xù)自由度大幅縮減。 對于β個陣元同時失效的情形， 可以根據(jù)是否對D產生影響來判斷失效陣元組成的子陣列在原陣列中的必要性。 于是定義： 當從P中移除某個陣元數(shù)為β的子陣B后， 使形成的差分共陣D發(fā)生改變時， 即若有P-=P＼B， 使得D-≠D， 稱B對陣列P是β-必要的［22］； 否則， B是β-非必要的。 這里， P-表示移除子陣后的受損陣列， D-為P-的差分共陣。

考慮到不同子陣失效對差分共陣D的影響一般是不同的， 以無孔的D為例， P中的必要陣元失效將導致原本均勻矩形虛擬陣列中產生新的孔洞Hf， 孔洞越大， 說明對D的結構破壞越嚴重， 進而使得在利用矩陣填充理論和最優(yōu)化算法重構數(shù)據(jù)時引入的誤差增大， 這種誤差積累將直接影響二維DOA估計的精度與可靠性。 因此， 定義子陣B在P中的重要度如下：

κD（B）1－D-D=HfD（12）

式中： κD（B）∈［0， 1］， 若B為非必要子陣， 有κD（B）=0。

2.2 陣元失效下虛擬陣列接收數(shù)據(jù)

首先考慮陣列P中非必要子陣失效的情況， 此時其差分共陣D不發(fā)生改變。 使用式（11）中對V的定義， 將虛擬陣列設為一大小為（2mxmax+1）×（2mymax+1）的均勻矩形陣列， D中陣元被完全包含其中， 下面從x和y兩個維度上討論虛擬陣列信號模型的構建。

基于陣列所在XOY平面， 沿X軸和Y軸方向將導向矢量與流形矩陣分為兩部分。 將式（11）中虛擬陣列的坐標代入， 則兩個維度上的陣列導向矢量為

aVx（～x， k）=［e-jπmxmax～x， k， …， 1， …， ejπmxmax～x， k］T（13）

aVy（～y， k）=［e-jπmymax～y， k， …， 1， …， ejπmymax～y， k］T（14）

之后， 將兩個陣列流形矩陣AVx和AVy定義為

AVx=［aVx（～x， 1）， aVx（～x， 2）， …， aVx（～x， K）］（15）

AVy=［aVy（～y， 1）， aVy（～y， 2）， …， aVy（～y， K）］（16）

于是， 在無噪聲的情況下， 令V表示的虛擬URA接收的單快拍數(shù)據(jù)矩陣表示為

Z=∑Kk=1pkaVx（～x， k）aTVy（～y， k）=AVxRSATVy（17）

式中： RS=diag（p1， …， pK）代表信號功率， 可知矩陣Z的大小為（2mxmax+1）×（2mymax+1）， 與虛擬陣列V的大小相對應。 為方便討論， 將其中的aVx（～x， k）aTVy（～y， k）展開觀察， 可以看出矩陣Z中各元素與V中元素有著一一對應的關系。

aVx（～x， k）aTVy（～y， k）=ejπ（－mxmax～x， k－mymax～y， k）…ejπ（－mxmax～x， k+mymax～y， k）

…

ejπ（mxmax～x， k－mymax～y， k）…ejπ（mxmax～x， k+mymax～y， k）（18）

結合式（13）～（14）， 設aV（～x， k， ～y， k）=aVx（～x， k）aVy（～y， k）， 該導向矢量按坐標順序對應V中各陣元， 所以有虛擬陣列的流形矩陣AV可表示為

AV=［aV（～x， 1， ～y， 1）， aV（～x， 2，～y， 2），? …，

aV（～x， K， ～y， K）］（19）

又已知D中元素在V中的索引， 對式（9）中向量zD在孔洞H處用0元素進行插值， 則虛擬陣列V接收的向量化單快拍數(shù)據(jù)可表示為

zV=AVp+σ2NιV （20）

式中： 噪聲項中ιV為式（9）中ιD的插值的結果。

然后， 將式（20）表示的單快拍列向量zV再次重排為（2mxmax+1）×（2mymax+1）的矩陣， 此時式（17）變?yōu)橛性肼暻闆r下的單快拍數(shù)據(jù)矩陣：

ZV=∑Kk=1pkaVx（～x， k）aTVy（～y， k）+σ2NΙV=

AVxRSATVy+σ2NΙV（21）

式中： σ2N為噪聲的方差； ΙV為向量ιV經矩陣重排后的結果， 是由0和1構成的二值矩陣。

現(xiàn)在， 考慮必要子陣失效時的情況， 根據(jù)2.1中的定義， 當P中有一β-必要子陣B失效時， 受損陣列為P-=P＼B， 新的差分共陣為D-≠D， 原虛擬陣列中新出現(xiàn)的孔洞坐標集合定義為Hf， 即Hf=D＼D-， 若原差分共陣不是無孔的， 則加上原本存在的孔洞H， 現(xiàn)虛擬陣列中所有孔洞集合記為H-=H∪Hf。 圖2展示了一個存在2個失效陣元的二維嵌套陣列與其受損后的差分虛擬陣列。 圖2 （a）中灰色圓點表示仍正常的物理陣元P-， 紅色圓點表示失效陣元B， 叉號以單位間隔表示陣元間的空白區(qū)域； 圖 2 （b）中灰色圓點表示有效的虛擬陣元D-， 紅色圓圈表示因陣元失效新增的孔洞Hf。

根據(jù)前文分析可知， 一方面， 對于虛擬陣列中原有孔洞H無對應的接收數(shù)據(jù)， 直接用0插值； 另一方面， 由于子陣B失效， 對應物理陣元接收數(shù)據(jù)為0， 使虛擬陣列中新增孔洞Hf。 為了定位矩陣ZV中需置零的全部元素， 且與集合H-中坐標對應， 定義映射矩陣G：

〈G〉i， j=1， （i－mxmax－1， j－mymax－1）∈D-0， （i－mxmax－1， j－mymax－1）∈H- （22）

式中： 〈·〉i， j表示矩陣的第i行第j列元素的索引。

于是在子陣B的陣元失效后， 式（21）的矩陣變?yōu)?/p>

Zf=ZVG（23）

式中： 表示Hadamard積， 即矩陣對應元素相乘。

為了在此條件下實現(xiàn)準確的二維DOA估計， 將采用基于解耦原子范數(shù)最小化理論的矩陣填充模型， 對矩陣Zf中的零元素進行數(shù)據(jù)填充， 重構出完整的虛擬陣列接收數(shù)據(jù)， 恢復虛擬陣列原有的連續(xù)性特征， 為后續(xù)應用子空間類DOA估計算法提供必要條件。

2.3 基于DANM的矩陣填充算法

根據(jù)文獻［15-16］中針對二維DOA估計問題給出的解耦原子范數(shù)最小化算法， 結合本文信號模型， 考慮以下問題。

回顧式（17）， 已知無噪聲的單快拍數(shù)據(jù)表示為Z=∑Kk=1pkaVx（～x， k）aTVy（～y， k）， 由此定義矩陣形式的原子集合為

A={aVx（～x）aTVy（～y）， ～x， ～y∈［－1， 1］}={Aatom（f）， f=（～x， ～y）∈［－1， 1］×［－1， 1］}（24）

式中： 矩陣Aatom（f）=aVx（～x）aTVy～y為集合中的一個原子， 則Z由原子集合A中K個原子組成。

在該原子集合上關于Z的原子范數(shù)表示為

ZA=infpk∈瘙綆{∑k|pk||Z=∑kpkAatom（f），

Aatom（f）∈A}（25）

式中：‖·‖A表示原子范數(shù)； inf表示下確界。

由于虛擬陣列D-中孔洞的存在導致其等價接收數(shù)據(jù)矩陣Zf中出現(xiàn)零元素， 所以需要利用式（22）定義的映射矩陣G， 將待優(yōu)化矩陣Z與矩陣Zf零元素相同位置的元素置零， 這樣可以避免較大的擬合誤差， 便于利用解耦原子范數(shù)最小化算法來恢復出滿秩的協(xié)方差數(shù)據(jù)矩陣Z^。

于是將最小化問題表述為

Z^=minZZA s.t. （ZG）－Zf2F≤ξ（26）

式中： ·F表示矩陣的Frobenius范數(shù)； ξ表示一個足夠小的正數(shù)， 作為閾值參數(shù)來約束擬合誤差。

為解決該原子范數(shù)最小化問題， 將其轉換為半正定規(guī)劃（Semidefinite Programming， SDP）問題。 在無噪聲情況下， 考慮數(shù)據(jù)矩陣Z∈瘙綇Lx×Ly， 其解耦原子范數(shù)如式（27）所示：

ZA=minux， uy12LxLy［tr（T （ux））+tr（T （uy））］s.t. T? （ux）ZZHT （uy）0（27）

式中： Lx=2mxmax+1， Ly=2mymax+1； T （ux）和T （uy）表示半正定Hermitian-Toeplitz矩陣， 以向量ux∈瘙綇Lx和uy∈瘙綇Ly分別作為兩矩陣的第1列進行構造； tr（·）表示矩陣的跡； 約束條件中“0”表示該矩陣是半正定的。

而在有噪聲條件下， SDP求解模型中將考慮式（26）中的約束項， 目標函數(shù)如式（28）所示：

minux， uy， Zμ2LxLy［tr（T （ux））+tr（T （uy））］+ZG－Zf2F? s.t.T （ux）ZZHT （uy）0（28）

式中： μ≥0為一加權因子， 作為正則化參數(shù)， 與噪聲方差和陣元總數(shù)相關， 其取值可參考文獻［24］。

式（28）中的SDP問題可使用CVX凸優(yōu)化工具箱求解， 優(yōu)化后的矩陣Z^不再有孔洞導致的零元素， 相當于無孔的虛擬URA接收的完整協(xié)方差矩陣數(shù)據(jù)。

為進一步削弱虛擬陣列模型帶來的信號相關性， 對矩陣Z^∈瘙綇Lx×Ly進行二維空間平滑操作［9， 25］， 可得到滿秩的數(shù)據(jù)矩陣RZ， 其大小為（Lx+1）（Ly+1）/4， 即將虛擬陣列的自由度縮減到了（mxmax+1）（mymax+1）， 但平滑操作可使二維DOA估計的精度顯著提高。 最后， 應用經典的二維MUSIC算法估計出全部K個信源的二維角度信息， 即（θk， φk）|k=1， 2， …， K， 使用二維譜峰搜索， 角度結果自動配對。

綜上所述， 所提出的基于解耦原子范數(shù)最小化的二維DOA估計算法的具體步驟如算法1所示。

算法1 陣元失效下基于DANM的二維DOA估計算法

輸入：

存在失效陣元的二維稀疏陣列的接收數(shù)據(jù)XP（t）， t=1， 2， …， T；

輸出： 二維DOA估計結果（θk， φk）， k=1， 2， …， K；

步驟1： 計算XP（t）的協(xié)方差矩陣R^P， 將其向量化為z；

步驟2：

根據(jù)式（10）求得稀疏陣列P-陣元失效前后的差分共陣D和D-， 并按式（11）設定均勻矩形陣列V， 根據(jù)式（22）設置映射矩陣G；

步驟3：

對應D的陣元坐標， 對向量z中元素進行去冗余操作并排序為zD， 按D在V中的索引將zD用0插值為zV；

步驟4：

將向量zV重排為與虛擬陣列V大小相當?shù)木仃嘮V， 其零元素位置由G確定， 對應陣元失效后的全部孔洞位置， 如式（23）所示；

步驟5：

構造SDP問題， 使用DANM算法進行矩陣填充， 利用CVX工具箱求解式（28）， 得到完整的數(shù)據(jù)矩陣Z^；

步驟6： 對Z^進行二維空間平滑操作得到協(xié)方差矩陣RZ；

步驟7： 對RZ使用經典的二維MUSIC算法進行二維DOA估計。

理論上， 陣列的自由度決定了其最大可估計信源數(shù)。 根據(jù)文獻［16］， 當使用Lx×Ly大小的URA接收數(shù)據(jù)， 應用DANM算法后直接使用范德蒙德分解估計二維角度信息時， 其最大可估計信源數(shù)受限于URA兩個方向上較短一邊的陣元數(shù)， 即min{Lx， Ly}， 這也是由于SDP的局限性。 而在本文方法中， 使用較少陣元的稀疏平面陣列， 由差分共陣和DANM算法構成Lx×Ly的虛擬URA， 再經過空間平滑操作后， 使自由度達到了（Lx+1）（Ly+1）/4。 另一方面， 陣元失效對虛擬陣列的連續(xù)性的破壞通常是較為嚴重的， 使用解耦原子范數(shù)最小化進行矩陣填充， 能夠避免虛擬陣列中連續(xù)URA部分縮小而浪費其他陣元數(shù)據(jù)， 有效維持了二維虛擬陣列的自由度。

3 仿真實驗

通過仿真實驗將本文所提出算法與現(xiàn)有方法進行比較， 驗證并評估算法在陣元失效情形下二維DOA估計性能。 為對比參數(shù)不同時對算法的影響， 使用均方根誤差（Root-Mean-Square Error， RMSE）來分析算法的DOA估計精度， 假設對K個入射信號進行Q次Monte-Carlo實驗， 則二維DOA估計的RMSE定義為

RMSE=1QK∑Qq=1∑Kk=1［（θ^k， q－θk）2+（φ^k， q－φk）2］（29）

式中： θk和φk分別為第k個入射信號的真實仰角和方位角； θ^k， q和φ^k， q為第k個入射信號在第q次實驗中的估計值。 特別地， 由于本仿真實驗中的入射信號較多， 在估計效果較差時常會出現(xiàn)多個估計值相近于同一個真實值附近的情況， 所以在計算RMSE時會選擇最準確的一個值， 而忽略另一些完全估計錯誤的角度值， 此時每次實驗中參與計算誤差的角度個數(shù)K′≤K。

另外， 在計算RMSE的同時記錄入射信號二維DOA的估計成功率ηSE， 其定義為符合要求的角度估計次數(shù)與估計總次數(shù)的比值：

ηSE=NsuccNall（30）

在本節(jié)的仿真實驗中， Nsucc均設為實驗中二維角度估計的平方根誤差≤2°的次數(shù)， Nall表示全部角度估計的總次數(shù)， 即在K個信號的Q次實驗中， Nall=QK。

下列實驗中， 本文方法基于圖2中的二維嵌套陣列， 假設在其42個物理陣元中， 位于（2， 15）和（11， 2）處的2個陣元失效， 原差分共陣D是無孔的， 陣元失效后D中產生38個孔洞， 經過矩陣填充， 使陣列達到最大自由度為192。 在不使用恢復孔洞數(shù)據(jù)的一般二維差分共陣（2D-DCA）方法中， 取圖2（b）虛擬陣列的陣元中從（－7， －11）到（7， 11）的中央連續(xù)矩形部分作為新的虛擬陣列， 同樣經二維空間平滑操作后進行DOA估計， 其陣列自由度可達96。 當基于大小為6×7， 同樣包含兩個失效陣元的物理均勻矩形陣列（URA）， 不經過數(shù)據(jù)恢復， 直接使用MUSIC算法進行二維DOA估計， 其自由度僅為42。

實驗1： 不同信源數(shù)對算法的估計性能的影響。

首先驗證本文方法的性能， 使用受損的二維嵌套陣列接收45個來自遠場的功率相等的獨立窄帶信號， 信號波長均滿足λ=2d， 所有信號振幅和相位均隨機。 信源的二維角度設置滿足文獻［26］中確保得到最優(yōu)解的最小角度間距。 信噪比（Signal to Noise Ratio， SNR）設置為20 dB， 快拍數(shù)為1 000， 二維譜峰搜索的掃描間隔設為（0.5°， 0.5°）。 結果如圖3所示， 所提算法可以估計全部45個信源， 超過了物理陣元的數(shù)量， 實現(xiàn)了穩(wěn)定的欠定DOA估計。

當信源數(shù)變化時， 將本文算法與另外兩種方法進行對比。 仿真實驗中信源數(shù)范圍為［1， 70］， 各信源在仰角-方位角空間按一定規(guī)律排布， 每10個角度為一行。 設置信噪比為0 dB， 快拍數(shù)為500， 每組均進行100次Monte-Carlo實驗， 計算三種方法對不同數(shù)量信源的估計成功率， 如圖4所示。 可以看出， 均勻陣列方法的最大可估計信源數(shù)被其陣元數(shù)所限定， 而本文方法利用稀疏陣列的差分共陣， 并使用DANM算法進行矩陣填充， 消除了陣元失效的影響， 保留了原虛擬陣列的最大自由度， 大大提高了可估計信源數(shù)， 且估計成功率更高， 表明該算法在面對數(shù)量大于物理陣元的信源時， 仍能維持一定的估計性能。

實驗2： 不同信噪比對算法的估計性能的影響。

使用的陣列形式和失效陣元同實驗1， 需估計信源個數(shù)為25， 角度固定不變， 快拍數(shù)均為500， 仿真中信噪比范圍為［－20， 20］dB， 對每組仿真實驗均進行100次Monte-Carlo實驗， 計算估計結果的估計成功率和均方根誤差， 其隨信噪比變化的曲線如圖5～6所示。 可見， 隨著信噪比的增大， 各算法的估計成功率不斷增大， 均方根誤差不斷變小， 但基于均勻陣列的算法結果在低信噪比時估計成功率一直較低， 其均方根誤差過高， 因此在圖6中只給出所提方法與2D-DCA方法的對比結果。 本文所提方法從-10 dB開始， 均方根誤差就已趨于平穩(wěn)， 估計成功率也均接近1， 且由于恢復了原虛擬陣列最大的自由度， 所得結果也優(yōu)于因孔洞縮減了自由度的一般差分共陣方法。 從實驗2可得， 所提方法在低信噪比的情況下有較強的魯棒性， 在多信源估計中有更高的準確性。

實驗3： 不同快拍數(shù)對算法的估計性能的影響。

使用的陣列形式和失效陣元同實驗1， 需估計信源個數(shù)仍為25， 其角度不變， 信噪比為0 dB， 仿真中快拍數(shù)最小為10， 最高至1 200， 對每組仿真實驗均進行100次Monte-Carlo實驗， 計算估計結果的估計成功率和均方根誤差， 其隨快拍數(shù)變化的曲線如圖7～8所示。 可以看出， 隨著快拍數(shù)的增大， 各算法的估計成功率同步增大， 均方根誤差均不斷減小。 本文所提方法在快拍數(shù)大于400時， 估計成功率均接近1， 當快拍數(shù)大于800時， 均方根誤差值開始趨于平穩(wěn)， 其性能明顯優(yōu)于其他方法， 精度也較高。 從實驗3可得， 當快拍數(shù)達到與陣列大小相適配的一定值后， 所提方法將充分發(fā)揮其虛擬陣列自由度大的優(yōu)勢， 在小快拍條件下表現(xiàn)出更高更準確的多信源角度估計性能。

實驗4： 不同重要度的子陣失效對本文所提算法的估計性能的影響。

根據(jù)前面的討論， 不同的子陣失效時對差分共陣的影響一般各不相同， 產生的新孔洞位置各異， 大小不一。 下面仍使用圖2中的二維陣列， 舉例關于該陣列的不同重要度κD的7個失效子陣B， 并與無失效陣元的情況進行對比， 驗證所提算法的有效性。? 設置信源數(shù)為45個， 信噪比為0 dB， 快拍數(shù)為500， Monte-Carlo實驗次數(shù)為100， 計算二維DOA估計的均方根誤差和估計成功率ηSE，? 結果如表1所示。 實驗結果顯示， 隨著失效子陣重要度的增加， 二維DOA估計的均方根誤差增大， 估計成功率同步降低， 說明所提算法對虛擬陣列的孔洞內插的效果受孔洞大小的影響， 孔洞越大， 越難以還原完整虛擬陣列的DOA估計效果。

4 結? 論

本文針對二維稀疏陣列中陣元失效下的DOA估計問題， 提出了一種有效的解決方案。 所提算法利用協(xié)方差矩陣數(shù)據(jù)的二階統(tǒng)計特性， 構建差分共陣模型， 在此基礎上引入解耦原子范數(shù)最小化算法進行矩陣填充恢復數(shù)據(jù)， 以實現(xiàn)對虛擬陣列中因物理陣元失效導致的孔洞進行內插。 此方法恢復了原始虛擬陣列的完整孔徑， 最大限度地提升了陣列的自由度， 確保高精度的DOA估計能力。 在成功恢復虛擬陣列數(shù)據(jù)后， 運用傳統(tǒng)的SS-MUSIC算法進行多信源的二維DOA估計， 實驗結果表明， 與已有的解決方案相比， 本文的方法能夠在相同陣元數(shù)量和陣元失效狀況下實現(xiàn)更多信源的成功估計， 并且在小快拍數(shù)、 低信噪比條件下展現(xiàn)出了更強的魯棒性。 盡管本文在該問題的探究上取得了一定效果， 但仍有若干挑戰(zhàn)亟待解決， 比如進一步提升欠定DOA估計精度， 優(yōu)化算法降低復雜度， 以及適應更廣泛的陣列構型， 這為后續(xù)深入研究提供了廣闊的探索空間。

參考文獻：

［1］ Mao Z H， Liu S H， Zhang Y D， et al. Joint DoA-Range Estimation Using Space-Frequency Virtual Difference Coarray［J］. IEEE Transactions on Signal Processing， 2022， 70： 2576-2592.

［2］ 曲明超， 司偉建， 袁雅芝. 基于不完全重合信號的單快拍DOA估計算法研究［J］. 通信學報， 2021， 42（12）： 88-95.

Qu Mingchao， Si Weijian， Yuan Yazhi. Research on Single Snapshot DOA Estimation Algorithm Based on Incompletely Overlapped Signal［J］. Journal on Communications， 2021， 42（12）： 88-95.（in Chinese）

［3］ 梁義魯， 司偉建， 李鎖蘭， 等. 基于混合陣列的DOA與極化信息聯(lián)合估計［J/OL］. 航空兵器， doi： 10.12132/ISSN.1673-5048.2023.0168.

Liang Yilu， Si Weijian， Li Suolan， et al. Joint Estimation of DOA and Polarization Information Based on Hybrid Arrays［J/OL］. Aero Weaponry， doi： 10.12132/ISSN.1673-5048.2023.0168. （in Chinese）

［4］ Schmidt R O. Multiple Emitter Location and Signal Parameter Estimation［J］. IEEE Transactions on Antennas and Propagation， 1986， 34（3）： 276-280.

［5］ Aboumahmoud I， Muqaibel A， Alhassoun M， et al. A Review of Sparse Sensor Arrays for Two-Dimensional Direction-of-Arrival Estimation［J］. IEEE Access， 2021， 9： 92999-93017.

［6］ Liu C L， Vaidyanathan P P. Cramér-Rao Bounds for Coprime and other Sparse Arrays， which Find more Sources than Sensors［J］. Digital Signal Processing， 2017， 61： 43-61.

［7］ Pal P， Vaidyanathan P P. Nested Arrays： A Novel Approach to Array Processing with Enhanced Degrees of Freedom［J］. IEEE Transactions on Signal Processing， 2010， 58（8）： 4167-4181.

［8］ Pal P， Vaidyanathan P P. Nested Arrays in Two Dimensions， Part Ⅰ： Geometrical Considerations［J］. IEEE Transactions on Signal Processing， 2012， 60（9）： 4694-4705.

［9］ Pal P， Vaidyanathan P P. Nested Arrays in Two Dimensions， Part Ⅱ： Application in Two Dimensional Array Processing［J］. IEEE Transactions on Signal Processing， 2012， 60（9）： 4706-4718.

［10］ Liu C L， Vaidyanathan P P. Remarks on the Spatial Smoothing Step in Coarray MUSIC［J］. IEEE Signal Processing Letters， 2015， 22（9）： 1438-1442.

［11］ Liu C L， Vaidyanathan P P， Pal P. Coprime Coarray Interpolation for DOA Estimation via Nuclear Norm Minimization［C］∥2016 IEEE International Symposium on Circuits and Systems （ISCAS）， 2016.

［12］ Zhou C W， Gu Y J， Shi Z G， et al. Off-Grid Direction-of-Arrival Estimation Using Coprime Array Interpolation［J］. IEEE Signal Processing Letters， 2018， 25（11）： 1710-1714.

［13］ Qiao H， Pal P. Unified Analysis of Co-Array Interpolation for Direction-of-Arrival Estimation［C］∥2017 IEEE InternationalConference on Acoustics， Speech and Signal Processing （ICASSP）， 2017.

［14］ Zhou C W， Gu Y J， Fan X， et al. Direction-of-Arrival Estimation for Coprime Array via Virtual Array Interpolation［J］. IEEE Transactions on Signal Processing， 2018， 66（22）： 5956-5971.

［15］ Tian Z， Zhang Z， Wang Y. Low-Complexity Optimization for Two-Dimensional Direction-of-Arrival Estimation via Decoupled Atomic Norm Minimization［C］∥2017 IEEE International Conference on Acoustics， Speech and Signal Processing （ICASSP）， 2017.

［16］ Zhang Z， Wang Y， Tian Z. Efficient Two-Dimensional Line Spectrum Estimation Based on Decoupled Atomic Norm Minimization［J］. Signal Processing， 2019， 163： 95-106.

［17］ 彭加強， 鄭桂妹. 基于解耦原子范數(shù)最小化的二維DOA估計［J］. 空軍工程大學學報： 自然科學版， 2021， 22（1）： 55-61.

Peng Jiaqiang， Zheng Guimei. A Two-Dimensional DOA Estimation Based on Decoupled Atomic Norm Minimization［J］. Journal of Air Force Engineering University： Natural Science Edition， 2021， 22（1）： 55-61.（in Chinese）

［18］ Lu A H， Guo Y， Li N， et al. Efficient Gridless 2-D Direction-of-Arrival Estimation for Coprime Array Based on Decoupled Atomic Norm Minimization［J］. IEEE Access， 2020， 8： 57786-57795.

［19］ Zhu C L， Wang W Q， Chen H， et al. Impaired Sensor Diagnosis， Beamforming， and DOA Estimation with Difference Co-Array Processing［J］. IEEE Sensors Journal， 2015， 15（7）： 3773-3780.

［20］ 孫兵， 吳晨曦， 阮懷林， 等. 陣元失效條件下的高精度DOA估計方法［J］. 電子學報， 2020， 48（9）： 1688-1694.

Sun Bing， Wu Chenxi， Ruan Huailin， et al. High Accuracy DOA Estimation Method with Array Sensor Failure［J］. Acta Electronica Sinica， 2020， 48（9）： 1688-1694.（in Chinese）

［21］ Chen J L， Zhang C， Fu S T， et al. Robust Reweighted l2，1-Norm Based Approach for DOA Estimation in MIMO Radar under Array Sensor Failures［J］. IEEE Sensors Journal， 2021， 21（24）： 27858-27867.

［22］ Liu C L， Vaidyanathan P P. Robustness of Difference Coarrays of Sparse Arrays to Sensor Failures-Part I： A Theory Motivated by Coarray MUSIC［J］. IEEE Transactions on Signal Processing， 2019， 67（12）： 3213-3226.

［23］ Liu C L， Vaidyanathan P P. Hourglass Arrays and other Novel 2D Sparse Arrays with Reduced Mutual Coupling［J］. IEEE Transactions on Signal Processing， 2017， 65（13）： 3369-3383.

［24］ Bhaskar B N， Tang G G， Recht B. Atomic Norm Denoising with Applications to Line Spectral Estimation［J］. IEEE Transactions on Signal Processing， 2013， 61（23）： 5987-5999.

［25］ Chen Y M. On Spatial Smoothing for Two-Dimensional Direction-of-Arrival Estimation of Coherent Signals［J］. IEEE Transactions on Signal Processing， 1997， 45（7）： 1689-1696.

［26］ Chi Y J， Chen Y X. Compressive Two-Dimensional Harmonic Retrieval via Atomic Norm Minimization［J］. IEEE Transactions on Signal Processing， 2015， 63（4）： 1030-1042.

Two-Dimensional DOA Estimation Algorithm for

Sparse Arrays under Sensor Failures

Si Weijian1， 2， Ma Wanyu1， 2*， Yao Lu3， Qu Mingchao1， 2， Liang Yilu1， 2

（1. College of Information and Communication Engineering， Harbin Engineering University， Harbin 150001， China;

2. Key Laboratory of Advanced Marine Communication and Information Technology， Ministry of Industry and

Information Technology， Harbin Engineering University， Harbin 150001， China;

3. Jiangsu Press and Publishing School， Nanjing 210012， China）

Abstract： To address the problem of the destruction of virtual array continuity and the degradation of degree of freedom due to missing data in two-dimensional sparse array under the conditions of sensor failures， a two-dimensional DOA estimation algorithm is proposed. Firstly， the virtual array is constructed based on the two-dimensional difference coarray， and then the covariance matrix data is recovered in the form of matrix completion by using decoupled atomic norm minimization to realize the virtual array interpolation. Finally， the SS-MUSIC algorithm is used for the two-dimensional DOA estimation of multiple sources. The proposed method compensates for the effects caused by the failure of physical sensors. It recovers the complete aperture characteristics of the original virtual array and maintains the degree of freedom of the virtual array， which ensures a higher-precision two-dimensional DOA estimation performance. Simulation results demonstrate that under the same number of physical elements and sensor failure， the proposed algorithm can effectively estimate more sources compared with the existing methods， and exhibits higher robustness under the conditions of a small number of snapshots and low signal-to-noise ratio. This approach maximally retains and utilizes the degree of freedom advantage of sparse arrays in the two-dimensional DOA estimation.

Key words：? two-dimensional DOA estimation; sparse array; difference coarray; sensor failure; decoupled atomic norm minimization; matrix completion