lcmv dbf中的lcmv算法

dbf中的lcmv算法，对应于《一种低旁瓣数字波束形成方法》中的面阵仿真

clear all;

%阵元大小
N=8;           
M=8;
T=512;%快拍数
c=3*10^8;%光速
f=3*10^8;%频率
d=c/f/2;
wave=c/f;%波长
fs1=0.5*10^6; %信号频率
f1=0.7*10^6;  %干扰1频率
f2=0.8*10^6;%干扰2频率
fs=4*10^6;
SNR=-10;%信噪比，信号功率为单位1
SIR=0;%信干比
pn= 1/10^（SNR/10）;  %噪声功率
noise=wgn（N*MTpn）;  %噪声信号
t=0:1:T-1
phi=45*pi/180;          %俯仰角
theta_s =200*pi/180;    %信号方向
theta_i1=60*pi/180;    %干扰方向1
theta_i2=230*pi/180;    %干扰方向2
for n=0:1:N-1
    for m=0:1:M-1
        a_theta_s（m+n*M+11）  = exp（j*2*pi*n*sin（phi）.*cos（theta_s）*d/wave+j*2*pi*m*sin（phi）.*sin（theta_s）*d/wave）;  %信号导向矢量 N*M * 1
        a_theta_i1（m+n*M+11） = exp（j*2*pi*n*sin（phi）.*cos（theta_i1）*d/wave+j*2*pi*m*sin（phi）.*sin（theta_i1）*d/wave）;%干扰导向矢量 N*M * 1
        a_theta_i2（m+n*M+11） = exp（j*2*pi*n*sin（phi）.*cos（theta_i2）*d/wave+j*2*pi*m*sin（phi）.*sin（theta_i2）*d/wave）;%干扰导向矢量 N*M * 1
    end
end
%信干比为0时是单位增益
S0=exp（j*2*pi*fs1*t/fs）;  %信号
%J1=exp（j*2*pi*f1*t/fs）; %干扰信号1
J1=8*exp（j*2*pi*f1*t/fs）; %改变干扰功率进行对比
%J2=exp（j*2*pi*f2*t/fs）; %干扰信号2
J2=8*exp（j*2*pi*f2*t/fs）;  %改变干扰功率进行对比
A=[a_theta_sa_theta_i1a_theta_i2];
f=[100]‘;%约束条件
Rxx（N*MN*M）=0;
S=a_theta_s*S0;     % N*M * T
I1=a_theta_i1*J1;   % N*M * T
I2=a_theta_i2*J2;   % N*M * T
%X=S+I1+noise;  %当只有一个干扰信号时的阵列接收信号
X=S+I1+I2+noise;  %当有两个干扰信号时的阵列接收信号
for i=1:T
   Rxx=Rxx+X（:i）*X（:i）‘/T;  %生成信号（包括干扰）的自相关矩阵（平均）
end
  invRxx=inv（Rxx）;
  Wopt=invRxx*A*inv（A‘*invRxx*A）*f;
theta_angle=0:1:360;
theta_angle_=theta_angle*pi/180;
Theta_angle_N=length（theta_angle）;
for n=0:1:N-1
    for m=0:1:M-1
        a（m+n*M+1:）=exp（j*2*pi*n*sin（phi）.*cos（theta_angle_）*d/wave+j*2*pi*m*sin（phi）.*sin（theta_angle_）*d/wave）;
    end;
end;

F（1:）=Wopt‘*a;
H（1:）=abs（F（1:））;
G（1:）=H（1:）./max（H（1:））;
Fa（1:）=20*log10（G（1:））;
%%%%%%%%%%%%%重复5次%%%%%%%%%%%%%%%%
for k=1:5
  %  求出旁瓣的平均值（包含栅瓣）
  G_d=zeros（Theta_angle_N1）;
  for i=2:1:（Theta_angle_N-1）
    if （（G（ki）>G（ki-1））&&（G（ki）>G（ki+1）））
        G_d（i）=G（ki）;
    end
  end
  [aabb]=max（G_d）;
  G_d（bb）=0;
  [aabb]=max（G_d）
  theta_q（k）=（（bb-1）*1）*pi/180
%%%%%%%%%%%增加约束%%%%%%%%%%%%%%%
  for n=0:1:N-1
    for m=0:1:M-1
        a_theta_q（m+n*M+1k）=exp（j*2*pi*n*sin（phi）.*cos（theta_q（k））*d/wave+j*2*pi*m*sin（phi）.*sin（theta_q（k））*d/wave）;%旁瓣方向导向矢量 N*M * 1
    end;
  end;
  A（:k+3）=a_theta_q（:k）;
  f（k+3）=0.0001%约束条件
  Wopt=invRxx*A*inv（A‘*invRxx*A）*f;%得到最优权值
  F（k+1:）=Wopt‘*a;
  H（k+1:）=abs（F（k+1:））;
  G（k+1:）=H（k+1:）./max（H（k+1:））;
  Fa（k+1:）=20*log10（G（k+1:））;
end
%%%%%%%%%画图%%%%%%%%%%%%%%%%
%plot（theta_angleFa（1:））;
hold on
plot（theta_angleFa（5:））;







xlabel（‘角度（deg）‘）;
ylabel（‘增益（dB）‘）;
axis（[0 360 -90 0]）;

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       2898  2012-07-05 22:00  lcmv.m

----------- ---------  ---------- -----  ----

                 2898                    1