阵元失效下稀疏阵列的二维DOA估计算法

司伟建 马万禹 姚璐 曲明超 梁义鲁

摘 要：本文针对二维稀疏阵列在阵元失效条件下， 因数据缺失导致虚拟阵列连续性被破坏及自由度下降的问题， 提出了一种二维DOA估计算法。 首先基于二维差分共阵构建虚拟阵列， 然后利用解耦原子范数最小化理論， 以矩阵填充的形式恢复协方差矩阵数据， 实现对虚拟阵列中丢失虚拟阵元的内插， 最后采用SS-MUSIC算法进行多信源的二维DOA估计。 所提方法弥补了物理阵元失效所造成的影响， 恢复了原始虚拟阵列的完整孔径特性， 保持了虚拟阵列的自由度， 从而确保了较高精度的二维DOA估计性能。 仿真实验结果表明， 在相同阵元数量及阵元失效情况下， 本文提出的算法相比已有方法能有效地估计更多信源， 并在小快拍数和低信噪比条件下表现出更高的稳健性， 最大限度地保留并利用了稀疏阵列在二维DOA估计中的自由度优势。

关键词：  二维DOA估计； 稀疏阵列； 差分共阵； 阵元失效； 解耦原子范数最小化； 矩阵填充

中图分类号： TJ765； TN911.7

文献标识码： A

文章编号：  1673-5048（2024）02-0114-09

DOI： 10.12132/ISSN.1673-5048.2024.0042

0 引  言

波达方向（Direction of Arrival， DOA）估计是阵列信号处理领域内的一个极其重要的分支， 长期以来受到学者们的广泛关注。 在过去的几十年间， DOA估计理论不断发展， 从雷达探测、 无源定位、 无线通信、 卫星导航， 乃至物联网环境中的智能感知， DOA估计技术在众多领域均展现出广泛而深远的应用价值［1-3］。 以多重信号分类（Multiple Signal Classification， MUSIC）算法［4］为代表的子空间分解类算法凭借高精度和超分辨的DOA估计能力为人所熟知， 并获得了广泛的应用与实践验证。 然而这类算法的信号接收模型中通常依赖于均匀阵列， 其理论上可估计的最大信源数目受限于阵元的个数。

因此， 在日益复杂的空间电磁环境下， 针对信号密集场景的精确测向问题， 基于稀疏阵列的欠定DOA估计成为新的研究热点［5］， 稀疏阵列突破了阵元间距必须为入射信号半波长整数倍的限制， 通过设计阵列构型以构建大孔径虚拟阵列， 从而显著提升了阵列自由度（Degree of Freedom， DOF）， 进而能够实现超过实际物理阵元数量的DOA估计［6］。 Pal等学者开创性地提出嵌套阵列的概念， 并基于此发明出一种欠定DOA估计算法［7］。 该算法首先通过对原始接收数据进行协方差矩阵的向量化处理， 形成与虚拟阵列相对应的单快拍数据， 随后执行空间平滑（Spatial Smoothing， SS）操作获得满秩的协方差矩阵， 使用MUSIC算法完成角度估计。 随后Pal等人进一步拓展， 又提出了适用于二维空间的平面嵌套阵列设计及相应DOA估计算法［8-9］。 尽管SS-MUSIC算法通常可以应用于大多数类型的稀疏阵列DOA估计中， 其有效性却高度依赖于虚拟阵列结构的连续性［10］， 一旦虚拟阵列呈现出不连续性， 即其间存在孔洞， 那么阵列的连续自由度扩展将受到显著限制。

通过使用协方差矩阵重构技术， 能够对虚拟孔洞位置的空缺数据信息进行近似化填充， 从而扩大虚拟阵列孔径。 一系列基于互质线阵内插的DOA估计算法被提出， 文献［11］通过核范数最小化构建凸优化模型， 填补虚拟阵列孔洞并提升DOA估计精度， 但核范数方法可能引入一定近似误差； 文献［12］则进一步将核范数最小化问题转化为迹最小化以求解， 实现离网格DOA估计； 文献［13］对虚拟阵列内插算法进行了完整分析， 优化模型限制了此过程中的误差积累。 近年来， 原子范数理论因其对信号空域稀疏特性的有效利用， 在DOA估计领域取得了显著进展， 原子范数最小化（Atomic Norm Minimization， ANM）作为一种无网格方法， 避免了一般压缩感知方法中的基失配问题。 文献［14］将ANM算法应用于互质线阵的虚拟阵列内插， 实现了高精度的无网格DOA估计。 文献［15-16］则针对二维DOA估计提出解耦原子范数最小化（Decoupled ANM， DANM）算法， 基于二维阵列接收数据模型， 将二维DOA估计问题转化为两个独立的一维DOA估计； 文献［17］中对DANM算法做出改进， 使之适用于多快拍情形下的二维DOA估计。 但这两种方法都局限于均匀矩形阵列， 且最大可估计信源数受限于较短一边的阵元数。 文献［18］基于互质平面阵列， 采用DANM算法进行二维虚拟阵列内插， 可实现欠定的二维DOA估计。

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

司伟建， 等： 阵元失效下稀疏阵列的二维DOA估计算法

在实际应用环境中， 由于硬件老化、 环境干扰等多种难以预测的因素导致性能退化乃至完全失效。 一旦阵元随机发生故障， 虚拟阵列的连续性很可能会因此受到显著破坏， 形成位置和大小不定的孔洞， 进而大幅降低连续自由度， 并对DOA估计精度产生严重影响［19］。 面对阵元失效情况下的DOA估计难题， 类似地， 可通过内插方法填补阵元失效后虚拟阵列中的孔洞， 近似恢复丢失虚拟阵元的数据信息， 以维持虚拟阵列的连续性， 从而最大限度地利用阵列原始的完整孔径优势。 文献［20］采用Toeplitz矩阵重构技术恢复协方差矩阵； 而文献［21］则基于重加权l2， 1范数最小化算法恢复了完整的阵列数据矩阵， 实现了一维DOA的有效估计。 然而这两种方法在将算法扩展到二维场景时仍面临较大困难。

鉴于当前研究现状， 阵元失效问题的研究重心主要集中在一维阵列上， 而在二维稀疏阵列DOA估计领域， 关于阵元失效下的虚拟阵列孔洞内插问题， 现有的研究尚不充分。 因此， 本文提出一种阵元失效下的二维稀疏阵列DOA估计算法， 该算法利用二维差分共阵形成虚拟阵列， 通过建立基于DANM理论的矩阵填充问题， 恢复协方差矩阵数据， 对阵元失效导致的孔洞内插， 以恢复原虚拟阵列的完整自由度， 并结合SS-MUSIC算法进行二维DOA估计。 对比传统算法， 本文方法能够避免舍弃部分虚拟阵元而造成的自由度损失， 有效维持了稀疏阵列二维DOA估计的优势性能。

1 信号模型

考虑空间中K个窄带远场独立信号入射到XOY平面中一个由M个标量天线构成的稀疏平面阵列上， 如图 1所示。 在平面直角坐标系中， 阵列中任一阵元的位置都可表示为（xid， yid）的形式， 其中： xi， yi∈瘙綄， i=1， 2， …， M； d=λ/2， 即半波长。 令ni=（xi， yi）∈瘙綄2， 代表单个阵元位置， 构成该稀疏平面阵列的所有物理阵元的二维位置整数集合记为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是虚数单位， k=1， 2， …， K。

假设噪声为加性高斯白噪声， 且噪声与信号相互统计独立， 则在t时刻时， 稀疏平面阵列的信号接收数据模型为

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个快拍数据时， 式（2）转化为

XP=APS+N （3）

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

计算阵列接收数据的协方差矩阵RP=E［XPXHP］， 其中E［·］表示数学期望， 将式（3）代入可得：

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

式中： RS为信号的自协方差矩阵， 其对角线元素表示信号功率； σ2N为输入噪声N的方差， I为单位阵。

由于采樣数量有限， 通常采用最大似然估计方法对式（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中元素个数一定小于M2， 但虚拟阵列的阵元数量较原物理阵列的大幅提升将使二维欠定DOA估计得以实现。

对应地， 式（7）向量中的元素也会有重复， 于是根据与D中各阵元对应关系， 对式（6）中向量z进行去冗余操作， 即将其中对应同一虚拟阵元的元素取平均值再按对应坐标大小重新排序［10］， 则对应D的单快拍虚拟接收数据为

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

式中： AD为虚拟阵列D接收数据的流形矩阵， 其大小为|D|×K（|·|表示集合的势， 即集合中元素的个数）， ιD表示式（6）中向量化单位阵vec（I）对应D的重排序。

2 阵元失效下二维DOA估计

2.1 二维差分共阵

根据式（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个物理阵元形成的差分共阵的阵元数最大可达|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）

但大多数二维稀疏阵列的D并不是一个完整的均匀矩形阵列（Uniform Rectangular Array， URA）， 所以将属于V但不属于D的位置坐标集合称为“孔洞”， 记为H=V＼D。 而当H=时， 有D=V， 则称D为无孔的虚拟阵列， 在这种情形下， 无需在虚拟阵列中央划取一个连续URA而舍弃掉外围的虚拟阵元， 能够最大化利用D的自由度。

当阵列中个别阵元出现故障以至于失效时， 阵元位置集合P中对应元素缺失， 其差分共阵D有可能发生改变， 破坏虚拟阵列中的连续结构， 使其连续自由度大幅缩减。 对于β个阵元同时失效的情形， 可以根据是否对D产生影响来判断失效阵元组成的子阵列在原阵列中的必要性。 于是定义： 当从P中移除某个阵元数为β的子阵B后， 使形成的差分共阵D发生改变时， 即若有P-=P＼B， 使得D-≠D， 称B对阵列P是β-必要的［22］； 否则， B是β-非必要的。 这里， P-表示移除子阵后的受损阵列， D-为P-的差分共阵。

考虑到不同子阵失效对差分共阵D的影响一般是不同的， 以无孔的D为例， P中的必要阵元失效将导致原本均匀矩形虚拟阵列中产生新的孔洞Hf， 孔洞越大， 说明对D的结构破坏越严重， 进而使得在利用矩阵填充理论和最优化算法重构数据时引入的误差增大， 这种误差积累将直接影响二维DOA估计的精度与可靠性。 因此， 定义子阵B在P中的重要度如下：

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

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

2.2 阵元失效下虚拟阵列接收数据

首先考虑阵列P中非必要子阵失效的情况， 此时其差分共阵D不发生改变。 使用式（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接收的单快拍数据矩阵表示为

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接收的向量化单快拍数据可表示为

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

式中： 噪声项中ιV为式（9）中ιD的插值的结果。

然后， 将式（20）表示的单快拍列向量zV再次重排为（2mxmax+1）×（2mymax+1）的矩阵， 此时式（17）变为有噪声情况下的单快拍数据矩阵：

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

AVxRSATVy+σ2NΙV（21）

式中： σ2N为噪声的方差； ΙV为向量ιV经矩阵重排后的结果， 是由0和1构成的二值矩阵。

现在， 考虑必要子阵失效时的情况， 根据2.1中的定义， 当P中有一β-必要子阵B失效时， 受损阵列为P-=P＼B， 新的差分共阵为D-≠D， 原虚拟阵列中新出现的孔洞坐标集合定义为Hf， 即Hf=D＼D-， 若原差分共阵不是无孔的， 则加上原本存在的孔洞H， 现虚拟阵列中所有孔洞集合记为H-=H∪Hf。 图2展示了一个存在2个失效阵元的二维嵌套阵列与其受损后的差分虚拟阵列。 图2 （a）中灰色圆点表示仍正常的物理阵元P-， 红色圆点表示失效阵元B， 叉号以单位间隔表示阵元间的空白区域； 图 2 （b）中灰色圆点表示有效的虚拟阵元D-， 红色圆圈表示因阵元失效新增的孔洞Hf。

根据前文分析可知， 一方面， 对于虚拟阵列中原有孔洞H无对应的接收数据， 直接用0插值； 另一方面， 由于子阵B失效， 对应物理阵元接收数据为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）的矩阵变为

Zf=ZVG（23）

式中： 表示Hadamard积， 即矩阵对应元素相乘。

为了在此条件下实现准确的二维DOA估计， 将采用基于解耦原子范数最小化理论的矩阵填充模型， 对矩阵Zf中的零元素进行数据填充， 重构出完整的虚拟阵列接收数据， 恢复虚拟阵列原有的连续性特征， 为后续应用子空间类DOA估计算法提供必要条件。

2.3 基于DANM的矩阵填充算法

根据文献［15-16］中针对二维DOA估计问题给出的解耦原子范数最小化算法， 结合本文信号模型， 考虑以下问题。

回顾式（17）， 已知无噪声的单快拍数據表示为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的原子范数表示为

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

Aatom（f）∈A}（25）

式中：‖·‖A表示原子范数； inf表示下确界。

由于虚拟阵列D-中孔洞的存在导致其等价接收数据矩阵Zf中出现零元素， 所以需要利用式（22）定义的映射矩阵G， 将待优化矩阵Z与矩阵Zf零元素相同位置的元素置零， 这样可以避免较大的拟合误差， 便于利用解耦原子范数最小化算法来恢复出满秩的协方差数据矩阵Z^。

于是将最小化问题表述为

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

式中： ·F表示矩阵的Frobenius范数； ξ表示一个足够小的正数， 作为阈值参数来约束拟合误差。

为解决该原子范数最小化问题， 将其转换为半正定规划（Semidefinite Programming， SDP）问题。 在无噪声情况下， 考虑数据矩阵Z∈瘙綇Lx×Ly， 其解耦原子范数如式（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）中的约束项， 目标函数如式（28）所示：

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

式中： μ≥0为一加权因子， 作为正则化参数， 与噪声方差和阵元总数相关， 其取值可参考文献［24］。

式（28）中的SDP问题可使用CVX凸优化工具箱求解， 优化后的矩阵Z^不再有孔洞导致的零元素， 相当于无孔的虚拟URA接收的完整协方差矩阵数据。

为进一步削弱虚拟阵列模型带来的信号相关性， 对矩阵Z^∈瘙綇Lx×Ly进行二维空间平滑操作［9， 25］， 可得到满秩的数据矩阵RZ， 其大小为（Lx+1）（Ly+1）/4， 即将虚拟阵列的自由度缩减到了（mxmax+1）（mymax+1）， 但平滑操作可使二维DOA估计的精度显著提高。 最后， 应用经典的二维MUSIC算法估计出全部K个信源的二维角度信息， 即（θk， φk）|k=1， 2， …， K， 使用二维谱峰搜索， 角度结果自动配对。

综上所述， 所提出的基于解耦原子范数最小化的二维DOA估计算法的具体步骤如算法1所示。

算法1 阵元失效下基于DANM的二维DOA估计算法

输入：

存在失效阵元的二维稀疏阵列的接收数据XP（t）， t=1， 2， …， T；

输出： 二维DOA估计结果（θk， φk）， k=1， 2， …， K；

步骤1： 计算XP（t）的协方差矩阵R^P， 将其向量化为z；

步骤2：

根据式（10）求得稀疏阵列P-阵元失效前后的差分共阵D和D-， 并按式（11）设定均匀矩形阵列V， 根据式（22）设置映射矩阵G；

步骤3：

对应D的阵元坐标， 对向量z中元素进行去冗余操作并排序为zD， 按D在V中的索引将zD用0插值为zV；

步骤4：

将向量zV重排为与虚拟阵列V大小相当的矩阵ZV， 其零元素位置由G确定， 对应阵元失效后的全部孔洞位置， 如式（23）所示；

步骤5：

构造SDP问题， 使用DANM算法进行矩阵填充， 利用CVX工具箱求解式（28）， 得到完整的数据矩阵Z^；

步骤6： 对Z^进行二维空间平滑操作得到协方差矩阵RZ；

步骤7： 对RZ使用经典的二维MUSIC算法进行二维DOA估计。

理论上， 阵列的自由度决定了其最大可估计信源数。 根据文献［16］， 当使用Lx×Ly大小的URA接收数据， 应用DANM算法后直接使用范德蒙德分解估计二维角度信息时， 其最大可估计信源数受限于URA两个方向上较短一边的阵元数， 即min{Lx， Ly}， 这也是由于SDP的局限性。 而在本文方法中， 使用较少阵元的稀疏平面阵列， 由差分共阵和DANM算法构成Lx×Ly的虚拟URA， 再经过空间平滑操作后， 使自由度达到了（Lx+1）（Ly+1）/4。 另一方面， 阵元失效对虚拟阵列的连续性的破坏通常是较为严重的， 使用解耦原子范数最小化进行矩阵填充， 能够避免虚拟阵列中连续URA部分缩小而浪费其他阵元数据， 有效维持了二维虚拟阵列的自由度。

3 仿真实验

通过仿真实验将本文所提出算法与现有方法进行比较， 验证并评估算法在阵元失效情形下二维DOA估计性能。 为对比参数不同时对算法的影响， 使用均方根误差（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次实验中的估计值。 特别地， 由于本仿真实验中的入射信号较多， 在估计效果较差时常会出现多个估计值相近于同一个真实值附近的情况， 所以在计算RMSE时会选择最准确的一个值， 而忽略另一些完全估计错误的角度值， 此时每次实验中参与计算误差的角度个数K′≤K。

另外， 在计算RMSE的同时记录入射信号二维DOA的估计成功率ηSE， 其定义为符合要求的角度估计次数与估计总次数的比值：

ηSE=NsuccNall（30）

在本节的仿真实验中， Nsucc均设为实验中二维角度估计的平方根误差≤2°的次数， Nall表示全部角度估计的总次数， 即在K个信号的Q次实验中， Nall=QK。

下列实验中， 本文方法基于图2中的二维嵌套阵列， 假设在其42个物理阵元中， 位于（2， 15）和（11， 2）处的2个阵元失效， 原差分共阵D是无孔的， 阵元失效后D中产生38个孔洞， 经过矩阵填充， 使阵列达到最大自由度为192。 在不使用恢复孔洞数据的一般二维差分共阵（2D-DCA）方法中， 取图2（b）虚拟阵列的阵元中从（－7， －11）到（7， 11）的中央连续矩形部分作为新的虚拟阵列， 同样经二维空间平滑操作后进行DOA估计， 其阵列自由度可达96。 当基于大小为6×7， 同样包含两个失效阵元的物理均匀矩形阵列（URA）， 不经过数据恢复， 直接使用MUSIC算法进行二维DOA估计， 其自由度仅为42。

实验1： 不同信源数对算法的估计性能的影响。

首先验证本文方法的性能， 使用受损的二维嵌套阵列接收45个来自远场的功率相等的独立窄带信号， 信号波长均满足λ=2d， 所有信号振幅和相位均随机。 信源的二维角度设置满足文献［26］中确保得到最优解的最小角度间距。 信噪比（Signal to Noise Ratio， SNR）设置为20 dB， 快拍数为1 000， 二维谱峰搜索的扫描间隔设为（0.5°， 0.5°）。 结果如图3所示， 所提算法可以估計全部45个信源， 超过了物理阵元的数量， 实现了稳定的欠定DOA估计。

当信源数变化时， 将本文算法与另外两种方法进行对比。 仿真实验中信源数范围为［1， 70］， 各信源在仰角-方位角空间按一定规律排布， 每10个角度为一行。 设置信噪比为0 dB， 快拍数为500， 每组均进行100次Monte-Carlo实验， 计算三种方法对不同数量信源的估计成功率， 如图4所示。 可以看出， 均匀阵列方法的最大可估计信源数被其阵元数所限定， 而本文方法利用稀疏阵列的差分共阵， 并使用DANM算法进行矩阵填充， 消除了阵元失效的影响， 保留了原虚拟阵列的最大自由度， 大大提高了可估计信源数， 且估计成功率更高， 表明该算法在面对数量大于物理阵元的信源时， 仍能维持一定的估计性能。

实验2： 不同信噪比对算法的估计性能的影响。

使用的阵列形式和失效阵元同实验1， 需估计信源个数为25， 角度固定不变， 快拍数均为500， 仿真中信噪比范围为［－20， 20］dB， 对每组仿真实验均进行100次Monte-Carlo实验， 计算估计结果的估计成功率和均方根误差， 其随信噪比变化的曲线如图5～6所示。 可见， 随着信噪比的增大， 各算法的估计成功率不断增大， 均方根误差不断变小， 但基于均匀阵列的算法结果在低信噪比时估计成功率一直较低， 其均方根误差过高， 因此在图6中只给出所提方法与2D-DCA方法的对比结果。 本文所提方法从-10 dB开始， 均方根误差就已趋于平稳， 估计成功率也均接近1， 且由于恢复了原虚拟阵列最大的自由度， 所得结果也优于因孔洞缩减了自由度的一般差分共阵方法。 从实验2可得， 所提方法在低信噪比的情况下有较强的鲁棒性， 在多信源估计中有更高的准确性。

实验3： 不同快拍数对算法的估计性能的影响。

使用的阵列形式和失效阵元同实验1， 需估计信源个数仍为25， 其角度不变， 信噪比为0 dB， 仿真中快拍数最小为10， 最高至1 200， 对每组仿真实验均进行100次Monte-Carlo实验， 计算估计结果的估计成功率和均方根误差， 其随快拍数变化的曲线如图7～8所示。 可以看出， 随着快拍数的增大， 各算法的估计成功率同步增大， 均方根误差均不断减小。 本文所提方法在快拍数大于400时， 估计成功率均接近1， 当快拍数大于800时， 均方根误差值开始趋于平稳， 其性能明显优于其他方法， 精度也较高。 从实验3可得， 当快拍数达到与阵列大小相适配的一定值后， 所提方法将充分发挥其虚拟阵列自由度大的优势， 在小快拍条件下表现出更高更准确的多信源角度估计性能。

實验4： 不同重要度的子阵失效对本文所提算法的估计性能的影响。

根据前面的讨论， 不同的子阵失效时对差分共阵的影响一般各不相同， 产生的新孔洞位置各异， 大小不一。 下面仍使用图2中的二维阵列， 举例关于该阵列的不同重要度κD的7个失效子阵B， 并与无失效阵元的情况进行对比， 验证所提算法的有效性。  设置信源数为45个， 信噪比为0 dB， 快拍数为500， Monte-Carlo实验次数为100， 计算二维DOA估计的均方根误差和估计成功率ηSE，  结果如表1所示。 实验结果显示， 随着失效子阵重要度的增加， 二维DOA估计的均方根误差增大， 估计成功率同步降低， 说明所提算法对虚拟阵列的孔洞内插的效果受孔洞大小的影响， 孔洞越大， 越难以还原完整虚拟阵列的DOA估计效果。

4 结  论

本文针对二维稀疏阵列中阵元失效下的DOA估计问题， 提出了一种有效的解决方案。 所提算法利用协方差矩阵数据的二阶统计特性， 构建差分共阵模型， 在此基础上引入解耦原子范数最小化算法进行矩阵填充恢复数据， 以实现对虚拟阵列中因物理阵元失效导致的孔洞进行内插。 此方法恢复了原始虚拟阵列的完整孔径， 最大限度地提升了阵列的自由度， 确保高精度的DOA估计能力。 在成功恢复虚拟阵列数据后， 运用传统的SS-MUSIC算法进行多信源的二维DOA估计， 实验结果表明， 与已有的解决方案相比， 本文的方法能够在相同阵元数量和阵元失效状况下实现更多信源的成功估计， 并且在小快拍数、 低信噪比条件下展现出了更强的鲁棒性。 尽管本文在该问题的探究上取得了一定效果， 但仍有若干挑战亟待解决， 比如进一步提升欠定DOA估计精度， 优化算法降低复杂度， 以及适应更广泛的阵列构型， 这为后续深入研究提供了广阔的探索空间。

参考文献：

［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与极化信息联合估计［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］ 彭加強， 郑桂妹. 基于解耦原子范数最小化的二维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