Fusion of Landsat 8 Operational Land Imager and Geostationary Ocean Color Imager for hourly monitoring surface morphology of lake ice with high resolution in Chagan Lake of Northeast China

The surface morphology of lake ice remarkably changes under the combined influence of thermal and mechanical forces. However, research on the surface morphology of lake ice and its interaction with climate is scarce. A large-scale linear structure has repeatedly appeared on satellite images of Chagan Lake in recent years. The Geostationary Ocean Color Imager (GOCI), with a 1 h revisit, and Landsat 8 Operational Land Imager (OLI), with a spatial resolution of 30 m, provide the possibility for the study of hourly changes in the large-scale linear structure. We merged the Landsat and GOCI images, using an Enhanced Spatial and Temporal Adaptive Reflectance Fusion Model (ESTARFM), and extracted the lengths and angles of the linear structure. We monitored the hourly changes in the surface morphology during the cold season from 2018 to 2019. The average length of the linear structure in the completely frozen period was 21 141.57 ± 68.36 m. The average azimuth angle was 335.48 ± 0.23^∘, nearly perpendicular to the domain wind in winter. Through two field investigations during the two recent cold seasons, we verified the linear structure as being ice fractures and ridges. The evolution of surface morphology is closely associated with air temperature, wind, and shoreline geometry.

1 Introduction

Lake ice is one of the essential climate variables in the cryosphere (Bojinski et al., 2014) and is closely associated with lake environments, ecological regulation, public transportation, and the safety of human activities (Hampton et al., 2017; Magnuson et al., 2000; Leppäranta, 2015; Brown and Duguay, 2010; Arp et al., 2020). The shortening of the ice cover duration and the thinning of ice thickness have been common trends throughout the world (IPCC, 2021; SROCC, 2019; Murfitt and Duguay, 2021). Recent work using remote sensing mainly focused on lake ice phenology (Weber et al., 2016; Zhang et al., 2021; Xie et al., 2020; Murfitt and Duguay, 2020; Du et al., 2017), lake ice classification (Hoekstra et al., 2020), ice thickness (Murfitt et al., 2018b; Kang et al., 2014; Gogineni and Yan, 2015), and ice albedo (Li et al., 2018; Lang et al., 2018). However, previous work on the surface morphology of lake ice is scarce. The surface morphology, i.e., ice ridges and fractures, is controlled by the dynamic processes of lake ice, which have attracted widespread concern in academia and the society. In this study, we monitored the surface morphology of Chagan Lake in Northeast China by combining high spatiotemporal remote sensing data and the results of field investigations and explored the potential influences of climate factors.

Satellite remote sensing is macroscopic, multi-source, and wide-ranging and has been successfully applied in the global remote sensing monitoring of lake ice (Murfitt and Duguay, 2021; Doernhoefer and Oppelt, 2016; Du et al., 2019). Visible light and multispectral data were first used to monitor lake ice based on the spectral difference between ice and water; however, the best period to monitor lake ice changes is prone to be missed due to the influences of clouds, fog, and light (Howell et al., 2009; Cai et al., 2019; Yang et al., 2019; Qi et al., 2020). Active microwave data are used to identify ice through the differences in the backscatter of water and ice, while passive microwave data are used to identify ice through the differences in brightness temperature (Cai et al., 2017). Microwave remote sensing can penetrate the dry snow on the surface of lake ice and observe the internal structure and stratification of lake ice (Jones et al., 2013) from the early stage of the qualitative differentiation between ground ice and floating ice to the quantitative inversion of the phenology and thickness of lake ice (Ke et al., 2013; Jeffries et al., 2013; Howell et al., 2009; Kang et al., 2014). Although the temporal resolution of active microwave remote sensing data has been improved from 30 d (ERS – European Remote Sensing) to daily return visits (Radarsat-2), the optimized technique is too costly and more suitable for case studies of small or medium lakes (Murfitt et al., 2018a; Geldsetzer et al., 2010). With a high temporal resolution and temporal coverage, passive remote sensing data can detect ice cover under all weather conditions and are limited by their low spatial resolution and significant mixed image effects. Thus, passive remote sensing is more suitable for monitoring large-scale lake ice (Du et al., 2017; Qiu et al., 2018). Although multi-source remote sensing is available to monitor lake ice processes, single-sensor remote sensing data cannot simultaneously achieve both accurate remote sensing monitoring and a high frequency.

The growth and decay processes of lake ice change very fast, requiring high temporal resolution to capture surface morphology. Satellite sensors with moderate spatial resolution, such as Visible Infrared Imaging Radiometer (VIIRS) and Moderate Resolution Imaging Spectroradiometer (MODIS), can monitor the temporal changes in lake ice daily but fail to reflect the spatial details of surface morphology. Satellite sensors with medium spatial resolution, such as Landsat and Sentinel, can provide fine-texture images but do not have frequent images to capture fast changes. Fusion methods include the unmixing method, weight function method, and dictionary–pair learning method (Sisheber et al., 2022; Zhu et al., 2016). The most common weight function method includes the Spatial and Temporal Adaptive Reflectance Fusion Model (STARFM; Feng et al., 2006), Spatial and Temporal Adaptive Algorithm (STAARCH; Hilker et al., 2009), Enhanced Spatial And Temporal Adaptive Reflectance Fusion Model (ESTARFM; Zhu et al., 2010), and Flexible Spatiotemporal Data Fusion (FSDAF; Zhu et al., 2016). Previous studies have proven that these spatiotemporal fusion methods can improve the monitoring abilities of remote sensing for specific applications, but that they fail to monitor the abrupt changes in landscapes and spectral differences. All these models are derived by pairs of coarse and fine resolution, e.g., one pair for the STARFM and two pairs for the ESTARFM. The ESTARFM method performs better than the STARFM in heterogeneous landscapes (Zhu et al., 2010; Knauer et al., 2016; Y. Wang et al., 2021; Jarihani et al., 2014). When monitoring spatial changes, the STAARCH method strictly requires two pairs of bases, including one before and one after the changes, which limits its wide application. The FSDAF is more robust than the other three methods but has limitations in detecting tiny changes (Zhu et al., 2016), making it difficult to monitor the changes in surface morphology. Therefore, we generated fusion images with a high spatiotemporal resolution based on the ESTARFM for further exploration.

The evolution of a lake ice season is mainly a thermodynamic process influenced by thermal and mechanical forces. Thermal forces enable the surface to melt and freeze at the turn of the day and night, and the mechanical strength of winds and currents makes the ice bulk move and collide, causing water courses, ice ridges, and ice fractures to appear, develop, and disappear. Therefore, the surface morphology of lake ice exhibits a periodical spatiotemporal difference, which differs significantly from the flat and smooth surface of lake ice. The horizontal and linear structures of lake ice are monitored by optical satellites for large lakes in Europe and have been explained by ice displacement (Leppäranta, 2015). High-resolution satellite–airborne synthetic-aperture radar (SAR) images have been used to monitor the surface deformation of sea ice, such as ice ridges (Dierking, 2010). Moreover, airborne aerial platforms, such as unmanned aerial vehicles (UAVs; Li et al., 2020), airborne radar (Jeffries et al., 2013), and ground-penetrating radar (Gusmeroli and Grosse, 2012), have effectively complemented remote sensing data sources to monitor the changes in the morphology of lake ice. The spatial distribution of the surface morphology of lake ice is complex, variable, and discontinuous, which is characterized by highlighting linear features on remote sensing images. Currently, there are few studies on the changes in the surface morphology of lake ice and the influence factors. It is meaningful to develop a quantitative method to describe the surface morphology and explore the potential influences.

This study proved the capability of high spatiotemporal remote sensing images for monitoring the surface morphology of lake ice in Chagan Lake, Northeast China. Our work aimed to (1) generate high spatiotemporal satellite images using Landsat and the Geostationary Ocean Color Imager (GOCI), (2) monitor the hourly spatial changes in the surface morphology, including the length and angle, of Chagan Lake, and (3) discuss the beneficial climate conditions during the formation of the surface morphology of lake ice.

Figure 1 The spatial distribution of Chagan Lake and field photographs. Panel (a) is provided by the Landsat 8 Operational Land Imager (OLI) on 10 February 2019, with the band composite of R(5) + G(4) + B(3). Panels (b)–(h) display the field photographs captured in field investigations.

2 Materials

2.1 Study area

As one of the 10 largest lakes in China, Chagan Lake (124^∘03^′–124^∘34^′ E, 45^∘09^′–45^∘30^′ N; Fig. 1) plays an essential role in fisheries, agricultural irrigation, and winter recreation in the surrounding areas (Wen et al., 2020). The average and maximum water depths are 2.5 and 4.5 m, respectively (Duan et al., 2007; Song et al., 2011). The lake has a water area of 329.72 km² and a perimeter of 201.03 km, according to the Landsat 8 Operational Land Imager (OLI), on 10 January 2019. The salinity of lake water ranges from 0.31 ‰ to 0.78 ‰ (Liu et al., 2020). The catchment of Chagan Lake is characterized by semi-arid and sub-humid continental monsoons, with air temperature, precipitation, and evaporation of 5.5 ^∘C, 430 mm, and 1496 mm, respectively (Song et al., 2011). The recharge sources mainly comprise precipitation, groundwater, and adjacent irrigation discharge (Liu et al., 2019). Salinized soil farmland and grassland pastures are widely distributed in the catchment area. Chagan Lake is a typical lake with seasonal ice cover, and the ice cover exists from November to April each cold season, with the maximum ice thickness ranging from 0.8 to 1.1 m (Liu et al., 2020; Hao et al., 2021). We conducted two field investigations on 30–31 December 2020 and 2–4 January 2022 to verify the results from remote sensing. We measured the ice thicknesses using an electronic digital calliper with a resolution of 0.01 mm. Moreover, we measured the water depths using a handheld sonar detector (Speedtech SM-5) with a resolution of 0.1 m and compared the depths in fall (17 September 2021) and winter (2–4 January 2022).

Table 1 The comparison of band range between GOCI and Landsat. The selected bands for merging Landsat and GOCI are marked in bold.

2.2 Materials

2.2.1 GOCI

The GOCI is the first satellite for detecting ocean color from a geostationary orbit and has been widely applied in deriving optical, biological, and biogeochemical properties (Ryu et al., 2012; Ryu and Ishizaka, 2012). The GOCI data have been available since April 2004, covering about 2500 km × 2500 km around the Korean Peninsula. It has six visible bands and two near-infrared bands, and Table 1 provides the details of the band information. The GOCI provides observations every hour from 08:30 to 15:30 local time (LT), with a spatial resolution of 500 m. The GOCI has the most significant advantage of a high temporal resolution, with eight images daily, which can offer the details of the freeze-up and break-up processes. A total of 96 GOCI images during the cold season from 2018 to 2019 were used in this study (Table 2). The atmospheric correction was performed by GOCI data processing software (GDPS).

Table 2 The usage of GOCI images from November 2018 to March 2019.

2.2.2 Landsat

The Landsat 8 satellite was launched on 11 February 2013. It carries the Operational Land Imager (OLI) and the Thermal Infrared Sensor (TIRS), which have been widely used to monitor lake and river ice (X. Wang et al., 2021; Yang et al., 2020). The OLI has nine bands, with a spatial resolution of 30 m for bands 1–7 and band 9 and a spatial resolution of 15 m for band 8. It has a temporal resolution of 16 d. Six Landsat images taken during the cold season from 2018 to 2019 were prepared for data fusion with the GOCI. The capture dates were 6 November 2018, 22 November 2018, 8 December 2018, 26 February 2019, 14 March 2019, and 15 April 2019. The path and row are 119 and 29, respectively. We downloaded the Landsat-calibrated surface reflectance Tier 1 collection for Landsat 8 OLI (LANDSAT/LC08/C01/T1_SR) from the Google Earth Engine (GEE) for further work. Table 1 compares the band information on Landsat and GOCI. From Table 1, we found that the band range of band 3 from GOCI (480–500 nm) overlapped with that of band 2 from Landsat 8 OLI (450–515 nm). The waterbody had relatively strong reflectance in the blue band (400–480 nm), and the blue band images clearly displayed the linear structure. Therefore, we merged band 3 of the GOCI and band 2 of the Landsat 8 OLI and generated the 96 fusion images for further work.

2.2.3 Auxiliary data

The lake ice phenology of the cold season from 2018 to 2019 was extracted from the combined time series of the surface temperature of lake water provided by the MOD11A1 and MYD11A1 products (Song et al., 2016; Hao et al., 2021). The freeze-up date is defined as the first day on which the surface temperature of lake water is below 0 ^∘C in winter; the break-up date is defined as the first day on which the surface temperature of lake water is above 0 ^∘C in spring. We utilized the daily air temperatures, precipitation, wind directions, and wind speeds of the Qian'an station (ID 50948) to explain the influence of climate on lake ice from 2010 to 2021. Qian'an has a longitude and latitude of 124.011^∘ E and 44.998^∘ N, respectively, with an elevation of 146.3 m. In total, 16 directions with an interval of 22.5^∘ were used to describe the wind directions, including north (N), north–northeast (NNE), northeast (NE), east–northeast (ENE), east (E), east–southeast (ESE), southeast (SE), south–southeast (SSE), south (S), south–southwest (SSW), southwest (SW), west–southwest (WSW), west (W), west–northwest (WNW), northwest (NW), and north–northwest (NNW). The climate records were used to explain the relationship between lake ice and climate.

Figure 2 The workflow of this study.

3 Methods

3.1 The framework of methodology

Figure 2 presents the flowchart of our work. We preprocessed the Landsat 8 OLI and GOCI and prepared the reflectance images of band 2 of Landsat 8 OLI and band 3 of GOCI. Then, we merged GOCI and Landsat 8 OLI using the ESTARFM and generated new fusion data with a spatial and temporal resolution of 30 m and 1 h, respectively. After that, the geographic location of the linear structure on the surface of the lake ice was identified, and the morphological parameters were extracted, including lengths and angles.

3.2 The ESTARFM fusion

The ESTARFM, in which two pairs of Landsat and GOCI images were used to generate spatiotemporal fusion data, was proposed by Zhu et al. (2010) and based on the STARFM. First, the coarse GOCI data were projected and resampled to a fine Landsat image at two known times t_m and t_n. Second, similar neighborhood pixels were searched with a moving window by setting spectral differences. Third, we calculated the normalized weight of each similar pixel by considering the spatial, spectral, and temporal differences. Then, the coarse GOCI values were transferred to fine Landsat data using the pixel-based conversion coefficients in the linear regression. Finally, the coarse GOCI data at the same time were used to calculate the fine fusion data at the predicted time (t_p), which is expressed as follows (Liu et al., 2021; Bai et al., 2017; Zhu et al., 2010):

$$\begin{array}{}\text{(1)}& \begin{array}{rl}{L}_b\left(x_{w/\mathrm{2}},y_{w/\mathrm{2}},t_{\mathrm{p}}\right)& =T_m\times L_{bm}\left(x_{w/\mathrm{2}},y_{w/\mathrm{2}},t_{\mathrm{p}}\right)\\ & +T_n\times L_{bn}\left(x_{w/\mathrm{2}},y_{w/\mathrm{2}},t_{\mathrm{p}}\right),\end{array}\end{array}$$

where $L_b(x_{w/\mathrm{2}},y_{w/\mathrm{2}},t_{\mathrm{p}})$ is the final predicted fine-resolution reflectance at the prediction time t_p. w represents the size of the moving window, and the corresponding center is $(x_{w/\mathrm{2}},y_{w/\mathrm{2}})$. $L_{bk}(x_{w/\mathrm{2}},y_{w/\mathrm{2}},t_{\mathrm{p}})$ is the fine-resolution reflectance at t_k (k=m or n) at the base date. T_k is the time weight, calculated from the magnitude of the detected change in the reflectance of the coarse spatial-resolution image between t_m and t_n and the prediction moment t_p.

$$\begin{array}{}\text{(2)}& T_k={\displaystyle \frac{\mathrm{1}/\left|\sum _{j=\mathrm{1}}^w\sum _{i=\mathrm{1}}^{w}C\left(x_j,y_i,t_k\right)-\sum _{j=\mathrm{1}}^w\sum _{i=\mathrm{1}}^{w}C\left(x_j,y_i,t_{\mathrm{p}}\right)\right|}{\sum _{k=m,n}\left(\mathrm{1}/\left|\sum _{j=\mathrm{1}}^w\sum _{i=\mathrm{1}}^{w}C\left(x_j,y_i,t_k\right)-\sum _{j=\mathrm{1}}^w\sum _{i=\mathrm{1}}^{w}C\left(x_j,y_i,t_{\mathrm{p}}\right)\right|\right)}},\left(k=m,n\right),\end{array}$$

where C(x_j, y_i, t_k) and $C(x_j,y_i,t_k)$ denote the image element values of similar image elements (x_i,y_j) within the moving window of the coarse spatial-resolution image at the reference moment t_k and prediction moment t_p, respectively.

3.3 The quantitative analysis of linear structures

In the beginning, the Landsat-GOCI fusion images were transformed into binary images. We extracted the original linear network with the Canny (1986) edge detection algorithm and then conducted edge detection to remove the outer boundaries. The morphological processing, including opening, filling, and eroding sequentially, was implemented for the inner part of the linear network. Then, the linear structure was derived from the largest connected domain of the linear network without boundaries, and the length is calculated by the shortest path of the largest connected domain. We connected the northernmost and southernmost ends into a straight line. The angle followed the definition of the wind direction above. We compared the auto-extraction and visual interpretation in our previous work (Hao et al., 2021). The R² values of the length and angle of 0.96 and 0.98 proved the good performance of the auto-extraction algorithm.

Figure 3 The actual image observed on 22 November 2018 (a) and its prediction images by the ESTARFM (b). Panels (c) and (d) display the magnified version of the red rectangle in panels (a) and (b), respectively.

Figure 4 Scatterplot of the real and the predicted reflectance by the ESTARFM for the blue band. The capture date was 22 November 2018.

4 Results

4.1 The performance of ESTARFM

We predicted the fine images from two pairs of fine Landsat and coarse GOCI data to fill the data gap caused by the low revisit frequency of Landsat. The two known pairs of data in the freeze-up process were captured on 6 November and 8 December 2018, and 53 fine ESTARFM fusion images were predicted from the coarse GOCI images. The two known pairs of data of the break-up process were captured on 26 February and 15 April 2019, and 43 fine ESTARFM fusion images were predicted. Figure 3 compares the spatial distribution of the original images and predicted images on 22 November 2018. In the predicted images, the texture of the ground objects was maintained, and the magnified versions of the figures in Fig. 3c and d clearly display the distribution of the linear structure. The predicted images were consistent with the original images and indicated the good fusion effect of ESTARFM. Figure 4 illustrates the scatterplots of the actual and predicted reflectance values along the 1:1 line. The R² value is 0.935, indicating that the predicted image was highly correlated with the actual image. The ranges of predicated and actual images were consistent; their mean reflectance values were both 0.10±0.03. The performance of the ESTARFM results was limited by (1) the limited image pairs available during the cold season from 2018 to 2019, (2) the time lag between the predicted and actual images, and (3) the inconsistency of the capture time between the predicted images and two pairs of input images (Lu et al., 2019; M. Liu et al., 2018). Therefore, the ESTARFM fusion images had a good performance and can provide reliable materials for further exploration.

Figure 5 The spatiotemporal changes in the linear structure for the fusion images of Lake Chagan during the cold season of 2018–2019.

Figure 6 The changes in the ice ridges during the cold season of 2018–2019. (a) Length changes during the growth process from 22 to 30 November 2018, as measured from 53 ESTARFM fusion images. (b) Length changes during the recession process from 15 to 24 March 2019, as measured from 43 ESTARFM fusion images. (c) The daily average length. (d) The angles of ice ridges in different stages.

4.2 The changes in surface morphology

We extracted the surface morphology of Chagan Lake from 96 fusion images during the cold season from 2018 to 2019. Figure 5 displays the spatial changes in the linear structure of the Landsat images in the freeze-up and break-up processes. Figures A1 and A2 present the original images of the GOCI, with a resolution of 500 m, the fusion images of Landsat and GOCI, with a spatial resolution of 30 m, and the network structure and the linear structure in the freeze-up and break-up processes, thus providing more details of the extraction process. The linear structure appeared on images from southeast to northwest, lasting from 22 to 30 November 2018. The linear structure disappeared from northwest to southeast, lasting from 15 to 24 March 2019. Figure 6 shows the average daily lengths of the ice ridges, based on the Landsat and GOCI remote sensing data. We monitored the growth and recession process of the linear structure via 96 Landsat-GOCI fusion images and monitored the stable process via four Landsat images. The growth stage lasted for 9 d, from 22 to 30 November 2018. The ice ridges had a length range of 5211.17–18 042.15 m and an average value of 12 680.32±4472.37 m, extending from southeast to northwest. The azimuth angles of the ice ridges in the growth stage ranged from 331.54 to 338.17^∘, with an average value of 334.38^∘ ± 2.08^∘. The lengths in the stable process ranged from 21 052.78 to 21 227.53 m, with an average value of 21 141.57±68.26 m, and the angles changed from 335.15 to 335.77^∘, with an average value of 335.48±0.20^∘. The recession stage lasted for 10 d, from 15 to 24 March 2019. The ice ridges had a length range of 19 178.18–5924.03 m and an average value of 13 288.59±4907.89 m, disappearing from northwest to southeast. The azimuth angles of the ice ridges in the recession stage ranged from 329.84 to 336.16^∘, with an average value of 332.90^∘ ± 2.54^∘. The changing rates in the growth and recession stages were 1425.66 and 1325.42 m per day, respectively, indicating that growth was slightly faster than recession. The large-scale structure extending from northwest to southeast has repeatedly appeared on Landsat images since 1986, which has been reported in our previous work (Hao et al., 2021).

4.3 The field investigation

Considering the safety of traveling on ice, we conducted two field investigations during the two recent cold seasons, from 30 to 31 December 2020 and 2 to 4 January 2022, respectively. We divided the lake area into three regions according to the surface morphology of lake ice. Region 1 was distributed along the linear structures. The surface of lake ice is uneven, and ice fractures and ice ridges were widely distributed. Region 2 was distributed along the northeastern coast, where the Ice and Snow Fishing and Hunting Cultural Tourism Festival of Chagan Lake has been held at the end of December each year since 2001. Region 3 covered the southern part of Chagan Lake. The lake ice in Regions 2 and 3 was flat and smooth, and snow cover was sporadically distributed. The color of the lake ice along the linear structure was distinguished from lake ice in the neighborhood. We infer that the frozen time of the linear structure was later than the lake ice in the neighborhood. The difference between ice fractures and ice ridges was the vertical height. Ice ridges were elevated sections formed on the upper and lower surfaces of lake ice, consisting mainly of ridge sails and keels. We located 10 sampling points along the linear structure on the satellite images and collected field photos of ice ridges and fractures (Fig. 1). We further verified that the large-scale fractures on the images were made up of ice fractures and ridges.

Figure 7 The ice thickness (mm) and water depth (m) of Chagan Lake was measured during the periods from 2 to 4 January 2022.

We also measured the ice thicknesses and water depths of 16 sampling points (Fig. 7). The ice thicknesses in the winter of 2021 ranged from 437.55 to 668.25 mm, with an average value of 582.24±58.14 mm. The average ice thicknesses of Regions 1, 2, and 3 were 551.58, 547.75, and 645.74 mm, respectively. The average water depths of Regions 1, 2, and 3 were 3.48, 2.99, and 3.00 m, respectively. Among the three regions, Region 2 had the smallest average values of ice thickness and water depth. The differences in water depth between the fall of 2021 and the winner of 2021 had an average value of 0.12±0.05 m and a maximum value of 0.2 m. The water depth in winter was lower than that in fall, and the decreasing water level also was a cause of lake ice fracturing in winter (Leppäranta, 2015). The ice features first formed in the nearshore area of the southeastern coast, where the water depth was relatively smaller than that in other regions in Fig. 7. The ice thicknesses and water depths showed spatial coherence with the surface morphology.

Figure 8 The wind field of Chagan Lake during the cold season of 2018. (a) Growth process from 22 to 30 November 2018. (b) Stable process. (c) Recession process from 15 to 24 March 2019. (d) Daily average and maximum wind speed from 1 November 2018 to 15 April 2019.

4.4 The climate condition

Lake ice processes are governed by the complex interaction of hydraulics, thermodynamics, and mechanics. The heat loss due to the decreasing air temperature exceeds the heat gained from surface water in late fall and early winter. When the water temperature falls below the freezing point, the cooled water provides a beneficial condition for ice crystals. Then, the volume of the lake ice expands, and the amount increases, which is followed by the formation of ice. We analyzed the wind roses of the daily average and maximum wind speeds from 1 November 2018 to 15 April 2019 (Fig. 8). The freeze-up and break-up dates of Chagan Lake in the cold season of 2019 derived from MODIS daily land surface temperature (LST) products were 14 November 2018 and 24 March 2019, respectively. The ice ridges appeared on 22 November 2018, 5 d after the freeze-up date; they disappeared on 24 March 2019 and were consistent with the break-up date.

The domain wind in the growth process was NW, WNW, W, and WSW, the domain wind in the stable process was WNW, WSW, and NNE, and the domain wind in the recession process was NW, WNW, and WSW. The WNW and WSW direction was the domain wind direction for all three stages. The frequency in the WSW direction in the growth, stability, and recession stages was 22.22 %, 30 %, and 14.23 %, respectively. The angles between the WSW direction and ice ridges were 87.98, 86.88, and 85.40^∘ for the three stages. The nearly perpendicular relationship was consistent with our previous study (Hao et al., 2021). Our previous used the yearly average values of wind from 2013 and 2020, and the work herein just exploited the wind in the cold season from 2018 to 2019. Besides, the wind speed also contributed to the formation of ice cracks and ridges. The daily average wind speed in the growth process and recession process was 3.88 and 3.63 m s⁻¹, and the average values for the whole cold season were 2.97 m s⁻¹. The daily maximum wind speeds for both the growth process and recession process were 3.88, 3.63, 7.05, and 6.86 m s⁻¹, and the average value for the whole cold season was 5.65 m s⁻¹. Not only the daily average value but also the maximum value was higher than the average level and revealed that the changing processes of ice fractures and ridges require the relatively strong action of wind. Therefore, wind speeds and directions played crucial roles in the development of ice ridges.

5 Discussion

Lake ice experiences different stages during the freezing and thawing cycle, including the phases of ice crystals, frazil ice, nails, pancake ice, and ice layers (Leppäranta, 2015). Lake ice expands and contracts as the air temperature rises and drops during cold seasons. The temperature difference between night and day results in the thermal expansion and contract of lake ice, which differ significantly within a given lake. Furthermore, long and narrow cracks are generated and likely to evolve into ice ridges under pressure when lake ice bulks, collides, and piles up. The definitions of lake ice are limited by the view ranges of field measurements, and satellite remote sensing provides a new perspective for surface morphology in a larger-scale observation. The large-scale linear structure has been found on remote sensing images during the cold season from 2018 to 2019. Similar phenomena have also been found in lakes and reservoirs in Northeast China (X. Liu et al., 2018). In our previous work, we used four Landsat 8 OLI images to monitor the monthly changes due to the limitation of temporal resolution (Hao et al., 2021). In this study, we took advantage of the hourly revisit of the GOCI and generated 53 and 43 ESTARFM fusion images in the freeze-up and break-up processes, respectively. This makes it possible to explore the linear structure in detail. The recurrent large-scale linear structure was further verified as being ice fractures and ice ridges in the fieldwork. Determining the spatial scales of ice fractures and ice ridges is challenging work, when considering the data source. A UAV is suitable to monitor the lake ice fractures at small scales (0–100 m), and the satellite sensors are suitable to monitor the ice ridges at large scales (10–100 km); both are suitable to monitor the horizon changes in the lake ice surface.

Besides the thermal forces, the lake ice fractures and ridges are also a dynamic process under the control of mechanical forces. The wind above the ice covers and water currents beneath the ice covers force the shift in ice bulk (Tan et al., 2012). Wang et al. (2006) compared the mechanical changes in the leads and ice covers, based on modeling results and satellite monitoring (Wang et al., 2006; Leppäranta, 2010), and revealed the influence of winds on the drift of ice. In the freeze-up process, the winds and water currents can push the ice toward the shore, preventing ice covers from freezing; in the break-up process, the wind can break ice covers and accelerate melting. The lake ice fractures were mainly controlled by thermal forces, and the ice ridges were mainly controlled by mechanical forces. The ice ridges underwent three stages during the cold season of 2018–2019, in which the wind directions and speeds exhibited remarkable differences. The ice ridges grew from southeast to northwest, with an average direction of 334.38^∘ and decayed from northwest to southeast with an average direction of 332.90^∘. The WSW direction frequently happened in all three stages, revealing the crucial role of winds in the development of ice ridges. The direction of the ice ridges was nearly perpendicular to the WSW direction (247.5^∘), which followed the principal laws of mechanics. The air temperature created a cold environment for the ice cover to freeze, the wind provided a mechanical force for the ice bulk to shift, and ice ridges and ice fractures formed. In addition, the direction of the ice ridges had a similar shape to the southwestern shoreline, and the stable shoreline geometry could explain the recurrent ice ridges with a specific direction, which has been reported in previous studies (Leppäranta, 2015).

Linear structures are common natural phenomena on the surfaces of sea ice and lake ice and profoundly influence light transfer and ice ecology. Lake ice ridges alter surface roughness and light transfer and then contribute to the thickness and volume of ice. People in cold regions have skillfully taken advantage of frozen ice covers for fishing, food storage, and commercial transportation. The capacity and stability of floating ice can be evaluated by the ice thickness and the spatial distribution of ice fractures and ridges (Tan et al., 2012). Generally, 30 cm is the thickness suggested for safe human activities on the ice (Leppäranta, 2015). Ice fractures and ice ridges potentially threaten human activities. In field investigations, we measured the ice thicknesses along the linear structure when the ice covers were steady. The ice thicknesses along the ice ridges were supposed to be thinner than in other areas (Leppäranta, 2015), but no significant difference in ice thickness had been found in our field measurements. Thus, the surface morphology of lake ice would be a reliable sign of travel on the ice being dangerous. Besides, we monitored the horizon changes in lake ice ridges using optical satellite images but ignored the vertical heights of ice ridges, which we need to consider in future work.

6 Conclusion

We generated high spatiotemporal remote sensing data of Landsat and GOCI, using ESTARFM to fill the gap in the fine monitoring of lake ice dynamics in Lake Chagan. We compared the reflectance of the fusion images and the original images on 22 November 2018. The R² value between the actual images and predicted images is up to 0.935, indicating that the predicted images were highly correlated with the actual images. Moreover, the consistency of the texture of the ground objects was maintained between the predicted images and the original images. Therefore, the ESTARFM fusion images provided reliable materials for further exploration.

We calculated the lengths and the angles of the linear structure on the fusion images in the freeze-up and break-up processes during the cold season from 2018 to 2019. Based on the satellite images, the linear structure experienced growth, stability, and recession stages. The growth stage lasted for 9 d, ranging from 22 to 30 November 2018. The recession stage lasted for 10 d, ranging from 15 to 22 March 2019. From southeast to northwest, the linear structure was 5211.17 to 18 042.15 m long during growth and from northwest to southwest, it disappeared. The average length of the ice ridges in the completely frozen period was 21 141.57±68.36 m. The average azimuth angle was 335.48^∘ ± 0.23^∘.

We performed field investigations and verified the linear structure as ice fractures and ridges. The direction of the linear structure was nearly perpendicular to the southwesterly wind direction, which is the dominant wind direction in winter. The deformation of the surface morphology was related to the meteorological conditions before the freeze-up process, including winds and air temperatures. This work demonstrated the capability of monitoring large-scale surface morphology using multi-source remote sensing and has profound implications for traveling safety on ice and ice engineering. We also plan to extend our findings to other large lakes in China and fill the knowledge gap of the surface morphology of lake ice.

Appendix A

Figure A1 The changes in ice ridges during the freeze-up process from 22 to 30 November 2019. (a–f) The original images from GOCI. (g–l) The fusion images from Landsat and GOCI. (m–r) The network structure of surface morphology. (s–x) The surface morphology.

Figure A2 The changes in ice ridges during the break-up process from 15 to 24 March 2019. (a–e) The original images from GOCI. (f–j) The fusion images from Landsat and GOCI. (k–o) The network structure of surface morphology. (p–t) The surface morphology.