Vector Tracking Loop Jitter

Abstract

The vector loop architecture is an advanced signal tracking algorithm for global navigation satellite system (GNSS) receivers. It has long been the subject of much research for its ability to address GNSS vulnerability. This paper presents a methodology to obtain the total jitter of the noncoherent vector tracking architecture using its state-space (SS) formulation, referred to as SS jitter, and relates the SS jitter to that determined with transfer function (TF) models. Thus, a theoretical framework for evaluating the performance of noncoherent vector loops is developed. In this regard, vector TFs are derived in the z-domain without decoupling the delay- and frequency-locked loops, unlike previous work. Two different noise equivalent bandwidths from the z-domain TF are also calculated. A close agreement between the SS and TF jitters is demonstrated, proving their equivalence. A metric to represent the collective performance of vector loops for loss-of-lock detection is identified using dynamically allocated bounds on the jitter. The steady-state jitter performance at different noise levels in the absence of dynamic stress error is analyzed in detail. Finally, using high-dynamic simulations of the Global Positioning System and Navigation with Indian Constellation signals, it is shown that the proposed metric holds promise for detecting loss of lock in the presence of jerk.

Keywords

1 INTRODUCTION

The vector tracking architecture is an advanced receiver algorithm developed to address the vulnerability of global navigation satellite system (GNSS) signals (Spilker Jr., 1995). Since its introduction in the 1980s (Copps et al., 1980), extensive research has been carried out to explore its potential for improving GNSS robustness against signal degradation, outages, and high platform dynamics (Lashley et al., 2009). Different variants of this tracking method with enhanced capabilities have also been proposed in the literature, including joint vector tracking (Liu et al., 2011), maximum likelihood vector architecture (Weill, 2010), differential vector phase-locked loops (Brewer & Raquet, 2016), and multi-receiver vector tracking (Ng & Gao, 2017). The superiority of the vector receiver over the traditional scalar one, however, comes at the price of increased computational complexity. Therefore, recent work has focused on reducing the computational load without sacrificing performance (Mu & Long, 2017; Abedi & Mosavi, 2022; Gómez et al., 2023), thus paving the way for real-time vector loop implementations. Multi-GNSS vector tracking receivers have also been studied in the literature for improved performance (Shytermeja, 2017; Tabatabaei et al., 2017).

Traditionally, the performance of a scalar tracking loop architecture is quantitatively evaluated by its total jitter (Ward et al., 2006). The jitter is determined by the sum of the 3−σ tracking error due to wideband and colored noises and dynamic stress error arising from abrupt changes in platform motion. Ward et al. (2006) studied different sources of jitter in detail for both conventional phase- and frequency-locked loops using appropriate transfer function (TF) models. With the help of the total jitter, Razavi et al. (2008) explained the relative importance of different tracking parameters in systematically designing a scalar loop. Jitter was also shown to be a reliable metric for loss-of-lock detection if the scalar tracking loop bandwidth (BW) is not too low. Yang et al. (2017) and Mute and Bhattacharyya (2023) analyzed how the tracking loop jitter varies in adaptive scalar loops whose parameters are adaptively chosen based on carrier-to-noise ratios (C/N₀) and platform dynamics (or line-of-sight jerk).

While the jitter of scalar tracking loops has been well studied, a similar treatment is lacking in the context of more complex vector architecture. Lashley et al. (2009) provided a methodology to calculate the vector jitter considering thermal noise, clock noise, and dynamic stress error for a noncoherent vector loop with 20-ms measurement update intervals. In their approach, the calculation of jitter due to noise involved updating the state error covariance using a two-step process. The dynamic stress error was shown to be derived from the state error difference between two separate Kalman filters, one of which had knowledge of the dynamics and, therefore, accurately modeled the change in acceleration and the other did not. As channels are coupled in vector tracking mode, the maximum jitter among all channels was compared against the threshold for loss-of-lock analysis. Simulation results showed that during an 8-g maneuver at 16 dB-Hz, the velocity estimates of the vector loops began to diverge. The jitters of the majority of channels also exceeded the threshold, thereby correctly predicting loss of lock. However, the performance was not validated at different C/N₀ ratios against specific platform jerk levels. To obtain the vector loop jitter with TF models similar to that of its scalar counterpart, appropriate TFs of the time-varying vector tracking architecture are needed. While Bhattacharyya (2018) described s-domain TF models and noise equivalent BWs of vector loops under certain assumptions, to the best of the author's knowledge, there is little work on vector jitter calculation using TFs.

This paper determines the total jitter of a noncoherent vector architecture with both state-space (SS) and TF formulations. The major contributions of this work are as follows. First, a theoretical framework for evaluating the performance of vector loops is presented in terms of SS- and TF-based tracking jitters. More precisely, a priori and a posteriori jitters based on the theoretical SS model, termed SS jitters, are obtained. Then, z-domain TF models of the vector loops are determined without decoupling the delay- and frequency-locked loops, unlike the work by Bhattacharyya (2018). The current TF formulation removes all other key assumptions of the previous work, is more generic in nature, and has a much simpler implementation. Two different noise equivalent BWs, one each for measurement noise and process noise, are also calculated. A close agreement between the SS and TF jitters is established to demonstrate their equivalence. Second, with the help of the SS jitter, a practical loss-of-lock detection algorithm is developed for the noncoherent vector loops. The SS jitters are initially obtained at the measurement update intervals. To reduce the computational load, the measurement updates of the noncoherent vector loops are often performed at an interval larger than the coherent integration time, which is usually 10 ms for carrier frequency tracking and 20 ms for code phase (Bhattacharyya & Gebre-Egziabher, 2010). Therefore, intermediate bounds on the jitter are calculated every 20 ms in order to apply existing thresholds for loss-of-lock detection. A suitable metric for representing the collective performance of all channels is also identified. The proposed metric is used to detect loss of lock in the presence of jerk in a high-dynamic simulation. It should be noted that the simulations are performed with the Global Positioning System (GPS) L1 and Navigation with Indian Constellation (NavIC) L5 signals. Thus, both geometric and frequency diversities are attained in a single scenario, which demonstrates the effectiveness of the proposed methodology for different constellations as well as different signals.

The remainder of this paper is organized as follows. First, the approach to calculating the SS jitter is described. Then, TF models are introduced, and noise equivalent BWs are determined. Using this, the TF jitter is obtained. Subsequently, an appropriate metric is proposed for loss-of-lock detection with the help of the SS jitter. Next, simulation studies show the equivalence between the SS and TF jitters, and the vector jitter performance is analyzed for different noise and jerk levels. Finally, conclusions are provided.

2 VECTOR LOOP JITTER

In this paper, a noncoherent vector tracking architecture is considered for GPS and NavIC constellations (see Figure 1) to determine the vector loop jitter. This architecture consists of vector delay-locked loops (VDLLs) for code phase tracking and vector frequency-locked loops (VFLLs) for carrier frequency tracking. The state transition and linearized measurement models of the vector navigation filter, which uses an extended Kalman filter (EKF) implementation, are as follows:

FIGURE 1

Noncoherent vector loop architecture

Appropriate variables at the channel and navigation levels are indicated in the figure. f_code and $f_{carr}^{G / N}$ are the code chipping rate and GPS/NavIC carrier frequency, respectively, and c is the speed of light.

${\underset{―}{ϰ}}_{ℓ} = F {\underset{―}{ϰ}}_{ℓ - 1} + {\underset{―}{w}}_{ℓ - 1}$ 1

${\underset{―}{z}}_{ℓ} = H_{ℓ} Δ {\underset{―}{ϰ}}_{ℓ} + {\underset{―}{ζ}}_{ℓ}$ 2

The vector $Δ {\underset{―}{ϰ}}_{ℓ}$ at t_ℓ is given by ( ${\underset{―}{ϰ}}_{ℓ} - {\hat{ϰ}}_{ℓ}^{-}$ ). The true state vector is $\underset{―}{ϰ}$ , and its a priori estimate is ${\underset{―}{\hat{x}}}^{-} . \underset{―}{ϰ}$ comprises the user position, velocity, and acceleration in the Earth-centered Earth-fixed (ECEF) coordinate frame as well as GPS and NavIC receiver clock biases, drifts, and drift rates; thus, this vector has 15 elements. It should be noted that the position estimation error in the ECEF frame can later be converted to the local navigation frame to compute the vertical and horizontal position errors. In the filter, ${\hat{\underset{―}{x}}}^{-}$ and the a posteriori estimate of $\underset{―}{ϰ}, {\hat{\underset{―}{x}}}^{+}$ , are propagated from a time instant $t_{k - 1}^{'}$ to $t_{k}^{'}$ by multiplying the matrix, $F^{'}$ , which is given as follows:

$F^{'} = I_{5 \times 5} \otimes [\begin{array}{c} 1 & Δ t^{'} & ({(Δ t^{'})}^{2}) / 2 \\ 0 & 1 & Δ t^{'} \\ 0 & 0 & 1 \end{array}]$ 3

where ⊗ denotes the Kronecker tensor product and I_5×5 is a 5×5 identity matrix. $t_{k}^{'}$ increments every 20 ms; thus, $Δ t^{'} = (t_{k}^{'} - t_{k - 1}^{'}) = 20 ms$ . In contrast, the subscript ℓ in Equation (1) indicates the time instants for measurement updates (i.e., t_ℓ). In the implementation of this study, t_ℓ increments every 0.1 s (i.e., $Δ t = (t_{ℓ + 1} - t_{ℓ}) = 0.1 s$ ). It is apparent that $t_{ℓ} = t_{k}^{'}$ when k = 5ℓ (with $ℓ = 1, 2, \dots$ ). Hence, $F = {(F^{'})}^{5}$ . The process noise vector, ${\underset{―}{w}}_{ℓ - 1}$ , is assumed to be added to obtain the true states every $t_{ℓ}$ , as shown in Equation (1). ${\underset{―}{w}}_{ℓ - 1}$ has zero mean and error covariance $Q . {\underset{―}{z}}_{ℓ}$ is the measurement innovation vector at $t_{ℓ}$ , and ${\underset{―}{ζ}}_{ℓ}$ is the measurement noise vector at $t_{ℓ}$ with zero mean and covariance $W_{ℓ} \cdot \underset{―}{w}$ and $\underset{―}{ζ}$ are unknown vectors; only their statistical properties (i.e., mean and covariance) can be predicted. ${\underset{―}{z}}_{ℓ}$ is obtained from five 20-ms code phase (early-minus-late envelope) and carrier frequency (atan2) discriminator outputs averaged over 0.1 s. It should be noted that in the beginning, prior to the transition from scalar to vector loops, the coherent integration time of the correlators is synchronized with data bit transitions. The data bit boundaries are continuously tracked afterward by counting 20 1-ms coarse acquisition (C/A)-code periods. Thus, the two 10-ms prompt correlator pair intervals of the atan2 carrier frequency discriminator are ensured not to fall across data bit boundaries. The measurement model matrix, $H_{ℓ}$ , relates ${\underset{―}{z}}_{ℓ}$ to $Δ {\underset{―}{ϰ}}_{ℓ}$ , as defined in Appendix A.1. Interested readers can also refer to the supplementary material of Mute and Bhattacharyya (2023) for a detailed discussion of the measurement model and code and carrier numerically controlled oscillator corrections of the dual-constellation GPS L1 and NavIC L5 vector tracking architecture.

If ${\underset{―}{e}}_{ℓ}^{+}$ is the state estimation error after a measurement update (i.e., ${\underset{―}{e}}_{ℓ}^{+} = ({\underset{―}{ϰ}}_{ℓ} - {\hat{\underset{―}{x}}}_{ℓ}^{+})$ ), then it can be propagated to the next time step as follows (Simon, 2006):

${\underset{―}{e}}_{ℓ + 1}^{+} = (I - K_{ℓ + 1} H_{ℓ + 1}) F {\underset{―}{e}}_{ℓ}^{+} + (I - K_{ℓ + 1} H_{ℓ + 1}) {\underset{―}{w}}_{ℓ} - K_{ℓ + 1} {\underset{―}{ζ}}_{ℓ + 1}$ 4

where $K_{ℓ + 1}$ is the Kalman gain at $t_{ℓ + 1}$ . From Equation (4), the propagation of the a priori state error before a measurement update, ${\underset{―}{e}}_{ℓ}^{-} (= Δ {\underset{―}{x}}_{ℓ})$ , is given as follows:

${\underset{―}{e}}_{ℓ + 1}^{-} = (F - K_{ℓ}^{'} H_{ℓ}) {\underset{―}{e}}_{ℓ}^{-} + {\underset{―}{w}}_{ℓ} - K_{ℓ}^{'} {\underset{―}{ζ}}_{ℓ}$ 5

where $K_{ℓ}^{'} = F K_{ℓ}$ . The vector loop jitter is defined for the VDLL and VFLL, using the a priori and a posteriori state estimation errors, ${\underset{―}{e}}_{ℓ + 1}^{-}$ and ${\underset{―}{e}}_{ℓ + 1}^{+}$ . This is because the code phase and carrier frequency of a channel are aligned with the incoming code phase and carrier frequency based on the estimated or predicted state vector of the vector navigation filter, as shown in Figure 1. Therefore, the state estimation/prediction error at the navigation level will reflect as code phase and carrier frequency alignment errors at the channel level. Hence, the a posteriori (and a priori) code phase and carrier Doppler errors together can be written as ${\underset{―}{ϵ}}_{e, ℓ + 1} = H_{ℓ + 1} {\underset{―}{e}}_{ℓ + 1}^{+}$ (and ${\underset{―}{ϵ}}_{p, ℓ + 1} = H_{ℓ + 1} {\underset{―}{e}}_{ℓ + 1}^{-}$ ). The vectors ${\underset{―}{ϵ}}_{e, ℓ + 1}$ and ${\underset{―}{ϵ}}_{p, ℓ + 1}$ are of dimension $2 n \times 1$ , where n is the number of visible satellites used for positioning. The code phase errors are given by the first n elements, and the carrier Doppler errors are obtained from the next n elements. The total a posteriori VDLL (or VFLL) jitter of a channel comprises the $3 - σ$ a posteriori code phase (or carrier Doppler) estimation error due to random effects or noise plus the transient/steady-state error due to abrupt changes in platform dynamics (i.e., the dynamic stress error, which is treated as a deterministic quantity). This is an extension of the definition provided by Ward et al. (2006) for scalar tracking loops. The total a priori jitter is defined in a similar manner. As the vector navigation filter estimates the user position, velocity, and acceleration, the dynamic stress error primarily arises from the user jerk. In this paper, jitter is studied after the initial transients of the navigation filter have died down.

To calculate the vector loop jitter for each channel, the effects of the process and measurement noises should be reflected in ${\underset{―}{ϵ}}_{e, ℓ + 1}$ through ${\underset{―}{e}}_{ℓ + 1}^{+}$ (and in ${\underset{―}{ϵ}}_{p, ℓ + 1}$ through ${\underset{―}{e}}_{ℓ + 1}^{-}$ ). For this purpose, the matrix pair $(H_{ℓ + 1}, (I - K_{ℓ + 1} H_{ℓ + 1}) F)$ must be observable in the case of $ϵ_{—e, ℓ + 1}$ , where $H_{ℓ + 1}$ is the measurement model matrix and $(I - K_{ℓ + 1} H_{ℓ + 1}) F$ is the system model matrix (Kailath, 1980). This is proved in Appendix A.2, assuming that $(H_{ℓ + 1}, F)$ is observable for a stable vector navigation filter and F is invertible (valid for the adopted model). It is also proved that the matrix pairs $((I - K_{ℓ + 1} H_{ℓ + 1}) F, (I - K_{ℓ + 1} H_{ℓ + 1}))$ and $((I - K_{ℓ + 1} H_{ℓ + 1}) F, K_{ℓ + 1})$ are reachable or stabilizable with regard to the noises ${\underset{―}{w}}_{ℓ}$ and ${\underset{―}{ζ}}_{ℓ + 1}$ for ${\underset{―}{e}}_{ℓ + 1}^{+}$ , respectively, where the first matrix of a pair is the system model matrix and the second one is the coefficient matrix for the input vector (see Appendix A.3–A.4 for derivations). It should be noted that the matrix pairs written as before are the standard representations of matrices in control theory for defining the system observability and reachability/stabilizability. While the latter is not essential in the usual sense for jitter calculation, for completeness, the proof in the Appendix justifies that some modes of ${\underset{―}{e}}_{ℓ + 1}^{+}$ are affected by noise. Hence, ${\underset{―}{ϵ}}_{e, ℓ + 1}$ would reflect the effects of ${\underset{―}{w}}_{ℓ}$ and ${\underset{―}{ζ}}_{ℓ + 1}$ on ${\underset{―}{e}}_{ℓ + 1}^{+}$ through observability. The observability of $ϵ_{p, ℓ + 1}$ is shown as part of the proof in Appendix A.2. The reachability criteria can be easily proved following the approach of the Appendix. Observability and stabilizability will also be assessed later in the context of the dynamic stress error. Next, the mathematical expression of the total SS jitter for vector loops is provided.

2.1 State-Space Jitter

Simon (2006) showed that for a given Kalman gain matrix, the steady-state estimation error covariance can be written as the sum of the covariance due to process noise only and the covariance due to measurement noise only. That is, after the initial transients die down, the estimation error covariance is determined by the process and measurement noises. Therefore, if $Σ_{p, ℓ}^{+}$ and $Σ_{m, ℓ}^{+}$ denote the a posteriori error covariances only for $\underset{―}{w}$ and $\underset{―}{ζ}$ , respectively, then the a posteriori error covariance, $Σ_{ℓ}^{+}$ , after the initial transients subside, is $(Σ_{p, ℓ}^{+} + Σ_{m, ℓ}^{+})$ . The same can also be written for the a priori error covariance, $Σ_{ℓ}^{-}$ . Thus, two jitters due to noise (called noise jitters later) can be calculated, one each with $Σ_{ℓ}^{-}$ and $Σ_{ℓ}^{+}$ . These terms are called the a priori and a posteriori SS noise jitters, respectively, in this paper. The a posteriori SS noise jitter for channel i in the VDLL is given by the following:

${\underset{―}{σ}}_{V D L L_{noise}, ℓ}^{+} (i) = \sqrt{G_{ℓ}^{+} (i, i)}$ 6

where $G_{ℓ}^{+} = H_{ℓ} Σ_{ℓ}^{+} H_{ℓ}^{T}$ and $i = 1, \dots, n$ . n is the number of satellites used for positioning at epoch ℓ and can change with time. ${\underset{―}{σ}}_{V D L L_{noise}, ℓ}^{+}$ represents the VDLL noise jitter vector for all channels, whereas its argument denotes the i-th channel. For the VFLL, the a posteriori noise jitter of channel η is given as follows:

${\underset{―}{σ}}_{V F L L_{noise}, ℓ}^{+} (η) = \sqrt{G_{ℓ}^{+} (η, η)}, η = n + 1, \dots, 2 n$ 7

Similarly, using $Σ_{ℓ}^{-}$ , the a priori noise jitters for the VDLL and VFLL, respectively, are as follows:

${\underset{―}{σ}}_{V D L L_{noise}, ℓ}^{-}^{(i)} = \sqrt{G_{ℓ}^{-} (i, i)} and {\underset{―}{σ}}_{V F L L_{noise}, ℓ}^{-} (η) = \sqrt{G_{ℓ}^{-} (η, η)}$ 8

where $G_{ℓ}^{-} = H_{ℓ} Σ_{ℓ}^{-} H_{ℓ}^{T}$ . Following the definition after Equation (5), the $1 - σ$ total jitter of a channel is given by the sum of the noise jitter and one-third of the dynamic stress error, which separately accounts for the effect of the nonzero mean of $\underset{―}{w}$ . It should be noted that $\underset{―}{w}$ is assumed to have zero mean, and its covariance, Q, is only considered in the previous noise jitter calculation. Because the user acceleration is modeled as constant in the vector navigation filter, user jerk will result in a dynamic stress error. In this paper, the final jitter for loss-of-lock detection is determined later using both a priori and a posteriori jitters. It should be noted that a separate consideration is made for the contributions of process noise (or Q) in the jitter and of dynamic stress error arising from jerk. The diagonal elements of Q corresponding to the user acceleration may have comparatively higher values than other elements so that the filter appropriately learns from its measurements during dynamics. However, the actual jerk need not be the same as what is expected from Q. Hence, the dynamic stress error must be included. While considering the full Q in the noise jitter along with the dynamic stress error may result in an overbound of the total jitter, simulation results will show that this does not lead to a false detection of loss of lock. This is partially due to the fact that the contribution of the diagonal acceleration terms of Q to $G_{ℓ}^{-}$ via the addition of Q in $Σ_{ℓ}^{-}$ is zero because the respective columns of $H_{ℓ}$ for the acceleration are zero (see Appendix A.1). Considering only the clock part of Q to separately update the noise jitter in a two-step process following the work of Lashley et al. (2009) would require fairly accurate knowledge of the clock parameters. This may not always be possible, and parameters may change with time. Thus, as a conservative approach, the full Q is preferred in this paper to overbound the noise jitter.

The dynamic stress error is calculated at each epoch for a given user jerk profile, instead of finding the steady-state error for a constant jerk level. This approach is used to account for real scenarios in which the jerk need not remain constant and can change with time. The calculation of the dynamic stress error is described next. The a priori state error due to nonzero jerk, $Δ {\underset{―}{ϰ}}_{d y n, ℓ}^{-}$ , before the measurement update at ℓ is as follows:

$Δ {\underset{―}{ϰ}}_{d y n, ℓ}^{-} = F Δ {\underset{―}{ϰ}}_{d y n, ℓ - 1}^{+} + B {\underset{―}{u}}_{j e r k, ℓ - 1}$ 9

where $Δ {\underset{―}{ϰ}}_{dyn, ℓ - 1}^{+}$ is the cumulative state error due to previous jerks after the measurement update at $ℓ - 1$ and is a zero vector when $ℓ = 1$ . $B {\underset{―}{u}}_{jerk, ℓ - 1}$ is the error introduced to the modeled user position, velocity, and acceleration at epoch ℓ due to user jerk, ${\underset{―}{u}}_{jerk, ℓ - 1}$ , at epoch $ℓ - 1 . {\underset{―}{u}}_{jerk, ℓ - 1}$ has three components given by ${[u x_{{jerk}_{ℓ - 1}} u y_{{jerk}_{ℓ - 1}} u z_{{jerk}_{ℓ - 1}}]}^{T}$ , which are the ECEF x-, y-, and z-components of jerk, respectively. The first nine rows of B are calculated over the time interval $Δ t (= 0.1 s)$ as follows:

$B (1 : 9, :) = I_{3 \times 3} \otimes {[\frac{Δ t^{3}}{6} \frac{Δ t^{2}}{2} Δ t]}^{T}$ 10

The last six rows of B have zero elements corresponding to the GPS and NavIC receiver clock errors. Following Equations (4) and (5), the state error, $Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{+}$ , after the measurement update and $Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{-}$ are as follows:

$Δ {\underset{―}{ϰ}}_{d y n, ℓ}^{+} = (I - K_{ℓ} H_{ℓ}) F Δ {\underset{―}{ϰ}}_{d y n, ℓ - 1}^{+} + (I - K_{ℓ} H_{ℓ}) B {\underset{―}{u}}_{j e r k, ℓ - 1}$ 11

$Δ {\underset{―}{ϰ}}_{d y n, ℓ}^{-} = F (I - K_{ℓ} H_{ℓ}) Δ {\underset{―}{ϰ}}_{d y n, ℓ - 1}^{-} + B {\underset{―}{u}}_{j e r k, ℓ - 1}$ 12

The a priori dynamic stress error is ${\underset{―}{ϵ}}_{dyn, ℓ}^{-} = H_{ℓ} Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{-}$ , and the a posteriori dynamic stress error is ${\underset{―}{ϵ}}_{dyn, ℓ}^{+} = H_{ℓ} Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{+}$ . Because $(H_{ℓ}, (I - K_{ℓ} H_{ℓ}) F)$ is observable and $((I - K_{ℓ} H_{ℓ}) F, (I - K_{ℓ} H_{ℓ}))$ is stabilizable (see Appendix A.2 and A.3), the effect of $B {\underset{―}{u}}_{jerk, ℓ - 1}$ is reflected in ${\underset{―}{ϵ}}_{dyn, ℓ}^{+}$ . A similar justification holds for ${\underset{―}{ϵ}}_{dyn, ℓ}^{-}$ . Using Equation (8) and the dynamic stress error, the total a priori $1 - σ$ SS jitters of a VDLL channel i and VFLL channel η before a measurement update are as follows:

${\underset{―}{σ}}_{V D L L, ℓ}^{-} (i) = {\underset{―}{σ}}_{V D L L_{n o i s e}, ℓ}^{-} (i) + \frac{| {\underset{―}{ϵ}}_{dyn, ℓ}^{-} (i) |}{3}$ 13

${\underset{―}{σ}}_{V F L L, ℓ}^{-} (η) = {\underset{―}{σ}}_{V F L L_{n o i s e}, ℓ}^{-} (η) + \frac{| {\underset{―}{ϵ}}_{dyn, ℓ}^{-} (η) |}{3}$ 14

where |.| represents the absolute value. The “−” in the superscript of the preceding equations is replaced with “+” after a measurement update to represent the total a posteriori $1 - σ$ SS jitters. It should be noted that these jitters are determined in this paper using theoretical models only. Next, the vector loop jitter is determined from its TF models.

2.2 Transfer Function Jitter

Traditionally, the jitter is calculated for a scalar tracking loop based on its TF model. Following this approach, the TF jitter is determined for the vector architecture and is later used to validate the SS jitter. Becaue TFs require constant parameters, first, piecewise constant models of vector loops are developed over time intervals, each of duration $T_{sub}$ , without decoupling the VDLL and VFLL. It should be noted that even when there is no dynamic stress error (i.e., for a static user or constant-velocity/acceleration user motion), the measurement model matrix, H, would change with time due to user-satellite relative motion, making the system model time-varying. K may also change owing to the adaptive calculation of W based on $C / N_{0}$ . Therefore, piecewise constant models are required for the TF formulations and subsequent jitter calculations. An algorithm to adaptively calculate $T_{sub}$ for a model is given in Appendix A.5. Over each interval of $T_{sub}$ duration, H, W, and K are assumed to be constant. It should be noted that the minimum value of $T_{sub}$ is 0.1 s, which is the measurement update interval. In this case, a TF model would last only 0.1 s. Using Equations (4) and (5), ${\underset{―}{e}}_{ℓ + 1}^{+}$ and ${\underset{―}{e}}_{ℓ + 1}^{-}$ for a given piecewise constant model q of duration $T_{sub}$ reduce to the following:

${\underset{―}{e}}_{ℓ + 1}^{+} = (I - K_{q} H_{q}) F {\underset{―}{e}}_{ℓ}^{+} + (I - K_{q} H_{q}) {\underset{―}{w}}_{ℓ} - K_{q} {\underset{―}{ζ}}_{ℓ + 1}$ 15

${\underset{―}{e}}_{ℓ + 1}^{-} = (F - K_{q}^{'} H_{q}) {\underset{―}{e}}_{ℓ}^{-} + {\underset{―}{w}}_{ℓ} - K_{q}^{'} {\underset{―}{ζ}}_{ℓ}$ 16

where $K_{q}$ and $H_{q}$ are the optimal Kalman gain and actual H matrices at the beginning of $T_{sub}$ for the model q, respectively. These matrices are kept constant until a new piecewise constant model $(q + 1)$ starts. The optimal K implies the Kalman gain obtained using the EKF equations with the actual time-varying model of the vector navigation filter. ℓ indicates the indices for the time instants of measurement updates, as before. Because the parameters of the preceding equations are constant, it is possible to find their z-transforms with respect to $\underset{―}{w}$ and $\underset{―}{ζ}$ as follows:

$\begin{matrix} {\underset{―}{e}}_{q}^{+} (z) = {\underset{―}{e}}_{1, q}^{+} (z) + {\underset{―}{e}}_{2, q}^{+} (z) and {\underset{―}{e}}_{q}^{-} (z) = {\underset{―}{e}}_{1, q}^{-} (z) + {\underset{―}{e}}_{2, q}^{-} (z) where \\ {\underset{―}{e}}_{1, q}^{+} (z) = {(z I - (I - K_{q} H_{q}) F)}^{- 1} (I - K_{q} H_{q}) \underset{―}{w} (z), and \\ {\underset{―}{e}}_{2, q}^{+} (z) = - {(z I - (I - K_{q} H_{q}) F)}^{- 1} K_{q} z \underset{―}{ζ} (z) \\ {\underset{―}{e}}_{1, q}^{-} (z) = {(z I - (F - K_{q}^{'} H_{q}))}^{- 1} \underset{―}{w} (z), and {\underset{―}{e}}_{2, q}^{-} (z) = - {(z I - (F - K_{q}^{'} H_{q}))}^{- 1} K_{q}^{'} \underset{―}{ζ} (z) \end{matrix}$ 17

The preceding two equations provide the z-transforms, ${\underset{―}{e}}_{q}^{+} (z)$ and ${\underset{―}{e}}_{q}^{-} (z)$ , of ${\underset{―}{e}}^{+}$ and ${\underset{―}{e}}^{-}$ , respectively, for the piecewise constant model $q \cdot {\underset{―}{e}}_{1, q} (z)$ and ${\underset{―}{e}}_{2, q} (z)$ represent the z-transform components with respect to $\underset{―}{w}$ and $\underset{―}{ζ}$ , respectively. The z-transforms of ${\underset{―}{ϵ}}_{e}$ and ${\underset{―}{ϵ}}_{p}$ (see after Equation (5) for a definition), ${\underset{―}{ϵ}}_{e} (z)$ and ${\underset{―}{ϵ}}_{p} (z)$ , are $H_{q} {\underset{―}{e}}_{q}^{+} (z)$ and $H_{q} {\underset{―}{e}}_{q}^{-} (z)$ , respectively. Because the z-transform is valid for a finite step size of 0.1 s, it is more appropriate for vector loops than the Laplace transform. The presented formulation does not require the decoupling of VDLL and VFLL, unlike the work of Bhattacharyya (2018), and also removes other key assumptions of the previous work. The z-domain TFs are used next to find the noise equivalent BWs and jitter.

2.2.1 Noise Equivalent Bandwidths

The TF matrices relating $\underset{―}{w} (z)$ and $\underset{―}{ζ} (z)$ to ${\underset{―}{ϵ}}_{e} (z)$ are denoted by $G_{e_{1, q}} (z)$ and $G_{e_{2, q}} (z)$ , respectively. Similarly, the matrices for $ϵ_{p} (z)$ are $G_{p_{1, q}} (z)$ and $G_{p_{2, q}} (z)$ . Using the previous two equations and the paragraph after them, the mathematical expressions of these matrices are summarized as follows:

$\begin{matrix} G_{e_{1, q}} (z) = H_{q} {(z I - (I - K_{q} H_{q}) F)}^{- 1} (I - K_{q} H_{q}), \\ G_{e_{2, q}} (z) = - z H_{q} {(z I - (I - K_{q} H_{q}) F)}^{- 1} K_{q} \end{matrix}$ 18

$\begin{matrix} G_{p_{1, q}} (z) = H_{q} {(z I - (F - K_{q}^{'} H_{q}))}^{- 1} \\ G_{p_{2, q}} (z) = - H_{q} {(z I - (F - K_{q}^{'} H_{q}))}^{- 1} K_{q}^{'} \end{matrix}$ 19

The dimensions of $G_{e_{1, q}} (z)$ and $G_{p_{1, q}} (z)$ are $2 n \times 15$ , and those of $G_{e_{2, q}} (z)$ and $G_{p_{2, q}} (z)$ are $2 n \times 2 n$ . Any arbitrary TF matrix, D(z), of dimension $m \times r$ can be written as follows:

$D (z) = [\begin{array}{l} D_{11} (z) & D_{12} (z) & \dots & D_{1 r} (z) \\ D_{21} (z) & D_{22} (z) & \dots & D_{2 r} (z) \\ ⋮ \\ D_{m 1} (z) & D_{m 2} (z) & \dots & D_{m r} (z) \end{array}]$ 20

where $D_{i η}$ is the TF element in the i-th row and η-th column. Its one-sided noise equivalent BW is given by the following:

${BW}^{i η} = \frac{1}{2 π} \int_{0}^{π} D_{i η} (e^{j ω}) {(D_{i η} (e^{j ω}))}^{*} d ω$ 21

where ${(D_{i η} (e^{j ω}))}^{*}$ is the complex conjugate transpose of $D_{i η} (e^{j ω})$ . BW ^iη is the BW of a unity gain boxcar filter that will pass the same noise power as that of $D_{i η} (e^{j ω})$ . Note that the direct current gain is not divided in the preceding equation, as different TF elements need not have the same gain (Bhattacharyya, 2018). Thus, ${BW}^{i η}$ is given by the square of the $H_{2}$ norm of $D_{i η} (e^{j ω})$ divided by 2 (i.e., ${BW}^{in} = {‖ D_{in} ‖}^{2} / 2$ ). In MATLAB, the noise equivalent BWs are calculated for a discrete step size of 0.1 s. Algorithm 1 calculates the BWs for a generic TF matrix, $D (z) = P {(z I - A_{1})}^{- 1} B_{1}$ , of dimension $m \times r$ , where $P, A_{1}$ , and $B_{1}$ are matrices of appropriate dimensions.

Because there are two input vectors (process noise, $\underset{―}{w}$ , and measurement noise, $\underset{―}{ζ}$ ), BWs are determined for each input vector using the corresponding TF matrix. These BWs are called the process noise and measurement noise BWs, respectively. The ${BW}^{in}$ terms for $G_{e_{1, q}} (z), G_{e_{2, q}} (z), G_{p_{1, q}} (z)$ , and $G_{p_{2, q}} (z)$ are denoted by ${BW}_{e_{1, q}}^{i η}, {BW}_{e_{2, q}}^{i η}, {BW}_{p_{1, q}}^{i η}$ , and ${BW}_{p_{2, q}}^{i η}$ , respectively. Next, the TF noise jitter and dynamic stress error are calculated.

ALGORITHM 1

Noise BW calculation

for $η = 1 : r %$ column index of the TF matrix, $D (z) = P {(z I - A_{1})}^{- 1} B_{1}$

% Find the zeros, poles, and gains of the TF matrix, D(z), from its state-space representation

% for the η-th input

[zero, pole, gain] = ss2zp(A₁, B₁, P, zero matrix of dimension $m \times r, η$ ) % see MATLAB help

% for details regarding “ss2zp"

for $i = 1 : m %$ row index of the TF matrix, $D (z)$

% Find the TF element, $D_{i η}$ , in the i-th row and η-th column in pole-zero form for a

% discrete step size of 0.1 s

$D_{i η} = zpk (zero (:, i)$ , pole, gain(i), 0.1); % zpk(.) finds the pole-zero form

% Calculate the noise BW of $D_{i η}$ as ${BW}^{i η} = {‖ D_{i η} ‖}^{2} / 2$

${BW}^{i η} = (norm (D_{i η}) \times norm (D_{i η})) / 2; %$ norm(.) calculates the $H_{2}$ norm of the TF, $D_{i η}$

end

2.2.2 Noise Jitter and Dynamic Stress Error

Considering both process and measurement noise BWs, the a posteriori noise jitters of a VDLL channel i and VFLL channel η for the q-th TF model are given as follows:

${\underset{―}{σ}}_{V D L L, {noise}_{t f}}^{+} (i, q) = {(2 \sum_{l = 1}^{15} {BW}_{e_{1, q}}^{i l} Q (l, l) + 2 \sum_{l = 1}^{2 n} {BW}_{e_{2, q}}^{i l} W_{q} (l, l))}^{1 / 2}, i = 1, \dots, n$ 22

${\underset{―}{σ}}_{V F L L, {noise}_{t f}}^{+} (η, q) = {(2 \sum_{l = 1}^{15} {BW}_{e_{1, q}}^{η l} Q (l, l) + 2 \sum_{l = 1}^{2 n} {BW}_{e_{2, q}}^{η l} W_{q} (l, l))}^{1 / 2}, η = (n + 1), \dots, 2 n$ 23

where Q is of dimension $15 \times 15 . W_{q}$ is the measurement error covariance of the q-th TF model obtained adaptively based on $C / N_{0}$ when the q-th model starts and is held constant until the start of the next model. It is evident that unlike scalar tracking, the noise jitter of a channel depends on the input noises of all channels. The a priori noise jitters of the same VDLL and VFLL channels are as follows:

${\underset{―}{σ}}_{V D L L, {noise}_{t f}}^{-} (i, q) = {(2 \sum_{l = 1}^{15} {BW}_{p_{1, q}}^{i l} Q (l, l) + 2 \sum_{l = 1}^{2 n} {BW}_{p_{2, q}}^{i l} W_{q} (l, l))}^{1 / 2}, i = 1, \dots, n$ 24

${\underset{―}{σ}}_{V F L L, {noise}_{t f}}^{-} (η, q) = {(2 \sum_{l = 1}^{15} {BW}_{p_{1, q}}^{η l} Q (l, l) + 2 \sum_{l = 1}^{2 n} {BW}_{p_{2, q}}^{η l} W_{q} (l, l))}^{1 / 2}, η = (n + 1), \dots, 2 n$ 25

These equations can be derived in the following way. If $D (e^{j ω})$ is a TF matrix and R is a diagonal input noise covariance, then the output noise variance of channel i is $C (i, i)$ , where C is as follows:

$C = \frac{1}{2 π} \int_{- π}^{π} D (e^{i ω}) R {(D (e^{j ω}))}^{*} d ω$ 26

The measurement noise covariance, $W_{q}$ , is assumed to be a diagonal matrix, determined based on $C / N_{0}$ and discriminator functions. The process noise covariance, Q, may have an off-diagonal element because of the receiver clock error, which is given by $2 π^{2} h_{- 2} c^{2} (Δ t)^{2} / 2$ (Brown & Hwang, 2012) $(\sim 1.77 \times 10^{- 4} m^{2} / s$ where, $h_{- 2}$ denotes the $h_{- 2}$ parameter of clock noise, having a numerical value equal to $2 \times 10^{- 20}$ for a low-quality temperature compensated crystal oscillator (TCXO) (Curran et al., 2012a), $Δ t = 0.1 s$ , and c is the speed of light). However, the effect of this term is not considered in the TF noise jitter, as it is small. The off-diagonal terms of Q for the user position, velocity, and acceleration are considered to be zero in this paper, thus allowing Equations (22)-(25) to be valid. A diagonal Q is also reasonable for an update interval of 0.1 s (Groves, 2013). Although neglecting the off-diagonal elements of Q does not compromise performance for the simulation studies in this paper, these elements will be accounted for in future work. Figure 2 depicts how input noises $\underset{―}{w}$ and $\underset{―}{ζ}$ propagate through the q-th TF model and contribute to the a priori noise jitter at the output of VFLL channel η. It should be noted that the TF noise jitter provides steady-state values. It is not possible to capture any transient effects in this method. The dynamic stress error is determined next. Unlike the noise jitter, this formulation considers transient effects. Transients can arise from changes in TF models and user motion.

FIGURE 2

Graphical representation of noise propagation from input to output through the q-th TF model for the a priori VFLL noise jitter in channel η

The superscripts of $G_{p_{1, q}}$ and $G_{p_{2, q}}$ denote the row and column indices. N(,) represents Gaussian distribution, with the mean and variance shown in parentheses.

To calculate the dynamic stress error, it is assumed that the user jerk profile (i.e., ${\underset{―}{u}}_{jerk}$ ) is known for every measurement update interval of 0.1 s and remains constant over 0.1 s. The state errors due to jerk before and after measurement updates, $Δ {\underset{―}{ϰ}}_{d y n, t f}^{-}$ and $Δ {\underset{―}{ϰ}}_{d y n, t f}^{+}$ , are propagated every 0.1 s through the q-th TF model of duration $T_{sub}$ using the following equations:

$Δ {\underset{―}{x}}_{d y n, t f, ℓ}^{+} = (I - K_{q} H_{q}) F Δ {\underset{―}{x}}_{d y n, t f, ℓ - 1}^{+} + (I - K_{q} H_{q}) B {\underset{―}{u}}_{j e r k, ℓ - 1}$ 27

$Δ {\underset{―}{x}}_{d y n, t f, ℓ}^{-} = (F - K_{q}^{'} H_{q}) Δ {\underset{―}{x}}_{d y n, t f, ℓ - 1}^{-} + B {\underset{―}{u}}_{j e r k, ℓ - 1}$ 28

Considering $B {\underset{―}{u}}_{jerk}$ to be the input vector, the TF matrices for $Δ {\underset{―}{ϰ}}_{dyn, tf}^{-}$ and $Δ {\underset{―}{ϰ}}_{dyn, tf}^{+}$ are $G_{dyn, q}^{-} (z) = {(z I - (F - K_{q}^{'} H_{q}))}^{- 1}$ and $G_{dyn, q}^{+} (z) = {(z I - (I - K_{q} H_{q}) F)}^{- 1}$ $(I - K_{q} H_{q})$ , respectively. $Δ {\underset{―}{x}}_{d y n, t f, ℓ}^{-}$ is calculated at each time instant ℓ as the sum of the impulse response of ${(z I - (F - K_{q}^{'} H_{q}))}^{- 1}$ due to $Δ {\underset{―}{x}}_{d y n, t f, ℓ - 1}^{-}$ (see the first term on the right-hand side of Equation (28)) and the step response of ${(z I - (F - K_{q}^{'} H_{q}))}^{- 1}$ over 0.1 s due to $B {\underset{―}{u}}_{j e r k, ℓ - 1}$ only. $Δ {\underset{―}{x}}_{d y n, t f, ℓ}^{+}$ is determined in a similar manner. Then, the corresponding dynamic stress errors, ${\underset{―}{ϵ}}_{dyn,tf, ℓ}^{-}$ and ${\underset{―}{ϵ}}_{dyn,tf, ℓ}^{+}$ , are obtained by multiplying $H_{q}$ by $Δ {\underset{―}{ϰ}}_{dyn,tf, ℓ}^{-}$ and $Δ {\underset{―}{ϰ}}_{dyn,tf, ℓ}^{+}$ , respectively.

The total a posteriori $1 - σ$ TF jitters of VDLL channel i and VFLL channel η after a measurement update for the q-th TF model are as follows:

${\underset{―}{σ}}_{V D L L, t f, ℓ}^{+} (i, q) = {\underset{―}{σ}}_{V D L L, {noise}_{t f}}^{+} (i, q) + \frac{| {\underset{―}{ϵ}}_{dyn,tf, ℓ}^{+} (i) |}{3}$ 29

${\underset{―}{σ}}_{V F L L, t f, ℓ}^{+} (η, q) = {\underset{―}{σ}}_{V F L L, {noise}_{t f}}^{+} (η, q) + \frac{| {\underset{―}{ϵ}}_{dyn,tf, ℓ}^{+} (η) |}{3}$ 30

Likewise, the a priori jitter can be calculated before measurement updates. Vector loop jitter thresholds and a collective metric to represent overall performance are described next. The thresholds and metric are then used in simulation studies to theoretically predict the conditions that would lead to loss of lock.

3 COLLECTIVE JITTER METRIC FOR LOSS-OF-LOCK DETECTION

In the previous formulations, the a priori and a posteriori vector loop jitters (both SS and TF) are calculated every 0.1 s, considering thermal noise and receiver clock and dynamic stress errors. In this section, a loss-of-lock detection mechanism is designed using the SS jitter described in Section 2.1. Following the work of Ward et al. (2006), the $3 - σ$ jitter threshold for loss-of-lock detection is $d / 2$ for code phase tracking, where d is the early-late correlator spacing (= 1 chip in this paper). For carrier frequency tracking, the threshold is $1 / (4 T)$ for the atan2 discriminator, where T is the coherent integration time. These values are proposed as conservative rule-of-thumb thresholds in the work by Ward et al. (2006), presumably to provide a safe margin from the noise-free pull-in range of the respective discriminators (Curran et al., 2012b; Curran, 2015). If the jitter exceeds its threshold, it is likely that the noisy discriminator output will also cross its pull-in range. It is important to note that the prescribed jitter threshold for carrier frequency tracking explicitly depends on T. In the vector loop implementation of this paper, T is 20 ms for the VDLL and 10 ms for the VFLL, as two 10-ms prompt correlator pairs are needed to form the carrier frequency discriminator outputs. Hence, for loss-of-lock detection, vector loop jitters of all channels are calculated every 20 ms so that they are obtained at the discriminator and code/carrier update rate, which is also the prediction step size of the navigation filter, and the aforementioned threshold can be applied in a manner similar to that of scalar tracking. In addition, the operations of channels are coupled in vector tracking. Thus, a collective metric accounting for the jitters of all channels is required. For this purpose, appropriate jitter bounds are proposed next.

Figure 3 illustrates the concept of jitter lower and upper bounds for a given channel. As shown in the figure, the time interval $Δ t (= 0.1 s)$ is divided into five sub-intervals of 20 ms each. The jitter is calculated every 20 ms using an approach similar to that for determining the SS jitter discussed earlier. The thick solid horizontal step represents the jitter upper bound for a channel over a 20-ms interval, whereas the dashed horizontal step indicates the lower bound. The predicted jitter is expected to be contained within the two bounds. Algorithm 2 describes the calculation of the jitter bounds every 20 ms over a 0.1-s interval.

FIGURE 3

Illustration of jitter lower and upper bounds for a channel over 20-ms intervals For a decreasing trend, the bounds will be interchanged. That is, there are two bounds: one corresponds to the beginning of a 20-ms interval and the other corresponds to the end.

The final jitter upper bound of a channel every 0.1 s is the average of five jitter upper bounds, each corresponding to the end of 20 ms, following the convention of Figure 3. Similarly, the final lower bound is the average of five lower bounds, each at the beginning of 20 ms. The reason for averaging the bounds is that if the majority of the five upper bounds over 0.1 s cross the threshold, their average is likely to cross the threshold as well. The mean of the averaged bounds across all VFLL channels and across all VDLL channels is also determined. First, the mean is separately calculated over all GPS channels and all NavIC channels, and then the two means are averaged. The averaged mean is called the mean jitter upper or lower bound for the VFLL (or VDLL). The mean bound is considered a collective jitter metric, reflecting overall vector loop performance, and is compared against the jitter threshold mentioned earlier for loss-of-lock detection every 0.1 s. The effectiveness of the mean jitter upper bound for detecting loss of lock will be shown in Section 4. It should be noted that for a decreasing jitter over 0.1 s (i.e., the reverse trend of Figure 3), the bound at the beginning of 20 ms becomes the upper bound. Because it would not be known if the jitter increases or decreases over a 0.1 s interval, both bounds at the beginning and end of 20 ms are calculated. The overall sequence of operations for calculating the SS and TF jitters and the jitter metric is described in Appendix A.6. In this paper, these terms are calculated offline. Simulation results are presented next.

ALGORITHM 2

Calculation of jitter bounds every 20 ms

% Calculation of noise jitter every 20 ms between time steps ℓ and $ℓ + 1$

Perform a measurement update of the a priori state error covariance:

% $Σ_{ℓ}^{+} = (I - K_{ℓ} H_{ℓ}) Σ_{ℓ}^{-} {(I - K_{ℓ} H_{ℓ})}^{T} + K_{ℓ} W_{ℓ} K_{ℓ}^{T}$

$Σ_{temp, 1} = Σ_{ℓ}^{+}; %$ Save $Σ_{ℓ}^{+}$ to a temporary variable for jitter upper bound calculation

$Σ_{temp, 2} = Σ_{e}^{+}; %$ Save $Σ_{e}^{+}$ to another temporary variable for jitter lower bound calculation

for index = 1:5

% Propagation of state error covariance over 20 ms for the upper bound (following the

% convention of Figure 3)

$Σ_{temp, 1} = F^{'} Σ_{temp, 1} {(F^{'})}^{T}$ ;

$Σ_{temp, 3} = Σ_{temp, 1} + (index)^{*} Q / 5 %$ One-fifth of the process noise covariance Q is added every

% 20 ms

$G_{temp} = H_{ℓ} Σ_{temp, 3} H_{ℓ}^{T}; %$ Project $Σ_{temp, 3}$ to code phase and Doppler domains

temp $= diag (G_{temp}); %$ Save the diagonal terms to an array, temp

noise_jitter_ub(:, index) = sqrt(temp); % 1-σ noise jitter upper bounds of all channels at

% the end of every 20 ms

% Propagation of state error covariance for the lower bound

$Σ_{temp, 3} = Σ_{temp, 2} + (index - 1)^{*} Q / 5;$

$Σ_{temp, 2} = F^{'} Σ_{temp, 2} {(F^{'})}^{T}; %$ For the next iteration

$G_{temp} = H_{ℓ} Σ_{temp, 3} H_{ℓ}^{T}; %$ Project $Σ_{temp, 3}$ to code phase and Doppler domains

temp $= diag (G_{temp}); %$ Save the diagonal terms to an array, temp

noise_jitter_lb(:, index) = sqrt(temp); % 1- σ noise jitter lower bounds of all channels at

% the beginning of every 20 ms

end

% Calculation of dynamic stress error every 20 ms

% Calculate $Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{+}$ and $Δ {\underset{―}{ϰ}}_{dyn, ℓ + 1}^{-}$ using Equations (11) and (12), respectively

% Calculate a posteriori and a priori dynamic stress errors at ℓ and $ℓ + 1$ , respectively

${\underset{―}{ϵ}}_{dyn, ℓ}^{+} = H_{ℓ} Δ {\underset{―}{ϰ}}_{dyn, ℓ}^{+} %$ A posteriori dynamic stress error

${\underset{―}{ϵ}}_{dyn, ℓ + 1}^{-} = H_{ℓ} Δ {\underset{―}{ϰ}}_{dyn, ℓ + 1}^{-} %$ A priori dynamic stress error; H does not change substantially

% over 0.1 s

% Assuming a linear variation of the error over 0.1 s, a slope is calculated. Note that the

% absolute value is taken:

slope $= \frac{| {\underset{―}{ϵ}}_{dyn, ℓ + 1}^{-} | / 3 - | {\underset{―}{ϵ}}_{dyn, ℓ}^{+} | / 3}{0.1};$

% Calculate the $1 - σ$ total jitters of both VDLL and VFLL channels at the beginning and end

% of every 20 ms

for index = 1:5

% 1 -σ upper bounds

% Note that the last iteration gives a jitter identical to the total a priori jitter after 0.1 s

% (i.e., Equations (13) and (14))

jitter_ub(:, index) = noise_jitter_ub(:, index) + $(| {\underset{―}{ϵ}}_{dyn, ℓ}^{+} | / 3 +$ slope ∗0.02∗ index);

% 1 -σlower bounds (Note that the first iteration gives the total a posteriori 1 -σ jitter at ℓ)

jitter_lb(:, index) = noise_jitter_lb(:, index) + $(| {\underset{―}{ϵ}}_{dyn, ℓ}^{+} | / 3 +$ slope ^∗0.02^∗( index - 1));

end

$Σ_{ℓ + 1}^{-} = Σ_{temp, 1} + Q %$ a priori state error covariance for time step $(ℓ + 1)$

4 SIMULATION STUDIES

GPS L1 and NavIC L5 signals were simulated for a high-dynamic unmanned aerial vehicle (UAV) trajectory using Accord's GNSS simulator (Accord Software and Systems, 2024). The generated signals were recorded by an IFEN intermediate-frequency (IF) data recorder with a 20-MHz sampling frequency. The recorded IF data were post-processed by a software vector receiver with the models discussed in the beginning of Section 2. The simulator introduced ionospheric and tropospheric delays using the Klobuchar and RTCA98 models, respectively. The receiver removed these errors using the Klobuchar and Hopfield models. The UAV trajectory, skyplot, and variations in jerk magnitudes are depicted in Figure 4. The jerk can be as high as 16 g/s, lasting approximately 20 s at two separate instants during the maneuver. Other jerk magnitudes are also labeled in the figure. While the signals were simulated at 42–44 dB-Hz, artificial white noise was injected into the recorded IF data to create low-C/N₀ environments. Thus, low $C / N_{0}$ along with high dynamics creates scenarios ideal for vector loop jitter analysis.

FIGURE 4

Simulated user trajectory, skyplot, and user jerk magnitude

“G” denotes GPS, and “N” denotes NavIC satellites in the skyplot. Circles indicate the end of a satellite trace.

The measurement noise covariance, W, is adaptively determined based on $C / N_{0}$ . Two sets of process noise covariance, denoted by $Q_{1}$ and $Q_{2}$ , are considered in this paper for loss-of-lock analysis with different jerks. For the user position, velocity, and acceleration states along the ECEF x-, y-, and z-axes, the covariances are $diag ([\begin{array}{l} 0.5 & 0.2 & 8 & 0.5 & 0.2 & 8 & 0.7 & 0.2 & 8 \end{array}])$ (used for all $C / N_{0}$ values) and $diag ([\begin{array}{l} 0.5 & 0.2 & 0.1 & 0.5 & 0.2 & 0.5 & 0.7 & 0.2 & 0.2 \end{array}]$ (used only for $C / N_{0}$ values of 30 dB-Hz or lower), respectively, where diag() is a diagonal matrix. The first three elements are used for the x-components of the user position, velocity, and acceleration in $m^{2}, (m / s)^{2}$ , and ${(m / s^{2})}^{2}$ , respectively; the next three are for the y-components, and the last three are for the z-components. Once set, these elements are not changed with $C / N_{0}$ or dynamics. For $Q_{1}$ , the clock bias and drift states are calculated using the $h_{0}$ and $h_{- 2}$ parameters of a TCXO (Brown & Hwang, 2012). The portion of the $Q_{1}$ matrix for the receiver clock bias (in $m^{2}$ ), drift (in $(m / s)^{2}$ ), and drift rate (in ${(m / s^{2})}^{2}$ ) (same for both GPS and NavIC) is $Q_{c l k} = [\begin{array}{c} q_{φ} Δ t + q_{ω} (Δ t)^{3} / 3 & q_{ω} (Δ t)^{2} / 2 & 0 \\ q_{ω} (Δ t)^{2} / 2 & q_{ω} (Δ t) & 0 \\ 0 & 0 & 1 \end{array}]$ . Here, $q_{φ} = \frac{h_{0}}{2} c^{2}$ and $q_{ω} = 2 π^{2} h_{- 2} c^{2}$ , the numerical values of the $h_{0}$ and $h_{- 2}$ clock parameters are $1 \times 10^{- 21}$ and $2 \times 10^{- 20}$ , respectively, and c is the speed of light. For $Q_{2}$ , at a mean $C / N_{0}$ (i.e., the average of all channel values) of 30 dB-Hz, $Q_{c l k}$ is taken as $diag ([32.51])$ to adequately track dynamics. It is not changed until 21 dB-Hz. For a mean $C / N_{0}$ at or below 20 dB-Hz, $Q_{clk}$ of $Q_{2}$ is considered $diag ([31.51])$ to stabilize the filter even without high dynamics. A larger $Q_{c l k}$ is chosen for $Q_{2}$ than that for $Q_{1}$ to compensate for small acceleration terms that are not adequate during high dynamics. It should be noted that $Q_{1}$ (except the clock bias and drift terms of $Q_{clk}$ ) and $Q_{2}$ are chosen by trial and error during the initial filter tuning process to enable the vector loops to track the user motion trajectory in Figure 4. Both sets of Q ensure that the filter learns from its measurements and that the a posteriori position, velocity, and acceleration errors are bounded by the respective $3 - σ$ bounds at steady state in the absence of user jerk across the different levels of $C / N_{0}$ considered in this paper. For $Q_{1}$ , the acceleration terms are kept high to account for propagation model inaccuracies during high dynamics with jerk, and the position and velocity terms have small, nominal values. For $Q_{2}$ , the acceleration terms are small and non-ideal. These non-ideal acceleration terms and those for the receiver clock in $Q_{2}$ result in a loss of lock at different jerks. While $Q_{1}$ is shown to better handle dynamics than $Q_{2}$ , the latter effectively demonstrates that jitter can detect loss of lock at different jerks in the same UAV simulation. The results also show how the total jitter changes with two different Q matrices, and for both cases, the loss-of-lock detection by the jitter is validated by the divergence of estimated states in the simulation.

Before loss-of-lock analysis, however, it is important to compare the SS and TF jitters and to show that jitter can predict the steady-state errors in the code and carrier channels in the presence of noise alone. To this end, the SS jitter is first validated with the TF jitter after initial transients of the filter have died down with $Q_{1}$ process noise covariance. The validation results of both the total a priori and a posteriori SS and TF jitters are shown in Table 1, where the mean $\pm 3 - σ$ values of the two jitter differences across all VDLL and VFLL channels and time steps are reported for different ranges of $C / N_{0}$ . It should be noted that the $C / N_{0}$ ratios are reduced from the default range of (42–44) dB-Hz at 126.1 s into the run and then restored at 1026.1 s. The scenario is run multiple times, each time with one set of $C / N_{0}$ noted in the table. The difference between the SS and TF jitters of a channel is expressed as a percentage of the TF jitter before the mean and standard deviation have been calculated. For a fair comparison, it is ensured that both the SS and TF jitters use the same piecewise constant models (see the beginning of Section 2.2); the jitters are then calculated theoretically but using the K, H, and W matrices obtained from the vector loop runs corresponding to each model. The jitter difference is calculated every $T_{sub}$ s (i.e., at the end of a constant model) so that the transient effects due to the model change captured by the SS jitter are reduced. The results in the table indicate that both the SS noise jitter and dynamic stress error are in good agreement with their TF counterparts, thus establishing an equivalence between two differently calculated jitters and validating the use of the SS jitter. The larger differences at lower $C / N_{0}$ can be attributed to the following factors: 1. The transition from high $C / N_{0}$ to relatively low $C / N_{0}$ at 126.1 s and vice versa at 1026.1 s results in transients, which take some time to subside. The TF noise jitter captures steady-state values, but the SS noise jitter reflects transients. 2. $T_{sub}$ is primarily 0.1 s at low $C / N_{0}$ , indicating that the TF models are changed just after 0.1 s. While changes can cause minor transients, these effects would not be evident in the TF noise jitter.

View this table:

Table 1

Mean ± 3-σ Value of the Difference between the SS and TF Jitters, Calculated over All Channels and Time Steps

The difference is expressed as a percentage of the corresponding TF jitter.

Having established an equivalence between the SS and TF jitters, the steady-state jitter performance is analyzed at different noise levels in the absence of dynamic stress error with 16 visible satellites. For this analysis, the duration of data between 330 s and 450 s into the run is considered when the UAV moves with a small constant acceleration (i.e., no jerk) and the filter has reached steady state. More precisely, it is evaluated if the jitter is able to predict the code phase and Doppler error variations. The selected duration is sufficiently long to study the steady-state performance after $C / N_{0}$ has been reduced.

For a stringent comparison, the code phase and carrier Doppler errors of a channel at the time of reception, $t_{r}$ , are determined differently from that defined after Equation (5), which forms the basis for the jitter calculation. If $ρ_{t r}$ is the true pseudorange of a channel at $t_{r}$ and $ρ_{p r}$ is the predicted pseudorange of the vector receiver, these quantities are given by the following:

$ρ_{t r} = c (t_{r} - t_{t}^{true})$ 31

$ρ_{p r} = c (t_{r} - t_{t}^{pred})$ 32

where $t_{t}^{true}$ is the true time of transmission and $t_{t}^{pred}$ is its receiver-predicted value. The difference between the two pseudoranges is as follows:

$ρ_{t r} - ρ_{p r} = c (t_{t}^{pred} - t_{t}^{true})$ 33

Because the transmitted code of a GNSS satellite gives a fine component of the signal transmission time, the code phase alignment error of the receiver results in an error in the calculation of the time. Thus, the a priori code phase error of the vector receiver in meters would be $(ρ_{t r} - ρ_{p r})$ . This error is expected to be bounded by the total a priori $3 - σ$ code phase jitter, whose $1 - σ$ value is given by Equation (13) after conversion to meters. It should be noted that the dynamic stress error component can be neglected in the absence of jerk. Similarly, the a priori version of Equation (29) should be considered for the TF jitter. The a priori carrier Doppler error of a channel is given as the difference between the receiver-predicted Doppler and true Doppler, which should be bounded by the $3 - σ$ value of the total jitter given by Equation (14). Likewise, the a priori TF jitter can be calculated following Section 2.2.

For a satellite, the true pseudorange measurement and its rate (needed for the Doppler calculation) are obtained from the simulator log files every 0.1 s, considering all error sources except the receiver clock bias and drift. Because the clock terms are not known, their estimates by the vector receiver at high $C / N_{0}$ ratios of (42-44) dB-Hz are assumed to be true and added to the measurements derived from the simulator files at appropriate times after extrapolation. The code phase and carrier Doppler errors of select GPS and NavIC satellites calculated as explained before, alongside the total a priori SS and TF jitters, are presented in Figures 5 and 6 for a mean $C / N_{0}$ of $30 dB - Hz$ . Both figures use the $Q_{1}$ process noise covariance. The code phase and carrier frequency discriminator outputs, which are the receiver's noisy estimates of the code and carrier errors, are also illustrated along with their respective $\pm 3 - σ$ bounds. Because the discriminator outputs are the innovation vector of the navigation filter, their bounds are determined from the diagonal elements of the innovation covariance matrix. In both figures, panel (c) shows the SS jitter whereas panel (d) shows the TF jitter.

Steady-state code phase jitter performance in two channels at a mean C/N0 of 30 dB-Hz Code phase discriminator outputs (O/P) are also shown along with ±3-σ bounds. Panel (c) depicts the SS jitter, and panel (d) presents the TF jitter, as obtained from their respective vector loop runs.

FIGURE 5

Steady-state code phase jitter performance in two channels at a mean C/N₀ of 30 dB-Hz

Code phase discriminator outputs (O/P) are also shown along with ±3-σ bounds. Panel (c) depicts the SS jitter, and panel (d) presents the TF jitter, as obtained from their respective vector loop runs.

Steady-state carrier Doppler jitter performance in two channels at a mean C/N0 of 30 dB-Hz Carrier frequency discriminator outputs are also shown along with ±3-σ innovation bounds. Panel (c) depicts the SS jitter, and panel (d) presents the TF jitter, as obtained from their respective vector loop runs.

FIGURE 6

Steady-state carrier Doppler jitter performance in two channels at a mean C/N₀ of 30 dB-Hz

Carrier frequency discriminator outputs are also shown along with ±3-σ innovation bounds. Panel (c) depicts the SS jitter, and panel (d) presents the TF jitter, as obtained from their respective vector loop runs.

For the results in panels (b) and (d) of Figures 5 and 6, the piecewise constant models discussed in the beginning of Section 2.2 are used because they are needed for calculating the TF jitter. For panels (a) and (c), the original, time-varying vector loop implementation is used, as calculating the SS jitter does not require piecewise constant models. Both figures show that the total a priori SS and TF jitters bound the corresponding code phase and carrier Doppler errors in their respective vector loop runs. Although not shown, other channels offer similar performance.

The metric proposed in Section 3 (i.e., the mean jitter upper bound of the VDLL/ VFLL) is expected to capture the average vector loop behavior over an interval of 0.1 s. Therefore, it is imperative to determine whether this metric can predict the code and carrier channel errors. For this purpose, Figure 7(a) presents the maximum and minimum code phase errors among all channels before and after the measurement updates (i.e., the maximum and minimum considering the a posteriori and a priori errors over a 0.1-s interval) at 30 dB-Hz and the corresponding $\pm 3 - σ$ VDLL jitter metric. Panel (b) of the same figure illustrates these errors for the carrier channels. It is evident that the jitter metric is a good indicator of the code and carrier channel performance. Finally, Figure 8 shows the carrier channel behavior for a mean $C / N_{0}$ of $23.1 dB - Hz$ for both $Q_{1}$ and $Q_{2}$ process noise covariances. While the carrier Doppler error sometimes exceeds the jitter metric in panel (a), the overall performance is acceptable. The code phase error plots, although omitted for brevity, are mostly bounded by their respective $3 - σ$ jitter metrics. It is important to note that the $3 - σ$ jitter metric represents the collective performance of all vector loop channels on average; hence, sporadic level crossings may occur at some time instants.

FIGURE 7

Comparison of ±3-σ jitter metric (or mean jitter upper bound) against code phase and carrier Doppler errors at 30 dB-Hz between 330 and 450 s into the run with Q1 process noise covariance

FIGURE 8

Comparison of ±3-σ VFLL jitter metric (or mean VFLL jitter upper bound) against carrier Doppler errors at 23 dB-Hz between 330 and 450 s into the run

The previous analysis indicates that jitter can predict the code and carrier channel errors in the presence of noise. Next, the loss-of-lock performance of the jitter is studied in the presence of jerk. For this purpose, the proposed metric with the SS jitter in Section 3 is shown in Figure 9 for the VFLL at different $C / N_{0}$ ratios. Artificial white noise is injected into the recorded IF data from 126.1 s into the run in panels (a) and (c) and from 76.1 s in panels (b) and (d) to reduce the $C / N_{0}$ from the default range of (42–44) dB-Hz. This results in an abrupt increase in the jitter before 200 s, as shown in the figure. The $C / N_{0}$ is restored after 1000 s. The mean of the reduced $C / N_{0}$ values over all channels and time steps before loss of lock is indicated at the top of each panel. It should be noted that each panel shows results for a different run executed independently of the others. The solid dark line in the figure represents the $1 - σ$ jitter threshold, $1 / (12 T)$ , for loss-of-lock detection, which is 8.33 Hz for $T = 10 ms$ . The label “State Divergence Onset” indicates the time instant when one of the vector loop states during the simulation starts to oscillate/diverge, meaning its error just crosses the $3 - σ$ bound of the filter. This is eventually followed by error growth in all states. The marker in the figure, on the other hand, indicates the time instant when the jitter metric crosses its threshold by a sufficient margin. Figures 9(a) and (b) show results for $Q_{1}$ process noise covariance whereas panels (c) and (d) show those for $Q_{2}$ . In (c), the mean jitter upper bound crosses its threshold just before the onset of the actual state error growth. In (b) and (d), the threshold is crossed at 5.2 s and 3 s after the onset, respectively. In contrast, panel (a) with $Q_{1}$ process noise covariance does not exhibit any loss of lock, as the vector loops work for the entire duration in this case. Thus, Figure 9 illustrates that the proposed jitter metric holds promise to represent the collective vector loop performance and could be a good indicator of loss of lock. Although the jitter lower bound is not shown, the time instants when the mean bound crosses the threshold are provided in the figure caption. A comparison with Figure 4(d) indicates that loss of lock occurs during 16.71 g/s in panels (b) and (c) of the figure and during 9.954 g/s in panel (d). It is evident that vector loops offer better dynamic performance with $Q_{1}$ than with $Q_{2}$ . However, the choice of $Q_{2}$ allows for the demonstration of loss of lock at two different jerk levels. It is important to note that the x-axis limits are less than 1000 s in panels (b) and (d) because the jitter is plotted until the vector loop algorithm has a minimum number of satellites for positioning. Jitter calculation is also not valid after loss of lock when the position and velocity errors are sufficiently large because the linearized model of the vector loops does not hold. Therefore, plotting the jitter after a few hundred seconds of continuous loss of lock does not provide any meaningful information, as is evident in panel (c) of the figure. In this case, after indicating loss of lock for 244.7 s, the jitter metric drops below the threshold when $C / N_{0}$ is restored to (42–44) dB-Hz, although the reduced value is larger than its previous level for the same $C / N_{0}$ . A careful inspection also reveals that after restoring $C / N_{0}$ , the individual channel jitters exhibit large variations over time in comparison with the earlier smooth profile. The actual vector loop simulation, on the other hand, has diverged significantly in 244.7 s, and in the remaining few seconds, the EKF navigation filter does not show any tendency to recover from the divergence.

Mean VFLL jitter upper bound for different mean C/N0 ratios Initial epochs with transients are not considered. Although not shown in the figure, the mean jitter lower bound (i.e., the mean of bounds calculated at the beginning of 20-ms intervals in Figure 3) exceeds the threshold for the first time at 789.3 s, 781.5 s, and 594.3 s in (b), (c), and (d), respectively. The abrupt change in jitter before 200 s and after 1000 s is due to the decrease in C/N0 from (42–44) dB-Hz and the increase afterward.

FIGURE 9

Mean VFLL jitter upper bound for different mean C/N₀ ratios

Initial epochs with transients are not considered. Although not shown in the figure, the mean jitter lower bound (i.e., the mean of bounds calculated at the beginning of 20-ms intervals in Figure 3) exceeds the threshold for the first time at 789.3 s, 781.5 s, and 594.3 s in (b), (c), and (d), respectively. The abrupt change in jitter before 200 s and after 1000 s is due to the decrease in C/N₀ from (42–44) dB-Hz and the increase afterward.

5 CONCLUSION

In this paper, the SS jitter and TF jitter of the noncoherent vector architecture were presented for a theoretical evaluation of the tracking loop performance. In this context, z-domain TF models of the time-varying vector architecture were developed without decoupling the VDLL and VFLL. Noise equivalent BWs were also calculated. The results show that the SS and TF jitters are in close agreement, thus proving their equivalence. Appropriate jitter bounds were proposed for measurement update intervals greater than 20 ms so that existing loss-of-lock thresholds could be applied. A metric to represent the collective performance of all channels was identified for loss-of-lock detection. The steady-state jitter performance in the absence of dynamic stress error was analyzed in detail. Then, high-dynamic simulations of GPS L1 and NavIC L5 signals were used to showcase the potential of the proposed metric to detect loss of lock at different jerk levels. In future work, extensive simulations with different user trajectories will be performed to assess the jitter performance. The feasibility of applying the user jerk calculated from the estimated user acceleration will be evaluated to obtain the total jitter in real time. The off-diagonal elements of the process noise covariance matrix will be accounted for in the TF jitter calculation. A better metric to represent the collective performance may also be identified. In addition, adaptive determination of the process noise covariance matrix may be implemented (Chen et al., 2024). Jitter is theoretically determined in this paper. A suitable methodology to estimate it will be developed in future work.

HOW TO CITE THIS ARTICLE:

Bhattacharyya, S. (2026). Vector tracking loop jitter. NAVIGATION, 73. https://doi.org/10.33012/navi.788

CONFLICT OF INTEREST

The author declares no conflicts of interest.

A | APPENDIX

In this section, the measurement model matrix of the vector navigation filter is defined. Then, the observability and reachability properties of some matrices are proven. Finally, the algorithms for T_sub and vector jitter calculations are described.

A.1 | Measurement Model Matrix of the Vector Navigation Filter

The measurement model matrix of the vector navigation filter at time t_l is $H_{ℓ} = [\begin{array}{l} H_{x, ℓ} & H_{y, ℓ} & H_{z, ℓ} & H_{c l k, ℓ}^{G} & H_{c l k, ℓ}^{N} \end{array}]$ . Its components are as follows:

$H_{ξ, ℓ} = {[\begin{matrix} {\tilde{u}}_{1, ℓ}^{ξ} \frac{f_{code}}{c} & \dots & {\tilde{u}}_{n_{G} + 1, ℓ}^{ξ} \frac{f_{code}}{c} & \dots & {\tilde{u}}_{n, ℓ}^{ξ} \frac{f_{code}}{c} & [0]_{1 \times n} \\ [0]_{1 \times n} & {\tilde{u}}_{1, ℓ}^{ξ} \frac{f_{carr}^{G}}{c} & \dots & {\tilde{u}}_{n, ℓ}^{ξ} \frac{f_{carr}^{N}}{c} \\ [0]_{1 \times 2 n} \end{matrix}]}_{3 \times 2 n}^{T}, ξ = x, y, z$ $H_{clk, ℓ}^{G} = {[\begin{array}{c} - \frac{f_{code}}{c} \dots - \frac{f_{code}}{c} & [0]_{1 \times n_{N}} & [0]_{1 \times n} \\ [0]_{1 \times n} & - \frac{f_{carr}^{G}}{c} \dots - \frac{f_{carr}^{G}}{c} & [0]_{1 \times n_{N}} \\ [0]_{1 \times 2 n} \end{array}]}_{3 \times 2 n}^{T}$ $H_{c l k, ℓ}^{N} = {[\begin{matrix} {[0]}_{1 \times n_{G}} & - \frac{f_{c o d e}}{c} \dots - \frac{f_{c o d e}}{c} & {[0]}_{1 \times n} \\ {[0]}_{1 \times n_{G}} & {[0]}_{1 \times n} & - \frac{f_{c a r r}^{N}}{c} \dots - \frac{f_{c a r r}^{N}}{c} \\ {[0]}_{1 \times 2 n} \end{matrix}]}_{3 \times 2 n}^{T}$

where [0] denotes a matrix or vector of zero elements. ${\tilde{u}}_{i, ℓ}^{ξ}$ is equal to the x−, y− or z−component of the line-of-sight vector in the ECEF frame at t_l for the i-th satellite, with $ξ = x, y, z . n_{G}$ is the number of visible GPS satellites, and n_N is the number of visible NavIC satellites at t_l (i.e., $n_{G} + n_{N} = n$ ). Columns 1–3 of $H_{c k, ℓ}^{G}$ and $H_{clk, ℓ}^{N}$ correspond to the receiver clock bias, drift, and drift rate for GPS and NavIC, respectively. f_code denotes the transmitted C/A-code chipping rate, which is the same for both constellations. $f_{carr}^{G}$ and $f_{carr}^{N}$ are the GPS and NavIC carrier frequencies, respectively. c is the speed of light.

Next, the observability and reachability properties of some matrices are proven. It should be noted that the subscript l denoting the time instants in the matrices is dropped for ease of representation.

A.2 | Proof that (H, (I − KH)F) is Observable

First, it is proved that (H, F(I − KH)) is observable. It has been assumed that (H, F) of the vector navigation filter is observable. Therefore, (F^T, H^T) is reachable. Using the reachability property (i.e., the Popov-Belevitch-Hautus eigenvector tests (Kailath, 1980)), it can be stated that, for any left eigenvector ${\underset{―}{v}}^{T}$ of F^T for the eigenvalue $λ, {\underset{―}{v}}^{T} F^{T} = λ {\underset{―}{v}}^{T}$ and ${\underset{―}{v}}^{T} H^{T} \neq 0$ . Let ${\underset{―}{u}}^{T}$ be a left eigenvector of (F (I − KH))^T such that ${\underset{―}{u}}^{T} H^{T} = 0$ . If the corresponding eigenvalue of (F (I − KH))^T is γ, then using ${\underset{―}{u}}^{T} H^{T} = 0$ , one can write the following:

${\underset{―}{u}}^{T} (F (I - K H))^{T} = {\underset{―}{u}}^{T} F^{T} = γ {\underset{―}{u}}^{T}$ 34

This relation implies that ${\underset{―}{u}}^{T}$ is a left eigenvector of (F^T) and ${\underset{―}{u}}^{T} H^{T} = 0$ , which contradicts the reachability criterion of (F^T, H^T) stated before. Thus, for a left eigenvector ${\underset{―}{u}}^{T}$ of (F (I − KH))^T, ${\underset{―}{u}}^{T} H^{T} \neq 0$ . That is, ((F(I-KH))^T, H^T) is reachable or (H, F(I − KH)) is observable.

Next, it is shown that (H,(I − KH)F) is observable (or (((I-KH)F)^T, H^T is reachable). Let ${\underset{―}{ω}}^{T}$ be a left eigenvector of ((I-KH)F)^T for the eigenvalue α. Then, the following relation holds:

${\underset{―}{ω}}^{T} ((I - K H) F)^{T} = α {\underset{―}{ω}}^{T}$ 35

By post-multiplying F^T on either side of the preceding equation and writing ${\underset{―}{ω}}^{T} F^{T}$ as ${\underset{―}{ω}}_{1}^{T}$ , one obtains ${\underset{―}{ω}}_{1}^{T} (I - K H)^{T} F^{T} = α {\underset{―}{ω}}_{1}^{T}$ . Thus, ${\underset{―}{ω}}_{1}^{T}$ is a left eigenvector of (F(I-KH))^T. Because ((F(I-KH))^T, H^T) has been proved to be reachable, ${\underset{―}{ω}}_{1}^{T} H^{T} \neq 0$ . Hence, $H {\underset{―}{ω}}_{1} \neq 0$ . Next, it is shown that $H \underset{―}{ω} \neq 0$ . Because ${\underset{―}{ω}}^{T}$ is a left eigenvector of ((I-KH)F)^T, proving $H \underset{―}{ω} \neq 0$ (or, ${\underset{―}{ω}}^{T} H^{T} \neq 0$ ) would imply that (((I-KH)F)^T, H^T) is reachable or (H, (I-KH)F) is observable.

From the expression of F after Equation (3), it is evident that F is invertible. Therefore, $H \underset{―}{ω} = H (F)^{- 1} {\underset{―}{ω}}_{1}$ . As $H {\underset{―}{ω}}_{1} \neq 0$ , ${\underset{―}{ω}}_{1} \in R_{H^{T}}$ (the range space of H^T), which is orthogonal to the null space of H, $N_{H}$ . Using the mathematical expression of F, it can be shown that ${\underset{―}{ω}}_{1}^{T} (F)^{- 1} {\underset{―}{ω}}_{1}$ cannot be zero for a nonzero ${\underset{―}{ω}}_{1}$ , meaning F⁻¹ will not project ${\underset{―}{ω}}_{1}$ from $R_{H^{T}}$ to $N_{H}$ . Therefore, $H \underset{―}{ω} \neq 0$ , and (H, (I-KH)F) is observable.

A.3 | Proof that ((I − KH)F, (I − KH)) is Stabilizable

Because (H, F) of the vector navigation filter is observable, there exists a K′ such that the eigenvalues of (F – K′H) are less than one (Anderson & Moore, 2005). Alternatively, it can be stated that by suitably selecting K′, one can arbitrarily place the eigenvalues of (F – K′H). Because F is invertible by definition, for each K′, there exists a unique K such that K = (F)⁻¹K′. Therefore, for any K′, (F – K′H) is identical to F(I – KH), where K = (F)⁻¹K′. It can be easily shown that the eigenvalues of (I-KH)F are identical to those of F(I – KH). Thus, the eigenvalues of F(I – KH) or (I – KH)F can be arbitrarily placed by selecting a suitable K.

Let ${\underset{―}{v}}^{T}$ be a left eigenvector of (I – KH)F such that ${\underset{―}{v}}^{T} (I - KH) = 0$ . Hence, ${\underset{―}{v}}^{T} (I - KH) F = 0$ , which implies that the corresponding eigenvalue of (I-KH)F is zero. Therefore, ${\underset{―}{v}}^{T} (I - KH) \neq 0$ for any left eigenvector, ${\underset{―}{v}}^{T}$ , of (I – KH)F for a nonzero eigenvalue. Thus, if K is selected such that the eigenvalues of (I – KH)F are nonzero, ((I-KH)F, (I – KH)) is reachable. If (I – KH)F has a zero eigenvalue, then ((I – KH)F, (I – KH)) is stabilizable, implying that any uncontrollable modes of (I – KH)F are stable. This is because the eigenvalues of (I – KH)F can be arbitrarily placed or made less than one by selecting a suitable K (see the preceding paragraph). This is true for a stable, working vector navigation filter. Although stabilizability is not required for the calculation of vector loop jitter, the proof in this section indicates that ${\underset{―}{e}}^{+}$ is affected by the process noise input in Equation (4), as the eigenvalues of (I – KH)F cannot all be zero; hence, some modes of ${\underset{―}{e}}^{+}$ would be controllable or affected by the process noise.

A.4 | Proof that ((I − KH)F, K) is Reachable

Let ${\underset{―}{v}}^{T}$ be a left eigenvector of (I − KH)F for the eigenvalue γ such that ${\underset{―}{v}}^{T} K = 0$ . Hence, ${\underset{―}{v}}^{T} F = γ {\underset{―}{v}}^{T}$ , which implies that γ is also an eigenvalue of F. Following the first paragraph of Appendix A.3, the eigenvalues of (I − KH)F can be arbitrarily placed by selecting a suitable K. For an asymptotically stable navigation filter (i.e., when vector loops are operational with a good number of satellites being tracked), K should be such that the eigenvalues of (I − KH)F are less than one. However, the eigenvalues of F are fixed and equal to one. Therefore, ${\underset{―}{v}}^{T} K \neq 0$ for any left eigenvector, ${\underset{―}{v}}^{T}$ , of (I − KH)F. Consequently, ((I − KH)F, K) is reachable.

A.5 | Algorithm for T_sub Calculation

The algorithm for calculating T_sub, over which H, W, and K matrices are assumed to be constant, is given in this section. A new model is started (i.e., the duration of T_sub for the existing TF model is determined) if one of the following happens: 1) the satellite visibility changes; 2) the 2-norm of the difference between the actual time-varying H matrix and that of the existing model exceeds a certain threshold; 3) the 2-norm of the difference between the optimal K (i.e., the K obtained using the EKF equations with the actual time-varying model of the filter) and that of the existing piecewise constant model exceeds a certain threshold. To determine whether the thrid condition arises, the optimal K is computed for each measurement update interval. The algorithm for the T_sub calculation is described next.

ALGORITHM 3

T_sub calculation

if(first measurement update)

tf_counter = 1; % TF model counter is one at the beginning

end

vis_flag =0; % Set visibility flag to zero at each epoch. If satellite visibility

% changes (i.e., a new satellite comes or a satellite sets/is blocked), a new model

% should be started. To indicate this, vis_flag is later set to one, when needed.

% Set the maximum 2-norm of the difference between the actual H and that of the

% existing TF model

delta_H_th = 0.01;

% Set the maximum 2-norm of the difference between the optimal K matrix and that of

% the existing TF model (delta_K_th) based on C/N0. The values of delta_K_th are selected

% by trial and error. It should be noted that lower values are also allowed, as the minimum

% T_sub is fixed to 0.1 s. A larger T_sub results in a TF model being valid for a longer duration,

% but performance should not degrade due to a larger T_sub.

if(average C/N0 over all channels > 30)

delta_K_th = 0.5;

else if(average C/N0 > 20)

delta_K_th = 0.1;

else

delta_K_th = 0.001;

end

if(a new satellite has come or one or more satellites have set/are blocked)

vis_flag = 1; % As a change in visibility may create transients, a new model should

% be started

end

if(vis_flag == 1 || number of measurement updates == 20) % The first 19 0.1-s intervals are

% skipped to avoid large transients

T_sub (tf_counter) = 0.1; % Start a new model

% Make K, H, and W of this model equal to the optimal K and actual H and W at the

% current time, respectively

tf_counter = tf_counter + 1; % Increment the model counter to point to the next model

else if(number of measurement updates > 20)

delta_K = Optimal K - constant K of existing model; % Find delta_K

delta_H = actual H - constant H of existing model; % Find delta_H

if(norm(delta_K) < delta_K_th and norm(delta_H) < delta_H_th)

T_sub (tf_counter-1) = T_sub (tf_counter-1) + 0.1; % Increment the duration of the

% existing model, and retain the previous K, H, and W

else

% T_sub (tf_counter) = 0.1; % Start a new model

% % Make K, H, and W of this model equal to the optimal K and actual H and W at

% % the current time, respectively

% tf_counter = tf_counter + 1; % Increment the model counter to point to the next model

end

A.6 | Algorithm for Overall Sequence of Operations for Jitter Calculation

ALGORITHM 4

Overall sequence of operations for the calculation of SS and TF jitters and jitter metric

% Obtain all H, W, and K matrices and $T_{s u b}$ values from the vector loop run

tf_index = 1; % TF model index is initialized to one

loop_count = 0.1; % Iteration counter has a step size of 0.1

for t = 0:0.1:total duration % A time loop iteration is performed every 0.1 s

% Calculation of the total a priori and a posteriori SS jitters and jitter metric for

% loss-of-lock detection

% Execute Algorithm 2 to calculate the five jitter lower and upper bounds over a 0.1-s

% interval between any two time steps l and (l+1). Note that the first of the five lower

% bounds is the total a posteriori SS jitter at l, and the last of the five upper bounds is the

% total a priori SS jitter at (l+1) for an increasing trend in jitter (see Figure 3).

% Calculate the jitter metric for loss-of-lock detection (i.e., the mean of bounds, as

% discussed in Section 3)

% Calculation of the total a priori and a posteriori TF jitters

if(loop_count == 0.1) % Execute only once for a given TF model

% Consider appropriate H, W, and K matrices for the TF model, which would last for

% T_sub(tf_index)

% Using Algorithm 1, calculate the process and measurement noise BWs for the TFs in

% Equations (18) and (19)

end

if(loop_count ∼= T_sub(tf_index)

loop_count = loop_count + 0.1; % Increment loop_count for the duration of

% T_sub)(tf_index) (≥ 0.1)

else

loop_count = 0.1; % Reset loop_count to 0.1 to re-calculate the BWs as a new TF model

% would start in the next iteration

tf_index = tf_index + 1; % Increment tf_index by one to point to the next TF model

end

% For the dynamic stress error at ℓ, determine the impulse response of

% $(z I - (I - K_{q} H_{q}) F)^{- 1}$ after 0.1 s due to $∆ {\underset{―}{x}}_{d y, t f, ℓ - 1}^{+}$ in Equation (27) (and similarly due

% to $∆ {\underset{―}{x}}_{d y, t f, ℓ - 1}^{-}$ in Equation (28)). Use the MATLAB function “impulse” for

% this purpose.

% Determine the step response of $(z I - (I - K_{q} H_{q}) F)^{- 1} (I - K_{q} H_{q})$ after 0.1 s due to

% $B u_{j e r k, ℓ - 1}$ in Equation (27) (and similarly in Equation (28)), and add this step response

% to the corresponding impulse response to determine the dynamic stress error

% Use Equations (22), (23), (29), and (30) to calculate the total a posteriori TF jitter.

% Similarly, calculate the a priori TF jitter.

end

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

REFERENCES

↵
Abedi, A. A., & Mosavi, M. R. (2022). Low computational-complexity vector tracking for low-cost GNSS receivers. Measurement, 192, 1–14. https://doi.org/10.1016/j.measurement.2022.111171
Google Scholar OpenURL
↵
Accord Software and Systems. (2024). Accord’s GNSS simulator. https://www.accord-soft.com/GNSS_simulator_simac6.php
Google Scholar
↵
Anderson, B. D. O., & Moore, J. B. (2005). Time-invariant filters. In Optimal filtering (pp. 62–89). Dover Publications.
Google Scholar
↵
Bhattacharyya, S. (2018). Vector loop transfer functions and noise bandwidths. NAVIGATION, 65(1), 55–72. https://doi.org/10.1002/navi.219
Google Scholar OpenURL
↵
Bhattacharyya, S., & Gebre-Egziabher, D. (2010). Development and validation of parametric models for vector tracking loops. NAVIGATION, 57(4), 275–295. https://doi.org/10.1002/j.2161-4296.2010.tb01783.x
Google Scholar OpenURL
↵
Brewer, J., & Raquet, J. (2016). Differential vector phase locked loop. IEEE Transactions on Aerospace and Electronic Systems, 52(3), 1046–1055. https://doi.org/10.1109/TAES.2016.140160
Google Scholar OpenURL
↵
Brown, R. G., & Hwang, P. Y. C. (2012). Kalman filter applications to the GPS and other navigation systems. In Introduction to random signals and applied Kalman filtering (pp. 324–328). John Wiley & Sons, Inc.
Google Scholar
↵
Chen, Z., Biggie, H., Ahmed, N., Julier, S., & Heckman, C. (2024). Kalman filter auto-tuning with consistent and robust Bayesian optimization. IEEE Transactions on Aerospace and Electronic Systems, 60(2), 2236–2250. https://doi.org/10.1109/TAES.2024.3350587
Google Scholar OpenURL CrossRef
↵
Copps, E. M., Geier, G. J., Fidler, W. C., & Grundy, P. A. (1980). Optimal processing of GPS signals. NAVIGATION, 27(3), 171–182. https://doi.org/10.1002/j.2161-4296.1980.tb01417.x
Google Scholar OpenURL
↵
Curran, J. T. (2015). Enhancing weak-signal carrier phase tracking in GNSS receivers. International Journal of Navigation and Observation, 2015(1), 1–15. https://doi.org/10.1155/2015/295029
Google Scholar OpenURL
↵
Curran, J. T., Lachapelle, G., & Murphy, C. C. (2012a). GNSS PLL design conditioned on thermal and oscillator phase noise. IEEE Transactions on Aerospace and Electronic Systems, 48(1), 180–196. https://doi.org/10.1109/TAES.2012.6129629
Google Scholar OpenURL
↵
Curran, J. T., Lachapelle, G., & Murphy, C. C. (2012b). Improving the design of frequency lock loops for GNSS receivers. IEEE Transactions on Aerospace and Electronic Systems, 48(1), 850–868. https://doi.org/10.1109/TAES.2012.6129674
Google Scholar OpenURL
↵
Gómez, M. A., Solera-Rico, A., Valero, E., Lázaro, J. A., & Fernández-Prades, C. (2023). Enhancing GNSS receiver performance with software-defined vector carrier tracking for rocket launching. Results in Engineering, 19, 1–15. https://doi.org/10.1016/j.rineng.2023.101310
Google Scholar OpenURL
↵
Groves, P. D. (2013). GNSS: User equipment processing and errors. In Principles of GNSS inertial and multisensor integrated navigation systems (pp. 399–417). Artech House.
Google Scholar
↵
Kailath, T. (1980). State space descriptions - some basic concepts. In Linear systems (pp. 79–145). Prentice-Hall, Inc.
Google Scholar
↵
Lashley, M., Bevly, D. M., & Hung, J. Y. (2009). Performance analysis of vector tracking algorithms for weak GPS signals in high dynamics. IEEE Journal of Selected Topics in Signal Processing, 3(4), 661–673. https://doi.org/10.1109/JSTSP.2009.2023341
Google Scholar OpenURL
↵
Liu, J., Cui, X., Chen, Q., & Lu, M. (2011, January). Joint vector tracking loop in a GNSS receiver. In Proceedings of the 2011 International Technical Meeting of the Institute of Navigation (pp. 1025–1032). https://www.ion.org/publications/pdf.cfm?articleID=9550
Google Scholar
↵
Mu, R., & Long, T. (2017). Design and implementation of vector tracking loop for high-dynamic GNSS receiver. Sensors, 21(16:5629), 1–19. https://doi.org/10.3390/s21165629
Google Scholar OpenURL
↵
Mute, D., & Bhattacharyya, S. (2023). A robust signal tracking algorithm for GPS and NAVIC constellations. GPS Solutions, 27(3), 1–20. https://doi.org/10.1007/s10291-023-01444-2
Google Scholar OpenURL
↵
Ng, Y., & Gao, G. X. (2017). GNSS multireceiver vector tracking. IEEE Transactions on Aerospace and Electronic Systems, 53(5), 2583–2593. https://doi.org/10.1109/TAES.2017.2705338
Google Scholar OpenURL
↵
Razavi, A., Gebre-Egziabher, D., & Akos, D. M. (2008). Carrier loop architectures for tracking weak GPS signals. IEEE Transactions on Aerospace and Electronic Systems, 44(2), 697–710. https://doi.org/10.1109/TAES.2008.4560215
Google Scholar OpenURL CrossRef
↵
Shytermeja, E. (2017, December). Design and performance of a GNSS single-frequency multi-constellation vector tracking architecture for urban environments [Doctoral dissertation, Université de Toulouse].
Google Scholar
↵
Simon, D. (2006). Additional topics in Kalman filtering. In Optimal state estimation: Kalman, H∞, and nonlinear approaches (pp. 298–329). Wiley-Interscience.
Google Scholar
↵
SpilkerJr., J. J. (1995). Fundamentals of signal tracking theory. In B. W.Parkinson & J. J.SpilkerJr. (Eds.), Global Positioning System: Theory and applications (Vol. 1, pp. 245–327). AIAA.
Google Scholar OpenURL
↵
Tabatabaei, A., Mosavi, M. R., Shahhoseini, H. S., & Borre, K. (2017). Vectorized and federated software receivers combining GLONASS and GPS. GPS Solutions, 21(3), 1331–1339. https://doi.org/10.1007/s10291-017-0615-8
Google Scholar OpenURL
↵
Ward, P., Betz, J. W., & Hegarty, C. J. (2006). Satellite signal acquisition, tracking, and data demodulation. In E. D.Kaplan & C. J.Hegarty (Eds.), Understanding GPS principles and applications (pp. 153–241). Artech House.
Google Scholar
↵
Weill, L. R. (2010, September). A high performance code and carrier tracking architecture for ground-based mobile GNSS receivers. In Proceedings of the 23rd International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS) (pp. 3054–3068). https://www.ion.org/publications/abstract.cfm?articleID=9413
Google Scholar
↵
Yang, R., Ling, K. V., Poh, E. K., & Morton, Y. (2017). Generalized GNSS signal carrier tracking-part II: Optimization and implementation. IEEE Transactions on Aerospace and Electronic Systems, 53(4), 1798–1811. https://doi.org/10.1109/TAES.2017.2674198
Google Scholar OpenURL

View Abstract

PDF

Previous Next

Abstract

1 INTRODUCTION

2 VECTOR LOOP JITTER

2.1 State-Space Jitter

2.2 Transfer Function Jitter

2.2.1 Noise Equivalent Bandwidths

ALGORITHM 1

2.2.2 Noise Jitter and Dynamic Stress Error

3 COLLECTIVE JITTER METRIC FOR LOSS-OF-LOCK DETECTION

ALGORITHM 2

4 SIMULATION STUDIES

5 CONCLUSION

HOW TO CITE THIS ARTICLE:

CONFLICT OF INTEREST

A | APPENDIX

A.1 | Measurement Model Matrix of the Vector Navigation Filter

A.2 | Proof that (H, (I − KH)F) is Observable

A.3 | Proof that ((I − KH)F, (I − KH)) is Stabilizable

A.4 | Proof that ((I − KH)F, K) is Reachable

A.5 | Algorithm for T_sub Calculation

ALGORITHM 3

A.6 | Algorithm for Overall Sequence of Operations for Jitter Calculation

ALGORITHM 4

REFERENCES

In this Issue

NAVIGATION: Journal of the Institute of Navigation

Jump to section

Keywords

Vector Tracking Loop Jitter

Abstract

1 INTRODUCTION

2 VECTOR LOOP JITTER

2.1 State-Space Jitter

2.2 Transfer Function Jitter

2.2.1 Noise Equivalent Bandwidths

ALGORITHM 1

2.2.2 Noise Jitter and Dynamic Stress Error

3 COLLECTIVE JITTER METRIC FOR LOSS-OF-LOCK DETECTION

ALGORITHM 2

4 SIMULATION STUDIES

5 CONCLUSION

HOW TO CITE THIS ARTICLE:

CONFLICT OF INTEREST

A | APPENDIX

A.1 | Measurement Model Matrix of the Vector Navigation Filter

A.2 | Proof that (H, (I − KH)F) is Observable

A.3 | Proof that ((I − KH)F, (I − KH)) is Stabilizable

A.4 | Proof that ((I − KH)F, K) is Reachable

A.5 | Algorithm for Tsub Calculation

ALGORITHM 3

A.6 | Algorithm for Overall Sequence of Operations for Jitter Calculation

ALGORITHM 4

REFERENCES

In this Issue

NAVIGATION: Journal of the Institute of Navigation

Alerts for this Article

Email This Article

Citation Tools

Vector Tracking Loop Jitter

Citation Manager Formats

Share This Article

Vector Tracking Loop Jitter

SAVE ITEM IN SELECTED FOLDER.

Jump to section

Related Articles

Cited By...

More in this TOC Section

Similar Articles

Keywords

A.5 | Algorithm for T_sub Calculation