锥面共形阵列天线实值盲极化DOA估计方法

张梓轩, 王 飞, 齐子森, 许 华, 梁 佳

(1.空军工程大学信息与导航学院,西安,710077;2.95133部队,武汉,430014)

来波信号极化特征与方位信息的耦合是共形阵列天线波达方向(direction-of-arrival,DOA)估计的难点[1-2]。共形阵列DOA估计大多采用旋转不变子空间算法(estimation of signal parameters via rotational invariance techniques ,ESPRIT)[3],将来波信号的二维角参数与极化参数解耦合,以达到不考虑极化参数而实现高分辨方位估计的效果[2]。但此类算法的子空间估计涉及复值矩阵的特征值分解,计算量庞大的问题较为突出,有必要寻求高效处理方法。

在超分辨空间谱估计30多年的发展历程中,学者们对经典线阵、面阵的高效超分辨算法的研究取得了丰硕的成果[4-9]。高效超分辨算法致力于减少系统的复杂度和算法的运算量,以满足应用系统对多信源实时测向的需求,实值/半实值类算法[10-14]就是其中的典型代表。实现超分辨算法实值化最具代表性的方法是基于阵列接收数据协方差矩阵的酉变换技术,文献[10]充分利用中心对称阵列(如均匀线阵)输出数据协方差矩阵为艾米特中心对称矩阵的特性,通过数据变换,将复数运算转换为实值运算,运算量可缩减75%。文献[12]提出了协方差矩阵分裂思想,对实值运算和阵列结构的任意性进行折中,进而发展了两种全新的基于半实值运算的超分辨算法,摆脱了算法对阵列结构的依赖性,但仍存在镜像方位的模糊问题,需要在算法应用过程中进一步处理。

本文将协方差矩阵分裂思想[15-18]引入到锥面共形阵列天线DOA估计,给出了一种高效的锥面共形阵列天线盲极化DOA估计方法。该方法将子空间估计过程实值化处理,通过子阵分割设计解决了协方差矩阵分裂引起的“镜像模糊”问题,并推导了参数估计的克拉美-罗界(Gramer-Rao bound,CRB),所提算法在保证DOA估计精度的前提下降低了算法复杂度。

1 锥面共形阵列天线快拍数据建模

对于锥面共形阵列,因其具有对称性,可将阵列划分为多个子阵,每个子阵的DOA估计算法均相同。因此,以下进行的算法分析,针对其中的一个子阵进行即可。

锥面共形阵列的子阵设置如图1和图2所示,以圆锥的顶点作为坐标轴的原点O,在圆锥面的母线上划分l1,l2,l33个子阵,每个子阵上有m个阵元均匀分布,并且处于同一条母线上的相邻阵元间距为λ/2(λ为入射信号波长)。

图1 锥面共形阵列(阵元全局分布)

图2 锥面共形阵列(局部阵元与全局分布关系)

子阵l1,l2,l3在全局坐标系XOY面投影如图3所示,3条母线与X轴顺时针的旋转方向的夹角分别为α1,α2,α3,且α1=π-α3,α2∈(α1,α3)。

图3 锥面共形阵列天线俯视图

以上子阵设置可实现方位角φ∈(α1,α3)的估计,所以将整个圆锥上的多个子阵综合起来,可实现φ∈(0,2π)的方位角估计。

由文献[17]给出的共形阵列天线快拍数据建模方法,可得锥面共形阵列的导向矢量为:

(1)

式中:ri为第i(i≤m)个阵元在全局坐标系中,对单位强度来波信号的响应;Gi为第i个阵元在全局坐标系中的位置矢量;u为来波信号的方向矢量;θ和φ分别为来波信号相对于全局坐标系的俯仰角和方位角。

当有n个来波信号时,阵列的流形矩阵为:

A=[a(θ1,φ1),a(θ2,φ2),…,a(θn,φn)]

(2)

所以,此锥面共形阵列的快拍数据模型为:

X=AS+N=(AθKθ+AφKφ)S+N

(3)

Kθ=diag(k1θ,k2θ,…,knθ)

(4)

Kφ=diag(k1φ,k2φ,…,knφ)

(5)

S=[s1,s2,…,sn]T

(6)

N=[n1,n2,…,nn]T

(7)

式中:K为对角阵;S为信号矢量;N为噪声矢量;θ和φ分别为入射信号在全局坐标系中的俯仰角和方位角。

2 实值盲极化DOA估计方法

由协方差矩阵R的公式有:

(8)

由此,可将R-AOCM和I-AOCM表示为:

Re(R)=

(9)

Im(R)=

(10)

式(8)表明接收数据矩阵X是由所有数据,也就是n个快拍数据组成的,它必须包含关于信源DOA的全部信息,而这些信息分别包含于Re(X)和Im(X)中。从这样的观点来看,Re(X)和Im(X)对于信源DOA的估计都是必要的,理论上必须同时用于DOA估计。

不过,根据式(9)和(10)可以直观地发现,它们都同时包含了Re(X)和Im(X),所以,R-AOCM和I-AOCM所包含的信息已经足够寻找所求信源的DOA。因此,单独通过R-AOCM或I-AOCM来求得信源DOA是有可能的。这两个部分单独使用时,也就不存在复数域的计算,即完成了原本算法的实值化,所以,在理论上,将显著降低原算法的计算复杂度。

原始的信号子空间span(S)及其共轭信号子空间span(S*)的交集为:

span(G)span(S)∩span(S*)

(11)

由式(11)可知,span(G)是span(S)的一个子集,它也自然包含了span(S)的一部分向量。因此,可用span(G)替代span(S)进行信源DOA的估计。

但是,span(G)与span(S)有一个根本性差别,即span(G)与span(A)(A为阵列流形)在信源的真实DOA方向上和它们的镜像方向上同时具有双重正交性[11],表示为:

(12)

所以,在对信号协方差矩阵R的实部或者虚部直接进行特征分解处理后,提取信号子空间时,需要提取的特征向量应该翻倍。因为,此时的信号子空间既包含了真实信源方位,还包含着镜像方向上的信源方位,也就是与真实信号等数量的虚假的信号。为了准确提取到真实信号方位的信息,必须将真实信号和虚假信号的特征信息同时提取,进行一系列处理直至求出方位信息后,再加以甄别,最终得出准确的信号方位。

如图1所示,每个子阵对均由同一条母线上分布的m个阵元构造,即1～m-1阵元,构成第1个阵列,2～m阵元构成另一个阵列,2个阵列之间的距离为λ/2。以此类推,可分别得到母线l1、l2、l3上距离矢量各不相同的3个子阵对。

分割母线l1上分布的阵元情况见图4。

图4 对母线上阵元的分割

如图4所示,子阵l11的接收数据为:

X11=A11S+N11=

(A11θKθ+A11φKφ)S+N11

(13)

子阵l12的接收数据为:

X12=A12S+N12=

(A11θKθ+A11φKφ)Φ1S+N12

(14)

式中:Φ1为相位差矩阵。

Φ1=diag[exp(-jw11),exp(-jw12),…,exp(-jw1n)]

(15)

w1i=(2π/λ)dv1ui=πv1ui

(16)

展开即得:

w1i=π[sin(θv1)cos(φv1)sin(θi)cos(φi)+

sin(θv1)sin(φv1)sin(θi)sin(φi)+cos(θv1)cos(θi)]

(17)

X11和X12的最后一行组成了母线l1的快拍数据X1:

(18)

式中:N11表示子阵l11接收的噪声矢量;N12表示子阵l12接收的噪声矢量;v1为子阵l11与子阵l12距离的方向矢量。同理,将l2、l3划分为l21、l22和l31、l32,即可获得其快拍数据X2和X3以及Φ2和Φ3。文献[3]中已证明,矩阵φ1,φ2,φ3与用对应的Φ1,Φ2,Φ3的特征值来分别构成的对角阵相等,此处不再详细阐述。接下来,我们便可通过Φ1,Φ2,Φ3这几个包含旋转不变关系的矩阵,结合式(16)求解信号的方向参数。

阵列的数据协方差矩阵为:

(19)

利用式(9)或(10)替换协方差矩阵,完成矩阵实值化处理,并进行特征分解:

(20)

式中:US为阵列信号子空间的特征矢量矩阵,根据式(18),可以知道US的第1～m-1行为子阵l11对应的信号子空间,将其记为US11;US的第2～m行为子阵l12对应的信号子空间,将其记为US12。同理还可得到US21,US22,US31,US32。

对应ESPRIT里的最小二乘法为:

(21)

式(15)中的exp(-jw1i)(i=1,2,…,n)即为Φ1的特征值t1i。同理可得Φ2的特征值t2i,Φ3的特征值t3i。

此时,结合α1=π-α3,α2∈(α1,α3)分别求解t1i,t2i,t3i,就可求出入射信号的方位角φi为:

φi=arctan

(22)

入射信号的俯仰角θi为:

(23)

zji=angle(tji)/(πλ),j=1,2,3

(24)

至此,在不考虑来波信号的极化参数的情况下,求解得到来波信号的方向信息,实现了盲极化DOA估计。

综上所述,锥面共形阵列天线的实值盲极化DOA估计算法的流程为:

1)在已求出快拍数据矩阵X的基础上,按照式(8)得到此阵列的数据协方差矩阵R;

2)采用实值化算法,在R-ACOM上进行特征分解,得到实值的信号子空间,从而分别得到母线l1、l2、l3对应的Φ1、Φ2、Φ3;

3)对Φ1、Φ2、Φ3进行特征分解,得到t1i、t2i、t3i;

4)根据式(22)～(24),得到对应的n个接收信号的方向信息;

5)通过限定子阵的空域覆盖范围,对估计出的信号俯仰角和方位角进行判定,排除虚假信号,从而得到最终准确的信号方位。

就运算复杂度而言,ESPRIT算法的运算复杂度主要集中在对数据协方差矩阵R的处理上,所以本节仅对协方差矩阵R特征分解部分的计算量进行讨论。表1给出了锥面共形阵列天线盲极化DOA估计算法和锥面共形阵列天线实值盲极化DOA估计算法的运算复杂度对比情况,其中,算法的运算复杂度均以各算法中包含的实数乘法运算次数来表示。表格中4×O(M2L)为给出的计算复数矩阵R的逆矩阵或EVD所包含的实数乘法次数[13]。其中,M为阵元个数,L为仿真实验的快拍次数。

表1 运算复杂度(实数乘法次数)对比

由表1可以看出,锥面共形阵列天线实值盲极化DOA估计算法在运算复杂度上相较于锥面共形阵列天线盲极化DOA估计算法,降低了约75%的运算量。按照此对比,节省时间约75%,但是在实际仿真中,由于本文所采用的实值化算法将真实信号和虚假信号的特征信息同时提取出来,进行一系列处理直至求出方位信息后,再加以甄别,最终得出准确的信号方位。这其中甄别真假信号的过程是原算法所没有的,所以,最终节省时间的百分比会小于75%。具体数值则与这部分操作所占时间比有关。

3 理论性能分析

ESPRIT算法可实现对信源方位角参数的无偏估计,与之对应的CRB给出了无偏参数估计方差的下限,本节将推导实值盲极化DOA估计的CRB。为了简化推导过程,假设信源的相关矩阵RS已知,且噪声方差归一化为1。则阵列协方差矩阵R中包含4n个未知参数,即n个俯仰角参数,n个方位角参数,2n个极化参数。估计参数可用矢量p表示为:

pT=[θ1,φ1,θ2,φ2,…,θn,φn,k1,k2,…,kn,

δ1,δ2,…,δn]

(25)

方位参数与极化参数联合估计的CRB设为BCR,则:

(26)

BCR=F-1

(27)

4n×4n的Fisher信息矩阵F可分块表示为:

(28)

式中:Fθθ为俯仰角估计块;Fφφ为方位角估计块;Fkk为极化矢量幅度比分量估计块;Fδδ为极化矢量相差的分量估计块;其它模块为相应参数估计的互相关块。Fisher矩阵的第i行j列元素Fij为[14]:

trace[DiRSAHR-1DjRSAHR-1]}

(29)

(30)

式中:trace(·)表示取矩阵(·)的迹;∂R/∂pi表示对矩阵R求pi的偏导,则有:

(31)

(32)

(33)

(34)

4 仿真实验

4.1 实验1:估计参数误差与信噪比分析

实验条件:阵列结构如图1所示,极化幅度比为5,相差为30°,实验执行次数5 000次;信噪比-5 dB～15 dB,步进1 dB;阵元数30个,快拍次数100次。波达方向1的俯仰角为65°,方位角为80°;波达方向2的俯仰角为80°,方位角为95°。

(a)俯仰角估计偏差与信噪比分析

由仿真实验可知,所提算法可以解决多辐射源方位估计问题,随着信噪比的增加,本文实值高效算法的DOA估计成功概率、估计方差以及估计偏差都逐步趋近文献[3]中所提算法,当信噪比大于5 dB时,对方位角的估计效果与文献[3]中算法基本相同,对俯仰角的估计效果与文献[3]中算法相比存在一定差距,但是本文算法的俯仰角估计偏差控制在0.5°以内,基本满足对多信源精细分辨的需求。

4.2 实验2:估计参数误差与快拍次数分析

实验条件:阵列结构如图1所示,极化幅度比为5,相差为30°,实验执行次数5 000次;信噪比5 dB;阵元数30个,快拍次数1 000次起,步长100次,结束于5 000次,波达方向1,俯仰角为65°,方位角为80°;波达方向2,俯仰角为:75°,方位角为90°。

(a)俯仰角估计偏差与快拍次数分析

由仿真实验可知,快拍次数的提升对本文算法以及文献[3]中算法的DOA估计效果影响十分明显,当快拍数在3 000次以上时,估计效果随着快拍的增加变化趋于平缓,但整体估计精度的变化并不如信噪比对估计结果的影响明显。本文算法对俯仰角与方位角的估计效果均未随快拍增加,与文献[3]中算法的实验结果趋势一致。这也说明协方差矩阵分裂造成了数据信息的丢失,影响了所采样本的信息含量,带来的差距不能单纯靠数据累计来补偿。

4.3 估计参数误差与阵元数分析

实验条件:阵列结构如图1所示,极化幅度比为5,相差为30°,实验执行次数5 000次;信噪比0 dB;快拍次数200次;阵元数150～600个,步进30个。波达方向1,俯仰角为65°,方位角为80°;波达方向2,俯仰角为75°,方位角为90°。

(a)俯仰角估计偏差与阵元个数分析

由仿真实验可知,阵元个数的增加对参数估计效果的提升十分明显,成功概率、估计偏差以及估计方差均随着阵元个数增加快速变化并与文献[3]中算法的估计效果趋于一致,这也说明参数估计效果对天线口径的变化极其敏感。

4.4 参数估计理论界仿真

实验条件:阵列结构如图1所示,极化幅度比为5,相差为30°,实验执行次数为5 000次;信噪比-5 dB～15 dB,步进1 dB;阵元数30个,快拍次数100次;波达方向1,俯仰角为65°,方位角为80°;波达方向2,俯仰角为75°,方位角为90°。文献[3]中算法的参数估计CRB推导详见该文献。

图8 估计参数的理论界限对比

由仿真结果可知:对协方差矩阵进行实值化处理后,参数的估计性能有损失,估计的理论线均大于实值化处理之前,但随着信噪比的增加,逐步趋于一致。仿真还呈现出俯仰角的估计效果好于方位角,主要原因是图1给出的阵列结构中母线上阵元数多于锥体横截面阵元数,所以对俯仰角分辨能力强,方位角分辨能力则差一些。

4.5 算法仿真时间分析

实验条件:阵列结构如图1所示,极化幅度比为5,相差为30°,信噪比10 dB;快拍次数为100次;运行次数100次;波达方向1,俯仰角为65°,方位角为80°;波达2,俯仰角为80°,方位角为95°;阵元数,150～600,步进150。

表2 算法仿真时间对比

由算法仿真时间对比可以看出,算法经过实值化处理后可显著降低算法复杂度,总体运行时间是文献[3]中所提算法的60%～70%之间,而且天线规模越大,阵元个数越多,复杂度下降比例越大,充分验证了算法的高效性。

5 结语

本文基于协方差矩阵分裂思想,将快拍数据协方差矩阵实值化处理,给出了一种高效的锥面共形阵列盲极化波达方向估计方法,所提算法将子空间估计过程进行实值化,并通过子阵分割解决了协方差矩阵分裂引起的“镜像模糊”问题,推导了实值盲极化DOA估计的克拉美-罗界,对参数估计的理论性能界进行了分析,讨论了算法复杂度改善情况,在保证DOA估计精度的前提下有效降低了算法复杂度.计算机仿真表明,新算法在信噪比为10 dB时,估计精度与原算法基本一致,而运算量仅为原算法的60%～70%,验证了算法的有效性。