基于协方差矩阵重构的互质阵列DOA估计方法

张 昊，吴晨曦

(1.安徽城市管理职业学院，安徽 合肥 230601；2.安徽省制约重点实验室，安徽 合肥 230037)

0 引言

波达方向估计(Direction of Arrival Estimation，DOA)是阵列信号处理领域研究的重要内容之一，在雷达、通信、电子对抗等诸多领域都有着广泛的应用。近几十年来，一系列高分辨方法相继被提出，但现有大部分方法仅适用于超定条件(阵元数多于信号数)， 随着现代电磁环境的日益复杂，欠定DOA估计问题越来越常见，如何利用较少的阵元数对更多的信号进行估计是目前研究的难点问题。

互质阵列[1]的提出为解决欠定DOA估计问题提供了一个新的思路。由于其具有易于构造、阵列扩展性好、物理阵元和虚拟阵元具有解析表达式等优点受到国内外学者的广泛关注[2-8]。现有成果可以分为两大类：一是基于子空间类方法[2-4]，该类方法的不足在于只能够利用差联合阵列中连续虚拟阵元部分而舍弃了非连续部分的虚拟阵元，导致阵列的虚拟孔径存在一定的损失同时估计精度受到扫描网格的影响；二是基于稀疏重构的方法[5-8]，该类方法虽然克服子空间方法的不足，但其前提条件是所有目标准确地位于预定的字典网格上，实际环境中并不一定满足即存在模型失配问题，此时估计性能会严重下降。虽然字典网格越密集，DOA信息与字典网格之间的误差越小，但计算量会随之增加同时过于密集的网格会造成字典原子之间的相关性增强，影响重构的性能。

阵列自由度与其相对应的差联合阵列中的虚拟阵元有关，而互质阵列的差联合阵列中由于存在空洞部分，若能够将空间部分进行有效填充，进而得到一个虚拟阵元更多、孔径更大的均匀阵列，从而提高可估计信号数和估计精度[9-11]。文献[9]通过阵列运动的方式来获得与空缺项对应的测量值，但该方法对运动状态有严格的要求，在实际应用中受限。文献[10]通过引入额外的工作频率来实现对空洞的填充，然而额外频率的使用会增大系统的工作带宽且会引起信号源特性起伏，在一定程度上增加了系统复杂性。文献[11]提出一种基于稀疏重构的外推技术，但该方法的估计性能同样受到模型失配的影响[12]。针对上述问题，本文提出基于协方差矩阵重构的互质阵列DOA估计方法。

1 信号模型

互质阵列由两个不同的均匀子阵组成，如图1所示，其中子阵1包含有N个物理阵元，分别位于{Mnd,0≤n≤N-1}，子阵2包含有2M-1个物理阵元，分别位于{Nmd,1≤m≤2M-1}，M,N是两个互质的整数，d=λ/2，λ为入射信号波长。

图1 互质阵列结构图Fig.1 Structure of coprime array

假设K个远场窄带平面波信号分别以[θ1,θ2,…,θK]入射到互质阵列上，t时刻阵列接收数据可表示为：

X(t)=AS(t)+n(t)

(1)

式(1)中，X(t)=[x1(t),x2(t),…,x2M+N-1(t)]T为阵列接收数据向量，S(t)=[s1(t),…,sK(t)]T为空间信号向量，n(t)=[n1(t),n2(t),…,n2M+N-1(t)]T为噪声向量，A=[a(θ1),…,a(θK)]为(2M+N-1)×K维阵列流型矩阵，其中a(θi)为入射角θi的导向矢量。

假设各入射信号之间互不相关，且噪声为高斯白噪声且与信号相互独立。则阵列接收数据的协方差矩阵可表示为：

RX=E[X(t)XH(t)]=ARsAH+n

(2)

在实际应用中，由于快拍数L有限，RX由采样协方差矩阵进行等效，即

(3)

2 基于协方差矩阵重构的DOA估计方法

首先，给出差联合阵列(difference coarray)的定义，定义如下集合

D={xi-xj},xi,xj∈C

(4)

式(4)中，xi为第i个阵元位置，C为所有阵元位置集合。D表示阵列中所有阵元位置差值构成的集合。由于集合D中存在冗余元素，定义集合D中所有不相同的元素组成的集合Du为差联合阵列。

图1所示的互质阵列的差联合阵列如图2所示，虚拟阵元的位置可具体表示为{±(Mnd-Nmd),0≤n≤N-1,0≤m≤2M-1}。

图2 差联合阵列示意图Fig.2 Structure of difference coarray

根据虚拟阵元分布可知，差联合阵列可以被划分为一个均匀阵列和两个非均匀阵列。其中，均匀阵列包含有2MN+2M-1个虚拟阵元以阵元间距d均匀分布且关于原点对称。而非均匀阵列由于存在空洞从而造成虚拟阵元呈现非均匀分布，两个非均匀阵列同样关于原点对称。进一步可知，非均匀阵列可再细分为单个空洞部分和多个空洞部分。现有的大部分方法仅利用了连续部分的虚拟阵元，而未得到有效利用的虚拟阵元数为MN+M+N-1，因此，若能将差联合阵列存在的空洞进行有效填充，则可显著提高可利用的虚拟阵元数，从而获得更大的阵列自由度，实现对更多信号的DOA估计。

实际上差联合阵列与阵列所能得到的波程差是一一对应的。对于图1所示的互质阵列，阵列接收数据协方差矩阵RX具体可表示为：

(5)

在实际应用中，由于受到信噪比和快拍数有限的影响，式协方差矩阵中同一波程差对应的元素并不完全相等，因此，为了充分利用样本信息提高估计的精度，对相同波程差的元素进行求平均运算作为该波程差对应的输出值，即

(6)

式(6)中，Ri(du)表示对应波程差R(du)的第i个元素。

进一步根据存在的波程差元素对协方差矩阵RX进行扩展，构成一个(2MN-N+1)×(2MN-N+1)维的Toeplitz矩阵RT：

(7)

矩阵RT可看作一个数据缺失的协方差矩阵，同时由于入射信号的稀疏性，因此，RT是低秩的。使得利用低秩矩阵重构理论[13]对矩阵中的零元素进行填充成为可能，具体可通过求解如下的优化问题：

(8)

式(8)中，R表示待重构的目标矩阵，η表示误差门限，该值的选取对算法性能的具有重要的影响。P是映射矩阵，一旦互质阵列结构确定，投影映射矩阵随之确定，具体可表示为：

(9)

然而，矩阵秩函数是非连续、非凸的，直接求解矩阵秩最小化问题是一个NP-hard问题。这里，利用迹范数最小化来对秩范数最小化进行凸松弛表示。因此，式(8)进一步可表示为：

(10)

vec(E)～AsN(0,W)

(11)

(12)

进一步可表示为：

(13)

式(13)中，Asχ2((2MN-N)2)表示自由度为(2MN-N)2的近似卡方分布。进一步在式(13)引入参数μ：

(14)

式(14)中，μ是一个决定噪声容错的参数，(在Matlab软件中可利用函数chi2inv(1-ρ,(2MN-N)2)计算)，μ与置信关联ρ相关联。为了获得较高的概率1-ρ，ρ的取值一般很小。例如，ρ=0.001。

通过将式(14)代入式(10)，DOA估计式可表示为：

mintr(R)

(15)

式(15)是一个凸优化问题，可通过CVX[15]优化工具箱进行求解，得到最优估计值R。

进一步对R进行特征分解，即

(16)

式(6)中，Us是大特征值对应的特征向量所构成的(2MN-N+1)×K维信号空间，UN是小特征值对应的特征向量所构成的(2MN-N+1)×(2MN-N+1-K)维噪声空间，Σs是大特征值组成的K×K维对角阵，ΣN为小特征值组成的(2MN-N+1-K)×(2MN-N+1-K)维对角阵。

为了避免角度搜索的复杂性，令P(z)=[1,z,…,z2MN-N]T，则DOA的求根多项式为：

(17)

(18)

3 仿真实验

本节中将通过仿真实验对所提方法的性能进行验证，并与SS-MUSIC方法、Lasso方法以及文献[11]中的EX-Lasso方法进行比较，其中，SS-MUSIC方法的角度搜索间隔为0.1°，Lasso方法和EX-Lasso方法的完备字典间隔为0.1°。以互质阵列作为接收阵列，最小阵元间距d为半波长。

实验1 为验证本文方法在欠定条件下的估计性能，假设9个等功率远场窄带信号入射到互质阵列上，阵元数为6，其中，M=2，N=3，入射角度分别为[-70.12°，-50.34°，-30.02°，-10.25°，10.65°，30.45°，50.12°，60.34°，80.23°]，信噪比SNR=10 dB，快拍数L=1 024，实验结果如图3所示。

图3 四种方法可行性比较Fig.3 Feasibility comparison of four methods

由图3的实验结果可知，本文方法利用6个阵元准确实现了对9个入射信号DOA的估计，而MUSIC方法理论上最多可估计信号数为7，Lasso方法的可估计信号数为8均小于入射信号数，均无法实现对所有信号DOA的准确估计。EX-Lasso方法虽然理论上能够实现对全部信号的DOA估计，但在重构过程中会受到模型失配的影响，从而对部分入射信号的DOA估计效果不佳。

实验2 为验证本文方法的估计精度与信噪比、快拍数之间的关系。假设2个等功率远场窄带信号入射到互质阵列上，阵元数为6，其中M=2，N=3，入射角度分别为[30.13°，60.55°]，分析信噪比、快拍数对四种方法估计精度的影响，每个信噪比、快拍数下进行500次蒙特卡罗实验。图4为快拍数L=1 024，信噪比SNR从-5～15 dB，步长为5 dB时，角度均方根误差随信噪比变化曲线；图5为信噪比SNR=10 dB，快拍数L从500～900，步长为100时，角度均方根误差随快拍数变化曲线。

图4 角度均方根误差随信噪比变化曲线Fig.4 RMSE versus SNR

图5 角度均方根误差随快拍数变化曲线Fig.5 RMSE versus snapshot number

由图4和图5的实验结果可知，在一定的信噪比和快拍数范围内，四种方法的估计均方根误差均随着信噪比的增加和快拍数的增大而减少。在相同的信噪比和快拍数条件下，本文方根的估计精度要明显高于其他三种方法，尤其是在低信噪比和少快拍数时，从而验证了本文方法具有较好的估计性能。这是由于一方面，本文方法对互质阵列的差联合阵列空洞部分的有效填充，提高了可利用的阵列自由度；另一方面，SS-MUSIC方法的估计精度受到角度搜索间隔的影响，Lasso方法和EX-Lasso方法在重构过程中均受到模型失配的影响，而本文文法有效避免了该问题，从而具有更高的估计精度。

实验3 为验证本文方法对相近信号的分辨能力，假设2个相近的等功率远场窄带信号入射到互质阵列上，阵元数为6，其中M=2，N=3，入射角度分别为[30.13°、32.26°]，信噪比SNR=5 dB，快拍数L=500。将四种方法的分辨性能进行比较，图6为四种方法具体的分辨结果。

图6 四种方法的分辨性能比较Fig.6 Comparison of the resolution of four methods

由图6的实验结果可知，当2个入射信号较近时，本文方法仍然能够准确估计出入射信号的DOA信息，而其他三种方法此时失效，体现了本文方法对相近信号较好的分辨能力。

4 结论

本文提出了基于协方差矩阵重构的互质阵列DOA估计方法，该方法首先根据差联合阵列与波程差一一对应的关系，将阵列协方差矩阵扩展为一个Toeplitz矩阵，进一步通过求解迹范数最小化问题实现对矩阵中零元素的恢复，最后利用求根MUSIC方法实现对DOA的估计。仿真实验结果表明，该方法对互质阵列的差联合阵列中空洞部分能够有效进行填充，从而获得更高的阵列自由度。在提高可估计信号数的同时，有效避免了传统稀疏重构中模型失配对估计性能的影响，提高了估计精度和分辨性能。