Topics in Dip-Moveout Correction and Prestack Time Migration

From SEG Wiki
Jump to: navigation, search
Seismic Data Analysis
Seismic-data-analysis.jpg
Series Investigations in Geophysics
Author Öz Yilmaz
DOI http://dx.doi.org/10.1190/1.9781560801580
ISBN ISBN 978-1-56080-094-1
Store SEG Online Store


E.1 Reflection point dispersal

We start the analysis using the geometry of Figure E-1. Place the source S at the origin (0, 0), and denote the coordinates of the points of interest as N : (x0 = SK, z0 = KN), H : (x1 = SB, z1 = BH), F : (x2 = SA, z2 = AF), and R : (x3 = SQ, z3 = QR). Our objective is to compute the coordinates for the normal-incidence reflection point N : (x0, z0) associated with midpoint M, and the reflection point R : (x3, z3) associated with the nonzero-offset source receiver separation x = SG.

First, we compute the coordinates for N : (x0, z0). From the geometry of Figure E-1, first, note that x0 = SK and z0 = KN. Then, from the triangle MKN, we have

which, by way of MK = x0x/2, yields


(E-1a)

Similarly, from the same triangle MKN, we have

which, by way of KN = z0, yields


(E-1b)

Finally, combine equations (E-1a) and (E-1b) to get a relation between x0 and z0


(E-2a)

Next, write the equation of the line that represents the dipping reflector in terms of the directional cosines and the normal d = SC to the line from the origin S


(E-2b)

Solve for z0 from equation (E-2a) and substitute into equation (E-2b)


(E-3)

Note that the directional cosines satisfy the relation cos2 α + cos2 ϕ = 1. By making use of this relation, simplify equation (E-3) to get the final expression for the coordinate x0


(E-4a)

Finally, substitute x0 back into equation (E-2a) to get the expression for the coordinate z0


(E-4b)
Figure E-1  Geometry of a dipping reflector to derive the expression for reflection-point dispersal.

The distance d = SC is given by

which can be rewritten using the geometry of Figure E-1 as follows:

Substitute the definitions SM = x/2 and D = MN, where D is the distance along the normal-incident raypath from the point N on the reflector to the midpoint location M, to obtain the relation between d and D


(E-5)

Make the further substitution of equation (E-5) into equations (E-4a) and (E-4b), and recall that cos2 α + cos2 ϕ = 1 to derive the expressions for the coordinates of N : (x0, z0) in terms of D


(E-6a)

and


(E-6b)

To derive the expressions for the coordinates of R : (x3, z3), we note that R is the intersection of the line defined by the reflecting interface, which yields the relation


(E-7a)

and the line SRF, which yields the relation


(E-7b)

where (x2, z2) are the coordinates of F. From the geometry of Figure E-1, first, note that x2 = SA, z2 = AF, and x3 = SQ and z3 = QR.

Solve equation (E-7b) for x3 and substitute into equation (E-7a) to get


(E-8)

The coordinates for F : (x2, z2) are derived in Section C.3 (equations C-18c,d)


(E-9a)

and


(E-9b)

Substitute equations (E-9a) and (E-9b) into equation (E-8), and replace d via equation (E-5) to obtain, after some involved algebra, the expression for z3


(E-10a)

Back substitution of equation (E-10a) into equation (E-7a) then yields the expression for x3


(E-10b)

Knowing the coordinates of N : (x0, z0) and R : (x3, z3), the distance Δ = NR can now be computed


(E-11a)

Substitute equations (E-10a) and (E-10b) for the coordinates of R : (x3, z3), and equations (E-6a) and (E-6b) for the coordinates of N : (x0, z0) into equation (E-11a) and carry out the algebraic steps to get the final expression


(E-11b)

Equation (E-11b) represents the reflection point smear, otherwise known as reflection point dispersal [1] — the distance between the normal-incidence reflection point N associated with the zero-offset raypath from midpoint M, and the reflection point R associated with the nonzero-offset raypath for a source-receiver pair separated by an offset x. For the case of zero-offset x = 0, reflection point smear Δ vanishes.

Define the two-way zero-offset time t0 at midpoint M as the time associated with the normal-incidence raypath MN by the relation vt0/2 = D. Then use this relation in equation (E-11b) to obtain the expression for reflection point dispersal Δ in terms of t0


(E-12)
Figure 3.1-14  Moveout for low-velocity event (a) is larger than for high-velocity event (b). Moveout for low-velocity dipping event (c) may not be distinguishable from high-velocity horizontal event (b). These observations are direct consequences of equation (7).

This equation is applicable to the general case of three-dimensional geometry (Figure 3.1-14), where α is the angle between the normal to the dipping reflector surface and the direction of the profile line [2].

For the 2-D geometry of the dipping reflector Figure E-1, note that


(E-13)

where ϕ is the dip angle of the reflector. Hence, equation (E-12) is written in terms of the reflector dip angle ϕ as follows:


(6)
(10)
(E-14)

Equation (E-14), by way of equation (6), is the same as equation (10) of the main text with offset x = 2h.

E.2 Equations for DMO correction

Figure 5.1-1  (a) The geometry of a nonzero-offset recording of reflections from a dipping layer boundary; (b) a sketch of the time section depicting the various traveltimes. NMO correction involves coordinate transformation from yn − t to yn − tn by mapping amplitude A at time t to B at time tn on the same trace. DMO correction involves coordinate transformation from yn − tn to y0τ0 by mapping amplitude B at time tn on the trace at midpoint location yn of the moveout-corrected common-offset section to amplitude C at time τ0 on the trace at midpoint location y0 of the zero-offset section. Zero-offset migration involves coordinate transformation from y0τ0 to ym − τ by mapping amplitude C at time τ0 on the trace at midpoint location y0 of the zero-offset section to amplitude D at time τ on the trace at midpoint location ym of the migrated section. Migration before stack involves direct mapping of amplitude A at time t on the trace at midpoint location yn of the common-offset section to amplitude D at time τ on the trace at midpoint location ym of the migrated section. See text for the relationships between the coordinate variables.

In this section, we shall derive the traveltime equation for dip-moveout (DMO) correction by using the geometry of Figure 5.1-1. Start with equation (C-22b) of Section C.3 which defines the traveltime t from source location S to the reflection point R to the receiver location G, written in prestack data coordinates


(E-15a)

where 2h is the offset, v is the medium velocity above the reflector, ϕ is the reflector dip, and t0 is the two-way zero-offset time at midpoint location yn. Dip-moveout correction is preceded by normal-moveout (NMO) correction using the dip-independent velocity v


(E-15b)

where tn is the event time after the NMO correction. To relate event time tn after the NMO correction and event time t0, first, write equation (E-15a) as


(E-16)

Set the right-hand sides of equations (E-15b) and (E-16) equal and simplify the result to obtain the dip-dependent moveout equation


(E-17)

Equations (E-15b) and (E-17) suggest that moveout correction can, in principle, be performed in two steps:

  1. Apply a dip-independent moveout correction using equation (E-15b) to map event time t to event time tn.
  2. Apply a dip-dependent moveout correction using equation (E-17) to map event time tn to event time t0.

This two-step moveout correction is equivalent to the one-step moveout correction using equation (E-15a) to map event time t directly to event time t0.

Our goal, however, is to map event time t to τ0 — two-way zero-offset time, not at midpoint location yn associated with the source-receiver pair S − G, but at midpoint location y0 associated with the normal-incidence reflection point R. This mapping requires coordinate transformation of moveout-corrected prestack data Pn(yn, tn; h) to P0(y0, τ0; h). Therefore, we need to compute the coordinates (y0, τ0) in terms of the coordinates (yn, tn).

From the geometry of Figure 5.1-1, note that


(E-18)

where Δ = NR is the distance along the reflector between the normal-incidence reflection point N associated with the zero-offset raypath from midpoint yn and the reflection point R associated with the nonzero-offset raypath for a source-receiver pair separated by an offset 2h.

Adapt equation (E-14) of Section E.1 for the prestack data coordinates x = 2h


(E-19)

and substitute into equation (E-18) to obtain


(E-20)

Solve equation (E-17) for t0 and write the result in the following form


(E-21)

or


(E-22)

where


(E-23)

Now, substitute equation (E-22) into equation (E-20) to get the final expression for y0


(E-24)

Also from the geometry of Figure 5.1-1, note that


(E-25)

Substitute equation (E-19) for Δ and carry out the required algebra to get the desired expression for τ0 in terms of t0


(E-26)

Now, substitute equation (E-22) into equation (E-26)


(E-27)

Finally, simplify equation (E-27) by way of equation (E-23) to obtain the desired expression for τ0 in terms of tn


(E-28)

We now remind ourselves of our objective to transform the normal-moveout-corrected prestack data Pn(yn, tn; h) from yn − tn coordinates to y0τ0 coordinates so as to obtain the dip-moveout-corrected data P0(y0, τ0; h). The transformation is done using equations (E-24) and (E-28) which relate the input data coordinates yn − tn to output data coordinates y0τ0. Note, however, equations (E-24) and (E-28) require knowledge of the reflector dip ϕ. To circumvent this requirement, we use the relation from Section D.1


(E-29)

which states that the reflector dip ϕ can be expressed in terms of wavenumber ky and frequency ω0, which are the Fourier duals of midpoint y0 and event time τ0, respectively. The variables y0 and τ0 defined by equations (E-24) and (E-28), respectively, are associated with the zero-offset section after DMO correction. By way of equation (E-29), these equations are recast explicitly independent of reflector dip as


(E-30)

and


(E-31)

where A of equation (E-23), by way of equation (E-29), is of the form


(E-32)

Since we have switched to the Fourier transform domain in our analysis, our objective now is to compute the dip-moveout-corrected data P0(ky, ω0; h) in the transform domain. Thus, start with the 2-D Fourier transform of P0(y0, τ0; h)


(E-33)

A wavefield is invariant under a coordinate transformation; hence, P0(y0, ω0; h) = Pn(yn, tn; h). Use the transform relations given by equations (E-30) and (E-31) and the invariance relation to write equation (E-33) in terms of the normal-moveout-corrected data Pn(yn, tn; h) [3]


(E-34)

From equation (E-30) it follows that


(E-35)

Square both sides of equation (E-31), substitute equation (E-32) for the quantity A and simplify to get


(E-36a)

Next, differentiate the result


(E-36b)

Rearrange the terms in equation (E-32) to the form


(E-36c)

and use in equation (E-36b) together with equation (E-31) for the τ0 in the denominator to obtain


(E-37)

Finally, substitute equations (E-35) and (E-37) into equation (E-34)


(E-38)

Now, we turn our attention to the phase term in equation (E-38)

Rearrange the terms, first, as follows:

then, by way of the expression for A as in equation (E-32)


(E-39)

Return to equation (E-38) and substitute equation (E-39) for the terms in the exponential


(E-40)

Use the Fourier transform integral


(E-41)

to rewrite equation (E-40) in the form


(E-42)

Equation (E-42) describes the dip-moveout correction process which transforms the normal-moveout-corrected prestack data with a specific offset 2h from yn − tn domain to y0τ0 domain. Referring to Figure 5.1-1, we remind ourselves that yn − tn coordinates are associated with event time tn at midpoint location yn, and y0τ0 coordinates are associated with event time τ0 at midpoint location y0 that corresponds to the normal-incidence reflection. The lateral excursion applied by dip-moveout correction is ΔyDMO = |yn − y0| and the vertical excursion is ΔtDMO = tn − τ0.

Figure 5.1-2  A flowchart for frequency-wavenumber dip-moveout correction algorithm. The scalar A is given by equation (13) and B = (2A2 − 1)/A3 as in equation (14a).

Once dip-moveout correction is applied, the data are inverse Fourier transformed


(E-43)

A flowchart of the dip-moveout correction in the frequency-wavenumber domain described above is presented in Figure 5.1-2.

E.3 Log-stretch DMO correction

A computationally efficient DMO correction can be formulated in the logarithmic time domain [4][5][6][3] [7]. The log-stretch time variable enables linearization of the coordinate transform equation (E-28), and as a result, the DMO correction is achieved by a simple multiplication of the input data with a phase-shift operator in the Fourier transform domain.

Define the following logarithmic variables that correspond to the time variables τ0 and tn of equation (E-28):


(E-44a)

and


(E-44b)

Hence, the inverse relationships are given by


(E-45a)

and


(E-45b)

In Section E.2, we derived the equations for DMO correction to transform the normal-moveout-corrected data from (yn, tn) coordinates to (y0, τ0) coordinates. The objective here is to derive equations for DMO correction in the log-stretch coordinates (y0, T0).

Square both sides of equation (E-28)


(E-46a)

and substitute equation (E-23) for A


(E-46b)

Here, we have to make the crucial assumption that A is close to unity. Refer to equation (E-23) and note that such an assumption implies one of the following: small h (offset), small ϕ (dip), large tn (time after moveout) or large v (velocity). As a result, equation (E-46b) can be approximated as


(E-46c)

The slope 0/dy0 measured on a zero-offset section in (y0, τ0) coordinates is given by


(E-46d)

Substitute equation (E-46d) into equation (E-46c)


(E-46e)

Apply the chain rule for differentiation


(E-47a)

and, by way of equation (E-45a), obtain the relationship


(E-47b)

Next, combine equations (E-29) and (E-46d) to get


(E-47c)

By analogy, we may write


(E-47d)

where Ω0 is the Fourier transform variable associated with the log-stretch time variable T0. Finally, substitute equation (E-47d) into equation (E-47b) to get


(E-47e)

Now, return to equation (E-46e), and use equations (E-45a,E-45b) and (E-47e) to obtain the expression


(E-48)

Take the logarithm of both sides and simplify to get the transform relation between the input log-stretch time variable Tn and the output log-stretch time variable T0


(E-49)

where


(E-50)

We now turn our attention to equation (E-24) and derive the corresponding equation in the log-stretch domain. First, substitute equation (E-23) into equation (E-24) and approximate the quantity A as before


(E-51a)

Substitute equation (E-46d) into equation (E-51a)


(E-51b)

Then, make the substitutions from equations (E-45b) and (E-47e), drop the high-order term, and simplify to get


(E-51c)

Note from equation (E-48)


(E-52)

Use this result in equation (E-51c) accompanied with the definition of equation (E-50) to obtain the final expression that relates y0 and yn in log-stretch coordinates


(E-53)

Equations (E-49), (E-50), and (E-53) correspond to equations (E-31), (E-32), and (E-30) in the log-stretch domain.

We remind ourselves of the objective to compute the dip-moveout-corrected data P0(ky, Ω0; h) in the transform domain. Thus, start with the 2-D Fourier transform of P0(y0, T0; h)


(E-54)

A wavefield is invariant under a coordinate transformation; hence, P0(y0, T0; h) = Pn(yn, Tn; h), the latter being the normal-moveout-corrected data in the log-stretch domain. Use the transform relations given by equations (E-49) and (E-53) and the invariance relation to write equation (E-54) in terms of the normal-moveout-corrected data Pn(yn, Tn; h) [3]


(E-55)

From equation (E-53) it follows that


(E-56a)

and from equation (E-49), we have


(E-56b)

Substitute equations (E-56a,E-56b) into equation (E-55) and rearrange the terms in the exponential to obtain


(E-57)

Assume that the Fourier transform variable Ω0 in the log-stretch domain is independent of Tn. Then, the double integral on the right-hand side represents the 2-D Fourier transform of the normal-moveout-corrected data in the log-stretch domain as


(E-58)

Therefore, by way of equation (E-58), equation (E-57) takes the form


(E-59)

Equation (E-59) describes the dip-moveout correction process which transforms the normal-moveout-corrected prestack data with a specific offset 2h from the log-stretch yn − Tn domain to y0T0 domain. Note that the relationship of input Pn(ky, Ω0; h) to output P0(ky, Ω0; h) given by equation (E-59) computationally is much simpler than that of equation (E-42). Once the phase-shift in the log-stretch domain given by equation (E-59) is applied to the input and the result is inverse Fourier transformed, the data are unstretched from (y0, T0) coordinates back to (y0, τ0) coordinates using the relationship expressed by equation (E-45a).

A variation of the phase-shift term in equation (E-59) is given by Notfors and Godfrey [6]. As in most log-stretch formulations of DMO correction, this reference assumes that under DMO correction the midpoint variable is invariant; hence, by way of equation (E-53), the first term in the exponential of equation (E-59) drops out. A further approximation, ln Ae = Ae − 1, and use of the definition for Ae given by equation (E-50) then lead to the following expression for DMO correction:


(E-60)

Various implementations of the log-stretch technique applied to shot records are described by Biondi and Ronen [8], Cabrera and Levy [9], and Zhou [7]. The latter also includes a method in double-log-stretch domain, which involves stretching both in time and midpoint coordinates. See principles of dip-moveout correction for a practical implementation of the log-stretch DMO correction described here.

E.4 The DMO ellipse

The objective is to derive the traveltime trajectory associated with the dip-moveout correction operator by using the method of stationary phase. Insert from equation (E-42) the expression for P0(ky, τ0; h) into equation (E-43)


(E-61)

where the total phase is given by


(E-62a)

Substitute equation (E-32) for A


(E-62b)

The main contribution to the integration in equation (E-61) occurs when the phase in equation (E-62b) stays nearly constant. We therefore determine the variation of Φ with respect to variables ω0 and ky


(E-63a)

and


(E-63b)

Then set each variation to zero. Rearranging the terms of the resulting expressions, we have


(E-64a)

and


(E-64b)

Then, square both sides of equations (E-64a,E-64b)


(E-65a)

and


(E-65b)

Finally, sum equations (E-65a,E-65b) to obtain


(E-66)

Equation (E-66) describes an ellipse with the following properties:

  1. Semi-major axis in midpoint y0 direction: a = h.
  2. Semi-minor axis in time τ0 direction: b = tn.

The ellipse of equation (E-66) in the y0τ0 plane describes the impulse response of a dip-moveout operator applied to a common-offset section with offset 2h.

E.5 Nonzero-offset traveltime equation

Prestack wave extrapolation is performed using the double-square operator which enables downward continuation of common-shot and common-receiver gathers (Section D.1). Stationary phase approximation to the double-square root operator yields the nonzero-offset traveltime equation (D-35) derived in Section D.2


(E-67)

where y, h, t, z are midpoint, offset, two-way traveltime and depth coordinates and v is the medium velocity. Here, we shall show that this equation represents an ellipse (Figure E-2) in the y − z plane for a constant t, h, and v, and derive the parameters of this ellipse. The nonzero-offset two-way time is associated with the raypath from source S to reflection point R to receiver G as sketched in Figure E-2. The origin of the y − z plane coincides with midpoint M.

Figure E-2  The prestack time migration ellipse. See Section E.5 for details.

Square both sides of equation (E-67) to get


(E-68)

Combine the second and third terms on the right-hand side and simplify the terms inside the square root


(E-69a)

Perform further algebraic manipulation to collect the terms in y and z


(E-69b)

Finally, normalize by the terms on the right-hand side and rearrange the terms in the denominators


(E-70)

Equation (E-70) represents an ellipse in the y − z plane for a constant t with the following parameters (Figure E-2):

  1. Semi-major axis in midpoint y direction: a = vt/2.
  2. Semi-minor axis in depth z direction:
  3. Distance from center to either focus:
  4. Distance from one focus to a point on the ellipse to the other focus: vt.

The ellipse of equation (E-70) in the y − z plane describes the impulse response of a nonzero-offset migration operator applied to prestack data.

When equation (E-70) is specialized to the zero-offset case, h = 0, we get


(E-71a)

which describes a circle in the y − z plane for a constant t with a radius vt/2. This circle represents the impulse response of a zero-offset migration operator applied to poststack data.

When equation (E-67) is specialized to the zero-offset case, h = 0, we get


(E-71b)

which describes the well-known diffraction hyperbola in the y − t plane for a constant z.

Reduction of the double-square-root equation (E-67) to a single-square-root equation (E-71b) can also be achieved by defining an equivalent offset he such that [10]


(E-72a)

Solving for he, we obtain [10] [11]


(E-72b)

where t is the two-way nonzero-offset traveltime of equation (E-67).

Poststack time migration can be conceptualized either by way of a semicircle superposition using equation (E-72a) or a diffraction summation along the hyperbolic traveltime trajectory using equation (E-72b). Similarly, prestack time migration can be conceptualized either by way of semi-elliptical superposition using equation (E-70) or diffraction summation over the traveltime surface described by equation (E-67). The traveltime surface is known as Cheops’ pyramid Claerbout [12] and is illustrated in Figure E-3a. The result of summation of amplitudes over the pyramidal surface is placed at its apex. The question that is of practical importance is how to define the summation paths over this surface.

Four possible choices of summation paths over the pyramidal surface of equation (E-67) to perform prestack time migration are:

  1. Summation curves of constant offset: Consider a set of vertical cross-sections of the traveltime pyramid illustrated in Figure E-3a parallel to the midpoint axis as shown in Figure E-3b. Sum the amplitudes along each of the constant-offset table-top traveltime curves, independently, and place the result for each at the apex of the summation curve. The summation collapses the pyramidal surface onto a hyperbolic traveltime curve, which is formed by combining the apex points of the constant-offset curves. This hyperbolic traveltime curve is orthogonal to the constant-offset summation curves.
  2. Summation curves of constant time: Consider a set of horizontal cross-sections of the traveltime pyramid as shown in Figure E-3c [13][14]. Sum the amplitudes along each of the constant-time curves, independently, and place the result for each at the maximum offset on the summation curve. The summation, again, collapses the pyramidal surface onto the hyperbolic traveltime curve, which is formed by combining the maximum-offset points of the constant-time curves. This hyperbolic traveltime curve is orthogonal to the constant-time summation curves. The event associated with the resulting hyperbolic moveout trajectory is on the he − t plane of a common-scatter-point (CSP) gather [10], where the equivalent offset he is given by equation (E-72b). The scatter point corresponds to the apex of the traveltime pyramid A0 in Figure E-3c.
  3. Summation curves of constant shot: Consider a set of vertical cross-sections of the traveltime pyramid as shown in Figure E-3d [15]. Sum the amplitudes along each of these constant-shot curves, independently, and place the result for each at the apex of the summation curve. The summation collapses the pyramidal surface onto the hyperbolic traveltime curve, which is formed by combining the apex points of the constant-shot curves. This hyperbolic traveltime curve is orthogonal to the constant-shot summation curves and is on the common-receiver plane that passes through the apex of the pyramid itself.
  4. Summation curves of constant angle: Consider a set of slanted cross-sections of the traveltime pyramid as shown in Figure E-3e [16]. These slanted traveltime curves are associated with constant angle of incidence (the slant-stack transform). Sum the amplitudes along each of the constant-angle curves, independently, and place the result for each at the apex of the summation curve. The summation collapses the pyramidal surface onto the hyperbolic traveltime curve which is formed by combining the apex points of the constant-offset curves. This hyperbolic traveltime curve is orthogonal to the constant-offset summation curves.
Figure E-3  The nonzero-offset traveltime surface associated with a point scatterer and the various summation trajectories for prestack time migration composed from Fowler Fowler [17]. See Section E.5 for details.

E.6 Prestack frequency-wavenumber migration

We start with prestack data in midpoint y and offset h coordinates, P(y, h, z = 0, t), where t is the event time in the unmigrated position, and perform 3-D Fourier transform


(E-73)

where ky, kh, and ω are the Fourier duals of the variables of the midpoint-offset coordinates, y, h, and t.

Extrapolate the prestack data from the surface z = 0 to a depth z by the following extrapolation equation which we borrow from Section D.1


(E-74)

The vertical wavenumber kz is given by


(E-75)

where, the double-square-root (DSR) operator takes the following form in midpoint-offset coordinates (equation D-22 of Section D.1)


(E-76)

The variables Y and H are the normalized midpoint and offset wavenumbers, respectively


(E-77a)

and


(E-77b)

Assume a horizontally layered earth model associated with a vertically varying velocity function v(z). By inverse Fourier transforming equation (E-74), we have


(E-78)

The imaging principle t = 0 then is applied to get the image volume P(y, h, z, t = 0),


(E-79)

This is the equation for the prestack phase-shift method. Equation (E-79) involves integration over frequency ω and 2-D inverse Fourier transformation along midpoint y and offset h axes.

We now consider the special case of constant velocity v. Stolt [18] devised a prestack migration technique that involves an efficient mapping in the 3-D Fourier transform domain from temporal frequency ω to vertical wavenumber kz.

First, combine equations (E-75) and (E-76)


(E-80)

Square both sides


(E-81)

and simplify


(E-82)

Define


(E-83)

and square both sides of equation (E-82). After some algebra, it follows that


(E-84)

which can be rewritten as


(E-85)

Now substitute equations (E-77a,E-77b) and (E-83) into equation (E-85), and simplify to obtain the final expression for the dispersion relation for prestack wave extrapolation


(E-86)

By setting the offset wavenumber kh = 0, we obtain the special case of zero-offset dispersion relation


(E-87)

as in equation (D-85) with x replaced by y.

By keeping the wavenumbers ky and kh unchanged and differentiating equation (E-86), we get


(E-88)

By setting the offset wavenumber kh = 0, we obtain the special case for zero-offset


(E-89)

as in equation (D-86) with x replaced by y.

When equations (E-86) and (E-88) are substituted into equation (E-79), we get


(E-90)

Finally, sum over kh to obtain the image at zero offset, h = 0


(E-91)

This is the equation for constant-velocity prestack Stolt migration. It involves two operations in the f − k domain. First, the temporal frequency ω is mapped onto the vertical wavenumber kz via equation (E-86). Second, the amplitudes are scaled by the quantity


(E-92)

The zero-offset image is then obtained by summing over the wavenumber kh (equation E-91), and inverse Fourier transforming in the midpoint y direction.

E.7 Velocity analysis by wavefield extrapolation

A method of migration velocity analysis based on wavefield extrapolation is described in migration velocity analysis. The main computational steps of this method [19] are outlined below. We work with seismic data in midpoint-(half) offset (y, h) coordinates. We want to obtain a volume of focused energy at zero offset in (y, v, τ) coordinates from a prestack data set in (y, h, t) coordinates. For a midpoint location y, migration velocity function then can be picked from the corresponding (v, τ) plane.

First, a 3-D Fourier transformation is applied to the upcoming wavefield P(y, h, τ = 0, t), which is recorded at the surface


(E-93)

where t is the two-way traveltime and


(E-94)

is the two-way vertical time equivalent of downward continuation depth z in a medium with velocity v(z). The variables (ky, kh, w) are the Fourier duals of (y, h, t).

The surface wavefield given by equation (E-93) then is extrapolated down to depth τ by


(E-95)

where


(E-96)

and Y and H are the normalized midpoint and offset wavenumbers given by equations (E-77a) and (E-77b). The −2 term puts the expression in retarded time form. (This term was not included in the previous definition of DSR given by equation E-76.) Equation (E-95) is used recursively to extrapolate the wavefield from one depth to another in steps of Δτ.

Next, we transform the extrapolated wave field P(ky, kh, τ, ω) into the space-time domain. In doing so, we only need to obtain the zero-offset information (h = 0). By summing the extrapolated wavefield over kh in equation (E-95), we get the wave field at zero offset, P(k, h = 0, τ, ω). By doing the 2-D inverse transform over (ky, ω), we obtain


(E-97)

Here, P(y, h = 0, τ, t) is the zero-offset section at various depth levels from which we want to extract velocity information.

Suppose that velocity ve were used to extrapolate the surface wavefield down to depth τ. Equation (E-95) is written with ve and τ as


(E-98a)

Now suppose that the true medium velocity v were used to extrapolate the surface wavefield down to depth τ = t. By rewriting equation (E-95) with v and t, we have


(E-98b)

Match the two extrapolated wavefields in equations (E-98a) and (E-98b) to get a relationship between ve, τ, v, and t


(E-99)

Because of the complexity of DSR [equation (E-96)], equation (E-99) does not provide an explicit expression for v in terms of the other three variables τ, t, and ve. However, we can get an approximate expression by expanding the square roots in equation (E-96) in the Taylor series and retaining the terms with Y2 and H2, only. By using this approximate form and the definitions given by equations (E-77a) and (E-77b), the following approximate relationship results


(E-100)

This expression suggests that downward continuation with the correct (medium) velocity to a wrong depth is equivalent to downward continuation to the correct depth with the wrong velocity [20].

The derivation of equation (E-100) assumes that ve is constant. When ve is depth-variable, then the relationship in equation (E-100) still holds because equation (E-95) is valid for a stratified earth model. However, quantity ve in that equation is replaced by the rms velocity.

Because the approximation made to equation (E-96) is best for small ratios of offset-to-reflector depth, the accuracy of the mapping procedure based on equation (E-100) degrades at very shallow depths. Refer to migration velocity analysis for the practical considerations of the migration velocity estimation technique described here.[21] [22] [23] [24] [25] [26]

References

  1. Deregowski and Rocca (1981), Deregowski, S. M. and Rocca, F., 1981, Geometrical optics and wave theory for constant-offset sections in layered media: Geophys. Prosp., 29, 374–406.
  2. Levin, 1971, Levin, F.K., 1971, Apparent velocity from dipping interface reflections: Geophysics, 36, 510–516.
  3. 3.0 3.1 3.2 Liner, 1990, Liner, C. L., 1990, General theory and comparative anatomy of dip-moveout: Geophysics, 55, 595–607.
  4. Bolondi et al., 1982, Bolondi, G., Loinger, E. and Rocca, F., 1982, Offset continuation of seismic sections: Geophys. Prosp., 30, 813–828.
  5. Bale and Jacubowicz, 1987, Bale, R. and Jacubowitz, H., 1987, Poststack prestack migration: 57th Ann. Internat. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 714–717.
  6. 6.0 6.1 Notfors and Godfrey, 1987, Marcoux, M. O., Godfrey, R. J., and Notfors, C. D., 1987, Migration for optimum velocity evaluation and stacking: Presented at the 49th Ann. Mtg. European Asn. Expl. Geophys.
  7. 7.0 7.1 Zhou et al. (1996), Zhou, B, Mason, I. M., and Greenalgh, S.A., 1996, Accurate and efficient shot-gather dip-moveout processing in the log-stretch domain: Geophys. Prosp., 43, 963–978.
  8. Biondi and Ronen (1987), Biondi, B. and Ronen, J., 1987, Dip moveout in shot profiles: Geophysics, 52, 1473–1482.
  9. Cabrera and Levy (1989), Cabrera, J. and Levy, S., 1989, Shot dip-moveout with logarithmic transformations: Geophysics, 54, 1038–1041.
  10. 10.0 10.1 10.2 Bancroft et al., 1998, Bancroft, J. C., Geiger, H. D., and Margrave, G., 1998, The equivalent offset method of prestack time migration: Geophysics, 63, 2042–2053.
  11. Margrave et al., 1999, Margrave, G., Bancroft, J. C., and Geiger, H. D., 1999, Fourier prestack migration by equivalent wavenumber: Geophysics, 64, 197–207.
  12. Claerbout (1985), Claerbout, J.F., 1985, Imaging the earth’s interior: Blackwell Scientific Publications.
  13. Bancroft and Geiger (1994), Bancroft, J. C. and Geiger, H. D., 1994, Equivalent-offset CRP gathers: 64th Ann. Internat. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 672–675.
  14. Bancroft et al., 1997, Bancroft, J. C., Margrave, G., and Geiger, H. D., 1997, A kinematic comparison of DMO-PSI and equivalent offset migration (EOM): 67th Ann. Internat. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 1575–1578.
  15. Berryhill, 1996, Berryhill, J. R., 1996, System and method of seismic shot-record migration: U. S. Patent 5,544,126.
  16. Ottolini (1982), Ottolini, R., 1982, Migration of seismic data in angle-midpoint coordinates: Ph.D. thesis, Stanford University.
  17. Fowler, 1997, Fowler, P., 1997, A comparative overview of prestack time migration methods: 67th Ann. Int. Soc. Explor. Geophys. Mtg., Soc. Expl. Geophys., Expanded Abstracts, 1571–1574.
  18. Stolt, 1978, Stolt, R.H., 1978, Migration by Fourier transform: Geophysics, 43, 23–48.
  19. Yilmaz and Chambers, 1984, Yilmaz, O. and Chambers, R., 1984, Migration velocity analysis by wave field extrapolation: Geophysics, 49, 1664–1674.
  20. Doherty and Claerbout (1974), Doherty, S.M. and Claerbout, J.F., 1974, Velocity analysis based on the wave equation: Stanford Expl. Proj., Rep. No. 1, Stanford University.
  21. Deregowski, S. M., 1982, Dip-moveout and reflector-point dispersal: Geophys. Prosp., 30, 318–322.
  22. Deregowski, S. M., 1986, What is DMO?: First Break, 4, 7–24.
  23. Hale, D., 1991, A nonaliased integral method for dip moveout: Geophysics, 56, 795–805.
  24. Hale, D., Hill, N. R., and Stefani, J., 1992, Imaging salt with turning seismic waves: Geophysics, 57, 1453–1462.
  25. Notfors, C. D. and Godfrey, R. J., 1987, Dip-moveout in the frequency-wavenumber domain: Geophysics, 52, 1718–1721.
  26. Rocca, F., Bolondi, G., and Loinger, E., 1982, Offset continuation of seismic sections: Geophys. Prosp., 30, 813–828.

See also

External links

find literature about
Topics in Dip-Moveout Correction and Prestack Time Migration
SEG button search.png Datapages button.png GeoScienceWorld button.png OnePetro button.png Schlumberger button.png Google button.png AGI button.png