Abstract

Wetlands provide many benefits for humans and the natural environment, but land-use changes have reduced their number and areal extent. Interest has grown in examining the landscape to determine those locations where, with minimal effort, it might be possible to develop a mitigation wetland – a location with sufficient water over a sufficient period of time to develop and maintain wetland functioning. This paper proposes a methodology to support the examination of the landscape for mitigation purposes through the application of open channel hydraulics principles to flow over a landscape. The methodology is part of a larger research effort ultimately combining hydrology and hydraulics, along with the landscape processes of infiltration and evapotranspiration, to perform a water balance assessment. Specifically, the methodology is implemented through readily available geographic information system tools along with Python scripts written for this study. The Python scripts automatically extract landscape characteristics from a digital elevation model and calculate hydraulic parameters that are used to determine water surface profiles using the Modified Euler's method. Multiple tests show that the script accurately produces profiles of flow between depressions over a landscape. Such determinations are the first step in understanding where water might exist on the surface to support mitigation wetland functions.

SOFTWARE AND DATA AVAILABILITY

The data, code and tools developed in this research are available via https://github.com/TrauthK/Wetlands/tree/master/Hydraulics/Spatial%20backwater%20Curves. The tools were developed in Python 2.7 with ArcGIS 10.3.

INTRODUCTION

As the name suggests, wetlands are locations that are associated with wetness. This wetness causes a transformation of the soil because of the absence of free oxygen, which in turn supports plants that thrive in this specific anaerobic environment. However, this wetness does not need to occur continuously, and may only be required for a short period of time, as during the growing season for the vegetation of interest.

Wetlands provide multiple beneficial services for the environment and human society. Different kinds of plants and animals are supported by wetlands, and the dynamic relationship among these creatures is very complex, which makes wetlands one of the most productive ecosystems on Earth, providing crucial habitats for many plants and animals, including birds, amphibians and reptiles (Snodgrass et al. 2000; Gibbons et al. 2006). Conversely, human activities, such as construction, pollutant discharges, deforestation and farming, can significantly impact the quality of water resources. For these reasons, interest has grown in the concept of mitigation wetlands – undertaking an effort to create or recreate wetland conditions. Informed decision-making regarding mitigation wetlands requires a thorough understanding of the water balance over time.

For some time, researchers have been interested in the connectivity between wetlands and how wetlands connect to downstream locations (Amezaga et al. 2002; Roe et al. 2009; Lang et al. 2012). Many geographically isolated wetlands can connect with other water bodies during wet seasons. During storms, wetlands receive water from their own contributing areas. After wetlands fill up, surface water may start to flow to neighboring wetlands at lower elevations or to downstream streams and rivers. Quantification of the processes governing hydrologic connectivity between wetlands is thus important for the scientific and decision-making communities (Golden et al. 2014). These wetlands that are connected, even if only periodically, can form a habitat that helps to support plant and animal populations.

Hydraulics background

The flow in a channel with a regular cross section (rectangular or trapezoidal) has been studied for many years (Manning 1891). Modeling the flow over a landscape is a difficult undertaking and introduces additional complexity. As water drains into a depression from the surrounding contributing area, it will exist at a certain depth that will vary based on the size of the contributing area and the type of storm generating the runoff. As water flows between depressions, there exists a water surface profile that connects the depth of water in one depression with the depth of water in the connecting upstream or downstream depression. A water surface profile depicts the changing flow depth in the longitudinal direction of flow. On the hydraulically mild slopes of interest for wetlands mitigation assessments, the profiles of interest are backwater (M1) and drawdown (M2) curves.

Water surface profiles play an important role in engineering practice. Infrastructure, such as channels, sewers or culverts, should be designed with the consideration of gradually varied flow. Different methods to calculate the water surface profiles have been developed, including the Standard Step method, the Direct Step method and the Modified Euler's method (Chow 1959; Sturm 2010). These methods are based on the energy, momentum and continuity equations.

Standard Step method

Currently, the most commonly used method is the Standard Step method which is utilized to compute a one-dimensional water surface profile for gradually varied flow. It is also embedded in HEC-RAS, developed by the U.S. Army Corps of Engineers (Brunner 1995). The flow characteristics at certain cross sections are calculated from the given channel information at another cross section and the reach length between them. The basic equations used in the Standard Step method to compute the water surface profile are the energy equation (Equation (1)) and the momentum equation (Equation (2)): 
formula
(1)
 
formula
(2)
where subscripts 1 and 2 represent cross sections 1 and 2, z is the bed elevation, y is the water depth, V is the velocity, is the velocity correction coefficient, h is the energy loss between sections 1 and 2, g is the gravitational acceleration, A is the cross-sectional area of flow, and is the distance below the free surface of the centroid of the area. The bed elevation, depth and velocity terms are collectively known as energy head, or simply, head. While energy is lost between sections (Equation (1)), momentum is conserved (Equation (2)).

The Standard Step computational procedure is a trial and error process. Based on the given discharge and the boundary and channel conditions, the method assumes an initial water surface elevation for a target cross section. It then recalculates a new water surface by calculating the headloss as a result of friction between the two cross sections. Iterative calculations with different assumed initial water surface elevations continue until the difference in the total heads is within an acceptable range.

Chow (1959) indicated that profile calculations should be carried out from downstream to upstream if the flow is subcritical and from upstream to downstream if the flow is supercritical. Performing calculations in the wrong direction could cause the water profile to diverge from the correct one. The Standard Step method can be used to perform calculations in either direction depending upon the flow conditions.

Direct Step method

The Direct Step method is used to calculate the distance associated with a specific depth change along the channel (Chow 1959). At any cross section of a channel, the total energy head, H, is calculated by 
formula
(3)
as H is differentiated with respect to x, which is the coordinate in the flow direction, Equation (3) can be converted to 
formula
(4)
and 
formula
(5)
in which is the slope of the energy grade line and is the bed slope. Using the finite difference form to approximate the derivative, , Equation (5) results in 
formula
(6)
where is the arithmetic mean slope of the energy grade line between cross sections i and . This method will calculate values step by step by assuming a depth y and a specific energy . Normally, the Direct Step method is recognized as the easiest approach to calculate a water surface, but it requires interpolation if the depth of flow is required at specified locations of interest.

Modified Euler's method

The Modified Euler's method is designed to compute the flow depth changes for specified changes in the distance along a channel. Sturm (2010) reported that this method is more appropriate for natural channels because the channel geometry can be collected from field surveys. It simply treats the water surface curve between two cross sections as a straight line to compute an approximate water depth for the next station. As a result, this method requires a small longitudinal increment to achieve accurate results. The Modified Euler's method moves forward and uses the arithmetic mean of the water surface slope values from two cross sections to obtain a more precise number. Because the slope value is unknown at the beginning, the method uses a predicted water depth and corrects it using the mean slope value from the predicted slope and ending slope value in each interval. This two-step method utilizes prediction (Equation (7)) and correction (Equation (8)) equations: 
formula
(7)
 
formula
(8)
where is the distance between the two calculation locations, and i indicates the current calculated depth and i + 1 indicates the depth at the next location. The method is straight forward enough so that once the channel parameters have been established, the water surface profile for a channel can be calculated using a spreadsheet (Sturm 2010).

Hydraulic modeling

HEC-GeoRAS (U.S. Army Corps of Engineers 2011) has been developed to perform hydraulic calculations utilizing landscape data input through ArcGIS to produce water surface profiles, among other outputs. This tool is able to derive cross sections at specified locations from a digital elevation model (DEM), once a user identifies the centerline of the flow path, the locations of cross-sections of interest and the locations of the channel banks. Manning's n values must be entered manually, as there is no capability to translate a land use/land cover layer to specific roughness values for a channel and the overbanks. While HEC-GeoRAS is strictly hydraulic in its capabilities, it is able to perform calculations based on the hydrologic output from HEC-GeoHMS (U.S. Army Corps of Engineers 2013). HEC-GeoHMS does incorporate various loss methods, such as the Green and Ampt infiltration method (Green & Ampt 1911), but the input parameters are held constant over time and do not represent the actual fluctuations that will occur. For example, the usage of the Green and Ampt Loss method requires the specification of a moisture deficit (the difference between the saturated moisture content and the initial moisture content) that is used for the entire period of calculation. In reality for a site, the soil moisture deficit will change over time with precipitation, infiltration and evapotranspiration (ET). ET itself is not incorporated into HEC-GeoHMS at all. Additionally, a set value for initial water loss must be entered, rather than being calculated as water that is infiltrated until the point of surface saturation.

Arc Hydro (Maidment & Morehouse 2002) provides users with a number of ArcGIS-based tools for landscape-based water resources investigations, which are similar to the spatial analysis techniques embedded in ArcGIS, including watershed delineation and generating flow paths. These capabilities include the ability to provide input data for HEC-GeoRAS. Also included within the category of hydrology and hydraulics is a Green and Ampt tool for the assessment of infiltration impacts. The tool requires the input of Green and Ampt parameters, which include the previously discussed moisture deficit and initial loss parameters that are treated as constants.

When modeling over a landscape, the capability to analyze a DEM and extract useful information for hydrologic and hydraulic purposes is quite mature. Various methods offer the ability to delineate a watershed boundary (Jenson & Domingue 1988) and generate the flow direction and flow path (Quinn et al. 1991). TauDEM (Tarboton 1997) is a set of tools capable of deriving landscape characteristics necessary for hydrologic analyses. Among its capabilities is a method to determine flow direction by proportioning flow between the two lowest downstream cells, recognizing that water flow is not limited to a single receiving cell. Being hydrologic in nature, the TauDEM suite does not incorporate infiltration or ET capabilities, nor does it incorporate the hydraulics necessary for the calculation of a water surface profile. In evaluating specific tools to utilize in a landscape analysis, the Height Above the Nearest Drainage (HAND) (Nobre et al. 2011) is a tool that can classify bare ground cells above a stream into different categories (e.g., 0–0.2, 0.2–0.5 and 0.5–1.0 m) above the stream bed or the water surface elevation). Such a classification could be used in developing inundation maps. However, HAND does not have the capability to work with the varying depths of flow associated with a backwater curve that may exist in slow-moving flood flows, nor is it integrated into either HEC-GeoRAS or Arc Hydro.

Objectives

Wetland mitigation efforts may want to focus on natural depressions that, with minimal effort, could be converted into locations with sufficient water for sufficient periods of time to actually function as wetlands. When investigating the potential for essentially creating new wetlands, the questions that must be answered include how much water will exist at a site, where it will exist and for how long. Answering these questions requires a water balance evaluation that includes precipitation as well as infiltration and ET. Infiltration and ET are processes that are intertwined based on soil moisture content and change throughout the year and with precipitation. Modeling of moisture content must be conducted with a very small time step (e.g., 1 h) and becomes very computationally intensive over the long time periods over which potential mitigation wetland sites must be evaluated (a year or more). Both infiltration and ET are impacted by the distribution of water over a landscape over time, where the distribution is impacted by both hydrologic and hydraulic processes. Thus, a thorough water balance assessment requires an integrated process incorporating not just hydrology and hydraulics but also the land cover and soil processes of infiltration and ET. Land managers engaged in decision-making regarding mitigation wetlands may need to investigate multiple sites, and thus need a methodology that can be implemented relatively quickly and easily.

The objectives of this study are to investigate and hydraulically characterize the flow that could exist between two depressions. The characterization is in terms of a water surface profile that could ultimately be used to determine the aerial extent of water flowing over a landscape. It is not necessary for a location to be continuously covered with water for wetland conditions to develop and be maintained – for plants, it is only necessary for surface inundation during the growing season. As can be seen with HEC-RAS output, the water surface profile in the context of a cross section shows the aerial extent of water on the ground surface. A water surface profile, in combination with an assessment of infiltration, could indicate locations with potential surface saturation.

This study is part of a larger effort to assess the spatial and temporal occurrence of water on a landscape for decision-making with respect to the creation of mitigation wetlands. The overall effort requires the incorporation of precipitation and the hydrology and hydraulics of water flowing over a landscape with the subsurface processes of infiltration and ET. As a landscape-based assessment, a geographic information system (GIS) analysis is required on a scale that is commensurate with the scale of depressions observed on the landscape (<1 ha). Because of the ultimate usage of the modeling results, this determination of water surface profiles must be undertaken within an overall system that incorporates the flow of water over a surface with the possibility of infiltration and ET.

The overall research intent is to combine hydrology and hydraulics on a landscape in a GIS format. It is directed at decision-making for mitigation wetlands, and thus it is necessary to track parameters within a water budget (namely infiltration and ET) over time. This set of capabilities is not currently available in a single suite of software. Thus, it was necessary to embark on a research effort to integrate all these capabilities. Because of the growing utilization of Python scripts within the GIS platform, it was decided to create a Python-based platform that would be able to accommodate the individual scripts addressing overland flow and water surface profiles as well as infiltration and evaporation.

METHODS

The overall project, of which this study is a part, requires the assessment of depressions on the landscape to identify locations that with minimal earthwork changes may have the potential to be converted to functioning wetlands. Such an analysis requires long-term continuous simulation (a year or more) of a water balance that incorporates hydrology and the hydraulics of water moving over a landscape but also the processes of infiltration and ET. Because existing tools have not been identified that meet all of these criteria, the overall plan became to develop a new Python-based tool. The study described in this paper regarding water surface profiles is only one part of the overall effort but must ultimately be able to be integrated into the larger model.

The Modified Euler's method is employed to calculate the water profile caused by backwater or drawdown effects between depressions on a mild slope. The Modified Euler's method was selected because it can directly calculate water surface profiles on the hydraulically mild slopes of interest here. In contrast to the iterative calculation process of the Standard Step method, the Modified Euler's method is more computationally efficient. This computational efficiency is essential given the fact that the overall model will perform a water balance over the period of at least a year, all the while calculating infiltration and ET on an hourly basis.

A work flow is proposed to simulate water surface profile development over a landscape (rather than in a defined channel) using a combination of readily available GIS tools and two Python scripts written for this study. Python is used because of its wide and growing usage in geospatial modeling. In this work flow, GIS tools extract detailed landscape information directly from a high-resolution DEM, which is very helpful, especially for the vast majority of depressions (potential wetland mitigation locations) where there are no field measurements. The two scripts written for this study are then used to derive hydraulic landscape information (i.e., cross-sectional area, wetted perimeter and hydraulic radius) and calculate the water surface profile for the landscape.

The consideration of locations for potential wetland mitigation sites may raise the issue of whether a DEM containing bathymetric information may be necessary for the analysis. Bathymetric DEMs are limited in availability and so could present a difficulty (U.S. Geological Survey 2019). Because the analysis is conducted on depressions, rather than on existing wetlands where there may be surface water, the scarcity of bathymetric DEMs is not an issue in this instance.

The computational process is shown in Figure 1. Once a user of the methodology identifies the general location of interest for analysis, a high-resolution DEM for the area is identified. Within the general area of interest, the user must identify two depressions for specific analysis. Tools readily available in ArcGIS are then used to delineate the contributing areas for each depression and identify the weir cross-section location (pour point) of the upper depression (which controls the outflow from the upstream depression). The starting point is normally recognized as the pour point of the upper depression which is the lowest point of the weir. Correspondingly, the lowest point of the downstream depression is considered as the ending point. The least-cost line (the steepest route) between these two points generated from the DEM is treated as the central line of the flow pathway.

Figure 1

Flow chart of the methodology.

Figure 1

Flow chart of the methodology.

A natural flow path is not a straight line and has variable characteristics at different locations, so the overall flow pathway is divided into several reaches which approximate straight lines and share similar characteristics. A user must identify the number of reaches into which the flow path is divided (fewer reaches for an initial assessment and a greater number of reaches for a more detailed assessment). The slope of each reach is calculated from the reach length and the starting and ending points. With the established reaches, a perpendicular line is created at the midpoint of each reach and is perpendicular to the invisible straight line of the starting and ending point of a reach. Then, a series of points is extracted from the original DEM using a perpendicular line, and the elevation data of these points are used to develop the cross section of each reach. In order to characterize the hydraulic characteristics of the landscape for the range of flow conditions, the hydraulic parameters of area, wetted perimeter and hydraulic radius are calculated for varying depths of flow (i.e., for different water surface elevations). The cross-sectional area of flow at a given depth is calculated as follows: 
formula
(9)
where is the given water surface elevation, h is the elevation of the land surface pixel i, l is the length of each pixel and m represents the number of pixels which are lower than the water surface elevation. The wetted perimeter is calculated as follows: 
formula
(10)
and the hydraulic radius is calculated as the cross-sectional area divided by the wetted perimeter. The last step is to implement the Modified Euler's method and compute the water surface profile for each reach. Based on an initial downstream water depth, the water surface profile is calculated from downstream to upstream. For each reach, the upstream water depth is used as the downstream water depth for the subsequent calculation in the next upper reach. Calculations continue until the upstream pour point is reached, and the set of calculated y values constitutes the water surface profile.

RESULTS AND DISCUSSION

The correct implementation of the Modified Euler's method using the scripts is demonstrated in two different ways. The first demonstration is conducted by using the scripts to extract the hydraulic input parameters for a trapezoidal channel from an artificial landscape and then calculating the backwater curve for an example provided by Sturm (2010). In this fashion, it is easy to ascertain that the Python script is carrying out the calculations properly. The second demonstration is performed for two different precipitation events on two connected depressions in the vicinity of Pershing State Park in Linn County, Missouri. This location is selected because of the historic presence of wetlands and current interest in wetland mitigation.

Demonstration with an established water surface profile

To verify the accuracy of the Python script, an artificial trapezoidal channel is created in ArcGIS using elevation data (Figure 2). It has the same channel geometry as an example provided by Sturm (2010). This trapezoidal channel has a channel slope of 0.001 m/m, a side slope ratio of 2:1 and a bottom width of 8 m. The Manning's n of the channel is 0.025. The Python script detects the DEM of the trapezoidal channel and extracts all geometric information directly from it. The script then calculates the hydraulic parameters (i.e., cross-sectional area, wetted perimeter and hydraulic radius), as they vary by water depth (Table 1).

Table 1

Hydraulic parameters calculated at given water surface elevations

Water surface elevation (m)Water depth (m)Area (m2)Wetted perimeter (m)Hydraulic radius (m)
100.49 
100.5 0.01 0.05 10.85 0.005 
100.51 0.02 0.19 16.74 0.011 
100.52 0.03 0.38 20.46 0.019 
100.53 0.04 0.6 24.49 0.025 
100.54 0.05 0.87 27.9 0.031 
100.55 0.06 1.16 30.38 0.038 
100.56 0.07 1.47 32.55 0.045 
… … … … … 
Water surface elevation (m)Water depth (m)Area (m2)Wetted perimeter (m)Hydraulic radius (m)
100.49 
100.5 0.01 0.05 10.85 0.005 
100.51 0.02 0.19 16.74 0.011 
100.52 0.03 0.38 20.46 0.019 
100.53 0.04 0.6 24.49 0.025 
100.54 0.05 0.87 27.9 0.031 
100.55 0.06 1.16 30.38 0.038 
100.56 0.07 1.47 32.55 0.045 
… … … … … 
Figure 2

Recreation of a trapezoidal channel in ArcGIS.

Figure 2

Recreation of a trapezoidal channel in ArcGIS.

Table 2 shows the results of the backwater calculations from Sturm (2010) and the Python script. Column 1 shows the distance from the starting point 0 to the upstream ending location 20 m away. At the beginning of the calculation, the x intervals are very small in order to obtain precise results where there may be significant curvature of the flow, as with an M2 water surface profile. Column 2 lists the water depths as given by Sturm (2010). Column 3 depicts the water depths calculated using a spreadsheet created for this purpose, while Column 4 displays the depths calculated using the Python script. A spreadsheet is used as an interim step in the coding to establish the calculational strategy. The very minor differences between these two processes indicate that the Modified Euler's method can be implemented using this Python script.

Table 2

A comparison of Excel Python calculations

x (m)Sturm (2010) y1 (m) By Excely1 (m) By Python
1.0300 1.0300 1.0300 
0.005 1.0575 1.0575 1.0576 
0.01 1.0578 1.0578 1.0580 
0.02 1.0585 1.0585 1.0586 
0.04 1.0596 1.0596 1.0598 
0.06 1.0608 1.0608 1.0610 
0.08 1.0619 1.0619 1.0621 
0.1 1.0629 1.0629 1.0630 
0.2 1.0677 1.0677 1.0676 
0.4 1.0757 1.0757 1.0751 
0.6 1.0824 1.0824 1.0815 
0.8 1.0883 1.0883 1.0869 
1.0936 1.0936 1.0922 
1.1147 1.1147 1.1142 
1.1446 1.1446 1.1434 
1.1669 1.1669 1.1650 
1.1852 1.1852 1.1826 
10 1.2010 1.2010 1.1978 
20 1.2596 1.2596 1.2545 
x (m)Sturm (2010) y1 (m) By Excely1 (m) By Python
1.0300 1.0300 1.0300 
0.005 1.0575 1.0575 1.0576 
0.01 1.0578 1.0578 1.0580 
0.02 1.0585 1.0585 1.0586 
0.04 1.0596 1.0596 1.0598 
0.06 1.0608 1.0608 1.0610 
0.08 1.0619 1.0619 1.0621 
0.1 1.0629 1.0629 1.0630 
0.2 1.0677 1.0677 1.0676 
0.4 1.0757 1.0757 1.0751 
0.6 1.0824 1.0824 1.0815 
0.8 1.0883 1.0883 1.0869 
1.0936 1.0936 1.0922 
1.1147 1.1147 1.1142 
1.1446 1.1446 1.1434 
1.1669 1.1669 1.1650 
1.1852 1.1852 1.1826 
10 1.2010 1.2010 1.1978 
20 1.2596 1.2596 1.2545 

Water surface profile demonstration on a landscape

The procedure of calculating the water surface profile at an actual landscape location is demonstrated in this section.

Process

Digital elevation model

A user of the methodology must identify a location of interest for potential mitigation wetlands. Pershing State Park in Linn County, Missouri currently contains wetlands. Because of the proximity to wetland locations and the potential for habitat connectivity with mitigation wetlands, this demonstration is undertaken in a nearby location. A 1-m resolution DEM is created based on the LiDAR data from the Missouri Spatial Data Information Service (MSDIS 2019). The scale of the demonstration location was specifically chosen to be representative of the scale of the areas over which this process would be implemented. Natural wetlands may be small, so this methodology is intended to be used at a very small scale (tens of hectares) to identify those locations that could possibly support small wetland mitigation sites. In addition, the simplicity of two depressions was specifically selected to mimic the assessments decision-makers may be conducting.

Identify upstream and downstream depressions of interest

A user of the methodology must identify specific depressions for analysis. The study area includes two depressions connected by a landscape with a mild natural slope (Figure 3). In the DEM, lighter areas represent higher elevations, while darker areas represent lower elevations. During a rainfall event, these two depressions would start to accumulate surface runoff from their respective contributing areas. Once the runoff would fill up the upper depression, excess runoff would flow out of it and contribute to the lower depression. The lower depression, with a larger contributing area, may have a deeper water depth which could cause drawdown or backwater effects.

Figure 3

Procedures to extract landscape information from a DEM. (a) Identify the contributing area; (b) create the central line of flow path; (c) divide the flow path into several reaches; and (d) generate a perpendicular line for each reach.

Figure 3

Procedures to extract landscape information from a DEM. (a) Identify the contributing area; (b) create the central line of flow path; (c) divide the flow path into several reaches; and (d) generate a perpendicular line for each reach.

Identify contributing basin boundaries

The contributing area for each subbasin can help to explain the flow accumulation. The first step is using the embedded spatial analysis tools in ArcGIS to generate boundaries for each depression (Figure 3(a)). In this example, depressions are identified from the DEM.

Create weir cross-section locations

With sufficient rain, surface water will flow out of the upper depression and contribute to the lower depression. The starting point for flow (Figure 3(a)) is located at the lowest point along the boundary of the upper depression (i.e., the pour point of the upstream watershed). The ending point is the lowest point in the lower depression (i.e., the sink).

Create the center line of the flow

Surface runoff always flows along the steepest slope over the landscape. The distance tools in ArcGIS can create a cost map based on the slope of the landscape and then generate the least-cost path (Figure 3(b)) which can be recognized as the centerline of the flow. For this study site, the elevation of the starting point is 228.65 m and the elevation of the ending point is 227.52 m. The total length of the flow is 229 m, and the slope of the overall flow pathway is 0.0049 m/m.

Split the channel into reaches, create a cross section for each reach and extract the cross-section parameters

In general, separate reaches are identified where the hydraulic characteristics of slope, land cover or geometry are consistent. For purposes of representing the different landscape characteristics in this demonstration, the flow path is divided into three reaches of equal length (Figure 3(c)). Each reach has different slope and cross-section parameters which are automatically extracted from the DEM using the Python script. First, the script creates perpendicular lines from the middle point of each reach (Figure 3(d)). The elevation data along each perpendicular line are extracted from the DEM. Second, the cross section for each reach is built from these elevation data. The width of each cross section is set at 100 m and can be modified by the user as necessary. In standard hydraulic convention, reach lengths are numbered starting from the downstream-most location.

Table 3 shows the elevation data of the starting and ending points for each reach, along with the length and slope. Figure 4 shows how the geometry of different cross sections can vary significantly, having a great impact on the hydraulic calculations. Compared with Reaches 1 and 3, the geometry of Reach 2 has steeper sides (impacting the cross-sectional area of flow), resulting in an increasing depth when water flows from Reach 3 to Reach 2 and a decreasing depth when the water flows from Reach 2 to Reach 1.

Table 3

Reach parameters

Length (m)
Ground surface (m)
Slope (m/m)
StartEndDifferenceStartEndDifference
Reach 1 76.35 76.35 227.52 227.80 0.28 0.0036 
Reach 2 76.35 152.7 76.35 227.80 228.16 0.36 0.0047 
Reach 3 152.7 229.05 76.35 228.16 228.65 0.50 0.0065 
Length (m)
Ground surface (m)
Slope (m/m)
StartEndDifferenceStartEndDifference
Reach 1 76.35 76.35 227.52 227.80 0.28 0.0036 
Reach 2 76.35 152.7 76.35 227.80 228.16 0.36 0.0047 
Reach 3 152.7 229.05 76.35 228.16 228.65 0.50 0.0065 
Figure 4

Cross-section geometry of each reach.

Figure 4

Cross-section geometry of each reach.

Apply Modified Euler's method to generate the water depth profile

The water surface profile for each reach is calculated from downstream to upstream in a continuous process. Starting with a given water depth, the script calculates the water depth at the end of Reach 1 and uses that depth as the starting water depth for Reach 2. Here, the profiles are calculated for two different flow circumstances (i.e., from a more intense Storm 1 generating greater flow and from a less intense Storm 2 generating less flow).

Results

Storm 1

Figure 5 shows the results from Storm 1 (with a larger discharge and a deeper downstream depth) as generated by the Python script. In this test, the discharge of the flow is 200 m3/s with an initial depth of 1.3 m at the bottom of the lower depression. The roughness coefficient for all three reaches is set to 0.025. The thicker solid line indicates the channel bottom, and the right vertical axis is elevation. The dotted line is the calculated uniform depth. The uniform depth changes between the different reaches because the landscape characteristics are different between the reaches. The thinner solid line represents the water surface profile, which always varies in terms of approaching uniform depth asymptotically. The 1.3 m depth in the downstream depression is greater than the calculated uniform depth of 0.90 m. Therefore, the water surface must decrease with an M1 curve in the process of approaching the equilibrium condition of uniform depth. Even as the uniform depth changes, as from Reach 1 to Reach 2, the water surface is shown to be continuous, as would actually occur on the landscape. The transition from Reach 2 to Reach 3 shows the continuous water surface, with the shape adjusting from the steeper sides of Reach 2 to the flatter sides of Reach 3. For most of Reach 3, the plot shows the water surface profile coinciding with the uniform depth line, indicating that the water surface profile is very close to the uniform depth in its asymptotic approach, as would be expected.

Figure 5

Water surface profile for Storm 1.

Figure 5

Water surface profile for Storm 1.

Storm 2

For Storm 2, the discharge is 5.0 m3/s with a given depth 0.3 m at the bottom of the lower depression. The roughness coefficient for Reach 2 is set to 0.05, and Reach 1 and 3 are unchanged at 0.025. In this scenario, the flow depth in Reach 1 decreases gradually following an M1 curve (Figure 6) because the water depth in the downstream wetland is larger than the uniform depth of 0.20 m in the region connecting the two wetlands. Unlike the plot from Storm 1, there is an M2 curve in the transition between Reach 1 and 2, as the water surface now needs to increase to approach uniform depth. This condition exists because a larger roughness coefficient is assumed for Reach 2 resulting in a deeper uniform depth of flow, even without changing any of the other parameters. Once in Reach 3, water depth again follows an M1 curve to decrease the depth of flow to approach the lower uniform depth.

Figure 6

Water surface profile for Storm 2.

Figure 6

Water surface profile for Storm 2.

CONCLUSIONS

This paper introduces an efficient methodology for using the Modified Euler's method to calculate water surface profiles over a natural landscape on the mild slopes connecting wetlands. Python scripts have been created to perform the calculations, including the extraction of landscape data directly from a DEM, eliminating the need for a field survey at every potential depression of interest. The efficiency of the methodology allows for widespread scoping calculations to assess locations for more detailed investigations for wetland protection or mitigation. Its simplicity and automation can help users with little hydraulics background to calculate profiles for decision-making.

Multiple tests have been performed to verify the accuracy of this methodology. The Python scripts created for this project automatically extract the necessary data and compute all necessary hydraulic parameters. The first test recreated a drawdown curve demonstration in a textbook. This was accomplished by generating an artificial landscape matching the trapezoidal channel of the example. The Python scripts extracted the necessary landscape parameters from the DEM, calculated the appropriate hydraulic parameters and developed the drawdown curve, exactly matching the published results. The second test was performed on a natural landscape with a discharge that would be associated with an intense storm. This set of conditions produced a sequence of backwater (M1) curves expected as landscape characteristics change, causing the water surface profile to adjust accordingly in moving from one reach to another. The third test was conducted on the same natural landscape (while changing one surface roughness value) but with a less intense storm. The lower flows associated with these more common, less intense storms, are of particular interest for continued wetland functioning. With the lower discharge, the script demonstrated the ability to calculate a water surface profile that incorporates both backwater (M1) and drawdown (M2) curves. For the development of the methodology, in order to focus on the surface processes, infiltration was assumed to be negligible.

This effort addresses the challenge of applying open channel principles to quantify the flow characteristics (i.e., a water surface profile) over a landscape in contrast to the current focus on flow in defined channels. It is able to accommodate the inherent variability in the landscape elevations and land covers. It is now possible to make assessments that could provide a first-step screening tool for wetland scientists to assist in the selection of locations that may warrant further investigation with respect to the preservation, restoration or mitigation of wetlands.

Now that the water surface profile methodology has been developed and demonstrated, it will be appropriate to extend the depiction of the water surface profile on the landscape to include the extent of the land surface that would be covered by water. Such water on the surface would be available for infiltration to support wetland function.

ACKNOWLEDGEMENT

This work was supported by US Environmental Protection Agency, Region 7 Wetland Program Development Grant CD97748001.

REFERENCES

REFERENCES
Amezaga
J.
Santamaría
L.
Green
A. J.
2002
Biotic wetland connectivity – supporting a new approach for wetland policy
.
Acta Oecologica
23
(
3
),
213
222
.
Brunner
G. W.
1995
HEC-RAS River Analysis System. Hydraulic Reference Manual. Version 1.0
.
Hydrologic Engineering Center
,
Davis, CA
.
Chow
V. T.
1959
Open Channel Hydraulics
.
McGraw-Hill Book Company, Inc.
,
New York
,
USA
.
Gibbons
J. W.
Winne
C. T.
Scott
D. E.
Willson
J. D.
Glaudas
X.
Andrews
K. M.
Todd
B. D.
Fedewa
L. A.
Wilkinson
L.
Tsaliagos
R. N.
2006
Remarkable amphibian biomass and abundance in an isolated wetland: implications for wetland conservation
.
Conservation Biology
20
(
5
),
1457
1465
.
Golden
H. E.
Lane
C. R.
Amatya
D. M.
Bandilla
K. W.
Kiperwas
H. R.
Knightes
C. D.
Ssegane
H.
2014
Hydrologic connectivity between geographically isolated wetlands and surface water systems: a review of select modeling methods
.
Environmental Modelling & Software
53
,
190
206
.
Green
W. H.
Ampt
G. A.
1911
Studies on soil physics, 1. The flow of air and water through soils
.
Journal of Agricultural Science
4
,
11
24
.
Jenson
S. K.
Domingue
J. O.
1988
Extracting topographic structure from digital elevation data for geographic information system analysis
.
Photogrammetric Engineering and Remote Sensing
54
(
11
),
1593
1600
.
Lang
M.
McDonough
O.
McCarty
G.
Oesterling
R.
Wilen
B.
2012
Enhanced detection of wetland-stream connectivity using LiDAR
.
Wetlands
32
(
3
),
461
473
.
Maidment
D. R.
Morehouse
S.
2002
Arc Hydro: GIS for Water Resources
.
ESRI, Inc.
,
Redlands, CA
.
Manning
R.
1891
On the flow of water in open channels and pipes
.
Transactions of the Institution of Civil Engineers of Ireland
20
,
161
207
.
Missouri Spatial Data Information Service
2019
MSDIS/NAIP2010. Available from: http://moimagery.missouri.edu/arcgis/rest/services/MSDIS/NAIP2010/ImageServer (accessed July 2019)
.
Nobre
A. D.
Cuartas
L. A.
Hodnett
M.
Rennó
C. D.
Rodrigues
G.
Silveira
A.
Waterloo
M.
Saleska
S.
2011
Height above the nearest drainage – a hydrologically relevant new terrain model
.
Journal of Hydrology
404
(
1–2
),
13
29
.
Quinn
P.
Beven
K.
Chevallier
P.
Planchon
O.
1991
The prediction of hillslope flow paths for distributed hydrological modelling using digital terrain models
.
Hydrological Processes
5
(
1
),
59
79
.
Sturm
T. W.
2010
Open Channel Hydraulics
.
McGraw-Hill
,
New York
.
U.S. Army Corps of Engineers
2011
HEC-GeoRAS: GIS Tools for Support of HEC-RAS Using ArcGIS, User's Manual, Version 4.3.93
.
U.S. Army Corps of Engineers Hydrologic Engineering Center
,
Davis, CA
.
U.S. Army Corps of Engineers
2013
HEC-GeoHMS: Geospatial Hydrologic Modeling Extension, User's Manual, Version 10.1
.
US Army Corps of Engineers Hydrologic Engineering Center
,
Davis, CA
.
U.S. Geological Survey
2019
Available from: https://wwwbrr.cr.usgs.gov/gstl/project-lidar.html (accessed December 2019)
.