Three-dimensional reverse-time migration of teleseismic receiver functions using the phase-shift-plus-interpolation method

Teleseismic receiver functions (RFs) are frequently used to determine depths of seismic discontinuities in the crust and upper mantle. We developed an efficient reverse-time migration (RTM) method that is applied to teleseismic receiver functions directly. Both the primary P-to-S converted phases and their crustal multiples in RFs can be used for imaging seismic discontinuities. The method uses the phase-shift-plus-interpolation algorithm to extrapolate both the source and receiver wavefields in a 3-D velocity model, which greatly reduces the computation costs compared with those using a full wave-equation numerical solver. Tests using synthetic data in various crustal models demonstrate the effectiveness of the method and its superiority over the common-conversion-point stacking method. In particular, the method handles diffraction caused by strong lateral structural variations correctly and there is no limitation on the maximum dip of the interface. We applied the method to real data of a linear array in the Wabash Valley Seismic Zone in the central USA and obtained a crustal structural image across a failed continental rift. We suggest that future passive-source seismic recording experiments for crustal scale imaging use station spacing less than 5 km, and a 2-D array with even smaller station spacing is desired for regions with strong lateral structural variations. With increasing numbers of sensors used in passive-source recording experiments nowadays, our RF-RTM method can be a useful tool for structural imaging on scales ranging from sedimentary basins, crust to lithosphere.