## 4 Satellite Synchronization  
As introduced before, by exploiting the dual one-way technique, it is possible to compute the clock bias between the two spacecraft at the time of the measurement. Unfortunately, computing only the time bias for a given specific time is not sufficient to model the deviation between two clocks. A more complex model that considers deterministic and random components of the clock deviation, such as bias, drift, and drift rate, is usually needed. A common solution is the use of stochastic differential equations. The clock time and frequency errors for a time epoch \(k\) can be derived from the previous one by exploiting, for example, the model proposed in [31] where two main clock noises are considered: White Frequency Modulation (WFM) and Random Walk Frequency Modulation (RWFM). These two noises generate on the clock's phase a Wiener noise (integral of the WFM) and an integrated Wiener noise (integral of the RWFM):  
\[X_{1}(t_{k+1})= X_{1}(t_{k})+(X_{2}(t_{k})+\mu_{1})\tau+\mu_{2}\frac{\tau^{2}}{2 }+J_{k,1} \tag{31}\] \[X_{2}(t_{k+1})= X_{2}(t_{k})+\mu_{2}\tau+J_{k,2}\] (32) \[\mathbf{J}_{k}\sim \mathcal{N}\left(0, \begin{bmatrix}\sigma_{1}^{2}\tau+\sigma_{2} ^{2}\frac{\tau^{3}}{2}&\sigma_{2}^{2}\frac{\tau^{2}}{2}\\ \sigma_{2}^{2}\frac{\tau^{2}}{2}&\sigma_{2}^{2}\tau\end{bmatrix} \right) \tag{33}\]  
where \(X_{1}\) represents the clock phase deviation (time Bias), \(X_{2}\) represents the clock frequency deviation, and \(\mu_{1}\) and \(\mu_{2}\) are the clock deterministic drift terms. The vector \(\mathbf{J}_{k}\) is a bivariate random variable with zero mean and variance that depends on two parameters: \(\sigma_{1}\), which is the diffusion coefficient of the phase noise, and \(\sigma_{2}\), which is the diffusion coefficient of the frequency noise; it represents the covariance matrix of the error related to the specific clock. By setting the \(\mu_{1}\), \(\mu_{2}\), \(\sigma_{1}\), and \(\sigma_{2}\), it is possible to set the noise model representing the desired type of clock to be represented for a long time. For short-time modeling, this model can be replaced by a simpler one that takes into account only the deterministic components of the clock deviation, which are the time bias, the drift, and the drift rate using the following polynomial representation for the bias of each satellite clock at given epoch (\(t_{0}\)) [32; 33]:  
\[\delta t=a_{0}+a_{1}(t-t_{0})+a_{2}(t-t_{0})^{2} \tag{34}\]  
where \(a_{0}\) is the clock offset, \(a_{1}\) is the clock drift, and \(a_{2}\) is the clock drift rate. This simpler model is commonly used to represent the clock error in the GNSS system and can be accurate enough to represent the clock error if it is frequently updated. For example, in the GPS system, the parameters of the clock are updated every two hours [34]. In this work, this second approach is proposed and, considering the main spacecraft clock as the reference one, the clock parameters to compute the bias with respect to the main spacecraft for the anchor satellite \(i\) can be collected in the vector \(\mathbf{a}_{i}=\left[a_{0}^{i},a_{1}^{i},a_{2}^{i}\right]\). Recalling that, for each time epoch \(t_{k}\), the TOA measurements are able to produce an estimation of the time bias \(\delta t(t_{k})\) between the spacecraft clocks as derived in Equation (4),for each spacecraft, it is possible to estimate \(a_{0}^{i},a_{1}^{i},a_{2}^{i}\) by linear regression or a Least Square solution, inverting the system of equations:  
\[\left\{ \begin{array}{rl}\delta t_{i}(t_{1})&=a_{0}^{i}+a_{1}^{i}(t_{1}-t_{0})+a _{2}^{i}(t_{1}-t_{0})^{2}+\epsilon(t_{1})\\ &\vdots\\ \delta t_{i}(t_{k})&=a_{0}^{i}+a_{1}^{i}(t_{k}-t_{0})+a_{2}^{i}(t_{k}-t_{0})^{2 }+\epsilon(t_{1})\\ &\vdots\\ \delta t_{i}(t_{K})&=a_{0}^{i}+a_{1}^{i}(t_{K}-t_{0})+a_{2}^{i}(t_{K}-t_{0})^{ 2}+\epsilon(t_{1})\\ \end{array} \right. \tag{35}\]  
where \(K\) is the number of epochs used for the estimation.