www.pudn.com > mvdr.zip > mvdr.m, change:2016-05-30,size:1265b


function [P_theta] =mvdr(M,s1,s2) 
% 定义基本变量 
theta_s1 = 0*pi/180;    
theta_s2 = 30*pi/180;   
sigma_s1 = 10^(s1/10);     
sigma_s2 = 10^(s2/10);    
sigma_n= 10^(0/10);  
f = 1000; 
c = 344; 
lamd = c/f; 
omega = 2*pi*f; 
     
 
% 表示坐标向量、速度、波数以及导向向量 
p_m = ones(2,1,M); 
p_m(1,1,:) = 0;  
p_m(2,1,:) = (0:M-1).*lamd/2;  
 
v_s1 = -[cos(theta_s1),sin(theta_s1)].'; 
k_s1 = omega/c*v_s1;  
 
v_s2 = -[cos(theta_s2),sin(theta_s2)].';  
k_s2 = omega/c*v_s2;  
 
a_theta0 = ones(M,2); 
for m=1:M 
    a_theta0(m,1) = exp(-1i*k_s1.'*p_m(:,:,m));  
    a_theta0(m,2) = exp(-1i*k_s2.'*p_m(:,:,m));  
end 
 
% 求解数据协方差 
rho_n = eye(M);  
Rs1 = sigma_s1*a_theta0(:,1)*a_theta0(:,1)'; 
Rs2 = sigma_s2*a_theta0(:,2)*a_theta0(:,2)'; 
Rn = sigma_n*rho_n; 
Rx = Rs1 + Rs2 + Rn; 
 
% 根据公式求解波束形成 
theta = (-90:.1:90); 
a_theta = ones(M,1); 
P_theta = ones(1,length(theta)); 
for n=1:length(theta) 
    v_theta = -[cos(theta(n)*pi/180),sin(theta(n)*pi/180)].'; 
    k_theta = omega/c*v_theta;  
    for m=1:M 
       a_theta(m,1) =  exp(-1i*k_theta.'*p_m(:,:,m));  
    end 
    w_theta = inv(Rx)*a_theta/(a_theta'*inv(Rx)*a_theta); 
    P_theta(1,n) = w_theta'*Rx*w_theta; % P(θ) = w(θ)'*Rx*w(θ) 
end 
 
P_theta = 10*log10(abs(P_theta));