CN105718660A - 临近空间大范围机动弹道三维包络计算方法 - Google Patents

Abstract

本发明公开了一种临近空间大范围机动弹道三维包络计算方法，包括以下步骤：建立弹道包络侧向宽度与纵程的数学模型；建立换极坐标系；建立换极坐标系中飞行器动力学模型；计算换极坐标系弹道参数；计算换极坐标系弹道侧向包络；计算一般坐标系弹道三维包络。本发明所述临近空间大范围机动弹道三维包络计算方法立足于飞行器机动能力分析，兼顾射向、纵程及建模误差等对弹道包络的影响，可广泛应用于具有相同动力学描述的飞行器，计算得到的包络可覆盖具有不同机动特性的弹道，能适应考虑禁飞区弹道及侧向最大范围机动弹道等极端情况，为飞行器机动能力分析和射前空间信息获取提供了重要的方法基础。

Description

临近空间大范围机动弹道三维包络计算方法

技术领域

本发明涉及飞行器弹道规划领域，特别涉及一种临近空间大范围机动弹道三维包络计算方法。

背景技术

本发明涉及的临近空间的大范围机动弹道特指高超声速滑翔飞行器包含初始下降段弹道和滑翔段弹道在内的弹道，其起点为再入点，终点为滑翔段弹道终点。该段弹道是高超声速滑翔飞行器全程弹道的主要组成部分，具有长时间处于临近空间并进行大范围侧向机动的特性，弹道特征显著区别于弹道式弹道。在再入点和目标点确定的条件下，临近空间大范围机动弹道的侧向包络可达几千公里，弹道形态因飞行任务不同而迥异。

进行弹道包络计算的意义包括但不限于如下两个方面。首先，计算弹道包络是分析飞行器机动能力的重要手段，为弹道规划和飞行任务制定提供先验知识。其次，计算弹道覆盖区域为射前高精度引力模型构建提供必要的空间位置信息，是修正引力高阶项影响的必要条件，对于快速发射和提升飞行器导航、制导精度具有重要意义。目前，尚无公开资料提出有关临近空间大范围机动弹道三维包络的计算方法或准则。

发明内容

本发明针对临近空间大范围机动弹道覆盖区域确定问题，首次提出一种临近空间大范围机动弹道三维包络计算方法。方法基于飞行器能力分析结果，能够根据再入点和目标点快速确定弹道三维包络。计算得到的包络可覆盖具有不同机动特性的弹道，为飞行器机动能力分析和射前空间信息获取提供了重要的方法基础。

本发明通过以下技术方案来实现上述目的：

一种临近空间大范围机动弹道三维包络计算方法，所述临近空间大范围机动弹道包括初始下降段弹道和滑翔段弹道，弹道起点为再入点，弹道终点为滑翔段弹道终点，所述计算方法包括以下步骤：

(1)建立弹道包络侧向宽度与纵程的数学模型；

(2)建立换极坐标系；

(3)建立换极坐标系中飞行器动力学模型；

(4)计算换极坐标系弹道参数；

(5)计算换极坐标系弹道侧向包络；

(6)计算一般坐标系弹道三维包络。

具体地，

(1)建立弹道包络侧向宽度与纵程的数学模型

令弹道再入点为I，其经度为λ_I、地心纬度为φ_I；弹道落点为T，其经度为λ_T、地心纬度为φ_T；以由再入点和落点确定的再入大圆弧平面为对称面，称沿射向方向对称面以左弹道为左侧机动弹道，沿射向方向对称面以右弹道为右侧机动弹道；记最大左侧机动弹道距对称面的最大距离为左边界B_l，记最大右侧机动弹道距对称面的最大距离为右边界B_r；

根据给定的飞行器再入点飞行状态参数、弹道终端飞行状态参数、飞行过程约束条件及终端约束条件，针对再入点为(λ_I,φ_I)＝(0°,0°)，落点为(λ_ki,0°)(i＝1,2,3,4,…)的低空大范围机动弹道进行弹道计算，记各落点对应的弹道纵程L_i，由式(1)计算L_i，

L_i＝R_earccos(sinφ_Isinφ_I+cosφ_Icosφ_Icos(λ_T-λ_T))(1)

式中：R_e为地球半径；

通过改变再入点初始速度方位角获得最大左侧机动弹道和最大右侧向机动弹道，记各落点对应的弹道侧向包络宽度W_io，由式(2)计算W_io，

W_io＝B_l-B_r(2)

考虑建模误差，取n倍W_io为实际侧向包络宽度W_i，对L_i(i＝1,2,3,4,…)和W_i(i＝1,2,3,4,…)进行拟合，建立如式(3)所示的弹道包络侧向宽度与纵程的数学模型，确定模型系数a_i，

$W=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{L}^i---\left(3\right)$

(2)建立换极坐标系

A.定义一个再入大圆弧平面作为换极赤道平面：

a)对目标点确定的情况，将再入点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；

b)对于目标点未确定的情况，根据再入点位置及速度方位角确定的再入大圆弧平面作为换极赤道面；

B.基于换极赤道平面定义换极坐标系O_E为地心，轴沿再入点地心矢径方向，轴在换极赤道面内垂直于轴指向目标点方向，轴与轴、轴构成右手系。

(3)建立换极坐标系中飞行器动力学模型

根据步骤(2)建立的换极坐标系建立换极坐标系中飞行器的动力模型，且用表示换极坐标系中各物理量，用x表示一般坐标系中各物理量：

在换极坐标系中建立以时间为自变量的滑翔飞行器动力学方程，其飞行状态量为换极后的经度地心纬度航迹偏航角速度速度倾角和地心距

$\left\{\begin{array}{c}\frac{d\hat{r}}{dt}=\hat{V}\sin \widehat{\mathrm{\θ}}\\ \frac{d\widehat{\mathrm{\λ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}}{\hat{r}\cos \widehat{\mathrm{\φ}}}\\ \frac{d\widehat{\mathrm{\φ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}}{\hat{r}}\\ \frac{d\widehat{\mathrm{\σ}}}{dt}=\frac{L\sin \mathrm{\υ}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}\tan \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}-{\hat{g}}_{{\mathrm{\ω}}_e}\frac{\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+C_{\mathrm{\σ}}+{\tilde{C}}_{\mathrm{\σ}}\\ \frac{d\widehat{\mathrm{\θ}}}{dt}=\frac{L}{\hat{V}}\cos \mathrm{\υ}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}+\frac{{\hat{g}}_r\cos \widehat{\mathrm{\θ}}}{\hat{V}}+\frac{{\hat{g}}_{{\mathrm{\ω}}_e}}{\hat{V}}(\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}}-\cos \widehat{\mathrm{\σ}}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}})+C_{\mathrm{\θ}}+{\tilde{C}}_{\mathrm{\θ}}\\ \frac{d\hat{V}}{dt}=-D+{\hat{g}}_r\sin \widehat{\mathrm{\θ}}+{\hat{g}}_{{\mathrm{\ω}}_e}(\cos \widehat{\mathrm{\σ}}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}}+\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}})+{\tilde{C}}_V\end{array}\right.---\left(4\right)$

其中，C_σ、C_θ为哥氏加速度项，和为牵连加速度项，

$\left\{\begin{array}{c}{C}_{\mathrm{\σ}}=(2{\mathrm{\ω}}_{ex}-2\tan \widehat{\mathrm{\θ}}({\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\σ}}\left)\right)\\ {\tilde{C}}_{\mathrm{\σ}}=-\frac{\hat{r}}{\hat{V}\cos \widehat{\mathrm{\θ}}}({\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}})\\ C_{\mathrm{\θ}}=2({\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}})\\ {\tilde{C}}_{\mathrm{\θ}}=\frac{\hat{r}}{\hat{V}}\[{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\cos \widehat{\mathrm{\θ}}\]\\ {\tilde{C}}_V=\hat{r}\[-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\sin \widehat{\mathrm{\θ}}\]\end{array}\right.---\left(5\right)$

其中，

$\left\{\begin{array}{c}{\mathrm{\ω}}_{ex}={\mathrm{\ω}}_e(\cos \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\sin \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\sin \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\\ {\mathrm{\ω}}_{ey}={\mathrm{\ω}}_e(-\sin \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\cos \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p)\\ {\mathrm{\ω}}_{ez}={\mathrm{\ω}}_e(-\cos \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p-\sin \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\cos \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\end{array}\right.---\left(6\right)$

其中，ω_e为地球旋转加速度矢量，λ_p和φ_p为换极后极点P的经度和地心纬度，A_p为P的方位角；

根据换极坐标系定义，一般坐标系与换极坐标系中地心距、当地速度倾角及速度的定义一致，

$r=\hat{r},\mathrm{\θ}=\widehat{\mathrm{\θ}},V=\hat{V}---\left(7\right)$

定义

$\left[\begin{array}{ccc}{\mathrm{cos\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{cos\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\φ}}_f\\ -{\mathrm{sin\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{cos\ψ}}_f{\mathrm{cos\φ}}_f\\ {\mathrm{cos\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& -{\mathrm{cos\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\φ}}_f\end{array}\right]---\left(8\right)$

$\stackrel{\mathrm{\Δ}}{=}\left[\begin{array}{ccc}{G}_{11}& G_{12}& G_{13}\\ G_{21}& G_{22}& G_{23}\\ G_{31}& G_{32}& G_{33}\end{array}\right]$

其中，ψ_f为点F的方位角。

$\{\begin{array}{c}{G}_{11}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{12}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{13}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_1\\ G_{21}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{22}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{23}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_2\\ G_{31}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{32}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{33}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_3\end{array}---\left(9\right)$

$\{\begin{array}{c}{G}_{11}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{21}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{31}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_1\\ G_{12}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{22}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{32}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_2\\ G_{13}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{32}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{33}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_3\end{array}---\left(10\right)$

由一般坐标系中λ和φ确定换极坐标系中和的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\\ \sin \widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\end{array}\right.& \left\{\begin{array}{c}\sin \widehat{\mathrm{\φ}}=k_3\\ \cos \widehat{\mathrm{\φ}}=\sqrt{k_1^2+k_2^2}\end{array}\right.\end{array}---\left(11\right)$

由换极坐标系中和确定一般坐标系中λ和φ的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \mathrm{\λ}={\tilde{k}}_1/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\\ \sin \mathrm{\λ}={\tilde{k}}_2/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.& \left\{\begin{array}{c}\sin \mathrm{\φ}={\tilde{k}}_3\\ \cos \mathrm{\φ}=\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.\end{array}---\left(12\right)$

由一般坐标系中σ确定换极坐标系中的表达式为

$\widehat{\mathrm{\σ}}=\mathrm{\σ}+\mathrm{\η}---\left(13\right)$

其中，

$\{\begin{array}{c}\sin \mathrm{\η}=\frac{\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{cos\φ}}_P}{\cos \widehat{\mathrm{\φ}}}\\ \cos \mathrm{\η}=-\cos (A_P-\widehat{\mathrm{\λ}})\cos (\mathrm{\λ}-{\mathrm{\λ}}_P)+\sin (A_P-\widehat{\mathrm{\λ}})\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{sin\φ}}_P\end{array}---\left(14\right)$

(4)计算换极坐标系弹道参数

根据弹道再入点经度λ_I、再入点地心纬度φ_I、落点经度λ_T、落点地心纬度φ_T计算得到换极弹道参数,由计算得记换极弹道落点经度为

(5)计算换极坐标系弹道侧向包络

由式(15)计算换极弹道纵程

$\hat{L}=R_{e}arccos(sin\widehat{\mathrm{\φ}}sin{\widehat{\mathrm{\φ}}}_T+cos{\widehat{\mathrm{\φ}}}_{I}cos{\widehat{\mathrm{\φ}}}_{T}cos({\widehat{\mathrm{\λ}}}_T-{\widehat{\mathrm{\λ}}}_I\left)\right)---\left(15\right)$

将代入式(16)，可计算得到换极弹道侧向包络宽度

$\hat{W}=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{\hat{L}}^i---\left(16\right)$

将换极弹道侧向包络描述为长度为宽度为的矩形，其边界由式(17)确定，

$\{\begin{array}{c}{\hat{E}}_{\mathrm{\λ}1}={\mathrm{\λ}}_I\\ {\hat{E}}_{\mathrm{\λ}2}=\frac{180\hat{L}}{{\mathrm{\πR}}_e}\end{array},\left\{\begin{array}{c}{\hat{E}}_{\mathrm{\φ}1}=-\frac{90\hat{W}}{{\mathrm{\πR}}_e}\\ {\hat{E}}_{\mathrm{\φ}2}=\frac{90\hat{W}}{{\mathrm{\πR}}_e}\end{array}\right.---\left(17\right)$

其中，为换极系侧向包络东向下边界，为换极系侧向包络东向上边界；为换极系侧向包络北向下边界，为换极系侧向包络北向上边界。

(6)计算一般坐标系弹道三维包络

按照东向dλ、北向dφ的间隔将由式(17)确定的换极坐标系弹道侧向包络剖分为互不重叠的若干网格，记网格节点坐标为

由步骤(3)所述坐标转换关系，计算一般坐标系网格节点坐标N(λ_n,φ_n)，由N(λ_n,φ_n)确定的区域即为一般坐标系弹道侧向包络区域。

优选地，考虑机动和偏差进行蒙特卡罗打靶，以弹道高度范围的m倍作为弹道包络的垂向范围ΔH。

本发明临近空间大范围机动弹道三维包络计算方法计算得到的弹道包络可覆盖具有各种机动特性的弹道，对考虑禁飞区的中程机动弹道和侧向最大范围机动弹道等极端情况均有较好的适应性，且具有以下优点：

1)首次提出了一种针对临近空间大范围机动弹道的三维包络计算方法。方法能够根据再入点和目标点快速确定弹道三维包络，计算得到的包络可覆盖具有不同机动特性的弹道，为飞行器机动能力分析和射前空间信息获取提供了重要的方法基础。

2)基于极点变换的思想，首次将一种换极坐标系引入弹道三维包络计算方法的建立中。通过将弹道再入大圆弧平面转移到换极赤道面附近，使得换极坐标系中的弹道包络关于换极赤道面对称，从而简化了确定包络边界的计算流程。

3)引入分片定位的思想，对换极坐标系弹道包络进行空域剖分，确定包络特征点位置。根据换极系特征点位置，通过坐标转换求解一般系特征点位置，进而确定一般系弹道包络。根据计算时间和计算量要求，可灵活调节用于特征点确定的空域剖分参数。

4)方法计算流程简洁、计算速度快，确定模型参数后即可完成弹道三维包络的快速计算。模型参数确定后，无需通过多次弹道规划和弹道计算来确定飞行器可达域，因此能够满足射前任务快速规划和空间信息快速获取的要求。

5)方法可广泛适用于具有相同动力学描述的飞行器，不受飞行器特征参数改变及飞行任务改变的影响。方法对考虑禁飞区的中程机动弹道以及最大侧向机动弹道等极端情况均有较好的适应性，且可根据模型偏差调整包络宽度系数，具有较强的鲁棒性。

附图说明

图1为实施例中不同射向弹道侧向包络示意图；

图2为实施例中东向纵程5560km弹道侧向包络示意图；

图3为实施例中东向纵程6671km弹道侧向包络示意图；

图4为实施例中东向纵程7784km弹道侧向包络示意图；

图5为实施例中东向纵程8896km弹道侧向包络示意图；

图6为实施例中东向不同纵程弹道侧向包络示意图；

图7为实施例中换极坐标系及飞行状态量定义示意图；

图8为实施例中弹道再入点I和换极后极点P的位置关系示意图；

图9为实施例中换极坐标系航迹偏航角的位置关系示意图；

图10为实施例中换极系弹道侧向包络节点示意图，其中n＝1.5；

图11为实施例中一般系弹道侧向包络节点示意图，其中n＝1.5；

图12为实施例中实施例中考虑禁飞区弹道与侧向包络位置关系示意图，其中n＝1.5；

图13为实施例中最大左侧机动弹道与侧向包络位置关系示意图，其中n＝1.5；

图14为实施例中最大右侧机动弹道与侧向包络位置关系示意图，其中n＝1.5；

图15为实施例中考虑禁飞区弹道与侧向包络位置关系示意图，其中n＝1.1；

图16为实施例中最大左侧机动弹道与侧向包络位置关系示意图其中n＝1.1；

图17为实施例中最大右侧机动弹道与侧向包络位置关系示意图其中n＝1.1。

具体实施方式

下面结合附图和具体实施例对本发明作进一步说明：

基于CAV-H飞行器模型进行仿真，基本仿真条件设置为：①质量m＝908kg，参考面积Sm＝0.48375m2；②再入点状态参数：速度V_I＝6500m/s，速度倾角θ_I＝0°，高度H_I＝80km；③滑翔段终端状态参数：速度V_k＝2500m/s，速度倾角θ_k＝0°，高度H_k＝30km，终端航迹偏航角偏差不大于±5°；④飞行约束条件：最大热流密度最大动压q_max＝100kPa，最大过载n_max＝3g；⑤终端结束时距目标点的待飞航程：S_togo＝100km。

仿真计算机配置为：Intel(R)Core(TM)i5-3470CPU3.20GHz，内存为3.46GB。软件环境为WindowXP操作系统，计算程序基于VC++6.0开发。其具体步骤如下：

第一步，建立弹道包络侧向宽度与纵程的数学模型

令弹道再入点为I，其经度为λ_I、地心纬度为φ_I；弹道落点为T，其经度为λ_T、地心纬度为φ_T；以由再入点和落点确定的再入大圆弧平面为对称面，称沿射向方向对称面以左弹道为左侧机动弹道，沿射向方向对称面以右弹道为右侧机动弹道；记最大左侧机动弹道距对称面的最大距离为左边界B_l，记最大右侧机动弹道距对称面的最大距离为右边界B_r；

针对不同射向的临近空间大范围机动弹道进行仿真，其再入点经度为λ_I＝0°、再入点地心纬度为φ_I＝0°，落点经度λ_T和落点地心纬度φ_T见表1。

表1不同射向弹道落点坐标

射向	正北	正东	正南	正西
					经度(°)	0	60	0	-60
地心纬度(°)	60	0	-60	0

计算得到不同射向弹道侧向包络示意图见图1，不同射向弹道侧向包络边界见表2。由图1和表2可见，由于地球自转和非球形因素的影响，东向弹道包络最大、南北向弹道包络次之，西向弹道包络最小；南北向弹道包络向东有小范围偏移。

表2不同射向弹道侧向包络边界

射向	正东	正南	正西	正北
					左边界Bl(km)	1597	1072	1072	1072
右边界Br(km)	1597	1513	1072	1528

针对射向为正东、纵程不同的临近空间大范围机动弹道进行仿真。弹道再入点经度为λ_I＝0°，再入点地心纬度为φ_I＝0°；落点地心纬度为0°，落点经度分别为50°、60°、70°和80°。记各落点对应的弹道纵程为L_i，由式(1)计算L_i，

L_i＝R_earccos(sinφ_Isinφ_I+cosφ_Icosφ_Icos(λ_T-λ_T))(1)

其中，R_e为地球半径。

计算得到不同落点对应的纵程L_i见表3。

表3不同落点对应的纵程

落点经度(°)	50	60	70	80
					纵程Li(km)	5560	6671	7784	8896

通过改变再入点初始速度方位角获得由飞行器能力决定的最大左侧机动弹道和最大右侧向机动弹道。不同纵程弹道侧向包络示意图见图2、图3、图4、图5和图6。

记各落点对应的弹道侧向包络宽度W_io，由式(2)计算W_io，

W_io＝B_l-B_r(2)

考虑建模误差，取n倍W_io为实际侧向包络宽度W_i，计算得到不同纵程弹道侧向包络宽度见表4。

表4不同纵程弹道侧向包络宽度

对L_i(i＝1,2,3,4,…)和W_i(i＝1,2,3,4,…)进行拟合，建立如式(3)所示的弹道包络侧向宽度与纵程的数学模型，

$W=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{L}^i---\left(18\right)$

具体为W＝-1.064813L+0.000879L²-1.225047e^-7L³+4.607571e^-12L⁴

第二步，建立换极坐标系

A.定义一个再入大圆弧平面作为换极赤道平面：

对目标点确定的情况，将再入点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；

对于目标点未确定的情况，根据再入点位置及速度方位角确定的再入大圆弧平面作为换极赤道面；

B.基于换极赤道平面定义换极坐标系O_E为地心，轴沿再入点地心矢径方向，轴在换极赤道面内垂直于轴指向目标点方向，轴与轴、轴构成右手系。见图7。

第三步，建立换极坐标系中飞行器动力学模型

在换极坐标系中建立以时间为自变量的滑翔飞行器动力学方程，其飞行状态量为换极后的经度地心纬度航迹偏航角速度速度倾角和地心距

$\left\{\begin{array}{c}\frac{d\hat{r}}{dt}=\hat{V}\sin \widehat{\mathrm{\θ}}\\ \frac{d\widehat{\mathrm{\λ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}}{\hat{r}\cos \widehat{\mathrm{\φ}}}\\ \frac{d\widehat{\mathrm{\φ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}}{\hat{r}}\\ \frac{d\widehat{\mathrm{\σ}}}{dt}=\frac{L\sin \mathrm{\υ}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}\tan \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}-{\hat{g}}_{{\mathrm{\ω}}_e}\frac{\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+C_{\mathrm{\σ}}+{\tilde{C}}_{\mathrm{\σ}}\\ \frac{d\widehat{\mathrm{\θ}}}{dt}=\frac{L}{\hat{V}}\cos \mathrm{\υ}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}+\frac{{\hat{g}}_r\cos \widehat{\mathrm{\θ}}}{\hat{V}}+\frac{{\hat{g}}_{{\mathrm{\ω}}_e}}{\hat{V}}(\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}}-\cos \widehat{\mathrm{\σ}}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}})+C_{\mathrm{\θ}}+{\tilde{C}}_{\mathrm{\θ}}\\ \frac{d\hat{V}}{dt}=-D+{\hat{g}}_r\sin \widehat{\mathrm{\θ}}+{\hat{g}}_{{\mathrm{\ω}}_e}(\cos \widehat{\mathrm{\σ}}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}}+\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}})+{\tilde{C}}_V\end{array}\right.---\left(4\right)$

其中，C_σ、C_θ为哥氏加速度项，和为牵连加速度项，

$\left\{\begin{array}{c}{C}_{\mathrm{\σ}}=(2{\mathrm{\ω}}_{ex}-2\tan \widehat{\mathrm{\θ}}({\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\σ}}\left)\right)\\ {\tilde{C}}_{\mathrm{\σ}}=-\frac{\hat{r}}{\hat{V}\cos \widehat{\mathrm{\θ}}}({\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}})\\ C_{\mathrm{\θ}}=2({\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}})\\ {\tilde{C}}_{\mathrm{\θ}}=\frac{\hat{r}}{\hat{V}}\[{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\cos \widehat{\mathrm{\θ}}\]\\ {\tilde{C}}_V=\hat{r}\[-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\sin \widehat{\mathrm{\θ}}\]\end{array}\right.---\left(5\right)$

其中，

$\left\{\begin{array}{c}{\mathrm{\ω}}_{ex}={\mathrm{\ω}}_e(\cos \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\sin \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\sin \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\\ {\mathrm{\ω}}_{ey}={\mathrm{\ω}}_e(-\sin \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\cos \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p)\\ {\mathrm{\ω}}_{ez}={\mathrm{\ω}}_e(-\cos \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p-\sin \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\cos \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\end{array}\right.---\left(6\right)$

其中，ω_e为地球旋转加速度矢量，λ_p和φ_p为换极后极点P的经度和地心纬度，A_P为P的方位角,见图8。

根据换极坐标系定义，一般坐标系与换极坐标系中地心距、当地速度倾角及速度的定义一致，

$r=\hat{r},\mathrm{\θ}=\widehat{\mathrm{\θ}},V=\hat{V}---\left(7\right)$

定义

$\left[\begin{array}{ccc}{\mathrm{cos\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{cos\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\φ}}_f\\ -{\mathrm{sin\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{cos\ψ}}_f{\mathrm{cos\φ}}_f\\ {\mathrm{cos\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& -{\mathrm{cos\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\φ}}_f\end{array}\right]---\left(8\right)$

$\stackrel{\mathrm{\Δ}}{=}\left[\begin{array}{ccc}{G}_{11}& G_{12}& G_{13}\\ G_{21}& G_{22}& G_{23}\\ G_{31}& G_{32}& G_{33}\end{array}\right]$

其中，ψ_f为点F的方位角。

$\{\begin{array}{c}{G}_{11}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{12}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{13}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_1\\ G_{21}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{22}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{23}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_2\\ G_{31}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{32}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{33}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_3\end{array}---\left(9\right)$

$\{\begin{array}{c}{G}_{11}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{21}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{31}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_1\\ G_{12}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{22}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{32}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_2\\ G_{13}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{32}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{33}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_3\end{array}---\left(10\right)$

由一般坐标系中λ和φ确定换极坐标系中和的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\\ \sin \widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\end{array}\right.& \left\{\begin{array}{c}\sin \widehat{\mathrm{\φ}}=k_3\\ \cos \widehat{\mathrm{\φ}}=\sqrt{k_1^2+k_2^2}\end{array}\right.\end{array}---\left(11\right)$

由换极坐标系中和确定一般坐标系中λ和φ的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \mathrm{\λ}={\tilde{k}}_1/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\\ \sin \mathrm{\λ}={\tilde{k}}_2/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.& \left\{\begin{array}{c}\sin \mathrm{\φ}={\tilde{k}}_3\\ \cos \mathrm{\φ}=\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.\end{array}---\left(12\right)$

由一般坐标系中σ确定换极坐标系中的表达式为

$\widehat{\mathrm{\σ}}=\mathrm{\σ}+\mathrm{\η}---\left(13\right)$

其中，

$\{\begin{array}{c}\sin \mathrm{\η}=\frac{\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{cos\φ}}_P}{\cos \widehat{\mathrm{\φ}}}\\ \cos \mathrm{\η}=-\cos (A_P-\widehat{\mathrm{\λ}})\cos (\mathrm{\λ}-{\mathrm{\λ}}_P)+\sin (A_P-\widehat{\mathrm{\λ}})\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{sin\φ}}_P\end{array}---\left(14\right)$

第四步，计算换极坐标系弹道参数

根据弹道再入点经度λ_I＝100°、再入点地心纬度φ_I＝20°、落点经度λ_T＝165°、落点地心纬度φ_T＝50°计算得到换极弹道参数,由计算得

第五步，计算换极坐标系弹道侧向包络

由式(15)计算换极弹道纵程

$\hat{L}=R_{e}arccos(sin\widehat{\mathrm{\φ}}sin{\widehat{\mathrm{\φ}}}_T+cos{\widehat{\mathrm{\φ}}}_{I}cos{\widehat{\mathrm{\φ}}}_{T}cos({\widehat{\mathrm{\λ}}}_T-{\widehat{\mathrm{\λ}}}_I\left)\right)---\left(15\right)$

计算带入式(16)，可计算得到换极弹道侧向包络宽度

$\hat{W}=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{\hat{L}}^i---\left(16\right)$

其中，模型系数a_i由第一步确定。

由计算的将换极弹道侧向包络描述为长度为宽度为的矩形，其边界由式(17)确定，

$\{\begin{array}{c}{\hat{E}}_{\mathrm{\λ}1}={\mathrm{\λ}}_I\\ {\hat{E}}_{\mathrm{\λ}2}=\frac{180\hat{L}}{{\mathrm{\πR}}_e}\end{array},\left\{\begin{array}{c}{\hat{E}}_{\mathrm{\φ}1}=-\frac{90\hat{W}}{{\mathrm{\πR}}_e}\\ {\hat{E}}_{\mathrm{\φ}2}=\frac{90\hat{W}}{{\mathrm{\πR}}_e}\end{array}\right.---\left(17\right)$

其中，为换极系侧向包络东向下边界，为换极系侧向包络东向上边界；为换极系侧向包络北向下边界，为换极系侧向包络北向上边界。

由计算得，

第六步，计算一般坐标系弹道三维包络

按照东向dλ＝1°、北向dφ＝1°的间隔将由式(17)确定的换极坐标系弹道侧向包络剖分为互不重叠的若干网格，记网格节点坐标为换极系弹道侧向包络节点见图10。

由第三步所述坐标转换关系，计算一般坐标系网格节点坐标N(λ_n,φ_n)，由N(λ_n,φ_n)确定的区域即为一般坐标系弹道侧向包络区域，一般系弹道侧向包络节点见图11。

考虑机动和偏差进行蒙特卡罗打靶，取弹道包络的垂向范围为沿参考弹道上下浮动10km，即ΔH＝20km。经过以上计算，最终可得给定再入点和落点弹道的三维弹道包络。

针对再入点与落点与第四步所述相同的临近空间大范围机动弹道进行仿真。图12、图13和图14分别给出了考虑禁飞区弹道、最大左侧机动弹道和最大右侧机动弹道与计算所得侧向包络的位置关系示意图。其中，禁飞区中心经度λ_nf＝138°、中心地心纬度φ_nf＝40°、高度H_nf＝100km、半径R_nf＝1000km。由图12、图13和图14可见，极端情况下的弹道均位于计算所得的包络内。

取n＝1.1，即以1.1倍W_io为实际侧向包络宽度W_i建立三维弹道包络数学模型。图15、图16和图17分别给出了考虑禁飞区弹道、最大左侧机动弹道和最大右侧机动弹道与计算所得包络的位置关系示意图，弹道仿真条件与图12、图13和图14中弹道相同。由图15、图16和图17可见，极端情况下的弹道均位于计算所得的包络内。

综合上述仿真结果可获得以下结论：

1)由本发明确定的方法进行临近空间大范围机动弹道三维包络计算，可使计算得到的包络覆盖具有不同机动特性的弹道，对考虑禁飞区的中程机动弹道以及最大侧向机动弹道等极端情况均有较好的适应性；

2)由本发明确定的方法进行弹道三维包络计算，模型参数一经固定，则无需再进行多次弹道计算来确定弹道可达域，避免了多次弹道积分带来的巨大运算量，因此具有计算速度快和轻量化的特征，能够满足射前任务快速规划和空间信息快速获取的要求；

3)建立模型所涉及的包络侧向宽度系数，确定包络特征点所涉及的空域剖分参数均可灵活调节，因此能够适应不同建模误差、不同计算时间要求和不同计算量要求下的弹道包络计算；

4)方法可广泛适用于具有相同动力学模型描述的飞行器，不受飞行器特征参数改变及飞行任务改变的影响，计算不出现奇点，具有适应范围广的特征。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (3)

1.一种临近空间大范围机动弹道三维包络计算方法，所述临近空间大范围机动弹道包括初始下降段弹道和滑翔段弹道，弹道起点为再入点，弹道终点为滑翔段弹道终点，其特征在于：所述计算方法包括以下步骤：

(1)建立弹道包络侧向宽度与纵程的数学模型；

(2)建立换极坐标系；

(3)建立换极坐标系中飞行器动力学模型；

(4)计算换极坐标系弹道参数；

(5)计算换极坐标系弹道侧向包络；

(6)计算一般坐标系弹道三维包络。

2.根据权利要求1所述的临近空间大范围机动弹道三维包络计算方法，其特征在于：

(1)建立弹道包络侧向宽度与纵程的数学模型

令弹道再入点为I，其经度为λ_I、地心纬度为φ_I；弹道落点为T，其经度为λ_T、地心纬度为φ_T；以由再入点和落点确定的再入大圆弧平面为对称面，称沿射向方向对称面以左弹道为左侧机动弹道，沿射向方向对称面以右弹道为右侧机动弹道；记最大左侧机动弹道距对称面的最大距离为左边界B_l，记最大右侧机动弹道距对称面的最大距离为右边界B_r；

根据给定的飞行器再入点飞行状态参数、弹道终端飞行状态参数、飞行过程约束条件及终端约束条件，针对再入点为(λ_I,φ_I)＝(0°,0°)，落点为(λ_ki,0°)(i＝1,2,3,4,…)的低空大范围机动弹道进行弹道计算，记各落点对应的弹道纵程L_i，由式(1)计算L_i，

L_i＝R_earccos(sinφ_Isinφ_I+cosφ_Icosφ_Icos(λ_T-λ_T))(1)

式中：R_e为地球半径；

通过改变再入点初始速度方位角获得最大左侧机动弹道和最大右侧向机动弹道，记各落点对应的弹道侧向包络宽度W_io，由式(2)计算W_io，

W_io＝B_l-B_r(2)

取n倍W_io为实际侧向包络宽度W_i，对L_i(i＝1,2,3,4,…)和W_i(i＝1,2,3,4,…)进行拟合，建立如式(3)所示的弹道包络侧向宽度与纵程的数学模型，确定模型系数a_i，

$W=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{L}^i---\left(3\right)$

(2)建立换极坐标系

A.定义一个再入大圆弧平面作为换极赤道平面：

a)对目标点确定的情况，将再入点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；

b)对于目标点未确定的情况，根据再入点位置及速度方位角确定的再入大圆弧平面作为换极赤道面；

B.基于换极赤道平面定义换极坐标系O_E为地心，轴沿再入点地心矢径方向，轴在换极赤道面内垂直于轴指向目标点方向，轴与轴、轴构成右手系。

(3)建立换极坐标系中飞行器动力学模型

根据步骤(2)建立的换极坐标系建立换极坐标系中飞行器的动力模型，且用表示换极坐标系中各物理量，用x表示一般坐标系中各物理量：

在换极坐标系中建立以时间为自变量的滑翔飞行器动力学方程，其飞行状态量为换极后的经度地心纬度航迹偏航角速度速度倾角和地心距

$\{\begin{array}{c}\frac{d\hat{r}}{dt}=\hat{V}\sin \widehat{\mathrm{\θ}}\\ \frac{d\widehat{\mathrm{\λ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}}{\hat{r}\cos \widehat{\mathrm{\φ}}}\\ \frac{d\widehat{\mathrm{\φ}}}{dt}=\frac{\hat{V}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}}{\hat{r}}\\ \frac{d\widehat{\mathrm{\σ}}}{dt}=\frac{L\sin \mathrm{\υ}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}\tan \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}-\hat{g}{\mathrm{\ω}}_e\frac{\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\σ}}}{\hat{V}\cos \widehat{\mathrm{\θ}}}+C_{\mathrm{\σ}}+{\tilde{C}}_{\mathrm{\σ}}\\ \frac{d\widehat{\mathrm{\θ}}}{dt}=\frac{L}{\hat{V}}\cos \mathrm{\υ}+\frac{\hat{V}}{\hat{r}}\cos \widehat{\mathrm{\θ}}+\frac{{\hat{g}}_r\cos \widehat{\mathrm{\θ}}}{\hat{V}}+\frac{\hat{g}{\mathrm{\ω}}_e}{\hat{V}}(\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}}-\cos \widehat{\mathrm{\σ}}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}})+C_{\mathrm{\θ}}+{\tilde{C}}_{\mathrm{\θ}}\\ \frac{d\hat{V}}{dt}=-D+{\hat{g}}_r\sin \widehat{\mathrm{\θ}}+\widehat{\mathrm{\θ}}{\mathrm{\ω}}_e(\cos \widehat{\mathrm{\σ}}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\φ}}+\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\φ}})+{\tilde{C}}_V\end{array}---\left(4\right)$

其中，C_σ、C_θ为哥氏加速度项，和为牵连加速度项，

$\left\{\begin{array}{c}{C}_{\mathrm{\σ}}=(2{\mathrm{\ω}}_{ex}-2\tan \widehat{\mathrm{\θ}}({\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\σ}}))\\ {\tilde{C}}_{\mathrm{\σ}}=-\frac{\hat{r}}{\hat{V}\cos \widehat{\mathrm{\θ}}}({\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}})\\ C_{\mathrm{\θ}}=2({\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\σ}})\\ {\tilde{C}}_{\mathrm{\θ}}=\frac{\hat{r}}{\hat{V}}\[{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\sin \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}+{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\sin \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\cos \widehat{\mathrm{\θ}}\]\\ {\tilde{C}}_V=\hat{r}\[-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ey}\cos \widehat{\mathrm{\θ}}\sin \widehat{\mathrm{\σ}}-{\mathrm{\ω}}_{ex}{\mathrm{\ω}}_{ez}\cos \widehat{\mathrm{\θ}}\cos \widehat{\mathrm{\σ}}+({\mathrm{\ω}}_{ey}^2+{\mathrm{\ω}}_{ez}^2)\sin \widehat{\mathrm{\θ}}\]\end{array}\right.---\left(5\right)$

其中，

$\left\{\begin{array}{c}{\mathrm{\ω}}_{ex}={\mathrm{\ω}}_e(\cos \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\sin \widehat{\mathrm{\λ}}\cos \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\sin \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\\ {\mathrm{\ω}}_{ey}={\mathrm{\ω}}_e(-\sin \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p+\cos \widehat{\mathrm{\λ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p)\\ {\mathrm{\ω}}_{ez}={\mathrm{\ω}}_e(-\cos \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{cosA}}_p-\sin \widehat{\mathrm{\λ}}\sin \widehat{\mathrm{\φ}}{\mathrm{cos\φ}}_p{\mathrm{sinA}}_p+\cos \widehat{\mathrm{\φ}}{\mathrm{sin\φ}}_p)\end{array}\right.---\left(6\right)$

其中，ω_e为地球旋转加速度矢量，λ_p和φ_p为换极后极点P的经度和地心纬度，A_p为P的方位角；

根据换极坐标系定义，一般坐标系与换极坐标系中地心距、当地速度倾角及速度的定义一致，

$r=\hat{r},\mathrm{\θ}=\widehat{\mathrm{\θ}},V=\hat{V}---\left(7\right)$

定义

$\begin{array}{c}\left[\begin{array}{ccc}{\mathrm{cos\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{cos\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\φ}}_f\\ -{\mathrm{sin\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{cos\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{cos\ψ}}_f{\mathrm{cos\φ}}_f\\ {\mathrm{cos\ψ}}_f{\mathrm{sin\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{cos\λ}}_f& -{\mathrm{cos\ψ}}_f{\mathrm{cos\λ}}_f-{\mathrm{sin\ψ}}_f{\mathrm{sin\φ}}_f{\mathrm{sin\λ}}_f& {\mathrm{sin\ψ}}_f{\mathrm{cos\φ}}_f\end{array}\right]\\ \stackrel{\mathrm{\Δ}}{=}\left[\begin{array}{ccc}{G}_{11}& G_{12}& G_{13}\\ G_{\mathrm{21}}& G_{\mathrm{22}}& G_{23}\\ G_{\mathrm{31}}& G_{\mathrm{32}}& G_{33}\end{array}\right]\end{array}---\left(8\right)$

其中，ψ_f为点F的方位角。

$\left\{\begin{array}{c}{G}_{11}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{12}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{13}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_1\\ G_{21}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{22}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{23}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_2\\ G_{31}\cos \mathrm{\φ}\cos \mathrm{\λ}+G_{32}\cos \mathrm{\φ}\sin \mathrm{\λ}+G_{33}\sin \mathrm{\φ}\stackrel{\mathrm{\Δ}}{=}{k}_3\end{array}\right.---\left(9\right)$

$\left\{\begin{array}{c}{G}_{11}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{21}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{31}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_1\\ G_{12}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{22}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{32}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_2\\ G_{13}\cos \widehat{\mathrm{\φ}}\cos \widehat{\mathrm{\λ}}+G_{23}\cos \widehat{\mathrm{\φ}}\sin \widehat{\mathrm{\λ}}+G_{33}\sin \widehat{\mathrm{\φ}}\stackrel{\mathrm{\Δ}}{=}{\tilde{k}}_3\end{array}\right.---\left(10\right)$

由一般坐标系中λ和φ确定换极坐标系中和的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\\ sin\widehat{\mathrm{\λ}}=k_1/\sqrt{k_1^2+k_2^2}\end{array}\right.& \left\{\begin{array}{c}\sin \widehat{\mathrm{\φ}}=k_3\\ \cos \widehat{\mathrm{\φ}}=\sqrt{k_1^2+k_2^2}\end{array}\right.\end{array}---\left(11\right)$

由换极坐标系中和确定一般坐标系中λ和φ的表达式为，

$\begin{array}{cc}\left\{\begin{array}{c}\cos \mathrm{\λ}={\tilde{k}}_1/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\\ sin\mathrm{\λ}={\tilde{k}}_2/\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.& \left\{\begin{array}{c}\sin \mathrm{\φ}={\tilde{k}}_3\\ \cos \mathrm{\φ}=\sqrt{{{\tilde{k}}_1}^2+{{\tilde{k}}_2}^2}\end{array}\right.\end{array}---\left(12\right)$

由一般坐标系中σ确定换极坐标系中的表达式为

$\widehat{\mathrm{\σ}}=\mathrm{\σ}+\mathrm{\η}---\left(13\right)$

其中，

$\left\{\begin{array}{c}\sin \mathrm{\η}=\frac{\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{cos\φ}}_P}{\cos \widehat{\mathrm{\φ}}}\\ \cos \mathrm{\η}=-\cos (A_P-\widehat{\mathrm{\λ}})\cos (\mathrm{\λ}-{\mathrm{\λ}}_P)+\sin (A_P-\widehat{\mathrm{\λ}})\sin (\mathrm{\λ}-{\mathrm{\λ}}_P){\mathrm{sin\φ}}_P\end{array}\right.---\left(14\right)$

(4)计算换极坐标系弹道参数

根据步骤(3)的坐标转换关系计算换极坐标系弹道参数，根据弹道再入点经度λ_I、再入点地心纬度φ_I、落点经度λ_T、落点地心纬度φ_T计算得到换极弹道参数,由计算得记换极弹道落点经度为

(5)计算换极坐标系弹道侧向包络

由式(15)计算换极弹道纵程

$\hat{L}=R_{e}arccos(\sin {\widehat{\mathrm{\φ}}}_I\sin {\widehat{\mathrm{\φ}}}_T+cos{\widehat{\mathrm{\φ}}}_{I}cos{\widehat{\mathrm{\φ}}}_{T}cos({\widehat{\mathrm{\λ}}}_T-{\widehat{\mathrm{\λ}}}_I\left)\right)---\left(15\right)$

其中，R_e为地球半径。

将代入式(16)，可计算得到换极弹道侧向包络宽度

$\hat{W}=\underset{i=0}{\overset{n-1}{\Σ}}{a}_i{\hat{L}}^i---\left(16\right)$

其中，模型系数a_i由第一步确定。

将换极弹道侧向包络描述为长度为宽度为的矩形，其边界由式(17)确定，

$\left\{\begin{array}{c}{\hat{E}}_{\mathrm{\λ}1}={\mathrm{\λ}}_I\\ {\hat{E}}_{\mathrm{\λ}2}=\frac{180\hat{L}}{{\mathrm{\πR}}_e}\end{array}\right.,\left\{\begin{array}{c}{\hat{E}}_{\mathrm{\φ}1}=-\frac{90\hat{W}}{{\mathrm{\πR}}_e}\\ {\hat{E}}_{\mathrm{\φ}2}=\frac{90\hat{W}}{{\mathrm{\πR}}_e}\end{array}\right.---\left(17\right)$

其中，为换极系侧向包络东向下边界，为换极系侧向包络东向上边界；为换极系侧向包络北向下边界，为换极系侧向包络北向上边界。

(6)计算一般坐标系弹道三维包络

按照东向dλ、北向dφ的间隔将由式(17)确定的换极坐标系弹道侧向包络剖分为互不重叠的若干网格，记网格节点坐标为

由步骤(3)所述坐标转换关系，计算一般坐标系网格节点坐标N(λ_n,φ_n)，由N(λ_n,φ_n)确定的区域即为一般坐标系弹道侧向包络区域。

3.根据权利要求2所述的临近空间大范围机动弹道三维包络计算方法，其特征在于：考虑机动和偏差进行蒙特卡罗打靶，以弹道高度范围的m倍作为弹道包络的垂向范围ΔH。