A multi-site intercomparison of integrated water vapour observations for climate change analysis

Water vapour plays a dominant role in the climate change debate. However, observing water vapour over a climatological time period in a consistent and homogeneous manner is challenging. On one hand, networks of groundbased instruments able to retrieve homogeneous integrated water vapour (IWV) data sets are being set up. Typical examples are Global Navigation Satellite System (GNSS) observation networks such as the International GNSS Service (IGS), with continuous GPS (Global Positioning System) observations spanning over the last 15 + years, and the AErosol RObotic NETwork (AERONET), providing long-term observations performed with standardized and well-calibrated sun photometers. On the other hand, satellite-based measurements of IWV already have a time span of over 10 years (e.g. AIRS) or are being merged to create long-term time series (e.g. GOME, SCIAMACHY, and GOME-2). This study performs an intercomparison of IWV measurements from satellite devices (in the visible, GOME/SCIAMACHY/GOME-2, and in the thermal infrared, AIRS), in situ measurements (radiosondes) and ground-based instruments (GPS, sun photometer), to assess their use in water vapour trends analysis. To this end, we selected 28 sites world-wide for which GPS observations can directly be compared with coincident satellite IWV observations, together with sun photometer and/or radiosonde measurements. The mean biases of the different techniques compared to the GPS estimates vary only between −0.3 to 0.5 mm of IWV. Nevertheless these small biases are accompanied by large standard deviations (SD), especially for the satellite instruments. In particular, we analysed the impact of clouds on the IWV agreement. The influence of specific issues for each instrument on the intercomparison is also investigated (e.g. the distance between the satellite ground pixel centre and the co-located ground-based station, the satellite scan angle, daytime/nighttime differences). Furthermore, we checked if the properties of the IWV scatter plots between these different instruments are dependent on the geography and/or altitude of the station. For all considered instruments, the only dependency clearly detected is with latitude: the SD of the IWV observations with respect to the GPS IWV retrievals decreases with increasing latitude and decreasing mean IWV.


Introduction
In climate research, the role of water vapour can hardly be overestimated. First of all, it is the most important natural greenhouse gas as it contributes to about 60 % of the natural greenhouse effect for clear skies (Kiehl and Trenberth, 1997). It also strongly influences atmospheric dynamics and the hydrological cycle through latent heat transport and diabatic heating, and is, in particular, a source of clouds and precipitation, directly affecting the climate. Unfortunately, clouds are the greatest source of uncertainty in climate models. In a more direct sense, water vapour also provides the largest known feedback mechanism for amplifying climate change (Soden and Held, 2006). The tropospheric water vapour

R. Van Malderen et al.: A multi-site IWV intercomparison
content increases in close association with warming, as dictated -under the condition of a constant tropospheric relative humidity -by the Clausius-Clapeyron equation which states that the water holding capacity of the atmosphere goes up at about 7 % K −1 increase in temperature (Trenberth et al., 2005;Wentz et al., 2007). Mears et al. (2007) demonstrated that, on a global scale, both climate models and satellite observations indicate an increase of the total amount of water in the atmosphere over tropical oceans with a rate of 5-7 % K −1 surface warming. Based on measurements of the total column water vapour from the Global Ozone Monitoring Experiment (GOME) for the period 1996-mid 2003, Wagner et al. (2006) found an increase of about 8 % K −1 for the tropics, and a comparable or even larger increase for the whole globe (8 and 12 % K −1 for the monthly and yearly averages, respectively).
However, finding observational evidence for this relationship over other or smaller regions is complicated by the observational constraints of water vapour as it has a very high temporal and spatial variability, in contrast to other greenhouse gases such as carbon dioxide or methane. For instance, there is a large gradient in the volume mixing ratio of water vapour from the ground (approx. 10 000 ppm) to the tropopause (approx. 4 ppm). The measurement of tropospheric water vapour is therefore a demanding task (Pałm et al., 2010). No standard instrument can measure it everywhere accurately in three dimensions. Luckily, about 45-65 % of the total column of water vapour is included in the surface-850 hPa layer (Ross and Elliot, 1996). In this paper, we will concentrate on this integrated water vapour (IWV) amount, also commonly known as precipitable water vapour (PW or PWV) or total column water vapour (TCWV), and we will express the IWV in units of mm. It is defined as the liquid equivalent of the total water vapour contained in an air column from the Earth's surface to the top of the atmosphere (Wang et al., 2007). This primary atmospheric variable can be measured by different devices, among which are ground-based, satellite-based and in situ instrumental techniques. The most commonly used ground-based instruments for water vapour monitoring are Fourier transform infrared (FTIR) spectrometers, microwave radiometers, lidars (for "light detection and ranging"), sun photometers and Global Positioning System (GPS) remote sensing (see e.g. Kämpfer, 2013, for an overview of most of these instruments). Recently,  developed an algorithm to retrieve the atmospheric water vapour column from Multi-AXis Differential Optical Absorption Spectroscopy (MAX-DOAS) based on ground-based observations in the yellow and red spectral range. In situ measurements by radiosondes have been used regularly since the 40s to measure the water vapour in the atmosphere and have widely been used to determine IWV trends in the literature (e.g. Ross and Elliott, 2001;Durre et al., 2009;Van Malderen and De Backer, 2010;Mattar et al., 2011). On meteorological satellites, the total precipitable water vapour is sounded using visible, near infrared, thermal infrared, passive microwave, and radio-occultation techniques (Urban, 2013). Each of these techniques has their strengths and weaknesses, and need to be inter-compared carefully under different conditions if progress is to be made on understanding the water vapour distribution and its time variability. Buehler et al. (2012) concluded that a literature survey reveals that reported systematic differences between different techniques are study-dependent and show no overall consistent pattern. Further improving the absolute accuracy of IWV measurements and providing climate-quality time series therefore remain challenging problems.
This study aims (1) to evaluate the quality and the consistency between different instruments measuring the IWV and (2) to assess their use for water vapour time series analysis and climate trend detection. Although an IWV multiinstrument intercomparison has been the subject of at least a dozen of papers in the literature (see Buehler et al., 2012, for a number of these, and Hocke and Martine, 2013, for tables of intercomparison studies), the added value of our analysis lies in focusing on homogeneous data sets or data sets taking part in a homogenisation procedure. Moreover, several studies in the literature concentrate on a multi-sensor intercomparison at a single site (e.g. Pałm et al., 2010;Schneider et al., 2010;Buehler et al., 2012). In this paper, we will first focus on the mid-latitude site Brussels (Belgium), as a case study (see Sect. 4), but subsequently extend our analysis to a set of 28 Northern Hemisphere (NH) stations to investigate the geographical dependency of the differences between the IWV data sets. Partly cloud-free skies are needed for water vapour observations by numerous instruments. A secondary aim of our study is therefore to study the impact of this "observation bias" on the IWV comparisons (this paper) and on the resulting IWV trends (paper in preparation). Finally, we identified in the literature some specific instrumental issues regarding the IWV retrievals, e.g. daytime/nighttime differences, sensor type change, (limiting) distance and/or height difference for co-location, and satellite scan angle. The last goal of this paper is to investigate the influence of these different issues on the IWV agreement between the different techniques.
The paper is organised as follows. In Sect. 2, we discuss the different instruments and data sets used for monitoring the water vapour, with a special emphasis on the (reported) uncertainties of the IWV retrievals. The general methodology used for the intercomparison is covered in Sect. 3. This methodology is then applied to the station at Brussels, Belgium,in Sect. 4, and to our selection of 28 stations in Sect. 5. The last section, Sect. 6, is reserved for drawing the conclusions and discussing the outlook.

R. Van Malderen et al.: A multi-site IWV intercomparison 2489 2 Instruments and data sets
The selection of instruments retrieving IWV presented in this paper (Table 1) is far from being complete but is based on the length of the time series and the data quality and homogeneity, so that they might be used for climate change analysis. This study focuses over land and includes some promising ground-based devices. Microwave satellite IWV measurements over ocean are often used as reference for global trend analysis and climate model intercomparisons and play a prominent role in the EUMETSAT (European Organisation for the Exploitation of Meteorological Satellites) operational climate monitoring from space initiative that will blend different satellite instruments. An intercomparison study over ocean including these microwave measurements has already been described in Mieruch et al. (2010) and was therefore not considered here.

GPS
Since 1997, the IGS (International GNSS Service, Dow et al., 2009) has provided operationally high-precision tropospheric zenith total delay (ZTD) estimates based on GPS observations recorded by continuously operating stations of its global network. Using surface measurements of pressure and temperature, these ZTD values can be turned into IWV values (Bevis et al., 1992(Bevis et al., , 1994Rocken et al., 1995;Ware et al., 1997) and used for atmospheric research.
In this study, the reprocessed IGS troposphere products will be used Bar-Sever, 2009, 2010;IGS ACC Website, 2013). They are computed with the GIPSY 1 software using the Precise Point Positioning (PPP) technique (Zumberge et al., 1997) and consist of the ZTDs and north and east gradient components (every 5 min) for 400 GPS stations from 1995 until the end of 2007. From 2008 to April 2011, these reprocessing results are complemented with the operational IGS final troposphere products which are fully consistent with the IGS reprocessing results. Consequently, we have access to 15+ years of continuous, worldwide and spatially well distributed ZTD values. However, even using consistent data analysis, the GPS ZTD can be affected by inhomogeneities due to changes at the stations, e.g. GPS equipment and/or operating procedures (e.g. elevation cut-off angle). Vey et al. (2009) concluded that only one third of 62 IGS stations -covering a period of at least 7 years with data gaps smaller than 3 months -could be assumed to be homogeneous.
After April 2011, IGS products are based on the new IGS08 terrestrial reference frame (Rebischung et al., 2012) and the new igs08.atx antenna models  and the consistency is not guaranteed anymore (Ray, 2011). This problem will be overcome when the next IGS back processing will take place 2 .

ZTD to IWV conversion scheme
The ZTD represents the total delay induced by the neutral atmosphere on the GPS signal propagation in the zenith direction. It includes the delay due to the whole density of the neutral atmosphere (named the zenith hydrostatic delay, ZHD) and an additional delay induced by the water vapour (called the zenith wet delay, ZWD). For climate applications, ZTDs are usually converted to IWVs to remove the hydrostatic effect: in a first step, the ZHD is subtracted from the ZTD to obtain the ZWD (all in m), (1) Assuming the hydrostatic equilibrium of the neutral atmosphere, the ZHD at a given GPS station can be modelled using the surface pressure (P s in hPa) and an estimation of the mean gravity (g m ) of the atmospheric column (Saastamoinen, 1972;Davis et al., 1985;Elgered et al., 1991), where g m depends on the latitude ϕ (in degrees) and the height h of the station above the Earth ellipsoid (in m): In a second step, the ZWD is converted into IWV (in kg m −2 or mm) using the following relationship (Hogg et al., 1981;Askne and Nordius, 1987;Bevis et al., 1992Bevis et al., , 1994: with κ(T m ) a proportionality factor containing constants like the specific gas constant of water vapour, the density of liquid water, atmospheric refraction constants and varying with the water vapour weighted mean temperature of the atmosphere T m (Bevis et al., 1992(Bevis et al., , 1994. T m can be either calculated from vertical profile data provided by radiosondes or global reanalyses, or estimated from surface temperature (T s ) observations using a linear empirical relationship (e.g. Bevis et al., 1992), the so-called T m − T s relationship: According to Wang et al. (2005), the estimation of T m using this latter relationship underestimates T m in the tropics and subtropics by up to 6 K and overestimates T m in the mid and high latitudes by up to 5 K. ≈ 10 % c ≈ 15 % for clear sky d ≈ 5 % e * for high latitudes, the broad swath orbits from GOME-2 overlap, i.e. for GOME-2 there can be two overpasses per day with a time difference of about 1.5 h. a Deblonde et al. (2005); Wang et al. (2007); Vey et al. (2010). b Miloshevich et al. (2009);Smit et al. (2013); Wang et al. (2013). c Alexandrov et al. (2009). d EUMETSAT (2010. e Fetzer et al. (2003). see Vey et al., 2009) we need to rely on the surface pressure and temperature observations from the large database of the World Meteorological Organization (WMO). The maximum separation distance between the IGS and the WMO stations was restricted to 50 km. The height difference is also taken into account by adjusting the surface values from the synoptic stations to the height of the IGS stations. For the adjustment of the surface temperature, we assume a standard lapse rate of −6.5 K km −1 , typical for wet adiabatic conditions. The hydrostatic and ideal gas equations are used to adjust the surface pressure. Hagemann et al. (2003) and Vey et al. (2009) showed that the pressure observations from neighbouring WMO stations precisely estimate the pressure at the GPS stations if the synoptic station is located within 50 km of the GPS site and for altitude differences less than 100 m (see also Gutman et al., 2003). We will come back to this point in Sect. 3.2.

Uncertainty of the GPS-based IWV
Three error sources have to be considered to estimate the uncertainty of GPS-based IWVs: (1) the uncertainty of the ZTD estimations, (2) the uncertainty of the ZHD modelling (Eq. 2) and (3) the uncertainty of the conversion from ZWD to IWV (Eqs. 3 and 4). Of these three, the main error source is generally the ZTD uncertainty, due to modelling limitations in the GPS analysis.
An estimate of the ZTD uncertainty is given by the formal error of the ZTD estimates. The analysis of the formal error for IGS stations over the whole period (1995-April 2011) is given in Fig. 1. The formal error ranges from 0.8 to 10 mm. 95 % of the ZTDs have a formal error below 2.6 mm and 99.9 % of them have a formal error below 4.4 mm. This is somewhat better than in Byun and Bar-Sever (2009), who  -4,19 4,4-4,59 4,8-4,99 5,2-5,39 5,6-5,79 6-6,19 6,4-6,59 6,8-6,99 7,2-7,39 7,6-7,79 8-8,19 8,4-8,59 8,8-8,99 9,2-9,39 9,6-9,79 10-10,19 Percentage of Observations ZTD Formal Error (mm) analysed a subset of 30 globally distributed IGS stations during 2003 and found a formal error between 1.5 and 5 mm of ZTD. Of course, this formal error does not take into account systematic errors (e.g. in GPS orbits and clocks) and is therefore slightly underestimating the actual error. With a claimed accuracy of the IGS ZTD product of 4 mm, Deblonde et al. (2005) derived a corresponding error of 0.6 mm in IWV. The ZHD is obtained from the surface pressure based on the hydrostatic equilibrium hypothesis (see Eq. 2). Its accuracy depends thus mainly on errors in the pressure measurements, causing offsets in the IWV estimates of at most 0.8 mm when considering an error in the pressure of 2 hPa (Vey et al., 2010). Other surface pressure uncertainties reported in the literature are 1 hPa (Deblonde et al., 2005) and 1.65 hPa (Wang et al., 2007), corresponding respectively to 0.4 and 0.6 mm uncertainty in IWV. Only during extreme (and rare) events, the atmosphere departs from hydrostatic equilibrium and consequently, this assumption can lead to significant errors in the IWV estimation (Brenot et al., 2006).

IGS ZTD Formal Error Histogram
During the conversion from ZWD tot IWV, the uncertainty comes from the estimation of the proportionality factor κ(T m ) in Eq. (3), which depends primarily on the mean atmospheric temperature. The uncertainty of T m (calculated from the surface temperature as in Eq. 4) is estimated to be around 5 K (or 1.8 % for T m = 273 K), which corresponds to an IWV error of 0.07 to 0.72 mm for respectively a dry or moist atmosphere (Deblonde et al., 2005).
Summed up these three error sources, the total uncertainty of GPS-based IWV measurements is generally less than 2 mm. This estimation agrees with previous studies where uncertainties are obtained by comparison with water vapour radiometers and radiosondes (Kuo et al., 1993;Rocken et al., 1993Rocken et al., , 1995Businger et al., 1996;Duan et al., 1996;Tregoning et al., 1998).

CIMEL sun photometer
A CIMEL sun photometer measures the sun and sky transmittance at selected wavelengths (filters) centred between 340 and 1020 or 1640 nm in an automated mode. Its fieldof-view is 1.2 • . Direct sun measurements are typically performed every 15 min between sunrise and sunset, under clear sky conditions. A CIMEL sun photometer is able to retrieve the IWV and the aerosol properties (e.g. the aerosol optical depth, Ångström exponent, single scattering albedo) using a combination of spectral filters and azimuth/zenith viewing controlled by a microprocessor. Using the 940 nm channel, centred on the 946 nm water vapour absorption line, the water vapour transmittance can be determined after first accounting for aerosol effects (using the 675 and 870 nm channels). Using a power law parameterisation, a conversion of the slant optical depth to the slant IWV can be obtained (e.g. Bruegge et al., 1992;Schmid et al., 2001). Uncertainties in this parameterisation and the Langley plot regression (due to variable atmospheric water vapour amounts) as well as deficits in the filter characterisation are the leading error sources. Alexandrov et al. (2009) estimate an IWV precision of about 10 % for this technique.
The CIMEL instrument is the standard instrument used in the AERONET (AErosol RObotic NETwork) international network (http://aeronet.gsfc.nasa.gov). AERONET provides a long-term, continuous and readily accessible public domain database of aerosol optical, microphysical and radiative properties for aerosol research and characterisation, validation of satellite retrievals, and synergism with other databases. The network is currently dense in Europe and the Americas. The data used in this paper is processed and archived with the version 2 of the AERONET IWV retrieval algorithm, defined as the quality-assured level (Holben et al., 2006). All the instruments in the AERONET are more or less annually calibrated with the Mauna Loa Observatory as the world standard reference.
The IWV retrieval with this instrument is limited to daytime and clear skies. This introduces a negative fair weather bias in the recorded IWVs, since cloudy conditions are often associated with higher IWV values. The reported IWV value in the AERONET database is the zenith value, but the solar slant value can easily be returned using the optical air mass table based on a standard atmosphere (Kasten and Young, 1989).

Radiosondes
Radiosondes (RS), launched on weather balloons for decades throughout the world, provide vertical profiles of pressure, temperature, relative humidity (or dew point temperature), wind speed and wind direction. The primary source of radiosonde data used in this paper was the Integrated Global Radiosonde Archive (IGRA, Durre et al., 2006). The IWV can be calculated from the IGRA sounding data if the pressure, temperature, and dew point depression are available at the surface and all mandatory above-surface levels up to and including 500 hPa. The calculation involves the conversion of dew point depression to specific humidity at each level followed by the integration of the specific humidity over all available levels between the surface and 500 hPa (Durre et al., 2009). The upper limit of 500 hPa for the integration is chosen because of the decreasing sensitivity of the radiosonde's humidity sensor with decreasing temperature (hence increasing height), so that we avoid the use of these data of lower quality. By this approach, we neglect the small contribution of the atmospheric layers above 500 hPa to the total column water vapour. Wang et al. (2007) pointed out that this introduces a dry bias of 2.44 % in IWV.
Global radiosonde climatic records suffer from three types of errors: (1) a systematic observational error, (2) spatial and temporal inhomogeneity, and (3) diurnal and spatial sampling errors (Wang and Zhang, 2008). For humidity measurements with the most commonly used capacitive polymer sensors (manufactured by Vaisala), the systematic observational errors give in general rise to a dry bias in the humidity measurements of the (lower) troposphere. The main sensor limitations are (i) chemical contamination: non-water molecules (e.g. from packaging material) occupy binding sites in the sensor polymer (dry bias), (ii) mis-calibration of the sensor, (iii) time lag, (iv) sensitivity of the daytime humidity observations to solar radiation (dry bias), (v) sensor ageing (dry bias). On the other side, a wet bias might be generated when a thin ice layer is formed on the humidity sensor after the passage through a cloud. After correction of all identified systematic biases and time-lag effects (see e.g. Wang et al., 2013, for a description of possible correction algorithms), Vaisala radiosondes may measure relative humidity with a relative uncertainty of about ±(3-5) % at ambient temperatures above −20 • C. However, at lower temperatures the relative uncertainty is increasing to R. Van Malderen et al.: A multi-site IWV intercomparison ±(5-10) % for the Vaisala radiosondes, although some differences exist between the different types (Smit et al., 2013).
These sensor-dependent errors and biases, together with other observational changes, often introduce non-climatic changes or inhomogeneities in historical records of humidity from radiosonde measurements. For that reason, radiosonde records need to be homogenised (e.g. Dai et al., 2011;Zhao et al., 2012) before they can be used to estimate long-term humidity and IWV trends (Wang et al., 2013). Such homogenisation algorithms can be neighbour-based procedures (e.g. Durre et al., 2009;McCarthy et al., 2009) and hence also partially anticipate for the spatial inhomogeneity of radiosonde measurements. The use or development of such a homogenisation procedure lies out of the scope of the present paper and the used IGRA database of radiosonde measurements therefore does not constitute a homogeneous data set. We nevertheless incorporate these RS measurements in our intercomparison, as they have been widely used in such analyses in the literature so far (see e.g. Buehler et al., 2012) and they provide the longest record (up to 60 years) of upper-air temperatures and humidity. Moreover, they have near-global coverage and high vertical resolution, are launched under all weather conditions and these daily daytime and nighttime observations take place at a large number of sites (see also Table 1).
Independent IWV uncertainty estimates derived from radiosonde profiles have not been reported in the literature so far. Nevertheless, Wang and Zhang (2008) concluded from comparisons with ground-based GPS measurements that the averaged systematic error for capacitive sensors, or their mean dry bias, is equal to −1.2 mm (−6.8 %) with a random error of 1.74 mm. Miloshevich et al. (2009) compared the IWV values calculated from the latest generation of Vaisala radiosondes (RS92) with microwave radiometer measurements and estimated a precision of 5 % for the IWV after applying an ultimate correction strategy.

GOME/SCIAMACHY/GOME-2
The "Global Ozone Monitoring Experiment" (GOME) on ERS-2, the "Scanning Imaging Absorption spectroMeter for Atmospheric CHartographY" (SCIAMACHY) on Envisat, and GOME-2 on MetOp-A, measure back-scattered light in the ultraviolet, visible and near-infrared parts of the spectrum. Total water vapour column amounts are retrieved from the spectral measurements in the visible wavelength range at 608-680 nm based on Differential Optical Absorption Spectroscopy (DOAS) providing a global data record since 1995 (GOME launch). Details on the MPI-C retrieval of water vapour columns can be found in Wagner et al. (2011).
The three instruments have different ground pixel sizes, see Table 1. All satellites are operated in a sun-synchronous orbit, with a constant equator crossing time (local) around 10:30 (GOME), 10:00 (SCIAMACHY) or 09:30 (GOME-2). Global coverage is provided after 3 (GOME), 6 (SCIAMACHY) and 1.5 (GOME-2) days. In the remainder of this paper, we will use the acronym GOMESCIA to denote the combined solution obtained from the three instruments.
Three main error sources contribute to the uncertainty of the GOMESCIA IWV retrievals (Wagner et al., 2011;EU-METSAT, 2010): (1) the errors of the spectroscopic data for the H 2 O and O 2 absorptions are estimated to about 5-10 %; (2) the uncertainty of the spectral retrieval as determined from the spectral residual is typically < 3 %; (3) the dominant error source is caused by uncertainties of the atmospheric radiative transfer, mainly due to effects of varying cloud cover and surface albedo. This error source is estimated to be about 15 % for clear sky observations and up to 100 % for individual observations in the presence of large cloud amounts. This range of uncertainty was largely confirmed by extensive validation exercises (Schröder and Bojkov, 2012) within the GLOBVAPOUR project (http://www.globvapour.info/index. html). Note that GOMESCIA measurements with the largest cloud contamination are filtered out using a threshold criterion for the simultaneously measured O 2 absorption (Wagner et al., 2006(Wagner et al., , 2011; see also Sect. 4.3. Recently, Schröder and Bojkov (2012) showed that systematic biases between the time series of the different sensors can occur, although they were merged using coincident observations during the respective overlap periods. This especially applies to the transition from SCIAMACHY to GOME-2 at the beginning of 2007. For observations over the continents, these IWV jumps can be up to 4 mm. Detailed investigations indicated that the jumps are caused by the differences in ground pixel size and swath width. Activities are ongoing to eliminate these jumps in the next version of the GOMESCIA data set by only using similar pixel sizes and viewing geometries.

AIRS
The Atmospheric Infrared Sounder (AIRS) is a cross-track scanning instrument aboard the polar orbiting Aqua satellite launched by NASA on May 2002 (Aumann et al., 2003). The scans from AIRS correspond to ellipsoidal foot prints with a major axis varying from 13.5 km (at nadir) to 31.5 km. A full global coverage is obtained in half a day, so AIRS provides two soundings per day at a given location on Earth between 45 • N and 45 • S. The instrument itself is a hyperspectral, scanning infrared sounder based on a grating spectrometer, that measures emitted radiation in the range from 3.7 to 15 µm, with a spectral resolution (λ/ λ) of 1200 (Aumann et al., 2003). Because of this high spectral resolution and excellent absolute accuracy, the AIRS retrieval algorithms are able to produce profiles of atmospheric temperature, moisture and trace gases, with a vertical resolution of a few kilometres throughout the troposphere. The AIRS IWV is obtained by integrating the vertical profile of water vapour mixing ratio retrieved from cloud-cleared radiances (Bedka et al., 2010, and references therein). The stated theoretical accuracy specification for the absolute AIRS total water vapour product is 5 % (Aumann et al., 2003). In a comparison with Vaisala operational radiosondes over sea/ocean, Fetzer et al. (2003) showed biases from −4 to 4 % for AIRS water vapour retrievals of the ∼ 2 km layer between the surface and 500 hPa. Also over (identical) ocean scenes, observational biases generally less than 5 % in IWV were found by Fetzer et al. (2006), who compared the IWV observations of AIRS and AMSR-E (Advanced Microwave Scanning Radiometer -EOS, on board Aqua). A validation over land of AIRS retrievals of IWVs (version 4) by measurements of more than 375 GPS receivers over the continental USA (from April to October 2004) is described in Rama Varma Raja et al. (2008). They concluded that the absolute biases between these two techniques range from 0.5 to 1.2 mm and root mean square differences from 3 to 4.5 mm, with consistently large monthly correlation coefficients, ranging from 0.91 to 0.98.
This study examines AIRS moisture products (level 2 version 5), which are available from September 2002 onwards (EOSDIS database; https://users.eosdis.nasa.gov). We took into account the flag Qual_H2O for quality control of the IWV product. This flag is based on the lowest good estimation of the pressure PBest, which depends on the cloud cover. The best data quality flag (Qual_H2O = 0) is assigned to the data if PBest equals the surface pressure, Qual_H2O = 1 for those data where PBest is lower than 300 hPa. If PBest is higher than 300 hPa, the value 2 is attributed to the flag Qual_H2O. In this paper, we will only deal with retrievals for which Qual_H2O = 0 or 1.

Methodology
As the IGS GPS network provides a long-term, world-wide (continental), homogeneously (re)processed and all-weather IWV database, it is used in this paper as the reference with which all other IWV data sets will be compared.

Co-location criteria and site selection
With a maximal horizontal distance of 30 km, we identified co-locations between IGS, AERONET sun photometer and IGRA radiosonde sites. At those selected sites, overpass measurements of GOME/SCIAMACHY/GOME-2 and AIRS are available. Because these different satellite devices have different ground pixel sizes, we apply different geometrical co-location criteria with the ground-based IGS stations: the best compromise with other criteria like the amount of co-locations and the presence of clouds was found for a maximum distance of 50 km for AIRS and the inclusion of the IGS station in the satellite ground pixel for GOME/SCIAMACHY/GOME-2. We will come back to this issue later in the paper, in Sects. 4.3 and 4.4 and Sects. 5.3 and 5.4. As addressed in Sect. 2.1.2, the WMO database of synoptic stations was consulted to find co-locations with IGS stations within a maximal horizontal distance of 50 km. With all these criteria, we end up with 28 locations where an IGS station is co-located with at least one CIMEL sun photometer or a radiosonde launch site. Figure 2 shows the locations of these 28 sites and the data availability for each instrument is given in Fig. 3. Unfortunately, our co-location sites are all Northern Hemisphere sites 3 (see Fig. 2), with the strongest Figure 3. IWV data availability over the last 15+ years for the different instruments at the selected sites. Note that, for high latitudes, no GOMESCIA observations are available in winter due to high SZA. concentration in Europe, but nevertheless representing more or less the whole longitude band. From Fig. 3 it should be noted that there is no single station in our selection at which all instruments are represented for the entire time period. Furthermore, it seems justified to use the GPS network as reference data set, as it provides a reasonable IWV time series for the majority of the selected stations.

Selection of synoptic stations for the ZTD-IWV conversion
To convert ZTDs to IWV values, coincident surface pressure and temperature observations are needed at the GPS site (see Sect. 2.1.1). As high-quality observations of these meteorological variables are rarely available at the GPS sites, we rely on WMO synoptic station measurements in the vicinity. For 13 out of the 28 selected IGS stations in our selection, more than one WMO station lies within the imposed maximum horizontal distance of 50 km, allowing different strategies to be applied for the ZTD to IWV conversion. For instance, Wang et al. (2007) propose to use the weighted average of the pressure observations (corrected for the altitude difference) of all selected synoptic stations to obtain the pressure at the GPS station height and location. However, different synoptic data sources of deviating quality and at sometimes diosonde and CIMEL sun photometer instruments is found in the Northern Hemisphere.
geographical distinct locations are then being mixed up. An alternative approach, used here because it better meets the aim of an "as purely observational as possible" intercomparison, consists of selecting the "best" WMO station for a given IGS station. To define the "best" WMO station, different criteria can be considered: (1) the data quality of the WMO station, i.e. no apparent inhomogeneities, small number of outliers, (2) the length and time frequency of the observation records, (3) the altitude difference with the GPS station, (4) the horizontal separation distance to the GPS station, and (5) the resulting GPS IWV correspondence with co-located IWV observations measured by other techniques. Therefore, we first analysed the impact of criteria as the height difference and horizontal separating distance on the calculated IWV values by creating, for a given IGS station, scatter plots between the IWV data sets calculated from the surface measurements provided by each of the co-located WMO SYNOP station. We will illustrate this procedure for the IGS station FFMJ. Three WMO stations are located within 50 km distance: 10532, 10635, and 10637, see Table 2. The altitude difference of almost 700 m between the stations 10635 and FFMJ is the largest value of all possible WMO-IGS station co-locations in this study and should therefore provide an upper limit. Of the three stations, the station 10637 has the smallest height difference and horizontal separating distance with FFMJ. So, we took this station as the reference surface data set for the comparison of the IWV values calculated with observations from those three WMO stations. The corresponding scatter plot properties are listed in Table 2, and illustrate that the large altitude difference between the IGS station and the WMO station 10635 cannot completely be corrected for, as a bias of 1.41 mm persists. On the other hand, the differences between the sets of FFMJ IWV values obtained from the surface measurements of stations 10532 and 10637 are really marginal. So, for this IGS station we selected the WMO station 10637 as surface data source. Additionally, we also analysed differences in scatter plot properties between GPS IWV values calculated from different co-located WMO stations on one hand, and IWV values observed with co-located instruments on the other hand. With this inter-technique analysis, it is possible to study the influence of the selection of the WMO station on the intercomparison. This is illustrated again for the IGS station FFMJ, which is now compared with the co-located CIMEL sun photometer at Mainz (located at a distance of 28 km from the IGS station and an altitude of only 20 m higher), see Table 3. The differences in scatter plot properties are of the same magnitude than the GPS IWV scatter plot differences and are almost insignificant for the IWV values generated from the WMO stations 10532 and 10637. In this case, the choice for the WMO station 10637, the closest station with minimum altitude difference, seems justified, as also the SD, correlation coefficient, and regression slope coefficient are slightly better.
To summarise, for a given IGS station, we chose the WMO surface station with minimum altitude difference and at minimum distance (in this order). If another WMO station has a much larger data record after removing spurious time periods and outliers, its surface data is used instead if the differences between the resulting IWV data sets are insignificant. For our selection of 28 IGS stations, the WMO surface stations finally selected, listed in Table 4, are located on average at a distance of 27.5 (±11.1) km and the average absolute height difference amounts to 35.7 (±42.2) m, with extreme values obtained from the IGS stations NISU and HLFX with the WMO stations located respectively 197 m below and 121 m above.

Correction for difference in altitude
An altitude difference between the different co-located ground-based or in situ instruments measuring IWV will introduce at least an artificial bias between the IWV data sets, because the device that is located at the lower altitude should logically measure a larger column of water vapour, and hence larger IWV values. Luckily, height differences between ground-based and in situ instruments at our selected sites are relatively small: the maximum height difference between the GPS antenna and the radiosonde launch site is at most 105 m (Munich, OBE2), and at most about 75 m (Munich, OBE2 and Paris, OPMT) between co-located GPS antenna and CIMEL sun photometers. For co-located IGS and radiosonde stations, an altitude correction is most easily applied to RS observations, because these of course provide vertical profile information. When the RS station is located below the corresponding IGS site, the vertical integration of RS specific humidity q is started at the pressure level corresponding to the IGS site altitude. This pressure level is computed from the heights of the IGS and RS sites based on the hydrostatic equilibrium hypothesis. At the starting pressure level, q is obtained by linear interpolation of the RS data. For RS stations located above IGS sites, RS surface data are extrapolated hydrostatically with an assumed temperature lapse rate of −6.5 K km −1 , while q is extrapolated by assuming that the dewpoint temperature depression (temperature minus dewpoint temperature) is constant with height, as explained by Deblonde et al. (2005). An alternative approach was proposed by Buehler et al. (2012) for the Arctic station Esrange, who found a relative bias IWV / IWV of −3.5 % per 100 m altitude difference. However, the actual scaling factors seem to depend on location

R. Van Malderen et al.: A multi-site IWV intercomparison
(for example, Bock et al., 2007, found a value of −4.0 % per 100 m for Africa) and might therefore not generally applicable to other IWV data. The altitude difference corrections applied here for the radiosonde IWV integrations improve the GPS-RS IWV comparisons, as the (mean) biases and SD are reduced, and higher (mean) scatter plot correlation coefficients are obtained. The only exception is that the regression line slope coefficients deviate more from unity after RS IWV correction.
Correcting the CIMEL IWV data for the altitude difference with a co-located GPS station is not as straightforward. The only possible altitude difference correction should then be applied to the GPS IWV data, using a correction applied to the GPS hydrostatic delay ( ZHD with a pressure difference applied to Eq. 2 instead of the surface pressure P s ). However, we cannot detect a real improvement in the coincident GPS-CIMEL IWV comparisons after adopting the altitude correction strategy to the GPS IWV retrievals. Furthermore, this correction introduces a dependence of the GPS-CIMEL IWV scatter plot properties on the GPS-CIMEL altitude differences. Finally, we prefer to have the same common reference GPS IWV data for the comparison with the data from the other instruments.

Coincidence criteria
In this paper, varying coincidence criteria are used to find a balance between the number of measurement matches and the different instruments measuring the same atmospheric conditions. For the ZTD to IWV conversion, "coincident" observations of surface temperature and pressure of a colocated WMO station were used. For this particular case, as the time frequency of the GPS-based ZTDs is very high (every 5 min), SYNOP observations at exactly the same time, i.e. within 1 min, were used. The time frequency of the resulting GPS IWV values is then mostly driven by the time frequency of the SYNOP observations, unless the GPS ZTD time series has significant gaps.
To create the IWV scatter plots, "coincident" measurements of two considered instruments are needed, for which we will not apply any time averaging or interpolation. Instead, all comparisons shown are point-by-point comparisons, which means that every IWV measurement of a given instrument (RS, CIMEL, GOMESCIA and AIRS) will be compared with the corresponding GPS IWV value, within a maximal time interval of 10 min for the CIMEL instrument, and 30 min for radiosondes, GOME(2)/SCIAMACHY and AIRS. When a measurement of an instrument coincides with several GPS IWV values within these time intervals, we use exclusively the coincidence with the minimal time difference. This strategy reflects the scope of the present paper to make an intercomparison based on purely observational data. We also compared hourly averages and the results obtained from that analysis are identical to the ones presented here. The difference between the applied time intervals for defining "coincident" data for the selected instruments stems from the different temporal resolution of these instruments (see Table 1).

Case study: a focus on Brussels, Belgium
For several reasons we first concentrate on the intercomparison at the IGS site BRUS located at Uccle, Brussels, Belgium (50 • 48 N, 4 • 21 E, 100 m a.s.l.). All ground-based instruments, as well as the weather station, are exactly colocated at the same site, so that we do not have to take into account any height difference nor separation distance between the different instruments. Also, several meteorological data (e.g. the cloud cover) are collected at Uccle providing additional information for the interpretation of IWV differences between different instruments. Finally, all metadata of each of the ground-based devices operated in Brussels are archived, so that we are aware of any instrumental change that might give rise to an inhomogeneity in the instrument's IWV time series.
An overview of all available IWV data sets in Brussels is shown in Fig. 4, also illustrating the presence of the wellestablished seasonal cycle in the IWV time series: maximum values are reached during the summers (June-August), when the surface temperatures are highest, and minimum values are attained during the winters (December-February). Since different instruments cover different observation periods, we selected a reference instrument against which the other techniques will be compared for the longest time interval possible, as already explained in Sect. 3. The choice of the GPS instrument as the reference for the Brussels station is justified by the following arguments: (1) BRUS ZTDs provided by IGS have a time frequency of 5 min, whereas the automatic weather station at Uccle reports measurements every 10 min, so that the resulting time frequency of GPS IWVs is 10 min, (2) data at this high temporal resolution are available since the launch of the automatic weather station at the end of 1999 (the IGS site in Uccle, BRUS, is operated continuously since 1993), which is sufficiently long compared to the other instruments, (3) the GPS ZTD time series only have minor data gaps, and (4) since GPS ZTDs are provided thanks to the reprocessing efforts carried out by the IGS, a homogeneously processed data set is available.

CIMEL sun photometer
The CIMEL sun photometer at Uccle, Brussels, was installed in July   Here, CIMEL IWVs are compared with IWVs of the IGS site "BRUS". The resulting scatter plot of more than 9000 coincident observations can be found in Fig. 5. A very good agreement is found between both techniques as shown by the high (above 0.99) linear Pearson correlation coefficient and an average bias of only 0.16 mm. Also, the SD has a very low value of only 0.96 mm, reflecting the negligible dispersion in the scatter plot. On the other hand, the slope of the linear regression line has a value significantly lower (0.913 ± 0.001) than the ideal one-to-one correlation. As can be derived from Fig. 5, the cause for this too low slope value is twofold: (1) at low IWV values (< 10 mm), the CIMEL IWVs are always higher than the corresponding IGS IWVs, and (2) at high IWV values (> 30 mm), the opposite is observed. The first cause can be attributed to the fact that, under dry conditions, the GPS data are known to be less precise (Wang et al., 2007). The second effect might reflect the weather observation bias present in the CIMEL data: the CIMEL instrument needs a clear sky in the direction of the sun. At high IWV values, there is a higher probability to have clouds. Water vapour  associated with these clouds might be captured by GPS signals, but never by the CIMEL sun photometer due to this weather observation bias, so that the CIMEL IWVs will always be lower than the corresponding IGS IWVs in these cases.

Influence of clouds
To test the impact of the cloud cover on our intercomparisons, we analysed the cloud meteorological data at the epochs of the coincident GPS-CIMEL IWV values. We end up with a total of about 1100 coincident IWV and cloud cover observations. The cloud data set consists of a classification of the cloud cover into eight "octas" classes, ranging from 0 (clear sky) to 8 (sky totally covered with clouds). To rule out the effect of the number of observations on the GPS-CIMEL IWV scatter plots, we combined several of these "octa" classes to have a representative number of observations for each "combined" class. Figure 6 shows the scatter plots and linear regression lines for IWV values under clear sky conditions, under very low cloud cover ("octas" equal to 1 or 2) and under moderate to heavy cloud cover ("octas" larger than 2). From this figure, it is obvious that the (mean) IWV values observed by both the CIMEL and the GPS increase with increasing cloud cover, also found by Gaffen and Elliott (1993) for radiosonde humidity data, especially over the continents in mid and high latitudes. Secondly, for increasing cloud cover, GPS is more frequently measuring higher IWV values than the CIMEL sun photometer does, as this latter displays a wet and dry bias with GPS for very low or no cloud cover and for high cloud cover, respectively. This is due to the fact that under such meteorological conditions, clouds contribute to the (zenith) IWV values observed by GPS in slant directions towards the different satellites. The CIMEL observations on the other hand, are always cloud-free measurements solely in the solar direction; with clouds contributing possibly only to the air mass measurement, needed to map the solar slant measurements to the zenith. As a consequence, the higher range of GPS IWVs for more cloudy skies, caused by the observation bias of the CIMEL instrument, give rise to lower regression line slope coefficients of the GPS-CIMEL scatter plots (see also

Solar slant integrated water vapour measurements
The hypothesis on the contribution of clouds to IWV measurements in different slants is based on differences observed between slant measurements mapped in the zenith direction only. As the CIMEL and GPS retrievals use different strategies to convert the measurements from slant directions to zenith, we evaluate the influence of these strategies on the GPS-CIMEL IWV comparison. One way forward is to compare coincident IWV values in the direction of the Sun (the so-called solar slant IWV). These solar slant IWVs are one step back in the CIMEL data process to obtain IWVs from the measured slant transmittance and are obtained by multiplying the (zenith) IWVs with the optical air mass as defined in Kasten and Young (1989). These CIMEL solar slant IWVs should then have no contributions from clouds, by definition.
To convert GPS IWVs into solar slant IWVs for Sun elevations higher than a cut-off angle of 7 • , we used the same wet mapping function as the one used during the IGS data analysis, i.e. the wet Global Mapping Function (GMF) of Boehm et al. (2006) for the ZTD contribution and the horizontal gradients developed by Chen and Herring (1997).
Comparing the GPS-CIMEL zenith and solar slant IWV scatter plots (Figs. 5 and 7), we find significantly higher bias and SD for the solar slant IWVs, illustrating the impact of the different mapping functions on the scatter plots.

Radiosondes
During the period of our intercomparison, three different radiosonde types were used at Uccle (see Fig. 4): the RS80 (until August 2007), the RS90 (November 2001-October 2003 and the RS92 (from September 2007), all produced by Vaisala. As the humidity sensors of the types RS90 and RS92 are identical besides some improvements in sensor design (e.g. to minimise the solar radiation heating), we will treat these two radiosonde types as one class, named RS9x. Both widely used radiosonde types RS80 and RS9x suffer from a well-known dry bias in their humidity measurements, caused by different error sources: chemical contamination and sensor ageing for the RS80 (see Sect. 2.3 and Van Malderen and De Backer, 2010, and references therein) and solar radiation for daytime RS9x observations (Vömel et al., 2007). The RS80 humidity measurements at Uccle were corrected by the method developed by Leiterer et al. (2005), which is currently the best available correction scheme for this radiosonde type (Suortti et al., 2008). Nighttime sounding data exist for the RS90, and for the RS80 until November 2001. For all these radiosonde types, measurements are taken/averaged every 10 s, so that the theoretical vertical resolution is about 100 m on average, which is much higher than the best resolution provided by the IGRA database (data points at the standard and significant levels only). A wet IWV bias with GPS is observed of 0.41 mm and 0.63 mm for RS80 and RS9x radiosondes, respectively. This has been attributed largely to a diurnal effect with a 0.06 mm (dry) and 0.32 mm (wet) bias for noon RS80 and RS9x observations respectively, relative to 1.53 mm (wet) and 1.11 mm (wet) for the two sensors' midnight observations. RS9x displays lower bias with respect to GPS, than RS80 for all the statistical analyses -demonstrating that this is a vastly improved IWV sensor. For both sensors daytime observations showed higher correlations with GPS than nighttime observations.
This daytime-nighttime difference might be explained by a different behaviour of the radiosonde humidity sensors and/or differences in the GPS IWV retrieval in daytime and nighttime conditions. For radiosondes, the heating of the humidity sensor by the solar radiation is likely to be at least partly responsible. Also the GPS data used here (IGS troposphere product) are not completely insensitive to the diurnal cycle, as the reprocessing is done with a time window of 24 h, changing at 00:00 UTC. As a consequence, IGS orbit discontinuities between adjacent days are detected Ray, 2009, 2013). The high temporal resolution of the IGS IWV data set at Brussels (10 min) enables us to investigate the differences between consecutive IGS processing cycles. Therefore, we compared both the GPS IWV retrievals from the end of a processing cycle (at times ranging between 23:30 and 23:50 UTC) and from the start of the next day processing cycle (between 00:00 and 00:30 UTC) with the radiosonde measurements at 00:00 UTC. From this test, we can conclude that the use of IWV values retrieved from the new processing cycle enhances the nighttime bias to the RS only up to 0.04 mm. Therefore, this effect could clearly be neglected for the remainder of the paper.
Comparing a mixture of daytime/nighttime RS observations, done with even different humidity sensor types, with

R. Van Malderen et al.: A multi-site IWV intercomparison
other devices, as we did in Fig. 5, should be handled with caution. However, as they are usually treated together in literature, we will do the same for this study.

GOMESCIA
Compared to the scatter plot properties obtained between the GPS IWV retrievals and the ground-based and in situ devices, the GPS-GOMESCIA scatter plot (see Fig. 5) exhibits a larger scatter and worse correlation, with a smaller linear regression slope coefficient. Of course, this worse agreement with the GPS IWV can be explained by the larger challenges that satellite IWV retrieval has to face compared to their ground-based or in situ counterparts: the cloud cover issue, the reduced sensitivity in the lower tropospheric layers and a gridded data set with issues on the pixel size. In the following subsections we investigate the spatial co-location criterion, the ground pixel size and the influence of clouds in more detail.

Influence of the spatial co-location criterion
We applied different selection criteria for the satellite measurements. For example we selected GOMESCIA measurements for which the ground pixel includes the IGS station. Alternatively, we can select measurements with minimum distance between the IGS station and the satellite ground pixel centre and additionally set an upper limit for this distance. This strategy has the advantage that it does not favour the satellite device with the largest pixel sizes (GOME), in contrast to the previous described method. In any case, the scatter plot properties of this selection criterion are very similar to the "IGS station in pixel" criterion, but of course depend on the limiting distance: the smaller the limiting distance, the lower the SD and the higher the correlation coefficient. The latter is a general finding, because it also applies to other selection criteria used in combination with this limiting distance (like simple or weighted averaging, and maximum value for the cloud flag). The effect of limiting the satellite scan angle in the measurement reduction process on the GPS-GOMESCIA scatter plot properties is less significant as compared to limiting the distance.

Influence of the satellite ground pixel size
Our GOMESCIA data set is composed of measurements with three different instruments, though using the same technique to retrieve the data from the measurements in the visible spectral range from 608 to 680 nm. Large differences in the ground pixel sizes exist between those three satellite devices, especially between GOME and SCIAMACHY/GOME-2. We therefore investigate whether the ground pixel size has a large impact on the agreement with the co-located IGS IWV retrievals. From Fig. 8, presenting the scatter plots with GPS for the three satellite devices separately, it can be derived that the GOME IWV measurements have a slight Figure 8. Scatter plot of coincident IWV measurements of either GOME, SCIAMACHY, or GOME-2 with the GPS device at Uccle, Brussels, Belgium. We selected satellite measurements for which the ground pixel includes the IGS station and with cloud flag above 1 (see Sect. 4.3.3).
inferior quality compared to SCIAMACHY/GOME-2 IWV measurements -the larger pixel size can be partly responsible for it -and that the SCIAMACHY measurements show the best agreement with the GPS IWV values. However, the inter-satellite IWV differences are not very important and it seems therefore justified to treat the three UV/VIS data sets together in the remainder of the paper. As the GOMES-CIA IWV retrieval used here (Wagner et al., 2011) applies instrument-dependent offsets to create homogeneous intersatellite time series, we also investigated the impact of this offset correction on the GPS-GOMESCIA scatter plot. For Brussels, these offsets have only a moderate impact on the GPS-GOME IWV bias (reduced by 0.5 mm), and a small impact on the GPS-SCIAMACHY bias (increased by 0.1 mm). The other scatter plot properties are not or insignificantly changed, but for the overall GPS-GOMESCIA comparison, applying the offsets brings a bias improvement of 0.22 mm and improvements in the correlation and regression slope coefficients by one thousandth.

Influence of clouds
To investigate the impact of cloud cover we selected measurements with different normalised O 2 column density for the correlation analyses versus GPS data. The normalised O 2 column density is used in the MPI-C to determine a cloud flag. A value of 1 for this flag is the applied threshold for cloud detection; higher values mean higher O 2 column densities, i.e. less cloud cover. Lower values can be caused by cloud shielding, but also by high mountains (not for the BRUS station). If all GOMESCIA observations (with the "IGS station in pixel" criterion) are used for the correlation analysis, we get a high negative bias of more than 2 mm, with a SD of almost 5 mm, and a correlation and regression slope coefficient around 0.75 and 0.65 respectively. If only observations with a cloud flag above 1 are considered much better agreement is found: the bias decreases now to −0.51 mm, the SD to 3.71 mm, while the correlation coefficient and regression slope increase to respectively 0.866 and 0.825 (see Fig. 5). If the O 2 column densities are further constrained, to e.g. the maximal O 2 column density (selecting the "best" cloud-free value), the slope of the GPS-GOMESCIA regression increases, to around 0.86 for the given example. However, this criterion might introduce for some stations systematic biases due to spatial sampling as observations for higher O 2 columns (caused by less clouds or lower surface elevation) are favoured. Finally, we perform correlation analyses for measurements with cloud flag values between 1 and 1.07 and above 1.07, respectively (Fig. 9). The value 1.07 is the median of the selected measurements, so that both data sets have more or less the same amount of observations. The best agreement, and in particular with the slope coefficient closest to one, is obtained from the GOMESCIA observations in the high range end of the O 2 absorption, hence the most cloud-free observations. Like for the CIMEL, measurements done at partly cloudy skies lead to lower regression slopes, higher SD, lower correlation coefficients and the turnover from a positive to a negative bias. However, in contrast to the CIMEL results, we assume that the worse agreement with GPS IWV observations under less clear sky conditions is caused by either an inferior data quality of the GOMESCIA IWV retrievals and/or the lower part of the atmosphere not completely sensed by the satellite device under such conditions. In this context, we mention that the mean of the GOMESCIA IWV values is higher (15.1 mm) for mostly cloud-free conditions (cloud flag above 1.07), compared to the mean for more cloudy scenes (13.3 mm), opposite to the means of the corresponding GPS-based values that decrease with decreasing cloudiness. This might actually point to the fact that the lowermost tropospheric layers are not sensed for less cloud-free scenes.

AIRS
The GPS-AIRS IWV scatter plot properties (see Fig. 5) are comparable to the results obtained for GOMESCIA and therefore reflect probably similar challenges for satellite IWV retrieval like the spatial co-location criterion and the cloud cover issue.

Influence of the spatial co-location criterion
The AIRS overpass measurements at the IGS station, shown in Fig. 5, were obtained by selecting the closest pixel to the BRUS station, passing the quality check (Qual_H2O flag values equal to 0 or 1), and with a maximum distance of 50 km between pixel centre and the station. The overall bias between AIRS and GPS is very small, only 0.01 mm, but at the cost of a rather high SD of 3.51 mm. The correlation coefficient is 0.883 and the slope coefficient is 0.842. Doubling the maximum distance between the IGS station and the AIRS pixel centre to 100 km has a negative impact on the GPS-AIRS scatter plot properties: higher SD, lower correlation coefficient, a lower slope, and a lower, even dry bias, between the AIRS and GPS IWV retrievals. This can be explained by a larger contribution from Qual_H2O = 1 data, as will become obvious in the next subsection.

Influence of clouds
Indeed, when comparing in Fig. 10 the AIRS data sets with different values for the quality flag Qual_H2O (which depend on clouds) with GPS, we found that the best quality AIRS data (Qual_H2O = 0) display a high wet bias with respect to GPS, while the AIRS data points with Qual_H2O = 1 have a small dry bias with GPS, likely due to an underestimation of the contribution of the lower tropospheric layers to the total column water vapour content. The best quality AIRS observations also show higher correlations and smaller SD with GPS. Contrary to the impact of the cloud flag on the GPS-GOMESCIA regression slope, no effect can be seen on the GPS-AIRS regression slope, which could pinpoint that the AIRS H 2 O quality flag is not very restrictive for the cloud cover.

Day-night differences
Bedka et al. (2010) reported a significant nighttime dry bias in the AIRS retrievals at IWV values above 20 mm for one of their sites and suggested that it exists on a significant spatial scale over the US Great Plains and desert Southwest, but not in the eastern United States or Canada. Fetzer et al. (2005) found an absolute bias of 0.5 mm in the IWVs retrieved by AIRS and AMSR-E during nighttime, but no bias during daytime observations. They attributed this daytime-nighttime difference to increased stratus clouds at night which have deleterious effects on the AIRS retrievals. As there are generally two AIRS overpasses a day at Brussels, one daytime and one nighttime, with identical sample sizes, we investigated the AIRS diurnal bias by studying the differences between the nighttime and daytime observations compared to the GPS IWV retrievals. The nighttime AIRS IWV data seem to suffer from a small wet bias with the GPS IWV retrievals (about 0.80 mm), whereas the daytime measurements show a dry bias (0.85 mm). Also the mean value of the coincident nighttime IWV retrievals is higher than its daytime counterpart. On the other hand, the AIRS nighttime observations show a better agreement with the GPS IWV retrievals (higher correlation coefficient, lower SD), although the regression slope between the daytime observations is closer to 1 (but with a higher sigma of this slope). Table 5 summarises the scatter plot properties of the different instrument intercomparisons at Brussels (with GPS as reference). The mean bias between the different techniques and GPS varies between −0.64 mm (SCIAMACHY) and +0.61 mm (RS9x). These are very small numbers, taking the instrumental and algorithm uncertainties of the different instruments into account. The best correlation and lowest dispersion is found for the GPS-CIMEL comparison. On the other hand, the slopes of the regression lines are closer to one for the all-weather device (RS) than for instruments demanding a (partly) clear sky (CIMEL, GOMESCIA, AIRS). When selecting only clear sky observations, these slopes increase for CIMEL and GOMESCIA, but hardly change for AIRS. Under dry conditions, GPS is less sensitive to low IWV values than the other devices (Fig. 5), which is in agreement with Wang et al. (2007). For low water vapour amounts, the ZTD is almost completely due to the ZHD. Therefore, small relative errors in these amounts produce a large relative error in their difference (ZWD) and consequently in the retrieved IWV (Schneider et al., 2010).

Summary for Brussels
The two satellite devices, although measuring water vapour in a different spectral window, agree similarly with GPS IWV. We also found only small differences between the GOME, SCIAMACHY and GOME-2 IWV comparisons with GPS, so that we will treat them as one data set for the remainder of this paper. We also confirm that Vaisala's state-of-the-art radiosonde type (RS9x) compares better with GPS IWV data than the preceding RS80 type. For both types, nighttime observations show a wet bias compared to daytime measurements and higher regression slopes, but the overall correlation with the GPS IWV values is better for daytime measurements despite the apparent radiation dry bias. For AIRS, nighttime observations also produce a wet bias with respect to GPS IWV values and AIRS daytime observations. The nighttime observations agree most accurately with the GPS IWV retrievals, despite the higher regression slopes for daytime AIRS measurements.
In addition, we compared the IWV measurements of all instruments with GPS separately for different seasons. We found a clear, identical seasonal dependence (for all instruments) of the bias and the SD: the bias is minimal (or driest) in summer, and maximal (wettest) in winter, whereas the SD is minimal in winter and maximal in summer. In this context, we remind that the largest mean IWV is obtained in Brussels in summertime and the smallest mean IWV in wintertime. Of course, this finding deserves further attention, but should also be analysed on a more global scale.

Exploitation of Northern Hemisphere IWV data sets
Before going into detail about the results obtained from the comparisons with the GPS instruments at the 28 selected sites, we draw the attention to the IGS station BRMU. As it will be clear from the forthcoming figures and discussion, IWV values retrieved from the ZTD data from this station, are exposed to a large mean bias of the order of 5 mm compared to the coincident data of the other co-located techniques. We believe that the origin lies in the ZTD to IWV conversion, and more precisely in the used surface pressure data of the co-located WMO station with code 78016. Indeed, we found large deviations between the surface pressure data of this station and reanalysis surface pressure data for the pixels surrounding this station. Using the latter data for the ZTD conversion, a better agreement between the CIMEL and the GPS IWV values is reached. As our intercomparison is meant to be purely observational, we nevertheless use the suspicious meteorological data and keep the BRMU station in the selection. This station is then illustrative for the deviations caused by using erroneous meteorological data in our analysis.

CIMEL sun photometer
A good agreement is achieved between the two ground-based devices for all the considered stations, see Fig. 11, apart from one exception (KSTU). As the GPS-based IWV values measured at the station KSTU compare well with the IWV data retrieved with GOMESCIA and AIRS, we presume that the CIMEL data of the station Krasnoyarsk has some data quality issues, although they are level 2 data. All regression slopes are inferior to 1, except for the scatter plot between BUCU and the CIMEL at Bucharest Inoe, for which a value of 1.01 is obtained. The regression slopes for the comparison of the same IGS station with two other co-located CIMEL instruments are below 1. The hypothesis that the weather observation bias (a clear sky needed in the direction to the sun) is responsible for regression slopes below 1, as postulated in Sect. 4.1.1 for Brussels, seems valid for all selected stations. Seven IGS stations are co-located with more than one CIMEL instrument within the imposed limiting distance of 30 km. This gives us the opportunity to analyse directly the data quality and uncertainties of individual CIMEL instruments and the geographical dependency of the comparisons for these co-locations. At four of these IGS stations (TLSE, VENE, OPMT, and BUCU), there is a CIMEL that has the largest differences in both the vertical and horizontal distances from the IGS station. In three of those cases (not BUCU), the data of this CIMEL show the weakest agreement with the GPS IWV values in all aspects (bias, SD, R 2 and slope); for BUCU this CIMEL has only the largest deviation from unity in the regression slope coefficient. Among the three IGS stations (TLSE, VENE, and BRMU) with a CIMEL at the minimal vertical or horizontal distance, TLSE and VENE exhibit the best GPS-CIMEL agreement in all considered aspects. Clearly, the geographical aspect plays an important role when comparing co-located CIMEL and GPS IWV values. Disentangling or even quantifying the effects of the vertical and horizontal distance on the GPS-CIMEL IWV scatter plots is not straightforward, as also other parameters (e.g. the number of coincident observations, the data quality of the measurements of each CIMEL) influence the quality of the scatter plots. In relative amounts, the bias and the SD vary most among different GPS-CIMEL scatter plots, and also the regression slope changes considerably when comparing data of different CIMELs to the same GPS instrument. The correlation coefficient seems less sensitive to the chosen CIMEL with which the GPS IWV values are compared.
Analysing the data set of GPS-CIMEL scatter plots of all selected IGS stations, we did not find any correlation between the scatter plot properties on one hand and the altitude difference or distance between the CIMEL and GPS stations on the other hand. Hence, it seems justified that we did not apply any altitude correction to the GPS IWV data when comparing with a co-located CIMEL station. Furthermore, we could not detect any dependency of any of the scatter plot properties shown in Fig. 11 on latitude, longitude, or mean observed IWV value. On the contrary, the SD decreases with latitude and increases with mean observed IWV value. These latter two dependencies seem coupled and we will come back to them in Sect. 5.5.

Radiosondes
An overview of the scatter plot properties of the GPS-RS co-locations is provided in Fig. 11. The IWV values given by co-located radiosondes and GPS devices compare fairly well, but the regression slope coefficients are smaller than 1 for all stations. The outlying station is now GLSV, with a large bias and SD and small correlation coefficient and regression slope. The poor agreement between both devices at this station is ascribed to the radiosonde data, as the GPS IWV data agree well with the IWV values retrieved from the CIMEL and the satellite instruments.
To calculate IWV values from radiosondes, we integrated the specific humidity profile starting from the co-located GPS station height, as discussed in Sect. 3.3. With this altitude correction applied, a small improvement of the bias, SD, and R 2 of the GPS-RS comparisons is obtained for the majority of the stations, but with lower regression slopes. We did not detect any correlation between the scatter plot properties and the height difference between the RS and the GPS station. Hence applying the altitude correction seems justified for our sample of stations. For completeness, note that we could also not observe a relationship between the distance between the RS and GPS station and any of the scatter plot properties.
There is neither a latitudinal (see Fig. 11) nor a longitudinal dependency of the GPS-RS bias, R 2 , and slope and we could not establish any relation between these properties and the altitude or mean IWV value of the station. On the other hand, the SD of the GPS-RS comparison decreases with latitude and increases with mean IWV. Figure 11. Column bar plots of scatter plot properties (count N , bias, R 2 and regression slope) of the different CIMEL and radiosonde instruments versus GPS for the selected sites. Sites are ordered with increasing latitude. The error bars represent the SD of the differences (bias panels) and the one-sigma uncertainty estimate of the regression slopes.

GOMESCIA
The most remarkable feature of the GPS-GOMESCIA scatter plot properties, shown in Fig. 12, is the larger variability among the different stations compared to the other two discussed instruments so far. The large SDs (error bars in the bias column bars) lead to a poorer determination of the other scatter plot properties.
As for the Brussels case study, we tried different strategies to reduce all satellite overpass measurements above the IGS station to a single overpass value. The same findings can be reached for the whole sample of stations -based on the weighted means of the scatter plot properties with weights equal to the number of coincident observations at a stationas for the Brussels case study. So, limiting the distance between the satellite ground pixel centre and the IGS station has a positive impact on the means of the SD and the correlation coefficient. Using other parameters to do the data reduction or filtering like the satellite scan angle affects to a lesser extent the means of the GPS-GOMESCIA scatter plot properties. Decreasing the satellite scan angle improves especially the SD and the R 2 of the GPS-GOMESCIA IWV scatter plots and the bias (not in absolute terms) with the GPS increases for almost all stations. The effect on the slope is more dependent on the stations considered and is therefore harder to catch on a more global scale.
As for the Brussels case study, we did not find any reason for treating the GOME, SCIAMACHY and GOME-2 data sets separately, as the differences in the scatter plots with the corresponding GPS IWV values are very minimal. We

R. Van Malderen et al.: A multi-site IWV intercomparison
compared the GPS-SCIAMACHY/GOME-2 scatter plots with the GPS-GOME scatter plots, given the similar pixel sizes of SCIAMACHY and GOME-2, but a much larger pixel size for GOME. Moreover, Noël et al. (2008) already concluded that the GOME-2 and SCIAMACHY water vapour total columns compare well on a global scale based on about 7 months of data, although an indication for a small scan angle dependency was reported. In the mean, we obtained the best agreement for the SCIAMACHY/GOME-2 and GPS IWV pairs, for all the discussed scatter plot properties, but the differences with the GPS-GOME scatter plot properties are minimal, except for the bias (−0.01 mm vs. 0.20 mm) and the slope (0.862 vs. 0.836). Applying instrument dependent offsets in the IWV retrieval algorithms to homogenise the time series has a small positive effect on the averaged overall bias of less than 0.1 mm, and practically no effect on the other scatter plot properties.
Imposing a minimum value of 1 for the GOMESCIA cloud flag is really a requisite to obtain a satisfactory agreement between GPS and GOMESCIA IWV data pairs. Also for the NH intercomparison, this selection criterion was used in combination with the requirement that the IGS station falls in the satellite ground pixel. When selecting the GOMES-CIA measurements with a maximal O 2 column density, the mean of the GPS-GOMESCIA regression slope coefficients increases to a large extent, but at the cost of introducing systematic biases due to the geographical dependence of the selection. For instance, in the case of the IGS station MARS (Marseille, France, a coastal station) this methodology resulted in selecting the bulk of the pixels over the Mediterranean Sea. To further analyse the impact of the cloud cover on the scatter plots, we split for each station the GPS-GOMESCIA IWV pairs in two samples according to their observed O 2 absorption with respect to a chosen overall O 2 absorption threshold above 1, and subsequently compare the GPS-GOMESCIA scatter plot properties for the two data set samples. We can conclude that selecting the most cloudfree observations in a sample decreases the SD and increases the correlation and slope coefficients substantially, but at the cost of an increasing (absolute) bias. This is true both for the weighted means of these parameters and for the large majority of the stations. For the most cloud-free observations, in general, the GOMESCIA instruments show a wet IWV bias with the co-located GPS device, and a dry bias for less cloudfree scenes.
In terms of geographical dependency, we cannot detect any correlation of the GPS-GOMESCIA bias, R 2 and regression slope with the latitude, longitude, altitude, and mean IWV of the GPS stations (see Fig. 12). Once again, the GPS-GOMESCIA SD decreases with the latitude and increases with the mean IWV. Nevertheless, we can conclude that the GOMESCIA IWV retrieval methodology seems consistent on a near-global (NH) scale.

AIRS
A summary of the GPS-AIRS scatter plot properties for the selected locations is given in Fig. 12. The biases between the GPS and AIRS coincident IWV values are small, but with relatively high SD. The correlation coefficients roughly range between 0.85 and 0.95, and all regression slope coefficients are smaller than 1, ranging between 0.70 and 0.95. This last finding is in qualitative agreement with the outcome of the GPS-AIRS comparison over the US undertaken by Rama Varma Raja et al. (2008) who concluded that, for mid-latitudes at least, the absolute values of AIRS derived total water vapour are dry biased in moist atmospheres (IWV > 40 mm) and wet biased in dry atmospheres (IWV < 10 mm). Comparing the GPS-satellite scatter plot parameters, it should be noted that there is less variability among the different stations for the AIRS device than for the GOMESCIA instruments.
The GPS-AIRS co-locations, on which the scatter plot properties in Fig. 12 are based, are obtained by selecting the AIRS measurements with Qual_H2O flag equal to 0 or 1 and with a maximum distance of 50 km between the ground pixel centre and the IGS station. As for the Brussels case study, we elaborate more on these selection criteria. Doubling the limiting distances between the AIRS ground pixel centres and the IGS stations leads to a decrease of the bias (or a drier bias), a higher SD, and lower correlation coefficients and regression slopes. The sample of data with the best quality (Qual_H2O = 0, or less affected by clouds), shows in the (weighted) mean, higher maximum IWV values, a higher bias to the coincident GPS IWV measurements, a lower SD and a higher correlation coefficient. On the other hand, the mean IWV value and regression slope coefficient are nearly identical between both samples of different data quality. Of course, it should be mentioned that we compare two samples of different sizes: the sample of the highest data quality only has one eighth of the amount of observations of the sample with Qual_H2O = 1, which might have a significant impact on the scatter plots. All the mentioned tendencies are hence in agreement with those found for the Brussels case study, although the absence of any impact of the Qual_H2O flag on the regression slope coefficient of the linear correlation with GPS is remarkable.
As AIRS has two overpasses a day above IGS stations with latitudes between 45 • N and 45 • S, we investigated the differences between daytime and nighttime IWV measurements. Observations with SZA lower than 90 • are classified as daytime, larger values of SZA are nighttime observations. The AIRS nighttime measurements are consistently showing a better agreement with the GPS IWV retrievals, as both in the mean and for the vast majority of the stations the (absolute) bias and SD are lower and the correlation coefficients higher. The daytime measurements of almost all stations have a negative (dry) bias compared with the GPS IWV values, while a positive (wet) bias is observed for the nighttime measurements at the majority of the stations. The reason for this apparent different behaviour of the AIRS instrument during the day or night is not clear to us. On the other hand, the daytime measurements have GPS-AIRS regression slope coefficients closer to one (but with higher standard deviations) than their nighttime counterparts.
Also for the GPS-AIRS biases, R 2 and slopes, no correlation with the latitude, longitude, altitude and mean IWV of the IGS stations could be detected (see Fig. 12). For the GPS-AIRS SDs, we observe again a decrease with increasing latitude and decreasing mean IWV. Once more, the AIRS IWV retrieval seems consistent geographically from our point of view.

Summary for the NH intercomparison
In Fig. 13, we present for each technique the resulting scatter plot parameters with respect to the GPS device, by averaging the scatter plot parameters for all stations with (normalised) weights equal to the (relative) number of observations at the station. This figure enables us to evaluate the different instruments, not only with the GPS instrument as reference, but also against each other.
The mean biases of the different instruments with the GPS device vary only between −0.3 to 0.5 mm, but these small biases are partially compensated by large SD values, especially for the satellite instruments. Based on the low scatter for the CIMEL instrument and the highest correlation coefficient, we conclude that the CIMEL instrument compares best with the GPS device for the IWV retrievals. According to the same criteria, the in situ radiosonde measurements of IWV also offer a good agreement with co-located GPS IWVs, although a world-wide radiosonde data set does not constitute a homogeneous database of IWV measurements. Radiosondes also have the highest mean regression slope coefficient with the GPS, although still considerably lower than 1. We attribute this to the fact that only radiosondes and GPS are all-weather devices, while the other three devices require a (partly) clear sky. Moreover, the GPS technique is known to be less sensitive for small amounts of IWV.
Comparing the two satellite instruments, we conclude that both reveal a similar qualitative and quantitative agreement with GPS, although the AIRS instrument shows less variability in the scatter plot properties among the different selected stations. The largest difference between both comparisons with GPS IWVs is in the regression slope coefficients: the mean slope of GOMESCIA is considerably higher than the AIRS mean slope. This can partly be explained by the fact that for GOMESCIA, a small number of the stations have regression slopes superior to 1, while the maximum value for AIRS is about 0.95. The minimum slope coefficient for both techniques lies in the range 0.70 to 0.75. It also seems that the AIRS slope coefficient is less sensitive to the selection of different samples of data with different fractions of cloud cover.
Another possibility is that the AIRS cloud cover information flags are not very restrictive.
As for the Brussels station, the IWV intercomparisons for the different instruments were done separately for the different seasons. The seasonal behaviour of the overall (weighted) averages of the different scatter plot properties is similar to the Brussels case: the bias is minimal (or driest) in summer (maximal mean IWV), and maximal (wettest) in winter (minimal mean IWV), whereas the SD is minimal in winter and maximal in summer. For the SD, this variability is linked to the detected dependency of the SD -also for all considered instruments -on latitude or mean IWV: the SD with GPS decreases with increasing latitude and decreasing mean IWV. Of course, in our sample there is an anti-correlation between the mean IWV and the latitude, so that this finding can be reduced to a strong correlation between SD and IWV and an apparent anti-correlation between bias and IWV. The anticorrelation of the bias with GPS and mean IWV value can be easily explained by the fact that GPS instruments seem to have different sensitivities to IWV at the IWV extremes than the other instruments: for larger IWV values, which are more frequent in NH summer season and at lower latitude stations, the GPS retrieval leads to higher IWV values than the other instruments, so that the IWV bias with GPS (instrument minus GPS) will be smaller (or negative). The reverse reasoning can then be done for the lower end IWV values. This variation of the bias as a function of the IWV value has also already been described in the literature (e.g. Ohtani and Naito, 2000;Deblonde et al., 2005;Prasad and Singh, 2009;Vey et al., 2010), but alternative explanations were formulated. For instance, Ohtani and Naito (2000) considered the effects of the seasonal variation of the GPS mapping function as one of the causes for their observed annual variation of the GPS-RS biases. The mapping function varies as the curvature of the atmosphere changes, which is determined basically by the changes in the ratio of the thickness of the atmosphere to the radius of the Earth. Hence, when the thickness changes according to the season, it results in the variation of the real mapping function. They found that the dependence of their GPS-RS bias could be reduced by tuning observed meteorological parameters like the surface temperature, the tropopause height, the temperature lapse rate, and the height of an isothermal layer in the mapping function. The tendency for the GPS-RS SD to increase with IWV was also found and discussed by other authors (Deblonde et al., 2005, and references therein). They attributed this feature in part to stronger humidity gradients that can exist between dry and moist air when moister air is involved. In the presence of strong gradients, the location and sampling differences between GPS and RS can be more significant than for lower IWV conditions. In addition, they claimed that the presence of strong horizontal gradients in atmospheric properties can have a negative impact on the ZTD accuracy due to a breakdown of the azimuthal symmetry assumption. Other authors point to the fact that uncertainties of some techniques are dependent on the absolute measured value (e.g. RS accuracy of 4 %, Ning et al., 2012). For the other scatter plot properties (bias, R 2 , and slope), a clear dependency on latitude, longitude, height and mean IWV of the IGS station could not be detected for any of the instruments.

Conclusions
This paper presents the results of a consistent, systematic, and uniform approach comparing IWV derived from groundbased GPS and sun photometers, in situ radiosondes, and satellite-based GOME, SCIAMACHY, GOME-2 and AIRS instruments at 28 (NH) sites.
For the majority of the sites, the best IWV agreement was obtained between the ground-based and in situ instruments, especially GPS and radiosondes, both all-weather devices. CIMEL and satellite IWV retrievals require skies with low cloud fractions. Under these conditions, both agree well with the GPS-derived values, while for instance limiting the distance between the satellite ground pixel centre and the GPS station is less crucial for a good GPS-satellite IWV agreement. For radiosondes and AIRS, the nighttime-daytime IWV differences show both a clear wet bias with GPS for nighttime observations, and a small bias (RS) or dry bias (AIRS) for daytime observations.
In general, the instruments suffering from an observation bias during cloud cover (CIMEL, GOMESCIA and AIRS) have lower correlation with GPS, higher SD and decreased biases (even going from wet to dry biases) with increasing cloud cover, and underestimate the IWV values at the high end IWV range. For the CIMEL instruments, this can be explained by the fact that for larger IWV values, there is a higher probability to have clouds, which contribute directly to the GPS IWVs, but not to the CIMEL IWVs as it requires a clear sky in the solar direction. In case of the satellite instruments, the lowermost tropospheric layers are probably not sensed correctly for less cloud-free scenes. On the other hand, similar to Wang et al. (2007) and Schneider et al. (2010), we confirm that the GPS technique seems less sensitive to low IWV values than the other devices.
These different sensitivities to both low and high IWV values of the GPS instrument and the other instruments, lead to seasonal and (less apparent here) latitudinal variations of the biases with GPS. In addition, there is a strong correlation between IWV and the SD with GPS, resulting also in a seasonal and latitudinal dependency of the SD. This feature might be attributed to more significant location and sampling differences between instruments for stronger humidity gradients or the dependency of the uncertainties of some techniques on the absolute measured value.
Apart from those general tendencies, the intercomparisons are site-and instrument-dependent (also different CIMEL sun photometers compare differently with an identical GPS at the same site), as shown in Figs. 11 and 12. This agrees with Buehler et al. (2012) who also concluded that "a literature survey reveals that reported systematic differences between different techniques are study-dependent and show no overall consistent pattern". Best agreement is obtained for sites such as Brussels, Belgium (with mean biases between ±0.65 mm), where all ground-based and in situ instruments are well characterised, have a high temporal resolution, and are located at the exact same site.
For climate applications, the homogeneity of the IWV time series is a major concern. Therefore, networks of the instruments considered here (like IGS, AERONET or GRUAN 4 ) put efforts in homogeneous data processing, regular instrument calibration, traceability, and measurement uncertainty estimations. On the other hand, measurements from different satellite instruments are processed by the same technique to build up climate data records. In a following paper, we aim at assessing the long-term homogeneity of the data sets used here and we will investigate the impact of the cloud cover observation bias for some instruments (CIMEL, GOMESCIA, AIRS) on the calculated IWV trends. Unfortunately, the mentioned networks of ground-based or in situ instruments suffer from a dearth of stations in the Southern Hemisphere and in polar regions, which are pollution-free areas that are yet susceptible to large climate model and cloud biases.