CN110824415A - 一种基于多发多收阵列的稀疏波达方向角度估计方法 - Google Patents

一种基于多发多收阵列的稀疏波达方向角度估计方法 Download PDF

Info

Publication number
CN110824415A
CN110824415A CN201911144835.3A CN201911144835A CN110824415A CN 110824415 A CN110824415 A CN 110824415A CN 201911144835 A CN201911144835 A CN 201911144835A CN 110824415 A CN110824415 A CN 110824415A
Authority
CN
China
Prior art keywords
array
matrix
target
iteration
value
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201911144835.3A
Other languages
English (en)
Other versions
CN110824415B (zh
Inventor
张新禹
刘永祥
姜卫东
霍凯
黎湘
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
National University of Defense Technology
Original Assignee
National University of Defense Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by National University of Defense Technology filed Critical National University of Defense Technology
Priority to CN201911144835.3A priority Critical patent/CN110824415B/zh
Publication of CN110824415A publication Critical patent/CN110824415A/zh
Application granted granted Critical
Publication of CN110824415B publication Critical patent/CN110824415B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/02—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using radio waves
    • G01S3/14—Systems for determining direction or deviation from predetermined direction
    • G01S3/143—Systems for determining direction or deviation from predetermined direction by vectorial combination of signals derived from differently oriented antennae

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明属于阵列信号处理领域，具体涉及一种基于多发多收阵列的稀疏波达方向角度估计方法。包括以下步骤：S1建立信号模型；S2信号协方差矩阵的估计与向量化；S3利用互质阵处理方法进行空间自由度拓展；S4采用近邻算子法来得到目标波达方向估计；S5确定目标来波方向和数量。本发明取得的有益效果为：通过本发明可提升低信噪比，小样本条件下多目标角度估计精度，并且本发明估计目标角度的迭代收敛速度快，运算效率高，可应用于实时多目标角度估计系统，具有重要的工程应用价值。

Description

一种基于多发多收阵列的稀疏波达方向角度估计方法
技术领域
本发明属于阵列信号处理领域，具体涉及一种基于多发多收(Multiple InputMultiple Output,MIMO)阵列的稀疏波达方向角度估计方法。
背景技术
目标角度估计(Direction of Arrival Estimation,DOAE)是阵列信号处理领域中非常重要的一个问题，相关的算法可广泛应用于雷达探测，声纳导航，多信道通信等领域。传统的目标角度估计算法大多是基于子空间信号流型分析或者属于最大似然估计类算法。但这两类算法往往分辨率有限，且受目标信号间相关性，采集样本数以及系统信噪比等因素影响较大，因此难以适应实际中各类复杂的应用环境。
随着稀疏恢复算法的引入，目标角度估计技术得到了很大的提升。稀疏恢复算法将传统的目标角度估计模型等效成字典学习模型，通过各类压缩感知算法得到更为精确的角度估计值。采用基于稀疏恢复算法的目标角度估计可以使得目标角度的获取对采集样本数以及信噪比的容忍度都大大提升，因此更有利于实际的应用。基于稀疏恢复算法的目标角度估计方法大致可以分为两类：一类采用了稀疏贝叶斯学习的框架，通过对超参数的估计来获得目标角度估计，这类方法依赖对模型变量分布的先验知识可以得到较好的恢复性能，但由于稀疏贝叶斯学习需要进行较大规模的矩阵运算，算法的运算量随着所需要的估计精度成指数倍增长，这一点又限制了该类算法的应用；另一类稀疏恢复方法通过引入基于lp,p＞0范数的稀疏性约束条件，并采用凸优化的方法得到信号的稀疏恢复结果，这一类方法的典型算法包括正交匹配法，基追寻法，LASSO算法等，这些传统算法的恢复结果要么容易收敛到较多的局部极值点，要么恢复结果难以达到要求的稀疏性。
传统的稀疏恢复算法依赖于字典矩阵的构建，好的字典矩阵往往可以得到较好的恢复性能。然而稀疏恢复的结果只能从字典矩阵所对应的网格中挑选，因此稀疏恢复中的一个重要问题是如何解决目标不在字典矩阵网格上的问题。对于稀疏贝叶斯学习类的方法，目标真实角度与相邻网格角度的偏移量被当成是待学习的参量，通过参数学习的方式来估计得到，但这类算法往往会大量提高稀疏贝叶斯学习方法的运算量，导致其难以在实际应用中施行。对于基于lp,p＞0范数的稀疏恢复方法，现有的方法往往采用缩小网格粒度的方法，这种做法也会大大提高算法的复杂度，导致无法在实际系统中应用。
随着阵列设计技术的发展，互质阵列应运而生。它采用两组不同间隔的线阵(间隔分别为m和n，m和n互质)，计算其接收信号的二阶统计量，再通过矢量化二阶统计矩阵的方法可以形成远超阵元个数的虚拟线阵，从而大大提高了阵列信号处理的自由度，为目标角度估计提供了更广阔的应用空间。
而传统的联合稀疏恢复方法(Zhao Tan and Arye Nehorai,“Sparse Directionof Arrival Estimation Using Co-Prime Arrays with Off-Grid Targets,”in IEEESIGNAL PROCESSING LETTERS,vol.21,No.1,2014,pp.26–29)无法利用MIMO阵列的多组发射信号，从而导致波达方向角度估计精度存在损失，在低信噪比，小样本数和较多目标数的条件下，难以满足工程实际的估计精度需求的问题。
发明内容
本发明要解决的技术问题是传统的稀疏恢复方法无法利用MIMO阵列的多组发射信号，从而导致波达方向角度估计精度存在损失，在低信噪比，小样本数和较多目标数的条件下，难以满足工程实际的估计精度需求的问题。
本发明的思路是针对MIMO阵列具有多组正交发射信号这一模型特点，提出基于多发多收阵列的新型稀疏波达方向角度估计方法。该方法首先通过互质阵孔径拓展方法达到阵列处理自由度拓展的目的，提高该算法可同时估计的目标信号数量；再通过近邻算子的推导得出适合于MIMO阵列角度估计模型的迭代估计算法；由于本算法采用了利用MIMO阵列中多组正交发射信号照射目标所得到的接收信号进行联合角度估计的方法，因此可以得到比传统方法更高的角度估计精度。
本发明采取的技术方案为，一种基于多发多收阵列的稀疏波达方向角度估计方法，包括以下步骤：
S1建立信号模型
考虑一个由G个发射天线和M个接收天线组成的MIMO阵列，收到来自不同方向θ＝[θ1 θ2 … θK]T的K个目标的回波信号，(·)T表示矩阵的转置。一般实际中用到的MIMO阵列，发射天线采用不同的正交编码，因此接收天线可以通过解码区分开来自不同发射天线的回波信号。于是，来自第i,i＝1,2…G个发射天线的阵列接收信号模型可以写为
yi(t)＝Axi(t)+ni(t) (1)
其中下标i用于指代某一序列或者矢量中的第i个元素，yi(t)表示t时刻接收阵列所接收到的来自第i个发射天线的信号，表示字典矩阵，
Figure BDA0002279448950000022
表示阵列的导向矢量，阵列的导向矢量为远场条件下的平面波入射到阵列所形成的接收信号矢量，它是一个复向量，维数等于阵列的阵元个数；字典矩阵A中的
Figure BDA0002279448950000023
组成了一个关于空间入射角的一种网格划分，N表示该划分网格点的总数，网格划分越细，目标角度估计的分辨率越高，一般的做法是将一个可能的空域区间做等间隔的划分，一般来说N的值要远大于接收天线个数M(例如N＝10M)。
接收天线阵列采用互质阵的排布方式，设阵元间最小间距为d，则该阵列在以下的位置处存在阵元：0:Pd:(Q-1)Pd，以及Qd:Qd:(2P-1)Qd，其中0:Pd:(Q-1)Pd表示从参考位置(即坐标为0的位置)每隔Pd放置一个阵元，直到(Q-1)P的位置，Q:Q:(2P-1)Q的定义与之相同，P和Q分别为互质的整数，且P<Q。对于这样的接收天线阵列，导向矢量
Figure BDA0002279448950000024
可以表示成：
Figure BDA0002279448950000025
上式的通项中h代表阵元的坐标位置，即h∈{0:P:(Q-1)P}∪{Q:Q:(2P-1)Q}，λ表示阵列发射的电磁波波长，π为圆周率常数，j表示虚数符号。
公式(1)中，xi(t)表示阵列接收到来自第i个发射天线的目标信号矢量，它的元素与字典矩阵中的每一列一一对应，当字典矩阵某一列所对应的入射角度存在目标时，与之相对应的xi(t)的该元素等于目标信号的复数值；当其所对应的角度不存在目标时，xi(t)的该元素值为0。ni(t)表示系统的加性噪声。这里我们假设在不同时刻的目标信号统计独立，这一假设在雷达信号处理中十分常见。在实际工程应用中，ni(t)通常为零均值的高斯白噪声，并且不同时刻间的噪声统计独立，即ni(t)～CN(0,σ0I)，这里ni(t)～CN(0,σ0I)表示ni(t)服从复数零均值的复高斯分布，其方差为σ0I，σ0表示噪声的功率，I表示单位矩阵。
S2信号协方差矩阵的估计与向量化
根据S1建立的信号模型，采用互质阵的信号处理方法来拓展阵列的自由度，更高的自由度意味着算法可以估计更多的目标角度。由于MIMO阵列具有G个不同的发射信号，因此一共可以形成G个不同的协方差矩阵和长矢量。按照下式计算得到第i个发射信号所对应的的协方差矩阵Ri
Figure BDA0002279448950000031
上面的公式中(·)H表示矩阵的共轭转置，L表示用来计算Ri所采用的信号时刻数，该值越大，根据上式所计算的协方差矩阵Ri越接近真实值。
通过对协方差矩阵向量化可以得到
ri＝vec(Ri) (4)
上式中vec(·)表示矩阵的向量化算子，它是将矩阵按照列由左至右的顺序依次堆叠成矢量。
S3利用互质阵处理方法进行空间自由度拓展
构造投影矩阵
Figure BDA0002279448950000032
(或者称为选择矩阵)(具体参照文献：P.Pal andP.Vaidyanathan,“Coprime sampling and the music algorithm,”in 2011IEEE DigitalSignal Processing Workshop and IEEE Signal Processing Education Workshop(DSP/SPE),2011,pp.289–294)，并按照下式计算得到投影后等效的线阵接收信号：
Figure BDA0002279448950000033
通过上面的投影运算，接收信号的维数实质上被降低了，即ri的维数小于
Figure BDA0002279448950000034
对于不同的发射阵元i，上式所采用的投影矩阵相同。
S4采用近邻算子法来得到目标波达方向估计
S4.1计算投影后等效的线性字典矩阵
Figure BDA0002279448950000041
不同于S1中的字典矩阵A，上式中
Figure BDA0002279448950000042
是投影后的字典矩阵，它等效于一个由均匀线阵所对应的字典矩阵。对于S1所建立信号模型中所描述的互质阵列排布，
Figure BDA0002279448950000043
所等效的均匀线阵维数等于2PQ+1，等效于一个从-PQd的位置到PQd的位置间隔为d的均匀线阵。因此，
Figure BDA0002279448950000044
中每一列是一个均匀线阵的导向矢量，其维数为2PQ+1，即
Figure BDA0002279448950000045
S1中的信号模型仅仅是针对当目标的实际波达方向落在字典矩阵的网格点上的情况，当实际波达方向不落在网格点上时，采用S1中的模型所得到的角度估计存在固定的系统误差，这一误差来源于字典矩阵所采用的角度网格
Figure BDA0002279448950000046
是离散的，因此S1中的模型所得到的角度估计只能是实际角度的一种近似。为此，本步骤通过泰勒展开公式来建立目标实际的导向矢量与由S1中网格化模型所得到的近似导向矢量之间的关系。对于字典矩阵
Figure BDA0002279448950000047
计算其一阶导数矩阵
Figure BDA0002279448950000048
Figure BDA0002279448950000049
此处引入一阶导数矩阵是为了解决目标不在字典矩阵网格上的问题。
Figure BDA00022794489500000411
的每一列可以按照下式计算得到
Figure BDA00022794489500000412
S4.2参数初始化：需要初始化的参数包括近邻算子系数μ，稀疏度约束系数
Figure BDA00022794489500000413
以及Lipschitz常量
Figure BDA00022794489500000414
近邻算子系数μ的选取需要结合实际的需要，由用户根据自己的系统环境来定。近邻算子系数μ越小，得到的波达方向角估计精度就越高，但越小的近邻算子系数就意味着算法需要迭代更多的次数，即更长的时间才能获得估计结果，一般的μ的值可以设置为
Figure BDA00022794489500000415
这里||·||F表示矩阵的Frobenius范数，表示由矩阵和矩阵
Figure BDA00022794489500000418
左右排列所形成的矩阵。
稀疏度约束系数
Figure BDA00022794489500000419
需要根据所估计目标的数量以及所设定网格点的疏密程度来定，一般来说，目标数量越小，
Figure BDA00022794489500000420
可以设的越大，在实际系统的应用中，该值可以设为10倍的噪声功率σ0。
Lipschitz常量的取值直接与μ有关，假设MIMO阵列发射了G组不同的正交信号，那么该常量的取值应当满足
Figure BDA00022794489500000422
值越大，算法收敛需要的时间越长。
S4.3通过迭代的方式得到目标的波达方向估计，需要迭代的中间变量包括zk和和zk是与信号直接相关的迭代变量，它们均是维数为2N的矢量，
Figure BDA0002279448950000053
zk和
Figure BDA0002279448950000054
中的下标k表示第k次迭代，
Figure BDA0002279448950000055
则代表第k次迭代更新前迭代步长，
Figure BDA0002279448950000056
和zk表示第k次迭代更新前的值，具体步骤如下：
S4.3.1初始化迭代变量：令初始的
Figure BDA0002279448950000057
这里12N表示一个元素全为1的向量，其维数为2N；令中间变量
Figure BDA0002279448950000058
迭代步长的初始值
Figure BDA0002279448950000059
迭代计数变量k＝1。
S4.3.2计算当前的误差梯度：令
Figure BDA00022794489500000510
它表示当前的z(信号估计值)所对应的模型与接收数据之间的误差。按照下式计算该函数的梯度：
Figure BDA00022794489500000511
这里的
Figure BDA00022794489500000512
(·)H表示矩阵的共轭转置。
S4.3.3计算稀疏度约束函数的梯度：定义稀疏度约束函数为h(zk)＝||zk||2,1，这里的
Figure BDA00022794489500000513
表示zk的联合范数，zk(i)和zk(i+N)分别表示zk的第i项和第i+N项，h(zk)的Moreau包络记作hμ(zk)，按照下面的表达式计算函数hμ(zk)的梯度
Figure BDA00022794489500000514
在上面的表达式中，proxμh(zk)表示函数hμ(zk)的近邻算子，它的第i项和第i+N项如下所示
Figure BDA00022794489500000516
上面的表达式中，proxμh(zk)i表示proxμh(zk)的第i项，proxμh(zk)i+N表示proxμh(zk)的第i+N项。依次令i＝1,2,…N，按照(12)进行计算就可以得到完整的proxμh(zk)。
S4.3.4更新通过下面的表达式来更新对
Figure BDA00022794489500000518
的估计
Figure BDA00022794489500000519
上式中的PF函数表示投影函数，其输入的自变量为
Figure BDA00022794489500000520
为了给出该函数的数学表达式，将其自变量拆分成前N个元素和后N个元素，并分别记为s和p，并且记s(i)和p(i)分别表示输入向量s和p的第i项，；那么根据s和p，函数PF输出结果的第i项和第i+N项可以写成：
Figure BDA0002279448950000061
上式中和
Figure BDA0002279448950000063
分别表示函数输出结果的第i项和第i+N项，r表示目标波达方向角度与最近网格点角度之间相差的最大绝对值，通常取值为网格点之间最小间隔d的一半，且单位为弧度(例如网格点之间最小间隔为1°，则
Figure BDA0002279448950000064
)，参数
Figure BDA0002279448950000065
|·|表示取绝对值。
S4.3.5更新
Figure BDA0002279448950000066
按照如下的表达式更新
Figure BDA0002279448950000068
上式中等号左边的表示第k次迭代更新后的值，将用于第k+1次迭代，等号右边的
Figure BDA00022794489500000610
表示更新前的值，在第一次迭代中
Figure BDA00022794489500000611
即为初始化的值。
S4.3.6更新zk+1：按照如下的表达式更新zk+1
Figure BDA00022794489500000612
这里等号左边的zk+1表示第k次迭代更新后的值，表示第k次迭代更新前的值，也即前一次迭代更新后的值。
S4.3.7收敛判断：在更新完zk+1之后，需要判断迭代是否需要继续，判断的条件如下：如果||zk+1-zk||1＞ε0并且k+1＜ε1，那么回到S4.3.2，令k＝k+1，重新开始迭代,；否则退出迭代，执行下一步；阈值参数ε0和ε1由用户决定，ε0越小，ε1越大，则得到的估计精度越高，但所花费的时间就越多；反之则估计精度降低，但花费的时间较少。通常ε0可以取为0.01，ε1可以取为2000。
S5确定目标来波方向和数量
目标来波方向和数量的确定通过设定门限ε2来进行判定得到。记算法输出结果zk+1的前N项为zf，后N项为zb。在算法的输出结果中，zf的物理含义是信号的幅度，它的各个元素与字典矩阵中的各列相对应，令zf(i)表示zf的第i项，那么zf(i)的值就代表了
Figure BDA0002279448950000072
方向上目标的幅度(如果没有目标，zf(i)的值会趋近于0)，zb是一个与目标实际波达方向相关的物理量，通过它可以求解出目标实际波达方向与网格点之间的差值。一般情况下，可以通过如下方法来确定来波方向和数量：如果zf(i)的值超过了门限ε2，则可判定zf(i)所对应的来波方向上存在目标，其波达方向角为
Figure BDA0002279448950000071
zb(i)表示zb的第i个元素；如果zf(i)的值低于门限ε2，则其所对应的来波方向上不存在目标。按照这种门限检测的方法在得到了目标来波方向的同时也就确定了目标的数量。通常情况下，门限ε2的取值为计算得到的zf的最大值一半。
当信噪比较低的情况下(如低于0dB的情况下)，为了获得准确的角度估计，则可以对向量zf进行峰值检测，认为有且只有在zf的峰值点处存在目标。
本发明取得的有益效果为：通过本发明可提升低信噪比，小样本条件下多目标角度估计精度，并且本发明估计目标角度的迭代收敛速度快，运算效率高，可应用于实时多目标角度估计系统，具有重要的工程应用价值。
附图说明
图1为本发明处理流程图；
图2为新方法与传统方法的空间功率谱图；
图3为较少采样样本数的条件下新方法与传统方法的性能随信噪比变化的比较；
图4为较多采样样本数的条件下新方法与传统方法的性能随信噪比变化的比较；
图5为较低信噪比的条件下新方法与传统方法的性能随采集样本数变化的比较；
图6为较高信噪比的条件下新方法与传统方法的性能随采集样本数变化的比较。
具体实施方式
下面结合附图对本发明进行进一步说明：
图1为本发明总处理流程。
本发明所述一种基于多发多收互质阵列的新型稀疏波达方向角度估计方法，包括以下步骤：
S1建立信号模型；
S2信号协方差矩阵的估计与向量化；
S3利用互质阵处理方法进行空间自由度拓展；
S4采用近邻算子法来得到目标波达方向估计；
S5确定目标来波方向和数量。
图2为本发明的方法与基于单一发射阵元的传统联合稀疏角度估计方法(ZhaoTan and Arye Nehorai,“Sparse Direction of Arrival Estimation Using Co-PrimeArrays with Off-Grid Targets,”in IEEE SIGNAL PROCESSING LETTERS,vol.21,No.1,2014,pp.26–29)以及MUSIC算法的空间功率谱图。该仿真基于一个10阵元的非均匀线阵，其阵元分布在0:Pd:(Q-1)Pd，以及Qd:Qd:(2P-1)Qd的位置上，其中
Figure BDA0002279448950000081
Q＝5，P＝3。考虑四个不相干的目标分别从20°，30°，-15°和-25°的位置处发射信号入射到阵列上。目标信号的信噪比均为5dB，阵列一共接收到500个采集样本。本发明所提出的方法以及传统的联合稀疏角度估计方法采用相同的字典矩阵，其网格点范围从-45°到45°，相邻网格点相差1°；MUSIC算法的搜索采用了穷举法，其搜索粒度为0.1°。由图可知，三种方法的空间谱都在目标入射的方向附近都出现了尖峰，因此通过极点检测的方法可以得到准确的角度估计结果；从图中可以看出，本发明所提出的算法以及传统的联合稀疏角度估计方法所得到的尖峰与目标实际方向较为匹配，而MUSIC算法的尖峰位置与实际目标的来波方向存在一定的偏差。
图3为三种算法的估计精度随着信噪比变化的情况，此时的采集样本数为100。估计精度由角度估计的均方根误差来表示。其表达式为
Figure BDA0002279448950000082
这里
Figure BDA0002279448950000083
表示第i次仿真得到的角度估计值，它是一个四个元素的矢量。由图可知，本发明算法在信噪比低至-5dB条件下仍然能够有效的得到目标角度的精确估计，而传统的联合稀疏角度估计方法在这一条件下则无法得到准确的估计，MUSIC算法虽然能够得到较为准确的估计，但估计误差较大。因此相比于传统的角度估计算法，本发明方法在信噪比较低时能得到更为精确的目标角度估计，估计精度提高较为显著。
图4为采集样本数改为500时，三种算法的估计精度随着信噪比变化的情况。由图可知，即使采样点数增加，本发明方法的估计精度相较于传统算法仍然有较为显著的提升。
图5为信噪比为0dB时，三种方法的估计精度随采集样本数变化的比较图。由图可知，三种算法的估计精度随着采样点数的增加而逐渐提高，但本发明所提出的方法相较于传统方法具有更好的估计精度，并且在采集样本数较少时，估计精度的提升更为显著。
图6为信噪比为5dB时，新方法与传统方法性能随采集样本数的变化情况。由图中可以看出，此时，新方法相较于传统方法在采集样本数不足时，估计精度有更加明显提升。
基于仿真的实验结果表明，本发明对噪声鲁棒性强，对小样本数据仍然有效，且角度估计精度明显高于传统方法，满足实际系统目标角度估计要求。本发明可在雷达回波数据质量受限条件下，实现多目标的入射角度精确估计，尤其为强对抗条件下导弹防御、空间目标监视中的空间目标识别提供了技术支持，工程应用价值高。

Claims (10)

1.一种基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于，该方法包括以下步骤：
S1建立信号模型
考虑一个由G个发射天线和M个接收天线组成的MIMO阵列，收到来自不同方向θ＝[θ1θ2…θK]T的K个目标的回波信号，(·)T表示矩阵的转置；一般实际中用到的MIMO阵列，发射天线采用不同的正交编码，因此接收天线可以通过解码区分开来自不同发射天线的回波信号；于是，来自第i,i＝1,2…G个发射天线的阵列接收信号模型可以写为
yi(t)＝Axi(t)+ni(t) (1)
其中下标i用于指代某一序列或者矢量中的第i个元素，yi(t)表示t时刻接收阵列所接收到的来自第i个发射天线的信号，表示字典矩阵，
Figure FDA0002279448940000012
表示阵列的导向矢量，阵列的导向矢量为远场条件下的平面波入射到阵列所形成的接收信号矢量，它是一个复向量，维数等于阵列的阵元个数；字典矩阵A中的
Figure FDA0002279448940000013
组成了一个关于空间入射角的一种网格划分，N表示该划分网格点的总数，网格划分越细，目标角度估计的分辨率越高；
接收天线阵列采用互质阵的排布方式，设阵元间最小间距为d，则该阵列在以下的位置处存在阵元：0:Pd:(Q-1)Pd，以及Qd:Qd:(2P-1)Qd，其中0:Pd:(Q-1)Pd表示从参考位置每隔Pd放置一个阵元，直到(Q-1)P的位置，Q:Q:(2P-1)Q的定义与之相同，P和Q分别为互质的整数，且P<Q；对于这样的接收天线阵列，导向矢量
Figure FDA0002279448940000014
可以表示成：
Figure FDA0002279448940000015
上式的通项中h代表阵元的坐标位置，即h∈{0:P:(Q-1)P}∪{Q:Q:(2P-1)Q}，λ表示阵列发射的电磁波波长，π为圆周率常数，j表示虚数符号；
公式(1)中，xi(t)表示阵列接收到来自第i个发射天线的目标信号矢量，它的元素与字典矩阵中的每一列一一对应，当字典矩阵某一列所对应的入射角度存在目标时，与之相对应的xi(t)的该元素等于目标信号的复数值；当其所对应的角度不存在目标时，xi(t)的该元素值为0；ni(t)表示系统的加性噪声；
S2信号协方差矩阵的估计与向量化
由于MIMO阵列具有G个不同的发射信号，因此一共可以形成G个不同的协方差矩阵和长矢量；按照下式计算得到第i个发射信号所对应的的协方差矩阵Ri
Figure FDA0002279448940000016
上面的公式中(·)H表示矩阵的共轭转置，L表示用来计算Ri所采用的信号时刻数，该值越大，根据上式所计算的协方差矩阵Ri越接近真实值；
通过对协方差矩阵向量化可以得到
ri＝vec(Ri) (4)
上式中vec(·)表示矩阵的向量化算子，它是将矩阵按照列由左至右的顺序依次堆叠成矢量；
S3利用互质阵处理方法进行空间自由度拓展
构造投影矩阵
Figure FDA0002279448940000021
并按照下式计算得到投影后等效的线阵接收信号：
Figure FDA0002279448940000022
通过上面的投影运算，接收信号的维数实质上被降低了，即ri的维数小于
Figure FDA00022794489400000217
对于不同的发射阵元i，上式所采用的投影矩阵
Figure FDA0002279448940000023
相同；
S4采用近邻算子法来得到目标波达方向估计
S4.1计算投影后等效的线性字典矩阵
不同于S1中的字典矩阵A，上式中
Figure FDA0002279448940000025
是投影后的字典矩阵，它等效于一个由均匀线阵所对应的字典矩阵；对于S1所建立信号模型中所描述的互质阵列排布，
Figure FDA0002279448940000026
所等效的均匀线阵维数等于2PQ+1，等效于一个从-PQd的位置到PQd的位置间隔为d的均匀线阵；因此，
Figure FDA0002279448940000027
中每一列是一个均匀线阵的导向矢量，其维数为2PQ+1，即
Figure FDA0002279448940000028
S1中的信号模型仅仅是针对当目标的实际波达方向落在字典矩阵的网格点上的情况，当实际波达方向不落在网格点上时，采用S1中的模型所得到的角度估计存在固定的系统误差，这一误差来源于字典矩阵所采用的角度网格
Figure FDA0002279448940000029
是离散的，因此S1中的模型所得到的角度估计只能是实际角度的一种近似；为此，本步骤通过泰勒展开公式来建立目标实际的导向矢量与由S1中网格化模型所得到的近似导向矢量之间的关系；对于字典矩阵
Figure FDA00022794489400000210
计算其一阶导数矩阵
Figure FDA00022794489400000211
Figure FDA00022794489400000212
此处引入一阶导数矩阵
Figure FDA00022794489400000213
是为了解决目标不在字典矩阵网格上的问题；
的每一列可以按照下式计算得到
Figure FDA00022794489400000215
S4.2参数初始化：需要初始化的参数包括近邻算子系数μ，稀疏度约束系数
Figure FDA00022794489400000216
以及Lipschitz常量
Figure FDA0002279448940000031
近邻算子系数μ的选取需要结合实际的需要，由用户根据自己的系统环境来定，近邻算子系数μ越小，得到的波达方向角估计精度就越高，但越小的近邻算子系数就意味着算法需要迭代更多的次数，即更长的时间才能获得估计结果；
稀疏度约束系数
Figure FDA0002279448940000032
需要根据所估计目标的数量以及所设定网格点的疏密程度来定，目标数量越小，
Figure FDA0002279448940000033
可以设的越大；
Lipschitz常量
Figure FDA0002279448940000034
的取值直接与μ有关，假设MIMO阵列发射了G组不同的正交信号，那么该常量的取值应当满足
Figure FDA0002279448940000035
Figure FDA0002279448940000036
值越大，算法收敛需要的时间越长；
S4.3通过迭代的方式得到目标的波达方向估计，需要迭代的中间变量包括
Figure FDA0002279448940000037
zk和和zk是与信号直接相关的迭代变量，它们均是维数为2N的矢量，
Figure FDA0002279448940000039
zk和中的下标k表示第k次迭代，
Figure FDA00022794489400000311
则代表第k次迭代更新前迭代步长，和zk表示第k次迭代更新前的值，具体步骤如下：
S4.3.1初始化迭代变量：令初始的这里12N表示一个元素全为1的向量，其维数为2N；令中间变量
Figure FDA00022794489400000314
迭代步长的初始值迭代计数变量k＝1；
S4.3.2计算当前的误差梯度：令它表示当前的z所对应的模型与接收数据之间的误差，按照下式计算该函数的梯度：
Figure FDA00022794489400000317
这里的
Figure FDA00022794489400000318
(·)H表示矩阵的共轭转置；
S4.3.3计算稀疏度约束函数的梯度：定义稀疏度约束函数为h(zk)＝||zk||2,1，这里的
Figure FDA00022794489400000319
表示zk的联合范数，zk(i)和zk(i+N)分别表示zk的第i项和第i+N项，h(zk)的Moreau包络记作hμ(zk)，按照下面的表达式计算函数hμ(zk)的梯度
Figure FDA00022794489400000320
在上面的表达式中，proxμh(zk)表示函数hμ(zk)的近邻算子，它的第i项和第i+N项如下所示
Figure FDA0002279448940000041
上面的表达式中，proxμh(zk)i表示proxμh(zk)的第i项，proxμh(zk)i+N表示proxμh(zk)的第i+N项；依次令i＝1,2,…N，按照(12)进行计算就可以得到完整的proxμh(zk)；
S4.3.4更新
Figure FDA0002279448940000042
通过下面的表达式来更新对
Figure FDA0002279448940000043
的估计
Figure FDA0002279448940000044
上式中的PF函数表示投影函数，其输入的自变量为为了给出该函数的数学表达式，将其自变量拆分成前N个元素和后N个元素，并分别记为s和p，并且记s(i)和p(i)分别表示输入向量s和p的第i项，；那么根据s和p，函数PF输出结果的第i项和第i+N项可以写成：
Figure FDA0002279448940000046
上式中
Figure FDA0002279448940000047
和
Figure FDA0002279448940000048
分别表示函数输出结果的第i项和第i+N项，r表示目标波达方向角度与最近网格点角度之间相差的最大绝对值，参数
Figure FDA0002279448940000049
|·|表示取绝对值；
S4.3.5更新按照如下的表达式更新
Figure FDA00022794489400000412
上式中等号左边的
Figure FDA00022794489400000413
表示第k次迭代更新后的值，将用于第k+1次迭代，等号右边的
Figure FDA00022794489400000414
表示更新前的值，在第一次迭代中即为初始化的值；
S4.3.6更新zk+1：按照如下的表达式更新zk+1
Figure FDA0002279448940000051
这里等号左边的zk+1表示第k次迭代更新后的值，
Figure FDA0002279448940000052
表示第k次迭代更新前的值，也即前一次迭代更新后的值；
S4.3.7收敛判断：在更新完zk+1之后，需要判断迭代是否需要继续，判断的条件如下：如果||zk+1-zk||1＞ε0并且k+1＜ε1，那么回到S4.3.2，令k＝k+1，重新开始迭代,；否则退出迭代，执行下一步；阈值参数ε0和ε1由用户决定，ε0越小，ε1越大，则得到的估计精度越高，但所花费的时间就越多；反之则估计精度降低，但花费的时间较少；
S5确定目标来波方向和数量
目标来波方向和数量的确定通过设定门限ε2来进行判定得到；记算法输出结果zk+1的前N项为zf，后N项为zb；在算法的输出结果中，zf的物理含义是信号的幅度，它的各个元素与字典矩阵中的各列相对应，令zf(i)表示zf的第i项，那么zf(i)的值就代表了
Figure FDA0002279448940000053
方向上目标的幅度，zb是一个与目标实际波达方向相关的物理量，通过它可以求解出目标实际波达方向与网格点之间的差值，通过如下方法来确定来波方向和数量：如果zf(i)的值超过了门限ε2，则可判定zf(i)所对应的来波方向上存在目标，其波达方向角为
Figure FDA0002279448940000054
zb(i)表示zb的第i个元素；如果zf(i)的值低于门限ε2，则不存在目标；按照这种门限检测的方法在得到了目标来波方向的同时也就确定了目标的数量。
2.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：S1做网格划分时，将一个可能的空域区间做等间隔的划分。
3.一种根据权利要求2所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：划分网格点的总数N的值要远大于接收天线个数M。
4.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：在实际工程应用中，ni(t)为零均值的高斯白噪声，并且不同时刻间的噪声统计独立，即ni(t)～CN(0,σ0I)，这里ni(t)～CN(0,σ0I)表示ni(t)服从复数零均值的复高斯分布，其方差为σ0I，σ0表示噪声的功率，I表示单位矩阵。
5.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：近邻算子系数μ的值设置为
Figure FDA0002279448940000055
这里||·||F表示矩阵的Frobenius范数，
Figure FDA0002279448940000056
表示由矩阵
Figure FDA0002279448940000057
和矩阵
Figure FDA0002279448940000058
左右排列所形成的矩阵。
6.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：稀疏度约束系数
Figure FDA0002279448940000061
的值设为10倍的噪声功率σ0。
7.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：目标波达方向角度与最近网格点角度之间相差的最大绝对值r取值为网格点之间最小间隔d的一半，单位为弧度。
8.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：门限ε0取为0.01，门限ε1取为2000。
9.一种根据权利要求1所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：门限ε2的取值为计算得到的zf的最大值一半。
10.一种根据权利要求1至9任一条所述基于多发多收阵列的稀疏波达方向角度估计方法，其特征在于：当信噪比较低的情况下，为了获得准确的角度估计，可以对向量zf进行峰值检测，认为有且只有在zf的峰值点处存在目标。
CN201911144835.3A 2019-11-19 2019-11-19 一种基于多发多收阵列的稀疏波达方向角度估计方法 Active CN110824415B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911144835.3A CN110824415B (zh) 2019-11-19 2019-11-19 一种基于多发多收阵列的稀疏波达方向角度估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911144835.3A CN110824415B (zh) 2019-11-19 2019-11-19 一种基于多发多收阵列的稀疏波达方向角度估计方法

Publications (2)

Publication Number Publication Date
CN110824415A true CN110824415A (zh) 2020-02-21
CN110824415B CN110824415B (zh) 2020-07-07

Family

ID=69557475

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911144835.3A Active CN110824415B (zh) 2019-11-19 2019-11-19 一种基于多发多收阵列的稀疏波达方向角度估计方法

Country Status (1)

Country Link
CN (1) CN110824415B (zh)

Cited By (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111431575A (zh) * 2020-03-04 2020-07-17 哈尔滨工程大学 基于常规波束形成的来波方向稀疏重构方法
CN111723686A (zh) * 2020-06-01 2020-09-29 中国电子科技集团公司第三十六研究所 一种信号个数估计方法、装置、电子设备和存储介质
CN112285647A (zh) * 2020-09-30 2021-01-29 中国船舶重工集团公司七五0试验场 一种基于稀疏表示与重构的信号方位高分辨估计方法
CN113050027A (zh) * 2021-03-08 2021-06-29 浙江大学 一种幅相误差情况下基于稀疏重构的波达方向估计方法
CN113093093A (zh) * 2021-04-07 2021-07-09 南京邮电大学 基于线性阵列波达方向估计的车辆定位方法
CN114019445A (zh) * 2021-09-22 2022-02-08 中国电子科技集团公司第二十九研究所 一种基于位置聚类动态稀疏重构的二维到达角度测量方法
CN114675257A (zh) * 2022-05-30 2022-06-28 中国人民解放军国防科技大学 一种非矩形阵列极化单脉冲雷达双目标分辨方法
CN115087881A (zh) * 2020-06-01 2022-09-20 华为技术有限公司 一种波达角aoa估计方法和装置
CN115267653A (zh) * 2022-07-29 2022-11-01 青岛理工大学 一种基于互质阵列的近场源位置估计方法
CN115656998A (zh) * 2022-11-22 2023-01-31 中国人民解放军空军预警学院 一种低样本数下阵列信号自适应检测方法与系统
CN115987337A (zh) * 2022-12-01 2023-04-18 东南大学 非自由空间中mimo系统的电磁有效自由度获取方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103399292A (zh) * 2013-07-22 2013-11-20 西安电子科技大学 一种基于软稀疏表示的doa估计方法
CN104977558A (zh) * 2015-06-16 2015-10-14 电子科技大学 一种基于贝叶斯压缩感知的分布源中心波达方向估计方法
CN109557503A (zh) * 2018-12-19 2019-04-02 成都理工大学 基于相关矩阵重建解相干的mimo互质阵列doa估计方法
CN109752710A (zh) * 2019-01-07 2019-05-14 中国人民解放军国防科技大学 一种基于稀疏贝叶斯学习的快速目标角度估计方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103399292A (zh) * 2013-07-22 2013-11-20 西安电子科技大学 一种基于软稀疏表示的doa估计方法
CN104977558A (zh) * 2015-06-16 2015-10-14 电子科技大学 一种基于贝叶斯压缩感知的分布源中心波达方向估计方法
CN109557503A (zh) * 2018-12-19 2019-04-02 成都理工大学 基于相关矩阵重建解相干的mimo互质阵列doa估计方法
CN109752710A (zh) * 2019-01-07 2019-05-14 中国人民解放军国防科技大学 一种基于稀疏贝叶斯学习的快速目标角度估计方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
崔琛 等: "稀疏阵列MIMO雷达高精度收发角度联合估计", 《应用科学学报》 *
陈蕾 等: "基于线性Bregman迭代的结构化噪声矩阵补全算法", 《计算机学报》 *

Cited By (18)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111431575B (zh) * 2020-03-04 2022-02-11 哈尔滨工程大学 基于常规波束形成的来波方向稀疏重构方法
CN111431575A (zh) * 2020-03-04 2020-07-17 哈尔滨工程大学 基于常规波束形成的来波方向稀疏重构方法
CN111723686A (zh) * 2020-06-01 2020-09-29 中国电子科技集团公司第三十六研究所 一种信号个数估计方法、装置、电子设备和存储介质
CN115087881A (zh) * 2020-06-01 2022-09-20 华为技术有限公司 一种波达角aoa估计方法和装置
CN111723686B (zh) * 2020-06-01 2023-09-22 中国电子科技集团公司第三十六研究所 一种信号个数估计方法、装置、电子设备和存储介质
CN112285647A (zh) * 2020-09-30 2021-01-29 中国船舶重工集团公司七五0试验场 一种基于稀疏表示与重构的信号方位高分辨估计方法
CN112285647B (zh) * 2020-09-30 2023-10-27 中国船舶重工集团公司七五0试验场 一种基于稀疏表示与重构的信号方位高分辨估计方法
CN113050027A (zh) * 2021-03-08 2021-06-29 浙江大学 一种幅相误差情况下基于稀疏重构的波达方向估计方法
CN113050027B (zh) * 2021-03-08 2023-09-19 浙江大学 一种幅相误差情况下基于稀疏重构的波达方向估计方法
CN113093093B (zh) * 2021-04-07 2023-08-18 南京邮电大学 基于线性阵列波达方向估计的车辆定位方法
CN113093093A (zh) * 2021-04-07 2021-07-09 南京邮电大学 基于线性阵列波达方向估计的车辆定位方法
CN114019445A (zh) * 2021-09-22 2022-02-08 中国电子科技集团公司第二十九研究所 一种基于位置聚类动态稀疏重构的二维到达角度测量方法
CN114675257B (zh) * 2022-05-30 2022-08-02 中国人民解放军国防科技大学 一种非矩形阵列极化单脉冲雷达双目标分辨方法
CN114675257A (zh) * 2022-05-30 2022-06-28 中国人民解放军国防科技大学 一种非矩形阵列极化单脉冲雷达双目标分辨方法
CN115267653A (zh) * 2022-07-29 2022-11-01 青岛理工大学 一种基于互质阵列的近场源位置估计方法
CN115267653B (zh) * 2022-07-29 2024-04-16 青岛理工大学 一种基于互质阵列的近场源位置估计方法
CN115656998A (zh) * 2022-11-22 2023-01-31 中国人民解放军空军预警学院 一种低样本数下阵列信号自适应检测方法与系统
CN115987337A (zh) * 2022-12-01 2023-04-18 东南大学 非自由空间中mimo系统的电磁有效自由度获取方法

Also Published As

Publication number Publication date
CN110824415B (zh) 2020-07-07

Similar Documents

Publication Publication Date Title
CN110824415B (zh) 一种基于多发多收阵列的稀疏波达方向角度估计方法
CN109932680B (zh) 一种基于平移互质阵列的非圆信号波达方向估计方法
CN106054123B (zh) 一种稀疏l阵及其二维doa估计方法
CN109738854B (zh) 一种天线阵列来波方向的到达角估计方法
CN106019213B (zh) 一种部分稀疏l阵及其二维doa估计方法
CN107450047B (zh) 嵌套阵下基于未知互耦信息的压缩感知doa估计方法
CN109655799A (zh) 基于iaa的协方差矩阵向量化的非均匀稀疏阵列测向方法
CN105467365A (zh) 一种改善mimo雷达doa估计性能的低旁瓣发射方向图设计方法
CN111257845B (zh) 一种基于近似消息传递的不在网格目标角度估计方法
CN104020440B (zh) 基于l型干涉式线性阵列的二维波达角估计方法
CN103260240A (zh) 大规模mimo中基于分布匹配的散射信源定位方法
CN110244273A (zh) 一种基于均匀分布式阵列的目标角度估计方法
CN113189592B (zh) 考虑幅相互耦误差的车载毫米波mimo雷达测角方法
CN108120953A (zh) 一种基于波达方向估计的无线电定位方法
CN110286350A (zh) 一种l型稀疏阵doa估计的精确配对方法及装置
CN112255629A (zh) 基于联合uca阵列的序贯esprit二维不相干分布源参数估计方法
CN113438682A (zh) 一种基于波束成形的sage-bem 5g无线信道参数提取方法
CN111308416B (zh) 基于四阶累积量的近场非圆信源参数估计方法
CN110231589B (zh) 一种大扩散角的多径信号波达方向估计方法
CN106125039B (zh) 基于局域联合处理的改进空时自适应单脉冲测角方法
CN113075649A (zh) 一种适用于分布式网络化雷达的信号级直接定位方法
CN117471396A (zh) 一种单基地对称嵌套阵mimo系统及近场定位方法
CN105158736B (zh) 一种mimo雷达发射方向图与波形设计方法
CN115236586B (zh) 一种基于数据预处理的极地冰下doa估计方法
CN113341371B (zh) 一种基于l阵和二维esprit算法的doa估计方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant