几种常用music算法

clear all; clc; close all radeg = 180/pi; derad=1/radeg; twpi = 2*pi; Melm= 7; kelm = 6; % dd = 0.5; %  d=0:dd:(Melm-1)*dd; %  iwave = 2; % number of DOA theta = [40 50]; % DOA n = 1024 % A=exp(-j*twpi*d.'*sin(theta*derad));%%%% direction matrix %构造相干信号源 S0=randn(iwave-1,n); S=[S0(1,:);S0]; X0=A*S; X=awgn(X0,10,'measured'); Rxxm=X*X'/n issp=2; %空间平滑 if issp==1 Rxx=ssp(Rxxm,kelm); elseif issp==2 Rxx=mssp(Rxxm,kelm); else Rxx=Rxxm; kelm=Melm; end [EV,D]=eig(Rxx); EVA=diag(D)'; [EVA,I]=sort(EVA); EVA=fliplr(EVA); EV=fliplr(EV(:,I)); for iang = 1:361 angle(iang)=(iang-181)/2; phim=derad*angle(iang); a=exp(-j*twpi*d(1:kelm)*sin(phim)).'; L=iwave; En=EV(:,L+1:kelm); SP(iang)=(a'*a)/(a'*En*En'*a); end SP=abs(SP); SPmax=max(SP); SP=10*log10(SP/SPmax); figure h=plot(angle,SP); set(h,'Linewidth',2) xlabel('angle (degree)') ylabel('magnitude (dB)') axis([-90 90 -60 0]) set(gca, 'XTick',[-90:10:90],'YTick',[-60:10:0]) grid on hold on legend('平滑MUSIC') function crs = mssp(cr, K) % modified spatial smoothing [M,MM]=size(cr); N=M-K+1; J = fliplr(eye(M)); crfb = (cr + J*cr.'*J)/2; crs = zeros(K,K); for in =1:N crs = crs + crfb(in:in+K-1,in:in+K-1) end crs = crs / N; end %End mssp.m function crs = ssp(cr, K) % spatial smoothing [M,MM]=size(cr); N=M-K+1; crs = zeros(K,K); for in =1:N crs = crs + cr(in:in+K-1,in:in+K-1) end crs = crs / N; end