River discharge estimation from radar altimetry: Assessment of satellite performance, river scales and methods

Journal Pre-proofs Research papers River discharge estimation from radar altimetry: assessment of satellite performance, river scales and methods E. Zakharova, K. Nielsen, G. Kamenev, A. Kouraev PII: DOI: Reference:

S0022-1694(20)30021-4 https://doi.org/10.1016/j.jhydrol.2020.124561 HYDROL 124561

To appear in:

Journal of Hydrology

Received Date: Revised Date: Accepted Date:

3 June 2019 16 December 2019 8 January 2020

Please cite this article as: Zakharova, E., Nielsen, K., Kamenev, G., Kouraev, A., River discharge estimation from radar altimetry: assessment of satellite performance, river scales and methods, Journal of Hydrology (2020), doi: https://doi.org/10.1016/j.jhydrol.2020.124561

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2020 Published by Elsevier B.V.

1

River discharge estimation from radar altimetry: assessment of satellite performance, river scales and methods

Zakharova1 E., Nielsen2 K., Kamenev3 G., Kouraev4,5 A. 1-

Institute of Water Problem, RAS, Moscow, Russia, [email protected] University of Denmark, Elektrovej, Building 328, 2800, Kgs. Lyngby, [email protected] 3-Dorodnicyn Computing Centre, Federal Research Centre “Computer Science and Control”, Russian Academy of Sciences, Vavilova str, 40, Moscow, 119333, Russia [email protected] 4-LEGOS, Université de Toulouse, CNES, CNRS, IRD, UPS Toulouse, France, [email protected] 5- Tomsk State University, Tomsk, Russia 2-Technical

Abstract The ENVISAT, Jason -2 and -3, Sentinel-3A, CryoSat-2, and AltiKa satellite altimeters were used to estimate the discharge of two Arctic rivers: the Ob and Pur in western Siberia. The accuracy of the water height estimates from ENVISAT varied with river size from 0.63 m for the Ob (~2 km in width) to 1.1 m for the Pur (~0.5 km in width). A new method for water level estimation combining the CryoSat-2, AltiKa, and Sentinel-3A altimeters provided better height accuracy (0.49 m for the Ob). The rating curve method of discharge estimation outperforms the Manning formulation, on the Ob, with a root mean square error (RMSE) of 13% vs 20% for daily discharge and 1% vs 5% for annual discharge. Sensitivity analysis performed on the Manning formulation showed high sensitivity to the river depth parameterisation. A functional relationship was identified between the depth of rivers, in the north of western Siberia, and their widths and deposits types. The combination, for the Ob, of the rating curve method with the multi-satellite water level retrieval algorithm provided higher discharge accuracy than previous studies relying on satellite measurements, on numerical modelling or on their combination. A synthesis of worldwide river altimetry studies supported our finding that the accuracy of altimetric discharge estimations decreases for narrower rivers. The rating curve method applied to Jason 2/3 height measurements allows estimation of daily discharges with RMSE of 18% for the wide Ob River but 38% for the narrow Pur River.

1. Introduction

2 Rivers are an important element of the Earth’s climate system. They ensure connection between atmosphere and ocean. In many regions of the world, rivers are the only water resource supporting local socio-economic development. Quantification and monitoring of river flow is important not only for sustainable management of this valuable resource, but also for prediction of future conditions. Observations of river discharge (Q) have significantly reduced over the last 30 years. Although many national services continue observations of key rivers at key stations, the availability of these observations for the scientific community and for climate research is limited. In 2010 the database of the Global River Data Centre contained information gathered from 3562 gauging stations, against 9544 known stations prior to 2010. By 2015 the number had reduced by 1000 stations. Among the rest, only 8% are located on big rivers (annual flow > 10 km3), which are generally the key water resource rivers and the most informative indicators of regional long-term hydroclimatic changes. Alternative methods of estimation and monitoring of river flow (modelling and satellite observations) have been rapidly developed to fill this gap. Hydrological modelling is a powerful tool giving insights at regional or basin scales; while satellites can provide worldwide observations. The main satellite instruments used for river discharge retrieval are optical and radar imaging sensors, providing measurements of river width [Smith and Pavelsky, 2008] and water velocity [Tarpanelli et al., 2015, Beltaos & Kaab, 2014], and altimetric radar, measuring water height [Zakharova et al., 2006, Papa et al., 2012, Tourian et al., 2013]. Through statistical relations or hydraulic equations, these observations may be converted to the river discharge. Decreasing numbers of in situ observations, needed for calibration and validation of hydrological models, has promoted the development of satellite-model coupling methods [Biancamaria et al., 2009, Getirana et al., 2013, Emery et al., 2018, Domeneghetti et al., 2014]. Altimetry provides reliable regular and weather-independent measurements of water height (H) at a river cross-section called a virtual station (VS). Over the last decade the number of on-going altimetric satellite missions has increased (CryoSat-2, AltiKa, Sentinel-3A and B, Jason 2 and 3) allowing for observation of different reaches along a river. Significant progress has been made during the last 25 years in altimetric instrumentation and in processing of altimetric measurements. This progress has resulted in important improvements in the accuracy of water height retrievals over inland waters and in diverse new applications of altimetry in continental hydrology. The size of rivers having

3 altimetry-built water level time series has decreased from 2–3 km [Birkett et al., 1998] to 300–400 m [Schwatke et al. 2015] and even to 80 m [Michailovsky et al., 2012]. The errors of level estimates has also reduced: for large rivers to 0.30 m [Schwatke et al. 2015] and for medium size rivers to 0.3–0.7 m using ENVISAT [Michailovsky et al., 2012] or to 0.20 m using Jason-2 and Jason-3 [Biancamaria et al., 2018]. The few studies that have applied altimetry to rivers covered by seasonal ice did not provide any estimation of the accuracy of water heights. Indeed, many studies excluded the winter period because of the negative impact of ice on H estimation. Nevertheless, several studies, dedicated to the rivers of the Arctic Ocean watershed, have demonstrated good potential for altimeters in high latitude hydrology [Kouraev et al., 2005, Birkinshaw et al., 2014, Sichangi et al, 2016, Zakharova et al., 2019]. Low accuracy of the water level retrievals, as a result of winter ice, is the main drawback of altimetrybased methods of discharge estimation [Tourian et al., 2013, Sichangi et al., 2016]. Low spatiotemporal altimeter resolution also constrains the application of this method in the Arctic. The polarorbiting altimeters (ERS, ENVISAT, SARAL/Altika, Cryosat-2 and Sentinel-3) have a poor temporal resolution (27 – 369 days) and the TOPEX/Jason series altimeters, with a 10 day repeat cycle, provide measurements only as far north as 66.7°N; while the mouths of many Arctic rivers are located at higher latitudes. A better understanding of the performance of the different altimetric missions, and different water height retrieval approaches, over rivers with seasonal ice could help to extend the application of altimeters in Arctic continental hydrology. The application of the satellite altimetry for estimation of the river discharge started with the large world rivers: the Ob [Kouraev et al., 2004], the Amazon [Zakharova et al., 2006], and the Ganges - Brahmaputra system [Papa et al., 2010]. As the accuracy of water level retrieval has improved, more rivers of smaller size and of challenging fluvial geomorphology have been studied. As a rule, smaller rivers have greater variability in their water regime, shorter duration of flood events and more complex and irregular fluvial geomorphology. These factors are important not only for the accuracy of the altimetry-retrieved water heights, but also for the accuracy of the derived and applied discharge conversions. The main methods of discharge estimation based on altimetry rely on establishment of statistical relations, rating curves (RC), between in situ observations, at the nearest gauging station, and remotely sensed water heights [Kouraev et al., 2004, Birkinshaw et al., 2010, Tarpanelli et al., 2013]; or on hydraulic equations and their parameters (width, slope, height) estimated from satellite images [LeFavour et al., 2005, Gleason et al., 2014]. A particular interest of use of hydraulic equations

4 is their application to ungauged rivers. Several hydraulics-based approaches have been developed recently [Durand et al., 2016] and the performance of these approaches has been evaluated on a number of rivers throughout the World. Such research seeks to understand the limits of their applicability and to investigate potential improvements [Bonnema et al., 2016, Michailovski et al., 2013, Sichangi et al., 2016]. This study evaluated the capacity of the various altimetric satellite missions, and associated methods, to estimate the discharge of large and medium size rivers in the complex Arctic environment. In section 2 we describe the study rivers (the Ob and Pur). Section 3 introduces the satellites, auxiliary data and methodology. Section 4 presents the results of the water level, water slope and discharge estimation using ENVISAT, Jason 2/3 and combination of CryoSat-2/AltiKa/Sentinel-3A (CS2/S3). Water discharge is estimated by two methods. For the Manning method we began with an “ungauged” experimental setup, where equation parameters were taken from auxiliary sources (topographical maps and manuals). Next, the roughness coefficients were optimised allowing comparison with the fully calibrated rating curve approach. Finally, in section 5 we explore the effect of river size on the accuracy of discharge estimates. We also discuss the performance of the various altimetric satellites in the context of Arctic rivers. After demonstration of the sensitivity of estimation accuracy to the parametrisation used and to the altimetric water slope retrievals, we reflect on the shortcomings of, and potential improvements to, the hydraulics-based method.

2. Region of study Two rivers, of different size, located in similar natural conditions on the Western Siberia plain were selected for this study: the wide Ob River and the narrow Pur River. The annual flow of the Pur is 14 times less than that of the Ob. Both rivers belong to the Arctic Ocean watershed. They drain boggy areas and have large floodplains. Natural regulation significantly extends the duration of the spring flood, which can last up to 2 months. Much of the year the watersheds are covered by seasonal snow and the river channels are icebound. The Ob River The Ob is the third largest river, by discharge, flowing to the Arctic Ocean. The catchment’s area is 2,926,000 km2 and the annual flow is 406 km3. The main source of water is snow. Snow cover installs in September (in the north) or October (in the south). The melt starts in the middle of May (in

5 the south) and progresses northward by 2° each 8 days [Zakharova et al., 2014]. The Western Siberian plain is very boggy. Wet areas occupy up to 70–80% of the area in the of the central parts. The bogs intercept melting water and rain, resulting in significant attenuation of the flood peak and an increase in the flood duration. At low reaches the flood begins at the end of May and lasts until the end of July. The seasonal magnitude of daily discharges (Q) is 35000–45000 m3/s, with a winter minimum of 2500–4000 m3/s and a spring flood maximum of 38000–47000 m3/s. From October to May, the river is covered by ice. In August–September several small rain floods of magnitude 2000–4000 m3/s can occur. In its lower reaches the Ob has numerous branches. Many of them are of the same width as the principal channel and can be separated by 10–15 km. The width of the flood plain can reach 50 km. A terminal station providing discharge measurements is located at Salekhard city: 66.6°N, 66.5°E. Here the water flow from numerous branches is concentrated in the principal channel and the single small secondary channel withdraws a negligible amount of water (Fig. 1a). This reach, of 175 km length, was selected for water height and water slope retrieval from altimetric satellites. A river cross-section located 47 km from Salekhard (line A on Fig. 3a), where information on river depth was available for both the main and the smaller channel, was chosen for discharge evaluation using the Manning approach. The width of the main channel at this cross-section is about 2000 m. The width of the secondary channel, located 22 km to the east, is 400 m. The floodplain has a width of 25 km here. Between the main and secondary channels, there are several third-order branches with widths of 50–100 m. As their connectivity is not clear, we considered that their contribution to the total discharge is low and can be neglected. The channels are surrounded by lowdensity riparian forest, while the majority of the floodplain is occupied by bogs or lakes of fluvial or periglacial origin. The Pur River Originating on the northern slopes of the Sibirskiye Uvaly hills located in western Siberia, the Pur flows northward and joins the Ob estuary near its centre. The territory is part of a discontinued permafrost zone. The watershed is covered by frozen bogs. Wetlands (bogs and thermokarst lakes) occupy about 55 % of the area [Zakharova et al., 2014]. The watershed area above the terminal station at Samburg is 95,000 km2 and the annual flow there is 28 km3 [Zakharova et al., 2011]. As limited data from Samburg station was available, we used an additional gauging station, Urengoy, located 140 km

6 upstream. The drainage area at Urengoy is 80,400 km2 and water flow is 24 km3. Similar to the Ob, daily discharges have high seasonal variation. The spring flood peak can be 20 times higher than the winter discharge. Peak discharge varies between dry and wet years ranging from 6390 to 9430 m3/s. In autumn rain events can produce floods of lower magnitude. The winter low-flow period starts at the beginning of October and lasts about 200 days. At the end of winter, river ice thickness can reach 120

cm. a

b

Figure 1. Studied reaches of the Ob (a) and Pur (b) Rivers and ENVISAT/AltiKa, Jason 2/3 and Sentinel-3A tracks.

3. Data and Methods 3.1. Data For estimation of the discharge, we used the data coming from several sources: 1) high-resolution Landsat 8 images, which were used for processing of the water masks and allowed an assessment of the channel width variation; 2) water level and water slope time series retrieved from ENVISAT, Jason 2 and Jason 3 satellites, as well as from combination of CryoSat-2, AltiKa and Sentinel-3A; 3) topographic map and 4) in situ water height and discharge at Salekhard and Urengoy stations.

7 3.1.1. Landsat 8 images High-resolution Landsat 8 multispectral optical images, drawn from the US Geological Service data portal, were used for geographical selection of altimetric measurement and for evaluation of the Ob River width changes during the flood. Precise selection of altimetric measurement locations is important for retrieval of high quality water level time series [Kouraev et al., 2004, Zakharova et al., 2005, Dubey et al., 2015]. A water mask was defined from Level 1 Landsat 8 images using the modified differential water index (MNDWI). In previous studies the algorithm has accurately delineated water bodies of different sizes and types at similar latitudes in Central Siberia and been used in the production of water level time series for both a big Arctic river [Zakharova et al., 2019] and small (<1 km2) Arctic ponds [Zakharova et al., 2017].

MNDWI= (Green- SWIR1)/( Green+ SWIR1),

eq. (1)

where, Green - green band (0.53-0.58 µm), SWIR1 - short wave infrared band (1.57 - 1.65µm). The pixels with MNDWI >0 were classified as water. The water masks derived from August Landsat 8 scenes corresponded to the low-level period, when islands and sand banks are generally easily identified. These masks were then used solely for selection of the altimetric measurements used in building water level and slope time series. Eight further unprocessed Landsat 8 images, covering different phases of the water regime, were also used for manual measurement of the channel width of the Ob at different seasonal stages and for determination of the coefficients in equations relating width and altimetric height (Halti).

3.1.2. Altimetry Five altimetric satellite missions were used for Ob River water level estimation: ENVISAT, Jason-2 and its successor Jason-3, AltiKa, CryoSat-2 and Sentinel-3A. All satellites, except for Jason2 and 3, are polar orbiting and cover the lower reaches of all Arctic rivers. The water level for the Pur River was estimated only from ENVISAT and Jason 2/3 satellites. The ENVISAT satellite, with onboard RA-2 altimeter, operated from 2002 to 2012. RA-2 measurements up to 2010 (before the satellite was moved to the new orbit) were used for water level

8 processing. The repeat period (cycle) for ENVISAT is 35 days. An 18 Hz sampling rate gives an alongtrack resolution of 370 m. The Jason-2 altimetric satellite operated with a 10-day repeat orbit during 2008–2016. Its Kuband Poseidon instrument provided measurements with a 20 Hz sampling rate. The orbit inclination was 66.033 °N, meaning that the northernmost mission track crossed the Ob River at 66.13 °N, 35 km south of the Salekhard gauging station. The Urengoy gauging station on the Pur River is located between two Jason-2 tracks at 65.97 °N. In February 2016 Jason-3 was launched onto the same orbit. In order to ensure continuity of measurements the missions flew in tandem, with an 80 s difference, for 8 months. The Jason-3 calibration period was not considered during the discharge analysis. Nevertheless, we used it to assess the continuity of height measurements generated by the instrument in the succeeding period. The SARAL/AltiKa altimeter was launched onto the prior ENVISAT orbit in 2013. The SARAL instrument operates at a frequency of 37 GHz (Ka band). A high sampling rate of 40 Hz provides measurements at ~180 m along-track spatial resolution. The enhanced bandwidth allows for better vertical resolution (30 cm instead of 47 cm for ENVISAT) and finer waveform sampling. Since June 2016 AltiKa has been operating in a drifting orbit. CryoSat-2, launched by ESA in 2010, is the first mission to operate in synthetic aperture radar (SAR) and SAR interferometric (SARIn) modes. Compared to the conventional low-resolution mode (LRM), the Delay/Doppler technique makes possible a finer along-track resolution of approximately 300 m [Raney, 1998]. Here we applied the ESA level 1b product from baseline C. This product contains waveforms with 256 and 1024 bins for SAR and SARIn, respectively, as well as the range and geophysical corrections needed to construct a surface elevation model. Sentinel-3A, launched in 2016, operates globally in SAR mode. The dual-frequency SAR radar altimeter (SRAL) is able to operate in open or closed loops. In the open loop the range window is positioned by an on-board digital elevation model (DEM), whereas in closed loop the range window is automatically positioned based on the previous measurements. After March 2019, the open loop has been applied between 60°S and 60°N. To improve the on board DEM, the user community is encouraged

to

provide

surface

elevations

at

virtual

stations

through

the

web

page

https://www.altimetry-hydro.eu/. Here we applied the Level 2 product “Enhanced measurements”. This product contains waveforms with 128 bins, ranges, and all necessary geophysical corrections.

9

3.1.3. Topographic map Topographic mapping of the Russian part of the Arctic regions at 1:200000 scale was sufficiently precise for our analysis. The maps of this series are the basis of the ESA DUE Permafrost DEM product [Strozzi et al., 2011]. They also map the relief of the flood plains with a vertical resolution of 1 m and contain information about the type of river deposits and about mean river width and depth during the low-level period. This mapping was used: 1) in our "ungauged-like" experimental setup providing information on initial river depth (D0) and 2) for evaluation of the critical water level, at which water begins to inundate the floodplain.

3.1.4. In situ data The gauging stations Salekhard and Urengoy operated by Russian Hydrological Service provided both water level and discharge data. Water level is measured daily, and discharge several times per year during the main phases of the hydrological regime. In hydrological analysis these instrumental discharge measurements are used for adjustment of the rating curves for each given year. Daily discharges are calculated by the Russian Hydrological Service from these re-adjusted curves. For the Ob River, daily discharge data were available only for two years 2002 and 2003, while the water level measurements covered the whole 2002–2016 period. For the Pur River, discharge at Urengoy station was available for 2008–2015 and water level measurements covered the ice-free seasons from 2002–2016.

3.2. Methods Two approaches to water level estimation were applied in this study. The first used the classical virtual station concept, where all measurements at a river-satellite crossing, and within a given satellite cycle, are averaged to produce the water height data for a given day. The Landsat 8 water mask was used to identify points within the river channel. As the mask was developed in lowlevel conditions, many sand banks and islands were exposed. These were excluded from the river mask ensuring less contamination of the altimetric signal by land. In the second approach the river level time series were derived from multi-mission altimetry observations gathered over a segment of the river.

10

3.2.1. Water level retrieval from ENVISAT and Jason 2/3 satellites Water levels were estimated from the range parameter, which are part of the Geophysical Data Records (GDR) of the ENVISAT RA-2 and Jason 2/3 satellites. GDRs include the range measurements retrieved from the returned altimetric signal using a number of algorithms (retrackers). Good performance in level estimation for inland waters has been reported for the Ice1 retracker, provided in both ENVISAT and Jason GDRs [Calmant et al., 2013]. Environmental and geophysical correction of the altimeter range measurements, relevant to our study region, included ionospheric (DORIS), dry (ECMWF model) and wet (CLS) tropospheric corrections, solid Earth tide correction and correction for the satellite's centre of gravity and geoid. These corrections were applied to the range estimates to calculate water levels. The Ob River water level time series were retrieved for two ENVISAT virtual stations (named for their corresponding satellite track) located 3 km south and 16 north of the Salekhard gauging station. Time series were also derived for two, more northerly, Jason 2/3 virtual stations situated 62 and 64 km south of Salekhard. For the Pur River water levels were retrieved for three ENVISAT track locations (326, 365 and 784) near to the Urengoy gauging station and for seven Jason tracks (238, 59,162, 135, 86, 211 and 10). Based on the satellite revisiting times, the ENVISAT water level measurements had a 35 days period, while the Jason data had a 10 day period.

3.2.2. Water level and water slope retrievals from multi-mission approach Water levels from CryoSat-2, AltiKa, and Sentinel-3A were derived in a similar manner, as from Jason and ENVISAT, except that the range parameter is based on a primary peak 80% threshold retracker. The primary peak retracker selects the first peak in the waveform, above a specified power threshold, as the primary peak. Retracking is then performed on the corresponding sub-waveform; a retracking point on the leading edge is identified as the decimal bin corresponding to 80 % of the OGOG amplitude. A detailed description of this retracker is given in Jain et al, (2015). CryoSat-2 operates in SARIn mode over the Ob River allowing for range correction from off-nadir reflections. We applied the off-nadir range correction described in Armitage and Davidson (2014). Water levels from CryoSat-2, AltiKa, and Sentinel-3 were derived from an approximately 175 km river segment within the

11 region 65.6–68.5 °E and 66.1–66.6°N. The satellite measurements over the river were extracted using the Landsat 8 water mask. To create the water level time series for the Ob River, we applied a state-space model, processed using an auto regressive model of order one (AR(1)) in the time direction. Water level observations were taken to be cubic spline functions of distance (eq.2).

(eq 2), where

η(ti) = ρ (ti-1) + ξ ,

-1 < ρ < 1,

ξi is an error term of normal distribution ~ N(0,ση2); β(sati) are

bias parameters related to the different satellite missions; εi is the observation noise. The splines based on points at equidistant locations respond to the underlying topography α(xi) and to the variations in the water level with distance τ(xi). Observational noise is assumed to follow a mixture of the Gaussian and heavy-tailed (Cauchy) distributions as described in Nielsen et al. (2015). We used 300 time steps in the AR(1), which is equivalent to a temporal resolution of approximately 7–8 days, with 5 nodes in the spline functions. The model allows for construction of a water level time series at any location along the considered river segment. This makes it possible to calculate the river slope at any location and time. Figure 2 displays the modelled river levels along the river for four dates as an example. For a given time, the slope is derived from the corresponding spline line as the water level difference between selected two locations divided by the distance between these locations. Here we selected 20 km segment near Salekhard station to derive the slope time series. The standard deviation of the river slope is derived using the delta method.

12

Figure 2: CS2/S3 modelled water levels (lines) as a function of distance displayed for 4 distinct times (2016.5, 2016.8, 2017.1 and 2017.3). All altimetric water level measurements extracted using water mask (symbols) are shown at the same time +/- 8 days. The distance is calculated from lower to upper reaches.

3.2.3. Discharge and water flow estimation In this study we used two different methods of discharge estimation from altimetric measurements: a rating curve method and the Manning method. The rating curve method was applied for both the Ob and Pur Rivers. As multi-mission water level retrieval for the Pur did not provide satisfactory results and, consequently, the water slope could not be estimated, the Manning-based estimation of Q was performed only for the Ob. Manning-based Q estimation relies on physical lows, which allow its application to ungauged river reaches. Satellite derived main variables are the inputs. To investigate the potential of the Manning method for ungauged Arctic rivers, initially, we used parameterisation found in manuals and the topographical map. Then, to compare its performance with that of the rating curve method, the parameters were calibrated. Rating curve model

13 Rating curves discharge estimation is a fundamental approach used to obtain daily discharges at gauging stations. Relations are established between simultaneously measured water heights and water discharges. The latter are calculated from instrumentally measured water area and velocity in a given river section. These measurements are then used for development of rating curves describing the H-Q relations. Daily discharges can then be calculated from the daily measurements of water level. The rating curves have a power form and can be approximated by equation (3).

Q=(H-a)b

eq. (3),

where H is water height, a and b are parameters of the equation. During periods of low variability in the water level, a polynomial function (equation (3)) can produce a better fit between Q and H [Zakharova et al., 2019].

Q = cnHn + cn−1Hn−1 + . . . +c0

eq. (4),

where H is water level, n and cn are parameters of the equation. For the winter period, we tested both equations, (3) and (4), and selected that which produced better accuracy. As water flow is not a stationary process, H-Q relations are rarely uniform, especially in the case of boreal or Arctic rivers with ice cover, or when rivers have large flood plains covered by high vegetation. Both ice and high vegetation significantly change the roughness conditions resulting in changes of water velocity and in point deviations from the main H-Q stationary line. In [Kouraev et al 2004 and Zakharova et al., 2019] we demonstrated, using in situ measurements for Salekhard station on the Ob River and for Kusur station on the Lena River, the existence and advantage of multiple H-Q relations. In this study we applied a similar approach and established a set of three rating curves (specifically for winter, flood rise and flood recession hydrological phases) for each virtual station or retrieved altimetric water level time series. Discharge estimation consisted of three steps. 1. Selection of the hydrological phase based on the altimetric water level and radar backscatter dynamics. 2. Establishment of rating curves between altimetric heights and in situ discharge using a calibration data set (see section 3.2.4) and fitting of equation (3) (or equation (4), if necessary).

14 3. Calculation of the altimetry derived discharge, for the validation period, from the developed Halti-Q rating curves. Manning method The river discharge is calculated as a product of cross-section area and mean water velocity:

Q = A × V,

eq. (5)

A - flow contributing section area, m2, V- water velocity m/s. For rivers, where the channel width is significantly larger than the depth, the area can be approximated by assuming a rectangular crosssection: A= B × h, where

h = D0 + ΔH

eq. (6), eq. (6a),

B - channel width, m; h - mean depth, m; D0 - initial depth at minimal water level, m; ΔH - water height correction, m. The Manning water velocity has following formulation:

where

V = 1/n ×R2/3 × S1/2

eq. (7),

R= A/P,

eq. (8),

n - Manning’s roughness coefficient, R - hydraulic radius, S - water surface slope m/m, P -wetted perimeter, m. As the channel width (B) varies with the water level, a relation B-H was established using the river width, manually measured using Landsat images, and the river height, retrieved from altimetry at the moment of the image acquisition.

Bm = am × Halti bm

and

Bs = as × Halti bs

eq. (9),

where Bm is the main channel width and Bs is the smaller channel width. The fitting of the points provided the following coefficients for equation (9): am =2042, bm =0.24, as =554 and bs =0.29. One of the important parameters in the Manning equation is the river mean water depth (h), which serves for calculation of the cross section area (A). The precise water depth can only be determined by in situ

15 profiling. Water depth information is provided as part of gauging station metadata and it is also recorded on ancient topographic and modern navigation maps. For our test aimed simulating an ungauged reach, we used the information provided by the topographic map. This depth (D0) refers to the low flow stage and needs to be corrected for water level changes (ΔH). As the shape of the channel section selected for this study is complex, we adopted a twoblock scheme for estimation of the effective area and hydraulic radius of the main branch (Fig. 3). Using Landsat images, altimetric heights and the topographic map (providing the topography of the floodplain) we estimated that water starts to flood the sandbank at Halti = 2.5 m, and at Halti = 4 m it exits onto the floodplain. The total width of the main branch Bm (including the main channel and the sandbank) was derived using equation (9). The area A at each altimetric water stage was estimated as the sum of the area of main channel and sandbank (if inundated), using equation (6), considering correction for depth ΔH in each sub-section. The hydraulic radius was then estimated using this area and the wetted perimeter P (black bold line on Fig. 3). A roughness coefficient reflects the resistance to flows in channels and over floodplains. Suggested values for Manning's n, tabulating the factors that affect roughness, can be found in Chow (1959). The detailed procedures for determining Manning's n values for natural channels and flood plains is described in Arcement (1989). The most important factors affecting channel n are the type and size of the bed/bank materials (sand, gravel, cobbles, rock, etc) and the shape of the channel

Main channel

Sandbank

Flood plain Halti=2.5 m abs

D0=12.2 m

Small channel

(depth, irregularities, obstructions, vegetation, meandering, etc).

Halti= 4 m abs D0s=6 m

Figure 3. Simplified schema of the Ob River cross-section adopted for estimation of the effective area and hydraulic radius for ice-free conditions. D0 – initial depth of the main channel taken from the topographic map and corresponding to the average low-flow season water height, D0s - the same for small channel, Halti =2.5 m abs – altimetric water level, at which the sandbank is submerged, Halti =4 m abs – altimetric water level, at which the floodplain inundation starts.

16

In winter, ice cover creates an additional boundary with associated hydraulic roughness. The net result is reduction of the hydraulic radius of the channel and changes in effective channel resistance. An existing formulation of composite (effective) roughness (nc) of an ice-covered river is based on weighting of the individual bed (nb) and ice (ni) roughness coefficients using the respective wetted perimeters of the channel bed and the ice cover. For this study, a more accurate formulation, taking into account the vertical distribution of water velocity under the ice, and known as the BelokonSabaneev equation, was adopted as proposed by Bruner (2016).

eq. (10)

Calibration of the parameters of the Manning equation The Manning equation is calibrated by selecting the combination of parameter values that leads to the lowest levels of a multicomponent series of errors when compared with observations. These component errors are aggregated into a number of objective functions. A solution to the combined objective functions is defined as efficient (Pareto optimal) if, by parameter variation, no single objective function can be improved without losses elsewhere in the solution set. Each possible combination of parameters values (potential calibration solution) leads to a particular set of objective functions values. This combination of parameters values and corresponding objective function outcomes can be Pareto efficient (Pareto optimal, Pareto dominant), or it can be within a Pareto efficient boundary and, thus, be sub-efficient. The set of all Pareto optimal objectives combinations is called the Pareto frontier (Pareto front) and corresponding calibration solutions form the efficient calibration solution set (Pareto optimal calibration solutions). To find the Pareto optimal calibration solutions, we used the multi-criteria identification sets method [Kamenev 2016, Kamenev, 2018]. This method is based on approximation and then visual examination of identification sets. These may be represented as a multidimensional graph of the vector objective function, its lower dimensional projections, or in other forms (projections, parts, sections, or slices), and may include sets of optimal, sub-optimal, efficient and sub-efficient parameter value combinations.

17 Both bed and ice roughness coefficients were simultaneously calibrated using the multicriteria identification sets method of Pareto domination with four objective functions. Two objective functions, Root Mean Square Error (RMSE) and bias, evaluated differences between simulations and observations during the ice-free period, and two similar functions were used for winter, when both ice and bed coefficients are relevant. The possible ranges of nb and ni were limited to values known to occur in natural conditions in large rivers. The nb range was 0.015–0.35, while for ni a larger band 0.01–0.04 was applied. Each combination of nb and ni values (calibration solution) leads to particular outcomes for the four objective functions.

3.2.4. Validation approach The altimetric datasets were divided into calibration and validation periods specifically for each mission (Table 1). Using the rating curves, established from simultaneous daily in situ H and Q for 2002–2003, the daily discharge at Salekhard station was estimated for the years 2004–2016. The estimation procedure was based on the set of the rating curves described in detail in Kouraev et al. (2004) and Zakharova et al., (2019). The reconstructed set of station data were then used for equation calibration and for validation of discharge estimates.

Table 1. Calibration and validation periods for different altimetric missions. Mission

Operational period

Calibration period

Validation period

ENVISAT

2002-2010

2003-2005

2006-2010

Jason2/3 (the Ob River)

2008-ongoing

2010-2012

2008-2009, 2013-2016

2009-2011

2012-2015

2011-2012

2013-2016

Jason2/3 (the Pur River) CS2/S3

2011-ongoing

The accuracy of the altimetric discharge and flow estimations was evaluated using bias, RMSE, and the Nash-Sutcliff coefficient (NS).

RMSE =  (Qalti-Qinsitu)2/N,

eq. (11),

NS= 1- ( Qalti-Qinsitu)2/  (Qinsitu -mean(Qinsitu))2,

eq.(12)

18 where Qalti- discharge estimated from satellite measurements, Qinsitu - in situ (or reconstructed) discharge at station, mean(Qinsitu) - average in situ (or reconstructed) discharge at station.

4. Results 4.1. Water level and water slope time series The water level data, from different satellites, for the Ob and Pur Rivers are presented in Fig. 4. The ENVISAT 51 track crosses the Ob River within 3 km of the gauging station. This proximity allows estimation of the error in Halti retrievals (0.63 m RMSE, see Table 2). During the low-level period, the width of the Ob, in the study reach, varies between 1.7 and 2.5 km. This width is large enough that the altimetric H average remains valid with low or moderate land contamination effects. ENVISAT water level time series for the Pur River were produced for three virtual stations. The virtual stations 326 and 365 are situated ~65 km north of the gauging station, while the virtual station 784 is within 2 km of the gauge. As the fluvial geometry at ENVISAT 784 virtual and at the gauging station are similar, we can assess the accuracy of the water height estimation, from ENVISAT, for a medium size Arctic river with complex fluvial morphology. The width of the Pur River varies from 0.5 to 1.2 km. Numerous sand banks and islands reduce the effective width significantly at the end of summer and in winter. Despite water mask application, signal pollution by the banks and islands results in their surface being detected as water surface. As a result, the altimetric time series had many anomalously high H measurements during flood recession and in the low-level period. We found that when ice-free, the RMSE between in situ observations and ENVISAT Halti is higher for the Pur than for the Ob River (0.79 m). As in situ measurements for winter are unavailable, we were unable to evaluate the accuracy of the height estimates over river ice. The continuity of measurements from Jason 2 and Jason 3 was verified during their 8 months of parallel operation on both wide and narrow rivers. An average bias between the instruments for the Pur River was 0.17 m and 0.01 m on tracks 59 and 162, respectively. The highest differences reach 2 m and occurred before or during the ice melt and during the summer minimum. At this time, the radar footprints are at their most heterogeneous. For the Ob River the average bias between missions was 0.10 and 0.26 m for tracks 112 and 187, respectively. The maximum difference was significantly less than for the narrow Pur River and did not exceed 1 m.

19 Among seven northernmost Jason tracks, only track 59 and track 162, with few anomalous Halti measurements, were selected for discharge estimation for the Pur River. At the crossover locations of these tracks with the river, the channel is 500–600 m wide and free of islands. A sand bank boarding right side of the channel is oriented such that only part of the return signals can be contaminated (suggesting ~200 m variation in the altimetric measurements around a mean nominal Jason orbit). As the tracks are located far from the gauging station, we could not determine the water level estimation errors, which can also be affected by changes in channel morphology. Water level time series retrieval from a combination of the CryoSat-2, AltiKa and Sentinel-3A instruments allowed altimetric measurements at the gauging station and, consequently, evaluation of the accuracy of estimated heights without any ambiguities caused by differences in hydraulic or morphometric conditions in the two compared reaches. The accuracy was high: 0.49 m RMSE and 0.18 m bias. The difference (bias) between in situ and satellite observations is low during the flood and in winter (0.04 m and 0.07 m, correspondingly). The bias increases up to 0.26 m during summer rain floods. The highest overestimation (0.79 m) occurred in April. On the Ob River at that time of year, water often exits over the ice through polynyas or pressure fractures. Wet ice or snow changes the return echoes (as could be seen in backscatter dynamics) and negatively affects the accuracy of the altimetric estimates. The CS2/S3 combination (similar to the ENVISAT) produced reliable estimates of water height during the winter over ice covered surfaces. The Jason Ice 1 algorithm behaved well at the beginning of the ice period, but overestimated the ice height during the second part of winter, especially before melting. The similar problem was noted for the Lena River in our previous study [Zakharova et al. 2019]. The low reaches of the Ob River are characterised by very low slopes, leading to development of a wide (up to 50 km) floodplain with multiple branches. The water slope, retrieved from the altimetric measurements, was 0.01–1.3 cm/km (Fig. 4c). A similar value, 2.7 cm/km, was found by Beltaos and Kaab (2014) for another big Arctic river, the Mackenzie, in September 2009. The water slope is strongly correlated (r=0.95) with the altimetric water height (not shown) and its response to short-term summer and main spring floods is clear. The error (2 σ) of the slope at a given time ranges between 0.14 cm/km and 1.0 cm/km. This error estimate might be underestimated since the correlation among along-track water level measurement is not accounted for. An average relative error is 22%. More than

20 80% of the S estimates have errors of less than 35%. The highest relative uncertainties were observed during the winters of 2012 and 2016, when the absolute slope values were close to zero.

Table 2. Accuracy of the ENVISAT and CS2/S3 water height retrievals at virtual stations located near in situ gauging stations. River, period

Satellite/ virtual station

RMSE, m

Ob River, annual

ENVISAT, 51

0.63

Ob River, annual

CS2/S3

0.49

Pur River, ice-free

ENVISAT, 784

0.79

21

Ob River

Water level,m

8

Q_insitu J2TR187 J2TR112 ENTR51 ENTR12 CS2/S3

4

0 a) 2004

Water level, m

16

2008

2012

2016

Pur River Q_insitu J2TR59 J2TR162 ENTR326 ENTR784 ENTR365

14

12

10

8 2004

2008

2012

2016

b)

2012

2016

c)

2 Ob River

Slope, cm/km

1.6 1.2 0.8 0.4 0 2004

2008 Time

Figure 4. Time series of the water level retrieved for the Ob (a) and the Pur (b) from different satellite missions, as well as the time series of the water slope estimated for the Ob from multi-satellite approach (c).

4.2. Discharge and water flow retrievals

22 The classical rating curve method was applied to the Jason altimeters for both rivers, and to ENVISAT and the CS2/AltiKa/S3 combination for the Ob River only. This method was not applied for the ENVISAT time series for the Pur River, as the number of Halti retrievals in the open water period was insufficient to establish reliable rating curves. The Manning method, involving both altimetric water height and water slope, was applied only to the Ob River.

4.2.1. Rating curve model The Q (discharge) values reconstructed from in situ water height data were used to evaluate the accuracy of the altimetry-simulated discharge for the Ob River. The accuracy of this reconstructed discharge (estimated as RMSE relative to the in situ Q for 2002–2003) depends on the hydrological phase. The lowest accuracy 13% (or 1917 m3/s) occurred in the flood rise period. For other periods the accuracy is very high, 6% RMSE (or 310 m3/s) for winter low flow and 7% RMSE (or 1790 m3/s) for summer recession. The accuracy of the Jason 2 and Jason 3 water level estimates deteriorated significantly during the winter critically affecting the Halti -Q relations during this period. Power fitting (eq. 3) did not produce reliable results and so we fitted the points by minimization of RMSE and bias using a first order polynomial function (eq. 4). This function also produced more accurate winter rating curves for other satellites. Table 3 shows the accuracy of the rating curve fitting. The number of ENVISAT and CS2/S3 water height estimates during the Ob flood rise was not sufficient to establish accurate rating curves. As the CS2/S3 water level estimates were relocated to the Salekhard station and the ENVISAT tracks are within 16 km of the gauge, we decided that the application of the station H-Q curve, instead of the Halti-Q curve, would provide better results for this period. Previous studies have demonstrated the advantage of densification of the accuracy of altimetric observations to discharge retrieval [Tourian et al., 2016, Zakharova et al., 2019]. Our final Qalti time series, presented in Table 4 for both rivers, were obtained by simple combination of the simulated discharges from virtual stations of the same satellite missions. This procedure did not take any additional mathematical treatment, as the virtual stations observations are not simultaneous and lateral water input is negligible.

23 Table 3. The accuracy of the rating curves built-up from altimetric height and in situ discharge and used for estimation of the Qalti Phase

winter RMSE,

flood R2

m3/s

RMSE

recession R2

m3/s

RMSE

R2

m3/s

J2_Tr_187

317(4)

0.05(4)

2932(3)

0.90(3)

1670(3)

0.95(3)

J2_Tr_112

209(4)

0.54(4)

4227(3)

0.92(3)

1670(3)

0.95(3)

J2_Tr_59

121(4)

0.43(4)

385(3)

0.93(3)

202(3)

0.89(3)

J2_Tr_162

204(4)

0.07(4)

461(3)

0.95(3)

357(3)

0.86(3)

EN_Tr_51

354(4)

0.91(4)

93*

0.99*

1852(3)

0.97(3)

EN_Tr_12

272(4)

0.93(4)

93*

0.99*

1632(3)

0.98(3)

CS2/S3

127(4)

0.88(4)

93*

0.99*

1684(3)

0.95(3)

*station rating curve is used (3) (4)

RC are fitted using power function (eq.3) - RC are fitted using first order polynomial function (eq.4)

On the Pur River, during the ice-free period, the number of the retrieved ENVISAT Halti values (only 1–3 per year) was insufficient for reliable rating curves, and the ENVISAT time series were excluded from discharge processing for this river. When enough simultaneous in situ discharge and altimetric water heights observations are available, and the morphology of the virtual and gauge reaches is similar, the rating curve approach allows for estimation of Q with very high accuracy (Fig. 5 and Table 4). The ENVISAT results had the highest Q errors, while the multi-satellite method had the lowest. For daily discharges, the error distribution among the missions is, in general, in agreement with the accuracy of the rating curves established for each virtual station (see Table 3). Integration over time reduces the errors by 3-6% for monthly and by 12–16% for annual flow. Also in agreement with the rating curves errors, the medium size Pur River had lower discharge retrieval accuracy than the Ob River. The difference between simulated and in situ discharge was higher during the dry years 2011–2013, when the annual flow halved. The largest differences were during flood recession and in winter. For the Ob River the frequency of water level measurement, by different altimeters, is not critical to monthly or annual flow estimation. The long flood duration permits adequate sampling, even by ENVISAT with its 35-day cycle. Combining discharges from several virtual stations, in the case of ENVISAT altimeter for the Ob River and Jason altimeters for the Pur River, increases the temporal

24 sampling frequency and moderately (1–2%) improves the accuracy of the monthly and annual flows. There was, however, no accuracy gain from combining Q estimates from the Jason 112 and 187 virtual stations, on the Ob.

Ob River

discharge, m3/s

50000 40000 30000 20000 10000

0 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016 2017 Time Qreconstructed

CS2/S3

J2_Tr187

J2_Tr112

EN_Tr12

EN_Tr51

a) Pur River

discharge, m3/s

8000 6000 4000 2000 0 2008

2009

2010

2011

2012 2013 Time

Q_Urengoy, insitu

J2_Tr59

2014

2015

2016

2017

J2_Tr162

b) Figure 5. Comparison of RC Q daily estimations with in situ/reconstructed observations for the Ob River (a) and for the Pur River (b).

4.2.2. Manning method To demonstrate the capacity of the Manning method to estimate water discharge for ungauged areas, we did not, initially, perform a calibration and instead used typical tabulated values of ice and bed roughness [Chow 1959, Bykov et al., 1965, Bruner, 2016]. For large rivers with no vegetation and sandy deposits, these sources propose an nb value of 0.025. As the freezing of Arctic rivers is often characterised by transport and accumulation of frazil ice under the ice edge, the lower ice boundary has a very rough surface [Beltaos et al., 2009]. We chose the highest value 0.03,

25 proposed by [Bruner, 2016] for ice of 30–90 cm thickness. In the studied reach, the Ob River valley is composed of two river branches and a large floodplain. The river discharge was estimated separately for each of these elements and then summed. Water discharge from the flood plain was estimated for conditions in which the water height exceeds a critical level. By analysing the topographical map, Landsat images and Halti values, a critical altimetric level of 4 m was chosen (see Fig. 3). The roughness coefficient of the flood plains (nf) is high compared to nb and typically varies between 0.1 and 0.2. We selected a value of 0.13, proposed by Bykov (1965), for conditions of dense brushes/lowdensity riparian forest. The simulated discharges, before calibration of parameters, were comparable with those obtained using the classical rating curve approach and gave the same high accuracy for annual water flow (Table 4). The discharge passes mainly through the primary river channel. The contribution of the small channel is about 2% (6% max) and of the floodplain is about 1.5 % (4% max). The Manning Qalti reproduced well the increased winter discharge of 2014–2016, although the method underestimated discharges in 2011–2013 and overestimated in 2015 and in 2016 (Fig. 6a). Flood discharge, when the altimetric observations caught the peak in 2013 and in 2016, had low errors of 5% and 12%, respectively. When the altimetric observations truncate the flood peak and discharge is interpolated, peak discharge estimation errors increase up to 35–45%.

a)

26

600

100

W insitu, km3

W insitu, km3

80 60 40 20

400

200

0

0 0

20

40

60

80

100

0

200

W alti, km3

400

600

W alti, km3

b)

c)

Figure 6. Comparison of Manning daily discharges (a), monthly (b) and annual (c) water flow (W) with in situ (reconstructed) observations at Salekhard station.

Table 4. Accuracy of the discharge estimations from the present and from other studies. For the ENVISAT and the Jason 2/3 missions, the discharge after densification of the virtual stations is used. The NS is not calculated for RC annual water flow as the validation period of 5 years is not long enough to well represent this statistic. Parameter (method)

RMSE, m3/s,

RMSE%

Bias, m3/s

Bias,%

NS

source

Ob River Q daily (altimetry/ RC) Q monthly (altimetry/ RC) Q annual (altimetry/ RC) Q daily (altimetry/ Manning) Q monthly (altimetry/ Manning) Q annual (altimetry/ Manning) Q daily (altimetry/ RC) Q monthly (altimetry/ RC) Q annual (altimetry/ RC) Q daily*** (altimetry, optic/Bjerklie eq.) Q monthly (altimetry/RC)

3290/2416/1957* 26/18/13*

581/187/32*

5/1/0.2* 0.91/0.95/0.97* This study

3091/1508/1278*

23/12/9*

-649/-385/84* -5/-1/1* 0.99/0.96/0.99* This study

1504/294/187*

11/3/1*

-882/-135/84* -7/-1/1*

3309/3033**

23/20**

-1051/-172**

2422/2343**

19/18**

883/642**

7/5**

na

This study

7/-1**

0.98/0.98**

This study

-662/531**

-5/4**

0.94/0.95**

This study

-662/531**

-5/4**

0.89/0.94**

This study

675

8

1440

11

400

3 0.98

880

16

Kouraev et al., 2004 Kouraev et al., 2004 Kouraev et al., 2004 Birkinshaw et al., 2018 Tourian et al., 2013

27 Q monthly (model- altimetry) Q daily (model) Q daily (model)

1917

14

0.99

Biancamaria et al., 2009 Decharme et al., 2012 van Vliet et al., 2012

2054 80

-10 Pur River

Q daily (altimetry/ RC) Q monthly (altimetry/ RC) Q annual (altimetry/ RC) Q daily (model)

341

38

-71

-4

0.93

This study

280

31

-71

-4

0.90

This study

157

17

-71

-4

0.77

This study

-10

0.78

Gusev 2015

et

al.,

* ENVISAT/Jason-2,3/CS2- AltiKa- S3; **using Manning equation before/after calibration *** for Kalpashevo station located in ~1200 km upperstream from Salekhard; na - not applicable.

In hydraulics models the roughness coefficient is usually calibrated against in situ discharge observations. Using the multicriteria identification sets method we found the ideal values of nb and ni to be 0.024 and 0.019 respectively. Interestingly, the optimisation suggested a lower ice roughness coefficient than the bed value. According to HEC-RAS model manual [Bruner, 2016], this calibrated value corresponds to rippled ice, while our initial thought was that Ob River ice ridges during formation and stays rough all winter; and we initially set ni to 0.03. As our initial choice for nb (0.025) was close to the calibrated value, the overall performance after optimisation did not change significantly. The annualized values of the accuracy metrics for daily discharges were 3033 m3/s (or 20 % from mean annual discharge), -172 m3/s and 0.98 for RMSE, bias and Nash-Sutcliffe efficiency, respectively. As the seasonal magnitude of the discharge is high (Qflood/Qwinter >10), it seems that high uncertainty in ni does not have a critical effect on annual accuracies. After calibration, the uncertainty in the peak discharge for years 2013 and 2016, with nearpeak Halti measurements, reduced to 1% and 8% respectively.

5. Discussion 5.1. Assessment of the performance of altimetric satellites for estimation of discharge of Arctic rivers. The accuracy of discharge estimation from altimetric satellite observations depends on the accuracy of the satellite-retrieved variables (water height and water slope), the accuracy of the

28 algorithm parameterisation and the degree of simplification (generalisation) of the hydrological processes within the applied algorithms. Water level retrievals for seasonally ice-covered rivers are rather rare [Kouraev et al., 2004, Frappart et al., 2010, Biancamaria et al., 2009, Tourian at al., 2013, Birkinshaw et al., 2016, Bjerklie et al., 2018, Zakharova et al., 2019, Sichangi et al., 2016]. The ice period is excluded from many studies because of high uncertainty in the water levels. On sea ice, it is thought that the signal of Kuband altimeters penetrates the snow and is reflected from the ice surface. In situ measurements of water level are made in ice pits where the water is close to the ice surface due to isostatic equilibrium and ice plasticity. A difference, <15 cm, between in situ and altimetric observations can be expected. A higher difference may be observed for the Ka-band AltiKa altimeter, as its signal is reflected primarily from the snow/air interface [2016, Remy et al., 1999, Guerreiro et al.]. In practice, the accuracy of height retrievals over river ice depends significantly on the instrument design (conventional/SAR altimeter) and on the retracker configuration. The present study, as well as data from the DAHITI (Database for Hydrological Time Series of Inland Waters) database1, demonstrated that the water level time series derived from ENVISAT, AltiKa, and Sentinel 3A, as well as from their combination, better describes the winter recession of Arctic rivers. In winter, the Jason satellites often detect the higher surfaces of the flood plain or even terraces [Zakharova et al., 2019]. This anomalous detection can last several cycles, be reproduced each year, be observed over several years, or occur only in the second part of the winter before flood. The reason of this unstable behaviour is unknown. Fortunately, for estimation of river discharge, this anomaly can be handled, in part, by establishing special winter rating curves (when the rating curve method is applied). The use of other methods based either on various modifications of the Manning equation or on coupling of numeric models with altimetric observations, can result in high winter discharge errors. In our case, the Jason 2 and 3 winter anomalies negatively affected the accuracy of the winter rating curves (see Table 3). Our river masks were defined during a low-level period when the channel width is minimal and sand islands and banks have maximal visibility. We expected the H estimates to be better during high waters, as there is less land at the mask border to affect the radar waveforms. However, analysis of the seasonal distribution of errors, from CS2/S3 and ENVISAT, found that three validated time series

1

https://dahiti.dgfi.tum.de/en/map/

29 had lower absolute errors in H during winter (respectively 0.44 m and 0.10–0.12 m) and higher errors during the flood (respectively 0.53 m and 0.56–0.87 m). The accuracy of the derived altimetric heights and the number of the valid observations used for rating curve definition, determined the accuracy of the fitted curves. The limited number of available flood measurements (due to low sampling frequency) did not permit a robust fitting of the flood Halti -Q curves for CS2/S3 and ENVISAT for the Ob River. We resolved this problem by applying the gauging station equations: the CS2/S3 measurements were relocated to the Salekhard reaches and the ENVISAT virtual stations are within 3 km and 16 km of the gauge. The few ENVISAT observations available for ice-free period constrained the application of this mission to the smaller Pur River. The overall errors in discharge estimation for the Ob River demonstrated an improvement in altimetric satellite performance from the earlier ENVISAT to the newer Cryosat-2/Sentinel-3 missions, and in their combination.

5.2. Effect of river size on the accuracy of altimetrically derived discharges. Two studies in the Amazon basin, which coupled numerical models with altimetric measurements [Getirana et al., 2013 and Emery et al., 2018], reported that the accuracy of their discharge simulations reduced with the river size. Our results also demonstrated lower discharge accuracy from the Jason satellites for the smaller Pur River (RMSE of 38%) compared to the Ob River (RMSE of 18%). The reduced accuracy is a product of the increasing errors in altimetric height and of the more complex relations between discharge and water level in narrower, or geometrically more complex channels, resulting in less accurate fits for the rating curves. A literature synthesis, not exhaustive, of the daily Q errors obtained using altimetry for rivers of different sizes is provided in Table 5. From these cases, one can conclude that using radar altimetry only, or in combination with other satellites or models, daily discharge accuracies of 5%–10% for the largest rivers (annual flow >1000 km3), 10–25% for large rivers (annual flow of 200–1000 km3) 30– 40% for medium size rivers (annual flow of 25–100 km3) can be anticipated. Global Earth systems models often combine monthly river discharges with climate projections to anticipate annual river flows. Our results for altimetric monthly and annual discharge demonstrated decreasing errors, comparing to daily values, if the altimetric observations coincide with the most

30 important flood events. Annually, the errors were halved for both large and small rivers (see Table 4), demonstrating the potential of radar altimetry data in climate studies.

Table 5. Accuracy of the daily discharges estimated using satellite altimetry for World Rivers of different size. River/Method

Normalizat ion method Qm dQ Qm Qm Qm Qm Qm Qm dQ Qm Qm Qm Qm Qm Qm Qm Qrms un Qm

W km3

Amazon, alti/model assimil. 5630 Amazonian rivers, RC <=5630 Amazon , RC(H+B) 5630 Amazon, Bjerklie eq. 5630 Congo, RC(H+B) 1250 Yangtze, RC(H+B) 790 Brahmaputra, RC 750 Lena R, RC 515 Zambezi, RC 504 Bjerklie eq.,S from DEM Mississippi, RC(H+B) 400 Ob, RC 406 Ob, Manning eq. 406 Ganges, RC 379 Volga, RC(H+B) 254 Niger-Benue, RC 272 Niger, RC quantile 272 Yukon, Bjerklie eq. (no winter) 202 Po, RC 46 Bjerklie eq. Pur, RC Qm 28 dQ - RMSE normalised on Q magnitude; Qm – Qrms – on RMSE of daily Q; un - unknown. From Sichangi, were selected.

RMSE, %

source

16 Emery et al., 2018 7-9 Getirana et al., 2013 1-6.5 Sichangi et al., 2016 13-31 Sichangi et al., 2016 5 Sichangi et al., 2016 7 Sichangi et al., 2016 10 Papa et al., 2012 46 Zakharova et al., 2019 6-20 Michailovsky et al., 20-73 2012 4.5 Sichangi et al., 2013 13 This study 20 This study 12 Papa et al., 2012 24 Sichangi et al., 2016 20 Tarpanelli et al., 2017 12 Tourian et al., 2017 2 Bjerklie et al., 2018 29-36 Tarpanelli et al., 2013 40-70 38 This study normalized on mean Q for study period, et al., (2016) only rivers where flood is not missed

5.3. Assessment of the methods of discharge estimation. The two methods used for estimating the discharge of the Ob River performed differently. The classical rating curve method produced lower errors (13% for daily discharge) comparing to the method based on the Manning equation (23%). The Manning formulation requires knowledge of the river depth and of the roughness parameters, which, in many cases, are difficult to determine. To estimate the river depth, different implicit methods have been tested [Leon, et al., 2006, Tourian et al., 2018], while the roughness evaluation remains largely an expert guess. To avoid guessing at the roughness, Bjerklie et al (2003) suggested a simplified equation, where only width, depth, and slope are used as predictors. Several studies compared the use of this equation with the rating curve approach [Tarpanelli et a., 2013, Michailovsky et al., 2012, Sichangi et al., 2016] and found better performance from the latter. Furthermore, Sichangi et al., (2016) and Tarpanelli et al. (2015) have

31 shown that the rating curve method is outperformed only when combined with satellite estimations of the river width. The greater accuracy of the rating curve method, based on the statistical relationship between just two parameters, can be expected when two key conditions are fulfilled: the water height satellite retrievals are sufficiently accurate and the river morphology at the gauging and the virtual stations is not radically different. The second condition ensures that the Q-H relation can be approximated by a unique rating curve (or by a limited set of rating curves) fitting the majority of cases, especially when the flow can be considered as quasi-stationary. Otherwise, dispersion of the points can either prevent an accurate fitting, result in erroneous rejection of satellite retrievals as outliers or introduce uncertainties due to oversight of a relationship branch (for example, the dumping effect of a tributary located downstream from a virtual station often produces the small loop in the curve). The Manning equation, by definition, necessarily produces less accurate discharge estimations than rating curves in cases when the above conditions are satisfied. As more variables (coming with their own errors) and more assumptions about parameterisation or the structure of equations are involved, error propagation leads to higher uncertainties in simulated discharge. Several assumptions were considered in our study. We treated the shape of the Ob river channel as rectangular. A trapezoidal approximation changes the simulated discharge by up to 17%. The use of the mean depth instead of integrating the bathymetry, in the calculation of the section areas, often produces an underestimation of the sectional area and, consequently, underestimation of the discharge. Another structural error can arise from the use of a fixed value for the roughness coefficient. We tried to overcome this problem by introducing ice roughness during the winter and floodplain roughness during the flood. This allows for seasonal variability of the effective Manning coefficient. Uncertainty can also be related to simplification of the floodplain flow. The maximal flood plain contribution to the discharge was evaluated as 4%. This contribution may have been underestimated as a result of setting a high roughness or by overestimating the water height threshold for inundation. The altimetry derived water slope for the 20 km Ob River reach represents average conditions along this reach, while the channel parameters (width and depth) used in the Manning equation pertain to a specific cross-section. Between the available cross-sections, and using the depth information provided on the topographical map, we selected the reach that suits better the requirements for instrumental discharge estimation, i.e. linear channel, absence of the islands and vegetation, regular form of the banks etc. We assumed that local

32 (high-resolution) variability of the water slope along the Ob reach was low and that our altimetric water slope was also representative of the cross-section. A comparison of the altimetric water slope retrieved over a 2 km and over 20 km reach centred at the same cross-section showed only a small difference (<1 %) between the two S retrievals (not shown). Nevertheless, one should be cautious of local slope variability when calculating discharge using the Manning equation and satellite measurements over ungauged rivers. The problem is less critical if simultaneous calibration of roughness and depth (which can absorb any potential errors arising from a slope discrepancy) is performed. Use of altimetric river heights in numerical hydrological models is a rapidly developing approach to discharge estimation. Table 5 demonstrates that this method shows good potential; generally being only slightly lower in accuracy to rating curves [Getirana et al., 2013, Biancamaria et al., 2009]. Preparations for the SWOT (Surface Water and Ocean Topography) satellite mission have included an extensive campaign for development of algorithms based on river hydraulic geometry or Manning equation simplification. An algorithm comparison study on 16 mainly medium size rivers [Durand et al., 2016] demonstrated that for 30% of cases these algorithms provide errors of <40% and for 8% of cases the errors are <20% (with the two lowest RMSE values being 5% and 6%). The best performer was a modified Manning algorithm solved for several river reaches (Metropolis-Manning). The authors argued that incorporation of various auxiliary data into the algorithms is highly beneficial and can result in significantly better discharge accuracy. In section 4.2.2 we demonstrated the positive effect of the incorporation of auxiliary data (Qinsitu) during calibration of the Manning equation. In the next section we discuss the impact of errors arising from uncertainties in variables as well as from parameterisation. We suggest a way to obtain auxiliary data on river depth in western Siberia.

5.4. Improving the accuracy of Manning-altimetry coupling discharge estimates. 5.4.1. Sensitivity of the Manning approach to errors in parametrisation. In equation (6) the mean river depth is important for estimation of section area, and consequently, through equations (7) and (8), for estimation of the water velocity. While combination of the altimetric and optical instruments permits evaluation of river width and water height, the initial section mean depth (D0) has to be obtained from other sources. The most reliable information on D0

33 is provided by in situ measurements at gauging stations. For reaches between stations, evaluation of the river depth is not a trivial task. The D0 used in this study, taken from the topographic map, is of unknown accuracy and based on measurements done 40–60 years ago. Our results indicate that even old cartography data is useful for discharge estimation. Nevertheless, the errors that could originate from uncertain river reach depths are of particular interest. In this study, river width variations obtained from satellite image analysis were compared to water height changes. The limited number of points used for definition of the B=f(Halti) function can introduce certain errors. The impact of inaccurate river widths on discharge estimates can be evaluated using a sensitivity test. The one-at-a-time approach to sensitivity analysis is an easy and established technique [Daniel, 1973]. This approach has been successfully applied in conjunction with complex hydrological models [van Griensven et al., 2006] and, in our opinion, is capable of demonstrating the impact of errors in the D0 and B approximations (eq.9) on discharge retrievals. An experiment was also performed to assess the sensitivity of the Manning method to variations in roughness coefficients and water slopes. The error range for altimetric slope retrieval was set to ± 30% reflecting an estimated mean error of 22% (which is seen as sufficiently accurate). This range, when applied to other parameters under investigation, accords with the expected natural variability of Manning coefficients in large rivers and is sufficient to cover the potential uncertainties (of ± 3.5 m) that might be introduced by using the cartographic value of river depth. As slope changes with time, each S measurement was independently randomized, within the given error range, and average scores were calculated over 2000 runs. Slope Depth,m a_m b_m n_ice n_bed

discharge errors 80% 60%

discharge errors 80%

60%

40% 20% 0% -30%

-20%

-10%

-20%

0%

10%

20%

parameter errors 30%

40%

20%

-40% -60%

0%

-80%

-30%

a)

-20%

-10%

b)

0%

10%

parameter errors 20% 30%

34

discharge errors 60%

RMSE

40%

20%

parameter errors

0% -60%

-40%

-20%

0%

20%

40%

60%

-20%

BIAS -40%

-60%

c) Figure 7. Discharge errors variation, expressed as bias (a) and RMSE (b), depending on introduced range of errors in slope, initial depth, coefficients related to river width (a_m and b_m) and bed and ice roughness in one-at-a-time sensitivity test and in overall sensitivity test (c).

Results of the one-at-a-time sensitivity test are presented in Fig. 7a and 7b. As already noted, the suggested algorithm is not sensitive to the roughness coefficient of river ice. Nevertheless, improvement of winter discharges estimations after optimisation of ni is possible. The algorithm has low sensitivity to errors in river width (lines corresponding to the coefficients am and bm have low slopes on Fig. 7a and b). This means that the time spent estimating river widths from optical images could be minimised by establishing B-H relations from a limited sample of images. The slope variable shows higher sensitivity, although slope errors within 30% do not significantly affect the daily discharges, staying within acceptable values (maximum 11% for bias and 30% for RMSE). The algorithm is more sensitive to the errors in the bed roughness parameter. Similar impacts on discharge estimation using Manning equation, from uncertainty in slope and bed roughness, have been reported for the Amazonian rivers [Sichangi et al., 2016]. The highest sensitivity found was for the initial river depth. Uncertainties of 30% in nb can lead to discharge errors of 40–60% for bias and 25–45% for RMSE, while the similar variation in initial river depth can increase Q errors up to 70–80% and 60–80% respectively. Discharge errors are lower when nb or D0 are overestimated than when underestimated. From equations 5 and 6 as well as from the Fig. 7a it is clear that bed roughness and the river depth

35 have opposite effects on discharge estimation; an underestimation of river depth could be compensated by underestimated bed resistance (or bed roughness) to the flow. In practice, if in situ discharge is available for calibration, optimisation of roughness and depth parameters at the same time produces a set of suitable pairs of nb and D0. In model optimisation routines, this is known as equifinality between parameters [Pianosi et al., 2016]. By performing this test, and using the multicriteria identification sets method, we found a direct quasi-linear relationship between these parameters in the ranges expected for this type of river (i.e., nb 0.018–0.028).

nb =0.0028×D0-0.008

eq (13)

Biancamaria et al. (2009) found that for the Lower Ob River the nb /D0 pair of 0.015/10 m produces the best river discharge simulation accuracy when using a coupled ISBA-LISFLOOD-FP model. This point lies close to the line defined by equation (13). According to this equation, at 10 m depth roughness equals 0.018. We hypothesise that if the calibration of LISFLOOD had been done with more discrete iteration steps, the optimal solution for nb /D0 in the LISFLOOD-FP model would have exactly matched our nb - D0 relation. To roughly evaluate the overall sensitivity of the Manning method, we simultaneously introduced errors in the relevant parameters using an error grid of enlarged space and reduced density (Fig. 7c). This demonstrated that to keep simulated discharge errors within acceptable values the parameters and slope errors should be less than 35–40%. If a calibration discharge set is not available, other methods of approximation of the river depth have to be implemented. Tourian et al, (2017) reviewed these methods, many of which rely on a combination of remote sensing data and empirical geomorphological relations. However, many of them, to a certain extent, use information on river discharge or river depth from alternative sources: e.g., in situ Q from another river reach, modelled discharge, in situ quintile discharge [Birkinshaw et al., 2014, Durand et al., 2008 , Tourian et al., 2017, Leon et al.,2006, Yoon et al., 2012]. Several attempts have been made to retrieve river depths using optical satellite sensors [Legleiter et al., 2009]. In seeking possible solutions, we explored the relationship between the initial river depth and the river width, hypothesising that some typical condition might apply on the West Siberian plain.

36 5.4.2. Estimation of river depth using a regional width-depth relation Leopold and Maddoc (1953) proposed a functional power relationship between river channel parameters and discharge. Since then, numerous studies have sought a connection between parameters of these functions at different hydrographical scales, from a single river reach to basin wide application [Miller, 2014]. Many works have reported reliable regional generalizations for rivers that have achieved quasi-equilibrium states [Phillips and Harlan, 1984; Singh, 2003]. In preparation for the SWOT mission, several river databases have been recently developed [Andreadis, 2013, Allen and Pavelsky, 2018]. The importance of auxiliary information on hydraulic geomorphology, which could be included in these databases, was clearly demonstrated [Durand et al., 2016]. By retrieving river width and depth from topographic maps for several lower reaches of the Ob River and adding data for other 8 rivers in north of western Siberia (62–68° N, 63–84° E), we found distinct clustering of the points on a B-D0 graph (Fig. 8). Several early studies have revealed effects on channel geometry adjustments and the B-D0 relationship from various parameters relevant to channel resistance (bed deposit diameter, vegetation) or to flow energy (slope) [Hey et al. 1986, Huang and Warner, 1995, Lee and Julien, 2006]. Examination of channel deposits made clear that the upper branch of the Fig. 8 graph, in the range of 90–500 m width, is mainly channels with higher resistance cobble-bolder deposit beds. Other points on this branch, in the width range 500–1200 m, correspond to western young secondary channel of the Ob River (called the Small Ob River and joining the study reach on the southern border) which has a bed of sand deposits; or to the upper reaches of small rivers of the region (i.e. potentially higher flow energy systems). Further analysis showed that the low branch in Fig. 8 corresponds to the condition of wide braided sand-bed reaches with important width variability, B<800 m, and negligible variation in depth. For wider braided channels, the depth variability is greater and the braided-channel line tend to converge with the main shallow sand-bedded meandering channel in a quasi-equilibrium state. Equations corresponding to each channel type are given in Table 6.

Table 6. Relation between channel width and depth at average water stage of low-level hydrological phase for rivers of the north of western Siberia. Channel type cobble channel, river upper reaches

Equation D0=0.18×B0.60

R2 0.89

RMSE,m 1.21

37 sand channels braiding channels

D0=0.08×B0.66

0.95

0.68

D0=(5.6×10-6) ×B1.88+0.83

0.98

0.70

16 14 12

D0,m

10 8 6 4 2 0 0

cobble channels

500

1000

width, m upper reaches, sand channels

1500

2000

braiding channels

2500

sand channels

Figure 8. Width – initial depth diagram for rivers of the north of western Siberia. Upper line corresponds to cobble-bolder channels and to sandy channel of the Small Ob River branch; middle line corresponds to sand-bed low meandering channels; lower line represents wide braided sand-bed reaches.

From hydraulic geometry theory, the river depth and river width could be related through an equation:

D0=c/(a (f/b))× B(f/b)

eq. (14)

where a and b are the coefficient and the power in relation B=f(Q), while c and f are the coefficient and the power in relation D0=f(Q). The hydraulic geometry coefficients found in Hey and Thorne (1986) for gravel-bed UK rivers produce following equation:

D0=0.12×B0.78

eq. (15)

38 In the 50–350 m width range the B–D0 points of equation (15) lie close to our cobble-bed line. For the Niger River, Neal et al. (2012) found similar coefficients, which best agreed with our findings at higher river widths. The measurements of B-D0 provided in [Bjerklie et al., 2018] for another big Arctic river, the Yukon, lie close to our sand-bed line and support the use of this type of relationship for estimation of initial depth. According to our sensitivity testing, the proposed relations for estimation of river discharge, in the northern part of the Western Siberian plain, can be justified for straight well-developed sand-bed river reaches wider than 800 m and for braided rivers wider than 1500 m. Under these conditions, the D0 error of 70 cm (see Table 6) is less than 10%. Consequently, the discharge errors due to depth uncertainties will not exceed 30% (see Fig. 7). For constrained river channels, for which the D0 can be estimated with RMSE of 1.21 m, acceptable Q errors of < 30% are expected for channels wider than 1000 m. These relations were developed for average water stage conditions during the low-level period; corrections might be needed if optically-derived river widths are used for D0 estimation.

6. Conclusion The accuracy of river discharge estimates from altimetric satellite observations depends on many factors. Among them are the accuracy of the water height estimates, the parametrisation of the algorithms and the availability of the auxiliary information. The accuracy of water height retrievals depends strongly on river fluvial morphology. Large river width, low banks, and the absence of the sand islands are the main favourable conditions for precise estimation of water level. The presence, at virtual stations, of sandy banks during the low-level period is less critical than islands. The accuracy (RMSE) of the H estimates from the ENVISAT mission for the wide Ob River was 0.63 m, but 1.1 m for the narrower Pur River. Recently developed approaches for water level estimation from a combination of CryoSat-2, AltiKa and Sentinel-3A altimeters provided better H accuracy (0.49 m RMSE). Two different methods of the water discharge estimation using different auxiliary information were evaluated in this study: the classical rating curves approach and the Manning formulation. The rating curve approach outperformed the Manning formulation with the errors of 13% vs 20% for daily Q and 1% vs 5% for annual Q. The combination of the rating curve method with multi-satellite water level

39 retrieval provides the most accurate discharge estimation of the Ob River; better than previous studies relying on satellite measurements, on numerical modelling or on their combination. The accuracy of the discharge estimates using the rating curve approach depends on the accuracy of the satellite water height and on the reliability of the built Q-H relationship. For the latter, two important factors can be reported: the similarity (or low difference) between the hydraulic geometry at the gauging and virtual stations and a sufficient number of observations used for curve fitting. The number of observations available for curve fitting could be constrained by the satellite repeat period. In the case of the Ob River the lowest discharge accuracy from the rating curve model was found for the ENVISAT altimeter which has the lowest sampling frequency. The Jason series altimeters, with 10-day cycles, had better performance. Finally, the multi-mission (CryoSat2/AltiKa/Sentinel-3A) provided the best estimates as a result of 1) better sampling frequency than ENVISAT, 2) the highest accuracy of the water level retrievals, and 3) the option to relocate the virtual station to the gauged reach and to use the station rating curves when the accuracy of the altimetric rating curves is low. The altimetry-based discharge accuracy decreases with reducing river size. Our results, and those of other published studies, show that, using different methods and approaches, the accuracy of altimetric daily Q estimates is 5–10% for the largest rivers (annual flow >1000 km3), 10–25% for large rivers (annual flow of 200–1000 km3) and 30–40% for medium size rivers (annual flow of 25–100 km3). The Manning formulation underpins many of the algorithms that have been recently developed in preparation for the SWOT satellite launch. Their potential advantage is in application to river reaches where simultaneous altimetric and in situ observations are absent. Our study demonstrated that the coupling of altimetry with the Manning approach permits estimation of Ob River daily Q with an accuracy of 20%. Sensitivity analysis showed that the Manning method is insensitive to the roughness parameter of river ice. Altimetry-based water slope retrievals with uncertainties less than 30% do not introduce significant errors in discharge estimations. The highest uncertainties come from the estimation of the initial river depth. Any error in parametrisation of initial depth, to achieve Q errors < 30%, has to be less than 10 %. However, they can be higher if compensated by the bed roughness parametrisation. The use of auxiliary in situ Q data for simultaneous calibration of D0 and nb can be beneficial. Initial depth information on Arctic rivers retrieved from topographic maps, and used in this

40 study, is sufficient for discharge estimations with uncertainty of 20%. We also found a functional relation between initial depth and river width. We hypothesize that this relation could be applied for rivers in the north of western Siberia for SWOT–based discharge estimation.

Author contribution: EZ designed the methods and estimated the discharge estimation KN designed the algorithm and retrieved water level and slope GK contributed to optimisation of the parameters of Manning equation AK contributed to analysis of the methods performance

References Allen G.H., and Pavelsky T.M., Global extent of rivers and streams, Science, 2018, Vol. 361, Issue 6402, pp. 585-588. Andreadis, K. M., G. J.-P. Schumann, and T. Pavelsky, 2013: A simple global river bankfull width and depth database. Water Resour. Res., 49, 7164–7168. Arcement G.J. and Schneider V.R., Guide for Selecting Manning's Roughness Coefficients for Natural Channels and Flood Plains. USGS Water Supply Paper 2339. 1989, 38p. Armitage, T. W. K., & Davidson, M. W. J. (2014). Using the interferometric capabilities of the ESA Cryosat-2 mission to improve the accuracy of sea ice freeboard retrievals. IEEE Transactions on Geoscience and Remote Sensing, 52(1), 529–536. Beltaos, S., and Kääb, A.: Estimating river discharge during ice breakup from nearsimultaneous

satellite

imagery,

Cold

Reg

Sci

Technol,

98,

35-46,

doi

10.1016/j.coldregions.2013.10.010, 2014. Beltaos, S., and T. Prowse. 2009. River-ice hydrology in a shrinking cryosphere. Hydrological Processes23: 122–144.Beltaos and Kaab (2013) Biancamaria, S., P. D. Bates, A. Boone, and N. M. Mognard, Largescale coupled hydrologic and hydraulic modelling of the Ob river in Siberia, J. Hydrol., (2009), 379(1–2), 136–150. Birkinshaw S. J., Moore. P., Kilsby. C. G., O'Donnell. G. M., Hardy. A. J., Berry. P. A. M. (2014). Daily discharge estimation at ungauged river sites using remote sensing. Hydrological Processes. 28(3). 1043–1054.

41 Bjerklie, D. M.; Birkett, C. M.; Jones, J. W.; Carabajal, C.; Rover, J. A.; Fulton, J W.; Garambois, P.A., Satellite remote sensing estimation of river discharge: Application to the Yukon River Alaska, Journal of Hydrology 561 (2018) 1000–1018. Bonnema MG, S Sikder, F Hossain, M Durand, CJ Gleason, DM Bjerklie, Benchmarking wide swath altimetry‐based river discharge estimation algorithms for the Ganges river system, Water Resources Research, 2016, 52 (4), 2439-2461 Bruner G.W., Hec-Ras, River Analysis System, Hydraulic Reference Manual, Ver.5, February, 2016. US Army Corps of Engineers, Hydrologic Engineering Center. Bykov V., Vasiliev A. Hydrometry. Hydrometeoizdat. 1965. Calmant S., Santos da Silva J; Medeiros Moreira D., Seyler F., Shum C.K., Crétaux J-F., Gabalda G., Detection of Envisat RA2/ICE-1 retracked radar altimetry bias over the Amazon basin rivers using GPS. Advances in Space Research, 2013, Volume 51, Issue 8, pp. 1551-1564 Chow V. T., Open-Channel Hydraulics, McGraw-Hill, New York, 1959. Daniel, C., One-at-a-time-plans. Journal of the American Statistical Association, 1973, 68, 353-360. Decharme B., Alkama R., Papa F., Faroux S., Douville H., Prigent C., Global off-line evaluation of the ISBA-TRIP flood model, Climate Dynamic, 2012, 38:1389–1412. Domeneghetti, A. (2016), On the use of SRTM and altimetry data for flood modelling in datasparse regions, Water Resources Res., 52, 2901–2918, Dubey A.K., Gupta P.K., Dutta S., Singh R.P., An improved methodology to estimate river stage and discharge using Jason-2 satellite data. Journal of Hydrology 529 (2015) 1776–1787. Durand, M., Gleason, C.J., Garambois, P.A., Bjerklie, D., Smith, L.C., Roux, H., Rodriguez, E., Bates, P.D., Pavelsky, T.M., Monnier, J., Chen, X., Di Baldassarre, G., Fiset, J.-M.,Flipo, N., Frasson, R.P.D.M., Fulton, J., Goutal, N., Hossain, F., Humphries, E., Minear, J.T., Mukolwe, M.M., Neal, J.C., Ricci, S., Sanders, B.F., Schumann, G., Schubert, J.E., Vilmin, L., 2016. An intercomparison of remote sensing river discharge estimation algorithms from measurements of river height, width, and slope. Water Resour.Res. 52, 4527–4549. Durand, M., K. M. Andreadis, D. E. Alsdorf, D. P. Lettenmaier, D. Moller, and M. Wilson (2008), Estimation of bathymetric depth and slope from data assimilation of swath altimetry into a hydrodynamic model, Geophys. Res. Lett., 35, L20401.

42 Emery C. M., Paris A., Biancamaria S., Boone A., Calmant S., Garambois P. A., da Silva J. S., Large-scale hydrological model river storage and discharge correction using a satellite altimetry-based discharge product, 2018, Hydrol. Earth Syst. Sci., 22, 2135–2162. Emery, C. M., Paris, A., Biancamaria, S., Boone, A., Calmant, S., Garambois, P.-A., and Santos da Silva, J.: Large-scale hydrological model river storage and discharge correction using a satellite altimetry-based discharge product, Hydrol. Earth Syst. Sci., 2018, 22, 2135-2162. Frappart, F., Papa, F., Güntner, A., Werth, S., Ramillien, G., Prigent, C., Rossow, W. B., and Bonnet, M.-P.: Interannual variations of the terrestrial water storage in the Lower Ob' Basin from a multisatellite approach, Hydrol. Earth Syst. Sci., 14, 2443-2453, https://doi.org/10.5194/hess-14-24432010, 2010. Getirana, A.C.V., Peters-Lidard, C., 2013. Estimating water discharge from large radar altimetry datasets. Hydrol. Earth Syst. Sci. 17 (3), 923–933. Gleason, C.J., Smith, L.C., Lee, J., 2014. Retrieval of river discharge solely from satellite imagery and at-many-stations hydraulic geometry: sensitivity to river form and optimization parameters. Water Resour. Res. 50 (12), 9604–9619. Guerreiro K, Fleury S., Zakharova E, Rémy F., Kouraev A., Potential for estimation of snow depth on Arctic sea ice from CryoSat-2 and SARAL/AltiKa missions. Remote Sensing of Environment, 2016., 186, 339-349. Hey, R. D., and C. R. Thorne (1986), Stable channels with mobile gravel beds, J. Hydraul. Eng., 112(8), 671–689. Huang H.Q. and Warner R.F., The multivariate controls of hydraulic geometry: a causal investigation in terms of boundary shear distribution. Earth Surface Processes and Landforms, 1995, 20: 115–130. Jain, M., Andersen, O. B., Dall, J., & Stenseng, L. (2015). Sea surface height determination in the Arctic using Cryosat-2 SAR data from primary peak empirical retrackers. Advances in Space Research, 55(1), 40–50. https://doi.org/10.1016/j.asr.2014.09.006 Kamenev G.K. Multicriteria Identification Sets Method // Computational Mathematics and Mathematical Physics, 2016, Vol. 56, No. 11, pp. 1843-1858. Kamenev G.K. A Multicriteria Method for Identification and Forecasting // Mathematical Models and Computer Simulations, 2018, Vol. 10, No. 2, pp. 154-163.

43 Kouraev. A.V., Zakharova. E.A., Samain. O., Mognard. N.M., Cazenave. A. 2004. Ob’ River discharge from TOPEX/Poseidon satellite altimetry (1992–2002). Rem. Sens. Environ. 93 (1). 238– 245. Lee, J.S. and Julien, P.Y. Downstream hydraulic geometry of alluvial channels. J. Hydraul. Eng., 2006), ASCE132(12), 1347–1352. Legleiter, C.J.; Roberts, D.A. A forward image model for passive optical remote sensing of river bathymetry. Remote Sens. Environ. 2009, 113, 1025–1045. LeFavour, G., Alsdorf, D., 2005. Water slope and discharge in the Amazon River estimated using the shuttle radar topography mission digital elevation model. Geophys. Res. Lett. 32 L17404– L17404. Leon, J.G., Calmant, S., Seyler, F., Bonnet, M.P., Cauhope, M., Frappart, F., 2006. Rating curves and estimation of average water depth at the upper Negro river based on satellite altimeter data and modeled discharges. J. Hydrol. 328 (3–4), 481–496. Leopold, L. B., and T. J. Maddock (1953), The hydraulic geometry of stream channels and some physiographic implications, U.S. Geol. Surv.Prof. Pap., 252, 56 pp. Michailovsky, C. I., S. McEnnis, P. A. M. Berry, R. Smith, and P. Bauer-Gottwein (2012), River monitoring from satellite radar altimetry in the Zambezi river basin, Hydrol. Earth Syst. Sci., 16(7), 2181–2192, doi:10.5194/hess-16-2181-2012. Miller Z. F., Pavelsky T. M., and Allen G. H., Quantifying river form variations in the Mississippi Basin using remotely sensed imagery, 2014, Hydrol. Earth Syst. Sci., 18, 4883–4895. Neal, J., G. Schumann, and P. Bates (2012), A subgrid channel model for simulating river hydraulics and floodplain inundation over large and data sparse areas, WaterResour.Res. , 48, W11506. Nielsen, K., Stenseng, L., Andersen, O. B., Villadsen, H., & Knudsen, P. (2015). Validation of CryoSat-2 SAR mode based lake levels. Remote Sensing of Environment, 171, 162–170. https://doi.org/10.1016/j.rse.2015.10.023. Nielsen, K., Zakharova E., Andersen, O. B. Compiling water level time series of non-repeat missions. ESA ArcFlux project publication, 2019. http://www.arcflux.eu/Publications/TS_model.pdf Papa, F., Bala, S.K., Kumar Pandey, R., Durand, F., Rahman, A., Rossow, W.B., 2012. Ganga-Brahmaputra river discharge. C11021. Geophysical Research 117.

44 Phillips, P. J. and Harlin, J. M., 1984, Spatial dependency of hydraulic geometry exponents in a subalpine stream. Journal of Hydrology, Vol. 71, pp. 277-283. Pianosi F., Beven K., Freer J., Hall JW, Rougier J., Stephenson DB, Wagener Th., Sensitivity analysis of environmental models: A systematic review with practical workflow, Environmental Modelling & Software, 2016, 79, 214-232. Raney, R. K. (1998). The delay/Doppler radar altimeter. Geoscience and Remote Sensing, IEEE Transactions On, 36(5), 1578–1588. Remy, F., Legresy, B., Vincent, P., 1999. New scientific opportunities from Ka-band altimetry. Geoscience and Remote Sensing Symposium. IGARSS '99 Proceedings. IEEE 1999,International vol. 1, pp. 506–508. Schwatke, C., Dettmering, D., Bosch, W., and Seitz, F.: DAHITI – an innovative approach for estimating water level time series over inland waters using multi-mission satellite altimetry, Hydrol. Earth Syst. Sci., 19, 4345-4364 Sichangi A.W., Wang L., Yang K., Chen D., Wang Zh., Li X.,, Zhou J., Liu W., Kuria D., Estimating continental river basin discharges using multiple remote sensing data sets. Remote Sensing of Environment 179 (2016) 36–53. Singh V. On the theories of hydraulic geometry, International Journal of Sediment Research, Vol. 18, No. 3, 2003, pp. 196-218. Smith, L.C., Pavelsky, T.M., 2008. Estimation of river discharge, propagation speed and hydraulic

geometry

from

space:

Lena

River

Siberia.

Water

Resour.

Res.

44.

https://doi.org/10.1029/2007WR006133. Strozzi, T. and Santoro, M. (2011): ESA DUE Permafrost - Pan-Arctic DEM V1 product guide , Vienna, Vienna University of Technology (Institute of Photogrammetry and Remote Sensing) . Tarpanelli, A., Barbetta, S., Brocca, L., Moramarco, T., 2013. River discharge estimation by using altimetry data and simplified flood routing modeling. Remote Sens. 5 (9), 4145–4162. Tarpanelli, A.; Amarnath, Giriraj; Brocca, Luca; Massari, Christian; Moramarco, Tommaso, Discharge estimation and forecasting by MODIS and altimetry data in Niger-Benue River, Remote Sensing of Environment 195 (2017) 96–106. Tarpanelli A. Brocca L., Barbetta S., Faruolo M., Lacava T., Moramarko T., Coupling MODIS and Radar Altimetry Data for Discharge Estimation in Poorly Gauged River Basins. IEEE Journal of

45 Selected Topics in Applied Earth Observations and Remote Sensing, 2015, vol. 8, issue 1, pp. 141148. Tourian, M.J.; Schwatke, C.; Sneeuw, N. River discharge estimation at daily resolution from satellite altimetry over an entire river basin. J. Hydrol. 2017, 546, 230–247. Tourian. M. J., Sneeuw. N., Bárdossy A. (2013). A quantile function approach to discharge estimation from satellite altimetry (ENVISAT). Water Resources Research. 49(7). 4174–4186. Tourian, M.J., Tarpanelli, A., Elmi, O., Qin, T., Brocca, L., Moramarco, T., Sneeuw, N., 2016. Spatiotemporal densification of river water level time series bymultimission satellite altimetry. Water Resour. Res. 52. van Griensven, A., Meixner, T., Grunwald, S., Bishop, T., Diluzio, M. and Srinivasan, R. (2006). A global sensitivity analysis tool for the parameters of multi-variable catchment models. Journal of Hydrology, 324(1-4), 10-23. van Vliet MTH, Yearsley JR, Franssen WHP, Ludwig F, Haddeland I, Lettenmaier DP, & Kabat P (2012). Coupled daily streamflow and water temperature modelling in large river basins. Hydrology and Earth System Sciences 16 (11): 4303-4321. Yoon, Y.; Durand, M.; Merry, C.J.; Clark, E.A.; Andreadis, K.M.; Alsdorf, D.E. Estimating river athymetry from data assimilation of synthetic {SWOT} measurements. J. Hydrol. 2012, 464–465, 363– 375. Zakharova. E.A., Kouraev. A.V., Cazenave. A., Seyler. F. 2006. Amazon River discharge estimated from TOPEX/Poseidon altimetry. C.R. Geosci. 338, (3). 188–196. Zakharova EA, Kouraev AV, Biancamaria S, Kolmakova MV, Mognard NM, Zemtsov VA, Kirpotin SN, Decharme B, Snow cover and spring flood flow in the Northern Part of Western Siberia (the Poluy, Nadym, Pur, and Taz Rivers), Journal of Hydrometeorology,2011, V12, 1498-1511. Zakharova, E.A., Kouraev, A.V., Rémy, F., Zemtsov, V.A., Kirpotin, S.N., Seasonal variability of the Western Siberia wetlands from satellite radar altimetry. J. Hydrol. 2014, 512, 366–378. Zakharova E.A., Kouraev A.V., Guillaso S., Garestier F., Desyatkin R.V., Desyatkin A.R. Recent dynamics of hydro-ecosystems in thermokarst depressions in Central Siberia from satellite and in situ observations: importance for agriculture and human life. Science of the total environment. 2018. 615. 290-1304.

46 Zakharova E.A., I.N. Krylenko, A.V. Kouraev, Use of non-polar orbiting satellite radar altimeters of the Jason series for estimation of river input to the Arctic Ocean, Journal of Hydrology, 2019, 568, 322-333.

Acknowledgements. This work was supported by the European Space Agency ArcFlux project of the Arctic+ ITT, by the RFBR project no. № 18-01-00465a and by RFBR project № 18-0560021-Arctic.

47 Sentinel-3 combination with CryoSat-2 and AltiKa gives the most accurate water height retrievals The rating curve method outperforms hydraulic- and model-based methods of discharge estimation River depth for Manning model can be estimated from regional depth= f(width,deposits) relations The accuracy of altimetric discharge retrievals decreases in parallel with the river size