next up previous [pdf]

Next: GPU implementation Up: Leader and Clapp: Accelerating Previous: GPU Implementation

Reverse Time Migration

Reverse time migration is a well known method, whereby reflection events are imaged in a kinematically accurate manner. Baysal and Sherwood (1983) and Claerbout (1985) pioneered the method, whereby an image is constructed by back propagating the receiver wavefield (recorded data) and forward propagating the source wavefield (modeled data, using assumed wavelet) and then cross correlating these at each imaging step. Mathematically this reduces to the imaging condition shown in Equation 2. It can be shown that this system is simply the adjoint of linearised modeling (Biondi, 2010).

$\displaystyle I(z,x,y) = \displaystyle\sum\limits_{i}^{nshots} \displaystyle\sum\limits_{t}^n P_s(t,z,x,y;\mathbf{s}_i) P_g(t,z,x,y;\mathbf{s}_i)$ (2)

$ I(x,y,z)$ is the image at point $ (x,y,z)$ , $ P_s(t,z,x,y;\bf s_i)$ is the source wavefield for shot $ \bf s_i$ and $ P_g(t,z,x,y;\bf s_i)$ is the corresponding recorded wavefield. This can be extended to a prestack imaging condition where subsurface offsets are included. This is necessary velocity estimation and many inversion methods. Equation 2 simply extracts the zero subsurface offset image.

$\displaystyle I(z,x,y,x_h,y_h) = \displaystyle\sum\limits_{i}^{nshots} \display...
...imits_{t}^n P_s(t,z,x+x_h,y+y_h;\mathbf{s}_i) P_g(t,z,x-x_h,y-y_h;\mathbf{s}_i)$ (3)

In Equation 3, $ x_h$ and $ y_h$ represent the subsurface half offsets. For each set of subsurface offsets the two wavefields are laterally shifted and the imaging condition reapplied, hence a 5D prestack image can be constructed (Biondi, 2006).

In 3D the computational intensity of RTM can become a hinderance relative to simpler migration operators. This is because the modeled shot, $ P_s(\bf s_i)$ , has to run forward in time, and the recorded data, $ P_g(\bf s_i)$ , is reversed in time. Since domain boundary damping makes the shot modeling non-reversable then practically we model and save the entirity of $ P_s(\bf s_i)$ , and hold this in memory as we reverse through $ P_g(\bf s_i)$ , applying the imaging condition at each time step. Of course, for large wavefields this is unpractical. One partial solution is referred in literature as checkpointing (Symes (2007); Dussaud and Cherrett (2008)). Here we effectively block the propagation and periodically reinject the source model. For example, let some intermediate position in the propagation be $ t_{check}$ , initially $ P_s(t_{max}-t_{check}; \bf s_i)$ is read, the receiver incrementally back propagated to $ t_{check}$ and the correlations applied. This can then be done $ t_{max}/t_{check}$ times (Clapp, 2008a) and can be parallelised. The memory problem is thus solved, however this process is fundamentally IO bound.

basic
basic
Figure 5.
A slice of the 3D wavefield after long propagation times, without and with random boundaries
[pdf] [png]

Clapp (2008b) proposed a different scheme - to use pseudo-random boundaries optimised to scatter these wavefields incoherently (Figure 5). The RTM imaging principle operates on the fact that incoherent data are either correlated out in imaging, or stacked out subsequently. Pseudo-random boundaries create a lot of noise in the propagated shot, however this will not correlate coherently with the receiver wavefield. Furthermore this propagation is now entirely time reversible, and by saving the velocity field the shot can be reverse propagated to its original state within machine precision (Figure 6). The benefits of this are twofold - computationally the source wavefield is time reversible, and so only two source wavefields are needed for any given time and no checkpointing is necessary; geophysically this algorithm provides accurate imaging, and any residual boundary scattered noise is well below the noise threshold. Thus a manner in which 3D RTM is both computationally and memory efficient can be constructed.

As seen in Figure 5, the low frequency component of the wavefield is not well scattered. For RTM this is not a significant problem, since low frequency artifacts are often prevalent and are removed by a low cut filter (or domain decomposition) before imaging (Taweesintananon (2011) provides an overview of these methods). However techniques such as waveform inversion methods rely on low frequency information. Shen and Clapp (2011) show how random boundaries can be further optimised, such that low and high frequencies become incoherent after scattering.

shotrand3
shotrand3
Figure 6.
An example of the time reversable wave propagation
[pdf] [png]



Subsections
next up previous [pdf]

Next: GPU implementation Up: Leader and Clapp: Accelerating Previous: GPU Implementation

2011-05-24