Design, implementation and application of the Marchenko Plane-wave algorithm

The Marchenko algorithm can suppress the disturbing effects of internal multiples that are present in seismic reflection data. To achieve this, a set of coupled equations with four unknowns is solved. These coupled equations are separated into a set of two equations with two unknowns using a time window. The two unknown focusing functions can be resolved by an iterative or direct method. These focusing functions, when applied to reflection data, create virtual point-sources inside the medium. Combining individual virtual point-sources into a plane-wave leads to an efficient computation of images without internal multiples. In this study the internal multiples are eliminated in a redatuming step which is part of the imaging algorithm. To use the Marchenko algorithm with plane-wave focusing functions, the time window that separates the unknowns must be adapted. The design of the plane-wave Marchenko algorithm is explained and illustrated with numerically modeled and measured reflection data.


INTRODUCTION
Seismic imaging is a technique to image geological structures in the subsurface of the earth, from reflected wavefields measured at the surface of the earth.The measured wavefields are usually originating from human activated and controlled sources like air-guns or vibrating plates.In passive seismic methods the source of the wavefield can originate from earthquakes, ocean waves, or uncoordinated human activities like traffic.The primary reflection of a geological structure, large and strong enough to get detected by a propagating wavefield, is of main interest and used to compute an image of the subsurface.At each geological structure wavefields are partly reflected upward and partly transmitted further downward.Between two strong reflecting structures the wavefield can bounce up-and down multiple times and generate so called internal multiples.These multiple reflections are also measured by the geophones and hard to distinguish from primary reflections.In the imaging step the reflections are migrated from time to depth and construct an image of the subsurface.If multiple reflections are not recognized as such, they will get imaged as being a primary reflection at the wrong depth.Those falsely imaged multiples distort the actual image of the subsurface; the distorted image will contain much more (ghost) structures that are positioned along with the primary reflections.In seismic imaging it is therefore important to Marchenko plane-wave recognize these multiple reflections, and if possible remove them from the computed image.This removal can be done before or after the imaging step.If the removal is done before the imaging step the multiples are removed from the measured data.For removal after the imaging step a computed prediction of imaged multiples is subtracted from the image to obtain an image without multiples.In this paper we discuss a method that removes internal multiples from the measured reflection data.
Besides internal multiples that are reflected between boundaries within the subsurface, there are also free-surface related multiples.These multiples are generated by reflections from upcoming waves that are bounced back into the subsurface by the surface of the earth.Free-surface multiples are not considered here.They are assumed to be removed prior to the removal of the internal multiples.
The Marchenko algorithm can eliminate internal multiples from seismic reflection data (Slob et al., 2014;Behura et al., 2014).In this algorithm the up-and down-going focusing functions, with a focal-point in the subsurface, are key to the method.The goal of the Marchenko method is to retrieve these up-and down-going parts of the focusing functions from the reflection data, by solving a coupled set of so-called Marchenko equations.This set of equations can be solved by iterative methods (Wapenaar et al. (2014a), Thorbecke et al. (2017)), or a direct method (van der Neut et al. (2015a), Ravasi (2017)).The Marchenko method has found many different applications that range from redatuming for time-lapse monitoring (van IJsseldijk et al., 2022), to adaptive subtraction of Marchenko estimated multiples (Staring et al., 2018), to homogeneous Green's function retrieval (Brackenhoff et al., 2019).In this paper a particular efficient application of the Marchenko method, using plane-waves, is highlighted and the implementation details of this method are discussed in more detail.Meles et al. (2018) show that besides focal-points, focal-planes can also be solved by the Marchenko equations.Meles et al. (2020) builds on the Marchenko Multiple Elimination (MME) method proposed by Zhang and Slob (2019) and introduces the plane-wave MME method.The big advantage of the plane-wave based Marchenko method is that with a minimal effort of one Marchenko run (for a single plane-wave), for each depth level (or time instant for MME), a multiple free image can be built up.Especially in 3D applications, the plane-wave Marchenko method is a computational efficient method to build an internal multiple-free image (Brackenhoff et al., 2022).One plane-wave can be sufficient to build an accurate image if the subsurface interfaces are near-flat.Multiple plane-waves, with different illumination angles, are needed for subsurface interfaces with varying dip, or geologically complex structures, to illuminate the subsurface properly (Almobarak, 2021).The Marchenko algorithm for the focal plane-wave method is similar to the focal-point algorithm.The main difference lies in the choice of the time-window to separate the Green's function from the focusing functions.In this paper we will discuss in detail how these timewindows are adapted for the Marchenko plane-wave method, discuss the implementation aspects and illustrate the method with applications on numerically modeled and field data.
The software accompanied with this paper contains scripts and source code to reproduce all the numerical examples presented in this paper.The code can also be found at its GitHub repository (Thorbecke et al., 2017;Thorbecke and Brackenhoff, 2022), where the most recent updated version and latest developments are available.To reproduce the figures and to carry out a few post-and pre-processing steps, Seismic Unix (Cohen and Stockwell, 2016) is required.

THEORY
The Marchenko method is introduced by two coupled equations that contain four unknown fields (up-and downgoing focusing functions and Green's functions) we would like to retrieve.These fields enable us to suppress internal multiples in seismic reflection data.The seismic reflection experiment is recorded at acquisition boundary S 0 .The reflection response R(x R , x S , t), a scaled version of the Green's function without the direct arrival (Wapenaar et al., 2012), is measured with sources and receivers positioned at x S and x R on this boundary.The recording time is denoted by t.This reflection data does not contain free-surface related multiple reflections neither a source wavelet.Pre-processing steps are required to remove free-surface multiples, the direct arrival, and deconvolve the wavelet from the measured reflection data.
The up-and downgoing parts of the focusing functions f − 1 and f + 1 , are used to define a relation between the decomposed Green's functions G −,+ and G −,− in the actual medium with the reflection response at the surface (Wapenaar et al., 2014b).The focusing functions have a focal point in the subsurface at x A .This focal point serves as a virtual source for the Green's function.The most rightward + and − superscripts of the decomposed Green's function refer to the direction of propagation (+ for down and − for up) from the virtual source at x A .The most left superscript indicate an up-(−) or downward(+) propagating field at the receiver locations.The relation between the two unknown focusing functions and the two unknown Green's function is given by the following two equations (Wapenaar et al., 2021); In the compact operator notation of van der Neut et al. (2015b) equations 1 and 2 are written as where ⋆ denotes the time-reverse.These two equations contain four unknowns: two Green's functions and two focusing functions.The only known in these equations is the measured reflection data R. Wapenaar et al. (2013) use the reasoning that the Green's function and focusing function can, under certain circumstances, be separated in time.Therefore, a time-window function (Wapenaar et al., 2014b) Θ(t) is defined that passes the focusing function and removes the Green's function from the left hand-side of equations 3 and 4.This function removes all events that arrive at later times than the direct wave, traveling from the virtual source position x A to the receiver position x R at surface S 0 , including the direct wave itself.This results in the following two equations that only have two unknowns (assuming the direct arrival f + 1,d is known) is the direct arrival of f + 1 , and f + 1,m the events that arrive before the direct arrival time t d .These time windows (Wapenaar et al., 2021) are defined as where θ(t) is a tapered Heaviside step function.Note that two time-windows are defined, one at time t a and one at t b .In the algorithm t b = t d − ε = t a and makes equation 8 equal to equation 7, with window function θ(t d − ε − t).The ε takes into account the finite length of the band-limited wavelet and ensures that the direct wave is removed from the right hand-side of equation 6 (Broggini et al., 2014).
To illustrate the application of the time-windows a virtual point source is defined at 800 meter depth in the laterally varying model of Figure 1 (after Meles et al., 2018).Figure 2 shows the focusing functions and Green's functions for the model of Figure 1, and the window functions (indicated with a dashed line) that separate these functions.Figure 2a represents the left-hand side of equation 3 and the time-window that separates f − 1 from G −,+ .In Figure 2b the time-window separates the time-reversal of f + 1 from G −,− and represents the left-hand side of equation 4. The correlation/convolution in equations 3 and 4 are carried out in the frequency domain and a discrete Fourier transform in time is used to transform the sampled reflection data to the frequency domain.The discrete Fourier transform makes the fields periodic in time with a periodicity equal to the number of time samples (n t ).Given this periodicity in time, reflection events can shift in time beyond n t and will end-up in negative times.The time windows defined in equations 7 and 8 pass all events earlier in time than t = t a and will also pass these time wrap-around events.To exclude these wrap-around events a time window is implemented for negative times as well.
In Figure 2a we can see that the focusing functions also include events at negative times and these events will not be present earlier than −t a = −t d + ε.Hence the cut-off point of the time window at negative times is chosen at −t a .The implemented time-windows Marchenko plane-wave become and the time-windows at negative times, to suppress time wrap-around, are indicated with the dotted lines in Figure 2.   To solve the unknown focusing functions in the coupled equations 5 and 6 different methods are developed.The iterative method described in Behura et al. (2014) and Wapenaar et al. (2014b) start with f + 1,d as the initial solution of f + 1 and solve equation 5 for f − 1 and substitute the solution in equation 6 to update f + 1 .This process is repeated until the updates to the focusing functions become very small.This iterative algorithm to solve the Marchenko equations is shown in Figure 3.The even iterations (starting the iteration count at 0 for the initial solution) in the scheme solve equation 5 and the odd iterations solve equation 6. Thorbecke et al. (2017) explains in more detail the implementation of this iterative algorithm.
Depending on the application it is not always needed to solve these equations until convergence.In Staring et al. (2018) the results of the first iteration are used to predict internal multiples and an adaptive subtraction method is used to suppress the predicted multiples.A direct least-squares inversion method to solve the coupled equations is discussed in van der Neut et al. (2015b) and Ravasi (2017).Marchenko plane-wave ) and up-going(f − 1 ) focusing functions are alternately updated.The scheme is finished after a pre-defined number of iterations and is usually chosen between 10-20 iterations.Meles et al. (2018) show that plane-wave focusing functions f + 1 and f − 1 and associated plane-wave Green's functions G−,+ and G−,− can be obtained by integrating an appropriate set of time-delayed focusing functions f + 1 and f − 1 , each involving the solution of a Marchenko equation.The plane-wave focusing functions (Wapenaar et al., 2021;Brackenhoff et al., 2022) are defined by the following integration with p = (p 1 , p 2 ) and p 1 and p 2 horizontal ray parameters.The surface S A is the depth level at which focusing takes place.The plane-wave Green's functions are defined by a similar integration.Here ), and (x 1,c , x 2,c ) is the rotation point of the plane-wave and p A = (p, x 3,A ) the ray parameter of the plane-wave at surface S A .The rotation point is chosen in the center of the lateral extent of the plane-wave.This rotation point also defines the time origin t = 0 of the plane-wave.By making this choice for t = 0 the time-axis in the computed plane-wave Green and focusing functions is the same as in the point source Marchenko solutions with a focal point at the rotation point of the plane-wave.
Note that the plane-wave integration in equation 11 for a time-reversed wave-field Marchenko plane-wave with p ′ A = (−p, x 3,A ), a plane-wave dipping with the opposite angle as a plane-wave with p.On a surface S A with a homogeneous velocity along the surface, p • x H,A are time shifts that are linearly proportional to the distance from the rotation point.
Applying the same integration as in equation 11 over all fields results in the plane-wave representations of equations 3 and 4 (Meles et al., 2018 Applying a time-window that separates the Green function from the focusing function gives again two equations with two unknowns (assuming , where f +⋆ 1,d is the direct arrival of the plane-wave with a propagation angle defined by (p A , x 3,A ), and f +⋆ 1,m contains the events that arrive before the direct arrival time td , where td is the first arrival time of a plane-wave with a propagation angle defined by (p A , x 3,A ).Note that the Green's functions in equations 13 and 14 contain plane-waves with opposite dipping angles defined by p A and p ′ A .To separate G −,+ (x R , p A , t) from f − 1 (x R , p A , t) the direct arrival td , that is the first arrival time of G −,+ (x R , p A , t), is required.
The notation for t d or t ′ d is a choice, our choice is the same as made in Wapenaar et al. (2021).This choice uses the ′ on t d the same as the ′ on the propagation angle p A in the upgoing Green's functions G −,+ (x R , p A , t) and G −,− (x R , p ′ A , t).Note that the defined arrival time td has the same arrival time of an upward propagating modeled plane at (x 3,A ) with propagation angle p ′ A .The time window Θb (t) removes G−,+ (x R , p A , t) from equation 13.The first non-zero ) from equation 14, hence all events at times later than t a = t′ d − ε are set to zero.Similar to the point-source scheme these time windows are implemented with an additional window at negative times to suppress time wrap-around and are given by Similar to Figure 2, Figure 4 is an example of the plane-wave focusing functions and Green's functions in the left hand-side of equations 13 and 14 and the window functions separating them.In this example the depth of the plane-waves is chosen at 800 m in the model shown in Figure 1.In the following section these two time-window functions are discussed in more detail and the differences with the point-source implementation of the Marchenko algorithm Marchenko plane-wave are explained.G−,+ (x R , p A , t) G−,− (x R , p ′ A , t)

BASIC ALGORITHM
The plane-wave method is illustrated with two numerical examples, a laterally invariant and laterally variant medium.In two dimensions the downward propagating plane-wave focusing function f + 1 (x, p A , t) in equation 11, defines a plane-wave at the focal plane S A with a dip angle α and p 1 = sin(α) c .For simplicity we assume a medium with a laterally invariant velocity c for each depth.Figure 6 shows f + 1 (x R , p A , t) for a focal plane at 900 meter depth at two different angles propagating through the model shown in Figure 5: Figure 6a with an angle of 0 degrees and Figure 6c with an angle of 3 degrees.In Figure 6b and 6d the focus function is shown for receivers at the focal level (900 m depth).Similar to the focal-point Marchenko method, the medium for the plane-wave focusing functions is chosen homogeneous below the focal level.The focus function f + 1 (x, p A , t) has a focus in time at t = p • x H,A .For a dipping plane-wave this time focus occurs at different positions in the model.Figure 6a and c shows the computed f + 1 (x R , p A , t) at the surface (0 m depth) and includes an extra event around -0.45 seconds that compensates for the multiples generated between the first and second reflector at 400 and 700 meter respectively (Zhang and Slob, 2019).The focus functions at focal level (Figure 6b and d) show only one event, the downgoing event present at the surface is compensated by the reflected event from the reflector at 700 meter depth.Marchenko plane-wave To illustrate this compensation effect, snapshots of the focusing function f + 1 (x, p A , t) (Figure 6a and 6c) propagating into the truncated medium (that is homogeneous below the last reflector at 700 m depth) are shown in Figure 7 for the same angles of 0 and 3 degrees.Figure 7a to 7d shows four different snapshots of the superposition of the down-going horizontal plane-wave f + 1 (x R , p A , t) and upgoing plane-wave f − 1 (x R , p A , t).The snapshots of a plane-wave with an angle of 3 degrees are shown in Figure 7e to 7h.At 0.05 seconds before t=0 (Figure 7a and 7e), there are two upward traveling reflected waves from the interfaces at 400 and 700 meter depth and two downgoing events from f + 1 (x R , p A , t).The next snapshots (Figure 7b and 7f at t = +0.05seconds) show that the second downward traveling event coincides at the first interface (400 m depth) with the upgoing reflection of the second interface (700 m depth) and these events compensate each other.This is indicated with an arrow in the pictures.The third and fourth snapshot show that after this compensation all the internal multiples between the reflectors at 400 and 700 m depth have vanished and only one downgoing direct wave and an upgoing reflected wavefield (from the reflector at 700 m depth) are remaining.The compensation of the first upward traveling multiple indicates that f + 1 (x R , p A , t) is a solution of the Marchenko equations (Zhang et al., 2019).The illustration in Figure 7 demonstrates that the same internal multiple compensation principle also holds for the plane-wave Marchenko method.
In Figure 8 the experiment is repeated in the laterally variant medium of Figure 1.The focal-plane is chosen at 800 m depth, just below the fourth reflector.The same observation is made; the downgoing events in f + 1 (x R , p A , t) compensate the upgoing events at interfaces and suppress the generation of internal multiples.Marchenko plane-wave

Horizontal plane-waves
To start the iterative Marchenko algorithm for plane-waves f + 1,d (x R , p A , t) is required; the first arrival of a plane-wave response from a focal-plane in the subsurface.This forward modeling step can be computed in a macro model estimated from the reflection data.The computed initial response is then muted below the first arrival times to get the time-reversal of the initial focusing field f + 1,d (x R , p A , t).A first example is made for a horizontal (zerodegree) plane-wave defined at 800 meter depth in the laterally varying model of Figure 1 (Meles et al., 2018).
Figure 9a shows the forward modeled plane-wave response of a horizontal plane-wave at 800 m depth and receivers at the surface.The first arrivals of that plane-wave, shown in Figure 9b, is the time-reversal of the input field f + 1,d (x R , p A , t) of the Marchenko scheme.The computed down-and up-going Green's functions, after 16 iterations of equations 15 and 16 and substituting the results in equations 13 and 14, are shown in 9c and 9d respectively and give the expected response.The horizontal plane-wave Marchenko mute-line (as defined by Θ′ a,b in equations 17 and 18), to separate the Green's function from the focusing functions, is symmetric around t = 0 because t d = t ′ d .This is the same time symmetry as in use for the Marchenko point-source algorithm.

Dipping plane-waves
The Marchenko algorithm for dipping plane-waves follows the same procedure as for horizontal plane-waves.As indicated by equations 17 and 18 the Marchenko time-windows have to be designed differently for dipping plane-waves.For horizontal plane-waves the implemented time-window is symmetric around t = 0 and Θ′ a = Θ′ b .This does not hold anymore for dipping plane-waves.Figure 10 is an example of time windows that are designed for a dipping plane-wave of +5 degrees.Two time-windows are designed: one for the even iterations (Figure 10a) in use by equation 15 and one for the odd iterations, that has a reverse angle (Figure 10b) and is used in equation 16.To design these windows for a positive angle the first arrival time for a plane-wave at the same depth level with a reverse angle, is needed as well.Marchenko plane-wave In the following we explain in more detail what is expressed in the plane-wave Marchenko equations 15 and 16.The plane-wave Marchenko scheme starts with forward modeling the response to a (dipping) plane-wave at depth.This modeled wave-field response is given in Figure 12a, where the depth of the plane-wave is chosen at 800 m.The direct arrival is selected from this field, it's time reverse will be f + 1,d (x R , p A , t), and is shown in Figure 12b, and together with the reflection data R, input of the Marchenko scheme.Marchenko plane-wave In the Marchenko algorithm the iterations are alternating between a convolution of the focusing functions with R or a correlation with R (Thorbecke et al., 2017).In the first step, correlation by equation 16, the wavefield is shifted backward in time related to the times of t d , and in the second step, convolution by equation 15, the wavefield is shifted forward in time related to the times of t ′ d .In Figure 11a the result of the first iteration (correlation) is shown and in Figure 11b the result of the second iteration (convolution) is shown.In Figure 11a we can see that the first event, that starts at negative time, shows an imprint of the undulation of the first reflector and has an opposite dip compared to the plane-wave in Figure 12a.
The result of the first iteration in Figure 11a is windowed in time (with the window in Figure 10a) to mute G−,+ (x R , p A , t), followed by time-reversal and convolution with R (equation 16).The result is shown in Figure 11b.In this second iteration the convolution step brought the data back to the arrival times corresponding to reflection times in the forward modeled plane-wave response of Figure 12a.Note that the arrival time in Figure 11b, starting from the left at time 0.2 s dipping to the right to 0.5 s, is the same as the first arrival time t d in Figure 12a and 12b.The result of this second iteration is muted in time (with the window in Figure 10b) to mute G−,− (x R , p ′ A , t) and the first arrival event at f + 1,d (x R , p A , t), as indicated in Figure 4.The odd iteration(s) are building up f − 1 (x R , p A , t) and G−,+ (x R , p A , t) and the even iteration(s) are building up f + 1 (x R , p A , t) and G−,− (x R , p ′ A , t). Figure 11c and 11d show f + 1 (x R , p A , t) and f − 1 (x R , p A , t) respectively after 16 iterations.
In Figure 12 three plane-wave responses are shown with angles of -5, 0 and 5 degrees after 16 iterations.Comparing these three plane-wave responses for G−,− and G−,+ shows that each angle illuminates different parts of the medium.This is clearly seen in the events that arrive later than 1.2 seconds.By combining different plane-wave responses into one image a fully illuminated image can be constructed by using only a few migrations (Meles et al., 2018).The migration strategy for plane-wave imaging without internal mul-Marchenko plane-wave tiples has the same possibilities as the with the point source Marchenko (double focusing van der Neut et al. (2017); Staring et al. (2018), MDD Ravasi et al. (2016)).Almobarak (2021) discusses the different plane-wave imaging methods and shows that the Marchenko Green's function plane-wave response is the basis for a computational efficient imaging method.
Figure 13 shows horizontal plane-wave images from the Troll field data-set located west of Norway and provided by Equinor.This data set is part of a time-lapse monitoring set.We have selected a small part of this data-set, that has source and receiver positions on the same locations.The receiver and source spacing is 12.5 meter and traces are recorded with a time sampling of 4 ms.The data-set has been pre-processed to remove free-surface multiples and deconvolve the wavelet (Qu and Verschuur, 2020).
The imaging is carried out according to the imaging method described in Meles et al. (2018).In the basic imaging method a forward modeled plane-wave response is computed at each depth level in an estimated smooth macro model of the data.This plane-wave depth response is correlated with the point-source responses of the recorded data and integrated over the receiver coordinate for each point-source response.This creates the plane-wave depth response of the data (Rietveld et al., 1992).This plane-wave response of the data is correlated with the same modeled plane-wave response at each depth level and the imaging condition t = 0 is used to construct the image for all depth levels.Combining these depth levels gives the left side picture in Figure 13 labeled 'standard'.
To compute the Marchenko based image, the first arrivals of the forward modeled planewave response at each depth level is input to the plane-wave Marchenko algorithm to create G−,+ .The computed G−,+ is, similar to the standard imaging method, correlated with the forward modeled plane-wave response for each depth level and the imaging condition at t = 0 constructs the image for all depth levels.The advantage of using G−,+ is that this field does not contain downgoing internal multiples from the layers above the focal/imaging depth.The middle picture of Figure 13 show the Marchenko-created image.The difference between the standard image and the Marchenko-based image is shown in Figure 13c.From this difference plot it is observed that the Marchenko method predicts and removes internal multiples.In this data set the effect of the multiple removal on the image is small.van IJsseldijk et al. (2022) shows that Marchenko multiple removal on this data set improves the confidence of the effects of small time-lapse changes.

CONCLUSIONS
The plane-wave Marchenko method is a straight-forward extension of the point-source Marchenko method.A counter-intuitive aspect of the plane-wave method is that the upand down-going Marchenko Green's functions have an opposite dipping angle.This is taken into account in the time-windows that separate the Green's function from the focusing function.In this paper, the use of these time-windows is illustrated with numerical and field data examples.The plane-wave Marchenko method can give a computational advantage over the point-source method.Specially for imaging applications with 3-dimensional data sets that has moderate lateral changes, in that case only a few plane-wave migrations are needed to compute a well illuminated image.

Figure 1 :
Figure 1: Multi layer model with velocity (a) and density (b) parameters.The location of the virtual point source is marked with a ⋆.

Figure 2 :
Figure 2: Illustration of the time-window function to separate the Green's function from the focusing function.The dashed black lines indicate the separation line of the timewindow and are indicated with white arrows.The dotted line indicates the time-window that suppresses time wrap-around.

Figure 3 :
Figure 3: Flow chart of the Marchenko algorithm.The down-(f + 1) and up-going(f − 1 ) focusing functions are alternately updated.The scheme is finished after a pre-defined number of iterations and is usually chosen between 10-20 iterations.

Figure 4 :
Figure 4: Illustration of the time-window function to separate the plane-wave Green's function from the focusing function.The dashed black lines indicate the separation line of the time-window and are indicated with white arrows.The dotted line indicates the timewindow that suppresses time wrap-around.
Figure 5: Two dimensional four layer model with velocity (a) and density (b) parameters.A common-source record, with source position x = (x 1 = 0, x 3 = 0) and receivers at x R = (x 1 = x R , x 3 = 0) (c).The source wavelet in R has a flat frequency spectrum from 5 to 90 Hz.

Figure 6 :Figure 7 :Figure 8 :
Figure 6: Time recordings of the plane-wave focusing function f + 1 (x R , p A , t) with a focal depth of 900 meter measured with receivers at x 3 = 0 and x 3 = 900 meter in the truncated medium for two different plane-wave propagation angles (0 and 3 degrees).At the endpoints of the plane-wave, diffraction curves are present due to the limited lateral extent of the constructed plane-wave.
Figure 9: Plane-wave Marchenko scheme for a horizontal plane-wave p A = 0.All figures are plotted with the same clipping factor.

Figure 10 :
Figure 10: The time windows for dipping plane-waves for even (a) and odd (b) iterations for an angle of +5 degrees.The wavefields in the black area of the windows pass, the fields in the white area are set zo zero.
Figure11: Basic plane-wave Marchenko results for a plane-wave with an angle of incidence of +5 degrees.Note, that the results of the first iteration (a) is dipping in the opposite direction as the second iteration (b) and the algorithm uses the time-windows, designed for dipping plane-waves as given in equations 17 and 18, to take this into account.

Figure 12 :Figure 13 :
Figure 12: Marchenko computed plane-wave responses for angles at 5 (a-d), 0 (e-h), and -5 (i-l) degrees.Note the difference in illumination in the decomposed Green's function for different dipping angles of the plane-wave.