Determining precipitable water vapour from upper‐air temperature, pressure and geopotential height

Radiosonde measurements of relative humidity (RH) are the main source of uncertainty in precipitable water vapour (PWV) calculation from pressure, temperature, and RH/dewpoint (PTU) data. This paper presents a formula expressing PWV in terms of pressure and temperature as functions of geopotential height (GPH), thereby allowing the PWV to be determined: (1) without any moisture‐related calculations other than those involved in measuring GPH (in radiosondes with a pressor sensor) or pressure (otherwise); (2) without relying on humidity measurements by using Global Positioning System (GPS)‐based GPH according to the gravity field, provided that pressure is directly measured. The numerical instability associated with random data errors or deviations from hydrostatic equilibrium makes the second approach unfeasible on short time scales, revealing discrepancies between the PTU‐ and GPS‐based GPHs; however, the estimation of long‐term average PWV above a location is not hindered. The estimation of PWV without humidity data was tested using high‐resolution data from 62 upper‐air stations operated by the NOAA National Weather Service. The seasonal mean {DJF, MAM, JJA, SON} PWV from the surface to 300 hPa calculated from PT and GPS data over the period 2016–2018, after rejecting individual estimates inconsistent with the 0%–100% RH range, showed a mean bias error of {−0.1, +0.1, −1.4, −0.9} kg·m−2 relative to the PTU‐based values across the stations, and a RMSE ranging from 2.4 (DJF) to 3.2 (JJA) kg·m−2. By restricting the analysis to observations with above‐average matching between the PTU‐ and GPS‐based GPH, the bias magnitude and RMSE reduced respectively to less than 0.5 and 1 kg·m−2 in all seasons. The results indicate that evaluating the long‐term agreement between the two PWV calculation methods at different sites could be useful in detecting systematic observation errors in GPS radiosonde systems using a pressure sensor.


INTRODUCTION
The precipitable water vapour (PWV) is an important climatic parameter for several reasons: it provides a bulk measure of water vapour in the atmosphere above a place; as it is generally expected to increase with the mean sea surface temperature, it is useful for monitoring global climate change; helps quantify atmospheric rivers and predict rainfall (e.g., Eiras-Barca et al., 2022;Gimeno et al., 2014;Trenberth et al., 2005;Tuller, 1968;Young et al., 2020).Radiosonde (hereinafter RS) observations probably remain the most reliable basis for determining PWV, although they are unevenly distributed over the globe and taken twice a day.Retrievals from satellite-borne thermal-infrared and microwave spectrometers can map the geographical PWV distribution better than RS-based measurements.Ground-based microwave radiometers and the processing of data from Global Navigation Satellite System (GNSS) receivers offer PWV time series with high temporal resolution.In the last two decades, the GNSS technique has proven to be particularly useful for PWV monitoring regardless of weather conditions (for a review, see Vaquero-Martínez & Antón, 2021).Nevertheless, the PWV derived from RS data has served as a baseline for the validation of different remote techniques.Leaving aside deviations from hydrostatic equilibrium and the spatial-temporal representativeness of balloon soundings, the accuracy of RS-based PWV measurement depends primarily on the accurate calculation of specific humidity, which in turn depends on the accuracy of pressure, temperature and humidity (PTU) observations.Given the uncertainty of relative humidity (RH) measurements, which depend on the performance of the humidity sensor under varying weather conditions, it would be interesting to have a method for deriving PWV based entirely on how atmospheric pressure and temperature change with altitude.The Global Positioning System (GPS) technology began to be applied in weather RSs in the mid-1990s, improving the tracking of the RS and thus the accuracy of upper-air wind measurements (Brettle & Galvin, 2003;CEPT/ERC, 1996;Dabberdt et al., 2002).Gradually, the so-called 'GPS radiosondes' have replaced the rawinsondes tracked by a radio direction-finding device (radio-theodolite or radar) or a ground-based navaid system (Loran-C, already in disuse at the time) (all tracking systems are described in WMO, 2021).According to the data reported in Ingleby's (2017) ECMWF Technical Memorandum, by the end of 2016 over two-thirds of the RS stations used GPS tracking; radar tracking, mostly used in the Russian Federation, Kazakhstan, Ukraine, and China, represented almost a fifth; and radio-theodolite plus other wind-finding systems stood for about one-sixth.A further advantage of GPS tracking is that it provides the altitude of the RS, allowing the way of measuring atmospheric pressure and geopotential height to be reversed.The geopotential height (hereafter GPH) can be obtained from the gravity field, instead of deriving it from the PTU data using the hydrostatic equation.As long as the surface pressure is accurately measured, the (hydrostatic) pressure aloft can be derived from the vertical data of GPS height, temperature and RH without the need to use an on-board pressure sensor -as envisaged by Brettle and Galvin (2003) for GPS RSs.The basic idea was not new, since radar RSs without a pressure sensor had been introduced in the 1960s and developed over time in the former Soviet Union (Zaitseva, 1993); nevertheless, radar-based height and derived pressure measurements do not achieve the same accuracy as GPS-based ones, especially at low radar elevation angles, nor equivalent consistency between stations (Nash, 2015).Once applied to GPS RSs, the new measuring technique proved to be a breakthrough in the accuracy and precision of pressure and GPH measurements at stratospheric levels (up to the balloon burst height, ∼30 km), as reported by experts from the World Meteorological Organization (WMO) Commission for Instruments and Methods of Observations (CIMO), who carried out international RS comparisons in the 2000s (Jeannet et al., 2008;Nash et al., 2006).The practical advantage of having to calibrate one less sensor, while reducing the cost of weather RSs, has led manufacturers to develop GPS RSs without a pressor sensor.This has brought a major instrumental change in upper-air observations.A significant example is the introduction of the RS41 by Vaisala in late 2013 as a replacement for the RS92(-SGP) -as this sonde was widely used at the time, mainly in Western Europe, South America, Canada, and Australia (Edwards, 2016;Edwards et al., 2014;Ingleby, 2017;Jensen et al., 2016).However, the low accuracy of the so-called 'GPS-based pressure' in the lower troposphere has led to the relaxation of the pressure accuracy requirements to accommodate the new measurements.In this regard, the GRUAN [Global Climate Observing System (GCOS) Reference Upper Air Network] report on the 2010 WMO intercomparison of radiosonde systems noted that 'it may still be desirable to retain an independent pressure measurement' (Miloshevich et al., 2012).In fact, RS manufacturers have designed hybrid solutions (examples taken from Ingleby, 2017 are in square brackets): versions with/without pressure sensor [Vaisala RS41;; option to measure 'barometric' pressure [Graw DM-09; Graw DFM-17]; direct and indirect pressure measurements made together .Although most of today's RS types no longer have a pressure sensor, a few years ago most national RS networks were still using RSs having the three P-T-U sensors (see Figure 1 in Ingleby, 2017).A reasonable explanation is that most radiosondes are 'single use', and so 'a balance is needed between cost and performance' (Ingleby & Edwards, 2015).Facing the end of production of the RS92 in August 2017, the impact of the RS92-RS41 transition on long-term RS data records was specifically accessed at GRUAN stations with twin sounding studies focused on temperature and humidity differences, not pressure (Dirksen et al., 2020;Sun et al., 2019).Nevertheless, sensor-based pressure measurements are at an advantage as an input to non-hydrostatic numerical weather prediction (NWP) models, as noted by Nash (2015).They will probably continue to be needed for less pressing purposes than NWP, such as gravity wave modelling and planetary boundary layer studies.
This paper explores the feasibility of determining the PWV based on data from GPS RS systems using a pressure sensor without relying on humidity data.This aim is grounded on the fact that moisture increases slightly the vertical distance between isobaric surfaces, because moist air is less dense than dry air.The GPH thickness ΔZ of a columnar layer bounded by two given pressure levels is proportional to the mean virtual temperature over the layer.There is no simple relationship between the change in mean layer virtual temperature and the change in PWV content in the layer, although they are roughly proportional (see Elliott et al., 1994).However, the thickness of an elemental layer neglecting humidity, Z dry can be deduced from its temperature and pressure difference p by the gas law and the hydrostatic equation; the ratio Z/Z dry (where Z is the true thickness) is equal to the ratio of virtual temperature to temperature, which bears a well-known one-to-one relationship with specific humidity; and the PWV is proportional to the mass-weighted vertical integral of specific humidity.Based on the above reasoning, the PWV at a given place can be formulated in terms of PT as functions of GPH over the local atmospheric column.Since the GPH can be inferred from the GPS position, in principle the PWV can be derived without humidity measurements (if pressure is measured directly, of course).This approach, although only applicable to a portion of the RS observations made so far, has the potential to circumvent the limitations of humidity measurements in cloudy environments or extremely dry or cold conditions (sensor biases or missing days), as well as inhomogeneities due to sensor differences among stations and sensor changes over time (Elliott, 1995;McCarthy et al., 2009).
The question is whether RS measurements of PT and GPS height are accurate enough to estimate the PWV from the increase (/decrease) in elemental layer thickness (/air density) related to moisture content.Predictably, the random data errors and flaws of the hydrostatic assumption can easily mask the small variations in layer thickness due to moisture changes.Deviations from the ideal gas law should also be considered.However, as long as PT and GPS height measurements have negligible systematic errors, it should be possible to determine the average PWV over long time periods and even capture seasonal variation.We will experiment with this hypothesis with 3 years of observations from the National Oceanic and Atmospheric Administration's National Weather Service (NOAA NWS) upper-air network.
The remainder of this paper is organised as follows.Section 2 aims to provide an understanding of the standard equations used for determining the PWV from RS observations, the approximations involved and the propagation of measuring errors, as well as of the methods for measuring the GPH and their limitations.In Section 3, the effect of moisture in the thickness of a pressure layer is studied in terms of the PWV within the layer.Section 4 introduces a general formula to determine the PWV from PT and GPH, regardless of how pressure and GPH are measured.Regarding its use without any information on humidity, requiring independent measurements of pressure and GPH, the numerical instability due to error propagation and/or deviations from hydrostatic equilibrium is discussed; conversely, the potential interest for climate studies is pointed out.Section 5 compares the seasonal and annual mean PWV calculated from PT and GPS data at 62 of the NOAA NWS upper-air stations over the period 2016-2018 with the corresponding values calculated as usual from PTU data.The agreement between the two mean PWV estimates is analysed in terms of the consistency between GPS-and PTU-based GPHs.In view of accurate daily calculations, the PWV time series calculated using the reported GPH instead of RH over half of the year 2018 at two stations with contrasting annual cycle are compared with the coincident time series obtained by the standard method.Conclusions follow in Section 6.

FUNDAMENTALS
The standard PWV formulation assumes hydrostatic equilibrium and involves relating specific humidity to the PTU variables.The formula presented later in this paper is grounded on the differential form of the equation used in the calculation of the geopotential difference between standard pressure levels, but accounting for deviations from ideal gas behaviour -as described in the sixth edition of the Smithsonian Meteorological Tables (hereinafter SMT) (List, 1949); i in addition, it uses the definition of GPH and the relation between virtual temperature and specific humidity.These topics must be carefully discussed in order to clarify the assumptions and uncertainties in routine calculations of PWV and GPH.

RS-based PWV calculation
In alphanumeric reports of upper-level PTU and wind data transmitted from fixed land stations, sea stations and mobile stations (TEMP, TEMP SHIP and TEMP MOBILE), temperature and dewpoint depression are reported (as available) at the standard pressure levels, significant levels, and tropopause level(s) (WMO, 2011, pp. 65-68).The pressures of non-standard levels from the surface up to and including the 100-hPa level are reported in whole hectopascals; lower pressures are reported in tenths of a hectopascal.Temperature is reported with a practical resolution of 0.2 • C; the coding process yields an average offset of −0.05 • C for temperatures measured to one decimal point (Ingleby & Edwards, 2015).Dewpoint depression is reported with a resolution of 0.1 for values <5 • C, and a resolution of 1 The standard method for estimating PWV from PTU data involves a certain number of approximations and uncertainties.(1) The acceleration of gravity is conveniently taken constant in the hydrostatic equilibrium equation.(2) Finding the specific humidity requires estimating the saturation vapour pressure of moist air and involves the propagation of temperature and RH/dewpoint depression errors.(3) For practical reasons, the vertical range of PWV integration must be limited.These aspects are addressed respectively in Sections 2.1.1-2.2.3.

Definition of PWV
The PWV above a location is the total mass of water vapour (WV) contained in a vertical column extending from the surface to a height representing the 'top' of the atmosphere, per unit cross-sectional area: where z is the height relative to mean sea level (MSL) (i.e., altitude) and  v is the mass of WV per unit volume of air.Although the integral extends virtually to the mesosphere where the highest and thinner clouds on Earth are found -and the designation 'integrated WV' (IWV) might be preferred -the contribution of the stratosphere and mesosphere to the total IWV is negligible and smaller than the uncertainty achieved by any measurement system.Alternatively, and in a literal sense, the PWV is often defined as the depth of water at the surface if all the WV contained in the column were to precipitate as rain: PWV * ≡ PVW ÷  w , where  w is the density of liquid water.For the air temperatures found at the Earth's surface, the value of PWV * in mm is nearly equal to that of PWV in kg⋅m −2 .At a temperature of 15 • C (≈ global average surface air temperature) 1∕ w = 0.001001 m 3 ⋅kg −1 and a PWV value of 1 kg⋅m −2 is equivalent to a PWV * value of 1.001 mm.Since the exact conversion factor depends on the temperature assumed for the precipitated water, the Equation (1) definition will be used in this paper.
Assuming hydrostatic equilibrium, where p is the atmospheric pressure,  is the air density and g is the acceleration due to gravity, Equation (1) becomes: where q =  v ∕ is the specific humidity (mass of WV per unit mass of air) and g 0 is a constant approximation for g allowing the PWV to be estimated without worrying about the local gravity.The specific humidity has to be derived from PTU data, as detailed later in Section 2.1.2.Equation (3) was first applied to upper-air observations taken by meteorographs (Solot, 1939) and has been in use ever since meteorological RSs began regularly measuring humidity.The value historically used for g 0 is 9.8 m⋅s −2 , which differs from the average gravity over the surface of the earth by ≈0.002 m⋅s −2 .But the standard gravity 9.80665 m⋅s −2 , broadly used in meteorology and differing from normal gravity at 45 • latitude by less than 0.001 m⋅s −2 (see Appendix B), can also be used in Equation (3) -making the PWV estimates smaller by 0.07%, which represents less than 0.05 kg⋅m −2 for PWV values as high as 70 kg⋅m −2 .At any rate, making g ≈ g 0 introduces negligible error in the PWV estimation, as shown below adopting standard gravity for g 0 and estimating g by theoretical gravity (see Appendix A).First, g increases from the equator to the poles and the relative difference between g 0 and g at sea level is within ±0.3%.Secondly, g decreases with increasing height by about 0.03% per km irrespective of latitude.The two variations combined give −0.3% ≤ (g 0 − g)∕g ≤ 0.6% in the layer 0-10 km, with −0.3% referring to the sea level at the poles and 0.6% referring to 10 km altitude at the equator.On average, the specific humidity of air decreases over height in an exponential-like manner; on zonal and annual average it decreases by about five times from the 1000-hPa to the 600-hPa level at latitudes below 70 degrees, regardless of the latitudinal differences in magnitude (see Figure 1.14 in Salby, 1996; note that WV mixing ratio ≈ specific humidity).Since the air density profile is nearly proportional to the pressure profile,  v ∝ qp decreases by about 5 × 1000 ÷ 600 = 8.3 times along the above-said layer; and since this is 4.1 km thick on global average, the e-folding scale of  v can be estimated as (4.1 km) ÷ ln 8.3 ≈ 2 km.From this result and Equation (1), one expects 50% of the PWV to be within approx.−(2 km) × ln (1 -0.5) = 1.4 km above the surface; 75% within approx.−(2 km) × ln (1 -0.75) = 2.8 km; and, likewise, 90% within 4.6 km.Consequently, taking g ≈ g 0 significantly affects the integral in Equation (3) only in the lower troposphere and so the resulting PWV error is mainly due the latitudinal dependence of g varying.A careful analysis indicates a relative PWV error varying from ≈ 0.2% at the poles to ≈ −0.4% at the equator.Considering the uncertainty in specific humidity due to uncertainties in RS measurements, changes in WV distribution during the balloon flight (which takes more than half an hour to rise 10 km) due to displacements of air masses (uplifts and wind advection) or to evaporation/condensation processes (cloud formation/dissipation and rain), local deviations from hydrostatic equilibrium, and finally representativeness issues related to the wind-driven balloon drift, it becomes clear that approximating g by standard gravity in Equation ( 3) is not only convenient for computational reasons but also completely reasonable.Note that this is equivalent to approximating geometric height by GPH in Equation ( 1) with the appropriate change of variable in Equation (2) (cf.Equation 2 in Section 2.2.2).

Specific humidity calculation
The equations to find specific humidity from PTU observations are described and discussed next.The molar mass of moist air, , where x v is the mole fraction of WV, M d is the mean molar mass of dry air and M v is the molar mass of WV, is usually formulated as where where e ′ denotes the WV pressure in air (the prime serves to distinguish from vapour pressure e of pure water at the same temperature and volume); the last expression uses the definition e ′ = x v p, which is based on Dalton's law of partial pressures.Using the definition of RH with respect to liquid water ii (WMO, 1953;WMO, 2021), where x vs is the mole fraction of WV which the air would have if it were saturated with respect to liquid water at the same PT, and e ′ s is the corresponding WV pressure, Equation (5) becomes: Note: the WV mixing ratio r (ratio of the mass of WV to the mass of dry air in a given volume) can be formulated by substituting Equation (7) into the relation r = q∕(1 − q); the approximation q ≈ r should be avoided in accurate calculations of PWV, particularly in the tropics (Elliott et al., 1991).The saturation pressure of WV in moist air, e ′ s , and in pure phase, e s , are not the same: where T is the air temperature.The enhancement factor f reflects the pressure exerted by air on water, the dissolution of air in water, and the interaction between air and WV molecules.By formally defining f = x vs p∕e s (as in Hyland andWexler, 1973 andHyland, 1975) and U = 100x v ∕x vs , it follows that x v = 0.01Ufe s ∕p; this approach is used by the BIPM (International Bureau of Weights and Measures) (Giacomo, 1982).
Regarding the temperature dependence of saturation vapour pressure over pure water (hereinafter SVP), Elliott andGaffen (1991), Alduchov andEskridge (1996) and Vömel (2016) provide information on the use of many formulations of e s (T), showing their differences over the atmospheric temperature range.Murphy and Koop (2005) give an insight on the difficulties of finding an accurate expression at temperatures below −30 • C. For the computation of different measures of humidity, the Goff and Gratch equation of 1946 adopted in the SMT (List, 1949) was recommended by the WMO in the past and remained widely used perhaps until the 1990s.Presently, the WMO ( 2021) CIMO Guide provides one of Sonntag's SVP equations (formula 'SA90' in Alduchov & Eskridge, 1996).Different SVP formulas diverge as temperature decreases below 0 • C; excluding the Magnus-Tetens formula, differences are generally less than 1%, 2% and 5% in the range (0, −20) • C, (−20, −40) • C and (−40, −60) • C respectively.Concerning the determination of specific humidity with the purpose of calculating PWV, outside the polar regions the choice of a particular SVP formula is irrelevant: the relative error in SVP is much smaller than the relative uncertainty in RH, notably in the coldest/highest tropospheric layers, and these layers contain a small fraction of the PWV.A related issue is the SVP equation used in the factory calibration of humidity sensors.By historical convention, RH is measured with respect to water irrespective of temperature (Equation 6).In modern RS designs, humidity sensors are designed to operate at temperatures from 40 • C down to −70 • C, aiming to provide useful RH measurements for NWP despite the uncertainties at the low end.However, the calibration of RH sensors under laboratory conditions at very low temperatures involves measuring saturation vapour pressures over ice, since saturation over water is physically impossible below −37 • C. Thus, the SVP with respect to water is needed to convert the WV pressures measured during calibration into a RH reading (Appendix A in Miloshevich et al., 2006;Nash et al., 2011 pp. 139-140).For perhaps two to three decades now, the well-known Wexler, Hyland and Wexler, and Sonntag equations (published in 1976, 1983and 1994 respectively) have been the most widely used SVP formulations among manufacturers in western countries for RS calibration and processing, as they give similar values in the entire temperature range of interest -not meaning that those three equations are inherently better than other for meteorological use; actually, their accuracy at very low temperatures, to which they must always be extrapolated, remains uncertain (Murphy & Koop, 2005).But older SVP formulas are still used by some RS manufacturers, mainly the Goff and Gratch equation.Ideally, when dealing with RS observations in cold regions such as the upper troposphere and the entire polar troposphere, the formula used by the RS manufacturer should be used again when converting the RH measured by the RS into WV pressure.The choice of a particular equation for the WV-water equilibrium vapour pressure is, however, irrelevant for the purpose of estimating the PWV outside of polar regions, given the tiny WV amounts residing in the upper troposphere.The Magnus formula updated by Alduchov and Eskridge (1996), which compares well to the best formulations known, will be used throughout this paper (see Appendix D).
Regarding the enhancement factor, Table 89 of the SMT (List, 1949) [reproduced as Table 4.10 in WMO (1966) International Meteorological Tables] contains reference values for f in the meteorological range based on data by Goff and Gratch.The wide-ranging estimates by Hyland (1975) include the atmospheric pressure values 1000, 500, and 250 hPa and cover all temperatures of meteorological interest.These determinations show that, in the meteorological range of PT: 1.000 ≤ f ≤ 1.005, except for near-sea-level pressures under extreme temperatures, roughly above 50 • C or below −40 • C when f > 1.005; and f − 1 drops from about 0.4% to 0.2% between the 1000and 500-hPa levels (see Appendix D).The approximation e ′ s (p, T) ≈ e s (T), with relative error ≈ 1 − f , is often used in practice.Since q is nearly proportional to e ′ s ∕p (cf.Equation 7) and roughly 90% of the PWV is found below the 500-hP level, the average percentage error introduced in the calculation of mid-latitude PWV by taking f = 1 should be between −0.2% and −0.4%.The relative error in q due to approximating WV pressure as e ′ ≈ 0.01Ue s is much smaller than that expected from the uncertainty of temperature and RH measurements, as shown below.
Taking the logarithm derivative of Equation (D1), using absolute temperature we can write: Assuming a constant lapse rate  in the troposphere, a pressure level p k below the tropopause corresponds to the temperature T k ≈ T 0 −  H p (ln p 0 − ln p k ), where T 0 and p 0 are respectively the air temperature and pressure at the MSL, and H p is the mean pressure scale height in the troposphere.Let T 0 ≈ 288 K, p 0 ≈ 1013 hPa and  H p ≈ 48 K, consistent with the International Standard Atmosphere (ISA).Then, for instance, an error T = ±0.5 K yields e s ∕e s ≈ ±3%, ±4% and ±5% at the 1000-, 500-, and 300-hPa levels respectively.In polar air, for which the air temperatures at successive pressure levels are much lower than the ISA values, |e s |∕e s ≤ 7% for Celsius temperatures ≥−70 • C and |T| ≤ 0.5 K. Regarding RH errors, |U|∕U can easily surpass |e s |∕e s .For instance, U = ±5% with U > 10% can introduce a percentage error between ±5% and ±50% in the product e s U. Considering the global vertical profile of RH from 1000 hPa (U ≅ 75%) to 300 hPa (U ≅ ≈ 35%) (Peixoto & Oort, 1996), U = ±5% represents an average relative error increasing from 7% at the surface to 14% at 300 hPa.Since specific humidity is nearly proportional to RH, its estimation will be affected by the same percentage error.Assuming that T and U are respectively measured with precisions of 0.2-0.5 • C and 2%-5% RH in the troposphere, from the above reasoning the RH errors should be the main factor of uncertainty of specific humidity calculations, and thus of RS-based PWV.Based on RS observations of NOAA's NWS, Connell and Miller (1995) showed that a consistent absolute error of 5% in RH measurements from the surface up to a height of 4 km would distort the PWV estimation in the same layer by ≈1 (/3) kg⋅m −2 in a dry (/tropical humid) environment -representing about 10% (/7%) of the measured PWV.
The amount of moisture in the air can be represented by the dewpoint temperature T d (temperature below which air at constant pressure becomes saturated with water).Indeed, for the exchange of RS observations over the WMO's Global Telecommunication System (GTS) the RH directly measured by RS is converted and reported as dewpoint depression, i.e.,: D ≡ T − T d (for current formulas for calculating T d from U and T, see OFCM, 1997 andWMO, 2021).Such conversion became common practice from around 1970, judging from the archived RS data from multiple sources (Durre et al., 2018).As remarked by Elliott and Gaffen (1991) referring to alphanumeric TEMP reports, in the backward conversion needed to recover RH the calculated RH is particularly influenced by the limited reporting precision of the observed D whenever D > 5 • C, roughly occurring with RH < 60%-75% depending on temperature.Noteworthy, the recovering of RH has been carried out whenever necessary in the construction of the Integrated Global Radiosonde Archive (Durre et al., 2006(Durre et al., , 2018)).
The specific humidity (q) can be readily calculated from p and T d = T − D. By definition of dewpoint, e ′ s (p, T d ) = e ′ .Recalling Equation ( 8), the WV pressure can thus be expressed as which can be inserted in Equation ( 5) to obtain q.In practice, Equation (10) can be simplified as e ′ ≈ e s (T d ).Elliott and Gaffen (1991) estimated that the standard deviation of specific humidity errors due to random errors in reported pressure and dewpoint depression represents about 6%-8% of the magnitude of q in different tropospheric layers.Ideally, e s (T d ) should be calculated using the same SVP formula that was the one used (in inverted form) by the RS system in deriving T d .The variety of formulas used by different countries over time to calculate T d (Gaffen, 1993, Tables 2-5) should be considered when using historical records of dewpoint, notably in cold, dry air (Elliott & Gaffen, 1993).However, this issue is not of importance for PWV estimation (including climate trends) outside the polar regions.

Upper bound of PWV calculation
The PWV that can be held in the atmosphere is constrained by its temperature distribution.From Equations ( 7) and ( 8) and since (1 − )e s ⋘ p, we have: q ≤ q (U=100%) ≈ fe s ∕p.
The mole fraction of WV in saturated air is x vs (p, T) = f e s (T)∕p.Thus, according to Equation (3): PWV∕p ≤ (∕g 0 )x vs , that is, the maximum PWV content of an elemental atmospheric column layer is proportional to x vs (p, T) for the PT conditions in the layer.Since e s varies roughly exponentially with T, its geometric mean over an interval (T 1 , T 2 ) is approximately e s (T) where T = 1 2 (T 1 + T 2 ).This rationale applies for instance to the range of temperature at a given latitude () and pressure level (p).Let T be the annual and zonal mean temperature as grossly modelled by the following equation for the troposphere ignoring hemispheric asymmetries (adapted from Salby, 1996, p. 514): ) , (11) with p 0 =1000 hPa. Figure 1a depicts the Celsius temperature (T( i , p)−273.15) for p between p 0 and 300 hPa at discrete latitudes  i with a 15 • interval; note that the upper level is selected such that Equation (11) can be applied to polar latitudes.Figure 1b shows x vs ≡ e s

[ T(𝜑 i , p)
] ∕p calculated using Equation (D1).It turns out that, on average over time and longitude, PWV∕p in saturated air decreases by two orders of magnitude from the equator to the poles in the lower and middle troposphere; at 60 • (/0 • ) latitude it drops by 10 (/20) times between the 1000-and 300-hPa levels, while at polar latitudes the vertical variation is moderate, being moreover reversed in the lower troposphere.Given the deviations of air temperature from T(, p), the values of x vs observed at a given time and position can deviate greatly from x vs .Nevertheless, Figure 1 indicates that, outside the polar regions, the lower troposphere is paramount to the amount of WV that can be held in the atmosphere while the upper troposphere is of little significance.
The early PWV calculations extended from the surface to 5-8 km height using vertical data recorded by balloon-borne meteorographs (Shands, 1949;Solot, 1939).Later climatological estimates using RS data from the time when RSs began employing a hygrometer, consisting mostly of monthly means, normally set the 500-hPa as the upper bound (Durre et al., 2009;Elliott et al., 1991; F I G U R E 1 (a) Schematic annual and zonal mean temperature plotted against pressure at different latitudes N/S. (b) Mole fraction of WV at saturation over water corresponding to the thermal structure of panel (a).The mean height of the 300-hPa level ranges from around 8.0 km at the poles (near the polar tropopause) to 9.5 km at the equator.Gaffen et al., 1992;Ho & Riedel, 1979;Ross & Elliott, 1996;Tuller, 1968).In the mid-latitudes, this pressure level is around 5.5 km above sea level.However, a higher top level was also used: 400 hPa (Gaffen & Elliott, 1993); 325 hPa (Reitan, 1960); 300 hPa (Elliott et al., 1991;Elliott et al., 1994;Zhang et al., 2018); 250 hPa (Rangarajan & Mani, 1982); and 200 hPa (Zhai & Eskridge, 1997).The vertical resolution of PTU data varied between standard pressure levels and pressure intervals with resolution up to 50 hPa.The wide acceptance of 500 hPa as the upper bound was for two reasons.First, only about 5%-7% of all the atmospheric WV is above the 500-hPa level (Elliott, 1995).Secondly, measuring accurately RH in the upper troposphere is difficult because of the slow response of RH sensors under extremely low WV concentrations at low temperatures and the contamination of the sensor from water or ice in clouds earlier in the ascent.iii Until the 1990s, humidity measurements by RS at temperature levels below −30 or −40 • C were generally considered unreliable.For decades, at levels below −40 • C, humidity was reported as missing in the US upper-air network (Elliott & Gaffen, 1991) and several other national networks worldwide (Gaffen, 1993).Thus, the usefulness of RS humidity observations was limited at most to heights up to the 400-hPa level on global average and 500 hPa over the polar regions, which fortunately is enough for NWP.(Cautiously, Ross & Elliott, 1996 confined their PWV climatology over North America to the levels up to 500 hPa arguing this was the region where annual mean temperatures were warmer than −30 • C at all RS sites from low latitudes to Alaska.).Furthermore, RH values below about 20% were considered very unlikely in the troposphere.At the US radiosonde stations, between 1973 and 1992, RH data below 20% (measured at the time by an hygristor) were flagged as '19%' and therefore reported as missing (Elliott & Gaffen, 1991).Other countries used cut-off values in the 10%-20% range depending on the hygrometer (Gaffen, 1993).This practice limited the frequency of humidity observations in the upper troposphere, making this region appear, on average, wetter than it really was (Gutzler, 1993).
The improvement of the humidity sensors used in weather RSs have led to a general increase in the vertical range of humidity (or complete PTU) observations over the years (Ferreira et al., 2019).Even so, RS observations of RH in the upper troposphere remain of relatively poor quality and therefore are scarcely assimilated into NWP systems (Lorenc et al., 1996).Moreover, dewpoint depression is often reported as missing up from a height much lower than the sounding termination level.Regarding the thin-film capacitor sensors widely used in modern RSs, and according to data from WMO RS intercomparisons gathered by Nash (2015), the uncertainty in RH measurements at temperatures higher than −50 • C mostly fall in the range of 5%-14%, depending on the temperature range, RH range and time of the day (day/night).On annual average, the temperature at the 300-hPa level ranges from about −35 to −55 • C between equatorial and polar latitudes (see Figure 1a).Since the upper troposphere contains a small amount of WV compared to the lower and middle troposphere, the relatively large RH errors in the former should not be critical to the PWV calculation; therefore, it is not unwise to set the upper bound at 300 hPa instead of 500 hPa.In PWV calculations based on reanalysis data (PTU gridded data), the assimilation of satellite observations allows the upper bound to be extended further to 100 hPa (as in Bock et al., 2005) thereby encompassing the tropical tropopause and the extratropical lowermost stratosphere.

Air density and layer thickness
Section 2.2.1 reviews the equation of state for moist air including the compressibility factor of air and the formulation of virtual temperature, comparing the effect of moisture versus deviations from ideal gas behaviour in air density.The related equations for hydrostatic equilibrium and pressure layer thickness are described in Section 2.2.2.Section 2.2.3 compares the alternative methods for measuring the GPH of isobaric surfaces by RS, discussing the overestimation of PTU-based GPH due to neglecting the compressibility factor of air (the theoretical gravity and GPH formula used in meteorology are examined in Appendix A).

Equation of state for real air
The virial expansion of the equation of state of a gas can be written as (Zemansky & Dittman, 1981): where p is the gas pressure, v is the molar volume, T is the absolute temperature, R * is the universal molar gas constant, and {A i } are the virial coefficients; the order i arises from the exponents of the expansion p = R * T ∑ ( A i ∕v i ) with A 1 = 1.The quantity C, which represents the ratio pv∕R * T, is called the compressibility factor.iv The difference C − 1 measures the deviation of the real gas from ideal gas behaviour due to molecular interactions.For a single-constituent gas, A 2 , A 3 , etc, are functions of T that depend on the nature of the gas; for a gas mixture, such as air, each mixture virial coefficient of order i ≥ 2 is determined by the single-gas virial coefficients and the mole fraction of each constituent (Jones, 1979;Sengers et al., 1971).In practice, it is sufficient to consider virial coefficients up to the third order, that is, A 2 and A 3 , which are related to interactions between pairs and triplets of molecules.In the case of moist air, whose composition varies mainly due to the variation of WV assuming a fixed mole fraction for CO 2 , these coefficients depend on temperature and WV mole fraction.Over the pressure range of the Earth's atmosphere, C < 1, which means that attractive intermolecular forces dominate over repulsive forces and thus tend to make the molar volume of air smaller than expected for an ideal gas at the same PT.Assuming small deviations from perfect gas behaviour, R * T ≈ pv with relative error 1 − C ≪ 1.Hence, substituting v accordingly in the first two higher-order terms of the right-hand side of Equation ( 12), where B 2 ≡ A 2 ∕R * and B 3 ≡ A 3 ∕(R * ) 2 .Appendix C quotes an empirical formula in the above form, with expressions for B 2 and B 3 in terms of Celsius temperature and WV mole fraction -as approved by the International Committee for Weights and Measures (CIPM) to accurately determine air density and thus eliminate buoyancy in weight measurements (Picard et al., 2008) -and discusses its extrapolation to the upper-air environment by comparison with the values of C tabulated in the SMT.
Defining virtual temperature as in Bohren and Albrecht (1998), substituting M = v, and combining with Equation ( 12) while omitting the virial expansion, the equation of state takes the form: where R d ≡ R * ∕M d is the specific gas constant of dry air and T ′ v ≡ CT v is the 'adjusted virtual temperature' -to use the terminology of the SMT (List, 1949).The ideal gas law for dry air (C = 1, M = M d ) reads: p = R d T. So, in the ideal gas approximation, the virtual temperature T v is the temperature that dry air would have if its pressure and density were equal to those of a given sample of moist air -this is the classical definition of virtual temperature.Considering real air as in Equation ( 15), the same definition applies to T ′ v , if instead of saying 'dry air' we say 'ideal dry air'.Eliminating x v between q = x v M d ∕M and Equation (4), solving in M d ∕M and substituting the resulting equation into Equation ( 14) yields which is a modern expression for the virtual temperature (alternative derivations can be found in Iribarne & Godson, 1981, Salby, 1996and Bohren & Albrecht, 1998).
The classical (and equivalent) mathematical definition of T v in terms of mixing ratio r rather than of specific humidity q can be obtained by substituting q = r∕(1 + r) into Equation ( 16).v Combining Equation ( 16) with Equation ( 7) and using e ′ s ≈ e s (T) (i.e., f ≈ 1), the virtual temperature of air at pressure p, temperature T and with RH U is given by: which is widely used in meteorology (Nash, 2015).The relative error introduced by the approximation e ′ s ≈ e s (T) in the above formula is less than 0.2% ×e s (T)∕p, since (f − 1)(1 − ) < 0.002; furthermore, e s is two to three orders of magnitude smaller than p in the lower to middle troposphere.Note that the WV pressure term 0.01U × e s (T) can be conveniently expressed as e s (T d ) when the dew point is known while the RH is unknown.
Using the defining Equations ( 4) and ( 14), the equation of state (Equation 15) can be explicitly stated as follows: This is the fundamental equation for the density of moist air, to be determined from directly measurable quantities (p, T, U) by substituting x v = 0.01Ufe s ∕p, provided that the functions C = C(p, T, x v ), e s = e s (T) and f = f (p, T) are formulated semi-empirically (Davis, 1992;Jones, 1978).Clearly, C is the factor by which the density calculated for the dry air-WV ideal mixture must be divided to obtain the density of real air, thus increasing density relative to the ideal gas value (C < 1 in the natural atmosphere).In order to see how much the deviations from the ideal gas law affects the air density in comparison with the moisture content, relative to ideal dry air at the same PT, let us write Equation ( 18) as where Note:  can be seen as the relative difference between T v and T (cf.Equation 17 with e s ⋘ p). Figure 2 describes the distribution of  in the troposphere up to the 300-hPa level considering the average distribution of temperature with pressure and latitude as modelled by Equation ( 11); C was calculated by Equations ( C1) and (C2) together with Equation (D1).Since (1 − ) ≈ 2 ∕ 5 and (typically) U ∼ 50%,  and  will be of the same order of magnitude if x vs ∼ 5 ×(1 − C).From this condition and comparing Figures 1 and 2 we obtain: (1) at polar latitudes, the deviations from ideal gas behaviour increase the density of the air much more than the moisture content decreases it; (2) at middle latitudes, the two effects are comparable in order of magnitude; and (3) at low latitudes moisture is the dominant factor.As will be shown, both factors influence the vertical distance between isobaric surfaces and thus their height, although this feature depends essentially on temperature.
F I G U R E 2 Vertical curves of 1 -C at different latitudes N/S for the round-year, zonal mean temperatures represented in panel (a) of Figure 1, with 0% RH (solid) and 100% RH (dashed).

Equations for hydrostatic equilibrium
Using the definition of geopotential altitude (Sissenwine et al., 1962) where  is the geopotential (potential energy of a unit mass relative to the MSL) and g 0 is the standard gravity, the fundamental hydrostatic equilibrium Equation (2) becomes:−p∕Z = g 0 .Combining this with the equation of state Equation ( 15), the hydrostatic equilibrium can be stated as: Note: in meteorology, Z is usually called 'geopotential height' (abbreviated in this paper as GPH), just as z is called height (or 'geometric height') instead of altitude, both being understood to be zero at MSL. Equation ( 21) differs from usual by having T ′ v instead of T v , allowing for deviations from ideal gas behaviour (C ≠ 1).Although g 0 = 9.80665 m ⋅ s −2 is an integral part of the GPH definition Equation (20) (see Appendix B), any other reference gravity acceleration would do for Equation ( 21), without changing it in any way; however, choosing it to be very close to the observed gravity at MSL and 45 • latitude has the obvious advantage of making Z similar to z in value.
Inverting the derivative in Equation ( 21) and integrating between any two pressure levels yields: This is a generalization of the commonly called 'hypsometric equation', as it includes the compressibility factor of air.It can be rewritten as (following basically from−p∕ = ).Equation ( 23) was published in the SMT, long before the concept of GPH was created, and was then referred to as 'hydrostatic equation' with the remark that 'for most meteorological purposes' it can be 'closely approximated' by replacing T ′ v with T v according to the ideal gas approach (List, 1949, p. 296).Along these lines, Equation ( 23) was restated with T v in place of T ′ v in the WMO International Meteorological Tables (WMO, 1966).
Equation ( 21) was used by Aparicio et al. (2009), who showed that the compressibility factor is important to realign the GPS-and PTU-based GPH scales when assimilating GPS radio occultation data into NWP systems.Herring and Quinn (1999) used it in an atmospheric model for reducing the upper-level atmospheric field of global weather analysis to a (topographic) surface pressure needed for atmospheric delay correction in satellite laser altimetry.As will be seen in Section 4, including the compressibility factor in the differential form of the hypsometric equation is also important for the purposes of this paper.

GPH in RS observations
Geopotential gradients on isobaric surfaces are related to the strength of synoptic-scale wind in the free atmosphere; moreover, the geopotential is a key quantity to accommodate the pressure as a vertical coordinate, allowing for simplification of the momentum equations in NWP and general circulation models (GCM).The value of g 0 used to scale geopotential for the purpose of weather analysis is immaterial as long as it is used consistently.The value 9.8 m⋅s −2 that sometimes appears in the literature to roughly define GPH (e.g.,: Wallace & Hobbs, 1977, Salby, 1996, WMO, 2011) is reminiscent of the fact that geopotential in upper-air reports and analysis were formerly expressed in 'geopotential metres' (gpm = 9.8 m 2 ⋅s −2 = 9.8 J⋅kg −1 ) -as decided in 1949 by the IMO (List, 1949, pp. 217 and 224) and subsequently adopted by the WMO (Introduction 3.1 of WMO, 1966).
The GPH definition was recommended by the WMO in the mid-1990s ('Suppl.No. 2 (IV.1996)' to WMO, 1988 Technical Regulations Vol.I).Since the GPH in metres is numerically equal to geopotential in 'standard geopotential metres' (std.gpm = 9.80665 m 2 ⋅s −2 ), eventually, the gpm unit was replaced by std.gpm to be in line with the GPH.For instance, the geopotential unit used in the US RS network changed from gpm to std.gpm in October 1993 (Elliott et al., 2002).
Traditionally, Z k at a pressure level p k above a RS station and at the observing time is obtained by adding Z 0 (station elevation converted to GPH according to local gravity) to the layer thickness Z k − Z 0 given by Equation ( 22) but with C = 1, starting the integration from the reading of the station's barometer.The integration is carried out making use of the RS data transmitted to the ground, although in alphanumerical TEMP reports the GPH is only given at the standard pressure levels.Since T v > T ′ v (C < 1), RS measurements overestimate GPH.In the mid-latitude troposphere 1 − C varies typically within 0.003-0.007(see Figure 2).Thus, setting C = 1 in Equation ( 22) should produce a bias between 3 and 7 m in the GPH at the 300-hPa level of mid-latitudes, typically at an altitude in the 9.0-9.5 km range, depending on latitude and month.Since C increases with increasing temperature at the same pressure, during winter/summer the overestimation of GPH should be larger/smaller than average.The largest GPH bias is expected at polar latitudes during polar winter.At any rate, the compressibility factor of air is disregarded for a simple reason: NWP models, as well as GCMs, use the ideal gas law; and atmospheric measurements need to be consistent with atmospheric models.So, the errors of PT measurements are the only concern regarding the limits of uncertainty required for the GPH at standard pressure levels according to the intended use of the data, NWP or climate studies (see Lenhard, 1970;Richner & Viatte, 1995;Inai et al., 2015;and Nash, 2015, p. 21 and annex B).Deviations of the atmosphere from hydrostatic equilibrium are not a major concern since the hydrostatic assumption is consistently used in large-scale NWP.
By definition, the GPH obeys the gravity field.Neglecting gravity anomalies and differences between geodetic height and height from MSL, Z(, z) is represented by Equations (A8), (A9) and (A1), with numerical expressions as given in Nash (2015).Directly measuring GPH from GPS data provides better accuracy and precision in the stratosphere, where PT errors yield relatively large errors in the hypsometric equation.Vertical profiles of GPS-based GPH with respect to pressure have been used as a reference to study discrepancies in PTU-based GPH obtained from different RS types and to quantify PT errors (Inai et al., 2009;Stauffer et al., 2014;Steinbrecht et al., 2008).Note, however, that the accuracy of GPS-based GPH is nearly the same as the accuracy of GPS height.The error due to neglecting geoid and gravity anomalies is estimated in Appendix A.
Getting the most out of GPS tracking, pressure tends to be measured indirectly.The method consists of using GPS-based Z k values to obtain recursively p k by using Equation ( 21) with C = 1 and Equation ( 17) for virtual temperature, initializing with the surface pressure (cf.Nash, 2015, p. 23).The surface pressure and the difference in altitude between the ground GPS antenna and the surface barometer are critical parameters to avoid a systematic error in the so derived 'GPS-based pressure'.Regarding random errors, in the lower troposphere the sensor-based pressure profiles are less noisy and more reproducible than the GPS-based ones; for a compelling campaign experiment with the RS92 and RS41, see Lehtinen et al. (2016).This feature reflects the fact that, according to the hydrostatic equation, random errors in GPS height produce errors in the vertical pressure change that must increase as the pressure increases.As part of GRUAN data processing for the RS92, Dirksen et al. (2014) developed a scheme to construct the GPH by combining the pressure and height data provided respectively by the pressure sensor and the GPS system.In recent years, the way how pressure and GPH are measured has probably varied among national RS networks, depending on the use of RSs with or without a pressure sensor.When pressure must be indirectly measured, in principle it would suffice to transmit GPH on the GTS leaving the determination of pressure to NWP centres (Ingleby & Edwards, 2015).The old way of measuring GPH is still eligible at stations using GPS RSs with a pressure sensor.This is the case with NOAA's NWS RS stations since the time they began using GPS tracking, to give a representative example that will be used in the experimental part of this paper (Sections 4 and 5).
In alphanumerical reports, the GPH is reported in whole metres at the standard pressure levels up to 700 hPa and in tens of metres at 500 hPa and higher (WMO, 2011, p. 195).The reporting precision of ±5 m at upper layers affects the global uncertainty of GPH observations; see discussion in Ingleby et al. (2016) for the 500-hPa level.In BUFR reports, the GPH is reported in metres at all reported levels (WMO, 2019).

EFFECT OF MOISTURE ON LAYER THICKNESS
It is of interest to investigate whether the effect of moisture on the distance between two given pressure levels can be measured in terms of PWV content within the corresponding layer.For simplicity, air will be assumed to be an ideal gas.Let ΔZ be the GPH thickness of an atmospheric layer bounded by pressures p 1 and p 2 (p 1 > p 2 ); and let ΔZ d be the thickness that the same layer would have in the absence of moisture if its temperature versus pressure profile were the same.From Equation ( 22) but with C = 1 and further using Equation ( 16) with  = 0.622, we have: Hence, the relative increase in layer thickness due to moisture is defined by: Or, according to Equation ( 21) with C = 1 and making proper use of T v ≈ T, quite simply: meaning that  is proportional to the vertical mean specific humidity over the layer.The vertical mean can be approximated by the mass-weighed mean as long as the layer is not too thick, i.e., with (⋅) m denoting mass-weighted vertical mean, and so, Equation ( 24) becomes: where qdp is the PWV content of the layer p 1 -p 2 .Equation ( 26) states that  is proportional to the ratio of WV mass to air mass within the layer.The factor 3 ∕ 5 corrects for the relative difference in molar mass between dry air and WV; 0.608 was rounded to 3 ∕ 5 ( = 0.6) because the approximation made above to ⟨qT⟩ leads to an overestimation of .Such error increases rapidly as the layer thickness increases but is negligible for the layers considered below.A trial calculation for 700-500 hPa under idealised midlatitude conditions [e-folding scale height of 3 km for q (≈2 km for the WV pressure) regardless of its magnitude; T decreasing 6.5 K⋅km −1 ; T = 270 K at 700 hPa] showed that Equation (26) differs from Equation ( 25) by a relative error of less than 1% (≈2% if the constant 0.608 were used in Equation 26).

TA B L E 1
Relative contribution of WV to layer thickness () corresponding to the average PWV in selected layers and latitude regions (cf.Table 2 in Elliott et al., 1994); dry layer thickness (ΔZ d ) assuming ideal behaviour and average temperature conditions, excluding polar regions; and thickness increase due to WV (ΔZ d ).

Layer (hPa)
PWV (kg⋅m Table 1 shows average values of  estimated by Equation ( 26) for four layers within the interval 1000-300 hPa, based on RS-based PWV estimates for the period 1973-1986 given in Elliott et al. (1994).The dry-layer thickness ΔZ d was computed assuming a lapse rate of 6.5 K⋅km −1 on average both over the globe and the equatorial region, while the average temperature at 1000 hPa was taken as 288 K globally and 300 K in the equatorial region.The product ΔZ d gives the average increase in layer thickness due to moisture.Two aspects stand out: (1) the PWV content in the layers 1000-850 and 850-700 hPa accounts respectively for about two-fifths and one-fourth (two-thirds in total) of the increase in thickness of the entire layer 1000-300 hPa due to moisture, both globally and in the equatorial region; and (2) the total increase in thickness in the equatorial region (32 m) is twice the global increase (16 m).For the polar regions, only  was calculated because the mean temperature profile is more difficult to describe; however, since the  values are about four times smaller than the corresponding global values and the height of the 300-hPa level is in the range 7-8 km at polar latitudes, the total increase in thickness due to moisture should be less than 4 m in the polar regions.As expected, Equation ( 26) does not hold for the too-thick 1000-300 hPa layer: from the summation values in Table 1, the corresponding error is 22% and 20% for the global and equatorial cases.
Inverting Equation ( 27) and considering a change from dry to moist air, PWV 1,2 = PWV 1,2 and (ΔZ This equation relates the PWV within a layer to its mean pressure, mean temperature, and geopotential thickness under the assumptions of ideal gas behaviour and hydrostatic equilibrium.However, it only applies with good approximation to relatively thin layers.In addition, the slight decrease of ΔZ expected from deviations from ideal gas behaviour (as seen in Section 2.2.3) has to be considered here if Z and p are measured independently, rather than one being artificially derived from the other using the ideal gas form of the hydrostatic equation.In the next section we will derive an integral equation similar to the above, but precise and including the compressibility factor of air.

DETERMINING PWV FROM PT AND GPH
Section 4.1 presents a formula expressing PWV in terms of pressure and temperature as functions of GPH which, by using independent measurements of pressure (sensor-measured) and GPH (GPS-based), permits, in theory, calculating PWV without relying on humidity measurements.The practical limitations of such a 'dry estimation' of PWV are detailed in Section 4.2.The potential for long-term PWV evaluation and its implications in accessing RS biases are outlined in Section 4.3.

Formulation
Using Equation ( 16), the basic equation for determining PWV, Equation (3), can be written as This equation assumes hydrostatic equilibrium and is independent of whether the air behaves as an ideal gas or not, just like Equation (3).Eliminating T v through the hypsometric equation in differential form for real air, Equation (21), and using g 0 dZ = gdz according to the definition in Equation (20), Equation (29) becomes: The final expression follows from taking g ≈ g 0 in Equation ( 29).The term p∕(R d CT) represents the density that air would have if it were dry and had the same pressure and temperature as the moist air.The first integral represents the columnar mass of air calculated under the above conditions, whereas the second integral is the actual column mass.The factor (1 accounts for the difference in molar mass between dry air and WV.(NB: Henceforth, we will use the approximate form of Equation ( 30) for the reason that it is consistent with the form of Equation ( 3) that is used in practice for calculating the PWV.) As most current RS systems use GPS tracking, the GPH can be calculated from GPS height and latitude using theoretical gravity, that is, approximating the geoid by an ellipsoid of revolution.Since the GPS refers to the World Geodetic System 1984 (WGS-84), it is adequate to use the WGS-84 Ellipsoidal Gravity Formula with free-air correction (see Appendix A).So, in theory, the PWV can be determined from the vertical data of PT and GPS height.Although C is slightly dependent on RH, notably for positive temperatures on the Celsius scale occurring in the lower troposphere (see Figure 2 and Table C2), it can be approximated so as to optimise the dependence on PT.
At this point it is worth comparing Equation (30) with Equation (3).They are physically equivalent but may lead to different results.Because of measuring errors or accelerated vertical motions, the measured values of p, T, U and Z do not exactly obey the hydrostatic equilibrium as given by Equation ( 21) along with T v = T v (p, T, U) given by Equation ( 17), unless one of the variables is derived from the others using the hydrostatic equation.Namely, if we use either i. sensor-based p and PTU-derived Z, or ii.GPS-based Z and ZTU-derived p, and the vertical resolution is the same for PTU and GPH data, then, regardless of any errors in RS measurements, Equation (30) should give the same result as Equation (3) except for minor differences between the SVP equation used in the calculations of T v by the RS system software and the one used in the calculation of specific humidity by the data user (apart from data uncertainties related to the reporting resolution).Because RS algorithms use the ideal gas approximation, for consistency we should take C = 1 in Equation (30).In either case (i) or (ii) the ground equipment needs RH data to derive Z or p -which means that Equation (30) provides an alternative method for finding the PWV without any advantage over Equation (3) other than computational: computing Equation ( 30) with C = 1 using the {p, T, Z} data as reported by a RS station is much simpler than computing specific humidity from the {p, T, U (/T d )} data and then PWV by Equation (3).(Note: Regarding alphanumeric reports, the only difficulty is that the GPH is not reported at significant levels; regarding historical data, g 0 should have exactly the same value as used in GPH measurements.)Contrastingly, if we apply Equation (30) to independent measurements of pressure and GPH, that is, iii.sensor-based p and GPS-based Z, measurement errors in any of the variables (p, T, U, Z) or departures from hydrostatic equilibrium in the layer of interest will lead inevitably to a different result from that obtained by Equation ( 3).On the one hand, Equation (30) assumes agreement between the GPS-based Z increments (given by Equations ( A10)-( A11)) and the hydrostatic Z increments (given by Equation ( 21) with C ≠ 1 duly included).On the other hand, GPS-and PTU-based GPH measurements are subject to errors; and their matching depends on the accuracy of the hydrostatic approximation, which in the present case is not artificially used in measurements.
Determining the PWV by applying Equation ( 30) to direct measurements of pressure and height is appealing, since the method does not depend on humidity measurements, and these are the main source of error in the usual PWV calculation.The drawback is the high sensitivity to small errors in the quantities involved.Taking p sfc ≈ 1000 hPa and p top = 300 hPa, 10 4 kg ⋅ m −2 is three orders of magnitude larger than the average PWV observed in mid-latitudes.Thus, a relative error of only about 0.1% in either of the main terms in brackets in Equation (30) will compromise the calculation of PWV due to catastrophic loss of accuracy.This is explained by the closeness between T v and T, leading to subtractive cancellation in Equation ( 29): 17), and x v varies from about 0.0001% in the coldest air to 5% in tropical near-surface air.
The relevance of the compressibility factor of air can be examined by analysing Equation (30) bearing in mind that the subtraction terms are nearly equal.Assuming that 1 − C ≈ n × 10 -4 with n being a positive integer ≪10 4 , it can easily be shown that C accounts for an increase in PWV 1000-300 hPa of n × 1.2 kg⋅m −2 .Considering the atmospheric temperature distribution over the year, n ≥ 2 throughout the troposphere (see Figure 2).For the annual average conditions of the mid-latitude troposphere, 1 − C is in the range 3-7 × 10 −4 , and so, the error in PWV 1000-300 hPa using Equation (30) with C = 1 must roughly lie between −4 and −8 kg⋅m −2 .During cold winters in mid-latitudes, and even more so in polar regions, the difference 1 − C is significantly larger.We must therefore include C in Equation ( 30), albeit approximated so as not to have to use humidity data: C ≈ C ≡ C(p, T, U), where U(p) represents the average RH above the place of interest.In practice, it suffices to use discrete values of U to distinguish the lower troposphere from the middle troposphere.The error ( C − C) when using U = 60% in the lower troposphere outside the polar regions, where RH affects C the most, is studied in Appendix C.
Upon substituting dZ = (g∕g 0 )dz in Equation ( 30), an accuracy of ±80 mGal (±8 × 10 −4 m⋅s −2 ) in g is required to have an accuracy of ±1 kg⋅m −2 in PWV 1000-300 hPa .Concerning the approximation of g by normal gravity with free-air correction, most free-air anomalies on the surface of the Earth are within a much narrower range than −80 to 80 mGal.Regarding the calculation of gravity, the formulation for meteorological use ignores the undulations of the geoid.As shown in Appendix A, this produces a relative error of less than 0.002% in most regions; the resulting uncertainty in PWV 1000-300 hPa is ±0.25 kg⋅m −2 .

Caveats
Assuming that C and g are estimated with enough accuracy at every point in the air column, any difference between calculating PWV by Equation (3) and by Equation (30) when using direct measurements of both pressure and height have four possible explanations: (i) RH errors affect only Equation ( 3); (ii) GPS height errors affect only Equation ( 30); (iii) PT data are accurate enough for Equation (3) but not for Equation (30); and (iv) deviations from hydrostatic equilibrium affect the result of Equation ( 30) much more than of Equation (3).The sources of error affecting Equation ( 30), (ii) to (iv), are discussed below.
It is questionable if GPS height measurements in RS systems, using civilian GPS receivers with only one GPS frequency, are accurate enough to find the PWV from only the temperature and height of the isobaric surfaces within the air column above a site.The number and the elevation angle of GPS satellites detectable by the station at the time of observation, GPS signal distortions caused by ionospheric disturbances, tropospheric signal delay caused by WV, electromagnetic noise from the sky and the atmosphere, thermal noise in the receiving system, and occasional loss of signal are all factors affecting GPS height measurements.Multipath propagation due to reflections and diffractions of GPS signals may cause particularly noisy or unstable measurements close to the ground in urban environments, within valleys and on mountain slopes.In addition, a systematic error in the measured height relative to MSL may occur if the height of the ground antenna relative to the surface barometer is not properly set due to a pressure error at the antenna.Finally, if the GPS lock is acquired after release of the weather balloon, the software processing system has to use the mean ascent rate to interpolate height between the time of launch and the time when measurements are reliable.GPS height uncertainties between 10 and 20 m (depending on the RS system) have been reported in WMO RS intercomparisons (Jeannet et al., 2008;Nash et al., 2006;Nash et al., 2011).As a general feature, in the troposphere the GPS height plotted against the time of flight shows much more noise than the GPH curves derived from PTU.The 2010 WMO intercomparison study of high-quality RS systems indicated that the systematic differences between simultaneous GPS heights of different GPS RS systems are much larger than the random error of the systems compared to the average (Nash et al., 2011).That said, the calculation of Equation ( 30) under the assumption of p and Z being independently measured will be affected by GPS height errors, unless these are vertically uniform.Noisy GPS height data will translate into large oscillations of Z k = (g k ∕g 0 )(z k+1 − z k ) in high-resolution data, contrasting with the regularity of the RS ascent speed inferred from the barometric change p k = p k+1 − p k .Unless such random oscillations cancel each other out in the numerical computation of the integral in Z of Equation ( 30), the resulting PWV will be strongly affected.
Concerning PT measurements, as noted qualitatively in Section 4.1, small errors introduce large errors in Equation ( 30) if pressure and GPH are independently measured.In the simplest case where the PT errors ( p ,  T ) are independent of height, the PWV error can be estimated as: (32) where  denotes the environmental lapse rate, assumed to be constant with height.Equations ( 31) and ( 32) are deduced from Equation (30) using the hypsometric (differential) Equation ( 21), with the following nuances as to Equation (32), permissible for an error estimate: regarding air as an ideal gas, neglecting moisture and substituting dT = − dZ, T ≈ T sfc (p∕p sfc ) A ; assuming  T ≪ T, 1∕(T +  T ) ≈ (1 −  T ∕T)∕T.In addition, ∕(1 − ) ≈ 5 ∕ 3 with 0.7% relative error.Let p sfc ≈ 1000 hPa, p top = 300 hPa, and T sfc = 288 K,  = 6.5 K⋅km −1 (A = 0.19) to represent average tropospheric conditions over the globe.Then, for instance,  p = ±0.2hPa yields E PWV ≈ ±4.1 kg⋅m −2 ; while  T = ±0.1 K yields E PWV ≈ ∓4.5 kg⋅m −2 .Clearly, a pressure (/temperature) error of a few tenths of a millibar (/degree) may produce a PWV error larger than the PWV being measured.The actual ranges of T sfc and p sfc are of little importance.The above figures are only indicative, firstly because the PT errors vary with height and secondly, with regard to Equation ( 32), the lapse rate may diverge greatly from 6.5 K⋅km −1 , notably in inversion layers and in the lower atmosphere of polar regions.However, they highlight the high sensitivity of Equation ( 30) to PT errors when p and Z are independently measured (Equations 31 and 32 do not apply otherwise).
The atmosphere departs significantly from hydrostatic equilibrium in convective storms, near frontal zones and during vertical mixing in the planetary boundary layer, while deviations from hydrostatic equilibrium are slight in most motions on a scale with meteorological significance (McIlveen, 1992).Given the near cancellation of the two subtraction terms in Equation ( 30) and since only the subtrahend relies on hydrostatic balance, large cancellation errors related to unbalanced vertical motions can occur.The hydrostatic balance within a column bounded by p sfc and p top implies the following equality: dz, p sfc − p top representing the column weight per unit surface area.An imbalance of, say, ±0.5 hPa between the left-and the right-hand sides will introduce an error of ∓8.5 kg⋅m −2 in PWV 1000-300 hPa given by Equation ( 30).According to Equation ( 31), it takes a pressure data error of ±0.4 hPa constant throughout the column to have about the same effect.The above example suggests that it is difficult to distinguish PWV errors due to non-hydrostatic conditions from those resulting from the propagation of pressure data errors.

Potential use in climate studies
Given the catastrophic cancellation described in Section 4.1 it is impractical to apply Equation (30) to find the PWV at a given location and time based on independent measurements of pressure and GPH.As shown in Section 4.2, the problem is ill-conditioned: small errors in PT or GPS height can result in large errors in the output; the formulation itself may lead to large errors because the hydrostatic assumption is inaccurate.This drawback can, however, be useful to detect systematic errors in PTU or GPS height measurements made by GPS RS systems with a pressure sensor by statistically comparing the PWV values obtained by Equations ( 3) and ( 7) with an adequate expression for e S (T) and by Equation ( 30) with adequate expressions for g and C -insofar as random measurement errors and deviations from hydrostatic equilibrium are expected to vanish on average for many observations.Moreover, to the extent that PT and GPS height do not suffer from systematic errors (except for height-independent errors in GPS height), averaging Equation (30) over several years at different stations should allow finding the seasonal mean PWV over the same years without resorting to humidity measurements.For demonstration purposes, three-year seasonal averages were attempted at various locations, as described in the next section.

Data
NOAA's NWS upper-air network covers the continental United States, plus some islands in the Pacific and Caribbean Sea (https://www.weather.gov/upperair/nws_upper).The NOAA NWS implemented the Radiosonde Replacement System (RRS) between 2005 (Facundo, 2004;Facundo & Fitzgibbon, 2005;Fitzgibbon & Bower, 2004) and 2015 to replace the Micro-ART system (https://www .weather.gov/upperair/rrs_overview).The RRS comprises a GPS tracking antenna, 1680-MHz GPS RSs, and its own PC workstation.Besides the standard WMO messages transmitted in real time, the RRS software builds have provided a data product in the BUFR format, including raw data with 1-s resolution save for occasional 1-s data gaps and processed data reported at the normalised interval of once a second.Until 2012 all RRS sites used the Mark IIA RS manufactured by Sippican; this sonde was phased out of NOAA's NWS in early 2015 apparently because of doubtful humidity measurements in dry or cold tropospheric air (Miloshevich et al. [2006] reported that the Mark IIA hygristor 'cannot reliably measure dry conditions (RH < 20%)' and 'becomes unresponsive at a temperature level that varies between −20 and −50 • C).Since then, the NOAA NWS has used two types of RSs vi : the LMS-6 built by Lockheed Martin Sippican (LMS) and the RS92-NGP built by Vaisala, which is the 1680-MHz version of the standard RS92-SGP (NGP stands for 'NWS GPS with Pressure').The LMS-6 uses a silicon sensor (capacitance aneroid cell), a chip thermistor, and a thin-film capacitance sensor; the RS92-NGP (like all RS92), uses a silicon sensor, a capacitive wire sensor, and a twin thin-film capacitance-heated sensor.
The RRS 1-s BUFR Dataset and the BUFR Decoder Program source code were downloaded from the former NOAA National Climatic Data Center (NCDC) ftp server [ftp://ftp.ncdc.noaa.gov/pub/data/ua/rrs-data/;last access: August 15, 2020].Data files with complete years after 2015 (when the RRS was fully implemented and the Mark IIA was no longer in use) were available until 2018, which explains the selection of the years 2016-2018 in this paper.Two sub-datasets were used: 'Processed u & v winds and position data' and 'Processed pressure, temperature and humidity (PTU) data', 'arrived at by applying normalisation, correction, smoothing, outlier removal, and data plausibility checks' (NCDC, 2008).The humidity variable is RH (dewpoint is given in a sub-dataset for standard and significant levels); pressure is the smoothed pressure; and temperature is the corrected temperature (solar correction by the LMS-6; solar and infrared correction by the RS92-NGP; for details and early testing with the Mark IIA, see Bower & Fitzgibbon, 2004).The NOAA NWS does not apply any correction to RH, though there may be some proprietary software that corrects the raw data sent to the RRS software (Ingleby, 2017).
The data were subsampled using a time step of 6 s starting from the balloon release time, for two reasons.First, a 1-s (≈5 m) resolution does not add useful information relative to a 6-s (≈30 m) resolution for estimating the PWV by the standard formula; trial calculations using 1-s versus 6-s resolution showed absolute differences less than 0.1 kg⋅m −2 in thousands of cases.Secondly, a sampling resolution better than 10 s can be regarded as high resolution (Ingleby et al., 2016).With the above in mind and taking into account the noisy character of GPS height data near the ground, roughly in the first minute of balloon flight (or the lowest 300 m), stations and soundings were primarily selected based on the following criteria: 1. Temporal completeness: stations having observations on at least 80% (/60%) of the days in each quarter (/month) of the year from 2016 through 2018.2. Vertical completeness: observations of PTU and GPS height with no missing values every 6 s between the time of launch and the time when the 300-hPa level is reached.3. Vertical consistency: pressure decreases with height at every 6 s after the release time; GPS height increases at every 6 s after the first minute of flight.4. Data quality: PTU and GPS data passed NCDC's automatic data quality checks designed for the dataset.5. Data integrity: station coordinates, barometer height, time of launch and flight time are all consistent between the wind-and-GPS-position and PTU subdatasets.
The BUFR dataset in question comprises 77 stations in the period 2016-2018.The nominal Zulu time (UTC) of observations is 0000Z and 1200Z (± 1 h), with extra launches at 0300Z and 1800Z (± 1 h) being present at some stations.The monthly data coverage, defined by the fraction of observation days in a month, ranges from 80.6% (September) to 98.4% (March) on average over 2016-2018 and across all stations.However, the month-by-month data coverage at individual stations varies considerably.Because of this, condition (1) removed 14 stations, including all eight Alaska sites that were part of the dataset.A further station was withdrawn (TUCSON, Arizona) because only 30% of the soundings meet the selection criteria (2)-( 5) together, despite the good data coverage.Of the 62 stations selected (mapped in Figure 3) 56 are located in the contiguous United States (90%); five are on tropical Pacific islands [one of the Hawaiian Islands (Big Island), one of American Samoa's Islands (Tutuila), two of the Caroline Islands (Yap Island, Moen Island), and one of thecMarshall Islands (Majuro)]; and one is on a Caribbean island (Puerto Rico).The aggregate of observations in the period 2016-2018 passing conditions (2)-( 5) represents ≈91% of the soundings in the data for the 62 sites; the data loss of 9% is almost entirely due to condition (3) of vertical consistency.

Computational formulas for PWV
The PWV from the surface to 300 hPa was alternatively calculated as follows.First, by Equation ( 3) with Equations ( 7) and ( 8) and f ≈ 1, written in discrete form: Secondly, by Equation ( 30) with g 0 dZ = d = gdz, written in discrete form: The quantities e s (T k ), ) were calculated respectively using: 1996) Magnus form approximation of SVP: Equation (D1).
• CIPM's formulation for the compressibility factor of air (Picard et al., 2008) but neglecting the enhancement factor and using a two-value RH to represent the low-middle troposphere over the contiguous United States (see Ross & Elliott, 1996, Figure 5): Equations ( C1) and (C2) with f = 1, U = 60% for p > 700 Pa, and U = 40% otherwise.
The lower bound of the summation in Equations ( 33) and ( 34), k = 1, refers to the lowest level in the sounding data (balloon release point); k 300 indicates the level whose pressure is closest to the 300-hPa level (according to data).Since the data were sampled to a resolution of 6 s [t k+1 − t k = 6 s, where t k is the flight time] and the ascent rate is ≈ 5 m⋅s −1 , the vertical distance between the upper bound of integration and the 300-hPa level is ≲15 m for individual observations, vanishing on average.In Equation ( 34) the midpoint approximation p∕T ≈ p∕T was used at each integration subinterval; since C and g are slowly varying functions of height, there is no need to include them in the trapezoidal rule.The change in latitude with the balloon drift is a minor detail, because g changes mainly in the vertical direction along the flight path.
Table 2 presents the values of the constants in Equations ( 33) and ( 34).R * , M v and M d are as recommended by the CIPM, with M d referring to 405 ppm of CO 2 by volume and calculated from Equation (4) in Picard et al. (2008); ϵ ≡ M v ∕M d and R d ≡ R * ∕M d are derived constants.(Note: R d = 287.05J⋅kg −1 ⋅K −1 is specified in WMO, 1988, Suppl. No. 2, April, 1996, while R d = 287.05J⋅kg −1 ⋅K −1 is specified in NOAA 's Federal Meteorological Handbook No. 3: OFCM, 1997.)While the rounded value ε ≈ 0.622 may be safely used in both PWV formulas, the value of R d requires care given the numerical instability of Equation (34).For CO 2 mole fractions in the interval 400-410 ppm, M d = 28.96552± 0.00006 g⋅mol −1 .That interval covers both the annual increase and the seasonal change of the global CO 2 within the period 2016-2108 according to monthly concentrations at Mauna Loa Observatory.vii From the original Equation ( 30) with p top = 300 hPa, it can be inferred that by rounding R d to one (two) decimal place(s) PWV 2 will be in error by about +2 (−0.1) kg⋅m −2 ; the uncertainty in R d due to global variations of CO 2 around the mean causes negligible error (< 0.03 kg⋅m −2 ).

TA B L E 2
Constants used in PWV calculations.3.

Consistency of PTU-and GPS-based GPH
Let ΔZ H 1,k denote the GPH thickness of a columnar layer in hydrostatic balance extending from the bottom level to an upper-level k ≥ 2, as given by Equation ( 22), where T v,i ≡ T v (p i , T i , U i ) is given by Equation ( 17).And let ΔZ G 1,k denote the corresponding GPH thickness based on the normal gravity field and geometric thickness: with ) can be calculated by the closed Equation (A8) for Z G ).The summations in Equations ( 35) and (36) were computed using RS data with a time step of 6 s, as in the PWV formulas Equations ( 33) and (34).To measure the fit between the PTU-and GPS-based GPHs we have defined two quantities.First, the absolute deviation between the alternative measures for the thickness of the layer surface-300 hPa: Second, and more informative, the RMS deviation between ΔZ H 1,k and ΔZ G 1,k over the vertical column: Under hydrostatic conditions and in the absence of measurement errors, both  300 and  col would be nearly zero and limited by the propagation of random errors related with the precision with which the PTU and GPS height data are reported: ±0.5 Pa, ±0.05 K, ±0.05% and ± 0.5 m.Values of  300 or  col of the order of metres are a clear indication of measurement errors or departures from hydrostatic balance.

5.2.3
Data sampling criteria for averaging PWV Given that PWV 2 (Equation 34) is extremely sensitive to pressure or temperature errors and to non-hydrostatic conditions, it may differ a lot from PWV 1 (Equation 33) or even result in a negative value when PWV 1 is small.However, a bias between PWV 1 and PWV 2 on average over many days at the same location can hardly be attributed to persistent non-hydrostatic conditions, but rather to mean biases in RS measurements.Hereafter, the seasonal mean over 2016-2018 is shortly referred to as 'seasonal mean' and indicated by an overbar irrespective of season.Since over 90% of RS stations are in the northern temperate zone [approx.two-thirds (one-fourth) in middle (subtropical) latitudes], except for six tropical stations, 'seasonal' refers to the four meteorological seasons, with December to February marking winter, etc. PWV 2 was compared to PWV 1 for each season and RS station, using three different sampling criteria: a.To ensure a fairly good agreement between the PTUand GPS-based GPHs, while avoiding excessive loss of data, individual observations (sounding profiles) were retained if they met the condition:  300 ≤ 10 m ∧  col ≤ 5 m. b.Observations were retained if they met the condition: min(PWV 1 ) < PWV 2 < max(PWV 1 ), where min and max are absolute extrema observed at the location for 1 year ( 2016), to keep information on moisture to a minimum.c.PWV 1 must be positive and lower than what one would have if the tropospheric column were fully saturated relative to water at the reported temperatures.Thus, observations were retained if they met the condition: 0 < PWV 2 < PWV 1 , with Since the computation of  300 and  col requires RH data, the comparison between PWV 1 and PWV 2 after sampling data by criterion (a) is only intended as proof of concept.Criterion (b) requires knowing the local range of PWV 1 to remove anomalous PWV 2 values due to the numerical instability of Equation (34).Criterion (c) discards unphysical PWV 2 values based on the range of PWV 1 for RH over the range 0%-100% but it does not use any information on humidity.Equation (39) results from substituting U = 100% in Equation ( 33) in addition to simplifying the expression for q and the numerical integration scheme.

An alternative robust formulation for PWV
To demonstrate how Equation (30) can be used without suffering from the ill-condition problem, the following finite-difference version was used with Υ k as in Equation ( 34), it being now understood that both p and Z take on values normally measured by RS -one variable being derived from the other using the ideal gas form of Equation ( 21) with C = 1 (differential hypsometric equation).The above equation differs from Equation (34) in that it does not involve any intermediate calculations.And it is simpler than the usual PWV calculation method represented by Equation ( 33) since it dispenses with humidity calculations, although it depends on the RH measurements implicated in the indirect measurement of either GPH or pressure.
With regard to the data used, {Z k } are the PTU-based GPHs as given in BUFR files but sampled at 6-s intervals as in the other PWV calculations.Equation ( 40) was tested against Equation (33) by computing time series for the year 2018 at two sites having entirely different PWV regimes.

5.3
Results and discussion

Seasonal and annual mean PWV
Figure 4 shows PWV 2 against PWV 1 for the 62 stations mapped in Figure 3 and identified by WBAN (Weather-Bureau-Army-Navy) number viii and geographic coordinates in Table 3, for the four temperate seasons and using the data sampling criteria (a), (b), (c) described in Section 5.2.3.The values of data points of each graph are given in Table S1, along with the mean bias error and mean absolute error of PWV 2 relative to PWV 1 for each station, season, and sample (note: the PWV 1 values differ slightly with the sample size and do not coincide exactly with the values shown in Table 3).For clarity, in the following, the deviation between concurrent station-season estimates (graphically, the vertical deviation between a data point and the line y = x shown in each scatter plot of Figure 4) is denoted as D S,i ≡ PWV S 2,i − PWV S 1,i where i indicates an individual station and S some season.The error metrics used in Table 4 are defined thus: MBE (/RMSE) = mean (/root mean square) of D S,i across the stations for a specified season; 1,i across the stations for a specified season.Since only 18 stations used the RS92 in the period studied and all tropical Pacific stations used the LMS6, no separate analysis was made between stations using the RS92 or the LMS6 (save for annual mean values at continental stations).
The statistics of the results depicted in

STATION
in 42% (70%); the MBE is in the range of 0.3 (MAM) to −0.8 (JJA) kg⋅m −2 ; and the RMSE increases to ≈2 kg⋅m −2 .In panel (c), in which observations leading to PWV 2 outside the acceptable range of PWV 1 (from dry to saturated air at the same temperatures) were excluded: 37% (61%) of the cases; the MBE is in the range of 0.1 (MAM) to −1.4 (JJA) kg⋅m −2 ; the RMSE increases further to ≈2-3 kg⋅m −2 .In all three trials, the largest MBE (negative in sign) is found in JJA.However, this is also when the PWV is largest for all stations -with the sole exception of the one located south of the equator, namely in Pago Pago (WBAN:61075).That is why the MAPE is found to be smallest in JJA and largest in DJF: imposing some consistence between PTU-and GPS-based GPH (case a), it varies from about 3% in JJA to 7% in DJF; without any data on moisture (case c), it varies from about 8% in JJA to 13% in DJF.
The condition  col < 5 m in data sampling (a) ensured an above-average consistency between Z G and Z H , as TA B L E 4 Coefficient of correlation (R) between PWV 1 and PWV 2 and error metrics of PWV 2 across the RS stations (Mean Bias Error, Root-Mean-Squared Error, Mean Absolute Percentage Error) by sampling (a, b, c) and season (DJF, etc.) corresponding to the panels of Figure 4. ) 300 hPa = 0.88 m (mean) (exceeding 4 m if Z H is calculated without the compressibility factor correction).It should be noted that the average moisture contribution to Z H at 300 hPa was calculated to be 14 m, which is in line with Table 1 given that most stations are in midlatitudes and on land.Regarding experiment (a) it was verified that, by imposing ever lower upper limits on  300 and  col , the general agreement between PWV 1 and PWV 2 increased progressively at the cost of severe to total loss of observations, with the critical factor being the threshold value for  col (not shown for the sake of brevity).The difference between panel (a) and panels (b) and (c) of Figure 4 simply demonstrates that PTU and GPS height data are to some extent inconsistent with the hydrostatic assumption underlying both PWV formulas used, Equations ( 33) and (34).

MBE (kg⋅m
To appreciate the systematic deviation of PWV 2 from PWV 1 at different sites, the annual means PWV 1 Moreover, for the group of 39 (/17) stations that use LMS6 (/RS92) sondes, MBE = −0.5 (/0.6) kg⋅m −2 and RMSE = 1.9 (/2.0) kg⋅m −2 respectively -which, if we take the 'dry estimates' as a reference, suggests that the LMS6 (/RS92) tends to overestimate (/underestimate) the PWV based on humidity data, although the RMS deviation between the alternative annual mean estimates is similar in both groups of stations.
Finding the reason why PWV 2 A − PWV 1 A varies significantly among stations is beyond the scope of this paper.However, it is possible to put forward some hypotheses.
On one hand, since the RH is around 60% on average over the tropospheric layers that contribute most to the mid-latitude PWV, a conjectural mean bias of ±5% RH would introduce a percentage bias of approx.±8% in PWV 1 A , representing an absolute bias of approx.± 2 kg⋅m −2 for PWV 1 A around 25.0 kg⋅m −2 (in the middle of the range shown in Figure 5).On the other hand, from Equations ( 31) and ( 32), a pressure (/temperature) bias of only ±0.1 hPa (/±0.05K) introduces a bias in PWV 2 A of ±2 kg⋅m −2 .The actual systematic errors in RH and temperature depend not only on the sensors used and the sonde design but also on the mean weather conditions above each station.For example, Figure 5 suggests an average overestimation of PWV 1 at island stations in the tropical Pacific, which may be due to a wet bias in RH measurements by the LMS6 in heavily clouded environments.As to the systematic pressure error, it depends on the correction between the station barometer height and the height of the balloon release point.In addition, some stations are prone to worse measurements of GPS height because of GPS signal obstruction or multipath interference caused by topography and buildings.Finally, the accuracy of the hydrostatic approximation, even if time-averaged, depends on local weather disturbances.Whatever the reason, evaluating the agreement station by station between long-term PWV estimates based on PTU data versus PT and GPS data, exemplified here for only 3 years, allows the detection of inconsistencies in GPS RS data vis-à-vis the hydrostatic balance, provided that a pressure sensor is used.

5.3.2
PWV time series at two example stations Alaska; second, the winter-summer increase in PWV by a factor of 3 (cf.Table 3), exemplifies the midlatitude continental PWV regime defined by Gaffen et al. (1992) in terms of monthly means, although the intermonthly variability seen here is striking even on a five-day time scale (centred moving average in Figure 6).Regarding the error PWV 2 − PWV 1 , the mean bias error (MBE) over 2018 is 0.4 kg⋅m −2 ; note that the negative trend of the error observed in the first half of 2018 (linear fit) is reversed in the second half (not shown).However, PWV 2 shows large and noisy deviations from PWV 1 in spite of the removal of unphysical PWV 2 values by using sampling (c).Clearly, Equation ( 34) is unsuitable for estimating the PWV at a given time and place, as predicted in Section 4. Reinforcing this result, the mean absolute error (MAE) of PWV 2 from PWV 1 as evaluated in each season over 2016-2018 varies in the range of 5-10 kg⋅m −2 for most of the other stations studied; for WILMINGTON, the winter (/summer) MAE is 4.7 (/9.0) kg⋅m −2 (see Table S1, column MAE, rows marked with 'c').
Figure 7 gives a second example, referring to PAGO PAGO INTL ARPT (WBAN:61075) located in American Samoa's main island (Tutuila).The large PWV and small fractional seasonal variation (cf.Table 3) is characteristic of a tropical oceanic PWV regime (Gaffen et al., 1992).The decrease in PWV from January to June is consistent with the slight decrease in surface air temperature between the wet and dry seasons in the South Pacific Ocean.The MBE of PWV 2 relative to PWV 1 over 2018 is −5.5 kg⋅m −2 , much larger than for the previous case.The range of daily deviations is considerable larger, although smaller in percentage.
Regarding the error PWV 3 − PWV 1 in the WILMING-TON series in Figure 6, the MBE is approx.−0.3 kg⋅m −2 both in the first and second half and throughout the year 2018; the linear-fit error trend of January-June is reversed in July-December (not shown).More importantly, the MAE is approx.0.3 kg⋅m −2 over the midyear shown and throughout the year.This small MAE is simply because Equation ( 40 that is, hydrostatically adjusted to PTU measurements.As for PAGO PAGO, the agreement between PWV 3 andPWV 1 is better than in the previous case: MBE ≈ 0.03 kg⋅m −2 and MAE ≈ 0.2 kg⋅m −2 .In both cases, the small systematic error of PWV 3 from PWV 1 (MAE) [linear fit line in bottom panels of Figures 6 and 7] is likely due to dissimilarities between the formulation for SVP of moist air in the LMS6 algorithm for computing the GPH (affecting PWV 3 ) and the one used here to compute the specific humidity (affecting PWV 1 ), being more pronounced in cold weather.The daily fluctuations of PWV 3 − PWV 1 can be explained by random errors related to processing of data in the BUFR dataset or roundoff errors in GPH (reported with precision in metres).
It was found that Equation (40) gives similar results if Z is calculated from the GPS position using normal gravity and p is then derived from Z, T and U, simulating GPS RS systems without a pressure sensor (not shown for the sake of brevity).So, Equation (40) serves as an alternative PWV calculation method to the standard method of Equation (33) whenever the pressure and GPH are not independently measured, no matter which of these variables is measured directly.If pressure is measured by sensor, the method is applicable irrespective of the tracking system.Although the accuracy of both Equations ( 33) and (40) depend on the accuracy of the PTU data, Equation ( 40) avoids moisture calculations by data users after RS observations.

CONCLUSIONS
Much of this study was concerned with discussing the equations widely used for calculating the PWV in a vertical column and the GPH of pressure levels, from which the following conclusions can be drawn: • The approximation of the acceleration of gravity to standard gravity, equivalent to approximating the altitude by GPH, introduces a relative error in PWV ranging from about −0.4% at the equator to 0.2% at the poles.
• Neglecting the WV enhancement factor for saturated air with respect to water underestimates the midlatitude PWV by typically 0.2% to 0.4% of the PWV magnitude.
The virtual temperature used in PTU-based GPH measurements is only affected by a relative error of less than 1/500 of the mole fraction of WV.
• A temperature error of ±0.5 • C yields at most a relative error of ±7% in the calculated specific humidity at temperatures above −70 • C; on global average, the error varies from ±3% at the surface to ±5% at the 300-hPa level.By comparison, a 5% RH error introduces on global average a relative error in specific humidity estimations ranging from about 7% at the surface to 14% at 300 hPa.
• The decrease in air density due to moisture content (relative to dry air density at the same PT conditions) will be greater or less than the increase in air density due to deviations from ideal gas behaviour (relative to the ideal gas density) depending on whether the relative difference between virtual temperature and temperature (approx.two-fifths of the mole fraction of WV) is greater or less than 1 − C where C is the compressibility factor of air.At low latitudes moisture content is the dominant factor, whereas at high latitudes deviations from ideal gas behaviour are more important than moisture.
• In the middle latitudes and under average temperature conditions, the ideal gas approximation (C = 1) entails a bias between 3 and 7 m in PTU-based GPH at the 300-hPa level.During cold winters in temperate latitudes and notably in the frigid zone, as 1 − C increases significantly, the overestimation of GPH should increase accordingly.
• At places where the geoid undulation (/free-air gravity anomaly) is lower than 60 m (/80 mGal) in magnitude, standard GPS-based GPH measurements suffer by an error magnitude of less than 0.02 m (/0.08 m) per km.
The accuracy of GPS-based GPH closely follows that of GPS height.
The hypothesis of finding the PWV from the observed pressure and temperature at different heights, seeking a relationship between the thickness of an atmospheric layer and its PWV content, is supported by the following evidence: • The fraction of the thickness of an elemental pressure layer due to humidity, taking as reference the layer thickness of dry air having the same temperature versus pressure profile, is proportional to the ratio between the PWV content within the layer and the pressure difference between its lower and upper bounds.
• Under ISA temperature conditions, an increase of 1 kg⋅m −2 in the PWV content of a pressure layer around the 925-Pa (/400-hPa) level should increase its thickness by approx.0.5 (/1.0) m.
• Based on layered estimates of PWV by Elliott et al. (1994) for the period 1973-1986, moisture accounts for approx.16 (/32) m of the global (/equatorial) average thickness of the 1000-300 hPa layer assuming a lapse rate of ≈6.5 K⋅km −1 , with the layers below the 700-hPa level accounting for two-thirds of the total moisture-related thickness.In the polar regions, the average moisture contribution to the 1000-300 hPa thickness reduces to less than 4 m.
A general formula for determining PWV from PT and GPH was presented (Equation 30).A 'dry estimation' of PWV on short time scales using sensor-based pressure and GPS-based GPH (turning Equation 30 into Equation 34) requires that the GPS-and PTU-based GPHs agree much beyond what is found in individual observations, because the method becomes ill-conditioned with respect to data errors and deviations from hydrostatic equilibrium.The feasibility of finding the average seasonal cycle of PWV above a station location without humidity data was studied by comparing three-month mean PWVs over 2016-2018 at 62 stations of NOAA's NWS RS network excluding Alaska (90% of them located in the contiguous United States), calculated from the surface to 300 hPa by the usual method (PWV 1 ≡ Equation 33) and the method proposed (PWV 2 ≡ Equation 34).The statistical measures of the deviations of PWV 2 with respect to PWV 1 across the stations are summarised thus: • Excluding (a) observations for which the PTU-based and GPS-based GPHs differ by more than 5 m in terms of RMS deviation across the column or more than 10 m in absolute deviation at 300 hPa: RMSE < 1 kg⋅m −2 in every season.
• Excluding (b) observation leading to PWV 2 values outside the range of PWV 1 as observed throughout the year 2016: MBE <1 kg⋅m −2 in magnitude and RMSE ≈ 2 kg⋅m −2 in all seasons.
The numerical instability of Equation (34) explains the large errors observed in individual PWV 2 estimates with respect to PWV 1 , with seasonal MAEs in the range of 5-10 kg⋅m −2 at most of the stations studied.Interestingly, by using as input the PTU-based GPH as measured and reported by RS, together with PT, Equation (30) with C = 1 (PWV 3 ≡ Equation 40) proved to be practically equivalent to the usual PWV formula, but simpler -although it depends on RH measurements, thereby not being a 'dry estimation'.This is attested by the following results for two stations with contrasting PWV regimes: The usefulness of a 'dry estimate' of PWV, taking advantage of having independent pressure and height measurements, should be seen in light of the time scale.Regarding long-term average PWV based on PTU data or PT and GPS data, the reason why some stations show better agreement than others requires further investigation, since it involves both systematic PTU errors and the average accuracy of the hydrostatic approximation.Understanding where the problem lies would be helpful in assessing the accuracy of most current RS observations, as needed to study long-term changes in temperature and humidity-related quantities.In order to ascertain whether the proposed method is capable of providing better estimates of average PWV than the usual method, comparisons with PWV measurements obtained with frost-point hygrometers may be necessary, albeit limited in time.It would be interesting to compare the two calculation methods in detail, by vertical layer (layered PWV) and latitude region.Regarding daily PWV estimation, the requirement for high accuracy of PTU and GPS height measurements to achieve reasonable agreement between the usual and 'dry' PWV estimates might be useful to find reliable baseline estimates for the validation of remote PWV sensing techniques.vation data used in this study, as archived by NOAA's former National Climatic Data Center.Three anonymous reviewers are thanked for their careful comments that improved the article.The research work was supported by the Consellería de Cultura, Educación e Universidade of the Regional Government of Galicia under the programme Consolidación e Estructuración de Unidades de Investigación Competitivas.Funding for open access publication was provided by the Universidade de Vigo/CISUG.

DATA AVAILABILITY STATEMENT
The radiosonde data supporting the findings of this study were obtained from NOAA's former National Climatic Data Center.The related raw data are currently available upon request from NOAA's National Centers for Environmental Information (https://www.ncei.noaa.gov/access/metadata/landing-page/bin/iso?id=gov.noaa.ncdc:C00966).

ENDNOTES
i The sixth revised edition of the SMT edited by Robert J. List reflected "the improvements in meteorological practice and the increase in accuracy of measurement of physical constants" which had "taken place since the publication of the fifth edition in 1931" (Bull, 1952); above all, it reflected the "increased interest and work with upper-air" (Sherman, 1952).This paper uses the fourth reprint issued in 1968 -with corrections to the original printing and the first reprint issued in 1958 (unchanged since the second reprint issued in 1963).
ii The classical definition of RH over water as the ratio of actual to saturated vapour pressure was used in the SMT from the first (1893) to the fifth (1931) edition.It was also adopted by the International Meteorological Organization (I.M.O.; 1873-1951) until 1946 and, subsequently, by the WMO in 1953 -removing a transitory definition proposed in 1947 by the I.M.O.Conference of Directors at Washington D.C., which had substituted mixing ratio for vapour pressure (Bilham (1950); WMO (1953), pp. 56-57;WMO (1966), Introduction 4.20-2).Unfortunately, the Washington definition was embodied in the sixth edition of the SMT published in the meantime (List, 1949), although the SMT moisture equations are flawless.Even today, it is sometimes mistakenly described as the WMO definition of RH.For the numerical difference between the two definitions see Bilham (1950) and Elliott and Gaffen (1991).
iii WV profiles in the upper troposphere/lower stratosphere are obtained by research RSs equipped with a frost-point hygrometer or, on a larger scale, derived from GPS radio-occultation.Despite their insignificance in relation to PWV, these measurements are critical for studying the coupled radiative-dynamic processes that shape the tropopause and their impact on the global climate system.iv In thermodynamics this quantity is always denoted by Z.Here we adopt the notation C from the SMT 6th ed.(List, 1949), since the letter Z was reserved there for altitude and has been used for geopotential altitude since the time this quantity was introduced in the US Standard Atmosphere, 1962 (Sissenwine et al., 1962).v The relation between water vapour mixing ratio and water vapour pressure and the expression for virtual temperature using mixing ratio instead of specific humidity were widely used in meteorology, as given in the SMT 6th edition (List, 1949) and International Meteorological Tables (WMO, 1966).vi apart from Vaisala's RS41 and the Graw DFM-17 developed by Graw/Noris Group GmbH, lately introduced at some US stations.a second-order Taylor expansion for the vertical gradient of normal gravity potential about h = 0, viz., with Where f e is the flattening of the ellipsoid,  is the angular velocity of the Earth and GM is the geocentric gravitational constant (Heiskanen & Moritz, 1967).The WGS84 Ellipsoid parameters required for Equation (A2) are all given in Table A1 (like a, 1∕f e is a defining constant; the other two defining constants,  and GM, are omitted).Equation (A2) deteriorates at high altitudes because the equipotential surfaces become less oblate with increasing distance from the Earth.But it is very accurate for meteorological purposes: the difference between  h and the accurate model gravity of the WGS84 Ellipsoid is less than 1 μGal (10 −8 m⋅s −2 ) for geodetic heights up to 20,000 m (NIMA, 2000).To first order, ) where 2∕a ≈0.31 mGal m −1 is known in gravimetry as free-air gradient.
The geoid -natural equipotential surface coincident with the MSL over the oceans after removing the effects of tides, winds, and currents -fluctuates around the reference ellipsoid by distances up to about 105 m.Over the oceans the deviation between the geoid and the ellipsoid (geoid undulation) is presently mapped by means of radar satellite altimetry.The global geoid undulation is modelled by a spherical harmonic representation of the Earth's gravitational potential, integrating data from satellite altimetry and gravity surveys.GPS height determinations are relative to the WGS84 ellipsoid, but the end user needs to know the height relative to the geoid, or MSL extended through the continents (orthometric height).GPS receivers convert the geodetic height information to orthometric height using a geoid model; RS systems use the EGM96 geoid (Ingleby, 2017;Miloshevich et al., 2012).Conceptually, the orthometric height (H) is measured along the vertical direction defined by the true for the purpose of defining the kilogram-force.It was later adopted in meteorology as the reference for the correction of mercury barometer readings to standard gravity, figuring for that purpose in the fourth (1918) and later editions of the SMT (from Marvin & Kimball, 1918to List, 1949).[In the first edition (Guyot et al., 1893), the 'standard gravity' (sic) assumed in the gravity correction of barometric readings was 'that prevailing at latitude 45 • and sea level', without risking a numerical value; this was needless because the gravity corrections for latitude and height were ingeniously determined by the ellipsoidal shape of the Earth.]As explained in International Meteorological Tables (WMO, 1966), a gravity reference for 'graduating the scales of mercury barometers or fixing the conversion factors for millibars and normal millimetres of mercury' 'cannot be related to the measured or theoretical value of the acceleration of gravity in specified conditions, for example sea-level at latitude 45 • , because such values are likely to change as new experimental data become available'.In many physical applications, CIPM's standard gravity has been assumed in the conversion from millimetres of mercury at standard temperature 0 • C to pressure units: 760 mm Hg = 1013.25 mb (mb ≡ hPa).
Standard gravity has been used in several 'Standard Atmosphere' models -stylized vertical descriptions of temperature, pressure, and density in a water-vapour-free atmosphere on an average day and over a wide latitude range -beginning with the US National Advisory Committee for Aeronautics (NACA) Standard Atmosphere officially adopted in 1925 (Diehl, 1926).In this model, replacing standard gravity for the acceleration of gravity in the hydrostatic equation was the basis of the formulae used to calculate atmospheric pressure at different altitudes assuming a linear decrease of temperature with height and standard values of temperature and pressure at MSL (15 • C; 760 mm Hg).The International Civil Aviation Organization (ICAO) Standard Atmosphere, 1952, used a geopotential unit based on standard gravity: the standard geopotential metre (9.80665 m 2 s −2 ); for practical reasons, the pressure-related geopotential expressed in these units was considered to represent altitude in metres.In the US Standard Atmosphere, 1962 -a joint product of Air Force Cambridge Research Laboratories, NASA and the US Weather Bureau -the 'geopotential altitude' was defined as geopotential divided by standard gravity, thus having the dimensions of a length (Sissenwine et al., 1962).Eventually, this quantity became conventional in meteorology (WMO, 1988, Appendix A inserted in 1996) in reporting the geopotential at standard pressure levels, as a replacement for the previous use of the geopotential metre unit (9.8 m 2 s −2 ) as in List (1949) and WMO (1966).It is now commonly called geopotential height, while in aeronautics the original name is preferred.
Standard gravity differs slightly from the zonal mean gravity acceleration at latitude ±45 • and at MSL as estimated by normal gravity on the reference ellipsoid (()) at the same latitude, viz.,  (±45 • ) ≈ 9.806 m⋅s −2 with an accuracy to three decimal places, according to the International Gravity Formula 1930 and subsequent revisions (described, e.g., in Li & Götze, 2001) up to the WGS84 Ellipsoidal Gravity Formula (NIMA, 2000).The global area average  ≈ 9.798 m⋅s −2 is lower than  (±45 • ) because half of the globe is located at latitudes below 30 • , where the centrifugal component of gravity is most pronounced.Interestingly, the value 9.8 m⋅s −2 is four times closer to  than 9.80665 m⋅s −2 .But the latter represents better the observed gravity around MSL in midlatitudes, being equal to WGS84  at a latitude of ±45.542 • (Mahoney, 2008).

APPENDIX C. COMPRESSIBILITY FACTOR OF AIR
In this paper the compressibility factor of air (denoted by C) was calculated by an interpolating formula that is part of the CIPM-81/91 and CIPM-2007 formulae for determining the density of moist air (Davis, 1992;Picard et al., 2008), viz., ) , (C1) with the constants given in Table C1, where p is the atmospheric pressure in Pa, T is the absolute air temperature, t = T − 273.15 is the corresponding Celsius temperature, and x v is the mole fraction of WV in air, given by where U is the RH in per cent, e s is the SVP of pure WV over a plane surface of water and f is the enhancement factor.The stated PT ranges for Equation (C1) with CIPM's formulas for f and e s to be used in Equation (C2) are: 600 hPa ≤ p ≤ 1100 hPa and 15  (Note: f is not a monotonic function with respect to t.).Buck (1981) presented a fitting equation to the Hyland data at 1000, 500 and 250 hPa over wide and meteorologically interesting temperature ranges, although it depends only on pressure since this is the main driving factor in the Earth's atmosphere.Likewise, AE96 developed an approximation dependent only on pressure but based on the Goff-Gratch data.From Buck's (1981), Equation [7c]) (quoted as Equation 19 in AE96): f − 1 ≈ 0.42%, 0.24%, and 0.16% at 1000, 500 and 250 hPa respectively.Using Hyland's (1975) Table B, similar values can be found for year-round, mid-latitude temperatures represented by the US Standard Atmosphere, 1976: f − 1 ≈ 0.39% at 1000 hPa (t ≈ 14 • C); 0.22% at 500 hPa (t ≈ −21 • C); and 0.14% at 250 hPa (t ≈ −52 • C).Concerning the calculation of WV pressure e ′ = f (p, t) × e s (t) × U from RS data, the relative uncertainty due to the propagation of systematic and random errors of temperature and RH measurements is much larger than f − 1.For this reason, the enhancement factor may be ignored in routine upper-air calculations of WV pressure and hence of WV mole fraction, WV mixing ratio, specific humidity and dewpoint.

g 0 =
9.80665 m⋅s −2 R * = 8.314472 J⋅mol −1 ⋅K −1 M d = 28.965(52)g⋅mol −1 R d = 287.04(72)J⋅kg −1 ⋅K −1 M v = 18.01528 g⋅mol −1 ε = 0.62195(5) Note: Figures in parentheses refer to 405 ppmv of CO 2 .F I G U R E 3 Locations of NOAA/NWS upper-air stations used in this study, and sondes used throughout 2016-2018.WBAN (Weather-Bureau-Army-Navy) numbers, geographic coordinates and elevations are given in Table Figure 4 are as follows: in panel (a), in which the averaging is restricted to the soundings meeting the condition  300 ≤ 10 m ∧  col ≤ 5 m: | | D S,i | | < 1 (2) kg⋅m −2 in 75% (96%) of the 62× 4 cases; the MBE varies between 0.2 kg⋅m −2 in DJF and −0.4 kg⋅m −2 in JJA; and the RMSE is less than .F I G U R E 4 Scatter diagrams of PWV 2 versus PWV 1 for the stations indicated in Figure 3 and Table 3 (one point per station), obtained TA B L E 3 Identification of the stations shown in Figure 3 and seasonally averaged PWVs over 2016-2018 calculated by Equation (33) on the basis of the observations meeting conditions (2)-(5) of Section 5.1.
2018 are plotted against each other in Figure 5; note that PWV 2 A was calculated without knowing anything about moisture content, similarly to panel (c) of Figure 4.At the five Pacific Island stations, all using LMS6 sondes, PWV 2 A is lower than PWV 1 A by 3-9 kg⋅m −2 depending on site.Among the 56 stations in the contiguous United States, the number of cases with positive or negative bias is about the same and | kg⋅m −2 in 25 (/40) sites.

Figure 6 F
Figure6shows a 180-day time series of PWV 1 (≡ Equation33), along with the deviations of PWV 2 and PWV 3 (≡ Equations 34 and 40 respectively) from PWV 1 , for WILMINGTON station (WBAN:13841), located at Wilmington Air Park, Ohio.This station was chosen for two reasons: first, its annual mean PWV is in the middle of the range of values observed at US stations excluding Figure7gives a second example, referring to PAGO PAGO INTL ARPT (WBAN:61075) located in American Samoa's main island (Tutuila).The large PWV and small fractional seasonal variation (cf.Table3) is characteristic of a tropical oceanic PWV regime(Gaffen et al., 1992).The decrease in PWV from January to June is consistent with the slight decrease in surface air temperature between the wet and dry seasons in the South Pacific Ocean.The MBE of PWV 2 relative to PWV 1 over 2018 is −5.5 kg⋅m −2 , much larger than for the previous case.The range of daily deviations is considerable larger, although smaller in percentage.Regarding the error PWV 3 − PWV 1 in the WILMING-TON series in Figure6, the MBE is approx.−0.3 kg⋅m −2 both in the first and second half and throughout the year 2018; the linear-fit error trend of January-June is reversed in July-December (not shown).More importantly, the MAE is approx.0.3 kg⋅m −2 over the midyear shown and throughout the year.This small MAE is simply because Equation (40) uses as input the GPH as measured by RS, As in Figure6, but for PAGO PAGO INTL ARPT (WBAN/WMO-ID: 61075/91765; Sonde: LMS6).
vii https://scrippsco2.ucsd.edu/graphics_gallery/mauna_loa_record/mauna_loa_record.html(last access: February 18, 2023).viii The WBAN number has been used in the United States since the 1950s to identify weather stations.WBAN numbers (alongside WMO index numbers) for NOAA/NWS upper-air stations, station names and country/state information can be found at https://www TA B L E A1 WGS84 Ellipsoid parameters (NIMA, 2000).a = 6378137.0m b = 6356752.3142m  e = 9.7803253359 m ⋅ s −2  p = 9.8321849378 m ⋅ s −2 1∕f e = 298.257223563e 2 = 6.69437999014 × 10 −3 k = 1.93185265241 × 10 −3 m = 3.44978650684 × 10 −3 Standard gravity, defined exactly as 9.80665 m⋅s −2 , was first introduced by the Comité International des Poids et Mesures (CIPM) in 1887 from gravimetric measurements made at Paris and corrected by that time to 45 • latitude based on the ideal figure of the Earth; it became accepted by the Conférence Générale des Poids et Mesures in 1901

−2 ) 𝜶 (cm⋅km −1 ) 𝚫Z d (m) 𝜶𝚫Z d (m)
The WBANs (Weather-Bureau-Army-Navy) underlined indicate tropical stations, one in the Caribbean (asterisked) and five on Pacific Islands.Station locations are within 7-8 km of the coordinates (LAT, LON) indicated.ELEV is the altitude of the station barometer.NOBS is the number of observations in each sample.SONDE is the RS type used throughout 2016-2018.1 kg⋅m −2 in all seasons.In panel (b), in which observations leading to PWV 2 values outside the observed range of PWV 1 over the year 2016 were excluded: | Note: Compressibility factor of moist air in the meteorological range of PTU using Equations (C1)-(C2) with f = 1 and Equation (D1).Relative differences with respect to SMT values are given in parentheses (r C in %).
• C ≤ t ≤ 27 • C.These narrow ranges are related to the requirement for high accuracy in laboratory determination of air density, designed to correctT A B L E C2