CN105612704A - 利用锁相环用于时间和频率同步的方法和装置 - Google Patents

Abstract

本发明涉及一种用于时间和频率同步的方法和装置，尤其在利用例如IEEE？1588精确时间协议(PTP)的分组网络上的用于时间和频率同步的方法和装置。定时协议消息暴露于诸如分组延迟变化(PDV)或者分组丢失的网络中的缺陷。本发明的实施方式提供了一种基于直接数字合成的数字锁相环(DPLL)，以提供在从装置(时间客户)使用的时间信号和频率信号二者。还提供了结合针对时钟偏移和偏斜估计的递归最小二乘方机制的该DPLL的示例。

Description

利用锁相环用于时间和频率同步的方法和装置

技术领域

本发明涉及用于时间和频率同步的方法和装置。具体但非排他地，涉及使用锁相环在使用例如IEEE1588精确时间协议(PTP)的分组网络上的时间和频率同步。

背景技术

移动网络分为两类：频分双工(FDD)，其使用两个分离的频带来发送/接收；以及时分双工(TDD)，其在单个频带上发送和接收。LTE-TDD、WiMAXTDD、CDMA网络(在北美普遍)、TD-CDMA和TD-SCDMA(除频率同步以外还)需要时间同步，而LTE-FDD、GSM(全球移动通信***)、W-CDMA和其它无线技术仅需要频率同步。即使随着LTE-FDD的使用，诸如网络MIMO和基于位置的服务这样的新的LTE移动服务也将需要精确时间同步。

本发明特别应用于分组网络上的时钟同步(时间和频率二者)，尤其是电信网络同步。与需要毫秒级准确性以良好操作的IT计算***和传感器网络不同，电信网络要求亚微秒(实际上，纳秒)级准确性。这种严格的时钟质量等级常规上由GPS、原子时钟和TDM定时链路来提供。因此，ITU-T、IEEE和其它标准体限定了特殊标准，以允许分组网络支持电信网络的特殊同步需求。现在被广泛接受且被电信业采用的一个这种最近的标准是IEEE1588精确定时协议(PTP)。甚至存在针对电信应用限定的特殊IEEE1588资料。即使用于IT计算***和传感器网络的时间同步具有相同的基本概念，这些***也具有不同的应用需求、协议、架构和实现目标，因此完全不将它们考虑在电信同步的范围内。

IEEE1588是现今工业接受的基于分组的方法/标准，其用于从主装置分配定时信息以使得分布式***的时钟能够以高精准(纳秒级准确性)同步。其基本原理是基于图1所示的同步消息的常规交换的主/从概念，从而通过在分组网络2上交换消息，将从装置3中的从时钟5同步到主装置1中的主时钟4。

IEEE1588通过将网络中的所有时钟调整到最高质量时钟(最优主时钟(grandmasterclock))来同步这些时钟。与仅支持频率转移的称作同步以太网的另一基于分组的技术不同，IEEE1588支持频率和时间转移二者转移。IEEE1588除时钟恢复机制以外还限定了宽范围的同步能力(伺服算法、PLL、定时器等)以在从(客户)装置使用以其本地时钟同步到主装置。

在需要以端对端方式实现时钟转移而没有分组网络辅助(例如，以不断跳转边界时钟(BC)[1][2]或者透明时钟(TC)的形式)的情况下，在分组网络的两端不存在可追踪到主要基准时钟(PRC)的基准时钟，或者在不存在网络GPS服务时，依赖定时的接收终端节点需要使用适应性定时技术以重建发送定时基准源的定时信号。接收终端将通常使用“基于分组的”时钟恢复机制，该时钟恢复机制使接收方时钟依从于发送方时钟。该时钟恢复机制能够处理在分组数据流中编码的所发送的时钟采样(时间戳)，以产生用于接收方的定时信号。该时钟恢复机制的目的是估计和补偿在发送方时钟和接收方时钟的振荡器之间发生的频率漂移。然而，分组延迟变化(PDV)和分组丢失的存在影响时钟估计/补偿处理的性能，使得发送方时钟显现得比其实际快或慢，最终，导致大部分为漫游(wander)的传播到接收方时钟信号。漫游是小于10Hz的时钟噪声，而抖动(jitter)是等于或大于10Hz的时钟噪声。

IEEE1588PTP基础

可使用时钟(物理或者虚拟)将事件与时间关联。单个事件的时间被称作时间戳，并是实数。在IEEE1588v2PTP中，消息分为事件消息和普通消息。所有IEEE1588PTP消息具有共同的报头。事件消息是时控消息，原因在于在各个消息的发送和接收时均产生准确时间戳。事件消息必须准确地加上时间戳，这是因为发送时间戳和接收时间戳的准确性直接影响时钟分布准确性。

普通消息不需要加上时间戳。在发送和接收任何事件消息时产生时间戳事件。事件消息的集合由Sync、Delay_Req、Pdelay_Req和Pdelay_Resp组成。普通消息的集合由Announce、Follow_Up、Delay_Resp、Pdelay_Resp_Follow_Up、Management和Singnaling组成。IEEE1588PTP允许两种不同的加时间戳方法：一步或者两步。一步时钟更新正在传输中的事件消息(Sync和Delay-Req)中的时间信息，而两步时钟在普通消息(Follow_Up和Delay-Resp)中传送分组的精确时间戳。

Sync、Delay_Req、Follow_Up和Delay_Resp消息用于产生和传递利用延迟请求-响应机制来同步普通和边界时钟所需的定时信息。Sync消息由主装置将发送至其从装置，并且包含其发送的确切时间或者随后跟随着包含该时间的Follow_Up消息。在两步普通或边界时钟中，Follow_Up消息传递针对特定Sync消息的离开时间戳的值。Delay_Req消息是用于接收节点利用Delay_Resp消息返回接收到Delay_Req消息的时间的请求。

在图1中例示了针对两步时钟的同步消息交换的基本模式。针对两步时钟的消息交换模式可解释如下。主装置1将Sync消息发送至从装置3，并且标记消息发送时间T1。从装置3接收Sync消息并且根据其本地时钟5标记接收时间T2。主装置1通过以下两种方式之一将时间戳T1传送至从装置：1)在Sync消息(未示出)中嵌入时间戳T1。这需要针对最高准确性和精确性的一些硬件处理(即，硬件加时间戳)。2)如图1所示在Follow_Up消息中嵌入时间戳T1。接着，从装置3将Delay_Req消息发送至主装置1并且标记消息发送时间T3。主装置1接收Delay_Req消息并且标记接收时间T4。主装置1通过将时间戳T4嵌入Delay_Resp消息中来将时间戳T4传送至从装置3。

在该PTP消息交换的最后，从装置3具有全部四个时间戳{T1、T2、T3、T4}。这些时间戳可用于计算从装置的时钟5相对于主时钟4的偏移和在两个时钟之间消息的平均传播时间。偏移和传播时间的计算假设主-从传播时间和从-主传播时间是相等的，即对称通信路径。

可按照端对端方式将时间/频率从主装置1转移至从装置3，而不涉及中间网络节点6，如图2所例示。在这种情景下，从装置3唯一负责正确地恢复主时钟信号。与其它方法(例如，利用不断跳转边界时钟或者透明时钟)相比，这里的时间/频率转移更有挑战性，这是因为从装置3暴露于由中间分组网络2(图3和图4)产生的所有PDV。网络装置(开关、路由器等)中的分组的处理和缓冲带来了经过分组网络2行进的分组的时间延迟的变化，这往往导致转移的时钟信号恶化。分组网络中固有的PDV是时钟噪声的主要来源。在从装置根据PTP定时信号中的恢复的时钟包含需要被去除或衰减的时钟噪声(大部分由PDV贡献)。在从装置使用滤波处理以滤掉时钟噪声，因此产生“平滑的”时钟输出。

发明内容

在其最宽的范围，本发明的从装置包括从时钟和锁相环，该从时钟和锁相环用于处理期望该从时钟与其同步的主时钟的时间和频率的估计。

本发明的第一方面提供了一种从装置，该从装置通过网络连接至具有主时钟的主装置，其中，该从装置包括：从时钟；以及数字锁相环，该数字锁相环包括相位检测器、环路滤波器、相位累加器和计数器，其中：该从装置设置成：与该主装置交换定时消息并且记录如下时间戳：根据主时钟从主装置发送所述定时消息的时间、根据从时钟接收所述定时消息的时间、根据从时钟发送所述定时消息的时间、以及根据主时钟接收所述定时消息的时间，根据所述时间戳估计从时钟相对于主时钟的偏斜和偏移；以及基于所估计的偏斜和偏移使所述从时钟同步到主时钟，以产生主时间估计。数字锁相环对主时间估计进行如下处理：在接收对主时间的第一估计时，计数器被初始化；相位检测器被设置成在接收到主时间的随后估计时，检测计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；由环路滤波器对误差信号进行滤波，以产生滤波后的误差信号；使用滤波后的误差信号来控制相位累加器的频率；以及相位累加器的输出使计数器计数递增并且还提供同步到主时钟的频率的从时钟的时钟频率。

可使用从装置中的数字锁相环使由网络(例如，可为分组网络)引入的时钟噪声衰减至等同于应用的时钟输出需求的水平。

优选地，定时消息是根据IEEE1588PTP的消息，例如，Sync、FollowUp、Delay_Req或者Delay_Response消息。

优选地，从装置利用计数器的输出作为同步到主时钟的时间的从时钟的时钟时间。

该方面的从装置因此可针对在分组网络上的端对端时间和频率分配提供基于时间戳的时钟恢复技术。该技术不要求来自分组网络的辅助，并且仍然能够提供亚微秒级时钟准确性。数字锁相环可提供在从装置(时间客户端)使用的时间和频率信号这二者。多数技术通常仅提供一种信号类型而非二者。

计数器输出(时间信号)随后可根据需要格式化成各种标准时间信号。

该方面的从装置还可包括直接数字合成器，该直接数字合成器根据相位累加器的输出产生模拟频率信号，该直接数字合成器包括：相位累加器；振荡器；映射装置；和数模转换器。从装置还可包括低通滤波器，该低通滤波器设置成对直接数字合成器的输出进行滤波。

从装置可按照其它方式处理或调整相位累加器的输出，以提供符合从装置中的终端应用的抖动需求的信号。具有接口形式因素的各种形式的信号(方波、正弦波等)可由相位累加器溢出输出构造。这些信号可由另一(通常是模拟)锁相环调整。

在一个实施方式中，数字锁相环的计数器还用于提供针对在从装置接收定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

在该实施方式中，计数器可在从装置从主装置接收到第一定时消息时被初始化，并且在接收第一主时间估计时，计数器可复位至所述第一主时间估计。

在另选实施方式中，从装置还包括自由运行第二的计数器，其中第二计数器用于提供针对在从装置接收到定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

从装置优选地使用用于估计(从装置相对于主从装置的)时钟偏移和偏斜的递归技术来产生主时间估计。该递归估计技术优选地为指数加权递归最小二乘方方法，其在结构上类似卡尔曼滤波器但不需像传统卡尔曼滤波器中那样需要要测量和处理噪声统计量的知识。这种方法的示例在以下实施方式中描述。

该方面的从装置优选地通过执行根据下面描述的第二方面的方法的相关步骤来操作。

本方面的从装置可包括上述优选和可选特征的一些、全部的任意组合，或者不需要上述优选和可选特征。

在其最宽的范围，本发明的方法提供使用锁相环处理从装置中的主时钟的时间和频率的初始估计来使从时钟的时间和频率同步到主时钟的方法。

本发明的第二方面提供一种将从装置中的从时钟的时间和频率同步到经网络连接至从装置的主装置中的主时钟的方法，该方法包括以下步骤：在主装置与从装置之间交换定时消息和以下时间戳：根据主时钟从主装置发送所述定时消息的时间；根据从时钟接收所述定时消息的时间；根据从时钟发送所述定时消息的时间；以及根据主时钟接收所述定时消息的时间，根据所述时间戳估计从时钟相对于主时钟的偏斜和偏移；以及基于所估计的偏斜和偏移使所述从时钟同步到主时钟，以产生主时间估计；使用包括相位检测器、环路滤波器、相位累加器和计数器的数字锁相环对主时间估计进行如下处理：在接收到主时间的第一估计时，初始化计数器；在接收到主时间的随后估计时，使用相位检测器来检测计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；使用环路滤波器对误差信号进行滤波，以产生滤波后的误差信号；使用滤波后的误差信号控制相位累加器的频率；以及使用相位累加器的输出使计数器计数递增并且获得同步到主时钟的频率的从时钟的时钟频率作为相位累加器的输出。

从装置中的数字锁相环可用于使由网络(例如，可为分组网络)引入的时钟噪声衰减至等同于应用的时钟输出需求的水平。

优选地，定时消息是根据IEEE1588PTP的消息，例如，Sync、FollowUp、Delay_Req或者Delay_Response消息。

优选地，所述方法还包括以下步骤：使用计数器的输出作为同步到主时钟的时间的从时钟的时钟时间。

该方面的方法因此可针对在分组网络上的端对端时间和频率分配提供基于时间戳的时钟恢复技术。该技术不要求来自分组网络的辅助，并且仍然能够提供亚微秒级时钟准确性。数字锁相环可提供在从装置(时间客户端)使用的时间和频率信号这二者。多数技术通常仅提供一种信号类型而非二者。

计数器输出(时间信号)随后可根据需要格式化成各种标准时间信号。

该方面的方法还可包括以下步骤：使用直接数字合成器根据相位累加器的输出产生模拟频率信号，该步骤包括以下步骤：映射相位累加器的输出以产生数字波形；以及使用数模转换器将所述数字波形转换为模拟波形。该方法还可包括以下步骤：对模拟波形进行低通滤波以产生平滑波形。

该方法可按照其它方式处理或调整相位累加器的输出，以提供符合从装置中的终端应用的抖动需求的信号。具有接口形式因素的各种形式的信号(方波、正弦波等)可由相位累加器溢出输出构造。这些信号可由另一(通常是模拟)锁相环调整。

在一个实施方式中，通过数字锁相环的计数器提供针对在从装置接到收定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

在该实施方式中，所述的方法还可包括以下步骤：在从装置从主装置接收到第一定时消息时初始化计数器；以及在接收到第一主时间估计时，将计数器复位至所述第一主时间估计。

在另选实施方式中，通过自由运行的第二计数器提供针对在从装置接收定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

该方法优选地利用用于估计(从装置相对于主从装置的)时钟偏移和偏斜的递归技术产生主时间估计。递归估计技术优选地为指数加权递归最小二乘方方法，其在结构上类似卡尔曼滤波器但不需像传统卡尔曼滤波器中那样要测量和处理噪声统计量的知识。这种方法的示例在以下实施方式中描述。

本方面的方法可包括上述优选和可选特征的一些、全部的任意组合，或者不需要上述优选和可选特征。

以上方面的方法优选地通过根据本发明的第一方面的从装置或者根据下面描述的本发明的第三方面的***来实现，但这不是要求。

本发明的其它方面包括在计算机***上运行的计算机程序，其可执行以上方面的方法，包括该方面的优选和可选特征的一些、全部的任意组合，或者不需要上述优选和可选特征

在其最宽的范围，本发明的***提供经网络在主装置与从装置之间的时间和频率同步，其中从装置使用锁相环处理从装置中的主时钟的时间和频率的初始估计。

本发明的第三方面提供了一种用于网络的时间和频率同步***，该***包括：具有主时钟的主装置；具有从时钟的从装置；以及连接主装置和从装置的网络，其中：从时钟包括：从时钟；和数字锁相环，该数字锁相环包括相位检测器、环路滤波器、相位累加器和计数器，其中：从装置设置成：与主装置交换定时消息并且记录如下时间戳：根据主时钟从主装置发送所述定时消息的时间、根据从时钟接收所述定时消息的时间、根据从时钟发送所述定时消息的时间、以及根据主时钟接收所述定时消息的时间，根据所述时间戳估计从时钟相对于主时钟的偏斜和偏移；以及基于估计的偏斜和偏移使所述从时钟同步到主时钟，以产生主时间估计；数字锁相环对主时间估计进行如下处理：在接收到对主时间的第一估计时，计数器被初始化；相位检测器被设置成在接收到对主时间的随后估计时，检测计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；由环路滤波器对误差信号进行滤波，以产生滤波后的误差信号；使用滤波后的误差信号来控制相位累加器的频率；以及相位累加器的输出使计数器计数递增并且还提供同步到主时钟的频率的从时钟的时钟频率。

可使用从装置中的数字锁相环使由网络(例如，可为分组网络)引入的时钟噪声衰减至等同于应用的时钟输出需求的水平。

优选地，定时消息是根据IEEE1588PTP的消息，例如，Sync、FollowUp、Delay_Req或者Delay_Response消息。

优选地，从装置利用计数器的输出作为同步到主时钟的时间的从时钟的时钟时间。

该方面的***因此可针对在分组网络上的端对端时间和频率分配提供基于时间戳的时钟恢复技术。该技术不要求来自分组网络的辅助，并且仍然能够提供亚微秒级时钟准确性。数字锁相环可提供在从装置(时间客户端)使用的时间和频率信号这二者。多数技术通常仅提供一种信号类型而非二者。

计数器输出(时间信号)随后可根据需要格式化成各种标准时间信号。

从装置还可包括直接数字合成器，该直接数字合成器从相位累加器的输出产生模拟频率信号，该直接数字合成器包括：相位累加器；振荡器；映射装置；和数模转换器。从装置还可包括低通滤波器，该低通滤波器被设置成对直接数字合成器的输出进行滤波。

从装置可按照其它方式处理或调整相位累加器的输出，以提供符合从装置中的终端应用的抖动需求的信号。具有接口形式因素的各种形式的信号(方波、正弦波等)可由相位累加器溢出输出构造。这些信号可由另一(通常是模拟)锁相环调整。

在一个实施方式中，数字锁相环的计数器还提供针对在从装置接收定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

在该实施方式中，计数器可在从装置从主装置接收到第一定时消息时被初始化；以及在接收到第一主时间估计时，计数器可复位至所述第一主时间估计。

在另选实施方式中，从装置还包括自由运行的第二计数器，其中第二计数器提供针对在从装置接收定时消息的时间的时间戳和从从装置发送定时消息的时间的时间戳。

从装置优选地利用用于估计(从装置相对于主从装置的)时钟偏移和偏斜的递归技术产生主时间估计。递归估计技术优选地为指数加权递归最小二乘方方法，其在结构上类似卡尔曼滤波器但不需像传统卡尔曼滤波器中那样要测量和处理噪声统计量的知识。这种方法的示例在以下实施方式中描述。

该方面的***优选地通过执行根据上面描述的第二方面的方法来操作。

本方面的***可包括上述优选和可选特征的一些、全部的任意组合，或者不需要上述优选和可选特征。

附图说明

现在将参照附图以举例的方式描述本发明的实施方式，其中：

图1示出了在IEEE1588PTP下两步时钟中的消息流，并且已经描述；

图2示出了通过网络的端对端时间/频率转移，并且已经描述；

图3按照示意的形式示出了分组网络中的分组延迟变化的影响，并且已经描述；

图4示出了分组延迟变化对PTP消息流的影响，并且已经描述；

图5示出了在有偏移和没有偏移的情况下，时钟偏斜对时钟的影响的简单模型；

图6示出了在PTP从装置估计服务器时间的计算处理；

图7示出了根据本发明的实施方式的从装置的架构；

图8示出了基于相位累加器的数字锁相环(DPLL)的操作；

图9示出了图8中的相位累加器的输出；

图10例示了在利用DPLL同步之后主装置与从装置之间的同步；

图11示出了在本发明的实施方式中使用的直接数字合成器(DDS)；

图12示出了在根据本发明的实施方式的从装置中操作的具有时间和频率输出的DPLL；

图13示出了具有直接数字合成器的PLL模型；

图14示出了DPLL测量和控制间隔；

图15示出了相位检测器的特征操作；

图16示出了用于为根据本发明的实施方式的时间客户端的操作建模的闭环控制模型；

图17示出了图16的闭环控制模型的进一步细节；

图18是模拟锁相环的基本框图；以及

图19示出了直接数字合成器的特征及其输入映射函数。

具体实施方式

基本时钟模型和测量和处理模型的发展

分别考虑时间服务器(主)时钟S和时间客户端(从)时钟C。时钟C相对于S在时间t≥0的差异或偏移为θ(t)＝(S(t)-C(t))。t用于指示真实或者基准(理想)时间。频率是时钟以其前进的速率，并且S在时间t的频率是s′(t)o偏斜α定义为时钟与另一时钟之间的归一化频率差。随后可针对同步问题定义一般化时钟偏移和偏斜方程。假设在任何特定时刻，具有时间线S(t)的主(服务器)时钟和具有时间线C(t)的从(客户)时钟之间的关系的瞬间示图可由图5所示的已知简单偏斜时钟模型来描述，并且通过以下方程来描述，

S(t)＝(1+α)C(t)+θ，(1)

其中，归一化的α是百万分之几数量级的非常小的量。该快照是两个时钟对齐(错位)的程度的瞬间示图。图5例示出了θ和α对时钟对齐的影响。

测量(观察)模型的发展

以上方程可扩展至将主时钟与从时钟通过有延迟的通信链路交换消息的情况考虑在内。假设从主装置行进至从装置的Sync消息经历了固定延迟d加上可变(随机)延迟ε。相似地，假设从主装置发送至从装置的Delay_Req消息经历了固定延迟d和可变延迟γ。还假设两个方向上的固定延迟分量相等(对称通信路径)，但是消息经历诸如排队延迟的可变延迟。主装置和从装置利用图1中描述的Delay_RequestDelay_Response机制来交换消息。认为***操作使得事件在离散时刻n＝0，1，2，3，...，载入记录，各个离散步长等于Δt秒的小间隔(采样)。

对于在时间戳T_1，n∈S(t)离开并且经历固定延迟d和可变延迟ε_n之后并且在时间戳T_2，n∈C(t)到达的第n个Sync消息(图1)，以上简单偏斜时钟模型可扩展以讲行进时间d+ε_n考虑在内从而获得以下表达式

(T_1，n+d+ε_n)＝(1+α_n)T_2，n+θ_n(2)

变量θ_n和α_n是在第n个Sync消息交换过程中的偏移和偏斜。

对于在时间戳T_3，n∈C(t)离开并且经历固定延迟d和可变延迟γ_n之后在时间戳T_4，n∈S(t)到达的第n个Delay_Req消息，获得以下表达式

T_4，n-d-γ_n＝(1+α_n)T_3，n+θ_n(3)

将(2)和(3)相加，得到

T_1，n+T_4，n+ε_n-γ_n＝(1+α_n)(T_2，n+T_3，n)+2θ_n

(T_1，n-T_2，n)+(T_4，n-T_3，n)＝2θ_n+α_n(T_2，n+T_3，n)+(γ_n-ε_n)(4)

以上表达式是以下形式的测量(观察)方程

$y_n={\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{X}}_n+v_n---\left(5\right)$

其中n是非负的时间索引，y_n＝(T_1，n-T_2，n)+(T_4，n-T_3，n)是标量，是矢量，是矢量，v_n＝(γ_n-ε_n)是测量噪声，并且T表示转置。诸如卡尔曼滤波器的滤波技术假设v_n是具有零平均值和方差为的白噪声。认为第n采样间隔是其中发生第n个Sync和第n个Delay_Req消息交换的时段。该测量方程将所有必要的协议(测量)变量({T_1，n，T_4，n}∈S(t)和{T_2，n，T_3，n}∈C(t))与时钟变量(θ_n和α_n)加上测量噪声v_n＝(γ_n-ε_n)关联。

可从(4)将时钟偏移表达如下

${\mathrm{\θ}}_n=-\[\frac{(T_{2,n}-T_{1,n})-(T_{4,n}-T_{3,n})}{2}\]-\frac{{\mathrm{\α}}_n(T_{2,n}+T_{3,n})}{2}-\frac{({\mathrm{\γ}}_n-{\mathrm{\ϵ}}_n)}{2}---\left(6\right)$

从该方程中可看出，对于具有零偏斜并且无噪声的***，时钟偏移简化为经典偏移方程

${\mathrm{\θ}}_n=-\[\frac{(T_{2,n}-T_{1,n})-(T_{4,n}-T_{3,n})}{2}\]---\left(7\right)$

当作为替代使用时钟模型C(t)＝(1+α)S(t)+θ时，可推导出没有负号的以上方程。无论按照哪种方式，可推导出相同的经典偏移方程。

状态(处理)方程的发展——时钟处理模型

通常，通过以下随机微分方程对时钟的时钟偏离建模

$\frac{dX\left(t\right)}{dt}=FX\left(t\right)+\mathrm{\ξ}\left(t\right),---\left(8\right)$

其中F是L×L真实矩阵并且X(t)是状态矢量

$X\left(t\right)=\left[\begin{array}{c}{x}_1\left(t\right)\\ .\\ .\\ .\\ x_L\left(t\right)\end{array}\right].---\left(9\right)$

量x₁(t)＝θ(t)表示时间偏差(或者时间偏移)，x₂(t)＝α(t)是时钟频率偏差，并且x₃(t)＝a(t)表示所谓的频率漂移或老化。

在该模型中，ξ(t)是统计学上独立的零均值白高斯噪声的矢量，

$\mathrm{\ξ}\left(t\right)=\left[\begin{array}{c}{\mathrm{\ξ}}_1\left(t\right)\\ .\\ .\\ .\\ {\mathrm{\ξ}}_L\left(t\right)\end{array}\right],---\left(10\right)$

其中对于k＝1，2，...，L，该矢量的各个元素具有自相关

E[ξ_k(t₁)ξ_k(t₂)]＝g_kδ(t₁-t₂)，(11)

并且其中δ(t)是狄拉克δ函数。ξ(t)的自相关矩阵因此为

E[ξ(t₁)ξ^T(t₂)]＝Gδ(t₁-t₂)，(12)

并且其中G是对角矩阵

可按照积分形式写随机微分方程(8)

$X\left(t\right)=\mathrm{\Φ}(t-t_0)X\left(t_0\right)+{\∫}_{t_0}^t\mathrm{\Φ}(t-\mathrm{\τ})\mathrm{\ξ}\left(\mathrm{\τ}\right)d\mathrm{\τ},---\left(14\right)$

其中t₀是初始时间，并且Φ(t)是过渡矩阵

Ф(t)＝e^Ft。(15)

时间轴可用采样时间Δt来离散化，从而t取值t_n＝nΔt，其中n＝0，1，2，...。(14)的积分形式在离散时间变成

$\begin{array}{c}X\left(t_n\right)=\mathrm{\Φ}\left(\mathrm{\Δt}\right)X\left(t_{n-1}\right)+w\left(t_{n-1}\right)\\ \⇒{\stackrel{\‾}{X}}_n=\stackrel{\‾}{A}{\stackrel{\‾}{X}}_{n-1}.+{\stackrel{\‾}{w}}_n\end{array}---\left(16\right)$

其中

$w\left(t_{n-1}\right)={\∫}_{t_{n-1}}^{t_n}\mathrm{\Φ}(t-\mathrm{\τ})\mathrm{\ξ}\left(\mathrm{\τ}\right)d\mathrm{\τ},---\left(17\right)$

是高斯随机变量，其均值为零，并且其协方差矩阵为

$\begin{array}{c}{Q}_{n-1}=Q\left(t_{n-1}\right)=E\[\left(w\left(t_{n-1}\right)-E\[w\left(t_{n-1}\right)\]\right){\left(w\left(t_{n-1}\right)-E\[w\left(t_{n-1}\right)\]\right)}^T\]\\ ={\∫}_{t_{n-1}}^{t_n}{\∫}_{t_{n-1}}^{t_n}\mathrm{\Φ}\left(t_n-\mathrm{\τ}\right)E\[\mathrm{\ξ}\left(\mathrm{\τ}\right){\mathrm{\ξ}}^T\left(s\right)\]{\mathrm{\Φ}}^T\left(t_n-s\right)dsd\mathrm{\τ}\\ ={\∫}_{t_{n-1}}^{t_n}{\∫}_{t_{n-1}}^{t_n}\mathrm{\Φ}\left(t_n-\mathrm{\τ}\right)G\mathrm{\δ}\left(s-\mathrm{\τ}\right){\mathrm{\Φ}}^T\left(t_n-s\right)dsd\mathrm{\τ}\\ ={\∫}_{t_{n-1}}^{t_n}\mathrm{\Φ}\left(t_n-\mathrm{\τ}\right){\mathrm{G\Φ}}^T\left(t_n-\mathrm{\τ}\right)d\mathrm{\τ}\end{array}---\left(18\right)$

对于L＝2，应该注意在时钟处理模型(16)中的矩阵A是

$\stackrel{\‾}{A}\≈\mathrm{\Φ}\left(\mathrm{\Δ}t\right)=I+F\mathrm{\Δ}t=\left[\begin{array}{cc}1& \mathrm{\Δ}t\\ 0& 1\end{array}\right],F=\left[\begin{array}{cc}0& 1\\ 0& 0\end{array}\right],---\left(19\right)$

其为(15)的幂级数展开的逼近，其为

$\mathrm{\Φ}\left(t\right)=e^{Ft}=I+F\left(t\right)+F^2\frac{t^2}{2}+F^3\frac{t^2}{6}+...+F^{q-1}\frac{t^{q-1}}{(q-1)!},---\left(20\right)$

其中I是单位矩阵。可通过以下双态动态模型来描述该***

${\stackrel{\‾}{X}}_n=\left[\begin{array}{c}{\mathrm{\θ}}_n\\ {\mathrm{\α}}_n\end{array}\right]=\left[\begin{array}{cc}1& \mathrm{\Δ}t\\ 0& 1\end{array}\right]\left[\begin{array}{c}{\mathrm{\θ}}_{n-1}\\ {\mathrm{\α}}_{n-1}\end{array}\right]+\left[\begin{array}{c}{w}_{\mathrm{\θ},n}\\ w_{\mathrm{\α},n}\end{array}\right]=\stackrel{\‾}{A}{\stackrel{\‾}{X}}_{n-1}+{\stackrel{\‾}{w}}_n---\left(21\right)$

其中A_n是已知的2×2状态过渡矩阵，并且是2维零均值白处理噪声矢量。如果Sync消息之间的时间固定并且被作为采样间隔，则得到ΔT_n＝(T_1，n-T_1，n-1)＝Δt。

处理噪声的方差与从装置处的本地PLL的振荡器有关。对于电信时钟同步问题，本地振荡器通常是温补晶体振荡器(TCXO)或者炉控晶体振荡器(OCXO)。诸如卡尔曼滤波器的滤波技术的实际实现要求获得测量噪声(v_n＝(γ_n-ε_n))的方差和处理噪声矢量的协方差的良好估计。下面描述的实施方式(其具有类似卡尔曼滤波器的结构)的指数加权递归最小二乘(RLS)技术[4]不要求这种噪声方差估计，但是像具有已知噪声统计量的卡尔曼滤波器那样得出准确的时钟参数估计。

指数加权递归最小二乘(RLS)法——类似卡尔曼滤波器的算法

在该实施方式中使用的模型基于所估计的变量的状态空间表示。状态空间暗含着在各个时间点，被建模的处理由对所有关注的有关量的状态变量进行概括的矢量来描述。在下面要描述的滤波算法使用***的时间行为的这种模型以及一些***变量的有噪声观察或测量结果，来产生对所有状态变量的最佳估计。这些估计随后在处理模型中使用以确定针对未来时间段的状态估计。

考虑由下面的方程对(即，分别为时钟测量方程和状态方程)描述的状态空间模型，

$y_n={\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{X}}_n+v_n,$

${\stackrel{\‾}{X}}_{n+1}={\stackrel{\‾}{A}}_n{\stackrel{\‾}{X}}_n+{\stackrel{\‾}{w}}_n,$

其中n是非负时间索引，是已知的M×M状态过渡矩阵，是M维状态(或参数)矢量，是M维零均值白处理噪声矢量，y_n是测量结果，是已知的M维测量矢量，v_n是具有零均值和方差为的白噪声，并且T表示转置。可按照的项表达测量结果y_n-i，0≤i≤n。从状态方程(16)，可写作

$\begin{array}{c}{\stackrel{\‾}{X}}_{n-i}=\left(\underset{j=n-i}{\overset{n-1}{\Π}}{\stackrel{\‾}{A}}_j^{-1}\right){\stackrel{\‾}{X}}_n-\underset{m=1}{\overset{i}{\mathrm{\Σ}}}\left(\underset{j=n-i}{\overset{n-m}{\Π}}{\stackrel{\‾}{A}}_j^{-1}\right){\stackrel{\‾}{w}}_{n-m}\\ ={\mathrm{\Ψ}}_n^{n-i}{\stackrel{\‾}{X}}_n-\underset{m=1}{\overset{i}{\mathrm{\Σ}}}{\mathrm{\Ψ}}_{n-m+1}^{n-i}{\stackrel{\‾}{w}}_{n-m}\end{array}---\left(22\right)$

其中

${\mathrm{\Ψ}}_a^b=\underset{l=b}{\overset{a-1}{\Π}}{A}_l^{-1},a>b---\left(23\right)$

是将和关联的后向过渡矩阵。在(5)中使用(22)，得到

$y_{n-i}={\stackrel{\‾}{D}}_{n-i}^T{\mathrm{\Ψ}}_n^{n-i}{\stackrel{\‾}{X}}_n-{\stackrel{\‾}{D}}_{n-i}^T\left(\underset{m=1}{\overset{i}{\mathrm{\Σ}}}{\mathrm{\Ψ}}_{n-m+1}^{n-i}{\stackrel{\‾}{w}}_{n-m}\right)+v_{n-i}.---\left(24\right)$

如果将(n+1)×M矩阵定义为

${\stackrel{\‾}{H}}_n={\[{\stackrel{\‾}{D}}_n,{\left({\stackrel{\‾}{D}}_{n-1}^T{\mathrm{\Ψ}}_n^{n-1}\right)}^T,{\left({\stackrel{\‾}{D}}_{n-2}^T{\mathrm{\Ψ}}_n^{n-2}\right)}^T,...,{\left({\stackrel{\‾}{D}}_0^T{\mathrm{\Ψ}}_n^0\right)}^T\]}^T,---\left(25\right)$

并且将该(n+1)维向量定义为

${\stackrel{\‾}{u}}_n=-{\[0,{\left({\stackrel{\‾}{D}}_{n-1}^T{\mathrm{\Ψ}}_n^{n-1}{\stackrel{\‾}{w}}_{n-1}\right)}^T,{\left({\stackrel{\‾}{D}}_{n-2}^T\underset{m=1}{\overset{2}{\mathrm{\Σ}}}{\mathrm{\Ψ}}_{n-m+1}^{n-2}{\stackrel{\‾}{w}}_{n-m}\right)}^T,...,{\left({\stackrel{\‾}{D}}_0^T\underset{m=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\Ψ}}_{n-m+1}^0{\stackrel{\‾}{w}}_{n-m}\right)}^T\]}^T,---\left(26\right)$

则测量结果的集合{y_n，y_n-1，...，y₀}可表达为

${\stackrel{\‾}{Y}}_n={\stackrel{\‾}{H}}_n{\stackrel{\‾}{X}}_n+{\stackrel{\‾}{n}}_n,---\left(27\right)$

其中(其具有零平均值)并且

考虑是待估计的M维向量。在最小二乘(LS)估计中，问题是找到矢量的估计作为测量结果的线性组合，使得估计将以下价值函数最小化：

$J\left(\hat{X}\right)={(\stackrel{\‾}{Y}-\stackrel{\‾}{H}\hat{X})}^T\stackrel{\‾}{W}(\stackrel{\‾}{Y}-\stackrel{\‾}{H}\hat{X}),---\left(28\right)$

其中是M×M权矩阵。该问题的解决方案可表达为

$\hat{X}={\left({\stackrel{\‾}{H}}^T\stackrel{\‾}{WH}\right)}^{-1}{\stackrel{\‾}{H}}^T\stackrel{\‾}{WY},---\left(29\right)$

其为的无偏估计。矩阵可设为等于对角矩阵

$\stackrel{\‾}{W}=diag\[1,\mathrm{\λ},{\mathrm{\λ}}^2,...,{\mathrm{\λ}}^{M-1}\]---\left(30\right)$

其中0＜λ≤1并且λⁿ(0≤n≤M-1)表示第(n+1)对角元素。这得到将价值函数(28)变为

$J\left(\hat{X}\right)=\underset{n=0}{\overset{M-1}{\Σ}}{\mathrm{\λ}}^{M-1-n}{(y_n-{\stackrel{\‾}{D}}_n^T\hat{X})}^2,---\left(31\right)$

其中y_n是的第(M-n)元素，并且是的第(M-n)行。该价值函数得到具有遗忘因子λ的指数加权RLS算法。

使用(29)，可将状态矢量估计为

${\hat{X}}_{n,n}={\left({\stackrel{\‾}{H}}_n^T{\stackrel{\‾}{W}}_n{\stackrel{\‾}{H}}_n\right)}^{-1}{\stackrel{\‾}{H}}_n^T{\stackrel{\‾}{W}}_n{\stackrel{\‾}{Y}}_n={\stackrel{\‾}{P}}_{n,n}{\stackrel{\‾}{H}}_n^T{\stackrel{\‾}{W}}_n{\stackrel{\‾}{Y}}_n,---\left(32\right)$

其中表示基于直至时间n为止的数据对的估计，是(n+1)×(n+1)权重矩阵，并且与卡尔曼滤波的情况不同，矩阵不是协方差矩阵。矩阵可表达为

${\stackrel{\‾}{W}}_n=diag\[1,{\mathrm{\λ}}_n,\underset{m=1}{\overset{2}{\Π}}{\mathrm{\λ}}_{n-m+1},...,\underset{m=1}{\overset{n}{\Π}}{\mathrm{\λ}}_{n-m+1}\].---\left(33\right)$

当λ₁＝λ₂＝…＝λ_n＝λ时，权重变为等于(30)。

为了减少(32)中针对的计算量，需要开发递归版本。将从(16)代入(27)中的得到

$\begin{array}{c}{\stackrel{\‾}{Y}}_{n-1}={\stackrel{\‾}{H}}_{n-1}{\stackrel{\‾}{X}}_{n-1}+{\stackrel{\‾}{n}}_{n-1}\\ ={\stackrel{\‾}{H}}_{n-1}{\stackrel{\‾}{A}}_{n-1}^{-1}{\stackrel{\‾}{X}}_n-{\stackrel{\‾}{H}}_{n-1}{\stackrel{\‾}{A}}_{n-1}^{-1}{\stackrel{\‾}{w}}_{n-1}+{\stackrel{\‾}{n}}_{n,n-1}\\ ={\stackrel{\‾}{H}}_{n,n-1}{\stackrel{\‾}{X}}_n+{\stackrel{\‾}{n}}_{n,n-1}\end{array}---\left(34\right)$

其中并且给定{y₀，y₁，…y_n-1}，对的一步预测可表达为

${\hat{X}}_{n,n-1}={\left({\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{H}}_{n,n-1}\right)}^{-1}{\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{Y}}_{n-1}={\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{Y}}_{n-1}---\left(35\right)$

其中

${\stackrel{\‾}{P}}_{n,n-1}={\left({\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{H}}_{n,n-1}\right)}^{-1}---\left(36\right)$

并且是权重矩阵。针对一步预测，对矩阵的合理选择为：

${\stackrel{\‾}{W}}_{n,n-1}={\mathrm{\λ}}_n{\stackrel{\‾}{W}}_{n-1}.---\left(37\right)$

通过该权重矩阵，包括最近的y_n-1在内的每个输入可针对一步预测被适当地加权。这里的目标是用和表达将(37)代入(36)，并且使用关系和给出

$\begin{array}{c}{\stackrel{\‾}{P}}_{n,n-1}=\frac{{\[{\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n-1}{\stackrel{\‾}{H}}_{n,n-1}\]}^{-1}}{{\mathrm{\λ}}_n}\\ =\frac{{\[{\left({\stackrel{\‾}{A}}_{n-1}^{-1}\right)}^T{\stackrel{\‾}{H}}_{n-1}^T{\stackrel{\‾}{W}}_{n-1}{\stackrel{\‾}{H}}_{n-1}{\stackrel{\‾}{A}}_{n-1}^{-1}\]}^{-1}}{{\mathrm{\λ}}_n}\\ =\frac{{\[{\left({\stackrel{\‾}{A}}_{n-1}^T\right)}^{-1}{\stackrel{\‾}{P}}_{n-1,n-1}^{-1}{\stackrel{\‾}{A}}_{n-1}^{-1}\]}^{-1}}{{\mathrm{\λ}}_n}\end{array}$

在以上结果中使用关系(abc)^-1＝c^-1b^-1a^-1，给出

${\stackrel{\‾}{P}}_{n,n-1}=\frac{{\stackrel{\‾}{A}}_{n-1}{\stackrel{\‾}{P}}_{n-1,n-1}{\stackrel{\‾}{A}}_{n-1}^T}{{\mathrm{\λ}}_n}.---\left(38\right)$

将(38)代入(35)并且使用关系得到下式

$\begin{array}{c}{\hat{X}}_{n,n-1}=\left(\frac{{\stackrel{\‾}{A}}_{n-1}{\stackrel{\‾}{P}}_{n-1,n-1}{\stackrel{\‾}{A}}_{n-1}^T}{{\mathrm{\λ}}_n}\right){\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{Y}}_{n-1}\\ ={\stackrel{\‾}{A}}_{n-1}{\stackrel{\‾}{P}}_{n-1,n-1}{\left({\stackrel{\‾}{H}}_{n,n-1}{\stackrel{\‾}{A}}_{n-1}\right)}^T\frac{{\stackrel{\‾}{W}}_{n,n-1}}{{\mathrm{\λ}}_n}{\stackrel{\‾}{Y}}_{n-1}\\ ={\stackrel{\‾}{A}}_{n-1}\left({\stackrel{\‾}{P}}_{n-1,n-1}{\stackrel{\‾}{H}}_{n-1}^T{\stackrel{\‾}{W}}_{n-1}{\stackrel{\‾}{Y}}_{n-1}\right)\end{array}$

在以上结果中使用(32)，得到

${\hat{X}}_{n,n-1}={\stackrel{\‾}{A}}_{n-1}{\hat{X}}_{n-1,n-1}.---\left(39\right)$

利用以下关系

${\stackrel{\‾}{H}}_n=\left[\begin{array}{c}{\stackrel{\‾}{D}}_n^T\\ {\stackrel{\‾}{H}}_{n,n-1}\end{array}\right],{\stackrel{\‾}{W}}_n=\left[\begin{array}{cc}1& 0\\ 0& {\stackrel{\‾}{W}}_{n,n-1}\end{array}\right],{\stackrel{\‾}{Y}}_n={\left[\begin{array}{cc}{y}_n& {\stackrel{\‾}{Y}}_{n-1}^T\end{array}\right]}^T,---\left(40\right)$

可导出针对和的递归。从(32)中的矩阵

$\begin{array}{c}{\stackrel{\‾}{P}}_{n,n}={\left(\left[\begin{array}{cc}{\stackrel{\‾}{D}}_n& {\stackrel{\‾}{H}}_{n,n-1}^T\end{array}\right]\left[\begin{array}{cc}1& 0\\ 0& {\stackrel{\‾}{W}}_{n,n-1}\end{array}\right]\left[\begin{array}{c}{\stackrel{\‾}{D}}_n^T\\ {\stackrel{\‾}{H}}_{n,n-1}\end{array}\right]\right)}^{-1}\\ ={\left({\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{H}}_{n,n-1}+{\stackrel{\‾}{D}}_n{\stackrel{\‾}{D}}_n^T\right)}^{-1}\end{array}$

将矩阵求逆引理(a+bcd)^-1＝a^-1-a^-1b(da^-1b+c^-1)^-1da^-1应用于上式，给出

${\stackrel{\‾}{P}}_{n,n}={\stackrel{\‾}{P}}_{n,n-1}-\frac{{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}}{{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n+1}.---\left(41\right)$

现在，期望的估计变为

$\begin{array}{c}{\hat{X}}_{n,n}={\stackrel{\‾}{P}}_{n,n}{\stackrel{\‾}{H}}_n^T{\stackrel{\‾}{W}}_n{\stackrel{\‾}{Y}}_n\\ ={\stackrel{\‾}{P}}_{n,n}\left(\left[\begin{array}{cc}{\stackrel{\‾}{D}}_n& {\stackrel{\‾}{H}}_{n,n-1}^T\end{array}\right]\left[\begin{array}{cc}1& 0\\ 0& {\stackrel{\‾}{W}}_{n,n-1}\end{array}\right]\left[\begin{array}{c}{y}_n\\ {\stackrel{\‾}{Y}}_{n-1}\end{array}\right]\right)\\ =\left({\stackrel{\‾}{P}}_{n,n-1}-\frac{{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}}{{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n+1}\right)\left({\stackrel{\‾}{H}}_{n,n-1}^T{\stackrel{\‾}{W}}_{n,n-1}{\stackrel{\‾}{Y}}_{n-1}+{\stackrel{\‾}{D}}_n{y}_n\right)\end{array}$

其在一些运算之后，缩减为

${\hat{X}}_{n,n}={\hat{X}}_{n,n-1}+{\stackrel{\‾}{K}}_n(y_n-{\stackrel{\‾}{D}}_n^T{\hat{X}}_{n,n-1}),---\left(42\right)$

其中

${\stackrel{\‾}{K}}_n=\frac{{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n}{{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n+1},---\left(43\right)$

是卡尔曼增益矢量。式(41)还可通过进一步表达为：

${\stackrel{\‾}{P}}_{n,n}=\left(\stackrel{\‾}{I}-{\stackrel{\‾}{K}}_n{\stackrel{\‾}{D}}_n^T\right){\stackrel{\‾}{P}}_{n,n-1},---\left(44\right)$

其中是单位矩阵。式(38)、(39)、(41)、(42)和(43)构成RLS算法，总结在如下步骤中：

1、通过以下设置来初始化算法：

${\hat{X}}_{0,-1}=\stackrel{\‾}{0}$

γ是小的正常数

n＝1

2、状态预测：

${\hat{X}}_{n,n-1}={\stackrel{\‾}{A}}_{n-1}{\hat{X}}_{n-1,n-1}$

3、计算

${\stackrel{\‾}{P}}_{n,n-1}=\frac{{\stackrel{\‾}{A}}_{n-1}{\stackrel{\‾}{P}}_{n-1,n-1}{\stackrel{\‾}{A}}_{n-1}^T}{{\mathrm{\λ}}_n}$

4、计算卡尔曼增益矢量：

${\stackrel{\‾}{K}}_n=\frac{{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n}{{\stackrel{\‾}{D}}_n^T{\stackrel{\‾}{P}}_{n,n-1}{\stackrel{\‾}{D}}_n+1}$

5、更新状态估计：

${\hat{X}}_{n,n}={\hat{X}}_{n,n-1}+{\stackrel{\‾}{K}}_n(y_n-{\stackrel{\‾}{D}}_n^T{\hat{X}}_{n,n-1})$

6、计算

${\stackrel{\‾}{P}}_{n,n}=(\stackrel{\‾}{I}-{\stackrel{\‾}{K}}_n{\stackrel{\‾}{D}}_n^T){\stackrel{\‾}{P}}_{n,n-1}$

7、n＝n+1

至步骤2

针对当前的时钟同步问题，初始化因子γ选为大约百万分之几的数量级的小的正常数(例如，γ＝10^-6)。遗忘因子λ_n是合适的常数(例如，λ_n＝0.9999)或者可以是动态变化的。

所提出的用于时间和频率恢复的DPLL架构

时钟偏移和偏斜可在由服务器广播各个Sync消息之后或者在多个时段的Sync消息广播之后由客户端估计。Sync消息之间的时段可用作***的采样时段。可利用如图6所示调整了偏斜和偏移的客户端的本地时钟C来计算服务器时间图7示出了所提出的在时间客户端处的同步机制的主要块。将自由运行的本地客户端振荡器25与估计的时钟偏移和偏斜一起使用，以得到服务器时钟估计

图7所示的架构的主要目的是在从装置重建对服务器时间和频率二者的准确估计。利用采用如图8所示的相位累加器23、环路滤波器22、相位检测器21和计数器24的数字锁相环(DPLL)20来实现这一点。本地辅助硬件时钟25(具有自由运行的计数器)仅用于为PTP消息加时间戳。自由运行的计数器可在从装置用首先到达的Sync消息初始化。相位累加器23的溢出脉冲按照如图9所示驱动计数器24的方式来控制DPLL20。在锁定模式下，在相位累加器溢出脉冲产生服务器频率信号的同时，DPLL计数器24产生连续服务器时间信号，。

在启动时，DPLL20等待首先计算的服务器时间估计该第一服务器时间估计用于初始化DPLL计数器从这一点起并且当在任何离散时刻n接收到后续服务器时间估计时，DPLL20开始按照闭环方式操作。在各个服务器时间估计由从装置标记DPLL计数器读数。然后，所计算的服务器时间估计与DPLL计数器读数之间的差提供误差信号该误差信号(e(n))被发送至输出控制相位累加器23的频率的环路滤波器21。相位累加器23的输出(溢出脉冲)继而提供从装置的时钟频率并且还驱动DPLL计数器24。一会儿过后，误差项期望收敛至零，这意味着DPLL20锁定至进入的主时基(频率和时间二者)。

如图8所示，在各个相位累加器溢出(输出)脉冲，DPLL计数器24由相位累加器溢出脉冲的标称周期(例如，125MHz标称相位累加器溢出输出频率为8ns)进行累加。按照如图10所例示DPLL计数器跟随计算出的服务器时间估计的方式来控制DPLL20。如图10所例示，控制递增的定时(由相位累加器溢出输入字控制)以使得增量沿着主时钟时基下降。通过正确控制的计数器增量的定时，DPLL20提供频率信号(相位累加器23溢出的输出)和时间信号(DPLL计数器24的输出)二者。

随后，可利用如图11和图12所示的其它技术调节频率信号，以提供符合终端应用的抖动需求的信号。可根据相位累加器溢出输出构造具有接口形式因素的各种形式的信号(方波、正弦波等)。可由另一模拟PLL(APLL)调节信号以由终端应用使用。也可将DPLL计数器输出(时间信号)格式化成各种时间标准信号。

可通过两个分离的计数器(自由运行的计数器和DPLL计数器)或者一个计数器(仅DPLL计数器)实现图7和图8所示的架构。两种架构均可用于合成时间和频率信号这二者：

·具有分离的计数器(自由运行的计数器25和DPLL计数器24)的架构：这些计数器的初始化如下实现：

o自由运行的计数器。该计数器仅用于对PTP消息加时间戳。其可在从装置由首先到来的Sync消息初始化。其初始值也可为用户配置的。

oDPLL计数器。该计数器针对锁定至主装置的从装置提供连续时间信号。由从装置用首先计算出的服务器时间估计进行初始化。

·具有一个计数器(仅DPLL计数器24)的架构：在该构造中，DPLL计数器用作时间戳计数器以及针对从装置提供连续时间信号的计数器这二者。初始化执行如下：

○DPLL按照开环模式开始。用首先到来的Sync消息来初始化DPLL计数器。其也可为用户构造的初始值。

○DPLL仍在开环下运行，但是DPLL计数器由相位累加器溢出脉冲驱动。DPLL计数器读数现在用于对后续PTP消息加时间戳。

○计算第一服务器时间估计并且将DPLL计数器复位至该第一时间估计。DPLL现在进入闭环模式。DPLL计数器读数用于为后续PTP消息加时间戳。

○后续服务器时间估计用作DPLL参考信号。DPLL继续在闭环模式操作。

相位累加器23是一种每当其接收到时钟脉冲时将存储在其中的数递增的可变模数计数器。当计数器溢出时其环绕(wraparound)，从而使得相位累加器的输出相连。添加的增量φ越大，累加器溢出越快，这导致更高的输出频率。相位累加器的输出频率f_DDS是***时钟频率f_o、相位累加器中的比特数q和相位增量值φ的函数。相位增量φ是无符号值。

$f_{DDS}=\frac{f_o}{2^q}\mathrm{\φ}---\left(45\right)$

从该方程中可以看出，相位累加器的频率分解是f_res＝f_o/2^q。假设相位累加器由对应于标称频率f_DDS＝f_nom的控制输入φ_nom操作，则可以从上面的讨论中看出，将量-φ_corr加至φ_nom(即，φ_DDS＝φ_nom-φ_corr)得到输出频率的减小，f_DDS＝f_nom-Δf，而将量+φ_corr加至(即，φ_DDS＝φ_nom+φ_corr)得到输出频率的增大，f_DDS＝f_nom+Δf。因此，通过适当地控制加至φ_nom的量φ_corr，可相应地控制相位累加器的输出频率f_DDS。

可根据需要采用相位累加器23作为如图11和图12所示的可从中产生模拟信号的直接数字合成器(DDS)30的一部分。尽管存在许多变化，但是可将常规DDS架构看作由仅三个普通数字组件构成的简单组件：相位累加器23(加法器/累加器)、映射装置26(诸如只读存储器(ROM)或随机存取存储器(RAM))和数模转换器(DAC)27。在许多情况下，低通滤波器28在DAC的输出实现，但是该组件通常不被看作是DDS的一部分。由于奈奎斯特定理，基准时钟f_o必须以高于合成时钟的更高的频率操作。

DPLL环路滤波器参数

为了分析和设计控制***，有必要获得***的量化数学描述或者模型。模型是***变量之间的一组数学关系。因为关注的***本质上是动态的，因此描述方程通常是差分(或差)方程。该模型或者差分方程组描述了***的动态行为。在建模中开发的差分方程通常是非线性的。因为它们比线性方程显著地更难求解，所以线性模型通常足够。线性化是找到接近非线性模型的线性模型的处理。与非线性模型相比，对于线性模型，分析和控制设计二者均容易得多。利用线性模型的理由是，如果小信号线性模型接***衡(equilibrium)(稳态)有效，并且稳定，则可存在包含该平衡的小区域，在该区域中非线性***是稳定的。换言之，假设从平衡的偏差是小的，从而非线性函数可由线性函数逼近。

此外，如果可将这些方程线性化，则可利用拉普拉斯变换来简化求解方法。实际上，***的复杂性和对所有相关因素的全面知识的缺乏使得有必要引入关于***操作的假设。因此，将发现的是，考虑物理***、勾画一些必要的假设和使***线性化是有用的。然后利用诸如拉普拉斯变换这样的数学工具，获得描述***的操作的解。

锁相环(PLL)主要是如图12所示的反馈控制***，并且其部件简化为图13中的“黑盒子”操作。因此，需要DDS30和相位检测器(PD)21的数学模型(例如，传递函数形式)来确定环路滤波器22的参数。在这一部分中，采用z变换技术来分析PLL的一般跟踪(即，备用状态)行为。在稳态假设下，相位误差采样小，并且一般非线性差方程可通过可由z-变换技术解开的线性方程逼近。应该注意，当PLL获得锁定，并且不由应用于其基准输入的大的相位步长、频率步长或者相位噪声拉出时，可通过线性模型分析其性能。

假设相位误差θ_e(即，θ_e(n)＝θ_s(n)-θ_OSC(n)是振荡器时钟相位θ_OSC(n)与基准时钟相位θ_s(n)之间的差)在有限范围内，作为反馈控制***的PLL可被简化为线性反馈控制***。该假设对于多数应用而言是合理的，这是因为真实PLL具有受约束和有限的范围(表达为标称操作频率的百万分之几，ppm)，在该范围外，不能确保锁定。因此，PLL的小信号线性分析可对于学***衡行为和稳定属性有用。在以下部分，针对相位检测器21、受控制的振荡器以及在给定的环路滤波器22的一些普通结构的情况下PLL作为整体开发控制模型。分析将进一步提供用于确定将符合预定设计和性能需要的环路滤波器22的参数的设计工序。

DDS的控制模型

在DDS30中，标称控制字φ_nom产生对应的标称频率f_nom。假设在离散时间n控制输入φ_nom改变量φ_corr。注意，改变在下一离散间隔起作用。这种改变得到输出频率

$f_{DDS}\left(n\right)=\frac{f_o}{2^q}({\mathrm{\φ}}_{nom}+{\mathrm{\φ}}_{corr}(n-1\left)\right)=f_{nom}+\mathrm{\Δ}f\left(n\right),---\left(46\right)$

或

$f_{DDS}\left(n\right)=f_{nom}+\frac{f_o}{2^q}{\mathrm{\φ}}_{corr}(n-1).---\left(47\right)$

其对应于以下角频率

${\mathrm{\ω}}_{DDS}\left(n\right)={\mathrm{\ω}}_{nom}+\frac{2{\mathrm{\πf}}_o}{2^q}{\mathrm{\φ}}_{corr}(n-1).---\left(48\right)$

以上方程还可写作

ω_DDS(n)＝ω_nom+K_DDSφ_corr(n-1)＝ω_nom+Δω(n)，(49)

其中

$K_{DDS}=K_{OSC}=\frac{2{\mathrm{\πf}}_o}{2^q}---\left(50\right)$

是DDS增益。通过定义，DDS的相位θ_DDS由频率变化Δω(n)＝ω_DDS(n)-ω_nom上的积分提供：

${\mathrm{\θ}}_{DDS}\left(n\right)=\underset{i=0}{\overset{n}{\Σ}}\mathrm{\Δ}\mathrm{\ω}\left(i\right)=K_{DDS}\underset{i=0}{\overset{n}{\Σ}}{\mathrm{\φ}}_{corr}\left(i\right).---\left(51\right)$

DDS在数字PLL中表现为数字积分器，如同VCO在模拟PLL中表现为模拟积分器。

假定ω_DDS＝2πf_oφ/2^q，DDS增益另选地可被获得为K_DDS＝dω_DDS/dφ＝2πf_o/2^q。

从以上积分式中，DDS的传递函数在z域中被给出为

$G_{DDS}\left(z\right)=\frac{{\mathrm{\Θ}}_{DDS}\left(z\right)}{{\mathrm{\Φ}}_{corr}\left(z\right)}=K_{DDS}.\frac{z^{-1}}{1-z^{-1}}=2^{1-q}{\mathrm{\πf}}_o.\frac{z^{-1}}{1-z^{-1}}---\left(52\right)$

其中z^-1表示延迟运算子(即，z^-1x(n)＝x(n-1))，并且Θ_DDS(z)＝Θ_OSC(z)和Φ_corr(z)分别是θ_DDS(n)＝θ_OSC(n)和φ_corr(n)的z变换。

相位检测器的控制模型

T_sp定义为估计之间的标称时间间隔，以及为在时间客户端处的DPLL的采样间隔。采样间隔T_sp由标称客户端时钟f_nom量化为M＝T_sp/t_nom个单位，如图14所示。也就是说，每M个DPLL时钟脉冲进行计算。相位检测器(PD)21在频率f_sp＝1/T_sp操作。作为在DPLL20测量和控制的基准操作间隔的间隔T_sp等于2π弧度，或者M个标称客户时钟周期(见图14)。在客户端在本地时钟时间L_clock的服务器时间的估计表示为

如果描绘PD输出e相对于相位误差θ_e，则得到图15所示的锯齿函数。这代表PD的特征曲线，并且覆盖大于2π或者小于-2π的相位误差。该曲线以时间段2π为周期。PD输出是针对输入相位差的从-2π至2π弧度的整个范围并且具有在M处的最大输出的理想直线，这是由于假设应用了线性控制***模型的PLL的稳态或锁定状态操作。注意，在锁定状态，所有频率接近它们的理想值。在这种状态下，相位误差范围[-2π，2π]映射至误差范围[-M，M]。

PD特征曲线的斜率等效于PD的增益。从图15中，斜率由以下给出

$K_{PD}=\frac{M}{2\mathrm{\π}}.---\left(53\right)$

当相位误差局限于范围-2π＜θ_e＜2π时，PD输出变为

$e\left({\mathrm{\θ}}_e\right)=\frac{M}{2\mathrm{\π}}{\mathrm{\θ}}_e=K_{PD}{\mathrm{\θ}}_e---\left(54\right)$

PD21测量由时间客户端DPLL控制的振荡器相位θ_OSC(n)与时间服务器(基准)时钟相位θ_s(n)之间的相位差θ_e(n)＝θ_s(n)-θ_OSC(n)，并且开发与该相位差θ_e(n)成比例的输出e(n)。这个运算可表达为

e(n)＝K_PD·θ_e(n)(55)

误差信号输出e(n)随后传送至环路滤波器G_LF(z)，以被处理为滤波后的误差PD的传递函数为

$G_{PD}\left(z\right)=\frac{E\left(z\right)}{{\mathrm{\Θ}}_e\left(z\right)}=K_{PD}=\frac{M}{2\mathrm{\π}},---\left(56\right)$

其中E(z)和Θ_e(z)分别是e(z)和θ_e(z)的z变换。

数字环路滤波器的控制模型

来自PD21误差信号e(n)传送至数字环路滤波器22，其输出用于调整振荡器的频率f_OSC＝f_DDS。许多形式的滤波器可用作环路滤波器22。例如，数字环路滤波器22可实现为具有以下传递函数G_LF(z)的比例加积分(PI)滤波器

$G_{LF}\left(z\right)=\frac{\tilde{E}\left(z\right)}{E\left(z\right)}=K_1+\frac{K_2}{1-z^{-1}},---\left(57\right)$

其中是滤波器输出的z变换，并且K₁和K₂分别是比例和积分路径增益。该传递函数等效于离散时间控制方程

$\tilde{e}\left(n\right)=\tilde{e}(n-1)+K_1\left(e\right(n)-e(n-1\left)\right)+K_{2}e\left(n\right)$

环路滤波器为产生二阶PLL的PI滤波器。比例增益K₁和积分增益K₂确定滤波器响应。滤波器增益K₁和K₂可在需要时在运行时动态调整，为了快速锁定(获取模式)在启动处理中具有较大增益，为了更好的稳定性和稳态误差(跟踪模式)在稳态下具有较小增益。

在时间客户端处的PLL的控制模型

具有设计良好的环路滤波器22的DPLL20最终可消除相位差，并且将所控制的振荡器输出相位和频率锁定至基准。图16和图17示出了作为闭环反馈控制***的DPLL20。

该部分描述了用于利用标准控制理论原理合成DPLL的方法。设计基于连续时间***的数字化，从而特定差分方程的s平面极点和零点利用匹配零极点(MPZ)方法映射至对应差分方程的z平面极点和零点。

PLL在s域的线性二阶模型

模拟或者连续时间PLL40(见图18)由相位检测器41、环路滤波器42和电压控制振荡器(VCO)43构成。相位检测器41可简单地由恒定增益K_PD表达。VCO43可在拉普拉斯域中被模型化为完整积分器，即G_VCO(s)＝K_VCO/s，其中K_VCO为其增益。环路滤波器42可在拉普拉斯域中指定为F(s)。

在没有噪声时，闭环传递函数和归一化相位误差响应在拉普拉斯域中分别指定为

$H\left(s\right)=\frac{{\mathrm{\Θ}}_{VCO}\left(s\right)}{{\mathrm{\Θ}}_s\left(s\right)}=\frac{K_{PD}{K}_{VCO}F\left(s\right)}{s+K_{PD}{K}_{VCO}F\left(s\right)},---\left(58\right)$

以及

$\frac{{\mathrm{\Θ}}_e\left(s\right)}{{\mathrm{\Θ}}_s\left(s\right)}=\frac{{\mathrm{\Θ}}_s\left(s\right)-{\mathrm{\Θ}}_{VCO}\left(s\right)}{{\mathrm{\Θ}}_s\left(s\right)}=\frac{s}{s+K_{PD}{K}_{VCO}F\left(s\right)}=1-H\left(s\right),---\left(59\right)$

其中Θ_VCO(s)、Θ_s(s)以及Θ_e(s)分别是VCO相位θ_VCO(t)、基准信号相位θ_s(t)和相位误差θ_e(t)的拉普拉斯变换。

环的阶等于环结构中的完整积分器的数量。由于VCO43模型化为完整积分器，环至少为1阶的。如果环路滤波器42包含一个完整积分器，则环为2阶。

环的阶可示为极大地影响环的稳态性能。稳态相位误差可通过终值定理从(59)中容易地确定，即，

$\underset{t\→\mathrm{\∞}}{\lim }{\mathrm{\θ}}_e\left(t\right)=\underset{s\→0}{\lim }{\mathrm{s\Θ}}_e\left(s\right)=\underset{s\→0}{\lim }\frac{s^2{\mathrm{\Θ}}_s\left(s\right)}{s+K_{PD}{K}_{VCO}F\left(s\right)}---\left(60\right)$

稳态误差被定义为在瞬态响应消失之后VCO相位从基准的偏差。稳态误差仅为θ_e(∞)。通过(60)可以示出一阶或更高阶环将以零稳态误差跟踪初始相位偏移。此外，需要二阶***来跟踪频率步长，而需要采用三阶环来以零稳态误差跟踪加速相位。

本文考虑具有以下传递函数的二阶滞后超前滤波器(还称作比例-积分(PI)滤波器)

$F\left(s\right)=\frac{1+{\mathrm{s\τ}}_2}{{\mathrm{s\τ}}_1}---\left(61\right)$

其中τ₁和τ₂是滤波器的时间常数。滤波器具有在s＝0的极点，因此作为积分器操作。其在零频率(至少在理论上)具有无限增益。获得具有该滤波器的PLL的闭环传递函数如下

$H\left(s\right)=\frac{2{\mathrm{\ζ\ω}}_{n}s+{\mathrm{\ω}}_n^2}{s^2+2{\mathrm{\ζ\ω}}_{n}s+{\mathrm{\ω}}_n^2}=\frac{2{\mathrm{\ζ\ω}}_{n}s+{\mathrm{\ω}}_n^2}{(s-s_0)(s-s_1)},---\left(62\right)$

其中ω_n和ζ分别是固有频率和阻尼因数，并且用K_PD、K_VCO、τ₁和τ₂指定为

${\mathrm{\ω}}_n=\sqrt{\frac{K_{PD}{K}_{VCO}}{{\mathrm{\τ}}_1}},---\left(63\right)$

$\mathrm{\ζ}=\frac{{\mathrm{\ω}}_n{\mathrm{\τ}}_2}{2}.---\left(64\right)$

这两个参数通常用于指定***的性能需求。闭环***的极点为

$s_{0,1}=-{\mathrm{\ζ\ω}}_n\±{\mathrm{j\ω}}_n\sqrt{1-{\mathrm{\ζ}}^2}.---\left(65\right)$

当ζ＞1时，极点是实的；并且当ζ＜1时，极点是复数和共轭的。当ζ＝1时，极点重复并且是实的，并且该条件称作临界阻尼。当ζ＜1时，响应欠阻尼并且极点为复数。

当ζ接近零时，闭环***的瞬态响应随着极点接近虚轴而振荡增大。以上模型可在连续时间域直接应用于PLL。但对于基于采样的数据的***，需要使用离散时间模型。

PLL在z域的线性二阶模型

可基于将***量化的非线性忽略的条件推导用于整个PLL40的线性化时不变的逼近传递函数。PD41、环路滤波器42和受控的振荡器43的z域表现分别如下

G_PD(z)＝K_PD，(66)

$G_{LF}\left(z\right)=K_1+\frac{K_2}{1-z^{-1}}=\frac{(K_1+K_2)z-K_1}{z-1},---\left(67\right)$

$G_{OSC}\left(z\right)=\frac{K_{OSC}{z}^{-1}}{1-z^{-1}}=\frac{K_{OSC}}{z-1}.---\left(68\right)$

利用这些传递函数，PLL40的闭环传递函数变为

$H\left(z\right)=\frac{G_{PD}\left(z\right)G_{LF}\left(z\right)G_{OSC}\left(z\right)}{1+G_{PD}\left(z\right)G_{LF}\left(z\right)G_{OSC}\left(z\right)},---\left(69\right)$

或者

$H\left(z\right)=\frac{K_{PD}{K}_{OSC}(K_1+K_2)z-K_{PD}{K}_{OSC}{K}_1}{z^2+\[K_{PD}{K}_{OSC}(K_1+K_2)-2\]z-(K_{PD}{K}_{OSC}{K}_1-1)}.---\left(70\right)$

现在，将匹配零极点(MPZ)方法应用于H(s)以获得具有H(z)形式(或与离散传递函数有关)的离散时间***H₂(z)。从该关系，将推导针对环路滤波器增益K₁和K₂的闭合形式表达。

匹配零极点(MPZ)方法

目的是将满足由ω_n和阻尼因数ζ指定的性能需求的***映射至z域中的对应模型。MPZ方法将模拟***的s平面极点和零点直接映射至离散时间***的对应的z平面极点和零点。这里，使用修改的MPZ(MMPZ)方法，其可以描述如下：

1、利用关系将s平面极点和零点映射至z平面，其中T_sp是采样间隔。

H(s)的在的极点将映射至H₂(z)的在的极点。H(s)的在的极点将映射至H₂(z)的在的极点。

在s＝-ω_n/2ζ的零点将被映射至H₂(z)的在的零点。

2、形成z平面中的离散时间传递函数，其具有先前步骤中确定的极点和零点。

$H_2\left(z\right)=\frac{K_{DC}(z-e^{-{\mathrm{\ω}}_n{T}_{sp}/2\mathrm{\ζ}})}{(z-e^{T_{sp}(-{\mathrm{\ω}}_n\mathrm{\ζ}+{\mathrm{j\ω}}_n\sqrt{1-{\mathrm{\ζ}}^2})})(z-e^{T_{sp}(-{\mathrm{\ω}}_n\mathrm{\ζ}-{\mathrm{j\ω}}_n\sqrt{1-{\mathrm{\ζ}}^2})})},---\left(71\right)$

其中K_DC是H₂(z)的直流或者低频增益。

3、将离散时间***H₂(z)的直流或低频增益设为等于连续时间***H(s)的那些。

通常使用终值定理来找到提供其拉普拉斯变换或z变换的时间函数的稳态值。假设具有函数x(t)，则s域中定理状态为

$\underset{t\→\mathrm{\∞}}{\lim }x\left(t\right)=\underset{s\→0}{\lim }sX\left(s\right),---\left(72\right)$

其中X(s)是x(t)的拉普拉斯变换，并且只要sX(s)的全部极点在s平面的左半平面(LHP)中。在z域中，定理状态为

$\underset{k\→\mathrm{\∞}}{\lim }x\left({\mathrm{kT}}_{sp}\right)=\underset{z\→1}{\lim }(1-z^{-1})X\left(z\right),---\left(73\right)$

其中X(z)是x(t)的z变换，并且如果(1-z^-1)X(z)的全部极点在单位圆内，则定理也可用于找到***的直流增益。在所有瞬变衰减之后，直流增益是***的输出与输入(输入推测常数)的比率。为了找到直流增益，假设存在单位步骤输入并且使用终值定理计算输出的稳态值。因此，对于传递函数为G(s)的***，直流增益(DCgain)定义为

$DCgain=\underset{s\→0}{\lim }sG\left(s\right)\frac{1}{s}=\underset{s\→0}{\lim }G\left(s\right),---\left(74\right)$

而对于传递函数为G(z)的***，直流增益定义为

$DCgain=\underset{z\→1}{\lim }(1-z^{-1})G\left(z\right)\frac{1}{1-z^{-1}}=\underset{z\→1}{\lim }G\left(z\right)---\left(75\right)$

获得H(s)的直流增益为

$\underset{s\→0}{\lim }H\left(s\right)=1.$

将H₂(z)的直流增益设为H(s)的直流增益，得到

K_DC＝1。

因此，传递函数H₂(z)简化为

$H_2\left(z\right)=\frac{(z-e^{-{\mathrm{\ω}}_n{T}_{sp}/2\mathrm{\ζ}})}{(z-e^{T_{sp}(-{\mathrm{\ω}}_n\mathrm{\ζ}+{\mathrm{j\ω}}_n\sqrt{1-{\mathrm{\ζ}}^2})})(z-e^{T_{sp}(-{\mathrm{\ω}}_n\mathrm{\ζ}-{\mathrm{j\ω}}_n\sqrt{1-{\mathrm{\ζ}}^2})})},---\left(76\right)$

数字环路滤波器增益

传递函数H₂(z)可进一步表达为

$H_2\left(z\right)=\frac{z-e^{-{\mathrm{\ω}}_n{T}_{sp}/2\mathrm{\ζ}}}{z^2-2e^{-2{\mathrm{\ζ\ω}}_n{T}_{sp}}cos\left({\mathrm{\ω}}_n{T}_{sp}\sqrt{1-{\mathrm{\ζ}}^2}\right)z+e^{-2{\mathrm{\ζ\ω}}_n{T}_{sp}}}.---\left(77\right)$

现在比较H(z)和H₂(z)的分母(或特征函数)，可以看出

$-K_{PD}{K}_{OSC}{K}_1+1=e^{-2{\mathrm{\ζ\ω}}_n{T}_{sp}},---\left(78\right)$

或者

$K_1=\frac{1}{K_{PD}{K}_{OSC}}\[1-e^{-2{\mathrm{\ζ\ω}}_n{T}_{sp}}\]---\left(79\right)$

以及

$K_{PD}{K}_{OSC}(K_1+K_2)-2=-2e^{-{\mathrm{\ζ\ω}}_n{T}_{sp}}cos\left({\mathrm{\ω}}_n{T}_{sp}\sqrt{1-{\mathrm{\ζ}}^2}\right)---\left(80\right)$

或者

$K_2=\frac{1}{K_{PD}{K}_{OSC}}\[1+e^{-2{\mathrm{\ζ\ω}}_n{T}_{sp}}-2e^{-{\mathrm{\ζ\ω}}_n{T}_{sp}}cos\left({\mathrm{\ω}}_n{T}_{sp}\sqrt{1-{\mathrm{\ζ}}^2}\right)\].---\left(81\right)$

通常，用于反馈控制***的性能规格通常涉及与***的时间响应关联的特定需求。建立时间t_set定义为***瞬变衰减的时间。对于PLL，t_set还指锁定时间。对于0≤ζ≤1的第二级***，针对安置在输入振幅的1％以内的***，建立时间为

$t_{set}=\frac{4.6}{{\mathrm{\ζ\ω}}_n}---\left(82\right)$

因此，对于二阶***，通过指定建立时间t_set以及阻尼因数(例如，ζ＝0.707)，无阻尼固有频率ω_n和滤波器增益K₁和K₂可从以上方程中容易确定。

电信PLL的设计人员通常采用以下逼近来确定PLL参数。PLL中的阻尼比率和电信***中的时钟通常分别具有0.1dB或者0.2dB的增益峰值(分别对应于阻尼比率4.3188和2.9585)。这得到电信同步过阻尼***。在二阶PLL中，环带宽B_w、阻尼比率ζ和固有频率ω_n通过下式关联

$B_w=\frac{{\mathrm{\ω}}_n}{2\mathrm{\π}}\sqrt{2{\mathrm{\ζ}}^2+1+\sqrt{{(2{\mathrm{\ζ}}^2+1)}^2+1}}.---\left(83\right)$

根据电信产业实际情况，具有PI滤波器的二阶PLL实现为具有B_w＝1Hz或者更小的闭环带宽和ζ≥3的阻尼比率。因此，对于本文的二阶PLL，通过指定B_w和ζ，可从以上方程确定固有频率ω_n和滤波器增益K₁和K₂。

通过利用环路滤波器增益K₁和K₂和误差e(n)，控制方程变为滤波后的误差随后可利用图19描绘的映射函数映射至对应的DDS输入控制字。简单的多映射函数是DDS校正因数如下计算

$\mathrm{\Δ}\mathrm{\φ}\left(n\right)=\tilde{e}\left(n\right)\·M,---\left(84\right)$

并且DDS控制字从下式得到

φ(n)＝φ_nom+Δφ(n)(85)

除描述的结构组件和用户交互以外，以上实施方式的***和方法可在计算机***中(具体地说，在计算机硬件或计算机软件中)实现。

术语“计算机***”包括根据上述实施方式实现***或执行方法的硬件、软件和数据存储装置。例如，计算机***可包括中央处理单元(CPU)、输入装置、输出装置和数据存储。优选地，计算机***具有用于提供视觉输出显示的监视器(例如，商业处理设计)。数据存储可包括RAM、盘驱动器或其它计算机可读介质。计算机***可包括通过网络连接的多个计算装置，并且能够经过该网络彼此通信。

以上实施方式的方法可设为计算机程序或者携带计算机程序的计算机程序产品或者计算机可读介质，所述计算机程序设置成在计算机上运行时执行上述方法。

术语“计算机可读介质”包括(但不限于)可直接由计算机或者计算机***读取或者访问的任何非临时媒体或介质。所述介质可包括(但不限于)诸如软盘、硬盘存储媒体和磁带的磁存储介质；诸如光盘或者CD-ROM的光学存储介质；诸如存储器的电存储介质，包括RAM、ROM和闪速存储器；以及以上的混合和组合，诸如磁/光学存储介质。

虽然已经结合上述示例性实施方式描述了本发明，但是当提供本公开时，许多等同修改形式和变形形式将对于本领域技术人员是清楚的。因此，认为上述本发明的示例性实施方式是示出性而非限制性的。在不脱离本发明的精神和范围的情况下，可对描述的实施方式进行各种改变。

具体地说，虽然以上实施方式的方法描述为在所述实施方式的***上实现，但是本发明的方法和***不需要彼此结合实现，而是可分别在替代性***上实现或者利用替代性方法实现。

参考文献

[1].M.Ouellette,KuiwenJi,SongLiu,andHanLi,“UsingIEEE1588andboundaryclocksforclocksynchronizationintelecomnetworks,”IEEECommun.Mag.,Feb.2011,pp.164-171.

[2].M.Ouellette，G.Garner，andS.Jobert，“SimulationsofachainofTelecomBoundaryClockscombinedwithSynchronousEthernetforphase/timetransfer,”2011InternationalIEEESymposiumonPrecisionClockSynchronizationforMeasurementControlandCommunication(ISPCS),12-16Sept.2011,pp.105-113.

[3].R.E.Kalman，“ANewApproachtoLinearFilteringandPredictionProblems,”TransactionoftheASME-JournalofBasicEngineering,March1960,pp.35-45.

[4].B.Chun，B.Kim；Y.H.Lee，“GeneralizationofexponentiallyweightedRLSalgorithmbasedonastate-spacemodel,”Proc.ofthe1998IEEEInter.Symp.onCircuitsandSystems,Vol.5,1998,pp.198-201.

以上全部参考文献以引用方式并入本文中。

Claims (21)

1.一种从装置，其通过网络连接至具有主时钟的主装置，其中，该从装置包括：

从时钟；以及

数字锁相环，该数字锁相环包括相位检测器、环路滤波器、相位累加器和计数器，其中：

所述从装置被设置成：

与所述主装置交换定时消息并且记录如下时间戳：根据所述主时钟从所述主装置发送所述定时消息的时间、根据所述从时钟接收所述定时消息的时间、根据所述从时钟发送所述定时消息的时间；以及根据所述主时钟接收所述定时消息的时间，

根据所述时间戳估计所述从时钟相对于所述主时钟的偏斜和偏移；以及

基于所估计的偏斜和偏移，使所述从时钟同步到所述主时钟，以产生主时间估计；

所述数字锁相环对所述主时间估计进行如下处理：

在接收到对主时间的第一估计时，所述计数器被初始化；

所述相位检测器被设置成在接收对主时间的随后估计时，检测所述计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；

由所述环路滤波器对所述误差信号进行滤波，以产生滤波后的误差信号；

使用滤波后的误差信号来控制所述相位累加器的频率；以及

所述相位累加器的输出使所述计数器计数递增并且还提供同步到所述主时钟的频率的所述从时钟的时钟频率。

2.根据权利要求1所述的从装置，其中，所述从装置利用所述计数器的输出作为所述从时钟的时钟时间，所述从时钟的所述时钟时间同步到所述主时钟的时间。

3.根据权利要求1或2所述的从装置，所述从装置还包括直接数字合成器，该直接数字合成器根据所述相位累加器的输出产生模拟频率信号，所述直接数字合成器包括：

相位累加器；

振荡器；

映射装置；以及

数模转换器。

4.根据权利要求3所述的从装置，所述从装置还包括低通滤波器，该低通滤波器被设置成对所述直接数字合成器的输出进行滤波。

5.根据权利要求1至4中任一项所述的从装置，其中，所述数字锁相环的所述计数器还用于提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。

6.根据权利要求5所述的从装置，其中，在所述从装置从所述主装置接收到第一定时消息时计数器被初始化，并且在接收到第一主时间估计时，所述计数器被复位至所述第一主时间估计。

7.根据权利要求1至4中任一项所述的从装置，所述从装置还包括自由运行的第二计数器，其中，该第二计数器用于提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。

8.一种使从装置中的从时钟的时间和频率同步到通过网络连接至该从装置的主装置中的主时钟的方法，该方法包括以下步骤：

在所述主装置与所述从装置之间交换定时消息和以下时间戳：根据所述主时钟从所述主装置发送所述定时消息的时间、根据所述从时钟接收所述定时消息的时间、根据所述从时钟发送所述定时消息的时间、以及根据所述主时钟接收所述定时消息的时间，

根据所述时间戳估计所述从时钟相对于所述主时钟的偏斜和偏移；以及

基于所估计的偏斜和偏移使所述从时钟同步到所述主时钟，以产生主时间估计；

使用包括相位检测器、环路滤波器、相位累加器和计数器在内的数字锁相环对上述主时间估计进行如下处理：

在接收到对主时间的第一估计时，初始化所述计数器；

在接收到对主时间的随后估计时，使用所述相位检测器来检测所述计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；

使用所述环路滤波器对所述误差信号进行滤波，以产生滤波后的误差信号；

使用滤波后的误差信号控制所述相位累加器的频率；以及

使用所述相位累加器的输出使计数器计数递增，并且

获得同步到所述主时钟的频率的所述从时钟的时钟频率作为所述相位累加器的输出。

9.根据权利要求8所述的方法，所述方法还包括以下步骤：使用所述计数器的输出作为同步到所述主时钟的时间的所述从时钟的所述时钟时间。

10.根据权利要求8或9所述的方法，所述方法还包括以下步骤：使用直接数字合成器根据所述相位累加器的输出产生模拟频率信号，该步骤包括以下步骤：

映射所述相位累加器的输出以产生数字波形；以及

使用数模转换器将所述数字波形转换为模拟波形。

11.根据权利要求10所述的方法，所述方法还包括以下步骤：对所述模拟波形进行低通滤波以产生平滑波形。

12.根据权利要求8至11中任一项所述的方法，其中，由所述数字锁相环的所述计数器提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。

13.根据权利要求12所述的方法，所述方法还包括以下步骤：

在所述从装置从所述主装置接收到第一定时消息时，初始化所述计数器；以及

在接收到第一主时间估计时，将所述计数器复位至所述第一主时间估计。

14.根据权利要求8至11中任一项所述的方法，其中，通过自由运行的第二计数器提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。

15.一种用于网络的时间和频率同步***，该***包括：

具有主时钟的主装置；

具有从时钟的从装置；以及

连接所述主装置和所述从装置的网络，其中：

所述从时钟包括：

从时钟；和

数字锁相环，该数字锁相环包括相位检测器、环路滤波器、相位累加器和计数器，其中：

所述从装置被设置成：

与所述主装置交换定时消息并且记录如下时间戳：根据所述主时钟从主装置发送所述定时消息的时间、根据所述从时钟接收所述定时消息的时间；根据所述从时钟发送所述定时消息的时间、以及根据所述主时钟接收所述定时消息的时间，

根据所述时间戳估计所述从时钟相对于所述主时钟的偏斜和偏移；以及

基于所估计的偏斜和偏移使所述从时钟同步到所述主时钟，以产生主时间估计；

所述数字锁相环对所述主时间估计进行如下处理：

在接收到对主时间的第一估计时，所述计数器被初始化；

所述相位检测器被设置成在接收到对主时间的随后估计时，检测计数器的输出与接收到的估计之间的相位差，并且产生表示该差的误差信号；

由环路滤波器对所述误差信号进行滤波，以产生滤波后的误差信号；

使用滤波后的误差信号来控制所述相位累加器的频率；以及

所述相位累加器的输出使所述计数器计数递增并且还提供同步到所述主时钟的频率的所述从时钟的时钟频率。

16.根据权利要求15所述的***，其中，所述从装置使用所述计数器的输出作为同步到所述主时钟的时间同步的所述从时钟的时钟时间。

17.根据权利要求15或16所述的***，其中，所述从装置还包括直接数字合成器，该直接数字合成器根据所述相位累加器的输出产生模拟频率信号，所述直接数字合成器包括：

相位累加器；

振荡器；

映射装置；以及

数模转换器。

18.根据权利要求17所述的***，其中，所述从装置还包括低通滤波器，该低通滤波器被设置成对所述直接数字合成器的输出进行滤波。

19.根据权利要求15至18中任一项所述的***，其中，所述数字锁相环的所述计数器还用于提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。

20.根据权利要求19所述的***，其中，在所述从装置从所述主装置接收到第一定时消息时所述计数器被初始化，在接收到第一主时间估计时，将所述计数器复位至所述第一主时间估计。

21.根据权利要求15至18中任一项所述的***，其中，所述从装置还包括自由运行的第二计数器，其中，该第二计数器用于提供针对在所述从装置接收定时消息的时间的时间戳和从所述从装置发送定时消息的时间的时间戳。