%% This file is automatically generated. Do not edit!
\documentclass[12pt]{geophysics}

\usepackage{color}
\usepackage{array}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsfonts}
\usepackage{algorithm}
\usepackage{algorithmic}
\usepackage{rplabels}
\usepackage{placeins}
\usepackage{morefloats}

\setfigdir{Fig}

\begin{document}
\chapter{Application to synthetic 2D data}

\section{Experiment setup}

In this section, I apply the JIRB method to the sedimentary section of the 2D Sigsbee 2A synthetic model (Figure \ref{fig:SigsA-orig-sq}). I display the background subsurface models in slowness squared, except when explicitly indicated. For the purpose of testing the validity of the method, I committed the inverse crime \citep[see discussion in][]{Schuster2017}.

\subsection{Synthetic data}

I smoothed the original model and added a negative Gaussian anomaly in slowness squared (Figure \ref{fig:SigsA-db}), corresponding to a positive velocity anomaly of approximately $850$ ft/s ($259$ m/s) at its apex, obtaining the model shown in Figure \ref{fig:SigsA-btrue}. Such a model extends horizontally for $20000$ ft ($6096$ m) spaced every $75$ ft ($22.86$ m), and $27000$ ft ($8230$ m) in depth spaced every $50$ ft ($15.24$ m), and constitutes the true background model, ${\bf b}_{\text{true}}$, from which I synthesized Born modeled data, ${\bf d}({\bf b}_{\text{true}})$, also using the reflectivity model shown in Figure \ref{fig:SigsA-refl}. The acquisition geometry consists of $54$ split-spread shots spaced every $500$ ft ($152.4$ m). Each shot gather contains $651$ receivers spaced every $75$ ft ($22.86$ m). The recording time is $11.92$ s with $4$ ms time step. Figure \ref{fig:SigsA-born} shows a shot gather that stands at the center of the model.

\multiplot[t]{0}{SigsA-orig-sq}{width=0.75\columnwidth}{Original Sigsbee A velocity model. The dashed rectangle indicates the sedimentary section employed for the numerical experiments. \NR}

\multiplot[t]{4}{SigsA-btrue,SigsA-refl,SigsA-b0,SigsA-db}{width=0.47\columnwidth}{a) True background model (${\bf b}_{\text{true}}$), b) true reflectivity model (${\bf r}_{\text{true}}$), c) wrong background model (${\bf b}_0$), and d) Gaussian anomaly representing the perturbation in the background (${\bf \Delta b}$). \ER}

\multiplot[t]{0}{SigsA-born}{width=0.75\columnwidth}{Single shot gather of the Sigsbee A model. \CR}

\multiplot[t]{0}{SigsA-velrbc-sq}{width=0.75\columnwidth}{Velocity framework for the RBC for the central shot. The dashed square indicates the imaging area. \NR}

For the inversion experiments, I prepared a model identical to ${\bf b}_{\text{true}}$, except for the Gaussian anomaly. This model without the anomaly became the wrong model, ${\bf b}_0$, shown in Figure \ref{fig:SigsA-b0}. In the numerical implementation, the subsurface models include an extension of $26250$ ft ($8001$ m) to the left and the right because the first and last shots stand at the limits of the image space. I obtain such an extension by replicating the last trace at each side. For the implementation of the random boundary conditions (RBC), I created velocity frames, initially within the spatial extension of each shot gather. Next, I included an extension of $1000$ ft ($304.8$ m) to the left and to the right to account for the migration aperture. Next, I incorporated a halo of $1000$ ft ($304.8$ m) that constitutes the support for the RBC. For the upper part of the model, I added extra space before the random halo to prevent artifacts that arise when the source stands at or close to the RBC. Figure \ref{fig:SigsA-velrbc-sq} illustrates the velocity frame for a single shot close to the center of the model space.

\subsection{Inversion preliminaries}

I pre-computed the Gauss-Newton Hessian of FWI, ${\bf H}({\bf b}_0)$, using the point-spread functions (PSFs) shown in Figure \ref{fig:SigsA-psf}. I produced this section by ``seeding" spikes every $15$ gridpoints and applying forward Born modeling followed by adjoint Born modeling. Gridpoints at the center of the PSFs correspond to zero-lag elements of the Hessian matrix, which are part of its main diagonal. Similarly, points picked at non-central locations of the PSFs correspond to off-diagonal elements of the Hessian matrix with the same non-zero lag.

For the Hessian application, I first pick samples corresponding to a determined lag. Next, I apply bilinear interpolation, followed by smoothing using a triangle filter, and then perform the Hessian multiplication for this lag ``on the fly". I move on to the next lag and repeat the operation.

\multiplot[t]{0}{SigsA-psf}{width=0.75\columnwidth}{Point-spread functions computed on the sedimentary section of Sigsbee A. I gained and clipped the section for visualization purposes. \CR}

I had previously performed conventional WEMVA tests to invert for the background model component independently. During these tests, I incorporated a gaining weighting factor of the form ${\bf E}=\text{diag}(z_0^n,z_1^n,z_2^n,...)$ (where $z_i$ represents depth levels, and $n$ is a positive integer), and applied it to the migration image to improve the inversion results (see Appendix \ref{appendix:prec}). I experimentally found $n=2.5$ as the best value for the power coefficient. I incorporated this weighting factor to the WEMVA operator member of the JIRB implementation. Likewise, I also incorporated B-splines to precondition the gradient.

After including both weighting and preconditioning, the JIRB objective function becomes

\begin{equation}
    \Phi(\mathbf{r}, \mathbf{p}) = \frac{1}{2} \left \Vert {\bf E} \big[\mathbf{H}(\mathbf{b}_0)\mathbf{r}
    - \mathbf{I}(\mathbf{Bp}) \big]\right \Vert _2^2
    - \frac{\lambda^2}{2} \left \Vert {\bf E} \mathbf{I}(\mathbf{Bp}) \right \Vert_2^2,
\label{eq:jirb_final_prec}
\end{equation}

where the inverted background component model can be recovered using ${\bf b}_{\text{inv}}={\bf B}{\bf p}_{\text{inv}}$. The JIRB experiments used the nonlinear steepest descent method. The inversion tests ran until the solver was unable to find an adequate step size to update the model.

I plotted the inverted reflectivity images displayed in the following section at the same color scale within the same figure set. For comparison purposes, I also computed conventional LWI (equation \eqref{eq:lwi}) using the correct background model. Such inversion ran for $30$ iterations.

For the background models, I compare perturbations in the background, ${\bf \Delta b}={\bf b}-{\bf b}_0$, rather than the absolute background models, for better visualization. I plotted the results at the same color scale within each figure set, except when specified.

\section{Numerical results}

\subsection{Comparing JIRB with LWI} \label{ssec:jirb-vs-lwi}

 Figure \ref{fig:SigsA-hessT-inv} shows the LWI reflectivity using the correct background model. Note how the amplitude values closely resemble those of the true reflectivity section in Figure \ref{fig:SigsA-true-refl}. Figure \ref{fig:SigsA-hessT-inv-zoom,SigsA-true-refl-zoom} shows a zoom into the depth interval of $12000-22000$ ft for a better visualization. Note that seismic features at LWI have the correct kinematics. I hereafter use this zoom level for the rest of this section. I also regarded the LWI reflectivity model as the desired result.

\multiplot[t]{2}{SigsA-hessT-inv,SigsA-true-refl}{width=0.71\columnwidth}{Comparison between a) conventional LWI with the true background model, and b) true reflectivity. \CR}

\multiplot[t]{2}{SigsA-hessT-inv-zoom,SigsA-true-refl-zoom}{width=0.71\columnwidth}{Zoom at sections in Figure \ref{fig:SigsA-hessT-inv,SigsA-true-refl}: a) LWI with true background model; b) true reflectivity. \CR}

Figure \ref{fig:SigsA-hess-inv-zoom,SigsA-true-refl-zoom} shows the LWI reflectivity using the wrong background model compared with the true reflectivity. As expected, reflections at the center of the image get pulled upwards because of the absence of the Gaussian anomaly (note in particular the central spike relative to the grid lines in both cases). There is also unfocused energy in the reflectors of the lower part of the section.

\multiplot[t]{2}{SigsA-hess-inv-zoom,SigsA-true-refl-zoom}{width=0.71\columnwidth}{Zoom at a) LWI with wrong background model, and b) true reflectivity. \CR}

Next, I display each JIRB result organized in sets of figures as follows: a) JIRB reflectivity vs. desired LWI reflectivity, both plotted with the same color scale and independent color scale, the latter to better appreciate reflectivity details otherwise obscured by the section whose amplitude dominates; b) JIRB perturbation in the background vs. true perturbation in the background; and c) trace comparisons of JIRB, LWI, and true reflectivities (I extracted the central trace for such comparisons).

Figure \ref{fig:SigsA-jirb-r-12-zoom,SigsA-lwi-b-12-zoom,SigsA-jirb-r-12-zoom-IG,SigsA-lwi-b-12-zoom-IG} shows the inverted reflectivity model obtained with JIRB using $\lambda=12$. Figure \ref{fig:SigsA-jirb-b-12,SigsA-dB} shows the comparison of the corresponding JIRB perturbation in the background with the true perturbation in the background. Here and hereafter I display the perturbations rather than the background models for the sake of better visualization. I chose this example as the baseline because I found that it approximately represents the minimum value of $\lambda$ at which we can obtain the desired correction. Note how the reflectors moved back to their correct position, and the reflection events are almost kinematically equivalent (Figures \ref{fig:SigsA-jirb-r-12-zoom-IG} and \ref{fig:SigsA-lwi-b-12-zoom-IG}). Compare with the LWI reflectivity. The correction in the reflectivity progresses at the same time as the background model gets corrected during the inversion. Nonetheless, the amplitudes are higher than the desired reflectivity (Figures \ref{fig:SigsA-jirb-r-12-zoom} and \ref{fig:SigsA-lwi-b-12-zoom}), as further corroborated by the trace comparison for JIRB, LWI, and the true reflectivity, shown in Figure \ref{fig:SigsA-jirb-12-traces}. One solution to this problem can be the incorporation of the minimum norm regularization of the reflectivity into the JIRB objective function. The main drawback is the introduction of an additional trade-off parameter that we ought to estimate.

%%%%%% lambda = 12

\multiplot[t]{4}{SigsA-jirb-r-12-zoom,SigsA-lwi-b-12-zoom,SigsA-jirb-r-12-zoom-IG,SigsA-lwi-b-12-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=12$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) but plotted with independent color scales. \CR}

\multiplot[t]{2}{SigsA-jirb-b-12,SigsA-dB}{width=0.71\columnwidth}{a) JIRB perturbation on the background using $\lambda=12$. b) True perturbation in the background. \CR}

\multiplot[t]{0}{SigsA-jirb-12-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=12$), LWI, and true reflectivity. \CR}

Elemental analysis of the objective function \eqref{eq:jirb_final_prec} promptly drives us to deduce that higher values of $\lambda$ underweight the reflectivity-related gradient contribution. Conversely, smaller values of $\lambda$ exacerbate the problem of higher amplitude, as Figures \ref{fig:SigsA-jirb-r-10-zoom,SigsA-lwi-b-10-zoom,SigsA-jirb-r-10-zoom-IG,SigsA-lwi-b-10-zoom-IG} and \ref{fig:SigsA-jirb-10-traces} exemplify for $\lambda=10$, notwithstanding that the JIRB perturbation in the background appears to be unaffected (Figure \ref{fig:SigsA-jirb-b-10,SigsA-dB}). Going for smaller values, e.g., $\lambda=5$, further underweights the background-related gradient to the point that it would not be able to correct the reflectivity inaccuracies. Figures \ref{fig:SigsA-jirb-r-5-zoom,SigsA-lwi-b-5-zoom,SigsA-jirb-r-5-zoom-IG,SigsA-lwi-b-5-zoom-IG} and \ref{fig:SigsA-jirb-5-traces} illustrate this problem. Note that the central spike does not reach the correct position. This issue is due to the comparatively poor recovery of the background component as observed in Figure \ref{fig:SigsA-jirb-b-5,SigsA-dB}, which the reader should compare to Figure \ref{fig:SigsA-jirb-b-12}).

%%%%%% lambda = 5

\multiplot[t]{4}{SigsA-jirb-r-5-zoom,SigsA-lwi-b-5-zoom,SigsA-jirb-r-5-zoom-IG,SigsA-lwi-b-5-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=5$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) but plotted with independent color scales. \CR}

\multiplot[t]{2}{SigsA-jirb-b-5,SigsA-dB}{width=0.75\columnwidth}{a) JIRB perturbation on the background using $\lambda=5$. b) True perturbation in the background. \CR}

\multiplot[t]{0}{SigsA-jirb-5-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=5$), LWI, and true reflectivity. \CR}

%%%%%% lambda = 10

\multiplot[t]{4}{SigsA-jirb-r-10-zoom,SigsA-lwi-b-10-zoom,SigsA-jirb-r-10-zoom-IG,SigsA-lwi-b-10-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=10$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) but plotted with independent color scales. \CR}

\multiplot[t]{2}{SigsA-jirb-b-10,SigsA-dB}{width=0.75\columnwidth}{a) JIRB perturbation on the background using $\lambda=10$. b) True perturbation in the background. \CR}

\multiplot[t]{0}{SigsA-jirb-10-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=10$), LWI, and true reflectivity. \CR}

Going for higher values of $\lambda$, I experimented using $15, 20,$ and $25$. I also used $\lambda=30$, but the sections are not very different from $\lambda=25$, hence they are not shown. Figures \ref{fig:SigsA-jirb-r-15-zoom,SigsA-lwi-b-15-zoom,SigsA-jirb-r-15-zoom-IG,SigsA-lwi-b-15-zoom-IG}, \ref{fig:SigsA-jirb-r-20-zoom,SigsA-lwi-b-20-zoom,SigsA-jirb-r-20-zoom-IG,SigsA-lwi-b-20-zoom-IG}, and \ref{fig:SigsA-jirb-r-25-zoom,SigsA-lwi-b-25-zoom,SigsA-jirb-r-25-zoom-IG,SigsA-lwi-b-25-zoom-IG} respectively show the corresponding reflectivity sections, while Figures \ref{fig:SigsA-jirb-15-traces}, \ref{fig:SigsA-jirb-20-traces}, and \ref{fig:SigsA-jirb-25-traces} show the corresponding trace comparisons. I omitted the perturbations in the background because they do not significantly differ from the case of $\lambda=12$. We observe that the reflectivity amplitudes of the JIRB results get reduced as $\lambda$ values increase, as anticipated. The best result is arguably obtained using $\lambda=25$. Nonetheless, note in trace comparisons of Figure \ref{fig:SigsA-jirb-25-traces} that, at shallow depths, amplitudes get underestimated, while at deep depths amplitudes still get overestimated.

One problem derived from increasing the ${\lambda}$ value is that the solver stops at earlier iterations. As a consequence, the JIRB images lose resolution in comparison to LWI. Conversely, keeping higher resolution by reducing the $\lambda$ value yields the higher-than-ideal amplitudes observed before. For a better appreciation of this issue, Figure \ref{fig:SigsA-jirb-5-traces-spec,SigsA-jirb-10-traces-spec,SigsA-jirb-15-traces-spec,SigsA-jirb-20-traces-spec,SigsA-jirb-25-traces-spec,SigsA-jirb-30-traces-spec} shows the comparison of the amplitude spectra of the JIRB and LWI traces plotted above. From these plots, we can arguably choose the $\lambda=20$ or the $\lambda=25$ reflectivity result as the best compromise between amplitude accuracy and seismic resolution.

%%%%%% lambda = 15

\multiplot[t]{4}{SigsA-jirb-r-15-zoom,SigsA-lwi-b-15-zoom,SigsA-jirb-r-15-zoom-IG,SigsA-lwi-b-15-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=15$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) with independent color scales. \CR}

\multiplot[t]{0}{SigsA-jirb-15-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=15$), LWI, and true reflectivity. \CR}

%%%%%% lambda = 20

\multiplot[t]{4}{SigsA-jirb-r-20-zoom,SigsA-lwi-b-20-zoom,SigsA-jirb-r-20-zoom-IG,SigsA-lwi-b-20-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=20$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) with independent color scales. \CR}

\multiplot[t]{0}{SigsA-jirb-20-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=20$), LWI, and true reflectivity. \CR}

%%%%%% lambda = 25

\multiplot[t]{4}{SigsA-jirb-r-25-zoom,SigsA-lwi-b-25-zoom,SigsA-jirb-r-25-zoom-IG,SigsA-lwi-b-25-zoom-IG}{width=0.47\columnwidth}{a) JIRB reflectivity using $\lambda=25$. b) LWI reflectivity with true background model. c) and d) Same as a) and b) with independent color scales. \CR}

\multiplot[t]{0}{SigsA-jirb-25-traces}{width=0.99\columnwidth}{Trace comparison between JIRB ($\lambda=25$), LWI, and true reflectivity. \CR}

%%%%%% Spectra of traces

\multiplot[t]{6}{SigsA-jirb-5-traces-spec,SigsA-jirb-10-traces-spec,SigsA-jirb-15-traces-spec,SigsA-jirb-20-traces-spec,SigsA-jirb-25-traces-spec,SigsA-jirb-30-traces-spec}{width=0.47\columnwidth}{Amplitude spectra of the trace comparisons for a) $\lambda=5$, b) $\lambda=10$, c) $\lambda=15$, d) $\lambda=20$, e) $\lambda=25$, and f) $\lambda=30$. \CR}

\FloatBarrier

\subsection{Using the JIRB background model for LWI}

Similar to (but less sophisticated than) reflection FWI \citep[e.g.][]{Xu2012}, JIRB incorporates the reflectivity into the inversion. Therefore, we can expect to recover a better estimation of the background model in comparison to WEMVA, where the background model is estimated alone. In this subsection, I performed the corresponding tests. The inverted background models using JIRB and WEMVA are employed to generate the migration images and the PSFs required by LWI in each case. I used the JIRB background model corresponding to $\lambda=25$. Both JIRB- and WEMVA-based LWI tests ran for 30 iterations.

Figures \ref{fig:SigsA-jirb-b-25} and \ref{fig:SigsA-wemva-b} show the perturbations in the background obtained using JIRB and WEMVA, respectively. The difference between both appears to be negligible, although the comparison of the color scales shows that the JIRB perturbation has higher amplitude.
Figures \ref{fig:SigsA-lwi-after-wemva-inv} and \ref{fig:SigsA-lwi-after-jirb-25-inv} show the reflectivity inversions for WEMVA and JIRB background models, respectively. I plotted the sections using the same color scale. Note how the JIRB result has more focused events, despite the small differences between the background models. The JIRB reflectivity also exhibits higher wavenumber content, as observed in the better resolved weak layers and confirmed by the traces comparison (Figure \ref{fig:SigsA-lwis-25-traces,SigsA-lwis-25-traces-spec}). However, noise arises in the form of high-dipping events crossing through the layers. In tests using the JIRB background model for other values of $\lambda$ (not shown here), there was no significant difference in the estimated reflectivity.

\multiplot[t]{2}{SigsA-jirb-b-25,SigsA-wemva-b}{width=0.71\columnwidth}{a) Perturbation in the background obtained using JIRB. b) Perturbation in the background obtained using WEMVA. \CR}

\multiplot[t]{2}{SigsA-lwi-after-wemva-inv,SigsA-lwi-after-jirb-25-inv}{width=0.71\columnwidth}{a) LWI reflectivity using the WEMVA background model. b) LWI reflectivity using the JIRB background model. \CR}

\multiplot[t]{2}{SigsA-lwis-25-traces,SigsA-lwis-25-traces-spec}{width=0.7\columnwidth}{a) Trace comparison between LWI obtained with the WEMVA background model, LWI obtained with the JIRB ($\lambda=25$) background model, and the true reflectivity. b) Spectra comparison of the LWI for WEMVA AND JIRB background models. \CR}

\section{Conclusions}

In this chapter, I showed the application of the JIRB method to a 2D synthetic dataset obtained from the sedimentary section of Sigsbee A.

The direct implementation of the joint inversion corrects the inaccuracies in the reflectivity at the same time that the background model is estimated. An adequate choice of the parameter $\lambda$ is crucial to have the better trade-off between obtaining higher-than-optimal amplitudes and seismic resolution. The main limitation is the possibility of not running enough iterations to obtain the desired resolution, even for small values of $\lambda$.

I also tested conventional LWI with the JIRB background model in comparison to the WEMVA background model. There is a small, yet noticeable improvement in crispness and seismic resolution despite the small difference between the JIRB and the WEMVA background models.
\end{document}
