CN113504505B - 一种适用于低信噪比环境下的一维doa估计方法 - Google Patents
一种适用于低信噪比环境下的一维doa估计方法 Download PDFInfo
- Publication number
- CN113504505B CN113504505B CN202110614512.7A CN202110614512A CN113504505B CN 113504505 B CN113504505 B CN 113504505B CN 202110614512 A CN202110614512 A CN 202110614512A CN 113504505 B CN113504505 B CN 113504505B
- Authority
- CN
- China
- Prior art keywords
- matrix
- signal
- vector
- convex optimization
- noise
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 51
- 239000011159 matrix material Substances 0.000 claims abstract description 160
- 239000013598 vector Substances 0.000 claims abstract description 78
- 238000005457 optimization Methods 0.000 claims abstract description 67
- 230000009467 reduction Effects 0.000 claims abstract description 28
- 238000005070 sampling Methods 0.000 claims abstract description 18
- 230000008569 process Effects 0.000 claims abstract description 15
- 238000012545 processing Methods 0.000 claims abstract description 7
- 238000000354 decomposition reaction Methods 0.000 claims description 14
- 238000001228 spectrum Methods 0.000 claims description 11
- 230000003595 spectral effect Effects 0.000 claims description 6
- 238000010276 construction Methods 0.000 claims description 3
- 230000007423 decrease Effects 0.000 claims description 3
- 230000000694 effects Effects 0.000 abstract description 4
- 230000007547 defect Effects 0.000 description 3
- 230000005540 biological transmission Effects 0.000 description 2
- 238000004891 communication Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000001629 suppression 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
- 238000013500 data storage Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 239000002699 waste material Substances 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/80—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 ultrasonic, sonic or infrasonic waves
- G01S3/86—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 ultrasonic, sonic or infrasonic waves with means for eliminating undesired waves, e.g. disturbing noises
-
- 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
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
Abstract
本发明公开了一种适用于低信噪比环境下的一维DOA估计方法,该方法通过对均匀线阵的接收信号进行处理,在对接收信号的协方差矩阵做降维处理的同时减去噪声信号的特征值达到降噪的效果,并在稀疏向量解的凸优化过程中迭代地更新加权矩阵的系数来逐步修正角度位置的偏差,减小了角度的估计误差,提高压缩感知类算法在低信噪比环境下的DOA估计精度。本发明基于压缩感知估计理论对接收信号进行采样,使得采样频率不受奈奎斯特定理的限制,通过对接收数据进行降噪预处理和迭代更新权重系数,提高了水下DOA的估计精度,在低信噪比环境下具有较好的降噪效果。
Description
技术领域
本发明涉及水下目标定位技术领域,具体涉及一种适用于低信噪比环境下的一维DOA估计方法。
背景技术
波达方向估计(DOA估计)是一种根据阵列接收信号获取目标信号的空域位置信息的技术,在水下目标定位等领域有着重要的应用和发展。传统的DOA估计算法在信号采样过程中采样频率需要高于奈奎斯特频率,这对于采样后信号的存储空间和传输过程带来一定的负担和资源浪费。近年来提出的基于压缩感知的DOA估计算法可以突破传统采样机制对于信号频率的限制,大大降低了信号采样的最低频率,从而弥补传统算法的不足。现有的平滑L0压缩感知算法利用平滑的L0范数优化稀疏重构过程,提高了信号稀疏重构的精度的同时具有较低的运算复杂度;还有一种加权的压缩感知算法利用噪声子空间与信号特征向量间的正交性构建加权矩阵,代入到凸优化过程中以提高角度的估计精度。但以上算法在快拍数较少的情况下对于多信源波达方向角的分辨率较低。
此外在水声通信领域中,由于环境因素的复杂性,水声通信往往处于信噪比较低的环境,压缩感知类算法在低信噪比环境下存在估计性能下降的问题,尤其对于色噪声信号的抑制效果较差。目前亟待提出一种DOA估计方法解决上述算法的缺陷以及在低信噪比环境下的估计精度问题。
发明内容
本发明的目的是为了解决现有技术中的上述缺陷,提供一种适用于低信噪比环境下的一维DOA估计方法。该方法通过对均匀线阵的接收信号进行处理,在对接收信号的协方差矩阵做降维处理的同时减去噪声信号的特征值达到降噪的效果,并在稀疏向量解的凸优化过程中迭代地更新加权矩阵的系数来逐步修正角度位置的偏差,最终减小角度的估计误差,提高压缩感知类算法在低信噪比环境下的DOA估计精度。
本发明的目的可以通过采取如下技术方案达到:
一种适用于低信噪比环境下的一维DOA估计方法,所述估计方法包括以下步骤:
S1、建立一维的均匀线阵的阵列信号模型,经过L次采样得到接收信号矩阵Y,并根据阵列流型矩阵A(θ)构建凸优化过程所需的观测矩阵
S2、对接收信号矩阵Y进行奇异值分解,得到从大到小排序的M个奇异值,对后M-K个奇异值取平均值作为噪声信号的平均功率;
S3、降低接收信号和稀疏信号的矩阵维数,同时对接收信号矩阵Y做降噪处理,并初始化加权矩阵为W(0)=IN,令i=1,i表示迭代轮数,开始迭代更新权重系数,IN为N阶单位矩阵,N为空域划分的角度数;
S4、对凸优化方程求解得到最优解并根据/>构造虚拟的接收信号矩阵/>对/>做奇异值分解,得到第i轮的噪声子空间/>
S5、利用噪声子空间与观测矩阵中目标信号对应矢量的正交关系更新第i轮迭代的加权矩阵W(i),令i=i+1,重复步骤S4到S5,直到收敛条件终止迭代;
S6、输出最后一次的凸优化方程解通过谱峰搜索获取前K个峰值相应的位置即为角度估计值。
进一步地,所述步骤S1的过程如下:假设均匀线阵上排布M个接收阵元,相邻阵元的间距为d,设有K个远场窄带声波信号入射到该均匀线阵上,阵元个数M和信源个数K满足关系K<M,每个信号对应的波长为λk,k=1,2,…,K,真实信号的波达方向角分别θ={θ1,θ2,…,θK},定义为信号入射方向与线阵所在直线的法线的夹角,接收阵元的采样信号y以列向量形式表示为:
其中,表示阵列流型矩阵,/>表示K个信源对应的目标信号矢量,/>表示M个阵元的噪声信号矢量,/>表示复数域,经过L次采样后得到的接收信号矩阵Y表达式如公式(2)所示:
Y=A(θ)S+G 公式(2)
其中,S表示包含K个信源信息的K×L维目标信号矩阵,G表示M×L维的噪声矩阵,A(θ)表示M×K维的阵列流型矩阵,其表达式为:
通过压缩感知类算法求解凸优化方程获得的N×K维稀疏解,该稀疏解中包含K个真实信源的角度信息,其中,凸优化方程的构建需要根据阵列流型矩阵中的元素构建M×N维的观测矩阵该观测矩阵/>的表示形式为:
其中,N表示接收阵元上空域划分的角度个数,满足M<<N,λ表示信号的中心频率对应的波长值;
根据压缩感知理论,通过对接收信号矩阵Y进行凸优化处理得到稀疏表示信号的最优解,该最优解中包含目标信号的入射角度,凸优化方程如公式(5)所示:
在凸优化方程中的观测矩阵基于阵列流型矩阵A(θ)中原子的参数构建;
针对接收信号模型中第一个阵元的接收信号,将目标信号矢量扩展成包含K个真实信号的N×1维虚拟原始信号矢量/>N表示空域划分的角度数,角度集合Θ所对应的虚拟原始信号矢量/>是一个K-稀疏信号,表明/>是可压缩的,在原始信号的稀疏解中包含来波方向角信息的标量位于虚拟原始信号矢量/>的非零元素位置。
进一步地,所述步骤S2过程如下:
首先对接收信号矩阵Y进行奇异值分解构建降噪的奇异值矩阵Σ′:
Y=UΣVH 公式(6)
Σ′=Σ-ΣN 公式(7)
其中,表示主对角线元素为奇异值构成的矩阵,/>表示实数域,表示Σ的主对角线元素中第K+1个到第M个奇异值之和的平均值作为噪声信号的平均功率。
进一步地,所述步骤S3过程如下:
对接收信号矩阵Y降维处理的同时,减去噪声平均功率得到处理后的接收信号矩阵YSVW:
YSVW=UΣ′DK 公式(9)
其中,凸优化方程中待估计信号的L1范数替代公式(5)中的L0范数作为衡量最优结果的参数,将约束条件转换为误差项并加入加权矩阵,凸优化方程变为:
其中,表示使凸优化方程解为最小值的XSV取值,/>和/>分别表示稀疏表示信号的初始值和凸优化方程解,/>表示矩阵的Frobenius范数的平方,||·||2,1表示矩阵中各列向量L2范数组成向量的L1范数值,/>表示加权矩阵,公式(11)对应的二阶锥规划问题表示为公式(12):
min p+ηq 公式(12)
min表示求解p+ηq的最小值,其约束条件为存在K个N×1维的矢量[Z1,Z2,…,Zk,…,ZK]和N×1维的矢量[γ1,γ2,…,γj,…,γN]T满足关系且||γ1,γ2,…,γj,…,γN||1≤q,Zk和γn的表达式为:
||WXSV(j,:)||2≤γj,j=1,2,…,N 公式(14)
其中,p、q为约束条件中两个范数和的临时变量,η表示正则化参数因子,用来平衡误差和稀疏性之间的关系,表示经过降噪和降维处理后的接收信号矩阵,YSVW(:,k)为YSVW的第k列列向量,XSV(:,k)表示XSV的第k列列向量,/>表示/>的第j行行向量,/>表示加权矩阵,/>表示的Zk共轭向量,||·||1、||·||2和/>分别表示向量的L1范数和L2范数以及Frobenius范数的平方,同时开始迭代前的加权矩阵W(0)初始化为N×N维的单位矩阵IN,重新令i=1,开始进入第一轮迭代。
进一步地,所述步骤S4过程如下:
将第i-1轮迭代的加权矩阵W(i-1)代入公式(11)中进行凸优化方程的求解,得到第i轮凸优化方程解根据其估计值构建出虚拟的接收信号矩阵/>并对/>进行奇异值分解:
其中,U(i)、Σ(i)、V(i)分别为第i轮迭代中的左奇异矩阵、奇异值矩阵和右奇异矩阵,U(i)矩阵中的前K列子矩阵/>和后M-K列子矩阵/>分别对应信号子空间和噪声子空间。
进一步地,所述步骤S5过程如下:
利用公式(16)得到的噪声子空间和包含真实信号角度信息的观测矩阵列向量,结合MUSIC算法的原理构造第i轮迭代的加权矩阵W(i):
其中,表示第i轮加权矩阵中的对角线上元素,/>表示第i轮迭代中噪声子空间的共轭矩阵,a(θj),j=1,2,…,N表示观测矩阵/>中第j列的列向量,根据噪声子空间与真实信号矢量的正交性,在真实信号角度对应列向量的加权系数会增大其谱峰,而其它位置的列向量对应谱峰值则会减小,令i=i+1,进入下一轮迭代,重复步骤S4和步骤S5,并将本轮迭代更新的加权矩阵W(i)代入下一轮迭代的凸优化方程中,直到满足迭代终止条件为止。
进一步地,所述步骤S5中迭代更新的终止条件包括以下两种:
(1)权重更新达到足够的迭代次数后终止;
(2)第i轮迭代中的凸优化方程解与第i-1轮迭代中的凸优化方程解/>之差的L2范数小于给定的误差参数时终止;
其中,ε为给定的误差参数,根据实际需要取值。
进一步地,所述步骤S6过程如下:
获取第i轮的凸优化方程解根据/>构建稀疏表示信号的空间谱密度函数Q:
其中,表示/>的第一列列向量,/>表示该向量的L1范数,通过对Q进行谱峰搜索获取前K个峰值的位置,即对应虚拟的原始信号矢量/>中真实信号入射角度的估计值/>
本发明相对于现有技术具有如下的优点及效果:
1、本发明公开的DOA估计方法基于压缩感知理论对接收信号进行压缩采样,相比传统的子空间类算法,该方法的采样频率不受奈奎斯特定理的限制,降低了采样后的数据存储和传输所消耗的资源,且通过稀疏重构即可恢复出原始信号。
2、本发明公开的DOA估计方法对接收信号矩阵进行了降维和降噪处理,在降低方法的运算复杂度的同时也提高了对噪声信号的抑制效果,在低信噪比的白噪声和色噪声环境下的DOA估计都具有更好的降噪性能。
3、本发明公开的DOA估计方法相比现有的加权凸优化算法具有精度更高的加权系数,通过迭代更新权重缩减了凸优化方程最优解的估计误差,进一步提高了算法的估计精度。
4、本发明公开的DOA估计方法对于多信源情况下的相邻入射角估计时的分辨力更好,相比传统的压缩感知类算法在抗模糊性上具有更高的性能。
附图说明
图1是本发明实施例提出的一种适用于低信噪比环境下的一维DOA估计方法中均匀线阵模型示意图;
图2是本发明实施例提出的一种适用于低信噪比环境下的一维DOA估计方法中接收信号模型示意图;
图3是本发明实施例提出的一种适用于低信噪比环境下的一维DOA估计方法的流程图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
实施例一
本实施例公开了一种适用于低信噪比环境下的一维DOA估计方法,所述估计方法包括以下步骤:
S1、建立一维的均匀线阵的阵列信号模型,经过L次采样得到接收信号矩阵Y,并根据阵列流型矩阵构建凸优化过程所需的观测矩阵
本实施例步骤S1的过程如下:
采用均匀线阵为接收阵列和基于远场窄带信号模型进行DOA估计,均匀线阵模型如附图1所示,假设均匀线阵上排布M个接收阵元,相邻阵元的间距为d,设有K个远场窄带声波信号入射到该均匀线阵上,阵元个数M和信源个数K满足关系K<M,每个信号对应的波长为λk,k=1,2,…,K,真实信号的波达方向角分别θ={θ1,θ2,…,θK},定义为信号入射方向与线阵所在直线的法线的夹角,接收阵元的采样信号y以列向量形式表示为:
其中,表示阵列流型矩阵,/>表示K个信源对应的目标信号矢量,/>表示M个阵元的噪声信号矢量,/>表示复数域,经过L次采样后得到的接收信号矩阵Y表达式如公式(2)所示:
Y=A(θ)S+G 公式(2)
其中,S表示包含K个信源信息的K×L维目标信号矩阵,G表示M×L维的噪声矩阵,A(θ)表示M×K维的阵列流型矩阵,其表达式为:
压缩感知类算法通过对接收信号矩阵Y求解凸优化方程获得的N×K维最优解中包含K个真实信源的角度信息,其中,凸优化方程的构建需要根据阵列流型矩阵A(θ)中的元素构建M×N维的观测矩阵/>该观测矩阵/>的表示形式为:
其中,N表示接收阵元上空域划分的角度个数,满足M<<N,且N越大,角度估计的位置越精确,λ表示信号的中心频率对应的波长值;
根据压缩感知理论,通过对接收信号矩阵Y进行凸优化处理得到稀疏表示信号的最优解,该最优解中包含目标信号的入射角度,凸优化方程如公式(5)所示:
在凸优化方程中的观测矩阵基于阵列流型矩阵A(θ)中原子的参数构建。
接收信号模型如附图2所示,针对该接收信号模型中各阵元的接收信号,将目标信号矢量扩展成包含K个真实信号的N×1维虚拟原始信号矢量/>角度集合Θ所对应的虚拟原始信号矢量/>是一个K-稀疏信号,表明/>是可压缩的,在原始信号的稀疏解中包含来波方向角信息的标量位于虚拟原始信号矢量/>的非零元素位置。
S2、对接收信号矩阵Y进行奇异值分解,得到从大到小排序的M个奇异值,对后M-K个奇异值取平均值作为噪声信号的平均功率;
本实施例的步骤S2过程如下:
在凸优化之前需要对接收信号矩阵Y做降维和降噪处理,目的是为了降低运算复杂度的同时抑制噪声的干扰。首先对接收信号矩阵Y进行奇异值分解构建降噪的奇异值矩阵Σ′:
Y=UΣVH 公式(6)
Σ′=Σ-ΣN 公式(7)
其中,表示主对角线元素为奇异值构成的矩阵,/>表示实数域,表示Σ的主对角线元素中第K+1个到第M个奇异值之和的平均值作为噪声信号的平均功率。
S3、降低接收信号和稀疏信号的矩阵维数,同时对接收信号矩阵Y做降噪处理,并在迭代开始前初始化加权矩阵为W(0)=IN,令i=1,开始迭代更新权重系数,IN为N阶单位矩阵,N表示空域划分的角度数;
本实施例的步骤S3过程如下:
对接收信号矩阵Y降维处理的同时,减去噪声平均功率得到处理后的接收信号矩阵YSVW:
YSVW=UΣ′DK 公式(9)
其中,通过右乘DK将列数为N的原接收信号矩阵压缩到只包含前K列信息,且保证包含信源信息的奇异值保留到降维后的矩阵中。相比现有的压缩感知估计方法,步骤S2在对接收信号矩阵Y作降维处理的同时,增加了降噪处理的步骤,目的是为了减少凸优化过程中由于噪声影响导致的估计偏差;另外,凸优化方程中待估计信号的L1范数替代了公式(5)中的L0范数作为衡量最优结果的参数,将约束条件转换为误差项并加入加权矩阵,凸优化方程变为:
其中,表示使凸优化方程解为最小值的XSV取值,/>和/>分别表示稀疏表示信号的初始值和凸优化方程解,/>表示矩阵的Frobenius范数的平方,||·||2,1表示矩阵中各列向量L2范数组成向量的L1范数值,/>表示加权矩阵。公式(11)对应的二阶锥规划问题表示为公式(12):
min p+ηq 公式(12)
min表示求解p+ηq的最小值,其约束条件为存在K个N×1维的矢量[Z1,Z2,…,Zk,…,ZK]和N×1维的矢量[γ1,γ2,…,γj,…,γN]T满足关系且||γ1,γ2,…,γj,…,γN||1≤q。Zk和γn的表达式为:
||WXSV(j,:)||2≤γj,j=1,2,…,N 公式(14)
其中,p、q为约束条件中两个范数和的临时变量,η表示正则化参数因子,用来平衡误差和稀疏性之间的关系,表示经过降噪和降维处理后的接收信号矩阵,YSVW(:,k)为YSVW的第k列列向量,XSV(:,k)表示XSV的第k列列向量,/>表示/>的第j行行向量,/>表示加权矩阵,/>表示的Zk共轭向量,||·||1、||·||2和/>分别表示向量的L1范数和L2范数以及Frobenius范数的平方。同时开始迭代前的加权矩阵W(0)初始化为N×N维的单位矩阵IN,令i=1,开始进入第一轮迭代。
S4、对凸优化方程求解得到最优解并根据/>构造虚拟的接收信号矩阵/>对/>做奇异值分解,得到第i轮的噪声子空间/>
本实施例的步骤S4过程如下:
将加权矩阵W(i-1)代入公式(11)中进行凸优化方程的求解,得到第i轮凸优化方程解根据其估计值构建出虚拟的接收信号矩阵/>并对/>进行奇异值分解:
其中,U(i)、Σ(i)、V(i)分别为第i轮迭代中的左奇异矩阵、奇异值矩阵和右奇异矩阵,U(i)矩阵中的前K列子矩阵/>和后M-K列子矩阵/>分别对应信号子空间和噪声子空间;由于加权系数构成需要由接收信号矩阵,因此通过迭代的方式利用上一轮的凸优化方程解构建虚拟的接收信号矩阵可以作为本轮的加权系数构成,与加权的L1-SVD估计方法相比,能够实现迭代地更新加权系数,进一步提高估计的精度。
S5、利用噪声子空间与观测矩阵中目标信号对应矢量的正交关系更新第i轮加权矩阵W(i),令i=i+1,重复步骤S4到S5,直到收敛条件终止迭代;
利用式(16)得到的噪声子空间和包含真实信号角度信息的观测矩阵列向量,结合MUSIC算法的原理构造第i轮迭代的加权矩阵W(i):
其中,表示第i轮加权矩阵中的对角线上元素,/>表示第i轮迭代中噪声子空间的共轭矩阵,a(θj),j=1,2,…,N表示观测矩阵/>中第j列的向量,根据噪声子空间与真实信号矢量的正交性,在真实信号角度对应列向量的加权系数会增大其谱峰,而其它位置的列向量对应谱峰值则会减小,令i=i+1,进入下一轮迭代,重复步骤S4和步骤S5,并将本轮迭代更新的加权矩阵W(i)代入下一轮迭代的凸优化方程中,直到满足迭代终止条件为止。
迭代更新的终止条件包括以下两种:
(1)权重更新达到足够的迭代次数后终止;
(2)第i轮迭代中的凸优化方程解与第i-1轮迭代中的凸优化方程解/>之差的L2范数小于给定的误差参数时终止;
其中,ε为给定的误差参数,根据实际需要取值。
S6、输出最后一次的凸优化方程解通过谱峰搜索获取前K个峰值相应的位置即为角度估计值。
本实施例的步骤S6过程如下:
获取最后一次的凸优化方程解根据/>构建稀疏表示信号的空间谱密度函数Q:
其中,表示/>的第一列列向量,/>表示该向量的L1范数。通过对Q进行谱峰搜索获取前K个峰值的位置,即对应原始信号矢量/>中真实信号入射角度的估计值/>与加权的L1-SVD算法相比,加权系数对角度位置加成更加准确,在步骤S1的条件下经过仿真对比,加权L1-SVD估计方法的估计偏差为0.18°,本算法的估计偏差为0.07°。
实施例二
本实施例提出一种基于迭代更新权重和一维均匀线阵的水下一维DOA估计方法,具体步骤如下:
T1、建立一维的均匀线阵的阵列信号模型。采用一维均匀线阵为接收阵列,阵元数M=8,设有K个非相关的远场窄带声波信号入射到该均匀线阵上,信源个数K=3,入射角度分别为[-30°,0°,30°],信号的传播速度为c=1500m/s,阵元间距d=0.05m,为信源信号的中心频率15kHz对应波长的一半,阵列单次采样的接收信号形式为:
y=[y1,y2,…,y8]T 公式(21)
经过快拍数L=100次采样得到的接收信号矩阵为:
根据阵列流型矩阵A(θ)表达式构建凸优化方程所需的观测矩阵其表达式如公式(6)所示。其中/>N表示在接收阵元上空域划分的角度数量,且满足关系M<<N。这里根据附图2规定角度搜索的空域范围为[-90°~90°],以均匀间隔为0.1°划分,共划分网格的角度数量N=1801。
T2、对接收信号矩阵Y做奇异值分解得到奇异值矩阵Σ,得到共8个奇异值,从大到小对所有奇异值排列,取后5个奇异值的平均值作为噪声信号的平均功率,并根据公式(8)、公式(9)、公式(10)构建降噪后的奇异值矩阵Σ′。
T3、对原接收信号矩阵Y做降噪处理,减去噪声平均功率并根据公式(11)对接收信号矩阵Y进行降维,得到处理后的接收信号矩阵YSVW。
T4、初始化加权矩阵W(0)=IN。令i=1,开始对加权矩阵进行更新。
T5、根据公式(13)求解凸优化方程得到第i轮凸优化方程解并依此构建第i轮迭代中虚拟的接收信号矩阵/>对/>做奇异值分解得到/>的噪声子空间/>
T6、根据公式(17)、(18)更新第i轮迭代的加权矩阵W(i)。令i=i+1,进入下一轮迭代,重复S4到S6的步骤,并将第i-1轮迭代更新的加权矩阵W(i-1)代入下一轮迭代的凸优化方程中,直到迭代轮数大于等于6。
T7、经过6次权重更新后结束迭代过程,根据第6轮的凸优化方程解构建稀疏表示信号的空间谱密度函数Q,并进行谱峰搜索,得到前3个谱峰值对应的角度为[-29.84°,0.09°,30.02°],对目标估计达到了预期精度,说明估计结果正确,本发明方法可行。
上述实施例为本发明较佳的实施方式,但本发明的实施方式并不受上述实施例的限制,其他的任何未背离本发明的精神实质与原理下所作的改变、修饰、替代、组合、简化,均应为等效的置换方式,都包含在本发明的保护范围之内。
Claims (7)
1.一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述估计方法包括以下步骤:
S1、建立一维的均匀线阵的阵列信号模型,经过L次采样得到接收信号矩阵Y,并根据阵列流型矩阵A(θ)构建凸优化过程所需的观测矩阵
其中,所述步骤S1的过程如下:假设均匀线阵上排布M个接收阵元,相邻阵元的间距为d,设有K个远场窄带声波信号入射到该均匀线阵上,阵元个数M和信源个数K满足关系K<M,每个信号对应的波长为λk,k=1,2,…,K,真实信号的波达方向角分别θ={θ1,θ2,…,θK},定义为信号入射方向与线阵所在直线的法线的夹角,接收阵元的采样信号y以列向量形式表示为:
其中,表示阵列流型矩阵,/>表示K个信源对应的目标信号矢量,/>表示M个阵元的噪声信号矢量,/>表示复数域,经过L次采样后得到的接收信号矩阵Y表达式如公式(2)所示:
Y=A(θ)S+G 公式(2)
其中,S表示包含K个信源信息的K×L维目标信号矩阵,G表示M×L维的噪声矩阵,A(θ)表示M×K维的阵列流型矩阵,其表达式为:
通过压缩感知类算法求解凸优化方程获得的N×K维稀疏解,该稀疏解中包含K个真实信源的角度信息,其中,凸优化方程的构建需要根据阵列流型矩阵中的元素构建M×N维的观测矩阵该观测矩阵/>的表示形式为:
其中,N表示接收阵元上空域划分的角度个数,满足M<<N,λ表示信号的中心频率对应的波长值;
根据压缩感知理论,通过对接收信号矩阵Y进行凸优化处理得到稀疏表示信号的最优解,该最优解中包含目标信号的入射角度,凸优化方程如公式(5)所示:
在凸优化方程中的观测矩阵
基于阵列流型矩阵A(θ)中原子的参数构建;
针对接收信号模型中第一个阵元的接收信号,将目标信号矢量扩展成包含K个真实信号的N×1维虚拟原始信号矢量/>N表示空域划分的角度数,角度集合Θ所对应的虚拟原始信号矢量/>是一个K-稀疏信号,表明/>是可压缩的,在原始信号的稀疏解中包含来波方向角信息的标量位于虚拟原始信号矢量/>的非零元素位置;
S2、对接收信号矩阵Y进行奇异值分解,得到从大到小排序的M个奇异值,对后M-K个奇异值取平均值作为噪声信号的平均功率;
S3、降低接收信号和稀疏信号的矩阵维数,同时对接收信号矩阵Y做降噪处理,并初始化加权矩阵为W(0)=IN,令i=1,i表示迭代轮数,开始迭代更新权重系数,IN为N阶单位矩阵,N为空域划分的角度数;
S4、对凸优化方程求解得到最优解并根据/>构造虚拟的接收信号矩阵/>对做奇异值分解,得到第i轮的噪声子空间/>
S5、利用噪声子空间与观测矩阵中目标信号对应矢量的正交关系更新第i轮迭代的加权矩阵W(i),令i=i+1,重复步骤S4到S5,直到收敛条件终止迭代;
S6、输出最后一次的凸优化方程解通过谱峰搜索获取前K个峰值相应的位置即为角度估计值。
2.根据权利要求1所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S2过程如下:
首先对接收信号矩阵Y进行奇异值分解构建降噪的奇异值矩阵Σ′:
Y=UΣVH 公式(6)
Σ′=Σ-ΣN 公式(7)
其中,表示主对角线元素为奇异值构成的矩阵,/>表示实数域,表示Σ的主对角线元素中第K+1个到第M个奇异值之和的平均值作为噪声信号的平均功率。
3.根据权利要求2所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S3过程如下:
对接收信号矩阵Y降维处理的同时,减去噪声平均功率得到处理后的接收信号矩阵YSVW:
YSVW=UΣ′DK 公式(9)
其中,凸优化方程中待估计信号的L1范数替代公式(5)中的L0范数作为衡量最优结果的参数,将约束条件转换为误差项并加入加权矩阵,凸优化方程变为:
其中,表示使凸优化方程解为最小值的XSV取值,/>和/>分别表示稀疏表示信号的初始值和凸优化方程解,/>表示矩阵的Frobenius范数的平方,||·||2,1表示矩阵中各列向量L2范数组成向量的L1范数值,/>表示加权矩阵,公式(11)对应的二阶锥规划问题表示为公式(12):
min p+ηq 公式(12)
min表示求解p+ηq的最小值,其约束条件为存在K个N×1维的矢量[Z1,Z2,…,Zk,,ZK]和N×1维的矢量[γ1,γ2,…,γj,…,γN]T满足关系且γ1,γ2,…,γj,…,γN||1≤q,Zk和γj的表达式为:
||WXSV(j,:)||2≤γj,j=1,2,…,N 公式(14)
其中,p、q为约束条件中两个范数和的临时变量,η表示正则化参数因子,用来平衡误差和稀疏性之间的关系,表示经过降噪和降维处理后的接收信号矩阵,YSVW(:,k)为YSVW的第k列列向量,XSV(:,k)表示XSV的第k列列向量,/>表示/>的第j行行向量,表示加权矩阵,/>表示的Zk共轭向量,||·||1、||·||2和/>分别表示向量的L1范数和L2范数以及Frobenius范数的平方,同时开始迭代前的加权矩阵W(0)初始化为N×N维的单位矩阵IN,重新令i=1,开始进入第一轮迭代。
4.根据权利要求3所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S4过程如下:
将第i-1轮迭代的加权矩阵W(i-1)代入公式(11)中进行凸优化方程的求解,得到第i轮凸优化方程解根据其估计值构建出虚拟的接收信号矩阵/>并对/>进行奇异值分解:
其中,U(i)、Σ(i)、V(i)分别为第i轮迭代中的左奇异矩阵、奇异值矩阵和右奇异矩阵,U(i)矩阵中的前K列子矩阵/>和后M-K列子矩阵/>分别对应信号子空间和噪声子空间。
5.根据权利要求4所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S5过程如下:
利用式(16)得到的噪声子空间和包含真实信号角度信息的观测矩阵列向量,结合MUSIC算法的原理构造第i轮迭代的加权矩阵W(i):
其中,表示第i轮加权矩阵中的对角线上元素,/>表示第i轮迭代中噪声子空间的共轭矩阵,a(θj),j=1,2,…,N表示观测矩阵/>中第j列的列向量,根据噪声子空间与真实信号矢量的正交性,在真实信号角度对应列向量的加权系数会增大其谱峰,而其它位置的列向量对应谱峰值则会减小,令i=i+1,进入下一轮迭代,重复步骤S4和步骤S5,并将本轮迭代更新的加权矩阵W(i)代入下一轮迭代的凸优化方程中,直到满足迭代终止条件为止。
6.根据权利要求5所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S5中迭代更新的终止条件包括以下两种:
(1)权重更新达到足够的迭代次数后终止;
(2)第i轮迭代中的凸优化方程解与第i-1轮迭代中的凸优化方程解/>之差的L2范数小于给定的误差参数时终止;
其中,ε为给定的误差参数,根据实际需要取值。
7.根据权利要求5所述的一种适用于低信噪比环境下的一维DOA估计方法,其特征在于,所述步骤S6过程如下:
获取第i轮的凸优化方程解根据/>构建稀疏表示信号的空间谱密度函数Q:
其中,表示/>的第一列列向量,/>表示该向量的L1范数,通过对Q进行谱峰搜索获取前K个峰值的位置,即对应虚拟的原始信号矢量/>中真实信号入射角度的估计值/>
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110614512.7A CN113504505B (zh) | 2021-06-02 | 2021-06-02 | 一种适用于低信噪比环境下的一维doa估计方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110614512.7A CN113504505B (zh) | 2021-06-02 | 2021-06-02 | 一种适用于低信噪比环境下的一维doa估计方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113504505A CN113504505A (zh) | 2021-10-15 |
CN113504505B true CN113504505B (zh) | 2023-11-03 |
Family
ID=78009296
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110614512.7A Active CN113504505B (zh) | 2021-06-02 | 2021-06-02 | 一种适用于低信噪比环境下的一维doa估计方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113504505B (zh) |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116009058B (zh) * | 2022-11-30 | 2024-04-19 | 杭州交大仪器设备有限公司 | 一种基于多探头传感器数据的地下管道定位方法 |
CN116879835B (zh) * | 2023-07-25 | 2024-06-25 | 安徽大学 | 一种投影最小最大凹函数波达方向估计方法和装置 |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106772225A (zh) * | 2017-01-20 | 2017-05-31 | 大连大学 | 基于压缩感知的波束域doa估计 |
CN107544052A (zh) * | 2017-08-07 | 2018-01-05 | 大连大学 | 一种基于矩阵补全的二阶统计量重构doa估计方法 |
WO2018094565A1 (zh) * | 2016-11-22 | 2018-05-31 | 深圳大学 | 脉冲噪声下的波束成形方法及装置 |
CN110109053A (zh) * | 2019-04-02 | 2019-08-09 | 华南理工大学 | 一种未知声速环境下快速doa估计方法 |
CN110308417A (zh) * | 2019-05-30 | 2019-10-08 | 电子科技大学 | 基于矩阵填充的嵌套阵阵元失效下的波达方向估计方法及装置 |
CN111190136A (zh) * | 2020-01-08 | 2020-05-22 | 华南理工大学 | 一种基于特定频率组合信号的一维doa估计方法 |
-
2021
- 2021-06-02 CN CN202110614512.7A patent/CN113504505B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2018094565A1 (zh) * | 2016-11-22 | 2018-05-31 | 深圳大学 | 脉冲噪声下的波束成形方法及装置 |
CN106772225A (zh) * | 2017-01-20 | 2017-05-31 | 大连大学 | 基于压缩感知的波束域doa估计 |
CN107544052A (zh) * | 2017-08-07 | 2018-01-05 | 大连大学 | 一种基于矩阵补全的二阶统计量重构doa估计方法 |
CN110109053A (zh) * | 2019-04-02 | 2019-08-09 | 华南理工大学 | 一种未知声速环境下快速doa估计方法 |
CN110308417A (zh) * | 2019-05-30 | 2019-10-08 | 电子科技大学 | 基于矩阵填充的嵌套阵阵元失效下的波达方向估计方法及装置 |
CN111190136A (zh) * | 2020-01-08 | 2020-05-22 | 华南理工大学 | 一种基于特定频率组合信号的一维doa估计方法 |
Also Published As
Publication number | Publication date |
---|---|
CN113504505A (zh) | 2021-10-15 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN111415676B (zh) | 一种基于分离矩阵初始化频点选择的盲源分离方法及系统 | |
CN113504505B (zh) | 一种适用于低信噪比环境下的一维doa估计方法 | |
CN107290709B (zh) | 基于范德蒙分解的互质阵列波达方向估计方法 | |
CN109298383B (zh) | 一种基于变分贝叶斯推断的互质阵波达方向角估计方法 | |
CN110007265A (zh) | 一种基于深度神经网络的波达方向估计方法 | |
CN110045323B (zh) | 一种基于矩阵填充的互质阵稳健自适应波束形成算法 | |
CN109375154B (zh) | 一种冲击噪声环境下基于均匀圆阵的相干信号参数估计方法 | |
CN108872926B (zh) | 一种基于凸优化的幅相误差校正及doa估计方法 | |
CN107576931B (zh) | 一种基于协方差低维度迭代稀疏重构的相关/相干信号波达方向估计方法 | |
CN107544051A (zh) | 嵌套阵列基于k‑r子空间的波达方向估计方法 | |
CN112731273A (zh) | 一种基于稀疏贝叶斯的低复杂度信号波达方向估计方法 | |
CN110895325B (zh) | 基于增强四元数多重信号分类的到达角估计方法 | |
CN110749855B (zh) | 一种基于协方差域零化的均匀线阵波达方向估计方法 | |
CN114726385A (zh) | 基于功率估计的卫星导航接收机空域抗干扰方法 | |
Liu et al. | Reduced‐dimension MVDR beamformer based on sub‐array optimization | |
CN112399366A (zh) | 基于Hankel矩阵及WKNN方差提取的室内定位法 | |
CN112087235A (zh) | 基于伪逆感知字典的稀疏度自适应doa估计方法及系统 | |
CN114415106B (zh) | 基于改进lamp网络的互耦阵列doa估计方法 | |
CN114996653A (zh) | 一种基于原子范数最小化的二维鲁棒自适应波束形成方法 | |
CN114184999B (zh) | 一种互耦小孔径阵列的生成式模型处理方法 | |
CN115014313A (zh) | 一种基于并行多尺度的偏振光罗盘航向误差处理方法 | |
CN110967664B (zh) | 基于cold阵列增强四元数esprit的doa估计方法 | |
CN110333477B (zh) | 杂波背景下天线阵列的信号波达方向估计方法 | |
CN107135026A (zh) | 未知互耦存在时基于矩阵重构的稳健波束形成方法 | |
CN110873866A (zh) | 一种互耦条件下单基地mimo雷达目标角度估计方法 |
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 |