Water-Resources Investigations Report 00-4158

U.S. Department of the Interior

U.S. Geological Survey

Delineation and Analysis of Uncertainty of 

Contributing Areas to Wells at the Southbury 

Training School, Southbury, Connecticut

Prepared in cooperation with the 

Connecticut Department of Environmental Protection and the

Connecticut Department of Mental Retardation

Y

Y

Y

Y

Y

Y
Water-Resources Investigations Report 00-4158

U.S. Department of the Interior

U.S. Geological Survey

Delineation and Analysis of Uncertainty of 

Contributing Areas to Wells at the Southbury 

Training School, Southbury, Connecticut

Prepared in cooperation with the 

Connecticut Department of Environmental Protection and the 

Connecticut Department of Mental Retardation 

By J. Jeffrey Starn, Janet Radway Stone, and John R. Mullaney

East Hartford, Connecticut

2000
Illustrations by: Stewart Crone

Manuscript preparation and editing by: B.A. Korzendorfer
U.S. DEPARTMENT OF THE INTERIOR

BRUCE BABBITT, Secretary

U.S. GEOLOGICAL SURVEY

Charles G. Groat, Director

For additional information write to: Copies of this report can be purchased from:

U.S. Geological Survey

Branch of Information Services

Box 25286, Federal Center

Denver, CO 80225

The use of firm, trade, and brand names in this report is for identification purposes only and 

does not constitute endorsement by the U.S. Geological Survey.

District Chief

U.S. Geological Survey

101 Pitkin Street

East Hartford, CT 06108
Contents III

CONTENTS

Abstract .................................................................................................................................................................................1

Introduction ...........................................................................................................................................................................2

Purpose and scope .......................................................................................................................................................2

Previous investigations ...............................................................................................................................................2

Acknowledgments ......................................................................................................................................................4

Data collection and analysis .................................................................................................................................................4

Geohydrology of the Transylvania Brook valley .................................................................................................................8

Geology .......................................................................................................................................................................8

Hydrology ...................................................................................................................................................................13

Hydraulic properties of the aquifer system .......................................................................................................13

Average ground-water recharge and relation to model calibration period .......................................................15

Water use ...........................................................................................................................................................21

Delineation of contributing areas ..........................................................................................................................................22

Watershed-scale simulation ........................................................................................................................................22

Model grid and layers .......................................................................................................................................22

Hydraulic properties of the aquifer ...................................................................................................................22

Boundary conditions and stresses .....................................................................................................................22

Model calibration ..............................................................................................................................................24

Results of the simulation ..................................................................................................................................25

Aquifer-scale simulation .............................................................................................................................................25

Model grid, layers, stresses, and time steps ......................................................................................................26

Boundary conditions and model parameter definition ......................................................................................29

Model calibration ..............................................................................................................................................31

Description of alternative models .....................................................................................................................32

Results of the simulation ..................................................................................................................................33

Diagnostic and inferential statistics ..................................................................................................................39

Analysis of uncertainty of contributing areas .......................................................................................................................40

Conditions used to estimate contributing areas ..........................................................................................................41

Deterministic contributing area ..................................................................................................................................42

Probabilistic contributing area based on variations in the defined parameters ..........................................................45

Summary and conclusions ....................................................................................................................................................50

References cited ....................................................................................................................................................................51

Appendix 1. Geohydrologic data for selected wells, Transylvania Brook drainage watershed, Connecticut ......................53
IV Contents

FIGURES

1-5. Maps showing:

 1. Location of the study area and long-term U.S. Geological Survey data-collection sites, 

Transylvania Brook watershed, Connecticut .............................................................................................................. 3

 2. Selected wells, Transylvania Brook watershed .................................................................................................... 5

 3. Observation wells, water-supply wells, and streamflow-gaging stations, Southbury Training School ............... 6

 4. Bedrock geology of the Transylvania Brook valley ............................................................................................. 9

 5. Stratified glacial deposits of the Transylvania Brook valley ............................................................................... 10

6-8. Diagrams showing:

 6. West-to-east geologic section through Transylvania Brook valley ...................................................................... 11

 7. North-to-south geologic section through Transylvania Brook valley .................................................................. 12

 8. Generalized water budget ..................................................................................................................................... 16

9-12. Graphs showing:

 9. Streamflow and ground-water recharge for the Pomperaug River, 1995-96 ....................................................... 17

10. Streamflow, precipitation, and ground-water recharge for the Pomperaug River, 

September to October 1996 ................................................................................................................................ 18

11. Water levels in four network wells and ground-water recharge in the Pomperaug River 

watershed, 1995-96 .............................................................................................................................................. 19

12. Water levels in observation well E2D, Southbury Training School ................................................................... 20

13-16. Maps showing:

13. Extent of watershed-scale model grid, Transylvania Brook watershed .............................................................. 23

14. Extent of aquifer-scale model grid and boundary conditions, Transylvania Brook valley ................................. 27

15. Altitude of top of model layer 3, aquifer-scale model, Transylvania Brook valley ............................................ 28

16. Aquifer zones, aquifer-scale model, Transylvania Brook valley ........................................................................ 30

17, 18. Graphs showing:

17. Simulated and observed water-level changes in coarse-grained deposits during the 

aquifer test, October 30 to November 4, 1996, Southbury Training School ...................................................... 35

18. Simulated and observed steady-state water levels, October 30 to November 4, 1996, 

Southbury Training School ................................................................................................................................. 36

19. Map showing altitude of simulated water table before and after the aquifer test, October 30 to 

November 4, 1996, Southbury Training School ......................................................................................................... 37

20. Graph showing observed and simulated water levels in observation well E2D, Southbury Training School ............... 38

21-25. Maps showing:

21. Deterministic contributing areas to wells PW-1 and PW-3, Southbury Training School ................................... 43

22. Selected pathlines to well PW-3, Southbury Training School ............................................................................ 44

23. Probabilistic contributing area to well PW-3, Southbury Training School ......................................................... 47

24. Probabilistic contributing area to well PW-1, Southbury Training School ......................................................... 48

25. Probabilistic contributing areas to wells PW-1 and PW-3, Southbury Training School .................................... 49
Contents V

TABLES

1. Data on observation wells drilled for the aquifer test, Southbury Training School, Connecticut .................................... 4

2. Instantaneous streamflow at Transylvania Brook, Connecticut........................................................................................ 7

3. Median values of hydraulic properties of hydrogeologic units ........................................................................................ 14

4. Summary table of transmissivity and storage from an aquifer test at well PW-3, Southbury Training School .............. 15

5. Ground-water levels in network wells in the Pomperaug River watershed ..................................................................... 19

6. Construction details of water-supply wells, Southbury Training School ......................................................................... 21

7. Recharge to valley from upland basins, Transylvania Brook drainage basin .................................................................. 25

8. Statistical measures of model fit of alternative models .................................................................................................... 33

9. Optimal parameter values for alternative models ............................................................................................................ 34

10. Correlation and variance/covariance for alternative model CAL0 ................................................................................ 40

CONVERSION FACTORS, VERTICAL DATUM, AND ABBREVIATIONS

Vertical datum: In this report, "sea level" refers to the National Geodetic Vertical Datum of 1929 (NGVD of 1929)—a geodetic

datum derived from a general adjustment of the first-order level nets of both the United States and Canada, formerly called

Sea Level Datum of 1929.

Altitude, as used in this report, refers to distance above or below sea level.

Transmissivity: The standard unit for transmissivity is cubic foot per day per square foot times foot of aquifer thickness 

[(ft3/d)/ft2]ft. In this report, the mathematically reduced form, foot squared per day (ft2/d), is used for convenience.

Multiply By To obtain

Length

inch (in.) 25.4 millimeter

foot (ft) 0.3048 meter

mile (mi) 1.609 kilometer

Area

acres 4,047 square meter

square foot (ft2) 0.09290 square meter

square mile (mi2) 2.590 square kilometer 

Volume

million gallons (Mgal) 3,785 cubic meter

Flow rate

cubic foot per second (ft3/s) 0.02832 cubic meter per second 

cubic foot per second per square mile (ft3/s/mi2) 0.01093 cubic meter per second per square kilometer

gallon per minute (gal/min) 0.06309 liter per second

gallon per day (gal/d) 0.003785 cubic meter per day

inch per year (in/yr) 25.4 millimeter per year

Hydraulic conductivity

foot per day (ft/d) 0.3048 meter per day

Transmissivity

foot squared per day (ft2

/d) 0.09290 meter squared per day
ABSTRACT 1

Delineation and Analysis of Uncertainty of Contributing 

Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

By J. Jeffrey Starn, Janet Radway Stone, and John R. Mullaney

ABSTRACT

Contributing areas to public-supply wells at the 

Southbury Training School in Southbury, Connecticut, 

were mapped by simulating ground-water flow in stratified glacial deposits in the lower Transylvania Brook 

watershed. The simulation used nonlinear regression 

methods and informational statistics to estimate parameters of a ground-water flow model using drawdown 

data from an aquifer test. The goodness of fit of the 

model and the uncertainty associated with model 

predictions were statistically measured.

A watershed-scale model, depicting large-scale 

ground-water flow in the Transylvania Brook watershed, was used to estimate the distribution of groundwater recharge. Estimates of recharge from 10 small 

basins in the watershed differed on the basis of the 

drainage characteristics of each basin. Small basins 

having well-defined stream channels contributed less 

ground-water recharge than basins having no defined 

channels because potential ground-water recharge was 

carried away in the stream channel. 

Estimates of ground-water recharge were used in 

an aquifer-scale parameter-estimation model. Seven 

variations of the ground-water-flow system were 

posed, each representing the ground-water-flow system 

in slightly different but realistic ways. The model that 

most closely reproduced measured hydraulic heads and 

flows with realistic parameter values was selected as 

the most representative of the ground-water-flow 

system and was used to delineate boundaries of the 

contributing areas. The model fit revealed no systematic model error, which indicates that the model is 

likely to represent the major characteristics of the 

actual system.

Parameter values estimated during the simulation are as follows: horizontal hydraulic conductivity of 

coarse-grained deposits, 154 feet per day; vertical 

hydraulic conductivity of coarse-grained deposits, 

0.83 feet per day; horizontal hydraulic conductivity of 

fine-grained deposits, 29 feet per day; specific yield, 

0.007; specific storage, 1.6E-05. Average annual 

recharge was estimated using the watershed-scale 

model with no parameter estimation and was determined to be 24 inches per year in the valley areas and 

9 inches per year in the upland areas.

The parameter estimates produced in the model 

are similar to expected values, with two exceptions. 

The estimated specific yield of the stratified glacial 

deposits is lower than expected, which could be caused 

by the layered nature of the deposits. The recharge estimate produced by the model was also lower—about 

32 percent of the average annual rate. This could be 

caused by the timing of the aquifer test with respect to 

the annual cycle of ground-water recharge, and by 

some of the expected recharge going to parts of the flow 

system that were not simulated. The data used in the 

calibration were collected during an aquifer test from 

October 30 to November 4, 1996. The model fit was 

very good, as indicated by the correlation coefficient 

(0.999) between the weighted simulated values and 

weighted observed values. The model also reproduced 

the general rise in ground-water levels caused by 

ground-water recharge and the cyclic fluctuations 

caused by pumping prior to the aquifer test.

Contributing areas were delineated using a 

particle-tracking procedure. Hypothetical particles of 

water were introduced at each model cell in the top 

layer and were tracked to determine whether or not 

they reached the pumped well. A deterministic contributing area was calculated using the calibrated model, 

and a probabilistic contributing area was calculated 

using a Monte Carlo approach along with the calibrated 

model.
2 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

The Monte Carlo simulation was done, using the 

parameter variance/covariance matrix generated by the 

regression model, to estimate probabilities associated 

with the contributing area to the wells. The probabilities arise from uncertainty in the estimated parameter 

values, which in turn arise from the adequacy of the 

data available to comprehensively describe the groundwater-flow system and the validity of the parameter 

definitions. The Monte Carlo data sets were conditioned to remove unrealistic parameter sets. Probabilities in the contributing area range from 1 to 100 

percent. The highest probabilities (greater than 50 

percent) are in the coarse-grained deposits that ring the 

head of the valley; this area is consistent with the deterministic contributing area defined using the estimated 

parameter values. The probabilities do not reflect 

subsurface variabilities within the defined parameter 

structure, but in this problem, the large-scale variations 

are expected to dominate contributing area uncertainty.

The contributing area shows that most water 

enters the stratified glacial deposits through the coarsegrained deposits that ring the head of the lower Transylvania Brook watershed. Some of these deposits are 

not saturated throughout the year and could not be 

simulated in the model; however, because the primary 

public-supply well receives most of its water from this 

area, the unsaturated deposits should be considered to 

be within the contributing area. The travel times for 

ground water in this area are less than 1 year, based on 

an assumed porosity of 0.20. Travel times for most of 

the rest of the contributing area are less than 2 years. 

Contributing areas for alternative models are similar, 

indicating that nonuniqueness in the design of the 

model does not seem to be a problem.

INTRODUCTION

The Connecticut Aquifer Protection Program 

(section 22a-354 of the Connecticut General Statutes) 

requires water suppliers to delineate aquifer-protection 

areas for all wells that obtain water from stratified 

glacial deposits and that provide water to more than 

1,000 people. State regulations specify the use of 

ground-water-flow models to delineate aquifer-protection areas. Ground-water-flow models commonly are 

calibrated using a trial-and-error method that is subjective, and the models do not provide a measure of the 

uncertainty of model predictions that is inherent in any 

type of model.

The Southbury Training School (STS), owned 

and operated by the Connecticut Department of Mental 

Retardation (DMR), maintains its own water-supply 

system from three wells that tap stratified glacial 

deposits in the Transylvania Brook watershed in central 

Connecticut (fig. 1). The STS supplies water to approximately 1,862 people; therefore, aquifer-protection 

areas must be mapped. In 1996, the U.S. Geological 

Survey (USGS), in cooperation with the Connecticut 

Department of Environmental Protection (DEP) and 

the Connecticut DMR, began a study to consider alternatives to traditional ground-water-flow modeling. The 

study, using the STS as a test site, was undertaken to 

develop and demonstrate methods that can be used to 

quantify the uncertainty in model predictions of 

contributing areas to wells, as well as meet the statutory 

requirements for delineating contributing areas. 

Modeling methods were enhanced in this study (1) to 

provide an objective model calibration that considers 

the amount of available data and (2) to estimate the 

uncertainty of model predictions of the contributing 

areas.

Purpose and scope

This report describes the hydrologic analyses 

and simulation modeling that were used to delineate 

contributing areas to public-water supply wells at the 

STS. The report presents an analysis of the uncertainty 

in the contributing areas, as determined by simulation 

models. It also presents geologic information for Transylvania Brook watershed and summarizes previous 

investigations and data collection in the study area.

Previous investigations

A previous investigation of the aquifer at the STS 

was done by Mazzaferro (D.L. Mazzaferro, Ground 

Water Inc., written commun., 1991). This study 

includes an estimate of the contributing area to the 

wells, details of water-supply and test-well construction, some geologic data, and estimates of the transmissivity of the aquifer. Details about the operation of the 

STS well field also are available from the Water-Supply 

Plan for the Southbury Training School (Corinne 

Fitting, Connecticut Department of Environmental 

Protection, oral commun., 1996). Meinzer and Stearns 

(1929) and Mazzaferro (1986) studied the hydrology of 

the Pomperaug River watershed, which is
INTRODUCTION 3

Figure 1. Location of the study area and long-term U.S. Geological Survey data-collection sites, Transylvania Brook 

watershed, Connecticut.
4 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

adjacent to the Transylvania Brook watershed (fig. 1). 

The climate, topography, and geology of the Pomperaug River watershed are similar to those of the Transylvania Brook watershed, and much of the hydrologic 

information is transferable, particularly with regard to 

recharge rates and water budget. Wilson and others 

(1974) described the hydrology of the lower Housatonic River watershed and presented information on the 

mean runoff of the Transylvania Brook watershed, the 

flow frequency duration for the Pomperaug River at 

Southbury, and the general yield of bedrock wells in the 

watershed.

Acknowledgments

The authors wish to thank the USGS personnel 

that contributed substantially to this report: George 

Casey provided the GIS analysis used in many aspects 

of the project, Timothy Frick provided the estimates of 

water use for the study area, Remo Mondazzi logged 

the boreholes and collected water-level and streamflow 

measurements throughout the project, and Steven 

Kiesman surveyed the well locations and instrumented 

the wells with pressure transducers for the aquifer test. 

The authors also thank Al Van Geersdaele, STS, for the 

high level of cooperation received in conducting the 

aquifer test using the STS water-supply infrastructure. 

DATA COLLECTION AND ANALYSIS

Sixteen boreholes were drilled at 11 locations for 

this study (table 1; figs. 2 and 3) to determine the depth 

to bedrock, altitude of the water level, and distribution 

of geologic units. Observation wells were completed in 

each borehole so that water levels could be monitored 

during an aquifer test. The observation wells are 

labelled N (north), E (east), or W (west), for direction 

from the main pumped well at STS, with a number

(1 being closest to the pumped well and increasing in 

distance from the pumped well). A letter at the end of 

the well label indicates a deep (D) or shallow (S) well. 

To supplement geologic information from the boreholes, 189 well-completion reports were obtained from 

files at the Connecticut DEP, inspected, and plotted on 

maps. Of these reports, 29 domestic wells were 

selected for inclusion in the USGS Ground-Water Site

Table 1. Data on observation wells drilled for the aquifer test, Southbury Training School (STS), Connecticut

[—, rock not encountered; all wells shown on fig. 3 except W3, which is shown on fig. 2]

Well 

number

Distance from 

well 

PW-3, 

in feet

Altitude of 

land surface, 

in feet above 

sea level

Altitude of 

rock (refusal), 

in feet above 

sea level

Altitude of top 

of screen, in 

feet above 

sea level

Altitude of 

bottom of 

screen, in feet 

above sea 

level

Altitude of 

water level, in 

feet above 

sea level

U.S. 

Geological 

Survey local 

well number1

1

Well number used in the U.S. Geological Survey Ground-Water Site Inventory Database

N1S 31.89 173.81 — 160.81 150.81 162.32 SB 101

N1D 34.57 173.77 96.77 114.77 112.77 162.44 SB 100

N2S 110.58 177.42 — 164.42 154.42 164.52 SB 103

N2D 114.13 177.59 90.59 118.59 116.59 162.83 SB 102

N3 260.81 175.23 99.23 130.23 128.23 163.38 SB 104

N5 1,500 175.55 111.55 127.55 125.55 178.22 SB 105

E1S 64.13 169.30 — l63.30 153.30 161.04 SB 107

E1D 60.94 169.33 93.33 110.33 108.33 162.39 SB 106

E2S 135.80 165.48 — 162.48 152.48 159.34 SB 109

E2D 136.50 165.38 85.38 105.38 103.38 162.33 SB 108

E3 470.00 152.10 72.10 87.10 85.10 159.84 SB 110

E4 1,050 215.30 — 121.30 119.30 161.52 SB 111

W1S 61.27 174.50 — 153.50 151.50 162.23 SB 113

W1D 62.52 174.56 112.56 128.56 126.56 162.22 SB 112

W2 227.32 181.52 131.52 151.52 149.52 164.74 SB 114

W3 1,750 185.69 160.69 162.69 160.69 170.02 SB 115
DATA COLLECTION AND ANALYSIS 5

Figure 2. Selected wells, Transylvania Brook watershed, Connecticut.
6 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 3. Observation wells, water-supply wells, and streamflow-gaging stations, Southbury Training School, Connecticut.
DATA COLLECTION AND ANALYSIS 7

Inventory (GWSI) database because of the quality and 

uniqueness of the information on the report (appendix 

1; fig. 2). Domestic wells are labelled with the GWSI 

site identification numbers (SB, for the town of Southbury, and a sequential well number for that town).

Streamflow was measured at two sites (fig. 3; 

table 2). The upstream site represents water flowing 

from Stibbs Lake into Transylvania Brook and into the 

area of the well field. The downstream site is just south 

of the well field and represents water that leaves the 

well field area. Historical measurements were available 

from the downstream site, as well as a third site at the 

mouth of Transylvania Brook (table 2; fig. 1). 

Table 2. Instantaneous streamflow at Transylvania Brook, Connecticut

[—, not measured; NA, not applicable]

Date

Streamflow at USGS gaging station

(locations shown on figs. 1 and 3)

Loss (-) or gain 

(+), in cubic 

feet per second between 

stations 

01204340 and 

01204350

Standard 

deviation1 of 

loss or gain, in 

cubic feet per 

second

1

Standard deviation is calculated as the square root of the sum of the variances of each measurement assuming each measurement to be 

accurate within 5 percent.

01204340

(upstream from

well field)

01204350

(downstream 

from well field)

01204400

(at mouth of 

Transylvania 

Brook)

June 10, 1965 — 1.27 3.73 — NA

July 15, 1965 — 0.29 0.98 — NA

September 27, 1965 — .15 .89 — NA

March 2, 1966 — 10.2 30. — NA

August 17, 1966 — .44 1.09 — NA

July 28, 1967 — — 1.38 — NA

October 14, 1966 — .97 — — NA

September 24, 1968 — — 1.54 — NA

June 27, 1969 — — 3.64 — NA

June 18, 1970 — — 4.07 — NA

July 20, 1971 — — 2.56 — NA

September 25, 1972 — — 1.14 — NA

August 28, 1973 — — 1.84 — NA

August 29, 1996 0.77 0.82 — +0.05 0.06

September 19, 1996 6.77 6.92 — +.15 .49

October 17, 1996 2.80 2.84 — +.04 .20

October 21, 1996 55.5 49 — -6.5 3.71

October 29, 1996 13.7 14 — +.03 .98

November 4, 1996 8.73 8.58 — -.15 .62

April 10, 1997 14.0 12.8 — -1.2 .95
8 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

An aquifer test was conducted at well PW-3 from 

October 30 to November 4, 1996 to determine the 

hydraulic properties of the aquifer (transmissivity, storativity, and boundary conditions). During the aquifer 

test, water levels were monitored in all observation 

wells, and inside and outside two streambed piezometers installed in Transylvania Brook (fig. 3); stage was 

measured at staff plates at each streamflow-gaging 

station. Streamflow was calculated from stage 

measurements based on stage/discharge relations 

determined at each streamflow-gaging station.

GEOHYDROLOGY OF TRANSYLVANIA 

BROOK VALLEY

The STS covers 1,600 acres, mostly in the watershed of Transylvania Brook, a tributary of the Pomperaug River. The Transylvania Brook watershed, like 

the rest of western Connecticut, is underlain by three 

principal hydrogeologic units—bedrock, glacial till, 

and stratified glacial deposits (often called “stratified 

drift”) (figs. 4 and 5). The watershed lies along the 

western edge of a bedrock structural basin (Gates, 

1954, 1959; Scott, 1974; Stanley and Caldwell, 1976; 

Rodgers, 1985). Stratified glacial deposits at the STS 

partially fill a bedrock valley in the southern half of the 

Transylvania Brook watershed. 

Geology

The bedrock of the structural basin (fig. 4) is 

similar in character and structure to that of the Hartford 

Basin in central Connecticut. Highlands to the east of 

Transylvania Brook are underlain by Mesozoic-age 

sedimentary bedrock, which is predominantly redbrown and includes arkosic sandstone, conglomerate, 

and shale. The highest ridges to the east reach 500 to 

650 ft in altitude and are composed of extrusive 

igneous basalt (also of Mesozoic age). Highlands to the 

west of Transylvania Brook are underlain by Paleozoic-age crystalline bedrock that is predominantly 

schist; these highlands reach altitudes greater than 

650 ft. The contact between sedimentary bedrock to the 

east and crystalline (metamorphic) bedrock to the west 

lies beneath glacial sediments in the Transylvania 

Brook watershed (fig. 4). Red-brown sandstone was 

encountered in the bottom of several boreholes (N2D, 

N3, N5, E1D, and E2D), indicating that the contact 

between sedimentary and crystalline bedrock lies 

farther west in the valley than shown on published 

geologic maps (Stanley and Caldwell, 1976; Rodgers, 

1985). This contact may be a fault zone or unconformity, although it is not shown as such on published 

geologic maps. 

Unconsolidated glacial deposits cover bedrock 

in most places in the Transylvania Brook watershed. 

These materials were deposited during the advance and 

retreat of Pleistocene continental glaciers, particularly 

the last (late-Wisconsinan) glaciation. Surficial materials have been mapped at a regional scale (Stone and 

others, 1992) and include (1) glacial till, which was laid 

down directly by ice on top of bedrock and is the surficial material on the valley sides and in the uplands; and 

(2) stratified glacial deposits, which partially fill the 

bedrock valley to an altitude of 200 to 250 ft (figs. 5, 6, 

and 7). 

Till consists of a nonsorted, generally nonstratified mixture of grain sizes ranging from clay to large 

boulders. The till matrix is composed dominantly of 

sand and silt, although boulders in and on the surface 

commonly are abundant. Till is typically a dense and 

compact material due its mode of deposition beneath 

the great weight of the ice sheet and to the presence of 

silt and clay in the matrix. This type of till commonly 

is identified as “hardpan” in well logs recorded by well 

drillers. A sandier and stonier, less dense and less 

compact facies of till (ablation till) may be present in 

some places. The color and lithology of till generally 

reflect the composition of the underlying bedrock to the 

north from which the till was derived. In the Transylvania Brook watershed, till is described in well logs as 

gray, having been derived predominantly from the 

crystalline bedrock underlying the highlands to the 

northwest. Red-brown till derived from the red-brown 

sedimentary rocks underlying the eastern side of Transylvania Brook and the highland to the east is likely to 

be present locally in the eastern part of the study area; 

thin red-brown till overlying red sandstone was 

encountered in boreholes E3 and N5 drilled at the STS 

during this study. Numerous well logs from the residential subdivision on the hillside northeast of the STS 

indicate the presence of thick gray till overlying redbrown sedimentary rock. Well logs for domestic wells 

inventoried indicate that till ranges from less than 3 ft 

to as much 140 ft in thickness in the study area. Several 

boreholes drilled at the site showed locations where till 

is absent and stratified glacial deposits directly overlie 

bedrock; other boreholes showed the presence of till.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 9

Figure 4. Bedrock geology of the Transylvania Brook valley, Connecticut.
10 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 5. Stratified glacial deposits of the Transylvania Brook valley, Connecticut.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 11

350

300

250

200

150

100

50

E4

E3 projected

Transylvania Brook

E1D, E1S

E2D, E2S

W2

W1D, W1S

PW-3 projected; SECTION B-B'

EXPLANATION

s

al

Fine sand and silt

Sand and gravel

Till

Sedimentary bedrock--lines show bedding

Alluvium

Sand

Crystalline bedrock

Well/well cluster and location(s)

E2D 

of screened interval

10-ft well screen

2-ft well screen

A A'

FEET

al

FEET

s

s

Datum is sea level Vertical exaggeration X 8

METERS

0 500 1000 1500 2000 FEET

0 150 300 450 600

350

300

250

200

150

100

50

Figure 6. East-to-west geologic section through Transylvania Brook valley, Connecticut.
12 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

EXPLANATION

Fine sand and silt

Sand and gravel

Till

Sedimentary bedrock--section

perpendicluar to bedding

Alluvium

s Sand

al

Crystalline bedrock

Well/well cluster and location(s)

E2D

of screened interval

10-ft well screen

2-ft well screen

Vertical exaggeration X 8 N5

Stibbs Lake

E3 projected

PW-3; SECTION A-A'

N1D, N1S

N2D, N2S

N3

Transylvania Brook

Sewage treatment plant

Transylvania Brook

al

al

B B'

FEET

50

100

150

200

250

0

300

FEET

Datum is sea level

0 500 1000 1500 2000

0 300 150 450 600 METERS

FEET

s s

50

100

150

200

250

0

300

Figure 7. North-to-south geologic section through Transylvania Brook valley, Connecticut.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 13

Stratified glacial deposits consist of layers of 

well to poorly sorted gravel, sand, silt, and clay laid 

down by meltwater in glacial lakes and streams that 

occupied valleys during retreat of the ice sheet. The 

stratified deposits in the Transylvania Brook watershed 

near the STS consist of both coarse-grained and finegrained sediments graded to or deposited in glacial 

lakes that occupied the Pomperaug and Transylvania 

Brook valleys in progressively northward positions in 

front of the retreating ice sheet. Coarse-grained 

deposits (gravel, sand and gravel, and sand) were laid 

down as ice-marginal deltas and fluviodeltaic deposits 

in close proximity to the retreating ice sheet. Finegrained deposits (very fine sand, silt, and clay) accumulated farther from the ice margin in the quiet water 

conditions of the glacial lake bottom. Near the STS, 

coarse-grained deposits, including sand and gravel and 

sand units, locally make up the entire thickness of the 

stratified section, particularly west of Rt. 172 and at the 

north end of the valley. In much of the area east of Rt. 

172, coarse-grained deposits underlie, overlie, or interfinger laterally with fine-grained deposits. The main 

production well (PW-3) at the STS and installed observation wells (W2, W1D, E1D, E2D, N1D, N2D, N3, 

and N5), are screened in a subsurface sand and gravel 

unit, which ranges from a few to about 40 ft in thickness. To the east, north, and south of PW-3, the sand 

and gravel unit is overlain by fine-grained deposits 

ranging from a few to 60 ft in thickness; to the west of 

PW-3, the sand and gravel unit is overlain by 30 to 

40 ft of sand (fig. 6). The subsurface sand and gravel 

unit likely extends west of Rt. 172 to connect with 

coarse-grained deposits on the valley side and north 

beneath Stibbs Lake at least as far as observation well 

N5. The unit pinches out to the east and south (figs. 6 

and 7). Split-spoon samples of the subsurface coarsegrained deposits indicate that the unit consists of interbedded layers of pebbly, medium to coarse sand, 

granule to pebble gravel, and coarse to very coarse sand 

and granules. Large pebbles and cobbles also may be 

present, although these were not sampled by the splitspoon sampler. The fine-grained deposits consist of 

thinly laminated very fine sand, silt, and clay. Individual layers are typically 0.1 to 1.0 in. in thickness and 

alternate between coarser (very fine sand to fine sand)

and finer layers (silt and/or clay). Coarser layers are 

typically thicker than finer layers.

Hydrology

Water in the Transylvania Brook watershed 

flows from till-covered uplands in stream channels and 

as subsurface flow over and into the stratified glacial 

deposits in the valley bottom. Some water from upland 

stream channels also seeps into the stratified glacial 

deposits. Most runoff (both ground-water and surfacewater) discharges to Transylvania Brook. Some 

ground-water runoff is intercepted by water-supply 

wells, both domestic and public; this water is returned 

after treatment to the aquifer or to Transylvania Brook. 

Hydraulic properties of the aquifer system

Three principal types of subsurface materials are 

found in the watershed—coarse-grained stratified 

glacial deposits; fine-grained deposits, which consist of 

till and fine-grained stratified deposits; and bedrock, 

which consists of sedimentary and crystalline rocks 

(figs. 4 and 5). The materials are grouped in this way 

for hydrogeologic purposes, which is somewhat 

different than the geologic grouping of units, because 

of the general similarity in their hydraulic properties. 

Melvin and others (1992) summarized the hydraulic 

properties of these materials in Connecticut on the 

basis of published aquifer test and laboratory test 

results (table 3), and these results are applicable to the 

materials at the STS site. 

Coarse-grained stratified deposits have the 

highest hydraulic conductivity because of the openness 

of their interconnected pore spaces, and thus form the 

most productive part of the ground-water system. Till 

and fine-grained stratified deposits have a lower 

hydraulic conductivity than the coarse-grained deposits 

because fine sand, silt, and clay occupy more of the 

open pore space within the deposits. Theoretically, the 

hydraulic conductivity of till, which is nonsorted and 

nonlayered, should be the same in all directions, but the 

hydraulic conductivity of the fine-grained deposits, 

composed of layered material, should be much higher 

in the horizontal direction (along the relatively coarser 

layers) than in the vertical direction (through both finer 

and coarser layers) (table 3). 
14 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Table 3. Median values of hydraulic properties of hydrogeologic units in Connecticut

[Modified from Melvin and others, 1992, table 1; —, no data or insufficient data]

Unit

Hydraulic 

conductivity, 

in feet per day

Storativity, 

dimensionless

Porosity, 

in percent Orientation

Coarse-grained stratified drift 170 .36 — —

Fine-grained stratified drift, silt .14 .29 — vertical

Fine-grained stratified drift, clay 0.0001

.82

—

—

—

—

vertical

horizontal

Loose surface till, crystalline-rock 

provenance

2.7 0.28 — —

Loose surface till, sedimentaryrock provenance.71 — 32 horizontal and 

vertical

Compact surface till .007 — — horizontal

Crystalline bedrock .6

Sedimentary bedrock 4.7

In bedrock, ground water flows primarily 

through fractures in the rock rather than through the 

pore spaces in the rock. Water-bearing pathways in 

fractured bedrock include nearly horizontal stress-relief 

fractures, layer-parallel fractures in layered crystalline 

rocks and sedimentary rocks, and high-angle to vertical 

fractures caused by movements of the earth’s crust. 

Bedrock wells typically penetrate many fractures, but 

only a few of these fractures may produce water. The 

passage of water from unconsolidated till and stratified 

deposits into or out of bedrock can be impeded by a 

compact till that overlies bedrock in most of Connecticut; this till, shown between stratified glacial deposits 

and bedrock in figs. 6 and 7, was compacted beneath 

the glaciers and thus has a lower hydraulic conductivity 

than loose surface till. The fractured bedrock is the 

source of water to all domestic wells inventoried in the 

area. The average yield of inventoried bedrock wells is 

less than 5 gal/min; however, the possibility of a highly 

fractured zone in the valley bottom is indicated by 

reported yields of 50 to 100 gal/min from several 

domestic wells in the area.

D.L. Mazzaferro (Ground Water, Inc. written 

commun., 1991) estimated the average transmissivity 

of the coarse-grained deposits at the STS well field to 

be 8,200 ft2

/d, using the specific capacity of the three 

production wells and grain-size analysis of two test 

borings. Mazzaferro estimated transmissivity at the 

production wells from 6,400 ft2

/d to 12,100 ft2/d and at 

the test borings from 6,400 to 9,500 ft2

/d. These estimates, however, which were based on specific capacity 

and grain-size analysis, are considered to be less accurate than estimates based on aquifer-test analysis.

An aquifer test was conducted by the USGS from 

October 30 to November 4, 1996 at STS well PW-3. 

During the test, well PW-1 was not pumping. Transmissivity and storativity were calculated using the semilogarithmic method described by Cooper and Jacob 

(1946) (table 4). The analyses of the aquifer-test data 

are on file at the Connecticut District office of the 

USGS. There were slight changes in the slope of the 

semi-logarithmic plot of drawdown with time that indicated an additional source of water to the aquifer, such 

as leakage from a surface-water body, recharge from 

precipitation, or leakage from fine-grained deposits. 

The drawdown curves were analyzed both at early and 

late times in the test.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 15

Table 4. Summary of transmissivity and storage from an aquifer test at well PW-3, Southbury Training School, 

Connecticut

[Well locations shown on fig. 3. S, shallow water-table well; D, deep well]

Well

Distance from 

pumped well,

in feet

Transmissivity, in 

feet squared per day

early/late time 

Storage

coefficient

early/late time

Drawdown at 

1,000 minutes into 

test, in feet

N1S (early)

N1S (late)

32 3,100

3,600

.007

.004

5.3

N2D (early)

N2D (late)

114 3,400

5,000

.0004

.00015

4.99

N3 (early)

N3 (late)

261 3,400

6,500

.0003

.00008

3.69

E1D (early)

E1D (late)

61 3,300

4,660

.0006

.00006

6.01

E2D (early)

E2D (late)

136.5 3,100

5,440

.0005

.00004

4.98

W1D (early)

W1D (late)

63 3,400

4,660

.0006

.0003

5.88

W1S (early)

W1S (late)

61.3 2,883

3,436

.005

.003

5.03

Average ground-water recharge and relation to 

model calibration period

To perform a transient ground-water simulation, 

as required by the aquifer protection mapping process, 

it is necessary to understand (1) the long-term average 

steady-state ground-water recharge in the basin of 

interest and (2) the relation between the simulated time 

period and average ground-water recharge. The relation of the study period to long-term average conditions 

will be discussed using long-term streamflow and 

ground-water level data from the adjacent Pomperaug 

River watershed. Limited data from the STS site also 

will be used to assess the hydrologic conditions at the 

time of the aquifer test.

Effective ground-water recharge is defined in 

this study as the amount of water that infiltrates from 

the land surface into the aquifer minus evapotranspiration from the aquifer (fig. 8). Average annual recharge 

was estimated using a regression-derived linear relation between ground-water outflow and the percentage 

of the drainage area underlain by coarse-grained 

deposits for Connecticut (Mazzaferro and others, 

1979). In Mazzaferro’s study, ground-water outflow 

was determined by hydrograph separation. Groundwater outflow is assumed to be a conservative estimate 

of recharge if changes in ground-water storage are 

small. The relation is 

, (1)

where:

is ground-water outflow as a percentage of total 

runoff, and

is the percentage of the watershed underlain by 

coarse-grained stratified glacial deposits.

Stratified glacial deposits occupy 17 percent of 

the land surface in the Transylvania Brook watershed, 

so that solution of equation 1 indicates that basin-wide 

recharge is 45.2 percent of total annual runoff. Mean 

annual runoff for the Transylvania Brook area is 1.71 

ft3

/s/mi2 (Wilson and others, 1974); therefore, groundwater outflow is calculated to be 3.89 ft3

/s, or 

10.5 in/yr over the 5.03-mi2

 area of the drainage basin.

Y= 35 0.6 + X

Y

X
16 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Precipitation

Land surface

Aquifer

Evapotranspiration

Transylvania Brook

Change in aquifer storage

Figure 8. Generalized water budget.

Water use

Mazzaferro and others (1979) observed that areas 

underlain by stratified glacial deposits have a recharge 

rate that is 2.7 times higher, on average, than areas 

underlain by till. Using this relation and the percentages 

of the basin underlain by till and stratified deposits, 

recharge to stratified glacial deposits is calculated to be 

about 22 in. (1.41 ft3/s) and recharge to till is about 

8 in. (2.48 ft3

/s). Although the rate of recharge to till is 

lower, the volume of water recharged is greater because 

till covers a larger area of the basin than do stratified 

deposits. 

Ground-water recharge occurs in annual cycles, 

and in order to simulate values of recharge in the model, 

it is necessary to understand how the time period of the 

aquifer test relates to the recharge cycle. Also, because 

steady-state hydrologic conditions are required as 

starting conditions for the transient model, it is necessary to assess whether the conditions at the time of the 

aquifer test approximate steady-state conditions, and if 

not, to determine what type of adjustment could be 

made to the data to make them comparable to steadystate values. Long-term data (more than 1 year) are 

needed to analyze the yearly cycle of recharge. No 

long-term data-collection sites are present in the Transylvania Brook watershed, but five USGS network 

wells and a streamflow-gaging station are in the adjacent Pomperaug River watershed (fig. 1). 

The Pomperaug River watershed has virtually 

the same surficial geology, climate, and topography as 

the Transylvania Brook watershed, so the hydrologic 

response of each basin to precipitation is expected to be 

similar. The main difference between the watersheds is 

size—the Pomperaug River watershed above the 

streamflow-gaging station is about 75 mi2

 and the 

simulated part of the Transylvania Brook watershed is 

about 5 mi2

. Long-term records are available for 

several streams of similar size to Transylvania Brook, 

but at more distant locations in Connecticut. Peaks on 

the streamflow hydrographs of these streams correspond closely with peaks on the hydrograph of the 

Pomperaug River; therefore the Pomperaug River is 

deemed to be an adequate surrogate for Transylvania 

Brook.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 17

Ground-water flow sustains streamflow between 

rains, and ground-water recharge can be estimated by 

analyzing streamflow records. Streamflow records can 

be used to calculate ground-water recharge using three 

computer programs—PART, RORA, and PULSE 

(Rutledge, 1997; 1998). These programs separate 

streamflow into ground-water and surface-water 

components and can be used to estimate ground-water 

discharge from a basin. PART is based on an empirical 

analysis of the streamflow record. RORA and PULSE 

are based on the drainage characteristics of a basin and 

require an estimate of the master recession curve 

(Rutledge, 1997; 1998). The master recession curve is 

fit to the recession segments of peaks; the difference in 

streamflow between adjacent recession curves is 

ground-water recharge. In this study, the results of all 

three programs are similar (fig. 9).

Streamflow and recharge data for the Pomperaug 

River from 1995-96 generally display the typical 

hydrologic cycle (fig. 9). Streamflow and recharge are 

high at the beginning of winter, because precipitation is 

abundant and plants are not using large amounts of 

water. As temperatures rise and plants begin to use 

water, streamflow and recharge decline and reach 

lowest levels in summer. After the first killing frost in 

fall, streamflow and recharge increase (Melvin, 1986). 

The decline in streamflow and recharge in spring 1995 

began earlier than normal because, as of June 1995, 

rainfall was 4.47 in. below normal for the year. As 

expected, less ground-water flowed to the stream in 

June to September during both years than in other 

months of the year, although streamflow and recharge 

were lower in summer 1995 than in summer 1996. 

Beginning in September 1995, a general rise in streamflow indicated increased ground-water discharge 

during this time period. Streamflow and recharge 

generally remained at this level until the summer of 

1996, when the yearly decline began. The aquifer test 

(October 30 to November 4, 1996) was conducted 

when the calculated ground-water component of 

streamflow was midway between the summer low and 

the rate at the end of the year.

1

10

100

1000

10000

J J MA M J J A SO N D J F MA M J J A SO ND

STREAMFLOW, IN CUBIC FEET 

PER SECOND 

0

2

4

6

8

10

12

14

MONTHLY GROUND-WATER 

RECHARGE, IN INCHES 

Streamflow at USGS station 01204000, Pomperaug River at Southbury

Ground-water recharge, estimated by application of

PULSE

RORA

PART

1995 1996

F

Figure 9. Streamflow and ground-water recharge for the Pomperaug River, 1995-96.
18 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

On October 20, 1996 (10 days prior to the aquifer 

test), 4.71 in. of precipitation fell in the Pomperaug 

River Basin (fig. 10). Stream stage in the Pomperaug 

River peaked on October 20 and subsequently declined 

throughout the aquifer test. An analysis of streamflow 

using the program PULSE shows that the 

2 months preceding the aquifer test was a period of 

generally increasing ground-water contribution to 

streamflow (fig. 10). From the lowest streamflow in 

September 1996 until the aquifer test, there were 5.1 in. 

of recharge, according to the program PULSE. Most of 

this recharge (3.1 in.) took place during the large rainfall 10 days before the aquifer test.

Ground-water recharge calculations presented in 

the previous paragraphs are corroborated by groundwater levels in USGS network wells (fig. 11). Water 

levels are measured by the USGS at four wells in the 

Pomperaug River watershed (biweekly October 1991 

to October 1996; monthly after October 1996) (figs. 1 

and 11; table 5). The wells are at various positions 

across the river valley—the valley bottom (SB 39), a 

stream terrace in stratified glacial deposits (SB 30), the 

side of a thick till deposit (SB 41), and on top of a thick 

till deposit (SB 42). The valley bottom well has a 

generally constant water level that is controlled by the 

stage of the Pomperaug River. Water levels in the other 

wells show a dampened fluctuation similar to the 

streamflow hydrograph. Water levels generally 

declined throughout 1995 until October when recharge 

increased. Rises in ground-water levels in the fall of 

1995 appear to lag behind calculated recharge events 

(fig. 11). By November, recharge was sufficient to 

maintain water levels (that is, recharge equaled groundwater outflow), so only small changes in water level in 

the aquifer took place until the summer of 1996, when 

ground-water levels began to decline. A period of 

rising water levels, indicating that recharge exceeded 

ground-water outflow, began in September 1996. 

1

10

100

1000

10000

September

1

September

11

September

21

STREAMFLOW AND GROUND-WATER 

October 1 October 11 October 21 October 31

RECHARGE,

 IN CUBIC FEET PER SECOND

0

1

2

3

4

5

6

7

8

9

10

PRECIPITATION, IN INCHES

Precipitation at Burlington, Connecticut

Streamflow at USGS station 01204000, Pomperaug River at Southbury

Ground-water recharge, estimated by application of PULSE

Aquifer test

Figure 10. Streamflow, precipitation, and ground-water recharge for the Pomperaug River, 

September to October 1996.
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 19

Table 5. Ground-water levels in network wells in the Pomperaug River watershed, Connecticut

Well Topographic 

position

Period of 

record

Depth to water, in feet below land surface

Median September 23, 

1996

October 30, 

1996

November 27, 

1996

SB 39 Valley bottom 10/24/91-

10/29/98

6.62 7.08 4.96 4.83

SB 30 Terrace 1/2/79-

10/29/98

18.97 20.05 16.9 16.93

SB 41 Hill side 10/24/91-

10/29/98

47.29 49.98 45.65 46.05

SB 42 Hill top 8/19/93-

10/29/98

13.60 16.30 12.1 12.3

F

0

10

20

30

40

50

60

70

80

J J MA M J J A SO ND J F MA M J J A SO ND

WATER LEVEL, IN FEET BELOW LAND SURFACE

0

2

4

6

8

10

12

14

16

18

20

MONTHLY GROUND-WATER RECHARGE, IN INCHES

Network wells

Valley bottom (SB39)

Terrace (SB30)

Hill side (SB 41)

Hill top (SB 42)

Ground-water recharge, estimated by application of

PULSE

RORA

PART

1995 1996

F

Figure 11. Water levels in four network wells and ground-water recharge in the 

Pomperaug River watershed, Connecticut, 1995-96.

F
20 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Ground-water levels in the four USGS network 

wells generally were close to median annual levels 

during the fall of 1996 (table 5). Water levels were 

below the median annual levels prior to the aquifer test 

on September 23, 1996 and rose from September 23 to 

October 30, probably in response to the precipitation 

on October 20. Water levels on November 27 were 

similar to the previous measurement, indicating that 

recharge was about equal to or greater than discharge. 

The frequency of water-level measurements 

(biweekly) was not sufficient to show when water 

levels reached their annual high level.

In addition to the long-term data available for the 

Pomperaug River, some limited data are available from 

the STS site that can be used to assess hydrologic 

conditions at the time of the aquifer test. Stream stage 

in Transylvania Brook declined throughout the aquifer 

test. In some situations, the amount of ground-water 

flow to or from the stream (referred to as streamflow 

gains and losses, respectively) can be calculated by 

subtracting the downstream flow from the upstream 

flow; in this situation, however, the difference between 

the streamflows was less than the standard deviation of 

the gain or loss (table 2). Water levels in two piezometers driven about 2 ft into the streambed near the 

pumped well were compared to water levels in the 

stream to see if water was leaking into or out of the 

stream. Throughout the aquifer test, water levels in the 

aquifer below the streambed were higher than water 

levels in the stream, indicating that the stream was 

gaining water throughout the test.This relation 

persisted as the stream stage declined throughout the 

test.

Water levels were recorded for several days prior 

to the aquifer test in observation wells used for the test; 

one water-level hydrograph is shown in figure 12. 

Large fluctuations caused by the cycling on and off of 

the pump in well PW-3 may obscure, to some extent, 

the natural fluctuations that would take place because 

of recharge; however, the high and low extremes in 

each pumping cycle increased before and after the 

aquifer test, indicating that ground-water recharge may 

have been taking place, as in the USGS network wells 

(table 5). 

Figure 12. Water levels in observation well E2D, Southbury Training School, Connecticut.

150

152

154

156

158

160

162

164

166

168

170

Oct 26 Oct 27 Oct 28 Oct 29 Oct 30

WATER LEVEL, 

IN FEET ABOVE SEA LEVE

L

Approximate trend line
GEOHYDROLOGY OF TRANSYLVANIA BROOK VALLEY 21

The hydrologic conditions during the aquifer test 

were complex. Ground water and streamflow relations 

during the aquifer test were typical of fall, but this relation was complicated by the 4.71 in. of precipitation on 

October 20. On the basis of existing data, however, it is 

reasonable to state that ground-water recharge was 

taking place, but at a rate below normal for fall and 

winter. To use water levels and streamflow as observations in model calibration, both sets of observations 

must reflect the same set of hydrologic conditions; in 

this case, however, ground-water levels were probably 

rising because of precipitation 10 days before the 

aquifer test whereas stream stage and streamflow were 

declining. Model calibration using these data should 

include the recharge on October 20, or the data should 

be adjusted before simulation to account for the 

recharge. Although the stream stage was declining 

throughout the aquifer test, ground-water levels may 

still have been rising because recharge is delayed by the 

time of travel through the unsaturated upper part of the 

aquifer system and from upland areas.

Water use

The use of water in the watershed falls into two 

categories—domestic water use, consisting of water 

withdrawn from small, privately owned wells 

completed in bedrock, and institutional water use, 

consisting of water withdrawn from wells completed in 

stratified glacial deposits. Residents of the STS are 

served by three wells identified as PW-1, PW-2, and 

PW-3 (table 6; fig. 2). Well PW-3 is the main production well and is in daily use at a rate of 0.338 Mgal/d 

(Al Van Geersdaele, Southbury Training School, oral 

commun., 1996). Well PW-1 is a standby well that 

automatically turns on and off during the day at a rate 

of 0.180 Mgal/d as demand dictates. Well PW-2 is classified as an emergency source of supply and is unused 

because of high concentrations of iron and manganese. 

Wells PW-1 and PW-2 were drilled in 1938 with 

reported yields at that time of 0.338 Mgal/d and 

0.72 Mgal/d, respectively. By 1967, the yields of these 

wells had declined to 0.144 Mgal/d for well PW-1 and 

0.41 Mgal/d for well PW-2. Well PW-3 was drilled in 

1970 and had a reported yield of 0.36 Mgal/d. Wateruse records were obtained from STS for January 1989 

to June 1996. The pumping rate during this period averaged 0.234 Mgal/d and ranged from 0.231 to 

0.263 Mgal/d, with slightly higher rates common 

during summer months. There was no discernible trend 

in water use (T.W. Frick, U.S.Geological Survey, 

written commun., 1996). Wastewater from the STS is 

collected at the sewage-treatment plant (fig. 2), treated, 

and discharged to Transylvania Brook downstream 

from USGS streamflow-gaging station 01204350.

A water-use survey was conducted by the USGS 

to determine the magnitude of domestic ground-water 

withdrawals in the Transylvania Brook watershed and 

what effect, if any, the withdrawals might have on 

ground-water flow. Annual domestic water use was 

estimated to be 1.1 Mgal during June to August 1995 

and 3.0 Mgal during September 1995 to May 1996 

(T.W. Frick, U.S.Geological Survey, written commun., 

1996). Domestic water is withdrawn from bedrock and 

returned through septic systems to the stratified 

deposits. About 15 percent of the domestic water is 

estimated to be lost through evapotranspiration. The 

net ground-water withdrawal from the aquifer system 

is about 0.0017 Mgal/d. This rate is small relative to 

pumping at the STS (0.338 Mgal/d) and is distributed 

over a large area; therefore, domestic ground-water 

withdrawals were not simulated in models in this study.

Table 6. Construction details of water-supply wells, Southbury Training School, Connecticut

[Data from D.L. Mazzaferro, Ground Water Inc., written commun., 1991; —, information not available]

Well Date 

drilled

Depth,

in feet

Diameter,

in inches

Water level, 

in feet below 

land surface 

(reported by 

driller)

Screen

length,

in feet

Screened 

interval,

in feet

Diversion 

permit, in 

million gallons 

per day

PW-1 1938 53 12 — 13 — 0.144

PW-2 1938 43 16 5 10 — .288

PW-3 1970 79 12 10.4 10 66-76 .36
22 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

DELINEATION OF CONTRIBUTING AREAS

Ground-water flow was simulated at two scales, 

a watershed scale and an aquifer scale. The watershedscale model, which encompasses upland and valley 

areas in the Transylvania Brook watershed, was used to 

understand how the topography, precipitation, and 

geology of the watershed affect the distribution of 

recharge to the stratified glacial deposits in the valley. 

The “distribution of recharge” in this report refers to 

the contribution of recharge from each upland basin, as 

a percentage of the total rate of recharge. Contributing 

areas were not estimated using the watershed-scale 

model. The aquifer-scale model encompassed the stratified glacial deposits in the valley and used the distribution of recharge estimated by the watershed-scale 

model. Contributing areas and the magnitude of 

recharge were estimated in the aquifer-scale model. 

The three-dimensional, finite-difference ground-water 

flow computer code known as MODFLOW-96 was 

used for both the simulations (McDonald and 

Harbaugh, 1988; Harbaugh and McDonald, 1996).

Watershed-scale simulation

The watershed-scale model was used to simulate 

processes in upland areas that affect the rates and distribution of ground-water recharge. The distribution of 

recharge indicated by this model was used in the 

aquifer-scale model. The simulation at the watershed 

scale is very coarse and should not be used to predict 

water levels in specific wells.

The Variable-Recharge package (Kontis, in 

press) in MODFLOW was used for the simulation. In 

this package, upland drainage basins (fig. 13) are simulated by specifying the land-surface altitude and the 

amount of water available for recharge (mean annual 

runoff). The model then calculates the head in the 

aquifer. If the head is above land surface, water is 

rerouted as channeled or unchanneled flow into the 

valley aquifer according to a user-specified ratio. 

Streamflow processes that might be simulated with the 

River or Stream packages in MODFLOW are not used 

in upland areas. The stratified glacial deposits in the 

valley receive water from upland areas and from direct 

recharge. As in upland areas, if the head is above land 

surface, water is rerouted. Streamflow losses from tributary streams in the valley bottom are calculated as in 

the Stream package (Prudic, 1989). The VariableRecharge Package can be used to estimate the amount 

of mean annual runoff that becomes ground-water 

recharge and the distribution, by drainage basin, of that 

recharge.

Model grid and layers

Ground-water processes were simulated at the 

watershed scale using a one-layer model (fig. 13). The 

bottom elevation of the aquifer was generated by 

subtracting 80 ft from land-surface altitude in the 

valley areas and 160 ft in the upland areas. The model 

grid has cells with a uniform spacing of 500 ft on each 

side. The grid is 65 rows by 27 columns and encompasses an active grid area of 5.03 mi2. The upland area 

is represented by 463 cells covering an area of 

4.15 mi2, and the valley is represented by 98 cells 

covering 0.88 mi2

.

Hydraulic properties of the aquifer

Simulated values for aquifer properties were 

based on previous studies (tables 3 and 4). Hydraulic 

conductivity of till and (or) bedrock, fine-grained stratified deposits, and coarse-grained stratified deposits 

was 0.5, 5, and 80 ft/d, respectively. The vertical 

hydraulic conductivity divided by the thickness of the 

riverbed deposits was 1 d-1. 

Boundary conditions and stresses

Water levels in shallow water-table aquifers 

often are a subdued replica of the land surface. No 

evidence to the contrary is present at the STS; therefore, the lateral boundaries of the aquifer system are 

assumed to be the watershed of Transylvania Brook 

(fig. 13) and are treated in the model as no-flow boundaries. In general, some differences could be present 

between the surface-water and ground-water divides, 

especially where there are large changes in hydraulic 

properties of the aquifer or large ground-water withdrawals near the boundary of the aquifer. The only 

known place in the study area where either of these 

conditions apply is along the possible fault between the 

arkose and crystalline bedrock. The fault could be a 

zone of hydraulic conductivity that is higher than the 

unfaulted bedrock that surrounds it, or it could be 

sealed shut by mineral deposits in the fault, by the 

smearing of rock debris into the fault openings, or by 

dense till forced into the top of the fault during glaciation. The fault was not simulated in the watershed-scale 

simulation, but the effects of the fault were tested in the 

aquifer-scale simulation.
DELINEATION OF CONTRIBUTING AREAS 23

Figure 13. Extent of watershed-scale model grid, Transylvania Brook watershed, Connecticut.
24 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

The lower boundary of the ground-water-flow 

system is formed by the imaginary surface that separates the water that discharges to Transylvania Brook 

from the water that flows beneath the watershed and 

discharges elsewhere or is relatively stagnant. Under 

this definition, the ground-water-flow system includes 

some water in the bedrock and all water in the glacial 

deposits. Some water in the bedrock may discharge to 

surface water in the watershed, but some water could 

flow out of the watershed through deep fractures or 

regionally extensive faults in the bedrock. No evidence 

indicates any significant quantities of water leaving the 

watershed through subsurface fractures or faults. 

The upper boundary of the ground-water-flow 

system is the water table, which is in coarse-grained 

deposits, fine-grained deposits, till, or bedrock. The 

water table fluctuates up and down as the amount of 

ground water in storage changes. 

The water available for recharge was applied 

over the entire modeled area. The watershed-scale 

model was intended to simulate natural conditions, and 

pumping was not simulated.

Model calibration

The watershed-scale model was checked for 

plausibility by comparing it to historical water levels 

recorded on well-completion reports and the amounts 

of recharge estimated using equation 1 for upland and 

valley areas. No changes were made to assigned values 

in the model except that mean annual runoff was 

lowered slightly from 8.60 ft3/s to 7.85 ft3/s to achieve 

a better match between simulated recharge rates and 

those predicted using the regression relation (equation 

1). This reduction was reasonable because the watershed-scale simulation considered only steady-state 

conditions; some of the mean annual runoff would be 

stormflow, a non-steady-state process not included in 

the model.

The watershed-scale model produced reasonable 

hydraulic heads within the limits of the data. Model 

heads were evaluated against measurements that were 

collected over a long period of time and a wide range 

of hydrologic conditions. Overall, different hydrologic 

conditions represented by the data probably average 

out, and the model and measured values are reasonably 

close. The mean error (observed head minus simulated 

head) was -1.16 ft. A value close to zero indicates that 

the model is unbiased and that the errors are random. In 

this case, the model is judged to be unbiased. The mean 

absolute error was 16.81 ft; this is the amount of deviation from the mean in the data. The measured heads in 

the watershed have a range of 349 ft from the lowest to 

highest measurement; therefore, the mean absolute 

error is about 5 percent of the total change in head. The 

root mean square error is 22.5 ft, or about 6 percent of 

the total change in head. This amount of error, relative 

to the total head loss in the system, was considered 

acceptable. 

Spatial bias also is important is assessing the 

plausibility of a model. A model should not consistently over- or under-predict water levels in any region 

of the modeled area. In the watershed-scale model, 

spatial bias was difficult to assess on the basis of 

measured water levels because the measurement locations are not evenly spread across the watershed (they 

are clustered around the valley), and the water levels 

were measured at different times of different years. 

Another way to assess spatial bias is to look at the 

predicted position of streams in the upland areas. The 

Variable-Recharge package uses the relation of landsurface and water-table altitudes to determine where a 

stream should be. In the watershed-scale model, 

streams were predicted in approximately the correct 

positions in basins 1, 2, 3, and 7 (fig. 13). Ponds near 

the headwater of the stream in basin 1 also were 

predicted by the model. These observations support the 

conclusion that predicted water levels were not consistently above or below the expected level near these 

streams.

The watershed-scale model yielded a net 

recharge rate of 8.8 in/yr to upland areas, slightly more 

than the 8 in/yr predicted by equation 1. Recharge from 

upland areas to stratified glacial deposits is equivalent 

to infiltration minus seepage losses to upland streams 

plus streamflow loss where upland tributaries cross 

over the coarse-grained deposits and unchanneled flow 

from upland basins. Net recharge from individual 

upland basins is given in table 7. 

The model yielded a net recharge rate of 

23.6 in/yr to the valley aquifer directly, which is higher 

than the 22 in/yr predicted by equation 1. Direct 

recharge from valley areas to stratified glacial deposits 

is equivalent to infiltration minus seepage losses to 

valley streams plus streamflow losses from the main 

channel of Transylvania Brook. Net direct recharge to 

individual basins in the valley is given in table 7.
DELINEATION OF CONTRIBUTING AREAS 25

Table 7. Recharge to stratified glacial deposits, Transylvania Brook drainage basin, Connecticut

Basin

Recharge from upland area Direct recharge to valley area Total recharge 

to valley, in 

cubic feet per 

second

in cubic feet 

per second

in inches per 

year

in cubic feet per 

second

in inches per 

year

1 0.18 1.7 0.05 23.6 0.23

2 .16 2.5 .03 23.6 .19

3 .50 20.9 .10 23.6 .60

4 .14 21.2 .04 23.6 .18

5 .31 11.7 .14 23.6 .45

6 .08 20.2 .10 23.6 .18

7 .27 12.2 .02 23.6 .29

8 .27 21.5 .15 23.6 .42

9 .32 21.1 .25 23.6 .57

10 .47 20.2 .02 23.6 .49

Results of the simulation

The main conclusion to be drawn from the watershed-scale simulation is that upland recharge is not 

evenly distributed (table 7). In the small, undrained 

basins, primarily those on the eastern side of Transylvania Brook (fig. 13), almost the entire amount of 

water available for recharge enters the aquifer and 

passes through the subsurface into the coarse-grained 

stratified deposits. In three basins on the western side 

of Transylvania Brook (fig. 13; basins 1, 2, and 7) that 

are drained by well-defined streams, most of the upland 

recharge enters the streams, which then flow out onto 

the valley aquifer. Where these streams flow over 

coarse-grained deposits, water flows from the stream 

into the aquifer and is counted as recharge. Tributary 

streams that do not have a large upland catchment area 

(basin 5 in fig. 13, for example) do not recharge the 

aquifer, rather they drain water from the aquifer.

Aquifer-scale simulation

The aquifer-scale model was used to calculate 

contributing areas to the two water-supply wells at the 

STS. Parameters of the model were estimated using 

nonlinear regression and drawdowns measured during 

an aquifer test. The distribution of recharge was 

obtained from the watershed-scale simulation, but the 

magnitude of recharge was estimated by the aquiferscale model to reflect hydrologic conditions at the time 

of the aquifer test. Contributing area simulations were 

done using hypothetical drought recharge conditions, 

which were based on information collected in the 

Pomperaug River Basin during the drought of record in 

the mid-1960s.

Ground-water-flow models typically are calibrated by a trial-and-error method, in which the model 

is adjusted by the modeler until a reasonable match 

between calculated and observed heads and flows is 

produced. Nonlinear regression makes calibration 

more efficient and objective because parameter values 

are adjusted to obtain automatically the best possible fit 

between simulated and observed values. The model fit 

is measured by the sum of squared weighted errors 

(SSE), where error equals the simulated minus 

observed values. The statistical framework of this 

process can be used to test the validity of the regression, the reliability of the parameter estimates, and the 

likelihood that a given model represents the system 

more accurately than an alternative model. The method 

is described by Cooley and Naff (1990), Hill (1992, 

1994, and 1998), and Poeter and Hill (1997). The 

computer program used is called UCODE, which 

stands for Universal Inverse Code (Poeter and Hill, 

1998). UCODE works with any model code by modifying model input files from user-defined templates, 

running the model, reading the model output, running 

the nonlinear regression, and modifying model parameters until the optimal parameter values are reached.
26 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

In this study, UCODE was used with 

MODFLOW-96. In UCODE/MODFLOW-96, the user 

defines parameters that represent the boundary conditions, stresses, and hydraulic properties of the aquifer 

system. Inputs to the transient model that apply to areas 

of the model grid, like aquifer properties and recharge, 

can be defined using zones (part of the model grid 

corresponding to a particular geologic unit). Model 

parameters also can be created using user-defined functions. For example, the vertical conductance between 

model layers, which is used as input by MODFLOW96, is calculated by dividing the parameter estimate by 

the distance between the centers of the overlying and 

underlying layers. 

Model grid, layers, stresses, and time steps

The aquifer-scale model grid was designed to 

simulate ground-water flow in the stratified glacial 

deposits and underlying bedrock in the Transylvania 

Brook watershed (fig. 14). The grid has 44 columns 

and 56 rows, variably spaced from about 10 ft near the 

pumped well (PW-3) to about 200 ft near the edges of 

the model grid. The active cells in the model cover 

about 1.459 x 107

 ft2, or about 0.5 mi2. 

The ground-water-flow system was simulated 

using three layers, although in some alternative 

models, described later in this report, bedrock was 

simulated as a fourth layer. All layers were simulated as 

convertible, meaning that MODFLOW-96 would 

determine if the layer was a water-table layer or a 

confined layer and then apply the appropriate hydraulic 

property values. In this study, all geologic units were 

simulated, so the bottom altitude of a given layer is the 

same as the top altitude of the underlying layer. Model 

layers are defined according to a percentage of the total 

thickness and to the geology, where the total thickness 

of the stratified deposits is the distance between land 

surface and the bedrock surface. Where the entire 

thickness of the deposits is the same geologic unit, 

layers 1 and 2 each represent one-quarter of the total 

thickness, and layer 3 is one-half of the total thickness. 

Where the deposits are in different geologic units, 

layers 1 and 2 are each one-half of the upper unit, and 

layer 3 is the lower unit. The top of layer 3 is shown on 

figure 15.

The stratified deposits are draped along the side 

of the bedrock valley wall (figs. 6 and 7), particularly 

on the western side of the valley. The water table probably is flatter than the slope of the bedrock surface 

along the lateral margins of the valley, and the saturated 

part of the aquifer gradually thins toward the edges of 

the model grid. This gradual thinning was incorporated 

into the model grid design by designating the minimum 

saturated combined thickness of layers 1 and 2 to be 

20 ft. The additional volume of aquifer material simulated in this way was offset by making model grid cells 

in layer 3 inactive where the total estimated saturated 

thickness was less than 40 ft. 

The pumped well (PW-3) was simulated as an 

aquifer stress using the Well package in MODFLOW96. The rate of pumping as measured during the test 

was used as the stress rate (185 gal/min), and the well 

was simulated in layer 3. Other wells at the STS were 

not in use during the test and therefore were not simulated. Domestic water wells also were not simulated, 

because the total volume of water pumped was small 

and spread over a large area. In addition, these wells 

were far from the pumped well at STS, thus minimizing 

any effect on the aquifer test.

The model consisted of 50 time steps that 

increase in size from 0.0001 to 0.8334 d by a factor of 

1.2. The total length of the stress period was 5 d. The 

initial conditions for the transient simulation were 

generated by a steady-state model. The parameter 

values in the model were estimated simultaneously for 

the steady-state and transient models.
DELINEATION OF CONTRIBUTING AREAS 27

Figure 14. Extent of aquifer-scale model grid and boundary conditions, Transylvania Brook valley, Connecticut.
28 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 15. Altitude of top of model layer 3, aquifer-scale model, Transylvania Brook valley, Connecticut.
DELINEATION OF CONTRIBUTING AREAS 29

Boundary conditions and model parameter 

definition

Nine model parameters were defined, although 

not all parameters were estimated in the final simulation. Six parameters define aquifer properties, and 

three define boundary conditions. The simultaneous 

estimation of all defined model parameters was not 

possible because the model was not sensitive to all 

parameters. If calibration data do not contain sufficient 

information about the parameter, the model is not 

sensitive to that parameter and the parameter cannot be 

reliably estimated. Parameters that cannot be estimated 

may nonetheless be important for model predictions. In 

this report, independent information on unestimated 

parameters is used to include their effect in model 

predictions.

The stratified deposits at the STS were divided 

into two zones related to the surficial geologic units 

(fig. 16). Coarse-grained deposits were designated 

zone 1, and fine-grained deposits were designated zone 

2. One hydraulic conductivity parameter is defined for 

zone 1 in layers 1, 2, and 3, and another for zone 2 in 

layers 1, 2, and 3. Layer 4 (bedrock), used only in an 

alternative model, also was divided into two zones, one 

representing bedrock (schist and arkose, zone 3), and 

one representing the possible fault contact between the 

two rock types (zone 4). 

The uppermost model layer contains the water 

table (the “water-table layer”), and the amount of water 

released from storage in the aquifer is proportional to 

the specific yield and the decline in water level during 

pumping. Although the water-table layer contains 

coarse and fine deposits, both were assigned the same 

specific yield, which was treated as a parameter. It is 

expected that the specific yield of the coarse and fine 

deposits represented by zones 1 and 2, respectively, are 

different, and the specific yield estimated by the model 

is more representative of zone 1 than zone 2. Treating 

them as the same is not expected to affect the model 

appreciably because drawdown was much less in zone 

2 than in zone 1. In layers below the water-table layer, 

water is not produced from de-watering of the aquifer 

as in the water-table layer. Instead, water is released 

from storage in the aquifer in proportion to the storage 

coefficient and the decline in pressure in the aquifer. 

One parameter is used for the specific storage in zones 

1 and 2. Specific storage is the storage coefficient 

divided by the thickness of the model layer, which 

UCODE multiplies by the layer thicknesses.

The stratification of the glacial deposits 

produces a difference in hydraulic conductivity in the 

vertical direction relative to the horizontal direction. A 

function is defined in UCODE to divide vertical 

hydraulic conductivity by the vertical distance between 

the center of layers to produce the vertical conductance 

term required by MODFLOW-96. The vertical 

hydraulic conductivity was represented by two parameters, one for each of zones 1 and 2. If vertically adjacent model cells contain different zones, the vertical 

hydraulic conductivity of zone 2 is used in the calibration. The effective hydraulic conductivity is the 

harmonic mean of the two, which is weighted heavily 

toward the lower of the two values (zone 2). 

Boundary conditions govern the flow of water 

into and out of the model grid during the simulation 

(fig. 14). The amount of this flow is governed by 

boundary condition terms that can be defined as model 

parameters. Two types of boundary conditions were 

used in the model, specified flow and head-dependent 

flow. Specified-flow boundaries were used to simulate 

the base and lateral boundaries of the model (flux = 0) 

and recharge. Recharge was applied using the Recharge 

package in MODFLOW-96. Simulated recharge 

included the combined ground-water recharge from 

upland basins and valley bottom (table 7), applied to 

the corresponding area in the valley bottom (fig. 14). A 

parameter was defined that proportionally changed this 

recharge. A parameter value of 1.0 produced the same 

amount of recharge as the watershed-scale simulation. 

Recharge was applied only to zone 1 areas (coarsegrained deposits) in the uppermost active model layer. 

Recharge was not applied to zone 2 because of the 

expected low vertical hydraulic conductivity of this 

unit, and because most of the recharge would be 

expected to occur near the valley wall, which corresponds to the location of zone 1 near the STS. This 

approach differs from that commonly used in modeling 

of valley aquifers in New England in which upland 

recharge is added through imaginary wells at the edge 

of the valley (Mazzaferro, 1986). 
30 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 16. Aquifer zones, aquifer-scale model, Transylvania Brook valley, Connecticut.
DELINEATION OF CONTRIBUTING AREAS 31

In this study, the parameter estimation procedure 

could change the simulated location of the edge of the 

valley, which made the use of imaginary wells problematic. Zone 1 constitutes a narrow band of cells 

(areas 1, 2, 3, 4, and 6 on fig. 14) in the area near the 

well field that is similar to the one-cell-width band that 

would receive recharge if imaginary wells were used. 

The approach used in this study spreads the source of 

recharge over a larger area than would be the case if 

wells were used, and thus slightly overestimates the 

contributing area.

Head-dependent flow boundaries are used to 

simulate flux to and from Transylvania Brook and 

Stibbs Lake. Transylvania Brook was simulated using 

the Stream Package (Prudic, 1989), and Stibbs Lake 

was simulated using the General-Head Package in 

MODFLOW-96. A head-dependent flow boundary 

requires the specification of a conductance, which 

controls the flow of water to (or from) an external 

source of water; flow equals the conductance times the 

difference in head between the boundary and the calculated value. In the case of Transylvania Brook and 

Stibbs Lake, the external hydraulic head is the average 

water level in the water body. The Stream Package also 

requires the altitudes of the top and bottom of the streambed, which were estimated as 1 and 2 ft below the 

average water level in the stream, respectively. The 

Stream Package accounts for water as it moves from 

one model cell to the next in the downstream direction 

and allows the stream to go dry if the water flowing out 

of the stream exceeds the amount of water flowing into 

the stream. 

In this implementation of UCODE and 

MODFLOW-96, parameters are defined for hydraulic 

conductivity of the stream and lake bed. For the stream 

boundary, a function is defined in UCODE that multiplies the hydraulic conductivity of the streambed by the 

area of the stream in each cell, divided by the thickness 

of the streambed. The stream areas were estimated 

using streamflow-measurement notes. For the lake 

boundary, the UCODE function multiplies the 

hydraulic conductivity of the lake bed by the area of the 

lake in each cell, determined using the 1:24,000-scale 

topographic map of the area. The thickness of the lake 

bed was assumed to be 1.0 ft. 

Model calibration

The model was calibrated to (1) steady-state 

water-level measurements made on October 30, 1996; 

(2) drawdown measurements made during the aquifer 

test from October 30-November 4, 1996; and (3) waterlevel rises after the rainfall on October 20, 1996. The 

calibration data consist of a set of measurements in 

each of 5 wells—one steady-state, pre-test water-level 

altitude in each well and 11 transient drawdown 

measurements in each well, for a total of 60 measurements. The 11 transient measurements that best represented drawdown in the aquifer were selected from 

many measurements made during the test. Streamflow 

measurements were not used in the calibration of the 

model because stream gains and (or) losses were within 

the margin of error for streamflow measurements.

Water-level data collected during the aquifer test 

were affected by factors external to the test that had to 

be removed from the data to calculate drawdown 

caused by pumping. External effects include (1) longterm water-level changes caused by years of pumping 

at the STS, (2) rise in water levels caused by the annual 

period of ground-water recharge that often takes place 

in Connecticut in the fall, and (3) rise in water levels 

caused by the pump being turned off for about 2 days 

prior to the start of the test. 

The effects of long-term water-level changes 

were assumed to be independent of time (in other 

words, average annual water levels are lower than they 

would be naturally, but have stabilized over decades of 

pumping at the STS); therefore, the long-term waterlevel correction depends only on distance from the 

pumped well. The radius of the area affected by 

pumping was calculated by assuming a circular 

recharge area around the well and a recharge rate of 

22 in/yr. The steady-state ground-water flow equation 

(Fetter, 1994, p. 218) was used to correct the steadystate water-level measurements made prior to the 

aquifer test by assuming no drawdown at the radius of 

the area affected by pumping. A hydraulic conductivity 

of 80 ft/d was used in this analysis, as determined by 

preliminary analysis of the aquifer test. This correction 

resulted in the addition of 0 to 3.50 ft to each of the 

steady-state water-level altitudes.

In Connecticut, water levels tend to rise from 

October to May (fig. 11); however, water levels in 1996 

probably rose faster than normal because of the 

extreme precipitation event on October 20 (fig. 10). 

This rise in water levels was treated as a linear trend 

defined by the increase in peaks in water levels 

measured several days before the aquifer test (fig. 12). 

A linear trend was fit to the rise of the lowest water 
32 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

level each day for each of the wells in zone 1 and was 

assumed to be steady during the aquifer test. The slope 

of each trend line was multiplied by the time, in days, 

since the beginning of the aquifer test and was added to 

drawdown. The trend correction resulted in a larger 

drawdown in all wells in zone 1 by 0.3 to 0.8 ft over the 

length of the test.

The effects of water-level recovery after 

pumping stopped, prior to the aquifer test, is shown 

following the last pumping cycle on figure 12. Ideally, 

the water levels should have reached a stable value by 

the time the aquifer test started, and corrections to 

drawdown would not be necessary; however, water 

levels were still rising at the start of the test. The nonsteady ground-water flow equation (the Theis equation; 

Fetter, 1994, p. 201) was used to correct the drawdowns 

for water-level recovery by assuming a transmissivity 

of 5,600 ft2

/d (from the preliminary aquifer-test analysis) and a storage coefficient of 0.0004, beginning 

when the pump was turned off. A constant pumping 

rate of 235 gal/min was used. The simulated waterlevel recovery at each time during the test, subtracted 

from the water level at the beginning of the test, was 

added to each measured drawdown. The result of this 

correction was an increase in drawdowns of 

0.09 ft over the duration of the test.

In UCODE, calibration data are assigned 

weights, thus allowing data that are known with a 

higher degree of accuracy to have a greater effect on 

the regression than data that are less well known. In this 

case, drawdown measurements made during the 

aquifer test (that help define the hydraulic properties of 

the aquifer and are accurately known) affect the regression more than the steady-state water-level altitudes 

(that help define the role of recharge and discharge in 

the aquifer and are less accurately known). According 

to regression theory (Hill, 1992, 1994, and 1998), the 

weights need to reflect possible error in the measurement of the data and are proportional to 1 divided by 

the variance of the measurement. Weights for measurements of water-level altitude are based on the assumption that 95 percent of the steady-state water level 

measurements are within 2.0 ft of their true value. This 

assumption is based on a qualitative assessment of how 

accurately the measured water level represents the 

horizontal and vertical heterogeneity of the aquifer 

within the entire model cell and on the estimated accuracy of the steady-state water-level correction. The 

measurement errors are assumed to be normally 

distributed, so that 2.0 equals the standard deviation 

(the square root of the variance) times the critical value 

at the 95-percent confidence interval (1.96). The variance for water-level altitudes is about 1.00. 

Water levels measured during the aquifer test 

were converted to drawdowns (the initial water level 

minus the water level at some time during the test). 

This subtraction cancels out many sources of error in 

the water-level altitude measurements. The weights for 

the drawdown measurements are based on the assumption that 95 percent of the measured drawdowns are 

within 0.1 ft of their true value. The actual measurements were made with a pressure transducer having an 

accuracy stated by the manufacturer of 0.03 ft; 

however, this accuracy may not have been realized 

under actual field conditions. The variance was calculated as above and equals 0.0026 for drawdown 

measurements.

Description of alternative models

Ground-water models are nonunique because 

hydrologic features in the flow system can be simulated 

in many ways. To evaluate the effect of this nonuniqueness, seven alternative models were posed in this study 

by simulating possible hydrologic features in various 

combinations. The alternative models are compared 

using various statistical measures of model fit and 

regression performance and, all else being equal, 

simpler models are given preference. In this study, one 

model was identified as best representing the flow 

system. The choice of alternative models is affected by 

the amount and type of data used for model evaluation. 

Some hydrologic features that might be important to 

simulating contributing areas may not be well represented in the model because the data used for model 

development did not contain enough information about 

those particular features. Composite-scaled sensitivities, which are computed by UCODE, can be used to 

evaluate if a model parameter can be estimated with the 

model. Composite-scaled sensitivities less than about 

0.01 times the largest composite-scaled sensitivity 

indicate that UCODE may not be able to estimate the 

parameter (Hill, 1998, p. 38). The alternative models 

are described below and are designated by the letters 

CAL (for “calibration”) followed by the numerical 

identifier of the simulation.

CAL0—The first alternative model was 

designed to include only the minimum number of 

hydrologic features needed to simulate the aquifer-test 

data. This model excludes Stibbs Lake (the hydraulic 
DELINEATION OF CONTRIBUTING AREAS 33

conductivity of the lake bed was set to 1x10-7 ft/d). 

This may be reasonable because the lake overlies finegrained deposits, which limit the lake as a source of 

water to the aquifer. The hydraulic conductivity of the 

streambed was set to 1.0 ft/d. Based on the compositescaled sensitivities computed at the optimal parameter 

values, the vertical hydraulic conductivity of zone 2 

could not be estimated and was fixed at 1x10-4 ft/d. Six 

parameters were estimated in the model—recharge, 

horizontal and vertical hydraulic conductivities of zone 

1, horizontal hydraulic conductivity of zone 2, specific 

yield, and specific storage.

CAL1—The second alternative model was identical to CAL0 except that the streambed hydraulic 

conductivity was fixed at 5.0 ft/d. A typical range for 

this parameter is 0.13 to 14 ft/d (Wilson and others, 

1974, p. 30); however, the model predicted significant 

negative recharge rates for streambed conductivities 

greater than 5.0 ft/d.

CAL2—This model was identical to CAL0 

except that Stibbs Lake was included by setting the 

hydraulic conductivity of the lake bed to 0.1 ft/d, based 

on values reported by Wilson and others (1974, p. 31). 

The sensitivities of this model at the optimal parameter 

values indicate that it would not be possible to estimate 

recharge while allowing leakage from the lake; therefore, recharge was fixed at the value estimated in 

CAL0.

CAL3—This also was identical to CAL0 except 

that a fourth layer was added to represent the bedrock. 

Previous alternative models were based on the assumption that the bedrock was impermeable; however, most 

domestic water wells in the area get their water from 

bedrock, so the bedrock is not completely impermeable. In this model, the hydraulic conductivity of the 

bedrock was set to 5 ft/d, based on table 3, and the 

thickness was assumed to be 200 ft. The specific 

storage of this layer was set to 1x10-5. 

CAL4—This model is identical to CAL3 except 

that the fault (fig. 4) in the bedrock is represented. The 

hydraulic conductivity of the fault was set to 200 ft/d to 

reflect the hypothesis that the fault might be a highly 

permeable zone within the bedrock. The thickness and 

specific storage of this layer are the same as for CAL3.

CAL5—This model is identical to CAL0, except 

that the vertical hydraulic conductivities of zone 2 and 

the lakebed were fixed at 0.1 ft/d. The model was 

designed to test the hypothesis that the vertical 

hydraulic conductivity of zone 2 is 0.1 ft/d (table 3). 

CAL6—This model was similar to CAL2, but 

bedrock was modeled as in CAL4. The model was 

designed to test the hypothesis that there might be 

interaction between Stibbs Lake and the fault. 

Results of the simulation

Alternative models CAL0, CAL1, and CAL5 

were very similar, based on their statistical properties 

(table 8), and all three fit the data more closely than the 

other alternatives. Model CAL1 was found to be sensitive to initial parameter estimates. In other words, 

different starting values resulted in different estimated 

parameters. In CAL1, recharge estimates tended to be 

lower than in CAL0 to offset the increased amount of 

water that could leak from the stream because of the 

higher streambed hydraulic conductivity. This was 

reflected in the model by a high correlation (0.94) 

between streambed hydraulic conductivity and 

recharge. In CAL1 and CAL5, the source of almost all 

water in the model is a few stream-boundary cells in the 

northern area of the model, and the optimal recharge 

rate was negative. This situation is considered unrealistic, and CAL0 was chosen as the most representative 

model.

Table 8. Statistical measures of model fit of alternative models

Statistical measure CAL0 CAL1 CAL2 CAL3 CAL4 CAL5 CAL6

Least squares objective function 278 264 343 313 482 278 460

Calculated error variance 5.15 4.89 6.23 5.80 8.93 5.15 8.52

Standard error of the regression 2.27 2.21 2.50 2.41 2.99 2.27 2.92

Correlation coefficient for observations .999 .999 .999 .999 .998 .999 .999

Correlation coefficient for normal 

residuals, R2

N

.955 .951 .983 .956 .974 .959 .974
34 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

The range of estimated parameter values among 

the alternative models was small (table 9). Estimated 

values generally were within the expected ranges based 

on typical values (table 3); however, recharge and 

specific-yield estimates were somewhat different than 

expected. The multiplication value for recharge was 

significantly smaller (0.24) than the value estimated by 

the watershed-scale model (1.0). Taken together with 

stream losses, this represents about 32 percent of the 

mean annual recharge estimated in the watershed-scale 

model. This is probably because the aquifer test was 

conducted at a time when ground-water recharge was 

just beginning for the year, and the full amount of 

annual recharge had not yet taken place. Also, there are 

some areas where coarse deposits overlie fine deposits 

(see, for example, the sand unit on the eastern edge of 

fig. 6). These deposits receive recharge but are not part 

of the coarse-grained aquifer and so were not simulated. This means that recharge estimated in the 

aquifer-scale model should be less than in the watershed-scale model.

Specific-yield estimates also were lower than 

expected. This could be because of the layered nature 

of the aquifer, which causes a delay in the time between 

when the head is lowered in the aquifer and when the 

water can physically drain downward through the 

various layers to the water table. This phenomenon has 

been noted in a study of a glacial aquifer with smallscale bedding in which short-term aquifer tests give 

specific yields that are much lower than values determined in laboratory studies (Nwankwor and others 

1984; 1992). In a review of Nwankwor’s studies, 

Moench (1994) found that a drainage delay can affect 

specific yield estimates, particularly in water-table 

wells and that unrealistically small specific yields can 

result from aquifer heterogeneity. 

Drawdowns in wells in deep, coarse-grained 

deposits were accurately simulated by model CAL0 

(fig. 17). Of particular significance is the fact that the 

simulated water levels followed the “s” shape of the 

data that is typical of water-table aquifers. This shape 

is caused by a change in the source of water from 

storage in the aquifer; the early part of the aquifer is 

dominated by specific storage and the later part is 

dominated by specific yield. The sensitivity of the 

model to each drawdown measurement shows the 

sensitivity to specific storage peaks early in the simulation and sensitivity to specific yield peaks at the end of 

the simulation period. 

Table 9. Optimal parameter values for alternative models

[Kh, horizontal hydraulic conductivity; Kv, vertical hydraulic conductivity; shading indicates parameter not estimated; —, not simulated]

Model parameter CAL0 CAL1 CAL2 CAL3 CAL4 CAL5 CAL6

Kh, zone 1 (ft/d) 154 155 157 148 132 165 133

Kh, zone 2 (ft/d) 29 26 21 18 12 13 11

Kv, zone 1 (ft/d) .83 0.87 .94 .66 .60 .73 .62

Kv, zone 2 (ft/d) 1.0E-04 1.0E-04 1.0E-04 1.0E-04 1.0E-04 1.0E-01 1.0E-04

Specific yield, unitless .0069 .0073 .006 .0082 0.010 .0048 .0095

Specific storage (/ft) 1.6E-05 1.6E-05 1.6E-05 1.6E-05 1.5E-05 1.6E-05 1.5E-05

Lakebed Kv (ft/d) 1.0E-07 1.0E-07 1.0E-01 1.0E-07 1.0E-07 1.0E-07 1.0E-01

Streambed Kv (ft/d) 1.0 5.0 1.0 1.0 1.0 1.0 1.0

Recharge multiplier, 

unitless

.23 -.002 .24 .23 .20 -.007 -.095

Bedrock Kh (ft/d) ——— 5 5 — 5

Fault Kh (ft/d) — — — — 200 — 200
DELINEATION OF CONTRIBUTING AREAS 35

0.1

1

10

0.001 0.01 0.1 1 10

Simulated

Observed

Observation well E2D

0.1

1

10

0.001 0.01 0.1 1 10

Simulated

Observed

0.1

1

10

0.001 0.01 0.1 1 10

Simulated

Observed

0.1

1

10

0.001 0.01 0.1 1 10

Simulated

Observed

0.1

1

10

0.001 0.01 0.1 1 10

Simulated

Observed

Observation well N2D

Observation well W1D Observation well N3

Observation well E1D

DRAWDOWN, IN FEET

TIME, IN DAYS

TIME, IN DAYS

DRAWDOWN, IN FEET

DRAWDOWN, IN FEET

TIME, IN DAYS

Figure 17. Simulated and observed water-level changes in coarse-grained deposits during the aquifer 

test, October 30 to November 4, 1996, Southbury Training School, Connecticut.

TIME, IN DAYS TIME, IN DAYS

DRAWDOWN, IN FEET

DRAWDOWN, IN FEET
36 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Steady-state water levels were accurately simulated in most, but not all wells (fig. 18). Observed water 

levels in coarse-grained deposits in two wells (W3 and 

N5) were significantly above simulated values. There 

are plausible explanations why these two wells were 

not simulated accurately. Well W3 is at the very edge of 

the model grid. In this area, two important features 

affect the water level in the well, and either could cause 

higher water levels than those predicted by the model: 

(1) the coarse deposits outside the model area may not 

have drained completely following the rainfall on 

October 20, 1996, and (2) the abrupt steepening of the 

bedrock valley in this area causes high vertical gradients over a short horizontal distance. The other well, 

N5, is a flowing well. The water level measured in that 

well was above the altitude of the top of the confining 

layer as it was mapped. The model could only simulate 

the flowing condition if the extent of the confining 

layer were increased to higher altitudes. The extent of 

the confining layer was not changed because no boring 

or other geologic data exist to support such an extension. The extent of the confining layer is very important 

in calculating contributing areas; however, an underestimate of the extent, as in this case, exposes more 

coarse-grained material at the land surface and probably slightly overestimates the contributing area.

Drawdowns in the fine-grained deposits were 

not as well simulated as levels in coarse-grained 

deposits, perhaps because of the relation among the 

layers in those deposits, the length of the well screen, 

and the thickness of the model layer. The fine-grained 

deposits are composed of many thin individual layers; 

this results in high vertical hydraulic gradients. This 

may cause observed water levels in wells with screens 

much shorter than the model layer thickness to be 

different from simulated water levels.

The simulated water-table map shows the 

expected pattern of ground-water flow. The shape of 

the 170-ft contour (fig.19) shows that streamflow 

recharges the aquifer in the northern part of the 

modeled area where it first comes in contact with the 

coarse-grained deposits. The 160-ft contour (fig.19) 

shows that flow is generally down valley toward the 

lower reaches of Transylvania Brook. At the end of the 

aquifer test, the 160-ft contour shows that streamflow 

was contributing flow to the aquifer. The 150-ft 

contour (fig. 19) shows little change caused by the 

aquifer test and that the stream is losing water to the 

aquifer. South of the area shown on figure 19, the 

water-table gradient is steeper where the brook is 

underlain by fine-grained deposits than where it is 

underlain by coarse-grained deposits. Ground water 

generally flows down valley until an area is reached 

where the brook is underlain by coarse-grained 

deposits. In these areas, ground water discharges to the 

stream at greater rates than where the stream is underlain by fine-grained deposits. 

Figure 18. Simulated and observed steady-state water levels, October 30 to November 4, 1996, 

Southbury Training School, Connecticut. [Error bars show an interval of plus or minus 2 feet.]

140

145

150

155

160

165

170

175

180

185

E1D E2D E2S E3 E4 N1S N2D N3 N5 W1D W1S W3

Observed

Simulated

STEADY-STATE WATER LEVEL,IN FEET ABOVE SEA LEVEL
DELINEATION OF CONTRIBUTING AREAS 37

Figure 19. Altitude of simulated water table before and after the aquifer test, October 30 to November 4, 1996, 

Southbury Training School, Connecticut.
38 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

The change from a gaining stream to a losing 

stream during the aquifer test is consistent with the 

streamflow measurements made before and after the 

aquifer test (table 3). The model simulates a gain in 

streamflow before the aquifer test from the upstream 

USGS streamflow-gaging station (01204340 on fig. 3) 

to the downstream USGS station (01204350 on fig. 3) 

of 0.02 ft3

/s. The measured rate of streamflow gain on 

this date (October 29, 1996) was 0.03 ft3

/s (table 2). At 

the end of the aquifer test, the model simulates a 

streamflow loss of 0.10 ft3

/s. The measured rate of 

streamflow loss on this date (November 4, 1996) was 

0.15 ft3/s. As mentioned previously, this relation was 

contradicted by data from the streambed piezometers, 

which indicated that the stream was gaining water 

throughout the aquifer test. This may be because (1) the 

stream really was gaining and the streamflow measurements are too imprecise to determine this small flow, or 

(2) the piezometers were not in places where the stream 

was losing water.

Model CAL0 was run and the results were 

compared to a second set of transient data to verify that 

the model is reasonable under other hydrologic conditions. Model response was compared to water levels 

collected in well E2D (fig. 20) for 4 days prior to the 

start of the aquifer test on October 30, 1996. These data 

show two hydrologic responses—an overall and 

gradual rise in response to the large amount of precipitation on October 20, 1996, and a daily cyclic fluctuation in response to the normal pumping at STS. To 

approximate this situation, a seasonal pattern of 

recharge was determined from the hydrograph-separation programs (fig. 9). The 140 days prior to October 

21, 1996 were a period of low ground-water recharge. 

Drought conditions (defined in next section) were 

assumed for this period. The beginning of the recharge 

period was simulated by assuming recharge at half the 

average annual rate. These conditions only approximate the true events, because the amount and distribution of recharge to the water table due to the storm are 

unknown. Although the simulated water levels were 

about 6 ft higher than the actual water levels, the model 

reproduces the rising trend and the cyclic fluctuations. 

Although this analysis is not quantitative, it supports 

model CAL0 as a reasonable representation of the flow 

system.

150

152

154

156

158

160

162

164

166

168

170

Oct 26 Oct 27 Oct 28 Oct 29 Oct 30

WATER LEVEL, 

IN FEET ABOVE SEA LEVE

L

Observed water level

Simulated water level

Figure 20. Observed and simulated water levels in observation well E2D, Southbury Training School, Connecticut.

Approximate trend line
DELINEATION OF CONTRIBUTING AREAS 39

Diagnostic and inferential statistics

The application of statistics to diagnose and 

analyze ground-water models is well documented 

(Cooley and Naff, 1990; Hill, 1992, 1994, 1998; 

Cooley, 1997; Poeter and Hill, 1997; Hill and others, 

1998). The discussion of model statistics in this report 

largely comes from that body of work and follows the 

methodology outlined by Hill (1998). The use of statistics with UCODE falls into two main categories: (1) the 

quality of the model calibration and (2) the accuracy of 

estimated parameters.

For the regression to produce a valid calibrated 

model, several conditions must be met. Assuming that 

the model is correct (the true geology and hydrology 

have been accurately represented by model parameters), the weighted residuals (the difference between 

the simulated and observed data divided by the variance) must come from a random distribution. Analysis 

of residuals, both numerical and graphical, can be used 

to determine how well the model is calibrated.

The calculated error variance (s2) and the standard error of the regression (the square root of the 

calculated error variance) (table 8), are quantitative 

indications of model fit. If the model fit is consistent 

with the data accuracy in the weighting matrix, these 

measures equal 1.0. Significant deviations of the calculated error variance or the standard error from 1.0 indicate that the fit is inconsistent with the weighting. 

Values of the calculated error variance and the standard 

error are typically greater than 1.0 in practice, 

reflecting the presence of model error as well as 

measurement error, which is not represented in the 

weight matrix.

The model fit also is indicated by the graphical 

relation of weighted residuals to weighted simulated 

values. Ideally, weighted residuals are scattered evenly 

about the line y=0, and their magnitude is not related to 

the simulated values. Plots were constructed for each 

alternative model, and weighted simulated residuals 

seem to be independent of weighted simulated values; 

therefore, the regression was judged to be valid 

according to this criterion.

Another measure of model calibration is that the 

observed values should be reasonably reproduced by 

the model, as reflected by the correlation coefficient 

between weighted simulated values and weighted 

observation that summarizes this relation. This correlation coefficient generally needs to be greater than 0.90 

(Hill, 1998). In all alternative models, the correlation 

coefficient is greater than 0.99 (table 8); therefore, the 

models are an adequate fit to the data by this criterion.

For a valid regression, the weighted residuals 

(simulated minus observed values, times the weight) 

need to be random and uncorrelated or correlated in a 

way that can be explained by the fitting imposed by the 

regression. If the weighted residuals are random, independent, and normally distributed, they fall on an 

approximately straight line in a normal probability 

graph. The summary statistic R2

N (table 8) is the correlation coefficient of this line and can be used to test the 

weighted residuals (Hill, 1992). For the number of 

observations in this model (60), the critical value for 

this statistic is 0.962 at the 0.05-percent significance 

level. Most alternative models have values close to or 

above the critical value (table 8), indicating that the 

residuals are nearly normal. The variations present are 

not considered to be important.

The estimated variance/covariance and correlation matrices (table 10) produced by UCODE can be 

used to measure the precision and correlation of the 

parameter estimates. The variance/covariance matrix is 

based on the optimal parameter values and includes the 

effects of parameters that are not estimated, such as 

lakebed hydraulic conductivity (Hill, 1998). In this 

model, streambed hydraulic conductivity and recharge 

are highly correlated (table 10), indicating that change 

in one can be offset by a change in the other to produce 

an identical model (correlation coefficient = 0.94). 

Uniqueness of a regression problem can be tested by 

varying the initial estimates of the parameters. If the 

regression converges to the same values, the estimated 

parameter values are likely to be unique. If the regression converges to different values, the parameters are 

too correlated and are not unique. It is possible that 

additional pairs of parameters are correlated, but that 

the correlation is obscured by inaccuracies in the sensitivities calculated by UCODE (Poeter and Hill, 1998). 

In this study, numerous initial estimates were tried in 

CAL0, and the final parameter values were judged to 

be unique.
40 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Table 10. Correlation and variance/covariance matrices for alternative model CAL0

[Correlations are in shaded part of table. Both matrices are symmetric, so the upper or lower diagonal of either matrix is the transpose of the part of the matrix 

shown. Correlations of parameters with themselves are equal to 1.0; italic type indicates parameter was fixed and not estimated. Kh, horizontal hydraulic 

conductivity; Kv, vertical hydraulic conductivity

]

Parameter

number

Lakebed

Kv

(1)

Streambed

Kv

(2)

Recharge

(3)

Specific

yield

(4)

Specific

storage

(5)

Kh,

zone 1

(6)

Kh,

zone 2

(7)

Kv,

zone 1

(8)

Kv,

zone 2

(9)

(1) 5.51E+0

0

0.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00 0.00E+00

(2) 0.00 3.42E+00 -7.48E-01 1.69E-01 -1.73E-02 4.13E-02 -4.67E-01 1.03E-01 -4.43E+00

(3) 0.00 -.94 .17 -5.05E-02 2.87E-03 -9.66E-03 1.10E-01 -2.16E-02 1.01E+00

(4) 0.00 .49 -.69 .04 1.75E-03 3.89E-03 -4.84E-02 3.40E-03 -5.14E-01

(5) 0.00 -.22 .14 .15 0.00 -1.03E-03 1.11E-03 -2.07E-04 7.57E-03

(6) 0.00 .55 -.57 .52 -.49 0.00 -8.06E-03 -3.83E-04 5.21E-03

(7) 0.00 -.84 .88 -.83 .11 -.66 .09 -1.42E-02 7.23E-01

(8) 0.00 .56 -.48 .10 -.04 -.18 -.45 .01 -2.74E-01

(9) 0.00 -.17 .19 -.22 .11 -.46 .25 .11 1.21E+02

Ground-water models are characteristically 

nonlinear; that is, the calculated sensitivities (here, 

derivatives of simulated hydraulic heads and drawdowns with respect to estimated parameters) of the 

model are related to the parameter value that is being 

estimated. Thus, changes in recharge, for example, 

could produce large changes in model head at one value 

of recharge but could produce small changes in model 

head at other values of recharge. The nonlinearity of 

ground-water models poses potential problems for the 

analysis of confidence intervals. Some models are 

approximately linear near the optimal parameter 

values, and Beale’s measure (Hill, 1994) can be used to 

quantify the degree of nonlinearity in a particular 

model. For model CAL0, if Beale’s measure is greater 

than 0.27, the model is nonlinear; if Beale’s measure is 

less than 0.024, the model is effectively linear. Beale’s 

measure for CAL0 is 0.28; therefore, the model is 

nonlinear.

ANALYSIS OF UNCERTAINTY OF 

CONTRIBUTING AREAS

In this study, deterministic and probabilistic 

contributing areas were estimated using a particletracking procedure. Deterministic contributing areas 

are estimated using the calibrated model and are correct 

if the geology has been mapped correctly, model 

parameters correctly represent the geology, groundwater recharge and discharge are well understood, and 

the calibration data have no errors. The effect of model 

construction on contributing areas can be analyzed by 

comparing the contributing areas produced by each 

alternative.

Estimation of probabilistic contributing areas 

uses additional information to better define the uncertainty in the calculated contributing area. Probabilistic 

contributing areas can be calculated using the calibrated model, the variance/covariance matrix generated 

by the nonlinear regression parameter estimation, and 
ANALYSIS OF UNCERTAINTY OF CONTRIBUTING AREAS 41

application of a probabilistic simulation technique, 

such as the Monte Carlo technique. In the type of probabilistic modeling done in this study, the calibration 

data are not assumed to be free from errors. Errors in 

the calibration data include (1) the inaccuracy of field 

measurements; (2) data that do not accurately represent 

the modeled feature (for example, water levels in wells 

with very short screens may not accurately represent 

water levels in very thick model layers; (3) factors 

external to the phenomena being modeled (such as 

steady-state water levels affected by a large rainfall that 

is not simulated); and (4) subsurface variation in 

hydraulic properties that is smaller in extent than the 

modeled feature (for example, water levels in a well 

screened in a clay layer of small extent may not represent water levels in the coarse-grained aquifer 

surrounding it). The latter type of error in the calibration data alternatively could be attributed to model 

error, because the model could be refined to accommodate the error (Hill, 1998). 

The main difference between the two types of 

contributing areas is that the deterministic method 

generates only one piece of information—either an 

area is inside the contributing area, or it is not. The 

probabilistic model generates more information: each 

area has a probability, ranging from 0 to 100 percent, of 

being within the true contributing area, based on errors 

in the model calibration data. In this way, the probabilistic method incorporates the uncertainty that is 

inherent, but often unstated, in deterministic contributing areas and offers water managers more information on which to base decisions. Probabilistic 

contributing areas generally are larger than deterministic areas because they include areas that have a low 

probability of contributing water to a well. This information can be used to make decisions, such as to err on 

the side of caution and manage land use over a greater 

area, or to collect additional data that could reduce the 

uncertainty, and the size, of the contributing area.

Particle tracking commonly is used to delineate 

areas that contribute water to wells (Franke and others, 

1998). In this study, a particle-tracking computer 

program known as MODPATH was used (Pollack, 

1994). Hypothetical “particles” are placed in the simulated system and moved in accordance with the groundwater velocities calculated by the model. In this study, 

one particle was placed in the center of each model cell 

at the water table. The particles were tracked forward 

until they discharged from the aquifer system to a 

pumped well or to Transylvania Brook. Each cell 

containing a particle that eventually discharged to a cell 

in which a pumped well was simulated was considered 

to be in the contributing area of the well. Contributing 

areas were calculated using calibrated model CAL0, 

with modifications to pumped well stresses and 

recharge rates. 

Particle tracking also requires an estimate of 

aquifer porosity. In this study, a uniform porosity was 

assumed; the porosity estimate does not affect the size, 

shape, or location of the contributing area, only the 

time-of-travel calculation. In this study, a uniform 

porosity of 0.20 was used, based on the low end of a 

range of typical values for a sand and gravel aquifer 

(Walton, 1984, p. 19). Using a low value of porosity 

produces a shorter travel times than if a higher value of 

porosity was used.

Conditions used to estimate contributing 

areas

The contributing area simulations were done 

using the registered diversion rate of 100 gal/min at 

well PW-1 and 250 gal/min at well PW-3. These wells 

were simulated in layer 3, which corresponds to the 

location of their screened interval. Well PW-2 was not 

simulated because it is classified as an emergency well 

only. 

Recharge was modified from CAL0 to simulate 

contributing areas under drought conditions to provide 

a conservative (larger area) estimate of the contributing 

area to the well. The drought is defined as (1) the 

maximum historical period with no direct recharge to 

the valley (180 days; R.L. Melvin, oral commun., 

2000), and (2) recharge from upland sources that would 

be typical of historical drought conditions, as determined from long-term USGS ground-water and 

surface-water records. 

Long-term USGS network well WY-1 (fig. 1) 

was used for drought analysis because it (1) has a long 

period of record (since 1944); (2) is in the same physiographic setting as the STS (Southwest Hills; Melvin, 

1986); and (3) is in a topographic setting that could be 

considered representative of the aquifer (valley terrace; 

Melvin, 1986). The lowest recorded water level for this 

well was 30.81 ft below land surface on September 24, 

1986. The mean water level for this well during 

1944-96 is 25.39 ft below land surface. The drought 
42 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

criterion based on WY-1 was defined to be an average 

5.42-ft drop in water level in the aquifer from average 

annual conditions. Upland recharge was varied, using a 

multiplier that proportionally changed recharge to 

produce the target drop in water level; 80 percent of the 

average annual upland recharge produces an average 

water level drop of 5.79 ft, and 70 percent of the 

average annual upland recharge produces an average 

water level drop of 7.08 ft. 

Several streamflow measurements were made on 

Transylvania Brook during a drought in the mid-1960s 

(table 2). The lowest recorded streamflow at USGS 

station 01204350 (fig. 3) was 0.15 ft3

/s on September 

27, 1965. At no time during the simulation did streamflow losses exceed this flow; therefore, the amount of 

available streamflow did not limit recharge to the 

aquifer. On this date, the gain in streamflow between 

USGS streamflow-gaging stations 01204350 and 

01204400 was 0.74 ft3

/s (table 2). A recharge multiplier of 80 percent produced 1.18 ft3/s of streamflow; a 

multiplier of 70 percent produced 1.01 ft3

/s. In simulations of contributing areas, 70 percent was chosen as 

the multiplier for the rate of recharge from uplands 

because it best met the combined criteria of groundwater level drop and streamflow gains. 

The flow system at the end of the 180-day transient simulation was used to calculate contributing 

areas. This assumption fixes the velocity of the ground 

water at the end of the drought period, and particles 

were tracked through this flow system as if the system 

were at steady state. In all likelihood, the flow system 

would change because of recharge at or before the end 

of 180 days. By assuming that drought conditions 

persist until all particles reach the pumped well, which 

could take several years, the contributing areas are 

overestimated. 

Deterministic contributing area

The contributing area calculated using CAL0 

includes the northern part of the Transylvania Brook 

valley within the model boundary (fig. 21). The gravel 

terrace (fig. 21) was outside the model boundary, but it 

should be considered to be in the contributing area to 

the wells. This terrace may not be saturated to any 

significant depth for most of the year, but rainfall infiltrates this area rapidly and then recharges the aquifer. 

The model is not detailed enough to simulate this 

ground-water runoff, but any surface contaminant 

could be introduced into the aquifer from this area. 

Travel times within the contributing area ranged from 

64 (0.18 year) to 1,808 days (4.95 years). The areas 

with the shortest travel times (less than 1 year) are at 

the base of the gravel terrace (fig. 21; PW-3) and the 

northeastern corner of Stibbs Lake (PW-1). Travel 

from most of the rest of the area is less than about 2 

years, with only a few small areas having longer travel 

times of up to 5 years. 

Selected pathlines of ground-water flow to well 

PW-3 (fig. 22) show the three-dimensional nature of 

the ground-water flow system at the STS. Pathline A 

begins on the gravel terrace, descends to layer 2, and 

remains in this layer until it comes close to the pumped 

well in layer 3. Pathline B travels in layer 1 perpendicular to the water-table contours through mostly fine 

deposits. The low vertical hydraulic conductivity of 

this material prevents the pathline from entering layer 

2. When pathline B passes onto the coarse-grained 

deposits, it descends over a short distance through layer 

2 into layer 3 and changes direction toward the pumped 

well. Pathline C travels through layer 1, descends over 

a short distance through layer 2 into layer 3, and then 

changes direction toward the pumped well. Other pathlines on figure 22 show that there is little vertical travel 

in the fine-grained deposits. Pathlines in layer 1 are 

perpendicular to the water-table contours until they 

pass over coarse-grained deposits, where they enter the 

underlying layer.

Contributing areas from the alternative models 

were computed, but are not published here; however, 

the contributing areas from the alternative models were 

very similar. The main difference was that contributing 

areas for two of the models (CAL2 and CAL6) were 

smaller than the contributing area calculated using 

CAL0. The contributing areas calculated using CAL2 

and CAL6 were more horseshoe-shaped than CAL0 

and included only the coarse-grained deposits that ring 

the head of the valley. The addition of bedrock in CAL3 

and CAL4 did not produce significantly different 

contributing areas than CAL0. None of contributing 

areas produced by alternative models extend south to 

the sewage-treatment facility (fig. 21). The general 

similarity among contributing areas indicates that 

nonuniqueness in model construction does not greatly 

affect the predicted contributing areas.
ANALYSIS OF UNCERTAINTY OF CONTRIBUTING AREAS 43

Figure 21. Deterministic contributing areas to wells PW-1 and PW-3, Southbury Training School, Connecticut.
44 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 22. Selected pathlines to well PW-3, Southbury Training School, Connecticut.
ANALYSIS OF UNCERTAINTY OF CONTRIBUTING AREAS 45

Probabilistic contributing areas based on 

variations in the defined parameters

The uncertainty in the estimated parameters as 

determined by UCODE is propagated through the analysis of the contributing areas through use of a Monte 

Carlo analysis. The Monte Carlo simulation involves 

generating sets of parameter values that, taken as a 

group, have the same statistical properties as the calibrated model parameters. Using the regression-derived 

variance/covariance matrix preserves the uncertainty 

with which parameters are estimated and preserves the 

correlation among the parameters. The variance/covariance matrix generated by UCODE for all parameters, 

whether they were estimated or not, is used in this study 

to quantify the effects of parameter uncertainty on the 

location and shape of the contributing area.

The Monte Carlo simulation is based on the definition of a random normal variable,

b = zσ + μ, (2)

where

b is a vector of model parameter values,

z is a vector of normally distributed random numbers,

σ is the square root of the variance/covariance matrix, 

and 

μ is a vector of optimal parameter values.

The Monte Carlo analysis was conducted as 

described by the following steps. 

Step 1. Flow model—Estimate optimal parameters of the ground-water-flow model and the variance/covariance matrix using UCODE (already 

discussed). 

Step 2. Deterministic contributing area—

Determine contributing areas by conducting particle 

tracking using MODPATH (Pollack, 1994). Place one 

particle in each model cell. Track the particles forward 

in time until all the particles have discharged to Transylvania Brook or the pumped wells (PW-3 and PW-1). 

Flag each cell according to the discharge point of its 

particle. Combine the flagged cells by discharge locations. Plot contributing areas on a map (fig. 21) 

(already discussed). 

Step 3. Generate parameter sets— Compute 

the square root of the variance/covariance matrix by 

taking its Cholesky decomposition matrix (Press and 

others, 1986). Multiply the resulting matrix by vectors 

of normally distributed random numbers to generate a 

large number of sets of model parameter values. 

Step 4. Monte Carlo simulation— Run the 

flow model, substituting one of the generated parameter sets for the optimal parameters. Run the particletracking simulation using this flow model. Save the 

particle discharge locations in a file.

Step 5. Compile results—Sum the number of 

times, over all simulations, that a particle reaches the 

pumped well(s) and divide by the number of simulations; thus, each cell will be associated with the number 

of times, in percent, a particle traveled from the cell to 

the pumped well. 

In this application of the Monte Carlo method, 

10,000 parameter sets were generated. The variance of 

two of the parameters, streambed hydraulic conductivity and hydraulic conductivity of zone 2, were very 

high. This led to many of those 10,000 parameter sets 

having unrealistically extreme (both high and low) 

values for those two parameters. The 95-percent linear 

confidence intervals for streambed hydraulic conductivity is 0.02 to 51 ft/d and for vertical hydraulic 

conductivity of zone 2 is 10-6 to 1,620 ft/d. These 

parameter ranges were not considered realistic, so the 

parameter sets were first conditioned to remove unrealistic sets. Parameter sets were excluded in which the 

vertical hydraulic conductivity of either the streambed 

or zone 2 was more than a factor of 10 different from 

the optimal value for CAL0 (table 9). The value for 

lakebed hydraulic conductivity was set to 107

 ft/d for 

these simulations. The number of reasonable data sets 

resulting from the conditioning was 490. No check was 

made for how the parameter sets affected model fit, and 

this would contribute to the uncertainty represented by 

the Monte Carlo analysis; thus, the actual uncertainty is 

overestimated.

Caution should be used in interpreting probabilistic contributing areas. A low probability at a given 

location does not mean that the location is likely to be 

outside the contributing area, rather it indicates that the 

model calibration data are not adequate to determine if 

the location is or is not outside the contributing area. 

The true contributing area is always unknown because 

one can never have a perfect description of the subsurface. Simulation is the best means of determining 

contributing areas (Franke and others, 1998) and model 

calibration data will always be limited; therefore, simulated contributing areas always will have some uncertainty associated with them. A location may be in the 

true contributing area but not in the simulated contributing area because a hydrologic feature, about which 

the calibration data do not contain much information, 

affects the simulation. The reverse is also true—a loca-
46 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

tion may be in the simulated contributing area and not 

in the true contributing area.

The uncertainty of the simulated contributing 

area could be reduced by collecting more data or by 

reformulating the model. If alternative models have 

been tested and ruled out, as they have in this study, 

reformulating the model without adding new data may 

not lead to great reductions in uncertainty. Parts of the 

probabilistic contributing areas at the STS are underlain by thinly laminated fine-grained deposits that have 

low probabilities. For reasons cited in this report, 

measuring representative hydraulic heads in thinly 

laminated fine-grained deposits may be problematic, 

and it may be difficult to improve the model calibration 

by collecting more data in these areas.

The number of times, expressed as a percentage 

of the total number of simulations (490), that a particle 

from a cell reached PW-3 ranged from 1 to 100 (fig. 

23). Particles from the coarse deposits along the 

western edge of the valley and from an area underlain 

by fine-grained deposits beneath and northwest of 

Stibbs Lake reached PW-3 more than 50 percent of the 

time. Particles from an area underlain by fine-grained 

deposits beneath and southwest of Stibbs Lake reached 

PW-3 less than 50 percent of the time. The probability 

is low here because of a combination of two factors that 

are particular to this model—this area is in the contributing area only if the vertical hydraulic conductivity of 

the fine-grained deposits is high, and most of the Monte 

Carlo parameter sets do not have high vertical 

hydraulic conductivity of the fine-grained deposits. If 

the calibration data better described the vertical 

hydraulic conductivity, there would be less variation of 

vertical hydraulic conductivity in the conditioned 

Monte Carlo data set, and the uncertainty that this area 

was either inside or outside the contributing area would 

be lower. 

The number of times, expressed as a percentage, 

that a particle from a cell reached PW-1 ranged from 1 

to 90 (fig. 24). Particles reached the pumped well more 

than 50 percent of the time from the area at the northernmost end of the modeled area. Other areas north of 

Stibbs Lake had lower probabilities (1 to 50 percent).

The combined probabilistic contributing areas 

depict the contributing area for both wells pumping 

simultaneously. Two lobe-shaped areas of high probability (in red on fig. 25) contribute water to wells PW1 and PW-3. Pathlines in the western lobe (pathline A 

on fig. 22) go directly from the contributing area to well 

PW-3. Pathlines in the eastern lobe area go directly 

from the contributing area to PW-1 (not shown on fig. 

22) and across the valley to PW-3 (pathlines B and C on 

fig. 22). Both lobes contain areas having short travel 

times and high probabilities. Short travel times are 

expected in these areas because of the high hydraulic 

conductivity of the coarse-grained deposits. High probabilities are expected in these areas because the calibration data were collected mainly in the coarse-grained 

deposits; therefore, the properties of the coarse-grained 

deposits are well known, and the probabilities are high. 

An area of long travel times and low probability 

separates the above mentioned lobes of the contributing area. Pathlines in this area go through finegrained deposits, which have lower hydraulic conductivity and less well known hydraulic properties. The 

deterministic contributing area (figs. 21) includes part 

of this area, but the probabilistic contributing area does 

not (fig. 25) because the probabilities are low and the 

Monte Carlo data set, being a random selection of 

likely parameter values, did not include the particular 

combination of parameter values used in the deterministic simulation. Any area in the deterministic contributing area should be considered also to be in the 

probabilistic contributing area. 

The Monte Carlo simulation shows that some 

areas outside the deterministic contributing area may 

actually be in the contributing area. The area having 

probabilities greater than 1 percent is larger than the 

deterministic contributing area. The deterministic 

contributing area for wells PW-1 and PW-3 (fig. 21) is 

roughly similar in size and shape to the area having 

probabilities greater than 25 percent (fig. 25). The two 

contributing areas (figs. 21 and 25) are significantly 

different in two places—southeast of Stibbs Lake and 

between well PW-3 and Transylvania Brook. These 

differences may have implications for management of 

the aquifer. For example, the area southeast of Stibbs 

Lake is underlain by coarse deposits and may appear to 

be appropriate for use as a septic tank/leach field site. 

In the deterministic simulation, such a decision might 

seem reasonable; however, the probabilistic simulation 

shows that water from this area may go to well PW-3. 

The area between well PW-3 and Transylvania Brook 

also is not in the deterministic area but is in the probabilistic area; therefore, water from Transylvania Brook 

may go to well PW-3. It is reasonable that these two 

areas are in the actual contributing area based on the 

geology of the site, but the calibration data were not 

adequate to define the contributing area with a higher 

probability. A greater probability might be achieved 

with more observation wells in these areas, a longer 

aquifer test at a higher pumping rate, and (or) an aquifer 

test within the area in question.
ANALYSIS OF UNCERTAINTY OF CONTRIBUTING AREAS 47

Figure 23. Probabilistic contributing area to well PW-3, Southbury Training School, Connecticut.
48 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Figure 24. Probabilistic contributing area to well PW-1, Southbury Training School, Connecticut.
ANALYSIS OF UNCERTAINTY OF CONTRIBUTING AREAS 49

Figure 25. Probabilistic contributing areas to wells PW-1 and PW-3, Southbury Training School, Connecticut.
50 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

SUMMARY AND CONCLUSIONS

The Southbury Training School (STS) lies 

mainly in the Transylvania Brook watershed in western 

Connecticut and withdraws ground water from stratified glacial deposits in the lower reaches of Transylvania Brook. The stratified deposits are found beneath 

the relatively (compared to the surrounding uplands) 

flat surfaces of the valley bottom. The deposits consist 

of coarse-grained deposits, which were deposited as 

ice-marginal deltas and fluviodeltaic deposits in close 

proximity to the retreating Pleistocene ice sheet, and 

fine-grained deposits, which were deposited in the 

quiet water conditions of a glacial lake bottom south of 

the retreating ice. The fine and coarse deposits interfinger and overlie one another in complex ways. The 

public-supply wells at STS withdraw most of their 

water from the coarse-grained deposits. 

Water flows into the stratified glacial deposits 

from upland runoff at the valley margins, through infiltration of stream water, and from direct precipitation on 

the valley bottom. Some water is captured by the watersupply wells at the STS, but most discharges to Transylvania Brook. Analysis of an aquifer test conducted 

at the primary public-supply well at the STS resulted in 

preliminary estimates of transmissivity for the stratified glacial deposits that ranged from about 2,900 to 

6,500 ft2

/d. The rate of ground-water recharge is estimated to be 24 in/yr over the valley bottom and 9 in/yr 

over the adjacent upland areas. 

Streamflow and ground-water-level data from 

nearby USGS network sites show the general hydrologic conditions during the aquifer test conducted from 

October 30 to November 4, 1996. During this period, 

the annual cycle of ground-water recharge was beginning, which led to generally rising ground-water levels. 

For the same period, stream levels were generally 

declining after a large amount of precipitation on 

October 20, 1996, 10 days before the aquifer test.

Ground-water flow was simulated for the entire 

watershed of lower Transylvania Brook to understand 

how topography, precipitation, and geology in the 

upland parts of the watershed interact to recharge the 

stratified glacial deposits in the valley. This simulation 

generated estimates of the distribution of recharge from 

10 small drainage basins in the Transylvania Brook 

watershed that were used in an aquifer-scale simulation. The estimates from each small drainage basin 

differed on the basis of the drainage characteristics of 

the basin. Small basins having well-defined stream 

channels contributed less recharge to the valley than 

basins having no defined channels because potential 

ground-water recharge was carried away in the stream 

channel.

Ground-water flow was simulated in the stratified glacial deposits to define the contributing areas to 

wells at STS. This simulation was done using a 

computer code that estimated the parameters of the 

ground-water-flow model and provided statistical 

measures of the goodness of fit of the model and the 

uncertainty associated with model predictions. The 

ground-water-flow simulation was done using 

MODFLOW-96, a general purpose, three-dimensional, 

finite-difference, ground-water-flow model. Contributing areas were computed using MODPATH, a 

particle-tracking code that works with MODFLOW-96 

models. The parameter estimation was done using 

UCODE, a universal nonlinear regression computer 

code. 

Parameters of the ground-water-flow model are 

defined in distinct layers and zones. Four zones are 

used in the simulation—zone 1 is the coarse-grained 

deposits, zone 2 is the fine-grained deposits, zone 3 is 

the bedrock, and zone 4 is the fault in the bedrock. Six 

hydraulic parameters are defined in the model—the 

horizontal and vertical hydraulic conductivities of 

zones 1 and 2, the specific yield of the uppermost layer 

and the specific storage of layers other than the uppermost. Parameters in zones 3 and 4 were not estimated 

in the simulation. Three additional parameters were 

used to estimate recharge and the hydraulic conductivities of the lake bed and streambed deposits. 

Regression modeling was done using aquifertest data collected at STS from October 30-November 

4, 1996. The data were first corrected for outside influences, which included the effects of long-term 

pumping at the STS, the effects of generally rising 

ground-water levels, and the effects of shutting off the 

pump at the STS prior to the test. Seven alternative 

models were posed, each representing the groundwater-flow system in slightly different but realistic 

ways. On the basis of the statistical measures of model 

fit and using the available data, one model was chosen 

as being the most representative of the ground-water 

flow system. The standard error of the regression for 

the chosen model is 2.27. The correlation coefficient 

between the weighted residuals from the regression and 

a normal distribution is 0.955, indicating that the residuals are nearly normal.
REFERENCES CITED 51

Parameter values estimated during the simulation are as follows: horizontal hydraulic conductivity 

of coarse-grained deposits, 154 feet per day; vertical 

hydraulic conductivity of coarse-grained deposits, 0.83 

feet per day; horizontal hydraulic conductivity of finegrained deposits, 29 feet per day; specific yield, 0.007; 

specific storage, 1.6E-05. Average annual recharge was 

estimated using the watershed-scale model with no 

parameter estimation and was determined to be 24 

inches per year in the valley areas and 9 inches per year 

in the upland areas.

The parameter estimates produced in the model 

are similar to expected values, with two exceptions. 

The estimated specific yield of the stratified glacial 

deposits (0.007) is lower than expected; this could be 

caused by the layered nature of the deposits. The 

recharge estimate produced by the model was also 

lower—about 32 percent of the annual average rate. 

This could be caused by the timing of the aquifer test 

with respect to the annual cycle of ground-water 

recharge and by some expected recharge going to parts 

of the flow system that were not simulated. The model 

reproduced the aquifer-test data accurately, as 

measured by the correlation coefficient (0.999) 

between the weighted simulated values and weighted 

observed values. The model also reproduced the 

general rise in ground-water levels caused by groundwater recharge and the cyclic fluctuations caused by 

pumping prior to the aquifer test. 

The simulated water table shows the effect of the 

fine-grained deposits on the shape of the contributing 

area to the well. The contributing area shows that most 

water enters the stratified glacial deposits through the 

coarse-grained deposits that ring the head of the lower 

Transylvania Brook valley. Some of these deposits are 

not saturated throughout the year and could not be 

simulated in the model; however, because the primary 

public-supply well receives most of its water from this 

area, the unsaturated deposits should be considered to 

be in the contributing area. Some travel times for 

ground water are less than 1 year; travel times for most 

of the contributing area are less than 2 years. 

A Monte Carlo simulation was done, using the 

variance/covariance matrix generated by the regression 

model, to estimate the probabilities associated with the 

contributing area to the wells. These probabilities arise 

from uncertainty in model, which in turn arises from 

the adequacy of the data available to comprehensively 

describe the ground-water-flow system. The Monte 

Carlo data sets were conditioned to remove unrealistic 

parameter sets. Probabilities in the contributing area 

range from 1 to 100 percent, and the highest probabilities (greater than 50 percent) are in the coarse-grained 

deposits that ring the head of the valley. The deterministic contributing area corresponds to the areas having 

probabilities of greater than 25 percent. 

REFERENCES CITED

Cooley, R.L., 1997, Confidence intervals for ground-water 

models using linearization, likelihood, and bootstrap 

methods: Ground Water, v. 35, no. 4, p. 869-879.

Cooley, R.L., and Naff, R.L., 1990, Regression modeling of 

ground-water flow: U.S. Geological Survey Techniques 

of Water-Resources Investigations, book 3, chap. B4, 

232 p.

Cooper, H.H., and Jacob, C.E., 1946, A generalized graphical method for evaluating formation constants and summarizing well field history: American Geophysical 

Union Transactions, v. 27, p. 526-534.

Fetter, C.W., 1994, Applied hydrogeology: Macmillan, 

691p.

Franke, O.L., Reilly, T.E., Pollock, D.W., and LaBaugh, 

J.W., 1998, Estimating areas contributing recharge to 

wells, lessons from previous studies: U.S. Geological 

Survey Circular 1174, 14 p.

Gates, R.M., 1954, The bedrock geology of the Woodbury 

quadrangle: Connecticut Geological and Natural History Survey Quadrangle Report 3, 23 p., scale, 

1:31,680.

———1959, Bedrock geology of the Roxbury quadrangle: 

U.S. Geological Survey Geologic Quadrangle Map 

GQ-121, scale 1:24,000.

Harbaugh, A.W., and McDonald, M.G., 1996, Users documentation for MODFLOW-96, an update to the U.S. 

Geological Survey modular finite-difference, groundwater-flow model: U.S. Geological Survey Open-File 

Report 96-485, 56 p.

Hill, M.C., 1992, MODFLOWP, A computer program estimating parameters of a transient, three-dimensional 

ground-water flow model using non-linear regression: 

U.S. Geological Survey Open-File Report 92-484, 

358 p.

———1994, Five computer programs for testing weighted 

residuals and calculating linear confidence intervals on 

results from the ground-water parameter-estimation 

computer program MODFLOWP: U.S. Geological Survey Open-File Report 93-481, 81 p.

———1998, Methods and guidelines for effective model 

calibration: U.S. Geological Survey Water-Resources 

Investigations Report 98-4005, 90 p.
52 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut

Hill, M.C., Cooley, R.L., and Pollock, D.W., 1998, A controlled experiment in ground water flow model calibration: Ground Water, v. 36, no. 3, p. 520-535.

Kontis, A.L., in press, Computer program for simulation of 

variable recharge with the U.S. Geological Survey 

modular finite difference ground-water flow model 

(MODFLOW): U.S. Geological Survey Open-File 

Report 00-173.

Mazzaferro, D.L., 1986, Ground-water availability and 

water quality in Southbury and Woodbury, Connecticut: 

U.S. Geological Survey Water- Resources Investigations Report 84-4221, 105 p. 

Mazzaferro, D.L., Handman, E.H., and Thomas, M.P., 1979, 

Water resources inventory of Connecticut, part 8, Quinnipiac River basin: Connecticut Water Resources Bulletin 27, 88 p. 

McDonald, M.G., and Harbaugh, A.W., 1988, A modular 

three-dimensional finite-difference ground-water flow 

model: U.S. Geological Survey Techniques of WaterResources Investigations, book 6, chap. A1, unpaginated.

Meinzer, O.E., and Stearns, N.D., 1929, A study of ground 

water in the Pomperaug basin, Connecticut, with special reference to intake and discharge: U.S. Geological 

Survey Water-Supply Paper 597-B, p. 73-146. 

Melvin, R.L., 1986, Connecticut observation wells—guidelines for network modification: U.S. Geological Survey 

Water-Resources Investigations Report 85-4079, 24 p.

Melvin, R.L., Stone, B.D., Stone, J.R., and Trask, N.J., 1992, 

Hydrogeology of thick till deposits in Connecticut: U.S. 

Geological Survey Open-File Report 92-43, 49 p.

Moench, A.F., 1994, Specific yield as determined by typecurve analysis of aquifer-test data: Ground Water, v. 32, 

no. 6, p. 949-957.

National Oceanic and Atmospheric Administration, 1991-

95, Climatological data, annual summaries, New 

England, 1991-95: Asheville, North Carolina, National 

Oceanic and Atmospheric Administration, v. 103-107, 

various pagination.

Nwankwor, G.I., Cherry, J.A., and Gillham, R.W., 1984, A 

comparative study of specific yield determinations for a 

shallow sand aquifer: Ground Water, v. 22, no. 6, 

p. 764-772.

Nwankwor, G.I., Gillham, R.W., van der Kamp, G., and 

Akindunni, F.F., 1992, Unsaturated and saturated flow 

in response to pumping of an unconfined aquifer—field 

evidence of delayed drainage: Ground Water, v. 30, 

no. 5, p. 690-700.

Poeter, E.P., and Hill, M.C., 1997, Inverse modeling, a necessary next step in groundwater modeling: Ground 

Water, v. 35, no. 2, p. 250-260.

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.

Pollack, 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.

Press, W.P., Flannery, B.P., Teukolosky, S.A., and Vetterling, 

W.T.,1986, Numerical recipes—the art of scientific 

computing: Cambridge University Press, 818 p.

Prudic, D.E., 1989, Documentation of a computer program 

to simulate stream-aquifer relations using a modular, 

finite-difference, ground-water flow model: U.S. Geological Survey Open-File Report 88-729, 113 p.

Rodgers, John (compiler), 1985, Bedrock geological map of 

Connecticut: Connecticut Geological and Natural History Survey, Natural Resources Atlas Series Map, 

2 sheets, scale 1:125,000. 

Rutledge, A.T., 1997, Model-estimated ground-water 

recharge and hydrograph of ground-water discharge to 

a stream: U.S. Geological Survey Water-Resources 

Investigations report 97-4253, 29 p.

———1998, Computer programs for describing the recession of ground-water discharge and for estimating mean 

ground-water recharge and discharge from streamflow 

records—update: U.S. Geological Survey WaterResources Investigations Report 98-4148, 43 p.

Scott, R.B., 1974, The bedrock geology of the Southbury 

quadrangle: Connecticut Geological and Natural History Survey Quadrangle Report 30, 63 p., scale 

1:24,000.

Stanley, R.S., and Caldwell, K.G., 1976, The bedrock geology of the Newtown quadrangle: Connecticut Geological and Natural History Survey, Quadrangle Report 33, 

43 p., scale 1:24,000.

Stone, J.R., Schafer, J.P., London, E.B., and Thompson, 

W.B., 1992, Surficial materials map of Connecticut: 

U.S. Geological Survey Special Map, 2 sheets, scale 

1:125,000.

Walton, W.C., 1984, Practical aspects of groundwater modeling—Flow, mass and heat transport, and subsidence; 

analytical and computer models: Dublin, Ohio, 

National Water Well Association, 566 p.

Wilson, W.E., Burke, E.L., and Thomas, C.E., Jr., 1974, 

Water resources inventory of Connecticut, part 5, lower 

Housatonic River basin: Connecticut Water Resources 

Bulletin 19, 79 p.
APPENDIX 1 53

Appendix 1. Geohydrologic data for selected wells, Transylvania Brook watershed, Connecticut 

[-, not available; 112DFSF, stratified glacial deposits; 112TILL, glacial till; 231 SDMN, sedimentary rock; NCBC, noncarbonate crystalline rock; 

SB, Southbury]

USGS identifier USGS aquifer 

code

Altitude of top of rock 

(feet above sea level)

Altitude of water level 

(feet above sea level)

Specific capacity 

(gallons per minute 

per foot)

SB 46 231SDMN 185.00 185.00 0.01

SB 47 300NCBC 45.00 135.00 0.05

SB 48 231SDMN 55.00 153.00 0.19

SB 50 300NCBC 295.00 285.00 0.01

SB 51 300NCBC 285.00 285.00 0.01

SB 52 300NCBC 125.00 149.00 0.01

SB 56 231SDMN 42.00 135.00 0.01

SB 59 300NCBC 35.00 126.00 0.23

SB 60 300NCBC 447.00 350.00 0.01

SB 61 300NCBC 400.00 430.00 0.00

SB 63 231SDMN 65.00 165.00 0.83

SB 69 300NCBC 315.00 275.00 0.01

SB 70 231SDMN 25.00 155.00 0.45

SB 74 231SDMN 139.00 137.00 0.28

SB 75 231SDMN 95.00 130.00 20.00

SB 76 300NCBC 37.00 181.00 0.11

SB 81 300NCBC 39.00 153.00 0.00

SB 84 300NCBC 15.00 153.00 0.20

SB 86 300NCBC 5.00 183.00 0.11

SB 88 300NCBC 289.00 279.00 0.01

SB 89 300NCBC 343.00 338.00 0.01

SB 90 231SDMN 335.00 379.00 0.20

SB 91 300NCBC 265.00 267.00 0.02

SB 92 231SDMN 358.00 345.00 0.04

SB 94 231SDMN 357.00 475.00 0.09

SB 95 231SDMN 280.00 437.00 0.01

SB 96 231SDMN 375.00 435.00 0.05

SB 97 300NCBC 250.00 420.00 0.01

SB 98 300NCBC 428.00 423.00 0.04
54 Delineation and Analysis of Uncertainty of Contributing Areas to Wells at the Southbury Training School, 

Southbury, Connecticut