Freely available online through the SEG open-access option.


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.


Estimating deep-reservoir properties in highly complex settings is a challenging task for exploration geophysicists. Data-misfit-based 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 redatuming 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 over- or underburden structures.


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)  
where p^B and G^A are shorthand notations for p^(x,xB) and G^(x,xA), 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 G^ 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^(xA,xB) 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 p^B and G^A (Wapenaar et al., 2011); i.e., because p^B are waves in the real medium, for the integrals to hold G^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 G^A(in)=0), equation 3 can be rewritten as  
where the contributions of the terms p^B(in)(nG^A(out))* and p^B(out)(nG^A(out)) 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^B(in) from the data, whereas reverse time (correlation) extrapolation uses only the outgoing waves p^B(out).


As mentioned above, it follows from the derivations of integral reciprocity that the interferometric identity in equation 3 only holds for the correct extrapolator G^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 p^B fields, using any trial model m that yields G^AG^A would result in  
As such, we can immediately define the frequency-domain metric as  
or its time-domain counterpart  
where F1 is the inverse Fourier transform. Here, we use I^=I^(p^B,G^A(m)) to denote the explicit dependency of the interferometric identity on the boundary data p^B and of the extrapolator on the model parameters G^A(m)). According to the exact interferometric relation in equation 3, these metrics will be minimum (i.e., equal to zero) when m is the model that produces the correct G^A(m))=G^A. Waveform inversions for the model m can then be pursued by either  
based on the frequency-domain equation 6, or  
based on the time-domain equation 7. In practice, in the presence of noisy finite-bandwidth data, we expect the frequency- and time-domain 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.


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 xB 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 panels — here, 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 Jt(xA,xB) in equation 7, stacked over time, for xA 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 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 (Figure 2a). 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.


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 between 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 multiple 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 target-volume 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 target-oriented 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 time-lapse 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.


We thank Schlumberger for supporting this research and for granting us permission to publish it. In addition, we thank A. Curtis, D.-Jan van Manen, J. Robertsson, and K. Wapenaar for discussions that contributed to the ideas herein. This paper has greatly benefited from the reviews by D.-J. van Manen and by an anonymous reviewer.