编程.docx
《编程.docx》由会员分享,可在线阅读,更多相关《编程.docx(19页珍藏版)》请在冰豆网上搜索。
编程
(1)在不考虑系统损耗、噪声杂波以及地球曲率等环境因素的条件下利用单脉冲根据最基本的雷达方程实现对静止目标的测距,雷达方程为:
Rmax^4=(Pt*G^2*nambda^2*sigma)/((4*pi)^3*Smin),G=4*pi*Ae/nambda^2
其中,Pt为0.15MW,G=20dB,f0=1.5GHz,sigma=1m^2,Smin=10^(-10)W
function[R]=radar_function(pt,g,f0,sigma,Smin)
c=3e8;
nambda=c/f0;
g_db=10^(g/10);
pg=pt*g_db^2*nambda^2*sigma;
pb=(4*pi)^3*Smin;
R=(pg/pb)^(1/4);
return
[R]=radar_function(1.5*10^5,20,1.5*10^9,1,10*10^(-10))
R=
2.3449e+003
(2)在考虑接收机噪声和系统损耗的条件下,在
(1)中的其他条件和参数不变的情况下根据新的雷达方程测距,雷达方程为:
Rmax^4=(Pt*G^2*nambda^2*sigma)/((4*pi)^3*k*T0*B*Fn*(S/N)min*Loss),G=4*pi*Ae/nambda^2
Pt为0.15MW,G=20dB,f0=1.5GHz,sigma=1m^2,k=1.38*10^(-23)J/K,T0=290K,B=5MHz,Fn=3dB,(S/N)min=2.5dB,loss=6dB
function[R]=radar_range(pt,g,f0,sigma,t0,b,fn,loss,snr)
c=3e8;
nambda=c/f0;
nambda_sqdb=10*log10(nambda^2);
pt_db=10*log10(pt);
sigma_db=10*log10(sigma);
four_pi_cub=10*log10((4*pi)^3);
k=1.38e-23;
kT0=k*t0;
kT0_db=10*log10(kT0);
b_db=10*log10(b);
pg=pt_db+2*g+nambda_sqdb+sigma_db;
pb=four_pi_cub+kT0_db+b_db+fn+loss+snr;
R=0.25*(pg-pb);
return
[R]=radar_range(0.15*10^6,20,1.5*10^9,1,290,5*10^6,3,6,2.5)
R=
42.5732%%用分贝表示的距离
%%用输出是SNR来表示的雷达方程,考虑了系统损耗、噪声系数
%%SNR=(Pt*G^2*nambda^2*sigma)/((4*pi)^3*k*T0*B*Fn*L*R^4)
function[snr]=radar_eq(pt,freq,g,sigma,te,b,nf,loss,range)
c=3e8;
lambda=c/freq;
p_peak=10*log10(pt);
lambda_sqdb=10*log10(lambda^2);
sigmadb=10*log10(sigma);
four_pi_cub=10*log10((4*pi)^3);
k_db=10*log10(1.38e-23);
te_db=10*log10(te);
b_db=10*log10(b);
range_pwr4_db=10*log10(range.^4);
num=p_peak+g+g+lambda_sqdb+sigmadb;
den=four_pi_cub+k_db+te_db+b_db+nf+loss+range_pwr4_db;
snr=num-den;
return
%%以L表示的雷达损耗
closeall;clearall
pt=1.5e6;
freq=5.6e9;
g=45;
sigma=0.1;
te=290;
b=5e6;
nf=3;
loss=6;
range=linspace(25e3,165e3,1000);
snr1=radar_eq(pt,freq,g,sigma,te,b,nf,loss,range);
snr2=radar_eq(pt,freq,g,sigma/10,te,b,nf,loss,range);
snr3=radar_eq(pt,freq,g,sigma*10,te,b,nf,loss,range);
figure
(1)
rangekm=range./1000;
plot(rangekm,snr3,'k',rangekm,snr1,'k-.',rangekm,snr2,'k:
')
grid
legend('\sigma=0dBsm','\sigma=-10dBsm','\sigma=-20dBsm')
xlabel('Detectionrange-Km');
ylabel('SNR-dB');
snr1=radar_eq(pt,freq,g,sigma,te,b,nf,loss,range);
snr2=radar_eq(pt*.4,freq,g,sigma,te,b,nf,loss,range);
snr3=radar_eq(pt*1.44,freq,g,sigma,te,b,nf,loss,range);
figure
(2)
plot(rangekm,snr3,'k',rangekm,snr1,'k-.',rangekm,snr2,'k:
')
grid
legend('Pt=2.16MW','Pt=1.5MW','Pt=0.6MW')
xlabel('Detectionrange-Km');
ylabel('SNR-dB');
对三种不同数值的RCS,SNR相对检测距离的曲线
对三种不同数值的雷达峰值功率,SNR相对检测距离的曲线
从上可以看出,对于一组给定的雷达参数,最大检测距离是所要求的最小可检测SNR的函数;
(3)考虑脉冲积累条件运用到雷达方程
closeall;clearall
pt=4;
freq=94e9;
g=47;%?
?
?
dB
sigma=20;%?
?
?
dB,RCS
te=290;
b=20e6;
nf=7;
loss=10;
range=linspace(1.e3,12e3,1000);
snr1=radar_eq(pt,freq,g,sigma,te,b,nf,loss,range);
Rnewci=(94^0.25).*range;
snrCI=snr1+10*log10(94);
figure
(1)
rangekm=range./1000;
plot(rangekm,snr1,'k',Rnewci./1000,snr1,'k-.')
axis([112-2045])
grid
legend('singlepulse','94pulseCI')
xlabel('Detectionrange-Km');
ylabel('SNR-dB');
snr_b10=10.^(snr1./10);
SNR_1=snr_b10./(2*94)+sqrt(((snr_b10.^2)./(4*94*94))+(snr_b10./94));
LNCI=(1+SNR_1)./SNR_1;
NCIgain=10*log10(94)-10*log10(LNCI);
Rnewnci=((10.^(0.1*NCIgain)).^0.25).*range;
snrnci=snr1+NCIgain;
figure
(2)
plot(rangekm,snr1,'k',Rnewnci./1000,snr1,'k-.',Rnewci./1000,snr1,'k:
')
axis([112-2045])
grid
legend('singlepulse','94pulseNCI','94pulseCI')
xlabel('Detectionrange-Km');
ylabel('SNR-dB');
前面考虑的都是在给定sigma即RCS恒定的前提下进行的,即是给定目标然后改变条件来测距,下面将进行的工作将是需要检测出目标并测距
考虑Swerling模型,将给定的RCS即sigma用Swerling模型来描述,即是认为RCS不再是恒定不变的,而是按照某种规则改变的,从而将引起起伏损耗,而对每一种Swerling模型来看,它计算出了作为检测概率、虚警概率和积累脉冲数的函数所要求的信噪比,也就是说这种规则即是按照Swerling模型给出的信噪比和检测概率、虚警概率等的函数,如果给定虚警概率和积累脉冲数就可以得到检测概率和信噪比的关系,
(4)以SwerlingI为例,认为目标满足SwerlingI型起伏,那么如何描述距离关系呢?
给定虚警概率和积累脉冲数,检测概率和信噪比相关,信噪比和距离相关,所以可以表示检测概率和距离的关系
现在考虑在某一距离范围内有一个RCS满足SwerlingI模型变化的静止目标,
设计一部地基雷达实现搜索和检测任务,威胁由具有平均RCS为4m^2的飞机和具有平均RCS为0.5m^2的导弹组成,飞机和导弹目标模型分别服从SwerlingI型和SwerlingIII型起伏,同时假设在最大距离处要求PD>=0.995,Pfa=10^(-7)或更好;导弹高度2km,飞机高度约7km,假定是具有360。
方位覆盖范围的扫描雷达,扫描速率小于等于2s/圈,假定L到X波段,要求150m的距离分辨率,假定噪声系数F=6dB,总的接收机损失L=8dB,目前使用方位波束宽度小于3。
的扇形波束,假定合理的检测门限是SNR为13dB,最后假设地球是平的。
%%计算改善因子,依据的是Peebles报道的一个根据经验推导的改善因子的表达式,精度在0.8dB之内
functionimpr_of_np=improv_fac(np,pfa,pd)
fact1=1+log10(1/pfa)/46.6;
fact2=6.79*(1+0.235*pd);
fact3=1-0.14*log10(np)+0.0183*(log10(np)^2);
impr_of_np=fact1*fact2*fact3*log10(np);
return
closeall;clearall;
np=4;
pfa=1e-7;
pdm=0.99945;
pda=0.99812;
Im=improv_fac(np,pfa,pdm);
Ia=improv_fac(np,pfa,pda);
Lm=10*log10(np)-Im;
La=10*log10(np)-Ia;
pt=114.7e3;
freq=3e9;
g=34.5139;
sigmam=0.5;
sigmaa=4;
te=290;
b=1e6;
nf=6;
loss=8;
losstm=loss+Lm;
lossta=loss+La;
range=linspace(20e3,120e3,1000);
snrmnci=radar_eq(np*pt,freq,g,sigmam,te,b,nf,losstm,range);
snrm=radar_eq(pt,freq,g,sigmam,te,b,nf,losstm,range);
snranci=radar_eq(np*pt,freq,g,sigmaa,te,b,nf,losstm,range);
snra=radar_eq(pt,freq,g,sigmaa,te,b,nf,losstm,range);
rangekm=range./1000;
figure
(1)
subplot(2,1,1)
plot(rangekm,snrmnci,'k',rangekm,snrm,'k-.')
grid
legend('With4-pulseNCI','Singlepulse')
ylabel('SNR-dB');
title('Missilecase')
subplot(2,1,2)
plot(rangekm,snranci,'k',rangekm,snra,'k-.')
grid
legend('With4-pulseNCI','Singlepulse')
ylabel('SNR-dB');
title('Aircraftcase')
xlabel('Detectionrange-Km')
(5)动目标的检测
要检测到目标并确定其距离,需要得出目标回波信号相对距离的波形(有目标就是有相应的尖峰,尖峰对应的距离就是目标距离)
上面讲到在最大可检测距离的时候得到最小可检测信噪比,从图中可以看到,当距离在减小的时候SNR在增加,即是说可以检测到目标(满足相应的检测概率和虚警率)
(6)给定一个发射信号,如脉冲信号Arect(t/tau),脉宽为tau,峰值功率为A^2,要检测到目标的存在并判断其距离就需要产生一个回波信号,对回波信号进行一系列的处理得到回波信号和距离的图形(如存在一个尖峰,对应尖峰的横坐标就是距离),可以根据图形判断是否存在目标以及雷达的探测威力,那么怎样来构造回波信号呢?
如回波信号的形式是K*A*rect((t-taum)/tau),也就是说要得到回波信号就需要知道信号的衰减和时延,时延可以用距离R来表示,R又可以通过满足一定条件(如Swerling模型)的雷达方程给出,即时延taum可以用雷达方程给出(雷达方程中还给出了噪声系数和系统衰减)那么回波信号中的衰减因子K就和传播过程中的衰减有关,前提是目标反射造成的损耗包括在雷达方程中,那么怎么来描述传播过程中的衰减呢?
而且这个衰减还是一个和距离有关的函数(变化量)。
(7)假设发射信号是线性调频信号,一个满足SwerlingI模型的目标从足够远的R0以恒定速度V0向雷达驶来,
目标的距离可以通过最小可检测信噪比得到,即是只要找到符合检测概率的最小可检测信噪比(绘出SNR和Pd的关系曲线)并将它带到雷达方程中就可以算出最大距离;
目标的检测即目标存在性的判定是通过门限实现的
对单个脉冲信号
因为虚警概率Pfa=exp(-VT2/(2*Ψ^2))
所以门限VT=(2*Ψ^2*ln(1/Pfa))^(1/2)
当使用非相干积累时
Pfa=1-gamma(VT/sqrt(np),np-1)可以得到VT
发射信号-Chirp信号
S(t)=rect(t/Tp)*exp(j*pi*beita*t^2)
Pt(峰值功率)
10KW
40dB
Tp(脉宽)
10us
-50dB
PRF(脉冲重复频率)
10000Hz
40dB
f(S波段-载频)
3GHz
94.77dB
B
10MHz
70dB
Ae
2m^2
3dB
Sigma(目标横截面积)
2m^2
3dB
G=4*pi*Ae/λ^2(天线增益)
800*pi
34dB
Fn(噪声系数)
6dB
6dB
L(系统损耗)
10dB
10dB
fs(采样频率)
50MHz
ThetaB(波束宽度)
1.5度
Wr(天线转速)
25转/min
N(积累脉冲数)ThetaB*fp/6/wr
10
10dB
KTe
4*10^(-21)
-203.98dB
V(目标速度)
500m/s
R0(目标初始距离)
脉压前
SNR=Pt*G^2*lamda^2*sigma/((4*pi)^3*k*Te*B*F*L*R^4)
=Pt+2*G+2*lamda+sigma-val1-val2-B-F-L-R^4(以dB为单位)
=
脉压增益
G1=1000
衰减因子k=sqrt(Pt*G*sigma*Ae/((4*pi)^2*R^4*K*Te*B*Fn))
=6.32*10^9
假定一速度v=200m/s的目标从距离R0向雷达靠近,目标服从Swerling0型,RCS=sigma=3m^2
回波信号Si(t)=k*rect((t-taum)/Tp)*exp(j*pi*beita*(t-taum)^2)
延迟时间taum=2*(R0-m*Tr*v)/c,R0为足够远的初始距离,m为脉冲积累数=10,由波束宽度和天线转速确定,Tr为PRT=1/1000,v=200m/s,c为电磁波在空中传播速度c=3e8m/s,k为衰减因子,由雷达方程算出
K=sqrt((Pt*G*Ae*sigma)/((4*pi)^2*Rm^4))
噪声为零均值、方差为sigma^2=k0*Te*B*Fn的高斯白噪声
y=wgn(1,N,10log(P))
门限的选择:
当只使用单个脉冲时,检测门限VT与虚警概率Pfa=exp(-VT^2/(2*fai^2))有关,fai^2为噪声功率;当使用非相干积累时,对于任意脉冲数,门限与虚警概率Pfa的一般关系式可以表示为:
Pfa=1-gamma(VT/sqrt(np),np-1),gamma函数形式为
对于np>1的情况,Marcum定义的虚警概率为Pfa=ln
(2)*(np/nfa),nfa为Marcum虚警数;
函数threshold.m使用在Newton-Raphson方法中使用的递归公式来计算门限
function[pfa,vt]=threshold(nfa,np)
delmax=.00001;
eps=0.000000001;
delta=10000.;
pfa=np*log
(2)/nfa;
sqrtpfa=sqrt(-log10(pfa));
sqrtnp=sqrt(np);
vt0=np-sqrtnp+2.3*sqrtpfa*(sqrtpfa+sqrtnp-1.0);
vt=vt0;
while(abs(delta)>=vt0)
igf=incomplete_gamma(vt0,np);
num=0.5^(np/nfa)-igf;
temp=(np-1)*log(vt0+eps)-vt0-factor(np-1);
deno=exp(temp);
vt=vt0+(num/(deno+eps));
delta=abs(vt-vt0)*10000.0;
vt0=vt;
End
[pfa,vt]=threshold(log
(2)*10^6,10)
pfa=
1.000000000000000e-005
vt=
29.522281*********
function[value]=incomplete_gamma(vt,np)
formatlong
eps=1.000000001;
if(np==1)
value1=vt*exp(-vt);
value=1-exp(-vt);
return
end
sumold=1;
sumnew=1;
calc1=1;
calc2=np;
xx=np*log(vt+0.0000000001)-vt-factor(calc2);%usefunctionfactor
temp1=exp(xx);
temp2=np/(vt+0.0000000001);
diff=.0;
ratio=1000;
if(vt>=np)
while(ratio>=eps)
diff=diff+1;
calc1=calc1*(calc2-diff)/vt;
sumnew=sumold+calc1;
ratio=sumnew/sumold;
sumold=sumnew;
end
value=1-temp1*sumnew*temp2;
return
else
diff=0;
sumold=1;
ratio=1000;
calc1=1;
while(ratio>=eps)
diff=diff+1;
calc1=calc1*vt/(calc2+diff);
sumnew=sumold+calc1;
ratio=sumnew/sumold;
sumold=sumnew;
end
value=temp1*sumnew;
end
得到当积累脉冲数np=10,虚警概率Pfa=10^(-5)时的门限值VT=29.522281*********
检测概率相对SNR的曲线