双极化面阵的快速波达方向估计算法
作者: 韩乐天 司伟建 曲明超
引用格式:韩乐天,司伟建,曲明超.双极化面阵的快速波达方向估计算法[J].航空兵器,2023,30(1):120-126.
HanLetian,SiWeijian,QuMingchao.AFastEstimationAlgorithm-for-theDirection-of-ArrivalofDoublePolarizedSurfaceArrays[J].AeroWeaponry,2023,30(1):120-126.(inChinese)
摘要:针对目前极化敏感面阵空域-极化域联合谱估计运算量大、耗时长的问题,提出一种降维求根MUSIC(MultipleSignalClassification)优化算法。通过对接收信号进行降维处理,提出新的求解模型将传统四维MUSIC转化为两个一维求根MUSIC求解空域波达方向和引用已求解出的空域信息结合拉格朗日乘子法解决来波信号极化信息估计问题。相比传统的4D-MUSIC和秩亏MUSIC,所提算法在不损失估计精度的前提下提高了运算速度,降低了运算复杂度,无需谱峰搜索过程,消除了因搜索步长而导致的量化误差。对日后大规模阵列计算及MIMO(MultipleInputMultipleOutput)雷达引入提供快速求解方法。仿真实验表明,所提算法在低信噪比0dB下空域误差约为0.85°,速度相比秩亏MUSIC提升了约64.7%,验证了该算法的有效性和高精度性。
关键词:极化敏感面阵;空域-极化域联合谱估计;MUSIC;MIMO;拉格朗日乘子法;信号处理
中图分类号:TJ760
文献标识码:A
文章编号:1673-5048(2023)01-0120-07
DOI:10.12132/ISSN.1673-5048.2022.0112
0引言
波达方向(Direction-of-Arrival,DOA)技术,是一门空域信号处理技术,利用天线阵列得到的空间电磁场信息来估计来波方向,随着电磁环境日益复杂,标量传感器阵列无法满足测向系统的要求。由电磁矢量传感器构成极化敏感阵列,以矢量形式接受电磁波信号,可以同时获得入射方向的空域信息和极化信息[1-3],有效地提高阵列信号处理的整体性能,改善了空间源信号多维参数估计、自适应波束形成[4-5]等性能,还可利用期望信号和噪声信号在极化域的不同,进行降噪滤波[6]。但引入极化信息提升信号处理性能的同时,也带来了较高的计算负担,除了需要估计二维空域信息外,还需进一步估计二维的极化信息,计算量大、运算复杂度高。
文献[7]将子空间旋转子不变参数估计方法(EstimationofSignalParametersviaInvarianceTechniques,ESPRIT)[8]引入极化敏感阵列测向中,实现了多信号DOA和极化参数的估计。文献[9-10]将MUSIC[11-12]算法推广到了极化敏感矢量阵列,但传统的极化MUSIC需进行四维谱峰搜索,不仅计算速度慢且复杂度高,已经不适用于极化敏感阵列测向技术。文献[13]提出了一种针对一维极化敏感线阵列的求根MUSIC算法,降低了复杂度,避免了谱峰搜索过程,但传统的测向阵列多为面阵,此算法无法满足相应要求。文献[14-15]提出了针对极化敏感阵列的降维秩亏MUSIC算法和极化模值约束降维算法,将传统的四维搜索缩减为二维,利用优化算法估计极化信息,但无法消除因谱峰搜索步长导致的量化误差,且运算过程复杂度较高。文献[16]提出一种针对普通标量阵的二维求根MUSIC,降低了运算复杂度且无需谱峰搜索。
本文提出一种基于双极化敏感面阵的二维求根优化MUSIC,首先提出一种新的求解模型,将传统的四维空域极化域信息通过化简转换为针对空域信息的二维信号估计,将求解空域二维根转换为求解两个一维根的方法解决空域信息的估计,采用拉格朗日乘子法估计极化信息,无需谱峰搜索,大大降低了运算复杂度,提高运算速度,消除了因搜索步长设置导致的量化误差。
1数据模型
1.1极化阵列摆放及其接收模型
本文接收单元采用双正交偶极子极化敏感阵元,组成M行,N列极化敏感面阵,排列方式如图1所示,相邻两正交偶极子阵列间距d=λ/2。
假设空间有K个远场点源,发射的窄带完全极化不相关电磁信号入射到接收阵列中。入射信源数应符合K≤2MN,K个信号的仰角和方位角分别为θk和φk,其中仰角θ为信号源入射信号与Z轴夹角,方位角φ为信号源入射信号与X轴夹角,取值范围分别为θk∈[0,π/2]和φk∈[0,2π],θ=0°时表示信号源正对天线阵入射。极化辅助角信息和极化相位差信息为γk和ηk,其中极化辅助角γ为两电场分量强度之间关系,极化相位差η为垂直分量超前水平分量的相对相位差,取值范围分别为γk∈[0,π/2]和ηk∈[-π,π],假设输入噪声为空-时-极化白噪声,且噪声与入射信号独立,入射信号互不相干。
设噪声信号矩阵为N(t),入射信号矩阵为S(t)=[s1(t),s2(t),…,sk(t)]T,则极化敏感阵列信号接收模型为
X=[x1(t),x2(t),…,xm(t)]T=
∑Kk=1a(θk,φk,γk,ηk)sk(t)+N(t)=
AS(t)+N(t)(1)
式中:A=[aθ1,φ1,γ1,η1,aθ2,φ2,γ2,η2,…,aθk,φk,γk,ηk]为空域-极化域阵列流形,aθk,φk,γk,ηk为第k个信源的空域-极化域导向矢量,表示为空域阵列流形与极化敏感阵元接收矢量的克罗内克积;aθ,φ,γ,η=sssp,ss=Ay⊙Ax=[ay(θ1,φ1)ax(θ1,φ1),…,ay(θk,φk)ax(θk,φk)]为空域导向矢量,ay(θk,φk)=[1,…,ej2π(M-1)dsinφksinθk/λ]T和ax(θk,φk)=[1,…,ej2π(N-1)dcosφksinθk/λ]T分别为Y轴和X轴方向阵元的导向矢量。
针对电磁信息不完备的极化敏感阵元,如双正交偶极子阵元,每个阵元只能接收到X和Y轴方向上的电场矢量,在入射波为完全极化信号的条件下其极化矢量为
sp=-sinφcosθcosφ
cosφcosθsinφcosγsinγejη(2)
针对不同种类的极化敏感阵元,极化矢量取决于对电磁的敏感程度。
1.2空域-极化域联合谱估计
阵列接收信号的协方差矩阵R^可以表示为
R^=1L∑Ll=1XHX(3)
式中:L为快拍数。
对协方差矩阵R^进行奇异值分解得到
R^=UsΣsUHs+UNΣNUHN(4)
式中:Us为K个较大特征值对应的特征向量张成的信号子空间;UN为2MN-K个较小特征值对应的特征向量张成的噪声子空间。
在理想条件下,接收空间的信号子空间与噪声子空间是相互正交的,即信号子空间中的导向矢量也与噪声子空间正交:
aHθ,φ,γ,ηUN=0(5)
但因为空间中噪声、干扰等情况的存在,得到的信号子空间并不能做到和噪声子空间完全的正交,因此,DOA是以最小优化搜索的条件实现的,即
θMUSIC=argθminaHθ,φ,γ,ηU^NU^HNaθ,φ,γ,η(6)
所以,传统MUSIC算法的谱峰估计公式为
PMUSIC=1aHθ,φ,γ,ηU^NU^HNaθ,φ,γ,η(7)
通过遍历(θ,φ,γ,η)四个参数,寻找出对应最大峰值的K个方位角、仰角、极化辅助角和极化相位差,则可得到K个来波信号的参数估计值。
但传统的MUSIC需要搜索四维信息,所需要运算的时间长、运算量大、算法复杂度高且难以实现。文献[14]提出的秩亏MUSIC算法,将四维谱峰搜索变为二维谱峰搜索。
利用下式替代式(6):
{θ,φ}=argmaxθ,φdet-1{Hθ,φ}(8)
式中:Hθk,φk=DHθk,φkU^NU^HNDθk,φk为经过秩亏变化后的谱峰函数,Dθk,φk=uθk,φkBθk,φk,uθk,φk=diag{ss}取空域阵列流形构成对角阵;sp=θk,φk·hγk,ηk,θk,φk=-sinφkcosθkcosφkcosφkcosθksinφk为空域极化域的坐标转换因子,与矢量阵元位置有关,hγk,ηk=cosγksinγkejηk为第K个信号的极化矢量;B=[b1x,b1y,…,bmx,bmy]T为广义极化敏感阵列,与极化敏感阵列类型与摆放位置有关,本文为电敏感的双正交偶极子bmx=[10],bmy=[01]。
文献[14]将空域-极化域四维搜索降维到二维空域角度搜索,后通过优化计算得到极化信息,一定程度上降低了运算量,提高了运算速度。但是因为需进行二维谱峰搜索,搜索步长导致的量化误差无法消除。本文提出一种二维求根约束MUSIC算法,通过构造多项式求根的方法求解空域信息,利用拉格朗日乘子法计算得出极化信息,无需谱峰搜索,进一步提高了运算速度,降低运算复杂度,且消除了因为搜索步长导致的量化误差。
2本文算法
2.1空域谱估计
利用空域导向矢量和极化域矢量信息的克罗内克积替换空域-极化域阵列流形,构造谱峰搜索函数:
PMUSIC=1aHθ,φ,γ,ηU^NU^HNaθ,φ,γ,η=
1[sssp]HU^NU^HN[sssp](9)
令
ψ(θ,φ,γ,η)=[sssp]HU^NU^HN[sssp](10)
通过克罗内克积的运算性质可得
ψ(θ,φ,γ,η)=[sssp]HU^NU^HN[sssp]=
sHp[ssI2]HU^NU^HN[ssI2]sp(11)
式中:I2为2×2维单位阵;矩阵类型、维数与阵列种类及摆放位置有关,本文应用为双正交偶极子且沿X,Y轴正方向摆放则该矩阵为2×2维单位阵。通过数学分析,sp只有在仰角θ=π/2,且极化辅助角γ=π/2时,矩阵不满秩。当仰角θ=π/2时,说明信号源入射方向来自于水平方向,在工程应用中此类情况常不做考虑,故可将矩阵sp近似看作行满秩,由矩阵运算的秩的关系中可知:
若A行满秩,则Rank(BA)=Rank(B);
若A列满秩,则Rank(AB)=Rank(B)。
故可知针对式(11)中的sp,因其满秩因此相乘任意矩阵不影响矩阵的秩,且可知,为使ψ=0必有解,可用下式取代式(11)进行求根计算:
U=[ssI]HU^NU^HN[ssI](12)
将式(12)代入(11):
ss=ay(θ,φ)ax(θ,φ)(13)
得到
P=[ay(θ,φ)axo(θ,φ)]H·