the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 31 Jul 2026
Spectral analysis of GNSS zenith wet delay and microwave radiometer precipitable water during FESSTVaL Lindenberg 2021: a multi-scale observational test of horizontal heterogeneity sensitivity
Abstract. Every Global Navigation Satellite System (GNSS) station estimates a zenith wet delay (ZWD), a direct measure of the water vapour integrated above it. The ZWD is a standard product of GNSS meteorology, yet its small-scale temporal structure is conventionally treated as noise. We ask whether the temporal spectrum of the ZWD reveals horizontal water vapour heterogeneity near the boundary layer top that a vertical sounding cannot see. Using FESSTVaL observations at Lindenberg in summer 2021, we compare the GNSS ZWD with the precipitable water from a microwave radiometer over 28 days. Both are described by a Matérn spectrum whose scale and amplitude are tracked on sliding windows. The radiometer scale follows the boundary layer height in calm conditions and shifts to a wind-driven, advective behaviour when the atmosphere becomes unstable, confirmed by tower turbulence. The GNSS scale stays stable across all conditions, because its oblique geometry ties it to horizontal structure at the boundary layer top, unchanged by the shifts in vertical stability. The strength of the GNSS signal grows with the horizontal heterogeneity measured independently at two scales (ρ = +0.433 and +0.525), while the vertical radiometer shows none and a scintillometer rules out a surface contribution. Accessible from any positioning solution, the GNSS scale measures horizontal coherence at the boundary layer top that is far more stable than the ERA5 boundary layer height (coefficient of variation 0.32 against 0.90) and persists through the night.
- Preprint
(3983 KB) - Metadata XML
- BibTeX
- EndNote
Status: open (extended)
-
RC1: 'Comment on egusphere-2026-4148', Anonymous Referee #1, 14 Sep 2026
reply
-
AC1: 'Reply on RC1', gael kermarrec, 15 Sep 2026
reply
RESPONSE TO REFEREE
Manuscript: Spectral analysis of GNSS zenith wet delay and microwave radiometer precipitable water during FESSTVaL Lindenberg 2021: a multi-scale observational test of horizontal heterogeneity sensitivity
We thank the referee for a careful and constructive report. The central request is well taken. The high-frequency variability of the ZWD is the fundamental observable of this study. Its sensitivity to the estimation strategy should therefore be established, not assumed. The constraint that governs it should be presented and motivated. We have documented the processing in full and carried out the reprocessing the referee suggested. A dedicated appendix reports it. Referee comments are reproduced in italics, followed by our reply. Line numbers refer to the revised manuscript.
GENERAL COMMENTS
================G1. Constraints in the ZWD estimation, and their impact on the variability
Although the temporal resolution of the ZWD estimates is 30 s, I assume that there are constraints in the estimation process ... Adjacent values are not independent estimates. The model and the constraint are of fundamental importance and shall be motivated and presented (see, e.g. Young et al., 2024). Furthermore, because the constraints used in the GNSS processing have a significant impact on the resulting variability it is also necessary to determine how significant it is.
We agree on both counts, and we state the point explicitly in the revised manuscript rather than leaving it implicit. Successive 30 s ZWD values are not independent estimates: the process noise of the random-walk model sets how much high-frequency variability the filter admits, and that variability is precisely what we analyse.
Model and constraint, now presented (Sect. 2.2). The ZWD is estimated as a random walk (the STO option of PRIDE PPP-AR) with an a priori standard deviation of 0.20 m and a process noise of 0.0300 m s^(-1/2); horizontal gradients are estimated jointly with a process noise of 0.005 m s^(-1/2). Section 2.2 now reports these values, together with the constellations used, the elevation cutoff, the satellite products and the ambiguity-resolution settings (see G4 and S2 below). We cite Young et al. (2024) at this point, whose variable random-walk approach makes the same argument that the process noise is a physical choice, not a numerical detail.
Motivation of the value. The constraint cannot be chosen neutrally for this application: it must admit the variability that the atmosphere actually produces in the band the Matérn fit targets, between roughly 1/(62.5 min) and the 30 s sampling. A constraint tightened well below that requirement no longer returns atmospheric structure, it returns the filter's own transfer function. Our value was calibrated against this requirement on an independent station. We also checked whether it is unusual among operational solutions. Applying our estimator to the output of a GFZ EPOS solution over June days gives an equivalent process noise between 0.0200 and 0.0271 m s^(-1/2), somewhat below our 0.0300 but of the same order. We claim no more from this than it supports: our solution admits slightly more high-frequency variability than that chain, not a categorically different amount. The argument that the ranks survive a change of constraint is made by the reprocessing below, not by this comparison.
Sensitivity, now quantified (new Appendix A). Following the referee's suggestion we reprocessed the full 28 days of LDB2 at three additional constraints. Two of them, 0.0150 and 0.0600 m s^(-1/2), are half and double the value used. The third, 0.0050 m s^(-1/2), is deliberately over-tightened. Every other setting was left unchanged. Each solution went through the same spectral chain. We state the prediction before the outcome. Because the constraint acts as a low-pass filter, the absolute value of the fitted scale must vary monotonically with it. Our conclusions, however, rest on rank correlations and on contrasts between regimes, never on absolute values: the four results are the CAPE-regime shift of the radiometer scale, the corresponding stability of the GNSS scale, the surface null test, and the correlations with independently measured heterogeneity. A monotone rescaling leaves a Spearman rank correlation unchanged, so the prediction was that the medians would shift while the rank structure would not.
Process noise (m s^(-1/2)) median λ (s^(-1)) relative to 0.0300 median σ (mm) relative to 0.0300 ρ against reference
========================== ================= ================== ============= ================== ===================
0.0050 0.00291 0.19 1.39 0.90 0.838
0.0150 0.01006 0.67 1.51 0.98 0.977
0.0300 (used here) 0.01504 1.00 1.54 1.00 1.000
0.0600 0.01727 1.15 1.56 1.01 0.986Over the range the referee proposed, from half to double the value used, the fitted scale changes by a factor of 1.72 while the fitted amplitude changes by 3 percent. The window-by-window rank correlation against the reference solution is 0.977 and 0.986. The rank ordering of the observable is thus preserved to better than 0.98 across a factor of four in constraint. All four solutions yield the same 8052 spectral windows. The constraint therefore does not alter the temporal segmentation. The comparisons are made on the same windows. The outlier screening removes no window at 0.0050, 0.0150 and 0.0300, and two at 0.0600. This also answers the referee's related question at S6 on the screened fraction. One point the referee may notice deserves a word. The count of 8052 is that of the continuous GNSS fit, while the main analysis reports 8001 windows. The 8001 follow from aligning the series on the 5 min grid common to the radiometric records. The difference is one of temporal matching. No window is discarded. The appendix now says so.
The over-tightened case points the same way. Take 0.0050 m s^(-1/2), six times tighter than the value used. The median scale falls to 0.19 of its reference value, as a filter suppressing variance inside the fitted band must do. Yet the amplitude decreases by only 10 percent. The rank correlations stay high: 0.838 for λ, 0.964 for σ. Tightening the constraint compresses the absolute values well before it disturbs their ordering, which is the quantity the analysis uses. The bottom line is short. Halving or doubling the constraint moves the fitted values. It does not move their ordering. Since every conclusion of the paper rests on orderings, none of them changes.
A test with an independently produced solution. This is the stronger of our two comparisons with outside processing. The EPOS comparison above only places our constraint among operational values. The comparison below tests something else: whether the rank structure survives a change of chain. A ZWD series for LDB2 over the same 28 days was produced independently for the FESSTVaL campaign, by a different operator and at a constraint of 0.0100 m s^(-1/2). Fitting it with the identical estimator and matching windows to ours leaves 7656 common windows, on which the Spearman rank correlation between the two fitted scales is +0.879, with a median ratio λ(0.0300)/λ(0.0100) of 2.33. The rank ordering of the observable therefore survives a change of processing chain at a factor of three in constraint, which is precisely the property the conclusions require, and it is a stronger statement than varying one parameter within a single setup.
The independent comparison corroborates this from outside our own processing. The median ratio it measures, 2.33 between the constraints 0.0300 and 0.0100, sits as expected on the monotone curve of the table above, so the constraint dependence is reproducible across processing chains as well as within ours.
G2. Treatment of observations during rain
Another general observation is that it is not clear how observations during rain are used, if they are used at all. The microwave radiometer produce more or less useless data during rain.
The referee's own observation resolves this. Our wording obscured it. No window is excluded from the analysis on the basis of precipitation, and none needs to be. Since the ZWD is essentially unaffected by rain at the intensities encountered here, as the referee notes, a rainy window is a valid ZWD window: the spectral characterisation of water vapour structure proceeds through precipitation exactly as it does in clear air. The sentence in Sect. 2.8 announcing that the RADOLAN events were "used to filter out windows with non-turbulent rain contamination in the GNSS PPP" described an operation that is neither performed nor required. It has been removed. RADOLAN is used as context, to document when precipitation occurred within the 28-day window.
This is the substance of the "all weather" statement in the introduction, which we have made explicit rather than left to be inferred. The claim is not that the delay measurement survives rain, which is well known, but that the characterisation of water vapour organisation can be carried out continuously, including during precipitation and through the night. That is the property which distinguishes the GNSS observable from every other instrument used here: the radiometer retrieval is compromised by liquid water on the radome, the scintillometer data are flagged out in rain. The radiosondes are four-times-daily snapshots. A continuous record of moisture structure that does not stop when it rains is precisely what the spectral approach offers. It is most valuable in the conditions where the alternatives fail.
G3. Overall context for future applications, and E-GVAP
the manuscript would benefit from a discussion on the overall context for future applications ... The IGS networks is mentioned but not the E-GVAP, run by EUMETNET.
Added to the Outlook (Sect. 5.2). We now note that the operational base in Europe is much larger than the IGS network. E-GVAP, run by EUMETNET, collects tropospheric delays from several thousand ground stations in near real time, for assimilation by the national weather services. That is the natural route by which a spectral product of this kind would reach forecasting practice. We have kept the discussion of a real-time forecasting application deliberately measured. The present study establishes an observational association at a single site; it does not demonstrate predictive skill. The test against convective initiation from radar is identified as the next step, not claimed here.
G4. Reference added
Young, Z. M., Blewitt, G., and Kreemer, C.: Improved GPS tropospheric path delay estimation using variable random walk process noise, J. Geod., 98, 89, https://doi.org/10.1007/s00190-024-01898-3, 2024, is now cited in Sect. 2.2.
SPECIFIC COMMENTS
=================S1. L2-3, "conventionally treated as noise"
does not agree with my experiences from the GNSS meteorology community ... This text can be deleted in order not to upset readers that are eager to improve the true temporal resolution in GNSS meteorology.
Point accepted, and the wording was careless on our part: there is indeed an active effort within GNSS meteorology towards higher effective temporal resolution. We have replaced the claim with the statement that the small-scale temporal structure is "rarely exploited as a signal in its own right", which is what we meant and which does not dismiss that work.
S2. L70+, GNSS observations and processing
Are all four GNSS used? What is the elevation cutoff angle? ... it would be a quick task to compare it with your time series for validation?
Section 2.2 now specifies the processing. All four global constellations enter the solution (GPS, GLONASS, Galileo and BeiDou, both BeiDou-2 and BeiDou-3), on the ionosphere-free combination of two frequencies per system. QZSS is enabled in the configuration but contributes no observations at the latitude of Lindenberg. We state this rather than list it among the systems used. Observations below an elevation of 7 degrees are removed in data editing. Orbits, clocks, biases and Earth rotation parameters are the CODE multi-GNSS final products.
The cutoff also corrects a number in the original manuscript. Section 4.2 stated elevation angles "between 5 and 90 degrees"; the editing cutoff is 7 degrees. The horizontal offset of a ray at the boundary layer top is correspondingly about 16 km at the cutoff rather than the "some twenty kilometres" originally written. The revised Discussion gives the offset at three elevations, 3.5 km at 30 degrees, 11 km at 10 degrees and 16 km at 7 degrees for z = 2000 m. The cutoff is now shown in Fig. 6 (see T7).
We have not added the comparison against the IGS product, because it cannot do the work the referee intends for it. The IGS combined ZTD is delivered at 5 min. Our Matérn fit is performed on 62.5 min windows of 30 s data. The spectral slope and saturation scale it estimates are set by content between roughly 1/(25 min) and the Nyquist frequency of the 30 s sampling. A 5 min product has no content above 1/(10 min): the entire band that carries our observable has been removed from it before the comparison starts. Agreement would therefore demonstrate only that the two solutions share the same slowly varying column, which no one disputes and which is not what the paper is about, while disagreement would be uninterpretable, since it could arise from the smoothing alone. The comparison that does bear on the high-frequency structure is the one against the independently produced FESSTVaL solution at 30 s and a different constraint, reported under G1, and we have relied on that instead.
S3. L77, citation for the LIN0 data quality statement
The statement was an unsupported appreciation. We have replaced it with the verifiable facts. LDB2 is a GRUAN reference station operating since 2005, with its antenna on a thermally protected 1.7 m steel mast. LIN0 was installed only in 2020. Over the 28 days analysed, LIN0 delivers about 2 percent fewer epochs than LDB2. The revised text gives these instead of a judgement on data quality. We also note that no conclusion of the paper depends on LIN0: it enters only the indicative bounding of the processing-related fraction of the excess variance, as already stated in Sect. 4.3.
S4. L82, why the Bevis coefficient rather than the ZWD provided by the radiometer
Why use the Bevis coefficients (which were derived from stations in North America) when the RPG instrument provides the zenith wet delay directly?
Because the comparison here is spectral rather than absolute, and the requirement is a conversion whose effect on the spectrum is known exactly. The Bevis mapping is a single multiplicative constant: it rescales the spectral amplitude by the square of the coefficient and leaves the fitted scale unchanged, as we now state explicitly after Eq. (1). The ZWD product delivered by the manufacturer (Radiometer Physics GmbH) is a separate retrieval. It may be site-tuned. It is not documented in a form that lets us establish its transfer function in the band of interest. Using it would introduce an unknown filter into the very quantity under study. We prefer a conversion whose spectral action is transparent, even if its coefficient is not optimal for the site, because an error in the constant affects only the amplitude scaling and not the shape of the spectrum or any of the rank correlations. Substituting it here would be counterproductive in a specific way. Retrieval algorithms of this kind are regression or optimal-estimation schemes whose weights are trained on radiosonde climatologies and which may include temporal smoothing or quality-dependent weighting; none of that is documented at a level that would let us state its frequency response. Feeding such a product into a spectral fit means fitting the convolution of the atmosphere with an undocumented filter. Any difference from the GNSS spectrum could then be attributed either to the atmosphere or to the retrieval, with no way to separate them. That is the one thing this analysis cannot afford, because the comparison of the two spectra is the whole measurement. The Bevis constant may well estimate the absolute column at Lindenberg less well than the manufacturer product. That is of no consequence here. A constant error in the coefficient shifts the amplitude and leaves the spectral shape, the fitted scale and every rank correlation untouched. Optimality of the retrieval and transparency of its spectral action are different requirements, and this study needs the second.
S5. L149, rain filtering in the GNSS PPP
Corrected, see G2. The GNSS estimates are not filtered for rain. The sentence implying otherwise has been removed.
S6. L173-175, percentage of windows where the MWR noise exceeds the GNSS noise
it will be informative to state the percentage of data when the MWR noise exceeds the GNSS noise.
This figure was already reported, and we have made it more prominent since the referee is right that it is the informative number. The excess variance is floored to zero on the windows where the converted radiometer amplitude exceeds the GNSS amplitude, which occurs on about 7 % of windows (Sect. 3.2). The floored windows are retained in the correlation analyses, so that no selection on the sign of the difference is introduced. The empirical ratio of the two amplitudes is about 3.7, which places most windows well away from the crossover. We agree with the referee that this ratio depends on the constraint used in the GNSS processing, and the new sensitivity assessment reports the floored fraction at each constraint (G1).
TECHNICAL CORRECTIONS
=====================Comment Action
================================== ==============================================================================================================================================================================================================================
"cadence" throughout Replaced by "temporal resolution" everywhere (six occurrences).
L10, define ρ ρ is now identified as the Spearman rank correlation at first use in the abstract. The collision the referee noticed is also resolved: ρ_w in Eq. (1) is now explicitly defined as the density of liquid water (see next row).
L84, define all symbols in Eq. (1) All symbols now defined: ρ_w the density of liquid water, R_v the specific gas constant of water vapour, k_2' and k_3 the refractivity constants, T_m the water-vapour-weighted mean temperature of the column.
L114, "data sets" Changed to "datasets".
L148, sentence beginning with "18" Reworded.
Figure 5, label positions Corrected, see T6 below.
Figure 6, MWR sensed volume Corrected, see T7 below.
Reference list, journal names All journal names given in full where the abbreviation is not in the Web of Science list.T6. Figure 5
The referee is right, and the figure was worse than a matter of placement: in the vector version the axis labels were broken up character by character and the exponents of the logarithmic tick labels were detached from their mantissas. Tracing this, the cause lies in the export, not the plotting. The axis text is typeset with the TeX interpreter in a font that the headless MATLAB used for the figures substitutes silently, so the vector PDF declares one font while the character positions were computed with the metrics of another. The figure has been regenerated. Figure 5 is now supplied as the 300 dpi raster produced by the same script, where the labels, the Greek symbols and the logarithmic ticks all render correctly. Three of the other result figures were already supplied in that form. We note that the panels themselves, and every number in them, are unchanged.
While checking the figure we also corrected its caption. Panel (b) was described as showing the systematic shift. That overstated what a reader sees. On an axis starting at zero the three tercile medians lie close together. The regime shift is carried by the change in the correlations reported in the text, not by the magnitude of the median. The caption now says so. We preferred this to tightening the axis, which would have made a modest effect look larger than it is.
T7. Figure 6
The referee is right on the geometry and we are grateful for the correction, which we had not considered. The MWR beam is in the far field over essentially the whole column. With a half-power beam width of about 4 degrees and an aperture of about 20 cm at K band, the near field ends within some metres of the instrument. The sensed volume is therefore a cone, not a cylinder, with a half-power diameter of about 140 m at a height of 2 km. Figure 6 has been redrawn with the radiometer volume as a cone of the correct opening. The GNSS elevation cutoff of 7 degrees is now marked as well, as suggested. The correction does not affect the analysis, since the argument rests on the contrast between a narrow vertical sounding and an oblique multi-ray integration spanning kilometres at the same height. A 140 m cone is narrow on that scale. It does make the figure an honest representation of the two geometries. It also strengthens the contrast rather than weakening it.
CHANGES NOT REQUESTED BY THE REFEREE
====================================One further change was made while preparing the revision. We list it for completeness. In Section 2.2 the description of the observable has been corrected to match the processing. The text stated that the corrected wet delay alone was used; the quantity actually analysed, and analysed throughout, is the full wet delay, formed as the sum of the a priori term and the estimated correction. The a priori term varies by about 1 mm over a day, so it carries no power in the fitted band and affects neither the fitted scale nor the fitted amplitude. No number in the paper changes; the description was simply incomplete.
We hope these revisions address the referee's concerns. We thank the referee again for a report that improved the documentation of the processing, the accuracy of the geometrical description, and the legibility and honesty of two of the figures.
Citation: https://doi.org/10.5194/egusphere-2026-4148-AC1
-
AC1: 'Reply on RC1', gael kermarrec, 15 Sep 2026
reply
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 171 | 88 | 38 | 297 | 34 | 36 |
- HTML: 171
- PDF: 88
- XML: 38
- Total: 297
- BibTeX: 34
- EndNote: 36
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
General Comments
The work assess the possibility to study the horizontal hetrogeneity in the atmosphere by using the short term variations in the atmospheric signal delay caused by water vapour. They are estimated from a ground-based GNSS receiver station in Lindenberg, Germany. In total. there are 28 days in June-July 2021 with GNSS estimates which are supported by and compared to observations by microwave radiometers and several other instruments.
Although the temporal resolution of the ZWD estimates is 30 s, I assume that there are constraints in the estimation process (often used is a random walk model with a corresponding constraint specified as a standard deviation). I have no hands on experience from using PRIDE but from what I understand the atmospheric estimates are obtained in a similar way as in the GIPSYX software package, also using PPP. Adjacent values are not independent estimates. The model and the constraint are of fundamental importance and shall be motivated and presented (see, e.g. Young et al., 2024). Furthermore, because the constraints used in the GNSS processing have a significant impact on the resulting variability it is also necessary to determine how significant it is. This is especially true in this case because the variability is the fundamental observable in the whole study. A possible method is to carry out a number of GNSS processings using different constraints. Given that it is only one station and 28 days of observations this is a quick task. If the constraint used so far has a reasonable value, given the expected atmospheric variability at Lindenberg, it will make sense to make solutions using say half of this value and the double value. Thereafter, the main question is if and how much these three time series will affect the results in the rest of the manuscript.
Another general observation is that it is not clear how observations during rain are used, if they are used at all. The microwave radiometer produce more or less useless data during rain. This is not discussed, see below.
When the above question about the impact of constraints used in the GSS processing is resolved the manuscript would benefit from a discussion on the overall context for future applications using this type of analysis. For example, will it have a real-time application for forecasting. The IGS networks is mentioned but not the E-GVAP, run by EUMETNET, including several thousands of GNSS ground stations in Europe (see https://egvap.dmi.dk/).
Young, Z.M., Blewitt, G., and Kreemer, C.: Improved GPS tropospheric path delay estimation using variable
random walk process noise, J. Geod., 98:89, https://doi.org/10.1007/s00190-024-01898-3, 2024.
Specific comments
Lines (L) 2-3: The statement " ... yet its small-scale temporal structure is conventionally treated as noise." does not agree with my experiences from the GNSS meteorology community (although it may be true in the geodetic community. focusing on daily estimates of 3D positions). This text can be deleted in order not to upset readers that are eager to improve the true temporal resolution in GNSS meteorology.
L 70+: In addition to the general comment above more information about the GNSS observations and processing is needed. Are all four GNSS used? What is the elevation cutoff angle? It determines the diameter vs. height of the sensed cone-shaped volume of air. One expects that a higher cutoff angle will increase the variations in the time series of the estimated zenith wet delay. Although the zenith total delay product from IGS does not have the high temporal resolution it would be a quick task to compare it with your time series for validation? For the future it may be a good idea to present the estimates time series from GNSS both those from using different constraints and the official IGS product.
L 77: Here you may want a citation to support the statement.
L82: Why use the Bevis coefficients (which were derived from stations in North America) when the RPG instrument provides the zenith wet delay directly? Most likely the RPG has optimized these algorithms for the Lindenberg site, depending on if it was stated at the time of placing the order.
L 149: You write "... used to filter out windows with non-turbulent rain contamination in the GNSS PPP." I do not understand this. GNSS estimates are more or less insensitive to rain unless it is an extrem rain intensity. In the introduction you refer to GNSS as a technique that " ... deliver, in all weather ...". There is a need for clarification.
L173-175: My experience is that the MWR noise in general exceeds the noise in the estimates based on GNSS data. Of course that depends on the constraint used in the GNSS process (discussed above). In any case it will be informative to state the percentage of data when the MWR noise exceeds the GNSS noise.
Technical Corrections
General: The word "cadence" is often used. It reminds me of music but I assume that the meaning is "temporal resolution". If so, why not use this well established term?
L 10: Define \rho (later in the manuscript it seems to represent both density and correlation coefficient)
L 84: Define all the symbols in Eq. (1)
L 114: "data sets" change to "datasets"
L148: Do not begin a sentence with "18"
Figure 5: Labels are not in the correct positions
Figure 6: The volume of air sensed by the MWR seems to be presented as a cylinder indicating that this part of the atmosphere is in the near field of the antenna. I think the half-power beam width at K-band of these RPG radiometers are around 4° and the aperture is around 20 cm, meaning that the near field ends at around 6 m. So the sensed air volume is a cone with a half-power diameter of some 140 m at the height of 2 km. It may not have an influence on your analyses but it would be nice with a figure that reflects reality. Additionally, also the elevation cutoff angle used in the GNSS processing could be included in the figure.
Reference list: The AMT guidelines are to use the full journal name when the abbreviation is not known. In this case I think most of the journals are found in the list published by "Web of Science: https://wos-help.webofscience.com/WOKRS535R111/help/WOS/A_abrvjt.html