In cooperation with the Bayfield County Land and Water Conservation Department and the U.S. Fish and Wildlife Service

Water-Resources Investigations Report 03–4130

Simulation of Ground-Water Flow and Rainfall Runoff with Emphasis on the Effects of Land Cover, Whittlesey Creek, Bayfield County, Wisconsin, 1999–2001

By Bernard N. Lenz, David A. Saad, and Faith A. Fitzpatrick

The PDF for the report is 7,399 kb

Table of Contents

Abstract

Introduction

Methods

Field Data

Ground-Water-Flow Simulations

Rainfall-Runoff Simulations

Summary and Conclusions

References Cited

Figures

Figure 1. Location map of Whittlesey Creek, Wis., study area and sampling a...

Figure 2. Land-cover data and subbasins used in the Whittlesey Creek, Wis.,...

Figure 3. Cross section showing generalized geology and conceptual model of...

Figure 4. Cross section showing the geology of the southern part of the Whi...

Figure 5. Hydrologic features simulated with analytic elements for the GFLO...

Figure 6. Boundary conditions and location of calibration targets used in t...

Figure 7. Block diagram showing the MODFLOW model hydraulic conductivity/re...

Figure 8. Digital Elevation Model of land surface for the Whittlesey Creek,...

Figure 9. Soils in the Whittlesey Creek, Wis., surface-water-contributing b...

Figure 10. Measured hourly stage data used in timing peak stages through Wh...

Figure 11. Land cover in the Whittlesey Creek, Wis., surface-water-contribu...

Figure 12. Land cover in the Whittlesey Creek, Wis., surface-water-contribu...

Figure 13. Base flows for Whittlesey Creek and tributaries, Wis., August 20...

Figure 14. Streambed heads for Whittlesey Creek and tributaries, Wis., Augu...

Figure 15. Streambed temperatures for Whittlesey Creek and tributaries, Wis...

Figure 16. MODFLOW ground-water-flow model parameter sensitivities from UCO...

Figure 17. Relation of measured and modeled (A) head and (B) stream flux (b...

Figure 18. Whittlesey Creek, Wis., ground-water-contributing area based on ...

Figure 19. Particle pathlines tracked backwards from Whittlesey Creek, Wis....

Figure 20. Hydrographs of daily mean flow from the Whittlesey Creek, Wis., ...

Figure 21. Flood hydrographs of measured daily mean and instantaneous (15-m...

Figure 22. Flood hydrographs of daily mean flow from the Whittlesey Creek, ...

Tables

Conversion Factors, Vertical Datum, and Abbreviated Units of Measurement

Table 1. Streamflow statistics for U.S. Geological Survey streamflow-gaging...

Table 2. Estimated range or value of aquifer properties, recharge, and stre...

Table 3. Physical characteristics of subbasins used in Whittlesey Creek, Wi...

Table 4. Soil and Water Assesment Tool (SWAT) model input parameter values ...

Table 5. Percentages of land-cover types used in the Whittlesey Creek, Wis....

Table 6. Results of field investigations of base flow and physical characte...

Table 7. Comparison of base flow and drainage area for Whittlesey Creek and...

Table 8. Preliminary estimate and UCODE-optimized parameter values used in ...

Table 9. Head statistics and simulated base flows at selected locations for...

Table 10. Simulated changes in base flow for Whittlesey Creek, Wis., as a r...

Table 11. Daily mean and intantaneous peak flows for U.S. Geological Survey...

Table 12. Average annual runoff and flood peaks (based on daily mean flows ...

Table 13. Percentage of total precipitation that becomes evapotranspiration...

Table 14. Percentage of streamflow at the Whittlesey Creek, Wis., U.S. Geol...

Table 15. Average ground-water infiltration rates by soil type in the Whitt...

Conversion Factors, Vertical Datum, and Abbreviated Units of Measurement

Temperature in degrees Celsius (°C) may be converted to degrees Fahrenheit (°F) as follows:

°F=(1.8°C)+32

Vertical coordinate information is referenced to the National Geodetic Vertical Datum of 1929 (NGVD 29).

Abbreviated units of measurement used in this report: Specific conductance is given in microsiemens per centimeter at 25°C (µS/cm). Permeability rate is given in inches per hour (in/h). Average annual runoff and recharge are given in inches per year (in/yr).

Abstract

The effects of land cover on flooding and base-flow characteristics of Whittlesey Creek, Bayfield County, Wis., were examined in a study that involved ground-water-flow and rainfall-runoff modeling. Field data were collected during 1999-2001 for synoptic base flow, streambed head and temperature, precipitation, continuous streamflow and stream stage, and other physical characteristics. Well logs provided data for potentiometric-surface altitudes and stratigraphic descriptions. Geologic, soil, hydrography, altitude, and historical land-cover data were compiled into a geographic information system and used in two ground-water-flow models (GFLOW and MODFLOW) and a rainfall-runoff model (SWAT). A deep ground-water system intersects Whittlesey Creek near the confluence with the North Fork, producing a steady base flow of 17–18 cubic feet per second. Upstream from the confluence, the creek has little or no base flow; flow is from surface runoff and a small amount of perched ground water. Most of the base flow to Whittlesey Creek originates as recharge through the permeable sands in the center of the Bayfield Peninsula to the northwest of the surface-water-contributing basin. Based on simulations, model-wide changes in recharge caused a proportional change in simulated base flow for Whittlesey Creek. Changing the simulated amount of recharge by 25 to 50 percent in only the ground-water-contributing area results in relatively small changes in base flow to Whittlesey Creek (about 2–11 percent). Simulated changes in land cover within the Whittlesey Creek surface-water-contributing basin would have minimal effects on base flow and average annual runoff, but flood peaks (based on daily mean flows on peak-flow days) could be affected. Based on the simulations, changing the basin land cover to a reforested condition results in a reduction in flood peaks of about 12 to 14 percent for up to a 100-yr flood. Changing the basin land cover to 25 percent urban land or returning basin land cover to the intensive row-crop agriculture of the 1920s results in flood peaks increasing by as much as 18 percent. The SWAT model is limited to a daily time step, which is adequate for describing the surface-water/ground-water interaction and percentage changes. It may not, however, be adequate in describing peak flow because the instantaneous peak flow in Whittlesey Creek during a flood can be more than twice the magnitude of the daily mean flow during that same flood. In addition, the storage and infiltration capacities of wetlands in the basin are not fully understood and need further study.

Introduction

Whittlesey Creek (fig. 1) is considered a regionally important stream for fish spawning and rearing (Gardner, 1994; Wisconsin Department of Natural Resources, 1996). The U.S. Fish and Wildlife Service established the Whittlesey Creek National Wildlife Refuge in 1999 and acquired 97 of the 570 acres proposed for the refuge that surrounds the mouth of Whittlesey Creek. The main objectives for the refuge are to protect and restore habitat in Whittlesey Creek and surrounding creeks for migration, spawning, and rearing of trout and salmon from Lake Superior and to protect important bird nesting areas. However, as in many Wisconsin tributaries to Lake Superior, aquatic habitat in Whittlesey Creek is possibly degraded because of accelerated runoff and associated sedimentation problems that have potentially resulted from agriculture, forestry practices, and development of roads (Red Clay Interagency Committee, 1967). In addition, the sources for abundant base flow in Whittlesey Creek are not well understood. This study was done in cooperation with the Bayfield County Land and Water Conservation Department and the U.S. Fish and Wildlife Service to identify land-cover effects on the flow characteristics of the creek.

Figure 1. Location map of Whittlesey Creek, Wis., study area and sampling and monitoring sites, 1999-2001.

Purpose and Scope

The purpose of this report is to present results of ground-water-flow and rainfall-runoff simulations as they related to historical, present, and potential future land-cover effects on flooding and base-flow characteristics of Whittlesey Creek. The objectives of the report include the following:

1. Quantify current streamflow characteristics—base flow, floods, and seasonal fluctuations,
2. quantify contribution of runoff and ground water to streamflow,
3. identify the ground-water contributing area¹ and gaining and losing reaches of Whittlesey Creek,
4. estimate presettlement hydrologic conditions
5. quantify the relation between present land-cover characteristics and base flows and floods in Whittlesey Creek, and
6. identify areas (if any) within the Whittlesey Creek surface-water contributing basin that have the greatest effect on runoff.

¹ In this report, "ground-water-contributing area" refers to the two-dimensional land-surface area projected from the three-dimensional subsurface volume of water discharging to the streambed of Whittlesey Creek.

The scope of work included conducting base-flow, streambed-head, and streambed-temperature surveys, monitoring stream stage and flow of Whittlesey Creek, and development of ground-water-flow and rainfall-runoff models for Whittlesey Creek. The field surveys were done in August 1999 and April and August 2000.

Description of Study Area

The Whittlesey Creek drainage basin, as delineated from topographic divides, covers an area of about 38 mi² (fig. 1). This basin includes surface-water-noncontributing and contributing areas (fig. 1). The noncontributing basin of Whittlesey Creek is in the Bayfield Highlands and is composed of sandy deposits with no surface-drainage features. The surface-water-contributing basin is only about 7.4 mi². Hydrologic features in the contributing basin include the main stem of Whittlesey Creek, the North Fork Whittlesey Creek, and numerous unnamed tributaries. Streamflow near the mouth of Whittlesey Creek is a consistent 17 to 18 ft³/s at base flow, with flood peaks of over 500 ft³/s. Just upstream from the confluence of the main stem and North Fork, the main stem becomes ephemeral or has base flow of less than 1 ft³/s.

The surface-water-contributing basin of Whittlesey Creek is 65.5 percent forested, 29.1 percent pasture/grassland/hay, 3.4 percent brushy rangeland, and 2 percent wetland based on 1992-93 land-cover data (fig. 2) (Reese and others, 2002). The basin generally has sloping plains in the uplands and deeply incised valleys. Surface water generally moves as sheet flow in the contributing basin until reaching a ditch or gully leading to the channel network. The gullies and the stream channels are generally steeply sloped, so water passes rapidly through the basin after reaching a ditch or gully. Near the mouth at Lake Superior, the channel slope flattens considerably. For example, the North Fork Whittlesey Creek has an average slope of 0.02, whereas the average slope of Whittlesey Creek below the confluence with the North Fork Whittlesey Creek is 0.005.

Figure 2. Land-cover data and subbasins used in the Whittlesey Creek, Wis., Soil and Water Assessment Tool (SWAT) model, 1992-93. Land cover from WISCLAND data from 1992-93 from Reese and others (2002). Drainage divides for subbasins were derived from 10-meter Digital Elevation Model data (Benchmark GIS, 2001).

Bedrock in the study area includes the Precambrian Bayfield Group, which consists dominantly of sandstone and siltsone and locally abundant shale and conglomerate (Morey and Ojakangas, 1982; Mudrey and others, 1982; Cannon and others, 1996). The Bayfield Group is estimated to be up to 6,900 feet thick (Morey and Ojakangas, 1982) and crops out at several sites along the Lake Superior shoreline near the study area. The Bayfield Group is overlain by glacial, glaciolacustrine, and fluvial deposits, including sandy till of the Copper Falls Formation in the Bayfield Highlands and the clayey Miller Creek Formation at lower altitudes towards Lake Superior (Clayton, 1984) (fig. 3). Unconsolidated deposits are several hundred feet thick in the Bayfield Highlands and generally thin toward Lake Superior. Where present, the Miller Creek Formation overlies the Bayfield Group and the Copper Falls Formation up to an altitude of about 1,100 ft. Locally, coarse-grained beach deposits overlie the Miller Creek Formation along relict early Holocene shorelines of ancient Lake Superior. Well logs commonly indicate alternating layers of sand and clay in the area underlain by the Miller Creek Formation (fig. 4).

Figure 3. Cross section showing generalized geology and conceptual model of the deep ground-water flow system of the Bayfield Peninsula, Wis.

Figure 4. Cross section showing the geology of the southern part of the Whittlesey Creek, Wis., study area. Location of section is shown on figure 1.

There are two ground-water flow systems in the study area: a deep flow system and a shallow flow system. In the deep system, which underlies the entire study area, ground water moves primarily through the Copper Falls Formation and the upper part of the Bayfield Group (fig. 3). The deep system receives a relatively large amount of recharge through the Copper Falls Formation in the Bayfield Highlands. The deep system is unconfined in the higher altitudes and confined at lower altitudes near Lake Superior, where the Miller Creek Formation overlies it. The deep system discharges to Lake Superior and deeply incised streams, such as Whittlesey Creek. The extent of the shallow system is difficult to delineate but probably includes much of the area underlain by the Miller Creek Formation. The alternating layers of clay and sand in this area result in isolated areas of perched water that may be separated from the deep system by 100 ft or more near the headwaters of Whittlesey Creek. The shallow system receives less recharge than the deep system because the Miller Creek Formation is less permeable than the Copper Falls Formation. Some ground water from the shallow system discharges to Whittlesey Creek, and some likely recharges the deep system. The two systems merge where the deep system intersects the stream channel near the confluence of the main stem and North Fork (fig. 4). The deep system and the nonperched part of the shallow system are included in simulations of ground-water flow, whereas perched flow is included in simulations of surface-water flow.

Climate

The regional climate is dominated by three air masses: the cold, dry polar continental air mass from the arctic and northwestern Canada; the warm, moist maritime tropical air masses from the Gulf of Mexico; and the mild, temperate maritime polar air mass from the Pacific (Eichenlaub, 1979). In addition, the moderating effects of Lake Superior influence the local climate. Average annual precipitation is 32 in. (Gebert, 1986), including an average annual snowfall of about 70 in. (Eichenlaub, 1979). Average annual evapotranspiration is estimated to be 17.7 in. (Young and Skinner, 1974). Daily mean July maximum temperature is 25.5°C, and daily mean January temperature is –4.5°C. There are an average of 120 frost-free days annually (Phillips and McCulloch, 1972). Because of the cool temperatures and short growing season, the ratio of precipitation to evapotranspiration is high compared to that in the southern part of Wisconsin.

Study Design

Ground-water and surface-water contributions to streamflow in Whittlesey Creek were simulated by use of separate models. Ground-water-flow models were used to understand the base-flow characteristics and delineate the ground-water-contributing area for the stream. A rainfall-runoff model was used to understand the surface-water characteristics of the basin and predict how land-cover changes would affect annual runoff and flood peaks.

Development of the ground-water flow models for the Whittlesey Creek study area followed the stepwise modeling approach described by Haitjema (1995). Stepwise modeling starts with a preliminary model (constructed from available data) that is used to guide additional data collection. Based on results of the preliminary model, additional data needs are identified, field data are collected, and the model is updated and improved. These steps are repeated until reasonable model results are achieved. Field data collected for the model included multiple measurements of base flow, streambed head and temperature, and other physical characteristics of Whittlesey Creek and adjacent streams.

A continuous rainfall-runoff model, called Soil and Water Assessment Tool (SWAT) (Neitsch and others, 2001), was used to simulate annual runoff and flood peaks in Whittlesey Creek. The model was calibrated to data collected at a USGS streamflow-gaging station near the mouth and to stage data collected at various upstream locations in the basin (fig. 1). After calibration, the SWAT model was run for five simulations of varied land cover: (1) present (1992-93), (2) land cover from 1928, representing the peak of agricultural land cover, (3) a forested basin, representing pre-European settlement prior to about 1870, (4) a reforested basin where forest is substituted for grassland and agriculture land, and (5) and a developed basin in which 25 percent of the basin is developed into urban residential land.

Acknowledgments

The project was made possible through the collaborative efforts of Maureen Gallagher, formerly with the U.S. Fish and Wildlife Service, Ashland, Wis., and Sandra Schultz, formerly with Bayfield County Land and Water Conservation Department. Ulf Gafvert, formerly with the Natural Resources Conservation Service, provided digital elevation model data for the rainfall-runoff model, assisted with field data collection, and provided much insight into the local soil and hydrologic characteristics of the basin. This study could not have happened without the cooperation of riparian landowners who provided access for base-flow surveys and for checking and downloading data from stage recorders.

In addition, the authors of this report acknowledge the late Jeff Carlson (trout conservationist and active member of Trout Unlimited) for his great efforts at preserving current fish habitat and his dedication toward reestablishing brook trout in the Whittlesey Creek National Wildlife Refuge. Jeff and other volunteers built the gage house for the streamflow-gaging station on Whittlesey Creek that shelters the stream monitoring equipment.

Methods

Field investigations included measurements of base flow, pH, specific conductance, water temperature, streambed head, and streambed temperature from 1999 through 2000 in Whittlesey Creek, Boyd Creek, Sioux River, and North Fish Creek (fig. 1). Surficial deposits were probed and core samples were collected at two locations in the study area by use of a hydraulically powered coring device. Continuous streamflow was monitored on Whittlesey Creek at one streamflow-gaging station. Stream stages were monitored at five remote sites on Whittlesey Creek from 2000 through 2001.

Ground-water flow for the study area was simulated using an analytic-element model (GFLOW) and a finite-difference model (MODFLOW). GFLOW (Haitjema, 1995) was used to construct a simplified one-layer screening model; MODFLOW (McDonald and Harbaugh, 1988) was used to construct a more detailed three-dimensional ground-water-flow model. The GFLOW model was used for hypothesis testing, estimating recharge, and improving the MODFLOW model. The analytic-element methodology is not as widely utilized for numerical modeling as finite-difference techniques (Hunt and Krohelski, 1996; Hunt and others, 2000). A complete description of analytic elements is beyond the scope of this report; a brief description is given below. Strack (1989) and Haitjema (1995) provide detailed discussions of this method.

An infinite aquifer is assumed in analytic element modeling. The problem domain does not require a grid or involve interpolation between cells. To construct an analytic-element model, features that affect ground-water flow (such as wells and surface-water bodies) are entered as mathematical elements or strings of elements. The amount of detail specified for the features depends on distance from the area of interest. Each element is represented by an analytic solution. The effects of these individual solutions are added together to arrive at a solution for the ground-water-flow system. Because the solution is not confined to a grid, heads and flows can be computed anywhere in the model domain without nodal averaging. In the GFLOW model, the analytic elements are two-dimensional and are used to simulate only steady-state conditions.

A continuous rainfall-runoff model, based on meteorological data, altitude, land cover, and soils, was created for the surface-water-contributing basin of Whittlesey Creek by use of SWAT (Neitsch and others, 2001). SWAT uses a geographic information system (GIS)-based interface and GIS coverages of digital elevation models to divide a drainage basin into similar-sized subbasins that are further divided, on the basis of soil type and land cover, into areas of like hydrologic characteristics. Meteorological data were used as inputs if available; if measured values were not available, then the model's estimating routines for nonexistent or missing periods of data were used. The SWAT model incorporates many commonly applied, widely accepted empirical formulas in its calculations of continuous runoff and flow routing. The model also accounts for antecedent conditions by tracking soil-moisture changes with continuous and seasonally adjusted mathematical routines representing soil and plants. The calibrated rainfall-runoff model was used to simulate the effects of five land-cover simulations and their effects on annual runoff and flood peaks in Whittlesey Creek.

Field Investigations

Streamflow measurements were made in Whittlesey Creek and surrounding streams during low-flow (base-flow) conditions by use of standard USGS techniques (Buchanan and Somers, 1984). Base-flow measurements were made in August 1999, April 2000, and August 2000. Stream-water temperature was measured near the center of flow with a hand-held thermometer or water-quality multimeter. In August 2000, pH and specific conductance also were measured with a water-quality multimeter (fig. 1).

Streambed head and temperature data were collected during August 2000 from a subset of the base-flow sites. Streambed head was measured by installing a minipiezometer to a depth of 3 ft below the stream bottom near the center of flow in the stream. The minipiezometers were constructed of a high-density polyethylene tubing attached to a 6-in. stainless steel mesh screen. The head in the minipiezometer was allowed to equilibrate and was then measured relative to the stream surface. Streambed temperature data were collected by use of a hand-held thermometer inserted 6-12 in. below the stream bottom. Temperature measurements were made near the center of flow and in shaded locations where possible.

Estimates of streambed leakance (the vertical hydraulic conductivity of the streambed sediments divided by its thickness) were calculated for three sites on Whittlesey Creek and one site on North Fish Creek. Estimates of vertical hydraulic conductivity (K_v) were based on Darcy's Law (Freeze and Cherry, 1979) and measurements of gradient between the surface-water and ground-water systems (minipiezometer head data) and unit ground-water discharge to the creek. Unit discharge to the creek was determined from the measured change in discharge (August, 2000 base-flow measurements) over a known stream area. The stream area was determined from stream-width measurements made in the field and stream-length measurements between base-flow sites made with a GIS. A streambed thickness of 3 ft was used to calculate leakance and was based on the depth of minipiezometer penetration. Streambed leakance was used as ground-water-flow model input to quantify the interaction between ground-water and surface-water features of the model.

Surficial deposits were probed and core samples collected at one location underlain by the Copper Falls Formation and another location underlain by the Miller Creek Formation during July 2000 (fig. 1). The Copper Falls site was approximately 1 mi west of the headwaters of Whittlesey Creek in the Chequamegon National Forest. The Miller Creek site was approximately 800 ft south of Whittlesey Creek near the confluence with the North Fork. Cores were collected to verify the texture of surficial deposits. Probing was done with a wire-wound stainless steel screen to determine the depth to the water table and determine whether there were perched ground-water conditions.

Continuous stream stage was recorded at 15-minute intervals at the streamflow-gaging station (fig. 1). Data collection started April 1, 1999, and continues to the present (2003); however, several periods of record are missing from times when the equipment was vandalized before a more secure structure was built. Periods of missing record are April 1-21, 1999; May 30-June 2, 1999; October 31, 1999-February 29, 2000; March 9-10, 2000; and April 19-August 28, 2000. During periods of missing data, stage at the monitoring site was measured weekly by a local observer, and peak stage was recorded by use of a crest-stage recorder (Buchanan and Somers, 1984). Flow measurements were made at varying stream stages throughout the year by use of standard USGS methods (Buchanan and Somers, 1984). Relations between stage and flow were determined from these data, and a rating table was developed to convert stage data into flow data. Flow was calculated for each 15-minute interval (when data were available) and used to determine the daily mean streamflow at the station (Kennedy, 1983). During periods of missing stage data, estimated flows were based on daily measurements of stage at the station and by comparison of recorded flow in the nearby North Fish Creek. Streamflow statistics for the station (USGS station number 040263205) for the period April 1, 1999, to September 30, 2001, are listed in table 1.

Table 1. Streamflow statistics for U.S. Geological Survey streamflow-gaging station 040263205, Whittlesey Creek near Ashland, Wis., for April 1, 1999, to September 30, 2001

[ft³/s, cubic feet per second]

Remote stage recorders were placed in the basin upstream from the streamflow-gaging station and used to determine the timing of peak flows through the channel network (fig. 1). The remote recorders collected stage data at 1-hour or 15-minute increments. Data from these recorders were used to determine the time of concentration (the amount of time it takes a drop of water falling at the hydraulically farthest point from the drainage basin outlet to reach the outlet used for peak runoff calculation) and to help determine the channel routing variables to use during model calibration. Damage from vandalism and the debris carried by floods resulted in a significant amount of lost data from these stage recorders and required moving the stage recorders between floods to get sufficient data for the SWAT model.

Simulation of Ground-Water Flow

The GFLOW model domain includes all major drainage basins in the vicinity of the Bayfield Peninsula, ranging from the Bois Brule River in the west to the Bad River in the east (fig. 5). Other surface-water features include Lake Superior to the north and several lakes and streams to the south. The geometry of the model layer includes a stepped bottom altitude ranging from 500 ft near Lake Superior to 750 ft beneath the Bayfield Highlands. The single layer represents all of the unconsolidated deposits and the shallow part of the Bayfield Group. Zones of horizontal hydraulic conductivity (K_h) and recharge to the ground-water system are represented in the GLFOW model as inhomogeneities and represent the sandy deposits of the Copper Falls Formation and the clayey Miller Creek Formation. The GFLOW model includes "far-field" and "near-field" linesinks. The far field surrounds the near field area of interest. In the far field, streams and lakes are simulated with coarse linesinks having little or no resistance between the surface-water feature and the ground-water system. The purpose of simulating the far field is to have the model explicitly define the regional ground-water-flow field in the vicinity of the area of interest. The near field represents the area of interest and includes several of the streams adjacent to and including Whittlesey Creek (fig. 5). Near-field streams are simulated using slightly more detailed linesinks with streambed resistance to control ground-water/surface-water interaction. In GFLOW, resistance is calculated by dividing streambed sediment thickness by the K_v of the sediments. Near-field streams are used for flux calibration; far-field streams are not. The linesinks representing streams and lakes were assigned altitudes based on USGS 7.5-minute quadrangle maps. Near-field linesinks were assigned stream widths (ranging from 2 to 75 ft) based on field measurements and stream order.

Figure 5. Hydrologic features simulated with analytic elements for the GFLOW model of the Whittlesey Creek, Wis., study area and vicinity. (The locations of flux and head targets are shown in figure 7.)

Fifty-nine head targets and eight flux targets were used for the GFLOW model calibration. Head targets were compiled using depth to water information from drillers' logs. Depth to water was converted to water-level altitude based on site location and land-surface altitudes on USGS 7.5-minute quadrangle maps. Flux targets included four sites on Whittlesey Creek, three sites on North Fish Creek, and one site on an unnamed tributary north of Whittlesey Creek (fig. 6). Recharge and hydraulic-conductivity values were adjusted until simulated heads and base flows were similar to measured values.

The MODFLOW model (figs. 6 and 7) includes all of the delineated surface-water drainages of Whittlesey Creek and much of North Fish Creek and the Sioux River, and part of Lake Superior. The MODFLOW model consists of two layers with layer 1 representing unconsolidated deposits and layer 2 representing the Bayfield Group. The model grid contains 200 rows and 200 columns per layer (40,000 model cells per layer). Each model cell is 462 ft on a side. Model-layer top and bottom altitudes were based on various data sources. Layer 1 top altitudes were based on 500-m (1,640-ft) and 30-m (98-ft) Digital Elevation Model (DEM) data (U.S. Geological Survey, 2001). Layer 1 bottom altitudes were determined from geologic and drillers' logs obtained from the Wisconsin Geological and Natural History Survey and from published bedrock altitudes (Wold, 1979) and unconsolidated deposit thickness maps of the area (Young and Skinner, 1974). Because of the extreme thickness of the Bayfield Group (as much as 6,900 ft), it was assumed the altitude of the base of the ground-water-flow system (bottom of layer 2) was much shallower than the bottom of the sandstone. It is possible that the ground-water-flow system may be bounded at depth by a fresh-water/saline-water interface. Information about the possible location of this interface is unavailable; therefore, the bottom of the ground-water flow system was arbitrarily set at an altitude of zero. Sensitivity of the model simulations to changes in bottom altitude was tested. Lateral boundaries for this model are represented by constant flux (MODFLOW well package), which was extracted from the corresponding area of the GFLOW model following the method described by Hunt and others (1998). This extraction feature is particularly important for this study area because the extent of the surface-water noncontributing basin is difficult to delineate, and surface-water topographic divides may not coincide with ground-water divides. The extent of the MODFLOW model includes several major hydrologic boundaries for Whittlesey Creek and therefore should have minimal impact on simulation of ground-water flow in the vicinity of the creek. The major boundaries include Lake Superior to the east, Sioux River to the north, North Fish Creek to the south, and a regional ground-water divide in the Bayfield Highlands to the west. No high-capacity wells are within the MODFLOW model area; therefore, the only wells included in the model are those that represent the constant flux boundaries.

Figure 6. Boundary conditions and location of calibration targets used in the MODFLOW model of the Whittlesey Creek, Wis., study area.

Figure 7. Block diagram showing the MODFLOW model hydraulic conductivity/recharge zonation, river nodes, and layer configuration. (Area shown is the same as MODFLOW model extent in figure 6.)

The MODFLOW model was divided into three zones of hydraulic conductivity and two zones of recharge (fig. 7). Model layer 1 includes two K zones and corresponding recharge zones representing areas underlain by the Copper Falls and Miller Creek Formations. Layer 2 has a third K zone that represents the Bayfield Group. Initial estimates of many of the hydrogeologic parameters used in the model are listed in table 2. Gaining reaches of streams and Lake Superior were simulated using the MODFLOW river package (McDonald and Harbaugh, 1988). Streambed leakance for most of the simulated streams was based on field measurements. Streambed leakance for Lake Superior was initially estimated and then adjusted during model calibration.

Table 2. Estimated range or value of aquifer properties, recharge, and stream characteristics in the study area

[K_h, horizontal hydraulic conductivity; K_v, vertical hydraulic conductivity; Fm, formation; ft, feet; d, day; in/yr, inches per year]

The same head and flux targets used in the GFLOW model were also used for calibration of the MODFLOW model (fig. 6). Of the 59 head targets, 53 are in unconsolidated deposits (model layer 1) and 6 are in bedrock (model layer 2).

Estimates of Hydraulic Conductivity and Recharge

Estimates of horizontal hydraulic conductivity (K_h) for the Bayfield Group and sandy Copper Falls Formation were calculated by use of the computer program TGUESS (Bradbury and Rothschild, 1985), which uses data collected during pumping tests conducted when a well is drilled. Data were available for 53 wells completed in the Copper Falls Formation and 6 wells completed in the Bayfield Group. The average value of K_h for wells in the Copper Falls Formation was 38 ft/d (values ranged from 0.48 to 230 ft/d) and 1.4 ft/d (ranging from 0.13 to 4.6 ft/d) for in wells in the Bayfield Group. K_h values calculated from a nearby pumping test, conducted by the USGS on Precambrian bedrock similar to the Bayfield Group, range from about 2 to 10 ft/d (Charles Dunning, USGS, oral commun., April 2000). No pumping-test data were available for the Miller Creek Formation. The Miller Creek Formation and underlying deposits were typically described as alternating layers of sand and clay in drillers' logs, and for this reason, K_h in this part of the study area was assumed to be lower than K_h of the Copper Falls Formation (a K_h of 30 ft/d was used in the preliminary models). Hydraulic conductivities for unconsolidated deposits underlying Lake Superior were assigned values similar to those used for the Miller Creek Formation (table 2). Because of the relatively homogeneous nature of the Copper Falls Formation, vertical hydraulic conductivity (K_v) was set equal to K_h. Because the Bayfield Group includes layers of sandstone, siltstone, and shale and the Miller Creek Formation includes layers of sand and clay, K_v was assumed to be lower than K_h for the Bayfield Group and Miller Creek Formation. K_v values of 0.14 and 0.03 ft/d for the Bayfield Group and Miller Creek Formations, respectively, were used in the preliminary model.

Ground-water recharge in the study area was estimated to be between 8 and 18 in/yr in the areas directly underlain by the Copper Falls Formation (mostly sand) and 1 to 3 in/yr in areas underlain by the Miller Creek Formation (mostly clay) (table 2). Initial estimates of recharge were based on base-flow estimates from streamflow data from nearby stations on North Fish Creek, and the Bois Brule, White, and Bad Rivers. Recharge was estimated from these stations by means of the computer program HYSEP (Sloto and Crouse, 1996) and ranged from 8.2 to 17.4 in/yr. Most of these rivers have ground-water-contributing areas directly underlain by both the Copper Falls and Miller Creek Formations. Therefore, the resulting estimates of recharge are a combination of higher recharge in areas underlain by sand and lower recharge in areas underlain by clay. Since there are no surface-drainage features in the Bayfield Highlands, another estimate of recharge through the Copper Falls Formation can be calculated as average annual precipitation (32 in.) minus average annual evapotranspiration (17.7 in.). The resulting estimate of 14.3 in/yr of recharge assumes negligible underflow and change in storage. An additional streamflow site on North Fish Creek near the headwaters has base flow that represents about 2 to 3 in/yr of ground-water recharge and is in an area mostly underlain by clay (Rose and Graczyk, 1996). The underlying assumptions in estimating recharge from streamflow records are that base flow represents ground water discharging to the stream and that the surface-water-contributing basin and the ground-water-contributing area are coincident. Although the second assumption is not necessarily true in the study area, it does provide a starting point for estimating recharge, which can be reevaluated during model calibration.

Parameter Sensitivity and Model Calibration

The computer program UCODE (Poeter and Hill, 1998) was used for parameter sensitivity testing and automatic calibration of the MODFLOW model. Observations used for sensitivity testing and model calibration include the 59 head targets and 8 flux targets. Observation weights were assigned as follows: all heads, standard deviation = 5 ft; all streamflow stations or sites measured during low-flow conditions, coefficient of variation = 0.1; and sites where streamflow was estimated from staff gages, coefficient of variation = 0.15. Parameters tested for sensitivity were: K_h and K_v of the Copper Falls and Miller Creek Formations and the Bayfield Group, recharge (two zones), Lake Superior streambed leakance, and the altitude of the bottom of layer 2. Parameters used in model calibration were: K_h of the Miller Creek and Copper Falls Formations, K_v of the Miller Creek Formation and the Bayfield Group, and streambed leakance for Lake Superior.

Pathline Analysis and Delineation of Ground-Water-Contributing Area

The ground-water-contributing area (from the deep flow system and nonperched part of the shallow system) for Whittlesey Creek was delineated using the computer program MODPATH (Pollock, 1994) and output from the calibrated MODFLOW model. Ten hypothetical particles of water were placed in each river cell that represents the gaining part of Whittlesey Creek. The particles were then tracked backward to the water table, and the corresponding pathlines were exported into a GIS coverage. The areal extent of the pathlines was manually digitized and represents the ground-water-contributing area. Values of porosity used in this analysis were 14 percent for the Bayfield Group (Ojakangas, 1986), 25 percent for the Copper Falls Formation, and 35 percent for the Miller Creek Formation. The porosities used for the unconsolidated deposits were based on ranges of values published by Freeze and Cherry (1979).

Simulation of Rainfall Runoff

For the SWAT model, the Whittlesey Creek surface-water-contributing basin was divided into 19 subbasins based on tributary junctions and similarities in soils, land cover, and slopes (table 3). Each subbasin was further divided into Hydrologic Response Units (HRUs) based on the significant land cover and soil type. In this study, land cover was considered significant if it made up 10 percent or more of the subbasin area and a soil type significant if it made up 15 percent or more of the subbasin area. An HRU was created for every combination of soil and land cover greater than those thresholds. Areas of land cover and soil type less than those thresholds were considered insignificant and were divided among the significant HRUs based on the ratio of HRU area to total area. Thus, if 60 percent of the basin was HRU1, then 60 percent of the insignificant area was treated the same as HRU1.

Table 3. Physical characteristics of subbasins used in Whittlesey Creek, Wis., Soil and Water Assesment Tool (SWAT) model

[mi, miles; mi², square miles; ft, feet; in/hr, inches per hour]

The Whittlesey Creek SWAT model used the Soil Conservation Service (SCS) Curve Number Method for runoff estimation. (Soil Conservation Service, 1972). This method uses soil and land-cover information to determine curve numbers, which in turn are used by the model to predict annual runoff. Time of concentration and the peak rain intensity are used to calculate the peak flood flows. SCS curve numbers used in the Whittlesey Creek SWAT model for the five land-cover simulations are listed in table 4. The SWAT model adjusts the SCS curve numbers (and thus the amount of runoff) between runoff events on the basis of antecedent soil moisture conditions. SWAT uses Manning's equation and channel roughness coefficients to define the rate and velocity of flow (Chow, 1959). Water is routed through the channel network by use of the Muskingum Method, which is a modification of the kinematic wave flood-routing model described by Chow and others (1988). The SWAT model uses a procedure developed by Lane (1983) for estimating transmission losses; that is, the loss of water through the bottom and sides of the channel.

Table 4. Soil and Water Assesment Tool (SWAT) model input parameter values for Whittlesey Creek, Wis., for five model simulations

[SCS, Soil Conservation Service (1972); °C, degrees Celsius.]

¹Represent a 25-percent increase in urban residential land having curve numbers of 68 for sand and 85 for clay.

Evapotranspiration was calculated in the model by the Priestly-Taylor Method (Priestly and Taylor, 1972) using solar radiation, air temperature, and relative humidity from the National Oceanic and Atmospheric Administration station in Gordon, Wis., about 50 mi southwest of the Whittlesey Creek Basin. Snowmelt was modeled by use of snowfall and temperature records and was based on equations developed by Anderson (1976). Hourly rainfall data at the streamflow-gaging station were limited because of repeated vandalizing of equipment at the site (no data from rain events greater than 1 in. were available); thus, precipitation data available from the Ashland Experimental Farm Weather Station of the National Weather Service were used in the model. Only daily data were available at the Ashland Experimental Farm; thus, a skewed normal distribution of daily rainfall was used in the model (Fiering, 1967, Nicks, 1974). Daily temperature means, maximums, and minimums used in the model also came from the Ashland Experimental Farm Weather Station.

DEM data were produced by use of softcopy photogrammetric software solutions as part of the Whittlesey Creek Orthophoto Production Project (fig. 8) (Benchmark GIS, 2001). GIS programs were used to develop high-resolution othophotography and spatially accurate terrain models. The DEM has 10-meter cell spacing and was derived by interpolation of vector contour lines from 1:24,000-scale USGS topographic maps (Benchmark GIS, 2001). The surface-water-contributing basin derived from DEM data and used in the SWAT model had an area of 7.2 mi², which is slightly smaller than the basin derived from USGS 7.5-minute topographic maps (7.4 mi²).

Figure 8. Digital Elevation Model of land surface for the Whittlesey Creek, Wis., surface-water-contributing basin. Drainage divides were derived from 10-meter Digital Elevation Model data (Benchmark GIS, 2001).

The WISCLAND (Wisconsin Initiative for Statewide Cooperation on Landscape Analysis and Data) land-cover data used in the model were derived from LANDSAT Thematic Mapper (TM) satellite imagery acquired from flyovers in August 1991; May, July, September, and October 1992; and May 1993 (Reese and others, 2002).

Copies of township-based maps of land cover in 1928 were obtained from the Wisconsin Historical Society Library (Wisconsin Department of Agriculture and Wisconsin Geological and Natural History Survey, 1928). The maps were digitized and spatially rectified by use of a GIS.

Soil-characteristic information for the study area was compiled from USSOILS, a national digital coverage of the State Soil Geographic Database (STATSGO) (Schwarz and Alexander, 1995). Soil characteristics examined in this study included clay content (percentage of soil with particles less than 2 µm in size), organic-matter content (percentage by weight), soil erodibility (K factor), permeability rates (in/hr), and slope (percent). The USSOILS coverage was intersected with the drainage subbasin outlines (fig. 9), and each of the soil characteristics of interest was calculated as an area-weighted average. The scale of the soils coverage is coarse (1:250,000) given the size of the Whittlesey Creek surface-water-contributing area, but it was the best source of digital soils data when the SWAT model was run.

Figure 9. Soils in the Whittlesey Creek, Wis., surface-water-contributing basin from STATSGO data (Schwartz and Alexander, 1995). Drainage divides were derived from 10-meter Digital Elevation Model data (Benchmark GIS, 2001).

SWAT Model Calibration

Whittlesey Creek had a consistent measured base flow of 17 to 18 ft³/s at the USGS streamflow-gaging station during the study. The consistency of this base flow and the results of the ground-water-flow modeling done as part of this study indicate that most of this base flow originates outside the surface-water-contributing basin for Whittlesey Creek, mostly from ground water in the deep regional flow system. The amount of base flow originating inside the surface-water-contributing basin is minimal. Without a regional source of base flow, all of Whittlesey Creek would likely be ephemeral.

The SWAT model is able to track water that falls within the surface-water-contributing basin. The model does not have the ability to account for water that crosses surface-water drainage divides, such as regional ground water. Because of this, the model was calibrated to the flow at the streamflow-gaging station (field number 2, fig. 1). This flow was first adjusted for the regional ground-water component by subtracting a constant 16.75 ft³/s from the station flow record. The value of 16.75 ft³/s was used because it was the lowest base flow recorded at the station during water year 1999 and the uncalibrated SWAT model indicated no base flow originated from the surface-water-contributing basin during extended periods of minimal precipitation. The lowest flow ever recorded at the streamflow-gaging station is 16 ft³/s, whereas the ground-water-flow model estimated regional base flow at the station to be 15 ft³/s.

In addition to data from the streamflow-gaging station, data for the timing of peak stage in the channel upstream of the station were gathered from continuous stage-recorder sites (fig. 1) and used to calculate the time of concentration for model calibration. Hourly stage data showing the shape and timing of a peak from the main stem of Whittlesey Creek above the confluence, the North Fork of Whittlesey Creek, and at the USGS streamflow-gaging station are shown in figure 10 for a runoff event that occurred on May 9-10, 2001.

Figure 10. Measured hourly stage data used in timing peak stages through Whittlesey Creek, Wis., for calibration of channel characteristics in the Soil and Water Assessment Tool (SWAT) model, May 9-10, 2001.

The period of record used in the SWAT model calibration was April 1, 1999, to October 15, 2001. The calibration targets were annual runoff for the period of calibration, the peak daily mean flows for three of the largest floods during the period of calibration, peak timing, base flow, snowmelt timing, and the regression rate of both peaks and base flow. Modeled annual runoff and flood peaks were calibrated by adjusting water-budget variables such as SCS curve numbers, evapotranspiration rates, overland-flow rates, and ground-water interaction parameters. It was determined during calibration that the USSOILS coverage, which lumps similar soil types in groups, was not adequate by itself for this scale model. Therefore, soil-profile parameters were adjusted for several subbasins on the basis of county-level soil survey maps (unpublished data, Ulf Gafvert, Natural Resources Conservation Service, written comm., 2002) that more accurately detail the complexity in the soils within the Whittlesey Creek Basin. Peak timing was calibrated by adjusting channel-routing characteristics and snowmelt parameters until peak timing matched the values calculated from actual stage data.

Use of a continuous model like SWAT can help account for varying antecedent conditions such as moisture content, seasonal variations in land cover, and even factors as detailed as the amount of detritus on the land surface. The inherent problems in modeling remain, however, in that empirical formulas incorporate a finite number of variables to represent a complex real system. Because of this, the model was calibrated with a best-fit approach. During model calibration, input parameters were adjusted until the modeled results best fit all calibration targets, though emphasis was placed on the average annual runoff for the period of calibration, the size and timing of the three largest peak flows, and base-flow recession rates.

Rainfall-Runoff Model Simulations

Using the calibrated SWAT model, conditions for 1992-93 land cover and four alternative land-covers were simulated (table 5). The first simulation (presettlement) represented conditions before European settlement (pre-1870), and the second (peak agriculture) represents historical conditions in the 1920s and 1930s, when row-crop agriculture in the area was at its peak. Both simulations use land-cover data from the 1928 Land Economic Inventory (Wisconsin Department of Agriculture and Wisconsin Geological and Natural History Survey, 1928) (fig. 11). The peak agriculture scenario uses the data as is, whereas all land cover other than wetland is converted to forest in the presettlement scenario because it is assumed that the forest was clear cut from about 1870 to 1900.

Table 5. Percentages of land-cover types used in the Whittlesey Creek, Wis., Soil and Water Assessment Tool (SWAT) model simulations

[All data are percentages. 1992-93 land cover from WISCLAND land cover (Reese and others, 2002). Peak agriculture data from 1928 Land Economic Inventory (Wisconsin Department of Agriculture and Wisconsin Geological and Natural History Survey, 1928)]

Figure 11. Land cover in the Whittlesey Creek, Wis., surface-water-contributing basin, from the 1928 Land Economic Survey (Wisconsin Department of Agriculture and Wisconsin Geologic and Natural History Survey, 1928). Drainage divides were derived from 10-meter Digital Elevation Model data (Benchmark GIS, 2001).

All land-cover categories for 1992-93 (other than wetland) were changed to forest in the third model simulation (reforested) (table 5). Additionally, overland and channel roughness coefficients in the model were increased to the maximum published values for forest to simulate the greatest amount of debris in the forest floors in the upland and the greatest amount of woody debris in the stream channel itself that would be reasonable (fig. 12). It should be noted that the reforested simulation has about 10 percent less wetland than the presettlement simulation.

Figure 12. Land cover in the Whittlesey Creek, Wis., surface-water-contributing basin used in the modeled reforested simulation. Drainage divides were derived from 10-meter Digital Elevation Model data (Benchmark GIS, 2001).

A fourth simulation (developed) was made by adjusting all curve numbers for each HRU in the calibrated model to represent a new mix of land cover, including 25 percent more urban residential land in each subbasin (table 5). Curve numbers from the 1992-93 land-cover categories HRU's were multiplied by 0.75; representing that 75 percent of the HRU stayed the same. The theoretical curve numbers for HRU's with urban-residential land use (1/4- to 1/8-acre lots and the same soil) were determined and multiplied by 0.25; representing development of 25 percent of the basin. These curve numbers were then added together to get the new HRU's curve numbers that represent 75 percent of the 1992-93 land cover and 25 percent developed land cover.

Field Data

Results of field measurements of streamflow, water temperature, pH, and specific conductance are shown in table 6. Base-flow data indicated ephemeral and/or losing reaches in the upper parts of the main stem and North Fork of Whittlesey Creek, where flows were typically zero to less than 0.10 ft³/s (fig. 13). Base flows increased substantially near the confluence of the main stem and North Fork. Specific conductance was highest and water temperatures lowest near the confluence, indicating that regional ground water was discharging to the stream near the confluence. Ratios of base flow to drainage areas are given in table 7 for Whittlesey Creek and several adjacent streams. The ratios in Whittlesey Creek increased significantly below the confluence and were greater than those in adjacent streams, including North Fish Creek, Sioux River, Boyd Creek, and a Lake Superior tributary to the south of Whittlesey Creek. Base flows were generally slightly higher in April than in August.

Table 6. Results of field investigations of base flow and physical characteristics of Whittlesey Creek and adjacent streams, Wis.

[USGS, U.S. Geological Survey; mi², square miles; ft³/s, cubic feet per second; µS/cm, microsiemens per centimeter at 25 °C; °C, degrees Celsius; in., inches; ns, not sampled; —, no USGS station number established; e, estimated. Drainage area reflects surface-water-contributing basin only.]

Figure 13. Base flows for Whittlesey Creek and tributaries, Wis., August 2000.

Table 7. Comparison of base flow and drainage area for Whittlesey Creek and adjacent streams, Wis.

[USGS, U.S. Geological Survey; e, estimated; mi², square miles; ft³/s, cubic feet per second; –, no USGS station number established. Average discharge reflects average of August 1999 and August 2000 discharge measurements. Drainage area reflects surface-water-contributing basin only.]

Streambed-head measurements ranged from –0.35 to 6.9 in. in Whittlesey Creek (fig. 14). Negative heads, mostly found near the headwaters, indicate losing reaches, whereas positive heads indicate gaining reaches. The highest measured heads were near the confluence. Streambed temperatures ranged from 7.5 to 19.5°C (fig. 15). Similar to water temperature, the highest streambed temperatures were measured in the headwater reaches, and the lowest were measured near the confluence. Slightly increasing streambed temperatures downstream from the confluence may indicate a larger proportion of contributions from warmer water of the shallow system or seasonal warming of the streambed sediments where the hydraulic gradient is more horizontal. Calculated streambed-leakance values, based on stream-channel area, discharge, and streambed head, ranged from 45 to 96 ft/d/ft for Whittlesey and North Fish Creeks. A leakance value of 45 ft/d/ft was used for the Sioux River. Leakance values for Boyd Creek and the other small tributaries in the area were not measured, but were assumed to be much lower than those for deeply incised, steep-gradient streams like Whittlesey Creek. A streambed leakance of 0.1 ft/d/ft for the unnamed tributary north of Whittlesey Creek resulted in simulated base flows that were close to measured values. For this reason, 0.1 ft/d/ft was used as the leakance value for Boyd Creek and the other small tributaries. Leakance values for Lake Superior sediments were assumed to be even lower than those of the small creeks and tributaries in the area because of greater thickness and lower permeability in the relatively undisturbed offshore depositional environment. A leakance value of 0.01 ft/d/ft was used for Lake Superior.

Figure 14. Streambed heads for Whittlesey Creek and tributaries, Wis., August 2000.

Figure 15. Streambed temperatures for Whittlesey Creek and tributaries, Wis., August 2000.

Core samples were collected to identify the soil types and determine the depth to the water table at two sites in the basin (fig. 1). The Copper Falls site (coring site 1) was probed to a depth of 46 ft, and core samples were collected from land surface to a depth of 12 ft. The surficial deposits associated with the Copper Falls Formation were easy to probe and consisted of loose, reddish-brown, medium-grained sand. The water table was deeper than 46 ft below land surface at this location. The Miller Creek site (coring site 2) was probed to a depth of 21 ft. Core samples were collected from land surface to a depth of 16 ft. The surficial deposits, consisting of layers of red sandy clay, red sand, and red clay, were very tight at this location, making probing difficult. At 21 ft below land surface, the surficial deposits were unsaturated; however, a perched water table was present at 11 ft below land surface.

Ground-Water-Flow Simulations

As a result of the stepwise modeling approach used in this study, several iterations of ground-water-flow models were developed for the study area. Each iteration included a GFLOW model simulation and a MODFLOW model simulation. Initially, preliminary models were constructed from available data and were used to help understand the ground-water-flow system, test hypotheses, and guide field-data collection. The second set of models incorporated updates to the preliminary models and included results of data and insights from field observations. The updated MODFLOW model was then calibrated by use of the computer program UCODE. The hydraulic conductivities from the calibrated MODFLOW model were incorporated into a final GFLOW model. The calibrated MODFLOW model was used to delineate the ground-water-contributing area to Whittlesey Creek.

A summary of the model results is given below. Optimized model parameters are listed in table 8. For each simulation, head statistics and simulated base flows at selected locations are listed in table 9.

Table 8. Preliminary estimate and UCODE-optimized parameter values used in the MODFLOW model

[K_h, horizontal hydraulic conductivity; K_v, vertical hydraulic conductivity; ft, feet; d, day; in/yr, inches per year]

¹Bayfield Group K_h was highly correlated with Miller Creek K_h; therefore, both could not be estimated effectively with UCODE. More data were available for estimating Bayfield Group K_h, so this parameter was excluded from optimization.

²Parameter not optimized because of model insensitivity. For Copper Falls Formation, K_v was set equal to optimized K_h.

³Parameter not optimized. Preliminary estimate from GFLOW model.

Table 9. Head statistics and simulated base flows at selected locations for preliminary, updated, and calibrated ground-water-flow models

[ft, feet; ft³/s, cubic feet per second]

¹Measured base flow target for Whittlesey Creek (field number 1) is 17 ft³/s.

²Measured base flow target for North Fish Creek (field number 27) is 86 ft³/s.

Preliminary Ground-Water-Flow Models

GFLOW Model

The preliminary GFLOW model includes the Bayfield Peninsula and vicinity (fig. 5) and was used as a screening tool to test hypotheses about the ground-water-flow system in the study area, estimate recharge, and to help improve the subsequent MODFLOW model. One question addressed with the GFLOW model was whether the expected distribution and range of recharge values could simulate sufficient base flow in Whittlesey Creek. Hydraulic conductivities of 30 and 100 ft/d and recharge rates of 2 and 17 in/yr (for the area underlain by the Miller Creek Formations and Copper Falls Formation, respectively) gave reasonable preliminary results and were within the expected range for those parameters. These values produced slightly low base flows in Whittlesey Creek and slightly high base flows for North Fish Creek (table 9). The absolute residual mean head error for this model was 22.3 ft.

MODFLOW Model

The preliminary MODFLOW model includes all of the Whittlesey and North Fish Creek drainages and much of the Sioux River drainage. The lateral boundaries of this model are represented by constant flux cells (simulated using the MODFLOW well package) that were extracted from the corresponding area of the preliminary GFLOW model. This extraction feature is particularly important for this study area because the extent of the surface-water drainage is difficult to delineate and, on the basis of GFLOW model, does not coincide with the ground-water divides. Horizontal hydraulic conductivities (K_h) of 100, 30 and 1.4 ft/d, for the Copper Falls and Miller Creek Formations and the Bayfield Group, respectively, and recharge values of 17 and 2 in/yr, for the area underlain by the Copper Falls and Miller Creek Formations, respectively, gave reasonable preliminary results (table 9). Again, the simulated base flows were slightly low for Whittlesey Creek and slightly high for North Fish Creek. The absolute residual mean head error for this model was 21.0 ft.

Modeled heads and fluxes for the preliminary GFLOW and MODFLOW models were reasonably close to measured values, with slightly better modeled head results from the MODFLOW model. Results of both ground-water-flow models indicate that much of the headwaters of Whittlesey Creek receive no ground-water contributions from the deep system. Both models also indicated that the deep ground-water-flow system intersected Whittlesey Creek just upstream from the confluence of the main stem and North Fork. In addition, both models indicated that the ground-water-contributing area to the stream was not coincident with the delineated surface-water-contributing basin. Because the preliminary models gave reasonable results, subsequent field investigations focused mainly on verifying preliminary model results and collecting additional base-flow data. In addition, field determination of streambed leakance was expected to improve modeled base-flow values.

Field measurements of base flow, streambed head, and streambed and water temperature indicated that the deep flow system intersected Whittlesey Creek near the confluence of the main stem and North Fork (table 6). Small streamflows, negative minipiezometer heads, and warmer streambed temperatures in headwater reaches indicated that those areas were receiving only a small amount of ground water from perched parts of the shallow system upstream from the confluence. Losing reaches and negative heads also indicated that some streamflow, derived from the perched part of the shallow system, recharged the deep system by percolating through the stream bottom.

Updated Models, Sensitivity, and Calibration

The GFLOW and MODFLOW models were updated with field-measured streambed leakance. In addition, the distribution of recharge in both models was altered to improve simulated base flows in North Fish Creek. The distribution of the high-recharge zone corresponding to the area underlain by the Copper Falls Formation was limited to the thick, permeable deposits of the Bayfield Highlands. Areas underlain by relatively thin deposits of the Copper Falls Formation, mainly south and west of the upper reaches of North Fish Creek, were assigned recharge rates similar to those used for areas underlain by the Miller Creek Formation. The updated models resulted in improved simulated heads and base flows for North Fish Creek. Simulated base flows for Whittlesey Creek were relatively unchanged. A range of recharge values was simulated in the GFLOW model—from 10 to 21 in/yr in the area underlain by the Copper Falls Formation and from 1 to 3 in/yr in the area underlain by the Miller Creek Formation. Values of 17 and 2 in/yr for the Copper Falls and Miller Creek areas, respectively, provided the best overall fit for heads and base flows simulated using GFLOW. In the MODFLOW model, the GFLOW-extracted boundary flux was based on the recharge values used in the GFLOW model. Because the boundary flux was fixed to the values used in the GFLOW model, aerial recharge in the MODFLOW model could not be properly estimated using UCODE. As a result, recharge values from the updated GFLOW model were used in the MODFLOW model instead of being estimated with UCODE.

Input parameters for the MODFLOW model were evaluated for model sensitivity by means of UCODE (fig. 16). The model is most sensitive to recharge and K_h of the Miller Creek and Copper Falls Formations. The model is least sensitive to K_v of the Copper Falls Formation and the bottom altitude of layer 2. The least sensitive parameters were excluded from UCODE estimation.

Figure 16. MODFLOW ground-water-flow model parameter sensitivities from UCODE for the Whittlesey Creek, Wis., study area. [Red asterisk (*) indicates parameters estimated with UCODE. K_h; horizontal hydraulic conductivity, K_v; vertical hydraulic conductivity.]

The MODFLOW model was calibrated using UCODE. The preliminary estimates and optimized values for model parameters are listed in table 8. The model parameters representing K_h of the Bayfield Group and Miller Creek Formations were highly correlated (r² = -0.99). Estimating highly correlated parameters with UCODE can result in nonunique results (Hill, 1998, p. 39); therefore, one of the highly correlated parameters should be excluded from estimation. Information on the K_h of the Bayfield Group was available, and the model was much less sensitive to this parameter; therefore, it was excluded from UCODE estimations. Most of the UCODE-estimated parameter values fell within the expected range or close to the estimated values listed in table 2. Using the parameter values from UCODE, the absolute residual mean head error for the calibrated MODFLOW model was 17.0 ft (table 9). Simulated base flows were improved, though still slightly low for Whittlesey Creek and slightly high for North Fish Creek. The relation of measured and modeled heads and stream flux for the calibrated MODFLOW model is shown in figure 17. Optimized parameter values for hydraulic conductivity from the calibrated MODFLOW model were then added to a final GFLOW model. The simulated base flows from the final updated GFLOW model were nearly unchanged; however, simulated heads improved slightly (table 9).

Figure 17. Relation of measured and modeled (A) head and (B) stream flux (base flow) for the calibrated MODFLOW model for the Whittlesey Creek, Wis., study area.

Simulated base flows in Whittlesey Creek were still lower than measured flows for several possible reasons. First, measured base flows included flow from the perched part of the shallow flow system. Simulations included only flow from the nonperched system. Second, measured base flows may have been slightly lower than 17 ft³/s under drier conditions, resulting in values closer to those modeled. Finally, the model is only an approximation of the real system. The model may require more detail and complexity (better geologic information, more layers, more recharge zones, more zones of hydraulic conductivity, accounting of perched flow, etc,) to provide a more accurate simulation of the system. With additional data, the model could be refined to address future issues or areas of interest.

Ground-Water-Contributing Area

Results from the final calibrated MODFLOW model were used to delineate the ground-water-contributing area to Whittlesey Creek from the deep flow system and the nonperched part of the shallow flow system. Ground-water pathlines were simulated with the computer program MODPATH. In all, 480 particles of water were tracked backwards, from the parts of Whittlesey Creek that received ground water from the deep flow system, to the water table. The 22-mi² land-surface area projected from the volume of aquifer encompassing those pathlines is the contributing area shown in figure 18. The average traveltime from the stream back to the water table was about 215 years; the median traveltime was about 94 years. Particle-tracking results indicated that the ground-water-contributing area did not coincide with the delineated surface-water-contributing basin. The Miller Creek Formation underlies approximately half of the ground-water-contributing area, and half is underlain by the Copper Falls Formation (fig. 19). Based on modeled recharge rates, about 90 percent of the base flow to Whittlesey Creek originated as recharge through the permeable deposits of the Copper Falls Formation in the Bayfield Highlands. On the basis of these results, concerns about land-cover effects on, and changes to, base flow in Whittlesey Creek should focus on this heavily forested part of the ground-water-contributing area.

Figure 18. Whittlesey Creek, Wis., ground-water-contributing area based on results from the calibrated MODFLOW model.

Figure 19. Particle pathlines tracked backwards from Whittlesey Creek, Wis., to the water table in the ground-water- contributing area.

Land-Cover Effects on Recharge and Base Flow

The relations among land-cover changes, recharge rates, and base flow in the vicinity of Whittlesey Creek are unclear. Changes in land cover could have an effect on ground-water recharge rates. Additionally, changing the recharge rate in an area can result in changes in base flow of nearby streams. Logging might increase ground-water recharge by reducing interception and evapotranspiration. Development on forested or agricultural land could reduce ground-water recharge by creating large impermeable areas (roofs, parking lots, roads, etc.). Tiling of agricultural land could reduce recharge by capturing and diverting water moving through the unsaturated zone. Land-cover changes and recharge effects on base flow in Whittlesey Creek were simulated using the calibrated MODFLOW model. Land-cover changes were represented by increasing recharge by 25 and 50 percent (a deforestation simulation) and decreasing it by 25 and 50 percent (an urbanization or tiling simulation). It is unknown exactly how much recharge in the vicinity of Whittlesey Creek would be affected by land-cover changes; however, the rates listed above should represent extreme cases.

Model-wide changes in recharge caused a proportional change in simulated base flow for Whittlesey Creek. Changes in base flows resulting from changes in recharge in only the ground-water-contributing area of Whittlesey Creek also were simulated (table 10). These localized changes in recharge were simulated for the entire ground-water-contributing area or for individual recharge zones (clayey zone or sandy zone) in the contributing area. Simulated changes in base flow are for Whittlesey Creek near USGS Station 040263205 (fig. 1).

Changing total recharge in the entire ground-water-contributing area by 25 percent resulted in an increase or decrease of about 6 percent in base flow (table 10). A 50-percent increase or decrease in recharge in the ground-water-contributing area resulted in about a 10 to 11 percent increase or decrease in base flow. Base flow did not change at the same rate as recharge because changes in recharge were applied only to the ground-water-contributing area of Whittlesey Creek. One outcome of the localized changes in recharge was that the size of the resulting contributing area to Whittlesey Creek changed. For example, increasing recharge by 50 percent in the original contributing area caused some mounding of the potentiometric surface there; the mounding changed the direction of ground-water flow and actually reduced the size of the resulting contributing area. Even though more recharge was applied, some of that recharge was lost to adjacent streams. Reducing recharge in the contributing area caused a depression in the potentiometric surface that resulted in a larger contributing area; when recharge was reduced, the larger contributing area made up for some of the decrease in recharge. If recharge was increased or decreased throughout the model area (not just the ground-water-contributing area), then localized mounding or depression of the potentiometric surface did not occur and the size of the contributing area hardly changed; the resulting percentage changes in base flow to Whittlesey Creek were, in this case, closely correlated with the percentage change in recharge.

Table 10. Simulated changes in base flow for Whittlesey Creek, Wis., as a result of changes in recharge in the ground- water-contributing area

[ft³/s, cubic feet per second; in/hr, inches per hour; %, percent]

¹ Simulated base flow for Whittlesey Creek at State Hwy. 13 is 14.8 ft³/s using the calibrated MODFLOW model

² Recharge rates used in the calibrated MODFLOW model are 17 and 2 in/hr for the sandy and clayey zones, respectively.

The two recharge zones within the ground-water-contributing area to Whittlesey Creek were about the same size; however, on the basis of the calibrated model, about 90 percent of base flow originated as recharge through the sandy Copper Falls Formation. Only about 10 percent of base flow was from recharge through the clayey Miller Creek Formation. Changing recharge in a single zone (sand or clay) would represent land-cover changes in part of the contributing area. A 25-percent change in recharge in the clayey zone resulted in about a 2-percent change in base flow (table 10). A 50-percent change in recharge in this zone resulted in about a 3-percent change in base flow. A likely land-cover change would be logging of forests in the sandy zone and a increase in recharge of possibly 25 percent, which would cause an increase in base flow of about 4 percent. A 50-percent change in recharge in the sandy zone would result in about a 7- to 8-percent change in base flow.

Rainfall-Runoff Simulations

Hydrographs of measured and simulated daily mean flows at the Whittlesey Creek streamflow station (field number 2) for the period of calibration are shown in figures 20 and 21. Daily mean flows were used in the simulation because of the lack of hourly rain data and because the SWAT model used was unable to calculate runoff in hourly time steps. Instantaneous peaks at the streamflow station from recorded 15-minute data were more than twice the simulated daily mean flows. Hydrographs of simulated daily mean flows for 1992-93 land cover, measured daily flow, and 15-minute data for the three largest peak flows during the period of model calibration are shown in figure 21. Daily mean and instantaneous peak flows and the calibrated model output are listed in table 11. Instantaneous peak flows were twice as high as daily peak flows during floods. Though modeling with a daily time step is adequate to determine average annual runoff and differentiate between surface runoff, evapotranspiration, and infiltration, it is not adequate to determine the absolute change in instantaneous peak flow. An hourly time step would improve the model results, especially for peak flows.

Figure 20. Hydrographs of daily mean flow from the Whittlesey Creek, Wis., streamflow-gaging station 040263205 and Soil and Water Assessment Tool (SWAT) model for 1992-93 land cover.

Figure 21. Flood hydrographs of measured daily mean and instantaneous (15-minute) flow from the Whittlesey Creek, Wis., streamflow-gaging station 040263205 and simulated daily mean flow from the Soil and Water Asessment Tool (SWAT) model for 1992-93 land cover, July 1999 and April 2001.

Table 11. Daily mean and intantaneous peak flows for U.S. Geological Survey streamflow-gaging station on Whittlesey Creek, Wis. and the calibrated SWAT model output using 1992-93 WISCLAND land-cover data (Reese and others, 2002)

[ft³/s; cubic feet per second]

For the modeled period, the average annual runoff from the surface-water- contributing area (excluding regional base flow) for the streamflow-gaging station was 10.9 in/yr (table 12). The SWAT model simulated the same value for 1992-93 conditions. The calibrated model underestimated the first large rainfall-derived daily mean flow for the flood peak on July 5, 1999, but it overestimated the daily mean flow for the peak that followed four days later (fig. 21). These runoff events were caused by thunderstorms, likely with rainfall intensity and amounts varying throughout the basin. The designed storm used in the simulation may not accurately represent this type of rain event. Hourly precipitation data from locations within the basin would have improved the model. The simulated daily mean flow for the April 23, 2001, rainfall- and snowmelt-derived peak closely matched actual values (fig. 21). The event on April 23, 2001, was the largest flood on record (3 years) at the streamflow-gaging station on Whittlesey Creek and data from nearby long-term streamflow-gaging stations indicate that it was an 80- to 100-year event.

Table 12. Average annual runoff and flood peaks (based on daily mean flows on peak-flow days) from the Whittlesey Creek, Wis., streamflow-gaging station 040263205 and the Soil and Water Assessment Tool (SWAT) model simulations

[in/yr, inches per year; ft³/s, cubic feet per second. Average annual runoff (April 4, 1999, through October 15, 2001) does not include the regional ground-water flow of 16.75 ft³/s and reflects the modeled surface-water-contributing basin of 7.2 square miles.]

Land-Cover Effects on Annual Runoff and Flood Peaks

Land-cover effects on average annual runoff in various model simulations for 1999 to 2001 were between 1 and 11 percent (table 12). The largest difference was between presettlement or reforested and peak agriculture. Flood peaks (based on daily mean flows on peak-flow days) increased in the developed and peak-agriculture simulations and decreased in the presettlement and reforested simulations (table 12). The increase in forest in the presettlement and reforested simulations resulted in smaller and broader flood peaks, whereas the developed and peak-agriculture simulations produced larger and steeper flood peaks compared to those peaks simulated with 1992-93 land use (fig. 22).

On the basis of data from North Fish Creek, an adjacent basin to the south of Whittlesey Creek, sediment loads were related to flow by a power function, indicating that decreasing peak flow, even without decreasing total runoff, can reduce erosion and resulting sediment loads (Fitzpatrick and others, 1999). It is assumed that Whittlesey Creek, which drains similar geologic deposits, would have a similar relation between sediment load and flow, and that small decreases in flood peaks could result in decreases in sediment loading.

Effects of the various land-cover simulations on flood peaks were affected by antecedent soil moisture; effects were less for subbasins with high runoff rates under nearly all conditions, such as urban development, intensive agriculture, and clay soils. Increased soil moisture has little effect on a land cover that can only absorb 10 percent of rainfall when dry, while a land cover that normally absorbs 50 percent of rainfall when dry will have much larger increases in runoff when saturated.

The storm of July 5, 1999, that occurred during dry conditions resulted in a greater difference in peaks between simulations than the storms of July 9, 1999, and April 23, 2001, which fell on nearly saturated soils (fig. 22). Simulations of developed and agricultural land cover resulted in larger increases in flood peaks during dry antecedent conditions (12 and 18 percent) than during wet conditions (0 and 7 percent). The entire basin, regardless of land cover, had nearly saturated soils during the storm of July 9, 1999, so the entire basin had high runoff for that storm in the calibrated model regardless of land cover. Differences in land cover between the simulations thus had a limited effect on the modeled results for that storm.

Figure 22. Flood hydrographs of daily mean flow from the Whittlesey Creek, Wis., Soil and Water Assessment Tool (SWAT) models of all land-cover simulations.

Flood peaks from the developed land-cover simulation were slightly smaller than those from the peak-agriculture simulation (fig. 22) because of the conversion of agriculture from predominately row crops in 1928 to predominately hay and grasslands in the developed simulation. SCS curve numbers for 25-percent developed grassland and forest are lower than those for row crops (table 4) and the area of urban land (25 percent) in the developed simulation is slightly smaller than the area of row crop agriculture (28.2 percent) in the peak-agriculture simulation (table 5).

Increases in forest in the reforested and presettlement simulations resulted in reduced daily mean flows on peak-flow days (12-14 percent and 1-4 percent, respectively) compared to 1992-93. It was assumed that the presettlement Whittlesey Creek Basin was 11.6 percent wetland as in the 1928 land economic survey (mostly in the upper basin), whereas 1.6 percent of the basin is wetland in the WISCLAND data used for the 1992-93 and reforested simulation. However, because of how the SWAT model simulates runoff from wetlands, modeled daily mean flows for peak-flow days in the presettlement simulation were higher than those for the reforested simulation. In the SWAT model, wetlands are assumed to be consistently wet or saturated, resulting in little storage and substantial runoff from the wetlands during storms. However, field observation of the wetlands in the Whittlesey Creek Basin show they likely do provide some degree of storage capacity, which would reduce rather than increase runoff peaks. It is also likely that more than 11.6 percent of the actual presettlement basin land cover was wetland. Even by 1928, some wetlands would have been drained. More research is needed to determine the extent and function of presettlement and 1992-93 wetlands in the Whittlesey Creek Basin before wetland storage can be accurately modeled. Because wetland storage could be an important part of the model that is missing, presettlement-simulation runoff and peak flows may be overestimated.

Land-Cover Effects on Water Budget

In all of the land-cover simulations, most of the precipitation falling in the Whittlesey Creek Basin (which averaged 37.3 in/yr for the modeled period) left the basin as evapotranspiration (greater than 65 percent; table 13). The remainder of the precipitation either ran off the land surface to become streamflow or percolated past the root zone of the soil profile to become ground water, which either reappeared as flow in the creek or became recharge to the deep aquifer system. Percentages of total precipitation that became evaporation, runoff, or ground water for the modeled simulations are listed in table 13. The percentages in table 13 do not add to 100 percent because of changes in soil-moisture storage in the basin during the 1999-2001 modeling period. Precipitation sources were 20.4 percent snowmelt and 79.6 percent rainfall for the 1999-2001 period.

Table 13. Percentage of total precipitation that becomes evapotranspiration, surface-water runoff, and ground water in the Whittlesey Creek, Wis., Soil and Water Assessment Tool (SWAT) model simulations

[Average annual precipitation for April 4, 1999, through October 15, 2001 was 37.3 inches per year. The percentages in this table do not add up to 100 percent because of changes in soil-moisture storage in the basin during the 1999-2001 modeling period.]

Evapotranspiration rates were slightly lower in the peak-agriculture simulation than in the other simulations. The peak-agriculture simulation also had the lowest percentage of forest. Forests have higher transpiration rates because the foliage of the trees has more surface area than that of crops and grasses and foliage has higher evaporation because it intercepts precipitation that then evaporates.

In the peak-agriculture and the developed simulations, much of the land cover had high runoff rates, so proportionately more precipitation became surface runoff and less water infiltrated past the root zone. In the presettlement and reforested simulations, land cover had lower runoff rates, so proportionately less of the precipitation became runoff.

Land-Cover Effects on Sources of Flow

Regional base flow, simulated with a constant 16.75 ft³/s, constituted nearly 75 percent of the flow in Whittlesey Creek in all land-cover simulations. The remaining flow to Whittlesey Creek simulated in the SWAT model came from either surface runoff or perched ground-water discharge. Simulated sources of all flow in Whittlesey Creek for the modeled period 1999-2001 are listed in table 14. The percentages of flow by source simulated with the SWAT model, excluding the regional base flow, are shown in parentheses in table 14.

Table 14. Percentage of streamflow at the Whittlesey Creek, Wis., U.S. Geological Survey streamflow-gaging station from surface runoff, soil-profile drainage, and regional ground water based on Whittlesey Creek Soil and Water Assessment Tool (SWAT) model simulations

[ft³/s, cubic feet per second. Numbers in parentheses are the percentages from just the SWAT model, excluding regional flow.]

The peak-agriculture and the developed simulations indicated that increased agriculture or development increased the percentage of Whittlesey Creek flow from surface runoff by 2.2 percent and 0.9 percent and decreased the percentage of the flow from perched ground water by 0.7 percent and 0.2 percent, respectively (table 14). In contrast, the presettlement and reforested simulations indicated that a more fully forested basin increased the percentage of shallow ground-water flow by 0.2 percent and 4.1 percent, respectively, and reduced surface runoff by 0.5 percent and 4.2 percent, respectively, compared to 1992-93 land-cover conditions.

Land-Cover Effects on Soil Infiltration Rates

Infiltration rates determined by the SWAT model varied by soil type and varied slightly between land-cover simulations within each soil type (table 15). Infiltration rates for sand were two to three times that for clay. Infiltration rates were lower for the peak-agriculture and developed simulations, whereas rates were highest for the presettlement simulation.

Table 15. Average ground-water infiltration rates by soil type in the Whittlesey Creek, Wis., Soil and Water Assessment Tool (SWAT) model simulations

Infiltration rates for clay soils in the SWAT model were similar to recharge rates used for the area underlain by the Miller Creek Formation in the ground-water- flow model (about 2 in/yr). Sandy-soil infiltration rates (4 to 5 in/yr) were much lower in the SWAT model than those used in the ground-water-flow model (17 in/yr) because the rates represented sands from two different areas. The SWAT model included only the surface-water-contributing basin of Whittlesey Creek, so infiltration rates for sand in the SWAT model represented mainly the sandy parts of the Miller Creek Formation, containing interbedded silts and clays. In contrast, the recharge rates for sand in the ground-water-flow model primarily represented the thick, sandy deposits of the Copper Falls Formation in the Bayfield Highlands, which have no surface-water drainage.

Infiltration rates in the developed and peak-agriculture simulations were lower than the other simulations for both sand and clay. Increased agriculture and increased development both caused a reduction in infiltration rates of soils and a decrease in the amount of time the water was on the land surface before reaching the stream; this in turn, resulted in more of the precipitation running off and reductions in ground-water recharge rates.

The reduction of infiltration due to changing land cover from forest to agriculture was more dramatic in the clay soils than sand because of the differing infiltration properties of the soils. Pure clay is relatively impervious, whereas sand is relatively permeable. Clay soils require additional structure in the soil profile (such as plant roots, organic matter, and such) to promote porosity; therefore, the infiltration rates of clay soils are susceptible to land-cover factors. Sand remains relatively porous regardless of soil structure, barring the overlayment of impervious surfaces (such as buildings or paving) or intense compaction. On the basis of generalized soil maps used in the SWAT model, subbasins 5, 6, 7, 8, 9, 10, 16, 17, 18, 19, and parts of subbasins 2, 4, and 15 are mostly clay (fig. 9) and therefore most dependent on land cover and soil structure for infiltration capacity. More detailed soils data would be needed to characterize the effects of localized land-cover changes.

The reforested and presettlement simulations both showed an increase in ground-water recharge rates from the 1992-93 simulation. Forests increase soil structure and the amount of time the runoff remains on the land surface because overland flow must pass around and through the forest floor, which is littered with detritus. Increased structure and runoff residence time on land surface resulted in increased infiltration and lower peaks. Forest cover increased the evapotranspiration rates, meaning less water is available for infiltration; however, the increase in infiltration on clay soils that change from grass to forest is greater than the effect of the increased evapotranspiration from the forest. Therefore, clay soils showed increases in local ground-water contributions to Whittlesey Creek in the forested simulations. The reforested simulation showed a slight drop in recharge rates on sand soils, a result of the increased evapotranspiration from the trees being greater than the increased infiltration caused by the changes in the soil structure from reforestation. The relation between increased infiltration and increased evapotranspiration when changing from grassland to a forest is not linear for all climates (Priestly and Taylor, 1972) and can be quite different and more dramatic in areas with different climate and vegetation than those in the Whittlesey Creek Basin. The increased presence of wetlands in the presettlement simulation resulted in increased infiltration on both soil types.

Summary and Conclusions

The U.S. Geological Survey studied the effects of land cover on flooding and base-flow characteristics of Whittlesey Creek, Bayfield County, Wis., by use of ground-water and rainfall-runoff models. The study was done in cooperation with the Bayfield County Land and Water Conservation Department and the U.S. Fish and Wildlife Service. These cooperators are interested in the protection and restoration of habitat in Whittlesey Creek and surrounding creeks for migration, spawning, and rearing of trout and salmon from Lake Superior and to protect important bird nesting areas.

The approach for the study involved collection of field data during 1999-2001 for base-flow, streambed head and temperature, meteorological data, continuous streamflow and stage, and other physical characteristics. Drillers' well logs provided data for potentiometric-surface altitudes and geologic descriptions. Geologic, soil, hydrographic, altitude, and historical land-cover data were compiled into a geographic information system. Streamflow characteristics for Whittlesey Creek—floods, base flow, seasonal fluctuations, relative contributions of runoff and ground water, and the relation between 1992-93 land-cover characteristics and base flows and flooding in Whittlesey Creek—were quantified. These data were used in two ground-water-flow models (GFLOW and MODFLOW) and a rainfall-runoff model (SWAT).

The synoptic surveys of base flow and streambed head and temperature for Whittlesey Creek and adjacent streams revealed that ground water from the deep flow system intersected Whittlesey Creek near the confluence of the main stem and North Fork. Small base flows, negative streambed heads, and warm streambed temperatures in headwater reaches of Whittlesey Creek indicated that those areas were receiving only a small amount of ground water from perched parts of the shallow ground-water-flow system upstream from the confluence. Losing reaches and negative heads above the confluence of the North Fork and main stem of Whittlesey Creek also indicated that some streamflow, derived from the perched part of the shallow system, recharged the deep system through the streambed.

Results from the MODFLOW and GFLOW ground-water-flow models were consistent with field data. Both models indicate that much of the headwaters of Whittlesey Creek do not receive ground-water contributions from the deep system. Results from both models also indicate that the deep ground-water-flow system intersected Whittlesey Creek just upstream from the confluence of the main stem and the North Fork and that the ground-water-contributing area did not coincide with the topographically delineated surface-water-contributing basin. Instead, about 90 percent of the base flow to Whittlesey Creek originated as recharge through the permeable sands of the Bayfield Highlands in the center of the Bayfield Peninsula, to the northwest of the topographically delineated Whittlesey Creek drainage basin. The land cover in the recharge area is heavily forested. Presently, base flow is 17-18 ft³/s at the streamflow-gaging station on Whittlesey Creek and through the Wildlife Refuge. Regional ground-water flow (baseflow) (74.4 percent), runoff from rainfall and snowmelt (15.7 percent), and perched ground-water flow (9.9 percent) constitute the annual mean streamflow in Whittlesey Creek (22 ft³/s) near the mouth.

The effects of changing land cover on ground-water recharge and base flow were simulated by changing recharge rates in the MODFLOW model. Deforestation could cause an increase in recharge due to decreased evapotranspiration, whereas decreased recharge could occur from development and resulting loss of pervious surfaces. Model-wide changes in recharge caused a proportional change in base flow to Whittlesey Creek. To test the sensitivity of base flow to more localized changes in recharge, recharge rates were increased by 25 and 50 percent and decreased by 25 and 50 percent within the ground-water-contributing area for Whittlesey Creek. These changes resulted in relatively small changes in base flow to Whittlesey Creek (about 2-11 percent, respectively) because the change in recharge also caused changes in the size of the resulting ground-water-contributing area. The most likely land-cover change for the ground-water-contributing area is logging of forests in the sandy zone. If this resulted in a 25 percent increase in recharge, the base flow of Whittlesey Creek would be expected to increase by about 4 percent.

The SWAT model simulated the effects of land-cover change on average annual runoff and flood peaks (daily mean flow on peak flow days) for Whittlesey Creek. The model was calibrated with data from continuous streamflow and stage data from Whittlesey Creek for the period April 1999 through October 2001. Conditions from 1992-93 WISCLAND satellite-derived land cover and four land-cover simulations were modeled. Simulations were: (1) presettlement (mainly forested before European settlement in the late 1800s); (2) peak agriculture, based on land-cover data from 1928; (3) reforested; and (4) developed (based on a hypothetical increase of 25-percent urban residential land). Presettlement and reforested simulations involved two different simulations because the 1928 land-cover map indicated 11.6 percent of the basin was wetland compared to 1.6 percent in 1992-93. In 1928, 28 percent of the basin was in row crop agriculture, compared to less than 1 percent in 1992-93. The amount of forested land was similar in 1928 and 1992-93.

Results from the SWAT model indicate that changes in land cover would have minimal effects on average annual runoff for Whittlesey Creek, but would affect flood peaks. An increase in forested land cover would result in a reduction of flood peaks by about 12-14 percent for floods with a recurrence interval of up to 100 years. The flood hydrographs from the more forested simulations had wider and broader peaks than the hydrographs for conditions in 1992-93. If land cover in 1992-93 were converted to forest, the daily mean flow for the day of the largest flood on record (370 ft³/s on April 23, 2001) would be reduced by about 12 percent. The data from nearby long-term USGS streamflow-gaging stations indicated that the event on April 23, 2001, was an 80- to 100-year event. The reduction of flood peaks under more forested conditions could potentially cause a reduction in sedimentation near the mouth of the creek, through the National Wildlife Refuge. If the basin were developed into 25-percent urban or returned to the intensive row-crop agriculture of the 1920s, daily mean flow on peak-flow days would be expected to increase up to 12 and 18 percent, respectively. These increases would be greatest for floods with dry antecedent conditions. The April 23, 2001, flood occurred under nearly saturated conditions, and the model showed an increase of only 1 percent for flood peaks for the developed and peak-agriculture simulations. The SWAT model is limited to a daily time step, which is adequate for describing the surface-water/ground-water interaction and percentage changes. It may not however, be adequate in describing instantaneous peak flow because the instantaneous peak flow in Whittlesey Creek during a flood can be more than twice the magnitude of the daily mean flow.

The models used in this study are inherently a simplistic approximation of a complex system. However, they adequately describe the effects of changing land cover on flood peaks, average annual runoff, and base flow in Whittlesey Creek. During the course of this study, several factors were identified that might improve models used in future studies. The SWAT model could be improved by the following: (1) use of a hourly time step (if this features became available in the model code), (2) a longer period of streamflow record for calibration, (3) use of more detailed county-level soil survey data, (4) increased understanding of the storage capacity of the large wetland in the upper basin, and (5) hourly precipitation data from within the drainage basin. The MODFLOW ground-water-flow model for the study could be improved by adding more detail and complexity such as additional layers, recharge zones, and zones of hydraulic conductivity. Flow from the perched part of the shallow flow system could also be incorporated into a future model.

References Cited

Anderson, E.A., 1976, A point energy and mass balance model of snow cover: U.S. Department of Commerce, National Weather Service, National Oceanic and Atmospheric Association Technical Review Report NWS 19, 150 p.

Benchmark GIS, 2001, Digital Elevation Model: Chapel Hill, N.C., scale 1:24,000.

Bradbury, K.R., and Rothschild, E.R., 1985, A computerized technique for estimating the hydraulic conductivity of aquifers from specific capacity data: Ground Water, v. 23, no. 2, p. 240–246.

Buchanan, T. J., and Somers, W. P., 1984, Discharge measurements at gaging stations: U.S. Geological Survey Techniques of Water-Resources Investigations, book 3, chap. A8, 65 p.

Cannon, W.F., Woodruff, L.G., Nicholson, S.W., and Hedgman, C.A., 1996, Bedrock geologic map of the Ashland and the northern part of the Ironwood, Wisconsin and Michigan: U.S. Geological Survey Miscellaneous Investigations Series Map I–2566, scale 1:100,000.

Clayton, Lee, 1984, Pleistocene geology of the Superior Region, Wisconsin: Wisconsin Geological and Natural History Survey Information Circular 46, 40 p.

Chow, V.T., 1959, Open-Channel Hydraulics: New York, McGraw-Hill, 680 p.

Chow, V.T., Maidment, D.R., and Mays, L.W., 1988, Applied Hydrology: New York, McGraw-Hill, 572 p.

Eichenlaub, V.L., 1979, Weather and Climate of the Great Lakes Region: Notre Dame, University of Notre Dame Press, 335 p.

Fitzpatrick, F.A., Knox, J.C., and Whitman, H.E., 1999, Effects of historical land-cover changes on flooding and sedimentation, North Fish Creek, Wisconsin: U.S. Geological Survey Water-Resources Investigations Report 99–4083, 12 p.

Freeze, R.A., and Cherry, J.A., 1979, Groundwater: Englewood Cliffs, N.J., Prentice-Hall, 604 p.

Fiering, M.B., 1967, Streamflow Synthesis: Cambridge, Mass., Harvard University Press, 159 p.

Gardner, M., 1994, Whittlesey Creek watershed news: Ashland, Wis., Bayfield County Land Conservation Dept., January, 1994, 4 p.

Gebert, W.A., 1986, Wisconsin—surface-water resources, in Moody, D.W., Chase, E.B., and Aronson, D.A., eds., National water summary 1985—hydrologic events and surface-water resources: U.S. Geological Survey Water-Supply Paper 2300, p. 485–490.

Haitjema, H.M., 1995, Analytic Element Modeling of Ground Water Flow: San Diego, Academic Press, 394 p.

Hill, M.C., 1998, Methods and guidelines for effective model calibration: U.S. Geological Survey Water-Resources Investigations Report 98–4005, 90 p.

Hunt, R.J., Lin, Y., Krohelski, J.T., and Juckem, P.F., 2000, Simulation of the shallow hydrologic system in the vicinity of Middle Genesee Lake, Wisconsin, using analytic elements and parameter estimation: U.S. Geological Survey Water-Resources Investigations Report 00–4136, 16 p.

Hunt, R.J., Anderson, M.P., and Kelson, V.A., 1998, Improving a complex finite-difference ground water flow model through the use of an analytic element screening model: Ground Water, v. 36, no. 6, p. 1011–1017.

Hunt, R.J., and Krohelski, J.T., 1996, The application of an analytic element model to investigate ground-water lake interactions at Pretty Lake, Wisconsin: Journal of Lake and Reservoir Management, v. 12, no. 4, p. 487–495.

Kennedy, E.J., 1983, Computation of continuous records of streamflow: U.S. Geological Survey Techniques of Water-Resources Investigations, book 3, chap. A13, 53 p.

Lane, L.J., 1983, Chapter 19, Transmission Losses, Soil Conservation Service, National Engineering Handbook, Section 4, Hydrology: U.S. Government Printing Office, Washington D.C., p. 19-1 to 19-21.

McDonald, M.G., and Harbaugh, A.W., 1988, A modular three-dimensional finite-difference ground-water flow model: U.S. Geological Survey Techniques of Water-Resources Investigations book 5, chap. A1 [variously paginated].

Morey, G.B., and Ojakangas, R.W., 1982, Keeweenawan sedimentary rocks of eastern Minnesota and northwestern Wisconsin, in Wold, R.J., and Hinze, W.J., eds., Geology and tectonics of the Lake Superior basin: Geological Society of America Memoir 156, p. 135–146.

Mudrey, M.G., Brown, B.A., and Greenberg, J.K., 1982, Bedrock geologic map of Wisconsin: Wisconsin Geological and Natural History Survey, scale 1:1,000,000.

Neitsch, S.L., Arnold, J.G., Kiniry, J.R., and Williams, J.R., 2001, Soil and Water Assessment Tool theoretical documentation: Version 2000, accessed February 2002 at URL: http://www.brc.tamus.edu/swat/swatdoc.html.

Nicks, A.D., 1974, Stochastic generation of the occurrence, pattern, and location of maximum amount of daily rainfall, Proceeding Symposium of Statistical Hydrology, Aug.–Sept. 1971, Tucson, Ariz.: U.S. Department of Agriculture, Miscellaneous Publication No. 1275, p. 154–171.

Ojakangas, R.W., 1986, Reservoir characteristics of the Keweenawan Supergroup, Lake Superior Region, in Mudrey, M.G., ed., Precambrian petroleum potential, Wisconsin and Michigan, 1986: Wisconsin Geological and Natural History Survey, Geoscience Wisconsin, v. 11, p. 25–31.

Phillips, D.W., and McCulloch, J.A.W., 1972, The climate of the Great Lakes Basin: Environment Canada, Climatological Studies, no. 20 [variously paginated].

Poeter, E.P., and Hill, M.C., 1998, Documentation of UCODE, a computer code for universal inverse modeling: U.S. Geological Survey Water-Resources Investigations Report 98–4080, 116 p.

Pollock, D.W., 1994, User's Guide for MODPATH/MODPATH-PLOT, Version 3—A particle tracking post-processing package for MODFLOW, the U.S. Geological Survey finite-difference ground-water flow model: U.S. Geological Survey Open-File Report 94–464 [variously paginated].

Priestly, C.H.B., and Taylor, R.J., 1972, On the assessment of surface heat flux and evaporation using large-scale parameters: Monthly Weather Review, v. 100, p. 81–92.

Red Clay Interagency Committee, 1967, Erosion and sedimentation control on the red clay soils of northwestern Wisconsin: Madison, Wis., Soil Conservation Board, 23 p.

Reese, H.M., Lillesand, T., Nagel, D.E., Stewart, J.S., Goldmann, R.A., Simmons, T.E., Chipman, J.W., and Tessar, P.A., 2002, Statewide land cover derived from multiseasonal Landstat TM data—A retrospective of the WISCLAND project: Remote Sensing of the Environment, v. 82, p. 224–237.

Rose, W.J., and Graczyk, D.J., 1996, Sediment transport, particle size, and loads in North Fish Creek in Bayfield County, Wisconsin, water years 1990–91: U.S. Geological Survey Water-Resources Investigations Report 95–4222, 18 p.

Schwarz, G.E. and Alexander, R.B., 1995, State Soil Geographic (STATSGO) Data Base for the conterminous United States (ed. 1.1): U.S. Geological Survey Open-File Report 95–449, geo-dataset on CD-ROM, scale 1:250,000.

Sloto, R.A., and Crouse, M.Y., 1996, HYSEP—A computer program for streamflow hydrograph separation and analysis: U.S. Geological Survey Water-Resources Investigations Report 96–4040, 46 p.

Soil Conservation Service, 1972, National Engineering Handbook, section 4, hydrology: Washington D.C., U.S. Government Printing Office [variously paginated].

Strack, O.D.L., 1989, Groundwater Mechanics: Englewood Cliffs, N.J., Prentice-Hall, 732 p.

U.S. Geological Survey, 2001, National Elevation Dataset: accessed May 14, 2001, at URL http://edcnts12.cr.usgs.gov/ned/default.htm.

Young, H.L., and Skinner, E.L., 1974, Water resources of Wisconsin-Lake Superior Basin: U.S. Geological Survey Hydrologic Investigations Atlas HA–524, 3 sheets.

Wisconsin Department of Agriculture and Wisconsin Geological and Natural History Survey, 1928, land economic survey: Bayfield County, Wis., scale 1:63,360.

Wisconsin Department of Natural Resources, 1996, A nonpoint source control plan for the Whittlesey Creek Priority Watershed Project: Wisconsin Department of Natural Resources Publication WR–474–96, 55 p.

Wold, R.J., 1979, Bedrock topography of Lake Superior: U.S. Geological Survey Miscellaneous Field Studies Map MF–1174, scale 1:600,000

Document Accessibility: Adobe Systems Incorporated has information about PDFs and the visually impaired. This information provides tools to help make PDF files accessible. These tools convert Adobe PDF documents into HTML or ASCII text, which then can be read by a number of common screen-reading programs that synthesize text as audible speech. In addition, an accessible version of Acrobat Reader 6.0, which contains support for screen readers, is available. These tools and the accessible reader may be obtained free from Adobe at Adobe Access.

AccessibilityFOIAPrivacyPolicies and Notices
U.S. Department of the Interior, U.S. Geological Survey Persistent URL: http://pubsdata.usgs.gov/pubs/wri/wrir-03-4130/ Page Contact Information: USGS Publishing Network Last modified: Wednesday, December 07 2016, 01:21:43 PM

---