CN114844544B - 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质 - Google Patents

一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质 Download PDF

Info

Publication number
CN114844544B
CN114844544B CN202210461073.5A CN202210461073A CN114844544B CN 114844544 B CN114844544 B CN 114844544B CN 202210461073 A CN202210461073 A CN 202210461073A CN 114844544 B CN114844544 B CN 114844544B
Authority
CN
China
Prior art keywords
tensor
matrix
signal
array
representing
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.)
Active
Application number
CN202210461073.5A
Other languages
English (en)
Other versions
CN114844544A (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 CN202210461073.5A priority Critical patent/CN114844544B/zh
Publication of CN114844544A publication Critical patent/CN114844544A/zh
Application granted granted Critical
Publication of CN114844544B publication Critical patent/CN114844544B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04B—TRANSMISSION
    • H04B7/00—Radio transmission systems, i.e. using radiation field
    • H04B7/02—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas
    • H04B7/04—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas
    • H04B7/06—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station
    • H04B7/0613—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission
    • H04B7/0615—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission of weighted versions of same signal
    • H04B7/0617—Diversity systems; Multi-antenna system, i.e. transmission or reception using multiple antennas using two or more spaced independent antennas at the transmitting station using simultaneous transmission of weighted versions of same signal for beam forming
    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04B—TRANSMISSION
    • H04B17/00—Monitoring; Testing
    • H04B17/30—Monitoring; Testing of propagation channels
    • H04B17/391—Modelling the propagation channel
    • 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
    • Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
    • Y02D30/00—Reducing energy consumption in communication networks
    • Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Signal Processing (AREA)
  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Variable-Direction Aerials And Aerial Arrays (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本发明公开了一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质，该方法包括步骤1：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并重排为多快拍虚拟阵元信号矩阵；步骤2：基于低管秩特性，对初始化张量进行补全；步骤3：对补全后的张量进行分解，得到三个因子矩阵；步骤4：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标角度和功率进行求解和匹配；步骤5：利用步骤4的参数估计值，计算波束成形权矢量。本发明所述技术方案利用张量的多维结构可提高参数估计的噪声鲁棒性；利用交替方向乘子法和张量分解的快速收敛性可显著减少计算量；目标方向通过参数反演得到，受先验误差影响较小。

Description

一种基于低管秩张量分解的互质阵列波束成形方法、系统及 介质
技术领域
本发明涉及阵列信号处理领域，尤其涉及基于低管秩张量分解的互质阵列波束成形方法、系统及介质。
背景技术
自适应波束成形是阵列信号处理领域的核心技术之一，由于其能够在保持目标信号的同时，抑制来自不同方向的干扰信号和背景噪声，该技术被广泛运用于雷达、声纳及无线通信等领域。自适应波束成形的性能通常随着阵列的自由度和孔径增大而改善。对于均匀线阵，其自由度和孔径与其阵元个数成线性关系，因此其阵列自由度与孔径的增加往往伴随着系统复杂度的提高和硬件成本的增加。互质阵列是按照互质准则设计的非均匀稀疏阵列，它通过增大阵元间距获得更大的孔径，又利用互质特性获得远大于阵元个数的自由度。与具有相同阵元数量的均匀线阵相比，互质阵列不仅具有更大的阵列孔径和自由度，还因为其更大的阵列间隔而具有更小的阵元互耦，因此互质阵列是一种更为有效的自适应波束成形结构。
目前，可以用于互质阵列自适应波束成形的算法主要分为两大类：通用类和特定类。通用类算法没有限定阵列的结构，因而可以直接用于互质阵列。然而，这类算法并没有利用互质阵列具有的高自由度特性，因而其性能受到物理阵元数量的限制。此外，部分通用类算法考虑的假设条件在互质阵列中不成立，也会导致该类算法运用在互质阵列中的性能严重下降。
特定类算法是指专门针对互质阵列所设计的自适应波束成形算法，其又可进一步细分为非填充类算法和填充类算法。互质阵列的高自由度来自于虚拟阵列，然而该虚拟阵列通常是非连续的，即有“孔洞”存在。非填充类算法和填充类算法的本质区别就是对孔洞的处理方式不同。
非填充类算法以利用连续虚拟阵元数据的方式来避免孔洞。这类算法的计算量较小，但是忽略了大量非连续虚拟阵元的携带的信息，因此算法的精度和最大可辨识信源数下降，尤其是在阵元个数较多时。
填充类算法通过对虚拟阵列中的孔洞进行填充，实现全虚拟阵元的利用。文件1“Robust adaptive beamforming via coprime coarray interpolation[J]”(SignalProcess,2020,169(107382))利用最小化核范数对多采样虚拟信号矩阵进行填充，然后利用多重信号分类(Multiple Signal Classification，MUSIC)方法估计目标角度，接着利用最小二乘(Least Square，LS)求解信号功率，进而重构干扰加噪声协方差，最后计算波束成形权矢量。然而，该方法利用内点法求解核范数优化问题，具有较高的运算量；MUSIC的使用致使该算法的噪声鲁棒性较弱。文件2“Coprime array-adaptive beamforming viaatomic-norm-based sparse recovery[J]”(IET Radar,Sonar&Navigation,2021,15:第1494页～第1507页)将多采样虚拟信号通过原子范数建模，以交替迭代的方式实现对目标角度和功率的估计。然而，该算法利用一阶近似迭代修正目标角度，具有较大的运算量和较小的角度误差容限。
发明内容
本发明提供一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质，其本质上是一种插值类算法，能够有效提升波束成形器的噪声鲁棒性和降低计算量
本发明提供的技术方案如下：
一方面，一种基于低管秩张量分解的互质阵列波束成形方法，包括：
步骤1：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并将初始化内插后的虚拟阵元信号重排为多快拍虚拟阵元信号矩阵；
步骤2：将多快拍虚拟阵元信号矩阵重排为初始化张量，并基于低管秩特性，对初始化张量进行补全；
对初始化张量通过张量核范数最小化得到补全后的张量数据；
步骤3：对补全后的张量进行分解，得到三个因子矩阵；
步骤4：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标来波信号角度和功率进行求解和匹配，分别获得目标信号和干扰信号的角度和功率的估计值；
步骤5：利用步骤4的参数估计值，重构波束成形器，计算波束成形权矢量，进行波束成形。
进一步的，所述互质阵列下的远场窄带非相干信号模型如下：
其中，x(l)是指互质阵列在第l个采样时刻的信号，a(θ)是互质阵列在θ方向的导向矢量，λ表示接收信号的波长，d表示单位阵元间距，且d≤λ/2，qid表示互质阵列中从左至右第i个阵元的位置，互质阵列由两个均匀子阵组成，第一个均匀子阵由2M个阵元组成，阵元间距为Nd；第二个均匀子阵由N个阵元组成，阵元间距为Md，M和N是一对互质的正整数；s0(l)表示第l个采样时刻的目标源信号，sk(l)表示l个采样时刻的第k个干扰的源信号，k＝1,2,…,K，K表示干扰个数；n(l)表示独立同分布的零均值复高斯噪声矢量，即/> 表示复高斯分布，pn表示噪声功率，I表示单位矩阵。
进一步的，基于低管秩特性，对初始化的张量进行补全的过程如下：
其中，表示补全后的张量，/>表示/>的张量核范数；ε表示填补缺失数据的张量，Ω是缺失数据的位置索引，Ωc为在/>张量结构下Ω的补集，/>表示ε向Ωc的投影；为初始化的张量。
进一步的，补全后的张量采用ADMM算法进行迭代求解，具体过程如下：
设定增广拉格朗日函数，并进行参数初始化设置：μ0＝10-4：
其中，表示拉格朗日乘子，μ表示惩罚系数，μ>0，以使得设定的增广拉格朗日函数最小化时，依次对/>ε,/>μ求解，求解过程中，令非求解变量视为常量：
步骤A1：令并求解得到
其中，表示/>的第i个正向切片；/> 表示从第三个维度对/>进行傅里叶变换；/>表示/>的第i个正向切片；表示从第三个维度对Wt进行傅里叶变换；i＝1,2,3,4； 和/>来自于/>的截断奇异值分解，/>表示/>的共轭转置矩阵；
由可以得到/>进而/>表示从第三个维度对/>进行傅里叶逆变换：
步骤A2：
步骤A3：
步骤A4：μt+1＝min(ρμt,μmax)，μmax设置为105，ρ表示每次更新的放大倍数，取值为1.1，迭代上述步骤直至收敛或者迭代次数小于500，获得/>的解算值。
进一步的，各因子矩阵和的正向切片关系如下：
其中，X(1)、X(2)、X(3)、X(4)分别表示的四个正向切片， 表示/>的共轭转置矩阵，θ0表示目标来波方向，θ1,…,θK表示K个干扰来波方向；/> 表示阵元位置为{qd|0≤q≤V-1,q是整数}的均匀线阵在θ方向上的导向矢量；P表示携带功率信息，Φ表示携带角度信息， 分别表示对补全后的张量进行分解得到的三个因子矩阵；/>表示/>的第i行，i＝1,2,3,4；/>表示将向量转化为对角矩阵的操作。
对补全后的张量进行分解，得到三个因子矩阵，其过程下：
S3.1：首先将张量通过Tucker分解压缩为/>其中F1＝K+1+P，F2＝min{F1,4}，P<V-K-2为非负整数，通常取1或2。
S3.2：对进行平行因子分解，简称CP分解。首先，随机生成20组CP分解的初始因子矩阵/>分别按以下更新规则迭代5次：
其中，和⊙分别代表Hadamard积和Khatri-Rao积。计算迭代5次以后的误差值
令将第Q组初始化因子矩阵迭代五次以后的结果记为{Ag,bg，Cg}，以此作为新的初始化因子矩阵，继续按照上述过程迭代更新，直到收敛。
S3.3：令A＝UAg，B＝VBg，C＝WCg作为对进行CP分解的初始值，对/>只能找上述更新步骤进行CP分解，直至/>即可得到/>的三个因子矩阵/> 和/>
进一步的，对目标来波信号角度和功率进行求解和匹配，获得干扰来波信号角度和功率的估计值的具体过程如下：
步骤B1：解缩放模糊；
基于中每个元素的模为1，令 除去缩放模糊，使得
其中，∏是同用表示置换模糊的未知列置换矩阵，是/>的逆算子，表示将矩阵的对角元素排列成列向量的操作；
则置换模糊的功率和角度的估计值为
其中，ik∈{0,1,…,K}且各不相同，k＝1,2,…,K+1；
步骤B2：解置换模糊，即目标信源匹配；
已知目标角度先验区域Ψ＝[θp-Δθ,θp+Δθ]，θp表示先验角度，Δθ表示角度置信半径；令以及/>分别表示信号源功率的估计值的集合以及对应的角度的估计值集合；
当时：
令则/>将/>和/>的其余元素作为干扰信号角度和功率的估计/>和/>重构干扰加噪声协方差矩阵，即 表示干扰加噪声协方差矩阵，I表示单位矩阵，/>表示噪声功率的估计值，/>是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵；
当时：
选取中较大的K个元素及其在/>中的对应元素作为干扰来波信号角度和功率的估计值，重构/>再计算目标信号协方差矩阵：/> 表示目标信号协方差矩阵，则目标来波方向上的导向矢量可估计为/>其中，/>x(l)表示阵列接收信号在第l时刻的采样值，xH(l)是x(l)的共轭转置矩阵，L表示总的采样时刻数，/>表示/>的主特征值。
进一步的，波束成形权矢量的计算公式如下：
其中， 和/>是对干扰信号功率和角度的估计值，/>是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵，/>是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵，/>表示噪声功率的估计值，I表示单位矩阵。
另一方面，一种基于低管秩张量分解的互质阵列波束成形系统，包括：
虚拟阵元信号初始化与重排单元：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并将初始化内插后的虚拟阵元信号重排为多快拍虚拟阵元信号矩阵；
张量初始化与补全单元：将多快拍虚拟阵元信号矩阵重排为初始化张量，并基于低管秩特性，对初始化张量进行补全；
张量分解单元：对补全后的张量进行分解，得到三个因子矩阵；
目标角度与功率估计单元：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标来波信号角度和功率进行求解和匹配，获得干扰来波信号角度和功率的估计值；
波束成形单元：利用目标角度与功率估计单元获得的参数估计值，重构波束成形器，计算波束成形权矢量。
进一步的，补全后的张量采用ADMM算法求解单元进行迭代求解。
再一方面，一种计算机存储介质，包括计算机程序，所述计算机程序指令被处理终端执行时使所述处理终端执行上述一种基于低管秩张量分解的互质阵列波束成形方法。
有益效果
相较于现有技术而言，本方案存在以下几点优点：
1.利用张量的多维结构可提高参数估计的噪声鲁棒性；
2.利用交替方向乘子法(Alternating Direction Method ofMultipliers，ADMM)和张量分解的快速收敛性可显著减少计算量；
3.在本算法中，目标方向是通过参数反演得到，不需要利用泰勒一阶近似进行修正，受先验误差影响较小。
附图说明
图1本发明提供的基于低管秩张量分解的互质阵列自适应波束成形方法流程示意图；
图2本发明涉及的互质阵列结构示意图；
图3张量分解过程示意图；
图4干信比为-5dB，2干扰源条件下的波束图，其中，(a)为无先验误差，(b)为存在2°先验误差；
图5干信比为-5dB，4干扰源条件下的波束图，其中，(a)为无先验误差，(b)为存在2°先验误差；
图6干信比为-40dB，2干扰源条件下的波束图，其中，(a)为无先验误差，(b)为存在2°先验误差；
图7干信比为-40dB，4干扰源条件下的波束图，其中，(a)为无先验误差，(b)为存在2°先验误差；
图8使用不同方法得到的输出SINR与输入的变化关系，其中，(a)为无先验误差，(b)为存在2°先验误差，(c)为存在-4°到4°随机误差；
图9使用不同方法得到的输出SINR与采样快拍数的变化关系，其中，(a)为无先验误差，(b)为存在2°先验误差，(c)为存在-4°到4°随机误差；
图10不同插值类方法的计算时间比较，其中，(a)为存在2个干扰，(b)为存在4个干扰。
具体实施方式
下面将结合附图和实施例对本发明做进一步的说明。
为了便于说明，先给出一些张量的有关定义：
对于张量ZH，ZL，ZF分别表示张量/>的水平展开，侧向展开和正向展开，具体地
如图1所示，基于低管秩张量分解的互质阵列自适应波束成形算法，具体包括以下步骤：
S1：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，通过重排得到初始化多采样虚拟信号矩阵。
互质阵列的结构如图2所示，包括两个稀疏均匀线阵：第一个由2M个阵元组成，阵元间距为Nd；第二个由M个阵元组成，阵元间距为Md。两个子阵共享第一个阵元，因此总阵元数为2M+N-1。综上，阵元位置可以表示为
其中，M和N是一对互质数，qid表示从左至右第i个阵元的位置。
假设K+1个远场窄带非相干信号从不同角度θ0,θ1,…,θK到达阵列。其中，θ0是目标信号方向，θk,k＝1,2,…,K分别是K个不同干扰信号的方向。由此，阵列在第l个采样时刻的信号模型为
其中，是阵列在θ方向的导向矢量，λ表示接受信号的波长，通常d≤λ。s0(l)表示目标源信号，和sk(l)，k＝1,2,…,K，分别表示第k个干扰的源信号。/>表示独立同分布的复高斯噪声矢量。
由接收信号计算采样协方差，作为协方差矩阵的近似：
其中，L表示总的信号采样数，[·]H表示转置共轭操作，p0和pk分别表示目标信号和第k个干扰信号的输入功率。pn表示噪声功率，I2M+N-1表示维度为2M+N-1的单位阵。
向量化得到一个虚拟信号长矢量，即
其中，[·]*表示取共轭操作，/>表示Kronecker积，1＝vec(I2M+N-1)。相应的虚拟阵元的分布为
将r中对应于同一虚拟阵元的元素求平均，以此替代这些冗余的元素，再按照对应的虚拟阵元位置由小到大进行排序，得到新的虚拟阵元分布表示，相应的虚拟信号为
其中，是虚拟阵元在θ角度下的导向矢量，/>表示中心元素为1，其余元素为0的矢量。
中拟阵元是不连续的，存在“孔洞”。若将这些空洞由虚拟阵元填满，可以得到一个连续的虚拟阵元，定义为
那么，相应的初始化虚拟信号可定义为：
其中，<·>q表示向量的第q个元素。
由于是单快拍信号且具有Hermitian对称性(即/>)，可以将其重排为一个Hermitian Toeplitz矩阵，称为多快拍虚拟信号矩阵：
其中， 表示/>的势，即/>中的元素个数。此时，对特征展开可以得到/>其中/>表示/>中最大的K+1的特征值构成的对角矩阵。相应地，/>是/>的其余较小特征值构成的对角矩阵。取这些较小特征值的平均作为噪声信号功率的估计值，即Tr(·)表示矩阵的迹。由此，去噪后的多快拍虚拟信号矩阵可以表示为/>
S2：将多快拍虚拟信号矩阵重排为初始化张量，利用其低管秩特性，通过张量核范数最小化得到补全后的张量数据。
初始化张量使得
其中，Y(i)称为张量的第i个正向切片，i＝1,2,3,4。显然，Y中的部分元素为0，意味着/>的相对应位置的元素也为0，代表着互相关信息的缺失，因此需要对/>进行补全。
利用张量数据的低管秩特性补全初始化张量该过程可以表示为
其中，表示补全后的张量，/>表示其张量核范数。ε表示填补缺失数据的张量，Ω是缺失数据的位置索引，Ωc为其在上述张量结构下Ω的补集，/>表示ε向Ωc的投影。
采用ADMM算法对上式进行求解，首先给出增广拉格朗日函数(忽略常数项)
其中，表示拉格朗日乘子，μ>0是惩罚系数。
按照如下步骤进行迭代：
1.
其等价为
其中，表示/>的第i个正向切片；/>表示从第三个维度对/>进行傅里叶变换，/>同理。由上式可以解得
其中，[(Z)+]ij＝max([Z]ij,0)，和/>来自于/>的截断奇异值分解，即 由/>可以得到/>进而/>表示从第三个维度对/>进行傅里叶逆变换。
2.
3.
4.μt+1＝min(ρμt,μmax)，μmax设置为105。
迭代上述步骤直至收敛或者迭代次数小于500。
S3：对补全后的张量数据进行分解，得到三个因子矩阵，其过程如图3所示。
S3.1：首先将张量通过Tucker分解压缩为/>其中F1＝K+1+P，F2＝min{F1,4}，P<V-K-2为非负整数，通常取1或2。
S3.2：对进行平行因子分解，简称CP分解。首先，随机生成20组CP分解的初始因子矩阵/>分别按以下更新规则迭代5次：
其中，和⊙分别代表Hadamard积和Khatri-Rao积。计算迭代5次以后的误差值
令将第Q组初始化因子矩阵迭代五次以后的结果记为{Ag,Bg,Cg}，以此作为新的初始化因子矩阵，继续按照上述过程迭代更新，直到收敛。
S3.3：令A＝UAg，B＝VBg，C＝WCg作为对进行CP分解的初始值，对/>只能找上述更新步骤进行CP分解，直至/>即可得到/>的三个因子矩阵/> 和/>
S4：利用各因子矩阵和各正向切片的对应关系，对目标角度和功率进行求解和匹配。
因子矩阵和的正向切片关系如下：
其中，表示阵元位置为{qd|0≤q≤V-1,q是整数}的均匀线阵在θ方向上的导向矢量；携带功率信息，/>携带角度信息， 表示/>的第i行；/>表示将向量转化为对角矩阵的操作。
S4.1：解缩放模糊。张量分解得到的存在缩放模糊和置换模糊，根据/>中每个元素的模为1，可以令 除去缩放模糊，使得/>
其中，∏是同用表示置换模糊的未知列置换矩阵，是/>的逆算子，表示将矩阵的对角元素排列成列向量的操作。则存在置换模糊的功率和角度估计为
其中，ik∈{0,1,…,K}且各不相同，k＝1,2,…,K+1。
S4.2：解置换模糊，即目标信源匹配。
已知目标角度先验区域Ψ＝[θp-Δθ,θp+Δθ]，以及/>
·当时：
1.则/>
2.其余元素作为干扰信号角度和功率的估计和/>重构干扰加噪声协方差矩阵/>
·当时：
1.选取中较大的K个元素及其在/>中的对应元素作为干扰信号功率和方向的估计重构/>
2.计算目标信号协方差矩阵则目标来波方向上的导向矢量可估计为/> 表示/>的主特征值。
S5：利用估计的参数，重构波束成形器，计算波束成形权矢量。
权矢量的计算公式为：
为了验证本发明提出的算法性能，通过四个仿真实验进行说明。实验中互质阵列的参数设置为M＝3，N＝5，d设置为半波长。设置目标的真实角度θ0＝5°，目标角度先验区域为Ψ＝[θp-4°,θp+4°]。设置两种干扰场景，其一包括角度为[-17°,25°]的两个干扰源，其二包括角度为[-25°,10°,15°,35°]的四个干扰源，两种场景的干噪比(Interference-to-NoiseRatio,INR)均为10dB。
仿真实验一从波束图的角度，将本发明所提算法得到的波束图与理想波束图相比，验证了本发明对目标信号的保持和对干扰信号的抑制能力。本实验设置快拍数为100，且分为两组，分别设置输入信噪比(Signal-to-NoiseRatio)为5dB和-30dB。在每一组中又细分为无先验误差和2°先验误差(|θ0-θp|＝2°)两种场景。图4给出了在SNR为5dB且存在2个干扰的情况下，无先验误差和2°先验误差分别对应的波束图。图5给出了存在4个干扰，其余条件相同时的波束图。从这两幅图可以看出，当目标信号功率不远远小于干扰信号时，不管是否存在先验误差，本发明的波束图和理想波束图在信号和干扰方向基本是一致的，说明本发明提出的算法的有效性和对先验误差的鲁棒性。图6给出了在SNR为-30dB且存在2个干扰的情况下，无先验误差和2°先验误差分别对应的波束图。图7给出了存在4个干扰，其余条件相同时的波束图。从图中可以看出，当不存在先验误差时，本发明的波束图仍与理想波束图在信号和干扰方向保持高度一致。当存在2°先验误差时，虽然本发明在信号方向的增益比理想波束图下降了约3dB，但是仍可以对干扰方向造成超过45dB的功率衰减，因此仍可以达到一定的保持信号抑制干扰的效果。
仿真实验二通过比较本发明所提波束成形算法与对比方法在不同输入信噪比下的输出SINR，对本发明性能的优越性进行了验证。具体地，对比方法包括文件1，文件3“Onrobust Capon beamforming and diagonal loading[J]”(IEEE Trans.SignalProcessing,2003,51(7):1702-1715)，文件4“Robustadaptivebeamformingusingworst-caseperformance:Asolutiontothe signal mismatch problem[J]”(IEEE Trans.SignalProcessing,2003,51(2):313-324)，文件5“Robust Adaptive Beamforming Based onInterference Covariance Matrix Reconstruction and SteeringVectorEstimation[J]”(IEEETrans.SignalProcessing,2012,60(7):3881-3885)，对比文件6“RobustAdaptive Beamforming Based on Steering Vector Estimation and CovarianceMatrix Reconstruction[J]”(IEEE Communications Letters,2015,19(9):1636-1639)，文件7“Coprime arrayadaptive beamforming with enhanced degrees-of-freedomcapability[J]”(2017 IEEE RadarConference(RadarConf),2017:1357-1361)所提的波束成形方法。本实验分为三个场景，包括无先验误差，2°先验误差以及-4°～4°随机先验误差。图8给出了在快拍数为100且存在4个干扰的情况下，三种场景所对应地输出SINR随着输入SNR变化(-20dB～20dB)的曲线。从图8可以看出，当无先验误差时，本发明基本能达到理想情况下的性能；对比文件3与对比文件4的方法在低SNR条件下与本发明性能接近，但是当SNR>0时，由于自相消效应，其性能下降明显；对比文件1和对比文件7在低SNR条件下性能较差，但随着SNR增加，其性能逐渐接近本发明。当存在2°先验误差以及-4°～4°随机先验误差时，所有算法的性能均有所下降。其中，对比方法6在SNR为-15dB～-3dB的范围内略优于本发明，在其余SNR范围，本发明均取得最好的性能。当SNR较高时，本发明与对比文件基本达到理想性能，验证了其在较高SNR条件下对先验误差的鲁棒性。
仿真实验三通过比较本发明所提波束成形算法与对比方法在不同采样快拍数下的输出SINR，对本发明性能的优越性进行了验证。其中包括的对比方法与实验二相同。本实验同样分为无先验误差，2°先验误差以及-4°～4°随机先验误差三个场景。图9给出了在SNR为0dB且存在4个干扰的情况下，三种场景所对应地输出SINR随着采样快拍数变化(25～250)的曲线。从图9可以看出，当无先验误差时，本发明性能明显由于对比方法；其中本发明与对比文件4所提方法较于其他方法在低信噪比条件下更为鲁棒，其余方法在低信噪比条件下性能明显下降。当存在2°先验误差时，所有算法的性能均有所下降，但本发明在整个快拍范围内均取得最佳性能，而对比文件2和对比文件3下降最为明显。当存在-4°～4°随机先验误差时，本发明与对比文件1的方法取得相近的性能，并优于其余对比方法。
仿真实验四通过比较本发明与其余插值类方法(即对比文件1和对比文件2)的计算时间来验证本发明在计算复杂度上的优越性。图10给出了在四种不同结构互质阵列下，三种方法分别在2干扰和4干扰情况下的计算时间。可以看到，本发明的计算量明显低于其余两种方法，且随着阵列的逐渐变大而变得更显著，证明了本发明在计算复杂度上的优势。此外，两种场景下各算法的计算时间相近，说明干扰个数对计算时间的影响较小。
应当理解，本发明各个实施例中的功能单元模块可以集中在一个处理单元中，也可以是各个单元模块单独物理存在，也可以是两个或两个以上的单元模块集成在一个单元模块中，可以采用硬件或软件的形式来实现。
本发明技术方案还提供了以下实例：
一种基于低管秩张量分解的互质阵列波束成形系统，包括：
虚拟阵元信号初始化与重排单元：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并将初始化内插后的虚拟阵元信号重排为多快拍虚拟阵元信号矩阵；
张量初始化与补全单元：将多快拍虚拟阵元信号矩阵重排为初始化张量，并基于低管秩特性，对初始化张量进行补全；
张量分解单元：对补全后的张量进行分解，得到三个因子矩阵；
目标角度与功率估计单元：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标来波信号角度和功率进行求解和匹配，获得干扰来波信号角度和功率的估计值；
波束成形单元：利用目标角度与功率估计单元获得的参数估计值，重构波束成形器，计算波束成形权矢量。
补全后的张量采用ADMM算法求解单元进行迭代求解。
本发明技术方案还提供了以下实例：
一种计算机存储介质，包括计算机程序，所述计算机程序指令被处理终端执行时使所述处理终端执行一种基于低管秩张量分解的互质阵列波束成形方法。
最后应当说明的是：以上实施例仅用以说明本发明的技术方案而非对其限制，尽管参照上述实施例对本发明进行了详细的说明，所属领域的普通技术人员应当理解：依然可以对本发明的具体实施方式进行修改或者等同替换，而未脱离本发明精神和范围的任何修改或者等同替换，其均应涵盖在本发明的权利要求范围当中。

Claims (8)

1.一种基于低管秩张量分解的互质阵列波束成形方法，其特征在于，包括：
步骤1：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并将初始化内插后的虚拟阵元信号重排为多快拍虚拟阵元信号矩阵；
所述互质阵列下的远场窄带非相干信号模型如下：
其中，x(l)是指互质阵列在第l个采样时刻的信号，a(θ)是互质阵列在θ方向的导向矢量，λ表示接收信号的波长，d表示单位阵元间距，且d≤λ/2，qid表示互质阵列中从左至右第i个阵元的位置，互质阵列由两个均匀子阵组成，第一个均匀子阵由2M个阵元组成，阵元间距为Nd；第二个均匀子阵由N个阵元组成，阵元间距为Md，M和N是一对互质的正整数；s0(l)表示第l个采样时刻的目标源信号，sk(l)表示l个采样时刻的第k个干扰的源信号，k＝1,2,…,K，K表示干扰个数；n(l)表示独立同分布的零均值复高斯噪声矢量，即/> 表示复高斯分布，pn表示噪声功率，I表示单位矩阵；
步骤2：将多快拍虚拟阵元信号矩阵重排为初始化张量，并基于低管秩特性，对初始化张量进行补全；
基于低管秩特性，对初始化的张量进行补全的过程如下：
其中，表示补全后的张量，/>表示/>的张量核范数；ε表示填补缺失数据的张量，Ω是缺失数据的位置索引，Ωc为在/>张量结构下Ω的补集，/>表示ε向Ωc的投影；/>为初始化的张量；
步骤3：对补全后的张量进行分解，得到三个因子矩阵；
步骤4：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标来波信号角度和功率进行求解和匹配，分别获得目标信号和干扰信号的角度和功率的估计值；
步骤5：利用步骤4的参数估计值，重构波束成形器，计算波束成形权矢量，进行波束成形。
2.根据权利要求1所述的方法，其特征在于，补全后的张量采用ADMM算法进行迭代求解，具体过程如下：
设定增广拉格朗日函数，并进行参数初始化设置：μ0＝10-4：
其中，表示拉格朗日乘子，μ表示惩罚系数，μ>0，以使得设定的增广拉格朗日函数最小化时，依次对/>ε,/>μ求解，求解过程中，令非求解变量视为常量：
步骤A1：令并求解得到
其中，表示/>的第i个正向切片；/> 表示从第三个维度对x进行傅里叶变换；/>表示/>的第i个正向切片；/>表示从第三个维度对Wt进行傅里叶变换；i＝1,2,3,4；/> 和/>来自于/>的截断奇异值分解，/>表示/>的共轭转置矩阵；
步骤A2：
步骤A3：
步骤A4：μt+1＝min(ρμt,μmax)，μmax设置为105，ρ表示每次更新的放大倍数，取值为1.1，迭代上述步骤直至x收敛或者迭代次数小于500，获得x的解算值。
3.根据权利要求1所述的方法，其特征在于，各因子矩阵和x的正向切片关系如下：
其中，X(1)、X(2)、X(3)、X(4)分别表示的四个正向切片，/> 表示/>的共轭转置矩阵，θ0表示目标来波方向，θ1,…,θK表示K个干扰来波方向；表示阵元位置为{qd|0≤q≤V-1,q是整数}的均匀线阵在θ方向上的导向矢量；P表示携带功率信息，/>Φ表示携带角度信息，/> 分别表示对补全后的张量进行分解得到的三个因子矩阵；/>表示/>的第i行，i＝1,2,3,4；/>表示将向量转化为对角矩阵的操作。
4.根据权利要求3所述的方法，其特征在于，对目标来波信号角度和功率进行求解和匹配，获得干扰来波信号角度和功率的估计值的具体过程如下：
步骤B1：解缩放模糊；
基于中每个元素的模为1，令/> 除去缩放模糊，使得
其中，Π是同用表示置换模糊的未知列置换矩阵，d(·)是的逆算子，表示将矩阵的对角元素排列成列向量的操作；
则置换模糊的功率和角度的估计值为
其中，ik∈{0,1,…,K}且各不相同，k＝1,2,…,K+1；
步骤B2：解置换模糊，即目标信源匹配；
已知目标角度先验区域Ψ＝[θp-Δθ,θp+Δθ]，θp表示先验角度，Δθ表示角度置信半径；令以及/>分别表示信号源功率的估计值的集合以及对应的角度的估计值集合；当/>时：
令则/>将/>和/>的其余元素作为干扰信号角度和功率的估计/>和/>重构干扰加噪声协方差矩阵，即 表示干扰加噪声协方差矩阵，I表示单位矩阵，/>表示噪声功率的估计值，/>是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵；
当时：
选取中较大的K个元素及其在/>中的对应元素作为干扰来波信号角度和功率的估计值，重构/>再计算目标信号协方差矩阵：/> 表示目标信号协方差矩阵，则目标来波方向上的导向矢量可估计为/>其中，/>x(l)表示阵列接收信号在第l时刻的采样值，xH(l)是x(l)的共轭转置矩阵，L表示总的采样时刻数，/>表示/>的主特征值。
5.根据权利要求1所述的方法，其特征在于，波束成形权矢量的计算公式如下：
其中， 和/>是对干扰信号功率和角度的估计值，是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵，/>是互质阵列在/>方向的导向矢量，/>是/>的共轭转置矩阵，/>表示噪声功率的估计值，I表示单位矩阵。
6.一种基于低管秩张量分解的互质阵列波束成形系统，其特征在于，包括：
虚拟阵元信号初始化与重排单元：建立互质阵列下的远场窄带非相干信号模型，对虚拟阵元信号进行初始化内插，并将初始化内插后的虚拟阵元信号重排为多快拍虚拟阵元信号矩阵；
所述互质阵列下的远场窄带非相干信号模型如下：
其中，x(l)是指互质阵列在第l个采样时刻的信号，a(θ)是互质阵列在θ方向的导向矢量，λ表示接收信号的波长，d表示单位阵元间距，且d≤λ/2，qid表示互质阵列中从左至右第i个阵元的位置，互质阵列由两个均匀子阵组成，第一个均匀子阵由2M个阵元组成，阵元间距为Nd；第二个均匀子阵由N个阵元组成，阵元间距为Md，M和N是一对互质的正整数；s0(l)表示第l个采样时刻的目标源信号，sk(l)表示l个采样时刻的第k个干扰的源信号，k＝1,2,…,K，K表示干扰个数；n(l)表示独立同分布的零均值复高斯噪声矢量，即/> 表示复高斯分布，pn表示噪声功率，I表示单位矩阵；
张量初始化与补全单元：将多快拍虚拟阵元信号矩阵重排为初始化张量，并基于低管秩特性，对初始化张量进行补全；
基于低管秩特性，对初始化的张量进行补全的过程如下：
其中，表示补全后的张量，/>表示/>的张量核范数；ε表示填补缺失数据的张量，Ω是缺失数据的位置索引，Ωc为在/>张量结构下Ω的补集，/>表示ε向Ωc的投影；/>为初始化的张量；
张量分解单元：对补全后的张量进行分解，得到三个因子矩阵；
目标角度与功率估计单元：利用各因子矩阵和补全后的张量各正向切片的对应关系，对目标来波信号角度和功率进行求解和匹配，获得干扰来波信号角度和功率的估计值；
波束成形单元：利用目标角度与功率估计单元获得的参数估计值，重构波束成形器，计算波束成形权矢量。
7.根据权利要求6所述的系统，其特征在于，补全后的张量采用ADMM算法求解单元进行迭代求解。
8.一种计算机存储介质，包括计算机程序，其特征在于，所述计算机程序指令被处理终端执行时使所述处理终端执行权利要求1至5任一项所述的方法。
CN202210461073.5A 2022-04-28 2022-04-28 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质 Active CN114844544B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210461073.5A CN114844544B (zh) 2022-04-28 2022-04-28 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210461073.5A CN114844544B (zh) 2022-04-28 2022-04-28 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质

Publications (2)

Publication Number Publication Date
CN114844544A CN114844544A (zh) 2022-08-02
CN114844544B true CN114844544B (zh) 2024-05-14

Family

ID=82567063

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210461073.5A Active CN114844544B (zh) 2022-04-28 2022-04-28 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质

Country Status (1)

Country Link
CN (1) CN114844544B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116184332B (zh) * 2023-04-21 2023-07-11 南京隼眼电子科技有限公司 一种雷达干扰的抑制方法、装置及存储介质

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110031794A (zh) * 2019-04-16 2019-07-19 中国人民解放军国防科技大学 一种基于差分共性阵重构的相干信源doa估计方法
CN111739551A (zh) * 2020-06-24 2020-10-02 广东工业大学 一种基于低秩与稀疏张量分解的多通道心肺音去噪系统
CN112731275A (zh) * 2020-12-28 2021-04-30 杭州电子科技大学 一种基于零化插值的互质阵部分极化信号参数估计方法
WO2021139208A1 (zh) * 2020-01-08 2021-07-15 华南理工大学 一种基于特定频率组合信号的一维doa估计方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2021068494A1 (zh) * 2020-05-03 2021-04-15 浙江大学 基于平面互质阵列虚拟域张量空间谱搜索的高分辨精确二维波达方向估计方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110031794A (zh) * 2019-04-16 2019-07-19 中国人民解放军国防科技大学 一种基于差分共性阵重构的相干信源doa估计方法
WO2021139208A1 (zh) * 2020-01-08 2021-07-15 华南理工大学 一种基于特定频率组合信号的一维doa估计方法
CN111739551A (zh) * 2020-06-24 2020-10-02 广东工业大学 一种基于低秩与稀疏张量分解的多通道心肺音去噪系统
CN112731275A (zh) * 2020-12-28 2021-04-30 杭州电子科技大学 一种基于零化插值的互质阵部分极化信号参数估计方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
High-order extended coprime array design for direction of arrival estimation;Shi Junpeng;《Journal of systems engineering and electronics》;20210831;全文 *
基于稀疏和低秩恢复的稳健DOA估计方法;王洪雁;于若男;;电子与信息学报;20200315(03);全文 *
新体制机载雷达波形优化设计研究综述;李志汇;《系统工程与电子技术》;20220330;全文 *

Also Published As

Publication number Publication date
CN114844544A (zh) 2022-08-02

Similar Documents

Publication Publication Date Title
CN108732549B (zh) 一种基于协方差矩阵重构的阵元缺损mimo雷达doa估计方法
CN108680891B (zh) 非均匀噪声条件下考虑互耦效应的doa估计方法
CN109254261B (zh) 基于均匀圆阵epuma的相干信号零陷加深方法
CN107576931B (zh) 一种基于协方差低维度迭代稀疏重构的相关/相干信号波达方向估计方法
CN111049556B (zh) 一种基于干扰协方差矩阵重构的互素阵稳健自适应波束形成方法
CN106788655B (zh) 互耦条件下未知互耦信息的干扰相干稳健波束形成方法
CN112731275B (zh) 一种基于零化插值的互质阵部分极化信号参数估计方法
CN111239678A (zh) 一种基于l型阵列的二维doa估计方法
CN114844544B (zh) 一种基于低管秩张量分解的互质阵列波束成形方法、系统及介质
CN109031227B (zh) 一种共轭梯度空时自适应处理方法及系统
CN111337873A (zh) 一种基于稀疏阵的doa估计方法
CN113219402A (zh) 基于Modified-ALM算法的稀疏阵列DOA估计方法
CN114884841A (zh) 基于高阶统计和非均匀阵列的欠定参数联合估计方法
CN104459635A (zh) 基于迭代收缩加权融合的自适应空滤滤波方法
CN113673419A (zh) 适用于强干扰环境的波束域快速稀疏贝叶斯方位估计方法
CN106685555B (zh) 基于低秩矩阵恢复的mimo水声系统信道状态信息反馈方法
CN115356678B (zh) 基于dpnalm算法的稀疏阵列doa估计方法
CN106877918B (zh) 互耦条件下稳健自适应波束形成方法
Blunt et al. A new framework for direction-of-arrival estimation
CN115587281A (zh) 基于因子矩阵先验的阵元失效mimo雷达角度估计方法
Chen et al. Finite data performance analysis of LCMV antenna array beamformers with and without signal blocking
CN114184999B (zh) 一种互耦小孔径阵列的生成式模型处理方法
CN114236471A (zh) 一种相关干扰源下的稳健自适应波束形成方法
CN114047481A (zh) 一种基于子空间正交性的稳健自适应波束形成方法
CN110958045A (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