Articles | Volume 11, issue 12
Research article
06 Dec 2018
Research article |  | 06 Dec 2018

Boundary-layer water vapor profiling using differential absorption radar

Richard J. Roy, Matthew Lebsock, Luis Millán, Robert Dengler, Raquel Rodriguez Monje, Jose V. Siles, and Ken B. Cooper

Remote sensing of water vapor in the presence of clouds and precipitation constitutes an important observational gap in the global observing system. We present ground-based measurements using a new radar instrument operating near the 183 GHz H2O line for profiling water vapor inside of planetary-boundary-layer clouds, and develop an error model and inversion algorithm for the profile retrieval. The measurement technique exploits the strong frequency dependence of the radar beam attenuation, or differential absorption, on the low-frequency flank of the water line in conjunction with the radar's ranging capability to acquire range-resolved humidity information. By comparing the measured differential absorption coefficient with a millimeter-wave propagation model, we retrieve humidity profiles with 200 m resolution and typical statistical uncertainty of 0.6 g m−3 out to around 2 km. This value for humidity uncertainty corresponds to measurements in the high-SNR (signal-to-noise ratio) limit, and is specific to the frequency band used. The measured spectral variation of the differential absorption coefficient shows good agreement with the model, supporting both the measurement method assumptions and the measurement error model. By performing the retrieval analysis on statistically independent data sets corresponding to the same observed scene, we demonstrate the reproducibility of the measurement. An important trade-off inherent to the measurement method between retrieved humidity precision and profile resolution is discussed.

1 Introduction

In this work, we discuss the implementation of differential absorption radar (DAR) for measuring humidity profiles inside of boundary-layer clouds (Lebsock et al.2015; Millán et al.2016). The DAR method, which is the microwave analog of the mature differential absorption lidar (DIAL) method (Browell et al.1979), combines the range-resolving capabilities of radar with the strong frequency dependence of atmospheric attenuation near a molecular rotational absorption line to retrieve density profiles of the absorbing gas along the line of sight. Recently, there was a demonstration of microwave integrated path differential absorption in airborne measurements of sea surface air pressure without range resolution (Lawrence et al.2011), utilizing the 60 GHz O2 line to measure the total oxygen column. More recently, our group demonstrated a ground-based DAR for humidity sounding operating between 183 and 193 GHz (Cooper et al.2018), with primary sensitivity to upper tropospheric water vapor due to significant attenuation in the lower troposphere at these frequencies. That work included a comparison of differential absorption measurements with a millimeter-wave propagation model showing good agreement, and left the topics of error analysis and profile inversion for future investigation. While the 183 to 193 GHz band is attractive for DAR measurements because of the large differential absorption values achievable, transmission at frequencies between 174.8 and 191.8 GHz is prohibited due to reservation for passive-only remote sensing (NTIA2015). On the other hand, the 167 to 174.8 GHz band offers fewer transmission restrictions, and features lower absolute absorption, thus enabling penetration into the boundary layer from an airborne or spaceborne platform. Of course, the smaller absolute absorption is accompanied by decreased differential absorption, making the profiling capabilities of this radar coarser than the 183 to 193 GHz DAR. Furthermore, the surface returns in both cloudy and clear-sky areas make a DAR measurement of the total water column possible.

The DAR approach has two unique aspects that complement existing methods for remotely sensing water vapor. First, because of its ranging capabilities it has precise height registration, unlike passive sounding whereby weighting functions can encompass broad swaths of the atmosphere. Second, in contrast with other methods the DAR signal increases with increasing cloud water content and precipitation, with the obvious caveat that the radar signal-to-noise ratio (SNR) will decrease from attenuation as the beam penetrates into the volume. The DAR therefore nicely complements the infrared and microwave sounding techniques, as well as differential absorption and Raman lidar techniques that are commonly used to remotely sense water vapor from the ground (Spuler et al.2015; Whiteman et al.1992; Wulfmeyer and Bösenberg1998), with a notable airborne DIAL system being the Lidar Atmospheric Sensing Experiment (LASE) (Browell et al.1998). Importantly, millimeter-wave transparency in clouds allows for airborne or spaceborne measurements of lower tropospheric humidity in cloudy scenes, while DIAL systems typically cannot measure inside boundary-layer clouds due to high optical thickness.

In addition to primary applications in profiling water vapor within clouds, the instrument architecture discussed here represents an important application of recent advances in solid-state G-band technology to meteorological radar. Indeed, there has been lingering interest within the atmospheric remote sensing community for decades in utilizing G-band radar for cloud and precipitation studies, with earlier attempts hampered by limited sensitivity due to available technology (Battaglia et al.2014). The addition of G-band reflectivity measurements to multi-frequency radar systems, for example a dual-frequency W- and Ka-band system, could provide significantly more information than additional measurements at a lower frequency because the scattering properties at G band for typical cloud particle sizes are not of Rayleigh character.

Here we present ground-based measurements using a 167 to 174.8 GHz DAR, provide in-depth measurement error analysis with emphasis on the role of background noise power, and develop a retrieval algorithm based on performing least squares fits of a spectroscopic model to the data. The retrieved profiles constitute the first active remote sensing measurements of water vapor profiles inside of clouds, and open up possibilities for a variety of scientific studies, including investigation of in-cloud humidity heterogeneity and the coupled relationship between boundary-layer clouds and thermodynamic profiles.

2 Measurement basis and method

2.1 Differential absorption radar

The DAR technique (Cooper et al.2018; Lawrence et al.2011; Millán et al.2014) utilizes range-resolved radar echoes at multiple carrier frequencies in the vicinity of a gaseous absorption line to probe the frequency-dependent optical depth between two points along the radar line-of-sight. The radar echoes, or returns, may originate from cloud hydrometeors or, in the case of an airborne system, from the Earth's surface as well, enabling total column optical depth measurements. For closely spaced transmission frequencies near the absorption line center, the hydrometeor scattering properties vary little, while the gaseous absorption exhibits strong frequency dependence. By comparing with a known propagation model, these measurements can be employed to retrieve range-resolved density profiles of the absorbing molecule. Furthermore, because of the differential nature of the measurement, one does not require absolute calibration of the radar receiver in order to obtain absolute density values for the absorbing molecule. In the case of a calibrated receiver, both range-resolved density profiles of the absorbing molecule and microphysical properties of the reflecting medium can be retrieved.

Assuming negligible multiple scattering, the radar echo power received from a collection of scatterers filling the beam at a distance r is

(1) P e ( r , f ) = C ( f ) Z ( r , f ) r - 2 e - 2 τ ( r , f ) ,

where C(f) includes the frequency dependence of the radar hardware (e.g., transmit power and gain), Z(r,f) is the (unattenuated) reflectivity, and τ(r,f) is the one-way optical depth including contributions from gaseous and particulate extinction. Taking the ratio of powers for two different ranges r1 and r2=r1+R and assuming frequency independence of the reflectivity and particulate extinction, we find

(2) P e ( r 2 , f ) P e ( r 1 , f ) = Z ( r 2 ) Z ( r 1 ) r 1 r 2 2 e - 2 β ( r 1 , r 2 , f ) R ,



is the average absorption coefficient between r1 and r2, ρj(r) is the density of the gas component with label j, κj(r,f) is the corresponding mass extinction cross section, which varies with r due to pressure and temperature, and βpart(r) is the particulate extinction coefficient integrated over local drop size distributions (DSDs).

Figure 1Gaseous absorption coefficient (one-way) calculated using the model from Read et al. (2004) and the parameters listed in the figure. The green shaded region and inset highlight the 167 to 174.8 GHz transmission band for this work.


Restricting our analysis to millimeter-wave propagation near the 183 GHz water vapor absorption line, the sum over gaseous absorption terms can be replaced by ρv(r)κv(r,f)+βgas,bg(r), where the subscript v corresponds to water vapor and βgas,bg is the background gas absorption coefficient due to all other components, which is assumed to be frequency-independent. Assuming that pressure and temperature vary slowly compared to the length scale R, we can therefore write Eq. (3) as


where the overbar symbol implies taking the mean value between r1 and r2. Thus, we see that measuring the frequency-dependent contribution to the optical depth between r1 and r2 reveals the average water vapor density given the known absorption line shape κv(f). Figure 1 shows the frequency dependence of the gaseous absorption coefficient ρvκv(f)+βgas,bg in the vicinity of the 183 GHz water vapor line for P=1000 mbar, T=285 K, and ρv=10 g m−3. For this work, we utilize the millimeter-wave propagation model from the EOS Microwave Limb Sounder (Read et al.2004). The 167 to 174.8 GHz transmission band is highlighted in green, as well as shown in the inset of Fig. 1, revealing a differential absorption coefficient of 3 dB km−1 for 10 g m−3 of water vapor.

Important to the validity of this DAR method is the dominance of gaseous differential absorption over particulate differential absorption, since we assume that βpart is frequency-independent. To investigate this for boundary-layer clouds, we perform Mie scattering calculations for liquid spheres and integrate the scattering parameters over DSDs corresponding to clouds and rain for a range of mean diameters. For the DSD, we use a modified gamma distribution of the form

(5) N ( D ) = N 0 Γ ( ν ) D D n ν - 1 1 D n e - D / D n ,

where N0 is a normalization factor with units of particle number per volume that fixes the total liquid water content , ν is the shape parameter, which is set to 4 for clouds and 1 for rain, and Dn is the characteristic diameter. For rain we enforce an additional constraint that N0=x1Dn1-x2, where x1=26.2 mx2-4 and x2=1.57 have been determined in previous studies by comparing to observations (Abel and Boutle2012). This allows the entire rain distribution to be determined by the liquid water content. The rain rate is calculated from this distribution by using the terminal velocity relation from Beard (1976).

Figure 2Differential particulate extinction coefficients for (a) rain and (b) cloud. In the case of rain, the DSD characteristic diameter Dn determines the liquid water content and hence the rain rate, while for clouds we fix ℒ=500 mg m−3. The legend in (a) gives the corresponding frequencies in GHz.


The results are shown in Fig. 2, where we plot the differential particulate extinction, Δβpart(f)=βpart(f)-βpart(f0), as a function of Dn. Here f0=167 GHz corresponds to the low-frequency end of the transmission band. In Fig. 2a, the corresponding rain rate is displayed on the upper horizontal axis. For the cloud species, the normalization parameter N0 is not fixed by any additional constraint, and is therefore determined at each Dn to fix , which is set here to 500 mg m−3. To find the differential particulate extinction for other values of , one can linearly scale the values in Fig. 2b. Clearly for precipitation scenarios, the differential extinction from rain is more than 2 orders of magnitude smaller than that from water vapor. For clouds in the limit of small diameter, the differential particulate extinction asymptotes to the Rayleigh value of ΔβpartRayleigh=6πLIm(Kw)(f-f0)/(ρwc)=0.2 dB km−1 for f=174.8 GHz, where ρw is the density of liquid water, Kw=(mw2-1)/(mw2+2), mw is the complex refractive index of water, and c is the speed of light. For larger values of Dn, the differential extinction is enhanced by a resonant feature characteristic of Mie scattering. Thus, for thick clouds with as large as 500 mg m−3, especially those that contain drizzle drops which tend to lie near this resonant size, there are important bias considerations that warrant future study in order to establish the application of DAR in these particular scenarios. Specifically, to mitigate the potential biases stemming from scattering by hydrometeors, the unattenuated reflectivity can be used to distinguish clouds from precipitation, and the frequency-dependent scattering effects can be modeled and incorporated in the retrieval.

2.2 FMCW radar basics and instrument details

Due to the lower transmit power as compared to conventional radar systems at lower frequencies, the 170 GHz radar is operated in a frequency-modulated continuous-wave (FMCW) mode, which can offer increased sensitivity relative to a pulsed system with the same power because the transmitter is always on. The basic principle of FMCW radar is outlined in Fig. 3. The transmitted signal is frequency-modulated with a linear chirp waveform of bandwidth ΔFchirp and duration Tchirp. After scattering off of a target at a distance r from the radar, the received chirp is delayed in time by an amount 2rc, leading to a fixed frequency offset of δf=2ΔFchirpr/(cTchirp) relative to the transmitted frequency chirp. By downconverting the received signal using the transmitted frequency f(t) shifted by 5 MHz for convenient amplification and detection, the fixed frequency offset between transmitted and received chirps is converted into a constant frequency signal in the intermediate frequency (IF) stage. Signal processing techniques are then used to convert the IF time-domain signal to a range-resolved power spectrum. In the IF power spectrum, the zero-range point is located at 5 MHz and the echo power from a range R is located at fIF(r)=5 MHz ±δf(r), where the positive(negative) sign applies for decreasing(increasing) frequency chirps.

Figure 3Basic FMCW radar schematic. See Sect. 2.2 for discussion.


Our system utilizes state-of-the-art millimeter-wave components designed at the Jet Propulsion Laboratory (JPL) and builds on years of FMCW radar development for security and planetary science applications (Cooper et al.2011, 2017). The architecture is similar to that presented in an earlier work (Cooper et al.2018) which demonstrated the DAR technique between 183 and 193 GHz, but modified to transmit in the 167 to 174.8 GHz band, in which transmission is not prohibited by international regulations (NTIA2015), to perform narrow-bandwidth frequency chirps, and to provide a 5 MHz offset of the zero-range radar signal from zero frequency within the IF band. The IF offset is helpful for future calibrated power measurements because of various effects that inhibit accurate power estimation near zero frequency. The radar has an average transmit power of 140 mW, is outfitted with a 6 cm primary aperture with corresponding gain of 40 dB, and uses a frequency chirp of bandwidth ΔFchirp=60 MHz and duration Tchirp=1 ms, resulting in a range resolution of Δr=c/2ΔFchirp=2.5 m. In general, the choice of radar range resolution involves a compromise between acquiring more statistically independent samples within a given target volume to reduce uncertainty for bright targets and having a longer integration time to reduce noise power and thus increase the SNR for weak targets. The choice of 2.5 m allows us to downsample the range dimension by a factor of 11 to realize our desired profile resolution of 27.5 m, with decreased uncertainty for the bright clouds measured in this work. A summary of relevant radar hardware parameters is given in Table 1.

Table 1Hardware and radar signal acquisition parameters used in this work. The noise figure reported is for a complex radar signal detected using a double-sideband front-end mixer.

Download Print Version | Download XLSX

To process the downconverted radar signal, we first sample it using an analog-to-digital converter (ADC) with a sampling frequency of 20 MHz for the 1 ms duration of the chirp. Then we apply a Hanning window in the time domain before performing a fast Fourier transform (FFT) to obtain the range-resolved power spectrum. Application of the Hanning window reduces side lobes from bright targets as well as the large transmit–receive leakage signal that is always present at zero range. For the radar parameters listed above, the corresponding conversion factor from IF frequency to the target range is δf(r)/r=400 kHz km−1.

2.3 Power measurement uncertainty

The starting point for assessing the achievable precision in humidity using DAR measurements is the statistical uncertainty of the radar power measurements themselves. Until this point, we have ignored the role of background noise power in the radar spectrum, which is an important factor in any realistic receiver. In general, the noise power within a given radar range bin Pn is proportional to the sum of the receiver noise temperature and the antenna temperature, which itself is proportional to the scene brightness temperature. By considering the simultaneous coherent detection of noise (Pn) and radar echo (Pe) power, one can show that the statistical uncertainty of the detected power, Pd=Pe+Pn, is given by (see Appendix A)

(6) σ d = 1 N p P e 2 + 2 P e P n + P n 2 1 / 2 ,

where Np is the number of radar pulses transmitted.

In order to accurately determine the frequency-dependent optical depth between two range bins, it is critical to obtain a separate measurement of the background noise power in the absence of radar echoes and subtract this off of Pd. To see why this is, consider Eq. (2) with the left-hand side replaced by Pd(r2,f)/Pd(r1,f), which is equivalent to interpreting the detected power as the true echo power, set Z(r2)=Z(r1) for simplicity, and consider the limit PePn (i.e., PdPn). In this case we would find that exp(-2β(r1,r2,f)R)1 regardless of the actual value of Pe, and thus would incorrectly estimate a vanishing water vapor density, when in fact it is the echo power which has vanished. Similarly, for modest values of the SNRPe/Pn, this would lead to a systematic underestimate of the true humidity. Therefore, after subtracting the separate noise power measurement from Pd we obtain a measurement of Pe with total uncertainty σe=(σd2+σn2)1/2, where σn=Pn/Np is the noise power measurement uncertainty (see Eq. 6 with Pe=0). The relative uncertainty in the measured echo power is therefore

(7) σ e P e = 1 N p 1 + 2 SNR + 2 SNR 2 1 / 2 .

As will be discussed in Sect. 3, the range dimension is purposefully oversampled in our measurements, allowing us to decrease the statistical power uncertainty at a given range by averaging Nb adjacent range bins. The resulting relative power uncertainty is given by

(8) σ e P e = ξ ( N b ) N p N b 1 + 2 SNR + 2 SNR 2 1 / 2 ,

where ξ(Nb)≥1 is a factor of order unity accounting for covariances between adjacent range bins that arise due to applying a window function to the time-domain radar signal before transforming to Fourier space. For the Hanning window used in this work, this function is given by ξ(Nb)=1+Nb-1Nb891/2.

2.4 Inversion algorithm for profile retrieval

Under the simplifying assumptions introduced in Sect. 2.1, and assuming that pressure and temperature are known as a function of range, the inverse problem to retrieve humidity can be solved directly. The implications of the latter assumption are explored in Appendix C. To invert the radar spectra, we consider a set of measured echo powers Pe(ri,fj) for ranges {r1,r2,,rm} and transmission frequencies {f1,f2,,fNf}, where ri+1-ri=Δr is the radar range resolution. We note that in most circumstances we employ a retrieval step size R that is larger than Δr, since, as we will show below, the precision in our retrieved humidity scales favorably with total optical depth and hence with increasing R. Then, given a step size such that R=ri+S-ri for some integer S, we form the frequency-dependent measured quantity

(9) γ i ( f j ) = - 1 2 R ln r i + S r i 2 P e ( r i + S , f j ) P e ( r i , f j )

for each starting range ri. From Eq. (2), we see that we can extract the average humidity between ri and ri+S by performing a least squares fit of the function

(10) γ ^ ( f ) = ρ κ v ( f ) + B

to the measurements for each i, where B is a frequency-independent offset containing information about dry air gaseous absorption, particulate extinction, and the relative reflectivity of the two ranges in question. We drop the v subscript on the water vapor density in the above equation for simplicity of notation. The resulting humidity estimates {ρ1,ρ2,,ρm-S} have a corresponding range axis {r1,r2,,rm-S}, where ri=(ri+ri+S)/2, and have associated uncertainties determined from the fitting procedure.

Using standard error propagation, the estimated uncertainty in the measured quantity γi(fj) defined in Eq. (9) is

(11) σ γ i ( f j ) = 1 2 R σ e P e r i + S , f j 2 + σ e P e r i , f j 2 1 / 2 .

In order to derive a simple analytical expression for the relative uncertainty in the retrieved humidity, we restrict ourselves for the moment to considering two transmission frequencies, f1 and f2. In this case, we can combine Eqs. (2), (4), and (9) to obtain the humidity directly,

(12) ρ ( r i ) = κ v ( f 2 ) - κ v ( f 1 ) - 1 γ i ( f 2 ) - γ i ( f 1 ) ,

with the associated relative uncertainty


where Δτ=κv(f2)-κv(f1)ρ(ri)R is the differential optical depth for f1 and f2 between range bins ri and ri+S. Equation (13) reveals that there are three linked quantities determining the sensitivity of the system: (1) the magnitude of the DAR signal quantified by Δτ, (2) the statistical uncertainty of the power measurements given by the quadrature sum of relative errors in Eq. (13), and (3) the relative uncertainty in the derived value for the humidity. Thus, given a set of measured echo powers and a specific value for the humidity, there is a trade-off between spatial resolution of the retrieval and relative uncertainty in the humidity estimate.

An important and subtle point regarding the uncertainty in the measured quantity γi(fj) is that Eq. (11) relies on a Taylor expansion in the relative error σePe, and therefore is only valid for measurements with SNR above some critical value that depends on the number of measurements Np. Because there is no closed-form expression for the probability distribution function (PDF) of γi(fj), we resort to a Monte Carlo analysis, which is described in Appendix B, to generate the relevant PDFs for the parameters used in this work numerically. From this analysis, we find that for Np=2000 pulses and Nb=11 averaged bins, the Taylor expansion method is accurate for measurements with SNR>-10 dB.

We note here that it is typical of differential absorption systems to utilize only two frequencies: one online and one offline. However, in this work we are concerned with validating both the spectroscopic model used and the radar hardware itself, which could be subject to unknown frequency-dependent systematic effects. The regression approach discussed above thus provides for a robust comparison of the measured frequency dependence γi(fj) with the model γ^(f), while a two-frequency approach would mask inconsistencies between measurements and model, or systematic hardware effects, since the two free parameters ρ and B are fully determined given two frequency points. Furthermore, a distributed set of frequencies allows for the possibility of extending retrievals deeper in range for moist atmospheres, as frequencies closer to the line center will be attenuated more strongly, and can be excluded from the fits described above when the critical SNR value is reached.

Figure 4DAR measurement spectra. (a) The bidirectional frequency chirp technique provides for accurate, real-time characterization of the background noise floor within the radar's IF band, with no loss of measurement duty cycle. Here the detected power spectrum for fj=167 GHz is shown. The IF frequency to range conversion factor is 400 kHz km−1. (b) Echo power spectra normalized to their value at 100 m for the 12 transmission frequencies. The large variability in the signals near 1.4 km indicates the system reaching the noise floor. (c) Echo power spectra after averaging Nb=11 adjacent bins, and filtered for points with SNR>-10 dB. (d) Measurement relative error (blue circles) for all traces in (c) compared with the statistical model (Eq. 8, dashed black line).


3 Boundary-layer measurements and analysis

3.1 Radar characteristics, spectra, and filtering

In this section we report on measurements performed at JPL on 15 March 2018 using the proof-of-concept differential absorption radar described in Sect. 2.2. For these measurements, we implement a new signal processing technique for real-time noise floor characterization, utilizing a triangle-wave frequency chirp (i.e., bidirectional) instead of a sawtooth-wave chirp (i.e., unidirectional). According to FMCW radar principles, the echo spectrum switches from residing on the low- to the high-frequency side of the zero-range signal (i.e., 5 MHz) for increasing and decreasing linear frequency chirps, respectively. As shown in Fig. 4a, this fast switching of the chirp direction alternately exposes the noise floor on each side of the zero-range point within the IF band, and provides accurate and nearly continuous estimation of the system noise power and the passive signal corresponding to the scene brightness temperature at each frequency bin. This technique is especially advantageous for airborne/spaceborne applications, as the brightness temperature of the observed scene can change on fast timescales due to different surface types (e.g., ocean versus land) and from the presence or absence of clouds.

Figure 4 showcases a few aspects of a single ground-based DAR measurement, for which the conditions were light drizzle and a cloud located a few hundred meters off the ground. For all the field measurements discussed in this work, we acquire Np=2000 pulses for each of 12 frequencies equally spaced between 167 and 174.8 GHz, with the radar positioned just inside a building, pointing at 30 elevation. The experimental sequence is as follows: first, we perform 40 frequency chirps at a given transmission frequency before switching to another frequency, which takes 1 ms. The received signal is downconverted to baseband, digitized in an ADC, and processed in real time as described in Sect. 2.2. We achieve a system duty cycle of >90 %, resulting in a total measurement/observation time of ≈25 s.

By subtracting the respective noise floors from the increasing and decreasing frequency chirp measurements (Fig. 4a), and subsequently combining the mirrored spectra, we obtain our estimate of the echo power spectra. In Fig. 4b, we plot the echo power spectra scaled by r2 for the 12 transmission frequencies before bin averaging, which reveals the range dependence of the quantity Z(r)exp(-2τ(r,f)). Each spectrum is normalized to its value at 100 m. Thus, we observe the differential absorption due to water vapor directly from the spreading of the spectra with increasing range, whereby for a particular range, the plotted values increase monotonically with decreasing transmit frequency. After averaging the quantity ri2Pe(ri,fj) within a swath of size Nb=11, we filter the spectra based on the Monte Carlo analysis in Appendix B, keeping only those points with SNR>-10 dB, and are left with the smoothed profiles shown in Fig. 4c. Figure 4d shows the relative error in the binned (Nb=11) echo power measurement (blue circles) plotted against the measured SNR for all 12 frequencies. The measured values agree very well with those predicted by Eq. (8) (black dashed line), indicating that our statistical model based on speckle noise, which underlies the Monte Carlo simulations implemented in this work, is accurate.

Figure 5Water vapor profile retrieval for DAR spectra from Fig. 4. (a) Three examples of least-squares fits of the millimeter-wave propagation model to DAR measurements. Artificial offsets are imposed in order to plot all three on the same graph. (b) The retrieved profile exhibits roughly constant absolute humidity error until SNR≈10 dB (1 km). See Sects. 2.4 and 3.2 for retrieval details. The green line shows the saturated water vapor density range dependence using a near-surface temperature of 11 C and lapse rate of 6 C km−1. The shaded regions correspond to deviations of ±2C.


3.2 Water vapor profile retrieval

Using the averaged, filtered spectra in Fig. 4c, we proceed towards retrieving the water vapor density profile using the procedure outlined in Sect. 2.4. For the profiles presented in this section, we utilize a retrieval step size of R=200 m. Beginning with an initial range of r1=100 m, we form the 12 quantities γi(fj) for each starting ri in the set {r1,r2,,rm-S}, and perform a least-squares fit of the function γ^(f) to the data at each range point. Note that the retrieved water vapor density ρi is related only to the difference between the value of the fitted function at 174.8 and 167 GHz, while the offset is related to particulate extinction and hydrometeor reflectivity, and is disregarded in this work. The pressure and temperature dependence of the absorption line shape is included in the fitting model using reported values at the surface from a nearby weather station, and assuming an exponential pressure profile with a scale height of 7.5 km and a temperature lapse rate of 6 C km−1. We note that for the relatively short vertical extents of the profiles from these ground measurements (e.g., 1.4sin30 km for Figs. 4 and 5), the retrieved ρ values are quite insensitive to the assumed thermodynamic profiles (see Appendix C).

An important element of the DAR technique in general is utilizing an accurate model for the absorption line shape. Examples of line shape fits to the data are shown in Fig. 5a for three different values of SNR, with arbitrary offsets imposed on the three traces to permit simultaneous plotting. To assign SNR values to these points, we compute the mean SNR for the 12 frequencies at ri and ri+S, and use the smaller of the two. Clearly the millimeter-wave model accurately captures the frequency dependence of the measurements, which is supported quantitatively by the typical reduced chi-square values of χred21 for these fits. The retrieved water vapor density profile is shown in Fig. 5b, where the range ri assigned to each fitted value ρi is the midpoint of ri and ri+S. Also plotted here is an estimate of the saturation vapor density given our lapse rate assumption. This profile is consistent with a cloud base between 400 and 600 m and shows qualitatively good agreement with the expectation that the relative humidity is approximately 100 % in liquid cloud layers. Note that because the retrieved values correspond to the mean humidity between ri and ri+S, we effectively retrieve the profile convolved with a box of size R (200 m here). For this retrieval, the absolute humidity errors lie between 0.55 and 0.60 g m−3 until around 1 km (SNR≈10 dB), where the error steadily increases until the final retrieval point at 1.25 km with σρ=2.9 g m−3. The value of σρ in the high-SNR regime (i.e., the first 1 km) remains roughly constant, even though ρ varies by a factor of 3, since the absolute humidity error is independent of the humidity itself, and depends only on the differential mass extinction cross section κv(174.8 GHz) κv(167 GHz), the retrieval step size R, and the power measurement uncertainty (see Eq. 13).

Figure 6Retrievals for different cloud and precipitation conditions. (a) Averaged DAR power spectra (Nb=11) for light rain near the surface, with a cloud extending from 1 to 2 km. (b) Retrieved humidity profiles for two independent data sets corresponding to the same scene from (a). (c) Averaged DAR power spectra (Nb=11) for heavy precipitation near the surface with strong particulate extinction. (d) Independent retrievals from two data sets for the scene in (c).


Though we do not have independent, coincident water vapor profile measurements with which to validate the accuracy of the retrieval, we can investigate repeatability of this DAR method by performing the retrieval on coincident, independent DAR measurements of the same exact scene. To do so, we acquire Np=4000 pulses at each frequency with a total measurement time of 50 s, and parse the data into two groups of Np=2000 pulses both spanning the full 50 s. The results are shown in Fig. 6, where we also present measurements of different cloud and precipitation scenarios than that presented in Figs. 4 and 5. In Fig. 6, panels (a) and (b) correspond to light rain at the surface with a cloud boundary at 1 km range, and panels (c) and (d) to heavy rain at the surface with strong particulate extinction. The retrievals from the two independent sample sets in both cases agree quite well, which showcases the reproducibility of the measurement and indicates that the estimated humidity error accurately captures the sample scatter.

Given a measured range-resolved echo power spectrum, what retrieval range resolution can we achieve for a specified minimum retrieval precision? As discussed briefly in Sect. 2.4, the relative error in the retrieved humidity σρ/ρ (see Eq. 13) for a given power measurement uncertainty varies inversely with the differential optical depth, and thus depends on both the retrieval step size R (i.e., retrieval resolution) used and the absolute value of the humidity ρ. Alternatively, one can look at the absolute error and rearrange Eq. (13) to find that, for a given pair of frequencies and power measurements at two ranges, the product of σρ and R is constant. Hence, reducing the retrieval step size by some factor increases the absolute humidity error by the same factor. In future work we will implement a retrieval algorithm that has adaptive range resolution based on both the inherent signal (i.e., humidity) and the measurement noise.

4 Conclusions

A proof-of-concept humidity-profiling DAR operating between 167 and 174.8 GHz has been constructed and tested from the ground. The instrument builds on progress made in an earlier version operating between 183 and 193 GHz (Cooper et al.2018), and employs a new signal processing technique for performing real-time noise power spectrum characterization and subtraction, providing for higher accuracy measurements of the radar echo power. A new direct inversion algorithm for retrieving humidity based on least squares fits to a spectroscopic model is applied to the measured echo power spectra, showing close agreement between the measurement and model frequency dependence. The humidity profiles retrieved from two statistically independent measurement sets of the exact same scene are in close agreement, highlighting the reproducibility of the method. The uncertainties in the power measurements, which in part determine the retrieved humidity uncertainty, agree very well with a statistical model based on radar speckle noise that incorporates the effects of background noise subtraction and downsampling, or binning, of the measured spectra.

Development of an operational airborne 167–174.8 GHz DAR is currently in progress, which will include an additional 20 dB of antenna gain and a factor of 4 increase in transmit power. Important future steps for this instrument include validation of the measurement accuracy using coincident measurements of humidity, pressure, and temperature (e.g., from radiosondes), and eventually testing from an airborne platform. Specifically, the surface returns while measuring from an airborne platform will be investigated for the retrieval of total column water within the boundary layer. A more significant augmentation of the system could include the addition of passive radiometric channels near the 183 GHz line. This would allow for continuous measurement of vertical humidity profiles when transitioning between clear-sky and cloudy areas, and opens the possibility to study biases in the humidity retrieved from radiometric measurements that are caused by scattering and emission from clouds.

Data availability

All of the data used in this paper are presented in the figures.

Appendix A: Error in detected power

In this appendix, we derive the expression for the detected power uncertainty within a single radar range bin in the presence of background noise (Eq. 6). To do so, we begin by assuming that all targets within the scattering volume are randomly distributed, leading to the well-known Rayleigh fading model for the received echo signal (Ulaby et al.1982). In the context of FMCW radar, we then consider the received complex electric field amplitude Ei corresponding to the ith frequency bin in the FFT spectrum, where we only consider the polarization direction that couples into the radar receiver. Within the Rayleigh fading model, it is shown that for Ei=E0,ieiϕi, the modulus of the field amplitude E0,i is normally distributed with zero mean and standard deviation σE, and the phase ϕi is uniformly distributed over the interval [0,2π]. Alternatively, we can write the corresponding voltage in the receiver as Ve,i=ai+ibi, where ai and bi are uncorrelated and are both normally distributed with zero mean. Then, from the expression converting electric field to power, Pe,i=|Ve,i|2=α|Ei|2, we find the probability distribution function for the received echo power

(A1) p ( P e , i 0 ) = 1 P e , i e - P e , i / P e , i ,

where the mean equals the variance and is given by Pe,i=2ασE2, α is a field-to-voltage conversion factor for the radar, and p(Pe,i<0)=0. Furthermore, we find that ai2=bi2=Pe,i/2. Though not proven here, the Rayleigh fading model also shows that the expectation value of the received power from N randomly distributed targets is the sum of the expectation values of the individual target echo powers.

Similarly, one can show that Gaussian white noise in the radar signal, which comes from both the scene brightness temperature and the radar electronics, results in a noise voltage within the ith frequency bin of the FFT spectrum with Fourier coefficient Vn,i=ci+idi, where ci=di=0 and ci2=di2=Pn,i/2. We proceed towards deriving Eq. (6) by considering the coherent detection of both the radar echo and noise signals. In this case, the detected voltage signal in the Fourier domain within the ith range bin is Vd,i=Ve,i+Vn,i, and the detected power is Pd,i=|Vd,i|2. Using the expectation values listed above, it is easy to show that

(A2) P d , i = P e , i + P n , i


(A3) Var ( P d , i ) = P e , i + P n , i 2 .

Therefore, we recover Eq. (6) by computing the standard error for N independent measurements, σd,i2=Var(Pd,i)/N.

Appendix B: Monte Carlo Analysis

As discussed in Sect. 2.3, subtracting off the noise power contribution to the detected power Pd is critical for accurate humidity estimation. However, for low values of SNR, we clearly expect the result Pe=Pd-Pn to be negative some of the time due to finite sampling, where “〈⋯〉” denotes the sample average. This is nonphysical. Therefore, in order to account for these finite sampling effects and the potential breakdown of standard error propagation when σePe is not small, we employ a Monte Carlo simulation of the DAR measurement. The PDF for the echo power received from randomly distributed hydrometeor targets within a single range bin is given by Eq. (A1). In order to generate random samples of the radar spectrum for transmission frequency fj, we begin with an idealized spectrum Pe(ri,fj)〉 for which we set C(f)=Z(r,f)=1 and ρv(r)=constant, and sample the distribution (Eq. A1) at each range bin ri. We then perform a fast Fourier transform (FFT) to obtain the corresponding time-domain radar signal, add the effects of background noise using Gaussian white noise, and apply a Hanning window. Taking the inverse FFT thus supplies a single random realization of a measured Pd spectrum. For the simulated spectra used in this work, we generate Np=2000 radar pulses with range resolution Δr=2.5 m and average them to realize a single radar measurement. These values for Np and Δr are the same parameters utilized in the field measurements presented in Sect. 3. We generate 10 000 averaged spectra for both Pd and Pn, giving 10 000 random realizations of the echo power measurement. For these simulations, we use fj=167 GHz and ρv=7.4 g m−3.

Figure B1Statistics of low-SNR estimation of two-way transmission. The blue data points represent the means and standard deviations of the Monte Carlo probability distributions, while the gray shaded area represents the error calculated using Eq. (8). For SNR<-10 dB, the systematic bias of the Monte Carlo mean and the underestimation of the true error by Eq. (8) imply that measurements in this region should be disregarded. In order to access lower values of SNR, one must increase the number of pulses Np for each measurement.


Our aim is to utilize the Monte Carlo simulations to inform where the Taylor expansion method for error propagation breaks down in our estimation of σγi(fj), and thus provide a criterion for filtering our measurements. To do so, we fix Nb=11 and the step size S=10 (i.e., 275 m) and compute the mean and standard deviation of the Monte Carlo probability distribution for the two-way transmission between ri and ri+S for each ri. Figure B1 shows the results, where we plot the Monte Carlo mean value divided by the a priori two-way transmission used to generate the Monte Carlo results, as a function of the SNR at ri+S. The gray shaded area represents the SNR-dependent errors predicted from Eq. (8). There are two notable deviations that arise for SNR values below 0.1: (1) the error estimated using the standard error propagation formalism begins underestimating the true standard deviation calculated using the Monte Carlo ensemble, and (2) the mean of the Monte Carlo-generated distribution systematically overestimates the true two-way transmission. We note here that this point of departure between the naive error propagation estimate and that from the Monte Carlo distributions does not depend on Nb or S, but is determined by the number of independent pulses Np used to realize a single radar measurement. From these simulations, we conclude that the standard error propagation model is sufficient for SNR>-10 dB. Therefore, after downsampling the measured spectra with Nb=11, we eliminate all measured values with SNR<-10 dB, as described in Sect. 3.1.

Appendix C: Retrieval dependence on assumed pressure and temperature values

To assess the dependence of the retrieved humidity on temperature and pressure, we will consider again the case of the two-frequency measurement, using transmission frequencies f1=167 GHz and f2=174.8 GHz. Then, for a given starting range ri and step size R, we use the measured quantities γi(f1) and γi(f2) to solve for the mean humidity between the two ranges,

(C1) ρ i = γ i ( f 2 ) - γ i ( f 1 ) κ v ( f 2 , P , T ) - κ v ( f 1 , P , T ) = γ i ( f 2 ) - γ i ( f 1 ) Δ κ v ( P , T ) ,

where we now explicitly write κv as a function of temperature T and pressure P, and we have defined the differential mass extinction cross section Δκv(P,T) for these two frequencies. Given reference values of P0=1000 mbar and T0=285 K, and corresponding retrieved humidity ρi,0, we calculate the error in our humidity estimate for different conditions P and T as

(C2) ρ i - ρ i , 0 ρ i , 0 = Δ κ v ( P 0 , T 0 ) Δ κ v ( P , T ) - 1 .

Figure C1 shows the humidity error for pressure deviations of ±20 mbar and temperature deviations of ±15 K. Here we see that the retrieved humidity is very weakly dependent on the assumed pressure, and only accrues an error of 10 % for a temperature deviation of about 8 K.

Figure C1Humidity error for assumed temperature T0=285 K and pressure P0=1000 mbar versus actual values T and P.


Competing interests

The authors declare that they have no conflict of interest.


This research was supported by NASA's Earth Science Technology Office under the Instrument Incubator Program, and was carried out at the Jet Propulsion Laboratory (JPL), California Institute of Technology, Pasadena, CA, USA, under contract with the National Aeronautics and Space Administration. Richard J. Roy's research was supported by an appointment to the NASA Postdoctoral Program at JPL, administered by Universities Space Research Association under contract with NASA.

Edited by: Murray Hamilton
Reviewed by: Gerald Mace and two anonymous referees


Abel, S. J. and Boutle, I. A.: An improved representation of the raindrop size distribution for single-moment microphysics schemes, Q. J. Roy. Meteor. Soc., 138, 2151–2162, 2012. a

Battaglia, A., Westbrook, C. D., Kneifel, S., Kollias, P., Humpage, N., Löhnert, U., Tyynelä, J., and Petty, G. W.: G band atmospheric radars: new frontiers in cloud physics, Atmos. Meas. Tech., 7, 1527–1546,, 2014. a

Beard, K. V.: Terminal Velocity and Shape of Cloud and Precipitation Drops Aloft, J. Atmos. Sci., 33, 851–864, 1976. a

Browell, E., Ismail, S., and Grant, W.: Differential absorption lidar (DIAL) measurements from air and space, Appl. Phys. B, 67, 399–410, 1998. a

Browell, E. V., Wilkerson, T. D., and Mcilrath, T. J.: Water vapor differential absorption lidar development and evaluation, Appl. Optics, 18, 3474–3483, 1979. a

Cooper, K. B., Dengler, R. J., Llombart, N., Thomas, B., Chattopadhyay, G., and Siegel, P. H.: THz Imaging Radar for Standoff Personnel Screening, IEEE T. Thz. Sci. Techn., 1, 169–182, 2011. a

Cooper, K. B., Durden, S. L., Cochrane, C. J., Monje, R. R., Dengler, R. J., and Baldi, C.: Using FMCW Doppler Radar to Detect Targets up to the Maximum Unambiguous Range, IEEE Geosci. Remote S., 14, 339–343, 2017.  a

Cooper, K. B., Monje, R. R., Millán, L., Lebsock, M., Tanelli, S., Siles, J. V., Lee, C., and Brown, A.: Atmospheric Humidity Sounding Using Differential Absorption Radar Near 183 GHz, IEEE Geosci. Remote S., 15, 163–167, 2018. a, b, c, d

Lawrence, R., Lin, B., Harrah, S., Hu, Y., Hunt, P., and Lipp, C.: Initial flight test results of differential absorption barometric radar for remote sensing of sea surface air pressure, J. Quant. Spectrosc. Ra., 112, 247–253, 2011. a, b

Lebsock, M. D., Suzuki, K., Millán, L. F., and Kalmus, P. M.: The feasibility of water vapor sounding of the cloudy boundary layer using a differential absorption radar technique, Atmos. Meas. Tech., 8, 3631–3645,, 2015. a

Millán, L., Lebsock, M., Livesey, N., Tanelli, S., and Stephens, G.: Differential absorption radar techniques: surface pressure, Atmos. Meas. Tech., 7, 3959–3970,, 2014. a

Millán, L., Lebsock, M., Livesey, N., and Tanelli, S.: Differential absorption radar techniques: water vapor retrievals, Atmos. Meas. Tech., 9, 2633–2646,, 2016. a

NTIA: Manual of Regulations and Procedures for Federal Radio Frequency Management, National Telecommunications & Information Administration (NTIA), Revision of the May 2013 Edn., 2015. a, b

Read, W. G., Shippony, Z., and Snyder, W.: EOS MLS forward model algorithm theoretical basis document, Jet Propulsion Laboratory, JPL D-18130/CL#04-2238, Pasadena, CA, USA, 2004. a, b

Spuler, S. M., Repasky, K. S., Morley, B., Moen, D., Hayman, M., and Nehrir, A. R.: Field-deployable diode-laser-based differential absorption lidar (DIAL) for profiling water vapor, Atmos. Meas. Tech., 8, 1073–1087,, 2015. a

Ulaby, F., Moore, R., and Fung, A.: Microwave Remote Sensing: Active and Passive, Vol. II, Addison-Wesley Publishing Company, Inc., Reading, Massachusetts, 1982. a

Whiteman, D. N., Melfi, S. H., and Ferrare, R. A.: Raman lidar system for the measurement of water vapor and aerosols in the Earth's atmosphere, Appl. Optics, 31, 3068–3082, 1992. a

Wulfmeyer, V. and Bösenberg, J.: Ground-based differential absorption lidar for water-vapor profiling: assessment of accuracy, resolution, and meteorological applications, Appl. Optics, 37, 3825–3844, 1998. a

Short summary
The measurement of water vapor profiles inside clouds with high spatial resolution represents an outstanding problem in atmospheric remote sensing. Here we present measurements from a proof-of-concept millimeter-wave (170 GHz) cloud radar aimed at filling this observational gap, and demonstrate the ability to retrieve in-cloud water vapor profiles with high precision and resolution. This technology could meaningfully impact future satellite-based measurements of water vapor.