A semi-analytical estimation of the effect of second-order ionospheric correction on the GPS positioning

doi: 10.1111/j.1365-246X.2005.02723.x A semi-analytical estimation of the effect of second-order ionospheric correction on the GPS positioning H. Munekane Geographical Survey Institute, 1 Kitasato, Tsukuba, Ibaraki 3050811, Japan. E-mail: munekane@gsi.go.jp Accepted 2005 June 27. Received 2005 June 15; in original form 2005 March 22 SUMMARY We developed a semi-analytical method to evaluate the effect of the second-order ionospheric correction on GPS positioning. This method is based on the semi-analytical positioning error simulation method developed by Geiger and Santerre in which, assuming the continuous distribution of the satellites, a normal equation is formed to estimate the positioning error taking all the contributions of the ranging error by the visible satellites into account. Our method successfully reproduced the averaged time-series of three IGS sites which is comparable to the rigorous simulation. We then evaluated the effect of the ionospheric error on the determination of the reference frame. We evaluated the additional Helmert parameters that are required for the ionospheric effect. We found that the ionospheric effect can lead to annual scale changes of 0.1 ppb, with an offset of 1.8 mm and a semi-annual oscillation of 1 mm in the z-direction. However, these values are too small to explain the current deviations between the GPS-derived reference frame and the ITRF reference frame. Next, we estimated the apparent scale changes due to the ionospheric error in the GEONET coordinate time-series in Japan. We could qualitatively reproduce the observed semi-annual scale changes peaking at the equinoxes and having asymmetrical amplitudes between the vernal and autumnal equinoxes. Key words: geodesy, GPS, ionosphere. 1 I N T RO D U C T I O N The GPS system is now being widely used for a variety of studies, such as crustal deformation (e.g. Miyazaki et al. 2004), local site deformation (e.g. Munekane et al. 2004) or loading deformation (e.g. Munekane & Matsuzaka 2004) among others. The progress in these studies has led to an increasing demand for more precise coordinate time-series. Correspondingly, many studies have been devoted to the evaluation and elimination of various sources of noise in GPS analysis such as tropospheric delay (Cucurull et al. 2002) or antenna phase centre and multipath (e.g. Park et al. 2004). Recently, among such noise sources, the effect of second-order ionospheric correction was investigated by Kedar et al. (2003). In the ordinary GPS dual-frequency analysis, the propagation group delay (phase advance) effect in the observed ranges due to the ionosphere is treated as inversely proportional to the square of the frequency of the propagating radiowaves (first order), and it is eliminated by taking a linear combination of the observed ranges in both frequencies. The study of Kedar et al. (2003) considered the second-order correction, which is inversely proportional to the cube of the frequency. They compared the coordinate time-series obtained by analysing the ranges that were corrected for the second-order ionospheric effect with that obtained by analysing the original ranges, and evaluated the impact of this correction. They found that the second-order effect is responsible for the positioning error at the subcentimetre level, 10 which is well above the accuracy limit specified in the current GPS analysis. Here, we present a simpler method in order to evaluate the effect of the second-order ionospheric correction on the coordinate time-series. This method is similar to the error estimation method proposed by Geiger (1988) and Santerre (1991). In this method, the observation equation which relates the ranging error from a satellite at a certain zenithal angle and azimuth to the positioning error is considered. Assuming that the satellites uniformly occupy a certain area in the observer’s sky, a normal equation is formed from the observation equation to estimate the positioning error taking all the contributions of the ranging errors by the visible satellites into account. The advantage of this method is that it provides us with better physical insights on how the ranging error affects the position estimates than can be obtained from a purely numerical evaluation. In this paper, we first reveal the general characteristics of the ionosphere-related positioning error for several analysing strategies using the semi-analytical evaluation method. Then, we specifically investigate the ionospheric effect in the Japanese GPS network. 2 SECOND-ORDER IONOSPHERIC CORRECTION In this section, we briefly summarize the characteristics of the second-order ionospheric correction. Readers should refer to
C 2005 RAS Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 GJI Geodesy, potential field and applied geophysics Geophys. J. Int. (2005) 163, 10–17 A semi-analytical estimation of the second-order ionospheric effect Brunner & Gu (1991), Bassiri & Hajj (1993) and Kedar et al. (2003) for detailed description. The ranges observed by the GPS are expressed as follows (eqs 1 and 2 in Kedar et al. 2003): s q Pi = ρ + 2 + 3 , (1) fi fi L i = ρ + n i λi − 1 s q − , 2 2 f i3 fi (2) where the index i denotes the GPS carrier frequency ( f 1 = 1.6 GHz and f 2 = 1.2 GHz); Pi and L i refer to the observed code and phase ranges, respectively; ρ represents the geometrical range between a satellite and a receiver; n i and λ i denote the integer ambiguity and wavelength, respectively. The coefficients q and s are given by s = 7527c
Nk · B dL, (3) (4) where N denotes the number density of electrons, and TEC denotes the total electron content along the propagation path of the GPS signal. c denotes the speed of light and k is the unit vector representing the direction of propagation of the GPS signal. B represents the magnetic field of the earth. The second and third terms in the right-hand side of eqs (1) and (2), respectively, correspond to the first-order ionospheric ranging error for the code and phase ranges while the third and fourth terms in the right-hand side of eqs (1) and (2), respectively, correspond to the second-order ionospheric ranging error for the code and phase ranges. The thin-shell approximation (e.g. Mannucchi et al. 1998) assumes that the ionosphere is localized in a ‘thin shell’. Considering this, we can take the dot product of k and B out of the integrand and eq. (4) is modified to s = 7527c(k · B)T EC. (5) Note that the coefficient for the second-order correction s is proportional to both TEC and the dot product of k and B. For the dual-frequency GPS analysis, we usually take a linear combination of eq. (2) in order to form the ionosphere-free combination, where the first-order ionospheric terms are eliminated (e.g. Hofmann-Wellenhof et al. 2001). The ionosphere-free combination L c is given by Lc = f 12 f2 f 12 L − 2 2 2 L2 2 1 − f2 f1 − f2 = ρ + d Lc, (6) where dLC , the ionospheric ranging error in the ionosphere-free combination, is expressed as   1 1 1 1 +  s.   d Lc = − (7) 2 f 1 f 12 − f 22 2 f 2 f 12 − f 22 3 E S T I M AT I O N O F P O S I T I O N I N G E R RO R 3.1 Positioning error when satellite positions are precisely known We first examine the positioning error when the positions of the satellites are kept fixed during the analysis, which is usually done when analysing a small GPS network.
C 2005 RAS, GJI, 163, 10–17 We estimated the positioning error by applying Geiger’s (1988) error estimation method in the case of single point positioning case. Henceforth, we assume that the satellite positions are precisely given. Let us consider the coordinate system whose axes point towards the east, north and the vertical direction, respectively. dLc , the ranging error for a satellite, is related to the positioning error as follows: d L c = Cdr, (8) where C is a design matrix and dr is the positioning error of the station. The components of C are given as C = (e1 , e2 , e3 ), (9) e1 = −sinθsinλ, (10) e2 = −sinθcosλ, (11) e3 = −cosθ. (12) Here, θ is the zenithal angle and λ denotes the azimuth angle measured clockwise from north. The observer’s clock error is not taken into account since it is effectively removed when taking the double difference of the observed ranges in the ordinary GPS analysis. We assumed that the distribution of the satellites is homogeneous within their area of visibility. This is a reasonable assumption since the current constellation of 29 satellites is considered. We integrate eq. (8) to form a normal equation, which is given by U = N dr. (13) The normal matrix N is expressed as   [e1 e1 ] [e1 e2 ] [e1 e3 ] N =  [e2 e1 ] [e2 e2 ] [e2 e3 ]  , [e3 e1 ] [e3 e2 ] [e3 e3 ] (14) and the vector U is given by   [e1 d L c ] U =  [e2 d L c ]  , [e3 d L c ] (15) where the bracket [f ] denotes the integration of the function f over the area where the satellites are visible. According to Santerre (1991), this area may be represented as in Fig. 1. In this case, the bracket is expressed as
λmax
θmax [f] = f sinθ dθ dλ, (16) λmin θmin where λ min and λ max are the azimuthal integration boundaries and θ min and θ max are the zenithal integration boundaries. To estimate dr using the normal equation, eq. (13), we need to evaluate the matrices N −1 and U. The calculation of matrix N −1 , or N, is straightforward given the definition in eqs (10)–(14). On the other hand, the calculation of U requires the knowledge of dLc for a certain zenithal angle and azimuth, which in turn requires the evaluation of the earth’s magnetic field B and the total electron content TEC at the ionosphere pierce point (IPP), where a radiowave from a satellite intersects the thin-shell ionosphere. The location of the IPP may be estimated as follows. We treat the earth as a sphere with radius Re and the ionosphere as a thin shell located at a distance h ion above the ground. Let b sta and l sta represent the latitude and longitude of the station, respectively, and θ and λ denote the zenithal angle and the azimuth of a satellite, respectively. The position vector of the IPP may be expressed using the position Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 q = 40.3 × T EC, 11 12 H. Munekane Table 1. Parameters to represent the satellite visibility used in this study. Station latitude −60◦ −20◦ < b sta < −20◦ < b sta < 20◦ 20◦ < b sta < 60◦ λ min λ max θ min θ max −135◦ 135◦ 0◦ 0◦ 45◦ 360◦ 315◦ 0◦ 0◦ 83◦ 83◦ 83◦ Figure 1. The ‘sky plot’ in the local coordinate system. The hatched area shows the area where the satellites are visible. The integration boundaries are represented by four parameters, that is, λ min and λ max for the azimuth and θ min and θ max for the zenithal angle. vector of the station and the distance between the IPP and the station as rIPP = rsta − r ′ k, (17) k = −sinθ cosλen −sinθ sinλee (18) −cosθeu , where r IPP , r sta and r ′ denote the position vector of the IPP, the position vector of the station, and the distance between the IPP and the station, respectively. k is the unit vector along the propagation direction of the radiowave, as previously defined in Section 2, and e e , e n and e u represent the unit vectors eastwards, northwards and upwards from the station, respectively. The distance r ′ is expressed as A = θ − sin−1 r′ = Re sinθ , Re + h ion  Re2 + (Re + h ion )2 − 2Re (Re + h ion ) cos A . (19) (20) We then evaluated the magnetic field at the IPP using the IGRF2000 geomagnetic model (Mandea et al. 2000). With respect to the TEC, we first evaluated the vertical total electron content (VTEC) using the Global Ionosphere Map (GIM) developed at the Centre for Orbit Determination in Europe (CODE). The VTEC is then mapped to the TEC using a regular mapping function as follows T EC = 1 V T EC, cosθ ′ (21) Re sinθ. Re + h (22) where sinθ ′ = The shape of the void coming from the GPS constellation is circular with its centre at the celestial pole. Therefore, the position of this void as seen from a ground station depends on the latitude. We selected the four parameters representing the satellite visibility as in Table 1 to take into account this latitudinal dependence of the void. 3.2 Positioning error when the satellite positions contain the ionospheric error Next, we consider the case where the satellite positions have errors due to the ionospheric effect. First, let us estimate the satellite position error due to the ionospheric effect. We assume that the satellite positions are estimated using the fiducial sites on the ground. In this case, the satellite position error may be estimated in the same manner as the positioning error of the ground stations. Let us also suppose that the GPS satellites are located above the ground at a height of h sat . We define the azimuth α and the zenithal angle δ as shown in Fig. 3, similar to the ground station case. The azimuth χ and ‘colatitude’ φ of the ground station are defined in the coordinate system such that its z-axis coincides with the vertical axis at the satellite, and the origin of the azimuth is parallel to the north axis at the satellite, as shown in Fig. 4. The ranging error dLc is related to the satellite position error r sat in the following equation: d L c = Csat drsat , (23) where the components of C sat are expressed as   Csat = e1sat , e2sat , e3sat , (24) e1sat = sinδsinα, (25) e2sat = sinδcosα, (26) e3sat = cosδ. (27) Assuming that the distribution of the ground stations is uniform, we integrate eq. (23) to obtain a normal equation as follows: Usat = Nsat drsat . (28) The normal matrix N sat is given by   {e1 e1 } {e1 e2 } {e1 e3 } Nsat =  {e2 e1 } {e2 e2 } {e2 e3 }  , {e3 e1 } {e3 e2 } {e3 e3 } (29) and the vector U sat is given by   {e1 d L c } Usat =  {e2 d L c }  , {e3 d L c } (30)
C 2005 RAS, GJI, 163, 10–17 Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 Fig. 2 shows the north component of the estimated station position error on 2000 April 2. It can be observed that our method reproduces a large northward error (>5 mm) near the magnetic equator that appears as large southward corrections in fig. 4 in the study of Kedar et al. (2003). In addition, one can observe that the large error in the subequatorial region around North India and China, reflecting the large north component anomaly of the Earth’s magnetic field around this area, which is apparent in the IGRF model. A semi-analytical estimation of the second-order ionospheric effect 13 Figure 2. North component of position error for precise satellite positions. Calculated on 2000 April 2. where the bracket f denotes the integration of the function f over the area on the earth from which the satellites are visible. The bracket is expressed as
2π
φmax {f} = f sinφ dφ dχ , (31) 0 0 where φ max represents the limit ‘colatitude’ of the stations from which the satellites are visible, and can be expressed using the maximum zenithal angle θ max as φmax = θmax − sin−1 Re sinθmax . Re + h sat (32) (33) and r ′′ = Re2 + (Re + h sat )2 − 2Re (Re + h sat )cosφ, δ = sin−1
C Re sinφ , r ′′ 2005 RAS, GJI, 163, 10–17 where r ′′ denotes the distance between the satellite and the station. The location of the IPP is expressed as rIPP = rsat + r ′′′ k (36) k = − sinδsinαesat e − sinδcosαesat n In order to complete the integration in N sat and U sat , we need to express the angles α and δ with χ and φ. These relationships may be expressed as α = χ − π, Figure 4. Definition of the azimuth and ‘colatitude’ of the ground station. Both the angles are defined in the coordinate system whose z-axis shares the vertical axis at the satellite, and the origin of the azimuth is parallel to the north axis of the satellite. (34) (35) − cosδesat u , (37) sat sat where esat e , en and eu represent the unit vectors directed eastwards, northwards and upwards from the satellite, respectively, and r ′′′ denotes the distance between the satellite and the IPP. r ′′′ is expressed as  r ′′′ = (Re + h sat )2 + (Re + h ion )2 − 2(Re + h sat )(Re + h ion ) cos φ . (38) The ionospheric mapping function similar to that in eqs (21) and (22) is used to convert the VTEC to the TEC. The zenithal angle at the Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 Figure 3. Definition of the azimuth and the zenithal angle of the radiowave propagating from a satellite. The azimuth is measured clockwise from the north direction. 14 H. Munekane Figure 5. North component of the estimated station position error when a satellite position error exists. The colour scale is narrower than that of Fig. 2. The error observed is much smaller than that in Fig. 2. Re + h sat sinφ . (39) r ′′′ If a satellite position error exists, the range error for the given satellite is modified to θ = sin−1 d L c = d L fix c + drsat · k, (40) where dLfix c is the range error when the exact satellite positions are known and dr sat denotes the satellite position error due to the second-order correction. Therefore, the station position error under the existence of the satellite position error may be obtained by replacing eq. (7) with eq. (40) when constructing the normal equation, eq. (13). The other conditions in the numerical evaluations are the same as those in the previous case. In this case, the north component of the estimated station position error is shown in Fig. 5. It can be seen that the large positioning error seen in Fig. 2 is significantly reduced. This is because the ranging error is mostly absorbed in the satellite positioning error, and dLc and dr sat · k in eq. (40) tend to cancel each other. 3.3 Fiducial-free approach In the practical GPS analysis, we often estimate the satellite orbits and station positions simultaneously to produce station coordinates in the centre of mass of the Earth system (CM) frame, and then transform the station positions in the CM frame into other frames such as the centre of figure (CF) frame via the Helmert transformation (e.g. Hofmann-Wellenhof et al. 2001). In this case, the station position error due to the ionospheric correction in the CM frame is the same as that in the first case. However, due to the ionospheric error, an additional Helmert transformation is required to align the station positions in other frames. This additional Helmert transformation is expressed as dr = Hadd r, (41) where H add denotes the matrix corresponding to the additional Helmert transformation; dr and r represent the ionospheric positioning error and station positions, respectively. We used virtual stations placed at every 5◦ from −60◦ to 60◦ in latitude and 0◦ to 360◦ in longitude and determined the parameters for the additional Helmert transformation (scale, shifts and rotations) in the least-squares method. First, let us take a look at the time-series of transformation parameters. The transformation parameters obtained from 2000 to 2003 are shown in Fig. 6. It can be seen that the ionospheric position Figure 6. The time-series for the estimated parameters of the Helmert transformation. The three rotation parameters are omitted since they are negligible. error mainly affects the scale and the shift in the z-direction. The scale changes exhibit an annual oscillation whose amplitude is approximately 0.1 ppb, and the network shifts show an offset of 1.8 mm and a semi-annual oscillation with an amplitude of 1 mm in the z-direction. Fig. 7 shows the station position error after the Helmert transformation was carried out. As compared with Fig. 2, it can be seen that the position error is slightly reduced since some portion of the error is absorbed in the transformation parameters. However, a position error exceeding 5 mm is still expected in the subequatorial region around North India and China.
C 2005 RAS, GJI, 163, 10–17 Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 station (θ) is essential for the calculation of the mapping function. It is expressed in terms of φ as A semi-analytical estimation of the second-order ionospheric effect 15 Figure 7. North component of the estimated station position error after the Helmert transformation. Calculated on 2000 April 2. In order to test the precision of the semi-analytic simulation that was developed here, we calculated the average time-series at three IGS sites (BAHR, COCO and GALA). This time-series is comparable to that shown in fig. 3 in the study of Kedar et al. (2003), which was estimated by direct simulation by comparing the estimated positions with the corrected/original GPS data. Fig. 8 shows the north component of the position error averaged for the three sites. The station position time-series estimated by the GIPSY-OASIS PPP analysis (Zumberge et al. 1997) is also shown (available at ftp://sideshow.jpl.nasa.gov/mbh/point). It should be noted that a remarkable correlation is observed between the two time-series. 4 A P P L I C AT I O N : I N V E S T I G AT I N G A N O M A LO U S S E A S O NA L S C A L E C H A N G E S I N T H E G E O N E T, J A PA N The GPS earth observation network (GEONET) is a dense GPS array deployed over the Japanese islands by the Geographical Survey Institute (GSI). It is well known that the coordinate time-series from the GEONET stations exhibit large systematic seasonal variations, and several studies have been devoted to determining the origins. Heki (2001, 2004) investigated the seasonal variations using the coordinate time-series in detail from the solutions of routine analyses and concluded that most of the variations originated from the snow load, excluding the anomalous annual scale changes of unknown origin whose amplitude measures 5–6 ppb (Hatanaka 2003). We reanalysed the scale changes from the GEONET stations. We first selected 77 stations that are evenly sampled in the GEONET network. Fig. 9 shows the stations that were used in the analysis. We estimated the coordinate time-series daily for the selected stations using the GIPSY-OASIS software with the PPP strategy (Zumberge
C 2005 RAS, GJI, 163, 10–17 Figure 9. GPS stations used to estimate the transformation parameters in the GEONET. et al. 1997). The network Helmert transformation for the network is defined as dr = H r, (42) where dr denotes the station coordinate variations after detrending, and r denotes the average station positions. The transformation parameters (scale, shifts and rotations) are estimated using the leastsquares method. We also estimated the transformation parameters for the time-series of the simulated positioning error in the same manner. Fig. 10 shows the scale changes estimated from the GPS analysis and the simulation. It can be observed that the scale changes obtained from the GPS analysis exhibit a semi-annual oscillation rather than an annual one, which were reported in the previous analyses (Heki 2001; Hatanaka 2003). The annual oscillation reported in the previous analyses is considered to be largely artificial due to a bug in the GPS analysis software which was used in the routine analysis (Hatanaka 2004). Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 Figure 8. The north component of the position error estimates averaged over three IGS sites (BAHR, COCO and GALA). Red dots represent the observed time-series of the site position and blue dots represent the time-series of the simulated position error. It should be noted that a good correlation is observed between the time-series of the position error and the position estimates obtained from the observed GPS data. 16 H. Munekane The peaks exhibited by the scale changes are maximum at the equinoxes and minimum at the solstices. It can also be observed that the peak amplitudes around the vernal equinoxes are slightly larger than those around the autumnal equinoxes. The simulated scale changes qualitatively reproduce the characteristics of the scale changes that are estimated from the GPS analysis, namely, the semi-annual oscillation and asymmetry in the peak amplitudes, accurately. 5 DISCUSSION It is interesting to note that while an annual oscillation is prominent in the global scale time-series, the semi-annual oscillation is prominent in the GEONET scale time-series. This aspect is explained as follows. At both the equinoxes, the northward shift peaks due to the peaking of the TEC values (e.g. Essex 1977). Therefore, an excess northward shift will be observed at each station after detrending the timeseries. The magnitude of the excess northward shift is greater if the stations are near the magnetic equator where the TEC values are higher. Therefore, after removing the network shifts, a north south contraction appears in the GEONET, which is observed as the scale contraction (Fig. 11). On the contrary, the excess southward shift is observed at both the solstices after detrending the time-series. Therefore, the north south expansion appears after removing the network shifts. At a global level, the distribution of the northward shift is symmetric across the magnetic equator at both the equinoxes since the TEC values are symmetric across the magnetic equator. As a result, the scale changes will not appear (Fig. 11) after eliminating the network shifts. On the other hand, at both the solstices, the northward shift is asymmetric across the magnetic equator, reflecting the asymmetric TEC distribution. This asymmetric distribution is interpreted as scale changes that have an opposite sign between solstices. We noted that the observed scale changes are larger at the vernal equinoxes than the autumnal equinoxes. This is because the northward shift observed at the vernal equinoxes is larger than that observed at the autumnal equinoxes since the TEC values are larger at the vernal equinoxes (e.g. Essex 1977). In Fig. 10, one can see that the present method underestimates the amplitudes of the semi-annual oscillation though it successfully reproduces the phase of the oscillation. For example, the estimated peak amplitudes at the solstices explain only 30–40 per cent of the observed amplitudes. We checked how the assumptions of an uniform sky distribution and of the sectoral representation of the Figure 11. Schematic figure that illustrates the positioning error at the solstices. (Top left) In the GEONET, the excess northward positioning error is larger at the southern part of the network near the magnetic equator. (Top right) After subtracting the network shifts, the contraction appears in the north south direction. (Bottom left) In the global network, the distribution of the northward shift is symmetrical across the magnetic equator. (Bottom right) After subtracting the network shifts, no contraction/expansion appears in the north–south direction. satellite visibility affect the amplitudes. For this purpose, we reestimated the positioning error and the Helmert parameters using the observed satellite positions rather than assuming the uniform satellite distribution, substituting the integration over the sphere in eq. (16) with the summation over the observed satellites. We found that the amplitudes of the scale change differ only by 3 per cent. Therefore, we conclude that the assumption of uniform satellite distribution and of the sectoral representation of the satellite visibility is not a main cause for the amplitude difference. We suppose that the factors which are not taken into account in this study such as the coupling between the position estimates and the zenithal tropospheric delay parameters (e.g. Hatanaka 2003) play important roles in the amplitude difference. Rigorous numerical simulation is needed to specify the causes for the amplitude difference. The additional Helmert parameters obtained in Section 3 may have contributed to the difference in the GPS-derived reference frame and the ITRF reference frame. The ionospheric effect can lead to annual scale changes of 0.1 ppb, an offset of 1.8 mm and a semi-annual oscillation of 1 mm in the z-direction. Heflin et al. (2002) examined the difference between these two frames. Their results show that the standard deviations are 1.0, 1.0 and 1.5 cm for the x, y and z geocentre components and 0.3 ppb for the scale. The additional Helmert transformations from this study are not sufficient to explain these large deviations. Therefore, we conclude that the contribution of the second-order ionospheric correction on the difference between the reference frames is negligible. We showed that the time-series of the position error estimated by the semi-analytical method shows a good correlation with the
C 2005 RAS, GJI, 163, 10–17 Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 Figure 10. Observed (red) and simulated (blue) scale changes in the GEONET. A semi-analytical estimation of the second-order ionospheric effect 6 C O N C LU S I O N We developed the semi-analytical method for evaluating the effect of second-order ionospheric corrections on the coordinate time-series. Our method successfully reproduced the time-series of the station positions of three IGS sites, which was comparable to the rigorous simulation carried out by Kedar et al. (2003). We also succeeded in qualitatively explaining the characteristics of the scale change timeseries observed in the GEONET in Japan. These results show the capability of the semi-analytical method in long-term estimations of the position error caused by various error sources such as water vapour, to obtain physical insights on how these error sources affect the position error. We estimated the additional Helmert parameters due to the ionospheric effect. Our results show that the scale has an annual oscillation whose amplitude is 0.1 ppb, network shifts have an offset of 1.8 mm and a semi-annual oscillation of 1 mm in the z-direction. However, the additional Helmert parameters are too small to explain the current deviations observed in the GPS-derived reference frame and the ITRF reference frame. AC K N OW L E D G M E N T S We are grateful to Jennifer Bonin for helpful comments on the English presentation. We would like to thank John Beavan, Kosuke Heki and an anonymous reviewer for their constructive suggestions. REFERENCES Bassiri, S. & Hajj, G.A., 1993. Higher-order ionospheric effects on the global positioning systems observables and means of modeling them, Manuscr. Geod., 18, 280–289.
C 2005 RAS, GJI, 163, 10–17 Brunner, F.K. & Gu, M., 1991. An improved model for dual frequency ionospheric correction of GPS observations, Manuscr. Geod., 16, 205– 214. Cucurull, L., Sedo, P., Behrend, E., Cardellach, E. & Ruis, A., 2002. Integrating NWP products into the analysis of GPS observables, Phys. Chem. Earth, 27, 377–383. Essex, E.A., 1977. Equinoctial variations in the total electron content of the ionosphere at northern and southern hemisphere stations, J. Atmos. Terr. Phys., 39, 645–650. Geiger, L., 1988. Simulating disturbances in GPS by continuous satellite distribution, J. Surv. Eng., 114, 182–194. Hatanaka, Y., 2003. Seasonal variation of scale of GEONET network and ZTD biases, paper presented at the Symposium JSG01, 23nd IUGG General Assembly, Sapporo, Japan. Hatanaka, Y., 2004. Re-evaluation of seasonal variations of scale in the GEONET solutions, paper presented at the 2004 fall meeting of Geodetic Society of Japan. Heflin, M., Argus, D., Jefferson, D., Webb, F. & Zumberge, J., 2002. Comparison of a GPS-defined global reference frame with ITRF 2000, GPS Solutions, 6, 72–75. Heki, K., 2001. Seasonal modulation of interseismic strain buildup in northeastern Japan driven by snow loads, Science, 293, 89–92. Heki, K., 2004. Dense GPS array as a new sensor of seasonal changes of surface loads, in The Sate of the Planet: Frontiers and Challenges in Geophysics, Vol. 150, pp.177–196, eds Sparks, R.S.J. & Hawkesworth, C.J., Geophys. Monogr., AGU, Washington DC. Hofmann-Wellenhof, B., Lichtenegger, H. & Collins, J., 2001. Global positioning system: theory and practice, Springer-Verlag, New York. Kedar, S., Hajj, A., Wilson, B.D. & Heflin, M.B., 2003. The effect of the second order GPS ionospheric correction on receiver positions, Geophys. Res. Lett., 30, 1829, doi:10.1029/2003GL017639. Mandea, M. et al., 2000. International geomagnetic reference field 2000, Geophys. J. Int., 141, 259–262. Mannucchi, A., Wilson, B.D., Yuan, D.N., Ho, U.J., Lindqwister, U.J. & Runge, T.F., 1998. A global mapping technique for GPS-derived ionospheric total electron content measurements, Radio Sci., 39, 1S10, 10.1029/2002RS002846. Miyazaki, S., Segall, P., Fukuda, J. & Kato, T., 2004. Space time distribution of afterslip following the 2003 Tokachi-oki earthquake: implications for variations in fault zone frictional properties, Geophys. Res. Lett., 31, L06623, doi:10.1029/2003GL019410. Munekane, H. & Matsuzaka, S., 2004. Nontidal ocean mass loading detected by GPS observation in the tropical Pacific region, Geophys. Res. Lett., 31, 08602, doi:10.1029/2004GL019773. Munekane, H., Tobita, M. & Takashima, K., 2004. Groundwater-induced vertical movements observed in Tsukuba, Japan, Geophys. Res. Lett., 31, 12 608, doi:10.1029/2004GL020158. Park, K.D. et al., 2004. Development of an antenna and multipath calibration system for Global Positioning System sites, Radio Sci., 39, 5002, 10.1029/2003RS002999. Santerre, R., 1991. Impact of GPS satellite sky distribution, Manuscr. Geod., 16, 28–53. Zumberge, J.F., Heflin, M.B., Jefferson, D.C., Watkins, M.M. & Webb, F.H., 1997. Precise point positioning for the efficient and robust analysis of GPS data from large networks, J. geophys. Res., 102, 5005–5017. Downloaded from https://academic.oup.com/gji/article/163/1/10/2087370 by guest on 17 July 2022 observed time-series. In addition, we showed that the semi-analytical method qualitatively reproduced the characteristics of the scale change time-series observed in the GEONET. These results show that it is possible to use the semi-analytical method as an effective tool to study the long-term (seasonal, interannual) effects of the second-order ionospheric correction or other error sources, such as water vapour, on the coordinate time-series to obtain physical insights on how these error sources affect the coordinate time-series. Finally, we would like to point out that the positioning error will change when the analytical strategies are changed. In particular, when the satellite position is fixed during the analysis, the positioning error will be small provided the satellite positions also contain the ionospheric error. Therefore, when making a correction for the ionospheric effect in the analysis with fixed satellite positions, we need to know whether the satellite positions are contaminated with the ionospheric effect. 17