Marchenko redatuming, imaging and multiple elimination, and their mutual relations

(1)

Marchenko redatuming, imaging and multiple elimination, and their mutual relations

Wapenaar, Kees; Brackenhoff, Joeri; Dukalski, Marcin; Meles, Giovanni; Slob, Evert; Staring, Myrna;

Thorbecke, Jan; Neut, Joost van der; Zhang, Lele; Reinicke Urruticoechea, C.

DOI

10.1190/geo2020-0854.1

Publication date

2021

Document Version

Accepted author manuscript

Published in

Geophysics

Citation (APA)

Wapenaar, K., Brackenhoff, J., Dukalski, M., Meles, G., Slob, E., Staring, M., Thorbecke, J., Neut, J. V. D.,

Zhang, L., & Reinicke Urruticoechea, C. (Accepted/In press). Marchenko redatuming, imaging and multiple

elimination, and their mutual relations. Geophysics, 1-103. https://doi.org/10.1190/geo2020-0854.1

Important note

To cite this publication, please use the final published version (if applicable).

Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent

of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons.

Takedown policy

Please contact us and provide details if you believe this document breaches copyrights.

We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

(2)

mutual relations

Kees Wapenaar

1 , Joeri Brackenhoff

1,2

, Marcin Dukalski

3 , Giovanni Meles

1,4

, Christian Reinicke

1,3

, Evert Slob

1 ,

Myrna Staring

1,5

, Jan Thorbecke

1 , Joost van der Neut

1 and Lele Zhang

1

1 _Delft

_University

_of

_Technology,

_Department

_of

_Geoscience

_and

Engineering,

Stevinweg

1,

2628

CN

Delft,

The

Netherlands

2 ETH Z¨

urich, Department of Earth Sciences, Sonneggstrasse 5, 8092 Z¨

urich, Switzerland

3 _{Aramco Overseas Company B.V., Informaticalaan 6-12, 2628 ZD Delft, The Netherlands}

4 Institute of Earth Sciences, University of Lausanne, Lausanne, 1015, Switzerland

5 Fugro Innovation & Technology, Prismastraat 4, Nootdorp, The Netherlands

Right-running head: Marchenko methods and mutual relations

(Dated: May 20, 2021)

With the Marchenko method it is possible to retrieve Green’s functions between virtual sources in

the subsurface and receivers at the surface from reflection data at the surface and focusing

func-tions. A macro model of the subsurface is needed to estimate the first arrival; the internal multiples

are retrieved entirely from the reflection data. The retrieved Green’s functions form the input

for redatuming by multidimensional deconvolution (MDD). The redatumed reflection response is

free of internal multiples related to the overburden. Alternatively, the redatumed response can be

obtained by applying a second focusing function to the retrieved Green’s functions. This process

is called Marchenko redatuming by double focusing. It is more stable and better suited for an

adaptive implementation than Marchenko redatuming by MDD, but it does not eliminate the

multiples between the target and the overburden. An attractive efficient alternative is plane-wave

Marchenko redatuming, which retrieves the responses to a limited number of plane-wave sources

at the redatuming level. In all cases, an image of the subsurface can be obtained from the

reda-tumed data, free of artefacts caused by internal multiples. Another class of Marchenko methods

aims at eliminating the internal multiples from the reflection data, while keeping the sources and

receivers at the surface. A specific characteristic of this form of multiple elimination is that it

predicts and subtracts all orders of internal multiples with the correct amplitude, without needing

a macro subsurface model. Like Marchenko redatuming, Marchenko multiple elimination can be

implemented as an MDD process, a double dereverberation process, or an efficient plane-wave

oriented process. We systematically discuss the different approaches to Marchenko redatuming,

imaging and multiple elimination, using a common mathematical framework.

INTRODUCTION

Building on the autofocusing method of Rose (2001,

2002), Broggini and Snieder (2012) showed how the

Marchenko method can be used to retrieve the 1D

Green’s function between a virtual source in the

sub-surface and a receiver at the sub-surface from the reflection

response at the surface. Unlike in seismic

interferome-try (Campillo and Paul, 2003; Wapenaar, 2003;

Schus-ter et al., 2004; Bakulin and Calvert, 2006; Gou´edard

et al., 2008), no physical receiver is needed at the

po-sition of the virtual source. The generalization of this

Green’s function retrieval method to 3D situations

(Wapenaar et al., 2014) formed the basis for the

devel-opment of Marchenko redatuming and imaging

meth-ods (Behura et al., 2014; Broggini et al., 2014). The

main characteristic of these methods is that internal

multiples are dealt with in a data-driven way. A

sub-surface image obtained with the Marchenko method

is free of artefacts related to internal multiples. The

required input consists of the reflection response at

the surface (deconvolved for the seismic wavelet and

free of surface-related multiples) and an estimate of

the direct arrivals of the Green’s functions. The

lat-ter can be obtained from a macro model of the

subsur-face. Hence, the required input is the same as that for

standard redatuming and imaging of primary

reflec-tions; the information needed to deal with the internal

multiples comes entirely from the reflection response

at the surface.

Since the introduction of the Marchenko method in

geophysics, many variants have been introduced. In

the initial approach, redatuming was achieved by

ap-plying multidimensional deconvolution (MDD) to the

downgoing and upgoing Green’s functions retrieved

with the Marchenko method (Broggini et al., 2014;

Ravasi et al., 2016). To obtain a more stable method,

suited for adaptive implementation, redatuming by

double focusing was developed (van der Neut et al.,

2015c; Staring et al., 2018). An important efficiency

gain was achieved with the plane-wave Marchenko

re-datuming approach (Meles et al., 2018). In all these

approaches, sources and receivers are redatumed from

the surface to virtual sources and receivers at one

or more depth levels in the subsurface.

This

re-quires a macro model of the overburden. To make the

Marchenko method less sensitive to the macro model,

it was proposed to extrapolate the virtual sources and

(3)

receivers upward to the acquisition surface (Meles et

al., 2016; van der Neut and Wapenaar, 2016). This led

to a class of Marchenko multiple elimination methods,

i.e., methods in which the sources and receivers stay

at the surface while the internal multiples are

elimi-nated from the data (Zhang et al., 2019a,b; Pereira

et al., 2019; Elison et al., 2020; Dukalski and de Vos,

2020; Meles et al., 2020; Staring et al., 2021).

In this paper we discuss the different Marchenko

methods in a systematic way, show their mutual

re-lations and discuss the specific properties of each

method. By using a consistent way of presenting these

methods, using a unified notation, we hope to convey

the systematics of the many Marchenko methods that

are currently around. The emphasis will be on

expla-nations with cartoon-like figures. Numerical examples

and field data applications can be found in the

refer-enced literature.

It is impossible to discuss all existing Marchenko

methods in a single paper. At various places we

in-clude references for variants that are not discussed

here. In particular, the discussion in this paper is

re-stricted to acoustic methods for lossless media. For

a discussion of the Marchenko method in dissipative

media we refer to Slob (2016) and for elastodynamic

Marchenko methods to Wapenaar and Slob (2014), da

Costa Filho et al. (2014) and Reinicke et al. (2020).

Throughout the paper we assume that the input data

are properly sampled. For Marchenko methods that

compensate for the effects of irregular sampling, see

Haindl et al. (2018), Peng et al. (2019) and van

IJs-seldijk and Wapenaar (2021). Recent developments

on the integration of the Marchenko method with full

waveform inversion are also beyond the scope of this

paper. For this subject we refer to Cui et al. (2020)

and Shoja et al. (2020). Finally, note that all

appli-cations indicated in this paper are restricted to the

seismic reflection method. For a discussion of the

‘ho-mogeneous Green’s function approach’ for

monitor-ing and forecastmonitor-ing the responses to induced seismic

sources, we refer to Brackenhoff et al. (2019a,b).

MARCHENKO REDATUMING AND IMAGING

Seismic redatuming is the process of virtually

mov-ing sources and/or receivers from the acquisition

sur-face to a new depth level (or ‘datum plane’) in the

subsurface. Traditionally this is done with one-way

wave field extrapolation operators (or ‘focusing

oper-ators’) which account for primaries only (Berkhout,

1982; Berryhill, 1984). Classical wave field

extrapo-lation and redatuming methods that account for

in-ternal multiples exist (Wapenaar et al., 1987;

Mul-der, 2005), but they require a very detailed subsurface

model. Redatuming methods that are based on

seis-mic interferometry (Schuster et al., 2004; Bakulin and

Calvert, 2006; van der Neut et al., 2011) do not need

a subsurface model, but they require the presence of

actual receivers at the depth level to which one wants

to redatum the sources.

Broggini and Snieder (2012) showed that with the

Marchenko method the same can be achieved as with

seismic interferometry, at least in 1D, without

requir-ing actual receivers in the subsurface. This was the

inspiration for the research into 3D Marchenko

reda-tuming and imaging, which is extensively discussed in

this section. Marchenko redatuming is a data-driven

method to create virtual sources and receivers in the

subsurface. It accounts for internal multiples in the

overburden, it only needs a macro model of the

over-burden and it does not require the presence of actual

receivers in the subsurface.

Focusing functions

The focusing function plays an essential role in

the Marchenko method. It is a 3D generalization of

the ‘fundamental solution’ in 1D scattering problems

(Lamb, 1980). On the other hand, from the seismic

perspective it can be seen as a generalization of the

3D focusing operator used in traditional redatuming,

accounting for internal multiples. Here we discuss its

basic properties. In the section “Retrieval of focusing

functions” we show how it can be retrieved from the

reflection response at the surface.

Consider a 3D inhomogeneous lossless acoustic

medium, with propagation velocity c(x) and mass

density ρ(x), where x = (x

1 , x

2 , x

3 ) is the Cartesian

coordinate vector. Here x

1 and x

2 are the

horizon-tal coordinates, in the following denoted by vector

x

H

= (x

1 , x

2 ); x

3 is the depth coordinate. For 2D

sit-uations, the coordinate vectors reduce to x = (x

1 , x

3 )

and x

H

= x

1 , respectively. The acquisition

bound-ary at x

3 = x

3,0

is denoted as

S

0 . Throughout this

paper we assume that

_S

0 is a transparent boundary

and that the upper half-space is homogeneous. We

choose a focal point x

A

= (x

H,A

, x

3,A

) in the

subsur-face, with x

H,A

= (x

1,A

, x

2,A

) (or x

H,A

= x

1,A

in the

2D situation) and x

3,A

> x

3,0

, and define a boundary

S

A

at the focal depth x

3,A

. We define a truncated

version of the medium, which is identical to the

ac-tual medium above

_S

A

and reflection free below

S

A

.

We introduce the focusing function f

1 (x, x

A

, t) (with

t denoting time) in this truncated medium. It

con-sists of a downgoing part f

1 +

(x, x

A

, t) and an

upgo-ing part f

₁

−

(x, x

A

, t) (the superscripts + and

−

re-fer to the propagation direction at the first

coordi-nate vector, here x). The downgoing focusing

func-tion f

1 +

(x

S

, x

A

, t), with x

S

at

S

0 , is defined such that

f

1 +

(x, x

A

, t) focuses at x = x

A

and t = 0 and

contin-ues as a diverging downgoing field into the

reflection-free half-space below

_S

A

. The upgoing focusing

func-tion f

₁

−

(x

R

, x

A

, t) is the response of the truncated

(4)

0

500 1000

1500

-2000

-1000

0 1000

2000

1900

2000

2100

2200

2300

2400

2500

De

pt

h

(k

m

)

Lateral position (km)

Ve

lo

ci

ty

(m

/s

)

0

0.5

1

1.5 -2 -1 0 1 2

a)

S

0

<latexit sha1_base64="j2hYY8psKNxn6vuJrxsTeRUAE0Q=">AAAB83icbVDLSgMxFL1TX7W+qi7dBIvgqsxUQRcuCm5cVrQP6Awlk2ba0EwmJBmhDP0NNy4UcevPuPNvzLSz0NYDgcM593JPTig508Z1v53S2vrG5lZ5u7Kzu7d/UD086ugkVYS2ScIT1QuxppwJ2jbMcNqTiuI45LQbTm5zv/tElWaJeDRTSYMYjwSLGMHGSr4fYzMOw+xhNnAH1Zpbd+dAq8QrSA0KtAbVL3+YkDSmwhCOte57rjRBhpVhhNNZxU81lZhM8Ij2LRU4pjrI5pln6MwqQxQlyj5h0Fz9vZHhWOtpHNrJPKNe9nLxP6+fmug6yJiQqaGCLA5FKUcmQXkBaMgUJYZPLcFEMZsVkTFWmBhbU8WW4C1/eZV0GnXvot64v6w1b4o6ynACp3AOHlxBE+6gBW0gIOEZXuHNSZ0X5935WIyWnGLnGP7A+fwB6aaRlA==</latexit>

S

A

<latexit sha1_base64="EHDi2TeBzlI2aXOJ77DEPyjNG5w=">AAAB83icbVDLSgMxFL3xWeur6tJNsAiuykwVdOGi4sZlRfuAzlAyaaYNzWSGJCOUob/hxoUibv0Zd/6NmXYW2nogcDjnXu7JCRLBtXGcb7Syura+sVnaKm/v7O7tVw4O2zpOFWUtGotYdQOimeCStQw3gnUTxUgUCNYJxre533liSvNYPppJwvyIDCUPOSXGSp4XETMKguxh2r/pV6pOzZkBLxO3IFUo0OxXvrxBTNOISUMF0brnOonxM6IMp4JNy16qWULomAxZz1JJIqb9bJZ5ik+tMsBhrOyTBs/U3xsZibSeRIGdzDPqRS8X//N6qQmv/IzLJDVM0vmhMBXYxDgvAA+4YtSIiSWEKm6zYjoiilBjayrbEtzFLy+Tdr3mntfq9xfVxnVRRwmO4QTOwIVLaMAdNKEFFBJ4hld4Qyl6Qe/oYz66goqdI/gD9PkDA3mRpQ==</latexit>

x

<latexit sha1_base64="iGe21gmbD6EMbggxzB3hoFRW/Uo=">AAAB8HicbVBNSwMxEJ2tX7V+VT16CRbBU9mtgj0WvHisaD+kXUo2zbahSXZJsmJZ+iu8eFDEqz/Hm//GbLsHbX0w8Hhvhpl5QcyZNq777RTW1jc2t4rbpZ3dvf2D8uFRW0eJIrRFIh6pboA15UzSlmGG026sKBYBp51gcp35nUeqNIvkvZnG1Bd4JFnICDZWekj7QYieZoO7QbniVt050CrxclKBHM1B+as/jEgiqDSEY617nhsbP8XKMMLprNRPNI0xmeAR7VkqsaDaT+cHz9CZVYYojJQtadBc/T2RYqH1VAS2U2Az1steJv7n9RIT1v2UyTgxVJLFojDhyEQo+x4NmaLE8KklmChmb0VkjBUmxmZUsiF4yy+vknat6l1Ua7eXlUY9j6MIJ3AK5+DBFTTgBprQAgICnuEV3hzlvDjvzseiteDkM8fwB87nD4rWkDQ=</latexit>

S

x

S

x

S

x

S

x

<latexit sha1_base64="lJicLHL/2KIVVtx6fzMZq1+VbKo=">AAAB8HicbVBNSwMxEJ2tX7V+VT16CRbBU9mtgj0WvHisYj+kXUo2zbahSXZJsmJZ+iu8eFDEqz/Hm//GbLsHbX0w8Hhvhpl5QcyZNq777RTW1jc2t4rbpZ3dvf2D8uFRW0eJIrRFIh6pboA15UzSlmGG026sKBYBp51gcp35nUeqNIvkvZnG1Bd4JFnICDZWekj7QYieZoO7QbniVt050CrxclKBHM1B+as/jEgiqDSEY617nhsbP8XKMMLprNRPNI0xmeAR7VkqsaDaT+cHz9CZVYYojJQtadBc/T2RYqH1VAS2U2Az1steJv7n9RIT1v2UyTgxVJLFojDhyEQo+x4NmaLE8KklmChmb0VkjBUmxmZUsiF4yy+vknat6l1Ua7eXlUY9j6MIJ3AK5+DBFTTgBprQAgICnuEV3hzlvDjvzseiteDkM8fwB87nD4lSkDM=</latexit>

R

x

R

x

R

x

A

<latexit sha1_base64="vQbknEF6b5gSqaC6WPLILGHEO7s=">AAAB8HicbVBNSwMxEJ2tX7V+VT16CRbBU9mtgp6k4sVjBfsh7VKyabYNTbJLkhXL0l/hxYMiXv053vw3pts9aOuDgcd7M8zMC2LOtHHdb6ewsrq2vlHcLG1t7+zulfcPWjpKFKFNEvFIdQKsKWeSNg0znHZiRbEIOG0H45uZ336kSrNI3ptJTH2Bh5KFjGBjpYe0F4Toadq/7pcrbtXNgJaJl5MK5Gj0y1+9QUQSQaUhHGvd9dzY+ClWhhFOp6VeommMyRgPaddSiQXVfpodPEUnVhmgMFK2pEGZ+nsixULriQhsp8BmpBe9mfif101MeOmnTMaJoZLMF4UJRyZCs+/RgClKDJ9Ygoli9lZERlhhYmxGJRuCt/jyMmnVqt5ZtXZ3Xqlf5XEU4QiO4RQ8uIA63EIDmkBAwDO8wpujnBfn3fmYtxacfOYQ/sD5/AFxXJAo</latexit>

Red rays : f

+

1 (x, x

A

, t)

<latexit sha1_base64="MfzsbJ1ErLopyz9ov1s57ldc/iw=">AAACHHicbVDLSgMxFM3UV62vqks3wSJUHMpMK1hcVdy4rGIf0KlDJpNpQzMPkoxYhvkQN/6KGxeKuHEh+Dem7Sy09YRwD+fcS3KPEzEqpGF8a7ml5ZXVtfx6YWNza3unuLvXFmHMMWnhkIW86yBBGA1IS1LJSDfiBPkOIx1ndDnxO/eECxoGt3Ickb6PBgH1KEZSSXaxlljchzfEtXTI0VicW/rkpJ6dmOndSTmxHA8+pHpW7QtdHtvFklExpoCLxMxICWRo2sVPyw1x7JNAYoaE6JlGJPsJ4pJiRtKCFQsSITxCA9JTNEA+Ef1kulwKj5TiQi/k6gYSTtXfEwnyhRj7jur0kRyKeW8i/uf1YunV+wkNoliSAM8e8mIGZQgnSUGXcoIlGyuCMKfqrxAPEUdYqjwLKgRzfuVF0q5WzFqlen1aatSzOPLgAByCMjDBGWiAK9AELYDBI3gGr+BNe9JetHftY9aa07KZffAH2tcPx2Of5g==</latexit>

Blue rays : f

₁

(x, x

A

, t)

<latexit sha1_base64="CceY6NLvVAshEkB96JbFaToPHUs=">AAACHXicbVDLSgMxFM3UV62vqks3wSJUqGWmFiyuqm5cVrAP6IxDJs20oZkHSUYchvkRN/6KGxeKuHAj/o1pOwttPSHcwzn3ktzjhIwKqevfWm5peWV1Lb9e2Njc2t4p7u51RBBxTNo4YAHvOUgQRn3SllQy0gs5QZ7DSNcZX0387j3hggb+rYxDYnlo6FOXYiSVZBfrick9eMkiYlYgR7E4NyuTk7p2YqR3J+XEdFz4kFayal9U5LFdLOlVfQq4SIyMlECGll38NAcBjjziS8yQEH1DD6WVIC4pZiQtmJEgIcJjNCR9RX3kEWEl0+1SeKSUAXQDrq4v4VT9PZEgT4jYc1Snh+RIzHsT8T+vH0m3YSXUDyNJfDx7yI0YlAGcRAUHlBMsWawIwpyqv0I8QhxhqQItqBCM+ZUXSadWNU6rtZt6qdnI4siDA3AIysAAZ6AJrkELtAEGj+AZvII37Ul70d61j1lrTstm9sEfaF8/qm6gXw==</latexit>

f

+

1,d

(x, x

A

, t)

<latexit sha1_base64="3aBKZ0qN5goMzMMVVDv0QNe2qP0=">AAACDXicbVDLSsNAFJ3UV62vqEs3g1WoWEpSBbusuHFZwT6giWEymbRDJw9mJmIJ+QE3/oobF4q4de/Ov3HaZqGtBy73cM69zNzjxowKaRjfWmFpeWV1rbhe2tjc2t7Rd/c6Iko4Jm0csYj3XCQIoyFpSyoZ6cWcoMBlpOuOriZ+955wQaPwVo5jYgdoEFKfYiSV5OhHvpOa1dTiAfSy7O60klquDx+yat6dy6o8cfSyUTOmgIvEzEkZ5Gg5+pflRTgJSCgxQ0L0TSOWdoq4pJiRrGQlgsQIj9CA9BUNUUCEnU6vyeCxUjzoR1xVKOFU/b2RokCIceCqyQDJoZj3JuJ/Xj+RfsNOaRgnkoR49pCfMCgjOIkGepQTLNlYEYQ5VX+FeIg4wlIFWFIhmPMnL5JOvWae1eo35+VmI4+jCA7AIagAE1yAJrgGLdAGGDyCZ/AK3rQn7UV71z5mowUt39kHf6B9/gCKC5qM</latexit>

-1.0

-0.5

0

0.5

1.0

1.5 -2000

-1000

0 1000

2000

b)

Lateral position (m)

c)

Lateral position (m)

Ti

m

e

(s

)

Ti

m

e

(s

)

f

+

1,d

(x

S

, x

A

, t)

<latexit sha1_base64="iA84TikRJGl/JCiJYRMLAJplOoU=">AAACD3icbZDLSsNAFIYn9VbrLerSzWBRKpaSVMEuK25cVrQXaGKYTCft0MmFmYlYQt7Aja/ixoUibt26822cthG09YeBj/+cw5nzuxGjQhrGl5ZbWFxaXsmvFtbWNza39O2dlghjjkkThyzkHRcJwmhAmpJKRjoRJ8h3GWm7w4txvX1HuKBhcCNHEbF91A+oRzGSynL0Q89JzHJicR/20vT2uJRYrgfvU+e6/EPnZXnk6EWjYkwE58HMoAgyNRz90+qFOPZJIDFDQnRNI5J2grikmJG0YMWCRAgPUZ90FQbIJ8JOJvek8EA5PeiFXL1Awon7eyJBvhAj31WdPpIDMVsbm//VurH0anZCgyiWJMDTRV7MoAzhOBzYo5xgyUYKEOZU/RXiAeIISxVhQYVgzp48D61qxTypVK9Oi/VaFkce7IF9UAImOAN1cAkaoAkweABP4AW8ao/as/amvU9bc1o2swv+SPv4BvxLm1I=</latexit>

-1.0

-0.5

0

0.5

1.0

1.5 -2000

-1000

0 1000

2000

FIG. 1 (a) Focusing functions f

+

1 and f

−

1 in the truncated

medium.

(b) The focusing function f

₁

+

(x

S

, x

A

, t) (fixed

x

A

, variable x

S

) with x

S

at

S

0 and x

A

at

S

A

. (c) The

focused field f

+

1 (x

0 A

, x

A

, t) (fixed x

A

, variable x

0 A

) with x

0 A

and x

A

at

S

A

.

the downgoing and upgoing functions are visualized in

Figure 1a. Note that downgoing and upgoing waves

meet each other at interfaces in such a way that only

the direct arrival of the focusing function, denoted by

f

_1,d

+

(x, x

A

, t), reaches the focal point.

The propagation of the focusing function through

the truncated medium, from

_S

0 to

S

A

, is formally

de-scribed by

f

1 +

(x

0 A

, x

A

, t) =

(1)

Z

S

0 dx

S

Z

∞

0 T (x

0 A

, x

S

, t

0 )f

1 +

(x

S

, x

A

, t

− t

0 )dt

0 ,

for x

0 A

at

S

A

, where T (x

0 A

, x

S

, t) is the transmission

response of the truncated medium. Throughout the

paper we assume that downgoing and upgoing fields

are power-flux normalized (Frasier, 1970; Kennett et

al., 1978; Ursin, 1983; Chapman, 1994), which

ex-plains why expressions like equation 1 do not contain

the vertical spatial derivative of one of the functions

under the integral. The formal focusing conditions are

f

1 +

(x

0 A

, x

A

, t) = δ(x

0 H,A

− x

H,A

)δ(t),

(2)

f

₁

−

(x

0 A

, x

A

, t) = 0,

(3)

for x

0 A

at

S

A

.

From equations 1 and 2 it

fol-lows that f

1 +

(x

S

, x

A

, t) is by definition the inverse

of T (x

A

, x

S

, t). For the truncated medium of

Fig-ure 1a, the downgoing function f

1 +

(x

S

, x

A

, t),

con-volved with a wavelet, is shown in gray-level display in

Figure 1b. Its direct contribution f

_1,d

+

(x

S

, x

A

, t) is a

hyperbolic-like event at negative time (actually this is

the traditional one-way focusing operator (Berkhout

and Wapenaar, 1993)). If no scattering occurred

be-tween

_S

0 and

S

A

, this would be the complete

focus-ing function. However, in a scatterfocus-ing medium,

ad-ditional events are present in f

1 +

(x

S

, x

A

, t) (as

visu-alized in Figures 1a and 1b), which avoid that

mul-tiply scattered waves reach

S

A

. Figure 1c shows the

focused field f

1 +

(x

0 A

, x

A

, t), with x

0 A

at

S

A

, obtained

by emitting f

₁

+

(x

S

, x

A

, t) (convolved with a wavelet)

into the truncated medium, according to equation

1. The amplitudes are clipped to emphasize the

de-tails.

Note that there are no artefacts related to

multiple scattering. Nevertheless, the focused field

deviates from the desired result, expressed by

equa-tion 2. The explanaequa-tion is that, in practical

situa-tions, the aperture is finite and focusing functions do

not compensate for evanescent waves, which implies a

spatial band-limitation (Berkhout and van Wulfften

Palthe, 1979, App.

C). Moreover, in practice the

seismic wavelet implies a temporal band-limitation.

Hence, the delta functions in equation 2 (and in the

re-mainder of the paper) should be interpreted as

band-limited delta functions. Consequently, in practice the

downgoing focusing function f

1 +

(x

S

, x

A

, t) is actually

a band-limited inverse of the transmission response

T (x

A

, x

S

, t).

Representations

We discuss two representations, which formulate

mutual relations between Green’s functions and the

focusing functions introduced in the previous section.

First we introduce the decomposed Green’s functions

G

−,+

(x

R

, x

A

, t) and G

−,−

(x

R

, x

A

, t), with x

A

at

S

A

inside the medium and x

R

at the acquisition boundary

S

0 , see Figure 2a. Unlike the focusing functions, the

Green’s functions are defined in the actual medium,

which in general is inhomogeneous also below

_S

A

.

Fol-lowing common conventions, the second coordinate

vector (here x

A

) denotes the position of the impulsive

source and the first coordinate vector (here x

R

) that

of the receiver. In the same order, the superscripts

de-note the propagation directions at the source and

re-ceiver. Since the half-space above

_S

0 is homogeneous,

(5)

0

500 1000

1500

-2000

-1000

0 1000

2000

1900

2000

2100

2200

2300

2400

2500

De

pt

h

(k

m

)

Lateral position (km)

Ve

lo

ci

ty

(m

/s

)

0

0.5

1

1.5 -2 -1 0 1 2

a)

S

0

<latexit sha1_base64="j2hYY8psKNxn6vuJrxsTeRUAE0Q=">AAAB83icbVDLSgMxFL1TX7W+qi7dBIvgqsxUQRcuCm5cVrQP6Awlk2ba0EwmJBmhDP0NNy4UcevPuPNvzLSz0NYDgcM593JPTig508Z1v53S2vrG5lZ5u7Kzu7d/UD086ugkVYS2ScIT1QuxppwJ2jbMcNqTiuI45LQbTm5zv/tElWaJeDRTSYMYjwSLGMHGSr4fYzMOw+xhNnAH1Zpbd+dAq8QrSA0KtAbVL3+YkDSmwhCOte57rjRBhpVhhNNZxU81lZhM8Ij2LRU4pjrI5pln6MwqQxQlyj5h0Fz9vZHhWOtpHNrJPKNe9nLxP6+fmug6yJiQqaGCLA5FKUcmQXkBaMgUJYZPLcFEMZsVkTFWmBhbU8WW4C1/eZV0GnXvot64v6w1b4o6ynACp3AOHlxBE+6gBW0gIOEZXuHNSZ0X5935WIyWnGLnGP7A+fwB6aaRlA==</latexit>

S

A

<latexit sha1_base64="EHDi2TeBzlI2aXOJ77DEPyjNG5w=">AAAB83icbVDLSgMxFL3xWeur6tJNsAiuykwVdOGi4sZlRfuAzlAyaaYNzWSGJCOUob/hxoUibv0Zd/6NmXYW2nogcDjnXu7JCRLBtXGcb7Syura+sVnaKm/v7O7tVw4O2zpOFWUtGotYdQOimeCStQw3gnUTxUgUCNYJxre533liSvNYPppJwvyIDCUPOSXGSp4XETMKguxh2r/pV6pOzZkBLxO3IFUo0OxXvrxBTNOISUMF0brnOonxM6IMp4JNy16qWULomAxZz1JJIqb9bJZ5ik+tMsBhrOyTBs/U3xsZibSeRIGdzDPqRS8X//N6qQmv/IzLJDVM0vmhMBXYxDgvAA+4YtSIiSWEKm6zYjoiilBjayrbEtzFLy+Tdr3mntfq9xfVxnVRRwmO4QTOwIVLaMAdNKEFFBJ4hld4Qyl6Qe/oYz66goqdI/gD9PkDA3mRpQ==</latexit>

Green rays : G

<latexit sha1_base64="mdDIHTem19LdKXEvM3zwCN2J9fg=">AAACIHicbZDLSgMxFIYzXmu9jbp0EyyC4lhmVGhxVXGhSxWrQqeWTHqmDWYyQ5IRyzCP4sZXceNCEd3p05jWCt5OCHz8/zkk5w8SzpR23TdrZHRsfGKyMFWcnpmdm7cXFs9UnEoKdRrzWF4ERAFnAuqaaQ4XiQQSBRzOg6v9vn9+DVKxWJzqXgLNiHQECxkl2kgtu5L5MsIHEkD4Dpakp3Z9p3/yg8ts09nI1zI/CPFN3jpxvmjP0estu+SW3UHhv+ANoYSGddSyX/12TNMIhKacKNXw3EQ3MyI1oxzyop8qSAi9Ih1oGBQkAtXMBgvmeNUobRzG0lyh8UD9PpGRSKleFJjOiOiu+u31xf+8RqrDajNjIkk1CPr5UJhyrGPcTwu3mQSqec8AoZKZv2LaJZJQbTItmhC83yv/hbOtsrdd3jreKdWqwzgKaBmtoDXkoQqqoUN0hOqIolt0jx7Rk3VnPVjP1stn64g1nFlCP8p6/wBxVKE/</latexit>

,+

(x

R

, x

A

, t)

Yellow rays : G

,

_(x

R

, x

A

, t)

<latexit sha1_base64="j/0xQD8bmUoznDx+6hd2NCsKCKY=">AAACIXicbZDLSgMxFIYz9VbrbdSlm2ARKkzLTBUsrioudFnFXqRTSybNtKGZC0lGLcO8ihtfxY0LRboTX8b0ImjrCYGP/z+H5PxOyKiQpvmppRYWl5ZX0quZtfWNzS19e6cmgohjUsUBC3jDQYIw6pOqpJKRRsgJ8hxG6k7/fOTX7wkXNPBv5CAkLQ91fepSjKSS2noptrkHbwljwYNtQI4G4tQ2Rie5uIvzRj7JxbbjwsekfW380JkhD9t61iyY44LzYE0hC6ZVaetDuxPgyCO+xAwJ0bTMULZixCXFjCQZOxIkRLiPuqSp0EceEa14vGECD5TSgW7A1fUlHKu/J2LkCTHwHNXpIdkTs95I/M9rRtIttWLqh5EkPp485EYMygCO4oIdygmWbKAAYU7VXyHuIY6wVKFmVAjW7MrzUCsWrKNC8eo4Wy5N40iDPbAPcsACJ6AMLkEFVAEGT+AFvIF37Vl71T604aQ1pU1ndsGf0r6+AYk2odY=</latexit>

x

R

· · · ·x

R

<latexit sha1_base64="ORw1yd+Cs6mVB4Sf0I4SCESjw1A=">AAACN3icnVDLSgMxFM3UV62vUZdugkVwVWaqoCspuHElVewD2mHIZDJtaCYZkoxYhv6VG3/DnW5cKOLWPzBtB9HWlQcSDufce5N7goRRpR3nySosLC4trxRXS2vrG5tb9vZOU4lUYtLAggnZDpAijHLS0FQz0k4kQXHASCsYnI/91i2Rigp+o4cJ8WLU4zSiGGkj+fZl1g0ieDfyr7s4FPqf1/cQ3y47FWcCOE/cnJRBjrpvP3ZDgdOYcI0ZUqrjOon2MiQ1xYyMSt1UkQThAeqRjqEcxUR52WTvETwwSggjIc3hGk7Unx0ZipUaxoGpjJHuq1lvLP7ldVIdnXoZ5UmqCcfTh6KUQS3gOEQYUkmwZkNDEJbU/BXiPpIIaxN1yYTgzq48T5rVintUqV4dl2tneRxFsAf2wSFwwQmogQtQBw2AwT14Bq/gzXqwXqx362NaWrDynl3wC9bnF6Q9rsk=</latexit>

x

A

<latexit sha1_base64="vQbknEF6b5gSqaC6WPLILGHEO7s=">AAAB8HicbVBNSwMxEJ2tX7V+VT16CRbBU9mtgp6k4sVjBfsh7VKyabYNTbJLkhXL0l/hxYMiXv053vw3pts9aOuDgcd7M8zMC2LOtHHdb6ewsrq2vlHcLG1t7+zulfcPWjpKFKFNEvFIdQKsKWeSNg0znHZiRbEIOG0H45uZ336kSrNI3ptJTH2Bh5KFjGBjpYe0F4Toadq/7pcrbtXNgJaJl5MK5Gj0y1+9QUQSQaUhHGvd9dzY+ClWhhFOp6VeommMyRgPaddSiQXVfpodPEUnVhmgMFK2pEGZ+nsixULriQhsp8BmpBe9mfif101MeOmnTMaJoZLMF4UJRyZCs+/RgClKDJ9Ygoli9lZERlhhYmxGJRuCt/jyMmnVqt5ZtXZ3Xqlf5XEU4QiO4RQ8uIA63EIDmkBAwDO8wpujnBfn3fmYtxacfOYQ/sD5/AFxXJAo</latexit>

T

d

(x

R

, x

A

, t)

<latexit sha1_base64="B81QHkGWUAl2ffW6tDqZsRI4rCE=">AAACCXicbZDLSsNAFIYn9VbrLerSzWARKpSSVMEuK25cVukNmhAmk0k7dHJhZiKWkK0bX8WNC0Xc+gbufBunbQRt/WHg4z/ncOb8bsyokIbxpRVWVtfWN4qbpa3tnd09ff+gK6KEY9LBEYt430WCMBqSjqSSkX7MCQpcRnru+Gpa790RLmgUtuUkJnaAhiH1KUZSWY4O205q8QB6WSW1XB/eZ85t9Ycuq/LU0ctGzZgJLoOZQxnkajn6p+VFOAlIKDFDQgxMI5Z2irikmJGsZCWCxAiP0ZAMFIYoIMJOZ5dk8EQ5HvQjrl4o4cz9PZGiQIhJ4KrOAMmRWKxNzf9qg0T6DTulYZxIEuL5Ij9hUEZwGgv0KCdYsokChDlVf4V4hDjCUoVXUiGYiycvQ7deM89q9ZvzcrORx1EER+AYVIAJLkATXIMW6AAMHsATeAGv2qP2rL1p7/PWgpbPHII/0j6+AeZImSU=</latexit>

-1.5

-1.0

-0.5

0

0.5

1.0

1.5 -2000

-1000

0 1000

2000

-1.5

-1.0

-0.5

0

0.5

1.0

1.5 -2000

-1000

0 1000

2000

b)

Lateral position (m)

c)

Lateral position (m)

Ti

m

e

(s

)

Ti

m

e

(s

)

G

,+

_(x

R

, x

A

, t)

<latexit sha1_base64="Wq0k163Ejq9svJK/QEC05wAMLSU=">AAACB3icbZDLSsNAFIYn9VbrLepSkMEiVKwlqYJdVlzosoq9QBvDZDpph04uzEzEErJz46u4caGIW1/BnW/jtI2grT8MfPznHM6c3wkZFdIwvrTM3PzC4lJ2Obeyura+oW9uNUQQcUzqOGABbzlIEEZ9UpdUMtIKOUGew0jTGZyP6s07wgUN/Bs5DInloZ5PXYqRVJat717cxkfFw6QQdxwX3if2dfGHzorywNbzRskYC86CmUIepKrZ+menG+DII77EDAnRNo1QWjHikmJGklwnEiREeIB6pK3QRx4RVjy+I4H7yulCN+Dq+RKO3d8TMfKEGHqO6vSQ7Ivp2sj8r9aOpFuxYuqHkSQ+nixyIwZlAEehwC7lBEs2VIAwp+qvEPcRR1iq6HIqBHP65FlolEvmcal8dZKvVtI4smAH7IECMMEpqIJLUAN1gMEDeAIv4FV71J61N+190prR0plt8EfaxzeEV5fI</latexit>

G

,

(x

R

, x

A

, t)

<latexit sha1_base64="Pm3KQGY+v7dtTd54HQBTROAmBhg=">AAACCHicbZDLSsNAFIYn9VbrLerShYNFqNCWpAp2WXGhyyr2Am0Mk+mkHTq5MDMRS8jSja/ixoUibn0Ed76N0zaCtv4w8PGfczhzfidkVEjD+NIyC4tLyyvZ1dza+sbmlr690xRBxDFp4IAFvO0gQRj1SUNSyUg75AR5DiMtZ3g+rrfuCBc08G/kKCSWh/o+dSlGUlm2vn9xG5eKpaQQdx0X3if2dfGHzooleWTreaNsTATnwUwhD1LVbf2z2wtw5BFfYoaE6JhGKK0YcUkxI0muGwkSIjxEfdJR6COPCCueHJLAQ+X0oBtw9XwJJ+7viRh5Qow8R3V6SA7EbG1s/lfrRNKtWjH1w0gSH08XuRGDMoDjVGCPcoIlGylAmFP1V4gHiCMsVXY5FYI5e/I8NCtl87hcuTrJ16ppHFmwBw5AAZjgFNTAJaiDBsDgATyBF/CqPWrP2pv2Pm3NaOnMLvgj7eMb+raYAQ==</latexit>

f

<latexit sha1_base64="JMwxplcAMpvt0bK33MOGkNTO2fM=">AAACBXicbZDLSsNAFIYn9VbrLepSF4NFqFBLUgW7rLhxWcVeoI1hMp20QyeTMDMRS+jGja/ixoUibn0Hd76N0zaCtv4w8PGfczhzfi9iVCrL+jIyC4tLyyvZ1dza+sbmlrm905BhLDCp45CFouUhSRjlpK6oYqQVCYICj5GmN7gY15t3REga8hs1jIgToB6nPsVIacs1933Xvj0uJB3Ph/cj97r4Q+dFdeSaeatkTQTnwU4hD1LVXPOz0w1xHBCuMENStm0rUk6ChKKYkVGuE0sSITxAPdLWyFFApJNMrhjBQ+10oR8K/biCE/f3RIICKYeBpzsDpPpytjY2/6u1Y+VXnITyKFaE4+kiP2ZQhXAcCexSQbBiQw0IC6r/CnEfCYSVDi6nQ7BnT56HRrlkn5TKV6f5aiWNIwv2wAEoABucgSq4BDVQBxg8gCfwAl6NR+PZeDPep60ZI53ZBX9kfHwDNfOXFA==</latexit>

1 (x

R

, x

A

, t)

f

<latexit sha1_base64="1k6qDPYZQK6RfuKx35q4xXlnKEo=">AAACBXicbZDLSsNAFIYn9VbrLepSF4NFqFhKUgW7rLhxWcVeoI1hMp20QyeTMDMRS+jGja/ixoUibn0Hd76N0zaCtv4w8PGfczhzfi9iVCrL+jIyC4tLyyvZ1dza+sbmlrm905BhLDCp45CFouUhSRjlpK6oYqQVCYICj5GmN7gY15t3REga8hs1jIgToB6nPsVIacs1933Xvj0uJB3Ph/cj97r4Q+dFdeSaeatkTQTnwU4hD1LVXPOz0w1xHBCuMENStm0rUk6ChKKYkVGuE0sSITxAPdLWyFFApJNMrhjBQ+10oR8K/biCE/f3RIICKYeBpzsDpPpytjY2/6u1Y+VXnITyKFaE4+kiP2ZQhXAcCexSQbBiQw0IC6r/CnEfCYSVDi6nQ7BnT56HRrlkn5TKV6f5aiWNIwv2wAEoABucgSq4BDVQBxg8gCfwAl6NR+PZeDPep60ZI53ZBX9kfHwDMr2XEg==</latexit>

1 +

(x

R

, x

A

, t)

f

+

1,d

(x

R

, x

A

, t)

<latexit sha1_base64="QK8J//WPi2e62xFTqWy3n7xsaq0=">AAACD3icbZDLSsNAFIYn9VbrLerSzWBRKpaS1IIuK25cVrEXaGKYTCft0MmFmYlYQt7Aja/ixoUibt26822cthG09YeBj/+cw5nzuxGjQhrGl5ZbWFxaXsmvFtbWNza39O2dlghjjkkThyzkHRcJwmhAmpJKRjoRJ8h3GWm7w4txvX1HuKBhcCNHEbF91A+oRzGSynL0Q89JzHJicR/20vT2uJRYrgfvU+e6/EPnZXnk6EWjYkwE58HMoAgyNRz90+qFOPZJIDFDQnRNI5J2grikmJG0YMWCRAgPUZ90FQbIJ8JOJvek8EA5PeiFXL1Awon7eyJBvhAj31WdPpIDMVsbm//VurH0zuyEBlEsSYCni7yYQRnCcTiwRznBko0UIMyp+ivEA8QRlirCggrBnD15HlrVinlSqV7VivVaFkce7IF9UAImOAV1cAkaoAkweABP4AW8ao/as/amvU9bc1o2swv+SPv4BvmGm00=</latexit>

t

b

= t

d

✏

<latexit sha1_base64="QIbR3l0NPdUTUw15iZE1+P0oBH0=">AAAB/3icbVDLSsNAFJ3UV62vquDGzWAR3FiSKtiNUHDjsoJ9QBPCZDJph04ezNwIJXbhr7hxoYhbf8Odf+OkzUJbD1w4nHMv997jJYIrMM1vo7Syura+Ud6sbG3v7O5V9w+6Kk4lZR0ai1j2PaKY4BHrAAfB+olkJPQE63njm9zvPTCpeBzdwyRhTkiGEQ84JaAlt3oErncNbmbLEPvTc5sliovcqJl1cwa8TKyC1FCBtlv9sv2YpiGLgAqi1MAyE3AyIoFTwaYVO1UsIXRMhmygaURCppxsdv8Un2rFx0EsdUWAZ+rviYyESk1CT3eGBEZq0cvF/7xBCkHTyXiUpMAiOl8UpAJDjPMwsM8loyAmmhAqub4V0xGRhIKOrKJDsBZfXibdRt26qDfuLmutZhFHGR2jE3SGLHSFWugWtVEHUfSIntErejOejBfj3fiYt5aMYuYQ/YHx+QPos5YA</latexit>

t

a

=

t

d

+ ✏

<latexit sha1_base64="25s/3teoxQAJ8+ean7oB4LzTyg0=">AAACAHicbVDLSsNAFJ3UV62vqgsXbgaLIIglqYLdCAU3LivYBzQhTCaTdujkwcyNUEI2/oobF4q49TPc+TdO2i609cCFwzn3cu89XiK4AtP8Nkorq2vrG+XNytb2zu5edf+gq+JUUtahsYhl3yOKCR6xDnAQrJ9IRkJPsJ43vi383iOTisfRA0wS5oRkGPGAUwJacqtH4JKbC3AzW4bYz89tliguCqdm1s0p8DKx5qSG5mi71S/bj2kasgioIEoNLDMBJyMSOBUsr9ipYgmhYzJkA00jEjLlZNMHcnyqFR8HsdQVAZ6qvycyEio1CT3dGRIYqUWvEP/zBikETSfjUZICi+hsUZAKDDEu0sA+l4yCmGhCqOT6VkxHRBIKOrOKDsFafHmZdBt167LeuL+qtZrzOMroGJ2gM2Sha9RCd6iNOoiiHD2jV/RmPBkvxrvxMWstGfOZQ/QHxucPUg+WNA==</latexit>

FIG. 2 (a) Green’s functions in the actual medium. (b)

The left-hand side of equation 4 (fixed x

A

, variable x

R

).

(c) The left-hand side of equation 5.

For a source at x

S

at the acquisition boundary

S

0 ,

the Green’s function G

−,+

_(x

R

, x

S

, t) is by definition

the reflection response of the medium. We denote this

by R(x

R

, x

S

, t). The following relations hold between

the power-flux normalized Green’s functions and

fo-cusing functions (Slob et al., 2014; Wapenaar et al.,

2014)

G

−,+

(x

R

, x

A

, t) + f

1 −

(x

R

, x

A

, t) =

(4)

Z

S

0 dx

S

Z

∞

0 R(x

R

, x

S

, t

0 )f

1 +

(x

S

, x

A

, t

− t

0 )dt

0 and

G

−,−

(x

R

, x

A

,

−t) + f

1 +

(x

R

, x

A

, t) =

(5)

Z

S

0 dx

S

Z

0 −∞

R(x

R

, x

S

,

−t

0 )f

1 −

(x

S

, x

A

, t

− t

0 )dt

0 .

The time integration boundaries acknowledge the fact

that the reflection response R(x

R

, x

S

, t) is a causal

function of time, i.e., R(x

R

, x

S

, t < 0) = 0. Both

equations account for internal multiple scattering.

Equation 4 is exact, whereas in equation 5

evanes-cent waves are neglected. The interpretation of

equa-tion 4 is as follows. The right-hand side quantifies

the reflection response of the actual medium to the

downgoing focusing function f

1 +

(x

S

, x

A

, t). The

left-hand side shows that this reflection response consists

of the upgoing focusing function f

₁

−

(x

R

, x

A

, t) (the

blue rays arriving at

S

0 in Figure 1a) and the Green’s

function G

−,+

(x

R

, x

A

, t). The latter can be

under-stood as follows. In Figure 1a it can be seen that

the focal point x

A

acts as a virtual source at t = 0

for downgoing waves. Figure 2a shows that the

re-sponse to this virtual source is the Green’s function

G

−,+

(x

R

, x

A

, t) (the green rays in this figure).

Equa-tion 5 is interpreted in a similar way. The right-hand

side quantifies the reflection response of the

time-reversed actual medium to the upgoing focusing

func-tion f

₁

−

(x

S

, x

A

, t). The left-hand side shows that this

reflection response consists of the downgoing focusing

function f

1 +

(x

R

, x

A

, t) and the time-reversed Green’s

function G

−,−

_(x

R

, x

A

,

−t). The functions on the

left-hand sides of equations 4 and 5, convolved with a

wavelet, are shown in gray-level display in Figures 2b

and 2c, respectively.

Recall that we assume that

S

0 is a transparent

boundary and that the half-space above

S

0 is

homoge-neous. Hence, the reflection response R in equations

4 and 5 contains no surface-related multiples, which

complies with the situation after surface-related

mul-tiple elimination (Verschuur et al., 1992; van

Groen-estijn and Verschuur, 2010). Alternatively, equations

4 and 5 can be modified to account for surface-related

multiples in R (see Ware and Aki (1969) for the 1D

situation and Singh et al. (2017) and Dukalski and de

Vos (2018) for 2D and 3D situations). A further

dis-cussion on the inclusion of surface-related multiples in

the representations is beyond the scope of this paper.

Retrieval of focusing functions

Assuming R is known, equations 4 and 5 form a

system of two equations for four unknowns (f

1 +

, f

1 −

,

G

−,+

_{and G}

−,−

_{). An inspection of Figures 2b and}

2c reveals that the Green’s functions reside in other

time intervals than the focusing functions. We discuss

window functions to suppress the Green’s functions

from equations 4 and 5, so that we are left with a

system of two equations for two unknowns. To define

a time window for equation 4, we need to know the

first possible arrival of G

−,+

_(x

R

, x

A

, t). This would

occur when there would be a reflector just below

_S

A

.

For the first possible arrival of this Green’s function

we write

{G

−,+

(x

R

, x

A

, t)

}

‘first’

∝ T

d

(x

R

, x

A

, t),

(6)

where ‘first’ stands for ‘first possible’, T

d

is the direct

arrival of the transmission response of the medium

between

_S

A

and

S

0 (which is also the direct arrival of

(6)

to’. Note that we ignored the reflection coefficient

of the hypothetical reflector (this is justified since

we only use equation 6 to derive the time window).

We denote the arrival time of the direct transmission

response as t

d

(x

R

, x

A

). Hence, G

−,+

(x

R

, x

A

, t) can

be suppressed from equation 4 by applying a time

window that removes everything beyond t = t

b

=

t

d

(x

R

, x

A

)

− (the dashed line in Figure 2b). Here

is a small positive time constant (typically half the

duration of a wavelet), to account for the fact that

in practice all terms in equations 4 and 5 are

band-limited.

To define a time window for equation 5,

we need to know the last arrival of the time-reversed

Green’s function G

−,−

_(x

R

, x

A

,

−t). This is given by

the time-reversed direct arrival, hence

{G

−,−

(x

R

, x

A

,

−t)}

last

=

−T

d

(x

R

, x

A

,

−t)

(7)

(the factor

_{−1 follows from the sign-convention for the}

source for upgoing waves in G

−,−

(Wapenaar, 1996),

but this sign is irrelevant for the derivation of the time

window). Hence, G

−,−

_(x

R

, x

A

,

−t) can be suppressed

from equation 5 by applying a time window that

re-moves everything before t = t

a

=

−t

d

(x

R

, x

A

)+ (the

dashed line in Figure 2c). Note that t

a

=

−t

b

. Based

on this analysis, we define two time windows as

Θ

a

(x

R

, x

A

, t) = θ(t

− t

a

) = θ(t + t

d

− ),

(8)

Θ

b

(x

R

, x

A

, t) = θ(t

b

− t) = θ(t

d

− − t),

(9)

where θ(t) is the Heaviside step function (or, in

prac-tice, a tapered version of the Heaviside step function).

These windows suppress the Green’s functions and

pass the focusing functions f

₁

+

and f

₁

−

, except the

direct arrival f

_1,d

+

(x

R

, x

A

, t), which coincides with the

last arrival of G

−,−

_(x

R

, x

A

,

−t), see Figure 2c.

A few words of caution are needed here. First,

equa-tion 6 is only correct for limited offsets: at large

off-sets, refracted waves in G

−,+

may arrive earlier than

T

d

. Second, in a laterally varying strongly scattering

medium, diffraction events in the focusing functions

may be unintentionally suppressed by the time

win-dows. Third, in practice the inherent band-limitation

may cause partial interference of focusing functions

and Green’s functions, particularly in the case of thin

layers (i.e., thin compared to the wavelength). In this

paper we assume that offsets are limited, lateral

vari-ations are mild, and layers are not thin. Dukalski et

al. (2019) discuss how to account for thin layering,

assuming the medium is horizontally layered.

Application of the window Θ

b

(x

R

, x

A

, t) to both

sides of equation 4 gives

f

1 −

(x

R

, x

A

, t) =

(10)

Θ

b

Z

S

0 dx

S

Z

∞

0 R(x

R

, x

S

, t

0 )f

1 +

(x

S

, x

A

, t

− t

0 )dt

0 .

Similarly, applying Θ

a

(x

R

, x

A

, t) to both sides of

equation 5 we obtain

f

1 +

(x

R

, x

A

, t)

− f

_1,d

+

(x

R

, x

A

, t) =

(11)

Θ

a

Z

S

0 dx

S

Z

0 −∞

R(x

R

, x

S

,

−t

0 )f

1 −

(x

S

, x

A

, t

− t

0 )dt

0 .

The term

_−f

_1,d

+

on the left-hand side of equation 11

accounts for the fact that f

_1,d

+

is not passed by the

window, see Figure 2c. Equations 10 and 11 form

a coupled system of Marchenko equations. We show

how f

1 +

and f

1 −

can be retrieved, assuming R and f

1,d

+

are known. We adopt the compact operator notation

introduced by van der Neut et al. (2015b). In this

notation, equations 10 and 11 read

f

₁

−

= Θ

b

Rf

1 +

,

(12)

f

1 +

= Θ

a

R

?

f

1 −

+ f

1,d

+

,

(13)

with superscript ? denoting time-reversal. For

sim-plicity we use the same fonts for operators as for

wave fields. Operations like Rf

₁

+

stand for a

(multidi-mensional) convolution process (see right-hand side of

equation 10), whereas operations containing a

time-reversal, like R

?

_f

−

1 , stand for a correlation process

(see right-hand side of equation 11). Window

func-tions are always applied in a multiplicative sense.

Sub-stitution of equation 12 into equation 13 gives

f

1 +

= Θ

a

R

?

Θ

b

Rf

1 +

+ f

1,d

+

.

(14)

The product notation Θ

a

R

?

Θ

b

Rf

1 +

should be

under-stood in the sense that operators and window

func-tions act on all terms to the right of it, hence it stands

for Θ

a

(R

?

(Θ

b

(Rf

1 +

))). For notational convenience we

will not use the brackets. We rewrite equation 14 as

{δ − Θ

a

R

?

Θ

b

R

}f

1 +

= f

+

1,d

,

(15)

where δ is the identity operator. This equation can

be solved by

f

₁

+

=

K

X

k=0

{Θ

a

R

?

Θ

b

R

}

k

f

1,d

+

,

(16)

where K is the number of iterations needed for the

scheme to converge with acceptable accuracy

(conver-gence is guaranteed for K

_{→ ∞ (Dukalski and de Vos,}

2018)). Other approaches to solve equation 15 are

pro-posed by van der Neut et al. (2015a), Dukalski and de

Vos (2018) and Becker et al. (2018). Once f

1 +

is found,

f

1 −

follows from equation 12. We call equation 16 the

Marchenko scheme. As input it requires the

reflec-tion response R(x

R

, x

S

, t) at the acquisition

bound-ary (i.e., the reflection measurements after

surface-related multiple elimination and deconvolution for the

wavelet) and the direct arrival f

_1,d

+

(x

R

, x

A

, t) of the

fo-cusing function. Analogous to equations 1 and 2, the

latter is related to the direct arrival of the

(7)

transmis-sion response via

δ(x

0 H,A

− x

H,A

)δ(t) =

(17)

Z

S

0 dx

R

Z

∞

0 T

d

(x

0 A

, x

R

, t

0 )f

1,d

+

(x

R

, x

A

, t

− t

0 )dt

0 ,

for x

A

and x

0 A

at

S

A

. Hence, f

1,d

+

(x

R

, x

A

, t) is the (in

practice band-limited) inverse of T

d

(x

A

, x

R

, t). When

a macro model of the medium between

_S

0 and

S

A

is available, T

d

can be derived from this model and

inverted to obtain f

_1,d

+

. For convenience, this

inver-sion is often approximated by time-reversal, according

to f

_1,d

+

(x

R

, x

A

, t)

≈ T

d

(x

A

, x

R

,

−t). The Marchenko

scheme appears to be quite robust with respect to

am-plitude and timing errors in the direct arrival of the

fo-cusing function (Broggini et al., 2014; Wapenaar et al.,

2014). Nevertheless, for horizontally layered media

the amplitude of the direct arrival can be corrected,

using the principle of energy conservation (Mildner et

al., 2019). For highly complex media it can be

ad-vantageous to account for wavefield complexity in the

initial estimate of the focusing function (Vasconcelos

et al., 2015; Vasconcelos and Sripanich, 2019). In the

section “Marchenko multiple elimination” we discuss

methods that are independent of the direct arrival of

the focusing function.

The essential expressions for the retrieval of the

fo-cusing functions, i.e., the Marchenko equations and

the Marchenko scheme in compact operator form, are

summarized in Box 1. This box also shows the main

expressions for the other methods discussed in the

cur-rent section “Marchenko redatuming and imaging”.

Retrieval of Green’s functions (source redatuming)

Once the focusing functions have been found, the

next step is the retrieval of the Green’s functions. We

define time windows Ψ

a

(x

R

, x

A

, t) and Ψ

b

(x

R

, x

A

, t)

via

Ψ

a,b

(x

R

, x

A

, t) = 1

− Θ

a,b

(x

R

, x

A

, t).

(18)

Note that the windows Ψ

a,b

(x

R

, x

A

, t) are

comple-mentary to Θ

a,b

(x

R

, x

A

, t), defined in equations 8 and

9. Hence, they pass the Green’s functions and

sup-press the focusing functions, except f

_1,d

+

, see Figure 2c.

Application of Ψ

b

(x

R

, x

A

, t) to both sides of equation

4 thus gives

G

−,+

(x

R

, x

A

, t) =

(19)

Ψ

b

Z

S

0 dx

S

Z

∞

0 R(x

R

, x

S

, t

0 )f

1 +

(x

S

, x

A

, t

− t

0 )dt

0 .

Similarly, applying Ψ

a

(x

R

, x

A

, t) to both sides of

equation 5 yields

G

−,−

(x

R

, x

A

,

−t) + f

1,d

+

(x

R

, x

A

, t) =

(20)

Ψ

a

Z

S

0 dx

S

Z

0 −∞

R(x

R

, x

S

,

−t

0 )f

1 −

(x

S

, x

A

, t

− t

0 )dt

0 .

We interpret these expressions as follows. The

reflec-tion response on the right-hand sides is the response

to an actual source at x

S

, observed by an actual

re-ceiver at x

R

, both at the acquisition surface

S

0 . This

is visualized in Figure 3a. The Green’s functions on

the left-hand sides are responses to a virtual source

for downgoing waves (equation 19) and upgoing waves

(equation 20) at x

A

in the subsurface, observed by the

actual receiver at x

R

at the surface. Hence, equations

19 and 20 accomplish source redatuming from x

S

at

the acquisition surface

S

0 to virtual-source position

x

A

in the subsurface. Figure 3b visualizes equation

19. Equation 19 and Figure 3b resemble the

virtual-source method proposed by Bakulin and Calvert

(2006), except that in their formulation the actual

re-ceivers are situated in a horizontal borehole and,

in-stead of using a Marchenko-derived focusing function,

they use a windowed time-reversed response between

sources at the surface and an actual receiver at x

A

in

the borehole. Hence, when measurements are

avail-able in a borehole, their method enavail-ables the retrieval

of the response to a virtual source at x

A

below a

com-plex overburden. Note, however, that their method

does not account for internal multiples.

In the compact operator notation, equations 19 and

20 for source redatuming become

G

−,+

_R,A

= Ψ

b

Rf

1 +

,

(21)

G

−,−?

_R,A

= Ψ

a

R

?

f

1 −

− f

+

1,d

,

(22)

where the subscripts R and A on the left-hand sides

refer to the actual receiver position x

R

at the surface

S

0 and the virtual-source position x

A

at the datum

plane

S

A

in the subsurface. In the following we discuss

different methods to redatum also the receivers from

the surface to a virtual-receiver position x

B

at

S

A

.

Receiver redatuming by MDD

We define the reflection response at datum plane

S

A

of the target below

S

A

as R

tar

(x

B

, x

A

, t) =

G

−,+

(x

B

, x

A

, t), with x

A

and x

B

both at

S

A

(sub-script ‘tar’ stands for ‘target’).

Note that, when

G

−,+

_(x

B

, x

A

, t) is defined in the actual medium, it

not only contains the response of the medium below

S

A

, but also multiples between reflectors below and

above

_S

A

(the dashed rays in Figure 4a). We define