matlab_基于功率倒置算法的GPS抗干扰天线仿真以及高斯噪声的产生函数(GPS antjamming)

clear; clc; N=1000; %yita=1; %deltafai=3*pi/180; fc=1575.42e+6; twodplmud=2; sita_d=15*pi/180; %卫星仰角 fai_d=30*pi/180; %卫星方位角 sita_j1=30*pi/180;%干扰1仰角 fai_j1=30*pi/180;%干扰1方位角 sita_j2=50*pi/180;%干扰2仰角 fai_j2=200*pi/180;%干扰2方位角 sita_j3=65*pi/180;%干扰3仰角 fai_j3=350*pi/180;%干扰3方位角 % T=1/(2*16e+8); sgma_n=1;%噪声均方差 sgma_d=0.1;%信号均方差 sgma_j1=100;%干扰1均方差 sgma_j2=100;%干扰2均方差 sgma_j3=100;%干扰3均方差 delsita=2*pi/180; delfai=2*pi/180; A=[1 -0.9]; B=1; %generate background noise for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_n); [X2(i),X2(i+1)]=gngauss(0,sgma_n); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N n0(i)=Xc(i)*cos(2*pi*fc*i*T)+Xs(i)*sin(2*pi*fc*i*T); end % for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_n); [X2(i),X2(i+1)]=gngauss(0,sgma_n); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N n1(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); end % for i=1:2:N; [X1(i),X1(i+1)]=gngauss(0,sgma_n); [X2(i),X2(i+1)]=gngauss(0,sgma_n); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N n2(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); end % for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_n); [X2(i),X2(i+1)]=gngauss(0,sgma_n); end Xc=filter(B,A,X1); X=filter(B,A,X2); for i=1:N n3(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); end % for i=1:2:N [x1(i),X1(i+1)]=gngauss(0,sgma_n); [X2(i),X2(i+1)]=gngauss(0,sgma_n); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N n4(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); end %generate satellitic singal for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_d); [X2(i),X2(i+1)]=gngauss(0,sgma_d); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N s10(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); s11(i)=s10(i)*exp(j*twodplmud*pi*sin(sita_d)*cos(pi/4-fai_d)); s12(i)=s10(i)*exp(j*twodplmud*pi*sin(sita_d)*cos(3*pi/4-fai_d)); s13(i)=s10(i)*exp(j*twodplmud*pi*sin(sita_d)*cos(5*pi/4-fai_d)); s14(i)=s10(i)*exp(j*twodplmud*pi*sin(sita_d)*cos(7*pi/4-fai_d)); end % generate jam for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_j1); [X2(i),X2(i+1)]=gngauss(0,sgma_j1); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N j10(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); j11(i)=j10(i)*exp(j*twodplmud*pi*sin(sita_j1)*cos(pi/4-fai_j1)); j12(i)=j10(i)*exp(j*twodplmud*pi*sin(sita_j1)*cos(3*pi/4-fai_j1)); j13(i)=j10(i)*exp(j*twodplmud*pi*sin(sita_j1)*cos(5*pi/4-fai_j1)); j14(i)=j10(i)*exp(j*twodplmud*pi*sin(sita_j1)*cos(7*pi/4-fai_j1)); end for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_j1); [X2(i),X2(i+1)]=gngauss(0,sgma_j1); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N j20(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); j21(i)=j20(i)*exp(j*twodplmud*pi*sin(sita_j2)*cos(pi/4-fai_j2)); j22(i)=j20(i)*exp(j*twodplmud*pi*sin(sita_j2)*cos(3*pi/4-fai_j2)); j23(i)=j20(i)*exp(j*twodplmud*pi*sin(sita_j2)*cos(5*pi/4-fai_j2)); j24(i)=j20(i)*exp(j*twodplmud*pi*sin(sita_j2)*cos(7*pi/4-fai_j2)); end % for i=1:2:N [X1(i),X1(i+1)]=gngauss(0,sgma_j1); [X2(i),X2(i+1)]=gngauss(0,sgma_j1); end Xc=filter(B,A,X1); Xs=filter(B,A,X2); for i=1:N j30(i)=Xc(i)*cos(2*pi*fc*i*T)+j*Xs(i)*sin(2*pi*fc*i*T); j31(i)=j30(i)*exp(j*twodplmud*pi*sin(sita_j3)*cos(pi/4-fai_j3)); j32(i)=j30(i)*exp(j*twodplmud*pi*sin(sita_j3)*cos(3*pi/4-fai_j3)); j33(i)=j30(i)*exp(j*twodplmud*pi*sin(sita_j3)*cos(5*pi/4-fai_j3)); j34(i)=j30(i)*exp(j*twodplmud*pi*sin(sita_j3)*cos(7*pi/4-fai_j3)); end %the receivers's signal x0=s10+j10+j20+j30+n0; x1=s11+j11+j21+j31+n1; x2=s12+j12+j22+j32+n2; x3=s13+j13+j23+j33+n3; x4=s14+j14+j24+j34+n4; x=[x1;x2;x3;x4]; %computing Rxx for i=1:4 for k=1:4 Rxx(i,k)=0; for m=1:N Rxx(i,k)=Rxx(i,k)+x(i,m)'*x(k,m)/N; end end end Rxx_1=eye(4)/Rxx; for i=1:4 rx0(i)=0; for m=1:N rx0(i)=rx0(i)+x(i,m)'*x0(m)/N; end end Wopt=Rxx_1*rx0'; sita=pi/180:pi/180:pi/2; for i=1:length(sita) y1(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita(i))*cos(pi/4-fai_j1))-Wopt(2)*exp(j*twodplmud*pi*sin(sita(i))*cos(3*pi/4-fai_j1))-Wopt(3)*exp(j*twodplmud*pi*sin(sita(i))*cos(5*pi/4-fai_j1))-Wopt(4)*exp(j*twodplmud*pi*sin(sita(i))*cos(7*pi/4-fai_j1)); y2(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita(i))*cos(pi/4-fai_j2))-Wopt(2)*exp(j*twodplmud*pi*sin(sita(i))*cos(3*pi/4-fai_j2))-Wopt(3)*exp(j*twodplmud*pi*sin(sita(i))*cos(5*pi/4-fai_j2))-Wopt(4)*exp(j*twodplmud*pi*sin(sita(i))*cos(7*pi/4-fai_j2)); y3(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita(i))*cos(pi/4-fai_j3))-Wopt(2)*exp(j*twodplmud*pi*sin(sita(i))*cos(3*pi/4-fai_j3))-Wopt(3)*exp(j*twodplmud*pi*sin(sita(i))*cos(5*pi/4-fai_j3))-Wopt(4)*exp(j*twodplmud*pi*sin(sita(i))*cos(7*pi/4-fai_j3)); end fai=pi/180:pi/180:2*pi; for i=1:length(fai) y4(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita_j1)*cos(pi/4-fai(i)))-Wopt(2)*exp(j*twodplmud*pi*sin(sita_j1)*cos(3*pi/4-fai(i)))-Wopt(3)*exp(j*twodplmud*pi*sin(sita_j1)*cos(5*pi/4-fai(i)))-Wopt(4)*exp(j*twodplmud*pi*sin(sita_j1)*cos(7*pi/4-fai(i))); y5(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita_j2)*cos(pi/4-fai(i)))-Wopt(2)*exp(j*twodplmud*pi*sin(sita_j2)*cos(3*pi/4-fai(i)))-Wopt(3)*exp(j*twodplmud*pi*sin(sita_j2)*cos(5*pi/4-fai(i)))-Wopt(4)*exp(j*twodplmud*pi*sin(sita_j2)*cos(7*pi/4-fai(i))); y6(i)=1-Wopt(1)*exp(j*twodplmud*pi*sin(sita_j3)*cos(pi/4-fai(i)))-Wopt(2)*exp(j*twodplmud*pi*sin(sita_j3)*cos(3*pi/4-fai(i)))-Wopt(3)*exp(j*twodplmud*pi*sin(sita_j3)*cos(5*pi/4-fai(i)))-Wopt(4)*exp(j*twodplmud*pi*sin(sita_j3)*cos(7*pi/4-fai(i))); end figure; plot(sita*180/pi,20*log10(abs(y1)/max(abs(y1))),'r',sita*180/pi,20*log10(abs(y2)/max(abs(y2))),'g',sita*180/pi,20*log10(abs(y3)/max(abs(y3))),'b'); Ylabel('GAIN/dB'); grid; figure %plot(fai*180/pi,20*log10(abs(y4)/max(abs(y4))),'r');grid; plot(fai*180/pi,20*log10(abs(y4)/max(abs(y4))),'r',fai*180/pi,20*log10(abs(y5)/max(abs(y5))),'g',fai*180/pi,20*log10(abs(y6)/max(abs(y6))),'b'); grid;