clc;
clear;
l=6;
lma=15;%??
%------入射角度设置-----
%ql=pi/2;
q2=pi/6;
q3=pi/2;
q4=3*pi/4;
%天线阵列因子
%A1=[exp(-2*pi*j*l/lma*[0:6]*cos(q1))]';
A2=[exp(-2*pi*j*l/lma*[0:6]*cos(q2))]';
A3=[exp(-2*pi*j*l/lma*[0:6]*cos(q3))]';
A4=[exp(-2*pi*j*l/lma*[0:6]*cos(q4))]';
A=[A2,A3,A4]; %求出A矩阵
R1=1e4;R2=5.0e4;R3=3e4;R4=8.12e4;
c=3.0e8;
taol=2*R1/c;
tao2=2*R2/c;
tao3=2*R3/c;
tao4=2*R4/c;
f0=20e6;
T=10e-6; %调制时间10us
B=3e6; %调频带宽30MHz
K=B/T; %调频斜率
Fs=10*B;Ts=l/Fs;
N=T/Ts;
%瑞利分布的热噪声
sigma=l;
m=rand(l,N);
noise=sqrt(2*log(1./m))*sigma;
t=linspace(-T/2,T/2,N);
%St1=exp(j*pi*K*(t-taol).^2+j*2*pi*f0*t)+0.2*rand(l,N);
St2=exp(j*pi*K*(t-tao2).^2+j*2*pi*f0*t);
St3=exp(j*pi*K*(t-tao3).^2+j*2*pi*f0*t);
St4=exp(j*pi*K*(t-tao4).^2+j*2*pi*f0*t);
st2=awgn(St2,20);
st3=awgn(St3,20);
st4=awgn(St4,20);
%-------输入信号矢量--------------------
d=[st2;st3;st4];
U=A*d;
figure(2);
plot(t,real(st2));
U1=(U)';
%-------输入信号协方差--------------------
c=cov(U*U1);
%------求协方差的特征矢量和特征值------
[s,h]=eig(c);
Vn=s(:,[1:3]);
qlb=[pi/180:pi/180:pi];
for i=1:length(qlb)
qla(i)=qlb(i);
Ala=[exp(-2*pi*j*1/lma*[0:6]*cos(qla(i)))]';
pmusic(i)=(Ala)'*Ala*inv(((Ala)'*Vn*(Vn)'*Ala));
T(i)=qla(i);
p1=abs(pmusic)/max(abs(pmusic));
end
figure(l);
T1=T*180/pi;
semilogy(T1,P1);
title('MUSIC Spectyum');
xlabel('Angle');
ylabel('Spectrum');
grid on;
评论0