Spatial and temporal correlations
of the GPS estimation errors
Borys Stoew and Gunnar Elgered
Onsala Space Observatory, Chalmers University of Technology, 439 92 Onsala, Sweden
Onsala
March 16, 2005
This report summarizes the work carried out within Work Package 3000 in the TOUGH
Project. It contains the final versions of Delivery D19 on spatial correlation of errors and
Delivery D22 on temporal correlation of errors.
TOUGH (Targeting Optimal Use of GPS Humidity Measurements in Meteorology) is a
project supported by the 5th Framework Programme of the European Commission.
Summary
We have studied the spatial correlation structure of the estimation errors when the GPS
technique is used to derive the atmospheric propagation delay in the zenith direction.
These error correlations are parameterized using an analytical function. For spatial scales
smaller than 200 km the GPS error correlations become significant; for larger spatial
scales, the error correlations are small and slowly decreasing with distance. These results
are based on a study of the differences between the propagation delays obtained from a
numerical weather model, and the corresponding GPS estimates.
The temporal correlations of the propagation delay estimation errors are also examined and the decorrelation times are of the order of 1–2 days. The temporal correlation
results are derived for the Swedish and Finnish GPS sites, using the relation between the
estimation errors for the vertical coordinate and the zenith propagation delay.
Contents
1 Introduction
1
2 Measuring and estimating the ZTD
2.1 Radiosonde data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 GPS data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 Microwave radiometer data . . . . . . . . . . . . . . . . . . . . . . . . . .
1
3
3
4
3 Extracting ZTD time series from HIRLAM
3.1 The non-linear observation operator for the ZWD . . . . . . . . . . . . . .
3.2 The tangent-linear observation operator . . . . . . . . . . . . . . . . . . .
3.3 ZTD Innovation vector . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4
4
5
5
4 Spatial correlations using analytic models
5
5 Spatial correlations based on NWP and GPS data
7
6 Spatial correlations using vertical coordinate residuals
12
7 Temporal correlations of slant delay residuals using ray tracing
7.1 Hydrostatic and wet mapping functions . . . . . . . . . . . . . . . . . . .
7.2 Ray tracing using radiosonde ascents . . . . . . . . . . . . . . . . . . . . .
13
13
14
8 Temporal correlations of ZTD from GPS position residuals
16
9 Discussion
18
References
19
ii
1
Introduction
In short term weather forecasting, various sources of data describing the atmospheric
systems and processes are used. The measurement techniques providing information
about the atmospheric parameters vary not only in sampling rate and spatial distribution,
but also in accuracy. The potential of using atmospheric delays estimated from GPS data
in near real time and the high temporal/spatial resolution of the measurements make the
GPS data assimilation worthwhile (Gutman et al., 2004).
Traditionally, the errors in the radiosonde measurements of humidity are assumed to
have no spatial correlation in the numerical weather prediction (NWP) model setup. Spatially correlated errors lead to a non-diagonal structure of the error covariance matrices,
whose inverting is computationally demanding. The use of a proper error-covariance
model ensures the positive definiteness of the respective covariance matrices, so that they
have no eigenvalues near zero and are non-singular. The errors in the estimates of atmospheric propagation delay, both from GPS processing and from NWP models, exhibit
spatial correlation structures which deserve special attention (Berre, 1997; Stoew et al.,
2001).
A statistical description of the temporal correlations of errors in the atmospheric delay
estimated from GPS data is needed for the proper operation of an NWP model — in particular, for the weighting of past zenith total delay (ZTD) estimates in an observation bias
reduction scheme (Lindskog, 2001; Daley, 1992). The spatial and temporal correlations
of the estimation errors are assessed in the following sections of this report.
Section 2 makes a brief overview of the relations between relevant atmospheric parameters and the ZTD. We introduce the radiosonde data which we later use to estimate
mapping function errors through the method of ray tracing. We also review the estimation of the ZTD using GPS and microwave radiometer data. Section 3 describes how
the ZTD time series are obtained from NWP models. Thereafter we focus on the spatial
(horizontal) correlation of ZTD errors. In Section 4, we discuss the analytic functions
used to model the horizontal correlations of the GPS estimation errors in the ZTD; we
provide the resulting parameters in Section 5. Using the residuals of the estimated vertical
coordinates, we assess the spatial error correlations in Section 6.
The temporal correlation of ZTD estimation errors is studied in Section 7 where we
discuss two different atmospheric delay mapping functions (MFs). The residual estimates
of the propagation delay, obtained by ray tracing of radiosonde data, are used to estimate
the temporal correlations of the errors contributed by the MFs used in the GPS processing.
The mapping function is only one component introducing a temporal correlation of the
ZTD estimation errors. For example, satellite orbits and clocks also contribute to the
effect. An attempt to estimate the total effect is made in Section 8 which presents
temporal correlations of the ZTD errors based on the assumption that these estimation
errors are strongly correlated with those of the estimated vertical positions of the GPS
sites. Section 9 ends the report with a discussion concerning the obtained results.
2
Measuring and estimating the ZTD
The refractivity in an air parcel is often written
N = k1
pd
e
e
+ k2 + k3 2 ,
T
T
T
(1)
where the parameters k1 , k2 and k3 are empirically determined from laboratory experiments, T is the absolute temperature, e represents the water vapour pressure, and p d is
the sum of the partial pressures of the “dry” gases. The refractivity is usually split into
1
a hydrostatic and a non-hydrostatic (wet) part:
N = N h + Nw .
(2)
Accordingly, the ZTD can be presented as a sum of a zenith hydrostatic delay (ZHD) and
a zenith wet delay (ZWD):
ℓ = ℓ h + ℓw .
(3)
The ZHD is due to the refractivity associated with the gasses in hydrostatic equilibrium,
obeying the hydrostatic equation for an ideal atmosphere:
dp
= −ρg .
(4)
dz
Above, p is the partial pressure of the gas at an altitude z, ρ is the gas density, and g
is the gravitational acceleration. The ZWD in (3) is mostly due to water vapour (e.g.,
Bevis et al., 1992).
Davis et al. (1985) proposed an approximate formula for estimating the ZHD at a
given site with latitude θ in degrees and height ho in km above the ellipsoid from pressure
measurements:
P
,
(5)
ℓh = (2.2768 ± 0.0024)
f (θ, ho )
where ℓh is in mm, the ground pressure P is in hPa, and
f (θ, ho ) = 1 − 0.00266 cos(2θ) − 0.00028 ho .
(6)
The amount of water vapour expressed as the height of the column equivalent liquid
water is termed integrated precipitable water vapour (IPWV) and is defined as
Z
1 ∞
I=
ρv dz ,
(7)
ρ l ho
where ρl and ρv are the densities of liquid water and water vapour. The IPWV is typically
measured in mm. The quotient
ℓw
,
(8)
Q=
I
where both ℓw and I are in mm, has been described by Askne and Nordius (1987):
µ
¶
k3
′
−8
Q = 10 ρl Rv k2 +
.
(9)
Tm
Above, Rv = 461.5 J kg−1 K−1 is the specific gas constant for water vapour; we assume
k2′ ≈ 22 K/hPa and k3 ≈ 3.7·105 K2 /hPa after Emardson and Derks (2000); the parameter
Tm is defined as the mean temperature of the water vapour:
Z ∞
ρv dz
ho
.
(10)
Tm = Z ∞
ρv
dz
ho T
Based on radiosonde profiles acquired over the US, Bevis et al. (1992) proposed the linear
regression Tm ≈ 70.2 + 0.72 To , where To is the surface temperature in Kelvin. Other
relations for Europe are reported by e.g. Emardson and Derks (2000).
The ZTD time series are dominated by the temporal variation of the ZWD, even though
the absolute values of the ZHD are about an order of magnitude larger those of the ZWD.
From the viewpoint of numerical weather modelling, the ZTD can easily be calculated
using integration over the different levels of the atmospheric model. Conversely, the ZTD
estimates based on the GPS measurements carry information about atmospheric pressure
and humidity. As the current amount of humidity measurements assimilated into NWP
models is insufficient for many applications, the GPS estimates of the ZTD are likely to
become valuable for the weather forecasting community.
2
2.1
Radiosonde data
Balloons carrying measurement sensors can be used to acquire vertical profiles of temperature, pressure, humidity, and wind. The regular launching of balloon-borne radiosondes (RS) is costly; this limitation has an impact on the amount of the acquired data and
on their quality. Recent studies over large regions of the northern hemisphere have shown
that a low RS launch frequency (e.g. only twice a day) may introduce bias of up to 3%
to the estimates of average amount of precipitable water (Dai et al., 2002).
The data obtained from the atmospheric soundings are often reduced to profiles containing the values at a set of standard pressure levels (surface, 1000, 850, 700, 500,
400, 300 hPa), along with the levels representing significant changes in the observed atmospheric parameters. The radiosondes used now in the Nordic countries are of the
Vaisala RS-80 type. Recent investigations on the performance of the Vaisala RS-80 relative humidity (RH) sensors suggest an age-related contamination by the packaging of
the sondes (Niell et al., 2001; Wang et al., 2002). As the sensors typically underestimate
the actual RH values by as much as 5% in very humid conditions, correction algorithms
should be applied (Lesht, 1999; Wade and Schwartz, 1993).
The vertical profiles of RH, and the air temperature and pressure can be used to derive
the partial pressure of water vapour and subsequently the wet refractivity. The RH
expressed in percent is defined as:
U% =
e
100 ,
esw
(11)
where e is the partial pressure of water vapour; esw is the saturation vapour pressure
over water surface, and characterizes the conditions for a given volume and temperature
when water can be found in equilibrium in both its liquid and gaseous phases. Empirical formulas for calculating esw (T ) have been derived by fitting a model to tabulated
measurements (e.g. Crane, 1976; List, 1958). Crane (1976) introduced the model
¶
µ
¶¸
·
µ
T
T − 273
− 5.31 ln
(12)
esw (T ) = 6.105 exp 25.22
T
273
using tabulated values of esw in hPa and T in Kelvin degrees. The model has an uncertainty within ±0.4% over the temperature range from 243 K to 303 K (−30 ◦ C to +30◦ C).
Note that the partial pressure of water vapour equals the saturation pressure at the dew
point temperature Td , i.e.
e = esw (Td ) .
(13)
2.2
GPS data
GPS is a navigation system using several satellite-borne radio transmitters. At present,
the constellation consists nominally of 24 active satellites in 6 orbital planes, allowing
for a simultaneous visibility of at least 4 satellites for any receiver on the surface of the
earth (Hofmann-Wellenhoff et al., 1997). For the purposes of meteorology (water vapour
in the atmosphere), high measurement accuracy is achieved by using state of the art GPS
receivers capable of measuring the carrier phase of the received signal (Blewitt, 1993).
Information about the moisture of the atmosphere can be extracted from the ZTD
fields, using a network of GPS receivers. The data processing can be carried out with the
GIPSY-OASIS software (Webb and Zumberge, 1993) using the Precise Point Positioning
solution (Zumberge et al., 1997). This processing is now routinely performed at the
analysis centre at the Onsala Space Observatory for a receiver network (Fig. 1) (e.g.,
Emardson et al., 1998).
3
KEVO
KEVO
KIRU
KIRU
SODA
SODA
OULU
OULU
ROMU
ROMU
SKEL
SKEL
VILH
VILH
KUUS
KUUS
OVER
OVER
ARJE
ARJE
UMEA
UMEA
OSTE
OSTE
JOEN
JOEN
VAAS
VAAS
SVEG
SVEG
SUND
SUND
OLKI
OLKI
KARL
KARL
VANE
VANE
GOTE
GOTE
RR
O
UO
TTU
LEKS MART
LEKS
MART
VIRO
VIRO
METS
METS
VAST
VAST
LOVO
LOVO
NORR
NORR
A
RA
OR
BBO SAA
S
NN
OO
JONK
JONK
VISB
VISB
OSKA
OSKA
10˚E
HASS
HASS
MALM
MALM
Fig. 1. Locations of the GPS sites in Sweden and Finland.
2.3
Microwave radiometer data
A microwave radiometer is used at the Onsala Space Observatory for independent assessment of atmospheric water vapour variability in general and the ZWD in particular.
In the following, it is referred to as a Water Vapour Radiometer (WVR). It operates at
21.0 and 31.4 GHz measuring the sky brightness temperatures with a long-term repeatability (order of years) of approximately 1 K. This corresponds to an absolute uncertainty
in the ZWD of 7 mm; the repeatability of the measured brightness temperatures becomes
better for shorter time scales. Typically, RMS differences obtained when WVR data are
compared to other independent measurements, such as radiosondes, are of the same order.
More details on the performance of the WVR have been documented elsewhere (Elgered
and Jarlemark, 1998).
3
3.1
Extracting ZTD time series from HIRLAM
The non-linear observation operator for the ZWD
The non-linear observation operator H(·) is used to project the background state onto
the observation space. The background IPWV estimate in the current HIRLAM version
is derived from (7):
N
X
1
I≈
qi (Pi+1/2 − Pi−1/2 ) ,
(14)
ρl gs (θ)
i=1
where gs (θ) is the gravity acceleration as a function of latitude of the site; qi is the
background specific humidity at model level i; N denotes the number of model levels. The
4
background pressure Pi±1/2 is evaluated at the model layer boundaries, and approximates
the partial pressure of the dry gases in the hydrostatic equation; the other parameters
are evaluated at the center of each grid cell. The ZWD part of the operator H(·) is based
on (9), and is defined by the numerically calculated IPWV for the model levels:
µ
¶
N
10−8 Rv X
k3
′
ℓw ≈
qi k2 +
(Pi+1/2 − Pi−1/2 ) ,
gs (θ)
Ti
(15)
i=1
where Ti is the background temperature at model level i.
3.2
The tangent-linear observation operator
The tangent-linear form H is the derivative of H(·) with respect to the state x of the
numerical weather model, evaluated at the background xb :
H = ∇x H(x) |x=xb .
(16)
The linear form of the observation operator for the IPWV over a single site, based on
(14), is:
" N
#
N
X
X
1
qi (P̄i+1/2 − P̄i−1/2 ) ,
(17)
q̄i (Pi+1/2 − Pi−1/2 ) +
I¯ ≈
ρl gs (θ)
i=1
i=1
where the tangent-linear pressure and specific humidity are denoted P̄i±1/2 and q̄i , respectively, and are the derivatives of the corresponding parameters at the background
state.
3.3
ZTD Innovation vector
ZTD estimates for one year were acquired, starting with 1 May 2000, at a sampling
rate of 5 minutes, acquired using 35 receivers on the territory of Sweden and Finland.
These time series were subsequently decimated at 6-hour sampling intervals in order to
match the HIRLAM forecast runs. These forecast runs were performed off-line at the
Swedish Meteorological and Hydrological Institute, and covered the same period. For
the calculation of the high resolution short term forecasts the HIRLAM 3DVAR was
used (Källén, 1996; Gustafsson et al., 2001; Lindskog et al., 2001).
The differences between the short term ZTD forecast fields and the estimates of ZTD
from GPS (the innovations) are shown in Fig. 2. Note that the HIRLAM–GPS innovation
time series (ℓH − ℓG) appear to have higher variability in the summer. For all the GPS
sites (Fig. 1) there is a persistent over estimation of ZTD by the GPS technique of about
1 cm. When a rectangular “sliding window” was applied to the differenced time series for
each of the GPS sites, there was no apparent seasonal dependence in the resulting slowly
varying bias between the HIRLAM and the GPS ZTDs. We do not correct the innovation
data for this time-varying component, and we only remove the mean difference when the
horizontal correlations of the innovations were calculated for each pair of GPS sites.
4
Spatial correlations using analytic models
When modelling atmospheric parameters, it is a common practice to assume horizontally
homogeneous, isotropic random field, locally time-stationary over a given period. Our
approach to spatial correlation modelling is based on the statistical analysis of turbulent
atmosphere (Tatarskii, 1971). We are interested in the spatial structure of the correlation
5
0.06
HIRLAM−GPS innovations in ZTD, m
0.04
0.02
0
−0.02
−0.04
−0.06
0
50
100
150
200
250
Days since May 1, 2000
300
350
400
Fig. 2. Time series of the differences between the HIRLAM and GPS zenith total delay estimates (i.e. the
innovations) for the site at Onsala. The corresponding time series of the innovations were produced for all GPS
sites, and then used to calculate the correlations of the GPS estimation errors.
coefficient ρ which is the normalized two-dimensional auto-covariance function of a random field. We use a sampled version (both in time and in space) of this field to calculate
ρ as a function of the site separation r (e.g., Stoew et al., 2001).
Different authors have proposed various analytical expressions which describe the spatial correlation of a random field as a function of distance (Daley, 1991; Tatarskii, 1971).
These expressions satisfy the requirement that any proper 2D correlation model should
yield a non-negative power spectral density (PSD) function. We consider the Gaussian
and the von Kármán analytical functions in order to model spatial correlations in polar
coordinates:
µ 2¶
r
(18)
̺1 (r) = a · exp − 2
r1
µ ¶ν µ ¶
b
r
r
̺2 (r) = ν−1
Kν
(19)
2 Γ(ν) r2
r2
where Kν (x) is the modified Bessel function of the second kind with parameter ν > 0; the
parameters a and b are scaling factors; r1 , r2 are characteristic lengths for the respective
function fitted to the data. Any linear combination of (18) and (19) forms a proper 2D
correlation function. To describe 2D spatial correlations, a second-order autoregressive
function is sometimes used in the literature (e.g., Daley, 1991); it is, however, a special
case of (19).
For a 2D isotropic field, the Fourier transform is equivalent to a zero-order Hankel
transform (Erdélyi et al., 1954). Thus we can write the normalized power spectral densities
6
corresponding to (18) and (19):
µ
¶
g1 (κ)
κ 2 r12
= exp −
g1 (0)
4
1
g2 (κ)
=¡
¢ν+1
g2 (0)
1 + κ 2 r22
(20)
(21)
where κ corresponds to a set of 2D spatial wave number pairs (k, l), for which κ 2 = k 2 +l2 .
This approach has its limitations which stem from the ignoring of the curvature of
the earth surface. Hollingsworth and Lönnberg (1986) described the isotropic component
of model-prediction error covariances using finite sums of Bessel functions, taking into
account the spherical harmonic expansion of the field. We refrain from using their representation of spatial correlations, as it yields large variations of the correlation model for
baselines where no data are available (Stoew, 2001).
5
Spatial correlations based on NWP and GPS data
It is possible to separate the contributions from two estimation error fields, assuming that
the spatial correlation structure is known for one of them. We denote the ZTD true values
as ℓi and ℓj for the locations i and j; the ZTD estimated using the HIRLAM and GPS are
G
H
H
G
G
ℓH
i and ℓi , respectively. The differences ei = ℓi − ℓi and ei = ℓi − ℓi represent the time
series of the estimation errors for the two methods. The covariance of the differenced time
series (the innovations) carries information about the spatial structure of the e H and eG
fields:
i
h
G
H
G
−
ℓ
cov ℓH
−
ℓ
,
ℓ
i
i
j
j
h
i
G
H
G
= cov ℓH
−
ℓ
+
ℓ
−
ℓ
,
ℓ
−
ℓ
+
ℓ
−
ℓ
i
i
j
j
i
i
j
j
h
i
G
H
G
= cov eH
i − ei , ej − ej
h
i
h
i
H
G
G
= cov eH
,
e
+
cov
e
,
e
(22)
i
j
i j ,
where in the last line we use the fact that there is no correlation between the HIRLAM
and the GPS estimation errors in ZTD (no GPS data are assimilated into HIRLAM).
H
An analytical model of the term cov[eH
i , ej ], derived using the NMC method (Parrish
and Derber, 1992), is used to fit a curve to the remaining term in order to obtain the
horizontal correlations of the GPS errors.
The NMC method was applied at SMHI for three different set-ups of HIRLAM, namely
for the 44-, 33-, and 22-km horizontal resolution. The resulting models of the HIRLAM
prediction error spatial correlations are presented in Fig.3. The spatial correlations depend on the chosen HIRLAM resolution, and as the model resolution increases, the correlation spatial scales become smaller (i.e., the correlation curves become steeper). In the
same plot we present also the spatial correlations derived using the ensemble assimilation method for a 44-km set-up of HIRLAM. The parameter values for the corresponding
models, ̺H(r), are summarized in Table 1.
The spatial correlations for the innovations ℓH − ℓG are shown in Fig. 4. These correlations can be expressed by the normalized pairwise covariances:
h
i
G
H
G
cov ℓH
i − ℓi , ℓj − ℓj
̺i,j(r) =
,
(23)
σi σj
7
Analytic function fit (NMC/44km)
Analytic function fit (NMC/33km)
Analytic function fit (NMC/22km)
Analytic function fit (ENS/44km)
0.9
0.8
Correlation coefficients
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−0.1
0
200
400
600
800
1000
Baseline length, km
1200
1400
1600
1800
Fig. 3. Analytic functions modelling the horizontal correlations of forecast errors for ZTD calculated for HIRLAM.
The solid line presents the model of the forecast error correlations, estimated using the NMC method for a 44-km
resolution NWP model. For the 33-km and the 22-km resolution model, the NMC method reveal even smaller
correlation scales (# and ▽, respectively). The results from the use of the ensemble assimilation method applied
to the 44-km resolution version of HIRLAM are also shown (×).
where r is the site separation; σi , σj are the standard deviations of the innovation time
series at sites i and j, respectively. They are of the order of 1.1–1.3 cm, and decrease
slightly with site latitude and height. Note that because spatially uncorrelated estimation/forecast errors contribute to the standard deviation of the innovations, ̺ i,j(0) ≤ 1 .
A composite model, based on a linear combination of three von Kármán curves (19),
Table 1. Parameters of the composite von Kármán correlation models (see text). The values of the parameters ν
and r2 for the ̺H(r) component were estimated using the NMC method at a horizontal resolution of 44, 33, and
22 km. The ensemble assimilation method was used to derive the parameters in the fourth group in the table.
Model
ai
Parameter
ν
r2
̺H(r)
̺G
1 (r)
̺G
2 (r)
0.27
0.53
0.07
1.05
3.33
2.98
192.80
37.02
685.33
NMC, 44-km resolution
̺H(r)
̺G
1 (r)
̺G
2 (r)
0.37
0.40
0.14
1.12
2.65
1.01
135.76
38.01
714.69
NMC, 33-km resolution
̺H(r)
̺G
1 (r)
̺G
2 (r)
0.36
0.42
0.13
0.90
2.21
1.01
127.10
42.31
726.27
NMC, 22-km resolution
̺H(r)
̺G
1 (r)
̺G
2 (r)
0.33
0.41
0.16
1.82
2.68
0.80
84.93
36.23
717.73
Ensemble assimilation,
44-km resolution
8
was fitted to the correlation values:
̺(r) = a1 ̺G
(r) + a2 ̺G
(r) + a3 ̺H(r) ,
1
2
(24)
where a1 , a2 , a3 are scaling factors; ̺H(r) represents the signature of the HIRLAM ZTD
forecast error correlations, and is plotted using the dash-dotted line in Fig. 4. The main
impact of the HIRLAM forecast errors is for baselines r < 1000 km. The other two
von Kármán curves, ̺G
(r) and ̺G
(r), model the contribution of GPS errors to the inno1
2
vation correlations (the solid line in Fig. 4). The composite model (24) is presented by
the dashed line, and has 7 degrees of freedom: three for each ̺G
(r) and ̺G
(r), and one
1
2
H
for ̺ (r). The parameter values for the model are shown in Table 1.
The slowly declining spatial correlations for baselines r > 400 km in the same graph
would have to be attributed to large scale spatial correlations in GPS estimation errors.
This behavior is in agreement with earlier simulation results (Jarlemark et al., 2001).
However, for small separations, the results depend strongly on the correctness of the
model ̺H(r), and on the availability of ZTD data from nearby GPS sites.
For the short baseline separations, the ZTD error correlations derived using the NMC
method decline at a slower rate compared to the HIRLAM–GPS difference data. This
result indicates that, in reality, the forecast error correlations may have shorter length
scales than those estimated by Berre (1997). It is further supported by the notion that
these correlations are calculated for the +6 hour forecasts, utilizing differences between
forecasts valid at the same time, but generated at different ranges (+24 and +48 hours
respectively). We speculate that the spatial scales of the forecast errors should decrease
if the forecast range is decreased. These scales also depend on the chosen horizontal
resolution of the NWP model, as demonstrated in Figures 4, 5, and 6.
A new method, used by the European Centre for Medium-Range Weather Forecasts,
based on ensemble assimilation indicates smaller length scales than those determined by
the NMC technique. The newer ensemble assimilation approach uses sets of randomly
perturbed observations and parameters of the numerical model to derive the error correlation structures in the forecasts (Houtekamer et al., 1996). The results from the ensemble
assimilation approach are shown in Fig. 7.
9
Horizontal correlations of the HIRLAM−GPS ZTD innovations
Innovation correlations
Model fit (innovations)
Model fit (Hirlam forecasts)
Model fit (GPS estimates)
0.9
0.8
Correlation coefficients
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−0.1
0
200
400
600
800
1000
Site separation, km
1200
1400
1600
1800
Fig. 4. Horizontal correlations of the HIRLAM–GPS zenith total delay innovations for a one year period (∗).
The dash-dotted line represents the contribution of the forecast errors correlations for ZTD, estimated using the
NMC method for HIRLAM at a resolution of 44 km. The dashed line shows the best fit to the innovation data of a
composite model of the correlations. The solid line represents the remaining term, corresponding to the correlations
of GPS errors.
Horizontal correlations of the HIRLAM−GPS ZTD innovations
Innovation correlations
Model fit (innovations)
Model fit (GPS estimates)
Model fit (Hirlam forecasts)
0.9
0.8
Correlation coefficients
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−0.1
0
200
400
600
800
1000
Site separation, km
1200
1400
1600
1800
Fig. 5. Horizontal correlations of the ZTD innovations for a one year period (∗). The dash-dotted line represents
the contribution of the forecast errors correlations for ZTD, estimated using the NMC method for HIRLAM at a
resolution of 33 km. The dashed line shows the best-fit composite model for the innovation correlations. The solid
line represents the remaining term, corresponding to the GPS error correlations.
10
Horizontal correlations of the HIRLAM−GPS ZTD innovations
Innovation correlations
Model fit (innovations)
Model fit (GPS estimates)
Model fit (Hirlam forecasts)
0.9
0.8
Correlation coefficients
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−0.1
0
200
400
600
800
1000
Site separation, km
1200
1400
1600
1800
Fig. 6. Horizontal correlations of the ZTD innovations for a one year period (∗). The dash-dotted line represents
the contribution of the forecast errors correlations for ZTD, estimated using the NMC method for HIRLAM at a
resolution of 22 km. The dashed line shows the best-fit composite model for the innovation correlations. The solid
line represents the remaining term, corresponding to the GPS error correlations.
Horizontal correlations of the HIRLAM−GPS ZTD innovations
Innovation correlations
Model fit (innovations)
Model fit (GPS estimates)
Model fit (Hirlam forecasts)
0.9
0.8
Correlation coefficients
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
−0.1
0
200
400
600
800
1000
Site separation, km
1200
1400
1600
1800
Fig. 7. Horizontal correlations of the ZTD innovations for a one year period (∗). The dash-dotted line represents
the contribution of the forecast errors correlations for ZTD, estimated using the ensemble assimilation method.
The dashed line shows the best fit to the innovation data of a composite model of the correlations. The solid line
represents the remaining term, corresponding to the correlations of GPS errors.
11
6
Spatial correlations using vertical coordinate residuals
It is well known that there is a strong correlation between the estimation errors in the
ZTD and the vertical coordinate. In order to assess the results in the previous section it
is worthwhile to study the spatial error correlations of the vertical coordinate. Vertical
coordinate residuals have been calculated using GPS data from 1996 to the beginning of
2001 from daily solutions (Johansson et al., 2002; Scherneck et al., 2003). The data have
been collected from the SWEPOS and FinnRef networks (see Fig. 1). We calculated the
cross-correlations of the time series for each pair of sites, and the results are shown in
Fig. 8.
The errors in the estimated station height have common sources with the errors in
the ZTD estimated using GPS. The plot in Fig. 8 supports the simulation results by
Jarlemark et al. (2001), suggesting a slow decrease of the horizontal correlations of the
GPS estimation errors. For the short baselines (∼200 km), the correlation values remain
below 0.9. This may indicate a random, uncorrelated component in the data, which has
a scaling effect on all correlation coefficient values.
Although the spatial correlations of these residuals are clearly present, the scatter of
the data points prevents us from providing a reliable fit of an analytical model of the type
discussed in Section 4.
1
0.9
0.8
Correlation coefficient
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0
200
400
600
800
1000
Site separation, km
1200
1400
1600
1800
Fig. 8. Horizontal correlations of the vertical coordinate residuals for the period 1996–2001, using data from
SWEPOS and FinnRef.
12
7
7.1
Temporal correlations of slant delay residuals using ray tracing
Hydrostatic and wet mapping functions
The slant propagation delay is usually converted to an equivalent zenith delay using
mapping functions (Davis et al., 1985). As an example, the wet-delay MF (WMF) m w (ε)
relates the equivalent zenith wet delay to the slant delay ℓw (ε):
ℓzw (ε) =
ℓw (ε)
,
mw (ε)
(25)
where ε is the elevation angle.
In GPS processing, the ZWD estimate is a weighted average of the equivalent zenith
delays ℓzw (εi ), estimated for the directions towards the satellites visible at a given epoch. It
is somewhat different from the actual value of the delay in the zenith direction, ℓ w (90◦ ),
due to (a) inhomogeneities of the atmosphere, and (b) imperfections in the mapping
function. The hydrostatic mapping functions (HMFs) and delays are related in a similar
way.
The ray-tracing software described below assumes the atmosphere to be horizontally
isotropic above the RS launch site. Most often, the total slant delays are expressed by
ℓ(ε) = ℓh mh (ε) + ℓw mw (ε) ,
(26)
where ℓh and ℓw denote the ZHD and the ZWD, respectively.
The GPS estimates utilize mapping functions which at present introduce an rms error
in total delay of less than 6 mm at satellite elevation angles larger than 10◦ (Bisnath et al.,
1997; Niell, 1996).
The Ifadis (1986) MFs are used for ε > 2◦ , and are in the form:
1
mI(ε) =
a
sin ε +
sin ε +
,
(27)
b
sin ε + c
The Ifadis HMF/WMF utilize surface measurements of the temperature, water vapour-,
and total pressure. Time series of the Ifadis HMF values for ε = 10◦ at the Landvetter
site are shown in Fig. 9.
The Niell (1996) HMF and WMF are based on statistics covering RS ascents from 26
locations and have the general form:
a
1+
1+
mN(ε) =
sin ε +
b
1+c
a
sin ε +
,
(28)
b
sin ε + c
where the parameters a, b, and c are obtained using only information about the site
location. These MFs are useful for ε > 3◦ . For the Niell WMF, these parameters depend
only on site latitude, while the behavior of the HMF is related also to the day of the year
and the station height (see Fig. 9). The differences between the parameterizations of the
Niell and the Ifadis mapping functions are summarized in Table 2.
13
5.564
Niell hydrostatic
Ifadis hydrostatic
5.562
Hydrostatic MF value at ε = 10°
5.56
5.558
5.556
5.554
5.552
5.55
5.548
5.546
0
50
100
150
200
250
Days since Jan 1, 2000
300
350
400
Fig. 9. Time series of the Niell and Ifadis hydrostatic mapping functions, calculated for ε = 10 ◦ at Landvetter. The
systematic differences between these two MFs can reach 0.1% in the winter for this site. The seasonal variations
in the Niell MF are only related to site latitude and the day of the year, while the Ifadis HMF utilizes surface
measurements.
7.2
Ray tracing using radiosonde ascents
The ray-tracing software uses vertical profiles of RH, air temperature and pressure to
calculate the wet and hydrostatic propagation delays, and the geometric delays (bending),
for a given elevation angle. The results from the ray tracing procedure are then used to
produce the normalized wet/hydrostatic atmosphere thickness (slant delays divided by
the corresponding zenith values) in order to obtain the residuals for the Niell and Ifadis
mapping functions, for ε = (5◦, 10◦, 15◦ ).
Ideally, the time correlations of these residuals should represent the temporal behavior
of the GPS slant delay estimation errors. The results for the Niell HMF/WMF are shown
in Fig. 10, spanning the year 2000. Note that the temporal correlations of the Niell WMF
residuals appear uncorrelated on scales larger than one day, while the HMF residuals are
correlated on scales of several days.
The results for the Ifadis HMF/WMF are shown in Fig. 11. Note the slowly-varying
bias in the HMF residuals (the HMF underestimates the ray-trace results in the winter);
this bias causes the correlations of the Ifadis HMF residual to be rather high on monthly
scales. Furthermore, there is an apparent diurnal signature in these correlations, which
should be attributed to the differences in the day-time and night-time variations of the
total pressure and temperature in the vertical direction.
The values of the temporal correlations for these structure function residuals change
insignificantly, only by about 1% for the different elevation angles.
Table 2. Summary of the properties of the Niell and the Ifadis mapping functions.
Hydrostatic MF
Wet MF
Niell
Depends on height, latitude and day of the year.
Depends only on station
latitude.
Ifadis
Parameterized using surface temperature, pressure
and humidity.
14
Slant delay residuals, m
0.02
Niell wet
Niell hydrostatic
a)
0.01
0
−0.01
−0.02
0
50
100
150
200
250
Days since Jan 1, 2000
300
350
400
1
Correlation coefficient
b)
0.8
0.6
0.4
0.2
0
−0.2
0
5
10
15
Time lag τ, days
20
25
30
Fig. 10. Residuals of the Niell HMF, calculated for ε = 10◦ at Landvetter (upper graph). The rms variations of
the HMF residuals are about 5 times larger than those of the Niell WMF. The temporal correlations of the Niell
mapping function residuals are shown in the lower graph. The WMF residuals (#) are practically uncorrelated on
scales larger than one day, while the Niell HMF residuals (×) are correlated on scales of several days.
Slant delay residuals, m
0.02
Ifadis wet
Ifadis hydrostatic
a)
0.01
0
−0.01
−0.02
0
50
100
150
200
250
Days since Jan 1, 2000
300
350
400
1
Correlation coefficient
b)
0.8
0.6
0.4
0.2
0
−0.2
0
5
10
15
Time lag τ, days
20
25
30
Fig. 11. Residuals of the Ifadis HMF, calculated for ε = 10◦ at Landvetter (upper graph). The rms variations of
the HMF residuals are about 6 times larger than those of the Ifadis WMF; there is also an apparent systematic
bias in the HMF for this site. The temporal correlations of the Ifadis MF residuals are shown in the lower graph.
The Ifadis WMF residuals (#) seem to be uncorrelated for time scales larger than one day. The Ifadis HMF
residuals (×) are correlated on scales of several days and exhibit a diurnal pattern.
15
8
Temporal correlations of ZTD from GPS position residuals
We calculated an estimate of the normalized auto-covariance function for a 1-year long
time series of differences between ZWDs, estimated using a WVR and GPS data at Onsala.
The results are shown in Fig. 12. These correlation values decrease to e −1 at a time lag of
about 1 day. The temporal structure of the WVR errors is unknown; therefore, the GPS
error temporal correlations cannot be deduced from the GPS–WVR differences without
reservations.
Correlation coefficient
1
0.8
0.6
0.4
−1
e
0.2
0
0
5
10
15
Time lag τ, days
20
25
30
Fig. 12. The temporal correlation values for the differences between GPS and WVR data at Onsala, for the year
2001. The time correlation of these differences decreases to about 1/e at a time lag τ ≈ 1 day.
Another way to analyze the temporal correlations of the errors in GPS-estimated
ZWD is to study the differences between time series of integrated water vapour derived
radiosonde- and GPS measurements. Although the errors in the RS data can be assumed temporally uncorrelated, the variance of the RS errors actually makes the timecorrelations of the differences appear smaller (see Fig. 13). Furthermore, there are few
colocated RS and GPS sites, and using similar differenced time series from nearby sites
would introduce a component due to the atmosphere itself.
Corr. coefficient
1
0.8
0.6
0.4
−1
e
0.2
0
0
5
10
15
Time lag τ, days
20
25
30
Fig. 13. The temporal correlation values for the differences between GPS and RS data at Sundsvall, for the year
1997. The time correlation of these differences decreases to about 1/e at a time lag τ ≈ 1/2 day, partly due to the
variance of the random measurement errors in the RS data.
The temporal behavior of the errors in GPS-estimated ZWD/ZTD can be described
indirectly, using the residual of the vertical coordinate of a permanent GPS site, as these
two parameters have common error sources. Their relation further depends on the data
processing strategy, the used zenith delay MFs, the antenna phase center modelling, and
the modelling of the ocean/atmosphere loading effects.
Time series of the vertical coordinate residuals, discussed earlier, are shown in Fig. 14,
together with the temporal correlations for the Onsala GPS site. The decrease in the
correlations of these residuals is consistently slow for all studied SWEPOS GPS sites; the
decorrelation time is about 1–2 days, and is summarized in Table 3. The variances of
the individual time series were found to be less than 100 mm2 . Apparently, the northern16
Vertical coordinate residuals, m
0.04
a)
0.02
0
−0.02
−0.04
−0.06
1996
1997
1998
1999
Calendar year
2000
2001
2002
Correlation coefficient
1
b)
0.8
0.6
0.4
−1
e
0.2
0
0
5
10
15
Time lag τ, days
20
25
30
Fig. 14. Time series of the vertical coordinate residuals for Onsala (upper graph), and the respective temporal
autocorrelation values (lower graph). The time correlation typically decreases to 1/e at a time lag τ ≈ 2 days.
most sites have slightly longer decorrelation times and slightly higher variances. This can
further be attributed to snow and ice, accumulating on the antenna radomes during the
winter (Jaldehag et al., 1996). Lastly, the data are generated at different rates (daily
coordinate estimates vs. 5-minute intervals for the ZTDs). These considerations imply
that the presented vertical residual correlations should be used as an approximation of
the ZWD/ZTD temporal error correlations cautiously, as the former over-estimate the
latter.
Correlation coefficient
1
0.8
0.6
0.4
−1
e−1
0.2
0
0
5
10
15
Time lag τ, days
20
25
30
Fig. 15. Temporal correlation values for the vertical coordinate residuals from 2001 at Onsala (#). The time
correlations of the GPS–WVR differences for 2001 are plotted for comparison (dashed line), and show consistent
correlation values.
Although the linear trends, annual and sub-annual oscillations (up to 3 cycles per year),
and air-pressure loading have been removed from these time series (Scherneck et al.,
2003), there are unmodelled components on time scales larger than one month, causing
the correlation estimates to remain rather high for time lags τ > 10 days for some GPS
sites, including the one at Onsala. In particular, the years 1996 and 2000 appear to have
short-term trends (see Fig. 14). Figure 15 presents the temporal correlations obtained
17
Table 3. Decorrelation of the residuals of the vertical coordinates for the SWEPOS sites. The data sets span the
years 1996–2001. The temporal correlation function is shown for time lag τ = 2 days. The GPS sites are sorted by
increasing latitude values.
Site
HASS
OSKA
ONSA
VISB
BORA
JONK
NORR
VANE
LOVO
KARL
MART
LEKS
SVEG
SUND
OSTE
UMEA
VILH
SKEL
OVER
ARJE
KIRU
Latitude
[◦ N]
Correlation ̺(τ )
at τ = 2 days
Variance
[mm2 ]
56.09
57.06
57.40
57.65
57.72
57.74
58.59
58.69
59.34
59.44
60.59
60.72
62.02
62.23
63.44
63.58
64.70
64.88
66.31
66.32
67.88
0.28
0.30
0.33
0.34
0.32
0.28
0.30
0.30
0.33
0.32
0.36
0.35
0.35
0.32
0.39
0.34
0.38
0.39
0.33
0.36
0.40
71
72
70
72
60
63
64
62
67
64
69
76
68
69
80
69
75
80
72
85
100
from a subset without an apparent trend in the vertical residual data for Onsala (the year
2001), together with a plot of the correlations of GPS–WVR differences from Fig. 12. We
also performed a study of smaller, seasonal subsets of these data and found no seasonal
dependence of the time correlation coefficients.
9
Discussion
The estimation of the atmospheric delay using GPS is associated with errors whose behavior is assumed to be random. The spatial structure of the error correlations in the ZTD
is derived most reliably using innovation time series. The spatial scales of the horizontal
correlations are of the order of several hundred kilometres, and the proposed analytical
models are described in Table 1 of Section 5. An implementation of the von Kármán
function for the purpose of a 3DVAR data assimilation should be simple.
Some contribution to the ZTD estimation errors can be systematic, due to modelling (imperfect mapping functions). The temporal correlations of the corresponding
errors have time-scales of days, and this is comparable to the time-scales of the atmospheric processes.
Furthermore, temporally correlated errors can be introduced into the ZTD estimates
by the GPS data processing strategy itself. In particular, a Kalman filter based data
processing software assumes the ZTD to be a random walk process with chosen parameters, which reflect the average behavior of the atmosphere. However, a rapid change in
the atmospheric state cannot be modelled well enough, resulting in short-term systematic differences. Such systematics have time scales of less than an hour, and are shorter
than the typical range of the numerical forecasts of today. These systematics will become
important in the future 4D data assimilation, when GPS data will be used at the precise acquisition epochs. A typical solution to the problem of time-correlated estimation
errors is “thinning” of the Kalman filter estimates, i.e. using results that are sufficiently
separated in time (Brown and Hwang, 1997).
The comparisons to WVR data in Section 8 suggest that the contribution of the atmosphere related systematic errors is small. The assumption that the daily vertical coordi18
nate residuals and the ZTD estimation errors have a similar time-correlation structure is
supported by the WVR-GPS comparison for the site at Onsala. This indicates that longterm, unmodelled trends in the satellite orbits and tidal motions of the sites propagate
into the ZTD estimation errors, with a variance of less than 100 mm2 .
Acknowledgements. The work on the spatial correlations of the zenith delay estimation errors was carried out in
close collaboration with Nils Gustafsson. We are grateful to Alan Rogers, Jim Davis, Tom Herring, and Arthur
Niell for providing us with the ray-tracing software.
References
Askne, J. and H. Nordius (1987). Estimation of tropospheric delay for microwaves from surface weather data.
Radio Sci., 22:379–386.
Berre, L. (1997). Non-separable structure functions for the HIRLAM 3DVAR. Tech. Rep. No. 30, Swed. Meteorol.
and Hydrol. Inst. Master Thesis.
Bevis, M., S. Businger, T.A. Herring, C. Rocken, R.A. Anthes and R.H. Ware (1992). GPS meteorology: Remote
sensing of atmospheric water vapor using the Global Positioning System. J. Geophys. Res., 97(D14):15,787–
15,801.
Bisnath, S.B., V.B. Mendes and R.B. Langley (1997). Effects of tropospheric mapping functions on space geodetic
data. In Proceedings of the IGS Analysis Center Workshop, Jet Propul. Lab., Pasadena, Calif., 12–14 March.
Blewitt, G. (1993). Advances in global positioning system technology for geodynamics investigations: 1978–1992.
In Smith, D.E. and D.L. Turcotte, eds., Contributions of Space Geodesy to Geodynamics: Technology, vol. 25,
pp. 195–213. AGU Geodynamics Series, Wiley & Sons, Inc.
Brown, R.G. and P. Hwang (1997). Introduction to Random Signals and Applied Kalman Filtering. John Wiley &
Sons.
Crane, R.K. (1976). Refraction effects in the neutral atmosphere. In Meeks, M.L., ed., Methods of Experimental
Physics, vol. 12B, pp. 186–200. Academic Press, N.Y.
Dai, A., J. Wang, R.H. Ware and T. van Hove (2002). Diurnal variation in water vapor over North America
and its implications for sampling errors in radiosonde humidity. J. Geophys. Res., 107(D10), 4090, 10.1029/
2001JD000642.
Daley, R. (1991). Atmospheric Data Analysis. Cambridge atmospheric and space science. Cambridge University
Press, Cambridge, U.K.
Daley, R. (1992). The effect of serially correlated observation and model error on atmospheric data assimilation.
Mon. Wea. Rev., 120(1):164–177.
Davis, J.L., T.A. Herring, I.I. Shapiro, A.E.E. Rogers and G. Elgered (1985). Geodesy by radio interferometry:
Effects of atmospheric modeling errors on estimates of baseline length. Radio Sci., 20(6):1593–1607.
Elgered, G. and P.O.J. Jarlemark (1998). Ground-based microwave radiometry and long term observations of
atmospheric water vapor. Radio Sci., 33(3):707–717, 98RS00488.
Emardson, T.R. and H.J.P. Derks (2000). On the relation between the wet delay and the integrated precipitable
water vapour in the European atmosphere. Meteorol. Appl., 7(1):61–68.
Emardson, T.R., G. Elgered and J.M. Johansson (1998). Three months of continuous monitoring of atmospheric
water vapor with a network of Global Positioning System receivers. J. Geophys. Res., 103(D2):1807–1820,
10.1029/97JD03015.
Erdélyi, A., W. Magnus, F. Oberhettinger and F.G. Tricomi (1954). Tables of Integral Transforms, vol. 2. McGrawHill Book Co.
Gustafsson, N., L. Berre, S. Hörnquist, X.Y. Huang, M. Lindskog, B. Navascués, K.S. Mogensen and S. Thorsteinsson (2001). Three-dimensional variational data assimilation for a limited area model. Part I: General formulation
and the background error constraint. Tellus, 53A(4):425–446.
Gutman, S.I., S.R. Sahm, S.G. Benjamin, B.E. Schwartz, K.L. Holub, J.Q. Stewart and T.L. Smith (2004). Rapid
retrieval and assimilation of ground based GPS precipitable water observations at the NOAA Forecast Systems
Laboratory: impact on weather forecasts. J. Met. Soc. Japan, 82(1B):351–360.
Hofmann-Wellenhoff, B., H. Lichtenegger and J. Collins (1997). Global Positioning System, Theory and Practice.
Springer, N.Y., 4th revised ed.
Hollingsworth, A. and P. Lönnberg (1986). The statistical structure of short-range forecast errors as determined
from radiosonde data. Part I: The wind field. Tellus, 38A(2):111–136.
Houtekamer, P.L., L. Lefaivre, J. Derome, H. Ritchie and H.L. Mitchell (1996). A system simulation approach to
ensemble prediction. Mon. Wea. Rev., 124(6):1225–1242.
Ifadis, I. (1986). The atmospheric delay of radio waves: Modeling the elevation dependence on a global scale. Tech.
Rep. No. 38L, Chalmers Univ. Tech., Göteborg, Sweden. ISBN: 99-0605353-4.
Jaldehag, R.T.K., J.M. Johansson, J.L. Davis and P. Elósegui (1996). Geodesy using the Swedish permanent GPS
network: Effects of snow accumulation on estimates of site positions. Geophys. Res. Lett., 23(13):1601–1604,
96GL00970.
Jarlemark, P., J. Johansson, B. Stoew, L. Gradinarsky and G. Elgered (2001). Spatial error correlation of GPS
atmospheres as determined from simulations. Phys. Chem. Earth, 26(6–8):451–456.
Johansson, J.M., J.L. Davis, H.G. Scherneck, G.A. Milne, M. Vermeer, J.X. Mitrovica, R.A. Bennett, B. Jonsson,
G. Elgered, P. Elósegui, H. Koivula, M. Poutanen, B.O. Rönnäng and I.I. Shapiro (2002). Continuous GPS
measurements of postglacial adjustment in Fennoscandia. 1. Geodetic results. J. Geophys. Res., 107(B8), 2157,
doi:10.1029/2001JB000400.
Källén, E. (1996). HIRLAM Documentation Manual. System 2.5. Tech. rep., SMHI, Norrköping, Sweden.
Lesht, B.M. (1999). Reanalysis of radiosonde data from the 1996 and 1997 water vapor intensive observation
periods: Application of the Vaisala RS-80H contamination correction algorithm to dual-sonde soundings. In
19
Proc. of the 9th Annual Atmospheric Radiation Measurements (ARM) Science Team Meeting, San Antonio,
Texas. [ Available online http://www.arm.gov/pub/docs/documents/technical/confp9903/lesht-99.pdf ].
Lindskog, M. (2001). Impact of observation processing on numerical weather prediction. Phil. Lic. Thesis, Dept.
Meteorology, Stockholm Univ.
Lindskog, M., N. Gustafsson, B. Navascués, K.S. Mogensen, X.Y. Huang, X. Yang, U. Andræ, L. Berre,
S. Thorsteinsson and J. Rantakokko (2001). Three-dimensional variational data assimilation for a limited area
model. Part II: Observation handling and assimilation experiments. Tellus, 53A:447–468.
List, R.J. (1958). Smithsonian Meteorological Tables. Smithsonian Inst., Washington, D.C.
Niell, A.E. (1996). Global mapping functions for the atmosphere delay at radio wavelengths. J. Geophys. Res.,
101(B2):3227–3246, 95JB03048.
Niell, A.E., A.J. Coster, F.S. Solheim, V.B. Mendes, P.C. Toor, R.B. Langley and C.A. Upham (2001). Comparison
of measurements of atmospheric wet delay by Radiosonde, Water Vapor Radiometer, GPS, and VLBI. J. Atmos.
Oceanic Technol., 18(6):830–850.
Parrish, D.F. and J.C. Derber (1992). The National Meteorological Centre’s spectral statistical interpolation
analysis system. Mon. Wea. Rev., 120(8):1747–1763.
Scherneck, H.G., J.M. Johansson, H. Koivula, T. van Dam and J.L. Davis (2003). Vertical crustal motion observed
in the BIFROST project. J. Geodynamics, 35:425–441.
Stoew, B. (2001). On the Temporal and Spatial Structure of the Estimation Errors in GPS Meteorology. Tech.
rep. No. 411L, Dept. of Radio and Space Science, Centre for Astrophysics and Space Science, Chalmers Univ. of
Technol., Göteborg, Sweden.
Stoew, B., J.M. Johansson and G. Elgered (2001). An assessment of estimates of IWV from ground-based GPS
data. Meteorol. Atmos. Phys., 77(1–4):99–107.
Tatarskii, V.I. (1971). The Effects of the Turbulent Atmosphere on Wave Propagation. Israel Program for Scientific
Translations, Jerusalem.
Wade, C.G. and B. Schwartz (1993). Radiosonde humidity observations near saturation. 8th Symposium on
meteorological observations and instrumentation, Anaheim, California, Jan 17–22.
Wang, J., H.L. Cole, D.J. Carlson, E.R. Miller, K. Beierle, A. Paukkunen and T.K. Laine (2002). Corrections
of humidity measurement errors from the Vaisala RS-80 radiosonde – Application to TOGA-COARE data. J.
Atmos. Oceanic Technol., 19:981–1002.
Webb, F.H. and J.F. Zumberge (1993). An introduction to the GIPSY/OASIS-II. JPL Publ. D-11088, Jet Propul.
Lab., Pasadena, California.
Zumberge, J.F., M.B. Heflin, D.C. Jefferson, M.M. Watkins and F.H. Webb (1997). Precise point positioning for
the efficient and robust analysis of GPS data from large networks. J. Geophys. Res., 102(B3):5005–5017.
20