See attached review.

This article applies a three-way error analysis or three-cornered hat method to investigate vertical correlations of ‘errors’ in refractivity profiles, using 3 sources of data that are as different as possible, namely from satellite, in situ, and model-based gridded dataset (respectively: radio occultation (RO), radiosondes (RS), and the ERA5 reanalysis).

I believe the readers will enjoy this article which raises a number of interesting questions - and systematically attempts to bring answers to them. I did not find any flaw in the work or approach, and I only have some questions to help improve the paper (and feed the discussion).

Detailed comments:

1. Given that the authors make an explicit attempt to tie the terminology to established documents like the GUM, other prior relevant publications may deserve to be cited, namely those that already considered the GUM and its applicability to Earth Observation data, e.g., from the FIDUCEO project: Merchant, C. , G. Holl, J. P. D. Mittaz, and E. R. Woolliams. 2019: Radiance Uncertainty Characterisation to Facilitate Climate Data Record Creation. Remote Sensing 11, no. 5: 474. doi:10.3390/rs11050474
2. Section 1.2 mentions “RO reference coordinates” and “RO reference time”: how are they defined in the present work?
3. Is it possible to indicate (or cite the appropriate reference) for the step where radiosondes and ERA5 data are projected into refractivity space?
4. The approach to calculate epsilon C and epsilon X needs to be detailed, and preferably with dedicated equations for clarity. This would also remove the need for forward references to sections 4.2 and 4.3 in section 1.2.
5. A map showing the locations (or a density map) of the 15,997 selected collocations may be a useful information for the readers.
6. Is it possible to comment on the possibility that the two steps of cubic spline interpolation, and of refractivity computations, both applied to ERA5 and radiosonde data, may (each) introduce correlations of uncertainties between these two datasets? Similarly, given that the assimilation of RO data improves the quality of a reanalysis, what are the prospects for some structural correlation between ERA5 forecast (even if, not analyses) and RO data?
7. The numbers of vertical levels in each dataset may need to be introduced in the data section.
8. As the work is using model forecast, whose quality decays as the integration time increases, may one expect a sensitivity of the results to the forecast integration time?
9. Relying on correlations to pick-up a signal exposes one to be sensitive to any transient or structural correlation that may exist in the input data that are correlated, whatever the reason (true signal, artefact of pre-processing, matching bias, …). In the present case, the step of bias removal seems to be limited to a mean profile subtraction, carried out at the scale of each entire dataset (RO, RS, ERA5), is this correct? If so, this would leave, present in the data, all the bias(es) that may exist within each subset of analysis. Would it be possible to consider applying the bias removal in each subset (i.e. at the step when expectation values are computed, modifying slightly equation (7) to introduce the removal of the means), and then display (or report on) how much this changes the results (or not).
10. Do the brackets in the right-hand side of in Equation (7) reflect the actual implementation? (i.e. averages are computed after adding all cross-products?)
11. Can you clarify how the data subsets (the sub-spaces in which expectation values and hence correlations are computed) are defined? This seems to be, at least initially, based solely by considering the vertical dimension, but then later in the paper other dimensions (for computing the correlations and presenting the results) are introduced. This may be done with various sets of subscripts (for the various dimensions: vertical, latitude band, …). In the ideal case where one would have many events, one could consider to compute these error estimates with subsets defined spatially (e.g. 5 deg x 5 deg). The resulting geographic patterns that may be obtained could be of interest.
12. Would it be possible to define early on, i.e., in the methodology section, the ‘raw uncertainty estimates’ mentioned in the results section? Also, the sigma symbol may deserve to be introduced numerically with an equation.
13. In figure 7, does the “0 m” line refer to no filter? If so, such a filter has an infinitely small width (Dirac), but is probably non-zero.
14. The results in figure 7 indicate that as one filters out small-scale variability in both other datasets, the dataset that appears to be most affected in its ‘error’ estimate is ERA5 (this one presents the largest spread, nearly 1%, between no-filter and 1800 m filter, in the lower troposphere). This would be consistent with that dataset containing the least small-scale vertical information, given that equation (7) suggests that for an error estimate to increase, there are two pathways: the cross-products of differences with respect to the two other datasets increase (first two terms in the sum), and/or the cross-product of the differences between the two other datasets decreases (third term, negative sign). The latter may be the mechanism by which removing small-scale information in RS and RO data (thus reducing the differences RO-RS) leads to ERA5 to appear of worse quality (when in fact its quality should be independent of that, but here the method uses the other data as references). Such hand-waving comments (for lack of a better expression) may be tried with a simple toy model. Similarly, the dataset whose ‘error’ estimate is the least affected by filtering the small-scale variability in the two other datasets seems to be the radiosondes, which is also consistent with that data source possibly containing the most small-scale vertical information in the lower troposphere (or is the figure 7(c), truncated at a maximum of 2.0%, showing something else?)
15. Figure 9 shows seemingly slightly different results because in this case one considers smoothing on only one (other) dataset a time. However, one finds consistency. When (only) RS or (only) RO data are filtered (in (a) and (c), respectively), then either one of the two may start to resemble more to ERA5 (but less to the other, i.e. RO or RS, respectively), so, in equation (7), the three terms that make up the total error estimate sees changes of different signs in its components (respectively, for the 3 terms in the right-hand-side of equation (7): decrease of differences ERA5-RS, no change ERA5-RO, and minus an increase of differences RS-RO – the net result is then a decrease of ERA5 estimated ‘errors’ when only RS is filtered). Similarly, this would explain that the ‘error’ estimate of RO increases in (b) when (only) RS data are filtered, making them resemble more ERA5 (respectively: increase of differences RO-RS, unchanged differences RO-ERA5, and minus a decrease of differences RS-ERA5). I note in passing that one missing piece of this puzzle would be to show what is happening to the error estimates of each dataset, when one filters that dataset only, and none of the other two datasets, as the results are not entirely predictable because they involve the sum of two terms moving in opposite directions, e.g., for RO error estimates, if filtering RO data: the differences RO-RS may increase, the differences RO-ERA5 may decrease, the differences ERA5-RS would be unchanged (so the net result is hard to predict -- but such a thought experiment may help shed light on the optimal ‘footprint’ to characterize each dataset).
16. With respect to the expected long-range correlations expected in RS data (but not picked-up as well as expected by the method), this may also be due to the choice of sub-setting considered here, analyzing together day-time and night-time data. The effects of the radiative corrections (leading to consistently positive or consistently negative differences in each profile) may, if not cancel out, possibly be reduced, when considered together. However, redoing the exercise by separating clearly night- and day- ascents (and possibly leaving aside those profiles ‘in between’), may show slightly different results. Such a separation for the RS data would somehow echo the efforts made to separate between rising and setting events for the RO part.
17. Figures 11-12, I fail to see the labels (a) to (f) (either add these or amend the figure caption?).

Paper Summary:

The authors develop and apply a variation of the three-cornered hat (3CH) method of analysis for three independent collocated (to within some limits) refractivity data sets: 1) from Metop-A, B, and COSMIC-1 radio occultation (RO); 2) from RS92 sondes as processed by GRUAN; and 3) from ERA5 forecast fields. The 3CH analysis is generalized to include off-diagonal elements such that certain error correlations are accounted for. Various sources of uncertainty are defined and analyzed. Derived RO uncertainties are reasonably consistent with previous “first principles” error analyses for RO. Altitude and latitude dependent uncertainties are presented, to gain insight into the RO data set and improve assimilation of RO into numerical weather models.

Review Summary:

The paper has the potential to be an original and useful advance in analysis of RO error characteristics. The generalization of 3CH (called “G3CH”) is potentially very valuable. However, the paper suffers from lack of clarity in 1) the definition of certain assumptions, and 2) how certain conclusions are reached. The paper requires major revision to improve and clarify the presentation. Clarifying certain assumptions might alter conclusions of the paper. It is difficult to be certain. In any case, the paper can be greatly improved with an altered presentation. Detailed comments follow.

Detailed Comments:

Line 36: this sentence is imprecise. One of the purposes of the paper is to take into account error correlations. Errors contain a random component. Therefore stating that the random error components are independent appears mis-leading. Random errors can be dependent and correlated. This should be re-phrased.

Line 39: while it is true that the authors focus on vertical error correlation, they have not made a convincing case that error correlation might not arise for other reasons. It seems that the current analysis could proceed at a particular vertical level, in which case it would seem incorrect to assume that all error correlation is from the vertical dimension. More on this point later.

Line 41: in light of the earlier sentences in this paragraph, are we to assume that the ERA5 and RS92 error covariance matrices contain off-diagonal terms only because of vertical error correlation? Please clarify.

Line 58: this sentence is not understood. Vertical footprint for RO profiles will be of order 10 km, which is much less than the distance between RO measurements for the data sets considered here.

Line 63: we suggest the authors add a figure to the paper that defines precisely what is meant by “footprint” for the data sets and for truth. This same figure should clarify the term “observation grid” (Line 69), since the observations are available at random places and times, and not on a grid.

Line 66: assuming that systematic errors are removed, i.e. errors have no bias, is a confusing aspect of this paper. While it is true that the authors remove global means from the data sets, this does not imply the errors as analyzed contain no bias. The reason is that the authors use data subsets in the analysis (e.g. latitude subsets, collocation-distance subsets, etc.) and these subsets may contain bias. An example would be lack of global bias arising because there are equal and opposite biases in the northern and southern hemispheres. The authors need to consider the possibility of biases in subsets of the data. If they take this into account, the analysis can proceed apace.

Line 69: we again recommend a figure be used to carefully define how the truth data set is “distorted” when mapped to the observation grid, and to carefully define what is meant by “representativeness error”.

Line 83: based on the writing so far, “error cross correlations” are error correlations between data sets. It is not immediately obvious to this reviewer how finite footprints of the data sets would lead to such correlations. (Again, the suggested figure might help here). For example, if the data sets are not overlapping in space, why would finite footprints lead to error correlation? Also, please clarify whether the footprints alluded to are horizontal, vertical or both. This question is posed because in several locations of the paper it is implied that vertical correlation is what leads to non-zero off-diagonal covariance, so vertical footprint would be relevant.

Line 120: Please clarify the notation. Do the different epsilon terms (x,y,z) each decompose into components I, R, C, X in equation (1)? We assume that the vector here represents different values along the vertical dimension. That could be stated explicitly.

Lines 127-128: we have remarked earlier how the bias free assumption may not apply. Can the authors verify they have removed bias from all data subsets they have worked on? Subtraction of one global bias will not guarantee there are no biases in subsets of the data. Also, the statement that randomness implies “bias free” is mis-leading. Please modify this statement.

Line 156: we raise again the concern that all data subsets might not be bias free. Has this been confirmed in the analysis?

Line 167: we ask the authors to define the “vertical footprint of truth”. The figure asked for earlier would again help here.

Line 167: why wouldn’t similarity of horizontal footprints also lead to cross-correlated errors?

Line 171: what is meant by “common grid”? What are the spacings of this grid?

Line 180: it would be useful to clarify the mathematical relationship between the uncertainty estimate and the footprint. What assumptions are made to derive this relationship?

Line 199: I believe what is being stated here is that G3CH is incorrectly assigning collocation error to RS92 error. We expect then, that if ERA5 is collocated to RS92 rather than RO, the RO would show the large uncertainty. Is that the case?

Line 218: is there a way to justify this interpretation using a mathematical model and showing it mathematically? Otherwise, it’s difficult to assess the validity of this interpretation.

Line 239: please refer back to the equations where this error covariance is defined. See the earlier comment about how the errors break down into the different components. One could, for example, insert those components into the covariance equations (7) and identify specific outcomes depending on the properties of these error components.

Line 270: see question raised earlier of how “footprint of the truth” is defined.

Line 275: the concept of “physical variability” is introduced here for the first time. How does it relate to the error components I, R, C, X defined earlier? Or is it a new component of error? In general, the statistical properties of the error distributions are not explicitly described (are they gaussian?) except that they are mean zero. If statistical error distribution is not relevant, and any mean-zero error distribution is acceptable, it should be stated.

Line 295: The Rieckh paper uses RO, radiosondes and analyses and forecasts. Please be more explicit why these data sets are not suitable for 3CH analysis, since they appear to be similar to the data sets used in this paper.