载人飞船返回舱外测高精度弹道计算与落点预报

2023-11-10段成林盛庆轩谢剑锋徐海涛慎千慧

载人航天 2023年5期

段成林，盛庆轩，谢剑锋，陈明，徐海涛，陈铭，慎千慧

（北京航天飞行控制中心，北京 100094）

1 引言

2022 年，中国已全面建成空间站，承担天地往返运输的载人飞船发射会更加频繁。为进一步提高返回任务执行效率，提升航天员舒适度，飞船系统实施了快速返回策略，载人飞船从正常撤离至轨返分离的飞行圈次从11 圈缩减为5 圈，返回时间缩短约9 h。为保障航天员安全着陆，赢得搜救时间，对返回舱进行准确的落点预报至关重要［1-2］。

返回舱再入后的高精度弹道计算结果可极大提高落点预报精度。再入后实时弹道计算有3 个数据源：返回舱GNC 惯性测量数据、GNSS 导航数据和地面光学、雷达等外测数据。根据GNC 数据估计的返回舱弹道在再入段由于无法对角速率陀螺的零位漂移进行实时校正，随着高度的下降，弹道估计偏差迅速增大［3-5］；GNSS 测量数据精度较高，但存在导航异常时卫星不定位、进入黑障后通讯中断等情况；地面光学、雷达测量数据精度高、数据可靠，是返回舱再入后弹道计算的主要数据源，是确保回收着陆任务成功的保底测控手段。实际返回舱高速再入过程中，经常存在速度变化快导致目标失锁、数据中断、跟踪设备多而造成数据融合困难等现象，外测弹道计算普遍存在精度低、稳定性差的问题。

目前，外测弹道计算主要采用单点定位或滤波方法。单点定位通过微分平滑求出速度信息，由于参与平滑点数少，导致精度较差；滤波方法在观测不充分时，收敛较慢，甚至可能出现发散现象。随着计算机性能的提高，将累积一段时间的数据用于事后处理的大数据融合方法实时化已逐渐成为可能［6-8］。北京航天飞行控制中心采用基于多项式拟合的实时弹道动态滑窗处理方法，将一段弹道数据联合处理，能够充分降低观测噪声，提高弹道计算精度。

对于使用不同测量设备和不同算法确定的多条测元弹道，需要采用弹道优选方法实时挑选出一条最符合客观实际的弹道。传统弹道优选方法通常根据经验事先指定某些测站组合和弹道为最优弹道［9］，但挑选的某组弹道通常并不是全程最优。优选弹道切换过程一般都伴随较大的跳变，而未被选择的弹道许多有效信息被丢弃。针对上述缺陷，本文提出了一种结合弹道置信度进行弹道加权融合计算的实时弹道自动优选方法，能够有效避免因某一测元弹道波动或间断引起的优选弹道波动，以提高优选弹道的可靠性和精度。

2 外测弹道计算方法

2.1 光学交会计算目标位置

利用光电经纬仪解算弹道至少需要2 台设备，采用方向余弦法求出位置信息，用微分平滑求出速度信息。方向余弦法将2 台光测设备的方位角和俯仰角的观测数据转换成各自的观测站与目标间向量的方向余弦，再根据正弦定理计算出各自到目标的斜距和弹道位置参数［10］。

步骤1：取2 台经纬仪的方位角Ai和俯仰角Ei（i＝1，2）的观测值，计算各自观测点与目标间向量在地固坐标系中的方向余弦li、mi、ni，得到式（1）。

式中，Ω为测站的地固坐标系到发射坐标系的转换矩阵，如式（2）所示。

式中，AT为导弹设计瞄准方向的天文射击方位角；λa、φa为发射坐标系原点的天文经度和纬度；λi、φi为测站i的天文经度和纬度。

步骤2：计算两测站到目标间夹角φ12的余弦，见式（3）。

步骤3：计算目标到测站连线与两测站间连线夹角的余弦，如式（4）所示。

式中，l12＝x02－x01，m12＝y02－y01，n12＝z02－z01，D12＝［（l12）2＋（m12）2＋（n12）2］1／2，x0i，y0i，z0i（i＝1，2）为第i个测站在地固坐标系下的站址坐标。

步骤4：计算目标到测站的斜距，如式（5）所示。

式中，sinφi＝（1－cos2φi）1／2，i＝1，2。

步骤5：计算目标在地固系中的位置参数，如式（6）所示。

多站交会计算目标位置，选择布站位置最好的一个测站，与其余测站配对，重复步骤1 ～5 的计算过程，得到每个组合计算的目标位置参数。

2.2 雷达计算目标位置

根据目标与测站几何关系及球面坐标关系式，可计算出目标在地固坐标系下的弹道位置参数，如式（7）所示。

式中，R，A，E为经数据预处理修正后的雷达测距、方位角和俯仰角观测量，Ω为该设备地固坐标系到发射坐标系的转换矩阵，x0，y0，z0为该观测站在地固坐标系下的站址坐标。

2.3 实时弹道动态滑窗确定目标位置速度

针对每条测元弹道，在时间序列t1，t2，…，tn上连续获取n个时刻的目标位置参数xi，yi，zi，（i＝1，…，n）后，对3 个方向的位置参数分别进行p（p≥3）阶正交多项式拟合。设拟合得到的多项式如式（8）所示。

对式（8）求导，可得式（9）。

假设输出的平滑时刻为t0，则t0时刻的目标位置速度参数为x（t0），y（t0），z（t0），x′（t0），y′（t0），z′（t0）。如果t0为窗口中心时刻，则有半个时间窗口的延迟；如果t0为窗口右端点tn时刻，则可以实现无滞后的弹道估计；如果t0为tn＋1时刻，则可以将多项式拟合值作为弹道预测值，与tn＋1时刻实测值比较，设置初始阈值作为野值剔除判断的一个标准。

在实时弹道计算中，通常对时间采样序列使用动态滑窗进行多项式拟合。当获取到tn＋1时刻目标位置参数后，对时间序列t2，…，tn＋1重复上述过程，可获得新的平滑时刻的目标位置速度参数。分析表明，选用3 阶多项式即可达到很好的逼近效果，在实际计算过程中可适度调整多项式阶数和滑动窗口宽度，以满足时间要求和精度要求。

3 实时弹道自动优选方法

为了确保对返回舱的持续、有效跟踪，在返回舱航区周围要布置各种类型的光学和雷达设备，为保证测量成功率，通常会组合形成多个典型测量方案，加之使用不同弹道计算方法，结果可能产生多达十几条测元弹道。为提高优选方法的自动性能，减少人工操作，降低操作风险，弹道选优使用自动选择方法。首先根据数据使用率计算弹道置信度，然后结合置信度进行弹道加权融合计算，每次迭代输出方差最小的弹道作为最优弹道。本文设计的弹道自动优选方法步骤如下：

1）将各测元弹道输出时间统一到同一时刻，根据测量数据使用率计算各测元弹道置信度。对于使用不同测量设备确定的弹道，由于各测量设备测量性能、观测几何均存在差异，导致计算得到的测元弹道可信度不尽一致。弹道计算中，将前一刻计算弹道外推到当前时刻，反算测量数据，并与实际接收到的测量数据比较。根据各测量设备的精度指标事先约定的最大允许偏差，剔除无效数据，计算测量数据使用率∂。如果∂≥80%，设置本条置信度ω＝1；如果∂≤40%，设置本条置信度ω＝0；其他情况ω＝∂／80%。弹道置信度信息1 为最优，0 为不可用。

2）依据先验信息，对不同方式解算的弹道赋不同权重。使用外测数据源，根据不同算法确定的弹道，由于算法计算机理不同，采用步骤1）计算各测元弹道的置信度，同样大小的置信度实际弹道精度可能会有较大差别。针对此种情况，需要根据事先仿真和实测数据验证的经验，在目标飞行正常情况下，对不同测元弹道采用不同权重。设需要融合的测元弹道数目为n，各弹道的置信度和权重分别为ωi，λi（i＝1，…，n）。结合各测元弹道的置信度对弹道权重进行归一化处理，设处理后各测元弹道的权重为λ－i，如式（10）所示。

3）对各测元弹道加权，计算平均弹道。设各测元弹道计算结果为yi（j）（i＝ 1，…，n；j＝ 1，2，3），其中j描述的是位置3 个方向的属性。由于速度信息的计算一般依赖于位置信息，计算平均弹道时可以不考虑速度影响。对各测元弹道加权求均值，设平均弹道为y＾（j），如式（11）所示。

5）根据各测元弹道距平均弹道的方差大小，动态调整各测元弹道的新权重。记各测元弹道的方差和为εa，如式（13）所示。

6）依据飞行时序，重复步骤1）～5），直至观测或弹道计算结束，全程可自动获取一条最优弹道。

实时自动弹道优选方法流程如图1 所示。

图1 实时弹道优选方法流程图Fig.1 Flow chart of real-time trajectory optimization method

本文实时弹道自动优选方法结合弹道置信度进行弹道加权融合优选计算，输出最优弹道时保留原始位置速度信息，有效避免了因某一方法弹道波动或间断引起优选弹道的波动，提高了优选弹道的可靠性和精度。

4 落点预报简化动力学模型

在地心惯性坐标系中，简化的返回舱落点预报动力学模型如式（16）、（17）所示。

在返回舱落点预报中，需要结合返回前1 h 着陆场实测气象风数据，制定落点预报气象风修正量插值表。实际落点计算时，将返回舱实时位置速度状态变量带入上述简化的动力学模型，通过积分得到后续时刻的目标位置参数（需先预估落点的高程，以控制外推结束时刻）。需要注意的是，积分过程的每一步都需要插值获取相应高程的气象风修正量，对目标位置速度进行修正［12-13］。

5 精度分析

5.1 初始条件设置

为验证实时弹道动态滑窗处理和弹道优选方法的精度，采用返回舱开伞前后一段700 s 弹道作为理论弹道，仿真测量数据并计算外测弹道。理论弹道的速度、高度曲线如图2 所示。

图2 返回舱速度和高度曲线Fig.2 Velocity and altitude curve of re-entry capsule

图2 中，返回舱速度从6 km／s 减少至0.01 km／s ，高度从53 km 下降至3 km。高度为10 km，时间300 s 左右返回舱开伞。

假设跟踪设备均为脉冲雷达，设备采样间隔为0.05 s，标称误差设置如表1 所示。

表1 雷达设备标称误差设置Table 1 Nominal error setting of radar equipment

5.2 实时弹道计算精度分析

根据以上条件，分别采用单点定位、滤波方法和实时弹道动态滑窗处理3 种方法计算弹道与理论弹道，比较分析弹道计算精度（图3～图5）。其中，实时弹道动态滑窗处理时，滑动窗口取200 个点，弹道输出时间取窗口第95%数据时刻。

图3 单点定位位置速度偏差Fig.3 Position and velocity deviation of single point positioning

图4 滤波方法位置速度偏差Fig.4 Position and velocity deviation of filtering method

图5 动态滑窗处理位置速度偏差Fig.5 Position and velocity deviation of dynamic sliding window processing

3 种弹道计算方法偏差统计结果，如表2 所示。可以看出，3 种弹道计算方法位置偏差量级一致，差别主要在速度方面。其中，单点定位解算的速度偏差波动较大，总的偏差为98.886 m／s，滤波方法速度偏差为36.570 m／s，实时弹道动态滑窗处理方法速度偏差为15.960 m／s。通过以上分析，滤波方法和动态滑窗处理方法能够大幅度提高单点定位的精度，动态滑窗方法精度最高。

表2 偏差统计分析结果Table 2 Statistical analysis of the deviation

5.3 弹道优选方法效果分析

仿真分析时，实时模拟4 条测元弹道。采用本文实时弹道自动优选方法，分析获取到的优选弹道是否为最优弹道，如图6 所示。从图6 可以看出，优选弹道有如下选取原则：

图6 实时弹道优选方法位置比较Fig.6 Position comparison of real-time trajectory optimization methods

1）只有1 条测元弹道时，输出该弹道为优选弹道，如图6 中第650 ～700 s 优选弹道输出测元弹道4；

2）只有2 条测元弹道时，输出2 条弹道的均值，如图6 中第560 ～650 s 优选弹道输出测元弹道3 和弹道4 的平均值；

3）有3 条或3 条以上测元弹道时，每一时刻均输出了方差最小的弹道。

图6 中除只有2 条测元弹道取均值外，其他时间优选弹道输出时不改变测元弹道的计算结果，本文方法最大限度地保留了目标的真实飞行情况，确保了优选弹道的可靠性和精度需求。

需要注意的是，方差最小弹道不一定是该时刻偏差最小的弹道，图6 中130 ～150 s 优选弹道选择了红色测元弹道1，而不是偏差最小的测元弹道2（绿色线），因为测元弹道2 明显偏离其他3 条弹道。实时弹道优选过程是一个逐步逼近最优弹道的过程，在选择过程中剔除距平均弹道方差最大的弹道，最后优选弹道一定是一组精度相当测元弹道中的某一条弹道。

5.4 落点预报精度分析

分析落点预报精度时，设定4 部非合作式雷达跟踪，1 倍设备标称测量误差，进行100 次Monte Carlo 仿真。每一时刻获取到返回舱高精度优选弹道后，采用返回舱落点预报动力学模型，结合返回前1 h 着陆场实测气象风数据，对10 km以下的着陆段，实时对落点预报结果进行气象风漂移修正，预报着陆时间和落点位置，并比较实际落点与理论落点的偏差，如图7 所示。

图7 返回舱着陆时间偏差和位置偏差Fig.7 Landing time error and position error of reentry module

图中，返回舱开始经过黑障段，目标特性发生剧烈变化，雷达难以连续平稳跟踪，此时间段落点预报偏差较大；开伞后，返回舱经过短暂的剧烈跳动，飞行趋于稳定；对实时弹道进行气象风漂移修正后，返回舱落点预报精度可达到百米量级，着陆时间偏差在1 min 以内。实际任务中，因为实时风向、风速与预报的高空风数据会有差别，加上伞降模型误差等因素，真实返回落点计算偏差与仿真会存在一些差异。