In the southern Gulf of California, the generation of new oceanic crust has resulted in linear magnetic anomalies and seafloor bathymetry that are characteristic of active seafloor-spreading systems. In the northern Gulf of California and the onshore (southeastern California, USA) Salton Trough region, a thick sedimentary package overlies the crystalline crust, masking its nature, and linear magnetic anomalies are absent. We use potential-field data and a geotherm analysis to constrain the composition of the crust (oceanic or continental) and develop a conceptual model for rifting. Gravity anomalies in the northern Gulf of California and Salton Trough are best fit with crustal densities that correspond to continental crust, and the fit is not as good if densities representative of mafic rocks, i.e., oceanic crust or mafic underplating, are assumed. Because extensive mafic underplated bodies would produce gravity anomalies that are not in agreement with observed gravity data, we propose, following earlier work, that the anomalies might be due to serpentinized peridotite bodies such as found at magma-poor rifted margins. The density and seismic velocities of such serpentinized peridotite bodies are in agreement with observed gravity and seismic velocities. Our conceptual model for the Salton Trough and northern Gulf of California shows that net crustal thinning here is limited because new crust is formed rapidly from sediment deposition. As a result, continental breakup may be delayed.
The Gulf of California and southern Salton Trough (northwestern Mexico and southeastern California, USA) form the oblique Gulf extensional province (Gastil et al., 1979) of generally right-lateral transform and strike-slip faults, short seafloor-spreading segments, and pull-apart and extensional basins, extending northwest from the East Pacific Rise spreading ridge to the San Andreas fault in the continental U.S. (Figs. 1, 2; e.g., Fuis and Kohler, 1984; Lonsdale, 1989, 1991; Stock and Hodges, 1989; Lizarralde et al., 2007; Umhoefer, 2011; Bennett and Oskin, 2014; Abera et al., 2016; van Wijk et al., 2017; Umhoefer et al., 2018). Extension in the northern Gulf extensional province may have started between ca. 19–17 Ma and 12.2 Ma (Bennett et al., 2016), progressively localizing into dextral transtension by ca. 8 Ma (Atwater and Stock, 1998, 2013; Bennett, 2009; Seiler, 2009; Seiler et al., 2010; Bennett et al., 2013, 2016; Bennett and Oskin, 2014; Darin et al., 2016). There is general agreement that extension in the gulf commenced sometime after the cessation of subduction of the Guadalupe and Magdalena microplates (fragments of the Farallon plate; e.g., Lonsdale, 1989, 1991; Stock and Lee, 1994; Umhoefer, 2011; Bennett and Oskin, 2014). Magnetic anomalies indicate that seafloor spreading started in the southern Gulf of California by ca. 3.5 Ma during the Pliocene (Larson, 1972; Klitgord et al., 1974; Lonsdale, 1989; DeMets, 1995), but north of ∼28° latitude, these anomalies are absent and the nature of the crust (continental, “transitional” [sometimes defined as thinned, with igneous intrusions, and underplating], or oceanic) is debated (Fuis and Kohler, 1984; Lachenbruch et al., 1985; Axen, 2008; Schmitt et al., 2013; Persaud et al., 2016; Han et al., 2016). The absence of seafloor-spreading magnetic anomalies in the northern Gulf may be explained either by the nature of the crust (if it is continental, no seafloor-spreading magnetic anomalies are expected), or by temperatures that are close to or above the Curie temperature so that underlying igneous (oceanic) crust would not be magnetized. We explore both ideas in this study with potential-field modeling and a one-dimensional geotherm analysis.
Potential-field data may discriminate between oceanic and continental crust based on differences in composition and thus density and magnetic susceptibility. Potential-field models are, however, nonunique and are therefore constrained in this study by seismic reflection and refraction profiles, receiver functions, and tomographic images (Fuis and Kohler, 1984; Parsons and McCarthy, 1996; Lewis et al., 2001; González-Escobar et al., 2009, 2014; Barak et al., 2015; Persaud et al., 2016; Han et al., 2016). Burial history (geohistory) and thermal models were developed for the Salton Trough to obtain the geotherm, the depth to the Curie isotherm, and the antigorite breakdown depth. We have developed gravity models for four transects: one across the Wagner and Consag Basins (here referred to as the Wagner basin transect) and three in the onshore Salton Trough region (Figs. 1, 2). We selected these transects because of the availability of geophysical data sets to constrain the potential-field models.
The nature of the crust in the northern Gulf of California and Salton Trough has been a topic of research for decades. The crust is covered with a thick layer of sediments that mask its nature, and the nonuniqueness of seismic velocities and gravity data, which have been used to understand the crustal structure in the northern Gulf extensional province, makes inferences on its nature difficult. Fuis and Kohler (1984) developed a gravity anomaly model across the Salton Trough with layers of sediments overlying a metasedimentary layer that, in turn, overlies a layer with a density of 3100 kg/m3. They obtained a reasonable fit between model and data, with an error between ∼3 and 8 mGal in the Salton Trough. From this model, they inferred the lower crust to be oceanic (mafic igneous) in nature. In their model, layer thicknesses and Moho depth were constrained only in a few locations. Han et al. (2016) noted that an upper crust composed of either stretched continental crust or a high-grade metamorphosed sedimentary layer would be in agreement with seismic velocities, and proposed that a gabbroic lower crust may be present below the Salton Trough region. Persaud et al. (2016) found that the seismic velocity structure in the Salton Trough agrees with a layer of possibly metasedimentary crust (Fuis et al., 1984; Dorsey, 2010) but is likely low-velocity (<∼7 km/s) stretched crystalline basement overlying a mafic underplating. Robinson et al. (1976) described basaltic rocks in the Salton Trough geothermal field as mantle-sourced (hydrothermally altered) low-potassium tholeiitic basalts with composition similar to that of the East Pacific Rise. They proposed a leaky transform model for emplacement of these basalts and rhyolites. Schmitt et al. (2013) found Nd isotopic ratios from young volcanic rocks in the northern Gulf of California and Salton Trough to overlap with those of East Pacific Rise basalts, and interpreted this as oceanic magmatism.
The Salton Trough and northern Gulf of California are characterized by locally very high (>800 mW/m2) surface heat-flow amplitudes (Lachenbruch et al., 1985; Elders and Sass, 1988; Prol-Ledesma et al., 2013; Neumann et al., 2017). The geothermal potential associated with these areas was recognized in the 1980s and was investigated with the Salton Sea Drilling Project (Ross and Forsgren, 1992). Geothermal studies in the Salton Trough (Muffler and White, 1969; Elders and Sass, 1988) recorded a bottom-hole temperature of 280 °C in the Chevron Wilson No. 1 well, drilled to a depth of 4.1 km (Muffler and White, 1969; Fig. 2). Even higher temperatures were found in the State 2-14 well, drilled to a depth of 3.22 km with bottom temperature of 355 °C (Elders and Sass, 1988; Figs. 2 and 8). The origin of these geothermal manifestations is debated, and might be from recent igneous activity, groundwater flow, and/or extreme crustal thinning (Han et al., 2016; Neumann et al., 2017). The geotherm outside of these geothermal regions is the focus of our geothermal analysis.
POTENTIAL FIELD STUDY
Potential-field models were constructed using Geosoft’s Oasis montaj software suite. The gravity and magnetic anomaly data for this region are from Bonvalot et al. (2012) and Maus et al. (2009), respectively. The pre-processed Bouguer gravity data (see Pavlis et al.  for processing details) were downloaded from the 2008 Earth Gravitational Model (EGM2008) from the Bureau Gravimétrique International (Toulouse, France; http://bgi.obs-mip.fr/). The regional grid is a 2.5 × 2.5 arc-minute grid, from 124°W to 104°W longitude and from 28°N to 38°N latitude. The gravity data were gridded to two different Universal Transverse Mercator (UTM) zones using the World Geodetic System 1984 standard (WGS84); the Wagner Basin data were mapped to UTM zone 12 and the Salton Trough data were mapped to UTM zone 11. Once the gravity data were gridded to the appropriate zone with the Oasis montaj software suite, a map of each zone was generated using the minimum curvature method with a cell size of 100 × 100 m. Hereby the software’s default iteration criteria were adopted. The resulting Bouguer gravity anomaly map for UTM zone 11 is shown in Figure 3A.
The total-field magnetic data used in this study were obtained from the U.S. National Oceanic and Atmospheric Administration (NOAA) National Centers for Environmental Information. These data are compiled from both land and sea surveys and processed into a single data set, EMAG2. EMAG2 is a 2–arc-minute–resolution data set spanning the whole Earth with the datum altitude set to 4 km above mean sea level. For areas with little available data, a least-squares co-location method was used. Additional information on the processing of the EMAG2 data set is described in Maus et al. (2009).
The EMAG2 data set was gridded to both UTM zones 11 and 12 using the Oasis montaj software suite for analyzing the Salton Trough region and Wagner Basin region, respectively. Once the magnetic data were loaded into the UTM grids, the grids were transferred into maps using the method of minimum of curvature with default iteration criteria in Oasis montaj. After the magnetic maps for UTM zones 11 and 12 were created, they were reduced to pole. Because each study area (Wagner Basin and Salton Trough) is relatively small, the magnetic data maps for UTM zones 11 and 12 were each reduced to pole at a single point. We chose for the Wagner Basin reduction 114.08°W longitude and 30.96°N latitude; this point is in the center of the Wagner Basin transect and results in an inclination angle of 56.90° and a declination angle of 11.37°. The location chosen for the Salton Trough region reduction to pole is 115.56°W longitude and 33.11°N latitude, which is a weighted average of the Salton Trough transects. This resulted in an inclination angle of 58.64° and a declination angle of 12.14°. Inclination angles and declination angles were determined using the International Geomagnetic Reference Field (IGRF) from 2008. May 2008 was used as the reference date because the total-field magnetic data were gathered around this time. Figure 3B shows the resulting magnetics map for UTM zone 12.
Potential-Field Model Setup and Constraints
Subsurface models of density structure were developed to calculate Bouguer gravity anomalies. These calculated anomalies were then compared to the measured values to determine how well the model fit the data. Misfits are expressed as errors in milligals, calculated by subtracting the predicted anomalies from the observed data set and averaging this difference over the transect. The Oasis montaj software suite uses the trends in calculated gravity anomaly values rather than the absolute values, which means that the software automatically creates a best-fit scenario by subtracting (or adding) a linear trend to the entire calculated anomaly values to match the measured field data. We did not attempt to fit calculated magnetic anomalies to measured values, as variations in magnetic anomaly magnitudes are small along the transects and locations of potential sources (such as intrusives) are not well constrained; instead, we provide a qualitative discussion. Table 1 lists the densities used for the gravity models.
The same densities were used for the Wagner Basin and Salton Trough transects. Whitmeyer and Karlstrom (2007) showed that east of the Salton Trough region, ca. 1.8 Ga Archean crust and ca. 1.7 Ga granitoids are present, while east of the northern Gulf of California, ca. 1.8 Ga Archean crust and ca. 1.65 Ga granitoids are found, among other formations. Because of a lack of more detailed evidence pointing at different overall basement densities across our transects, the same densities were used for all transects. Uppermost crust is composed of a layer of unconsolidated to poorly consolidated sediments (here referred to as “sediments”), and consists of sand and sandstone, silt and siltstone, and clay and claystone (e.g., Herzig and Elders, 1988; Dorsey et al., 2011; Table 2). The layer of metamorphosed sedimentary rocks consists mainly of greenschist-facies metasedimentary rocks (Herzig and Jacobs, 1994). Densities of the crustal rocks must actually vary with depth and temperature, and the layered structure of the continental and juvenile oceanic (or gabbroic) rocks in our potential-field models is therefore a simplification. We selected average densities for these layers from a range of densities used in studies by Nafe and Drake (1961), Schubert and Sandwell (1989), Carlson and Herrick (1990), Maystrenko and Scheck-Wenderoth (2009), van Avendonk et al. (2009), and Tenzer et al. (2011). Lower oceanic crust is also referred to as gabbroic material in this study. Crustal layers are named according to standard nomenclature (upper, middle, and lower continental or oceanic crust) for ease of discussion, but we note that these names do not necessarily reflect their origin, as discussed later. Upper mantle densities are quite low (3200 kg/m3), reflecting the possible presence of fluids and/or partial melt (Wang et al., 2009; Zhang et al., 2009; Barak et al., 2015; Persaud et al., 2015).
Three transects are modeled in the Salton Trough (Fig. 2). They extend from 116.68°W, 33.98°N to 114.56°W, 32.56°N (line 1), from 116.39°W, 33.02°N to 115.31°W, 33.52°N (line 2), and from 116.02°W, 32.77°N to 114.73°W, 33.02°N (line 3). These transects are close to the transects shown in Han et al. (2016): line 1 of this study corresponds with line 1 of Han et al. (2016), but has an endpoint shifted east; line 2 of this study corresponds with line 3 of Han et al. (2016), but the entire transect is shifted somewhat northwest; and line 3 of this study corresponds with line 2 of Han et al. (2016) and is in approximately the same location. Line 3 of this study also corresponds with the transect in Persaud et al. (2016), which is used to constrain some of the models. None of the lines cross the magnetic highs associated with the Peninsular Ranges batholith (Langenheim and Jachens, 2003).
Five models, labeled based on their layer properties (cc is continental crust, oc is oceanic crust, ms is metasediments, gb is gabbro, sp is serpentinized peridotite), are developed for line 1 in the Salton Trough (Figs. 3, 4). Models shown in Figures 4A, 4B, and 4C follow the interpretations by Han et al. (2016); models shown in Figures 4D and 4E follow the interpretations by Persaud et al. (2016). Line 1 model cc (Fig. 4A) assumes that there is continental crust below a layer of sediments in the Salton Trough. Line 1 model cc + ms (Fig. 4B) is the same, but now a layer of metasediments underlies the sedimentary layer. Comparison between these models demonstrates the effect of the metasedimentary layer. Line 1 model cc + ms + gb (Fig. 4C) assumes that both metamorphic sediments and gabbroic underplating are present in the Salton Trough region. This model follows the subsurface structure suggested by Han et al. (2016), whereas the models shown in Figures 4A and 4B follow only the general geometry of the structure suggested by Han et al. (2016). Model cc for line 1 (Fig. 4D) assumes the geometry from Persaud et al. (2016) without a layer of metasediments; model cc + ms (Fig. 4E) is the same, but includes a layer of metasediments.
Three models are developed for the line 2 transect (Figs. 3, 5). Line 2 model cc (Fig. 5A) assumes a stretched continental crust with no metasediments or gabbroic layer. Line 2 model cc + ms (Fig. 5B) assumes that the Salton Trough is underlain by stretched continental crust and a layer of metasediments. Line 2 model cc + ms + gb (Fig. 5C) is the same, except for a layer of gabbro that underlies the stretched crust. This model is inferred from Persaud et al. (2016).
Five models are developed for line 3 (Figs. 3, 6). Line 3 model cc (Fig. 6A) consists of a layer of sediments overlying stretched continental crust. Line 3 model cc + ms (Fig. 6B) is the same, but includes a layer of metasediments. Line 3 models cc + ms + sp and cc + ms + gb (Figs. 6C, 6D) consists of stretched continental crust beneath a metasedimentary layer and overlying a layer of serpentinized peridotite or a gabbroic layer, and a crustal structure inferred from Persaud et al. (2016). Line 3 model oc + ms (Fig. 6E) uses the general Persaud et al. (2016) geometry, but assumes that all crust along the axis of the Salton Trough below the metasedimentary layer is oceanic in nature. Indicated in Figure 6B are depth constraints used for this model from Ichinose et al. (1996), Zhu and Kanamori (2000), and Han et al. (2016) and the 5.65 km/s P-wave isovelocity contour from Persaud et al. (2016) which indicates the base of the sediment layer beneath the Imperial Valley. Extrapolation of Moho depths on the west side of lines 2 and 3 is aided by data from Lewis et al. (2001). All geophysical model constraints have uncertainties within which boundaries between our model layers were allowed to be moved up or down. The 1 km gridded velocity structure in Persaud et al. (2016) for example is smoothed over 5 × 5 × 2 km, and we allow for 2–5 km of precision in picking our top of basement. Han et al. (2016) used a smoothing ratio of 20:1 in the lower crust; this provides allowable precision in picking our Moho. Zhu and Kanamori (2000) published a Moho depth map of southern California using receiver functions. We have used this map to constrain the Moho depths outside of the rift zone. Uncertainties in Moho depths from Lewis et al. (2001) are estimated to be ±3 km. Wells drilled in the region provide some information on the metasediments (Palmer et al., 1975; Sass et al., 1988; Rose et al., 2013). The gabbroic layer in Barak et al. (2015) and Han et al. (2016) constrains our gabbroic layer.
The Wagner Basin transect (line 4 in Figs. 2 and 3) used for the potential-field models extends from 115.43°W, 31.64°N to 112.73°W, 30.27°N; this orientation is chosen so that it is nearly perpendicular to the trend of the basin axis. The transect crosses both the Consag fault located on the western side of the Wagner Basin and the Wagner fault on the eastern side of the basin. We developed three models for this transect. In line 4 model cc (Fig. 7A), a layer of sediments overlies stretched continental crust. In line 4 model cc + ms (Fig. 7B), a layer of sediments overlies a layer of metasediments, which in turn overlies a layer of stretched continental crust. Comparison between models cc and cc + ms shows the effect of the layer of high-density metasediments on Bouguer gravity anomalies. Model oc + ms (Fig. 7C) is different from the model in Figure 7B in that the basement is now of oceanic composition. Comparison of line 4 models cc + ms and oc + ms thus shows the effect of oceanic versus continental crust on Bouguer anomalies. All three models are constrained using the same geophysical data sets, then adjusted to match the potential-field data as closely as possible. As with the Salton Trough transects, adjustments are made to layer thicknesses in places that are not constrained by available data, or within uncertainties.
Basement depths in the Wagner Basin transect are constrained using González-Escobar et al. (2009) and González-Escobar et al. (2014). Bathymetry is constrained by Bischoff and Niemitz (1980). The locations of the Wagner and Consag faults are obtained from González-Escobar et al. (2009). Two-way travel times from seismic line 5023 from González-Escobar et al. (2009) are used to estimate basement depth. For this estimation, a seismic velocity of 4000 m/s is assumed for sediments, in agreement with densities found in the State 2-14 well sediments located in the Salton Trough (Sass et al., 1988) and the Nafe and Drake (1961) curve for density-velocity relationships. This seismic velocity is, however, somewhat higher than typical for northern Gulf of California sediments. In the northern Gulf of California, Phillips (1964), González-Fernández et al. (2005), and Persaud et al. (2003) found velocities <4 km/s in the upper ∼3 km of the sediment package. Our depth to basement is therefore considered a minimum depth, and may locally be ∼500–1000 m deeper than our initial estimate. We have adjusted this boundary to larger depths where needed to fit the gravity data.
Transect C-C′ from Lewis et al. (2001) is used for Moho depth constraints and crosses the Wagner Basin transect used in this study. In the Lewis et al. (2001) study, the closest station to our transect is the ROKO station (located at 114.45°W, 31.1°N), which is <0.2° south of the Wagner Basin transect. The Moho depth is ∼15 ± 3 km at the ROKO station. González-Fernández et al. (2005) found a ∼15-km-deep Moho in the nearby Delfin Basin.
Line 1. Model cc for line 1 in the Salton Trough (Fig. 4A) is a simplified model based on the crustal and mantle geometry of Han et al. (2016). This model should be considered a base model, as there is strong evidence of metasediments from wells located in the region (Elders and Sass, 1988; Sass et al., 1988; Rose et al., 2013). Adjustments (<∼2 km) to layer thicknesses did not result in a reasonable fit to Bouguer anomalies; thus we show the initial model in Figure 4A to illustrate where most of the misfit occurs. The average misfit between model-predicted and observed gravity anomalies is 37.87 mGal. The fit on the southeastern side of the transect is particularly poor; the crust thickens here in the model, but this geometry was poorly constrained according to Han et al. (2016). The calculated Bouguer gravity is too high between 60 and 205 km along the transect. This indicates that the Moho here is likely deeper than assumed in this model.
Model cc + ms (Fig. 4B) is similar to model cc, but includes a metasedimentary layer. This model also is not adjusted to match the Bouguer gravity data; the average misfit is 34.73 mGal. Including the metasedimentary layer reduces the misfit between 60 and 205 km along the transect relative to model cc. Model cc + ms + gb (Fig. 4C) follows Han et al. (2016) and includes an underplated gabbroic layer below stretched continental crust. This increases the average misfit along the transect to 39.50 mGal.
Model cc in Figure 4D is set up following the general velocity model in Persaud et al. (2016), through the intersection point with line 3 (line 3 corresponds to the transect in Persaud et al. ; see Fig. 2). It should be considered a base model as it does not include a metasedimentary layer. This model has been adjusted to fit the observed Bouguer gravity anomaly data within the uncertainties of the geophysical constraints; the fit is quite good (2.01 mGal misfit). Model cc + ms in Figure 4E includes a metasedimentary layer. The fit also is good (1.75 mGal misfit). We also developed the same model as shown in Figure 4E with a gabbroic underplated layer (not shown); this reduced the fit significantly. Of all models tested, model cc + ms in Figure 4E provides the best fit.
The magnetic anomalies along line 1 (Fig. 4A, top panel) could not be related to known subsurface structures. The horizontal resolution of the magnetic data set is ∼4 km, and magnetic anomalies typically vary ∼60 nT over a distance of ∼60 km. There is no relation between the magnetic anomalies and a possible gabbroic underplated layer in model cc + ms +gb (Fig. 4C).
Line 2. Model cc for line 2 in the Salton Trough (Fig. 5A) includes only stretched continental crust, and is the base model. Where it intersects line 1, the correlated velocity structure of Persaud et al. (2016) from line 3 is used to constrain layer depths. The model produces an average difference between the calculated anomalies and observed anomalies of 3.15 mGal, and errors are mostly from the ends of the transect, where it crosses mountainous terrain that we did not attempt to fit, as it is outside of our region of interest. Below the Salton Trough, the Moho is at ∼17.5 km depth in the model. Model cc + ms (Fig. 5B) includes a layer of metasediments. This improves the fit to the observed gravity data somewhat (average error is 3.09 mGal).
In model cc + ms + gb (Fig. 5C), we included a layer of gabbro. This worsens the fit to observed Bouguer anomaly data, and the average misfit is 5.69 mGal. Addition of the high-density gabbroic layer increases the calculated gravity trend over the entire transect. This demonstrates nicely the effect of subtracting regional trends from gravity data; the software package automatically generates a best match for the observed gravity by removing regional trends. The introduction of higher-density material in the center of the transect causes a decrease in the magnitude of the regional trend that is removed from the entire transect, resulting in higher predicted values across the entire transect. The error in this model is reduced by reducing the thickness and lateral extent of the gabbro layer.
The magnetic anomaly data for line 2 (Fig. 5A, upper panel) provide little additional information. Magnetic anomalies typically vary ∼40 nT over a distance of ∼35 km. We interpret the higher magnetic anomalies within the Salton Trough as likely associated with shallow igneous intrusions in the area.
Line 3. Model cc for line 3 in the Salton Trough (Fig. 6A) is based on the velocity structure of Persaud et al. (2016). This is the base model, consisting of stretched continental crust. The Moho follows the 7 km/s velocity contour from Persaud et al. (2016); although the ray coverage at those depths is poor (Persaud et al., 2016), this Moho depth provides a good fit to the observed gravity data. The shallowest Moho depth in the gravity model is ∼16 km. Model cc is a good fit to Bouguer anomalies, with an average misfit between model-predicted and observed gravity anomalies of only 0.73 mGal; this misfit is distributed evenly along the transect. Model cc + ms (Fig. 6B) consists of stretched continental crust with a metasedimentary layer. Including this layer reduces the misfit slightly to 0.7 mGal. Constraints for the Moho and top of basement are from Ichinose et al. (1996), Zhu and Kanamori (2000), Han et al. (2016), and Persaud et al. (2016). Figure 6C shows model cc + ms + sp with serpentinized peridotite instead of gabbroic material (see Fig. 6D) below the stretched crust. The density of the serpentinized peridotite was assumed to be the same as the density of continental crust in this model, so the predicted Bouguer anomaly does not change from that shown in Figure 6B.
Model cc + ms + gb (Fig. 6D) includes a gabbro layer below stretched continental crust. The model is based on the velocity structure of Persaud et al. (2016) with the assumption that their mafic crust is a gabbro layer, similar to the interpretation of Han et al. (2016). Introducing the gabbroic layer increases the misfit to 2.5 mGal between the calculated and observed gravity anomaly. At several locations along the transect, we could not reduce the misfit without moving the Moho depth to significantly greater depths (i.e., outside the range of geophysical uncertainties) than estimated from other studies. The lateral extent of the gabbro layer is somewhat smaller than in Persaud et al. (2016); the misfit increased with a wider gabbroic layer.
Model oc + ms (Fig. 6E) assumes that oceanic crust is present below the metasediments. Constraints for the metasediments and sediments are from Persaud et al. (2016); the subsurface structure was then adjusted within geophysical constraints (Fig. 6C) to match the gravity data as closely as possible. The model has an average error of 2.1 mGal between the observed and calculated gravity anomaly. The model fits the observed gravity data well between 63 km and the end of the transect, but we were not able to obtain a better fit between 0 and 63 km along the transect.
In summary, models in which the lower crust consists of low-density (2800 kg/m3) material provide the best fit with gravity data. The magnetic anomaly data for line 3 (Fig. 6A, upper panel) are characterized by long-wavelength variations that are not clearly related to the subsurface structure (gabbroic body and/or oceanic crust).
Three potential-field models are developed for line 4 across the Wagner Basin (Fig. 7): one with only stretched continental crust (model cc), another with a layer of metasediments overlying stretched continental crust (model cc + ms), and a third with oceanic crust (model oc + ms). Model cc (Fig. 7A) fits the observed Bouguer gravity anomaly measurements well, with an error between the modeled and observed anomalies of 2.58 mGal. A noticeable misfit between data and the model exists between 115 km and 140 km along the transect, where the calculated gravity is lower than the observed values. This is at the location of the Wagner Basin, indicating that the density of sediments in the deep Wagner Basin must be higher than assumed in the model; we attribute this to the absence of a metasedimentary layer in the model, or the presence of sills (Persaud et al., 2003).
Model cc + ms is the best-fit model for the Wagner Basin transect (Fig. 7B). The Bouguer gravity anomaly data match the model predicted anomalies well throughout the entire transect, including the Wagner Basin. In this model, the crust is very thin: a Moho depth of ∼14 km underneath the Wagner Basin and metasediments to a depth of 8.7 km result in predicted gravity anomalies that are in agreement with data. This is in reasonable agreement with González-Escobar et al. (2009), who reported the top of acoustic basement at ∼7 km depth. Lewis et al. (2001) found the Moho to be ∼15 km deep just south of the Consag fault; González-Fernández et al. (2005) reported a ∼15-km-deep Moho in the nearby Delfin Basin, in general agreement with our model prediction. The misfit between the modeled and observed gravity is 2.32 mGal.
Wagner Basin model oc + ms (Fig. 7C) results in an error between the observed and calculated gravity anomalies of 3.31 mGal. This error is located mainly between 121 km and 153 km along the transect, where the basin is underlain by mafic crust. This is a clear indication that mafic crust below the Wagner Basin would result in somewhat higher gravity anomalies than observed, despite our choice of conservative density values for the mafic crust. Moho depth beneath the basin in this model is ∼14 km, consistent with the findings of Lewis et al. (2001). The mafic crustal layer in this model is ∼5.8 km thick.
Magnetic anomaly data along the Wagner Basin transect are shown in the upper panel in Figure 7A. Amplitudes gradually increase from ∼−140 nT to ∼120 nT from northwest to southeast (Fig. 3, 7A). High amplitudes in the southeast may be caused by the Comondù Group andesites. These data give no additional insight into whether the Wagner Basin is underlain by mafic crust.
Abundant evidence supports the presence of small-scale igneous intrusions in the study area (Persaud et al., 2003; Schmitt and Vazquez, 2006; Schmitt and Hulen, 2008; Barak et al., 2015; Karakas et al., 2017). Because the locations and dimensions of sills and dikes are not precisely known along our Wagner Basin transect, we have not attempted to build a model with these small-scale intrusions. We did build general models with shallow small-scale intrusions to test if they would fit the observed gravity anomalies. We specifically modified the Wagner Basin models (not shown here) to include small intrusions within the sedimentary and metasedimentary layers, and were able to obtain a good fit with observations after slightly adjusting layer thicknesses.
In summary, Wagner Basin models with a lower-density (2800 kg/m3) lower crust fit the gravity anomalies best. We note that with the nonuniqueness inherent to gravity data, the limited geophysical constraints, and the lack of constraints on exact layer densities, it is not possible to confidently constrain layer depths any further.
The magnetic field in the study area is characterized by a lack of linear anomalies (Fig. 3, top panels in Figs. 4–7) that could be correlated with juvenile oceanic crust or gabbroic bodies. This might be explained by these structures being absent, or, alternatively, gabbroic or other material with high magnetic susceptibility being present but located at depths close to or below the Curie isotherm. There are few constraints on the geotherm, and thus the depth of the Curie isotherm in the Salton Trough and northern Gulf of California, and these are mainly from geothermal areas.
Wellbore bottom temperatures in geothermal areas show very high temperatures at shallow depths (for example, 280 °C at 4.1 km depth in the Chevron Wilson No.1 well [Figs. 2, 8; Muffler and White, 1969] and 355 °C at a depth of 3.22 km in the State 2-14 well [Figs. 2, 8; Elders and Sass, 1988]). Heat flow studies in the Wagner Basin have found anomalies >800 mW/m2 in geothermal areas (Neumann et al., 2017). Petrological constraints and thermal modeling suggest a partially molten lower crust below the Salton Trough (Karakas et al., 2017).
The geotherm in a sedimentary basin is affected by many processes, including the rate of sediment deposition, crustal and mantle-lithosphere thinning, processes such as groundwater circulation, and magmatic intrusions. Rüpke et al. (2013) found that the sediment cover in a basin with low to normal sedimentation rates (<∼0.4 mm/yr) acts as an insulating blanket, maintaining high temperatures in the stretched basement. To understand how rapid sedimentation and crustal thinning affect the geotherm, we calculated a one-dimensional thermal evolution of the Fish Creek–Vallecito Basin, located west of the floor of the Salton Trough (Dorsey et al., 2011). We used the stratigraphic column from Dorsey et al. (2011). It consists of a continuous section of late Miocene–Pleistocene rocks in the Fish Creek–Vallecito Basin, deposited in the hanging walls of the Vallecito fault (early part of section; Winker and Kidwell, 1996; Kerr, 1984; Dorsey et al., 2007; Shirvell et al., 2009) and the West Salton detachment fault (Axen and Fletcher, 1998; location in Fig. 2). We consider this an analog to the stratigraphic buildup of the Salton Trough from ca. 8 to ca. 1 Ma. The measured section shows ∼5.5 km of sediment deposition during these ∼7 m.y., with average sedimentation rates of ∼0.78 mm/yr, and exceeding 2 mm/yr during rapid subsidence. The phase of subsidence is followed by uplift caused by initiation of dextral faulting in the region during the past ∼1 m.y. (Lutz et al., 2006; Dorsey et al., 2011, 2012). This uplift event did not occur in the present Salton Trough, and we have ignored it in the thermal history analysis. In these thermal models, groundwater circulation and magmatic intrusions are not included.
The Schlumberger PetroMod 1D software suite was used to construct the burial history (also called geohistory) and thermal models (Hantschel and Kauerauf, 2009). The stratigraphic column was backstripped, and as sediments were deposited in a forward model, temperature was calculated in the stratigraphic column and underlying lithosphere. Input parameters for the Petromod 1D model are listed in Table 2. Lithology and ages are from Dorsey et al. (2011), and we refer to that study for a discussion on possible sources of error. The PetroMod lithology library was used (Table 2; Hantschel and Kauerauf, 2009); where geologic formations consist of several lithologies, the formation was modeled as multiple layers with each representing the relative abundance and thickness of the different lithologies. Ages of the tops of these layers were estimated assuming sedimentation rates from Dorsey et al. (2011). Constant temperatures at the surface (20 °C) and lithosphere-asthenosphere boundary (1300 °C) and a time-varying heat flow at the base of the column were used to calculate temperatures in the lithosphere. We developed four different temperature models, in which the basal heat flow into the basin is varied. This basal heat flow was calculated according to the amount of crustal and mantle lithosphere thinning that occurs during rifting (the stretching factors) and thus varies over time as thinning continues. In model A (Figs. 8A, 8B), both crustal and mantle have stretching factors of two. Model B (Fig. 8B) has a stretching factor of two for the crust and three for the mantle lithosphere, and thus produces a higher basal heat flow, and model C (Fig. 8B) is the warmest of our models, with stretching factors of two for the crust and four for the mantle lithosphere. Model D is the coldest model, with a stretching factor of 1.5 for the crust, and two for the mantle. We assume that the crustal thickness in the Salton Trough was 30 km before the onset of extension (based on crustal thickness estimates at the flanks of the Gulf of California; see potential field models). The thickness of the mantle lithosphere before the onset of extension was assumed to be 70 km. During rifting, the mantle lithosphere layer thins to 35–17.5 km in the models, dependent on the mantle stretching factor. With these assumptions, the basal heat flow into the basin increases from ∼40 mW/m2 at 7 Ma to 63–72 mW/m2 today, with higher stretching factors corresponding to higher basal heat-flow values. Vitrinite reflectance (a measure of thermal maturity) values are calculated for each model, and compared with reported values from well Chevron Wilson No. 1.
Sedimentation rate in the Fish Creek–Vallecito Basin was very high from ca. 4.6 Ma to 3.1 Ma (Dorsey et al., 2011; Fig. 8A), with total subsidence occurring at a rate of ∼2.1 mm/yr. During this 1.5 m.y. time span, the rate of sedimentation in the basin was so fast that it caused cooling of the basin. This is observed in Figure 8A by the deepening of isotherms starting ca. 4.6 Ma. In the deeper part of the stratigraphic column, isotherms move to shallower depths starting ca. 3.1 Ma; this indicates slow conductive heating of the deeper basin starting around this time. In shallow parts of the stratigraphic column, isotherms move to shallower depths somewhat earlier, ca. 3.3 Ma. An important result of these thermal analyses is that very rapid sedimentation cools the uppermost crust. Earlier work by Rüpke et al. (2013) did not report this relation between sedimentation rate and geotherm of the upper crust, probably because sedimentation rates in the Salton Trough are much higher than tested in that study. When sedimentation rates are low to average, the sediment cover acts as a thermal blanket.
Present-day model temperatures at the bottom of the stratigraphic column strongly depend on the stretching factors; model C has the highest temperature (∼200 °C; Fig. 8) and model D has the lowest temperature (∼160 °C) at the base of the stratigraphic column. Even temperatures in our warmest model (model C) are much lower than observed in geothermal wells (Muffler and White, 1969; Elders and Sass, 1988): 280 °C at 3.22 km depth in the Chevron Wilson No. 1 well (Muffler and White, 1969) and 355 °C at a depth of 3.22 km in the State 2-14 well (Elders and Sass, 1988).
Our models also predict vitrinite reflectance (Hantschel and Kauerauf, 2009). The Chevron Wilson No. 1 well (see Fig. 2 for location) was found to have a mean vitrinite reflectance of 3.0 at a depth of ∼4.1 km (13,433 ft) (Barker, 1983). Model A reaches a mean vitrinite reflectance of only 0.65 at 4 km depth. Model C, the model with the highest basal heat flow, calculates a mean vitrinite reflectance of ∼0.85 at 4 km depth.
The rate and amount of lithospheric thinning influence the subbasin geotherm; in model C, this geotherm is 33 °C/km, in model B it is 29 °C/km, and in model A it is 22.5 °C/km (Fig. 8). Assuming a Curie temperature of 540 °C for basalt (Zablocki and Tilling, 1976) or a Curie temperature of 575 °C for gabbro (Kent et al., 1978), model A corresponds to a Curie temperature depth range of ∼20.0–20.8 km, model B corresponds to a Curie temperature depth range of ∼17.3–18 km, and model C corresponds to a Curie temperature depth range of ∼15–15.6 km. Localized intrusive igneous activity in the Salton Trough (Schmitt and Hulen, 2008) would have increased the geotherm through addition of heat in the sedimentary layers, or through promoting fluid circulation in the subsurface (Neumann et al., 2017). These effects were ignored in our models, and our Curie isotherm depths may locally be maxima.
The depth of the lithosphere-asthenosphere boundary is ∼40–50 km beneath the northern Gulf of California (Fernández and Pérez-Campos, 2017) and 40 km below the Salton Trough (Lekic et al., 2011), but varies somewhat. The boundary in our model B approaches this value best, and indicates a present-day lithosphere-asthenosphere boundary at 38.8 km. This suggests that the depth of the Curie isotherm below the Salton Trough is ≤17.3 km, and that any gabbroic lower crust (e.g., in our potential-field models) is close to the Curie temperature.
LITHOSPHERE EVOLUTION IN THE NORTHERN GULF OF CALIFORNIA AND SALTON TROUGH
For all transects, models in which continental crust densities are used produce gravity anomalies that are a better fit to gravity measurements than models in which bodies of higher-density (corresponding to mafic [oceanic] crust) are used. In models that include high-density gabbroic bodies or juvenile oceanic crust, the fit with observations is better when the igneous body is of smaller size. In the Salton Trough, we were not able to adjust line 1 models (northwest-southeast oriented) that include igneous crust within the limits of other geophysical constraints so that the misfit between modeled and observed gravity data became acceptable. We therefore conclude that such an extensive body of high-density material likely does not exist in the Salton Trough transect. Along the shorter lines 2 and 3 (northeast-southwest oriented), high-density material corresponding to oceanic crust is permitted, in the sense that the misfit between data and models is not large. We note, however, that the fit with the gravity data along these transects improves when the size of the high-density bodies is decreased. From these results, we conclude that densities in the Salton Trough at depths between ∼10 and 20 km are likely closer to 2800 kg/m3 than to 3000 kg/m3, and higher densities do not agree with the gravity data. In the Wagner Basin, high-density (oceanic) crust is permitted; the fit between models and data is however slightly better when the crust has lower (continental) densities.
In summary, the potential-field models favor relatively low densities (∼2800 kg/m3) between ∼10 and 20 km depth in the Salton Trough. We labeled these densities “continental crust” in our models (Figs. 4–7). If these layers consisted of silica-rich continental crust material, their seismic velocities are expected to be lower than observed velocities; seismic studies (e.g., Fuis et al., 1984; Parsons and McCarthy, 1996; Han et al., 2016; Barak et al., 2015; Persaud et al., 2016) consequently find high P-wave velocities (∼6.9 km/s: Parsons and McCarthy, 1996; ∼6.8 km/s: Han et al., 2016) and low S-wave velocities (<4 km/s: Barak et al., 2015) at ∼10–20 km depth beneath the Salton Trough. This velocity structure has been interpreted to be in agreement with a gabbroic lithology, as a result of decompression melting below the rifted zone (Han et al., 2016; Barak et al., 2015). However, gabbroic rocks have higher densities than continental crust material (3100 versus 2800 kg/m3).
This discrepancy between potential-field models and seismic studies can be explained if the high–P-wave, low–S-wave velocity bodies between ∼15 and 20 km depth consist of serpentinized peridotite (Nicolas, 1985; Mével, 2003), and are thus mantle rocks. Serpentinized mantle is found at (magma-poor) rifted margins (for example, Iberia: Chian et al., 1999; Newfoundland: Tucholke et al., 2004), and may have started to form there before it was exhumed (Skelton and Valley, 2000). In this interpretation, the Moho (Fig. 9) is located above the serpentinized peridotite body. Serpentine lowers the density of peridotite, and at 500 °C, ∼35% serpentine in peridotite corresponds to a density of 2800 kg/m3 (Mével, 2003). This percentage of serpentine lowers P-wave velocities to ∼6.9 km/s and S-wave velocities to ∼3.7 km/s (Mével, 2003), similar to observed values (Parsons and McCarthy, 1996; Barak et al., 2015). The presence of serpentinized peridotite below the Salton Trough was considered previously (Nicolas, 1985), but deemed unlikely because of inferred high temperatures at 10–20 km depth (Barak et al., 2015). Our thermal models show that the very rapid subsidence and sediment infill in the Salton Trough suppress the crustal geotherm (Fig. 8B), so that serpentinized peridotite is stable (Ulmer and Trommsdorff, 1995; Wunder and Schreyer, 1997) at the temperature range at 10–20 km depth below the Salton Trough (Fig. 8). This layer of serpentinized peridotite would be located ∼3 km below the brittle-ductile transition (Fig. 9), which implies that faults and shear zones act as fluid conduits to transport water to these depths.
Serpentinization of peridotite increases the magnetic susceptibility of the rock as a result of the formation of magnetite. When they acquire magnetic remanence, such serpentinized bodies could produce strong magnetic anomalies (Dyment et al., 1997; Maffione et al., 2014). We do not find evidence for linear magnetic anomalies along our transects, while temperatures below the Salton Trough at 10–20 km depth are somewhat below the Curie temperature of magnetite (Fig. 8). The lack of linear magnetic anomalies in the northern Gulf of California and Salton Trough region is thus enigmatic. When temperatures are very close to the Curie temperature locally, close to igneous intrusions, magnetization may not have occurred.
The velocity structure of the upper mantle below the Gulf of California and Salton Trough is characterized by a series of low-velocity regions, extending from ∼40 km to 150 km depth (Fig. 1; Wang et al., 2009; Di Luccio et al., 2014; Barak et al., 2015). Wang et al. (2009) noted that they are spaced ∼250 km apart, but Barak et al. (2015) noted that this spacing is disrupted in the Cerro Prieto region, and Di Luccio et al. (2014) noted slightly different locations for the low-velocity anomalies and that their spacing is 300–400 km in the southern gulf. The low-velocity zones have been explained as zones of asthenosphere upwelling (Wang et al., 2009; Di Luccio et al., 2014; Barak et al., 2015). None of these low-velocity zones are centered below rift or ridge axes (Fig. 1); e.g., Wang et al. (2009) showed that they underlie the Ballenas transform fault, the region south and west of the Wagner Basin axis, and northeast of the Carmen spreading center (Fig. 1). The low-velocity zone beneath the Salton Trough (not shown in Fig. 1) is west of the rift axis (Barak et al., 2015).
Because the low-velocity regions are all off axis and do not spatially coincide with potential seafloor-spreading segments, we propose that, instead of them being Gulf of California rift-related upwelling zones, they are relics of the spreading segments that separated the Pacific plate from the microplates that are now partially subducted beneath the North American plate and Baja California (Fig. 1). Seafloor spreading formed asthenospheric upwelling systems below the spreading ridges that separated the Farallon plate (and later, the microplates) from the Pacific plate. Below the spreading ridges, decompression melting occurred. After subduction of the microplates ceased, the spreading centers became inactive, but depleted asthenosphere remained. Cessation of Guadalupe-Pacific spreading is estimated at ca. 11 Ma (Michaud et al., 2006). Since that time, the North American plate has moved ∼550 km in a west-southwestward direction with respect to the sublithospheric mantle, assuming a relative velocity of 5 cm/yr between the North American plate and underlying mantle (Silver and Holt, 2002; Fig. 1), and the Gulf of California is now located above the fossil upwelling systems of depleted mantle that are imaged as low-velocity zones. When extension-related melting occurs in the gulf, the underlying fossil upwelling systems melt, and the depleted mantle produces mid-ocean ridge basalt (MORB)–like melts (Robinson et al., 1976; Schmitt et al., 2013) similar in chemistry to those of the East Pacific Rise further south (Schmitt et al., 2013). Off-axis melting anomalies and volcanism have been found elsewhere; in the Gulf of Aden (Leroy et al., 2010), they have been explained by channeling of plume material away from the Afar upwelling. However, the upper mantle low-velocity regions in the Gulf of California seem to not be correlated with reported off-axis volcanism.
In our preferred conceptual model for the Salton Trough (Fig. 9, constructed after line 3), the sedimentary and metasedimentary layers are underlain by a thin layer of stretched continental crust. The (diffuse) Moho is at ∼15 km depth. At sub-Moho depths where P-wave velocities are high, S-wave velocities are low, and densities are low, the uppermost mantle is serpentinized, with a decreasing percentage of serpentine with increasing depth, so that densities and seismic velocities gradually increase. Magmatic processes transport melts derived from depleted, now-extinct spreading-ridge systems upward. Our thermal model calculates that the brittle-ductile transition is at ∼11 km depth below the rift, in agreement with seismicity (Lin et al., 2007; Hauksson et al., 2012; Lin, 2013). In the conceptual model, the rapid deposition of sediments creates new crust as crustal material is removed by thinning. The rift crust thus undergoes limited net thinning during extension, possibly delaying continental rupture. Intense crustal deformation forms a locally diffuse Moho, as described at magma-poor rifted margins elsewhere (e.g., Enachescu, 1988).
More constraints on the temperatures and seismic velocities beneath the Wagner Basin are needed to determine whether more stretched continental crust, gabbroic juvenile oceanic crust, or serpentinized peridotite are present below the sedimentary, metasedimentary, and/or stretched continental crust layers there. The potential-field models allow for small, high-density gabbroic bodies to be present in the Wagner Basin, but continental crust densities provide better fits with observed gravity anomalies.
Gravity models were developed to obtain a better understanding of the subsurface density structure in the Salton Trough and Wagner Basin, which are areas of rapid stretching of the lithosphere in an oblique extensional setting. The models show that densities beneath ∼10–15 km depth need to be low (∼2800 kg/m3) in the Salton Trough to match gravity data. Because P-wave velocities are high at these depths (∼6.9 km/s or so), this layer likely does not consist of normal continental crust material, while its low density excludes the presence of extensive gabbroic bodies. We propose (following Nicolas ) that a layer of serpentinized peridotite could be present beneath the thinly stretched crust. A simple one-dimensional geothermal analysis shows that temperatures at these depths are within the antigorite stability field. Our gravity models do not constrain the nature of the crust in the Wagner Basin.
If no juvenile oceanic crust is present, igneous rocks with East Pacific Rise chemical signature in the study area may be explained by melting of relict, depleted upwelling systems. These are the fossil upwelling systems that were present below the now-extinct Farallon plate–Pacific plate spreading ridges. They have been overridden by the North American plate, and are now located below the Gulf of California. They may have been imaged seismically as low-velocity zones.
Rapid sedimentation creates new continental crust as continental crust is removed by extension and thinning in the Salton Trough. Net crustal thinning is thus limited, possibly delaying continental rupture.
We would like to thank Arturo Martin-Barajas for sharing a seismic profile that crosses the Wagner Basin, and Sascha Brune and an anonymous reviewer for constructive comments.