CN105354380B - 面向摄动因素影响补偿的滑翔弹道快速修正方法 - Google Patents

Abstract

本发明提供一种面向摄动因素影响补偿的滑翔弹道快速修正方法，包括：建立换极坐标系；换极坐标系中飞行器动力学模型建立；计算纵程、横程偏差；建立纵程、横程与滑翔弹道终端经纬度的关系式；计算纵程、横程关于落点经纬度的偏导数；计算弹道终端经纬度偏差；计算弹道终端经纬度偏差关于阻力加速度的偏导数；计算弹道终端经纬度偏差关于弹道诸元的偏导数；计算修正后弹道诸元；计算修正后弹道。由本发明方法进行滑翔弹道快速修正，可使修正后弹道横程、纵程偏差控制在1km左右，弹道精度提高约20倍，满足滑翔弹道快速修正的精度要求。且本发明能够针对多种摄动因素、适应各种条件下滑翔弹道的快速修正，计算不出现奇点，适应范围广。

Description

面向摄动因素影响补偿的滑翔弹道快速修正方法

技术领域

本发明属于飞行器弹道规划领域，特别涉及面向摄动因素影响补偿的高超声速滑翔飞行器滑翔弹道快速修正方法。

背景技术

滑翔弹道是高超声速滑翔飞行器全程弹道的主要组成部分，具有显著区别于其他飞行器弹道的复杂弹道特性。由于飞行器本体及地球物理环境因素的建模复杂性以及对地球物理摄动因素的物理认知不完善，通常基于简化模型建立弹道规划方法，由模型近似引入的误差统称为摄动因素。摄动因素将导致飞行器实际运动状态偏离设计状态，是弹道规划及制导控制设计中方法误差的主要来源。摄动因素主要包括飞行器本体建模误差和地球物理环境因素误差两类。其中，飞行器本体建模误差主要包括弹体结构误差、发动机特性误差和气动系数误差等，地球物理环境因素误差主要包括大气模型不确定性及扰动引力等。摄动因素的影响不可忽略。

以飞行过程中扰动引力为例，其引起的滑翔段终端位置偏差可达几十公里，远超末制导的有效半径，如不考虑补偿则可能使飞行器不能完成既定飞行任务。

从现有工作来看，鲜有滑翔弹道快速规划方法纳入对摄动因素影响的修正，相近领域的理论和方法难以直接应用。总体上，现有理论和方法与待解决问题存在如下差距：1)针对不确定因素的分析手段多通过对不确定项的冗余设计和对提高弹道规划算法的鲁棒性来削弱弹道规划方法对不确定因素的敏感度，不仅加重了对弹道规划算法鲁棒性和可靠性设计的负担，并且当影响超出方法能力时将导致规划任务不能完成；2)不确定因素和扰动引力等不可测量量有本质区别，针对前者建立的方法不适用于后者；3)尚无工作采用针对滑翔弹道的弹道规划参数进行摄动因素影响补偿的思路，因此未能将对摄动因素影响的补偿嵌入到滑翔飞行器的弹道规划方法中。

针对上述问题，基于D-E剖面跟踪的滑翔弹道规划方法，提出一种面向摄动因素补偿的弹道快速修正方法。方法的基本思想为采用等效补偿策略将摄动因素的影响量归算到弹道诸元中，通过调整弹道诸元实现对摄动因素影响的补偿。为计算弹道诸元补偿量，方法首先基于状态空间摄动理论建立了摄动因素沿滑翔弹道的误差传播模型，获得了滑翔弹道终端经纬度偏差关于摄动因素的半解析表达式；其次建立弹道横程、纵程偏差与弹道诸元的数学关系，进而确定诸元修正量。方法具有计算速度快、补偿精度高的特征，是首个考虑了摄动因素影响补偿的滑翔弹道快速修正方法，由于能够解析地表征摄动因素与弹道诸元的关系，因此从根本上补偿了摄动因素对弹道规划精度的影响，为深化滑翔弹道误差传播机理认识、建立高精度弹道快速规划方法提供了理论基础和方法支撑。

发明内容

本发明针对高超声速滑翔飞行器滑翔弹道快速修正问题，首次提出一种面向摄动因素影响补偿的滑翔弹道快速修正方法。该方法综合利用摄动理论、坐标系转换、误差传播分析和飞行器动力学建模等实现弹道快速修正方法的建立。建立的方法立足于传统的D-E剖面跟踪弹道规划方法，基于等效补偿的思路将摄动因素的影响量归算到弹道诸元修正量中，可广泛适用于具有相同动力学描述的飞行器弹道，可用于补偿任意摄动因素对滑翔弹道的影响。

本发明的弹道快速修正方法构建思路具体包括以下步骤：第一步，建立换极坐标系；第二步，换极坐标系中飞行器动力学模型的建立；第三步，计算纵程、横程偏差；第四步，建立纵程、横程与滑翔弹道终端经纬度的关系式；第五步，计算纵程、横程关于落点经纬度的偏导数；第六步，计算弹道终端经纬度偏差；第七步，计算弹道终端经纬度偏差关于阻力加速度的偏导数；第八步，计算弹道终端经纬度偏差关于弹道诸元的偏导数；第九步，计算修正后的弹道诸元；第十步，计算修正后的弹道。

本发明中，换极是极点变换的意思，换极坐标系是指基于极点变换思想以重新定义的换极赤道面为基准建立的坐标系，而一般坐标系是指以地球赤道面为基准建立的坐标系。

具体地，其中第三步为由参考弹道和干扰弹道计算纵程、横程偏差；第九步为建立纵程、横程偏差和弹道诸元补偿量的关系，计算修正后的弹道诸元；第十步为基于修正后弹道诸元计算修正后的弹道，并通过坐标变换获得一般坐标系中的修正弹道。

本发明的技术方案主要包括以下步骤：

第一步，建立换极坐标系

首先引入一个换极坐标系。为表述方便，用X表示换极坐标系中各物理量，用表示一般坐标系中各物理量。用于构建弹道快速修正方法的各步骤均在换极坐标系下完成，以下不再赘述。

按如下方式建立换极坐标系：

①定义一个再入大圆弧平面作为换极赤道平面：1)对目标点确定的情况，将滑翔起点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；2)对于目标点未确定的情况，根据滑翔起点位置及方位角确定的再入大圆弧平面作为换极赤道面。

②基于换极赤道平面定义换极坐标系O_E-XYZ：O_E为地心，X轴沿滑翔起点地心矢径方向，Y轴在换极赤道面内垂直于X轴指向目标点方向，Z轴与X轴、Y轴构成右手系。

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

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

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

其中，

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

第三步，计算纵程、横程偏差

不考虑摄动因素计算一条滑翔弹道，称为参考弹道。其由滑翔弹道起点到滑翔弹道终点的纵程为L^*，弹道终端偏离目标终端的横程为H^*。

考虑摄动因素计算一条滑翔弹道，称为干扰弹道。其由滑翔弹道起点到滑翔弹道终点的纵程为L，弹道终端偏离目标终端的横程为H。

由摄动因素引起的滑翔弹道纵程偏差ΔL和横程偏差ΔH分别为，

第四步，建立纵程、横程与滑翔弹道终端经纬度的关系式

定义弹道起点为F(经度λ_f，地心纬度φ_f)，参考弹道终点为M(经度λ_m，地心纬度φ_m)，干扰弹道终点为C(经度λ_c，地心纬度φ_c)，参考航程角为β₀，实际航程角为β′，纵程角为β，横程角为ζ，则

干扰弹道终点与参考弹道终点之间的角距为，

cosζ′＝sinφ_m sinφ_c+cosφ_m cosφ_c cos(λ_c-λ_m)

(60)

干扰弹道终点与参考弹道终点相对弹道起点F的张角为，

纵程角β和横程角ζ为

纵程L和横程H为，

第五步，计算纵程、横程关于落点经纬度的偏导数

根据第四步建立的纵程、横程与滑翔弹道终端经纬度的关系式，可推导得纵程、横程关于落点经纬度的偏导数，

其中，

第六步，计算弹道终端经纬度偏差

在换极坐标系中建立以飞行状态偏差量为状态变量的摄动方程如下，

对式(16)进行一次积分即可求解状态转移矩阵Φ(t_k,t)的伴随矩阵G(t,t_k)，进而通过式(17)求解Φ(t_k,t)。

Φ(t_k,t)＝G^T(t,t_k) (71)

令

求解式(70)得，

G₁₁＝1,G₁₂＝G₁₃＝G₁₄＝G₁₅＝G₁₆＝0 (73)

G₂₂＝1,G₂₁＝G₂₃＝G₂₄＝G₂₅＝G₂₆＝0 (74)

摄动方程(69)的通解为，

不考虑初态误差，则由式(79)可得换极坐标系下终端经纬度偏差为，

第七步，计算弹道终端经纬度偏差关于阻力加速度的偏导数

求(80)关于阻力加速度D的偏导数，

对式(81)做如下变换，

高超声速滑翔飞行器能量E关于时间t的偏导数为，

阻力加速度D关于能量E的偏导数为，

其中，ρ为大气密度，S_r为飞行器参考面积，C_D为阻力系数，M为飞行器质量。

将式(83)和式(84)代入式(82)可得，

式(80)求关于时间t的偏导数为，

基于瞬时平衡假设，可推导得，

将式(87)～式(94)代入式(86)，再代入式(85)，可得终端经纬度偏差关于阻力加速度的偏导数。

第八步，计算弹道终端经纬度偏差关于弹道诸元的偏导数

以分段函数形式描述D-E剖面：

其中C₁和C₂为待修正的弹道诸元。

求式(95)关于C₁和C₂的偏导数，

将式(96)拟合为函数Eⁿ(n＝0,1,2,…)的线性组合，使其具有形式统一的数学描述，

其中，a_i(i＝0,1,2)和b_i(i＝1,3,5,7,9,11)为拟合系数。

联立式(97)和式(85)，可得弹道终端经纬度偏差关于弹道诸元的偏导数，

第九步，计算修正后弹道诸元

滑翔弹道纵程、横程偏差与弹道诸元补偿量的关系式为，

则基于第三步计算得到的纵程、横程偏差和J₁和J₂的推导结果，可得

记不考虑摄动因素进行弹道规划得到的参考弹道诸元为C₁ ^*和C₂ ^*。则考虑摄动因素影响，修正后弹道诸元C₁ ^c和C₂ ^c为，

第十步，计算修正后弹道

将修正后诸元代入式(95)，可得修正后D-E剖面，

通过跟踪修正后D-E剖面，即可得到换极坐标系下的修正后弹道。

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

定义

其中，和分别为点F在一般坐标系中的经度、纬度和方位角。

由换极系下经纬度λ和φ确定一般系下经纬度和的表达式为

由换极系下航迹偏航角σ确定一般系下航迹偏航角的表达式为

其中，

经过以上变换，最终可得到一般坐标系中修正后滑翔弹道。

至此，经过上述十步可最终建立一种面向摄动因素影响补偿的滑翔弹道快速修正方法，该方法计算速度快、补偿精度高，能够满足滑翔弹道快速修正的计算轻量化和补偿精度要求。与现有研究基础相比，本发明提出的方法具有以下优点：

1)首次在高超声速滑翔飞行器滑翔弹道计算中考虑了对摄动因素影响的修正，开创性地提出了一种考虑摄动因素影响补偿的滑翔弹道快速修正方法，方法能够解析地表征摄动因素与弹道诸元的关系，因此可从根本上补偿摄动因素对弹道规划精度的影响，为深化滑翔弹道误差传播机理认识、建立高精度弹道快速规划方法提供了理论基础和方法支撑。

2)是首个针对传统D-E剖面跟踪弹道规划方法进行摄动因素影响补偿的弹道快速修正方法。方法立足于D-E剖面跟踪方法的本质，首次以剖面调节参数作为待修正的弹道诸元，通过调整弹道诸元来调整D-E剖面，最终实现对滑翔弹道的修正。

3)首次将等效补偿的思想引入滑翔弹道修正方法的建立中，将摄动因素对关键指标的影响量归算到待修正的弹道诸元中，通过计算弹道诸元补偿量实现对弹道的解析修正，避免了多次弹道规划中涉及的大量弹道积分导致的不可容忍的计算量。

4)基于极点变换的思想，首次将一种换极坐标系引入弹道快速修正方法的建立中，并推导了换极坐标系中的飞行器动力学方程，将滑翔弹道再入平面转移到换极赤道面附近，从而简化了弹道规划过程以及摄动方程状态转移矩阵的求解。

5)首次借助几何辅助原理，在由滑翔弹道起点和终点确定的球面三角形内建立弹道纵程、横程与弹道起点和终点经纬度的关系式。基于弹道摄动理论和换极坐标系建立滑翔弹道的误差传播模型，从而建立摄动因素和弹道终端经纬度偏差的半解析关系式。

6)方法可广泛适用于具有相同动力学描述的飞行器，不受飞行器特征参数改变的影响。方法所针对的摄动因素基于相同的补偿方法，涵盖飞行器本体建模误差和地球物理环境摄动因素，可适用于任意单一摄动因素或复杂摄动因素对弹道特性的影响补偿。

7)该方法计算速度快、补偿精度高，可使修正后滑翔弹道的横程、纵程偏差控制在1km以内，弹道精度提高约20倍。修正后弹道满足热流密度、过载、动压等过程约束条件。方法适用于任意滑翔弹道，不出现奇点。

附图说明

图1为换极坐标系及飞行状态量定义示意图；

图2为滑翔弹道起点F和换极后极点P的位置关系示意图；

图3为纵程与横程计算示意图；

图4为D-E剖面示意图；

图5为换极坐标系航迹偏航角的位置关系示意图；

图6为典型滑翔弹道参考D-E剖面示意图；

图7为不考虑扰动引力的典型滑翔弹道曲线及实际D-E剖面示意图；

图8为沿典型滑翔弹道的扰动引力示意图；

图9为考虑扰动引力的典型滑翔弹道曲线及实际D-E剖面示意图；

图10为典型滑翔弹道修正后D-E剖面示意图；

图11为修正后的典型滑翔弹道曲线及实际D-E剖面示意图；

图12为平原地区沿滑翔弹道扰动引力示意图；

图13为丘陵地区沿滑翔弹道扰动引力示意图；

图14为特大山区沿滑翔弹道扰动引力示意图。

具体实施方式

下面结合具体实例，对本发明作进一步的说明：

基于CAV-H飞行器模型，以飞行过程中扰动引力为单一摄动因素进行仿真。采用1080阶球谐函数法计算沿滑翔弹道的扰动引力。

仿真初始条件设置为：

(1)针对某典型滑翔弹道进行方法精度分析：①初始点状态参数：速度V₀＝6500m/s，速度倾角θ₀＝0°，高度H₀＝80km，经度λ₀＝0°，地心纬度φ₀＝0°；②结束点状态参数：速度V_k＝2500m/s，速度倾角θ_k＝0°，高度H_k＝30km，经度λ_k＝55°，地心纬度φ_k＝30°；③飞行约束条件：最大热流密度最大动压q_max＝100kPa，最大过载n_max＝3g；④终端结束时距目标点的待飞航程：S_togo＝100km。

(2)针对位于不同区域、具有不同射程、东向发射的系列滑翔弹道进行方法适应性分析：①弹道起点位于平原地区(经度λ_p0＝115°E，地心纬度φ_p0＝35°N)：射程为5000km弹道(终端经度λ_pk1＝166°E，终端地心纬度φ_pk1＝24°N)；射程为6000km弹道(终端经度λ_pk2＝176°E，终端地心纬度φ_pk2＝20°N)；射程为7000km弹道(终端经度λ_pk3＝183°E，终端地心纬度φ_pk3＝15°N)；②弹道起点位于丘陵地区(经度为λ_h0＝110°E，地心纬度为φ_h0＝27°N)：射程为5000km弹道(终端经度λ_hk1＝160°E，终端地心纬度φ_hk1＝18°N)；射程为6000km弹道(终端经度λ_hk2＝168°E，终端地心纬度φ_hk2＝15°N)；射程为7000km弹道(终端经度λ_hk3＝175°E，终端地心纬度φ_hk3＝12°N)；③弹道起点位于特大山区(经度为λ_m0＝80°E，地心纬度为φ_m0＝42°N)：射程为5000km弹道(终端经度λ_mk1＝134°E，终端地心纬度φ_mk1＝28°N)；射程为6000km弹道(终端经度λ_mk2＝142°E，终端地心纬度φ_mk2＝23°N)；射程为7000km弹道(终端经度λ_mk3＝149°E，终端地心纬度φ_mk3＝17°N)。其他仿真条件不变。

仿真计算机配置为：Intel(R)Core(TM)i5-3470CPU@3.20GHz，内存为3.46GB。软件环境为Window XP操作系统，计算程序基于VC++6.0开发。

其具体步骤如下：

第一步，建立换极坐标系

首先引入一个换极坐标系，见图1。为表述方便，用X表示换极坐标系中各物理量，用表示一般坐标系中各物理量。用于构建弹道快速修正方法的各步骤均在换极坐标系下完成，以下不再赘述。

按如下方式建立换极坐标系：

①定义一个再入大圆弧平面作为换极赤道平面：1)对目标点确定的情况，将滑翔起点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；2)对于目标点未确定的情况，根据滑翔起点位置及方位角确定的再入大圆弧平面作为换极赤道面。

②基于换极赤道平面定义换极坐标系O_E-XYZ：O_E为地心，X轴沿滑翔起点地心矢径方向，Y轴在换极赤道面内垂直于X轴指向目标点方向，Z轴与X轴、Y轴构成右手系。

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

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

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

其中，

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

第三步，计算纵程、横程偏差

不考虑摄动因素计算一条滑翔弹道，称为参考弹道。其由滑翔弹道起点到滑翔弹道终点的纵程为L^*，弹道终端偏离目标终端的横程为H^*。

考虑摄动因素计算一条滑翔弹道，称为干扰弹道。其由滑翔弹道起点到滑翔弹道终点的纵程为L，弹道终端偏离目标终端的横程为H。

由摄动因素引起的滑翔弹道纵程偏差ΔL和横程偏差ΔH分别为，

第四步，建立纵程、横程与滑翔弹道终端经纬度的关系式

纵程横程计算示意图见图3。定义弹道起点为F(经度λ_f，地心纬度φ_f)，参考弹道终点为M(经度λ_m，地心纬度φ_m)，干扰弹道终点为C(经度λ_c，地心纬度φ_c)，参考航程角为β₀，实际航程角为β′，纵程角为β，横程角为ζ，则

干扰弹道终点与参考弹道终点之间的角距为，

cosζ′＝sinφ_m sinφ_c+cosφ_m cosφ_c cos(λ_c-λ_m) (114)

干扰弹道终点与参考弹道终点相对弹道起点F的张角为，

纵程角β和横程角ζ为

纵程L和横程H为，

第五步，计算纵程、横程关于落点经纬度的偏导数

根据第四步建立的纵程、横程与滑翔弹道终端经纬度的关系式，可推导得纵程、横程关于落点经纬度的偏导数，

其中，

第六步，计算弹道终端经纬度偏差

在换极坐标系中建立以飞行状态偏差量为状态变量的摄动方程如下，

对式(70)进行一次积分即可求解状态转移矩阵Φ(t_k,t)的伴随矩阵G(t,t_k)，进而通过式(71)求解Φ(t_k,t)。

Φ(t_k,t)＝G^T(t,t_k) (125)

令

求解式(124)得，

G₁₁＝1,G₁₂＝G₁₃＝G₁₄＝G₁₅＝G₁₆＝0 (127)

G₂₂＝1,G₂₁＝G₂₃＝G₂₄＝G₂₅＝G₂₆＝0 (128)

摄动方程(123)的通解为，

不考虑初态误差，则由式(133)可得换极坐标系下终端经纬度偏差为，

以扰动引力为单一摄动因素，则摄动项为，

其中，δg_r为沿地球矢径方向的扰动引力分量，为沿地球自转方向的扰动引力分量。

则式(134)可化为，

第七步，计算弹道终端经纬度偏差关于阻力加速度的偏导数

求式(136)关于阻力加速度D的偏导数，

对式(137)做如下变换，

高超声速滑翔飞行器能量E关于时间t的偏导数为，

阻力加速度D关于能量E的偏导数为，

其中，ρ为大气密度，S_r为飞行器参考面积，C_D为阻力系数，M为飞行器质量。

将式(139)和式(140)代入式(138)可得，

式(136)求关于时间t的偏导数为，

基于瞬时平衡假设，可推导得，

将式(143)～式(148)代入式(142)，再代入式(141)，可得终端经纬度偏差关于阻力加速度的偏导数。

第八步，计算弹道终端经纬度偏差关于弹道诸元的偏导数

以分段函数形式描述D-E剖面，见图4：

其中C₁和C₂为待修正的弹道诸元。

对能量进行归一化处理，并取E₀＝0，E₁＝0.3，E₂＝0.5，E₃＝0.7，E₄＝0.8，E_f＝1，可得D(E)关于C₁和C₂的偏导数，

将式(150)拟合为函数Eⁿ(n＝0,1,2,…)的线性组合，使其具有形式统一的数学描述，

联立式(151)和式(141)，可得弹道终端经纬度偏差关于弹道诸元的偏导数，

第九步，计算修正后弹道诸元

滑翔弹道纵程、横程偏差与弹道诸元补偿量的关系式为，

则基于第三步计算得到的纵程、横程偏差和J₁和J₂的推导结果，可得

记不考虑摄动因素进行弹道规划得到的参考弹道诸元为C₁ ^*和C₂ ^*。则考虑摄动因素影响，修正后弹道诸元C₁ ^c和C₂ ^c为，

第十步，计算修正后弹道

将修正后诸元代入式(149)，可得修正后D-E剖面，

通过跟踪修正后D-E剖面，即可得到换极坐标系下的修正后弹道。

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

定义

其中，和分别为点F在一般坐标系中的经度、纬度和方位角。

由换极系下经纬度λ和φ确定一般系下经纬度和的表达式为

由换极系下航迹偏航角σ确定一般系下航迹偏航角的表达式为(见图5)，

其中，

经过以上变换，最终可得到一般坐标系中修正后滑翔弹道。

图6给出了某典型滑翔弹道不考虑扰动引力的参考D-E剖面。由图6可见，参考剖面满足分段函数的数学描述，其对应的参考弹道诸元为图7给出了不考虑扰动引力的影响，根据弹道规划结果进行开环弹道仿真得到的参考滑翔弹道曲线及其实际阻力加速度剖面。由图7可见，通过跟踪图6所示的D-E剖面得到的实际飞行弹道能够满足终端位置约束，其实际阻力加速度剖面能较好地跟踪设计剖面，且满足热流密度、过载和动压等过程约束条件。

图8给出了利用1080阶球谐函数方法计算得到的沿某典型滑翔弹道的扰动引力。由图8可见，由于飞行高度随飞行时间逐渐下降，扰动引力随飞行时间呈总体上升的趋势，且具有较为复杂的变化特性。其中，最大值在接近滑翔弹道终端的低空区域，约为72.849mgal；最小值在接近弹道起始点的高空区域，约为1.562mgal。

图9给出了考虑扰动引力的影响，根据弹道规划结果进行开环弹道仿真计算得到的干扰滑翔弹道曲线及其实际阻力加速度剖面。由图9可见，由于扰动引力的影响，干扰弹道逐渐偏离参考弹道，两者之间的偏差随飞行时间累积。

表1 弹道规划诸元

	C<sub>1</sub>	C<sub>2</sub>
			参考诸元	2.449936	3.733374
诸元修正量	-0.128264	0.309661
			修正后诸元	2.578200	3.423713

表1给出了参考弹道诸元以及根据弹道快速修正方法计算得到弹道规划诸元修正量和修正后诸元。图10给出了修正后D-E飞行剖面。

图11给出了根据修正后的D-E飞行剖面进行弹道仿真得到的修正后滑翔弹道曲线及其实际阻力加速度剖面。由图11可见，修正弹道终端较干扰弹道更为接近参考弹道，其实际飞行的阻力加速度剖面能够满足热流密度、过载和动压等过程约束条件。

表2给出了参考弹道、干扰弹道和修正弹道的终端位置，以及干扰弹道和修正弹道的横程、纵程偏差。由表2可见，由扰动引力引起的滑翔弹道终端纵程偏差为-27662m，终端横程偏差为2861m，终端高度偏差为-2.9m，偏差较大而不容忽略。基于弹道快速修正方法进行修正后，修正弹道的终端纵程偏差为-1564m，终端横程偏差为382m，终端高度偏差为-11.6m。修正后弹道终端位置精度提高了约20倍。

表2 终端状态对比

	经度(deg)	纬度(deg)	高度(km)	纵程偏差(m)	横程偏差(m)
						参考弹道	55.2840	29.0019	29.9489	-	-
干扰弹道	55.0245	28.8982	29.9460	-27662	2861
						修正弹道	55.2701	28.9942	29.9373	-1564	382

图12～图14分别给出了沿弹道起点分别位于平原地区、丘陵地区和特大山区，射程分别为5000km、6000km和7000km滑翔弹道的扰动引力示意图。总体而言，特大山区平均扰动引力最大，平原地区平均扰动引力最小。

表3～表8分别给出了滑翔弹道起点位于平原地区、丘陵地区和特大山区，射程为5000km、6000km和7000km滑翔弹道的弹道诸元和弹道终端状态。计算结果表明：1)扰动引力对弹道的影响随飞行时间累积，因此射程越长，终端偏差越大；2)扰动引力受地形影响，特大山区平均扰动引力最大，平原地区平均扰动引力最小，因此位于特大山区的滑翔弹道具有较大的终端状态偏差；3)本发明的方法能较好地修正扰动引力对弹道终端偏差的影响，修正后弹道终端状态偏差控制在1km以内，弹道精度提高约20倍。

表3 平原地区弹道规划诸元

表4 平原地区滑翔弹道终端状态对比

表5 丘陵地区弹道规划诸元

表6 丘陵地区滑翔弹道终端状态对比

表7 特大山区弹道规划诸元

表8 特大山区滑翔弹道终端状态对比

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

1)由本发明确定的方法进行滑翔弹道快速修正，可使修正后弹道横程、纵程偏差控制在1km左右，弹道精度提高约20倍，满足滑翔弹道快速修正的精度要求；

2)由本发明确定的方法进行滑翔弹道快速修正，避免了精确模型的引入，也避免了多次弹道规划带来的大量弹道积分运算，因此具有计算速度快和轻量化的特征，具有广阔的工程应用前景；

3)由本发明建立的弹道快速修正方法，可解析地表征摄动因素与待修正弹道诸元之间的数学关系，有助于加深摄动因素误差传播的机理认识，从本质上抑制摄动因素对弹道精度的影响；

本发明能够针对多种摄动因素、可适应各种条件下滑翔弹道的快速修正，计算不出现奇点，具有适应范围广的特征。

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

Claims (2)

1.一种面向摄动因素影响补偿的滑翔弹道快速修正方法，包括设置仿真初始条件，(1)针对典型滑翔弹道进行方法精度分析，(2)针对位于不同区域、具有不同射程、东向发射的系列滑翔弹道进行方法适应性分析；且具体包括以下步骤：第一步，建立换极坐标系；第二步，换极坐标系中飞行器动力学模型的建立；第三步，由参考弹道和干扰弹道计算纵程、横程偏差；第四步，建立纵程、横程与滑翔弹道终端经纬度的关系式；第五步，计算纵程、横程关于落点经纬度的偏导数；第六步，计算弹道终端经纬度偏差；第七步，计算弹道终端经纬度偏差关于阻力加速度的偏导数；第八步，计算弹道终端经纬度偏差关于弹道诸元的偏导数；第九步，建立纵程、横程偏差和弹道诸元补偿量的关系，计算修正后的弹道诸元；第十步，基于修正后弹道诸元计算修正后的弹道，并通过坐标变换获得一般坐标系中的修正弹道；其中，

第一步，建立换极坐标系

按如下方式建立换极坐标系：

①定义一个再入大圆弧平面作为换极赤道平面：1)对目标点确定的情况，将滑翔起点和目标点地心矢径构成的再入大圆弧平面作为换极赤道平面；2)对于目标点未确定的情况，根据滑翔起点位置及方位角确定的再入大圆弧平面作为换极赤道面；

②基于换极赤道平面定义换极坐标系O_E-XYZ：O_E为地心，X轴沿滑翔起点地心矢径方向，Y轴在换极赤道面内垂直于X轴指向目标点方向，Z轴与X轴、Y轴构成右手系；

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

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

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

其中，

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

第三步，计算纵程、横程偏差

不考虑摄动因素计算一条滑翔弹道，称为参考弹道；其由滑翔弹道起点到滑翔弹道终点的纵程为L^*，弹道终端偏离目标终端的横程为H^*；

考虑摄动因素计算一条滑翔弹道，称为干扰弹道；其由滑翔弹道起点到滑翔弹道终点的纵程为L，弹道终端偏离目标终端的横程为H；

由摄动因素引起的滑翔弹道纵程偏差ΔL和横程偏差ΔH分别为，

第四步，建立纵程、横程与滑翔弹道终端经纬度的关系式

定义弹道起点为F，其经度λ_f，地心纬度φ_f；参考弹道终点为M，其经度λ_m，地心纬度φ_m；干扰弹道终点为C，其经度λ_ζ，地心纬度φ_c；参考航程角为β₀，实际航程角为β′，纵程角为β，横程角为Z，则

干扰弹道终点与参考弹道终点之间的角距为，

cosζ′＝sinφ_msin_φ+cosφ_mcos_φcos(λ_φ-λ_m) (6)

干扰弹道终点与参考弹道终点相对弹道起点F的张角为，

纵程角β和横程角ζ为

纵程L和横程H为，

；第五步，计算纵程、横程关于落点经纬度的偏导数

根据第四步建立的纵程、横程与滑翔弹道终端经纬度的关系式，可推导得纵程、横程关于落点经纬度的偏导数，

其中，

第六步，计算弹道终端经纬度偏差

在换极坐标系中建立以飞行状态偏差量为状态变量的摄动方程如下，

对式(16)进行一次积分即可求解状态转移矩阵Φ(t_k，t)的伴随矩阵G(t，t_k)，进而通过式(17)求解Φ(t_k，t)；

Φ(t_k，t)＝G^T(t，r_k) (17)

令

求解式(16)得，

G₁₁＝1，G₁₂＝G₁₃＝G₁₄＝G₁₅＝G₁₆＝0 (19)

G₂₂＝1：G₂₁＝G₂₃＝G₂₄＝G₂₅＝G₂₆＝0 (20)

摄动方程(15)的通解为，

不考虑初态误差，则由式(25)可得换极坐标系下终端经纬度偏差为，

第七步，计算弹道终端经纬度偏差关于阻力加速度的偏导数

求(26)关于阻力加速度D的偏导数，

对式(27)做如下变换，

高超声速滑翔飞行器能量E关于时间t的偏导数为，

阻力加速度D关于能量E的偏导数为，

其中，ρ为大气密度，S_r为飞行器参考面积，C_D为阻力系数，M为飞行器质量；

将式(29)和式(30)代入式(28)可得，

式(26)求关于时间t的偏导数为，

基于瞬时平衡假设，可推导得，

将式(33)～式(40)代入式(32)，再代入式(31)，可得终端经纬度偏差关于阻力加速度的偏导数；

第八步，计算弹道终端经纬度偏差关于弹道诸元的偏导数

以分段函数形式描述D-E剖面：

其中C₁和C₂为待修正的弹道诸元；

求式(41)关于C₁和C₂的偏导数，

将式(42)拟合为函数Eⁿ(n＝0，1，2，…)的线性组合，使其具有形式统一的数学描述，

其中，a_i(i＝0，1，2)和b_i(i＝1，3，5，7，9，11)为拟合系数；

联立式(43)和式(31)，可得弹道终端经纬度偏差关于弹道诸元的偏导数，

第九步，计算修正后的弹道诸元

滑翔弹道纵程、横程偏差与弹道诸元补偿量的关系式为，

则基于第三步计算得到的纵程、横程偏差和J₁和J₂的推导结果，可得

记不考虑摄动因素进行弹道规划得到的参考弹道诸元为C₁ ^*和C₂ ^*；则考虑摄动因素影响，修正后弹道诸元C₁ ^c和C₂ ^c为，

第十步，计算修正后的弹道

将修正后诸元代入式(41)，可得修正后D-E剖面，

通过跟踪修正后D-E剖面，即可得到换极坐标系下的修正后弹道；

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

定义

其中，和分别为点F在一般坐标系中的经度、纬度和方位角；

由换极系下经纬度λ和φ确定一般系下经纬度和的表达式为

由换极系下航迹偏航角σ确定一般系下航迹偏航角的表达式为

其中，

经过以上变换，最终可得到一般坐标系中修正后滑翔弹道。

2.根据权利要求1所述的修正方法，其特征在于：所述摄动因素包括飞行器本体建模误差和地球物理环境因素建模误差。