capon.m

doa估计经典capon算法, 简单波达方向估计算法,空间角度分辨率有限,有待进一步改进算法

%七单元天线阵DOA估计
clear;clc;
d=1;               %天线阵元的间距
lma=2;             %信号中心波长
q1=1*pi/4;q2=1*pi/3;q3=1*pi/6;q4=3*pi/4;%四输入信号的方向
A1=[exp(-2*pi*j*d*[0:6]*cos(q1)/lma)]';%求阵因子
A2=[exp(-2*pi*j*d*[0:6]*cos(q2)/lma)]';
A3=[exp(-2*pi*j*d*[0:6]*cos(q3)/lma)]';
A4=[exp(-2*pi*j*d*[0:6]*cos(q4)/lma)]';
A=[A1,A2,A3,A4];                       %得出A矩阵
n=1:1900;
v1=.015;                               %四信号的频率
v2=.05;
v3=.02;
v4=.035;
d=[1.3*cos(v1*n);1*sin(v2*n);1*sin(v3*n);1*sin(v4*n)];%构造输入信号矢量
U=A*d;                             %总的输入信号
U1=(U)';
c=cov(U*U1);                       %总输入信号的协方差矩阵
[s,h]=eig(c);                      %求协方差的特征矢量及特征值
Vn=s(:,[5:7]);                     %取出与零特征值对应的特征矢量
ci=inv(c);                         %求协方差矩阵的逆矩阵
q1b=[pi/180:pi/180:pi];
for n=1:length(q1b)
    q1a(n)=q1b(n);
    A1a=[exp(-2*pi*j*1*[0:6]*cos(q1a(n))/lma)]';
    Pcapon(n)=(A1a)'*A1a*(inv((A1a)'*Vn*(Vn)'*A1a));%应用Cpon法估计输出
    T(n)=q1a(n);
    P=abs(Pcapon);
end
figure(1)                                          %绘出各波达方向图
   polar(T,P);title('Capon-1')
  figure(2);
  T1=T*180/pi;
  semilogy(T1,P);grid;title('Capon-2')