the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Simultaneous assimilation of water levels from river gauges and satellite flood maps for nearrealtime flood mapping
Fernando Nardi
Fabio Castelli
Hydrometeo hazard early warning systems (EWSs) are operating in many regions of the world to mitigate nuisance effects of floods. EWS performances are majorly impacted by the computational burden and complexity affecting flood prediction tools, especially for ungauged catchments that lack adequate river flow gauging stations. Earth observation (EO) systems may integrate the lack of fluvial monitoring systems supporting the setting up of affordable EWSs. But, EO data, constrained by spatial and temporal resolution limitations, are not sufficient alone, especially at medium–small scales. Multiple sources of distributed flood observations need to be used for managing uncertainties of flood models, but this is not a trivial task for EWSs. In this work, a nearrealtime flood modelling approach is developed and tested for the simultaneous assimilation of both water level observations and EOderived flood extents. An integrated physically based flood wave generation and propagation modelling approach, that implements an ensemble Kalman filter, a parsimonious geomorphic rainfall–runoff algorithm (width function instantaneous unit hydrograph, WFIUH) and a quasi2D hydraulic algorithm, is proposed. An approach for assimilating multiple stage gauge observations is proposed to overcome stability issues related to the updating of the quasi2D hydraulic model states. Furthermore, a methodology to retrieve distributed observed water depths from satellite images to update 2D hydraulic modelling state variables is implemented. Performances of the proposed approach are tested on a flood event for the Tiber River basin in central Italy. The selected case study shows varying performances depending on whether local and distributed observations are separately or simultaneously assimilated. Results suggest that the injection of multiple data sources into a flexible data assimilation framework constitutes an effective and viable advancement for flood mitigation to tackle EWS uncertainty and numerical stability issues. Specifically, our findings reveal that the simultaneous assimilation of observations from static sensors and satellite images led to an overall improvement of the Nash–Sutcliffe efficiency (NSE) between 5 % and 40 %, the Pearson correlation up to 12 % and bias reduction up to 80 % with respect to the openloop simulation. Moreover, this combined assimilation allows us to reduce the flood extent uncertainty with respect to the disjoint assimilation simulations for several hours after the satellite image acquisition.
 Article
(9564 KB)  Fulltext XML
 BibTeX
 EndNote
Floods represent one of the most costly and deadly natural disasters (EMDAT, 2016), affecting annually on average more than 21 million people and producing economic loss greater than USD 100 billion (Desai et al., 2015). The ability to understand and predict floods represents a crucial aspect of river basin management strategies (Knight and Shamseldin, 2005). Numerical simulations of flood scenarios are used for proper design of structural (e.g. levees, diversion channels and dams) and nonstructural mitigation measures (e.g. land use regulations, flood zoning, flood proofing, flood forecasting and warning, disaster prevention, preparedness and response; Thampapillai and Musgrave, 1985). Early warning systems (EWSs) are nowadays increasingly used for the timely detection of flood events (Kundzewicz, 2013).
EWSs generally require integrated geospatial modelling of floodplain domains supporting integrated topographic–hydrologic–hydraulic modelling chains to produce inundation predictions (Krzhizhanovskaya et al., 2011). Digital terrain models (DTMs) and rainfall and runoff observations are required by EWSs for flood nowcasting and forecasting (Grimaldi et al., 2016). In the case of mediumterm forecast (i.e. days/weeks ahead), rainfall and runoff observations are not sufficient, and numerical weather prediction (NWP) models are required, especially for basins whose concentration time is limited so that emergency measures, such as evacuation, cannot be properly applied on time (Hopson and Webster, 2010). In this regard, recent advances in NWP models in weather forecasting were developed to adopt ensemble prediction systems (EPSs) (Buizza et al., 2005) as inputs of hydrological and hydraulic models. Flood models are computer and dataintensive applications with input data requirements (i.e. accuracy and distribution) that are often unmet, especially river flow observations (Wing et al., 2020). As a result, EWSs suffer from several occurring uncertainties associated with boundary conditions, numerical parametrizations and discretizations of floodplain features and processes (Demeritt et al., 2007). The calibration and validation of flood models for datascarce regions constitute, thus, a significant challenge for flood modellers who are often compelled to understand and manage parameter and output uncertainties (Moradkhani et al., 2005; Hostache et al., 2011; Liu and Gupta, 2007).
Data assimilation (DA) methods represent effective means of reducing these uncertainties (Cloke and Pappenberger, 2009). DA improves EWS performances by adjusting flood model parameters, input, output or state variables using available observations. DA models are used both in NWP and hydrologic–hydraulic modelling.
Advances in EPS approaches and increasing of computational power allowed the accuracy of NWP models to be improved as inputs of flood forecasting systems (Yu et al., 2016). Successful examples of advanced EPS approaches in NWP models for flood forecasting services at large scale are the EPSECMWF, from the European Centre for MediumRange Weather Forecasts (De Roo et al., 2003), and the COSMOLEPS, from the Consortium for SmallScale Modelling – Limitedarea Ensemble Prediction System (Marsigli et al., 2005).
Flood models can be updated in DA approaches by ingesting outputs of NWP models or direct rainfall–runoff observations. Stream gauge observations are the most used for updating hydrologic (McLaughlin, 2002; Moradkhani et al., 2005; Liu and Gupta, 2007) and hydraulic (Madsen and Skotner, 2005; Neal et al., 2007) model variables. However, single (or sparse) gauging stations generally fail to provide accurate flow observations during extreme events due to the distributed complex nature of flood processes (e.g. split flows, tributary junctions, overbank flow conditions and bridge overtopping). This is particularly critical for EWSs covering secondary ungauged river networks (Biancamaria et al., 2011; Mason et al., 2012). Moreover, we are observing a global decrease of the gauging stations in the river network (Stokstad, 1999; Sivapalan et al., 2003; Schumann et al., 2015).
To tackle these issues, in the last 10 years, Earth observation (EO) data have been used to inject water altimetry observations in DA frameworks for updating flood models, usually adopting radar synthetic aperture radar (SAR) technologies and 1D (Matgen et al., 2007; Neal et al., 2009; Matgen et al., 2010; Giustarini et al., 2011) or 2D (Andreadis et al., 2007; Hostache et al., 2010; Mason et al., 2012; GarcíaPintado et al., 2013; Andreadis and Schumann, 2014) hydraulic routing algorithms. One of the critical issues of the model state updating is the persistence of the improvements of the model performances. Regardless of the DA algorithm (e.g. direct insertion, particle filter (PF), EnKF), the assimilation of the model states in real and synthetic scenarios caused more accurate predictions immediately after the updating step, and they quickly decrease, depending on the specific case study, a few hours or even a few minutes after the state updating, going back to the same performances of the openloop model realization (Andreadis et al., 2007; Matgen et al., 2010; GarcíaPintado et al., 2013; Andreadis and Schumann, 2014). Some of these studies demonstrated that the updating of inflow boundaries can increase the persistence of the error reductions between the observations in both 1D (Matgen et al., 2010) and 2D (Andreadis et al., 2007; GarcíaPintado et al., 2013) hydraulic models. Other studies investigated the spatial weighting of remotesensingderived water level observations in DA approaches (Grimaldi et al., 2016). For example, Giustarini et al. (2011) found significant benefits in a local weighting procedure of assimilating unbiased, very precise water levels observations, while a global weighting procedure is recommended for water level observations in ungauged basins. However, if the local weighting is combined with poorly spatially distributed field data, the model updating can lead to an overcorrection that could even decrease the overall model performances. In fact, the frequency of the model corrections seems to be effective mostly during the rising limb of the flow hydrograph, while it seems not to be significant efficient during the recession limb (Giustarini et al., 2011; GarcíaPintado et al., 2013). GarcíaPintado et al. (2015) proposed a novel methodology to test the performance of a global formulation, a traditional local formulation and their own novel local formulation of the EnKF model to improve the forecast of a 2D hydraulic model assimilating SARderived water levels. Their novel local formulation of the EnKF was able to remove the unphysical relationships and spurious correlations that characterized the global filter. The authors also proved that the updating of the 2D hydraulic model friction and channel bathymetry seems to have a secondorder effect, with respect to the inflow updating, in flood inundation models applied to gradually varied flow in large rivers. Andreadis and Schumann (2014) applied a local EnKF for assimilating synthetic SARderived water levels, inundation width and flood extent in a 2D hydraulic model, partitioning the Ohio River (516 km) in reaches of equal lengths. The authors obtained similar results for reach lengths varying from 5 to 50 km.
Beside satellitederived altimetry, in the last years SARderived inundation extent mapping techniques have been tested to provide spatially distributed information to support nearrealtime flood detection services (Martinis et al., 2015; Pierdicca et al., 2009).
There are recent examples of DA research proving the value of assimilating satellite images for diverse purposes. In this regard, several aspects have been investigated to assimilate flood extent observations in a flood forecasting framework, such as the relation between the flood extent and the model state variable, the updating of the model inflows and parameters, the impact of the typology, the timing, the location and the frequency of the satellitederived flood extent observations on the performances of the DA performances. Lai et al. (2014) applied a variational data assimilation (4DVar) method for updating the friction (i.e. Manning the values) of a 2D hydraulic model based on shallow water equations proposing a novel cost function able to relate the satellitederived flood extent to indirectly observed flood depths. RevillaRomero et al. (2016) proposed an EnKF approach for updating the streamflow values and the parameters of a global rainfall–runoff (LISFLOOD) model using flood extent observations gathered from the Global Flood Detection System (GFDS). Hostache et al. (2018) proposed a PF approach for updating a highresolution hydraulic model directly using ENVISAT ASARderived water extents for nearrealtime flood forecasting. The authors analysed improved performances of EWSs in reducing water level estimation errors when compared to openloop (OL) simulations (i.e. not updating flood state variables with observations). Hostache et al. (2018) underlined opportunities of SAR images, overcoming visibility issues of optical sensors due to clouds, but also stressed some limitations of water altimetry approaches. In particular, the need of highresolution topographic data, challenging preprocessing and hydraulic modelling development make SARderived DA approaches hard to replicate and to be applied at varying scales (Mason et al., 2012; Wood et al., 2016). Dasgupta et al. (2021b) proposed a novel mutualinformationbased likelihood function for assimilating SARderived flood extents in a highresolution 2D hydraulic model adopting a PF approach. Dasgupta et al. (2021a) investigated the timing, the positioning and the frequency of the SARderived flood extents, in the performances of the PF assimilation of a 2D hydraulic model, finding that the optimal strategy for the image acquisition depends on the river morphology and flood wave arrival timing. Moreover, it was found that the number of observations to significantly improve the performances of the DA model increases with the narrowing of the floodplain valley.
Despite the remarkable progress in the integration of remotely sensed observations in DA frameworks, there are still major challenges (Grimaldi et al., 2016). For example, an approach able to assimilate heterogeneous observations from both local and distributed datasets coming from different sources (i.e. traditional stage gauges and remotely sensed flood extents) is still missing. Moreover, quasi2D and 2D hydraulic models can be sensitive to different simultaneous local state updating (i.e. water level corrections at specific time steps) because contiguous channel–floodplain cells can be characterized by different elevations, geometry and roughness; therefore instability issues can rise during the model state corrections with standard localization techniques. Another critical issue is that largescale flood forecasting models need to provide timely predictions, but their spatial resolution can limit the effectiveness of the assimilation of satellitederived flood extents if limited changes of water depths do not imply significant changes in flood extension and if the model does not have a sufficient resolution (Hostache et al., 2018).
In this work, a DA framework supported by heterogeneous observations coming from both local water level observations (i.e. stage gauges) and spatially distributed information gathered from satellite images is proposed and tested. This research seeks to develop a more flexible DA scheme that may value all available sources of observations for distributed flood modelling updates. The aim of this work is to mitigate flood prediction uncertainties by combining heterogeneous data and an integrated topographic–hydrologic–hydraulic modelling approach while maintaining inundation forecasting robustness, scalability and numerical stability. In achieving this goal, novel scientific advances and technical challenges of EOdriven DA approaches for flood prediction are investigated, in particular, a methodology for updating the state variable from multiple local stage gauge (SG) observations, propagating the state variable corrections instead of applying localization in a hydraulic model for distributed flood routing in floodplain domains, and the gathering of spatially distributed water level observations by means of flood extension processing and detection from satellite images, also adopting GISbased algorithms for overcoming the issues of the different resolutions between the ensembles of the flood extents retrieved from the satellitederived images and the ones generated from the hydraulic model simulations. This work conceptualizes and tests a framework for updating state variables of a quasi2D hydraulic model adopting the ensemble Kalman filter (EnKF) method to take advantage of observations gathered from heterogeneous sources. The Tiber River basin in central Italy is selected as a case study that was recently the subject of flood events at the mesoscale level (approximately 100 km^{2} of floodprone domain) to investigate on improved flood modelling performances.
The paper is organized as follows: Sect. 2 describes the adopted hydrologic, hydraulic and DA modelling methodologies. Section 3 illustrates the case study, the available data and the proposed DA implementation procedure. Section 4 discusses case study results, while Sect. 5 provides conclusive remarks underlining advantages, limitations and suggested future developments of this research.
2.1 Hydrologic and 2D hydraulic modelling
The physically based quasi2D hydraulic model (FLO2D; O'brien et al., 1993) was selected for flood wave routing and propagation. In fact, the regular grid mesh of its computational domain and the open format of the input and output files mean that the model can be integrated into a data assimilation framework. The model solves the differential form of the dynamic wave approximation of the de Saint Venant equations with a central, finite difference numerical scheme. The numerical solution is applied along the river flow path for inchannel 1D flood wave routing and for outofchannel unconfined flood propagation considering eight potential flow directions in a bidimensional (2D) domain. Channel and floodplain grid cells are assigned an absolute elevation, defining the floodplain and channel top surface topography. Channel conveyance capacity is considered by assigning a cross section to each cell with top banks associated with the corresponding floodplain cell elevation. The channel–floodplain flow exchange is simulated to take into account overbank and return flows within the riverine system. River and floodplain bridges, culverts, levees and any obstruction within the simulation domain are simulated in FLO2D by means of rating curve, width and areal reduction factors. Two main boundary conditions were defined for the application of the quasi2D model: (a) the floodplain domain extent and (b) the hydrologic forcing from upstream for both the main river stem (i.e. source node) and for the tributaries.
The upstream and tributary flow boundary conditions (b) are derived by stage gauge observations (where available) or from the application of a rainfall–runoff model (considering rainfall observations are available at river basin scale). A parsimonious hydrological model tailored for ungauged basins was selected to simulate the hydrologic forcing following Grimaldi et al. (2012) approach. This rainfall–runoff approach is based on the application of the geomorphic characterization of the instantaneous unit hydrograph (IUH) adopting the WFIUH (width function IUH) method (Mesa and Mifflin, 1986). In the WFIUH, the shape of the river basin response to the rainfall forcing is associated with rainfall drop residency time distribution. The width function (WF) distribution may be expressed by estimating the flow paths and, associated with each flow path, the travel time to reach the outlet (RodriguezIturbe and Rinaldo, 1997). The WFIUH distribution is, thus, estimated by applying channel (v_{c}(x)) and hillslope (v_{h}(x)) velocities to their corresponding flow paths L_{c}(x) and L_{h}(x):
The WFIUH distribution can be estimated using the DTM as main input information and applying terrain analysis algorithms for river basin hydrologic processing (pit removal, Jenson and Domingue, 1988; flow direction and flow accumulation; Tarboton et al., 1991) to estimate flow paths at the basin scale. Channel velocities are considered constant according to Grimaldi et al. (2012). The hillslope velocity distributions v_{h} are calculated according to NRCS (1997) as a function of the local slope and land use (Haan et al., 1994; McCuen, 2009). The adopted runoff modelling approach also considers distributed rainfall input and related infiltration losses using the SCSCN method (Cronshey, 1986). Input rain gauge observations are interpolated using the Thiessen polygon methodology to properly assess distributed rainfall input for the hydrologic model (Thiessen, 1911).
2.2 Data assimilation (DA) framework
A scheme of the whole DA framework with the reference of the related sections is illustrated in Fig. 1. The ensemble Kalman filter method (EnKF; Evensen, 2003) was selected for DA application on the proposed 2D hydraulic modelling approach. EnKF, widely used in literature for DA, was selected for its efficiency in dealing with the significant nonlinear flood dynamics (Reichle et al., 2002). The EnKF model is a sequential DA method that estimates the model state at time t+1 (h_{t+1}) based on the observations at the time steps in which they are available. The DA process is characterized by two steps: the forecast step and the updated step, whose variables will be represented by the superscript − for forecasting and + for updating. The method is based on ensemble generations: the forecast (a priori) state error covariance matrix ${\mathbf{P}}_{t+\mathrm{1}}^{}$ is approximated propagating the ensemble of the model states, according to the model errors expressed as a noise term w_{t+1}, from the previous time step; at the same time, an ensemble of observations y_{t+1} at each update time is generated according to their error distribution introducing the noise term η_{t+1}. The updated probability density function (pdf) of the model states is given by a combination of data likelihood and forecast pdf of the model states by means of Bayesian update. Specifically, the posterior estimate of the i element of the ensemble ${h}_{t+\mathrm{1}}^{i+}$ is calculated using the observation ${y}_{t+\mathrm{1}}^{i}$, performing a linear correction with the Kalman filter to the forecast state ensemble members:
where K_{t+1} is the Kalman gain matrix, H(…)_{t+1} is a propagator relating the state variables to the measured variables and providing the expected value of the output given the model state, ${v}_{t+\mathrm{1}}^{i}$ is the sample of the observation errors and ${R}_{t+\mathrm{1}}^{y}$ is the variance of the observation error.
The performance of the ensemble forecast is influenced by the spread of the ensemble Murphy (1988) and Anderson (2001) but also by the ensemble size. The size has to be sufficiently large to represent a statistically significant sample, but at the same time it has to be computational efficient considering the purpose of the application (e.g. nearrealtime updating of a flood model). In this work, the approach proposed by Anderson (2001) was selected. Therefore, the optimal ensemble size was selected to reach a normalized RMSE ratio (NRR) equal to 1.
The EnKF method application for the proposed quasi2D distributed hydraulic model was developed as follows. The state variable h_{t+1} is associated with the water depth in a specific point of the computational floodplain domain. In case the observation is a stage gauge measurement, the state variable position is determined by identifying the closest channel cell. The correction is then also applied to the closest floodplain cells and propagated upstream and downstream as illustrated in Sect. 2.2.2. In the case that the observation is gathered from a satellite image, the EnKF method is applied to both the channel and the floodplain cells for the entire computational domain as illustrated in Sect. 2.2.3, “Errors related to satellite image observations”. The model error w_{t+1} is estimated considering the uncertainties related to the input forcing I_{t+1} and the model parameters as explained in Sect. 2.2.3. The observation y_{t+1} is a water depth value gathered indirectly by the sensor. For this reason, the observation transition operation H introduced in Eq. (2) is an identity matrix, showing a direct relationship between state variables and observations.
2.2.1 Model updating
Stage gauge observations
The application of the EnKF adopting one or more point measurements as observations for updating the state variables of a physical model has been studied in depth in scientific literature. Localization is a widespread method aimed to both reduce/avoid spurious (unphysical) forecast error correlations and reduce the dimension of the state vector (thus reducing the computational time) in computationally heavy models (Anderson, 2007; Hunt et al., 2007). Many localization methods include distancebased approaches for specifying the area of influence of an observation with a userspecified distance (Ott et al., 2004; Sakov and Bertino, 2011). Alternatively, there are adaptive localization methods (Anderson, 2007; Bishop and Hodyss, 2009) aimed to remove spurious correlations if distancebased localization is critical because of the model structure (Rasmussen et al., 2015), for example using a “hierarchical ensemble filter”, where a portion of ensemble filters is used to detect sampling error (Anderson, 2007). Localization can be applied by assimilating independent local subdomains (domain localization; Ott et al., 2004) or multiplying the covariance of the forecast error by a Gaussian shaped correlation, which can be developed with a distancebased method (covariance localization; Houtekamer and Mitchell, 1998), allowing observations to be assimilated “serially” (Tippett et al., 2003). Alternatively, observation localization (Hunt et al., 2007) applies the inverse of the abovementioned distancebased correlation to the covariance of the observation error and is demonstrated to have a better performance than the covariance localization in some applications (Whitaker and Hamill, 2002). GarcíaPintado et al. (2015) are among the few cases in scientific literature in which localization (specifically observation localization) is applied to a 2D hydraulic model. The authors took into account the physical connectivity of flows proposing a novel distancebased metric approach that considers channel network distance (instead of only Euclidean distance) for weighting the covariance of the observation error and improving their model forecast skill with respect to the global filter and the traditional observation localization filter based on the Euclidean distance weight. Observation localization can be applied considering absolute water levels (with respect to the average sea level) or water depths (with respect to the terrain elevation). However, river water depths can dramatically change among contiguous cells of the hydraulic domain, for example, moving from a channel to a floodplain cell or because of changes of the local geometry (e.g. cross section shape). In fact, usually, stage gauge measurements are located under hydraulic structures such as bridges, where the geometry of the cross section (that can be reshaped to be adapted to the bridge geometry) can have large differences with respect to the surrounding natural cross sections. Therefore, the localization techniques should be better applied to absolute water levels. For this work, an observation localization technique (GarcíaPintado et al., 2015) applying a weight to the error covariance (i.e. a distance metric also based on a channel network distance) was implemented. A fifthorder polynomial alongchannel distancebased weighting function (Gaspari and Cohn, 1999) to correct the observation error covariance matrix corresponding to a local analysis domain was applied. However, even by changing the scale length of the correlation function, instability issues were encountered when updating the water levels far from the observation location, especially in those areas with a higher channel slope, because the of the changes of terrain elevation from upstream to downstream. Therefore we proposed a simplified methodology that aimed to assimilate observations at stage gauge locations and propagate water depths correction (the difference between the posterior and the forecast state variables) for the surrounding channel and floodplain cells. The alongchannel upstream and downstream water level correction is performed, applying a distancebased gain function by adopting an approach similar to Madsen and Skotner (2005):
where g(x_{k}) is the gain assigned to the k cell, A is the gain amplitude (assumed equal to 1), and the g^{′}(x_{k}) term is expressed as
where x_{obs}, x_{k}, x_{uc} and x_{dc} are the linear coordinates along the channel of the cell with the observation, the k cell to be updated, and the upstream and downstream bounds for the gain function, respectively. The last two terms allow us to consider how far the updating could be inferred to correct the flood water profile. The gain function also allows more than one observation to be injected into the DA for the same time step. The bounds of the gain for a k cell are limited by the position of the closest stage gauge cells. Figure 2 provides a scheme of the adopted channel and floodplain model updating, depicting the propagation of the gain function upstream and downstream with respect to the observation point. The same correction at the k cell is then assigned to the floodplain cells closest according to a distance measure along the flow path. Furthermore, in order to properly assimilate more than one stage gauge observation, the channel segment (and its floodplain) that falls between two different simultaneous stage observations is updated, weighting the observation values by a multiplying factor, expressed as the inverse of the distance between the observation and target channel cells. The water level correction for the k cell ΔH(x_{k}) is given by the following expression:
where Δh(x_{obs,u}) and Δh(x_{obs,d}) are the water level updates in the upstream and downstream stage gauges respectively, g(x_{k,u}) and g(x_{k,d}) are the gains relative to the upstream and downstream observation respectively, and x_{obs,u} and x_{obs,d} are the linear coordinates along the channel of the upstream and downstream cells of the observation respectively. When the gain function is propagated upstream, and the water level correction is positive, a water profile counterslope may be inferred, causing a numerical instability issue in the hydraulic model. To avoid this issue, a further condition was imposed: the absolute water level in the channel cell ${h}_{\mathrm{abs}}^{+}\left({x}_{k}\right)$, cannot be lower than the adjacent downstream channel cell ${h}_{\mathrm{abs}}^{+}\left({x}_{k+\mathrm{1}}\right)$:
The proposed simplified updating approach allows the dimension of the state vector to be remarkably reduced, considering only the locations with observations; therefore it is expected to avoid filter convergence issues while pursuing acceptable computational efficiency. The model updating procedure is invoked at each time step when one or more observations become available. The hydraulic simulation is stopped, saving distributed floodplain water levels and volume conservation to binary files. Then, the EnKF is applied, and the water depth corrections are applied to update model states in the binary files.
Satellite image observations
The assimilation of flow depths derived from satellite image processing is developed following three main steps:

flood detection from satellite image(s);

comparison of the flood extent detected from the satellite image with the ensemble of flood extents simulated by the hydraulic model;

derivation of the water depth distribution related to the satellite image starting from the ensemble of the water elevation distributions of the hydraulic model.

The proposed methodology aims to be applicable to both multispectral and SAR images to take advantage of all available satellite observations of a flood event. Considering multispectral images suffer from significant limitations due to cloud cover and light conditions, the modified normalized difference water index (MNDWI) proposed by Xu (2006) is applied. The MNDWI is expressed as
$$\begin{array}{}\text{(7)}& \mathrm{MNDWI}={\displaystyle \frac{{\mathit{\rho}}_{\mathrm{bg}}{\mathit{\rho}}_{\mathrm{bm}}}{{\mathit{\rho}}_{\mathrm{bg}}+{\mathit{\rho}}_{\mathrm{bm}}}},\end{array}$$where ρ_{bg} and ρ_{bm} are the reflectance indices of the green and midinfrared (MIR) bands respectively. For SAR images, the image histogram thresholding methodology is implemented following Brivio et al. (2002).

The satellitedetected water extension is, then, compared with the flood extension ensemble simulated by the hydraulic model (HM) at the time step of the satellite image's acquisition date. In order to avoid the impact of resolution issue on the comparison, the simulated flood raster is downscaled at the same resolution of the satellite image by following this procedure:

The water surface elevation is interpolated at the satellite image resolution by applying the Kriging method (Matheron, 1969; Oliver and Webster, 1990) using the maximum floodplain extent polygon.

The interpolated water surface elevation (WSE) is intersected with a highresolution DTM to flag positive values as potentially flooded.
The two raster flood maps (water extension from SI and from HM) are, then, quantitatively compared by applying the measureoffit F index (Horritt and Bates, 2001):
$$\begin{array}{}\text{(8)}& F={\displaystyle \frac{A}{A+B+C}}.\end{array}$$For the generic k cell pertaining to the hydraulic modelling domain, the satellitederived indirectly observed water depth ${h}_{\mathrm{o},t}^{k}$ at the time t is expressed as
$$\begin{array}{}\text{(9)}& {h}_{\mathrm{o},t}^{k}={h}_{\mathrm{m},{i}_{\mathrm{1}},t}^{k}\cdot {\displaystyle \frac{{F}_{{i}_{\mathrm{1}}}}{{F}_{{i}_{\mathrm{1}}}+{F}_{{i}_{\mathrm{2}}}}}+{h}_{\mathrm{m},{i}_{\mathrm{2}},t}^{k}\cdot {\displaystyle \frac{{F}_{{i}_{\mathrm{2}}}}{{F}_{{i}_{\mathrm{1}}}+{F}_{{i}_{\mathrm{2}}}}},\end{array}$$where ${F}_{{i}_{\mathrm{1}}}$ and ${F}_{{i}_{\mathrm{2}}}$ are the two bestfitting flood maps from the ensemble of the HM compared to the flood extent from the SI, and ${h}_{\mathrm{m},{i}_{\mathrm{1}},t}^{k}$ and ${h}_{\mathrm{m},{i}_{\mathrm{2}},t}^{k}$ are their related the flow depths of the k cell at time t.

2.2.2 Model errors
The uncertainty related to model errors is numerically managed within the proposed DA by perturbing the following:

the hydrologic forcing input given by the upstream static sensors and the rainfall–runoff modelling output;

the hydraulic model parametrization associated with channel roughness, expressed by the distributed Manning coefficients.
In both cases, the flow discharge values at time t of the s input for the i element of the ensemble are expressed using a similar approach to GarcíaPintado et al. (2013). The matrix of inflow errors is expressed as the composition of a temporally correlated error and a heteroscedastic error, whose variance is proportional to the flow value at time t. The temporally correlated error for the i ensemble member, ${q}_{\mathrm{s},t}^{i}$, evolves individually according to the expression proposed by Evensen (2003):
where ρ_{s,t} is a temporal autocorrelation coefficient, and N(0,R_{s}) is a white noise with a given variance R_{s}. The temporal autocorrelation coefficient between two time steps t and t+Δt is imposed as a function of Δt and a specified time decorrelation length τ (Evensen, 2003) as follows:
The variance R_{s} of the white noise is imposed equal to 1 (Evensen, 2003). This spatially independent value is reasonable for SGderived flows which should not include spatial correlation errors, since errors and uncertainties in SG measurements should not depend on the gauge position. On the other hand, the variance R_{I} related to an input derived from the hydrological model should depend on the distance between the locations of the other inflows, considering that the precipitation field is the main input forcing of the hydrologic model and one of the most impacting factors in flood mapping uncertainties for hydrologic–hydraulic modelling (Annis et al., 2020). Therefore, the spatial correlation ${R}_{\mathrm{I}}^{x,y}$ between two inflow errors x and y, i.e. expressed with a Gaussiandecay correlation model (GarcíaPintado et al., 2009), is
where d_{x,y} is the distance between the x and y locations, and θ is a spatial correlation coefficient.
As proposed by GarcíaPintado et al. (2013), the matrix of the heteroscedastic error is obtained as the elementwise product of the abovementioned temporally correlated error with the following factor:
where ${Q}_{\mathrm{s},t}^{\mathrm{os}}$ is the SGderived or simulated streamflow value by the s input at time t, α_{s} is the coefficient of variation related to the uncertainty of the discharge and b is a heteroscedasticity factor. Equation (13) infers the intuitive principle that high discharge values are more uncertain than low values. The resulting heteroscedastic error is then applied to the term $\mathit{\gamma}{Q}_{\mathrm{s},t}^{\mathrm{os}}$, where γ is a multiplicative bias factor.
The uncertainty related to discharge observations gathered from static sensors (SG) is the sum of two components (Clark et al., 2008): the estimation of the water level from the static sensor reading (EWL) and the conversion of the water level into discharge using the fluvial cross section rating curve (ERC). In this work, the coefficient of variation related to the static sensor was set to α_{SG}=0.1, where α_{EWL}=0.1 (Weerts and El Serafy, 2006; Clark et al., 2008; Rakovec et al., 2012b), and α_{ERC} is considered negligible with respect to the α_{EWL} (Di Baldassarre and Montanari, 2009). The coefficient of variation related to the input provided by the hydrologic model (α_{I}) can be derived from a validation analysis of the hydrologic model that calculates the distribution of the simulated flow errors. For both uncertainties related to SG and I, γ and b values were set equal to 1.
In addition to the uncertainty due to the hydrologic forcing, the uncertainty related to the channel roughness is also considered as follows (Clark et al., 2008; McMillan et al., 2013):
where ${p}_{\mathrm{s}}^{i}$ is the perturbed model parameter for the i element of the ensemble, p_{s} is the calibrated model parameter and ϵ_{p} is the fractional parameter error.
To avoid potential systematic underestimation of the model covariance error due to the limited ensemble size (overconfidence in prior estimates), the inflation method (Anderson, 2001) was implemented. The percentage of increment of the ensemble forecast anomalies can be considered as constant (Anderson and Anderson, 1999; Whitaker and Hamill, 2002) or timevariant (Ott et al., 2004), such as a ratio between a userdefined standard deviation with respect to the forecast standard deviation error (Rasmussen et al., 2015), or between the forecast and the updated forecast standard deviation (GarcíaPintado et al., 2015). In this work we adopted the first approach where, according to Evensen (2003), each i element of the forecast state variable ${h}_{t}^{i}$ is expressed as follows:
where λ is an input inflation parameter imposed equal to 1.01 (Evensen, 2003), and $\stackrel{\mathrm{\u203e}}{{h}_{t}^{}}$ is the average value of the state variable ensemble at time t.
2.2.3 Observation errors
Errors related to stage gauge observations
The errors associated with observations of the SG within the floodplain domain are considered by performing a perturbation of the observed value using a similar approach adopted for perturbing the input flow from stage gauges, as a combination of a temporally correlation error (Eq. 10) and a heteroscedastic error expressed as
where ${h}_{\mathrm{SG},t}^{\mathrm{obs}}$ is the observed water level value by the static sensor at time t. In this case, there is no error due to the rating curve application, considering the water level observations are directly compared to the simulated ones; therefore the coefficient of variation α_{SG} is assumed to be equal to 0.02 m (Schmidt, 2002; Pappenberger et al., 2006).
Errors related to satellite image observations
The procedure adopted for deriving water depth distributions from satellite images is affected by a series of errors that must be taken into account and in particular the following:

Error in the water surface detection from satellite images. This error is due to the water detection technique that could overestimate or underestimate the flood extension. Both multispectral and SAR image processing for water extent mapping require a threshold to apply in the water index and the backscatter coefficient respectively. Literature values of these thresholding values could lead to inaccuracies considering optimal threshold values are usually casestudy or eventspecific; therefore, a perturbation of the threshold value is performed by adopting a normal distribution with zero mean and a standard deviation derived from literature values (Pierdicca et al., 2009).

Error of the water surface extraction from the simulated WSE of the hydraulic model. This error is mainly due to the vertical error of the DTM that is used in the water surface elevation interpolation procedure. The generic iDTM of the ensemble is perturbed by generating a vertical error with a normal distribution characterized by a zero mean and a variance that is uniformly distributed between 0 and 0.3 m, U(0,0.3), according to literature values (Hodgson and Bresnahan, 2004; Leon et al., 2014; Brouwer et al., 2017). Considering the proposed normally distributed independent errors do not take the spatial continuity of the elevation data into account (Raaflaub and Collins, 2006; Heuvelink et al., 2007), a GIS algorithm for inferring spatially autocorrelated errors is applied. A correlation distance error (CDE) equal to 100 m is applied according to Li et al. (2011), Livne and Svoray (2011), Mudron et al. (2013) and Leon et al. (2014). The abovementioned GIS algorithm includes the following steps: (1) generation of a raster (NR) of random values with a normal Gaussian distribution (μ=0, s=1) for the entire extension of the DTM, (2) generation of a raster (SR) with the average of the NR values within a neighbourhood equal to CDE, (3) creation of a error distribution raster (Err) that divides the SR raster by its spatially averaged standard deviation and multiplies the result for the adopted variance U(0,0.3) and (4) addition of the Err raster to the original DTM.

Error of the water depths derived from the ensemble of hydraulic modelling. Equation (9) assumes a linear relationship between water depth values of two hydraulic profiles and the weight associated with their relative F indexes, expressing the comparison with the observed water extension from SI. The application of this weighted mean of the simulated water depths could lead to an inaccuracy in the vertical estimation of the water depths, especially at the boundaries of the two different simulated flood extents. The perturbation error due to the profile derivation for the i element of the ensemble and the generic k cell is expressed as a random uniform noise:
$$\begin{array}{}\text{(17)}& {\mathrm{err}}_{\mathrm{PD}}^{i,k}=U\left(c\cdot \mathrm{\Delta}{h}_{\mathrm{12}}^{k},+c\cdot \mathrm{\Delta}{h}_{\mathrm{12}}^{k}\right),\end{array}$$where $\mathrm{\Delta}{h}_{\mathrm{12}}^{k}$ is the water level difference at the k cell of the two bestfitting hydraulic simulations (see Eq. 9), and c is a coefficient ranging between 0 and 1, considering that the gentle terrain slopes in floodplains limit the error of water depths derivation in an interval smaller than $\mathrm{\Delta}{h}_{\mathrm{12}}^{k}$.
3.1 The Tiber River in central Italy
The selected case study is represented by the Tiber River upstream of the city of Rome (Fig. 3). The fluvial transect goes from the village of Orte Scalo to the northern boundary of the city of Rome corresponding to the Castel Giubileo dam. The entire floodplain domain of the Tiber Orte–Castel Giubileo transect has an extension of 5881 km^{2}, with a main tributary represented by the Nera River (drainage area of 4180 km^{2}) and 15 minor ungauged tributaries draining into the selected fluvial domain. The Tiber river at the upstream Orte boundary section has a drainage area of 8400 km^{2}, while at the downstream end of Castel Giubileo the drainage area is 14 850 km^{2} (total Tiber River basin catchment area at the Tyrrhenian sea outlet is approximately 17 400 km^{2}). The floodplain domain is mostly characterized by agricultural use, but major road and railway infrastructures were developed to connect several urbanized areas along the Tiber floodplain with the four main towns of Orte Scalo, Fiano Romano, Monterotondo and the northern part of the city of Rome that have been subject to frequent floods in January 2014, November 2012, November 2010 and November 2005, causing damage to buildings, roads and bridges. This floodplain domain also has a strategic importance for the flood risk mitigation of the city of Rome, considering flood volume accumulation in this domain determines a significant flood peak attenuation that propagates through the historical city centre. Understanding, monitoring and predicting flood scenarios in this fluvial domain is crucial for protection of the socioeconomic and cultural assets of the Italian capital city. The city of Rome EWS strictly relies on the flood modelling predictions of the selected area.
3.2 Parametrization of the flood forecasting model
3.2.1 Topography and hydrologic modelling
Topographic data to represent the morphology of the selected Tiber river subbasin domain were gathered from the TINITALY 10 m resolution DTM (Tarquini et al., 2012) for supporting the hydrologic modelling. Rainfall time series for rainfall–runoff modelling were gathered from 94 rain gauges with a temporal frequency ranging from 1 to 15 min. SCS infiltration method parametrization used fourthlevel CORINE Land Cover dataset gathered from the Istituto Superiore per la Protezione e la Ricerca Ambientale (ISPRA) repository, with ancillary data for the lithology and permeability layers gathered from Autorità di Bacino Distrettuale dell'Appennino Centrale. The river basin terrain analysis procedure needed to provide WFIUH hydrologic modelling input parameters used a value of 1 km^{2} to define stream network source cells, a constant parameter that adequately matched the fluvial network extension observed from aerial images of the basin. WFIUH kinematic parameters were calibrated using four small gauged basins (Naja, Niccone, Puglia, Sovara) estimating channel flow velocities equal to 2 m s^{−1} and distributed hillslope flow velocities in the range 0.01 to 0.1 m s^{−1}.
3.2.2 2D hydraulic modelling
The bathymetry needed to represent the channel conveyance capacity in the hydraulic model was derived by interpolating available surveyed cross sections. The continuity of the channel–floodplain topographic domain was obtained using the available highresolution lidar DTM (1 m resolution), gathered from Ministero dell'Ambiente e della Tutela del Territorio e del Mare, that covers most of the selected Tiber river floodplain area. Where lidar was not available, the 5 m resolution DTM from Regione Lazio was used. The floodplain grid cell resolution used for 2D flood unconfined flow routing was set equal to 150 m to produce computationally efficient model runs. The consistency of the coarseresolution hydraulic model was validated by matching flood simulations with available observations from three recent inundation events (November 2005, November 2010, November 2012). Available flood observations were also used to calibrate the roughness parametrization. A constant Manning coefficient for the channel equal to 0.04 m${}^{\mathrm{1}/\mathrm{3}}$ and variable Manning parameters for the floodplain, classified using CORINE Land Cover classes, were calibrated. The 2D hydraulic model validation was performed by comparing simulations with available flood stage time series from seven gauges that were used as floodplain control sections. This also allowed the consistency of the numerical representation of major urban features of the study domain to be verified. Flow and stage rating tables available for the main bridges and weirs gathered from the Centro Funzionale regionale del Lazio were used to finalize the 2D hydraulic model validation.
3.3 DA implementation
The ensemble size for the application of the EnKF model was defined, according to Anderson (2001), by analysing the available stage gauges and the selected flood events. The optimal ensemble size was set equal to 40.
To define the temporal autocorrelation error (Eq. 11), a time decorrelation length τ equal to 3 d was imposed, while a spatial correlation coefficient θ was assumed equal to 60 km. GarcíaPintado et al. (2013) considered these values as representative of a spatially distributed or semidistributed model that ingests continuous rainfall field inputs after being calibrated with previous flood events.
The hydrologic model is affected by different sources of uncertainties: the structural uncertainties, given by the simplification of the modelled physical processes (e.g. we adopted a WFIUH approach, neglecting groundwater flow, mud and debris flow), the input uncertainties (given by the rainfall values and antecedent soil moisture conditions) and parametric uncertainties due to the inaccuracy of the model calibration). These sources of uncertainty should be considered separately. For example, input rainfall uncertainty from rain gauges can be estimated considering quantitative precipitation ensembles (Clark and Slater, 2006), such as sequential Gaussian simulations (Goovaerts, 1997; Rakovec et al., 2012a). Precipitation ensembles generated with NWPs can then be coupled with hydrologic models to improve flood forecasting (Jasper et al., 2002; Sorooshian et al., 2008). In this work we adopted a simplified procedure taking into account all the modelling uncertainties considering the frequency distribution of the errors between the observed and simulated flow values obtained by the calibration and validation of four small tributaries of the Tiber River basin in past flood events.
From the validation of the hydrologic model, the frequency distribution of the relative flow errors was characterized by a mean equal to zero and a standard deviation equal to 0.28⋅Q^{os}; thus α_{I}=0.28 (Eq. 13).
The application of Eq. (14) to consider channel roughness parametrization uncertainties resulted in a p_{s}=0.04 m${}^{\mathrm{1}/\mathrm{3}}$ s^{−1} (according to the hydraulic model calibration), with ϵ_{p} assumed equal to 0.125. This limits the Manning channel value variations between a minimum of 0.035 and a maximum of 0.45 m${}^{\mathrm{1}/\mathrm{3}}$ s^{−1}. The floodplain roughness uncertainty was considered less significant considering that the governing factor for the Tiber floods in the selected domain is the volume, while minor specific urban features and singularities characterized the selected flood events. It is also noted that, for the selected events, the flow is conveyed by the channel for most of the simulation time.
A Landsat 7 image (acquisition date: 14 November 2012 – 09:43 LT) was processed using Eq. (7) to extract the observed flood extension. Landsat 7 products are affected by minor corruptions due to a failure of the satellite scan line corrector (Scaramuzza and Barsi, 2005) that creates a data void (i.e. empty stripes) in the water mask. These irregularities were analysed and interpolated, allowing us to overcome this issue and define a correct delineation of the flood extension, clearly visible from the available Landsat 7 image. Figure 4 shows the resulting flood detection map. The Landsat 7 image extension covers only a portion of the selected study domain. The satellite image acquisition time is consistent with peak flow conditions within these two gauging stations (Fig. 4). A threshold value equal to 0 for MNDWI was chosen (Xu, 2006). To apply the EnKF method, the variance of the random noise related to the threshold value of MNDWI was set equal to 0.2.
Three sets of DA scenarios were selected: (1) stage gauge observations, (2) satellite image observations, and (3) both SG and SI. In the case of assimilation from SGderived observations, two subscenarios were implemented: the simultaneous assimilation of four stage gauges (ASS 4SG) and the assimilation of one upstream stage gauge, Ponte Felice (ASS 1SG). These latter scenarios were introduced to analyse how the model performance varies downstream far from the observation location.
The 2005, 2010 and 2012 floods were selected for applying the DA framework using stage gauge observations. In the case of scenarios 2 and 3, the Landsat 7 image, described in Sect. 3, was assimilated for the 2012 flood event.
4.1 Assimilation of stage gauge observations
Figure 5 shows the comparison of observed and simulated water level time series at Stimigliano, Nazzano and Ponte del Grillo stations (16.6, 48.2 km and 59.1 km away from the upstream Ponte Felice station respectively) for the 2005, 2010 and 2012 flood events. The first two flood events are characterized by multiple peaks, whose rising and recession curves are not properly represented by the openloop (OL) simulation. This limitation is probably due to the coarse resolution of the flood model. In fact, wetting and drying phenomena along preferential flow pathways are usually influenced by the microtopography of the domain and are better represented in higher resolution models (Nicholas and Mitchell, 2003; Neal et al., 2011). The simultaneous assimilation of four stage gauges (ASS 4SG) is able to overcome this issue, and the spread of assimilated ensemble at each stage gauges is much lower than the one of the OL simulation because of the small error associated with the stage gauge observation as illustrated in Sect. 2.2.3, “Errors related to stage gauge observations”. On the other hand, the assimilation of the single upstream stage gauge (ASS 1SG) provides a slight benefit, with respect to the OL simulation, at the closest stage gauge (Stimigliano), mostly for high values of flow depths (e.g. at the peak of the 2012 flood event), while it does not imply any substantial changes downstream, where local inflow conditions and terrain geometry completely attenuate the upstream water level corrections. Table 1 shows the performance indexes for OL, ASS 4SG and ASS 1SG simulations considering the stage gauge observations as true values, given their low uncertainty. Bias is expressed as the ratio of the sum of the observed and simulated water levels. The scenario ASS 4SG improves the prediction of the water levels in terms of Nash–Sutcliffe efficiency (NSE), R and bias performance indexes (Table 1), while the performances related to the scenario ASS 1SG decrease gradually downstream with respect to the assimilated gauge observation. Note that the ASS 1SG performances are even slightly worse than OL in Nazzano and Ponte del Grillo. This is due to the fact that propagation of corrections related to faraway observation locations can be counterproductive (Giustarini et al., 2011). Bias for OL and both ASS simulations tends to increase above 1 in the downstream stage gauges because of their overestimation of the water levels, especially in the recession limb (that is less important for EWS). For the same reason, for all the simulations, the further the flow is from the upstream inflow, the more the R coefficient tends to decay, but in the case of the ASS 4SG simulation, this decay is mitigated.
Figure 6 shows the channel water depth profiles for three different time steps (correspondent to three different position of the peak flow along the channel) for ASS 4SG (left panels), ASS 1SG (right panels) and OL (both right and left panels) simulations for the November 2012 flood event. The ASS 4SG scenario shows improvements in terms of reduction of ensemble spread in most of the channel domain with respect to the OL simulation, even if the adopted gain function illustrated in Eq. (3) attenuates the correction mostly in the downstream part of the domain, far from the flow observations. Note the assimilated water depths right upstream of the Stimigliano stage gauge (second black dot from left) because of the propagation of the water depth correction of the upstream stage gauge (Ponte Felice). In the ASS 1SG scenario, the corrections of the water depths a respect to the OL simulation are almost negligible in the downstream part of the domain.
Spatial covariances of the assimilated flow depths (represented in Fig. 7 for a portion of the computational domain) increase with simulation time during the rising limb for both assimilation scenarios (ASS 4SG and ASS 1SG), but ASS 4SG shows lower values of covariance with respect to the ASS 1SG simulation, especially close to the downstream SG location. This confirm that the quasi2D hydraulic model is able to improve its performances mostly if multiple stage gauge observations are assimilated, since propagation of corrections is attenuated by downstream inflows and local geometry.
Figure 8 shows the impact of water level updating model ASS 4SG on the simulated flood extent. In Fig. 8b a simulated flood extension subset is shown considering the mean levels of the OL and ASS 4SG simulations for the November 2012 event. The ASS SG flood extent simulations are on average 6.5 km^{2} larger than the OL simulations. Major differences are located in the flat areas where flood extents are very sensitive to flow depth variations. In this figure, the resolution of the flood extension of the hydraulic model was increased interpolating the water surface elevation at the 5 m resolution DTM by applying the Kriging method and filtering out the cells with lower resolution a respect to the 5 m resolution DTM. This method can be, thus, accepted in the case of long persistence (e.g. several hours) of water levels at the same location, as in the selected case study. Inaccurate flood extent mapping is expected for small basins with low flood persistence. The application of this methodology for smaller basins should require a higher flood model and DTM resolution.
4.2 Assimilation of the satellitederived flood extent observations
Figure 9 shows the observed and simulated flood hydrographs at Stimigliano, Nazzano and Ponte del Grillo stations for both the OL and ASS models for the 2012 flood event. The updated mean water levels at the SI acquisition time are slightly higher and stay higher than the OL simulation for a few hours. The spread of the ensemble of the ASS simulation is significantly reduced in correspondence of the SI observation. This reduction is gradually damped until it completely disappears in approximately 8 h. The improvement of simulation NSE and bias for the SI assimilation case is not significant (Table 2), considering that the updating persists for only few hours, as shown in Fig. 9.
Figure 10a shows the reduction of the water level uncertainty and the correction of the mean value along the channel profile at the SI acquisition time. Note that the proposed methodology for gathering the indirectly observed satellitederived water depths allowed observations to be made at each cell; therefore the EnKF is applied serially to the whole domain with positive depths values instead of only at the SG locations, thus avoiding increase of the water depth ensemble spread in cells far from the SG locations as shown in Fig. 6. Reduction of the spatial covariances right after the SI acquisition time with respect to the OL simulation are shown in Fig. 10b and c. Covariances for OL are much lower than ASS SI, mostly along the channel cells. Covariances of ASS SI simulation are also lower than the ones generated by the ASS 4SG and ASS 1SG observations at the same step (see the right panels of Fig. 7), especially far from the SG locations. The adopted updating procedure allows the flood extent of 4 km^{2} to be increased at the time of the SI acquisition (see an inset of the flood extension in a flat area of the floodplain domain in Fig. 11), leading the false negative rate to be reduced by 7 % and the false positive rate to be increased by only 1 %. Note that the satellitederived flood extent is considered a true flood map. Despite the smaller water level observation errors given by the stage gauge observations with respect to the SI, the overall change of the mean flood extent between the OL and ASS SI simulations (14 %) does not show significant differences as respect the stage gauge DA (9 %).
Figure 12 shows the variability of bias, root mean square errors (RMSEs) and standard deviation of the ensembles calculated, starting from the time of the acquisition of the satellite image and comparing OL and ASS modelling results at Stimigliano, Nazzano and Ponte del Grillo stations. Improvements in terms of bias and RMSE are significant for 20 h after the SI acquisition, while the uncertainty reduction (i.e. the difference in the ensemble amplitude between the OL and ASS SI simulations) persists for 8 h. This behaviour suggests that observations right after the SI assimilation step (e.g. the SG observations) could benefit from the reduction of the model uncertainty for the whole computational domain.
4.3 Assimilation of both stage gauge and satellite observations
At the positions where the channel flow gauges are located, it is reasonable to expect that stage gauge observations support better model performances as compared to the model updated based on satellite observations (see Figs. 5 and 9). However satellite image observations provide spatially distributed information that is important for flood wave routing in a complex domain, especially when flooding impacts large unconfined domains. Obtaining multiple stage gauge observations and satellite images constitutes the optimal support for DA application for EWSs.
To investigate the potential benefit of taking advantage of multiple heterogeneous distributed data sources of flood observations, we simulated the simultaneous assimilation of the four SG observations and the SIderived flood extent (ASS TOT). As expected, the water level time series at the stage gauge locations are similar to the ones of the disjoint ASS 4SG simulation (bottom panels of Fig. 5) in terms of mean and spread of the ensemble (therefore plots of hydrographs at SG are not shown for the joint assimilation simulation). Specifically, improvements in performances with respect to the OL simulation are observed in terms of increase of NSE (5 %–40 %), increase of Pearson correlation (up to 12 %) and bias reduction up to 80 % (see Table 3).
Conversely, differences with respect to the 4SG assimilation are found far from the SG locations in terms of ensemble spread of channel water levels and spatial covariances (Fig. 13) at the SI the acquisition time (Fig. 13a and c) and 2 h later (Fig. 13b and d). The spread of the ensemble of the water levels along the channel profile is reduced with respect to both the ASS 4SG assimilation and the ASS SI assimilation. Moreover, covariance values of ASS TOT are smaller with respect to the ones of ASS SI and ASS SG simulations for the same time step (see Figs. 7c, 10c and 13c). In fact, the combination of the model pdf right before the SI acquisition (that is already narrowed by the previous SG assimilation) with the pdf of the SIderived observed water depths allows a water depth ensemble to be generated with a lower spread with respect to both the disjoint SI and SG assimilations. The joint assimilation of SG and SI observations (ASS TOT) also has a positive impact on the flood extent uncertainty (see Fig. 14) when compared with the OL and disjoint SG and SI assimilation. Specifically, the ASS TOT simulation led to a reduction of the maximum flood extent variability for several hours after the SI acquisition with respect to the ASS SG and ASS SI simulations. Figure 14 also shows a consistent difference of the flood extent ensemble spread between the ASS 1SG and the ASS 4SG simulations, confirming that the propagation of the assimilation of a single stage gauge has a spatially limited effect (only for the cells around the SG location), while the joint SG and SI assimilation improves the performances after the SI acquisition time for almost 15 h.
4.4 Pros and limitations
The proposed modelling chain approach is affected by some limitations but also advantages that are summarized in this section. Firstly, the adopted hydraulic model has a coarse spatial resolution (150 m cell size), and its performance can be considered acceptable for highmagnitude flood events; however some limitations could arise in representing the flow patterns for lowmagnitude events, where microtopography can have an important role in the flow propagation (Bates, 2012). Further tests considering a smaller domain with a higher resolution 2D hydraulic model are needed to verify the stability of the model updating, when water level corrections are applied to cells characterized with smaller dimensions. Moreover, the model uncertainty can be considered quite large during the peak flows of the selected flood events with an amplitude of the ensemble of the water levels even equal to 2 m (see Figs. 9 and 10). More accurate models could not benefit from the assimilation of the satellite image if the indirect water level distribution derivation is affected by uncertainties larger than the ones related to the forecasting model. Tables 1–3 also show some bias in the hydraulic model, mostly prominent in the recession limb of the hydrographs. This is a limitation, since the optimality of the DA techniques is realized if the observations and the models are not biased (Dee, 2005; Liu et al., 2012). Bias reductions and further improvement of the simulation can be done by updating the model inputs (i.e. the inflow hydrographs) and parameters (e.g. the roughness) using an augmented state vector approach in the EnKF framework (Montzka et al., 2012).
The simplified rainfall–runoff modelling allowed input flow hydrographs to be generated very quickly according to the needs of a nearrealtime flood modelling purpose. However, the model can only be considered appropriate for small basins characterized by an impulsive response, for which the groundwater component can be neglected and complex topography and flow control structures are absent to avoid equifinality issues during the calibration/validation analysis (Beven, 2006). Furthermore the application of the SCSCN model at subdaily timescale (Grimaldi et al., 2013) is a strong limitation, and more advanced models should be preferred to reduce the model uncertainties.
Despite the several measures adopted to prevent instability issues, instability can occur during the updating of the water levels from the stage gauges. For this reason, the model is tailored to remove the critical elements causing instability from the ensemble and to generate new elements in order to keep the sample size constant. This measure can slow down the model, that should be as fast as possible for a proper nearrealtime application. The instability issues that sometimes can occur could be due to the fact that the FLO2D model does not allow the flow velocities to be updated, rather only the flow depths. This limitation also does not enable control of the volume conservation, which is an important factor to verify the accuracy of the model simulation. The need of updating the model states each time the stage gauge observations are available can affect the efficiency of the model in terms of computational performance. The proposed model requires averagely 3.7 min for each simulation hour with a standard laptop (processor: 4 core with 180 GHz each and 8 GB RAM). Alternative approaches such as the asynchronous ensemble Kalman filter (AEnKF; Sakov et al., 2010), allowing past observations to be ingested over a time window to update model state at a specific time step, could help to reduce the time when the model is updated, even if each model updating requires slightly higher costs concerning the computational time and storage requirements.
A methodology for indirectly deriving the distribution of the water depths from the water footprint gathered from a satellite image was proposed; this methodology is affected by a series of errors that were taken into account for assigning a proper reliability to the observation related to the satellite image. This reliability is numerically lower than the one related to stage gauge observations but, at the same time, provides distributed information instead of the observations given by the static sensors. The derivation of the water depths from the flood extent gathered by the satellite image was performed with a linear combination of the values given by the ensemble of the results provided by the hydraulic forecasting model. Since the latter has to be updated by an observation that is indirectly derived by the model itself, this approach can be considered disputable; however, practically it was demonstrated not to cause instability issues during the model updating, since the distribution of the flow depths is coherent with the model state variable, and it has been demonstrated to slightly improve the model performance. This approach can be considered as a hybrid methodology of two literature approaches that consider prognostic and diagnostic variables for assimilating satellitederived information (Hostache et al., 2018), taking advantage of the rapid flood detection algorithms adopted for direct assimilating the flood extent and, at the same time, directly retrieving the water levels that are prognostic variables, thus more straightforward to assimilate than the flood extent (Hostache et al., 2018).
During a flood event, the adoption, as an observation, of a multispectral image, potentially corruptible by the cloud covering and sunlight, is much less likely than a SAR image; however, the proposed approach for the model updating can be applicable, regardless of the type of image as input observation. The satellite revisit time of the current satellite missions is still a strong limitation, since it can be much higher than the travel time of small basins. Furthermore, usually multispectral SAR images require time to be processed. However, new satellite missions and also the combination of more constellations will considerably reduce the revisit time in the near future, allowing different images to be made for the same area with a temporal frequency higher than the persistence of the model correction (8 h). Moreover, recent automatic satellite image techniques for extracting the flood extension have been implemented in realtime services for flood mapping (Martinis et al., 2015).
The proposed DA framework investigated the opportunities and challenges of assimilating multiple sources of observations for improving the performances of nearrealtime flood predictions in the case of some real flood events selected as case studies. Specifically, stage gauge water level readings and satellitederived flood extents were used for testing the proposed DA framework. We infer the following main conclusions:

The assimilation of multiple (four) stage gauges significantly improved the flood model performances in terms of NSE and R bias and also reduced simulated inundation extent uncertainties; however the spatial influence of assimilation is limited if only one SG observation is adopted.

The assimilation of distributed flood depths, indirectly retrieved from satellitederived flood extent, provided slightly better modelling performances only a few hours after the SI acquisition and allowed the water level and flood extent uncertainties to be reduced for the whole computational domain.

The simultaneous assimilation of SG and SI observation enabled a reduction of the water level and flood extent uncertainty for several hours with respect to the disjoint SI and SG assimilation after the SI acquisition time.
Future tests are needed that take advantage of the increasing availability of satellitederived flood extents at higher temporal and spatial resolution, to discover the effective capacity of the proposed flexible multisource DA framework and to value a larger EO data availability.
Moreover, the flexibility of the proposed model to assimilate local and distributed observations suggests the suitability of using other data sources, also gathered from informal observation systems. In particular, future work will allow us to investigate the use of crowdsourced observations to apply the proposed flood prediction framework in ungauged basins. Crowdsourcing has already proved to be an effective means to improve hydrologic (Mazzoleni et al., 2017) and hydraulic (Annis and Nardi, 2019) modelling performances, but further work is needed to test the use of crowdsourced data in realtime flood modelling approaches. Significant improvements are expected in the near future for improved weather predictions by valuing all data sources, especially affordable citizendriven observations for developing countries that are the most vulnerable to hydroextremes (Alley et al., 2019).
FLO2D modelling software is available at https://flo2d.com/ (FLO2D, 2022).

The TINITALY DEM (10 m resolution) data were downloaded from the INGV website (https://doi.org/10.13127/TINITALY/1.0, Tarquini et al., 2007).

The regional DTM (5 m resolution) was directly provided by Regione Lazio.

The lidar (1 m resolution) was directly provided by the Ministero dell'Ambiente e della Tutela del Territorio e del Mare.

The CORINE Land Cover data were downloaded from the Istituto Superiore per la Protezione e la Ricerca Ambientale.

Historical rainfall and runoff data from the Tiber River basin were directly provided by Regione Lazio.

Landsat images were downloaded from the USGS Earth Explorer website: https://earthexplorer.usgs.gov/ (USGS, 2022).
AA conceptualized the framework with the supervision of FN and FC. AA developed the model code and performed the simulations. AA prepared the manuscript with contributions from all coauthors.
The contact author has declared that neither they nor their coauthors have any competing interests.
Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Fernando Nardi and Antonio Annis acknowledge the support received by the WARREDOC center of University for Foreigners of Perugia through the WARREDOC–Fondazione ENI Enrico Mattei (FEEM) research agreement. Fernando Nardi acknowledges the support received by the Southeast Environmental Research Center in the Institute of Environment at Florida International University.
This research has been supported by the Ministero dell'Ambiente e della Tutela del Territorio e del Mare (SimPro).
This paper was edited by Patricia Saco and reviewed by Maurizio Mazzoleni and two anonymous referees.
Alley, R. B., Emanuel, K. A., and Zhang, F.: Advances in weather prediction, Science, 363, 342–344, 2019. a
Anderson, J. L.: An ensemble adjustment Kalman filter for data assimilation, Mon. Weather Rev., 129, 2884–2903, 2001. a, b, c, d
Anderson, J. L.: Exploring the need for localization in ensemble data assimilation using a hierarchical ensemble filter, Physica D, 230, 99–111, 2007. a, b, c
Anderson, J. L. and Anderson, S. L.: A Monte Carlo implementation of the nonlinear filtering problem to produce ensemble assimilations and forecasts, Mon. Weather Rev., 127, 2741–2758, 1999. a
Andreadis, K. M. and Schumann, G. J.: Estimating the impact of satellite observations on the predictability of largescale hydraulic models, Adv. Water Resour., 73, 44–54, 2014. a, b, c
Andreadis, K. M., Clark, E. A., Lettenmaier, D. P., and Alsdorf, D. E.: Prospects for river discharge and depth estimation through assimilation of swathaltimetry into a rasterbased hydrodynamics model, Geophys. Res. Lett., 34, 1–5, 2007. a, b, c
Annis, A. and Nardi, F.: Integrating VGI and 2D hydraulic models into a data assimilation framework for real time flood forecasting and mapping, GeoSpat. Inform. Sci., 22, 223–236, 2019. a
Annis, A., Nardi, F., Volpi, E., and Fiori, A.: Quantifying the relative impact of hydrological and hydraulic modelling parameterizations on uncertainty of inundation maps, Hydrolog. Sci. J., 65, 507–523, 2020. a
Bates, P. D.: Integrating remote sensing data with flood inundation models: how far have we got?, Hydrol. Process., 26, 2515–2521, 2012. a
Beven, K.: A manifesto for the equifinality thesis, J. Hydrol., 320, 18–36, 2006. a
Biancamaria, S., Durand, M., Andreadis, K., Bates, P., Boone, A., Mognard, N., Rodriguez, E., Alsdorf, D., Lettenmaier, D., and Clark, E.: Assimilation of virtual wide swath altimetry to improve Arctic river modeling, Remote Sens. Environ., 115, 373–381, 2011. a
Bishop, C. and Hodyss, D.: Ensemble covariances adaptively localized with ECORAP. Part 1: Tests on simple error models, Tellus A, 61, 84–96, 2009. a
Brivio, P., Colombo, R., Maggi, M., and Tomasoni, R.: Integration of remote sensing data and GIS for accurate mapping of flooded areas, Int. J. Remote Sens., 23, 429–441, 2002. a
Brouwer, T., Eilander, D., van Loenen, A., Booij, M. J., Wijnberg, K. M., Verkade, J. S., and Wagemaker, J.: Probabilistic flood extent estimates from social media flood observations, Nat. Hazards Earth Syst. Sci., 17, 735–747, https://doi.org/10.5194/nhess177352017, 2017. a
Buizza, R., Houtekamer, P., Pellerin, G., Toth, Z., Zhu, Y., and Wei, M.: A comparison of the ECMWF, MSC, and NCEP global ensemble prediction systems, Mon. Weather Rev., 133, 1076–1097, 2005. a
Clark, M. P. and Slater, A. G.: Probabilistic quantitative precipitation estimation in complex terrain, J. Hydrometeorol., 7, 3–22, 2006. a
Clark, M. P., Rupp, D. E., Woods, R. A., Zheng, X., Ibbitt, R. P., Slater, A. G., Schmidt, J., and Uddstrom, M. J.: Hydrological data assimilation with the ensemble Kalman filter: Use of streamflow observations to update states in a distributed hydrological model, Adv. Water Resour., 31, 1309–1324, 2008. a, b, c
Cloke, H. and Pappenberger, F.: Ensemble flood forecasting: A review, J. Hydrol., 375, 613–626, 2009. a
Cronshey, R.: Urban hydrology for small watersheds, Tech. rep., US Dept. of Agriculture, Soil Conservation Service, Engineering Division, https://www.nrcs.usda.gov/Internet/FSE_DOCUMENTS/stelprdb1044171.pdf (last access: 21 February 2022), 1986. a
Dasgupta, A., Hostache, R., Ramsankaran, R., Schumann, G. J.P., Grimaldi, S., Pauwels, V. R., and Walker, J. P.: On the Impacts of Observation Location, Timing, and Frequency on Flood Extent Assimilation Performance, Water Resour. Res., 57, e2020WR028238, https://doi.org/10.1029/2020WR028238, 2021a. a
Dasgupta, A., Hostache, R., Ramsankaran, R., Schumann, G. J.P., Grimaldi, S., Pauwels, V. R., and Walker, J. P.: A mutual informationbased likelihood function for particle filter flood extent assimilation, Water Resour. Res., 57, e2020WR027859, https://doi.org/10.1029/2020WR027859, 2021b. a
Dee, D. P.: Bias and data assimilation, Q. J. Roy. Meteorol. Soc., 131, 3323–3343, 2005. a
Demeritt, D., Cloke, H., Pappenberger, F., Thielen, J., Bartholmes, J., and Ramos, M.H.: Ensemble predictions and perceptions of risk, uncertainty, and error in flood forecasting, Environ. Hazards, 7, 115–127, 2007. a
De Roo, A. P., Gouweleeuw, B., Thielen, J., Bartholmes, J., BongioanniniCerlini, P., Todini, E., Bates, P. D., Horritt, M., Hunter, N., Beven, K., Pappenberger F., Heise, E., Rivin, G., Hils, M., Hollingsworth, A., Holst, B., Kwadijk, J., Reggiani, P., Van Dijk, M., Sattler K., and Sprokkereef, E.: Development of a European flood forecasting system, Int. J. River Basin Manage., 1, 49–59, 2003. a
Desai, B., Maskrey, A., Peduzzi, P., De Bono, A., and Herold, C.: Making development sustainable: the future of disaster risk management, global assessment report on disaster risk reduction, https://www.preventionweb.net/english/hyogo/gar/2015/en/garpdf/GAR2015_EN.pdf (last access: 21 February 2022), 2015. a
Di Baldassarre, G. and Montanari, A.: Uncertainty in river discharge observations: a quantitative analysis, Hydrol. Earth Syst. Sci., 13, 913–921, https://doi.org/10.5194/hess139132009, 2009. a
EMDAT: The OFDA/CRED International Disaster Database, Universite Catholique de Louvain, Brussels, Belgium, https://emdat.be/frontpage (last access: 21 FEbruary 2022), 2016. a
Evensen, G.: The ensemble Kalman filter: Theoretical formulation and practical implementation, Ocean Dynam., 53, 343–367, 2003. a, b, c, d, e, f
FLO2D: TwoDimensional Flood Routing Model, FLO2D [code], https://flo2d.com/, last access: 21 February 2022. a
GarcíaPintado, J., Barberá, G. G., Erena, M., and Castillo, V. M.: Rainfall estimation by rain gaugeradar combination: A concurrent multiplicativeadditive approach, Water Resour. Res., 45, W01415, https://doi.org/10.1029/2008WR007011, 2009. a
GarcíaPintado, J., Neal, J. C., Mason, D. C., Dance, S. L., and Bates, P. D.: Scheduling satellitebased SAR acquisition for sequential assimilation of water level observations into flood modelling, J. Hydrol., 495, 252–266, 2013. a, b, c, d, e, f, g
GarcíaPintado, J., Mason, D. C., Dance, S. L., Cloke, H. L., Neal, J. C., Freer, J., and Bates, P. D.: Satellitesupported flood forecasting in river networks: A real case study, J. Hydrol., 523, 706–724, 2015. a, b, c, d
Gaspari, G. and Cohn, S. E.: Construction of correlation functions in two and three dimensions, Q. J. Roy. Meteorol. Soc., 125, 723–757, 1999. a
Giustarini, L., Matgen, P., Hostache, R., Montanari, M., Plaza, D., Pauwels, V. R. N., De Lannoy, G. J. M., De Keyser, R., Pfister, L., Hoffmann, L., and Savenije, H. H. G.: Assimilating SARderived water level data into a hydraulic model: a case study, Hydrol. Earth Syst. Sci., 15, 2349–2365, https://doi.org/10.5194/hess1523492011, 2011. a, b, c, d
Goovaerts, P.: Geostatistics for natural resources evaluation, Oxford University Press on Demand, ISBN 9780195115383, 1997. a
Grimaldi, S., Petroselli, A., and Nardi, F.: A parsimonious geomorphological unit hydrograph for rainfall–runoff modelling in small ungauged basins, Hydrolog. Sci. J., 57, 73–83, 2012. a, b
Grimaldi, S., Petroselli, A., and Romano, N.: GreenAmpt CurveNumber mixed procedure as an empirical tool for rainfall–runoff modelling in small and ungauged basins, Hydrol. Process., 27, 1253–1264, 2013. a
Grimaldi, S., Li, Y., Pauwels, V. R., and Walker, J. P.: Remote sensingderived water extent and level to constrain hydraulic flood forecasting models: Opportunities and challenges, Surv. Geophys., 37, 977–1034, 2016. a, b, c
Haan, C. T., Barfield, B. J., and Hayes, J. C.: Design hydrology and sedimentology for small catchments, Elsevier, ISBN 9780080571645, 1994. a
Heuvelink, G. B., Brown, J. D., and van Loon, E. E.: A probabilistic framework for representing and simulating uncertain environmental variables, Int. J. Geogr. Inform. Sci., 21, 497–513, 2007. a
Hodgson, M. E. and Bresnahan, P.: Accuracy of airborne LiDARderived elevation, Photogram. Eng. Remote Sens., 70, 331–339, 2004. a
Hopson, T. M. and Webster, P. J.: A 1–10day ensemble forecasting scheme for the major river basins of Bangladesh: Forecasting severe floods of 2003–07, J. Hydrometeorol., 11, 618–641, 2010. a
Horritt, M. and Bates, P.: Predicting floodplain inundation: rasterbased modelling versus the finiteelement approach, Hydrol. Process., 15, 825–842, 2001. a
Hostache, R., Lai, X., Monnier, J., and Puech, C.: Assimilation of spatially distributed water levels into a shallowwater flood model. Part II: Use of a remote sensing image of Mosel River, J. Hydrol., 390, 257–268, 2010. a
Hostache, R., Matgen, P., Montanari, A., Montanari, M., Hoffmann, L., and Pfister, L.: Propagation of uncertainties in coupled hydrometeorological forecasting systems: A stochastic approach for the assessment of the total predictive uncertainty, Atmos. Res., 100, 263–274, 2011. a
Hostache, R., Chini, M., Giustarini, L., Neal, J., Kavetski, D., Wood, M., Corato, G., Pelich, R.M., and Matgen, P.: NearRealTime Assimilation of SARDerived Flood Maps for Improving Flood Forecasts, Water Resour. Res., 54, 5516–5535, 2018. a, b, c, d, e
Houtekamer, P. L. and Mitchell, H. L.: Data assimilation using an ensemble Kalman filter technique, Mon. Weather Rev., 126, 796–811, 1998. a
Hunt, B. R., Kostelich, E. J., and Szunyogh, I.: Efficient data assimilation for spatiotemporal chaos: A local ensemble transform Kalman filter, Physica D, 230, 112–126, 2007. a, b
Jasper, K., Gurtz, J., and Lang, H.: Advanced flood forecasting in Alpine watersheds by coupling meteorological observations and forecasts with a distributed hydrological model, J. Hydrol., 267, 40–52, 2002. a
Jenson, S. K. and Domingue, J. O.: Extracting topographic structure from digital elevation data for geographic information system analysis, Photogram. Eng. Remote Sens., 54, 1593–1600, 1988. a
Knight, D. and Shamseldin, A.: River basin modelling for flood risk mitigation, CRC Press, ISBN 9780415383448, 2005. a
Krzhizhanovskaya, V. V., Shirshov, G., Melnikova, N., Belleman, R. G., Rusadi, F., Broekhuijsen, B., Gouldby, B., Lhomme, J., Balis, B., Bubak, M., Pyayt, A., Mokhov, I., Ozhigin, A., Lang, B., and Meijer, R.: Flood early warning system: design, implementation and computational modules, Proced. Comput. Sci., 4, 106–115, 2011. a
Kundzewicz, Z. W.: 15 Floods: lessons about early warning systems, Late lessons from early warnings: science, precaution, innovation, p. 25, https://www.eea.europa.eu/publications/latelessons2/latelessonschapters/latelessonsiichapter15/view (last access: 21 February 2022), 2013. a
Lai, X., Liang, Q., Yesou, H., and Daillet, S.: Variational assimilation of remotely sensed flood extents using a 2D flood model, Hydrol. Earth Syst. Sci., 18, 4325–4339, https://doi.org/10.5194/hess1843252014, 2014. a
Leon, J. X., Heuvelink, G. B., and Phinn, S. R.: Incorporating DEM uncertainty in coastal inundation mapping, PLoS One, 9, e108727, https://doi.org/0.1371/journal.pone.0108727, 2014. a, b
Li, S., MacMillan, R., Lobb, D. A., McConkey, B. G., Moulin, A., and Fraser, W. R.: Lidar DEM error analyses and topographic depression identification in a hummocky landscape in the prairie region of Canada, Geomorphology, 129, 263–275, 2011. a
Liu, Y. and Gupta, H. V.: Uncertainty in hydrologic modeling: Toward an integrated data assimilation framework, Water Resour. Res., 43, 1–18, 2007. a, b
Liu, Y., Weerts, A. H., Clark, M., Hendricks Franssen, H.J., Kumar, S., Moradkhani, H., Seo, D.J., Schwanenberg, D., Smith, P., van Dijk, A. I. J. M., van Velzen, N., He, M., Lee, H., Noh, S. J., Rakovec, O., and Restrepo, P.: Advancing data assimilation in operational hydrologic forecasting: progresses, challenges, and emerging opportunities, Hydrol. Earth Syst. Sci., 16, 3863–3887, https://doi.org/10.5194/hess1638632012, 2012. a
Livne, E. and Svoray, T.: Components of uncertainty in primary production model: the study of DEM, classification and location error, Int. J. Geogr. Inform. Sci., 25, 473–488, 2011. a
Madsen, H. and Skotner, C.: Adaptive state updating in realtime river flow forecasting – A combined filtering and error forecasting procedure, J. Hydrol., 308, 302–312, 2005. a, b
Marsigli, C., Boccanera, F., Montani, A., and Paccagnella, T.: The COSMOLEPS mesoscale ensemble system: validation of the methodology and verification, Nonlin. Processes Geophys., 12, 527–536, https://doi.org/10.5194/npg125272005, 2005. a
Martinis, S., Kersten, J., and Twele, A.: A fully automated TerraSARX based flood service, ISPRS J. Photogram. Remote Sen., 104, 203–212, 2015. a, b
Mason, D., Schumann, G.P., Neal, J., GarciaPintado, J., and Bates, P.: Automatic near realtime selection of flood water levels from high resolution Synthetic Aperture Radar images for assimilation into hydraulic models: a case study, Remote Sens. Environ., 124, 705–716, 2012. a, b, c
Matgen, P., Schumann, G., Pappenberger, F., and Pfisterz, L.: Sequential assimilation of remotely sensed water stages in flood inundation models, IAHS Publ., 316, 78–88, 2007. a
Matgen, P., Montanari, M., Hostache, R., Pfister, L., Hoffmann, L., Plaza, D., Pauwels, V. R. N., De Lannoy, G. J. M., De Keyser, R., and Savenije, H. H. G.: Towards the sequential assimilation of SARderived water stages into hydraulic models using the Particle Filter: proof of concept, Hydrol. Earth Syst. Sci., 14, 1773–1785, https://doi.org/10.5194/hess1417732010, 2010. a, b, c
Matheron, G.: Universal kriging, in: Matheron's Theory of Regionalised Variables, Oxford University Press, 123–180, ISBN13 9780198835660, https://doi.org/10.1093/oso/9780198835660.001.0001, 1969. a
Mazzoleni, M., Verlaan, M., Alfonso, L., Monego, M., Norbiato, D., Ferri, M., and Solomatine, D. P.: Can assimilation of crowdsourced data in hydrological modelling improve flood prediction?, Hydrol. Earth Syst. Sci., 21, 839–861, https://doi.org/10.5194/hess218392017, 2017. a
McCuen, R. H.: Uncertainty analyses of watershed time parameters, J. Hydrol. Eng., 14, 490–498, 2009. a
McLaughlin, D.: An integrated approach to hydrologic data assimilation: interpolation, smoothing, and filtering, Adv. Water Resour., 25, 1275–1286, 2002. a
McMillan, H. K., Hreinsson, E. Ö., Clark, M. P., Singh, S. K., Zammit, C., and Uddstrom, M. J.: Operational hydrological data assimilation with the recursive ensemble Kalman filter, Hydrol. Earth Syst. Sci., 17, 21–38, https://doi.org/10.5194/hess17212013, 2013. a
Mesa, O. J. and Mifflin, E. R.: On the relative role of hillslope and network geometry in hydrologic response, in: Scale problems in hydrology, Springer, 1–17, ISBN 9789400946781, 1986. a
Montzka, C., Pauwels, V., Franssen, H.J., Han, X., and Vereecken, H.: Multivariate and multiscale data assimilation in terrestrial systems: A review, Sensors, 12, 16291–16333, 2012. a
Moradkhani, H., Hsu, K.L., Gupta, H., and Sorooshian, S.: Uncertainty assessment of hydrologic model states and parameters: Sequential data assimilation using the particle filter, Water Resour. Res., 41, 1–17, 2005. a, b
Mudron, I., Podhoranyi, M., Cirbus, J., Devečka, B., and Bakay, L.: Modelling the Uncertainty of Slope Estimation from a LidarDerived Dem: a Case Study from a LargeScale Area in the Czech Republic/Modelovanie Neistoty Vo Vỳpočte Sklonov Z Lidarovỳch Dmr; Prípadová Štúdia Vybraného Malého Územia V Čr, GeoSci. Eng., 59, 25–39, 2013. a
Murphy, J.: The impact of ensemble forecasts on predictability, Q. J. Roy. Meteorol. Soc., 114, 463–493, 1988. a
Neal, J., Schumann, G., Bates, P., Buytaert, W., Matgen, P., and Pappenberger, F: A data assimilation approach to discharge estimation from space, Hydrol. Process., 23, 3641–3649, 2009. a
Neal, J., Schumann, G., Fewtrell, T., Budimir, M., Bates, P., and Mason, D.: Evaluating a new LISFLOODFP formulation with data from the summer 2007 floods in Tewkesbury, UK, J. Flood Risk Manage., 4, 88–95, 2011. a
Neal, J. C., Atkinson, P. M., and Hutton, C. W.: Flood inundation model updating using an ensemble Kalman filter and spatially distributed measurements, J. Hydrol., 336, 401–415, 2007. a
Nicholas, A. and Mitchell, C.: Numerical simulation of overbank processes in topographically complex floodplain environments, Hydrol. Process., 17, 727–746, 2003. a
NRCS: Ponds – Planning, design, construction, Agriculture Handbook no. 590, Natural Resources Conservation Service Washington, DC, USA, https://nrcspad.sc.egov.usda.gov/distributioncenter/pdf.aspx?productID=115 (last access: 21 February 2022), 1997. a
O'brien, J., Julien, P., and Fullerton, W.: Twodimensional water flood and mudflow simulation, J. Hydraul. Eng., 119, 244–261, 1993. a
Oliver, M. A. and Webster, R.: Kriging: a method of interpolation for geographical information systems, Int. J. Geogr. Inform. Syst., 4, 313–332, 1990. a
Ott, E., Hunt, B. R., Szunyogh, I., Zimin, A. V., Kostelich, E. J., Corazza, M., Kalnay, E., Patil, D., and Yorke, J. A.: A local ensemble Kalman filter for atmospheric data assimilation, Tellus A, 56, 415–428, 2004. a, b, c
Pappenberger, F., Matgen, P., Beven, K. J., Henry, J.B., Pfister, L., and de Fraipont, P.: Influence of uncertain boundary conditions and model structure on flood inundation predictions, Adv. Water Resour., 29, 1430–1449, 2006. a
Pierdicca, N., Chini, M., Pulvirenti, L., Candela, L., Ferrazzoli, P., Guerriero, L., Boni, G., Siccardi, F., and Castelli, F.: Using COSMOSkyMed data for flood mapping: Some casestudies, in: 2009 IEEE International Geoscience and Remote Sensing Symposium, vol. 2, II933, 2009. a, b
Raaflaub, L. D. and Collins, M. J.: The effect of error in gridded digital elevation models on the estimation of topographic parameters, Environ. Model. Softw., 21, 710–732, 2006. a
Rakovec, O., Hazenberg, P., Torfs, P. J. J. F., Weerts, A. H., and Uijlenhoet, R.: Generating spatial precipitation ensembles: impact of temporal correlation structure, Hydrol. Earth Syst. Sci., 16, 3419–3434, https://doi.org/10.5194/hess1634192012, 2012a. a
Rakovec, O., Weerts, A. H., Hazenberg, P., Torfs, P. J. J. F., and Uijlenhoet, R.: State updating of a distributed hydrological model with Ensemble Kalman Filtering: effects of updating frequency and observation network density on forecast accuracy, Hydrol. Earth Syst. Sci., 16, 3435–3449, https://doi.org/10.5194/hess1634352012, 2012b. a
Rasmussen, J., Madsen, H., Jensen, K. H., and Refsgaard, J. C.: Data assimilation in integrated hydrological modeling using ensemble Kalman filtering: evaluating the effect of ensemble size and localization on filter performance, Hydrol. Earth Syst. Sci., 19, 2999–3013, https://doi.org/10.5194/hess1929992015, 2015. a, b
Reichle, R. H., McLaughlin, D. B., and Entekhabi, D.: Hydrologic data assimilation with the ensemble Kalman filter, Mon. Weather Rev., 130, 103–114, 2002. a
RevillaRomero, B., Wanders, N., Burek, P., Salamon, P., and de Roo, A.: Integrating remotely sensed surface water extent into continental scale hydrology, J. Hydrol., 543, 659–670, 2016. a
RodriguezIturbe, I. and Rinaldo, A.: Fractal river networks: chance and selforganization, Cambridge University Press, ISBN 9780521004053, 1997. a
Sakov, P. and Bertino, L.: Relation between two common localisation methods for the EnKF, Comput. Geosci., 15, 225–237, 2011. a
Sakov, P., Evensen, G., and Bertino, L.: Asynchronous data assimilation with the EnKF, Tellus A, 62, 24–29, 2010. a
Scaramuzza, P. and Barsi, J.: Landsat 7 scan line correctoroff gapfilled product development, in: vol. 16, Proceeding of Pecora, 23–27, http://www.asprs.org/a/publications/proceedings/pecora16/Storey_J.pdf (last access: 21 February 2022), 2005. a
Schmidt, A. R.: Analysis of stagedischarge relations for openchannel flows and their associated uncertainties, PhD thesis, University of Illinois at UrbanaChampaign, http://hdl.handle.net/2142/83191 (last access: 21 February 2022), 2002. a
Schumann, G. J.P., Bates, P. D., Neal, J. C., and Andreadis, K. M.: Measuring and mapping flood processes, in: Hydrometeorological hazards, risks and disasters, Elsevier, 35–64, https://doi.org/10.1016/B9780123948465.000023, 2015. a
Sivapalan, M., Takeuchi, K., Franks, S., Gupta, V., Karambiri, H., Lakshmi, V., Liang, X., McDonnell, J., Mendiondo, E., O'connell, P., et al.: IAHS Decade on Predictions in Ungauged Basins (PUB), 2003–2012: Shaping an exciting future for the hydrological sciences, Hydrolog. Sci. J., 48, 857–880, 2003. a
Sorooshian, S., Hsu, K.L., Coppola, E., Tomassetti, B., Verdecchia, M., and Visconti, G.: Hydrological modelling and the water cycle: coupling the atmospheric and hydrological models, in: vol. 63, Springer Science & Business Media, ISBN 9783540778431, 2008. a
Stokstad, E.: Scarcity of rain, stream gages threatens forecasts, Science, 285, 1199–1200, https://doi.org/10.1126/science.285.5431.1199, 1999. a
Tarboton, D. G., Bras, R. L., and RodriguezIturbe, I.: On the extraction of channel networks from digital elevation data, Hydrol. Process., 5, 81–100, 1991. a
Tarquini, S., Isola, I., Favalli, M., and Battistini, A.: TINITALY, a digital elevation model of Italy with a 10 meters cell size (Version 1.0), Istituto Nazionale di Geofisica e Vulcanologia (INGV) [data set], https://doi.org/10.13127/TINITALY/1.0, 2007. a
Tarquini, S., Vinci, S., Favalli, M., Doumaz, F., Fornaciai, A., and Nannipieri, L.: Release of a 10mresolution DEM for the Italian territory: Comparison with globalcoverage DEMs and anaglyphmode exploration via the web, Comput. Geosci., 38, 168–170, 2012. a
Thampapillai, D. J. and Musgrave, W. .: Flood damage mitigation: A review of structural and nonstructural measures and alternative decision frameworks, Water Resour. Res., 21, 411–424, 1985. a
Thiessen, A. H.: Precipitation averages for large areas, Mon. Weather Rev., 39, 1082–1089, 1911. a
Tippett, M. K., Anderson, J. L., Bishop, C. H., Hamill, T. M., and Whitaker, J. S.: Ensemble square root filters, Mon. Weather Rev., 131, 1485–1490, 2003. a
USGS: EarthExplorer, USGS [data set], https://earthexplorer.usgs.gov/, last access: 22 February 2022. a
Weerts, A. H. and El Serafy, G. Y.: Particle filtering and ensemble Kalman filtering for state updating with hydrological conceptual rainfallrunoff models, Water Resour. Res., 42, 1–17, 2006. a
Whitaker, J. S. and Hamill, T. M.: Ensemble data assimilation without perturbed observations, Mon. Weather Rev., 130, 1913–1924, 2002. a, b
Wing, O. E., Quinn, N., Bates, P. D., Neal, J. C., Smith, A. M., Sampson, C. C., Coxon, G., Yamazaki, D., Sutanudjaja, E. H., and Alfieri, L.: Toward Global Stochastic River Flood Modeling, Water Resour. Res., 56, e2020WR027692, https://doi.org/10.1029/2020WR027692, 2020. a
Wood, M., Hostache, R., Neal, J., Wagener, T., Giustarini, L., Chini, M., Corato, G., Matgen, P., and Bates, P.: Calibration of channel depth and friction parameters in the LISFLOODFP hydraulic model using mediumresolution SAR data and identifiability techniques, Hydrol. Earth Syst. Sci., 20, 4983–4997, https://doi.org/10.5194/hess2049832016, 2016. a
Xu, H.: Modification of normalised difference water index (NDWI) to enhance open water features in remotely sensed imagery, Int. J. Remote Sens., 27, 3025–3033, 2006. a, b
Yu, W., Nakakita, E., and Jung, K.: Flood forecast and early warning with highresolution ensemble rainfall from numerical weather prediction model, Proced. Eng., 154, 498–503, 2016. a