ism.m

非相干子空间算法源程序（CSM）－－－这是一种宽带聚焦矩阵算法

%非相干信号子空间（ISM）算法
clc
clear all
close all
M=8;                                       %阵元数
N=1024;                                    %快拍数
ts=0.0125;                                 %时域采样间隔
f0=100;                                    %入射信号中心频率 
f1=80;                                     %入射信号最低频率
f2=120;                                    %入射信号最高频率
c=1500;                                    %声速
lambda=c/f0;                               %波长
d=lambda/2;                                %阵元间距
SNR=5;                                     %信噪比
b=pi/180;  
theat1=12*b;                               %入射信号波束角1
theat2=8*b;                                %入射信号波束角2
n=ts:ts:N*ts;                            
theat=[theat1 theat2]';

%%%%%%%%%%%%%%%% produce signal %%%%%%%%%%%%%%%%

s1=chirp(n,80,1,120);                      %生成线形调频信号1
sa=fft(s1,N);                              %进行FFT变换
s2=chirp(n+0.125,80,1,120);                %生成线形调频信号2
sb=fft(s2,N);                              %进行FFT变换

%%%%%%%%%%%%%%%%%%%%% ISM算法 %%%%%%%%%%%%%%%%%%
P=1:2;
a=zeros(M,2);
sump=zeros(1,181);
for i=1:33
    f=80+(i-1)*1.25;
    s=[sa(i) sb(i)]';
    for m=1:M
        a(m,P)=exp(-j*2*pi*f*d/c*sin(theat(P))*(m-1))';
    end
    R=a*(s*s')*a';
    [em,zm]=eig(R);
    [zm1,pos1]=max(zm);
    for l=1:2
        [zm2,pos2]=max(zm1);
        zm1(:,pos2)=[];
        em(:,pos2)=[];
    end
    k=1;
    for ii=-90:1:90
        arfa=sin(ii*b)*d/c;
        for iii=1:M
            tao(1,iii)=(iii-1)*arfa;
        end
        A=[exp(-j*2*pi*f*tao)]';
        p(k)=A'*em*em'*A;
        k=k+1;
    end
    sump=sump+abs(p);
end
pmusic=1/33*sump;
pm=1./pmusic;
thetaesti=-90:1:90;
plot(thetaesti,20*log(abs(pm)));
xlabel('入射角/度');
ylabel('空间谱/dB');
title('ISM算法');
grid on