GOFS16: an operational global ocean analysis and forecasting system at eddy-resolving resolution

Masina, Simona; Cipollone, Andrea; Iovino, Doroteaciro; Ciliberti, Stefania; Lecci, Rita; Cretí, Sergio; Lyubartsev, Vladyslav; Coppini, Giovanni; Clementi, Emanuela

doi:10.5194/egusphere-2026-887

Preprints

https://doi.org/10.5194/egusphere-2026-887

Preprints

09 Apr 2026

| 09 Apr 2026

GOFS16: an operational global ocean analysis and forecasting system at eddy-resolving resolution

Simona Masina, Andrea Cipollone, Doroteaciro Iovino, Stefania Ciliberti, Rita Lecci, Sergio Cretí, Vladyslav Lyubartsev, Giovanni Coppini, and Emanuela Clementi

Abstract. The Global Ocean Forecast System GOFS16 is an operational ocean analysis and forecast system that has been running daily at the Euro-Mediterranean Center on Climate Change since early 2017. GOFS16 produces 6-day forecasts of the state of the global ocean and sea ice (three-dimensional ocean temperature, salinity, and currents, as well as sea level, sea ice thickness, concentration and drift) with a system based on the NEMO model configured in a global eddying ocean at 1/16° horizontal resolution and 98 vertical levels. To compute the initial conditions for the forecasting, in situ observations of temperature and salinity, altimeter data of sea level anomaly and satellite sea surface temperature fields are jointly assimilated each day over a 1-day observation window with a 3DVar scheme, OceanVar, adapted to handle the global high-resolution grid. This paper introduces the first version of the GOFS16 system and describes its components and the validation procedure, which routinely produces statistics of the system skill. The scientific assessment is presented for the period January 2022–December 2023, at global and regional scales. Results indicate that GOFS16 performs within the expected range of skill for current global systems. However, the impact of the eddy-resolving resolution is hidden either by system weaknesses or traditional verification metrics not adequate at such resolution.

Received: 14 Feb 2026 – Discussion started: 09 Apr 2026

Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims made in the text, published maps, institutional affiliations, or any other geographical representation in this paper. While Copernicus Publications makes every effort to include appropriate place names, the final responsibility lies with the authors. Views expressed in the text are those of the authors and do not necessarily reflect the views of the publisher.

Download & links

Simona Masina, Andrea Cipollone, Doroteaciro Iovino, Stefania Ciliberti, Rita Lecci, Sergio Cretí, Vladyslav Lyubartsev, Giovanni Coppini, and Emanuela Clementi

Status: final response (author comments only)

RC1:
'Comment on egusphere-2026-887', Anonymous Referee #1, 08 May 2026
The paper presents the GOFS16 short-range ocean forecasting system run operationally at CMCC. The system is run at 1/16° resolution using the NEMO model and OceanVar 3DVar data assimilation scheme. The model, data assimilation and observation components of the system are described, as well as operational implementation details. Assessment of the operational system is presented for the 2-year period January 2022 to December 2023 using bias, RMSE and anomaly correlations of the analyses and forecasts compared against observations. A section on the impact of a tropical cyclone on the ocean state in GOFS16 is also included.
The manuscript is generally clearly written with a useful amount of detail in the system description, though some additional details could be included to strengthen the paper, as suggested below. The assessment of results follows standard procedures used in the ocean forecasting community and gives a clear picture of the quality of the GOFS16 analysis and forecasts. However, the assessment results are often compared with results from other systems published more than 10 years ago, and do not necessarily give a fair impression of the quality of GOFS16 compared to other operational systems of today. Conclusions about the impact of the high resolution of GOS16 also seem rather strong given that there isn’t a direct comparison with a lower resolution version of the same system, and the fact that so many aspects of the system (apart from the resolution) can affect the forecast skill.
Major comments
In the abstract, it is stated that “Results indicate that GOFS16 performs within the expected range of skill for current global systems”. Instead, since there isn’t a direct comparison with the other current global systems included in the paper, it would be better to give some details of the actual skill of the GOFS16 system based on the assessments presented in the paper.

The description of the data assimilation set up would benefit from some clarifications and more information. Specific comments are made in the minor comments section below.

In numerous places, (e.g. L345, L363, L381, L401, L425, L433, L466) results from GOFS16 are compared with those published over 10 years ago in Ryan et al., 2015 and one single recent system description in Liu et al., 2021. The other systems referred to in Ryan et al., 2015 have been developed since that publication with eddy resolving configurations implemented by many centres which is not clear from your descriptions. More recent publications with skill information from some of the other systems is available, e.g. Mignac et al., 2025; Lellouche et al., 2018. Apart from the issue with the old reference and updates to those other systems, it is not a like-with-like comparison of the statistics in the sense that results in the older references are for a different period with different observation sampling, different quality atmospheric forcing etc.

Conclusions are drawn on L616 about the impact of model resolution on the results (and are implied in other places such as L329, L345, L432). However, there is no comparison of GOFS16 to a lower resolution version of the same system. The comparisons to other, much older, lower resolution systems has major flaws, particularly in that those older systems were run for much older periods with different observing systems and different atmospheric forcings available.

It is interesting to see the response of the GOFS16 system to the tropical cyclone in Section 4.4., but this section would be much stronger if there were comparisons of the model forecasts to observations, particularly if there were Argo profiles in the region of the cold wake to show how well the GOFS16 system was representing the true ocean response.

Minor comments
L44-45: “…the eddying 1/16° configuration (GLOB16) improves the representation of the boundary currents…” – improves compared to what? Do you mean compared to ¼ and 1-degree configurations as described in the Iovino 2023 paper.

L47. Maybe change “prone” to “designed”?

L56. You could state here when the system was made operational.

L57. It isn’t clear from this sentence which fields are available hourly or daily.

L58. When I visited the website (https://gofs.cmcc.it) the display of the outputs was not working and it said “invalid request” over the map.

L85. 3.4 is a very old version of NEMO – are there plans to update it?

L85. “…that integrates numerically primitive equations …” this could be reworded.

L112-113. Why is the temporal frequency of the surface forcing so low, especially for radiative fields, requiring then an imposed diurnal cycle described on L120?

L120-121. If the mean sea level can drift in the model, how is this dealt with in the data assimilation?

L150-151. You talk here about non-linear observation operators. Is this an FGAT scheme then in which the model is used to generate the innovations at the correct time? If not, then which non-linearities are included in H? What time is the background field used in the calculation of the innovations valid for?

L151. “construct” -> “constructing”

L154. Does the vertical operator define the magnitude of the background error variances? Please describe the magnitude and spatial structure of the background error variances in more detail and perhaps show a figure of them.

L158. Why are no SSH increments produced? How does the model react to T/S increments without balancing increments to SSH or velocity?

L161. Is the recursive filter performed in 1Dx1D or 2D? How are the normalisation factors calculated to ensure unit correlations at zero separation?

L162. Are the length-scales are coming from a ½ degree resolution model? How would that affect the results of DA using 1/12 degree model?

Table 1. For SLA you list both L3 and L4 products – why is that?

L170-176. I found it hard to understand the reasoning behind the SLA bias correction scheme. Are the innovation bias and other terms calculated at each model grid point or some lower resolution grid? Is the resulting bias estimate d_overbar then removed from the SLA obs on the current day? Why is the bias adjusted in the way it is depending on the observation and model variability? If the obs or model variability is large then the bias is reduced to zero, but in that case do you then assimilate the raw (un-bias-corrected) observations?

L176. Suggest removing “instead”.

L178. How is the SST diurnal signal treated when calculating the innovations? Presumably you are not aiming to correct the diurnal signal since the increments are added over a whole day, so how do you avoid diurnal cycle errors from contaminating the innovations?

L183. Why do you thin in this way rather than averaging to create superobs?

L187-189. This description of the error variance would fit better with the earlier error covariance descriptions, e.g. around L162.

L189. What do you mean here by the horizontal length-scales being scaled to the high-resolution grid?

L191-192. Why do you need to include a nudging scheme for SST during the analysis when you are assimilating SST observations? Is the SST nudging applied in the forecast?

L193. Is the SSS nudging scheme included in the forecast?

L201. “…reduce the artificial spatial correlation among observations that are inevitably introduced at the length scale of grid spacing”. Including observations on the model grid at this high resolution will surely give you observations which do contain spatially correlated errors, and more so than if you thinned to lower resolution model grids.

L202. Typo “explicitly”.

L213. By “inspected” do you mean they are manually inspected as part of the operational production/dissemination? Is this done every day including weekends?

Fig 1. Since every day the pattern is the same, perhaps you could zoom in on only two or three of the days so that it is easier to see the details of each cycle.

L221-222. I don’t agree that this procedure is analogous to nowcasting in the atmosphere. In nowcasting, high resolution observations are used with the latest analysis to provide a very short-range estimate, whereas here there are no observations included.

L231. Why did you choose to use GFS atmospheric forcing rather than ECMWF? The spatial and temporal resolution of the forcing is very low compared to the ocean model resolution.

L238. The EN4 objective analysis is not available close to real time so with what delay is the SSS field available and how do you deal with the time offset?

L247. Is the first guess field a daily average, instantaneous in the middle of the day or something else? How is the interpolation of the model values done in the vertical for profile observations?

L252. Do you mean that “the increments are added to the model fields with an …” and not added to the restart fields?

L255-258. You discuss in situ observations being available close to real time, but satellite SST observations are also available within a few hours.

L259. Which product data are hourly and which daily? Or do you provide all data at the full 3D resolution at hourly and daily frequencies?

L262. “0.33 wall clock hours”. Do you mean that you run the two 1-day NEMO runs (background trajectory and IAU) at 1/16 degree resolution and the OceanVar step all within 20 minutes elapsed time?

L274. Is the assessment against observations carried out using the model best estimate fields? In that case, please clarify that the observations have already been assimilated and are not independent.

L319. Did any of the bug fixes or continuous improvements have a large enough impact to describe what they are here?

L325. What is the cold bias in SST due to?

L331. What is happening during 2023 in Fig. 3a? There's a different behaviour in the bias during May-Aug 2023 compared to 2022 with a positive rather than negative bias evolution, particularly in the 5-day forecast. Also, why is there a gap in Oct 2023?

L340. The bias values are negative, not positive as stated in the text.

L352. Any idea why there are degradations in the N. Pacific?

L358. Why is there a bias in the SLA differences, despite the bias removal process? Could this be related to the drift in the global mean SSH you mentioned earlier?

L366. “…with barely any differences among them.” Not sure what you mean here.

Fig. 7. The climatology SST RMSEs on this figure are very different in magnitude to the ones in Fig. 6 of Blockley et al., 2014, especially for the N. Pacific where you have a value of over 3 degrees whereas Blockley et al. have it at about 1.2 degrees. Can you explain this large discrepancy?

L409. Why do you do the assessment down to 1500 m depth and not the lowest depth of Argo floats at 2000 m?

L413. What happens to the RMSE in T and S at the start of 2023 where there is a discontinuity in the RMSE time-series plots of Fig. 9, including also a jump in the salinity bias?

Fig 9 caption. Part b is for T profiles not SST. Also a typo “nowvast”.

L433-435. You claim here that it is the high vertical resolution is giving improvements, but the maximum of the temperature profile RMSE is not any lower than more recently published results from other systems with fewer vertical levels.

L435-436 and also L454. Do you have evidence that it is the vertical mixing in the model causing the large degradation in temperature and salinity? Other explanations could include an unbalanced initial state or overfitting of the profile observations in the analysis.

L446. I think you mean Fig. 10b?

L474. Argo profiles are now provided close to real time so it would be useful to provide evidence if you think this is not the case.

L496. The observations of drifter-derived velocities from CMEMS do have some temporal filtering which allow a fairer comparison with models.

L497. The NEMO model should include near-inertial oscillations (see for example Waters et al., 2024; https://doi.org/10.3389/fmars.2024.1383522).

L506. “Data are vertically averaged”. Do you mean the GOFS16 data?

L508. Is the correlation an anomaly correlation as stated in the caption of Fig. 12? If so, what climatology is used to calculate the anomalies?

Fig. 13. Is this the GOFS16 best estimate SST field? Perhaps you could also include the field from before the TC (e.g. 30^th August) to show the change in SST more clearly.

L660. You could also mention here more sophisticated data assimilation algorithms such as 4DVar and hybrid ensemble/variational methods which could use the limited observations in a better way.

L696. Multi-scale methods don’t have to be done in two steps, see e.g. Mirouze et al., 2016 (https://doi.org/10.3402/tellusa.v68.29744).
Citation: https://doi.org/10.5194/egusphere-2026-887-RC1
RC2:
'Comment on egusphere-2026-887', Anonymous Referee #2, 20 May 2026
The authors present a global implementation of a 3DVAR assimilation scheme with a high resolution oceanographic model. I think that the manuscript is well written and the evaluation of the impact of data assimilation within the new global operational system is very detailed and informative. I do not fully agree with some conclusions of the study and I have several small comments that may be addressed by the authors before the final publication.

Comment:
The authors discuss the modest improvement of the data assimilation performance with four times higher model resolution. They suspect that this is mainly due to the inadequate forecast evaluation methodology developed for low resolution models. They also suspect that arbitrary reducing the background horizontal covariance scales in the data assimilation scheme might negatively impact its performance. I think that there may be other equally or even more important reasons for lower than expected system’s performance. An important reason might be that this data assimilation implementation produces only temperature and salinity increments that are dynamically unbalanced with respect to sea level and velocity increments. The IUA methodology is used to add increments to the initial state, but this may result in excessively smoothing and eventually losing a significant part of information introduced by data assimilation. This may be a less important issue with low resolution models characterised by simple dynamics, but in high resolution models the dynamical consistency of increments like those in Oddo et al (2026) (doi:10.5194/gmd-19-423-2026) might be a crucial ingredient for maintaining the impact of data assimilation in short term forecasts.

Minor comments:
Line 111: Does “fields” mean “fields of atmospheric variables”

Line 112: What are “turbulent variables”?

Line 113: What are “daily fields”?

Line 116: I guess that it should be “concentration” instead of “value”.

Line 166: I could not find information on how mean dynamic topography is included in the calculation of SLA misfits.

Lines 170-176: Is there a reference on the debias methodology? On line 171 “d” is defined as innovation, and on line 172 as bias. How are model and observational standard deviations (sigma_n and sigma_o) calculated? Are they prescribed arbitrary or updated dynamically? Are they estimated along tracks, regions or have a high spatial resolution? What is exactly the physical reasoning behind this scheme? I think that this all is important, because “unbiased” SLA observations that might be dynamically inconsistent are used to propagate temperature and salinity increments over deep water columns, and debiasing might significantly influence the forecasting performance.

Line 176: It should be “unbiased” instead of “unbias”.

Line 178: Does “the shallowest” mean “the top”?

Line 179: What does the “lower bound for the error” mean? Are misfits with magnitude larger than this excluded? This value is much lower than the average RMS for SST misfits shown in Fig. 3a.

Line 185: How is an observation flagged bad? It could be expected that misfits are large at the thermocline and might be flagged as bad. All observations below are excluded? On the other hand, in Fig. 10 the number of assimilated observations seems to increase with depth.

Line 188: I guess that both background and observational error covariances change significantly with four times higher resolution. Should they be recalculated for the new model?

Line 189: I do not understand why horizontal correlation scales of model errors should linearly scale with the model resolution.

Lines 190-194: Is the SST and salinity nudging in conflict with assimilation of other observations, in particular, SST observations described in one of the previous paragraphs?

Line 198: The sentence is not clear. What shares the same grid?

Line 218: It should be “increments” instead of “corrections”. Can you provide a reference for IAU?

Lines 435-439: Is this discussion related to Fig 10? This result may also indicate the presence of dynamical inconsistency between initial fields of temperature, salinity, sea level and velocity.

Line 446: I think it should be Fig 10b instead of 8b. Arbitrary filtering out observations with large misfits may influence the system evaluation statistics when comparing with other data assimilation systems.

Line 670: What does it mean to correct scales? What are “coarse observations”? Is one scale for the spatial bias and the other is for the unbiased errors? How can these coarse and fine scales impact model dynamics?
Citation: https://doi.org/10.5194/egusphere-2026-887-RC2

Simona Masina, Andrea Cipollone, Doroteaciro Iovino, Stefania Ciliberti, Rita Lecci, Sergio Cretí, Vladyslav Lyubartsev, Giovanni Coppini, and Emanuela Clementi

Viewed

Total article views: 447 (including HTML, PDF, and XML)

HTML	PDF	XML	Total	BibTeX	EndNote
289	135	23	447	17	21

HTML: 289
PDF: 135
XML: 23
Total: 447
BibTeX: 17
EndNote: 21

Views and downloads (calculated since 09 Apr 2026)

Month	HTML	PDF	XML	Total
Apr 2026	190	82	15	287
May 2026	93	50	5	148
Jun 2026	6	3	3	12

Cumulative views and downloads (calculated since 09 Apr 2026)

Month	HTML	PDF	XML	Total
Apr 2026	190	82	15	287
May 2026	93	50	5	148
Jun 2026	6	3	3	12

Viewed (geographical distribution)

Total article views: 445 (including HTML, PDF, and XML) Thereof 445 with geography defined and 0 with unknown origin.

Country	#	Views	%

Latest update: 08 Jun 2026

Short summary

The paper presents GOFS16, an eddy-resolving global operational ocean and sea ice forecasting system which provides 6-day forecasts of three-dimensional temperature, salinity, currents, sea level, and sea ice properties. The system assimilates satellite and in situ observations using a 3D variational data assimilation scheme. Validation is conducted routinely using global and regional metrics. Results indicate that GOFS16 performs within the expected range of skill for current global systems.


Total:	0
HTML:	0
PDF:	0
XML:	0