高精度晶体振荡器的守时补偿方法
【专利摘要】本发明提供一种高精度晶体振荡器守时补偿方法，包括：通过FPGA采集板卡对晶体振荡器的频率进行处理,得到晶体振荡器与GPS秒脉冲之间的对时误差,测量晶体振荡器的工作温度，构造温度特性曲线，得到晶体振荡器温度特性模型,对温度补偿之后的对时误差数据进行卡尔曼滤波，得到构造的基准时钟，从而获得GPS误差特性曲线,在分离了温度误差及GPS误差之后，采用最小二乘法预测守时补偿值。本发明根据GPS接收机的误差特性、晶体振荡器的相位随机误差特性、温度频率偏移误差特性以及老化特性等，对晶振守时误差进行晶振温度?频率变化特性模型的重构，并使用卡尔曼滤波的方法构建基准时钟，以此分离各误差源的特性，得到高精度的守时补偿值。
【专利说明】
高精度晶体振荡器的守时补偿方法
技术领域
[0001] 本发明涉及电力系统继电保护及同步采样技术领域，具体而言涉及基于实时温度 测量更新晶体振荡器温度模型参数，以此构造基准时钟模型并生成守时补偿值的方法。
【背景技术】
[0002] 基于GPS的晶振驯服技术是利用GPS接收机输出的秒脉冲信号锁定本地时钟，当 GPS模块失锁或者秒脉冲信号丢失之后，则主要依靠晶振自身的稳定度及预设的补偿量来 保证晶振在守时阶段的准确度。因此，当晶振自身的稳定度一定时，对于守时阶段的补偿算 法的研究具有重要意义。
[0003] 目前对于温度补偿晶体振荡器的守时补偿算法主要是针对算法本身的改进，如采 用改进kalman滤波或者滑动平均等算法来实现守时。由于kalman滤波对于已知误差模型的 预测较为有效，滑动平均算法则主要应用于具有正态分布特性的误差中，而守时误差中包 含各类误差特性，若将其统一处理，将会导致由于模型不准确而造成的预测偏差。另外还有 利用历史数据采用拟合的方式预测守时补偿值。预测需基于准确的模型，若采用各误差特 性相叠加的对时误差进行分析处理，会导致由于误差特性混叠从而加大预测算法的复杂 度，同时降低预测准确度。因此迫切需要提出一种分离各个误差特性的守时方法。
【发明内容】
[0004] 本发明的目的在于解决目前守时算法将误差统一建模处理而未考虑各误差源的 分布特性导致守时精度降低的问题。
[0005] 为解决上述技术问题，本发明所采用的技术方案如下：
[0006] -种高精度晶体振荡器守时补偿方法，包括以下步骤：
[0007] 通过FPGA采集板卡对晶体振荡器的频率进行处理，得到晶体振荡器与GPS秒脉冲 之间的对时误差，采用高精度温度变送器测量晶体振荡器的工作温度，构造温度特性曲线， 得到晶体振荡器温度特性模型，对温度补偿之后的对时误差数据进行卡尔曼滤波，得到构 造的基准时钟，从而获得GPS误差特性曲线，在分离了温度误差及GPS误差之后，采用最小二 乘法预测守时补偿值。
[0008] 在一些实施例中，具体步骤为：
[0009] 第一步，以晶体振荡器的频率作为基准，对振荡器自身时钟进行计数，得到相对时 标(程序运行开始累计的时间值）。同时，判断GPS接收机的时钟是否有效，若有效进一步判 断GPS秒脉冲是否失准，得到以GPS秒脉冲为基准的晶体振荡器对时误差。
[0010] 第二步，得到晶体振荡器的初始温度频率特性曲线。
[0011]选取测量的温度数据，将其中的最大值及最小值去除（避免温度传感器的突变量 对测量造成较大影响），其余值取算术平均。同时对振荡器的频率处理方式也相同，将最大 与最小温度值对应的频率去除，其余频率值取算术平均：
[0013]其中，化为第i秒的晶体振荡器频率计数值，η为计数次数。将频率与温度的采样值 代入如下方程中进行求解：
[0015] 式中，△ fi为对应温度Ti时测得的频率fi和fQ的差值；△ Ti为实时温度Ti与To的差 值；Δ E频率修正量。通过引入频率修正量实现对三次函数的ao、bo、co三个温度特性曲线的 参数进行更新修正，得到晶体振荡器的初始温度频率特性曲线。
[0016] 第三步，得到晶体振荡器温度频率特性曲线。
[0017] 在滤除GPS接收机及振荡器频率偏差的综合误差曲线中的高频随机分量之后，采 用滤波后的数据取负值作为基准时钟信号，再次采集温度及振荡器频率，进而计算出新的 温度频率特性的各次参数。重复上述步骤，直到温度频率特性参数更新值小于设定范围，则 停止更新。最终得到的温度频率特性曲线。
[0018] 第四步，采用卡尔曼滤波得到基准时钟。
[0019] 实际测得的误差量包含了 GPS的随机抖动及晶体振荡器的各误差特性曲线，其表 达式为
[0020] XGPS = ErrGPS+X( 3)式中，XGPS表示总的误差曲线;ErrGPS为GPS接收机的随机抖动;X 为晶体振荡器初始温度特性曲线造成的频率偏移。
[0021] 采用卡尔曼滤波算法对XCPS进行预处理，。卡尔曼滤波的核心由5个基本公式组成， 主要可以分为时间更新方程和测量更新方程。其中时间更新方程负责及时向前推算当前状 态变量和误差协方差估计值，以便为下一个时间状态构造先验估计，测量更新方程将先验 估计和新的测量变量结合以构造改进的后验估计。
[0022] 时间更新方程：
[0024]状态更新方程：
[0025] Κ, =ΡΚ1!' {IIP, IIr + R):
[0026] i = Α + Α%(5)式中，喊为k时刻系统观测值，pL为k_l时刻对系统 pL =(\-Κ,Π )Ρκ
[0027] 的控制量，Α、Β为系统参数，ΡιΓ为误差协方差，Q为观测随机误差的方差;Kk为卡尔 曼增益，Η为测量系统参数，ZkSk时刻的测量值，时刻的预测值。
[0028]经过卡尔曼滤波处理可以滤除GPS接收机随机抖动中的高频分量，从而可以得到 一个较为准确的时间基准。
[0029] 第五步，构造 GPS误差特性曲线。
[0030] 对晶体振荡器计数模块进行两路并行频率计数，其中一路采用GPS为基准的频率 计数得到Errl，另一路采用之前经过卡尔曼滤波处理之后的数据为基准进行频率计数得到 Err2,其表达式为：
[0031] Erri = S-Xgps = S-Errcps-X (6)
[0032] Err2 = S-X(7)式中S为第四步构建出的较为准确的基准时钟，Em为GPS接收机有 效时的守时误差，Err2表不GPS丢失之后，利用S作为基准时钟而产生的守时误差，Xgps表不 GPS接收机与晶体振荡器的综合误差，ErrCPS表示GPS接收机的误差。
[0033] 由于在守时系统中采用两路并行频率计数的形式，Em可以直接得到，Em经过滤 波及求差之后得到的。因此，由式(6)、式⑴可以得到GPS的误差Em PS为：
[0034] Errcps = Εγγ2~Εγπ (8)
[0035] 第六步，在重构了GPS误差特性曲线的基础上，根据式
[0036] X=Xcps-ErrGPs(9)
[0037] 得到晶体振荡器总的误差特性曲线X，利用温度传感器及构建的晶体振荡器温度 特性模型，对晶体振荡器温度特性进行补偿，得到温度补偿之后的晶体振荡器误差曲线为 Xte_。当经过温度补偿之后，该晶体振荡器误差模型就可以表示为:fo'zfB+fKxoHAT+ft (υ+δ(1〇)式中，ΔΤ为考虑到温度传感器的误差;f\(X)为晶体振荡器的随机相位偏差;f3 (X)为晶体振荡器的老化偏差，S为其他系统误差;fV为进行了温度补偿之后的晶体振荡器 输出频率。
[0038] 对式(10)中的fV进行最小二乘预测。最小二乘法的表达式为：
[0039]
式中，Xpre 是采用最小二乘预测的守时数据。N为进行最小二乘算法采用的数据量。拟合得到的直线反 映晶体振荡器频率的变化，同时也为生成补偿值提供依据。
[0040] 应当理解，前述构思以及在下面更加详细地描述的额外构思的所有组合只要在这 样的构思不相互矛盾的情况下都可以被视为本公开的发明主题的一部分。另外，所要求保 护的主题的所有组合都被视为本公开的发明主题的一部分。
[0041] 结合附图从下面的描述中可以更加全面地理解本发明教导的前述和其他方面、实 施例和特征。本发明的其他附加方面例如示例性实施方式的特征和/或有益效果将在下面 的描述中显见，或通过根据本发明教导的【具体实施方式】的实践中得知。
【附图说明】
[0042] 附图不意在按比例绘制。在附图中，在各个图中示出的每个相同或近似相同的组 成部分可以用相同的标号表示。为了清晰起见，在每个图中，并非每个组成部分均被标记。 现在，将通过例子并参考附图来描述本发明的各个方面的实施例，其中：
[0043]图1是分频计数的流程图。
[0044]图2是基于GPS的晶体振荡器历史数据的示意图。
[0045] 图3是最小二乘拟合结果的示意图。
[0046] 图4是基于误差分离的改进守时算法仿真曲线的示意图。
[0047]图5是总仿真曲线的示意图。
[0048]图6是晶体振荡器守时结构图。
【具体实施方式】
[0049] 为了更了解本发明的技术内容，特举具体实施例并配合所附图式说明如下。
[0050] 在本公开中参照附图来描述本发明的各方面，附图中示出了许多说明的实施例。 本公开的实施例不必定意在包括本发明的所有方面。应当理解，上面介绍的多种构思和实 施例，以及下面更加详细地描述的那些构思和实施方式可以以很多方式中任意一种来实 施，这是因为本发明所公开的构思和实施例并不限于任何实施方式。另外，本发明公开的一 些方面可以单独使用，或者与本发明公开的其他方面的任何适当组合来使用。
[0051 ]根据本发明的实施例，提出一种高精度温度补偿晶体振荡器的守时方法。该方法 综合考虑了 GPS接收机的误差、晶振振荡器的相位随机误差、温度频率偏移误差以及老化特 性误差，利用对时误差的历史数据，通过误差反馈更新基准时钟源，得到精确的晶体振荡器 温度特性模型，进一步以更新生成的时钟源为基准得到晶体振荡器误差特性，通过最小二 乘法预测出守时补偿值。
[0052]在一些具体的实施例中，我们下面将结合附图1-附图6来更加详细地说明前述方 法的实现。
[0053]第一步，以晶体振荡器的频率作为基准，对振荡器自身时钟进行计数，得到相对时 标(程序运行开始累计的时间值）。同时，判断GPS接收机的时钟是否有效，若有效进一步判 断GPS秒脉冲是否失准，得到以GPS秒脉冲为基准的晶体振荡器对时误差。程序流程图如图1 所示。
[0054]第二步，得到晶体振荡器的初始温度频率特性曲线：
[0055]采用标称频率为25MHz的温度补偿晶体振荡器，其晶体为AT切型，温度频率偏移特 性为三次函数：
[0057]式中，fo为参考温度点的频率，可以设置零度或者室温为参考温度点；Af为实际 温度与参考温度点之间的偏差;T为实测温度;To是参考温度;AQ、BQ、C()分别为温度特性曲线 的一次、二次、三次温度系数。因此，要得到AT切型晶体振荡器的温度特性曲线只需确定其 三个温度系数即可。GPS接收机标准差为25ns。历史数据是基于GPS接收机的晶体振荡器输 出频率偏差数据，其24小时内的历史数据如图2所示(1秒为间隔进行采样）。
[0058]选取测量的温度数据，采取连续取5s的数据，将其中的最大值及最小值去除(避免 温度传感器的突变量对测量造成较大影响），其余值取算术平均。同时对振荡器的频率处理 方式也相同，将最大与最小温度值对应的频率去除，其余频率值取算术平均：
[0059]
其中，&为第1秒的晶体振荡器频率计数 值，η为计数次数。将频率与温度的采样值代入如下方程中进行求解：
[0061]式中，△ fi为对应温度Ti时测得的频率fi和fo的差值；△ Ti为实时温度Ti与To的差 值；Δ E频率修正量。通过引入频率修正量实现对三次函数的ao、bo、co三个温度特性曲线的 参数进行更新修正，得到晶体振荡器的初始温度频率特性曲线。
[0062]第三步，得到晶体振荡器温度频率特性曲线：
[0063] 采用卡尔曼滤波滤除GPS接收机及振荡器频率偏差的综合误差曲线中的高频随机 分量，采用滤波后的数据取负值作为基准时钟信号，再次采集温度及振荡器频率，计算出新 的温度频率特性的各次参数。重复上述步骤，直到温度频率特性参数更新值小于设定范围， 则停止更新，最终得到温度频率特性曲线。为了避免得到病态方程，设置温度变化在〇. 1度 以上时再进行数据更新。
[0064] 第四步，采用卡尔曼滤波得到基准时钟：
[0065] 实际测得的误差包含了 GPS的随机抖动及晶体振荡器的各误差特性曲线，其表达 式为：
[0066] XGPS = Err gps+X (4)式中，Xgps表示总的误差曲线;Errcps为GPS接收机的随机抖动;X 为晶体振荡器使用初始温度特性曲线所造成的频率偏移。
[0067]采用卡尔曼滤波算法对Xgps进行预处理。卡尔曼滤波的核心由5个基本公式组成， 主要可以分为时间更新方程和测量更新方程。其中时间更新方程负责及时向前推算当前状 态变量和误差协方差估计值，以便为下一个时间状态构造先验估计，测量更新方程将先验 估计和新的测量变量结合以构造改进的后验估计。
[0068] 时间更新方程：
[0070]状态更新方程：
[0071] Kk^ i\llj {IiPk 11' +R) ζ
[0072] 1 =XA:+[A,(Z,f-(6)式中，4为k时刻系统观测值，?/:-#k-l时刻对系统
[0073] /1=(1-^//)/-
[0074] 的控制量，A、B为系统参数，ΡιΓ为误差协方差，Q为观测随机误差的方差;Kk为卡尔 曼增益，Η为测量系统参数，ZkSk时刻的测量值，1,为1^时刻的预测值。
[0075]经过卡尔曼滤波处理滤除GPS接收机随机抖动中的高频分量，从而得到一个较为 准确的时间基准。
[0076] 第五步，构造 GPS误差特性曲线：
[0077] 对晶体振荡器计数模块进行两路并行频率计数，其中一路采用GPS为基准的频率 计数得到Errl，另一路采用之前经过卡尔曼滤波处理之后的数据为基准进行频率计数得到 Err2,其表达式为：
[0078] Erri = S-Xgps = S-Errcps-X (7)
[0079] Err2 = S-X(8)式中S为第四步构建出的较为准确的基准时钟，Em为GPS接收机有 效时的守时误差，Err2表不GPS丢失之后，利用S作为基准时钟而产生的守时误差，Xgps表不 GPS接收机与晶体振荡器的综合误差，ErrCPS表示GPS接收机的误差。
[0080] 由于在守时系统中采用两路并行频率计数的形式，Em可以直接得到，Em经过滤 波及求差之后得到的。因此，由式(7)、式⑶可以得到GPS的误差Em PS为：
[0081 ] Errcps = Err2-Erri (9)
[0082] 第六步，在重构了GPS误差特性曲线的基础上，根据式
[0083] X=Xops-Errcps (10)
[0084] 得到晶体振荡器总的误差特性曲线X，利用温度传感器及构建的晶体振荡器温度 特性模型，对晶体振荡器温度特性进行补偿。得到温度补偿之后的晶体振荡器误差曲线为 Xte_。当经过温度补偿之后，该晶体振荡器误差模型就可以表示为:fo'zfB+fKxoHAT+ft (υ+δ(11)式中，ΔΤ为考虑到温度传感器的误差;f\(X)为晶体振荡器的随机相位偏差;f3 (X)为晶体振荡器的老化偏差，S为其他系统误差;fV为进行了温度补偿之后的晶体振荡器 输出频率。
[0085] 对于直线拟合y = a*x+b
[0091] 这里E(x)为求均值。
[0092] 根据上述公式，对式(11)中的fV进行最小二乘预测。最小二乘法的表达式为：
[0094] 式中XlSti时刻对应的晶体振荡器误差值，Xpre是采用最小二乘预测的守时数据，N 为进行最小二乘算法采用的数据量t
[0095] 拟合得到的直线反映晶体振荡器频率的变化，同时也为生成补偿值提供依据。通 过多次实验反复调整式中的参数，以N=200为例对上述频率特性曲线进行最小二乘拟合， 结果如图3所示。基于误差分离的守时算法的守时曲线如图4所示。
[0096] 通过仿真来对比本发明提出的基于误差分离的改进守时算法与仅基于温度补偿 的算法以及基于预测的守时算法，得到10小时守时曲线如图5所示。
[0097] 由图5曲线(a)可见，由于传统的守时补偿方法没有考虑晶体振荡器的温度特性以 及GPS的特性，其10小时守时误差为165微秒。在图5曲线(b)中是基于温度补偿之后的晶体 振荡器守时误差曲线，由图可知，其10小时守时误差在45微秒以内。这是由于对晶体振荡器 进行了建模，考虑了温度的影响，但是由于模型本身的误差，以及忽略了GPS的随机抖动，其 守时精度仍然具有改进空间。图5曲线(c)采用本发明的守时算法，可见在基准信号GPS丢失 10小时后，守时误差呈现一个抖动的下降趋势，基本保证10小时守时误差在2微秒以内。
[0098] 虽然本发明已以较佳实施例揭露如上，然其并非用以限定本发明。本发明所属技 术领域中具有通常知识者，在不脱离本发明的精神和范围内，当可作各种的更动与润饰。因 此，本发明的保护范围当视权利要求书所界定者为准。
【主权项】
1. 一种晶体振荡器的守时补偿方法，其特征在于，包括以下步骤： 1) 采用温度变送器测量晶体振荡器的工作温度，通过FPGA采集板卡对晶体振荡器的频 率进行处理得到对时误差； 2) 根据采集的工作温度数据，构造晶体振荡器的温度特性曲线，更新晶体振荡器的温 度特性曲线参数，得到晶体振荡器温度特性模型并进行温度-频率误差补偿； 3) 将温度-频率误差补偿之后的对时误差数据进行卡尔曼滤波，得到构造的基准时钟， 从而获得GPS误差特性； 4) 在分离了温度误差及GPS误差之后，采用最小二乘法预测守时补偿值。2. 根据权利要求1所述的晶体振荡器守时补偿方法，其特征在于，在所述步骤1)中，通 过温度传感器获取晶体振荡器的工作温度，并在FPGA采集板卡设置分频电路，得到晶体振 荡器的频率与GPS秒脉冲之间的对时误差。3. 根据权利要求2所述的晶体振荡器守时补偿方法，其特征在于，在根据前述步骤1)得 到晶体振荡器的实时温度及其相应频率值后，在所述步骤2)中，按照下述方式处理进行温 度-频率误差补偿： 先设定初始温度特性参数，然后对其进行修正以获取温度特性模型；其中在选取测量 数据构造温度特性模型时，采用连续若干点的数据，将其中的最大值及最小值去除，然后将 其余值取算术平均，同时对晶体振荡器的频率处理方式也相同，将最大与最小温度值对应 的频率去除，然后将其会颇率倌取筧术平询，具体朴理如下：(1) 其中，N1为第i秒的晶体振荡器频率计数值，η为计数次数； 将频率与温度的采样值代入如下方程中进行求解：⑵ 式中，Δ fi为对应温度Ti时测得的频率fi和fo的差值;fo为基准频率，Δ Ti为实时温度Ti 与To的差值；△ E频率修正量，通过引入该修正量实现对三次函数的ao、bo、Co三个温度特性 曲线的参数进行更新修正，实现温度-频率误差补偿修正。4. 根据权利要求1所述的高精度晶体振荡器守时补偿方法，其特征在于，在所述步骤3) 中，采用卡尔曼滤波得到基准时钟，具体如下： 测得的数据为基于GPS时钟的晶体振荡器频率计数值，即测得的数据为XCPS，其中包含 了 GPS的随机抖动及晶体振荡器的各误差特性曲线： Xcps = ErrGPs+X (3) 式中，Xgps表不总的误差曲线，也即实际测得的误差量;Errcps为GPS接收机的随机抖动； X为晶体振荡器初始温度特性曲线造成的频率偏移； 采用卡尔曼滤波算法对Xcps进行预处理，卡尔曼滤波的核心由5个基本公式组成，包括 时间更新方程和测量更新方程，其中时间更新方程负责及时向前推算当前状态变量和误差 协方差估计值，以便为下一个时间状态构造先验估计，测量更新方程将先验估计和新的测 量变量结合以构造改进的后验估计； 时间更新方程：Cl) ㈤ 式中，S为k时刻系统观测值，(6/^为k-Ι时刻对系统的控制量，A、B为系统参数，/T为误 差协方差，Q为观测随机误差的方差;Kk为卡尔曼增益，H为测量系统参数，ZkSk时刻的测量 值，;^为k时刻的预测值。 经过卡尔曼滤波处理之后的数据，滤除GPS接收机随机抖动中的高频分量，从而可以得 到一个准确的时间基准。5.根据权利要求4所述的高精度晶体振荡器的守时补偿方法，其特征在于： 在重构了 GPS误差特性曲线的基础上，便可得到晶体振荡器总的误差特性曲线X: X = XGPS-Errcps (6) 在得到晶体振荡器总的误差曲线之后，基于温度传感器及构建的晶体振荡器温度特性 模型，对晶体振荡器温度特性进行补偿； 得到温度补偿之后的晶体振荡器误差曲线为Xtemp，当经过温度补偿之后，该晶体振荡器 误差模型就可以表示为： fo' =fB+fi(x〇)+ Δ T+f3(t)+5 (7) 式中，A T为考虑到温度传感器的误差，本文取温度变送器的分辨率对应的误差;A(X0) 为晶体振荡器的随机相位偏差;f3(t)为晶体振荡器的老化偏差，δ为其他系统误差;为进 行了温度补偿之后的晶体振荡器输出频率； 然后，对式（7)中的fη'讲行最小二乘预测，最小二乘法的表汰式为：(S) 式中，XAt1时刻对应的晶体振荡器误差值，Xpre3是采用最小二乘预测的守时数据，N为 进行最小二乘算法采用的数据量，其M
【文档编号】H03B5/04GK106026919SQ201610324056
【公开日】2016年10月12日
【申请日】2016年5月16日
【发明人】卜京, 贾鹏辉, 游步新, 殷明慧, 谢云云, 姚娟, 杨阳
【申请人】南京理工大学