## Abstract

Traditional approaches of using dispersion curves for S-wave velocity reconstruction have limitations, principally, the 1D-layered model assumption and the automatic/manual picking of dispersion curves. At the same time, conventional full-waveform inversion (FWI) can easily converge to a non-global minimum when applied directly to complicated surface waves. Alternatively, the recently introduced wave equation dispersion spectrum inversion method can avoid these limitations, by applying the adjoint state method on the dispersion spectra of the observed and predicted data and utilizing the local similarity objective function to depress cycle skipping. We apply the wave equation dispersion spectrum inversion to three real datasets of different scales: tens of meters scale active-source data for estimating shallow targets, tens of kilometers scale ambient noise data for reservoir characterization and a continental-scale seismic array data for imaging the crust and uppermost mantle. We use these three open datasets from exploration to crustal scale seismology to demonstrate the effectiveness of the inversion method. The dispersion spectrum inversion method adapts well to the different-scale data without any special tuning. The main benefits of the proposed method over traditional methods are that (1) it can handle lateral variations; (2) it avoids direct picking dispersion curves; (3) it utilizes both the fundamental and higher modes of Rayleigh waves, and (4) the inversion can be solved using gradient-based local optimizations. Compared to the conventional 1D inversion, the dispersion spectrum inversion requires more computational cost since it requires solving the 2D/3D elastic wave equation in each iteration. A good match between the observed and predicted dispersion spectra also leads to a reasonably good match between the observed and predicted waveforms, though the inversion does not aim to match the waveforms.

This is a preview of subscription content, access via your institution.

## References

Beaty K, Schmitt D, Sacchi M (2002) Simulated annealing inversion of multimode Rayleigh wave dispersion curves for geological structure. Geophys J Int 151(2):622–631

Bensen GD, Ritzwoller MH, Barmin MP, Levshin AL, Lin F, Moschetti MP, Shapiro NM, Yang Y (2007) Processing seismic ambient noise data to obtain reliable broad-band surface wave dispersion measurements. Geophys J Int 169:1239–1260. https://doi.org/10.1111/j.1365-246X.2007.03374.x

Boore DM, Asten MW (2008) Comparisons of shear-wave slowness in the Santa Clara Valley, California, using blind interpretations of data from invasive and noninvasive methods. Bull Seismol Soc Am 98(4):1983–2003

Borisov D, Modrak R, Gao F, Tromp J (2018) 3D elastic full-waveform inversion of surface waves in the presence of irregular topography using an envelope-based misfit function. Geophysics 83(1):R1–R11

Cheng F, Xia J, Behm M, Hu Y, Pang J (2019) Automated data selection in the Tau-p domain: application to passive surface wave imaging. Surv Geophys 40(5):1211–1228

Choi Y, Alkhalifah T (2012) Application of multi-source waveform inversion to marine streamer data using the global correlation norm. Geophys Prospect 60(4):748–758

Çubuk-Sabuncu Y, Taymaz T, Fichtner A (2017) 3-D crustal velocity structure of western Turkey: constraints from full-waveform tomography. Phys Earth Planet Int 270:90–112

Dai T, Xia J, Ning L, Xi C, Liu Y, Xing H (2021) Deep learning for extracting dispersion curves. Surv Geophys 42:69–95. https://doi.org/10.1007/s10712-020-09615-3

Dziewonski AM, Anderson DL (1981) Preliminary reference Earth model. Phys Earth Planet Inter 25(4):297–356

Fang H, Yao H, Zhang H, Huang Y-C, van der Hilst RD (2015) Direct inversion of surface wave dispersion for three-dimensional shallow crustal structure based on ray tracing: methodology and application. Geophys J Int 201(3):1251–1263

Fomel S (2007) Local seismic attributes. Geophysics 72(3):A29–A33

Foti S, Hollender F, Garofalo F, Albarello D, Asten M, Bard P-Y et al (2018) Guidelines for the good practice of surface wave analysis: a product of the InterPACIFIC project. Bull Earthq Eng 16(6):2367–2420

Garofalo F, Foti S, Hollender F, Bard P, Cornou C, Cox B et al (2016a) InterPACIFIC project: comparison of invasive and non-invasive methods for seismic site characterization. Part II: inter-comparison between surface-wave and borehole methods. Soil Dyn Earthq Eng 82:241–254

Garofalo F, Foti S, Hollender F, Bard P, Cornou C, Cox BR et al (2016b) InterPACIFIC project: comparison of invasive and non-invasive methods for seismic site characterization. Part I: intra-comparison of surface wave methods. Soil Dyn Earthq Eng 82:222–240

Gu YJ, Webb SC, Lerner-Lam A, Gaherty JB (2005) Upper mantle structure beneath the eastern Pacific Ocean ridges. J Geophys Res: Solid Earth 110(B6)

IRIS D (2012) Data services products: Ancc-ciei, western us ambient noise cross-correlations

Ivanov J, Tsoflias G, Miller RD, Peterie S, Morton S, Xia J (2016) Impact of density information on Rayleigh surface wave inversion results. J Appl Geophys 135:43–54

Komatitsch D, Tromp J (2002) Spectral-element simulations of global seismic wave propagation—I. Validation. Geophys J Int 149(2):390–412

Lei W, Ruan Y, Bozdağ E, Peter D, Lefebvre M, Komatitsch D et al (2020) Global adjoint tomography-model GLAD-M25. Geophys J Int 223(1):1–21

Li J, Hanafy S, Schuster GT (2018). Dispersion inversion of guided P-waves in a waveguide of arbitrary geometry. In: SEG technical program expanded abstracts 2018. Society of Exploration Geophysicists, pp 2526–2530

Li Z, Chen X (2020) An effective method to extract overtones of surface wave from Array Seismic Records of earthquake events. J Geophys Res: Solid, Earth 125(3):e2019JB018511

Lin F-C, Moschetti MP, Ritzwoller MH (2008) Surface wave tomography of the western United States from ambient seismic noise: Rayleigh and Love wave phase velocity maps. Geophys J Int 173(1):281–298

Liu DC, Nocedal J (1989) On the limited memory BFGS method for large scale optimization. Math Program 45(1–3):503–528

Liu Q, Peter D, Tape C (2019) Square-root variable metric based elastic full-waveform inversion—Part 1: theory and validation. Geophys J Int 218(2):1121–1135

Luo Y, Xia J, Miller RD, Xu Y, Liu J, Liu Q (2008) Rayleigh-wave dispersive energy imaging using a high-resolution linear Radon transform. Pure Appl Geophys 165(5):903–922

Masoni I, Boelle J-L, Brossier R, Virieux J (2016). Layer stripping FWI for surface waves. In: SEG technical program expanded abstracts 2016. Society of Exploration Geophysicists, pp 1369–1373

Métivier L, Brossier R (2016) The SEISCOPE optimization toolbox: a large-scale nonlinear optimization library based on reverse communication. Geophysics 81(2):F1–F15

Mi B, Xia J, Shen C, Wang L (2018) Dispersion energy analysis of Rayleigh and love waves in the presence of low-velocity layers in near-surface seismic surveys. Surv Geophys 39(2):271–288

Mikesell TD, Malcolm AE, Yang D, Haney MM (2015) A comparison of methods to estimate seismic phase delays: numerical examples for coda wave interferometry. Geophys J Int 202(1):347–360

Nakata N, Snieder R, Tsuji T, Larner K, Matsuoka T (2011) Shear wave imaging from traffic noise using seismic interferometry by cross-coherence Shear wave imaging from traffic noise. Geophysics 76(6):SA97–SA106

Nocedal J, Wright S (2006) Numerical optimization. Springer

Pan Y, Gao L, Bohlen T (2019). From multichannel analysis to waveform inversion of shallow-seismic surface waves. In: SEG technical program expanded abstracts 2019. Society of Exploration Geophysicists, pp 5020–5024

Park CB, Miller RD, Xia J (1999) Multichannel analysis of surface waves. Geophysics 64(3):800–808

Pica A, Diet J, Tarantola A (1990) Nonlinear inversion of seismic reflection data in a laterally invariant medium. Geophysics 55(3):284–292

Plessix R-E (2006) A review of the adjoint-state method for computing the gradient of a functional with geophysical applications. Geophys J Int 167(2):495–503

Ren L, McMechan GA, Guo P (2019) Concurrent elastic inversion of Rayleigh and body waves with interleaved envelope-based and waveform-based misfit functions. In: SEG international exposition and annual meeting

Rix GJ, Leipski EA (1995) Accuracy and resolutions of surface wave inversion. In Recent advances in instrumentation. ASCE, Data Acquisition and Testing in Soil Dynamics, pp. 17–32

Savage M, Sheehan A (2000) Seismic anisotropy and mantle flow from the Great Basin to the Great Plains, western United States. J Geophys Res: Solid Earth 105(B6):13715–13734

Saygin E, Kennett BLN (2012) Crustal structure of Australia from ambient seismic noise tomography. J Geophys Res 117:B01304. https://doi.org/10.1029/2011JB008403

Shapiro NM, Campillo M, Stehly L, Ritzwoller MH (2005) High-resolution surface-wave tomography from ambient seismic noise. Science 307(5715):1615–1618

Shen W, Ritzwoller MH (2016) Crustal and uppermost mantle structure beneath the United States. J Geophys Res: Solid Earth 121(6):4306–4342

Shragge J, Yang J, Issa N, Roelens M, Dentith M, Schediwy S (2021) Low-frequency ambient distributed acoustic sensing (DAS): case study from Perth, Australia. Geophys J Int 226(1):564–581

Shragge J, Yang J, Issa NA, Roelens M, Dentith M, Schediwy S (2019). Low-frequency ambient Distributed Acoustic Sensing (DAS): useful for subsurface investigation? In: SEG technical program expanded abstracts 2019. Society of Exploration Geophysicists, pp 963–967

Smith JA, Borisov D, Cudney H, Miller RD, Modrak R, Moran M, et al. (2019) Tunnel detection at Yuma Proving Ground, Arizona, USA—Part 2: 3D full-waveform inversion experiment stunnel detection at ypg-part 2: 3d fwi. Geophysics 84(1):B95–B108

Snieder R (2004) Extracting the Green’s function from the correlation of coda waves: a derivation based on stationary phase. Phys Rev E 69(4):046610

Socco LV, Foti S, Boiero D (2010) Surface-wave analysis for building near-surface velocity models-established approaches and new perspectives. Geophysics 75(5):75A83–75A102

Solano CP, Donno D, Chauris H (2013) 2D surface wave inversion in the FK domain. In: 75th eage conference and exhibition incorporating SPE Europec

Strobbia C, Andreas L, Vermeer P, Glushchenko A (2011) Surface waves: use them then lose them. Surface-wave analysis, inversion and attenuation in land reflection seismic surveying. Near Surface. Geophysics 9(6):503–513

Tanimoto T (2008) Normal-mode solution for the seismic noise cross-correlation method. Geophys J Int 175(3):1169–1175

Tape C, Liu Q, Maggi A, Tromp J (2010) Seismic tomography of the southern California crust based on spectral-element and adjoint methods. Geophys J Int 180(1):433–462

Tromp J (2020) Seismic wavefield imaging of earth’s interior across scales. Nat Rev Earth Environ 1(1):40–53

Virieux J, Operto S (2009) An overview of full-waveform inversion in exploration geophysics. Geophysics 74(6):WCC1–WCC26

Wang J, Gu YJ, Chen Y (2020) Shear velocity and radial anisotropy beneath southwestern Canada: evidence for crustal extension and thick-skinned tectonics. J Geophys Res: Solid, Earth 125(2):e2019JB018310

Wang Y, Miller RD, Peterie SL, Sloan SD, Moran ML, Cudney HH, et al. (2019) Tunnel detection at Yuma proving ground, Arizona, USA—Part 1: 2D full-waveform inversion experiment. Geophysics 84(1):B95–B105

Wolfe P (1969) Convergence conditions for ascent methods. SIAM Rev 11(2):226–235

Xia J, Miller RD, Park CB (1999) Estimation of near-surface shear-wave velocity by inversion of Rayleigh waves. Geophysics 64(3):691–700

Xia J, Miller RD, Park CB, Tian G (2003) Inversion of high frequency surface waves with fundamental and higher modes. J Appl Geophys 52(1):45–57

Yi J, Liu Y, Yang Z, Lu H, He B, Zhang Z (2019) A least-squares correlation-based full traveltime inversion for shallow subsurface velocity reconstruction. Geophysics 84(4):R613–R624

Yoon K, Suh S, Cai J, Wang B (2012) Improvements in time domain FWI and its applications. In: 2012 SEG annual meeting

Yuan YO, Simons FJ, Bozdağ E (2015) Multiscale adjoint waveform tomography for surface and body waves. Geophysics 80(5):R281–R302

Zhang Z-d, Alajami M, Alkhalifah T (2020) Wave-equation dispersion spectrum inversion for near-surface characterization using fibre-optics acquisition. Geophys J Int 222(2):907–918

Zhang Z-d, Alkhalifah T (2019a) Local-crosscorrelation elastic full-waveform inversion. Geophysics 84(6):R897–R908

Zhang Z-D, Alkhalifah T (2019b) Wave-equation Rayleigh-wave dispersion inversion using fundamental and higher modes. Geophysics 84(4):EN57–EN65

Zhang Z-D, Liu Y, Schuster G (2015). Wave equation inversion of skeletonized surface waves. In: SEG technical program expanded abstracts 2015. Society of Exploration Geophysicists, pp 2391–2395

Zhang Z-d, Schuster G, Liu Y, Hanafy SM, Li J (2016) Wave equation dispersion inversion using a difference approximation to the dispersion-curve misfit gradient. J Appl Geophys 133:9–15

Zheng Y, Hu H (2017) Nonlinear signal comparison and high-resolution measurement of surface-wave dispersion. Bull Seismol Soc Am 107(3):1551–1556

Zhu H, Bozdağ E, Peter D, Tromp J (2012) Structure of the European upper mantle revealed by adjoint tomography. Nat Geosci 5(7):493–498

Zigone D, Ben-Zion Y, Campillo M, Roux P (2015) Seismic tomography of the southern California plate boundary region from noise-based Rayleigh and love waves. Pure Appl Geophys 172(5):1007–1032

## Acknowledgements

We appreciate Caroline Johnson’s help in improving the manuscript. E. Saygin and L. He wish to acknowledge financial assistance provided through Australian National Low Emissions Coal Research and Development (ANLEC R&D). E. Saygin was supported by CSIRO’s Deep Earth Imaging Future Science Platform. L. He was supported by China Scholarship Council. ANLEC R&D is supported by Australian Coal Association Low Emissions Technology Limited and the Australian Government through the Clean Energy Initiative. We thank the InterPACIFIC project for providing the field data (data available here: http://interpacific.geopsy.org/). We thank David Lumley for his contribution to the retrieval of the continuous part of the SW HUB dataset (data available here: https://wapims.dmp.wa.gov.au/WAPIMS/Search/SwHubCarbonStorage). We thank Weisen Shen for providing the reference S-wave velocity model (data available here: http://ciei.colorado.edu/Models/). We extend our thanks to Roman Pevzner for providing the reference S-wave velocity of Harvey-1 well. Z. Zhang and T. Alkhalifah thank KAUST for its support and specifically the seismic wave analysis group members for their valuable insights. For computer time, this research used the resources of the Supercomputing Laboratory at King Abdullah University of Science and Technology (KAUST) in Thuwal, Saudi Arabia.

## Author information

### Affiliations

### Corresponding author

## Additional information

### Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

## Appendix A: The wave equation dispersion spectrum inversion

### Appendix A: The wave equation dispersion spectrum inversion

In full-waveform inversion, objective functions intend to measure the mismatch between the predicted and observed data. One of the most intuitive measurements is the \(L_2\) norm distance, which is given by

where \(d^{p}(s,r,t,m)\) and \(d^{o}(s,r,t)\) are predicted and observed data, respectively. *m* denotes the model parameter, which is, in this paper, shear wave velocity. \(\phi (m)\) measures the waveform difference. *s*, *r* and *t* are indexes of time-domain common-shot-gathers, which are for source, receiver and time, respectively.

The widely used \(L_2\) norm objective function suffers from the known cycle-skipping problem, which is even worse for dispersive surface waves. Inspired by the success of surface wave dispersion curve inversions (e.g., MASW), we found that the dispersion spectrum (e.g., f-v spectrum) is less influenced by the cycle-skipping than waveform. The shape of the dispersion spectrum retains the kinematic information needed for S-wave recovery in the near-surface, while the amplitude is influenced by many factors, making it difficult to predict by numerical solutions. As a result, we propose to measure the curvature difference between observed and predicted dispersion spectra. In particular, we use the local crosscorrelation to quantify the difference between the predicted and observed dispersion spectra. The objective function is given by Zhang and Alkhalifah (2019b)

with

where \(S(s,f,v,f',m)\) is a local crosscorrelation function, which quantifies the local similarity (range 0 and 1) between the predicted and observed f-v spectra. Note that, the predicted and observed data are dispersion spectra in frequency-velocity domain. *f* and *v* are the available frequency and selected phase-velocity. \(|C^p(s,f,v,m)|\) and \(|C^o(s,f,v,f')|\) are the predicted and observed Rayleigh wave dispersion spectra (the amplitude spectrum), respectively. The calculation of such spectra will be explained later. \(f'\) denotes the frequency extension of the observed spectrum. Here, we use trace stretching to calculate the phase shift (Mikesell et al. 2015). A typical extension range is \(f'=0.8-1.2\). *F* and *w* denote the central point and the length of the sliding window, respectively. In practice, we use the shaping regularization introduced by Fomel (2007) to perform the crosscorrelation of the predicted and observed dispersion spectra locally. Unlike the global correlation, the local measurements can suppress cross talk between neighboring modes (the fundamental and higher modes). The penalty function, \(W(f')\), forces the predicted spectrum to be similar to the observed one. Instead of the widely used linear and Gaussian penalty functions, we use a polynomial function to satisfy the following boundary conditions: \(W|_{ f'_{\min }} = 0; W|_{ f'_{\max }} = 0; W'|_{ f'_{\min }} = 0; W'|_{ f'_{\max }} = 0 ; W|_1 = 1; W'|_1=0\). The proposed penalty function is decided by the selected extension range (\(f'_{\min }\), \(f'_{\max }\)), thus there are fewer parameters to manipulate. More importantly, it reduces to zero smoothly at the extension boundaries and thus eliminating discontinuities from the objective function.

The f-v spectrum is calculated using a linear Radon transform (Luo et al. 2008). After a temporal Fourier transform of the shot gather, the linear Radon transform can be calculated for each temporal frequency component *f* as:

and its adjoint form is given by

*D*(*f*, *x*) denotes the Fourier transform of the common-shot-gathers in the time axis. *C*(*f*, *v*) is the complex f-v spectrum. Notice that, we use the amplitude spectrum (|*C*(*f*, *v*)|) in the objective function [Eq. (A.3)] and the corresponding adjoint source [Eq. (A.7)]. \(x_{\min }\) and \(x_{\max }\) are the minimum and the maximum offsets used for the integration, respectively. \(v_{min}\) and \(v_{max}\) are the boundaries of the selected phase velocities for integral, which control the stacking slopes of Rayleigh waves at different frequencies. We can set a physical range for them that is dependent on the target area.

The predicted seismic data are calculated by solving a first order elastic wave equation, which is given by

where \(\Psi = (v_1, v_2, v_3, \sigma _1, \sigma _2, \sigma _3, \sigma _4, \sigma _5, \sigma _6)^T\) is a vector containing three particle velocities and six stresses components, *E* denotes the spatial-differentiation operator, \(\mathbf{I} _3\) is a 3 by 3 identity matrix. \(\mathbf{C}\) represents the stiffness matrix and \(\mathbf{s}\) denotes the point source used for modeling.

The gradient for S-wave velocity updating with respect to the objective function can be calculated by the chain rule, which is \(\frac{\partial \phi (m)}{\partial m}=\frac{\partial \phi (m)}{\partial d^p(m)}\frac{\partial d^p(m)}{\partial m}\), where \(\frac{\partial \phi (m)}{\partial d^p(m)}\) is the adjoint source and \(\frac{\partial d^p(m)}{\partial m}\) is the Fréchet derivative (Plessix 2006). The adjoint source of the proposed objective function is given by Zhang and Alkhalifah (2019b)

and

By crosscorrelating of the forward-propagated source wavefield and the backward-propagated adjoint-source [Eqs. (A.7) and (A.8)] wavefield, we can calculate the gradient of the objective function. S-wave velocity is iteratively updated using the Limited-memory Broyden–Fletcher–Goldfarb–Shanno (l-BFGS) optimization scheme (Liu and Nocedal 1989), which can be written as

where \(H^{-1}\) denotes the inverse Hessian which is approximated by the l-BFGS and \(\alpha\) is the step length used for S-wave velocity updates, which is determined by Wolfe conditions (Wolfe 1969):

with \(p_k =-g_k\) for gradient descent and \(0<c_1<c_2<1\).

The Wolfe conditions have been successfully used in solving waveform inversion problems (Métivier and Brossier 2016; Liu et al. 2019). An acceptable step length should decrease the data misfit (the Armijo rule) and the slope (the curvature condition) sufficiently. Nocedal and Wright (2006) suggested \(c_1=10^{-4}\) and \(c_2=0.9\) to obtain sufficient reductions for quasi-Newton methods. The iteration is terminated when there are ten fails in searching for \(\alpha\).

## Rights and permissions

## About this article

### Cite this article

Zhang, Zd., Saygin, E., He, L. *et al.* Rayleigh Wave Dispersion Spectrum Inversion Across Scales.
*Surv Geophys* **42, **1281–1303 (2021). https://doi.org/10.1007/s10712-021-09667-z

Received:

Accepted:

Published:

Issue Date:

### Keywords

- Rayleigh waves
- Near-surface
- Elastic inversion
- High modes