Skip Nav Destination

Research Article| October 18, 2019

Velocity-independent Marchenko focusing in time- and depth-imaging domains for media with mild lateral heterogeneity

Yanadet Sripanich;

Yanadet Sripanich

Formerly

Utrecht University

, Department of Earth Sciences, Utrecht,

The Netherlands

; presently

PTT Exploration and Production Public Company Limited

, Bangkok,

Thailand

. E-mail: [email protected] (corresponding author).

Search for other works by this author on:

GSW

Google Scholar

Ivan Vasconcelos;

Ivan Vasconcelos

Utrecht University

, Department of Earth Sciences, Utrecht,

The Netherlands

. E-mail: [email protected].

Search for other works by this author on:

GSW

Google Scholar

Kees Wapenaar

Delft University of Technology

, Department of Geoscience and Engineering, Delft,

The Netherlands

. E-mail: [email protected].

Search for other works by this author on:

GSW

Google Scholar

Author and Article Information

Yanadet Sripanich

Formerly

Utrecht University

, Department of Earth Sciences, Utrecht,

The Netherlands

; presently

PTT Exploration and Production Public Company Limited

, Bangkok,

Thailand

. E-mail: [email protected] (corresponding author).

Ivan Vasconcelos

Utrecht University

, Department of Earth Sciences, Utrecht,

The Netherlands

. E-mail: [email protected].

Kees Wapenaar

Delft University of Technology

, Department of Geoscience and Engineering, Delft,

The Netherlands

. E-mail: [email protected].

Publisher: Society of Exploration Geophysicists

Received: 28 Jan 2019

Revision Received: 11 Jun 2019

Accepted: 03 Aug 2019

First Online: 23 Jan 2020

Online ISSN: 1942-2156

Print ISSN: 0016-8033

Funding

Funder(s):
- Award Id(s):
  320639
Funder(s):
- Award Id(s):
  742703

© 2019 The Authors. Published by the Society of Exploration Geophysicists. All article content, except where otherwise noted (including republished material), is licensed under a Creative Commons Attribution 4.0 International License (CC BY). See https://creativecommons.org/licenses/by/4.0/. Distribution or reproduction of this work in whole or in part commercially or noncommercially requires full attribution of the original publication, including its digital object identifier (DOI).

The Authors

Freely available online through the SEG open-access option.

Geophysics (2019) 84 (6): Q57–Q72.

https://doi.org/10.1190/geo2019-0059.1

ABSTRACT

The Marchenko method retrieves Green’s functions between the acquisition surface and any arbitrary point in the medium. The process generally involves solving an inversion starting with an initial focusing function, e.g., a direct-wave Green’s function from the desired subsurface position, typically obtained using an approximate velocity model. We have formulated the Marchenko method in the time-imaging domain. In that domain, we recognize that the traveltime of the direct-wave Green’s function is related to the Cheop’s traveltime pyramid commonly used in time-domain processing, which in turn can be readily obtained from the local slopes of the common-midpoint gathers. This observation allows us to substitute the velocity-model-based initial focusing operator with that from a data-driven slope estimation process. Moreover, we found that working in the time-imaging domain allows for the specification of the desired subsurface position in terms of vertical time, which is connected to the Cartesian depth position via the time-to-depth conversion. Our results suggest that the prior velocity model is only required when specifying the position in depth, but this requirement can be circumvented by making use of the time-imaging domain within its usual assumptions (e.g., mild lateral heterogeneity). Provided that those assumptions are satisfied, the estimated Green’s functions from the proposed method have comparable quality to those obtained with the knowledge of a prior velocity model.

time-domain, acoustic, reciprocity, multiples, common midpoint

INTRODUCTION

Green’s functions between the surface and any subsurface point are the main ingredient in seismic redatuming and imaging. The Marchenko method is a framework to obtain such information using solely the reflection data at the surface and an initial estimate of the direct-wave Green’s function from the desired subsurface position (Broggini and Snieder, 2012; Broggini et al., 2012; da Costa Filho et al., 2014; Slob et al., 2014; Wapenaar et al., 2014a, 2014b, 2017; Ravasi, 2017; Singh et al., 2017). Given a prior approximate (smooth) velocity model of the subsurface in Cartesian coordinates, one can specify a subsurface position and obtain the direct-wave Green’s function from this position to the surface using conventional forward extrapolation (Wapenaar et al., 2014b; Thorbecke et al., 2017). An alternative strategy involves a separate inversion for the direct-wave Green’s function from the common focal point technology based on the same starting velocity model (Berkhout, 1997; Thorbecke, 1997). Therefore, much like conventional seismic migration, a caveat to the current Marchenko imaging implementation is the requirement of a priori velocity knowledge.

Conventional seismic imaging can be accomplished in either the time or the depth domain. The former generally performs with higher computational efficiency but becomes less accurate than the latter when dealing with geologically complex areas such as subsalt regions (Yilmaz, 2001). The shortcomings of time-imaging methods are largely due to the following (Fomel, 2013, 2014):

1)
Approximate direct-wave Green’s functions are used for imaging, which typically depend on hyperbolic or nonhyperbolic traveltime approximations.
2)
Each time-domain image point is associated with its own approximate effective velocity under the assumption of straight-ray geometry relative to the surface.
3)
When lateral heterogeneity is present, the final images are generated in distorted coordinates defined by image rays (Hubral, 1977) as shown in Figure ¹.

However, in areas with moderately complex geology where such assumptions are approximately valid, we can turn these limitations into advantages. In particular, recent research on time-domain imaging has led to an alternative data-driven time-imaging workflow for improved efficiency and accuracy with local event slopes from common-midpoint (CMP) gathers instead of velocity (Fomel, 2007). This development leads to an opportunity of relying on these velocity-independent data-driven techniques from time-domain imaging to estimate the direct-wave Green’s functions that can be used in a new form of the Marchenko method.

In this paper, we first study the Marchenko method in the time-imaging domain and establish relationships between the focusing functions obtained from the Marchenko methods in the time- and depth-imaging domains. Making use of the slope-based time-domain processing workflow, we subsequently propose a scheme to obtain the direct-wave Green’s function from the desired subsurface position on a reflector to the surface. We show that the newly estimated direct-wave Green’s function can be used in the Marchenko method and leads to comparable results to those that rely on the prior knowledge of a velocity model. We rely on two synthetic models with and without lateral heterogeneity and a field data set from the North Sea to demonstrate the proposed approach.

MARCHENKO METHOD IN THE TIME-IMAGING DOMAIN

Reciprocity theorems on curvilinear surfaces

The key components to deriving the single-sided Green’s function representations for the Marchenko method are the one-way reciprocity theorems of convolution and correlation type (Wapenaar and Grimbergen, 1996; Wapenaar et al., 2014a; van der Neut et al., 2015). Assuming that the image rays are well defined with no caustics, we first recognize that a constant depth level in the Cartesian coordinates generally corresponds to a curved level in the time-imaging domain and vice versa (Figure ²). In other words, the current Marchenko method has already been implemented with respect to a curved level in the time-imaging domain. To show that a converse relationship exists, we need to find the single-sided representations for a curvilinear level in depth that corresponds to a constant time surface.

Because the time-imaging domain is defined by image rays (Figure ¹), it represents a special curvilinear coordinate system of the semiorthogonal type (Sava and Fomel, 2005) due to the orthogonality between the ray direction and the wavefront (Figure ²). In such curvilinear systems, we denote coordinates along the curved boundary surfaces as $ξ = (ξ_{1}, ξ_{2}, ξ_{3})$ and assume that there exists a one-to-one reversible map between the ray coordinates $ξ$ and the physical Cartesian coordinates. We use $ξ$ in our subsequent derivation to ensure its generality and will specify the meaning to each time-domain coordinate $ξ_{i}$ in a later section.

The one-way reciprocity theorems for semiorthogonal curvilinear systems in the Fourier

ω

domain are given by (Frijlink and Wapenaar, 2010):

\int_{S_{a}} (p_{A}^{+} p_{B}^{-} - p_{A}^{-} p_{B}^{+}) d ξ_{a} = \int_{S_{f}} (p_{A}^{+} p_{B}^{-} - p_{A}^{-} p_{B}^{+}) d ξ_{f},

(1)

\int_{S_{a}} (p_{A}^{+} p_{B}^{+ *} - p_{A}^{-} p_{B}^{- *}) d ξ_{a} = \int_{S_{f}} (p_{A}^{+} p_{B}^{+ *} - p_{A}^{-} p_{B}^{- *}) d ξ_{f},

(2)

where equations and represent the one-way reciprocity theorems of convolution and correlation type, respectively. Here,

p

denotes the flux-normalized wavefields in the frequency domain decomposed into upgoing

(-)

and downgoing

(+)

constituents with respect to

ξ_{3}

at the acquisition surface

S_{a}

and the focusing surface

S_{f}

⁠. The superscript

*

denotes the complex conjugation. The integrations are done along the coordinates

ξ_{a}

and

ξ_{f}

at surfaces

S_{a}

and

S_{f}

⁠, which no longer need to represent constant-depth surfaces. Similar to the current derivation of the Marchenko method, the subscripts

A

and

B

denote the two acoustic states — the truncated medium and the true medium, respectively (Figure ³). Moreover, the considered volume between

S_{a}

and

S_{f}

is assumed to have equal medium parameters in both states and is source free.

Marchenko equations in image-ray coordinates

We can clearly observe that equations 1 and 2 are similar to the one-way reciprocity theorems for the case of flat datum levels except that the integrations are now done over curvilinear surfaces. Therefore, we can follow a similar procedure as in the previous works and derive the single-sided Green’s function representations (Wapenaar et al., 2014a, 2014b; van der Neut et al., 2015). We first consider two acoustic states shown in Figure ³ and specify the pertaining parameters in equations 1 and 2 as follows:

1)
State A at $S_{a}$ ⁠:
- •
  $p_{A}^{+} = δ (ξ_{a} - ξ_{a}^{A})$ representing a point source for a downgoing wavefield just above $S_{a}$ ⁠.
- •
  $p_{A}^{-} = R_{A} (ξ_{a}^{A}, ξ_{a}, ω)$ representing the reflected response in the truncated medium from $ξ_{a}$ to $ξ_{a}^{A}$ ⁠.
2)
State A at $S_{f}$ ⁠:
- •
  $p_{A}^{+} = T (ξ_{f}, ξ_{a}^{A}, ω)$ representing the transmission response from the surface $S_{a}$ to $S_{f}$ ⁠.
- •
  $p_{A}^{-} = 0$ indicating no reflection from below $S_{f}$ ⁠.
3)
State B at $S_{a}$ ⁠:
- •
  $p_{B}^{+} = δ (ξ_{a} - ξ_{a}^{B})$ representing a point source for a downgoing wavefield just above $S_{a}$ ⁠.
- •
  $p_{B}^{-} = R (ξ_{a}^{B}, ξ_{a}, ω)$ representing the total reflected response in the true medium from $ξ_{a}$ to $ξ_{a}^{B}$ ⁠.
4)
State B at $S_{f}$ ⁠:
- •
  $p_{B}^{+} = g^{+} (ξ_{f}, ξ_{a}^{B}, ω)$ representing the downgoing part at $S_{f}$ of the Green’s function between $S_{a}$ and $S_{f}$ ⁠.
- •
  $p_{B}^{-} = g^{-} (ξ_{f}, ξ_{a}^{B}, ω)$ representing the upgoing part at $S_{f}$ of the Green’s function between $S_{a}$ and $S_{f}$ ⁠.

The superscripts in

ξ_{a}^{A}

and

ξ_{a}^{B}

are used to distinguish between the coordinates at

S_{a}

of state A and B. The

ξ_{a}

in the argument of

R

and

R_{A}

is used to denote a general position on

S_{a}

⁠. From this point onward, we will also omit the argument

ω

for conciseness. The above substitution leads to

R (ξ_{a}^{B}, ξ_{a}^{A}) - R_{A} (ξ_{a}^{A}, ξ_{a}^{B}) = \int_{S_{f}} T (ξ_{f}, ξ_{a}^{A}) g^{-} (ξ_{f}, ξ_{a}^{B}) d ξ_{f}

(3)

and

δ (ξ_{a}^{B} - ξ_{a}^{A}) - \int_{S_{a}} R_{A} (ξ_{a}^{A}, ξ_{a}) R^{*} (ξ_{a}^{B}, ξ_{a}) d ξ_{a} = \int_{S_{f}} T (ξ_{f}, ξ_{a}^{A}) g^{+ *} (ξ_{f}, ξ_{a}^{B}) d ξ_{f} .

(4)

To retrieve the Green’s function from the two above equations, we define the focusing function

f_{1}^{+}

as the inverse of transmission

T

in the truncated medium (state A) as follows:

δ (ξ_{f} - ξ_{f}^{A}) = \int_{S_{a}} T (ξ_{f}, ξ_{a}^{A}) f_{1}^{+} (ξ_{a}^{A}, ξ_{f}^{A}) d ξ_{a}^{A},

(5)

where

ξ_{f}^{A}

denotes the desired focusing position on the surface

S_{f}

in the truncated medium (state A). As a result,

f_{1}^{+} (ξ_{a}^{A}, ξ_{f}^{A})

represents a purpose-built downgoing field that travels from

ξ_{a}^{A}

on the surface

S_{a}

to focus at

ξ_{f}^{A}

S_{f}

in the truncated medium. We apply the focusing function (equation ) to equations and , and we integrate over

S_{a}

⁠, which leads to the desired one-sided integral representations:

g^{-} (ξ_{f}^{A}, ξ_{a}^{B}) = \int_{S_{a}} R (ξ_{a}^{B}, ξ_{a}^{A}) f_{1}^{+} (ξ_{a}^{A}, ξ_{f}^{A}) d ξ_{a}^{A} - f_{1}^{-} (ξ_{a}^{B}, ξ_{f}^{A}),

(6)

and

- g^{+ *} (ξ_{f}^{A}, ξ_{a}^{B}) = \int_{S_{a}} R^{*} (ξ_{a}^{B}, ξ_{a}) f_{1}^{-} (ξ_{a}, ξ_{f}^{A}) d ξ_{a} - f_{1}^{+} (ξ_{a}^{B}, ξ_{f}^{A}),

(7)

where we define the reflected response

f_{1}^{-}

of the truncated medium to

f_{1}^{+}

f_{1}^{-} (ξ_{a}, ξ_{f}^{A}) = \int_{S_{a}} R_{A} (ξ_{a}^{A}, ξ_{a}) f_{1}^{+} (ξ_{a}^{A}, ξ_{f}^{A}) d ξ_{a}^{A} .

(8)

In other words,

f_{1}^{-} (ξ_{a}, ξ_{f}^{A})

represents the upgoing reflected response of the truncated medium at

ξ_{a}

on the surface

S_{a}

to the propagation of

f_{1}^{+} (ξ_{a}^{A}, ξ_{f}^{A})

that is designed to focus at

ξ_{f}^{A}

S_{f}

⁠.

From equations –, we can see how different

ξ

variables in the single-sided Green’s function representations have emerged. For example,

g^{-} (ξ_{f}^{A}, ξ_{a}^{B})

corresponds to the Green’s function from

ξ_{a}^{B}

in the true medium (state B) to the focusing position

ξ_{f}^{A}

defined in state A according to equation . For simplicity, we recast equations and by switching the notation according to

ξ_{a}^{B} = ξ_{a}

⁠,

ξ_{f}^{A} = ξ_{f}

and changing the dummy integrating variable to

ξ_{a}^{'}

⁠. This leads to the following final system of equations in the frequency domain:

g^{-} (ξ_{f}, ξ_{a}) = \int_{S_{a}} R (ξ_{a}, ξ_{a}^{'}) f_{1}^{+} (ξ_{a}^{'}, ξ_{f}) d ξ_{a}^{'} - f_{1}^{-} (ξ_{a}, ξ_{f}), - g^{+ *} (ξ_{f}, ξ_{a}) = \int_{S_{a}} R^{*} (ξ_{a}, ξ_{a}^{'}) f_{1}^{-} (ξ_{a}^{'}, ξ_{f}) d ξ_{a}^{'} - f_{1}^{+} (ξ_{a}, ξ_{f}),

(9)

which has the same integral form as the original single-sided Green’s function representations in Cartesian depth coordinates, except for the integration over the curvilinear boundary

S_{a}

⁠. The variables

ξ_{a}

and

ξ_{a}^{'}

are defined on the surface

S_{a}

⁠, whereas

ξ_{f}

is defined at the focusing level

S_{f}

⁠. In case that the acquisition surface

S_{a}

is flat at depth

x_{3}^{'}

⁠, we have

ξ_{a}^{'} = x_{a}^{'} = (x_{1}^{'}, x_{2}^{'}, x_{3}^{'})

⁠. A similar reflection response

R

to the constant-depth case can be used, and we can recast the fully curvilinear system in equation as

g^{-} (ξ_{f}, x_{a}) = \int_{S_{a}} R (x_{a}, x_{a}^{'}) f_{1}^{+} (x_{a}^{'}, ξ_{f}) d x_{a}^{'} - f_{1}^{-} (x_{a}, ξ_{f}), - g^{+ *} (ξ_{f}, x_{a}) = \int_{S_{a}} R^{*} (x_{a}, x_{a}^{'}) f_{1}^{-} (x_{a}^{'}, ξ_{f}) d x_{a}^{'} - f_{1}^{+} (x_{a}, ξ_{f}),

(10)

which has a similar form to the Cartesian-space one-sided representations of the conventional Marchenko approach (e.g., Wapenaar et al., 2014b; van der Neut et al., 2015). We do, however, point out two key points that differentiate this system from its original counterpart:

•
In equation 10, the focusing functions and Green’s functions retain the conventional Cartesian $x_{a}$ as their (data-related) surface argument, but the focal-point argument is $ξ_{f}$ which here lies on the curvilinear depth datum $S_{f}$ connected to a fixed time-domain datum (Figure ²; more on this point below).
•
As a result, the $\pm$ superscripts in this theoretical construct correspond to propagation directions normal to the curvilinear surface local to $ξ_{f}$ ⁠, meaning that the ideal $g^{\pm}$ differs in radiation from their conventional Marchenko counterparts at some fixed $ξ_{f} = x_{f}$ ⁠. In practice, whether this radiation distinction expresses itself on the retrieved fields is entirely dependent on how one chooses to solve the focusing problem numerically.

Though these are notable differences, it is also very convenient that in the integral kernels of equation 10, we see the conventional Cartesian reflection response $R (x_{a}, x_{a}^{'})$ ⁠. Because of it, this final system of coupled Marchenko equations allows us to solve for the time-domain-coordinate focusing functions from the usual, data-related, reflection response $R$ ⁠, by applying causality constraints to equation 10 as described by Wapenaar et al. (2014b) and van der Neut et al. (2015).

Because equation 10 is similar to that in the case of a constant-depth focusing level, we can argue that the form of the corresponding Marchenko equations remains the same as long as there exists a transformation between the Cartesian coordinates and some semiorthogonal curvilinear coordinates, whose level curve matches the desired curvilinear datum level. In the case of the time-imaging domain, the coordinate transformation is defined by the mapping of image rays for the time-to-depth conversion (Cameron et al., 2007; Iversen and Tygel, 2008; Sripanich and Fomel, 2018). In this paper, we define the image-ray coordinates as $ξ_{1} = x_{0}$ ⁠, $ξ_{2} = y_{0}$ ⁠, and $ξ_{3} = t_{0}$ (Figure ¹). The first two coordinates $x_{0}$ and $y_{0}$ define the escape location of the image rays at the acquisition surface $S_{a}$ and $t_{0}$ is their one-way traveltime. The curved datum level at depth then corresponds to an image-wavefront surface tied to some constant time $t_{0}$ (Figure ²). Given the same focusing functions at the acquisition surface, the focusing position defined in the time-imaging domain (⁠ $x_{0}$ ⁠, $y_{0}$ ⁠, and $t_{0}$ ⁠) can be translated to its corresponding Cartesian position through the same mapping. Equipped with these results, we can proceed with making use of efficient time-domain techniques to solve the Marchenko equations and obtain focusing functions associated with some specified position in the time-imaging domain.

SLOPE-BASED TIME-DOMAIN PROCESSING

The process of time-domain imaging can be conceptually summarized as shown in Figure ⁴ (Fowler, 1997). The recorded CMP data, which are dependent on midpoint

m

⁠, half-offset

h

⁠, and time

t

are first migrated to zero-offset (⁠

x_{n}

⁠,

t_{n}

⁠) through a combination of normal and dip moveout operations. Poststack time migration subsequently maps the result to (⁠

x_{0}

⁠,

t_{0}

⁠) to correct subsurface reflection-point positions. The entire process constitutes the time-imaging routine, and it is equivalent to the prestack time migration process — stacking along the Cheops traveltime pyramid (Claerbout, 1996). Fomel (2007) shows that under the regular assumptions of hyperbolic traveltime and straight-ray geometry of time imaging, prestack time migration (mapping) can be done with local event slopes as follows:

t_{0}^{2} = \frac{t p_{h} [(t - h p_{h})^{2} - h^{2} p_{m}^{2}]^{2}}{4 (t - h p_{h})^{2} [t p_{h} + h (p_{m}^{2} - p_{h}^{2})]},

(11)

x_{0} = m - \frac{h t p_{m}}{t p_{h} + h (p_{m}^{2} - p_{h}^{2})},

(12)

where

t

is the two-way reflection traveltime. The terms

p_{m} = \partial t / \partial m

and

p_{h} = \partial t / \partial h

are estimated local event slopes from the CMP gathers in the midpoint and half-offset directions, respectively.

From Figure ⁴, we can deduce that the traveltime of direct-wave Green’s function from the same subsurface position is represented by the same value of $t_{0}$ in this mapping (migration) process. Consequently, the desired traveltime of the direct-wave Green’s function from that location is a contour of time $t_{0}$ of the one-way traveltime map (equation 11). Therefore, we can summarize the steps to obtain slope-based direct-wave Green’s function and solve the Marchenko system as follows:

1)
Given the CMP gathers, measure the slopes of primaries using methods such as plane-wave destruction (Fomel, 2002), which iteratively solves a partial differential equation that governs a local plane wave.
2)
Generate the traveltime $t_{0} (m, h, t)$ and distance $x_{0} (m, h, t)$ maps according to equations 11 and 12.
3)
Remap the $t_{0}$ data according to the $x_{0}$ data to correctly position the traveltime across the midpoint coordinate and obtain $t_{0} (x_{0}, h, t)$ ⁠.
4)
Specify a desired focusing position in (⁠ $x_{0}$ ⁠, $t_{0}$ ⁠) and obtain the traveltime of the direct-wave Green’s function from the corresponding contour of $t_{0} (x_{0}, h, t)$ ⁠.
5)
Run a hyperbolic regression using picked $t_{0}$ to ensure a smooth hyperbolic curve on discretized grids and obtain a best-fitted velocity.
6)
Generate the hyperbolic traveltime curve from the obtained velocity and convolve it with a zero-phase wavelet.
7)
Scale the amplitudes with appropriate geometric spreading factor for straight-ray geometry to obtain an approximate direct-wave Green’s function.
8)
Use the time-reversed direct-wave Green’s function as the initial focusing function and solve the time-constrained Marchenko system that follows from equation 10 (Wapenaar et al., 2014b; van der Neut et al., 2015).

EXAMPLES

In this section, we first apply the proposed workflow to two cases of synthetic horizontally layered and laterally heterogeneous media. The primary-only CMP gathers to be used for slope estimation are generated with Kirchhoff modeling based on an accurate two-point ray tracing scheme (Sripanich and Fomel, 2014), whereas the full data set with all orders of multiples used in the construction of the reflection operator $R$ in equation 10 is obtained from finite-difference modeling. We subsequently show a field-data example from the North Sea previously studied by Szydlik et al. (2007) and Ravasi et al. (2015, 2016). In practice, we note that the required CMP gathers for slope estimation can be obtained from the results of prior demultiple processes using, for example, Radon-based filtering of multiples. These demultipled gathers are the same inputs that are generally used in prestack time/depth migrations.

Horizontally layered model

We first consider a horizontally layered model shown in Figure ⁵ and look at the Green’s function from (0 m, 1000 m) on the third reflector. The input CMP corresponding to this position is shown in Figure ⁶. Its estimated dip and the associated one-way traveltime map $t_{0} (x_{0}, h, t)$ generated according to the proposed workflow are shown in Figure ⁶ and ⁶, respectively. The contour of $t_{0} = 0.4753 s$ is chosen in association with the specified point (0 m, 1000 m), and it is shown in Figure ⁶ as a dotted-dashed black line. The contour is subsequently convolved with a choice of zero-phase wavelet and multiplied by the weight to properly account for geometric spreading. For simplicity, we express this weight along the desired traveltime contour as a panel (Figure ⁶). The weighting can then be easily accomplished by a point-wise multiplication to the panel of the convolved traveltime contour. This weight is determined from the ratio $1 / \sqrt{r}$ in a 2D medium (or $1 / r$ in 3D), where $r$ is the distance from the specified point to different positions on the surface in the time-imaging domain. The use of this ratio is permissible due to the straight-ray and effective velocity assumptions used in the time-imaging domain.

The true Green’s function from forward modeling is shown in Figure ⁷ overlain by the traveltime prediction (the dotted-dashed magenta line) of the direct wave using the proposed slope-based workflow. We use this prediction to generate the initial focusing function and solve the Marchenko system (equation 10). The result from the usual velocity-based workflow (van der Neut et al., 2015), with the initial focusing function generated by numerical wave propagation in a background velocity model, is shown in Figure ⁷ in comparison with that from the proposed slope-based method in Figure ⁷. A comparison between the central traces of the true and estimated Green’s functions from Figure ⁷ is shown in Figure ⁸. We observe that the result from the slope-based workflow is slightly inferior due to the assumption of zero-phase wavelet, which leads to small errors in the predicted coda. That aside, the results are comparable in quality indicating the validity of the proposed method. We also investigate another potential choice of wavelet that corresponds to the characteristic phase behavior of 2D far-field Green’s functions in the “Discussion” section.

Laterally heterogeneous model

Next, we turn to a laterally heterogeneous model with lateral heterogeneity (Figure ⁹). In this example, the image rays are no longer vertical and the focusing positions in the time- and depth-imaging domains are related through the mapping defined by image rays. We consider the Green’s function from (⁠ $- 35 m$ ⁠, 1000 m) on the third reflector because this is the position at which the image ray originating from (0 m, 0 m) will pass through. Figure ¹⁰ shows the inputs and intermediate products of the proposed slope-based workflow. Note that there is a nonzero dip in the midpoint direction (Figure ¹⁰). In this example, the contour of $t_{0} = 0.4165 s$ is chosen. The final comparison of Green’s functions is shown in Figure ¹¹, and a comparison of central traces is shown in Figure ¹². The results again are comparable in quality, and similar conclusions can be drawn.

We emphasize that the focusing position is defined in terms of (⁠ $x_{0}$ ⁠, $t_{0}$ ⁠) as opposed to (⁠ $x$ ⁠, $z$ ⁠) in the usual Marchenko workflow. To confirm that the specified focusing (⁠ $x_{0}$ ⁠, $t_{0}$ ⁠) translate to the Cartesian (⁠ $x$ ⁠, $z$ ⁠) according to the image ray mapping, we back propagate the computed focusing function, and the result is shown in Figure ¹³. The dashed line is vertical, whereas the solid line is the image ray originating from (0 m, 0 m). We can observe that the response is now positioned along the image ray and is at the reflector that we chose to compute our slope for the traveltime prediction in the first place.

Field-data example from the North Sea

For the last example, we test our proposed approach against an ocean-bottom cable data set from the North Sea. The reference velocity model in Cartesian coordinates of this data set is shown in Figure ¹⁴. The receiver line placed at the sea bottom of 92 m depth spans a length of 6 km (from 3147 to 9098 m) and contains 235 sensors with an approximate interval of 25 m. The shot interval is 50 m along the 12 km sail line across the extent of the model. We specify two focusing points at a shallow position (6000 m, 1306 m) and at a deep position (6000 m, 2167 m), where we will compare the results of our approach with those from the velocity-based workflow.

To facilitate the comparison between results from our proposed approach and those from the previous velocity-based workflow, we first convert the reference velocity model (Figure ¹⁴) from depth to time according to a 1D-medium assumption. This process is approximately valid in this case due to the nearly lateral homogeneous characteristic of the model close to the desired focusing positions around the midpoint 6000 m. A comparison between the velocity models before and after this process is shown in Figure ¹⁵. We will use the velocity in Figure ¹⁵ to obtain a CMP gather for our slope-based workflow that in principle, should allow a direct comparison between our results and those from the velocity-based workflow in depth.

An example raw CMP gather of this data set at the location 6000 m is shown in Figure ¹⁶. To obtain an accurate estimation of slopes that corresponds to primary reflections, we apply a demultiple process with the use of a hyperbolic Radon transform based on the converted reference velocity in Figure ¹⁵. This process amounts to transforming the input raw CMP gather using the hyperbolic Radon transform with given reference velocity, selective windowing of nearly flat events, and inverse transforming these selected data. This requirement of primary-only gathers in our proposed approach substitutes the previous need of a depth velocity model in the conventional Marchenko workflow. As described in the previous section, the slope of primary-only gathers is used to obtain the initial focusing function, but the reflection operator $R$ in equation 10 still makes use of the full data set containing all orders of multiples. We will discuss more about this point in the next section. The resulting demultipled gather at a similar location is shown in Figure ¹⁷ together with the estimated slope and chosen $t_{0}$ contours at 0.6 and 1.0 s associated with the two focusing positions. Figure ¹⁸ shows the direct-wave Green’s functions obtained from the conventional workflow with wave extrapolation in the reference velocity model (Figure ¹⁴) overlain by the traveltime prediction (the dotted-dashed magenta line) using the proposed slope-based workflow. We can observe a good agreement between the results from both approaches, which corroborates our proposition. In our experiment, we observe that the processing step on hyperbolic regression plays a crucial role in controlling the final quality of the traveltime predictions. The number of samples of picked traveltime contours used for the regression must be chosen in accordance with the offset range and the quality of measured slopes at a large offset. Here, we choose to run a regression in the offset range of $- 350$ to 350 and $- 600$ to 600 for the shallow and deep focusing positions, respectively.

Using the obtained traveltime predictions, we follow the proposed workflow and obtain the final estimates of $f_{1}^{+}$ ⁠, $f_{1}^{-}$ ⁠, $g^{+}$ ⁠, $g^{-}$ for the shallow focusing position shown in Figure ¹⁹. Their counterparts from the conventional velocity-based workflow are shown in Figure ²⁰. The results appear comparable in quality to similar prominent events present and no apparent artifacts. To further compare the results, we investigate the difference between Figures ¹⁹ and ²⁰ as shown in Figure ²¹. We can observe that the only difference in estimated $f_{1}^{+}$ and $g^{+}$ corresponds to the estimate of the direct-wave Green’s function. This is not surprising because we arbitrarily choose a zero-phase wavelet as an approximate waveform in our workflow. However, the estimated $f_{1}^{-}$ and $g^{-}$ appear quite different between the results from the two methods. We emphasize that this observation does not invalidate our slope-based workflow because in the Marchenko system (equation 10), as discussed by van der Neut et al. (2015), all of the physical events in the upgoing Green’s function $g^{-}$ are retrieved with the correct amplitudes at the first iteration by the action of the reflection operator $R$ on the initial focusing function (the time-reversed direct-wave Green’s function), whereas later iterations amount to eliminating artifacts. However, physical events in the downgoing Green’s function $g^{+}$ are updated throughout all the iterations. Therefore, the fact that we do not observe any difference in the $g^{+}$ panel apart from the direct-wave Green’s function, but do observe some small difference throughout the $g^{-}$ panel — while using the same $R$ in both experiments — is due entirely to the difference in our choice of the direct-wave Green’s function (i.e., the zero-phase waveform). The absence of any prominent artifacts between both results also suggests that our Marchenko solver is functioning properly. This finding is also supported by a comparison of the central traces of the estimated $f_{1}^{-}$ and $g^{-}$ as shown in Figure ²². There, we observe that the coda events captured by the velocity-based and the proposed workflows are similar, except for minor amplitude differences and slight shifts due to the use of different direct-wave Green’s functions. These relatively small mismatches lead to the errors shown in the difference panels in Figure ²¹ and ²¹. Similar results for the deep focusing position are shown in Figures ²³, ²⁴, ²⁵, and ²⁶, from which the same conclusions can be drawn.

Finally, it is important to emphasize that in practice, we do not have access to such a reference velocity model (Figure ¹⁴) and its corresponding velocity in time (Figures ¹⁵ and ²⁷) at the early stage of processing. Given raw CMP gathers, we would have to follow the conventional preprocessing procedure (e.g., surface-consistent correction, deghosting, and demultiple) and obtain primary-only gathers before implementing the proposed workflow. These processes may also involve a choice of background velocity trend that requires quality control. Therefore, a more realistic use of this data set can be achieved by first picking a background velocity trend from the velocity analysis (Figure ²⁷). We subsequently use this picked velocity in our simple demultiple process with the hyperbolic Radon transform. Figure ²⁸ shows the resulting gather together with the estimated slope and chosen $t_{0}$ contours associated with the two focusing positions. The final traveltime predictions after hyperbolic regression are shown in Figure ²⁹ indicating a good performance. In this example, these predictions lead to almost similar results to those in Figures ¹⁹ and ²³.

DISCUSSION

We emphasize first that using the proposed method, a prior knowledge on a smooth velocity model is not required and a forward wavefield extrapolation is no longer needed to generate the direct-wave Green’s function (the time-reversed initial focusing function). By using the knowledge of local slopes, we show how to obtain the direct-wave Green’s function directly from the primary-only CMP gathers under the usual assumptions of time-domain imaging. The original recorded gathers are still required to construct the reflection operator $R$ used in the Marchenko system.

In view of time-imaging assumptions, we emphasize that in the presence of strong lateral heterogeneity, image rays that define the time-imaging map may cross-cut each other and lead to caustics. In this study, we restrict ourselves to cases in which the lateral heterogeneity is mild, and the image rays are uniquely defined, with no complications from caustics. Future investigations are required to properly formulate the Marchenko equations appropriate for time-domain imaging under such circumstances. Moreover, we rely on the hyperbolic traveltime approximation, which is exact for plane reflectors in a constant-velocity (homogeneous isotropic) or elliptically anisotropic overburden and approximately valid for small-offset data in other cases. For example, we can already observe the decrease in accuracy in Figure ¹¹, where the estimated traveltime denoted by the dotted-dashed magenta line deviates slightly at larger distances. In such cases, other (more accurate) traveltime approximations can, in principle, be used and additional measurements (more parameters) in addition to local slopes may be needed (Khoshnavaz et al., 2016; Stovas and Fomel, 2016). The topic of traveltime approximation is one of the primary research subjects in time-domain imaging, and there are many choices of traveltime approximation with different ranges of accuracy and applicability (Tsvankin and Thomsen, 1994; Alkhalifah and Tsvankin, 1995; Gelchinsky et al., 1999; Jäger et al., 2001; Fomel and Stovas, 2010; Santos et al., 2011; Blias, 2013; Dell et al., 2013; Farra et al., 2016; Ravve and Koren, 2017; Sripanich et al., 2017).

Apart from the limitations that come with those time-imaging assumptions, we observe that the choice of wavelet also plays a notable role in controlling the final quality of the results from the Marchenko method and in ensuring consistent results from both setups in the time- and depth-imaging domains. In the current approach, the amplitude information is only partially handled by weighting according to the geometric spreading factor corresponding to the straight-ray assumption from the focusing position. An alternative is to also consider the behavior of the far-field Green’s function in a 2D constant-velocity medium. This involves an additional 45° phase shift and a division by $\sqrt{ω}$ ⁠. The results from our proposed slope-based workflow with this choice of wavelet in the two synthetic models are shown in Figure ³⁰. Similar to the case of a zero-phase wavelet, we can still observe minor differences between the exact Green’s functions and the estimated ones. Finding a dynamically appropriate waveform with nonzero phase from slopes using, for example, wavelet estimation methods and far-field signatures is the subject of future research. In a recent notable example, the augmented Marchenko approach by Dukalski et al. (2019) relies on energy conservation for the focusing fields to yield amplitude and phase corrections to the initial focusing functions — the approach currently applies only to laterally invariant media, but the framework could probably be extended to more complex scenarios. In our case, we note that our current choice of a zero-phase wavelet will still lead to a focus at the desired position, but the focused response will contain a slightly different amplitude signature than that computed from the velocity-based workflow, which intrinsically involves a nonzero-phase wavelet. Nonetheless, as shown by our experiments, the Marchenko method can still provide appropriate coda waveforms for focusing functions and Green’s functions while accommodating various wavelet choices for the initial focusing functions used as an input to the system. In a companion paper (Sripanich and Vasconcelos, 2019), we discuss additional insights on Marchenko focusing, by considering the method in the time-imaging domain and analyzing the amplitude and directionality of the focused responses with respect to the surface data aperture.

Accurate knowledge of local slopes also plays a crucial role in the traveltime approximation of the direct-wave Green’s functions. In this study, we use the dip estimator based on the plane-wave destruction filter described by Fomel (2002), which requires the computational cost of $O (N_{d} N_{f} N_{i})$ ⁠, where $N_{d}$ is the number of data points, $N_{f}$ is the length of the filter, and $N_{i}$ is the total number of iterations. Several other methods exist for this purpose (Harlan et al., 1984; Marfurt, 2006; Hale, 2007; Schleicher et al., 2009; Chen et al., 2013) and should be chosen on a case-by-case basis depending on the desired level of accuracy and efficiency.

Moreover, in our scheme, it is crucial that only the slopes of target primary reflection events are used for traveltime prediction. This particular requirement is similar to most processing/imaging techniques that work with primary reflections. In view of the Marchenko method, this replaces the need of an approximate depth velocity model, which is relatively much harder to obtain and requires costly tomographic or full-waveform inversion. Provided the primary-only gathers, our method can readily provide estimates of the initial focusing functions based on data-driven slope measurements without directly computing for subsurface velocity. However, if the gathers still contain residual multiples, the traveltime predictions in equations 11 and 12 are expected to be erroneous because both equations were derived using the geometry of primary reflections (Fomel, 2007). For example, in our horizontally layered model, if the raw gather with primaries and multiples is used in the slope estimation process, the results would look like those shown in Figure ³¹. In comparison to those in Figure ⁶, we observe that the local slopes of primaries and multiples are estimated, leading to inaccurate traveltime estimates. One should expect that such errors will be higher in more complex models, in which primaries and multiples may overlap and lead to conflicting dips. An expeditious multiple removal based on local slopes (Figure ³¹) using simple velocity filtering exists (Cooke et al., 2009) and may be incorporated in the workflow.

We note that the recent development of the plane-wave Marchenko method (Meles et al., 2018) can also potentially benefit from a consideration based on the time-imaging domain. In this scheme, a different focusing condition is proposed to solve for

F_{1}^{+} (x_{a}, z_{f}) = \int_{S_{f}} f_{1}^{+} (x_{a}, x_{f}) d x_{f},

(13)

which represents a summation of conventional focusing functions along some focusing surface

S_{f}

at depth

z_{f}

⁠. Using the one-way reciprocity theorem for curvilinear coordinates (equations and ) and considering the time-imaging domain, one can follow a similar procedure and formulate a system of equations to solve for

F_{1}^{+} (ξ_{a}, t_{f}) = \int_{S_{f}} f_{1}^{+} (ξ_{a}, ξ_{f}) d ξ_{f},

(14)

which represents a summation of focusing functions along some focusing surface

S_{f}

(image wavefront) at vertical time

t_{f}

⁠. The use of time-imaging tools (such as local slopes) in this context requires further investigation, but we note that the focusing position defined by

t_{f}

in this case, again, does not require the knowledge of velocity as

t_{f}

is directly related to the coordinates of the recorded CMP data.

Finally, we highlight that this work can also be deemed complementary to that of van der Neut and Wapenaar (2016). In their study, it was shown that the single-sided Green’s function representations in the Cartesian coordinates could be modified by redatuming with a choice of direct-wave Green’s function related to some fictitious reflector. By doing so, revised single-sided representations can be developed, and the initial focusing function (strictly the inverse of the direct-wave Green’s function) for the corresponding Marchenko method becomes a band-limited delta function. The direct-wave Green’s function that was used to write these revised representations in the first place does not need to be specified explicitly, conveniently obviating the requirement for a smooth velocity model. Van der Neut and Wapenaar (2016) subsequently show that in this redatumed domain, it is possible to separate and eliminate the effects of primary and multiple reflections from the overburden above the chosen fictitious reflector, and an effective internal multiple elimination technique can be developed. However, this regime is not applicable to Green’s function retrieval because the knowledge of the direct-wave Green’s function is never explicitly computed. On the other hand, in this study, we approach the Marchenko method from the perspective of time-domain imaging and show that it is possible to obtain an approximate direct-wave Green’s function straightforwardly from primary-only CMP gathers using local slopes without needing any velocity model. Therefore, it is conceptually possible to adopt the method of van der Neut and Wapenaar (2016) or its variant (Zhang and Staring, 2018) to obtain primary-only gathers needed by our proposed approach, which then leads directly into our velocity-independent Marchenko-based focusing and redatuming workflow.

CONCLUSION

In this paper, we formulate a new form of the Marchenko system in the time-imaging domain defined by image rays. We show that the resulting Marchenko equations have the same form as in the previous case of a constant-depth datum level, except that the focal-point coordinates and radiation directions are now determined along curved space boundaries of constant time (image-wavefront surfaces) that define the time-imaging domain. A priori knowledge of subsurface velocity is no longer needed to generate the initial focusing function and the Marchenko scheme can be carried out by making use of local event slopes measured directly from recorded CMP gathers. The focusing positions are now defined in terms of image-ray coordinates related to the surface location and vertical traveltime, easily obtained from the data without any prior knowledge of the subsurface model. The resulting focusing functions and Green’s functions obtained using the proposed method are of comparable quality to those conventionally obtained with prior knowledge of the subsurface velocity.

ACKNOWLEDGMENTS

We thank J. Brackenhoff, F. Broggini, S. Fomel, G. Meles, J. Thorbecke, J. Trampert, and C. Urruticoechea for helpful discussions. This work is supported by the European Research Council (ERC) under the European Union’s Seventh Framework Programme (FP/2007-2013) grant agreement number 320639 (iGEO) and under the European Union’s Horizon 2020 research and innovation programme grant agreement number 742703.

DATA AND MATERIALS AVAILABILITY

Data associated with this research are available and can be obtained by contacting the corresponding author.

Figure 1.

View large Download slide

The relationship between depth- and time-imaging coordinates. An example image ray originating from $x_{a}$ with an orthogonal slowness vector to the surface is shown. Every point along this ray is mapped to the same lateral distance in the time domain with a different traveltime $t_{a}$ ⁠.

Figure 2.

View large Download slide

The relationship between the variable-depth surface in Cartesian coordinates and the constant time surface in the time-imaging (image-ray) coordinates. The solid arrow lines denote image rays that originate from the acquisition surface $S_{a}$ at right angles. The focusing surface of $S_{f}$ is chosen as the surface of equal image-ray traveltime (image wavefront) that gets mapped to a flat surface in the time-imaging domain.

Figure 3.

View large Download slide

The two acoustic states A and B used in the derivation of the system of Marchenko equations (equation 9). The setup is similar to that for the conventional Marchenko derivation in the Cartesian coordinates, except that the surfaces $S_{a}$ and $S_{f}$ are flat in some choice of semiorthogonal curvilinear coordinates, which are the image-ray coordinates in our case. State A is identical to state B above $S_{f}$ but has a reflection-free region (homogeneous) below. The reflection and transmission responses (operators) in state A are denoted by $R_{A}$ and $T$ ⁠, respectively. State B is the true medium, where the region below $S_{f}$ can be heterogeneous and cause upward scattering. The term $R$ denotes the total reflection response related to the originally recorded data at the surface $S_{a}$ ⁠. The desired down- and upgoing Green’s functions at $S_{f}$ are denoted by $g^{+}$ and $g^{-}$ ⁠, which can be obtained from solving the Marchenko system.

Figure 4.

View large Download slide

A schematic summarizing the concept behind prestack time migration after Fowler (1997).

Figure 5.

View large Download slide

A horizontally layered model with vertical image rays. The star denotes the specified focusing position at (0 m, 1000 m) on the third reflector. Sources and receivers are placed at the surface with lateral spacing of 10 m.

Figure 6.

View large Download slide

Inputs and products from the slope-based workflow in the horizontally layered model: (a) CMP gather from Kirchhoff modeling, (b) estimated dip with a plane-wave destruction filter, (c) one-way traveltime map $t_{0} (x_{0}, h, t)$ ⁠, and (d) the geometric spreading weight calculated based on the straight-ray assumption in the time-imaging domain.

Figure 7.

View large Download slide

A comparison in the horizontally layered model of the true Green’s function (a) and estimated ones using the true velocity to generate initial focusing functions (b) and that from the proposed slope-based workflow (c). The dotted-dashed magenta line in (a) denotes the estimated traveltime of the direct wave using local slopes.

Figure 8.

View large Download slide

A comparison in the horizontally layered model of the central trace of the true Green’s function (blue) and estimated ones using (a) the true velocity to generate initial focusing function and (b) that from the proposed slope-based workflow.

Figure 9.

View large Download slide

A laterally heterogeneous model with bending image rays. The star denotes the specified focusing position at (⁠ $- 35 m$ ⁠, 1000 m) on the third reflector. Sources and receivers are placed at the surface with a lateral spacing of 10 m.

Figure 10.

View large Download slide

Inputs and products from the slope-based workflow in the laterally heterogeneous model: (a) CMP gather from Kirchhoff modeling, (b) estimated dip in the offset direction, (c) estimated dip in the midpoint direction, and (d) one-way traveltime map $t_{0} (x_{0}, h, t)$ ⁠. The geometric spreading weight can be similarly computed based on the straight-ray assumption as before.

Figure 11.

View large Download slide

A comparison in the laterally heterogeneous model of (a) the true Green’s function, (b) estimated ones using the true velocity to generate initial focusing functions, and (c) that from the proposed slope-based workflow. The dotted-dashed magenta line in (a) denotes the estimated traveltime of the direct wave using local slopes.

Figure 12.

View large Download slide

A comparison in the laterally heterogeneous model of the central trace of the true Green’s function (blue) and estimated ones using (a) the true velocity to generate the initial focusing function and (b) that from the proposed slope-based workflow.

Figure 13.

View large Download slide

The focused response in the laterally heterogeneous model. Due to the specification of focusing position in the time-imaging domain (⁠ $x_{0}$ ⁠, $t_{0}$ ⁠), the focused position lies along the image ray originating from (0 m, 0 m) (solid) as opposed to the vertical line (dashed).

Figure 14.

View large Download slide

Volve migration velocity model with the stars denoting the chosen focusing positions at (6000 m, 1306 m) and (6000 m, 2167 m). The source positions are denoted by a green box and cover the entire lateral extent of the model, whereas the receivers are denoted by white upside-down triangles.

Figure 15.

View large Download slide

(a) The original velocity in depth at 6000 m and (b) its corresponding effective velocity expressed in two-way time.

Figure 16.

View large Download slide

An example raw CMP gather of the Volve data set at 6000 m.

Figure 17.

View large Download slide

Inputs and products from the slope-based workflow in the Volve example using the known prior reference velocity model: (a) demultipled CMP gather, (b) estimated dip with plane-wave destruction filter, (c) one-way traveltime map $t_{0} (x_{0}, h, t)$ with a contour related the shallow focusing point (6000 m, 1306 m), and (d) one-way traveltime map with a contour related the deep focusing point (6000 m, 2167 m). Note the different offset range in (c and d).

Figure 18.

The direct-wave Green’s functions obtain from wave extrapolation using the reference velocity model (Figure 14) overlain by the traveltime prediction (the dotted-dashed magenta line) using the proposed slope-based workflow with converted reference velocity (Figure 27a) for the (a) shallow focusing point and the (b) deep focusing point.

View large Download slide

The direct-wave Green’s functions obtain from wave extrapolation using the reference velocity model (Figure ¹⁴) overlain by the traveltime prediction (the dotted-dashed magenta line) using the proposed slope-based workflow with converted reference velocity (Figure ²⁷) for the (a) shallow focusing point and the (b) deep focusing point.

Figure 27.

A comparison of (a) the effective velocity at approximately 6000 m obtained from the depth velocity model (Figure 14) that is approximately laterally homogeneous (1D) and (b) picked velocity from velocity analysis, which is done following the conventional processing procedure.

View large Download slide

A comparison of (a) the effective velocity at approximately 6000 m obtained from the depth velocity model (Figure ¹⁴) that is approximately laterally homogeneous (1D) and (b) picked velocity from velocity analysis, which is done following the conventional processing procedure.

Figure 19.

View large Download slide

Estimated results using the proposed workflow for the shallow focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠. Clipping is adjusted to give a better visualization of the coda.

Figure 20.

View large Download slide

Estimated results using velocity-based workflow for the shallow focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠. Clipping is adjusted to give a better visualization of the coda.

Figure 21.

Difference between the estimated results in Figures 19 and 20 for the shallow focusing position: (a) f1+, (b) f1−, (c) g+, and (d) g−.

View large Download slide

Difference between the estimated results in Figures ¹⁹ and ²⁰ for the shallow focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠.

Figure 22.

View large Download slide

A comparison of the central traces of estimated (a) $f_{1}^{-}$ and (b) $g^{-}$ for the shallow focusing position from the velocity-based workflow (the solid blue lines) in Figure ²⁰ and from the proposed workflow (the dashed red lines) in Figure ¹⁹. We can observe that similar codas (events) are captured by both workflows despite minor amplitude differences and slight shifts in some areas.

Figure 23.

View large Download slide

Estimated results using the proposed workflow for the deep focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠. The clipping is adjusted to give a better visualization of the coda.

Figure 24.

View large Download slide

Estimated results using a velocity-based workflow for the deep focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠. Clipping is adjusted to give a better visualization of the coda.

Figure 25.

The difference between estimated results in Figures 23 and 24 for the deep focusing position: (a) f1+, (b) f1−, (c) g+, and (d) g−.

View large Download slide

The difference between estimated results in Figures ²³ and ²⁴ for the deep focusing position: (a) $f_{1}^{+}$ ⁠, (b) $f_{1}^{-}$ ⁠, (c) $g^{+}$ ⁠, and (d) $g^{-}$ ⁠.

Figure 26.

View large Download slide

A comparison of the central traces of estimated (a) $f_{1}^{-}$ and (b) $g^{-}$ for the deep focusing position from the velocity-based workflow (the solid blue lines) in Figure ²⁴ and from the proposed workflow (the dashed red lines) in Figure ²³. Again, we can observe that similar codas (events) are captured by both workflows despite minor amplitude differences and slight shifts in some areas.

Figure 28.

View large Download slide

Inputs and products from the slope-based workflow in the Volve example as should be implemented in the conventional processing procedure: (a) demultipled CMP gather, (b) estimated dip with plane-wave destruction filter, (c) one-way traveltime map $t_{0} (x_{0}, h, t)$ with a contour related the shallow focusing point (6000 m, 1306 m), and (d) one-way traveltime map with a contour related the deep focusing point (6000 m, 2167 m). Note the different offset range in (c and d).

Figure 29.

View large Download slide

The direct-wave Green’s functions obtain from wave extrapolation using the reference velocity model (Figure ¹⁴) overlain by the traveltime prediction (the dotted-dashed magenta line) using the proposed slope-based workflow with picked velocity (Figure ²⁷) for the (a) shallow focusing point and the (b) deep focusing point.

Figure 30.

View large Download slide

(a) Estimated Green’s function and (b) its middle trace from the slope-based workflow with modified waveform according to the behavior of the 2D far-field Green’s function in the horizontally layered model. The trace of the true Green’s function is plotted in the solid blue. (c and d) Similar plots for the laterally heterogeneous model. Similar to the results with zero-phase wavelet, we observe minor differences between the estimated results and the true values due to an approximate choice of waveform used in our slope-based workflow.

Figure 31.

View large Download slide

Inputs and products from the slope-based workflow in the horizontally layered model: (a) raw CMP gather that contains both the primaries and multiples, (b) estimated dip with a plane-wave destruction filter, (c) one-way traveltime map $t_{0} (x_{0}, h, t)$ ⁠, and (d) the difference between the estimated one-way time in the presence of multiple (c) and that from the primary-only gather in Figure ⁶. We can observe an error in traveltime predictions when the gather does not contain only primary reflections.

Velocity-independent Marchenko focusing in time- and depth-imaging domains for media with mild lateral heterogeneity

ABSTRACT

INTRODUCTION

MARCHENKO METHOD IN THE TIME-IMAGING DOMAIN

Reciprocity theorems on curvilinear surfaces

Marchenko equations in image-ray coordinates

SLOPE-BASED TIME-DOMAIN PROCESSING