CN105974777B - 一种利用Algos和Kalman组合产生原子时标的方法 - Google Patents

Abstract

本发明公开一种利用Algos和Kalman组合产生原子时标的方法，首先对原始钟差数据进行离群值检测；计算前六个月的TA‑K哈德玛方差以确定每台钟的权重，作为计算本月TA‑K权重；其次用改进的Algos算法得到的结果，作为Kalman算法的噪声矩阵和初值的输入，利用哈德玛方差确定每台钟的权重；最后预测下一时刻的钟差值。采用本发明的技术方案，既可以降低Kalman算法的发散性，又能使时间尺度稳定性好。

Description

一种利用Algos和Kalman组合产生原子时标的方法

技术领域

本发明属于原子时标的技术领域，尤其涉及一种利用原子钟数据通过Algos和Kalman组合产生原子时标的方法。

背景技术

时间频率的研究是基础科学研究中一个重要的分支，时间频率在科研、定位系统、电力系统、军事、国家安全等方面具有举足轻重的地位。为此世界各国大都建有自己的守时实验室，产生自己国家的原子时标并参与国际比对。

中国计量院承担着向全国发送标准的时间和频率，因此精确的时间尺度具有很大的意义。在中国计量院，通过原子时算法控制钟组来产生时间和频率的输出。在原子钟组稳定的前提下，原子时算法直接影响原子时标水平。

当前计量院守时实验室采用的原子时算法是Algos算法，得到的时间尺度具有较高的长期稳定性，但是属于滞后性算法。并且该算法采用标准方差来衡量原子钟噪声的大小，这已经不能精确的描述原子钟噪声的大小。

发明内容

本发明要解决的技术问题是，提供一种利用Algos和Kalman组合产生原子时标的方法，可以既保证了原子时标的稳定性，又保证了原子时标的实时性。

为解决上述问题，本发明采用如下的技术方案：

一种利用Algos和Kalman组合产生原子时标的方法包括以下步骤：

步骤1、对实验室测得的原子钟钟差数据进行离群值检测；

步骤2、根据步骤1得到的原子钟钟差数据计算每台原子钟分配权重；

步骤3、计算实验室单个原子钟与实验室原子钟组平均时间尺度之差x_i(t)，即TA-K数据，

若每台原子钟的权重大于最大权重上限，则该台原子钟的权重为最大权重上限；设h'_i(t)为t时刻原子钟H_i时间改正量，则将权重、原子钟数据、h'_i(t)通过如下公式计算：

其中，x_ij(t)为原子钟H_j和原子钟H_i在t时刻的钟差值；x_i(t)为实验室单个原子钟与实验室原子钟组平均时间尺度之差，x_i(t)与x_j(t)不相同，N为1…i…j…N；

步骤4、根据TA-K数据x_i(t)和每台钟的权重计算钟组的原子时标；

将x_i(t)数据作为Kalman滤波器模型的输入，计算t时刻原子钟频差值的预测值为：其中，A(t)为t时刻原子钟频差值，为A(t)的临时原子钟频差值预测值，K为Kalman滤波器的增益，H为测量矩阵，

计算t时刻钟组的原子时标TA(t)：

其中，ω_i为每台原子钟的权重，A_i(t)为每台原子钟的当前原子钟频差值，为每台原子钟的预测原子钟频差值。

作为优选，所述离群值检测利用莱特准则，检测的原子钟频差数据为{x₁,x₂,x₃,…x_N-1,x_N}，样本均值为标准偏差为莱特准则是若测量值则该值是离群值，剔出所述离群值得到了离群值检测之后的原子钟频差数据。

作为优选，步骤2具体为：初始化时，根据实际情况第一个月对每台原子钟分配权重ω₁为N为钟组钟的个数；处理数据小于六个月时，计算前几个月的哈德玛方差作为下个月每台钟的权重；在处理钟差数据时间大于六个月的时候，采用前六个月钟差数据计算权重，即首先计算六个月的哈德玛方差σ²，即然后对每台原子钟分配权重

作为优选，步骤3中，最大权重上限ω_max为a为经验常数，N为原子钟的个数。

作为优选，步骤4中，计算A(t)的临时原子钟频差值预测值和Kalman滤波器的增益K公式如下：K(t)＝P^-(t)H^T(HP^-(t)H^T)^-1，其中，P^-(t)为误差协方差矩阵的预测值，为t-τ时刻原子钟频差值的预测值，H为测量矩阵，τ为数据采样间隔。

本发明采用技术方案为：首先对原始钟差数据进行预处理，包括离群值检测；计算前六个月的TA-K(K表示每台钟)哈德玛方差以确定每台钟的权重，作为计算本月TA-K权重；其次用改进的Algos算法得到的结果，作为Kalman算法的噪声矩阵和初值的输入，利用哈德玛方差确定每台钟的权重；最后预测下一时刻的钟差值。采用本发明的技术方案，既可以降低Kalman算法的发散性，又能使时间尺度稳定性好。

附图说明

图1为本发明产生原子时标的方法的流程示意图。

具体实施方式

以下结合具体实施例，并参照附图，对本发明进一步详细说明。

本发明实施例提供一种利用Algos-Kalman组合产生原子时标的方法，由于卡尔曼滤波器的发散性问题，会导致计算的时间尺度发散，而改进的Algos算法得到的时间尺度稳定性好。因此，利用改进的Algos算法计算钟的相位、频率结果及系统噪声作为卡尔曼算法的初值和系统噪声方差，防止卡尔曼算法估值结果误差随时间无穷增大；本发明Algos-Kalman组合产生的时间尺度，用于频率预报和主钟驾驭。

如图1所示，一种产生原子时标的方法包括以下步骤：

步骤1、对实验室测得的原始钟钟差数据进行离群值检测。

所述离群值检测是利用莱特准则，检测的原子钟频差数据为{x₁,x₂,x₃,…x_N-1,x_N}，样本均值是标准偏差为莱特准则是若测量值则该值为离群值，应剔出，得到了离群值检测之后的原子钟频差数据。

步骤2、根据步骤1得到的原子钟钟差数据计算每台原子钟分配权重。

初始化时，根据实际情况第一个月对每台原子钟分配权重ω₁为N为钟组钟的个数；处理数据小于六个月时，计算前几个月的哈德玛方差作为下个月每台钟的权重；在处理钟差数据时间大于六个月的时候，采用前六个月钟差数据计算权重，即首先计算六个月的哈德玛方差σ²，即然后对每台原子钟分配权重

步骤3、计算实验室单个原子钟与实验室原子钟组平均时间尺度之差x_i(t)，即TA-K数据。

利用步骤2得到的每台原子钟的权重，判断每台原子钟的权重是否大于最大权重上限，若大于，则该台原子钟的权重为最大权重上限；其中，最大权重上限ω_max为a为经验常数，N为原子钟的个数，根据经验ω_max一般设置为30％。

设h'_i(t)为t时刻原子钟H_i时间改正量，目的是为了确保在原子钟的权发生变化时或原子钟的数目变化时时间尺度的连续性，则将权重、原子钟数据、h'_i(t)通过如下方程组进行运算：

其中，x_ij(t)为原子钟H_j和原子钟H_i在t时刻的钟差值；x_i(t)为实验室单个原子钟与实验室原子钟组平均时间尺度之差，即TA-K数据(TA为地方原子时，K为实验室每台钟)；x_i(t)与x_j(t)不相同，N取值为1…i…j…N。

步骤4、根据TA-K数据x_i(t)和每台钟的权重计算钟组的原子时标。

首先，计算噪声方差矩阵C，其计算公式如下：

其中，σ_k ²为不同数据采样次数下原子钟的哈德玛方差，k＝1,2,3；τ₀为该次数据采样间隔。

Kalman滤波器的模型为A(t)＝Φ(τ)A(t-τ)+W，W为过程噪声，A(t-τ)为t-τ时刻原子钟频差值，A(t)为t时刻原子钟频差值，τ为数据采样间隔。

设为A(t)无偏估计，则P(t)为误差协方差矩阵，

将x_i(t)数据作为Kalman滤波器模型的输入，计算t时刻原子钟频差值的预测值如下：

计算误差协方差矩阵的预测值P^-(t)：

P^-(t)＝ΦP(t-τ)Φ^T+C

计算Kalman滤波器的增益K：

K(t)＝P^-(t)H^T(HP^-(t)H^T)^-1，其中，H为测量矩阵。

计算A(t)的临时预测值

其中，为t-τ时刻原子钟频差值的预测值。

计算A(t)的预测值以及t时刻的误差协方差矩阵P(t)：

P(t)＝P^-1(t)-P^-1(t)H^TK^T(t)

计算t时刻钟组的原子时标TA(t)：

其中，ω_i为每台原子钟的权重，A(t)为当前原子钟频差值，为预测原子钟频差值。

以上实施例仅为本发明的示例性实施例，不用于限制本发明，本发明的保护范围由权利要求书限定。本领域技术人员可以在本发明的实质和保护范围内，对本发明做出各种修改或等同替换，这种修改或等同替换也应视为落在本发明的保护范围内。

Claims (5)

1.一种利用Algos和Kalman组合产生原子时标的方法，其特征在于，包括以下步骤：

步骤1、对实验室测得的原子钟钟差数据进行离群值检测；

步骤2、根据步骤1得到的原子钟钟差数据计算每台原子钟分配权重；

步骤3、计算实验室单个原子钟与实验室原子钟组平均时间尺度之差x_i(t)，即TA-K数据，

若每台原子钟的权重大于最大权重上限，则该台原子钟的权重为最大权重上限；设h'_i(t)为t时刻原子钟H_i时间改正量，则将权重、原子钟数据、h'_i(t)通过如下公式计算：

其中，x_ij(t)为原子钟H_j和原子钟H_i在t时刻的钟差值；x_i(t)为实验室单个原子钟与实验室原子钟组平均时间尺度之差，x_i(t)与x_j(t)不相同，N为1…i…j…N，N为原子钟组中原子钟的数量；

步骤4、根据TA-K数据x_i(t)和每台钟的权重计算钟组的原子时标；

将x_i(t)数据作为Kalman滤波器模型的输入，计算t时刻原子钟频差值的预测值为：其中，A(t)为t时刻原子钟频差值，为A(t)的临时原子钟频差值预测值，K为Kalman滤波器的增益，H为测量矩阵，

计算t时刻钟组的原子时标TA(t)：

其中，ω_i为每台原子钟的权重，A_i(t)为每台原子钟的当前原子钟频差值，为每台原子钟的预测原子钟频差值。

2.如权利要求1所述的利用Algos和Kalman组合产生原子时标的方法，其特征在于，步骤1中，所述离群值检测利用莱特准则，检测的原子钟频差数据为{x₁,x₂,x₃,…x_N-1,x_N}，样本均值为标准偏差为莱特准则是若测量值则该值是离群值，剔出所述离群值得到了离群值检测之后的原子钟频差数据。

3.如权利要求1所述的利用Algos和Kalman组合产生原子时标的方法，其特征在于，步骤2具体为：初始化时，根据实际情况第一个月对每台原子钟分配权重ω₁为N为钟组钟的个数；处理数据小于六个月时，计算前几个月的哈德玛方差作为下个月每台钟的权重；在处理钟差数据时间大于六个月的时候，采用前六个月钟差数据计算权重，即首先计算六个月的哈德玛方差σ²，即然后对每台原子钟分配权重

4.如权利要求1所述的利用Algos和Kalman组合产生原子时标的方法，其特征在于，步骤3中，最大权重上限ω_max为a为经验常数，N为原子钟的个数。

5.如权利要求1所述的利用Algos和Kalman组合产生原子时标的方法，其特征在于，步骤4中，计算A(t)的临时原子钟频差值预测值和Kalman滤波器的增益K公式如下：其中，P-(t)为误差协方差矩阵的预测值，为t-τ时刻原子钟频差值的预测值，H为测量矩阵，τ为数据采样间隔。