CN110208733A

CN110208733A - 基于四阶累量的非圆信号阵列波达方向角估计方法

Info

Publication number: CN110208733A
Application number: CN201910285680.9A
Authority: CN
Inventors: 蔡晶晶; 刘应盼; 李鹏
Original assignee: Xidian University
Current assignee: Xidian University
Priority date: 2019-04-10
Filing date: 2019-04-10
Publication date: 2019-09-06
Anticipated expiration: 2039-04-10
Also published as: CN110208733B

Abstract

本发明公开了一种基于四阶累量的非圆信号均匀阵列波达方向角估计方法，主要解决现有技术中阵元利用率低，信号识别数量少的问题，其实现步骤是：1)构建一个非圆信号均匀线性阵列；2)获取非圆信号均匀线性阵列的输出信号；3)根据均匀线性阵列的输出信号计算所有四阶累量；4)根据所有四阶累量构造四阶累量协方差矩阵；5)计算四阶累量矩阵的噪声子空间；6)构造阵列导向矢量；7)根据四阶累量的噪声子空间和构造阵列导向矢量计算空间谱；8)根据空间谱绘制幅度谱图，得到波达方向角。本发明在阵元数量有限的情况下大大提高了阵列可识别的信源数目，可用于对飞机、舰船运动目标的侦察与无源定位。

Description

基于四阶累量的非圆信号阵列波达方向角估计方法

技术领域

本发明属于信号处理技术领域，特别涉及一种电磁信号的阵列信号波达方向角估计方法，可用于对飞机、舰船运动目标的侦察与无源定位。

背景技术

信号的波达方向角DOA估计是阵列信号处理领域的一个重要分支，它是指利用天线阵列对空间声学信号、电磁信号进行感应接收，再运用现代信号处理方法快速准确的估计出信号源的方向，在雷达、声纳、无线通信等领域具有重要应用价值。随着科技的不断进步，对阵列在进行信号波达方向估计时达到的自由度也有越来越高的要求。

针对该问题的研究中，出现较早、应用较为广泛的是多重信号分类MUSIC子空间的模型，对于一个L阵元的典型线性均匀阵列，传统的MUSIC类计算方法可检测的信源数目是L-1个。之后的大部分算法都是利用该模型生成的，例如信号参数估计旋转不变技术ESPRIT。这些算法由于采用典型的线性均匀阵列，造成估计的信号数目低于阵元数目，当目标个数很多时甚至无法识别，导致目标捕获失败。

目前实际使用的信号多为非圆信号，利用非圆信号的性质可在有限阵元条件下，提高角度自由度，检测更多的信源。Feifei Gao等人在其发表的论文“Improved MUSICUnder the Coexistence of Both Circular and Noncircular Sources”(《IEEETRANSACTIONS ON SIGNAL PROCESSING,VOL.56,NO.7,JULY 2008)中公开了一种基于二阶累量非圆信号阵列的DOA估计方法，该方法能够使用M个阵元，利用M个阵元输出信号之间的二阶累积量，可得到多个虚拟阵元。该方法具有估计多于阵元数目的信号数的能力，但是，该方法仍然存在的不足之处是，在阵元数量一定的情况下，能估计的信号数量有时仍不能满足实际应用。

为解决上述问题，近年来，出现了使用四阶累量的波达角度估计方法，可使信号估计数量得到较大提升，这种方法利用输出信号之间共轭乘积和相乘的方式获得不同的组合来实现更多的信号估计数量，但这种方法对于非圆信号均匀阵列而言，由于没有合理估计阵元输出信号而不能使阵元达到最大利用率。而实际应用中，在给定一定数量阵元的情况下，如果不能合理利用这些阵元获得足够多的虚拟阵元，就不能估计足够多的信号，将造成侦察和定位资源的浪费。

发明内容

本发明的目的在于针对上述现有技术存在的不足，提出一种基于四阶累量的非圆信号均匀线性阵列波达方向角估计方法，为了在有限的阵元数量条件下，提高能够进行估计的信号数量，减小因不能合理利用阵元造成的资源浪费。

为实现上述目的，本发明技术方案包括如下：

1.一种基于四阶累量的非圆信号波达方向角估计方法，其特征在于，包括如下：

(1)用M个天线接收机形成非圆信号均匀线性阵列L，将每个天线接收机称为一个阵元，其阵元间距为d，定义非圆信号均匀线性阵列L的第一个阵元为起始阵元，其中，M≥1，0＜d≤λ/2，λ为入射到阵列的窄带信号波长；

(2)由非圆信号均匀线性阵列L中天线接收机对空间目标信号进行采样，得到均匀线性阵列输出信号：Y(t)＝[y₁(t)，y₂(t)，…，y_i(t)，…，y_M(t)]，其中y_i(t)表示第i个阵元的输出信号，1≤i≤M：

(3)计算非圆信号均匀线性阵列输出信号Y(t)的两组四阶累量；

(3a)对非圆信号均匀线性阵列输出信号Y(t)求解第一组四阶累量：

c₁(k1，k2，k3，k4)＝cum(y_k1(t)，y_k2(t)^*，y_k3(t)，y_k4(t)^*)，

(3b)对非圆信号均匀线性阵列输出信号Y(t)求解第二组四阶累量：

c₂(k1，k2，k3，k4)＝cum(y_k1(t)，y_k2(t)^*，y_k3(t)^*，y_k4(t)^*)，

其中，k1，k2，k3，k4是非圆信号均匀线性阵列中阵元的位置，y_k1(t)，y_k2(t)，y_k3(t)，y_k4(t)是非圆信号均匀线性阵列输出信号Y(t)的输出信号，1≤k1，k2，k3，k4≤M，cum表示求解四阶累量操作，(·)^*表示向量的共轭运算；

(4)根据(3)中的非圆信号均匀线性阵列输出信号Y(t)的两组四阶累量，分别构造出对应的四阶累量矩阵；

(4a)定义中间变量D＝M-1，从第一组c₁(k1，k2，k3，k4)所有四阶累量中依次找出k1，k2，k3，k4满足下述条件的四阶累量：

k1+k3-k2-k4＝-2D，...，-1，0，1，...，2D，

将这些四阶累量依次定义为四阶累量元素g(-2D)，...，g(-1)，g(0)，g(1)，...，g(2D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₁；

(4b)从第二组c₂(k1，k2，k3，k4)所有四阶累量中依次找出k1，k2，k3，k4满足下述条件的四阶累量：

k1-k2-k3-k4＝-3D，...，-1，0，1，...，D，

将这些四阶累量依次定义为四阶累量元素g(-3D)，...，g(-1)，g(0)，g(1)，...，g(D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₂：

(5)利用(4)中得到的两组四阶累量矩阵G₁、G₂，生成最终四阶累量协方差矩阵G：

其中G₂ ^H为G₂的共轭转置矩阵；

(6)计算四阶累量矩阵G的噪声子空间U_n；

(7)构造阵列导向矢量α(θ)；

(8)根据阵列导向矢量α(θ)和四阶累量矩阵G的噪声子空间U_n，计算空间谱P(θ)；

(9)以波达方向角范围θ的值为x轴坐标，以空间谱P(θ)的幅度值为y轴坐标，绘制幅度谱图，从该幅度谱图中按照从高到低的顺序寻找幅值前K个谱峰，K表示均匀线性阵列的空间目标信号个数，K≥1，这些谱峰的峰值点所对应的x轴坐标即为目标的波达方向角度值。

本发明与现有技术相比具有以下优点：

1)本发明采用了四阶累积阵列模型进行波达方向角度估计，克服了现有技术中采用典型的线性均匀阵列造成估计的信号数目低于阵元数目的缺点，提高了在阵元数目相同的条件下的阵列可识别信源数目。

2)本发明将四阶累量应用到非圆信号均匀线性阵列的DOA估计中，通过使用四阶累量，设置阵元位置，可使用M个阵元获得2(4M-3)个连续虚拟阵元，大大提高了阵列利用率，同时通过获得更多的虚拟阵元，进一步增加了阵列可识别的信源数目。

3)本发明采用了非圆信号均匀线性阵列模型进行波达方向角度估计，相比于使用其他新型阵列模型对阵元的数目要求更低，提高了阵元数目使用的灵活性。

附图说明

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

具体实施方式

以下参照附图，对本发明的技术方案和效果作进一步的详细说明。

参附图1，本发明的实现步骤如下：

步骤1：构造非圆信号均匀线性阵列。

将M个天线接收机作为M个阵元，每个阵元间距为d，形成非圆信号均匀线性阵列L，定义非圆信号均匀线性阵列L的第一个阵元为起始阵元，位置为0d，最后一个阵元的位置为(M-1)d，其中，M≥1，0＜d≤λ/2，λ为入射到阵列的窄带信号波长。

步骤2：获得非圆信号均匀线性阵列输出信号。

将波长为λ的窄带信号入射到非圆信号均匀线性阵列L中，由非圆信号均匀线性阵列L中天线接收机对空间目标信号进行采样，得到非圆信号均匀线性阵列输出信号：

Y(t)＝[y₁(t)，…，y_i(t)，…，y_M(t)]，

其中y_i(t)表示第i个阵元的输出信号，1≤i≤M；

步骤3：计算均匀线性阵列输出信号的两组四阶累量。

(3a)对非圆信号均匀线性阵列输出信号Y(t)求解两组四阶累量c₁(k1，k2，k3，k4)、c₂(k1，k2，k3，k4)：

c₁(k1，k2，k3，k4)＝cum(y_k1(t)，y_k2 ^*(t)，y_k3(t)，y_k4 ^*(t))，

c₂(k1，k2，k3，k4)＝cum(y_k1(t)，y_k2 ^*(t)，y_k3 ^*(t)，y_k4 ^*(t))，

(3b)对(3a)中计算非圆信号均匀线性阵列输出信号Y(t)的四阶累量操作cum按如下方式进行：

假设A(t)，B(t)，C(t)，D(t)为所要求的四阶累量信号，则使用cum计算A(t)，B(t)，C(t)，D(t)得到的四阶累量为：

其中，T为信号采集快拍数，·表示乘积运算，(g)^*表示向量的共轭运算，表示求和运算。

步骤4：根据上述求得的两组四阶累量c₁(k1，k2，k3，k4)、c₂(k1，k2，k3，k4)构造四阶累量矩阵G。

(4a)形成四阶累量矩阵G₁：

定义中间变量D＝M-1，从第一组c₁(k1，k2，k3，k4)所有四阶累量中依次找出k1，k2，k3，k4满足下述条件的四阶累量：

k1+k3-k2-k4＝-2D，...，-1，0，1，...，2D，

将这些四阶累量依次定义为如下四阶累量元素：g(-2D)，...，g(-1)，g(0)，g(1)，...，g(2D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₁：

(4b)形成四阶累量矩阵G₂：

从第二组c₂(k1，k2，k3，k4所有四阶累量中依次找出k1，k2，k3，k4满足下述条件的四阶累量：

k1-k2-k3-k4＝-3D，...，-1，0，1，...，D，

将这些四阶累量依次定义为四阶累量元素：g(-3D)，...，g(-1)，g(0)，g(1)，...，g(D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₂：

(4c)利用得到的两组四阶累量矩阵G₁、G₂生成最终四阶累量协方差矩阵G：

其中G₂ ^H为G₂的共轭转置矩阵。

步骤5：计算四阶累量矩阵的噪声子空间U_n。

(5a)对四阶累量矩阵G进行如下特征分解：

G＝U·Λ.U^H，

其中，Λ为四阶累量矩阵(·)^H的特征值矩阵，U为四阶累量矩阵G的特征值所对应的特征向量矩阵，(·)^H表示矩阵的共轭转置运算；

(5b)将特征值矩阵Λ中的特征值按从大到小排序，取后2(2D+1)-K个较小特征值对应的特征向量矩阵作为噪声子空间：U_n＝[U_K+1，U_K+2，…，U_2(2D+1)]，其中U_K+1到U_2(2D+1)分别为第K+1到第2(2D+1)位置特征值对应的特征向量矩阵，K表示均匀线性阵列的空间目标信号个数。

步骤6：构造阵列导向矢量α(θ)。

构造的阵列导向矢量α(θ)，按如下公式进行：

其中，Z是(2D+1)×1维元素全部为0的向量，α₁(θ)和α₂(θ)按如下公式计算：

(·)^T表示矩阵转置运算，j为虚数单位；

步骤7：计算空间谱P(θ)。

(8a)将观测空域[-90°，90°]等间隔划分成E个角度，定义为波达方向角范围：

θ＝[θ₁，θ₂，…，θ_e，…θ_E]，

其中，θ_e表示第e个目标波达方向角，e＝1，2，...，E，E＞＞M；

(8b)计算空间谱P(θ)：

当波达方向角范围θ＝[θ₁，θ₂，…，θ_e，…θ_E]时，对应的空间谱：

P(θ)＝[P(θ₁)，P(θ₂)，…，P(θ_e)，…P(θ_E)]

其中，α(θ)表示阵列导向矢量，P(θ_e)表示第e个目标空间谱，(g)^H表示矩阵的共轭转置运算，U_n表示四阶累量矩阵G的噪声子空间，det表示求矩阵的行列式运算。

步骤8：绘制幅度谱图。

(8a)以波达方向角范围θ的值为x轴坐标，以空间谱P(θ)的幅度值为y轴坐标，绘制幅度谱图；

(8b)从绘制的幅度谱图中，按照从高到低的顺序将所有谱峰幅度值进行排序，将前K个谱峰的峰值点所对应的x轴坐标做为目标的波达方向角度值，其中，K表示均匀线性阵列的空间目标信号个数，且假设空间目标信号在传播过程中加入了均值为零的复高斯白噪声，K≥1。

下面结合仿真实例对本发明的效果做进一步的描述。

仿真实例

将M个天线接收机作为M个阵元，设M＝3，每个阵元间距为d。构建非圆信号均匀阵列四阶累量,仿真步骤如下：

(1.1)非圆信号均匀线性阵列L阵元位置为[0,1,2]d。

(1.2)利用非圆信号均匀线性阵列L计算非圆信号均匀线性阵列所有二阶累量的虚拟阵元位置为[-2,-1,0,1,2]d、[0,1,2,3,4]d、[-4,-3,-2,-1,0]d。

(1.3)根据(1.2)中的得到的二阶累量，计算得到所有四阶累量的虚拟阵元位置，结果如表1，

表1基于四阶累量的非圆信号均匀线性阵列虚拟阵元

由表1可知，基于四阶累量的非圆信号均匀线性阵列结构在阵元数量一定的情况下，可以得到更多的阵元位置信息，从而增加阵列可识别信源数目，同时，基于四阶累量的非圆信号均匀线性阵列结构相比于使用其他新型阵列模型对阵元的数目要求更低，提高了阵元数目使用的灵活性。

综上，本发明解决了现有技术阵元利用率低，识别信源数目少，无源定位估计误差大的问题，降低了对阵元数目的要求，保证了阵元数目使用的灵活性，提高了阵列可识别的信源数目以及低信噪比下对信号方向角的估计性能。

Claims

(1)用M个天线接收机形成非圆信号均匀线性阵列L，将每个天线接收机称为一个阵元，其阵元间距为d，定义非圆信号均匀线性阵列L的第一个阵元为起始阵元，其中，M≥1,0＜d≤λ/2，λ为入射到阵列的窄带信号波长；

(2)由非圆信号均匀线性阵列L中天线接收机对空间目标信号进行采样，得到均匀线性阵列输出信号：Y(t)＝[y₁(t),y₂(t),…,y_i(t),…,y_M(t)]，其中y_i(t)表示第i个阵元的输出信号，1≤i≤M；

(3)计算非圆信号均匀线性阵列输出信号Y(t)的两组四阶累量；

c₁(k1,k2,k3,k4)＝cum(y_k1(t),y_k2 ^*(t),y_k3(t),y_k4 ^*(t))，

c₂(k1,k2,k3,k4)＝cum(y_k1(t),y_k2 ^*(t),y_k3 ^*(t),y_k4 ^*(t))，

其中，k1,k2,k3,k4是非圆信号均匀线性阵列中阵元的位置，y_k1(t),y_k2(t),y_k3(t),y_k4(t)是非圆信号均匀线性阵列输出信号Y(t)的输出信号，1≤k1,k2,k3,k4≤M，cum表示求解四阶累量操作，(·)^*表示向量的共轭运算；

(4a)定义中间变量D＝M-1，从第一组c₁(k1,k2,k3,k4)所有四阶累量中依次找出k1,k2,k3,k4满足下述条件的四阶累量：

k1+k3-k2-k4＝-2D,...,-1,0,1,...,2D，

将这些四阶累量依次定义为四阶累量元素g(-2D),...,g(-1),g(0),g(1),...,g(2D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₁；

(4b)从第二组c₂(k1,k2,k3,k4)所有四阶累量中依次找出k1,k2,k3,k4满足下述条件的四阶累量：

k1-k2-k3-k4＝-3D,...,-1,0,1,...,D，

将这些四阶累量依次定义为四阶累量元素g(-3D),...,g(-1),g(0),g(1),...,g(D)，再将这些四阶累量元素重新排列形成四阶累量矩阵G₂：

其中G₂ ^H为G₂的共轭转置矩阵；

(6)计算四阶累量矩阵G的噪声子空间U_n；

(7)构造阵列导向矢量α(θ)；

2.根据权利要求1所述的方法，其中(3)中计算非圆信号均匀线性阵列输出信号Y(t)的四阶累量操作cum按如下方式进行：

假设A(t),B(t),C(t),D(t)为所要求的四阶累量信号，则使用cum计算A(t),B(t),C(t),D(t)得到的四阶累量为：

其中，T为信号采集快拍数，·表示乘积运算，(·)^*表示向量的共轭运算，表示求和运算。

3.根据权利要求1所述的方法，其中(6)中计算四阶累量矩阵G的噪声子空间U_n，按如下步骤进行：

(6a)对四阶累量矩阵G进行如下特征分解：

G＝U·Λ·U^H，

(6b)将特征值矩阵Λ中的特征值按从大到小排序，取后2(2D+1)-K个较小特征值对应的特征向量矩阵作为噪声子空间U_n,K表示均匀线性阵列的空间目标信号个数。

4.根据权利要求1所述的方法，其中(7)中构造的阵列导向矢量α(θ)，按如下公式进行：

其中,Z是(2D+1)×1维的元素全部为0的向量，α₁(θ)和α₂(θ)按如下公式计算：

(·)^T表示矩阵转置运算，j为虚数单位。

5.根据权利要求1所述的方法，其中(8)中计算空间谱P(θ)，按如下步骤进行：

(8a)将观测空域[-90°,90°]等间隔划分成E个角度,定义为波达方向角范围：

θ＝[θ₁,θ₂,…,θ_e,…θ_E]，

其中，θ_e表示第e个目标波达方向角，e＝1,2,...,E，E＞＞M；

(8b)计算空间谱P(θ)：

当波达方向角范围θ＝[θ₁,θ₂,…,θ_e,…θ_E]时，对应的空间谱：

P(θ)＝[P(θ₁),P(θ₂),…,P(θ_e),…P(θ_E)]

其中，α(θ)表示阵列导向矢量，P(θ_e)示第e个目标空间谱，(·)^H表示矩阵的共轭转置运算，U_n表示四阶累量矩阵G的噪声子空间，det表示求矩阵的行列式运算。