MVDR频率估计方法及其Matlab代码实现

1.前言

$MVDR$频率估计方法的理论知识此处略过,直接介绍最后得到的$MVDR$谱估计公式。

2.$MVDR$谱估计公式

$$\hat{P}_{MVDR} \left( w \right )= \frac{1}{a^{H} \left( w \right )\hat{R}^{-1} a\left( w \right )}$$

其中$a_{w} = a_{M} \left( w \right ) = \begin{bmatrix} 1 \\  e^{-jw} \\*** \\  e^{-j \left(M-1 \right)w}\end{bmatrix}$,$\hat{R} = E \left \{ x(n) \hat{x}(n) \right \} = \begin{bmatrix}
r(0) & r(1) & \cdots & r(M-1)\\
r(-1) & r(0) & \cdots & r(M-2)\\
\cdots & \cdots & \ddots & \vdots \\
r(1-M) & r(2-M) & \cdots & r(0)
\end{bmatrix}$

$M$为自相关矩阵的阶数。

3.算法步骤

  1. 由$N$个观测样本$x\left( 0 \right), x\left( 1 \right),\cdots,x\left( N*1 \right)$估计样本相关矩阵$\hat{R}$。
  2. 在$\left [ -\pi , \pi\right ]$内改变$w$,画出$\hat{P}_{MVDR} \left( w \right)$,峰值位置就是信号角频率的估计值。

4.算例及代码实现

4.1 算例

设随机过程$u \left( n \right)$为$u \left( n \right) = e^{j0.5\pi n + j \phi _{1}} + e^{-j0.3\pi n + j \phi _{2}} + v_{n}$,其中,$v_{n}$是0均值,方差为1的白噪声,$\phi _{1}$、$\phi _{2}$是相互独立并在$\left[ 0,2\pi \right]$上服从均匀分布的随机相位,使用$MVDR$方法进行信号频率估计,画出频率估计谱线。(要求:信号样本数取1000,估计的自相关矩阵为8阶。)

4.2 $Matlab$实现

$MVDR$方法的输出功率为信号功率加上一个数,因此利用$MVDR$方法进行谱估计时,若在某个频点有信号,则该点会出现一个比它真实功率值大一点的数值,若是没有信号,信号和噪声都被滤波器抑制,按公式估计的值会很小。

4.3 代码及下载

具体代码如下,点此下载源码

N = 1000;%样本点数
noise = (randn(1,N) + 1i * randn(1,N))/sqrt(2);%复噪声
%产生输入信号
n=0:N-1;
signal1 = exp(1i*0.5*pi*n + 1i*2*pi*rand);
signal2 = exp(-1i*0.3*pi*n + 1i*2*pi*rand);
un = signal1 + signal2 + noise;
%计算自相关矩阵
M = 8;%自相关矩阵阶数
for k = 1:N-M
xs(:,k) = un(k+M-1:-1:k).'; %样本矩阵
end
R = xs * xs'/(N-M);%用时间平均代替统计平均
%计算MVDR谱
NF = 2048;%MVDR扫描点数
for m = 1:NF
Aq = exp(-1i*2*pi*(m-1)/NF*(0:M-1)');
Pmvdr(m) = 1/(Aq' * inv(R) * Aq); %MVDR谱
end
Pmvdrshift = [Pmvdr(NF/2+1:NF) Pmvdr(1:NF/2)];%将0:NF移到-NF/2:NF/2
figure(1)
plot(-0.5:1/(NF-1):0.5,10*log10(abs(Pmvdrshift)/max(abs(Pmvdrshift))));
xlabel('w/2Π')
ylabel('归一化功率谱 (dB)')
title('MVDR频率估计方法');

你可能也喜欢

5 1 vote
Article Rating
Subscribe
提醒
0 评论
Inline Feedbacks
View all comments
0
Would love your thoughts, please comment.x
()
x
Scroll to Top