With a given GPS receiver you get a certain positioning accuracy. So, the first question you may ask is: can I do anything to improve it?
The answer most often is: yes, go and buy a better receiver! Most ranging errors are determined by physics and you can do little to improve the situation.
With a given GPS receiver you get a certain positioning accuracy. So, the first question you may ask is: can I do anything to improve it?
The answer most often is: yes, go and buy a better receiver! Most ranging errors are determined by physics and you can do little to improve the situation.
However, it is interesting to analyze the information you may get from a single “one-way” range, that is, a single set of observables (P1, Φ1, P2, Φ2) related to an individual satellite seen by a single dual-frequency receiver. By so doing you also obtain an insight into the nature of the error sources — the ones that do matter and the ones that can be modeled or eliminated.
From the ephemeris of the tracked satellite we can compute the elevation angle, the ionospheric delay is estimable because the data originate from a dual-frequency receiver, the tropospheric delay is computed from a standard model, and the remaining errors we call multipath errors.
Finally, we estimate the orbital errors by comparing satellite positions computed from broadcast and precise ephemerides. The data sample illustrates these different error contributions very well, their variation over time, and, therefore, their dependency on elevation angle.
We select a specific one-way observation set that was acquired at a station south of Aalborg on January 19, 2007. (This is from the data set already used in easy15, published in the January/February 2010 issue of Inside GNSS.) We investigated PRN14 in the full session length, that is, about 42 minutes — including the 95 epochs that we omitted in the easy15 discussion because of missing data at either the rover or the master, or repeated estimation of ambiguities.
The subplots in Figure 1 depict the satellite elevation angle, which varies from 18 to 33 degrees, and one-way errors in ionosphere, troposphere, and multipath.
When dealing with the ionospheric delay I, we assume that we have both pseudorange and carrier phase observations on L1 and L2 at our disposal. We estimate the delay I according to the following matrix equation:
As usual we define the constant α = (f1 / f2)2 = I2 / I1. First we estimate the ambiguities on L1 and L2 — N1 and N2, respectively — as reals; incorrect values for N1 and N2 just change the level for I. Next we estimate I alone as
For the tropospheric delay we use the Goad-Goodman model. The delay T ranges from 7.4 m to 4.2 meters. If the satellite passed zenith, the corresponding delay would have been about 2.5 meters.
Multipath describes the situation where signals coming from the satellite propagate along several paths to the receiver antenna. The main part of the signal radiates directly from the satellite, but part of the signal is reflected from surfaces near the receiver.
Multipath occurs when the signal arrives at the antenna from these reflected surfaces in addition to the line-of-sight source. The reflected signal is phase-shifted with respect to the original transmission and appears as additive noise at the antenna.
Multipath depends on satellite geometry and the antenna environment, which makes multipath difficult to model. For long observation periods — 24 hours or more — the multipath effects are partly reduced by averaging. However, observation periods usually last for only a few hours and often much less; this is why multipath is a problem.
Because the antenna locations are different, the multipath signature at each antenna is unique and the error is not common mode.
Let the geometrical distance between satellite and receiver be ρ, let I be the ionospheric delay, and M the code multipath including receiver noise. The pseudoranges observed on L1 and L2 for PRN14 may then be expressed as
Remaining errors are identified as multipath and receiver noise — similarly for the phase observations:
The multipaths m1 and m2 on phase observations are so small that we subsequently put mi = 0. We want to find an expression for M1. We start by subtracting (5) from (3):
and subtracting (6) from (5) yields
We insert (9) into (7) and obtain
We can reasonably assume that E{M1}= 0. The second term on the left side of (11) is a constant; so, it is possible to reduce M1 to M1* such that E{M1*] = 0:
Analogously, by exchanging subscripts, we have for the multipath on L2:
For all epochs with all four observations P1, P2, Φ1, and Φ2 available, we compute multipath according to Equation (12). In the present case, the average error is below two meters and noisier at low satellite elevation angles. The noisy character of the plot reflects the actual noise in any pseudorange observation.
Finally we investigate the difference between precise and broadcast ephemerides. Any file containing broadcast ephemerides is specific to the site from which the satellites were tracked. Therefore, we need to identify the tracked satellites in the time period during which we want to compare the precise and broadcast ephemerides.
This can be done by inspecting the corresponding observation file in the PRN14 navigation message and getting the date and the seconds of the week; next we find the corresponding Julian Day Number (jd). We find an epoch (modulo 15 minutes) 60 minutes ahead of jd, which is counted in unit of day. (The variable jd is a real number where the integer part is number of whole days and the decimal part is made up of hours, minutes, seconds and fractions of seconds all convert into fractions of a day). Finally the GPS week number is computed by the Matlab function gps_time.
The precise ephemerides that we used can be found at igscb.jpl.nasa.gov/igscb/product/1410, the number 1410 being the GPS week number. The file extracts to igr14105.sp3, which contains precise ephemerides for January 19, 2007, in the SP3 format. The IGS SP3 files typically contain satellite positions (X,Y,Z) and clock offsets for each satellite every 15 minutes, from 0:00 hours to 23:45. (The SP3 file starts with some header lines, which for our purposes we can skip.)
We read precise orbits for the interpolation period plus three times 15 minutes ahead of the observation period and three times 15 minutes after the period, i.e., 13:15 hours, 13:30 hours, … , 15:15 hours — in total nine sets. The precise coordinates for all tracked satellites are stored in the matrix Xp.
In the following discussion, we introduce a time scale with units of 15 minutes. We could have chosen whole minutes instead. The first observation is from 13.85 hours, and the last observation is taken at 14.55 hours. The entire observation period is 2,535 seconds = 42.25 minutes = 2.8 quarters of an hour.
A Lagrange interpolation, one for each coordinate, of at least seventh order is used to compute the actual position. The Matlab function intp does an (n – 1)th order polynomial Lagrange interpolation. Let y be an n × 1vector of data given at the discrete times x, (x is as well an n × 1vector), (See the relevant discussion in Matrix Analysis by R. A. Horn and C. R. Johnson, Volume 1, pages 29–30, cited in Additional Resources at the end of this article.)
We want an interpolated precise position each minute; hence, we divide by 15 to keep units in quarters of an hour. The first parameter in the intp procedure describes the abscissae for the points at which we know the satellite coordinates. The second parameter contains these coordinates, and the third parameter describes the point set at which we want interpolated values — all in units of quarter of an hour. Each coordinate is interpolated separately!
Figure 2 shows the difference in satellite position as computed from broadcast ephemerides given via the navigation file, and the postprocessed satellite positions for PRN14. When comparing with broadcast ephemerides we assume the precise positions to be the true ones.
The actual influence of the orbit error on the receiver position is given by the projection of the difference vector onto the line between the receiver and satellite. Let vector b be the difference vector between the satellite positions computed from the broadcast and the precise ephemerides, and a to be the vector between the receiver and the satellite position as
computed from broadcast ephemerides. Then the projected difference vector p onto vector a is given as
Figure 3 shows vector p for PRN14. The computed ephemeris error varies in the range of ±2 meters. The discontinuities in the error reflect the routine ephemeris updates at two-hour intervals.
Additional Resources
Horn, R. A., and C. R. Johnson, Matrix Analysis.
Cambridge University Press, Cambridge, 1985
