CN108872928A - 基于四线性分解的平面阵列角度的估算方法 - Google Patents

Abstract

本发明公开了一种基于四线性分解的平面阵列角度的估算方法，其根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；通过最小二乘方法估计目标的DOD与DOA；所述估算方法能有效的对抗非均匀白噪声，其在小快拍和大快拍条件下均能获得较好的估计效果；其参数估计的精度高于ESPRIT算法和HOSVD算法，但复杂度低于PARAFA C算法。

Description

基于四线性分解的平面阵列角度的估算方法

技术领域

本发明涉及一种雷达信号处理技术，具体的说涉及一种基于四线性分解的平面阵列角度的估算方法。

背景技术

阵列信号处理作为信号处理的一个重要分支，其本质是利用空间分散排列的传感器阵列和多通道接收来获取信号的时域和空域等多维信息，达到检测信号和提取其参数的目的。它在通信、雷达、声呐、地震勘探和射电天文等领域获得了广泛的应用和迅速的发展，有很大的应用前景。近年来，利用阵列信号处理进行多参数估计得到了越来越多的关注，尤其是空间波达角(direction of arrival，DOA)估计已成为研究的热点问题之一。二维参数估计的问题是一个非常贴合实际的问题，传统的DOA估计一般是指在一个平面内对信号源的方位角进行估计，属于一维参数估计问题，存在许多局限性。在实际通信***中，用户并非位于同一平面，因此需要二维波达角来表示。

二维DOA估计在过去的数十年间取得了广泛的关注。至今已涌现出大量优秀的估计算法，典型的主要代表有多重谱峰分类(multiple signal classification algorithm，MUSIC)算法、基于旋转不变技术的参数估计 (estimation of signal parameters viarotational invariance techniques， ESPRIT)算法、传播算子(propogator method，PM)算法。其中，MUSIC算法是利用噪声子空间和信号子空间的正交特性，通过谱峰搜索来估计DOA，计算量大，复杂度高；ESPRIT算法不需要进行谱峰搜索，但仍然需要奇异值分解或者特征值分解。以上两种方法都是特征值子空间类算法，对于面阵而言，特征值分解运算量很大。而PM算法不需要特征值分解，因而复杂度较低。以上算法往往是利用矩阵分解技术进行参数估计，其往往忽略了阵列数据固有的多维结构。张量分解技术能够有效的挖掘阵列数据内部的多维结构，因而近年来收到越来越多的关注。典型的张量分解技术主要有高阶子空间分解(high-order singular value decomposition，HOSVD)算法和平行因子(parallel factor，PARAFAC) 算法，其中PARAFAC算法采用迭代的方法获得参数估计，相比HOSVD算法，其估计精度往往更高，且复杂度更低，因而更适合大规模阵列信号处理。

然而，上述算法都是假设在均匀的高斯白噪声条件下进行的，在该假设下，噪声的协方差矩阵与单位矩阵是一个线性关系，因此能够比较容易的将信号和噪声进行分离，但是这种假设往往与实际不符，实际工程中，受阵元增益相位不一致等因素的影响，噪声往往是非均匀的。实际上，高斯白噪声可被视为非均匀噪声的一种特例。在非均匀噪声条件下，噪声协方差矩阵为一个未知的对角矩阵。此时，信号和噪声无法有效分离，从而引起角度估计性能的下降，特别是在低信噪比条件下，现有算法有可能会失真乃至完全失效，不适宜于工程运用。

发明内容

鉴于以上原因，有必要提供一种能够降低大快拍背景下PARAFAC复杂度，并可以消除非均匀噪声条件下的影响，而且适用于实际工程中的硬件运算的基于四线性分解的平面阵列角度的估算方法。

本发明提供一种基于四线性分解的平面阵列角度的估算方法，所述基于四线性分解的平面阵列角度的估算方法包括如下步骤：

S1、根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；

S2、对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；

S3、利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；

S4、通过最小二乘方法估计目标的DOD与DOA；

其中，所述非均匀噪声协方差矩阵的具体公式如下：

Q_est＝R_est-diag{R_est}

上式中，diag{·}表示对角化操作；R_est＝XX^H/L；L表示快拍数；D_m(A_y)表示将A_y的第m行的元素按顺序构成一个对角矩阵。

本发明所述基于四线性分解的平面阵列角度的估算方法，其利用阵列信号的多维结构，根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；通过最小二乘方法估计目标的DOD与DOA；本发明所述双基地MIMO雷达角度估算方法其与传统算法相比，能有效的对抗非均匀白噪声，其在小快拍和大快拍条件下均能获得较好的估计效果；还能够正确估计目标的方位信息，并能自动配对所估计的角度；其参数估计的精度高于ESPRIT算法和HOSVD算法，但复杂度低于PARAFA C 算法，且由于降噪的协方差矩阵能抑制非均匀白噪声的影响，因而其性能会优于其对应的非降噪的算法，同其没有孔径损失，因而其有更高的估计精度。

附图说明

图1是平面阵列角度估计示意图；

图2是所述估算方法在SNR＝-10dB时的散点图结果；

图3是所述估算方法与现有其他算法在不同SNR条件下RMSE性能的比较结果；

图4是所述估算方法与现有其他算法在不同L条件下RMSE性能的比较结果。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明，应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

本发明提供一种基于四线性分解的平面阵列角度的估算方法，其特征在于，所述基于四线性分解的平面阵列角度的估算方法包括如下步骤：

S1、根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；

S2、对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；

S3、利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；

S4、通过最小二乘方法估计目标的DOD与DOA。

具体的，如图1所示，该平面矩阵有N×M个阵元，N和M分别表示x和y 方向的阵元个数。阵元间距均为λ/2，λ表示发射信号波长。假设阵列远场中有K个非相干元，第k(1≤k≤K)个阵源的方位是其中为目标相对阵元的方位角，θ_k为目标相对于阵元的俯仰角，则空间第i个阵元与参考阵元0 之间的波程差为

面阵一般在x-y面内，所以z_i一般为0，(x_i,y_i)为第i个阵元的坐标。假设x轴上N个接收阵元的方向矩阵为A_x，y轴上M个接收阵元的方向矩阵为A_y，则如图1所示的子阵1的方向矩阵A_x，而子阵2的方向矩阵就得考虑沿y轴的偏移，每个阵元相对于参考阵元的波程差就等于子阵1的阵元的波程差加上所以可得每个子阵列的方向矩阵分别为

其中，D_m(A_y)表示将A_y的第m行的元素按顺序构成一个对角矩阵。A_x和A_y分别如表达式(3)和表达式(4)所示。

上式中表示x轴上的N个接收阵元的导引矢量，表示y轴上M个接收阵元的导引矢量。则阵列输出可写成

其中，为信源矩阵，表示非均匀噪声矩阵。根据表达式6，方向矩阵也可以表示为⊙表示Khatri-Rao乘积。定义如下协方差矩阵:

上式中E{.}和{.}^H分别表示数学期望和共轭转置。假设N与S是非相关的，根据表达式(5)和表达式(7)得

R＝E{XX^H}＝APA^H+Q 表达式(8)

在工程中，信号协方差矩阵可通过L个快拍的数据估计R_est＝XX^H/L，L表示快拍数。

考虑信源是非相关的，则

其中p_k表示第k个信号的功率。当噪声为非均匀白噪声时，噪声协方差矩阵为

其中，q_i表示第i个阵元的噪声功。令a_t为A的第t(t＝1,2,…,MN)行。通过分析发现，矩阵APA^H的第t个对角元素为即其为一个常数。通过如下简化的协方差矩阵可有效的抑制非均匀噪声，且信号特征不受影响，则简化的非均匀噪声协方差矩阵为

Q_est＝R_est-diag{R_est}＝APA^H-λI 表达式(11)

上式中，diag{·}表示对角化操作，I为一个单位矩阵。

表达式(11)也可以表示为如下一个四阶张量

式中，分别为A_x，A_y的共轭，表达式(12)定义了一个四线性模型，R也称为的对称多模展开形式。也可以表示成如下四线性分解的形式

其中，为一个对角张量，其对称多模展开为P，为了估计K个目标的DOA，需要估计方向矩阵A_x和A_y，协方差张量可以通过计算一个最佳近似来实现,即

||·||_F表示Frobenius范数。ALS(alternating least squares)的核心是使一个矩阵变化而其他矩阵固定不变，然后将张量优化问题转化为线性最小二乘问题。根据表达式(14)，将展开成n(n＝1,2,3,4)模可以得到四个矩阵：

式中，表示的模-n展开形式[10]。表达式(15)可以进一步表述成如下联合优化问题：

利用最小二乘准则可知A_x，A_y，的最小二乘估计分别为

其中A₁,A₂,A₃,A₄分别代表最小二乘估计中的A_x,A_y,根据表达式 (8)，P可由如下表达式求得

则可以直接求得P的最小二乘估计

ALS算法通过迭代的方法获取矩阵的估计，其核心是假设构成张量的所有矩阵中仅一个未知，而其余矩阵均认为已知，利用最小二乘准则获取该矩阵的估计。例如，固定P,A_y,通过最小二乘准则求A_x；固定P,A_x,求A_y；然后固定P,A_x,A_y,求固定P,A_x,A_y,估计最后假设A_x,A_y, 固定，更新P。ALS算法重复上述步骤直到算法收敛，即迭代次数到达一个阈值，或或|ε^new-ε^old|/ε^old≤10^-10，ε^new和ε^old分别表示上次迭代和当前迭代的残差。

A_x,A_y,P可以被随机初始化，然而这种初始化会导致算法收敛速度缓慢。在本发明专利中，方向矩阵的初值可以通过PM算法获取，以加速算法收敛。由于A_x，A_y分别和是共轭的，因此在ALS算法中只需要使用最小二乘准则更新A_x，A_y，P即可。

下文中的定理1给出了本发明所述估算方法参数估计唯一性的条件。

定理1：对于本发明中的四线性模型，假设的k-秩分别为和如果其满足

那么对于ALS所估计的A₁，A₂，A₃和A₄，除了列模糊和尺度模糊外，参数估计是唯一的。

四线性分解的参数估计的模糊性可以表述为如下分解矩阵的形式

其中∏是置换矩阵，N₁，N₂，N₃，N₄分别表示对应的拟合误差矩阵,Δ₁，Δ₂，Δ₃，Δ₄为四个对角矩阵，其对角元素分别表示相应的尺度因子，且其满足Δ₁Δ₂Δ₃Δ₄＝I_K。注意A₁，A₂的相位仍然是线性的，此时，仍可以利用LS方法估计 DOA。利用LS的方法对第k个目标的DOA进行拟合的方法为

其中， phase(·)表示求复数的相位，和分别表示A₁，A₂的第k列。由此得到c₁，c₂的LS解为

显然的第二个元素分别为和的估计值。因此可求出θ_k和的估计值为

针对本发明提出的非均匀噪声背景下均匀矩形阵列二维DOA估计方法，进行了部分数值仿真实验。仿真中假设K＝3个信源处于远场，其二维DOA分别为 和信源矩阵为高斯矩阵。仿真实验中X轴上阵元的个数N＝6，Y轴上阵元数的个数M＝8，二者均为均匀线性阵列，阵元间距均为发射信号半波长。非均匀噪声协方差矩阵对角元素为通过 MATLAB中randn函数随机产生的序列。仿真中的信噪比(signal-to-noise ratio， SNR)定义为

附图2是本发明所提算在SNR＝-10dB、快拍数L＝100时，方向角与仰角估计的散点图。由附图2可知，本发明所述估算方法能够正确估计目标的方位信息，并能自动配对所估计的角度。

为比较本发明所述估算方法同现有算法的估计精度，对算法进行500次蒙特卡洛仿真，角度估计的精度用均方根误差(Root Mean Squared Error，RMSE) 进行性能评价，其中RMSE定义为

式中和分别为第i次蒙特卡洛仿真中获得的对θ_k与的估计。为进一步说明算法性能，将本发明的算法同降噪的ESPRIT算法(标记为ESPRIT)、直接张量法(标记为HOSVD)、降噪的协方差张量法(标记为HOSVD-2)、及PARAFAC 算法(标记为PARAFAC)的性能相比较，其中，前缀降噪表示利用降噪的协方差矩阵进行参数估计。

附图3给出了本发明所述估算方法(标记为proposed)同现有算法在不同 SNR条件下性能比较的结果，其中L＝500。由仿真结果可知，随着信噪比的增加，所有算法的估计精度均有所提升。由于非均匀噪声的影响，在信噪比较低时， HOSVD和PARAFAC的性能都会有所下降。由于降噪的协方差矩阵能抑制非均匀白噪声的影响，因而其性能会优于其对应的非降噪的算法。同时，由于在最终角度估计时本发明所采用的是LS拟合的方法，因而其没有孔径损失，因而其有更高的估计精度。

附图4给出了本发明所述估算方法(标记为proposed)和同现有算法在不同L条件下性能比较的结果，其中SNR＝0dB。由仿真结果可知，随着快拍数的增加，所有算法的估计精度均有所提升。但是受孔径损失或非均匀噪声的影响， ESPRIT算法、HOSVD算法以及PARAFAC算法均会带来估计性能的下降。而本发明所述估算方法能有效的对抗非均匀白噪声，其在小快拍和大快拍条件下均能获得较好的估计效果。

本发明所述基于四线性分解的平面阵列角度的估算方法，其利用阵列信号的多维结构，根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；通过最小二乘方法估计目标的DOD与DOA；本发明所述双基地MIMO雷达角度估算方法其与传统算法相比，能有效的对抗非均匀白噪声，其在小快拍和大快拍条件下均能获得较好的估计效果；还能够正确估计目标的方位信息，并能自动配对所估计的角度；其参数估计的精度高于ESPRIT算法和HOSVD算法，但复杂度低于PARAFAC 算法，且由于降噪的协方差矩阵能抑制非均匀白噪声的影响，因而其性能会优于其对应的非降噪的算法，同其没有孔径损失，因而其有更高的估计精度。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (5)

1.一种基于四线性分解的平面阵列角度的估算方法，其特征在于，所述基于四线性分解的平面阵列角度的估算方法包括如下步骤：

S1、根据接收阵元的阵列信号模型构建非均匀噪声协方差矩阵，并进一步构建非均匀噪声协方差矩阵的四阶张量；

S2、对接收方向矩阵、发射方向矩阵、目标协方差矩阵进行初始化，将接收方向矩阵、发射方向矩阵、目标协方差矩阵初始化为一个单位矩阵；

S3、利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算，直至满足收敛条件；

S4、通过最小二乘方法估计目标的DOD与DOA；

其中，所述非均匀噪声协方差矩阵的具体公式如下：

Q_est＝R_est-diag{R_est}

上式中，diag{·}表示对角化操作；R_est＝XX^H/L；L表示快拍数；D_m(A_y)表示将A_y的第m行的元素按顺序构成一个对角矩阵。

2.根据权利要求1所述基于四线性分解的平面阵列角度的估算方法，其特征在于，所述交替最小二乘法通过随机矩阵分别对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行初始化。

3.根据权利要求1所述基于四线性分解的平面阵列角度的估算方法，其特征在于，所述协方差矩阵的四阶张量模型如下：

上式中，为协方差矩阵的四阶张量，为一个对角张量，其对称多模展开为P；A_x为发射方向矩阵；A_y为接收方向矩阵，分别为矩阵A_x，A_y的共轭。

4.根据权利要求3所述基于四线性分解的平面阵列角度的估算方法，其特征在于，步骤S3中所述利用交替最小二乘法对接收方向矩阵、发射方向矩阵以及目标协方差矩阵进行迭代计算的具体公式如下：

5.根据权利要求4所述基于四线性分解的平面阵列角度的估算方法，其特征在于，步骤S4中所述通过最小二乘方法估计目标的DOD与DOA的具体方法如下：

利用LS方法对DOD和DOA的拟合分别是

上式中， phase(·)表示求复数的相位，和分别表示A₁，A₂的第k列，由此得到c₁，c₂的LS解为

显然的第二个元素分别为和的估计值，因此可求出θ_k和的估计值为