Operational Service for Mapping the Baltic Sea Landfast Ice Properties

Abstract

The Baltic Sea is partly covered by sea ice in every winter season. Landfast ice (LFI) on the Baltic Sea is a place for recreational activities such as skiing and ice fishing. Over thick LFI ice roads can be established between mainland and islands to speed up transportation compared to the use of ferries. LFI also allows transportation of material to or from islands without piers for large ships. For all these activities, information on LFI extent and sea ice thickness, snow thickness and degree of ice deformation on LFI is very important. We generated new operational products for these LFI parameters based on synthetic aperture radar (SAR) imagery and existing products and prediction models on the Baltic Sea ice properties. The products are generated daily and have a 500 m pixel size. They are visualized in a web-portal titled “Baltic Sea landfast ice extent and thickness (BALFI)” which has free access. The BALFI service was started in February 2019. Before the BALFI service, information on the LFI properties in fine scale (<1 km) was not available from any single source or product. We studied the accuracy and quality of the BALFI products for the ice season 2019–2020 using ice charts and in-situ coastal ice station data. We suggest that the current products give usable information on the Baltic LFI properties for various end-users. We also identify some topics for the further development of the BALFI products.

1. Introduction

The Baltic Sea is a semi-enclosed brackish sea water basin in Northern Europe. The ice cover in the Baltic Sea usually begins to form in November, and has its largest extent between January and March [1,2,3]. The normal ice break-up starts in April and the ice melts completely by the beginning of June. The maximum annual ice cover ranges from 9% to 100% of the whole Baltic Sea area, and the average is 50% [3,4]. The Baltic Sea ice extent and duration of the ice season depend on the indices of the Arctic Oscillation (AO) and North Atlantic Oscillation (NOA)—e.g., the annual maximum ice extent generally decreases with increasing indices of AO and NAO [2]. The ice in the Baltic Sea occurs as landfast ice (LFI) and drift ice. LFI occurs in the coastal and archipelago areas, and usually extends to the 5–15 m isobath. Drift ice has a dynamic nature due to forcing by winds and currents. The motion of drift ice results in an uneven and broken ice field with distinct floes up to several kilometers in diameter, leads and cracks, brash ice barriers, rafted ice and ice ridges. In the Bay of Bothnia, the annual maximum level ice thickness is typically 0.65–0.80 m, and it reaches 0.3–0.5 m even in mild winters [1,2]. The measured all time maximum is around 1.2 m. In the Southern Baltic Sea, the coastal areas of Germany and Poland and the Danish Straits, the annual maximum level ice thickness seldom exceeds 0.5 m [5]. The thickness of ice ridges (sail height plus keel depth) is typically 5 to 15 m [6]. The salinity of the Baltic Sea ice is typically from 0.2 to 2‰ depending on the location, time and weather history [7].

Wintertime ship navigation and icebreaker operations in the Baltic Sea rely heavily on sea ice information provided by national ice services, such as the Finnish Ice Service (FIS) at Finnish Meteorological Institute (FMI). The most important source of sea ice information satellite is synthetic aperture radar (SAR) imagery, which can be acquired regardless of the amount of solar illumination or cloud cover. Other sources include optical and thermal infrared imagery, sea ice model data and in-situ sea ice observations. The most important sea ice parameters comprise the location of the sea ice edge, sea ice and snow thickness (SIT and SNT), degree of ice deformation (DoD) as well as sea ice concentration (SIC) and location of larger leads [8].

Every winter, LFI exists at least along the coasts of the Bay of Bothnia (coastal states are Finland and Sweden) and the Gulfs of Finland and Riga (Finland, Estonia, Latvia and Russia), and in the Finnish Archipelago Sea [9]. During an average winter, LFI also covers the coast of the Sea of Bothnia (Finland and Sweden). In a severe winter, it is also present in the Southern Baltic Sea (e.g., Danish Straits). The average extent of LFI in late winter, 1st March [9], is shown in Figure 1. LFI cover is not always restricted to coastal areas and skerries but occurs also sometimes as “fast ice bridges” in the narrow parts of the Baltic Sea such as the Quarken and the Åland Sea. Statistics on the average duration of the Baltic Sea ice season, and the average and the largest measured SIT on LFI, are shown in

Table 1. Statistics on the Baltic Sea ice season duration and the thickness of the landfast ice. Compiled by the Finnish Ice Service.

Average Duration of the Ice Season		Average Ice Thickness	Largest Ice Thickness
Bay of Bothnia	4–6 months	0.4–0.7 m	>1.0 m (1.22 m/21.4.1985)
Sea of Bothnia	1–2 months	0.2–0.4 m	~0.7 m
Archipelago Sea	3–4 months
Gulf of Finland, western parts	1–3 months	0.2–0.5 m	>0.7 m (0.85 m/2003)
Gulf of Finland, eastern parts	>4 months
Northern Baltic Proper	<1 month	0.05–0.2 m	~0.5 m

. Along the Finnish and Swedish coasts, there is LFI cover from one month up to a half a year.

The Baltic Sea LFI zone may extend from a few hundred meters to several tens of kilometers from the coast and is the place for many recreational activities such as skiing, skating, snowmobiling and ice fishing. Ice roads (for cars and trucks) can be established over thick LFI between mainland and populated islands (e.g., Hailuoto in the Bay of Bothnia) to speed up transportation compared to the use of ferries. LFI also allows transport of material to or from islands without piers for large ships or usage of slow and costly barges. For all these activities, the extent and SIT of LFI are very important. Other important parameters include SNT and structure of ice (e.g., columnar ice vs. granular snow-ice), as granular sea ice may be mechanically weak and break under heavy machinery, or even under a single person. In addition, for the safety of the activities on LFI, the location of ship tracks in the sea ice needs to be known. Close to the coast ships sail only along fixed shipping lanes to and from harbors, but at the LFI edge, when water is deep enough, routes may vary depending on ice conditions (SIT and deformation). Historical LFI information is necessary for the planning of the location of off-shore infrastructure, such as wind power turbines, and lighthouses. Economic impacts on having available information on the Baltic Sea LFI properties are related to cost and time savings, increased safety (less accidents and search and rescue operations) and improvements of services.

Currently, Baltic Sea ice information is freely available by ice charts prepared by national ice services such as FIS and at Swedish Meteorological and Hydrological Institute (SMHI). The ice charts are prepared manually using all available information, the SAR imagery being the most important one (more detailed description will be in Section 2.5) [8].

There are also automated Baltic Sea ice products produced and distributed by Copernicus Marine Environment Monitoring Service (CMEMS) [10] and various national institutes. The products by CMEMS are freely available, whereas some by national institutes are not. Here, we focus only on the CMEMS Baltic Sea ice products. They include: (1) Baltic Sea—Sea Ice Concentration and Thickness Charts [11]; (2) Baltic Sea—SAR Sea Ice Thickness and Drift, Multisensor Sea Ice Concentration [11]; and (3) Baltic Sea Physics Analysis and Forecast [12].

The first CMEMS product gives SIC and SIT information on a 1 km spatial scale. These two parameters are based on the gridded digitized Baltic ice chart produced by the ice analysts at FIS or SMHI. The ice chart is produced daily in the afternoon. Although the grid size is acceptable, the ice chart segments of equal SIT and SIC can be very large—e.g., up to 100–200 km in diameter.

The second product includes four datasets: (a) gridded SIT chart based on the digitized ice chart refined by a SENTINEL-1 (S-1) Extra Wide (EW) swath or RADARSAT-2 (RS-2) ScanSAR Wide (SCW) SAR image [11]; (b) gridded ice drift based on the S-1 or RS-2 SAR image pairs [13]; (c) gridded daily SIT mosaic covering the whole Baltic Sea; and (d) gridded SIC daily mosaic over the Baltic Sea based on the SAR data and the Advanced Microwave Scanning Radiometer 2 (AMSR2) radiometer data [14].

The gridded SIT chart is based on a combination of ice chart SIT data and SAR image analysis [15]. First, the SAR image is preprocessed by performing an incidence angle correction and speckle filtering. Next, an image segmentation is performed and segmentwise image autocorrelation feature median is computed. This textural feature is used to locate the open water areas, typically having low autocorrelation. After this, the SAR segments are mapped to corresponding ice chart polygons. A SAR segment corresponds to an ice chart polygon which covers the largest fraction of the segment. For all SAR segments mapped to a certain ice chart polygon, SIT is interpolated between the minimum and maximum SIT of the polygon. The mapping of SIT is performed linearly such that the lowest segment-wise backscattering coefficient ($\mathsf{\sigma }°$) median corresponds to the ice chart SIT minimum and the highest $\mathsf{\sigma }°$ median to the ice chart SIT maximum. Finally, the open water areas are updated according to the SAR open water classification based on thresholding the segment-wise autocorrelation computed earlier. The pixel size of the gridded SIT chart is 500 m. The SIT mosaic is produced twice a day by overlaying the most recent swath SIT data over the older one.

Sea ice drift estimation is based on the computing phase correlation of pairwise SAR data windows sampled from two images in two resolutions [13]. The nominal resolution for the sea ice drift data is about 800 m. SIC estimation is based on SAR image segmentation and a multilayer perceptron (MLP) neural network having as input segment-wise SAR textural features and AMSR2 polarization ratio and spectral gradient ratio signatures [14]. The MLP has three layers—i.e., one hidden layer, and the structure of 29 × 60 × 1 neurons. The number of inputs consists of 12 SAR features, 16 AMSR2 features, and a constant value. The activation function of the hidden layer is a sigmoid function, and for the output layer a linear one. The single output is the SIC estimate. The SIC field of the digitized ice chart was used as training data. The daily SIC chart has a 500 m pixel size.

The Baltic Sea Physics Analysis and Forecast product contains model calculations for the physical condition in the Baltic Sea [12,16]. It includes hourly and daily mean data for SIC and SIT. The product is based on simulations with the HBM ocean model code (HIROMB-BOOS-Model; HIROMB denotes High Resolution Operational Model for the Baltic sea, and BOOS is Baltic Operational Oceanographic System). The horizontal resolution of the model data is around 1 nautical mile (NM) (i.e., 1852 m).

None of the aforementioned CMEMS products explicitly give the extent of LFI. An end-user may estimate the extent from time series analysis of the sea ice drift data (LFI has zero drift speed), but this is time consuming and requires some knowledge on the SAR-based ice drift estimation. The SIT for LFI is in the products, but without any ice type information, such as LFI vs. drift ice, and thus, only at very close to coastline (<1–2 km) one can assume SIT data to correspond typically to that of LFI.

LFI extent and SIT information are available in the Baltic Sea ice charts by FMI and SMHI, see an example in Figure 4. However, the extent information is based only on visual analysis of satellite imagery (mainly SAR) and is not always updated on a daily basis. Additionally, as the operational ice chart is mainly aimed to serve the winter navigation, the manual LFI analysis gets less attention in areas where shipping is less active—i.e., shallow or large archipelago areas far away from the main shipping routes. The SIT information is mainly based on in-situ thickness measurements at coastal ice stations. The ice charts are available to the general public only as image data (e.g., as pdf-files), and not as WMO SIGRID-3 vector archive format data. We are not aware of any products by national Ice Services or institutes which give information on the Baltic Sea LFI properties in fine scale (<1 km).

In this study, we developed methods and products for operational monitoring of the Baltic Sea LFI properties, including LFI extent and thickness, snow thickness (SNT) on LFI, and degree of ice deformation (DoD) (qualitative index). The LFI products are based on C-band SAR imagery, a sea ice thermodynamic model run at FMI, and CMEMS Baltic Sea products. Our LFI products have around 500 m pixel size (the actual resolution is coarser for some LFI parameters). They are distributed through a service portal titled “Baltic Sea landfast ice extent and thickness” (BALFI) [17]. The BALFI portal has free access and it is targeted for people living on the Baltic Sea coast and islands who require information on LFI properties for recreational activities such as skiing, skating and ice fishing, and for transporting people and goods along ice roads to/from islands, governmental and local authorities/institutions such as national Ice Services, icebreaker management, harbor authorities and commercial activities such as tourism. The BALFI service was started in February 2019, and it is active during the Baltic Sea ice season, typically from November/December to May.

The Baltic Sea LFI extent estimation is based on time series of S-1 and RS-2 SAR imagery. The qualitative index on DoD is also estimated with the SAR imagery. SNT on LFI is estimated with a high-resolution one-dimensional thermodynamic snow and ice model HIGHTSI [18,19,20] run at FMI. SIT is taken from the CMEMS Baltic Sea Physics Analysis and Forecast product [12], which does not include SNT. The CMEMS model includes sea ice dynamics, but HIGHTSI only ice thermodynamics. SIT, SNT and DoD data are mapped over the LFI area derived from the SAR imagery.

The BALFI products and service were developed in 2017–2019 within a CMEMS User Uptake program titled “Development and promotion of demonstrations of CMEMS downstream services”. A CMEMS downstream service is targeted to provide either (a) value-added ocean products, or (b) information or derived products developed from information provided by CMEMS. Therefore, CMEMS Baltic Sea products were used as input data for the BALFI products.

In the following, Section 2 first describes the used datasets and their processing. Next, our algorithms and procedures for calculating the BALFI products are described in Section 3. This section also shows examples of the products and includes a short description of the BALFI web-portal. In Section 4 (Results), we investigate accuracy and quality of the BALFI products with available in-situ and ice chart data for the ice season 2019–2020. In Section 5 (Discussion), we compare the BALFI service to existing sea ice services, summarize the results on the product validation and identify topics for further developments. Finally, we give our conclusions in Section 6.

2. Datasets and Processing

In this section, the used data sets used and their pre-processing methods are presented. The datasets are for the Baltic Sea ice season in 2019–2020.

2.1. SAR Imagery

RADARSAT-2 ScanSAR and SENTINEL-1 EW mode Ground Range Detected Medium Resolution (GRDM) imagery with the HH/HV polarization combination are used for the BALFI products. The SAR imagery are freely available to FMI through CMEMS services. The preprocessing of the SAR images includes calibration; calculation of HH- and HV-polarization backscattering coefficients (${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ and ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$); thermal noise removal in the HV-polarization images; geo-rectification into a Mercator projection with 100 m pixel size; incidence angle correction in the HH-polarization images, land masking and SAR mosaicking [14,21]. This SAR image processing method is also used for the FMI’s CMEMS Baltic Sea ice products.

The Mercator projection is the same one used in the FMI and SMHI Baltic Sea ice charts, and it has a reference latitude of 61 2/3 degrees and the WGS84 datum. As $\mathsf{\sigma }°$ in the RS-2 ScanSAR and S-1 EW HH-polarization imagery is heavily dependent on the incidence angle (${\mathsf{\theta }}_0$), a correction is necessary before SAR image classifications. Here, an incidence angle correction method described in [22] was applied. The correction maps the ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ values using an empirical linear dependence for the ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ in dB-scale to a predefined ${\mathsf{\theta }}_0$ of 30°.

The ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$ values for sea ice are typically close to the RS-2 and S-1 SAR noise floors which are not constant, but have a variation as function of ${\mathsf{\theta }}_0$ (i.e., in the range direction). Noise floor modulation of the ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$ values in the SAR imagery complicates their visual and automated classification. The noise floor variation in the range direction was corrected by applying the noise floor data available in the metadata of each SAR image.

A mosaic of the SAR imagery is built twice a day by always overlaying newer S-1 and RS-2 imagery, such that at each mosaic grid cell the most recent SAR data are available. The SAR mosaic has Mercator easting-northing coordinates with a 500 m pixel size, and it covers the whole Baltic Sea—see Figure 1. The SAR mosaics have eight bit linear scaling of ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ and ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$ values in dB-scale—from −30 to 0 dB for ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$, and from −40 to 0 dB for ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$. A SAR mosaic at HH-polarization is depicted in Figure 2. Time series of daily SAR mosaics are used to estimate the LFI extent and DoD.

The land mask for the SAR imagery and mosaic was based on the Global Self-consistent Hierarchical, High-resolution Geography (GSHHG) database from the National Oceanic and Atmospheric Administration (NOAA) coastline data [23]. The same land mask is also applied to the BALFI products.

2.2. HIGHTSI Model

A high-resolution one-dimensional thermodynamic snow and ice model HIGHTSI [18,19,20] is applied in the BALFI for SNT estimation. The physical configuration of HIGHTSI is schematically illustrated in Figure 3. A surface layer is defined, and a snow and ice surface heat balance equation is used to calculate surface temperature and melt:

$$\left(1-\alpha \right)Q_s-I_0+Q_d-Q_b+Q_h+Q_{le}+F_c+F_m={\rho }_{i,s}{\left(L_f\right)}_{i,s}d{h}_{i,s}/dt,$$

(1)

where $\alpha$ is surface albedo; $Q_s$ is the incoming short-wave radiative flux at surface; $I_0$ is the solar radiation penetrating below the surface layer; $Q_d$ is the incoming atmospheric long-wave radiation; $Q_b$ is the long-wave radiation emitted by the surface; $Q_h$ and $Q_{le}$ are the turbulent fluxes of sensible heat and latent heat, respectively; $F_c$ is the surface conductive heat flux; $F_m$ is the heat flux due to surface melting; $L_f$ is the latent heat of fusion; $\rho$ is the density of ice or snow; $h$ is the thickness of ice or snow. Subscripts $i$ and $s$ denote ice and snow, respectively. $Q_h$ and $Q_{le}$ are calculated taking the atmospheric stratification into account on the basis of the observed wind speed ($V_s$), air temperature ($T_a$) and relative humidity ($R_h$) as well as the modelled ice surface temperature ($T_{sfc}$). A partial differential heat conduction equation solves the snow and ice temperature field:

$${\left(\rho c\right)}_{s,i}\frac{\partial T_{s,i}\left(z,t\right)}{\partial t}=\frac{\partial }{\partial z}\left(k_{s,i}\frac{\partial T_{s,i}\left(z,t\right)}{\partial z}\right)-\frac{\partial q_{s,i}\left(z,t\right)}{\partial z},$$

(2)

where $T$ is the temperature, $c$ is specific heat, $k$ is thermal conductivity, $q$ is the amount of solar radiation penetrating below the surface and $\left(z,t\right)$ is the vertical coordinate and time. The penetration of solar radiation within snow and ice is carefully parameterized which ensures the quantitative calculation of sub-surface snow and ice melting.

The freezing and melting at ice bottom are determined by the heat conductive flux at lower part of ice floe and upward oceanic heat flux below ice bottom:

$$-{(k_i\partial T_i/\partial z)}_{bot}+F_w=-{\rho }_i(L_f)dH/dt,$$

(3)

where ${\rho }_i$ is the sea ice density at the basal layer, $dH/dt$ is the ice growth rate at bottom and $F_w$ is the oceanic heat flux.

The thermal properties ($\rho$, $c$, $k$) of ice are not constants but parameterized according to [24,25]. The snow insulation effect was considered by taking into account a time dependent parameterization of ${\rho }_s$ [26] and $k_s$ [27]. The refreezing of snow flooding slush (snow-ice formation) and refreezing of surface and sub-surface melted snow (superimposed ice formation) at the snow-ice interface are considered in the model. The flooding slush is a result of isostatic imbalance of overload snow on top of total ice floe. For sea conditions, the snow-ice is different from the superimposed ice because slush is a mixture of snow and salty sea water. The details of model physics can be found in [18,28].

HIGHTSI was initially designed to investigate snow and ice mass balance for the brackish seas such as Baltic Sea and Bohai Sea [18,19,20]. The model has been further developed and applied for snow and sea ice studies in the Arctic [29,30,31] and Antarctic [32,33].

In the HIGHTSI model runs for the BALFI, the ice drift in the Baltic Sea is taken into account by incorporating information on SIC changes into the model. The SIC data are extracted from the CMEMS Baltic Sea—Sea Ice Concentration and Thickness Charts product, which has daily SIC information on a 1 km spatial scale [11].

For the BALFI SNT product, the HIGHTSI model uses atmosphere parameters from the ECMWF numerical weather prediction (NWP) forecasts. On each day, ECMWF forecasts made at 00:00 are applied. The forecasts last for 10 days. The forecasts are interpolated into 1 h time interval, and then a HIGHTSI run is then carried out for a 10-day period on a daily basis. When a grid cell is at least partly covered by sea ice (SIC > 20%), SIT and SNT are calculated at that particular grid cell. If SIC at a certain grid cell is reduced below 20%, SIT and SNT will remain as the values of the previous time step. The calculation of SIT and SNT are resumed once SIC is again larger than 20%. The HIGHTSI forecast at 24:00 is used as the initial condition for the following day’s forecasts. The HIGHTSI model grid over the Baltic Sea has 0.0665° × 0.133° resolution in latitude–longitude (196 × 166 grid points). In the Mercator projection used for the BALFI products, the grid resolution is 7.1 km in the easting and 5.9–8.6 km in the northing. The HIGHTSI SNT data at the morning 06:00 LT are used as daily SNT values over the LFI zone (from SAR mosaic analysis).

2.3. CMEMS Baltic Sea Products

2.3.1. Baltic Sea Physics Analysis and Forecast

The Baltic Sea Physics Analysis and Forecast product provides forecasts for the physical conditions in the Baltic Sea [12]. The Baltic forecast is updated twice a day providing a new six-day forecast with hourly data for sea level variations, SIC and SIT at the surface, and temperature, salinity and horizontal velocities for a 3D field. The product is produced by the 3D ocean model code HBM (HIROMB-BOOS-Model) developed within the Baltic ocean community. For details on the HBM ocean code, see [16] and references in [12,34]. Recent model developments include improvement of the thermodynamic routine which should improve especially the description of sea ice during growing and melting season, and the sea ice module by including fast ice [12].

The forecasting system used in CMEMS is denoted as HBM-CMEMS-V4. The product grid has a resolution of 1 NM in the horizontal (delta longitude is 1′40″ and delta latitude is 1′), and up to 25 vertical depth levels. The area covers the Baltic Sea including the transition area towards the North Sea (i.e., the Danish Belts, the Kattegat and Skagerrak). The HBM-CMEMS-V4 forecast is produced by the Danish Meteorological Institute (DMI) with forcing from the DMI’s HARMONIE 2.5 km system for the first 2.5 days, and the available ECMWF forcing for the next 3.5 days. No data assimilation is used in the ocean model. The HBM-CMEMS-V4 model is run twice per day for the 00 UTC or 12 UTC analyses; resulting 12 h best estimate for the hindcast period (before 00 UTC or 12 UTC), which substitutes the first 12 h of the previous analysis forecast, and the latest available 144 h forecast [12]. The Baltic ocean data are updated at the CMEMS server at 10 UTC and 22 UTC at the latest for the 00 UTC and 12 UTC analysis, respectively. For the BALFI service, SIT hourly data are downloaded twice a day at 06:00 and 18:00. The most recent SIT data are used as daily SIT values over the LFI zone.

The quality of the HBM-CMEMS-V4 forecast system for the physical condition in the Baltic Sea was assessed by using a two year hindcast for the time periods of 1 October 2013 to 30 September 2014, and 1 October 2014 to 30 September 2015 [34]. The total Baltic Sea ice extent (SIE; area with SIC ≥ 15%) and sea ice volume (SIV) were validated against the FIS ice charts. These two variables are derived from SIC and SIT data. During the first year long period, both SIE and SIV were close to the ice chart values during the freeze-up season lasting to around mid-January 2013. After that, they were both overestimated by the HBM-CMEMS-V4 up to mid-March 2014. For the rest of the ice season (up to mid-May) SIV was close to the ice chart data, but SIE was slightly overestimated. In 2014–2015, both SIE and SIV were continuously overestimated by the HBM-CMEMS-V4. The daily SIV and SIE overestimations were up to around 5 × 10³ km³ and 35 × 10³ km², respectively. The correlation coefficient ($r$) between the ice chart and HBM-CMEMS-V4 SIEs was very high, 0.95, for 2013–2014, but smaller, 0.90, for 2014–2015. SIE centered pattern root mean square difference (cRMSD) between the two datasets was around 10 × 10³ and 7 × 10³ km² for the two ice seasons (for SIT and SIV, statistics were not given). Sea ice data by the HBM-CMEMS-V4 have not yet been validated with in-situ measurements.

2.3.2. Baltic Sea—Sea Ice Concentration and Thickness Charts

The SIC dataset from the Baltic Sea—Sea Ice Concentration and Thickness Charts product is used in the HIGHTSI model runs, see Section 2.2. The SIC data are based on the gridded digitized Baltic ice chart produced daily by the ice analysts at FIS or SMHI, and tit has a 1 km pixel size [11].

2.4. In-Situ Data from Coastal Ice Stations

The Finnish Ice Service (FIS) conducts weekly in-situ SIT and SNT measurements on coastal ice stations along the Finnish coast, from the Bay of Bothnia to the Gulf of Finland. The stations are all located in the fast-ice regime, typically within the shelter of the archipelago. The total number of stations is 21, but not all of them are active during an ice season, depending mainly on the ice conditions. For the ice season 2019–2020, long time series of data (from 11 November 2019 to 11 May 2020) are available only from three stations in northern part of the Bay of Bothnia as the ice season was very mild. In addition, some data (from four to eight weekly measurements) are available from four other stations in the Bay of Bothnia. The locations of these stations are shown in Figure 7a. The in-situ weekly data include maximum and minimum SIT and SNT at a station location, and also max and min thickness of snow-ice. The difference between max and min SIT and SNT was up to 20 and 10 cm, respectively. Here, the max and min thickness values were averaged, but in the case the min SNT was 0 cm, the max value is used directly. The data are used to evaluate the BALFI SIT (from the CMEMS model) and SNT (from HIGHTSI) products accuracies.

2.5. Baltic Sea Ice Chart

The Baltic Sea ice chart is during the winter-season issued on a daily basis by the ice services at FMI or SMHI. The daily chart is published in the early afternoon, depicting the ice situation in the morning of the issuing date. The chart’s time stamp is 12:00 UTC, albeit indicative. An ice analyst uses the latest available Earth observation (EO) data, mainly SAR imagery, and in-situ data for producing the ice chart in a GIS production system [8]. The in-situ data originate from Finnish and Swedish icebreakers, coastal ice observation stations of the Baltic Sea ice services, drift buoys and other occasional sources such as merchant vessels, harbors and citizen observations.

The primary publication format of the ice chart is a PDF product, issued in co-production by FMI and SMHI on both institutes’ webpages. In addition to ice information, the chart includes information to winter navigation issued by the winter navigation authorities, such as assistance restrictions to all Baltic Sea ports and notifications to mariners. Figure 4 shows an example of the ice charts publicly available.

Ice chart information is also published in data format, both in the gridded and vector form. The parameters in the numerical outputs are Ice Type, Ice Concentration, Level Ice Thickness (minimum, mean, maximum), Degree of Deformation (DoD) and sea surface temperature. Symbol features in the visual output (PDF) are only saved in the Shape-like SIGRID-3 (Sea Ice Grid) data format.

The DoD parameter quantifies deformation-induced navigational difficulty and includes ridging information. It is a six-digit scale numeral: 0 (undeformed), 1 (rafted ice), 2 (slightly ridged ice), 3 (moderately ridged ice), 4 (heavily ridged ice) and 5 (brash, typically indicating a compacted ice edge zone of brash ice, connected to the brash ice barrier symbol). As all parameters are applied to one polygon, the DoD gives only one class for the whole polygon area. DoD also represent approximately the degree of ice ridging.

LFI extent from the ice chart is used here to validate the BALFI automated LFI extent from the SAR imagery. The ice chart LFI extent is not available in any CMEMS Baltic Sea product. The BALFI DoD product will be compared against the ice chart DoD data.

3. Baltic Sea Landfast Ice Products and BALFI Service

Here, algorithms and procedures for calculating the BALFI products are described. The products are: (a) LFI extent and SIT, (b) SNT on LFI and (c) degree of sea ice deformation (DoD) (qualitative index) on LFI. The format of the products is also presented. A schematic production chain for the BALFI products is shown in Figure 5. The production is fully automatic and operational. In the processing chain, the inputs (on the left in Figure 5) are the CMEMS Baltic Sea products, ECMWF NWP model data for the HIGHTSI model and S-1 and RS-2 SAR imagery. The starting point for the BALFI products is extraction of the LFI extent. It is defined as a combination of the 15 day temporal SAR image cross-correlation second smallest value (minimum proved to be more sensitive to random errors) between adjacent day SAR mosaics. SIT from the HBM-CMEMS-V4 model and SNT from the HIGHTSI model are mapped over the detected LFI. Using SAR imagery DoD is estimated for LFI. Details of the algorithms and procedures are given in Section 3.1, Section 3.2, Section 3.3 and Section 3.4 below.

The LFI products are operationally generated at FMI and distributed via the BALFI service web-portal [17]. On the web-portal, the products are visualized over a map layer, and users can zoom in and out in the map and see product data values at a selected geographical location. The BALFI service is further described in Section 3.5.

3.1. Landfast Ice Extent

Remote sensing methods for the LFI extent estimation were recently reviewed in [35]. In brief, SAR imagery allows the estimation of the LFI extent with high resolution, from ~10 m to a few hundred meters, regardless of cloud and daylight conditions which hamper the usage of the optical imagery. Mahoney et al. [36,37] detected LFI based on vector grayscale gradient fields of three subsequent SAR images. LFI grown all the way to the ocean bottom within brackish or freshwater Arctic river deltas (sea ice has very low salinity) can be identified based on the SAR backscatter magnitude because if there is no ice–water interface, the dielectric contrast at the ice bottom is significantly reduced [38]. Time series of sea ice drift data from SAR image pairs can be used to derive static ice areas, which are then be interpreted as LFI, assuming the time series of ice drift at a certain location is long enough. SAR ice drift estimation algorithms are typically based on either temporal cross-correlation—that is, cross-correlation between co-registered spatially (partly) overlapping SAR images acquired at different time instants [39]—temporal phase correlation [40] or optical flow [41]. In the phase correlation, only the phases of the two signals are correlated, thus making the approach less sensitive to differences in local signal magnitude variations. Optical flow by definition is the apparent pattern of motion of a scene caused by the relative motion between an observer and the scene. Assuming the motion is reasonably small, it can be estimated between two images over the same scene at two time instants based on the first-order Taylor expansion of the optical flow equation for image intensity at the two time instants. In [13], the cumulative Baltic Sea ice drift estimated from multi-temporal SAR imagery was used for locating LFI by indicating the areas where no ice motion had occurred within a predefined and long-enough time period (typically around two weeks). Temporal cross-correlation minimum was used to locate the Baltic LFI to aid SAR-based SIC estimation in [42]. SAR interferometry (InSAR) can also be used for LFI detection [43,44], as the phase difference is random for drift ice and coherent for the static ice fields. However, selection of the suitable temporal baseline is not straightforward as in order to have a high coherence over LFI, a short temporal baseline is needed, but on the other hand, high coherence dynamic by a long baseline are needed for separating LFI from the drift ice [44].

The LFI extent (i.e., area) along the Baltic Sea coast is estimated here by thresholding a temporal second smallest value of the cross-correlation between each pair of adjacent day HH and HV SAR mosaics for a 15 day period. The cross-correlation values at the same location (grid cell) corresponding to different time instants (the 14 adjacent day cross-correlations computed for the 15 day period) are ordered by their magnitude for the HH and HV channels separately, and then at each location the second smallest cross-correlation values are selected and channel-wise thresholds are applied. An example of the thresholding is depicted in Figure 6. This procedure follows [35] where the average cross-correlation was used instead. The results with the HH- and HV-polarization mosaics are combined with the logic OR function. Additionally, some morphological filtering is applied to smooth the LFI boundary [35]. The pixel size for the estimated LFI extent image is the same as for the SAR mosaic—500 m.

3.2. Sea Ice and Snow Thickness on Landfast Ice

For SIT and SNT estimation in the Arctic, many studies have been carried out using microwave radiometer data [45,46,47,48,49,50,51], radar and laser altimeter data [52,53,54,55] or with a combination of radar altimeter and radiometer data [56]. There are also operational products available. In addition, SAR imagery and scatterometry data have been used for estimating SNT or snow water equivalent (or their changes in time) on smooth landfast FYI in the Canadian Arctic Archipelago (CAA) [57,58,59]. This is based on a thermodynamic-brine-dielectric effect within the FYI snow cover and its effect on radar signatures [59,60]. In general, thicker snow insulates the basal snow layer from $T_a$ changes keeping the brine volume and dielectric properties constant, but thin snow is more thermally diffusive which leads to a larger change in brine volume and dielectric constant, and further in $\mathsf{\sigma }°$ [59].

In the Baltic Sea, SIT estimation with Soil Moisture and Ocean Salinity (SMOS) mission L-band radiometer was investigated in [61]. In this case study, the root-mean-square difference (RMSD) between the SMOS SIT and modal SIT from an airborne electromagnetic (EM) sounding instrument SIT data for 12 large circular areas (diameter from 36 to 66 km) was only 11 cm. However, it was concluded that the interpretation of SMOS-retrieved SIT ideally requires some knowledge on the shape of the SIT distribution within the SMOS footprint. Unfortunately, the spatial resolutions of various radiometers—e.g., AMSR2 and SMOS—are too coarse (~5 to 50 km) to map properties of the Baltic Sea LFI. In addition, the mixed-pixel effect—a mixture of brightness temperature signatures from land, islands, sea ice and ocean—would significantly deteriorate the quality of the retrievals. There are not yet any operational SIT products tailored for the Baltic Sea from radar or laser altimeter data. The Baltic Sea SIT estimation using radar or laser altimeter data is very challenging due to rather thin ice, and therefore, small freeboards, compared to the Arctic Ocean. SAR-based estimation of SNT or SWE for the Baltic landfast ice has not yet been studied. The Baltic Sea ice has very low salinity—from 0.2 to 2‰ depending on the location, time, and weather history [7]—and thus, there are not likely to be the same kind of observable changes on the landfast ice $\mathsf{\sigma }°$ due to the basal snow layer brine volume variation as in the CAA. The LFI extent and SIT and SNT products are obtained by mapping SIT from the HBM-CMEMS-V4 model and SNT from the HIGHTSI model over the detected LFI. This requires resampling of the coarser resolution SIT and SNT data to 500 m pixel size of the LFI extent image. The resolution of the SIT data is around 1 NM (i.e., 1852 m), and that of the SNT is 7.1 km by 5.9–8.6 km, depending on latitude. The resampling is conducted with bilinear interpolation, and in addition, along the coastline nearest neighbor extrapolation is needed to cover all the details of the coastline.

Figure 7 shows the LFI extent and SIT and SNT charts on 5 March 2020 when the Baltic Sea ice had largest extent during the ice season 2019–2020. This date was determined by FIS. The charts cover only the Bay of Bothnia as there was very little LFI elsewhere. The LFI extent can be visually compared to the that in the manual ice chart in Figure 4, and some differences in the extent can be seen—e.g., more LFT in the northeast part of the Bay of Bothnia in the SAR-based chart. Differences between the ice chart and SAR-based LFI extents will be studied more in Section 4.1.

3.3. Degree of Sea Ice Deformation on Landfast Ice

Currently, there are no automated sea ice products for the Arctic or Baltic Sea showing explicitly information on sea ice deformation or roughness. SAR-imagery-based sea ice classification systems typically perform classification to open water and different ice types, such as new ice, first-year-ice, multi-year-ice, or just open water vs. sea ice classification [62]. Some ice types may have also qualitative deformation categories, such as “smooth/level” and “rough/deformed”. In addition, algorithms for SIC estimation using SAR imagery have been developed [42,62]. Roughness of wintertime FYI in the CAA has been quantitatively linked to $\mathsf{\sigma }°$ in the high resolution (~3 to 10 m) C-band RS-2 and L-band ALOS-2 PALSAR-2 SAR imagery [63]. The sea ice surface roughness was measured with an airborne laser scanner. The $\mathsf{\sigma }°$ at all four linear polarizations from a C-band image at high ${\mathsf{\theta }}_0$ (45°) and from a L-band image at moderate ${\mathsf{\theta }}_0$ (28°) were strongly and positively correlated with the surface roughness ($r>$0.75). Retrieval models for the surface roughness were developed with the best case root-mean-square error (RMSE) of 5 mm. Studies with polarimetric SAR (POLSAR) and/or interferometric SAR (InSAR) data have also shown good results in estimation of sea ice deformation or roughness [64,65,66]. For any operational product showing quantitative information on ice deformation or roughness, either high-resolution SAR, POLSAR or InSAR data are very likely needed. A method for estimation of the degree of sea ice ridging (DIR) in the Baltic Sea using dual-polarized RS-2 SAR imagery was developed in [21]. DIR is characterized with four classes which follow the DoD classes in the FIS ice charts (see Section 2.5): DIR 1 is a combination of DoD classes 0 and 1, and DIR 2, 4 and 3 are equal to DoD 2, 3 and 4. The HH- and HV-polarization SAR image are first segmented. The sizes of the resulting segments exhibited high variation: the segment sizes varied from 20 to over 11,000 km²; the mean of them was about 1000 km². The goal is for the resulting segments to be mainly composed of one DIR class. For each segment, a feature vector was computed—that is, a set of SAR image texture features (in total 13) related to the ice ridging. The features include backscattering coefficients, entropy, autocorrelation, coefficient of variation, edge density, segment size and kurtosis. As an additional feature, the mean SIC for each segment was extracted from the FIS ice charts. The second stage classifies every segment using the feature vector and assigns one DIR class to each segment. The classification used a random forest classifier that also assesses the efficiency of each feature in the discrimination. The eight most informative segment-wise features were used in the classification. The DIR class for a feature vector was selected using majority vote within an ensemble of 200 classification trees. Training and testing DIR data were extracted from the FIS ice charts. The overall agreement between the ice-chart-based manual DIR and the automated DIR results varied monthly, from 63% to 83%. The SAR-based DIR estimation is only possible under cold conditions (dry snow cover).

We decided not to apply the SAR DIR estimation by Gegiuc et al. [21] in here as it is a segment-based approach, and over LFI the resulting segments are typically very large. We would like to obtain DoD or DIR estimates for LFI in fine scale, not greater than 1 km, so that data are well usable by the end-users. This DIR estimation method is also quite complex and slow to perform for operational usage, and it is currently limited to dry snow conditions.

Here, the DoD estimation was conducted using the cross-correlation between the HH- and HV-polarization SAR images. Use of this feature is based on a visual comparison between visual interpretation of DoD for a set of SAR images in January–April 2018 and several SAR texture features, including the cross-correlation, and ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ and ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$ images. Based on this visual analysis, we found that the cross-correlation corresponded the best to the visual DoD observations. The more deformed ice, the more distinguishable features (above the noise level) the SAR imagery has at both polarizations, resulting in an increasing inter-channel cross-correlation as a function of increasing DoD. The physical interpretation is that deformation structures (pressure ridges, rubble fields) relevant for people moving on LFI are characterized by moderate or large-scale surface roughness. These are typically present on both SAR channels as local fluctuations in $\mathsf{\sigma }°$ due to varying ice ridge structure (width, height, etc.) and orientation of ice blocks forming the ridges and other deformations. Higher cross-correlation values then indicate the areas where the ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ and ${\mathsf{\sigma }}_{HV}^{\mathrm{o}}$ fluctuation patterns are similar indicating deformed structures, but neglecting the areas with large small-scale surface roughness, likely to produce rather high ${\mathsf{\sigma }}_{HH}^{\mathrm{o}}$ without any significant spatial pattern. The cross-correlation value corresponding to 10% of ordered 15-day time series values at each grid cell was used as DoD—i.e., a temporal filtering was performed for each grid cell. The qualitative scale for DoD is from one (no deformation) to 249 (highly deformed). Like in the case of SIT and SNT, the estimated DoD was mapped over the detected LFI extent, see an example in Figure 7c.

Ship tracks were also visible as deformed (high pixel values) linear structures in the DoD product. However, ship tracks within deformed LFI typically get mixed with the deformed ice background.

It is noted that the DoD product is still experimental, and linking it to quantitative figures of sea ice deformation, such as height of ice ridges, would require a combined dataset of SAR imagery and sea ice roughness from an airborne 3-D laser scanner data.

3.4. Format of the BALFI Products

The BALFI products cover the whole Baltic Sea, including Kattegat. They are generated daily, and areas with new S-1 and RS-2 SAR image coverage are updated. The products are 8-bit GeoTiff-files, and give SIT and SNT with 1 cm resolution in the range of 0–255 cm, and DoD scaled from 1 (no deformation) to 249 (highly deformed ice). However, SIT values over 100 cm and SNT over 20–30 cm are rare. In all GeoTiff-files, pixel value 250 is land mask, and in the SIT and DoD charts, pixel value zero is open ocean, but in the SNT chart, it shows both open ocean and no snow on LFI. Thus, the SNT chart does not show the full extent of LFI.

The map data are given in a latitude–longitude coordinate system (converted from the Mercator projection of the FIS ice charts). The resolution of the product is about 500 m, and the pixel size in latitude–longitude is 0.0047° by 0.0094°. The latitude–longitude grids are included in the GeoTiff-files.

The formats and standards of the BALFI products follow those used by national Ice Services and CMEMS. The BALFI products in GeoTiff-format can be easily read, processed and visualized using freely available GIS tools (e.g., QGIS, Python with suitable GIS libraries).

3.5. BALFI Service

The BALFI service [17] is hosted at the Finnish National Satellite Data Centre (NSDC) operated and managed by FMI (see http://nsdc.fmi.fi/). The service portal is targeted for end-users with high-speed internet connection, and using the portal with desktop PCs, laptops, tablets and smartphones. Figure 8 shows a screenshot of the BALFI portal. The main section of the portal has a map and a selection tool for the BALFI products. One BALFI product is visualized at a time. Users can zoom in and out of the map, move around the map using a mouse, and see latitude–longitude coordinates at the mouse cursor. Clicking on the map produces a popup label with information about the product value at that particular location. The map window also shows fixed colormaps for the products—e.g., SIT product has a colormap for the 0–50 cm SIT range. On the left side of the map window are sub-sections for general information on the BALFI service and BALFI products, conditions for product use, etc. The BALFI products can be downloaded using OGC Web Coverage Service, and can also be visualized using WMS services.

In the map window, the raster format BALFI maps are shown over the vector format geographic dataset (i.e., map). This geographic dataset has much better spatial resolution than the 500-m BALFI charts, and therefore, near the coast there may be missing pixels (i.e., no-data areas) in the BALFI maps, or pixels slightly overlapping land. For small coastal areas without data in the BALFI charts, data from the nearest chart pixels can be used.

4. Results

We investigate here the accuracy and quality of the BALFI products for the ice season 2019–2020 using coastal ice station data (SIT and SNT) and LFI extent from the ice chart. Usage statistics of the BALFI service since its start in February 2019 are also presented.

4.1. Landfast Ice Extent

The SAR-based LFI extent (i.e., area) for the whole Baltic Sea was compared to manually derived extent in the ice chart in the time period from 10 December 2019 to 11 May 2020 (in total 154 days). Time series of the daily LFI extent from the two datasets is shown in Figure 9. The BALFI LFI extent was always larger than the ice chart LFI extent, except in late melting season in May 2020, but their temporal changes matched each other. The mean difference between the BALFI and ice chart LFI extents was around 1006 km², RMSD was 1323 km² and the correlation was 0.94. The mean difference and RMSD were rather large compared to the ice chart LFI extent variation in Figure 9, from roughly 3000 to 9000 km². The high correlation between the two datasets suggests that our automated LFI extent mapping retrieves LFI changes as good as the manual estimation for the ice charts.

The difference between the ice chart and BALFI LFI extents can be explained by the ice chart preparation process where the manual LFI analysis is based on a broader source of information, such as imagery other than SAR, but also on continuous monitoring and empirical knowledge of typical ice break-up zones. In the typical evolution of LFI, it tends to get expanded by consolidating drift ice to its existing edge. Depending on the weather development, the edge line evolves according to freezing and wind-induced break-up from the edge. The ice analysts tend to avoid overestimating the expansion of LFI in order to indicate probable break-off zones and that way to diminish risks caused by rapid day-to-day changes. We acknowledge this difference against the ice chart but also identify benefits in the automated Baltic Sea LFI extent mapping, especially in areas where ship navigation is rare and in shallow or detailed parts of the archipelagos. The BALFI LFI extent underestimation in late melting season is due to wet snow/ice which reduces greatly contrast and features in SAR imagery, leading to small cross-correlations between SAR image pairs.

Next, we looked into relative differences between the two LFI extent datasets on a daily basis. The following two statistics were calculated: a fraction of the pixels in the ice chart LFI extent covered by the BALFI LFI extent: $f_p=N_p^B/N_p^{IC}$, where $N_p^B$ is the number of the BALFI LFI pixels covering the ice chart LFI pixels $N_p^{IC}$, and a relative difference between the whole LFI extents ($(A_{LFI}^B-A_{LFI}^{IC})/A_{LFI}^{IC}$). The BALFI LFI extent covers nearly all of the ice chart LFI extent from the start of the ice season to mid-Apr 2020; $f_p$ was always over 95%, and the mean was 99%. During the late melting season, $f_p$ decreased continuously down to 62%. This disagreement can be again explained with deteriorating effects of the wet snow and ice conditions on the SAR imagery quality. However, good quality manual mapping of the LFI extent from SAR imagery in wet conditions is possible by experienced ice analysts. The relative difference of the LFI extents varied daily from −28% to +62%, and the average was +17%. Negative differences occurred only in late melting season, as Figure 9 shows. Here, how this difference varies from day-to-day is more important. There were sometimes large absolute changes, up to 40%, but 90% of them are at maximum only 10%. This shows that the difference in the total LFI extent between the two products does not typically have large day-to-day variations.

An example of the daily difference in the LFI extents in shown in Figure 10. The date selected for this example was 15 March 2020 when the LFI extent was the largest in the BALFI chart—around 10,850 km². The BALFI LFI extent was noticeably larger in the Bay of Bothnia than that in the ice chart. This figure demonstrates an observed deficiency in the automated LFI extent mapping; it does not sometimes detect LFI in smaller bays or when there is sea ice only up to few kilometers from the coast. In the latter case, sea ice is typically thin (<30 cm) and featureless in SAR imagery. However, in manual mapping for the ice chart, it can be better detected.

In summary, the BALFI LFI extent is always larger than ice chart LFI extent, except in late melting conditions. This difference can be explained by the ice chart preparation process. Typically, the BALFI LFI extent covers nearly all (>95%) of the LFI in the ice chart. The BALFI SAR-based LFI detection has a rather poor accuracy in the late melting conditions, and for detecting very narrow band (a few kms) LFI adjacent to open water.

4.2. Sea Ice and Snow Thickness on Landfast Ice

Next, we compared BALFI SIT (from the HBM-CMEMS-V4 model) and SNT (from the HIGHTSI model) data (denoted here as $h_i^B$ and $h_s^B$) to the in-situ data from the FIS coastal ice stations ($h_i^{IS}$ and $h_s^{IS}$). For the comparison, the $h_i^B$ and $h_s^B$ data (500 m pixel size) at the ice station (IS) locations were averaged within a 3 × 3 pixel block (1500 by 1500 m) to match the original resolution of the HBM-CMEMS-V4 SIT data, around 1 NM (1852 m).

First, we looked into long time series of SIT and SNT at the three coastal ice stations—see Figure 7a. Time series of SIT and SNT at the northernmost IS ”Röyttä” (65°45.8′N, 24°13.5′E) are shown in Figure 11. $h_i^{IS}$ was always significantly larger than $h_i^B$: $h_i^{IS}$ increased from 10 cm on 11 November 2019 at the start of the ice season (first in-situ measurement was possible) to 75 cm at the beginning of April 2020, and then to 42 cm on 11 May 2020 in the late melting season (the last measurement). The start of the $h_i^B$ time series was later, on 12 December 2019, with $h_i^B$ of 18 cm. The maximum $h_i^B$ also occurred in April 2020, but it was only 33 cm. The $h_i^{IS}$ data show the start of the LFI melting in late April, but the $h_i^B$ data around 1.5 weeks later. The difference between $h_i^B$ and $h_i^{IS}$ (i.e., bias) varied from −43 to −9 cm, and the average was −27 cm. Despite the large $h_i^B-h_i^{IS}$ $r$ was rather high, 0.88, showing similar temporal behavior between the two datasets. We suggest that there are two possible reasons for the large $h_i^B-h_i^{IS}$. First, the Röyttä IS location did not represent average ice conditions in the area, or the HBM-CMEMS-V4 just failed to predict the correct SIT evolution at the Röyttä IS.

The $h_s^B$ time series showed only small temporal variation, and the max $h_s^B$ was only 8 cm. $h_s^{IS}$ was typically slightly larger, and had the max value of 20 cm. The difference $h_s^B-h_s^{IS}$ had a high variation, from −16 to +3 cm, and $r$ between $h_s^B$ and $h_s^{IS}$ was only 0.50. However, when $h_s^{IS}$ was very small, below 5 cm, then $h_s^B$ also showed equally small values. Differences between $h_s^B$ and $h_s^{IS}$ can be possibly explained as in the case of $h_i^B$ and $h_i^{IS}$, and in addition, the effects of snowdrift dynamics are currently missing in the HIGHTSI model.

At the Kemi Ajos IS (65°38.9′N, 24°30.3′E), which is also very close to the land and only 21 km from the Röyttä IS, a similar underestimation of $h_i^{IS}$ by $h_i^B$ was observed—e.g., max $h_i^{IS}$ was 65 cm, but the corresponding $h_i^B$ was only 26 cm. Here, the $h_i^{IS}$ data showed ice-free conditions in May 2020, but $h_i^B$ was still 10–20 cm. For the $h_s^B$ and $h_s^{IS}$ datasets, the same kind of differences and equivalences were present as for the Röyttä case.

The time series data at the third IS “Hailuoto” (64°56.0′N, 24°40.0′E) are shown in Figure 12. This station was much farther away from the land, around 4 km, than the other two stations. Here, $h_i^B$ and $h_i^{IS}$ were close to each other up to mid-February 2020, and after that, $h_i^{IS}$ was significantly larger. When $h_i^{IS}$ had its max value of 80 cm, then $h_i^B$ was only 20 cm. The mean $h_i^B-h_i^{IS}$ was then smaller, −14 cm, than in the Röyttä IS data, as here, the datasets matched over half of the ice season. $r$ was only 0.50 as $h_i^B$ was nearly constant in March to mid-April when $h_i^{IS}$ had large changes. A noticeable disagreement between $h_i^B$ and $h_i^{IS}$ occurred on 24 February 2020 when $h_i^B$ was 0 cm—i.e., showed open water, but $h_i^{IS}$ was 23 cm. Around this date, the outer boundary of the SAR detected LFI extent was very close to the Hailuoto IS, and it had some daily variations. $h_s^{IS}$ showed only snow-free ice or very thin snow cover—max $h_s^{IS}$ was only 3 cm. $h_s^B$ was also mostly ≤3 cm. The largest $h_s^B$ (7 cm) occurred in the beginning of February 2020, but $h_s^{IS}$ showed snow-free ice. If the IS location was on a snow-free patch of sea ice, but on average, on a larger spatial scale, there was some snow on the ice, then this would explain the difference. Unfortunately, visual observations of the snow and sea ice conditions at the ISs were not conducted.

Next, all data from the seven ISs (see Section 2.4) were compared to the corresponding BALFI data. A scatterplot between $h_i^{IS}$ and $h_i^B$ is shown in Figure 13. The plot shows clearly the overall underestimation of $h_i^{IS}$ by $h_i^B$. A linear regression model between $h_i^{IS}$ by $h_i^B$ is:

$$h_i^B=0.40·h_i^{IS}+3.4 \mathrm{cm},$$

(4)

with the coefficient of determination ($r^2$) of 0.77. The number of samples for (4) was 73. On average, $h_i^B$ is only 0.4 times of $h_i^{IS}$, and the mean $h_i^B-h_i^{IS}$ is −16 cm. $h_s^B$ is typically smaller than $h_s^{IS}$; the slope coefficient of the linear regression is only 0.21, and the mean $h_s^B-h_s^{IS}$ −2 cm. The scatter between $h_s^B$ and $h_s^{IS}$ is very large; $r^2$ is only 0.29.

In general, the comparison between the IS and BALFI data showed that $h_i^B$ typically underestimates $h_i^{IS}$ significantly. Additionally, $h_s^B$ was smaller than $h_s^{IS}$ when $h_s^{IS}\ge$5 cm. These underestimations could be due to the ISs being mostly very close to land, which location may pose difficulties for the HBM-CMEMS-V4 and HIGHTSI models to predict correct $h_i^B$ and $h_s^B$ data, or the IS locations did not represent average ice conditions predicted by the models: 1 NM (1852 m) scale for the HBM-CMEMS-V4, and 6–9 km scale for the HIGHTSI. On the contrary to the results here, the Baltic Sea total SIE and SIV, and thus also SIT, were found to be typically overestimated by the HBM-CMEMS-V4 compared to the ice charts [34]. More in-situ LFI data, especially at locations away from the coastline, are needed to further investigate the BALFI SIT and SNT accuracies, and what kind of improvements the HBM-CMEMS-V4 and HIGHTSI models would possibly need. Feed-back and observations by end-users would also be valuable in further work.

4.3. Degree of Ice Deformation on Landfast Ice

The BALFI DoD product was compared to the DoD information in the Baltic ice chart. The ice chart DoD is determined manually, and it quantifies ice deformation-induced navigational difficulty with a numeral on a six-digit scale—for details, see Section 2.5. For many LFI ice chart polygons, DoD was zero here. The BALFI DoD is a qualitative figure from 1 to 249. In the comparison, only those ice chart polygons with DoD ≥1 were used. In total, there was 1424 such polygons. For these polygons, median BALFI DoD was calculated. The results are in

Table 2. Comparison between the degree of deformation (DoD) assigned manually to the ice chart polygons and the BALFI DoD. For each ice chart, polygon median BALFI DoD was calculated, and then their averages for each ice chart DoD class derived.

Ice Chart DoD	Number of Polygons	Average Median BALFI DoD
1	304	20
2	385	26
3	265	47
4	37	111
5	433	55

. The average median DoD increased with increasing ice chart DoD, if the DoD class 4, which has very small number of median DoD samples, was excluded. $r$ between the ice chart DoD and the average median BALFI DoD was rather high, 0.68, but $r$ between the ice chart DoD and single median DoD values was only 0.22. The small $r$ can be explained with a large difference in spatial resolutions—500 m for the BALFI DoD and typically from tens of to several hundred km² for the ice chart polygons. For further validation of the BALFI DoD, we would need quantitative data on sea ice ridging or roughness—e.g., from an airborne laser scanner.

4.4. Usage Statistics of the BALFI Service

During December 2019–May 2020, the BALFI service had in total 2572 users. The total number of sessions was 3589. In February–May 2019, when the service was operated for the first time, these figures were much higher—20,360 and 25,600, respectively. We assume that this was due to a very mild Baltic Sea ice season in 2019–2020, and thus, there was much less need for the BALFI service.

In 2019–2020, around 75% of the users were from Finland, 15.5% from Sweden and 4.7% from Estonia. In these Baltic Sea states, and also in Russia, the BALFI service is the most useful as their coasts always have some LFI during a Baltic Sea ice season. Interest in the BALFI service in Southern Baltic States (Denmark, Germany, Poland, Lithuania) depends on the occurrence of LFI—i.e., on the severity of the ice season. In the Southern Baltic Sea, there was at the last time a significant amount of LFI in 2012. Over half of all users (~60%) accessed the BALFI portal with mobile phones and tablets.

5. Discussion

Operational monitoring and forecasting (both in short and seasonal time scales) of sea-ice conditions is very important for winter navigation and offshore activities, and for monitoring of sea ice processes and trends. Sea ice information is currently made available by several international and national institutes—e.g., CMEMS, EUMETSAT OSI SAF and Ice Services at FMI and SMHI, and in many different formats. This information is typically not easy to access and use by the general public interested in landfast ice (LFI) conditions for recreational activities, such as skiing and ice fishing, or for operations to transport people or cargo over LFI (ice roads). Ice charts in image-format (see Figure 4) which are available to the public show easily decipherable overviews of sea ice conditions, but provide quite generalized levels of detail for LFI, and typically no data on snow thickness on sea ice (SNT). We developed here three operational products on the properties of the Baltic Sea LFI (extent and sea ice thickness (SIT), SNT and degree of ice deformation (DoD)) based on SAR imagery and existing products and models, and a web-portal (BALFI) for their visualization and usage. The BALFI portal is targeted to the general public, and for companies and authorities needing sea ice information in their operations. It is the first service of its kind for the Baltic Sea. It seems that there are not many similar services for the Arctic; one example is the Greenland Community Ice Information Service [67] which gives sea ice information around Qaanaaq (in the north-west coast of Greenland) for safe and efficient navigation when hunting or travelling. Another one is SIKU in Canada which is a mobile application and web platform by and for Inuit providing tools and services for ice safety, language preservation and weather [68]. The sea ice information includes satellite imagery, an ice roughness product and SIT and SNT data from in-situ sensors.

The BALFI service started operation in the middle of the Baltic Sea ice season in February 2019, and for 2019–2020, it was available for the whole ice season. The BALFI products for the 2019–2020 ice season were compared against the Baltic ice charts and in-situ SIT and SNT data from the FIS coastal ice stations, and their quality and accuracy were investigated. We did not have any other data available for the validation of the products.

The BALFI LFI extent was found to always be larger than ice chart LFI extent, except in late melting conditions. This difference can be explained by the ice chart preparation process where the ice analysts tend to avoid overestimating the expansion of the LFI by consolidating drift ice in order to indicate probable wind-induced break-off zones and that way to diminish risks caused by rapid day-to-day changes. It also identifies benefits in the automated BALFI LFI extent mapping, especially in areas where ship navigation is rare, and in shallow or detailed parts of the archipelagos. Unfortunately, the SAR-based LFI detection has currently a rather poor accuracy in the late ice melting conditions, and it sometimes does not detect very narrow band (a few kms) LFI adjacent to open water. The BALFI SIT and SNT data typically underestimated the in-situ data from the ice stations very close to the land. It is possible that this location is problematic for the HBM-CMEMS-V4 and HIGHTSI models to predict correct SIT and SNT, or the measurements at the stations did not represent average ice conditions predicted by the models. Against the ice chart DoD information, the BALFI DoD on the average increased with the increasing ice chart DoD, but the correlation was small. However, the ice chart DoD from polygons had much coarser spatial resolution than the BALFI DoD.

We have not yet received any feed-back by the BALFI end-users on the quality and accuracy of the BALFI products. Based on the validation studies of the BALFI products, we assume them to give usable information on the Baltic LFI properties for various end-users. However, further development and validation of the BALFI are needed. For example, we would need at least LFI extent estimations from high-resolution optical imagery, such as SENTINEL-2, and more in-situ LFI data, especially at locations away from the coast line. In addition, data on sea ice ridging or roughness from an airborne laser scanner, and SIT data from an airborne EM instrument would be highly valuable. A severe Baltic Sea ice season when LFI covers large areas, also in the Southern Baltic Sea, would significantly increase opportunities for product development and validation. SIT data from other models could be also used as input to the BALFI products—e.g., SIT from the HIGHTSI or the FMI’s Helsinki multicategory sea ice model (HELMI) [69,70]. Snowdrift dynamics is still missing in HIGHTSI, which affects the temporal variation of modelled SNT during the winter season. For development of the HIGHTSI model, SIT and SNT data will be available from a thermistor string-based snow and ice mass balance buoy (SIMBA) program initiated a few years ago at FMI. Two SIMBA buoys were deployed on LFI outside Hailuoto Island and Kalajoki on the Finnish coast. The measurement data are under quality control, and an automatic SIMBA algorithm [71] will be adapted for the Baltic Sea, and we expect SIMBA data will eventually be useful for the improvement of the HIGHTSI model operations.

Finally, we identified following topics for further development of the BALFI products and web-portal: (a) Identifying ice convergence (deformation) or divergence (opening) along the LFI boundary using SAR-based ice drift data. (b) Combination of the HBM-CMEMS-V4 and HIGHTSI SIT data, and possible data from other models, for a more accurate SIT estimate for LFI. (c) Showing shipping lane information in the BALFI products or in the vector map layer of the BALFI portal for safety of the people traveling on LFI. (d) Optimization of the BALFI web-portal for smartphones and using GPS-data for automatic geolocation when viewing the products. (e) Development of the SIMBA observation network (in-situ SIT and SNT data).

6. Conclusions

We developed operational products for the monitoring of the Baltic Sea LFI extent, sea ice and snow thickness and degree of deformation. The products are based on the S-1 and RS-2 SAR imagery, the HBM-CMEMS-V4 and HIGHTSI models and the SAR-based LFI extent mapping algorithm developed in [35]. They are generated daily and have a 500 m pixel size (the actual resolution is coarser for SIT and SNT). The products are visualized in the BALFI web-portal [17]. The BALFI service (products and portal) is targeted for the general public interested in the Baltic LFI conditions for recreational activities, such as skiing and ice fishing, for operations to transport people or goods over LFI (ice roads), and for authorities/institutions such as Ice Services and icebreaker management. The BALFI service was started in February 2019, and it is active during the Baltic Sea ice season, typically from November/December to May. It is the first service of its kind for the Baltic Sea, and there are also only few similar ones for the Arctic.

In December 2019–May 2020, the BALFI service had in total 2572 users. In February–May 2019, when it was operated for the first time, the number of users was much higher—20360. This difference was likely due to the very mild Baltic Sea ice season in 2019–2020 when there was less need for the BALFI service.

In this study, we also validated the BALFI products with available independent data (ice charts and in-situ coastal ice station data). We suggest that the current products give usable information on the Baltic LFI properties for various end-users, but their further development would be beneficial. Especially, the SAR-based LFI extent mapping has poor accuracy in late ice melting conditions, and sometimes narrow bands of LFI adjacent to open water are not detected. The BALFI SIT and SNT data typically underestimated the in-situ data from the ice stations. This could be possibly corrected by using also SIT data from other models, and by improving FMI’s HIGHTSI model. We also identified some topics for the further development of the BALFI service, and data requirements to enhance development and validation activities. A severe Baltic Sea ice season when LFI covers large areas would enhance these activities.

A similar service can be set up for local communities in the Arctic, if the required input data are available.