用MVDR方法估计相干函数.doc

用MVDR方法估计相干函数 一、MVDR谱 用MVDR方法做信号的谱估计是基于滤波器分解。 设输入X(n)是一组均值为0的随机数，作为长度为N，的K阶滤波器的输入。 滤波器的输出为,输出信号的能量为： 其中，。 考虑这样一个L*K的矩阵 并且有，k=0,1,2…K-1，K=L。矩阵叫做傅里叶矩阵并且是酉矩阵，则。在MVDR谱中为了使滤波器的输出方差最小，所以滤波器系数必须满足下面的约束条件： (1) 在这个约束条件下，输入X(n)在频率Wk处无失真，在其他频率处被衰减。这也相当于是使下面的公式最小： (2) 为拉格朗日乘数。使(2)式最小的解为： (3) 定义X(n)在Wk的谱为 (4) 把3式代入4式可得 (5) 由3式和5式可得出： 考虑所有的向量fk，总的表达式变为： ， ，为对角阵。 二、MVDR交叉谱 这里假设有两个0均值的稳态随机信号和，他们各自的谱分别为和。由前面的公式可得出两个滤波器为 （6） 和在Wk出的谱为： (7) 其中： (8) 设和分别为和的输出。定义和的交叉谱： (9) 由公式9可以的出 (10) 其中 是和的互相关矩阵。把6式代入10式中可得： (11) 三、MVDR相关函数 定义X1(n)和X2(n)幅值均方相关函数(MSC)为： (12) 通过公式11可以得出： (13) 把公式13和11代入公式12中，MSC变为： (14) 并且有 四、代码和仿真结果： function [MSC]=coherence_MVDR(x1,x2,L,K); xx1 = zeros(L,1); xx2 = zeros(L,1); r11 = zeros(L,1); r22 = zeros(L,1); r12 = zeros(L,1); r21 = zeros(L,1); F= zeros(L,K); l= [0:L-1]; f= exp(2*pi*l*j/K); for k = 0:K-1 F(:,k+1) = f.^k; end F = F/sqrt(L); n = length(x1); for i = 1:n xx1 = [x1(i);xx1(1:L-1)]; xx2 = [x2(i);xx2(1:L-1)]; r11 = r11 + xx1*conj(xx1(1)); r22 = r22 + xx2*conj(xx2(1)); r12 = r12 + xx1*conj(xx2(1)); r21 = r21 + xx2*conj(xx1(1)); end r11 = r11/n; r22 = r22/n; r12 = r12/n; r21 = r21/n; R11 = toeplitz(r11); R22 = toeplitz(r22); R12 = toeplitz(r12,conj(r21)); Dt1 = 0.01*r11(1)*diag(diag(ones(L))); Dt2 = 0.01*r22(1)*diag(diag(ones(L))); Ri11 = inv(R11 + Dt1); Ri22 = inv(R22 + Dt2); Rn12 = Ri11*R12*Ri22; Si11 = real(diag(F*Ri11*F)); Si22 = real(diag(F*Ri22*F)); S12 = diag(F*Rn12*F); MSC = real(S12.*conj(S12))./(Si11.*Si22); clc clear all close all n = 1024; nT= [0:n-1]; Nf= 5; f = zeros(Nf,1); f(1) = 0.05; f(2)