Using Remote Sensing Data to Parameterize Ice Jam Modeling for a Northern Inland Delta

The Slave River is a northern river in Canada, with ice being an important component of its flow regime for at least half of the year. During the spring breakup period, ice jams and ice-jam flooding can occur in the Slave River Delta, which is of benefit for the replenishment of moisture and sediment required to maintain the ecological integrity of the delta. To better understand the ice jam processes that lead to flooding, as well as the replenishment of the delta, the one-dimensional hydraulic river ice model RIVICE was implemented to simulate and explore ice jam formation in the Slave River Delta. Incoming ice volume, a crucial input parameter for RIVICE, was determined by the novel approach of using MODIS space-born remote sensing imagery. Space-borne and air-borne remote sensing data were used to parameterize the upstream ice volume available for ice jamming. Gauged data was used to complement modeling calibration and validation. HEC-RAS, another one-dimensional hydrodynamic model, was used to determine ice volumes required for equilibrium jams and the upper limit of ice volume that a jam can sustain, as well as being used as a threshold for the volumes estimated by the dynamic ice jam simulations using RIVICE. Parameter sensitivity analysis shows that morphological and hydraulic properties have great impacts on the ice jam length and water depth in the Slave River Delta.


Introduction
In cold regions, the ice regime of rivers can be divided into three main phases during the winter: river freeze-up, continuous solid ice cover, and ice cover breakup.With ice cover breakup, rubble ice in rivers can pose threats to infrastructure adjacent to the river via jamming and subsequent flooding.However, some ecosystems, particularly inland deltas (e.g., the Slave River Delta), rely on ice jam flooding for replenishment of moisture and sediment.One of the main contributions of flow in the Slave River is from the Peace River.Ice jam development can be separated into three phases: formation, extension, and release.Studying ice jam processes can help extend our understanding of flooding and drying behaviors in river deltas such as the Slave River Delta.Ice jams have rough undersides, which increase flow resistance leading to backwater effects.Within the ice jam, cohesion is a force due to additional freezing which can uphold the structural integrity of the ice jam.Beltaos applied the Rising Limb Analysis Method (RLAM) to analyze hydrodynamic forces in the lower Mackenzie River [1].He found that pre-breakup conditions, such as Cumulative Melting Degree Days, as well as the rate of flow increases, are important factors in ice jam formation and release.He substantiated these findings with a hydrodynamic model for simulating ice jam waveforms.
River ice jamming and subsequent flooding can significantly impact local aquatic ecosystems [2].On the other hand, the reduced frequency of ice-jam flooding can lead to drying trends and subsequent shifts in flora and fauna composition and succession.Beltaos stated the effects of flow regulation by the dam construction on the river and climate change factors with respect to the drying trends in the Peace-Athabasca Delta [3].Based on calculations under different flow conditions, he concluded that regulation is the dominant contributor to the lower frequency and severity of ice jam flooding.Prowse and Conly pointed out that ice-jam-induced backwater flowing into the Peace-Athabasca Delta [4] is a major source of that delta's water supply.Although the Slave River Delta is approximated 400 km downstream of the Peace-Athabasca Delta, and its flows are augmented by the Lake Athabasca and the Athabasca River systems, flow regulation may still affect the Slave River Delta [5].
Inland deltas are particularly susceptible to ice jamming and subsequent flooding.Examples include flooding along the Peace-Athabasca Delta [6], and the lower Red River in Manitoba [7].To understand the formation and development of ice jams and the flood hazard they induce, the river ice hydraulic model RIVICE was implemented and calibrated to determine the conditions required for ice jamming.RIVICE has proven to be an effective tool for predicting ice jam formation and its impacts on flooding [8].
Due to safety issues, it is difficult to carry out field measurements on ice jams.To obtain the information necessary for better understanding river ice processes, researchers have taken advantage of satellite remote sensing technology, such as MODIS and RADARSAT-2, to monitor and characterize river ice behavior [9].RADARSAT-2, which is a commercial Earth observation radar satellite, was launched in December 2007 to monitor the environment for natural resource management.It provides high-quality data of the earth's surface features with a finer resolution (2 m/8 m) but at longer revisit intervals (24 days) than MODIS [10].Due to the longer wavelength emitted and received by the sensor, an advantage of RADARSAT-2 is that image acquisition is not hindered by weather conditions.Lindenschmidt et al. used RADARSAT-2 to measure ice thickness along the Red River, and determined the relationship between snow depth, ice thickness, and backscattering signal [11].Due to the high penetration and resolution of RADARSAT-2 signals, it is possible to analyze ice types based on RADARSAT-2 backscatter.Chu and Lindenschmidt classified river ice types using RADARSAT-2 imagery [12].In another study of the Slave River, different breakup patterns were tracked along the river, using RADARSAT-2 satellite imagery [13].These researchers used space-borne remote sensing technology (RADARSAT-2 and MODIS) to determine ice extent and ice types.However, it is difficult to acquire ice cover data along the Slave River in consecutive days by RADARSAT-2.Hence, MODIS can be utilized to monitor daily ice changes instead of RADARSAT-2.MODIS (MODerate-resolution Imaging Spectroradiometer) regularly scans the surface of the earth, albeit with a coarse resolution, permitting the monitoring of daily changes in river ice, which is essential for ice cover breakup research [14].MODIS has been used to facilitate various river-ice-related research, including the determination of the extent of ice and prediction of ice-related hazards along the Susquehanna River [15,16], and ice breakup monitoring on the Mackenzie River [17].Chaouch et al. designed an automated algorithm to monitor river ice using MODIS data [15].In addition to river ice data, researchers have also obtained climatic, snow cover and surface temperature data from MODIS products [18][19][20].One drawback of MODIS is that it is an optical sensor dependent on cloud cover and daylight conditions.
Studies focused on monitoring and modeling ice jams in the Slave River Delta are sparse.Our research goal was to obtain more understanding of river ice formation and deformation processes in the Slave River Delta by focusing on the breakup period during which ice jam flooding is most prevalent.

Study Site
The Slave River flows approximately 434 km northward from the Peace River and Peace-Athabasca Delta in Alberta to Great Slave Lake.This study focuses mainly on the delta of the river, the Slave River Delta, which is located at the mouth of the river at Great Slave Lake.The study area of this research extends from the Jean River to Great Slave Lake, approximately 25 km along the most downstream part of the main channel of the river (Figure 1).The drop in elevation at this part of the river is very low, with a gradient of 0.00003 m/m.The Resdelta Channel is the main channel, from which four other tributary channels, the Middle, Steamboat, and Nagle channels, and the Jean River, branch through the Slave River Delta.The mean annual flow of the Slave River is 3400 m 3 /s [5].Regulation of flow by the W.A.C Bennett Dam in the upper reaches of the Peace River moderates the seasonal difference in the flow, typically resulting in increasingly low natural winter flows and decreasingly high natural summer flows.The ice-cover season in the Slave River Delta extends from November to May.Before 1980, the frequency of high-magnitude discharge was once every two years; since 1980, this frequency has been once every five years [21].This decreasing trend in flood frequency the past several decades shows the drier conditions in the Slave River delta.
Water 2017, 9, 306 3 of 18 most downstream part of the main channel of the river (Figure 1).The drop in elevation at this part of the river is very low, with a gradient of 0.00003 m/m.The Resdelta Channel is the main channel, from which four other tributary channels, the Middle, Steamboat, and Nagle channels, and the Jean River, branch through the Slave River Delta.The mean annual flow of the Slave River is 3400 m 3 /s [5].Regulation of flow by the W.A.C Bennett Dam in the upper reaches of the Peace River moderates the seasonal difference in the flow, typically resulting in increasingly low natural winter flows and decreasingly high natural summer flows.The ice-cover season in the Slave River Delta extends from November to May.Before 1980, the frequency of high-magnitude discharge was once every two years; since 1980, this frequency has been once every five years [21].This decreasing trend in flood frequency the past several decades shows the drier conditions in the Slave River delta.Geomorphological factors of the river system, such as the sinuosity, channel width, and slope, have a great influence on the river ice regime and, due to its high flows and low gradient, make the Slave River Delta more conducive to ice jam flooding than the rest of the river [22].

RIVICE Model
The hydraulic model, RIVICE, was used to simulate ice jamming in the Slave River Delta.RIVICE is a one-dimensional model based on the continuity and full dynamic wave equations of flow and momentum.Although there are some limitations in one-dimensional models, they are widely used and have proven to be reasonably accurate.RIVICE can be used to simulate ice processes such as frazil ice generation, ice lodgment, ice cover progression through juxtaposition, ice deposition, hanging dam formation, and telescoping leading to ice thickening [23].Geomorphological factors of the river system, such as the sinuosity, channel width, and slope, have a great influence on the river ice regime and, due to its high flows and low gradient, make the Slave River Delta more conducive to ice jam flooding than the rest of the river [22].

RIVICE Model
The hydraulic model, RIVICE, was used to simulate ice jamming in the Slave River Delta.RIVICE is a one-dimensional model based on the continuity and full dynamic wave equations of flow and momentum.Although there are some limitations in one-dimensional models, they are widely used and have proven to be reasonably accurate.RIVICE can be used to simulate ice processes such as frazil ice generation, ice lodgment, ice cover progression through juxtaposition, ice deposition, hanging dam formation, and telescoping leading to ice thickening [23].
The model setup was as follows: Fifty cross-sections of the Slave River spaced 500 m apart (Figure 1) were extracted from a bathymetric model, different parts of which were surveyed between 2013 and 2015 and used to set up the RIVICE model.Cross-sections were also available for the complete reach of the Slave River, approximately 10 km apart, as part of an extensive water surveying campaign carried out in 1980.This campaign included gauging water levels and flows along the main channel and side-channels of the delta.These data served as a basis for calibration and validation of open-water and ice-covered conditions of the flow regime.The gauging stations included the following: Jean River (07NC001), Nagle River (07NC002), and Resdelta Channel (07NC009) (Figure 1).Only water level and flow data from 1980 are available at these stations in the Delta.
MODIS data were used to determine ice volumes approximated by: Ice Volume = Average width of ice cover × Ice cover length × Average ice thickness (1) All three factors on the right hand of Equation ( 1) can cause errors in ice volume estimation, which is an input parameter in the RIVICE model.Overestimated ice volume may lead to higher backwater levels.To improve the accuracy of the ice thickness calculation and enhance the reliability of the RIVICE model, aerial photos were used to validate the MODIS ice cover estimation (shown in the Results section).
When simulating an ice jam, both the upstream discharge boundary and incoming ice volume data are necessary.For other years, only the discharge data from the Fitzgerald station (07NB001) (Figure 1) can be drawn upon for the ice jam modeling.This station, located in Alberta, has provided continuous daily flow data for the past 68 years [24].
Before running the model, boundary conditions need to be set, particularly the downstream water levels and the upstream discharge.After running RIVICE, simulated water level profiles of the river at different time steps are generated and compared to the gauge data.The fit can be fine-tuned by adjusting the parameter settings.For example, Manning's coefficient of bed roughness n bed (Table 1) is adjusted until the open water simulations are reasonably in line with recorded water levels.The n 8m (Table 1) coefficient can be adjusted to calibrate the roughness along the ice cover underside and match simulated ice-cover conditions with back-water staging observations.Parameters that define the condition and state of the ice regime were also included, such as water temperature, ice cover strength parameters, ice transport parameters and ice cover and flow characteristics.Initial values from the literature and other studies were used to set up the model, but were then adjusted and "fine-tuned" during the calibration process.The model run parameter settings from the Monte-Carlo analysis were chosen for the calibration, as this yielded the best fit between the simulation and observed data.
To quantify which factor influences the formation of ice jamming the most, a local parameter sensitivity analysis was carried out to calculate the sensitivity of different parameters on ice jam flooding.Local sensitivity (e) is defined as in [25]: in which P refers to the parameter, such as n bed , n 8m , etc. (Table 1); O 0 is the initial calibrated output variable; and O i is the output from a change in a parameter value.Sensitive parameters have higher impacts on the output variable and generate higher ∆O values.Some parameters such as inflowing ice volume (v ice ), downstream water level, and discharge (Q) are determined by measurements including remote sensing and direct gauge data.Thus, errors from remote sensing data processing may lead to the inaccuracy of the local sensitivity calculation.To improve the accuracy of the analysis, more than one remote sensing technology was utilized in this research.

Remote Sensing Dataset
RADARSAT-2 uses microwaves (radar) for ranging and distance measuring.The waves are transmitted obliquely to the normal surface of the earth and the sensor then receives the backscatter signals from features on the surface of the earth which reveal different properties and characteristics of the surface features [26].After processing these signals, a two-dimensional image can be constructed from which further analyses of the surface or sub-surface features can be made, which, in our case, is the ice cover of the river.Radar data were used to help determining the breakup time of the Slave River, and more importantly, to be compared with MODIS image, improving our understanding of the ice phenology.
MODIS is an optical remote sensing satellite sensor.In this study, MOD09GQ (MODIS Surface Reflectance) was applied to assess river ice progression.Since MODIS satellites provide products of daily scans of the earth's surface, it is possible to acquire a continuous dataset of fine temporal resolution, allowing better tracking of ice cover breakup behavior and potentially identifying factors that cause ice cover retraction [22].In the MODIS images of the river, brighter pixels represent ice whereas darker ones correspond to open water.The contrast between the light and dark pixels permits the determination of the length of the ice-covered river channel from which ice volume can be calculated.MODIS datasets from 2000 to 2015 were downloaded from the NASA website [27].Through pixel analyses of the MODIS images, ice cover information can be attained, which are crucial input data for the RIVICE model.ArcGIS [28] was used to pre-process the raw MODIS imagery data and obtain the ice-covered channel length.In ArcMap, the centerline and cross-sections were marked on the MODIS image from which the length of the ice-covered portion was measured.Using these, together with the average ice cover width and estimated ice thickness, ice volume was determined.

Ice Volume Calculation by Using MODIS
Various methods for ice volume evaluation have been used, including modeling, direct measurements and the combination of these two [19,29].We present a novel technique in which remote sensing datasets are applied to calculating the extent of the ice cover and, using an empirical equation to estimate ice thicknesses (h) as a function of both Cumulative Freezing Degree Days (CFDD) and the average channel width in the Slave River Delta, derive an ice volume for downstream ice jamming.After comparing many equations for ice thickness (details about two other equations can be found in [30]), Stefan's equation showed the highest consistency with gauged data.Although it is more accurate and reliable to measure the ice thickness using direct measurement or remote sensing technology, the Stefan equation can still reasonably estimate thermal ice growth, even if it does not consider snow depth, which can affect the accuracy of ice thickness calculations.However, while snow depth data were not available for the Slave River Delta, the Stefan equation still provided reasonable ice thickness estimates.RADARSAT-2 will be utilized for ice thickness estimation in future research to compliment the ice thickness calculation.Plots of the CFDD for all winters from 2000 to 2015 and the winter of the year 1980 is shown in Figure 2. The relationship between ice thickness and CFDD is h = c √ CFDD, where c is a constant.Thawing effects on ice volume can be accounted for by determining the cumulative average air temperature degrees surpassing −5 • C. Thus, thawing effects were used to attenuate CFDD in spring.
Water 2017, 9, 306 6 of 18 sensing technology, the Stefan equation can still reasonably estimate thermal ice growth, even if it does not consider snow depth, which can affect the accuracy of ice thickness calculations.However, while snow depth data were not available for the Slave River Delta, the Stefan equation still provided reasonable ice thickness estimates.RADARSAT-2 will be utilized for ice thickness estimation in future research to compliment the ice thickness calculation.Plots of the CFDD for all winters from 2000 to 2015 and the winter of the year 1980 is shown in Figure 2. The relationship between ice thickness and CFDD is h = c√CFDD, where c is a constant.Thawing effects on ice volume can be accounted for by determining the cumulative average air temperature degrees surpassing −5 °C.Thus, thawing effects were used to attenuate CFDD in spring.MOD09GQ band2 (wavelength: 841-876 nm, resolution: 250 m) datasets were used to establish time series to estimate the ice cover situation of the Slave River [31].Different pixel values in the MODIS images represent different backscatters to which threshold values are applied to distinguish water from ice.By creating the histogram, the distribution of pixel values for each image can be determined.Clear water provides little backscatter signal to the satellite since it absorbs longer wavelengths and the reflected shorter wavelengths tend to scatter in the air.Hence, longer wavelengths appear more sensitive to the difference between water and ice.Since MOD09GQ band2 is longer in wavelength than band1, it has been chosen for our ice phenology research.A heat balance calculation based on upstream river temperature may help determine how much ice is melting and moving respectively.This calculation requires more extensive field sampling and can be addressed in future work.Due to this characteristic of clear water, it appears dark on MODIS imagery.On the other hand, ice provides more reflectance and appears lighter than clear water on the MODIS imagery [32].This difference leads to the significantly different pixel values in the MODIS imagery.Pixel values, with only one feature variable (a pixel is either ice or water), are normally distributed.Theoretically, a mixture of two normal distributions forms a bimodal distribution.Figure 3 shows an example of the MODIS pixel-value distribution from 17 November 2002.On that day, both ice and water appear in the river, so that two distinct peaks (local maxima) representing ice and water can be found in the histogram, which is consistent with the bimodal distribution theory.There is a big difference between the pixel values of the ice and water, and a value separating the two establishes the threshold.On cloudy days, backscatter values are difficult or even impossible to attain.Hence ice volumes for those days have to be interpolated between clear-day images [12].The algorithm Harmonic Analysis of Time Series (HANTS) was used to replace cloud-contaminated pixels with MOD09GQ band2 (wavelength: 841-876 nm, resolution: 250 m) datasets were used to establish time series to estimate the ice cover situation of the Slave River [31].Different pixel values in the MODIS images represent different backscatters to which threshold values are applied to distinguish water from ice.By creating the histogram, the distribution of pixel values for each image can be determined.Clear water provides little backscatter signal to the satellite since it absorbs longer wavelengths and the reflected shorter wavelengths tend to scatter in the air.Hence, longer wavelengths appear more sensitive to the difference between water and ice.Since MOD09GQ band2 is longer in wavelength than band1, it has been chosen for our ice phenology research.A heat balance calculation based on upstream river temperature may help determine how much ice is melting and moving respectively.This calculation requires more extensive field sampling and can be addressed in future work.Due to this characteristic of clear water, it appears dark on MODIS imagery.On the other hand, ice provides more reflectance and appears lighter than clear water on the MODIS imagery [32].This difference leads to the significantly different pixel values in the MODIS imagery.Pixel values, with only one feature variable (a pixel is either ice or water), are normally distributed.Theoretically, a mixture of two normal distributions forms a bimodal distribution.Figure 3 shows an example of the MODIS pixel-value distribution from 17 November 2002.On that day, both ice and water appear in the river, so that two distinct peaks (local maxima) representing ice and water can be found in the histogram, which is consistent with the bimodal distribution theory.There is a big difference between the pixel values of the ice and water, and a value separating the two establishes the threshold.On cloudy days, backscatter values are difficult or even impossible to attain.Hence ice volumes for those days have to be interpolated between clear-day images [12].The algorithm Harmonic Analysis of Time Series (HANTS) was used to replace cloud-contaminated pixels with Fourier series values from MODIS images [33].After interpolation, pixel values change so that another threshold must be applied to the interpolated images.Applying the threshold to the MODIS image pixel-value distributions, percentages of water and ice can be measured as well as ice volumes.During the ice breakup period, the diminishing ice of the upstream reach is considered to be inflowing ice at the downstream ice jam.We do make allowance for the recession of the ice front by tracking the location of the ice front at a daily time resolution using the MODIS imagery.Ice motion velocities in the year 2000-2015 were approximately 48-100 km/day.Based on the ice front movement, inflowing ice can be tracked.One point of uncertainty is that not all of the ice downstream of the ice front may be moving; rather, some could be melting while stationary, due to warm inflow river water.However, it is very difficult or perhaps impossible to distinguish how much ice has moved from upstream to the delta, as opposed to just melting before the next image was acquired.A shorter time interval between image acquisitions would be necessary.Fourier series values from MODIS images [33].After interpolation, pixel values change so that another threshold must be applied to the interpolated images.Applying the threshold to the MODIS image pixel-value distributions, percentages of water and ice can be measured as well as ice volumes.
During the ice breakup period, the diminishing ice of the upstream reach is considered to be inflowing ice at the downstream ice jam.We do make allowance for the recession of the ice front by tracking the location of the ice front at a daily time resolution using the MODIS imagery.Ice motion velocities in the year 2000-2015 were approximately 48-100 km/day.Based on the ice front movement, inflowing ice can be tracked.One point of uncertainty is that not all of the ice downstream of the ice front may be moving; rather, some could be melting while stationary, due to warm inflow river water.However, it is very difficult or perhaps impossible to distinguish how much ice has moved from upstream to the delta, as opposed to just melting before the next image was acquired.A shorter time interval between image acquisitions would be necessary.Ice volume calculation processes are shown in Figure 4 Using the Geographic Information System (GIS) [28], the water surface area of the Slave River Delta was measured, which is 12.8 km 2 .Percentages of the ice cover on the Slave River Delta are calculated from MODIS datasets.An example of calculating ice percentages is shown in Figure 4.

Data Preparation for Model
Essential inputs to the RIVICE control files include upstream discharge boundaries, downstream water level boundaries, incoming ice volume, the location of ice lodgments, and cross-sectional data.The water level data can be obtained from gauge stations (Figure 5).Flow data of Resdelta Channel (07NC009) (Figure 1) and Fitzgerald station are shown in Figure 6 and were used as the upstream discharge boundary.The downstream water levels were acquired at Great Slave Lake, which is used as the boundary condition for modelling according to the daily water level data at Yellowknife Bay station (07SB001).Ice volume calculation processes are shown in Figure 4 Using the Geographic Information System (GIS) [28], the water surface area of the Slave River Delta was measured, which is 12.8 km 2 .Percentages of the ice cover on the Slave River Delta are calculated from MODIS datasets.An example of calculating ice percentages is shown in Figure 4.

Data Preparation for Model
Essential inputs to the RIVICE control files include upstream discharge boundaries, downstream water level boundaries, incoming ice volume, the location of ice lodgments, and cross-sectional data.The water level data can be obtained from gauge stations (Figure 5).Flow data of Resdelta Channel (07NC009) (Figure 1) and Fitzgerald station are shown in Figure 6 and were used as the upstream discharge boundary.The downstream water levels were acquired at Great Slave Lake, which is used as the boundary condition for modelling according to the daily water level data at Yellowknife Bay station (07SB001).

Ice Volume
The CFDD was determined for each winter from mean daily air temperature records between 1971 and 2015.In the equation, the constant c depends on many factors such as meteorological conditions, geomorphological parameters, and water quality.It is derived empirically as the slope of a linear relationship between measured ice thicknesses and √CFDD (Figure 7).Based on the field survey and relationship between ice thicknesses and CFDD at the Slave River Delta and its tributaries

Ice Volume
The CFDD was determined for each winter from mean daily air temperature records between 1971 and 2015.In the Stefan equation, the constant c depends on many factors such as meteorological conditions, geomorphological parameters, and water quality.It is derived empirically as the slope of a linear relationship between measured ice thicknesses and √CFDD (Figure 7).Based on the field survey and relationship between ice thicknesses and CFDD at the Slave River Delta and its tributaries

Ice Volume
The CFDD was determined for each winter from mean daily air temperature records between 1971 and 2015.In the Stefan equation, the constant c depends on many factors such as meteorological conditions, geomorphological parameters, and water quality.It is derived empirically as the slope of a linear relationship between measured ice thicknesses and √ CFDD (Figure 7).Based on the field survey and relationship between ice thicknesses and CFDD at the Slave River Delta and its tributaries [13], the constant c was estimated to be 0.021.The blue trendline shows the relationship between the ice thickness and √ CFDD of the Slave River near Fort Smith.The constant c equals 0.018.Hence, 0.02 has been chosen as the constant, and the trendline with a slope of 0.02 m/m is shown by the orange line in Figure 7.
Water 2017, 9, 306 10 of 18 [13], the constant c was estimated to be 0.021.The blue trendline shows the relationship between the ice thickness and √CFDD of the Slave River near Fort Smith.The constant c equals 0.018.Hence, 0.02 has been chosen as the constant, and the trendline with a slope of 0.02 m/m is shown by the orange line in Figure 7.After determining ice percentage, ice thickness, and the ice-covered water surface area, ice volumes were calculated.
It is difficult to determine the ice cover extent from MODIS images during cloud-covered periods, particularly for the breakup events of 2007 and 2010.The ice cover thresholds for these two years cannot be determined.Thus, ice volumes cannot be calculated for these two events.Aerial photos acquired along the Slave River on 29 April 2015 were used to validate the ice volume estimation result of MODIS (Figure 8).The RADARSAT-2 image acquired on 30 April was compared with MODIS data for the same day (Figure 9).The noticeable difference between images from these two satellites is partly due to the different acquisition times, other reasons can be found in [12].
In Figure 10, three ice volume clouds (scatters showing the relationship between ice volumes and water levels) are shown.The cluster of blue dots show the modeling results of HEC-RAS, and the red ones and green ones show RIVICE modeling results with and without thawing effects.Thawing comes into effect for days on which average temperatures are above −5 °C.Clearly, the ice volumes are lower when thawing is incorporated in the calculation.The RIVICE clusters yield lower backwater elevations than HEC-RAS, which is to be expected, because not every jam that is formed is expected to be a fully-evolved equilibrium jam as premised by HEC-RAS.In HEC-RAS, stresses acting on an ice jam are calculated based on the ice jam force balance equation, while the deposition and erosion of ice covers are not taken into account [23,34].In RIVICE, ice jam formation and development were considered in a dynamic balance process.
The incoming ice volume of 1.1 million m 3 will lead to backwater staging effects causing ice jam flooding, which replenishes moisture in the delta during ice-cover breakup.When the incoming ice volume is higher than 1.5 million m 3 , the backwater effects are substantial.After determining ice percentage, ice thickness, and the ice-covered water surface area, ice volumes were calculated.
It is difficult to determine the ice cover extent from MODIS images during cloud-covered periods, particularly for the breakup events of 2007 and 2010.The ice cover thresholds for these two years cannot be determined.Thus, ice volumes cannot be calculated for these two events.Aerial photos acquired along the Slave River on 29 April 2015 were used to validate the ice volume estimation result of MODIS (Figure 8).The RADARSAT-2 image acquired on 30 April was compared with MODIS data for the same day (Figure 9).The noticeable difference between images from these two satellites is partly due to the different acquisition times, other reasons can be found in [12].
In Figure 10, three ice volume clouds (scatters showing the relationship between ice volumes and water levels) are shown.The cluster of blue dots show the modeling results of HEC-RAS, and the red ones and green ones show RIVICE modeling results with and without thawing effects.Thawing comes into effect for days on which average temperatures are above −5 • C. Clearly, the ice volumes are lower when thawing is incorporated in the calculation.The RIVICE clusters yield lower backwater elevations than HEC-RAS, which is to be expected, because not every jam that is formed is expected to be a fully-evolved equilibrium jam as premised by HEC-RAS.In HEC-RAS, stresses acting on an ice jam are calculated based on the ice jam force balance equation, while the deposition and erosion of ice covers are not taken into account [23,34].In RIVICE, ice jam formation and development were considered in a dynamic balance process.
The incoming ice volume of 1.1 million m 3 will lead to backwater staging effects causing ice jam flooding, which replenishes moisture in the delta during ice-cover breakup.When the incoming ice volume is higher than 1.5 million m 3 , the backwater effects are substantial.

Model Calibration
Water levels from 1980 recorded during the open-water and ice-covered period were used to calibrate the RIVICE model.
Open water conditions: In 1980, during the open-water period, the daily flow reached a maximum of 4870 m 3 /s, on 13 June 1980 (Figure 6), when the water level elevation of Great Slave Lake was 156.43 m [24].These values were used as the boundary conditions in RIVICE.The change of water levels at the Jean River gauge (07NC001) occurred approximately one day earlier than at the Nagle Channel gauge (07NC002) and two days earlier than at the Resdelta River gauge (07NC009) (Figure 1).
Ice-covered conditions: RIVICE was also implemented to simulate ice-covered conditions of the main channel in the Slave River Delta.The water level of the Resdelta River gauge (07NC009) (Figure 1) on 31 January 1980, was used to calibrate the model.
Figure 12 shows good agreement between the modeling result and the gauge recording of Resdelta (07NC009).

Model Calibration
Water levels from 1980 recorded during the open-water and ice-covered period were used to calibrate the RIVICE model.
Open water conditions: In 1980, during the open-water period, the daily flow reached a maximum of 4870 m 3 /s, on 13 June 1980 (Figure 6), when the water level elevation of Great Slave Lake was 156.43 m [24].These values were used as the boundary conditions in RIVICE.The change of water levels at the Jean River gauge (07NC001) occurred approximately one day earlier than at the Nagle Channel gauge (07NC002) and two days earlier than at the Resdelta River gauge (07NC009) (Figure 1).
Ice-covered conditions: RIVICE was also implemented to simulate ice-covered conditions of the main channel in the Slave River Delta.The water level of the Resdelta River gauge (07NC009) (Figure 1) on 31 January 1980, was used to calibrate the model.
Figure 12 shows good agreement between the modeling result and the gauge recording of Resdelta (07NC009).The change of water levels at the Jean River gauge (07NC001) occurred approximately one day earlier than at the Nagle Channel gauge (07NC002) and two days earlier than at the Resdelta River gauge (07NC009) (Figure 1).
Ice-covered conditions: RIVICE was also implemented to simulate ice-covered conditions of the main channel in the Slave River Delta.The water level of the Resdelta River gauge (07NC009) (Figure 1) on 31 January 1980, was used to calibrate the model.
Figure 12 shows good agreement between the modeling result and the gauge recording of Resdelta (07NC009).Ice-jammed conditions: Simulating ice jams is more complicated than simulating open water and ice-covered conditions since more factors need to be used in RIVICE calculations to simulate ice jam formation and progression.Water level data at the Great Slave Lake (156.557m) and the Resdelta Channel (156.884m) on 27 April 1980 were used to set up and calibrate the model, respectively.After defining the incoming ice volume on ice jam, the water level profile can be determined (Figure 13).Gauged water levels and HEC-RAS modeling results are applied to calibrate the ice-jam modeling results.HEC-RAS is a one-dimensional hydrodynamic model which is mainly used for simulating steady and unsteady flow routing based on an implicit four-point finite difference scheme.The HEC-RAS model simulates only equilibrium ice jams, i.e., in steady state, not dynamically, as is the case in RIVICE.Hence, the HEC-RAS result can be seen as the maximum backwater level that can be determined if the ice jam formation is not limited by the supply of incoming rubble ice.Some of the gauges in the delta were installed during the winter, prior to ice-cover breakup.Ice jams regularly occur along this reach during the final breakup of the Slave River.The error between simulated values and observed ones is shown in Figure 14, which is acceptable.Ice-jammed conditions: Simulating ice jams is more complicated than simulating open water and ice-covered conditions since more factors need to be used in RIVICE calculations to simulate ice jam formation and progression.Water level data at the Great Slave Lake (156.557m) and the Resdelta Channel (156.884m) on 27 April 1980 were used to set up and calibrate the model, respectively.After defining the incoming ice volume on ice jam, the water level profile can be determined (Figure 13).Gauged water levels and HEC-RAS modeling results are applied to calibrate the ice-jam modeling results.HEC-RAS is a one-dimensional hydrodynamic model which is mainly used for simulating steady and unsteady flow routing based on an implicit four-point finite difference scheme.The HEC-RAS model simulates only equilibrium ice jams, i.e., in steady state, not dynamically, as is the case in RIVICE.Hence, the HEC-RAS result can be seen as the maximum backwater level that can be determined if the ice jam formation is not limited by the supply of incoming rubble ice.Ice-jammed conditions: Simulating ice jams is more complicated than simulating open water and ice-covered conditions since more factors need to be used in RIVICE calculations to simulate ice jam formation and progression.Water level data at the Great Slave Lake (156.557m) and the Resdelta Channel (156.884m) on 27 April 1980 were used to set up and calibrate the model, respectively.After defining the incoming ice volume on ice jam, the water level profile can be determined (Figure 13).Gauged water levels and HEC-RAS modeling results are applied to calibrate the ice-jam modeling results.HEC-RAS is a one-dimensional hydrodynamic model which is mainly used for simulating steady and unsteady flow routing based on an implicit four-point finite difference scheme.The HEC-RAS model simulates only equilibrium ice jams, i.e., in steady state, not dynamically, as is the case in RIVICE.Hence, the HEC-RAS result can be seen as the maximum backwater level that can be determined if the ice jam formation is not limited by the supply of incoming rubble ice.Some of the gauges in the delta were installed during the winter, prior to ice-cover breakup.Ice jams regularly occur along this reach during the final breakup of the Slave River.The error between simulated values and observed ones is shown in Figure 14, which is acceptable.Some of the gauges in the delta were installed during the winter, prior to ice-cover breakup.Ice jams regularly occur along this reach during the final breakup of the Slave River.The error between simulated values and observed ones is shown in Figure 14, which is acceptable.

Ice Jam Flooding
Figure 15 shows the floodwater extensions, calculated with HEC-RAS, within the floodplain for the ponding cover water surface at elevations of 156, 157, and 158 m.For instance, when the backwater level exceeds 156 m, the area to the left-hand side of the red line will be flooded.To determine the backwater effects of ice jamming, which causes flooding, the RIVICE model was run three times under the following conditions: open water, ice-covered, and ice-jammed conditions.These three scenarios were simulated using the same boundary downstream water level elevation (156.45 m), but different upstream water flows.Associated with the HEC-RAS simulated flood extension under different backwater levels, the RIVICE modeling results can help determine the maximum possible flooding extent in the floodplain under different situations.

Ice Jam Flooding
Figure 15 shows the floodwater extensions, calculated with HEC-RAS, within the floodplain for the ponding cover water surface at elevations of 156, 157, and 158 m.For instance, when the backwater level exceeds 156 m, the area to the left-hand side of the red line will be flooded.To determine the backwater effects of ice jamming, which causes flooding, the RIVICE model was run three times under the following conditions: open water, ice-covered, and ice-jammed conditions.These three scenarios were simulated using the same boundary downstream water level elevation (156.45 m), but different upstream water flows.Associated with the HEC-RAS simulated flood extension under different backwater levels, the RIVICE modeling results can help determine the maximum possible flooding extent in the floodplain under different situations.

Local Sensitivity Analysis
Applying equation ( 2) to the local sensitivity analysis, e-values are positively proportional to parameter sensitivity.Given that, parameter sensitivity can be evaluated based on e-values.
Hydraulic properties, including upstream discharge (Q), downstream water level (W), and the riverbed roughness coefficient (n bed ), have the greatest impact on the backwater levels.The upstream discharge (Q) has the greatest negative impact on the ice jam length.
Ice transporting properties (v dep , v e , Table 1) have little impact on the ice jamming since the main process of ice jam formation is ice thickening as forces shove and retract the ice cover.Surprisingly, ice volume (V ice ) has some impact on the ice jamming length but not on backwater staging (Figure 16).

Local Sensitivity Analysis
Applying equation ( 2) to the local sensitivity analysis, e-values are positively proportional to parameter sensitivity.Given that, parameter sensitivity can be evaluated based on e-values.
Hydraulic properties, including upstream discharge (Q), downstream water level (W), and the riverbed roughness coefficient (n ), have the greatest impact on the backwater levels.The upstream discharge (Q) has the greatest negative impact on the ice jam length.
Ice transporting properties (v , v , Table 1) have little impact on the ice jamming since the main process of ice jam formation is ice thickening as forces shove and retract the ice cover.Surprisingly, ice volume (V ) has some impact on the ice jamming length but not on backwater staging (Figure 16).Physical properties of the ice (FT, PC, PS, ST) have some influence on the ice jam morphology.Increasing porosity will increase compaction between ice blocks when shoved, leading to a reduction in ice jam length.An increase in ice front thickness means more surface area is available for water shoving on the ice front to retract the ice jam cover, which corresponds to the shortening of the ice jam.
The location of the lodgment X has the largest positive impact on the ice jam length, potentially due to varying cross-section flow dimensions along the river.Moving the location further downstream leads to a decrease in backwater staging and shortening of the ice jam cover.
Mechanical properties, K1 and K2, have a significant impact on the ice jam length and partly on the water level, which is consistent with the research in [35].Ice cover internal resistance will increase with an increase in K1, which provides more resistance to shoving.Therefore, more ice can stay in place and the ice jam length is extended.Greater internal resistance corresponds to less thickening of the ice cover, allowing more flow underneath the ice cover and decreasing backwater effects.

Discussion
Both numerical modeling and remote sensing images are useful tools for studying river ice processes.These two methods were combined successfully to characterize ice jamming within the Slave River Delta.The excellent agreement between modeling results and measurement readings endorses the reliability of the river ice model RIVICE.Under open water conditions, this agreement verified the hydraulic properties set for the Slave River Delta.Under ice-covered conditions, with the limitation of sparse gauged data (07NC009), model calibration cannot be firmly supported by the field measurements.HEC-RAS estimated a greater ice volume than RIVICE, but lower than when no thawing effects are considered, which is reasonable.Physical properties of the ice (FT, PC, PS, ST) have some influence on the ice jam morphology.Increasing porosity will increase compaction between ice blocks when shoved, leading to a reduction in ice jam length.An increase in ice front thickness means more surface area is available for water shoving on the ice front to retract the ice jam cover, which corresponds to the shortening of the ice jam.
The location of the lodgment X has the largest positive impact on the ice jam length, potentially due to varying cross-section flow dimensions along the river.Moving the location further downstream leads to a decrease in backwater staging and shortening of the ice jam cover.
Mechanical properties, K1 and K2, have a significant impact on the ice jam length and partly on the water level, which is consistent with the research in [35].Ice cover internal resistance will increase with an increase in K1, which provides more resistance to shoving.Therefore, more ice can stay in place and the ice jam length is extended.Greater internal resistance corresponds to less thickening of the ice cover, allowing more flow underneath the ice cover and decreasing backwater effects.

Discussion
Both numerical modeling and remote sensing images are useful tools for studying river ice processes.These two methods were combined successfully to characterize ice jamming within the Slave River Delta.The excellent agreement between modeling results and measurement readings endorses the reliability of the river ice model RIVICE.Under open water conditions, this agreement verified the hydraulic properties set for the Slave River Delta.Under ice-covered conditions, with the limitation of sparse gauged data (07NC009), model calibration cannot be firmly supported by the

Figure 1 .
Figure 1.The Slave River Delta with important channels.(The dark lines in the bottom panel represent cross-sections spaced 500 m apart.Widths of cross-sections vary between 500 and 1200 m).

Figure 1 .
Figure 1.The Slave River Delta with important channels.(The dark lines in the bottom panel represent cross-sections spaced 500 m apart.Widths of cross-sections vary between 500 and 1200 m).

Figure 2 .
Figure 2. Plots of cumulative freezing degree days for winters (1 October-19 April) between the years 1980 and 2015 at Fort Smith, Northwest Territories.

Figure 2 .
Figure 2. Plots of cumulative freezing degree days for winters (1 October-19 April) between the years 1980 and 2015 at Fort Smith, Northwest Territories.

Figure 3 .
Figure 3.The histogram of pixel values of MODIS imagery.

Figure 3 .
Figure 3.The histogram of pixel values of MODIS imagery.

Figure 4 .
Figure 4. Steps taken in determining ice cover percentages.Figure 4. Steps taken in determining ice cover percentages.

Figure 4 .
Figure 4. Steps taken in determining ice cover percentages.Figure 4. Steps taken in determining ice cover percentages.

Figure 5 .
Figure 5. Water level data for open water case of the three gauge stations.

Figure 6 .
Figure 6.Measured flow data at Resdelta Channel and Fitzgerald Station.

Figure 5 .
Figure 5. Water level data for open water case of the three gauge stations.

Figure 5 .
Figure 5. Water level data for open water case of the three gauge stations.

Figure 6 .
Figure 6.Measured flow data at Resdelta Channel and Fitzgerald Station.

Figure 6 .
Figure 6.Measured flow data at Resdelta Channel and Fitzgerald Station.

Figure 7 .
Figure 7. Relationship between the cumulative freezing degree days and ice thicknesses.

Figure 7 .
Figure 7. Relationship between the cumulative freezing degree days and ice thicknesses.

Figure 8 .
Figure 8. Examples of aerial photographs used to validate the MODIS ice cover estimation (the number shows the location of the aerial photos.Data were acquired on 29 April 2015).

Figure 8 .
Figure 8. Examples of aerial photographs used to validate the MODIS ice cover estimation (the number shows the location of the aerial photos.Data were acquired on 29 April 2015).

Figure 8 .
Figure 8. Examples of aerial photographs used to validate the MODIS ice cover estimation (the number shows the location of the aerial photos.Data were acquired on 29 April 2015).

Figure 10 .
Figure 10.Ice volume calculation results from HEC-RAS and RIVICE.

Figure 11 .
Figure 11.Calibration for RIVICE simulation under open water conditions: (a) water level profile under minimum daily flow; (b) water level profile under maximum daily flow.

Figure 10 .
Figure 10.Ice volume calculation results from HEC-RAS and RIVICE.

3. 2 .
Model Calibration Water levels from 1980 recorded during the open-water and ice-covered period were used to calibrate the RIVICE model.Open water conditions: In 1980, during the open-water period, the daily flow reached a maximum of 4870 m 3 /s, on 13 June 1980 (Figure 6), when the water level elevation of Great Slave Lake was 156.43 m [24].These values were used as the boundary conditions in RIVICE.The minimum daily flow during the open water period in 1980 was 2660 m 3 /s, which occurred on 21 May 1980.Comparisons of gauged and simulated data are provided in Figure 11, showing good agreement between the two.

Figure 10 .
Figure 10.Ice volume calculation results from HEC-RAS and RIVICE.

Figure 11 .
Figure 11.Calibration for RIVICE simulation under open water conditions: (a) water level profile under minimum daily flow; (b) water level profile under maximum daily flow.

Figure 11 .
Figure 11.Calibration for RIVICE simulation under open water conditions: (a) water level profile under minimum daily flow; (b) water level profile under maximum daily flow.

Figure 12 .
Figure 12.Calibrated ice-covered water profile for the Slave River Delta.

Figure 13 .
Figure 13.Simulated ice-jammed water profile for the Slave River Delta.

Figure 12 .
Figure 12.Calibrated ice-covered water profile for the Slave River Delta.

Figure 12 .
Figure 12.Calibrated ice-covered water profile for the Slave River Delta.

Figure 13 .
Figure 13.Simulated ice-jammed water profile for the Slave River Delta.

Figure 13 .
Figure 13.Simulated ice-jammed water profile for the Slave River Delta.

Figure 14 .
Figure 14.Error between simulated values and observed data.

Figure 14 .
Figure 14.Error between simulated values and observed data.

Figure 14 .
Figure 14.Error between simulated values and observed data.

Figure 15
Figure15shows the floodwater extensions, calculated with HEC-RAS, within the floodplain for the ponding cover water surface at elevations of 156, 157, and 158 m.For instance, when the backwater level exceeds 156 m, the area to the left-hand side of the red line will be flooded.To determine the backwater effects of ice jamming, which causes flooding, the RIVICE model was run three times under the following conditions: open water, ice-covered, and ice-jammed conditions.These three scenarios were simulated using the same boundary downstream water level elevation (156.45 m), but different upstream water flows.Associated with the HEC-RAS simulated flood extension under different backwater levels, the RIVICE modeling results can help determine the maximum possible flooding extent in the floodplain under different situations.

Table 1 .
Parameters used in RIVICE model.