VESPA-22: a ground-based microwave spectrometer for long-term measurements of polar stratospheric water vapor

The new ground-based 22 GHz spectrometer, VESPA-22 (water Vapor Emission Spectrometer for Polar Atmosphere at 22 GHz) measures the 22.23 GHz water vapor emission line with a bandwidth of 500 MHz and a frequency resolution of 31 kHz. The integration time for a measurement ranges from 6 to 24 h, depending on season and weather conditions. Water vapor spectra are collected using the beam-switching technique. VESPA-22 is designed to operate automatically with little maintenance; it employs an uncooled front-end characterized by a receiver temperature of about 180 K and its quasi-optical system presents a full width at half maximum of 3.5^∘. Every 30 min VESPA-22 measures also the sky opacity using the tipping curve technique. The instrument calibration is performed automatically by a noise diode; the emission temperature of this element is estimated twice an hour by observing alternatively a black body at ambient temperature and the sky at an elevation of 60^∘. The retrieved profiles obtained inverting 24 h integration spectra present a sensitivity larger than 0.8 from about 25 to 75 km of altitude during winter and from about 30 to 65 km during summer, a vertical resolution from about 12 to 23 km (depending on altitude), and an overall 1σ uncertainty lower than 7 % up to 60 km altitude and rapidly increasing to 20 % at 75 km.

In July 2016, VESPA-22 was installed at the Thule High Arctic Atmospheric Observatory located at Thule Air Base (76.5^∘ N, 68.8^∘ W), Greenland, and it has been operating almost continuously since then. The VESPA-22 water vapor mixing ratio vertical profiles discussed in this work are obtained from 24 h averaged spectra and are compared with version 4.2 of concurrent Aura/Microwave Limb Sounder (MLS) water vapor vertical profiles. In the sensitivity range of VESPA-22 retrievals, the intercomparison from July 2016 to July 2017 between VESPA-22 dataset and Aura/MLS dataset convolved with VESPA-22 averaging kernels shows an average difference within 1.4 % up to 60 km altitude and increasing to about 6 % (0.2 ppmv) at 72 km.

1 Introduction

The polar atmosphere is a very complex system in which water vapor plays an important role. Water vapor has a major impact on the radiative balance affecting both infrared radiation by greenhouse effect and visible radiation by clouds coverage. The importance of water vapor in the Arctic region is enhanced by the so-called Arctic amplification effect (Serreze and Francis, 2006), a positive feedback that links the quantity of water vapor in the atmosphere, the presence of clouds, the ice coverage, and the surface temperature. Additionally, about 10 % of the surface warming measured during the last two decades can be ascribed to stratospheric water vapor, as shown by Solomon et al. (2010).

The characterization of the water vapor profile, particularly in the stratosphere and mesosphere, is important to understand many chemical processes. Polar water vapor is directly involved in the ozone chemistry as the main source of the OH radical in reactions that cause the ozone destruction (Solomon, 1999). It is also related to the formation of polar stratospheric clouds (PSCs) that grant the catalytic surfaces on which heterogeneous reactions take place.

The main sources of middle atmosphere water vapor are transport through the tropical tropopause and methane oxidation, whereas the main sink is photolysis in the Lyman band. In the middle atmosphere, the lifetime of water vapor varies from several days to weeks, making it a valuable tracer for the investigation of dynamical processes that characterize the polar regions especially during winter and spring.

The processes that lead to long-term variations in stratospheric and mesospheric water vapor are not completely understood. Positive trends were observed during the 1980–2000 period (Nedoluha et al., 1999; Rosenlof et al., 2001) and Oltmans et al. (2000) suggested that only one-half of these changes are related to anthropogenic activities.

The comprehension of the change in climate that affects the Arctic and the peculiar characteristics of this region's atmosphere calls for long-term measurements. For this task, ground-based microwave remote sensing is a powerful tool to measure water vapor profiles and total amount of precipitable water vapor (PWV).

In this paper, VESPA-22 (water Vapor Emission Spectrometer for Polar Atmosphere at 22 GHz), a microwave spectrometer developed at the INGV (Istituto Nazionale di Geofisica e Vulcanologia) is presented. In July 2016, during the Study of the water VApor in the polar AtmosPhere (SVAAP) measurements campaign, a key part of a research effort devoted to the study of the impact of water vapor and clouds on the radiation budget at the ground, VESPA-22 was installed at the Thule High Arctic Atmospheric Observatory (THAAO) located at Thule Air Base (76.5^∘ N, 68.8^∘ W), Greenland. The instrument general features, the measurement physics and technique, a comparison between VESPA-22 and the Aura/Microwave Limb Sounder (MLS) (Waters et al., 2006) version 4.2 datasets are discussed in this work.

2 Instrumental setup

VESPA-22 collects the microwave radiation emitted by the water vapor transition at 22.235 GHz with a spectral resolution of 31 kHz and a bandwidth of 500 MHz. From the spectral measurements water vapor profiles can be retrieved with a temporal resolution of 1–4 profiles a day, depending on season and weather conditions. The instrument also measures the sky opacity in few minutes by performing a scan of the sky at various angles called tipping curve.

In Fig. 1, a photo and the layout of the spectrometer are displayed. VESPA-22 collects the sky signal coming from different elevation angles using a parabolic mirror. The signal is reflected by the mirror to a Gaussian choked horn antenna (Teniente et al., 2002). This antenna has a far field directivity of 23.5 dB with a pattern that can be approximated at 99.85 % to a Gaussian beam with a beam waist of 22.4 mm (Bertagnolio et al., 2012). The quasi-optical system, consisting of an antenna and parabolic mirror, has a full width at half maximum (FWHM) ${\mathit{\theta }}_{\mathrm{3}\mathrm{dB}}=\mathrm{3.5}{}^{\circ }$ (Bertagnolio et al., 2012). This high directivity of the system allows the observation at angles as low as 12^∘ above the horizon, therefore maximizing the number of emitting molecules and the signal-to-noise ratio of the measurement.

Figure 1 A photo of VESPA-22 installed at the THAAO, Thule Air Base, Greenland (left), and the layout of the instrument (right).

A first uncooled low noise amplification stage amplifies the incoming signal with a gain of 35 dB. Two noise diodes manufactured by Noisecom are inserted in the radio frequency (RF) chain after the choke antenna by means of two 20 dB broadwall direction couplers. These two diodes produce a signal which is measured to average at about 119 and 78 K and are used to calibrate the sky signal observed by VESPA-22, as described in Sect. 3.3.

After the 20 dB couplers the signal is amplified and down-converted to lower frequencies and eventually sampled at 1 Gs s⁻¹ by a fast Fourier transform spectrometer (FFTS; Acquiris AC240/Agilent U1080 A-001) with 16 384 channels, resulting in a spectral resolution of 31 kHz for a total bandwidth of 500 MHz.

During measurements, the antenna is moved back and forth to change its distance from the mirror by $\mathit{\lambda }/\mathrm{4}=\mathrm{3.34}\mathrm{mm}$ in order to avoid multiple internal reflections of the signal which would cause the formation of standing waves affecting the observed spectrum. One full spectrum is obtained by averaging together two 6 min spectra collected with the antenna positioned at these two different distances from the mirror. Once a day, two photodiodes are employed to check and reset the correct distance between the antenna and the mirror.

VESPA-22 is installed indoors to protect it from the strong winds and storms which can occur during most of the year. The indoor installation prevents the deposition on the quasi-optical system of snow in winter and dust in summer, as well as the condensation of water droplets, therefore improving the durability of the equipment. Additionally, the parabolic mirror and its driving motor are not exposed to strong winds, and VESPA-22 is therefore characterized by a pointing offset that is very stable with time. The instrument is located in a small wooden annex (Fig. 2) to the main observatory in order to minimize the presence of metal surfaces which could also yield standing waves in the observed spectrum.

Figure 2 A photo of the exterior of the wooden annex hosting VESPA-22. The observing window of the signal beam is visible on the side of the annex.

The spectrometer observes the sky through two 5 cm thick Plastazote LD15 windows, one covering the observation at angles from 10 to 60^∘ above the horizon and a smaller one covering the zenith direction (see Figs. 1 and 2). This material was tested in the laboratory and proved to have a negligible absorption in the microwave region of interest. Furthermore, the windows are not perpendicular to the antenna beam in order to minimize any potential formation of standing waves. The opacity of the LD15 window sheets was estimated during VESPA-22 installation and it was measured to be less than 0.0005 Np. Three powerful fans are employed to blow off the snow from the observing windows, preventing deposition or ice formation on the external side of the windows. VESPA-22 compares the sky emission from the zenith direction (also called reference direction) to the sky emission from an angle close to the horizon (also called signal direction) in order to perform a stratospheric measurement, as will be explained in Sect. 3.2. The zenith is observed through a white Delrin^® acetal homopolymer resin sheet, hereafter simply delrin (see Fig. 1), which adds a grey-body emission to the reference beam; this sheet is set at Brewster's angle with respect to the incident beam in order to minimize reflections and the formation of standing waves.

During an hour of standard instrument operations, about 3 min is dedicated to tipping curves (see Sect. 3.4) used to measure tropospheric opacity, about 35 min is dedicated to measure the emissions from the signal and reference beams, and the rest of the hour is dedicated to instrumental operations, such as the rotation of the parabolic mirror or the instrument calibration by means of the noise diodes.

3 Measurement technique

VESPA-22 collects the 22.235 GHz radiation emitted by water vapor molecules. At this frequency, the Rayleigh–Jeans approximation for the Plank's law can be used and the radiation intensity can be expressed in terms of brightness temperature. The line shape of the emission from a stratospheric altitude z is mainly a function of atmospheric pressure:

$$\begin{array}{}\text{(1)}& {\displaystyle }\mathrm{\Delta }\mathit{\nu }\propto P\left(z\right){\left({\displaystyle \frac{T_{\mathrm{0}}}{T\left(z\right)}}\right)}^n,\end{array}$$

where Δν is the line width, n<1 is a constant coefficient and T₀=300 K (see Table 3, Sect. 4.1, for the spectroscopic values used for VESPA-22 retrievals). Using this dependence and knowing the pressure and temperature atmospheric profiles, the measured spectrum can be inverted using the reverse problem theory to obtain the vertical water vapor mixing ratio profile. In the mesosphere, the Doppler broadening overcomes the pressure broadening of Eq. (1) determining an upper limit to the altitude range in which the water vapor profile can be retrieved. In the lower stratosphere, the broadening makes the line width comparable to the VESPA-22 spectral bandpass, therefore setting a lower altitude limit for the deconvolved profile.

3.1 Atmospheric opacity

The majority of the incoming radiation measured by VESPA-22 comes from the troposphere. The tropospheric emission needs to be filtered in order to study the stratospheric signal. In order to simplify the radiative transfer equation, the following approximations can be adopted (e.g., Nedoluha et al., 1995).

• The troposphere is represented as an isothermal layer absorbing the signal to be measured. The first kilometers of the atmosphere produce the greatest part of this layer emission which has a line width larger than VESPA-22 bandwidth. This contribution is treated as an emission constant in frequency.

• The contribution of the stratospheric water vapor absorption to the opacity τ is small. This approximation was tested by calculating the atmospheric opacity by means of the radiative transfer simulation software ARTS (Eriksson et al., 2011) and using water vapor vertical profiles with and without their stratospheric component. At the frequency of maximum absorption, the contribution of the stratospheric profile to the overall atmospheric opacity is between 2 and 5 % depending on the season. When averaged over the VESPA-22 frequency range, the difference between the opacity calculated using a full water vapor vertical profile and a profile with no water vapor above the tropopause is between 0.4 and 1.5 %.

• The opacity τ_ν can be substituted by its mean value τ. The maximum difference between τ_ν and its mean value τ is between 1.6 and 3.6 % depending on season.

• The only signal coming from outside the atmosphere, T₀, is the cosmic background radiation with a constant brightness temperature of 2.73 K.

With these approximations the radiative transfer equation describing the radiation received at the ground from an elevation angle θ can be written as

$$\begin{array}{}\text{(2)}& {\displaystyle }{T}_{\mathrm{S}}\left(\mathit{\nu }\right)=T_{\mathrm{0}}{e}^{-\mathit{\mu }\mathit{\tau }}+T_{\text{trop}}\left(\mathrm{1}-e^{-\mathit{\mu }\mathit{\tau }}\right)+\mathit{\mu }T\left(\mathit{\nu }\right)e^{-\mathit{\mu }\mathit{\tau }}.\end{array}$$

The left term is the radiation received by the spectrometer when pointing at the elevation angle θ. The first term on the right side is the extra-atmospheric emission, the second term is the solution of the radiative transfer equation for an isothermal domain where the tropospheric emission is indicated with T_trop and represents the mean tropospheric temperature weighted with the water vapor concentration. The third term on the right is the emission coming from the stratosphere and mesosphere, T(ν), proportional to the air mass factor μ and attenuated in the troposphere by a factor e^−μτ. The air mass factor at the elevation angle θ is computed according to the formula presented in the work of de Zafra (1995).

3.2 Beam-switching technique

The technique that VESPA-22 employs to measure the stratospheric signal is called the beam-switching technique or Dicke switching technique (e.g., Parrish et al., 1988). The instrument compares the emission coming from an observation angle θ (signal beam) with a reference signal with the same mean power over the passband. The observation angle depends on the atmospheric opacity and for VESPA-22 it varies from 12 to 25^∘ above the horizon. The reference signal used by VESPA-22 is the sky emission at the zenith (reference beam). In clear-sky conditions, the emission at the zenith is smaller than the emission at a much larger zenith angle. Therefore, in order to ensure that the reference beam has the same mean power of the signal beam, a thin sheet of delrin is inserted in the reference beam. The delrin sheet acts as a grey body so that

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }{T}_{\mathrm{R}}\left(\mathit{\nu }\right)=T_{\mathrm{0}}{e}^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}+T_{\text{trop}}\left(\mathrm{1}-e^{-\mathit{\tau }}\right)e^{-{\mathit{\tau }}_{\mathrm{d}}}\\ \text{(3)}& {\displaystyle }& {\displaystyle }+T\left(\mathit{\nu }\right)e^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}+T_{\mathrm{d}}\left(\mathrm{1}-e^{-{\mathit{\tau }}_{\mathrm{d}}}\right),\end{array}$$

where T_R is the radiation observed by VESPA-22 coming from the zenith, partially absorbed by the delrin sheet, T_d is the physical temperature of the sheet and τ_d its mean opacity value over the spectral passband. In this equation, the air mass factor μ is equal to 1.

During data-taking operations, VESPA-22 alternates reference and signal observations. The instrument constantly checks if the two beams have the same mean power and continuously changes the signal angle to minimize the difference between them. When the two beams have the same intensity, the frequency independent terms of Eqs. (3) and (2) can be equated, obtaining

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }{T}_{\mathrm{0}}{e}^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}+T_{\text{trop}}\left(\mathrm{1}-e^{-\mathit{\tau }}\right)e^{-{\mathit{\tau }}_{\mathrm{d}}}+T_{\mathrm{d}}\left(\mathrm{1}-e^{-{\mathit{\tau }}_{\mathrm{d}}}\right)\\ \text{(4)}& {\displaystyle }& {\displaystyle }\approx T_{\mathrm{0}}{e}^{-\mathit{\mu }\mathit{\tau }}+T_{\text{trop}}\left(\mathrm{1}-e^{-\mathit{\mu }\mathit{\tau }}\right),\end{array}$$

where the stratospheric contribution to the mean beam intensity (about 1 %) is neglected (de Zafra, 1995).

The stratospheric signal T(ν) is obtained by subtracting reference from signal. Using Eqs. (2), (3) and (4) it is possible to write

$$\begin{array}{ll}\text{(5)}& {\displaystyle }& {\displaystyle }{T}_{\mathrm{S}}\left(\mathit{\nu }\right)-T_{\mathrm{R}}\left(\mathit{\nu }\right)\approx T\left(\mathit{\nu }\right)\left(\mathit{\mu }{e}^{-\mathit{\mu }\mathit{\tau }}-e^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}\right){\displaystyle }& {\displaystyle }\text{and}\\ \text{(6)}& {\displaystyle }& {\displaystyle }T\left(\mathit{\nu }\right)={\displaystyle \frac{T_{\mathrm{S}}\left(\mathit{\nu }\right)-T_{\mathrm{R}}\left(\mathit{\nu }\right)}{\mathit{\mu }{e}^{-\mathit{\mu }\mathit{\tau }}-e^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}}}.\end{array}$$

Three delrin sheets with different thicknesses (3, 5, and 9 mm) and opacities can be employed, depending on the season, in order to maintain the signal angle between 12 and 25^∘ above the horizon. Spectra collected are smoothed using a 50-channel moving average. This smoothing process is, however, not performed in a 6 MHz interval centered around the emission line to maintain the maximum frequency resolution near the line peak.

In order to estimate τ_d, the mean opacity value of delrin over the spectral passband, during normal data-taking operations the signal angle is locked to its balanced position and the delrin sheet is removed from the reference beam. τ_d can then be calculated using

$$\begin{array}{}\text{(7)}& {\displaystyle }{\mathit{\tau }}_{\mathrm{d}}=-\ln \left({\displaystyle \frac{T_{\mathrm{d}}-{\stackrel{\mathrm{‾}}{T}}_{\mathrm{S}}}{T_{\mathrm{d}}-{\stackrel{\mathrm{‾}}{T}}_{{\mathrm{R}}_{\text{nod}}}}}\right),\end{array}$$

where ${\stackrel{\mathrm{‾}}{T}}_{{\mathrm{R}}_{\text{nod}}}$ is the mean intensity of the reference beam without the sheet and ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{S}}$ is the mean value of the signal beam. The brightness temperatures ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{S}}$ and ${\stackrel{\mathrm{‾}}{T}}_{{\mathrm{R}}_{\text{nod}}}$ are calculated using parameters obtained by means of a calibration performed using liquid nitrogen (LN₂; see Sect. 3.3) which is always carried out before estimating τ_d. In order to run these operations, qualified personnel must be at the observatory and therefore τ_d was measured only in July and November 2016 and February 2017 (Table 1). A total of seven measurements of the compensating sheet opacity were performed, two in July, three in November, and two in February 2017. Half of the difference between the minimum and maximum τ_d  values obtained during the same period is used as measurement uncertainty. The mean value of the delrin opacity changes with time, possibly due to a certain level of degassing, i.e., the property of absorbing/releasing water vapor molecules from/to the environment, which depends on atmospheric humidity and is often noticed in plastic materials. During winter, as the air is drier, the compensating sheets appear to release some water vapor and lower their opacities. VESPA-22 spectral data are calibrated using a linear interpolation between the τ_d values measured over time.

3.3 Calibration scheme

The broad-band response of VESPA-22 to the signal coming from the sky can be written as

$$\begin{array}{}\text{(8)}& {\displaystyle }V\left(\mathit{\nu }\right)=\mathit{\alpha }\left(\mathit{\nu }\right)\left(T\left(\mathit{\nu }\right)+T_{\text{rec}}\left(\mathit{\nu }\right)\right)+V_{\mathrm{0}}\left(\mathit{\nu }\right).\end{array}$$

In this equation, V is the incoming signal in count numbers of the FFT back-end spectrometer, T is the signal brightness temperature, α the gain of the instrument, T_rec the receiver noise temperature and V₀ the “zero” signal of the FFTS. All the quantities represented in Eq. (8) are frequency dependent. The “zero” signal, which amounts to approximately 0.5 % of the incoming signal to the FFTS, is measured approximately every 15 min and it is subtracted to each signal and reference 15 min integration spectra which are eventually saved on the control and acquisition PC (see Fig. 1). Using Eq. (8), Eq. (6) can be written as

$$\begin{array}{}\text{(9)}& {\displaystyle }T\left(\mathit{\nu }\right)={\displaystyle \frac{\mathrm{1}}{\mathit{\alpha }\left(\mathit{\nu }\right)}}{\displaystyle \frac{V_{\mathrm{S}}\left(\mathit{\nu }\right)-V_{\mathrm{R}}\left(\mathit{\nu }\right)}{\mathit{\mu }{e}^{-\mathit{\mu }\mathit{\tau }}-e^{-\mathit{\tau }-{\mathit{\tau }}_{\mathrm{d}}}}}.\end{array}$$

VESPA-22 measures the gain parameter α by means of the noise diode. As mentioned before, VESPA-22 has two different noise diodes, one used to perform regular calibrations and the second inserted to ensure the stability over time of the first one, as described by Gomez et al. (2012). During regular measurements, the calibration noise diode is switched on and its emission is added to the reference signal. The noise diode is then switched off and VESPA-22 measures only the radiation coming from the zenith. The gain parameter can be obtained by subtracting these two measurements according to

$$\begin{array}{}\text{(10)}& {\displaystyle }\mathit{\alpha }\left(\mathit{\nu }\right)={\displaystyle \frac{V_{\mathrm{R}+\mathrm{nd}}\left(\mathit{\nu }\right)-V_{\mathrm{R}}\left(\mathit{\nu }\right)}{T_{\mathrm{nd}}\left(\mathit{\nu }\right)}},\end{array}$$

where V_R+nd and V_R are the signals expressed in count numbers measured with the noise diode turned on and off, respectively, and T_nd is the noise diode emission temperature at the various frequencies, estimated during a calibration procedure. The calibration consists in measuring the emission from two sources at two different known emission temperatures. The first source is a black body at ambient temperature made with an Eccosorb CV-3 panel by Emerson and Cuming; the second source is mainly the sky at an angle of 60^∘ above the horizon. Every 3 or 4 months, approximately, the sky emission can be replaced by a second CV-3 panel immersed in LN₂. The use of LN₂ likely grants more accurate results, as the physical temperature of the emitting body has a smaller uncertainty with respect to the estimated sky temperature at 60^∘, but the LN₂ calibration cannot be performed automatically by the instrument and has so far been carried in July and November 2016 and in February 2017.

Table 1 Mean values and SDs (standard deviations) of the measured opacity for the two used delrin sheets during different seasons.

The general calibration equations are described in what follows, whereas the tipping curve procedure that allows the use of the sky as calibration source is explained in Sect. 3.4. Knowing the emission temperature of two sources α, T_rec, and T_nd can be obtained by using the following equations:

$$\begin{array}{ll}\text{(11)}& {\displaystyle }& {\displaystyle }\mathit{\alpha }\left(\mathit{\nu }\right)={\displaystyle \frac{V_{\text{hot}}\left(\mathit{\nu }\right)-V_{\text{cold}}\left(\mathit{\nu }\right)}{T_{\text{hot}}-T_{\text{cold}}}},{\displaystyle }& {\displaystyle }{T}_{\text{rec}}\left(\mathit{\nu }\right)\\ \text{(12)}& {\displaystyle }& {\displaystyle }={\displaystyle \frac{T_{\text{hot}}\left(V_{\text{cold}}\left(\mathit{\nu }\right)-V_{\mathrm{0}}\left(\mathit{\nu }\right)\right)-T_{\text{cold}}\left(V_{\text{hot}}\left(\mathit{\nu }\right)-V_{\mathrm{0}}\left(\mathit{\nu }\right)\right)}{V_{\text{hot}}\left(\mathit{\nu }\right)-V_{\text{cold}}\left(\mathit{\nu }\right)}},{\displaystyle }& {\displaystyle }\text{and}\\ \text{(13)}& {\displaystyle }& {\displaystyle }{T}_{\mathrm{nd}}\left(\mathit{\nu }\right)={\displaystyle \frac{V_{\text{cold + nd}}\left(\mathit{\nu }\right)-V_{\text{cold}}\left(\mathit{\nu }\right)}{\mathit{\alpha }\left(\mathit{\nu }\right)}}.\end{array}$$

T_hot and T_cold and V_hot and V_cold are the emission temperatures and the recorded count numbers of the two sources, respectively. For the black body, the emission temperature is assumed to be equal to its physical temperature.

The noise diode produces a signal that is measured to be quite stable in frequency. In fact, single-channel T_nd values are always within 1.5 % of the spectral mean of the diode temperature brightness ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$. The spectra originated from black-body measurements (especially those from the CV-3 immersed in LN₂) can be affected by standing waves and in order to avoid to transfer them in the calibrated sky spectral measurements, T_nd is averaged over the central 11 000 channels of the FFTS, as suggested by Gomez et al. (2012). Therefore, it can be written

$$\begin{array}{}\text{(14)}& {\displaystyle }{\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}=\left(T_{\text{hot}}-T_{\text{cold}}\right)\sum _i{\displaystyle \frac{V_{\text{cold + nd}}\left({\mathit{\nu }}_i\right)-V_{\text{cold}}\left({\mathit{\nu }}_i\right)}{V_{\text{hot}}\left({\mathit{\nu }}_i\right)-V_{\text{cold}}\left({\mathit{\nu }}_i\right)}},\end{array}$$

where V_{cold + nd}(ν_i) is the value of channel i with the noise diode turned on, while V_cold(ν_i) is the value of the same channel with the noise diode turned off.

3.4 Tipping curve procedure

In this section, the tipping curve calibration technique employed to calibrate the noise diodes using the sky signal is described. The tipping curve procedure is performed twice every hour as it is also used to measure the atmospheric opacity needed in Eq. (9). During a tipping curve, VESPA-22 collects the radiation coming from different elevation angles, approximately every 5^∘ from 35 to 60^∘ above the horizon. The measured spectra are averaged using the 11 000 central channels of the spectrometer. Radiation from the stratosphere contributes less than 1 % and can be neglected, so the signal intensity can be described by means of the following equation:

$$\begin{array}{}\text{(15)}& {\displaystyle }\stackrel{\mathrm{‾}}{T}\left({\mathit{\theta }}_i\right)\cong T_{\mathrm{0}}{e}^{-\mathit{\mu }\left({\mathit{\theta }}_i\right)\mathit{\tau }}+T_{\text{trop}}\left(\mathrm{1}-e^{-\mathit{\mu }\left({\mathit{\theta }}_i\right)\mathit{\tau }}\right).\end{array}$$

The atmospheric temperature T_trop can be estimated from the surface temperature T_surf:

$$\begin{array}{}\text{(16)}& {\displaystyle }{T}_{\text{trop}}=T_{\text{surf}}-d,\end{array}$$

where the value of d can be affected by seasonal variations. In order to characterize this parameter several radiosoundings were launched during July, November, and December 2016 and February 2017. The value of the parameter d as a function of time is obtained from a linear interpolation between the mean values of T_surf−T_trop measured during these four periods, as described in Table 2. T_trop is the average temperature obtained from radiosonde data by weighting the tropospheric temperature vertical profile with the water vapor concentration profile. Using Eq. (15), it is possible to explicitly give the relation between the opacity and the mean brightness temperature of the received signal:

$$\begin{array}{}\text{(17)}& {\displaystyle }\mathit{\mu }\left({\mathit{\theta }}_i\right)\mathit{\tau }=\ln \left({\displaystyle \frac{T_{\mathrm{0}}-T_{\text{trop}}}{\stackrel{\mathrm{‾}}{T}\left({\mathit{\theta }}_i\right)-T_{\text{trop}}}}\right).\end{array}$$

A linear regression of the opacities $\ln \left(\frac{T_{\mathrm{0}}-T_{\text{trop}}}{\stackrel{\mathrm{‾}}{T}\left({\mathit{\theta }}_i\right)-T_{\text{trop}}}\right)$ observed at θ_i and the air mass factors μ(θ_i) allows us to retrieve the opacity at the zenith, τ (Nedoluha et al., 1995). Substituting for $\stackrel{\mathrm{‾}}{T}\left({\mathit{\theta }}_i\right)$ using Eq. (8) and approximating all the spectral quantities with their mean values (indicated by a bar) over the central 11 000 channels, Eq. (17) can be written as

$$\begin{array}{}\text{(18)}& {\displaystyle }\mathit{\mu }\left({\mathit{\theta }}_i\right)\mathit{\tau }=\ln \left({\displaystyle \frac{T_{\mathrm{0}}-T_{\text{trop}}}{\frac{\left(\stackrel{\mathrm{‾}}{V}\left({\mathit{\theta }}_i\right)-{\stackrel{\mathrm{‾}}{V}}_{\mathrm{0}}\right)}{\stackrel{\mathrm{‾}}{\mathit{\alpha }}}-{\stackrel{\mathrm{‾}}{T}}_{\text{rec}}-T_{\text{trop}}}}\right),\end{array}$$

where it appears that the mean values ${\stackrel{\mathrm{‾}}{T}}_{\text{rec}}$ and $\stackrel{\mathrm{‾}}{\mathit{\alpha }}$ are needed to perform the calculation.

Table 2 Mean values and SDs of T_surf−T_trop obtained from the radiosoundings.

The emission from two calibration sources at different temperatures is used to calculate $\stackrel{\mathrm{‾}}{\mathit{\alpha }}$ and ${\stackrel{\mathrm{‾}}{T}}_{\text{rec}}$ according to Eqs. (11) and (12). However, since the sky emission temperature $T_{\text{cold}}^{\text{sky}}$ is not known, an iterative procedure is used to obtain both τ and $T_{\text{cold}}^{\text{sky}}$. An initial opacity value, τ₀, is adopted as a first guess to obtain $T_{\text{cold},\mathrm{0}}^{\text{sky}}$ using Eq. (21) with $\mathit{\theta }=\mathrm{60}{}^{\circ }$. $T_{\text{cold},\mathrm{0}}^{\text{sky}}$ is then used to obtain α₀ and $T_{{\text{rec}}_{\mathrm{0}}}$ by means of Eqs. (11) and (12); μ(θ_i)τ is calculated for different elevation angles using Eq. (18), and ultimately a linear fit allows us to calculate a new estimate for τ, τ₁. The iterative procedure goes on until the intercept value is minimized. Figure 3 shows the results of a tipping curve measurement carried out on 10 December 2016.

Figure 3 (a) The values of $\ln \left(\frac{T_{\mathrm{0}}-T_{\text{trop}}}{\stackrel{\mathrm{‾}}{T}\left({\mathit{\theta }}_i\right)-T_{\text{trop}}}\right)$ (indicated with τ_sig on the y axis of panel a) as function of the air mass factor μ(θ_i) (blue points) measured during a tipping curve on 10 December 2016, and the fit result (green line). (b) The residuals (measurements minus fit) as a function of the elevation angle.

Figure 4 Time series of the noise diodes mean emission temperature calculated by means of the tipping curve procedure (blue and cyan solid circles) compared with values obtained by means of LN₂ calibrations (orange and red stars). The green dots and the right y axis display the ratio of the emissions from two noise diodes.

The value of $T_{\text{cold}}^{\text{sky}}$ measured with this procedure is used in Eqs. (11) and (13) to estimate ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$. In order to avoid the use of data acquired during inhomogeneous sky conditions, all measurements producing linear fits (see Fig. 3a) with a root mean square larger than 0.4 are discarded. Figure 4 shows the time series of both noise diodes' mean emission temperatures (blue and cyan) and their ratio (green dots and right y axis) from July 2016 to July 2017. The noise diode in blue is the one used as calibration diode. In the same plot, ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$ values obtained using a LN₂ cooled CV-3 as the cold load are also depicted (orange and red stars). The mean relative difference between ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$ values calculated with LN₂ and with tipping curves carried out immediately before or after LN₂ calibrations is (0.4±0.4) and (0.2±0.3) % for the calibration and the backup diodes, respectively. The ratio between the two noise diodes can provide insights into potential drifts of one diode with respect to the other. Since the estimated uncertainty on such a ratio is 0.05 (3.4 %), there appears to be no drift in the time frame discussed in this work.

4 Retrieval process

VESPA-22 water vapor vertical profiles are obtained using the optimal estimation theory (Rodger, 2000). In what follows the retrieved water vapor profile is indicated as x, whereas x_a indicates the a priori profile and $\stackrel{\mathrm{̃}}{\mathit{x}}$ the real water vapor atmospheric profile. The quantities y and y_a represent the measured and a priori spectra, respectively. The matrix K is the weighting functions matrix, whereas S_e and S_a are the covariance matrices associated to the spectral measurements and the a priori vertical profile, respectively. The profile x can be used to calculate the synthetic spectrum y_fit according to

$$\begin{array}{}\text{(19)}& {\displaystyle }{\mathit{y}}_{\text{fit}}={\mathit{y}}_{\mathrm{a}}+\mathbf{K}\left(\mathit{x}-{\mathit{x}}_{\mathrm{a}}\right).\end{array}$$

The altitude grid used for VESPA-22 retrievals starts from 10 km and goes up to 110 km altitude, at steps of 1 km. This range is much larger than the sensitivity interval of the instrument, which is in fact limited by the Doppler broadening at high altitudes and by the FFTS bandwidth and the tropospheric influence on the lower stratosphere at low altitudes. Only the central 400 MHz of the measured spectrum are used in the retrieval. The S_a matrix is computed according to

$$\begin{array}{}\text{(20)}& {\displaystyle }{S}_{\mathrm{a},ij}={\mathit{\sigma }}_i{\mathit{\sigma }}_j{e}^{\frac{-\left|z_i-z_j\right|}{h}},\end{array}$$

where σ_i is the root mean square of the variance of the a priori profile (expressed in volume mixing ratio, or vmr; see Fig. 5) at the altitude z_i, while h is a correlation altitude set to be 5 km. Values of σ_i were empirically chosen in order to optimize the characteristics of the retrieval (i.e., maximize sensitivity range and vertical resolution without introducing unphysical oscillations in the retrieved vertical profile).

Figure 5 The value of σ_i used to compute the a priori covariance matrix as a function of altitude (Eq. 20).

In the inversion process for VESPA-22 spectra, S_e is a diagonal matrix with its diagonal elements all equal and calculated using a two-step process. A first retrieval of the original, not smoothed (see Sect. 3.2), spectral measurement averaged over 24 h is performed using a fixed value ($\mathrm{1}\times {\mathrm{10}}^{-\mathrm{5}}$) for the S_e diagonal elements and the obtained profile, x₀, is used to calculate a synthetic spectrum y_0,fit by means of Eq. (19). In order to consider the spectral measurement noise in the retrieval process, a second and final inversion is then performed, this time with the S_e diagonal elements set to the ${\left({\mathit{y}}_{\text{unsmoothed}}-{\mathit{y}}_{\mathrm{0},\text{fit}}\right)}^{\mathrm{2}}$ mean value and the spectrum smoothed as indicated in Sect. 3.2. Note that also in this second inversion the diagonal elements of S_e are held constant. For a 24 h integration time the S_e diagonal values range from $\mathrm{3}\times {\mathrm{10}}^{-\mathrm{4}}$ K² (maximum value during summer) to $\mathrm{8}\times {\mathrm{10}}^{-\mathrm{6}}$ K² (minimum value obtained during winter).

A second-order polynomial (light blue curve in Fig. 7a) is also added to the retrieval in order to take into account the spectral emission from the upper troposphere and lower stratosphere, as well as a potential contribution to the baseline from the delrin sheet. The polynomial is calculated independently for each retrieved profile. The addition of this extra degree of freedom to the retrieval process reduces the altitude interval in which VESPA-22 retrievals can be considered reliable, raising the lower limit of the sensitivity range (see definition below and Fig. 8b) by approximately 6 km altitude. The use of a second-order term also introduces an additional source of uncertainty into the retrieved mixing ratio profile. Such a contribution is taken into account in the error analysis discussed in Sect. 4.3 (see Fig. 11).

An important quantity used to characterize the retrieval quality is the averaging kernels matrix A (Rodgers, 2000). The rows of A are called averaging kernels (AKs) and can be used to characterize the sensitivity of the water vapor retrieval at a given altitude to variations in the water vapor concentration profile at all altitudes (Rodger, 2000). If the AKs are well-peaked functions, centered at their nominal altitude, a perturbation in the atmospheric water vapor concentration at a specific altitude is transferred by the algorithm to the correct altitude layer of the retrieved profile. Furthermore, the area enclosed under each AK is an indication of the total sensitivity of the retrieved value at that altitude to atmospheric variations in water vapor concentration. A sensitivity value close to 1 at a certain altitude indicates that the major contribution to the retrieved value at that altitude comes from the spectral measurements rather than from the a priori water vapor profile. Following what was indicated by Tschanz et al. (2013), VESPA-22 retrievals are also considered valid for scientific use in the altitude range where the sensitivity is above 0.8. The AKs' FWHM can be used as a rough estimate of the local vertical resolution of the obtained water vapor mixing ratio vertical profile.

Figure 6 The a priori profiles employed by VESPA-22 for the October–April period (indicated as winter a priori, blue line) and for the June–August period (indicated as summer a priori, red line) obtained from climatology.

Table 3 Spectroscopic parameters (“reference” model) used for VESPA-22 retrievals. Indicated parameters, from left to right, are emission line frequency, intensity, lower state energy, pressure broadening coefficient, pressure broadening temperature dependence, self-broadening coefficient, and self-broadening temperature dependence. The line intensity is given for a reference temperature of 296 K.

4.1 Forward model and a priori profiles

Figure 6 displays the a priori profiles used by the VESPA-22 retrieval algorithm during 10 out of the 12 months of the year and based on local climatology (3 years' worth of data of Aura/MLS v4.2). These profiles are identical below 48 km and diverge above, due to the large difference in polar regions between summer and winter water vapor mesospheric profiles. The summer a priori (red line) is used during the period from 1 June to 15 September, while the winter a priori (blue line) is used during the period from 16 October to 30 April. During the month of May and from 16 September to 15 October, there are transition periods in which a linear daily interpolation from one a priori to the other was used for continuity.

The matrix K and the y_a spectrum are calculated using the radiative transfer simulation software ARTS (Eriksson et al., 2011), adopting a Voigt–Kuntz lineshape and the line intensity provided by the JPL 2012 catalogue (Pickett et al., 1998, and https://spec.jpl.nasa.gov/). Following the work of Seele (1999) and Tschanz et al. (2013), the line described by the JPL 2012 catalogue is divided into three emission lines indicating the hyperfine splitting of the 22.235 GHz water vapor line. The employed pressure broadening and self-broadening parameters are those reported by Liebe (1989). Table 3 summarizes the spectroscopic parameters used for the analysis of VESPA-22 spectral measurements which, in what follows, are indicated as the “reference” model.

The profile x_ARTS is the profile used in the forward model calculations to compute the a priori spectrum and the weighting functions. x_ARTS matches the a priori profile from 12 km of altitude upward and, below 9 km, it is consistent with the measurements of precipitable water vapor (PWV) collected by the HATPRO radiometer (Rose and Czekala, 2009; Pace et al., 2015) installed at the THAAO. This lower part of x_ARTS is calculated according to

$$\begin{array}{}\text{(21)}& {\displaystyle }{\mathit{x}}_{\text{ARTS}}\left(\text{from the ground to 9\hspace{0.17em}km}\right)={\displaystyle \frac{{\text{PWV}}_{\text{Hatpro}}}{{\text{PWV}}_{\text{Eu}}}}{\mathit{x}}_{\text{Eu}},\end{array}$$

where x_Eu is a water vapor mixing ratio profile obtained by daily smoothing monthly averages calculated from the radiosoundings launched from the Eureka station, PWV_Eu is the associated water vapor column content, and PWV_Hatpro is the column content measured by the HATPRO located at the THAAO. Data from Eureka were chosen (instead of those from Alert, Canada, for example) because they show the closest resemblance to the tropospheric profiles measured at Thule by local radiosoundings, when the latter are available. In order to avoid discontinuities in x_ARTS, values at altitudes between 9 and 12 km are obtained with a linear interpolation between x_ARTS(9 km) and x_ARTS(12 km).

The pressure and temperature profiles needed to run the forward calculation are built merging NASA Goddard Space Flight Center (GSFC), Aura/MLS and climatological temperature and pressure profiles. The NASA GSFC profiles obtained through the Goddard Automailer Service (Lait et al., 2005) are used to build the tropospheric meteorological state, from the ground up to 9 km of altitude. Between 10 and 87 km altitude, the MLS temperature and pressure profiles collected during VESPA-22 observations, in a radius of 300 km from the observation point of VESPA-22, are averaged together to produce a single set of daily meteorological vertical profiles. The VESPA-22 observation point coordinates are chosen to be 74.8^∘ N and 73.5^∘ W, and represent an estimate of the geographical coordinates of the air mass that is observed by VESPA-22 (which points southwest, at about 220^∘) at 60 km altitude when the instrument aims at an elevation of 15^∘ above the horizon. Daily temperature and pressure profiles from 97 to 110 km of altitude are obtained by daily smoothing zonal monthly averages from the COSPAR International Reference Atmosphere (Fleming et al., 1990).

The absence of vertical discontinuities in the temperature and pressure daily profiles is assured by a smoothing process performed at the altitudes where the three different datasets (GSFC, MLS, and COSPAR) are stitched together, between 9 and 12 km and between 87 and 97 km.

The software ARTS is employed to simulate the emission from the zenith, ${\stackrel{\mathrm{̃}}{\mathit{y}}}_{\mathrm{r}}$, and from an angle close to the horizon, y_s. The emission ${\stackrel{\mathrm{̃}}{\mathit{y}}}_{\mathrm{r}}$ is then rescaled using Eq. (22), therefore simulating the effect of the delrin compensating sheet:

$$\begin{array}{}\text{(22)}& {\displaystyle }{\mathit{y}}_{\mathrm{r}}={\stackrel{\mathrm{̃}}{\mathit{y}}}_{\mathrm{r}}+T_{\mathrm{d}}\left(\mathrm{1}-e^{-{\mathit{\tau }}_{\mathrm{d}}}\right),\end{array}$$

whereas the mean difference y_s−y_r is

$$\begin{array}{}\text{(23)}& {\displaystyle }{\mathit{y}}_{\mathrm{s}}-{\mathit{y}}_{\mathrm{r}}={\mathit{y}}_{\mathrm{a}}\left(\stackrel{\mathrm{̃}}{\mathit{\mu }}{e}^{-\stackrel{\mathrm{̃}}{\mathit{\mu }}\stackrel{\mathrm{̃}}{\mathit{\tau }}}\right)-e^{-\stackrel{\mathrm{̃}}{\mathit{\tau }}-{\mathit{\tau }}_{\mathrm{d}}}.\end{array}$$

The opacity τ_d used in Eq. (22) is the opacity of the compensating sheet; the temperature T_d is the temperature of the sheet, measured by a sensor installed next to it; and the signal beam observing angle, $\stackrel{\mathrm{̃}}{\mathit{\theta }}$, is chosen in order to minimize the mean difference y_s−y_r, as it is in fact attained by VESPA-22 in its data-taking process. The a priori spectrum is calculated according to the same equation used for the measured signal (Eqs. 5 and 6),

$$\begin{array}{}\text{(24)}& {\displaystyle }{\mathit{y}}_{\mathrm{a}}={\displaystyle \frac{{\mathit{y}}_{\mathrm{s}}-{\mathit{y}}_{\mathrm{r}}}{\stackrel{\mathrm{̃}}{\mathit{\mu }}{e}^{-\stackrel{\mathrm{̃}}{\mathit{\mu }}\stackrel{\mathrm{̃}}{\mathit{\tau }}}-e^{-\stackrel{\mathrm{̃}}{\mathit{\tau }}-{\mathit{\tau }}_{\mathrm{d}}}}},\end{array}$$

where $\stackrel{\mathrm{̃}}{\mathit{\mu }}$ and $\stackrel{\mathrm{̃}}{\mathit{\tau }}$ are the air mass factor associated to the simulated signal beam and the zenith opacity calculated from the x_ARTS profile.

Figure 7 (a) An example of VESPA-22 measured spectrum (blue) collected on 23 December 2016, with the a priori spectrum in green and the fit spectrum in red. The second-order polynomial is indicated with a cyan solid line. (b) The residual y_fit−y. The central part of the spectrum is unsmoothed in order to maintain the maximum spectral resolution near the peak and its residual is larger.

Deriving Eqs. (24) and (22), and using the K matrix definition (Rodgers, 2000), the retrieval weighting function matrix can be written as

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }\mathbf{K}={\displaystyle \frac{\partial }{\partial \mathit{x}}}{\left({\displaystyle \frac{{\mathit{y}}_{\mathrm{s}}-{\mathit{y}}_{\mathrm{r}}}{\stackrel{\mathrm{̃}}{\mathit{\mu }}{e}^{-\stackrel{\mathrm{̃}}{\mathit{\mu }}\stackrel{\mathrm{̃}}{\mathit{\tau }}}-e^{-\stackrel{\mathrm{̃}}{\mathit{\tau }}-{\mathit{\tau }}_{\mathrm{d}}}}}\right)}_{{\mathit{x}}_{\mathrm{a}}}\\ \text{(25)}& {\displaystyle }& {\displaystyle }\cong {\displaystyle \frac{{\left[\frac{\partial {\mathit{y}}_{\mathrm{s}}}{\partial \mathit{x}}\right]}_{{\mathit{x}}_{\mathrm{a}}}-{\left[\frac{\partial {\mathit{y}}_{\mathrm{r}}}{\partial \mathit{x}}\right]}_{{\mathit{x}}_{\mathrm{a}}}}{\stackrel{\mathrm{̃}}{\mathit{\mu }}{e}^{-\stackrel{\mathrm{̃}}{\mathit{\mu }}\stackrel{\mathrm{̃}}{\mathit{\tau }}}-e^{-\stackrel{\mathrm{̃}}{\mathit{\tau }}-{\mathit{\tau }}_{\mathrm{d}}}}}={\displaystyle \frac{{\mathbf{K}}_{\mathrm{s}}-{\mathbf{K}}_{\mathrm{r}}{e}^{-{\mathit{\tau }}_{\mathrm{d}}}}{\stackrel{\mathrm{̃}}{\mathit{\mu }}{e}^{-\stackrel{\mathrm{̃}}{\mathit{\mu }}\stackrel{\mathrm{̃}}{\mathit{\tau }}}-e^{-\stackrel{\mathrm{̃}}{\mathit{\tau }}-{\mathit{\tau }}_{\mathrm{d}}}}},\end{array}$$

where K_s and K_r are the weighting function matrices that ARTS calculates for the simulated signal and reference beams. As a first approximation, the dependence of $\stackrel{\mathrm{̃}}{\mathit{\tau }}$ on the stratospheric water vapor profile in Eq. (25) is neglected.

4.2 Retrieval example

Figure 7a shows a VESPA-22 spectrum integrated for 24 h (blue line) on 23 December 2016, its corresponding synthetic spectrum y_fit (red line) and the a priori spectrum (green line), while the residual (defined as the difference between fit and measured spectrum, y_fit−y) is plotted in Fig. 7b. The cyan line is the second-degree polynomial retrieved by the inversion algorithm. Figure 8a shows the result of the inversion of the measured spectrum depicted in Fig. 7 with the a priori profile and retrieval 1σ uncertainty. The details on the uncertainty calculation are discussed in Sect. 4.3.

4.3 Retrieval uncertainty

The uncertainty characterizing VESPA-22 retrieved profiles can be divided into four major contributions: (1) the uncertainty due to the linear approximation used in the optimal estimation, (2) the uncertainty due to the various parameters used in the spectra calibration and pre-processing, (3) the uncertainty due to spectral noise and potential artifacts, and (4) the uncertainty introduced by the use of the second-order polynomial in the retrieval process. One additional error source is the limited vertical resolution inherent to concentration vertical profiles obtained by means of this ground-based observing technique. This leads to solution profiles that can be considered a smoothed version of the real atmospheric concentration profiles. In discussing the optimal estimation method, Rodgers (2000) suggests that this error, called “smoothing error”, should be estimated only if accurate knowledge of the variability of the atmospheric fine structure is available. This approach is used here and the smoothing error is not included in the error estimate.

The first contribution can be evaluated observing the difference Δy_lin between the fit spectrum without the addition of the second-order polynomial, ${\mathit{y}}_{\text{fit}}^{\ast }$, and the spectrum obtained using ARTS to calculate the emission expected from the retrieved profile x (in the calculation, ARTS does not perform any linear approximation):

$$\begin{array}{}\text{(26)}& {\displaystyle }\mathrm{\Delta }{y}_{\text{lin}}={\mathit{y}}_{\text{fit}}^{\ast }-f\left(x\right).\end{array}$$

The uncertainty Δx_lin that Δy_lin causes on the retrieved profile can be calculated with

$$\begin{array}{}\text{(27)}& {\displaystyle }\mathrm{\Delta }{x}_{\text{lin}}=\mathbf{G}\mathrm{\Delta }{y}_{\text{lin}},\end{array}$$

where G is the gain matrix (Rodgers, 2000). The uncertainty Δx_lin has a negligible contribution to the total uncertainty (with a maximum of 0.1 % at 70 km altitude).

Figure 8 (a) The retrieved VESPA-22 profile (blue solid line) correspondent to the spectrum showed in Fig. 7. The a priori profile is indicated with a green solid line. The two red dashed lines describe the uncertainty of this VESPA-22 retrieval (for details on the estimated uncertainty of VESPA-22 mixing ratio vertical profiles see Sect. 4.3). (b) Rows of the A matrix multiplied by a factor of 10 as a function of altitude (some A functions are highlighted in colors). The vertical profile of the sensitivity is shown in red.

Table 4 Uncertainties on the various parameters used in the calibration process. When the uncertainty is a function of altitude the minimum and maximum values of the uncertainty are reported.

Figure 8b shows the corresponding averaging kernels (AKs, black and colored solid lines). AKs are multiplied by a factor of 10 and the sensitivity in indicated with a red solid line.

In order to evaluate the second contribution listed above, the effects on the retrieved profile due to the variation of each single parameter used in the measurements calibration and pre-processing was investigated. The difference between the profile retrieved using the “correct” value of a specific parameter and the retrieval obtained by changing such a value by the estimated relative uncertainty of the parameter is considered to be the contribution σ_i of this parameter to the total calibration and pre-processing uncertainty. The total uncertainty from these sources is defined as calibration uncertainty. Table 4 summarizes the uncertainties on the various parameters involved in the calibration and pre-processing of VESPA-22 spectra. When the uncertainty is a function of altitude the minimum and maximum values of the uncertainty are reported. The total calibration uncertainty σ_cal is given by

$$\begin{array}{}\text{(28)}& {\displaystyle }{\mathit{\sigma }}_{\text{cal}}=\sqrt{\sum _i{\mathit{\sigma }}_i^{\mathrm{2}}}.\end{array}$$

Figure 9 Relative contributions to the calibration uncertainty defined by Eq. (28) and indicated with a red solid curve. The contributions are the signal beam angle (yellow), the noise diode temperature (green), the sky opacity (blue), the MLS meteorological profiles (orange and brown), the compensating sheet opacity (magenta), and the spectroscopic parameters (cyan).

Figure 10 (a) Mean retrieved water vapor profiles obtained using different spectroscopic models (the blue line is the model eventually adopted for the analysis of VESPA-22 spectra, called reference model, or RM). The spectroscopic parameters for the additional models tested are obtained from different versions of the JPL catalogue (years 1985, 2001 and 2012; Pickett et al., 1998) with the pressure broadening parameters taken from Liebe (1989) and Cazzoli et al. (2007). Two models include values from HITRAN 2004 and 2012 (Rothman et al., 2013). Models are indicated in the legend with their abbreviations: JPL2012Liebe1989 (JPL12L89), JPL1985Liebe, 1989 (JPL85L89), JPL1985Cazzoli, 2007 (JPL85C07), JPL2001Liebe, 1989 (JPL01L89), HITRAN2012 (H12), and HITRAN2004 (H04). (b) Vertical profiles of mean relative differences between VESPA-22 profiles obtained using the reference model and the different spectroscopic models. Data represented here range from 4 October 2016 to 22 May 2017. The average retrieval sensitivity range is marked by dashed horizontal green lines.

In Fig. 9, the relative contributions of the various parameters needed in the calibration and pre-processing procedures are shown. The yellow line shows the σ_i contribution due to the uncertainty on the signal beam angle. In order to minimize this contribution, the choked antenna and the parabolic mirror are aligned using a He-Ne laser. During the period from April to October, the VESPA-22 pointing offset can be verified by scanning rapidly at around the sun elevation angle and comparing the measured position of the center of the sun against the known ephemerides. For the measurements discussed here an offset of $\mathrm{\Delta }\mathit{\theta }=\mathrm{0.2}{}^{\circ }\pm \mathrm{0.1}{}^{\circ }$ was estimated. The green solid line shows the potential relative error on the water vapor mixing ratio vertical profile due to the uncertainty on the estimated noise diode temperature ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$. In order to evaluate the uncertainty on ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$, the 0.4 % difference between the values obtained by means of the LN₂ calibrations and by means of the skydips immediately following or preceding the LN₂ calibrations (see Sect. 3.4 and Fig. 4) is considered. On top of this, the fluctuations of the signal produced by the calibration diodes must also be taken into account. These can be evaluated by using the standard deviation (SD) of the difference over time between the ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$ values of the two noise diodes, measured to be 1.1 %. The tipping curve procedure allows for calculation of the noise diodes brightness temperature averaged over the 11 000 central channels of the spectrometer. Since the noise diodes brightness temperatures do have a small frequency dependence, using their mean values introduces a source of uncertainty which is estimated to be 0.9 %. An additional source of uncertainty is due to sky inhomogeneities, and amounts approximately to 1 %. All these sources of uncertainty for ${\stackrel{\mathrm{‾}}{T}}_{\mathrm{nd}}$ are added in quadrature to obtain a total uncertainty of 1.8 %. The blue line shows the contribution due to the uncertainty Δτ on the sky zenith opacity τ. In order to estimate this uncertainty, the uncertainties introduced by sky inhomogeneity and by the estimation of the effective tropospheric temperature T_trop are considered. The first contributes to about 2 % and is estimated using the uncertainty on the slope of the linear fit of the tipping curve (see Sect. 3.4 and Fig. 3). The second contribution can be evaluated observing the daily fluctuations of T_trop measured at Eureka (Canada), Aasiaat (southern Greenland), and Alert (Canada). The natural day-to-day temperature fluctuations and the lack of tropospheric meteorological data at Thule during long intervals of time led us to estimate an uncertainty on T_trop of about 6 K, which then makes the total uncertainty Δτ∼5 %.

The brown and orange lines show the σ_i's due to temperature and geopotential height profiles used in the forward calculation. The uncertainties on these parameters are obtained from the MLS data quality and description document (Livesey et al., 2015). In Fig. 9, the magenta line shows the contribution of the uncertainty on the compensating sheet opacity, Δτ_d. Although the uncertainty on the measurements of τ_d is about 2 % (see Table 1), the use of a linear interpolation suggests the use of a more conservative estimate, eventually set at 10 %. The cyan line shows the contribution due to uncertainties in the employed spectroscopic parameters, and it has the largest impact on the calibration uncertainty. Following Straub et al. (2010), the emission line intensity and the pressure broadening coefficient were assigned an uncertainty of $\mathrm{8.7}\times {\mathrm{10}}^{-\mathrm{22}}$ m² Hz and 1014 Hz Pa⁻¹, respectively.

In addition to these uncertainties, Fig. 10 shows the results of retrieving VESPA-22 data from October 2016 to May 2017 using different sets of values for the spectroscopic parameters (see also Haefele et al., 2009). In the figure, the mean profiles retrieved using the different models listed in the legend are depicted. Such models include values taken from HITRAN (Rothman et al., 2013) and JPL catalogues (Pickett et al., 1998), for different years, and from Liebe (1989) and Cazzoli et al. (2007). The blue line represents the mean profile obtained using VESPA-22 reference model (RM). In Fig. 10b, the mean relative difference profiles between retrievals obtained using the reference model and the other spectroscopic models are drawn (see legend). The water vapor mixing ratio retrieved profiles show a strong dependence on the spectroscopic model of choice, with a relative difference between profiles obtained using different models reaching a maximum of 8 % at the top of the sensitivity range.

Figure 11 shows the uncertainty calculated for the retrieval of 23 December 2016, shown in Fig. 8. The overall calibration uncertainty is shown with a red curve. The uncertainty due to spectral noise and artifacts can be evaluated using the S uncertainty matrix (Rodger, 2000) obtained as

$$\begin{array}{}\text{(29)}& {\displaystyle }\mathbf{S}={\left({\mathbf{K}}^T{\mathbf{S}}_{\mathrm{e}}^{-\mathrm{1}}\mathbf{K}+{\mathbf{S}}_{\mathrm{a}}^{-\mathrm{1}}\right)}^{-\mathrm{1}}.\end{array}$$

The square root of the diagonal elements of S represents the spectral uncertainty of the retrieved profile at different altitudes and is indicated in Fig. 11 with a blue line. This term depends on the S_e value of the retrieval and its value increases during summer and with poor weather conditions. As mentioned in Sect. 4, the diagonal values of the S_e matrix are held constant and do not take into account the spectral noise reduction obtained with the 50-channel smoothing process which is applied to most of the spectrum (Sect. 3.2). This may lead to a small overestimation of the spectral uncertainty depicted in Fig. 11 in the altitude range 20–50 km (with a maximum potential overestimation of 0.3 % at 25 km altitude). The green curve displays the estimated uncertainty due to the use of the second-order polynomial in the retrieval process. The total uncertainty is calculated as

$$\begin{array}{}\text{(30)}& {\displaystyle }{\mathit{\sigma }}_{\text{tot}}=\sqrt{{\mathit{\sigma }}_{\text{cal}}^{\mathrm{2}}+{\mathit{\sigma }}_{\text{spec}}^{\mathrm{2}}+{\mathit{\sigma }}_{\text{poly}}^{\mathrm{2}}}\end{array}$$

and is represented with a purple line in Fig. 11.

5 VESPA-22 and Aura/MLS water vapor profile intercomparison

VESPA-22 water vapor vertical retrievals obtained from 24 h integration spectra (from 00:00 to 23:59 UT) have been compared with version 4.2 of Aura/MLS water vapor vertical profiles. The MLS profiles used for this intercomparison are daily mean profiles obtained averaging all MLS profiles collected within a radius of 300 km centered around VESPA-22 observation point (74.8^∘ N, 73.5^∘ W; see Sect. 4.1).

Figure 11 Vertical profiles of the calibration uncertainty (red), the retrieval uncertainty (blue), and the total uncertainty (magenta) of VESPA-22 water vapor mixing ratio vertical profiles obtained inverting a spectrum collected on 23 December 2016 and integrated for 24 h.

The intercomparison is carried out using data from 15 July 2016 to 2 July 2017. In general, spectra collected during July, August, and September are less continuous and noisier than spectra observed during the rest of the selected period. This is due to extensive testing of the equipment, poor weather conditions, and snow covering the zenith observing window (in November 2016 a second powerful fan was installed outside the roof window in order to reduce snow deposition). Furthermore, there are no stratospheric measurements carried out by VESPA-22 between 4 and 22 November 2016, due to snow covering the lower portion of the reference beam window. A few isolated days in which large sky inhomogeneities do not allow the correct balance of signal and reference beams have also been removed from the intercomparison. It is worth recalling that the signal-to-noise ratio of VESPA-22 spectra, and therefore the quality of the retrievals, depends on the sky opacity and, consequently, on the season, being noticeably better during winter and poorer in summer. This is particularly true for microwave spectrometers operating in polar regions, where seasonal fluctuations of the tropospheric water vapor column content are significant.

Table 5 summarizes the characteristics of the Aura/MLS water vapor version 4.2 retrievals (Livesey et al., 2015). In order to compare the two datasets, MLS vertical profiles are convolved with VESPA-22 averaging kernels in order to match the resolution of VESPA-22 profiles according to Rodgers (2000):

$$\begin{array}{}\text{(31)}& {\displaystyle }{\mathit{x}}_{\text{MLS}}={\mathit{x}}_{\mathrm{a}}+\mathbf{A}({\stackrel{\mathrm{̃}}{\mathit{x}}}_{\text{MLS}}-{\mathit{x}}_{\mathrm{a}}),\end{array}$$

where ${\stackrel{\mathrm{̃}}{\mathit{x}}}_{\text{MLS}}$ is the raw (high resolution or HR) MLS water vapor vertical profile, x_a is the a priori profile, A is the averaging kernel matrix, and x_MLS is the convolved MLS profile. Additionally, VESPA-22 profiles were also compared with MLS profiles smoothed in the vertical by using a 10 km moving average. This second set of degraded MLS profiles was generated in order to study the correlation between VESPA-22 and MLS datasets without introducing the dependency from one another brought by the convolution process (affecting MLS convolved profiles).

Figure 12 The MLS profile measured on 23 December 2016, with its high original vertical resolution (HR, magenta line) compared with the a priori profile (green line), the MLS convolved profile obtained from Eq. (31) (red line), and the VESPA-22 retrieved profile (blue line).

Figure 13 (a) VESPA-22 averaged water vapor mixing ratio vertical profile with its SD (blue line and blue dashed lines) and the mean MLS convolved profile with its SD (red line and red dashed lines); (b) vertical profiles of the mean relative difference between VESPA-22 and MLS convolved (red line), and between VESPA-22 and MLS smoothed (blue line) with their SDs (red and blue dashed lines, respectively); (c) vertical profiles of the correlation coefficient between VESPA-22 and MLS convolved profiles (red line), and between VESPA-22 and MLS smoothed profiles (blue line); (d) vertical profile of the mean FWHM of VESPA-22 averaging kernels with its SD. The data used for the intercomparison range from 12 July 2016 to 2 July 2017. The green dashed horizontal lines display the mean altitudes of the sensitivity range extremes.

Table 5 The main specifications of Aura/MLS water vapor mixing ratio profiles version 4.2 (Livesey et al., 2015).

Figure 12 shows the MLS measured profile on 23 December 2016 (magenta line), the MLS convolved profile obtained from Eq. (31) (red line), the a priori profile (green line), and the VESPA-22 retrieved profile (blue line) for the same day. The MLS convolved profile tends to the a priori profile below 26 km and above 72 km, where the retrieval sensitivity drops.

Figure 13a shows the mean VESPA-22 retrieved profile (in blue) and the mean MLS convolved profile (in red), with their SDs indicated with dashed lines. The mean sensitivity of VESPA-22 retrieved profiles is larger than 0.8 from about 26 to 72 km altitude and is indicated by the two green horizontal dashed lines. This interval can vary from day to day depending on the noise level of the 24 h integrated spectra (see white solid lines in Fig. 15). Figure 13b displays the relative difference of VESPA-22 water vapor mixing ratio mean vertical profile with respect to the MLS mean convolved profile (red line) and with respect to the MLS smoothed mean vertical profile (blue line) with their SDs (blue and red dashed lines, respectively). The largest relative and absolute (not shown) differences between the two datasets occur at 72 km, the upper limit of the VESPA-22 sensitivity range, and it is about −6 %, with Aura/MLS convolved mean profile being larger than VESPA-22 mean retrieval. The SD of the differences increases in the mesosphere as a result of the larger atmospheric variability and the larger relative uncertainty of both instruments at these altitudes (see Fig. 11 and Table 5) with respect to lower levels.

Figure 14 Time series at different altitudes of water vapor mixing ratio values obtained with VESPA-22 (blue), by convolving MLS HR vertical profiles with VESPA-22 averaging kernels (red), and by smoothing MLS HR vertical profiles with a 10 km running average (yellow). Please note that the time series at the five different altitude levels have different scales on the y axis in order to better show the differences between datasets.

Figure 13c shows the vertical profiles of the correlation coefficient between VESPA-22 and MLS convolved data (red), and between VESPA-22 and MLS smoothed data (blue). The latter shows values of 0.8 or higher over the entire sensitivity range. The good correlation of VESPA-22 with MLS smoothed profiles down to 26 km altitude suggests that VESPA-22 retrievals are reliable at these stratospheric altitudes. Such a high correlation was obtained by using for the VESPA-22 retrieval algorithm an a priori profile held constant in time up to 48 km, therefore providing no contribution to the correlation. Below 25 km the correlation with MLS smoothed quickly deteriorates, whereas the correlation of VESPA-22 with MLS convolved profiles remains high, emphasizing the dependency of both datasets on VESPA-22 a priori profile and averaging kernels. Figure 13d shows the vertical profile of the mean values of the FWHM of VESPA-22 averaging kernels (blue solid line) calculated over the comparison period. The FWHM of the averaging kernels is a measure of the retrieval vertical resolution (Rodgers, 2000). The single-profile vertical resolution can vary depending on the level of noise affecting the measured spectrum.

Figure 14 shows the time series at different altitudes of water vapor mixing ratio values from VESPA-22 (blue dots), MLS convolved (red dots), and MLS smoothed (yellow dots) datasets. The 25 and 35 km water vapor time series display rapid and intense variations in winter caused by the polar vortex moving over Thule in mid-December and then away from Thule at the end of January. The 25 km time series is also useful to evaluate the quality of the VESPA-22 measurements at the bottom limit of the sensitivity range. Figure 14 shows that, at this altitude, VESPA-22 is capable of depicting the rapid variations in water vapor but does not have the necessary sensitivity to match the most intense variations observed by MLS (mid-January and early February). Both MLS and VESPA-22 time series at 65 km show a peak in mid-April, although VESPA-22 underestimates the peak intensity with respect to MLS possibly due to the reduced sensitivity of VESPA-22 retrievals at this altitude during spring.

In order to provide a more complete, albeit less quantitative, overview of the VESPA-22 and convolved MLS time series, Fig. 15 shows contour maps of VESPA-22 (in color) and MLS convolved (black lines) water vapor profiles. The data collected by both instruments reveal an absolute maximum in August 2016 at 50 km of height of about 8.3 ppmv. During fall and winter the maximum of the water vapor mixing ratio profile lowers its altitude to reach about 35 km in January 2017, due also to the downward motion of air inside the polar vortex. The two steep gradients in water vapor mixing ratio due to the polar vortex moving away from Thule are clearly visible between 25 and 40 km altitude.

Figure 15 A map showing the VESPA-22 retrieved profiles (colored areas) compared to MLS convolved profiles (black lines). The blank areas indicate the lack of VESPA-22 data for more than 3 days. White solid lines indicate the time series of the sensitivity interval.

Both instruments also observe the return of the water vapor mixing ratio to pre-winter values in mid-April, possibly indicating the occurrence of the vortex final warming, with the re-establishment of a maximum in the mixing ratio profile at about 50 km altitude.

6 Summary

VESPA-22 was installed at the Thule High Arctic Atmospheric Observatory (THAAO; http://www.thuleatmos-it.it/) located at Thule Air Base, Greenland, in July 2016 for long-term observations of the polar middle atmospheric water vapor and column content. The instrument is characterized by a FWHM encompassing 3.5^∘ (Bertagnolio et al., 2012), granting the observation of the signal beam at angles as low as about 12^∘ above the horizon. The instrument is installed indoor in a wooden annex to the main laboratory and observes the sky emission through windows made of 5 cm thick Plastazote LD15 sheets. During clear-sky conditions, there are no evident artifacts larger than 2 mK affecting the measured spectra. VESPA-22 operated automatically with minimum need for maintenance for about 1 year, proving the robustness of hardware and acquiring system.

The instrument calibration is regularly performed using two noise diodes. The emission temperature of these elements is measured during tipping curve calibrations and checked periodically with liquid nitrogen calibrations. The calibrating noise diode temperature is estimated with an uncertainty of 1.8 %.

VESPA-22 retrieval algorithm is based on the optimal estimation technique (Rodgers, 2000); the retrieved profiles from 24 h integration spectra have an average sensitivity larger than 0.8 from about 26 to 72 km of altitude and a vertical resolution from about 12 to 23 km. The forward model is provided by the ARTS software (Eriksson et al., 2011), which is tuned to reproduce the VESPA-22 measurement technique. In order to perform a realistic simulation of the troposphere, the PWV measured by the HATPRO radiometer (Rose and Czekala, 2009; Pace et al., 2015) operating side by side with VESPA-22 is taken into account in the forward model.

The uncertainty on VESPA-22 retrieved water vapor mixing ratio vertical profiles is evaluated as the sum of the contributions from calibration, pre-processing, spectroscopic parameters, measurement noise, and use of a second-order polynomial baseline in the retrieval process. In the sensitivity range of VESPA-22 retrievals, the total uncertainty is estimated to be about 5–6 % from 26 to 60 km, increasing to about 18 % at 72 km.

VESPA-22 and MLS water vapor profiles are compared during a period from July 2016 to July 2017. The VESPA-22 data used in the comparison are the results of retrievals from 24 h integration spectra, while MLS data are the daily mean vertical profiles collected by MLS in a radius of 300 km around VESPA-22 observation point, convolved with VESPA-22 averaging kernels. No significant biases were observed between the two datasets in the altitude range from 25 to 60 km, whereas from 60 to 72 km the value of the mean difference (VESPA-22 – MLS) increases, reaching −6 % (−0.2 ppmv) at 72 km. A good correlation is also found between VESPA-22 and MLS smoothed (see Sect. 5) vertical profiles in the entire vertical range (from 26 to 72 km altitude, on average) where VESPA-22 profiles are considered valid.

The results described in this paper proved that VESPA-22 is capable of carrying out reliable middle atmospheric water vapor measurements during different seasons and weather conditions, although the vertical range where the dataset is recommended for scientific use is reduced in summer with respect to winter. VESPA-22 retrievals are capable of correctly representing the rapid variations in water vapor concentrations that can occur in the stratosphere, as revealed by the large water vapor gradients measured both by VESPA-22 and Aura/MLS in mid-December 2016, late January/early February 2017, and then again during April and May (see Fig. 15). At the same time, VESPA-22 spectral measurements and corresponding mixing ratio retrievals appear to have the necessary stability and consistency to observe slow seasonal variations, such as the water vapor subsidence occurring inside and at the edge of the polar vortex.