Atterholt and Zhan
,
Sci. Adv.
10
, eadr3327 (2024) 27 November 2024
ScieNce Adv
ANceS
|
ReSeAR
ch AR
ticle
1 of 11
GEOPHYSICS
Fine-
scale Southern California Moho structure
uncovered with distributed acoustic sensing
James Atterholt*† and Zhongwen Zhan
Moho topography yields insights into the evolution of the lithosphere and the strength of the lower crust. The
Moho reflected phase (PmP) samples this key boundary and may be used in concert with the first arriving P phase
to infer crustal thickness. The densely sampled station coverage of distributed acoustic sensing arrays allows for
the observation of PmP at fine-
scale intervals over many kilometers with individual events. We use PmP recorded
by a 100-
km-
long fiber that traverses a path between Ridgecrest, CA and Barstow, CA to explore Moho variability
in Southern California. With hundreds of well-
recorded events, we verify that PmP is observable and develop a
technique to identify and pick P-
PmP differential times with high confidence. We use these observations to con-
strain Moho depth throughout Southern California, and we find that short-
wavelength variability in crustal thick
-
ness is abundant, with sharp changes across the Garlock Fault and Coso Volcanic Field.
INTRODUCTION
Observations of Moho structure are important because they provide
evidence of processes that deform the lithosphere and yield con-
straints on the rheology of the lower crust (
1
,
2
). Crustal thickness
varies over many scales, and fine-
scale variability (on the order of a
few kilometers) is important because it constrains the depth-
extent
and depth-
dependent behavior of localized tectonic processes ob-
served at the surface (
3
–
6
). The most popular methods for uncovering
Moho structure are receiver functions and controlled-
source seismic
surveys. Receiver functions, which leverage secondary phases in tele-
seismic data, are a powerful tool for robustly determining Moho depth
(
6
–
9
), but are a low-
frequency measurement that requires seismome-
ters near the measurement points, limiting their spatial resolution and
geographic coverage. Controlled-
source surveys are performed by re-
cording numerous sources with shorter transitory arrays to produce
seismic reflection profiles (
10
–
12
) or by measuring source wavefields
over large distances with longer arrays in so-
called wide-
angle refrac-
tion and reflection surveys (
13
–
17
). These surveys can produce Moho
phases that can be used to uncover crustal thickness variability at im-
pressive resolution, but they are expensive, are logistically challenging,
and sometimes fail to penetrate Moho depths.
An alternative to these techniques is to use arrival times of the
Moho reflected phase (PmP) in regional earthquake wavefields, which
can be used to constrain Moho depth (
7
,
18
–
21
) and lower crustal
velocity structure (
22
–
25
). The advantage of using this observation
to constrain crustal thickness over other methods is that measure-
ments are made passively using higher-
frequency regional earthquakes
and, with the right seismicity distribution, can resolve structure out-
side a seismic network. This technique has not been used as frequently
as others because PmP is difficult to confidently identify on individual
seismograms (
20
). Recent studies have found success using multi-
event estimates of PmP moveout and machine learning techniques
to generate expanded PmP catalogs (
7
–
18
). These techniques are not
readily transferrable to or ideal for dense array datasets that generally
have lower data quality but high spatial sampling. The high-
frequency
content of PmP observations combined with dense array datasets may
allow
for very high-
resolution sampling of short-
wavelength features
near the crust-
mantle boundary.
Distributed acoustic sensing (DAS) is an increasingly popular tech-
nique in seismology that transforms fiber optic cables into dense arrays
of strainmeters using phase interferometry of backscattered light
(
26
). This technique is powerful because it enables the deployment of
extensive, dense seismic networks for long periods of time at low logisti-
cal burden. DAS is routinely used for shallow crustal imaging (
27
–
29
)
and has recently been used to resolve structure in the middle crust using
travel-
time tomography (
30
). The dense recordings enabled by DAS
facilitate short-
wavelength observations over large distances, and
the long deployment times allow for the passive recording of high-
quality earthquake wavefields with diverse source locations. While
receiver functions have been successfully computed by combining DAS
data with nearby broadband data (
31
), wider adoption of deep crustal
imaging techniques with DAS has been slow, because DAS has high
noise levels at low frequencies (
32
). PmP is an especially promising
avenue for realizing the potential of DAS for deep crustal imaging for a
few reasons. First, PmP is a high-
frequency phase, and DAS is most sen-
sitive to higher-
frequency wavefields. Second, the high spatial density
and abundance of channels in DAS arrays allows for both spatial coher
-
ence and array-
side moveout to be used to identify PmP. Third, although
DAS arrays have many channels, the cable geometry often imposes a
narrow geographic footprint; with PmP, Moho depth measurements
may be made over a much larger area than is covered by the fiber.
Here, we introduce a method with which to identify and pick
relative arrival times of secondary phases in DAS data. We apply this
method to a DAS array in the Mojave Desert to measure PmP differ
-
ential times and use these observations to invert for Moho depth over a
wide area. We find good correspondence with previous results, and
we observe short-
wavelength features in the Moho across the Garlock
Fault and the Coso Volcanic Field (CVF). This technique offers a prom-
ising outlook for using DAS arrays to make fine-
scale observations of
lithospheric structure over broad geographic areas.
RESULTS
Autocorrelation for secondary phase retrieval
In August of 2021, a DAS array was deployed on a 100-
km segment of
dark fiber between Ridgecrest, CA and Barstow, CA (Fig. 1). This array
california i
nstitute of
technology, Pasadena, c
A, USA.
*c
orresponding author. email: jatterholt@
usgs.
gov
†Present address: United States Geological Survey, Golden, c
O, USA.
copyright © 2024
the
Authors, some rights
reserved; exclusive
licensee American
Association for the
Advancement of
Science. No claim to
original U.S.
Government Works.
distributed under a
creative c
ommons
Attribution license 4.0
(cc BY).
Downloaded from https://www.science.org on December 03, 2024
Atterholt and Zhan
,
Sci. Adv.
10
, eadr3327 (2024) 27 November 2024
ScieNce Adv
ANceS
|
ReSeAR
ch AR
ticle
2 of 11
has 10,000 channels with 10-
m channel spacing and 100-
m gauge
length. Over the course of 2 years, the array recorded 440 M2.5
+
regional earthquakes. Many of these earthquake wavefields exhibit
strong secondary phase arrivals that coincide with the expected onset
time of the Moho reflected phase (PmP) for a standard earth model. The
relative arrival time between the first arriving P phase and PmP provides
a strong constraint on crustal thickness (Fig. 2). Coupled with the spatial
density of DAS, P-
PmP differential times can yield very short wave-
length profiles of Moho depth.
Identifying and measuring the relative arrival times of PmP on
DAS data can be a challenging problem. First, although PmP is often
the strongest phase in
P
- wave coda, it has historically been very chal-
lenging to identify in broadband data, with only
~
1% of records
showing readily identifiable PmP waveforms (
18
). DAS presents an
advantage toward identifying this phase; the spatial density of DAS
arrays allows for the evaluation of candidate phases across a broad
spatial window. This spatial window allows for the use of spatial co-
herence and array-
side moveout to assist in identifying secondary
phases like PmP. Once identified, however, it is impractical to pick the
thousands of phase arrivals recorded by DAS for every earthquake
wavefield by hand. Traditional automated phase picking methods are
usually inadequate for these data because P-
phase onsets are often
weak and complicated in DAS due to broadside insensitivity of axial
strainmeters and strong surface wave scattering due to local heteroge-
neities (
33
). Machine learning methods, such as PhaseNet DAS (
34
),
have shown promise toward addressing this problem, but as of now,
they have not been adapted to secondary phases.
We develop a simple and effective approach to both identify and
pick relative arrival times of secondary phases like PmP. Because the
scattered surface waves that are nearly ubiquitous in DAS data are
generated locally (
33
–
35
), these waves are common to both the first
arrival and secondary arrivals. These phases consequently generate
complicated but highly correlated wavefields. The autocorrelation of
DAS-
recorded earthquake wavefields is thus a potentially powerful
tool for making highly accurate relative phase arrival time measure-
ments. We apply a straightforward and semi-
automated framework
for identifying PmP and obtaining P-
to-
PmP differential times for
DAS data (see Materials and Methods). This methodology is outlined
using an example in Fig. 3. Of the initial set of events, 229 events
had visually identifiable P-
phase onsets and some high-
quality first
arrival picks. Of this smaller subset of events, we can observe and pick
PmP arrivals over large spatial windows for 72 events. Additionally,
as shown in Fig. 2, by evaluating the change in differential time
with source-
to-
receiver distance, we determine for each observation
whether PmP is trailing a direct wave (Pg) or a mantle head wave (Pn)
first arrival. From this dataset, our workflow yields over 200,000 accu
-
rate and precise P-
PmP differential times that may be used to con-
strain the velocity structure and crustal thickness throughout
Southern California. These picks and the corresponding correlo-
grams are shown in Fig. 4.
Fine-
scale crustal thickness variability from dense Moho
depth profiles
As shown in Fig. 2, the differential time between Pg and PmP is a
function of both the difference in path length between these two
phases, which is dependent on the Moho depth, and the difference in
the velocity structure to which these two phases are sensitive. Pn has
a longer ray path than PmP, but it arrives earlier than PmP and even-
tually overtakes Pg because the upper mantle velocity is much higher
Fig. 1.
Experiment setting and data.
the d
AS array geometry (blue curve) plotted
against the kept (green) and discarded (red) M2.5
+
events recorded by the array.
Also plotted are the fault traces included in the USGS Quaternary Fault database
(2018) and locations of relevant tectonic provinces and features (Gv
, Great
valley;
SN, Sierra Nevada; BR, Basin and Range; Md
, Mojave d
esert; cvF, c
oso
volcanic
Field; GF, Garlock Fault; SGM, San Gabriel Mountains; SBM, Santa Barbara Moun-
tains; iA, isabella Anomaly).
the green star corresponds to the event that produced
the wavefield shown in Fig. 3.
A
BC
Fig. 2.
Important phases and relative arrival times.
(
A
) Ray paths for the three phas-
es of interest in this study for a fixed source-
receiver distance and representative veloc
-
ity model: the direct phase (Pg), the Moho head wave (Pn), and the Moho reflected
wave (PmP) (S not pictured for simplicity). (
B
) Reduced (8.1 km/s) arrival times of phases
plotted against source-
receiver distance. (
C
) differential times of phases relative to the
first arriving P phase plotted against source-
receiver distance.
Downloaded from https://www.science.org on December 03, 2024
Atterholt and Zhan
,
Sci. Adv.
10
, eadr3327 (2024) 27 November 2024
ScieNce Adv
ANceS
|
ReSeAR
ch AR
ticle
3 of 11
than the crustal velocity. When using P-
PmP differential times, we
thus need to consider three things: crustal structure, upper mantle
velocity, and crustal thickness. Crustal structure and thickness trade
off with each other, and considering the two simultaneously would
result in non-
uniqueness, particularly because parameterizing the
crust requires both seismic velocities and a layering structure. Short-
wavelength variability in PmP differential times is much more sensi-
tive to changes in crustal thickness. This is because sensitivity due to
crustal structure depends on integrated velocities along the ray path,
not detailed velocity structure; we subsequently illustrate this with
synthetic tests. Since we are most interested in sharp changes in
the Moho, we choose to fix a representative crustal model using an
ensemble of velocity profiles (see Materials and Methods) and invert
for Moho topography. The distance at which Pn supersedes Pg as the
first arrival and the rate at which the Pn-
PmP differential time in-
creases with distance are functions of the upper mantle velocity inte-
grated along the ray path. While this value also trades off with
absolute Moho depth, we can resolve a robust estimate of a represen-
tative upper mantle velocity for the entire dataset in a straightforward
way. We can then incorporate this velocity into our model when in-
verting for crustal thickness.
Using a representative crustal model, we determine a best-
fitting
combination of Moho depth and upper mantle velocity for our entire
dataset (see Materials and Methods). We find that the optimal pair
is a Moho depth of 31.5 km and an upper mantle velocity of 8.1 km/s.
The ensemble of crustal models is shown in fig. S1. This Moho depth
is very similar to that of the Hadley-
Kanamori (HK) model for
Southern California (32 km) (
36
), and the upper mantle velocity
Fig. 3.
Example of autocorrelation for phase retrieval.
(
A
) earthquake wavefield with representative data quality. i
ncluded are the PhaseNet d
AS picks and computed
P and PmP times for the representative velocity model used in this study with a Moho depth of 32 km. Green bars are a reference for the channel bounds within which
differential time picks were made. (
B
) Autocorrelation image created for the corresponding earthquake wavefield using the framework outlined in the text. i
ncluded are
the computed P-
PmP differential times and the picks made on the positive correlation peak associated with the PmP arrival. l
ocation of this event is indicated as a green
star in Fig. 1.
Fig. 4.
Summary of all PmP picks made in this study.
Plotted are the autocorrelation correlograms for each channel for which estimates of P-
PmP differential time could
be made. Black dotted lines correspond to pick times, and correlograms are organized by decreasing pick times when following Pg and increasing pick times when fol-
lowing Pn. Gray dotted line marks the transition from Pg to Pn as first arrival.
Downloaded from https://www.science.org on December 03, 2024
Atterholt and Zhan
,
Sci. Adv.
10
, eadr3327 (2024) 27 November 2024
ScieNce Adv
ANceS
|
ReSeAR
ch AR
ticle
4 of 11
is consistent with Pn tomography results from this region (
20
,
37
).
We then use our fixed crustal and optimal upper mantle velocities
combined with our P-
PmP differential times and first arrival classi-
fications to invert for crustal thickness along each of our bounce-
point profiles (see Materials and Methods). The resultant profiles
provide estimates of Moho depth at sub-
kilometer intervals, which
is smaller than the minimum wavelength expected in the earth-
quake wavefields used in this experiment. These Moho depth esti-
mates and select cross sections are shown in Fig. 5.
These profiles yield Moho topography that is in close agreement
with previously resolved long-
wavelength features in Southern
California (
38
) and provide remarkably high spatial frequency reso-
lution along profiles that illuminate highly localized changes in crustal
thickness. Broadly, we observe a Moho depth of approximately 30 km
throughout most of the Mojave block that shallows to the east. We also
resolve a deep Moho over and around the Isabella anomaly (
39
) that has
been regularly observed in this region (
40
). Consistent with other
studies, the crust thickens sharply at the transition to the San Bernadino
Mountains, corresponding to the mountain root (
11
,
40
). The Moho
is relatively deep just to the northwest of these mountains, which is con-
sistent with some recent results (
7
). In general, the absolute values of
Moho depth across the region agree with previous localized mea-
surements of Moho depth made using receiver functions (
9
), as is illus-
trated in Fig. 5 and fig. S2.
There are trade-
offs between variability in the velocity model and
crustal thickness, and we evaluate and discuss these trade-
offs later.
However, the primary advantage of combining a high-
frequency
measurement like P-
PmP differential times and the high spatial
sampling of DAS is toward resolving fine-
scale variability in crustal
thickness, rather than broad-
scale absolute values. Velocity model
variability will generally exert a longer wavelength effect on P-
PmP
differential times because the velocity contribution to the differential
time is path integrated. Crustal thickness variability, however, can
exert shorter wavelength effects on P-
PmP differential times, because
a sharp change in the Moho depth abruptly changes the path length
of PmP. We demonstrate this using a set of synthetic tests (see Mate-
rials and Methods) shown in Fig. 6. These tests show that even very
sharp and high-
amplitude local velocity perturbations have a smooth
and muted effect on P-
PmP differential time with source-
to-
receiver
distance. Comparatively, a Moho step produces an abrupt jump in
the P-
PmP differential time that is representative of the amplitude of
the step in the Moho.
Sharp Moho changes across the Garlock Fault and CVF
Our results show that there is considerable variation in Moho depth
along our crustal thickness profiles. We choose to focus on two short-
wavelength features that are observable only because of the unique
combination of the station density of DAS and the regional sensitivity
of the PmP measurements. The correlogram profiles that inform
these observations are shown in Figs. 7 and 8. One of these interesting
features shown in Fig. 5 is a sharp step (3 to 4 km) in the Moho that is
very close to the surface trace of the western segment of the Garlock
Fault. As shown in Fig. 7, this step is clear in the data and occurs over
only a few kilometers. Because this profile is recovered using indi-
vidual event wavefields, uncertainties related to source location do
not affect the relative change. A less precise but more objective way
to confirm the presence of the step is to sum the envelope of these
correlograms along the expected P-
PmP differential times for a range
Fig. 5.
Summary of resolved Moho depth variability.
(
A
) Moho depth inversion results along corresponding bounce-
point locations for all events included in this study.
included as squares are Moho depth estimates reported in previous study (
9
). Purple dotted line marks the
dAS array location. (
B
)
cross sections of Moho depth corre
-
sponding to the profiles mapped in (A). i
ncluded are all estimates of Moho depth within 10 km of the profile reported in this study (black dots) and all estimates of Moho
depth within 20 km of the profile reported in (
9
) (orange squares). l
ettering indicates different structural features where they cross the profiles (llFZ, little lake Fault
Zone).
the cvF location marks the approximate southern extent of the mid-
crustal low-
velocity zone reported in (
55
). lateral dotted line corresponds to smoothed Moho
depth estimates from the c
ommunity Moho Model (
38
). (
C
) Schematic diagrams of the cross sections in (B).
Downloaded from https://www.science.org on December 03, 2024
Atterholt and Zhan
,
Sci. Adv.
10
, eadr3327 (2024) 27 November 2024
ScieNce Adv
ANceS
|
ReSeAR
ch AR
ticle
5 of 11
of Moho depths. When we do this along the spatial windows to the
north and south of the expected step location, we find that there are
peaks at distinct Moho depths (Fig. 9). Another shorter profile to the
east, also shown in Fig. 7, shows a step-
like feature of the same polar
-
ity, providing additional evidence for this Moho discontinuity.
This observation may be attributed to a discontinuity in the crust
imposed by the Garlock Fault and would suggest that the western seg-
ment of the fault penetrates the Moho. This observation is distinct from
some early results suggesting that the Garlock Fault is truncated (
10
) or
approaches a horizontal angle in the middle crust (
41
). The results pre-
sented in this study have deeper penetration depths and require fewer
assumptions than these earlier studies. The deeper penetration depths
resolved in this study add increased importance to the Garlock Fault as
a physical boundary between the Mojave and the Sierra Nevada and
Basin and Range terrains to the north (Fig. 1). This agrees with some
geologic studies that suggest the Garlock Fault delineates, and slip on
the fault may be driven by, a difference in extensional behavior between
the Basin and Range and the Mojave (
42
). Additionally, this feature
suggests that along this segment of the fault the Garlock is nearly verti-
cal through the entire crust. This extends the results of earlier studies
that used focal mechanisms (
43
) and imaging techniques (
44
,
45
) to
infer a near-
vertical dipping fault at seismogenic depths. The width
of this step is intriguing, because the behavior of continental faults at
depth can shed light on the strength of the lower crust and upper man-
tle (
6
,
46
). Since the Garlock Fault consists of a narrow step, this obser
-
vation supports the possibility of a narrow shear zone at depth, rather
than a broad deformational zone. The wavelength of this step observed
at the Garlock Fault (
<
3 km) suggests an upper bound for the width of
the deformation zone at the crust-
mantle boundary. The PmP observa-
tion becomes faint in the narrow spatial interval at the step, and we can-
not further characterize the properties of the step within this interval.
That such a narrow zone with steep offset is maintained on a relatively
slow slipping fault (
47
) may help constrain the strength contrast at the
crust mantle boundary. However, future study is required to quantify
the rheological conditions necessary to preserve this step. This observa-
tion is comparable to that of the Denali fault, another intracontinental
AB
Fig. 6.
Synthetic tests.
(
A
)
velocity models and source-
array geometry (shown as blue line at surface) used in these synthetic tests. d
ouble-
sided arrows correspond to
the approximate extent and locations of bounce points in these tests. Solid black line marks the Moho, and dotted black lines outline velocity anomalies. (
B
) differential
times from autocorrelation wavefields for each synthetic test. Shown are reference relative arrival times of PmP for a flat Moho at 25 and 30 km depth. Black arrow indicates
the diffracted phase.
Downloaded from https://www.science.org on December 03, 2024