CN103217671B - 色噪声环境下的多输入多输出雷达收发角度快速估计方法 - Google Patents

Abstract

本发明提供的是一种色噪声环境下的多输入多输出雷达收发角度快速估计方法。1）将A/D采样模块采集的各接收阵元的回波信号进行匹配滤波，得到虚拟阵列的接收数据；2）将虚拟阵列的接收数据进行重排；3）对接收数据和重排后的矩阵分别进行子阵划分，满足发射端和接收端的旋转不变特性；4）构造出两个四阶累积量矩阵对；5）分别通过特征值分解获得包含发射和接收旋转不变特性的信号子空间；6）利用ESPRIT算法估计出目标的收发角度；7）利用最大似然估计方法实现收发角度的正确配对。本发明可有效抑制高斯色噪声的影响，对发射/接收阵元数无特殊要求，其估计精度高、计算复杂度低，可用于跟踪和制导中对海面或低空目标的定位。

Description

色噪声环境下的多输入多输出雷达收发角度快速估计方法

技术领域

本发明涉及的是一种双基地MIMO雷达多目标收发角度估计方法。

背景技术

多输入多输出（MIMO）雷达是借鉴通信领域的多输入多输出技术而提出一种具有高探测性能的新体制雷达。与传统的相控阵雷达相比，MIMO雷达在发射端发射相互正交的波形，在接收端可以通过匹配滤波形成一个很大的虚拟阵列孔径，提高了探测目标的空间分辨率同时也增加了雷达的自由度。目前MIMO雷达主要分为以下两种，一种是共址MIMO雷达，也即单/双基地MIMO雷达，该雷达的发射端和接收端的阵元均为集中式的分布，通过不同发射天线发射正交信号在接收端匹配形成一个很长的虚拟阵列，这些虚拟阵元可以形成一个低旁瓣的窄带波束，提高目标角分辨率和角度估计的精度，获得很好的参数估计性能；同时增加雷达的自由度，可以探测更多的目标。另一种为统计MIMO雷达，该MIMO雷达的发射端和接收端的阵元均为宽分布式的排列方式，可以获得空间分集特性抑制目标闪烁的缺点。同时，MIMO雷达在抗干扰、载波抑制和低截获率方面均具有很大的优势。

在MIMO雷达系统中，目标相对于发射阵列的波离方向(DOD)和相对于接收阵列的波达方向(DOA)为多目标定位时的重要估计参数。目前对多目标的收发角度估计主要集中于假设背景噪声为白噪声的情况进行研究的：C Duofang等在Electronics Letters期刊的2008年第44卷第12期的770页至771页提出了一种ESPRIT算法，该方法将ESPRIT算法同时应用在双基地MIMO雷达的发射端和接收端，利用信号子空间的旋转不变性实现了收发角度的正确配对，但需要进行额外配对。陈金立等在电子与信息学报的2009年第31卷第7期的1664页至1668页提出了一种双基地MIMO雷达快速多目标定位方法，该方法利用两次一维ESPRIT算法之间的关系实现了收发角度的自动配对。但是背景噪声为高斯色噪声时，以上这些算法的角度估计性能会严重下降，有时甚至会发生错误估计。针对色噪声的环境，Jin Ming等在SignalProcessing期刊的2009年第89卷第2期的244页至251页提出了一种在发射端利用三天线的ESPRIT算法消除色噪声的影响，刘晓莉等在电子与信息学报的2010年第32卷第9期的2179页至2183页提出了一种在接收端利用双天线的MUSIC-ESPRIT算法消除色噪声的影响。但上述方法分别对发射阵元和接收阵元的数目有特殊要求，不具有实用性。

发明内容

本发明的目的在于提供一种可以有效抑制高斯色噪声的影响，提高目标收发角度的估计性能，降低运算复杂度，同时对发射/接收阵元的数目没有特殊要求，更具有实用性的色噪声环境下的多输入多输出雷达收发角度快速估计方法。

本发明的目的是这样实现的：

（1）双基地MIMO雷达系统由M个发射阵元和N个接收阵元组成。M个发射阵元发射相互正交的波形，在接收端通过A/D采样模块采集每个接收阵元的回波信号，将每个接收阵元的采集信号通过一个匹配滤波器组，分离出相应的各个发射通道的信号，得到一个阵元数为NM的虚拟阵列的接收数据Y_r(t)；

（2）对虚拟阵列的接收数据Y_r(t)进行重排，得到一新矩阵Y_t(t)；具体方法为：

将Y_r(t)＝[y₁,y₂,...,y_M,y_M+1,y_M+2,...,y_2M,...,y_(N-1)M+1,y(_N-1)M+2,...,y_NM]^T进行有限次的行对调变换，得到一个与Y_r(t)等价的新矩阵Y_t(t)，

Y_t(t)＝[y₁,y_M+1,...,y_(N-1)M+1,y₂,y_M+2,...,y_(N-1)M+2,...,y_M,y_2M,...,y_NM]^T     (1)

＝C_tS_t(t)+N_t(t)

其中：y_i为接收数据Y_r(t)＝C_rS_r(t)+N_r(t)的第i个列矢量(i＝1,2,...,NM)，C_r为收发角度的联合导向矢量矩阵，S_r(t)是P个目标散射回来的信号矩阵，N_r(t)为噪声矩阵；C_t、S_t(t)和N_t(t)分别为重排后的联合导向矢量矩阵、信号矩阵和噪声矩阵，(·)^T表示向量或矩阵的转置运算。

同时矩阵Y_r(t)和Y_t(t)两者之间的等价关系表示为

Y_t(t)＝JY_r(t)=JC_rS_r(t)+JN_r(t)                     (2)

＝C_tS_t(t)+N_t(t)

其中：J为MN×MN维转换矩阵，对应的有限次对调变换可表示为其中 和a(θ_p)分别为第p个目标的发射和接收导向矢量，和θ分别为第p个目标的发射角和接收角，N_t(t)=JN_r(t)。

（3）将Y_r(t)和Y_t(t)分别进行划分，获得发射端和接收端的旋转不变特性；具体方法为：

分别取Y_r(t)的前M(N-1)行和后M(N-1)行的输出数据构造成M(N-1)×P维矩阵，记为Y_r1(t)和Y_r2(t)；分别取Y_t(t)的前N(M-1)行和后N(M-1)行的输出数据构造成N(M-1)×P维矩阵，记为Y_t1(t)和Y_t2(t)：

$\left\{\begin{array}{c}{Y}_{r1}\left(t\right)=C_{r1}{S}_{r1}\left(t\right)+N_{r1}\left(t\right)\\ Y_{r2}\left(t\right)=C_{r1}{\mathrm{\Φ}}_r{S}_{r2}\left(t\right)+N_{r2}\left(t\right)\end{array}\right.---\left(3\right)$

$\left\{\begin{array}{c}{Y}_{t1}\left(t\right)=C_{t1}{S}_{t1}\left(t\right)+N_{t1}\left(t\right)\\ Y_{t2}\left(t\right)=C_{t1}{\mathrm{\Φ}}_t{S}_{t2}\left(t\right)+N_{t2}\left(t\right)\end{array}\right.---\left(4\right)$

其中：为包含了所有发射方向信息的旋转对角矩阵，Φ_r=diag[exp(-j2πd_rsinθ₁/λ),exp(-j2πd_rsinθ₂/λ),...,exp(-j2πd_rsinθ_P/λ)]为包含了所有接收方向信息的旋转对角矩阵，其中λ为载波波长，d_t和d_r分别为发射和接收阵元间距；S_r1(t)和S_r2(t)分别为S_r(t)的前M(N-1)行和后M(N-1)行构成的矩阵；S_t1(t)和S_t2(t)分别为重排后的信号矢量S_t(t)的前N(M-1)行和后N(M-1)行构成的矩阵；N_r1(t)、N_r2(t)∈C^M(N-1)×T和N_t1(t)、N_t2(t)∈C^N(M-1)×T为噪声矩阵，T为采样拍数。

（4）根据四阶累积量的基本定义，构造出四阶累积量矩阵对：

假设空间中存在N维零均值复平稳随机过程X，四阶累积量定义式为

$C_{4x}(p_1,p_2,p_3,p_4)=\mathrm{cum}\{x_{p_1}^{*},x_{p_2},x_{p_3},x_{p_4}^{*}\}$

$=E\{x_{p_1}^{*},x_{p_2},x_{p_3},x_{p_4}^{*}\}-E\{x_{p_1}^{*},x_{p_3}\}E\{x_{p_2},x_{p_4}^{*}\}----\left(5\right)$

$E\{x_{p_1}^{*},x_{p_2}\}E\{x_{p_3},x_{p_4}^{*}\}-E(x_{p_1}^{*},x_{p_4}^{*})E\{x_{p_2},x_{p_3}\}$

式中： $X={[x_{p_1},x_{p_2},...,x_{p_N}]}^T,$ $E\{x_{p_1},x_{p_2},x_{p_3},x_{p_4}\}$ 和E{x_i,x_j}分别为X的四阶矩和二阶矩，其中1≤p₁,p₂,p₃,p₄≤N，(·)^*表示向量或矩阵的共轭运算。

根据四阶累积量的基本定义，构造出四阶累积量矩阵对：

$\left\{\begin{array}{c}{R}_{r1}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r1}^H\}\\ R_{r2}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r2}^H\}\end{array}\right.---\left(6\right)$

$\left\{\begin{array}{c}{R}_{t1}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{t1}_1},y_{t1}^H\}\\ R_{t2}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{\mathrm{tr}1}_1},y_{t2}^H\}\end{array}\right.---\left(7\right)$

式中：y_r1和y_t1分别为矩阵Y_r1和Y_t1的第1个列矢量，y_r2和y_t2分别为矩阵Y_r1和Y_t1的第2个列矢量，和分别为y_r1和y_t1的第一个元素；为四阶累积量矩阵的计算表达式，(·)^H表示向量或矩阵的共轭转置运算。

将式(6)中矩阵R_r1的第(i,j)个元素进一步展开，则有

由四阶累积量的性质可知，高斯色噪声的四阶累积量恒等于零，即 从而证明了四阶累积量可有效抑制高斯色噪声的影响。又由于P个目标相互不相关，则有以下等式成立：

将式(9)带入到式(8)可得：

同理，式(6)中的R_r2可整理为：

最后将式(10)和式(11)经过化简写成矩阵的形式：

$R_{r1}=C_{r1}{R}_{4s}{C}_{r1}^H---\left(12\right)$

$R_{r2}=C_{r1}{R}_{4s}{\mathrm{\Φ}}_r^H{C}_{r1}^H---\left(13\right)$

同理，将式(7)展开与化简可有：

$R_{t1}=C_{t1}{R}_{4s}{C}_{t1}^H---\left(14\right)$

$R_{t2}=C_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H{C}_{t1}^H---\left(15\right)$

（5）对四阶累积量矩阵R_t1和R_r1进行特征值分解：

$R_{t1}=U_{\mathrm{ts}}{\mathrm{\Σ}}_{\mathrm{ts}}{U}_{\mathrm{ts}}^H+U_{\mathrm{tn}}{\mathrm{\Σ}}_{\mathrm{tn}}{U}_{\mathrm{tn}}^H---\left(16\right)$

$R_{r1}=U_{\mathrm{rs}}{\mathrm{\Σ}}_{\mathrm{rs}}{U}_{\mathrm{rs}}^H+U_{\mathrm{rn}}{\mathrm{\Σ}}_{\mathrm{rn}}{U}_{\mathrm{rn}}^H---\left(17\right)$

式中：Σ_ts和Σ_tn分别为R_t1进行特征值分解后得到的P个较大特征值构成的P×P维对角矩阵和其余NM-P个特征值构成的对角矩阵，Σ_rs和Σ_rn分别为R_r1进行特征值分解后得到的P个较大特征值构成的P×P维对角矩阵和其余NM-P个特征值构成的对角矩阵；U_ts和U_rs分别是与R_t1和R_r1的P个较大特征值对应特征矢量构成的NM×P维矩阵，U_tn和U_rn分别为是对应剩下的特征矢量构成的矩阵，其中U_ts和U_rs分别为包含了发射和接收旋转不变特性的信号子空间，P为目标数。

（6）根据信号子空间，利用ESPRIT算法分别得到目标的收发角度。

根据ESPRIT算法的基本原理可知，存在一个满值矩阵T_t，使得

U_st=C_t1T_t                   (18)

那么就是P×P维的满秩矩阵，则存在唯一的P×P维线性算子Ψ_t，使得下式成立：

${\mathrm{\Ψ}}_t{U}_{\mathrm{st}}^H{C}_{t1}{R}_{4s}=U_{\mathrm{st}}^H{C}_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H---\left(19\right)$

将上式等号左右两侧同时乘以C_t1U_st，则有

${\mathrm{\Ψ}}_t{U}_{\mathrm{st}}^H{C}_{t1}{R}_{4s}{C}_{t1}{U}_{\mathrm{st}}=U_{\mathrm{st}}^H{C}_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H{C}_{t1}{U}_{\mathrm{st}}---\left(20\right)$

即

${\mathrm{\Ψ}}_t{U}_{\mathrm{st}}^H{R}_{t1}{U}_{\mathrm{st}}=U_{\mathrm{st}}^H{R}_{t2}{U}_{\mathrm{st}}---\left(21\right)$

则可以得到线性算子Ψ_t：

${\mathrm{\Ψ}}_t=\left(U_{\mathrm{st}}^H{R}_{t2}{U}_{\mathrm{st}}\right){\left(U_{\mathrm{st}}^H{R}_{t1}{U}_{\mathrm{st}}\right)}^{-1}---\left(22\right)$

对上式进行特征值分解，即可估计出第p个目标的发射角：

式中：为对角矩阵Ψ_t的第p个对角元素。

同理，存在一个满秩矩阵T_r，使得U_sr=C_r1T_r成立，最后整理得到

${\mathrm{\Ψ}}_r=\left(U_{\mathrm{sr}}^H{R}_{r2}{U}_{\mathrm{sr}}\right){\left(U_{\mathrm{sr}}^H{R}_{r1}{U}_{\mathrm{sr}}\right)}^{-1}---\left(24\right)$

对其进行特征值分解，即可估计出第p个目标的接收角：

${\widehat{\mathrm{\θ}}}_p=-\arcsin (\mathrm{\λ}\∠\left({\widehat{\mathrm{\β}}}_p^{*}\right)/2\mathrm{\π}{d}_r)---\left(25\right)$

式中：为对角矩阵Ψ_r的第p个对角元素。

（7）利用最大似然估计方法，实现收发角度的正确配对。

利用最大似然估计算法对角度进行匹配，将P²种不同收发角度组合方式代入到下式中，其中P个使得代价函数Θ取值最小值的角度组合方式，即为正确的收发角度配对方式。其中代价函数为

式中：为与张成的子空间相互正交的投影矩阵，其中为联合导向矢量矩阵的第p个列矢量(p=1,2,...,P)，I_K为K×K维的单位阵，(·)^#表示向量或矩阵的违逆运算。

本发明提供了一种色噪声环境下的MIMO雷达收发角度快速估计方法，可以有效抑制高斯色噪声的影响，提高目标收发角度的估计性能，进一步降低了运算复杂度，同时对发射/接收阵元的数目没有特殊要求，更具有实用性。

本发明的技术思路为：双基地MIMO雷达的收发阵列均为半波长的线性等距阵列，各个发射阵元发射相互正交的波形；在接收端通过A/D采样模块采集每个接收阵元的回波信号，将每个接收阵元的采集信号通过一个匹配滤波器组，分离出相应的各个发射通道的信号，得到一个阵元数为MN的虚拟阵列的接收数据；对虚拟数据矩阵进行重排和划分，分别满足发射端和接收端的旋转不变特性；从四阶累积量的基本定义出发，构造出四阶累积量矩阵对；对四阶累积量矩阵R_t1和R_r1进行特征值分解，分别获得包含发射和接收旋转不变特性的信号子空间；利用ESPRIT算法分别估计出目标的收发角度；最后利用最大似然估计方法，实现收发角度的正确配对。

本发明与现有技术相比较有以下特点：

（1）现有的技术大部分都是针对噪声背景为白噪声的情况对目标的收发角度进行估计的，而当噪声为高斯色噪声时，其角度估计性能会严重下降，有时甚至会发生错误估计；而本发明通过利用四阶累积量的盲高斯特性，有效抑制了高斯色噪声的影响，实现了收发角度的正确估计，提高了收发角度的估计精度。

（2）本发明避免了高斯色噪声情况下Jin Ming等提出的在发射端利用三天线的ESPRIT算法和刘晓莉等提出的在接收端利用双天线的MUSIC-ESPRIT算法对发射和接收阵元的数目具有特殊要求的缺点，从而本发明更有利于应用到实际应用中。

（3）本发明避免了MUSIC算法的空间谱搜索过程，降低了运算复杂度，更有利于实时处理和硬件上的实现。

附图说明

图1是本发明的实现流程图；

图2是本发明目标定位的整理流程图；

图3是本发明的多目标收发角度配对结果图；

图4是本发明与其他算法的角度估计均方根误差与信噪比的变化曲线图；

图5是本发明与其他算法的角度估计均方根误差与采样拍数的变化曲线图；

图6是本发明在接收阵元数N不变时的角度估计均方根误差与信噪比的变化曲线图；

图7是本发明在发射阵元数M不变时的角度估计均方根误差与信噪比的变化曲线图。

具体实施方式

下面结合实现流程图，举例对本发明做更加具体的描述。

步骤一、接收端通过A/D采样模块采集每个接收阵元的回波信号

假设双基地MIMO雷达是由M个发射阵元和N个接收阵元组成，并且发射和接收阵列都是半波长的均匀等距线阵，所有阵元都是全向的，d_t和d_r分别为发射和接收阵元间距，则在接收端通过A/D采样模块采集每个接收阵元的回波信号x(l,t)可表示为

式中：和θ分别为第p个目标的发射角和接收角；a(θ)=[1,exp(-j2π(d_r/λ)sinθ),...,-j2π(d_r/λ)(N-1)sinθ)]^T为N×1维的接收导向矢量，λ为载波波长，为M×1维的发射导向矢量，β为目标散射系数，其与目标雷达的横截面积有关，f_d为目标信号的多普勒频率；S=[s₁(l)…s_M(l)]是发射的正交窄带信号，l代表发射信号的脉冲时间，t表示雷达脉冲数，w(l,t)为高斯色噪声，(·)^T表示矢量或者矩阵的转置运算。

步骤二、将每个接收阵元的采集信号通过一个匹配滤波器组，分离出相应的各个发射通道的信号，得到一个阵元数为MN的虚拟阵列的接收数据Y_r(t)

Y_r(t)=C_rS_r(t)+N_r(t)                  (32)

式中：是MN×P维的收发角度联合导向矢量矩阵，其中 代表克罗内克积； $S_r\left(t\right)={[{\mathrm{\β}}_1{e}^{j2\mathrm{\π}{f}_{d_1}t},...,{\mathrm{\β}}_P{e}^{j2\mathrm{\π}{f}_{d_P}t}]}^T$ 是P个目标散射回来的信号，N_r(t)为通过匹配滤波器组后的噪声矢量。

步骤三、对虚拟接收数据Y_r(t)进行重排，得到一新矩阵Y_t(t)

将Y_r(t)=[y₁,y₂,...,y_M,y_M+1,y_M+2,...,y_2M,...,y_(N-1)M+1,y_(N-1)M+2,...,y_NM]^T进行有限次的行对调变换，得到一个与Y_r(t)等价的新矩阵Y_t(t)：

Y_t(t)=[y₁,y_M+1,...,y_(N-1)M+1,y₂,y_M+2,...,y_(N-1)M+2,...,y_M,y_2M,...,y_NM]^T      (33)

式中：y_i为接收数据Y_r(t)的第i个列矢量(i=1,2,...,NM)，(·)^T表示向量或矩阵的转置运算。

步骤四、将虚拟阵列的接收数据Y_r(t)和重排后的矩阵Y_t(t)分别进行子阵划分，获得发射端和接收端的旋转不变特性；

分别取Y_r(t)的前M(N-1)行和后M(N-1)行的输出数据构造成M(N-1)×P维矩阵，记为Y_r1(t)和Y_r2(t)；分别取Y_t(t)的前N(M-1)行和后N(M-1)行的输出数据构造成N(M-1)×P维矩阵，记为Y_t1(t)和Y_t2(t)：

$\left\{\begin{array}{c}{Y}_{r1}\left(t\right)=C_{r1}{S}_{r1}\left(t\right)+N_{r1}\left(t\right)\\ Y_{r2}\left(t\right)=C_{r1}{\mathrm{\Φ}}_r{S}_{r2}\left(t\right)+N_{r2}\left(t\right)\end{array}\right.---\left(34\right)$

$\left\{\begin{array}{c}{Y}_{t1}\left(t\right)=C_{t1}{S}_{t1}\left(t\right)+N_{t1}\left(t\right)\\ Y_{t2}\left(t\right)=C_{t1}{\mathrm{\Φ}}_t{S}_{t2}\left(t\right)+N_{t2}\left(t\right)\end{array}\right.---\left(35\right)$

式中：为包含了所有发射方向信息的旋转对角矩阵，Φ_r=diag[exp(-j2πd_rsinθ₁/λ),exp(-j2πd_rsinθ₂/λ),...,exp(-j2πd_rsinθ_P/λ)]为包含了所有接收方向信息的旋转对角矩阵；S_r1(t)和S_r2(t)分别为S_r(t)的前M(N-1)行和后M(N-1)行构成的矩阵；S_t1(t)和S_t2(t)分别为重排后的信号矢量S_t(t)的前N(M-1)行和后N(M-1)行构成的矩阵；N_r1(t)、N_r2(t)∈C^M(N-1)×T和N_t1(t)、N_t2(t)∈C^N(M-1)×T为噪声矩阵，T为采样拍数。

步骤五、根据四阶累积量的基本定义，构造出两个四阶累积量矩阵对

$\left\{\begin{array}{c}{R}_{r1}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r1}^H\}=C_{r1}{R}_{4s}{C}_{r1}^H\\ R_{r2}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r2}^H\}=C_{r1}{R}_{4s}{\mathrm{\Φ}}_r^H{C}_{r1}^H\end{array}\right.$

$\left\{\begin{array}{c}{R}_{t1}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{t1}_1},y_{t1}^H\}=C_{t1}{R}_{4s}{C}_{t1}^H\\ R_{t2}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{\mathrm{tr}1}_1},y_{t2}^H\}=C_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H{C}_{t1}^H\end{array}\right.---\left(37\right)$

式中：y_r1和y_t1分别为矩阵Y_r1和Y_t1的第1个列矢量，y_r2和y_t2分别为矩阵Y_r1和Y_t1的第2个列矢量，和分别为y_r1和y_t1的第一个元素；为四阶累积量矩阵的计算表达式，(·)^*表示向量或矩阵的共轭运算，(·)^H表示向量或矩阵的共轭转置运算。

步骤六、对四阶累积量矩阵对R_t1和R_r1进行特征值分解，分别获得包含发射和接收旋转不变特性的信号子空间

$R_{t_1}=U_{\mathrm{ts}}{\mathrm{\Σ}}_{\mathrm{ts}}{U}_{\mathrm{ts}}^H+U_{\mathrm{tn}}{\mathrm{\Σ}}_{\mathrm{tn}}{U}_{\mathrm{tn}}^H---\left(38\right)$

$R_{r_1}=U_{\mathrm{rs}}{\mathrm{\Σ}}_{\mathrm{rs}}{U}_{\mathrm{rs}}^H+U_{\mathrm{rn}}{\mathrm{\Σ}}_{\mathrm{rn}}{U}_{\mathrm{rn}}^H---\left(39\right)$

式中：Σ_ts和Σ_tn分别为R_t1进行特征值分解后得到的P个较大特征值构成的P×P维对角矩阵和其余NM-P个特征值构成的对角矩阵，Σ_rs和Σ_rn分别为R_r1进行特征值分解后得到的P个较大特征值构成的P×P维对角矩阵和其余NM-P个特征值构成的对角矩阵；U_ts和U_rs分别是与R_t1和R_r1的P个较大特征值对应特征矢量构成的NM×P维矩阵，U_tn和U_rn分别为是对应剩下的特征矢量构成的矩阵，其中U_ts和U_rs分别为包含了发射和接收旋转不变特性的信号子空间，P为目标数。

步骤七、根据信号子空间，利用ESPRIT算法分别得到目标的收发角度

根据ESPRIT算法的基本原理可知，存在一个满值矩阵T_t，使得

U_st=C_t1T_t                             (40)

那么就是P×P维的满秩矩阵，则存在唯一的P×P维的线性算子Ψ_t，使得 ${\mathrm{\ψ}}_t{U}_{\mathrm{st}}^H{C}_{t1}{R}_{4s}=U_{\mathrm{st}}^H{C}_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H$ 成立，然后经整理与化简可有：

${\mathrm{\Ψ}}_t=\left(U_{\mathrm{st}}^H{R}_{t2}{U}_{\mathrm{st}}\right){\left(U_{\mathrm{st}}^H{R}_{t1}{U}_{\mathrm{st}}\right)}^{-1}---\left(41\right)$

对其进行特征值分解，即可估计出第p个目标的发射角：

式中：为对角矩阵Ψ_t的第p个对角元素。

同理，存在一个满秩矩阵T_r，使得U_sr=C_r1T_r成立，最后得到

${\mathrm{\Ψ}}_r=\left(U_{\mathrm{sr}}^H{R}_{r2}{U}_{\mathrm{sr}}\right){\left(U_{\mathrm{sr}}^H{R}_{r1}{U}_{\mathrm{sr}}\right)}^{-1}---\left(43\right)$

对其进行特征值分解，即可估计出第p个目标的接收角：

${\widehat{\mathrm{\θ}}}_p=-\arcsin (\mathrm{\λ}\∠\left({\widehat{\mathrm{\β}}}_p^{*}\right)/2\mathrm{\π}{d}_r)---\left(44\right)$

式中：为对角矩阵Ψ_r的第p个对角元素。

步骤七、利用最大似然估计方法，实现收发角度的正确配对

将P²种不同的角度组合方式代入到下式中，其中P种使得代价函数Θ取值最小值的收发角度组合方式，即为正确的收发角度配对方式。其代价函数为

式中：为与张成的子空间相互正交的投影矩阵，其中为联合导向矢量矩阵的第p个列矢量(p=1,2,...,P)，I_K为K×K维的单位阵，(·)^#表示向量或矩阵的违逆运算。

本发明的效果可通过以下仿真说明：

（一）仿真条件与内容：

试验一、本发明对多目标收发角度的配对性能

假设双基地MIMO雷达的发射阵元数M=4，接收阵元数N=6，收发阵元间距均为信号的半个波长。空间中存在3个互不相关的目标，且各目标的位置分别为 三个目标的信噪比均为SNR=10dB，采样拍数T=500，对三个目标的发射角和接收角的配对问题进行测试，独立进行200次Monte-Carlo试验。

试验二、多目标角度估计均方根误差与信噪比的变化关系

为了与刘晓莉提出的在接收端利用双天线的MUSIC-ESPRIT算法相比较，我们假设双基地MIMO雷达的发射阵元数M=8，接收阵元数N=2，收发阵元间距均为信号的半个波长。空间中存在3个互不相关的目标，且各目标的位置分别为 采样拍数T=500。定义目标的均方根误差为其中和分别为第l次Monte-Carlo仿真实验的第i个目标的发射角度和接收角度估计值，和θ_i分别为第i个目标的发射角和接收角的真实值。P和L分别为目标数目和Monte-Carlo仿真实验的次数，这里独立进行200次Monte-Carlo试验。这里采用本发明、传统的ESPRIT算法和刘晓莉提出的MUSIC-ESPRIT算法进行仿真比较。

试验三、多目标角度估计均方根误差与采样拍数的变化关系

三个目标的信噪比均为SNR=10dB，其他条件与试验二相同，这里采用本发明、传统的ESPRIT算法和刘晓莉的MUSIC-ESPRIT算法进行仿真比较。

试验四、不同阵元数的MIMO雷达对多目标的角度估计性能

三个目标的采样拍数均为T=500，其他条件与试验二相同，这里选择不同的阵元数目对多个目标进行角度估计，独立进行200次Monte-Carlo试验。

(二)仿真结果

试验一、本发明对多目标收发角度的配对性能

图3为本发明的多目标收发角度配对结果图，从图中可知，本发明可以正确估计出双基地MIMO雷达的发射角和接收角，并且可实现收发角度的正确配对，从而验证了本发明对多个目标收发角度估计的有效性。

试验二、多目标角度估计均方根误差与信噪比的变化关系

图4为本发明与其他两种算法的角度估计均方根误差与信噪比的变化曲线图。从图中可知，本发明的收发角度估计性能比其他两种算法更优越，这是因为MUSIC-ESPRIT算法首先利用一维谱峰搜索估计出双基地MIMO雷达的发射角，在根据收发角度的关系利用ESPRIT算法估计出目标的接收角；但该算法角度的估计精度受搜索步长的影响，如果希望估计出的角度越精确，就要求搜索步长越小，但运算量越大，这里进行仿真时设置角度的估计精度为0.01。而本发明提出的算法分别构造出四阶累积量矩阵，然后利用ESPRIT算法估计出MIMO雷达的收发角度，本发明利用四阶累积量的盲高斯特性，有效地抑制了高斯色噪声的影响，提高了目标的角度估计性能；并且本发明所提算法对发射或接收阵元的数目没有特殊的要求，更具有实用性，同时避免了谱峰搜索过程，降低了运算复杂度，算法的运算时间更短，更加有利于信号的实时处理和硬件实现。

试验三、多目标角度估计均方根误差与采样拍数的变化关系

图5为本发明与其他两种算法的角度估计均方根误差与采样拍数的变化曲线图。从图中可知，随着采样拍数的变化，本发明比传统ESPRIT算法和MUSIC-ESPRIT算法具有更好的角度估计性能，且随着采样拍数的增加，均方根误差呈递减的趋势。

试验四、不同阵元数的MIMO雷达对多目标的估计性能

图6和图7分别为本发明在不同阵元数的情况下角度估计均方根误差与信噪比的变化曲线图。从图中可知，无论是发射阵元数还是接收阵元数的增加，本发明的角度估计性能都会有所提高，可更好的抑制高斯色噪声的影响，这是因为MIMO雷达的阵元数越多，其获得的分集增益越大。因此，对目标进行探测时，选择合适的采样拍数和合理的阵元数目，都能够进一步提高本发明的角度估计性能。

Claims (7)

1.一种色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是包括如下过程：

(1)M个发射阵元同时发射相互正交的窄带信号，N个接收阵元接收回波信号，并且发射和接收阵列都是半波长均匀线阵，所有阵元都是全向的；

(2)通过A/D采样模块采集每个接收阵元的回波信号，将每个接收阵元的采集信号通过一个匹配滤波器组，分离出相应的各个发射通道的信号，得到一个阵元数为NM的虚拟阵列的接收数据Y_r(t)；

(3)将虚拟阵列的接收数据进行重排，得到一新矩阵Y_t(t)；

(4)将虚拟阵列的接收数据Y_r(t)和重排后的矩阵Y_t(t)分别进行子阵划分，获得发射端和接收端的旋转不变特性；

(5)根据四阶累积量的基本定义，构造出两个四阶累积量矩阵对；

(6)对四阶累积量矩阵进行特征值分解，分别获得包含发射和接收旋转不变特性的信号子空间；

(7)根据信号子空间，利用ESPRIT算法分别估计出目标的收发角度；

(8)利用最大似然估计方法，实现收发角度的正确配对。

2.根据权利要求1所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的对虚拟阵列的接收数据Y_r(t)进行重排，得到一新矩阵Y_t(t)，按如下方式进行：

将Y_r(t)＝[y₁,y₂,...,y_M,y_M+1,y_M+2,...,y_2M,...,y_(N-1)M+1,y_(N-1)M+2,...,y_NM]^T进行有限次的行对调变换，得到一个与Y_r(t)等价的新矩阵Y_t(t):

Y_t(t)＝[y₁,y_M+1,...,y_(N-1)M+1,y₂,y_M+2,...,y_(N-1)M+2,...,y_M,y_2M,...,y_NM]^T

＝C_tS_t(t)+N_t(t)

其中：y_i为接收数据Y_r(t)＝C_rS_r(t)+N_r(t)的第i个列矢量，i＝1,2,...,NM，C_r为收发角度的联合导向矢量矩阵，S_r(t)是P个目标散射回来的信号矩阵，N_r(t)为噪声矩阵；C_t、S_t(t)和N_t(t)分别为重排后的联合导向矢量矩阵、信号矩阵和噪声矩阵，(·)^T表示向量或矩阵的转置运算。

3.根据权利要求2所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的将虚拟阵列的接收数据Y_r(t)和重排后的矩阵Y_t(t)分别进行子阵划分，获得发射端和接收端的旋转不变特性，按如下方式进行：

1)分别取Y_r(t)的前M(N-1)行和后M(N-1)行的输出数据构造成M(N-1)×P维矩阵，记为Y_r1(t)和Y_r2(t)，

2)分别取Y_t(t)的前N(M-1)行和后N(M-1)行的输出数据构造成N(M-1)×P维矩阵，记为Y_t1(t)和Y_t2(t)：

$\left\{\begin{array}{c}{Y}_{r1}\left(t\right)=C_{r1}{S}_{r1}\left(t\right)+N_{r1}\left(t\right)\\ Y_{r2}\left(t\right)=C_{r1}{\mathrm{\Φ}}_r{S}_{r2}\left(t\right)+N_{r2}\left(t\right)\end{array}\right.$

$\left\{\begin{array}{c}{Y}_{t1}\left(t\right)=C_{t1}{S}_{t1}\left(t\right)+N_{t1}\left(t\right)\\ Y_{t2}\left(t\right)=C_{t1}{\mathrm{\Φ}}_t{S}_{t2}\left(t\right)+N_{t2}\left(t\right)\end{array}\right.$

其中：为包含了所有发射方向信息的旋转对角矩阵，Φ_r＝diag[exp(-j2πd_rsinθ₁/λ),exp(-j2πd_rsinθ₂/λ),...,exp(-j2πd_rsinθ_P/λ)]为包含了所有接收方向信息的旋转对角矩阵，其中和θ分别为第p个目标的发射角和接收角，λ为载波波长，d_t和d_r分别为发射和接收阵元间距；S_r1(t)和S_r2(t)分别为S_r(t)的前M(N-1)行和后M(N-1)行构成的矩阵，S_t1(t)和S_t2(t)分别为重排后的信号矢量S_t(t)的前N(M-1)行和后N(M-1)行构成的矩阵；N_r1(t)、N_r2(t)∈C^M(N-1)×T和N_t1(t)、N_t2(t)∈C^N(M-1)×T为噪声矩阵，其中T为采样拍数。

4.根据权利要求3所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的根据四阶累积量的基本定义，构造出两个四阶累积量矩阵对，按如下方式进行：

$\left\{\begin{array}{c}{R}_{r1}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r1}^H\}=C_{r1}{R}_{4s}{C}_{r1}^H\\ R_{r2}=\mathrm{cum}\{y_{{r1}_1}^{*},y_{r1},y_{{r1}_1},y_{r2}^H\}=C_{r1}{R}_{4s}{\mathrm{\Φ}}_r^H{C}_{r1}^H\end{array}\right.$

$\left\{\begin{array}{c}{R}_{t1}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{t1}_1},y_{t1}^H\}=C_{t1}{R}_{4s}{C}_{t1}^H\\ R_{t2}=\mathrm{cum}\{y_{{t1}_1}^{*},y_{t1},y_{{t1}_1},y_{t2}^H\}=C_{t1}{R}_{4s}{\mathrm{\Φ}}_t^H{C}_{t1}^H\end{array}\right.$

其中：y_r1和y_t1分别为矩阵Y_r1和Y_t1的第1个列矢量，y_r2和y_t2分别为矩阵Y_r1和Y_t1的第2个列矢量，和分别为y_r1和y_t1的第一个元素；为四阶累积量矩阵的计算表达式，(·)^*表示向量或矩阵的共轭运算，(·)^H表示向量或矩阵的共轭转置运算。

5.根据权利要求4所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的对四阶累积量矩阵R_t1和R_r1进行特征值分解，分别获得包含发射和接收旋转不变特性的信号子空间，按如下步骤进行：

1)对四阶累积量矩阵R_t1和R_r1进行特征值分解，分别得到NM个特征值和NM个特征值所对应的特征矢量；

2)取R_t1分解得到的P个大特征值所对应的特征矢量构成信号子空间U_st；取R_r1分解得到的P个大特征值所对应的特征矢量构成信号子空间U_sr，其中P为目标数。

6.根据权利要求5所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的求解目标的收发角度，按如下步骤进行：

1)利用ESPRIT算法求解旋转发射不变关系矩阵：然后对Ψ_t进行特征值分解，由P个特征值可以求解出P个目标的发射角度；

2)利用ESPRIT算法求解接收旋转不变关系矩阵：对Ψ_r进行特征值分解，由P个特征值求解出P个目标的接收角度。

7.根据权利要求6所述的色噪声环境下的多输入多输出雷达收发角度快速估计方法，其特征是所述的利用最大似然估计方法实现收发角度的正确配对，按如下方式进行：

利用最大似然估计算法对角度进行匹配，其中P个使得代价函数Θ取值最小值的收发角度组合方式，即为正确的收发角度配对方式；所述代价函数为

其中：张成的子空间相互正交的投影矩阵，其中为联合导向矢量矩阵的第p个列矢量，p＝1,2,...,P，I_K为K×K维的单位阵，(·)^#表示向量或矩阵的违逆运算。