Free-surface elastic wave injection and wavefield reconstruction with applications to elastic reverse time migration

Accurately injecting and extrapolating multicomponent seismic data is a major step in elastic reverse time migration (RTM) and full waveform inversion (FWI). The theoretical basis of this process is the seismic representation theorem. We propose a numerical approach based on the collocated-grid MacCormack finite-difference (FD) scheme to implement the injection and extrapolation. The MacCormack FD scheme defines wavefield parameters on the same grid. When dealing with the boundary condition in the representation theorem, it is more advantageous than the commonly used staggered-grid FD scheme which defines the stress and velocity on offset grids; thus interpolations are required to link the boundary values to nearby nodes. On the contrast, with the MacCormack scheme, the traction and particle velocity data can be directly injected to the boundary nodes as equivalent concentrated forces, force dipoles and double couples. The accuracy of this approach is validated by comparing the regenerated wavefield with the original wavefield. Our primary purpose is to develop an injection method for elastic RTM, but this approach can also be used to numerically investigate other interesting issues in the elastic wave injection, including errors and artifacts caused by injection from unclosed boundaries, incomplete boundary values, and uncertain velocity models. In particular, seismic data are usually acquired on the ground surface where the traction-free condition is naturally maintained. We simulate a land acquisition scenario by generating a synthetic multicomponent dataset on a traction-free surface. The reconstructed wavefields from velocity data are compared to the original forward propagated waves to evaluate the accuracy of the boundary treatments. Two multicomponent elastic RTM examples are presented to demonstrate the applications of the proposed injection method. By combining the regenerated receiver-side wavefield and the vectorized image condition, the PP, PS, SP and SS images are calculated.