musicalgorithm.txt

music的matlab程序 可测定波达方向

%MUSIC ALGORITHM
%DOA Estimation by ULA
clear all; 
close all;
clc
N_x=1024; % Length of Signal 
N=8; % Size of Rx Matrix 

l=1.8;%波长
d=0.5*l;%阵元间距
M=3; % Number of Signals 
w=[pi/6 pi/10 pi/8]';%信号频率 
xx=pi/4;yy=0*pi/4;zz=-pi/4; %三个信号的入射角,小于5度就难以辨别了
B=[1 exp(-j*2*pi*d*sin(xx)/l) exp(-j*2*2*pi*d*sin(xx)/l) exp(-j*2*3*pi*d*sin(xx)/l) exp(-j*2*4*pi*d*sin(xx)/l) exp(-j*2*5*pi*d*sin(xx)/l) exp(-j*2*6*pi*d*sin(xx)/l) exp(-j*2*7*pi*d*sin(xx)/l);
1 exp(-j*2*pi*d*sin(yy)/l) exp(-j*2*2*pi*d*sin(yy)/l) exp(-j*2*3*pi*d*sin(yy)/l) exp(-j*2*4*pi*d*sin(yy)/l) exp(-j*2*5*pi*d*sin(yy)/l) exp(-j*2*6*pi*d*sin(yy)/l) exp(-j*2*7*pi*d*sin(yy)/l);
1 exp(-j*2*pi*d*sin(zz)/l) exp(-j*2*2*pi*d*sin(zz)/l) exp(-j*2*3*pi*d*sin(zz)/l) exp(-j*2*4*pi*d*sin(zz)/l) exp(-j*2*5*pi*d*sin(zz)/l) exp(-j*2*6*pi*d*sin(zz)/l) exp(-j*2*7*pi*d*sin(zz)/l)]'; %阵列流型,信号源决定行数，阵元数决定列数
xxx=2*exp(j*w*[0:N_x-1]);%仿真信号
x=B*xxx+randn(8,N_x)+j*randn(8,N_x); %加噪声
R=x*x';
[V D]=eig(R);
[lambda,index] = sort((diag(D)));
UU=V(:,index(1:5));
theta=-90:0.1:90;%ULA估计角度变化的范围和频率选择
for i = 1:length(theta)
AA=[1 exp(-j*2*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*2*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*3*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*4*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*5*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*6*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*7*pi*d*sin(theta(i)/180*pi)/l)];%方向矢量
WW=AA*UU*UU'*AA';
Pmusic(i)=abs(1/WW);%角谱
% PCap(i) = abs(1/(AA*inv(R)*AA'));
end
Pmusic = 20*log10(Pmusic);
sita=-90:0.1:90;
plot(sita,Pmusic);
grid
xlabel('rad') 
------------------------------


------------------------------
[parse]王文生 wrote:[/parse]
%MUSIC ALGORITHM
%DOA Estimation by ULA
clear all; 
close all;
clc
N_x=1024; % Length of Signal 
N=8; % Size of Rx Matrix 

l=1.8;%波长
d=0.5*l;%阵元间距
M=3; % Number of Signals 
w=[pi/6 pi/10 pi/8]';%信号频率 
xx=pi/4;yy=0*pi/4;zz=-pi/4; %三个信号的入射角,小于5度就难以辨别了
B=[1 exp(-j*2*pi*d*sin(xx)/l) exp(-j*2*2*pi*d*sin(xx)/l) exp(-j*2*3*pi*d*sin(xx)/l) exp(-j*2*4*pi*d*sin(xx)/l) exp(-j*2*5*pi*d*sin(xx)/l) exp(-j*2*6*pi*d*sin(xx)/l) exp(-j*2*7*pi*d*sin(xx)/l);
1 exp(-j*2*pi*d*sin(yy)/l) exp(-j*2*2*pi*d*sin(yy)/l) exp(-j*2*3*pi*d*sin(yy)/l) exp(-j*2*4*pi*d*sin(yy)/l) exp(-j*2*5*pi*d*sin(yy)/l) exp(-j*2*6*pi*d*sin(yy)/l) exp(-j*2*7*pi*d*sin(yy)/l);
1 exp(-j*2*pi*d*sin(zz)/l) exp(-j*2*2*pi*d*sin(zz)/l) exp(-j*2*3*pi*d*sin(zz)/l) exp(-j*2*4*pi*d*sin(zz)/l) exp(-j*2*5*pi*d*sin(zz)/l) exp(-j*2*6*pi*d*sin(zz)/l) exp(-j*2*7*pi*d*sin(zz)/l)]'; %阵列流型,信号源决定行数，阵元数决定列数
xxx=2*exp(j*w*[0:N_x-1]);%仿真信号
x=B*xxx+randn(8,N_x)+j*randn(8,N_x); %加噪声
R=x*x';
[V D]=eig(R);
[lambda,index] = sort((diag(D)));
UU=V(:,index(1:5));
theta=-90:0.1:90;%ULA估计角度变化的范围和频率选择
for i = 1:length(theta)
AA=[1 exp(-j*2*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*2*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*3*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*4*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*5*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*6*pi*d*sin(theta(i)/180*pi)/l) exp(-j*2*7*pi*d*sin(theta(i)/180*pi)/l)];%方向矢量
WW=AA*UU*UU'*AA';
Pmusic(i)=abs(1/WW);%角谱
% PCap(i) = abs(1/(AA*inv(R)*AA'));
end
Pmusic = 20*log10(Pmusic);
sita=-90:0.1:90;
plot(sita,Pmusic);
grid
xlabel('rad')