CN109407045B - 一种非均匀传感器阵列宽带信号波达方向估计方法 - Google Patents

一种非均匀传感器阵列宽带信号波达方向估计方法 Download PDF

Info

Publication number
CN109407045B
CN109407045B CN201811178910.3A CN201811178910A CN109407045B CN 109407045 B CN109407045 B CN 109407045B CN 201811178910 A CN201811178910 A CN 201811178910A CN 109407045 B CN109407045 B CN 109407045B
Authority
CN
China
Prior art keywords
vector
gamma
array
sparse representation
iteration
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
CN201811178910.3A
Other languages
English (en)
Other versions
CN109407045A (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.)
Suzhou University
Original Assignee
Suzhou 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 Suzhou University filed Critical Suzhou University
Priority to CN201811178910.3A priority Critical patent/CN109407045B/zh
Publication of CN109407045A publication Critical patent/CN109407045A/zh
Application granted granted Critical
Publication of CN109407045B publication Critical patent/CN109407045B/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
    • 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/802—Systems for determining direction or deviation from predetermined direction

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)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种非均匀传感器阵列宽带信号波达方向估计方法，包括以下步骤：S1、对非均匀传感器阵列中的样本序列作滑窗快速傅里叶变换，并构造实值加权样本协方差矢量yl；S2、构造对应的过完备阵列流型矩阵Φl；S3、利用实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl设置迭代程序中空间联合稀疏表示向量γ的初始值γinit和正则化参数λ的初始值λinit；S4、利用实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl更新空间联合稀疏表示向量γ和正则化参数λ；S5、当第i次迭代得到的γi与第i+1次迭代得到的γi+1之间满足||γi‑γi+1||2/||γi||2＜α或者达到了最大迭代次数Pmax，终止迭代；否则返回S4；S6、空间联合稀疏表示向量γ的所有峰值，并输出宽带信号的波达方向估计值。

Description

一种非均匀传感器阵列宽带信号波达方向估计方法
技术领域
本发明涉及阵列信号处理领域，更具体的说涉及一种非均匀传感器阵列宽带信号波达方向估计方法。
背景技术
宽带信号具有抗干扰能力强、调制样式复杂、携带信息量大等特点，在雷达、通信、电子对抗侦察等领域具有十分广泛的应用，同时自然界中可接收的自然信号绝大多数都是宽带信号，因此对宽带信号的波达方向(Direction-of-Arrival，DOA)进行估计的问题正日益受到关注。
近年来，日益复杂多变的电磁环境导致来波信号数量增多，来波信号个数大于阵元数的欠定情况时有发生，而传统的传感器位置均匀排列的均匀阵列无法解决欠定DOA估计问题。对于欠定DOA估计问题，利用嵌套阵列或互质阵列等非均匀阵列进行DOA估计的研究正受到越来越多的关注。
利用非均匀阵列进行宽带信号的DOA估计，目前已有的研究都是先做多个窄频带划分，然后基于差分阵列的概念来构造虚拟扩展孔径阵列，而最终的DOA估计主要从两个不同的角度出发：1、合成各频带上虚拟阵列下的增广协方差矩阵，然后利用宽带信号DOA估计中经典的非相干子空间方法(Incoherent Signal Subspace Method,ISSM)实现DOA估计；2、基于稀疏表示思想，将该问题转化为多个频带上共享相同稀疏支撑集的非负稀疏信号重构问题并进行求解。
文献N.Hu,et al.Underdetermined DOA estimation method for widebandsignals using joint nonnegative sparse Bayesian learning,IEEE SPL,2017中公开了一种联合非负稀疏贝叶斯学习方法，该文献中利用各个频带上数据的联合非负稀疏性，建立了层级贝叶斯模型。利用期望最大化算法实现这种联合非负稀疏贝叶斯学习可以在迭代中给出波达方向的估值；该方法算法复杂度低、收敛速度快。然而在其层级贝叶斯模型中并未对基于凸优化的方法中常用于多频带联合稀疏约束的
Figure GDA0002450397090000011
范数进行建模，因此其在低信噪比或者采样点数较少的情况下性能会受到影响。
发明内容
本发明要解决的技术问题是提供一种非均匀传感器阵列宽带信号波达方向估计方法，其能够提高波达方向的估计精度。
为了解决上述技术问题，本发明提供了一种非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，包括以下步骤：
S1、对非均匀传感器阵列中各传感器接收通道采集的样本序列作滑窗快速傅里叶变换，并构造L个窄频带上的实值加权样本协方差矢量yl；
S2、构造与各窄频带上的实值加权样本协方差矢量yl相对应的过完备阵列流型矩阵Φl；
S3、利用L个窄频带上的实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl设置迭代程序中空间联合稀疏表示向量γ的初始值γinit，并根据实值加权样本协方差矢量yl、过完备阵列流型矩阵Φl和空间联合稀疏表示向量的初始值γinit计算正则化参数λ的初始值λinit；
S4、利用L个窄频带上的实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl更新空间联合稀疏表示向量γ和正则化参数λ；
S5、当第i次迭代得到的γi与第i+1次迭代得到的γi+1之间满足||γi-γi+1||2/||γi||2＜α或者达到了最大迭代次数Pmax，终止迭代；否则继续返回S4迭代；其中α为一个任意趋于0的小数；
S6、当满足迭代终止条件时，搜索当前更新所得的空间联合稀疏表示向量γ的所有峰值，并根据空间联合稀疏表示向量γ的峰值位置输出宽带信号的波达方向估计值。
进一步的，步骤S1包括：
S11、对非均匀传感器阵列中各传感器接收通道采集的样本序列作滑窗快速傅里叶变换，得到输出快拍矢量序列xl(t),l＝1,2,…,L；t＝1,2,…,T；
S12、利用输出快拍矢量序列构造实值加权样本协方差矢量yl,l＝1,2,…,L。
进一步的，步骤S4包括：
S41、利用实值加权样本协方差矢量yl、过完备阵列流型矩阵Φl以及上一轮迭代更新的空间联合稀疏表示向量γ，更新隐变量向量wl所服从的非负高斯分布的均值向量μl与协方差矩阵Σl，其中，若当前为第一轮迭代则空间联合稀疏表示向量γ＝γinit；
S42、利用隐变量向量wl所服从的非负高斯分布的均值μl、协方差矩阵Σl和上一轮迭代更新的正则化参数λ，更新空间联合稀疏表示向量γ的值；其中，若当前为第一轮迭代则正则化参数λ＝λinit；利用迭代更新的空间联合稀疏表示向量γ，更新正则化参数λ的值。
进一步的，步骤S11包括：
S111、利用窗函数对非均匀传感器阵列中各传感器接收通道采集的样本序列进行滑动数据选取，其中阵元个数为M；
S112、对各窗口中的加权非均匀阵列采样数据做快速傅里叶变换，其中，所加窗函数的窗长以及快速傅里叶变换点数为NFFT，相邻窗交迭Noverlap点；
S113、根据目标宽带信号所覆盖的频域范围，在第t个快速傅里叶变换窗口所得到的NFFT个频域值中选择L个作为待使用频带上的值，由于有M个传感器接收通道，可构造L个子频带上的每个子频带含有M个元素的向量：xl(t),l＝1,2,…,L；t＝1,2,…,T，其中每个向量是该频带上的“频域快拍”；
S114、对非均匀阵列采样数据共做T个滑窗快速傅里叶变换，可构造L个子频带上每个频带上T个“频域快拍”的数据向量：xl(t),l＝1,2,…,L；t＝1,2,…,T。
进一步的，步骤S12包括：
S121、计算L个窄频带上样本协方差矩阵
Figure GDA0002450397090000031
其中(·)H表示共轭转置；
S122、对矩阵
Figure GDA0002450397090000032
作向量化得到样本协方差矢量
Figure GDA0002450397090000033
S123、用加权矩阵J对
Figure GDA0002450397090000034
进行加权处理以去除噪声项，得到加权样本协方差矢量
Figure GDA0002450397090000035
其中权重矩阵J＝[J1,J2,…JM-1]Τ，Jm＝[em(M+1)-M+1,em(M+1)-M+2,…,em(M+1)]，含M2个元素的向量em除第m个元素为1之外其他所有元素都为0，其中(·)T表示矩阵转置；
S124、利用J与
Figure GDA0002450397090000036
构造加权样本协方差矢量的协方差矩阵
Figure GDA0002450397090000037
S125、构造实值加权样本协方差矢量
Figure GDA0002450397090000038
其中
Figure GDA0002450397090000039
Re(·)与Im(·)分别表示取实部与取虚部。
进一步的，步骤S2包括：
S21、以一定的角度网格间距构造可以覆盖阵列探测角度范围的密集角度抽样网格
Figure GDA00024503970900000310
S22、对于角度
Figure GDA00024503970900000311
构造L个窄频带上与实值加权样本协方差矢量相对应的阵列导向矢量
Figure GDA00024503970900000312
其中
Figure GDA00024503970900000313
Figure GDA00024503970900000314
为Kronecker积，原始阵列导向矢量
Figure GDA00024503970900000315
fl为第l个频带的频率，dm为第m个传感器相对于参考传感器的间距，v为波速，(·)*表示求共轭；
S23、对于网格
Figure GDA0002450397090000041
中的每一个角度
Figure GDA0002450397090000042
计算
Figure GDA0002450397090000043
并由此构造L个窄频带上与各实值加权样本协方差矢量相对应的过完备阵列流型矩阵Φl：
Figure GDA0002450397090000044
进一步的，步骤S3包括：
S31、计算空间联合稀疏表示向量γ的初值
Figure GDA0002450397090000045
其中⊙为Hadamard积，||·||F为Frobenius范数；
S32、计算
Figure GDA0002450397090000046
其中
Figure GDA0002450397090000047
进而计算正则化参数λ的初值
Figure GDA0002450397090000048
其中||·||2为
Figure GDA00024503970900000418
范数；
S33、设置正则化参数λ的无信息超先验分布参数固定参数值c,d，其中c,d→0。
进一步的，步骤S41包括：
S421、计算L个窄频带上的隐变量向量wl中各元素的二阶矩：
Figure GDA0002450397090000049
其中
Figure GDA00024503970900000410
特别当
Figure GDA00024503970900000411
时，直接令
Figure GDA00024503970900000412
S422、利用计算所得的
Figure GDA00024503970900000417
以及上一轮迭代所更新的正则化参数λ，更新空间联合稀疏表示向量γ的第n个元素
Figure GDA00024503970900000414
S423、利用上一轮迭代所更新的空间联合稀疏表示向量γ，更新正则化参数λ：
Figure GDA00024503970900000415
进一步的，步骤S6包括，当满足迭代终止条件时，搜索当前更新所得的空间联合稀疏表示向量γ的所有峰值，其峰值位置与网格
Figure GDA00024503970900000416
中对应的角度即为宽带信号的波达方向估计值。
本发明的有益效果：
1、本发明将非均匀传感器阵列的宽带信号波达方向估计问题转化为多个窄频带上的具有相同非负稀疏支撑集的多向量联合非负稀疏重构问题，基于稀疏贝叶斯学习以及贝叶斯生成网络理论，建立层级贝叶斯模型并以此为基础进行求解，可实现欠定情况下的宽带信号波达方向估计；
2、本发明基于目标求解向量的非负性，将非均匀传感器阵列的宽带信号波达方向估计所对应的复数问题转化为实数问题，降低了算法的复杂度；
3、可确保欠定条件下信号源个数较多时的DOA估计精度；
4、可确保信噪比低、采样点数少的情况下DOA估计精度；
5、本发明所提出的非均匀传感器阵列的宽带信号波达方向估计方法，不需要预设任何经验参数，且其迭代算法的初值完全由传感器阵列输出数据得来，从而提高了DOA估计精度。
附图说明
图1是非均匀传感器阵列宽带信号波达方向估计方法的流程图；
图2是在SNR＝-16dB，频域快拍数T＝100条件下的空间谱图；
图3是在频域快拍数T＝100条件下的RMSE随着SNR变化的趋势图；
图4是在SNR＝-16dB条件下的RMSE随着频域快拍数变化的趋势图。
具体实施方式
下面结合附图和具体实施例对本发明作进一步说明，以使本领域的技术人员可以更好地理解本发明并能予以实施，但所举实施例不作为对本发明的限定。
实施例
一种非均匀传感器阵列宽带信号波达方向估计方法，参照图1，仿真产生中心频率为100Hz，带宽40Hz宽带声信号源(声信号的传播速度为340m/s)。采用二级嵌套声传感器阵列作为非均匀阵列的例子，包括6个传感器，分别位于[0,d,2d,3d,7d,11d]，其中d＝1.7m。共生成9个覆盖相同频域范围且不相关的宽带信号源，从远场-69°、-45°、-30°、-15°、0°、12°、24°、42°、65°这9个不同的方向入射到该非均匀阵列上，各传感器接收信号时叠加与信号不相关的高斯白噪声，最后阵列输出采样率为240Hz的数据。
对各传感器通道输出数据样本序列作滑窗FFT以得到L＝8个窄频带上的阵列输出快拍矢量序列：对一小段非均匀阵列采样数据做滑窗FFT，所加汉明窗的窗长以及FFT点数为NFFT＝64，相邻窗交迭Noverlap＝32点；根据目标宽带信号所覆盖的频域范围80～120Hz，在每个FFT窗口所得到的NFFT＝64个频域值中选择86.25Hz、90Hz、93.75Hz、97.5Hz、101.25Hz、105Hz、108.75Hz和112.5Hz共L＝8个频带作为目标子频带数据以供后续使用，对该段数据共做T个滑窗FFT，从而构造L＝8个子频带上每个子频带T个“频域快拍”的数据向量：xl(t),l＝1,2,…,L；t＝1,2,…,T。
利用L个窄频带上的阵列输出快拍矢量序列构造L个窄频带上的实值加权样本协方差矢量：首先，计算每个频点对应的样本协方差矩阵
Figure GDA0002450397090000061
其中(·)H表示共轭转置；然后，对
Figure GDA0002450397090000062
作向量化得到样本协方差矢量
Figure GDA0002450397090000063
其次，用加权矩阵J对
Figure GDA0002450397090000064
进行加权处理以去除噪声项，得到加权样本协方差矢量
Figure GDA0002450397090000065
其中权重矩阵J＝[J1,J2,…JM-1]Τ，Jm＝[em(M+1)-M+1,em(M+1)-M+2,…,em(M+1)]，含M2个元素的向量em除第m个元素为1之外其他所有元素都为0，其中(·)T表示矩阵转置；利用J与
Figure GDA0002450397090000066
构造加权样本协方差矢量的协方差矩阵
Figure GDA0002450397090000067
最后，构造实值加权样本协方差矢量
Figure GDA0002450397090000068
其中
Figure GDA0002450397090000069
Re(·)与Im(·)分别表示取实部与取虚部。
构造L个窄频带上与各实值加权样本协方差矢量相对应的过完备阵列流型矩阵：首先，以1°的角度网格间距构造可以覆盖阵列探测角度范围的密集角度抽样网格
Figure GDA00024503970900000610
然后，对于角度
Figure GDA00024503970900000611
构造L个窄频带上与实值加权样本协方差矢量相对应的阵列导向矢量
Figure GDA00024503970900000612
其中
Figure GDA00024503970900000613
Figure GDA00024503970900000614
为Kronecker积，
Figure GDA00024503970900000615
中原始阵列导向矢量
Figure GDA00024503970900000616
fl为第l个频带的频率，dm为第m个传感器相对于参考传感器的间距，v为波速，(·)*表示求共轭；最后，对于网格
Figure GDA00024503970900000617
中的每一个角度
Figure GDA00024503970900000618
计算
Figure GDA00024503970900000619
并由此构造L个窄频带上与各实值加权样本协方差矢量相对应的过完备阵列流型矩阵Φl：
Figure GDA00024503970900000620
设置迭代程序的初值以及固定参数：首先，计算空间联合稀疏表示向量γ的初值
Figure GDA0002450397090000071
其中⊙为Hadamard积，||·||F为Frobenius范数；然后，计算
Figure GDA0002450397090000072
其中
Figure GDA0002450397090000073
进而计算空间联合稀疏表示向量λ的初值
Figure GDA0002450397090000074
其中||·||2为l2范数；设置固定参数值c＝d＝10-6。
更新L个窄频带上的隐变量向量wl所服从的非负高斯分布的均值与协方差矩阵：首先，利用上一轮迭代所得的γ，更新L个窄频带上的隐变量向量wl所服从的非负高斯分布的协方差矩阵
Figure GDA0002450397090000075
其中Γ＝diag(γ)；然后，更新L个窄频带上的隐变量向量wl所服从的非负高斯分布的均值
Figure GDA0002450397090000076
迭代更新空间联合稀疏表示向量γ与正则化参数λ的值：首先，计算L个窄频带上的隐变量向量wl中各元素的二阶矩：
Figure GDA0002450397090000077
其中
Figure GDA0002450397090000078
特别当
Figure GDA0002450397090000079
时，直接令
Figure GDA00024503970900000710
然后，利用计算所得的
Figure GDA00024503970900000715
以及上一轮迭代所更新的正则化参数λ，更新空间联合稀疏表示向量γ的第n个元素
Figure GDA00024503970900000712
最后，利用上一轮迭代所更新的空间联合稀疏表示向量γ，更新正则化参数λ：
Figure GDA00024503970900000713
计算是否满足本实施例设置的迭代终止条件(最大迭代次数Pmax达到50)，若不满足继续迭代；若满足，终止迭代，搜索当前更新所得的空间联合稀疏表示向量γ的所有峰值，其峰值位置与网格
Figure GDA00024503970900000714
中对应的角度即为宽带信号的波达方向估计值。
参照图2，在信噪比(signal-to-noise ratio,SNR)为-16dB，频域快拍数T＝100时本发明所提方法的空间谱图，其中各垂直虚线标出了真实DOA的位置，可见在欠定条件(信号源个数大于传感器个数)以及信噪比非常低的情况下该方法仍能通过搜索空间谱的峰值较为精确地对各信号进行定位，同时所提方法中所隐含的层级稀疏贝叶斯分布很好地抑制了空间谱伪峰的出现。
参照图3和图4，以信噪比和快拍数为变量，以均方根差(root-mean-suqare-error,RMSE)评估本发明所提宽带DOA估计方法的精度，图3为T＝100条件下RMSE随着SNR变化的趋势图，图4为SNR＝-5dB条件下RMSE随着频域快拍数变化的趋势图。在欠定条件下信号源个数较多时，以及信噪比低、采样点数少的情况下，本发明仍能达到接近于克拉美罗界(所有估计方法的误差都不会小于该下界)的估计效果。
此外，本发明所提方法所涉及到的每次迭代的计算复杂度约为O(L·min{M3(M-1)3,N3})，同时该方法能快速收敛，从而相比于其他稀疏重构宽带DOA估计方法在计算复杂度上要小得多。
本发明所涉及迭代计算的理论推导过程：
以估计一维波达方向的线型阵列为例，假设某非均匀分布线型传感器阵列共有M个全向传感器，分别位于0,d1,…,dM-1位置，其中dm-1为第m个传感器距参考传感器的距离。假设有K个不相关的远场宽带信号源，分别从K个方向θk,k＝1,2,…,K入射。利用FFT或者窄带滤波器将阵列输出采样数据划分为L个窄带，每个窄带上的“频域快拍”数为T
xl(t)＝Alsl(t)+nl(t),l＝1,2,…,L；t＝1,2,…,T, (1)
其中xl(t),sl(t),nl(t),分别表示第l个子频带上的阵列输出向量、源信号向量和加性噪声向量。Al为第l个子频带(fl频点)上的导向矢量矩阵：Al＝[a(fl,θ1),a(fl,θ2),…,a(fl,θK)]，其中阵列导向矢量
Figure GDA0002450397090000081
v为波速。假设信号源与噪声不相关，则xl(t)的协方差矩阵可以写作
Figure GDA0002450397090000082
其中E{·}表示期望，
Figure GDA0002450397090000083
是K个信号源在fl频点上的方差，Ql是fl频点上的对角的噪声协方差矩阵。
对Rl进行矢量化以构建协方差矢量，得到
Figure GDA0002450397090000084
式中vec(·)为矢量化运算，符号
Figure GDA0002450397090000085
表示Khatri-Rao积，
Figure GDA0002450397090000086
νl＝vec(Ql)。
Figure GDA0002450397090000087
的列数反映了阵列的自由度。由于νl的非零项仅存在于特定位置，加权矩阵J可用来消除噪声，并得到
Figure GDA0002450397090000088
其中J＝[J1,J2,…JM-1]Τ，Jm＝[em(M+1)-M+1,em(M+1)-M+2,…,em(M+1)]。
在实际应用中由于快拍数有限，我们只能得到样本协方差矩阵
Figure GDA0002450397090000091
当源信号服从圆对称高斯分布时，样本协方差矩阵分布满足渐进复高斯分布，即
Figure GDA0002450397090000092
这里
Figure GDA0002450397090000093
表示Kronecker积，利用加权矩阵J对
Figure GDA0002450397090000094
进行加权，得到
Figure GDA0002450397090000095
其中
Figure GDA0002450397090000096
又因为pl为非负实向量，(6)式可转换为如下的实值高斯分布形式
Figure GDA0002450397090000097
且
Figure GDA0002450397090000098
若以稀疏贝叶斯学习(SBL)思想进行DOA估计，设置密集角度抽样网格
Figure GDA0002450397090000099
覆盖整个入射方向的范围，(7)式中的条件分布借助抽样网格可转化为稀疏表示形式
Figure GDA00024503970900000910
其中Φl是
Figure GDA00024503970900000911
的过完备形式，列为
Figure GDA00024503970900000912
其中
Figure GDA00024503970900000913
wl是非负稀疏方差向量，它的非零项的位置与抽样网格
Figure GDA00024503970900000914
中实际的DOA相对应。当l≠m，yl|wl与ym|wm相互独立，故可进一步得
Figure GDA00024503970900000915
其中wl,l＝1,2,…,L具有相同的稀疏支撑集，因此假设这些向量独立同分布且满足n维非负高斯分布，其联合分布为
Figure GDA00024503970900000916
其中
Figure GDA00024503970900000917
表示非负高斯概率密度函数。
空间联合稀疏表示向量γ决定了与向量wl,l＝1,2,…,L相同的稀疏支撑集，这里我们认为空间联合稀疏表示向量γ是一个随机向量，并假定所有的γ[n]独立同分布且满足Gamma先验分布
Figure GDA0002450397090000101
将(11)式与(10)式结合，并利用广义逆高斯(GIG)概率密度函数的可积性和第二类修正贝塞尔函数
Figure GDA0002450397090000102
可以得出：
Figure GDA0002450397090000103
由(12)式可知，该模型进行稀疏约束的方式与广泛应用的BPDN方法中的
Figure GDA0002450397090000108
范数的作用相似。
最后假设λ具有伽马超先验分布
Figure GDA0002450397090000104
其中超参数c,d被设置为常量，且应有c,d→0。
贝叶斯准则下，可推导wl的后验分布为非负高斯分布，且m≠l时wm|(ym,γ,λ)与wl|(yl,γ,λ)无关，则有
Figure GDA0002450397090000105
其中
Figure GDA0002450397090000106
Γ＝diag(γ)。
基于期望最大化(EM)思想，将wl,l＝1,2,,…,L视为隐藏变量，可推导空间联合稀疏表示向量γ和正则化参数λ的更新公式：
(1)推导γ[n]的更新公式，忽略p(y1,y2,…,yL,w1,w2,…,wL,γ,λ)中的γ[n]，从而最大化
Figure GDA0002450397090000107
对(15)式中的γ[n]求导，令导数为0时可求得：
Figure GDA0002450397090000111
其中服从后验分布p(wl|yl,γ,λ)的wl的第n个元素的二阶矩
Figure GDA0002450397090000112
为
Figure GDA0002450397090000113
(2)推导正则化参数λ的更新公式，忽略其他无关项，只需要将logp(γ|λ)p(λ)最大化，可得
Figure GDA0002450397090000114
以上所述实施例仅是为充分说明本发明而所举的较佳的实施例，本发明的保护范围不限于此。本技术领域的技术人员在本发明基础上所作的等同替代或变换，均在本发明的保护范围之内。本发明的保护范围以权利要求书为准。

Claims (9)

1.一种非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，包括以下步骤：
S1、对非均匀传感器阵列中各传感器接收通道采集的样本序列作滑窗快速傅里叶变换，并构造L个窄频带上的实值加权样本协方差矢量yl；
S2、构造与各窄频带上的实值加权样本协方差矢量yl相对应的过完备阵列流型矩阵Φl；
S3、利用L个窄频带上的实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl设置迭代程序中空间联合稀疏表示向量γ的初始值γinit，并根据实值加权样本协方差矢量yl、过完备阵列流型矩阵Φl和空间联合稀疏表示向量的初始值γinit计算正则化参数λ的初始值λinit；
S4、利用L个窄频带上的实值加权样本协方差矢量yl和过完备阵列流型矩阵Φl更新空间联合稀疏表示向量γ和正则化参数λ；
S5、当第i次迭代得到的空间联合稀疏表示向量γi与第i+1次迭代得到的空间联合稀疏表示向量γi+1之间满足||γi-γi+1||2/||γi||2＜α或者达到了最大迭代次数Pmax，终止迭代；否则继续返回S4迭代；其中α为一个任意趋于0的小数；
S6、当满足迭代终止条件时，搜索当前更新所得的空间联合稀疏表示向量γ的所有峰值，并根据空间联合稀疏表示向量γ的峰值位置输出宽带信号的波达方向估计值；
步骤S3包括：
S31、计算空间联合稀疏表示向量γ的初值
Figure FDA0002424416180000011
其中⊙为Hadamard积，||·||F为Frobenius范数；
S32、计算
Figure FDA0002424416180000012
其中
Figure FDA0002424416180000013
进而计算正则化参数λ的初值
Figure FDA0002424416180000014
其中||·||2为l2范数；
S33、设置正则化参数λ的无信息超先验分布参数固定参数值c,d，其中c,d→0。
2.如权利要求1所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S1包括：
S11、对非均匀传感器阵列中各传感器接收通道采集的样本序列作滑窗快速傅里叶变换，得到输出快拍矢量序列xl(t),l＝1,2,…,L；t＝1,2,…,T；
S12、利用输出快拍矢量序列构造实值加权样本协方差矢量yl,l＝1,2,…,L。
3.如权利要求1所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S4包括：
S41、利用实值加权样本协方差矢量yl、过完备阵列流型矩阵Φl以及上一轮迭代更新的空间联合稀疏表示向量γ，更新隐变量向量wl所服从的非负高斯分布的均值向量μl与协方差矩阵Σl，其中，若当前为第一轮迭代则空间联合稀疏表示向量γ＝γinit；
S42、利用隐变量向量wl所服从的非负高斯分布的均值μl、协方差矩阵Σl和上一轮迭代更新的正则化参数λ，更新空间联合稀疏表示向量γ的值；其中，若当前为第一轮迭代则正则化参数λ＝λinit；利用迭代更新的空间联合稀疏表示向量γ，更新正则化参数λ的值。
4.如权利要求2所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S11包括：
S111、利用窗函数对非均匀传感器阵列中各传感器接收通道采集的样本序列进行滑动数据选取，其中阵元个数为M；
S112、对各窗口中的加权非均匀阵列采样数据做快速傅里叶变换，其中，所加窗函数的窗长以及快速傅里叶变换点数为NFFT，相邻窗交迭Noverlap点；
S113、根据目标宽带信号所覆盖的频域范围，在第t个快速傅里叶变换窗口所得到的NFFT个频域值中选择L个作为待使用频带上的值，由于有M个传感器接收通道，可构造L个子频带上的每个子频带含有M个元素的向量：xl(t),l＝1,2,…,L；t＝1,2,…,T，其中每个向量是该频带上的“频域快拍”；
S114、对非均匀阵列采样数据共做T个滑窗快速傅里叶变换，可构造L个子频带上每个频带上T个“频域快拍”的数据向量：xl(t),l＝1,2,…,L；t＝1,2,…,T。
5.如权利要求2所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S12包括：
S121、计算L个窄频带上样本协方差矩阵
Figure FDA0002424416180000021
其中(·)H表示共轭转置；
S122、对矩阵
Figure FDA0002424416180000022
作向量化得到样本协方差矢量
Figure FDA0002424416180000023
S123、用加权矩阵J对
Figure FDA0002424416180000024
进行加权处理以去除噪声项，得到加权样本协方差矢量
Figure FDA0002424416180000025
其中权重矩阵J＝[J1,J2,…JM-1]T，Jm＝[em(M+1)-M+1,em(M+1)-M+2,…,em(M+1)]，含M2个元素的向量em除第m个元素为1之外其他所有元素都为0，其中(·)T表示矩阵转置；
S124、利用J与
Figure FDA0002424416180000031
构造加权样本协方差矢量的协方差矩阵
Figure FDA0002424416180000032
S125、构造实值加权样本协方差矢量
Figure FDA0002424416180000033
其中Re(·)与Im(·)分别表示取实部与取虚部。
6.如权利要求1所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S2包括：
S21、以一定的角度网格间距构造可以覆盖阵列探测角度范围的密集角度抽样网格
Figure FDA0002424416180000035
S22、对于角度
Figure FDA0002424416180000036
构造L个窄频带上与实值加权样本协方差矢量相对应的阵列导向矢量
Figure FDA0002424416180000037
其中
Figure FDA0002424416180000038
Figure FDA0002424416180000039
为Kronecker积，原始阵列导向矢量
Figure FDA00024244161800000310
fl为第l个频带的频率，dm为第m个传感器相对于参考传感器的间距，v为波速，(·)*表示求共轭；
S23、对于网格
Figure FDA00024244161800000311
中的每一个角度
Figure FDA00024244161800000312
计算
Figure FDA00024244161800000313
并由此构造L个窄频带上与各实值加权样本协方差矢量相对应的过完备阵列流型矩阵Φl：
Figure FDA00024244161800000314
7.如权利要求3所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S41包括：
S411、利用上一轮迭代所得的空间联合稀疏表示向量γ，更新L个窄频带上的隐变量向量wl所服从的非负高斯分布的协方差矩阵
Figure FDA00024244161800000315
其中Γ＝diag(γ)；
S412、更新L个窄频带上的隐变量向量wl所服从的非负高斯分布的均值
Figure FDA00024244161800000316
8.如权利要求3所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S42包括：
S421、计算L个窄频带上的隐变量向量wl中各元素的二阶矩：
Figure FDA0002424416180000041
其中
Figure FDA0002424416180000042
特别当
Figure FDA0002424416180000043
时，直接令
Figure FDA0002424416180000044
S422、利用计算所得的
Figure FDA0002424416180000045
以及上一轮迭代所更新的正则化参数λ，更新空间联合稀疏表示向量γ的第n个元素
Figure FDA0002424416180000046
S423、利用上一轮迭代所更新的空间联合稀疏表示向量γ，更新正则化参数λ：
Figure FDA0002424416180000047
9.如权利要求6所述的非均匀传感器阵列宽带信号波达方向估计方法，其特征在于，步骤S6包括，当满足迭代终止条件时，搜索当前更新所得的空间联合稀疏表示向量γ的所有峰值，其峰值位置与网格
Figure FDA0002424416180000048
中对应的角度即为宽带信号的波达方向估计值。
CN201811178910.3A 2018-10-10 2018-10-10 一种非均匀传感器阵列宽带信号波达方向估计方法 Active CN109407045B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811178910.3A CN109407045B (zh) 2018-10-10 2018-10-10 一种非均匀传感器阵列宽带信号波达方向估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811178910.3A CN109407045B (zh) 2018-10-10 2018-10-10 一种非均匀传感器阵列宽带信号波达方向估计方法

Publications (2)

Publication Number Publication Date
CN109407045A CN109407045A (zh) 2019-03-01
CN109407045B true CN109407045B (zh) 2020-05-22

Family

ID=65466930

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811178910.3A Active CN109407045B (zh) 2018-10-10 2018-10-10 一种非均匀传感器阵列宽带信号波达方向估计方法

Country Status (1)

Country Link
CN (1) CN109407045B (zh)

Families Citing this family (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110045323B (zh) * 2019-03-14 2022-02-18 电子科技大学 一种基于矩阵填充的互质阵稳健自适应波束形成算法
CN110109050B (zh) * 2019-04-08 2022-05-03 电子科技大学 嵌套阵列下基于稀疏贝叶斯的未知互耦的doa估计方法
CN110457748B (zh) * 2019-07-04 2022-10-11 中国人民解放军63892部队 一种等水平的二部覆盖阵的试验设计方法
CN110579737B (zh) * 2019-07-17 2022-10-11 电子科技大学 一种杂波环境中基于稀疏阵列的mimo雷达宽带doa计算方法
CN111505567B (zh) * 2020-03-25 2022-09-06 山东省科学院海洋仪器仪表研究所 基于单矢量水听器波达方向估计的多目标跟踪方法
CN111983553B (zh) * 2020-08-20 2024-02-20 上海无线电设备研究所 一种基于互质多载频稀疏阵列的无网格doa估计方法
CN111929638A (zh) * 2020-09-24 2020-11-13 成都启英泰伦科技有限公司 一种语音波达方向估计方法及装置
CN112198508B (zh) * 2020-10-29 2024-07-26 中国人民武装警察部队工程大学 一种基于支撑集约束的雷达目标成像与识别方法
CN112487703B (zh) * 2020-11-09 2024-05-28 南京信息工程大学滨江学院 基于稀疏贝叶斯在未知噪声场的欠定宽带信号doa估计方法
CN113406560B (zh) * 2021-05-20 2024-07-16 华南理工大学 一种非相干分布宽带源的角度和频率参数估计方法
CN113411142B (zh) * 2021-06-29 2022-09-06 齐鲁工业大学 宽带确定性最大似然测向方法、系统、终端及存储介质
CN113589255B (zh) * 2021-08-23 2023-08-01 武汉大学 一种基于多频联合稀疏贝叶斯学习的到达角估计方法
CN114660535A (zh) * 2022-03-15 2022-06-24 中国人民解放军火箭军工程大学 一种基于离格的稀疏阵列宽带信号欠定doa估计算法
CN115656998B (zh) * 2022-11-22 2023-06-13 中国人民解放军空军预警学院 一种低样本数下阵列信号自适应检测方法与系统
CN115792792B (zh) * 2022-11-22 2023-10-27 中国人民解放军空军预警学院 极低样本下基于协方差矩阵迭代估计的测向方法与系统

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103954950B (zh) * 2014-04-25 2016-09-07 西安电子科技大学 一种基于样本协方差矩阵稀疏性的波达方向估计方法

Also Published As

Publication number Publication date
CN109407045A (zh) 2019-03-01

Similar Documents

Publication Publication Date Title
CN109407045B (zh) 一种非均匀传感器阵列宽带信号波达方向估计方法
Wen et al. Direction finding in MIMO radar with large antenna arrays and nonorthogonal waveforms
WO2021139208A1 (zh) 一种基于特定频率组合信号的一维doa估计方法
US9357293B2 (en) Methods and systems for Doppler recognition aided method (DREAM) for source localization and separation
CN110113085B (zh) 一种基于协方差矩阵重构的波束形成方法及系统
Zheng et al. Classification and localization of mixed near-field and far-field sources using mixed-order statistics
CN107015190A (zh) 基于虚拟阵列协方差矩阵稀疏重建的互质阵列波达方向估计方法
CN108845325B (zh) 拖曳线列阵声纳子阵误差失配估计方法
CN107576931B (zh) 一种基于协方差低维度迭代稀疏重构的相关/相干信号波达方向估计方法
CN104408278A (zh) 一种基于干扰噪声协方差矩阵估计的稳健波束形成方法
CN103901395A (zh) 一种冲击噪声环境下相干信号波达方向动态跟踪方法
CN112487703B (zh) 基于稀疏贝叶斯在未知噪声场的欠定宽带信号doa估计方法
Ahmad et al. Wideband DOA estimation based on incoherent signal subspace method
Tian et al. Passive localization of mixed sources jointly using MUSIC and sparse signal reconstruction
Steinwandt et al. Beamspace direction finding based on the conjugate gradient and the auxiliary vector filtering algorithms
Zhang et al. FFT-based DOA estimation for coprime MIMO radar: A hardware-friendly approach
Ang et al. Multiband grid-free compressive beamforming
Zhang et al. An extended TOPS algorithm based on incoherent signal subspace method
CN108594165B (zh) 一种基于期望最大化算法的窄带信号波达方向估计方法
CN109116338A (zh) 一种基于四阶累积量的宽带凸优化doa估计方法
Gan et al. DOA estimation of wideband signals based on slice-sparse representation
Hassanien et al. Single-snapshot beamforming using fast iterative adaptive techniques
Tian et al. A novel DOA estimation for distributed sources in an impulsive noise environment
Dong et al. Experimental study on the performance of DOA estimation algorithm using a coprime acoustic sensor array without a priori knowledge of the source number
Pan Spherical harmonic atomic norm and its application to DOA estimation

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