Matlab波束形成程序.doc

(word完整版)Matlab波束形成程序波束形成与智能天线1.均匀线阵方向图%8阵元均匀线阵方向图，来波方向为0度clc；clear all；close all;imag=sqrt（-1）;element_num=8;%阵元数为8d_lamda=1/2；%阵元间距d与波长lamda的关系theta=linspace(pi/2,pi/2,200）;theta0=0；%来波方向w=exp（imag2*pi*d_lamda*sin（theta0)0:element_num1）；for j=1:length(theta) a=exp(imag*2*pid_lamda*sin（theta(j)*0:element_num-1）； p(j)=wa；endfigure；plot（theta,abs（p）,grid onxlabel(theta/radian)ylabel（amplitude)title（8阵元均匀线阵方向图）当来波方向为45度时，仿真图如下：8阵元均匀线阵方向图如下,来波方向为0度，20log（dB）随着阵元数的增加，波束宽度变窄，分辨力提高：仿真图如下：2。 波束宽度与波达方向及阵元数的关系clcclear all close allima=sqrt(-1）；element_num1=16； %阵元数element_num2=128; element_num3=1024； lamda=0。03; 波长为0.03米d=1/2*lamda; 阵元间距与波长的关系theta=0:0。5：90；for j=1：length（theta)； fai(j）=theta（j）*pi/180asin(sin（theta(j）pi/180）-lamda/(element_num1d）； psi（j)=theta(j）pi/180-asin（sin(theta(j)*pi/180）-lamda/（element_num2d）; beta（j）=theta（j）*pi/180asin(sin(theta(j)*pi/180）-lamda/（element_num3d)；endfigure；plot（theta,fai，r，theta,psi，b，theta,beta,g），grid onxlabel（theta）;ylabel（Width in radians)title(波束宽度与波达方向及阵元数的关系）仿真图如下：3. 当阵元间距 时,会出现栅瓣，导致空间模糊。仿真图如下：4. 类似于时域滤波，天线方向图是最优权的傅立叶变换仿真程序和仿真图如下：clcclear all close allima=sqrt(-1）;element_num=32; %阵元数source_num=1； %信源数d_lamda=1/2; %阵元间距与波长的关系theta=linspace(-pi/2,pi/2，200)；theta0=0; %来波方向w=exp(ima*2*pi*d_lamdasin（theta0）0:element_num1)；for j=1：length（theta); a=exp（ima*2*pi*d_lamdasin（theta（j）)0:element_num-1); p（j)=w*a； end figure；subplot(1,2，1）plot(theta,abs（p)）,grid onxlabel（theta/radian)ylabel（amplitude)title（按定义的方向图)pfft=fftshift(fft(w，128）；subplot(1，2,2)plot(linspace(-pi/2,pi/2，128），abs（pfft)）,grid onxlabel（theta/radian）ylabel（FFT_amplitude)title(最优权的傅里叶变换)5。 最大信噪比准则方向图和功率谱clc；clear all;close all;ima=sqrt(1）;element_num=8; 阵元数为8d_lamda=1/2; 间距为半波长theta=-90：0。5:90; %范围theta0=0; %来波方向theta1=20; %干扰方向L=512; %采样单元数for i=1：L amp0=10*randn（1); amp1=200randn（1）； ampn=1； s(：，i）=amp0exp（ima*2*pi*1/2*sin(theta0pi/180）*0:element_num1）； j(:，i)=amp1exp（ima2*pi*1/2sin（theta1pi/180）*0：element_num-1); n(：，i)=ampn*(randn(element_num，1）+ima*randn（element_num，1）;endRs=1/L*ss； %信号自相关矩阵Rnj=1/L（jj+nn）; 干扰+噪声的自相关矩阵V，D=eig(Rs，Rnj）; （Rs，Rnj）的广义特征值和特征向量D，I=sort(diag(D）； 排序Wopt=V（:,I（8）)； 最优权矢量for j=1：length（theta) a=exp(ima*2*pi*d_lamdasin(theta(j）pi/180)0：element_num1）; f(j)=Wopta； p（j)=a*Rsa+aRnj*a；endF=20*log10（abs(f)/max(max（abs(f）);P=20*log10(abs（p）/max(max(abs(p)）)；subplot（121)plot（theta,F)；grid on；hold onplot（theta0，-50：0，。);plot(theta1,50：0,。）xlabel（theta/0)；ylabel(F in dB）;title(maxSNR 方向图);axis（90 90 -50 0);hold onsubplot(122）plot(theta，P，r）；grid onxlabel(theta/0)；ylabel(功率 in dB)；title(maxSNR功率谱)仿真图如下：6. ASC旁瓣相消MSE准则clc;close all;clear allima=sqrt（1）;M=32; 辅助天线的数目d_lamda=.5;theta0=-30; 来波方向theta1=60； 干扰方向 L=512； 采样单元数s=zeros（1,512)； 预划分一个区域for ii=1:L； amp0=1randn(1); 信号的幅度随机产生,保证信号之间是不相关的 amp1=200randn(1); ampn=1； jam(:,ii)=amp1*exp(ima*2*pi*0.5sin(theta1pi/180）*0:M-1)+ampn*(randn（M，1）+imarandn(M，1)； %干扰+噪声 s（ii）=amp0*exp(ima*2*pi0。5*sin(theta0*pi/180)）+amp1*exp(ima2pi*0。5sin(theta1*pi/180）+ampn（randn(1，1）+ima*randn(1，1）);接收信号（信号+干扰+噪声） s0（ii）=amp0*exp(ima*2*pi*0。5sin(theta0*pi/180)）；endRx=1/Ljamjam； 噪声自相关矩阵，相当于X（t）r_xd=1/Ljams;Wopt=pinv(Rx)*r_xd；delta=s0-(s-Wopt*jam）;delta1=abs（mean(delta.2）（mean（delta)）.2) 方差theta=linspace（-pi/2，pi/2，200)；for jj=1：length（theta） a=exp(ima*2*pi*。5sin(theta(jj）)*0:M-1); f(jj)=Wopta;endF=20log10（abs（f)/（max(max（abs(f)）)；figure（1)plot（theta180/pi,F),grid on，hold onplot（theta0,-50：0,。）plot(theta1，50：0，。）xlabel(theta/o）；ylabel(F/dB)；title(MSE准则下的方向图)axis（-90 90 -50 0)；可为x轴和y轴设置一个极限范围仿真图如下:7。 %线性约束最小方差(LCMV）准则 clc；clear all ；close all;ima=sqrt(-1）；element_num=8; %阵元数d_lamda=1/2; %阵元间距与波长的关系theta=90:0。5:900。3; %搜索范围theta0=0； %三个信号源的来波方向theta1=30; theta2=60;L=512； 采样单元数for i=1：L； amp0=10randn(1）; amp1=100*randn(1）; amp2=10*randn(1); ampn=10； x（：,i)=amp0*exp（ima2*pi1/2sin（theta0pi/180）0：element_num-1）+.。 amp1exp(ima*2pi1/2sin（theta1*pi/180）0:element_num1)+.。. amp2*exp(ima2*pi1/2sin（theta2pi/180）*0：element_num-1)+.。 ampn*（randn(element_num，1）+ima*randn(element_num,1）);endRx=1/L*x x;steer1=exp（ima*2pi1/2*sin（theta0pi/180)0:element_num-1)steer2=exp(ima2pi*1/2*sin（theta1*pi/180）0:element_num1)steer3=exp(ima2*pi1/2*sin（theta2*pi/180)*0:element_num-1)C=steer1 steer2 steer3;F=1 0 1; 把三个方向都作为来波方向w=inv(Rx)*C(inv(C*inv（Rx)C)F;for j=1：length（theta）; a=exp(ima*2pi*d_lamda*sin(theta(j）pi/180）0：element_num1）; f（j)=w*a; p(j)=1/(a*inv（Rx）*a);endf=10log10(abs(f)/（max（max(abs(f）)）);figure(1)subplot（121)plot(theta,f）,grid on,hold onplot(theta0，20:0，。）plot(theta1,20：0,.）plot(theta2，20:0,.）xlabel(theta/o）；ylabel(F/dB)；title（Capon beamforming方向图）axis（90 90 -20 0)；可为x轴和y轴设置一个极限范围P=10*log10(abs(p)/(max（max(abs(p)）);subplot(122）plot（theta,P），grid on，hold onplot（theta0，-20:0,.)plot(theta1，20:0,。)plot(theta2,20:0,。）xlabel(theta/o)；ylabel（功率/dB)；title(Capon beamforming功率谱)仿真图如下：8. %Capon beamformingClc;clear all ；close all；ima=sqrt(1）;element_num=8; 阵元数d_lamda=1/2; %阵元间距与波长的关系theta=90：0.5:90； 范围theta0=0; 来波方向theta1=20; 干扰方向theta2=60； %干扰方向L=1000； 采样单元数for i=1:L; amp0=10*randn（1）;信号的幅度随机产生,保证信号之间是不相关的 amp1=200randn(1); amp2=200randn（1); ampn=3； x(:，i)=amp0*exp(ima*2*pi*1/2sin(theta0pi/180)0:element_num1)+。. amp1exp(ima*2*pi*1/2*sin(theta1*pi/180)*0:element_num-1）+。 amp2exp(ima*2pi*1/2*sin(theta2*pi/180）0:element_num-1）+.。 ampn（randn（element_num，1)+ima*randn(element_num，1);endRx=1/Lx x;R=inv（Rx）；steer=exp(ima2pi1/2sin（theta0pi/180）0:element_num1）;w=R*steer/(steerR*steer）;Capon最优权矢量for j=1:length(theta)； a=exp(ima*2*pid_lamda*sin（theta(j）pi/180)*0：element_num1）； f(j)=wa； p（j)=1/（aRa)；endF=20log10(abs(f)/(max(max(abs(f)）)；P=20log10（abs（p）/（max（max（abs（p）)；此处是功率的对数形式Figure;subplot(121)plot（theta,F)，grid on,hold onplot（theta0，50:0，。）plot（theta1，50:0,。）plot(theta2,50:0,.)xlabel(theta/o）;ylabel(F/dB);title（Capon beamforming方向图）axis（-90 90 -50 0）;subplot(122)plot（theta，P）,grid onxlabel（theta/o);ylabel(功率/dB);title(Capon beamforming功率谱)仿真图如下：9。 不同方法估计协方差矩阵的Capon波束形成clc;clear all;close all;ima=sqrt（-1);element_num=8； %阵元数为8d_lamda=1/2; 间距为半波长theta=-90:0.5:90； 范围theta0=0; 来波方向theta1=50； %干扰方向L=512;%采样单元数for i=1:L amp0=10randn（1)； amp1=50randn(1）; ampn=0。5; s(:，i）=amp0*exp(ima*2*pi*1/2sin(theta0pi/180)*0:element_num-1）； j（:,i）=amp1exp(ima2pi*1/2sin(theta1pi/180)0:element_num-1)； n(:，i）=ampn*exp（ima*2pi*randn(1)*0：element_num1）;endRx=1/L*(s+j+n)(s+j+n);%接收信号自相关矩阵Rnj=1/L*(j+n）*（j+n)；干拢+噪声的自相关矩阵e=exp（ima2*pid_lamda*sin（theta0*pi/180)0：element_num-1）;来波方向信号Wopt_Rx=inv(Rx）*e/(e*inv(Rx)e);%采用接收信号的权矢量Wopt_Rnj=inv（Rnj)*e/(e*inv（Rnj)*e）；采用干拢+噪声信号的权矢量for j=1:length(theta) a=exp(ima*2*pi*d_lamda*sin(theta（j）pi/180)0:element_num-1）; f1(j）=Wopt_Rx*a； f2(j)=Wopt_Rnja；endF1=20*log10(abs（f1）/max（max(abs(f1)）)；F2=20*log10（abs(f2）/max（max(abs（f2）)；figure;plot(theta，F1,theta，F2,r）,grid on;hold onplot(theta0,50:0，.);plot(theta1,50:0，。)xlabel(theta/0)；ylabel(F（1,2）/dB）；title(不同方法估计协方差矩阵的Capon波束形成）；axis(-90 90 50 0)；仿真图如下：10%多点约束的Capon波束形成和方向图clc；clear all ;close all;ima=sqrt（1）；element_num=8； %阵元数d_lamda=1/2; %阵元间距与波长的关系theta=-90:0.3:90； %搜索范围theta0=0； %来波方向theta1=20; 干扰方向theta2=50; 干扰方向L=512； %采样单元数Rx=zeros（element_num，element_num）；产生协方差矩阵for i=1:L； amp0=10*randn(1)； amp1=10*randn(1); amp2=50randn(1）； namp=0.5randn（1);%噪声的幅度随机产生，保证噪声与信号之间是不相关的 J（:，i)=amp1exp（ima2*pi1/2sin（theta1pi/180）*0：element_num-1)+.。 amp2*exp(ima2pi1/2*sin（theta2*pi/180)*0:element_num-1)+.。 nampexp（ima2*pirandn（1）0:element_num1); x（:,i）=amp0exp（ima2*pi1/2*sin(theta0*pi/180)0：element_num-1）+.。 J(：，i）;表示接收信号endRx=Rx+1/Lx*x；R=inv（Rx);w=amp0*exp(ima2*pi*1/2sin(theta0*pi/180）0：element_num1）+。amp2*exp（ima2pi*1/2*sin（theta2pi/180）*0:element_num1）+.amp1exp（ima*2*pi*1/2sin(theta1pi/180）*0：element_num-1）;for j=1：length（theta); a=exp(ima2pi*d_lamda*sin(theta（j）*pi/180)*0:element_num1)； f（j)=wa； p（j）=1/（aR*a)；endF=10*log10（abs(f）/max(max(f);取对数的方向图P=10log10(abs（p)/max(max（p）；%取对数的功率谱估计figure;subplot(121)plot（theta，F)，grid on，hold onplot(theta0，-50：0，。)；plot(theta1,50：0，。）;plot（theta2,-50：0，。）;xlabel(theta/o）；ylabel(F/dB)；title(capon beamforming方向图)subplot(122)plot（theta,P）,grid on,hold onplot（theta0，50：0,。);plot（theta1，-50：0，.）；plot(theta2,-50:0，。);xlabel(theta/o）;ylabel（P/dB）；title(capon beamforming功率谱）仿真图如下：大作业%自适应波束形成方向图clc；clear all; close all;ima=sqrt（-1）；element_num=8； %阵元数c=3e8;f=500e6；lamda=c/f；d=1/2*lamda; 阵元间距与波长的关系theta=90:0.5:90； 范围theta0=0; 来波方向theta1=45； 干扰方向theta2=60; %干扰方向L=512; 采样单元数for i=1:L； amp0=10randn（1）；%信号的幅度随机产生，保证信号之间是不相关的 amp1=100*randn（1); amp2=100randn(1）； ampn=10; x(:,i)=amp0exp（ima2pi*1/2sin（theta0*pi/180)*0:element_num-1）+.。 amp1*exp（ima2pi1/2sin(theta1pi/180)*0:element_num1）+。.。 amp2*exp（ima*2pi*1/2*sin（theta2*pi/180）0:element_num-1)+.。. ampn*(randn(element_num，1）+imarandn(element_num,1）；产生阵列信号endRx=1/Lx x;steer1=exp（ima*2pi1/2sin（theta0pi/180）0:element_num1）steer2=exp（ima*2pi1/2sin（theta1*pi/180）*0：element_num-1)steer3=exp(ima2*pi*1/2sin（theta2pi/180)*0：element_num-1）C=steer1 steer2 steer3;F=1 0 0;%把0度方向作为来波方向w=inv(Rx)C*（inv（C*inv(Rx)C）)*F;for j=1:length(theta)； a=exp（ima2pi1/2sin（theta（j）*pi/180)0:element_num-1)； f（j)=wa； p（j)=1/（ainv（Rx)*a）;endF=20*log10(abs（f)/(max(max（abs（f)）);P=20*log10（abs（p）/(max(max（abs(p)）);此处是功率的对数形式figuresubplot（121）plot（theta，P）,grid on,hold onxlabel（theta/o）;ylabel（F/dB);title(阵列信号方向图）axis（-90 90 -50 0)；subplot（122）plot（theta，F)，grid on,hold on;plot(theta0,50:0,。)plot(theta1，-50:0,。)plot（theta2,50:0，。）xlabel(theta/o）；ylabel(功率/dB)；axis(90 90 50 0)；title（自适应波束形成方向图)仿真图如下：当采样数L=2048时仿真图如下：当采样数L=2048时仿真图如下：当L=512，幅度为amp0=10randn（1）; amp1=200*randn（1)；amp2=200randn(1）；ampn=10；时仿真图如下：当L=512,幅度为amp0=100*randn(1)； amp1=200*randn(1）；amp2=200*randn(1);ampn=10；时仿真图如下：当L=512，幅度为amp0=100randn（1)； amp1=10randn（1）;amp2=10randn(1）；ampn=10;时仿真图如下：当L=2048，幅度为amp0=10randn(1）； amp1=50randn（1）；amp2=50randn(1)；ampn=10;时仿真图如下:当L=2048,幅度amp0=10randn(1）；amp1=500randn（1）;amp2=500randn（1)；ampn=10；时仿真图如下：