Nội dung

Global Positioning Systems, Inertial Navigation, and Integration, Mohinder S. Grewal, Lawrence R. Weill, Angus P. Andrews Copyright # 2001 John Wiley & Sons, Inc. Print ISBN 0-471-35032-X Electronic ISBN 0-471-20071-9 5 GPS Data Errors 5.1 SELECTIVE AVAILABILITY ERRORS Prior to May 1, 2000, Selective Availability (SA) was a mechanism adopted by the Department of Defense (DoD) to control the achievable navigation accuracy by nonmilitary GPS receivers. In the GPS SPS mode, the SA errors were speci®ed to degrade navigation solution accuracy to 100 m (2D RMS) horizontally and 156 m (RMS) vertically. In a press release on May 1, 2000, the President of the United States announced the decision to discontinue this intentional degradation of GPS signals available to the public. The decision to discontinue SA was coupled with continuing efforts to upgrade the military utility of systems using GPS and supported by threat assessments which concluded that setting SA to zero would have minimal impact on United States national security. The decision was part of an ongoing effort to make GPS more responsive to civil and commercial users worldwide. The transition as seen from Colorado Springs, CO., U.S.A. at the GPS Support Center is shown in Figure 5.1. The ®gure shows the horizontal and vertical errors with SA, and after SA was suspended, midnight GMT (8 PM EDT), May 1, 2000. Figure 5.2 shows mean errors with and without SA, with satellite PRN numbers. Aviation applications will probably be the most visible user group to bene®t from the discontinuance of SA. However, precision approach will still require some form of augmentation to ensure that integrity requirements are met. Even though setting SA to zero reduces measurement errors, it does not reduce the need for and design of WAAS and LAAS ground systems and avionics. 103 104 GPS DATA ERRORS Fig. 5.1 Pseudorange residuals with SA. Time and frequency users may see greater effects in the long term via communication systems that can realize signi®cant future increases in effective bandwidth use due to tighter synchronization tolerances. The effect on vehicle tracking applications will vary. Tracking in the trucking industry requires accuracy only good enough to locate in which city the truck is, whereas public safety applications can require the precise location of the vehicle. Maritime applications have the potential for signi®cant bene®ts. The personal navigation consumer will bene®t from the availability of simpler and less expensive products, resulting in more extensive use of GPS worldwide. Because SA could be resumed at any time, for example, in time of military alert, one needs to be aware of how to minimize these errors. There are at least two mechanisms to implement SA. Mechanisms involve the manipulation of GPS ephemeris data and dithering the satellite clock (carrier frequency). The ®rst is referred to as epsilon-SA (e-SA), and the second as clockdither SA. The clock-dither SA may be implemented by physically dithering the frequency of the GPS signal carrier or by manipulating the satellite clock correction data or both. Though the mechanisms to implement SA and the true SA waveform are classi®ed, a variety of SA models exist in the literature (e.g., [4, 15, 21, 134]). These references show various models. One proposed by Braasch [15] appears to be the most promising and suitable. Another used with some success for predicting SA is a Levinson predictor [6]. 5.1 HDOP < =1.0 HDOP < =2.0 HDOP > 2.0 Display All Logging All SELECTIVE AVAILABILITY ERRORS SA Watch (SA) 150m 0.9 100m 11 05 50m 21 15 30 25 60ϒ 23 30ϒ Horizon Weighted Mean HDOP < =1.0 HDOP < =2.0 HDOP > 2.0 Display All Logging All SA Watch (no SA) 150m 100m 05 50m 29 11 30 25 60ϒ 21 15 23 22 Weighted Mean Fig. 5.2 Pseudorange residuals with and without SA. 30ϒ Horizon 105 106 GPS DATA ERRORS The Braasch model assumes that all SA waveforms are driven by normal white noise through linear system (autoregressive moving-average, ARMA) models (see Chapter 3 of [46]). Using the standard techniques developed in system and parameter identi®cation theory, it is then possible to determine the structure and parameters of the optimal linear system that best describes the statistical characteristics of SA. The problem of modeling SA is estimating the model of a random process (SA waveform) based on the input=output data. The technique to ®nd an SA model involves three basic elements:  the observed SA,  a model structure, and  a criterion to determine the best model from the set of candidate models. There are three choices of model structures: 1. an ARMA model of order ( p; q), which is represented as ARMA( p; q); 2. an ARMA model of order ( p; 0† known as the moving-average MA( p) model; and 3. an ARMA model of order q; 0†, the auto regression AR(q) model. Selection from these three models is performed with physical laws and past experience. 5.1.1 Time Domain Description Given observed SA data, the identi®cation process repeatedly selects a model structure and then calculates its parameters. The process is terminated when a satisfactory model, according to a certain criterion, is found. We start with the general ARMA model. Both the AR and MA models can be viewed as special cases of an ARMA model. An ARMA p; q† model is mathematically described by a1 y k ‡ a 2 y k 1 ‡


‡ aq y k q‡1 ˆ b1 x k ‡ b2 x k 1 ‡


‡ bp x k p‡1 ‡ ek 5:1† or in a concise form by q P iˆ1 ai y k i‡1 ˆ p P jˆ1 bj x k j‡1 ‡ ek ; 5:2† where ai ; i ˆ 1; 2; . . . ; q, and bj ; j ˆ 1; 2; . . . ; p, are the sets of parameters that describe the model structure, xk and yk are the input and output to the model at any time k for k ˆ 1; 2; . . . , and ek is the noise value at time k. Without loss of generality, it is always assumed that al ˆ 1. 5.1 SELECTIVE AVAILABILITY ERRORS 107 Once the model parameters ai and bj are known, the calculation of yk for an arbitrary k can be accomplished by yk ˆ q P iˆ2 a i yk i‡1 ‡ p P jˆ1 bj x k j‡1 ‡ ek : 5:3† It is noted that when all of the ai in Eq. 5.3 take the value of 0, the model is reduced to the MA p; 0† model or simply MA( p). When all of the bj take the value of 0, the model is reduced to the AR(0; q) model or AR(q). In the latter case, yk is calculated by yk ˆ q P iˆ2 ai y k i‡1 ‡ ek : 5:4† 5.1.1.1 Model Structure Selection Criteria Two techniques, known as Akaike's ®nal prediction error (FPE) criterion and the closely related Akaike's information theoretic criterion (AIC), may be used to aid in the selection of model structure. According to Akaike's theory, in the set of candidate models, the one with the smallest values of FPE or AIC should be chosen. The FPE is calculated as FPE ˆ 1 ‡ n=N V; 1 n=N 5:5† where n is the total number of parameters of the model to be estimated, N is the length of the data record, and V is the loss function for the model under consideration. Here, V is de®ned as V ˆ n P iˆ1 e2i ; 5:6† where e is as de®ned in Eq. 5.2. The AIC is calculated as AIC ˆ log‰ 1 ‡ 2n=N †V Š 5:7† In the following, an AR(12) model was chosen to characterize SA. This selection was based primarily on Braasch's recommendation [14]. As such, the resulting model should be used with caution before the validity of this model structure assumption is further studied using the above criteria. 5.1.1.2 Frequency Domain Description The ARMA models can be equivalently described in the frequency domain, which provides further insight 108 GPS DATA ERRORS into model behavior. Introducing a one-step delay operator Z 1, Eq. 5.2 can be rewritten as AZ 1 †yk ˆ B Z l †xk ‡ ek ; 5:8† where AZ 1 BZ 1 q P †ˆ iˆ1 †ˆ ai Z p P iˆ1 i‡l ; 5:9† bi Z i‡1 ; 5:10† and Z 1 yk ˆ yk 1 : 5:11† It is noted that A Z 1 † and B Z 1 † are polynomials of the time-shift operator Z 1 and normal arithmetic operations may be carried out under certain conditions. De®ning a new function H Z 1 † as B Z 1 † divided by A Z 1 † and expanding the resulting H Z 1 † in terms of operator Z 1, we have H Z 1 †ˆ BZ AZ 1 † 1† ˆ 1 P iˆ1 hi Z i‡1 : 5:12† The numbers of {hi } are the impulse responses of the model. It can be shown that hi is the output of the ARMA model at time i ˆ 1; 2; . . . when the model input xi takes the value of zero at all times except for i ˆ 1. The function H Z 1 † is called the frequency function of the system. By evaluating its value for Z 1 ˆ e jo, the frequency response of the model can be calculated directly. Note that this process is a direct application of the de®nition of the discrete Fourier transform (DFT) of hi . 5.1.1.3 AR Model Parameter Estimation with structure AZ 1 The parameters of an AR model †yk ˆ ek 5:13† may be estimated using the least-squares (LS) method. If we rewrite Eq. 5.13 in matrix format for k ˆ q; q ‡ 1; . . . ; n, we get 2 yn 6 yn l 6 6 .. 4 . yq yn yn .. . yq l 2 1 32 3 2 3


yn q‡1 en a1


y n q 7 6 a2 7 6 e n 1 7 76 7 6 7 .. 76 .. 7 ˆ 6 .. 7 .. 5 4 5 4 5 . . . . a e


y1 q q 5:14† 5.1 SELECTIVE AVAILABILITY ERRORS 109 or H
A ˆ E; 5:15† where 2 yn 6 6 yn l 6 H ˆ6 . 6 .. 4 yq yn yn .. . yq 2 3


yn q‡1 7


yn q 7 7 ; .. 7 .. . 7 . 5 1


l 5:16† y1 T A ˆ ‰a1 a2 a3


aq Š ; 5:17† and E ˆ ‰en en 1 en 2


eq Š T : 5:18† The LS estimation of the parameter matrix A can then be obtained by A ˆ H T H† 1 H T E: 5.1.2 5:19† Collection of SA Data To build effective SA models, samples of true SA data must be available. This requirement cannot be met directly as the mechanism of SA generation and the actual SA waveform are classi®ed. The approach we take is to extract SA from ¯ight test data. National Satellite Test Bed (NSTB) ¯ight tests recorded the pseudorange measurements at all 10 RMS (Reference Monitoring Station) locations. These pseudorange measurements contain various clock, propagation, and receiver measurement errors, and they can, in general, be described as PRM ˆ r ‡ DTsat ‡ DTrcvr ‡ DTiono ‡ DTtrop ‡ DTmultipath ‡ SA ‡ Dtnoise 5:20† where r is the true distance between the GPS satellite and the RMS receiver; DTsat and DTrcvr are the satellite and receiver clock errors; DTiono and DTtrop are the ionosphere and troposphere propagation delays, DTmultipath is the multipath error; SA is the SA error; and Dtnoise is the receiver measurement noise. To best extract SA from PRM , values of the other terms were estimated. The true distance r is calculated by knowing the RMS receiver location and the precise orbit data available from the National Geodetic Survey (NGS) bulletin board. GIPSY= OASIS analysis (GOA) was used for this calculation, which re-created the precise orbit and converted all relevant data into the same coordinate system. Models for propagation and satellite clock errors have been built into GOA, and these were used 110 GPS DATA ERRORS to estimate DTsat ; DTiono , and DTtrop . The receiver clock errors were estimated by the NSTB algorithm using data generated from GOA for the given ¯ight test conditions. From these, a simulated pseudorange PRsim was formed.: PRsim ˆ rsim ‡ DTsatsim ‡ DTrcvrsim ‡ DTionosim ‡ DTtropsim 5:21† where DTsatsim ; DTrcvrsim ; DTionosim , and DTtropsim are, respectively, the estimated values of DTsat ; DTrcvr ; DTiono , and DTtrop in the simulation. From Eqs. 5.20 and 5.21, pseudorange residuals are calculated: DPR ˆ PRM PRsim ˆ SA ‡ DTmultipath ‡ Dtnoise ‡ DTmodels ; 5:22† where DTmodels stands for the total modeling error, given by   
 DTmodels ˆ r rsim ‡ DTsat DTsatsim ‡ DTrcvr DTrcvrsim     ‡ DTiono DTionosim ‡ DTtrop DTtropsim 5:23† It is noted that the terms DTmultipath and Dtnoise should be signi®cantly smaller than SA, though it is not possible to estimate their values precisely. The term DTmodels should also be negligible compared to SA. It is, therefore, reasonable to use DPR as an approximation to the actual SA term to estimate SA models. Examination of all available data show that their values vary between 80 m. These are consistent with previous reports on observed SA and with the DoD's speci®cation of SPS accuracy. 5.2 IONOSPHERIC PROPAGATION ERRORS The ionosphere, which extends from approximately 50 to 1000 km above the surface of the earth, consists of gases that have been ionized by solar radiation. The ionization produces clouds of free electrons that act as a dispersive medium for GPS signals in which propagation velocity is a function of frequency. A particular location within the ionosphere is alternately illuminated by the sun and shadowed from the sun by the earth in a daily cycle; consequently the characteristics of the ionosphere exhibit a diurnal variation in which the ionization is usually maximum late in midafternoon and minimum a few hours after midnight. Additional variations result from changes in solar activity. The primary effect of the ionosphere on GPS signals is to change the signal propagation speed as compared to that of free space. A curious fact is that the signal modulation (the code and data stream) is delayed, while the carrier phase is advanced by the same amount. Thus the measured pseudorange using the code is larger than the correct value, while that using the carrier phase is equally smaller. The magnitude of either error is directly proportional to the total electron count (TEC) in a tube of 1 m2 cross section along the propagation path. The TEC varies spatially 5.2 IONOSPHERIC PROPAGATION ERRORS 111 due to spatial nonhomogeneity of the ionosphere. Temporal variations are caused not only by ionospheric dynamics, but also by rapid changes in the propagation path due to satellite motion. The path delay for a satellite at zenith typically varies from about 1 m at night to 5±15 m during late afternoon. At low elevation angles the propagation path through the ionosphere is much longer, so the corresponding delays can increase to several meters at night and as much as 50 m during the day. Since ionospheric error is usually greater at low elevation angles, the impact of these errors could be reduced by not using measurements from satellites below a certain elevation mask angle. However, in dif®cult signal environments, including blockage of some satellites by obstacles, the user may be forced to use low-elevation satellites. Mask angles of 5 7:5 offer a good compromise between the loss of measurements and the likelihood of large ionospheric errors. The L1 -only receivers in nondifferential operation can reduce ionospheric pseudorange error by using a model of the ionosphere broadcast by the satellites, which reduces the uncompensated ionospheric delay by about 50% on the average. During the day errors as large as 10 m at midlatitudes can still exist after compensation with this model and can be much worse with increased solar activity. Other recently developed models offer somewhat better performance. However, they still do not handle adequately the daily variability of the TEC, which can depart from the modeled value by 25% or more. The L1 =L2 receivers in nondifferential operation can take advantage of the dependency of delay on frequency to remove most of the ionospheric error. A relatively simple analysis shows that the group delay varies inversely as the square of the carrier frequency. This can be seen from the following model of the code pseudorange measurements at the L1 and L2 frequencies: ri ˆ r  k ; fi2 5:24† where r is the error-free pseudorange, ri is the measured pseudorange, and k is a constant that depends on the TEC along the propagation path. The subscript i ˆ 1, 2 identi®es the measurement at the L1 or L2 frequencies, respectively, and the plus or minus sign is identi®ed with respective code and carrier phase pseudorange measurements. The two equations can be solved for both r and k. The solution for r for code pseudorange measurements is rˆ f12 f12 f22 r1 f12 f22 f22 r2 ; 5:25† where f1 and f2 are the L1 and L2 carrier frequencies, respectively, and r1 and r2 are the corresponding pseudorange measurements. An equation similar to the above can be obtained for carrier phase pseudorange measurements. However, in nondifferential operation the residual carrier phase 112 GPS DATA ERRORS pseudorange error can be greater than either an L1 or L2 carrier wavelength, making ambiguity resolution dif®cult. With differential operation ionospheric errors can be nearly eliminated in many applications, because ionospheric errors tend to be highly correlated when the base and roving stations are in suf®ciently close proximity. With two L1 -only receivers separated by 25 km, the unmodeled differential ionospheric error is typically at the 10±20-cm level. At 100 km separation this can increase to as much as a meter. Additional error reduction using an ionospheric model can further reduce these errors by 25±50%. 5.2.1 Ionospheric Delay Model J. A. Klobuchar's model [37,76] for ionospheric delay in seconds is given by  Tg ˆ DC ‡ A 1 where x ˆ 2p t x2 x4 ‡ 2 24  for jxj  p 2 5:26† Tp † (rad) P DC ˆ 5 ns (constant offset) Tp ˆ phase ˆ 50,400 s A ˆ amplitude P ˆ period t ˆ local time of the earth subpoint of the signal intersection with mean ionospheric height (s) The algorithm assumes this latter height to be 350 km. The DC and phasing Tp are held constant at 5 ns and 14 h (50,400 s) local time. Amplitude (A) and period (P) are modeled as third-order polynomials: Aˆ 3 P nˆ0 an fnm s†; Pˆ 3 P nˆ0 bn fnm s†; where fm is the geomagnetic latitude of the ionospheric subpoint and an ; bn are coef®cients selected (from 370 such sets of constants) by the GPS master control station and placed in the satellite navigation upload message for downlink to the user. For Southbury, Connecticut, an ˆ ‰0:8382  10 8 ; bn ˆ ‰0:8806  105 ; 0:745  10 8 ; 0:3277  105 ; The parameter fm is calculated as follows: 0:596  10 7 ; 0:596  10 7 Š; 0:1966  106 ; 0:1966  106 Š