Hydrology of Mountain Blocks in Arizona and New Mexico as Revealed by Isotopes in Groundwater and Precipitation

Abstract

Mountain-block groundwater in the Southern Basin-and-Range Province shows a variety of patterns of δ¹⁸O and δ²H that indicate multiple recharge mechanisms. At 2420 m above sea level (masl) in Tucson Basin, seasonal amount-weighted means of δ¹⁸O and δ²H for summer are −8.3, −53‰, and for winter, −10.8 and −70‰, respectively. Elevation-effect coefficients for δ¹⁸O and δ²H are as follows: summer, −1.6 and −7.7 ‰ per km and winter, −1.1 and −8.9 ‰ per km. Little altitude effect exists in 25% of seasons studied. At 2420 masl, amount-weighted monthly averages of δ¹⁸O and δ²H decrease in summer but increase in winter as precipitation intensity increases. In snow-banks, δ¹⁸O and δ²H commonly plots close to the winter local meteoric water line (LMWL). Four principal patterns of (δ¹⁸O, δ²H) data have been identified: (1) data plotting along LMWLs for all precipitation at >1800 masl; (2) data plotting along modified LMWLs for the wettest 30% of months at <1700 masl; (3) evaporation trends at all elevations; (4) other patterns, including those affected by ancient groundwater. Young, tritiated groundwater predominates in studied mountain blocks. Ancient groundwater forms separate systems and mixes with young groundwater. Recharge mechanisms reflect a complex interplay of precipitation season, altitude, precipitation intensity, groundwater age and geology. Tucson Basin alluvium receives mountain-front recharge containing 50%–90% winter precipitation.

1. Introduction

In the Basin and Range Province [1] of Southwestern North America, mountain blocks play an important role in regional hydrology. Water from the highest mountain ranges is conveyed to adjacent basin-fill aquifers either as surface flow or as sub-surface mountain-block recharge [2]. Since the 1880s, the alluvial aquifers have been a crucial water source for agriculture, mining and urban centers in the region. Groundwater hydrology studies in mountain ranges are potentially difficult because of the challenges of data collection in remote and inaccessible terrain with few sampling points [2].

Studies of mountain-system hydrology have been conducted at a variety of spatial scales. Smaller-scale studies have focused on water budgets of individual hillslopes or small watersheds, using parameters such as water chemistry, temperature and isotopes to examine the distribution of available water between runoff, soil water and bedrock groundwater [3,4,5,6,7,8]. At a larger scale, isotope and noble gas ratio data have identified groundwater moving from a mountain block into basin alluvium in Salt Lake Valley, Utah, USA [9,10]. Stable isotope studies have addressed the residence time of mountain-block groundwater of Table Mountain, South Africa, and the Cascade Range, Oregon, USA [11,12].

A potential approach to mountain-block hydrology consists in comparing isotope data of groundwater from springs and wells with isotopes in mountain precipitation, taking into account the altitude and seasonal effects, with a view to constraining recharge mechanisms and flow paths. Such studies have been undertaken in the Basin-and-Range Province [13,14,15,16,17]. Figure 1 is a conceptual sketch of a mountain block, showing possibilities for water movement that might be evaluated using stable isotope data.

In this article, stable oxygen and hydrogen isotope data from environmental waters are reviewed in an assessment of groundwater characteristics in six mountain blocks (Santa Catalina, Rincon, Tucson, Santa Rita, Galiuro, and Sacramento) of the Basin and Range Province in Southern Arizona and New Mexico, and in the Mogollon Highlands at Payson, Arizona (Figure 2). Data for the Huachuca and Chiricahua Mountains (Mts.) (Figures S1 and S2), which illustrate less clearly the phenomena to be described in the other locations, are reviewed in Supplementary Information. The aims of the study are: (1) detailed evaluation of isotope character of high-altitude precipitation in the Santa Catalina (Mts); (2) comparison of isotope data for groundwater in multiple mountain blocks; (3) interpretation of stable O and H isotope data in the context of isotope altitude effects, local geology and groundwater ages; (4) interpretation of isotope data patterns in terms of recharge mechanisms; and (5) re-evaluation of seasonality of mountain-front recharge to the alluvial aquifer in Tucson Basin.

2. Materials and Methods

2.1. Study Area

The mountain blocks considered here are hard-rock ranges of the Basin-and-Range province, plus part of the southern rim of the Colorado Plateau. The Basin-and-Range Province formed by Neogene tectonic extension, part of an extended episode of taphrogenesis reflecting the interaction of mid-ocean ridges and triple junctions with trenches along the western edge of the continent [18]. Extension split a hard-rock plateau into a series of mountain ranges separated by basins filled with alluvium, a process that began 15 Ma ago in Southeastern Arizona [18]. The ranges are fault-bounded horsts or tilted blocks of consolidated Early Proterozoic to Neogene rocks [19]. The complex interior structure of the ranges reflects evolution from a Late Cretaceous compressive stress regime to a Neogene extensional regime [19].

The region has an arid to semiarid climate, except at the crests of the higher mountain ranges. Annual average precipitation ranges from 75 to over 1000 mm annually, with the higher totals occurring at the range crests [20]. Potential evaporation greatly exceeds rainfall [20,21]. Rainfall typically occurs in two seasons or modes: in summer, localized, potentially intense convective precipitation results from the North American monsoon, and in winter, widespread orographic precipitation originates from Pacific frontal systems. These precipitation modes are consistently separated in the spring by a hot period with little or no rain during May and June. Dry, lower mountain slopes typically support grasses, forbs, succulents and shrubs, succeeded upwards by stunted evergreen oaks and conifers. Well-watered range-crests support forests of conifers and deciduous trees. Most woodland parts of the study areas have little human settlement and are protected by National Forest or National Park designation; in such areas, fire has been suppressed for many years. Since the onset of long-term drought in the late 1990s, the vegetation and surface hydrology of large areas have been greatly affected by forest fires. Higher-density human settlement in these mountain ranges is present locally in the Southern Sacramento, Santa Catalina and Tucson Mts. and the Mogollon Highlands.

When sampling for this study began in the 1990s, springs of modest discharge were widely distributed in most of the ranges considered here, except for the Tucson Mts. Attempts at re-sampling since 2012 have shown that many springs no longer discharge in the Santa Catalina and Galiuro Mts. The change in hydrology reflects drought conditions that have prevailed since about 1998 [22].

The hydrogeology of each range is strongly influenced by its lithologies and structures. In the following summary, representative references are given for the geology. Detailed geological maps and cross-sections can be found in those references.

2.2. Santa Catalina and Rincon Mts

These ranges form the northern and northeastern boundary of Tucson Basin, and constitute a single structural block of the Basin-and-Range Province. In the Santa Catalina Mts. (SCM), the highest peak is Mt. Lemmon, 2791 masl. Steep slopes and V-shaped drainages are characteristic. In the Rincon Mts. (RM), the highest peak is Mica Mountain, 2641 masl. The central massifs of both ranges consist mainly of granitoid intrusions, including thick sheets of Paleogene equigranular peraluminous leucogranite and associated pegmatite intruding tilted Proterozoic and Paleozoic strata to the northeast, and intensely deformed Proterozoic granite to the southwest [23,24]. The SCM-RM block formed the lower plate of a Neogene low-angle detachment fault marked by a band of impermeable mylonite [25]. A pronounced set of sub-horizontal leucogranite sheets on Mica Mountain (RM) creates a stepped topography that is absent in the SCM.

Previous stable isotope studies in the SCM addressed the origin of spring water [26], the origin of groundwater associated with the detachment fault [27], the nature of precipitation at a station at 2420 masl [28] and relationships among precipitation, critical zone groundwater, soil water and stream flow [6,8,29,30].

2.3. Tucson Mts

The Tucson Mountains form a terrain of steep slopes and sharp peaks interspersed with broad valleys and pediments, and bounding Tucson Basin to the west. Wassen Peak, 1428 masl, is the highest point. Late Cretaceous–Paleogene granitoids intrude coeval felsic volcanic rocks and Permian to Jurassic sedimentary rocks [24]. Fracture-hosted groundwater is pumped from as much as 170 m below valley-floors and pediments; water levels are declining rapidly. Detailed isotope (stable O, H in water, stable S, O in sulfate, ³H, ¹⁴C, 87Sr/86Sr) and geochemical data were presented for fracture-hosted groundwater from a 9 km² area near the Yuma Mine, 32.19° N, 111.73° W [31].

2.4. Santa Rita Mts

In the steep-sided massif of Mt. Wrightson (2882 masl), Triassic and Jurassic granitoids intrude Triassic felsic volcanic rocks with interbedded eolian sandstone. Paleozoic limestone and Mesozoic clastic sedimentary rocks underlie the eastern slope of the massif [32]. Published groundwater isotope data (stable O, H and tritium) are available for the Rosemont area at 1500–1700 masl on the northern ridge of the range [33], where complex faulting juxtaposes Proterozoic granitoid, Paleozoic marine shelf strata and Mesozoic clastic sedimentary rocks [32,34,35].

2.5. Galiuro Mts

Steep slopes and V-shaped topography are characteristic. The highest point is Bassett Peak, 2332 masl. Late Paleogene felsic and intermediate volcanic rocks unconformably overlie discontinuous outcrops of Late Cretaceous–Paleogene andesite and granitoids [36]. At the south end of the range, a group of hot springs (Hooker and nearby hot springs) discharges at 53 °C. Stable O, H isotope data were measured for a few springs at the north end of the range [37].

2.6. Mogollon Highlands

Paleozoic strata, predominantly limestone, sandstone and siltstone, unconformably overlie Proterozoic granitoids [38]. Highest elevations, 2300–2400 masl, occur at the Mogollon Rim (the southern edge of the Colorado Plateau) where the strata are flat-lying. South of the Mogollon Rim, at elevations of 1500–1700 masl, the Paleozoic strata and Proterozoic granitoids are juxtaposed in complex fashion by normal faults. Available reports address the hydrogeology of the Mogollon Highlands [39] and isotopes in groundwater [40].

2.7. Sacramento Mts

The highest point in the range is Cathey Peak, 2923 masl. Landforms at high elevation include V-shaped canyons, broad valleys with thick soil, and local karst features. Paleozoic strata, predominantly limestone, are bounded by steep normal faults to the west, and dip eastward. Small springs discharge where shale beds form aquitards, but regional hydrology is dominated by long-distance karst permeability [16]. Groundwater isotope (stable O, H in water, ³H, ¹⁴C) and geochemical data amassed during wet years, 2006–2009, [16] were compared with isotope data collected in 2003, following a dry period [17].

2.8. Recharge Seasonality

The difference in stable isotope composition between summer and winter rain in Tucson Basin has been used with isotope mass-balance to determine that winter recharge is dominant in the basin aquifer [41,42,43,44]. Both O and H isotopes should be used in the mass-balance approach [45] because of an isotope effect related to precipitation intensity like that in tropical regions [46]. Two principal recharge mechanisms involving local precipitation (not river water of remote derivation) operate in alluvial basins of Arizona [45]. Basins in Northwestern Arizona receive recharge from winter precipitation, with strong evaporation effects. Basins in Southeastern and Central Arizona receive both summer and winter recharge, but only from precipitation falling in the wettest months. The seasonality of recharge in mountain blocks may not correspond to that in adjacent basins, and may vary with altitude [47].

2.9. Altitude Effects

The altitude effect on stable isotopes in precipitation [48] consists in decreases in δ¹⁸O and δ²H with increasing sampling elevation. The altitude-δ¹⁸O relationship has been used to determine altitudes of recharge [49,50] and paleoelevations in ancient mountain ranges [51,52]. For worldwide stations at latitudes <70° and altitudes <5000 m, the relationship between δ¹⁸O and altitude is linear, and the range of isotopic lapse rate for δ¹⁸O is −1.0 to −5.1‰/km with an average of −2.8 ± 1.3 (1σ) ‰/km [53]. Snow sublimation may be superimposed on isotope altitude effects, decreasing lapse rates of δ¹⁸O and leading to underestimates of paleo-elevation [54].

Isotope altitude effects are commonly absent on mountain ranges in continental interiors and in snow [55]. Altitude effects in interior ranges of Arizona and New Mexico appear to be complex [45]. Altitude effects are present in summer but not in winter on a north-facing transect near Prescott, Arizona, outside the present study area; on a neighboring south-facing transect, altitude effects are present in both seasons, but are imperfectly correlated with altitude [56]. On the east-facing slope of the Sacramento Mts., New Mexico, normal altitude effects were observed in summer, but reverse effects (δ¹⁸O increasing with altitude) occurred in two winters of three studied [16].

In Tucson Basin, data collected at 2420 masl between 1995 and 2004 at Palisades Ranger Station (PRS) in the Santa Catalina Mts. [28] (Supplementary Table S1) have been compared with a 32-year dataset collected at 747 masl at the University of Arizona (UA) station in Tucson [22,28]. PRS has an average precipitation of 750 mm/year; precipitation at PRS is more frequent than at UA, where average precipitation is 360 mm/yr. The long-term estimates of altitude coefficients (Figure 3a) are based on all precipitation events for each station.

2.10. Sampling and Analyses

Spring water was sampled as close as possible to the discharge point. Water supply wells in continual use were sampled without further purging. Samples were stored in robust, tightly sealed plastic bottles. Measurements of δ¹⁸O and δ²H were made on an automated gas-source isotope ratio mass spectrometer (Finnigan Delta S) at the Environmental Isotope Laboratory of the University of Arizona, Tucson, Arizona. For H isotopes, samples were reacted at 750 °C with Cr metal using a Finnigan H/Device coupled to the mass spectrometer. For O isotopes, samples were equilibrated with CO₂ at approximately 15 °C in an automated equilibration device coupled to the mass spectrometer. International reference materials VSMOW and SLAP were used for standardization [57]. On the basis of repeated internal standards, analytical precision (1σ) was 0.9‰ or better for δ²H and 0.08 ‰ or better for δ¹⁸O. Prior to 1992, water samples were prepared manually by reduction to H₂ gas using Zn metal for H isotopes, and CO₂ equilibration for O isotopes. Analytical precision (1σ) was poorer, 1.5‰ to 2‰ for δ²H, and 0.15‰ for δ¹⁸O. Tritium was measured on 0.18 L samples following electrolytic enrichment, using Quantulus 1220 spectrophotometers with a detection limit of 0.6 tritium units (TU). Calibration was based on National Institute of Standards and Technology standard reference materials 4361 B and C. Carbon-14 measurements were made on 50 L water samples at the Radiocarbon Laboratory, University of Arizona. Dissolved inorganic carbon was precipitated as BaCO₃, processed to benzene, and measured using Quantulus 1220 liquid-scintillation spectrophotometers. The detection limit for full-size samples was 0.8 percent modern carbon (pMC) or lower. Calibration was based on international standard Oxalic Acid I.

In this article, the presence of finite tritium in a groundwater sample is used as evidence that the sample contains at least some recharge from precipitation that fell since about 1953. Finite tritium does not preclude the presence of older recharge mixed with the younger water.

3. Results

Previously unpublished isotope data are listed with sample location data in Supplementary Tables S1–S3. Sources of published groundwater data are cited in each section, below.

3.1. Precipitation, Santa Catalina Mts

Figure 3a shows all available data for PRS [28]. Gaps exist, and many data represent composites of rain events over multiple weeks (Table S1). Data for summer (defined here as June to October) and winter (November to May) span broad, largely identical ranges; most data plot on a linear trend above the global meteoric water line (GMWL) [58] between δ¹⁸O values of −19‰ and −3‰. Summer data predominate at δ¹⁸O > −5‰. Seasonal local meteoric water lines (LMWL), weighted for precipitation amount [59] differ in slope and intercept (Figure 3a). Long-term amount-weighted averages (Figure 3b) are (δ¹⁸O, δ²H) = (−10.8‰, −70‰) for winter and (−8.3‰, −53‰) for summer.

3.2. Altitude Effects at Annual Time Scale

Differences in seasonal amount-weighted mean δ¹⁸O between PRS and UA are distributed evenly between 0, signifying no altitude effect, and 4.9, corresponding to −2.9 ‰/km (Figure 3c). In Southern Arizona, estimation of average isotope altitude effects clearly requires data for multiple years. Available data indicate little to no altitude effect (< 1‰ difference in Figure 3c) in four of sixteen individual seasons.

3.3. Isotopes in Snow

Winter precipitation above 2000 masl in the Santa Catalina Mts. is commonly snow, locally forming snow-banks that may persist until the spring equinox. Vertical sections through a snow-bank on Mt. Bigelow Road (32.4158° N, 110.7302° W) were sampled between 2007 and 2012 to investigate the occurrence of isotope shifts due to sublimation or other snow-metamorphic processes that are known elsewhere in the Basin and Range Province [54,60]. Except in 2007, the data form a scattered linear trend (Figure 3d), essentially conforming to the winter LMWL. In 2007, much of the snow plotted far above the GMWL, reflecting snow-bank isotope fractionation opposite in sense to that observed in Southeast California and Nevada [54]. Three meltwaters plot with snow-bank data, and a fourth plots to the right of the GMWL, as predicted for snow that has partly sublimed [14]. Large isotope fractionation within snow-banks does not occur in most years in the Santa Catalina Mts. When present, it may result from sublimation and recondensation of vapor along a thermal gradient [61]. Snow metamorphism studies in Alaska [61] and the Southwest USA [54,60] do not show the large isotope shifts to the left of the GMWL that were identified in this study (Figure 3d).

3.4. Precipitation Intensity Effects

Recharge in subtropical and tropical regions has the isotope composition of precipitation falling in the wettest months [46]. This effect is also present in Southern Arizona, where an isotope intensity effect exists in the UA precipitation isotope dataset, and recharge appears to be dominated by precipitation falling during approximately the wettest 30% of months (Figure 4a) [45]. The PRS precipitation dataset is not ideal for such analysis, because of composite samples that are not easily arranged into calendar-month brackets, and because of the relatively small number of monthly data (30 for winter, 37 for summer, compared with 159 and 134, respectively, at UA). Where necessary, data were combined into approximate month-long non-calendrical brackets. Limitations notwithstanding, the wettest summer months (plotted at x ≥ 60 in Figure 4b) at PRS have lower δ¹⁸O and δ²H than the long-term amount-weighted averages (plotted at x = 0). In winter, the opposite occurs: the wettest months have higher values of δ¹⁸O and δ²H than the long-term averages. Means representing the wettest 30% of months at PRS are nearly collinear with the LMWL as defined by long-term means of all data (Figure 5a).

Therefore, isotope data cannot reveal a relationship between precipitation intensity and recharge at the elevation of PRS; such relationships rely on a distinction between the LMWL (long-term means) and the line including wettest-month means [45,46]. Recharge may nonetheless occur mainly during the wettest summer months at high altitudes. In winter, this is less likely because snow-banks and their meltwater can represent all winter precipitation events. The isotope composition of high-altitude winter recharge is, therefore, better estimated by the weighted mean for all winter events.

3.5. Groundwater, Santa Catalina Mts

Figure 5a shows data for the main range of the Santa Catalina Mts., sampled in 1994–1996 [26], omitting data for samples that had apparently evaporated between the discharge and sampling points. Most samples plot close to the GMWL. Springs discharging at elevations above 2300 masl have a restricted range of isotope composition, plotting near or just below the amount-weighted mean for winter precipitation. Springs discharging below 2300 masl have a much broader range that extends from mean summer precipitation to values far below mean winter precipitation. For springs discharging below 1400 masl, the range of isotope composition is more restricted and possibly related (with one exception) to mean precipitation for the wettest months near 1400 masl.

Plotting δ¹⁸O as a function of discharge altitude (Figure 5b), it is possible to constrain recharge seasonality. Reasoning that recharge must occur at elevations above the discharge points, winter recharge is predominant for samples plotting within the green rectangle, regardless of recharge altitude. Similarly, samples plotting in the blue rectangle have approximately equal summer and winter recharge, or predominant summer recharge. Points to the left of the green rectangle can only be explained by recharge of water from unusually low-δ¹⁸O events, recorded at PRS in both summer and winter (Supplementary Table S1).

All spring discharge in the main massif of the range contained tritium, and, therefore, contains some post-1953 recharge [26]. Groundwater flow appears to occur on a time scale of decades, most likely through shallow fracture systems. The drying of many springs since the beginning of drought conditions in about 1998 is additional evidence of the short time scale of flow. Winter recharge at the range crest could supply discharge to springs with data plotting in the green rectangle (Figure 5b). All springs considered could potentially have catchments extending to the range crest. However, more than half of the springs below 2300 masl were supplied with water of other isotope composition. Recharge was, therefore, occurring at lower elevations, with a wide range of combinations of winter and summer precipitation. No explanation is possible at present for apparent preferential recharge after unusually low-δ¹⁸O events; in winter, such events are not the largest precipitation events (Figure 4b; Supplementary Table S1).

3.6. Springs, Rincon Mts

Two sets of samples are considered (Figure 6). Spring water was collected in 2010 on Mica Mountain, the main massif of the range, by National Parks Service staff [62], mostly from small pools into the bottom of which spring water was seeping; some of the samples showed strong evidence of evaporation and others may be slightly evaporated. A second set of spring and well water samples was collected close to outcrop of the Catalina detachment fault by one of us (CJE), Joy Gillick, and Alex Leonard between 1996 and 2014.

Data for springs discharging above 2180 masl form an array between mean winter and summer precipitation at PRS, and show no evidence of evaporation, with one exception. Springs discharging below 2180 masl on Mica Mountain form an evaporation trend of slope near 3.5, the trend apparently also including two of the springs discharging above 2180 masl, and the mean isotope composition of the wettest 30% of summer months. Two intensely evaporated samples with (δ¹⁸O, δ²H) = (+7.3, −2‰) and (+15.0, +23‰) are omitted from the plot, but lie on the evaporation trend and are included in the regression calculation. Most of the spring waters contained finite tritium when sampled in 2010.

Groundwater drawn from wells in hard rock and from Agua Caliente spring near the Catalina detachment fault at the base of the Rincon Mts. plots near (δ¹⁸O, δ²H) = (−10, −76‰). A similar set of data, predating automated analytical techniques, was collected in the early 1980s [27]. The groundwater is of sodium-sulfate type; sulfate-oxygen isotopes indicate that the sulfate source is oxidized sulfide [65]. Tritium is below detection, and ¹⁴C content ranges from 16 to 47 percent modern carbon (pMC). Corrected dates are not offered, because sulfide oxidation has probably enhanced the dissolution of accessory carbonate of unknown δ¹³C value from the aquifer.

Three recharge mechanisms are present in the Rincon Mts. (1) At high-altitudes, various combinations of winter and summer recharge plot on the high-altitude LMWL. Most groundwater of this kind bears detectable tritium (Supplementary Table S2) and appears to flow with short residence time through shallow fractures. (2) High-altitude recharge that is either equal parts high-altitude winter and summer precipitation (considering seasonal amount-weighted means for all precipitation data) or entirely precipitation from the wettest 30% of summer months. This kind of water appears to recharge and discharge cyclically in stream beds, undergoing evaporation in the process. (3) Recharge of slightly evaporated water with (δ¹⁸O, δ²H) lower than the values for present-day weighted-mean winter precipitation at PRS. The groundwater is of long residence time, probably thousands of years, apparently flowing through a deep fracture system to an aquifer confined beneath mylonite of the Catalina detachment fault. Such groundwater may be a mixture between ancient precipitation (represented in Figure 6 by eleven data points for Southeastern Arizona groundwater containing <8 pMC [44,63,64] with more recent high-altitude recharge. Note that the data from [64] include four samples without ¹⁴C data, but from the same deep, saline aquifer as dated samples.

3.7. Groundwater, Tucson Mts

The Tucson Mts. are of lower elevation and have lower precipitation than the other ranges surrounding Tucson Basin. Samples were taken from domestic wells and from one spring, all in hard rock. Values of (δ¹⁸O, δ²H) (Supplementary Table S2) form a linear array approximately parallel to and to the right of the GMWL (Figure 7). Samples that contain finite tritium plot at the upper end of the trend, while the oldest sample, containing 7.8 pMC, plots at the lower end. This suggests a relationship between stable O and H isotopes and groundwater residence time, although the age distribution (indicated approximately by uncorrected pMC) is complex in the middle of the array. Values of (δ¹⁸O, δ²H) at the upper end of the trend fit a modified LMWL for the wettest 30% of months at 1000 masl, rather than the LMWL for all data (both lines interpolated between UA and PRS altitudes). A precipitation intensity effect operates on recent recharge, which arises equally from winter and summer precipitation, or predominantly from summer precipitation. The lower end of the array indicates ancient recharge at a time of lower (δ¹⁸O, δ²H) in average precipitation.

3.8. Groundwater, Santa Rita Mts

Data for springs higher than 1800 masl on the Mt. Wrightson massif (Figure 8) plot on the LMWL for precipitation at PRS, between summer and winter means, with no indication of evaporation. The range of δ¹⁸O is similar over about 800 m of altitude. All samples in which tritium was measured contain finite tritium (Supplementary Table S2). In this case, all recharge may occur at the range crest, with approximately equal winter and summer contributions.

Data from Rosemont, for wells and springs below 1600 masl, plot either along a line passing through weighted-mean averages of (δ¹⁸O, δ²H) for the wettest 30% of months at 1600 masl (which is close to collinear with the GMWL), or to the right of that line. Most data form a cluster with a corresponding evaporation trend (line A in Figure 8, slope 4), and are consistent with about-equal contributions of summer and winter recharge. A few data points suggest predominant winter recharge but may also reflect mixing with ancient recharge. The field of data is delimited downward by an evaporation trend (line B), which appears to be related to the data for ancient water in Figure 6. No ¹⁴C measurements are available for this data set.

3.9. Groundwater, Galiuro Mts

Most samples taken in 1994 plot along an evaporation trend of slope near 4 (Figure 9). Care was taken to avoid water that might have evaporated since discharge. A few samples, including Hooker Hot Springs, plot along the GMWL. Samples taken in 1982 [37] plot close to the evaporation trend. Samples taken in 2012 at Upper and Lower Ash Springs plot closer to the GMWL. Other springs that might have been re-sampled proved to be dry in 2012. Most samples taken in 1994 contained finite tritium, between 1 and 6 TU, at the time of sampling. A sample from Hooker Hot Spring contained tritium below detection, and 17 pMC.

The data are consistent with two main recharge mechanisms in the range. Cold groundwater of decadal residence time, evaporated prior to infiltration, originated as recharge plotting below present mean winter precipitation at 2000 masl (based on Tucson Basin data) but resembling mean winter precipitation that has no altitude effect in Northern Arizona [45,56]. By 2012, this recharge regime may no longer have operated in the catchments of Upper and Lower Ash Springs. At Hooker Hot Springs at the south end of the range, hot groundwater of long residence time, probably thousands of years, was not evaporated prior to recharge and flows through a deeper fracture system.

3.10. Groundwater Mogollon Highlands

Groundwater from an area between the lower flanks of the Mogollon Rim and the Town of Payson yielded (δ¹⁸O, δ²H) data forming a linear evaporation trend of slope 5 (Figure 10), and show no relationship to altitude [40]. Samples with δ¹⁸O > −10‰ are affected by recharge of treated water from the Payson Town Lake. The data trend is closely similar to trends for alluvial basins in Northwestern Arizona and for major rivers that rise along the Mogollon Rim. Such trends originate in regional mean winter precipitation that has no altitude effect in Northern Arizona [45,56].

3.11. Groundwater, Sacramento Mts

The data (Figure 11) represent springs and precipitation sampled in 2006–2009, which were wetter than average years [16], and springs and wells sampled in 2003, following several dry years [17]. All samples are from carbonate-dominated strata.

Data collected in 2003 from springs and wells at 1750 to 2550 masl form a linear trend (SMT) of slope 5.6, consistent with evaporated winter precipitation (measured later, in 2007 and 2008) as the principal source of recharge. The slope is high for an evaporation trend and may be influenced by summer recharge. Summer recharge is clearly present during 2006–2009, when the data for springs discharging between 1750 and 2400 masl form an array stretching between the SMT and summer precipitation. Discharge from springs between 2400 and 2800 masl form a second array stretching between the SMT and a mixture of summer and winter precipitation. The groundwater system is dynamic, susceptible to change in isotope content at annual time scales. Winter recharge with evaporation dominates at all altitudes above 1750 masl in dry years, but in years when above-average summer rainfall occurs, summer recharge appears in groundwater.

3.12. Seasonality of Mountain-front Recharge, Tucson Basin

Water from the Santa Catalina and Rincon Mts. is supplied to Tucson basin by three major catchments: Rillito Creek, Cañada del Oro and Rincon Creek (Figure 12a). Mountain-front recharge occurs at distances up to 10 km from range fronts [44]. Recharge occurs mainly during the wettest 30% of months [45].

Assumptions for the re-evaluation are: (1) that precipitation amount and isotope lapse rates are linear functions of altitude, established by data for the UA and PRS sites, and (2) that runoff from any altitude has the same likelihood of reaching the basin floor. Therefore, the bulk isotope content of runoff from the mountain block is the isotope content of precipitation at some altitude A, above and below which the volume of precipitation is the same. Altitude A is determined by considering each mountain range as a cone divided into n equally spaced layers (Figure 12b). The results for a pyramid, or for intersecting cones and pyramids are identical to those for the cone, and do not depend on apical angle. Each layer has an interpolated annual precipitation amount and an annular footprint area; the volume of precipitation for the layer is the product annular footprint x precipitation (mm). Volumes are summed for all layers above and below a chosen layer. For the Santa Catalina Mts., using 100 layers from 1000 to 2800 masl, altitude A is about 1660 masl (Supplementary Table S4). Interpolated, amount-weighted, seasonal values of (δ¹⁸O, δ²H) in precipitation at 1660 masl are compared with groundwater isotope data for the parts of Tucson Basin recharged from mountain runoff (Figure 13).

The groundwater data do not conform to the LMWL defined by the seasonal means of all precipitation data at 1660 masl (Figure 5b), but instead fall on a modified LMWL corresponding to the wettest 30% of months at 1660 masl. The modified LMWL was calculated by applying the long-term isotope lapse rates for Tucson Basin (Figure 3a) to the amount-weighted mean isotope compositions for the wettest 30% of months at the UA site (Figure 4). In order to avoid scatter of isotope data in younger, less integrated groundwater, only samples containing less than 1.5 TU were used [44]. In areas recharged from Rincon Creek and Cañada del Oro, winter recharge makes up 75% to 90% of the total. In areas recharged from Rillito Creek, winter recharge constitutes 50% to 75% of the total. Water with about 50% winter recharge (samples in Figure 13 with δ¹⁸O ≤ −9.3‰) occurs in discrete areas at some distance from recharge sources (Figure 12a).

4. Discussion

4.1. Residence Times

Where tritium data are available, they indicate the prominence of young groundwater in each of the mountain blocks. Residence times of multiple decades are likely. Ancient groundwater is recognized in certain blocks, either from ¹⁴C data or from distinctive stable O and H isotope data. Available ¹⁴C data without correction suggest residence times of thousands of years. The shortest residence times are found in the carbonate-dominated aquifers of the Sacramento Mountains, where groundwater isotopes change at a time-scale of 5 years or less.

Contrasting systems of young and ancient groundwater have been distinguished in the Rincon and Galiuro Mts, and possibly in the Rosemont area of the Santa Rita Mts. Except in the Galiuro Mts, there is isotope evidence of mixing between old and young water. Preservation of ancient water in a mountain block is likely to be related to the depth and connectivity of fracture systems. In general, the depth of penetration of groundwater in mountain block aquifers is uncertain. Studies from regions of active tectonics have demonstrated groundwater flow in fractured rock to depths reaching 300 m below the surface [66] and references therein. In the Mont Blanc massif between France and Italy, 1950s or younger recharge has penetrated to a depth exceeding 2000 m on the basis of tritium measurements of water dripping into a highway tunnel [67]. Measurements of ³⁶Cl and tritium in a tunnel within the unsaturated zone at Yucca Mountain, Nevada, USA indicated post-bomb groundwater infiltration along fractures to depths as great as 500 m [68]. The presence of warm water seems consistent with deep flow in the Galiuro Mts. In the Tucson Mts., ancient water is produced from fractured rock up to 170 m below the surface. In the Rosemont area, close juxtaposition of wells producing young and ancient groundwater [33] suggests limited connectivity between sets of fractures.

Climate is also implicated in the persistence of ancient groundwater. O and H isotope data like those for ancient water in this study are interpreted as evidence that ancient water was recharged in the late Pleistocene, at a time of cooler and wetter climate in the Southwest USA [69]. Recharge was more abundant than at present under such conditions. In the low-elevation Tucson Mts., present rainfall may be insufficient to bring about recharge, leading to the dominance of ancient groundwater in much of that mountain block.

4.2. Multiple Recharge Mechanisms

The patterns of stable isotope data described above result from multiple recharge mechanisms affecting young groundwater, in combination with mixing between young and ancient groundwater. The data patterns may be grouped as follows. (i) Data plotting on a modified present-day LMWL defined by seasonal amount-weighted means for the wettest months; (ii) data falling along the present-day LMWL drawn through amount-weighted means for all summer and winter precipitation, and plotting between the means; (iii) data falling along evaporation trends; (iv) data with no relationship to present-day LMWLs or evaporation trends.

Group (i) matches Pattern 2 of groundwater in alluvial basins [45], and process B (Figure 1). It occurs at Rosemont in the Santa Rita Mts., and in recent recharge in the Tucson Mts., in both cases at altitudes below 1700 masl. This data pattern arises when recharge occurs during the wettest 30% (approximately) of months and reflects the loss of surface water to evaporation during drier months. At Rosemont, equal amounts of summer and winter recharge are typical (Figure 8), and in the Tucson Mts., more summer than winter recharge (Figure 7). Mixing of recent recharge with ancient groundwater is indicated in both cases. Groundwater from springs discharging at elevations below 1400 masl in the Santa Catalina Mts. also largely fits this pattern, relative to mean isotope compositions for precipitation at 1400 masl (Figure 5a).

Group (ii) is characteristic of groundwater from elevations above 1800 masl in the Santa Rita and Rincon Mts. (Figure 6 and Figure 8), and matches process A (Figure 1). Precipitation from wettest months at such altitudes is collinear with the LMWL for all precipitation (Figure 6). Both summer and winter recharge are indicated in both ranges, but the proportions remain uncertain because of the collinearity of the LMWL and the modified LMWL, and the possibility of summer recharge mainly during the wettest months.

Group (iii) corresponds to Pattern 1 recharge in alluvial basins of Central and Northwestern Arizona [45] and process C (Figure 1). The pattern reflects recharge of winter precipitation only, typically evaporated, in a region with no altitude dependence of O and H isotopes in winter precipitation. Group (iii) is observed in the Mogollon Highlands (Figure 10) and the Galiuro Mts. until 1996, but possibly not since 2012 (Figure 9). In both cases, the evaporation trend passes through the point (δ¹⁸O, δ²H) = (−11.2, −77‰), corresponding to mean winter precipitation in the Pattern 1 region [45]. The persistent evaporation trend in the Sacramento Mts. [17] also originates in winter precipitation, in this case with (δ¹⁸O, δ²H) near (−13.5, −87‰); after wet summers, however, isotope compositions move towards the LMWL (Figure 11). The evaporation trend at Rosemont in the Santa Rita Mts. originates in a mixture of summer and winter recharge (Figure 8). The evaporation trend in the Rincon Mts. originates mainly as high-elevation summer precipitation (Figure 7). In both cases, the isotope data are consistent with multiple cycles of infiltration, discharge into streambeds, evaporation and renewed infiltration.

Group (iv) comprises other patterns, and includes ancient groundwaters, as established by ¹⁴C data in the Rincon, Tucson and Galiuro Mts. (Figure 6, Figure 7 and Figure 9; Supplementary Table S2). Transmission of such water into basin alluvium (process D of Figure 1) has been documented near the Rincon Mts. [65].

In the Santa Catalina Mts., Group (iv) also includes spring waters from 1400 to 2300 masl (Figure 5a), many of which cannot be derived from long-term mean winter precipitation, regardless of the altitude of recharge (Figure 5b). Precipitation from certain events or seasons (Figure 3b,c) might account for the very low values of (δ¹⁸O, δ²H) in springs between 1400 and 2300 masl, but an explanation of that nature is problematic because it requires preferential infiltration of water of extreme isotope composition. A similar infiltration phenomenon was identified, but not satisfactorily explained, in soil water in the Marshall Gulch catchment, at 2300 to 2600 masl in the Santa Catalina Mts. [70]. Infiltration may be favored by progressive wetting of soil during a precipitation season, later precipitation events being more likely to survive evapotranspiration and lead to recharge. However, precipitation events with low (δ¹⁸O, δ²H) values occur randomly in both summer or winter rainy seasons in the Santa Catalina Mts. (Supplementary Table S2). Infiltration may also be favored by precipitation intensity (Figure 4b), as explained above. The highest precipitation intensity in summer corresponds with bulk monthly values of (δ¹⁸O, δ²H) as low as (−10.3, −67‰) in rainwater, insufficient to explain the observed values in spring water. In winter, precipitation intensity fails as an explanation because highest precipitation intensity corresponds with monthly (δ¹⁸O, δ²H) values that are higher than the long-term winter mean.

4.3. Seasonality of Recharge

As in surrounding alluvial basins [45], recharge in mountain blocks is locally dominated by winter precipitation, by summer precipitation, or by subequal mixtures. A particular recharge seasonality may be characteristic of an entire mountain range (e.g., Mogollon Highlands, Galiuro Mts. in the 1990s) or part of a range (Rincon Mts). Alternatively, recharge seasonality may depend upon altitude (Santa Rita Mts., Santa Catalina Mts.). Altitude-dependence may reflect the tendency of winter snow to accumulate in banks and melt in bulk at winter’s end at high altitudes, whereas snow from individual events melts separately throughout winter at lower levels. Such an effect may account for the greater importance of winter infiltration above 2400 masl in the Sacramento Mts. (Figure 11). It may also explain the relative isotope uniformity of spring water (near mean winter precipitation) above 2300 masl in the Santa Catalina Mts., contrasting with the broad range of δ¹⁸O at 1400–2300 masl (Figure 5a), but leaves unexplained the implied preferential recharge from low− δ¹⁸O events.

4.4. Influence of Geology

All ranges studied consist predominantly of silicate rocks, except for the carbonate-dominated Sacramento Mts. Unhindered infiltration into the highly permeable limestones of the Sacramento Mts. is reflected in the short-term changes in isotope data in that range. Thick soil profiles in typical high-elevation valleys appear to play a role in generating the dry-period evaporation trend in groundwater of the range [17]. A contrast in permeability is present at the crest and on the north slope of the Santa Catalina Mts, where tilted strata overlie less fractured, less permeable granite. Most of the springs of the range discharge from these strata [26]. Whether this geological environment is implicated in the unusual distribution of groundwater isotopes with altitude (Figure 5b) is not known.

The very different isotope patterns in groundwater from the Santa Catalina and Rincon Mts. (Figure 5a and Figure 6) occur in proximate ranges with closely similar lithology. The greater influence of summer recharge in the Rincon Mts. may arise from structural geology, the stepped topography governed by intrusion style at the range crest enhancing retention of summer runoff more effectively than the V-shaped canyons of the Santa Catalina Mts.

In the Tucson Mts., the Galiuro Mts., and the Mt Wrightson massif of the Santa Rita Mts., groundwater is hosted by fractured volcanic rock. In these cases, geology does not appear to impose a single pattern on isotope data in groundwater, and altitude, as discussed above, may be the primary control. The evaporation trend of groundwater in the Galiuro Mts. is like that in the Mogollon Highlands. Both mountain blocks are contiguous with, and behave as part of, a regional zone in North-Central Arizona in which winter recharge is dominant, and no isotope altitude effect is present [45].

4.5. Future Research

The complexity of (δ¹⁸O, δ²H) data patterns and recharge mechanisms in the mountain blocks reviewed here stands in sharp contrast to the simpler zonation of neighboring alluvial basins [45]. From the discussion above, it emerges that hydrologic processes in the mountain blocks reflect a complicated interplay of factors including precipitation seasonality, altitude, precipitation intensity, groundwater age and geology. The recharge mechanisms described in this article operate within a small number of mountain blocks that lie within a small portion of the North American Cordillera. In that area, the isotopic distinction between winter and summer precipitation provides a useful means of evaluating groundwater isotope data. Detailed studies of groundwater in other mountain blocks, both in the region of subequal winter and summer precipitation and beyond, will add to the understanding of mountain block hydrogeology. It is impossible to predict, given the complexity of results from the present study, the recharge mechanisms that will operate in neighboring mountain blocks, let alone in more distant examples. Detailed study of each mountain block will be required. Such research will benefit greatly from additional detailed information about isotope lapse rates in precipitation. The two detailed datasets [28,31] for sites in Arizona reveal intricacies that can be resolved only through long-term commitment to the gathering of data.

5. Conclusions

• In 1995–2004, precipitation at 2420 masl at PRS in the Santa Catalina Mts. had seasonal amount-weighted averages of (δ¹⁸O, δ²H) and LMWLs as follows: in summer (−8.3, −53‰), and δ²H = 7.0δ¹⁸O + 4.2; in winter, (−10.8, −70‰) and δ²H = 8.3δ¹⁸O + 20.1. Wetter months have lower (δ¹⁸O, δ²H) in summer, but higher (δ¹⁸O, δ²H) in winter.
• In 2008–2013, snow-bank samples defined a line δ²H = 9.0δ¹⁸O + 26.3; in 2007, isotope metamorphism gave rise to (δ¹⁸O, δ²H) values deviating strongly from this line.
• Long-term coefficients for isotope altitude effect (δ¹⁸O, δ²H) between Tucson and PRS are (−1.6, −7.7 ‰ per 1000 m) in summer, and (−1.1, −8.9 ‰ per 1000 m) in winter. Altitude effect was absent in 25% of seasons examined.
• Groundwater containing finite tritium predominates in all ranges studied except the Tucson Mts., and is interpreted as young groundwater, or mixtures containing young groundwater. Ancient groundwater, indicated by ¹⁴C and distinctive (δ¹⁸O, δ²H) data, is present in the Galiuro, Tucson, Rincon and Santa Rita Mts.
• Multiple recharge mechanisms of young groundwater are indicated by patterns of (δ¹⁸O, δ²H) data. The patterns are influenced by altitude, geology and mixing with ancient groundwater.
• Regional, zoned recharge mechanisms of neighboring alluvial basins [45] are present in the silicate mountain blocks. Winter-only recharge with evaporation occurs in the Mogollon Highlands and the Galiuro Mts. Infiltration of winter and summer infiltration from the wettest months, in various proportions, occurs below 1700 masl in the Santa Rita and Tucson Mts, and below 1400 masl in the Santa Catalina Mts.
• Additional recharge mechanisms include: (1) infiltration of both winter and summer precipitation in various proportions at high elevations in the Rincon and Santa Rita Mts., yielding (δ¹⁸O, δ²H) data plotting along high-altitude LMWLs; (2) Infiltration of high-elevation summer-dominant precipitation, with strong evaporation, in the Rincon Mountains; (3) recharge of water with (δ¹⁸O, δ²H) values below those of long-term mean winter precipitation above 1400 masl in the Santa Catalina Mts.; and (4) alternation of winter recharge with evaporation (dry years) and recharge of mixed winter and summer precipitation without evaporation (wet years) in the carbonate rocks of the Sacramento Mts.
• A complex interplay of isotope effects related to altitude, precipitation seasonality, precipitation intensity and groundwater age is responsible for the multiple patterns of isotope data.
• Assuming that mountain runoff from any altitude has an equal likelihood of reaching the basin floor, mountain-front recharge to alluvium in Tucson Basin consists of 50%–90% winter runoff from surrounding mountain blocks.