Abstract
Fractures in Earth's critical zone influence groundwater flow and storage and promote chemical weathering. Fractured materials are difficult to characterize on large spatial scales because they contain fractures that span a range of sizes, have complex spatial distributions, and are often inaccessible. Therefore, geophysical characterizations of the critical zone depend on the scale of measurements and on the response of the medium to impulses at that scale. Using P-wave velocities collected at two scales, we show that seismic velocities in the fractured bedrock layer of the critical zone are scale-dependent. The smaller-scale velocities, derived from sonic logs with a dominant wavelength of ~0.3 m, show substantial vertical and lateral heterogeneity in the fractured rock, with sonic velocities varying by 2,000 m/s over short lateral distances (~20 m), indicating strong spatial variations in fracture density. In contrast, the larger-scale velocities, derived from seismic refraction surveys with a dominant wavelength of ~50 m, are notably slower than the sonic velocities (a difference of ~3,000 m/s) and lack lateral heterogeneity. We show that this discrepancy is a consequence of contrasting measurement scales between the two methods; in other words, the contrast is not an artifact but rather information—the signature of a fractured medium (weathered/fractured bedrock) when probed at vastly different scales. We explore the sample volumes of each measurement and show that surface refraction velocities provide reliable estimates of critical zone thickness but are relatively insensitive to lateral changes in fracture density at scales of a few tens of meters. At depth, converging refraction and sonic velocities likely indicate the top of unweathered bedrock, indicative of material with similar fracture density across scales.
Introduction
Interactions among Earth's atmosphere, biosphere, and hydrosphere drive chemical and physical processes that, over geologic time, convert bedrock into soil. These processes create and maintain a structure that extends from Earth's surface vegetation to subsurface depths where bedrock is no longer physically or chemically altered by near-surface processes. This thin region, which supports virtually all terrestrial life, is referred to as Earth's critical zone (CZ) (Richter et al., 2009; Brantley et al., 2017; Riebe et al., 2017). In eroding landscapes, the CZ is conceptualized as layers of soil, saprolite, fractured bedrock, and less altered bedrock (Figure 1). The spatial distribution and characteristics of each of these layers provide valuable information about the processes involved in the weathering of rocks (Anderson et al., 2007; St. Clair et al., 2015; Rempe and Dietrich, 2018), the generation of porosity (Graham et al., 1997; Navarre-Sitchler et al., 2015; Hayes et al., 2019), and water movement in landscapes (Salve et al., 2012; Brooks et al., 2015; Lovill et al., 2018).
Figure 1
Many hypotheses that explain CZ evolution have grown out of observing and understanding the existing architecture of the CZ (Figure 1) (e.g., Riebe et al., 2017). Fundamental questions about CZ architecture remain incompletely answered, including: (1) How thick is the CZ? (2) Are there distinct layers and boundaries within the CZ? (3) If so, how deep are they, how well-defined are they, and how do they vary spatially across a landscape? The answers to these questions will provide crucial tests of hypotheses about the processes that form and maintain the CZ. For example, the interaction between topographic and compressive regional tectonic stresses may control the opening of subsurface fractures, leading to inverted bedrock topography under some hillslopes (Slim et al., 2015; St. Clair et al., 2015; Moon et al., 2017). CZ structure may depend strongly on strong contrasts in hydrogeologic properties such as porosity or hydraulic conductivity (Dralle et al., 2018; Rempe and Dietrich, 2018; Hahm et al., 2019). In some locations, lithologic fabric or foliation (Novitsky et al., 2018; Leone et al., 2020), composition (Hahm et al., 2014; Buss et al., 2017) or geochemical reaction rates (Dixon and Thorn, 2005; Gabet and Mudd, 2009; Hilley et al., 2010; Braun et al., 2016) may control CZ structure and evolution. Testing any of these process-based hypotheses about CZ evolution, or developing new hypotheses, requires quantifying CZ architecture over large areas.
Quantifying CZ architecture requires both drilling and geophysics and must balance inevitable tradeoffs between detail and scale. Soil pits provide detailed, high-resolution information about the top 1–3 m of the CZ (Richter and Markewitz, 1995; Fletcher et al., 2006; Graham et al., 2010), including root distributions (Canadell et al., 1996; Pawlik et al., 2016; Hasenmueller et al., 2017) and water infiltration (Vereecken et al., 2007; Bales et al., 2011; Carey et al., 2019; Kotikian et al., 2019). Road cuts and outcrops provide a detailed look at the CZ (Dethier and Lazarus, 2006; Dewandel et al., 2006; Goodfellow et al., 2014) over larger scales than soil pits, but their locations are not always convenient for testing hypotheses. In contrast, boreholes, although spatially limited, provide in situ windows into deep regions of CZ at high vertical resolution. Boreholes create unique and detailed measurements of the physical structure and chemical composition of the CZ, including fracture density and orientations (Anderson et al., 2002; Kuntz et al., 2011; Holbrook et al., 2019; Gu et al., 2020b), geochemical reactions (Brantley and White, 2009; Brantley et al., 2013; Holbrook et al., 2019), and hydrogeologic properties like porosity and hydraulic conductivity (Walsh et al., 2013; Flinchum et al., 2018a; Ren et al., 2019). A major challenge in CZ science is to obtain measurements that have both high spatial resolution and cover large areas.
Geophysical methods are a tool for exploring deep CZ structure at landscape (1,000s of m) and watershed (10s−100s of m) scales. Near-surface geophysical methods have variable resolution and depth penetration and have been successfully used to quantify lateral and vertical variability in CZ structure (Befus et al., 2011; Leopold et al., 2013; Orlando et al., 2016), estimate bedrock fracture density (Clarke and Burbank, 2011; Golebiowski et al., 2016; Novitsky et al., 2018), and estimate the depth to bedrock (Nyquist et al., 1996; Yaede et al., 2015; Flinchum et al., 2018b; Moravec et al., 2020). Given the diversity of near-surface geophysical methods (e.g., Parsekian et al., 2015), it is critical to select the appropriate geophysical tool to test a given scientific hypothesis. It is also critical to understand how the geophysical properties are related to the underlying physical properties that define CZ structure (Figure 1).
Seismic refraction is a geophysical method that is being increasingly used to study the deeper regions of the CZ (Befus et al., 2011; Olyphant et al., 2016; Callahan et al., 2020; Moravec et al., 2020). Seismic refraction uses elastic energy to produce 2D images of P-wave velocity on the order of 10s−100s of meters per deployment. P-wave velocity is a geophysical property that describes how fast compressional waves travel through a material. P-waves propagate by successively compressing and dilating the material through which they pass (e.g., Pride, 2005). Those compressions and dilatations are strains that occur at the scale of the seismic wavelength: that is, a long-wavelength wave strains a larger volume than a short-wavelength wave, as governed by the universal wave equation, v = λf where v is velocity, λ is wavelength, and f is frequency.
Seismic refraction methods are advantageous for CZ studies because (1) data can be collected quickly with a small team over large scales in rugged terrain; (2) the method quantifies P-wave velocity over 100s of meters; (3) P-wave velocity near Earth's surface is primarily controlled by porosity, a parameter highly relevant to CZ science; (4) the technique provides P-wave velocity estimates to depths of ~50 m or deeper, covering the full depth of the CZ in most areas; and (5) inversions can be run on laptop and desktop computers (Holbrook et al., 2014; Olyphant et al., 2016; Flinchum et al., 2018b; Seyfried et al., 2018; West et al., 2019).
The growing number of studies using seismic velocities to explore the deep CZ make it important for CZ scientists to understand what P-wave velocities can tell us about the CZ. In this paper, we focus on an apparent discrepancy between P-wave velocities collected via seismic refraction using a sledgehammer source (wavelengths of 10s of meters) and P-wave velocities collected from a sonic logging tool at multiple borehole locations (wavelengths of 10s of cm). Although both seismic refraction and sonic logging measure P-wave velocity, the methods have important differences in the scale over which they quantify velocity. To aid the reader in distinguishing the two measurements, in this paper we refer to P-wave velocities measured by the sonic log as “sonic velocities” and those measured by seismic refraction as “refraction velocities.”
This study is motivated, in part, by recent results reported by Holbrook et al. (2019), who showed a detailed comparison of geophysical and geochemical properties in the deep CZ in the Piedmont of South Carolina. They showed, based on both P-wave velocities and geochemical data, that a ~20-m-thick zone of fractured and weathered bedrock exists between saprolite and intact bedrock. This layer contains fractures that reduce seismic velocity and allow chemical weathering (Holbrook et al., 2019), perhaps through the opening of fractures by tectonic and topographic stresses (St. Clair et al., 2015; Moon et al., 2017). A layer of fractured and weathered bedrock is likely to be a ubiquitous feature of the CZ, but its properties and potential roles in geochemical and hydrological processes are poorly understood.
At the South Carolina site, both seismic refraction and sonic velocities were measured (Figure 2). Two contrasts between the sonic velocities and refraction velocities are immediately apparent (Figure 2). First, the sonic velocities are highly variable (changes > 1,500 m/s) above ~40 m depth (Figure 2B), where refraction velocities (and the entire model in Figure 2A) lack sharp changes. Second, the refraction velocity in the fractured bedrock layer is substantially slower than the average sonic velocity. This latter contrast is fundamental, as it cannot be explained as a simple or harmonic average of the sonic velocities. These apparent discrepancies immediately raise several questions: (1) Are refraction velocities accurate in the deep CZ? (2) If so, what accounts for the discrepancy between the seismic and sonic velocities? and (3) What scales of velocity heterogeneity in the deep CZ can be detected by seismic refraction methods?
Figure 2
We use a detailed set of borehole and seismic refraction observations in a weathered and fractured granite in the Laramie Range, Wyoming, to show the scales over which seismic refraction data provide valid information about CZ structure. Our observations in the Laramie Range are consistent with those at the South Carolina site (Holbrook et al., 2019), in that they show a distinct discrepancy between sonic (10s of cm scale) and seismic refraction velocities (10s of m scale) in the fractured rock region of the CZ. However, at this site we have another observation, 500 m away in a similar rock type, where the refraction and sonic velocities agree. Here we demonstrate that seismic velocity in the CZ is scale-dependent, and that the discrepancy between sonic and refraction velocities is not an error, but rather information about fracturing. This discrepancy, attributed to the contrasting measurement scales, is observed in other disciplines as well and is similar to measuring hydraulic conductivity on a small core sample and comparing this value to a large-scale pump test, or comparing the geochemical composition of a small piece of core against the composition of cuttings averaged over 10s of meters. We hypothesize that the magnitude of the velocity discrepancy derives from the scale of heterogeneities present in the fractured rock layer of the CZ. Although refraction images are unable to resolve small-scale heterogeneity (10s of cm to a few meters), the resulting structure is still clearly influenced by this heterogeneity. When interpreted in this context, we argue that seismic refraction remains a crucial geophysical tool to characterize the deep CZ structure across landscapes.
Site Description
Our study site is in the Laramie Range, Wyoming, where Sherman granite weathers to a thick grus layer (Eggler et al., 1969; Evanoff, 1990; Frost et al., 1999). The field site sits on gently undulating topography, which has been interpreted to be an ancient weathering surface that was uplifted during the Laramide Orogeny and is still being eroded (Gregory and Chase, 1994; Chapin and Kelley, 1997). Outcrops of Sherman granite exist in the northwest and southwest corners of the study area but are otherwise absent; one notable outcrop is located ~100 m southeast of hole BW-3 (Figure 3). There is no evidence of Quaternary glaciation, but it is likely that the site was much colder than the modern average temperature of 5.4°C. The site receives 620 mm of precipitation annually, most of which is in the form of snow (NRCS, 2021).
Figure 3
Over the course of 3 years (2013–2016) 30 seismic refraction profiles, totaling more than 7 km in length, were collected by the Wyoming Center for Hydrology and Geophysics (Flinchum et al., 2018a,b, 2019; Novitsky et al., 2018; Hayes et al., 2019). The deep CZ structure has also been imaged with passive seismic methods (Keifer et al., 2019; Wang et al., 2019a,b). The site contains nine boreholes ranging in depth from 16 to 73 m (Figure 3) (Flinchum et al., 2018a; Ren et al., 2019). Downhole optical televiewer logs and casing depths revealed that intact, unweathered bedrock, defined seismically by velocities >4,000 m/s, lies at depths up to 60 m below the surface (Flinchum et al., 2018b). The bedrock velocity of 4,000 m/s was determined from a seismic refraction survey conducted over a relatively unfractured outcrop nearby (Flinchum et al., 2018b). The soft and friable saprolite above the weathered bedrock interface was sampled through a combination of push coring and seismic refraction analysis and shown to be porous (~34%) and have velocities >1.1–1.2 km/s (Flinchum et al., 2018a, 2019).
Methods
Measuring P-Wave Velocity With Shallow Seismic Refraction and Sonic Logging
P-wave velocities presented in this manuscript were acquired by two methods: shallow seismic refraction (SSR) and sonic logging (SL). Both methods rely on the refraction of elastic energy at elastic boundaries. The elastic properties are dependent on the density and bulk and shear moduli of a given material. Refraction is advantageous because the velocity of a material can be determined simply by the first-arrival time of the seismic energy—the phase and amplitude of the arrival can be ignored.
During a SSR survey, the travel times between a seismic source to an array of geophones are measured and inverted to generate a 2D tomogram of seismic velocity (Woodward, 1992; Sheehan et al., 2005; Gance et al., 2012; Zelt et al., 2013). This type of SSR surveying is a robust method that has been used on crustal scales (Mooney and Brocher, 1987; Korenaga et al., 2000; Hayes et al., 2013), hillslope scales (Whiteley and Eccleston, 2006; McClymont et al., 2012; Thayer et al., 2018; Leone et al., 2020), and over outcrops (Holbrook et al., 2014; Flinchum et al., 2018b). SSR sources utilize a low-frequency (10s of Hz) source that produces relatively large wavelengths (10s−100s of meters). In an SSR survey seismic energy can travel 10s−100s of meters both vertically and laterally before being recorded at a geophone.
SL is less common in the CZ because it requires an uncased, large-diameter (2” or greater) borehole. In addition, sonic logs only work beneath the water table, since water couples the sonic source and receivers to the surrounding rock. Due to these constraints, only a limited number of studies have used sonic velocities to study the top 100 m of Earth's surface (Watson et al., 2005; Martí et al., 2006; Holbrook et al., 2014; Raef et al., 2015; Flinchum et al., 2018a; Gu et al., 2020a). The acquisition of a sonic velocity via SL works on the same principles as the SSR survey; that is, velocity is obtained by measuring the time taken to travel between a source and multiple receivers. In the SL case, however, most of the energy is traveling vertically, at right angles to the majority of travel paths in SSR.
There are a few important differences between the SSR and SL. First, during a SL survey, the maximum offset between the source and receivers will be limited by the tool length, often to 10′s of centimeters to a few meters, in contrast to the 10′s−100's of meters that might separate a source and geophone during SSR. The tool is moved up and down the hole in centimeter increments, providing estimates of velocity every few centimeters. Second, the sonic tool uses a source frequency that is 2–3 orders of magnitude higher than those used in a seismic refraction survey (our study is comparing ~80–15,000 Hz). The higher frequencies produce wavelengths 2–3 orders smaller than refraction. Finally, given that water couples the source and receiver to the borehole wall, the first arrival from energy will be critically refracted along the water/rock interface making it sensitive to a thin region along the borehole wall.
Shallow Seismic Refraction (SSR)
We reprocessed two SSR refraction lines that were part of a larger study presented by Flinchum et al. (2018a,b). L29, which extends from the top of the gently undulating ridge into the valley to the east, intersects boreholes BW-1 and BW-4 and comes within 10 m of holes BW-6 and BW-7 (Figure 3). A sledgehammer source was used by striking a stainless steel plate; hammer swings were stacked at least 8 times to improve the signal-to-noise. L29 was collected with 240 channels at 1 m spacing with shots every 12 m. L25 was collected with 96 channels at 2.5 m spacing with shots every 10 m. On L25, three off-end shots extending 30 m before the first geophone and 30 m past the last geophone were also collected. There were no off-end shots collected on L29.
For the SSR surveys, arrival times were picked manually. To pick the arrivals, shot gathers were trace normalized. We made as many picks as possible on the trace-normalized data without aid of any filtering. After picking on the unfiltered data, we applied a band-pass filter between 10 and 150 Hz to help identify far-offset arrivals. After picking, L29 contained 3935 first arrival picks and L25 contained 2039 picks (Supplementary Figure S1). In cases where the first arrival was noisy or could not be determined, no pick was made, to avoid adding errors into the final models.
To invert the first-arrival times and produce a 2D image of seismic velocities we use the Python Geophysical Inversion Modeling Library (PyGIMLi) (Rücker et al., 2017). PyGIMLi utilizes a shortest path algorithm (Dijkstra, 1959; Moser, 1991; Moser et al., 1992), which models seismic energy propagation as rays, not as waves. The shortest path algorithm is used in many other travel-time tomography codes (e.g., Zelt et al., 2013). The shortest path algorithm calculates a path through the velocity model—referred to as a ray path. The ray paths are used to mask areas in the final velocity that are not constrained by the data. In other words, if a ray path does not pass through a given cell in the model it will not be plotted. Ray paths tend to be curved in appearance because they refract toward the surface when velocity increases with depth. As a result, the resolved depths of the models vary spatially as driven by the data.
PyGIMLi utilizes a deterministic Gauss-Newton inversion scheme (Gance et al., 2012) which requires a data weight for each travel-time pick. These data weights are used to drive the error weight matrix and ultimately influence the final model fits. Here, instead of assigning a constant error value to each pick we used a linear equation as a function of offset (distance from source) line length, where offset is the distance [L] between the source and the geophone. Using Equation 1 means that the errors start at 0.5 ms and increase linearly to a maximum of 3 ms at the spread length. We chose a maximum error of 3 ms because of the high quality of the data. Each associated error is read into the inversion and used to weight the travel-times, such that travel times with the smallest errors are given the most weight.
To compare P-wave velocities derived from SSR with those from SL, it is important to quantify uncertainty on the velocities produced by the tomography. To quantify uncertainty, we randomly resampled the travel-time data, to utilize the large amount of input travel-times and randomly remove influences for each one. To achieve this, both L29 and L25 were inverted 20 times on the same mesh with the same inversion parameters, with 40% of the travel-time picks randomly removed before each inversion. Velocities at each pixel in the mesh were averaged together from all inversion realizations to produce an average 2D velocity model as well-standard deviations expressed as a percentage of velocity (Supplementary Figures S2, S3). We then use the average velocity profile to trace all the picked arrival times, so that the final fitting criteria for the tomograms shown here use all of the first arrival picks. We also compute the vertical velocity gradient (dv/dz) for each averaged tomographic model to highlight vertical changes in velocity.
Sonic Logging (SL)
For the SL acquisition, we used a Mount Sopris 2AFWS (Colorado, USA). The downhole instrument is ~ 3.0 m long with an insulator separating the source from the three receivers. The distance between the source and receivers is fixed at 0.91, 1.22, and 1.52 m. The source has a center frequency of 15 kHz. To process the downhole data, we recorded the elastic wavefield response at each of the three receivers. The elastic wavefield recording allows us to determine P-wave velocities by picking the first arrival times. As the tool travels through the borehole a source is generated at 5 cm intervals. To aid in the picking of arrivals, a semblance plot is created by summing the amplitudes of the waveforms along fifteen different slowness values (1/velocity) using a 50 μs window for every depth sampled. The best-fitting velocity from each source location (every 5 cm) was auto-picked based on the highest semblance value. After the automatic picks were made across the entire borehole depth, the picks are visually inspected, and erroneous picks are manually adjusted or removed. Because water in the borehole is required to couple the seismic source and geophones, sonic velocity logs do not exist for the unsaturated zone near the surface. Our data set consists of P-wave velocity profiles below the casing and below the water table for all five boreholes in this study (Figure 3).
The Relationship Between Velocity and Porosity
In areas with crystalline bedrock, we can assume that the porosity of intact bedrock is near zero, so that changes in the seismic velocity must be due to porosity generated through chemical and physical weathering processes (including fractures). The strong relationship between porosity and velocity has allowed for the estimation of porosity in the CZ from seismic refraction data (Holbrook et al., 2014; Flinchum et al., 2018a; Hayes et al., 2019; Callahan et al., 2020). Several rock physics models link porosity to velocity, including porous media models (Hashin and Shtrikman, 1963; Dvorkin and Nur, 1996; Dvorkin et al., 1999) and differential effective media (DEM) (Budiansky and Oconnell, 1976; Berryman et al., 2002a,b). A DEM model is appropriate when the pore spaces do not form a single connected network (Berryman et al., 2002b), consistent with the conceptual framework of the fractured rock region (Berryman et al., 2002b).
Regardless of which rock physics model is used, an increase in porosity will always result in a decrease in P-wave velocity if lithology and saturation remain constant. The rate of velocity decrease depends on the pore structure and mineralogy. In fractured media models, small increments in porosity can significantly reduce velocity. Following work by Flinchum et al. (2018a), we modelled porosities and velocities measured on granitic rocks around the world using a DEM derived by (Berryman et al., 2002b), which assumes interconnected penny-shaped inclusions throughout the material. We obtain the bulk and shear moduli of the unfractured medium using the Voigt-Reuss-Hill average (Hill, 1963; Mavko et al., 2009). The effective elastic moduli and of the fractured material are then obtained by solving a system of differential equations from (Berryman et al., 2002b) (Supplementary Material S1). Flinchum et al. (2018a) showed that a DEM model with a crack aspect ratio of 0.016 fit porosity and velocity measurements of weathered granites from multiple studies (Figure 4) (Begonha and Braga, 2002; Sousa et al., 2005; Novakova et al., 2012). The porosity-velocity relationship defined previously by Flinchum et al. (2018a) was calculated using a DEM rock physics model with an approximation for the Sherman granite (Table S1) (Flinchum et al., 2018a). The DEM model produces an exponential decrease in seismic velocity with increasing porosity (Figure 4). Based on the model and observed data (Figure 4) a porosity decrease of 0.10 m3/m3 can result in a velocity increase of ~5,000 m/s. As a result, we expect seismic velocity to be highly sensitive to variations in porosity or fracture density in the CZ.
Figure 4
Results
Heterogeneity
The CZ structure in the weathered and fractured granite at the study site is laterally heterogeneous. A Mount Sopris QL40-OBI-2G Optical Televiewer was deployed to provide a 360° unwrapped optical borehole image (OBI) at each of the boreholes. We converted the images to black and white and boosted the contrast to emphasize fractures and heterogeneity in the rock. The OBI logs and sonic logs based on both seismic velocities from the four deep (> 50 m) boreholes, all located within 70 m of each other, show strong lateral heterogeneity in the fractured rock directly below the casing (Figure 5). Downhole hydraulic conductivity logging revealed that hydraulic conductivity varies over two orders of magnitude but generally decreases with depth (Ren et al., 2019). Other geophysical studies have shown that the influence of heterogeneous remnant fabric from fractures extends upward into the saprolite (Novitsky et al., 2018) and that the depth to the intact bedrock is also laterally heterogeneous on the 100s−1,000s of m scale, extending as deep as 60 m below the surface (Flinchum et al., 2018b; Keifer et al., 2019; Wang et al., 2019a). In addition, nuclear magnetic resonance (NMR) logging was completed at all of the boreholes in the study area and show low but variable water contents (<0.15 m3/m3) in the fractured rock region of the CZ; the water is focused in visibly fractured regions (Flinchum et al., 2019; Ren et al., 2019). The fracture density in the fractured rock layer inferred from NMR is highly heterogeneous, with porosity varying from 0 to 0.15 m3/m3.
Figure 5
Comparing P-Wave Velocities Derived From SSR and SL
At L29, the SSR profile is generally divided into two domains. The first domain is defined by velocities < ~1,800 m/s (warmer colors in Figure 6A). The second domain is defined by velocities > ~4,000 m/s (cooler colors in Figure 6A). These two domains are separated by a thin transitional region where velocity increases rapidly over a short vertical distance (~5–10 m) with gradients of 100–300 m/s/m (Figure 6B). The vertical velocity gradient image (Figure 6B) highlights this transition or boundary separating these two domains. This boundary is deepest (> 50 m) under the ridge (between 25 and 100 m along the profile) (Figures 6A,B) and comes within 10 m of the surface under the valley (between 210 and 240 m along the profile) (Figures 6A,B). The strong increase in velocity can also be observed in the travel-time picks and represents the most significant seismological boundary revealed by the surface refraction data (Supplementary Figure S1).
Figure 6
Sonic velocity logs from BW1, BW-4, BW-6, and BW-7, which are on or near seismic profile L29 (Figure 3), show significant variability above 40 meters depth (Figure 6C). This variability is surprising, given that BW-4, BW-6, and BW-7 lie within 30 m of each other and BW-1 is only 70 m from BW-4 (Figure 3). Between 20 and 25 meters depth, for example, the velocity at BW-4 is ~3,000 m/s, while at BW-6 it is ~5,000 m/s. This variability disappears below about 40 m depth, where the sonic velocities converge to a velocity of ~4,500 m/s (Figure 6C).
Extracted 1D velocity profiles from the surface refraction data are smooth and lack both the detailed vertical structure and the lateral heterogeneity that is observed in the sonic velocities (Figure 6C). Above 40 m depth, where strong variability in sonic velocities occurs, refraction velocities are slower by as much as 4,000 m/s (Figure 6C). Despite the large differences in the shallow parts of the borehole, the tomographic velocity and sonic velocities converge to a common velocity between 4,000 and 5,000 m/s (Figure 6C), close to the 4,000–4,500 m/s measured by surface refraction on granite outcrops in the Southern Sierras (Holbrook et al., 2014) and nearby in the Laramie Range (Flinchum et al., 2018b). The discrepancy between the sonic velocities and travel-time velocities above 40 m is not a processing or inversion artifact (Figure 6C). We conducted a forward model study using the sonic velocities to predict SSR travel times for L29. The forward modeling results show clearly that the sonic velocities are too fast to explain the observed travel times (see Supplementary Materials S2, S3). Furthermore, similar discrepancies have been observed in the Piedmont of South Carolina (Figure 2) (Holbrook et al., 2019) and at a minimum of three other sites worldwide (Watson et al., 2005; Martí et al., 2006; Raef et al., 2015).
Just 500 m northwest of profile L29, the sonic and refraction velocities on line L25 show a much better match at all depths with overlapping observations (Figure 7). As with L29, the SSR profile can be divided into two domains: one with velocities < ~1,800 m/s, and one with velocities > ~4,000 m/s (Figure 7A). In this second refraction profile (Figure 7A), the seismic refraction velocity shows a high velocity zone that comes to the surface and extends 100s of meters laterally (Figure 7A).
Figure 7
The L25 profile crossed a depression in the landscape, so the water table is at ~7 m depth. This shallow water table allowed us to acquire sonic velocities closer to the surface than in any of the other boreholes (Figure 3). At this location, the sonic velocities and the seismic refraction velocities are in much better agreement (Figure 7C). The match is not perfect, however: the SSR velocity at L25 is about ~10% higher (~4,300 compared to ~3,900 m/s) than the sonic velocities observed in the borehole (Figure 7C). This is most likely caused by the strong vertical velocity gradient near the surface, which is a difficult structure for travel-time tomography to resolve. This area of the model has our highest uncertainty estimates (Supplementary Figure S3e). The uncertainty in the strong vertical velocity gradient near the surface would have the tendency to make the velocities near the surface slightly slower and the velocities at depth slightly faster to fit the observed travel times. Thus, we believe that the 1D velocities at L25 (Figure 7C) agree within error. Additionally, the surface refraction velocity profile still lacks the velocity detail observed in the sonic profile, and both velocities converge to an average velocity of ~4,000–4,500 m/s, consistent with velocities collected over outcrop in the Southern Sierras (Holbrook et al., 2014) and the Laramie Range (Flinchum et al., 2018b).
Discussion
The Seismic/Sonic Velocity Discrepancy
One of the challenges of studying CZ processes is integrating a wide variety of observational data across different spatial scales (e.g., mineralogical vs. hillslope vs. landscape scales). The discrepancy between refraction and sonic velocities (Figures 2B, 6C) might appear to beg an obvious question: which velocity is “correct”? Here we argue that the answer is that both velocities are correct, at the wavelengths over which they are measured. Neither velocity is in error; rather, the velocities differ because they sense fracture densities that are scale-dependent. Before presenting our argument, we must first address whether the observed sonic/refraction velocity discrepancy is due to vertical anisotropy. Anisotropy is the ratio between seismic velocities measured in different directions. In general, seismic velocities travel slower across foliation or fractures (Novitsky et al., 2018; Eppinger et al., 2021). The seismic waves produced by SSR travel largely in the horizontal direction, whereas the seismic waves produced by SL largely travel in the vertical direction, raising the idea that pervasive, near-vertical fractures would cause slower refraction velocities than sonic velocities. However, vertical anisotropy cannot explain the velocity discrepancy we observe for several reasons. First, the anisotropy values required to explain the SSR/SL discrepancy would need to be as high as ~100%, which is far higher than observed anisotropy due to fracturing elsewhere (e.g., Crampin et al., 1980; Maultzsch et al., 2003). Second, faster velocities in the vertical direction would imply pervasive near-vertical fracturing. While fracturing exists at our site, the amount of anisotropy observed in the saprolite (Novitsky et al., 2018) is modest (maximum 12%), and we are unaware of processes that would explain 100% anisotropy in the fractured bedrock with only ~10% anisotropy in the saprolite. Third, numerous near-horizontal fractures are observed in the borehole optical logs (Figures 5A–D)—a pattern that should create faster velocities in the horizontal direction. Finally, decreases in sonic velocity between the boreholes at Blair correlate strongly with fracture density observed in the boreholes (Figures 5A–D). This suggests that sonic velocities are controlled by the number of fractures in a given volume (i.e., fracture density), not the orientation of the fractures.
In the case of our sonic and seismic waves passing through the fractured bedrock layer, the wavelengths differ by two orders of magnitude (v = fλ): a seismic P-wave with a dominant frequency of 40 Hz traveling at 2,000 m/s would have a dominant wavelength of 50 m, while a sonic P-wave at 6 kHz and 3,000 m/s has a wavelength of 50 cm (Figure 8). P-waves propagate by compressing and dilating the material at the scale of the seismic wavelength. This, in turn, means that the seismic wave senses all the fractures present in a volume ~50 m in diameter, while the sonic wave only senses fractures within a 0.5 m diameter. The volume that is deformed by this compression is related to the wavelength and is known as the representative elementary volume. The volume can change with scale depending on the heterogeneity of the material. The representative elementary volume is an important concept that can used to explain the variation of a physical property with upscaling (Matonti et al., 2015; Bailly et al., 2019).
Figure 8
To use a hydrogeology analogy, we postulate that comparing sonic and refraction velocities is similar to comparing hydraulic conductivity measurements of a fractured aquifer obtained with a pump (i.e., transmissivity) test with measurements made in the lab on a small sample. Apparent hydraulic conductivity is well-known to scale with the measurement, where small scale measuremetns tend to show lower hydraulic conducitvites (Neuman, 1988; Clauser, 1992; Rovey and Cherkauer, 1995; Schulze-Makuch and Cherkauer, 1998). In most cases, depending on how the sample was obtained, the lab sample will have a lower hydraulic conductivity because the pump test is an average of high and low hydraulic conductivity fractures along the entire length of the borehole (Clauser, 1992) (equivalent to the large yellow area in Figure 8). It is not uncommon in fractured rock aquifers to have one or two highly productive fractures while others do not contribute to flow (Sausse and Genter, 2005; Lee and Kim, 2015).
The analogy of hydraulic conductivity is useful but imperfect, in one respect: in the case of hydraulic conductivity, the length of the screened area or the location of the sample are known, while in the case of seismic refraction, the wavefield is complex (Supplementary Figure S5), and wavelength varies as a function of velocity and source frequency (Supplementary Figures S4, S5). The velocity is affected by the porosity distribution of the CZ, which is, by definition, site specific. (Liu et al., 1976; Xia et al., 2003; Tisato and Marelli, 2013). Fluid in pore space can cause frequency dependent attenuation in the form of Wave Induced Fluid Flow (WIFF) (Müller et al., 2010). In seismic refraction, the elastic energy can travel 100s of meters from source to receiver (Figure 8). The combination of these physical phenomena makes it difficult to explicitly define the sample region represented by the travel-time observation.
We hypothesize that velocities measured across such vastly different scales in a porous medium would only be equal if the fraction of pore space was equal across those scales—essentially, if the distribution of pore space is uniform. Refraction and sonic velocities would match if we collected a seismic refraction profile over a borehole with sonic velocities in a location where the underlying rock has a homogenous porosity distribution—that is, if porosity is equal regardless of the sample volume. While such a fairly uniform distribution of porosity can occur in a sedimentary rock (e.g., a clean sandstone or a gravel aquifer), it is unlikely to be true when porosity is due to fractures. Fractured materials can have highly non-uniform fracture distributions across a wide variety of scales. We illustrate this point with a simple example: a layer of bedrock pervasively broken by evenly spaced fractures (Figure 8).
If fractures are 1 m apart and each have an aperture of 4 cm, then within the 50-m wavelength of a seismic P-wave, there will be 50 × 0.04 = 2 m of open fracture space, or a fracture porosity of 4% (2/50). Given the strong sensitivity of P-velocities to fracture porosity, adding 4% porosity to the bedrock could easily reduce seismic velocities from, say, 4,500–2,000 m/s (Figure 4). The seismic wave travels at this relatively slow velocity throughout the fractured bedrock layer, and its travel times in the refraction data reflect that slow propagation. The sonic tool, in contrast, with its wavelength of around 50 cm, would measure non-weathered and unfractured rock (velocity 4,000–4,500 m/s) except when crossing an individual fracture (Figure 8), when the velocity would drop precipitously—exactly as observed in the sonic logs of Figures 6B, 7B.
For this reason, we hypothesize that the persistent discrepancy between seismic and sonic velocities in the fractured bedrock layer of the CZ in South Carolina (Holbrook et al., 2019) (Figure 2) and the Laramie Range (Figure 6) should be viewed not as an error, but rather as valuable information. Seismic velocities that consistently underestimate sonic velocities are a sign of pervasive fracturing at scales between the respective seismic and sonic wavelengths. Where the two velocities agree, in contrast, then porosity is likely uniformly distributed. This is the case in Figures 2, 6B at depths >40 m, where the relatively high velocities suggest that the meter-scale fractures are largely closed due to ambient pressure (e.g., St. Clair et al., 2015).
What Do P-Wave Velocities Tell Us About CZ Structure?
Under our proposed hypothesis, velocity models produced by seismic refraction tomography will be blind to small-scale (10s of meters) heterogeneity but will still be able to quantify the thickness of the CZ across landscapes over 100s of m. Seismic refraction still remains a particularly powerful geophysical tool to explore the deeper regions of the CZ, especially if the CZ boundaries of interest are controlled by changes in porosity or fracture density. The velocity profiles derived from seismic refraction represent large-scale averages that reflect the dominant seismic wavelength. It remains difficult to quantify the volume that this average represents, since quantifying the volumes would require a fully 3D seismic waveform solution. Full waveform inversions, even in 2D, are still in their infancy when working in the strong velocity gradients associated with near-surface materials (Gao et al., 2006, 2007; Virieux and Operto, 2009; Romdhane et al., 2011).
Despite the lack of sensitivity to small-scale structures, seismic refraction remains one of the only geophysical tools that can quantify the general location of boundaries defined by changes in porosity in the CZ over landscape to watershed scales. The seismic refraction method suffers from the classic tradeoff between velocity detail and spatial coverage, and users need to be aware of the relevant wavelength scales to interpret the images correctly. Although high-frequency sonic velocity logs in the CZ remain rare, our data show that P-wave seismic velocity is scale-dependent. Given the interdisciplinary nature of CZ science, it is likely that the 2D profiles of velocity will be used to inform or test conceptual models of CZ evolution (West et al., 2019; Leone et al., 2020). It is important to realize that the velocities in seismic refraction tomography represent large-scale averages over the scale of their dominant wavelength. While smaller-scale heterogeneity certainly exists beneath the lateral and vertical resolution of the refraction inversion, the velocity discrepancy itself is not a consequence of the inversion, but rather of the physics of wave propagation at different wavelengths. Understanding this large-scale averaging is particularly important at sites where velocities are used quantitatively to estimate porosity (Mota and Santos, 2010; Flinchum et al., 2018a; Callahan et al., 2020; Gu et al., 2020b) or strain (Hayes et al., 2019).
Our results show that the refraction and sonic velocities converge at depths where the fracture density is uniform, consistent with the idea of a protolith that is defined by a small, uniform porosity. The sample volume is not straightforward, as wave propagation physics is complicated, but a good rule of thumb is that the wavelengths, which depend on the velocity and frequency content of the source, provide a rough estimate of the averaging scale that the velocities represent.
Conclusions
P-wave velocities estimated from surface refraction tell us about the CZ structure on scales of 10s−100s of m. Sonic velocities are more sensitive to fracture density and CZ heterogeneity. The velocity profiles generated by seismic refraction tomography are blind to small-scale (< 10 m) heterogeneity but successfully image the thickness of the CZ across landscapes. The lack of sensitivity to smaller-scale structures means that refraction velocities must be interpreted in the context of their wavelengths, as they will not be able to determine if a boundary is sharp, gradational, or laterally heterogeneous. Despite the lack of sensitivity to smaller-scale structures, surface seismic velocity remains one of the most useful physical properties for detecting the bulk transitions from weathered to less weathered bedrock.
There is still much to learn from seismic velocities of CZ structure. There still exists a spatial gap between the cm scale (sonic velocities) and 10s−100s of m scale (seismic refraction). The lack of high-frequency velocity measurements on samples, either by sonic logs or laboratory tests, leaves a significant gap in our understanding of the scalability of seismic velocities in the CZ. Our work shows that mismatches between sonic and refraction velocities are properly viewed not as artifacts or errors, but rather as information about the spacing and density of fractures.
Funding
This work was supported by NSF (EPS-1208909, 1531313, and EAR- 2012353).
Publisher's Note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Statements
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author. Geophysical data used in this article is hosted on Hydroshare and can be downloaded here: http://www.hydroshare.org/resource/5c6dd575e27f44159f826e4c1c20bf19 (Flinchum, 2021).
Author contributions
All authors have been actively involved in data collection, processing, and development of ideas presented in this manuscript.
Acknowledgments
We are grateful to many others for thoughtful discussions and constructive criticism to improve the manuscript. We are also grateful for Dr. Ty Ferre for very insightful comments during review of this manuscript.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/frwa.2021.772185/full#supplementary-material
References
1
AndersonS. P.BlanckenburgF.von WhiteA. F. (2007). Physical and chemical controls on the critical zone. Elements3, 315–319. 10.2113/gselements.3.5.315
2
AndersonS. P.DietrichW. E.BrimhallG. H. (2002). Weathering profiles, mass-balance analysis, and rates of solute loss: Linkages between weathering and erosion in a small, steep catchment. Gsa Bull.114, 1143–1158. 10.1130/0016-7606(2002)114<1143:WPMBAA>2.0.CO;2
3
BaillyC.FortinJ.AdelinetM.HamonY. (2019). Upscaling of elastic properties in carbonates: a modeling approach based on a multiscale geophysical data set. J. Geophys. Res. Solid Earth124, 13021–13038. 10.1029/2019JB018391
4
BalesR. C.HopmansJ. W.O'GeenA. T.MeadowsM.HartsoughP. C.KirchnerP.et al. (2011). Soil moisture response to snowmelt and rainfall in a sierra nevada mixed-conifer forest. Vadose Zone J.10, 786–799. 10.2136/vzj2011.0001
5
BefusK. M.SheehanA. F.LeopoldM.AndersonS. P.AndersonR. S. (2011). Seismic constraints on critical zone architecture, boulder creek watershed, front range, Colorado. Vadose Zone J.10, 1342–1342. 10.2136/vzj2010.0108er
6
BegonhaA.BragaM. A. S. (2002). Weathering of the Oporto granite: geotechnical and physical properties. Catena49, 57–76. 10.1016/S0341-8162(02)00016-4
7
BerrymanJ. G.BergeP. A.BonnerB. P. (2002a). Estimating rock porosity and fluid saturation using only seismic velocities. Geophysics67, 391–404. 10.1190/1.1468599
8
BerrymanJ. G.PrideS. R.WangH. F. (2002b). A differential scheme for elastic properties of rocks with dry or saturated cracks. Geophys. J. Int.151, 597–611. 10.1046/j.1365-246X.2002.01801.x
9
BrantleyS. L.HolleranM. E.JinL.BazilevskayaE. (2013). Probing deep weathering in the Shale Hills Critical Zone Observatory, Pennsylvania (USA): the hypothesis of nested chemical reaction fronts in the subsurface. Earth Surf. Processes38, 1280–1298. 10.1002/esp.3415
10
BrantleyS. L.McDowellW. H.DietrichW. E.WhiteT. S.KumarP.AndersonS. P.et al. (2017). Designing a network of critical zone observatories to explore the living skin of the terrestrial Earth. Earth Surf. Dynam.5, 841–860. 10.5194/esurf-5-841-2017
11
BrantleyS. L.WhiteA. F. (2009). Approaches to modeling weathered regolith. Rev. Mineral. Geochem.70, 435–484. 10.2138/rmg.2009.70.10
12
BraunJ.MercierJ.GuillocheauF.RobinC. (2016). A simple model for regolith formation by chemical weathering. J. Geophys. Res. Earth Surf.121, 2140–2171. 10.1002/2016JF003914
13
BrooksP. D.ChoroverJ.FanY.GodseyS. E.MaxwellR. M.McNamaraJ. P.et al. (2015). Hydrological partitioning in the critical zone: Recent advances and opportunities for developing transferable understanding of water cycle dynamics. Water Resour. Res.51, 6973–6987. 10.1002/2015WR017039
14
BudianskyB.OconnellR. (1976). Elastic-moduli of a cracked solid. Int. J. Solids Struct.12, 81–97. 10.1016/0020-7683(76)90044-5
15
BussH. L.LaraM. C.MooreO. W.KurtzA. C.SchulzM. S.WhiteA. F. (2017). Lithological influences on contemporary and long-term regolith weathering at the Luquillo Critical Zone Observatory. Geochim. Cosmochim Ac196, 224–251. 10.1016/j.gca.2016.09.038
16
CallahanR. P.RiebeC. S.PasquetS.FerrierK. L.GranaD.SklarL. S.et al. (2020). Subsurface weathering revealed in hillslope-integrated porosity distributions. Geophys. Res. Lett.47. 10.1029/2020GL088322
17
CanadellJ.JacksonR. B.EhleringerJ. B.MooneyH. A.SalaO. E.SchulzeE.-D. (1996). Maximum rooting depth of vegetation types at the global scale. Oecologia108, 583–595. 10.1007/BF00329030
18
CareyA. M.PaigeG. B.CarrB. J.HolbrookW. S.MillerS. N. (2019). Characterizing hydrological processes in a semiarid rangeland watershed: a hydrogeophysical approach. Hydrol. Process.33, 759–774. 10.1002/hyp.13361
19
ChapinC. E.KelleyS. A. (1997). The Rocky Mountain erosion surface in the Front Range of Colorado, in Colorado Front Range guidebook: Denver, eds BolyardD. W.SonnenbergS. A. (Cham: Rocky Mountain Association of Geologists), 101–114.
20
ClarkeB. A.BurbankD. W. (2011). Quantifying bedrock-fracture patterns within the shallow subsurface: Implications for rock mass strength, bedrock landslides, and erodibility. J. Geophys. Res. Earth Surf.116:F04009. 10.1029/2011JF001987
21
ClauserC. (1992). Permeability of crystalline rocks. Eos. Trans. Am. Geophys. Union73, 233–238. 10.1029/91EO00190
22
CrampinS.McGonigleR.BamfordD. (1980). Estimating crack parameters from observations of P-wave velocity anisotropy. Geophysics45, 345–360. 10.1190/1.1441086
23
DethierD. P.LazarusE. D. (2006). Geomorphic inferences from regolith thickness, chemical denudation and CRN erosion rates near the glacial limit, Boulder Creek catchment and vicinity, Colorado. Geomorphology75, 384–399. 10.1016/j.geomorph.2005.07.029
24
DewandelB.LachassagneP.WynsR.MaréchalJ. C.KrishnamurthyN. S. (2006). A generalized 3-D geological and hydrogeological conceptual model of granite aquifers controlled by single or multiphase weathering. J. Hydrol.330, 260–284. 10.1016/j.jhydrol.2006.03.026
25
DijkstraE. W. (1959). A note on two problems in connexion with graphs. Numer. Math.1, 269–271. 10.1007/BF01386390
26
DixonJ. C.ThornC. E. (2005). Chemical weathering and landscape development in mid-latitude alpine environments. Geomorphology67, 127–145. 10.1016/j.geomorph.2004.07.009
27
DralleD. N.HahmW. J.RempeD. M.KarstN. J.ThompsonS. E.DietrichW. E. (2018). Quantification of the seasonal hillslope water storage that does not drive streamflow. Hydrol. Process.32, 1978–1992. 10.1002/hyp.11627
28
DvorkinJ.NurA. (1996). Elasticity of high-porosity sandstones; theory for two North Sea data sets. Geophysics61, 1363–1370. 10.1190/1.1444059
29
DvorkinJ.PrasadM.SakaiA.LavoieD. (1999). Elasticity of marine sediments: rock physics modeling. Geophys. Res. Lett.26, 1781–1784. 10.1029/1999GL900332
30
EgglerD. H.LarsonE. E.BradleyW. C. (1969). Granites grusses and sherman erosion surface southern laramie range colorado-wyoming. Am. J. Sci.267:510. 10.2475/ajs.267.4.510
31
EppingerB. J.HayesJ. L.CarrB. J.MoonS.CosansC. L.HolbrookW. S.et al. (2021). Quantifying depth-dependent seismic anisotropy in the critical zone enhanced by weathering of a piedmont schist. J. Geophys. Res. Earth Surf.126:e2021JF006289. 10.1029/2021JF006289
32
EvanoffE. (1990). Early Oligocene paleovalleys in southern and central Wyoming: Evidence of high local relief on the late Eocene unconformity. Geology18, 443–446. 10.1130/0091-7613(1990)018<0443:EOPISA>2.3.CO;2
33
FletcherR. C.BussH. L.BrantleyS. L. (2006). A spheroidal weathering model coupling porewater chemistry to soil thicknesses during steady-state denudation. Earth Planet. Sc. Lett.244, 444–457. 10.1016/j.epsl.2006.01.055
34
FlinchumB. A. (2021). Supporting Data Repository for paper for What do P-wave Velocites Tell Us about the Critical Zone, HydroShare. Available online at: http://www.hydroshare.org/resource/5c6dd575e27f44159f826e4c1c20bf19
35
FlinchumB. A.HolbrookW. S.GranaD.ParsekianA. D.CarrB. J.HayesJ. L.et al. (2018a). Estimating the water holding capacity of the critical zone using near-surface geophysics. Hydrol. Process.32, 3308–3326. 10.1002/hyp.13260
36
FlinchumB. A.HolbrookW. S.ParsekianA. D.CarrB. J. (2019). Characterizing the critical zone using borehole and surface nuclear magnetic resonance. Vadose Zone J.18, 1–18. 10.2136/vzj2018.12.0209
37
FlinchumB. A.HolbrookW. S.RempeD.MoonS.RiebeC. S.CarrB. J.et al. (2018b). Critical zone structure under a granite ridge inferred from drilling and three-dimensional seismic refraction data. J. Geophys. Res. Earth Surf.123, 1317–1343. 10.1029/2017JF004280
38
FrostC. D.FrostB. R.ChamberlainK. R.EdwardsB. R. (1999). Petrogenesis of the 1.43 Ga Sherman Batholith, SE Wyoming, USA: a reduced, rapakivi-type anorogenic granite. J. Petrol40, 1771–1802. 10.1093/petroj/40.12.1771
39
GabetE. J.MuddS. M. (2009). A theoretical model coupling chemical weathering rates with denudation rates. Geology37, 151–154. 10.1130/G25270A.1
40
GanceJ.GrandjeanG.SamynK.MaletJ.-P. (2012). Quasi-Newton inversion of seismic first arrivals using source finite bandwidth assumption: application to subsurface characterization of landslides. J. Appl. Geophys.87, 94–106. 10.1016/j.jappgeo.2012.09.008
41
GaoF.LevanderA.PrattR. G.ZeltC. A.FradelizioG.-L. (2007). Waveform tomography at a groundwater contamination site: surface reflection dataWaveform tomography: reflection data. Geophysics72, G45–G55. 10.1190/1.2752744
42
GaoF.LevanderA. R.PrattR. G.ZeltC. A.FradelizioG. L. (2006). Waveform tomography at a groundwater contamination site: VSP-surface data setHigh-Resolution Velocity Modeling. Geophysics71, H1–H11. 10.1190/1.2159049
43
GolebiowskiT.PorzucekS.PasierbB. (2016). Ambiguities in geophysical interpretation during fracture detection—case study from a limestone quarry (Lower Silesia Region, Poland). Near Surf. Geophys.14, 371–384. 10.3997/1873-0604.2016017
44
GoodfellowB. W.ChadwickO. A.HilleyG. E. (2014). Depth and character of rock weathering across a basaltic-hosted climosequence on Hawai‘i. Earth Surf. Processes39, 381–398. 10.1002/esp.3505
45
GrahamR.RossiA.HubbertR. (2010). Rock to regolith conversion: producing hospitable substrates for terrestrial ecosystems. Gsa Today20, 4–9. 10.1130/GSAT57A.1
46
GrahamR. C.AndersonM. A.SternbergP. D.TiceK. R.SchoenebergerP. J. (1997). Morphology, porosity, and hydraulic conductivity of weathered granitic bedrock and overlying soils. Soil Sci. Soc. Am. J.61:516. 10.2136/sssaj1997.03615995006100020021x
47
GregoryK. M.ChaseC. G. (1994). Tectonic and climatic significance of a late Eocene low-relief, high-level geomorphic surface, Colorado. J. Geophys. Res. Solid Earth99, 20141–20160. 10.1029/94JB00132
48
GuX.HeaneyP. J.ReisF. D. A. A.BrantleyS. L. (2020a). Deep abiotic weathering of pyrite. Science370:eabb8092. 10.1126/science.abb8092
49
GuX.MavkoG.MaL.OakleyD.AccardoN.CarrB. J.et al. (2020b). Seismic refraction tracks porosity generation and possible CO2 production at depth under a headwater catchment. Proc. Natl. Acad. Sci.U.SA.117, 18991–18997. 10.1073/pnas.2003451117
50
HahmW. J.RempeD. M.DralleD. N.DawsonT. E.LovillS. M.BrykA. B.et al. (2019). Lithologically controlled subsurface critical zone thickness and water storage capacity determine regional plant community composition. Water Resour. Res.55, 3028–3055. 10.1029/2018WR023760
51
HahmW. J.RiebeC. S.LukensC. E.ArakiS. (2014). Bedrock composition regulates mountain ecosystems and landscape evolution. Proc. Natl. Acad. Sci. U.S.A. 111, 3338–3343. 10.1073/pnas.1315667111
52
HasenmuellerE. A.GuX.WeitzmanJ. N.AdamsT. S.StinchcombG. E.EissenstatD. M.et al. (2017). Weathering of rock to regolith: the activity of deep roots in bedrock fractures. Geoderma300, 11–31. 10.1016/j.geoderma.2017.03.020
53
HashinZ.ShtrikmanS. (1963). A variational approach to the theory of the elastic behaviour of multiphase materials. J. Mech. Phys. Solids11, 127–140. 10.1016/0022-5096(63)90060-7
54
HayesJ. L.HolbrookW. S.LizarraldeD.AvendonkH. J. A.BullockA. D.MoraM.et al. (2013). Crustal structure across the Costa Rican Volcanic Arc. Geochem. Geophys. Geosystems.14, 1087–1103. 10.1002/ggge.20079
55
HayesJ. L.RiebeC. S.HolbrookW. S.FlinchumB. A.HartsoughP. C. (2019). Porosity production in weathered rock: where volumetric strain dominates over chemical mass loss. Sci. Adv.5:eaao0834. 10.1126/sciadv.aao0834
56
HillR. (1963). Elastic properties of reinforced solids: some theoretical principles. J. Mechan. Phys. Solids11, 357–372. 10.1016/0022-5096(63)90036-X
57
HilleyG. E.ChamberlainC. P.MoonS.PorderS.WillettS. D. (2010). Competition between erosion and reaction kinetics in controlling silicate-weathering rates. Earth Planet. Sci. Lett.293, 191–199. 10.1016/j.epsl.2010.01.008
58
HolbrookW. S.MarconV.BaconA. R.BrantleyS. L.CarrB. J.FlinchumB. A.et al. (2019). Links between physical and chemical weathering inferred from a 65-m-deep borehole through Earth's critical zone. Sci. Rep.9:4495. 10.1038/s41598-019-40819-9
59
HolbrookW. S.RiebeC. S.ElwaseifM.HayesJ. L.Basler-ReederK.HarryD. L.et al. (2014). Geophysical constraints on deep weathering and water storage potential in the Southern Sierra Critical Zone Observatory. Earth Surf. Processes39, 366–380. 10.1002/esp.3502
60
KeiferI.DuekerK.ChenP. (2019). Ambient Rayleigh wave field imaging of the critical zone in a weathered granite terrane. Earth Planet. Sc Lett.510, 198–208. 10.1016/j.epsl.2019.01.015
61
KorenagaJ.HolbrookW. S.KentG. M.KelemenP. B.DetrickR. S.LarsenH. -C.et al. (2000). Crustal structure of the southeast Greenland margin from joint refraction and reflection seismic tomography. J. Geophys. Res. Solid Earth105, 21591–21614. 10.1029/2000JB900188
62
KotikianM.ParsekianA. D.PaigeG.CareyA. (2019). Observing heterogeneous unsaturated flow at the hillslope scale using time-lapse electrical resistivity tomography. Vadose Zone J.18, 1–16. 10.2136/vzj2018.07.0138
63
KuntzB. W.RubinS.BerkowitzB.SinghaK. (2011). Quantifying solute transport at the shale hills critical zone observatory. Vadose Zone J.10, 843–857. 10.2136/vzj2010.0130
64
LeeH.KimB. (2015). Characterisation of hydraulically-active fractures in a fractured granite aquifer. Water Sa41, 139–148. 10.4314/wsa.v41i1.17
65
LeoneJ. D.HolbrookW. S.RiebeC. S.ChoroverJ.FerréT. P. A.CarrB. J.et al. (2020). Strong slope-aspect control of regolith thickness by bedrock foliation. Earth Surf Processes. 30, 2998–3010. 10.1002/esp.4947
66
LeopoldM.VölkelJ.HuberJ.DethierD. (2013). Subsurface architecture of the Boulder Creek Critical Zone Observatory from electrical resistivity tomography: critical zone architecture at boulder creek critical zone observatory. Earth Surf. Processes38, 1417–1431. 10.1002/esp.3420
67
LiuH.AndersonD. L.KanamoriH. (1976). Velocity dispersion due to anelasticity; implications for seismology and mantle composition. Geophys. J. Int.47, 41–58. 10.1111/j.1365-246X.1976.tb01261.x
68
LovillS. M.HahmW. J.DietrichW. E. (2018). Drainage from the critical zone: lithologic controls on the persistence and spatial extent of wetted channels during the summer dry season. Water Resour. Res.54, 5702–5726. 10.1029/2017WR021903
69
MartíD.CarbonellR.Escuder-VirueteJ.Pérez-EstaúnA. (2006). Characterization of a fractured granitic pluton: P- and S-waves' seismic tomography and uncertainty analysis. Tectonophysics422, 99–114. 10.1016/j.tecto.2006.05.012
70
MatontiC.GuglielmiY.ViseurS.BrunaP. O.BorgomanoJ.DahlC.et al. (2015). Heterogeneities and diagenetic control on the spatial distribution of carbonate rocks acoustic properties at the outcrop scale. Tectonophysics638, 94–111. 10.1016/j.tecto.2014.10.020
71
MaultzschS.ChapmanM.LiuE.LiX. Y. (2003). Modelling frequency-dependent seismic anisotropy in fluid-saturated rock with aligned fractures: implication of fracture size estimation from anisotropic measurements. Geophys. Prospect.51, 381–392. 10.1046/j.1365-2478.2003.00386.x
72
MavkoG.MukerjiT.DvorkinJ. (2009). The Rock Physics Handbook, (New York, NY: Cambridge Press), 21–80. 10.1017/CBO9780511626753
73
McClymontA. F.HayashiM.BentleyL. R.LiardJ. (2012). Locating and characterising groundwater storage areas within an alpine watershed using time-lapse gravity, GPR and seismic refraction methods. Hydrol. Process.26, 1792–1804. 10.1002/hyp.9316
74
MoonS.PerronJ. T.MartelS. J.HolbrookW. S.ClairJ.St. (2017). A model of three-dimensional topographic stresses with implications for bedrock fractures, surface processes, and landscape evolution. J. Geophys. Res. Earth Surf.122, 823–846. 10.1002/2016JF004155
75
MooneyW. D.BrocherT. M. (1987). Coincident seismic reflection/refraction studies of the continental lithosphere: a global review. Rev. Geophys.25, 723–742. 10.1029/RG025i004p00723
76
MoravecB. G.WhiteA. M.RootR. A.SanchezA.OlshanskyY.ParasB. K.et al. (2020). Resolving deep critical zone architecture in complex volcanic terrain. J. Geophys. Res.125:e2019JF005189. 10.1029/2019JF005189
77
MoserT. (1991). Shortest-path calculation of seismic rays. Geophysics56, 59–67. 10.1190/1.1442958
78
MoserT. J.NoletG.SniederR. (1992). Ray bending revisited. Bull. Seismol. Soc. Am.82, 259–288.
79
MotaR.SantosF. A. M. (2010). 2D sections of porosity and water saturation from integrated resistivity and seismic surveys. Near Surf. Geophys.8, 575–584. 10.3997/1873-0604.2010042
80
MüllerT. M.GurevichB.LebedevM. (2010). Seismic wave attenuation and dispersion resulting from wave-induced flow in porous rocks — a reviewSeismic wave attenuation. Geophysics75, 75A147–75A164. 10.1190/1.3463417
81
Navarre-SitchlerA.BrantleyS. L.RotherG. (2015). How porosity increases during incipient weathering of crystalline silicate rocks. Rev. Mineral. Geochem.80, 331–354. 10.2138/rmg.2015.80.10
82
NeumanS. P. (1988). Groundwater Flow and QualityModelling, (Norwell, MA: Kluwer Academic), 331–362. 10.1007/978-94-009-2889-3_19
83
NovakovaL.SosnaK.Bro,ŽM.NajserJ.NovakP. (2012). The matrix porosity and related properties of a leucocratic granite from the Krudum Massif, West Bohemia. Acta Geodyn. Geomater.9, 521–540.
84
NovitskyC. G.HolbrookW. S.CarrB. J.PasquetS.OkayaD.FlinchumB. A. (2018). Mapping inherited fractures in the critical zone using seismic anisotropy from circular surveys. Geophys. Res. Lett.45, 3126–3135. 10.1002/2017G.L.075976
85
NRCS SSS (2021). Natural Resources Conservation Service, United States Department of Agriculture. Soil Climate Research Station Data Crow Creek, Wyoming. Available online at: https://wcc.sc.egov.usda.gov/nwcc/site?sitenum=1045
86
NyquistJ. E.DollW. E.DavisR. K.HopkinsR. A. (1996). Cokriging surface topography and seismic refraction data for bedrock topography. J. Environ. Eng. Geoph.1, 67–74. 10.4133/JEEG1.1.67
87
OlyphantJ.PelletierJ. D.JohnsonR. (2016). Topographic correlations with soil and regolith thickness from shallow-seismic refraction constraints across upland hillslopes in the Valles Caldera, New Mexico. Earth Surf. Processes41, 1684–1696. 10.1002/esp.3941
88
OrlandoJ.ComasX.HynekS. A.BussH. L.BrantleyS. L. (2016). Architecture of the deep critical zone in the Río Icacos watershed (Luquillo Critical Zone Observatory, Puerto Rico) inferred from drilling and ground penetrating radar (GPR). Earth Surf. Processes41, 1826–1840. 10.1002/esp.3948
89
ParsekianA. D.SinghaK.MinsleyB. J.HolbrookW. S.SlaterL. (2015). Multiscale geophysical imaging of the critical zone. Rev. Geophys.53, 1–26. 10.1002/2014RG000465
90
PawlikŁ.PhillipsJ. D.ŠamonilP. (2016). Roots, rock, and regolith: Biomechanical and biochemical weathering by trees and its impact on hillslopes—a critical literature review. Earth-Sci. Rev.159, 142–159. 10.1016/j.earscirev.2016.06.002
91
PrideS. R. (2005). Relationships between seismic and hydrological properties, in Water Science and Technology Library, eds RubinY.HubbardS. S. (Dordrecht: Springer Netherlands), 253–290. 10.1007/1-4020-3102-5_9
92
RaefA.GadS.Tucker-KuleszaS. (2015). Multichannel analysis of surface-waves and integration of downhole acoustic televiewer imaging, ultrasonic Vs and Vp, and vertical seismic profiling in an NEHRP-standard classification, South of Concordia, Kansas, USA. J. Appl. Geophys.121, 149–161. 10.1016/j.jappgeo.2015.08.005
93
RempeD. M.DietrichW. E. (2018). Direct observations of rock moisture, a hidden component of the hydrologic cycle. Proc. Natl. Acad. Sci. U.S.A.115, 2664–2669. 10.1073/pnas.1800141115
94
RenS.ParsekianA. D.ZhangY.CarrB. J. (2019). Hydraulic conductivity calibration of logging NMR in a granite aquifer, Laramie Range, Wyoming. Groundwater57, 303–319. 10.1111/gwat.12798
95
RichterD.deB.MobleyM. L. (2009). Monitoring earth's critical zone. Science326, 1067–1068. 10.1126/science.1179117
96
RichterD. D.MarkewitzD. (1995). How deep is soil?Bioscience45, 600–609. 10.2307/1312764
97
RiebeC. S.HahmW. J.BrantleyS. L. (2017). Controls on deep critical zone architecture: a historical review and four testable hypotheses. Earth Surf. Processes42, 128–156. 10.1002/esp.4052
98
RomdhaneA.GrandjeanG.BrossierR.RejibaF.OpertoS.VirieuxJ. (2011). Shallow-structure characterization by 2D elastic full-waveform inversionRomdhane et al.Elastic FWI for shallow structures. Geophysics76, R81–R93. 10.1190/1.3569798
99
RoveyC. W.CherkauerD. S. (1995). Scale dependency of hydraulic conductivity measurements. Groundwater33, 769–780. 10.1111/j.1745-6584.1995.tb00023.x
100
RückerC.GüntherT.WagnerF. M. (2017). pyGIMLi: an open-source library for modelling and inversion in geophysics. Comput. Geosci.109, 106–123. 10.1016/j.cageo.2017.07.011
101
SalveR.RempeD. M.DietrichW. E. (2012). Rain, rock moisture dynamics, and the rapid response of perched groundwater in weathered, fractured argillite underlying a steep hillslope. Water Resour. Res.48:W11528. 10.1029/2012WR012583
102
SausseJ.GenterA. (2005). Types of permeable fractures in granite. Geol. Soc. Lond. Spec. Publ.240, 1–14. 10.1144/GSL.SP.2005.240.01.01
103
Schulze-MakuchD.CherkauerD. S. (1998). Variations in hydraulic conductivity with scale of measurement during aquifer tests in heterogeneous, porous carbonate rocks. Hydrogeol. J.6, 204–215. 10.1007/s100400050145
104
SeyfriedM.LohseK.MarksD.FlerchingerG.PiersonF.HolbrookW. S. (2018). Reynolds creek experimental watershed and critical zone observatory. Vadose Zone J.17, 1–20. 10.2136/vzj2018.07.0129
105
SheehanJ. R.DollW. E.MandellW. A. (2005). An evaluation of methods and available software for seismic refraction tomography analysis. J. Environ. Eng. Geophys.10, 21–34. 10.2113/JEEG10.1.21
106
SlimM.PerronJ. T.MartelS. J.SinghaK. (2015). Topographic stress and rock fracture: a two-dimensional numerical model for arbitrary topography and preliminary comparison with borehole observations. Earth Surf. Processes40, 512–529. 10.1002/esp.3646
107
SousaL. M. O.RíoL. M. S.del CallejaL.ArgandoñaV. G. R.de ReyA. R. (2005). Influence of microfractures and porosity on the physico-mechanical properties and weathering of ornamental granites. Eng. Geol.77, 153–168. 10.1016/j.enggeo.2004.10.001
108
St. ClairJ.MoonS.HolbrookW. S.PerronJ. T.RiebeC. S.MartelS. J.et al. (2015). Geophysical imaging reveals topographic stress control of bedrock weathering. Science350, 534–538. 10.1126/science.aab2210
109
ThayerD.ParsekianA. D.HydeK.SpeckmanH.BeverlyD.EwersB.et al. (2018). Geophysical measurements to determine the hydrologic partitioning of snowmelt on a snow-dominated subalpine hillslope. Water Resour. Res.54, 3788–3808. 10.1029/2017WR021324
110
TisatoN.MarelliS. (2013). Laboratory measurements of the longitudinal and transverse wave velocities of compacted bentonite as a function of water content, temperature, and confining pressure. J. Geophys. Res.118, 3380–3393. 10.1002/jgrb.50252
111
VereeckenH.KasteelR.VanderborghtJ.HarterT. (2007). Upscaling hydraulic properties and soil water flow processes in heterogeneous soils: a review. Vadose Zone J.6, 1–28. 10.2136/vzj2006.0055
112
VirieuxJ.OpertoS. (2009). An overview of full-waveform inversion in exploration geophysics. Geophysics74, WCC1–WCC26. 10.1190/1.3238367
113
WalshD.TurnerP.GrunewaldE.ZhangH.ButlerJ. J.RebouletE.et al. (2013). A small-diameter NMR logging tool for groundwater investigations. Groundwater51, 914–926. 10.1111/gwat.12024
114
WangW.ChenP.KeiferI.DuekerK.LeeE.-J.MuD.et al. (2019a). Weathering front under a granite ridge revealed through full-3D seismic ambient-noise tomography. Earth Planet. Sc. Lett.509, 66–77. 10.1016/j.epsl.2018.12.038
115
WangW.ChenP.LeeE.-J.MuD. (2019b). Earthquake and Disaster Risk: Decade Retrospective of the Wenchuan Earthquake, (Singapore, Springer), 203–231. 10.1007/978-981-13-8015-0_8
116
WatsonD. B.DollW. E.GameyT. J.SheehanJ. R.JardineP. M. (2005). Plume and lithologic profiling with surface resistivity and seismic tomography. Ground Water43, 169–177. 10.1111/j.1745-6584.2005.0017.x
117
WestN.KirbyE.NybladeA. A.BrantleyS. L. (2019). Climate preconditions the Critical Zone: elucidating the role of subsurface fractures in the evolution of asymmetric topography. Earth Planet. Sc. Lett.513, 197–205. 10.1016/j.epsl.2019.01.039
118
WhiteleyR. J.EcclestonP. J. (2006). Comparison of shallow seismic refraction interpretation methods for regolith mapping. Explor. Geophys.37, 340–347. 10.1071/EG06340
119
WoodwardM. J. (1992). Wave-equation tomography. Geophysics57, 15–26. 10.1190/1.1443179
120
XiaJ.MillerR. D.ParkC. B.TianG. (2003). Inversion of high frequency surface waves with fundamental and higher modes. J. Appl. Geophys.52, 45–57. 10.1016/S0926-9851(02)00239-2
121
YaedeJ. R.McBrideJ. H.NelsonS. T.ParkC. B.FloresJ. A.TurnbullS. J.et al. (2015). A geophysical strategy for measuring the thickness of the critical zone developed over basalt lavas. Geosphere11, 514–532. 10.1130/GES01142.1
122
ZeltC. A.HainesS.PowersM. H.SheehanJ.RohdewaldS.LinkC.et al. (2013). Blind test of methods for obtaining 2-D near-surface seismic velocity models from first-arrival traveltimes. J. Environ. Eng. Geophys.18, 183–194. 10.2113/JEEG18.3.183
Summary
Keywords
critical zone, seismic refraction, velocities, critical zone architecture, geophysics, sonic velocities, multiple scales
Citation
Flinchum BA, Holbrook WS and Carr BJ (2022) What Do P-Wave Velocities Tell Us About the Critical Zone?. Front. Water 3:772185. doi: 10.3389/frwa.2021.772185
Received
07 September 2021
Accepted
29 November 2021
Published
10 January 2022
Volume
3 - 2021
Edited by
Marc Dumont, Sorbonne Universités, France
Reviewed by
Cédric Bailly, UMR8148 Géosciences Paris Sud (GEOPS), France; Ty Ferre, University of Arizona, United States
Updates
Copyright
© 2022 Flinchum, Holbrook and Carr.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Brady A. Flinchum bflinch@clemson.edu
This article was submitted to Water and Critical Zone, a section of the journal Frontiers in Water
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.