The Origins of Sulfate in Cenozoic Non-Marine Evaporites in the Basin and-Range Province, Southwestern North America

: Cenozoic evaporites (gypsum and anhydrite) in southwestern North America have wide ranges of δ 34 S ( − 30 to +22‰; most +4 to +10‰) and δ 18 O SO4 (+3 to +19‰). New data are presented for ﬁve basins in southern Arizona. The evaporites were deposited in playas or perennial saline lakes in closed basins of Oligocene or younger age. Very large accumulations in Picacho, Safford and Tucson Basins have isotope compositions plotting close to a linear δ 34 S-δ 18 O SO4 relationship corresponding to mixing of two sources of sulfur: (1) sulfate recycled from Permian marine gypsum and (2) sulfate from weathering of Laramide-age igneous rocks that include porphyry copper deposits. In the large evaporites, sulfate with δ 34 S > +10‰ is dominantly of Permian or Early Cretaceous marine origin, but has locally evolved to higher values as a result of bacterial sulfate reduction (BSR). Sulfate with δ 34 S < − 10‰ formed following exposure of sulﬁdes, possibly formed during supergene enrichment of a porphyry copper deposit by BSR, and have values of δ 18 O SO4 higher than those of local acid rock drainage because of participation of evaporated water in BSR. Accumulations of 30 to 100 km 3 of gypsum in Picacho and Safford Basins are too large to explain as products of contemporaneous erosion of Permian and Laramide source materials, but may represent recycling of Late Cretaceous to Miocene lacustrine sulfate.


Study Context and Scope
Sulfate minerals, gypsum and anhydrite, commonly occur in terrestrial and marine evaporites as bedded deposits in association with limestone, shale, clay and halite. Environmental factors, including climate, location, topography, and geology, control the supply and enrichment of solutes and determine the type of precipitates formed in closed terrestrial basins. Gypsum formation is commonly limited by the availability of sulfur. Hence, one needs to understand the geochemical cycling of sulfur in order to account for the origins of evaporites [1].
An important aspect of evaporite research centers on their potential use as records of seawater chemistry, tectonic setting, paleoclimate, and palaeohydrology. Sulfur and oxygen isotope compositions in sulfate are useful for identifying sulfate sources and inferring the geochemical environment of sulfate-bearing water [2,3]. The δ 34 S and δ 18 O SO4 age curves in marine evaporites [4] have been extensively used to date marine evaporites and to determine the marine or non-marine origin of evaporite deposits [5][6][7]. For instance, δ 34 S and δ 18 O SO4 have been used to determine the origins of sulfur in the evaporites in Paris Basin, France [8], southern Australia [9] and the Namib Desert, Namibia [10].
Non-marine evaporites are largely limited to the Cenozoic; nearly all older evaporites are of marine origin, reflecting the poorer preservation potential of lacustrine sedimentary rocks [11]. Non-marine evaporites generally have lower δ 34 S values than marine evaporites because of addition of low-δ 34 S sulfur from rainwater, sulfide minerals, volcanic sulfur,

Picacho Basin
Picacho Basin is one of several basins in central Arizona known collectively as the Gila Low, a region of internal drainage early in the formation of the Basin-and-Range Province [28]. Picacho Basin accumulated over 3000 m of Neogene continental sediments [29]. Fine sand, silt, clay, and evaporites were deposited in playas or ephemeral lakes [30]. The Exxon State (74)-1 well penetrated about 1800 m of massive anhydrite containing minor interbeds of shale, tuff, halite, and limestone nodules. The evaporites are of middle to late Miocene age, 14.9 to 10.5 Ma [20].

Tucson Basin
Basin fill consisting of conglomerate, sandstone and gypsiferous mudstone of Oligocene to Holocene age, is up to 3000 meters thick in the central graben of the basin [21].

Picacho Basin
Picacho Basin is one of several basins in central Arizona known collectively as the Gila Low, a region of internal drainage early in the formation of the Basin-and-Range Province [28]. Picacho Basin accumulated over 3000 m of Neogene continental sediments [29]. Fine sand, silt, clay, and evaporites were deposited in playas or ephemeral lakes [30]. The Exxon State (74)-1 well penetrated about 1800 m of massive anhydrite containing minor interbeds of shale, tuff, halite, and limestone nodules. The evaporites are of middle to late Miocene age, 14.9 to 10.5 Ma [20].

Tucson Basin
Basin fill consisting of conglomerate, sandstone and gypsiferous mudstone of Oligocene to Holocene age, is up to 3000 meters thick in the central graben of the basin [21]. Gypsum or anhydrite evaporites occur in the Oligocene (Pantano Formation) to Miocene (Tinaja Beds) strata, both at depth and at the basin margins. In the Exxon State (32)-1 well the interval 546 to 686 m contained an abundance of gypsum and anhydrite crystals [20,21]. In the Pantano Formation, small nodules of gypsum were sampled from outcrop at Marsh Station at the southeastern margin of the basin, and massive gypsum from a disused gypsum quarry at the northern margin. Gypsum samples from the Tinaja beds were obtained from cuttings from three drill holes near the center of the basin.

San Pedro Basin
Commercial gypsum deposits are quarried at Mammoth in the Miocene-Pliocene Quiburis Formation, which accumulated in the closed San Pedro Trough between 11 and 5 Ma ago [31]. The northern part of the Formation, the Redington Member, consists of clayrich sediments with local massive gypsum beds. Three additional sites of different age were sampled. In the upper San Pedro Basin, the alluvial-lacustrine St. David Formation [32] contains minor gypsum of probable Pliocene or Pleistocene age in outcrop about 8 km southwest of St. David. At Cascabel, a small lacustrine clay deposit of probable Holocene age next to the San Pedro River contains abundant, dispersed crystals of a sodium sulfate mineral. On the east flanks of the San Pedro valley at Cascabel, the Mineta Formation, a volcaniclastic redbed conglomerate of probable Oligocene age [33] contains a lens of recrystallized gypsum up to 1 m thick in a lacustrine member.

Safford Basin
Seismic reflection profiling [34,35] revealed the basin to be a half-graben, with the southwest side downfaulted along a secondary breakaway fault. Basin-fill consists of conglomerate, silt, clay and evaporites with minor lacustrine limestone and diatomite [36,37]. In the basin center, the Sanchez Beds consist of 250 m of clay, gypsum/anhydrite, and halite, and are overlain by 520 m of Pliocene lacustrine units. Fine-grained Pliocene lacustrine sediments and alluvial fan deposits, the Bear Springs Wash Beds (BSWB), interfinger with the Sanchez beds and older basin fill in the southwest part of the basin.
A 1971 drill hole (Tenney #3 State well) about 32 km south of Safford, encountered 700 m of gypsum and anhydrite below 366 m and terminated in evaporites at 1068 m. The full thickness is unknown [22]. Gypsum/anhydrite of the Sanchez Beds was taken from three sets of well cuttings, from the Tenney #3 State, Claridge and City of Safford wells. In addition, gypsum was collected from the BSWB in Frye, Hall and Riggs Mesas along the southwestern basin margin, where outcrop includes red and green shales with scattered thin beds of gypsum, pyritic shale and white tuff. They remain exposed along the mountain front in an area where breakaway faulting has migrated toward the basin center. Sulfate was also extracted from lacustrine clay and ocher outcrops along the Gila River between Safford and Geronimo.

Unnamed Basin, Douglas
The basin is located about 8 km east of Douglas, Arizona. Two disused quarries expose irregular, near-surface lenses of light buff to white gypsum and gypsite up to 2 m thick. The lenses overlie Cenozoic, possibly Pleistocene, clay. Coarse yellow crystals of gypsum occur locally at the base of the soil. An active sulfate source is apparent: in D Hill, 1 km south of the quarries, altered Paleozoic limestone contains goethite casts of coarse pyrite crystals, some with sulfide cores. Small gullies drain the north slope of D Hill and lead to the quarries. The gypsum appears to have formed where sulfate-bearing water from D Hill has evaporated on infiltrating basin-fill sediment in a small basin that is open to the north. Gypsum, pyrite from the limestone and sulfate from dry stream sediment (representing integrated sulfate supplied in runoff) were sampled.

Modern Playas
Animas Playa, New Mexico, Willcox Playa, Arizona and Playa Guzmán, Chihuahua are ephemeral lakes within internal drainage basins. The Salton Sea, California, has been a permanent lake since the early 20th Century.

Environment of Evaporite Formation
The terrestrial Cenozoic setting of the province indicates that most of the evaporites accumulated in isolated interior-drainage basins, mostly likely in environmental settings like modern playas and saline lakes [14,16]. Lacustrine conditions ended in most areas between 10.5 and 5 Ma [20,31].
Evaporitic basins are hydrologically closed; that is, evaporation exceeds inflow under an arid to semi-arid climate regime. Evaporite deposits can form only where sufficient supply of solutes is available in inflow. These seemingly conflicting requirements are reconciled in orographic deserts, where high mountains catch precipitation and generate rain shadows on neighboring basin floors [12]. In this regard, the Basin-and-Range Province has provided ideal conditions in the past for the formation of lacustrine evaporites. The thicker known evaporites occur near the northeastern margin of the Province, proximal to the Transition Zone bordering the Colorado Plateau.
Deposition of great thicknesses of evaporites in short time intervals is well documented. High evaporation in arid regions permits accumulation of sulfates at rates of 1-40 m/1000 years and halite at 10-100 m/1000 years [38,39]. At such rates, the thick anhydrite evaporites in Picacho basin could form in 45,000 to 1.8 million years.

Origins of Sulfur in Non-Marine Evaporites
Sulfur cycling and the behavior of sulfate isotopes in surface environments have been reviewed in [40] and references therein. In non-marine basins, numerous sources of sulfur and numerous pathways of sulfur supply are possible. Rainwater, dust, surface water, and groundwater can all transport the sulfate from source areas to the lacustrine environments where sulfate evaporites form. Possible sulfate sources are: (1) marine aerosols (cyclic salts) and dimethysulfide (DMS); (2) dissolution of pre-existing evaporites within the drainage basin; (3) sulfate in connate water or dispersed in non-evaporitic marine sediments; (4) weathering of sulfide-bearing rock; (5) biogenic reduced sulfur; and (6) anthropogenic sulfur, including soil amendments, fertilizers and combustion products. Compiled ranges of δ 34 S from studies of non-marine evaporites published prior to 1994 have been compared with the ranges for marine sulfate and DMS [40]; the authors noted that single or mixed multiple sources of sulfur may apply in a single basin.
Sulfate generated by oxidation of sulfide commonly has δ 34 S values close to those the sulfide and low values of δ 18 O SO4 [40]. Values of δ 18 O SO4 result from the incorporation into sulfate of O from water in addition to atmospheric O 2 , e.g., according to the following reactions for oxidation of pyrite [41]: The second reaction is rate-limited under acid conditions by the reaction which is greatly accelerated by Thiobacillus ferrooxidans in submersed and anaerobic environments [42,43]. Atmospheric O 2 , with δ 18 O +23.5‰ [44] differs greatly in δ 18 O from meteoric water; therefore Equations (1) and (2) should yield sulfate with distinctive ranges of δ 18 O SO4 . Both reactions may occur concurrently, but in acid mine drainage in the western USA, Equation (2) typically dominates [43]. In Miocene gypsum from southern Poland, Equation (1) appears to dominate; sulfate derived by oxidation of coal-bed reduced sulfur has δ 18 O SO4 values of +18 to +20‰, approaching that of atmospheric O 2 [45].

Evolution of δ 34 S and δ 18 O in Sulfate Evaporites
Under lacustrine conditions, no oxygen isotope exchange occurs between sulfate ions and water [46]. Aqueous sulfate in lake water preserves the δ 34 S and δ 18 O SO4 of its sources unless the sulfate is subjected to one of the following processes.
Isotope fractionation occurs between gypsum and coexisting aqueous sulfate. Under evaporitic conditions, ∆ 34 S and ∆ 18 O SO4 (gypsum-aqueous sulfate) are about +1.65‰ and +3.5‰, respectively [47][48][49]. In a small sulfate reservoir, early-deposited gypsum has higher values of δ 34 S and δ 18 O SO4 than later-crystallized gypsum. If cyclic dissolution and recrystallization of gypsum are accompanied by removal of some aqueous sulfate, fractionation may accumulate in the residual solid. Such cumulative fractionation may occur in vadose playas, as suggested for Cl isotopes [50]. Gypsum dehydration and diagenetic redistribution of sulfate do not generate significant isotope fractionation [51].
Bacterial sulfate reduction (BSR) and thermochemical sulfate reduction (TSR) produce sulfide species from dissolved sulfate by reaction with organic matter [52]. BSR occurs in near-surface environments and water columns at temperatures up to 80 • C [53] or possibly higher [52]. TSR, which occurs at temperatures above 100 • C [53] is improbable in the depositional environment of evaporites or in weathering profiles, but is possible where buried gypsum recrystallizes as anhydrite at temperatures above 90.5 • C [54].
Partial BSR causes increases in both δ 34 S and δ 18 O SO4 in residual sulfate in experiments utilizing single bacterial species and single organic-compound substrates, provided that ambient water is not greatly depleted in 18 O. Change in δ 18 O SO4 appears to arise from isotope exchange between water and S 4+ reaction intermediaries, and produces scattered results, even in simple experimental situations; it does not show a simple relationship with change in δ 34 S [55,56]. Natural BSR by multiple microbial species exploiting complex organic substrates is less likely to produce predictable relationships between δ 34 S and δ 18 O SO4 in residual sulfate. While BSR can readily fractionate sulfate isotopes in aqueous sulfate, it is unlikely to shift δ 34 S of pre-existing massive gypsum or anhydrite. Minute pyrite crystals surrounding gypsum crystals buried beneath Bristol Dry Lake, California [57] represent BSR of sulfate dissolved from the gypsum, but the BSR is unlikely to have affected the δ 34 S of the remaining gypsum.

Previous Geochemical and Isotope Studies of Nonmarine Evaporites in Southwest North America
Low Br content (2-6 ppm) suggested a nonmarine origin for Miocene halite in Luke Basin, Arizona [23]. Modern perennial saline lakes in the western USA have a δ 18 O SO4 range of +13.3 to +23.2‰ [58]. In the Bristol Dry Lake basin, California, the range of δ 34 S values of near-surface (<2 m) gypsum, anhydrite, and celestite was +6 to +9‰, while 12 core samples (286 to 517 m) had lower δ 34 S values in the range 0 to +4‰ [57]. In the Hualapai basin halite body, Arizona, δ 34 S and δ 18 O SO4 values were +6.8‰ and +5.9‰, respectively, in an intercalated anhydrite layer, and +11.3‰ and +11.5‰ in the anhydrite cap [16]. The δ 34 S value of the intercalated anhydrite is significantly lower than those of Phanerozoic marine sulfate [4], indicating a non-marine origin. A 200,000-year record of gypsum, glauberite and mirabilite in Death Valley, California shows δ 34 S values from +5 to +22‰, and δ 18 O SO4 from +9 to +23‰ [59]. The δ 18 O SO4 values probably reflect changes in lake level, reduction of sulfate and re-oxidation of sulfide and/or differences in crystallization temperature of the sulfate minerals. Ephemeral lake conditions (warm, shallow water) resulted in precipitation of sulfate with lower δ 18 O SO4 values. Wet, cold, stratified, and perennial lake conditions resulted in precipitation of sulfate with high δ 18 O SO4 values because of re-oxidation of sulfide formed by BSR in the lower water or bottom mud. At the White Sands and Estancia gypsum dune fields, New Mexico, δ 34 S values imply recycling of Permian marine gypsum, with superimposed BSR at Estancia [60].

Sulfate Isotopes in Precipitation
Anthropogenic emissions have accounted for 75-84% of sulfur in the atmosphere of the Northern Hemisphere since 1980 [61,62]. In the southwest USA, major sources have been copper smelters and coal-fired electric power stations. Rainwater sulfate collected in Tucson Basin between 1996 and 2003 [63,64] had a δ 34 S range of +2.1 to +8.0‰ with a mean, weighted for precipitation amount, of +4.0‰. A small dataset for δ 18 O SO4 ranged from +8.4 to +12.0‰, with an average of +10.2‰. Sulfate in pre-industrial precipitation, lacking anthropogenic input, contained mostly sulfur of volcanic and marine biogenic sources [61], except in desert regions where wind-blown dust may contribute sulfate. In the study area, detectable marine sulfate aerosol appears to penetrate only to windward slopes of mountain ranges bordering the Gulf of California, on the evidence of SO 4 /Cl ratios in groundwater [65]. Desert varnish sampled as a proxy for pre-industrial atmospheric sulfate had average values of +7.9 ± 0.9‰ for δ 34 S and +6.4 ± 1.1‰ for δ 18 O SO4 in sulfate [66]; the authors determined means of +5.8 ± 1.4‰ (δ 34 S) and +11.2 ± 1.1‰ (δ 18 O SO4 ) for modern atmospheric sulfate. In Tucson Basin, young groundwater (containing detectable tritium) had lower values of δ 34 S than pre-bomb groundwater [67], indicating a sense of change like that of the desert varnish data.

Bedrock Sulfate and Sulfide
Available δ 34 S and δ 18 O SO4 data for bedrock, including PCDs, in the study area are summarized in Table 1. Sulfate collected from springs or wells drilled in bedrock is assumed to represent the sulfur disseminated in the bedrock. Both δ 34 S and δ 18 O SO4 in the bedrock have wide ranges, from −6 to +25‰ and −2 to +13‰, respectively. Bedrock types important in this study area are: (a) Late Paleogene felsic volcanic rocks and associated granitoids, low in sulfur content. No S isotope data are available.
(c) Laramide (Late Cretaceous-Paleocene) intermediate volcanic rocks and associated granitoids with numerous large PCDs. PCDs and overlying cap rock have high concentrations of sulfur. Wells drilled in these strata generally produce groundwater with high sulfate concentration. In the Patagonia Mountains, wells, springs and acid rock drainage generally have high concentrations of sulfate (60-3000 ppm); δ 34 S ranges from −10.7 to +0.1‰, and δ 18 O SO4 from −2.2 to +3.4‰ [63,69,70]. In the Tucson Mountains, δ 34 S in groundwater hosted by Laramide igneous rock ranges from +1 to +3‰, and δ 18 O SO4 from +2 to +7‰ [66]. Low δ 18 O SO4 values in both cases result from oxidation of sulfide.
The mass ratio of hypogene sulfate to hypogene sulfide may be about 1:2 in the PCDs of southern Arizona (E. Seedorff, personal communication). Assuming mean δ 34 S values of 0‰ in sulfide and +10‰ in hypogene sulfate (averages of available data), then the δ 34 S of bulk sulfate weathered from PCDs would be +3.3‰. The Mineral Park deposit has more positive mean hypogene δ 34 S values than other PCDs in the region, +2.4‰ for sulfide and +14.7‰ for sulfate [73]; in this case the bulk δ 34 S would be +6.5‰. It is difficult to estimate the total sulfur content of an entire PCD system, including ore, capping rock and lateral disseminations. Giant PCD orebodies may contain up to 10 9 t of sulfide-S, and 1.5 × 10 9 t total S [74,75]; for comparison, the Morenci deposit in Arizona may have contained 0.24 × 10 9 t S in supergene-enriched ore [72]. Supergene enrichment (and therefore erosion) of PCDs in southern Arizona began in the Mid-Miocene at 13 Ma and was most intense from 11 to 7 Ma [72,76] contemporaneous with the formation of deep evaporite basins and deposition of massive evaporites.
(d) Mesozoic sedimentary rocks. These generally have low sulfur contents; e.g., 30 to 285 ppm in Cretaceous sedimentary rocks in southern Arizona [77]. The last marine incursion in southern Arizona resulted in the deposition of the Bisbee Group of Early Cretaceous age [78]. Corresponding marine sulfate (δ 34 S + 15 ± 2‰, δ 18 O SO4 + 13.5 ± 1.5‰) [4] may have survived as traces of evaporitic gypsum in a coastal lagoonal facies [79] or in connate water. In a lacustrine unit of the terrestrial Late Cretaceous Cabullona Group [80] near Huépac, northern Sonora, Mexico, δ 34 S ranges from +10 to +16‰, and δ 18 O SO4 from +6 to +15‰ [81]. of bulk sulfate weathered from PCDs would be +3.3‰. The Mineral Park deposit has more positive mean hypogene δ 34 S values than other PCDs in the region, +2.4‰ for sulfide and +14.7‰ for sulfate [73]; in this case the bulk δ 34 S would be +6.5‰. It is difficult to estimate the total sulfur content of an entire PCD system, including ore, capping rock and lateral disseminations. Giant PCD orebodies may contain up to 10 9 t of sulfide-S, and 1.5 × 10 9 t total S [74,75]; for comparison, the Morenci deposit in Arizona may have contained 0.24 × 10 9 t S in supergene-enriched ore [72]. Supergene enrichment (and therefore erosion) of PCDs in southern Arizona began in the Mid-Miocene at 13 Ma and was most intense from 11 to 7 Ma [72,76] contemporaneous with the formation of deep evaporite basins and deposition of massive evaporites.
(e) Paleozoic marine strata. The Permian Epitaph Formation, containing the oldest known evaporitic gypsum-anhydrite in southern Arizona, formed in a restricted basin transitional between marine and continental environments [82]. A gypsum sample from the Whetstone Mountains gave δ 34 S and δ 18 0SO4 values of +11.7‰ and +12.8‰, respectively. Gypsum from the Permian Kaibab Formation in northern Arizona yielded a δ 34 S value of +13‰ [73]. These data are typical of Permian marine evaporites [4]. Pyrite from Paleozoic sulfidic limestone on D Hill near Douglas has a δ 34 S range of +3 to +21‰. This limestone is heavily altered; therefore, the pyrite may be secondary.
(e) Paleozoic marine strata. The Permian Epitaph Formation, containing the oldest known evaporitic gypsum-anhydrite in southern Arizona, formed in a restricted basin transitional between marine and continental environments [82]. A gypsum sample from the Whetstone Mountains gave δ 34 S and δ 18 O SO4 values of +11.7‰ and +12.8‰, respectively. Gypsum from the Permian Kaibab Formation in northern Arizona yielded a δ 34 S value of +13‰ [73]. These data are typical of Permian marine evaporites [4]. Pyrite from Paleozoic sulfidic limestone on D Hill near Douglas has a δ 34 S range of +3 to +21‰. This limestone is heavily altered; therefore, the pyrite may be secondary.
(g) Lower Proterozoic strata and intrusions. A sample of groundwater from Paleoproterozoic granite at Payson, Arizona gave a δ 34 S value of +0.6‰, which may represent igneous sulfur; higher δ 34 S values in other groundwater samples from granite near Payson appear to reflect dissolved Permian gypsum [83]. In central Arizona, lower Proterozoic volcanogenic massive sulfide deposits have a δ 34 S range of 0 ± 3‰, probably representing contemporaneous igneous sulfide from the host felsic volcanic suite [84].

Materials and Methods
Sulfate minerals and sulfate-rich sediment samples were collected from outcrop, well cuttings and drill core. Single water samples were collected from the modern playas, except for Willcox Playa, where gypsiferous dune sand was collected. Sample locations are shown in Figure 1. Detailed sample locations for Tucson and Safford Basins and the site near Douglas are given as supplementary Figures S1-S3. Approximately 300 mg of gypsum/anhydrite was dissolved in warm 1N HCl. Water samples were acidified with 1N HCl to remove carbonate species. The warm solutions were filtered to remove any solid residues, and sufficient 10% BaC1 2 solution was added to precipitate all sulfate as BaSO 4 ; solutions were allowed to stand for an hour to complete the precipitation. The BaSO 4 was filtered, washed thoroughly with deionized water and dried. Both δ 34 S and δ 18 O were analyzed in the Environmental Isotope Laboratory, University of Arizona. For S isotopes, 15 mg of BaSO 4 was mixed with 60 mg each of Cu 2 O and SiO 2 , combusted under vacuum at 1100 • C, and the resultant SO 2 was cryogenically purified [86]. Isotope ratios in SO 2 gas were measured on a modified VG602C gas source mass spectrometer. For O isotopes, 1 mg BaSO 4 was analyzed using a continuous-flow isotope ration mass spectrometer (Finnigan Delta X Plus) coupled with a thermal combustion elemental analyzer. Standardization is based on international standards NBS123 and OGS-1 (equivalent to NBS-127) [87]. Values of δ 34 S and δ 18 O SO4 are reported relative to Vienna Canyon Diablo Triolite (VCDT) and Vienna Standard Mean Ocean Water (VSMOW), with analytical precisions of 0.13‰ (1σ) and 0.4‰ (1σ), respectively, according to repeated analysis of laboratory standards.

Results
The δ 34 S and δ 18 O SO4 values of sulfate minerals, sulfate-rich sediment and water from playa lakes are listed in Supplementary Table S1. Figure 4 is a histogram of δ 34 S values, and Figure 5 is a plot of δ 34  range in gypsum sand is +12 to +14‰ [60], Estancia, +19 to +21‰ [60] and Death Valley, +14 to +19‰ [59]. Figure 4. Histogram of δ 34 S data for evaporites and playas, and data for rain and dust in Tucson Basin. Data are from this study, except: ranges for White Sands and Estancia Basin [60], Hualapai Basin [16], Bristol Dry Lake [57], Death Valley [59].
Many of the ancient evaporite sites have δ 34 S values between +4 and +10‰, including large-volume evaporites and gypsiferous sediment in Tucson, Picacho and Safford basins. Outliers with δ 34 S > +10‰ include large-volume occurrences in Safford Basin (Claridge well and part of the Tenney well) and sulfate-rich sediment in Death Valley; small-volume gypsum bodies in this bracket occur at Douglas (+19 to +22‰), Marsh Station in Tucson Basin (+12 to +13‰) and Cascabel in the San Pedro Basin (+11‰ and +16 to +17‰) Outliers with δ 34 S < +4‰ include large gypsum bodies at depth below Bristol Dry Lake, California (0 to +4‰) and in Safford Basin, gypsum intersected by the City of Safford well (mainly −20 to −10‰) and minor occurrences in the Bear Springs Wash Beds, (−30 to −10‰). Many of the ancient evaporite sites have δ 34 S values between +4 and +10‰, including large-volume evaporites and gypsiferous sediment in Tucson, Picacho and Safford basins. Outliers with δ 34 S > +10‰ include large-volume occurrences in Safford Basin (Claridge well and part of the Tenney well) and sulfate-rich sediment in Death Valley; small-volume gypsum bodies in this bracket occur at Douglas (+19 to +22‰), Marsh Station in Tucson Basin (+12 to +13‰) and Cascabel in the San Pedro Basin (+11‰ and +16 to +17‰) Outliers with δ 34 S < +4‰ include large gypsum bodies at depth below Bristol Dry Lake, California (0 to +4‰) and in Safford Basin, gypsum intersected by the City of Safford well (mainly −20 to −10‰) and minor occurrences in the Bear Springs Wash Beds, (−30 to −10‰). Geosciences 2021, 11, x FOR PEER REVIEW 11 of (a) (b) Vertical variations of δ 34 S and δ 18 0SO4 are present in the Tenney #3 State well, Safford Bas in which two separate gypsum layers appear to be present ( Figure 6). Although depths of w cuttings cannot be correlated from one part of Safford Basin to another, occurrences of gyps in the basin appear to indicate a complex evolution of isotopes over time. Vertical variatio of δ 34 S and δ 18 0SO4 are documented in Death Valley [59] (Figure 6) and of δ 34 S at Bristol D Lake [57]. For gypsum with δ 34 S > 0‰, δ 34 S and δ 18 0SO4 show a positive, linear correlation Safford Basin (Figure 5b). Elsewhere, δ 34 S and δ 18 0SO4 are scattered around the δ 18 0SO4 = δ line of the plot, much of the scatter arising in the minor gypsum occurrences at Ma Station and Douglas. In the large gypsum accumulations of Picacho and Tucson Basi some data plot near the line δ 18 0SO4 = δ 34 S, and some above the line.  Vertical variations of δ 34 S and δ 18 0SO4 are present in the Tenney #3 State well, Safford Basin, in which two separate gypsum layers appear to be present ( Figure 6). Although depths of well cuttings cannot be correlated from one part of Safford Basin to another, occurrences of gypsum in the basin appear to indicate a complex evolution of isotopes over time. Vertical variations of δ 34 S and δ 18 0SO4 are documented in Death Valley [59] (Figure 6) and of δ 34 S at Bristol Dry Lake [57].

Sources of Sulfur, Principal Controls
With two exceptions, the bulk values of δ 34 S in all evaporite deposits in the study area are lower than those of sulfate in Late Paleogene to Neogene marine evaporites, which have δ 34 S values ranging from +21 to +22‰ [4]. These values indicate that the evaporites are of non-marine origin. The exceptions, at Douglas and in the Estancia Basin [60], are forming at present in non-marine basins. No spatial zonation of δ 34 S is present at regional scale.
In Figure 5, the data plot close to or above the line δ 18 O SO4 = δ 34 S, except for the Douglas site where the sulfur source is local and unusual for the region (Section 4.4). The relationship with the line δ 18 O SO4 = δ 34 S can be compared with the mixing trend, also close to δ 18 O SO4 = δ 34 S, of sulfate near the Patagonia Mountains [88]. Acid rock drainage with low values of δ 34 S and δ 18 O SO4 mixes with basin groundwater containing reworked Permian marine sulfate in which δ 34 S and δ 18 O SO4 are both near +14‰ (Figure 5a). A similar trend is present in Tucson Basin groundwater [63]. Therefore, a principal regional control on sulfate isotopes appears to be the mixing of sulfur of Permian marine and bedrock sulfide sources. Mixing with marine sulfate is required, because the sulfur flux from PCDs has δ 34 S and δ 18 O SO4 less than about 6.5‰ (Section 1.8). Cyclic recrystallization (Section 4.4) may explain the data above the line δ 18 O SO4 = δ 34 S.

Sources of Sulfur, Picacho Basin
The 50 to 100 km 3 anhydrite body in Picacho Basin requires a very large source of sulfur. The range of δ 34 S in vertical profile, 4‰ (Supplementary Table S1), indicates relatively little change in the ratio of sulfur sources during deposition. Table 2 summarizes mass-balance requirements for sulfur, with different assumptions about sources. The average of available δ 34 S data, +7.8‰, might be explained by mixing about-equal amounts of sulfur from PCDs (bulk δ 34 S = +3.3‰) and Permian marine gypsum (+13‰). The quantities required, multiple giant PCDs and 50 to 100 km 3 of Permian gypsum, are unlikely to have been available in the immediate region around Picacho Basin as the evaporites accumulated, although larger volumes of Permian gypsum might have been available in the Colorado Plateau. If Laramide andesite containing 60 ppm S with δ 34 S = 0‰ supplied the low-δ 34 S endmember, about 40% of the total sulfur according to isotope balance, volumes of andesite of about 10 5 km 3 would be required if all S were extracted. In this case, a very large amount of clastic detritus would have to be accounted for.
Present-day rainwater in Tucson Basin contains 0.7 ppm S as sulfate [66]. If such rainwater collected over a catchment of 10 5 km 3 , comparable to the area of the Basin-and-Range province in southern Arizona, it would take 15 to 30 Ma to supply enough sulfur. The δ 34 S value of ancient rain, 7.9‰, matches the bulk δ 34 S of the evaporite body, but the chloride that would also be supplied in rainwater is missing. The large differences in δ 34 S among sites for this study indicate that local geology is more important than regional atmospheric deposition in generating the isotope characteristics of the evaporite deposits.
Singly or in combination, the sulfur sources do not satisfactorily explain the evaporite accumulation in Picacho Basin or in Safford Basin where 33 to 67 km 3 of gypsum may be present [22].

Sources of Sulfur, Safford Basin
Drill-core samples show stratification of δ 34 S and δ 18 O SO4 (Figures 5b and 6). In the Tenney #3 State well, δ 34 S values are +10 to +17‰ (type A) from 347 to 670 m below surface, and +8 to +10‰ (type B) below 670 m. Although sediment depths cannot be correlated with certainty between drill sites, type A gypsum in the Claridge and Tenney wells may be of the same age. Types A and B plot close to the line δ 18 O SO4 = δ 34 S (Figure 6), consistent with the mixing of sulfur of PCD and Permian or Early Cretaceous marine origins. In this case, the sulfur sources appear to have evolved over time, from marine + PCD in type B to marine alone in type A.
The City of Safford well (δ 34 S from −18 to −5‰), which penetrates down-faulted strata close to the basin center, and outcrop of the Bear Springs Wash Beds, BSWB (δ 34 S from −30 to +10‰) contain gypsum of a third type, C. The δ 18 O SO4 values range from +2.4 to +7.5‰ (Figure 5b), higher than values that arise from oxidation of sulfide elsewhere in the region, e.g., 0 ± 1‰ in Tucson Basin ( Table 1). The higher δ 18 O SO4 values in Safford Basin may reflect the participation of evaporated lake water [40]. In the BSWB, type C occurs as thin gypsum beds and as dispersed sulfate in clastic sediments and ocher along the northeastern flank of the Pinaleno Mountains (Supplementary Figure S2). Thin (a few cm) sulfide-bearing beds weathering to Fe sulfates are present in outcrop. Oxidation of similar sulfide appears to be the source of the gypsum with δ 34 S values < −25‰ in the BSWB. The range of δ 34 S in the City of Safford well cuttings can be explained by mixing of sulfate derived from sulfides with sulfate of type A and B. Gypsum in the well cuttings may include 25 to 30% of sulfide-sulfur according to isotope balance.
The amount of sulfide visible in outcrop of the BSWB is small and appears insufficient to account for the 267 m of low-δ 34 S massive gypsum in the City of Safford well. A possible explanation for the low-δ 34 S sulfate in the well cuttings may be found in the formation  and subsequent erosion of the supergene enrichment blanket overlying the giant Morenci PCD in the Gila River catchment about 40 km upstream of Safford Basin. Supergene enrichment is a near-surface weathering process occurring at about 25 • C near the water table [89]; therefore BSR is possible but not TSR [53]. The few measurements available for the deposit are consistent with primary sulfide δ 34 S near 0‰, and supergene sulfide with a δ 34 S range of −31 to 0‰ as a result of BSR, for which textural evidence is present [72]. At the PCD, supergene enrichment may initially have released sulfate with δ 34 S > 0‰ (contributing to types A and B) during the formation stage, leaving residual sulfide with low δ 34 S values in the supergene ore of the deposit. Subsequent erosion of the supergene ore may have released sulfate with δ 34 S < 0‰ (contributing to type C). The occurrence of type C gypsum at 95−367 m below surface, in contrast to types A and B at 329−1037 m (Supplementary Table S1) appears consistent with this explanation. Safford Basin is the only example in which a large-scale contribution of BSR to evaporite sulfate in the study area is probable.

Sources of Sulfur, Gypsum with δ 34 S > 10‰
Gypsum in the White Sands playa and dune-field, New Mexico, has δ 34 S values of 12 to 13‰ [60], readily explained by reworking of Permian marine gypsum in nearby strata [90]. In the Estancia gypsum dune-field, New Mexico, subjacent Permian gypsum is also the main source of sulfate [60], but dune sand has a δ 34 S range of +19 to +21‰. The authors attributed the increase in δ 34 S to BSR in the groundwater system feeding the Estancia playa.
The sulfur source for the small gypsum occurrence near Douglas, Arizona, is pyritic limestone in D Hill. Three samples of the pyrite had δ 34 S values of +3.8, +10.9 and +21.2‰. Dry stream sediment on the hill slopes preserves sulfate that provides a better estimate, +15.5 to +16.8‰, of average δ 34 S in runoff from the hill. The gypsum has δ 34 S values of +20.8 to +21.8‰. Values of δ 18 O SO4 increase as δ 34 S increases (Supplementary Table S1). BSR is unlikely, because gypsum is forming in near-surface sediment that is usually dry.   Table S1) appears consistent with this explanation. Safford Basin is the only example in which a large-scale contribution of BSR to evaporite sulfate in the study area is probable.

Sources of Sulfur, Gypsum with δ 34 S > 10‰
Gypsum in the White Sands playa and dune-field, New Mexico, has δ 34 S values of 12 to 13‰ [60], readily explained by reworking of Permian marine gypsum in nearby strata [90]. In the Estancia gypsum dune-field, New Mexico, subjacent Permian gypsum is also the main source of sulfate [60], but dune sand has a δ 34 S range of +19 to +21‰. The authors attributed the increase in δ 34 S to BSR in the groundwater system feeding the Estancia playa.
The sulfur source for the small gypsum occurrence near Douglas, Arizona, is pyritic limestone in D Hill. Three samples of the pyrite had δ 34 S values of +3.8, +10.9 and +21.2‰. Dry stream sediment on the hill slopes preserves sulfate that provides a better estimate, +15.5 to +16.8‰, of average δ 34 S in runoff from the hill. The gypsum has δ 34 S values of +20.8 to +21.8‰. Values of δ 18 0SO4 increase as δ 34 S increases (Supplementary Table S1). BSR is unlikely, because gypsum is forming in near-surface sediment that is usually dry. The occurrence of coarse-grained, clear gypsum crystals indicates recrystallization, possibly affecting some of the gypsum each wet season. A stream-sediment sample with δ 34 S and δ 18 0SO4 values of +15.5‰ and +3.3‰, respectively, may be close to unevolved sulfate formed by weathering of pyrite because of its low δ 18 0SO4 value. This sample provides the starting point for Rayleigh fractionation modeling of isotope fractionation arising from cyclic recrystallization (Figure 7). The array of samples from Douglas could be explained by repeated dissolution and crystallization of sulfate with δ 18 0SO4 near +3‰ and δ 34 S of +12 to +15.5‰. The residual gypsum samples would represent only about 2% of the original sulfate. Repeated loss of sulfate from the system is plausible because it is open to groundwater flow to the northwest. Figure 7. Plots of δ 18 0 (sulfate) vs. δ 34 S in gypsum, where fractionation is due to cyclic dissolution and recrystallization with loss of some aqueous sulfate in each cycle. The theoretical values (colored points) were calculated from chosen starting points, using fractionation factors from [47,48]. "f" represents the fraction of original sulfate remaining as gypsum.
At Marsh Station, most of the δ 34 S values plot close to those of Permian sulfate (+12 to +13‰), but values of δ 18 0SO4 are about +15‰. Inorganic recrystallization starting from (+12, +12‰) does not explain the δ 18 0SO4 values (Figure 7). In this case, gypsum occurs as nodules in gray mud, and BSR may have affected aqueous sulfate prior to gypsum formation, but from a starting material that was not Permian marine sulfate.  (Figure 7). In this case, gypsum occurs as nodules in gray mud, and BSR may have affected aqueous sulfate prior to gypsum formation, but from a starting material that was not Permian marine sulfate.
At Cascabel, data for gypsum from the Mineta Formation appear to form a recrystallization trend (Figure 7). The gypsum retains bedding structure but is strongly recrystallized and occurs in strata with evidence of hydrothermal circulation (bleaching of Fe 3+ oxides, local deposition of pyrite). The lowest values of δ 34 S and δ 18 O SO4 fall within the field of Early Cretaceous marine sulfate [4], a possible source. An Early Cretaceous marine transgression represented by Mural Limestone, part of the Bisbee Group, extended to a coastline at most 15 km south of the gypsum outcrop [91]. In the immediate area of the gypsum, the depositional environment of the Bisbee Group was subaerial [92]. Textural evidence of sulfate minerals in the Mural Limestone is rare [93], but survival of marine sulfate in connate water is not precluded.

Bristol Dry Lake and Death Valley
These occurrences differ from those in the Basin and Range Province in Arizona in the abundance of halite, and in the case of Death Valley [59], the presence of highly-soluble sulfate minerals glauberite and thenardite, exceeding gypsum in abundance. At Bristol Dry Lake, BSR was suggested to explain low values of δ 34 S in gypsum from drill core beneath the present playa, implicitly assuming constant δ 34 S in sulfate supplied to the basin [57]. BSR is implausible because it is difficult to decrease the δ 34 S of residual gypsum by that means. Gypsum in the sample set ranges from 3.5 Ma old to Holocene; a change in sulfate source over time seems a more likely explanation.
At Death Valley, δ 34 S of sulfate ranges from +5 to +22‰ (Figure 8) in sediment dating from the past 200 Ka [59]. The authors explained the range of δ 34 S by changes in the relative importance of three sources, atmospheric deposition (+7‰), local groundwater (+16 to +18‰) and regional groundwater (+21‰), the changes mediated by long-term variations in rainfall amount. Such an explanation may hold for long core intervals in which δ 34 S has small ranges, but is less satisfactory for other core intervals with large ranges of δ 34 S and fluctuations observed at a scale of 5 to 20 cm. BSR acting in small, closed-system cells may explain the latter core intervals. The sulfate occurs in soluble minerals and is or was available to bacteria. At 110 to 120 m depth, persistent δ 34 S values of +5 to +7‰ may correspond to the dominant sulfate source at the time (a period of halite deposition), and higher values of δ 34 S to sulfate remaining after BSR.

Evaporite Evolution
The near absence of halite in large sulfate accumulations in Picacho and Tucson Basins, was noted as a problem in Section 4.2. A related problem applies in the Hualapai and Luke Basins, where large halite accumulations with little sulfate occur [16,23]. Even if the

Episodes of Evaporite Formation
It has been assumed above that all sulfur in the larger evaporite bodies originated from bedrock weathered at the time of Basin-and-Range extension that began about 15 Ma ago in southern Arizona. Older gypsum exists in the region: in the Mineta Formation of Oligocene age, in the Magdalena and Tubutama Formations of Early Miocene (about 22 Ma) age [17] in Sonora, Mexico, in the mid-Miocene or older Pantano Formation [33] of Tucson Basin, Arizona and dispersed in lacustrine sediments of the syn-volcanic Late Cretaceous Cabullona Group [80] at Huépac, Sonora [81]. These have isotope characteristics like those of the Basin-and-Range evaporites (Figures 4 and 5a). This suggests mixing of igneous or PCD and marine sulfur in non-marine sediment of the Cordillera over a protracted interval. Sulfate accumulated in the deep troughs of the Basin and Range Province may partly represent re-working of earlier non-marine evaporites, the salts progressively moving towards the deepest depressions. The need for very large sources of sulfur in the Picacho and Safford Basins stands, but the sources may be distributed widely in space and time.

Evaporite Evolution
The near absence of halite in large sulfate accumulations in Picacho and Tucson Basins, was noted as a problem in Section 4.2. A related problem applies in the Hualapai and Luke Basins, where large halite accumulations with little sulfate occur [16,23]. Even if the sulfate evaporites mainly result from weathering of sulfur-rich, chloride-poor bedrock, some chloride is supplied in dust and rainwater, as observed in present-day halite enrichment in Tucson Basin soil profiles [64]. Cyclic reworking of evaporites could separate chloride from sulfate. An evaporite deposit about to be eroded may lose halite to groundwater before gypsum becomes mobile because of the greater solubility of halite; such a separation is occurring at present where Permian marine halite of the Holbrook Basin, northeastern Arizona, is dissolving into superjacent groundwater [94]. Permanent separation of chloride and sulfate remains difficult, the solubility difference notwithstanding. Chloride-bearing groundwater would need to flow into a downgradient basin before substantial transport of sulfate occurred. Dewatering of gypsum to form anhydrite would have provided abundant hot water to leach halite, and could have transported the chloride into the Luke Basin from the Picacho and Chandler Basins, 115 km and 50 km to the southeast, respectively.
The Death Valley lacustrine deposits in California provide an example of rift-basin evaporite that has undergone little modification through interaction with groundwater. The strata retain halite, soluble Na sulfates and a complex, fine-scale layering of δ 34 S and δ 18 O SO4. Original large evaporite deposits in southern Arizona were most likely of similar nature, prior to reworking by groundwater.

Study Significance and Future Research
This study highlights the complexity of sulfur geochemistry in a tectonically-active continental region in which interior basins have formed and eroded over 100 million years, and in which basins have at times been closed and at others hydraulically connected. Similar complex evolution is likely in other active orogenic zones. The large sizes and non-marine isotope compositions of several of the evaporite deposits indicate that their potential contribution to surface water and groundwater should be considered in addition to regional marine evaporites when δ 34 S is used in weathering budget calculations.
Several questions remain inadequately answered, e.g.,: (1) What are the bulk value of δ 34 S and the total sulfur content in PCDs of the study area? (2) What mechanisms separate sulfate and chloride as evaporites are deposited and reworked? (3) Why is there an apparent discrepancy between the size of sulfate evaporites and local sources of sulfur? (4) What is the sulfur content and sulfur isotope composition in certain major bedrock units in this region? Larger datasets are required in order to address such questions. This will be possible for numbers (1), (3) and (4). Question (2) requires detailed studies of basin development, which may prove difficult because collections of twentieth-century drill-core cuttings have been discarded since the initiation of this study.

Conclusions
1. Cenozoic gypsum and anhydrite evaporites in southwest North America have δ 34 S values ranging from 0 to +22‰, except in a limited area of Safford Basin, for which the range is −30 to −10‰. Most δ 34 S values fall within the range +4 to +10‰. The range of δ 18 O SO4 is +3 to +19‰.
2. The evaporites are non-marine in origin, on the basis of their isotope values and geological setting.
3. A linear relationship between δ 18 O SO4 and δ 34 S, evident in the Sonoita Creek Basin [89], is a principal control of sulfate isotope composition in most of the evaporites, and is consistent with mixing between sulfate from Permian or Cretaceous marine gypsum with sulfate from weathered bedrock sources such as Laramide igneous rocks and associated porphyry copper deposits.
4. The isotope composition of sulfate from which gypsum formed remained about constant in certain large deposits in Tucson and Picacho Basins but changed greatly over time in Safford Basin.
5. Evaporites with δ 34 S greater than +10‰ mostly derived sulfur from Permian marine gypsum, e.g., at White Sands and Estancia, New Mexico [60] and Safford and Marsh Station (this study). Re-worked Early Cretaceous marine sulfate may be present in the Safford and San Pedro Basins. A minor evaporite at Douglas, Arizona, is deriving sulfur from pyritic limestone, and may have evolved to high δ 18 O SO4 and δ 34 S values through cyclic recrystallization.
6. Bacterial sulfate reduction influenced sulfate isotopes strongly in few cases: Estancia, New Mexico [60] and possibly in limited stratigraphic intervals at Death Valley, California In Safford Basin, BSR led to deposition of pyrite in an anoxic lake, but of insufficient amount to account for a large volume of low-δ 34 S gypsum in the basin. BSR operating during formation of a supergene enrichment blanket at a nearby porphyry copper deposit may have provided sufficient low-δ 34 S sulfate.
7. The amount of sulfur in large evaporite bodies in southern Arizona is larger than can be explained by weathering of source rock at the time of deposition of the bodies. Long-term reworking of earlier non-marine evaporites may have occurred.