Research ArticleOriginal Article
Open Access

Implementation and Accuracy of Doppler Navigation with LEO Satellites

Ariel Baron, Pini Gurfil and Hector Rotstein
NAVIGATION: Journal of the Institute of Navigation June 2024, 71 (2) navi.649; DOI: https://doi.org/10.33012/navi.649

Abstract

Utilizing broadband low Earth orbit satellite signals in an opportunistic manner for navigation is becoming increasingly popular. This paper deals with a particularly useful approach for navigation based on satellite signals of opportunity, which uses carrier Doppler-shift observables. We provide analytically derived and simplified formulas for the Jacobian involved in the numerical computation of the navigation solution and derive a global navigation satellite system-like dilution-of-precision metric that can be used to assess accuracy. A numerical study provides preliminary computational results.

Keywords

1 INTRODUCTION

Using signals from low Earth orbit (LEO) satellites as signals of opportunity (SoOps) has recently been proposed as a resilient alternative to global navigation satellite systems (GNSSs), motivated by the fact that signals from LEO satellites have desirable characteristics for positioning, navigation, and timing (PNT) applications. LEO satellites are 25–40 times closer to Earth than their counterparts in GNSS constellations, which reside in medium Earth orbits (MEOs). Reduced orbit radii give rise to two important advantages from a user viewpoint: first, signals received on or close to the Earth’s surface are potentially significantly stronger; second, the orbital periods of LEO satellites are much shorter than those of MEO satellites. Stronger signals imply stronger immunity to jamming, particularly when combined with high spectral distribution. Smaller orbital radii result in higher angular velocities with respect to the receiver, thereby providing rich Doppler measurements that can be used for navigation.

Therefore, it is not surprising that interest in LEO satellite navigation has been increasing, fueled by companies such as OneWeb, SpaceX, Boeing, Samsung, Kepler, Telesat, and LeoSat, planning to launch thousands of LEO satellites in order to provide broadband internet and communication services globally. As a consequence, signals from LEO satellites will become abundant and diverse in frequency and direction. Indeed, SoOps that might be useful for navigation will come in different frequency bands, including the L-band, S-band, Ku-band, and more. The position and velocity of LEO satellites can be approximately calculated for a given epoch by using the two-line-element (TLE) files publicly available from the North American Aerospace Defense Command (Spiridonov et al., 2020) although precise applications may require epehemerides provided by the constellation operator or dedicated ground-control segment. Consequently, exploiting broadband LEO satellite signals in an opportunistic manner for navigation is becoming an increasingly popular approach that will become mature in the near future.

The current work revisits a particularly interesting approach for LEO-satellite SoOp navigation that uses carrier Doppler-shift observables, originally proposed by Psiaki (2021a). The main results include analytically derived and simplified formulas for the Jacobian involved in the numerical computation of the navigation solution, the specification of a GNSS-like dilution-of-precision (DOP) metric that can be used to assess accuracy, and a numerical study showing preliminary results. These developments expand upon the work of Psiaki (2021a), with the DOP formulations and numerical simulations developed via a newly formulated analytical Jacobian matrix. Additionally, the current work uses numerical simulations to examine navigation performance for constellations not previously studied.

This paper is organized as follows. Section 2 contains a review of existing literature. Although LEO navigation is relatively new, several alternative approaches have already been considered. Section 3 concentrates on the algorithm for using carrier Doppler-shift measurements and shows how computations can be improved by using explicit formulas for the Jacobian instead of numerical approximations. Section 4 shows how the DOP can be formulated for the LEO navigation case. Section 5 provides results of numerical simulations of the DOP formulations, the navigation algorithm performance, and a Monte-Carlo comparison of differing options for Jacobian matrices. Section 6 presents conclusions.

2 PREVIOUS WORK

The literature on LEO navigation is already quite extensive. Some reports, such as those by Reid et al. (2018), Reid et al. (2016), and Reid et al. (2020), as well as studies on the commercial Xona Ravens constellation, have proposed using special-purpose LEO satellites by tailoring the broadband protocol and infrastructure of the constellation to support precise PNT. Others have favored the opportunistic use of LEO satellite signals designed for communication to compute a navigation solution (Ardito et al., 2019; Khalife & Kassas, 2019; Morales et al., 2018; Pesyna et al., 2012). The latter approach is followed in this paper.

Multiple SoOp methodologies exist, including time of arrival (TOA), time difference of arrival (TDOA), angle of arrival (AOA), Doppler, and differential Doppler methods (Jafarnia-Jahromi et al., 2012; Neinavaie et al., 2022; Raquet, 2013). One-way TOA methods generally require precise transmitter clocks. LEO satellite clocks, unless disciplined by GNSS, do not possess the needed clock accuracy. Moreover, satellites in LEO constellations need not be tightly time-synchronized with each other (El-Mowafy et al., 2022), at least not to the level required for positioning trilateration. We note that range estimation is easier and more accurate when additional knowledge of the transmitted signal is available, whereas in opportunistic LEO measurements, the code of the transmitted signals may not be public knowledge. “Uncooperative” signals may be successfully decoded (see, e.g., work by Humphreys et al. (2022)) and code-phase ranging made possible, under the assumption that constellation operators will not change parts of the code in the foreseeable future.

Phase ranging and TOA require integer ambiguity resolution via either a more advanced least-squares method or the use of a base station (Subirana et al., 2013; Teunissen, 2017; Zhang et al., 2020). TDOA and differential Doppler methods rely on the use of a base station, whereas AOA methods require receivers with multiple-element antennas (Al Aziz, 2015). In contrast, Doppler methods do not require integer ambiguity resolution (in the absence of a cycle slip), the use of a base station, or multiple-element antennas. Therefore, Doppler-based positioning from LEO satellites appears to be an excellent alternative to GNSS.

Humphreys et al. (2022) explored pseudorange LEO navigation by decoding the broadcasted code from STARLINK and performing code-phase ranging. While not specifically intended for LEO, the work by Morales et al. (2016) explored the use of a pseudorange SoOp to aid an inertial navigation system (INS).

Carrier phase ranging from SoOps of LEO satellites was introduced in the work by Pesyna et al. (2012), in which a standalone carrier phase ambiguity resolution of a time division multiple access signal was attempted. Alternatively, Kassas (2020) explored the use of a base station to resolve the carrier phase ambiguity from a LEO SoOp.

The combination of SoOp pseudorange and Doppler-shift measurements has also been addressed via various methodologies. For instance, Levanon (1999) explored the use of range and range rate for positioning by intersecting a sphere of constant range, a cone of constant range rate, and the Earth’s surface. This method was successfully implemented by Satelles, Inc., in their navigation platform, which claimed 20-m accuracy for a static receiver.1

Kassas (2020), Kassas et al. (2019), Morales et al. (2019), and Ardito et al. (2019) developed simultaneous tracking and navigation (STAN) with a LEO satellite signal framework, which simultaneously navigates and estimates the position of LEO satellites via the use of Doppler measurements and pseudoranges. In the STAN framework, the vehicle utilizes GNSS signals and uncertain information about LEO satellite states for navigation. The vehicle combines GNSS pseudorange mea-surements with the onboard INS while GNSS signals are available and tracks LEO satellites to refine estimates of their states. In the absence of GNSS signals, the vehicle switches to a mode wherein it simultaneously tracks LEO satellites and navigates by integrating SoOp pseudorange and Doppler measurements with its onboard INS.

The DOP associated with LEO SoOPs was initiated by McLemore & Psiaki (2022) under both Doppler shifts and pseudorange measurements. This issue was explored in order to improve the clock bias DOP (CBDOP) obtained when only Doppler shifts are used. Neinavaie et al. (2022) developed a method for using frequency difference-of-arrival SoOP signals from LEO satellites.

Although especially attractive for LEO navigation given the relatively large angular velocity between satellites and a receiver, Doppler shifts have also been proposed to enhance GNSS positioning. For instance, van Diggelen (2009) studied the use of Doppler shifts in an assisted Global Positioning System (A-GPS) context in order to jump-start an initial guess for GNSS positioning and hence reduce the time-to-first-fix under low signal-to-noise scenarios. Lehtinen (2002) and Progri et al. (2020) examined GNSS positioning with Doppler shifts for conditions in which code-pseudoranging was unavailable.

The use of Doppler navigation information as a SoOP from LEO satellites has also been explored. For instance, Khalife & Kassas (2019) examined the receiver architecture needed for Doppler information extraction and developed an extended Kalman filter for SoOp navigation from LEO Doppler shifts, whereas Benzerrouk et al. (2019) examined different estimation techniques for standalone Doppler navigation from the IRIDIUM-next constellation.

In this paper, we follow the alternative approach proposed by Psiaki (2021a), which uses a high-fidelity model of carrier Doppler shifts, to formulate a batch filter for simultaneously computing the eight unknowns of the underlying problem: three-dimensional position and velocity, clock bias, and clock drift.

3 LEO DOPPLER-SHIFT BATCH FILTER

GNSS and SoOp Doppler navigation are based on using variations of the Doppler-shift model (see, e.g., works by van Diggelen (2009), Lehtinen (2002), Khalife & Kassas (2019), Benzerrouk et al. (2019), and Morales-Ferre et al. (2020)). This model can be expressed as follows:

Embedded Image 1

Here, Dj denotes the Doppler shift from satellite j, λ is the wavelength, Embedded Image is the user-satellite unit line-of-sight (LOS) vector, vj is the satellite velocity vector, c is the speed of light, Embedded Image is the user clock drift calculated by differentiating the user clock bias δR with respect to the true received signal time TR, Embedded Image is the satellite clock drift calculated by differentiating the satellite clock bias δj with respect to signal transmission time Tj, and v is the user velocity. In this paper, unless otherwise mentioned, all position, velocity, and acceleration vectors are resolved in Earth-centered Earth-fixed (ECEF) coordinates. The standard use of this model (Angrisano et al., 2022; Petovello, 2015; Psiaki, 2021b) assumes that the position of the user is either known or has been estimated from pseudorange measurements, and the velocity is then computed. In contrast, in this paper, we will follow the approach of Psiaki (2021a), with the following eight equations in this section emanating from that reference. The current work utilizes the simultaneous dependence of the Doppler-shift observable on user position, velocity, clock bias, and clock bias rate. This approach is practical because of the large number of satellites that can be simultaneously observed in LEO constellations, together with the previously mentioned relatively large angular velocity with respect to the observer. Thus, we consider the eight-dimensional state vector:

Embedded Image 2

where r and v are the respective user position and velocity vectors in ECEF coordinates and δR and Embedded Image are the user clock bias and drift, respectively. Now, let us consider the connection between the range-rate equivalent Doppler shift (RREDS) and the derivative of the accumulated delta range (ADR):

Embedded Image 3

The ADR can be modeled by the following equation:

Embedded Image 4

where rj is the satellite position vector computed from ephemeris or TLE data, ωE is the Earth’s rotation rate, Embedded Image is the propagation time from the satellite to the receiver, Embedded Image is the signal delay due to the troposphere, Embedded Image is the signal advance due to the ionosphere, and β is the bias on the measured beat carrier phase. The time propagation Embedded Image is given as follows:

Embedded Image 5

The matrix Embedded Image, which corrects the satellite’s position at the time of reception for Earth’s rotation during signal propagation, is given by the following:

Embedded Image 6

By substituting Equation (4) into Equation (3) and neglecting the ionospheric and tropospheric terms, Psiaki (2021a) proposed a more accurate model:

Embedded Image 7

where:

Embedded Image

is the velocity of the satellite at the time of signal reception and:

Embedded Image

is a term arising from the time rate of change of the propagation delay. ωE is Earth’s rotational angular velocity vector in ECEF coordinates.

Because Equation (7) includes the effect of the satellite clock bias rate and propagation time on the relative LOS velocity vector, the model in Equation (7) is more accurate than the model in Equation (1), although Psiaki (2021a) noted that this model neglects the time derivatives of the tropospheric and ionospheric terms. These ionospheric and tropospheric error terms are functions of the receiver time-dependent user positions and clock biases, and therefore, they must be differentiated as such. As an alternative, Psiaki (2021a) proposed using a finite-difference approximation:

Embedded Image 8

and then implementing a batch filter minimizing the performance index:

Embedded Image 9

This approach is referred to herein as the Psiaki batch filter (PBF) approach. Because the cost function is nonlinear in the variables, the problem can be solved numerically, e.g., by using a recursive least-squares approach.

3.1 PBF Implementation

The PBF algorithm can be used to minimize Equation (9) by applying Algorithm 1. In Algorithm 1, a descent direction may be computed by using the Jacobian of the vector f(xk). The most straightforward way to do this is by using perturbations to compute a numerical approximation, as shown in Equation (10). Although relatively easy to implement, a numerical approximation has two main disadvantages. First, a numerical approximation may become time-consuming, especially when a relatively large number of measurements are available. Second, a careful tuning of the perturbation parameters is required to compute a reasonable approximation:

Embedded Image 10

ALGORITHM 1
ALGORITHM 1

PBF Algorithm

A more sophisticated approach would be to compute an analytic formula for the Jacobian by direct differentiation. Although this procedure requires lengthy calculations, once this step has been completed, this approach becomes efficient from a computational viewpoint, and, moreover, is highly accurate. Our subsequent numerical studies also show this approach to be superior in terms of convergence to a solution. The interested reader is referred to Appendix A for a statement and derivation of the analytical Jacobian matrix.

In practical applications, one may trade-off accuracy for computational complexity, with a preference for computing a modified descent direction by using a simplified Doppler-shift model. For instance, the following approximation to the Jacobian can be used:

Embedded Image 11

This Jacobian was derived by Psiaki (2021a) as a Jacobian matrix for the simplified Doppler-shift model given in Equation (1). Implementing this approximation is relatively straightforward when compared with the fully analytical formula or the numerical approximation. The derivative of the unit vector of geometric range Embedded Image between the satellite and receiver can be written as follows:

Embedded Image 12

We note that in Equation (12), the time derivative of the propagation delay is not considered because, during the derivation of the analytical Jacobian matrix in Appendix A.1.1.1, it was found that the magnitude of Embedded Image is on the order of 10−9 s/m. Therefore, the time propagation is considered constant for the purposes of this derivation. The acceleration Embedded Image of the j-th satellite in ECEF coordinates can be obtained as follows:

Embedded Image 13

where Embedded Image is the skew-symmetric matrix equivalent of the Earth angular velocity vector transformed from Earth-centered inertial to ECEF coordinates (Groves, 2013). Alternatively, more detailed analytical formulas for Embedded Image can be obtained from the TLE using formulas by Zhang et al. (2006) or via simple numerical differentiation. All three alternatives for Embedded Image, along with Equation (12) as an expression for Embedded Image, are possibilities for implementation of the simplified Jacobian matrix PBF.

4 DOP ANALYSIS

Given uncertain measurements, the DOP quantifies the expected accuracy of each computed variable as a function of the measurement variance. Assuming measurements with zero-mean Gaussian distributions with differing variances, the uncertainty on each PBF measurement can be modeled by λσdopp, i.e., the uncertainty in the RREDS measurements. A physically meaningful DOP can then be computed by scaling the output variables so that the matrix relating measurements to variables becomes non-dimensional (Psiaki, 2021a), thereby allowing an “apples and oranges” comparison of the different errors (e.g., position vs. velocity). To that end, the analytical expression for the Jacobian shown in Appendix A will be used. Given the state in Equation (2), the analytical Jacobian from Equation (A26) can be restated as follows:

Embedded Image 14

Once the optimization algorithm has converged, the RREDS measurement noise and state errors are related via the following:

Embedded Image 15

Now, let S be a diagonal matrix of scaling factors. We can then write the following:

Embedded Image 16

The rules for a geometric DOP (GDOP) analysis are stated below.

  1. A GDOP analysis requires that all elements of the Jacobian matrix be non-dimensional. Examination of Equation (14) leads to the following observations: the first three columns have units of Embedded Image; the fourth column has units of Embedded Image; the fifth through seventh columns are dimensionless; the eighth column has units of m/s.

  2. A GDOP analysis requires that all estimation errors have the same units and that those units be identical to the units of the measurement errors. Examination of the error state in Equation (15) shows that Δr has units of meters, ΔδR has units of seconds, Δv has units of Embedded Image, and Embedded Image has units of Embedded Image.

The above rules indicate that all error states must be scaled to have units of Embedded Image. Additionally, we must simultaneously divide the corresponding column in the Jacobian matrix in Equation (14) so that it is dimensionless. The following manipulations are identical to those performed by Psiaki (2021a).

4.1 Scaling of Embedded Image

The simplest example is the case in which Embedded Image, as it has Embedded Image units. Therefore, the last column of H in Equation (14) must be divided by c, and the error state Embedded Image must be multiplied by c.

4.2 Scaling of Δr

One approach for rescaling Δr is to multiply Δr by the upper bound on the corresponding column in the H matrix in Equation (14). If we examine the first three columns in H, we see that they are the time rate of change of the unit vector from the satellite to the receiver. The maximum of this value is called the maximum LOS sweep rate, which is given by the following (Psiaki, 2021a):

Embedded Image 17

where RE is the radius of the Earth, aorb is the semimajor axis of the satellite’s orbit, and μ is Earth’s gravitational constant. Thus, we can divide the first three columns of H by γ and multiply the error Δr by γ.

4.3 Scaling of δR

The fourth column of H in Equation (14) is the negative acceleration in the direction of the unit vector between the satellite and the receiver. This value is referred to as the range acceleration, and its maximum is given as follows (Psiaki, 2021a):

Embedded Image 18

4.4 GDOP Formulas

Consequently, we can synthesize the following relation:

Embedded Image 19

where:

Embedded Image 20

We can calculate the GDOP as follows:

Embedded Image 21

When N satellites are used, a value must be chosen for aorb; the average or minimum Embedded Image value can be chosen.

4.5 Additional DOP Definitions

In this subsection, we provide additional DOP definitions for individual components of the state. For the DOP associated with position and velocity, the values should be related to the expected accuracy in the north–east–down (NED) reference frame. Therefore, we slightly modify AGDOP from Equation (20):

Embedded Image 22

where Embedded Image is the rotation matrix from the ECEF reference frame to the NED reference frame. All DOP quantities depend on the following quantity:

Embedded Image 23

The different DOP quantities are listed below.

  1. Position DOP (PDOP):

    Embedded Image 24

  2. Horizontal Position DOP (HDOP):

    Embedded Image 25

  3. Clock Bias DOP (CBDOP):

    Embedded Image 26

  4. Velocity DOP (VDOP):

    Embedded Image 27

  5. Clock bias rate DOP (CBRDOP):

    Embedded Image 28

To calculate the expected output accuracy of a state element, the DOP value is utilized according to Table 1. Based on the relation in Table 1, there are several issues to consider with regard to the ability of obtaining a good expected accuracy of a state element (a smaller DOP is desired for better accuracy). Additionally, having lower noise for λσDopp improves accuracy. Lastly, a large η and a large γ also improve accuracy. LEO satellites have much larger values for η and γ, and therefore, their use with the PBF yields a better expected accuracy relative to higher-orbit satellites.

View this table:
TABLE 1

GDOP Expected Accuracy

5 NUMERICAL SIMULATIONS

In this section, we provide results from three simulations. The first is a DOP simulation. This simulation calculates the expected DOP and, thereby, the expected accuracy of the PBF for any constellation for any time and receiver location. The second simulation is a navigation simulation, which is also capable of propagating satellite constellations to a desired time and location. However, the navigation simulation uses the receiver time and location as well as satellite information to simulate RREDS measurements subject to noise. The simulation then uses the measurements along with the satellite orbits to perform PNT processing and obtain a navigation solution via the PBF. These first two simulations complement each other; the first calculates the expected accuracy using DOP formulations, whereas the second examines the actual PBF accuracy when RREDS measurements are subject to noise. Lastly, a Monte-Carlo performance comparison for different formulations of the PBF Jacobian matrices is performed.

5.1 DOP Simulation

In the DOP simulation, we simulated the expected DOP and accuracy for a stationary receiver located in Ramat Aviv, Israel, on June 14, 2022, at 14:59:41 for a duration of 10 min. Our code used a TLE file from the same day to propagate a multitude of constellations via the MATLAB satellite communications toolbox. Our simulation subsequently used the satellite positions and velocities and the receiver location as inputs to the AGDOP matrix from Equation (22). Subsequently, DOP values were calculated via Equations (24)–(28). We obtained the expected accuracy using equations from Table 1. We simulated a RREDS measurement error of 0.1 m/s. Each constellation had its own simulation, and we obtained DOP values for every time point during the simulation for each state member. In Table 2, we show the results of the expected position accuracy for all analyzed constellations.

View this table:
TABLE 2

DOP Simulation Results

From Table 2, it is evident that, with a measurement noise level of 0.1 m/s, STARLINK has the best expected horizontal accuracy. GPS and GNSS constellations will have poorer accuracy because of their small η and γ values. LEO satellites have much larger values for η and γ, and therefore, their expected accuracy is better.

5.2 Navigation Simulation

For our navigation simulation, we used TLEs and the MATLAB satellite communications toolbox to propagate satellite constellations to a desired time, and from those satellites that were in view of a ground-truth location, we generated RREDS measurements subject to Gaussian noise using Equation (8). Our satellite data, from a multitude of constellation and receiver ground-truth data, were used as inputs into the ADR (Equation (4)). We also provided ground-truth clock biases δR and clock bias rates Embedded Image. It should be noted that atmospheric effects were not modeled in the generated measurements. We then finite-differenced the ADR values according to Equation (8) and added 0.1-m/s random Gaussian noise to the simulated measurements. When at least eight satellites were in view of the receiver, we had the receiver use Algorithm 1 to estimate its own location, velocity, and clock data. We then compared the estimated state to the ground truth used to generate measurements. We again simulated a receiver located in Ramat Aviv, Israel, on June 14, 2022, at 14:59:41 for a duration of 10 min. This simulation helps determine how well a LEO SoOp receiver would perform.

The results of the simulation are shown in Table 3, which shows that, again, STARLINK has the best performance. These results agree with the order of magnitude of the DOP results for the expected accuracy in position shown in Table 2. We define accuracy as the root mean square (RMS) of the magnitude of the error vector. Additionally, in the absence of real signals from LEO constellations, this simulation gives a reasonable approximation of the real expected accuracy of the PBF for any constellation at any time and receiver location.

View this table:
TABLE 3

Navigation Simulation Results

5.3 Comparison of Jacobian Matrices

The numerical Jacobian in Equation (10), the analytical Jacobian in Equation (A26), and the simplified Jacobian in Equation (11) were tested using Monte-Carlo conditions as input into Algorithm 1 for the PBF for 100 trials according to Sections 5.3.1, 5.3.2, and 5.3.2. The purpose of this analysis was to investigate the feasibility of using alternative Jacobian matrices in the PBF.

5.3.1 Truth Parameter Sampling

We selected truth locations according to the same distributions used by Psiaki (2021a). Truth locations were randomly selected across the Earth’s surface, with an equal probability distribution over the entire globe. Truth altitudes were chosen randomly from a flat distribution within the range of 0–9,144 m above sea level (30,000 ft). Truth receiver clock offsets were sampled from a flat distribution between −0.25 s and +0.25 s. The individual ECEF components of truth velocities were randomly sampled from a Gaussian distribution with a mean of 0 and a standard deviation of 137 m/s, which results in a 10% probability that the truth velocity magnitude exceeds the speed of sound. Truth receiver clock offset rates were sampled from a Gaussian distribution with a mean value of 0 and a standard deviation of 3.336 × 10−9 s/s, which corresponds to a standard deviation of 1 m/s for the truth range-rate-equivalent receiver clock offset rates Embedded Image (where Embedded Image represents the offset rate).

5.3.2 Initial Solution Sampling

We selected initial solution guesses for the PBF according to the following distributions. ECEF initial position errors were sampled from a uniform distribution between 143 and 151 km. The filter was initialized with zero velocity for each ECEF component. The filter was initialized with a receiver clock offset of 0 s. The receiver was initialized with a receiver clock offset rate of 0 s/s for all 100 cases.

5.3.3 Miscellaneous Conditions

This simulation was conducted with the following miscellaneous conditions. The elevation mask used in the simulation was set to 7.5°. To model the satellite orbits, the simulation utilized TLE data obtained from June 14, 2022. The specific time point covered by the simulation was set to 14:59:41. Lastly, because the purpose of this analysis is to compare the performance of the Jacobian matrices, a noise signal with 0.1-meter-per-second (mps) 1-σ was added to the generated RREDS measurements used in the simulation.

5.3.4 Results

From Table 4, it is clear that the analytical Jacobian matrix given in Equation (A26) has the best performance in some categories. All of the Jacobian matrices had an equivalent number of unconverged trials, and the analytical and simplified Jacobian matrices required similar numbers of iterations to reach convergence. This analysis also shows that the simplified Jacobian matrix in Equation (11) can feasibly be used with the PBF in Algorithm 1, with only a small degradation in performance, in some categories, relative to the analytical Jacobian. In some cases, the simplified Jacobian matrix even outperformed the analytical Jacobian matrix. The results from Table 4 are significant because they show that there are no downsides to using the simplified Jacobian. There are no outliers in any of the simulation cases with a much greater number of required iterations, outright failure to converge, or convergence to a much less accurate solution compared with the other two Jacobians. We note that in both the navigation and DOP simulations detailed above, a measurement error of 0.1 m/s was chosen. It remains unclear whether RREDS measurements can be obtained with this level of accuracy from the main lobe or side lobe of LEO satellite transmissions at the specified 7.5° elevation mask.

View this table:
TABLE 4

Results of Jacobian Comparison

6 CONCLUSIONS

This paper provides a comprehensive analysis of Doppler navigation with LEO satellites, presenting a useful approach for navigation using satellite SoOps and demonstrating a method for assessing accuracy using a GNSS-like DOP metric. Results on numerical simulations of the DOP formulations, the navigation algorithm performance, and a Monte-Carlo comparison of differing options for Jacobian matrices are presented, which provide valuable insights into the implementation and accuracy of this approach. This paper also provides a review of existing literature and discusses the algorithm for using carrier Doppler-shift measurements. Finally, the paper concludes that the simplified Jacobian is a viable option for implementing the PBF. This work provides closer insight for researchers and practitioners into Doppler navigation with LEO satellites.

HOW TO CITE THIS ARTICLE

Baron, A., Gurfil, P., & Rotstein, H. (2024). Implementation and accuracy of Doppler navigation with LEO satellites. NAVIGATION, 71(2). https://doi.org/10.33012/navi.649

APPENDIX A DERIVATION OF AN ANALYTICAL JACOBIAN

A.1 Derivation

The Jacobian matrix is defined as the derivative of f(x) in Algorithm 1 with respect to each element of the state defined in Equation (2). In practice, because RREDS measurements do not depend on elements of the state, the Jacobian can be defined as the derivative of the measurement model defined in Equation (8) with respect to each state element. We denote the measurement model as follows:

Embedded Image A1

Embedded Image A2

Embedded Image A3

Embedded Image A4

The ADR is given in Equation (4).

A.1.1 Term 1

Substituting the finite differences from Term 1 of Equation (A1) into Equation (4), we have the following:

Embedded Image A5

Here, we define the following:

Embedded Image A6

and denote:

Embedded Image A7

A.1.1.1 Derivative of Term 1 with Respect to r

We differentiate Equation (A5) with respect to r:

Embedded Image A8

The term Embedded Image requires a numerical approximation. However, Embedded Image is on the order of 10−9 s/m whereas Embedded Image is on the order of 105 m/s; thus, their product is on the order of 10−4. Consequently, the term can be considered negligible. Moreover, Embedded Image. We have the following for the derivative:

Embedded Image

T the z-coordinate of Embedded Image will be zero. Additionally, experimental tests show that in the x and y coordinates, the product Embedded Image is on the order 102 m/s. When Embedded Image is multiplied by Embedded Image, which is on the order of 10−9 s/m, the product is on the order of 10−7. This term is deemed negligible. With these simplifications, Equation (A8) reduces as follows:

Embedded Image A9

A.1.1.2 Derivative of Term 1 with Respect to δR

We start with the following:

Embedded Image A10

The following term:

Embedded Image

must be evaluated numerically, but again, Embedded Image is on the order of 10−2 s/s whereas Embedded Image is on the order of 10−17 s/s; thus, their product is on the order of 10−19 s/s. Therefore, the term is negligible when compared with the −1 term. Embedded Image is on the order of 10−17 s/s and, when multiplied by Embedded Image, which is on the order of 102 m/s, the product is considered negligible. Noting that Embedded Image is the velocity of the satellite at time a1 and using the simplifications above, we can rewrite Equation (A10) as follows:

Embedded Image A11

We realize that the dot product of the LOS unit vector and the satellite velocity is the component of RREDS arising from the satellite velocity. We denote this as follows:

Embedded Image A12

With this notation, Equation (A10) can be simplified as follows:

Embedded Image A13

A.1.1.3 Derivative of Term 1 with Respect to v

Here, we have the following:

Embedded Image A14

Following Psiaki (2021a), we assume here that Embedded Image.

A.1.1.4 Derivative of Term 1 with Respect to Embedded Image

Again, following Psiaki (2021a), we assume that Embedded Image and approximate Embedded Image as (r – 2vΔtR) in the first term in the unit vector below. Therefore, we can write the following:

Embedded Image A15

We can see experimentally that Embedded Image is on the order of 105 m/s and Embedded Image. Their product leaves the left-hand term in Equation (A15) on the order of 105 m. Thus, the right-hand term in Equation (A15) is clearly dominant. Therefore, we can approximate Embedded Image.

A.2 Jacobian Statement

The above process is repeated for Terms 2, 3, and 4, and a five-point finite difference of results is performed:

Embedded Image A16

Consequently, the analytical Jacobian matrix for the PBF is stated in MATLAB notation:

Embedded Image A17

Embedded Image A18

Embedded Image A19

Embedded Image A20

where:

Embedded Image A21

Embedded Image A22

Embedded Image A23

Embedded Image A24

and:

Embedded Image A25

In Equations (A16)–(A24), as mentioned by Psiaki (2021a), a finite-difference interval of 0.1 s ≤ ΔtR ≤ 0.25 s can be applied. A compact version of the analytical Jacobian can be written as follows:

Embedded Image A26

Abbreviations

ADR
accumulated delta range
AOA
angle of arrival
CBDOP
clock bias dilution of precision
CBRDOP
clock bias rate dilution of precision
DOP
dilution of precision
ECEF
Earth-centered Earth-fixed
GDOP
geometric dilution of precision
GNSS
global navigation satellite system
GPS
Global Positioning System
HDOP
horizontal position dilution of precision
INS
inertial navigation system
LEO
low Earth orbit
LOS
line of sight
MEO
medium Earth orbit
mps
meters-per-second
NED
north–east–down
PDOP
position dilution of precision
PBF
Psiaki batch filter
PNT
positioning, navigation, and timing
RREDS
range-rate equivalent Doppler shift
SoOp
signal of opportunity
STAN
simultaneous tracking and navigation
TDOA
time difference of arrival
TLE
two-line element
TOA
time of arrival
VDOP
velocity dilution of precision

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

  1. ↵
    1. Al Aziz, M.
    (2015). Navigation for UAVs using signals of opportunity [PhD Thesis, University of Kansas]. https://kuscholarworks.ku.edu/handle/1808/20929?show=full
  2. ↵
    1. Angrisano, A.,
    2. Cappello, G.,
    3. Del Pizzo, S., &
    4. Gaglione, S.
    (2022). Time-differenced carrier phase technique for precise velocity estimation on an Android smartphone. Sensors, 22(21), 8514. https://doi.org/10.3390/s22218514
  3. ↵
    1. Ardito, C.,
    2. Morales, J.,
    3. Khalife, J.,
    4. Abdallah, A., &
    5. Kassas, Z.
    (2019). Performance evaluation of navigation using LEO satellite signals with periodically transmitted satellite positions. Proc. of the 2019 International Technical Meeting of the Institute of Navigation, Reston, VA, 306–318. https://doi.org/10.33012/2019.16743
  4. ↵
    1. Benzerrouk, H.,
    2. Nguyen, Q.,
    3. Xiaoxing, F.,
    4. Amrhar, A.,
    5. Rasaee, H., &
    6. Landry, R., Jr.
    (2019). LEO satellites based Doppler positioning using distributed nonlinear estimation. International Federation of Automatic Control-PapersOnLine, 52(12), 496–501. https://doi.org/10.1016/j.ifacol.2019.11.292
  5. ↵
    1. El-Mowafy, A.,
    2. Wang, K., &
    3. Allahvirdizadeh, A.
    (2022). The potential of LEO mega-constellations in aiding GNSS to enable positioning in challenging environments, Proc. of the 27th International Federation of Surveyors Congress 2022 (FIG Congress 2022), Warsaw, Poland, 1–11. https://www.fig.net/resources/proceedings/fig_proceedings/fig2022/papers/ts04g/TS04G_el-mowafy_wang_et_al_11396.pdf
  6. ↵
    1. Groves, P.
    (2013). Principles of GNSS, Inertial, and Multisensor Integrated Navigation Systems (2nd Ed.), Artech. https://ieeexplore.ieee.org/document/9101092
  7. ↵
    1. Humphreys, T. E.,
    2. Iannucci, P. A.,
    3. Komodromos, Z., &
    4. Graff, A. M.
    (2022). Signal structure of the Starlink Ku-band downlink. IEEE Transaction on Aerospace and Electronic Systems, 59(5), 6016–6030. https://doi.org/10.1109/TAES.2023.3268610
  8. ↵
    1. Jafarnia-Jahromi, A.,
    2. Broumandan, A.,
    3. Nielsen, J., &
    4. Lachapelle, G.
    (2012). GPS vulnerability to spoofing threats and a review of antispoofing techniques. International Journal of Navigation and Observation, 2012, 127072. https://doi.org/10.1155/2012/127072
  9. ↵
    1. Kassas, Z.
    (2020). Navigation from low-Earth orbit. In Y. J. Morton, F. van Diggelen, J. J. Spilker, Jr., & B. Parkinson (Eds.), Position, navigation, and timing technologies in the 21st century: Integrated satellite navigation, sensor systems, and civil applications, 1381–1412. John Wiley & Sons, Ltd. https://doi.org/10.1002/9781119458555.ch43b
  10. ↵
    1. Kassas, Z.,
    2. Morales, J., &
    3. Khalife, J.
    (2019). New-age satellite-based navigation – STAN: Simultaneous tracking and navigation with LEO satellite signals. Inside GNSS, 14(4), 56–65. https://insidegnss.com/new-age-satellite-based-navigation-stan-simultaneous-tracking-and-navigation-with-leo-satellite-signals/
  11. ↵
    1. Khalife, J., &
    2. Kassas, Z.
    (2019). Receiver design for Doppler positioning with LEO satellites. Proc. of the IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), Brighton, UK, 5506–5510. https://doi.org/10.1109/ICASSP.2019.8682554
  12. ↵
    1. Lehtinen, A.
    (2002). Doppler Positioning with GPS [Master’s Thesis, Tampere University of Technology]. Yumpu. https://www.yumpu.com/en/document/view/18205420/antti-lehtinen-doppler-positioning-with-gps-matematiikan-laitos
  13. ↵
    1. Levanon, N.
    (1999). Instant active positioning with one LEO satellite. NAVIGATION, 46(2), 87–95. https://doi.org/10.1002/j.2161-4296.1999.tb02397.x
  14. ↵
    1. McLemore, B., &
    2. Psiaki, M. L.
    (2022). DOP analysis for a LEO navigation constellation that relies on Doppler shift and pseudorange. Proc. of the 35th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2022), Denver, CO, 2503–2524. https://doi.org/10.33012/2022.18430
  15. ↵
    1. Morales, J.,
    2. Khalife, J.,
    3. Abdallah, A.,
    4. Ardito, C., &
    5. Kassas, Z.
    (2018). Inertial navigation system aiding with Orbcomm LEO satellite Doppler measurements. Proc. of the 31st International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2018), Miami, FL, 2718–2725). https://doi.org/10.33012/2018.16059
  16. ↵
    1. Morales, J.,
    2. Khalife, J.,
    3. Cruz, U., &
    4. Kassas, Z.
    (2019). Orbit modeling for simultaneous tracking and navigation using LEO satellite signals. Proc. of the 32nd International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS + 2019), Miami, FL, 2090–2099. https://doi.org/10.33012/2019.17029
  17. ↵
    1. Morales, J.,
    2. Roysdon, P. F., &
    3. Kassas, Z.
    (2016). Signals of opportunity aided inertial navigation. Proc. of the 29th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2016), Portland, OR, 1492–1501. https://doi.org/10.33012/2016.14652
  18. ↵
    1. Morales-Ferre, R.,
    2. Lohan, E. S.,
    3. Falco, G., &
    4. Falletti, E.
    (2020). GDOP-based analysis of suitability of LEO constellations for future satellite-based positioning. Proc. of the 2020 IEEE International Conference on Wireless for Space and Extreme Environments (WiSEE), Vicenza, Italy, 147–152. https://doi.org/10.1109/WiSEE44079.2020.9262624
  19. ↵
    1. Neinavaie, M.,
    2. Shadram, Z.,
    3. Kozhaya, S., &
    4. Kassas, Z. M.
    (2022). First results of differential Doppler positioning with unknown Starlink satellite signals. Proc. of the 2022 IEEE Aerospace Conference (AERO), Big Sky, MT, 1–14. https://doi.org/10.1109/AERO53065.2022.9843493
  20. ↵
    1. Pesyna, K. M.,
    2. Kassas, Z. M., &
    3. Humphreys, T. E.
    (2012). Constructing a continuous phase time history from TDMA signals for opportunistic navigation. Proc. of the 2012 IEEE/ION Position, Location and Navigation Symposium (PLANS), Myrtle Beach, SC, 1209–1220. https://doi.org/10.1109/PLANS.2012.6236977
  21. ↵
    1. Petovello, M.
    (2015). How does a GNSS receiver estimate velocity? Inside GNSS. https://insidegnss.com/how-does-a-gnss-receiver-estimate-velocity/
  22. ↵
    1. Progri, I.,
    2. Hill, J., &
    3. Michalson, W.
    (2020). A GNSS Doppler navigation algorithm. Journal of Geolocation Geoinformation and Geo-intelligence, 2020, 2020071603, 22–29. https://doi.org/10.18610/JG3.2020.071603
  23. ↵
    1. Psiaki, M. L.
    (2021a). Navigation using carrier Doppler shift from a LEO constellation: TRANSIT on steroids. NAVIGATION, 68(3), 621–641. https://doi.org/10.1002/navi.438
  24. ↵
    1. Psiaki, M. L.
    (2021b). Navigation using carrier Doppler shift from a LEO constellation: TRANSIT on steroids [Webinar]. The Institute of Navigation. https://www.ion.org/publications/webinar-psiaki.cfm
  25. ↵
    1. Raquet, J. F.
    (2013). Navigation using pseudolites beacons and signals of opportunity. NATO STO lecture series SET-197, navigation sensors and systems in GNSS degraded and denied environments, 1–18. https://api.semanticscholar.org/CorpusID:160030256
  26. ↵
    1. Reid, T.,
    2. Gunning, K.,
    3. Perkins, A., &
    4. Neish, A.
    (2020). System and method for navigation (No. US 20200379118 A1). Retrieved from https://patents.google.com/patent/US20200379118A1/en
  27. ↵
    1. Reid, T.,
    2. Neish, A.,
    3. Walter, T., &
    4. Enge, P.
    (2016). Leveraging commercial broadband LEO constellations for navigating. Proc. of the 29th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2016), Portland, OR, 2300–2314. https://doi.org/10.33012/2016.14729
  28. ↵
    1. Reid, T.,
    2. Neish, A.,
    3. Walter, T., &
    4. Enge, P.
    (2018). Broadband LEO constellations for navigation. NAVIGATION, 65(2), 205–220. https://doi.org/10.1002/navi.234
  29. ↵
    1. Spiridonov, A. A.,
    2. Saetchnikov, V. A.,
    3. Ushakov, D. V.,
    4. Cherny, V. E., &
    5. Kezik, A. G.
    (2020). Small satellite orbit determination methods based on the Doppler measurements by Belarusian State University ground station. IEEE Journal on Miniaturization for Air and Space Systems, 2(2), 59–66. https://doi.org/10.1109/JMASS.2020.3047456
  30. ↵
    1. Subirana, J.,
    2. Zornoza, J.,
    3. Hernandez-Pajares, M.,
    4. Agency, E. S., &
    5. Fletcher, K.
    (2013). GNSS data processing (No. v. 1). ESA Communications. https://books.google.co.il/books?id=RO8xngEACAAJ
  31. ↵
    1. Teunissen, P. J.
    (2017). Carrier phase integer ambiguity resolution. In P. J. Teunissen & O. Montenbruck (Eds.), Springer handbook of global navigation satellite systems, 661–685. Springer International Publishing. http://doi.org/10.1007/978-3-319-42928-1_23
  32. ↵
    1. van Diggelen, F.
    (2009). A-GPS: Assisted GPS, GNSS, and SBAS (Vol. 1). Artech. https://ieeexplore.ieee.org/document/9100796
  33. ↵
    1. Zhang, J.,
    2. Zhang, K.,
    3. Grenfell, R., &
    4. Deakin, R.
    (2006). GPS satellite velocity and acceleration determination using the broadcast ephemeris. The Journal of Navigation, 59(2), 293–305. https://doi.org/10.1017/S0373463306003638
  34. ↵
    1. Zhang, X.,
    2. Zhang, Y., &
    3. Zhu, F.
    (2020). A method of improving ambiguity fixing rate for post-processing kinematic GNSS data. Satellite Navigation, 1, 20. https://doi.org/10.1186/s43020-020-00022-y
Back to top

In this issue

Print
Download PDF
Article Alerts
Email Article
Citation Tools
Implementation and Accuracy of Doppler Navigation with LEO Satellites
Ariel Baron, Pini Gurfil, Hector Rotstein
NAVIGATION: Journal of the Institute of Navigation Jun 2024, 71 (2) navi.649; DOI: 10.33012/navi.649
Twitter logo Facebook logo Mendeley logo
Bookmark this article

Jump to section

Keywords