Front. Phys. Frontiers in Physics Front. Phys. 2296-424X Frontiers Media S.A. 867861 10.3389/fphy.2022.867861 Physics Original Research Two Novel One-Way Delay Clock Skew Estimators and Their Performances for the Fractional Gaussian Noise/Generalized Fractional Gaussian Noise Environment Applicable for the IEEE 1588v2 (PTP) Case Avraham and Pinchas OWD Clock Skew Estimator - PTP Avraham Yehonatan Pinchas Monika * Department of Electrical and Electronic Engineering, Ariel University, Ariel, Israel

Edited by: Ming Li, Zhejiang University, China

Reviewed by: Hu Shaoxiang, University of Electronic Science and Technology of China, China

Nan Mu, Michigan Technological University, United States

 *Correspondence: Monika Pinchas, monika.pinchas@gmail.com

This article was submitted to Interdisciplinary Physics, a section of the journal Frontiers in Physics

 22 03 2022 2022 10 867861 01 02 2022 14 02 2022 Copyright © 2022 Avraham and Pinchas. 2022 Avraham and Pinchas

This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

Papers in the literature dealing with the Ethernet network characterize packet delay variation (PDV) as a long-range dependence (LRD) process. The fractional Gaussian noise (fGn) or the generalized fractional Gaussian noise (gfGn) belong to the LRD process. The IEEE1588v2 is a two-way delay (TWD) protocol that uses the messages from the Forward (Master to Slave) and the Reverse (Slave to Master) paths. Suppose we have a significant difference between the PDV variances of the Forward and the Reverse paths. Thus, if we can use only the path with the lowest PDV variance (namely, only the one-way delay (OWD) technique), we might get a better clock skew performance from the mean square error (MSE) point of view compared with the traditional TWD method. This paper proposes two OWD clock skew estimators, one for the Forward path and one for the Reverse path applicable for the white-Gaussian, fGn and gfGn environment. Those OWD estimators do not depend on the unknown asymmetry between the fixed delays in the Forward and Reverse paths and nor on the clock offset between the Master and Slave. We also supply two closed-form approximated expressions for the MSE related to our new proposed OWD clock skew estimators. In addition, we supply some conditions, summarized in a table, guiding us whether we should use the OWD clock skew estimator for the Forward path or for the Reverse path, or just use the TWD algorithm. Simulation results confirm that our new proposed OWD clock skew estimators achieve better clock skew performances from the MSE point of view, compared with the TWD clock skew estimator recently proposed by the same authors and compared with two literature known OWD methods (the maximum likelihood and Kalman clock skew estimators).

 PTP PDV LRD fGn gfGn TWD OWD IEEE1588V2 1 Introduction

The Precision Time Protocol (PTP), named also as the IEEE 1588v2 standard [1] is a TWD exchange scheme where the Slave exchanges a series of synchronization packets with its Master so the packet timestamps can be employed to estimate the clock skew relative to the Master. In other words, the PTP as a TWD algorithm uses messages received from the Forward (Master to Slave) and from the Reverse (Slave to Master) paths in order to estimate the clock skew relative to the Master. The synchronization packets can encounter several intermediate switches and routers along the network path between the Master and the Slave [2]. Networks often suffer large unpredictable queuing delays at switches and routers (thus having heavy PDV in the network) due to the presence of background traffic [3]. This background traffic may be a real traffic or one caused by a cyber attack [2] where a malicious intermediate node deliberately delays the transmission of synchronization messages. A heavy PDV can be seen in the Forward path, in the Reverse path or in both paths. Usually, the PDV in the Forward path is different from the PDV in the Reverse path. The problem is that the PDV of the Forward and Reverse paths can significantly hamper the accuracy of the clock skew estimation [3]. A lower PDV will lead to a more accurate clock skew estimation compared to a higher PDV. Thus, if the difference in the PDVs encountered in Forward and Reverse paths is high, the clock skew estimation accuracy obtained with the TWD technique is mainly decreased due to the path with the higher PDV. Thus, if we could use for the clock skew estimation task only the path with the lowest PDV (namely, using the OWD technique), the clock skew performance from the MSE point of view might be improved compared with the case where we have also to consider the path with the higher PDV (the TWD approach). Since the lower PDV path may occur in the Forward path as well as in the Reverse path, two different OWD clock skew estimators are needed (one for the Forward path and the other one associated with the Reverse path). So far we have seen that for the clock skew estimation task, the OWD technique may be more useful compared with the TWD approach in cases where the Forward PDV variance is very different from the Reverse PDV variance. According to [4], PTP (which is a TWD exchange scheme) has more unknown parameters than available equations. Thus, in order to solve the problem, a symmetric path is usually assumed. Namely, the fixed delay in the Forward path is usually assumed to be the same as the fixed delay for the Reverse path. But, in practical scenarios, this is not the case. Thus, for an asymmetrical path, when the symmetric path assumption is applied, a degradation in the clock skew estimation may be obtained when the TWD approach is applied. Now, the OWD technique relies only on one path, on the Reverse path or on the Forward path. Thus, the symmetric assumption is not needed in the OWD technique which can be considered here as an advantage compared with the TWD approach. Suppose for a moment that we have three clock skew estimators applicable for the PTP case. Namely, we have one OWD clock skew estimator for the Forward path, one OWD clock skew estimator for the Reverse path and a TWD clock skew estimator. Next we wish to know which of the three clock skew estimators should be taken for the clock skew estimation task given a network where different PDV variances are seen on both Forward and Reverse paths but the fixed delay of the Forward path is equal to the fixed delay of the Reverse path. It is quite reasonable to think that when the Forward path PDV variance is equal or close to equal to the Reverse path PDV variance, the TWD clock skew estimator is preferable over the OWD clock skew estimator due to the “averaging” effect of the variances in the TWD clock skew estimator. But, when the PDV variances of the Forward and Reverse paths are different and on the same time the difference in the variances is not very high, it is not clear if the OWD clock skew estimator for the Forward path or the OWD clock skew estimator for the Reverse path or maybe the TWD clock skew estimator should be applied for the clock skew estimation task. Thus, some guiding lines (closed-form expressions, conditions) are needed here, telling us which approach should be applied in order to get the best clock skew performance in the MSE point of view. Namely, which clock skew estimator should be taken: the OWD clock skew estimator for the Forward path or the OWD clock skew estimator for the Reverse path or perhaps the TWD clock skew estimation approach. Recently [4], we proposed a new TWD clock skew estimator for the PTP case that has the best clock skew performance in the MSE point of view compared to the relevant literature known estimators [5–7]. This clock skew estimator [4] is suitable for the white-Gaussian and fGn/gfGn cases and does not depend on the asymmetric fixed delay between the Forward and the Reverse paths, nor on the offset between the Master and the Slave clocks. This paper is a direct continuation of our previous work [4]. Thus, please refer to [4] in order to find a detailed overview of the existing TWD and recently proposed OWD approaches for the PTP case. Please note that the two recently proposed OWD clock skew estimators [6,7], are both OWD clock skew estimators associated with the Forward path. Thus, if the Forward path PDV variance is much higher compared with the Reverse path PDV variance, [6,7] may not get better clock skew performance from the MSE point of view compared with the TWD approach and compared with the OWD clock skew estimator associated with the Reverse path. As already was mentioned, this paper is a direct continuation of our previous work [4] where we proposed a novel TWD clock skew estimator applicable for the PTP case. In this paper we propose:

1. A novel OWD clock skew estimator for the Forward path based on [4], applicable for the white-Gaussian and fGn/gfGn environment.

2. A novel OWD clock skew estimator for the Reverse path based on [4], applicable for the white-Gaussian and fGn/gfGn environment.

3. A closed-form-approximated expression for the clock skew performance (MSE) related to our OWD proposed clock skew estimator for the Forward path.

4. A closed-form-approximated expression for the clock skew performance (MSE) related to our OWD proposed clock skew estimator for the Reverse path.

5. Guiding lines (closed-form expressions, conditions), summarized in a table (please refer to Table 1), telling us if we should use the OWD clock skew estimator for the Forward path or the OWD clock skew estimator for the Reverse path or perhaps the TWD clock skew estimator proposed by [4] in order to get the best clock skew performance from the MSE point of view.

Summary of the conditions where the suggested estimator have the possible lower MSE.

Z σ ω 1 2 = σ ω 2 2
Z = 1 ( σ ω 1 2 = σ ω 2 2 = σ 2 ) Z ≠ 1
white-Gaussian if σ 2 ≥ σ ̃ 2 when Z > 1: if Z ≥ Z ̃
then we use OWD in Eq. 62 then we use OWD in Eq. 61
otherwise, we use TWD in Eq. 60 otherwise, we use TWD in Eq. 60
σ ̃ 2 = 2 A B T s y n 2 when Z < 1: if Z ≤ Z ̄
then we use OWD in Eq. 62
otherwise, we use TWD in Eq. 60
Z ̃ = 3 1 + B A σ ω 1 2 T s y n 2
Z ̄ = 1 3 1 + B A σ ω 1 2 T s y n 2
fGn/gfGn if σ 2 ≥ σ ̃ f G n / g f G n 2 when Z > 1: if Z ≥ Z ̃ f G n / g f G n
then we use OWD in Eq. 62 then we use OWD in Eq. 61
otherwise, we use TWD in Eq. 60 otherwise, we use TWD in Eq. 60
σ ̃ f G n / g f G n 2 = 2 A B C + D C T s y n 2 when Z < 1: if Z ≤ Z ̄ f G n / g f G n
then we use OWD in Eq. 62
otherwise, we use TWD in Eq. 60
Z ̃ f G n / g f G n = 3 1 + C C + D B A σ ω 1 2 T s y n 2
Z ̄ f G n / g f G n = 1 3 1 + C C + D B A σ ω 1 2 T s y n 2

The clock skew performances (MSE) of our new proposed OWD clock skew estimators were compared via simulation with the clock skew performances (MSE) obtained with two TWD clock skew estimators [4,5] and with the literature known OWD clock skew estimators [6,7]. Simulation results will show the advantage in performance (MSE) of our new proposed OWD clock skew estimators compared to [4–7]. Simulation results will also show the effectiveness of our closed-form-approximated expressions for the clock skew performance (MSE) associated to the Forward and Reverse paths as well as the effectiveness of our proposed guiding lines, leading us to the right choice of the clock skew estimator from the MSE point of view.

The paper is organized as follows. Section 2 briefly introduces the system under consideration and the assumptions we applied for our algorithm. Section 3 proposes the OWD clock skew estimators for the Forward and Reverse paths. Section 4 suggests the closed-form approximated expressions for the MSE related to our new proposed OWD clock skew estimators where the PDV is a white-Gaussian process. Section 5 suggests the closed-form approximated expressions for the MSE related to our new proposed OWD clock skew estimators where the PDV is an fGn/gfGn process. In Section 6, we derive some guiding lines (conditions), summarized in Table 1, telling us under what condition should we prefer the OWD for the Forward path over the OWD for the Reverse path or should just prefer the TWD clock skew estimator obtained in [4]. Section 7 presents simulation results, and in Section 8, a conclusion is given.

 2 System Description

As already was mentioned earlier, this paper is a direct continuation of our previous work [4]. Thus, the system description is the same as in [4]. Please refer to [4], for having a detailed description of the message exchange flow between the Master and the Slave. Let us recall Figure 1 from [4] where based on [8–10] we may write: t 1 j + d m s + ω 1 j = t 2 j 1 + α + Q t 4 j − d s m − ω 2 j = t 3 j 1 + α + Q where Q is the time difference between the Master and the Slave clocks (offset) and α is the clock skew. The Forward and the Reverse fixed delays are denoted as d ms , d sm respectively. The Forward PDV is denoted as ω 1[j] and the Reverse PDV is denoted as ω 2[j]. The total number of the Sync messages periods is denoted as J, where j = 1, 2, 3, … , J. At timestamp t 1, the Master sends a Sync message to the Slave. The Slave receives this Sync message at timestamp t 2 and sends back to the Master a Delay Req message at timestamp t 3. The Master receives this Delay Req message at timestamp t 4. Please note that t 1[j], t 2[j], t 3[j], t 4[j] are the timestamps of t 1, t 2, t 3 and t 4 respectively at the jth Sync message period.

PTP messaging timing diagram.

We consider two different models for the PDV, as was done in [4]:

1. The PDV is modeled as a white-Gaussian noise with zero mean and the variance E ω n [ j ] , ω n [ m ] is σ ω n 2 when j = m and zero when j ≠ m

where E [.] denotes the expectation operator on (.) and n = 1, 2.

2. The PDV is modeled as an fGn/gfGn process with zero mean. Based on [11–13] we have:

a. When j = m: E ω n [ j ] , ω n [ m ] = σ ω n 2 .

b. When j ≠ m: E ω n [ j ] , ω n [ m ] = σ ω n 2 2 ( j − m ) a | − 1 2 H − 2 ( | ( j − m ) a | ) 2 H + ( | ( j − m ) a | + 1 ) 2 H .

 where a = 1 is the fGn case.

In addition, we use also the same assumptions as were made in [4]:

1. The Forward and the Reverse PDVs are independent. Thus, we have: E ω 1 [ j ] , ω 2 [ m ] = 0 ∀ j , m .

2. In the Slave clock the time between t 2 [j] to t 3 [j] is constant and is denoted as X. Thus, we have: t 3 [j] − t 2 [j] = X.

In this paper we propose two novel OWD clock skew estimators, one for the Forward path and one for the Reverse path. Both OWD clock skew estimators are based on our previously TWD clock skew estimator [4] given by: α ̂ = 1 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 1 , j i T 2 , j i + T 4 , j i T 2 , j i − 1 where T 1 , j i = t 1 j + i − t 1 j , T 2 , j i = t 2 j + i − t 2 j , T 4 , j i = t 4 j + i − t 4 j

3 The OWD Clock Skew Estimators

In the following, we present in Theorem 1 our new proposed OWD clock skew estimator for the Forward path and in Theorem 2 our new proposed OWD clock skew estimator for the Reverse path.

3.1 Theorem 1

The clock skew estimator for the OWD in the Forward path (Master to Slave) can be written as: α ̂ F = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 1 , j i T 2 , j i − 1 where α ̂ F is the is the clock skew estimator that based only on timestamps from the Forward path.

Proof of Theorem 1

In order to avoid the fixed delay, we can subtract between two timestamps from different Sync periods. Therefore, based on Eq. 1 we have: T 1 , j i + Ω 1 , j i = T 2 , j i 1 + α j , i F where α j , i F is the clock skew between the (j + i)-th and ith Sync period, and Ω1,j (i) is: Ω 1 , j i = ω 1 j + i − ω 1 j

Based on Eq. 6 the clock skew can be written as: α j , i F = T 1 , j i T 2 , j i + Ω 1 , j i T 2 , j i − 1

The OWD clock skew in the Forward path can be defined as: α F = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i α j , i F

By putting Eq. 8 into Eq. 9 we define the clock skew in the Forward path as: α F = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 1 , j i T 2 , j i + Ω 1 , j i T 2 , j i − 1 = α ̂ F + 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i Ω 1 , j i T 2 , j i

This completes our proof.

 3.2 Theorem 2

The skew clock estimator for the OWD in the Reverse path (Slave to Master) can be written as: α ̂ R = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 4 , j i T 2 , j i − 1 where α ̂ R is the is the clock skew estimator that based only on timestamps from the Reverse path.

Proof of Theorem 2

In order to avoid the fixed delay, we can subtract between two timestamps from different Sync periods. Therefore, based on Eq. 2 we have: T 4 , j i − Ω 2 , j i = T 3 , j i 1 + α j , i R where α j , i R is the clock skew between the (j + i)-th and ith Sync period, and T 3 , j i = t 3 j + i − t 3 j ; Ω 2 , j i = ω 2 j + i − ω 2 j

Based on the definition that t 3 [j] − t 2 [j] = X, as mentioned in the section of System Description (assumption 2), we can write: T 3 , j i = t 3 j + i − t 3 j = t 2 j + i + X − t 2 j + X = T 2 , j i

By using Eq. 14, we can write Eq. 12 as: T 4 , j i − Ω 2 , j i = T 2 , j i 1 + α j , i R

Therefore, based on Eq. 15 the clock skew can be written as: α j , i R = T 4 , j i T 2 , j i − Ω 2 , j i T 2 , j i − 1

The OWD clock skew in the Reverse path can be defined as: α R = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i α j , i R

By putting Eq. 16 into Eq. 17 we define the clock skew in the Reverse path as: α R = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 4 , j i T 2 , j i − Ω 2 , j i T 2 , j i − 1 = α ̂ R − 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i Ω 2 , j i T 2 , j i

This completes our proof.

 4 The Clock Skew Performance for the White-Gaussian Case

In the following, Theorems 3 and 4 present the closed-form approximated expressions for the MSE related to our new proposed OWD clock skew estimator for the Forward path and Reverse path, respectively. According to [4], the MSE for the TWD clock skew estimator for the Gaussian case is given by: E e 2 ≈ 1 J J − 1 2 σ ω 1 2 + σ ω 2 2 T s y n 2 A 1 + 1 P where P is: P = A B σ ω 1 2 + σ ω 2 2 σ ω 1 4 T s y n 2 and A, B are given by: A = 2 ∑ i = 1 J − 1 J − i i 2 + ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 k ≠ i J − 1 ∑ m = 1 m = j m = j + i − k J − k 1 i k − ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 m = j + i m = j − k J − k 1 i k B = 12 ∑ i = 1 J − 1 J − i i 4 + 6 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 k ≠ i J − 1 ∑ m = 1 m = j m = j + i − k m = j + i m = j − k J − k 1 i k 2 + 4 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 m ≠ j m ≠ j + i − k m ≠ j + i m ≠ j − k J − k 1 i k 2

4.1 Theorem 3

For Ω 1 , j ( i ) T 1 , j ( i ) ≪ 1 , the closed-form-approximated expression for the MSE related to the OWD clock skew estimator (for the Forward path), can be defined as: E e F 2 ≈ 2 J J − 1 2 σ ω 1 2 T s y n 2 A 1 + 1 P F where P F is: P F = A B T s y n 2 σ ω 1 2 where A and B are defined in Eq. 21, 22 respectively.

Proof of Theorem 3

Based on Eq. 10 the error of the OWD clock skew estimator (for the Forward path) is: e F = α F − α ̂ F = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i Ω 1 , j i T 2 , j i

Now, according to Eq. 6 we can write T 2,j (i) as following: T 2 , j i = T 1 , j i + Ω 1 , j i 1 + α j , i F

Based on Eq. 26 we may write the expectation of Eq. 25 as: E e F = 2 1 + α F J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i E a j , i 1 + a j , i where a j , i = Ω 1 , j i T 1 , j i

For Ω 1 , j ( i ) T 1 , j ( i ) ≪ 1 we can write Eq. 27 as: E e F ≈ 2 1 + α F J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i E a j , i 1 − a j , i

Based on Eq. 29 the approximated MSE related to the OWD clock skew estimator (for the Forward path) can be define as: E e F 2 ≈ 4 1 + α F 2 J J − 1 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ j = m J − k E a j , i a m , k − E a j , i 2 a m , k − E a j , i a m , k 2 + E a j , i 2 a m , k 2

According to [4], we can write T 1,j (i) and T 1,m (k) as: T 1 , j i = t 1 j + i − t 1 j = i T s y n T 1 , m k = t 1 m + k − t 1 m = k T s y n where T syn is the Sync messages period.

Based on Eq. 28 and on Eq. 31 we can simplify the expressions in Eq. 30: E a j , i a m , k = E Ω 1 , j i T 1 , j i Ω 1 , m k T 1 , m k = E Ω 1 , j i Ω 1 , m k i k T s y n 2 E a j , i 2 a m , k = E Ω 1 , j 2 i T 1 , j 2 i Ω 1 , m k T 1 , m k = E Ω 1 , j 2 i Ω 1 , m k i 2 k T s y n 3 E a j , i a m , k 2 = E Ω 1 , j i T 1 , j i Ω 1 , m 2 k T 1 , m 2 k = E Ω 1 , j i Ω 1 , m 2 k i k 2 T s y n 3 E a j , i 2 a m , k 2 = E Ω 1 , j 2 i T 1 , j 2 i Ω 1 , m 2 k T 1 , m 2 k = E Ω 1 , j 2 i Ω 1 , m 2 k i 2 k 2 T s y n 4

Since the PDV has zero mean Eqs 33 and 34 can be set to zero. Now, by putting Eqs 32 and 35 into Eq. 30, we can write the following expression: E e F 2 ≈ 2 1 + α F J J − 1 T s y n 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 1 , j i Ω 1 , m k i k + 1 T s y n 2 E Ω 1 , j 2 i Ω 1 , m 2 k i k 2

Based on [4] we can write the two summation parts in Eq. 36 as: ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 1 , j i Ω 1 , m k i k = σ ω 1 2 A ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k 1 T s y n 2 E Ω 1 , j 2 i Ω 1 , m 2 k i k 2 = σ ω 1 4 T s y n 2 B where A and B are defined in Eq. 21, 22 respectively.

Based on Eq. 37 we may write (36) as: E e F 2 ≈ 2 1 + α F J J − 1 T s y n 2 σ ω 1 2 A + σ ω 1 4 T s y n 2 B

For practical systems, the two clocks (Master and Slave) operate at almost the same frequency. Therefore we can write that (1 + α F ) ≈ 1.

After rearranging Eq. 38 we can write: E e F 2 ≈ 2 J J − 1 T s y n 2 σ ω 1 2 A 1 + σ ω 1 2 T s y n 2 B A

Now, it can be easily seen that based on Eq. 39 we can write Eq. 23, and this completes our proof.

 4.2 Theorem 4

For Ω 1 , j ( i ) T 1 , j ( i ) ≪ 1 , the closed-form-approximated expression for the MSE related to the OWD clock skew estimator (for the Reverse path), can be defined as: E e R 2 ≈ 2 J J − 1 2 σ ω 2 2 T s y n 2 A where A is defined in Eq. 21.

Proof of Theorem 4

Based on Eq. 18 the error of the OWD clock skew estimator (for the Reverse path) is: e R = α R − α ̂ R = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i − Ω 2 , j i T 2 , j i

Let us recall (26): T 2 , j i = T 1 , j i + Ω 1 , j i 1 + α j , i F

Based on Eq. 26, we may write the expectation of Eq. 41 as: E e R ≈ 2 1 + α F J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i − E b j , i 1 + a j , i where a j , i = Ω 1 , j i T 1 , j i b j , i = Ω 2 , j i T 1 , j i

As mentioned before for practical systems, the two clocks (Master and Slave) operate at almost the same frequency. Therefore, we can write (1 + α F ) ≈ 1.

For Ω 1 , j ( i ) T 1 , j ( i ) ≪ 1 we can write Eq. 42 as: E e R ≈ 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i − E b j , i − a j , i b j , i

Based on the assumption of independence of the Forward and the Reverse messages (as mentioned in assumption 1 in section 2), we may write Eq. 44 as: E e R ≈ 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i − E b j , i

Based on Eq. 45 the closed-form-approximated expression for the MSE related to the OWD clock skew estimator (for the Reverse path) can be define as: E e R 2 ≈ 4 J J − 1 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ j = m J − k E b j , i b m , k

Based on Eq. 43 and on Eq. 31 we can rewrite the expression E b j , i b m , k as: E b j , i b m , k = E Ω 2 , j i T 1 , j i Ω 2 , m k T 1 , m k = E Ω 2 , j i Ω 2 , m k i k T s y n 2

Now, based on Eq. 47 the closed-form-approximated expression for the MSE related to the OWD clock skew estimator (for the Reverse path) is: E e R 2 ≈ 2 J J − 1 T s y n 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 2 , j i Ω 2 , m k i k

According to [4], we can write the summation part in Eq. 48 as: ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 2 , j i Ω 2 , m k i k = σ ω 2 2 A where A is defined in Eq. 21.

Now, it can be easily seen that by putting Eq. 49 into Eq. 48 we get Eq. 40, and this completes our proof.

 5 The Clock Skew Performance for the LRD Case

In this section we applied the fGn/gfGn model for the LRD process. This model has the Hurst exponent in the range of 0.5 ≤ H < 1 and the a parameter in the range of 0 < a ≤ 1, where for a = 1 we have the fGn case. In the following, Theorems 5 and 6 present the closed-form approximated expressions for the MSE related to our new proposed OWD clock skew estimator for the Forward path and Reverse path, respectively. According to [4] the closed-form-approximated expression for the MSE related to the TWD clock skew estimator [4] is given by: E e 2 ≈ 1 J J − 1 2 σ ω 1 2 + σ ω 2 2 T s y n 2 1 + 1 P C + D where C and D are given by: C = ∑ i = 1 J − 1 J − i i 2 2 − f G H i , H + ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 k ≠ i J − 1 ∑ m = 1 m = j m = j + i − k J − k 1 i k 1 + 1 2 f G H ∗ i − k , H , a − f G H ∗ i , H , a − f G H ∗ k , H , a − ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 m = j + i m = j − k J − k 1 i k 1 − 1 2 f G H ∗ i , H , a − f G H ∗ k , H , a + f G H ∗ i + k , H , a D = ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 m ≠ j m ≠ j + i m ≠ j − k m ≠ j + i − k J − k 1 2 i k f G H ∗ j − m , H , a − f G H ∗ j + i − m , H , a − f G H ∗ j − m − k , H , a + f G H ∗ j + i − m − k , H , a the function f G H ∗ ( . ) is: f G H ∗ x , H , a = ‖ x a | − 1 | 2 H − 2 | x a | 2 H + | x a | + 1 2 H and P is defined in Eq. 20.

5.1 Theorem 5

The closed-form approximated expression for the MSE related to our new proposed OWD clock skew estimator (for the Forward path) can be defined as: E e F 2 ≈ 2 J J − 1 2 σ ω 1 2 T s y n 2 1 + 1 P F C + D where C, D and P F are defined in Eq. 51, 52 and 24 respectively.

Proof of Theorem 5

The MSE for the fGn/gfGn case is based on the MSE of the OWD clock skew estimator for the Forward path defined in Eq. 36. Based on the fact that 1 + α F ≈ 1, we can write Eq. 36 as: E e F 2 ≈ 2 J J − 1 T s y n 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 1 , j i Ω 1 , m k i k + 1 T s y n 2 E Ω 1 , j 2 i Ω 1 , m 2 k i k 2

According to [4], we can write the first part in Eq. 55 as: ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 1 , j i Ω 1 , m k i k = σ ω 1 2 C + D where C and D are defined in Eq. 51, 52 respectively.

The calculation of the second expression in Eq. 55 is quite difficult to carry out for the fGn/gfGn case. Therefore, following [4] we can write: ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 1 , j i Ω 1 , m k i k + 1 T s y n 2 E Ω 1 , j 2 i Ω 1 , m 2 k i k 2 = σ ω 1 2 C 1 + 1 P F + D

Now, by putting Eq. 57 into Eq. 55, we get Eq. 54 and this completes our proof.

 5.2 Theorem 6

The closed-form approximated expression for the MSE related to our new proposed OWD clock skew estimator (for the Reverse path) can be defined as: E e R 2 ≈ 2 J J − 1 2 σ ω 2 2 T s y n 2 C + D where C and D are defined in Eq. 51, 52.

Proof of Theorem 6

The MSE for the fGn/gfGn case is based on the MSE of the OWD clock skew estimator for the Reverse path defined in Eq. 48. Let as recall Eq. 48: E e R 2 ≈ 2 J J − 1 T s y n 2 ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 2 , j i Ω 2 , m k i k

We can write the summation part in Eq. 48 as was done in [4] as: ∑ i = 1 J − 1 ∑ j = 1 J − i ∑ k = 1 J − 1 ∑ m = 1 J − k E Ω 2 , j i Ω 2 , m k i k = σ ω 2 2 C + D

Now, by putting Eq. 59 into Eq. 48, we get Eq. 58 and this completes our proof.

 6 The Preferred Clock Skew Estimator for Each Scenario

In this paper we proposed two OWD clock skew estimators (Eqs 5, 11). Thus, we can consider now two OWD clock skew estimators and one TWD clock estimator proposed by [4]. Let us recall the three estimators.

At first, we recall the TWD clock skew estimator from [4]: α ̂ = 1 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 1 , j i T 2 , j i + T 4 , j i T 2 , j i − 1

The OWD clock skew estimator for the Forward path is given by Eq. 5: α ̂ F = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 1 , j i T 2 , j i − 1

The OWD clock skew estimator for the Reverse path is given by Eq. 11: α ̂ R = 2 J J − 1 ∑ i = 1 J − 1 ∑ j = 1 J − i T 4 , j i T 2 , j i − 1

It would be very helpful for the system designer if he could know which of the above listed clock skew estimators (Eqs 60–62) he should choose in order to achieve the best clock skew performance in the MSE point of view. In this section we will give the system designer guidelines (conditions) that will help him to choose wisely the best clock skew estimator in order to achieve the best clock skew performance from the MSE point of view.

Please note, in this section we define Z as: Z = σ ω 2 2 σ ω 1 2

Thus we have:

a. when Z > 1 → σ ω 1 2 < σ ω 2 2 ,

b. when Z < 1 → σ ω 1 2 > σ ω 2 2 ,

c. when Z = 1 → σ ω 1 2 = σ ω 2 2 .

In the following we have Theorem seven supplying us guidelines for choosing the preferable clock skew estimator from the above listed clock skew estimators (Eqs 60–62) leading to the best clock skew performance from the MSE point of view for the white-Gaussian case. Theorem 8 supplies us guidelines for choosing the preferable clock skew estimator from the above listed clock skew estimators (Eqs 60–62) leading to the best clock skew performance from the MSE point of view for the fGn/gfGn case.

6.1 Theorem 7

For the white-Gaussian process we have the following conditions:

Case a: For Z > 1:

If Z ≥ Z ̃ , then we use the OWD clock skew estimator in the Forward path (Eq. 61).

where Z ̃ = 3 1 + B A σ ω 1 2 T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Case b: For Z < 1:

If Z ≤ Z ̄ , then we use the OWD clock skew estimator in the Reverse path (Eq. 62),

where Z ̄ = 1 3 1 + B A σ ω 1 2 T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Case c: For Z = 1 ( σ ω 1 2 = σ ω 2 2 = σ 2 ) :

If σ 2 ≥ σ ̃ 2 , then we use the OWD clock skew estimator in the Reverse path (Eq. 62),

where σ ̃ 2 = 2 A B T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Proof of Theorem 7

We rewrite the closed-form-approximated expression for the MSE related to the TWD clock skew estimator (Eq. 19), and the closed-form-approximated expressions for the MSE related to our new proposed OWD clock skew estimators (Eqs 23, 40) with the help of Z (Eq. 63).

The closed-form approximated expressions for the MSE related to the TWD clock skew estimator (Eq. 19) can be written as: E e 2 ≈ 1 J J − 1 2 σ ω 1 2 1 + Z T s y n 2 A 1 + 1 P

The closed-form approximated expression for the MSE related to the OWD clock skew estimator for the Forward path (Eq. 23) can be written as: E e F 2 ≈ 2 J J − 1 2 σ ω 1 2 T s y n 2 A 1 + 1 P F

The closed-form approximated expression for the MSE related to the OWD clock skew estimator for the Reverse path (Eq. 40) can be written as: E e R 2 ≈ 2 J J − 1 2 Z σ ω 1 2 T s y n 2 A where A, B, P and P F are defined in Eqs 21, 22, 20 and 24 respectively.

In the following, we define M S E G T , M S E G F and M S E G R as the MSE of the TWD case (Eq. 67), OWD in the Forward path case (Eq. 68) and OWD in the Reverse path case (Eq. 69) respectively.

For Z > 1, we wish to find the value for Z where we have M S E G T ≥ M S E G F . Thus, we may write with the help of Eqs 20, 24: 1 J J − 1 T s y n 2 σ ω 1 2 1 + Z A + B σ ω 1 4 T s y n 2 ≥ 1 J J − 1 T s y n 2 4 σ ω 1 2 A + B σ ω 1 4 T s y n 2 where Eq. 70 can be written also as: σ ω 1 2 1 + Z A + B σ ω 1 4 T s y n 2 ≥ 4 σ ω 1 2 A + B σ ω 1 4 T s y n 2

After rearranging Eq. 71 we can write: Z σ ω 1 2 A ≥ 3 σ ω 1 2 A + B σ ω 1 4 T s y n 2

We can divide Eq. 72 by A σ ω 1 2 ( A σ ω 1 2 > 0 ) and then we can write: Z ≥ Z ̃ where Z ̃ is defined in Eq. 64.

This completes our proof of Theorem 7, case a.

Now, we continue the proof of case b. For Z < 1, we wish to find the value for Z where we have M S E G T ≥ M S E G R . Thus, with the help of Eq. 20 we may write: 1 J J − 1 T s y n 2 σ ω 1 2 1 + Z A + B σ ω 1 4 T s y n 2 ≥ 1 J J − 1 T s y n 2 4 Z σ ω 1 2 A where Eq. 74 can be written also as: σ ω 1 2 1 + Z A + B σ ω 1 4 T s y n 2 ≥ 4 Z σ ω 1 2 A

After rearranging Eq. 75 we can write: σ ω 1 2 A + B σ ω 1 2 T s y n 2 ≥ 3 Z σ ω 1 2 A

We can divide Eq. 76 by 3 A σ ω 1 2 ( 3 A σ ω 1 2 > 0 ) and then we can write: Z ≤ Z ̄ where Z ̄ is defined in Eq. 65.

This completes our proof of Theorem 7, case b.

Now we can continue the proof of case c of Theorem 7. For Z = 1, we defined that σ ω 1 2 = σ ω 2 2 = σ 2 . We wish to find the value for Z where we have M S E G T ≥ M S E G R (please note that for Z = 1: M S E G T < M S E G F ∀ j , since B has only positive values). Thus, with the help of Eq. 20 we can write: 1 J J − 1 T s y n 2 2 σ 2 A + B σ 4 T s y n 2 ≥ 1 J J − 1 T s y n 2 4 σ 2 A where Eq. 78 can be written also as: 2 σ 2 A + B σ 4 T s y n 2 ≥ 4 σ 2 A

After rearranging Eq. 79 we can write: σ 2 ≥ σ ̃ 2 where σ ̃ 2 is defined in Eq. 66. Now, we have completed our proof of Theorem 7.

 6.2 Theorem 8

For the fGn/gfGn process we have the following conditions:

Case a: For Z > 1:

If Z ≥ Z ̃ f G n / g f G n , then we use the OWD clock skew estimator in the Forward path (Eq. 61),

where Z ̃ f G n / g f G n = 3 1 + C C + D B A σ ω 1 2 T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Case b: For Z < 1:

If Z ≤ Z ̄ f G n / g f G n , then we use the OWD clock skew estimator in the Reverse path (Eq. 62),

where Z ̄ f G n / g f G n = 1 3 1 + C C + D B A σ ω 1 2 T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Case c: For Z = 1 ( σ ω 1 2 = σ ω 2 2 = σ 2 ) :

If σ 2 ≥ σ ̃ f G n / g f G n 2 , then we use the OWD clock skew estimator in the Reverse path (Eq. 62),

where σ ̃ f G n / g f G n 2 = 2 C + D C A B T s y n 2

Otherwise, we use the TWD clock skew estimator (Eq. 60).

Proof of Theorem 8

We rewrite the closed-form-approximated expression for the MSE related to the TWD clock skew estimator (Eq. 50), and the closed-form-approximated expressions for the MSE related to our new proposed OWD clock skew estimators (Eqs 54, 58 with the help of Z (Eq. 63).

The closed-form approximated expressions for the MSE related to the TWD clock skew estimator (Eq. 50) can be written as: E e 2 ≈ 1 J J − 1 2 σ ω 1 2 1 + Z T s y n 2 C 1 + 1 P + D

The closed-form approximated expressions for the MSE related to the OWD clock skew estimator for the Forward path Eq. 54 can be written as: E e F 2 ≈ 2 J J − 1 2 σ ω 1 2 T s y n 2 C 1 + 1 P F + D

The closed-form approximated expressions for the MSE related to the OWD clock skew estimator for the Reverse path Eq. 58 can be written as: E e R 2 ≈ 2 J J − 1 2 Z σ ω 1 2 T s y n 2 C + D where C, D, P and P F are defined in Eqs 51, 52, 20 and 24 respectively.

In the following, we define: M S E f G n / g f G n T , M S E f G n / g f G n F and M S E f G n / g f G n R as the MSE of the TWD case (Eq. 84), OWD in the Forward path case (Eq. 85) and OWD in the Reverse path case (Eq. 86) respectively.

For Z > 1, we wish to find the value for Z where we have M S E f G n / g f G n T ≥ M S E f G n / g f G n F . Thus, we may write: 1 J J − 1 T s y n 2 σ ω 1 2 1 + Z 1 + 1 P C + D ≥ 1 J J − 1 T s y n 2 4 σ ω 1 2 1 + 1 P F C + D where Eq. 87 can be written also as: 1 + Z 1 + 1 P C + D ≥ 4 1 + 1 P F C + D

Based on the definition of P and P F in Eqs 20, 24 respectively, we can write: P = A B 1 + Z T s y n 2 σ ω 1 2 , P F = A B T s y n 2 σ ω 1 2 → P = P F 1 + Z

Based on Eq. 89 we can write Eq. 88 as: 1 + Z 1 + 1 P F 1 + Z C + D ≥ 4 1 + 1 P F C + D

After rearranging Eq. 90 we can write: Z ≥ 3 1 + 1 P F C + D C + D

By putting the definition of P F Eq. 24 into Eq. 91 we can write: Z ≥ 3 C + D + B C A σ ω 1 2 T s y n 2 C + D → Z ≥ Z ̃ f G n / g f G n where Z ̃ f G n / g f G n is defined in Eq. 81.

This completes our proof of Theorem 8, case a.

Now we continue the proof of case b. For Z < 1, we wish to find the value for Z where we have M S E f G n / g f G n T ≥ M S E f G n / g f G n R . Thus, we may write: 1 J J − 1 T s y n 2 σ ω 1 2 1 + Z 1 + 1 P C + D ≤ 1 J J − 1 T s y n 2 4 Z σ ω 1 2 C + D where Eq. 93 can be written also as: 1 + Z 1 + 1 P C + D ≤ 4 Z C + D

By using the definition of P in Eq. 20 into Eq. 95 we may write: 1 + Z 1 + B A σ ω 1 2 1 + Z T s y n 2 C + D ≤ 4 Z C + D

After rearranging Eq. 95 we can write: Z ≤ 1 3 C + D + C B A σ ω 1 2 T s y n 2 C + D → Z ≤ Z ̄ f G n / g f G n where Z ̄ f G n / g f G n is defined in Eq. 82.

This completes our proof of Theorem 8, case b.

Now we can continue the proof of case c of Theorem 8. For Z = 1, we defined that σ ω 1 2 = σ ω 2 2 = σ 2 . We wish to find the value for Z where we have M S E f G n / g f G n T ≥ M S E f G n / g f G n R . Note that for Z = 1 M S E f G n / g f G n T < M S E f G n / g f G n F ∀ j since, A, B and C have only positive values. Thus, we can write: 1 J J − 1 T s y n 2 2 σ 2 C 1 + 1 P + D ≥ 1 J J − 1 T s y n 2 4 σ 2 C + D where Eq. 97 can be written also as: 2 σ 2 C 1 + 1 P + D ≥ 4 σ 2 C + D

By using the definition of P in Eq. 20, we may write Eq. 98 as: 2 σ 2 C 1 + B A σ 2 2 T s y n 2 + D ≥ 4 σ 2 C + D

After rearranging Eq. 99 we can write: σ 2 ≥ σ ̃ f G n / g f G n 2 where σ ̃ f G n / g f G n 2 is defined in Eq. 83. Now, we have completed our proof of Theorem 8.

 7 Simulation Results

In this section, we first start to test our guidelines (conditions) from Theorem 8, summarized in Table 1. Figures 2, 3 show the simulated clock skew performance (MSE) comparison between the OWD clock skew estimator for the Forward path (Eq. 61) with the TWD clock skew estimator proposed by Avraham and Pinchas [4] for the fGn case. Figures 4–6 show the simulated clock skew performance (MSE) comparison between the OWD clock skew estimator for the Reverse path (Eq. 62) with the TWD clock skew estimator proposed by Avraham and Pinchas [4] for the fGn case. Figure 7 shows the simulated clock skew performance (MSE) comparison between the OWD clock skew estimator for the Forward path (Eq. 61) with the TWD clock skew estimator proposed by Avraham and Pinchas [4] for the gfGn case. Figure 8 shows the simulated clock skew performance (MSE) comparison between the OWD clock skew estimator for the Reverse path (Eq. 62) with the TWD clock skew estimator proposed by Avraham and Pinchas [4] for the gfGn case. In Figures 2, 3 and Figure 7, the Forward PDV variance was set lower than the Reverse PDV variance (Z > 1). Now, according to test case a of Theorem 8, if Z ≥ Z ̃ f G n / g f G n we should choose the OWD clock skew estimator for the Forward path over the TWD clock skew estimator, else we should choose the TWD clock skew estimator. Namely, if Z ≥ Z ̃ f G n / g f G n , a better clock skew performance from the MSE point of view can be obtained with the OWD clock skew estimator for the Forward path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4]. In Figure 2 and Figure 7, we have Z ≥ Z ̃ f G n / g f G n while in Figure 3 we have Z < Z ̃ f G n / g f G n . Indeed we can see that in Figure 2 and Figure 7 a lower MSE is obtained with the OWD clock skew estimator for the Forward path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4] which clearly demonstrates that the OWD clock skew estimator for the Forward path should be chosen for the clock skew estimation task. In Figure 3 a lower MSE is obtained with the TWD clock skew estimator compared with the OWD clock skew estimator for the Forward path which means that the TWD clock skew estimator should be chosen for the clock skew estimation task. Thus, we may say that according to Figures 2, 3 and Figure 7, test case a of Theorem 8 works correctly. In Figures 4, 5 and Figure 8, the Reverse PDV variance was set lower than the Forward PDV variance (Z < 1). Now, according to test case b of Theorem 8, if Z ≤ Z ̄ f G n / g f G n we should choose the OWD clock skew estimator for the Reverse path over the TWD clock skew estimator, else we should choose the TWD clock skew estimator. Namely, if Z ≤ Z ̄ f G n / g f G n a better clock skew performance from the MSE point of view can be obtained with the OWD clock skew estimator for the Reverse path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4]. In Figure 4 and Figure 8, we have Z ≤ Z ̄ f G n / g f G n while in Figure 5 we have Z > Z ̄ f G n / g f G n . Indeed we can see that in Figure 4 and Figure 8 a lower MSE is obtained with the OWD clock skew estimator for the Reverse path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4] which clearly demonstrates that the OWD clock skew estimator for the Reverse path should be chosen for the clock skew estimation task. In Figure 5 a lower MSE is obtained with the TWD clock skew estimator compared with the OWD clock skew estimator for the Reverse path which means that the TWD clock skew estimator should be chosen for the clock skew estimation task. Thus, we may say that according to Figures 4, 5 and Figure 8, test case b of Theorem 8 works correctly. In Figure 6, the Reverse PDV variance was set equal to the Forward PDV variance (Z = 1). Now, according to test case c of Theorem 8, if σ 2 ≥ σ ̃ f G n / g f G n 2 we should choose the OWD clock skew estimator for the Reverse path over the TWD clock skew estimator, else we should choose the TWD clock skew estimator. Namely, if σ 2 ≥ σ ̃ f G n / g f G n 2 a better clock skew performance from the MSE point of view can be obtained with the OWD clock skew estimator for the Reverse path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4]. According to Figure 6, up to approximately J = 80 we have that σ 2 < σ ̃ f G n / g f G n 2 . Thus, up to approximately J = 80 we see according to Figure 6 that a lower MSE is obtained with the TWD clock skew estimator proposed by Avraham and Pinchas [4] compared with the OWD clock skew estimator for the Reverse path. But, for J > 80, σ 2 > σ ̃ f G n / g f G n 2 thus a lower MSE is obtained with the OWD clock skew estimator for the Reverse path compared with the TWD clock skew estimator proposed by Avraham and Pinchas [4]. Figure 6 clearly demonstrates the effectiveness of test case c of Theorem 8. According to Figures 2–8, our guidelines (conditions) from Theorem 8 indeed may help the system designer to choose wisely the preferred approach among Eqs 60–62 that should be applied for the clock skew estimation task in order to get the best clock skew performance from the MSE point of view. According to Figures 2–8 we can also see the advantage of having two OWD clock skew estimators (one for the Forward path and one associated with the Reverse path) that can supply better clock skew performance from the MSE point of view compared to our recently proposed TWD clock skew estimator Avraham and Pinchas [4], when complying with our proposed guidelines (conditions) from Theorem 8.

Test case a of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Forward path (Eq. 61) and the TWD clock skew estimator (Eq. 60). The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 4 e − 8 [ sec 2 ] , σ ω 2 2 = 6.4 e − 7 [ sec 2 ] , Z = 16, Z ̃ f G n / g f G n ( J = 500 ) = 3.233 (for J < 500, Z ̃ f G n / g f G n < 3.233 ). The results were obtained for 100 Monte-Carlo trails.

Test case a of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Forward path (Eq. 61) and the TWD clock skew estimator (Eq. 60). The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 1 e − 8 [ sec 2 ] , σ ω 2 2 = 1.44 e − 8 [ sec 2 ] , Z = 1.44, Z ̃ f G n / g f G n ( J = 500 ) = 3.058 (for J < 500, Z ̃ f G n / g f G n < 3.058 , Z ̃ f G n / g f G n ( J = 10 ) = 3.004 ). The results were obtained for 100 Monte-Carlo trails.

Test case b of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Reverse path (Eq. 62) and the TWD clock skew estimator (Eq. 60). The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 6.4 e − 7 [ sec 2 ] , σ ω 2 2 = 4 e − 8 [ sec 2 ] , Z = 0.0625, Z ̄ f G n / g f G n ( J = 500 ) = 0.748 (for J < 500, Z ̄ f G n / g f G n < 0.748 , Z ̄ f G n / g f G n ( J = 10 ) = 0.36 ). The results were obtained for 100 Monte-Carlo trails.

Test case b of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Reverse path (Eq. 62) and the TWD clock skew estimator (Eq. 60). The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 1.44 e − 8 [ sec 2 ] , σ ω 2 2 = 1 e − 8 [ sec 2 ] , Z = 0.69, Z ̄ f G n / g f G n ( J = 500 ) = 0.342 7 (for J < 500, Z ̄ f G n / g f G n < 0.342 7 ). The results were obtained for 100 Monte-Carlo trails.

Test case c of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Reverse path (Eq. 62) and the TWD clock skew estimator (Eq. 60). The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.6, d ms = 2.5 ms, d sm = 3 ms, σ ω 1 2 = σ ω 2 2 = 4 e − 6 [ sec 2 ] , Z = 1, σ ̃ f G n / g f G n 2 ( J = 500 ) = 1.03 e − 6 (for J < 500: σ ̃ f G n / g f G n 2 > 1.03 e − 6 , σ ̃ f G n / g f G n 2 ( J = 100 ) = 3 e − 6 , σ ̃ f G n / g f G n 2 ( J = 60 ) = 4.3 e − 6 , σ ̃ f G n / g f G n 2 ( J = 10 ) = 1.72 e − 5 ). The results were obtained for 100 Monte-Carlo trails.

Test case a of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Forward path (Eq. 61) and the TWD clock skew estimator (Eq. 60). The PDV is an gfGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.95, a = 0.08, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 1 e − 8 [ sec 2 ] , σ ω 2 2 = 9 e − 8 [ sec 2 ] , Z = 9, Z ̃ f G n / g f G n ( J = 500 ) = 3.02 (for J < 500, Z ̃ f G n / g f G n < 3.02 ). The results were obtained for 100 Monte-Carlo trails.

Test case b of Theorem 8. Performance (MSE) comparison between the OWD clock skew estimator in the Reverse path (Eq. 62) and the TWD clock skew estimator (Eq. 60). The PDV is an gfGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.95, a = 0.08, d ms = 5 ms, d sm = 5.5 ms, σ ω 1 2 = 9 e − 8 [ sec 2 ] , σ ω 2 2 = 1 e − 8 [ sec 2 ] , Z = 0.111, Z ̄ f G n / g f G n ( J = 500 ) = 0.356 (for J < 500, Z ̄ f G n / g f G n < 0.356 , Z ̄ f G n / g f G n ( J = 10 ) = 0.337 ). The results were obtained for 100 Monte-Carlo trails.

In Figures 9–12 we compared the clock skew performance (MSE) of our new proposed OWD clock skew estimators (Eqs 61, 62) with the clock skew performance (MSE) that is obtained from three clock skew algorithms: a.) TWD clock skew estimator, named as the ML-like estimator (MLLE) proposed by Noh et al. [5], b.) OWD clock skew estimator, named as the maximum likelihood estimator proposed by Levy and Pinchas [7], c.) OWD clock skew estimator, named as the Kalman estimator proposed by Chaloupka et al. [6]. According to Noh et al. [5] we have: β ̂ = T 2,1 J − 1 2 + T 3,1 J − 1 2 T 1,1 J − 1 T 2,1 J − 1 + T 3,1 J − 1 T 4,1 J − 1 − 1 where β ̂ = 1 α ̂ + 1 − 1 T 2,1 (J − 1) = t 2 [J] − t 2 [1], T 1,j (i), T 2,j (i), T 3,j (i) and T 4,j (i) are defined in Eqs 4, 13.

Performance (MSE) comparison between the OWD clock skew estimator for the Forward path (Eq. 61), likelihood clock skew estimator proposed by Noh et al. [5] and maximum likelihood clock skew estimator proposed by Levy and Pinchas [7]. The PDV is an gfGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.95, a = 0.08, d ms = 5.5 ms, d sm = 5 ms, σ ω 1 2 = 4 e − 8 [ sec 2 ] , σ ω 2 2 = 1.6 e − 7 [ sec 2 ] . The results were obtained for 100 Monte-Carlo trails.

Performance (MSE) comparison between the OWD clock skew estimator for the Reverse path (Eq. 62), likelihood clock skew estimator proposed by Noh et al. [5] and maximum likelihood clock skew estimator proposed by Levy and Pinchas [7]. The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5.5 ms, d sm = 5 ms, σ ω 1 2 = 1.6 e − 7 [ sec 2 ] , σ ω 2 2 = 4 e − 8 [ sec 2 ] . The results were obtained for 100 Monte-Carlo trails.

Performance (MSE) comparison between the OWD clock skew estimator for the Forward path (Eq. 61), likelihood clock skew estimator proposed by Noh et al. [5] and Kalman clock skew estimator proposed by Chaloupka et al. [6]. The PDV is an gfGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.95, a = 0.08, d ms = 5.5 ms, d sm = 5 ms, σ ω 1 2 = 4 e − 8 [ sec 2 ] , σ ω 2 2 = 1.3 e − 7 [ sec 2 ] , L = 100, Q KAL = 0, δ σ = δ μ = 1e − 4, μ ̂ [ 1 ] = 0 . The results were obtained for 100 Monte-Carlo trails.

Performance (MSE) comparison between the OWD clock skew estimator for the Reverse path (Eq. 62), likelihood clock skew estimator proposed by Noh et al. [5] and Kalman clock skew estimator proposed by Chaloupka et al. [6]. The PDV is an fGn process. α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.7, d ms = 5.5 ms, d sm = 5 ms, σ ω 1 2 = 3.6 e − 7 [ sec 2 ] , σ ω 2 2 = 4 e − 8 [ sec 2 ] , L = 100, Q KAL = 0, δ σ = δ μ = 1e − 4, μ ̂ [ 1 ] = 0 . The results were obtained for 100 Monte-Carlo trails.

According to Levy and Pinchas [7] we have: ε ̂ = 2 ∑ i = 1 J ∑ j = 1 J ∑ k = 1 min i , j − 1 A max J , i , j , k , H △ T s 1 . i − △ T s 2 . i △ T s 1 . j − △ T s 2 . j ∑ i = 1 J ∑ j = 1 J ∑ k = 1 min i , j − 1 A max J , i , j , k , H △ T m . j △ T s 1 . i − △ T s 2 . i + △ T m . i △ T s 1 . j − △ T s 2 . j − 1 − 1 where A max (J, i, j, k, H) is: A max J , i , j , k , H = Γ J Γ − H + 0.5 Γ J − H + 0.5 2 Γ i − k − H + 0.5 Γ J − H + 0.5 − i + k Γ j − k − H + 0.5 Γ J − H + 0.5 − j + k Γ J − i + k Γ J − j + k Γ i − k Γ j − k − 1 1 i − k 1 j − k − 1 J − i + k 1 J − j + k and ε ̂ is: ε ̂ = 1 α ̂ + 1 − 1 Γ(.) denotes the Gamma function, △ denotes the difference between two consecutive timestamps. T m. i is the timestamp in the ith period when the Master sends the Sync message. T s1. i is the timestamp in the ith period when the dual-Slave receives the Sync message. T s2. i is the timestamp in the ith period when the Slave receives the Sync message.

The Kalman estimator from Chaloupka et al. [6] depends on a predefined parameter L that defines the sliding window’s length in the algorithm. The L parameter impacts the performance (MSE). As we increase L, it reduces the MSE. However, L also depends on the total number of sync periods, which we set for the frequency synchronization task as 500. Therefore, L must be smaller than 500.

According to Chaloupka et al. [6] the Kalman’s measurement equation is: T 1 , j L − T 2 , j L = T 2 , j L α j + Ω 1 , j L

The Kalman’s state equation is: α ̂ j + 1 = α ̂ j + u j . where the variance of u[j] is Q KAL . The estimate of the noise measurement variance is given by Chaloupka et al. [6]: R ̂ j = 1 − δ σ R ̂ j − 1 + δ σ x j − μ ̂ j 2 where μ ̂ j = 1 − δ μ μ ̂ j − 1 + δ μ x j ; x j = T 1 , j L − T 2 , j L δ μ and δ σ are smoothing factors which are between zero and one.

According to Figures 9, 10, our new proposed OWD clock skew estimator for the Forward path (Eq. 61) (Figure 9) or for the Reverse path (Eq. 62) (Figure 10) achieves a lower MSE compared to the clock skew estimators proposed by Noh et al. [5] and Levy and Pinchas [7]. Please note that for the simulation results presented in Figure 10, the PDV for the Reverse path was set lower than the PDV for the Forward path. Since the OWD clock skew estimator proposed by Levy and Pinchas [7] is based on the Forward path only, the clock skew accuracy with this estimator Levy and Pinchas [7] is indeed decreased (Figure 10).

According to Figures 11, 12, our new proposed OWD clock skew estimator for the Forward path (Eq. 61) (Figure 11) or for the Reverse path (Eq. 62) (Figure 12) achieves a lower MSE compared to the clock skew estimators proposed by Noh et al. [5] and Chaloupka et al. [6]. Please note that for the simulation results presented in Figure 12, the PDV for the Reverse path was set lower than the PDV for the Forward path. Since the OWD clock skew estimator proposed by Chaloupka et al. [6] is based on the Forward path only, the clock skew accuracy with this estimator Chaloupka et al. [6] is indeed decreased (Figure 12).

Next we tested our new proposed OWD clock skew estimators’ performances (MSE) for the Forward and Reverse paths (Eqs 61, 62) with our closed-form-approximated expressions for the MSE for the LRD case (Eqs 54, 58) for the Forward and Reverse path, respectively. Case 1 in Figures 13, 14 presents the clock skew performance of our new proposed OWD clock skew estimator for the Forward path (Eq. 61) compared with our closed-form-approximated expression for the MSE in Eq. 54. Case 2 in Figures 13, 14 presents the clock skew performance of our new proposed OWD clock skew estimator for the Reverse path (Eq. 62) compared with our closed-form-approximated expression for the MSE in Eq. 58. According to Figures 13, 14 it can be clearly seen that our new closed-form-approximated expressions for the MSE in Eqs 54, 58 supply very close results to the simulated one.

case 1: Performance comparison for the fGn case, between our new proposed clock skew estimator for the Forward path (Eq. 61) with the performance results for our new proposed expression for the MSE (Eq. 54). σ ω 1 = 0.8 e − 3 [ sec ] . Case 2: Performance comparison for the fGn case, between our new proposed clock skew estimator for the Reverse path (Eq. 62) with the performance results for our new proposed expression for the MSE (Eq. 58). σ ω 2 = 0.1 e − 3 [ sec ] . For both cases: α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.6, d ms = 5.5 ms, d sm = 5 ms. The results were obtained for 100 Monte-Carlo trails.

case 1: Performance comparison for the gfGn case, between our new proposed clock skew estimator for the Forward path (Eq. 61) with the performance results for our new proposed expression for the MSE (Eq. 54). σ ω 1 = 0.8 e − 3 [ sec ] . Case 2: Performance comparison for the gfGn case, between our new proposed clock skew estimator for the Reverse path (Eq. 62) with the performance results for our new proposed expression for the MSE (Eq. 58). σ ω 2 = 0.1 e − 3 [ sec ] . For both cases: α = 50ppm, Q = 5 ms, T s y n = 15.6 m s ( 64 p a c k e t sec ) , H = 0.95, a = 0.95, d ms = 5.5 ms, d sm = 5 ms. The results were obtained for 100 Monte-Carlo trails.

 8 Conclusion

In this paper we derived two novels OWD clock skew estimators for the Forward and Reverse paths applicable for white-Gaussian process and for the fGn/gfGn environment. Those estimators do not depend on the unknown fixed paths nor on the clock offset between the Master and Slave. In addition, we derived also closed-form-approximated expressions for the clock skew performance (MSE) for the new proposed OWD clock skew estimators for the Forward and Reverse paths. In order to help the system designer to choose the right clock skew estimator that may get the best clock skew performance from the MSE point of view, some guidelines (conditions) were derived, helping choosing the right clock skew estimator wisely. Simulation results has confirmed that our new OWD clock skew estimators indeed achieve better clock skew performance from the MSE point of view compared to the literature known clock skew estimators. Simulation results have also confirmed that our closed-form-approximated expressions for the MSE related to our new proposed OWD Forward and Reverse estimators are indeed efficient.

 Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

 Author Contributions

All authors listed have made a substantial, direct, and intellectual contribution to the work and approved it for publication.

 Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

 Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

 References [Dataset] Arnold D . IEEE 1588-2019-IEEE Standard for a Precision Clock Synchronization Protocol for Networked Measurement and Control Systems (2019). Available from: https://standards.ieee.org/standard/1588-2019.html (Accessed June 17, 2021). K.Karthik A S.Blum R . Estimation Theory Based Robust Phase Offset Estimation in the Presence of Delay Attacks (2016). Available from: https://arxiv.org/pdf/1611.05117.pdf (Accessed August 10, 2021). Guruswamy A Blum RS Kishore S Bordogna M . On the Optimum Design of L-Estimators for Phase Offset Estimation in IEEE 1588. IEEE Trans Commun (2015) 63:5101–15. 10.1109/TCOMM.2015.2493534 Avraham Y Pinchas M . A Novel Clock Skew Estimator and its Performance for the IEEE 1588v2 (PTP) Case in Fractional Gaussian Noise/Generalized Fractional Gaussian Noise Environment. Front Phys (2021) 9:1–21. 10.3389/fphy.2021.796811 Noh K-L Chaudhari QM Serpedin E Suter BW . Novel Clock Phase Offset and Skew Estimation Using Two-Way Timing Message Exchanges for Wireless Sensor Networks. IEEE Trans Commun (2007) 55:766–77. 10.1109/TCOMM.2007.894102 Chaloupka Z Alsindi N Aweya J . Clock Skew Estimation Using Kalman Filter and IEEE 1588v2 PTP for Telecom Networks. IEEE Commun Lett (2015) 19:1181–4. 10.1109/LCOMM.2015.2427158 Levy C Pinchas M . Maximum Likelihood Estimation of Clock Skew in IEEE 1588 with Fractional Gaussian Noise. Math Probl Eng (2015) 2015:1. 10.1155/2015/174289 Karthik AK Blum RS . Robust Clock Skew and Offset Estimation for IEEE 1588 in the Presence of Unexpected Deterministic Path Delay Asymmetries. IEEE Trans Commun (2020) 68:5102–19. 10.1109/TCOMM.2020.2991212 Karthik AK Blum RS . Robust Phase Offset Estimation for IEEE 1588 PTP in Electrical Grid Networks. In: 2018 IEEE Power & Energy Society General Meeting (2018). 10.1109/PESGM.2018.8586488 Karthik AK Blum RS . Optimum Full Information, Unlimited Complexity, Invariant, and Minimax Clock Skew and Offset Estimators for IEEE 1588. IEEE Trans Commun (2019) 67:3624–37. 10.1109/TCOMM.2019.2900317 Li M Zhao W . On Bandlimitedness and Lag-Limitedness of Fractional Gaussian Noise. Physica A: Stat Mech its Appl (2013) 392:1955–61. 10.1016/j.physa.2012.12.035 Cappe O Moulines E Pesquet J-C Petropulu AP Xueshi Y . Long-range Dependence and Heavy-Tail Modeling for Teletraffic Data. IEEE Signal Process Mag (2002) 19:14–27. 10.1109/79.998079 Li M . Generalized Fractional Gaussian Noise and its Application to Traffic Modeling. Physica A: Stat Mech its Appl (2021) 579:126138. 10.1016/j.physa.2021.126138