MUSIC 和 ESPRIT 程序(仅供参考学习) 本人最近做 DOA 估计，MUSIC 和 ESPRIT 算法是最经典的算法，看到很多初学者都在求 助 MUSIC 和 ESPRIT 的程序，其实大家只要把原理搞明白了，很简单的。下面我把我编的程序 奉献给大家，仅供大家参考，如果有错误和不足之处请指正。 鉴于有些人下了程序不看为什么，就直接复制、粘贴了，我把程序转化为 PDF 文档，进 行了加密处理，不能复制，只能看，相信你看懂了自热就明白了，程序只是算法的实现而已， 最重要的是算法。 另外 MUSIC 和 ESPRIT 算法都需要进行特征分解，在雷达上应用受到限制，建议大家研究 下 DOA 的并行运算、非特征分解的 DOA 估计和一次快拍的 DOA 估计，有兴趣的同学可以和 我交流，大家共同进步。 程序我是经过测试过的，没有任何问题，如果你在仿真时出了问题请检查下你敲的代码 是不是有问题。 clc clear all format long N=200;%%快拍数 doa=[20 40]/180*pi; %%信号到达角， w=[pi/4 pi/3]';%%信号频率 M=8;%%阵元数 P=length(w); %%信号个数，也可以用特征分解的大特征值来决定 lambda=150;%波长 d=lambda/2;%阵元间距 snr=15;%%信噪比 %%%导向向量 B=zeros(P,M); for k=1:P B(k,:)=exp(-j*2*pi*d*sin(doa(k))/lambda*[0:M-1]); end B=B'; %%%导向向量 xx=2*exp(j*(w*[1:N])); %生成信号 x=B*xx; %%%%噪声平均 因为matlab产生的噪声不太好，不是严格意义上的白噪声， %%%平均后结果更好 [pp,ppp]=size(x); xx=zeros(pp,ppp); cycle=200; for kkk=1:cycle xx=xx+awgn(x,snr); end x=xx/cycle; %%%%噪声平均结束 R=x*x'; %数据协方差矩阵 %针对相干源的时候进行平衡 J=fliplr(eye(M)); R=R+J*conj(R)*J; By han in HIT 2010-4-14 

MUSIC 和 ESPRIT 程序(仅供参考学习) %%%%以下是MUSIC的程序 [U,V]=eig(R); UU=U(:,1:M-P); theta=-90:0.5:90; %%谱峰搜索 for ii=1:length(theta) AA=zeros(1,length(M)); for jj=0:M-1 AA(1+jj)=exp(-j*2*jj*pi*d*sin(theta(ii)/180*pi)/lambda); end WW=AA*UU*UU'*AA'; Pmusic(ii)=abs(1/WW); end Pmusic=10*log10(Pmusic/max(Pmusic)+eps); figure(1) plot(theta,Pmusic) xlabel('Angle \theta/degree') ylabel('P(\theta) /dB') title('Spectrum') grid on %%%%MUSIC程序结束 %%以下是ESPRIT程序 Rxx=R(1:M-1,1:M-1);%%%M-1维的自相关函数 Rxy=R(1:M-1,2:M);%%%M-1维的互相关函数 b=[zeros(1,M-2);eye(M-2)]; b=[b zeros(M-1,1)]; Cxx=Rxx-min(eig(Rxx))*eye(M-1); Cxy=Rxy-min(eig(Rxx))*b; a=eig(Cxx,Cxy); %找出最接近1的a值其对应的角度即为φ a1=abs(abs(a)-1); for i=1:P [c,d]=min(a1); a1(d)=1000; bb(i)=a(d); a(d)=1000; end if P>1 disp('The angles of signals are') else disp('The angle of signal is') end DOA=asin(angle(bb)/pi)/pi*180 %%%%ESPRIT程序结束 By han in HIT 2010-4-14