untitled.m

关于高斯白噪声下信源数目的估计程序

clear all;
clc;
%参数赋值
c=3.0e+8;                                           %光速
lmda=0.24;                                          %L波段
fc=c/lmda;                                          %载频
N=16;                                               %阵元数
DOA_sigA=0;                                  %信号源1来波方向与阵法线方向夹角 
DOA_sigB=20;                                  %信号源2来波方向与阵法线方向夹角 
d=lmda/2;                                          %阵元间距
fs=1.0e+3;                                            %采样率
chirp_rateA=20; 
chirp_rateB=10;                                                %调频率
P=2;                                                  %信源数 
%计算信源复包络
t=[0:1/fs:1]';                                            %采样时间
s1=exp(j*pi*chirp_rateA*t.^2);                                     
s2=exp(j*pi*chirp_rateB*t.^2); 
carw=exp(j*2*pi*fc*t);                                %载波
sc1=s1.*carw;                                        %信源复包络                                       
sc2=s2.*carw;

%产生方向矢量1(steering vector)
sv_sigA=[];
for i=1:N
    sv_sigA=[sv_sigA;exp(j*2*pi*d*(i-1)*sin(DOA_sigA*pi/180)/lmda)];
end
%产生方向矢量2(steering vector)
sv_sigB=[];
for i=1:N
    sv_sigB=[sv_sigB;exp(j*2*pi*d*(i-1)*sin(DOA_sigB*pi/180)/lmda)];
end

x=sv_sigA*sc1.'+sv_sigB*sc2.'+randn(N,length(t));  %N元等距线阵接收信号矢量

Rx=zeros(N,N);%阵列输出协方差矩阵
for i=1:length(t)
    Rx=Rx+x(:,i)*x(:,i)';
end
Rx=Rx/length(t);

d=eig(Rx);%对R做特征值分解
%用AIC准则进行估计
AIC=[];

for k=0:N-1
    
    total1=1;
    total2=0;
    
    for i=k+1:N
        total1=total1*(d(i)^(1/(N-k)));
        total2=total2+d(i);
    end
    
     total3=(1/(N-k))*total2;
     
    AIC=[AIC,(-2)*log((total1/total3)^((N-k)*P))+2*k*(2*N-k)];
end