New Insights on Flood Mapping Procedure: Two Case Studies in Poland

: The use of the Mike11 one-dimensional (1D) hydraulic model, together with o ﬃ cial hydrology, represents a standard approach of the National Water Management Authority (NWMA) in Poland for ﬂood mapping procedures. A di ﬀ erent approach, based on the hydrological Event-Based Approach for Small and Ungauged Basins (EBA4SUB) model and the Flood-2 Dimensional (FLO-2D) hydraulic model has here been investigated as an alternative procedure. For the analysis, two mountainous rivers in Poland were selected: Kamienica Nawojowska is characterized by a narrow valley, while Skawinka has a broad valley. It was found that the ﬂood zones can enormously di ﬀ er locally, with larger zones generated by the Mike11 / NWMA model in some cases and by the EBA4SUB / FLO-2D model in other situations. The beneﬁts of using the two-dimensional (2D) model are consistent in areas without drainage and where the connection to the main channel is insu ﬃ cient. The use of 1D modeling is preferred for the possibility of mapping the entire river network in a short computational time.


Introduction
Floods are a natural phenomenon that occurred in the past and that will be repeated in the future. For a natural ecosystem, periodic flooding from the river plays a positive role. In situations of intensive urban development, when an increase in urbanization of areas located in the proximity of rivers is observed, floods can cause considerable damages. Adequate flood protection can only be implemented by combining technical measures with proper urban development planning [1,2]. According to the provisions of the Floods Directive [3], flood risk management plans, and their derived flood-prone areas, are crucial elements of flood hazard reduction. In Poland, actions have been carried out for many years in order to delineate such flood-prone areas. The National Security System (ISOK) has the necessary task to provide an effective system for protecting the country from extraordinary hazards, which is particularly essential due to the growing number of such events and the increasing scale of their economic and social impact. Following this task, flood hazard management plans were developed for the majority of catchments in Poland, mainly the ones having discharge gauging stations, allowing local authorities to make decisions about the identification of flood hazard zones.
When determining the flood hazard zones, hydraulic models are used to forecast the water flow surface during the flood wave [4]. When using one-dimensional (1D) models to determine the thanks to the Curve Number for Green-Ampt (CN4GA) internal routine and its automatic calibration, we benefit from the accuracy of a physically based infiltration scheme (the Green-Ampt equation) mixed with the simplicity of an empirical approach (the CN method). Second, for excess rainfall-direct runoff transformation, we determine the instantaneous unit hydrograph (IUH) shape using detailed geomorphological information pixel by pixel, avoiding the use of synthetic methods [20].
Bearing in mind the limitations associated with the hydrological rainfall-runoff models used in Poland for determining flood zones, this study aimed to analyze the possibilities of using the EBA4SUB model to determine the extent of flood zones. The EBA4SUB model was here coupled with FLO-2D, and the results have been compared with the official procedure used in Poland by the National Water Management Authority (i.e., the use of Mike11/NWMA). The novelty of this work is the use of a combined approach based on the hydrological EBA4SUB model with a 2D hydraulic model, FLO-2D, to generate flood zones in two mountainous catchments in Poland. It should be emphasized that there has been no analysis so far regarding the possibility of using this approach to determine the design hydrographs in this region. The study was carried out for selected catchments of the Upper Vistula river basin. Therefore, the conducted study will allow us to determine whether this model may be an alternative to commonly used methodological approaches. Moreover, comparing the official procedure with the proposed approach based on the combination of EBA4SUB and FLO-2D will allow us to understand if an update of the existing flood-prone areas maps is possible, minimizing the user subjectivity.

Study Area
The study regarded two catchments of the Upper Vistula basin: Skawinka and Kamienica Nawojowska. The Upper Vistula basin accounts for about 25% of the entire Vistula basin and about 15% of Poland's area and covers the south-eastern region of the country. It covers part of the Carpathians, Subcarpathian valleys, and Małopolska Uplands. The Upper Vistula basin area shows large variations in elevation [21]. Physiographic parameters of the studied catchments (Skawinka and Kamienica Nawojowska, respectively) are: catchment area (A) 316.0 and 238.0 km 2 , main river length (L) 34.0 and 33.1 km, average main river slope (I) 10.3 and 17.3% , average catchment slope (Ψ) 18.6 and 31.0% , soil imperviousness index (N) 34.0 and 33.0%, and runoff coefficient (Φ) 0.62 and 0.82. The Skawinka catchment area includes agricultural areas (75.6%), forests (20.8%), and artificial surfaces (3.6%). The Kamienica Nawojowska catchment is characterized by forests and seminatural areas (58.6%), agricultural areas (36.4%), and artificial surfaces (5.0%). The studied catchments have a different hydrographic system and the analyzed zones include estuary river sections with a length of about 9 km ( Figure 1). The layout of the river network indicates differences between Kamienica Nawojowska and Skawinka. The tributaries of Kamienica are concentrated in the central part of the basin and the supply is mainly on the left bank. Skawinka has tributaries distributed more evenly along the entire length.

Materials
The physiographic characteristics for the analyzed catchments were determined on the basis of the Hydrographic Division Map for Poland, the DEM with a resolution of 100 m, the Corine Land Cover database, and the Soil Map for Poland. The source of DEM, aerial photos, and images is the Polish Central Office of Geodesy and Cartography [GUGiK], which provides open access services. In

Materials
The physiographic characteristics for the analyzed catchments were determined on the basis of the Hydrographic Division Map for Poland, the DEM with a resolution of 100 m, the Corine Land Cover database, and the Soil Map for Poland. The source of DEM, aerial photos, and images is the Polish Central Office of Geodesy and Cartography [GUGiK], which provides open access services. In particular, the open access services used in the present study are: http://www.gugik.gov.pl/ pzgik/zamow-dane/ortofotomapa, http://www.gugik.gov.pl/pzgik/dane-bez-oplat/dane-dotyczacenumerycznego-modelu-terenu-o-interwale-siatki-co-najmniej-100-m-nmt_100, and http://www.gugik. gov.pl/pzgik/zamow-dane/numeryczny-model-terenu. The source of river channel data was the Mike11 model, which was developed in 2015 in connection with the implementation of the flood hazard maps and flood risk maps. The Polish National Water Management Authority implemented the project. The daily rainfall and flow data for the calculation of the 100-years return period of annual maximum daily precipitation and annual maximum flow were obtained from the Institute of Meteorology and Water Management in Warsaw. The obtained data covered the period from 1981 to 2014. Starting from the available data, hydrologic calculations were first performed employing the EBA4SUB rainfall-runoff model. Then, hydraulic calculations were performed employing one 1D model (Mike11) and one 2D model (FLO-2D), as explained in the following sections. Flood-prone areas were determined using both the hydraulic models assuming a return period of 100 years.

The Proposed Approach
The data needed to apply the EBA4SUB hydrologic model include the selection of a gross rainfall design hyetograph, the DEM of the area, and data on the catchment land use and soil properties [20]. The gross rainfall cumulative value was here assessed on the based on the maximum daily precipitation for a 100-year return period, calculated by log-normal distribution on observed data. It was here assumed that rainfall causes the highest runoff from catchment with a duration equal to concentration time (Tc), so the calculated maximum daily precipitation was rescaled to a rainfall duration equal to Tc. In this approach we assume that the critical rainfall duration is equal to the catchment concentration time, following the theoretic assumption that this choice causes the maximum peak discharge at the outlet compared to shorter or longer rainfall durations [22]. It is noteworthy, anyway, that this hypothesis is debated in the literature. In many practical applications, rainfall durations 2-3 times larger than the concentration time are used in order to maximize the peak discharge [23].
The gross rainfall design hyetographs in this study were determined distributing, within Tc, the gross rainfall cumulative value using the Chicago, the beta distribution, and the Deutscher Verband für Wasserwirtschaft und Kulturbau (DVWK) methods, respectively. For beta distribution, the following parameters were assumed: 5 and 2; 2 and 2; 2 and 3. The choice of these parameters resulted from the verification of the impact of changing the shape of the hyetograph of precipitation on the peak culmination and, in consequence, on the extent of the flood zones. The excess rainfall hyetograph was determined according to the Curve Number for Green-Ampt (CN4GA) procedure proposed by Grimaldi et al. [24]. CN4GA is based on the Curve Number (CN) method and the Green-Ampt (GA) equation. The CN method first estimates the excess rainfall cumulative value with the following formula: where: P n -excess rainfall cumulative value [mm], P-gross rainfall cumulative value [mm], S-maximum potential catchment retention [mm] that is determined based on the CN value of the investigated area. The Green-Ampt equation is then used to determine the excess rainfall hyetograph [25]: for t > t pon (2) where: q 0(t) -infiltration rate, t pon -ponding time, I(t)-cumulative infiltration, K s -saturated hydraulic conductivity, ∆θ-change in soil-water content between the initial value and the field saturated soil-water content, ∆H-the difference between the pressure head at the soil surface and the matrix pressure head at the moving wetting front. The CN4GA automatically assures that the calculated excess rainfall event has the same cumulative value and the same initial abstraction value derived with the CN method. However, it presents a physically based time distribution.
The CN parameter is the only parameter governing the relationship among the gross rainfall cumulative value and the excess rainfall cumulative value. Its value is determined based on soil type, land cover, and land management condition, according to the official tables provided in the the United States Department of Agriculture (USDA) [26]. The The Antecedent Moisture Condition (AMC) II (average condition for soil moisture) was here assumed. The soil types were established based on the Soil Map for Poland. The land cover and land management were established from the Corine Land Cover base.
Runoff hydrograph is then determined using the DEM that was adequately preprocessed, using the width function based instantaneous unit hydrograph (WFIUH). In particular, the DEM preprocessing analysis was performed as follows. First, pits and flat areas were removed. Then, flow paths were defined using an optimized flow direction algorithm [27]. Third, the river network was extracted using the drop analysis and the WFIUH was calculated using the following equation [28]: where: L c , L h -length of the path for the channel and hill slope cell of the DEM, V c , V h -surface flow velocity for the channel and hillslope cell.
Regarding V c and V h , the hillslope surface flow velocity is linked according to empirical formulas to the local slope and land cover data. In contrast to the channel surface flow, velocity is calibrated so that the projection on the time axis of the WFIUH centre of mass is equal to the basin lag time(T L ). The lag time is expressed as the 60% of T c , that is calculated using the Giandotti formula [29]. After having defined the WFIUH, the design hydrograph Q(t) can be calculated with the following equation: In order to characterize the design hydrograph shape, the values of wave slenderness coefficients were determined from the relationship [30]: Sustainability 2020, 12, 8454 6 of 17 The quality of the simulations obtained using the EBA4SUB model was assessed with the mean absolute percentage error (MAPE), described by the relationship [30]: where: Q m , max -maximum flow with a specified return period, calculated using EBA4SUB [m 3 ·s −1 ], Q s , max -maximum flow with the same return period, calculated using the log-normal distribution [m 3 ·s −1 ]. The 2D FLO-2D model was used to determine flood-prone areas. FLO-2D employs the dynamic wave momentum equation solved on a numerical grid of square cells, whose resolution depends on the adopted hydrograph peak discharge [31]. In the present study, aerial photos and images allowed the estimation of hydraulic parameters such as floodplain roughness coefficients. In particular, Manning's roughness coefficients of 0.04 and 0.045 were assigned for the Kamienica Nawojowska main channel and the Skawinka main channel, respectively, with the value of 0.08 assigned for the floodplain areas of both case studies. The channel cross sections were approximated, starting from the existing Mike11 data, to a trapezoidal shape (45 • side slope) with a maximum width of 45 m and a maximum depth of 2 m for Kamienica Nawojowska, and to a triangular shape with a maximum width of 20 m and a maximum depth of 4 m for Skawinka. The selected hydraulic domains are represented by the last 9 km of the Skawinka main channel (the total hydraulic modeling area is equal to 27 km 2 ) and by the last 9 km for the Kamienica Nawojowska main channel (the total hydraulic modeling area is equal to 18 km 2 ). The duration of computer simulations can take a few hours for each one. In the FLO-2D model, the flood-prone areas are expressed in a gridded way, and they are the result of the maximum flow depth occurring in the specific pixel in the whole simulation. In both investigated areas, the key parameters affecting the extension of the zones and their shape are the design hydrograph peak discharge and its total volume. The spatial resolution for the EBA4SUB data was 100 m, for FLO-2D input data 8m, for FLO-2D output data 50 m.

The Standard Approach
The process of generating flood hazard maps is usually carried out in Poland using several stages (Mike11/NWMA procedure). (1) Identification of input data regarding the catchment, river system, floodplains, and hydraulics structures.  [32].
Transformation of the excess rainfall into direct runoff is performed using the Soil Conservation Service Unit Hydrograph (SCS-UH). The excess rainfall hyetograph is determined according the Curve Number (SCS-CN) method described above. The shape of rainfall hyetograph is assumed based on beta distribution [21]. In regards to the SCS-UH method, in order to determine the Synthetic Unit Hydrograph (SUH) shape from the nondimensional q/q p versus the t/tp hydrograph, the peak discharge q p and the time to peak t p are computed, as: where: A-the catchment area, t r -the duration of rainfall T l -the lag time from centroid of rainfall to peak discharge. The T l can be calculated from watershed characteristics using main stream length L, watershed slope s, and curve number (CN): 1900·s 0.5 (9) In this approach, the 1D Mike11 model was used for calculations of unsteady flow in watercourses based on continuity equation and the Saint-Venant's equation of momentum conservation [33]: where: Q-discharge, A-cross-section of the bed, q-lateral tributary, h-ordinate of the water surface, x-longitudinal coordinate measured alongside the river bed, R-hydraulic radius, α c -Coriolis coefficient, t-time, g-gravitational acceleration.
The last equation considers incompressibility and homogeneity of water, small bottom declines, flux parallel to the bottom, and steady motion. For a rapid flow, a reduced equation is used: The exactness of water surface reproduction in the Mike11 software results from the exactness of iterative discharge calculation (10 −4 ) and the area (10 −3 ) and subsequent recalculation for the water surface elevation [34]. The duration of computer simulations last a few seconds (from 4 to 17). The Mike11 model was used for computation of water surface levels at cross-sections in main channel and flood zones. River models accounted for a river of the length of significance in terms of flood protection, i.e., for Kamienica Nawojowska 27.6 km and for Skawinka 35.0 km. Boundary conditions covered all required tributaries concentrated and distributed as provided for in the guidelines [35]. In Mike11 flood-prone areas are determined from the intersection of the numerical water surface model (WSM) with the DEM. WSM was generated using ordinates in cross-sections. Moreover, areas without a hydraulic connection in the riverbed were omitted, as were areas where the water depth is less than the accuracy of DEM. Regarding the input data, for Mike 11 flood zones generation procedure the resolution is less than 1 m.
For both the proposed approach and the standard approach, the source of DEM, aerial photos and images was the Polish Central Office of Geodesy and Cartography [36], which provides open access services. The source of river channel data was a Mike11 model, which was developed in 2015 in connection with the implementation of the flood hazard maps and flood risk maps. The Polish National Water Management Authority implemented the project. Regarding the land cover data and FLO-2D, references are [37,38]. In both the approaches, the final modelled grids were resampled at 50 m × 50 m resolution to provide a direct comparison.
Regarding the calibration of both the approaches, in the proposed approach EBA4SUB is characterized by two advantages. First, the excess rainfall calculation is automatically performed matching the cumulate excess rainfall values computed by applying Equation (2) and Equation (1). It is noteworthy that this approach combines the accuracy of a physically based infiltration scheme (the Green-Ampt equation) with the simplicity of an empirical approach (the SCS-CN method) employing only one parameter (the CN). Second, the WFIUH is calculated based on the basin concentration time that is automatically estimated by the model. Regarding the hydraulic modeling, the proposed approach calibration was not performed (note: for calibration we mean the fitting of parameters of a model to achieve results of simulation more similar to observation). In case of the standard approach, the calibration process was performed during project "Flood protection programme in the Upper Vistula (2014-2015)". The results of calibration of Mike11/NWMA procedure follow the National Water Management Authority rules, where the used metrics are the correlation coefficient (that were higher than 0.70), special correlation coefficient (≥0.70), total square error (≤10%) and errors of culmination level (≤0.15 m), culmination flow (≤10%), culmination dislocation (≤1 1 2 h), and flood wave volume (≤10%). In case of calibration, the model error must qualify it to the range of excellent-good, and for verification to the range of excellent-quite good.

Initial Analysis
The initial analysis of hydrological data was conducted for annual maximum flows time series for Kamienica Nawojowska and Skawinka. The analysis included the determination of some descriptive statistics: positional measures: min, average, and max; measures of dispersion: standard deviation (S) and coefficient of variation (C s ), and measures of shape of the studied variate distribution: coefficient of skewness (Ske) and kurtosis (K). The results are presented in Table 1. Based on the results summarized in Table 1, it was found that the differences between minimum and maximum observed annual maximum flow was 94% for Kamienica Nawojowska and 96% for Skawinka. The dynamics of changes in maximum annual flows remained at a high level, which was evidenced by the coefficient of variation. The coefficients of variations for analyzed catchments are strongly related to the high variability of maximum precipitation, which has an enormous impact on flooding. Flood size in this region depends on multiple factors such as geology, soils, geomorphological evolution, landscape structure, or land use. These factors turn precipitation transformation into the runoff. Analyzing the values of the skewness coefficients, for Kamienica Nawojowska and Skawinka, it was shown that it is greater than 0. Therefore, right-sided asymmetry of the empirical distributions of random variables was found. This is because in the analyzed time-series most observations are smaller than their average value. Therefore, the average values are greater than the medians of time series. In turn, the kurtosis values indicated the platykurtic (negative value) and leptokurtic (positive value) empirical distributions for Kamienica Nawojowska and Skawinka respectively.

EBA4SUB Design Hydrographs
The EBA4SUB model was used to determine the design hydrographs to be propagated at the beginning of the investigated river sections. Hydrographs were determined for the following gross rainfall distributions: beta, DVWK, and Chicago. The characteristics of the modelled design hydrographs are presented in Table 2. The hydrographs shapes are shown in Figure 2a dislocation (≤1½ h), and flood wave volume (≤10%). In case of calibration, the model error must qualify it to the range of excellent-good, and for verification to the range of excellent-quite good.

Initial Analysis
The initial analysis of hydrological data was conducted for annual maximum flows time series for Kamienica Nawojowska and Skawinka. The analysis included the determination of some descriptive statistics: positional measures: min, average, and max; measures of dispersion: standard deviation (S) and coefficient of variation (Cs), and measures of shape of the studied variate distribution: coefficient of skewness (Ske) and kurtosis (K). The results are presented in Table 1. Based on the results summarized in Table 1, it was found that the differences between minimum and maximum observed annual maximum flow was 94% for Kamienica Nawojowska and 96% for Skawinka. The dynamics of changes in maximum annual flows remained at a high level, which was evidenced by the coefficient of variation. The coefficients of variations for analyzed catchments are strongly related to the high variability of maximum precipitation, which has an enormous impact on flooding. Flood size in this region depends on multiple factors such as geology, soils, geomorphological evolution, landscape structure, or land use. These factors turn precipitation transformation into the runoff. Analyzing the values of the skewness coefficients, for Kamienica Nawojowska and Skawinka, it was shown that it is greater than 0. Therefore, right-sided asymmetry of the empirical distributions of random variables was found. This is because in the analyzed timeseries most observations are smaller than their average value. Therefore, the average values are greater than the medians of time series. In turn, the kurtosis values indicated the platykurtic (negative value) and leptokurtic (positive value) empirical distributions for Kamienica Nawojowska and Skawinka respectively.

EBA4SUB Design Hydrographs
The EBA4SUB model was used to determine the design hydrographs to be propagated at the beginning of the investigated river sections. Hydrographs were determined for the following gross rainfall distributions: beta, DVWK, and Chicago. The characteristics of the modelled design hydrographs are presented in Table 2. The hydrographs shapes are shown in Figure 2a    Based on the results, it was found that the highest peak discharge value for the Kamienica Nawojowska river was obtained when a Chicago hyetograph is assumed. The lowest peak discharge value was obtained for the beta hyetograph with parameters 2 and 2. The difference between the highest and lowest value in the peak discharge value is 12%. The values of the wave slenderness coefficient for the obtained hydrographs were close to 1.000. This means that the volumes of the rising and falling parts are similar to each other. The MAPE error values indicate that the peak discharge determined by the EBA4SUB model is close to the value from the statistical method for the hydrograph where precipitation distribution was assumed according to the Chicago method. In the case of the Skawinka River, the highest peak discharge value from the EBA4SUB model was obtained for the hydrograph where the hyetograph was assumed according to the DVWK method. The lowest peak discharge value, similarly to the Kamienica Nawojowska river, was obtained for the beta distribution with parameters 2 and 2. The difference between the highest and the smallest size of the peak discharge is 27%. The slenderness coefficient values are below 1.000. It means that the volumes of the rising part are larger than the corresponding of the falling part. The MAPE error values indicate that the peak discharge determined by the EBA4SUB model is close to the value from the statistical method if precipitation distribution is assumed according to the DVWK method. This statement is close to early work of the authors [22], where it was showed that peak discharges from the EBA4SUB model is closer to quantiles of peak discharges achieved from statistical method than Qp from other rainfall-runoff models, like Snyder or SCS-UH. A second work [21] showed that Qp from the EBA4SUB model where gross rainfall hyetograph was determined using beta distribution was closer to Qp from Pearson III type distribution. A third work [39] shows that the EBA4SUB model can be used as part of new method to generate design peak discharges in mountainous catchments in Poland. The presented results confirm that it is necessary to calibrate and verify hydrological models to eliminate errors that affect the quality of the obtained results [37].

Hydraulic Modeling and Flood-Prone Areas Determination
Flood-prone areas modeled employing the previously described methodologies are shown in Figure 3.   Table 3 summarizes the total areas of flood hazard map covering the investigated river sections with a length of about 9 km. In the case of Skawinka, the flood-prone area obtained with Mike11/NWMA procedure is the largest, with a total area A = 2.380 km 2 and the smallest (FLO-2D with B(2;2)) has the value A = 1.973 km 2 . FLO-2D with B (2;2) flood-prone area is 17.1% smaller than the corresponding obtained with Mike11/NWMA. The largest zone for Kamienica Nawojowska was also generated based on Mike11/NWMA, A = 1.156 km 2 , and the smallest based on the FLO-2D application (again with B(2;2)) presents A = 0.666 km 2 , with a difference of 42.4%. The smallest differences with Mike 11/NWMA were found in the case of Skawinka with the DVWK model and with the B(2;3) model in the case of Kamienica Nawojowska. As for the differences between the EBA4SUB models, the differences are up to 8.2% for Skawinka and 26.8% for Kamienica Nawojowska. A detailed analysis was carried out by identifying areas and focusing on those locations where the zones generated by the FLO-2D model are the same, or underestimated (Figure 4a) or overestimated (Figure 4b), compared to the Mike11/NWMA application. On average, the difference in the areas of FLO-2D-B(2;2) and Mike11/NWMA zones is 17.1%, with local differences that may not occur or be much larger in specific locations.

Discussion
The final result of the procedures, i.e., the map of flood-prone areas, consists of the work of many people, including specialists in spatial analysis, topography, hydrology and, finally, hydraulic modeling. The selection of the hydroinformatics tools affects not only the accuracy of results but, most of all, the computational time of the calculations needed for the map generation. In the case of 1D models, the short calculation time is counterbalanced by the need to obtain the most accurate information from various databases, like a detailed survey of the river cross-section. In the case of 2D models, the increased computation time is counterbalanced by a simplified model preparation and by a limited postprocessing operational time.
In the selected case studies, the flood-prone areas obtained routing on the topography with FLO-2D and the different design hydrographs obtained with the EBA4SUB model slightly differ from each other (Figure 3). In the case of the Mike11/NWMA application, the modelled zones ( Figure 3) significantly differ from the corresponding obtained with the EBA4SUB/FLO-2D application. The nature of these differences is not the same in the case of Skawinka and Kamienica Nawojowska rivers. This circumstance is due to the different structure of the river valley and floodplains. Differences between individual zones are more evident at the local scale. They can be found in both rivers due to the use of different calculation methods.
In order to present the diversity of the flood-prone area over the entire section and, thus, the differences in the use of particular methods, the entire river section was divided into sectors

Discussion
The final result of the procedures, i.e., the map of flood-prone areas, consists of the work of many people, including specialists in spatial analysis, topography, hydrology and, finally, hydraulic modeling. The selection of the hydroinformatics tools affects not only the accuracy of results but, most of all, the computational time of the calculations needed for the map generation. In the case of 1D models, the short calculation time is counterbalanced by the need to obtain the most accurate information from various databases, like a detailed survey of the river cross-section. In the case of 2D models, the increased computation time is counterbalanced by a simplified model preparation and by a limited postprocessing operational time.
In the selected case studies, the flood-prone areas obtained routing on the topography with FLO-2D and the different design hydrographs obtained with the EBA4SUB model slightly differ from each other (Figure 3). In the case of the Mike11/NWMA application, the modelled zones ( Figure 3) significantly differ from the corresponding obtained with the EBA4SUB/FLO-2D application. The nature of these differences is not the same in the case of Skawinka and Kamienica Nawojowska rivers. This circumstance is due to the different structure of the river valley and floodplains. Differences between individual zones are more evident at the local scale. They can be found in both rivers due to the use of different calculation methods.
In order to present the diversity of the flood-prone area over the entire section and, thus, the differences in the use of particular methods, the entire river section was divided into sectors (Kamienica Nawojowska into 23 sectors, Skawinka into 22 sectors, see Figure 5), and the areas of flood-prone zones in the single sectors were calculated. This allowed providing the analysis of zone extension variability along the length of the longitudinal profile of both streams ( Figure 6). The flood-prone area map generated based on the EBA4SUB/FLO-2D modeling approach takes on a larger or smaller area compared to the Mike11/NWMA approach, for the specific sector. The size of the zones in the different sectors can be influenced by both the hydrological and hydraulic model. If the primary source of differences between the zones is due to the hydraulic model, in which a more realistic flood propagation can be expected in the 2D model, two separate trends should result. The first is the difference in flow capacity of the section. The second is created when some parts of the Mike 11 floodplain are always removed in the same way, i.e., in the case that sufficient depth is reached to ensure contact of the zone with the rest of the flood risk area. The observed tendency presenting something different, i.e., variability, indicates the inclusion of other factors found both in the Mike11/NWMA and FLO-2D methods.     In the case of Kamienica Nawojowska, the dimensions of the Mike11/NWMA zones are consistent with the zones obtained from the FLO-2D-B(5;2) application, where the zones are small, almost exclusively narrow, and the section relatively compact and under the influence of pungent and throw pressure manifested by floodplain terraces. In the upper part of the investigated section, where the zones are created above hydrotechnical structures (bridges and water steps), the differences between the zones obtained using the B(2;2), B(2,3) and B(5;2) hydrographs are small and the Mike11/NWMA zone is substantially larger. This is due to the application of insufficient accuracy of the model grid in FLO-2D, which under these conditions generalizes the actual layout of the terrain. This approach is necessary due to the practical usability of the software in modeling of significant fragments of the riverbed and allowing avoiding increasing computation time compared to the 1D model over 1000 times [40]. At the bottom of Kamienica Nawojowska, the regulated riverbed in the FLO-2D model has been generalized, making water access to floodplains possible. The different size of zones due to the characteristics of the flood wave is also visible here: the FLO-2D-B(2;3) application presents the largest zone.
In the case of Skawinka, we are dealing with a watercourse that creates larger zones and the relative differences between the distributions are not reflected in the simulations with FLO-2D, which present similar values. In the outlet section, the differences considering the Mike11/NWMA zones result from merging with the receiver zone following the procedure and in other places with insufficient representation of the terrain layout by generalization of embankments Their absence allows for the creation of zones that the standard procedure does not produce. Regarding the spatial distribution of the flood zone areas B5.2, B2.3, and B2.2, the differences are due to the larger flow capacity in the B2.3 zone. Another aspect for the Mike11/NWMA approach zone is to consider water depth in the grid in relation to the accuracy of DEM. In this case, since the DEM is characterized by an accuracy inside the range up to 0.2 m, if we remove B2.3 grid nodes having the depth less than 0.2 m, the zones would be similar to the Mike11/NWMA approach.
In the traditional approach of generating flood-prone areas, i.e., intersecting the WSM obtained from a 1D model with DEM on floodplain terraces, areas that require expert interpretation emerge [39]. This applies in particular to areas without drainage and those where the connection in the main channel is too small to fill with water during the peak of the wave. The FLO-2D model has a functionality that allows the analysis of details of the generated floodplains [41]. This is particularly important when verifying the physical connection of part of the zone in the main channel with the zone in the floodplain. The advantage of the FLO-2D model is the possibility of performing the procedure in an automated manner and without a great subjective assessment [42]. However, it depends mainly on the accuracy of the data. Moreover, the 2D model allows the user to check the results of hydraulic calculations, e.g., flow depth and local speed, at each pixel of the investigate area. It is noteworthy that this possibility can be time consuming due to the high number of pixels contained in the gridded structure of the model.
Final considerations can be done concerning the calibration issue. Commonly known hydrological approaches used to design flood-prone areas need calibration of parameters in both hydrological and hydraulic models. Calibration is a critical component of model development because parameters generally cannot be determined directly from measurement but are instead inferred indirectly by calibrating the model to observed responses [43,44]. In the last years in Poland, a new approach was applied taking into account hydrological, hydraulic, social, environmental, economic, and implementation aspects. In order to choose a solution, a multicriteria analysis was used, which is based on 13 indicators that use the results of mathematical modeling of flood waves for the formulated variants compared to the current state. This method was applied for the main tributaries of the Vistula River and its results were used in the Flood Risk Management [45][46][47]. The EBA4SUB and FLO-2D model can alternatively be used to determine flood-prone areas in ungauged catchments. The main advantage of the approach is the simplicity of its use. Concerning EBA4SUB, the simplicity means that it can be efficiently and successfully used by practitioners to determine the size of credible flows for estimating flood-prone areas or for designing hydrotechnical constructions. This approach also resolves the problem of having a high uncertainty in hydrological and hydraulic calculations, because parameters of both models have physical representation in real catchments characteristics and can be easily inferred from topography and land use, and in doing so minimizing the subjectivity of the user during the calibration procedure. With the introduction of Water Framework Directive of UE (WFD) in Poland over ten years ago, all water bodies are subject now to the same procedure. Standardized, simplified with use of 1D modeling (Mike 11 was commonly used) rivers and some streams were modelled and hazard zones were generated (and later risk zones). The authors of this study participated in many projects where flood zones were designed in Polish conditions [44][45][46][47][48]. This had two significant consequences. First is that in a short period of time (about 3 years) all rivers were to some degree measured and analyzed. Second is that the attention of researchers was focused elsewhere with the job done thanks to the country offices financial involvement. Now with the perspective of obligatory repeating every several years, the procedure is repelling many attempts of scientific involvement. Years pass and new capabilities in the field of remote survey and modeling are accessible. If we look at the Mike11/NWMA procedure itself, the growing labor costs in Poland, and on the opposite side the better 2D modeling capabilities, we cannot expect this to be the official procedure in the future for much longer. We believe any attempt to use alternative procedures needs to go under the careful eye of the science world. Moreover, practitioners around the world need modern and quick tools to assess flood risk because they are responsible for protecting the infrastructure and maintaining its functionality.

Conclusions
The differences in the flood-prone area maps created using different hydrological and hydraulic models have been investigated in the present manuscript. The obtained results showed that the proposed modeling approach, characterized by the combined use of the hydrological EBA4SUB model and the hydraulic FLO-2D model, represent an accurate alternative, although much more demanding in terms of computational time, to the standard approach used in Poland. Finding a solution capable of the same results without the 2D modeling effort is the main issue of the research today. While FLO-2D requires a high computational effort, the same way as any 2D modeling tool, it also offers fast access to detailed results. Pixel by pixel verification provides evidence of flow capacity, with FLO-2D modeled flood zones that were here generally narrower in respect to the official procedure. One question arises: is the simplifying of the source data in FLO-2D leading to such differences or is the 1D tool generally inaccurate?
Is FLO-2D discredited? In our opinion: no. Extensions of flood-prone areas appear in places where the official procedure results would not be created without manual intervention (for instance the cutting of parts of the flood zones). An essential aspect of any modeling tool is time effectiveness, which also applies to the FLO-2D model. The procedure using FLO-2D can be utilized as a valid alternative to the official approach. Moreover, in this work a second focus concerned the use of a recently developed rainfall-runoff model, named EBA4SUB, for estimating the design hydrograph to be propagated with the hydraulic model. The obtained results are in line with the official hydrological approach, so the use of the combined approach EBA4SUB-FLO-2D in our opinion can pave the way for a fast update of the flood-prone areas at a national level. Further research is needed, in particular regarding the application of the proposed procedure to a large number of catchments, eventually assessing the quality of the results comparing the modeled flood area maps with occurred ones.