基于扩展互质阵列，使用空间平滑，MUSIC
clear all
close all
derad = pi/180; % deg -> rad
-
radeg = 180/pi;
twpi=2*pi;
f=1*10^8;
c=3*10^8;
lamda=c/f;
d=lamda;
M=3; %阵元数
N=4;
iwave=2;
theta=[-70 80];
snr=100;
n = 500;
S1=exp(j*pi/6*[1:n]); %信源
S2=exp(j*pi/4*[1:n]);
S=[S1;S2];
A1=exp(-j*twpi*[0:2*M-1]'*N*d.*sin(theta*derad)/lamda);
A2=exp(-j*twpi*[0:N-1]'*M*d.*sin(theta*derad)/lamda);
A=[A1;A2];
X=awgn(A*S,snr,'measured');
%协方差
Rxx=X*X'/n;
z=Rxx(:);
z1=roundn(z,-4);
z2=unique(z1);
%z2=reshape(z1,1,length(z1(:)));
z3=angle(z2);
[z4,I]=sort(z3);
Z=(z2(I)); %虚拟阵列接收信号
for i=1:iwave;
B(:,i)=kron(conj(A(:,i)),A(:,i));
end
%-M*N*d:d:M*N*d
r=-M*N*d:d:M*N*d;
MM=(length(r)+1)/2;
Z1=Z(29:53); %取中间2*M*N+1个数据
%Rzz=Z1*Z1'/n;
R=ssp(Z1,M*N+1);
% for i=1:M*N+1
% for n=0:M*N
% sub(i)=(-i+1+n)*d;
% end
% end
%特征值
[EV,D]=eig(R);%%%%
EVA=diag(D)';
[EVA,Ix]=sort(EVA);
EVA=fliplr(EVA);
EV=fliplr(EV(:,Ix));
% MUSIC
for iang = 1:361
angle(iang)=(iang-181)/2;
phim=derad*angle(iang);
% aa=exp(1j*2*pi*sin(thetatest(t))*(0:MM-1)'*d/lambda);
a=exp(-j*twpi*[0:MM-1]'*d.*sin(phim)/lamda);
% a=exp(-j*twpi*d*sin(phim)).';
Un=EV(:,iwave+1:end); % 噪声子空间
% En=EV(:,iwave:end);
SP(iang)=1/(a'*Un*Un'*a);
end
%
SP=abs(SP);
SPmax=max(SP);
SP=10*log10(SP/SPmax);
h=plot(angle,SP);
set(h,'Linewidth',2)
xlabel('angle (degree)')
ylabel('magnitude (dB)')
axis([-90 90 min(SP) inf])
set(gca, 'XTick',[-90:30:90])
grid on
运行结果不正确
我已经检查修改了很多次，但不清楚哪个代码步骤存在问题，希望您能指正，自己实在找不出来了，thanks！
粗略看了一眼