运行结果有误，但不知道哪部分代码存在问题，希望您能够指正帮助！谢！

Question

基于扩展互质阵列，使用空间平滑，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！

Answer 1

粗略看了一眼

1. 信源代码生成不对，Line17和18，改成S1 = randn(1, n) ;
2. 阵列流型矩阵生成不对，Line20~22
   indexA = sort(unique([(0:1:2*M-1).N (0:1:N-1).M]));
   A = exp(-1jtwpiindexA'*d.sin(thetaderad)/lamda);
3. B生成的不对，kron代码写的不对
4. ssp是什么？
5. 感觉你没有理解互质阵的概念，建议重新把论文看一遍，再推导一下公式就懂了