CN115587281A - 基于因子矩阵先验的阵元失效mimo雷达角度估计方法 - Google Patents

Abstract

本发明公开了一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，步骤如下：对存在阵元失效的双基地MIMO雷达接收信号进行匹配滤波，得到虚拟阵列的输出信号

并计算获得虚拟阵列协方差矩阵

构造四阶协方差张量

并将其表示成张量CANDECOMP/PARAFAC分解模型；建立具有因子矩阵先验约束的张量填充模型；将张量填充模型转换无约束的增广拉格朗日函数形式；利用ADMM算法迭代求解增广拉格朗日函数，迭代结束时获得因子矩阵U⁽¹⁾，U⁽²⁾和对角矩阵Δ⁽³⁾，Δ⁽⁴⁾；构造出完整协方差张量

再通过对称厄米特展开还原为完整的协方差矩阵R，最后采用ESPRIT算法进行目标角度估计。本发明能够有效地恢复出四阶协方差张量中多个切片缺失数据，提高阵元失效下MIMO雷达角度估计性能。

Description

基于因子矩阵先验的阵元失效MIMO雷达角度估计方法

技术领域

本发明涉及MIMO雷达角度估计，特别是一种基于因子矩阵先验的阵元失效 MIMO雷达角度估计方法。

背景技术

与传统相控阵雷达的机制不同，多输入多输出(Multiple Input MultipleOutput， MIMO)雷达利用多个发射天线同时发射多路正交信号波形，多个接收天线接收 目标回波信号并通过匹配滤波器分离各通道目标信息，从而提高最大可识别目标 数、抗干扰能力以及参数估计精度。

目标方位角估计是阵列信号处理中的一个主要研究方向。双基地MIMO雷达 采用收发分置的阵列配置，可以在接收端同时估计出目标的发射角(direction of departure，DOD)和接收角(direction of arrival，DOA)。目前，针对双基地MIMO 雷达目标角度估计问题，学者们提出了利用虚拟阵列协方差矩阵进行目标角度估 计的MUSIC、Capon、ESPRIT等算法，其中ESPRIT算法因避免二维谱峰搜索， 大幅减少了运算量而得到广泛关注。在实际应用中，由于工作环境恶劣、元器件 老化等原因，在长期使用的MIMO雷达中收发阵列可能会出现失效阵元，这破坏 了雷达阵列的几何结构，使得虚拟阵列中出现大量失效阵元，从而导致虚拟阵列 协方差矩阵中出现大量整行整列的数据缺失，导致雷达角度估计性能不可避免地 下降。因此，有效补偿失效阵元所致的缺失数据对MIMO雷达角度估计影响是极其重要的。

Chen等在论文《Joint Sensor Failure Detection and Corrupted CovarianceMatrix Recovery in BistaticMIMO Radar With Impaired Arrays》(IEEE SensorsJournal,2019, 19(14):5834-5842)中，提出一种基于块Hankel矩阵重构的MIMO雷达角度估计 算法，该方法将协方差矩阵构造为四重Hankel矩阵，使得重构后的矩阵中每一行 和每一列均存在非零元素，然后采用矩阵填充算法(Matrix Completion，MC)恢 复出缺失数据，但由于Hankel操作会使得矩阵维度扩大，因此该方法计算复杂度 高，运算时间长。现有MIMO雷达阵元失效下角度估计方法大多都是在矩阵框架 下恢复故障阵元数据，但矩阵形式不能很好地反映数据间的多维约束关系，因此 此类方法的失效阵元缺失数据恢复性能有待进一步的提升。为了利用张量数据的 多维结构，可将MIMO雷达的虚拟阵列协方差数据矩阵构造出四阶协方差张量。 但由于阵元失效，MIMO雷达四阶协方差张量中存在部分切片完全缺失的情况， 此时利用现有张量填充算法无法有效恢复这种结构性缺失数据。Ying等人在论文 《Hankel Matrix Nuclear Norm Regularized Tensor Completion for N-dimensional Exponential Signals》(IEEETransactions on Signal Processing,2017,65 (14):3702-3717)中，提出一种基于Hankel矩阵核范数正则化的张量填充算法 (HankelMatrix Nuclear Norm Regularized Tensor Completion，HMRTC)。该方法 联合利用张量的低CANDECOMP/PARAFAC(CP)秩和因子矩阵的范德蒙德结构， 对由每个因子矩阵的列向量所构造Hankel矩阵施加低秩约束，能从少量采集数据 中重建多维完整信号，从而能恢复张量中结构性缺失的数据。由于MIMO雷达协 方差张量经CP分解后的因子矩阵具有范德蒙德结构，因此HMRTC算法能适用 于MIMO雷达失效阵元缺失数据的重建问题。但是HMRTC算法在恢复缺失数据 时未考虑张量分解后的因子矩阵之间的关系，其失效阵元缺失数据的恢复性能不 能达到最优。在阵元失效MIMO雷达中，利用回波信号的多维结构特征以及张量 因子矩阵的先验信息来提高失效阵元数据的恢复精度，对提高MIMO雷达抗阵元 故障能力是非常有必要的。

发明内容

发明目的：本发明的目的是提供一种基于因子矩阵先验的阵元失效MIMO雷 达角度估计方法，以快速恢复阵元失效下MIMO雷达中的缺失数据，从而提高阵 元失效下的角度估计精度。

技术方案：本发明所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度 估计方法，包括以下步骤：

(1)对存在阵元失效的双基地MIMO雷达接收信号进行匹配滤波，得到虚 拟阵列的输出信号

并计算获得虚拟阵列协方差矩阵

(1.1)双基地MIMO雷达的发射阵列和接收阵列皆为均匀线性阵列，分别 由M个发射阵元和N个接收阵元组成；对阵元失效下双基地MIMO雷达所采集 到的回波信号进行匹配滤波，经过脉冲积累后，获得MN个虚拟阵元输出的回波 信号矩阵为：

其中，

为在Q个快拍下虚拟阵列回波信号矩阵；

和

分别为存在失效阵元时的发射阵列和接收阵列流形矩阵；

为高 斯白噪声矩阵；S∈^P×Q为目标系数矩阵，其中P为相互独立的远场目标；⊙表示 Khatri-Rao积；

(1.2)Ω_T和Ω_R分别为失效发射阵元和失效接收阵元的位置集合，当发射阵 列中第p_t∈Ω_T个阵元出现故障时，其流形矩阵

中的第p_t行全为零；当接收阵列 中第p_r∈Ω_R个阵元出现故障时，流形矩阵

中的第p_r行全为零；虚拟阵列协方 差矩阵在Q脉冲周期下的最大似然估计为：

其中，(·)^H表示共轭转置，R_s＝SS^H/Q＝diag(ρ)表示信号协方差矩阵，diag(ρ) 表示由向量ρ生成的对角矩阵，ρ＝[ρ₁,ρ₂，…,ρ_P]，ρ_p,p＝1,2,…,P，表示第p个目标 的反射系数，

表示噪声协方差矩阵。

(2)利用虚拟阵列协方差矩阵

构造四阶协方差张量

将协方差张量

表 示成张量CANDECOMP/PARAFAC(CP)分解模型。

(2.1)双基地MIMO雷达虚拟阵列协方差矩阵

能够重新排列成一个四阶 协方差张量

其中，(·)^*表示共轭；

表示Tucker运算。

(2.2)重新调整因子矩阵

使得g中的元素都为1，此时将g消 除得到协方差张量的简化表达式

(2.3)由于CP分解的置换和缩放效应的扰动，协方差张量

的因子矩阵表 示为

和

其中

和

表示四 个因子矩阵；

表示置换矩阵；Λ₁、Λ₂、Λ₃和Λ₄为P×P的实值对角矩阵， 对角元素对应着缩放系数；N₁、N₂、N₃和N₄表示拟合误差。

(2.4)由于MIMO雷达的发射和接收导向矩阵A_t、A_r均为范德蒙德结构， 而因子矩阵

和

的每一个列向量是用不同的常数对A_t和A_r的对应列向量分别 缩放得到的，因此当忽略拟合误差项影响时，矩阵

和

也具有范德蒙德结构， 即因子矩阵U⁽ⁿ⁾,n＝1,…,4具有范德蒙德结构。

(2.5)忽略拟合误差项的影响，得到因子矩阵之间的关系为

和

其中对角矩阵

对角 矩阵

(3)利用因子矩阵的范德蒙德结构特点、因子矩阵之间的相互关系以及张量 的低CP秩，建立具有因子矩阵先验约束的张量填充模型。

为有效恢复出失效阵元的缺失数据，利用因子矩阵的范德蒙德结构特点、因 子矩阵之间的相互关系以及张量的低CP秩，建立用于失效阵元缺失数据恢复的 张量填充模型：

其中，M⁽ⁿ⁾,n＝1,2,3,4，表示辅助矩阵，Ω表示不完整协方差张量

中已知非 零元素的集合；

表示投影到集合Ω的投影算子；||·||_*表示核范数；||·||_F表示Frobenius范数；λ表示正则化参数；

表示块Hankel矩阵操作，定义为将矩阵

变换到块Hankel矩阵

的操作，其中，u_p表示U的第p列，p＝1,2,…,P；S₁+S₂＝S+1；

表示Hankel 操作，定义为将向量

变换为Hankel矩阵

的操作，其中， S₁+S₂＝S+1。

(4)为求解具有约束的张量填充模型，将步骤(3)中的张量填充模型转换 无约束的增广拉格朗日函数形式。

将步骤(3)中的优化模型表示成无约束的增广拉格朗日函数形式：

其中，为了简洁表示，定义

和

β＞0表示惩罚系数；D⁽ⁿ⁾表示拉 格朗日乘子矩阵,n＝1,2,3,4；<·>表示内积。

(5)利用ADMM算法迭代求解步骤(4)中的增广拉格朗日函数，迭代结 束时获得因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾。

利用交替方向乘子(Alternating Direction Method of Multipliers，ADMM)算法迭代求解步骤(4)中的优化问题，其中第k次迭代更新变量

和

的步 骤如下：

(5.1)在第k次迭代中，通过求解优化问题

来更新变量

将该子问题转化为关于求解U⁽¹⁾和U⁽²⁾的凸优化问题：

式中，

(·)^T表示转置；R_(n)和

分 别为张量

和

的模式n展开；凸优化问题(3)的解满足：

其中，

表示

的逆变换；通过求解U⁽ⁿ⁾每一行的闭 式解来获得U⁽ⁿ⁾的解；令

和W_i＝diag(Ω⁽ⁿ⁾(i,:))，i∈{1,2,…,I_n}； 当n＝1时，I₁＝M；当n＝2时，I₂＝N；则式(4)能够转化为：

式中，U⁽ⁿ⁾(i,:)表示矩阵U⁽ⁿ⁾第i行；C(i,:)表示矩阵C中第i行。因此，行向量U⁽ⁿ⁾(i,:) 的闭式解表示为：

然后，通过求解U⁽ⁿ⁾的每一行

得到因子矩阵 U⁽ⁿ⁾＝[(U⁽ⁿ⁾(i,:))^T…(U⁽ⁿ⁾(I_n,:))^T]^T，n＝1,2。

(5.2)关于求解Δ＝{Δ⁽³⁾,Δ⁽⁴⁾}的凸优化问题表示为：

类似地，该凸优化问题的解满足：

式(8)能够进一步表示为：

式中，

为求解第i个方程组所得到的对角矩阵，n＝3,4,i＝1,2,…,I_n且当n＝3时，I₃＝M，当n＝4时，I₄＝N；因此

的闭式解能够表示为：

其中，

表示广义逆；为减小估计误差，对

取平均作为Δ⁽ⁿ⁾的解：

(5.3)关于变量

的优化问题表示为：

利用SVT算法求解矩阵核范数最小化问题，上述优化问题的解为：

其中，

为奇异值软阈值收缩算子，1/β_k为阈值。

(5.4)变量

的更新表达式为：

(5.5)更新β的表达式为β_k+1＝ρβ_k。

(6)根据因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算构造出完整 协方差张量

再通过对称厄米特展开还原为完整的协方差矩阵R，最后采用 ESPRIT算法进行目标角度估计。

根据步骤(5)得到的因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算 构造出完整协方差张量

并将

通过对称厄米特展开 为协方差矩阵，最后再利用ESPRIT算法估计出目标角度。

一种计算机存储介质，其上存储有计算机程序，该计算机程序被处理器执行 时实现上述的基于因子矩阵先验的阵元失效MIMO雷达角度估计方法。

一种计算机设备，包括储存器、处理器及存储在存储器上并可在处理器上运 行的计算机程序，所述处理器执行所述计算机程序时实现上述的基于因子矩阵先 验的阵元失效MIMO雷达角度估计方法。

有益效果：与现有技术相比，本发明具有如下优点：

1、阵元失效下双基地MIMO雷达虚拟阵列协方差矩阵具有整行和整列的数 据缺失。为利用回波数据的多维结构来恢复失效阵元的缺失数据，本发明将虚拟 阵列协方差矩阵重新排列为协方差张量，但协方差张量存在多个切片数据完全缺 失的问题，而传统的张量填充算法无法重建具有切片缺失的协方差张量。本发明 提出一种基于因子矩阵先验的张量填充算法，能够有效地恢复出四阶协方差张量 中多个切片缺失数据，提高阵元失效下MIMO雷达角度估计性能。

2、本发明挖掘了MIMO雷达的四阶协方差张量的多维结构特征、因子矩阵 的范德蒙德结构特点以及因子矩阵之间的关系，能更精确地恢复协方差张量中的 缺失数据，因此能更好地补偿因阵元失效而带来的角度估计性能损失。

3、本发明能有效解决MIMO雷达协方差张量中整个切片缺失的数据恢复问 题，具有较好的失效阵元缺失数据恢复精度以及较高的实时性。

附图说明

图1为本发明所述方法的步骤流程图；

图2为双基地MIMO雷达的目标角度估计星座图；其中，图2(a)中阵元失效 时采用ESPRIT算法，图2(b)中阵元失效时采用现有技术二，图2(c)中阵元失效 时采用本发明所述方法；

图3为本发明方法与现有技术二的收敛性变化曲线图；

图4为目标角度估计RMSE随信噪比的变化曲线图；

图5为目标角度估计RMSE随快拍数的变化曲线图；

图6为失效发射阵元数为1时角度估计RMSE随接收故障阵元数的变化曲线 图；

图7为失效发射阵元数为2时RMSE随接收故障阵元数的变化曲线图。

具体实施方式

下面结合附图对本发明的技术方案作进一步说明。

如图1所示，基于因子矩阵先验的阵元失效MIMO雷达角度估计方法包括以下步骤：

步骤1：对存在阵元失效的双基地MIMO雷达接收信号进行匹配滤波，得到 虚拟阵列的输出信号

并计算获得虚拟阵列协方差矩阵

双基地MIMO雷达的发射阵列具有M个发射阵元，其阵元间距为d_t，接收 阵列具有N个接收阵元，其阵元间距为d_r，且均为均匀线性阵列(ULA)。假设空 间远场存在P个非相干目标，第p(p＝1,2,…,P)个目标的波离方位角(DOD)为θ_p， 波达方位角(DOA)为

各发射阵元同时发射相互正交的编码信号 W＝[w₁,w₂,…,w_M]^T(即(1/K)WW^H＝I_M×M)，其中，

表示第m(m＝1,2,…,M)个 发射阵元在每个脉冲周期内的基带编码信号，K表示每个脉冲周期内的相位编码 个数，(·)^H表示共轭转置。接收阵列中的每个阵元接收目标反射回波，第q(q＝1,2,…,Q)个脉冲周期时的接收信号X_q为：

式中，

其中

为接收阵列导向矢量；

其中，

为 发射阵列导向矢量；

表示由矢量s_q构成的对角矩阵，其中

ρ_p表示第p个目标的反射系数，f_dp表示第p 个目标的多普勒频率，f_s为脉冲重复频率，(·)^T表示进行转置运算；

为 噪声矩阵。

利用发射波形的相互正交性，将(1/K)W^H右乘式(1)得到匹配滤波后的输出为：

其中，Y_q＝X_qW^H/K；

令y_q＝vec(Y_q)，Q个脉冲周期的回波 信号经匹配滤波后，虚拟阵列的输出为：

Y＝[y₁，y₂，…，y_Q]＝(A_r⊙A_t)S+Z (3)

式中，

则虚拟阵列协方差矩阵可表示为：

R＝YY^H/Q＝(A_r⊙A_t)R_s(A_r⊙A_t)^H+R_Z (4)

式中，

为信号协方差矩阵；

为噪声协方差 矩阵。对于非相干目标，R_s可表示为一个对角矩阵，即R_s＝diag([ρ₁,ρ₂，…,ρ_P])。

假设Ω_T和Ω_R分别为故障发射和接收阵元的位置集合，向量

和

指示在发射阵列和接收阵列中正常阵元的位置，其中，在故障阵元的索 引上元素为0，在其余位置上为1。因此，故障的发射和接收阵列导向矩阵可表示 为

将失效阵元对应的虚拟阵元输出数据置零，则阵元失效下的双基地MIMO雷 达虚拟阵列输出数据矩阵可表示为：

式中，

为经失效阵元对应的虚拟阵元输出数据置零后的噪声矩阵。

那么阵元失效下双基地MIMO雷达的虚拟阵列协方差矩阵可表示为：

式中，

为噪声协方差矩阵。

步骤2：利用虚拟阵列协方差矩阵

构造四阶协方差张量

将协方差张量

表示成张量CANDECOMP/PARAFAC(CP)分解模型。

双基地MIMO雷达的协方差矩阵R可重新排列成一个四阶协方差张量

可表示成：

式中，

为对角核心张量，其第(p,p,p,p)个元素等于ρ_p，其余为零；

为噪声协方差张量。

同样，阵元失效下双基地MIMO雷达的协方差矩阵

可重新排列成一个四阶 协方差张量

可表示成：

式中，

为阵元失效下的噪声协方差张量。分析式(11)可知，协方差张量

中 的切片数据缺失可表示为：

由(12)和(13)可以看出，在MIMO雷达中，失效阵元导致四阶协方差张量中一 些切片完全缺失。区别于大多数张量补全方法中的数据随机缺失的问题，我们将 协方差张量中的这些数据缺失形式称为结构性缺失。

忽略式(10)中的噪声项，协方差张量

的CP分解形式可表示为：

其中，

可以通过缩放因子矩阵使得g中的所有元素变 为1。此时可将g消除得到协方差张量的简化表达式：

MIMO雷达的发射和接收导向矩阵具有范德蒙德结构特点，而因子矩阵

和

的每一个列向量是用不同的常数对A_t和A_r的对应列向量进行缩放得到的，因 此，因子矩阵

和

也具有范德蒙德结构。

考虑到CP分解的置换和缩放效应，协方差张量

的因子矩阵可估计为：

式中，

和

表示因子矩阵；

表 示置换矩阵；Λ₁、Λ₂、Λ₃和Λ₄为P×P的实值对角矩阵，对角元素对应着缩放系 数；N₁、N₂、N₃和N₄表示拟合误差。忽略公式(16)中的拟合误差项，可以从(16a) 中得到：

将(17)代入到(16c)中得到：

其中，

为对角矩阵。类似地，U⁽²⁾、U⁽⁴⁾之间的关系可表示为：

步骤3：利用因子矩阵的范德蒙德结构特点、因子矩阵之间的相互关系以及 张量的低CP秩，建立具有因子矩阵先验约束的张量填充模型。

为有效恢复出失效阵元的缺失数据，利用因子矩阵的范德蒙德结构特点、因 子矩阵之间的相互关系以及张量的低CP秩，建立用于失效阵元缺失数据恢复的 张量填充模型：

其中，Ω表示不完整协方差张量

中已知非零元素的集合；

表示投影到 集合Ω的投影算子；||·||_*表示核范数；||·||_F表示Frobenius范数；λ表示正则化参 数；

表示块Hankel矩阵操作，定义为将矩阵

变换到块 Hankel矩阵

其中，S₁+S₂＝S+1；

表示 Hankel操作，定义为将向量

变换为Hankel矩阵

其中， S₁+S₂＝S+1，可表达为：

同时，

为

的逆操作，具体表示为：

其中，

为

的逆操作，可表示为：

根据

和

的定义，可以得到：

步骤4：为求解具有约束的张量填充模型，将步骤(3)中的张量填充模型转 换无约束的增广拉格朗日函数形式。

将式(20)的优化模型表示成增广拉格朗日函数形式：

其中，为了简洁表示，定义

和

β＞0表示惩罚系数；D⁽ⁿ⁾(n＝1,2,3,4)表示拉格朗日乘子矩 阵；<·>表示内积。

步骤5：利用ADMM算法迭代求解步骤4中的增广拉格朗日函数，迭代结束 时获得因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾。

利用ADMM(交替方向乘子，Alternating Direction Method of Multipliers)算法求解公式(25)，当迭代优化一个变量时，固定其他变量的最新值。其中，第k 次迭代更新变量

和

的步骤如下：

步骤5.1：在第k次迭代中，通过求解优化问题

来更新变 量

，将该子问题转化为关于U⁽ⁿ⁾的凸优化问题：

其中，

R_(n)和

分别为张量

和

的模式n展开。该凸优化问题的解满足：

式中，

表示

的逆变换。通过求解U⁽ⁿ⁾每一行的闭 式解来获得U⁽ⁿ⁾的解。定义

和W_i＝diag(Ω⁽ⁿ⁾(i,:))， i∈{1,2,…,I_n}，且当n＝1时，I₁＝M，当n＝2时，I₂＝N。则式(27)可以转化为：

式中，U⁽ⁿ⁾(i,:)表示矩阵U⁽ⁿ⁾第i行；C(i,:)表示矩阵C中第i行。因此，行向 量U⁽ⁿ⁾(i,:)的闭式解可表示为：

然后，可以通过求解U⁽ⁿ⁾的每一行

得到因子矩阵 U⁽ⁿ⁾＝[(U⁽ⁿ⁾(i,:))^T…(U⁽ⁿ⁾(I_n,:))^T]^T，n＝1,2。

步骤5.2：关于Δ的凸优化问题表示为：

类似地，该凸优化问题的解满足：

式(8)可进一步表示为：

其中，

为求解第i个方程组所得到的对角矩阵，n＝3,4；i＝1,2,…,I_n，且当n＝3时，I₃＝M；当n＝4时，I₄＝N。因此

的闭式解可表示为：

其中，

表示广义逆。为减小估计误差，对

取平均作为Δ⁽ⁿ⁾的解：

步骤5.3：关于变量

的优化问题表示为：

利用SVT算法求解矩阵核范数最小化问题，上述优化问题的闭式解为：

其中，

为奇异值软阈值收缩算子，1/β_k为阈值。

步骤5.4：变量D⁽ⁿ⁾的更新表达式为：

步骤5.5：更新β的表达式为β_k+1＝ρβ_k。

步骤6：根据因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算构造出 完整协方差张量

再通过对称厄米特展开还原为完整的协方差矩阵R，最后采 用ESPRIT算法进行目标角度估计。

根据步骤5得到的因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算构 造出完整协方差张量

并将

通过对称厄米特展开为 协方差矩阵，最后再利用ESPRIT算法估计出目标角度。

本发明的技术效果可通过以下仿真实验说明。将本发明方法与现有技术一(Jinli Chen,Tingxiao Zhang,Jiaqiang Li,et al.Joint Sensor Failure Detectionand Corrupted Covariance Matrix Recovery in Bistatic MIMO Radar With ImpairedArrays[J].IEEE Sensors Journal,2019,19(14):5834-5842)和现有技术二(J.Ying etal. Hankel Matrix Nuclear Norm Regularized Tensor Completion for N-dimensional Exponential Signals,IEEETransactions on Signal Processing,2017,65(14):3702-3717) 方法的DOA估计性能进行对比，并以阵元正常时直接采用ESPRIT算法进行角度 估计性能作为参照。需要说明的是，为了公平比较各种方法恢复失效阵元缺失数据的性能，均采用ESPRIT算法从不同方法重构出的完整数据中估计目标角度。 假设MIMO雷达的发射和接收阵元数分别为M＝7和N＝15，发射与接收阵元间距 均为信号波长的一半，各发射阵列发射相互正交的Hadamard编码信号，每个脉 冲重复周期内的编码个数为256，存在3个远场非相干目标

角度估计均方根误差定义为

其中，P为目标个数；M_T为蒙特卡 洛实验次数；

为第p个目标在第m_t次蒙特卡洛实验中的DOD估计值；

为 第p个目标在第m_t次蒙特卡洛实验中的DOA估计值；信噪比的定义为

以下仿真实验中取M_T＝100。本发明方法中，设置 ρ＝1.1，正则化参数λ＝20，β₀＝10^-1，当迭代条件满足

或者到达最大迭代次数k＝1000时，迭代停止，其中，

表示第k次迭代时协方 差张量估计值。

仿真实验1：阵元失效下不同方法的目标角度估计星座图

本实验设置信噪比SNR＝-5dB，快拍数Q＝100，蒙特卡洛实验次数M_T＝100， 星号表示角度估计值，交叉号表示角度真实值，发射阵列中第2和第5个阵元失 效，接收阵元第1，4,7,12，14个阵元失效。由图2(a)可以看出，在没有对失效阵 元进行补偿的情况下，ESPRIT算法失效。本发明方法和现有技术二都利用阵列 数据的多维结构来恢复失效阵元的缺失数据，如图2(b)所示，对于大多数蒙特卡 洛试验，现有技术二可以正确估计目标角度，但是某些试验中出现不正确的角度 估计结果。这是因为现有技术二在重构失效阵元缺失数据时只考虑了因子矩阵的 范德蒙德结构，而忽略因子矩阵之间的关系，导致在重构缺失数据时出现不稳定 的情况。从图2(c)可以看出，本发明方法在所有蒙特卡洛试验中能够准确地估计 目标角度，这主要得益于在张量填充中充分利用了因子矩阵的范德蒙德结构以及 因子矩阵之间的关系。

仿真实验2：本发明方法和现有技术二的收敛性分析

本次实验对本发明方法和现有技术二的收敛性进行了分析。图3为收敛速率 与迭代次数的关系曲线图，其中，仿真参数与仿真实验1相同。纵坐标为残差值

表示为对数形式。从图3可以看出，随着迭代次数的增 加，本发明方法的残差值相对变化下降的更快，收敛速率也更快，而现有技术二 需要更多的迭代次数才能达到收敛。现有技术二对因子矩阵进行随机初始化，收 敛速率较慢，而本发明方法通过对阵元失效下MIMO雷达的协方差张量进行CP 分解得到的因子矩阵作为初始化值，可以加快收敛速度，且可以降低重构误差。

仿真实验3：不同方法目标角度估计误差随信噪比的变化关系

本实验设置信噪比变化范围为-30～0dB，其余仿真参数与仿真实验1相同。 由图4可知，由于阵元失效破坏了协方差矩阵结构的完整性，因此阵元失效下直 接采用ESPRIT算法时角度估计精度较差，即无法有效估计目标的角度；由于对 失效阵元的缺失数据进行了恢复，现有技术一和现有技术二以及本发明方法具有 良好的估计性能，尤其在低信噪比区域，现有技术一和现有技术二以及本发明方 法与ESPRIT算法在正常情况下(即无阵元失效情况)几乎具有相同的性能，但总体 来看，本发明方法的性能明显优于现有技术一和现有技术二，特别是在高信噪比 的情况下。

仿真实验4：不同方法目标角度估计误差随快拍数的变化关系

本实验设置快拍数的变化范围为50～350，信噪比为-5dB，其余仿真参数与仿 真实验1相同。从图5可以看出，随着快拍数的增加，所有方法的角度估计精度 均有所提升，但本发明方法明显优于基于矩阵填充的方法即现有技术一，这表明 MIMO雷达的阵列数据多维结构对提高阵元失效的缺失数据的恢复精度具有非常 重要的作用。此外，由于在张量填充中增加了因子矩阵之间的约束关系，因此提 高了缺失数据恢复的稳定性，本发明方法比现有技术二具有更好的角度估计性能。

仿真实验5：不同方法角度估计误差随失效阵元数的变化关系

在本次实验中，考虑两种阵元失效情形:

情形1：发射阵列的第四个阵元失效，接收阵列的失效阵元数由1～8依次增 加，且每次失效阵元的位置均随机变化。

情形2：发射阵列的第二个和第五个阵元失效，接收阵列的失效阵元数由1～8 依次增加，且每次失效阵元的位置均随机变化。

在上述两种情形下，角度估计性能随接收失效阵元数的变化关系如图6和图 7所示，其中信噪比SNR＝-5dB，快拍数Q＝100。针对每种情形分别进行了100 次实验。从图6和图7可以看出，随着接收失效阵元数量的增加，每种方法的角 度估计性能都有不同程度的恶化，但是本发明方法在两种情况下均具有最好的角 度估计精度。

仿真实验6：不同方法的运行时间对比

表1不同DOA估计方法的运行时间

本实验仿真设置与仿真实验1相同，运行软件为MATLAB2018a，CPU为 IntelCorei5-4570，内存为8GB。由表1可知，相较于现有技术一和现有技术二， 本发明方法的运行时间更短，同时角度估计性能更优。

Claims (9)

1.一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，包括以下步骤：

(1)对存在阵元失效的双基地MIMO雷达接收信号进行匹配滤波，得到虚拟阵列的输出信号

并计算获得虚拟阵列协方差矩阵

(2)利用虚拟阵列协方差矩阵

构造四阶协方差张量

将协方差张量

表示成张量CANDECOMP/PARAFAC分解模型；

(3)利用因子矩阵的范德蒙德结构特点、因子矩阵之间的相互关系以及张量的低CP秩，建立具有因子矩阵先验约束的张量填充模型；

(4)为求解具有约束的张量填充模型，将步骤(3)中的张量填充模型转换无约束的增广拉格朗日函数形式；

(5)利用ADMM算法迭代求解步骤(4)中的增广拉格朗日函数，迭代结束时获得因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾；

(6)根据因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算构造出完整协方差张量

再通过对称厄米特展开还原为完整的协方差矩阵R，最后采用ESPRIT算法进行目标角度估计。

2.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(1)具体为：

(1.1)双基地MIMO雷达的发射阵列和接收阵列皆为均匀线性阵列，分别由M个发射阵元和N个接收阵元组成；对阵元失效下双基地MIMO雷达所采集到的回波信号进行匹配滤波，经过脉冲积累后，获得MN个虚拟阵元输出的回波信号矩阵为：

其中，

为在Q个快拍下虚拟阵列回波信号矩阵；

和

分别为存在失效阵元时的发射阵列和接收阵列流形矩阵；

为高斯白噪声矩阵；S∈^P×Q为目标系数矩阵，其中P为相互独立的远场目标；⊙表示Khatri-Rao积；

(1.2)Ω_T和Ω_R分别为失效发射阵元和失效接收阵元的位置集合，当发射阵列中第p_t∈Ω_T个阵元出现故障时，其流形矩阵

中的第p_t行全为零；当接收阵列中第p_r∈Ω_R个阵元出现故障时，流形矩阵

中的第p_r行全为零；虚拟阵列协方差矩阵在Q脉冲周期下的最大似然估计为：

其中，(·)^H表示共轭转置，R_s＝SS^H/Q＝diag(ρ)表示信号协方差矩阵，diag(ρ)表示由向量ρ生成的对角矩阵，ρ＝[ρ₁,ρ₂，…,ρ_P]，ρ_p,p＝1,2,…,P，表示第p个目标的反射系数，

表示噪声协方差矩阵。

3.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(2)具体为：

(2.1)双基地MIMO雷达虚拟阵列协方差矩阵

能够重新排列成一个四阶协方差张量

其中，(·)^*表示共轭；

表示Tucker运算；

(2.2)重新调整因子矩阵

使得g中的元素都为1，此时将g消除得到协方差张量的简化表达式

(2.3)由于CP分解的置换和缩放效应的扰动，协方差张量

的因子矩阵表示为

和

其中

和

表示四个因子矩阵；

表示置换矩阵；Λ₁、Λ₂、Λ₃和Λ₄为P×P的实值对角矩阵，对角元素对应着缩放系数；N₁、N₂、N₃和N₄表示拟合误差；

(2.4)由于MIMO雷达的发射和接收导向矩阵A_t、A_r均为范德蒙德结构，而因子矩阵

和

的每一个列向量是用不同的常数对A_t和A_r的对应列向量分别缩放得到的，因此当忽略拟合误差项影响时，矩阵

和

也具有范德蒙德结构，即因子矩阵U⁽ⁿ⁾,n＝1,…,4，具有范德蒙德结构；

(2.5)忽略拟合误差项的影响，得到因子矩阵之间的关系为

和

其中对角矩阵

对角矩阵

4.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(3)具体为：

为有效恢复出失效阵元的缺失数据，利用因子矩阵的范德蒙德结构特点、因子矩阵之间的相互关系以及张量的低CP秩，建立用于失效阵元缺失数据恢复的张量填充模型：

其中，M⁽ⁿ⁾,n＝1,2,3,4，表示辅助矩阵，Ω表示不完整协方差张量

中已知非零元素的集合；

表示投影到集合Ω的投影算子；||·||_*表示核范数；||·||_F表示Frobenius范数；λ表示正则化参数；

表示块Hankel矩阵操作，定义为将矩阵

变换到块Hankel矩阵

的操作，其中，u_p表示U的第p列，p＝1,2,…,P；S₁+S₂＝S+1；

表示Hankel操作，定义为将向量

变换为Hankel矩阵

的操作，其中，S₁+S₂＝S+1。

5.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(4)具体为：

将步骤(3)中的优化模型表示成无约束的增广拉格朗日函数形式：

其中，为了简洁表示，定义

Δ＝{Δ⁽³⁾,Δ⁽⁴⁾}、

和

β＞0表示惩罚系数；D⁽ⁿ⁾表示拉格朗日乘子矩阵,n＝1,2,3,4；<·>表示内积。

6.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(5)具体为：

利用交替方向乘子算法迭代求解步骤(4)中的优化问题，其中第k次迭代更新变量

Δ，

和

的步骤如下：

(5.1)在第k次迭代中，通过求解优化问题

来更新变量

将该子问题转化为关于求解U⁽¹⁾和U⁽²⁾的凸优化问题：

式中，

(·)^T表示转置；R_(n)和

分别为张量

和

的模式n展开；凸优化问题(3)的解满足：

其中，

表示

的逆变换；通过求解U⁽ⁿ⁾每一行的闭式解来获得U⁽ⁿ⁾的解；令

和W_i＝diag(Ω⁽ⁿ⁾(i,:))，i∈{1,2,…,I_n}；当n＝1时，I₁＝M；当n＝2时，I₂＝N；则式(4)能够转化为：

式中，U⁽ⁿ⁾(i,:)表示矩阵U⁽ⁿ⁾第i行；C(i,:)表示矩阵C中第i行。因此，行向量U⁽ⁿ⁾(i,:)的闭式解表示为：

然后，通过求解U⁽ⁿ⁾的每一行

得到因子矩阵

(5.2)关于求解Δ＝{Δ⁽³⁾,Δ⁽⁴⁾}的凸优化问题表示为：

类似地，该凸优化问题的解满足：

式(8)能够进一步表示为：

式中，

为求解第i个方程组所得到的对角矩阵，n＝3,4,i＝1,2,…,I_n且当n＝3时，I₃＝M，当n＝4时，I₄＝N；因此

的闭式解能够表示为：

其中，

表示广义逆；为减小估计误差，对

取平均作为Δ⁽ⁿ⁾的解：

(5.3)关于变量

的优化问题表示为：

利用SVT算法求解矩阵核范数最小化问题，上述优化问题的解为：

其中，

为奇异值软阈值收缩算子，1/β_k为阈值；

(5.4)变量

的更新表达式为：

(5.5)更新β的表达式为β_k+1＝ρβ_k。

7.根据权利要求1所述的一种基于因子矩阵先验的阵元失效MIMO雷达角度估计方法，其特征在于，所述步骤(6)具体为：

根据步骤(5)得到的因子矩阵U⁽¹⁾,U⁽²⁾和对角矩阵Δ⁽³⁾,Δ⁽⁴⁾，利用Tucker运算构造出完整协方差张量

并将

通过对称厄米特展开为协方差矩阵，最后再利用ESPRIT算法估计出目标角度。

8.一种计算机存储介质，其上存储有计算机程序，其特征在于，该计算机程序被处理器执行时实现如权利要求1-7中任一项所述的基于因子矩阵先验的阵元失效MIMO雷达角度估计方法。

9.一种计算机设备，包括储存器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述计算机程序时实现如权利要求1-7中任一项所述的基于因子矩阵先验的阵元失效MIMO雷达角度估计方法。

Images