EAS591T – Space Geodetic Measurements of Active Crustal
Motions
LAB 7
The pseudorange measurements ^{j}R_{i}(t) can be modeled as:
_{} (1)
t = time of epoch
^{j}R_{i} = pseudorange measurement
^{j}r_{i} = satellitereceiver geometric distance
c = speed of light
^{j}d = satellite clock bias
d_{i} = receiver clock bias
DI = ionospheric propagation error
DT = tropospheric propagation error
MP = multipath
e = receiver noise
(ranges in meters, time in seconds)
Neglecting the propagation, multipath, and receiver errors, eq.(1) becomes:
_{} (2)
The geometric distance between satellite j and receiver i is given by:
_{} (3)
with [^{j}X, ^{j}Y, ^{j}Z] = satellite position, [X_{i}, Y_{i}, Z_{i}] = receiver position in an ECEF coordinate system.
Our mission, if we accept it, is to solve for [X_{i}, Y_{i}, Z_{i}, d_{i}], assuming that we know [^{j}X, ^{j}Y, ^{j}Z, ^{j}d]. A major problem here is that the unknowns [X_{i}, Y_{i}, Z_{i}] are not linearly related to the observables…
Assuming that we now the approximate coordinates of the receiver [X_{o}, Y_{o}, Z_{o}], one can write that the actual coordinates equal the approximate coordinates plus a slight adjustment:
_{} (4)
DX_{i},
DY_{i}, DZ_{i} are our new unknowns. We can now write:
_{} (5)
Since we know the approximate point [X_{o}, Y_{o}, Z_{o}], we can now expand f(X_{o}+DX_{i}, Y_{o}+DY_{i}
, Z_{o}+DZ_{i}) using a Taylor’s series with respect to that point:
_{} (6)
We intentionally truncate the Taylor’s expansion after the linear
terms. Recall from eq.(3) that:
_{} (7)
The partial derivatives in eq.(6) are therefore given by:
_{} (8)
We can now substitute eq.(8) into eq.(6):
_{} (9)
We now have an equation that is linear with respect to the unknowns DX_{i}, DY_{i}, DZ_{i}.
Now let us go back to our pseudorange measurements ^{j}R_{i}(t) and rewrite eq.(2):
_{} (10)
We can rearrange eq.(10) by separating the known and unknown terms of
each side (recall that the satellite clock correction ^{j}d(t)
is provided in the navigation
message):
_{} (11)
We can simplify the notation by assigning:
_{} (12)
Let us assume that we have 4 satellites visible simultaneously. We use
eq.(11) and write it for the 4 satellites::
_{} (13)
Tired of carrying along all these terms, subscripts, and superscripts? Me too. Let us introduce:
_{} (14)
L = vector of n
observations. Must have at
least 4 elements (i.e. 4 satellites), but in reality will have from 4 to 12
elements depending on the satellite constellation geometry.
X = vector of u
unknowns. Four elements in
our case.
A = matrix of linear
functions of the unknowns (= design matrix), n rows by u columns.
Now we can write our problem (eq.13) in a matrixvector form:
_{} (15)
In
general, n > u, leading to an overdetermined system.
Because actual data contain observational errors and noise, this system is in
nonconsistent. In order to make it consistent, one must account for a noise
vector r. Eq.(15) becomes:
_{} (16)
The
“noise vector” r represents residuals, i.e. observations (L)
minus model (AX). The least squares solution to
eq.(16) is:
_{} (17)
P is the weight matrix,
defined by:
_{} (18)
s_{o}^{2} = a priori variance
S_{L} = covariance matrix of the observations.
The law of covariance propagation gives the covariance matrix of the
unknowns S_{X}:
_{} (19)
In the case of pseudoranges, the observations are independant and have
equal variance s_{o}^{2}. Therefore S_{L }is the diagonal matrix:
_{} (20)
Assuming that the weight matrix is I, eq.(17) can be
simplified to:
_{} (21)
Now that DX_{i}, DY_{i}, DZ_{i} are found, the antenna coordinates [X_{i}, Y_{i},
Z_{i}] are obtained
using eq.(4).
The associated covariance matrix of the unknowns S_{x} is:
_{} (22)
We can transform S_{x }from an ECEF frame to a local topocentric frame using the law of variance propagation (disregarding the timecorrelated components of S_{x}:
_{} (23)
where R is
the rotation matrix (cf. lab 1):
_{} (24)
with j = geodetic latitude of the site, l = geodetic longitude of the site.
The DOP factors (Dilution Of Precision) are given by:
_{} (25)
Assignment:
Write a MATLAB program to compute the position and clock bias of a GPS
receiver and the GDOP using:
Compare solutions using C1, P1, and P2
The a priori position of the receiver in ECEF frame (in meters) is:
Xo=4433470.0
Yo=362670.0
Zo=4556210.0
Possible program structure:
I find:

C1 
P1 
P2 
DX 
37.448 
36.926 
36.136 
DY 
52.132 
52.521 
53.487 
DZ 
60.883 
60.628 
59.357 
DT
(nsec) 
187.16 
186.30 
172.53 
Xa 
4433432.552 
4433433.074 
4433433.864 
Ya 
362722.131 
362722.521 
362723.487 
Za 
4556149.117 
4556149.372 
4556150.643 
GDOP 
5.2 
5.2 
5.2 