CN110174659A - 基于迭代近端投影的mimo雷达多测量矢量doa估计方法 - Google Patents
基于迭代近端投影的mimo雷达多测量矢量doa估计方法 Download PDFInfo
- Publication number
- CN110174659A CN110174659A CN201910541040.XA CN201910541040A CN110174659A CN 110174659 A CN110174659 A CN 110174659A CN 201910541040 A CN201910541040 A CN 201910541040A CN 110174659 A CN110174659 A CN 110174659A
- Authority
- CN
- China
- Prior art keywords
- matrix
- mimo radar
- sparse
- vector
- doa estimation
- 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
Links
- 239000013598 vector Substances 0.000 title claims abstract description 76
- 238000000034 method Methods 0.000 title claims abstract description 40
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 58
- 238000005457 optimization Methods 0.000 claims abstract description 38
- 230000009467 reduction Effects 0.000 claims abstract description 24
- 230000001427 coherent effect Effects 0.000 claims abstract description 13
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 11
- 238000013213 extrapolation Methods 0.000 claims abstract description 8
- 241000932075 Priacanthus hamrur Species 0.000 claims abstract 3
- 239000011159 matrix material Substances 0.000 claims description 112
- 238000005259 measurement Methods 0.000 claims description 24
- 230000009466 transformation Effects 0.000 claims description 9
- 238000001914 filtration Methods 0.000 claims description 5
- 238000001228 spectrum Methods 0.000 claims description 5
- 230000005540 biological transmission Effects 0.000 claims description 4
- 230000008569 process Effects 0.000 claims description 4
- 230000017105 transposition Effects 0.000 claims description 4
- 208000001992 Autosomal Dominant Optic Atrophy Diseases 0.000 claims 1
- 206010011906 Death Diseases 0.000 claims 1
- 238000012804 iterative process Methods 0.000 abstract 1
- 238000000342 Monte Carlo simulation Methods 0.000 description 5
- 238000012545 processing Methods 0.000 description 5
- 238000009499 grossing Methods 0.000 description 4
- NUHSROFQTUXZQQ-UHFFFAOYSA-N isopentenyl diphosphate Chemical compound CC(=C)CCO[P@](O)(=O)OP(O)(O)=O NUHSROFQTUXZQQ-UHFFFAOYSA-N 0.000 description 4
- 230000035945 sensitivity Effects 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 238000003491 array Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 206010057239 Post laminectomy syndrome Diseases 0.000 description 2
- 238000013459 approach Methods 0.000 description 2
- 230000008602 contraction Effects 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- OAICVXFJPJFONN-UHFFFAOYSA-N Phosphorus Chemical compound [P] OAICVXFJPJFONN-UHFFFAOYSA-N 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000013499 data model Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 230000004807 localization Effects 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000011084 recovery Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 238000002945 steepest descent method Methods 0.000 description 1
- 230000001131 transforming effect Effects 0.000 description 1
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
-
- 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
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
- G01S7/418—Theoretical aspects
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T10/00—Road transport of goods or passengers
- Y02T10/10—Internal combustion engine [ICE] based vehicles
- Y02T10/40—Engine management systems
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Computer Networks & Wireless Communication (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，通过高维回波数据转换至低维空间以降低空域维度，并对降维后的数据进行奇异值分解，提取信号子空间以降低时域维度，利用近端函数优化模型来表示MIMO雷达多测量矢量DOA估计中的非凸非平滑稀疏优化问题，然后在迭代过程中通过外推步骤和SCAD函数获得近端算子以求解该优化问题。本发明方法在低快拍和低信噪比下相干信源的DOA估计性能优于现有算法。
Description
技术领域
本发明涉及多输入多输出雷达目标参数估计方法，特别是涉及一种基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法。
背景技术
多输入多输出(Multiple Input and Multiple Output,MIMO)雷达是基于MIMO通信技术发展起来的一种新体制雷达。相比于传统相控阵雷达，MIMO雷达在目标探测、抗干扰、目标参数估计和目标识别等方面显示出诸多优势，因此已成为学术界研究的热点。
波达方向(Direction of Arrival,DOA)估计是MIMO雷达参数估计的一个重要内容。如今关于MIMO雷达DOA估计的算法已经多不胜数，诸如Capon算法，多重信号分类(Multiple Signal Classification,MUSIC)算法，基于旋转不变技术的信号参数估计(Estimation of Signal Parameters via Rotational Invariance Technique,ESPRIT)算法等。
通常，目标相对整个观测空间高度稀疏，因此许多学者将压缩感知理论应用到MIMO雷达DOA估计中。然而，求解l0范数最小化问题属于NP-hard问题，需要组合搜索、当维度增加时难以实现。Malioutov等人在论文“A sparse signal reconstructionperspective for source localization with sensor arrays”(IEEE Transactions onSignal Processing,2005,53(8):3010-3022)中针对传统阵列提出一种基于l1-SVD(l1norm-Singular Value Decomposition)的DOA估计方法，该方法通过对阵列接收数据进行奇异值分解，并建立信号子空间构建l2,1范数联合稀疏模型，然后采用二阶锥规划(Second Order Cone Programming,SOCP)求解该模型。Liu等人在论文“Direction ofarrival estimation via reweighted l1norm penalty algorithm for monostaticMIMO radar”(Multidimensional Systems and Signal Processing,2018,29(2):733-744)中提出一种基于加权l1范数算法的单基地MIMO雷达DOA估计方法，该算法利用降维Capon(Reduced-Dimensional Capon,RD-Capon)空间谱的系数构造加权矩阵，对小稀疏向量加以较大权重惩罚，对大稀疏向量加以较小权重惩罚以促进解的稀疏性，但由于采用SOCP求解方法导致其计算复杂度较高。Mohimani等人在论文“A Fast Approach forOvercomplete Sparse Decomposition Based on Smoothed L0Norm”(IEEE Transactionson Signal Processing,2009,57(1):289-301)中提出一种平滑l0范数(Smoothed l0norm，SL0)算法，该方法利用梯度投影原理以及最速下降法，采用一个平滑的高斯函数近似l0范数，从而将离散函数的优化问题转化为连续函数的优化问题。该算法的运算效率较高，在保证相同精度的条件下，能够比基追踪算法的重构速度快2～3倍。为了进一步提高目标的重构性能，Jing等人在论文“Joint Smoothed-Norm DOA Estimation Algorithm forMultiple Measurement Vectors in MIMO Radar”(Sensors,2017,17(5):1068)中提出一种基于平滑l0范数的MIMO雷达多测量矢量(Multiple Measurement Vector，MMV)DOA估计方法，利用接收数据的四阶累积量矩阵构建稀疏表示模型以有效抑制噪声，并将MMV问题转化为一个联合平滑函数的求解问题，且该方法相比l1-SVD方法快约两个数量级。Sadeghi等人在论文“Iterative Sparsification-Projection:Fast and Robust Sparse SignalApproximation”(IEEE Transactions on Signal Processing,2016,64(21):5536-5548)中基于近端方法提出一种相对SL0算法更为广泛的稀疏重构算法以解决非凸平滑问题，称为迭代稀疏投影(Iterative Sparsification-Projection,ISP)算法。针对非凸非平滑优化问题，Ghayem等人在论文“Sparse Signal Recovery Using Iterative ProximalProjection”(IEEE Transactions on Signal Processing,2018,66(4):879-894)中提出一种迭代近端投影(Iterative Proximal Projection,IPP)算法，该算法通过构造近端投影优化模型，并利用SCAD(Smoothly Clipped Absolute Deviation Penalty)阈值函数获得近端算子来求解该模型以促进解的稀疏度，同时引用外推步骤改善该算法的收敛性能。传统基于压缩感知的MIMO雷达DOA估计算法将非凸非平滑稀疏表示问题近似成凸或平滑函数问题进行求解，因此在稀疏表示模型上存在一定程度上的误差，从而导致传统基于压缩感知的DOA估计算法的性能不能达到最优。
发明内容
发明目的：本发明解决的技术问题是针对基于压缩感知的MIMO雷达DOA估计问题，提供一种能有效解决非凸非平滑稀疏表示问题的方法以提高DOA估计性能。
技术方案：本发明提供了一种基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，包括如下步骤：
步骤(1)：对MIMO雷达接收阵列信号进行匹配滤波，取多个快拍下虚拟阵列的输出信号矩阵；
步骤(2)：对MIMO雷达虚拟阵列输出信号矩阵进行降维变换，得到降维后的接收数据矩阵；
步骤(3)：对MIMO雷达降维后的输出信号进行奇异值分解，得到由多测量矢量构成的矩阵；
步骤(4)：根据稀疏重构理论，将搜索空域按等角度间隔划分，将步骤(3)得到的矩阵转换成稀疏表示模型；
步骤(5)：利用近端函数模型建立MIMO雷达多测量矢量DOA估计的稀疏优化问题，在迭代过程中通过外推步骤和SCAD函数获得近端算子以求解该优化问题；
步骤(6)：步骤(5)获得稀疏解之后，通过搜索其谱峰所在位置得到真实目标DOA估计值。
进一步地，步骤(1)具体为：对具有M个发射阵元和N个接收阵元的MIMO雷达接收阵列信号进行匹配滤波，取J个快拍下MIMO雷达虚拟阵列的输出信号X＝AS+N，其中为接收信号矩阵；为信号矩阵，其中表示复数域；为高斯噪声矩阵；为发射接收联合导向矩阵，其中为对应第p(p＝1,2,...,P)个目标的发射阵列的导向向量，为对应第p(p＝1,2,...,P)个目标的接收阵列的导向向量，(·)T表示矩阵转置，θp为第p(p＝1,2,...,P)个目标的方位角，表示Kronecker积，P是相干目标的数目。
步骤(2)具体为：对MIMO雷达虚拟阵列输出信号矩阵X进行降维变换，得到降维后的接收数据矩阵其中，为降维矩阵， 为转换矩阵， 0N×M为N×M维的零矩阵；(·)H表示共轭转置运算
步骤(3)具体为：为了进一步降低计算复杂度和对噪声的敏感性，对MIMO雷达降维后的输出信号矩阵进行奇异值分解，得到其中USV由P个大特征值对应的左奇异值矢量组成的信号子空间矩阵，UNV由其余M+N-1-P个小特征值对应的左奇异特征值矢量组成的噪声子空间矩阵，V为右奇异特征值矢量组成的矩阵，Λ为的特征值构成的对角矩阵。令则为由多测量矢量构成的矩阵，可表示为其中为降维后的阵列流形矩阵，为虚拟均匀线阵导向矩阵，SSV＝SVDP，ΛP×P由P个大特征值组成的对角矩阵，0P×(J-P)为P×(J-P)维的零矩阵。
步骤(4)具体为：根据稀疏重构理论，将搜索空域[-90°,90°]按等角度间隔划分为L个单元，且L＞＞P，表示空域内所有可能的入射方向，定义冗余字典其中则又可转换成稀疏表示模型：其中,与SSV具有相同的行支撑，即Sθ是P行稀疏矩阵，Sθ中的非零行元素对应冗余字典中目标的DOA，因此DOA估计问题可以转化为求解稀疏矩阵Sθ中非零行元素位置的问题。
步骤(5)具体为：利用近端函数优化模型建立MIMO雷达多测量矢量DOA估计的稀疏优化问题：其中，是非凸非平滑函数，z为辅助变量，是由矩阵z的每一行向量的l2范数构成的列向量，是由矩阵Sθ的每一行向量的l2范数构成的列向量，定义为可行集的指示函数。该稀疏优化问题需要多次迭代来求解，其中在第k次迭代中的稀疏解可表示为其中为非凸非平滑函数的近端算子，代表可行集的投影。
步骤(6)：在获得稀疏解之后，通过搜索其谱峰所在位置得到真实目标DOA估计值。
进一步地：步骤(5)中的通过多次迭代来求解稀疏优化问题的具体步骤为：
步骤(a)：为了避免迭代近端投影算法在求解上述非凸非平滑稀疏优化问题时易陷入局部最小值的问题，利用外推步骤改善迭代近端投影算法性能，则在第k次迭代中的稀疏解表示为其中w≥0为权重常数。该模型在第k次迭代中的稀疏解可进一步表示为其中为非凸非平滑函数的近端算子。
步骤(b)：由于非平滑函数的近端算子可通过SCAD(Smoothly Clipped Absolute Deviation Penalty)惩罚函数产生相应的SCAD阈值函数来计算，其中，为矢量中第l(l＝1,2,…,L)个元素，λ为调整参数，a为常量，一般取值为a＞2，sign(·)为符号函数，(α)+＝max(α，0)。
发明原理：由于传统基于压缩感知的MIMO雷达DOA估计算法通常将非凸非平滑稀疏优化问题近似成凸或平滑函数优化问题进行求解，由于稀疏优化问题模型误差的存在而导致DOA估计性能不理想。因此，本发明提出一种基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法。该方法首先将高维回波数据转换至低维空间以降低空域维度，并对降维后的数据进行奇异值分解，提取信号子空间以降低时域维度，利用近端函数优化模型来表示MIMO雷达多测量矢量DOA估计中的非凸非平滑稀疏优化问题，然后在迭代过程中通过外推步骤和SCAD函数获得近端算子以求解该模型。本发明方法在低快拍和低信噪比下相干信源的DOA估计性能优于现有算法。
有益效果：与现有技术相比
(1)基于压缩感知的MIMO雷达DOA估计优化问题是一种非凸非平滑稀疏优化问题，由于该问题求解难度较大，因此传统方法通常将此问题近似成凸或平滑函数优化问题进行求解，而本发明利用了近端函数优化模型较好地表示了该非凸非平滑稀疏优化问题，尽可能地降低稀疏优化问题模型的误差，从而提高了MIMO雷达DOA估计性能；
(2)本发明采用SCAD函数获得近端算子来求解以近端函数模型所表示的非凸非平滑稀疏优化问题，不仅克服了硬阈值收缩函数对数据中微小波动的敏感性，而且避免了软阈值收缩函数带来的偏差，因此能够进一步促进解的稀疏性，有效提高了本发明算法对MIMO雷达DOA估计的精度；
(3)本发明方法将高维回波数据转换至低维空间以降低空域维度，并对降维后的数据进行奇异值分解(SVD)，提取信号子空间以降低时域维度，这样不仅提高了算法的实时性，而且降低了算法对噪声的敏感性。
附图说明
图1是本发明的流程图；
图2是不同算法的DOA估计均方根误差随信噪比变化的曲线图；
图3是不同算法的DOA估计均方根误差随快拍数的变化关系；
图4是不同算法的DOA估计均方根误差随角度间隔变化关系。
具体实施方式
下面结合实施例对本发明作更进一步的说明。
如图1所示，本发明的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法包括如下步骤：
步骤1：对具有M个发射阵元和N个接收阵元的MIMO雷达接收阵列信号进行匹配滤波，取J个快拍下MIMO雷达虚拟阵列的输出信号X＝AS+N，其中为接收信号矩阵；为信号矩阵，其中表示复数域；为高斯噪声矩阵；为发射接收联合导向矩阵，其中为对应第p(p＝1,2,...,P)个目标的发射阵列的导向向量，为对应第p(p＝1,2,...,P)个目标的接收阵列的导向向量，(·)T表示矩阵转置，θp为第p(p＝1,2,...,P)个目标的方位角，表示Kronecker积，P是相干目标的数目。
假设窄带单基地MIMO雷达系统具有M个发射阵元和N个接收阵元，发射和接收阵列均为均匀线阵，其阵元间隔分别为dt＝λ/2和dr＝λ/2，λ为接收信号波长。假设存在P个远场窄带相干目标，其入射角度分别为θ1,θ2,...,θP，MIMO雷达接收阵列信号经匹配滤波后可表示为：
x(t)＝As(t)+n(t) (1)
式中，为发射接收联合导向矩阵；为对应第p(p＝1,2,...,P)个目标的发射阵列的导向向量；为对应第p(p＝1,2,...,P)个目标的接收阵列的导向向量；(·)T表示矩阵转置；表示复数域；表示Kronecker积；为相干信源信号矢量，其中，sp(t)＝αps1(t)，复常数αp表示sp(t)相对于s1(t)的相干系数；为接收阵列的噪声矩阵，服从零均值，方差为σ2的高斯分布，即n(t)～Nc(0,σ2IMN)，IMN表示MN×MN维的单位矩阵。取J个快拍下MIMO雷达虚拟阵列的输出信号，即
X＝AS+N (2)
式中，为接收信号矩阵；为信号矩阵；为高斯噪声矩阵。
步骤2：对MIMO雷达虚拟阵列输出信号矩阵X进行降维变换，得到降维后的接收数据矩阵其中，为降维矩阵， 为转换矩阵， 0N×M为N×M维的零矩阵；(·)H表示共轭转置运算。
对于发射阵列和接收阵列均为阵元间距等于半波长的均匀线阵的单基地MIMO雷达来说，其有效虚拟阵元个数为M+N-1，因此，MN×1维的目标均匀线阵导向矢量可由(M+N-1)×1维的虚拟均匀线阵导向矢量通过线性变换来表示，即
式中，表示(M+N-1)×1维的虚拟均匀线阵导向矢量，为转换矩阵，其中， 0N×M为N×M维的零矩阵。根据式(3)，阵列流形矩阵A可进一步表示为
A＝GB (4)
式中，B＝[b(θ1),b(θ2),...,b(θP)]为(M+N-1)×P维的虚拟均匀线阵导向矩阵。为了降低算法的计算复杂度，可定义一个大小为(M+N-1)×MN的降维矩阵对接收数据进行降维预处理。为了使得降维后的噪声是服从N～(0,σ2IM+N-1)的高斯白噪声，降维矩阵T需满足TTH＝IM+N-1，因此降维矩阵可选取为其中
式中，min(·)表示取最小的元素；diag(·)表示对角化操作；(·)H表示共轭转置运算。经过降维变换后，式(2)的接收数据矩阵可表示为：
式中，为降维后的接收数据矩阵；为降维后的阵列流形矩阵；为降维后的高斯白噪声矩阵。
步骤3：为了进一步降低计算复杂度和对噪声的敏感性，对MIMO雷达降维后的输出信号矩阵进行奇异值分解，得到其中USV由P个大特征值对应的左奇异值矢量组成的信号子空间矩阵，UNV由其余M+N-1-P个小特征值对应的左奇异特征值矢量组成的噪声子空间矩阵，V为右奇异特征值矢量组成的矩阵，Λ为的特征值构成的对角矩阵。令则为由多测量矢量构成的矩阵，可表示为其中为降维后的阵列流形矩阵，为虚拟均匀线阵导向矩阵，SSV＝SVDP， ΛP×P由P个大特征值组成的对角矩阵，0P×(J-P)为P×(J-P)维的零矩阵。
为了进一步降低计算复杂度和对噪声的敏感性，对接收数据矩阵进行奇异值分解(SVD)，
其中，USV由P个大特征值对应的左奇异值矢量组成的信号子空间矩阵，UNV由其余M+N-1-P个小特征值对应的左奇异特征值矢量组成的噪声子空间矩阵，V为右奇异特征值矢量组成的矩阵，Λ为X的特征值组成的对角矩阵。令进一步推导式(7)得：
其中，ΛP×P由P个大特征值组成的对角矩阵，0P×(J-P)为P×(J-P)维的零矩阵。令SSV＝SVDP，则降维后的接收数据模型(8)可简化为：
对比式(6)和式(9)可知，阵列接收数据矩阵的维数从(M+N-1)×J降到(M+N-1)×P，而实际应用中，目标个数P远远小于快拍数J，即P＜＜J，因此，式(9)的信号矩阵维度得到显著降低，有利于目标DOA估计的快速实现。
步骤4：根据稀疏重构理论，将搜索空域[-90°,90°]按等角度间隔划分为L个单元，且L＞＞P，表示空域内所有可能的入射方向，定义冗余字典其中则又可转换成稀疏表示模型：其中,与SSV具有相同的行支撑，即Sθ是P行稀疏矩阵，Sθ中的非零行元素对应冗余字典中目标的DOA，因此DOA估计问题可以转化为求解稀疏矩阵Sθ中非零行元素位置的问题。
根据稀疏重构理论，将搜索空域[-90°,90°]按等角度间隔划分为L个单元，且L＞＞P，表示空域内所有可能的入射方向，定义冗余字典其中则式(9)可转换成稀疏表示模型：
其中,与SSV具有相同的行支撑，即Sθ是P行稀疏矩阵，Sθ中的非零行元素对应冗余字典中目标的DOA，因此DOA估计问题可以转化为求解稀疏矩阵Sθ中非零行元素位置的问题。从多个测量矢量中恢复多个未知的稀疏信号源Sθ，该重建模型(10)可以称为多测量矢量(MMV)模型。该模型采用多快拍下的阵列接收数据，提高了在低信噪比下DOA估计性能。
步骤5：利用近端函数优化模型建立MIMO雷达多测量矢量DOA估计的稀疏优化问题：其中，是非凸非平滑函数，z为辅助变量，是由矩阵z的每一行向量的l2范数构成的列向量，是由矩阵Sθ的每一行向量的l2范数构成的列向量，定义为可行集的指示函数。该稀疏优化问题需要多次迭代来求解，其中在第k次迭代中的稀疏解可表示为其中为非凸非平滑函数的近端算子，代表可行集的投影。
为获得式(10)的稀疏解，可将MIMO雷达的多测量矢量DOA估计问题转化为近端函数优化问题，其稀疏优化问题为：
其中，是非凸非平滑函数，z为辅助变量，是由矩阵z的每一行向量的l2范数构成的列向量，是由矩阵Sθ的每一行向量的l2范数构成的列向量。定义为可行集的指示函数，即：
在式(11)的稀疏优化模型中引入惩罚函数方法，即：
式中，λ＞0为惩罚因子，||·||2表示l2范数。通过交替最小化方法将式(13)转化为关于和的两个子问题的迭代求解，即：
根据近端投影的定义，式(14)进一步化简为：
其中，为非平滑函数的近端算子。该模型第k次迭代中的稀疏解可表示为：
其中，代表可行集的投影。
为了避免在求解非凸非平滑问题时易陷入局部最小值的问题，利用外推步骤改善迭代近端投影(IPP)算法的收敛性能，则在第k次迭代中的稀疏解表示为：
其中，w≥0为权重常数。则式(13)在第k次迭代中的的稀疏解可进一步表示为：
其中，为非凸非平滑函数的近端算子。
非平滑函数的近端算子可通过不同的惩罚函数产生相应的阈值函数获得。采用l0惩罚函数来表示即则非平滑函数的近端算子可表示为
式中，硬阈值函数可表示为
式中，λ为调整参数。采用l1惩罚函数来表示即则非平滑函数的近端算子可表示为：
式中，软阈值函数可表示为：
当采用SCAD惩罚函数来表示即
式中，a为常量，一般取值为a＞2。则近端算子可采用SCAD阈值函数来表示，即
式中，SCAD阈值函数可表示为
式中，sign(·)为符号函数，(α)+＝max(α，0)。该SCAD阈值函数不仅克服了硬阈值收缩函数对数据中微小波动的敏感性，而且避免了软阈值收缩函数带来的偏差，因此能够进一步促进解的稀疏性。
在步骤5中采用外推步骤和SCAD阈值函数来求解稀疏优化问题的具体步骤如下：
初始化：
(a)设初始值
(b)选取一组合适的序列[λ1,λ2,...,λK]，且λk+1＝cλk，0＜c＜1，λ1＝max{|u0|}，a＝2000，w＝0.9
算法迭代：
Fork＝1,2,...,K,其中K为外循环次数
(a)令λ＝λk
(b)进行Q次迭代求解全局最小值，并将该最小值投影到可行集上
(1)令
(2)Forq＝1,2,...,Q,其中Q为内循环次数
a)计算
b)计算
c)将投影到可行集， 代表可行集的投影
(3)令
最终解
步骤6：在获得稀疏解之后，通过搜索其谱峰所在位置得到真实目标DOA估计值。
本发明的技术效果可通过以下仿真进一步说明，为了验证本发明方法在MIMO雷达DOA估计方面的优势，选取MUSIC算法、前后向空域平滑MUSIC(Forward and BackwardSpatial Smoothing_Music,FBSS_MUSIC)算法、SL0(Smoothedl0norm,SL0)算法与l1-SVD(l1-norm Singular Value Decomposition,l1-SVD)算法进行对比，将本发明称为IPP_SCAD_SVD算法。假设均匀线阵单基地MIMO雷达的发射和接收阵元数均为5，收发阵元间隔为dt＝dr＝λ/2，在空间角度范围[-90°,90°]以角度间隔0.05°等分。假设存在3个远场窄带相干目标，设置各目标的DOA分别为θ1＝-19.8°，θ2＝0°，θ3＝20.6°。DOA估计的均方根误差定义为其中，表示第p个目标在第mt次蒙特卡罗实验中的目标DOA估计值，MT为蒙特卡罗实验次数。信噪比定义为在本发明方法中，设置SCAD阈值函数中的参数a＝2000，权重常数w＝0.9，外循环次数K＝3，衰减因子c＝0.8。
仿真实验1：图1为各种算法的DOA估计均方根误差随信噪比变化的曲线图。设置信噪比在-10～15dB之间变化，快拍数J＝100，进行100次蒙特卡罗实验。从图1中可以看出，由于信源相干导致协方差矩阵的秩亏缺，信号特征向量发散到噪声子空间，造成MUSIC算法的DOA估计方法失效；FBSS_MUSIC算法利用前后向空间平滑技术实现信号的解相干，能够有效估计DOA，但在高信噪比下该方法的DOA估计性能较差。SL0_SVD算法、L1_SVD算法以及IPP_SCAD_SVD算法均属于基于压缩感知的DOA估计方法，这些方法能够对相干信源的DOA进行有效估计，且DOA估计的均方根误差随信噪比的增加而减小，而本发明算法的DOA估计精度明显优于其他算法。
仿真实验2：图2为各种算法的DOA估计均方根误差随快拍数的变化关系曲线。设置信噪比为-5dB，进行100次蒙特卡罗实验，快拍数J在50～350之间变化。从图2中可以看出，MUSIC算法无法处理相干信源，其他各种算法的角度估计精度均随着快拍数的增加而提高，而本发明算法相比其他算法具有更高的DOA估计精度。
仿真实验3：图3为各种算法的DOA估计均方根误差随两个相干目标的入射角度间隔Δθ的变化关系。两个相干目标的入射角度分别为θ1＝6°，θ2＝6°+Δθ，其中Δθ∈[14°22°]，信噪比为5dB，快拍数J＝100，进行100次蒙特卡罗实验。从图3中可以看出，MUSIC算法无法分辨相干目标，FBSS_MUSIC算法、SL0_SVD算法、L1_SVD算法和本发明算法的DOA估计精度均随着目标角度间隔的增大而提高，而本发明算法的DOA估计精度始终高于其他算法，表明本发明算法相比其他算法能获得更高的空间角度分辨率。
Claims (7)
1.一种基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，包括如下步骤：
步骤(1)：对MIMO雷达接收阵列信号进行匹配滤波，取多个快拍下虚拟阵列的输出信号矩阵；
步骤(2)：对MIMO雷达虚拟阵列输出信号矩阵进行降维变换，得到降维后的接收数据矩阵；
步骤(3)：对MIMO雷达降维后的输出信号进行奇异值分解，得到由多测量矢量构成的矩阵；
步骤(4)：根据稀疏重构理论，将搜索空域按等角度间隔划分，将步骤(3)得到的矩阵转换成稀疏表示模型；
步骤(5)：利用近端函数模型建立MIMO雷达多测量矢量DOA估计的稀疏优化问题，在迭代过程中通过外推步骤和SCAD函数获得近端算子以求解该优化问题；
步骤(6)：步骤(5)获得稀疏解之后，通过搜索其谱峰所在位置得到真实目标DOA估计值。
2.根据权利要求1所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于：步骤(1)具体为：
对具有M个发射阵元和N个接收阵元的MIMO雷达接收阵列信号进行匹配滤波，取J个快拍下虚拟阵列的输出信号矩阵X＝AS+N；其中为接收信号矩阵；为信号矩阵，其中表示复数域；为高斯噪声矩阵；为发射接收联合导向矩阵，其中为对应第p(p＝1,2,...,P)个目标的发射阵列的导向向量，为对应第p(p＝1,2,...,P)个目标的接收阵列的导向向量，(·)T表示矩阵转置，θp为第p(p＝1,2,...,P)个目标的方位角，表示Kronecker积，P是相干目标的数目。
3.根据权利要求1所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，步骤(2)具体为：
对MIMO雷达虚拟阵列输出信号矩阵X进行降维变换，得到降维后的接收数据矩阵其中，为降维矩阵， 为转换矩阵，m＝0,1,...,M-1；0N×M为N×M维的零矩阵；(·)H表示共轭转置运算。
4.根据权利要求1所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，步骤(3)具体为：
对MIMO雷达降维后的输出信号进行奇异值分解，得到 其中USV由P个大特征值对应的左奇异值矢量组成的信号子空间矩阵，UNV由其余M+N-1-P个小特征值对应的左奇异特征值矢量组成的噪声子空间矩阵，V为右奇异特征值矢量组成的矩阵，Λ为的特征值构成的对角矩阵；令则为由多测量矢量构成的矩阵，可表示为其中为降维后的阵列流形矩阵，为虚拟均匀线阵导向矩阵，SSV＝SVDP，ΛP×P由P个大特征值组成的对角矩阵，0P×(J-P)为P×(J-P)维的零矩阵。
5.根据权利要求1所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，步骤(4)具体为：
根据稀疏重构理论，将搜索空域[-90°,90°]按等角度间隔划分为L个单元，且L＞＞P，表示空域内所有可能的入射方向，定义冗余字典其中则又可转换成稀疏表示模型：其中,与SSV具有相同的行支撑，即Sθ是P行稀疏矩阵，Sθ中的非零行元素对应冗余字典中目标的DOA。
6.根据权利要求1所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，步骤(5)具体为：
利用近端函数优化模型建立MIMO雷达多测量矢量DOA估计的稀疏优化问题：其中，是非凸非平滑函数，z为辅助变量，是由矩阵z的每一行向量的l2范数构成的列向量，是由矩阵Sθ的每一行向量的l2范数构成的列向量，定义为可行集的指示函数；该稀疏优化问题需要多次迭代求解，其中在第k次迭代中的稀疏解可表示为其中为非凸非平滑函数的近端算子，代表可行集的投影。
7.根据权利要求6所述的基于迭代近端投影的MIMO雷达多测量矢量DOA估计方法，其特征在于，步骤(5)中的通过多次迭代求解稀疏优化问题的具体步骤为：
步骤(a)：利用外推步骤改善算法性能，在第k次迭代中的稀疏解表示为其中w≥0为权重常数；该模型在第k次迭代中的稀疏解可进一步表示为其中为非凸非平滑函数的近端算子；
步骤(b)：非平滑函数的近端算子可通过SCAD惩罚函数产生相应的SCAD阈值函数来计算，其中，为矢量中第l(l＝1,2,…,L)个元素，λ为调整参数，a为常量，一般取值为a＞2，sign(·)为符号函数，(α)+＝max(α，0)。
Applications Claiming Priority (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910382095 | 2019-05-08 | ||
| CN2019103820950 | 2019-05-08 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110174659A true CN110174659A (zh) | 2019-08-27 |
| CN110174659B CN110174659B (zh) | 2023-05-23 |
Family
ID=67698839
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910541040.XA Active CN110174659B (zh) | 2019-05-08 | 2019-06-21 | 基于迭代近端投影的mimo雷达多测量矢量doa估计方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110174659B (zh) |
Cited By (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111060885A (zh) * | 2019-12-13 | 2020-04-24 | 河海大学 | 一种mimo雷达的参数估计方法 |
| CN111107023A (zh) * | 2019-09-23 | 2020-05-05 | 南京邮电大学 | 在大规模mimo中基于光滑范数的压缩感知信道估计方法 |
| CN111142063A (zh) * | 2020-01-06 | 2020-05-12 | 西安邮电大学 | 一种基于降维优化的快速压缩感知低空目标测角方法 |
| CN111680262A (zh) * | 2020-05-07 | 2020-09-18 | 北京控制工程研究所 | 一种自适应增益梯度投影辨识方法 |
| CN112462352A (zh) * | 2020-10-30 | 2021-03-09 | 哈尔滨工程大学 | 一种适用于低信噪比条件下的线谱增强方法 |
| CN113835085A (zh) * | 2021-09-30 | 2021-12-24 | 南京信息工程大学 | 一种基于复杂地形补偿的雷达快速测高方法 |
| CN116362316A (zh) * | 2023-05-29 | 2023-06-30 | 成都阿加犀智能科技有限公司 | 一种模型转换方法、装置、存储介质及电子设备 |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080231505A1 (en) * | 2007-03-23 | 2008-09-25 | Weiqing Zhu | Method of Source Number Estimation and Its Application in Method of Direction of Arrival Estimation |
| US20120259590A1 (en) * | 2011-04-11 | 2012-10-11 | Jong Chul Ye | Method and apparatus for compressed sensing with joint sparsity |
| CN103954950A (zh) * | 2014-04-25 | 2014-07-30 | 西安电子科技大学 | 一种基于样本协方差矩阵稀疏性的波达方向估计方法 |
| CN103983958A (zh) * | 2014-05-16 | 2014-08-13 | 哈尔滨工程大学 | 基于多测量矢量稀疏表示的mimo雷达连续目标角度估计方法 |
| JP2017036990A (ja) * | 2015-08-10 | 2017-02-16 | 三菱電機株式会社 | 到来方向推定装置 |
| CN106772225A (zh) * | 2017-01-20 | 2017-05-31 | 大连大学 | 基于压缩感知的波束域doa估计 |
| CN108957388A (zh) * | 2018-05-21 | 2018-12-07 | 南京信息工程大学 | 一种基于协方差匹配sl0算法的mimo雷达相干信源doa估计方法 |
-
2019
- 2019-06-21 CN CN201910541040.XA patent/CN110174659B/zh active Active
Patent Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20080231505A1 (en) * | 2007-03-23 | 2008-09-25 | Weiqing Zhu | Method of Source Number Estimation and Its Application in Method of Direction of Arrival Estimation |
| US20120259590A1 (en) * | 2011-04-11 | 2012-10-11 | Jong Chul Ye | Method and apparatus for compressed sensing with joint sparsity |
| CN103954950A (zh) * | 2014-04-25 | 2014-07-30 | 西安电子科技大学 | 一种基于样本协方差矩阵稀疏性的波达方向估计方法 |
| CN103983958A (zh) * | 2014-05-16 | 2014-08-13 | 哈尔滨工程大学 | 基于多测量矢量稀疏表示的mimo雷达连续目标角度估计方法 |
| JP2017036990A (ja) * | 2015-08-10 | 2017-02-16 | 三菱電機株式会社 | 到来方向推定装置 |
| CN106772225A (zh) * | 2017-01-20 | 2017-05-31 | 大连大学 | 基于压缩感知的波束域doa估计 |
| CN108957388A (zh) * | 2018-05-21 | 2018-12-07 | 南京信息工程大学 | 一种基于协方差匹配sl0算法的mimo雷达相干信源doa估计方法 |
Non-Patent Citations (3)
| Title |
|---|
| FATEME GHAYEM 等: "Sparse Signal Recovery Using Iterative Proximal Projection", 《IEEE TRANSACTIONS ON SIGNAL PROCESSING 》 * |
| 谭伟杰 等: "基于稀疏表示的平行互素阵二维测向方法", 《系统工程与电子技术》 * |
| 陈金立 等: "基于协方差匹配SL0算法的MIMO雷达DOA估计", 《雷达科学与技术》 * |
Cited By (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111107023A (zh) * | 2019-09-23 | 2020-05-05 | 南京邮电大学 | 在大规模mimo中基于光滑范数的压缩感知信道估计方法 |
| CN111107023B (zh) * | 2019-09-23 | 2024-02-02 | 南京邮电大学 | 在大规模mimo中基于光滑范数的压缩感知信道估计方法 |
| CN111060885A (zh) * | 2019-12-13 | 2020-04-24 | 河海大学 | 一种mimo雷达的参数估计方法 |
| CN111142063A (zh) * | 2020-01-06 | 2020-05-12 | 西安邮电大学 | 一种基于降维优化的快速压缩感知低空目标测角方法 |
| CN111680262A (zh) * | 2020-05-07 | 2020-09-18 | 北京控制工程研究所 | 一种自适应增益梯度投影辨识方法 |
| CN111680262B (zh) * | 2020-05-07 | 2023-07-14 | 北京控制工程研究所 | 一种自适应增益梯度投影辨识方法 |
| CN112462352A (zh) * | 2020-10-30 | 2021-03-09 | 哈尔滨工程大学 | 一种适用于低信噪比条件下的线谱增强方法 |
| CN113835085A (zh) * | 2021-09-30 | 2021-12-24 | 南京信息工程大学 | 一种基于复杂地形补偿的雷达快速测高方法 |
| CN113835085B (zh) * | 2021-09-30 | 2023-07-25 | 南京信息工程大学 | 一种基于复杂地形补偿的雷达快速测高方法 |
| CN116362316A (zh) * | 2023-05-29 | 2023-06-30 | 成都阿加犀智能科技有限公司 | 一种模型转换方法、装置、存储介质及电子设备 |
| CN116362316B (zh) * | 2023-05-29 | 2023-12-12 | 成都阿加犀智能科技有限公司 | 一种模型转换方法、装置、存储介质及电子设备 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110174659B (zh) | 2023-05-23 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110261841B (zh) | 基于迭代加权近端投影的mimo雷达单测量矢量doa估计方法 | |
| CN110174659A (zh) | 基于迭代近端投影的mimo雷达多测量矢量doa估计方法 | |
| CN104749553B (zh) | 基于快速稀疏贝叶斯学习的波达方向角估计方法 | |
| Wang et al. | A sparse representation scheme for angle estimation in monostatic MIMO radar | |
| CN104020439B (zh) | 基于空间平滑协方差矩阵稀疏表示的波达方向角估计方法 | |
| CN104251989B (zh) | 基于压缩空间谱的单基地mimo雷达目标波达方向估计方法 | |
| CN109375154B (zh) | 一种冲击噪声环境下基于均匀圆阵的相干信号参数估计方法 | |
| CN107870315B (zh) | 一种利用迭代相位补偿技术估计任意阵列波达方向方法 | |
| CN108303683B (zh) | 单基地mimo雷达实值esprit非圆信号角度估计方法 | |
| CN107991659B (zh) | 基于字典学习的米波雷达低仰角目标测高方法 | |
| CN112782695B (zh) | 基于isar图像和参数优化的卫星姿态和尺寸估计方法 | |
| CN112379327A (zh) | 一种基于秩损估计的二维doa估计与互耦校正方法 | |
| CN107121665B (zh) | 一种基于稀疏阵的近场相干源的无源定位方法 | |
| CN113567913A (zh) | 基于迭代重加权可降维的二维平面doa估计方法 | |
| CN107450046B (zh) | 低仰角多径环境下的波达角估计方法 | |
| CN105242236B (zh) | 宽带信号超分辨测向中的阵元位置误差校正方法 | |
| CN106483193A (zh) | 一种基于高阶累计量的波达快速估计方法 | |
| CN109507634B (zh) | 一种任意传感器阵列下的基于传播算子的盲远场信号波达方向估计方法 | |
| CN110133641A (zh) | 一种尺度自适应的穿墙成像雷达目标跟踪方法 | |
| CN115329261A (zh) | 一种基于空间平滑稀疏重构的mimo雷达低仰角估计方法 | |
| CN112858995B (zh) | 一种基于分布式阵列的联合测角方法和装置 | |
| CN105303009A (zh) | 基于压缩感知与正则mfocuss的超分辨谱估计方法 | |
| CN109633635B (zh) | 基于结构化递归最小二乘的米波雷达测高方法 | |
| CN103792509B (zh) | 电磁信号的二维波达方向角估计方法 | |
| CN110412535A (zh) | 一种序贯的空时自适应处理参数估计方法 |
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 |