Root-MUSIC_pinghua.rar_ROOT_music_root MUSIC_smoothing music_平滑R

clear all; close all; jj = sqrt(-1); N = 8; %阵元数为8 K = 4; %子阵数 P = N-K+1; %子阵阵元 d_lamda = 0.5; %阵元间距为半波长 M = 2; %信号源数为2 f0 = 12*10e3; %信号频率12k fs = 4*f0; %采样频率 theta = [22 30]/180*pi; %两信号源角度 SNR = 20; %信噪比为10dB L = 64; %快拍数256 t = (1:L)*1/fs; s1 = 1*rand(1,L).*exp(jj*2*pi*f0*t); % s2 = 1*rand(1,L).*exp(j*2*pi*f0*t); s2 = 1.5*s1; A = exp(-jj*2*pi*d_lamda*[0:(N-1)].'*sin(theta)); S = A*[s1;s2]; n = awgn(S,SNR,'measured'); X = S+n; %接收信号 Rx = X*X'/L; Ry = zeros(K,P,P); RRy = zeros(K,P,P); for k = 1:K Z = [zeros(P,k-1) eye(P,P) zeros(P,K-k)]; Q = [zeros(P,k-1) fliplr(eye(P,P)) zeros(P,K-k)]; Ry(k,:,:) = Z*Rx*Z'; RRy(k,:,:) = Q*Rx'*Q'; end RY = squeeze(sum(Ry,1)/K); RRY = squeeze(sum(RRy,1)/K); R = (RY+RRY)/2; [U,D] = eig(R); [D,index] = sort(diag(D),'descend'); UU = U(:,index); Us = UU(:,[1:M]); Un = UU(:,[M+1:P]); syms z; pz = z.^(0:P-1).'; pz1 = z.^(0:-1:-(P-1)).'; fz = z^(P-1)*pz1.'*Un*Un'*pz; a = sym2poly(fz); r = roots(a); rr = abs(r); [rr,index] = sort(abs(rr-1)); if (abs(angle(r(index(1)))-angle(r(index(3))))<0.05) root = r(index([1 2])); else root = r(index([1 3])); end % root = r(index([1,3])); th = asin(-angle(root)/pi)*180/pi