-
PDF
- Split View
-
Views
-
Cite
Cite
Aydin Shoja, Joost van der Neut, Kees Wapenaar, Target-oriented least-squares reverse-time migration using Marchenko double-focusing: reducing the artefacts caused by overburden multiples, Geophysical Journal International, Volume 233, Issue 1, April 2023, Pages 13–32, https://doi.org/10.1093/gji/ggac438
- Share Icon Share
SUMMARY
Geophysicists have widely used least-squares reverse-time migration (LSRTM) to obtain high-resolution images of the subsurface. However, LSRTM is computationally expensive and it can suffer from multiple reflections. Recently, a target-oriented approach to LSRTM has been proposed, which focuses the wavefield above the target of interest. Remarkably, this approach can be helpful for imaging below complex overburdens and subsalt domains. Moreover, this approach can significantly reduce the computational burden of the problem by limiting the computational domain to a smaller area. Nevertheless, target-oriented LSRTM still needs an accurate velocity model of the overburden to focus the wavefield accurately and predict internal multiple reflections correctly. A viable alternative to an accurate velocity model for internal multiple prediction is Marchenko redatuming. This method is a novel data-driven method that can predict Green’s functions at any arbitrary depth, including all orders of multiples. The only requirement for this method is a smooth background velocity model of the overburden. Moreover, with Marchenko double-focusing, one can make virtual sources and receivers at a boundary above the target and bypass the overburden. This paper proposes a new algorithm for target-oriented LSRTM, which fits the Marchenko double-focused data with predicted data. The predicted data of the proposed method is modelled by a virtual source term created by Marchenko double-focusing on a boundary above the target of interest. This virtual source term includes all the interactions between the target and the overburden. Moreover, the Marchenko double-focused data and the virtual source term are free of multiples generated in the overburden. Consequently, our target-oriented LSRTM algorithm suppresses the multiples purely generated inside the overburden. Our algorithm correctly accounts for all orders of multiples caused by the interactions between the target and the overburden, resulting in a significant reduction of the artefacts caused by the overburden internal multiple reflections and improves amplitude recovery in the target image compared to conventional LSRTM.
1 INTRODUCTION
In seismic exploration, geophysicists estimate the Earth’s wave propagation parameters, such as the subsurface velocity and density. To evaluate these parameters, an array of sources and receivers is deployed at the surface of the Earth to generate and record seismic waves. A common approach to retrieve the subsurface parameters is to consider the model of the subsurface as a superposition of a background model (m0), which corresponds to the long wavelengths, and a short-wavelength reflectivity model (δm). This scale separation is based on a weak-scattering assumption (Schuster 2017). In seismic migration, we aim to obtain a structural image of δm.
Migration is an imaging process based on a linear relation between the recorded data and the reflectivity model (δm). One of the common migration techniques is called reverse time migration (RTM; Baysal et al. 1983; Zhang et al. 2018; Zhou et al. 2018). In RTM, the underlying linear relation is the Born integral, so the RTM image is produced by applying the adjoint of the Born operator to the observed data. However, this approach may suffer from many errors such as noise (Dutta et al. 2017), temporal band limitations, limited recording aperture (Tang 2009), spatial aliasing and multiple reflections (Liu et al. 2016).
One way to reduce some of these errors is to resolve δm by least-squares inversion. This approach, known as least-squares reverse-time migration (LSRTM), is the most common algorithm for producing a high-resolution reflectivity model of the subsurface. Nevertheless, LSRTM is computationally expensive (Dai et al. 2012) as it requires solving the wave equation and its adjoint for a large model in many iterations and storing large data volumes (Tang 2009; Herrmann & Li 2012; Farshad & Chauris 2021). Reducing the computational burden of LSRTM is possible by reducing the computational domain. In other words, by limiting ourselves to a relatively small target region of interest. This so-called ‘target-oriented’ approach (Valenciano et al. 2006; Haffinger et al. 2013; Willemsen et al. 2016; Yuan et al. 2017; Zhao & Sen 2018; Guo & Alkhalifah 2020) divides the medium into two parts: an overburden and a target. Target-oriented LSRTM (TOLSRTM) achieves this goal by bypassing the overburden and focusing the wavefield above the target of interest. This process, in which a wavefield is extrapolated from its current acquisition surface to another surface, is called redatuming (Berryhill 1984; Schuster & Zhou 2006; Luiken & van Leeuwen 2020; Barrera et al. 2021).
Target-oriented LSRTM, in addition to its lower computational costs, can also be helpful in cases the overburden produces strong internal multiple reflections that mask the reflections of the target region. Nonetheless, for this method to work correctly, a model of the overburden, which can produce all orders of reflections, is needed. Researchers have shown a strong interest in solving the problems caused by inaccurate overburden models and the artefacts that internal multiple reflections can create in LSRTM images. For example, Guo & Alkhalifah (2019,2020) introduced a least-squares waveform redatuming method to redatum the data to the boundary of the target and combined it with full waveform inversion.
Among different redatuming methods, Marchenko redatuming showed itself as a promising data-driven tool for predicting the Green’s functions at the boundary of the target region with all orders of internal multiples (Wapenaar et al. 2014; Ravasi et al. 2016; Dukalski & de Vos 2018; Vargas et al. 2021). The Marchenko redatuming can predict the internal multiples by only having the direct arrival between the surface and the target boundary in a smooth background velocity model of the overburden. Once we solve the Marchenko equations, we get the response of the medium stimulated by virtual sources at the boundary of the target and recorded by receivers at the surface, including all orders of multiple reflections. Hence, the first step of Marchenko redatuming can be interpreted as source-side redatuming in the physical medium (including reflections in the overburden that are not encoded in the background model). Diekmann et al. (2021) used Marchenko equations to predict the in-volume Green’s function of a target of interest and used this wavefield as the incident wavefield inside the scattering integral. Cui et al. (2020) used the reciprocity theorem to derive a forward modelling equation based on Marchenko wavefields for Full waveform inversion. Moreover, one of the byproducts of the Marchenko integrals is the downgoing part of the focusing function, which is the inverse of the transmission response of the overburden. This focusing function can also be applied at the receiver side to achieve double-sided (meaning source and receiver) redatuming (Staring et al. 2018). The process of double-sided redatuming with the help of the Marchenko focusing function is called Marchenko double-focusing.
In this paper we consider the upgoing Marchenko double-focused data as the ‘observed data’ that we feed to an LSRTM algorithm. To complete our TOLSRTM algorithm, we formulate a forward modelling to match the Marchenko double-focused data with the predicted data of the target region. This forward modelling is based on the downgoing wavefields from the Marchenko equation and a Born model of the target. This Marchenko downgoing wavefield acts as virtual source for the modelling and the inversion.
The organization of this paper is as follows. First, we briefly review LSRTM and target-oriented LSRTM. Secondly, Marchenko redatuming and double-focusing are explained. Then, our method of combining target-oriented LSRTM and Marchenko double-focusing is introduced. In the results section, we show numerical examples of our method and we exhibit the consequences of using an inaccurate velocity model for redatuming and TOLSRTM. Finally, results and potential future directions are discussed.
2 THEORY
2.1 LSRTM
2.2 Target-oriented least-squares reverse time migration
Researchers have proposed to implement the imaging and inversion in a target-oriented manner (Valenciano et al. 2006; Haffinger et al. 2013; Willemsen et al. 2016; Yuan et al. 2017; Zhao & Sen 2018; Guo & Alkhalifah 2020). In this way, the wavefield can be focused below a complex overburden to produce virtual sources and receivers right above the target (Berryhill 1984; Schuster & Zhou 2006; Luiken & van Leeuwen 2020; Barrera et al. 2021). These redatumed wavefields can be used for imaging and inversion. The advantages of this approach are as follows: (1) The wavefield is focused below specific structures, such as salt domes and (2) the computational costs are lower because the computational domain can be reduced significantly. However, current methods need a reasonable estimation of the velocity model of the overburden (i.e. the part of the model above the target of interest), which can predict internal multiple reflections. Otherwise, these overburden-generated multiple reflections produce artefacts inside the retrieved image of the target.
Marchenko redatuming enables us to take these overburden-generated multiples into account and improve the quality of the target-oriented LSRTM images. The following section gives a brief review of Marchenko redatuming by double-focusing.
2.3 Marchenko redatuming by double-focusing
Marchenko redatuming is a novel data-driven method that can retrieve the Green’s function at a surface above the target area with all orders of its multiple-scattered events. The only requirements for this method are the reflection response at the surface and a smooth background velocity model of the overburden that can predict the direct arrival from the surface to the redatuming level.

The Green’s functions resulting from Marchenko redatuming and double-focusing. (a) Downgoing part of the Marchenko Green’s function, (b) upgoing part of the Marchenko Green’s function, (c) downgoing Marchenko double-focused Green’s function and (d) upgoing Marchenko double-focused Green’s function.
Through Marchenko double-focusing, we retrieve a downgoing Green’s function (|$G^{+,+}_{\rm df}$|) which contains a band-limited delta function and interactions between target and overburden, and an upgoing Green’s function (|$G^{-,+}_{\rm df}$|), which is the upgoing response to |$G^{+,+}_{\rm df}$| from the target at the redatuming level, still including interactions between target and overburden on the source side (Figs 1c and d). In comparison, ‘conventional double-focusing’ is the act of using the inverse of the direct arrival of the transmission response of the overburden instead of the downgoing Marchenko focusing function for double-focusing. Considering that the inverse of the direct arrival of the transmission response of the overburden does not contain multiple reflections, it is unable to predict and remove the multiples generated by the overburden. In the next section, ‘double-focusing’ is a general term referring to both methods, and it is mentioned wherever a distinction is necessary.
2.4 Target-oriented LSRTM by Marchenko double-focusing
It is possible to compute eqs (27) and (28) in two ways. The first way is to compute the background Green’s function (G0) with an algorithm of choice, for example finite-difference, and then insert it in eqs (27) and (28). The second way is to inject |$G^{+,+}_{\rm df}$| as a downward radiating dipole line source for each |${\bf x}^{\prime }_{\rm vs}$| location in the target background medium with a numerical PDE solver to compute |$P^{\rm inc}_{\rm df}$|. Then inject |$\frac{\omega ^2}{\rho _0}\delta m({\bf x}) P^{\rm inc}_{\rm df}$| as a contrast source to compute |$\hat{P}^{\rm scat}_{\rm pred}$|. The decision of choosing the appropriate computation method is based on the dimensions of the target and the number of virtual sources.
3 NUMERICAL RESULTS
3.1 Marchenko double-focusing versus conventional double-focusing
To validate our theory, we use two different synthetic models. We compute the reflection data sets at the acquisition surface using a finite-difference algorithm (Thorbecke & Draganov 2011) and a dipole (vertical force) source with a flat spectrum wavelet (0–100 Hz). The injection rate is represented as an iω factor in eq. (2). Further, we use a spatial sampling of 2.5 m in both directions and a time sampling of 0.4 ms. To reduce the data size, we set the receiver temporal sampling to 4 ms. For the sake of better visualization, we convolve all time–space records with a 30 Hz Ricker wavelet.
Note that we do not commit an inverse crime in these examples. As we mentioned earlier, to compute |$\hat{P}^{\rm scat}_{\rm obs}$|, we compute the reflection response at the acquisition surface by a finite-difference algorithm. Then we apply Marchenko double-focusing to this data set. On the other hand, to compute |$\hat{P}^{\rm scat}_{\rm pred}$| we use eqs (30), (28) and 27. In all of the mentioned steps, we deal with convolutions that lead to a loss of resolution. To avoid it, we compute everything with a flat spectrum wavelet (0–100 Hz), and then we convolve both |$\hat{P}^{\rm scat}_{\rm obs}$| and |$\hat{P}^{\rm scat}_{\rm pred}$| with a 30 Hz Ricker wavelet.
For comparison, we compute the redatuming operator in two different ways for both models. First, we isolate the direct arrival of the transmission response of the overburden with a time window and use it as an estimation of the redatuming operator and redatum the surface data with it. This approach is equivalent to conventional redatuming of primary reflections (Berryhill 1984). Secondly, we solve the coupled Marchenko equations with help of this direct arrival to compute the redatumed Green’s and focusing functions (Thorbecke et al. 2017). We calculate the direct arrival in this way to avoid kinematic errors. We discuss the issue of kinematic errors later in Section 3.3.
3.1.1 Simple model
As a first example, we define a model with a single-layered overburden and a single small diffractor with dimensions of 5 m by 5 m inside a homogeneous background as the target (Fig. 2). The background velocity is set to 1500 m s–1, and the velocity of the overburden layer is 3000 m s–1. The density is constant and equal to 2500 kg m–3 everywhere. For this model, 101 sources and receivers with a spacing of 10 m are used. The number and spacing of real and virtual sources and receivers are the same. In Fig. 2 blue crosses are the actual source and receiver locations, and the red dots are the virtual source and receivers locations. Moreover, the overburden is designed such that a higher-order multiple reflection coincides with the primary event of the diffractor. We design it in this way to show that our method can handle the overburden-related multiple reflections correctly.

Velocity model with a single-layered overburden and a small diffractor in the target.
In Fig. 3, the common-source record at the acquisition surface is shown for this model. The diffractor response of the simple model is visible around 0.5 s in this figure and it partly coincides with a higher-order multiple reflection.

Common-source record at the surface of a source located at lateral distance = 0 in the simple model. A linear time-varying (|$\frac{t}{{\rm d}t}$|) gain is applied to this record to amplify weaker events.
In Fig. 4, we show the double-focused data of this model. From here onward, we use the subscript ‘Cdf’ for the data produced by ‘conventional double-focusing’ and ‘Mdf’ for ‘Marchenko double focused’ data. In Fig. 4(b), it is visible that the overburden-related multiple reflections have been suppressed by applying Marchenko double focusing, yet the interactions between the overburden and the target are still present. Moreover, in Fig. 5, we show the double-focused downgoing Green’s function (eq. 12). The downgoing wavefield of the Marchenko double-focusing (Fig. 5b) contains the higher-order interactions between the target and the overburden, whereas the downgoing wavefield of the conventional double-focusing does not include these interactions. Injecting the downgoing Marchenko double-focused wavefield into the target from the focusing depth reproduces these high-order interactions in the predicted data of the TOLSRTM algorithm.

Double-focused (upgoing) data of model with a single-layered overburden corresponding to a virtual source at the middle. (a) Conventional double-focused data (|$G^{-,+}_{\rm Cdf}$|) by applying conventional redatuming operators and (b) Marchenko double-focused data (|$G^{-,+}_{\rm Mdf}$|). The maximum and the minimum value of the grey-level scale of both figures are the same and a linear time-varying (|$\frac{t}{{\rm d}t}$|) gain is applied to both records to amplify weaker events.

Virtual source wavefield of the model with a single-layered overburden at the middle virtual source location. (a) The conventional double-focused source (|$G^{+,+}_{\rm Cdf}$|) by applying conventional redatuming operators and (b) Marchenko double-focused source (|$G^{+,+}_{\rm Mdf}$|). The maximum and the minimum value of the grey-level scale of both figures are the same and a linear time-varying (|$\frac{t}{{\rm d}t}$|) gain is applied to both records to amplify weaker events.
In Fig. 6 the predicted data of both approaches after 20 iterations are shown. The predicted data of the conventional double focusing is not a good fit for the conventional double-focused data. In contrast, the predicted data of the Marchenko double-focusing approach is a better fit for Marchenko double-focused data. In Fig. 6(b), the weak multiple reflections generated by interactions between target and overburden can be seen inside the red box. This confirms that our method is not only able to suppress multiple reflections generated inside the overburden significantly but also able to predict the interactions between the target and the overburden.

Predicted data with a virtual source located at the middle of the model with a single-layered overburden after 20 iterations. (a) Conventional double-focusing and (b) Marchenko double-focusing. Weak multiple reflections that are predicted by our method can be seen inside the red box. The maximum and the minimum value of the grey-level scale of both figures are the same and a linear time-varying (|$\frac{t}{{\rm d}t}$|) gain is applied to both records to amplify weaker events.
Figs 7 and 8 show the imaging result of these two approaches for the first iteration and after 20 iterations for the model with single-layered overburden. Note the horizontal event in Fig. 6, coming from the overburden multiple, which overlaps the image of the target. In Fig. 7 this horizontal event is strongly suppressed. To demonstrate the resolution improvement of our method, we compare the horizontal and vertical sections of the image of different approaches in Figs 9(a) and (b). As these figures are showing, our approach has a superior resolution. A comparison of the convergence rate of the cost functions of these approaches is also shown in Fig. 10.

Retrieved target image of the model with a single-layered overburden by conventional double-focusing operators. (a) RTM and (b) LSRTM after 20 iterations.

Retrieved target image of the model with a single-layered overburden by Marchenko double-focusing. (a) RTM and (b) LSRTM after 20 iterations.

A Comparison of the Horizontal and vertical cross-section of the retrieved image of the simple model with different approaches. (a) A horizontal cross-section at depth of 450 m and (b) a vertical cross section at lateral distance zero.

A comparison between convergence rate of cost functions of different approaches to TOLSRM for the simple model.
3.1.2 Syncline model
As a second example, we design a model with more layers and a syncline-shaped layer in its overburden (Fig. 11a) and a line reflector with a width of 5 m in the target region. In this model, the background velocity of the target area is 2400 m s–1, and the velocity of the reflector is 2700 m s–1. Moreover, density variations are introduced in the overburden (Fig. 11b) to generate substantial internal multiples, but the density of the target region is constant and equal to 1000 kg m–3. For this model, 301 sources and receivers with the same spacing as in the first model are deployed. The number and spacing of real and virtual sources and receivers are the same. In Figs 11(a) and (b) blue crosses are the actual source and receiver locations, and the red dots are the virtual source and receivers locations.

Velocity and density of the model with the syncline overburden and a line reflector in the target.
Fig. 12 shows the middle common-source record at the acquisition surface for the model with the syncline overburden. In Fig. 13, we show the conventional double-focused data (Fig. 13), the predicted data in case of conventional double-focusing (Fig. 13a), the Marchenko double-focused data (Fig. 13c) and the predicted data in case of Marchenko double-focusing (Fig. 13b) for the model with the syncline overburden. In Fig. 13(d), it can be seen, compared to conventional double-focused data (Fig. 13a), that the multiples purely generated in overburden are suppressed by applying Marchenko double focusing. In addition, by comparing Figs 13(c) and (d) we see that the predicted data of Marchenko double-focusing not only predicts the primary reflection, but also predicts the multiples that are generated by interactions between the target and the overburden.

Common-source record at the surface of a source located at lateral distance = 0 in the model with the syncline overburden. A linear time-varying (|$\frac{t}{{\rm d}t}$|) gain is applied to this record to amplify weaker events.

A comparison between Double-focused data and predicted data of the model with the syncline overburden corresponding to the middle virtual source. (a) Conventional double-focused data, (b) Marchenko double-focused data, (c) predicted data of conventional double-focusing and (d) predicted data of Marchenko double-focusing. The first few time samples are set to zero to mute focusing artefacts. The maximum and the minimum value of the grey-level scale of all figures are the same.
Moreover, in Fig. 14, we show the double-focused downgoing Green’s function (eq. 12). We see that the virtual source wavefield for the conventional double-focused predicted data contains only a bandlimited delta function, whereas the virtual source wavefield for the Marchenko double-focused predicted data contains additional events to predict the target and overburden interactions in the Marchenko double-focused data in Fig. 13(b).

Virtual source wavefield of the model with the syncline overburden for the middle virtual source. (a) The conventional double-focused source (|$G^{+,+}_{\rm Cdf}$|) by applying conventional redatuming operators and (b) Marchenko double-focused source (|$G^{+,+}_{\rm Mdf}$|). The maximum and the minimum value of the grey-level scale of both figures are the same.
Consequently, the conventional double-focusing approach is unable to fit the predicted data, whereas the Marchenko double-focusing approach can produce predicted data with an acceptable fit to the Marchenko double-focused data. Note that the predicted data of our proposed method can predict the interactions between the target and the overburden, without the knowledge of the overburden reflectivity model, by taking |$G^{+,+}_{\rm Mdf}$| as the virtual source wavefield. This source contains these interactions, as can be seen in Fig. 14(b).
Figs 15 and 16 show the imaging result of these two approaches for the first iteration and after 20 iterations for the model with the syncline overburden. To demonstrate the improvement of vertical resolution and lateral continuity, we compare the horizontal and vertical sections of the image of different approaches in Figs 17(a) and (b). As these figures are showing, our approach has a superior vertical resolution and lateral continuity. A comparison of the cost functions of these different approaches to TOLSRTM is also shown in Fig. 18. TOLSRTM converges faster when Mdf (instead of Cdf) is used for redatuming, since the multiple reflections are correctly taken into account.

Retrieved target image of the model with the syncline overburden by conventional double-focusing operators. (a) RTM and (b) LSRTM after 20 iterations.

Retrieved target image of the model with the syncline overburden by Marchenko double-focusing. (a) RTM and (b) LSRTM after 20 iterations.

A Comparison of the Horizontal and vertical cross-section of the retrieved image of the model with the syncline overburden with different approaches. (a) A horizontal cross-section at depth of 1200 m and (b) a vertical cross-section at lateral distance zero.

A comparison between convergence rate of cost functions of different approaches to TOLSRM for the model with the syncline overburden.
3.2 Investigating the effect of the point-spread function
In this section, we investigate the effect of the PSF and the necessity of applying it to the reflection response of the target area. In Fig. 19, we show the effect of the PSF (eq. 16) on the true reflection response (as modelled during each iteration of TOLSRTM) of the target area of the simple model and compare it with Marchenko double-focused data. Fig. 19(a) shows the modelled data (eq. 28) in the true perturbation model without PSF. Fig. 19(b) shows the same after the PSF has been applied. Fig. 19(c) shows the Marchenko double-focused data. We can see the effect of the PSF by analysing Figs 19(b) and (c). The PSF not only imposes a band limitation on the data but also dampens the far-offsets such that the predicted data fits the observed data better. Thus, applying the PSF to the reflection response is necessary for better convergence. Moreover, in Fig. 20 the cost function of TOLSRTM with and without PSF is shown. From this figure, we conclude that including the PSF in our formulation successfully improves the convergence.

Comparison of the reflection response of target area of the model with a single-layered overburden. (a) Without PSF, and (b) With PSF and (c) double-focused data. The maximum and the minimum value of the grey-level scale of figures are set to the maximum and the minimum of the input.

Comparison of the cost function of the TOLSRTM with and without PSF.
3.3 Inaccurate velocity model
To determine the sensitivity of our algorithm to kinematic errors, we define a model similar to the model with a syncline in the overburden but with a high-velocity inclusion in the overburden and a more complex target (Fig. 21a). We use a smooth version of this model (Fig. 21b) to estimate the direct arrival between the surface and the focusing level. In addition, we use the target part of this smooth velocity model to compute the background Green’s functions with a finite difference algorithm for the forward and adjoint kernels. The density model is shown in Fig. 21(c). In Fig. 22, we show images as obtained by applying our methodology with the inaccurate velocity model. The yellow arrows and ellipses in Fig. 22 indicate discontinuities and depth estimation errors caused by velocity inaccuracies. Moreover, the red arrows in this figure designate a couple of artefacts that are removed by our method, and the white arrows and ellipses show the amplitude improvements resulting from our method. Despite the fact that we use a wrong velocity model, our algorithm which is based on Marchenko double-focusing is still able to correctly predict and reduce the artefacts caused by the overburden-related multiple reflections and improves the quality and amplitude recovery of the LSRTM image. The observation that Marchenko imaging is not sensitive to errors in the velocity model is not new, see for instance Broggini et al. (2014).

Velocity model: (a) true model, (b) smooth model and (c) density model.

Target-oriented LSRTM with the Erroneous velocity model. (a) Conventional double-focusing and (b) Marchenko double-focusing. The red arrows show some overburden artefacts that are removed, the white arrows and ellipses show amplitude improvement and the yellow ones show some of the effects of velocity errors. The value of the grey-level scale is set to 25 per cent of the maximum value.
4 DISCUSSION
Our numerical results suggest that our Marchenko target-oriented LSRTM can recover a clearer image of the target of interest than conventional TOLSRTM. The retrieved image of the model with single-layered overburden using Marchenko double-focusing (Fig. 8b) has no (or at least a reduced) imprint of the overburden-related spurious reflector. In contrast, the recovered image using conventional redatuming operators (Fig. 7b) could not handle this spurious reflector even after 20 iterations. A further investigation into the resolution improvement of our LSRTM algorithm compared to the conventional double-focused LSRTM, shows that our algorithm is able to improve the resolution in both horizontal and vertical directions (Fig. 9). The model with the syncline overburden demonstrates the power of using Marchenko double-focusing for target-oriented LSRTM in removing the strong overburden-related multiple reflection artefacts. The image that is obtained with the conventional redatuming operator (Fig. 15b) is overwhelmed by artefacts from the overburden. On the other hand, our algorithm is able to remove these artefacts and retrieves a clear image of the target area (Fig. 16b). In addition, our TOLSRTM algorithm can model data that can predict the interactions between the target area and the overburden without the knowledge of the overburden reflectivity model. Moreover, the Marchenko equation has been derived for elastic media too (da Costa Filho et al. 2014; Wapenaar & Slob 2014; Reinicke et al. 2020). Hence, we expect that the current methodology can also be extended to elastic LSRTM.
Furthermore, a comparison between Figs 4(b) and 19 shows the necessity of applying the PSF to the modelled reflection data of the target. Without imposing this PSF on the modelled data, the predicted data (Fig. 19a) do not match the observed data (Fig. 19c). Hence, the TOLSRTM without PSF converges slowly compared to the case where the PSF is applied to it (Fig. 20).
In our last example, we addressed the issue of velocity inaccuracies in the overburden. As Fig. 22 shows our method is able to improve the quality of the image by suppressing the overburden-generated multiples and predicting the interactions between target and overburden, even if the overburden velocity model is inaccurate.
The Marchenko redatuming algorithm used in this paper requires surface-related multiples to be removed from the acquired reflection data. However, researchers exploited the possibility of including the surface-related multiple reflections within the Marchenko framework (Ravasi 2017; Singh et al. 2017; Dukalski & de Vos 2018). Consequently, it is quite straightforward to include surface-related multiples in our algorithm.
5 CONCLUSION
We have introduced a target-oriented LSRTM algorithm that can correctly handle the internal multiple reflections generated in the overburden and interactions between target and overburden. To achieve this, we solved the Marchenko equations to obtain the Marchenko double-focused data and the downgoing part of the focusing function. By having these Green’s functions and focusing operators, we can formulate a Born integral that can be used as the forward modelling operator and an adjoint modelling operator can also be constructed for target-oriented LSRTM. Importantly, we have to apply a point-spread function to the reflection response of the target to account for the finite spatial bandwidth caused by the overburden and the finite recording aperture at the acquisition surface.
Since LSRTM and full waveform inversion are closely related, the question may arise: ‘Can we modify our algorithm for target-oriented full waveform inversion?’. To address this question, we have to mention that the double-focused data does not include the diving waves, so it is dominated by the medium-to-short wavelengths and does not include long wavelengths. These long wavelengths are crucial for velocity recovery in the full waveform inversion process. Moreover, the evanescent waves, such as refracted waves, are not handled by the current redatuming methodologies. The aforementioned issues put a limit on modifying our algorithm to be applied for target-oriented full waveform inversion, which requires long-to-medium wavelengths. However, recently some advancement has been made towards including evanescent waves within the Marchenko framework (Wapenaar 2020; Diekmann & Vasconcelos 2021; Wapenaar et al. 2021). Nevertheless, the double-focused data contains all target reflections correctly (but not the other waveforms). Therefore, another possible future direction of this research is to investigate the potential of our method in improving inversion methods that rely on reflections only, such as joint migration inversion (Verschuur et al. 2016) and reflection waveform inversion (Brossier et al. 2015).
In recent years, the focus of the seismic exploration community has shifted toward fast and high-resolution imaging and inversion methods. Our proposed target-oriented LSRTM is able to produce such images in a relatively small target region and has a relatively faster convergence rate while correctly accounting for multiple reflections caused by the overburden.
ACKNOWLEDGEMENTS
We thank the reviewers, Katrin Löer and Ole Edvard Aaker, for their constructive comments, which helped us to improve the paper. This work has received funding from the European Union’s research and innovation programme: European Research Council (grant agreement 742703).
DATA AVAILABILITY
This study uses numerical data only. All of the models and data sets are produced by the authors, and the codes used will be shared on reasonable request to the corresponding author.