Capon-波束形成matlab仿真(附源代码)教学内容
波束形成Matlab程序
1?均匀线阵方向图%8阵元均匀线阵方向图,来波方向为clc; clear all; close all;0度imag=sqrt(_1);element_num=8;% 阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=li nspace(-pi/2,pi/2,200);theta0=0;% 来波方向w=exp(imag*2*pi*dl_lamda*sin(theta0)*[0:eleme nt_nu m-1]');for j=1:le ngth(theta)a=exp(imag*2*pi*dd_l amda*si n(theta(j))*[0:eleme nt_nu m-1]'); p(j)=w'*a;end figure;plot(theta,abs(p)),grid on xlabel('theta/radia n')ylabel('amplitude')title('8 阵元均匀线阵方向图')°2 8阵元均匀线阵方向图7654321-15 -1 -0 5 0 06thetaradian1 158当来波方向为45度时,仿真图如下8阵元均匀线阵方向图如下,来波方向为0 度,20log (dB)8阵元均苛銭阵方向图来波方向为0度随着阵元数的增加,波束宽度变窄,分辨力提高:仿真图如下Q d pE =ro 二2. 波束宽度与波达方向及阵元数的关系clcclear allclose 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; forj=1:length(theta);fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um1*d));psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um2*d));beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_ num3*d));endfigure; plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid on xlabel('theta'); ylabel('Width in radians')title(' 波束宽度与波达方向及阵元数的关系') 仿真图如下:。
Capon 波束形成matlab仿真(附源代码)
Capon波束形成阵列N=16,信号0-30θ︒=,干扰为160θ︒=,219θ︒=,345θ︒=,干扰功率分别为:40dB,35dB,50dB。
Capon波束形成后的方向图和功率谱如下:为了比较接收数据直接估计噪声协方差矩阵和利用干扰+噪声估计协方差矩阵的Capon 波束形成的差异,进行如下仿真:可以看出利用干扰+噪声估计协方差矩阵的方向图性能较优于接收数据直接估计噪声协方差矩阵的方向图。
代码:clc;clear all ;close all;ima=sqrt(-1);element_num=8; %阵元数d_lamda=1/2; %阵元间距与波长的关系theta=-90:0.5:90; %范围theta0=-30; %来波方向theta1=60; %干扰方向1theta2=19; %干扰方向2theta3=45; %干扰方向3L=1000; %采样单元数for i=1:L;amp0=10*randn(1);%信号的幅度随机产生,保证信号之间是不相关的amp1=100*randn(1);%输入阵列的噪声amp2=sqrt(10^3.5)*randn(1);%输入阵列的噪声amp3=sqrt(10^5)*randn(1);%输入阵列的噪声ampn=3;%噪声x(:,i)=amp0*exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]')+...amp1*exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]')+...amp2*exp(ima*2*pi*1/2*sin(theta2*pi/180)*[0:element_num-1]')+...amp3*exp(ima*2*pi*1/2*sin(theta3*pi/180)*[0:element_num-1]')+...ampn*(randn(element_num,1)+ima*randn(element_num,1));endRx=1/L*x* x';R=inv(Rx);steer=exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]');w=R*steer/(steer'*R*steer);%Capon最优权矢量for j=1:length(theta);a=exp(ima*2*pi*d_lamda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/(max(abs(f))));P=20*log10(abs(p)/(max(abs(p))));%此处是功率的对数形式figure;% subplot(121)plot(theta,F),grid on,hold onplot(theta0,-80:0,'.')plot(theta1,-80:0,'.')plot(theta2,-80:0,'.')plot(theta3,-80:0,'.')xlabel('theta/o');ylabel('F/dB');title('Capon beamforming方向图') % axis([-90 90 -50 0]);% subplot(122)figure;plot(theta,P),grid onxlabel('theta/o');ylabel('功率/dB');title('Capon beamforming功率谱')。
波束形成Matlab程序
1•均匀线阵方向图%8阵元均匀线阵方向图,来波方向为0度clc;clear all;close all;imag=sqrt(_1);element_num=8;% 阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=li nspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*dl_lamda*si n(theta0)*[0:eleme nt_nu m-1]');for j=1:le ngth(theta)a=exp(imag*2*pi*dd_l amda*si n(theta(j))*[0:eleme nt_nu m-1]'); p(j)=w'*a;endfigure;plot(theta,abs(p)),grid onxlabel('theta/radia n')ylabel('amplitude')title('8阵元均匀线阵方向图')8阵元均匀线阵方向图当来波方向为45度时,仿真图如下:8阵元均匀线阵方向图如下,来波方向为0度,20log (dB )8阵元均苛銭阵方向图来波方向为0度 S =-s ==d E B随着阵元数的增加,波束宽度变窄,分辨力提高:仿真图如下: Qp=二 d E ro2. 波束宽度与波达方向及阵元数的关系clcclear allclose 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/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um1*d));psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_n um2*d));beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_ num3*d));endfigure; plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid on xlabel('theta');ylabel('Width in radians')title(' 波束宽度与波达方向及阵元数的关系') 仿真图如下:。
自适应波束形成与Matlab程序代码注解
1.均匀线阵方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;%阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=linspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag);patterndBnorm=20*log10(patternmagnorm);figure(1)plot(theta*180/pi,patternmag);grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);hold on;figure(2)plot(theta,patterndBnorm,'r');grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);axis([-1.5 1.5 -50 0]);(2)仿真结果A.来波方向为0°归一化B.来波方向为45°C.随着阵元数的增加,波束宽度变窄,分辨力提高,仿真图如下:2.波束宽度与波达方向及阵元数的关系(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num1=16;element_num2=128;element_num3=1024;lambda=0.1;d=0.5*lambda;theta=0:0.5:90;for j=1:length(theta)fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num1*d)); psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num2*d)); beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num3*d)); endfigureplot(theta,fai,'r',theta,psi,'b',theta,beta,'g');grid on;xlabel('theta');ylabel('width in radians');title('波束宽度与达波方向及阵元数目的关系');legend('N=16','N=128','N=1024');(2)仿真结果结果3. 当阵元间距/2d λ>时,会出现栅瓣,导致空间模糊(1)仿真结果4. 类似于时域滤波,天线方向图是最优权的傅立叶变换(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;source_num=1;d_lambda=0.5;theta=linspace(-pi/2,pi/2,200);theta0=0;w=exp(imag*2*pi*d_lambda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j))*[0:element_num-1]'); p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag);patterndBnorm=20*log10(patternmagnorm);figure(1)subplot(1,2,1);plot(theta,patterndBnorm);grid on;xlabel('theta/radian');ylabel('amplitude/dB');axis([-2.0 2.0 -50 0]);subplot(1,2,2);pfft=fftshift(fft(w,256));pfftmag=abs(pfft);pfftmagnorm=pfftmag/max(max(pfftmag));pfftdB=20*log10(pfftmagnorm);pfftdBnorm=20*log10(pfftmagnorm);plot(linspace(-pi/2,pi/2,256),pfftdBnorm);grid on;xlabel('theta/radian');ylabel('FFT_amplitude/dB');axis([-2.0 2.0 -50 0]);(2)仿真结果5.最大信噪比准则方向图和功率谱(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=0.5;%间距为半波长theta=-90:0.5:90;%扫描围theta0=0;%来波方位theta1=20;%干扰方向L=512;%采样点数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);ampn=1;s(:,i)=amp0*exp(imag*2*pi*0.5*sin(theta0*pi/180)*[0:element_num-1]');j(:,i)=amp1*exp(imag*2*pi*0.5*sin(theta1*pi/180)*[0:element_num-1]');n(:,i)=ampn*exp(randn(element_num,1)+imag*randn(element_num,1)); endRs=1/L*s*s';%信号自相关矩阵Rnj=1/L*(j*j'+n*n'); %干扰+噪声的自相关矩阵[V,D]=eig(Rs,Rnj); %(Rs,Rnj)的广义特征值和特征向量[D,I]=sort(diag(D)); %特征向量排序Wopt=V(:,I(8));%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=Wopt'*a;p(j)=a'*Rs*a+a'*Rnj*a;endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-80:0,'.');plot(theta1,-80:0,'.');xlabel('theta/0');ylabel('F in dB');title('max-SNR 方向图');axis([-90 90 -80 0]);hold on;subplot(1,2,2);plot(theta,P,'r');grid on;xlabel('theta/0'); ylabel('功率 in dB'); title('max-SNR 功率谱'); grid on;axis([-90 90 -80 0]); (2)仿真结果6.ASC旁瓣相消----MSE准则(1) matlab 程序clc;close all;clear all;imag=sqrt(-1);M=32;%辅助天线数目d_lambda=0.5;%阵元间距theta0=-30;%来波方向theta1=60;%干扰方向L=512;%采样单元数s=zeros(1,512); %预划分一个区域for ii=1:Lamp0=1*randn(1);%信号的幅度随机产生,保证信号之间是不相关的amp1=200*randn(1);ampn=1;jam(:,ii)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:M-1]' )+ampn*(randn(M,1)+imag*randn(M,1)); %干扰+噪声s(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180))+amp1*exp(im ag*2*pi*d_lambda*sin(theta1*pi/180))+ampn*(randn(1,1)+imag*randn( 1,1));%接收信号(信号+干扰+噪声)s0(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180));endRx=1/L*jam*jam';r_xd=1/L*jam*s';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(imag*2*pi*d_lambda*sin(theta(jj))*[0:M-1]');f(jj)=Wopt'*a;endF=20*log10(abs(f)/max(max(abs(f))));figure(1)plot(theta*180/pi,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.'); xlabel('theta/°');ylabel('F/dB');title('MSE准则下的方向图'); axis([-90 90 -50 0]);(2)仿真结果7.线性约束最小方差(LCMV)准则(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=0.5;%阵元间距与波长的关系theta=-90:0.5:90; %搜索围theta0=0; %三个信号源的来波方向theta1=30;theta2=60;L=512;%采样单元数for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=10*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -20 0]);P=20*log10(abs(p)/(max(max(abs(p))))); subplot(1,2,2)plot(theta,P);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -20 0]);(2)仿真结果8.Capon beamforming(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=0.5;%阵元间距与波长的关系theta=-90:0.5:90; %搜索围theta0=0; %三个信号源的来波方向theta1=20;theta2=60;L=1000;%采样单元数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);amp2=200*randn(1);ampn=3;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';R=inv(Rx);steer=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');w=R*steer/(steer'*R*steer);%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p))))); subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -90 0]);(2)仿真结果9.不同方法估计协方差矩阵的Capon波束形成(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=0.5;%间距为半波长theta=-90:0.5:90;%扫描围theta0=0;%来波方向theta1=50;%干扰方向L=1024;%采样单元数for i=1:Lamp0=10*randn(1);amp1=50*randn(1);ampn=0.5;s(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]'); n(:,i)=ampn*exp(imag*2*pi*randn(1)*[0:element_num-1]');endRx=1/L*(s+j+n)*(s+j+n)';%接收信号自相关矩阵Rnj=1/L*(j+n)*(j+n)';%%干拢+噪声的自相关矩阵e=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f1(j)=Wopt_Rx'*a;f2(j)=Wopt_Rnj'*a;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 on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');xlabel('theta/°');ylabel('F(1,2)/dB');title('不同方法估计协方差矩阵的Capon波束形成');axis([-90 90 -60 0]);(2)仿真结果10.多点约束的Capon波束形成和方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;d_lambda=0.5;theta=-90:0.3:90;theta0=0;theta1=20;theta2=50;L=512;Rx=zeros(element_num,element_num);%产生协方差矩阵for i=1:Lamp0=10*randn(1);amp1=10*randn(1);amp2=50*randn(1);ampn=0.5*randn(1);%噪声的幅度随机产生,保证噪声与信号之间是不相关的j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]') +amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*exp (imag*2*pi*randn(1)*[0:element_num-1]');x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]') +j(:,i);%表示接收信号endRx=1/L*x*x';R=inv(Rx);w=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*ex p(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2* pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));figure;subplot(1,2,1);plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming方向图'); axis([-90 90 -50 0]);subplot(1,2,2);plot(theta,P);hold on;grid on;plot(theta0,-90:0,'.');plot(theta1,-90:0,'.');plot(theta2,-90:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming功率谱');(2)仿真结果11.自适应波束形成方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;c=3e8;f=5e8;lambda=c/f;d_lambda=0.5;theta=-90:0.5:90;theta0=0;theta1=45;theta2=60;L=2048;for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=100*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]');C=[steer1 steer2 steer3];F=[1 0 0]';w=inv(Rx)*C*(inv(C'*inv(Rx)*C))*F;for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('自适应波束形成方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p))))); subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('功率谱');axis([-90 90 -50 0]);(2)仿真结果(3)GUI界面。
自适应波束形成与Matlab程序代码注解
xlabel('theta/radian')
ylabel('amplitude/dB')
title([num2str(element_num)'阵元均匀线阵方向图','来波方向为'num2str(theta0*180/pi)'度']);
axis([-1.5 1.5 -50 0]);
(2)仿真结果
A.来波方向为0°
不归一化
归一化
B.来波方向为45°
不归一化
归一化
C.随着阵元数的增加,波束宽度变窄,分辨力提高,仿真图如下:
非归一化
归一化
不归一化
归一化
2.波束宽度与波达方向及阵元数的关系
(1)matlab 程序
clc;
clearall;
closeall;
imag=sqrt(-1);
element_num1=16;
psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num2*d));
beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num3*d));
end
gridon;
xlabel('theta/radian')
ylabel('amplitude/dB')
title([num2str(element_num)'阵元均匀线阵方向图','来波方向为'num2str(theta0*180/pi)'度']);
Capon-波束形成matlab仿真(附源代码)教学内容
Ca p o n - 波束形成ma t l a b 仿真(附源代码)Capon波束形成阵列N=16,信号o -30 ,干扰为! 60 , 2 19 , 3 45 ,干扰功率分别为:40dB,35dB,50dB。
Capon波束形成后的方向图和功率谱如下:-10-20-30E -40B0 1D0Capon beamForming \门;[7Capon beftinfonning 勺牛吃fo o o -7-6-9Co o o o o o12 3 4 5-BO-80-100 3 -60 40 '20 0 20 4Q &0tTieta>o-60 40 -20 0 20 40 60 00 100ttiefa/o为了比较接收数据直接估计噪声协方差矩阵和利用干扰+噪声估计协方差矩阵的Capon波束形成的差异,进行如下仿真:可以看出利用干扰+噪声估计协方差矩阵的方向图性能较优于接收数据直接估计噪声协方差矩阵的方向图。
代码:clc;clear all ;close all;ima=sqrt(-1);eleme nt_num=8; %阵元数d_lamda=1/2; %阵元间距与波长的关系theta=-90:0.5:90; %范围theta0=-30; %来波方向theta仁60; %干扰方向1theta2=19; %干扰方向2theta3=45; %干扰方向3■30買E-O-so-90-80 ^0 -40 -20 0 20 40 GO B0theta/O不冋方也佶计协方垫年阵的Capon泼束形成5_-L=1000; %采样单元数for i=1:L;amp0=10*ra ndn( 1);%信号的幅度随机产生,保证信号之间是不相关的amp1=100*ra ndn (1);%输入阵列的噪声amp2=sqrt(10A3.5)*ra ndn( 1);%输入阵列的噪声amp3=sqrt(10A5)*ra ndn⑴;%输入阵列的噪声amp n=3;% 噪声x(:,i)=amp0*exp(ima*2*pi*1/2*si n(theta0*pi/180)*[0:eleme nt_num-1]')+...amp1*exp(ima*2*pi*1/2*si n(theta1*pi/180)*[0:eleme nt_num-1]')+...amp2*exp(ima*2*pi*1/2*si n(theta2*pi/180)*[0:eleme nt_num-1]')+...amp3*exp(ima*2*pi*1/2*si n(theta3*pi/180)*[0:eleme nt_num-1]')+...amp n*(ra ndn (eleme nt_nu m,1)+ima*ra ndn (eleme nt_nu m,1));endRx=1/L*x* x';R=i nv(Rx);steer=exp(ima*2*pi*1/2*si n(theta0*pi/180)*[0:eleme nt_num-1]');w=R*steer/(steer'*R*steer);%Capon 最优权矢量for j=1:le ngth(theta);a=exp(ima*2*pi*d_lamda*si n( theta(j)*pi/180)*[0:eleme nt_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/(max(abs (f))));P=20*log10(abs(p)/(max(abs(p))));%此处是功率的对数形式figure;% subplot(121)plot(theta,F),grid on,hold onplot(theta0,-80:0,'.')plot(theta1,-80:0,'.')plot(theta2,-80:0,'.')plot(theta3,-80:0,'.') xlabel('theta/o'); ylabel('F/dB');title('Capon beamforming方向图') % axis([-90 90 -50 0]); % subplot(122)figure;plot(theta,P),grid on xlabel('theta/o'); ylabel('功率/dB'); title('C apon beamformi ng 功率谱')。
波束形成 Matlab程序教学文稿
波束形成M a t l a b程序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(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]'); for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]'); p(j)=w'*a;endfigure;plot(theta,abs(p)),grid onxlabel('theta/radian')ylabel('amplitude')title('8阵元均匀线阵方向图')当来波方向为45度时,仿真图如下:8阵元均匀线阵方向图如下,来波方向为0度,20log(dB)随着阵元数的增加,波束宽度变窄,分辨力提高:仿真图如下:2.波束宽度与波达方向及阵元数的关系clcclear allclose 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/180-asin(sin(theta(j)*pi/180)-lamda/(element_num1*d)); psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num2*d)); beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num3*d)); endfigure;plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid onxlabel('theta');ylabel('Width in radians')title('波束宽度与波达方向及阵元数的关系')仿真图如下:3. 当阵元间距 /2dλ> 时,会出现栅瓣,导致空间模糊。
Capon波束形成器的仿真
Capon 波束形成器的仿真一、原理1、波束形成的定义波束形成就是从传感器阵列重构源信号,然后通过增加期望信源的贡献和抑制掉干扰源来实现的。
实际上波束形成可以看成空域(或时域)滤波器。
2、波束形成器的最佳权向量的推导首先,假定一天线阵列中个阵元的接收信号向量为()n x ,权向量为12[,,,]M w w w =w ,由于阵列的输出是对阵元的接收信号向量在个阵元上个分量的加权和,所以阵列的输出可写作:*1()()()MHmm m n n w x n ===∑y w x再令空间远场期望信号为()d t ,J 个干扰信号为()j i t ,1,2,,j J =,每个阵元上的加性白噪声为()k n t ,他们都具有相同的方差2n σ。
在这些假设条件下,第k 个阵元上的接收信号可以表示为1()()()()()()j Jk k d k i j k j x t a d t a i t n t θθ==++∑111212()()()()()()(),(),()()()()J d i i J M M x t d t n t x t i t n t a a a i t n t x t θθθ⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎡⎤=+⎣⎦⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦ 矩阵形式如下:1()()()()()()()()j Jd i j j t t t a d t a i t n t θθ==+=++∑x As n由以上可知:N 个快拍的波束形成器的输出为()()(1,2,,)H y t w x t t N ==输出的平均功率为2211222222111111()()()111()()()()()j N N H t t N J NNH Hd j i t j t t P w y t w x t N N w a d t i t w a w n t N N N θθ========⎡⎤=++⎢⎥⎣⎦∑∑∑∑∑∑当N →∞时,上式变为{}{}{}{}22222221()()()()()()()()jHHHJHHdj i nj P w E y t w E x t x t w wRwE d t w a E i t w a wθθσ=====++∑为了保证来自方向d θ的期望信号的正确接收,并完全抑制其它J 个干扰,关于权向量的约束条件应为:()1H d w a θ=和()0j H i w a θ=加上约束条件后{}222()()n P w E d t w σ=+波束形成器最佳权向量的确定以叙述为:在置零约束条件的约束下,求满足{}{}2min ()min HwwE y t w Rw =的权向量w 。
自适应波束形成与Matlab程序代码注解
1.均匀线阵方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;%阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=linspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag);patterndBnorm=20*log10(patternmagnorm);figure(1)plot(theta*180/pi,patternmag);grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);hold on;figure(2)plot(theta,patterndBnorm,'r');grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);axis([ -50 0]);(2)仿真结果A.来波方向为0°不归一化归一化B.来波方向为45°不归一化归一化C.随着阵元数的增加,波束宽度变窄,分辨力提高,仿真图如下:非归一化归一化不归一化归一化2.波束宽度与波达方向及阵元数的关系(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num1=16;element_num2=128;element_num3=1024;lambda=;d=*lambda;theta=0::90;for j=1:length(theta)fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num1*d)); psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num2*d)); beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num3*d)); endfigureplot(theta,fai,'r',theta,psi,'b',theta,beta,'g');grid on;xlabel('theta');ylabel('width in radians');title('波束宽度与达波方向及阵元数目的关系');legend('N=16','N=128','N=1024');(2)仿真结果结果3. 当阵元间距时,会出现栅瓣,导致空间模糊(1)仿真结果非归一化归一化4. 类似于时域滤波,天线方向图是最优权的傅立叶变换(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;source_num=1;d_lambda=;theta=linspace(-pi/2,pi/2,200);theta0=0;w=exp(imag*2*pi*d_lambda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag);patterndBnorm=20*log10(patternmagnorm);figure(1)subplot(1,2,1);plot(theta,patterndBnorm);grid on;xlabel('theta/radian');ylabel('amplitude/dB');axis([ -50 0]);subplot(1,2,2);pfft=fftshift(fft(w,256));pfftmag=abs(pfft);pfftmagnorm=pfftmag/max(max(pfftmag)); pfftdB=20*log10(pfftmagnorm);pfftdBnorm=20*log10(pfftmagnorm);plot(linspace(-pi/2,pi/2,256),pfftdBnorm); grid on;xlabel('theta/radian');ylabel('FFT_amplitude/dB');axis([ -50 0]);(2)仿真结果5.最大信噪比准则方向图和功率谱(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=;%间距为半波长theta=-90::90;%扫描范围theta0=0;%来波方位theta1=20;%干扰方向L=512;%采样点数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);ampn=1; s(:,i)=amp0*exp(imag*2*pi**sin(theta0*pi/180)*[0:element_num-1]'); j(:,i)=amp1*exp(imag*2*pi**sin(theta1*pi/180)*[0:element_num-1]');n(:,i)=ampn*exp(randn(element_num,1)+imag*randn(element_num,1));endRs=1/L*s*s';%信号自相关矩阵Rnj=1/L*(j*j'+n*n'); %干扰+噪声的自相关矩阵[V,D]=eig(Rs,Rnj); %(Rs,Rnj)的广义特征值和特征向量[D,I]=sort(diag(D)); %特征向量排序Wopt=V(:,I(8));%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=Wopt'*a;p(j)=a'*Rs*a+a'*Rnj*a;endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-80:0,'.');plot(theta1,-80:0,'.');xlabel('theta/0');ylabel('F in dB');title('max-SNR 方向图');axis([-90 90 -80 0]);hold on;subplot(1,2,2);plot(theta,P,'r');grid on;xlabel('theta/0');ylabel('功率 in dB');title('max-SNR 功率谱');grid on;axis([-90 90 -80 0]);(2)仿真结果6.ASC旁瓣相消----MSE准则(1) matlab 程序clc;close all;clear all;imag=sqrt(-1);M=32;%辅助天线数目d_lambda=;%阵元间距theta0=-30;%来波方向theta1=60;%干扰方向L=512;%采样单元数s=zeros(1,512); %预划分一个区域for ii=1:Lamp0=1*randn(1);%信号的幅度随机产生,保证信号之间是不相关的amp1=200*randn(1);ampn=1;jam(:,ii)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:M-1]')+ampn*(randn(M,1)+imag*randn(M,1)); %干扰+噪声s(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180))+amp1*exp(imag*2*pi*d_lamb da*sin(theta1*pi/180))+ampn*(randn(1,1)+imag*randn(1,1));%接收信号(信号+干扰+噪声)s0(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180));endRx=1/L*jam*jam';r_xd=1/L*jam*s';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(imag*2*pi*d_lambda*sin(theta(jj))*[0:M-1]');f(jj)=Wopt'*a;endF=20*log10(abs(f)/max(max(abs(f))));figure(1)plot(theta*180/pi,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('MSE准则下的方向图');axis([-90 90 -50 0]);(2)仿真结果7.线性约束最小方差(LCMV)准则(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=;%阵元间距与波长的关系theta=-90::90; %搜索范围theta0=0; %三个信号源的来波方向theta1=30;theta2=60;L=512;%采样单元数for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=10*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*exp(imag*2*pi *d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/1 80)*[0:element_num-1]')+ampn*(randn(element_num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -20 0]);P=20*log10(abs(p)/(max(max(abs(p)))));subplot(1,2,2)plot(theta,P);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -20 0]);(2)仿真结果8.Capon beamforming(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=;%阵元间距与波长的关系theta=-90::90; %搜索范围theta0=0; %三个信号源的来波方向theta1=20;theta2=60;L=1000;%采样单元数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);amp2=200*randn(1);ampn=3;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*exp(imag*2*pi *d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/1 80)*[0:element_num-1]')+ampn*(randn(element_num,1)+imag*randn(element_num,1));endRx=1/L*x*x';R=inv(Rx);steer=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');w=R*steer/(steer'*R*steer);%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]'); f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p)))));subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -90 0]);(2)仿真结果9.不同方法估计协方差矩阵的Capon波束形成(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=;%间距为半波长theta=-90::90;%扫描范围theta0=0;%来波方向theta1=50;%干扰方向L=1024;%采样单元数for i=1:Lamp0=10*randn(1);amp1=50*randn(1);ampn=;s(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]'); n(:,i)=ampn*exp(imag*2*pi*randn(1)*[0:element_num-1]');endRx=1/L*(s+j+n)*(s+j+n)';%接收信号自相关矩阵Rnj=1/L*(j+n)*(j+n)';%%干拢+噪声的自相关矩阵e=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f1(j)=Wopt_Rx'*a;f2(j)=Wopt_Rnj'*a;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 on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');xlabel('theta/°');ylabel('F(1,2)/dB');title('不同方法估计协方差矩阵的Capon波束形成'); axis([-90 90 -60 0]);(2)仿真结果10.多点约束的Capon波束形成和方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;theta=-90::90;theta0=0;theta1=20;theta2=50;L=512;Rx=zeros(element_num,element_num);%产生协方差矩阵for i=1:Lamp0=10*randn(1);amp1=10*randn(1);amp2=50*randn(1);ampn=*randn(1);%噪声的幅度随机产生,保证噪声与信号之间是不相关的j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag* 2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*exp(imag*2*pi*randn(1)*[0:element_ num-1]');x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+j(:,i);%表示接收信号endRx=1/L*x*x';R=inv(Rx);w=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*exp(imag*2*pi*d_la mbda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[ 0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));figure;subplot(1,2,1);plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming方向图');axis([-90 90 -50 0]);subplot(1,2,2);plot(theta,P);grid on;plot(theta0,-90:0,'.');plot(theta1,-90:0,'.');plot(theta2,-90:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming功率谱');(2)仿真结果11.自适应波束形成方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;c=3e8;f=5e8;lambda=c/f;d_lambda=;theta=-90::90;theta0=0;theta1=45;theta2=60;L=2048;for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=100*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*exp(imag*2*pi *d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/1 80)*[0:element_num-1]')+ampn*(randn(element_num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]');C=[steer1 steer2 steer3];F=[1 0 0]';w=inv(Rx)*C*(inv(C'*inv(Rx)*C))*F;for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('自适应波束形成方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p)))));subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('功率谱');axis([-90 90 -50 0]);(2)仿真结果(3)GUI界面。
Matlab仿真相控阵波束扫描的“动图”(附代码)
Matlab仿真相控阵波束扫描的“动图”(附代码)在《无源/有源相控阵雷达、电扫阵列及典型系统》中有下面这样一幅动图,有人问是怎么画出来的。
今天我们就来分析一下,并且手把手教大家动手画出来。
动图分析:我们观察到这是一个极坐标下线阵的波束扫描图,其中外圈标注的是度数0~360°,0°表示垂直于线阵的方向。
其中,一圈一圈的环代表的是幅度,标出的100,80,60...是dB值。
波束增益最大的是指向0°的方向,最大值100dB,旁瓣比主瓣低接近40dB,说明进行了降低旁瓣的加权处理,并且波束在扫描到30多度后出现了删瓣。
动图仿真:这里我们做了一定的简化,首先,没有进行降低旁瓣的加权。
幅度进行了归一化,也就是指向0°时是0dB。
仿真时设置:N=32,d=λ/2,波束指向θB从0~60°,再从60~0°,得到下面的动图:从该动图可以看出:1. 随着波束指向的增大,波束宽度变胖,当达到60°时,变胖了一倍。
2. 随着波束指向的增大,增益也有降低。
3. 没有进行降低旁瓣的加权,旁瓣较高。
几个知识点科普:波束宽度与扫描角θB的关系:当扫描的最大角度为θmax时,为了不出现删瓣,阵元间距d和波长λ需要满足关系:也就是说当阵元间距小于半波长时,即使扫描到90°都不会出现删瓣。
如果你想看看出现删瓣的情况,这里也有(设置d=0.7λ)。
有一点需要注意,这里仿真图中90~270°内出现波束是因为这里仿真的每个阵元都是全向天线。
仿真程序如下:Matlab中有函数polarplot可以直接使用,但不是dB表示的。
需要polarplotdb函数以及上面仿真程序的可以发送“180724”直接下载。
自适应波束形成与Matlab程序代码注解讲诉
1.均匀线阵方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;%阵元数为8d_lamda=1/2;%阵元间距d与波长lamda的关系theta=linspace(-pi/2,pi/2,200);theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag);patterndBnorm=20*log10(patternmagnorm);figure(1)plot(theta*180/pi,patternmag);grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);hold on;figure(2)plot(theta,patterndBnorm,'r');grid on;xlabel('theta/radian')ylabel('amplitude/dB')title([num2str(element_num) '阵元均匀线阵方向图','来波方向为' num2str(theta0*180/pi) '度']);axis([-1.5 1.5 -50 0]);(2)仿真结果A.来波方向为0°不归一化归一化B.来波方向为45°不归一化归一化C.随着阵元数的增加,波束宽度变窄,分辨力提高,仿真图如下:非归一化归一化不归一化归一化2.波束宽度与波达方向及阵元数的关系(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num1=16;element_num2=128;element_num3=1024;lambda=0.1;d=0.5*lambda;theta=0:0.5:90;for j=1:length(theta)fai(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num1*d)); psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num2*d)); beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lambda/(element_num3*d)); endfigureplot(theta,fai,'r',theta,psi,'b',theta,beta,'g');grid on;xlabel('theta');ylabel('width in radians');title('波束宽度与达波方向及阵元数目的关系');legend('N=16','N=128','N=1024');(2)仿真结果结果3. 当阵元间距/2d λ>时,会出现栅瓣,导致空间模糊(1)仿真结果非归一化归一化4. 类似于时域滤波,天线方向图是最优权的傅立叶变换(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=32;source_num=1;d_lambda=0.5;theta=linspace(-pi/2,pi/2,200);theta0=0;w=exp(imag*2*pi*d_lambda*sin(theta0)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j))*[0:element_num-1]'); p(j)=w'*a;endpatternmag=abs(p);patternmagnorm=patternmag/max(max(patternmag));patterndB=20*log10(patternmag); patterndBnorm=20*log10(patternmagnorm); figure(1)subplot(1,2,1);plot(theta,patterndBnorm);grid on;xlabel('theta/radian');ylabel('amplitude/dB');axis([-2.0 2.0 -50 0]);subplot(1,2,2);pfft=fftshift(fft(w,256));pfftmag=abs(pfft);pfftmagnorm=pfftmag/max(max(pfftmag)); pfftdB=20*log10(pfftmagnorm);pfftdBnorm=20*log10(pfftmagnorm);plot(linspace(-pi/2,pi/2,256),pfftdBnorm); grid on;xlabel('theta/radian');ylabel('FFT_amplitude/dB');axis([-2.0 2.0 -50 0]);(2)仿真结果5.最大信噪比准则方向图和功率谱(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=0.5;%间距为半波长theta=-90:0.5:90;%扫描范围theta0=0;%来波方位theta1=20;%干扰方向L=512;%采样点数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);ampn=1;s(:,i)=amp0*exp(imag*2*pi*0.5*sin(theta0*pi/180)*[0:element_num-1]');j(:,i)=amp1*exp(imag*2*pi*0.5*sin(theta1*pi/180)*[0:element_num-1]');n(:,i)=ampn*exp(randn(element_num,1)+imag*randn(element_num,1)); endRs=1/L*s*s';%信号自相关矩阵Rnj=1/L*(j*j'+n*n'); %干扰+噪声的自相关矩阵[V,D]=eig(Rs,Rnj); %(Rs,Rnj)的广义特征值和特征向量[D,I]=sort(diag(D)); %特征向量排序Wopt=V(:,I(8));%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=Wopt'*a;p(j)=a'*Rs*a+a'*Rnj*a;endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-80:0,'.'); plot(theta1,-80:0,'.'); xlabel('theta/0'); ylabel('F in dB');title('max-SNR 方向图'); axis([-90 90 -80 0]); hold on;subplot(1,2,2);plot(theta,P,'r');grid on;xlabel('theta/0'); ylabel('功率 in dB'); title('max-SNR 功率谱'); grid on;axis([-90 90 -80 0]); (2)仿真结果6.ASC旁瓣相消----MSE准则(1) matlab 程序clc;close all;clear all;imag=sqrt(-1);M=32;%辅助天线数目d_lambda=0.5;%阵元间距theta0=-30;%来波方向theta1=60;%干扰方向L=512;%采样单元数s=zeros(1,512); %预划分一个区域for ii=1:Lamp0=1*randn(1);%信号的幅度随机产生,保证信号之间是不相关的amp1=200*randn(1);ampn=1;jam(:,ii)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:M-1]' )+ampn*(randn(M,1)+imag*randn(M,1)); %干扰+噪声s(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180))+amp1*exp(im ag*2*pi*d_lambda*sin(theta1*pi/180))+ampn*(randn(1,1)+imag*randn( 1,1));%接收信号(信号+干扰+噪声)s0(ii)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180));endRx=1/L*jam*jam';r_xd=1/L*jam*s';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(imag*2*pi*d_lambda*sin(theta(jj))*[0:M-1]'); f(jj)=Wopt'*a;endF=20*log10(abs(f)/max(max(abs(f))));figure(1)plot(theta*180/pi,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('MSE准则下的方向图');axis([-90 90 -50 0]);(2)仿真结果7.线性约束最小方差(LCMV)准则(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=0.5;%阵元间距与波长的关系theta=-90:0.5:90; %搜索范围theta0=0; %三个信号源的来波方向theta1=30;theta2=60;L=512;%采样单元数for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=10*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -20 0]);P=20*log10(abs(p)/(max(max(abs(p)))));subplot(1,2,2)plot(theta,P);grid on;hold on;plot(theta0,-20:0,'.');plot(theta1,-20:0,'.');plot(theta2,-20:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -20 0]);(2)仿真结果8.Capon beamforming(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数d_lambda=0.5;%阵元间距与波长的关系theta=-90:0.5:90; %搜索范围theta0=0; %三个信号源的来波方向theta1=20;theta2=60;L=1000;%采样单元数for i=1:Lamp0=10*randn(1);amp1=200*randn(1);amp2=200*randn(1);ampn=3;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';R=inv(Rx);steer=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');w=R*steer/(steer'*R*steer);%最优权矢量for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming 方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p)))));subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('Capon beamforming 功率谱');axis([-90 90 -90 0]);(2)仿真结果9.不同方法估计协方差矩阵的Capon波束形成(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;%阵元数为8d_lambda=0.5;%间距为半波长theta=-90:0.5:90;%扫描范围theta0=0;%来波方向theta1=50;%干扰方向L=1024;%采样单元数for i=1:Lamp0=10*randn(1);amp1=50*randn(1);ampn=0.5;s(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]'); n(:,i)=ampn*exp(imag*2*pi*randn(1)*[0:element_num-1]');endRx=1/L*(s+j+n)*(s+j+n)';%接收信号自相关矩阵Rnj=1/L*(j+n)*(j+n)';%%干拢+噪声的自相关矩阵e=exp(imag*2*pi*d_lambda*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(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f1(j)=Wopt_Rx'*a;f2(j)=Wopt_Rnj'*a;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 on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');xlabel('theta/°');ylabel('F(1,2)/dB');title('不同方法估计协方差矩阵的Capon波束形成');axis([-90 90 -60 0]);(2)仿真结果10.多点约束的Capon波束形成和方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;d_lambda=0.5;theta=-90:0.3:90;theta0=0;theta1=20;theta2=50;L=512;Rx=zeros(element_num,element_num);%产生协方差矩阵for i=1:Lamp0=10*randn(1);amp1=10*randn(1);amp2=50*randn(1);ampn=0.5*randn(1);%噪声的幅度随机产生,保证噪声与信号之间是不相关的j(:,i)=amp1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]') +amp2*exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*exp (imag*2*pi*randn(1)*[0:element_num-1]');x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]') +j(:,i);%表示接收信号endRx=1/L*x*x';R=inv(Rx);w=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+amp1*ex p(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(imag*2* pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]');for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*R*a);endF=20*log10(abs(f)/max(max(abs(f))));P=20*log10(abs(p)/max(max(abs(p))));figure;subplot(1,2,1);plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('Capon beamforming方向图'); axis([-90 90 -50 0]);subplot(1,2,2);plot(theta,P);hold on;grid on;plot(theta0,-90:0,'.');plot(theta1,-90:0,'.');plot(theta2,-90:0,'.');xlabel('theta/°');ylabel('P/dB');title('Capon beamforming功率谱');(2)仿真结果11.自适应波束形成方向图(1)matlab 程序clc;clear all;close all;imag=sqrt(-1);element_num=8;c=3e8;f=5e8;lambda=c/f;d_lambda=0.5;theta=-90:0.5:90;theta0=0;theta1=45;theta2=60;L=2048;for i=1:Lamp0=10*randn(1);amp1=100*randn(1);amp2=100*randn(1);ampn=10;x(:,i)=amp0*exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]')+am p1*exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]')+amp2*exp(im ag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]')+ampn*(randn(element_ num,1)+imag*randn(element_num,1));endRx=1/L*x*x';steer1=exp(imag*2*pi*d_lambda*sin(theta0*pi/180)*[0:element_num-1]');steer2=exp(imag*2*pi*d_lambda*sin(theta1*pi/180)*[0:element_num-1]');steer3=exp(imag*2*pi*d_lambda*sin(theta2*pi/180)*[0:element_num-1]');C=[steer1 steer2 steer3];F=[1 0 0]';w=inv(Rx)*C*(inv(C'*inv(Rx)*C))*F;for j=1:length(theta)a=exp(imag*2*pi*d_lambda*sin(theta(j)*pi/180)*[0:element_num-1]');f(j)=w'*a;p(j)=1/(a'*inv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs(f)))));subplot(1,2,1)plot(theta,F);grid on;hold on;plot(theta0,-50:0,'.');plot(theta1,-50:0,'.');plot(theta2,-50:0,'.');xlabel('theta/°');ylabel('F/dB');title('自适应波束形成方向图');axis([-90 90 -50 0]);P=20*log10(abs(p)/(max(max(abs(p))))); subplot(1,2,2)plot(theta,P);grid on;hold on;xlabel('theta/°');ylabel('P/dB');title('功率谱');axis([-90 90 -50 0]);(2)仿真结果(3)GUI界面。
常规波束形成matlab
常规的波束形成方法通常使用MATLAB实现。
以下是一个简单的MATLAB代码示例,演示了如何使用常规的波束形成方法:```matlab设置参数c = 343; 光速,单位为m/sf = 1000; 信号频率,单位为HzN = 128; 阵元数目d = c/f; 阵元间距,单位为mtheta = 1:1:180; 方向角,单位为度t = 0:1/f:1-1/f; 时间,单位为sA = zeros(1,length(theta)); 信号幅度A(5) = 5; 在第5个方向上有目标W = zeros(N,length(theta)); 阵列流形计算阵列流形for i = 1:1:Nfor j = 1:length(theta)W(i,j) = exp(-2*pi*f*(i-1)*d*sind(theta(j))/c);endend计算波束形成输出S = zeros(length(theta),length(t));for j = 1:length(theta)S(j,:) = sin(2*pi*f.*t).*A(j);endS = W*S;画图figure;plot(theta,abs(S));xlabel('方向角(度)');ylabel('幅度');title('波束形成输出');```这段代码使用了常规的波束形成方法,首先设定了光速、信号频率、阵元数目、阵元间距、方向角和时间等参数。
然后,根据给定的信号幅度和阵列流形计算波束形成输出。
最后,使用MATLAB的绘图功能将波束形成输出绘制成图形。
波束形成-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=45/180*pi;%来波方向(我觉得应该是天线阵的指向)%theta0=0;%来波方向w=exp(imag*2*pi*d_lamda*sin(theta0)*[0:element_num-1]');for j=1:length(theta) %(我认为是入射角度,即来波方向,计算阵列流形矩阵A) a=exp(imag*2*pi*d_lamda*sin(theta(j))*[0:element_num-1]');p(j)=w'*a; %(matlab中的'默认为共轭转置,如果要计算转置为w.'*a)endfigure;plot(theta,abs(p)),grid onxlabel('theta/radian')ylabel('amplitude')title('8阵元均匀线阵方向图')见张小飞的书《阵列信号处理的理论和应用2.3.4节阵列的方向图》当来波方向为45度时,仿真图如下:8阵元均匀线阵方向图如下,来波方向为0度,20log(dB)随着阵元数的增加,波束宽度变窄,分辨力提高:仿真图如下:2.波束宽度与波达方向及阵元数的关系clcclear allclose 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/180-asin(sin(theta(j)*pi/180)-lamda/(element_num1*d));psi(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num2*d));beta(j)=theta(j)*pi/180-asin(sin(theta(j)*pi/180)-lamda/(element_num3*d)); endfigure;plot(theta,fai,'r',theta,psi,'b',theta,beta,'g'),grid onxlabel('theta');ylabel('Width in radians')title('波束宽度与波达方向及阵元数的关系')仿真图如下:3. 当阵元间距 /2d λ> 时,会出现栅瓣,导致空间模糊。
CAPON波束形成-Matlab程序教学文稿
C A P O N波束形成-M a t l a b程序CAPON波束形成器仿真1.实验原理波束形成就是从传感器阵列重构源信号。
(1)、通过增加期望信源的贡献来实现;(2)、通过抑制掉干扰源来实现。
经典的波束形成需要观测方向(期望信源的方向)的知识。
盲波束形成试图在没有期望信源方向信息的情况下进行信源的恢复。
波束形成技术的基本思想是:通过将各阵元输出进行加权求和,在一时间内将天线阵列波束“导向”到一个方向上,对期望信号得到最大输出功率的导向位置即给出DOA估计。
虽然阵列天线的方向图是全方向的,但阵列的输出经过加权求和后,却可以被调整到阵列接收的方向增益聚集在一个方向上,相当于形成了一个”波束”。
这就是波束形成的物理意义所在。
在智能天线中,波束形成是关键技术之一,是提高信噪比、增加用户容量的保证,能够成倍地提高通信系统的容量,有效地抑制各种干扰,并改善通信质量。
波束形成器的最佳权向量w取决于阵列方向向量)(ak,而在移动通信里用户的方向向量一般未知,需要估计(称之为DOA估计)。
因此,在计算波束形成的最佳权向量之前,必须在已知阵列几何结构的前提下先估计期望信号的波达方向。
Capon 波束形成器求解的优化问题可表述为w arg min P(w)θ=其约束条件为 1)(a w H =θCapon 波束形成器在使噪声和干扰所贡献的功率为最小的同时,保持了期 望信号的功率不变。
因此,它可以看作是一个尖锐的空间带通滤波器。
最优加 权向量w 可以利用Lagrange 乘子法求解,其结果为 )(a R ˆ)(a )(a R ˆw 1H 1CAP θθθ--= 当μ不取常数,而取作)(a R ˆ)(a 11H θθμ-=时,最佳权向量就转变成Capon 波束形成器的权向量。
空间谱为)(a R ˆ)(a 1)(P 1-H CAP θθθ=2.变量定义 M :均匀线阵列数目P :信号源个数nn :快拍数angle1、angle2、angle3:信号来波角度u :复高斯噪声Ps :信号能量refp :信噪比(实值)X:接收信号Rxx:接收信号的相关矩阵doa:波达方向估计3.仿真结果采用上述算法进行仿真,结果如图所示。
基于Matlab的Capon波束形成和空域谱算法仿真分析
基于Matlab的Capon波束形成和空域谱算法仿真分析毛海波
【期刊名称】《计算机仿真》
【年(卷),期】2005(000)0z1
【摘要】采用Matlab语言分别仿真了Capon波束形成,Bartlett和Music空域谱估计算法,仿真结果展示了阵列天线参数与天线波束图之间的关系,使复杂的天线参数直观化,算法的理论结果可视化.通过改变信号和噪声的方向,幅度以及调整阵列天线的参数,对比前后的天线波束图仿真结果的变化.通过分析仿真的结果,指出了这几种算法上的优缺点,为将来算法的改进和适用范围的扩大以及新的算法的提出打下了铺垫.
【总页数】4页(P491-494)
【作者】毛海波
【作者单位】北京航天自动控制研究所,北京,100854
【正文语种】中文
【中图分类】TP391.9
【相关文献】
1.基于协方差矩阵估计的稳健Capon波束形成算法 [J], 陈明建;罗景青;龙国庆
2.基于鲁棒Capon波束形成的探地雷达成像算法 [J], 周琳;黄春琳;粟毅
3.基于协方差矩阵重构的稳健Capon波束形成算法 [J], 季科;王洪雁;裴炳南
4.基于迭代梯度方法的线性约束稳健Capon波束形成快速算法 [J], 郭翔宇;鄢社锋;王文侠
5.基于干扰噪声矩阵估计的稳健Capon波束形成算法 [J], 徐定杰;刘向锋;沈锋;韩浩
因版权原因,仅展示原文概要,查看原文内容请购买。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
C a p o n-波束形成m a t l a b仿真(附源代
码)
Capon波束形成
阵列N=16,信号
0-30
θ︒
=,干扰为
160
θ︒=,
219
θ︒
=,
345
θ︒
=,干扰功率分别为:40dB,35dB,50dB。
Capon波束形成后的方向图和功率谱如下:
为了比较接收数据直接估计噪声协方差矩阵和利用干扰+噪声估计协方差矩阵的Capon波束形成的差异,进行如下仿真:
可以看出利用干扰+噪声估计协方差矩阵的方向图性能较优于接收数据直接估计噪声协方差矩阵的方向图。
代码:
clc;
clear all ;
close all;
ima=sqrt(-1);
element_num=8; %阵元数
d_lamda=1/2; %阵元间距与波长的关系
theta=-90:0.5:90; %范围
theta0=-30; %来波方向
theta1=60; %干扰方向1
theta2=19; %干扰方向2
theta3=45; %干扰方向3
L=1000; %采样单元数
for i=1:L;
amp0=10*randn(1);%信号的幅度随机产生,保证信号之间是不相关的 amp1=100*randn(1);%输入阵列的噪声
amp2=sqrt(10^3.5)*randn(1);%输入阵列的噪声
amp3=sqrt(10^5)*randn(1);%输入阵列的噪声
ampn=3;%噪声
x(:,i)=amp0*exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]')+... amp1*exp(ima*2*pi*1/2*sin(theta1*pi/180)*[0:element_num-1]')+...
amp2*exp(ima*2*pi*1/2*sin(theta2*pi/180)*[0:element_num-1]')+...
amp3*exp(ima*2*pi*1/2*sin(theta3*pi/180)*[0:element_num-1]')+...
ampn*(randn(element_num,1)+ima*randn(element_num,1));
end
Rx=1/L*x* x';
R=inv(Rx);
steer=exp(ima*2*pi*1/2*sin(theta0*pi/180)*[0:element_num-1]');
w=R*steer/(steer'*R*steer);%Capon最优权矢量
for j=1:length(theta);
a=exp(ima*2*pi*d_lamda*sin(theta(j)*pi/180)*[0:element_num-1]');
f(j)=w'*a;
p(j)=1/(a'*R*a);
end
F=20*log10(abs(f)/(max(abs(f))));
P=20*log10(abs(p)/(max(abs(p))));%此处是功率的对数形式
figure;
% subplot(121)
plot(theta,F),grid on,hold on
plot(theta0,-80:0,'.')
plot(theta1,-80:0,'.')
plot(theta2,-80:0,'.')
plot(theta3,-80:0,'.')
xlabel('theta/o');
ylabel('F/dB');
title('Capon beamforming方向图') % axis([-90 90 -50 0]);
% subplot(122)
figure;
plot(theta,P),grid on
xlabel('theta/o');
ylabel('功率/dB');
title('Capon beamforming功率谱')。