chaper3_3.asv

空间谱估计很好的学习例子

%%%%chaper3_3%%%
%利用空域最大熵方法求解信号的估计值
%sig_nums为信号源数目
%M为
clc;clear;
M=6;
sig_nums=2;
array_nums=8;
si=1;%为0号阵元的信号；
fai=[-30,30];%信号的来波角度
fai=fai*pi/180;
theta=pi*sin(fai);
y=zeros(1,(M+2));% 用来放8个信号
for k=0:7
    y(k+1)=(exp(j*k*theta(1))+exp(j*k*theta(2)))+randn(1);%求出8个阵元的信号值； 
end
r_corr=zeros(1,M);
for i=1:M
 r_corr(i)=1/2*(y(7)*conj(y(7-i))+y(8)*conj(y(8-i)));
end
 r_corr_0=1/2*(y(7)*conj(y(7))+y(8)*conj(y(8)));%求出r0的相关函数
r_corr_new=[r_corr_0,r_corr];%成为8列的矩阵
%%%%对R_corr进行赋值%%%%
R_corr=zeros(M+1,M+1);
for i=1:M+1
    R_corr(i:(M+1),i)=r_corr_new(1,1:(M+1)- (i-1));
    R_corr(i,i:(M+1))=conj(r_corr_new(1,1:(M+1)-(i-1)));
end
%%%求g(m)的值%%%%%
g_fmin=1;%先让g_fmin为1
G=zeros((M+1),1);
G(1,1)=g_fmin;
R=rank(R_corr)
c=[R_corr,G];
gm=rref(c);
gm1=gm(1,8);
gm=gm(2:7,8).';
 gm=gm/(gm1);
%%%%求功率空间谱%%%%
theta=linspace(-100,100,1000);
theta1=theta*pi/180;
Y_theta=zeros(1,1000);
for k=1:1000
    sums=0;
    for n=1:M
         sums=sums+gm(n)*exp(-j*n*pi*sin(theta1(k)));
    end
    Y_theta(k)=1/(abs(1+sums).^2);%这里白噪声的方差为1；
end
plot(theta,Y_theta)