Seismic Data Crunching for Fun and Profit

1
Seismic Data Crunching for Fun and Profit

• Peter Shearer
• IGPP, U.C. San Diego
• AGU Gutenberg Lecture
• December 5, 2005

2
Seismic Data Crunching for Fun and Profit

• Guy Masters, Paul Earle, Luciana Astiz, Megan
  Flanagan, Michael Hedlin, Egill Hauksson,
  Guoqing Lin, German Prieto, Jesse Lawrence,
  Heidi Houston, John Vidale, Miaki Ishii, Kris
  Walker

3
Beno Gutenberg

• Core radius properties
• Global travel-time curves
• Ocean/continent differences
• Mantle low-velocity zone
• Gutenberg-Richter magnitude
• Magnitude/energy relations
• b-value
• Microseisms
• Four major books, 300 articles

1889 1960
4
Gutenbergs 1959 book, Physics of the Earths
Interior, concludes with
THE DATA MUST BE GREATLY AMPLIFIED AND
STRENGTHENED.
This might be called the motto of this book.
5
Thats one down and only 3,827 to go!
6
1970s Digital seismograms
7
1989 CD-ROM data distribution

• 650 MB capacity
• Cheap to produce
• Selected events released by NEIC
• First reasonably practical access to large global
  datasets for individual seismologist

8
(No Transcript)
9
Seismic Data Crunching Rules

• Analyze entire dataset whenever possible
• Use simple methods to get sense of data before
  doing complicated inversions
• Consider reflection seismology methods like
  stacking and back-projection
• Avoid any hand-processing of seismograms!

10
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image earthquake rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

11
Global Stacking using Automatic Gain Control (AGC)

• Calculate average absolute value in 5 s bins
• Divide each bin by average of previous 24 bins.
  This normalizes the amplitude of each trace.
• Stack in 0.5 distance bins

12
AGC Stack Long-period vertical
90
60
Time (minutes)
30
0
0
180
90
Distance (degrees)
from Shearer (1991)
13
AGC Stack Long-period transverse
90
60
Time (minutes)
30
0
0
180
270
360
90
Distance (degrees)
from Shearer (1991)
14
SH is faster than SV
SH (solid) vs. SV (dashed) from stacked
seismograms
S (17 23)
SS (37 43)
SSS (57 63)
Distance (degrees)
from Shearer (1991)
15
AGC Stack Short-period vertical
from Astiz et al. (1996)
16
19881994 IRIS Farm archive 834
earthquakes 27,000 seismograms
Vertical
Transverse
Radial
from Astiz et al. (1996)
17
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image earthquake rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

18
Stacking using a reference phase
Unaligned SH waves
Aligned SH waves
1 minute
Stack
Reference pulse stacks for 20 different range bins
19
CD-ROM stacks (1991)
P wave (vertical)
Topside reflections
PP
P
660-km discontinuity
S wave (transverse)
SS
410-km discontinuity
No global 220-km discontinuity
S
20
CD-ROM stack SS precursors
SS-wave stack (transverse)
4
SS
2
S660S
SS
0
Time (minutes)
-2
410-km discontinuity
520
-4
660
-6
-8
80
100
120
140
160
180
Range (degrees)
Sdiff
from Shearer (1991)
21
No coherent reflectors above 410 or below 660
22
SS precursors are ideal for global mantle
discontinuity studies
Source
Bounce point
Receiver
Good global distribution of bounce points
from Flanagan Shearer (1998)
23
Depression in 660 in NW Pacific
from Shearer (1991)
24
660 topography from SS precursors
Shearer Masters (1992)
Flanagan Shearer (1998)
Gu et al. (2002)
blue depressed (1020 km)
red elevated
25
CD-ROM stacks (1991)
SV/P discontinuity conversions (Faber Muller,
1984)
P-wave stack (radial)
PcSdiff
P/SV discontinuity conversions (Vinnik, 1977)
26
Receiver functions at GSN stations
Shearer (1991)
Lawrence Shearer (2005)
27
A discrepancy in receiver function versus SS
precursor results

• Easiest to compare transition zone thickness
  (WTZ d660 - d410) because it is insensitive to
  structure above 410 km
• Chevrot et al. (1999) performed global receiver
  function analysis, found a difference in WTZ with
  SS precursor results
• Why is mean WTZ 10 km greater in Chevrot et al.
  study compared to SS precursor results?

28
New global receiver function study

• 118 stations
• Mean WTZ 242 2 km
• Agrees with SS precursors

Bias in prior studies due to
Plane wave approximation Non-uniform
station coverage
Lawrence Shearer (2005)
29
Transition Zone Thickness Models
Receiver functions
SS precursors
Gu et al. (1998)
Lawrence Shearer (2005)
Flanagan Shearer (1998)
30
Slabs in the transition zone
P-wave tomography
660 topography
Flanagan Shearer (1998)
from Karson and van der Hilst (2000)
31
Slabs in the transition zone
Response of 660-km discontinuity to slab
50100 km deflection in vicinity of slab
Lesser deflection in large region beneath slab
32
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image earthquake rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

33
Stacking complications at short periods
34
Stacking incoherent waves
1 Hz seismograms are incoherent
35
Uniquely valuable PKP precursors

• Core P velocity drop bends rays so that scattered
  waves can arrive before direct phases

36
PKP Precursor Examples
PKP
Precursors
First observed by Gutenberg and Richter in 1934!
37
Average PKP Precursor Wavefield
10
0
PKP(DF)
Precursors
Time (s)
-10
-20
120
145
125
130
140
135
Range (degrees)
38
Predicted PKP precursor envelopes for scattering
at different depths
39
Stacks suggest whole mantle scattering
120
125
130
Range (degrees)
135
140
145
-20
0
10
-10
Time (s)
40
PKP Precursor Interpretation

• 0.5 to 1 RMS velocity perturbations at 10 km
  scale length
• Recent analyses show scattering extends at least
  1000 km above CMB (Hedlin et al., 1997 Cormier,
  1999 Margerin Nolet, 2003).
• Early studies put scattering near CMB (e.g.,
  Cleary Haddon, 1972)

Mantle
CMB
Core
41
Pdiff coda provides more evidence for mid-mantle
scattering
Good fit to data stack obtained with 1 RMS
velocity heterogeneity throughout the mantle
coda
Pdiff
from Earle Shearer (2001)
42
Mantle mixing calculations
Heterogeneity is likely at all scales
Davies (2002)
Xie Tackley (2004)
43
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

44
Relative event location using waveform
cross-correlation
Nakamura (1978) Poupinet (1984)
Ito (1985) Fremont Malone
(1987) Xie et al. (1991) Got
et al. (1994) Haase et al. (1995)
Dodge et al. (1995) Nadeau et al. (1995)
Phillips et al. (1997) Rubin et al. (1999)
Waldhauser et al. (1999)
45
Cross-correlation of pair of similar events
P-waves
Quake 1
Quake 2
Time (s)
Time (s)
46
Caltech/UCSD Waveform Cross-Correlation Project

• 1984 to 2003 waveforms now online at Caltech
• Cross-correlation completed for 17 million event
  pairs
• Two relocated catalogs now available at SCEDC
• Hauksson et al. double-difference
• Shearer et al. cluster analysis (SHLK_1.0)

47
340,000 events in SHLK_1.0. Similar event
clusters, shown in black, are relocated using
cross-correlation data. Other events are colored
by year and are located using phase pick data
alone.
48
Catalog Locations
SHLK_1.0 locations
49
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

50
UCSD/Caltech spectral analysis
Online database of seismograms, 19842003 gt
300,000 earthquakes P and S multi-taper
spectra computed for all records 60 GB in
special binary format
51
Isolating Spectral Contributions
52
Receiver response
Observed spectrum
Source spectrum
Distance term to account for Q
logu(f)



log(f)

• gt 60,000 earthquakes, gt350 stations
• 1.38 million P-wave spectra (STN gt 5, 5-20 Hz)
• Iterative least squares approach with outlier
  suppression

53
Stacked P-wave source spectra

• Solid lines are stacks of source spectra in bins
  of equal moment, corrected for Q using empirical
  Greens function (EGF) method
• Dashed lines are best fit for w-2 source model
  with constant stress drop
• Ds 1.60 MPa (constant) for Madariaga (1976)
  source model

Mw 3.1
Relative Log Moment
Mw 1.9
10
20
50
5
2
Frequency (Hz)
54
Source-specific EGF method
For each event, find 500 neighboring events
Fit moment binned spectra to Ds and EGF
Then subtract EGF from target event spectrum and
compute Ds for this event
55
Observed source Ds using spatially varying EGF
method
65,070 quakes gt300,000 spectra
19892001 gt 4 spectra/event 5 - 20 Hz
band EGF corrected
56
Stress drop versus depth

• Average Ds increases from 0.6 to 2 MPa from 0 to
  8 km
• Nearly constant from 8 to 18 km
• Suggestion of falloff below 18 km
• Large scatter at all depths

90
Median
10
from Shearer et al. (2005)
57
Seismic data crunching examples

• Teleseismic
• AGC stacks
• Reference phase stacks
• Envelope function stacks for scattered energy
• Back-projection to image earthquake rupture
• Local
• Waveform cross-correlation for locations
• Spectra stacks for source properties

58
Migration in Reflection Seismology
For each pixel in image, sum values from each
trace at time of predicted source-to-scatterer-to-
receiver travel time
59
Back-projection to image earthquake rupture
Japanese Hi-Net array of 700 stations
2004 Sumatra-Andaman earthquake
from Ishii et al. (2005)
60
Direct Back-projection
Stack along predicted P-wave travel time curves
Assume grid of possible source locations
61
Problem Incoherent stacking from time shifts
from 3-D structure
Unmodeled 3-D velocity perturbations cause time
shifts in wavefront
62
Sumatra earthquake P-waves
Aligned on theoretical (iasp91) P-wave travel
times
63
Migration in Reflection Seismology
Problem Time shifts from 3-D structure can
destroy stack coherence Solution Statics
corrections (station terms)
slow
slow
fast

64
Align P-waves with cross-correlation
P onset
65
Method forces coherent stack at hypocenter
Cross-correlation times correct for
perturbations along each hypocenter-station ray
path
66
But coherence not guaranteed for sources offset
from hypocenter
Time shifts here not identical to hypocenter
shifts
Calibrated time corrections at hypocenter
67
Stacks at different source points
Time (seconds)
from Ishii et al. (2005)
68
Stacks and Time Slice (60 seconds)
Time (seconds)
69
Stacks and Time Slice (300 seconds)
Time (seconds)
70
(No Transcript)
71
Short-period radiation from Hi-net backprojection
(Ishii et al., 2005)
Harvard multiple CMT solution (Tsai et al., 2005)
72
2001 Mw 7.8 Kunlun earthquake
from Walker et al. (2005)
73
2005 M 7.6 Pakistan earthquake
1 2.5 s 2 7.5 s 3 12.5 s 4 17.5
s 5 22.5 s
from Walker et al. (2005)
74
P-wave Back-Projection Method

• Suited for a global near-real-time system
• Will give rapid information about rupture extent
  and likely areas of strong ground motion, which
  go beyond hypocenter and magnitude
• We are working with U.S. Geological Survey
  scientists to implement this method

75
THE DATA MUST BE GREATLY AMPLIFIED AND
STRENGTHENED.