空间碎片激光测距的距离预报偏差实时修正算法

高健，梁智鹏，韩兴伟，2，董雪，安宁，温冠宇，刘承志，2

（1 中国科学院国家天文台长春人造卫星观测站，长春 130117）

（2 中国科学院空间目标与碎片观测重点实验室，南京 210008）

0 引言

卫星激光测距技术（Satellite Laser Ranging，SLR）作为单点精度最高的空间大地测量技术之一，其测距精度已接近毫米级，将SLR 技术应用在空间碎片的精密测定轨中能够使其轨道精度得到大幅提升［3-6］。自2002 年澳大利亚EOS 公司GREENE B 等［7］首次报告实现了空间碎片的漫反射激光测距后，国内外多个激光测距台站也相继开展了空间碎片测距实验，并成功获取到空间碎片的激光观测数据［8-11］。

空间碎片测距相比合作目标测距的主要难点在于极低的漫反射回波接收信噪比和较大的初轨预报偏差。近年来针对空间碎片测距的方法创新和技术改进主要集中于提高测距系统的微弱回波探测能力，如近红外波长测距［12］、多脉冲测距［13-14］、单站发射多站接收［15-17］、超导纳米线单光子探测［18-19］等多种技术手段的组合与叠加。然而测距成功率由探测概率和虚警概率共同决定，为了抑制天空背景噪声，使虚警概率保持在较低水平，普遍采用时间滤波的距离门控技术，其依赖于高精度的距离预报值来计算精确的探测器开门时刻，而较大的距离预报偏差将加大回波捕获难度，降低回波搜索效率，甚至无法在目标过境弧段内完成测距过程。

在提升空间碎片的初轨预报精度以及预报偏差的实时改进方面，通过利用已产生的短弧长稀疏激光数据或融合光学测角数据能够提升有限数量目标的测站预报精度［20-21］，但“后处理”的改进方式以预先获取激光数据为前提，无法满足激光测距对于高精度预报的高实时性需求。通过解算并应用沿迹方向的时间偏差能够在一定程度上间接减小距离偏差［22］，奥地利Graz 观测站在空间碎片白天测距中即应用时间偏差实时解算技术来减小距离搜索范围［23］，但仍无法对视线方向的距离预报偏差直接进行修正。

本文借鉴卫星编队飞行［24］的概念，将空间碎片的预报轨道与真实轨道看作邻近轨道，从空间目标轨道相对运动的角度来解决TLE 预报偏差的改进问题。以光学测角信息作为观测量，通过预报轨道上虚拟预报目标与真实目标的视位置角度偏差来反演距离预报偏差的演化过程。以卡尔曼滤波为状态估计手段，构建了近圆轨道目标的距离预报偏差实时修正算法，并通过数值仿真对算法的有效性进行了验证。

1 理论模型建立

空间碎片预报目标与真实目标分别沿各自轨道绕地心运行，其相互间的轨道运动状态服从相对运动动力学规律。空间碎片相对于观测站的距离预报偏差作为预报偏差的视线方向分量，无法直接通过光学测角手段进行测量和修正，而处在邻近轨道的空间碎片预报目标与真实目标的视位置角度偏差则可通过光学测角手段实时获取，且相对运动动力学方程为预报偏差的角度分量与距离分量之间的相互转化提供了数学基础。因此，以相对运动动力学方程为理论建模依据，通过建立相应的空间几何模型、状态演化模型与角度观测模型，即可实现距离偏差真值的求解。

1.1 空间几何模型

将空间碎片预报目标作为主星，真实目标作为从星，构建图1 所示的空间几何模型，其包含了主星、从星、地心、观测站之间的空间位置几何关系。图1 中，O点为预报目标质心，P点为真实目标质心，E点为地心，S点为观测站，Oxyz为固连在主星O上的星基轨道坐标系，x轴和y轴在轨道平面内，x轴方向由地心E指向主星O，y轴垂直于x轴并指向主星O的运动方向，z轴指向轨道平面法向，与x轴和y轴构成右手直角坐标系。γ为由观测站S指向主星O的单位矢量，同时定义γ1为单位矢量γ在轨道平面内的投影与x轴的夹角，γ2为单位矢量γ与z轴的夹角。

(1)庐山交通发展快速。2017年7月，庐山索道正式开通，并在开通后的首月共接待游客超过20万人次，受到了广大游客的喜爱。庐山索道项目也是九江市旅游的重点工程，仅8分钟就可以上下庐山，它的建造得到了江西省委、省政府、市委、市政府和庐山管理局的高度重视和支持，有效缓解了庐山旅游旺季的交通压力，为游客上庐山提供了极大的便利。

1.2 状态演化模型

在相对运动动力学方程中，Clohessy-Wiltshire（C-W）方程［25］是描述近圆轨道目标相对运动的简洁形式。考虑距离偏差修正过程仅针对空间碎片的单次短时过境，其摄动项影响远小于测量误差，因此，将C-W方程作为求解预报偏差矢量演化过程的数学模型，其解析表达式如式（1），其中(x，y，z)为真实目标在预报目标质心坐标系下的相对位置空间坐标，对应于预报偏差矢量的空间状态，n为其绕地心运行的平角速度。

由式（2）和式（3）的C-W 方程解的形式可知，预报目标与真实目标在任意时刻t的相对位置和相对速度状态仅由相对运动初始状态决定，即预报偏差矢量随时间的演化过程是确定性的。只要获取多组不同时刻的相对运动状态量或其线性组合，即可对未知量相对运动初始状态进行求解，进而得到任意时刻t的相对运动状态

1.3 角度观测模型

预报偏差矢量在空间中的演化过程对地面观测站而言表现为预报目标与真实目标在光学测角平面上投影点的相对运动。为了量化这种映射关系，构建如图2 所示的角度观测模型。其中，图2 所示平面为光学测角平面，O'与P'分别为空间几何模型中预报目标O与真实目标P在光学测角平面上的投影点，v是预报目标O的空间飞行速度方向在光学测角平面上的投影。定义α轴方向为与v平行的沿迹方向，β轴方向为垂直于v的垂迹方向，若将角度偏差观测量按所定义的沿迹方向和垂迹方向进行分解和量化，则P'在以O'为原点的O'αβ直角坐标平面内的坐标(α，β)即为角度观测模型的角度偏差观测量。

联合所构建的空间几何模型与角度观测模型，可得到角度偏差观测量(α，β)与Oxyz坐标系下相对位置状态量(x，y，z)的几何映射关系为

式中，c1=cosγ1，s1=sinγ1，c2=cosγ2，s2=sinγ2，ρ为距离预报值，对应于预报矢量的模长。同时可得到距离预报偏差Δρ与相对位置状态量(x，y，z)的几何映射关系为

式（4）和式（5）所描述的是某一时刻的静态几何映射关系，而其中c1，s1，c2，s2，ρ均隐含为时间t的函数，与角度偏差的观测时刻相对应，即几何映射关系随时间动态更新。联合式（4）的几何映射关系与状态演化模型可知，只要获取多组不同时刻的角度偏差观测量，即可映射得到多组相对运动初始状态的线性组合，进而实现对相对运动初始状态的求解。进一步地，只要得到相对运动初始状态的近似解或最优估计值，即可通过式（2）得到任意时刻t的相对位置状态量(xt，yt，zt)，并通过式（5）的几何映射关系最终得到各时刻所对应的距离预报偏差值。

2 算法设计

距离预报偏差实时修正算法设计的关键是通过角度偏差观测量的累积使距离预报偏差估计值迅速向真值收敛。从状态演化模型中C-W 方程解的形式可知，其描述的是由相对运动初始状态所定义的线性动态系统，可使用经典的卡尔曼滤波算法对其进行状态估计。定义X=为卡尔曼滤波器的状态向量，Z=(ρα ρβ)T为观测向量，其无控状态方程和观测方程分别为

式中，wk-1为过程噪声向量，vk为观测噪声向量，wk-1和vk均为白噪声向量且相互独立。式（6）的状态方程中Xk和Xk-1均对应于相对运动初始状态，因此状态转移矩阵Ak取单位矩阵E。定义与式（4）对应的几何映射矩阵H为

与式（2）和式（3）对应的相对运动状态转移矩阵M为

则式（7）的观测方程中观测矩阵Ck可表示为

即通过几何映射矩阵H和相对运动状态转移矩阵M，角度偏差观测量与相对运动初始状态直接建立了映射关系，其中矩阵H隐含有时间变量t，且H和M中的时刻t与角度偏差量的观测时刻t=tk相对应。

卡尔曼滤波对状态量的迭代求解包括预测和更新两个步骤。其预测过程如式（11），包括状态向量的先验估计与协方差矩阵计算。算法设定状态向量估计值的初值为零向量，同时设定状态协方差矩阵Pk的初值P0=diag((Δx)2(Δy)2(Δz)2，其中相对位置标准差Δx，Δy，Δz均取100 m，相对速度标准差均取0.1 m/s，过程噪声协方差矩阵Qk取零矩阵。

卡尔曼滤波的更新过程如式（12），包括卡尔曼增益Kk的计算以及状态向量估计值与协方差矩阵Pk的更新过程。算法设定观测噪声协方差矩阵Rk=diag((ρΔα)2(ρΔβ)2)，其中Δα，Δβ均与角度观测量标准差σA相等，对应于望远镜系统的角度测量精度指标，设置为5″以内，ρ对应于t=tk时刻的距离预报值。

式（11）和式（12）的卡尔曼滤波算法作为最优状态求解器实时迭代产生相对运动初始状态估计值，在每一次迭代计算的同时，将产生的代入式（2）可得到t=tk时刻的相对位置状态量估计值而后更新式（5）的几何映射关系并将代入，即得到了本次迭代对应的距离预报偏差估计值。

在算法的实际应用中，卡尔曼滤波器的角度偏差观测量输入由空间碎片真实位置观测值和预报位置之间的角度偏差依所建立的角度观测模型实时转换得到，其中真实位置通过方位俯仰或赤经赤纬形式由观测系统实时提供，与望远镜的指向方向无关，而各时刻的预报位置在观测开始之前即已计算得到，且作为已知量在观测过程中不再改变。即算法运行的距离预报偏差修正过程与目标跟踪过程中的望远镜指向偏差修正动作相互独立，互不影响。

在观测弧段的起始，系统开始接收光学测角数据并实时转换为卡尔曼滤波器的角度偏差观测量输入，算法迭代产生的距离预报偏差估计值被实时传输至测距系统。测距系统依据接收到的距离预报偏差修正量实时对TLE 距离预报进行修正并搜索回波，直至弧段结束，由此即实现了距离预报偏差实时修正技术辅助的空间碎片测距过程。

3 数值仿真与结果分析

为了验证距离预报偏差实时修正算法的有效性，选取轨道高度在1 500 km 以下的9 颗ILRS 近地联测激光卫星作为数值仿真对象，基本涵盖了近地目标不同的轨道高度和轨道倾角。所选目标同时具有TLE 两行根数预报和ILRS 提供的高精度CPF 预报，其中无修正TLE 预报轨道与CPF 预报轨道分别对应算法模型中的预报目标轨道和真实目标轨道。本文研究的是距离预报偏差改进问题，将CPF 预报轨道作为参考轨道，其轨道精度相比TLE 预报能够满足算法验证的轨道精度需求［2］。

为了模拟真实测距场景，在数值仿真中将CPF 角度预报作为卡尔曼滤波器的角度观测量输入，同时将CPF 距离预报作为真实距离观测值参考。定义CPF 距离预报与无修正TLE 距离预报之间的距离偏差为参考距离偏差ΔρREF，则算法实时产生的距离偏差估计值向参考距离偏差ΔρREF的收敛速度以及与其接近程度将作为数值仿真结果的评价标准。

获取9 颗近地目标于2021-12-01～2021-12-20 日间在长春站过境的所有弧段，在筛除最大俯仰角小于10°以及初始距离偏差小于100 m 的弧段后共计633 圈，各弧段取10°作为观测的起止俯仰角并做预报。从实用化的角度，将TLE 距离预报值与真实距离观测值之间的偏差修正至小于100 m 即可满足空间碎片观测时段常规探测器的距离门宽范围，从而有效增加回波捕获概率。因此，数值仿真以从弧段起始至|ΔρREF|＜100 所用时长作为修正时间，将修正时间与修正时间的全弧段时长占比作为评价指标。

仿真并统计各目标20 天内所有弧段的平均弧段时长、平均初始距离偏差、平均最大距离偏差、平均修正时间与平均修正时间弧段时长占比，各目标仿真统计结果如表1 所示，其中弧段时长和修正时间单位为分钟（min）。为了仿真对比角度测量精度对算法修正精度和收敛速度的影响，并验证算法对于角度测量精度指标的宽容度，将卡尔曼滤波器中角度观测量标准差σA分别设置为2″和5″。

从仿真统计结果可以看出，所选近地合作目标的平均弧段时长在16 min 以内，平均初始距离偏差量均小于500 m，且弧段最大距离偏差略大于初始距离偏差，两者相比差别不大。角度观测量标准差σA为2″时，算法能够在大部分弧段起始的1 min 以内将TLE 预报的距离偏差量修正至小于100 m，各目标平均修正时间弧段时长占比最大值小于15%，即大部分弧段能够在弧段起始的前15%完成距离偏差修正过程。对比σA为5″的仿真结果，平均修正时间较σA为2″时略有增加，表明较高的角度测量精度有利于加快距离偏差的修正速度，各目标平均修正时间弧段时长占比最大值小于20%，仍在可接受范围内，即算法对于角度测量精度指标具有较好的宽容度，能够适用于具备不同角度测量精度的望远镜系统。

以所选目标中平均初始距离偏差最大的GRACE-FO-1 为例，图3 显示了GRACE-FO-1 于协调世界时（Universal Time Coordinated，UTC）2021 年12 月18 日01：29：17～01：36：21 在长春站过境弧段的完整距离预报偏差修正过程（σA=2″）。

图3 中弧段时长7.1 min，初始距离偏差接近1 000 m。在弧段的起始，距离偏差估计值的初值为0，随着光学角度偏差观测量的累积，修正算法使得距离偏差估计值迅速向参考距离偏差ΔρREF收敛；在弧段开始1 min 后，ΔρREF与的偏差量减小至100 m 以内，也即TLE 距离预报值与真实距离观测值之间的偏差量被修正至100 m 以内。此时已满足了常规探测器的距离门宽搜索范围要求，激光测距系统可稳定地获取有效回波，随后距离偏差量进一步减小，一直稳定在ΔρREF附近，直至弧段结束。

本文算法主要在于通过减小距离预报偏差来缩减距离门回波搜索范围，提高有效回波信号落入距离门宽内的概率，进而提升距离门控回波搜索捕获效率和整体测距成功概率。其对于测距成功概率的提升程度则主要与不同观测时段不同天空背景噪声强度下所允许的探测器最大距离门宽，以及距离预报偏差修正前后的距离偏差量有关。空间碎片测距能力的提升是多种技术手段综合运用的结果，在实际应用中，将本文提出的距离预报偏差实时修正技术与其它测距技术手段相结合，能够发挥出更大的作用。

4 结论

非合作目标空间碎片的距离预报偏差修正问题一直是领域内的难题，将视位置角度偏差作为观测量的距离预报偏差实时修正方法为其提供了新的思路和解决方案，仿真统计结果表明方法具有很高的实用价值，而能够持续稳定地获取角度偏差观测量是方法成功的关键。当前空间碎片白天测距仍是领域内的难题，受限于白天强天空背景噪声下的探测器最大距离门宽限制，其对于距离预报精度的要求更加严苛，若能够很好地解决白天空间目标的可视问题，方法将在空间碎片白天测距中也发挥重要的作用，有利于空间碎片白天测距能力的提升。算法模型中指向偏差、蒙气差、角度采样延时、大气延迟误差等观测量和状态量噪声对于算法稳定性的影响将通过进一步的数值仿真和实测实验进行研究和验证。