专利摘要

本发明公开了一种基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法。本发明利用准平稳宽带信号频谱的帧内相关性和帧间独立性为信号设置相应的块稀疏先验分布模型，通过块稀疏贝叶斯模型对稀疏信号进行估计，从而获得精度更高的估计结果。即首先对阵列接收信号进行适当的分帧处理，对每帧信号分别做傅里叶变换并在频域建立其各自的块稀疏贝叶斯模型，然后在各帧信号独立的假设下，联合所有帧的信息建立一个总的贝叶斯模型，并通过超参数向量来控制所有帧待重构信号的稀疏性，最后利用期望最大化准则获得该超参数向量的迭代更新式。本发明充分利用准平稳宽带阵列信号的短时平稳特性来建立块稀疏模型，能够获得更高的波达方向估计精度。

权利要求

1.基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法，其特征在于，当利用N个阵元组成的均匀线性阵列对K个远场宽带信号进行DOA估计时，包括下列步骤：

步骤1：将阵列接收信号均匀划分成L帧，并对每帧进行F点离散时间傅里叶变换，得到每帧的阵列接收信号向量yl＝[y1,l,y2,l,…,yF,l]，其中l＝1,2,...,L，f＝1,2,...,F；

将空间角度平均划分成M个网格，每个网格的角度表示为θm，其中m＝1,2,...,M；

步骤2：分帧构建块稀疏贝叶斯学习模型：yl＝Φxl+nl；

其中，xl表示每帧待重构信号源向量，nl表示每帧的噪声向量，Φ表示N×F的稀疏贝叶斯学习模型的感知矩阵，N为阵元数；

各阵列接收信号向量yl的条件概率密度函数为复高斯分布p(yl|xl；λ)～CN(Φxl,λINF)，λ表示噪声协方差，每帧待重构信号源向量xl的先验分布假设为复高斯分布p(xl；Γl,Bl)～CN(0,Σl)，其中Γl＝diag(γ1,l,...,γM,l)控制各xl的稀疏性，Bl表示第l帧信号频谱的帧内相关性矩阵，符号表示Kronecker积运算；因为每帧待重构信号源向量xl具有相同的稀疏性，即Γl非零元素的位置相同，将每帧的稀疏性Γl约束成同一矩阵Γ来控制这个相同的稀疏性，即Γ＝diag(γ1,γ2,...,γM)，且

步骤3：联合每帧的阵列接收信号向量yl，构建总稀疏贝叶斯学习模型：y＝Ψx+n，其中IL表示L×L的单位矩阵，阵列接收信号向量y＝vec([y1,y2,...,yL])，信号源向量x＝vec([x1,x2,...,xL])，噪声向量n＝vec([n1,n2,...,nL])，其中vec(·)为向量化算符；

信号源向量x的先验分布为复高斯分布Σ0＝diag(Σ1,...,ΣL)，阵列接收信号向量y的条件概率密度函数为复高斯分布INFL是维数为NFL×NFL的单位矩阵；

利用条件概率密度函数和先验分布函数得到信号源向量x的后验分布

p(x|y；Θ)～CN(μx,Σx)，则信号源向量x的后验均值μx＝Σ0ΨH(λINFL+ΨΣ0ΨH)-1y，后验协方差矩阵(·)H表示矩阵的共轭转置操作；

从而得到超参数集合Θ＝{Γ,B1,...,BL,λ}，并基于迭代更新公式对γm、Bl和噪声方差λ进行迭代更新，其中γm、Bl和λ的初始值为预设值：

其中，符号Tr(·)表示求矩阵的迹运算，(·)H表示矩阵的共轭转置操作，J＝M(l-1)+m；Φm＝Φ(:,(m-1)F+1:mF)；

当满足迭代停止条件时，停止对γm、Bl和λ的迭代更新，并执行步骤4；

步骤4：由M个γm中大于0的索引下标m所对应的角度θm获得信源波达方向估计结果。

2.如权利要求1所述的方法，其特征在于，迭代停止条件为迭代次数达到指定次数。

3.如权利要求1所述的方法，其特征在于，迭代停止条件为连续两次迭代估计信号的误差小于指定门限，其中表示基于当前迭代更新后的参数γm、Bl和λ计算得到的信号源向量的后验均值，表示基于当前迭代的上一次迭代后的参数γm、Bl和λ计算得到的信号源向量的后验均值。

4.如权利要求1所述的方法，其特征在于，各γm的初始值为1。

5.如权利要求1所述的方法，其特征在于，各Bl的初始值为单位矩阵。

说明书

技术领域

本发明属于阵列信号处理领域，具体涉及一种基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法。

背景技术

波达方向估计是阵列信号处理的一个重要研究方向，随着宽带信号的广泛应用，如何有效地实现对宽带信号波达方向的高精度、高分辨率估计成为当前的一个研究热点。现实生活中准平稳宽带阵列信号的DOA估计具有广泛的应用，如麦克风阵列对声源的定位、航空站系统对视频信号的定位等。在已有的宽带阵列信号DOA估计方法中，子空间类方法如双边相关变换(Two-side Correlation Transformation，TCT)算法，需要精确地预估信源方位，在预估不准确的情况下算法性能会急剧下降；稀疏重构类方法如分布式压缩感知同步正交匹配追踪(Distributed Compressive Sensing-SimultaneousOrthogonalMatching Pursuit，DCS-SOMP)算法虽然有较低的运算复杂度，但其重构性能较差，另外，将L1-SVD(Singular Value Decomposition)方法扩展到宽带信号处理的系列算法因为采用了奇异值分解，需要已知信源个数，而且在低信噪比的情况下性能不理想。

发明内容

本发明的发明目的在于：针对上述存在的问题，提供一种高估计精度的基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法。

本发明的基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法包括：(1)对阵列接收信号进行分帧处理，对每帧信号分别做离散时间傅里叶变换，利用空域稀疏性构建出稀疏模型；(2)分帧构建块稀疏贝叶斯学习模型；(3)在各帧信号独立的假设下，联合所有帧的信息建立一个总的贝叶斯模型，并将控制每帧待重构信号稀疏性的超参数向量设成一致，利用期望最大化准则得到超参数集合的迭代更新公式；(4)通过迭代求解超参数向量从而获得DOA估计结果。

当利用N个阵元组成的均匀线性阵列对K个远场宽带信号进行DOA估计时，本发明的具体步骤为：

步骤1：将阵列接收信号均匀划分成L帧，并对每帧进行F点离散时间傅里叶变换，得到每帧的阵列接收信号向量yl＝[y1,l,y2,l,...,yF,l]，其中l＝1,2,...,L，f＝1,2,...,F；

将空间角度平均划分成M个网格，每个网格的角度表示为θm，其中m＝1,2,...,M。

基于空域稀疏性构建关于yf,l的稀疏模型：

yf,l＝Af(θ)xf,l+nf,l f＝1,2,...,F,l＝1,2,...,L

其中Af(θ)是方向矩阵，xf,l是待重构信号源向量，nf,l是噪声(高斯白噪声)向量的傅里叶变换。

步骤2：分帧构建块稀疏贝叶斯学习模型：

yl＝Φxl+nl

其中xl表示每帧待重构信号源向量，nl表示每帧的噪声向量，Φ表示N×F的稀疏贝叶斯学习模型的感知矩阵，每帧的阵列接收信号向量yl的条件概率密度函数为复高斯分布p(yl|xl；λ)～CN(Φxl,λINF)，λ是噪声协方差，每帧待重构信号源向量xl的先验分布假设为复高斯分布p(xl；Γl,Bl)～CN(0,Σl)，其中Γl＝diag(γ1,l,...,γM,l)控制着xl的稀疏性，Bl揭示第l帧信号频谱的帧内相关性， 符号 表示Kronecker积运算。因为每帧的xl具有相同的稀疏性，即Γl非零元素的位置是一样的，所以将每帧的Γl约束成同一矩阵Γ来控制这个相同的稀疏性。即Γ＝diag(γ1,γ2,...,γM)，对应的有

步骤3：联合每帧的阵列接收信号向量yl，构建总稀疏贝叶斯学习模型：y＝Ψx+n，其中 其中IL表示L×L的单位矩阵，阵列接收信号向量y＝vec([y1,y2,...,yL])，信号源向量x＝vec([x1,x2,...,xL])，噪声向量n＝vec([n1,n2,...,nL])，其中vec(·)为向量化算符。

因信号源向量x的先验分布为复高斯分布 Σ0＝diag(Σ1,...,ΣL)，阵列接收信号向量y的条件概率密度函数为复高斯分布 INFL是维数为NFL×NFL的单位矩阵。

根据贝叶斯公式，利用条件概率密度函数和先验分布函数可以得到信号源向量的后验分布p(x|y；Θ)～CN(μx,Σx)，则其后验均值μx＝Σ0ΨH(λINFL+ΨΣ0ΨH)-1y，后验协方差矩阵 (·)H表示矩阵的共轭转置操作。

因此，可得超参数集合Θ＝{Γ,B1,...,BL,λ}，根据期望最大化准则得到各超参数的迭代更新公式：

符号Tr表示求矩阵的迹运算，J＝M(l-1)+m，采用MATLAB的记法定义：

Φm＝Φ(:,(m-1)F+1:mF)；从矩阵Φ的第(m-1)F+1列到第mF列的这个子矩阵

从向量μx的第(J-1)F+1个元素到第JF个元素的这个子向量

从矩阵Σx的第(J-1)F+1行到第JF行，第(J-1)F+1列到第JF列的这个子矩阵。

基于γm、Bl和λ预设的初始值进行迭代更新，当满足迭代停止条件时，停止对γm、Bl和λ的迭代更新，并执行步骤4。

其中，优选的初始化方式为：γm初始化为1，Bl初始化为单位矩阵。而λ则取决于应用的场景。

步骤4：由M个γm中大于0的索引下标m所对应的角度θm获得信源波达方向估计结果。

其中迭代停止条件可以是：迭代次数达到指定次数；或则连续两次迭代估计信号的误差 小于指定门限，其中 表示基于当前迭代更新后的参数γm、Bl和λ计算得到的信号源向量的后验均值， 表示基于当前迭代的上一次迭代后的参数γm、Bl和λ计算得到的信号源向量的后验均值。

综上所述，由于采用了上述技术方案，本发明的有益效果是：利用准平稳宽带信号功率谱的短时平稳特性，建立一个具有分块对角协方差矩阵的块稀疏贝叶斯学习模型，探索信号短时谱的帧内相关性，并联合利用各帧信息来估计揭示信源DOA的超参数向量，从而在不需要知道信源个数的前提下，获得信源波达方向估计结果，且在低信噪比、低快拍数、来波方向较为密集的情况下也有较好的估计精度。

附图说明

图1是本发明的处理流程图；

图2是本发明(提出算法)与其他两种算法(SBL-SML(stochastic maximumlikelihood，SML)算法、忽略相关法)的均方根误差随信噪比变化的比较图；

图3是本发明(提出算法)与忽略信号频谱相关性的算法、SBL-SML算法的均方根误差随帧数变化的比较图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面结合实施方式和附图，对本发明作进一步地详细描述。

参见图1，在将本发明的估计方法用在由N个阵元组成的均匀线性阵列对K个方向为 的远场宽带信号进行DOA估计时，具体实现步骤如下：

步骤1、将阵列接收信号划分成L帧，并对每帧进行F点离散时间傅里叶变换，则阵列接收信号可以表示成

其中方向矩阵 sf,l是信号源向量的傅里叶变换，nf,l是高斯白噪声向量的傅里叶变换， 是阵列对 方向上的信号在频率点f处的导向矢量，v表示信号传播速率。

为构建稀疏模型，将空间角度平均划分成M个网格，构造新的方向矩阵Af(θ)＝[af(θ1),...,af(θM)]，则阵列接收信号可以重新表示成：

yf,l＝Af(θ)xf,l+nf,l f＝1,2,...,F,l＝1,2,...,L

因为M＞＞K，所以每个待重构向量xf,l均是稀疏的，且具有相同的支撑集，支撑集代表着信号源方向。

步骤2、分帧构建块稀疏贝叶斯学习模型。

为综合利用每个频谱的信息，探索其相关性，分别给每帧信号构造联合矩阵Yl＝[y1,l,y2,l,…,yF,l]，Xl＝[x1,l,x2,l,...,xF,l]，Nl＝[n1,l,n2,l,...,nF,l]，则有

vec(Yl)＝diag(A1,A2,...,AF)vec(Xl),l＝1,...,L

为构造块稀疏结构向量，令 vec(Xl)和 之间实际上是行变换的关系，所以令D∈RMF×MF是初等行变换矩阵，且满足：

令yl＝vec(Yl)∈CNF×1，nl＝vec(Nl)∈CNF×1，G＝diag(A1,A2,...,AF)，则有

yl＝Gvec(Xl)+nl＝GDxl+nl＝Φxl+nl l＝1,...,L

假设所有噪声向量相互独立，且服从复高斯分布，方差均为λ，则阵列接收信号的条件概率密度函数也服从复高斯分布：

p(yl|xl；λ)～CN(Φxl,λINF),l＝1,...,L

INF是维数为NF×NF的单位矩阵。

假设不同方向上的信号源是相互独立的，即Xl的行向量 是相互独立的，也服从复高斯分布：

超参数γm,l≥0控制着Xl的行稀疏性，当γm,l＝0时， 为零，当γm,l＞0时，意味着对应方向θm有来波信号。Bm,l是正定矩阵，用来表示信号频谱之间的相关性，为避免过拟合，采取T-SBL(temporal extension of SBL)算法中的措施，将不同行的 约束成相同矩阵Bl。

根据 的分布以及高斯分布的特性，可以得到向量xl的先验分布也服从复高斯分布：

其中Γl＝diag(γ1,l,...,γM,l)， 是协方差矩阵，符号 表示Kronecker积运算。

步骤3、联合利用每帧信号，构建总的稀疏贝叶斯模型。

因为每帧待重构向量xl具有相同的稀疏结构，利用这一先验知识，将所有 约束成同一个Γ来控制所有xl的块稀疏性，在此，假设 是相互独立的，构造矩阵X＝[x1,x2,...,xL]，则X的概率分布函数为：

其中Σ0＝diag(Σ1,...,ΣL)是协方差矩阵，估计矩阵X等价于估计其向量化算符x＝vec(X)，

p(x；Γ,B1,...,Bl)＝p(X；Γ,B1,...,Bl)～CN(0,Σ0)

定义Y＝[y1,y2,...,yL]，令y＝vec(Y)，则y的条件概率密度函数是：

其中 IL是维数为L×L的单位矩阵。

根据贝叶斯公式，利用条件概率密度函数和先验分布函数可以得到信号源向量的后验分布：

p(x|y；Θ)～CN(μx,Σx)

其中，超参数集合Θ＝{Γ,B1,...,BL,λ}，后验均值μx＝Σ0ΨH(λINFL+ΨΣ0ΨH)-1y，后验协方差矩阵

要得到超参数集合Θ的估计值，利用EM(期望最大化准则)方法，将x当作隐藏变量，最大化代价函数：

Θ(old)表示前一次迭代中超参数的估计值，E表示对x的后验分布logp(x|y；Θold)求期望。

在估计Γ和Bl时，只需要最大化代价函数中含有Γ和Bl的项，因此Q(Θ)可以简化为

然后对Q(Γ,B1,...,BL)分别求γm(m＝1,...,M)的偏导数，得

令上式为0，即可获得γm的迭代更新公式为：

符号Tr表示求矩阵的迹运算，其中J＝M(l-1)+m，采用MATLAB的记法定义：

对Q(Γ,B1,...,BL)分别求Bl(l＝1,...,L)的偏导数

同理可得到Bl的迭代更新公式：

估计噪声方差λ时，只需要最大化代价函数中含有λ的项，因此Q(Θ)可以简化为

利用T-SBL算法的简化步骤，可以得到噪声方差λ的估计式：

采用MATLAB的记法定义：Φm＝Φ(:,(m-1)F+1:mF)。

步骤4、利用上述各超参数的迭代更新公式进行迭代，直至连续两次迭代估计信号的误差 小于指定门限或者迭代次数达到指定次数，根据迭代停止时的输出结果Γ即可获得信源波达方向估计结果。

本发明的效果通过仿真图2和图3进行说明：

仿真实验中，信号源为语音信号，取其700Hz-1200Hz的频谱，四个信号源的真实方向是[-10°,-5°,3°,10°]，均匀线性阵列的阵元个数是28，阵元间距是半波长，分帧处理之后，每帧信号取了9个频点。蒙特卡洛实验次数为W＝200，用均方根误差(RMSE)来比较每个算法的性能，其定义为 是第w次实验估计的第k个信源的方向。

图2是帧数取5时，本发明与忽略信号频谱相关性的算法、SBL-SML算法的均方根误差随信噪比SNR变化的比较图，SNR取值范围是[4:12]dB，间隔为2dB。

图3是信噪比取7dB时，本发明与忽略信号频谱相关性的算法、SBL-SML算法的均方根误差随帧数变化的比较图，帧数的取值范围是[2:14]，间隔为3。

从实验结果可以看出，针对准平稳信号，通过分帧考虑信号频谱之间的相关性，本发明获得了更高的方向估计精度。

以上所述，仅为本发明的具体实施方式，本说明书中所公开的任一特征，除非特别叙述，均可被其他等效或具有类似目的的替代特征加以替换；所公开的所有特征、或所有方法或过程中的步骤，除了互相排斥的特征和/或步骤以外，均可以任何方式组合。

基于块稀疏贝叶斯学习的准平稳宽带阵列信号波达方向估计方法专利购买费用说明

Q：办理专利转让的流程及所需资料

A：专利权人变更需要办理著录项目变更手续，有代理机构的，变更手续应当由代理机构办理。

1：专利变更应当使用专利局统一制作的“著录项目变更申报书”提出。

2：按规定缴纳著录项目变更手续费。

3：同时提交相关证明文件原件。

4：专利权转移的，变更后的专利权人委托新专利代理机构的，应当提交变更后的全体专利申请人签字或者盖章的委托书。

Q:专利著录项目变更费用如何缴交

A：（1）直接到国家知识产权局受理大厅收费窗口缴纳,（2）通过代办处缴纳，（3）通过邮局或者银行汇款,更多缴纳方式

Q：专利转让变更，多久能出结果

A：著录项目变更请求书递交后，一般1-2个月左右就会收到通知，国家知识产权局会下达《转让手续合格通知书》。