CAPON和MUSIC MATLAB 算法

J=sqrt(-1); 
source_number=4; 
source_doa=[30 45 60 135]; 
sensor_number=7; 
snapshot_number=2000; 
snr=10; 
  
A=exp(-J*(0:sensor_number-1)'*pi*sin(source_doa*pi/180)); 
s=(randn(source_number,snapshot_number)+J*randn(source_number,snapshot_number))/sqrt(2); 
x=A*s; 
y=awgn(x,snr); 
R=y*y'/snapshot_number; 
  
[V,D]=eig(R); 
Un=V(:,1:sensor_number-source_number); 
Gn=Un*Un'; 
  
searching_doa=0:0.1:90; 
for i=1:length(searching_doa) 
    a_theta=exp(-J*(0:sensor_number-1)'*pi*sin(pi*searching_doa(i)/180))
    P_con(i)=abs(a_theta'*R*a_theta);
    P_BF(i)=abs((a_theta'*R*a_theta)./(a_theta'*a_theta)); 
    P_capon(i)=1./abs((a_theta'*inv(R)*a_theta)); 
    P_music(i)=1./abs((a_theta'*Gn*a_theta)); 
end 
plot(searching_doa,P_con/max(P_con),'k');hold on;
plot(searching_doa,P_BF/max(P_BF),'r'); hold on;  
plot(searching_doa,P_capon/max(P_capon),'g'); hold on;
plot(searching_doa,P_music/max(P_music),'b'); hold off;grid on;
xlabel('ang');
ylabel('功率谱估计');
legend('conditional spectrum','Bartlett spectrum','Capon spectrum','Music spectrum');