CN105974366A - 互耦条件下基于四阶累积量稀疏表示的mimo雷达波达方向估计方法 - Google Patents

Abstract

本发明属于单基地MIMO雷达***技术领域，具体涉及互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法。本发明包括：发射阵列发射相互正交的相位编码信号，接收端进行匹配滤波处理后获得接收数据，并利用发射和接收阵列都具有的互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响；构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵。本发明由于四阶累积量技术和加权稀疏表示框架的应用，成功地抑制了色噪声，在高斯色噪声情况下，本发明提供精确的波达方向估计，具有更高的角度分辨率和更好的角度估计性能。

Description

互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估 计方法

技术领域

本发明属于单基地MIMO雷达***技术领域，具体涉及互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法。

背景技术

多输入多输出(MIMO)雷达由新型的阵列结构构成，发射相互正交的波形，由于比传统相控阵雷达具有更高分辨率和更好的参数可识别性等很多优势，而吸引越来越多的关注。MIMO雷达可以分为两类，统计MIMO雷达和相干MIMO雷达(IEEESignalProcessingMagazine，2007，24(5)：106-114)。统计MIMO雷达通过大天线距获得空间增益，而相干MIMO雷达通过紧密的天线距形成虚拟阵列大孔径，获得更高的空间分辨率和更多的自由度。相干MIMO雷达包括双基地和单基地MIMO雷达，前者发射阵列和接收阵列是互相分开的，后者是紧密共置的。本发明中，我们研究的是单基地MIMO雷达中的波达方向估计问题。

角度估计在信号处理和雷达应用中具有重要作用，利用MIMO雷达和传统阵列信号模型的相似之处，大量关于MIMO雷达DOA估计的文献涌现出来，其中最具有代表性的是基于角度搜索或旋转不变特性的子空间角度估计方法，具体包括MUSIC，Capon和ESPRIT方法等。MUSIC和Capon具有很大的计算复杂度而不利于信号实时处理。ESPRIT算法具有高运算效率和更好的角度估计性能。派生的方法包括降维Capon，降维ESPRIT(RD-ESPRIT，ElectronicsLetters，2011，47(4)：283-284)和酉ESPRIT，RD-ESPRIT和酉ESPRIT比ESPRIT方法具有更低的计算量。发射波束域能量集中技术(IEEETransactionsonSignalProcessing，2011，59(6)：2669-2682)，使信噪比(SNR)增益得到最大化。然而，当考虑互耦时，以上子空间方法性能下降甚至失效。为了消除互耦影响，一些方法被提出来，相似MUSIC(SignalProcessing，2012，92(2)：517-522)通过角度搜索对角度进行估计，求根MUSIC具有更低的计算复杂度，相似ESPRIT方法(SignalProcessing，2012，92(12)：3039-3048)具有良好的角度估计性能，同时避免了空间峰值搜索。

在阵列信号处理中，经常遇到空间高斯色噪声的情况而不是白噪声，以上提到方法的波达方向估计性能进一步下降。高阶累积量在DOA估计上的应用能够解决这个问题，因为四阶累积量(FOC)能够抑制任意方差的高斯噪声。基于这个特性，FOC-MUSIC，FOC-ESPRIT等方法被提出来。互耦条件下，一种改进的FOC-MUSIC方法(SignalProcessing，2009，89：1839-1843)在高斯白噪声和色噪声条件下都具有良好的估计性能。

新兴的稀疏表示在信号分析领域吸引越来越多的关注，通过找到数据最稀疏的表示能够将其应用到DOA估计中。所有仿真结果证实，与传统子空间方法相比，稀疏表示方法具有显著优点，能够更好地适应特殊情况，此外提供更高角度分辨率并且更少地依赖入射信号数目的先验信息。为了估计波达方向一些方法被提出来，例如l₁-SVD，FOCUSS，l₁-SRACV和W-l₁-SRACV等。W-l₁-SRACV(IEEESignalProcessingLetters，2012，19(3)：155-158)通过设计阵列协方差矢量的加权l₁范数最小化框架而具有很好的分辨率。实值稀疏表示算法(IEEEAntennasandWirelessPropagationLetters，2013，12：376-379)降低了计算量。互耦条件下，改进的l₁-SVD(IEEEAntennasandWirelessPropagationLetters，2012，11：1210-1213)角度估计方法被提出来，未知互耦影响得以消除。然而，当考虑互耦或者噪声为高斯色噪声时，以上方法的估计性能下降。

发明内容

本发明的目的在于提供一种互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法。

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

互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，包括如下步骤：

(1)发射阵列发射相互正交的相位编码信号，接收端进行匹配滤波处理后获得接收数据，并利用发射和接收阵列都具有的互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响；

(2)构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵；

(3)对四阶累积量观测矩阵进行降维处理，获得稀疏表示框架下的相应模型，并利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解；

(4)设计加权l₁范数约束最小化的稀疏表示框架，利用编程软件包SOC计算方法，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声和色噪声条件下MIMO雷达***中目标DOA的精确估计。

步骤(1)中利用单基地MIMO雷达互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响按如下步骤：

(1.1)在快拍下，接收数据为

x(t)＝CAs(t)+n(t)

其中C_t和C_r分别是发射和接收阵列的互耦矩阵，是Kronecker积操作。A＝[a(θ₁),...,a(θ_P)]且p＝1,2,...P，P是远场目标总数。和分别是发射和接收导向矢量，θ_p是第p个目标的DOA，M和N分别是发射和接收天线数；s(t)＝[s₁(t),...,s_P(t)]^T是非高斯信号，n(t)是零均值高斯白噪声或色噪声；

(1.2)消除未知互耦影响，构造选择矩阵J，新的数据矢量为

其中I_Q表示Q×Q维单位矩阵，k+1为非零互耦系数个数。是新噪声， D是包含未知互耦信息的对角矩阵。

所述的步骤(2)中构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵按如下步骤：

(2.1)对消除互耦后的数据矩阵进行降维处理，构造降维转换矩阵J₃如下

其中且

(2.2)将降维转换矩阵J₃与数据矢量y(t)相乘，新的数据矢量为

其中B＝[b(θ₁),...,b(θ_P)]且 是去除互耦和冗杂行后的新噪声向量；

(2.3)构造具有特殊形式的四阶累积量矩阵：

其中表示中第k₂行且第k₁列的元素，通过收集的J个快拍，获得的四阶累积量矩阵的样本估计中每一项

其中是的第k_i个元素，在收集J个快拍条件下，

所述的步骤(3)中对四阶累积量观测矩阵进行降维处理，获得稀疏表示框架下的相应模型，并利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解按如下步骤：

(3.1)对四阶累积量观测矩阵进行降维处理：

其中V_s由对应P个大奇异值的右奇异向量组成，U和V分别由奇异值分解的左和右奇异向量组成，D_mc和C_s分别为包含奇异值，未知互耦和信号四阶累积量的对角矩阵，S_FB中第i行第j列元素为S_FB(i,j)＝|[Fb(θ_j)]_i|²[Fb(θ_j)]_i，[Fb(θ_j)]_i是Fb(θ_j)第i个值。

(3.2)利用导向矢量和噪声子空间的正交性，设计权值矩阵加强稀疏解：

其中是稀疏表示框架下的完备字典，为所有感兴趣DOA的离散样本网格，L＞＞P，包含完备字典中可能对应目标真实DOAs的P个导向矢量。 由中剩下导向矢量组成；U_n为噪声子空间，由U中第(P+1)到列组成；潜在地包含和 的第i个值是W中第i行2-范数，

所述步骤(4)中设计加权l₁范数约束最小化的稀疏表示框架，利用编程软件包SOC计算方法，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声和色噪声条件下MIMO雷达***中目标DOA的精确估计按如下步骤：

设计加权l₁范数约束最小化的稀疏表示框架，

其中是真实值的样本估计；是L×1维向量，稀疏矩阵满足η正则化参数选择为在四阶累积量估计误差的方差均一化基础上，具有自由度高概率1-ε置信区间的卡方分布的上限值，ε＝0.001足够；利用公式计算η；最后，通过SOC程序计算，测绘寻找P个峰值获得目标的DOAs。

本发明的有益效果在于：

1、本发明由于四阶累积量技术和加权稀疏表示框架的应用，成功地抑制了色噪声，在高斯色噪声情况下，本发明提供精确的波达方向估计，比相似ESPRIT，FOC-MUSIC和l₁-SVD 方法具有更高的角度分辨率和更好的角度估计性能。

2、本发明由于加权稀疏表示技术的应用，权值矩阵强化了稀疏解，无论噪声是高斯白噪声还是色噪声，本发明估计性能均优于相似ESPRIT，FOC-MUSIC和l₁-SVD方法，具有更低的SNR阈值，提供更高的角度分辨率；

3、本发明由于降维转换技术的应用和所设计的四阶累积量矩阵的特殊形式，计算量大大降低，本发明的计算复杂度比FOC-MUSIC方法更合理。

附图说明

图1本发明的整体框架图；

图2本发明互耦K＝2和K＝3时，SNR＝-10dB和SNR＝0dB不同情况下三个目标角度估计的空间谱；

图3不同方法互耦白噪声时三个目标角度估计均方根误差和信噪比关系；

图4不同方法互耦色噪声时三个目标角度估计均方根误差和信噪比关系；

图5不同方法互耦白噪声时两个目标角度估计均方根误差和角度间隔关系；

图6不同方法互耦色噪声时两个目标角度估计均方根误差和角度间隔关系；

图7不同方法互耦白噪声时三个目标角度估计均方根误差和快拍数关系；

图8不同方法互耦色噪声时三个目标角度估计均方根误差和快拍数关系；

图9不同方法互耦白噪声时三个目标角度估计分辨概率和信噪比关系；

图10不同方法互耦色噪声时三个目标角度估计分辨概率和信噪比关系。

具体实施方式

下面结合附图对本发明做进一步描述。

本发明提供一种未知互耦误差条件下基于四阶累积量稀疏表示的单基地多输入多输出(multiple-input multiple-output，MIMO)雷达波达方向(direction of arrival，简称DOA)估计方法，主要是为了解决目前在稀疏表示框架下，单基地MIMO雷达***中的角度估计方法存在不适用于高斯色噪声和互耦误差情况，估计精度不理想等缺点。首先通过线性转换消除互耦影响，构造降维转换矩阵减少计算复杂度。然后利用四阶累积量矩阵的有利形式，设计加权l₁范数约束最小化稀疏表示框架获得目标的DOA。其过程为：获得单基地MIMO雷达***的接收信号，利用发射和接收阵列都具有的互耦矩阵(MCM)的特殊结构，通过线性变换消除未知互耦的影响；然后构造降维转换矩阵，对消除互耦后的数据进行降维，进而基于新的数据矩阵构造四阶累积量观测矩阵，并进一步进行降维处理；获得稀疏表示模型，利用导向矢量和噪声子空间的正交性，设计权值矩阵；最后设计加权l₁范数约束最小化稀疏表示框架，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声或色噪声情况下MIMO雷达***中目标DOA的精确估计。与互耦误差条件下单基地MIMO雷达波达方向估计传统的相似ESPRIT，FOC-MUSIC和l₁-SVD方法相比，本发明在高斯白噪声和色噪声条件下都提供精确的DOA估计，同时具有更高的角度分辨率和更好的角度估计性能。

本发明的目的在于克服上述方法的缺陷，提供一种新的未知互耦误差条件下基于四阶累积量稀疏表示的单基地MIMO雷达波达方向估计方法。本发明利用发射和接收阵列都具有的互耦矩阵(MCM)带状对称Toeplitz结构，消除未知互耦的影响；然后基于新的接收数据，构造降维转换矩阵，进而构造具有特殊形式的四阶累积量矩阵以降低稀疏信号恢复的计算复杂度；最后设计加权l₁范数约束最小化的稀疏表示框架，通过找到恢复矩阵中的非零行获得目标的DOA。与未知互耦误差条件下单基地MIMO雷达角度估计的传统方法相比，由于利用了四阶累积量技术和加权稀疏表示框架，本发明在高斯白噪声和色噪声条件下都提供精确的DOA估计，同时具有更高的角度分辨率和更好的角度估计性能。本发明波达方向估计主要包括以下几个方面：

1、根据互耦条件下单基地MIMO雷达接收数据的结构和发射、接收阵列都具有的互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响。

考虑一个窄带单基地MIMO雷达***，发射阵列和接收阵列分别由M和N个天线组成，阵列均为半波长空间距的均匀线性阵列(ULAs)。在单基地MIMO雷达中发射和接收阵列紧密共置，因此对于一个远场目标可以看作它们具有相同的角度(即波达方向(DOA))。发射阵列利用M个天线发射M个具有相同带宽和中心频率的不同窄带正交波，假设目标数目为P，θ_p表示第p个目标的DOA。在有K＝k+1个非零互耦系数和min{M,N}＞2k的情况下，发射阵列和接收阵列都考虑互耦的影响。于是，在某个快拍下接收端匹配滤波器的输出表示为

其中是非高斯信号，是具有零均值的高斯白噪声或色噪声向量。和分别是发射和接收导向矢量，p＝1,2,...P。其中C_t和C_r分别是发射和接收阵列的互耦矩阵，互耦矩阵为带状对称Toeplitz矩阵并且表示为

其中c_ij(i＝t,r；j＝0,1,...,k)是互耦系数，与两个天线间的距离有关，并且满足0<|c_ik|<...<|c_i1|<|c_i0|＝1。基于的结构特点和Kronecker积操作的特性，未知互耦条件下MIMO雷达的接收数据为

其中是导向矩阵，并且关于信号和噪声引入如下假设：

1)信号是零均值平稳非高斯信号，彼此相互独立；

2)噪声是零均值高斯噪声，可能为白噪声或空间相关色噪声；

3)噪声和信号互相独立。

稀疏表示方法由于未知互耦矩阵的存在而失效。为构造有效的四阶累积量稀疏表示框架，首先需要消除发射和接收阵列的互耦影响。利用公式(2)中互耦矩阵带状对称Toeplitz结构特点，对于C_ta_t(θ_p)，p＝1,2,...,P，定义维选择矩阵以选择C_t的中间行。令有

其中是标量，并且定义另一个选择矩阵类似公式(4)，则

其中因此进一步推导出

其中令将J乘到x(t)的左边，得到

其中y(t)是新接收数据矢量，对角矩阵D＝diag[|c(θ₁)|,...,|c(θ_P)|]包含未知互耦的信息。和分别是新导向矩阵和噪声。由于互耦矩阵C转变成为了对角矩阵，因此互耦不再对导向矩阵存在影响，即公式(7)中互耦误差影响得到补偿。

2、构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵。

构造降维转换矩阵，对消除互耦后的数据进行降维处理，如下所示

在消除未知互耦影响后，为了降低下面四阶累积量计算的复杂度，基于单基地MIMO雷达的特殊结构构造降维转换。的第p列为

可以发现中包含很多重复项，并且满足

其中G和b(θ_p)的具体表达式如下所示

其中且 重复项在下面的四阶累积量计算中将导致大量的冗杂信息，进而在稀疏信号恢复中导致很大的计算量。为了解决这个问题，构造降维转换如下

因此，将降维转换矩阵J₃与数据矢量y(t)相乘，有

其中 和分别是消除互耦影响和去除冗杂行后的新数据向量和噪声。同时，F计算为

其中降维转换不仅将接收数据矩阵的维数从显著减少到而且避免了附加空间色噪声，J为快拍数。

基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵，如下所示

在补偿了未知互耦误差影响和降低了接收数据维数之后，通过收集的J个快拍，数据的四阶累积量定义为

其中是的第k_i个元素，在J个快拍条件下，的估计中的每一项按如下公式获得

根据公式(14)，设计四阶累积量矩阵如下所示

其中表示中第k₂行且第k₁列的元素。基于四阶累积量的特性和信号s_i(i＝1,2,...,P)之间相互独立等假设，推导出

其中和分别是b(θ_p)中第k₂和第k₁个元素，是信号s_p的四阶累积量。构造四阶累积量矩阵后，个值被合理地排列在中。由于中包含k₁和k₂的所有可能的排列置换，推导出如下所示

其中D_mc＝diag(|c(θ₁)|⁴,...,|c(θ_P)|⁴)和C_s＝diag(β₁,...,β_P)为对角矩阵。同时，S_FB中的第i行且第j列元素为

S_FB(i,j)＝|[Fb(θ_j)]_i|²[Fb(θ_j)]_i(19)

其中1≤j≤P，[Fb(θ_j)]_i是Fb(θ_j)中的第i个值。

3、对四阶累积量观测矩阵进行降维处理，获得稀疏表示框架下的相应模型，并利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解。

对四阶累积量观测矩阵进行降维处理，如下所示

为了进一步降低DOA估计稀疏表示方法的计算复杂度，四阶累积量观测矩阵能够从维转换成为维，如下所示

其中包含绝大部分的信号能量，可以用来代替对DOAs进行精确的估计。V_s由对应着大奇异值的右奇异向量组成，此外，U和V分别由左奇异向量和右奇异向量组成，是对角矩阵，且是奇异值。令D_mc和C_s均为对角矩阵。当不考虑互耦时原始的恢复矩阵模型中的某些行是零向量时，T中的对应行也必然为零，即使D_mc和的组成是完全未知的。这表明D_mc，和V_s的引入不会改变稀疏解，因此通过公式(7)的线性转换消除互耦影响后，在公式(20)现有的形式中，稀疏表示方法能够成功地对波达方向进行估计。

获得稀疏表示框架下的相应模型，如下所示

为了应用稀疏表示方法估计DOA，令为所有感兴趣DOAs的离散样本网格。可能的DOA数目远远大于目标数目，即L＞＞P，进而为可能DOAs构造完备字典令满足当空域某点有真实目标存在时，中对应行为非零向量而其他为零向量。因此，与原始矩阵T具有相同的行支持，我们能够将DOA估计问题转换成为找到完备字典中对应非零行的位置。为估计非零行的最小数目，直接的稀疏度量是l₀-范数处罚，但l₀-范数最小化是非凸优化且NP困难而无法解决。l₁-范数处罚能够用来合理地解决这个问题。为获得稀疏矩阵构造如下l₁-范数约束最小化框架

其中且η是正则化参数，通过收集J个快拍是真实值的样本估计，从公式(14)(15)(16)和(20)获得。

利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解，如下所示

为了加强l₁-范数处罚以更好地接近l₀-范数处罚，我们提出如下加权l₁-范数最小化稀疏表示方法。将完备字典按列分为两部分，一部分包含可能对应目标真实DOAs的P个导向矢量，另一部分由剩下的对应着恢复矩阵中零行的导向矢量组成。因此，能够表示为于是选择U中从(P+1)到的共列作为噪声子空间U_n，表示为

其中当J→∞时，已经被证明W₁中的每一项趋于零，因此基于此特性，设计一个加权矩阵如下所示

其中小，且同时当J→∞时， 和分别是和的第i项和第j项。潜在地包含两部分，故潜在地包含和W_r能够实现大权值处罚稀疏向量中可能为零的元素，小权值保留值更大的项。换言之，加权矩阵的引入能够加强稀疏解，提高DOA估计的精确度。

4、设计加权l₁范数约束最小化的稀疏表示框架，利用编程软件包SOC(二阶锥)计算方法，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声和色噪声条件下MIMO雷达***中目标DOA的精确估计。

构造加权稀疏表示框架如下所示

其中η是设置误差量的正则化参数，η数值的选择在最后的DOA估计中起到很重要的作用。令为四阶累积量估计误差，为了获得的分布特性，首先推导的分布，并且 是真实值的样本估计，因此基于的元素结构，可以认为是一系列的四阶累积量估计值。根据文献(Signals，SystemsandComputers，1993，2：1186-1190)，对于统计过程四阶累积量估计中所有项具有满足渐近正态分布的特性。此外，因此多快拍情况下，是具有零均值的渐近正态分布，即

其中AsN(μ,V)表示均值为μ，协方差矩阵为的渐近正态分布。并且根据文献可知四阶累积量估计误差的协方差如下所示

其中τ＝k₁,k₂,k₃,k₄，ρ＝k₁',k'₂,k₃',k'₄，缩写表示Q₄₄＝Q₄₄(τ,ρ)，m₂(λ)＝m₂(k₁,λ)，且m₂(λ')＝m₂(k₁',λ')。根据可以获得估计值，并且k₁₁＝k₂,k₁₂＝k₃,k₁₃＝k₄；k₂₁＝k₃,k₂₂＝k₂,k₂₃＝k₄；k₃₁＝k₄,k₃₂＝k₃,k₃₃＝k₂，同理k'具有相同定义。同时，利用公式(15)获得k阶样本矩

当k₂＝k₃＝k₄且k'₂＝k₃'＝k'₄，V由构成，来自于推导出i和k₁，k₂的关系如下所示

其中[·]_rem表示值的余数，[·]_ceil表示不小于值的最小的整数。j和k₁'，k'₂的关系与公式(27)相同。由于利用矢量化操作的特性，可以推导出

其中基于公式(28)和中渐近正态分布的线性转换不变特性，推导出中关于均值和协方差矩阵的结论如下所示

其中即满足均值为零且协方差矩阵为E_vV(E_v)^H的渐近正态分布。因此，根据文献(IEEETransactionsonSignalProcessing，2005，53(8)：3010-3022)，本发明中正则化参数η选择为在四阶累积量估计误差的方差均一化基础上，具有自由度高概率1-ε置信区间的卡方分布的上限值，ε＝0.001足够。利用Matlab软件，公式可以用来计算η。另一方面，公式(24)可以通过SOC(二阶锥)程序软件包计算，例如SeDuMi和CVX。最后，通过解决公式(24)，测绘获得目标的DOAs。

实施例1

步骤一、建立未知互耦误差条件下单基地MIMO雷达的接收信号模型，根据互耦矩阵结构特点，利用线性变换消除未知互耦影响。

考虑一个窄带单基地MIMO雷达***，发射阵列和接收阵列分别由M和N个天线组成，阵列均为半波长空间距的均匀线性阵列。在单基地MIMO雷达中发射和接收阵列紧密共置，因此对于一个远场目标可以看作它们具有相同的角度(即波达方向(DOA))。发射阵列利用M个天线发射M个具有相同带宽和中心频率的不同窄带正交波，假设目标数目为P，θ_p表示第p个目标的DOA。在有K＝k+1个非零互耦系数和min{M,N}＞2k的情况下，发射阵列和接收阵列都考虑互耦的影响。于是，在某个快拍下接收端匹配滤波器的输出表示为

其中是非高斯信号，是具有零均值的高斯白噪声或色噪声向量。和分别是发射和接收导向矢量，p＝1,2,...P。其中C_t和C_r分别是发射和接收阵列的互耦矩阵，互耦矩阵为带状对称Toeplitz矩阵并且表示为

其中c_ij(i＝t,r；j＝0,1,...,k)是互耦系数，0<|c_ik|<...<|c_i1|<|c_i0|＝1，与天线间距离有关。基于结构和Kronecker积特性，互耦条件下MIMO雷达接收数据为

其中是导向矩阵，并且关于信号和噪声作如下假设：信号是零均值平稳非高斯信号，彼此相互独立；噪声是零均值高斯噪声，可能为白噪声或空间相关色噪声；噪声和信号互相独立。

稀疏表示方法由于未知互耦矩阵的存在而失效。为构造有效的四阶累积量稀疏表示框架，首先需要消除发射和接收阵列的互耦影响。利用公式(31)中互耦矩阵带状对称Toeplitz结构特点，对于C_ta_t(θ_p)，p＝1,2,...,P，定义维选择矩阵以选择C_t的中间行。令有

其中是标量，且定义另一个选择矩阵类似公式(33)，可以得到

其中因此进一步推导出

其中c(θ_p)＝c_t(θ_p)c_r(θ_p)，令将J乘到x(t)的左边，得到

其中y(t)是新接收数据矢量，对角矩阵D＝diag[|c(θ₁)|,...,|c(θ_P)|]包含未知互耦信息。和是新导向矩阵和噪声。互耦矩阵转变成了对角矩阵，因此互耦不再对导向矩阵存在影响，公式(36)中互耦误差得到补偿。

步骤二、构造降维转换矩阵，对消除互耦后的数据进行降维处理。

在消除未知互耦影响后，为了降低下面四阶累积量计算的复杂度，基于单基地MIMO雷达的特殊结构构造降维转换。的第p列为

可以发现中包含很多重复项，并且满足

其中G和b(θ_p)的具体表达式如下所示

其中且 重复项在下面的四阶累积量计算中将导致大量的冗杂信息，进而在稀疏信号恢复中导致很大的计算量。为了解决这个问题，构造降维转换如下

因此，将降维转换矩阵J₃与数据矢量y(t)相乘，得到新数据如下

其中 是新噪声向量。同时，F计算为

其中降维转换不仅将接收数据矩阵的维数从显著减少到而且避免了附加空间色噪声，J为快拍数。

步骤三、基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵。

在补偿了未知互耦误差影响和降低了接收数据维数之后，通过收集的J个快拍，数据的四阶累积量定义为

其中是的第k_i个元素，在收集J个快拍条件下，的样本估计中的每一项按如下公式获得

根据公式(43)，设计四阶累积量矩阵如下所示

其中表示中第k₂行且k₁列的元素。基于信号和噪声的假设，推导出

其中和分别是b(θ_p)中第k₂和第k₁个元素，是信号s_p的四阶累积量。构造四阶累积量矩阵后，个值被合理地排列在中。由于中包含k₁和k₂的所有可能的排列置换，推导出如下所示

其中D_mc＝diag(|c(θ₁)|⁴,...,|c(θ_P)|⁴)和C_s＝diag(β₁,...,β_P)为对角矩阵。同时，S_FB中的第i行且第j列元素为

S_FB(i,j)＝|[Fb(θ_j)]_i|²[Fb(θ_j)]_i(48)

其中1≤j≤P，[Fb(θ_j)]_i是Fb(θ_j)中的第i个值。

步骤四、对四阶累积量观测矩阵进行降维处理。

为了进一步降低下面DOA估计稀疏表示方法的计算复杂度，四阶累积量观测矩阵能够从维转换成为维，如下所示

其中包含绝大部分信号能量，可以用来代替对DOA进行精确估计。V_s由对应着大奇异值的右奇异向量组成，此外，U和V分别由左和右奇异向量组成，是对角矩阵，是奇异值。令D_mc和C_s均为对角矩阵。当不考虑互耦时原始恢复矩阵模型中的某些行是零向量时，T中的对应行也必然为零，即使D_mc和的组成完全未知。这表明D_mc，和V_s的引入不会改变稀疏解，故通过公式(36)消除互耦后，在公式(49)现有形式中，稀疏表示方法能够成功地对DOA进行估计。

步骤五、获得稀疏表示框架下的相应模型，利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解。

为了应用稀疏表示方法估计DOA，令为所有感兴趣DOAs的离散样本网格。可能的DOAs数目远远大于目标数目，即L＞＞P，进而为可能DOAs构造完备字典令满足当空域某点有真实目标存在时，中对应行为非零向量而其他为零向量。因此，与原始矩阵T具有相同行支持，我们能够将DOA估计问题转换成为找到完备字典中对应非零行的位置。为估计非零行的最小数目，直接的稀疏度量是l₀-范数处罚，但l₀-范数最小化是非凸优化且NP困难而无法解决。l₁-范数处罚能够用来合理地解决这个问题。为了加强l₁-范数处罚以更好地接近l₀-范数处罚，我们提出如下加权l₁-范数最小化稀疏表示方法。将完备字典按列分为两部分，一部分包含可能对应目标真实DOAs的P个导向矢量，另一部分由剩下的对应着恢复矩阵中零行的导向矢量组成。因此，能够表示为于是选择U中从(P+1)到的共列作为噪声子空间U_n，表示为

其中当J→∞时，已经被证明W₁中的每一项趋于零，因此基于此特性，设计一个权值矩阵如下所示

其中小，且同时当J→∞时， 和分别是和的第i项和第j项。潜在地包含两部分，故潜在地包含和W_r实现了大权值处罚那些稀疏向量中可能为零的元素，小权值保留值更大的项。换言之，加权矩阵的引入能够加强稀疏解，提高DOA估计精确度。

步骤六、设计加权l₁范数约束最小化稀疏表示框架，获得恢复矩阵，寻找恢复矩阵中的非零行，获得目标精确的DOA。

构造加权稀疏表示框架如下所示

其中且η是正则化参数，通过收集J个快拍是真实值的样本估计，从公式(43)(44)(45)和(49)获得。η是设置误差量的正则化参数，η数值的选择在最后的DOA估计中起到很重要的作用。

令为四阶累积量估计误差，为了获得的分布特性，我们首先推导的分布，并且 是真实值的样本估计，因此基于的元素结构，可以认为是一系列的四阶累积量估计值。根据文献(Signals，SystemsandComputers，1993，2：1186-1190)，对于统计过程 且四阶累积量估计中的所有项具有满足渐近正态分布的特性，因此多快拍情况下，是具有零均值的渐近正态分布，即

其中AsN(μ,V)表示均值为μ，协方差矩阵为的渐近正态分布。并且根据文献可知四阶累积量估计误差的协方差如下所示

其中τ＝k₁,k₂,k₃,k₄，ρ＝k₁',k'₂,k₃',k'₄，缩写表示Q₄₄＝Q₄₄(τ,ρ)，m₂(λ)＝m₂(k₁,λ)，且m₂(λ')＝m₂(k₁',λ')。根据可以获得估计值，并且k₁₁＝k₂,k₁₂＝k₃,k₁₃＝k₄；k₂₁＝k₃,k₂₂＝k₂,k₂₃＝k₄；k₃₁＝k₄,k₃₂＝k₃,k₃₃＝k₂，同理k'具有相同定义。同时，利用公式(44)获得k阶样本矩

当k₂＝k₃＝k₄且k'₂＝k₃'＝k'₄，V由构成，来自于推导出i和k₁，k₂的关系如下所示

其中[·]_rem表示值的余数，[·]_ceil表示不小于值的最小的整数。j和k₁'，k'₂的关系与公式(55)相同。由于利用矢量化操作的特性，可以推导出

其中基于公式(56)和中渐近正态分布的线性转换不变特性，推导出中关于均值和协方差矩阵的结论如下所示

其中即满足均值为零且协方差矩阵为E_vV(E_v)^H的渐近正态分布。因此，根据文献(IEEETransactionsonSignalProcessing，2005，53(8)：3010-3022)，本发明中正则化参数η选择为在四阶累积量估计误差的方差均一化基础上，具有自由度高概率1-ε置信区间的卡方分布的上限值，ε＝0.001足够。利用Matlab软件，公式可以用来计算η。另一方面，公式(52)可以通过SOC(二阶锥)程序软件包计算，例如SeDuMi和CVX。最后，通过解决公式(52)，测绘获得目标的DOAs。

本发明的效果可通过以下内容说明：

(一)运算复杂度分析：

本发明的主要计算量是四阶累积量的计算和通过SOC程序包解决公式(52)获得恢复矩阵的稀疏解。前者由于所设计的四阶累积量矩阵的特殊形式只需要后者需要其中L₁是离散样本网格总数，因此本发明的主要计算复杂度是然而，在FOC-MUSIC方法中，由于所有的(k₁,k₂,k₃,k₄)在矩阵C_x中都需要被任意地排列置换，因而需要的计算量为降维转换的应用和具有优势的四阶累积量矩阵形式使得计算量大大减少，故本发明的计算复杂度比FOC-MUSIC方法更合理。

(二)仿真条件与内容：

我们呈现一些仿真结果来证实本发明的有效性和优点。用相似ESPRIT，FOC-MUSIC和l₁-SVD方法与本发明进行对比。角度估计的均方根误差(RMSE)用来评估角度估计性能，其定义为

其中是第i次MonteCarlo试验DOA真实值θ_p的估计，Q是MonteCarlo试验总次数，下面的仿真中选择Q＝500。考虑具有M个发射天线和N个接收天线的窄带单基地MIMO雷达***，发射和接收阵列均为半波长空间距的均匀线性阵列。此外，仿真中的接收数据是混有高斯噪声和阵元间互耦影响的二进制相移键控(BPSK)信号，大多数情况下，假设不相关目标的数目P是已知的，SNR定义为S'和N'是J快拍下信号和噪声矩阵。此外，对于本发明以及FOC-MUSIC方法和l₁-SVD方法，离散样本网格均从-90°到90°变换，均匀间隔为0.05°。

(三)仿真结果：

1、本发明互耦K＝2和K＝3时，SNR＝-10dB和SNR＝0dB不同情况下三个目标角度估计的空间谱

图二展示了本发明对于非零互耦系数K＝2和K＝3时，SNR＝-10dB和SNR＝0dB不同情况下的空间谱，其中M＝N＝7，J＝2000，三个不相关目标的DOAs分别为θ₁＝-20°，θ₂＝0°和θ₃＝20°。当K＝3时，发射阵列的非零互耦系数为[c_t0,c_t1,c_t2]＝[1,0.6+j0.2,0.02+j0.1]，接收阵列为[c_r0,c_r1,c_r2]＝[1,0.5+j0.3,0.01+j0.2]。当K＝2时，它们为[c_t0,c_t1]＝[1,0.0174+j0.0377]和[c_r0,c_r1]＝[1,0.0521-j0.1029]，并且将这种情况应用到下面的仿真中。如图二所示，本发明的空间峰值很尖锐，旁瓣压制很低，即使在k增加，SNR降低的情况下。这表明本发明能够提供极好的估计性能。

2、不同方法互耦白噪声时三个目标角度估计均方根误差和信噪比关系

图三展示了高斯白噪声条件下，不同方法DOA估计的RMSE和SNR的关系，其中三个不相关目标为θ₁＝-20°，θ₂＝0°和θ₃＝20°，同时M＝N＝7，J＝4000。从图三可以看出，在低SNR区域l₁-SVD方法比其他方法具有更好的性能。然而，本发明由于利用了所设计的四阶累积量矩阵的特殊形式和加权稀疏表示框架，当SNR超过大约0dB时，本发明与相似ESPRIT，FOC-MUSIC和l₁-SVD方法相比，提供最好的角度估计性能。

3、不同方法互耦色噪声时三个目标角度估计均方根误差和信噪比关系

图四展示了高斯色噪声情况下，DOA估计的RMSE和SNR的关系，其他条件与图三保持一致。与l₁-SVD，相似ESPRIT和FOC-MUSIC方法相比，图四很明显说明了本发明的估计性能在所有SNR区域都是最好的。因为利用了四阶累积量的特殊性质，本发明成功地抑制了色噪声，然而在l₁-SVD和相似ESPRIT方法中，色噪声导致了性能的下降，尤其在低SNR区域。

4、不同方法互耦白或色噪声时两个目标角度估计均方根误差和角度间隔关系

图五和图六展示了RMSE和角度间隔的关系，其中M＝N＝7，J＝4000，图五中考虑高斯白噪声以及SNR＝5dB，而图六考虑的是高斯色噪声且SNR＝0dB。每种方法中有两个不相关目标，DOA分别为θ₁＝0°和θ₂＝0°+Δθ，并且Δθ从2°到14°变化。由图五可知，当Δθ在某个值之内时，本发明的估计性能和相似ESPRIT方法相近，并且都比FOC-MUSIC和l₁-SVD方法性能好。当角度间隔超过大约6°时，本发明在所有方法中性能最好。如图六所示，色噪声下本发明对于空间相近目标提供最好的角度估计性能，即本发明具有最高的空间角度分辨率。

5、不同方法互耦白或色噪声时三个目标角度估计均方根误差和快拍数关系

图七和图八展示了不同方法DOA估计的RMSE和快拍数的关系，M＝N＝7。图七中噪声是高斯白噪声且SNR＝5dB，而图八为高斯色噪声且SNR＝0dB。考虑三个具有不同DOA的不相关目标为θ₁＝-11.5°，θ₂＝0°和θ₃＝11.5°。从图七和图八可以看出，对于高斯白噪声和色噪声，本发明在所有快拍区域都提供很好的估计性能，此外，本发明相比于相似ESPRIT，FOC-MUSIC和l₁-SVD方法的估计优势随着J的增加而越来越突出。

6、不同方法互耦白或色噪声时三个目标角度估计分辨概率和信噪比关系

图九和图十分别展示了高斯白噪声和高斯色噪声条件下，不同方法的目标分辨概率和SNR的关系，其中M＝N＝7，J＝2000。三个不相关目标的DOAs分别为θ₁＝-11°，θ₂＝0°和θ₃＝11°。此外，当三个目标所有DOA估计误差的绝对值都在0.1°之内时，目标被认为成功分辨。如图所示，当SNR足够高时，所有方法提供100％的目标分辨概率，然而，每种方法的分辨概率在某一点开始下降，定义为SNR临界值。从图九和图十可以清楚看到，本发明与相似ESPRIT，FOC-MUSIC和l₁-SVD方法相比具有最低的SNR临界值，同时它们表明无论噪声是高斯白噪声还是色噪声，本发明都提供更高的目标分辨概率。

Claims (5)

1.互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，其特征在于：包括如下步骤：

(1)发射阵列发射相互正交的相位编码信号，接收端进行匹配滤波处理后获得接收数据，并利用发射和接收阵列都具有的互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响；

(2)构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵；

(3)对四阶累积量观测矩阵进行降维处理，获得稀疏表示框架下的相应模型，并利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解；

(4)设计加权l₁范数约束最小化的稀疏表示框架，利用编程软件包SOC计算方法，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声和色噪声条件下MIMO雷达***中目标DOA的精确估计。

2.根据权利要求1所述的互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，其特征在于：步骤(1)中利用单基地MIMO雷达互耦矩阵带状对称Toeplitz结构特点，通过线性变换消除未知互耦的影响按如下步骤：

(1.1)在快拍下，接收数据为

x(t)＝CAs(t)+n(t)

其中C_t和C_r分别是发射和接收阵列的互耦矩阵，是Kronecker积操作。A＝[a(θ₁),...,a(θ_P)]且P是远场目标总数。和分别是发射和接收导向矢量，θ_p是第p个目标的DOA，M和N分别是发射和接收天线数；s(t)＝[s₁(t),...,s_P(t)]^T是非高斯信号，n(t)是零均值高斯白噪声或色噪声；

(1.2)消除未知互耦影响，构造选择矩阵J，新的数据矢量为

$y\left(t\right)=Jx\left(t\right)=\tilde{A}Ds\left(t\right)+\tilde{n}\left(t\right)$

其中I_Q表示Q×Q维单位矩阵，k+1为非零互耦系数个数。是新噪声， D是包含未知互耦信息的对角矩阵。

3.根据权利要求1所述的互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，其特征在于：所述的步骤(2)中构造降维转换矩阵，对消除互耦后的数据进行降维处理，进而基于新的数据矩阵构造具有特殊形式的四阶累积量矩阵按如下步骤：

(2.1)对消除互耦后的数据矩阵进行降维处理，构造降维转换矩阵J₃如下

$J_3={\left(G^{H}G\right)}^{(-\frac{1}{2})}{G}^H$

其中且

(2.2)将降维转换矩阵J₃与数据矢量y(t)相乘，新的数据矢量为

$\tilde{y}\left(t\right)=J_{3}y\left(t\right)=FBDs\left(t\right)+\stackrel{\&OverBar;}{\tilde{n}}\left(t\right)$

其中B＝[b(θ₁),...,b(θ_P)]且 是去除互耦和冗杂行后的新噪声向量；

(2.3)构造具有特殊形式的四阶累积量矩阵：

$C_{\tilde{y}}(k_2,k_1)=cum\{{\tilde{y}}_{k_2},{\tilde{y}}_{k_1}^{*},{\tilde{y}}_{k_1},{\tilde{y}}_{k_1}^{*}\}$

其中表示中第k₂行且第k₁列的元素，通过收集的J个快拍，获得的四阶累积量矩阵的样本估计中每一项

$\begin{array}{c}cum({\tilde{y}}_{k_1},{\tilde{y}}_{k_2}^{*},{\tilde{y}}_{k_3},{\tilde{y}}_{k_4}^{*})=E\left({\tilde{y}}_{k_1}{\tilde{y}}_{k_2}^{*}{\tilde{y}}_{k_3}{\tilde{y}}_{k_4}^{*}\right)-E\left({\tilde{y}}_{k_1}{\tilde{y}}_{k_2}^{*}\right)E\left({\tilde{y}}_{k_3}{\tilde{y}}_{k_4}^{*}\right)\\ -E\left({\tilde{y}}_{k_1}{\tilde{y}}_{k_3}\right)E\left({\tilde{y}}_{k_2}^{*}{\tilde{y}}_{k_4}^{*}\right)-E\left({\tilde{y}}_{k_1}{\tilde{y}}_{k_4}^{*}\right)E\left({\tilde{y}}_{k_2}^{*}{\tilde{y}}_{k_3}\right)\end{array}$

其中是的第k_i个元素，在收集J个快拍条件下，

4.根据权利要求1所述的互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，其特征在于：所述的步骤(3)中对四阶累积量观测矩阵进行降维处理，获得稀疏表示框架下的相应模型，并利用导向矢量和相应噪声子空间的正交性，设计权值矩阵以加强稀疏解按如下步骤：

(3.1)对四阶累积量观测矩阵进行降维处理：

${\tilde{C}}_{\tilde{y}}=C_{\tilde{y}}{V}_s={\mathrm{FBD}}_{mc}{C}_s{S}_{FB}^H{V}_s$

其中V_s由对应P个大奇异值的右奇异向量组成，U和V分别由奇异值分解的左和右奇异向量组成，D_mc和C_s分别为包含奇异值，未知互耦和信号四阶累积量的对角矩阵，S_FB中第i行第j列元素为S_FB(i,j)＝|[Fb(θ_j)]_i|²[Fb(θ_j)]_i，[Fb(θ_j)]_i是Fb(θ_j)第i个值。

(3.2)利用导向矢量和噪声子空间的正交性，设计权值矩阵加强稀疏解：

$W_r=diag\&lsqb;{\left(W_1^{\left(l_2\right)}\right)}^T,{\left(W_2^{\left(l_2\right)}\right)}^T\&rsqb;/max\left(W_2^{\left(l_2\right)}\right)$

其中是稀疏表示框架下的完备字典，为所有感兴趣DOA的离散样本网格，L＞＞P，包含完备字典中可能对应目标真实DOAs的P个导向矢量。 由中剩下导向矢量组成；U_n为噪声子空间，由U中第(P+1)到列组成；潜在地包含和 的第i个值是W中第i行2-范数，

5.根据权利要求1所述的互耦条件下基于四阶累积量稀疏表示的MIMO雷达波达方向估计方法，其特征在于：所述步骤(4)中设计加权l₁范数约束最小化的稀疏表示框架，利用编程软件包SOC计算方法，获得恢复矩阵，寻找恢复矩阵中的非零行，实现对高斯白噪声和色噪声条件下MIMO雷达***中目标DOA的精确估计按如下步骤：

设计加权l₁范数约束最小化的稀疏表示框架，

$\begin{array}{cc}\underset{T_{\widehat{\mathrm{\&theta;}}}}{\min }\left|\right|W_r{T}_{\widehat{\mathrm{\&theta;}}}^{\left(l_2\right)}|{|}_1,& s.t.\left|\right|{\hat{\tilde{C}}}_{\tilde{y}}-{\tilde{B}}_{\widehat{\mathrm{\&theta;}}}{T}_{\widehat{\mathrm{\&theta;}}}|{|}_F\&le;\sqrt{\mathrm{\&eta;}}\end{array}$

其中是真实值的样本估计；是L×1维向量，稀疏矩阵满足η正则化参数选择为在四阶累积量估计误差的方差均一化基础上，具有自由度高概率1-ε置信区间的卡方分布的上限值，ε＝0.001足够；利用公式计算η；最后，通过SOC程序计算，测绘寻找P个峰值获得目标的DOAs。