CN110208735A - 一种基于稀疏贝叶斯学习的相干信号doa估计方法 - Google Patents

一种基于稀疏贝叶斯学习的相干信号doa估计方法 Download PDF

Info

Publication number
CN110208735A
CN110208735A CN201910506316.0A CN201910506316A CN110208735A CN 110208735 A CN110208735 A CN 110208735A CN 201910506316 A CN201910506316 A CN 201910506316A CN 110208735 A CN110208735 A CN 110208735A
Authority
CN
China
Prior art keywords
signal
distribution
indicate
matrix
gamma
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
CN201910506316.0A
Other languages
English (en)
Other versions
CN110208735B (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.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical University
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 Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN201910506316.0A priority Critical patent/CN110208735B/zh
Publication of CN110208735A publication Critical patent/CN110208735A/zh
Application granted granted Critical
Publication of CN110208735B publication Critical patent/CN110208735B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Direction-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/02Direction-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/14Systems for determining direction or deviation from predetermined direction
    • G01S3/143Systems for determining direction or deviation from predetermined direction by vectorial combination of signals derived from differently oriented antennae
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02DCLIMATE 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/00Reducing energy consumption in communication networks
    • Y02D30/70Reducing energy consumption in communication networks in wireless communication networks

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Analysis (AREA)
  • Algebra (AREA)
  • Pure & Applied Mathematics (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明提供了一种基于稀疏贝叶斯学习的相干信号DOA估计方法,获取接收阵列的输出信号,网格化观测空间,构造超完备阵列流形,结合稀疏表示的思想,将DOA估计问题转化为稀疏信号重构问题,求解如下稀疏矩阵方程,建立稀疏贝叶斯概率模型,采用变分贝叶斯推断方法计算各隐变量的近似后验分布,计算入射信号的DOA估计值。本发明无需预先估计入射信号个数,且不涉及解相干操作,有效地实现相干信号的DOA估计,针对相干信号,所提DOA估计方法的角度分辨能力好,与现有的其它DOA估计方法相比,该方法具有更好的测向精度。

Description

一种基于稀疏贝叶斯学习的相干信号DOA估计方法
技术领域
本发明涉及信号处理技术领域,尤其是一种相干信号DOA的估计方法。
背景技术
波达方向(DOA,Direction-of-Arrival)估计是阵列信号处理领域中的核心任务之一,该技术对由空间中按照某种特定方式放置的传感器阵列接收到的来波信号进行处理,从而估计出目标信号的个数、来波方向等信号参数,被广泛应用于雷达、声纳、无线通信等领域。DOA估计问题与谐波恢复有着密切联系,而谐波恢复问题可以通过参数化技术来解决,特别是以多重信号分类(MUSIC,Multiple Signal Classification)算法为代表的具有高角度分辨性能的子空间类算法。在利用MUSIC算法进行DOA估计时,对阵列输出数据的协方差矩阵进行特征值分解,得到信号子空间和噪声子空间,利用二者的正交性构造空间方位谱函数,以角度搜索的方式得到DOA估计结果。然而,子空间算法需要预知目标信号个数,在快拍数较少和信噪比较低的情况下,算法的估计性能出现严重下降,并且不能用于处理相干信号(协方差矩阵出现秩亏损,无法正确分离信号子空间和噪声子空间)。若使用子空间类算法处理相干信号,需要进行解相干操作。
不同于子空间类算法,基于稀疏贝叶斯学习的DOA估计算法将DOA估计问题转换为稀疏信号重构问题,在贝叶斯估计框架下对入射信号的统计特性进行估计,进而获得DOA估计结果。基于稀疏贝叶斯学习的DOA估计算法不需要提前预知目标信号个数,且不涉及能够处理相干信号。在实际应用中,由于信号的多径传播特性,阵列接收到的信号不再是单一的不相关信号,因此,研究基于稀疏贝叶斯学习的高分辨、稳健的相干信号DOA估计具有重要的应用价值。
发明内容
为了克服现有技术的不足,本发明提供一种稀疏贝叶斯学习框架下的相干信号DOA估计方法。本发明可以解决现有子空间类DOA估计算法在处理相干信号时需要预知目标个数以及解相干操作这一问题。
本发明解决其技术问题所采用的技术方案是包括如下步骤:
步骤一:获取接收阵列的输出信号Y;
设定使用M个全指向性传感器组成接收阵列,且假定空间中有N个远场窄带相干信号,分别以角度θn入射到接收阵列上,其中,n=1,2,…,N,利用接收阵列对入射信号进行接收采样,阵列的输出信号为:
Y=[y(1),y(2),…y(L)] (1)
其中,y(t)(t=1,…,L)表示阵列在t时刻的输出信号,L表示快拍数;
步骤二:网格化观测空间,构造超完备阵列流形A;
将观测空间角度在范围[-90°,90°]内以1°的角度间隔均匀划分,得到角度网格点集Θ={θ1,…,θK},其中K为网格点总数,且K>>N;根据角度网格点集Θ,构造阵列流形:
A=[a(θ1),a(θ2),…,a(θk),…,a(θK)] (2)
其中,
表示在网格点θk上的导向矢量,简记为ak,dm为接收阵列中第m个传感器的位置坐标,m=1,…,M,λ为入射相干信号的波长,j为虚数单位;
步骤三:结合稀疏表示的思想,将DOA估计问题转化为稀疏信号重构问题,求解如下稀疏矩阵方程:
Y=AX+n (4)
其中,X表示K×L维的信号矩阵,n表示M×L维的加性高斯白噪声矩阵;由于X中仅有N个非零行向量,X为稀疏矩阵;
步骤四:建立稀疏贝叶斯概率模型;
首先,对阵列输出信号Y的每一列向量进行复高斯分布假设,则Y的似然函数表示为:
其中,Y·i,X·i分别表示矩阵Y,X的第i列向量,IM表示单位矩阵,β>0表示噪声精确(precision),噪声精确是指噪声方差的倒数,对β进行伽马先验分布假设,即:
p(β)=Gamma(β|c,d) (6)
其中,c,d为伽马分布的参数;
接下来,对信号矩阵X构造分层稀疏先验;
在第一层先验中,对X的每一列进行复高斯先验假设,则X的概率分布为:
其中,Z·i表示D×L维矩阵Z的第i列向量,W是一个K×D维的权矩阵,μ=[μ1,…,μK]T,Λ=diag(γ),γ=[γ1,…,γK]T,diag(·)表示生成对角矩阵运算;超参数γ含有预设网格点方向上入射信号的功率信息;
在第二层先验中,分别对超参数W,Z,μ,γ进行先验假设,假设W,Z的每一列以及μ均服从零均值复高斯分布,γ的每个元素独立同分布,均服从伽马分布,即:
其中,W·i表示W的第i列向量,γi表示γ的第i个元素,α=[α1,…,αD]T,a,b,δ为分布的参数,(·)-1表示矩阵求逆运算;
在第三层先验中,对超参数α的每一个元素进行伽马先验假设:
其中,g,h为伽马分布的参数;
设置参数初值a0=b0=c0=d0=g0=h0=10-6,δ0=10-3;记ξ={X,W,Z,μ,α,β,γ},称作隐变量集;
步骤五:采用变分贝叶斯推断方法计算各隐变量ξi的近似后验分布q(ξi),得到:
即X中第i列向量X·i的后验分布为复高斯分布,其均值和方差ΣX分别为:
ΣX=[<β>AHA+<Ψ>-1]-1 (14)
其中,Ψ=Λ-1=(diag(γ))-1,Y·i,Z·i分别表示Y,Z的第i列向量,<·>表示求期望运算,(·)H表示矩阵的共轭转置运算;
即W中第j列向量的后验分布为复高斯分布,其均值和方差分别为:
其中,Xji表示X的第j行第i列元素,γjj分别表示γ,μ的第j个元素;
即Z中第i列向量Z·i的后验分布为复高斯分布,其均值和方差ΣZ分别为:
ΣZ=<I+WHΛW>-1 (18)
即μ的后验分布为复高斯分布,其均值μμ和方差Σμ分别为:
Σμ=<LΛ+δI>-1 (20)
即α中第i个元素αi的后验分布为伽马分布,该分布的参数g,h为:
g=g0+K (21)
其中,||·||2表示向量的2范数运算;
⑥q(β)=Gamma(β|c,d),即β的后验分布为伽马分布,该分布的参数c,d分别为:
c=c0+LM (23)
其中,tr(·)表示矩阵的迹;
即γ中第j个元素γj的后验分布为伽马分布,该分布的参数a,b分别为:
a=a0+L (25)
其中,W表示W的第j列向量,|·|表示绝对值;
根据{X,W,Z,μ,α,β,γ}的更新公式(13)~(26),在设置隐变量初值后,进行隐变量的迭代更新直至满足收敛条件,则停止迭代;隐变量初值X(0)=AH(AAH)-1Y,W(0)=1K×M,Z(0)=1M×L(0)=1K×1(0)=1M×1(0)=1K×1,其中(·)(r)表示第r步迭代中的变量,收敛条件为:||γ(r)(r-1)||2≤10-4
步骤六:计算入射信号的DOA估计值;
根据步骤五中得到的各隐变量估计值,以γ为纵坐标Θ为横坐标绘制空间方位谱图,在谱图中峰值点所在的网格点角度区域内,通过一维搜索的方式计算得到入射信号的DOA估计值计算公式如下:
其中,表示导向矢量aj对θj的偏导值,Re(·)表示取复数的实部,θi表示空间方位谱图中第i个峰值点所在的网格点角度区域。
由于超参数γ含有入射信号的功率信息,空间方位谱图中峰值点所对应的网格点即为入射信号DOA估计值;然而,当入射信号没有落在预设网格点上,即网格失配时,根据谱图峰值点所确定的DOA估计值存在较大的误差,需要进行进一步的精确计算。构造目标函数:
其中,表示中与γj无关的部分,表示中与γj有关的部分;将对γj求偏导并令结果为0,得到:
其中,对θj求偏导并令结果为0,即:
将式(29)代入式(30)中得到:
根据式(31)即可得到DOA估计值的计算公式(27)。
本发明的有益效果在于由于本发明无需预先估计入射信号个数,且不涉及解相干操作,有效地实现相干信号的DOA估计。本发明采用了多层先验模型,使用非零均值复高斯分布对目标信号的空域稀疏性和信号之间的相干性进行描述,有利于相干信号的稀疏重构。针对相干信号,所提DOA估计方法的角度分辨能力好,与现有的其它DOA估计方法相比,该方法具有更好的测向精度。本发明也能应用于独立信号和混合信号(相干+不相关)的DOA估计问题。
附图说明
图1是本发明DOA估计结果空间方位谱图。
图2是本发明与现有五种DOA估计方法在不同信噪比条件下对两组相干信号DOA估计结果的均方根误差对比图。
图3是本发明与现有五种DOA估计方法在不同快拍数条件下对两组相干信号的DOA估计结果的均方根误差对比图。
具体实施方式
下面结合附图和实施例对本发明进一步说明。
本发明的步骤为:
步骤一:获取接收阵列的输出信号Y;
设定使用M个全指向性传感器组成接收阵列,且假定空间中有N个远场窄带相干信号,分别以角度θn入射到接收阵列上,其中,n=1,2,…,N,利用接收阵列对入射信号进行接收采样,阵列的输出信号为:
Y=[y(1),y(2),…y(L)] (1)
其中,y(t)(t=1,…,L)表示阵列在t时刻的输出信号,L表示快拍数;
步骤二:网格化观测空间,构造超完备阵列流形A;
将观测空间角度在范围[-90°,90°]内以1°的角度间隔均匀划分,得到角度网格点集Θ={θ1,…,θK},其中K为网格点总数,且K>>N;根据角度网格点集Θ,构造阵列流形:
A=[a(θ1),a(θ2),…,a(θk),…,a(θK)] (2)
其中,
表示在网格点θk上的导向矢量,简记为ak,dm为接收阵列中第m个传感器的位置坐标,m=1,…,M,λ为入射相干信号的波长,j为虚数单位;
步骤三:结合稀疏表示的思想,将DOA估计问题转化为稀疏信号重构问题,求解如下稀疏矩阵方程:
Y=AX+n (4)
其中,X表示K×L维的信号矩阵,n表示M×L维的加性高斯白噪声矩阵;由于X中仅有N个非零行向量,X为稀疏矩阵;
步骤四:建立稀疏贝叶斯概率模型;
首先,对阵列输出信号Y的每一列向量进行复高斯分布假设,则Y的似然函数表示为:
其中,Y·i,X·i分别表示矩阵Y,X的第i列向量,IM表示单位矩阵,β>0表示噪声精确(precision),噪声精确是指噪声方差的倒数,对β进行伽马先验分布假设,即:
p(β)=Gamma(β|c,d) (6)
其中,c,d为伽马分布的参数;
接下来,对信号矩阵X构造分层稀疏先验;
在第一层先验中,对X的每一列进行复高斯先验假设,则X的概率分布为:
其中,Z·i表示D×L维矩阵Z的第i列向量,W是一个K×D维的权矩阵,μ=[μ1,…,μK]T,Λ=diag(γ),γ=[γ1,…,γK]T,diag(·)表示生成对角矩阵运算;超参数γ含有预设网格点方向上入射信号的功率信息;
在第二层先验中,分别对超参数W,Z,μ,γ进行先验假设,假设W,Z的每一列以及μ均服从零均值复高斯分布,γ的每个元素独立同分布,均服从伽马分布,即:
其中,W·i表示W的第i列向量,γi表示γ的第i个元素,α=[α1,…,αD]T,a,b,δ为分布的参数,(·)-1表示矩阵求逆运算;
在第三层先验中,对超参数α的每一个元素进行伽马先验假设:
其中,g,h为伽马分布的参数;
设置参数初值a0=b0=c0=d0=g0=h0=10-6,δ0=10-3;记ξ={X,W,Z,μ,α,β,γ},称作隐变量集;
步骤五:采用变分贝叶斯推断方法计算各隐变量ξi的近似后验分布q(ξi),得到:
即X中第i列向量X·i的后验分布为复高斯分布,其均值和方差ΣX分别为:
ΣX=[<β>AHA+<Ψ>-1]-1 (14)
其中,Ψ=Λ-1=(diag(γ))-1,Y·i,Z·i分别表示Y,Z的第i列向量,<·>表示求期望运算,(·)H表示矩阵的共轭转置运算;
即W中第j列向量的后验分布为复高斯分布,其均值和方差分别为:
其中,Xji表示X的第j行第i列元素,γjj分别表示γ,μ的第j个元素;
即Z中第i列向量Z·i的后验分布为复高斯分布,其均值和方差ΣZ分别为:
ΣZ=<I+WHΛW>-1 (18)
即μ的后验分布为复高斯分布,其均值μμ和方差Σμ分别为:
Σμ=<LΛ+δI>-1 (20)
即α中第i个元素αi的后验分布为伽马分布,该分布的参数g,h为:
g=g0+K (21)
其中,||·||2表示向量的2范数运算;
⑥q(β)=Gamma(β|c,d),即β的后验分布为伽马分布,该分布的参数c,d分别为:
c=c0+LM (23)
其中,tr(·)表示矩阵的迹;
即γ中第j个元素γj的后验分布为伽马分布,该分布的参数a,b分别为:
a=a0+L (25)
其中,W表示W的第j列向量,|·|表示绝对值;
根据{X,W,Z,μ,α,β,γ}的更新公式(13)~(26),在设置隐变量初值后,进行隐变量的迭代更新直至满足收敛条件,则停止迭代;隐变量初值X(0)=AH(AAH)-1Y,W(0)=1K×M,Z(0)=1M×L(0)=1K×1(0)=1M×1(0)=1K×1,其中(·)(r)表示第r步迭代中的变量,收敛条件为:||γ(r)(r-1)||2≤10-4
步骤六:计算入射信号的DOA估计值;
根据步骤五中得到的各隐变量估计值,以γ为纵坐标Θ为横坐标绘制空间方位谱图,在谱图中峰值点所在的网格点角度区域内,通过一维搜索的方式计算得到入射信号的DOA估计值计算公式如下:
其中,表示导向矢量aj对θj的偏导值,Re(·)表示取复数的实部,θi表示空间方位谱图中第i个峰值点所在的网格点角度区域。
由于超参数γ含有入射信号的功率信息,空间方位谱图中峰值点所对应的网格点即为入射信号DOA估计值;然而,当入射信号没有落在预设网格点上,即网格失配时,根据谱图峰值点所确定的DOA估计值存在较大的误差,需要进行进一步的精确计算。构造目标函数:
其中,表示中与γj无关的部分,表示中与γj有关的部分;将对γj求偏导并令结果为0,得到:
其中,对θj求偏导并令结果为0,即:
将式(29)代入式(30)中得到:
根据式(31)即可得到DOA估计值的计算公式(27)。
本发明的实施例的步骤如下:
步骤一:得到接收阵列的输出信号Y。
假定有N个远场相干窄带信号以角度入射到阵元数为M的接收阵列上,设阵列在t(t=1,2,…,L)时刻的输出信号为y(t),则阵列的输出信号Y=[y(1),y(2),…y(L)],其中,L为快拍数。
步骤二:网格化观测空间,构造超完备阵列流形A。
基于入射信号在空域角度方位内是有限的、稀疏分布的这一特性,将观测空间角度范围[-90°,90°]以某一角度间隔均匀划分成K(K>>N)个网格点,得到网格点集Θ={θ1,…,θK},构造超完备阵列流行A=[a(θ1),a(θ2),…,a(θK)],其中,表示在网格点θk(k=1,…,K)上的导向矢量,简记为ak,dm(m=1,…,M)为接收阵列中第m个传感器的位置坐标,λ为入射相干信号的波长,j为虚数单位,(·)T表示矩阵转置运算。
步骤三:结合稀疏表示的思想,将DOA估计问题转化为稀疏信号重构问题,即求解稀疏矩阵方程Y=AX+n,其中,X表示K×L维的信号矩阵,其为仅有N个非零行向量的稀疏矩阵,n表示M×L维的加性高斯白噪声矩阵。
步骤四:根据稀疏贝叶斯学习理论,建立稀疏贝叶斯概率模型。
指定阵列输出信号Y的每一列向量服从复高斯分布,引入超参数β>0表示高斯白噪声的精度,则Y的似然函数表示为其中,Y·i,X·i分别表示矩阵Y,X的第i列向量,I表示单位矩阵;对β进行伽马先验分布假设,即p(β)=Gamma(β|c,d)。
构建分层稀疏先验分布,在第一层先验中,引入超参数W,Z,μ,γ,W,Z分别为K×D维和D×L维的矩阵,μ=[μ1,…,μK]T,γ=[γ1,…,γK]T,指定信号矩阵X的每一列向量服从复高斯分布,其中,Z·i表示Z的第i列向量,Λ=diag(γ),diag(·)表示生成对角矩阵运算;在第二层先验中,对超参数W,Z,μ,γ进行先验假设,再引入超参数α=[α1,…,αD]T,指定W的每一行向量均服从零均值、方差为(diag(α))-1的复高斯分布,指定Z的每一列向量均服从零均值、方差为ID的复高斯分布,指定μ服从零均值、方差为δ-1IK的复高斯分布,指定γ的每一个元素均服从参数为a,b的伽马分布;在第三层先验中,对超参数α的每一个元素指定参数为g,h的伽马分布。
步骤五:根据稀疏贝叶斯学习理论,采用变分贝叶斯推断方法求解隐变量X,W,Z,μ,α,β,γ的后验分布,得到各自统计特性(均值、方差等)的更新公式;设置合适的隐变量初始值,迭代更新隐变量直至收敛,得到隐变量的最优估计值,各隐变量的更新公式如下:
①X中第i(i=1,…,L)列向量X·i的均值和方差ΣX
ΣX=[<β>AHA+<Ψ>-1]-1
其中,Ψ=Λ-1=(diag(γ))-1,Y·i,Z·i分别表示Y,Z的第i列向量,<·>表示求期望运算,(·)H表示矩阵的共轭转置运算,(·)-1表示矩阵求逆运算;
②W中第j(j=1,…,D)列向量的均值和方差为:
其中,Xji表示X的第j行第i列元素,γjj分别表示γ,μ的第j个元素;
③Z中第i(i=1,…,L)列向量Z·i的均值和方差ΣZ
ΣZ=<I+WHΛW>-1
④μ的均值μμ和方差Σμ为:
Σμ=<LΛ+δI>-1
⑤α中第i(i=1,…,D)个元素αi的后验分布为伽马分布,分布的参数g,hi
g=g0+K
其中,||·||2表示向量的2范数运算;
⑥β的后验分布为伽马分布,分布的参数c,d:
c=c0+LM
其中,tr(·)表示矩阵的迹;
⑦γ中第j(j=1,…,K)个元素γj的后验分布为伽马分布,分布的参数a,b:
a=a0+L
其中,W表示W的第j行向量,|·|表示绝对值。
步骤6:计算目标信号的DOA估计值。
以γ最优估计值取10倍的以10为底的对数值为纵坐标(单位为分贝dB),以网格点集Θ为横坐标,绘制空间方位谱图;由于γ是一个稀疏向量,其大多数元素趋于0,仅包含N个非零值,即谱图中的峰值,取第i(i=1,…,N)个峰值对应网格点所在的角度区域θi为搜索区间,计算该搜索区间内存在的目标信号DOA估计值:
其中,表示导向矢量aj对θj的偏导值,Re(·)表示取复数的实部。
本发明的效果可通过以下仿真说明:
1.仿真条件:
假设有2组共4个远场相干窄带信号(每组由2个信号组成,同一组信号相干,不同组信号不相关)入射到接收阵列上。观测空间角度为[-90°,90°],空间网格化划分间隔为1°。将所提方法记为CASBL,与MFOCUSS、SS-MUSIC、OGSBI、cRVM、L1-SVD和SS-MUSIC这五种方法进行性能比较,以蒙特卡罗曲线为参考,比较各方法在不同条件下的均方根误差(RMSE,average root mean-square error)曲线。
RMSE的计算公式表示为:
其中,Q表示实验次数,Q=200,为第q次实验中第n个入射信号的DOA估计值,为第n个入射信号的真实DOA值。
2.仿真内容与结果:
仿真1:假设这两组相干信号的入射角度分别为[-20.8°,-12.6°]和[12.3°,19.5°],相关系数分别为[-0.0349+j0.9994,-0.6490+j0.2622]T和[0.7092+j0.5541,0.7999+j0.0140]T。使用12元非均匀线阵进行DOA估计,阵列的阵元位置坐标为[-9,-6,-5,-3,-2,-1,1,2,3,5,6,9]λ/2,其中λ为入射窄带信号的波长。快拍数为100,信噪比(SNR,signal-to-noise ratio)为0dB。采用本发明进行DOA估计,得到的空间方位谱图如图1所示,图中横坐标为角度值,纵坐标为归一化幅度谱值(各元素值除以最大元素值)。
从图1可以看出,本发明能够正确分辨这两组相干信号,且谱图中的峰值较为尖锐,表明所提方法的角度分辨能力好;由于这两组信号之间是不相关的,仿真结果说明该方法不仅能够处理相干信号,也能处理独立信号和混合信号(相干+不相关)的DOA估计问题。
仿真2:假设两组相干信号的入射方向分别为随机变量在[-0.5°,-0.5°]上均匀分布。使用15元均匀线阵进行DOA估计,阵元间隔为入射窄带信号的半波长。采样快拍数为50,信噪比从-5dB增大到20dB,分别采用本发明和其它五种方法进行200次独立的DOA估计实验,计算不同信噪比条件下各方法估计结果的RMSE,得到均方根误差-信噪比曲线如图2所示。图2中横坐标为信噪比,纵坐标为DOA估计结果的RMSE。
从图2可以看出,CASBL在不同的信噪比条件下均具有最小的RMSE,且最贴近蒙特卡罗曲线,即最接近理想估计结果;特别是在低信噪比的情况下,本发明具有最佳的DOA估计性能。
仿真3:在仿真2的基础上,固定信噪比为0dB,将快拍数从30增大到120,分别采用本发明和其它五种方法进行200次独立的DOA估计实验,计算各方法在不同快拍数条件下RMSE,得到均方根误差-快拍数曲线如图3所示,图中横坐标为快拍数,纵坐标为DOA估计结果的RMSE。
从图3中可以看出,本发明在不同的快拍数条件下均具有小的RMSE,即具有最优的估计性能。

Claims (2)

1.一种基于稀疏贝叶斯学习的相干信号DOA估计方法,其特征在于包括下述步骤:
步骤一:获取接收阵列的输出信号Y;
设定使用M个全指向性传感器组成接收阵列,且假定空间中有N个远场窄带相干信号,分别以角度θn入射到接收阵列上,其中,n=1,2,…,N,利用接收阵列对入射信号进行接收采样,阵列的输出信号为:
Y=[y(1),y(2),…y(L)] (1)
其中,y(t)(t=1,…,L)表示阵列在t时刻的输出信号,L表示快拍数;
步骤二:网格化观测空间,构造超完备阵列流形A;
将观测空间角度在范围[-90°,90°]内以1°的角度间隔均匀划分,得到角度网格点集Θ={θ1,…,θK},其中K为网格点总数,且K>>N;根据角度网格点集Θ,构造阵列流形:
A=[a(θ1),a(θ2),…,a(θk),…,a(θK)] (2)
其中,
表示在网格点θk上的导向矢量,简记为ak,dm为接收阵列中第m个传感器的位置坐标,m=1,…,M,λ为入射相干信号的波长,j为虚数单位;
步骤三:结合稀疏表示的思想,将DOA估计问题转化为稀疏信号重构问题,求解如下稀疏矩阵方程:
Y=AX+n (4)
其中,X表示K×L维的信号矩阵,n表示M×L维的加性高斯白噪声矩阵;由于X中仅有N个非零行向量,X为稀疏矩阵;
步骤四:建立稀疏贝叶斯概率模型;
首先,对阵列输出信号Y的每一列向量进行复高斯分布假设,则Y的似然函数表示为:
其中,Y·i,X·i分别表示矩阵Y,X的第i列向量,IM表示单位矩阵,β>0表示噪声精确(precision),噪声精确是指噪声方差的倒数,对β进行伽马先验分布假设,即:
p(β)=Gamma(β|c,d) (6)
其中,c,d为伽马分布的参数;
接下来,对信号矩阵X构造分层稀疏先验;
在第一层先验中,对X的每一列进行复高斯先验假设,则X的概率分布为:
其中,Z·i表示D×L维矩阵Z的第i列向量,W是一个K×D维的权矩阵,μ=[μ1,…,μK]T,Λ=diag(γ),γ=[γ1,…,γK]T,diag(·)表示生成对角矩阵运算;超参数γ含有预设网格点方向上入射信号的功率信息;
在第二层先验中,分别对超参数W,Z,μ,γ进行先验假设,假设W,Z的每一列以及μ均服从零均值复高斯分布,γ的每个元素独立同分布,均服从伽马分布,即:
其中,W·i表示W的第i列向量,γi表示γ的第i个元素,α=[α1,…,αD]T,a,b,δ为分布的参数,(·)-1表示矩阵求逆运算;
在第三层先验中,对超参数α的每一个元素进行伽马先验假设:
其中,g,h为伽马分布的参数;
设置参数初值a0=b0=c0=d0=g0=h0=10-6,δ0=10-3;记ξ={X,W,Z,μ,α,β,γ},称作隐变量集;
步骤五:采用变分贝叶斯推断方法计算各隐变量ξi的近似后验分布q(ξi),得到:
即X中第i列向量X·i的后验分布为复高斯分布,其均值和方差ΣX分别为:
x=[〈β>AHA+<Ψ>-1]-1 (14)
其中,Ψ=Λ-1=(diag(γ))-1,Y·i,Z·i分别表示Y,Z的第i列向量,<·>表示求期望运算,(·)H表示矩阵的共轭转置运算;
即W中第j列向量的后验分布为复高斯分布,其均值和方差分别为:
其中,Xji表示X的第j行第i列元素,γjj分别表示γ,μ的第j个元素;
即Z中第i列向量Z·i的后验分布为复高斯分布,其均值和方差ΣZ分别为:
ΣZ=<I+WHΛW>-1 (18)
即μ的后验分布为复高斯分布,其均值μμ和方差Σμ分别为:
Σμ=<LΛ+δI>-1 (20)
即α中第i个元素αi的后验分布为伽马分布,该分布的参数g,h为:
g=g0+K (21)
其中,||·||2表示向量的2范数运算;
⑥q(β)=Gamma(β|c,d),即β的后验分布为伽马分布,该分布的参数c,d分别为:
c=c0+LM (23)
其中,tr(·)表示矩阵的迹;
即γ中第j个元素γj的后验分布为伽马分布,该分布的参数a,b分别为:
a=a0+L (25)
其中,Wj.表示W的第j列向量,|·|表示绝对值;
根据{X,W,Z,μ,α,β,γ}的更新公式(13)~(26),在设置隐变量初值后,进行隐变量的迭代更新直至满足收敛条件,则停止迭代;隐变量初值X(0)=AH(AAH)-1Y,W(0)=1K×M,Z(0)=1M×L(0)=1K×1(0)=1M×1(0)=1K×1,其中(·)(r)表示第r步迭代中的变量,收敛条件为:||γ(r)(r-1)||2≤10-4
步骤六:计算入射信号的DOA估计值;
根据步骤五中得到的各隐变量估计值,以γ为纵坐标Θ为横坐标绘制空间方位谱图,在谱图中峰值点所在的网格点角度区域内,通过一维搜索的方式计算得到入射信号的DOA估计值计算公式如下:
其中, 表示导向矢量aj对θj的偏导值,Re(·)表示取复数的实部,θi表示空间方位谱图中第i个峰值点所在的网格点角度区域。
2.根据权利要求1所述的一种基于稀疏贝叶斯学习的相干信号DOA估计方法,其特征在于:
由于超参数γ含有入射信号的功率信息,空间方位谱图中峰值点所对应的网格点即为入射信号DOA估计值;然而,当入射信号没有落在预设网格点上,即网格失配时,根据谱图峰值点所确定的DOA估计值存在较大的误差,需要进行进一步的精确计算;构造目标函数:
其中,表示中与γj无关的部分,表示中与γj有关的部分;将对γj求偏导并令结果为0,得到:
其中,对θj求偏导并令结果为0,即:
将式(29)代入式(30)中得到:
根据式(31)即可得到DOA估计值的计算公式(27)。
CN201910506316.0A 2019-06-12 2019-06-12 一种基于稀疏贝叶斯学习的相干信号doa估计方法 Active CN110208735B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910506316.0A CN110208735B (zh) 2019-06-12 2019-06-12 一种基于稀疏贝叶斯学习的相干信号doa估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910506316.0A CN110208735B (zh) 2019-06-12 2019-06-12 一种基于稀疏贝叶斯学习的相干信号doa估计方法

Publications (2)

Publication Number Publication Date
CN110208735A true CN110208735A (zh) 2019-09-06
CN110208735B CN110208735B (zh) 2022-11-11

Family

ID=67792241

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910506316.0A Active CN110208735B (zh) 2019-06-12 2019-06-12 一种基于稀疏贝叶斯学习的相干信号doa估计方法

Country Status (1)

Country Link
CN (1) CN110208735B (zh)

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111175692A (zh) * 2020-01-17 2020-05-19 西北工业大学 基于分层合成Lasso先验模型的离格稀疏贝叶斯DOA估计方法
CN111257845A (zh) * 2020-02-11 2020-06-09 中国人民解放军国防科技大学 一种基于近似消息传递的不在网格目标角度估计方法
CN111273301A (zh) * 2020-02-18 2020-06-12 西北工业大学 水声目标辐射噪声线阵波束输出信号频谱重构方法
CN111610512A (zh) * 2020-06-01 2020-09-01 桂林电子科技大学 一种基于稀疏贝叶斯学习的频控阵雷达离网目标定位方法
CN112731273A (zh) * 2020-12-09 2021-04-30 南京邮电大学 一种基于稀疏贝叶斯的低复杂度信号波达方向估计方法
CN112948606A (zh) * 2020-12-14 2021-06-11 西南交通大学 一种基于自适应网格的信号估计方法及装置
CN113406570A (zh) * 2021-06-05 2021-09-17 西北工业大学 一种平稳干扰环境下的贝叶斯稳健波束形成方法
CN113406571A (zh) * 2021-06-05 2021-09-17 西北工业大学 一种运动干扰环境下的贝叶斯稳健波束形成方法
CN113673158A (zh) * 2021-08-19 2021-11-19 西北工业大学 适用于强干扰环境下的波束域变分贝叶斯方位估计方法
CN113985348A (zh) * 2021-10-25 2022-01-28 合肥工业大学 基于多任务学习的单快拍相干超分辨doa估计技术
CN114157538A (zh) * 2021-11-22 2022-03-08 清华大学 一种基于双通道接收机的无线信号到达角估计方法及系统
CN114415105A (zh) * 2021-12-31 2022-04-29 西北工业大学 一种阵列互耦情况下波达方向估计方法
CN114415109A (zh) * 2022-01-10 2022-04-29 西北工业大学 一种稀疏贝叶斯学习的直接定位方法
CN115407008A (zh) * 2021-05-26 2022-11-29 株式会社岛津制作所 分析方法和诊断辅助方法
CN114415105B (zh) * 2021-12-31 2024-05-24 西北工业大学 一种阵列互耦情况下波达方向估计方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR3046250A1 (fr) * 2015-12-23 2017-06-30 Thales Sa Procede de determination de la direction d'arrivee en presence de repliement spectral et dispositif associe
CN109116293A (zh) * 2018-08-22 2019-01-01 上海师范大学 一种基于离格稀疏贝叶斯的波达方向估计方法
CN109298382A (zh) * 2018-09-10 2019-02-01 西北工业大学 一种基于期望极大算法的非均匀直线阵波达方向角估计方法
CN109444810A (zh) * 2018-12-24 2019-03-08 哈尔滨工程大学 一种非负稀疏贝叶斯学习框架下的互质阵列非网格doa估计方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR3046250A1 (fr) * 2015-12-23 2017-06-30 Thales Sa Procede de determination de la direction d'arrivee en presence de repliement spectral et dispositif associe
CN109116293A (zh) * 2018-08-22 2019-01-01 上海师范大学 一种基于离格稀疏贝叶斯的波达方向估计方法
CN109298382A (zh) * 2018-09-10 2019-02-01 西北工业大学 一种基于期望极大算法的非均匀直线阵波达方向角估计方法
CN109444810A (zh) * 2018-12-24 2019-03-08 哈尔滨工程大学 一种非负稀疏贝叶斯学习框架下的互质阵列非网格doa估计方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
JIE YANG: "An Efficient Compressed Sensing-based DOA Estimation Method in Nested MIMO Sonar", 《OCEANS 2017 - ABERDEEN》 *
冯明月: "基于Bessel 先验快速稀疏贝叶斯学习的互质阵列DOA估计", 《电子与信息学报》 *

Cited By (25)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111175692A (zh) * 2020-01-17 2020-05-19 西北工业大学 基于分层合成Lasso先验模型的离格稀疏贝叶斯DOA估计方法
CN111175692B (zh) * 2020-01-17 2022-09-02 西北工业大学 基于分层合成Lasso先验模型的离格稀疏贝叶斯DOA估计方法
CN111257845A (zh) * 2020-02-11 2020-06-09 中国人民解放军国防科技大学 一种基于近似消息传递的不在网格目标角度估计方法
CN111257845B (zh) * 2020-02-11 2020-09-22 中国人民解放军国防科技大学 一种基于近似消息传递的不在网格目标角度估计方法
CN111273301A (zh) * 2020-02-18 2020-06-12 西北工业大学 水声目标辐射噪声线阵波束输出信号频谱重构方法
CN111610512B (zh) * 2020-06-01 2022-08-09 桂林电子科技大学 一种基于稀疏贝叶斯学习的频控阵雷达离网目标定位方法
CN111610512A (zh) * 2020-06-01 2020-09-01 桂林电子科技大学 一种基于稀疏贝叶斯学习的频控阵雷达离网目标定位方法
CN112731273A (zh) * 2020-12-09 2021-04-30 南京邮电大学 一种基于稀疏贝叶斯的低复杂度信号波达方向估计方法
CN112731273B (zh) * 2020-12-09 2023-06-23 南京邮电大学 一种基于稀疏贝叶斯的低复杂度信号波达方向估计方法
CN112948606A (zh) * 2020-12-14 2021-06-11 西南交通大学 一种基于自适应网格的信号估计方法及装置
CN112948606B (zh) * 2020-12-14 2022-10-21 西南交通大学 一种基于自适应网格的信号估计方法及装置
CN115407008A (zh) * 2021-05-26 2022-11-29 株式会社岛津制作所 分析方法和诊断辅助方法
CN113406571A (zh) * 2021-06-05 2021-09-17 西北工业大学 一种运动干扰环境下的贝叶斯稳健波束形成方法
CN113406571B (zh) * 2021-06-05 2024-04-09 西北工业大学 一种运动干扰环境下的贝叶斯稳健波束形成方法
CN113406570B (zh) * 2021-06-05 2024-04-12 西北工业大学 一种平稳干扰环境下的贝叶斯稳健波束形成方法
CN113406570A (zh) * 2021-06-05 2021-09-17 西北工业大学 一种平稳干扰环境下的贝叶斯稳健波束形成方法
CN113673158A (zh) * 2021-08-19 2021-11-19 西北工业大学 适用于强干扰环境下的波束域变分贝叶斯方位估计方法
CN113985348B (zh) * 2021-10-25 2024-05-07 合肥工业大学 基于多任务学习的单快拍相干超分辨doa估计技术
CN113985348A (zh) * 2021-10-25 2022-01-28 合肥工业大学 基于多任务学习的单快拍相干超分辨doa估计技术
CN114157538B (zh) * 2021-11-22 2023-06-06 清华大学 一种基于双通道接收机的无线信号到达角估计方法及系统
CN114157538A (zh) * 2021-11-22 2022-03-08 清华大学 一种基于双通道接收机的无线信号到达角估计方法及系统
CN114415105A (zh) * 2021-12-31 2022-04-29 西北工业大学 一种阵列互耦情况下波达方向估计方法
CN114415105B (zh) * 2021-12-31 2024-05-24 西北工业大学 一种阵列互耦情况下波达方向估计方法
CN114415109A (zh) * 2022-01-10 2022-04-29 西北工业大学 一种稀疏贝叶斯学习的直接定位方法
CN114415109B (zh) * 2022-01-10 2024-04-26 西北工业大学 一种稀疏贝叶斯学习的直接定位方法

Also Published As

Publication number Publication date
CN110208735B (zh) 2022-11-11

Similar Documents

Publication Publication Date Title
CN110208735A (zh) 一种基于稀疏贝叶斯学习的相干信号doa估计方法
CN106772226B (zh) 基于压缩感知时间调制阵列的doa估计方法
Liu et al. An efficient maximum likelihood method for direction-of-arrival estimation via sparse Bayesian learning
CN109298383B (zh) 一种基于变分贝叶斯推断的互质阵波达方向角估计方法
CN104749553B (zh) 基于快速稀疏贝叶斯学习的波达方向角估计方法
CN110109051B (zh) 基于频控阵的互耦阵列doa估计方法
WO2018094565A1 (zh) 脉冲噪声下的波束成形方法及装置
CN108375763B (zh) 一种应用于多声源环境的分频定位方法
CN109450499B (zh) 一种基于导向矢量和空间功率估计的鲁棒波束形成方法
CN110045321B (zh) 基于稀疏和低秩恢复的稳健doa估计方法
CN105259550B (zh) 基于压缩感知的多输入多输出雷达二维角度估计方法
CN106646344A (zh) 一种利用互质阵的波达方向估计方法
CN109407046A (zh) 一种基于变分贝叶斯推断的嵌套阵列波达方向角估计方法
CN105306123A (zh) 一种抗阵列系统误差的稳健波束形成方法
CN107576931A (zh) 一种基于协方差低维度迭代稀疏重构的相关/相干信号波达方向估计方法
CN109600152A (zh) 一种基于子空间基变换的自适应波束形成方法
CN111580042B (zh) 一种基于相位优化的深度学习测向方法
CN109298382A (zh) 一种基于期望极大算法的非均匀直线阵波达方向角估计方法
CN112255629A (zh) 基于联合uca阵列的序贯esprit二维不相干分布源参数估计方法
CN106680779B (zh) 脉冲噪声下的波束成形方法及装置
CN106788655A (zh) 互耦条件下未知互耦信息的干扰相干稳健波束形成方法
CN104215939A (zh) 一种融合广义对称结构信息的知识辅助空时自适应处理方法
CN111880143B (zh) 改进稀疏贝叶斯学习的高精度定位方法、存储介质及设备
CN113671439A (zh) 基于非均匀智能超表面阵列的无人机集群测向系统及方法
CN113805139A (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