Estimating Spatial Variations in Atmospheric Delays using GPS
Thomas A. Herring, Massachusetts Institute of Technology, Cambridge, MA and
Seiichi Shimada, Solid Earth Science Division, National Research Institute of Earth Science Disaster Prevention, Tsukuba-Shi, Japan.
1. Introduction
In the analysis of Global Positioning System (GPS) data, it is often noticed that during conditions of disturbed atmospheric conditions, the carrier phase residuals increase. Often the station position and other geodetic parameters are also effected in these conditions. In this paper, we address the origin of these increases in the residuals. Are the variations due to inhomogeneous atmospheric propagation delays that can not be removed with the normal estimates of time varying spherically symmetric zenith delays? Or, are they due to changes in the electromagnetic (EM) environment near the GPS antenna?
We examine evidence for changes in the phase residuals during conditions of disturbed atmospheric conditions and we discuss methods that might be used to reduce the effects of atmospheric conditions on geodetic parameter estimates. We will consider several data sets that are summarized in Table 1. The choice of data sets is to show several characteristics of phase residuals and their changes under different weather conditions.
Table 1. Data sets analyzed.
|
Region |
Date |
Comments |
|
Central Asia |
Jan, 1998 |
Winter, stable weather conditions |
|
California |
Feb, 1998 |
Winter, rain on some days |
|
Japan |
Dec, 1997 |
Fall, calm conditions |
|
Japan |
Aug, 1997 |
Summer, normal day 213, Typhoon south of array days 217-218 |
2. Comparison of data sets characteristics
In our analysis of the GPS data in the data sets listed above we have computed postfit one-way phase residuals by explicitly estimating the satellite and station clocks at each 30 second epoch. Conditions are applied that make the average phase residual to each satellite, and at each site, zero at every epoch. The initial phase biases are also estimated based on pseudorange data and are refined in an iterative analysis with the estimation of the clock corrections. The estimation strategy assumes white-noise clocks. In Figure 1, we show examples of ionosphere-free phase residuals plotted as projections onto the sky. For the data sites in Central Asia and Japan the phase residuals are small. However, for the Californian site, HOLC, the residuals are large. The large residuals at HOLC appear to be due to a fiberglass box which houses the antenna [pictures of this and other sites can be found at http://www-socal.wr.usgs.gov/scign/dgga/log_sheet.html].
2.1 Environmental effects on phase residuals
The station HOLC provides an example of likely electromagnetic (EM) effects on the phase residuals. When it rains at this site the fiberglass box appears to change its EM characteristics increasing the phase residual RMS from typical values of 10 mm to nearly double this amount. Associated with these increases in residuals are changes in the estimates of the position of the site of up to 15 mm in latitude. Surprisingly, the effects on the height estimates are less than 10 mm. Figure 2 shows the change in residuals for one satellite (PRN 16) between a dry day (day 026) and a wet day (day 034). This satellite tracks from horizon-to-horizon. In the latter part of the track, the wet conditions appear to increase the size of oscillations on the residuals suggesting that these increases are due to the EM effects of the fiberglass box.
Figure 1
. Plots of GPS ionosphere-free phase residuals for selected sites. The phase residuals are plotted for a 4-hour interval and are shown as deviations from the tracks of the satellites across the sky. The sites shown are (a) SELE in Central Asia, 1998/029, RMS 4.1 mm; (b) 3048 in Japan (GSI station), 1997/347, RMS 4.4; and (c) HOLC in California, 1998/026, RMS 8.4 mm. The date is given as year and day-of-year. The darker shaded regions are positive phase residuals. The scales are such that a 190-mm phase residual would appear as a 19-mm deviation from the satellite trajectory. [Examples of these types of plots can be found at http://www-gpsg.mit.edu/~tah/cont98g/sky.html].
2.2 Atmospheric effects on residuals
Figure 2. Phase residuals at HOLC in California on dry (026) and wet (034) days. The residuals have been aligned to the same sidereal time.
In many cases changes to the phase residuals during unsettled conditions seem to reflect variations in the atmospheric delays rather than changes in the EM environment of the site. In Figure 3, we show the phase residuals for three days in August 1997, at the same time each day, for the Geographic Survey Institute (GSI) site 3048 on the Izu peninsular. On days 217 and 218, a typhoon was passing to the South of this region. In this case, the increase in the residuals is thought to be due to inhomogeneities in the atmospheric delays because the increases do not appear to be enhancements of pre-existing patterns in the residuals.
Figure 3.
Sky plots of the phase residuals at GSI site 3048 on days (a) 213, (b) 217, and (c) 218 for the interval 16-20 hrs UT. The RMS scatters of the residuals on the three days are 5.7, 7.9, and 11.1 mm, respectively. The scales are the same as in Figure 1.
All of the GPS sites in the Izu peninsular region of Japan show increased residuals during this time. In Table 2 we give the RMS scatter of the phase residuals on these days and on day 346 (December 1). The RMS scatter of the residuals increase by up to 14 mm on day 218 and show increases of up to 11 mm between the day in Fall and the summer day (213) before the typhoon.
3. Modeling the atmospheric delay
The one-way GPS phase residuals present an opportunity to map the atmospheric delays as a function of azimuth and elevation around a station. With sufficient density of sites it should by possible relate the atmospheric delays in different directions between sites at different locations. Coherence between atmospheric delay estimates from different locations derived from phase residuals to different satellites is one method that could be used to assess the reliability of the estimation procedures. We show a schematic model of this type in Figure 4. The assumptions in the model are that a thin layer at height H above the GPS sites can represent the atmospheric delay perturbations. We can then express the distance, D, to the perturbation as
(1)
where e is the elevation angle and we have assumed that D is sufficiently small that a flat-Earth approximation is valid. The North and East difference in position can be computed from the azimuth and distance. The effective atmospheric delay at this location is computed from the phase residual divided by the atmospheric delay mapping function (i.e., approximately 1/sin e in the flat-Earth approximation).
Figure 4.
Schematic model to represent inhomogeneities in the atmospheric propagation delay. The perturbations are assumed to be in a layer height, H, above surface. The position of the perturbation can be computed from the azimuth and elevation angle to a satellite and the assumed value of H. In this model, different GPS stations could see the same perturbations with measurements to different satellites (as shown in the schematic).However, to obtain reliable estimates of the atmospheric delay it is necessary to remove the effects of multipath and other non-atmospheric delay contributions to the phase residuals. In Figure 5 we show an example of large multipath that is largely eliminated by subtracting the phase residuals from a previous day shifted by 4 min-per-day. (This shift is possible because the GPS satellites are in orbits with periods very close to 12 sidereal hours and thus repeat at the same sidereal time each day). For finding the time offset which best can align longer separations of days, the time offset that brings the azimuth and elevation angles to the satellites in close agreement can be used. A consequence of removing the residuals from a previous day is that any perturbations in atmospheric delays on the reference day are transferred to the day being analyzed. When GPS data are collected continuously it should be possible to find the average multipath signal that can be removed from the data. The caveat to this approach is that the multipath may be effected by the local moisture conditions (see Figure 2).
Figure 5. Correction for multipath by removing the phase residuals from a previous day shifted by 4 min-per-day. The raw phase residuals on day 347 for PRN 17 at station KWN0 are shown with a line and symbol X; the corrected phase residuals using data from day 346 are shown with the open squares. The large oscillations between 8.5 hrs and 9.5 hrs UT in the raw phase residuals are almost certainly due multipath and not due to atmospheric delay perturbations.
4. Example application
In Figure 6, we show an example of the application of the model discussed in the previous section. We show results from day 218, 1997. Day 213 was used as the reference day to remove the effects of multipath. To generate these results, we shifted the epochs of the phase residuals on day 213 by 40 minutes to align the sidereal times of the phase measurements. The shifted residuals were subtracted from the corresponding residuals on day 218 and the changes in the zenith-delay estimates from the analysis that generated the phase residuals were added back into the residuals. The multipath-corrected residuals were converted to effective zenith delays by dividing by the atmospheric delay mapping function and their horizontal position determined assuming an effective atmospheric layer height of 5 km (see parts (a) and (b) of Figure 6). For comparison purposes we performed the same computation for the 347-346 pair of days in the fall of 1997. The effective zenith delays are shown part (c) of Figure 6.
Of significance in Figure 6 is the rapid increase in zenith delay seen with PRNs 9 and 10 seen between 8 and 9 hrs UTC. The peaks of these increased have been marked with numerals 2 and 3. The horizontal positions of these events (see Figure 6b) suggest that they were due to the motion of an area in increased atmospheric delay that propagated NE at a velocity of ~30 km/hr. If we trace this motion back in time, we find that measurements to PRN 24 at site 0297, approximately 20 km SW of 3046, made around 8 hrs UTC show a similar rapid increase in zenith delay. This event is labeled 1 in Figure 6. These results suggest that we were able to track an anomalous variation in zenith delay with a spatial scale of only a few kilometers over a time period of almost an hour. Analysis of the tracks of other satellites from other stations in the regions suggests that none of the other stations were able to see satellites in correct direction and at the correct time to detect this feature. There is some indication at site 3085, 37-km SW of 3048 may have seen the edge of the feature at 7.8 hrs UTC but in this case the increase in zenith delay was not as pronounced as it was at sites 3048 and 0297.
5. Conclusions
The effects of changes in weather conditions can be readily seen in GPS carrier phase residuals. In this paper we have examined several examples of such changes. In some cases, the changes seem to be due to variations in the electromagnetic environment around the GPS antenna rather than changes in atmospheric propagation delays. In other cases, the changes do appear to be due to spatial inhomogeneities in the atmospheric propagation delay. In one example shown, it appears that we were able to track an anomalous atmospheric region between two stations and over an interval of time. We showed also that for correct interpretation of atmospheric delay variations, multipath needed to be accounted for at individual stations. In the method outlined here, correcting for multipath results in our zenith delay estimates being relative measures of changes in zenith delay between pairs of days. With sufficient duration of data available, it would be possible to determine average the multipath values provided these values did not depend too much on local environmental conditions. The techniques used here could lead to more refined models of the variations of a spatially inhomogeneous atmosphere that ultimately could be used to improve the accuracy of geodetic parameters estimates.
Figure 6
. Example of conversion of one-way phase residuals to effective spatially distributed zenith delay estimates. In (a), we show phase residuals mapped to effective zenith delays for a four-hour period at GSI site 3048, except for the line labeled "PRN 24 0297" which is from GSI site 0297 located 8.2 km South and 19.1 West of 3048. The data are from days 218 and 213, 1997. We have marked the peaks of three rapid increases in the zenith delays at times 8.0, 8.4 and 8.7 hrs UTC, labeled with the numerals 1, 2 and 3 respectively. In (b), we show the effective horizontal positions of the zenith delay estimates computed assuming H=5 km (see Figure 4). For PRN 24 at site 0297 we show the position of the track relative to site 3048. The numerals (1, 2, and 3) mark the horizontal positions of the peaks in (a). For comparison, in (c) we show a similar analysis except the days used were 347 and 346 from the fall of 1997. In this latter case, the variations between the satellites are much smaller than those variations seen on days 218-213.Table 2
. Root-mean-square (RMS) scatter of postfit phase residuals for sites in the NIED and GSI networks. The columns labeled DRMS give the change in RMS between day 213 and the day number given. For day 346, the difference is in the sense of RMS on day 213 minus day 346.|
Sites |
RMS 213 (mm) |
RMS 217 (mm) |
RMS 218 (mm) |
RMS 346 (mm) |
D RMS 217 (mm) |
D RMS 218 (mm) |
D RMS 346 (mm) |
|
0229 |
7.2 |
11.0 |
10.9 |
4.9 |
8.3 |
8.2 |
5.3 |
|
0230 |
7.9 |
11.6 |
14.0 |
7.4 |
8.5 |
11.6 |
2.8 |
|
0297 |
7.1 |
11.4 |
13.2 |
5.5 |
8.9 |
11.1 |
4.5 |
|
2106 |
8.2 |
11.3 |
12.2 |
7.3 |
7.8 |
9.0 |
3.7 |
|
2107 |
6.4 |
9.6 |
11.1 |
4.9 |
7.2 |
9.1 |
4.1 |
|
2108 |
6.3 |
9.1 |
10.7 |
4.9 |
6.6 |
8.6 |
4.0 |
|
2109 |
6.4 |
8.3 |
10.0 |
5.6 |
5.3 |
7.7 |
3.1 |
|
3031 |
7.6 |
11.6 |
10.6 |
5.6 |
8.8 |
7.4 |
5.1 |
|
3035 |
8.5 |
12.2 |
12.8 |
7.1 |
8.8 |
9.6 |
4.7 |
|
3042 |
8.1 |
12.7 |
15.2 |
5.4 |
9.8 |
12.9 |
6.0 |
|
3048 |
5.7 |
7.9 |
11.1 |
4.4 |
5.5 |
9.5 |
3.6 |
|
3058 |
6.7 |
10.1 |
9.8 |
4.7 |
7.6 |
7.2 |
4.8 |
|
3066 |
8.2 |
10.3 |
12.0 |
5.6 |
6.2 |
8.8 |
6.0 |
|
3068 |
8.9 |
13.1 |
12.8 |
6.0 |
9.6 |
9.2 |
6.6 |
|
3077 |
7.0 |
9.1 |
9.6 |
6.6 |
5.8 |
6.6 |
2.3 |
|
3080 |
5.9 |
9.0 |
10.6 |
5.1 |
6.8 |
8.8 |
3.0 |
|
3081 |
6.7 |
9.7 |
10.9 |
6.0 |
7.0 |
8.6 |
3.0 |
|
3085 |
7.1 |
9.8 |
10.4 |
5.3 |
6.8 |
7.6 |
4.7 |
|
3087 |
6.0 |
8.8 |
10.8 |
5.4 |
6.4 |
9.0 |
2.6 |
|
3088 |
6.4 |
11.2 |
13.7 |
5.6 |
9.2 |
12.1 |
3.1 |
|
3092 |
6.1 |
8.9 |
9.7 |
5.5 |
6.5 |
7.5 |
2.6 |
|
3094 |
6.2 |
9.3 |
10.6 |
5.9 |
6.9 |
8.6 |
1.9 |
|
3101 |
5.6 |
8.2 |
9.7 |
5.2 |
6.0 |
7.9 |
2.1 |
|
4111 |
7.0 |
9.6 |
11.3 |
5.6 |
6.6 |
8.9 |
4.2 |
|
5105 |
6.1 |
9.3 |
11.5 |
4.6 |
7.0 |
9.7 |
4.0 |
|
HMO0 |
14.9 |
16.8 |
17.1 |
12.5 |
7.8 |
8.4 |
8.1 |
|
HRT0 |
11.2 |
13.2 |
14.6 |
0.0 |
7.0 |
9.4 |
11.2 |
|
HTK0 |
10.0 |
12.2 |
14.0 |
7.7 |
7.0 |
9.8 |
6.4 |
|
HTS0 |
12.4 |
14.7 |
16.3 |
12.0 |
7.9 |
10.6 |
3.1 |
|
KRK0 |
9.7 |
11.1 |
10.7 |
8.5 |
5.4 |
4.5 |
4.7 |
|
KWN0 |
11.0 |
12.9 |
12.5 |
10.5 |
6.7 |
5.9 |
3.3 |
|
NBK0 |
10.6 |
13.9 |
17.4 |
8.8 |
9.0 |
13.8 |
5.9 |
|
NRY0 |
11.4 |
14.2 |
15.9 |
10.2 |
8.5 |
11.1 |
5.1 |
|
OYM0 |
8.9 |
13.1 |
13.9 |
7.5 |
9.6 |
10.7 |
4.8 |
|
SIB0 |
8.7 |
11.0 |
11.6 |
7.6 |
6.7 |
7.7 |
4.2 |
|
SMD0 |
12.3 |
14.4 |
14.6 |
11.5 |
7.5 |
7.9 |
4.4 |
|
TSKB |
6.5 |
10.9 |
11.5 |
6.0 |
8.7 |
9.5 |
2.5 |
|
USUD |
10.3 |
15.2 |
10.4 |
5.8 |
11.2 |
1.4 |
8.5 |
|
ALL |
8.4 |
11.4 |
12.5 |
7.1 |
7.7 |
9.3 |
4.5 |
Acknowledgements:
This works was supported by through the Special Coordination Funds of the Science and Technology Agency of the Japanese Government. We also thank Mr. Yuki Hatanaka for the GSI data available for this study. This work was also supported by NASA Grant NAG5-3550.