Pseudoranges contain a number of errors (see e.g. the theory page). Some of these errors can be estimated and used as corrections to the pseudoranges. In this page the following corrections are described:
Again, the GPS SPS
Signal Specification is followed closely.
A fourth correction, the receiver clock correction, results from the final
stand alone position calculation and is emphasized
in that page.
1. SV Clock Correction
The parameters, required to calculate the SV clock correction are present in subframe 1 of the navigation message (see e.g. the conversion page). They are:
The SV clock correction for the C/A code pseudorange observation is given by the expression:
dTclck = af0 + af1 * (Ttr - Toc) + af2 * (Ttr - Toc)^2 + F * ec * rootA * sin(E) - Tgd
with:
The implementation of this expression is straightforward.
The term: F * ec * rootA * sin(E) is called the relativistic correction term
Trel [sec]. It may be a good idea to calculate this term within the
'satpos' procedure since the required parameters
are available there. The term: T = Ttr - Toc should be corrected for end
of week crossovers: if T > 302400 then subtract 604800 from T, if T <
-302400 then add 604800 to T.
Below follow fragments of the satpos procedure with the additions for implementation of Trel in italics.
procedure satpos (eph : vec16;
{ephemeris}
Ttr : real;
{satellite GPS time}
var Trel: real;
{satellite relativistic correction term}
var X : vec3); {satellite
position}
const
pi = 3.1415926535898; {WGS 84 value
of pi}
mu = 3.986005E+14; {WGS
84 value of earth's univ. grav. par.}
Wedot = 7.2921151467E-5; {WGS 84 value of earth's rotation
rate}
F = -4.442807633E-10; {rel cor
term constant}
var
M0, dn, ec, A, W0, i0, w Wdot, Cuc, Cus, Crc, Crs, Cic, Cis,
Toe, Idot: real;
n0, M, snu, cnu, nu, phi, du, dr, di, u, r, i, Xdash, Ydash,
Wc: real;
i: integer;
begin
code .....
..........
code ....
X[1] := Xdash*cos(Wc) - Ydash*cos(i)*sin(Wc);
X[2] := Xdash*sin(Wc) + Ydash*cos(i)*cos(Wc);
X[3] := Ydash*sin(i)
Trel := F * ec * rootA * sin(e); {relativistic correction term}
end; {procedure satpos}
An example. Take the clock correction parameters from the example subframe 1 in the conversion page, and take an SV transmission time of 588480. When testing the satpos procedure with the same transmit time, I printed the parameter E = -1.5747877 rad. Entering all parameters in the expression above gives for dTcl: -3.649408E-6 sec, or, multiplied by the speed of light: -1094.065 m.
2. Ionospheric Delay.
Subframe 4 page 18 contains the ionospheric delay coefficients. Furthermore, the approximate user position, and the SV's azimuth and elevation are required. And here we run into a problem: the user position is the final result of the calculation, and is not yet known at this stage. Some solutions are:
In the following I assume, that the approximate user position is known. The
following calculation steps are extracted from 'Ionospheric Effects on GPS'
by J.A. Klobuchar. All angular units are in semicircles (!), El is the SV's
elevation [sc], Az [sc] the SV's azimuth, Latu [sc] the appr. user latitude,
and Lonu [sc] the appr. user longitude(the calculation of El and Az is given
in the Calculation of SV's Elevation and Azimuth
page).
Remarks:
An example of a Turbo Pascal procedure has been incorporated in a Turbo Pascal program to calculate receiver position from raw data.This program can be viewed or downloaded by clicking here.
3. Tropospheric Delay.
Two strategies can be followed to determine the tropo delay:
There are numerous tropo models availble, I picked Hopfield's model which splits the tropo delay into two parts: a contribution of the 'dry' atmosphere and a contribution of the 'wet' atmosphere. For a detailed explanation, see any textbook on GPS.
The zenith delay of the dry component is given by:
Kd = 1.55208E-4 * Pamb * (40136.0 + 148.72 * Tamb) / (Tamb + 273.16)
[m]
The zenith delay of the wet component is given by:
Kw = -0.282 * Pvap / (Tamb + 273.16) + 8307.2 * Pvap/ (Tamb + 273.16)^2
[m]
Multiply the zenith delay's with their 'mapping functions' to correct for elevations lower than 90 deg and add the components to obtain the SV's tropospheric delay correction:
dRtrop = Kd / sin(sqrt(El * El + 1.904E-3)) + Kw / sin(sqrt(El * El + 0.6854E-3)) [m]
with El the SV's elevation in [rad].
When we use the standard atmosphere values, the expression reduces to
dRtrop = 2.312 / sin(sqrt(El * El + 1.904E-3)) + 0.084 / sin(sqrt(El * El + 0.6854E-3)) [m]
Example. Using standard atmosphere values for an SV with an elevation of 12.86 deg the Tropo delay is 10.575 m.
4. The corrected pseudorange.
The last step in this page is the determination of the corrected pseudorange:
PRcorr = PR + dTclck * c - dTiono * c - dRtrop [m], with c the speed of light.
Suppose the SV 30's raw pseudorange is 24401509.104 m, then the corrected
range is:
24401509.104 - 1094.065 - 3.811 - 10.575 = 24400400.653 m