Flood hazard of major river segments, Benevento Province, Southern Italy

ABSTRACT On 15 October 2015, a storm-induced flood hit the central sector of Benevento Province (southern Italy) causing two deaths and severe damage to infrastructure, buildings and local agriculture. This area has a long history of similar events and since 1924 its major river segments have been monitored with several hydrometric stations. We used data from two of these stations and a LiDAR derived high-resolution topography to develop a flood hazard map. For map computation, we first derived a flood inundation map from topography. Subsequently we estimated the probability of exceedance of each specific fluvial stage from the combination of a Generalized Extreme Value and a Gamma fits of available hydrometric data. As boundary condition, we considered a reference scenario corresponding to an estimated 500 year flood. The hazard maps provide an overview of the flood hazard in the central sector of Benevento Province and floodplains zonation in flood perspective.


Introduction
In the last three decades, ongoing climate change has been responsible for a global increase in the frequency and magnitude of heavy precipitation (e.g. Groisman et al., 2005). In the Mediterranean regions, such an increase is influencing the recurrence and magnitude of damaging flooding events and consequently the exposure of urban and rural communities (e.g. Cevasco et al., 2015;Diodato, 2007;Knox, 1993). The exposure to flooding is considered as a proxy of the risk, so that an evaluation of the extent of the potential inundation areas is crucial for hazard and risk assessment (e.g. Jongman, Ward, & Aerts, 2012;Peduzzi, Dao, Herold, & Mouton, 2009). This evaluation is often based on the outcome of deterministic hydrodynamic models that simulate the water movement across the floodplain (e.g. Bates & De Roo, 2000;Di Baldassarre, Schumann, Bates, Freer, & Beven, 2010;Gonçalves, Marafuz, & Gomes, 2015;Teng et al., 2017;Toda, Yokingco, Paringit, & Lascoad, 2017) or, in presence of monitoring (fluvial stage or discharge) and topographic data, it can be completed using statistical models associated with GIS processing (e.g. Alfonso, Mukolwem, & Di Baldassarre, 2016). Statistical models have been widely applied for explaining earth surface and hydrologic processes (e.g. Diodato et al., 2014;Diodato et al., 2015;Guerriero et al., 2015) and have the potential to support natural hazard evaluation (e.g. Di Martire, De Rosa, Pesce, Santangelo, & Calcaterra, 2012;Grelle et al., 2014). The combination of statistical models and GIS processing is particularly suitable for flood hazard map production (e.g. flood directive 2007/60/CE; Barredo, 2009).
Evaluating flood hazard requires the production of flood inundation maps showing the water depth during flood events of specific magnitudes and the definition of a reference extreme flood scenario (i.e. design event, e.g. Nuswantoro, Diermanse, & Molkenthin, 2016;Woo & Waylen, 1986). Such scenario represents the boundary condition of the analysis. For regions where long-records of hydrometric data are available, hazard maps can be obtained from some sort of inundation models and probability analysis of available time series (e.g. Woo & Waylen, 1986). If no data exists, the hazard evaluation, as well as the production of a flood inundation map, can be based on historical information provided by locals and geomorphological observations (Furdada, Calderón, & Marqués, 2008;Hosking & Wallis, 1986;Magliulo & Valente, 2014;Montané, Buffin-Bélanger, Vinet, & Vento, 2017;Piacentini, Urbano, Sciarra, Schipani, & Miccadei, 2016). Ideally, a flood inundation map indicating the floodable area in relation to a specific fluvial stage, should be produced using numerical flood plain models developed on the basis of LiDAR data (Montané et al., 2017). This would make the floodplain zonation as accurate as expected for land planning purposes, especially in the presence of urban areas.
Similarly, the reference extreme flood scenario can be identified on the basis of historical data (e.g. Sutcliffe, 1987) but it can be a very challenging task, that depends on the availability, the length, the reliability and the continuity of the records (e.g. Hosking & Wallis, 1986). In the absence of historical data, it can be estimated on the basis of geomorphological field observations (e.g. Furdada et al., 2008;Magliulo & Cusano, 2016;Montané et al., 2017). If the evaluation is based on discontinuous short records and poorly constrained historical information, the extreme scenario could correspond to a well-known event selected because society wants to protect itself from against an event of corresponding magnitude (Garry & Graszk, 1999). Additionally, if this event has a return period shorter than 100 years, an event of magnitude corresponding to that of a 100 year flood should be selected. In Italy, the guidelines of the Ministry of Environment and Land Protection indicate as reference extreme event, or in other words the lower boundary of the hazard zonation (lower boundary of the low hazard zone), an event of magnitude corresponding to that of a 300 to 500 year flood.
On 15 October 2015, a destructive overflow of the Tammaro and Calore rivers hit the town of Benevento in southern Italy, and the central sector of its province, already known for other natural hazards like landslides and earthquakes (e.g. Guerriero et al., 2016;Guerriero et al., 2017;Maresca, Castellano, De Matteis, Saccorotti, & Vaccariello, 2003;Revellino, Grelle, Donnarumma, & Guadagno, 2010). The event caused two deaths and severe damage to infrastructure, buildings and local agriculture. Most of the protection embankments, built along the Calore river in the twentieth century, were surmounted and/or damaged inundating otherwise protected areas. It was caused by a storm characterized by a maximum intensity of 27.4 mm/10 min and a maximum cumulative amount of rainfall in 19 h of 415.6 mm (registered in Paupisi, Benevento Province). As reported in Table 1, this event has a number of historical precedents in this area (data from Andriola, Delmonaco, Margottini, Serafini, & Trocciola, 1996;Cardinali et al., 1998;Diodato, 1999;Rossi & Villani, 1994;Zazo, 1949). Among these events, the destructive flood of October 1949 is remembered for its impact on the society and territory. It was induced by a storm and was responsible for 47 fatalities and severe damage to buildings, infrastructure and service lines. During this event, the hydrometric station of the 'Ponte di Annibale' located along the lower course of the Volturno river indicated a fluvial stage of about 10 meters that corresponds to an estimated discharge of 3200 m 3 /s (for further details see the AVI project, http://avi.gndci.cnr.it/). Despite its impact no data about its magnitude (i.e. water depth) is available for the town of Benevento and the surrounding area.
On this basis, and considering the need for an updated tool to guide future land planning in the Benevento Province, we use LiDAR data acquired in 2012 across the Benevento Province and discontinuous hydrometric data from 1924 to 2016 to produce a flood hazard map of its two major river segments. The first is the eastern, upper, segment of the Calore river starting approximately at the confluence with the Tammaro river and covering all of the industrial The second is the western, lower, segment of the Calore river corresponding to the Telesina Valley, known for it's high value agricultural production (Figure 1(b)). During the flood, the hydrometric stations located within the eastern and the western segments were destroyed, therefore, we made a comprehensive analysis of the flood in order to estimate its magnitude (i.e. fluvial stage) and associated return period.

Methods
Below, we describe the methods used to study the flood event of October 2015 and consequently develop the flood hazard map using high spatial resolution topography and probability analysis of hydrometric time series.

October 2015 event characterization
We used a post event Top of the Atmosphere (TOA) Sentinel 2 multispectral image taken on 15 November 2015 to map the area inundated by the flood. To recognize this area from the multispectral image, we used false color image composition (i.e. 11-8-3 bands, SWIR-NIR-Green, Figures 2(a,b)) and, since most of the floodplain is not covered by high vegetation, we computed the Normalized Difference Water Index (NDWI; e.g. McFeeters, 1996;Figures 2(c,d)).
The NDWI is a satellite-derived index from the reflectance in green and Near-Infrared (NIR) channels that can be applied to delineate open water features and has the potential to distinguish flooded and nonflooded areas (Hudson & Colditz, 2003;Wang, Colby, & Mulcahy, 2002). Since it might not work properly under high vegetation that typically characterize fluvial shores, it is often used in combination with the NDVI index (e.g. Powell, Jakemam, & Croke, 2014). In this case, zones characterized by high vegetation are limited because most of the floodplain is used for cropping and due to the effect of the flood that removed most of the trees that colonize river embankments. Therefore, in this case the NDVI was not used.
Once recognized through visual interpretation of the images, the inundated area was field validated considering the position of the geomorphologic signature of the flood on inactive fluvial scarps bounding river terraces and considering the photographic documentation acquired immediately after the event (e.g. Figures 3(a,b)). After validation, the mapped flooded area was used as boundary condition (target water level) to estimate water depth during the event next to the monitoring stations of Benevento and Solopaca. This estimation was completed through the reconstruction of the cross-sectional geometry of the floodplain and intersection of a line at the target water level and using lake analysis. Lake analysis was completed through the r.lake module implemented in GRASS GIS (e.g. Cannata & Marzocchi, 2012) over an area surrounding the monitoring stations. The r.lake module fills a lake to a target water level from a given starting point using as input data a raster map with at least one cell value greater than zero (e.g. https://grass. osgeo.org/grass72/manuals/r.lake.html). In this case, the input data for lake analysis was the LiDAR derived flood inundation map (See below for flood inundation map generation). The maximum relative elevation (i.e. water depth) above the river, estimated from cross sections (Figures 4(a,b)) and lakes (Figures 4(c,d)), was used as reference value for this event in the annual maxima hydrometric time series. We used probability analysis to estimate the theoretical return period of this event on the basis of available data and understand if it can be considered as the reference flood scenario for the Benevento Province according to the Garry and Graszk (1999) method. The theoretical return period was calculated as the inverse of the probability of exceedance.

Flood inundation and hazard mapping
As a first step in deriving the flood hazard map, we compute a flood inundation map (i.e. relative elevation above the river) as a difference between a bare-earth, airborne LiDAR derived Digital Surface Model and a digital elevation model of the water surface (wDEM). The LiDAR DSM has a 2 × 2 m cell size and is resampled from a 1 × 1 m cell size using cubic convolution (e.g. Wu, Li, & Huang, 2008). The wDEM is computed through the interpolation of a triangulated irregular network of the water elevation data at cross-sections constructed along the river segments. Cross-section positions are consistently selected as function of the As a second step, we use a type III Generalized Extreme Value distribution function (GEV, ξ < 0; e.g. Hosking & Wallis, 1986) where ξ is the shape parameter, σ is the scale parameter and µ is the location parameter. As indicated by Coles (2001), the Generalized Extreme Value distribution function is suitable to describe the statistical behavior of annual maxima but has the disadvantage of underestimating the intensity of very high return period events.
To overcome this drawback, we use it in combination with a Gamma distribution function (e.g. Yue, 2001) that has the potential to better estimate the intensity of very high return period events and has the following form: where α is the shape parameter and β is the scale parameter and Γ(α) is the gamma function that is calculated as follow: GEV and Gamma parameterization is automatically completed using the distribution fitter tool implemented in MATLAB (GEV Benevento: ξ = −0.168, σ = 1.740, µ = 3.415; GEV Solopaca: ξ = −0.153, σ = 1.765, µ = 3.752; Gamma Benevento α = 2.9812, β = 1.39506; Gamma Solopaca: α = 4.61855, β = 0.979929). Function combination is completed considering the univariate statistics of each time series. Especially, we use the GEV functions for inundation stages lower or equal to the time series 75 th percentile (5.23 m for the Benevento hydrometric station and 5.53 m for the Solopaca hydrometric station; Figures (b and e) of the Main Map) and the Gamma function for higher inundation stages. Since we have two time-series, the statistical analysis provides two specific GEV distribution functions and two specific Gamma distribution functions, a first pair for the eastern segment crossing the Benevento town (Benevento station) and a second pair for the western segment comprising most of the Telesina Valley (Solopaca station). Each pair is used to derive a composite fit for each station that is shown in the graphs c and f of the Main Map using the dashed red line.
In both cases, residual analysis indicates a very high goodness of fit for both functions (R 2 ≥ 0.99 for the GEV functions and R 2 ≥ 0.95 for the Gamma functions, see Figs. c and f in the Main Map). Even if the goodness of fit of the Gamma function is lower that the GEV function, the visual analysis of the graphs of figures c and f of the Main Map suggests the higher performances of the Gamma function at higher water stages. Additionally, a calculation of the RMSE between the measured and modeled fluvial stages (for the two stations and as function of the theoretical probability of exceedance) indicates that for higher stages, the RMSE of the GEV function is consistently higher than the RMSE of the Gamma function (e.g. considering the highest three stages: Benevento, GEV RMSE = 0.83, Gamma RMSE = 0.77, Solopaca, GEV RMSE = 1.04, Gamma RMSE = 0.69; considering the highest two stages: Benevento, GEV RMSE = 1.0, Gamma RMSE = 0.73, Solopaca, GEV RMSE = 1.26, Gamma RMSE = 0.65). The combined functions (GEV for x ≤ 75th percentile and Gamma for x > 75th percentile) are used to derive the flood hazard maps from inundation maps as annual probability of exceedance maps (Main Map) and floodplain zonation maps (inset maps of the Main Map). Floodplain zonation is completed considering four zones, each corresponding to a specific hazard level. These zones are selected accordingly to the local Basin Authority guideline and considering the morphology of the floodplain. Especially, we identified a i) very high hazard zone that could be flooded as a consequence of 1 to 5 year floods, a ii) high hazard zone that could be flooded as a consequence of 5 to 30 year floods, a iii) medium hazard zone that could be flooded as a consequence of 30 to 100 year floods and a iv) low hazard zone that could be flooded as a consequence of 100 to 500 year floods. Flood magnitude for each specific return period is derived from probability analysis considering that the theoretical return period is the inverse of the probability of exceedance.
It is important to note that flood hazard estimation is completed using a simple elevation-based approach (e.g. Teng et al., 2017) that has the disadvantage of classifying, areas topographically disconnected from the river course and/or the flood plain, as floodable. This drawback would suggest that a preliminary identification of the area connected or connectable (e.g. area protected by river embankments) to the river is needed. Once identified, these areas can be manually removed from the model. Having said that, we do not exclude these areas from the analysis because even if they might not be directly flooded by the river outflow, they have the potential to be flooded due to secondary pluvial process like concentrated runoff or could be flooded in the case of embankment damage as in the event of October 2015. Thus, for the study area this represents a cautionary condition that might be considered as acceptable. Additionally, the flood inundation map is developed considering as an upper boundary the maximum stage of the river estimated for a 500 year flood. This is according to the Italian Basin Authorities guidelines.
Moreover, the flood inundation maps on which the hazard computation is based are not fully comparable (from the spatial viewpoint and considering its total extent) with any event-map. However, they are able to simulate or be reliable for different flooding scenarios of same magnitude but generated by different rainfall spatial distribution. This because the presented method is simpler than hydrodynamic models and only based on the relation between floodplain topography and fluvial stage. This makes the method able to work at both local and regional scale. A further simplification is that we do not consider the effect of climate change and channel adjustment that might influence the reliability of the analysis in the longterm (Magliulo, Valente, & Cartojan, 2013;Scorpio & Rosskopf, 2016;Whitfield, 2012).

The October 15th event and the extreme flood scenario
The flood was triggered by a storm that started on the night of 14 October and stopped on 15 October 2015, early in the morning. Data from the monitoring network of the Campania Region (http://www. agricoltura.regione.campania.it/meteo/rete.htm) indicated that a maximum cumulative rainfall value of ∼415 mm was recorded at the Paupisi meteorological station. More than 100 mm was recorded at 20 meteorological stations located within the Benevento Province in an area of about 2000 km². Approximately 50% of this area corresponds to the Tammaro river basin (∼800 km²). The destructive overflow hit the middle and lower course of the Tammaro river downstream from the Campolattaro dam and the middle and lower course of the Calore river that crosses the industrial area, the town of Benevento and the Telesina Valley ( Figure 2). The flood involved an area across rivers floodplains of about 30 km 2 and was characterized by a maximum river stage of ∼9.8 m estimated next to the Solopaca station (Figure 4(b,d)). A water stage of ∼9.5 m was estimated within the town of Benevento (Figure 4(a,c)). In the industrial area of Benevento, next to the confluence between the Tammaro and the Calore rivers, all of the factories were damaged or completely destroyed. It was estimated that the flood caused about €1 billion worth of damage to local agriculture (Magliulo & Cusano, 2016).
Despite its socioeconomic impact, the probability analysis of available data indicates that this event has a theoretical return period of between 25 and 50 years (Solopaca and Benevento hydrometric stations, respectively; see gray stars in graphs c and f of the Main Map). Flooding events with return periods shorter than 100 years should not be considered as reference for hazard analysis even if they represent a well-known event that society wants to protect itself from (Garry & Graszk, 1999). Additionally, the guideline of the Italian Ministry of Environment and Land Protection underline that for hazard evaluation a 300 to 500 year event should be considered as extreme flood reference. For the analysis, we consider the 500 year flood and estimated its magnitude in terms of water depth using the combined function that for high intensity event uses the Gamma distribution function. It indicates that the 500 year event has an estimated magnitude of ∼13 m for the Solopaca hydrometric station (i.e. the western segment) and ∼14.5 m for the Benevento hydrometric station (i.e. the eastern segment; red stars in graphs c and f of the Main Map).

The flood hazard map
The Main Map reports the hazard maps derived from GIS processing and probability analysis. It is composed by the annual probability of exceedance maps and the hazard zonation maps (inset maps of the Main Map). The annual probability of exceedance maps for the major river segments of the Benevento Province are an interpretation of the hazard level across the Calore, Tammaro and Ufita floodplains connected to the occurrence of storm-induced floods of different magnitudes. These maps are simple tools suitable for identifying the hazard level of each human structure located within the floodplains on annual basis and is the first step toward a flood risk evaluation for this area. Being estimated on the basis of a flood inundation map (see methods section), the flood hazard tends to assume higher values for low values of the relative elevation above the river. It also varies as function of the distance from the river and the cross-sectional geometry of the river valley. In fact, where the river banks are sufficiently high, the hazard level evolves form very high (around 100% of probability), within or next to the river bed, to low values, across the adjoining floodplain. An example is the upper Telesina Valley where the river banks has a relative elevation up to 7.5 m.
The zonation maps (inset maps of the Main Map) show four zones of the floodplain that can be flooded by events of specific return period (see methods section). As indicated by the maps, most of the analysed floodplain is classified from very high to high hazard zones and the zones of medium and low hazard are poorly developed. It is important to notice that a part of the industrial zone, the western part of the Benevento town and the area located along the Sabato river are classified as very high to high hazard zones. In this context, the maps have the potential to support future land planning, guiding recommendation of uses, occupation restrictions and flood insurance plan development of the analysed floodplain segments (e.g. Linnenluecke & Griffiths, 2010). Additionally, since most of the existing protection embankments were surmounted and/or damaged by the October 2015 event, they should support the designing of new flood mitigation measures in terms of typology, position and geometry and the identification of zones of flood expansion for lamination purposes.

Conclusions
The hazard maps provide an overview of the flood hazard in the central sector of the Benevento Province which was hit by the flood event of October 2015. Representing the first step toward an evaluation of flood risk in the central sector of Benevento Province, they allow to estimate the specific hazard level for each human structure located across the floodplain of Calore (lower sector), Tammaro (a small part of its lower sector) and Ufita (a limited part of its lower sector) rivers and provide a zonation in flood perspective. Conversely from products that are based on hydrodynamic models of flood, these maps derive from statistical analysis of measured fluvial stage and LiDAR derived high-resolution topography which make them suitable to guide future land planning in an urban area like that of Benevento. In this way, the presented method has the potential to be applied in many fluvial contexts for which fluvial stage and topographic data are available and has no limitation to work at different scales, from local to regional. The resulting hazard maps, developed in direct collaboration with the Benevento Province authority, should be considered as a tool to support policy decision oriented to flood mitigation in the analysed floodplain segments (i.e. recommended uses, occupation restrictions and flood insurance plan development). Being a tool for preventing and reducing losses caused by floods, they also represent a contribution toward a better understanding of the environmental risk in this sector of the Benevento Province.

Software
The map was made using Golden Software Map Viewer 8. Hydrometric data were processed using Microsoft excel and Matlab tm 2017 (https://it.mathworks.com/ products/matlab.html). DEMs were generated and or modified using Quantum GIS, Open Source Geographic Information System, licensed under the GNU General Public License. Quantum GIS is an official project of the Open Source Geospatial Foundation (OSGeo, public domain software).

Funding
This work was supported by the Benevento Province Authority.