User Tools

Site Tools


public:vlbi_fundamantals:theoretical_delays

Theoretical Delays

In order to calculate observed minus computed values for the least-squares adjustment, several models need to be applied. At first, the station coordinates at the observation epoch have to be determined. Then the station coordinates are rotated from the terrestrial to the celestial system (see section about Earth orientation) where the computed delays between two stations forming a baseline and a radio source are built, taking into account relativistic corrections and applying troposphere delay, and other corrections. For all the modeling details the reader is referred to the Conventions of the International Earth Rotation and Reference Systems Service (Petit et al., 2010) and its online updates at http://tai.bipm.org/iers/convupdt/convupdt.html, as well as to special IVS Conventions such as those for the treatment of the thermal expansion of radio telescopes (Nothnagel, 2009).

Station Coordinates at the Time of Observation

At first, coordinates (valid at a reference epoch; e.g. J2000.0 = 1 January 2000 at 12 h Terrestrial Time TT) and velocities of a specific realization of the International Terrestrial Reference System (ITRS) are taken to determine the mean coordinates at the time of observation. It should be mentioned already here that typically these realizations are TT frames and that the coordinates are provided in a conventional tide-free system. Examples are the International Terrestrial Reference Frame 2008 (Altamimi et al., 2011) or specific VLBI realizations like the VTRF2008 (Böckmann et al., 2010), the VLBI contribution to ITRF2008.

Then several corrections are added to get closer to the true station coordinates at the observation epoch. These corrections include periodic and aperiodic deformations of the Earth's crust. The largest periodic corrections are for the solid Earth tides, ocean tide loading, and pole tide loading. Solid Earth tides show mainly diurnal and semi diurnal oscillations which cause vertical deformations in a range of $\pm$20 cm and horizontal displacements of about 30\% of the vertical effect (Mathews et al., 1997). More difficult to model is the loading by the water masses of ocean tides and currents (ocean loading), which amounts to as much as a decimeter on some coastal or island sites (e.g. \citeauthor{scherneck1991},~\citeyear{scherneck1991}, \citeauthor{schuh2000},~\citeyear{schuh2000}). Additionally, there is also a periodic deformation at the S1 (24 hours) and S2 (12 hours) periods caused by pressure tides due to thermal heating of the atmosphere. All these effects, which should be corrected at the observation level, are described in detail in the IERS Conventions 2010 \citep{iersconv2010} and its electronic updates. The analysis strategy is not so clear with aperiodic deformations, e.g. with non-tidal atmosphere (but also non-tidal ocean and hydrological) loading, although their significance in VLBI analysis has been shown repeatedly (e.g. Rabbel & Schuh, 1986, vanDam et al., 1994). So far, there is no general agreement within the international space geodesy community whether these corrections should be applied at the observation level: Arguments against the application of non-tidal atmosphere loading at the observation level are that there is no consensus model available, that the accuracy of existing models (Petrov & Boy, 2004) is still not sufficient, and that geophysicists are interested in station coordinate time series which would show the loading signals. On the other hand – and this holds in particular for VLBI with a small number of stations (6 to 8) taking part in typical 24 hour sessions – neglected a priori atmosphere loading is partly absorbed by no-net-rotation (NNR) and no-net-translation (NNT) conditions and does not show up in the estimated station coordinate time series (Böhm et al., 2009b). Furthermore, there is significant variation in non-tidal atmospheric loading corrections within 24 hours which would be neglected if the correction was applied at a later stage of data analysis.

Earth orientation

In the next step, we need to transform the station coordinates from the International Terrestrial Reference System (ITRS) into the Geocentric Celestial Reference System (GCRS) at the epoch of the observation $t$. The transformation matrix can be written as

$$ \left[GCRS\right]=Q(t)\cdot R(t)\cdot W(t)\cdot\left[ITRS\right], $$

where $Q(t)$, $R(t)$, and $W(t)$ are the transformation matrices arising from the motion of the celestial pole in the celestial reference system, from the rotation of the Earth around the axis associated with the pole, and from polar motion respectively (Petit & Luzum, 2010). Matrix $W$ ('wobble') includes as parameters the coordinates $x_p$ and $y_p$ (polar motion) of the Celestial Intermediate Pole (CIP) in the Earth-fixed frame and the correction angle $s'$ which locates the position of the Terrestrial Intermediate Origin (TIO) on the equator of the CIP{}. Terrestrial (TIO) and Celestial Intermediate Origin (CIO) realize reference meridians in the respective systems. These terms are part of the 'CIO-based' transformation concept following the Non-Rotating Origin (NRO), which replaces the former 'equinox-based' transformation. Matrix $R(t)$ is the Earth rotation matrix with the angle $\theta$ between TIO and CIO. The conventional relationship defining UT1 from the Earth rotation angle $\theta$ is given by Capitaine, 2000 as

$$ \theta(T_u)=2\pi\cdot(0.7790572732640+1.00273781191135448\cdot T_u), $$

where $T_u$ = (Julian UT1 date - 2451545.0), and UT1 = UTC + (UT1-UTC). The difference between Universal Time (UT1) and Universal Time Coordinated (UTC, which differs by a known integer number of SI-seconds from TAI, the International Atomic Time, realized as a weighted mean of signals provided by atomic clocks located all over the world) can be uniquely observed by VLBI. All satellite techniques like the Global Navigation Satellite Systems (GNSS) or Satellite Laser Ranging (SLR) – due to the direct dependence between UT1 and the right ascension of the ascending node of the satellite orbit – can only observe length-of-day, which is the negative time derivative of (UT1-UTC), but need external information about UT1-UTC every few days.

The precession/nutation matrix, denoted $Q(t)$, includes the rotations around the angles $X$ and $Y$ (which are the coordinates of the CIP in the celestial frame) and the correction angle $s$ which positions the CIO on the equator of the CIP. The CIP is the reference pole for space geodetic measurements, i.e. it defines the observed axis. This is a pure convention realized by an appropriate theory of precession and nutation as will be described below. The orientation of the CIP does not coincide with that of a physical axis like the rotation axis, the figure axis, or the angular momentum axis, but it can be related to all of them. By definition the CIP is an intermediate pole which divides the motion of the pole of the ITRS w.r.t. the GCRS into a celestial and a terrestrial part. The celestial part (precession and nutation, [$X, Y$]) includes all motions with periods >2 days observed in the celestial frame, and this corresponds to all frequencies between 0.5 (retrograde) and +0.5 (prograde) cycles per sidereal day in the GCRS. The terrestrial part (polar motion, [$x_p$, $y_p$]) includes all motions outside of the retrograde daily band in the ITRS, i.e. it includes frequencies below 1.5 and above 0.5 cycles per sidereal day in the ITRS. The largest part of the celestial motion of the CIP can be calculated with a conventional precession/nutation model. Presently the model IAU 2006/2000A is recommended by the IERS Conventions 2010. However, remaining unmodeled parts of the celestial motion can be observed with VLBI and are provided by the IERS as so-called celestial pole offsets [$dX$, $dY$]. These offsets stem from residual errors of the a priori precession/nutation model and the phenomenon of the Free Core Nutation (FCN) which is a resonance mode due to the deviation of the rotation axis of the mantle from the rotation axis of the core (Dehant & Mathews, 2009). This retrograde motion with a period of about 430 days in the GCRS and a varying amplitude of up to 0.3 mas ($\sim10$ mm on the Earth surface) (Herring et al., 2002) is not predictable and cannot be neglected if someone wants to achieve highest positioning accuracies. Thus, neither precession/nutation nor polar motion and UT1-UTC can be predicted accurately enough with models but have to be observed by space geodetic techniques. A combination of these estimates is provided by the IERS, e.g. in the IERS 05 C04 series (Bizouard and Gambis, 2009), which can be used by space geodetic techniques as a priori information. If not estimated from the observations, the standard pole coordinates to be used are those published by the IERS $(x, y)_{IERS}$ with additional terms to account for the diurnal and semi-diurnal variations caused by ocean tides $(\Delta x, \Delta y)_{ocean\_tides}$ (Englich et al., 2008) and for libration $(\Delta x, \Delta y)_{libration}$:

$$ (x_p,y_p)=(x,y)_{IERS}+(\Delta x,\Delta y)_{ocean\_tides}+(\Delta x,\Delta y)_{libration}. $$

Here $(\Delta x,\Delta y)_{libration}$ are the forced variations in pole coordinates corresponding to motions with periods less than two days in space that are not part of the IAU 2000A nutation model. The IERS Earth Orientation Parameters (EOP) Product Center provides a subroutine to interpolate ('Lagrange' interpolation) in the $(x, y)_{IERS}$ pole coordinates which are typically released at midnight. However, this kind of interpolation (unlike linear interpolation) does not allow the rigorous estimation of estimated polar motion rates to be used for Earth rotation excitation studies. The situation is similar for the Earth rotation angle $\theta$ with models for the effects of ocean tides and libration that have to be added to the IERS 05 C04 values of UT1-UTC. In the case of the Earth rotation angle, tidal terms (with periods from 5 to 35 days) are usually removed before the Lagrange interpolation of the IERS values, and are restored afterwards.

General relativistic model for the VLBI time delay

The general relativistic model for the time delay is developed in the frame of the IAU Resolutions, i.e. general relativity using the Barycentric Celestial Reference System (BCRS) and the Geocentric Celestial Reference System (GCRS). The procedure to compute the VLBI time delay according to the so-called consensus model is taken from \citet{eubanks1991}, and it is summarized in the IERS Conventions 2010 (Petit & Luzum, 2010 , Chapter 11). The model has been developed for VLBI observations to extragalactic radio sources taken from the Earth surface, but not for observations to objects in our solar system like Earth- or Moon-orbiting satellites. In the model, which is accurate to the picosecond-level, it is assumed that the ionospheric delays have already been removed.

Theoretically, the VLBI time delays are measured in proper times of the station clocks, whereas the VLBI model is expressed in terms of coordinate time in a given reference system. We consider the VLBI time delay from the correlator to be equal to the Terrestrial Time (TT, agrees with SI second on the geoid and is equal to TAI, apart from a constant offset: TT = TAI + 32.184 sec) coordinate time interval $d_{TT}$ between the arrival of a radio signal at the reference point of the first station and the arrival of the same signal at the reference point of the second station. From a TT coordinate interval, $d_{TT}$, the Geocentric Coordinate Time (TCG) coordinate interval, $d_{TCG}$, can be determined by scaling: $d_{TCG} = d_{TT} /(1-L_G)$ with $L_G = 6.969290134\times 10^{-10}$ (Petit & Luzum, 2010).

The Terrestrial Reference System (TRS) space coordinates from the analysis of VLBI observations, $x_{VLBI}$, are termed 'consistent with TT' if derived from $d_{TT}$ intervals, and the TRS coordinates recommended by the IAU and IUGG resolutions, $x_{TCG}$, can be derived with $x_{TCG} = x_{VLBI}/(1-L_G)$ (Petit, 2000). Presently, all VLBI analysis centers provide their coordinate solutions as consistent with TT. Since also SLR analysis centers submit their solutions consistent with TT, it was decided that the coordinates should not be re-scaled to $x_{TCG}$ for the computation of ITRF2008 so that the scale of ITRF2008 (and earlier realizations) in this respect does not fully comply with IAU and IUGG resolutions (Petit & Luzum, 2010).

In the remaining part of this section, we provide the equations for the general relativistic VLBI time delay model and follow the IERS Conventions 2010 (Petit & Luzum, 2010). Although the delay to be calculated is by convention the time of arrival at station 2 minus the time of arrival at station 1, it is the time of arrival at station 1, $t_1$, that serves as the time reference for the measurement. Thus, all quantities are assumed to be calculated at $t_1$, including the effects of the troposphere. Assuming that the reference time is the UTC time of arrival of the VLBI signal at receiver 1, and that it is transformed to the appropriate timescale (e.g., UT1 or TT) to be used to compute each element of the geometric model, the following steps are carried out to calculate the VLBI time delay. First, the barycentric station vectors $X_i(t_1)$ for the receivers are determined with

$$ X_i(t_1)= X_\oplus(t_1)+ x_i(t_1) $$

where $t_1$ is the TCG time of arrival of the radio signal at the first receiver, $X_\oplus(t_1)$ is the barycentric radius vector of the geocenter, and $x_i(t_1)$ the GCRS radius vector of the i-th receiver. Then, we calculate the vectors $R_{iJ}$ from the Sun, the Moon, and each planet (with the barycentric coordinates $X_{J}$) except the Earth to receivers 1 and 2. The time $t_{1J}$ of the closest approach of the signal to the gravitating body $J$ can be determined with

$$ t_{1J}=min\left[t_1,t_1-\frac{K\cdot(X_J(t_1)-X_1(t_1))}{c}\right] $$

so that

$$ R_{1J}(t_1)=X_1(t_1)-X_{J}(t_{1J}) $$

and

$$ R_{2J}(t_1)=X_2(t_1)-\frac{V_\oplus}{c}(K\cdot{b})-X_{J}(t_{1J}) $$

with $V_\oplus$ the barycentric velocity of the geocenter, $b=x_2(t_1)-x_1(t_1)$ the GCRS baseline vector at the TCG time of arrival $t_1$, $K$ the unit vector from the barycenter to the source in the absence of gravitational or aberrational effects, and $c$ the velocity of light. The differential TCB gravitational delay for each of those bodies $J$ can then be calculated by

$$ \Delta T_{grav,J}=(1+\gamma)\frac{GM_J}{c^3}ln\frac{\left|R_{1J}\right|+K\cdot R_{1J}}{\left|R_{2J}\right|+K\cdot R_{2J}}, $$

where $M_J$ is the rest mass of the $J^{th}$ gravitating body and $G$ is the gravitational constant (Petit & Luzum, 2010). According to the General Theory of Relativity (GRT) the so-called light deflection parameter $\gamma$ is usually set to unity, but $\gamma$ can also be estimated as a solve-for parameter in a VLBI global solution comprising many or all VLBI sessions that have ever been carried out. See (Lambert & Le Poncin-Lafitte, 2009), (Heinkelmann & Schuh, 2010), or (Lambert & Le Poncin-Lafitte, 2011) for further details.

Analogously, we determine the TCB gravitational delay due to the Earth with

$$ \Delta T_{grav\oplus}=(1+\gamma)\frac{GM_\oplus}{c^3}ln\frac{\left| x_{1}\right|+K\cdot x_{1}}{\left|x_{2}\right|+K\cdot x_{2}}, $$

where $M_\oplus$ is the rest mass of the Earth, and we sum up the effects of all gravitating bodies to find the total differential TCB gravitational delay with

$$ \Delta T_{grav}=\sum_{J}\Delta T_{grav,J}. $$

We need to consider the Sun, Earth, Jupiter, the Earth's moon, and the other planets. If observations pass close to them, the major satellites of Jupiter, Saturn, and Neptune should also be added. If the ray path passes very close to some massive bodies, extra terms need to be included for accuracies better than 1 ps (see Klioner, 1991). For observations made very close to the Sun, higher order relativistic time delay effects become increasingly important. The largest correction is due to the change in delay caused by the bending of the ray path by the gravitating body $J$ described in (Richter & Matzner, 1983) and (Bibliography). The correction is

$$ \delta T_{grav,J}=\frac{4G^2M^2_J}{c^5}\frac{b\cdot(N_{1J}+K)}{(\left|R\right|_{1J}+R_{1J}\cdot K)^2}. $$

which should be added to $\Delta T_{grav}$ in the Equ. from above. $N_{1J}$ is the unit vector from the $J^{th}$ gravitating body to the first receiver, and $M_J$ is the rest mass of the $J^{th}$ gravitating body. Next, we compute the vacuum delay between $t_{vi}$, which are the 'vacuum' TCG times of arrival of a radio signal at the $i^{th}$ VLBI receiver including the gravitational delay but neglecting the troposphere propagation delay and the change in the geometric delay caused by the existence of the troposphere propagation delay:

$$ t_{v2}-t_{v1}=\frac{\Delta T_{grav}-\frac{K\cdot b}{c}\left[1-\frac{(1+\gamma)\cdot U}{c^2}-\frac{\left|V_\oplus\right|^2}{2c^2}-\frac{V_\oplus\\omega_2}{c^2}\right]-\frac{V_\oplus}{c^2}(1+K\cdot V_\oplus/2c)}{1+\frac{K\cdot(V_\oplus+\omega)}{c}}. $$

In the Equ. above $\omega_i$ is the geocentric velocity of the $i^{th}$ receiver. The aberrated radio source vectors $k_i$ (the unit vector from the $i^{th}$ station to the radio source after aberration) for use in the determination of the troposphere propagation delays are calculated with

$$ k_i=K+\frac{V_\oplus+\omega_i}{c} - K \frac{ K\cdot(V_\oplus+\omega_i)}{c}. $$

Thus, we can add the geometric part of the troposphere propagation delay to the vacuum delay with

$$ t_{g2}-t_{g1}=t_{v2}-t_{v1}+\Delta L_1\cdot\frac{K\cdot(\omega_2-\omega_1)}{c}. $$

where $\Delta L_i$ is the troposphere propagation TCG delay for the $i^{th}$ receiver ($=t_i-t_{gi}$). The total delay can be found by adding the best estimate of the troposphere propagation delay:

$$ t_2-t_1=t_{g2}-t_{g1}+(\Delta L_2-\Delta L_1). $$

The troposphere propagation delays in the Equations above need not be from the same model. The estimate in the Equ. from directly above ($t_2-t_1$) should be as accurate as possible (see Troposphere delay modeling), while the $\Delta L_1$ model in the Equ. ($t_{g2}-t_{g1}$) from above need only be accurate to about 10 ns (Petit & Luzum, 2010). Sections Earth orientation and Relativistic delay model are symbolically summarized in the Fig. below.

Symbolic representation of the derivation of the time delay $\tau$ in the TT frame starting with station coordinates in the ITRS and source coordinates in the ICRS (by courtesy of Lucia Plank)

Troposphere delay modeling

The troposphere path delay $\Delta L(e)$ at the elevation angle $e$ is usually represented as the product of the zenith delay $\Delta L^z$ and an elevation-dependent mapping function $mf(e)$:

$ \Delta L(e)=\Delta L^z\cdot mf(e). $

This concept is not only used to determine a priori slant delays for the observations, but the mapping function is also the partial derivative to estimate residual zenith delays, typically every 20 to 60 minutes. In the analysis of space geodetic observations, not only zenith delays are estimated, but also other parameters like stations heights and clocks. Whereas the partials for the clocks ($= 1$) and the station heights ($={\rm sin}(e)$) are exactly known, the partial derivative for the zenith delays (i.e., the mapping function) is only known with a limited accuracy. Via the correlations between zenith delays, station heights, and clocks, any imperfection of the mapping function is also manifested as station height error (and clock error) (Nothnagel et al., 2002). To reduce these correlations, observations at low elevations need to be included in the analysis; however, care has to be taken because mapping function errors increase rapidly at very low elevations, i.e. at $5^{\circ}$ or below. Presently, the best trade-off between reduced correlations and increasing mapping function errors is found for cutoff elevation angles of about $7^{\circ}$ (MacMillan & Ma, 1994, Teke et al., 2008) or by appropriate down-weighting (Gipson & MacMillan, 2009). Simulations of VLBI2010 observing scenarios indicate that in the future with faster antennas and significantly more observations the cutoff elevation angle can be increased so that the mapping function is less critical for the accuracy of VLBI analysis (Petrachenko, 2009).

Considering the Equ. from above we find the following relationship: If the erroneous mapping function was too large, the estimated zenith delay $\Delta L^z$ becomes too small, because the observed troposphere delay $\Delta L(e)$ stays the same. Consequently, the estimated station height moves up to account for the reduced zenith delay. MacMillan & Ma, 1994 set up a rule of thumb specifying that the error in the station height is approximately 0.22 of the delay error at the lowest elevation angle included in the analysis. Böhm (2004) confirmed this rule of thumb for VLBI analysis (and a cutoff elevation angle of $5^{\circ}$) specifying that the station height error is about one fifth of the delay error at $5^{\circ}$ elevation angle. The corresponding decrease of the zenith delay becomes about one half of the station height increase. Assuming azimuthal symmetry of the neutral atmosphere around the station (i.e., at a constant elevation angle the delay is not dependent on the azimuth of the observation), the approach as described in the Eq. below (e.g. Davis et al., 1985) is generally applied:

$ \Delta L(e)=\Delta L_h^z\cdot mf_h(e)+\Delta L_w^z\cdot mf_w(e), $

where $\Delta L(e)$ is the total path delay of the microwaves in the neutral atmosphere and $e$ is the elevation angle of the observation to the quasar (vacuum or geometric elevation angle). $\Delta L_h^z$ and $\Delta L_w^z$ are the a priori zenith hydrostatic and the estimated zenith wet delays, and $mf_h(e)$ and $mf_w(e)$ are the mapping functions which provide the ratio of the slant delay to the delay in zenith direction. The input to both mapping functions is the vacuum elevation angle $e$, because the bending effect is accounted for by the hydrostatic mapping function. The underlying continued fraction form to all mapping functions is that proposed by Herring (1992):

$ mf(e)=\frac{1+\frac{a}{1+\frac{b}{1+c}}}{{\rm sin}(e)+\frac{a}{{\rm sin}(e)+\frac{b}{{\rm sin}(e)+c}}}. $

Presently, the most accurate mapping functions globally available are the Vienna Mapping Functions 1 (VMF1, Böhm et al., 2006). Whereas the $b$ and $c$ coefficients in the Equ. above are provided as analytical functions depending on day of the year and station latitude, the $a$ coefficients (hydrostatic and wet) are provided as time series with a 6 hour time resolution on global ($2.0^{\circ}$ in latitude times $2.5^{\circ}$ in longitude) grids as well as for all VLBI sites for the complete history of VLBI observations. The coefficients are available from the website http://ggosatm.hg.tuwien.ac.at/ of the Vienna University of Technology as derived from operational analysis data as well as from forecast data of the European Centre for Medium-Range Weather Forecasts (ECMWF). Vienna Mapping Functions 1 coefficients from forecast data can be used for VLBI real-time applications without significant loss of accuracy as shown by Böhm et al. (2009b). Böhm et al. (2006) tested the concept of a 'total' Vienna Mapping Function 1, i.e. using the same function for mapping the a priori zenith delay to the elevation of the observation and for the estimation of the residual zenith delays. This concept was not as successful (in terms of baseline length repeatability, which is the standard deviation after removing a linear trend from a time series of baseline lengths) as the separation into a hydrostatic and a wet part, because the variation of the zenith wet delay is faster than can be described by coefficients with a 6 hour time resolution. A total mapping function (which is close to the hydrostatic mapping function) cannot account for this variation whereas the wet mapping function is able to account for it by estimating zenith wet delays every hour or even faster.

The a priori zenith hydrostatic delays in the Equ. for $\Delta L(e)$ above can be determined very accurately from the atmosphere pressure at the site (see Saastamoinen, 1972, Davis et al., 1985). If locally recorded pressure values are not available, it is recommended to take values retrieved from numerical weather models as e.g. provided by the Vienna University of Technology (http://ggosatm.hg.tuwien.ac.at/) together with the coefficients of the VMF1. If those are also not accessible, it is recommended to use an analytical expression like the Global Pressure and Temperature model (GPT, Böhm et al., 2007). Zenith wet delays are usually fully estimated in the VLBI analysis, although Gipson & MacMillan (2009) found a slight improvement if zenith wet delays were already added to the a priori delays and only residual zenith delays were estimated. Approximate values of zenith wet delays are provided with the VMF1 files. Errors in the zenith hydrostatic delays or the mapping functions have an influence on station heights as can be described with the rule of thumb by Böhm (2004) mentioned above: Exemplarily, the zenith hydrostatic and wet delays shall be 2000 mm and 200 mm, respectively, the minimum elevation angle is $5^{\circ}$, and the corresponding values for the hydrostatic and wet mapping functions are 10.15 ($mf_h(5^{\circ})$) and 10.75 ($mf_w(5^{\circ})$).

  1. We consider an error in the wet mapping function of 0.1 ($mf_w(5^{\circ}) = 10.85$ instead of 10.75) or in the hydrostatic mapping function of 0.01 ($mf_h(5^{\circ}) = 10.16$ instead of 10.15). The error at $5^{\circ}$ elevation in both cases is +20 mm, i.e. the error in the station height is approximately +4 mm.
  2. We assume an error in the total pressure measured at the station of -10 hPa. Let us assume that we take the 'mean' pressure from GPT (Böhm et al., 2007) although the real pressure at the site is larger by 10 hPa. -10 hPa correspond to about -20 mm zenith hydrostatic delay (Saastamoinen, 1972), which is then mapped with the wrong mapping function (factor $0.6 = 10.75 - 10.15$). At $5^{\circ}$ elevation the mapping function error causes +12 mm delay error, and one fifth of it, i.e. +2.4 mm, is the resulting station height error (Böhm et al., 2006).

On the other hand - due to atmosphere loading - the station height is decreased by about 4 mm (assuming a regression coefficient of 0.4 mm/hPa) if the actual pressure is larger than the mean pressure (from GPT) by 10 hPa. This implies that the application of GPT (or any other mean pressure model) for the determination of the a priori zenith hydrostatic delays of VLBI observables (or any other technique using microwave signals) partly corrects for atmosphere loading (Steigenberger, et al., 2009). In our example 2.4 mm out of 4.0 mm are compensated, because the application of GPT causes an error that goes into the same direction as the atmosphere loading corrections. In addition to the azimuthal symmetry of the troposphere delays according to the Eq. for $\Delta L(e)$ from above, we need to take azimuthal asymmetry into account, which in North-South direction is caused by the larger extension of the troposphere above the equator compared to the poles and which can also be due to local weather phenomena, e.g. if a VLBI site is located close to the coast. Typically, North and East gradients are estimated in VLBI analysis with a resolution of 2 to 24 hours. MacMillan (1995) proposed to use the equation

$ \Delta L(a,e)=\Delta L_0(e)+mf_h(e)\cdot {\rm cot}(e)\cdot\left[G_n {\rm cos}(a)+G_e {\rm sin}(a)\right] $

which goes back to Davis (1993) with $a$ denoting the azimuth and $\Delta L_0$ the symmetric delay. Chen & Herring (1997) published the formula

$ \Delta L(a,e)=\Delta L_0(e)+\frac{1}{{\rm sin}(e){\rm tan}(e)+C}\cdot\left[G_n {\rm cos}(a) + G_e {\rm sin}(a)\right] $

and recommended to use $C = 0.0032$ if the total (hydrostatic plus wet) gradients are estimated. Böhm & Schuh (2007) tested the application of a priori gradients derived from data of the ECMWF and found that estimating gradients from VLBI observations provides better baseline length repeatabilities than keeping those values fixed to non-zero a priori values derived from numerical weather models.

At present, many investigations are carried out on 'direct ray-tracing', i.e. deriving the slant delay for every single observation from numerical weather models (Hobiger et al., 2008). Gipson & MacMillan (2009) found better baseline length repeatabilities for CONT08 (a 14-days VLBI campaign in August 2008) with slant delays derived from data of the Goddard Modelling and Assimilation Office (GMAO) compared to using VMF1. Böhm et al. (2010) obtained improvement for UT1 estimates from Intensive sessions on the baseline Wettzell-Tsukuba if using ray-traced delays at station Tsukuba.

Antenna deformation

It has already been mentioned that in addition to the troposphere delays, there are other effects which depend on azimuth and elevation of the observations and need to be taken into account. For instance, deformations of the radio telescope structure which occur during the 24 hours of an observing session or between the observing sessions have to be considered if we want to achieve highest accuracies. They can be caused by snow and ice loading of the antenna (Haas, 1999) or thermal expansion of the radio telescopes (Nothnagel, 2009). Wresnik et al. (2007) compared thermal deformation measurements with invar rods at Wettzell (Germany) and Onsala (Sweden) to thermal deformation models applying local air temperatures as well as structure temperatures measured at those sites. They found an improvement if the air temperatures are first converted to structure temperatures before using them for modeling the thermal deformation of the radio telescopes (see Fig. below). Recent investigations also deal with gravitational deformations (Abbondanza & Sarti, 2010) which can significantly influence the estimated station heights if neglected in VLBI analysis (Sarti et al., 2011).

Measured vertical deformation (invar rod, upper curve), calculated antenna deformation using the air temperature to model the structure temperature (central curve), and the thermal deformation using the measured air temperature directly at Wettzell. The curves are offset by 2 mm for clarity (wresnik2007).

Axis offsets

At the radio telescope the distance between the feed horn and the axis intersection, which constitutes the baseline reference point, is assumed to be constant at the millimeter-level. In this case the corresponding time offset becomes part of the clock offset parameter. However, an axis offset model is applied to each antenna where the pointing axes do not intersect (Sovers, 1998). The axis offsets are provided by the IVS Analysis Coordinator together with the coefficients of the thermal expansion models (Nothnagel, 2009). Large radio telescopes such as the Effelsberg 100 m antenna exhibit elevation-dependent changes in the focal distance which can however be modeled to a level of a few millimeters (Rius et al., 1987). In the case of radio telescopes of the type alt-azimuth, the correction $\Delta\tau$ due to an axis offset $AO$ can be calculated with

$ \Delta \tau=-AO\cdot {\rm{sin}}(zd') $

if $zd'$ is the zenith distance corrected for refraction.

Source structure

A major problem is constituted by the fact that most of the observed radio sources tend to show structure at the level of a few milliarcseconds which often varies with time. These effects, in particular the changes in the source structure, pose a limit on the accuracy of the radio reference frame. For highest accuracies, a regular monitoring of the structure, which is also accomplished by analyzing VLBI data, can be done in parallel to the geodetic analysis, thus providing a means to correct for the source structure effects (see Collioud & Charlot, 2009, Schuh, 2000, and other references therein).

A few examples of constituents to the delay

The various effects that have been described in the previous sections on the calculated time delay $\tau$ for simulated observations with stations Westford (U.S.A) and Wettzell (Germany) to source 0642+449 are evident in the Fig. below. These effects are ordered by magnitude from top to bottom. The upper plot displays the total delay which is mainly caused by the Earth rotation and can be as large as one Earth radius. The next plot depicts the contribution of the hydrostatic delays at both sites to the total delay, i.e., the delay at Wettzell minus the delay at Westford. At each site, the hydrostatic delay at low elevation angles can be larger than 20 m, and it increases rapidly as the elevation angle decreases (Section Troposphere delay modeling). The third sub-plot illustrates the influence of the axis offsets at Wettzell and Westford (Section Axis offsets), of the gravitational effect of the Sun (Section General relativistic model for the VLBI time delay), and of the solid Earth tides (Section Station coordinates)). The bottom plot contains the influence of the high-frequency Earth rotation model accounting for the ocean tides (Section Earth orientation), as well as the effect of ocean loading (Section Station coordinates) at the sites on the calculated delay $\tau$.

Influence on the calculated time delay $\tau$ for simulated observations with stations Westford (U.S.A) and Wettzell (Germany) to radio source 0642+449. From top to bottom: Total delay (mainly caused by Earth rotation), troposphere, axis offsets, gravitational effect of the Sun, solid Earth tides, high-frequency model for the Earth rotation caused by ocean tides, and ocean loading (by courtesy of Lucia Plank)


| <= back |

public/vlbi_fundamantals/theoretical_delays.txt · Last modified: 2014/07/22 12:02 by admin