The authors aim to examine the effect of assimilating GRACE data, along with other datasets, on state estimates produced by the project for Estimating the Circulation and Climate of the Ocean (ECCO). By way of comparing the results against a reference run that does not include any data constraints the authors claim that the ECCO optimization leads to large adjustments in bottom pressure (pb) fields at monthly and longer timescales largely upon assimilation of GRACE data. Another conclusion drawn is that a mean ocean mass constraint is essential for mitigating large imbalances in freshwater fluxes derived from atmospheric reanalyses (used as prior forcing) and for producing a realistic barystatic sea level curve. Inspecting the residuals, the authors also point to problems with the GRAC data that appear to be inconsistent with other information about the ocean circulation and its variability.

The subject of the paper is important, and ultimately the manuscript should be published. However, some mayor shortcomings must be remedied first. First and foremost, the experimental setup is flawed in that only one assimilation experiment constraint simultaneously by many data sets is being used to pinpoint the influence of a specific data set – the GRACE data - by comparing the results against a reference run in which no data were assimilated at all. Obviously, this cannot work and conclusions drawn are not backed up by the results shown. This holds even more so as the GRACE data kick in at the same time when Argo data become available and so the solution will be impacted in all its aspect by both data sets and all other data as well. Ideally, exactly the same set-up should be run twice with and without GARCE data involved. At a minimum, what the authors need to do here is compare results against a previous optimization run (Release 4 described in detail by Fukumori et al. (2019) which included almost all data as constraints, but not the GRACE data. Moreover, many important details are missing in the paper regarding the approach but also the assimilation experiment itself that need to be added. More detailed comments are provided below. Once all those have been addressed satisfactorily I belief that the paper can become a significant contribution.

Major Comments:

• The authors use a GRACE and GARCE-Follow on data product. But no information is provided as to what they do to this data set by themselves. This relates especially to the question if a trend was removed, with which temporal resolution the data are available and being used as constraint. Also, the source of the uncertainty shown in Fig. 1a is not specified. In Section 2.2 a different GRACE data set is being introduced that apparently is being used for comparison but differs from the data set assimilated. Would be useful to expand on this; e.g., why was the GRD aspects in the assimilated although the model cannot deal with this?
• What is the effective spatial resolution of the GRACE data and how is this incorporated into the cost function? Are coastal regions excluded in the cost function or is the model constrain by the GARCE data also over shallow regions?
• How is the Antarctic ice-shelf melt incorporated into the model. The text refers only to a heat transfer coefficient which is being adjusted. But what happens to the freshwater? Is this being use as freshwater flux as well?
• The handling of spurious changes resulting from the Boussinesq model physics mentioned only in passing in Figure caption 10 needs to be explained in detail in the text as this is an important aspect that needs to be understood to appreciate the results, e.g., with respect to the global barystatic sea level curve.
• The existence of earthquake imprints in the GRACE results is mentioned many time throughout the paper. I suggest discussing this early on once to avoid frequent repetition.

Minor Comments:

• Line 114: suggest referring to later discussions of rms difference to make this text here a bit more quantitative.
• Line 117: eddy variability in GARCE data? What is the effective resolution of GARCE data (see above) and which eddies are resolved in this data set? I guess what is meant here is that GRACE observes a field that also contains the effect of eddies; but that individual mesoscale eddies are not resolved in the data set. In contrast the model does not resolve eddies and does not show the effect of these unresolved eddies.
• Line 118: “ … as discussed in relation to Fig. 1a …” I am not sure where this text points to; at least I could not find it.
• 3 and related text: the expression “cost” is being used here without introducing it. I find this term misleading and suggest using a different word, such as normalized quadratic misfit. Cost specifically points to the cost function that was not even mentioned before.
• Line 139: change “converging” to “converged”
• Related to Fig. 2d: How realistic is the GRACE error provided?
• Line 143: “Peak value of distributions”? I suggest rephrasing the text.
• Line 152: I suggest changing “observations” to GRACE data, as these are also estimates and not observations.
• 4: Details hard to see, specifically those situated in shallow regions; holds also for Fig. 6 and possibly others. Which period is considered here and in other figures?
• Line 164: Seasonal cycle main component? Is this consistent with the figure shown here? What is the percentage of the seasonal cycle relative to the total variance?
• Line 185: Contrary to what is stated I still see a significant seasonal cycle in these principal components 3 - 5.
• 7: you detrended over a longer period; yet the GRACE trend is still prominently present in mode 3!
• Caption of Fig. 7 and others: Change ECCO fields into ECCO pb fields.
• 8: bowls shape probably due to model drift/adjustment. Should be discussed.
• 9: The figure is misleading; largest adjustments happen over all shallow region and ACC, not in the Arctic as stated. All shallow regions will be affected by Argo data. So, what the impact of GRACE is as opposed to other information is entirely unclear. I suspect there in no impact of GARCE at all.
• Lines 253 and Fig. 11: what about the jumps in PC1 and in PC5? I suggest putting in the times of the earthquakes as I suspect that they are the reason for those. What happens in mode 5 after 2017? Impact of GRACE-FO?

Like Professor Stammer, I also feel this paper is worth publishing after some revisions are made. My major concerns (as indicated in the point-by-point comments below) concern the lack of some quantified statistics in the discussion. The discussion throughout could be improved by giving explicit values of differences (e.g., X% of the oceans have a change > 1 cm RMS).

It is a shame that the authors do not have an assimilation run on hand that was done in the same framework without assimilating GRACE data. I understand the reason why (computation costs are expensive), but I do encourage the authors to consider Dr. Stammer’s suggestion on a potential way forward by considering an older run (in a slightly different assimilation framework) and comparing for that specific period around 2002 when GRACE data is introduced. This would provide the best evidence that it’s not Argo profile measurements that cause the change at that period.

With that said, I do find Figure 8 and 9 compelling evidence that it is GRACE that drives the improvement. However, I note the evidence that suggests this is the gap between the end of GRACE and the start of GRACE-FO, but this is not explicitly discussed in the section (see Point 9 below). I encourage the authors to take another look at this particular time-period and perhaps do some additional analysis on it to tease out the GRACE improvements further with the data they have.

Point-by-Point Comments

======================

1. Although the figure captions clarify that global means are removed before computing Pb variance, this should be discussed more in the text as to why this is done – i.e., there is a ± 1-2 cm variation in global mass (barystatic sea level), primarily with seasonal periods but also a trend. By removing monthly global means, the authors are focusing on the dynamic ocean mass variations, not geodetic. This information is glossed over in the introduction and data description but should be discussed.
2. Lines 115-120: “In areas of western boundary currents (e.g., Agulhas retroflection, Gulf Stream extension, Argentine Basin), GRACE fields look more energetic than both ECCOv4r5 and ECCOctrl, consistent with the presence in GRACE of eddy variability that cannot be represented in either ECCOv4r5 or ECCOctrl, as discussed in relation to Figure 1a”. Some more discussion is warranted here. For example, the authors should reference previous works (e.g., Fu’s studies in the Argentine Basin) indicating there are known and large variation is barotropic sea level (Pb) at mesoscales in these regions. While the ECCO resolution can’t capture the magnitude of Pb variance due to the resolution, GRACE does capture more of it. I doubt GRACE is capturing the full magnitude or resolving the mesoscale signals, but it does likely have information there that ECCO is not capturing. It should also be noted that the GRACE data is explicitly down-weighted in these regions (see Figure 1a), which may explain why it has limited influence on ECCO there.
3. Some quantification should be added in the discussion around Figure 2, especially Figure 2d. For example, what are the values of RMSD in certain areas (e.g., the Weddell Sea). The color bars make it difficult to see a difference between 2.5 and 1.5 cm, for instance. It may be that even though the change looks small, it still may be significant. This is shown to some extent in Figure 2d, but the authors need to guide the reader a little more on what the negative vs. positive values mean – e.g., negative means the assimilation pulled it closer to GRACE while positive pulled it away. Also, please add some statistics: what percentage of the ocean improved? By how much? What was the average RMSD increase in smaller areas where there is a larger change.
4. I understand the “Model Cost” concept, but I wonder if expressing it int erms of percent variance reduced/explained (or a similar metric) would be understood by a broader audience. Remember, this is not just being written for modelers or experts on data assimilation.
5. I notice there is no discussion of correlation. My first thought reading this is does the temporal correlation in each grid improve with assimilation? And by how much? I suggest adding some analysis on this.
6. One aspect of the trends in the control run that are not mentioned but are quite obvious to me is the negative trend throughout the Atlantic much of the Southern Ocean and positive trends throughout the Indian and Pacific Oceans north of north of about 60°S. This indicates a relatively large exchange of mass between those basins, enough to cause a change of ± 6 cm over the 20+ year period! This is indicative of a large-scale drift in the control run which is fixed via the assimilation. While there are some mentions of changes in the Southern Ocean, I feel this should be more explicitly discussed. In some respects, it’s a bigger change than some of the more localized differences. That “mode” does seem to be captured in EOF3 (Figure 7). One can especially see the change in the system in that mode when GRACE is introduced. I don’t see this is even discussed in that section! I think it is one of the more obvious impacts starting around 2002!
7. For the annual analysis, I would be interested in seeing a phase analysis. I know the authors use a climatology, so a phase value is not a direct estimate (unlike with a sinusoid fit). However, they could look at the month of the peak value and see if there is any evidence of a phase shift from the control to the assimilation run. Looking at the EOF analysis, I think there will be and it would be an interesting thing to show.
8. For Figure 8. Please add some statistics of the actual values with GRACE obs and without. It will be a small change, maybe 1 cm RMS, but worth noting with a quantified metric. Same for Figure 9 – What’s the percentage of the ocean where changes are > 0.5 cm? Between 0 and 0.5 cm? etc.
9. One of the important periods for assessing the impact of GRACE is for the year when GRACE observations ended and before GRACE-FO came on-line. No other change in the observing system occurred then (unlike before 2002 with the introduction of Argo profiles into the assimilation as well). It is clear from Figure 8 that the updates to the control run Pb return to the pre-2002 levels, then move back up in 2018 when GRACE-FO is used. I see no discussion of this at all, but to me, it is the clearest indicator in the analysis that GRACE is responsible and not another data set.