Validation of an ocean shelf model for the prediction of mixed-layer properties in the Mediterranean Sea west of Sardinia

The Regional Ocean Modeling System (ROMS) has been employed to explore the sensitivity of the forecast skill of mixed-layer properties to initial conditions, boundary conditions, and vertical mixing parameterisations. The initial and lateral boundary conditions were provided by the Mediterranean Forecasting System (MFS) or by the MERCATOR global ocean circulation model via one-way nesting; the initial conditions were additionally updated through the assimilation of observations. Nowcasts and forecasts from the weather forecast models COSMO-ME and COSMO-IT, partly melded with observations, served as surface boundary conditions. The vertical mixing was parameterised by the GLS (generic length scale) scheme (Umlauf and Burchard, 2003) in four different set-ups. All ROMS forecasts were validated against the observations which were taken during the REP14-MED survey to the west of Sardinia. Nesting ROMS in MERCATOR and updating the initial conditions through data assimilation provided the best agreement of the predicted mixed-layer properties with the time series from a moored thermistor chain. Further improvement was obtained by the usage of COSMO-ME atmospheric forcing, which was melded with real observations, and by the application of the k-ω vertical mixing scheme with increased vertical eddy diffusivity. The predicted temporal variability of the mixed-layer temperature was reasonably well correlated with the observed variability, while the modelled variability of the mixed-layer depth exhibited only agreement with the observations near the diurnal frequency peak. For the forecasted horizontal variability, reasonable agreement was found with observations from a ScanFish section, but only for the mesoscale wave number band; the observed submesoscale variability was not reproduced by ROMS.


Introduction
In ocean acoustics research, the diagnostics and prediction of selected mixed-layer properties, such as the mixed-layer depth and the mixed-layer temperature, are of primary interest because they have a profound impact on the propagation of sound in the ocean.In this article, a high-resolution ocean circulation numerical model is presented which provides nowcasts and forecasts of these properties.The objectives are (i) to evaluate the sensitivity of the properties to different set-ups of the initial conditions, lateral boundary conditions, atmospheric forcing patterns, vertical grid, and vertical mixing parameterisations and (ii) to find a set-up which reproduces and best predicts the depth and the temperature of the mixed layer and the associated spatio-temporal variabilities obtained from observations.By definition, temperature and salinity are constant in the mixed layer, and the sound speed increases slightly with depth due to the pressure effect (Dietrich et al., 1975).Therefore, sound rays in the mixed layer are refracted upwards and reflected at the sea surface.Hence, the mixed layer acts as a surface duct (Katsnelson et al., 2012).On the other hand, at a depth greater than the mixed-layer depth and because of the decreasing temperature, the rays are refracted in the other direction, i.e. towards greater depths.Consequently, in terms of passive acoustic monitoring, if a sound source is within the mixed layer, the sound cannot be "heard" at depths greater than the mixed-layer depth.If the sound source is located below the mixed layer, it cannot be heard in the mixed layer.The equivalent is true for the location of objects by active sonar: if the sonar is within the mixed layer, the acoustic signal can hardly reach an object at a greater depth, and vice versa.This is, of course, an idealised model based on ray Published by Copernicus Publications on behalf of the European Geosciences Union.theory which does not take account of the non-linear and frequency-dependent effects, but it clearly emphasises that knowledge about the depth of the mixed layer is mandatory for the planning and conduction of acoustic experiments.
The sound speed c in seawater is a function of temperature T , salinity S, and pressure p.Hence, small changes dc in the sound speed can be described by the total differential dc = ∂c ∂T S 0 ,p 0 dT + ∂c ∂S T 0 ,p 0 dS + ∂c ∂p T 0 ,S 0 dp, where the subscripts T 0 , S 0 , and p 0 indicate that T , S, and p, respectively, are held constant during the execution of the partial differential.For the mid-latitudes and close to the sea surface (T 0 = 15 • C, S 0 = 35, p 0 = 0 dbar), the partial differentials in Eq. ( 1) yield ∂c/∂T ≈ 3.2 m s −1 • C −1 and ∂c/∂S ≈ 1.2 m s −1 , which means that the fractional change in the sound speed with temperature is about 3 times larger than the change with salinity (Chen and Millero, 1977).Moreover, as typical spatio-temporal variations in temperature are O(10 • C), but those of salinity are only O(1) at best, the first two terms in Eq. (1) yield 31.2 and 1.2 m s −1 .Hence, changes in the sound speed are largely controlled by changes in the temperature, and the impact of salinity variations in the mixed layer can confidently be ignored for the calculation of the sound speed.However, one may note that this is only true for the open ocean.In coastal areas, estuaries, and in polar regions, the salinity variations are frequently larger and the concurrent variations in temperature smaller.
Besides the temperature and the depth of the mixed layer, the temporal and horizontal variability of these two quantities require special attention (Pace and Jensen, 2002).The temporal variability at a fixed location is affected by temporal changes in the following: air-sea fluxes in momentum, heat, and fresh water; sea state conditions and internal waves; horizontal advection; vertical motion; and optical properties of the seawater.
The horizontal variability is due to spatial differences of the same quantities and, in addition, to the presence of mesoscale and sub-mesoscale features like fronts, meanders, eddies, and filaments (e.g.Medwin and Clay, 1998).Both the temporal and the horizontal variability impact the sound speed and the underwater sound propagation.
The objective of this article is to find a model set-up which predicts in the best possible way the mixed-layer properties and their spatio-temporal variabilities.While for the temporal variabilities the main focus of attention is directed at timescales between O(1 h) and O(10 days), the intention is to resolve the horizontal variabilities on scales of 10 km and below.This requires a circulation model with a built-in vertical mixing scheme that accurately reproduces the diurnal cycle.A state-of-the-art scheme has recently been published by Ling et al. (2015).It is an enhancement of the turbulence closure model of Noh et al. (2011), which is similar to Mellor and Yamada (1982) but additionally takes into account the effects of wave breaking and Langmuir circulation.Ling et al. (2015) developed new numerical techniques and improved the schemes for the physics in Noh's model, which amongst others intensified the diurnal amplitude of the simulated sea surface temperature.Noh et al. (2016) incorporated Noh's model into a global ocean general circulation model, and they could show that the new mixing scheme helped to correct too-high mixed-layer temperatures and too-shallow mixed-layer depths in the high-latitude ocean.A similar approach was pursued by a series of papers by Bernie et al. (2005Bernie et al. ( , 2007Bernie et al. ( , 2008)).In the first publication, a one-dimensional mixed-layer model was developed based on the K-profile parameterisation of Large et al. (1994).The model was forced with observed fluxes from a mooring in the tropical Pacific Ocean, and it qualitatively reproduced the observed diurnal variability in the sea surface temperature over a period of about 4 months.However, most of the time, the modelled temperature was higher than the observed one by up to 1 • C. Bernie et al. (2007) implemented the turbulent kinetic energy scheme of Blanke and Delecluse (1993) in an ocean circulation model, and this circulation model was coupled with an atmospheric circulation model (Bernie et al., 2008).The major outcome of the latter publication was that the inclusion of the diurnal cycle leads to a tropic-wide increase in the mean sea surface temperature, and, in addition, the authors could demonstrate that the modelled diurnal cycle was modulated by intraseasonal variations.The vertical mixing in all papers mentioned above was accomplished by turbulence closure models.By contrast, Gentemann et al. (2009) improved the parameterisation of the absorption of solar radiation in the diurnal heating bulk model of Fairall et al. (1996a).This change, combined with a reduction in accumulated heat and momentum, increased the model's responsiveness to changes in the surface heat flux and surface stress.Amongst others, the improved model predicted the vertical temperature profile within the diurnal thermocline, increased warming at low wind speeds, and decreased warming at high wind speeds.
The experimental area addressed in this study is situated in the Mediterranean Sea to the west of Sardinia (see Figs. 3,4,5,and 7).Here, the observational data from a June 2014 oceanographic survey are used to drive the aforementioned ocean circulation model and to validate the model results (Onken et al., 2014(Onken et al., , 2016)).
The reader may note that within this article all dates refer to the year 2014 and all times are in UTC (Coordinated Universal Time ) except where otherwise stated.

Observational data
The observational data originate from the REP14-MED experiment, which took place in the eastern Sardo-Balearic Sea in the period of 6-25 June 2014.The collection of in situ data was accomplished by the NATO Research Vessel Alliance, the Research Vessel Planet of the German Ministry of Defence, a fleet of underwater gliders, surface drifters, one subsurface float, and six oceanographic moorings.Throughout this article, however, the author will refer only to the data of one mooring, denoted as M1, and to the CTD (conductivitytemperature-depth) data collected by the survey vessels and the gliders.For more details of the observations, see Onken et al. (2016).

Mooring M1
M1 (Fig. 1) was launched on 8 June at 07:14 at 8 • 12.98 E, 39 • 30.80 N (Fig. 3) and recovered on 20 June at 13:55.The water depth at the launch position was ≈ 150 m.M1 consisted of the central mooring M1 CTR and a sidewaysextending appendix M1 APP floating at the sea surface.M1 CTR was equipped with an upward-looking ADCP (acoustic doppler current profiler) mounted at a nominal depth of 100 m below the sea surface, a CTD probe at 1 m of depth, and a meteorological buoy at the surface.The appendix was connected by a 50 m long rope to M1 CTR and extended to about 40 m in the vertical direction.Forty Starmon thermistors (Star-Oddi, Gardabaer, Iceland) were mounted along the vertical cable to record the temperature with high vertical resolution.In addition, four RBR data loggers (RBR, Ottawa, Canada) determined the actual depth of the Starmons.The nominal and actual vertical positions and the recorded parameters of the sensors are summarised in Table 1.
The Starmons recorded time t (10 s) and temperature T in intervals of 10 s.The RBRs additionally recorded pressure p.The depth z was calculated internally by the RBR software.As the Starmons did not have a pressure gauge, their actual vertical position was evaluated thereafter from the depth records of the RBRs because the positions of the Starmons relative to the RBRs was known.This was accomplished in the following way: the Starmon records at 3.0 and 7.0 m of depth were rejected because the recorders failed (sensor positions 7 and 15 in Table 1); all records before 8 June at 07:18 and after 20 June at 13:30 were clipped and then sub-sampled in 5 min intervals; at each instant, a second-order polynomial fit was applied to the actual depth of the RBRs versus their nominal depth; and the actual depths of the Starmon recorders were calculated with the polynomial using the previously calculated coefficients.
This procedure was advisable in order to correct for potential depth changes of the Starmons due to the actions of waves, internal waves, horizontal advection, and vertical shear.However, it turned out that these corrections were rather small: right at the sea surface at sensor position 1, the actual depth of the Starmon (Table 1) varied between 0 and 1.85 m, and the depth of the sensors at position 41 varied between 39.89 and 41.81 m.Hence, the applied corrections were around ±1 m.The time series of the recorded temperature and the vertical temperature gradient are shown in Fig. 2, several properties of which may be challenging to reproduce for the circulation model.
-The near-surface temperature varied between about 22 • C and more than 24 • C. At the beginning of the recording period, it was around 23 -The mixed-layer depth may be defined approximately by means of the depth at which the vertical temperature gradient is maximal.Between about 15 and 20 June, there is clear evidence for such a signal: the mixedlayer depth varied between about 4 m on 15 June and about 13 m on 18 and 19 June.However, between 8 and 14 June, the signal is rather indistinct: the nighttime mixed-layer depth ranged from about 2 m on 9 June to about 6 m on 10 June, but during daylight hours the maximum gradient was sometimes found right at the surface.Thus, a mixed layer in the "classical" sense did not exist.
-There are clear signals for temporal variability for both temperature and the depth of the mixed layer.The variability occurred on all scales between about 2 weeks and the Nyquist period of 10 min (twice the sampling interval; see above).2.2 Data collected by lowered CTD, gliders, and towed measuring systems On both survey vessels, casts with lowered CTD probes were conducted during the entire survey, but only the casts taken during the 7-11 June period were used here (for a more detailed description of the probes and their calibration, see Onken et al., 2016).These casts belonged to the initialisation survey, the purpose of which was to provide realistic temperature and salinity data for the initialisation of the numerical models.In total, 108 casts were taken on a regular horizontal grid with a mesh size of ≈ 10 km (Fig. 3), resolving the internal Rossby radius of deformation, the first mode of which is around 13 km (Grilli and Pinardi, 1998).Eleven gliders (for their payloads, see Onken et al., 2016) were deployed on 8 and 9 June at the positions marked "L" in Fig. 4 and directed to their nominal zonal tracks.The scheduled tracks were arranged parallel to the CTD sections (Fig. 3) but offset by 5 km in the meridional direction.For the validation of the model forecasts, Alliance conducted a survey during 21-24 June with a ScanFish (EIVA, Skanderborg, Denmark).The tracks are shown in Fig. 5.
3 The circulation model

ROMS
The employed numerical ocean circulation model is ROMS, the Regional Ocean Modeling System.ROMS is a hydrostatic, free-surface, primitive equations ocean model, the algorithms of which are described in detail in Shchepetkin andMcWilliams (2003, 2005).In the vertical, the primitive equations are discretised over variable topography using stretched terrain-following coordinates, or so-called s coordinates (Song and Haidvogel, 1994).In the version used for this article, spherical coordinates on a staggered Arakawa C grid are applied in the horizontal, and the horizontal mixing of the momentum and the tracers is along isopycnic surfaces.
The vertical mixing is parameterised by means of the GLS (generic length scale) scheme (Umlauf and Burchard, 2003) using the stability function of Kantha and Clayson (1994).
The air-sea interaction boundary layer in ROMS is based on the bulk parameterisation of Fairall et al. (1996b).

The model domain and discretisation
The model domain is situated to the west of Sardinia and it is identical to the area shown in Fig. 3.The west and east boundaries are at 6 • 30.5 and 8 • 35.5 E, while in the south and north the domain is limited by the 38 • 36.4 and 40 • 59.6 N latitude circles, respectively.In the east-west direction, the domain is separated into 120 grid cells and in the southnorth direction into 178 cells, which yields an average grid spacing of x ≈ y ≈ 1500 m in the zonal and meridional direction, respectively.A comparison with Fig. 3 reveals that the domain boundaries are kept away from the observations by about 30 arcmin; this was intended in order to mitigate a deterioration of the model solution at the observational sites due to false advection from the open boundaries.Bathymetry data from the General Bathymetric Chart of the Oceans (GEBCO) with a spatial resolution of 1 arcmin were provided by the British Oceanographic Data Centre (BODC), and coastline data were obtained from NOAA (National Oceanic and Atmospheric Administration).In order to avoid the crowding of the s coordinates in shallow water regions, the bathymetry was clipped at 20 m, which is the minimum allowed water depth.For the smoothing of the bathymetry, a second-order Shapiro filter was applied.After smoothing, the so-called rx0 parameter resulted in 0.31, which is about 50 % higher than the maximum value  (2000).However, rx0 is still lower than 0.4 as suggested in the ROMS forum (https://www.myroms.org/forum).
In the vertical direction, the model domain was separated into 70 s layers, the position of which is controlled by three parameters (θ s , θ b , h c ) and two functions, V tr and V str .Here, V tr is the transformation equation, V str is the vertical stretching function, θ s and θ b are the surface and bottom control parameters, and h c is the critical depth controlling the stretching (for more details, see https://www.myroms.org/wiki/).For all ROMS runs shown below, the functions and parameters were selected as V tr = 2, V str = 1, θ s = 5, and θ b = 0.4, while h c was kept a variable.

Initialisation
ROMS was initialised from nowcasts of the coarser Mediterranean Forecasting System (MFS; Dobrowsky et al., 2009;Tonani et al., 2014) or the MERCATOR global ocean circulation model (Drévillon et al., 2008).In either case, the downscaling from the parent to the child was accomplished first Ocean Sci., 13, 235-257, 2017 by the linear horizontal interpolation of the prognostic fields on the ROMS grid.As the maximum horizontal resolution of MFS was close to 7 km (1/16 • ) and that of MERCATOR was 9.25 km (1/12 • ), the scale factors were around 4.7 and 6.2, respectively.Thereafter, all fields were interpolated vertically from the horizontal depth levels to the s coordinates.A special issue was the alignment of the land masks of the parent and the child; if any wet grid cell of the child was covered by a dry grid cell of the parent, a smooth transition of all variables was created by taking the average of the surrounding parent cells.However, as this may lead to a violation of continuity by non-zero horizontal velocities normal to the land mask, all horizontal velocities next to the ROMS land mask were set to zero.

Lateral boundary conditions and nesting
The ROMS code includes various methods for the treatment of open boundaries.After extensive sensitivity studies, it was found that the following algorithms served best for the posed problem: for the sea surface elevation, the Chapman condition was selected (Chapman, 1985), and for all other quantities (i.e.barotropic and baroclinic momentum, turbulent kinetic energy, temperature, and salinity), the mixed radiation-nudging conditions after Marchesiello et al. (2001) were applied.
The lateral time-dependent boundary conditions were provided by the parent by means of one-way nesting.However, the information from the parent was not instantaneously superimposed to the ROMS solution; additional nudging was applied to all prognostic variables (except for the sea surface elevation), which allowed these fields to adjust slowly to the parent values at the boundaries within an e-folding timescale of 2 days.In addition, a factor of 5 was used for the nudging timescales, which caused a stronger nudging on the inflow.

Surface boundary conditions
At the sea surface, the boundary conditions for the air-sea exchange of fresh water, momentum, and heat were evaluated from the outputs of two numerical weather prediction models and from the measurements of the meteorological buoy on top of M1 (see Fig. 1).This was accomplished by means of the wind field at 10 m (2 m for M1) of height, air temperature and relative humidity at 2 m, air pressure at sea level, total cloud cover (not available from M1), net shortwave radiation, and precipitation.The output of the weather prediction models was made available by the Italian weather service CNMCA (Centro Nazionale di Meteorologia e Climatologia Aeronautica) in two different set-ups, COSMO-ME and COSMO-IT.COSMO-ME covers the entire Mediterranean Sea with a horizontal resolution of 7 km and provided 72 h forecasts, while COSMO-IT encompasses Italy and the adjacent waters at the very high resolution of 2.2 km; however, the forecast range was only 24 h.The temporal resolution of both models was 1 h.The time series of all available variables from COSMO-ME, COSMO-IT, and the meteorological buoy are shown in Fig. 6 at the M1 position.

Data assimilation
In most of the model runs presented below, the temperature and salinity data from the shipborne CTD probes and gliders were assimilated using objective analysis (OA; see Bretherton et al., 1976;Carter and Robinson, 1987;Thomson and Emery, 2014).Namely, ROMS includes a module which enables data assimilation with the 4D-Var method.However, as 4D-Var is based on variational methods, it is rather expensive in terms of computer resources; according to parallel ROMS runs using 4D-Var (A.Funk, personal communication, 2016), the CPU time increases by about a factor of 10 compared to OA.During the integration of ROMS, the engine conducting the data assimilation was invoked every day at midnight, and it was controlled by six parameters: -W : the width of the time window which determines what data are assimilated.In all ROMS runs below, W = 48 h; this setting was found to provide the best forecast skill (Onken, 2017).Hence, all temperature and www.ocean-sci.net/13/235/2017/Ocean Sci., 13, 235-257, 2017 R. Onken: Mixed-layer prediction salinity data of the previous and the following 24 h were selected for assimilation.
-C: the isotropic correlation length scale.C = 15 km was used throughout, which is approximately the internal Rossby radius of the Western Mediterranean in summer (Grilli and Pinardi, 1998).Isotropic correlation is a strong assumption, especially close to the coast.However, according to the observations from the ADCP measurements (I.Borrione, personal communication, 2016), predominantly meridional currents prevailed only in a 10 km wide strip along the Sardinian coast, while the rest of the 180 km wide model domain was characterised by an eddy field with alternating currents.Here, the usage of a non-isotropic correlation scale would deteriorate the results.
-δT obs , δS obs : the observational errors of temperature and salinity, respectively.δT obs = 0.5 • C and δS obs = 0.16 were used throughout.These values were obtained from the variance of all CTD casts in the upper thermocline.

Integration and output
All ROMS runs presented below were initialised on 1 June at 00:00 and integrated forward for 24 days until 25 June at 00:00.To satisfy the horizontal and the vertical CFL criterion, a baroclinic time step t = 108 s (800 steps per day) was chosen, and the number of barotropic time steps between each baroclinic time step was 40.Harmonic mixing along the isopycnals with an eddy diffusivity coefficient of 5 m 2 s −1 was used for the horizontal diffusion of the tracers T and S, and a horizontal viscosity coefficient of 1 m 2 s −1 was selected for the diffusion of momentum.Further on, a quadratic law using a coefficient of 0.003 was applied for the bottom friction, and the pressure gradient term was computed using the standard density Jacobian algorithm of Shchepetkin and Williams (2001).
The three-dimensional volume of all prognostic fields was written to an output file at 6 h intervals.For comparison of the ROMS results with the observed records at mooring M1, the time series of the vertical temperature profiles right at the position of M1 were written to an extra file at the full temporal resolution.

Sensitivity of near-surface temperature and mixed-layer depth
The purpose of this section is to investigate the impacts of the following on the temperature between the surface and about 42 m of depth (which was the vertical range of the M1 observations) and the depth of the mixed layer: initialising ROMS from different data sets; the set-up of the vertical grid; different atmospheric forcing patterns; different vertical mixing schemes; and the background eddy diffusivity.
This was achieved with 5 series of ROMS runs named A-E (see Table 2 for the parameter settings and the results of each model run) with a total of 28 runs.The task of each series was to assess the sensitivity of the ROMS forecast skill to variations in the mechanisms mentioned in the bullets above.For each run, the ability of ROMS to predict the temperature was assessed by means of the root mean square (rms) difference between the observed temperature T obs and the predicted temperature T ROMS at each depth level of the observations.T was evaluated for the period of 15 June at 00:00 to 20 June at 13:55 where N observations were available in 5 min intervals (Sect.2.1).This interval was selected because it enabled the comparison of all runs with those which were forced by data assimilation until 12 June at 00:00.The 3-day lag between the last assimilation on 12 June and the start of the evaluation period on 15 June was granted to ROMS in order to recover from "assimilation shocks" which frequently become noticeable in the form of strong inertial oscillations.The experience from the precursor model runs has shown that such oscillations die off after about three to four inertial periods (18.7 h at 40 • N).In order to synchronise the modelled and the observed temperature, T ROMS was linearly interpolated in space and time on the observations.The equivalent method was also applied to the mixed-layer depth, D, which due to the lack of salinity observations at the M1 position was defined as the depth at which the temperature was 1 • C colder than the temperature at the surface for the first time (Lamb, 1984;Wagner, 1996).Hence, is the rms difference of the mixed-layer depths.

Series A: initialising ROMS from different data sets
In this series, h c = 10 m was selected for the critical depth.
In the first run, referred to as A1, ROMS was initialised from MFS, while in A2 the initial conditions were provided by MERCATOR.A3 was initialised from MERCATOR as well, but temperature and salinity data from the CTD casts and 10 gliders taken during 7-12 June at 00:00 were additionally Figure 8a shows the time series of the near-surface temperatures at 0.81 m of depth from runs A1-A3 in comparison with the corresponding observations of the uppermost Starmon sensor in M1 at the same depth level.In A1 and A2, the predicted temperatures agree reasonably well with the observations after 15 June, but before then the temperature exceeds the observations by several degrees.Extreme differences are visible during 12-14 June with differences of close to 3 • C. Figure 6 shows that during this period the predicted and observed wind speeds were close to 0 m s −1 and the shortwave radiation flux reached maximum values of more than 800 W m −2 .Hence, as these quantities are the major drivers of the mixed-layer temperature, it is concluded that the selected parameterisation of the vertical mixing in ROMS is not adequate for such calm situations.By contrast, as soon as the wind became stronger after 14 June, the maximum difference between the predicted and measured temperature is less than 1 • C. In A3 before 12 June, there is better agreement between the modelled and the observed temperature.However, as can be seen from the sudden drop in the modelled temperature at midnight on 10-12 June, the data assimilation led to an underestimation of the surface temperature.The reasons for this are twofold: first, some of the assimilated profiles started at 2 or even 3 m of depth because the measurements close to the surface were not reliable.In such cases, the uppermost measurements were extended to the surface and led to an underestimation of the near-surface temperature, which was sometimes significant because of the extremely shallow or even non-existent mixed layer.Second, the OA "advected" properties from the neighbouring casts which were not representative for the M1 position.On 13 June, the modelled temperature again exceeds the observations by almost 2 • C, but the difference is less than in A1 and A2.After 15 June, the A3 temperature is rather close to the temperatures in A1 and A2.As a skill measure for the forecasted near-surface temperature, T was evaluated for all runs and resulted in T = 0.30 • C in A1, T = 0.53 • C in A2, and T = 0.51 • C in A3 (see also the legend box in Fig. 8 and Table 2).The temporal evolution of the mixed-layer depth is displayed in Fig. 8b.As revealed by the decreasing rms differ-ences D between the modelled and observed mixed-layer depths, the forecast skill increases from A1 to A2 and from A2 to A3.The close agreement between the observed and modelled mixed-layer depth in A3 before 12 June, which was forced by the assimilation, is noteworthy.The mismatch between the model and the observations during 12-15 June is also remarkable as another indication of an inadequate parameterisation of the mixed-layer dynamics at low wind speeds.The vertical distribution of the rms temperature differences T of all runs in the A series is shown in Fig. 9.It is demonstrated that at most depth levels, T is lower or equal in A2 compared to A1.The assimilation in A3 led to a further significant decrease between about 4 and 35 m of depth; only above 4 m and below 35 m of depth is T higher in A3.The generally better forecast skill of A3 is also supported by T , the vertical mean of T , which is greater than 1 • C in A1 and A2 but only 0.90 • C in A3 (see also Table 2).In summary, nesting ROMS in MERCATOR and assimilating CTD profiles provided the best forecasts for the temperature and the depth of the mixed layer and the thermocline temperature below about 4 m of depth.Therefore, all runs in the B series will be based on A3.
The temporal evolution of the modelled temperature in A3 at the position of mooring M1 is shown in Fig. 10b.In comparison with Fig. 10a, the modelled temperature close to the sea surface is too high on 13 and 14 June, while at depths greater than about 3-10 m, T A3 appears too low.This is confirmed by Fig. 10c, which exhibits the temperature difference T A3 −T M1 : in approximately the upper 2 m depth range, T A3 partly exceeds T M1 by about 2 • C on these days, and just below, T A3 is up to more than 3 • C lower than the observed temperature.This aberrant cold layer can be identified during the whole model run.Apparently, the modelled mixed-layer depth is shallower than the observed one.This is illustrated by the vertical temperature gradient in Fig. 10e.Namely, a comparison with Fig. 10d reveals that the generally descending trend of the maximum gradient is similar, but the depth of the modelled maximum is always less than the observed one.Moreover, the observed variability is significantly higher than the modelled one.While for the entire period there is clear evidence of a strong diurnal variability in the observations (e.g. the deep mixed layer in the early morning and the shallow mixed layer in the afternoon), the modelled variability is much less pronounced.Another feature worth mentioning is that the thermocline is too warm during the assimilation phase before 12 June (Fig. 10c).It has been verified that this was caused by the assimilation of the glider data because this feature is not present in a run where only casts from lowered CTD were assimilated (not shown).As can be seen from Figs. 3 and 7, two CTD casts were taken exactly at the M1 position, while numerous glider casts are close to M1 (note that the meridional offset of the glider tracks with respect to the CTD meridional sections was 5 km).Thus, as the correlation scale of the OA was 15 km, the modelled temperature at M1 was primarily determined by the glider measurements because the large number of glider profiles reduced the statistical weight of the two CTD casts.

Series B: sensitivity to the set-up of the vertical grid
If the transformation equation, the vertical stretching function, and the total number of layers are held constant, the layer thicknesses of the ROMS vertical grid are controlled by the surface and bottom control parameters, θ s and θ B , and the critical depth, h c .For mixed-layer modelling in shelf areas, it would be desirable to have a high vertical resolution close to the surface, which can be achieved by either increasing θ s or decreasing h c .However, as increasing θ s would make the vertical transformation more non-linear, it was decided to keep θ s = 5 constant and vary only h c .In this series, the sensitivity of the ROMS results to five different settings of the critical depth is investigated using h c ∈ {10, 20, 50, 100, 200} m.For each of these choices, the impact on the layer thicknesses at the position of mooring M1 is illustrated in Fig. 11.A minimum layer thickness of 0.27 m right at the sea surface is achieved by h c = 10 m in run B1, while the thickness of that layer gradually increases in B2-B5.In the latter (h c = 200 m), the thickness is close to 1.3 m.B1, because it is identical to A3, is the control run.
For all runs in this series, the temporal evolution of the mixed-layer properties is displayed in Fig. 12.Although still too high around 14 June, the near-surface temperatures in all runs of this series most resemble the observations during the entire integration period, which is also expressed by the corresponding low values for T ; the minimum T = 0.44 • C  is obtained from B5, while the highest is found in B1 ( T = 0.51 • C).For the mixed-layer depth, there is no clear evidence of which run might do best.D varies only in a rather narrow range between 2.62 m in B1 and 2.75 m in B5.The vertical distributions of T (Fig. 13) and the vertical averages T are almost identical for all runs.However, right at the surface, T is minimal in B5 as shown in Fig. 12a.As the above results did not reveal a clear tendency of which choice for h c yielded the best results, it was decided to continue with B1 (h c = 10 m) as the control run in series C below.This decision was guided by Bernie et al. (2008), who asserted that a minimum vertical resolution of 1 m is mandatory to resolve the diurnal cycle of the sea surface temperature.Another criterion for this decision was the rx1 grid parameter (i.e. the Haney condition, after Haney, 1991) being at a minimum in B1 (see Table 2).

Series C: sensitivity to atmospheric forcing
Series C consists of three model runs, C1, C2, and C3.C1 is identical to B2; in C2, the surface boundary conditions were provided by COSMO-IT instead of COSMO-ME.In C3, the atmospheric forcing was defined by means of the observations of the meteorological buoy on top of mooring M1.Here, the observations were spread uniformly across the entire model domain whenever available.If no observations were available, i.e. before 8 June and after 20 June, the atmospheric fields of COSMO-ME were used.As observations of cloudiness were not available from M1, the corresponding fields from COSMO-ME were used throughout.
According to Fig. 14a, the predicted near-surface temperature from C2 closely resembles that of C1, except for 14-  17 June when the temperatures in C2 are about 1 • C higher.Apparently, this was driven by the different wind forecasts of the weather prediction models (Fig. 6).Before 14 June, the wind forecasts of both models were almost identical, but for the following 2 days during a period of stronger winds, the forecasts differ from each other.Overall, the nearsurface temperature does not appear to be very sensitive to the choice of the weather forecast models.This is also expressed by T , which attains similar values of 0.51 and 0.42 • C. The signature of the temperature changes considerably when ROMS was driven by the weather observed at M1; this is already evident during 8-10 June when the modelled temperature in C3 is different from C1 and C2.After 15 June, it is mostly higher than both the observations and the predictions of C1 and C2, which correspondingly leads to a higher T of 0.80 • C. With respect to the modelled mixedlayer depth (Fig. 14b) and based on the D criterion, C1  is superior to C2 and C3, but the large discrepancies during 12-15 June between the predictions and the observation are still present in all three runs.This corroborates the above hy-pothesis that the mismatch is not caused by the atmospheric forcing because the most appropriate forcing was applied in C3.
A surprising result was obtained from the vertical structure of the rms temperature difference (Fig. 15).Below about 3 m of depth, T C1 is about 0.1 • C lower than T C2 , but a considerable improvement in the predicted stratification is provided by C3.In the entire vertical range below about 5 m, T C3 is up to 0.4 • C lower than T C1 .Only right at the surface is T C3 approximately 0.3 • C higher than the corresponding values from C1 and C2, which is obviously due to the above-mentioned mismatch after 15 June.C3 provides the best results for the temperature stratification in the thermocline.As the temperature in this depth range was definitely not affected by the heat exchange at the sea surface (≈ 90 % of the shortwave radiation is absorbed in the uppermost 1 m depth range), its improvement could only be achieved by lateral advection, which is controlled by the wind; apparently, the wind is better represented in the observations than in the weather forecasts.To summarise, the objective skill measure T for the near-surface temperature and D for the mixed-layer depths indicate that C1 provides the best forecast, while T C3 is clearly superior to T C1 and T C2 in the thermocline.The latter confirms the decision to use C3 as a control run in series D because the advective processes are obviously reproduced best.
The decision for C3 is supported by Fig. 16.By visual inspection, the evolution of the predicted temperature pattern in C3 (Fig. 16c) resembles the observations (Fig. 16a) more than in C1 (Fig. 16b).Namely, the near-surface temperature is too high, but the thickness of the warm layer during 16- 20 June is roughly the same as in the observations, close to 10 m.Moreover, the depth and the variability of the maximum vertical temperature gradient in C3 resembles the observed pattern to a larger degree (Fig. 16d, e, f), although the vertical temperature gradient is still too weak.

Series D: sensitivity to the vertical mixing parameterisation
The GLS scheme (Umlauf and Burchard, 2003) provides a generalisation of a class of differential length-scale equations used in turbulence models for oceanic flows.Commonly used models, like the k-kl model of Mellor and Yamada (1982), the kmodel (Rodi, 1987), and the k-ω model (Wilcox, 1988), are recovered as special cases of the generic scheme.
Here, k is the turbulent kinetic energy, l is the length scale of the turbulence, is the dissipation rate, and ω is the specific dissipation rate.In series A-C, the GLS vertical mixing scheme was applied using its generic parameters as formulated by Umlauf and Burchard (2003).In the following, D1 is identical to C3 serving as the control run, the GLS scheme with the k-kl parameterisation is applied in D2, the kparameters are applied in D3, and the k-ω parameterisation, which was adjusted to oceanic conditions by Umlauf et al. (2003), is applied in D4.
After 12 June, the near-surface temperature of all runs is correlated with the observations (Fig. 17a), but is mostly still too high.Moreover, the graphs indicate that the temperatures from D2, D3, and D4 are closer to the observed ones, which is also expressed by T D2 = 0.50 • C, T D3 = 0.51 • C, and T D4 = 0.41 • C, while T D1 = 0.80 • C. For the mixed-layer depth (Fig. 17b), the best agreement with the observations was obtained from D4 with D D4 = 2.71 m.However, the mixed layer was mostly too shallow in all runs in this series.Hence, based on the T and D criteria, the k-ω mixing scheme in D4 definitely performs the best.This is also supported by the vertical structure of T displayed in Fig. 18.There is clear evidence that the k-kl scheme (D2), the kscheme (D3), and the k-ω scheme (D4) do better than the generic GLS (D1).Between the surface and about 5 m of depth, the best result was obtained from D4.Therefore, D4 will serve as the control run in the following E series.
An indicator of why the k-ω parameterisation performed better than the other closure schemes is possibly found in the publication of Reffray et al. (2015).Here, a one-dimensional model implemented in a three-dimensional circulation model was used to investigate physical and numerical turbulentmixing behaviour.Amongst others, the k-kl, the k-, and the k-ω scheme were compared to each other.It turned out that the k-ω scheme was the most sensitive to the vertical resolution.In a coarse (about 10 m) resolution model, k-kl and kclearly did better than k-ω, while at a high (about 1 m) resolution, all three schemes yielded suitable results.In the D series, the vertical resolution close to the sea surface is 0.27 m (see Fig. 11 and Sect.4.2 above).Hence, one may speculate that the k-ω formulation becomes superior to the other schemes when the vertical resolution is increased.

Series E: sensitivity to the background vertical eddy diffusivity
The shortcoming of all the model runs conducted so far was that the mixed layer was too warm and too shallow, and the thermocline was too cold with respect to the observational data.This is also in agreement with the findings of Reffray et al. (2015).Hence, it was conjectured that the parameterisation of the vertical transport of heat and/or momentum was not adequate.Several attempts were undertaken to finetune the D4 results by varying the vertical eddy viscosity coefficient and the turbulent closure parameters, but the outcomes were sobering; a significant improvement in the forecast skills for the mixed-layer properties was not achieved.Hence, in this series, the background vertical eddy diffusiv- ≤ 3 × 10 −5 m 2 s −1 in E4 and E5, which is somewhat higher than (1.7±0.2)×10−5 m 2 s −1 obtained from the tracer measurements in the thermocline during the North Atlantic Tracer Release Experiment (Ledwell et al., 1998;Thorpe, 2007).By contrast, the minimum of D = 2.05 m is found in E9 for A VT = 7 × 10 −5 m 2 s −1 .
Figure 20 shows the near-surface temperature in E4 and the mixed-layer depth in E9 together with the corresponding quantities of the control run E1 and the observations.After 15 June, the increase in A VT from 1×10 −6 to 2×10 −5 m 2 s −1 shifted the near-surface temperature by about 0.1 • C closer to the observations.Most of the time, the modelled signal is correlated with the observations, although the modelled maximum and minimum temperatures are frequently lagged a few hours behind the observed extreme values.Similar features were also described by Gentemann et al. (2009) when comparing time series of observed sea surface temperatures with those generated by the model of Fairall et al. (1996a).In their improved model (see the Introduction), they demonstrated that the peak warming in the afternoon was shifted earlier.For the mixed-layer depth, the increase in the eddy diffusivity to 7 × 10 −5 m 2 s −1 caused a significant reduction in D from 2.71 m in E1 to 2.05 m in E9.While in the precursor series the mixed layer was always too shallow, it now agrees remarkably well with the observations, except for T between the modelled temperature T ROMS and the observed temperature T obs evaluated at the actual depths of the observations.
T was computed only for the period after 15 June at 00:00.
large discrepancies on 19 and 20 June where the predicted mixed layer is up to 4 m shallower than the observed one.As the M1 wind speed was very low on these days (Fig. 6), other processes leading to a deepening of the mixed layer were probably inadequately parameterised, such as Langmuir circulation and wave breaking ( Noh et al., 2011( Noh et al., , 2016)).

Temporal variability
In order to assess the modelled temporal variability of the temperature and the depth of the mixed layer, the normalised spectra of the near-surface temperature amplitude T at 0.81 m of depth and of the mixed-layer depth amplitude D were computed by the Fourier transform, both from the observations and the ROMS outputs of runs E4 and E9, respectively.To enable a sufficient spectral resolution for the cycle periods of around 1 day, the entire time series between 8 and 20 June was used as input for the Fourier transform.At first glance, the modelled spectrum of the near-surface temperature (Fig. 21a1) resembles the observations in the cycle period range between about 0.1 and 1 days, but significant differences are evident in the bands between about 0.1 and 0.4 days where the modelled amplitude is up to 1 order of magnitude different from the observed one.This mismatch is not surprising, because here the temporal variability is controlled mainly by internal waves that are either not reproduced or are only partially reproduced by the model.By contrast, the range of 0.4-0.8days (10-19 h) is dominated by tides and inertial motions.Theoretically, at 40 • N, the inertial peak is at 18.7 h (0.78 days), but a correspondingly small peak is only visible in the modelled spectrum; no such peak is noticeable in the observations.Probably, the modelled peak is a leftover of the assimilation shock on 12 June.Additional peaks are found in both the modelled and the observed spectrum at about 0.4, 0.5, and 0.6-0.7 days (≈ 10, ≈ 12, and ≈ 14-17 h).While the sources of the first and the latter are unknown, the 12 h peak might be related to a semi-diurnal tidal component.However, as there was no tidal forcing in the ROMS version utilised in this study and the MERCATOR forcing at the lateral boundaries was defined by means of daily averages, the semi-diurnal variability could only be caused by tides in the atmosphere.Both the modelled and the observed spectrum are dominated by the diurnal variability represented by the peak at 1 day.In the red part of the spectrum between 1 and about 10 days, the modelled and observed amplitudes exhibit some weak correlation, and they are of about the same order of magnitude.This matter is not discussed here because it is potentially impacted by long-period fluctuations in the forcing at the surface and at the lateral boundaries.More detailed information on the correlation corr r TROMS , Tobs between the modelled and the observed temperature amplitudes is shown in Fig. 21a2.The correlation coefficient r = 0.74 together with the p value p = 3.05 × 10 −22 proves a high significant correlation, and the regression coefficients a 0 = 0 and a 1 = 1.74 indicate that, in general, the modelled amplitudes are overestimated.By contrast, there is less but still significant correlation between the modelled and the observed mixed-layer amplitudes DROMS and Dobs (Fig. 21b2), which is indicated by corr r DROMS , Dobs = 0.50, and p = 4.52 × 10 −9 .This finding is also supported by the spectrum (Fig. 21b1) in which a slight correlation of the amplitudes is only found for the diurnal and semi-diurnal cycles.

Horizontal variability
In order to assess the capability of ROMS to reproduce and predict the horizontal variability of mixed-layer properties, the results of run E9 were analysed along the ScanFish tracks A03, A05, A07, A09, and A10 (see Fig. 5) and compared with the data collected by the towed device.E9, using A VT = 7×10 −5 m 2 s −1 , was selected for this comparison because both T and D were acceptable.The details of the ScanFish tracks are summarised in Table 3.As ROMS output was only available in 6 h intervals starting at midnight, in each case the output cycle was used which fell within the time window when the tracks were conducted.This assumed synopticity of the ScanFish tracks is justified by the fact that the maximum duration of the tracks was 5 h 28 min for A03.
To make the ScanFish observations and the ROMS products comparable, the ScanFish temperature was interpolated vertically on 1 dbar standard levels, and the ROMS temperature was mapped on the same levels.As the upper inflection point of the ScanFish varied between about 5 and 10 dbar, there was frequently no information on the near-surface temperature available.In such cases, the temperature at the inflection level was extended to the surface.The same method was applied to the ROMS temperature, which was not defined right at the surface but in the centre of the first s layer below the surface.In deep-water regions, this was located at about 3 m of depth.
Figure 22a shows a temperature section from the Scan-Fish measurements along the central track of A05, and the corresponding section from ROMS is displayed in Fig. 22b.The overall features of both sections resemble each other, but  Table 3.The timing and nominal positions of the ScanFish tracks considered in this study (cf.Fig. 5).The ROMS analysis determines the instant of the model output which was used for comparison.the small-scale horizontal variability of the ScanFish temperature was not reproduced by ROMS.This is probably due to the smoothing effect of the OA, the combined action of the horizontal eddy diffusivity, and the numerical diffusion.However, as the last assimilation cycle was conducted on 12 June, 10 days prior to the ScanFish observations, one may exclude the possibility that the OA removed the small-scale features.Moreover, the vertical temperature gradient is much weaker in ROMS, which was already noted above.Hence, this is apparently not caused by the increased vertical diffusivity but by the vertical resolution of ROMS.The sea surface temperatures and the mixed-layer depths from the ScanFish and ROMS are displayed in Fig. 22c and d.For the surface temperature, the observed large-scale west-east trend is re-produced by ROMS, but there are differences of up to 0.5 • C in the central portion of the section.The maximum differences between the modelled and the observed mixed-layer depth in the 0-20 km range are close to 5 m at 13 km of distance, while in the eastern half of the section, the modelled and the observed mixed-layer depths approach each other.However, the smaller-scale O(1 km) observed variability was not reproduced by ROMS for both the sea surface temperature and the mixed-layer depth.
To investigate why the small-scale variability was not predicted correctly, run E9 was repeated using a smaller horizontal eddy diffusivity coefficient of 1 m 2 s −1 instead of 5 m 2 s −1 , which was used for all the model runs so far.However, no significant changes were noticeable.Thus, one has to settle for the fact that the present set-up of ROMS is only able to reproduce the horizontal variability of mixed-layer properties on scales which are comparable to the Rossby radius.

Conclusions
ROMS has been utilised to diagnose and predict the properties of the ocean mixed layer.The sensitivity of the model results to the choice of the initial and boundary conditions, the set-up of the vertical grid, and the vertical mixing schemes were investigated.The initial and lateral boundary conditions for ROMS were taken from two different parent models through one-way nesting.At the surface, ROMS was forced by two different weather forecasts or by observations.All ROMS nowcasts and forecasts were validated against observations which were taken in June 2014 to the west of Sardinia in the Mediterranean Sea.
To explore the sensitivity of the near-surface temperature and the mixed-layer depth to the choice of the initial conditions, ROMS was alternatively initialised by the Mediterranean Forecasting System (MFS) and the global MERCA-TOR model.In addition, observed temperature and salinity data were assimilated.For validation, the time series of temperature were compared with the observations from a mooring.Initialising ROMS from MERCATOR instead of MFS provided better agreement between the model and the observations, but significant improvement was obtained from a ROMS run initialised from MERCATOR and updated with assimilated data from CTD casts and gliders.This applied both to the near-surface temperature and the mixed-layer depth as well as to the temperature distribution in the upper thermocline.
To investigate the impact of the surface boundary conditions, atmospheric forcing fields were taken from the weather prediction models COSMO-ME and COSMO-IT and from the observations of a meteorological buoy acting as a point source.With respect to the mixed-layer depth, the best agreement with the observations was obtained from a model run forced with COSMO-ME, while the near-surface temperature exhibited the best match when ROMS was forced by COSMO-IT.However, the stratification in the upper thermocline was best represented when the point source was applied.The obvious reason for this surprising result is that the momentum forcing was overestimated by both COSMO-ME and COSMO-IT.
For the vertical mixing, four different configurations of the GLS scheme of Umlauf and Burchard (2003) were applied, representing the generic version: the k-kl model of Mellor and Yamada (1982), the kmodel (Rodi, 1987), and the k-ω model (Wilcox, 1988).The best performance was obtained from the k-ω model.
Regardless of which initial conditions or surface boundary conditions were applied, the modelled mixed layer was always too shallow and too warm.Therefore, the background vertical eddy diffusivity coefficient, A VT , was varied over more than 1 order of magnitude.The best agreement of the mixed-layer temperature was obtained for A VT ≈ 2 × 10 −5 m 2 s −1 , while A VT = 7 × 10 −5 m 2 s −1 provided the best match of the mixed-layer depth with the observations.A positive and significant correlation was found between the modelled and the observed temporal variability in the mixed-layer temperature.The modelled variability resembled the observed variability predominantly for cycle periods in the spectral ranges between about 0.5 and 1 days.By contrast, less correlation was found between the modelled and the observed variability in the mixed-layer depth.Slight agreement was only found for the diurnal period.
The horizontal variability was validated against measurements from a high-resolution zonal ScanFish section.Both the modelled mixed-layer temperature and the mixed-layer depth closely resembled the observations, but only on the larger scales of O(10 km).Hence, the mesoscale variability was rather well reproduced, but the sub-mesocale variability was not.

Figure 1 .
Figure 1.The design drawing of mooring M1.For explanations, see the text.

Figure 3 .
Figure 3.The CTD casts taken by Planet and Alliance during 7-11 June and the position of mooring M1.The colour bar indicates the water depth [m].

Figure 4 .
Figure 4.The actual glider tracks during 8-23 June.The small circles along the tracks show the surfacing positions.Each glider is marked by a different colour.The colour code for the bathymetry is the same as in Fig. 3.

Figure 5 .
Figure 5.The ScanFish tracks of Alliance during 21-24 June.The colour code for the bathymetry is the same as in Fig. 3.

Figure 6 .
Figure6.The time series of the measured and predicted atmospheric parameters at the site of mooring M1 from the observations of the meteorological buoy on top of M1, COSMO-ME, and COSMO-IT.U -wind and V -wind denominate the zonal and meridional wind components, respectively.The cloud cover was not recorded at M1.The precipitation is not shown because no precipitation was predicted or measured during the entire period.

Figure 7 .
Figure 7.The actual surfacing positions of all assimilated gliders during 7-11 June.Each glider is marked by a different colour.The colour code for the bathymetry is the same as in Fig. 3.

Figure 8 .
Figure 8. ROMS runs A1, A2, and A3: the time series of the (a) near-surface temperature at 0.81 m of depth, (b) the mixed-layer depth (MLD), and the corresponding observations at mooring M1.The numbers on the abscissae indicate June dates.The period for which the data are assimilated is highlighted with grey shading.

Figure 9 .
Figure 9. ROMS runs A1, A2, and A3: the rms temperature differences T [ • C] between the modelled temperature T ROMS and the observed temperature T obs evaluated at the actual depths of the observations.The vertical mean T is written in the second column of the legend box.T was computed only for the period after 15 June at 00:00.

Figure 10 .
Figure 10.(a) The observed temperature at mooring M1, (b) the modelled temperature from ROMS run A3, and (c) the difference between the modelled and the observed temperature.The vertical temperature gradient from (d) the observations and (e) from A3.The instant of the last data assimilation is indicated by the the grey dashed vertical line.

Figure 11 .
Figure 11.The layer thicknesses at the position of mooring M1 for various assumptions of the critical depth h c .

Figure 12 .
Figure 12.ROMS runs B1-B5: the time series of (a) the near-surface temperature at 0.81 m of depth, (b) the mixed-layer depth (MLD), and the corresponding observations at mooring M1.The numbers on the abscissae indicate June dates.The period for which the data are assimilated is highlighted with grey shading.

Figure 13 .
Figure13.ROMS runs B1-B5: the rms temperature differences T [ • C] between the modelled temperature T ROMS and the observed temperature T obs evaluated at the actual depths of the observations.The vertical mean T is written in the second column of the legend box.T was computed only for the period after 15 June at 00:00.

Figure 14 .
Figure 14.ROMS runs C1, C2, and C3: the time series of (a) the near-surface temperature at 0.81 m of depth, (b) the mixed-layer depth (MLD), and the corresponding observations at mooring M1.The numbers on the abscissae indicate June dates.The period for which the data are assimilated is highlighted with grey shading.

Figure 15 .
Figure 15.ROMS runs C1, C2, and C3: the rms temperature differences T between the modelled temperature T ROMS and the observed temperature T obs evaluated at the actual depths of the observations.T was computed only for the period after 15 June at 00:00.

Figure 16 .
Figure 16.(a) The observed temperature at mooring M1, (b) the modelled temperature from ROMS runs C1 and (c) C3, (d) and the vertical temperature gradient from M1, (e) C1, and (f) C3.The instant of the last data assimilation is indicated by the the grey dashed vertical line.

Figure 17 .
Figure 17.ROMS runs D1-D4: the time series of (a) the near-surface temperature at 0.81 m of depth, (b) the mixed-layer depth (MLD), and the corresponding observations at mooring M1.The numbers on the abscissae indicate June dates.The period for which the data are assimilated is highlighted with grey shading.

Figure 18 .
Figure 18.ROMS runs D1-D4: the rms temperature differencesT between the modelled temperature T ROMS and the observed temperature T obs evaluated at the actual depths of the observations.T was computed only for the period after 15 June at 00:00.

Figure 19 .
Figure19.Series E: T and D for ROMS runs E1-E13.Both quantities were computed only for the period after 15 June at 00:00.

Figure 20 .
Figure 20.The time series of (a) the near-surface temperature at 0.81 m of depth from E1 and E4, (b) the mixed-layer depth (MLD) from E1 and E9, and the corresponding observations at mooring M1.The numbers on the abscissae indicate June dates.The period for which the data are assimilated is highlighted with grey shading.

Figure 21 .Figure 22 .
Figure 21.The spectra and correlation parameters of the modelled and observed amplitudes (a) T of the near-surface temperature at 0.81 m of depth from run E4 and (b) D of the mixed-layer depth from run E9.The spectra were evaluated for the entire time series during 8-20 June where observations were available.

Table 1 .
The nominal and actual depths of the Starmon and RBR sensors mounted on the appendix of mooring M1.For the meaning of the recorded variables, see the text.

Table 2 .
The parameter settings and results of the ROMS runs in series A-E.The bold text indicates the parameters or boundary forcing patterns which are varied within the respective series.The best run in each series is marked by an asterisk and serves as the control run for the successive series.