Numeric implementation of wave-equation migration velocity analysis operators |

Wavefield reconstruction for zero-offset migration based on the
one-way wave-equation is performed by recursive phase-shift in depth
starting with data recorded on the surface as boundary conditions. In
this configuration, the imaging condition extracts the image as time
from the reconstructed wavefield at every location in
space. Thus, the surface data need to be extrapolated backward in time
which is achieved by selecting the appropriate sign of the phase-shift
operation (which depends on the sign convention for temporal Fourier
transforms):

where represents the spatially-variable slowness at depth level . Equations 1-2 describe wavefield extrapolation using a pseudo-differential operator, which justifies our use of laterally-varying slowness . As indicated earlier, the image is obtained from this extrapolated wavefield by selection of time , which is typically implemented as summation of the extrapolated wavefield over frequencies:

Phase-shift extrapolation using wavenumbers computed using equations 1 and 2 is not feasible in media with lateral variation. Instead, implementation of such operators is done using approximations implemented in a mixed space-wavenumber domain (Huang et al., 1999; Stoffa et al., 1990; Ristow and Ruhl, 1994). A brief summary the mixed-domain implementation of the split-step Fourier (SSF) operator is presented in Appendix A.

For velocity analysis, we assume that we can separate the *total
slowness*
into a known background component
and an
unknown component
. With this convention, we can linearize
the depth wavenumber relative to the background slowness
using a truncated Taylor series expansion

(5) |

(6) |

A seismic image is produced using migration by wavefield extrapolation as follows: for each frequency, read data at all spatial locations ; then, for each depth, sum the wavefield into the image at that depth level (i.e. apply the imaging condition) and then reconstruct the wavefield to the next depth level (i.e. perform wavefield extrapolation). The ``-'' sign in this algorithm indicates that extrapolation is anti-causal (backward in time), and the factor ``2'' indicates that we are imaging data recorded in two-way traveltime with an algorithm designed under the exploding reflector model. Wavefield extrapolation is usually implemented in a mixed domain (space-wavenumber), as briefly summarized in Appendix A.

The wavefield perturbation
caused at depth by a
slowness perturbation
at depth is obtained by
subtraction of the wavefields extrapolated from to in the
true and background models:

Here, and correspond to a given depth level and frequency . A similar relation can be applied at all depths and all frequencies.

Equation 7 establishes a non-linear relation between
the wavefield perturbation
and the slowness perturbation
. Given the complexity and cost of numeric optimization based
on non-linear relations between model and wavefield parameters, it is
desirable to simplify this relation to a linear relation between model
and data parameters. Assuming small slowness perturbations, i.e. small
phase perturbations, the exponential function
can be linearized using the approximation
which is valid for small values of the phase
. Therefore the wavefield perturbation
at depth
can be written as

Equation 8 defines the zero-offset forward scattering operator , producing the scattered wavefield from the slowness perturbation , based on the background slowness and background wavefield at a given frequency . The image perturbation at depth is obtained from the scattered wavefield using the time imaging condition, similar to the procedure used for imaging in the background medium:

(9) |

Given an image perturbation
, we can construct the scattered
wavefield by the adjoint of the imaging condition

(10) |

Equation 11 defines the adjoint scattering operator , producing the slowness perturbation from the scattered wavefield , based on the background slowness and background wavefield . A typical implementation of zero-offset forward and adjoint scattering uses the following algorithms:

The forward zero-offset wave-equation MVA operator follows a similar pattern to the implementation of the downward continuation operator: for each frequency and for each depth, read the slowness perturbation at all spatial locations , then apply the scattering operator ( ) given equation 11 to the slowness perturbation and accumulate the scattered wavefield for downward continuation; then, apply the imaging condition ( ) producing the image perturbation at depth and then reconstruct the scattered wavefield backward in time using the wavefield extrapolation operator ( ) to the next depth level. The adjoint zero-offset wave-equation MVA operator also follows a similar pattern to the implementation of the downward continuation operator: for each frequency and for each depth, reconstruct the scattered wavefield forward in time using the wavefield extrapolation operator ( ) to the next depth level, then apply the adjoint of the imaging condition ( ) by adding the image to the scattered wavefield and then apply the adjoint wavefield scattering operator ( ) to obtain the slowness perturbation . Both forward and adjoint scattering algorithms require the background wavefield, , to be precomputed at all depth levels, although more efficient implementations using optimal checkpointing are possible (Symes, 2007). Scattering and wavefield extrapolation are implemented in the mixed space-wavenumber domain, as briefly explained in Appendix A.

Numeric implementation of wave-equation migration velocity analysis operators |

2013-08-29