基于root‑MUSIC算法的目标波达方向获取方法与流程

本发明属于雷达技术领域，特别涉及一种基于root-MUSIC算法的目标波达方向获取方法，即基于求根的多重信号分类(root-MUSIC)算法的目标波达方向获取方法，适用于确定任意波束的目标波达方向。

背景技术：

目标的波达方向DOA估计在阵列信号处理中具有重要的意义，在过去的近40年里，人们提出了一系列的目标的波达方向DOA估计算法；其中，MUSIC算法由于其自身的高分辨能力得到了广泛的应用。然而，传统的MUSIC算法实现需要在整个角度范围内进行一维或多维搜索，计算量很大。为了降低计算量，人们提出了root-MUSIC算法，然而它的计算量是阵元数的三次方，并且阵元数的增加会带来更大的计算量；与此同时，考虑到传统的波束空间处理的优势，如计算量低，分辨目标需要的SNR低，人们提出了波束空间MUSIC算法，研究表明，波束空间MUSIC算法在波束形成共同的带外零点特性时，如波束形成矩阵为离散傅里叶变换DFT矩阵的情况，波束空间root-MUSIC算法的实现的计算量仅为波束数的三次方，共同的带外零点特性通常是不现实的。

技术实现要素：

针对传统的波束空间root-MUSIC算法计算量大，且对波束形成矩阵有限制的问题，本发明的目的在于提出一种基于root-MUSIC算法的目标波达方向获取方法，该种基于root-MUSIC算法的目标波达方向获取方法利用多窗低阶实多项式求根来实现波束空间root-MUSIC，相比传统算法计算量更低，估计的结果更稳定，且可适用于任意波束形成矩阵中目标波达方向的确定。

为实现上述目的，本发明采用如下技术方案予以实现。

一种基于root-MUSIC算法的目标波达方向获取方法，包括以下步骤：；

步骤1，确定雷达，所述雷达包含N个阵元，并确定雷达检测范围内的波束形成矩阵B，所述雷达检测范围内的波束形成矩阵B中包含N_s个目标信号；

然后确定雷达检测范围内的波束形成矩阵的空间频谱区域范围，并将所述空间频谱区域范围以间隔ε均匀划分为M个空间频谱区域，并顺序排列，其中M＝2N_b-1，N_b表示波束形成矩阵B包含的波束个数，N、N_b、M、ε分别为大于0的整数；

步骤2，选取顺序排列后M个空间频谱区域中的第1个空间频谱区域，并计算得到第1个空间频谱区域的窗矢量矩阵T，进而计算得到第1个空间频谱区域的窗矢量矩阵估计

步骤3，确定雷达接收到的回波数据为X，然后计算得到雷达回波数据的采样协方差矩阵估计对雷达回波数据的采样协方差矩阵估计进行特征分解，并计算得到雷达回波数据的采样协方差矩阵估计的噪声子空间矩阵

步骤4，利用雷达检测范围内的波束形成矩阵B和雷达回波数据的采样协方差矩阵估计的噪声子空间矩阵计算得到雷达检测范围内的波束形成矩阵的噪声协方差矩阵V；

步骤5，将雷达检测范围内的波束形成矩阵的噪声协方差矩阵V中包含的所有的对角元素构造列矢量v；

步骤6，确定M个空间频谱区域的变换矩阵W，并根据列矢量v和第1个空间频谱区域的窗矢量矩阵估计计算得到中间矩阵D；

然后，利用M个空间频谱区域的变换矩阵W对中间矩阵D进行变换，得到M个空间频谱区域对应的M组低阶实多项式系数；

步骤7，根据M个空间频谱区域对应的M组低阶实多项式系数，对应得到M个方程，其中第m个方程为f^l-1为第l-1次幂的自变量，f⁰＝1，γ_ml为第m组低阶实多项式系数中的第l阶系数，求解第m个方程后得到L个共轭根，每个共轭根包含实部和虚部，然后在L个共轭根中选择虚部满足设定要求的N'_s个共轭根，所述设定要求为虚部模值小于10^-6的共轭根，该N'_s个共轭根各自对应的实部，分别对应为雷达检测范围内的波束形成矩阵B内N_s个目标信号各自空间角度的正弦值，并且所述N_s个目标信号各自角度的正弦值，分别对应为雷达检测范围内的波束形成矩阵B内N_s个目标信号各自的波达方向；

其中，雷达检测范围内的波束形成矩阵B中包含的目标信号个数与L个共轭根中虚部满足设定要求的共轭根个数相等且一一对应。

本发明具有以下优点：

第一，本发明将高阶多项式分解为多组低阶多项式，通过低阶实多项式求根来降低计算复杂度，提高实时处理能力；另外，相比传统的root-MUSIC仅适用于DFT波束形成矩阵的情况，本发明方法能够适用于任意波束。

第二，本发明方法利用多窗实多项式求根来降低波束空间root-MUSIC算法的复杂度，并针对传统的波束空间root-MUSIC算法在大阵列情况下的计算量问题，将高阶多项式分解为多组低阶多项式，通过低阶实多项式求根来降低计算复杂度，提高实时处理能力。

第三，相比传统的root-MUSIC在实际操作中通常需要利用加窗等锥削方法来降低旁瓣，进而抑制干扰，因此传统的root-MUSIC仅适用于DFT波束形成矩阵；相比传统的root-MUSIC，本发明方法考虑到在小角度区域内(如-3dB波束宽度内)导向矢量的相关性，该小角度区域内任意导向矢量可以利用其对应空间频率的低阶实多项式来近似，使得本发明方法能够利用多组低阶实多项式求根来实现，相比传统算法计算量更低，估计的结果更稳定，且可适用于任意波束形成矩阵。

附图说明

下面结合附图和具体实施方式对本发明作进一步详细说明。

图1是本发明的一种基于root-MUSIC算法的目标波达方向获取方法流程图；

图2是不同阶数下导向矢量近似误差图；

图3是理想情况下多种算法估计目标角度RMSE对比图；

图4是锥削情况下多种算法估计目标角度RMSE对比图。

具体实施方式

参照图1，为本发明的一种基于root-MUSIC算法的目标波达方向获取方法流程图；所述基于root-MUSIC算法的目标波达方向获取方法，包括以下步骤：

步骤1，确定雷达，所述雷达包含N个阵元，并确定雷达检测范围内的波束形成矩阵B，所述雷达检测范围内的波束形成矩阵B中包含N_s个目标信号，下标s代表信号signal，且所述雷达检测范围内的波束形成矩阵B包含N'_b个列矢量，其中每一个列矢量都是为了将雷达中对应波束指向的空间角度而设置的权值矢量，表示N×N_b阶的复数域矩阵，∈表示属于，N_b表示波束形成矩阵B包含的波束个数，下标b代表波束beam；N≥N_b≥N_s。

然后确定雷达检测范围内的波束形成矩阵的空间频谱区域范围为

[-N_b/N N_b/N]，并将该空间频谱区域范围[-N_b/N N_b/N]以间隔ε均匀划分为M个空间频谱区域，并顺序排列，其中M＝2N_b-1，N、N_b、N_s、M、ε分别为大于0的整数。

步骤2，选取顺序排列后M个空间频谱区域中的第1个空间频谱区域，并利用最小二乘方法计算得到第1个空间频谱区域的窗矢量矩阵T，

T＝[t₁,t₂,…t_L+1]，T∈C^{(2N-1)×(L+1)}，C^{(2N-1)×(L+1)}表示(2N-1)×(L+1)阶的复数域矩阵，L表示设定的多项式最高阶数，L≥3；进而计算得到第1个空间频谱区域的窗矢量矩阵估计

步骤2的具体实现子步骤如下：

2.1选取顺序排列后M个空间频谱区域中的第1个空间频谱区域，所述第1个空间频谱区域包含若干个频点，选取其中Q个频点，并分别将第q个频点的频率记为f_q，将第q个频点的频率对应的扩展导向矢量记为上标T表示转置，N表示雷达包含的阵元个数，e表示指数函数。

2.2利用最小二乘法计算得到第q个频点的频率对应的扩展导向矢量的低阶实多项式近似，即T∈C^{(2N-1)×(L+1)}，C^{(2N-1)×(L+1)}表示(2N-1)×(L+1)阶的复数域矩阵，T为第1个空间频谱区域的窗矢量矩阵，是∑ζ(f_q)取最小值时对应的矢量矩阵；∑ζ(f_q)为第q个频点的频率对应的扩展导向矢量的低阶实多项式的近似误差，f_q为第q个频点的频率对应的扩展导向矢量的L阶多项式矢量，L表示设定的多项式最高阶数，L≥3，f_q为第q个频点的频率，q∈{1,2,…,Q}，Q为第1个空间频谱区域中选取的频点个数。

2.3计算得到第1个空间频谱区域的窗矢量矩阵估计F＝[f₁,…,f_q,…,f_Q]，q∈{1,2,…,Q}，Q为第1个空间频谱区域中选取的频点个数，上标T表示转置，上标-1表示求逆操作。

所述第1个空间频谱区域的窗矢量矩阵估计包含L+1个列矢量，

分别记为t₁,t₂,…,t_h,…,t_L+1，h∈{1,2,…,L+1}，t_h表示第1个空间频谱区域的窗矢量矩阵估计中第h个列矢量，L表示设定的多项式最高阶数。

步骤3，确定雷达接收到的回波数据为X，将雷达接收到的回波数据中的第k个采样数据记为X(k)，该雷达接收到的回波数据中的第k个采样数据X(k)为N×1维列矢量，k∈{1,2,3,…,K}，K为雷达接收到的回波数据X包含的采样点数；然后计算得到雷达回波数据的采样协方差矩阵估计上标H表示共轭转置操作。

对雷达回波数据的采样协方差矩阵估计进行特征分解，分别得到N”个特征值和N”个特征值各自对应的特征矢量；将N”个特征值按从大到小的顺序排列，并分别将前N_s个大的特征值作为雷达回波数据的采样协方差矩阵估计的对角矩阵将前N_s个大的特征值各自对应的特征矢量作为雷达回波数据的采样协方差矩阵估计的信号子空间矩阵将N”-N_s个特征值作为雷达回波数据的采样协方差矩阵估计的噪声子空间对角矩阵将N”-N_s个特征值各自对应的特征矢量作为雷达回波数据的采样协方差矩阵估计的噪声子空间矩阵即雷达回波数据的采样协方差矩阵估计进行特征分解后满足：→表示对雷达回波数据的采样协方差矩阵估计进行特征值分解；所述M个空间频谱间隔中包含的目标信号总个数与N”个特征值按从大到小的顺序排列后选取的大的特征值个数相等。

步骤4，利用雷达检测范围内的波束形成矩阵B和雷达回波数据的采样协方差矩阵估计的噪声子空间矩阵计算得到雷达检测范围内的波束形成矩阵的噪声协方差矩阵V，为雷达回波数据的采样协方差矩阵估计的噪声子空间矩阵，上标H表示共轭转置操作。

步骤5，将雷达检测范围内的波束形成矩阵的噪声协方差矩阵V中包含的所有的对角元素构造列矢量v，所述列矢量v为中心对称的厄尔米特Hermitian列矢量，且将列矢量v中第k个元素记为v_k，k∈{1,…,N'}，N'为列矢量v包含的元素个数；列矢量v包含的元素个数与雷达包含的阵元个数相等。

步骤6，确定M个空间频谱区域的变换矩阵W，W＝[w₁,w₂,...,w_m,...w_M]，w_m为第m个空间频谱区域的变换矢量，

w_m＝[e^{-j2π(N-1)(m-1)ε},e^{-j2π(N-2)(m-1)ε},…1,…e^{j2π(N-2)(m-1)ε},e^{j2π(N-1)(m-1)ε}]^T。

根据M个空间频谱区域的变换矩阵W、列矢量v和第1个空间频谱区域的窗矢量矩阵估计计算得到中间矩阵D，

l∈{1,2,…,L+1}，t_l表示第1个空间频谱区域的窗矢量矩阵估计中第l个列矢量，L表示设定的多项式最高阶数；1_L+1表示由L+1个1构成的行矢量，表示克罗内科积，⊙表示阿达玛积。

然后，利用M个空间频谱区域的变换矩阵W对中间矩阵D进行变换，得到M个空间频谱区域对应的M组低阶实多项式系数，其中第m组低阶实多项式系数为γ_m，

γ_m＝[γ_m1,γ_m2,…,γ_ml,…γ_m(L+1)]，γ_ml为第m组低阶实多项式系数中的第l阶系数，m∈{1,2,3,…,M}，

l∈{1,2,3,…,L+1}，L表示设定的多项式最高阶数，L≥3；上标T表示转置。

步骤7，根据M个空间频谱区域对应的M组低阶实多项式系数，对应得到M个方程，其中第m个方程为f^l-1为第l-1次幂的自变量，f⁰＝1，γ_ml为第m组低阶实多项式系数中的第l阶系数，求解第m个方程后得到L个共轭根，每个共轭根包含实部和虚部，然后在L个共轭根中选择虚部满足设定要求的N'_s个共轭根，所述设定要求为虚部接近于0，即虚部模值小于10^-6或更小的共轭根，该N'_s个共轭根各自对应的实部，分别对应为雷达检测范围内的波束形成矩阵B内N_s个目标信号各自角度的正弦值，并且所述N_s个目标信号各自空间角度的正弦值，分别对应为雷达检测范围内的波束形成矩阵B内N_s个目标信号各自的波达方向。

其中，雷达检测范围内的波束形成矩阵B中包含的目标信号个数与L个共轭根中虚部满足设定要求的共轭根个数相等且一一对应。

本发明的效果可通过以下仿真实验进一步说明。

(一)仿真实验数据说明

仿真中采用雷达阵元数N＝32的均匀线阵，理想情况下，波束形成利用N_b＝3的DFT矩阵，在波束域频带内有两个信源，角度分别为-2度和1度，单个信源SNR＝15dB，多项式阶数L＝5；锥削情况下采用-30dB的泰勒窗，其他参数不变

(二)仿真结果及分析

本发明的仿真结果由图2，图3和图4所示，图2是不同阶数下导向矢量近似误差图；图3是理想情况下多种算法估计目标角度RMSE对比图；图4是锥削情况下多种算法估计目标角度RMSE对比图。

图2表明误差随着阶数的增加而降低，当阶数大于4时，误差在可接受范围内；图3表明RMSE随着SNR增加而降低，本发明方法相比另外两种算法计算量更低，但是在所有SNR可以获得与另外两种算法几乎同样的性能。由图4可以看出，由于传统的波束域root-MUSIC在锥削情况下是不可实现的，此处不进行对比，图4表明本发明方法相比波束域频谱MUSIC算法计算量更低，但是在所有SNR都可以获得几乎同样的性能。

综上所述，仿真实验验证了本发明的正确性，有效性和可靠性。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围；这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。