稀疏双圆阵二维波达方向估计

晏 辉, 司伟建, 夏新凡

(1.哈尔滨工程大学信息与通信工程学院,黑龙江 哈尔滨 150001;2.哈尔滨工程大学先进船舶通信与信息技术工业和信息化部重点实验室,黑龙江 哈尔滨 150001;3.上海无线电设备研究所,上海 201109)

0 引言

波达方向(direction of arrival,DOA)估计作为阵列信号处理的关键技术,已经广泛应用于通信、声纳、雷达等领域[1]。均匀圆阵(uniform circular array,UCA)作为测向系统中常用的阵列之一,因其具备360°的方位角测向能力、容易共形和导向矢量共轭对称等优势受到广泛关注[2]。针对其复杂的导向矢量,文献[3]提出模式空间变换的概念,将圆阵变为虚拟均匀线阵(uniform linear array,ULA),使其导向矢量具备范德蒙结构。文献[4]通过构造实值波束空间变换矩阵,结合多重信号分类(MUSIC)算法和旋转不变子空间(ESPRIT)算法提出基于均匀圆阵的实值波束空间多重信号分类(UCA-RB-MUSIC)算法和基于均匀圆阵的旋转不变子空间(UCA-ESPRIT)算法。但是波束空间变换会引入估计误差。文献[5-6]通过对误差定量分析,采用迭代的方法消除误差主要部分,算法性能得以提高,但是计算量也随之增大。在阵列稀疏情况下,文献[7]针对波束变换带来的误差问题,对导向矢量进行误差补偿,减少了阵元数需求。文献[8]利用阵列流形分离技术,结合传播因子(PM)和求根MUSIC算法估计角度信息,但是会引入映射误差。文献[9]引入互质稀疏阵列,不仅能够减小互耦误差,还能估计角度间隔很小的信源。但上述三种稀疏阵列算法都只能估计一维角度信息,无法估计二维角度信息。文献[10-13]利用双圆阵估计二维角度信息,但所需阵元数较多且没有考虑阵列稀疏的情况。文献[14]在稀疏阵列波束域误差补偿的基础上进行相位校正,可以估计二维相干信号,但是估计的信源数有限。

本文采用稀疏双圆阵列,首先根据两个子阵的旋转不变性,用子阵接收数据协方差矩阵构造波达矩阵,估计信号的俯仰角[15];其次利用波束变换理论,将UCA变为虚拟的ULA,使导向矢量具备范德蒙结构;再次对稀疏阵列的导向矢量进行误差补偿,减少稀疏阵元下波束空间变换带来的误差影响;最后利用求根MUSIC算法估计信号的方位角。所提方法可以应用于阵列稀疏的实际工程场景中。

1 阵列结构和信号模型

假设K个波长为λ的远场窄带非相干信号入射到各向同性的双圆形稀疏阵列,如图1所示。阵列子阵阵元数为M,相邻阵元间距大于λ/2,阵列半径为r,两个子阵间距d=λ/2,入射信号俯仰角φ∈[0,π],方位角θ∈ [0,2π]。

图1 双圆形稀疏阵列结构图

在不考虑阵列误差、通道不一致和互耦影响的前提下,子阵1和子阵2的接收信号矢量表达式分别为

式中:A为M×K维阵列导向矢量矩阵;S(t)=[s1(t),s2(t),...,sK(t)]T为K×1维入射信号矢量,T 表示转置运算;n(t)=[n0(t),n1(t),…,nM-1(t)]T为M×1维噪声矩阵;Φ为相位差矩阵。本文假设噪声为加性高斯白噪声,彼此独立且与信号不相关,噪声功率为δ2。导向矢量矩阵A可以表示为

其中

式中:a(θi,φi)为信号i的导向矢量,其中,φi和θi分别为信号i的俯仰角和方位角;k=2π/λ为信号波数;γn=2πn/M为阵元n的位置。

相位差矩阵Φ可表示为

式中:diag(·)为矩阵对角化函数。

2 基于稀疏阵列的二维DOA估计

2.1 俯仰角估计

文献[15]首次提出波达方向矩阵法,用于双线性平行阵列二维参量估计,文献[10]将其扩展到双圆阵中。该方法首先根据接收信号矢量得到两个子阵的自协方差矩阵RXX和互协方差矩阵RYX,表达式为

其中

定义波达方向矩阵

文献[15]已证明,当导向矢量矩阵A和信号协方差矩阵RSS满秩时,波达方向矩阵R满足

即通过对R的特征分解得到特征值和特征向量,其中非零特征值和对应的特征向量分别与对角阵Φ和导向矢量A对应相等,通过此对等关系即可求出俯仰角和方位角。本文采用文献[10]中的方法对俯仰角求解,计算公式为

式中:Arg(·)为求取复数辐角函数;ηi为R特征分解得到的非零特征值。但是在求解方位角时,首先特征分解操作并不能保证特征向量唯一(存在比例关系),导向矢量矩阵A的对等关系不一定成立;其次由于圆阵方位角θ∈[0,2π],容易出现测角模糊,想要得到准确的角度信息必须进行额外处理,增大了计算难度。

2.2 方位角估计

针对2.1节方位角估计所遇到的问题,引入稀疏阵列下的波束空间变换法。首先利用波束空间变换,将圆阵转化为导向矢量具备范德蒙结构的虚拟线阵;然后进行导向矢量补偿,消除稀疏阵列情况下的误差影响;最后通过求根MUSIC算法估计出信号的方位角。

在得到俯仰角后,利用波束空间变换理论[3],直接用波束变换矩阵W乘以去噪后的接收信号协方差矩阵RXX,0,使UCA的导向矢量具备范德蒙结构。波束变换矩阵的表达式为

其中

式中:ε=ceil(kr)为最大相位模式数,其中ceil(·)表示向上取整。当M＞2ε+1时,波束域协方差矩阵

其中

式 中:Γ(φi)=diag(V-ε(φi),…,Vm(φi),…,Vε(φi))为波束域导向矢量的俯仰角部分,其中Vm(φi)=jmJm(kr sinφi),Jm·（）表示阶数为m的第一类贝塞尔函数;av(θi)为虚拟线阵导向矢量。由此看出经过波束空间变换后,导向矢量中的方位角和俯仰角已经分离。

当M＜2ε+1(阵元间距大于λ/2)时,阵列变为稀疏阵列,波束空间变换带来的残差不可忽略。为了消除误差影响,本文采用文献[7]中的方法,重新定义稀疏下的阵元数为2N+1(N＜ε),则新的波束变换矩阵可表示为

波束域虚拟线阵长度变为2N+1,稀疏阵列第i个信号源的导向矢量可表示为

其中

式中:a0(θi,φi)为稀疏波束域中不含误差的导向矢量部分;Δa(θi,φi)为相位模式数在 [N,ε]之间的波束变换估计误差;Γs(φi)为稀疏阵列导向矢量的俯仰角部分;avs(θi)为稀疏情况下虚拟线阵的导向矢量。新的波束域协方差矩阵表达式为

式中:I为(2N+1)×(2N+1)维单位矩阵;Jl是单位矩阵I的后ε-N列(2N+1)×(ε-N)维矩阵(ε≤3N+1);Jr是单位矩阵I的前ε-N列(2N+1)×(ε-N)维矩阵。则补偿后的波束域导向矢量

式中:Δa′(θi,φi)为补偿后的残差项,可以通过选择合适的模式数达到任意小值。误差消除后利用求根MUSIC算法可以估计出方位角,第i个信号的求根多项式为

则本文所提方法的二维DOA估计具体步骤为:

a)获取两个阵列的接收数据矢量X(t),Y(t);

b)分别求取子阵1的自协方差矩阵RXX和两个子阵的互协方差矩阵RYX,然后根据式(11)构造波达方向矩阵R;

c)对R进行特征分解,然后利用式(13)估计俯仰角φi;

f)根据式(24),利用求根MUSIC算法估计方位角θi。

3 仿真实验

为了验证所提方法的正确性,做如下仿真实验:实验1验证本文方法的有效性;实验2验证本文方法随信噪比变化的情况,对比算法为文献[4]中的UCA-ESPRIT算法和UCA-RB-MUSIC算法;实验3验证高度相关情况下信噪比对本文方法的影响;实验4验证快拍数对本文方法的影响。在构造补偿矩阵P时,本文通过大量仿真实验证明,当ε=3N+1时,估计性能较好。所以本文所有仿真实验均取N=(ε-1)/3。

(1)实验1

假设4个独立信号分别以入射角(100°,35°),(110°,38°),(190°,40°),(220°,45°)入射到图1所示阵列。取r=1.5λ,则N=3,M=7,信噪比为25 dB。仿真的快拍数为128,独立进行100次蒙特卡罗实验。仿真结果如图2所示。对于间隔较近的信号,本文所提方法能够准确估计出结果。

图2 4个信号二维DOA估计结果

针对3个信源情况,可以进一步减少阵元数。取r=λ,则N=2,M=5,其他条件不变,得到入射角分别为(50°,20°),(100°,80°),(150°,50°)的3个独立信号的估计结果,如图3所示。

图3 3个信号二维DOA估计结果

(2)实验2

假设4个独立信号分别以入射角(60°,20°),(120°,30°),(32°,35°),(200°,40°)入射到图1所示阵列,其他仿真条件与实验1相同。为了定量分析,本文算法阵元数M取7,其他两种算法阵元数均取14。均方根误差定义为

式中:n为蒙特卡罗实验次数;(θi,φi)为信号i入射角的真实值;(θij,φij)为信号i入射角的第j次估计值。

信噪比从5 dB增加到29 dB,步长3 d B,进行100次蒙特卡罗实验,仿真得到本文算法、文献[4]中UCA-ESPRIT算法和UCA-RB-MUSIC算法估计结果的均方根误差随信噪比变化情况,如图4所示。

图4 实验2均方根误差随信噪比变化情况

从图4中可以看出,本文算法性能介于UCA-RB-MUSIC算法和UCA-ESPRIT算法之间。相比二维谱峰搜索的UCA-RB-MUSIC算法,本文算法计算量大大减少,运行时间显著提高;相比UCA-ESPRIT算法,本文算法由于多次特征分解,计算量稍大,但是分辨率更高,估计性能更好。

(3)实验3

由于多径传播和障碍物遮挡等原因,实际环境中充斥着大量高度相关信号,需验证高度相关情况下信噪比对本文方法的影响。假设4个相关信号分别以入射角(120°,20°),(60°,40°),(320°,60°),(200°,80°)入射到图1所示阵列,各信号间的相关系数为0.9,其他条件与实验1相同。信噪比从5 d B增加到29 dB,步长3 d B,进行100次蒙特卡罗实验,仿真得到本文算法估计结果的均方根误差随信噪比变化情况,如图5所示。可知所提方法对于高度相关信号具备较好的分辨性能,当信噪比高于15 dB时,均方根误差显著减小。

图5 实验3均方根误差随信噪比变化情况

(4)实验4

假设4个独立信号分别以入射角(120°,20°),(50°,40°),(200°,60°),(270°,50°)入射到图1所示阵列,其他仿真条件与实验1相同。进行100次蒙特卡罗实验,改变快拍数,仿真得到本文方法估计结果的均方根误差随快拍数变化情况,如图6所示。当快拍数大于50时,所提方法估计结果的均方根误差已然很小,算法实时性较好,基本符合实际工程应用中小快拍数要求。

图6 实验4均方根误差随快拍数变化情况

4 结束语

本文利用稀疏双圆阵接收数据矢量的协方差矩阵构造波达方向矩阵,通过波达方向矩阵特征分解估计入射信号俯仰角;然后对子阵1进行波束空间变换,使圆阵变为导向矢量具备范德蒙结构的虚拟线阵;再对稀疏情况下的导向矢量进行补偿,消除了波束空间变换引起的误差;最后利用求根MUSIC算法估计入射信号的方位角。所提方法无需进行复杂的二维谱峰搜索,方位角和俯仰角自动配对,可以应用在阵列稀疏情况的二维测向系统中。虽然双圆阵列损失了一个子阵的阵元,但实际上稀疏情况下总的阵元数并不多,满足工程实际的要求。未来主要的工作是将本文方法应用在信号相干情况下。