空地导弹全弹道建模
空地导弹全弹道建模一、自由下滑段：
此时发动机未工作，P=0.
二、固定弹道倾角下滑段
三、转平飞过渡段
四、平飞段
五、比例导引段
% 程序名称：空地导弹弹道设计
% 参考文献：某型空地导弹总体设计，西工大硕士论文
% 编写：appe1943，XJTU
% zhistar@https://www.wendangku.net/doc/d011274042.html,
% 2009-10-10
%==================================================================== ============
% 弹道设计说明：
% 整个弹道分为五部分：自由下滑段、固定弹道倾角下落、转平飞过渡段、平飞段、导引
段
% 1、自由下滑段：导弹自6000m高度从载机发射，初始速度270m/s，弹道倾角0°，迎角0°。
% 自由下滑时间2s，呈现抛物线运动。
% 2、固定弹道倾角下落：自由下滑段结束后，发动机工作，在自动驾驶仪的控制下，导弹以固定
% 弹道倾角-20°快速下落。
% 3、转平飞过渡段：在到达3000m高度时，过渡段采用指数拉平关系。
% 4、平飞段：巡航高度500m
% 5、导引段：飞至目标10km时，按比例导引攻击目标
%==================================================================== =============
close all
clear
clc
Xm=0; %导弹坐标X-axis
Ym=6000; %导弹坐标Y-axis
m=1000; %导弹初始质量
v=270; %导弹初始速度
theta=-0.1; %导弹初始弹道倾角
alpha=0; %导弹初始攻角
T=0; %俯仰角
g=9.81; %当地重力加速度
rou=1.1673; %大气密度
s=0.2552; %特征面积
C_X=0.25023; %阻力系数
C_Y=0.60242; %升力系数
Xt=1.0e+5; %地面目标坐标X-axis
Yt=0; %地面目标坐标Y-axis
t=0; %飞行时间
Dt=0.01; %积分长度
i=1;
%-------------------------------------------------------------------------------------
%自由下滑段
while (t<=2)
P=0; %自由下滑段发动机不工作
q=rou*v^2/2; %动压头
X=C_X*q*s; %阻力
Y=C_Y*q*s; %升力
v=v+(P*cos(alpha)-X-m*g*sin(theta))/m*Dt;
theta=theta+(P*sin(alpha)+Y-m*g*sin(theta))/(m*v)*Dt;
Xm=Xm+v*cos(theta)*Dt;
Ym=Ym+v*sin(theta)*Dt;
t=t+Dt;
%存储中间过程，便于以后作图
Xstore(i)=Xm;
Ystore(i)=Ym;
thetaStore(i)=theta;
vstore(i)=v;
alphastore(i)=alpha;
i=i+1;
end
%--------------------------------------------------------------------------
%固定倾角下滑段
%注意：在固定弹道倾角下滑阶段，发动机一定保持在工作状态，从而通过调节攻角变化使倾角为定值；
while Ym>=3000
P=30000;
deltaM=-1; %秒流量
theta=-20*pi/180; %固定弹道倾角
q=rou*v^2/2; %动压头
X=C_X*q*s; %阻力
Y=C_Y*q*s; %升力
v=v+(P*cos(alpha)-X-m*g*sin(theta))/m*Dt;
Xm=Xm+v*cos(theta)*Dt;
Ym=Ym+v*sin(theta)*Dt;
alpha=asin((m*g*cos(theta)-Y)/P);
m=m+deltaM*Dt;
t=t+Dt;
%存储中间过程，便于以后作图
Xstore(i)=Xm;
Ystore(i)=Ym;
thetaStore(i)=theta;
vstore(i)=v;
alphastore(i)=alpha;
i=i+1;
end
%--------------------------------------------------------------------------
%转平飞过渡段
t0=t;
k=log(3000)-log(500);
while Ym>500
P=30000;
deltaM=-1; %秒流量
q=rou*v^2/2; %动压头
X=C_X*q*s; %阻力
Y=C_Y*q*s; %升力
v=v+(P*cos(alpha)-X-m*g*sin(theta))/m*Dt;
theta=theta+(P*sin(alpha)+Y-m*g*sin(theta))/(m*v)*Dt; Xm=Xm+v*cos(theta)*Dt;
Ym=Ym-3000*k*exp(-k*(t-t0))*Dt;
m=m+deltaM*Dt;
%alpha=asin((m*g*cos(theta)-Y)/P);
t=t+Dt;
%存储中间过程，便于以后作图
Xstore(i)=Xm;
Ystore(i)=Ym;
thetaStore(i)=theta;
vstore(i)=v;
alphastore(i)=alpha;
i=i+1;
end
%-------------------------------------------------------------------------- %平飞段弹道
R=sqrt((Xt-Xm)^2+(Yt-Ym)^2);%弹目距离
while R>10000
P=30000;
deltaM=-1; %秒流量
q=rou*v^2/2; %动压头
X=C_X*q*s; %阻力
Y=C_Y*q*s; %升力
v=v+(P*cos(alpha)-X)/m*Dt;
theta=0;
Xm=Xm+v*Dt;
m=m+deltaM*Dt;
alpha=asin((m*g-Y)/P);
t=t+Dt;
R=sqrt((Xt-Xm)^2+(Yt-Ym)^2);
%存储中间过程，便于以后作图
Xstore(i)=Xm;
Ystore(i)=Ym;
thetaStore(i)=theta;
vstore(i)=v;
alphastore(i)=alpha;
i=i+1;
end
%-------------------------------------------------------------------------- %导引弹道
N=3; %比例导引常数
Dt1=Dt;
while R>20
P=30000;
deltaM=-1;
%%Dr=-v*(Xt-Xm)/R;
ee=atan(-Ym/(Xt-Xm))-theta;
Dq=N*v/R*sin(ee);
%A=N*Dr*Dq; %法向加速度，此处没有用到
q=rou*v^2/2; %动压头
X=C_X*q*s; %阻力
Y=C_Y*q*s; %升力
v=v+(P*cos(alpha)-X-m*g*sin(theta))/m*Dt1;
theta=theta+Dq*Dt1;
Xm=Xm+v*cos(theta)*Dt1;
Ym=Ym+v*sin(theta)*Dt1;
if Ym<0 || Xm>Xt
stop
end
m=m+deltaM*Dt1;
R=sqrt((Xt-Xm)^2+(Yt-Ym)^2)
t=t+Dt1;
%存储中间过程，便于以后作图
Xstore(i)=Xm;
Ystore(i)=Ym;
thetaStore(i)=theta;
vstore(i)=v;
alphastore(i)=alpha;
i=i+1;
end
%-------------------------------------------------------------------------- %输出结果
figure(1)
plot((1:i-1)*Dt,thetaStore*180/pi)
grid on
title('弹道倾角变化曲线')
xlabel('时间/s')
ylabel('弹道倾角/Deg')
figure(2)
plot(Xstore,Ystore,Xt,Yt,'r.')
grid on
title('导弹运动轨迹')
xlabel('水平')
ylabel('垂直')
figure(3)
plot((1:i-1)*Dt,vstore)
grid on
xlabel('时间/s')
ylabel('速度/(m/s)')
figure(4)
plot((1:i-1)*Dt,alphastore*180/pi)
grid on
xlabel('时间/s')
ylabel('攻角/Deg')
参考文献：
某型空地导弹总体设计。西北工业大学硕士论文【2002】。王勇。谷良贤