CN108387915B - 一种融合钟差预报方法 - Google Patents

Abstract

一种融合钟差预报方法，对原始钟差序列进行一次差分处理得到钟差差分序列；根据钟差差分序列中学习段数据采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行拟合段的粗预报，将钟差差分序列中拟合段序列的钟差差分值与得到的拟合段粗预报值作差得到拟合误差序列；根据拟合误差序列用灰色预报模型进行精调整预报得到精调整预报值；根据钟差差分序列采用二次多项式模型进行粗预报获得钟差差分序列的粗预报值；将钟差差分序列的粗预报值和精调整预报值相加得到融合预报序列；将融合预报序列进行差分逆运算获得最终钟差预报序列；将原始钟差序列与最终钟差预报序列作差得到预报误差。本发明方法可以提高预报精度。

Description

一种融合钟差预报方法

技术领域

本发明属于原子钟数据处理技术领域，尤其涉及一种卫星星载原子钟的钟差预报方法。

背景技术

由于卫星与地球的相对位置及其运行关系，导致卫星在部分时段无法与地面站点进行双向时间比对。在无法比对的时段中，根据已有的钟差数据进行钟差预报是保持卫星高精度时间同步的重要手段。同时钟差预报也是广义多基站高精度时间同步系统(如多基地雷达系统、组网雷达系统等)各基站在无法进行时间比对的情况下保持系统时间同步的一种有效方法，是提高系统时间同步抗风险能力的重要举措，具有重要军事意义和推广应用价值。

卫星钟差预报一般是由于卫星与地球的相对位置情况导致无法比对而采取的自主预报运行手段，预报时长一般较短，也存在由于比对设备失效、钟源故障或者受到干扰等因素导致无法比对而需要进行较长时长预报的情况，总体来说，卫星钟差预报一般为中短期预报，长时间无法比对的情况比较少见。

对于中短期钟差预报问题，王宇谱等人针对北斗钟差数据的特点，研究了线性模型、二次多项式模型、灰色模型和Kalman滤波模型的单独预报效果(王宇谱，陈正生，李伟杰，等.BDS卫星钟差短期预报性能分析.大地测量与地球动力学，2017，37(5):450-456)。王志浩等人研究了灰色模型和二次多项式模型在中短期钟差预报中的应用，得出了灰色模型预报精度高于二次多项式模型的结论(王志浩，胡永辉，侯雷，等.二次多项式模型和灰色理论模型在接收机钟差预报中的应用和比较.时间频率学报，2014，37(3):157-163)。梁月吉等人提出了一种基于钟差一次差的灰色模型预报算法，提高了传统灰色预报精度(梁月吉，任超，杨秀发，等.基于一次差的灰色模型在卫星钟差预报中的应用.天文学报，2015，56(3):P264-277)。但以上方法都只研究了单一预报模型的钟差预报。采用单一模型进行钟差预报，存在预报精度不高和较明显的钟差累积效应等缺点，例如二次多项式模型算法虽然简单，但误差累积效应明显，预报误差会随预报时间的增加而显著增大；灰色模型需要样本少，存在短期预报精度较低的问题。

发明内容

本发明的目的在于提供一种可以提高钟差预报精度，同时可以降低钟差预报中的误差累积效应的融合钟差预报方法。

为了实现上述目的，本发明采取如下的技术解决方案：

一种融合钟差预报方法，包括以下步骤：

步骤一、获取原始钟差序列；

步骤二、对原始钟差序列进行一次差分处理，得到钟差差分序列，将钟差差分序列分为学习段序列、拟合段序列和预报段序列；

步骤三、根据钟差差分序列中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行拟合段的粗预报，将钟差差分序列中拟合段序列的钟差差分值与得到的拟合段粗预报值作差，得到拟合误差序列；

步骤四、根据拟合误差序列，采用灰色预报模型进行精调整预报，得到精调整预报值：

步骤五、根据钟差差分序列的钟差差分值，采用二次多项式模型进行钟差差分序列的粗预报，获得钟差差分序列的粗预报值；

步骤六、将钟差差分序列的粗预报值和精调整预报值相加，得到融合预报序列；

步骤七、将融合预报序列进行差分逆运算，获得最终钟差预报序列；

步骤八、将原始钟差序列与最终钟差预报序列作差，得到预报误差。

优选的，还包括滑动更新步骤，滑动更新包括以下子步骤：

将步骤二中钟差差分序列的学习段序列和拟合段序列组成滑动窗序列；

根据滑动窗序列中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行滑动窗序列拟合段的粗预报，将滑动窗序列中拟合段序列的钟差差分值与得到的滑动窗序列拟合段粗预报值作差，得到滑动窗序列的拟合误差序列；

根据滑动窗序列的拟合误差序列，采用灰色预报模型进行精调整预报，得到滑动窗序列的精调整预报值；

根据滑动窗序列的钟差差分值，采用二次多项式模型进行滑动窗序列的粗预报，获得滑动窗序列的粗预报值；

将滑动窗序列的粗预报值和滑动窗序列的精调整预报值相加，得到滑动更新前融合预报序列；

用滑动更新前融合预报序列对钟差差分序列的学习段序列和拟合段序列进行更新，得到滑动更新后的钟差差分序列，滑动更新结束；

基于滑动更新后的钟差差分序列执行步骤三至步骤六，得到滑动更新后融合预报序列，将滑动更新前融合预报序列与滑动更新后融合预报序列合并得到完整融合预报序列；

基于完整融合预报序列执行步骤七至步骤八，获得最终钟差预报值及预报误差。

由以上技术方案可知，本发明在对原始钟差序列进行差分处理的基础上，利用多项式模型对钟差数据进行粗预报，采用灰色模型对多项式模型拟合误差进行精调整，将多项式模型粗预报值与灰色模型精调整预报值相加构成融合钟差预报值，此时钟差预报值为原始钟差差分序列的预报值，通过差分逆运算后获得最终的钟差预报值。本发明方法相较于单一模型的预报方法，大幅提高了预报精度。同时在优选方案中，利用经过一定预报时长的融合预报后，引入滑动窗对学习样本进行更新，利用更新后的数据进行钟差预报，可以减小中期预报中钟差累积效应对钟差预报的影响，进一步提高预报的准确性。

附图说明

图1为本发明实施例1方法的流程图；

图2为钟差差分序列分段的示意图；

图3为本发明实施例2方法的流程图；

图4为钟差差分序列滑动更新的示意图

图5为基于多项式模型的短期预报误差曲线图；

图6为基于灰色模型的短期预报误差曲线；

图7为实施例1的短期预报误差曲线图；

图8为实施例1、基于多项式模型及基于灰色模型的三种预报方法的中期预报误差对比图；

图9为实施例1和实施例2的中期预报误差对比图。

下面结合附图和各实施例对本发明进一步详细说明。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下获得的所有其他实施例，都属于本发明保护的范围。

实施例1

图1为本发明实施例1的流程图，本实施例用于对短期钟差进行预报。如图1所示，本实施例融合钟差预报方法的步骤如下：

步骤一、获取某一时间段的原始钟差序列X，X＝{x⁽⁰⁾(1),x⁽⁰⁾(2),…,x⁽⁰⁾(n)}，x⁽⁰⁾(1)、x⁽⁰⁾(2)、x⁽⁰⁾(n)分别表示不同时刻的原始钟差值；

步骤二、对原始钟差序列X进行一次差分处理，得到钟差差分序列ΔX，ΔX＝{Δx⁽⁰⁾(1),Δx⁽⁰⁾(2),…,Δx⁽⁰⁾(n-1)}，将钟差差分序列ΔX分为学习段序列、拟合段序列和预报段序列，如图2所示，n为原始钟差序列的时长，l为学习段序列的时长，m-l为拟合段序列的时长，n-m为预报段序列的时长；由于原子钟相邻钟差数值相差很小，一般为大致递增或者递减的曲线，通过基于差分处理的数据预处理步骤可以得到有效数字位减少的数据序列，从而简化数据结构的复杂性，便于后续的预报处理；钟差差分序列的分段可根据预报时长进行相应调整；

步骤三、根据钟差差分序列ΔX中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行钟差差分序列ΔX拟合段的粗预报，将钟差差分序列ΔX中拟合段序列的钟差差分值与得到的拟合段粗预报值作差，得到拟合误差序列；

根据钟差差分序列中学习段序列的钟差差分值采用二次多项式模型进行学习预报，即利用学习段序列进行学习，将学习段序列的钟差差分值带入二次多项式模型进行学习，获得预报参数的估计值，得到预报参数的估计值后再次采用二次多项式模型进行拟合段的粗预报；

二次多项式模型为y(t)＝α₁+α₂t+α₃t²，其中，α₁、α₂和α₃表示预报参数，t为预报时刻，利用学习段序列进行学习时，t＝1,…,l；将二次多项式模型写成Y＝Zα，采用最小二乘法计算预报参数α₁、α₂、α₃的估计值，得到

其中，

其中，

分别为α₁、α₂和α₃的估计值，t(1)…t(l)表示各具体时刻(1,…,l时刻)，y(1)…y(l)表示学习段序列中各具体时刻(1,…,l时刻)对应的钟差差分值；

根据预报参数的估计值，采用二次多项式模型进行拟合段的粗预报，得到拟合段的粗预报值

t＝l+1,…,m；

将钟差差分序列ΔX中拟合段序列的钟差差分值与拟合段的粗预报值作差，得到拟合误差值f⁽⁰⁾(t)，

t＝l+1,…,m，Δx⁽⁰⁾(t)为t时刻对应的钟差差分值，拟合误差值构成拟合误差序列F⁽⁰⁾；

步骤四、根据步骤三得到的拟合误差序列F⁽⁰⁾，采用灰色预报模型进行精调整预报，得到精调整预报值；

采用灰色预报模型进行预报的步骤如下：

将拟合误差序列F⁽⁰⁾＝{f⁽⁰⁾(l+1),f⁽⁰⁾(l+2),…,f⁽⁰⁾(m)}进行一阶累加得到一阶累加序列F⁽¹⁾＝{f⁽¹⁾(l+1),f⁽¹⁾(l+2),…,f⁽¹⁾(m)}，其中，f⁽⁰⁾(l+1)表示t＝l+1时的拟合误差值，f⁽¹⁾(l+1)为f⁽⁰⁾(l+1)对应的一阶累加值；

根据一阶累加序列的值，采用灰色预报模型进行精调整预报，得到精调整预报值

t＝m+1,…,n，其中，

(·)^T表示矩阵转置，f⁽¹⁾(1)、f⁽¹⁾(2)、f⁽¹⁾(m-l)分别表示一阶累加序列中的第1个、第2个、第m-l个元素；

采用灰色预报模型进行预报前，先检查拟合误差序列中各元素(拟合误差值)的正负号是否一致，若一致则直接计算精调整预报值，若不一致，则给每个元素加一个常数，该常数与拟合误差序列中绝对值最大的元素的符号一致，且该常数的绝对值大于拟合误差序列中最大元素绝对值，通过灰色预报模型获得预报值后减去该常数获得精调整预报值；

步骤五、在步骤三进行二次多项式模型学习的基础上，根据钟差差分序列ΔX的钟差差分值，采用二次多项式模型进行钟差差分序列的粗预报，获得钟差差分序列的粗预报值

t＝m+1,…,n；

步骤六、将步骤五得到的钟差差分序列的粗预报值和步骤四得到的精调整预报值相加，得到融合预报序列

t＝m+1,…,n；

步骤七、将由步骤六获得的融合预报序列进行差分逆运算，获得最终钟差预报序列

x⁽⁰⁾(m)表示m时刻的原始钟差值；

步骤八、将原始钟差序列与最终钟差预报序列作差，得到预报误差error(t)＝x⁽⁰⁾(t)-x′(t)，t＝m+1,…,n。

实施例2

实施例1的方法用于对短期钟差进行预报，当用于中期钟差预报时，为了减小钟差预报误差累积效应的影响，本实施例利用融合预报值对钟差差分序列进行滑动更新后，再基于滑动更新后的数据进行预报。本实施例方法与实施例1的方法不同的地方在于：增加了一次滑动更新的步骤，滑动更新的六个子步骤与实施例1的步骤一至六的执行过程相同，只是将钟差差分序列的学习段和拟合段作为滑动窗序列来计算融合预报序列，然后再用由步骤六获得的融合预报序列替代原钟差差分序列中的学习段序列和拟合段序列，得到新的学习样本后，再基于新的学习样本进行后续的学习、粗预报、精调整预报及融合等步骤，最后将滑动更新前得到的融合预报序列和滑动更新后得到的融合预报序列合并，得到完整融合预报序列，基于完整融合预报序列计算最终钟差预报值和预报误差。

下面结合图3和图4对本实施例方法进行详细说明，本实施例以一周内的钟差数据为原始钟差序列，预报第三至第七天的钟差。将原始钟差序列一次差分处理后得到差分钟差序列，差分钟差序列分为学习段序列(第一天)、拟合段序列(第二天)和预报段序列(第三至第七天)。第一天和第二天的钟差差分数据组成滑动窗序列，用滑动窗序列进行融合预报，获得第三天和第四天的融合预报数据，该数据为滑动更新前融合预报序列，用该序列更新原钟差差分序列的学习段序列和拟合段序列后，基于更新后的钟差差分数据预报第五天至第七天的数据，获得滑动更新后融合预报序列，将滑动更新前融合预报序列与滑动更新后融合预报序列合并得到完整的融合预报序列，将完整融合预报序列差分逆运算后得到最终钟差预报值及预报误差。本实施例方法的具体步骤如下：

步骤一、获取原始钟差序列X；

步骤二、对原始钟差序列X进行一次差分处理得到钟差差分序列ΔX，将钟差差分序列ΔX分为学习段序列、拟合段序列和预报段序列，钟差差分序列ΔX的学习段序列和拟合段序列组成的滑动窗序列ΔX’，即滑动窗序列ΔX’的学习段序列和拟合段序列亦为钟差差分序列ΔX的学习段序列和拟合段序列；后续步骤中基于滑动窗序列ΔX’进行学习、拟合和预报，获得滑动更新前融合预报序列；

步骤三、根据滑动窗序列ΔX’中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行滑动窗序列拟合段的粗预报，将滑动窗序列ΔX’中拟合段序列的钟差差分值与得到的滑动窗序列拟合段粗预报值作差，得到滑动窗序列的拟合误差序列F⁽⁰⁾’；

步骤四、根据步骤三得到的滑动窗序列的拟合误差序列F⁽⁰⁾’，采用灰色预报模型进行精调整预报，得到滑动窗序列的精调整预报值；

步骤五、在步骤三进行二次多项式模型学习的基础上，根据滑动窗序列ΔX’的钟差差分值，采用二次多项式模型进行滑动窗序列ΔX’的粗预报，获得滑动窗序列的粗预报值

步骤六、将步骤五得到的滑动窗序列的粗预报值和步骤四得到的滑动窗序列的精调整预报值相加，得到滑动更新前融合预报序列r(t)’；

步骤七、用滑动更新前融合预报序列替换钟差差分序列ΔX的学习段序列和拟合段序列，即学习样本(包括学习段序列和拟合段序列)的更新，用第三、第四天的融合预报结果替换原钟差差分序列中第一、第二天的钟差差分数据，得到滑动更新后的钟差差分序列ΔX”，滑动更新结束；

步骤八、根据滑动更新后的钟差差分序列ΔX”中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行滑动更新后的钟差差分序列拟合段的粗预报，将滑动更新后的钟差差分序列中拟合段序列的钟差差分值与得到的滑动更新后的钟差差分序列拟合段粗预报值作差，得到滑动更新后的拟合误差序列F⁽⁰⁾”；

步骤九、根据步骤八得到的滑动更新后的拟合误差序列F⁽⁰⁾”，采用灰色预报模型进行精调整预报，得到滑动更新后的精调整预报值；

步骤十、在步骤八进行二次多项式模型学习的基础上，根据滑动更新后的钟差差分序列ΔX”的钟差差分值，采用二次多项式模型进行更新后的钟差差分序列ΔX”的粗预报，获得滑动更新后的粗预报值

步骤十一、将步骤十得到的滑动更新后的粗预报值和步骤九得到滑动更新后的精调整预报值相加，得到滑动更新后融合预报序列r(t)”；

步骤十二、将步骤六得到的滑动更新前融合预报序列与步骤十一得到的滑动更新后融合预报序列合并，得到完整融合预报序列；

步骤十三、将由步骤十二获得的完整融合预报序列进行差分逆运算，获得最终钟差预报序列；

步骤十四、将原始钟差序列与最终钟差预报序列作差，得到预报误差。

为了验证本发明方法的预报效果，随机选用第1970GPS周(2017.10.8-2017.10.14)的GPS部分星载原子钟精密卫星的钟差数据进行预报，同时采用现有技术的基于多项式模型的预报方法和基于灰色模型的预报方法进行对比。采用均方根误差(MSE)来表示预报误差，

β(i)、

分别为校验钟差和预报钟差。在统计中对预报误差最大最小值进行统计。

每种类型原子钟随机选择两颗星钟，具体如表1所示。

表1

表1中08号卫星星钟为铯钟，其余卫星星钟为铷钟。

选用2017.10.8一天的钟差数据，进行预报时长为20session(5h)的预报，1session代表15分钟。由于预报时长较短，采用实施例1的预报方法进行预报。

图5至图7分别为多项式模型(PM)、灰色模型(GR)和实施例1的预报误差曲线图，统计性能如表2所示。从图5至图7和表2可以看出，短期预报方面，相对于多项式模型/灰色模型单独预报，在短期预报上，实施例1的方法可将预报精度提升50％以上，预报MSE均值优于0.5ns(5h)。图5至图7中的“PGx”表示星载原子钟编号。

选用一周的钟差数据，采用实施例1的预报方法进行预报时长为5天的预报，将前两天的钟差数据作为学习样本，同样采用现有技术的基于多项式模型的预报方法和基于灰色模型的预报方法进行对比。三种方法的预报误差统计对比图如图8所示，统计性能如表2所示。从图8并结合表2统计可看出，中期预报方面，实施例1的方法获得了预报MSE均值小于2ns的统计精度(5d)，预报性能优于多项式模型/灰色模型单独预报。

选用一周的钟差数据，采用实施例2的预报方法进行预报时长为5天的预报，滑动更新前用第一天和第二天的钟差差分数据为学习样本，预报第三天和第四天的融合预报数据，然后用该融合预报数据更新学习样本，预报第五天至第七天的钟差差分序列。实施例1和实施例2的预报误差统计对比图如图9所示，统计性能如表2所示。从图9并结合表2可以看出，通过增加滑动窗性能更新的步骤，降低了钟差预报中的误差累积效应，较不添加滑动窗的预报方法的预报精度提高了约31.67％。

表2

注：Avg表示均值，Std表示标准差。(单位：ns)

以上所述，仅是本发明的较佳实施例而已，并非对本发明做任何形式上的限制，虽然本发明已以较佳实施例揭露如上，然而并非用以限定本发明，任何熟悉本专业的技术人员，在不脱离本发明技术方案范围内，当可利用上述揭示的技术内容做出些许更动或修饰为等同变化的等效实施例，但凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属于本发明技术方案的范围内。

Claims (1)

1.一种融合钟差预报方法，其特征在于，包括以下步骤：

步骤一、获取原始钟差序列；

步骤二、对原始钟差序列进行一次差分处理，得到钟差差分序列，将钟差差分序列分为学习段序列、拟合段序列和预报段序列，钟差差分序列的学习段序列和拟合段序列组成滑动窗序列；

步骤三、根据滑动窗序列中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行滑动窗序列拟合段的粗预报，将滑动窗序列中拟合段序列的钟差差分值与得到的滑动窗序列拟合段粗预报值作差，得到滑动窗序列的拟合误差序列；

步骤四、根据滑动窗序列的拟合误差序列，采用灰色预报模型进行精调整预报，得到滑动窗序列的精调整预报值：

步骤五、根据滑动窗序列的钟差差分值，采用二次多项式模型进行滑动窗序列的粗预报，获得滑动窗序列的粗预报值；

步骤六、将滑动窗序列的粗预报值和滑动窗序列的精调整预报值相加，得到滑动更新前融合预报序列；

步骤七、用滑动更新前融合预报序列替换钟差差分序列的学习段序列和拟合段序列，得到滑动更新后的钟差差分序列；

步骤八、根据滑动更新后的钟差差分序列中学习段序列的钟差差分值，采用二次多项式模型进行学习预报，根据学习结果采用二次多项式模型进行滑动更新后的钟差差分序列拟合段的粗预报，将滑动更新后的钟差差分序列中拟合段序列的钟差差分值与得到的滑动更新后的钟差差分序列拟合段粗预报值作差，得到滑动更新后的拟合误差序列；

步骤九、根据滑动更新后的拟合误差序列，采用灰色预报模型进行精调整预报，得到滑动更新后的精调整预报值；

步骤十、根据滑动更新后的钟差差分序列的钟差差分值，采用二次多项式模型进行滑动更新后的钟差差分序列的粗预报，获得滑动更新后的粗预报值；

步骤十一、将滑动更新后的粗预报值和滑动更新后的精调整预报值相加，得到滑动更新后融合预报序列；

步骤十二、将滑动更新前融合预报序列与滑动更新后融合预报序列合并，得到完整融合预报序列；

步骤十三、将完整融合预报序列进行差分逆运算，获得最终钟差预报序列；

步骤十四、将原始钟差序列与最终钟差预报序列作差，得到预报误差。

Images