1 Introduction

Surface currents affect surface waves through many different mechanisms that induce changes in the wave frequency, wave height, and the propagation direction of the waves. Spatially varying currents cause refraction-related changes in the propagation direction of the waves. Furthermore, waves encountering opposing currents will steepen, whereas waves affected by aligned currents will increase in wavelength.

Even though in many practical applications waves can be predicted with sufficient accuracy while neglecting wave-current interactions, the effects can be decisive in certain conditions. Opposing currents can lead to sudden increases in SWH and cause changes in other wave parameters, which can affect navigation safety (Chen 2018). Accounting for currents can also be important when estimating local wave power potential in coastal areas (Barbariol et al. 2013).

In areas with strong currents, such as western boundary currents, the wave-current interaction has been quite evident. These areas have thus been a focus of many modelling and measurements studies (e.g. Sugimori 1973, McKee 1977, Hayes1980, Holthuijsen and Tolman 1991, Wang et al. 1994). These strong shear currents are able to refract waves, inducing significant changes to the wave direction and the wave length and energy. In addition, small scale variability in the current fields, e.g. eddies and fronts of order of 10–100 km, has been shown to cause variability patterns in the SWH field on the same scales (Ardhuin et al. 2017; Quilfen et al. 2018). Also, Wang et al. (2020) have shown that wave-current interactions are significant when compiling wave statistics. According to them, including surface currents in a wave model changed the monthly mean SWH up to 0.2 m and Tp up to 0.5 s (relative changes varied up to 20%).

In coastal regions, current-induced changes in the wave field are typically caused by strong, tide- or storm-induced currents. For example, Ardhuin et al. (2012) showed that including the effect of tidal currents in the wave model reduced the errors in SWH by over 30%. However, wave-current interaction is also important during storm conditions, when strong currents occasionally occur. In situations of strong wind–driven currents, the currents and simultaneous waves have roughly the same propagation direction. Thus, in areas with strong wind–driven currents, wave-current interaction typically causes waves to grow longer and the SWH to be reduced. For example, Benetazzo et al. (2013) reported a reduction in SWH of up to 0.6 m during two typical high wind events in the Adriatic Sea. Wang and Sheng (2016) similarly showed a general reduction of SWH results during slowly moving storms. However, they also noted that during very fast moving storms, SWH was reduced by over 11% one side of the storm track and increased by around 5% on the other side of the track.

Wave-current interaction has been studied for a long time (e.g. Johnson 1947, Longuet-Higgins and Stewart 1960, Kenyon 1971) and even though the formulations have been included in the third-generation wave models (e.g. Komen et al., 1994), these effects have only recently been accounted for in operational wave forecasting systems (e.g. Staneva et al. 2015).

In the Baltic Sea, the mean surface currents are relatively weak, c. 0.05–0.10 m s− 1. During storms and high wind conditions, surface currents can be considerably stronger, up to 0.6 m s− 1 and in narrow straits even higher (Stanev et al. 2018). Surface currents have shown to affect the wave fields in the Baltic Sea in storm conditions. e.g. Viitak et al. (2016) found that surface currents induced up to 50-cm changes in SWH in the Gulf of Finland.

Until recently, none of the Baltic Sea operational wave forecasts have included the influence of surface currents. In July 2019, Copernicus Marine and Environmental Monitoring Service’s (CMEMS) Baltic Monitoring and Forecasting Centre (BAL MFC) launched a new wave forecast suite (BALTICSEA_ANALYSIS_FORECAST_WAV_003_010, (Vähä-Piikkiö et al. 2019)) that utilises surface currents calculated by the BAL MFC 3D ocean model suite (BALTICSEA_ANALYSIS_FORECAST_PHY_003_006, Golbeck et al. ()).

The validation of this new product revealed that the current refraction induced small changes in the SWH, especially near the coastal areas. This study expands on these findings by analysing a 2-year model simulation from the third-generation wave model WAM that was forced by precalculated surface currents from the 3D hydrodynamical model NEMO; as a reference, we ran the wave model without currents. We identify which areas, and under which conditions, the surface currents affect the wave field in the Baltic Sea the most. The results and current-induced changes are also evaluated against measured wave data from wave buoys and satellite altimeters using both integrated parameters and wave spectra.

2 Materials and methods

2.1 The WAM model

We conducted two different model runs with the third-generation wave model WAM (Komen et al. 1994) from 7 June 2016 to 1 July 2018. One run (WCUR) used currents U(x,t) varying in space x and time t, and the other run (hence reference run, WREF) did not (i.e. U(x,t) ≡ 0).

The wave model solves the action balance equation in spherical coordinates as defined by Komen et al. (1994):

$$ \begin{aligned} \frac{\partial N}{\partial t} + (cos\phi)^{-1}\frac{\partial}{\partial \phi}(\dot{\phi}\cos\phi N) + \frac{\partial}{\partial\lambda}(\dot{\lambda}N) + \frac{\partial}{\partial\omega}(\dot{\omega}N) \\ + \frac{\partial}{\partial\theta}(\dot{\theta}N) = S_{in} + S_{nl}+ S_{ds} + S_{bt} + S_{br}. \end{aligned} $$
(1)

Here, N = F/σ is the action density, σ is the intrinsic frequency, F(t,ϕ,λ,ω,𝜃) represents spectral density as a function of time t, latitude ϕ, longitude λ, angular frequency ω and wave direction 𝜃. Superscript dots denote partial time derivatives. The terms on the right side of the equation represent parameterised source terms: wind input Sin (Janssen 1991), wave–wave interactions Snl (Hasselmann et al. 1985), white-capping dissipation Sds (Komen et al. 1994), bottom friction Sbt (Hasselmann et al. 1973) and depth-induced wave breaking Sbr (Battjes and Janssen 1978).

On the left side of the equation, the partial time derivates are:

$$ \dot{\phi} = (c_{g}\cos \theta - u)R^{-1}, $$
(2)
$$ \dot{\lambda} = (c_{g}\sin \theta - v)(R \cos\phi)^{-1}, $$
(3)
$$ \dot{\omega} = \partial\boldsymbol{\varOmega}/\partial t, $$
(4)
$$ \dot{\theta} = c_{g}\sin\theta\tan\phi R^{-1} + \dot{\theta}_{D}, $$
(5)

where

$$ \dot{\theta}_{D} = \left( \sin\theta\frac{\partial}{\partial\phi}\boldsymbol{{\varOmega}} - \frac{\cos\theta}{\cos\phi}\frac{\partial}{\partial\lambda}\boldsymbol{{\varOmega}}\right)(kR)^{-1}, $$
(6)

where R is the radius of the earth, cg the group velocity, k the wave number vector and k its magnitude, and u and v are the north and east components, respectively, of the surface current field U. Ω stands for the dispersion relation:

$$ \omega = \boldsymbol{\varOmega} (\mathbf{k},\mathbf{x},t) = \sigma + \mathbf{k} \cdot \mathbf{U} = \sqrt{gk \tanh{(kh)}} + \mathbf{k} \cdot \mathbf{U}. $$
(7)

The model grid resolution was approximately 1 nmi (1.852 km) with 36 spectral directions and 35 frequencies (0.04177–1.06719 Hz). The open boundary was located in Skagerrak at longitude 8.0139 E. The bathymetry, grid and spectral resolution are the same as used in the CMEMS BAL MFC wave analysis and forecast system operational since July 2019 (marine.copernicus.eu, BALTICSEA_ANALYSIS_ FORECAST_WAV_003_010). The model was run in finite depth mode, implying that the presence of the sea bed affects the waves. The attenuation of wave energy due to subgrid-scale islands and islets that were smaller than grid resolution was also accounted for via grid obstruction fields (Tuomi et al. 2014). As the study focuses on the northern Baltic Sea, the incoming wave energy was defined as zero at the open boundary at Skagerrak. This will affect the wave model results in the Skagerrak, Kattegat and in some cases also in the Danish Straits. Therefore, we limit our analysis to the area east of the Bornholm islands (east from longitude 13 E). The study area is presented in Fig. 1.

Fig. 1
figure 1

Study area and locations of wave buoys. Legend shows abbreviations of the buoy names and the depth of the buoy location. BB represents the Bothnian Bay, BS the Bothnian Sea, FG the Finngrundet, NBP the Northern Baltic Proper and GoF the Gulf of Finland

In the presence of sea ice, surface waves are attenuated. We used a 30% ice concentration as the threshold for setting wave energy to zero (Tuomi et al. 2011; Tuomi et al. 2019). Ice information was taken hourly from the NEMO model, together with currents.

2.2 The NEMO model

Circulation was modelled with the 3D hydrodynamical NEMO model (Nucleus for European Modelling of the Ocean; Madec et al. (2017)) version 4.0. The model setup is based on the Nemo-Nordic configuration (Hordoir et al. 2019). The domain covers the English Channel, the North Sea and the Baltic Sea. A regular, evenly spaced latitude–longitude grid is used with a resolution of approximately 1 nmi (Hordoir et al. 2019). In the vertical direction, 56 equipotential z levels were used. Vertical resolution ranges from 1 m at the surface to about 24 m at a depth of 700 m. In the Baltic region, the horizontal NEMO grid coincides with the WAM grid. Sea ice evolution was simulated with the coupled SI3 model using 5 ice categories and one snow layer. The NEMO-SI3 simulation was carried out offline. Horizontal surface currents and sea-ice concentration were saved with a resolution of 1 h and used as forcing fields in WAM.

2.3 Model forcing

The WAM and NEMO models were forced with the ECMWF atmospheric forecast model medium range HRES product (Owens and Hewson 2018). Data was collected from 00 UTC runs using the first 24 h. In NEMO, hourly atmospheric pressure, 10 m wind, long and short wave radiation, 2 m air temperature and specific humidity, precipitation and showfall fields were used. WAM was forced with 10 m winds.

At the open boundaries, NEMO was forced with 3D water temperature and salinity fields, and sub-tidal water elevation and depth-averaged velocity fields from CMEMS 1/12 global ocean forecast (GLOBAL_ANALYSIS_ FORECAST_PHY_001_024, Lellouche et al. (2018)). In addition, tides were imposed by interpolating TPXO global tidal model data (Egbert and Erofeeva 2002) to the boundary. Initial conditions for water temperature and salinity were generated from a previous hindcast simulation and observation data sets.

The model bathymetry is based on the global 30 arc-second GEBCO bathymetric data set (GEBCO_2014 Grid, version 20141103). In NEMO, the bathymetry has been modified in certain shallow regions, e.g. the coastal region of Denmark, to improve water elevation and salinity inflow modelling skill with respect to tide gauge and mooring observations. WAM used the NEMO bathymetry with some differences as the wave model grid in the Baltic Sea accounts for the effects of irregular shoreline by estimating the effective fetch (Tuomi et al. 2012), and therefore the wave model shoreline is typically further away from the actual shoreline than in the 3D model grid. When providing the surface current fields to the wave model, only sea points that existed in both models were included, meaning that U = 0 for the ocean model land cells.

2.4 Wave observations

2.4.1 Wave buoys

Data from five wave buoys were used to evaluate model accuracy (Fig. 1). All buoys were Datawell Directional Waveriders, which use 27-min time series to calculate the wave spectrum, including the mean direction and directional spreading for each spectral component.

At the Gulf of Finland and Bothnian Sea sites, both Mk-III and DWR4 buoys have been used during 2016–2018. For the other three locations, the Mk-III has been used exclusively. The model Mk-III has a sampling time of 0.78 s and the wave spectrum is determined for the interval 0.025–0.580 Hz with a resolution of 0.010 Hz (0.005 Hz below 0.1 Hz). The newer DWR4 has a sampling time of 0.39 s and calculates the spectrum up to 1 Hz, with the higher 0.005-Hz resolution extending up to 0.25 Hz.

The Finngrundet buoy is operated by the Swedish Meteorological and Hydrological Institute (SMHI). Hourly wave data consisting of the peak period and significant wave height was downloaded from the CMEMS data portal (INSITU_BAL_NRT_OBSERVATIONS_013_032, https://marine.copernicus.eu/). The other four wave buoys are operated by the Finnish Meteorogical Institute (FMI); for these buoys, wave spectra were available every 30 min.

In the Baltic Sea, wave buoys are typically removed before the arrival of the seasonal ice cover (Tuomi et al. 2011). Nonetheless, year-round measurements can be possible during mild ice winters, which was the case for the wave buoy in Northern Baltic Proper for the model period June 2016 to June 2018. In the Gulf of Finland, the buoy was removed during both ice winters (13 January to 18 April 2017 and 29 January to 14 April 2018). The Bothnian Sea buoy was removed during the 2017/2018 ice winter (22 February to 8 April 2018), but data is also unavailable for 30 March to 9 May 2017 because of equipment failure. Finngrundet data was unavailable during 2 October 2016 to 6 May 2017 and 1 March to 4 May 2018. In the Bothnian Bay, the buoy was removed during both ice winters (8 February to 6 June 2017 and 30 November 2017 to 14 June 2018).

Shorter gaps, lasting less than a day, are present in the data during times of device maintenance or replacement.

2.4.2 Satellite altimeters

Satellite altimeter Ku/Ka-band SWH data was retrieved from the Radar Altimeter Database System (RADS) database (http://rads.tudelft.nl/rads/rads.shtml, Scharroo et al. 2013, Scharroo 2012). We used data from three altimeters (JASON-2, JASON-3 and SARAL) operating during our study period. Measurements flagged as good data and as sea (land/sea flag4 = 0) within longitudes of 13 to 30 and latitudes of 53.5 to 67 were used. We also included an additional quality assessment that was proposed by Kudryavtseva and Soomere (2016) when they evaluated the quality of JASON-2 and SARAL data in the Baltic Sea. According to their suggestions, data points with Ku/Ka-band backscatter coefficient > 13.5 dB and SWH normalised standard deviation > 0.5 m were removed. Values of SWH = 0 were also removed as unreliable.

2.5 Error metrics

We used standard statistical metrics to evaluate the performance of the model against measurements and to compare the two model runs with each other. These same measures were used in CMEMS BAL MFC wave product validation (Vähä-Piikkiö et al. 2019).

Bias,

$$ \text{bias} = \frac{1}{N}\sum\limits_{i=1}^{N}(m_{i}-r_{i}), $$
(8)

root mean square error:

$$ \text{RMSE} = \sqrt{\frac{1}{N}\sum\limits_{i=1}^{N}(m_{i}- r_{i})^{2}}, $$
(9)

and correlation coefficient:

$$ R = \frac{\text{cov}(m,r)}{\sqrt{\text{var}(m)\text{var}(r)}}, $$
(10)

where cov denotes covariance and var variance. Here, m refers to model and r to the reference (observation) time series. In model comparisons, WREF run was used as the reference series. In this case, RMSE was called root mean square deviance (RMSD).

3 Results

3.1 Validation

Validation of the model results was conducted against five wave buoys located in different areas of the Baltic Sea (Fig. 1, Section 2.4.1) and satellite altimeter data (Section 2.4.2). The significant wave height (SWH) and peak period (Tp, only wave buoys) were validated following the same methods as used in Vähä-Piikkiö et al. (2019). Bias, root mean square error (RMSE) and Pearson correlation (R) are presented in Table 1.

Table 1 Comparison of modelled significant wave height (SWH) and peak period (Tp) against wave buoys (locations shown in Fig. 1) and three satellite altimeters (Jason-2, Jason-3 SARAL)

Table 1 shows that there were only small differences in the quality of the model runs. The wave buoys are located offshore (Fig. 1), where currents are typically weaker in the Baltic Sea, so it was expected to not see many changes in the wave parameters at the buoy locations. The most notable differences can be seen at the Gulf of Finland (GoF) wave buoy, where the bias of the peak period changes from − 0.505 s (reference run WREF) to − 0.624 s (run with currents WCUR), and the RMSE of the same variable increases from 1.320 to 1.354 s. However, the bias of SWH at the GoF buoy decreased from 0.014 m (WREF) to 0.004 m (WCUR). As this buoy location seems to be most affected by current refraction, it is studied in more detail in Section 3.5.

Data from satellite altimeters was used to validate larger areas, including near-shore locations. The measurements were compared to the closest model grid point and allowing for a maximum 30-min time difference. A comparison with the altimeter measurements showed that the current refraction had only a small effect on the quality of the model results: RMSE decreased from 0.269 (WREF) to 0.268 m (WCUR) and correlation coefficient increased from 0.946 (WREF) to 0.947 (WCUR). Although altimeters have larger spatial coverage than wave buoys, their temporal coverage is rather limited, and most cases where the model runs showed significant wave-current interactions were not captured by the altimeters.

The validation shows that both model runs performed well in the Baltic Sea. The inclusion of the current refraction to the wave model caused no major changes to the model quality at the studied locations.

3.2 Surface currents

The mean current field over the study period (7 June 2016–30 June 2018) has an average speed of 0.02 m s− 1 (Fig. 2). However, there was large seasonal variation, and the monthly mean currents were occasionally as high as 0.06 m s− 1. The mean currents were strongest from October to December (0.04–0.06 m s− 1) and weakest during March–May (0.02–0.03 m s− 1). The variability in meteorological conditions also led to large year-to-year variation in the monthly mean values.

Fig. 2
figure 2

Mean surface current fields calculated from hourly NEMO output from June 7, 2016, to July 1, 2018. Colours indicate current magnitude and arrows direction. Directions are shown for every 15 longitude and latitude grids

The strongest currents were mostly limited to narrow areas by the coasts. The currents occasionally exceeded 0.5 m s− 1 on the western coast on Åland Sea (Å S), the western coast of the Northern Baltic Proper (NBP), the eastern coast of the Western Gotland Basin (WGB), Gdańsk Bay (GB) and some smaller areas on the coasts of the eastern Gotland basin (EGB). Changes in the wave fields by the coastal areas are discussed separately in Section 3.4.

The weakest mean currents were found in the Bothnian Bay (BB) and the Gulf of Riga (GoR). In most areas of these basins, current magnitudes exceeded 0.1 m s− 1 less than 25% of the time over the study period. Similarly, weak currents also occurred in the eastern part of the Gulf of Finland (GoF) and central areas of the Bothnian Sea (BS).

Although earlier studies have shown that NEMO simulates the general features of the Baltic Sea current field fairly well (e.g. Hordoir et al. 2019, Westerlund et al. 2019), it has not been extensively validated against current measurements. Long-term current measurements are rare in the Baltic Sea, especially in the northern parts, and data is mostly restricted to specific measurement campaigns. Nevertheless, as NEMO is one of the state-of-the-art models used in the Baltic Sea, we consider the current fields to be sufficiently well represented in order to study the wave-current interaction in the Baltic Sea.

3.3 Wave statistics

We calculated wave statistics from the reference run (WREF) and the run with currents (WCUR) to identify areas in which there were larger, or more frequent, differences caused by currents. The differences between two model runs were calculated by subtracting the reference run from the run with currents (WCUR−WREF). This means that positive differences indicate larger values in the WCUR and vice versa. Changes in SWH and Tp were evaluated independently. Statistics were calculated over the whole study period of 7 June 2016 to 1 July 2018 with metrics presented in Section 2.5. Ice was treated with type F statistics (Tuomi et al. 2011), meaning that statistics for each grid point were calculated only for times with no ice. According to FMI’s Ice service, ice winter 2016/2017 was mild, whereas 2017/2018 was average. During 2016/2017, there was significant ice cover only in the BB (Bothnian Bay), GoF (Gulf of Finland) and GoR (Gulf of Riga). In 2017/2018 also BS (Bothnian Sea) and NBP (Northern Baltic Proper) had an ice cover. Generally, the NEMO model tends to overestimate the ice cover to some extent (Tuomi et al. 2019).

3.3.1 Mean and maximum values

Figure 3a shows the mean SWH field for the reference run (WREF) averaged over the study period. The highest mean values of SWH were up to 1.31 m and mean Tp was up to 5.3 s (Fig. 4a) in both of the model runs. Overall, the mean values of this period showed similar characteristics to earlier studies (e.g. Tuomi et al. 2011, Björkqvist et al. 2018). The longest and highest waves were located in the the southern part of EGB. The BS, GoF and the GoR have significantly milder wave climates than the Gotland Basins (EGB and WGB) and the NBP.

Fig. 3
figure 3

Statistics of the significant wave height (SWH): mean SWH of the reference run (WREF, panel a), the difference of the mean SWH between the run with current refraction (WCUR) and WREF, calculated WCUR-WREF (b), root mean square deviation (RMSD) between WREF and WCUR (c), maximum SWH (WREF, d), and difference of maximum SWH (WCUR-WREF, e). Values are calculated over the whole study period 7 June 2016–1 July 2018

Fig. 4
figure 4

Statistics of Tp over the whole study period (7 June 2016–1 July 2018.): mean Tp (WREF, panel a), the difference of mean Tp (WCUR-WREF, b), and RMSD between WCUR and WREF (c)

The mean difference (WCUR - WREF) of SWH (Fig. 3b) varied between − 2.3 and 3.4 cm. A decrease in SWH due to currents was limited to narrow areas near the coast, whereas the increase in the mean values was more widespread and located mostly in the open sea areas. Nonetheless, there were a few coastal areas where the mean SWH increased, such as at the southern coast of the GoF and near the Gotland and Öland islands in the WGB. Areas where the currents are generally weaker, such as the GoF, had clearly smaller mean differences.

Similar to the mean difference, the root mean square deviation (RMSD) of SWH between the two model runs was also small in the GoR, the BS and the BB, together with easternmost areas of the GoF (Fig. 3c). RSMD values highlight the areas in WGB and the eastern coast of EGB. The maximum RMSD of 0.05 m was on the western coast of Gotland. As the mean difference between model runs was significantly smaller (0.02 m) than the RMSD (0.05 m), these larger RMSD values are most likely caused by the variation that the temporally varying currents induce to the wave field.

The difference in the maximum values showed large values ranging from − 0.51 m to 0.61 m (Fig. 3e). In many areas, the run with currents (WCUR) had larger maximum values of SWH than the reference run (WREF). The largest differences are not necessarily linked to the extreme storm conditions, as the maximum values could have occurred at different times. However, our study period includes storm Toini on 11 January 2017, during which the highest measured SWH at the NBP wave buoy was 8.0 m (Björkqvist et al. 2017). The highest modelled SWH during this event was 8.78 m located west from the NBP buoy (Fig. 3d). During this storm, current refraction induced the largest difference in maximum values of SWH (0.61 m), occurring after the highest peak of the storm on 12 January in the western part of the NBP. Maximum instantaneous differences of the modelled values are studied separately later in Section 3.3.3.

The mean difference in Tp (WCUR-WREF) was negative in most parts of the Baltic Sea (Fig. 4b), varying between − 0.25 and 0.12 s. Changes were mostly caused by spatial displacement of the wave field due to refraction. An increase in mean Tp was only seen in limited areas in the southern coast of GoF, ÅS and WGB. The mean Tp also increased slightly in the NBP. In general, the implementation of current refraction to the wave model shortened the peak wave period in the Baltic Sea, and it was most notable in GoF and the eastern coast of EGB.

The RMSD of Tp was highest in the entire GoF, the NBP and the WGB (Fig. 4c). The RMSD of Tp between the two runs was smallest by the BS, BB and the GoR, similar to differences in SWH. The largest RMSD is, however, in the ÅS with the maximum value being 1.0 s. As with SWH, the time-varying currents induce variation in the Tp field, which shows in higher RMSD values.

3.3.2 Changes in wave direction

Because the Baltic Sea is small, the wave field is typically wind wave dominated. Swell, however, occurs under certain conditions. We therefore evaluated the effect of currents on the wind sea and swell directions separately.

Currents steered swell more strongly than wind waves: the absolute difference in mean swell direction was over 2 in most areas, whereas mean changes in wind wave direction were mostly below 1 (Fig. 5). Differences in swell direction were largest in the Archipelago and Åland Sea, with mean differences of over 10. Wind wave direction differences of over 1 are mostly located close to the coast where the strongest currents occur (Fig. 2).

Fig. 5
figure 5

Absolute differences of mean direction for swell (a) and wind waves (b) calculated over (7 June 2016–1 July 2018). Note the different ranges of the colour bars

3.3.3 Maximal instantaneous differences

Maximum values of SWH presented in Section 3.3.1 could have occurred at different times in different parts of the sub-basins. To further evaluate the effect of currents on waves, we calculated the differences in SWH and Tp (WCUR-WREF) for each time step and evaluated the largest differences.

In terms of SWH, the maximum differences (Fig. 6) were seen throughout the Baltic Proper and in the Quark (Fig. 1). Including currents generally increased the wave height, but the opposite was true in some areas near the coast with the maximum increase and decrease of SWH (WCUR-WREF) ranging from − 0.26 m to 0.76 m.

Fig. 6
figure 6

Maximum differences of the SWH at the simultaneous timestamps (WCUR-WREF) during the study period (7 June 2016–1 July 2018). The maximum increase (a) and maximum decrease (b) are presented separately

Differences in the maximum values of Tp were large (Fig. 7), which was due to the shift from one spectral peak to another in situations with two peaked spectra. These spectral peaks refer to wind sea and swell and are further discussed in Section 3.5. These differences varied from − 12.28 to 12.28 s. The largest increase in Tp was seen in the ÅS and the northern coast of NBP.

Fig. 7
figure 7

The maximum differences of the Tp calculated at the same timestamp (WCUR - WREF) during the study period (7 June 2016–1 July 2018). The maximum increase (a) and maximum decrease (b) are shown separately

3.3.4 Probability of occurrence

Including currents rarely led to large changes in SWH or Tp, and the changes were typically restricted to a specific area. We calculated the fraction of the time steps when the difference in SWH exceeded 0.1 m. Similar calculations were done for Tp with a threshold of 1 s. The calculations were performed separately for positive and negative differences.

In general, the increase in SWH was more frequent than the decrease. There were certain areas, such as the NBP, western GoF and WGB, where an increase of over 0.1 m was found in over 0.5% of the time (Fig. 8a). A decrease due to current interactions occurring more than 0.5% of the time was seen mainly along the coasts of the EGB and WGB (Fig. 8b). Most often, the decrease in the SWH was seen on the eastern coast of EGB (more specifically by the westcoast of Saaremaa island). A decrease in SWH occurred in the same regions where currents over 0.3 m s− 1 were modelled. However, not all currents of this magnitude caused strong SWH differences of more than − 0.1m. There were large areas in the Baltic Sea where the difference in SWH between model runs was always within ± 0.1 m, e.g. most of BS and BB, GoR and eastern GoF.

Fig. 8
figure 8

Percentage of time the difference in SWH (WCUR-WREF) was over or equal to 0.1 m (a) and under/equal − 0.1 m (b) over the whole study period calculated with 1-h temporal resolution. Occurrence frequencies are presented as percentages. Markers show the locations of wave buoys (also presented in Fig. 1) together with locations A and B studied in Section 3.4

Differences larger than 0.1 m were, naturally, even rarer, for example differences exceeding 0.2 m occurred less than 0.2% of time (not shown).

Figure 9 shows the fraction of the time when current refraction increased Tp over 1 s (a) and decreased it over 1 s (b). The use of current refraction led to prominent increase of Tp in the ÅS, the WGB and southern coast of the GoF (Fig. 9a). The NBP and WGB had large areas where an increase of 1 s in Tp was seen over 0.5% of time. An increase of 1 s occurred over 1% of the time only in coastal areas.

Fig. 9
figure 9

Percentage of times the difference in Tp (WCUR-WREF) was over or equal to 1 s (a) and under (or equal to) 1 s (b) over the whole study period calculated with 1-h temporal resolution. Occurrence frequencies are presented as percentages. Markers show the locations of wave buoys (also presented in Fig. 1) together with locations A and B studied in Section 3.4

A decrease in Tp resulting from wave-current interactions occurred in much larger areas than increases (Fig. 9b). There were large areas in the Bornholm Basin, the WGB, the northern half of NBP, ÅS, the Quark and the GoF where the current forcing decreased Tp by more than 1 s at more than 0.5% of the time. A decrease was also prominent in the eastern coast of EGB and the Gdańsk Bay (GB). A decrease by 1 s was seen over 2% of the time only in narrow areas near the coast.

3.3.5 Seasonal variation

To further evaluate when currents typically caused larger changes in the wave fields, we looked into the seasonality of the changes presented above in Section 3.3.4.

Large changes in SWH, of over 0.1 m, were most common during winter, especially in December and January (not shown). Spring had the fewest cases of notable changes in SWH. This is in accordance with the typical seasonal variability in the Baltic Sea region. Winds, surface currents and waves are typically stronger in autumn and winter compared to conditions during spring and summer.

The seasonal variation in the Tp differences was considerably smaller. However, strong changes were less frequent in the spring compared to other seasons.

3.4 Coastal regions

Wave-current interactions typically caused a decrease in SWH in the coastal areas (Figs. 3b and 6b) where the strongest currents also occurred (Fig. 2). Strong currents are topographically steered and flow in an along-shore direction, whereas waves typically propagate towards the coast at some angle. A varying current field with a strong gradient along the coast causes the waves to be refracted away from these coasts, thus increasing SWH in the open sea or at the opposite coast (Figs. 3b and 6a).

The increase in SWH due to refraction from strong coastal currents was clearly seen, for example, during storm Toini (Section 3.3.1), currents along both the eastern and western coasts refracted the waves, thus increasing SWH by over 0.2 m in large areas of central NBP. However, the largest changes of SWH (up to 0.61 m) in the western part of the NBP were a combination of refraction in the NBP and strong opposing currents in the Åland Sea. This more direct interaction with concurrent wave and current propagation directions is presented in more detail in this section.

The magnitude of the wave-current interaction depends on the angle between the wave and current propagations, based on dispersion relation (Eq. 7). A perpendicular angle between the wave and current directions causes no changes to the wave fields. In all other cases, currents have some component that is either along or opposing the direction of wave propagation, thus affecting the wave field. When the difference in direction is below 90, currents have an along-wave component that decreases angular frequency (ω in Eq. 7), whereas currents with opposing direction respectively have component that increases ω.

We selected a few of the locations with the most frequent changes of SWH (Fig. 8) and Tp (Fig. 9) for a closer analysis. Modelled data from each location was divided into two cases: currents that were either along, or opposing, to the wave direction. Components were identified by calculating the absolute difference between current and mean wave direction at spectral peak (WREF). Cases with opposing and following currents are evaluated separately and the magnitude of the current component is represented by the intensity of the colour in the following Figs. 10 and 11. The locations of these two areas, which we decided to present as examples, are marked in Figs. 8 and 9 as A and B (correspondingly).

Fig. 10
figure 10

Difference of SWH (WCUR-WREF) by the coast of Stockholm archipelago (point A in Fig. 8). Red represents situations when currents had opposing component to the waves and blue when the component was following. Colours indicate the absolute direction difference between WREF peak direction and current flow direction

Fig. 11
figure 11

Difference of SWH (WCUR-WREF) by the NW coast of Gotland (point B in Fig. 8). Red represents situations when currents had opposing component to the waves and blue when the component was following. Colours indicate the absolute direction difference between WREF peak direction and current flow direction

The areas where we most often see a larger decrease in SWH (Fig. 8b) are the same locations where the strongest currents occur. In these coastal areas, the current direction was mainly along the wave direction (per the definition given two paragraph prior). Location A (shown in Fig. 8) near the Stockholm archipelago is a good example of this type of coastal area, where strong currents aligned with the wave direction led to lower values of SWH compared to the reference run (Fig. 10, blue). At this location, opposing currents occurred rarely (Fig. 10, red). This same type of situation was also very evident in the eastern coast of EGB.

The largest changes in the SWH were seen when the current speed exceeded 0.3 m s− 1, having a strong alongside current component. When looking at the opposing currents (shown in red in Fig. 10), the importance of the relative angle between waves and currents is very evident. For example, currents with speed of 0.45 m s− 1 hardly affected SWH as the opposing fraction of this current was very small. Figure 10 also shows some large changes in the SWH with small current component (light coloured dots). These changes are more likely due to refraction of the waves caused by currents in the other regions of the basin. Generally, current speeds of less than 0.1 m s− 1 caused only small differences to the SWH that varied between both positive and negative directions.

In most of the coastal areas, currents were mainly aligned with the waves, like by the Stockholm archipelago (Fig. 10), where currents were opposing the waves on only 30% of the occasions. However, there were a few areas that more frequently had currents opposing the wave field direction. One of these areas was the western coast of Gotland (location B in Figs. 8 and 9). There, currents of over 0.5 m s− 1 increased SWH close to 0.4 m (Fig. 11, in red). A decrease of SWH was also prominent in this region, where the currents are aligned (Fig. 11, in blue), similar to the situation by Stockholm archipelago (Fig. 10, in blue). This behaviour explains well the mean increase in SWH and Tp around the WGB. However, current speeds did not correlate as clearly with the increases in these parameters in the ÅS and the GoF, where current refraction seems to play a larger role in the chances of the wave field.

3.5 Gulf of Finland

The Gulf of Finland (GoF) is a special area with respect to the wave field, since the narrowness of the gulf limits wave growth under high wind conditions and consequently the wave direction is typically steered towards the axis of the gulf (e.g. Pettersson et al. 2010).

Although the GoF wave buoy is not in an optimal location considering the largest differences seen in this area, the spectral data available from this buoy allows us to demonstrate how surface currents change the modelled wave field and how it compares to the measurements. We selected a case from August 2016, when there was a large difference in the predicted peak period together with some small differences in the SWH between WREF and WCUR.

The time series of Tp (Fig. 12) showed that on 10–11 August, there were up to 4 s differences in the simulated peak period between WREF and WCUR, with WREF typically showing higher values than WCUR. There were also a few instances where the differences in the peak period were accompanied by differences in the mean direction at spectral peak.

Fig. 12
figure 12

Timeseries from GoF buoy location, shown in Fig. 1. Panel a with Δ SWH and c with ΔTp show the difference between the two model runs (WCUR-WREF) together with the current speed at the same location (green dashed line and right axis). Panels b and d show the modelled values compared to the buoy-measured ones (black lines). Panel e shows the directions of each value

The buoy spectra (Fig. 13) showed that the old wind sea energy started to decrease and turn into swell in the evening of 10 August, and a new wind sea starts to develop, exceeding the energy of swell peak at 18 UTC on 10 August. The two peaked spectra with swell towards NE and wind sea towards N prevailed until midday on 11 August. The decrease in swell energy at the buoy location was due to swell being refracted by the currents and shifting towards the southern part of the gulf, increasing the SWH and Tp on the SE coast. The differences in the SWH and Tp between the southern and northern part of the gulf, seen in the earlier analysis of wave statistics (Figs. 3b4b) and the frequency of the Tp differences (Fig. 9), are mainly caused by these types of events. Large shifts of Tp occurred most often in the GoF region, yet still with frequencies of only 3% (Section 3.3.4). The inclusion of wave-current interaction improves the modelled two-peaked (swell and wind sea) spectra, as without currents the energy at the swell peak is much higher than in the measurements. When currents are used, the swell wave energy decreases faster and the wind sea peak starts to have more energy than the swell peak, shifting the peak period and peak wave direction to wind waves and showing better agreement with the measured spectra.

Fig. 13
figure 13

Development of wave spectrum between between 12:00 on 10 August 2016 and 6:00 on 11 August 2016 by GoF buoy (location shown in Fig. 1). Left axis shows spectral power and right axis corresponding direction (dots)

4 Discussion

Overall, the effects induced by surface currents on the simulated wave fields were relatively small. Even near the coasts, where currents are strongest, magnitudes rarely exceeded 0.3 m s− 1, which was shown to be the threshold for considerable changes in wave fields. The current effects were typically seen in areas from which very little measured wave data was available. The wave buoys are located in the central parts of the basins, where the effect was small or negligible and altimeter data, although covering wider spatial areas, was temporally limited and concurrent measurements during significant current-induced changes were rare. As the tuning of the wind–wave coupling in the wave model has been performed for the WREF, the validation of the SWH only presents the differences between the two setups, not necessarily the actual accuracy of that of an operational product.

A detailed evaluation of the performance of simulated wave-current interaction would require more measurements from coastal locations where the changes were most prominent, e.g. the Gulf of Finland and the Eastern Gotland basin. Availability of both wave (with spectrum) and current measurements from such locations would enable the evaluation of the quality of NEMO surface currents and their impact on the simulated wave fields.

The cases in which the current forcing induced large changes to the wave fields required both high current speeds, of over 0.3 m, and a substantial difference between the wave and current directions. Our results showed that in most areas, where current effects were prominent, there was a decrease in the SWH. In a few areas, opposing currents induced an increase in SWH. The largest differences were during autumn and winter, when the winds and waves are typically highest in the Baltic. For example, during storm Toini, currents refracted some of the wave energy towards the central basin of the Northern Baltic Proper, thus increasing the SWH in the central locations. As our simulation period is short, only 2 years, the number of high wind and wave events is quite small. To more comprehensively evaluate the wave-current interaction in high wind and storm conditions, a longer simulation period is needed. The present climate reanalysis, such as ERA5 (Hersbach et al. 2020), provides sufficiently good forcing to perform longer simulation, even in relatively small sea areas. However, a long time-series of surface current forcing with a high enough spatial and temporal resolution to solve the small-scale features of the Baltic Sea is not readily available.

Although the current-induced effects on the wave field are relatively small in the Baltic Sea, except for some high wind and storm events, the occasional, considerably large changes in the peakwave period and also changes in the gradients of SWH encourage the use of current forcing in operational systems. When planning for the coastal structures, the typical differences in the model runs done without current refraction should be noted, especially when planning for structures that are vulnerable to extreme events.

So far, avoiding the inter-dependencies of different operational systems in addition to the relatively small effect on the model results has hindered the use of current refraction in the operational wave forecast systems. In the CMEMS BAL MFC wave production system, the surface current effects have been accounted for since the July 2019 upgrade. In December 2020, the system was upgraded to utilise the surface current and ice concentration forcing from the new BAL MFC physical analysis and forecasting system based on NEMO. The results presented in this study reflect the improvements in the BAL MFC wave forecast system, although the model setups used in this study have some differences to the operational versions regarding the model setups, bathymetry, meteorological forcing and boundary conditions.

5 Conclusion

We performed two different wave model runs to see how much, where and how often sea surface currents affected the wave fields in the Baltic Sea. This was done by comparing two WAM model runs, one with and one without hourly forcing from sea surface currents (from NEMO) for the period of 7 June 2016–1 July 2018. Both WAM and NEMO were forced with winds from the ECMWF medium range atmospheric forecast system.

Validation of the model runs against five operational wave buoys and satellite altimeters showed only small differences between model runs. The buoys were located in the central parts of the basins, where the currents had no notable effects on the wave field and altimeters did not have a sufficiently high temporal resolution to catch these rarely occurring changes between the two model runs. Only one of the buoys was located in an area with notable changes, namely the Gulf of Finland.

Currents had only a small effect on the mean statistics of the Baltic Sea wave climate. Mean SWH varied up to 2 cm; however, during storms, the local maximums of SWH differed up to 60 cm. Instantaneous differences between model runs showed even higher variation (up to 76 cm), indicating that wave-current interaction can occasionally affect wave fields significantly. On average, currents induced a decrease in Tp in most areas of the Baltic Sea (up to 0.2 s). This decrease in Tp was most pronounced in the Gulf of Finland and the eastern coast of the Eastern Gotland Basin. Changes in Tp were mostly caused by spatial displacement of the wave field due to refraction. However, larger changes in Tp were connected to changes in swell energy, which occasionally shifted the peak of the spectral power between swell and wind sea. A comparison of mean directions of the swell and wind sea showed that currents turned swell more strongly than wind sea.

Calculation of the occurrence differences showed that, even though maximum differences between model runs are large, they happen rarely. Current-induced changes of over ± 0.1 m in SWH or ± 1 s in Tp occurred only in few areas with frequencies of around 3%. Many areas, like the Bothnian Sea and the Bothnian Bay, experienced hardly any changes of this magnitude. There was also clear seasonal separation in the occurrences of these current-induced changes. Currents changed SWH most often during the winter and least during the spring. Seasonal variation was not as notable for changes in Tp; however, they were less frequent during the spring when when the weakest mean currents were simulated.

In general, wave model runs with currents tend to predict weaker values of SWH by the coasts and larger values in open sea areas. The areas of more frequent decreases in the SWH were confined to the areas where currents of over 0.3 m s− 1 occurred. However, there were few coastal areas, e.g. the coasts of the Western Gotland basin, the Åland Sea and the southern Gulf of Finland, which were locally affected with more severe wave conditions (longer and higher waves) when current refraction was included in the simulations. In these basins, more severe wave conditions were due to currents opposing waves and refracting long waves to the opposite coast. We showed an example of the opposing currents in the Western Gotland Basin and studied the refraction of the longer waves in the Gulf of Finland, where we had buoy measurements. As currents refracted the swell, they decreased the swell energy in the northern and central part of the gulf and increased it in the southern part of the gulf, thus changing the peak period, we saw a clear improvement in the prediction of swell in the Gulf of Finland.