Category: Research & Development

Innovation: Doppler-Aided Positioning
Improving Single-Frequency RTK in the Urban Enviornment

By Mojtaba Bahrami and Marek Ziebart

A look at how Doppler measurements can be used to smooth noisy code-based pseudoranges to improve the precision of autonomous positioning as well as to improve the availability of single-frequency real-time kinematic positioning, especially in urban environments.

INNOVATION INSIGHTS by Richard Langley

WHAT DO A GPS RECEIVER, a policeman’s speed gun, a weather radar, and some medical diagnostic equipment have in common? Give up? They all make use of the Doppler effect. First proposed in 1842 by the Austrian mathematician and physicist, Christian Doppler, it is the change in the perceived frequency of a wave when the transmitter and the receiver are in relative motion.

Doppler introduced the concept in an attempt to explain the shift in the color of light from certain binary stars. Three years later, the effect was tested for sound waves by the Dutch scientist Christophorus Buys Ballot. We have all heard the Doppler shift of a train whistle or a siren with their descending tones as the train or emergency vehicle passes us. Doctors use Doppler sonography — also known as Doppler ultrasound — to provide information about the flow of blood and the movement of inner areas of the body with the moving reflectors changing the received ultrasound frequencies. Similarly, some speed guns use the Doppler effect to measure the speed of vehicles or baseballs and Doppler weather radar measures the relative velocity of particles in the air.

The beginning of the space age heralded a new application of the Doppler effect. By measuring the shift in the received frequency of the radio beacon signals transmitted by Sputnik I from a known location, scientists were able to determine the orbit of the satellite. And shortly thereafter, they determined that if the orbit of a satellite was known, then the position of a receiver could be determined from the shift. That realization led to the development of the United States Navy Navigation Satellite System, commonly known as Transit, with the first satellite being launched in 1961. Initially classified, the system was made available to civilians in 1967 and was widely used for navigation and precise positioning until it was shut down in 1996. The Soviet Union developed a similar system called Tsikada and a special military version called Parus. These systems are also assumed to be no longer in use — at least for navigation.

GPS and other global navigation satellite systems use the Doppler shift of the received carrier frequencies to determine the velocity of a moving receiver. Doppler-derived velocity is far more accurate than that obtained by simply differencing two position estimates. But GPS Doppler measurements can be used in other ways, too. In this month’s column, we look at how Doppler measurements can be used to smooth noisy code-based pseudoranges to improve the precision of autonomous positioning as well as to improve the availability of single-frequency real-time kinematic positioning, especially in urban environments.

Correction and Further Details

The first experimental Transit satellite was launched in 1959. A brief summary of subsequent launches follows:
- Transit 1A launched 17 September 1959 failed to reach orbit
- Transit 1B launched 13 April 1960 successfully
- Transit 2A launched 22 June 1960 successfully
- Transit 3A launched 30 November 1960 failed to reach orbit
- Transit 3B launched 22 February 1961 failed to deploy in correct orbit
- Transit 4A launched 29 June 1961 successfully
- Transit 4B launched 15 November 1961 successfully
- Transits 4A and 4B used the 150/400 MHz pair of frequencies and provided geodetically useful results.
- A series of Transit prototype and research satellites was launched between 1962 and 1964 with the first fully operational satellite, Transit 5-BN-2, launched on 5 December 1963.
- The first operational or Oscar-class Transit satellite, NNS O-1, was launched on 6 October 1964.
- The last pair of Transit satellites, NNS O-25 and O-31, was launched on 25 August 1988.
“Innovation” is a regular column that features discussions about recent advances in GPS technology and its applications as well as the fundamentals of GPS positioning. The column is coordinated by Richard Langley of the Department of Geodesy and Geomatics Engineering at the University of New Brunswick, who welcomes your comments and topic ideas. To contact him, email lang @ unb.ca.

Real-time kinematic (RTK) techniques enable centimeter-level, relative positioning. The technology requires expensive, dedicated, dual-frequency, geodetic-quality receivers. However, myriad industrial and engineering applications would benefit from small-size, cost-effective, single-frequency, low-power, and high-accuracy RTK satellite positioning. Can such a sensor be developed and will it deliver? If feasible, such an instrument would find many applications within urban environments — but here the barriers to success are higher. In this article, we show how some of the problems can be overcome.

Single-Frequency RTK

Low-cost single-frequency (L1) GPS receivers have attained mass-market status in the consumer industry. Notwithstanding current levels of maturity in GPS hardware and algorithms, these receivers still suffer from large positioning errors. Any positioning accuracy improvement for mass-market receivers is of great practical importance, especially for many applications demanding small size, cost-effectiveness, low power consumption, and highly accurate GPS positioning and navigation. Examples include mobile mapping technology; machine control; agriculture fertilization and yield monitoring; forestry; utility services; intelligent transportation systems; civil engineering projects; unmanned aerial vehicles; automated continuous monitoring of landslides, avalanches, ground subsidence, and river level; and monitoring deformation of built structures. Moreover, today an ever-increasing number of smartphones and handsets come equipped with a GPS receiver. In those devices, the increasing sophistication of end-user applications and refinement of map databases are continually tightening the accuracy requirements for GPS positioning.

For single-frequency users, the RTK method does appear to offer the promise of highly precise position estimates for stationary and moving receivers and can even be considered a candidate for integration within mobile handhelds. Moreover, the RTK approach is attractive because the potential of the existing national infrastructures such as Great Britain’s Ordnance Survey National GNSS Network-RTK (OSNet) infrastructure, as well as enabling technologies such as the Internet and the cellular networks, can be exploited to deliver RTK corrections and provide high-precision positioning and navigation.

The basic premise of relative (differential) positioning techniques such as RTK is that many of the sources of GNSS measurement errors including the frequency-dependent error (the ionospheric delay) are spatially correlated. By performing relative positioning between receivers, the correlated measurement errors are completely cancelled or greatly reduced, resulting in a significant increase in the positioning accuracy and precision.

Single-Frequency Challenges. Although RTK positioning is a well-established and routine technology, its effective implementation for low-cost, single-frequency L1 receivers poses many serious challenges, especially in difficult and degraded signal environments for GNSS such as urban canyons. The most serious challenge is the use of only the L1 frequency for carrier-phase integer ambiguity resolution and validation. Unfortunately, users with single-frequency capability do not have frequency diversity and many options in forming useful functions and combinations for pseudorange and carrier-phase observables. Moreover, observations from a single-frequency, low-cost receiver are typically “biased” due to the high level of multipath and/or receiver signal-tracking anomalies and also the low-cost patch antenna design that is typically used. In addition, in those receivers, measurements are typically contaminated with high levels of noise due to the low-cost hardware design compared to the high-end receivers. This makes the reliable fixing of the phase ambiguities to their correct integer values, for single-frequency users, a non-trivial problem. As a consequence, the reliability of single-frequency observations to resolve ambiguities on the fly in an operational environ
ment is relatively low compared to the use of dual-frequency observations from geodetic-quality receivers. Improving performance will be difficult, unless high-level noise and multipath can be dealt with effectively or unless ambiguity resolution techniques can be devised that are more robust and are less sensitive to the presence of biases and/or high levels of noise in the observations.

Traditionally, single-frequency RTK positioning requires long uninterrupted initialization times to obtain reliable results, and hence have a time-to-fix ambiguities constraint. Times of 10 to 25 minutes are common. Observations made at tens of continuous epochs are used to determine reliable estimates of the integer phase ambiguities. In addition, these continuous epochs must be free from cycle slips, loss of lock, and interruptions to the carrier-phase signals for enough satellites in view during the ambiguity fixing procedure. Otherwise, the ambiguity resolution will fail to fix the phase ambiguities to correct integer values. To overcome these drawbacks and be able to determine the integer phase ambiguities and thus the precise relative positions, observations made at only one epoch (single-epoch) can be used in resolving the integer phase ambiguities. This allows instantaneous RTK positioning and navigation for single-frequency users such that the problem of cycle slips, discontinuities, and loss of lock is eliminated. However, for single-frequency users, the fixing of the phase ambiguities to their correct integer values using a single epoch of observations is a non-trivial problem; indeed, it is considered the most challenging scenario for ambiguity resolution at the present time.

Instantaneous RTK positioning relies fundamentally upon the inversion of both carrier-phase measurements and code measurements (pseudoranges) and successful instantaneous ambiguity resolution. However, in this approach, the probability of fixing ambiguities to correct integer values is dominated by the relatively imprecise pseudorange measurements. This is more severe in urban areas and difficult environments where the level of noise and multipath on pseudoranges is high. This problem may be overcome partially by carrier smoothing of pseudoranges in the range/measurement domain using, for example, the Hatch filter. While carrier-phase tracking is continuous and free from cycle slips, the carrier smoothing of pseudoranges with an optimal smoothing filter window-width can effectively suppress receiver noise and short-term multipath noise on pseudoranges. However, the effectiveness of the conventional range-domain carrier-smoothing filters is limited in urban areas and difficult GNSS environments because carrier-phase measurements deteriorate easily and substantially due to blockages and foliage and suffer from phase discontinuities, cycle-slip contamination, and other measurement anomalies. This is illustrated in Figure 1. The figure shows that in a kinematic urban environment, frequent carrier-phase outages and anomalies occur, which cause frequent resets of the carrier-smoothing filter and hence carrier smoothing of pseudoranges suffers in robustness and effective continuous smoothing.

Figure 1. Satellite tracking and carrier-phase anomaly summary during the observation time-span. These data were collected in a dense urban environment in both static and kinematic mode. The superimposed red-points show epochs where carrier-phase observables are either missing or contaminated with cycle slips, loss of locks, and/or other measurement anomalies.

Doppler Frequency Shift. While carrier-phase tracking can be discontinuous in the presence of continuous pseudoranges, a receiver generates continuous Doppler-frequency-shift measurements. The Doppler measurements are immune to cycle slips. Moreover, the precision of the Doppler measurements is better than the precision of pseudoranges because the absolute multipath error of the Doppler observable is only a few centimeters. Thus, devising methods that utilize the precision of raw Doppler measurements to reduce the receiver noise and high-frequency multipath on pseudoranges may prove valuable especially in GNSS-challenged environments. Figure 2 shows an example of the availability and the precision of the receiver-generated Doppler measurements alongside the delta-range values derived from the C/A-code pseudoranges and from the L1 carrier-phase measurements. This figure also shows that frequent carrier-phase outages and anomalies occur while for every C/A-code pseudorange measurement there is a corresponding Doppler measurement available.

Figure 2. Plots of C/A-code-pseudorange-derived delta-ranges (top), L1 carrier-phase-derived delta-ranges (middle), and L1 raw receiver-generated Doppler shifts that are transformed into delta-ranges for the satellite PRN G18 during the observation time-span when it was tracked by the receiver (bottom).

Smoothing. A rich body of literature has been published exploring aspects of carrier smoothing of pseudoranges. One factor that has not received sufficient study in the literature is utilization of Doppler measurements to smooth pseudoranges and to investigate the influence of improved pseudorange accuracy on both positioning and the integer-ambiguity resolution. Utilizing the Doppler measurements to smooth pseudoranges could be a good example of an algorithm that maximally utilizes the information redundancy and diversity provided by a GPS/GNSS receiver to improve positioning accuracy. Moreover, utilizing the Doppler measurements does not require any hardware modifications to the receiver. In fact, receivers measure Doppler frequency shifts all the time as a by-product of satellite tracking.

GNSS Doppler Measurement Overview

The Doppler effect is the apparent change in the transmission frequency of the received signal and is experienced whenever there is any relative motion between the emitter and receiver of wave signals. Theoretically, the observed Doppler frequency shift, under Einstein’s Special Theory of Relativity, is approximately equal to the difference between the received and transmitted signal frequencies, which is approximately proportional to the receiver-satellite topocentric range rate.

Beat Frequency. However, the transmitted frequency is replicated locally in a GNSS receiver. Therefore, strictly speaking, the difference of the received frequency and the receiver locally generated replica of the transmitted frequency is the Doppler frequency shift that is also termed the beat frequency. If the receiver oscillator frequency is the same as the satellite oscillator frequency, the beat frequency represents the Doppler frequency shift due to the relative, line-of-sight motion between the satellite and the receiver. However, the receiver internal oscillator is far from being perfect and therefore, the receiver Doppler measurement output is the apparent Doppler frequency shift (that includes local oscillator effects). The Doppler frequency shift is also subject to satellite-oscillator frequency bias and other disturbing effects such as atmospheric effects on the signal propagation.

To estimate the range rate, a receiver typically forms an average of the delta-range by simply integrating the Doppler over a very short period of time (for example, 0.1 second) and then dividing it by the duration of the integration interval. Since the integration of frequency over time gives the phase of the signal over that time interval, the procedure continuously forms the carrier-phase observable that is the integrated Doppler over time. Therefore, Doppler frequency shift can also be estimated by time differencing carrier-phase measurements. The carrier-phase-derived Doppler is com
puted over a longer time span, leading to smoother Doppler measurements, whereas direct loop filter output is an instantaneous measure produced over a short time interval.

Doppler frequency shift is routinely used to determine the satellite or user velocity vector. Apart from velocity determination, it is worth mentioning that Doppler frequency shifts are also exploited for coarse GPS positioning. Moreover, the user velocity vector obtained from the raw Doppler frequency shift can be and has been applied by a number of researchers to instantaneous RTK applications to constrain the float solution and hence improve the integer-ambiguity-resolution success rates in kinematic surveying. In this article, a simple combination procedure of the noisy pseudorange measurements and the receiver-generated Doppler measurements is suggested and its benefits are examined.

Doppler-Smoothing Algorithm Description

Motivated both by the continual availability and the centimeter-level precision of receiver-generated (raw) Doppler measurements, even in urban canyons, a method has been introduced by the authors that utilizes the precision of raw Doppler measurements to reduce the receiver noise and high-frequency multipath on code pseudoranges. For more detail on the Doppler-smoothing technique, see Further Reading. The objective is to smooth the pseudoranges and push the accuracy of the code-based or both code- and carrier-based positioning applications in GNSS-challenged environments.

Previous work on Doppler-aided velocity/position algorithms is mainly in the position domain. In those approaches, the improvement in the quality of positioning is gained mainly by integrating the kinematic velocities and accelerations derived from the Doppler measurement in a loosely coupled extended Kalman filter or its variations such as the complementary Kalman filter. Essentially, these techniques utilize the well-known ability of the Kalman filter to use independent velocity estimates to reduce the noise of positioning solutions and improve positioning accuracy. The main difference among these position-domain filters is that different receiver dynamic models are used.

The proposed method combines centimeter-level precision receiver-generated Doppler measurements with pseudorange measurements in a combined pseudorange measurement that retains the significant information content of each.

Two-Stage Process. The proposed Doppler-smoothing process has two stages: (1) the prediction or initialization stage and (2) the filtering stage. In the prediction stage, a new estimated smoothed value of the pseudorange measurement for the Doppler-smoothing starting epoch is obtained. In this stage, for a fixed number of epochs, a set of estimated pseudoranges for the starting epoch is obtained from the subsequent pseudorange and Doppler measurements. The estimated pseudoranges are then averaged to obtain a good estimated starting point for the smoothing process. The number of epochs used in the prediction stage is the averaging window-width or Doppler-smoothing-filter length. In the filtering stage, the smoothed pseudorange profile is constructed using the estimated smoothed starting pseudorange and the integrated Doppler measurements over time. The Doppler-smoothing procedures outlined here can be performed successively epoch-by-epoch (that is, in a moving filter), where the estimated initial pseudorange (the averaged pseudorange) is updated from epoch to epoch. Alternatively, an efficient and elegant implementation of the measurement-domain Doppler-smoothing method is in terms of a Kalman filter, where it can run as a continuous process in the receiver from the first epoch (or in post-processing software, but then without the real-time advantage). This filter allows real-time operation of the Doppler-smoothing approach.

In the experiments described in this article, a short filter window-width is used. The larger the window width used in the averaging filter process, the more precise the averaged pseudorange becomes. However, this filter is also susceptible to the ionospheric divergence phenomenon because of the opposite signs of the ionospheric contribution in the pseudorange and Doppler observables. Therefore, the ionospheric divergence effect between pseudoranges and Doppler observables increases with averaging window-width and the introduced bias in the averaged pseudoranges become apparent for longer filter lengths.

Using the propagation of variance law, it can be shown that the precision of the delta-range calculated with the integrated Doppler measurements over time depends on both the Doppler-measurement epoch interval and the precision of the Doppler measurements, assuming that noise/errors on the measurements are uncorrelated.

Experimental Results

To validate the improvement in the performance and availability of single-frequency instantaneous RTK in urban areas, the proposed Doppler-aided instantaneous RTK technique has been investigated using actual GPS data collected in both static and kinematic pedestrian trials in central London. In this article, we only focus on the static results and the kinematic trial results are omitted. It is remarked, however, that the data collected in the static mode were post-processed in an epoch-by-epoch approach to simulate RTK processing.

In the static testing, GPS test data were collected with a measurement rate of 1 Hz. At the rover station, a consumer-grade receiver with a patch antenna was used. This is a single-frequency 16-channel receiver that, in addition to the C/A-code pseudoranges, is capable of logging carrier-phase measurements and raw Doppler measurements. Reference station data were obtained from the Ordnance Survey continuously operating GNSS network. Three nearby reference stations were selected that give different baseline lengths: Amersham (AMER) ≈ 38.3 kilometers away, Teddington (TEDD) ≈ 20.8 kilometers away, and Stratford (STRA) ≈ 7.1 kilometers away. In addition, a virtual reference station (VRS) was also generated in the vicinity (60 meters away) of the rover receiver.

Doppler-Smoothing. Before we present the improvement in the performance of instantaneous RTK positioning, the effect of the Doppler-smoothing of the pseudoranges in the measurement domain and comparison with carrier-phase smoothing of pseudoranges is given. To do this, we computed the C/A-code measurement errors or observed range deviations (the differences between the expected and measured pseudoranges) in the static mode (with surveyed known coordinates) using raw, Doppler-smoothed and carrier-smoothed pseudoranges. FIGURE 3a illustrates the effect of 100-second Hatch-filter carrier smoothing and FIGURE 3b shows a 100-second Doppler-smoothing of the pseudoranges for satellite PRN G28 (RINEX satellite designator) with medium-to-high elevation angle. The raw observed pseudorange deviations (in blue) are also given as reference. The quasi-sinusoidal oscillations are characteristic of multipath. Comparing the Doppler-smoothing in Figure 3b to the Hatch carrier-smoothing in Figure 3a, it can be seen that Doppler-smoothing of pseudoranges offers a modest improvement and is more robust and effective than that of the traditional Hatch filter in difficult environments.

Figure 3. Smoothed pseudorange errors (observed range deviations) using the traditional Hatch carrier-smoothing filter. Smoothing filter length in the experiments for both filters was set to 100 seconds. Satellite PRN G28 was chosen to represent a satellite at medium-to-high elevation angle.

Figure 3. Smoothed pseudorange errors (observed range deviations) using the Doppler-smoothing filter. Smoothing filter length in the experiments for both filters was set to 100 seconds. Satellite PRN G28 was chosen to represent a satellite at medium-to-high elevation angle.

Figure 4a illustrates carrier-phase Hatch-filter smoothing for low-elevation angle satellite PRN G18. In this figure, the Hatch carrier-smoothing filter reset is indicated. It can be seen that due to the frequent carrier-phase discontinuities and cycle slips, the smoothing has to be reset and restarted from the beginning and hardly reaches its full potential. In contrast, Doppler smoothing for PRN G18 shown in FIGURE 4b had few filter resets and managed effectively to smooth the very noisy pseudorange in some sections of the data.

Figure 4. Smoothed pseudorange errors (observed range deviations) and filter resets and filter length (window width) using the traditional Hatch carrier-smoothing filter. Smoothing filter length in the experiments for both filters was set to 100 seconds. Satellite PRN G18 was chosen to represent a satellite at low elevation angle as it rises from 10 to 30 degrees.

Figure 4. Smoothed pseudorange errors (observed range deviations) and filter resets and filter length (window width) using the Doppler-smoothing. Smoothing filter length in the experiments for both filters was set to 100 seconds. Satellite PRN G18 was chosen to represent a satellite at low elevation angle as it rises from 10 to 30 degrees.

Considering RTK in this analysis, we can demonstrate the increase in the success rate of the Doppler-aided integer ambiguity resolution (and hence the RTK availability) by comparison of the obtained integer ambiguity vectors from the conventional LAMBDA (Least-squares AMBiguity Decorrelation Adjustment) ambiguity resolution method using Doppler-smoothed pseudoranges with those obtained without Doppler-aiding in post-processed mode. The performance of ambiguity resolution was evaluated based on the number of epochs where the ambiguity validation passed the discrimination/ratio test. The ambiguity validation ratio test was set to the fixed critical threshold of 2.5 in all the experiments. In addition to the ratio test, the fixed solutions obtained using the fixed integer ambiguity vectors that passed the ratio test were compared against the true position of the surveyed point to make sure that indeed the correct set of integer ambiguities were estimated.

The overall performance of the single-epoch single-frequency integer ambiguity resolution obtained by the conventional LAMBDA ambiguity resolution method without Doppler-aiding is shown in Figure 5 for baselines from 60 meters up to 38 kilometers in length. In comparison, the performance of the single-epoch single-frequency integer ambiguity resolution from the LAMBDA method using Doppler-smoothed pseudoranges are shown in Figure 6 for those baselines and they are compared with integer ambiguity resolution success rates of the conventional LAMBDA ambiguity resolution method without Doppler-aiding. Figure 6 shows that using Doppler-smoothed pseudoranges enhances the probability of identifying the correct set of integer ambiguities and hence increases the success rate of the integer ambiguity resolution process in instantaneous RTK, providing higher availability. This is more evident for shorter baselines. For long baselines, the residual of satellite-ephemeris error and atmospheric-delay residuals that do not cancel in double differencing potentially limits the effectiveness of the Doppler-smoothing approach. It is well understood that those residuals for long baselines strongly degrade the performance of ambiguity resolution. Relative kinematic positioning with single frequency mass-market receivers in urban areas using VRS has also shown improvement.

Figure 5. Single-epoch single-frequency integer ambiguity resolution success rate obtained by the conventional LAMBDA ambiguity resolution method without Doppler-aiding.

Figure 6. Plots of integer ambiguity resolution success rates: single-epoch single-frequency integer ambiguity resolution success rate obtained by the conventional LAMBDA ambiguity resolution method without Doppler-aiding (in blue) and using Doppler-smoothed pseudoranges (in green).

Conclusion

In urban areas, the proposed Doppler-smoothing technique is more robust and effective than traditional carrier smoothing of pseudoranges. Static and kinematic trials confirm this technique improves the accuracy of the pseudorange-based absolute and relative positioning in urban areas characteristically by the order of 40 to 50 percent.

Doppler-smoothed pseudoranges are then used to aid the integer ambiguity resolution process to enhance the probability of identifying the correct set of integer ambiguities. This approach shows modest improvement in the ambiguity resolution success rate in instantaneous RTK where the probability of fixing ambiguities to correct integer values is dominated by the relatively imprecise pseudorange measurements.

The importance of resolving the integer ambiguities correctly must be emphasized. Therefore, devising innovative and robust methods to maximize the success rate and hence reliability and availability of single-frequency, single-epoch integer ambiguity resolution in the presence of biased and noisy observations is of great practical importance especially in GNSS-challenged environments.

Acknowledgments

The study reported in this article was funded through a United Kingdom Engineering and Physical Sciences Research Council Engineering Doctorate studentship in collaboration with the Ordnance Survey. M. Bahrami would like to thank his industrial supervisor Chris Phillips from the Ordnance Survey for his continuous encouragement and support. Professor Paul Cross is acknowledged for his valuable comments. The Ordnance Survey is acknowledged for sponsoring the project and providing detailed GIS data.

Manufacturer

The data for the trial discussed in this article were obtained from a u-blox AG AEK-4T receiver with a u-blox ANN-MS-0-005 patch antenna.

Mojtaba Bahrami is a research fellow in the Space Geodesy and Navigation Laboratory (SGNL) at University College London (UCL). He holds an engineering doctorate in space geodesy and navigation from UCL.

Marek Ziebart is a professor of space geodesy at UCL. He is the director of SGNL and vice dean for research in the Faculty of Engineering Sciences at UCL.

FURTHER READING

• Carrier Smoothing of Pseudoranges

“Optimal Hatch Filter with an Adaptive Smoothing Window Width” by B. Park, K. Sohn, and C. Kee in Journal of Navigation, Vol. 61, 2008, pp. 435–454, doi: 10.1017/S0373463308004694.

“Optimal Recursive Least-Squares Filtering of GPS Pseudorange Measurements” by A. Q. Le and P. J. G. Teunissen in VI Hotine-Marussi Symposium on Theoretical and Computational Geodesy, Wuhan, China, May 29 – June 2, 2006, Vol. 132 of the International Association of Geodesy Symposia, Springer-Verlag, Berlin and Heidelberg, 2008, Part II, pp. 166–172, doi: 10.1007/978-3-540-74584-6_26.

“The Synergism of GPS Code and Carrier Measurements” by R. Hatch in Proceedings of the 3rdInternational Geodetic Symposium on Satellite Doppler Positioning, Las Cruces, New Mexico, February 8-12, 1982, Vol. 2, pp. 1213–1231.

• Combining Pseudoranges and Carrier-phase Measurements in the Position Domain

“Position Domain Filtering and Range Domain Filtering for Carrier-smoothed-code DGNSS: An Analytical Comparison” by H. Lee, C. Rizos, and G.-I. Jee in IEE Proceedings Radar, Sonar and Navigation, Vol. 152, No. 4, August 2005, pp. 271–276, doi:10.1049/ip-rsn:20059008.

“Complementary Kalman Filter for Smoothing GPS Position with GPS Velocity” by H. Leppakoski, J. Syrjarinne, and J. Takala in Proceedings of ION GPS/GNSS 2003, the 16th International Technical Meeting of the Satellite Division of The Institute of Navigation, Portland, Oregon, September 9–
12, 2003, pp. 1201–1210.

“Precise Platform Positioning with a Single GPS Receiver” by S. B. Bisnath, T. Beran, and R. B. Langley in GPS World, Vol. 13, No. 4, April 2002, pp. 42–49.

“GPS Navigation: Combining Pseudorange with Continuous Carrier Phase Using a Kalman Filter” by P. Y. C. Hwang and R. G. Brown in Navigation, Journal of The Institute of Navigation, Vol. 37, No. 2, 1990, pp. 181–196.

• Doppler-derived Velocity Information and RTK Positioning

“Advantage of Velocity Measurements on Instantaneous RTK Positioning” by N. Kubo in GPS Solutions, Vol. 13, No. 4, 2009, pp. 271–280, doi: 10.1007/s10291-009-0120-9.

• Doppler Smoothing of Pseudoranges and RTK Positioning

Doppler-Aided Single-Frequency Real-Time Kinematic Satellite Positioning in the Urban Environment by M. Bahrami, Ph.D. dissertation, Space Geodesy and Navigation Laboratory, University College London, U.K., 2011.

“Instantaneous Doppler-Aided RTK Positioning with Single Frequency Receivers” by M. Bahrami and M. Ziebart in Proceedings of PLANS 2010, IEEE/ION Position Location and Navigation Symposium, Indian Wells, California, May 4–6, 2010, pp. 70–78, doi: 10.1109/PLANS.2010.5507202.

“Getting Back on the Sidewalk: Doppler-Aided Autonomous Positioning with Single-Frequency Mass Market Receivers in Urban Areas” by M. Bahrami in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, 22–25 September 2009, pp. 1716–1725.

• Integer Ambiguity Resolution

“GPS Ambiguity Resolution and Validation: Methodologies, Trends and Issues” by D. Kim and R. B. Langley in Proceedings of the 7th GNSS Workshop – International Symposium on GPS/GNSS, Seoul, Korea, 30 November – 2 December 2000, Tutorial/Domestic Session, pp. 213–221.

The LAMBDA Method for Integer Ambiguity Estimation: Implementation Aspects by P. de Jong and C. Tiberius. Publications of the Delft Geodetic Computing Centre, No. 12, Delft University of Technology, Delft, The Netherlands, August 1996.

“A New Way to Fix Carrier-phase Ambiguities” by P.J.G. Teunissen, P.J. de Jonge, and C.C.J.M. Tiberius in GPS World, Vol. 6, No. 4, April 1995, pp. 58–61.

“The Least-Squares Ambiguity Decorrelation Adjustment: a Method for Fast GPS Integer Ambiguity Estimation” by P.J.G. Teunissen in Journal of Geodesy, Vol. 70, No. 1–2, 1995, pp. 65–82, doi: 10.1007/BF00863419.
May 1, 2011
Innovation: GLONASS

Developing Strategies for the Future

By Yuri Urlichich, Valeriy Subbotin, Grigory Stupak, Vyacheslav Dvorkin, Alexander Povalyaev, and Sergey Karutin

A team of authors from Russian Space Systems, a key developer of navigation and geospatial technologies in the Russian aerospace industry, describes the new L3 CDMA signal to be broadcast by GLONASS-K satellites and the progress to date in developing the SDCM augmentation system.

INNOVATION INSIGHTS by Richard Langley

IT’S NO LONGER JUST A GPS WORLD. Russia’s GLONASS, or Global′naya Navigatsionaya Sputnikova Sistema, will soon have a full complement of satellites in orbit providing positioning, navigation, and timing worldwide.

The Soviet Union began development of GLONASS in 1976 just a few years after work started on GPS. The first satellite was launched in 1982 and a fully populated constellation of 24 functioning satellites was achieved in early 1996. However, due to economic difficulties following the dismantling of the Soviet Union, by 2002 the constellation had dropped to as few as seven satellites. But the Russian economy improved, and restoration of GLONASS was given high priority by the Russian government. The satellite constellation was gradually rejuvenated using primarily a new modernized spacecraft, GLONASS-M. The new design offered many improvements, including better onboard electronics, a longer lifetime, an L2 civil signal, and an improved navigation message. The GLONASS-M spacecraft still used a pressurized, hermetically sealed cylinder for the electronics, as had the earlier versions. Today, 26 functional GLONASS-M satellites are on orbit, 22 of them in service and providing usable signals, with four more having reserve status. A full constellation of 24 satellites should be available later this year with launches of several GLONASS-M satellites and the latest variant, the GLONASS-K satellite.

GLONASS-K satellites are markedly different from their predecessors. They are lighter, use an unpressurized housing (similar to that of GPS satellites), have improved clock stability, and a longer, 10-year design life. They also include, for the first time, code-division-multiple-access (CDMA) signals accompanying the legacy frequency-division-multiple-access signals. There will be two versions: GLONASS-K1 will transmit a CDMA signal on a new L3 frequency, and GLONASS-K2, in addition, will feature CDMA signals on L1 and L2 frequencies. The first GLONASS-K1 satellite was launched on February 26 and is now undergoing tests.

GLONASS is being further improved with a satellite-based augmentation system. Called the System for Differential Correction and Monitoring or SDCM, it will use a ground network of monitoring stations and Luch geostationary communication satellites to transmit correction and integrity data using the GPS L1 frequency. The first of these satellites, Luch-5A, will be launched this year.

In this month’s column, a team of authors from Russian Space Systems, a key developer of navigation and geospatial technologies in the Russian aerospace industry, describes the new L3 CDMA signal to be broadcast by GLONASS-K satellites and the progress to date in developing the SDCM augmentation system.

The Russian Global Navigation Satellite System (GLONASS) is once again approaching full operation. As of March, 22 satellites are in service, providing nearly continuous global coverage. These satellites are modernized GLONASS or GLONASS-M satellites, transmitting the legacy frequency-domain-multiple-access (FDMA) navigation signals in the L1 and L2 frequency bands.

The structure of the navigation signals transmitted by the satellites determines the accuracy of the pseudorange measurements, which, in turn, affects a user’s position accuracy. Evolution of the GLONASS navigation signals is a top priority for the overall system development. A new version of the satellites, GLONASS-K, will broadcast a code-division-multiple-access (CDMA) signal in the L3 band for the first time in the system’s history. In addition to the change in signal parameters, new navigation information will be transmitted to users through this signal. Further GLONASS navigation signal development assumes that a new CDMA civil signal will also become available in the L1 and L2 bands.

The evolution of GNSS augmentation is also an important task in the development of satellite navigation in Russia. The Russian satellite-based augmentation system (SBAS), the System for Differential Correction and Monitoring (SDCM), is entering into the deployment phase and that is why some aspects of interoperability and compatibility with other SBASs become important. Taking into account the fact that satellite channels are the most efficient and universal tool to supply GNSS users with precise ephemeris and clock parameters and the positive experience of regional systems (such as the Quasi-Zenith Satellite System), we can see the potential for the development of a precise positioning service.

In this article, we will discuss plans for modernizing GLONASS, covering both the new signals and the augmentation service.

Navigation Signals

The main task for GLONASS development is an extension of the ensemble of navigation signals. This extension means that new CDMA signals in the L1, L2, and L3 bands will be added to the existing FDMA signals. The GLONASS satellites will keep broadcasting the legacy signals until the last receiver stops working.

The first phase in the implementation of CDMA technology on GLONASS-K satellites includes a new signal in the L3 band on a carrier frequency of 1202.025 MHz. The first GLONASS-K satellite was launched on February 26, 2011, and is undergoing tests. The ranging code chipping rate for the CDMA signal is 10.23 megachips per second with a period of 1 milliseconds. It is modulated onto the carrier using quadrature phase-shift keying (QPSK), with an in-phase data channel and a quadrature pilot channel. The signal spectrum is shown in Figure 1.

Figure 1. L3 CDMA signal spectrum (frequencies in MHz).

A block diagram of how the GLONASS L3 signal is formed is presented in Figure 2. The set of possible ranging codes consists of 31 truncated Kasami sequences. (Kasami sequences are binary sequences of length 2^m – 1 where m is an even integer. These sequences have good cross-correlation values approaching a theoretical lower bound. The Gold codes used in GPS are a special case of Kasami codes.) The full length of these sequences is 2¹⁴ – 1 = 16,383 symbols, but the ranging code is truncated to a length of N = 10,230 with a period of 1 milliseconds and with the following initial state (IS) in the generator (G) registers: G2 – IS = 00110100111000, G1 IS = n, G3 IS = n + 32. It these equations, n is the system number of the satellite in the orbit constellation. For these codes, inter-channel jamming is about –40 dB.

Figure 2. Formulation of L3 CDMA signal.

The navigation message symbols (NSs) are transmitted at a rate of 100 bits per second with half-rate convolution coding (CC) with a memory of 6. This means that the duration of an NS is 10 milliseconds and the duration of the CC symbols is 5 milliseconds. The CC switch (see Figure 2) should be in the lower position for the first half of each NS.

The pseudorandom sequence of the L3 data signal, PRS-D, is modulo-2 summed with a periodic 5-bit Barker code (BC = 00010) b
efore phase modulation. Barker code symbols have a duration of 1 millisecond and are synchronized with the pseudorandom code symbols. The pseudorandom sequence of the L3 pilot signal, PRS-P, is modulo-2 summed with a 10-bit Neuman-Hoffman code (NH = 0000110101). The Neuman-Hofman code symbols have a duration of 1 millisecond and are synchronized with the information symbols. The Barker and Neuman-Hoffman codes are used for CC synchronization in the L3 user’s receiver (see Further Reading for background details).

The navigation message superframe (2 minutes long) will consist of 8 navigation frames (NFs) for 24 regular satellites in the GLONASS first modernization stage and 10 NFs (lasting 2.5 minutes) for 30 satellites in the future. Each NF (15 seconds long) includes 5 strings (3 seconds each). Every NF has a full set of ephemerides for the current satellite and part of the system almanac for three satellites. The full system almanac is broadcast in one superframe. A time marker is located at the beginning of a string and given as a number of a string within the current day in the satellite time scale.

The GLONASS system and the satellites’ time scales are coordinated with the Russian national time scale, UTC(SU), which is periodically adjusted for a leap seconds. A special flag, A, is used in each frame to inform users about an anomalous fifth string of this frame. If А = 0, the fifth string will be normal with a 3-second duration; if А = 1, the fifth string will be either 2 seconds or 4 seconds. The correction value (+1 second or –1 second) is also transmitted in the special NF flag, KP. If KP = 11, the fifth string will be shorter due to a correction of –1 second; if KP = 01, it will be longer due to a correction of +1 second. A user should not use the short string. A string is lengthened by adding “0” to the normal string. This algorithm is implemented with the objective of simplifying the time scale correction process in user equipment.

Modulation and Multiplexing. There are intensive studies being carried out for developing new CDMA signals in the L1 and L2 bands in addition to the L3 signal described above. The main difficulties to be overcome in these studies are to ensure a low-power spectral density (PSD) of –238 dBW/m²/Hz in the 1610.6–1613.8 MHz radio astronomy band and the multiplexing of more than two signal components, providing a constant signal level.

The first task could be solved by using a modulation with a low PSD level in the radio astronomy band, such as a binary offset carrier (BOC) modulation with a subcarrier frequency of 5.115 MHz and a spreading code chipping rate of 2.5575 megachips per second (BOC(5, 2.5)) as shown in Figure 3.

Figure 3. BOC(5, 2.5) signal spectrum (frequencies in MHz).

There are two well-known methods of signal multiplexing — time multiplexing and amplitude equalizing. The time multiplexing technique is used for the GPS L2C signal, while the amplitude equalizing method is used for the composite BOC (CBOC) signals in the Galileo E1/L1 band and the alternative BOC (AltBOC) signals in E5a-E5b bands. This method has the disadvantage of about 10–16 percent loss of the transmitter power on the equalization. However, it has an advantage: simple user equipment architecture and, more importantly, the possibility of step-wise implementation of the multicomponent signal. The step-wise approach is compatible with older receivers. New user equipment will be able to track both old and new signal components, as well as a combined signal consisting of old and new components. Vector and phase diagrams for two-, four-, and six-component AltBOC signals are shown in Figure 4. Even with six components, losses are lower than about 16 percent, but it is possible to avoid any loss using time multiplexing. That is why the final decision about future GLONASS signals has not yet been made.

Figure 4. Vector and phase relationships for AltBOC signals with (a) 2-, (b) 4-, and (c) 6-components, with losses of 0, and about 15 and 16 percent respectively.

There have been extensive studies on the definition of the ensemble of code sequences with a minimum level of interchannel jamming. It was found that the jamming level for random shifts does not depend on the code type, but rather depends on the number of symbols, N, in a period. Cross-correlation functions for Kasami 4095 and Weil 10230 codes are shown in Figures 5 and 6. (Kasami codes, as previously mentioned, are being used for the GLONASS L3 CDMA signal. Weil codes are prime length sequences constructed from the well-known Legendre sequences and are used for the GPS L1C signal.) For comparison, we show cross-correlation functions for random codes with equal lengths on the same figures. It is obvious that the histograms of predefined and random codes are close to being equal. Sidelobe dispersion levels are lower than 0.1 dB.

The results obtained from the studies allow us to draw a conclusion about the invariance of the stochastic characteristics of inter-channel interference using a code structure with a fixed length of N symbols. That is why it is possible to choose an ensemble of binary code sequences on the basis of generation simplicity.

Figure 5. Kasami and random code cross-correlation functions (4,095 symbols).

Figure 6. Weil and random code cross-correlation functions (10,230 symbols).

GLONASS Augmentation Development

SDCM has been under development since 2002. The main elements of the system, including the network of reference stations in Russia and abroad, the central processing facility (CPF), and the SDCM information distribution channel, have been designed.

Ground Stations. The SDCM uses 14 monitor stations in Russia and two in Antarctica at Russia’s Bellingshausen and Novolazarevskaya research stations. Eight more monitor stations will be added in Russia and several more outside Russia. The additional overseas stations may include sites in Latin America and the Asia-Pacific region.

Central Processing. Raw measurements (GLONASS and GPS L1 and L2 pseudorange and carrier-phase measurements) from the ground stations come to the SDCM CPF. The CPF calculates the precise satellite ephemerides and clocks, controls integrity, and generates the SBAS messages. The format of these messages is compliant with the international standard also used by the Wide Area Augmentation System (WAAS), the European Geostationary Navigation Overlay Service (EGNOS), and the Japanese Multi-functional Transport Satellite (MTSAT) Satellite-based Augmentation System (MSAS).

Format Limitations. The current SBAS format has a limited capability for broadcasting corrections for GLONASS and GPS satellites combined. There is space for only 51 satellites, insufficient for the current number of satellites on orbit. As a result, studies are looking into the efficiency of SDCM data broadcasting in an attempt to resolve this contradiction. The three main options are: use a dynamic satellite mask, use two CDMA signals, or provide an additional SBAS message.

Under the first option, SDCM satellites would only broadcast corrections and integrity data for those GLONASS and other GNSS satellites in view of users in the territory of the Russian Federation. For the second option, SDCM satellites would transmit two CDMA signals with independent sets of correc
tions and integrity data on each signal. The third option assumes that the SDCM data stream would have additional messages with information about satellites not included in the initial list of 51.

The first scenario is possible with the current version of the SBAS format. The other two options require some changes in the format of SBAS messages and international coordination. But the SDCM CPF is ready to operate in all of these modes.

Distribution. The main advantage of SBAS is its universal space channel to users. The SDCM orbit constellation will consist of three geostationary satellites from the multifunctional space relay system Luch (see Figure 7). Luch, which means “ray” or “beam” in Russian, will be used to relay communications between low Earth-orbiting spacecraft and ground facilities in Russia in a similar fashion to that of NASA’s Tracking and Data Relay Satellite System. The satellites will also include transponders for relaying SDCM signals from the CPF to users. The first satellite, Luch-5A, will be launched this year and will occupy an orbital slot at 16° west longitude. Luch-5B will be launched in 2012 to a slot at 95° east longitude. The full constellation will be deployed by 2014 with the launch of Luch-4 into a slot at 167° east longitude.

Figure 7. Multifunctional relay system Luch.

Wideband transponders (22 MHz) will be installed on board the Luch-5A and Luch-5B satellites. These transponders will transmit signals on a carrier frequency of 1575.42 MHz. As the SDCM service area is Russian territory, the main beam will be directed to the north with an angle of 7 degrees relative to the direction to the equator. The transmitted power will be 60 watts and will give a signal power level at the Earth’s surface roughly equal to that of GLONASS and GPS signals, about –158 dBW.

SDCM will also provide service through the Internet. A system website (www.sdcm.ru) already gives users information about real-time and a posteriori GLONASS and GPS monitoring (see Figure 8). An SDCM data-broadcasting ground system has been developed and is being tested now. It will help to verify SDCM data before the Luch satellites are launched. SDCM SBAS messages will be transmitted through the Internet in real time using the SISNeT (Signal in Space through the Internet) approach. The SISNeT protocol was developed for relaying EGNOS messages over the Internet.

Figure 8. SDCM website, www.sdcm.ru.(Click to enlarge.)

A set of experiments was carried out to evaluate SDCM performance. In one experiment, 130 hours of raw pseudorange data was processed to generate the results shown in Figure 9. The upper plot shows the positioning results of a stand-alone receiver working only with the GLONASS and GPS signals. The lower plot presents results of GLONASS/GPS/SDCM navigation. It is clear that the SDCM ephemeris and clock corrections improve user accuracy by more than a factor of two.

Figure 9. SDCM tests results; (a) without and (b) with SDCM corrections.

However, precise point positioning (PPP) technology, based on post-processing dual-frequency carrier-phase measurements with precise satellite ephemeris and clock data, expands the areas of practical use of satellite positioning without complex user ground infrastructure of reference stations and wireless communication channels. Studies have already demonstrated that decimeter-level PPP is possible using GLONASS data or GLONASS data in combination with GPS data. Tests are under way to deliver the precise satellite ephemeris and clock data over the Internet to allow real-time PPP. We can envisage that some time in the future, the ephemeris and clock data could be provided to users in real time using satellite signals.

Future SDCM Satellites. The first SDCM satellites will provide service over the main part of Russia, excluding northern regions. To cover those regions, the SDCM orbit constellation could be enlarged using satellites in circular, inclined geosynchronous orbit (GSO); inclined, elliptical geosynchronous orbit (IGSO); or Molniya-type highly elliptical orbit (HEO) with an orbital period of precisely one-half of a sidereal day.

A comparative availability analysis for satellites with different orbits shows that using four GSO/IGSO/HEO satellites in two planes allows a user anywhere in Russia to continuously receive a signal from two satellites with a minimum elevation angle of 5 degrees. If the elevation mask angle is 30 degrees, availability will fall to 0.9 for IGSO satellites and 0.8 for HEO satellites. An orbit constellation of GSO satellites provides an availability of 0.8 and 0.3 for 5- and 30-degree mask angles respectively.

It is important to point out that the development of satellite orbit and clock prediction technology allows us to consider the possibility of using GSO, IGSO, or HEO satellites for ranging signal broadcasting. In that case, the navigation message could include precise ephemerides and clock data for all GNSS satellites to provide the data for a PPP service as mentioned earlier.

Conclusion

GLONASS development is entering a new historical phase. New CDMA navigation signals and deployment of a national SBAS system will provide not only a new quality of navigation service, but the basis for a regional precise navigation system with an accuracy of a few decimeters for users in Russia and neighboring countries.

Acknowledgment

This article is based on the paper “GLONASS Developing Strategy” presented at ION GNSS 2010, the 23rd International Technical Meeting of The Institute of Navigation, Portland, Oregon, September 21–24, 2010.

Yuri Urlichich is the general director and general designer of the Joint Stock Company (JSC) Russian Space Systems, formerly the Russian Institute of Space Device Engineering, headquartered in Moscow. He is a GLONASS general designer, doctor of science, professor, and author of more than 150 papers and 20 patents.

Valeriy Subbotin is a first deputy general director and general designer of JSC Russian Space Systems and a doctor of science. He has worked in the space industry for more than 40 years and has published more than 45 papers.

Grigory Stupak is a deputy general director and general designer of JSC Russian Space Systems, a GLONASS deputy general designer, and a professor of Bauman Moscow State Technical University (BMSTU). He has worked in the space industry for 35 years and has published more than 150 papers.

Vyacheslav Dvorkin is a deputy general designer of JSC Russian Space Systems and a doctor of science. Dvorkin has been developing GLONASS, GNSS augmentations, and user equipment for more than 35 years. He is an author of 50 papers in the satellite navigation field.

Alexander Povalyaev is a deputy head of division in JSC Russian Space Systems and a professor of Moscow Aviation Institute. He has been developing methods and algorithms for processing GNSS carrier-phase measurements for 30 years and has published more than 40 papers.

Sergey Karutin is a deputy head of division in JSC Russian Space Systems and an assistant professor at BMSTU. Karutin has been on the GLONASS team since 1998, developing GNSS augmentations and user equipment. He received a Ph.D. degree in 2004.

FURTHER READING

• GLONASS Background and Use

“GPS, GLONASS, and More: Multiple Constellation Processing in the International GNSS Service” by T. Springer and R. Dach in GPS World, Vol. 21, No. 6, June 2010, pp. 48–58.

“The Future is Now: GPS + GLONASS + SBAS = GNSS” by L. Wanninger in GPS World, Vol. 19, No. 7, July 2008, pp. 42–48.

“GLONASS: Review and Update” by R.B. Langley in GPS World, Vol. 8, No. 7, July 1997, pp. 46–51.

• GLONASS Current and Future Signal Structures

GLONASS Interface Control Document, Edition 5.1, Russian Institute of Space Device Engineering, Moscow, 2008.

“The Spreading and Overlay Codes for the L1C Signal” by J.J. Rushanan in Navigation, Vol. 54, No. 1, Spring 2007, pp. 43–51.

Spread Spectrum Systems for GNSS and Wireless Communications by J.K. Holmes, Artech House, Inc., Norwood, Massachusetts, 2007.

“The Galileo Code and Others” by G.W. Hein, J.-A. Avila-Rodriguez, and S. Wallner in Inside GNSS, Vol. 1, No. 6, September 2006, pp. 62–74.

• System for Differential Correction and Monitoring

“Russian System for Differential Correction and Monitoring: A Concept, Present Status, and Prospects for Future” by S.V. Averin, V.V. Dvorkin, and S.N. Karutin in Proceedings of ION GNSS 2006, the 19th International Technical Meeting of the Satellite Division of The Institute of Navigation, Fort Worth, Texas, September 26–29, 2006, pp. 3037–3044.

Minimum Operational Performance Standards for Global Positioning/Wide Area Augmentation System Airborne Equipment, RTCA/DO-229D, prepared by SC-159, RTCA Inc., Washington, D.C., December 13, 2006.

“Appendix B. Technical Specifications for the Global Navigation Satellite System (GNSS)” in Aeronautical Telecommunications: International Standards and Recommended Practices, Annex 10 to the Convention on International Civil Aviation, Vol. I. Radio Navigation Aids, (6th ed.), International Civil Aviation Organization, Montreal, Quebec, Canada, 2006.

• SISNeT

“Proposal of an Internet-Based EGNOS Receiver Architecture and Demonstration of the SISNeT Concept” by E. González, M. Toledo, A. Catalina, C. Barredo, F. Torán, and A. Salonico in Proceedings of ION GPS/GNSS 2003, the 16th International Technical Meeting of the Satellite Division of The Institute of Navigation, Portland, Oregon, September 9-12, 2003, pp. 1628–1641.

• Precise Point Positioning

“An Evaluation of OmniStar XP and PPP as a Replacement for DGPS in Airborne Applications” by J.S. Booth, and R.N. Snow in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 22–25, 2009, pp. 1188–1194.

“Precise Point Positioning for Real-Time Determination of Co-Seismic Crustal Motion” by P. Collins, J. Henton, Y. Mireault, P. Héroux, M. Schmidt, H. Dragert, and S. Bisnath in Proceedings of ION GNSS 2009, Savannah, Georgia, September 22–25, 2009, pp. 2479–2488.

“Orbits and Clocks for GLONASS Precise-Point-Positioning” by R. Píriz, D. Calle, A. Mozo, P. Navarro, D. Rodríguez, and G. Tobías in Proceedings of ION GNSS 2009, Savannah, Georgia, September 22–25, 2009, pp. 2415–2424.

“Study on Precise Point Positioning Based on Combined GPS and GLONASS” by X. Li, X. Zhang, and F. Guo in Proceedings of ION GNSS 2009, Savannah, Georgia, September 22–25, 2009, pp. 2449–2459.

April 1, 2011
Innovation: Realistic Randomization

A New Way to Test GNSS Receivers

By Alexander Mitelman

INNOVATION INSIGHTS by Richard Langley

GNSS RECEIVER TESTING SHOULD NEVER BE LEFT TO CHANCE. Or should it? There are two common approaches to testing GNSS receivers: synthetic and realistic. In synthetic testing, a signal simulator is programmed with specific satellite orbits, receiver positions, and signal propagation conditions such as atmospheric effects, signal blockage, and multipath. A disadvantage of such testing is that the models used to generate the synthetic signals are not always consistent with the behavior of receivers processing real GNSS signals. Realistic testing, on the other hand, endeavors to assess receiver performance directly using the signals actually transmitted by satellites. The signals may be recorded digitally and played back to receivers any number of times. While no modeling is used, the testing is specific to the particular observing scenario under which the data was recorded including the satellite geometry, atmospheric conditions, multipath behavior, and so on. To fully examine the performance of a receiver using data collected under a wide variety of scenarios would likely be prohibitive. So, neither testing approach is ideal. Is there a practical alternative? The roulette tables in Monte Carlo suggest an answer.

Both of the commonly used testing procedures lack a certain characteristic that would better assess receiver performance: randomness. What is needed is an approach that would easily provide a random selection of realistic observing conditions. Scientists and engineers often use repeated random samples when studying systems with a large number of inputs especially when those inputs have a high degree of uncertainty or variability. And mathematicians use such methods to obtain solutions when it is impossible or difficult to calculate an exact result as in the integration of some complicated functions. The approach is called the Monte Carlo method after the principality’s famous casino. Although the method had been used earlier, its name was introduced by physicists studying random neutron diffusion in fissile material at the Los Alamos National Laboratory during the Second World War.

In this month’s article, we look at an approach to GNSS receiver testing that uses realistic randomization of signal amplitudes based on histograms of carrier-to-noise-density ratios observed in real-world environments. It can be applied to any simulator scenario, independent of scenario details (position, date, time, motion trajectory, and so on), making it possible to control relevant parameters such as the number of satellites in view and the resulting dilution of precision independent of signal-strength distribution. The method is amenable to standardization and could help the industry to improve the testing methodology for positioning devices — to one that is more meaningfully related to real-world performance and user experience.

Virtually all GNSS receiver testing can be classified into one of two broad categories: synthetic or realistic. The former typically involves simulator-based trials, using a pre-defined collection of satellite orbits, receiver positions, and signal propagation models (ionosphere, multipath, and so on). Examples of this type of testing include the 3rd Generation Partnership Project (3GPP) mobile phone performance specifications for assisted GPS, as well as the “apples-to-apples” methodology described in an earlier GPS World article (see Further Reading).

The primary advantage of synthetic testing is that it is tightly controllable and completely repeatable; where a high degree of statistical confidence is required, the same scenario can be run many times until sufficient data has been collected. Also, this type of testing is inherently self-contained, and thus amenable to testing facilities with modest equipment and resources.

Synthetic approaches have significant limitations, however, particularly when it comes to predicting receiver performance in challenging real-world environments. Experience shows that tests in which signal levels are fixed at predetermined levels are not always predictive of actual receiver behavior. For example, a receiver’s coherent integration time could in principle be tuned to optimize acquisition at those levels, resulting in a device that passes the required tests but whose performance may degrade in other cases. More generally, it is useful to observe that the real world is full of randomness, whereas apart from intentional variations in receiver initialization, the primary source of randomness in most synthetic tests is simply thermal noise.

By comparison, most realistic testing approaches are designed to measure real-world performance directly. Examples include conventional drive testing and so-called “RF playback” systems, both of which have also been described in recent literature (see Further Reading). Here, no modeling or approximation is involved; the receiver or recording instrument is physically operated within the signal environment of interest, and its performance in that environment is observed directly. The accuracy and fidelity of such tests come with a price, however. All measurements of this type are inherently literal: the results of a given test are inseparably linked to the specific multipath profile, satellite geometry, atmospheric conditions, and antenna profile under which the raw data was gathered. In this respect, the direct approach resembles the synthetic methods outlined above — little randomness exists within the test setup to fully explore a given receiver’s performance space.

Designing a practical alternative to the existing GNSS tests, particularly one intended to be easy to standardize, represents a challenging balancing act. If a proposed test is too simple, it can be easily standardized, but it may fall well short of capturing the complexities of real-world signals. On the other hand, a test laden with many special corner cases, or one that requires users to deploy significant additional data storage or non-standard hardware, may yield realistic results for a wide variety of signal conditions, but it may also be impractically difficult to standardize.

With those constraints in mind, this article attempts to bridge the gap between the two approaches described above. It describes a novel method for generating synthetic scenarios in which the distribution of signal levels closely approximates that observed in real-world data sets, but with an element of randomness that can be leveraged to significantly expand testing coverage through Monte Carlo methods. Also, the test setup requires only modest data storage and is easily implemented on existing, widely deployed hardware, making it attractive as a potential candidate for standardization.

The approach consists of several steps. First, signal data is gathered in an environment of interest and used to generate a histogram of carrier-to-noise-density (C/N₀) ratios as reported by a reference receiver, paying particular attention to satellite masking to ensure that the probability of signal blockage is calculated accurately. The histogram is then combined with a randomized timing model to create a synthetic scenario for a conventional GNSS simulator, whose output is fed into the receiver(s) under test (RUTs). The performance of the RUTs in response to live and simulated signals is compared in order to validate the fidelity and usefulness of the histogram-based simulation. This hybrid approach combines the benefits of synthetic testing (repeatability, full control, and compactness) with those of live testing (realistic, non-static distribution of signal levels), while avoiding many of the drawbacks of each.

Histograms

The method explored in this article relies on cumulative histograms of C/N₀ values reported by a receiver in a homogeneous signal environment. This representation is compact and easy to implement with existing simulator-based test setups, and provides information that can be particularly useful in tuning acquisition algorithms.

Motivation and Theoretical Considerations. To motivate the proposed approach, consider an example histogram constructed from real-world data, gathered in an environment (urban canyon) where A-GPS would typically be required. This is shown in FIGURE 1, together with a representative histogram of a standard “coarse-time assistance” test case (as described in the 3GPP Technical Standard 34.171, Section 5.2.1) for comparison. (Note that the x-axis is actually discontinuous toward the left side of each plot: the “B” column designates blocked signals, and thus corresponds to C/N₀ = –∞.)

From the standpoint of signal distributions, it is evident that existing test standards may not always model the real world very accurately.

FIGURE 1a. Example histogram of a real-world urban canyon, the San Francisco financial district;.

Figure 1b. Example histograms of 3GPP TS 34.171 “coarse-time assistance” test case).

The histogram is useful in other ways as well. Since the data set is normalized (the sum of all bin heights is 1.0), it represents a proper probability mass function (PMF) of signal levels for the environment in question. As such, several potentially useful parameters can be extracted directly from the plot: the probability of a given signal being blocked (simply the height of the leftmost bin); upper and lower limits of observed signal levels (the heights of the leftmost and rightmost non-zero bins, respectively, excluding the “blocked” bin); and the center of mass, defined here as

(1)

where y[n] is the height of the nth bin (dimensionless), x[n] is the corresponding C/N₀ value (in dB-Hz), and x[“B”] = –∞ by definition.

Finally, representing environmental data as a PMF enables one additional theoretical calculation. The design of the 3GPP “coarse-time assistance” test case illustrated above assumes that a receiver will be able to acquire the one relatively strong signal (the so-called “lead space vehicle (SV)” at -142 dBm) using only the assistance provided, and will subsequently use information derivable from the acquired signal (such as the approximate local clock offset) to find the rest of the satellites and compute a fix. Suppose that for a given receiver, the threshold for acquisition of such a lead signal given coarse assistance is P_i (expressed in dB-Hz). Then the probability of finding a lead satellite on a given acquisition attempt can be estimated directly from the histogram:

(2)

where is the average number of satellites in view over the course of the data set. A similar combinatorial calculation can be made for the conditional probability of finding at least three “follower” satellites (that is, those whose signals are above the receiver’s threshold for acquisition when a lead satellite is already available).

The product of these two values represents the approximate probability that a receiver will be able to get a fix in a given signal environment, expressed solely as a function of the receiver’s design parameters and the histogram itself. When combined with empirical data on acquisition yield from a large number of start attempts in an environment of interest, this calculation provides a useful way of checking whether a particular histogram properly captures the essential features of that environment. This validation may prove especially useful during future standardization efforts.

Application to Acquisition Tuning. In addition to the calculations based on the parameters discussed above, histograms also provide useful information for designing acquisition algorithms, as follows.

Conventionally, the acquisition problem for GNSS is framed as a search over a three-dimensional space: SV pseudorandom noise code, Doppler frequency offset, and code phase. But in weak signal environments, a fourth parameter, dwell time – the predetection integration period, plays a significant role in determining acquisition performance. Regardless of how a given receiver’s acquisition algorithm is designed, dwell time (or, equivalently, search depth) and the associated signal detection threshold represent a compromise between acquisition speed and performance (specifically, the probabilities of false lock and missed detection on a given search). To this end, any acquisition routine designed to adjust its default search depth as a function of extant environmental conditions may be optimized by making use of the a priori signal level PMF provided by the corresponding histogram(s).

Data Collection

The hardware used to collect reference data for histogram generation is simple, but care must be taken to ensure that the data is processed correctly. The basic setup is shown in FIGURE 2.

Figure 2. Data collection setup with a reference receiver generating NMEA 0183 sentences or in-phase and quadrature (I/Q) raw data and one or more test receivers performing multiple time-to-first-fix (TTFF) measurements.

It is important to note that the individual components in the data-collection setup are deliberately drawn here as generic receivers, to emphasize that the procedure itself is fundamentally generic. Indeed, as noted below, future efforts toward standardizing this testing methodology will require that it generate sensible results for a wide variety of RUTs, ideally from different manufacturers. Thus, the intention is that multiple receivers should eventually be used for the time-to-first-fix (TTFF) measurements at bottom right in the figure. For simplicity, however, a single test receiver is considered in this article.

Procedure. The experiment begins with a test walk or drive through an environment of interest. Since an open sky environment is unlikely to present a significant challenge to almost any modern receiver, a moderately difficult urban canyon route through the narrow alleyways of Stockholm’s Gamla Stan (Old Town) was chosen for the initial results presented in this article. The route, approximately 5 kilometers long, is shown in FIGURE 3 (top). For the TTFF trials gathered along this route, assisted starts with coarse-time aiding (±2 seconds) were used to generate a large number of start attempts during the walk, ensuring reasonable statistical significance in the results (115 attempts in approximately 60 minutes, including randomized idle intervals between successive starts).

Once the data collection is complete, the reference data set is processed with a current almanac and an assumed elevation angle mask (typically 5 degrees) to produce an individual histogram for each satellite in view, along with a cumulative histogram for the entire set, as shown in Figure 3 (bottom). The masking calculation is particularly important in properly classifying which non-reported C/N₀ values should be ignored because the satellite in question is below the elevation angle mask at that location and time, and which should be counted as blocked signals.

Figure 3a. Data collection, Gamla Stan (Old Town), Stockholm (route and street view).

Figure 4. Fluctuation timing models (top: “Multi SV” variant; bottom: “Indiv SV” variant).

In addition to proper accounting for satellite masking, the raw source data should also be manually trimmed to ensure that all data points used to build the histogram are taken homogeneously from the environment in question. Thus the file used to generate the histogram in Figure 3 was truncated to exclude the section of “open sky” conditions between the start of the file and the southeast corner of the test area, and similarly between the exit from the test area and the end of the file.

Finally, the resulting histogram is combined with a randomized timing model to create a simulator scenario, which is used to re-test the same RUTs shown in Figure 2.

Reference Receiver Considerations. The accuracy of the data collection described above is fundamentally limited by the performance of the reference receiver in several ways.

First, the default output format for GNSS data in many receivers is that of the National Marine Electronics Association (NMEA) 0183 standard (the histograms presented in this article were derived from NMEA data). This is imperfect in that the NMEA standard non-proprietary GSV sentence requires C/N₀ values to be quantized to the nearest whole dB-Hz, which introduces small rounding errors to the bin heights in the histograms. (In this study, this effect was addressed by applying a uniformly distributed ±0.5 dB-Hz dither to all values in the corresponding simulated scenario, as discussed below.) If finer-grained histogram plots are required, an alternative data format must be used instead.

Second, many receivers produce data outputs at 1 Hz, limiting the ability to model temporal variations in C/N₀ to frequencies less than 0.5 Hz, owing to simple Nyquist considerations. While the raw data for this study was obtained at walking speeds (1 to 2 meters per second), and thus unlikely to significantly misrepresent rapid C/N₀ fading, studies done at higher speeds (such as test drives) may require a reference receiver capable of producing C/N₀ measurements at a higher rate.

A third limitation is the sensitivity of the reference receiver. Ideally, the reference device would be able to track all signals present during data gathering regardless of signal strength, and would instantaneously reacquire any blocked signals as soon as they became visible again. Such a receiver would fully explore the space of all available signals present in the test environment. Unfortunately, no receiver is infinitely sensitive, so a conventional commercial-grade high sensitivity receiver was used in this context. Thus the resulting histogram is, at best, a reasonable but imperfect approximation of the true signal environment.

Finally, a potentially significant error source may be introduced if the net effects of the reference receiver’s noise figure plus implementation loss (NF+IL) are not properly accounted for in preparing the histograms. (If an active antenna is used, the NF of the antenna’s low-noise amplifier essentially determines the first term.) The effect of incorrectly modeling these losses is that the entire histogram, with the exception of the “blocked” column, is shifted sideways by a constant offset.

The correction applied to the histogram to account for this effect must be verified prior to further acquisition testing. This can be done by generating a simulator scenario from the histogram of interest, as described below, and recording a sufficiently long continuous data set using this scenario and the reference receiver. A corresponding histogram is then built from the reference receiver’s output, as before, and compared to the histogram of the original source data. The amplitude of the “blocked” column and the center of mass are two simple metrics to check; a more general way of comparing histograms is the two-sided Kolmogorov-Smirnov test (see “Results”).

Timing Models

The histograms described in the preceding section specify the amplitude distribution of satellite signals in a given environment, but they contain no information about the temporal characteristics of those signals. This section briefly describes the timing models used in the current study, as well as alternatives that may merit further investigation.

In real-world conditions, the temporal characteristics of a given satellite signal depend on many factors, including the physical features of the test environment, multipath fading, and the velocity of the user during data collection. Various timing models can be used to simulate those temporal characteristics in laboratory scenarios.

Perhaps the simplest model is one in which signal levels are changed at fixed intervals. This is trivial to implement on the simulator side, but it is clearly unlikely to resemble the real-world conditions mentioned above. A second alternative would be to generate timing intervals based on the Allan (or two-sample) variance of individual C/N₀ readings observed during data collection as a measure of the stability of the readings. While this is more physically realistic than an arbitrarily chosen interval as described above, it is still a fixed interval. These observations suggest that a timing model including some measure of randomness may represent a more realistic approach.

One statistical function commonly used for real-world modeling of discrete events (radioactive decay, customers arriving at a restaurant, and so on) is the Poisson arrival process. This process is completely described with a single non-negative parameter, λ, which characterizes the rate at which random events occur. Equivalently, the time between successive events in such a process is itself a random variable described by the exponential probability distribution function:

(3 )

The resulting inter-event timings described by this function are strictly non-negative, which is at least physically reasonable, and directly controllable by varying the timing parameter λ. For simplicity, then, the Poisson/exponential timing model was chosen as an initial attempt at temporal modeling, and used to generate the results presented in this article.

Two variants of the Poisson/exponential timing model are considered. In the first, defined herein as the “Multi SV” case, a single thread determines the timing of fluctuation events, and the power levels of one or more satellites are adjusted at each event. In the second variant, defined as the “Indiv SV” case, each simulator channel receives its own individual timing thread, and all fluctuation events are interleaved in constructing the timing file for the simulator. These two variants are shown schematically in FIGURE 4.

Figure 4. Fluctuation timing models (top: “Multi SV” variant; bottom: “Indiv SV” variant).

Constructing Scenarios

Once a target histogram is available, it is necessary to generate random signal amplitudes for use with a simulator scenario. This is done by means of a technique known as the probability integral transform (PIT). This approach uses the c
umulative distribution function (or, in the discrete case considered here, a modified formulation based on the cumulative mass function) of a probability distribution to transform a sequence of uniformly distributed random numbers into a sequence whose distribution matches the target function.

Finally, the random signal levels generated by the PIT process are assigned to individual simulator channels according to a set of timed events as described in the preceding section, completing the randomized scenario to be used for testing.

Results

Given a simulator scenario constructed as described above, the RUTs originally included in the data collection campaign are again used to conduct acquisition tests, this time driven from the simulator.

To validate that a particular fluctuating scenario properly represents the live data, it is necessary to quantify two things: how well a generated histogram matches the source data, and how well a receiver’s acquisition performance under simulated signals matches its behavior in the field. At first these may appear to be two qualitatively different problems, but a mathematical tool known as the two-sided Kolmogorov-Smirnov (K-S) test can be used for both tasks.

Validation of Experimental Setup. As a first step toward validating that the C/N₀ profile of the simulated signals matches that of the reference data, TABLE 1 gives the values of the two-sided K-S test statistic, D (a measure of the greatest discrepancy between a sample and the reference distribution), for histograms generated with the reference receiver for the two timing-thread models described above and several values of the Poisson/exponential parameter, λ. The reference cumulative mass function (CMF) for each test was derived from the histogram generated for the raw (empirically collected) data set.

These results illustrate good agreement (D < 0.05) between the overall signal distribution profile in the empirical data set and that in each of the six simulated fluctuating scenarios.

As a further check, TABLE 2 shows the same K-S statistic for the histogram generated from the “Multi SV” timing model as a function of several NF+IL values. As before, the reference CMF comes from the raw (empirically collected) data set, and the same reference receiver was used to generate data from the simulator scenario. Evidently, an NF+IL value of 4 dB gives good agreement between empirical and simulated data sets.

Validation of Receiver Performance. Finally, TTFF tests with the simulated scenarios described above are conducted with the same receiver(s) used in the original data gathering session. Here, the K-S test is used to compare the live and simulated TTFF results rather than signal distributions. An example result, illustrating cumulative distribution functions of TTFF, is shown in FIGURE 5 for the live data set collected during the original data gathering session, alongside three results from the “Multi SV” fluctuating model, generated with NF+IL = 4 dB and several different values of the Poisson/exponential timing parameter, λ. While agreement with live data is not exact for any of the simulated scenarios, the λ^-1 = 3.0 seconds case appears to correspond reasonably well (D < 0.10).

FIGURE 5 Time-to-first-fix cumulative distribution functions from live and simulated data (“Multi SV” variant with NF+IL = 4 dB).

Conclusions and Future Work

This article has introduced a novel approach to testing GNSS receivers based on histograms of C/N₀ values observed in real-world environments.

Much additional work remains. For the proposed method to be amenable to standardization, it is obviously necessary to gather data from many additional environments. Indeed, it appears likely that no one histogram will encapsulate all environments of a particular type (such as urban canyons), so significant additional experimentation and data collection will be required here. Also, as mentioned at the beginning of the article, the proposed method will need to be tested with multiple receivers to verify that a particular result is not unique to any specific brand or architecture. Finally, higher rate C/N₀ source data may also be necessary to capture the rapid fades that may be encountered in dynamic scenarios, such as drive tests, and the fluctuation timing models will need to be revisited once such data becomes available.

Acknowledgments

The author gratefully acknowledges the assistance of Jakob Almqvist, David Karlsson, James Tidd, and Christer Weinigel in conducting the experiments described in this article. Thanks also to Ronald Walken for valuable insights on the accurate treatment of the source environment in calculating target histograms. This article is based on the paper “Fluctuation: A Novel Approach to GNSS Receiver Testing” presented at ION GNSS 2010.

Alexander Mitelman is the GNSS research manager at Cambridge Silicon Radio, headquartered in Cambridge, U.K. He earned his S.B. degree from the Massachusetts Institute of Technology and M.S. and Ph.D. degrees from Stanford University, all in electrical engineering. His research interests include signal-quality monitoring and the development of algorithms and testing methodologies for GNSS.

FURTHER READING

• GNSS Receiver Testing in General
GPS Receiver Testing, Application Note by Agilent Technologies. Available online at http://cp.literature.agilent.com/litweb/pdf/5990-4943EN.pdf.

• Synthetic GNSS Receiver Testing
“Apples to Apples: Standardized Testing for High-Sensitivity Receivers” by A. Mitelman, P.-L. Normark, M. Reidevall, and S. Strickland in GPS World, Vol. 19, No. 1, January 2008, pp. 16–33.

Universal Mobile Telecommunications System (UMTS); Terminal conformance specification; Assisted Global Positioning System (A-GPS); Frequency Division Duplex (FDD), 3GPP Technical Specification 34.171, Release 7, Version 7.0.1, July 2007, published by the European Telecommunications Standards Institute, Sophia Antipolis, France. Available online at http://www.3gpp.org/.

• Realistic GNSS Receiver Testing
“Record, Replay, Rewind: Testing GNSS Receivers with Record and Playback Techniques” by D.A. Hall in GPS World, Vol. 21, No. 10, October 2010, pp. 28–34.

“Proper GPS/GNSS Receiver Testing” by E. Vinande, B. Weinstein, and D. Akos in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 22–25, 2009, pp. 2251–2258.

“Advanced GPS Hybrid Simulator Architecture” by A. Brown and N. Gerein in Proceedings of The Institute of Navigation 57th Annual Meeting/CIGTF 20th Guidance Test Symposium, Albuquerque, New Mexico, June 11–13, 2001, pp. 564–571.

• Receiver Noise
“Measuring GNSS Signal Strength: What is the Difference Between SNR and C/N₀?” by A. Joseph in Inside GNSS, Vol. 5, No. 8, November/December 2010, pp. 20–25.

“GPS Receiver System Noise” by R.B. Langley in GPS World, Vol. 8, No. 6, June 1997, pp. 40–45.

Global Positioning System: Theory and Applications, Vol. I, edited by B.W. Parkinson and J.J. Spliker Jr., published by the American Institute of Aeronautics and Astronautics, Inc., Washington, D.C., 1996.

• Test Statistics
“The Probability Integral Transform and Related Results” by J. Agnus in SIAM Review (a publication of the Society for Industrial and Applied Mathematics), Vol. 36, No. 4, December 1994, pp. 652–654, doi:10.1137/1036146

“Kolmogorov-Smirnov Test” by T.W. Kirkman on the College of Saint Benedict and Saint John’s University Statistics to Use website: http://www.physics.csbsju.edu/stats/KS-test.html.

• NMEA 0183
NMEA 0183, The Standard for Interfacing Marine Electronic Devices, Ver. 4.00, published by the National Marine Electronics Association, Severna Park, Maryland, November 2008.

“NMEA 0183: A GPS Receiver Interface Standard” by R.B. Langley in GPS World, Vol. 6, No. 7, July 1995, pp. 54–57.

Unofficial online NMEA 0183 descriptions: NMEA data; NMEA Revealed by E.S. Raymond, Ver. 2.3, March 2010.

March 1, 2011
Innovation: GNSS and the Ionosphere

What’s in Store for the Next Solar Maximum?

By Anna B.O. Jensen and Cathryn Mitchell

Although the sun can become disturbed at any time, solar activity is correlated with the approximately 11-year cycle of spots on the sun’s surface. We are just coming out of a minimum in the solar cycle and headed for the next maximum, predicted to occur around the middle of 2013. How significantly will GNSS users be affected? In this month’s column, two ionosphere experts tell us what might be in store.

INNOVATION INSIGHTS by Richard Langley

“HERE COMES THE SUN / here comes the sun / And I say / it’s all right.”

Is it? Of course, George Harrison was referring to the welcome return of the sun after a long dreary English winter. But can GNSS users sing the same refrain?

The signals from global navigation satellites must transit the ionosphere on their way to receivers on or near the Earth’s surface. The passage exacts a toll in the form of an added delay of the pseudorandom-noise-code signals and an advance of the phase of the signals’ carriers, due to the presence of the ionosphere’s free electrons. These perturbations must be ameliorated in some way to achieve high accuracy in GNSS positioning, navigation, and timing applications.

Where do the ionosphere’s electrons come from? For the most part, they are valence electrons, stripped from upper atmosphere atoms and molecules by the extreme ultraviolet light continuously emitted by the sun. On the Earth’s night-side, the electrons and the ionized atoms and molecules tend to recombine. This ionization and recombination process, along with the interactions of the particles with the Earth’s magnetic field, governs the density of the electrons at a particular location and time. The ionosphere is also affected by the solar wind, and its associated magnetic field, but the cocoon established by the Earth’s magnetic field (the magnetosphere) tends to deflect the solar wind so that it usually has little influence on the ionosphere.

Normally, the sun is quiescent: its electromagnetic and particle radiation is fairly constant, and its effects on the ionosphere benign. The delay in GNSS code observations and the advance in phase observations can be readily estimated and removed from the observations using a variety of models and methods. However, the sun can become disturbed, giving rise to occasional violent outbursts with large increases in electromagnetic and particle radiation. These outbursts can radically change the distribution of the electrons in the ionosphere, reducing the effectives of some amelioration methods. The electron density variability can become so rapid that a GNSS receiver can lose lock on satellite signals. And an increase in the sun’s radio emissions can become so large as to drown out GNSS signals on the sunlight side of the Earth.

Although the sun can become disturbed at any time, solar activity is correlated with the approximately 11-year cycle of spots on the sun’s surface. We are just coming out of a minimum in the solar cycle and headed for the next maximum, predicted to occur around the middle of 2013. How significantly will GNSS users be affected? In this month’s column, two ionosphere experts tell us what might be in store.

GNSS satellite signals are affected by the space environment and the Earth’s atmosphere as they travel from satellites at an altitude of about 20,000 kilometers above the surface of the Earth to receivers located at, or close to, the surface.

In the upper part of the Earth’s atmosphere, the ionosphere, which is located from about 80 to 1,000 kilometers above the surface of the Earth, satellite signals are affected by the free electrons stripped from atoms and molecules by ionization. The signals are refracted by this plasma, which changes their speed of travel. The effect is mainly a function of the number of free electrons present, the electron density.

In the lower parts of Earth’s atmosphere, in the troposphere and the stratosphere — where the atoms and molecules are electrically neutral — the satellite signals experience additional refraction. Here the effect is a function of pressure, temperature, and humidity. The effect of the troposphere and stratosphere is often just referred to as the “tropospheric effect” in GNSS positioning as it is in the troposphere where most of the neutral atmosphere refraction occurs.

The ionospheric and tropospheric effects on satellite signals must be accounted for in the GNSS positioning process in order to obtain reliable and accurate position solutions. In this article, we look at the ionospheric effect on satellite signals. Although the variation in signal speed is the largest direct ionospheric effect on the GNSS satellite signals, scintillation is another important effect. Scintillation occurs when irregularities in the electron density of the ionosphere cause rapid changes in the phase and amplitude of the transmitted signals. These changes might cause a GNSS receiver to lose lock on a satellite signal. This means in practice that satellite signals are lost, or signal tracking can be rather difficult, during scintillation events. However, we restrict our article to the subject of the propagation speed of the signals and do not consider scintillation further.

In the following, we review characteristics of the ionospheric effect on GNSS satellite signals as well as the predictions of increased ionospheric activity for the coming years and the consequences for GNSS users.

Signals

The ionosphere as a whole is electrically neutral, but it contains a significant number of free electrons and ions. The negatively charged free electrons affect the electromagnetic satellite signals in various ways. Most important is the signal delay affecting code (pseudorange) measurements, also called the “ionospheric delay” (and the associated advance of carrier-phase measurements), which is caused by a change in the refractive index along the signal path. The refractive index changes continuously as a function of the composition of the transmission media all the way from the satellites to the GNSS receivers.

For the majority of the signal path — that is, from the satellite at an altitude of about 20,000 kilometers down to approximately 1,000 kilometers above the surface of the Earth — the change in the refractive index is usually sufficiently small to ignore when the GNSS satellite signals are used for positioning at the surface of the Earth (although, at times, the region above the ionosphere — the plasmasphere — can affect GNSS signals). We therefore use the approximation that the first part of the signal path is in a vacuum where the propagation of GNSS satellite signals is not affected.

Then, when the signals enter the ionosphere, we must consider the signal delay, and even though the density of electrons is largest at an altitude around 300 kilometers, we must consider the total number of electrons experienced by a satellite signal all the way through the ionosphere.

The size of the so-called first order effect of the signal delay, d, given in meters, can be modeled by the expression in Equation (1),

(1)

where f is the GNSS signal frequency, for instance 1.57542 x 10⁹ Hz for the GPS L1 frequency. The constant 40.3 is derived from the values of the electron charge, the electron mass, and the permittivity of free space. Finally, TEC is an abbreviation for total electron content and this value is given by integrating the number of free electrons along the signal path in a cross section of one square meter.

It turns out that the “delay” affecting carrier-phase measurements has exactly the same magnitude as the signal delay but is negative. In other words, the phase is advanced.

In practice, for single-frequency receivers, it is not possible to obtain the actual number of electrons along the signal path for every satellite signal, and we therefore need other models to predict or estimate the electron density or the signal delay.

A large number of models and methods for estimating the ionospheric signal delay have been developed. A comparison of some of them is given in a paper by Allain and Mitchell (see Further Reading). The most widely used model is probably the Klobuchar model, named after John Klobuchar, its developer. Coefficients for the Klobuchar model are determined by the GPS control segment and distributed with the GPS navigation message to GPS receivers where the coefficients are inserted into the model equation and used by receivers for estimation of the signal delay caused by the ionosphere.

Dispersion. The ionosphere is dispersive for radio waves, which means that the GNSS ionospheric signal delay is a function of the frequency of the signal. If pseudorange measurements from more than one frequency are available, for instance from dual-frequency GPS receivers, this can be used for enhanced modeling of the ionospheric effect by using combinations of the measurements made on both frequencies.

The basic expression for estimation of the ionospheric delay for dual-frequency code-based positioning is shown in Equation (2),

(2)

where d is the ionosphere delay, P denotes pseudorange, and f denotes frequency. The subscript notation L1 and L2 refers to the GPS L1 and L2 frequencies, respectively.

For high-accuracy carrier-phase-based positioning, an ionosphere-free combination of carrier-phase observations of the L1 and L2 frequencies is often used to reduce the effect of the ionospheric phase advance in the positioning process.

Estimating the ionosphere delay with Equation (2) for code observations or utilizing the ionosphere-free combination of the phase observations compensates for the first order ionospheric effect. This is the major part of the effect, but higher order effects are present, and the size of the residual higher order effects is increased (up to some centimeters) when the ionospheric activity is increasing.

For high-accuracy applications, the difference in the time of transmission and reception of the satellite signals of the various frequencies also must be considered as the signals on various frequencies are not transmitted from the satellites (nor received at a GNSS receiver) at exactly the same time epochs. These differences are normally referred to as the satellite and receiver differential code biases.

It is important also to note in this context that the noise level on the pseudorange corrected for the ionosphere and on the ionosphere-free carrier-phase observation is increased compared to using the pure single-frequency observations for positioning, but nevertheless these first-order approaches are used successfully in most software and receiver firmware for dual-frequency positioning.

Further developments of ionosphere-free combinations will evolve in the future as the new GPS L5 frequency and the new Galileo and GLONASS frequencies become fully available for multi-frequency ionosphere-free combinations. These more advanced combinations have the potential to further reduce the residual effect of the ionospheric delay in the positioning process.

Summing up, the GNSS signal delay caused by the ionosphere is a function of the electron density of the ionosphere. But what is driving the variation in electron density, and how do we know if it is changing?

Solar Activity and Sunspots

Equation (1) shows that the ionospheric signal delay is a direct function of the total electron content. The number of free electrons in the ionosphere is not constant; it varies significantly with time and space. The number of free electrons is driven by the ionization and recombination processes of the ionosphere, and these processes are in turn driven mainly by extreme ultraviolet radiation from the sun. Radiation from other cosmic sources also has an influence but it is minor compared to the effect of the solar radiation. There are also significant short-term (minutes to hours) changes caused by wave activity from the neutral atmosphere. The ionosphere itself is embedded in the neutral atmosphere — at these altitudes this is known as the thermosphere. The thermosphere is in constant movement due to waves and tides that are generated in situ or ascending from the underlying atmosphere. This thermosphere activity affects the ionosphere and causes some of the short-term variability in the electron density. However, the term “ionospheric activity” generally refers to the variability in electron density as driven by solar activity.

The fact that ionospheric activity is mainly driven by solar activity implies that the temporal variation of the electron content of the ionosphere follows a daily cycle, with the largest TEC values in the early afternoon local time, when the effect of the solar radiation has reached a maximum. Consequently, we see the lowest activity late at night just before sunrise.

There is also a geographic variability in the electron content with the highest electron density in the equatorial region and the lowest density in the high latitude regions. The latter, however, is affected by a larger variability, correlated with auroral activity.

The geographic variation of TEC is illustrated with a global ionosphere map from the Center for Orbit Determination in Europe (CODE) shown in Figure 1. Global ionosphere maps are generated at CODE on a daily basis, and the maps are available on the CODE website (see Further Reading).

Figure 1. Global ionosphere map for November 22, 2010, at 14:00 UTC. (Map generated by CODE, University of Bern.)

The TEC is provided in TEC units (TECU), where one TECU equals 1016 electrons per meter squared.

The sun also emits a constant flow of charged particles called the solar wind. The particles, mostly electrons and protons with energies between about 10 and 100 kilo-electron-volts, travel at an average speed of about 450 kilometers per second, but varying from 200 to 900 kilometers per second depending on solar activity. Although the Earth’s magnetosphere deflects most of the solar wind, the interplanetary magnetic field, which is associated with the solar wind, can cause disturbances in the geomagnetic field. When this happens, particles of the solar wind enter the geomagnetic field and cause increased ionization in the ionosphere. The solar wind therefore also has a large influence on the variability of ionospheric activity. Also, sudden eruptions of the sun such as solar flares and coronal mass ejections (CMEs) cause increased ionization and thereby a larger ionospheric variability.

Figure 2 shows a CME blast and subsequent impact at the Earth.

Figure 2. Coronal mass ejection (CME) and subsequent impact at the Earth. The left part of the illustration is composed of an image from NASA’s Solar Dynamics Observatory spacecraft superimposed on an image from the Solar and Heliospheric Observatory spacecraft jointly operated by NASA and the European Space Agency. The CME cloud arrives at the Earth about two to four days later and is shown being mostly deflected around the Earth’s magnetosphere. The blue paths emanating from the Earth’s poles represent some of its magnetic field lines. (Image: NASA/Goddard Space Flight Center.)

Solar activity and the quantity of emissions from the sun are highly correlated with the number of sunspots on its surface. A sunspot looks like a dark spot because the temperature in a sunspot is lower than that in its surroundings. The generation of sunspots is not well understood, but it is related to anomalies in the solar magnetic field. What is well known, however, is the history of the number of sunspots, because these have been observed since the early 1600s.

The number of sunspots generally follows a cycle of about 11 years. During the last few years (2007–2009), we have experienced a time period with a low number of sunspots. In fact, there were many days in a row without any sunspots visible (see Figure 3). During the next three to four years, the number of sunspots is expected to increase, and this will be followed by a decrease until we reach a new period of low solar activity in 2019–2020.

Figure 3. Images of the sun taken by the Solar and Heliospheric Observatory spacecraft. On the left is an image taken on March 27, 2001, at the peak of the last sunspot cycle. The daily sunspot count was 241. On the right is an image taken on December 15, 2008, near the minimum of the last sunspot cycle, showing no sunspots. (Image: Solar and Heliospheric Observatory)

Numerous investigations of time series of sunspot numbers have been carried out, and even though the cycles generally last 11 years, cycles of 9 and 13 years’ duration have been observed. Also, the cycles vary with respect to the maximum number of sunspots observed during a cycle, and various “cycles of cycles” appear to be present with respect to the strength of the sunspot cycles. For instance, a cycle with a period of about 420 years has been identified in the historic listings of sunspot numbers combined with other observations contributing to the knowledge of solar activity. A very low number of sunspots was observed for a number of years between 1645 and 1715 when the sun was especially calm. This period is often referred to as the Maunder Minimum after the solar astronomer Edward W. Maunder. If the theory of the 420-year cycle is correct, then we will see a period with lower solar activity and fewer sunspot numbers by the end of this century.

But let’s turn our attention to the previous and current sunspot cycles referred to as cycles number 23 and 24 (The 1755–1766 cycle is traditionally numbered “1.”). A new cycle begins with the first observed high-latitude, reversed-polarity sunspot. Reversed polarity means a sunspot with opposite magnetic polarity compared to sunspots from the previous solar cycle. Sunspots from the new and previous cycles initially coexist. Eventually, only the new-cycle sunspots are present. Cycle 24 began on January 4, 2008, when the first reversed-polarity sunspot appeared.

Analyses of observations of solar activity show that the density of the solar wind increases with increasing sunspot number. Also, with a large sunspot number, solar flares and CMEs happen more frequently. Ionospheric storm activity is more common when the sunspot number is high, and this activity increases the variability in ionospheric delays. This all adds up to an increased number of free electrons in the ionosphere and a larger variability, which provides a larger and more variable signal delay for all types of GNSS-based positioning, navigation, and timing during periods with high sunspot numbers.

We know that the sunspot number is expected to increase during the next three to four years. What can be expected and what can we do to minimize the effects of the increased ionospheric activity on positioning, navigation, and timing applications?

The Last Solar High

As mentioned earlier, the current solar cycle is referred to as cycle 24. During the last solar cycle, cycle 23, the GNSS community was alert and aware of what could happen, and therefore many events were observed and analyzed. Among the most well-known events is a sequence of storms during October and November 2003, commonly referred to as the Halloween Storms. The most extreme was the storm on October 30, 2003, which resulted from a CME on October 29 at 20:49 UTC, which subsequently impacted Earth’s magnetic field at 16:20 UTC on October 30 and produced a great geomagnetic storm, which lasted for many hours.

Effects on GPS positioning of this storm have been documented by the GNSS research group of the Royal Observatory of Belgium, where kinematic analyses of data from 36 GNSS stations in Europe showed position errors of more than 10 centimeters in the horizontal and up to 26 centimeters in the vertical between 21:00 and 22:00 UTC on October 30. The position errors were largest for locations in northern Europe including Sweden and Norway. The data analysis was carried out using high-quality carrier-phase data, and the processing was based on using an ionosphere-free linear combination of observations from the L1 and L2 frequencies, whereby the first-order effect of the ionosphere is removed from the results. The position errors are thus caused by mainly higher order ionospheric effects.

For navigation-grade GPS positioning, a U.S. National Atmospheric and Oceanic Administration technical memorandum (see Further Reading) reported that the Wide Area Augmentation System (WAAS) vertical error limit of 50 meters was exceeded for a period of about 11 hours on October 30, 2003. This means that, in practice, WAAS was not available for precision aircraft approaches during that time. The European Geostationary Navigation Overlay Service (EGNOS) was not transmitting during the storm, but simulations carried out later by ESA showed that the boundary regions of the EGNOS coverage area would have been especially affected by a reduction in service availability of about 20–60 percent during that day. The simulations also showed, however, that in the center of the EGNOS coverage area (in the vicinity of northern Italy), the effect would have been much smaller with a reduction in service availability of only 5–6 percent over the day.

Such large storms are also often accompanied by displays of aurora (aurora borealis and aurora australis) at lower latitudes than normal. Figure 4 shows full-sky aurora observed near Fredericton, New Brunswick, Canada (46 degrees north latitude) on October 31, 2003

Figure 4. Photo of red and green auroras observed near Fredericton, New Brunswick, Canada (46 degrees north latitude) early on October 31, 2003. (Courtesy of Richard and Marg Langley.)

During a storm event on November 20, 2003, auroral activity was visible at mid-latitudes over most of North America as far south as Florida and in southern Europe including Italy and Greece.

Eruptions of the sun, often occurring in connection with high sunspot numbers, can have other effects besides the influence on GNSS-based positioning, navigation, and timing. Power-grid blackouts are known to have happened because of geomagnetic storms in connection with the sunspot peaks of both cycles 22 and 23 in 1989 and in 2003, respectively. For instance, the southern part of Sweden experienced a power blackout for several hours during the evening of October 30, 2003.

Also, orbiting satellites can experience problems with the increased radiation and solar wind density. Solar panels are, for instance,
susceptible to increased aging. And many types of satellite communication can be affected by increased ionospheric activity, not only GNSS satellite signals. Signals used for satellite phones, satellite TV, and so on can be affected.

Another phenomenon that can affect GNSS positioning is solar radio storms (also referred to as solar radio bursts) caused by events on the sun, often a solar flare, which creates radio waves that are emitted from the solar atmosphere and can propagate to the Earth where they cause an increased noise level in radio signals. Solar radio storms can cover a wide range of frequencies, including the frequencies used for GNSS. One such storm occurring on December 6, 2006, did affect GNSS positioning. With an increased noise level on the satellite signals, GNSS performance is reduced. If the noise level becomes too large, as a consequence of, for instance, a solar radio storm, GNSS receivers will lose lock on the GNSS signals, whereby positioning performance is further reduced or positioning might even be impossible. Solar radio storms are expected to happen more frequently during the peak of a solar cycle, but the event in December 2006 happened during a period with low solar activity, highlighting the fact that GNSS performance can be affected at any time, even when the sunspot number is low.

Predictions for the Next Solar High

Many predictions for the present solar cycle have been made. Because of the very long period with low solar activity during 2007–2009, some predictions expected a sudden outburst of activity and a very large cycle maximum, while other predictions foretold another increase in solar activity might not occur for many years.

However, with a general increase in the number of sunspots during 2010, it looks like we are now well into solar cycle number 24. Things can still change, but the current predictions say the maximum of the current solar cycle will be lower than the maximum of the last cycle encountered in 2001.

Predictions of sunspot numbers are based on history, logged information on sunspot numbers, and on observations of related geomagnetic activity.

The latest prediction for the current cycle as generated by NASA is shown in Figure 5.

Figure 5. Sunspot cycle 23 and predictions for cycle 24 from NASA’s Marshall Space Flight Center. (Image: NASA)

The curves in Figure 5 show the observed smoothed sunspot number, with smoothing over a period of a year or so, and the predicted value for the remainder of cycle number 24. The dotted lines indicate the observed or expected range of the monthly-averaged sunspot numbers. The plot is updated every month as new data is obtained.

The current prediction for cycle 24 gives a smoothed sunspot number maximum of about 59 in June/July of 2013. This peak is much lower than that of the previous cycle. We are currently two years into cycle 24 and the predicted size continues to fall. According to forecasters, predicting the behavior of a sunspot cycle is fairly reliable once the cycle is well under way (about three years after the minimum in sunspot number occurs). Prior to that time, the predictions are less reliable but nonetheless equally as important.

Even though the maximum of the current solar cycle is expected to be lower than the last peak, it is important for GNSS users to be aware of the effects to be expected during the coming years.

Consequences for GNSS Users

As discussed earlier in this article, GNSS users experience a general satellite signal delay caused by the ionosphere. This signal delay is always present but varies in size. The delay is generally well modeled by most receivers and software to an extent that makes GNSS useable for all of the purposes we know today.

During enhanced ionospheric activity, GNSS users can experience residual ionospheric effects, which can cause reduced positioning, navigation, and timing performance. In such cases, dual-frequency receivers might improve the situation because of the enhanced possibilities for handling the ionospheric effect with dual-frequency data.

During enhanced ionospheric or geomagnetic storm activity caused by sudden eruptions of the sun, increased ionospheric variability will occur. Apart from causing an increased ionospheric signal delay, and thereby increased residual effects in the positioning process, this will also cause increased scintillation effects. These might cause GNSS receivers to lose lock on some or all GNSS satellite signals, reducing performance of the GNSS receiver. In the few very worst cases, GNSS-based positioning, navigation, and timing might not be possible at all for a short interval of time during very high ionospheric activity.

These worst-case scenarios are more prone to happen close to the peak of a solar cycle, which we will meet next during 2013–2014.

However, it is worth noting that for the next peak of the solar cycle, we are much better prepared for the consequences than during the last cycle. GNSS software and receiver technology has been improved to better resist the challenges of increased ionospheric activity during this solar cycle. The improvements are based on experiences gained during the last solar cycles and are to the benefit of many GNSS users. For example, users of wide area augmentation systems such as WAAS and EGNOS have correction and integrity information available, which can be a great help in identifying time epochs when positioning and navigation solutions might not be trustable because of increased ionospheric activity. The integrity information is transmitted from geostationary satellites, and during time periods with extremely high ionospheric activity, the signals with integrity information might be disrupted. This should, however, be detected by the GNSS receiver, so warning messages will be displayed for navigators.

High-accuracy real-time kinematic (RTK) positioning is today often carried out with RTK correction data from a service provider generated using a network of reference stations. Here, indications of increased ionospheric activity can be detected by the software operated by the service provider, and warnings can be distributed to the RTK users.

Warning systems have been improved, and a number of sites on the Internet provide information on current and predicted ionospheric activity (see Further Reading).

Also, in the future, GNSS users will be able to benefit from the increased number of GNSS frequencies available. These frequencies open up opportunities for new and improved methods for correction of the ionospheric delay to the benefit of users who will experience more stable and reliable GNSS performance.

Summary and Conclusion

In this article we have reviewed the ionospheric effects on GNSS satellite signals, how these can be modeled and mitigated, and how they are related to solar activity and the number of sunspots. We have also described how sudden eruptions of the sun can cause increased ionospheric activity and how these events are often correlated with a high sunspot number. Some examples of consequences for GNSS users during the last solar high have been provided, and we have evaluated the predictions for the next solar high and possible consequences for GNSS users.

We are heading towards a period of increased solar activity. GNSS users must expect more disturbances compared to what we have seen for the last four to five years. The peak of the current solar cycle is expected to be lower than the last peak, and therefore consequences for GNSS users should also be less significant. Most of the time GNSS will work very well. But we will likely see a few days with major effects, and since the number of GNSS users is increasing, the overall consequences might also be more severe, not because the ionospheric activity is worse, but simply because more people will be affected.

ANNA B.O. JENSEN is the owner of AJ Geomatics in Copenhagen and a part-time associate professor of the National Space Institute at the Technical University of Denmark (DTU Space). She has a Ph.D. from the University of Copenhagen with co-supervision from the University of Calgary, and has worked in research and development within GNSS and geodesy for more than 15 years. Her current research interests include ionospheric modeling, high accuracy positioning, and navigation in the Arctic.

CATHRYN MITCHELL is a professor in the Department of Electronic and Electrical Engineering at the University of Bath in the United Kingdom and heads the INVERT Centre, which studies inverse problems and tomography over a range of scientific fields, including navigation, space science, and medical imaging. She has a Ph.D. from the University of Wales in Aberystwyth. Mitchell has a particular interest in the use of GNSS measurements to characterize and map the ionosphere.

FURTHER READING

• Introduction to the Ionosphere and Its Effects on GNSS
“The Perfect Solar Storm” by D.N. Baker and J.L. Green in Sky & Telescope, Vol. 121, No. 2, February 2011, pp. 28–34.

Severe Space Weather Events–Understanding Societal and Economic Impacts: A Workshop Report by the National Research Council Committee on the Societal and Economic Impacts of Severe Space Weather Events, published by National Academies Press, Washington, D.C., 2008; available on line: http://www.nap.edu/openbook.php?record_id=12507.

“A Beginner’s Guide to Space Weather and GPS” by P.M. Kintner, Jr., October 31, 2006; available on line: http://gps.ece.cornell.edu/SpaceWeatherIntro_ed2_10-31-06_ed.pdf.

“Combating the Perfect Storm: Improving Marine Differential GPS Accuracy with a Wide-Area Network” by S. Skone, R. Yousuf, and A. Coster in GPS World, Vol. 15, No. 10, October 2004, pp. 31–38.

“Space Weather: Monitoring the Ionosphere with GPS” by A. Coster, J. Foster, and P. Erickson in GPS World, Vol. 14, No. 5, May 2003, pp. 42–49.

The High-Latitude Ionosphere and its Effects on Radio Propagation by R.D. Hunsucker and J.K. Hargreaves, published by Cambridge University Press, Cambridge, U.K., 2002.

“GPS, the Ionosphere, and the Solar Maximum” by R.B. Langley in GPS World, Vol. 11, No. 7, July 2000, pp. 44–49.

• The Effects of the Halloween Storms on GNSS
“Impact of the Halloween 2003 Ionospheric Storm on Kinematic GPS Positioning in Europe” by N. Bergeot, C. Bruyninx, P. Defraigne, S. Pireaux, J. Legrand, E. Pottiaux, and Q. Baire in GPS Solutions, Online First, 2010, doi: 10.1007/s10291-010-0181-9.

“Assessment of EGNOS Performance Under Worst-Case Ionospheric Conditions (Solar Storm of October/November 2003)” by C. Montefusco, J. Ventura-Traveset, B. Arbesser-Rastburg, F. Froment, D. Flament, E. Tapias, S. Radicella, and R. Leitinger in EGNOS – The European Geostationary Navigation Overlay System – A Cornerstone of Galileo, ESA SP-1303, published by the European Space Agency Publications Division, Noorwijk, The Netherlands, 2006, pp. 259–268.

Halloween Space Weather Storms of 2003 by M. Weaver, W. Murtagh, C. Balch, D. Biesecker, L. Combs, M. Crown, K. Doggett, J. Kunches, H. Singer, and D. Zezula, NOAA Technical Memorandum OAR SEC-88, published by the Space Environment Center, National Oceanic and Atmospheric Administration, Office of Oceanic and Atmospheric Research, Boulder, Colorado, June 2004; available on line: http://www.swpc.noaa.gov/Services/HalloweenStorms_assessment.pdf

• Ionospheric Models and Corrections
“Ionospheric Delay Corrections for Single-Frequency GPS Receivers over Europe Using Tomographic Mapping” by D.J. Allain and C.N. Mitchell in GPS Solutions, Vol. 13, No. 2, 2009, pp. 141–151, doi: 10.1007/s10291-008-0107-y.

“Good, Better, Best: Expanding the Wide Area Augmentation System” by T.R. Schempp in GPS World, Vol. 19, No. 1, January 2008, pp. 62–67.

“Ionospheric Time-Delay Algorithm for Single-Frequency GPS Users” by J.A. Klobuchar in IEEE Transactions on Aerospace and Electronic Systems, Vol. AES-23, No. 3, May 1987, pp. 325–331, doi: 10.1109/TAES.1987.310829.

“Global Ionosphere Maps Produced by CODE” on the website of the Astronomical Institute of the University of Bern, Bern, Switzerland: http://aiuws.unibe.ch/ionosphere/.

• Solar Cycle and Solar Weather Predictions:
“Solar Weather Event Modelling and Prediction” by M. Messerotti, F. Zuccarello, S.L. Guglielmino, V. Bothmer, J. Lilensten, G. Noci, M. Storini, and H. Lundstedt in Space Science Reviews, Vol. 147, 2009, pp. 121–185, doi: 10.1007/s11214-009-9574-x.

“Predicting Solar Cycle 24 and Beyond” by M.A. Clilverd, E. Clarke, T. Ulich, H. Rishbeth, and M.J. Jarvis in Space Weather, Vol. 4, S09005, 2006, doi: 10.1029/2005SW000207.

• Current Space Weather and Warnings
European Space Agency Space Weather Web Server, http://www.esa-spaceweather.net/spweather/current_sw/index.html

National Weather Service Space Weather Prediction Center, http://www.swpc.noaa.gov/

Swedish Institute of Space Physics (Institutet för rymdfysik) “Today’s and Recent Space Weather,” http://www.lund.irf.se/HeliosHome/spwfo.html

SpaceWeather.com – News and Information About the Sun-Earth Environment, http://www.spaceweather.com/

February 1, 2011
Innovation: The Distress Alerting Satellite System
Taking the Search out of Search and Rescue

By David W. Affens, Roy Dreibelbis, James E. Mentall, and George Theodorakos

In 1997, a Canadian government study determined that an improved search and rescue system would be one based on medium-Earth orbit satellites, which can provide full global coverage, can determine beacon location, and would need fewer ground stations. This month’s column examines the architecture of the GPS-based Distress Alerting Satellite System and takes a look at early test results.

INNOVATION INSIGHTS by Richard Langley

IT IS NOT COMMONLY KNOWN that the GPS satellites carry more than just navigation payloads. Beginning with the launch of the sixth Block I satellite in 1980, GPS satellites have carried sensors for the detection of nuclear weapons detonations to help monitor compliance with the Non-Proliferation Treaty. The payload is known as the Nuclear Detonation (NUDET) Detection System (NDS) and is jointly supported by the U.S. Air Force and the Department of Energy.

And now a third task is being assigned to the GPS satellites — that of search and rescue. Since the mid-1980s, a combination of low Earth orbit (LEO) and geostationary orbit (GEO) satellites have been used to detect and locate radio beacons activated by mariners, aviators, and others in distress virtually anywhere in the world and at any time. Some 28,000 lives have been saved worldwide since the search and rescue satellite-aided tracking, or SARSAT, system was implemented.

But the current system has some drawbacks. LEO satellites can determine a beacon’s position using the Doppler effect but their field-of-view is limited and one of them may not be in range when a beacon is activated. Furthermore, a large number of ground stations is needed to relay data from these satellites to search and rescue authorities. GEO satellites, on the other hand, have a large field of view (although missing parts of the Arctic and Antarctic), but they cannot position a beacon unless its signal contains location information provided by an integral satellite navigation receiver.

In 1997, a Canadian government study determined that a better SARSAT system would be one based on medium Earth orbit (MEO) satellites. A MEO system can provide full global coverage, determine beacon location, and do this with fewer ground stations. GPS was identified as the ideal MEO constellation.
And so was born the Distress Alerting Satellite System (DASS) that will become fully operational on Block III satellites. But already nine GPS satellites are hosting prototype hardware that is being used for proof-of-concept testing.

In this month’s column, we examine the architecture of DASS (including its relationship with the NDS), and take a look at some of the very positive test results already obtained — results that support the claim that DASS will take the search out of search and rescue.

NASA, which pioneered the technology used for the satellite-aided search and rescue capability that has saved thousands of lives worldwide since its inception nearly three decades ago, has developed new technology that will more quickly identify the locations of people in distress and reduce the risk to rescuers.

The Search and Rescue (SAR) Mission Office at the NASA Goddard Space Flight Center, in collaboration with several government agencies, has developed a next-generation satellite-aided search and rescue system, called the Distress Alerting Satellite System (DASS). NASA, the National Oceanic and Atmospheric Administration (NOAA), the U.S. Air Force, the U.S. Coast Guard, and other agencies are now completing the development and testing of the new system and expect to make it operational in the coming years after a complete constellation of DASS-equipped satellites is launched.

When completed, DASS will be able to almost instantaneously detect and locate distress signals generated by emergency beacons installed on aircraft and maritime vessels or carried by individuals, greatly enhancing the international community’s ability to rescue people in distress, This improved capability is made possible because the satellite-based instruments used to relay the emergency signals will be installed on the GPS satellites.

A recent satellite-aided rescue started on June 10, 2010, when 16-year-old Abby Sunderland on her 40-foot (12.2-meter) sailboat “Wild Eyes” encountered heavy seas approximately 2,000 miles (3,200 kilometers) west of Australia in the Indian Ocean. Her sailboat was dismasted and an emergency situation resulted. Ms. Sunderland activated her two emergency beacons whose signals were picked up by orbiting satellites. Using coordinates derived from the signals, a search plane spotted Ms. Sunderland the next day, and a day later she was rescued by a fishing boat directed to the scene. This highly publicized event is one of thousands of successful rescues made possible by years of NASA research and development.

Background

The beginnings of satellite-aided search and rescue date back to 1970, when a plane carrying two U.S congressmen crashed in a remote region of Alaska. A massive search and rescue effort was mounted, but to this day, no trace of them or their aircraft has ever been found. At the time, search for missing aircraft was conducted by search aircraft flying over thousands of square kilometers hoping to sight the missing aircraft. As a result of this tragedy, Congress recognized this inefficient search method and passed an amendment to the Occupational Safety and Health Act of 1970 requiring most aircraft flying in the United States to carry emergency locator beacons (ELTs) to provide a local homing capability. NASA then developed the technology to detect and locate an ELT from ground stations using the beacon signal relayed by satellites to provide more global coverage. This concept evolved into a highly successful international search and rescue system called COSPAS-SARSAT (COSPAS is an acronym for the Russian words “Cosmicheskaya Sistema Poiska Avariynyh Sudov,” which translates to “Space System for the Search of Vessels in Distress;” SARSAT is an acronym for Search and Rescue Satellite-Aided Tracking). Established by Canada, France, the United States, and the former Soviet Union in 1979, the system has 43 participating countries and has been instrumental in saving more than 28,000 lives worldwide, including 6,400 in the U.S. — all as a result of NASA’s innovations.

Since this auspicious beginning, NASA has continued to perform SAR research and development as a member of the National Search and Rescue Committee, and supports the National Search and Rescue Plan through an interagency memorandum of understanding with the Coast Guard, the Air Force, and NOAA. NOAA is responsible for operation of the U.S. portion of current COSPAS-SARSAT system that relies on SAR payloads on weather satellites in low-earth and geostationary orbits. As shown in Figure 1, the satellites relay distress signals from emergency beacons to a network of ground stations and ultimately to the U.S. Mission Control Center (USMCC) operated by NOAA. The USMCC distributes the alerts to the appropriate search and rescue authorities: the U.S. Air Force or the Coast Guard. The Air Force coordinates search and rescue for the mainland U.S. SAR region and operates the Air Force Rescue Coordination Center. The Coast Guard performs maritime search and rescue and oversees the U.S. national SAR policy.

FIGURE 1. Overall concept of search and rescue system. (Image: Cospas-Sarsat)

Beacons

Three types of distress emergency locator beacons are in use that are compatible with the COSPAS-SARSAT system:
- EPIRBs (emergency position-indicating radio beacons) designed for maritime use.
- ELTs (emergency locator transmitters) for use on aircraft.
- PLBs (personal locator beacons) for personal use. These can be used by persons engaged in high-risk activities such as mountain climbing and backcountry skiing.
Originally, emergency locator beacons transmitted an analog signal on two frequencies: 121.5 MHz and 243 MHz in the civil and military aeronautical communications bands, respectively, so that they would be audible over aircraft radios. Later, a signal that was encoded with a digital message and transmitted at 406 MHz was added. Since February 1, 2009, only the 406-MHz-encoded signals are relayed by satellites supporting the international COSPAS-SARSAT system. Therefore, older beacons that only transmit the 121.5/243-MHz signals are now only detectable by ground-based receivers and aircraft overflying a crash site.

The 406-MHz beacons transmit an approximately half-second message, or burst, approximately every 50 seconds, beginning 50 seconds after being activated. The actual time of burst transmission is dithered in time so that no two beacons will have all of their bursts coincident. A 406-MHz beacon may also have an integral global navigation satellite system (GNSS) receiver. Such a beacon uses the GNSS receiver to attempt to determine its location for inclusion in the transmitted digital message. In this way, the beacon will be located once it is detected by a low-Earth-orbit (LEO) or geostationary orbit (GEO) satellite.

Distress messages contain information such as:
- The beacon’s country of origin.
- A unique 15-digit hexadecimal beacon ID.
- Location, when equipped with an integrated GNSS receiver.
- Whether or not the beacon contains a 121.5-MHz homing signal.
Room for Improvement

SARSAT first became operational in the mid-1980s. The current system uses instruments placed on LEO and GEO weather satellites to detect and locate mariners, aviators, and recreational enthusiasts in distress almost anywhere in the world at anytime and in almost any condition. Previously, dedicated Russian LEO satellites were also implemented but the use of these satellites was discontinued in 2007.

Although it has proven its effectiveness, as evidenced by the number of persons rescued over the system’s lifetime, the current capability does have limitations. LEO spacecraft orbit the Earth 14 times a day and use the Doppler effect with satellite orbital ephemeris data to calculate the position of a beacon. However, a satellite may not be in a position to pick up a distress signal the moment a user activates the beacon. Time is critical in responding to an emergency situation. Unfortunately, delays of two hours or longer are possible, especially near the equator.

LEO spacecraft carry two instruments: a Search and Rescue Repeater (SARR) supplied by the Canadian Department of National Defence, and a Search and Rescue Processor (SARP) provided by the French Centre National d’Etudes Spatiales (CNES). The SARR is a pure repeater, which relays the beacon signal to a local ground station where the data is analyzed to obtain a location. The SARP processes the received beacon signal by measuring the Doppler shift as a function of time, and decoding the digital message included in the 406-MHz signal. This information is stored until it can be transmitted to a ground station using the SARR’s downlink transmitter. Under most conditions beacon locations can be determined to within a radius of 5 kilometers.

Geostationary weather satellites, on the other hand, orbit above the Earth in a fixed location over the equator. Although they do provide continuous visibility of much of the Earth, they cannot independently locate a beacon unless it contains a GNSS receiver that determines its position and includes it in the beacon’s digital message. Currently, not all beacons contain integral GNSS receivers. Furthermore, even if a beacon contains a GNSS receiver, the navigation signal may be obstructed by terrain or thick foliage.

The next-generation system, DASS, overcomes these limitations and will improve accuracy and response time to provide an even more capable life-saving system.

Distress Alerting Satellite System

A 1997 Canadian government study of possible alternative satellite systems for SARSAT, including commercial sources, determined that the ideal system is based on medium Earth orbit (MEO) satellites. A MEO system will be able to provide superior global detection and location data with fewer ground stations than the existing COSPAS-SARSAT system. The GPS constellation was identified as an ideal MEO platform.

The concept of the DASS system is straightforward. Three or more antennas track different GPS satellites equipped with search and rescue repeaters that receive the distress signal and retransmit the signal to the ground. Since each satellite is in a different orbit, each received signal has a different Doppler-shifted arrival frequency and time of arrival. Knowing the position and orbit of each satellite, it is possible to determine the position of the distress beacon.

Future improvement in location accuracy is made possible by one of the strengths of the DASS space segment. That is, the DASS location algorithm optimizes location accuracy utilizing time and frequency measurements of beacon signals that were not designed for that purpose. The DASS space segment allows for the beacon signal to be modified in the future, enhancing the performance of this type of location process.

Other advantages of DASS over the existing system are fairly obvious. Reception of the emergency signal is immediate. Locations can be determined after receiving a single beacon burst since it does not rely on measuring the Doppler shift over time to determine position, as in the current LEO system. A full constellation of DASS-equipped GPS satellites in orbit will ensure that four or more satellites are in view of the transmitting emergency beacon anywhere in the world while requiring fewer ground stations.

Another key strength of the DASS system is the promise of SARSAT transponders on each satellite in the large and well-managed GPS constellation. There are at least 24 GPS active satellites in orbit at any given time (currently, 31 are active). When the GPS constellation is fully populated by satellites with DASS transponders, it will provide global coverage for satellite-supported search and rescue and provide capabilities for rapid detection and location of distress beacons.

Efforts are ongoing to integrate a satellite beacon repeater instrument, to be provided by the Canadian government, onto the GPS Block III B and C satellites to provide the DASS space segment for operational use.

DASS Development

DASS development will proceed in phases referred to as the definition and development, proof of concept, demonstration and evaluation, initial operating capability, and final operating capability. The proof of concept (POC) phase was completed in January 2009. The POC testing and results are summarized in this article. At the time of this writing, preparations are ongoing to initiate the demonstration and evaluation phase.

Definition and Development. In 2000, as part of the definition and development phase, the NASA GSFC SAR Mission Office began discussions with the Department of Energy’s Sandia National Laboratories (SNL) to determine if it would be feasible to add a SAR repeater function to a Department of Energy (DOE) instrument on GPS satellites. Sandia representatives thought it possible, and NASA agreed to fund a study to determine if, with minor modification, one could include a search and rescue repeater function to their instrument. The SNL feasibility study concluded that the GPS DOE package could, with minor modifications, perform the SAR mission. The study also determined that accurate locations could be calculated after a single beacon transmission and improved with each subsequent beacon transmission. Based on this information, NASA, with the cooperation of the U.S. Air Force Space Command and SNL, proceeded with the development of the new space-based search and rescue system, which was named the Distress Alerting Satellite System.

Proof of Concept. In 2003, a memorandum of agreement (MOA) between NASA, NOAA, the Air Force, the Coast Guard, and the Department of Energy tasked NASA to perform a POC program for DASS. The MOA included the development of a POC space segment and a prototype ground station to perform post-launch checkout, performance testing, and implementation planning of an operational DASS system. It stressed the need for DASS, gave authority to each participating agency to participate in the POC demonstration, and defined the roles of each.

The Air Force Space Command approved the addition of modified equipment on GPS satellites. The DASS POC space segment operates as a subcomponent of GPS Block IIR and IIF satellites. Nine GPS Block IIR satellites carry experimental DASS payloads, and all 12 IIF satellites are scheduled to. Therefore, the final POC space segment will consist of 21 DASS-equipped GPS satellites. Each payload receives 406-MHz SAR signals on an extant GPS UHF antenna and relays the signals at a GPS S-band frequency on a second extant antenna.

It is important to note that the performance of the DASS POC space segment will be exceeded by the performance of the operational space segment being designed specifically for DASS and planned for launch on GPS Block III satellites.

A prototype DASS ground station (Figure 2) was funded by NASA and installed at GSFC. The DASS prototype ground system consists of four antennas, four receivers, and the workstations and servers necessary to process the received data, command and control the operation of the ground station, and display and analyze the results. The antennas are located on the corners of the roof of a building connected by fiber-optic cable to signal processing equipment located in another building two kilometers away.

FIGURE 2. Prototype ground station at NASA GSFC. (Images: NASA)

Proof of Concept Testing

The overall objectives of the POC tests were to demonstrate the effectiveness of the DASS concept and to define its technical and operational characteristics. The primary technical objective was to demonstrate the system’s ability to detect and locate 406-MHz emergency beacons under various controlled conditions. This is the most important measure of the system’s ability to perform as expected.

The specific objectives of the DASS POC demonstration were to
- Confirm the expected performance of the DASS concept.
- Determine if new or enhanced requirements needed to be established.
- Define preliminary performance levels that will be used to establish the scope and content of the next phase of development, referred to as the demonstration and evaluation phase.
Therefore, during POC testing, performance measurements were taken for the probability of detection, probability of location, and location accuracy, defined as follows.
- Probability of detection is the probability of detecting the transmission of a 406-MHz beacon and recovering a valid beacon message from any available satellite.
- Probability of location is the probability of obtaining a location solution within a given time after beacon activation, independently of any encoded position data in the 406-MHz beacon message.
- Location accuracy is the distance from the location solution obtained within 5 minutes after beacon activation, to the actual beacon location. The required performance is specified as the probability that a given solution is within a given distance of the actual location.
It is important to note that the predicted performance of DASS assumes a full constellation of DASS-equipped GPS satellites. In fact, one of the key strengths of DASS is the promise of DASS transponders on each satellite in the GPS constellation. When a full constellation is equipped with DASS transponders, there will typically be between seven and 13 GPS satellites visible at the NASA ground station. Thus, it will be possible to schedule the ground-station antennas to receive data from the best satellites in terms of geometry, signal strength, processing capability, and other factors.

However, at the time of the POC testing, there were only eight GPS satellites equipped with DASS transponders. A maximum of three DASS-equipped GPS satellites were visible at the same time at the NASA ground station (above a 15-degree elevation angle), and there were times when only one DASS-equipped GPS satellite was visible. Thus, it was impossible to optimize satellite selection since there was never an opportunity to select from an excess of satellites that a full constellation would provide.

In particular, satellite geometry and its effect on performance is never as optimal as what would be obtained from a full constellation of GPS satellites. To predict the results of a full constellation using the results from a severely reduced constellation, a calculation based on “dilution of precision” was used.

Dilution of precision (DOP) or geometric dilution of precision, to be specific, is used to describe the geometric strength of satellite configuration on GPS accuracy. When visible satellites are close together in the sky, the geometry is said to be weak and the DOP value is high; when far apart, the geometry is strong and the DOP value is low. Thus a low DOP value gives rise to a better GPS positional accuracy due to the wider angular separation between the satellites used to calculate a beacon’s position.

Location accuracy results can be scaled to reflect the true DOP that would be obtained by a satellite constellation of 24 GPS satellites. The DOP error caused by uncertainty in time and frequency measurements is used for scaling. The DOP of the satellites actually used to calculate a location solution, denoted by ftDOP_ACT, is always bigger than the DOP that would have been available from a constellation of 24 GPS satellites, ftDOP₂₄. The raw location errors need to be multiplied by the ratio ftDOP₂₄ / ftDOP_ACT to reflect the results that would have been obtained if all 24 satellites were present.

The raw average location error, err_avg, is given by the following:

err(b) = err(lat(b),lon(b))= distance from the known location to (lat(b),lon(b))

err_avg(b₀) = err(lat_avg(b₀),lon_avg(b₀))

where Ω(b₀) is the set of seven or fewer consecutive burst locations within 5 minutes, starting with burst b₀.

The scaled location error is the location error scaled by the DOP ratio:

Since DOP changes little over 5 minutes, the error of the average is approximately

where ftDOP_ACT(b) is the time-frequency DOP of burst b calculated with either three or four satellite geometries depending on
the number of measurements used in the location calculation.

Test Source

A custom-designed beacon simulator was used to generate the transmissions of multiple COSPAS-SARSAT 406-MHz beacons over an extended period of time. To represent expected operational realism in the tests, the beacon simulator was used to transmit beacons at the limits of the five major beacon parameters specified by COSPAS-SARSAT as well as the nominal values. The five major beacon parameters are transmit power, modulation index, bit rate, un-modulated carrier duration, and modulation rise and fall times (see TABLE 1).

Table 1. Cospas-Sarsat beacon specifications. (Data: Cospas-Sarsat)

During POC testing, five beacons were transmitted using three scenarios: maximum beacon parameter values, minimum beacon parameter values, and variable power. The parameter values changed in each test scenario and are highlighted in TABLE 2. Beacon detection and location performance is measured for periods when there are three or more satellites visible at the same time, and for durations sufficient to collect a statistically significant amount of data.

Table 2. Beacon parameter values for each test scenario. (Data: Authors)

Two characteristics of the test source that affect system performance are the beacon antenna pattern and ground mask. To simulate beacons, the beacon simulator has a monopole antenna with the gain pattern shown in Figure 3. There is a substantial reduction in the transmitted signal at high-elevation angles (above 60°). DASS-equipped GPS satellites are often at high-elevation angles during a typical day. As expected, the effect of the pattern on test results can clearly be seen upon close inspection of the data. However, the beacon antenna pattern is an unavoidable reality and is, therefore, fully represented in the data used to generate the results presented here. Additionally, there were significant ground obstructions of the beacon signal in certain directions. The effect of beacon antenna pattern is fully included in the results presented in this article, but ground mask is taken into account by limiting satellite visibility to an elevation cut-off angle of 15 degrees.

FIGURE 3. Beacon simulator transmit antenna gain pattern.

POC Test Results

In this section, we discuss the POC test results in terms of probability of detection, probability of location, and location accuracy.

Probability of Detection. As previously mentioned, probability of detection is the probability of detecting the transmission of a 406-MHz beacon and recovering a valid beacon message from any available satellite. The requirement is that 95 percent of individual transmitted messages are detected.

Test results are given in TABLE 3 and show that the probability of detection is approximately 99 percent for all scenarios, even though only three satellites were in view at a time. Obviously, the probability of detection is dependent on the number of available satellites and performance would improve with continuous coverage by four or more satellites.

Table 3. Probability of detection test results. (Data: Authors)

Probability of Location. Again, the probability of location is the probability of obtaining a location solution within a given time after beacon activation, independently of any encoded position data in the 406-MHz beacon message. The requirement is that the probability of calculating a beacon location is 98 percent within 5 minutes.

Since the probability of location is dependent on the number of visible satellites, our performance was limited by the reduced constellation of DASS-equipped satellites. Results from periods of three-satellite coverage were 85 percent within 5 minutes, 92 percent within 10 minutes, and 94 percent within 15 minutes.

Again, the probability of location is dependent on the number of visible satellites, and performance would improve with continuous coverage by four or more satellites. To investigate the possible improvement with enhanced satellite coverage, we reduced the minimum satellite elevation angle from 15 to 10 degrees. This allowed a fourth satellite to become visible for a limited time at very low elevation angles. Even though the signal quality from such a satellite was poor, the probability of location during this period of four-satellite coverage improved as follows: 91 percent within 5 minutes, 96 percent within 10 minutes, and 97 percent within 15 minutes.

As can be seen from these results, even adding a satellite with a very low elevation-angle pass significantly improves performance. The expectation is that having a full constellation of satellites available would improve performance even more. Furthermore, the increase in satellite performance expected in the operational system will also improve probabilities of detection and location.

Location Accuracy. Recall that location accuracy is measured as the percentage of location solutions obtained within five minutes after beacon activation that are within five kilometers of the actual beacon location.

The requirement is to obtain 95 percent of the locations to within 5 kilometers of the actual location and 98 percent within 10 kilometers within five minutes after beacon activation.

As mentioned earlier, the requirements included in the performance specification assume a constellation of 24 DASS-equipped GPS satellites. POC testing was done with a system that had only eight DASS-equipped GPS satellites available. However, location errors can be scaled to reflect what the DOP would be if the satellite constellation contained all 24 GPS satellites. Therefore, it is the scaled results that can be used to determine whether performance will meet the requirement.

TABLE 4, therefore, presents the location accuracy results as measured, and after being scaled by DOP.

Table 4. Location accuracy for 5-minute periods. (Data: Authors)

Another important performance metric for DASS is location accuracy obtained after a single beacon burst is received. Even though there is not currently a requirement for single burst location accuracy, it is a very desirable feature of DASS since an emergency situation does not guarantee that more than a single burst will be received. Single burst location accuracy was, therefore, measured with the results shown in TABLE 5. Once again, the results are scaled by DOP values to remove the effect of non-optimal satellite geometry.

Table 5. Single burst location accuracy. (Data: Authors)

More insight into this performance can be gained by examining the single burst location accuracy distribution as a function of distance error, as shown in TABLE 6. It can be seen that, for these beacons, computed locations are within 9 kilometers of the actual location 95 percent of the time. Again, the expectation is that having a full constellation of satellites available would improve this performance. For instance, having more satellites to choose from might allow the system to select data from satellites with stronger or less noisy links.

Table 6. Single burst location accuracy by distance error. (Data authors)

Conclusion

The promise of search and rescue instruments on each satellite in the large and well-managed GPS constellation will provide a significant advancement in the capabilities of the already highly successful COSPAS-SARSAT system. The new system will provide global coverage for satellite-supported search and rescue and provide capabilities for rapid detection and location of distress beacons while requiring fewer ground stations.

The DASS POC system has validated, by test, the predictions made by analysis during the definition and development phase. The DASS POC testing has demonstrated reliable detection and accurate location of beacons within five minutes of activation. Accurate locations are also produced after even a single burst of a newly activated beacon, which is a desirable feature of DASS, since an emergency situation does not guarantee that more than a single burst will be received.

The performance obtained using a reduced constellation of satellites equipped with a modified, existing instrument not only demonstrates the existing capability, but also confirms the improvements to come with the operational system. In fact, the success of DASS is being emulated by the European Union in the design of their future Galileo GNSS constellation and the Russians in an upgraded GLONASS GNSS constellation, all of which will be interoperable by international agreement.

DASS will contribute to NASA’s goal of taking the search out of search and rescue. Achieving this goal will not only improve the chances of rescuing people in distress quickly, which is critical to their survival; it will also reduce the risk to rescuers who often put themselves in dangerous situations to affect a rescue. That is why the motto of the Search and Rescue Office is “Saving more lives, reducing risks to search personnel, and saving resources.”

David W. Affens is the manager of the NASA Search and Rescue (SAR) Mission Office at the Goddard Space Flight Center (GSFC) in Greenbelt, Maryland, where he began working in 1990. He holds a degree in electronic engineering. Before joining NASA, he worked in various aspects of submarine warfare and intelligence gathering for the U.S. Navy over a span of 21 years.

Roy Dreibelbis is a consultant who has worked in rescue-related jobs since 1957, including helicopter rescue missions in Vietnam. As an officer in the U.S. Air Force, he was the director of Inland SAR at rescue headquarters for the coterminous 48 states, was commander of the 33rd Air Rescue Squadron, and served as deputy chief of staff for rescue operations at rescue headquarters from 1979 until 1981. Upon retirement from the Air Force, he was employed by the State of Louisiana as flight operations director and chief pilot. In 1987, he accepted employment with contractors in the District of Columbia area that supported NASA and NOAA SARSAT activities.

James E. Mentall is the NASA/GSFC Search and Rescue Instrument Manager. He has a Ph.D. in physics and has spent more than 42 years of his professional life at GSFC. For 15 of those years, he has been responsible for the integration and test of the Search and Rescue Repeater and the Search and Rescue Processor on the NOAA Polar-orbiting Operational Weather Satellites. He has also served as the deputy mission manager for the Search and Rescue Mission Office and played a significant role in the procurement of the DASS antenna system and ground station.

George Theodorakos is the chief staff engineer for MEI Technologies, Inc. He received his B.S. summa cum laude and M.S. degrees in electrical engineering from the University of Maryland, College Park, Maryland, in 1978 and 1987, respectively. Since 2002, in his role as chief staff engineer at MEI, he has provided technical management support to the Search and Rescue Mission Office at GSFC.

FURTHER READING

• Distress Alerting Satellite System (DASS)
“Distress Alerting Satellite System (DASS)” on the NASA Search and Rescue Mission Office website, Goddard Space Flight Center, Greenbelt, Maryland.

• Search and Rescue Satellite-Aided Tracking (SARSAT)
“Search and Rescue,” Chapter 6 in Review of the Space Communications Program of NASA’s Space Operations Mission Directorate by the Committee to Review NASA’s Space Communications Program, Aeronautics and Space Engineering Board, Division on Engineering and Physical Sciences, National Research Council, published by the National Academies Press, Washington, D.C., 2007.

National Search and Rescue Plan of the United States, authored on behalf of the National Search and Rescue Committee by the United States Coast Guard, Washington, D.C.

• Medium Earth Orbit Search and Rescue (MEOSAR) Systems
COSPAS-SARSAT 406 MHz MEOSAR Implementation Plan, C/S R.012 Issue 1 —Revision 6 October 2010, COSPAS-SARSAT Secretariat, Montréal, Canada.

“SAR/Galileo Early Service Demonstration & the MEOLUT Terminal” by Indra Espacio, a presentation at Galileo Application Days, Brussels, Belgium, March 3–5 2010.

“Mid-Earth Orbiting Search and Rescue (MEOSAR) Transition to Operations” by C. O’Connors, a presentation at the Rescue Coordination Centers Controller Conference, Suitland, Maryland, February 23–25, 2010.

“Overview of MEOSAR System Status” by J. King, a presentation at BMW-2009, Beacon Manufacturers Workshop, St. Pete Beach, May 8, 2009.

“MEOSAR to the Rescue” by J. King in Channels, the EMS SATCOM Quarterly, published by EMS Technologies, Inc., January 31, 2007.

• Nuclear Detonation (NUDET) Detection System
“Detecting Nuclear Detonations with GPS” by P.R. Higbie and N.K. Blocker in GPS World, Vol. 5, No. 2, February 1994, pp. 48–50.
January 1, 2011
Innovation: Record, Replay, Rewind

Testing GNSS Receivers with Record and Playback Techniques

By David A. Hall

Is there a way to perform repeatable tests on GNSS receivers using real signals? This month’s column looks at how to use an RF vector signal analyzer to digitize and record live signals, and then play them back to a GNSS receiver with an RF vector signal generator.

INNOVATION INSIGHTS by Richard Langley

AS A PROFESSOR, I’m quite familiar with testing — of students, that is. It’s how we check their performance — how well they have mastered the course material. Outside academia, testing is also quite common. We have to pass a driving test before we can get a license. We might have to pass a physical fitness test before starting a job. And manufacturers have to test or stress their products to make sure they are fit for purpose. As David Ogilvy, the father of advertising once quipped, “Never stop testing, and your advertising will never stop improving.” But it’s not just manufacturers who should test products. Consumers, or their representatives, should test products on offer — not only to corroborate (or dispute) manufacturers’ claims but also to compare one manufacturer’s product against another. There’s a whole slew of magazines, television programs, and web resources devoted to testing and comparing everything from laundry detergent to automobiles. And GNSS receivers are no exception.

When we conduct tests, we are usually trying to get answers to certain questions — just like those posed to students on their exams. In testing GNSS receivers, what are some appropriate questions? When a receiver is turned on, how long does it take until the position of the receiver is determined? When a weak signal area is encountered, can the receiver still determine its position? If the signal is interrupted and then restored, how long does it take for the receiver to recover and resume calculating its position? And what is the position accuracy under different situations?

While we can certainly hook up an antenna to a receiver to get answers to these questions in a certain environment on a certain day at a certain time with certain signals, the scenario cannot be repeated — not exactly. If we tweak a receiver operating parameter, for example, we don’t know for certain whether any observed change is due to the tweaking or a change in the scenario. We could use a radio-frequency (RF) simulator — a device for mimicking the radio signals generated by the satellites. This would allow us to define scenarios, including receiver trajectories, and to replay them as many times as necessary while varying the operating parameters of the receiver. Or we could modify the scenario from run to run. Such test scenarios could include those difficult to carry out with live signals such as determining how a receiver would perform in low Earth orbit. While extremely useful, these are tests with simulated signals.

Is there a way to perform repeatable tests on GNSS receivers using real signals? In this month’s column, we learn how to use an RF vector signal analyzer to digitize and record live signals, and then play them back to a GNSS receiver with an RF vector signal generator — a procedure we can repeat as often as we like.

While GNSS simulators have long provided the de facto technique for testing GPS receivers, radio frequency (RF) record and playback has emerged as an innovative method to introduce real-world impairments to GNSS receivers. In this article, we will provide a hands-on tutorial on how to test a navigation device using the record and playback technique.

The premise of RF record and playback is to capture GNSS signals off the air with a vector signal analyzer (VSA) and then replay them to a receiver with an RF vector signal generator (VSG). With recorded GNSS signals, one is able to introduce a signal that contains natural impairments — instead of an ideal signal — to the GNSS receiver. As a result, one can observe how a receiver will behave in a real-world environment where interference, multipath fading, and other impairments are present.

A VSA combines traditional superheterodyne radio receiver technology with high-speed analog-to-digital converters and digital signal processors to perform a variety of measurements on complex modulated signals. It is widely used in the telecommunications industry as a test instrument. Digitized signals can be recorded for future analysis. A VSG reverses the process, taking a digital representation of a complex waveform and, using digital-to-analog converters, generating an appropriately modulated RF signal.

Recording GPS or GLONASS signals off the air can be done in a fairly straightforward manner. An RF recording system combines appropriate antennas, amplifiers, and an RF signal recorder (usually a VSA) to capture many hours of continuous RF signal. In such a system, the basic components include the RF front end, the RF signal-acquisition device, and high-volume storage media. A block diagram of a typical recording system is shown in Figure 1.

Figure 1. GPS receivers implement cascaded low-noise amplifiers. The RF signal acquisition block includes analog-to- digital conversion (ADC) and digital down conversion (DDC) to select the data of interest.

In the figure, the RF front end is designed to condition the GNSS signal in such a way that it can be captured — with maximum dynamic range — by the recording device. The recording device digitizes a given signal bandwidth, and then stores in-phase and quadrature (IQ) waveforms to disk.

In general, RF recording devices are designed to tune to a broad range of frequencies and can thereby record many different types of signals. Thus, selecting the signal to record is as simple as setting the center frequency and bandwidth of the recording device. For example, to record the GPS C/A-code L1 signal, the center frequency should be set to 1575.42 MHz. Because each satellite generates the same carrier frequency, one can capture C/A-code signals from all satellites simply by capturing all signals within a 2.046 MHz (twice the code chipping rate) band around the carrier frequency.

By contrast, recording GLONASS signals requires slightly different settings. Because the GLONASS constellation uses frequency division multiplexing, every satellite generates the same code, but each pair of antipodal satellites transmits at a unique center frequency. Thus, recording L1 signal information for the entire GLONASS constellation requires a recorder to capture signals that range from 1598.0625 MHz (channel -7) to 1605.375 MHz (channel 6). In order to capture the entire bandwidth of each satellite, a recorder is actually required to capture 1.022 MHz of signal for each carrier (again, twice the code chipping rate). Therefore, the total recording bandwidth is actually 1597.5515 MHz to 1605.886 MHz, a span of 10.3345 MHz. On the RF signal analyzer, one can record GLONASS signals simply by setting the center frequency to 1601.71875 MHz, and the bandwidth to ≥ 10.3345 MHz.

Modern RF signal recorders are capable of recording both GPS and GLONASS C/A-code signals on a single wideband recording channel. For example, one of our RF signal analyzers is capable of recording up to 50 MHz of signal bandwidth. With this instrument, one can simultaneously record both GPS and GLONASS by setting the center frequency to 1590.1415 MHz and the bandwidth to ≥ 31.489 MHz. However, while RF recording systems can be used to capture a wide range of GNSS signals including GPS L1/L2/L5, GLONASS L1/L2, Galileo, and others, this article focuses primarily on the GPS C/A-code signal.

Setting up the RF Front End

The trickiest aspect of recording GPS signals is the selection and configuration of the appropriate antenna and low noise amplifier (LNA). When connecting a typical off-the-shelf GPS passive patch antenna to a signal analyzer, the peak power in the GPS L1 band ranges from -120 to -110 dBm. Because the power level of GPS signals is small, significant amplification is required to ensure that the VSA can capture the full dynamic range of the signal.

The simplest method to amplify an off-the-air GPS signal so that it can be captured by an RF signal recorder is the combination of an active GPS antenna and one or more external LNAs. Note that many professional GPS antennas offer the best performance because they combine high element gain with an LNA and even pre-selection filtering, which improves the dynamic range of the RF recorder.

With the RF front end appropriately configured, one can verify system performance using a simple spectrum analyzer demonstration panel. The demo panel allows one to visualize the RF spectrum in the GPS L1 band. If all is set up correctly, the C/A-code GPS signal should be visually present on the display. Figure 2 illustrates a screenshot of the spectrum on a virtual spectrum analyzer display.

Note that visualizing the GPS signal in the frequency domain with an RF signal recorder (or spectrum analyzer) requires careful attention to settings such as resolution bandwidth and averaging. Because the signal-to-noise ratio (SNR) of the GPS signal is so small, the settings shown in Figure 2 require a narrow resolution bandwidth (10 Hz) and significant averaging (20 averages per measurement record, so a 20-second interval for 1 Hz data). With these settings applied, one can easily visualize a modulated signal above the noise floor with approximately 1 MHz of bandwidth and centered at 1575.42 MHz. This signal is the GPS C/A-code. In Figure 2, the reference level of the signal analyzer was set to -50 dBm to reduce the noise floor of the instrument to the lowest possible level. Note that setting the signal analyzer’s reference level provides a simple mechanism to adjust the front-end attenuation or amplification. In general, RF signal analyzers provide the greatest dynamic range when the reference level of the instrument matches closely with the average power of the signal connected to the front end. In this case, setting the reference level of our signal analyzer to -50 dBm removes all front-end attenuation, giving the analyzer a more optimal noise figure for signal recording.

Figure 2. GPS is visible in the spectrum only if a narrow resolution bandwidth is used. This spectrum was obtained with a center frequency of 1575.42 MHz, a frequency span of 4 MHz, a resolution bandwidth of 10 Hz, root-mean-square averaging with 20 averages, and a reference level of 250 dBm.

Hardware Connections

With the reference level appropriately set, it is important to properly configure the RF front end of the recording device. As previously mentioned, one can achieve the best RF recording results by using an active GPS antenna. The active antenna used in our experiment utilized a built-in LNA to provide up to 30 dB of gain with a 1.5 dB noise figure. (Recall that the noise figure is the difference in dB between the noise output of a device and the noise output of an “ideal” device with the same gain and bandwidth when it is connected to sources at the standard noise temperature — usually 290 K.) However, the LNA must be powered by supplying a DC bias to the RF connection. While there are several methods to supply the DC bias, we will look at two of the easiest methods.

Method 1: Active Antenna Powered by GPS Receiver. The first method to power an active antenna is with a bias tee or DC power injector. Using this three-port component, a DC voltage (3.3 V in this case) is fed to its DC port, which applies the appropriate DC offset to the active antenna connected to the RF-in port while blocking it on the RF-out port. The device gets its name from the fact that the three ports are often arranged in the shape of a “T.” Note that the precise DC voltage one should apply depends on the DC power requirements of the active antenna. A diagram illustrating the connections is shown in Figure 3.

Observe in Figure 3 that one can use off-the-shelf hardware such as a programmable DC power supply to supply the DC bias signal. Also, one can use a generic off-the-shelf bias tee as long as it has bandwidth up to 1.58 GHz.

Figure 3. This set-up shows the use of a DC bias tee to power an active GPS antenna.

Method 2: Active GPS Antenna Powered by Receiver. A second method of powering the active GPS antenna is with the receiver itself. Most off-the-shelf GPS receivers use a single port to power and receive signals from an active GPS antenna, and this port is already biased with an appropriate DC voltage. Combining an active GPS receiver, a power splitter, and a DC blocker, one can power an active LNA and simply record essentially the same signal as that observed by the GPS receiver. A diagram of the appropriate connections is shown in Figure 4. Some splitters incorporate a DC block on all but one of the output ports.

As Figure 4 illustrates, the DC bias from the GPS receiver is used to power the LNA. This method is particularly useful for drive tests because one can observe the receiver’s characteristics, such as velocity and dilution of precision, while recording.

Figure 4. With a DC blocker, one can record and analyze the same GPS signals being tracked by a GPS receiver.

Selecting the Right LNA

Recording GPS signals with generic RF signal recorders is possible largely because external LNAs can be used to reduce the effective noise floor of the receiver. Today, one can find off-the-shelf spectrum analyzers with noise figures ranging from 15 dB to 20 dB. One of our analyzers, for example, has a 15 dB noise figure while applying up to 60 dB of gain. By applying external amplification to the front of an RF signal analyzer, however, one can substantially reduce the noise figure of the RF recording system.

To calculate the total noise that will be added to the recorded GPS signal, one must calculate the noise figure for the entire RF front end. As a matter of principle, the noise figure of the entire system is always dominated by the first amplifier in the system. Thus, careful selection of the first and second stage LNAs is crucial for a successful signal recording.

We can calculate the noise figure of the RF recording system by using the Friis formula for noise figure, named for engineer Harald Friis, a Danish-American radio engineer who worked at Bell Telephone Laboratories. To use this formula, first convert the gain and noise figure of each component to its linear equivalent; the latter is called the “noise factor.” For cascaded systems such as our RF recording system, the Friis formula provides us with the noise factor of the entire system:

(1)

Note that both noise factor (nf) and gain (g) are shown in lowercase to distinguish them as linear measures rather than logarithmic measures. The conversion from linear to logarithmic gain and noise figure (and vice v
ersa) is shown in the following equations:

An active GPS antenna using a built-in LNA typically provides 30 dB of gain while introducing a noise figure that is typically on the order of 1.5 dB. The second part of the recording instrumentation provides 30 dB of additional gain as well. Though its noise figure is higher (5 dB), the second amplifier actually introduces very little noise into the system. As an academic exercise, one can use the Friis formula to calculate the noise factor for the entire RF front end of the recording instrumentation. Gain and noise figure values are shown in Table 1.

Table 1. Noise figures and factors of the first two components of the RF front end.

According to the calculations above, one can determine the overall noise factor for the receiver:

(6)

To convert noise factor into a noise figure (in dB), apply Equation 2, which yields the following results:

(7)

As Equation 7 illustrates, the noise figure of the first LNA (1.5 dB) dominates the noise figure of the entire RF recording system. Thus, with the VSA configured such that the noise floor of the instrument is less than that of the input stimulus, one’s recording introduces only 1.507 dB of noise to the off-the-air signal.

Saving Data to Disk

Each GNSS produces slightly varying requirements for an RF recorder’s signal bandwidth and center frequency. For the GPS C/A-codes, the essential requirement is to record 2.046 MHz of RF bandwidth at a center frequency of 1575.42 MHz.

In the tests described here, we set the IQ sample rate of our RF recorder at 5 megasamples per second (Ms/s). Since each 16-bit I and Q sample is 32 bits (or 4 bytes each), the actual recording data rate is 20 megabytes per second (MB/s) to ensure the entire bandwidth was captured. Capturing more than 4 MHz of bandwidth is sufficient to record the 2.046 MHz C/A-code signals.

Because one can achieve data rates of 20 MB/s or more with standard PXI controller hard drives (PXI is the open, PC-based platform for test, measurement, and control), one does not need to use an external redundant array of independent disks (RAID) volume to stream GPS signals to disk when using a PXI recording system. In general, data rates exceeding 20 MB/s require the use of an external RAID volume. External RAID systems are capable of storing more than 600 MB/s of data and can be used to support wide bandwidth channels or even multi-channel recording applications. For example, the recording system shown in Figure 5 uses an external RAID volume for high-speed signal recording. This system combines PXI RF signal generators and analyzers with external amplifiers and filter banks for a ready-to-use GNSS record and playback solution.

Figure 5. Two-channel record and playback system from Averna.

In our tests, we decided to use a 320 GB USB drive for better portability. With a disk speed of 5400 revolutions per minute, we were able to benchmark it ahead of time and observed that we were able to achieve read and write speeds exceeding 25 MB/s. Thus, we were easily able to use this disk drive and still record IQ samples at 5 MS/s (20 MB/s) when recording off-the-air signals. With the existing hard-drive setup, we could record more than 4 hours of continuous IQ signal. Note that capturing longer recordings simply requires a larger hard disk. By using a 2 terabyte RAID volume (the largest addressable disk size in the Windows XP operating system), we can extend our recording time to 25 hours. With this setup, we could also reduce the IQ sample rate to 2.5 MS/s (still sufficient to capture the GPS C/A-code signals) and extend the recording time to 50 hours.

Receiver Performance

Once the off-the-air signal of a GNSS band is recorded to disk, it can be re-generated and fed to a receiver using an RF signal generator. With an RF signal generator that is able to reproduce the real-world GNSS signal, engineers are able to test a wide range of receiver characteristics. Because recorded signals contain a rich set of channel impairments such as ionosphere distortion and interference from other transmitters, design engineers often use recorded signals to prototype the baseband processing algorithms on a GNSS receiver.

In our case, we used a VSG directly connected to a GPS evaluation board. In the experiments described below, the receiver’s latitude, longitude, and velocity were tracked over time. Data was read from the receiver using a serial port, which read NMEA 0183 sentences at a rate of one per second. NMEA 0183 is a standard protocol developed by the National Marine Electronics Association for communications between marine electronic devices. NMEA 0183 has been adopted by virtually all GPS receiver manufacturers. In our LabVIEW graphical development environment, one can parse all sentences to return satellite and position-fix information.

For practical testing purposes, GPS dilution of precision and active satellites (GSA), GPS satellites in view (GSV), course over ground and ground speed (VTG), and GPS fix data (GGA) sentences are the most useful. More specifically, one can use information from the GSA sentence to determine whether the receiver has achieved a position fix and is used in time-to-first-fix measurements. When performing sensitivity measurements in this example, the GSV sentence was used to return carrier-to-noise-density ratios (C/N0) for each satellite being tracked. In addition, the VTG sentence allows us to observe the velocity of the receiver. Finally, the GGA sentence provides the receiver’s precise position by returning latitude and longitude information. See the references in Further Reading for in-depth information on the NMEA 0183 protocol.

Using the receiver’s reported latitude and longitude information, we are able to test its ability to report a repeatable position when the recorded signal is played back to the receiver. In this experiment, we tracked the receiver position over 10 minutes. For the best results, the command interface of the receiver should be tightly synchronized with the start trigger of the RF signal generator. The results in Figure 6 show that the RF vector signal generator in this experiment was synchronized with the GPS receiver by using the data line of the serial communications (COM) port (RxD, pin 2) as a start trigger. Using this synchronization method, the vector signal generator and GPS receiver were synchronized to within one clock cycle of the VSG’s arbitrary waveform generator (100 MS/s). Thus, the maximum skew should be limited to 10 microseconds. Given our receiver’s maximum velocity of 15 meters per second (our maximum speed on the drive test), we can determine that the maximum error induced by clock offset of the signal generator is 10 microseconds x 15 meters per second, or 0.15 millimeters.

Using the configuration described above, one is able to report the receiver’s latitude and longitude over time, as shown in Figure 6.

Figure 6A. Receiver latitude over a four-minute span.

Figure 6B. Receiver longitude over a four-minute span.

As the data from Figure 6 illustrate, a recorded test-drive signal reports static, position, and velocity information. In addition, one can observe that this information is relatively repeatable from one trial to the next, as evidenced by the difficulty in graphically observing each individual trace. To better characterize the deviation between each trace, one can also compute the standard deviation between each sample in the waveforms. Figure 7 illustrates the standard deviation between each of the 10 trials, calculated for every one-second interval, versus time.

Figure 7. Standard deviation of both latitude and longitude over time.

When observing the horizontal standard deviation, it is interesting to note that the standard deviation appears to rapidly increase at time = 120 seconds. To investigate this phenomenon further, we can plot the total horizontal standard deviation against the receiver’s velocity and a proxy for C/N₀. In this case, we simply averaged the C/N₀ values for the four highest satellites reported by the receiver. Since four satellites are required to achieve a three-dimensional position fix, our assumption was that position accuracy would closely correlate with the signal strength of these important satellite signals.

One simple method to evaluate the horizontal repeatability of the receiver position versus time is to calculate the standard deviation on a per-sample basis of each recorded latitude and longitude (in degrees). Once the standard deviation is measured in degrees, we can roughly convert this to meters with the following equation:

Note that Equation 8 represents a highly simplified error calculation method, which assumes that the Earth is a perfect sphere. For a more precise calculation of repeatability, the geodesic formula (which presumes that the Earth is ellipsoidal) should be used. In our simple experiment, the goal is merely to correlate repeatability with other factors that we can measure from the receiver. Figure 8 illustrates the standard deviation of horizontal position repeatability over 10 trials and at one-second time intervals.

Figure 8. Correlation of position accuracy and C/N₀.

As one can observe in Figure 8, the peak horizontal error (measured by standard deviation) occurring at time = 120 seconds is directly correlated with satellite C/N₀ and not correlated with receiver velocity. At this sample, the standard deviation is nearly 2 meters while it is less than 1 meter during most other times. Concurrently, the top four C/N₀ averages drop from nearly 45 dB-Hz to 41 dB-Hz.

The exercise above illustrates not only the effect of C/N₀ on position accuracy but also the types of analysis that one can conduct using recorded GPS data. For this experiment, the drive recording of the GPS signal was conducted in Huizhou, China (a city north of Shenzhen), but the actual receiver was tested at a later date in Austin, Texas.

Conclusion

In this article, we’ve illustrated how to use commercially available off-the-shelf products to record GPS signals with an RF recorder, and then play the signal back to a receiver. As the results illustrate, recorded GPS signals can be used to measure a wide range of receiver characteristics. Not only can receiver designers use these test techniques to better prototype a receiver baseband processor, but also to measure system-level performance such as position repeatability.

Manufacturers

The tests discussed in this article used a National Instruments PXIe-5663E, 6.6 GHz, RF signal analyzer; a National Instruments PXI-5690, 100 kHz to 3 GHz, two-channel programmable amplifier and attenuator; a National Instruments PXIe-5672, 2.7 GHz, RF vector signal generator with quadrature digital upconversion; a 320 GB USB Passport hard drive from Western Digital Corp.; a National Instruments PXI-4110 programmable, triple-output, precision DC power supply; and a ZX85-12G-S+ bias tee manufactured by Mini-Circuits. The article also mentioned the RP-3200 2-channel record and playback system manufactured by Averna, which incorporates National Instruments modules.

David Hall is an RF product manager for National Instruments. He holds a bachelor’s of science with honors in computer engineering from Pennsylvania State University.

FURTHER READING

More on GNSS Receiver Record and Playback Testing

GPS Receiver Testing, tutorial published by National Instruments, Austin, Texas.

Friis Formula and Receiver Performance

RF System Design of Transceivers for Wireless Communications by Q. Gu, published by Springer, New York, 2005.

Global Positioning System: Signals, Measurements, and Performance, 2nd edition, by P. Misra and P. Enge, published by Ganga-Jamuna Press, Lincoln, Massachusetts, 2006.

“Measuring GPS Receiver Performance: A New Approach” by S. Gourevitch in GPS World, Vol. 7, No. 10, October 1997, pp. 56-62.

“GPS Receiver System Noise” by R.B. Langley in GPS World, Vol. 8, No. 6, June 1997, pp. 40–45.

Global Positioning System: Theory and Applications, Vol. I, edited by B.W. Parkinson and J.J. Spliker Jr., published by the American Institute of Aeronautics and Astronautics, Inc., Washington, D.C., 1996.

GNSS Receiver Testing Using Simulators

“Testing Multi-GNSS Equipment: Systems, Simulators, and the Production Pyramid” by I. Petrovski, B. Townsend, and T. Ebinuma in Inside GNSS, Vol. 5, No. 5, July/August 2010, pp. 52–61.

“GPS Simulation” by M.B. May in GPS World, Vol. 5, No. 10, October 1994, pp. 51–56.

GNSS Receiver Testing Using Software

“GPS MATLAB Toolbox Review” by A.K. Tetewsky and A. Soltz in GPS World, Vol. 9, No. 10, October 1998, pp. 50–56.

GNSS L1 Signal Descriptions

Navstar GPS Space Segment / Navigation User Interfaces, Interface Specification, IS-GPS-200 Revision E, prepared by Science Applications International Corporation, El Segundo, California, for Global Positioning System Wing, June 2010.

Global Navigation Satellite System GLONASS, Interface Control Document, Navigational Radio Signal in Bands L1, L2 (Edition 5.1), prepared by Russian Institute of Space Device Engineering, Moscow, 2008.

NMEA 0183

NMEA 0183, The Standard for Interfacing Marine Electronic Devices, Ver. 4.00, published by the National Marine Electronics Association, Severna Park, Maryland, November 2008.

“NMEA 0183: A GPS Receiver Interface Standard” by R.B. Langley in GPS World, Vol. 6, No. 7, July 1995, pp. 54–57.

Unofficial online NMEA 0183 descriptions: NMEA data; NMEA Revealed by E.S. Raymond, Ver. 2.3, March 2010.

October 1, 2010
Innovation: Friendly Reflections

Monitoring Water Level with GNSS

A receiver can selectively acquire scattered signals and the resulting measurements can be interpreted to reveal certain characteristics of the source of the scattering. This article discusses the design and application of a GNSS instrument that uses scattered signals for monitoring the level and roughness of inland and coastal water surfaces for the betterment of planet Earth.

By Alejandro Egido and Marco Caparrini

INNOVATION INSIGHTS by Richard Langley

WHY IS THE SKY BLUE? This is an age-old question, interesting to anyone with a curiosity about his or her surroundings. But what has it got to do with global navigation satellite systems? Believe it or not, there is a connection.

Some of you might remember the explanation of the sky’s color from your Physics 101 course but to bring everyone up to the same level, let’s review. Everything we see is the result of the interaction of light and matter. And by matter, we mean the atoms, molecules, and particles making up matter. Light causes matter to vibrate. And vibrating matter (due to its electrical charges) in turn emits light, which combines with the original light. But matter not only re-emits light in the forward direction, it re-emits light in all other directions. This is called scattering.

Now, the light from the sun includes all colors and so if look directly at the sun when it is high in the sky (don’t try this at home), it looks white or slightly yellowish. We are seeing the light propagating directly toward our eyes. When we look at the sky away from the sun, we are seeing scattered light. And this scattered light is predominantly blue. Why? It turns out that scattering is proportional to the fourth power of frequency. Light that is of a higher frequency, say a factor of two, is sixteen times more intensely scattered. So, blue light, which has about twice the frequency of light from the red end of the visible spectrum, is scattered much more than red light. Violet light is scattered even more but our eyes are not as sensitive to violet light as they are to blue light. Hence the sky looks blue.

So what has this got to do with GNSS? As we know, for the best positioning and navigation results, we need the satellite signals to travel along a direct path to the receiver’s antenna. There may be slight changes in the speed and direction of propagation of these direct-path signals caused by the interaction of the electromagnetic waves with the matter making up the ionosphere and the neutral atmosphere, but these are readily accounted for in the position fixes.

However, once they reach the Earth’s surface, the signals can be reflected by buildings, vegetation, the ground, water surfaces, and so on. The signals are actually being scattered by the matter they encounter. A receiver can selectively acquire the scattered signals and the resulting measurements can be interpreted to reveal certain characteristics of the source of the scattering.

In this month’s column, we learn about the design and application of a GNSS instrument that uses scattered signals for monitoring the level and roughness of inland and coastal water surfaces–yet one more use of GNSS signals for the betterment of planet Earth.

Lakes and water reservoirs are the world’s most important sources of accessible fresh water. Despite its paramount importance — not only for a large variety of human activities, but also for the sustainability of ecosystems — fresh water is already scarce in many regions. The problem is envisaged to become worse in the coming decade. In addition, in climatological studies surface water storage is a critical element of the water cycle since the analyses integrate all hydrologic processes (precipitation, runoff, evapotranspiration, and so on) over a given basin; and for hydroelectric companies, it is the main parameter to be kept under observation for efficient energy production. All of these concerns make the monitoring of fresh water resources a prime activity for a wide variety of stakeholders including governments, climate research organizations, and hydroelectric production companies.

Coastal management is also a wide-ranging issue with large social and economic impacts. Care of our coasts includes dealing with threats such as storm surges and flooding, coastal erosion, and conflicting land-use issues. Coastal areas support the greatest concentration of living resources and people on the planet. In the past few decades, these regions have experienced a population density increase, which is envisioned to grow steadily. Furthermore, conflicts between commercial interests, recreational activities, infrastructure development, environment conservation, and exploitation of natural resources will become increasingly important and contentious. In fact, the coastal zone is a peculiar environment in which terrestrial, oceanic, atmospheric, and human inputs of energy and matter converge. Storm surges and coastal flooding events have caused considerable damage and economic loss on European coasts in particular. Such events, possibly linked to the world climate change, are expected to get worse in the near future, due to sea level rise and storm activity.

So, close monitoring of both inland waters and coastal regions is necessary for the well being of the planet. And since the need is so pervasive, monitoring systems should be characterized by a relatively low cost, low maintenance, and easy deployment, to serve the widest possible user community. We have developed a patent-pending solution using signals from global navigation satellite systems (GNSS).

Called Oceanpal, our monitoring system exploits reflected GNSS signals as signals of opportunity for passive remote sensing of the Earth’s water surfaces. These multipath signals are usually considered to be nuisance signals since they reduce the accuracy of GNSS positioning applications. But for monitoring various processes affecting the Earth’s surface, they are very beneficial. The technique is known as GNSS reflectometry (GNSS-R), and during the past decade, its use as a technique for Earth observation purposes has taken root.

GNSS-R is basically a bistatic radar technique. While most radar systems, such as those used for monitoring air space and harbor approaches and for weather forecasting, combine the radar transmitter and receiver at the same site — so-called monostatic radar — bistatic systems use transmitters and receivers separated by a considerable distance. Such systems have been used for studying certain atmospheric phenomena and for military applications where simple line-of-sight reflections from the target of interest are inadequate or insufficient.

The concept of bistatic radar can be extended to satellite signals. Since some of the signal transmitted by a satellite gets reflected off the Earth’s surface, detecting this reflected signal by a separate passive receiver would provide some information about the reflecting surface. While any satellite signal could be used in principle, GPS (and other GNSS) turn out to be particularly useful. The concept of using GPS signal reflections was initially proposed in 1993 by Manuel Martín-Neira, working at the European Space Agency’s European Space Research and Technology Centre in Noordwijk, The Netherlands. Since then, the technique has been successfully implemented by an increasing number of researchers.

We could list several reasons for the continuous growing interest in GNSS as a remote sensing tool, but two main ones stand out: first, the global availability and stability of GNSS signals enables their use as reliable signals of opportunity; and second, GNSS makes use of L-band radiation, which is highly interactive with the natural scattering medium but relatively impervious to atmospheric conditions. Moreover, the passive nature of this concept allows for the production of cost- and resource-effective instruments.

Navigation signals are sensitive to a wide variety of geophysical parameters including topography, surface roughness, surface moisture, ionospheric electron content, tropospheric water vapor, water salinity, and vegetation. Research targeting related geophysical applications has been ongoing for many years, and the first pre-operational services exploiting reflected GNSS signals are now available. In fact, while the scientific community is waiting for a dedicated GNSS-R space mission to confirm the theoretical predictions about the characteristics of reflected signals observed from space, ground-based and airborne sensors have already been developed and validated for a number of applications.

The GNSS-R research area that has been most thoroughly investigated concerns the reflection of navigation signals from water surfaces, given the highly reflective nature of water. However, from water the interest has now moved towards ice and land applications, more specifically to the detection of sea ice and the monitoring of soil moisture. Recently, GNSS-R has also been proposed as a possible tool to monitor vegetation. This article focuses on the presentation of the Oceanpal sensor, and the description of the altimetry algorithms for monitoring the levels of sea (coastal) and inland waters.

Our Instrument

As mentioned above, Oceanpal is a GNSS-R-based sensor designed for operational monitoring of coastal and inland waters. The instrument comprises three subsystems: a radio frequency (RF) section, an intermediate frequency (IF) section, and a data-processing section. The RF section features a pair of low gain L-band antennas. A right-hand circularly polarized (RHCP) zenith-facing antenna collects the direct GNSS signals while a left-hand circularly polarized (LHCP) nadir-facing antenna collects the sea- or lake-surface reflected GNSS signals. (On reflection, the signals become predominantly LHCP.) Data bursts of some minutes’ duration are acquired from each antenna using two GPS L1 receivers (front ends) that down-convert the signals to IF. Within the IF sections, the signals are one-bit sampled and stored on a hard disk.

These direct and reflected raw data are then fed into the processing section of the instrument, where a pair of software GNSS receivers detects and tracks the available signals in the direct channel (which works as a master) and blindly despreads the reflected signals in the reflected (slave) channel. The result of this processing is a set of direct and reflected electromagnetic field time series for each satellite in view, plus some ancillary information, such as the satellite pseudorandom noise code (PRN) numbers and GPS time references, among others. The architecture described above is shown in FIGURE 1.

Figure 1. Basic operation of Oceanpal and the principle of GNSS-R-based sea-surface monitoring. Right-hand and left-hand circularly polarized antennas feed signals to radio frequency (RF) receiver front-ends that, in turn, feed software (SW) receiver back-ends and subsequent processing algorithms.

The data products provided by Oceanpal are so-called “Level-2” or derived products, namely the significant wave height (a statistical measure of trough-to-crest wave height), and the height of the nadir antenna over the mean level of the water surface under observation. To make this data available for the user in a friendly way, the observations are uploaded to a web server and displayed on a web page.

Oceanpal requires low maintenance compared to its competitors. Standard oceanographic buoys, which use accelerometers and a magnetic compass, or GPS buoys, featuring a conventional GPS receiver, are in contact with water, which implies costly infrastructures and frequent maintenance operations. Pressure sensors and air bubblers, commonly used to monitor the level of water reservoirs, also require frequent maintenance because of sediment accumulation. Compared to the alternatives, our sensor is a less costly and lower maintenance solution.

GNSS-R Altimetry Algorithms

The inland-water/sea-level monitoring is based on the estimation of the height of the Oceanpal antennas above the water/sea surface. This height is retrieved by the comparison of the delay (in time or distance) between the reflected and the direct signals. The reflection geometry is shown in FIGURE 2. Such a delay can be estimated using either the PRN code or the carrier phase of the incoming signals. The phase-based estimation provides more precise values, but it is only available for calm water surfaces where coherent constructive scattering (specular reflection) is predominant. In the case of rougher surfaces, the reflected signal’s coherency is lost, and therefore the code-based algorithm must be used.

Figure 2. The geometry of GNSS signal reflections for altimetry applications.

The basic equation that links the delay of arrival of both signals with the height of the antennas over the surface as a function of time (t) can be written as equation (1):

   (1)

where τ represents the lapse between the time of arrival of the reflected and the direct signals (as determined using either phase or code measurements), h is the height to be estimated, e is the elevation angle of the satellite considered, and b is the system bias, which is considered unknown but constant during every estimation. Solving a linear system with many such equations for different satellites over, say, one minute provides the sought estimation of h (and b).

Measuring the Level of a Water Reservoir

As mentioned before, when the water surface is sufficiently flat, the coherency of the reflected signal is maintained, thus its phase can be used to retrieve estimates of the height of the antennas over the surface. This algorithm is the so-called phase altimetry algorithm. The basic observable for this algorithm is the interferometric complex field (ICF), defined as the ratio between the reflected and direct complex correlation waveform peaks:

     (2)

where P_R and P_D represent the time series of waveform peaks for the reflected and direct signals, respectively. In computing this ratio, adverse propagation effects such as the extra delay induced by the ionosphere and troposphere cancel out. Measuring the phase of the ICF, , one is then considering the phase single difference, , between the reflected and direct signals as given in equation (3):

     (3)

where k is the wave number of the GPS carrier frequency (the reciprocal of the wavelength), noise_φ is the noise present in the ICF phase and is the unknown integer cycle ambiguity. D is the excess path of the reflected with respect to the direct signal, which can be directly linked to the height of the antennas over the surface. In order to solve for the cycle ambiguities, phase double differences among satellites are calculated, and by means of an ambiguity resolution algorithm (we use the null-space method developed by Manuel Martín-Neira and colleagues) the unknown phase-cycle ambiguities can be determined. It is then a straightforward procedure to work out the excess path of the reflected signals to finally deduce the height of the antennas over the water surface.

La Baells Experiment. An experimental campaign was carried out with an Oceanpal instrument at the La Baells water reservoir (near Berga in Catalonia, Spain) in cooperation with the Catalan Water Agency. This experiment was designed to study the feasibility of accurate altimetry measurements at lakes and reservoirs using our technique.

Within this campaign, one week of data was gathered early in March 2008 to compare the Oceanpal GNSS-R phase-altimetry measurements with those from the La Baells in-situ sensor (a water bubbler known to have centimeter-level accuracy). The results from this campaign are shown in FIGURE 3. After referencing the measurements to the antennas’ position with respect to the mean water level, the accuracy obtained from the Oceanpal measurements with respect to the ground truth (the water bubbler) was better than 2 centimeters (after a five-minute integration time).

Figure 3. Results of a one-week campaign at the La Baells water reservoir near Berga, Spain, in March 2008. Lake height in meters with respect to mean sea level.

Despite the fact that the phase altimetry algorithm is precise, it requires the simultaneous observation of several reflections from different satellites to converge and accurately solve for the phase ambiguities. However, this cannot be done for all scenarios, and in these situations the conventional phase altimetry algorithm cannot be applied.

Lake Laja Experiment. A case where we couldn’t use the phase approach was project Hydro. This was an initiative developed by our organization in collaboration with Pontificia Universidad Católica de Chile (the Pontifical Catholic University of Chile) and funded by ENDESA (Empresa Nacional de Electricidad S.A.), one of the world’s largest electricity companies. An Oceanpal instrument was installed at Lake Laja, in the Biobío Region, Chile, a water reservoir managed by ENDESA Chile. The Hydro project aims to use remote sensing assets to predict and monitor water flow in the Laja River basin. For that, having precise measurements of Lake Laja’s level is a must.

The instrument was installed on the shore of the lake as seen in FIGURE 4. However, the high variability of the lake’s level, more than 10 meters in one year, and the abruptness of the terrain, results in the number of observed reflections from the water surface being quite low. This is especially the case when the level of the lake is low. In this situation, the number of different GPS satellites observed per hour was calculated to be fewer than two for more than 45 percent of the time, and fewer than three for more than 85 percent of the time. Given this scarcity of reflections, we could not use the phase altimetry algorithm as described above.

Figure 4. The Oceanpal installation at Lake Laja, Chile.

We developed a new phase altimetry algorithm, which considers the interferometric phase evolution over time. The resulting relative phase parameter can be linked to the height of the antennas over the water surface by means of the same geometrical relationship as before. Despite the fact that measuring a relative phase increases the measurement noise with respect to the case in which an absolute phase is used, the phase ambiguity and the bias between the direct and reflected receiving channels do not need to be calculated, thus reducing the complexity of the algorithm and its convergence requirements. A Kalman filter is used to smooth the inherently noisy behavior of the relative phase.

The Oceanpal measurements were compared to those of a sensor operated by the Dirección General de Aguas (DGA), the Chilean water management agency. An accuracy better than 9 centimeters was achieved in determining the lake’s level during the austral winter, when the lake is at its minimum level and therefore the satellites’ reflections from the water surface are scarce. The lake level has its maximum during the summer after the melting season. During this period of time, the achieved accuracy of Oceanpal with the new phase algorithm was better than 5 centimeters. A comparison of Oceanpal and DGA’s sensor measurements of the water level is shown in FIGURE 5.

Figure 5A. A comparison of measurements of Lake Laja’s water level by Oceanpal and a water bubbler sensor operated by Dirección General de Aguas (DGA) for two periods of time corresponding to the austral winter (from late April 2009 until early August 2009).

Figure 5B. A comparison of measurements of Lake Laja’s water level by Oceanpal and a water bubbler sensor operated by Dirección General de Aguas (DGA) for two periods of time corresponding to the austral summer (from late November 2009 until late January 2010).

Measuring Sea Level

Sea level is obtained from Oceanpal measurements by means of the code altimetry algorithm due to the inherent roughness of the sea surface. This technique derives altimetric information from the displacement of reflected waveforms with respect to the direct ones. Such a displacement can be directly related to the delay between the direct and reflected signals (the so-called lapse), and is used in a similar way to the phase-based method to extract the altimetry information of the water surface being monitored.

Despite the fact that the code altimetry algorithm is not as precise as the phase altimetry algorithm, it is not subject to the coherence requirement for the reflected signal. Therefore, it can be applied to rough, dynamic surfaces such as the open ocean and coastal areas. The use of code altimetry in rough water conditions results in a clear observation of tide dynamics but, as expected, with a higher error range compared to situations where phase altimetry can be applied.

Scheveningen Pier Experiment. The performance of the code-based algorithm was tested during an experimental campaign carried out on Scheveningen Pier in Den Haag (The Hague), The Netherlands. An Oceanpal instrument was installed close to a Radac X-band radar tide gauge. FIGURE 6 shows the tide variation at the installation site estimated by the Radac instrument and by Oceanpal. As can be seen, a good agreement between both estimates is achieved with a standard deviation of the difference of 12 centimeters.

Figure 6. Daily tidal variation at Scheveningen Pier, The Hague, The Netherlands, on May 3-4, 2008, measured by X-band radar and Oceanpal.

To improve this result, a combination of code and phase estimation is being investigated, involving the alignment of the phase using the code information. The combination of these two parameters may provide the best of both worlds. However, with the signals from modernized GPS and those of the forthcoming Galileo system, the code-ranging precision is envisioned to increase by a factor of four or five, which is expected to impact directly on the precision of the code altimetry algorithm.

Conclusion and Outlook

During the past decade, the scientific community’s interest in GNSS-R has grown, leading to the continuous development of new applications and to an increasing relevance in specific market niches. Some of these applications, especially those related to the monitoring of water surfaces, have reached an operational level of maturity, and provide end users with valuable information.

In this brief article, we have described the Oceanpal instrument and outlined its use in altimetric measurements of water surfaces. It was shown that using the phase of reflected signals with respect to that of direct signals, accurate measurements of a lake’s level could be obtained. In addition, we overviewed a new algorithm that incorporates the evolution of this phase in time. This algorithm is suitable for low satellite visibility scenarios. For example, using this algorithm, the level of Lake Laja in Chile was determined with an overall accuracy better than 7 centimeters. Such a level of accuracy meets the monitoring requirements necessary for improving the stream-flow prediction in the Laja River basin. We also showed that code altimetry can be successfully used to monitor sea level variations associated with tides, with a demonstrated accuracy of 12 centimeters.

These encouraging results are expected to be further improved with the evolution of GPS, the refurbishment of the Russian GLONASS system, and the deployment of the European Galileo system. First of all, when all three navigation systems are fully deployed, it is calculated that at least 20 navigation satellites will be visible at the same time. A GNSS-R instrument could take advantage of this large number of available signals. In addition, the quality of these signals is expected to be largely improved in terms of signal-to-noise ratio, bandwidth, and ranging precisions, which will in turn improve the performance of GNSS-R altimetry algorithms. As a result, the prospects for GNSS-R altimetry over water surfaces, not only for ground-based systems, but also airborne and even spaceborne systems, are extremely promising.

Manufacturers

The Oceanpal instrument was developed by Starlab, Barcelona, Spain. The Scheveningen Pier experiment used a Radac, Haarlem, The Netherlands, WaveGuide radar level gauge.

ALEJANDRO EGIDO has a B.Sc. degree in electrical engineering from the University of Zaragoza, Spain. After his studies, he worked on the Sentinel-1 remote sensing satellite project at the European Space Agency (ESA), where he performed the interference analysis of the synthetic aperture radar instrument. Since 2007, he has been a research engineer at Starlab, Barcelona, while pursuing a Ph.D. at the Polytechnic University of Catalonia. His main research field is the use of GNSS signals as sources of opportunity for remote sensing applications, with special interest in estimating terrestrial bio-geophysical parameters.

MARCO CAPARRINI received the “Laurea” degree in electronic engineering — remote sensing from the University “La Sapienza” in Rome. He has worked as a research engineer at ESA’s European Space Research and Technology Centre in Noordwijk, The Netherlands; at the German Aerospace Center in Oberpfaffenhofen, Germany; and at the Swiss Federal Institute of Technology in Zurich. His main research field is the use of GNSS signals as sources of opportunity for remote sensing of planet Earth, and he is the Starlab manager for the space research and development area.

FURTHER READING

• Principles of GNSS Reflectometry (GNSS-R)
“The PARIS Concept: An Experimental Demonstration of Sea Surface Altimetry Using GPS Reflected Signals” by M. Martín-Neira, M. Caparrini, J. Font-Rossello, S. Lannelongue, and C. Serra Vallmitjana in IEEE Transactions on Geoscience and Remote Sensing, Vol. 39, No. 1, January 2001, pp. 142–150, doi: 10.1109/36.898676.

• Overview of GNSS-R Applications
“GNSS Reflectometry and Remote Sensing: New Objectives and Results” by J. Shuanggen and A. Komjathy in Advances in Space Research, Vol. 46, 2010, pp. 111–117, doi:10.1016/j.asr.2010.01.014.

• GNSS-R Experimental Campaigns

“Oceanpal: Monitoring Sea State with a GNSS-R Coastal Instrument” by M. Caparrini, A. Egido, F. Soulat, O. Germain, E. Farrès, S. Dunne, and G. Ruffini in Proceedings of the 2007 International Geoscience and Remote Sensing Symposium, Barcelona, Spain, July 23–28, 2007, pp. 5080–5083.

“The Eddy Experiment: Accurate GNSS-R Ocean Altimetry from Low Altitude Aircraft” by G. Ruffini, F. Soulat, M. Caparrini, O. Germain, M. Martín-Neira in Geophysical Research Letters, Vol. 31, L12306, 4 pp., 2004, doi:10.1029/2004GL019994.

“The Eddy Experiment: GNSS-R Speculometry for Directional Sea- Roughness Retrieval from Low Aircraft” by O. Germain, G. Ruffini, F. Soulat, M. Caparrini, B. Chapron, and P. Silvestrin in Geophysical Research Letters, Vol. 31, L21307, 4 pp., 2004, doi: 10.1029/2004GL020991.

“Wind Speed Measurement Using Forward Scattered GPS Signals” by V. Zavorotny, J. Garrison, A. Komjathy, and S. Katzberg in IEEE Transactions on Geoscience and Remote Sensing, Vol. 40, No. 1, January 2002, pp. 50–65, doi: 10.1109/36.981349.

• GNSS-R for Monitoring Soil Moisture

“The SAM Sensor: An Innovative GNSS-R System for Soil Moisture Retrieval” by A. Egido, C. Martin-Puig, D. Felip, M. Garcia, M. Caparrini, E. Farrés, and G. Ruffini in Proceedings of NAVITEC 2008, the 4th ESA Workshop on Satellite Navigation User Equipment Technologies, Noordwijk, The Netherlands, December 10–12, 2008.

• GNSS-R for Ice Detection

“Reflecting on GPS: Sensing Land and Ice from Low Earth Orbit” by S.T. Gleason in GPS World, Vol. 18, No. 10, October 2007, pp. 44–49.

• GNSS-R for Ocean Surface Monitoring

“GPS: A New Tool for Ocean Science” by A. Komjathy, J.L. Garrison, and V. Zavorotny in GPS World, Vol. 10, No. 4, April 1999, pp. 50–56.

• Using Signal-to-Noise Ratio as a Multipath Observable

“It’s Not All Bad: Understanding and Using GNSS Multipath” by Andria Bilich and Kristine M. Larson in GPS World, Vol. 20, No. 10, October 2009, pp. 31–39.

• Carrier-Phase Ambiguity Resolution Techniques

“GPS Ambiguity Resolution and Validation: Methodologies, Trends and Issues” by D. Kim and R.B. Langley in Proceedings of the 7th GNSS Workshop – International Symposium on GPS/GNSS, Seoul, Korea, Nov. 30 – Dec. 2, 2000, Tutorial/Domestic Session, pp. 213–221.

September 1, 2010
Innovation: Better Weather Prediction Using GPS

Water Vapor Tomography in the Swiss Alps

By Simon Lutz, Marc Troller, Donat Perler, Alain Geiger, and Hans-Gert Kahle

A team of Swiss researchers is using data from a network of GPS receivers and the technique of tomography to obtain profiles of how moisture is distributed with height, which might lead to better weather forecasts.

INNOVATION INSIGHTS by Richard Langley

WEATHER FORECASTING IS STILL AN IMPERFECT ART. Humankind has been trying to predict the weather for millennia. Early attempts were based on general observations such as “Red sky at night, shepherd’s delight; Red sky in morning, sailor’s warning.” But it was only with advances in scientific thought and the invention of measuring devices, such as the mercury barometer, that more specific predictions could be made.

Towards the end of the 18th century, the father of modern chemistry, Antoine Laurent Lavoisier, said “It is almost possible to predict one or two days in advance, within a rather broad range of probability, what the weather is going to be; it is even thought that it will not be impossible to publish daily forecasts, which would be very useful to society.”

Forecasting ability has improved over the years as measurement technology, communications, and the understanding of atmospheric processes have improved. Meteorologists use measurements from various types of sensors together with mathematical models describing the physics of the atmosphere to predict its future state. The temporal and spatial density of the measurements and the sophistication of the models have a direct impact on the success of the forecasts. Weather stations on the Earth’s surface, radar installations, radiosondes, and satellite sensors all provide data for modern forecasts. Yet better sampling of the current state of the atmosphere, particularly the distribution of water vapor, is required to produce more accurate and more timely forecasts of its future state. GPS can help.

The signals from the GPS satellites must transit the atmosphere on their way to a receiver on the Earth’s surface. The atmosphere’s atoms and molecules slow down the signals so that they arrive slightly later than they would if the Earth was surrounded by a vacuum, and this effect shows up in the GPS receiver measurements. The receiver or measurement processing software needs to remove or model the effect to obtain accurate receiver positions. On the other hand, if all parameters affecting GPS measurements such as satellite and receiver coordinates are well known, then the delay imparted by the atmosphere can be estimated. It is possible to separate the effect of water vapor from that of the dry gases such as nitrogen, oxygen, and carbon dioxide and to provide a measure of the atmosphere’s moisture content. Several national weather agencies are ingesting such estimates from networks of GPS receivers into experimental or operational numerical weather forecast models. But these values represent an integrated measure of moisture above a receiver. Profiles of how moisture is distributed with height would be more useful and might lead to better weather forecasts. In this month’s column, a team of Swiss researchers discuss how they use data from a network of GPS receivers and the technique of tomography to obtain such profiles.

“Innovation” features discussions about advances in GPS technology, its applications, and the fundamentals of GPS positioning. The column is coordinated by Richard Langley, Department of Geodesy and Geomatics Engineering, University of New Brunswick.

Water vapor plays an essential role in the dynamics and thermodynamics of the atmosphere — especially storm systems — on local, regional, and global scales. It is a precursor of precipitation. Furthermore, a signiﬁcant fraction of the energy released to the atmosphere comes from water vapor via latent heat. And much of the “greenhouse effect” is caused by the presence of water vapor in the atmosphere.

Beginning in 1992, a number of researchers successfully tested the hypothesis that the Global Positioning System (GPS) could be used to detect long- and short-term global and regional air-mass changes by estimating the amount of water vapor in the air above a GPS receiver. The arrival of GPS signals at a receiver is delayed by the presence of the Earth’s atmosphere. The satellite signals slow down when they encounter the atmosphere’s electrons, atoms, and molecules. In particular, the signals are affected by the presence of water vapor. Through a careful analysis of the GPS receiver’s measurements, the amount of water vapor along the signal path can be estimated. This is an integrated value that depends on the density of the water vapor molecules, or alternatively, the associated humidity at each point along the signal path. But from a single integrated value, there is no way to determine the profile of humidity — how the humidity varies with height above the surface. However, if a network of GPS receivers is deployed over a region, it is possible to determine the three-dimensional structure of humidity in the atmosphere above the receivers using tomography in a similar way to that used for medical imaging — albeit with radio waves rather than X-rays.

At the Swiss Federal Institute of Technology in Zurich (familiarly known by its German abbreviation ETH), we have developed the Atmospheric Water Vapor Tomography Software (AWATOS) for estimating humidity profiles. We have tested it with data from various measurement campaigns, including one in Hawaii. We have also used it to determine 40 humidity profiles over Switzerland with data from the Automated GNSS Network for Switzerland (AGNES) of the Swiss Federal Office of Topography, Swisstopo. And recently, we have implemented it in an operational testbed analyzing AGNES data together with observations from the Automated Swiss Weather Station Network (ANETZ) of the Swiss Federal Office of Meteorology and Climatology, MeteoSwiss.

To assess the potential of ground-based GPS water vapor tomography to support meteorological forecasting systems, the tomographic results must be available within near real-time and must be produced with an accuracy comparable to that of existing meteorological measurement techniques and numerical weather prediction models. With those goals in mind, we have carried out a project to determine humidity profiles in a region of the Swiss Alps. In this article, we outline the project, including the background theory, and discuss how we validated the results by comparing them to radiosonde measurements and weather prediction models.

Theoretical Background

Before looking at the project, we will briefly describe the theory behind our tomographic technique.

Radio Wave Refractivity. The propagation of radio waves through the Earth’s ionosphere and the electrically neutral atmosphere (the air) is accompanied by phase and amplitude variations caused by the varying refractive index of the media. Since the effect of the ionosphere on GPS signals can be removed almost completely by processing measurements on both the L1 and L2 frequencies, we are only concerned with the effect of the neutral atmosphere here. In 1951, Essen and Froome published a general formula for the refractive index of air, n, and the corresponding atmospheric refractivity, N, using the three meteorological parameters: total (barometric) air pressure, p, measured in hectopascals; air temperature, T, in kelvins; and the partial pressure of water vapor, e, in hectopascals (see Equation 1). The associated empirically determined constants k1, k2, and k3 have been continuously improved over the years.

In the weighted mean formula for non-dispersive radio wave refractivity for air with 0.0375 percent carbon dioxide content, k1 is set to 77.6890 kelvins per hectopascal, k2 to 71.2952 kelvins per hectopascal, and k3 to 375463 kelvins-squared per hectopascal. The k1 term of Equation 1 can be associated with dry refractivity (N^dry), the refractivity of the dry constituents of air, and the second and third terms with the wet part (N^wet), which is proportional to the partial water vapor pressure.

Tropospheric Refraction. The speed of propagation of a radio wave is governed by the refractivity or index of refraction along the signal (slant) path. The path itself is determined by Snell’s Law relating angle of incidence to angle of refraction at the boundary of two media with differing refractive indices. As mentioned previously, GPS measurements include the additional or excess delay due to the presence of the neutral atmosphere. Since the bulk of the effect occurs in the lower, denser part of the atmosphere — the troposphere — we commonly refer to it as the tropospheric delay. The tropospheric slant path delay, , between station p and satellite r is deﬁned by the following integral along the signal ray path, s:

By integrating the individual components of N in Equation 2 and applying Equation 1, the tropospheric slant path delay can be written as a function of the meteorological parameters p, T, and e.

Tropospheric delay as a function of the observation zenith angle, , (90° minus the elevation angle) is calculated using appropriate mapping functions. The mapping function, , is deﬁned as the ratio of the electrical path length through the troposphere at a particular geometrical zenith angle to the electrical path length in the zenith direction. Typically, separate mapping functions are used for the dry and wet components. Furthermore, the slant wet delay, , for elevation angles down to 3 degrees can be represented as the sum of the isotropic term, ZWD_p (zenith wet delay at station p) with its corresponding mapping function, and a non-isotropic component, :

The Tomographic Voxel Model. Separate slant delays only provide integrated measures of the tropospheric refractivity — a one-dimensional view, if you like. To get the three-dimensional structure of refractivity, we need a different approach. We divide the troposphere into small volume elements or voxels (short for volumetric pixel). With multiple, simultaneous raypaths criss-crossing the model volume, it is possible, in principle, to estimate the refractivity of each voxel and hence get a height profile of refractivity.

The tomographic voxel model is a three-dimensional geometrical structure with ellipsoidal borders. The grid spacing deﬁnes the resulting resolution of the tomographic analysis. In the horizontal plane, the voxel model covers the whole catchment area. For each voxel, an unknown but constant refractivity is introduced. Figure 1 illustrates the principle by means of one single observation.

$FIGURE 1 Principle of GPS tomography. The refractivity in the atmosphere along the raypath of a GPS satellite signal to a ground-based receiver is discretized by a three- dimensional voxel model.$
FIGURE 1 Principle of GPS tomography. The refractivity in the atmosphere along the raypath of a GPS satellite signal to a ground-based receiver
is discretized by a three- dimensional voxel model.

According to Equation 2, the wet part of the slant path delay ( ) for one observation between station p and satellite r can be expressed as a summation over each individual voxel i of the voxel model with a total of k voxels, through which the GPS signal passes:

The refractivity value, N_i, of each voxel is determined by performing a least-squares adjustment. A priori models and inter-voxel constraints can be introduced into the tomographic inversion system. The a priori tomographic model consists of selected voxels, which have externally estimated refractivity values. Inter-voxel constraints provide a spatially smoothing characteristic, as the actual state (or the refractivity) of the atmosphere is smoothly changing from point to point.

Double-Difference GPS Tomography. The software package AWATOS is based on double-difference GPS observations; that is, the difference of measurements made by a pair of receivers between a pair of satellites. Common errors such as those of satellite and receiver clocks difference out. The remaining errors in the observation equation are primarily just those due to atmospheric refraction. The inﬂuence of the ionosphere can be corrected to ﬁrst order by using a linear combination of dual-frequency data as previously mentioned.

Therefore, in double-difference processing, the tropospheric slant path delay, , can be reconstructed by combining four observations (between two stations p and q and two satellites r and s). Similar to Equation 3, the total double-difference path delay, , can be written as a function of the GPS processing output (the zenith path delay, ZPD, and the double-difference phase residual ):

Usually, the dry and wet path delays are treated individually with appropriate models and mapping functions. This separation is carried out within the software package AWATOS for both the path delays and the phase residuals.

Introducing the double-difference slant path delays, , as well as the estimation of the zenith total delays, ZTD, for each station, a priori refractivity values, N₀, and inter-voxel constraints , (with the scalar product condition ), into the tomographic equation system, the ﬁnal form of the inversion equation for the unknown refractivity, N, according to Equation 4 including the design matrix A of the observations is:

To obtain only the wet part of the resulting refractivity ﬁeld (values of refractivity and their gradients, , the individual components of the tomographic observation vector (the left-hand side of Equation 7) have to be correspondingly preprocessed. This is done by introducing additional meteorological observations or numerical data as well as tropospheric mapping functions and models.

Data Description

We recently carried out two measurement campaigns to study the feasibility of our method on a non-permanent densification network in the Swiss Alps. We were interested in investigating such a small-scale high-resolution conﬁguration to see how it can help to determine and model water vapor over a local, mountainous catchment area. We also carried out these campaigns with an eye towards the development of a near real-time analysis procedure with a high update rate of less than one hour and the potential to support short- and medium-range weather forecasts and hydrological hazard assessment.

The Project Area. Two ﬁeld campaigns, each lasting seven days, were carried out in an area of about 50 kilometers by 50 kilometers in the eastern part of the mountainous canton of Valais in the southwest of Switzerland (see Figure 2) in July and October 2005. This region was selected because of its high degree of exposure to hydrological hazards such as flooding of river valleys.

Figure 2. Project area (identified by the rectangle) in the Swiss Alps in the southwest of Switzerland. The elevation of the topography varies from 500 meters to over 4000 meters above mean sea level.

Besides the impact of the fast-changing meteorological situation in the project area, the rough topography is also a challenge for high-precision GPS analysis because of limited ﬁelds of view.

GPS Network and Meteorological Data. Ground-based geodetic GPS stations with dual-frequency receivers were deployed for continuous measurement during the period of the two campaigns. The network was complemented by permanent GPS stations from the national network. The ensemble of all stations used in July 2005 is portrayed in Figure 3.

Figure 3. The 21 GPS stations in the project area in the mountainous canton of Valais (see also Figure 2) used during the measurement campaign in July 2005. The stations’ altitudes vary between 527 meters (SION) and 3119 meters (ZER2).

In October 2005, the non-permanent three-dimensional geodetic Turtmann network was operated with six additional stations in the vicinity of the stations BRAE, SUST, and EMSH (see Figure 3). Furthermore, for this second campaign, data from three stations of the permanent geodynamics/tectonics network in Valais, TECVAL, in the northwestern part of the project area was available.

Several GPS stations were collocated with non-permanent meteorological measurement systems collecting surface temperature, humidity, and air pressure data. Also, rainfall was recorded for validation purposes at ﬁve ANETZ stations within the project area. The temperature, humidity, and air pressure observations were processed with the software package Collocation of Meteorological Data for Interpolation and Estimation of Tropospheric Path Delays (COMEDIE) developed at the Geodesy and Geodynamics Lab, ETH Zurich. COMEDIE provides a four-dimensional modeling of meteorological data in space and time. It is based on the method of least-squares collocation and interpolation, meaning that the model is described by a functional and a stochastic part. The interpolated data was used for the separation of the total delay (and refractivity) into a dry and a wet part and to obtain a priori values, N₀, for the tomographic analysis (see Equation 7).

To compare the results from GPS processing and tomography, independent measurement techniques were used during the measurement campaigns: solar spectrometry, using the Geodetic Mobile Solar Spectrometer (GEMOSS), for integrated path delays as well as weather balloon soundings up to the tropopause for meteorological proﬁle data.

The Numerical Weather Model COSMO-7. MeteoSwiss uses the COSMO-7 model, developed by the Consortium for Small-scale Modelling, for its operational numerical weather forecasts. The model domain is covered by a grid of 383 × 325 points over western and central Europe with a horizontal resolution of 7 kilometers. The model consists of 45 levels vertically distributed between the ﬁltered orography (or mountain topography) and an altitude of 22.5 kilometers.

For comparison and validation, a subset of the reanalyzed COSMO-7 vertical grid point proﬁle data was processed in order to obtain refractivity proﬁles as well as integrated and interpolated time series of zenith path delays using another of our software packages, Collocation and Interpolation of Tropospheric Path Delays (COITROPA).

Results

We processed the GPS data from the two measurement campaigns and have compared the results with those from GEMOSS, COMEDIE, radiosonde data, and COSMO-7.

GPS Data Processing. The GPS processing yields high-quality receiver coordinates, tropospheric parameter estimates (ZPD), and ionosphere-free double-difference residuals to reconstruct the slant path delays (see Equation 5). International GNSS Service (IGS) precise products, including satellite orbits, were used to analyze the data, and for ray tracing in AWATOS.

Bernese GPS Software, version 5.0, was chosen for the processing of the GPS data due to its ﬂexibility, modular design, and state-of-the-art modeling characteristics. The network solution was obtained by using minimally constrained coordinates of selected stations of the IGS reference frame with baseline lengths of up to 1,000 kilometers. The mean repeatabilities for the north, east, and up components of the daily coordinate solutions for all stations within the project area are given in Table 1. Final as well as ultra-rapid orbits and broadcast ephemerides were used to compare the best possible results with those that could be expected in real time.

The larger number of stations during the October campaign has a positive inﬂuence on the mean coordinate repeatabilities in the horizontal plane, whereas the up component remains at the same order of magnitude. Depending on the antenna and receiver types, there was a slightly positive or negative correlation discovered between the trend of the daily coordinates and the ZTD estimates.

Comparison of ZTD Time Series. The time series of zenith total delays (ZTD) from GPS, GEMOSS, integrated ground meteorological data (COMEDIE), and radiosondes coincide well. In particular, the passage of a cold front with heavy rainfall in the middle of the October 2005 campaign is reflected in the two local delay maximums on October 23 (see Figure 4).

Figure 4. Comparison of zenith total delay (ZTD) at station SUST obtained with COMEDIE, GPS, the local radiosondes (RS) and solar spectrometry (GEMOSS) for the October campaign in 2005. The mean values of the ZTD time series and the standard deviation are given for each technique for comparison purposes in parentheses.

The ZTD from the balloon soundings show an almost systematic overestimation. This may be due to an inaccurate self-calibration of the sensors or a lack of data in the upper atmosphere, and the related mismodeling of the zenith path delay. Differences in the COMEDIE time series are due to meteorological inhomogeneities in the lowest part of the troposphere and the inﬂuence of distant radiosondes, which were added to get the vertical information in the upper part. The interpolated ZTD values derived from the numerical weather model COSMO-7 are on average smaller than the GPS estimates (see Table 2).

The ZTD time series of both methods, GPS and the numerical weather model, correlate well with rainfall data. There is a slow increase of the zenith path delay before the precipitation event due to the accumulation of atmospheric water vapor and an abrupt decrease afterwards. Usually, the impact of short periods of localized precipitation is more clearly represented in the GPS data of the dense observation network than in the data of the weather model. The COSMO-7 time series seem to be too smooth.

Effect of Voxel Model Resolution. In order to assess the quality of the results obtained by applying the high-resolution GPS tomographic technique, special time series contour plots were created. They consist mainly of the wet refractivity proﬁles for each voxel model column between mean sea level and an altitude of 10 kilometers. The height of the nearest GPS station is given by a dashed line.

Figure 5 and Figure 6 give two examples from station SUST in the northwest of the project area (see Figure 3) during the October 2005 campaign. Figure 5 shows the wet refractivity variation from a 16-layer voxel model with 5-kilometer horizontal grid spacing, whereas Figure 6 was calculated with 32 layers with the same horizontal resolution.

$Figure 5. Vertical wet refractivity distribution (in parts per million) from a 16-layer voxel model (the increasing vertical distances with height of the voxels are given by black tick marks on the left-hand side) in October 2005. The time series of integrated wet delays (ZWD) from 15 radiosondes (RS), the interpolated proﬁles from the numerical weather model COSMO-7, and the GPS tomographic results (AWATOS) are shown for comparison purposes with a corresponding scale on the right-hand side. Mean values and their standard deviations are shown in parentheses.$
Figure 5. Vertical wet refractivity distribution (in parts per million) from a 16-layer voxel model (the increasing vertical distances with height of the voxels are given by black tick marks on the left-hand side) in October 2005. The time series of integrated wet delays (ZWD) from 15 radiosondes (RS), the interpolated proﬁles from the numerical weather model COSMO-7, and the GPS tomographic results (AWATOS) are shown for comparison purposes with a corresponding scale on the right-hand side. Mean values and their standard deviations are shown in parentheses

$Figure 6. Vertical wet refractivity distribution (in parts per million) from a 32-layer voxel model over the timespan of the October campaign in 2005 and ZWD time series of integrated AWATOS, COSMO-7, and the corresponding radiosonde (RS) proﬁles for comparison purposes.$
Figure 6. Vertical wet refractivity distribution (in parts per million) from a 32-layer voxel model over the timespan of the October campaign in 2005 and ZWD time series of integrated AWATOS, COSMO-7, and the corresponding radiosonde (RS) proﬁles for comparison purposes.

The integrated wet refractivity proﬁles and the reference radiosonde measurements agree better the more layers that are introduced into the tomographic voxel model. The largest differences between the results with different numbers of layers can be detected in the middle troposphere between 4- and 6-kilometers altitude (see Figure 7). It is also recognizable that in the lower troposphere, voxel models with a large number of layers are even able to resolve refractivity inversions.

$Figure 7. Tomographic wet refractivity proﬁles (in parts per million) from 16-, 26-, and 43-layer voxel models, and that of the corresponding radiosonde (RS) launched at station GRUB at 1844 meters above mean sea level on July 13, 2005, 17:04 UTC.$
Figure 7. Tomographic wet refractivity proﬁles (in parts per million) from 16-, 26-, and 43-layer voxel models, and that of the corresponding radiosonde (RS) launched at station GRUB at 1844 meters above mean sea level on July 13, 2005, 17:04 UTC.

We analyzed tomographic voxel models with horizontal resolutions of 15, 10, 7.5, 5, 3.75, and 3 kilometers. Increasing the horizontal resolution of the model leads to an increase in the estimated wet refractivity above an altitude of 6 kilometers compared to both the radiosondes and the numerical weather model. Due to the mean distance of about10 kilometers between the ground-based GPS stations in the project area, the best results were obtained with a similar resolution.

Table 3 gives the results of the comparison between the wet refractivity proﬁles of the tomographic analysis and the radiosondes launched within the project area.

Effect of Temporal Resolution. The tomographic results shown here are based on one-hour time windows for the GPS double-difference data. Higher update rates are also possible without changing the input options of AWATOS. Figure 8 shows the wet refractivity variation based on a 10-minute window together with rainfall data at the ANETZ station Evolène in the western part of the project area.

$Figure 8. Wet refractivity distribution at station Evolène (EVOL) in the western part of the project area from a 26-layer tomographic voxel model with an update interval of 10 minutes. Rainfall data in millimeters per 10 minutes is shown as vertical bars with the corresponding scale on the right-hand side.$
Figure 8. Wet refractivity distribution at station Evolène (EVOL) in the western part of the project area from a 26-layer tomographic voxel model with an update interval of 10 minutes. Rainfall data in millimeters per 10 minutes is shown as vertical bars with the corresponding scale on the right-hand side.

Even though the wet refractivity proﬁles are affected by higher-frequency variations in the upper troposphere, precipitation and weather changes are still recognizable in the 10-minute time series.

Although the Bernese GPS Software is not designed for real-time parameter estimation, near real-time conditions can be simulated by introducing speciﬁc input ﬁles. Thus, the sensitivity of AWATOS to real-time conditions can be assessed. In terms of coordinate repeatability, the results of the horizontal components degrade by about 30 percent when using the predicted part of the ultra-rapid products (see also Table 1). Using broadcast ephemerides, the three-dimensional accuracy suffers even more.

Implications. We collected input data for both the spatially and temporally high-resolution GPS tomographic analysis and the validation of the results. The inhomogeneous distribution of rainfall in the local project area would necessitate even more rain gauges in the meteorological measurement network to perform a hydrological hazard assessment.

The comparison of independent techniques showed that the ZTD time series agree within 2 centimeters on average; that is, to better than 1 percent. The correlation of the GPS data and the data derived from the numerical weather model is greater than 70 percent. However, local rain showers are sometimes more clearly represented by the data of the dense GPS network than by COSMO-7.

It is possible to increase the spatial and temporal resolution in GPS tomography, so it can enhance numerical weather models. The better agreement of the tomographic proﬁles with radiosonde data, compared to the COSMO-7 estimates, indicates that the numerical weather prediction models will beneﬁt from additional information on the vertical distribution of water vapor provided by high-resolution GPS tomography.

To assess the potential of near real-time GPS tomography, IGS satellite products with short latency and fast update rates were tested in the GPS processing. With ultra-rapid orbits, we obtained satisfactory results for the tropospheric parameters in almost real-time mode. The use of predicted orbits in the tomographic processing degrades the results of the wet refractivity proﬁles by 20 percent compared to using ﬁnal (that is, best available) products.

Conclusions

In this brief article, we have shown that high-resolution GPS tomography is well suited for application in mountainous regions, especially in view of its potential to contribute to hydrological hazard assessment. We have been able to estimate the wet refractivity ﬁeld with a spatial
and temporal resolution comparable with the current and the next generation of numerical weather models (COSMO-2 with 2-kilometer horizontal resolution).

We have been successful in illustrating several beneﬁcial aspects of GPS tomography in supporting high-resolution numerical weather prediction models. We would also point out that tomographically determined wet refractivity ﬁelds may also be used in conjunction with directly estimated integrated slant path delays to adjust the GPS observations, especially those obtained at low elevation angles. Implemented in GPS processing software, GPS tomography could provide completely anisotropic tropospheric corrections for very high-accuracy positioning applications.

Acknowledgments

The research discussed in this article was ﬁnancially supported by the Swiss National Science Foundation and the Swiss Geodetic Commission.

The Swiss Federal Office of Meteorology and Climatology, MeteoSwiss, and the Swiss Federal Office of Topography, Swisstopo, provided necessary data sets for processing, analyzing, and validating the results.

Furthermore, O. Heller and Dr. A. Somieski supported the ﬁeld work and several residents or public organizations in the canton of Valais offered their premises for temporary mounting of the campaign measurement systems.

SIMON LUTZ is a research fellow at the Astronomical Institute of the University of Bern, Switzerland, and a member of the Center for Orbit Determination in Europe analysis center team. He received M.S. and Ph.D. degrees in geodesy and geodynamics from the Swiss Federal Institute of Technology (ETH), Zurich, Switzerland.

MARC TROLLER is a communications, navigation, and surveillance (CNS) expert at Swiss Air Navigation Services Ltd., Skyguide, Switzerland. He received M.S. and Ph.D. degrees in geodesy and geodynamics from ETH Zurich.

DONAT PERLER is a Ph.D. candidate at ETH Zurich. He received an M.S. degree in computer science from ETH Zurich.

ALAIN GEIGER is a professor in the Geodesy and Geodynamics Lab of the Institute of Geodesy and Photogrammetry at ETH Zurich. He received an M.S. degree in physics and a Ph.D. in geodesy and geodynamics, both from ETH Zurich.

HANS-GERT KAHLE is professor emeritus of geodesy and geodynamics at ETH Zurich and was leader of the Geodesy and Geodynamics Lab from 1979 to 2009. He received a Ph.D. degree from the University of Kiel, Germany, and the habilitation in geophysics from ETH Zurich.

FURTHER READING

• Seminal Paper on Use of GPS for Meteorology

“GPS Meteorology: Remote Sensing of Atmospheric Water Vapor Using the Global Positioning System” by M. Bevis, S. Businger, T.A. Herring, C. Rocken, R.A. Anthese, and R.H. Ware in Journal of Geophysical Research, Vol. 97, No. D14, 1992, pp. 15787–15801, doi:10.1029/92JD01517.

• Other Studies on Using GPS to Monitor the Atmosphere

“Using the Global Positioning System to Study the Atmosphere of the Earth: Overview and Prospects” by J.L. Davis, M.L. Cosmo, and G. Elgered in GPS Trends in Precise Terrestrial, Airborne, and Spaceborne Applications edited by G. Beutler, G.W. Hein, W.G. Melbourne, and G. Seeber, editors, Volume 115 of the International Association of Geodesy Symposia, Springer Verlag, Berlin, 1996, pp. 233–242.

“GPS Meteorology: Direct Estimation of the Absolute Value of Precipitable Water” by J. Duan, M. Bevis, P., Fang, Y. Bock, S. Chiswell, S. Businger, C. Rocken, F. Solheim, T. van Hove, R. Ware, S. McClusky, T.A. Herring, and R.W. King in Journal of Applied Meteorology, Vol. 35, No. 6, 1996, pp. 830–838, doi:10.1175/1520-0450(1996)035<0830:GMDEOT>2.0.CO;2.

• Effect of the Atmosphere on GPS Positioning

“Atmospheric Modelling in GPS Analysis and Its Effect on the Estimated Geodetic Parameters” by T.R. Emardson and P.O.J. Jarlemark in Journal of Geodesy, Vol. 73, No. 6, 1999, pp. 322–331, doi:10.1007/s001900050249.

• GPS Tomography

High-Resolution GPS Tomography in View of Hydrological Hazard Assessment by S.L. Lutz, Ph.D. dissertation, Eidgenössische Technische Hochschule (ETH) Zürich, Nr. 17675, Zürich, Switzerland, 2008, doi:10.3929/ethz-a-005648120.

“Determination of the Spatial and Temporal Variation of Tropospheric Water Vapour Using CGPS Networks” by M. Troller, A. Geiger, E. Brockmann, and H.-G. Kahle in Geophysical Journal International, Vol. 167, No. 24, 2006, pp. 509–520, doi:10.1111/j.1365-246X.2006.03101.x.

“Diagnosis of Three-Dimensional Water Vapor Using a GPS Network” by A.E. MacDonald, Y. Xie, and R.H. Ware in Monthly Weather Review, Vol. 130, No. 2, 2002, pp. 386–397, doi:10.1175/1520-0493(2002)130<0386:DOTDWV>2.0.CO;2.

“3-D Refractivity Field from GPS Double Difference Tomography” by M. Troller, B. Bürki, M. Cocard, A. Geiger, and H.-G. Kahle in Geophysical Research Letters, Vol. 29, No. 24, 2149, 2002, 4 pp. doi:10.1029/2002GL015982.

• Radio Wave Refractivity of Air

“Refractive Index Formulae for Radio Waves” by J.M. Rüeger in Proceedings of the XXII International Federation of Surveyors (FIG) International Congress, Washington, D.C., April 19–26, 2002.

• Previous GPS World Articles on Tropospheric Propagation Delay

“Tropospheric Delay: Prediction for the WAAS User” by P. Collins and R.B. Langley in GPS World, Vol. 10, No. 7, July 1999, pp. 52–58.

“The Effect of Weather Fronts on GPS Measurements” by T. Gregorius and G. Blewitt, in GPS World, Vol. 9, No. 5, May 1998, pp. 52–60.

“Effect of the Troposphere on GPS Measurements” by F.K. Brunner and W.M. Welsch, in GPS World, Vol. 4, No. 1, January 1993, pp. 42–51.

July 1, 2010
Innovation: GPS, GLONASS and More
Multiple Constellation Processing in the International GNSS Service

By Tim Springer and Rolf Dach

Does combining GPS and GLONASS observations make a difference? The International GNSS Service (IGS) has been providing such data for several years. Representatives from two IGS analysis centers discuss the past, present, and future of IGS GNSS monitoring and product development.

INNOVATION INSIGHTS by Richard Langley

ARE WE THERE YET — at a multiple-constellation GNSS world? The European Galileo system only has two test satellites in orbit, with constellation completion not scheduled until 2014. The Chinese Beidou/Compass system has launched some test satellites, but global coverage is not promised until 2020. And the first Japanese Quasi-Zenith Satellite System space vehicle is scheduled for launch this year with the system not fully operational until 2013. So, does this mean GPS is still the only game in town? No, not by a long shot. We have overlooked Russia’s GLONASS.

Standing for Global’naya Navigatsionnaya Sputnikova Sistema, GLONASS was conceived by the former Soviet Ministry of Defence in the 1970s, perhaps as a response to the announced development of GPS. The first satellite was launched on October 12, 1982. But because of launch failures and the characteristically brief lives of the satellites, a further 70 satellites were launched before a fully populated constellation of 24 functioning satellites was achieved in early 1996. Unfortunately, the full constellation was short-lived. Russia’s economic difficulties following the dismantling of the Soviet Union hurt GLONASS. Funds were not available, and by 2002 the constellation had dropped to as few as seven satellites, with only six available during maintenance operations! But Russia’s fortunes turned around, and with support from the Russian hierarchy, GLONASS was reborn. Longer-lived satellites were launched, as many as six per year, and slowly but surely the constellation has grown to 21, with two in-orbit spares.

But are there any users outside Russia? Although dual-system GPS/GLONASS receivers have been around for at least a decade, manufacturers have taken notice of GLONASS’s recent phoenix-like rebirth. All of the high-end manufacturers now offer receivers with GLONASS capability. Does combining GPS and GLONASS observations make a difference? You bet — just ask any surveyor who uses both systems in the real-time kinematic (RTK) approach. Scientific applications requiring high-accuracy satellite orbit and clock data also benefit. The International GNSS Service (IGS) has been providing such data for several years, and in this month’s article representatives from two IGS analysis centers discuss the past, present, and future of IGS GNSS monitoring and product development.

So, getting back to our question, are we there yet? Many early adopters of GPS plus GLONASS data and products would reply with a resounding “yes.”

“Innovation” features discussions about advances in GPS technology, its applications, and the fundamentals of GPS positioning. The column is coordinated by Richard Langley, Department of Geodesy and Geomatics Engineering, University of New Brunswick.

In 2005, the International GPS Service (IGS) was renamed the International GNSS Service. With this change, the IGS governing board and the IGS community expressed their expectation to extend activities from the well-established GPS to other active and planned global navigation satellite systems such as GLONASS, Galileo, and Compass. Meanwhile, the GLONASS satellite constellation, as well as the IGS GNSS tracking network, have evolved significantly. Since 2003, the GLONASS satellite constellation has been improving steadily, leading to the current, May 2010, constellation with 21 operational satellites and two in-orbit spares. And starting in 2008, the GNSS capabilities of the IGS tracking network have been greatly enhanced giving rise to a truly global GNSS tracking system with more than 100 GNSS (GPS plus GLONASS) receivers. The almost-complete GLONASS satellite constellation, coupled with a readily available global tracking network with high-quality receivers, have greatly increased the interest in and need for GNSS products such as precise satellite orbit ephemerides. However, the IGS analysis center products are still mainly GPS-only. Only two analysis centers provide true multi-GNSS solutions. Two analysis centers provide GLONASS-only solutions (a GLONASS combined IGS product is available but without accurate clocks). No combined IGS GNSS product exists. In view of the large interest from the user community, this is a really disappointing situation. In particular, because experiences gathered with handling GPS plus GLONASS will make the incorporation of other GNSS such as Galileo, Compass, and the Quasi-Zenith Satellite System (QZSS) that much easier.

However, during a meeting of the IGS analysis centers in December 2009, it became clear that many of the centers had started to implement and enhance the GLONASS processing capabilities in their software. This was happening as a direct consequence of the improvements in the GLONASS constellation, the IGS GNSS tracking network, and increased user interest (if not demand). Throughout 2010 and 2011, we will therefore see a significant increase in the number of true GNSS solutions within the IGS. A very positive development for the GNSS world.

In this article, we give an overview of the recent developments in the area of multi-GNSS processing within the IGS in general, but with a focus on the activities of the two analysis centers in the IGS that are leading the GNSS efforts: the Center for Orbit Determination in Europe (CODE) and the European Space Operations Center (ESOC) of the European Space Agency.

Why GNSS?

Within the IGS, we are often confronted with the question: Why GNSS? Why should I go through the burden of adding GNSS capabilities to my software, having larger processing loads, and so on, for little or no benefit? Well, from an IGS analysis center point of view, this question is valid. The accuracies achieved with GPS alone are so good that there will be little visible benefit of including another system. Nevertheless, there are indeed benefits.

There is a large number of users worldwide who would see benefits of using GNSS products compared to GPS-only products. Clearly, all real-time users will benefit enormously from the increased number of satellites. Figure 1, showing the so-called position dilution of precision (PDOP), demonstrates this very clearly. The two panels in Figure 1 show the GPS-only PDOP and the GPS-plus-GLONASS PDOP using the satellite constellation of May 3, 2009.

FIGURE 1A. Effect of GLONASS on position dilution of precision.

FIGURE 1B. Effect of GLONASS on position dilution of precision.

Figure 2 shows the PDOP improvement in percentage when comparing the GPS-only to the GPS-plus-GLONASS PDOP values. At high latitudes, that is, above 55 degrees, the improvement is at the 30 percent level. At mid-latitudes, the improvements are still well above 15 percent, demonstrating the significant improvements real-time GNSS users may expect compared to real-time GPS-only users.

Figure 2. Position dilution of precision improvement using GLONASS.

With the current GPS constellation, daily solutions are not limited by the number of available satellites, but rather by the analysis models (such as that for the troposphere), calibration uncertainties (such as models for antenna phase-center variation), and environmental effects (such as multipath). For these reasons, IGS-like processing strategies, in which data from reference stations are processed in 24-hour batches, will not show clear benefits from adding data from more satellites and other systems.

However, besides real-time users, users at high latitudes (including the whole of Canada and most of Europe) will see improvements. Recently, several researchers have noticed that for latitudes higher than 50 degrees, the addition of GLONASS brings benefit. This is, of course, thanks to the higher orbital inclination of the GLONASS satellites (about 64 degrees) compared to the inclination of the GPS satellites (about 55 degrees), which is also very nicely demonstrated in the PDOP (see Figure 1). So, from a service point of view — the “S” in IGS — there is a clear need to provide GNSS solutions to the user community. Besides offering significant benefits in terms of accuracy, the increased number of satellites will also make solutions more reliable and robust. The completely different repeat cycle of the GLONASS satellite orbits is especially important as it changes the sensitivity to multipath completely. Multipath effects in GPS-only data repeat almost perfectly from day to day with a 4-minute time shift giving rise to spurious, near yearly signals in GPS time series. Satellites from other constellations, such as GLONASS, introduce other system-related frequencies, which results in a general reduction of such GNSS-induced frequencies in a multi-GNSS solution.

Because of the constellation design, each GPS satellite follows its own ground track in each orbit cycle. That means that at a ground station, each GPS satellite is observed on one and the same track each day so that a systematic influence of a satellite (such as a mismodeling of the satellite antenna position with respect to the satellite’s center of mass) has a systematic effect on the obtained (daily) station positions. This systematic translation of satellite-related errors into station-related parameters doesn’t happen for any other GNSS constellation.

IGS GNSS Analysis Centers

A detailed description of the IGS is beyond the scope of this article; an excellent overview was provided in an earlier Innovation column. We simply point out here that it is important to know that the IGS serves as the reference in many GNSS applications by providing data and products of the highest possible quality. Very well known and widely used are the tracking data from the IGS station network — the raw pseudorange and carrier-phase measurements — and the orbit and clock products of the GPS satellites. The IGS generates these products by combining the orbit and clock solutions of the individual analysis centers that contribute to the IGS. For the GPS-only products, 10 different analysis centers contribute to three different product series called the ultra-rapid, rapid, and final products. The final products deliver the highest possible quality but have the longest delay, as they become available 12 days after the end of the observation week. The rapid products are roughly comparable in quality to the IGS final products, but they are delivered daily with a delay of only 17 hours after the end of the observation day. The ultra-rapid products are delivered four times per day 3 hours after the end of the last used observation. For example, at 03:00 UTC, an ultra-rapid product is delivered that used data up to 00:00 UTC. It consists of two parts: an estimated part and a predicted part that may be used for real-time purposes. The quality of the estimated part is very similar to that of the rapid products. The predicted part is, of course, significantly less accurate, although the orbits have an astonishing precision of well below 30 millimeters — much better than that of the orbits in the satellites’ broadcast navigation messages.

In addition to these GPS-only products, there is also a GLONASS product. However, contrary to the GPS side of things, for GLONASS, only a final product is generated. Four analysis centers provide products for the IGS GLONASS combination: the Bundesamt für Kartographie und Geodäsie (BKG), Frankfurt am Main, Germany; CODE, based at the Astronomical Institute of the University of Bern, Switzerland; ESOC, Darmstadt, Germany; and the Information-Analytical Center (IAC) of Roscosmos, Moscow, Russia.

The analysis centers BKG and IAC determine the GLONASS satellite orbits, introducing the information for the GPS satellites from the IGS solution without further estimation. The analysis center CODE provides, since May 2003, orbits for GPS and GLONASS based on a rigorously combined analysis of the data of both GNSS, that is, a true multi-GNSS solution. Since January 2008, ESOC follows this strategy as well. From these four analysis centers, only two, ESOC and IAC, provide satellite clock estimates for the GLONASS satellites. This situation prevents the IGS from making a robust and reliable combined GLONASS clock product. With four analysis centers contributing to the orbits, the IGS can and does make an excellent GLONASS combined orbit product.

In our definition of true multi-GNSS solutions, the measurements from each system contribute to all relevant parameters to the same extent. This can only be achieved by the rigorous combined processing of the data from all available GNSS. The two-step approach, introducing the GPS solution when solving for the GLONASS orbits and satellite clocks, is regarded as an extension of a GPS-only solution to GLONASS. As the contributions from BKG and IAC in the IGS GLONASS product demonstrate, this two-step procedure provides excellent results.

From a user point of view, a big disadvantage is the fact that the IGS does not provide a real GNSS product. The IGS provides a high-quality GPS product and a high-quality GLONASS orbit product, but there is no combined GNSS product. Also, the IGS is only capable of generating final GLONASS products because only two analysis centers, CODE and ESOC, submit GNSS products for the rapid and ultra–rapid products. IGS policy requires contributions from at least three analysis centers for a meaningful and robust combined product.

Users of GNSS orbits and/or clocks therefore have to use the products of one of the individual analysis centers or combine the GPS-only and GLONASS-only products from the IGS. Here, the GNSS products of the CODE and ESOC analysis centers are clearly preferable over those of the IGS and other analysis centers since these are the only two true GNSS products that guarantee the full consistency between the two GNSS.

GLONASS Tracking Network

Until 2003, the IGS had established a GLONASS tracking network of merely 20 stations. In 2003, this number grew rapidly from 20 to 30, but after 2003 the number of stations remained stable for quite a long time with a very inhomogeneous distribution. For example, there were only a few stations in the whole western hemisphere. In 2006/2007, a new generation of combined GPS/GLONASS receivers became available, produced by several well–known GPS receiver manufacturers. With this new equipment available, the number of GLONASS tracking stations in the IGS network started to increase steadily. In 2008, the increase rate went up significantly (see Figure 3) and, more importantly, the global distribution of the receivers improved as, finally, significant numbers of stations started to emerge in both North and South America. Orbits and clocks of the GLONASS satellites are, since ear
ly 2009, determined from the data of more than 100 globally well-distributed tracking stations in the IGS network (see Figure 4). A good global distribution of observing sites is extremely important for orbit determination and even more so for the clock determination. Until early in 2008, the GLONASS clock determination suffered from gaps in the global tracking network, which had severe impact on the clock estimates. If tracking gaps cause an interruption of the carrier-phase tracking of a GNSS satellite, the clock estimates are basically reset and a jump will occur. The size of the jump is delimited by the accuracy of the code (pseudorange) observations, that is, at the 1-meter level, or 3 nanoseconds in clock terms.

We may state that today orbit and clock determination for the GLONASS satellites may be based on a truly global tracking network of high-quality geodetic–type receivers. This significant improvement is due to the efforts of many IGS station managers and their institutions.

Figure 3. Number of sites in the IGS network providing GLONASS data, used for orbit determination at CODE.

Figure 4. Current distribution of IGS combined GPS and GLONASS tracking stations.

GLONASS Constellation

After reaching a full orbit constellation of 24 satellites in early 1996, the GLONASS constellation degraded rapidly due to Russia’s economic difficulties following the break-up of the Soviet Union coupled with the short lifetime of the GLONASS satellites. Since 2002, the GLONASS constellation has slowly but surely been rebuilt (see Figure 5). Currently, there are 21 active modernized GLONASS (GLONASS-M) satellites, which have a significantly longer lifespan compared to the original satellites. Additionally, there are two reserve satellites on orbit.

Figure 5. Development of the GLONASS satellite constellation since 1982.

Russia intends to have a full 24-satellite constellation in place by the end of 2010. To achieve this goal, two more triple-satellite launches are planned, one in August and one in November. The November launch could include a new type of GLONASS satellite, GLONASS-K. The GLONASS-K version is a lighter, unpressurized spacecraft, with a design lifetime of 10 years. In addition to the legacy frequency-division-multiple-access signals, it will transmit code-division-multiple-access signals and use an additional frequency band overlapping with the GPS L5 band.

Orbit and Clock Accuracy

The developments of both the GLONASS tracking capabilities of the IGS station network as well as the steady increase in the number of GLONASS satellites has had a positive influence on the accuracy of the GLONASS orbits and clocks. It also has significantly increased the interest in the GLONASS system. The enhancement of the IGS GNSS tracking network from an almost purely European network to a truly global network between 2008 and now has had a significant impact on the quality of the GLONASS orbits and clocks. To show the effect on the quality of the GLONASS orbit estimates, we look at the difference between two independent consecutive solutions spanning 24 hours from 0 to 24 hours GPS Time. We compare the “midnight point” of both solutions, that is, the solution at the end of one day (or arc) and the beginning of the next day (or arc). This will give us a worst-case estimate for the orbit quality because typically the orbit is less accurate at the boundary of the orbital arc compared to the middle of the orbital arc. We have analyzed these orbit differences for all GPS and GLONASS satellites separately for four half-year time spans using the routine IGS GNSS solutions from ESOC. The differences are computed in three different satellite-orbit-related directions: radial, along-track, and cross-track. The times spans are:
- January to June 2008 (6 months)
- July to December 2008 (6 months)
- January to June 2009 (6 months)
- July to December 2009 (6 months)
The results are shown in Figure 6. For the GPS satellites, we cannot see any improvement over time. The quality of the GPS orbits is excellent at the 25- to 35-millimeter level for all three components.

Figure 6. Evolution of GPS and GLONASS orbit quality from January 2008 to December 2009.

Remember, we are looking at the worst-case differences here. For GLONASS, we can see a significant improvement over the four time spans. Early in 2008, the orbit quality was at the 120-millimeter level (cross-track), which has improved significantly to the 85-millimeter level. It is important to note that no processing changes were made during this time interval, and that the improvements are thanks to the improvements in the station tracking network and the GLONASS satellite constellation.

The clock quality is more difficult to assess, but over the timeframe of 2008 to 2009 we have noticed that the clock estimates of the GLONASS satellites have become complete. In 2008, with the still-far-from-global tracking network, there were many gaps in the tracking of the GLONASS satellites. This means that at some epochs no stations were tracking a GLONASS satellite. Such gaps cause jumps in the satellite clock estimates, because the carrier-phase observations become discontinuous, and these jumps are at the 1-meter (3-nanosecond) level. With the improvements of the IGS GNSS tracking network, the GLONASS tracking is now complete and clocks for all epochs are estimated. A comparison of the clocks of the two analysis centers that provide estimated clocks for the GLONASS satellites shows an agreement at the 80-picosecond level, which is only slightly worse than the agreement between the GPS clocks. Significant biases at the few-hundred-nanosecond level exist only in the GLONASS clocks because of receiver internal frequency-dependent delays. The ESOC GNSS orbit and clock products are, however, perfectly suited for precise point positioning using either GPS, GLONASS or, even better, both GNSS simultaneously. It should be noted that since February 2010, the ESOC IGS clock products are now sampled at 30 rather than 300 seconds, which further enhances their suitability.

Conclusions and Outlook

The IGS has promised to become a GNSS service by changing its name in 2005, more than four years ago. Meanwhile, the GLONASS satellite constellation as well as the IGS GNSS tracking network have matured and are practically complete. For the IGS to become a true GNSS service, a substantial number of the analysis centers should provide GNSS contributions to all IGS products: final, rapid, ultra-rapid, and real-time. These products should come from performing a rigorous combined analysis of the observations of all active GNSS satellites. It is expected that over the next two years, we will see a significant increase in the number of true GNSS solutions within the IGS, a very positive development for the GNSS world.

Within the IGS, the analysis centers CODE and ESOC are leading the GNSS efforts. CODE has provided fully consistent GPS/GLONASS products from a rigorously combined processing approach for all IGS products (final, rapid, and ultra-rapid) since May 2003, or for seven years. Since the beginning of 2008, ESOC has followed this good practice for its final products, and in February 2010 ESOC started to produce rapid and ultra-rapid GNSS products. A unique feature of the ESOC products is that they include the clocks for the GLONASS satellites, even with a sampling rate of 30 seconds for the final products. CODE will add GLONASS clocks to its IGS products very soon, during the fi
rst half of 2010. The GLONASS orbit and clock product quality has become comparable to that of the GPS products within the IGS. However, because GLONASS carrier-phase integer ambiguity resolution is difficult, the GLONASS products are and will remain somewhat less accurate than the GPS products.

The experiences gathered at CODE and ESOC by fully combining the observations from the GPS and GLONASS systems will pave the way for the integration of additional systems and signals within the IGS. Hence, IGS will retain its leading position in providing the reference, in the broadest sense of the word, for all GNSS. In the near future, this means the integration of QZSS and Galileo observations as well as the integration of the new triple-frequency signals from the latest generation of GPS satellites, Block IIF, the first of which was scheduled for launch last month.

The positive GNSS developments within the IGS will require an update of the IGS combination software to enable a true GNSS combination. The CODE and ESOC analysis centers have indicated that they are interested in taking on this important task of rewriting and enhancing the IGS orbit and clock combination software to make the IGS a true GNSS service.

Acknowledgments

CODE is a collaboration among the Astronomical Institute, University of Bern (Bern, Switzerland), the Swiss Federal Office for Topography (Wabern, Switzerland), the Bundesamt für Kartographie und Geodäsie (Frankfurt am Main, Germany), and the Institut für Astronomische und Physikalische Geodäsie of the Technische Universität München (Munich, Germany).

The authors are very grateful to the IGS and its numerous contributors for providing the global GNSS tracking data network.

TIM SPRINGER received his Ph.D. in physics from the Astronomical Institute of the University of Bern (AIUB) in 1999. He has been a key person in the development of the Center for Orbit Determination in Europe (CODE), one of the IGS analysis centers, located at AIUB. Since 2004, he has been working for the Navigation Support Office (OPS-GN) at the European Space Operations Centre (ESOC) of the European Space Agency (ESA) in Darmstadt, Germany. In this group, he has led the development of the new ESOC GNSS software, which is used for most GNSS activities at OPS-GN, including GIOVE-A and -B analyses.

ROLF DACH received his Ph.D. in geodesy at the Institut für Planetare Geodäsie of the University of Technology in Dresden, Germany. Since 1999, he has been working as a scientist at AIUB, where he is head of the GNSS research group. He oversees the development of the Bernese GPS Software, used at CODE for activities in the frame of the AIUB IGS analysis center and elsewhere.

FURTHER READING

• GLONASS Status and History

Russian Space Agency’s Information–Analytical Center website: www.glonass-ianc.rsa.ru.

“Renovated GLONASS: Improved Performances of GNSS Receivers” by A.E. Zinoviev, A.V. Veitsel, and D.A. Dolgin in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 22–25, 2009, pp. 3271–3277.

“Other Satellite Navigation Systems” by S. Feairheller and R. Clark, Chapter 11 in Understanding GPS: Principles and Applications, 2nd edition, edited by E.D. Kaplan and C.J. Hegarty, published by Artech House, Boston, 2006.

“GLONASS Performance, 1995–1997, and GPS-GLONASS Interoperability Issues” by G.L. Cook in Navigation, Vol. 44, No. 3, Fall 1997, pp. 291–300.

“GLONASS Review and Update” by R.B. Langley in GPS World, Vol. 8, No. 7, July 1997, pp. 46–51.

• The International GNSS Service

“The International GNSS Service in a Changing Landscape of Global Navigation Satellite Systems” by J.M. Dow, R.E. Neilan, and C. Rizos in Journal of Geodesy, Vol. 83, No. 3-4, March 2009, pp. 191–198, doi:10.1007/s00190-008-0300-3; erratum: Vol. 83, No. 7, July 2009, p. 689, doi: 10.1007/s00190-009-0315-4.

“GNSS Processing at CODE: Status Report” by R. Dach, E. Brockmann, S. Schaer, G. Beutler, M. Meindl, L. Prange, H. Bock, A. Jäggi, and L. Ostini in Journal of Geodesy, Vol. 83, No. 3-4, March 2009, pp. 353–365, doi:10.1007/s00190-008-0281-2.

“The International GNSS Service: Any Questions?” by A.W. Moore in GPS World, Vol. 18, No. 1, January 2007, pp. 58–64.

IGS Central Bureau website. IGS FAQ, Site Guidelines, data and product access information, and network details are available: http://igscb.jpl.nasa.gov

• Benefits of Multi-GNSS

“The Benefits of Multi-constellation GNSS: Reaching up Even to Single Constellation GNSS Users” by B. Bonet, I. Alcantarilla, D. Flament, C. Rodriguez, and N. Zarraoa in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 22–25, 2009, pp. 1268–1280.

“Assessment of GPS/GLONASS RTK Under Various Operational Conditions” by R.B. Ong, M.G. Petovello, and G. Lachapelle in Proceedings of ION GNSS 2009, the 22nd International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 22–25, 2009, pp. 3297–3308.

“The Future is Now: GPS + GLONASS + SBAS = GNSS” by L. Wanninger in GPS World, Vol. 19, No. 7, July 2008, pp. 42–48.

• GNSS Signal Anomalies

“Anomalous Harmonics in the Spectra of GPS Position Estimates” by J. Ray, Z. Altamimi, X. Collilieux, and T. van Dam in GPS Solutions, Vol. 12, No. 1, January 2008, pp. 55–64, doi:10.1007/s10291-007-0067-7.
June 1, 2010
Innovation: Accuracy versus Precision

A Primer on GPS Truth

By David Rutledge

True to its word origins, accuracy demands careful and thoughtful work. This article provides a close look at the differences between the precision and accuracy of GPS-determined positions, and should alleviate the confusion between the terms — making abuse of the truth perhaps less likely in the business of GPS positioning.

INNOVATION INSIGHTS by Richard Langley

JACQUES-BÉNIGNE BOSSUET, the 17th century French bishop and pulpit orator, once said “Every error is truth abused.” He was referring to man’s foibles, of course, but this statement is much more general and equally well applies to measurements of all kinds. As I am fond of telling the students in my introduction to adjustment calculus course, there is no such thing as a perfect measurement. All measurements contain errors. To extract the most useful amount of information from the measurements, the errors must be properly analyzed.

Errors can be broadly grouped into two major categories: biases, which are systematic and which can be modeled in an equation describing the measurements, thereby removing or significantly reducing their effect; and noise or random error, each value of which cannot be modeled but whose statistical properties can be used to optimize the analysis results.

Take GPS carrier-phase measurements, for example. It is a standard approach to collect measurements at a reference station and a target station and to form the double differences of the measurements between pairs of satellites and the pair of receivers. By so doing, the biases in the modeled measurements that are common to both receivers, such as residual satellite clock error, are canceled or significantly reduced. However, the random error in the measurements due to receiver thermal noise and the quasi-random effect of multipath cannot be differenced away. If we estimate the coordinates of the target receiver at each epoch of the measurements, how far will they be from the true coordinates?

That depends on how well the biases were removed and the effects of random error. By comparing the results from many epochs of data, we might see that the coordinate values agree amongst themselves quite closely; they have high precision. But, due to some remaining bias, they are offset from the true value; their accuracy is low. Two different but complementary measures for assessing the quality of the results.

In this month’s column, we will examine the differences between the precision and accuracy of GPS-determined positions and, armed with a better understanding of these often confused terms, perhaps be less likely to abuse the truth in the business of GPS positioning.

“Innovation” features discussions about advances in GPS technology, its applications, and the fundamentals of GPS positioning. The column is coordinated by Richard Langley, Department of Geodesy and Geomatics Engineering, University of New Brunswick.

For many, Global Positioning System (GPS) measurement errors are a mystery. The standard literature rarely does justice to the complexity of the subject. A basic premise of this article is that despite this, most practical techniques to evaluate differential GPS measurement errors can be learned without great difficulty, and without the use of advanced mathematics. Modern statistics, a basic signal-processing framework, and the careful use of language allow these disruptive errors to be easily measured, categorized, and discussed.

The tools that we use today were developed over the last 350 years as mathematicians struggled to combine measurements and to quantify error, and to generally understand the natural patterns. A distinguished group of scientists carried out this work, including Adrien-Marie Legendre, Abraham de Moivre, and Carl Friedrich Gauss. These luminaries developed potent techniques to answer numerous and difficult questions about measurements.

We use two special terms to describe systems and methods that measure or estimate error. These terms are precision and accuracy. They are terms used to describe the relationship between measurements, and to underlying truth. Unfortunately, these two terms are often used loosely (or worse used interchangeably), in spite of their specific definitions. Adding to the confusion, accuracy is only properly understood when divided into its two natural components: internal accuracy and external accuracy.

GPS measurements are like many other signals in that with enough samples the probability distribution for each of the three components is typically bell-shaped, allowing us to use a particularly powerful error model. This bell-shaped distribution is often called a Gaussian distribution (after Carl Friedrich Gauss, the great German mathematician) or a normal distribution. Once enough GPS signal is accumulated, a normal distribution forms. Then, potent tools like Gauss’s normal curve error model and the associated square-root law can be brought to bear to estimate the measurement error.

An interesting aspect of GPS, however, is that over short periods of time, data are not normally distributed. This is of great importance because many applications are based upon small datasets. This results in a fundamental division in terms of how measurement error is evaluated. For short periods of time, the gain from averaging is difficult to quantify, and it may or may not improve accuracy. For longer periods of time the gain from averaging is significant, a normal distribution forms, and the square-root law is used to estimate the gain. The absence of a Gaussian distribution in these datasets (1 hour or less) is one source of the confusion surrounding measurement error. Another source of confusion is the richly nuanced concept of accuracy. By closely looking at each of these, a clear picture emerges about how to effectively analyze and describe differential GPS measurement error.

The GPS Signal

It is helpful to consider consecutive differential GPS measurements as a signal, and thus from the vantage of signal processing. Here, we use the term measurement to refer to position solutions rather than the raw carrier-phase and pseudorange measurements a receiver makes. Sequential position measurements from a GPS system are discrete signals, the result of quantization, transformation, and other processing of the code and carrier data into more meaningful digital output. In comparison, a continuous signal is usually analog based and assumes a continuous range of values, like a DC voltage. A signal is a way of describing how one value is related to another.

Figure 1 shows a time series consisting of a discrete signal from a typical GPS dataset (height component). These data are based on processing carrier-phase data from a pair of GPS receivers, in double-difference mode, holding the position of one fixed while estimating that of the other. The vertical axis is often called the dependent variable and can be assigned many labels. Here it is labeled GPS height. The horizontal axis is typically called the independent variable, or the domain. This axis could be labeled either time or sample number, depending on how we want this variable to be represented. Here it is labeled sample number. The data in Figure 1 are in the time domain because each GPS measurement was sampled at equal intervals of time (1 second). We’ll refer to a particular data value (height) as x_i.

Figure 1. A 10-minute sample of GPS height data.

Ten minutes of GPS data are displayed in Figure 1. These data are the first 600 measurements from a larger 96-hour dataset that forms the basis of this paper. The mean (or average) is the first number to calculate in any error-assessment work. The mean is indicated by . There is nothing fancy in computing the mean; simply add all of the measurements together and divide by the total sample number, or N. Equation 1 is its mathematical form:

[1]

The mean for these data is 474.2927 meters, and gives us the average value or “center” of the signal. By itself, the mean provides no information on the overall measurement error, so we start our investigation by calculating how far each GPS height determination is located away from the mean, or how the measurements spread or disperse away from the center. In mathematical form, the expression denotes how far the ith sample differs from the mean.

As an example, the first sample deviates by 0.0038 meters (note that we always take the absolute value). The average deviation (or average error) is found by simply summing the deviations of all of the samples and dividing by N. The average deviation quantifies the spreading of the data away from the mean, and is a way of calculating precision. When the average deviation is small, we say the data are precise. For these data, the average deviation is 0.0044 meters.

For most GPS error studies, however, the average deviation is not used. Instead, we use the standard deviation where the averaging is done with power rather than amplitude. Each deviation from the mean, , is squared, , before taking the average. Then the square root is taken to adjust for the initial squaring. Equation 2 is the mathematical form of the standard deviation (SD):

[2]

The standard deviation for the data in Figure 1 is 0.0052 meters.

But note that these data have a changing mean (as indicated by the slowly varying trend). The statistical or random noise remains fairly constant, while the mean varies with time. Signals that change in this manner are called nonstationary. In this 10-minute dataset, the changing mean interferes with the calculation of the standard deviation. The standard deviation of this dataset is inflated to 0.0052 meters by the shifting mean, whereas if we broke the signal into one-minute pieces to compensate, it would be only 0.0026 meters.

To highlight this, Figure 2 is presented as an artificially created (or synthetic) dataset with a stationary mean equal to the first data point in Figure 1, and with the standard deviation set to 0.0026 meters. This figure, with its stable mean and consistent random noise, displays a Gaussian distribution (as we will soon see graphically), and illustrates what our dataset is not.

Figure 2. A 10-minute sample of synthetic data.

Contrasting these two datasets helps us to understand a critical aspect of differential GPS data. Analyzing a one-minute segment of GPS data from Figure 1 would provide a correct estimate of the standard deviation of the higher frequency random component, but would likely provide an incorrect estimate of the mean. This is because of its wandering nature; a priori we do not know which of the 10 one-minute segments is closer to the truth. It is tempting then to think that by calculating the statistics on the full 10 minutes we will conclusively have a better estimate of the mean, but this is not true.

The mean might be moving toward or away from truth over the time period. It is not yet centered over any one value because its distribution is not Gaussian. What’s more, when we calculate the statistics on the full 10 minutes of data, we will distort the standard deviation of the higher frequency random component upwards (from 0.0026 meters to 0.0052 meters).

This situation results in a great deal of confusion with respect to the study of GPS measurement error. When we look at Figures 1 and 2 side by side we see the complication. Figure 2 is a straightforward signal with stationary mean and Gaussian noise. Averaging a consecutive series of data points will improve the accuracy. Figure 1 is composed of a higher frequency random component (shown by the circle), plus a lower frequency non-random component. It is the superimposition of these two that causes the trouble. We cannot reliably calculate the increase in accuracy as we accumulate more data until the non-random component converges to a random process. This results in a very interesting situation; in numerous cases gathering more data can actually move the location parameter (the mean, ) away from truth rather than toward it.

To fully understand the implications of this, consider its effect on estimating accuracy. If the mean is stationary, statistical methods developed by Gauss and others could be used to estimate the measurement error of an average for any set of N samples. For example, the so-called standard error of the average (SE) can be computed by taking the square root of the sample number, multiplying it by the standard deviation, and then dividing by the sample number (a method to provide an estimate of the error for any average that is randomly distributed). Equation 3 is its mathematical form:

[3]

which simplifies to S/√N . This model can only be used if the data have a Gaussian distribution. Clearly this model cannot be used for the data in Figure 1, but can be used for the data in Figure 2. The implications are significant. The data from Figure 1 are not Gaussian because of the nonstationary mean, so we do not know if the gain from 10 minutes of averaging is better or worse than the first measurement. By contrast, the data in Figure 2 are Gaussian, so we know that the average of the series is more accurate than any individual measurement by a factor equal to the square root of the measurements.

By shifting these data into another domain we can see this more clearly. Figure 3 shows the 10 minutes of GPS data from Figure 1 plotted as a histogram or distribution of the number of data values falling within particular ranges of values. We call each range a bin. The histogram shows the frequencies with which given ranges of values occur. Hence it is also known as a frequency distribution. The frequency distribution can be converted to a probability distribution by dividing the bin totals by the total number of data values to give the relative frequency. If the number of observations is increased indefinitely and simultaneously the bin size is made smaller and smaller, the histogram will tend to a smooth continuous curve called a probability distribution or, more technically, a probability density function. A normal probability distribution curve is overlain in Figure 3 for perspective. This curve simultaneously demonstrates what a normal distribution looks like, and serves to graphically display the underlying truth (by showing the correct frequency distribution, mean, and standard deviation). It was generated by calculating the statistics of the 96-hour dataset, then using a random-number generator with adjustable mean and standard deviation (this is an example of internal accuracy, and will be discussed at length in an upcoming section). We can see that our Figure 1 dataset is not Gaussian because it does not have a credible bell shape. By contrast, when we convert the synthetic data from Figure 2 into a frequency distribution, we see the effect of the stationary mean — the data are distributed in a normal fashion because the mean is not wandering.

Figure 3. Frequency distribution of a 10-minute sample of GPS height data.

Recall that all that is needed to use the Gauss model of measurement error is the presence of a random process. Mathematically, the measurement accuracy for the average of the data in Figures 1 and 3 is the overall standard deviation, or 0.0052 meters, because there is no gain per the square-root law. In comparison, the measurement accuracy for the average in Figure 4 is SE = (√ 600•0.0026) / 600 = 0.0001 meters. The standard deviation from the mean is still 0.0026 meters, but the accuracy of the averaged 600 samples is 0.0001 meters. Recall that precision is the spreading away from the mean, whereas accuracy is closeness to truth. When a process is normally distributed, the more data we collect the closer we come to underlying truth. The difference between the two is remarkable. Measurement error can be quickly beaten down when the frequency distribution is normal. This has significant implications for people who collect more than an hour of data, and raises the following question: At what point can we use the standard error model?

Figure 4. Frequency distribution of a 10-minute sample of synthetic data.

Frequency Distribution

In an ideal world, GPS data would display a Gaussian distribution over both short and long time intervals. This is not the case because of the combination of frequencies that we saw earlier (random + non-random). As an aside, this combination is a good example of why power is used rather than amplitude to calculate the deviation from the mean. When two signals combine, the resultant noise is equal to the combined power, and not amplitude.

Interesting things happen as we accumulate more data and continue our analysis of the 96-hour dataset. Earlier we discussed calculating the SD and the mean, and we looked at short intervals of GPS data in the time domain and the frequency-distribution domain. Moving forward, we will continue to look at the data in the frequency-distribution domain because it is far easier to recognize a Gaussian distribution there. The goal is to discover the approximate point at which GPS data behave in a Gaussian fashion as revealed by the appearance of a true bell curve distribution.

Figure 5 shows one minute of GPS data along with the “truth” curve for perspective. This normal curve, as discussed above, was generated using a random number generator with programmable SD and mean variables. The left axis shows the probability distribution for the GPS data, and the right axis shows the probability distribution function for the normal curve. This figure reinforces what we already know: one minute of GPS data are typically not Gaussian (Figure 3 shows the same thing for 10 minutes of data).

Figure 5. Frequency distribution of a 1-minute sample of GPS height data.

Figure 6 shows 1 hour of GPS data. The data in Figure 6 show the beginnings of a clear normal distribution. Note that the mean of the GPS data is still shifted from the mean of the overall dataset. The appearance of a normal distribution at around 1 hour of data indicates that we can begin use of the standard error model, or the Gaussian error model. Recall that this states that the average of the collection of measurements is more accurate that any individual measurement by a factor equal to the square root of the number of measurements, provided the data follow the Gauss model and are normally distributed. For one hour of data, the gain is square root of 1 times the SD divided by N. In effect, no gain. But from this point forward each hour of data provides √N gain. Figure 7 shows 12 hours of data with a gain of √12. By calculating the standard error for the average of 12 hours of data, SE = (√12•0.0069)/12, or 0.0020 meters, we see a clear gain in accuracy. Notice also that at 12 hours the normal curve and the GPS data are close to being one and the same.

Figure 6. Frequency distribution of a 1-hour sample of GPS height data.

Figure 7. Frequency distribution of a 12-hour sample of GPS height data.

Several things are worth pointing out here. The non-stationary mean converts to a Gaussian process after approximately 1 hour. There is nothing magical about this; conversion at some point is a necessary condition for the system to successfully operate. If it did not, the continually wandering mean would render it of little use as a commercial positioning system. Because it is non-stationary over the shorter occupations considered normal for many applications, it is confusing. Collecting more data in some instances can contribute to less accuracy. This situation also creates a gulf between those who collect an hour or two, and those who collect continuously. It is worth emphasizing that the distribution of data under our “truth” curve fills out nicely after 12 hours. This coincides with one pass of the GPS constellation, suggesting (as we already know) that a significant fraction of the wandering mean is affected by the geometrical error between the observer and the space vehicles overhead.

By looking at the 12 one-hour Gaussian distributions that comprise a 12-hour dataset, we see clearly what Francis Galton discovered in the 1800s. A normal mixture of normal distributions is itself normal, as Figure 8 shows. This sounds simple, but in fact it has significant implications. The unity between consecutive 1-hour segments of our dataset is the normal outline, reinforcing the increasing accuracy of the location parameter, , as more and more normal curves are summed together.

Figure 8. (a) Frequency distribution of 12 1-hour samples of GPS height data; (b) the 12 1-hour samples combined.

Internal vs. External Accuracy

Figure 9 shows the relationship between precision and accuracy. The dashed vertical line indicates the mean of the dataset (the inflection point at which the histogram balances). The red arrows bracket the spread of the dataset at 1 standard deviation from the mean (precision), while the black arrows bracket the offset of the mean from truth (accuracy). Notice that the mean ( ) is a location parameter, while the standard deviation (<e
m>s) is a spread parameter. What we do with the mean is accuracy related; what we do with the standard deviation is precision related.

Figure 9. Relationship between precision and accuracy.

Accuracy is the difference between the true value and our best estimate of it. While the definition may be clear, the practice is not. Earlier we discussed two techniques used to calculate precision — the average deviation, and the standard deviation. We also discussed the square-root law that estimates the measurement error of a series of random measurements. As we saw, it was not possible to calculate this until roughly 1 hour of data had been collected. Furthermore, the data were said to be accurate when a good correlation appeared between the overlain curve and the GPS data at 12 hours.

But here is the interesting thing; the truth curve was derived internally. As previously discussed, data were accumulated for 96 hours, and then statistics were calculated on the overall dataset. Then a random number generator with programmable mean and standard deviation was used to generate a perfectly random distribution curve with the same location parameter and spread. This was declared as truth, and then smaller subsets of the same dataset were essentially compared with a perfect version of itself! This is an example of what is called internal accuracy.

By contrast, external accuracy is when a standard, another instrument, or some other reference system is brought to bear to gauge accuracy. A simple example is when a physical standard is used to confirm a length measurement. For instance, a laser measurement of 1 meter might be checked or calibrated against a 1-meter platinum iridium bar that is accepted as a standard. The important point here is that truth does not just appear — it has to be established through an internal or external process.

Accuracy can be evaluated in two ways: by using information internal to the data, and by using information external to the data. The historical development of measurement error is mostly about internal accuracy. Suppose that a set of astronomical measurements is subjected to mathematical analysis, without explicit reference to underlying truth. This is internal accuracy, and was famously expressed by Isaac Newton in Book Three of his Principia: “For all of this it is plain that these observations agree with theory, so far as they agree with one another.”

Internal accuracy constrains and simplifies the problem. It eliminates the need to bring other instruments or systems to bear. It makes the problem manageable by allowing us to use what we already have. Most importantly, it eliminates the need to consider point of view. Because we are not venturing outside of the dataset, it becomes the reference frame. By contrast, when you ponder bringing an external source of accuracy to bear it gets complicated, especially with GPS.

For example, is it sufficient to use one GPS receiver to check the accuracy of another, or should an entirely different instrument be used? Is it suitable to use the Earth-centered, Earth-fixed GPS frame to check itself, or should another frame be used? If we use another frame, should it extend beyond the Earth, or is it sufficient to consider accuracy from an Earth perspective? When we say a GPS measurement is accurate, what we are really saying is that it is accurate with respect to our reference frame. But what if you were an observer located on the Sun? An Earth-centric frame no longer makes sense when the point that you wish to measure is located on a planet that is rotating in an orbit around you. For an observer on the Sun, a Sun-centered, Sun-fixed reference frame would probably make more sense, and would result in easier to understand measurements. But we are not on the Sun, so a reference frame that rotates with the Earth — making fixed points appear static — makes the most sense. The difference between the two is that of perspective, and it can color our perception of accuracy.

Internal accuracy assessments sidestep these complications, but make it difficult to detect systematic errors or biases. Keep in mind that any given GPS measurement can be represented by the following equation: measurement = exact value + bias + random error. The random-error component presents roughly the same problem for both internal and external assessments. The bias however, requires external truth for detection. There is no easy way to detect a constant shift from truth in a dataset by studying only the shifted dataset.

In practice, people generally look for internal consistency, as Newton did. We look for consistency within a continuous dataset, or we collect multiple datasets at different times and then look for consistency between datasets. It is not uncommon to use the method taken in this article: let data accumulate until one is confident that the mean has revealed truth, and then use this for all further analysis. For this approach, accuracy implies how the measurements mathematically “agree with one another.”

All of this shows that accuracy is a very malleable term. Internal accuracy assumes that the process is centered over truth. It is implicitly understood that more measurements will increase the accuracy once the distribution is normal. The standard error is calculated by taking the square root of the sample number, multiplying it by the standard deviation, and then dividing by the sample number. With more samples, the standard error of the average decreases, and we say that the accuracy is increasing. Internal accuracy is a function of the standard deviation and the frequency distribution.

External accuracy derives truth from a source outside the dataset. Accuracy is the offset between this truth and the measurement, and not a function of the standard deviation of the dataset. The concept is simple, but in practice establishing an external standard for GPS can be quite challenging. For counterpoint, consider the convenient relationship between a carpenter and a tape measure. He is in the privileged position of carrying a replica of the truth standard. GPS users have no such tool. It is impossible to bring a surrogate of the GPS system to bear to check a measurement. Fortunately, new global navigation satellite systems are coming on line to help, but a formal analysis of how to externally check GPS accuracy leads one into a morass of difficult questions.

Accuracy is not a fundamental characteristic of a dataset like precision. This is why accuracy lacks a formal mathematical symbol. One needs to look no further than internal accuracy for the proof. For a dataset that is shifted away from truth, or biased, no amount of averaging will improve its accuracy. Because it is possible to be unaware of a bias using internal accuracy assessments, it follows that accuracy cannot be inherent to a dataset.

Looking at the interplay between mathematical notation and language provides more insight. For example, we describe the mathematical symbol with the word mean. We don’t stop there, however; we also sometimes call it the average. Likewise, the mathematical symbol s is described by the words standard deviation, but we also know s as precision, sigma, repeatability, and sometimes spread. English has a wealth of synonyms, giving it an ability to describe that is unparalleled. In fact, it is one of only a few languages that require a thesaurus. This is why it is important to make a clear distinction between the relatively clear world of mathematical notation and the more free-form world of words. Language gives us flexibility and power, but can also confound with its ability to provide subtle differences in meaning.

When we look at the etymology of the word accuracy, we can see that it is aptly named. It comes from the Latin word accuro, which means to take care of, to prepare with care, to trouble about, and to do painstakingly. Accuro is itself derived from the root cura, which means roughly the same thing and is familiar to us today in the form of the word curator. It is fitting language for a process that requires so much care.

When we discuss measurement error we seldom use mathematical symbols; we use language that is every bit as important as the symbols. The word error itself derives from the Latin erro, which means to wander, or to stray, and suitably describes the random tendency of measurements.

Whether we describe it with mathematics or language, error describes a fundamental pattern we see in nature; independent measurements tend to randomly wander around a mean. When the frequency distribution is normal, accuracy from the underlying truth occurs in multiples of √N. Error is the umbrella covering the other terms because it is the natural starting point for any discussion. Because of this, precision and accuracy are naturally subsumed under error, with accuracy further split into internal and external accuracy. By contemplating all of this, we expose the healthy tension between words and mathematical notation. Neither is perfect. Mathematics establishes natural patterns and provides excellent approximation tools, but is not readily available to everyone. Language opens the door to perspective and point of view, and invites questions in a way that mathematical notation does not.

Final Notes

Making sense of GPS error requires that we take a close look at the intricacies of the GPS signal, with particular attention to the ramp up to a normal distribution. It also requires a good hard look at the language of error. Shifting the GPS data back and forth between the frequency-distribution and time domains nicely illustrates the complications imposed by a non-stationary mean. Datasets that are an hour or less in duration do not always increase in accuracy when the measurements are averaged. Averaging may provide a gain, but it is not a certainty. When the non-stationary mean converges to a Gaussian process after an hour or so, we begin to see what De Moivre discovered almost 275 hundred years ago: accuracy increases as the square root of the sample size.

The GPS system is so good that the division of accuracy into its proper internal and external accuracy components is shimmering beneath the surface for most users. It is rare that a set of GPS measurements has a persistent bias, so internal accuracy assessments are usually appropriate. This should not stop us from being careful with how we discuss accuracy, however. Some attempt should be made to distinguish between the two types, and neither should be used interchangeably with precision. What’s more, while accuracy is not something intrinsic to a dataset like precision, it is still much more than just a descriptive word. Accuracy is the hinge between the formal world of mathematics and point of view. Its derivation from N and s in internal assessments stands in stark contrast to the more perspective-driven derivation often found in external assessments. When carrying out internal assessments, we must be aware that we are assuming that the measurements are centered over truth. When carrying out external assessments, we must be mindful of what outside mechanism we are using to provide truth. True to its word origins, accuracy demands careful and thoughtful work.

David Rutledge is the director for infrastructure monitoring at Leica Geosystems in the Americas. He has been involved in the GPS industry since 1995, and has overseen numerous high-accuracy GPS projects around the world.

FURTHER READING

• Highly Readable Texts on Basic Statistics and Probability
The Drunkard’s Walk: How Randomness Rules Our Lives by L. Mlodinow, Pantheon Books, New York, 2008.

Noise by B. Kosko,Viking Penguin, New York, 2006.

• Basic Texts on Statistics and Probability Theory
A Practical Guide to Data Analysis for Physical Science Students by Louis Lyons, Cambridge University Press, Cambridge, U.K., 1991.

Principles of Statistics by M.G. Bulmer, Dover Publications, Inc., New York, 1979.

• Relevant GPS World Articles
“Stochastic Models for GPS Positioning: An Empirical Approach” by R.F. Leandro and M.C. Santos in GPS World, Vol. 18, No. 2, February 2007, pp. 50–56.

“GNSS Accuracy: Lies, Damn Lies, and Statistics” by F. van Diggelen in GPS World, Vol. 18, No. 1, January 2007, pp. 26–32.

“Dam Stability: Assessing the Performance of a GPS Monitoring System” by D.R. Rutledge, S.Z. Meyerholtz, N.E. Brown, and C.S. Baldwin in GPS World, Vol. 17, No. 10, October 2006, pp. 26–33.

“Standard Positioning Service: Handheld GPS Receiver Accuracy” by C. Tiberius in GPS World, Vol. 14, No. 2, February 2003, pp. 44–51.

“The Stochastics of GPS Observables” by C. Tiberius, N. Jonkman, and F. Kenselaar in GPS World, Vol. 10, No. 2, February 1999, pp. 49–54.

“The GPS Observables” by R.B. Langley in GPS World, Vol. 4, No 4, April 1993, pp. 52–59.

“The Mathematics of GPS” by R.B. Langley in GPS World, Vol. 2, No. 7, July/August 1991, pp. 45–50.

May 1, 2010
Innovation: GPS by the Numbers

A Sideways Look at How the Global Positioning System Works

In his 200th Innovation column, Contributing Editor Richard Langley takes a look at GPS by the numbers, getting a sense of how GPS works by examining the key numbers that govern its remarkable capabilities, from zero to pi and beyond.

INNOVATION INSIGHTS by Richard Langley

WELCOME TO INNOVATION COLUMN NUMBER 200. I have managed this column continuously since the first issue of GPS World magazine, which appeared back in 1990. From the outset, we established that the column should deal with issues that have broad application and interest and are presented in terms that are accessible to as wide a range of readers as possible. Since 1990, we have covered a wide range of topics, some of them at the leading edge of GPS development and some of them reviewing the basics of GPS operation in tutorial fashion. The column has appeared 199 times and now we come to number 200.

So clearly 200 is an important number for me and, I hope, for you. But the number 200 is interesting for other reasons, too. It is the smallest base 10 unprimeable number — you can’t turn it into a prime number by changing just one of its digits to any other digit. It’s how many dollars you get when you pass Go in Monopoly. And in 2012, it will be how many years have elapsed since The War of 1812 — the last time Canada and the United States had a serious quarrel (other than in hockey). But more to the point of this column, it is the designation of the basic reference document that describes how GPS works: IS-GPS-200. Formerly known as an Interface Control Document or ICD, it has gone through several revisions since its first public release in July 1991. It is full of numbers. Numbers that tell us how the GPS signals are generated and how a receiver is to interpret the signals to provide a position fix.

If you are a regular reader of the Innovation column, then likely you have an inquisitive bent. You like to know how things work — GPS in particular. And you don’t have to be convinced about the importance of numbers and their role in understanding the world around us. As Sir William Thomson, a.k.a. Lord Kelvin, said in one of his lectures,

“I often say that when you can measure what you are speaking about, and express it in numbers, you know something about it; but when you cannot measure it, when you cannot express it in numbers, your knowledge is of a meagre and unsatisfactory kind.”

So in this column, the 200th, we’re going to look at GPS by the numbers, getting a sense of how GPS works by examining some of the key numbers that govern its remarkable capabilities, from the smallest to the largest. I’ll draw heavily on material from the past 199 columns.

Let’s get started.

“Innovation” features discussions about advances in GPS technology, its applications, and the fundamentals of GPS positioning. The column is coordinated by Richard Langley, Department of Geodesy and Geomatics Engineering, University of New Brunswick.

Numbers. We use them for counting and measuring, for labeling and ordering, and for codes and calculations. The number of numbers is infinite. However, there are some special numbers that characterize how GPS works. Some of these are peculiar to GPS; others are more common, finding utility in other global navigation satellite systems or even in our everyday lives. In this article, we’ll take a look at some of these special numbers and their importance to GPS.

We’ll begin with the smallest non-negative number and work our way up to one of the largest GPS-relevant numbers, concluding with an imaginary but very important number.

0

Zero. The smallest cardinal number and the smallest non-negative integer. While zero is a pure real number (a number on an infinitely long number line), it is also a purely imaginary number (see the last entry in this article) because it lies on both the real and imaginary axes on the complex plane. It is used to indicate a null amount. The English mathematician, Alfred North Whitehead, wrote in his 1911 book An Introduction to Mathematics, “The point about zero is that we do not need to use it in the operations of daily life. No one goes out to buy zero fish. It is in a way the most civilized of all the cardinals, and its use is only forced on us by the needs of cultivated modes of thought.” Perhaps it was not needed for daily operations in 1911 but it is indispensible in our modern world. For zero is also one of the two binary digits (the other is one, of course) used in the binary or base-2 number system that is fundamental to how computers, digital electronics, and communications systems operate. For example, we represent the GPS pseudorandom noise (PRN) ranging codes and the navigation message as sequences of zeros and ones and the zeros are just as important as the ones.

The C/A- and P(Y)-codes (see entries 1023 and 235,469,592,765,000), along with the navigation message, are modulated onto the signal carriers using binary phase-shift keying or BPSK. BPSK is a digital modulation scheme that conveys a signal by changing, or modulating, the phase of the carrier wave between two values separated by 180°. The spectrum of a BPSK-modulated signal is a sinc function, with most of the power concentrated around the carrier frequency. An alternative modulation technique is binary offset carrier (BOC) modulation. BOC modulation uses a square-wave subcarrier to offset the spectral power from the carrier frequency and thus allows a BOC-modulated signal to share the same bandwidth as a BPSK signal. The new GPS M-code on L1 and L2 uses a BOC(10.23,5.115) — abbreviated as BOC(10,5) — modulation, which specifies a subcarrier frequency of 1021.023 MHz and a spreading-code chipping rate of 5.115 megachips per second. The spreading code is a pseudorandom bit stream from a signal protection algorithm, having no apparent structure or period. The future L1C signal, the new civil signal to be implemented on L1 by Japan’s Quasi-Zenith Satellite System and the GPS III satellites, will also use BOC modulation. And Europe’s Galileo system, now in development, will also use this modulation technique, which has already been tested in space by the forerunner GIOVE test satellites.

0.00000000000001

(Or 1 x 10^-14 in scientific notation). The approximate frequency stability of the rubidium atomic frequency standards in the GPS Block IIR satellites. These devices are used to control the frequency and timing of all aspects of the navigation signals, including the generation of the carrier frequencies and the pseudorandom noise modulation codes. Given their role in controlling the timing of the signals, they are also referred to as clocks.

Each Block IIR satellite contains three rubidium clocks, only one of which is active at any time. The others are spares and the GPS Control Segment carries out a “clock swap” when the performance of an active clock begins to deteriorate, cycling through the remaining units. Many of the Block IIR satellites are still on their first clock.

The Block IIF satellites will also contain three clocks, however, only two will be rubidium clocks. The third clock will be a cesium clock. This mixture of clock types is patterned after the arrangement used on the Bloc
k II and IIA satellites, which used two rubidium clocks and two cesiums.

0.77922077922…

The rational number 60/77. A rational number is any number that can be expressed as a fraction or quotient a/b of two integers, with the denominator b not equal to zero. The decimal expansion of a rational number always either terminates after a finite number of digits or begins to repeat the same sequence of digits over and over again. The digits 077922 of this particular rational number repeat ad infinitum. And why should we be interested in this particular number? It is the ratio of the L2 and L1 carrier frequencies. This number and its inverse ( — the dot above the 3 indicates it repeats indefinitely) are used in various combinations of GPS measurements. For example, if we let η = 60/77, then the ionosphere-free pseudorange combination is

where P₁ is a pseudorange measurement on L1 and P₂ is the corresponding pseudorange measurement on L2.

1

The loneliest number according to the American rock band Three Dog Night and the Google calculator (try typing “loneliest number” into the Google search engine). It was also the space vehicle number (SVN) of the first Block I GPS satellite, which was launched on February 22, 1978. The satellite did not stay lonely for long. By the end of the year, three more Block I satellites were launched. In total, 10 Block I satellites were successfully orbited between 1978 and 1985 to demonstrate the feasibility of GPS. SVN1 continued in operation until July 17, 1985.

The first satellite of the Block II operational constellation was launched in February 1989. The four-year hiatus in launches was due, in part, to the Space Shuttle Challenger disaster as it had been planned to launch the operational satellites using the Shuttle. Following the accident, it was decided to continue with expendable rockets for GPS launches but to switch to the newly designed Delta II rocket.

The pace of Block II launches was rapid, with five launches of the original Block II design in 1989 and four in 1990. A modified version of the Block II satellite — the IIA — was developed, and between 1990 and 1997 19 Block IIAs were launched. The Block II and IIA satellites established the operational GPS constellation. Full operational capability was declared on April 27, 1995.

A new variant of the Block II satellite was developed for replenishing the constellation as the earlier satellites were retired. Following an initial launch failure, 12 of the Block IIR satellites were launched between 1997 and 2004.

Under the GPS modernization program, the remaining eight Block IIR satellites were modernized with a new navigation payload that included the L2C and M-code signals as well as a new antenna panel (also included on the last four of the classic Block IIRs). The IIR-M satellites were launched between 2005 and 2009, bringing the total number of GPS satellites ever placed in orbit to 58.

One is also the PRN number of SVN49, the Block IIR-M satellite that was modified to transmit the first L5 GPS signals (see 1176.45).

2.4

The approximate delay, in meters, experienced by a GPS signal propagating vertically (from the zenith) through the neutral atmosphere to a receiver at mean sea level. Although the electrically neutral, or unionized, atmosphere extends from ground level up to 50 kilometers and more, the bulk of it is in the lowest most part we call the troposphere. Consequently, the neutral atmosphere delay is often termed the tropospheric delay. The delay varies with actual atmospheric conditions and the elevation angle at which a GPS signal arrives at the receiver’s antenna. If unaccounted for, tropospheric delay would result in position errors of several meters in the horizontal plane and two to three times these values in the vertical. Predictive or “blind” tropospheric models based on climatology attempt to significantly reduce the effect of the troposphere on GPS position fixes.

One such model is UNB3m, developed at the University of New Brunswick. Using a look-up table of surface meteorological parameter values from standard atmospheric models, it can compute the tropospheric delay for a given day of year, latitude, and station height. For example, the UNB3m zenith delay for a sea-level site at a latitude of 60° on day-of-year 201 is 2.435 meters. UNB3m is able to predict zenith delays with an average root-mean-square error of 4.9 centimeters. Better and more consistent performance has been obtained with a wide-area model developed specifically for North America, UNBw.na.

A version of an earlier UNB model became the basis of the RTCA Minimum Operational Performance Standards (MOPS) troposphere model that is included in the firmware of most GPS receivers.

3.1415926…. π

Every nerd’s favorite number. It is the ratio of a circle’s circumference to its diameter in conventional or Euclidean space. We use it, for example, to convert angles measured in radians to degrees (π radians 4 180 degrees). π is an irrational number, which means that its value cannot be expressed exactly as a fraction m/n, where m and n are integers. Consequently, its decimal representation never ends or repeats. But we sometimes use an easily remembered fraction, such as 22/7, to get an approximate value. In this case, 3.14. But, if we compute more digits with this fraction, we get 3.1428571…, clearly an incorrect result. A better way to remember π to eight digits is to count the number of letters in each word of the mnemonic “May I have a large container of coffee?”

In computations related to GPS, how many digits of π should be used? It depends. If you are developing your own algorithms and software for modeling GPS observations or determining precise orbits for the satellites, you’ll likely need π to 16 digits for double-precision floating-point calculations. But it would be a mistake to use π to this precision in computing the position of a satellite from the broadcast ephemeris. The GPS interface specification document, IS-GPS-200, specifies a 14-digit value for π (3.1415926535898) in the satellite coordinate computation. Use fewer or more digits, and the resulting satellite coordinates will not be as accurate.

4

This is the minimum number of satellites that a receiver needs to track and generate a pseudorange measurement to produce a three-dimensional “instantaneous” position fix. The receiver solves a system of four nonlinear equations to obtain the three receiver coordinates and the offset of the receiver’s clock from GPS (System) Time. It is possible to use fewer than four satellites for positioning or navigation, but then additional information must come from elsewhere. For example, if we are navigating and we know our height accurately or can safely assume a value, say, so many meters above the sea surface, then only three pseudoranges would be needed to determine the horizontal coordinates. If the number of satellites drops to two, then another assumption must be made to continue navigation (for example, holding the receiver clock offset constant or assuming a constant driving direction). If the clock offset is held constant, then position accuracy deteriorates quickly since the actual receiver clock offset will diverge from the assumed value. On the other hand, if the direction of travel is held constant, the GPS receiver can at least compute the position along the assumed trajectory. In reality, the vehicle will likely not travel along a perfectly straight path and navigation fails after the first turn.

For single-satellite navigation, three assumptions must be made concerning height, trajectory, and clock offset, but the navigation results are, at best, educated guesses.

If a receiver can acquire and track more than four satellites, it typically uses an optimizing Kalman filter procedure to obtain its position.

10.23

The frequency of a GPS satellite atomic frequency standard in megahertz and the fundamental frequency for signal generation in the satellite. The carrier frequencies and the code chipping rates are harmonically related to this frequency. The P-code chipping rate is identical to the fundamental frequency, while the C/A-code rate is 1/10 of it.

12.5

The length of the full navigation message in minutes. To convert the measured signal delays or pseudoranges between the receiver and the satellites, the receiver must know where the satellites are. To do this in real time requires that the satellites broadcast this information. Accordingly, there is a message superimposed on both the L1 and L2 carriers along with the PRN codes. Each satellite broadcasts its own message, which consists of orbital information (the ephemeris) to be used in the position computation, the offset of its clock from GPS Time, and information on the health of the satellite and the expected accuracy of the range measurements.

The message also contains almanac data for all of the satellites in the GPS constellation, as well as their health status and other information. The almanac data, a crude description of the satellite orbits, is used by the receiver to determine the location of each satellite. The receiver uses this information to quickly acquire the signals from satellites that are above the horizon but are not yet being tracked. So, once one satellite is tracked and its message decoded, acquisition of the signal from other satellites is quite rapid. A receiver will store a copy of the almanac to speed up initial acquisition of satellites when it is switched on.

The GPS navigation message is sent at a relatively slow rate of 50 bits per second, taking 12.5 minutes for all of the information to be transmitted. To minimize the time it takes for a receiver to obtain an initial position, the ephemeris and satellite clock offset data is repeated every 30 seconds.

24

The number of satellites in the current GPS baseline constellation. The GPS constellation went through a number of design alternatives even after the first satellites were launched with different numbers of orbit planes, satellites per plane, and orbit inclinations. The current design has four satellites, irregularly spaced, in six orbit planes. The orbit planes, labeled A through F, are spaced at 60° intervals around the equator with a nominal inclination to the equator of 55°. However, we have typically had a surfeit of satellites with more than 24 in operation since the mid-1990s. In fact, during 2008, as many as 31 satellites were transmitting healthy signals at the same time. However, although a modern GPS receiver should be able to handle a 32-satellite active constellation, there are limits imposed by the GPS Control Segment and some legacy military equipment that currently imposes a 30-satellite active constellation limit.

Although the number of active satellites is well in excess of 24, the constellation has been operated as a 24-satellite constellation without optimizing the orbit locations of the “bonus” satellites. In fact, several pairs of satellites are bunched together minimizing geometrical performance. This is in the process of being changed. The GPS Wing recently announced the transition to a 24+3 or “Expandable 24” baseline constellation. Taking about 24 months to complete, six on-orbit satellites are being rephased within their respective orbit planes to improve the overall geometry of the active constellation so that the number of GPS satellites in view from anywhere on Earth will increase, enhancing the possibility of getting a position fix in partially obscured environments, and potentially improving the accuracy of fixes.

Of course, 24 is also the number of hours in the day during which GPS is available at any point on the Earth’s surface with good sky visibility. It is also the title of a popular American TV series, whose protagonist, Jack Bauer, frequently makes use of imaginary GPS tracking capabilities.

40.3

The scaling factor, which together with the signal frequency and the total electron content, is used to compute the delay experienced by a GPS signal as it propagates through the ionized part of the Earth’s atmosphere. The total electron content is the integrated electron density along the signal’s path. Basically, it is the total number of electrons in a tube with a cross-sectional area of one square meter centered on the signal path. To a very good approximation, the delay, in meters, is computed as

where TEC is the total electron content in so-called TEC units or TECUs (10¹⁶ electrons per square meter) and f is the signal carrier frequency in MHz. The scaling factor is a function of the electron’s charge and mass and a constant of electromagnetism theory called the permittivity of free space also known as the electric constant. The scaling factor is actually 40.308193 but this much precision is not generally needed in GPS calculations.hile the code signals are delayed, making pseudoranges longer than they would be in the absence of the ionosphere, the phases of the signal carriers are advanced, make carrier-phase measurements shorter — but by exactly the same magnitude as the code delays.

TEC is highly variable both temporally and spatially. The dominant variability is diurnal following the variation in incident solar radiation. Maximum ionization occurs at approximately 1400 local time. On the ionosphere’s nighttime side, in the absence of solar radiation, free electrons and ions tend to recombine, thereby reducing the TEC. The protonosphere, or uppermost region of the ionosphere, may contribute up to 50 percent of the electron content during the nighttime hours. Typical nighttime values of vertical TEC for mid-latitude sites are of the order of 10 TECU or less with corresponding daytime values of the order of 100 TECU. However, such typical daytime values can be exceeded by a factor of two or more, especially in near-equatorial regions. TEC also varies seasonally with higher values during equinoxes.

1023

This is the number of chips in the C/A-code. The C/A-, or coarse/acquisition-, code is one of the two legacy PRN ranging codes that have been transmitted by all GPS satellites. These PRN codes consist of sequences of binary values (zeros and ones) that, at first sight, appear to have been randomly chosen. But a truly random sequence can only arise from unpredictable causes over which, of course, we would have no control, and which we could not duplicate. However, using a mathematical algorithm or special hardware devices called tapped feedback registers, we can generate sequences that do not repeat until after some chosen interval of time. Such sequences are termed pseudorandom. The apparent randomness of these sequences makes them indistinguishable from certain kinds of noise such as the hiss heard when a conventional AM radio is tuned between stations.

The C/A-code is a sequence of 1,023 binary digits, or chips, which is repeated every millisecond. This means that the chips are generated at a rate of 1,023 million per second and that one chip has a duration of about 1 microsecond.

The C/A-code is generated by two 10-cell feedback registers referred to as G1 and G2. A delayed version of the G2 sequence is obtained by binary adding the contents of a pair of tapped G2 cells and binary adding that result to the output of G1. That becomes the C/A-code. The various alternative pairs of G2 taps (delays) are used to generate the complete set of 36 unique PRN C/A-codes. There are actually 37 PRN C/A-codes, but two of them (34 and 37) are identical. The first 32 codes are assigned to satellites. Codes 33 through 37 are reserved for other uses such as for ground transmitters. This family of codes is a subset known as Gold codes, which have the property that any two have a very low cross correlation (are nearly orthogonal). The term for the codes comes from the inventor, Robert Gold, not from their lustrous properties.

The C/A-code is modulated only onto the L1 carrier, unlike the P(Y)-code (see 235,469,592,765,000) which appears on both L1 and L2. However, beginning with the first Block IIR-M or modernized Block IIR satellite, a new civil code, L2C or L2 Civil, has been transmitted on the L2 frequency. The future Block IIF satellites will also transmit L2C.

1023 is also the maximum value of the GPS week. This is the number of full weeks that have elapsed since the GPS Time zero point of midnight UTC beginning January 6, 1980 — but with a special counting procedure. GPS weeks are numbered consecutively with week zero starting on January 6 and ending on January 12, 1980. The GPS week, together with the Z-count (see 403199), specifies an epoch or event related to GPS signals or measurements. The current GPS week is included in subframe one of the navigation message, which — along with other subframes containing satellite clock, ephemeris data, and other user-required information — is transmitted every 30 seconds. Only 10 bits are used to represent the GPS week, and so the largest possible week number is 1023 (2¹⁰–1). In other words, the GPS week number is modulo 1024 (see the “Modular Arithmetic” sidebar). At the end of week number 1023, the week number rolls over to zero. This first occurred on August 21/22, 1999, and caused difficulties for some GPS receivers as their manufacturers had failed to account for the “end-of-week rollover” in receiver firmware. The next occurrence will be in April 2019. By that time, the new Civil Navigation (CNAV) message will be in use, in which the GPS week number is represented as a 13-bit value, meaning it rolls over after 8192 weeks, or about every 157 years.

Although officially the GPS week number is still modulo 1024, some agencies, such as the International GNSS Service, prefer to use a running count of the GPS week, ignoring the rollover. The number of the week beginning April 4, 2010, is then alternatively given as 1578 or 554. Or, mathematically speaking (see sidebar), 1578 ≡ 554 (mod 1024).

1176.45

The L5 carrier frequency in megahertz. The L5 carrier frequency is obtained by electronically multiplying the satellite 10.23 MHz standard frequency by 115. It is the lowest and the newest of the GPS frequencies and is used for the new civil-only GPS signal.

The addition of the L5, or Link 5, civil signal to the suite of signals transmitted by the satellites is a key feature of GPS modernization. The introduction of such a signal on a different carrier frequency than that used by the legacy L1 GPS signal was proposed in the 1995 reports by the U.S. National Research Council and the National Academy of Public Administration on the future of GPS. The reports argued that an unencrypted signal on a second frequency would offer civil users the benefit of ionospheric delay correction, wide-lane carrier-phase ambiguity resolution, improved interference rejection, and faster accuracy recovery in multipath environments. The frequency is in a protected aeronautical radionavigtion services band and, unlike L2, means that L5 can be used for safety-of-life services. The L5 signal will be standard on all Block IIF and future satellites. An L5 demonstration payload was included on Block IIR-M satellite SVN49 to secure the L5 frequency under the rules of the International Telecommunication Union.

1227.60

The L2 carrier frequency in megahertz. The L2, or Link 2, carrier is modulated with the P(Y)-code. Additionally, starting with the Block IIR-M satellites, a new civil ranging code, L2C, is being transmitted on L2 along with the new military M-code. These new signals are also part of the GPS modernization effort. The L2 carrier frequency is obtained by electronically multiplying the satellite 10.23 MHz standard frequency by 120.

1381.05

The L3 carrier frequency in megahertz. This frequency is used in conjunction with the GPS satellites’ secondary purpose, which is to detect nuclear detonations. The L3 carrier frequency is obtained by electronically multiplying the satellite 10.23 MHz standard frequency by 135.

1575.42

The L1 carrier frequency in megahertz. The L1, or Link 1, carrier is modulated with the C/A-code and the P(Y)-code. Starting with the Block IIR-M satellites, the new military M-code is also transmitted on L1. The L1 carrier frequency is obtained by electronically multiplying the satellite 10.23 MHz standard frequency by 154. If you’ve been counting, you’ll have noticed that we didn’t list an L4 frequency. L4 has never been implemented but it has been studied. For example, a frequency of 1841.40 MHz (10.23×180) was once considered for ionospheric correction.

403199

The maximum value of the GPS time of week count. The GPS satellites count and communicate GPS Time in a unique manner that is ultimately related to how they generate the PRN ranging codes. As described below, the P-code is generated by combining two shorter PRN codes, X1 and X2, which are clocked in phase at a chipping rate equal to the satellite’s 10.23-MHz oscillator frequency. X1 has a repetition interval, or period, of 1.5 seconds — a fundamental GPS timing unit. The start of each 1.5-second interval identifies an epoch. The number of X1 epochs since the beginning of the week is called the time of week (TOW) count, which runs from zero to 403,199 at the end of week. The TOW count returns to zero coincident with the resetting of the PRN codes.

The TOW count can be represented as a 19-bit binary number, a truncated version (the 17 most significant bits) of which is part of the handover word (HOW) that a satellite transmits every six seconds. The HOW appears as the second word in each data subframe of the navigation message. These 17 bits correspond to the TOW count at the X1 epoch that occurs at the start of the immediately following subframe, and so effectively preannounces the arrival of a time marker, just like telephone “speaking clocks” and shortwave radio time and frequency stations.

The TOW count by itself cannot be used to unambiguously establish the date of an event. It can only time an event at modulo 604,800 seconds [(403199+1)x1.5] because it is reset every week. This time ambiguity is reduced by noting the number of full weeks that have elapsed since January 6, 1980 modulo 1024 — the GPS week number (see 1023). The TOW count and the GPS week number combine to form the 29-bit Z count. The 19 least-significant bits are the TOW count and the 10 most-significant bits are the GPS week number.

299,792,458

The speed of light in meters per second. This is the speed with which all electromagnetic radiation propagates in a vacuum. Until 1983, the speed of light was measured experimentally using adopted standards for the length of the second and the length of the meter. However, compared to the second, the definition of the meter had a large uncertainty. So in 1983, the 17th General Assembly of Weights and Measures defined the meter as the distance travelled by light in a vacuum during 1/299,792,458 of a second, fixing the speed of light at 299,792,458 meters per second — exactly. This constant is used by a GPS receiver, for example, to convert the measured signal propagation time in seconds to a pseudorange in meters.

235,469,592,765,000

(Or 2.35469592765000 x 10¹⁴ in scientific notation). This is the number of chips in the P-code if it were allowed to continue without being reset. The P-, or precision code, is one of the two legacy PRN ranging codes that have been transmitted by all GPS satellites. The other is the C/A-code, already discussed.

The P-code is actually the product of two PRN codes, each of which is generated with a pair of feedback registers. The X1 code has a length of 15,345,000 chips while the X2 code has a length of 15,345,037 or 37 chips longer. So the complete P-code has a length equal to the product of the lengths of the X1 and X2 codes, or 235,469,592,765,000. The codes are clocked at a rate of 10.23 MHz so that each chip has a length of about 0.097752 microseconds. This means the pattern of chips in the full P-code would not repeat for almost 266 days. Each satellite is assigned a unique one-week segment of the P-code, which is reset at Saturday/Sunday midnight each week. The individual P-codes have low cross-correlations with each other. In other words, no significant segments of the P-code of one satellite matches that of another.

Before transmission, a P-code chip sequence is encrypted to form a new sequence called the Y-code. The combined sequence is usually referred to as P(Y). Although civil GPS receivers cannot use conventional correlation procedures to acquire and track the P(Y) code, they can use knowledge of the underlying P-code sequence and C/A-code tracking on L1 to produce pseudorange and carrier-phase measurements on both the L1 and L2 frequencies.

√-1

The square root of -1. This is the unit of imaginary numbers. The concept of imaginary numbers, actually known to the ancient Greeks, was introduced in the effort to solve algebraic equations. Not all equations can be solved using real numbers. In particular, x²=-1 has no real-valued solution. But we can say some solution exists and represent it by the symbol i. (Mathematicians and physicists use this symbol whereas electrical engineers prefer j, since i usually describes a varying electrical current.) Then i has the property — by definition — that its square is -1. Of course, that equation would also permit the solution -i. An imaginary number — sometimes called pure imaginary — is any number of the form bi, where b is a non-zero, real number. A real number, a, and an imaginary number, bi, can be combined into a complex number, a+bi, or a+ib, the more usual notation. Using complex numbers and a set of rules governing their manipulations, any algebraic equation can be solved.

It is useful to consider the real and imaginary parts of a complex number to be orthogonal so that we can represent a complex number geometrically on a plane — the complex plane — where the real component is plotted on the x-axis and the imaginary component on the y-axis. We can then represent a 2-dimensional vector as a complex number, with one component considered real and the other imaginary. The magnitude or modulus of the vector, r, is the positive square root of the sum of the squares of the real and imaginary components with the vector making an angle, Φ , with respect to the positive real axis.

It can be easily shown that

This is Euler’s famous formula, which provides an enlightening connection between plane geometry and algebra.
And, we may also write any complex number in the form

or even more compactly as

If the vector rotates counterclockwise with angular speed ω, its projection onto the real axis generates a sine wave. The modulus of this vector is the amplitude of the oscillations, while its argument is the total phase,

where t is time. The phase constant θ represents the angle that the vector forms with the real axis at t=0. This representation of a sine wave as a phase vector, or phasor, finds great utility in signal theory including descriptions of the propagation of radio waves such as those emitted by GPS satellites.

Modular Arithmetic

GPS Time, like all time systems, is based on modular arithmetic. This arithmetic is a little different from conventional arithmetic in that numbers, typically restricted to integers, have a finite maximum value. Adding one to that number doesn’t get you a larger number — it gets you a smaller one, a much smaller one: zero.

Modular arithmetic is known to us all as clock arithmetic. Take the 24-hour time system as an example. If it’s currently 1800, then 8 hours later we say it’s 0200, not 2600. Similarly, if it’s currently 0400, then 6 hours earlier it was not -0200 but 2200. The idea here is that if two numbers differ by 24 or a multiple of 24, then they are “equal.” We could simply write 26 = 2 but this could be confusing. So we write 26 ≡ 2 (mod 24), and -2 ≡ 22 (mod 24), or in words, 26 is congruent (or somehow “equal”) to 2 (modulo 24) and -2 is congruent to 22 (modulo 24). In arithmetic modulo 24, any number larger than 24 is congruent to some number less than 24 because we can always subtract a multiple of 24 from the larger number to get the smaller one. Similarly, any negative number is congruent to some positive number less than 24, and 24 is congruent to 0. This means that in arithmetic modulo 24, we need deal only with integers from 0 to 23.

We can choose any positive integer for the modulus and carry out arithmetic operations accordingly. Using a modulus of 4, for example, we would have 2+2=0 in our loose notation — a disturbing result if interpreted as conventional arithmetic. But when written 2+2 ≡ 0 (mod 4), the meaning is clear.

As another example of modular arithmetic, consider this question: If today is Monday, what day of the week is it 185 days from now? The modulus here of course is 7, the number of days in the week. So, mathematically stated: 1+185 ≡ ? (mod 7). The answer: 4 or Thursday. The answer is obtained by dividing the sum on the left side of the congruency by 7, using “long division,” and noting the remainder. Or, alternatively, the sum is divided by 7, and the decimal part of the result is then multiplied by 7.

An interesting quirk of modular arithmetic is that a number and the sum of its digits are congruent, modulo 9. This property is the basis for a formerly well-known procedure (before the days of calculators and computers) for checking the correctness of hand multiplication — the rule for casting out nines, which states that the product of two numbers and the product of the sums of their digits must have the same remainder on division by 9.

Many computer languages have a built-in modular arithmetic function or operator. Typically called MOD, it returns the remainder from an integer division operation. In BASIC, for example, if we enter 5 MOD 2, we get 1 because 5 divided by 2 is 2 with a remainder of 1. The same computation is coded 5 2 MOD in Forth, 5 % 2 in Python, and MOD (5,2) in FORTRAN.

The following line of FORTRAN code by Henry Fliegel of The Aerospace Corporation inherently uses modular arithmetic by way of integer division to determine the Julian day (JD) number from the year, month, and day of an AD Gregorian calendar date, incorporating all leap year rules:

JD=367*Y-7*(Y+(M+9)/12)/4-3*((Y+(M-9)/7)/100+1)/4+275*M/9+D+1721029

And just how can the GPS end-of-week rollover be described using modular arithmetic? Very simply: 1023+1 ≡ 0 (mod 1024).

FURTHER READING

• GPS Interface Control Documents

Navstar GPS Space Segment / Navigation User Interfaces, Interface Specification, IS-GPS-200 Revision D (IRN-200D-001), prepared by ARINC Engineering Services LLC, El Segundo, California, March 2006.

Navstar GPS Space Segment / User Segment L5 Interfaces, Interface Specification, IS-GPS-705 Revision IRN-705-003, prepared by ARINC Engineering Services LLC, El Segundo, California, 22 September 2005.

Navstar GPS Space Segment / User Segment L1C Interfaces, Interface Specification, IS-GPS-800,
prepared by Science Applications International Corporation, El Segundo, California, 4 September 2008.

• Imaginary Numbers

An Imaginary Tale: The Story of √-1 by Paul J. Nahin, Princeton University Press, Princeton, New Jersey, paperback edition, 2007.

• Innovation Column 1

“GPS: A Multipurpose System” by D. Wells and A. Kleusberg in GPS World, Vol. 1, No. 1, January/February 1990, pp. 60–63. (Scanned versions of the first 15 Innovation columns are available.)

• Innovation Column 100

“Smaller and Smaller: The Evolution of the GPS Receiver” by R.B. Langley in GPS World, Vol. 11, No. 4, April 2000, pp. 54–58.

• Some Other Numerically Relevant Innovation Columns

“Why is the GPS Signal so Complex?” by R.B. Langley in GPS World, Vol. 1, No. 3, May/June 1990, pp. 56–59.

“The Orbits of GPS Satellites” by R.B. Langley, in GPS World, Vol. 2, No. 3, March 1991, pp. 50–53.

“Time, Clocks, and GPS” by R.B. Langley in GPS World, Vol. 2, No. 10, November/December 1991, pp. 38–42.

“Detecting Nuclear Detonations with GPS” by P. Highie and N. K. Blocker in GPS World, Vol. 5, No. 2, February 1994, pp. 48–50.

“The Promise of a Third Frequency” by R.R. Hatch in GPS World, Vol. 7, No. 5, May 1996, pp. 55–58.

“The GPS End-of-Week Rollover” by R.B. Langley in GPS World, Vol. 9, No. 11, November 1998, pp. 40–47.

“Tropospheric Delay: Prediction for the WAAS User” by P. Collins and R.B. Langley in GPS World, Vol. 10, No. 7, July 1999, pp. 52–58.

“GPS, the Ionosphere, and the Solar Maximum” by R.B. Langley in GPS World, Vol. 11, No. 7, July 2000, pp. 44–49.

“The New L5 Civil GPS Signal” by A.J. Van Dierendonck and C. Hegarty in GPS World, Vol. 11, No. 9, September 2000, pp. 64–71.

“Time for a Better Receiver: Chip-Scale Atomic Frequency References” by J. Kitching in GPS World, Vol. 18, No. 11, November 2007, pp. 52–57.

“The GPS L2C Signal: A Preliminary Analysis of Data Quality” by R.F. Leandro, R.B. Langley, L. Sükeová, T. Thirumurthi, and M.C. Santos in GPS World, Vol. 19, No. 10, October 2008, pp. 42–47.

“GPS L5 First Light: A Preliminary Analysis of SVN49’s Demonstration Signal” by M. Meurer, S. Erker, S. Thölert, O. Montenbruck, A. Hauschild, and R.B. Langley in GPS World, Vol. 20, No. 6, June 2009, pp. 49–58.

April 1, 2010
Innovation: Hybrid Positioning

A Prototype System for Navigation in GPS-Challenged Environments

By Chris Rizos, Dorota A. Grejner-Brzezinska, Charles K. Toth, Andrew G. Dempster, Yong Li, Nonie Politi, Joel Barnes, Hongxing Sun, and Leilei Li

A team of Australian and U.S. researchers have integrated a ground-based system with GPS and INS to create a hybrid system that provides precise and accurate position information continuously in a variety of environments where GPS alone comes up short.

INNOVATION INSIGHTS by Richard Langley

GPS HAS ITS LIMITATIONS. Although it is a 24/7 global system, it doesn’t work everywhere. The microwave radio signals transmitted by the satellites are rather weak, and although they can provide excellent positioning performance when a receiver’s antenna has a direct line-of-sight view of a sufficient number of satellites well spread out in the sky, positioning accuracy degrades or becomes impossible when the signals are effectively blocked by obstacles such as trees, rock faces, and buildings outdoors and by roofs, ceilings, and walls indoors.

In many obstructed environments, the signals aren’t completely blocked but rather their power is severely attenuated so that they are no longer strong enough to be acquired and tracked by a conventional GPS receiver. Remarkable progress has been made in the development of super-sensitive receivers that, in conjunction with an appropriate antenna and assistance information provided over a mobile phone network, can provide position fixes in such environments. However, the precisions and accuracies of these pseudorange-based positions are often very poor — perhaps as low as 100 meters or more.

So, is it possible to obtain precise and accurate positions in obstructed environments? Well, we could add measurements from GLONASS (or other satellites) to GPS measurements, but GLONASS suffers the same problem as GPS, and while the additional satellites could be an advantage in some partially obscured areas there are many places where we won’t be any better off. We could use an inertial navigation system (INS), but such devices have their own weaknesses such as the requirement of initial calibration and the accumulation of position error with time. Are there any other technologies available?

We know GPS works very well when there is a direct line-of-sight view between the satellite transmitters and the receivers and carrier-phase measurements can provide decimeter- and even centimeter-accuracies. So why not develop a ground-based system that works in a similar way to GPS, which would allow you to place the transmitters wherever you like? Well, such a system has indeed been developed and in this month’s column, a team of Australian and U.S. researchers describes how they integrated the ground-based system together with GPS and INS to create a hybrid system that provides precise and accurate position information continuously in a variety of environments where GPS alone comes up short.

“Innovation” features discussions about advances in GPS technology, its applications, and the fundamentals of GPS positioning. The column is coordinated by Richard Langley, Department of Geodesy and Geomatics Engineering, University of New Brunswick.

The determination of the position and orientation (or “pointing direction”) of a device (or platform to which it is attached), to high accuracy, in all outdoor environments, using reliable and cost-effective technologies is something of a “holy grail” quest for navigation researchers and engineers.

However, ongoing research has identified two classes of applications that place stringent demands on the positioning/orientation device: (a) man-portable mapping and imaging systems that operate in a range of difficult urban and rural environments, often used for the detection of underground utility assets (such as pipelines, cables, conduits), unexploded ordnances and buried objects, and (b) the guidance/control of construction or mining equipment in environments where good “sky view” is not guaranteed.

The solution to this positioning/orientation problem is increasingly seen as being based on an integration of several technologies: satellite (GNSS including GPS) and terrestrial ranging systems, inertial navigation systems (INSs), laser guidance/scanning systems, and even electro-optical devices such as surveyors’ total stations or laser scanners. Each has its shortcomings, but within an integrated system, advantage can be taken of the complementary characteristics of several of these sensor technologies.

Centimeter-level accuracy positioning systems for outdoor use typically have at their core the GPS technology. GPS is, in fact, the most effective general-purpose navigation tool ever developed because of its ability to address a wide variety of applications: air, sea, land, and space navigation; precise timing; geodesy; surveying and mapping; machine guidance/control; military and emergency services operations; hiking and other leisure activities; personal location; and location-based services. The varied applications use different and appropriate receiver instrumentation, operational procedures, and data processing techniques. But all require signal availability from a minimum of four GPS satellites for three-dimensional fixes.

However, one of the usual limiting factors in using GPS is the need for direct line-of-sight between the satellites and the ground receiver. In particular, the robustness of positioning is compromised when GPS receivers are near or under trees, in urban/suburban areas, or in deep open-pit mines and construction sites, where there is partial sky view obstruction by buildings or walls. The traditional means of overcoming the gaps in navigation coverage due to satellite signal blockages is to use an INS. An INS (with its inertial measurement unit or IMU) is also the most convenient means of determining the orientation of the device or platform. The integration of GPS and INS can, in principle, overcome the defects of standalone INS (sensor errors that grow unbounded with time) and GPS (signal availability requirement). But navigation accuracy degrades rapidly if there are no GPS measurements to calibrate the INS sensor errors.

A new terrestrial RF-based distance measurement technology offers promise of continuous signal coverage, even in difficult urban/rural environments. This technology is known as “Locata.”

The Locata approach is to deploy a network of ground-based transceivers that cover an area with strong time-synchronized ranging signals. When a Locata receiver uses four or more ranging signals it can compute a high-accuracy position entirely independent of GPS or INS. However, a standalone Locata receiver has its own shortcomings: (a) in some situations it may be difficult to achieve good vertical dilution of precision due to logistical constraints of placing transmitters (to give a variation in elevation angle between the terrestrial transmitters and the receiver whose positions are to be determined), and (b) as with GPS, multiple receivers/antennas are required to derive orientation information.

What is therefore required is several carefully selected navigation sensor technologies, integrated within a single hardware package, the measurements from which are simultaneously processed to provide continuous, reliable, and accurate navigation solutions (that is, both position and orientation information).

In cooperation with Locata Corporation, the SNAP Laboratory within the School of Surveying and Spatial Information Systems at the University of New South Wales (UNSW) and the SPIN Laboratory at The Ohio State University have assembled a working prototype of a hybrid system based on GPS, inertial navigation, and Locata receiver technology to provide seamless and reliable navigation aimed at supporting vehicle guidance and control, open-pit mining, mobile and GIS mapping, and industrial applications.

Locata Technology

The SNAP Lab has been conducting pseudolite research for many years, and has experimented with pseudolites in nonsynchronous and synchronized modes for a variety of applications, using both the GPS L1 frequency as well as the 2.4 GHz ISM band frequencies. Locata Corporation has developed state-of-the-art RF terrestrial positioning technology (“Locata”), which consists of a network (“LocataNet”) of time-synchronized pseudolite-like transceivers (“LocataLites”). UNSW has assisted in the development of the technology through experimental testing and benchmarking. In a relatively open outdoor environment, the LocataNet can provide real-time stand-alone kinematic positioning (without a base station) at centimeter-level accuracy. Even in an indoor environment where LocataLite signals arrive at a Locata receiver via non-line-of-sight paths (penetrating the walls of buildings), the static positioning quality can be at the sub-centimeter level, and also at the sub-meter level for kinematic positioning.

Locata has several advanced features that have been developed over a period of about 10 years through several technology generations, including a time-synchronized positioning network, network propagation to many LocataLites, improved signal penetration, change of transmitting frequency and signal structure, and spatial and frequency diversity.

In TABLE 1, the key characteristics of the two generations of Locata technology are listed. Using 2.4 GHz not only means the frequency is license-free, but also permits transceiver output power of up to 1 watt, which means greater operating distances (up to 10 kilometers). Using dual-frequency signals changes the initial phase-bias resolution from known-point initialization to on-the-fly (OTF), where the initial phase bias is resolved while the receiver is moving. The higher chipping rate (10 MHz) results in less pseudorange multipath error, because the delay in a reflected signal will rarely be more than two chips. The 10-Hz measurement rate allows relatively high velocities of the receiver.

Table 1. Specification summary of Locata’s first- and second- generation systems.

In terrestrial-based RF-based positioning, multipath error is more severe than with GPS, because the terrestrially transmitted signal arrives at the receiver at a very low (typically less than 10 degrees) or even a negative elevation angle, which can result in severe multipath signal fading. In the second-generation Locata system, spatial and frequency diversity techniques are employed. Spatial and frequency diversity are two of the three types of diversity principles (the other being polarization) that are common practices in terrestrial RF communications to mitigate against signal fading. The LocataLite transceiver uses two spatially separated (usually in the vertical) antennas, which transmit two signals at different frequencies. This gives a cluster of four diverse signals transmitted from one LocataLite. With this diversity technology, Locata kinematic positioning in moderately obstructed environments can provide centimeter-level quality with 100-percent coverage, as well as consistent geometry and high reliability. The Locata’s multipath mitigation technology is very important and relevant to this project, because the operational environments are often vegetated or wooded.

Triple Integration

As discussed in the preceding sections, there are both advantages and disadvantages to every navigation sensor. GPS and Locata have high positioning accuracy in open or moderately obstructed environments, but they are sensitive to signal blockage such as the case in dense forests, urban canyons, deep mine pits, and indoors. In contrast, INS is totally autonomous — that is, independent of external signal sources — and has high output rate for position, velocity, and attitude, but its unaided navigation error grows rapidly with time.

The most common data-processing tool to integrate GPS and INS is the Kalman filter, which forms the basis for multi-sensor integration in this research. The basic Kalman filter applies to linear system models. Therefore, several variations were developed to cope with the non-linear navigation model, such as the extended Kalman filter and the unscented Kalman filter.

The following discussion of the integration of the GPS/INS/Locata sensors is focused on two aspects: 1) the system state selection, and 2) the measurement model or integration model that decides which information to pass to the filter.

The error state vector consists of a nine-dimensional navigation error state sub-vector (three for the position, three for the velocity, and three for the orientation), an accelerometer error state sub-vector, a gyroscope error state sub-vector, and a three-dimensional gravity disturbance state sub-vector. Of course, other sensor error models can be considered for the gyroscope and accelerometer sensors, such as a combination of random constants, first-order Gauss-Markov variables, scale factors, and so on. In this case, the state space could have a dimension of more than 30. The objective is to adjust the sensor error model later based on experimental results (if needed). However, because of the limitations of observability, it is not yet known whether an augmented error state vector would give better results.

When integrating INS hardware with other sensors, the sensors cannot share the same physical location, which would be ideal from a theoretical point of view. Knowing the spatial relationship among the sensors is important to ensure the highest possible navigation performance. The displacement vectors or mounting biases are offsets, also referred to as lever arms, from the center of the IMU to the centers of the other sensors. These lever-arm parameters may be included in the Kalman filter and thus can be estimated. However, if the lever arms are precisely measured during the assembly of the system, they do not need to be included in the filter as estimable parameters.

For multiple sensor integration in a Kalman filter, there are essentially two types of general models: loosely coupled and tightly coupled. The loosely-coupled model uses a decentralized filter that has several sub-filters to process the sub-systems independently. In other words, the Kalman filter solutions from the sub-systems are combined in an overall Kalman filter that provides the integrated navigation solution. In contrast, the tightly-coupled model uses a single main filter to process the output of all sensors. In GPS/INS integration, tightly-coupled systems have obvious advantages in environments where GPS signals are frequently lost, because they can rely on the other sensor(s) when GPS positioning becomes impossible.

In the tightly-coupled model, the raw observations of all sensors will be input to the main filter. For GPS and Locata, the primary observations will be the carrier phase measurements, as code (pseudorange) observations cannot provide the required accuracy. High-accuracy GPS positioning needs to address the issue of carrier-phase ambiguity. The ambiguity can be treated as an unknown in the Kalman filter, but it may take several minutes to resolve the ambiguity using GPS alone. Using certain ambiguity resolution techniques, however, the ambiguity can be resolved outside the main filter in the GPS/INS high-precision (carrier-phase) integration filter. Note that if the ambiguity were to be resolved within the filter, this would increase the number of states of the filter. For the GPS component, ionospheric delay should be included for applications that cover a large area. Ionospheric delay can be resolved using network-based differential techniques,
but it will affect the ambiguity resolution for single baseline differential positioning if it is not included in the local solution. The filter is designed either to use, or not to use, ionospheric delay, which can ensure flexibility to accommodate network-based and single-baseline differential positioning.

As mentioned above, the measurement model in the tightly-coupled model is based on the raw observations. For GPS and Locata, the observations will be the carrier-phase observations. The approximate values for the linearization of the GPS and Locata measurement equations are provided by the INS navigation solution.

The GPS carrier-phase ambiguity is solved independently outside the Kalman filter with OTF techniques. The GPS differential positioning coefficient matrix remains the same regardless of whether or not a network-based differential technique is used. For velocity determination, the double-differenced Doppler observation is used to eliminate the clock error rate as an unknown (because it is difficult to model this in the filter). The initial carrier-phase bias of the Locata is also not included in the filter, because it can be resolved instantaneously with dual-frequency data in the Locata second-generation system.

The implementation of the filter will be flexible, so adjustments can be made to account for actual environmental conditions. The filter is designed with an open interface and is modular in structure, so that components can be added (or removed) from the model. In open-sky areas, GPS is sufficient for system positioning, so only its observations need to be processed. In moderately obstructed environments, GPS and Locata observations will be processed. In this case the number of GPS observation equations is limited and sometimes will be less than four. FIGURE 1 illustrates the flowchart of the triple-integration of GPS, INS, and Locata.

Figure 1. Workflow of the integrated GPS/ INS/Locata system.

Field Tests

For experimental purposes, we used a dual INS, based on a navigation grade unit and a tactical grade unit. In addition, a Locata receiver and a dual-frequency GPS receiver were placed on a vehicle at Locata’s Numeralla Test Facility (NTF) near Canberra, Australia. This test site features both open-sky and obscured environments, allowing for testing the system’s performance under truly challenging scenarios. The test was repeated by mounting the devices on an autonomous electrical car, driven on the UNSW campus. In both cases, the separation between the rover and the terrestrial transmitters was between a few tens of meters to several kilometers. The GPS and Locata data were processed separately (for testing the internal consistency) as well in a hybrid solution, resulting in few-centimeter-level accuracy per coordinate, depending primarily on GPS availability and the geometry between the rover and Locata devices, as well as the level of multipath fading.

Test 1: NTF. The first integration test was conducted at the NTF on March 17, 2008. The NTF covers an area of approximately three hundred acres (2.5 kilometers × 0.6 kilometers) and is ideally suited to real-world system testing over a wide area. At the NTF, a number of LocataNet configurations are possible through the installation of permanent antenna towers. The network configuration used for this experiment is illustrated in FIGURE 2.

Figure 2. NTF: LocataLite network.

Before the test, a special mounting platform was designed and built. The platform, shown in FIGURE 3, consists of a two-level metal frame. The bottom level can accommodate two inertial measurement units, while the top level can hold up to four antennas. The platform can be easily attached to either the roof of the NTF test vehicle or to the body of UNSW’s small electric car (described later).

Figure 3. Devices setup for the NTF test.

The devices used in the test include two dual-frequency GPS receivers (one used as the rover receiver and the other as the base station), one navigation grade INS, and one Locata rover unit. The GPS antenna and the Locata antenna were mounted with the INS together on the top of a truck. The GPS data rates were set to 1 Hz. The average length of the GPS differential baselines was about 1.2 kilometers. The GPS observation conditions were good during the testing period. The Locata data rate was set to 10 Hz, while INS data rate was 256 Hz, and both were synchronized with the GPS time using SNAP-Lab-developed time synchronization devices based on field-programmable gate array (FPGA) technology.

The GPS/INS data were first processed in tightly-coupled mode. The trajectory is depicted in FIGURE 4. The standard deviation of position, velocity, and attitude are shown in FIGURES 5-7 respectively.

Figure 4. The trajectory of the vehicle in the NTF test

Figure 5. The standard deviation of position in the test.

Figure 6. The standard deviation of velocity in the test.

Figure 7. The standard deviation of attitude in the test.

In Figures 5-7, it can be seen that the standard deviations of position and velocity are less than 0.02 meters and 0.01 meters per second respectively. The standard deviations of pitch and roll angles are less than 0.001 degrees as well as that of yaw, which is less than 0.01 degrees after the vehicle starts to move, at about the 1500th second.

The Locata data was post-processed using Locata’s Integrated Navigation Engine (LINE). It provides an unsmoothed single point position using carrier-phase measurements. The initial ambiguity bias was resolved using the data from the GPS carrier-phase position. Following this initialization, the Locata solution was computed independently of GPS. A 15-meter tower LocataLite location in the vicinity of the start and end of the test (indicated by the “figure eight” pattern in FIGURE 8) allowed sufficient geometry for 3D positioning using Locata. For the rest of the data where there was insufficient vertical geometry, GPS height aiding was used. Figures 8 and 9 show the independent Locata and GPS solutions (without lever arm correction) for the section of the trajectory in the vicinity and the end of the test, respectively. The Locata solution compared to the GPS solution to within a few centimeters for the entire trajectory.

Figure 8. Section of trajectory showing independent Locata solution (black) vs. GPS (blue) with no lever-arm correction.

Figure 9. End of trajectory showing independent Locata solution (black) vs. GPS (blue) with no lever-arm correction.

To test the GPS/INS/Locata integration, some GPS observation epochs were deleted to simulate two GPS blockages from seconds of week 94100 to 94250 and from 94500 to 94600. The INS standalone navigation errors with this deleted GPS data were about 8 meters and 2.6 meters, respectively.

In the final GPS/INS/Locata integration test, Locata compensated for the missing GPS data. The integration result was almost identical to the GPS/INS integration result obtained with the original GPS observed data clearly showing that the Locata system could seamlessly replace GPS in this scenario.

Test 2: Electric Car. Early in 2007, UNSW researchers established a permanent LocataNet on the university campus to provide a research and test facility at UNSW devoted to the Locata technology. The LocataNet setup at UNSW is illustrated in FIGURE 10. It consists of four dual-frequency LocataLites situated on tops of four buildings surrounding a lawn test area. The master LocataLite is on the Civil Engineering building and the other three LocataLites are synchronized to it.

Figure 10. LocataLites on the UNSW campus.

Currently, to be able to obtain a carrier-phase position solution with Locata, the initial ambiguities need to be resolved by initializing the rover receiver on a known position. For this purpose, a point in the middle of the test area was surveyed, and the coordinates were used to initialize the Locata receiver.

SNAP Lab has developed a small electric car that can be driven using an attached handheld controller (see FIGURE 11). The controller enables the car to move in both forward and reverse and to steer the front wheels.

Figure 11. The electronic car used in the test.

For these tests, the same mounting platform as the one used in the previous experiment allowed all the sensors and ancillary equipment to be attached to the car. For this experiment, we used the following equipment: a Locata receiver, two GPS receivers, a tactical grade INS, a 360-degree prism (tracked by a robotic total station), and two time-sync FPGA data-logging devices.

The starting position was the known point in the middle of the Locata network. The car was then driven in a circular path three times before finishing back at the starting position.

During the test the raw data stream from the Locata receiver, the GPS receivers, and the INS were recorded using the time-sync data-logging devices. In addition, a robotic total station (RTS), which was set up at the edge of the test area, automatically tracked the prism position (the data was recorded internally).

The Locata data was post-processed using LINE to give a single point unsmoothed carrier-phase solution. The initial ambiguity bias was resolved using the data from the GPS carrier-phase position. Following this initialization, the Locata solution was computed independently of GPS. Where there was insufficient vertical geometry (at the very west end of the trajectory shown in FIGURE 12), GPS height aiding was used. The Locata-only solution and the RTS result are shown in Figure 12. The two solutions compare to within a few centimeters of each other.

Figure 12. The trajectory from the Locata-only and robotic total station solutions.

We then carried out the integrated GPS/INS processing. To test the GPS/INS/Locata integration, two GPS outages were simulated by simply removing the data from the GPS file, between seconds of week 103703 and 103713 and 103834 and 103844, respectively.

We then carried out the integrated GPS/INS processing. To test the GPS/INS/Locata integration, two GPS outages were simulated by simply removing the data from the GPS file, between seconds of week 103703 and 103713 and 103834 and 103844, respectively.

In comparison to the original GPS/INS integration, the standalone INS solution has errors of about 35 meters and 12 meters during the first and second outages, respectively.

The Locata/INS integration significantly reduced the navigation error during the GPS outages, as summarized in TABLE 2.

Table 2. The difference between the Locata/INS solution and the original GPS/ INS solution

From Table 2 it can be seen that 3D position differences between the Locata/INS and the original GPS/INS integration result have been reduced to 1.143 meters and 0.053 meters during the two GPS outages, respectively. However, the improvement in the accuracy of the attitude angles is not obvious because a 10-second GPS outage is not long enough to cause a significant INS drift.

Concluding Remarks

The test experiments described here are a demonstration of the proof-of-concept of a triple-integration GPS/INS/Locata system. The navigation results indicate that this sensor combination may support navigation in GPS-denied environments, as long as some partial view of the LocataLites within the network is available. Further development of this triple integration system is being undertaken.

Acknowledgments

The research is funded by the Australian Research Council. This article is based on the paper “A Hybrid System for Navigation in GPS-challenged Environments: A Case Study,” presented at ION GNSS 2008, the 21st International Technical Meeting of the Satellite Division of The Institute of Navigation, Savannah, Georgia, September 16-19, 2008.

Manufacturers

The Numerella test equipment included Locata devices, a Honeywell H-764G navigation-grade INS, a Boeing (now Systron Donner) C-MIGITS II tactical grade INS, and a Leica System 1200 dual-frequency GPS receiver. The UNSW campus test equipment included Locata devices, an Omnistar GPS receiver, a Leica MC500 GPS receiver, a Boeing C-MIGITS II INS, a Leica GRZ4 360-degree prism, and a Leica robotic total station TCRP 1203+.

CHRIS RIZOS is a graduate of the University of New South Wales (UNSW), Sydney, Australia, where he obtained a Ph.D. in satellite geodesy. He is head of the School of Surveying and Spatial Information Systems at UNSW.

DOROTA BRZEZINSKA is a professor and leader of the Satellite Positioning and Inertial Navigation (SPIN) Laboratory at The Ohio State University (OSU) in Columbus, Ohio. She received her M.S. and Ph.D. in geodetic science from OSU.

CHARLES TOTH is a senior research scientist at OSU’s Center for Mapping. He received a Ph.D. in electrical engineering and geo-information sciences from the Technical University of Budapest, Hungary.

ANDREW G. DEMPSTER is the director of research in the School of Surveying and Spatial Information Systems at UNSW.

YONG LI is a senior research fellow at the SNAP Lab. He obtained a Ph.D. in aerospace engineering.

NONIE POLITI is a graduate of the School of Electrical Engineering and Telecommunications at UNSW. He obtained a Bachelor’s degree in Telecommunication Engineering and an M.Eng.Sc. in electronics.

JOEL BARNES is director of navigation R&D for Locata Corporation and is also a senior visiting research fellow at the SNAP Lab.

HONGXING SUN is a post-doctoral researcher in the SPIN Lab. He received a bachelor’s degree in geodesy and M.S. and Ph.D. degrees in photogrammetry from Wuhan University, China.

LEILEI LI is a Ph.D. candidate at Chongqing University, China. He is also a visiting Ph.D. student in the SPIN Lab. He received an M.S. degree in instrument science and technology from Chongqing University.

FURTHER READING

• Locata

“Locata: A New Technology for High Precision Positioning” by N. Politi, Y. Li, F. Khan, M. Choudhury, J. Bertsch, J.W. Cheong, A. Dempster, and C. Rizos in Proceedings of ENC-GNSS 2009, the European Navigation Conference, Naples, Italy, May 3-6, 2009.

“Deploying a Locata Network to Enable Precise Positioning in Urban Canyons” by J.-P. Montillet, G.W. Roberts, C. Hancock, X. Meng, O. Ogundipe, and J. Barnes in Journal of Geodesy, Vol. 83, 2009, pp. 91–103 (doi: 10.1007/s00190-008-0236-7).

“LocataLites as a Solution to Open-cut Mining Applications” by J. Barnes in GPS World’s online TechTalk blog, posted February 21, 2008.

“High Accuracy Positioning Using Locata’s Next Generation Technology” by J. Barnes, C. Rizos, M. Kanli, A. Pahwa, D. Small, G. Voigt, N. Gambale, and J. Lamance in Proceedings of ION GNSS 2005, the 18th International Technical Meeting of the Satellite Division of The Institute of Navigation, Long Beach, California, September 13–16, 2005, pp. 2049–2056.

“A Positioning Technology for Classically Difficult GNSS Environments from Locata” by J. Barnes, C. Rizos, M. Kanli, and A. Pahwa in Proceedings of IEEE/ION PLANS 2006, the Position, Location, and Navigation Symposium, San Diego, California, April 25–27, 2006, pp. 715–721.

• Integrated Positioning

“Seamless Navigation Through GPS Outages – A Low-cost GPS/INS Solution” by Y. Li, P. Mumford, and C. Rizos in Inside GNSS, Vol. 3, No. 5, July/August 2008, pp. 39–45.

“Ubiquitous Positioning: Anyone, Anything: Anytime, Anywhere” by X. Meng, A. Dodson, T. Moore, and G. Roberts in GPS World, Vol. 18, No. 6, June 2007, pp. 60–65.

“Photogrammetry for Mobile Mapping: Bridging Degraded GPS/INS Performance in Urban Centers” by T. Hassan, C. Ellum, S. Nassar, W. Cheng, and N. El-Sheimy in GPS World, Vol. 18, No. 3, March 2007, pp. 44–48.

“Development of a GPS/INS Integrated System on the Field Programmable Gate Array Platform” by Y. Li, P. Mumford, J. Wang, and C. Rizos in Proceedings of ION GNSS 2006, the 19th International Technical Meeting of the Satellite Division of The Institute of Navigation, Fort Worth, Texas, September 26–30, 2006, pp. 2222–2231.

“An Integrated Positioning System: GPS + INS + Pseudolites” by Y. Yi, D. Grejner-Brzezinska, C. Toth, J. Wang, and C. Rizos in GPS World, Vol. 14, No. 7, July 2003, pp. 42–49.

• Kalman Filtering for Integrated Systems

“Tightly-coupled GPS/INS Integration Using Unscented Kalman Filter and Particle Filter” by Y. Yi and D.A. Grejner-Brzezinska in Proceedings of ION GNSS 2006, the 19th International Technical Meeting of the Satellite Division of The Institute of Navigation, Fort Worth, Texas, September 26–30, 2006, pp. 2182–2191.

“Low-cost Tightly Coupled GPS/INS Integration Based on a Nonlinear Kalman Filtering Design” by Y. Li, J. Wang, C. Rizos, P. Mumford, and W. Ding in Proceedings of NTM 2006, the National Technical Meeting of The Institute of Navigation, Monterey, California, January 18–20, 2006, pp. 958–966.

• Data Time Synchronization

“A Time-synchronisation Device for Tightly Coupled GPS/INS Integration” by P. Mumford, Y. Li, J. Wang, C. Rizos, and W. Ding in Proceedings of IGNSS Symposium 2006, International Global Navigation Satellite Systems Society, Gold Coast, Australia, July 17–21, 2006.

March 1, 2010