Subsurface-domain , interferometric objective functions for target-oriented waveform inversion

Resolving details in subsurface reservoir parameters from surface waveform data is a challenging problem, particularly when a reservoir is beneath a complex overburden whose properties are also poorly known. We have developed new metrics for detecting and quantifying errors in subsurface models that are well-suited for target-oriented inversion. Our metric, when combined with state-of-the art redatuming, aims at enabling waveform inversion for target volume parameters with no need to resolve model features elsewhere (e.g., in the overburden). We refer to these metrics as “interferometric objective functions” because they rely on extrapolation from reciprocity integrals commonly used in seismic interferometry. As in seismic interferometry, wavefield extrapolation retrieves the wave response between two points by combining observed data with an extrapolator that describes the response between the subsurface and the data boundary. When the source point is outside a target volume, either forward time or reverse time extrapolation produces the same field. However, because they physically rely on different components from the boundary data, the forward time and reverse time extrapolated fields are only equal when the model used is consistent with the real subsurface within the target volume. As such, we use the difference between the forward time and reverse time extrapolated fields to define subsurface-domain metrics that quantify model errors. Our approach thus provides a new metric for target-oriented nonlinear inversion in the subsurface domain, one that fundamentally differs from other subsurface domain metrics based on, e.g., focusing or image extensions. INTRODUCTION Estimating deep-reservoir properties in highly complex settings is a challenging task for exploration geophysicists. Data-misfitbased full-waveform inversion methods (Plessix, 2006; Virieux and Operto, 2009; Warner et al., 2013) can handle nonlinearities in the waveforms, but they struggle to achieve high-resolution models from reflection data. On the other hand, image-domain tomography (Symes, 2008; Yang and Sava, 2011; Fleury and Perrone, 2012; Biondi and Almomin, 2013) relies on finite-frequency information from deep reflectors but, in turn, cannot easily handle nonlinear effects and it requires the generation of depth-image gathers. Building on the concepts of interferometry-based imaging (Halliday and Curtis, 2010; Vasconcelos et al., 2010; Vasconcelos, 2013) and exact boundary extrapolation (van Manen et al., 2007; Vasmel et al., 2013), here, we propose a new subsurface-domain waveform inversion metric that can fully handle nonlinear effects, does not require an imaging step or image gathers, and is aimed at inversion in a target-oriented fashion. Our approach requires data on a boundary of receivers enclosing a volume of interest, with sources anywhere outside the bounded volume. Even though this at first may seem as an impediment to the method because it would not usually be available as directly observed data, we point out that the needed subsurface fields can in principle be reconstructed by redatuming procedures. One approach that may enable this is the data-driven Marchenko redatumManuscript received by the Editor 30 November 2016; revised manuscript received 10 February 2017; published ahead of production 02 May 2017; published online 19 June 2017. Formerly Schlumberger Gould Research, Cambrige, UK; presently Utrecht University, Utrecht, the Netherlands. E-mail: i.vasconcelos@uu.nl. Formerly University of Edinburgh, Edinburgh, UK; presently Statoil, Bergen, Norway. E-mail: mrava@statoil.com. Delft Technical University, Delft, the Netherlands. E-mail: j.r.vanderneut@tudelft.nl. © The Authors.Published by the Society of Exploration Geophysicists. All article content, except where otherwise noted (including republished material), is licensed under a Creative Commons Attribution 4.0 Unported License (CC BY-NC-ND). See https://creativecommons.org/licenses/by-nc-nd/4.0/ Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its digital object identifier (DOI). Commercial reuse and derivatives are not permitted. A37 GEOPHYSICS, VOL. 82, NO. 4 (JULY-AUGUST 2017); P. A37–A41, 4 FIGS. 10.1190/GEO2016-0608.1 Downloaded from https://pubs.geoscienceworld.org/geophysics/article-pdf/82/4/A37/4075963/geo-2016-0608.1.pdf by guest on 15 January 2019 ing method as proposed by Broggini et al. (2012) and Wapenaar et al. (2013, 2014a). These approaches are shown to robustly retrieve full-waveform redatumed fields, with balanced amplitudes and reconstructed internal multiples, using little detailed knowledge of the subsurface (van der Neut et al., 2015), even in the presence of overburden model errors (Thorbecke et al., 2013; Wapenaar et al., 2014b). Here, we point out that the combination of these new redatuming schemes with our proposed metric is the main motivation of this work: It allows for waveform inversion to act only over a target volume, with no need to estimate model details from overor underburden structures. AN INTERFEROMETRIC IDENTITY FOR A TARGET VOLUME All acoustic waves propagating between the points xB and xA, shown in the geometric configuration of Figure 1, can be retrieved by crosscorrelations by means of the integral (Wapenaar et al., 2011)


INTRODUCTION
Estimating deep-reservoir properties in highly complex settings is a challenging task for exploration geophysicists.Data-misfitbased full-waveform inversion methods (Plessix, 2006;Virieux and Operto, 2009;Warner et al., 2013) can handle nonlinearities in the waveforms, but they struggle to achieve high-resolution models from reflection data.On the other hand, image-domain tomography (Symes, 2008;Yang and Sava, 2011;Fleury and Perrone, 2012;Biondi and Almomin, 2013) relies on finite-frequency information from deep reflectors but, in turn, cannot easily handle nonlinear effects and it requires the generation of depth-image gathers.Building on the concepts of interferometry-based imaging (Halliday and Curtis, 2010;Vasconcelos et al., 2010;Vasconcelos, 2013) and exact boundary extrapolation (van Manen et al., 2007;Vasmel et al., 2013), here, we propose a new subsurface-domain waveform inversion metric that can fully handle nonlinear effects, does not require an imaging step or image gathers, and is aimed at inversion in a target-oriented fashion.
Our approach requires data on a boundary of receivers enclosing a volume of interest, with sources anywhere outside the bounded volume.Even though this at first may seem as an impediment to the method because it would not usually be available as directly observed data, we point out that the needed subsurface fields can in principle be reconstructed by redatuming procedures.One approach that may enable this is the data-driven Marchenko redatum-ing method as proposed by Broggini et al. (2012) and Wapenaar et al. (2013Wapenaar et al. ( , 2014a)).These approaches are shown to robustly retrieve full-waveform redatumed fields, with balanced amplitudes and reconstructed internal multiples, using little detailed knowledge of the subsurface (van der Neut et al., 2015), even in the presence of overburden model errors (Thorbecke et al., 2013;Wapenaar et al., 2014b).Here, we point out that the combination of these new redatuming schemes with our proposed metric is the main motivation of this work: It allows for waveform inversion to act only over a target volume, with no need to estimate model details from overor underburden structures.

AN INTERFEROMETRIC IDENTITY FOR A TARGET VOLUME
All acoustic waves propagating between the points x B and x A , shown in the geometric configuration of Figure 1, can be retrieved by crosscorrelations by means of the integral (Wapenaar et al., 2011) where pB and ĜA are shorthand notations for pðx; x B Þ and Ĝðx; x A Þ, respectively.These fields in the integrand are the superposition of all the ingoing and outgoing waves shown in Figure 1.Throughout this paper, with wavefield extrapolation/redatuming in mind, we assume that all p fields correspond to "boundary data," which are either observed directly or are indirectly reconstructed (Wapenaar et al., 2013).All Ĝ fields are full-wave "wavefield extrapolators," and they denote the impulse response of a wavefield modeling operator of choice (e.g., finite differences), given a subsurface model.The notation ∂ n indicates local normal derivatives.The Ã superscript denotes complex conjugation in the frequency domain, corresponding to a time reversal in the time domain.The geometric configuration in Figure 1 also allows for retrieving pðx A ; x B Þ by means of crossconvolution (Wapenaar et al., 2011): The main and most general integral identity we discuss in this study is based on recognizing that the correlation-and convolution-based integrals in equations 1 and 2 retrieve the same response (with opposite polarity).Thus, by adding equations 1 and 2, we obtain the identity Physically, this identity states that forward time (convolution-based) and reverse time (correlation-based) extrapolations yield the same fields.A key underlying point in all these integral results is that the medium parameters must be the same for pB and ĜA (Wapenaar et al., 2011); i.e., because pB are waves in the real medium, for the integrals to hold ĜA must be the "true medium" extrapolator within the target integration domain.
For our purposes, it is important to note that although forward time and reverse time extrapolations yield the same final result, they each rely on different physical contributions from the boundary data.To illustrate this point, by giving the data with ingoing and outgoing wave components (Figure 1a); e.g., if we choose an extrapolation model that is homogeneous outside of the bounded domain (setting the ingoing waves ĜðinÞ A ¼ 0), equation 3 can be rewritten as where the contributions of the terms pðinÞ B ð∂ n ĜðoutÞ A Þ Ã and pðoutÞ

B ð∂ n ĜðoutÞ
A Þ vanish at their stationary regions (Wapenaar and Fokkema, 2006), respectively, for reverse and forward time extrapolation.From the contributions in the integrand, we see that the forward time (convolution) extrapolation uses only the ingoing waves pðinÞ B from the data, whereas reverse time (correlation) extrapolation uses only the outgoing waves pðoutÞ B .

SUBSURFACE-DOMAIN OBJECTIVE FUNCTION
As mentioned above, it follows from the derivations of integral reciprocity that the interferometric identity in equation 3 only holds for the correct extrapolator ĜA (i.e., corresponding to the true medium properties) that satisfies both integral relations in equations 1 and 2. Therefore, for a given set of observed pB fields, using any trial model m that yields Ĝ 0 A ≠ ĜA would result in As such, we can immediately define the frequency-domain metric as or its time-domain counterpart where F −1 is the inverse Fourier transform.Here, we use Î ¼ Îð pB ; Ĝ 0 A ðmÞÞ to denote the explicit dependency of the interferometric identity on the boundary data pB and of the extrapolator on the model parameters Ĝ 0 A ðmÞÞ.According to the exact interferometric relation in equation 3, these metrics will be minimum (i.e., equal to zero) when

NUMERICAL EXAMPLES
The following numerical experiment is designed to study the effect of choosing different media outside of the bounded domain for extrapolation, as well as the dependence of the medium within the target volume.With these objectives, we use the boundary acquisition and model configurations shown in Figure 2. In this case, the boundary data used in all of the extrapolation examples are acquired using the model in Figure 2a; i.e., the data contain all scattering interactions with the inhomogeneities inside and outside of the bounded domain.We use a single, fixed-source position at x B in all of our examples, with receivers at discrete positions x, uniformly sampled at 7.6 m along the boundary (Figure 2), recording pressure, and both components of particle velocity.The source excitation function is a zero-phase Ricker wavelet with 15 Hz peak frequency.The wavefield extrapolation follows the scheme by Vasconcelos (2013).For the purpose of this paper, we use the directly modeled boundary fields as input for testing our metric because these direct observations would be the best-case scenario.In the more realistic case in which only surface reflection data are available, these fields would be retrieved by means of redatuming (van der Neut et al., 2015): The use of our metric with such redatuming approaches, along with an analysis of the assumptions and errors associated to redatuming, is the subject of the ongoing research.
Given that the input boundary data for all of the panels are the same (i.e., acquired with the model in Figure 2a), in Figure 3, we show the forward time and reverse time (i.e., equations 2 and 1) extrapolated fields and the corresponding differences used in the interferometric identity (equation 3).When comparing the extrapolated fields in Figure 3a and 3b, it is impossible to visually distinguish which field was extrapolated forward in time and which is time reversed: Indeed the wavefield differences are negligible (Figure 3c).Note that for the top panels of Figure 3, the extrapolation is done with a homogeneous overburden (Figure 2b), and the difference between these fields and those using the true model (not shown) are also negligible.This changes for the bottom panelshere, the forward time (Figure 3d) and reverse time (Figure 3e) fields differ substantially as shown by Figure 3f.This discrepancy is caused by the model errors only within the target volume, and it contains the information that is used by our interferometric identities.
The equivalence of choosing either of the models of Figure 2 for extrapolation is further supported by evaluating the interferometric misfit J t ðx A ; x B Þ in equation 7, stacked over time, for x A varying over a target volume (Figure 4).This observation verifies a key aspect of the interferometric misfit in equation 7: If the model properties inside the target domain are correct, our J metric is also minimum when using an extrapolation model that is completely homogeneous outside the bounded domain.Because of this key feature, the interferometric misfit can be used to detect model errors in the inside of the boundary for data collected using the model in

CONCLUSION
In this paper, we present a metric to quantify errors in subsurface models using wavefield extrapolation schemes based on convolution-and correlation-type reciprocity.Because they rely on different physical contributions from the boundary wavefields, the fields extrapolated in forward and reverse time are the same only when the subsurface model used for extrapolation is consistent with that of the real parameters that generate the observed data.We define a subsurface-domain metric for inversion based on the difference be-tween forward time and reverse time extrapolated fields.This metric relies on full-waveform data and extrapolators; therefore, it allows for nonlinear inversion of the waveform misfit.
Our proposed subsurface-domain misfit has features that distinguish it from existing tomography/inversion metrics.Although this misfit is defined in the subsurface domain, the interferometric misfit is not a depth-image-based metric because it does not require evaluating an imaging condition or depth-domain image gathers.We show that our misfit can in fact be evaluated in the subsurface domain even for a single source point, although it can also take advantage of multi- ple source locations.In addition, being well-defined in the frequency and time domain, our metric, in principle, allows for multiscale inversion approaches in the subsurface domain.In this paper, we validate our metric as a reliable measure for model errors and demonstrate that it is sensitive only to the chosen volume; the calculation of model gradients and detailed quantitative benchmarking against conventional FWI of surface data are topics of ongoing research.In terms of computations, although this has not yet been tested, we expect the cost of inversion with this metric to be primarily controlled by the cost of wave-equation modeling only in the target volume (assuming the redatuming cost to be negligible compared with modeling); for targetvolume inversion, this method should remain as costly as, but more likely less costly than, conventional FWI.
The main objective of our metric is for it to be used in a targetoriented fashion; i.e., given sufficient enclosing boundary data acquired in an arbitrarily heterogeneous medium, the interferometric misfit can be used to quantify subsurface errors only in the inside of a target-bounded domain, with little knowledge of the medium properties outside it.This follows from the fundamentals of integral reciprocity, which allow for different choices of extrapolation models outside the boundary of interest.When combined with redatuming approaches that can retrieve the necessary target-bounding wavefields while using conventional migration velocity models (e.g., redatuming by the Marchenko or inverse scattering series), our metric provides a means for waveform inversion for model details only in the target volume.This sidesteps the need for resolving model details everywhere in the medium (not only in the target) to fit surface data, as would be the case in conventional FWI.This feature, which holds for arbitrarily complex models, may prove advantageous for applications such as targeted reservoir characterization and localized timelapse monitoring using local waveform information.In practice, the reliance of this approach on redatuming (such as the Marchenko method) will impose limitations for inversion because of the assumptions and errors associated with the redatuming process.We are currently investigating these issues and will report on them in the future.Finally, because of the general nature of the reciprocity integrals used in our extrapolation-based metrics, the extension to arbitrary elastic media is straightforward, including anisotropy and attenuation effects.

Figure 1 .
Figure 1.(a) Illustration of the data and (b) extrapolator fields.The sources are at x A and x B , and the receivers are on the boundary.The ingoing and outgoing components of the fields are also shown.
Figure 2a, without knowledge of parameters outside of the boundary.This case is shown in Figure 4c, in which we use a smoothed version of the model in Figure 2b for extrapolation.In Figure 4c, warm colors indicate misfit increases due to model changes only in the interior of the bounded domain, with no actual need for the exterior model parameters present in the actual configuration in which the data are acquired (Figure2a).Although the panels in this figure feature the unknown parts of the model, we note that they are the objective function values at each position x in the target volume and are not to be confused with, e.g., gradient-based back projections as commonly seen in conventional full-waveform inversion approaches.

Figure 2 .
Figure 2. The configuration used for the numerical experiment.A star denotes the source position, and triangles denote the receivers.The receivers used for recording the data are placed on the top-and-bottom enclosing boundary ∂D ¼ ∂D t ∪ ∂D b (green and pink lines), whereas all points x A in the interior of the bounded domain eventually act as virtual receiver locations after wavefield extrapolation.(a) A model that is heterogeneous outside of the bounded domain and contains the "true model" inside of the domain.Panel (b) is identical to that in (a) inside of the bounded domain, but it is homogeneous outside.

Figure 3 .Figure 4 .
Figure 3. Snapshots of target-volume extrapolated wavefields and wavefield differences.Panels (a-d) correspond to the model in Figure 2b, whereas (d-f) pertain to an incorrect, smoothed version of the model in Figure 2b.Panels (a and d) are the forward time extrapolated fields, (b and e) are the reverse time fields, and (c and f) are the difference between the corresponding forward time and reverse time snapshots.All plots use the same amplitude scale, set according to the true subsurface field (not shown).
m is the model that produces the correct Ĝ 0 A ðmÞÞ ¼ ĜA .Waveform inversions for the model m can then be pursued by either based on the time-domain equation 7.In practice, in the presence of noisy finite-bandwidth data, we expect the frequency-and timedomain metrics to behave differently in the inverse problem.One important practical difference is that the frequency-domain metric allows for natural multiscale approach by separating different frequencies or frequency bands.In addition, it is likely that each of these metrics may require different approaches to preconditioning and model regularization.These issues are the subject of ongoing research.