波束形成Matlab程序

1、實(shí)用文檔1.均勻線陣方向圖%8 陣元均勻線陣方向圖，來波方向?yàn)? 度clc;clear all;close all;imag=sqrt(_1);element_num=8;% 陣元數(shù)為 8d_lamda=1/2;%陣元間距 d 與波長 lamda 的關(guān)系theta=li nspace(-pi/2,pi/2,200);theta0=45/180*pi;%來波方向(我覺得應(yīng)該是天線陣的指向)%theta0=0;%來波方向w=exp(imag*2*pi*dl_lamda*si n(theta0)*0:eleme nt_nu m-1);for j=1:le ngth(theta) %(我認(rèn)為是入射角度

2、，即來波方向，計(jì)算陣列流形矩陣 A)a=exp(imag*2*pi*damda*si n(theta(j)*0:eleme nt_nu m-1);p(j)=w*a; %(matlab中的默認(rèn)為共軛轉(zhuǎn)置，如果要計(jì)算轉(zhuǎn)置為w.*a)endfigure;plot(theta,abs(p),grid onxlabel(theta/radia n)ylabel(amplitude)title(8陣元均勻線陣方向圖)見小飛的書陣列信號(hào)處理的理論和應(yīng)用2.3.4 節(jié)陣列的方向圖 -T-1-0.500511.52thetaJradian實(shí)用文檔swsmmEE當(dāng)來波方向?yàn)?5度時(shí)，仿真圖如下:8陣元均勻絨陣才向

3、團(tuán).來波方向拘4G度8 陣元均勻線陣萬向圖如下，來波萬向?yàn)? 度，20log (dB)日陣元均勻線陣方向囹，束披方向?yàn)?度31-1 -a s n 0斤iheld/i adidi i11桿iheta/ractian?oJI一實(shí)用文檔隨著陣元數(shù)的增加，波束寬度變窄，分辨力提高：仿真圖如下:實(shí)用文檔實(shí)用文檔012.波束寬度與波達(dá)萬向及陣元數(shù)的關(guān)系clcclear allclose allima=sqrt(_1);eleme nt_nu m1=16; % 陣元數(shù)eleme nt_n um2=128;eleme nt_nu m3=1024;lamda=0.03;% 波長為 0.03 米d=1/2*lam

4、da; %陣元間距與波長的關(guān)系theta=0:0.5:90;for j=1:le ngth(theta);fai(j)=theta(j)*pi/180-asi n(si n( theta(j)*pi/180)-lamda/(eleme nt_num1*d);psi(j)=theta(j)*pi/180-asi n(s in (theta(j)*pi/180)-lamda/(element_num2*d);beta(j)=theta(j)*pi/180-asi n(si n( theta(j)*pi/180)-lamda/(eleme nt_num3*d); end figure;plot(th

5、eta,fai,r,theta,psi,b,theta,beta,g),grid onxlabel(theta);ylabel(Width in radia ns)title(波束寬度與波達(dá)方向及陣元數(shù)的關(guān)系)仿真圖如下：液束寬度與浪達(dá)方向及陣元載的關(guān)系07芒逹翱=-1實(shí)用文檔010.690908070705Q5Q4040he實(shí)用文檔3.當(dāng)陣元間距d/2時(shí)，會(huì)出現(xiàn)柵瓣，導(dǎo)致空間模糊。仿真圖如下:4.類似于時(shí)域?yàn)V波，天線方向圖是最優(yōu)權(quán)的傅立葉變換 仿真程序和仿真圖如下：clcclear allclose allima=sqrt(-1);eleme nt_num=32;%陣元數(shù)實(shí)用文檔source

6、_ num=1;%信源數(shù)d_lamda=1/2; %陣元間距與波長的關(guān)系theta=li nspace(-pi/2,pi/2,200);thetaO=O; %來波方向(ayy應(yīng)該是陣列指向方向)w=exp(ima*2*pi*d_lamda*s in (thetaO)*O:eleme nt_nu m-1); for j=1:length(theta);a=exp(ima*2*pi*d_lamda*s in (theta(j)*O:eleme nt_nu m-1); p(j)=w*a;endfigure;subplot(1,2,1)plot(theta,abs(p),grid on xlabel(

7、theta/radia n) ylabel(amplitude) title(按定義的方向圖) pfft=fftshift(fft(w,128);subplot(1,2,2)plot(li nspace(-pi/2,pi/2,128),abs(pfft),grid on xlabel(theta/radia n)ylabel(FFT_amplitude)title(最優(yōu)權(quán)的傅里葉變換)實(shí)用文檔5. %最大信噪比準(zhǔn)則方向圖和功率譜clc;clear all;close all;ima=sqrt(-1);element_num=8;%陣元數(shù)為8d_lamda=1/2; %間距為半波長theta=-

8、90:0.5:90;%圍thetaO=O; %來波方向theta1=20; %干擾方向L=512;%采樣單元數(shù)for i=1:Lamp0=10*ra ndn (1);amp1=200*ra ndn (1);amp n=1;s(:,i)=amp0*exp(ima*2*pi*1/2*si n(theta0*pi/180)*0:element_num-1);j(:,i)=amp1*exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1);n (:,i)=amp n*(ra ndn (eleme nt_nu m,1)+ima*ra ndn (elem

9、e nt_num,1);endRs=1/L*s*s;%信號(hào)自相關(guān)矩陣Rnj=1/L*(j*j+n*n);%干擾+噪聲的自相關(guān)矩陣V,D=eig(Rs,Rnj); %（Rs,Rnj）的廣義特征值和特征向量按定裂的方甸劉theta-radifsnthet/radian實(shí)用文檔D,l=sort(diag(D); %排序Wopt=V(:,l(8);%最優(yōu)權(quán)矢量實(shí)用文檔-36-00for j=1:le ngth(theta)a=exp(ima*2*pi*d_lamda*si n( theta(j)*pi/180)*0:element_num-1); f(j)=Wopt*a;p(j)=a*Rs*a+a*R

10、nj*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 indB);title(max-SNR方向圖);axis(-90 90 -50 0);hold onsubplot(122)plot(theta,P,r);grid onxlabel(theta/0);ylabel(功率in dB);

11、title(max-SNR功率譜)仿真圖如下：-20-30-4C-4f) )CD70% 編赳-100-50O50100ttietEI/0實(shí)用文檔6. %AS（旁瓣相消-MSE準(zhǔn)貝Ufor ii=1:L;amp0=1*randn(1); %amp1=200*ra ndn (1);amp n=1;jam(:,ii)=amp1*exp(ima*2*pi*0.5*si n(theta1*pi/180)*0:M-1)+amp n*(randn( M,1)+ima*ra ndn( M,1); %干擾+噪聲s(ii)=amp0*exp(ima*2*pi*0.5*si n(theta0*pi/180)+amp

12、1*exp(ima*2*pi*0.5*sin(theta1*pi/180)+ampn*(randn(1,1)+ima*randn(1,1);%接收信號(hào)(信號(hào)+干擾+噪聲)s0(ii)=amp0*exp(ima*2*pi*0.5*si n(theta0*pi/180);endRx=1/L*jam*jam; %噪聲自相關(guān)矩陣，相當(dāng)于X(t)r_xd=1/L*jam*s;Wopt=p in v(Rx)*r_xd;delta=s0-(s-Wopt*jam);delta仁abs(mea n(delta42)-(mea n( delta).A2) %方差theta=li nspace(-pi/2,pi/2

13、,200);for jj=1:le ngth(theta)a=exp(ima*2*pi*.5*s in (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 onplot(theta0,-50:0,.)plot(theta1,-50:0,.)xlabel(theta/o);ylabel(F/dB);clc;close all;clear allima=sqrt(-1);M=32;%d_lamda=.5;theta0=-30;%

14、theta1=60;%L=512;%輔助天線的數(shù)目s=zeros(1,512);%來波方向干擾方向采樣單元數(shù)預(yù)劃分一個(gè)區(qū)域信號(hào)的幅度隨機(jī)產(chǎn)生，保證信號(hào)之間是不相關(guān)的實(shí)用文檔title(MSE準(zhǔn)則下的方向圖)實(shí)用文檔axis(-90 90 -50 0);%仿真圖如下：可為x軸和y軸設(shè)置一個(gè)極限圍實(shí)用文檔7. %線性約束最小方差(LCMV準(zhǔn)則clc;clear all ;close all;ima=sqrt(-1);eleme nt_num=8;%陣元數(shù)d_lamda=1/2; %陣元間距與波長的關(guān)系theta=-90:0.5:90-0.3;%搜索圍thetaO=O;%三個(gè)信號(hào)源的來波方向thet

15、a1=30;theta2=60;L=512;%采樣單元數(shù)for i=1:L;amp0=10*ra ndn (1);amp 1=100*ra ndn (1);amp2=10*ra ndn (1);amp n=10;x(:,i)=amp0*exp(ima*2*pi*1/2*si n(theta0*pi/180)*0:element_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)+.amp n

16、*(ra ndn (eleme nt_nu m,1)+ima*ra ndn (eleme nt_nu m,1);endRx=1/L*x* x;steer仁exp(ima*2*pi*1/2*si n(theta0*pi/180)*0:eleme nt_num-1)steer2=exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1)steer3=exp(ima*2*pi*1/2*si n(theta2*pi/180)*0:eleme nt_num-1)C=steer1 steer2 steer3;F=1 0 1;%把三個(gè)方向都作為來波方向w=i

17、 nv(Rx)*C*(i nv(C* in v(Rx)*C)*F;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*i nv(Rx)*a);endf=10*log10(abs(f)/(max(max(abs (f);figure(1)subplot(121)plot(theta,f),grid on,hold on實(shí)用文檔plot(theta0,-20:0,.)plot(theta1,-20:0,.)plot(theta2,-20:0

18、,.) xlabel(theta/o);ylabel(F/dB);title(Capon beamforming方向圖)axis(-90 90 -20 0);%可為x軸和y軸設(shè)置一個(gè)極限圍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(C apon beamformi ng功率譜)仿真圖如下：實(shí)用文檔8

19、.%Capon beamformingClc;clear all ; close all;ima=sqrt(-1);eleme nt_num=8;%damda=1/2;%theta=-90:0.5:90;%thetaO=O;%theta1=20;%theta2=60;%L=1000;%for i=1:L;陣元數(shù)陣元間距與波長的關(guān)系 圍來波方向干擾方向干擾方向采樣單元數(shù)4-S-W-12Capon beamforming向SHP.博-14-ie-203%oC叩口口beamfiorming功率譜5-60050theta.-Q100實(shí)用文檔ampO=1O*ra ndn;信號(hào)的幅度隨機(jī)產(chǎn)生，保證信號(hào)之間

20、是不相關(guān)的amp1=200*ra ndn (1);amp2=200*ra ndn (1);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)+. amp n*(ra ndn(eleme nt_nu m,1)+ima*ra ndn (eleme nt_nu m

21、,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);%Capo n最優(yōu)權(quán)矢量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(max(abs (f);P=20*log10(abs(p)/(ma

22、x(max(abs(p);%此處是功率的對(duì)數(shù)形式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(C apon beamformi ng功率譜)

23、仿真圖如下：實(shí)用文檔9. %不同方法估計(jì)協(xié)方差矩陣的Capon波束形成clc;clear all;close all;ima=sqrt(-1);Capon b-eamformingT向圖theta/oCapon b日mfoimirig功率譜實(shí)用文檔element_num=8;%陣元數(shù)為8d_lamda=1/2; %間距為半波長theta=-90:0.5:90;%圍thetaO=O;%來波方向theta1=50;%干擾方向L=512;%采樣單元數(shù)for i=1:Lamp0=10*ra ndn (1);amp1=50*ra ndn (1);amp n=0.5;s(:,i)=amp0*exp(ima

24、*2*pi*1/2*si n(theta0*pi/180)*0:eleme nt_num-1);j(:,i)=amp1*exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1);n (:,i)=amp n*exp(ima*2*pi*ra ndn (1)*0:eleme nt_nu m-1);endRx=1/L*(s+j+n)*(s+j+n);%接收信號(hào)自相關(guān)矩陣Rnj=1/L*(j+n)*(j+n);%干攏+噪聲的自相關(guān)矩陣e=exp(ima*2*pi*d_lamda*sin(theta0*pi/180)*0:eleme nt_num-1);

25、%向信號(hào)Wopt_Rx=i nv(Rx)*e/(e* in v(Rx)*e);%采用接收信號(hào)的權(quán)矢量Wopt_Rnj=i nv(Rnj)*e/(e* in v(Rnj)*e);%采用干攏+噪聲信號(hào)的權(quán)矢量for j=1:le ngth(theta)a=exp(ima*2*pi*d_lamda*si n( theta(j)*pi/180)*0:eleme nt_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);figu

26、re;plot(theta,F1,theta,F2,r),grid on;hold on plot(theta0,-50:0,.);plot(theta1,-50:0,.)xlabel(theta/0);ylabel(F(1,2)/dB);title(不同方法估計(jì)協(xié)方差矩陣的Cap on波束形成);axis(-90 90 -50 0);仿真圖如下：來波方實(shí)用文檔10.%多點(diǎn)約束的Capon波束形成和方向圖clc;clear all ;close all;ima=sqrt(-1);eleme nt_num=8;%陣元數(shù)damda=1/2;%陣元間距與波長的關(guān)系theta=-90:0.3:90;%

27、搜索圍theta0=0;%來波方向theta1=20;%干擾方向theta2=50;%干擾方向L=512;%米樣單兀數(shù)Rx=zeros(element_num,element_num);%產(chǎn)生協(xié)方差矩陣for i=1:L;amp0=10*ra ndn (1);amp1=10*ra ndn (1);amp2=50*ra ndn (1);namp=0.5*ra ndn( 1);%噪聲的幅度隨機(jī)產(chǎn)生，保證噪聲與信號(hào)之間是不相關(guān)的J(:,i)=amp1*exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1)+.amp2*exp(ima*2*pi*1

28、/2*si n(theta2*pi/180)*0:eleme nt_num-1)+.n amp*exp(ima*2*pi*ra ndn (1)*0:eleme nt_nu m-1);x(:,i)=amp0*exp(ima*2*pi*1/2*si n(theta0*pi/180)*0:eleme nt_num-1)+.*20LUET -2530&0-40-2002040 GO 30th&tS 0T同方法佔(zhàn)計(jì)協(xié)方差巨陣的門液束形成o實(shí)用文檔J(:,i);%表示接收信號(hào)endRx=Rx+1/L*x*x;R=i nv(Rx);w=amp0*exp(ima*2*pi*1/2*si n(theta0*pi

29、/180)*0:eleme nt_num-1)+.amp2*exp(ima*2*pi*1/2*si n(theta2*pi/180)*0:eleme nt_num-1)+.amp1*exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1);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=10*log10(abs(f)/max(max (f);%取對(duì)數(shù)的方向圖P

30、=10*log10(abs(p)/max(max(p);%取對(duì)數(shù)的功率譜估計(jì)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

31、,-50:0,.);xlabel(theta/o);ylabel(P/dB);title(cap on beamform ing功率譜)仿真圖如下：實(shí)用文檔cspon beanrf( (jrming功率譜0-5-10-15-20-25-3035-40-50大作業(yè)%自適應(yīng)波束形成方向圖clc;clear all;close all;ima=sqrt(-1);cap口口beamfaiTTiing?!向圖-100-50050100theta+ * #-1000050100thetfiQ實(shí)用文檔eleme nt_num=8;%陣元數(shù)c=3e8;f=500e6;lamda=c/f;d=1/2*lamda

32、;%陣元間距與波長的關(guān)系theta=-90:0.5:90;%圍theta0=0;%來波方向theta1=45;%干擾方向theta2=60;%干擾方向L=512;%米樣單兀數(shù)for i=1:L;ampO=1O*ra ndn;信號(hào)的幅度隨機(jī)產(chǎn)生，保證信號(hào)之間是不相關(guān)的amp 1=100*ra ndn (1);amp2=100*ra ndn (1);amp n=10;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

33、 nt_num-1)+.amp2*exp(ima*2*pi*1/2*si n(theta2*pi/180)*0:eleme nt_num-1)+.amp n*(ra ndn (eleme nt_nu m,1)+ima*ra ndn (eleme nt_nu m,1);%產(chǎn) 生陣歹U信號(hào)_endRx=1/L*x* x;steer仁exp(ima*2*pi*1/2*si n(theta0*pi/180)*0:eleme nt_num-1)steer2=exp(ima*2*pi*1/2*si n(theta1*pi/180)*0:eleme nt_num-1)steer3=exp(ima*2*pi*

34、1/2*si n(theta2*pi/180)*0:eleme nt_num-1)C=steer1 steer2 steer3;F=1 0 0;%把0度方向作為來波方向w=i nv(Rx)*C*(i nv(C* in v(Rx)*C)*F;for j=1:le ngth(theta);a=exp(ima*2*pi*1/2*s in (theta(j)*pi/180)*0:eleme nt_num-1);f(j)=w*a;p(j)=1/(a*i nv(Rx)*a);endF=20*log10(abs(f)/(max(max(abs (f);P=20*log10(abs(p)/(max(max(abs(p);%此處是功率的對(duì)數(shù)形式figuresubplot(12