
\chapter{Wave-equation inversion in practice}
\label{regu}

\section{Introduction}
Wave-equation inversion can be done either by explicit computation of the Hessian matrix (chapters \ref{theo} and \ref{hess}) or by an implicit method \cite[]{nemeth99,myseg02,KuhlSacchi03,Mclapp122}. The main advantage of the explicit method is that after the Hessian matrix is computed (chapter \ref{hess}), the regularization, weights, preconditioners, or solvers can be changed without having to recompute Green's functions. For example, 
regularization  always requires parameter adjustment, which can demand several iterations of the inversion process. This can be prohibitively expensive in the implicit method.

In this chapter, I compare two different regularization schemes for the wave-equation inversion in the poststack (equation \ref{eqn:damping}), and in the prestack  (equations \ref{eqn:regugene}) image domains. To solve the linear systems in equations \ref{eqn:damping} and \ref{eqn:regugene},  I use conjugate gradients.

I choose the Sigsbee2B model \cite[]{paffenholz:2122} to illustrate the action of the different regularizations. The migration and the Hessian matrix are computed using the exact (stratigraphic) velocity model  (Figure \ref{fig:Sis-vel}), so that the migration is as kinematically correct as possible, creating the ideal conditions for the linearized inversion to work. 


\section{Inversion of the point-scatterers image }

In chapter \ref{hess}, I created a model with 32 point scatterers and convolved it with the Hessian matrix for the Sigsbee2B model. The resulting image (Figure \ref{fig:conv-Sis-new}) is good for illustrating how the inversion works when the Hessian matches the physics used in migration (one-way wave-equation), and in the absence of noise. Ideally the inversion should recover the scatterers; in practice, the performance of the inversion depends on the condition number of the Hessian matrix. A large condition number both increases the number of iterations required, and limits the accuracy to which the solution can be obtained \cite[]{numerical-recipes,practical-optimization}. 

\subsection{Wave-equation inversion in the poststack-image domain }

In this subsection I illustrate how the wave-equation inversion works  when the subsurface offset dimension is not computed. In this situation, the regularization is limited to damping (equation \ref{eqn:damping}), and as was shown in chapter \ref{theo}, no {\it a priori} knowledge of the properties of the reflectivity can be included besides constraining its size.

Figure \ref{fig:conv-Sis-inv} shows the image of the 32 point scatterers obtained by solving the linear system in equation \ref{eqn:damping} (Figure  \ref{fig:conv-env-Sis-inv} shows the envelope). The result is shown after seven iterations of a conjugate-gradient iterative solver. No regularization was applied.  In the well-illuminated areas, the Point Spread Function (PSF) has been reduced in size (compare Figure \ref{fig:conv-Sis-inv} with \ref{fig:conv-Sis-new} and Figure \ref{fig:conv-env-Sis-inv} with \ref{fig:conv-env-Sis-new}). This is the case for all  points located in the sediments to the left of the salt body such as the point with coordinates $ {\bf x}=(31050,11000)$. 

\plot{conv-Sis-inv}{width=3.5in}{Inversion of the image of the 32 point scatterers after seven iterations  of a conjugate-gradient iterative solver with no regularization.  Compare with Figure \ref{fig:conv-Sis-new}, and notice how the amplitudes have been balanced, and the PSF size has been reduced.}

As the point scatters enter a shadow zone, the inversion increases their amplitudes. Compare the point scatterer at position  $ {\bf x}=(34500,16000) $ in Figures \ref{fig:conv-Sis-new} and \ref{fig:conv-Sis-inv}. Unfortunately the performance of the inversion degrades most in this area, where we are most interested in recovering the reflectors. This makes clear the need for regularization. The point located at $ {\bf x}=(35625,11000) $ is a good example; here the inversion has increased the amplitudes, but the noise has also been increased. As we see later with the Sigsbee2B image, there are almost no primaries in this area of the shadow zone.

\plot{conv-env-Sis-inv}{width=3.5in}{Envelope of the inversion of the image of the 32 point scatterers after 7 iterations  of a conjugate-gradient iterative solver and no regularization. Compare with the migration in Figure \ref{fig:conv-env-Sis-new}, and notice how the amplitudes have been balanced, and the PSF size has been reduced.}

Figure \ref{fig:res-ir} shows how the norm of the residuals (${\bf r \cdot r}$) as they change with the number of iterations for different values of regularization parameter $\varepsilon$. The choice of stooping the inversion for $\varepsilon=0.0$ at iteration seven is not arbitrary, since iterative methods have an inherent regularization property when applied to solve the system of equations (like equation \ref{eqn:lse}), and are stopped at the optimal iteration number. \cite{hanke} discuss regularization techniques alternative to Tikhonov regularization. They are base on the {\it semiconvergence} property of the norm of the residuals --convergence in the beginning of the iteration, divergence eventually. In this thesis I only make use of  Tikhonov regularization, but other regularization strategies could be the topic for future research. 


\plot{res-ir}{width=3.5in}{Norm of the residuals as they change with number of iterations  of a conjugate-gradient iterative solver for different regularization parameters $\varepsilon$ in the poststack image domain. The parameter choice $\varepsilon=0.5$ decreases the norm of the residuals, stabilizing the inversion. }

The convergence of the inversion in the poststack-image domain can be improved by adding a damping factor (equation \ref{eqn:damping}).  In this case, $\varepsilon=0.5$ decreases the the norm of the residuals, stabilizing the inversion
The inversion result for the regularization parameter value $\varepsilon=0.5$ is shown in Figures \ref{fig:conv-Sis-inv05} and \ref{fig:conv-env-Sis-inv05}. Compared to that in Figures \ref{fig:conv-Sis-inv} and \ref{fig:conv-env-Sis-inv}, this result is more stable in the shadow zones. Adding a ``global" regularization, however could degrade the performance of the inversion in the well-illuminated areas. This suggests that the use of a spatially variable regularization parameter that would be applied only in the shadow zones should improve the regularize inversion results.% Unfortunately a variable regularization is usually difficult to implement, making the solution model-dependent. 

\plot{conv-Sis-inv05}{width=3.5in}{Inversion of the image of the 32 point scatterers after 200 iterations  of a conjugate-gradient iterative solver with regularization parameter $\varepsilon=0.5$.  Compare with Figures \ref{fig:conv-Sis-new},  notice how the amplitudes have been balanced, and the PSF size has been reduced. Compared with \ref{fig:conv-Sis-inv} the amplitudes in the shadow zone are weaker, but the result is cleaner.}

\plot{conv-env-Sis-inv05}{width=3.5in}{Envelope of the inversion of the image of the 32 point scatterers after 200 iterations  of a conjugate-gradient iterative solver with regularization parameter $\varepsilon=0.5$.  Compare with Figures \ref{fig:conv-Sis-new},  notice how the amplitudes have been balanced, and the PSF size has been reduced. Compared with \ref{fig:conv-Sis-inv} the amplitudes in the shadow zone are weaker, but the result is cleaner.}

\subsection{Wave-equation inversion in the prestack-image domain}

As explained in chapter \ref{theo}, a more sophisticated regularization could be to force the model to fulfill the semblance principle. In the subsurface-offset dimension, this consists of penalizing the energy in the image that is not focused at zero subsurface offset (equation \ref{eqn:regusemb}). To apply this regularization, the model must be expanded to include the subsurface offset dimension.

Adding the new dimension increases the size and the number of nonzero coefficients of the Hessian matrix, increasing the cost of the inversion by an order of magnitude. The regularization becomes more important now since, as the size of the matrix increases, the inversion becomes more susceptible to numerical errors. 

\plot{res-ir-off}{width=3.5in}{Norm of the residuals as they change with the number of iterations  of a conjugate-gradient iterative solver, for different regularization parameters $\varepsilon$  in the prestack image domain. The parameter choice $\varepsilon=3.0$ makes the norm of the residuals decrease, stabilizing the inversion. }

As with inversion in the poststack image domain, the result of the inversion is unstable in the absence of regularization (Figure \ref{fig:res-ir-off}). Figure \ref{fig:conv-Sis-inv-offzero} shows the image obtained in the subsurface-offset domain by solving the linear system in equation \ref{eqn:regugene}, using seven iterations of the conjugate-gradients solver with no regularization. The front panel (to the left) shows the zero subsurface offset, and the side panel (to the right) shows the subsurface-offset gather at position $x=35550$ ft. Compared to Figure \ref{fig:conv-Sis-off} notice how the amplitudes have been balanced everywhere except in the shadow zone. Also notice how the size of the three dimensional PSF has been reduced, mainly in the $x$ and $z$ dimensions.

\plot{conv-Sis-inv-offzero}{width=5.0in}{Inversion of the image of the 32 point scatterers in the subsurface-offset domain after seven iterations of the conjugate-gradient solver and  no regularization.  Compare with the migration (Figure \ref{fig:conv-Sis-off}); notice how the amplitudes have been balanced, and the size of the three dimensional PSF has been reduced, but the result is unstable in the shadow zone. }

The reduction in the size of the PSF in the subsurface-offset dimension can be better judged after transforming the inverted image to the reflection-angle domain by slant stacking \cite[]{SavaFomel:03} (Figure \ref{fig:conv-Sis-inv-angzero}). The front panel shows a constant-angle section (at $18^\circ$), and the side panel shows the reflection-angle gather at position $x=35550$ ft. Comparing this figure with the migration (Figure \ref{fig:conv-Sis-ang}), we see that the gaps in the angle coverage have been partially filled by the inversion. We also see again that the inversion becomes unstable inside the shadow zone. 

\plot{conv-Sis-inv-angzero}{width=5.0in}{Inversion of the  image of the 32 point scatterers in the reflection-angle domain (offset-to-angle transformation), after seven iterations of the conjugate-gradient solver and  no regularization. Compare with Figure \ref{fig:conv-Sis-ang}; notice how the gaps in the angle coverage has been partially filled by the inversion, the amplitudes have been balanced, and the PSF size has been reduced. }

Adding the regularization in the subsurface offset makes the inversion result stable (from the norm of the residuals in Figure \ref{fig:res-ir-off}). Figure \ref{fig:conv-Sis-inv-off}  shows the image of the 32 point scatterers in the subsurface-offset domain, obtained by solving the linear system in equation \ref{eqn:regugene} with the differential semblance as the regularization operator and $\varepsilon=3.0$ (equation \ref{eqn:regusemb}). The front panel shows the zero subsurface-offset image, and the side panel shows the subsurface-offset gather at position $x=35550$ ft. Compared to the migration shown in Figure \ref{fig:conv-Sis-off}, besides the usual benefit of balancing the amplitudes, the width of the PSF has been reduced in the all three dimensions. Also, if we compare the inversion results with and without regularization (Figure  \ref{fig:conv-Sis-inv-off} versus  \ref{fig:conv-Sis-inv-offzero}), we can see a reduction of the PSF in the subsurface-offset dimension. This shows that the regularization is working.

\plot{conv-Sis-inv-off}{width=5.0in}{Inversion of the image of the 32 point scatterers in the subsurface-offset domain after 200 iterations of the conjugate-gradient solver and  regularization parameter $\varepsilon=3.0$.  Compare with the migration in Figure \ref{fig:conv-Sis-off} and the inversion without regularization in Figure \ref{fig:conv-Sis-inv-offzero}; notice how the amplitudes have been balanced, and the size of the three dimensional PSF has been reduced.}

 Figure \ref{fig:conv-Sis-inv-ang} shows the inverted image in the reflection-angle domain. Comparing this result with the migration (Figure \ref{fig:conv-Sis-ang}) and with the inversion with no regularization (Figure \ref{fig:conv-Sis-inv-angzero}), we see that the inversion is stable in the shadow zone. Also, the inversion plus regularization equalizes the amplitudes in the angle gathers, filling the illumination gaps.   

%\plot{conv-env-Sis-inv-off}{width=6in}{Envelope of the inversion of the convolution of the image of the 32 point scatterers in the  subsurfe-offset  domain after 200 iterations and  regularization parameter $\varepsilon=3.0$. Compare with Figure \ref{fig:conv-env-Sis-off}, notice how the amplitueds had been balance and the size of the three dimensional PSF has been reduced.}

\plot{conv-Sis-inv-ang}{width=5.0in}{ Inversion of the image of the 32 point scatterers in the reflection-angle domain (offset-to-angle transformation) after 200 iterations of the conjugate-gradient solver and  regularization parameter $\varepsilon=3.0$.  Compare with Figures \ref{fig:conv-Sis-ang} and \ref{fig:conv-Sis-inv-angzero}; notice that the amplitudes have been balanced, and the PSF size has been reduced. Also, the inversion with regularization equalizes the amplitude in the angle gathers, filling the illumination gaps.}

%\plot{conv-env-Sis-inv-ang}{width=6in}{Envelope of the inversion of convolution of the image of the 32 point scatterers with Hessian matrix of the Sigsbee2B model after 200 iterations and  regularization parameter $\varepsilon=3.0$. Compare with Figures \ref{fig:conv-env-Sis-new} and \ref{fig:conv-env-Sis-inv}, notice how the amplitueds had been balance and the PSF size has been reduced.}

Another interesting comparison can be made between the regularized inversion in the poststack image-domain and the zero subsurface-offset of the regularized inversion in the prestack image-domain. 
Figure \ref{fig:conv-Sis-inv-zoff} shows the zero subsurface offset (front panel of Figure \ref{fig:conv-Sis-inv-off}) at the same scale as Figure \ref{fig:conv-Sis-inv} (poststack image-domain inversion) to allow direct comparison. The regularized inversion in the prestack image-domain shows the same level of amplitudes inside and outside the shadow zone, whereas the regularized inversion in the poststack-image domain shows lower amplitudes inside the shadow zone than outside. This comparison shows that the considerable cost of adding one dimension to the model space allows a regularization that has the potential of yielding better results. 

\plot{conv-Sis-inv-zoff}{width=3.5in}{Zero subsurface-offset of the regularized inversion in the prestack image-domain (front panel of Figure \ref{fig:conv-Sis-inv-off}), at the same scale as Figures \ref{fig:conv-Sis-new}, \ref{fig:conv-Sis-inv}, and \ref{fig:conv-Sis-inv05}.  Notice how the amplitudes have been balanced and the PSF size has been reduced.}

\section{Inversion of the Sigsbee2B data}

The salt body in the Sigsbee2B model produces irregular illumination of the subsurface. Two main shadow zones can be seen by analyzing the diagonal of the Hessian matrix for this model (Figure \ref{fig:diag-Sis}). Inside the shadow zones, the primaries are weak and are contaminated with salt-related multiples. This is a problem we often encounter in field data. I address the multiples problem in Appendix \ref{app-mult}.  

\subsection{Irregular angle illumination}

Figure \ref{fig:mig-stack} shows the angle stack of the shot-profile migration (using a cross-correlation imaging condition) corresponding to the portion of Sigsbee2B model shown inside the target box in Figure \ref{fig:Sis-vel}. Note that the shadow zones in Figure \ref{fig:diag-Sis} coincide with the areas where the amplitudes of the reflectors fade in the migration. From the reflection coefficients shown in Figure \ref{fig:refle-new}, we know that the amplitude of the events should not change as they get closer to the salt. It is this weakening of the amplitude in this region of  Figure \ref{fig:mig-stack} that we would like to address.

\plot{mig-stack}{width=3.5in}{Angle stack of the Sigsbee2B shot-profile migration (Figure \ref{fig:migang}).}
 
\plot{refle-new}{width=3.5in}{Reflection coefficients of the Sigsbee2B model.}

The explanation for this discrepancy can found by analyzing the migration image before stacking. Figure \ref{fig:migang} shows the migration in the reflection-angle domain. The front face corresponds to the image at a $25^\circ$ reflection angle, and the side face corresponds to the reflection-angle gather at a position below the tip of the salt ($x=32300$ ft). Figure \ref{fig:migoff} shows the migration result in the subsurface-offset domain. The front face corresponds to the image at zero subsurface offset, and the side face corresponds to the subsurface-offset gather at position $x=32300$ ft. 

In the shadow zone below the tip of the salt, the reflectors are imaged at a few near reflection angles (not focused at zero subsurface-offset), leaving gaps in the angle illumination. The irregular illumination with angle translates into a biased stack image that is difficult to use for interpretation.

\plot{migang}{width=5.0in}{Sigsbee shot-profile migration (reflection angle) using a cross-correlation imaging condition. The front face corresponds to the image at a $25^\circ$ reflection angle, and the side face corresponds to the reflection-angle gather at  $x=32300$ ft. }

\plot{migoff}{width=5.0in}{Sigsbee shot-profile migration (subsurface-offset) using cross-correlation imaging condition. The front face corresponds to the image at zero subsurface-offset, and the side face corresponds to the subsurface-offset gather at  $x=32300$ ft. }

Another well-known problem in shadow zones are the ``migration smiles" that are created in their proximity. They give the misleading impression that the dip of the reflectors is curving up as they enter the shadow zones. This increases the risk of misinterpreting the geometry of the reservoirs.

\subsection{Processing before inversion}
\label{filt}

By comparing the migrated image in Figure \ref{fig:mig-stack} with the reflection coefficients shown in Figure \ref{fig:refle-new} we see events in the migration that do not correspond to the geology. These are well-known artifacts that correspond to multiple reflections inside the salt body that mask the weak primaries signal that we wish to increase by the inversion. The Sigsbee2B model Hessian (Figure \ref{fig:hmatrix-Sis-new}), and the migration (Figure \ref{fig:mig-stack})  were computed using one-way wave-equation propagators. 
As the result, the multiples were neither included in the modeling nor propagated with the correct velocity/direction during migration. For this reason, the multiples should be treated as noise for the inversion. Later, I will show how they can make the inversion fail to recover the primaries signal in the shadow zones. 

One solution to this problem is to suppress the multiples (events that do not fit the physics of the operator) from the migration. This requires the multiples to be separated from the primaries in such a way, that during removal, the primaries are minimally altered. In Appendix \ref{app-mult}, I give a detailed explanation of the demultiple processing. 
The resulting angle-stack of the demultiple migration can be seen in Figure \ref{fig:mig-stack-filt}. Comparing it with the unfiltered stack (Figure \ref{fig:mig-stack}), we see the effectiveness of the process.

\plot{mig-stack-filt}{width=3.5in}{Angle stack of the Sigsbee2B shot-profile migration after filtering the salt-related multiple reflections (Figure \ref{fig:migang1-filt}).}

\subsection{Inversion in the poststack image-domain}

In this subsection, I solve the linear system in equation \ref{eqn:damping} for two different right-hand-side vectors. One is the  migrated image at zero subsurface offset without filtering (front panel, Figure \ref{fig:migoff1}), and the other is the  migrated image at zero subsurface offset after filtering (front panel, Figure \ref{fig:migoff1-filt}). 

The results of the inversion are shown in Figures \ref{fig:invz} and \ref{fig:invz-filt} (unfiltered and filtered input), and should be compared with the migrations shown in Figures \ref{fig:mig-stack} and \ref{fig:mig-stack-filt} respectively. The inversion with the unfiltered migration is more unstable, needing a higher value of the regularization parameter, which reduces the effectiveness of the inversion outside the shadow zones. Thus, I choose to limit the number of iterations instead of adding a regularization parameter for this case. The result in Figure  \ref{fig:invz} shows the inversion of the unfiltered migration after seven iterations,  without regularization. 

\plot{invz}{width=3.5in}{Inversion in the poststack image domain with the unfiltered migration as the input (equation \ref{eqn:damping})  after seven iterations  of a conjugate-gradient iterative solver and no regularization.}

 The inversion of the filtered migration shows a convergence behavior similar to that of the image of the 32 scatterers (Figure \ref{fig:res-ir}). The choice of regularization parameter, $\varepsilon=2.0$, decreases the norm of the residuals and stabilizes the inversion.

\plot{invz-filt}{width=3.5in}{Inversion in the poststack image domain with a filtered migration as the input (equation \ref{eqn:damping}) with a regularization parameter $\varepsilon=2.0$. }

Overall, the inversion results have more balanced amplitudes, allowing the continuation of the reflector inside the shadow zones with the much-improved kinematics (reducing the migration smiles). The reflectors  also gain vertical and horizontal resolution, particularly seen at the two faults present in the reflectivity model (Figure \ref{fig:refle-new}). They can be better interpreted in the inversion than in the migration. 

Figures \ref{fig:resinvz}, and \ref{fig:resinvz-filt} show the residuals of the inversion with the unfiltered migration, and with the filtered migration, after seven iteration of the conjugate-gradient solver.  As expected, the residuals in the unfiltered inversion become large in the shadow zone, making the conjugate-gradient solver increase the amplitude of the multiples that were not included in the computation of the Hessian. 
On the other hand, the filtering of the multiples before inversion allows the inversion process to boost the amplitude of the primaries in the shadow zone. Unfortunately the multiples were not completely removed by filtering, so the inversion procedure also increases their amplitudes.  

\plot{resinvz}{width=3.5in}{Residuals of the inversion in the poststack image domain with the unfiltered migration as the input (equation \ref{eqn:damping}) after seven iterations, with a regularization parameter $\varepsilon=0.0$. }
\plot{resinvz-filt}{width=3.5in}{Residuals of the inversion in the poststack image domain with a filtered migration as the input (equation \ref{eqn:damping}) after seven iterations, with a regularization parameter $\varepsilon=0.0$. }

\subsection{Inversion in the  prestack image domain}

When imaging subsalt sediments, it is often desirable  to compute prestack image domain gathers. \cite{winbow1430} show that different angles contribute to different portions of the image space depending on the geometry of the salt, so that stacking different angle ranges improve interpretation.
 
To obtain inverted images in the prestack image-domain a more general inversion problem has to be solved as a result of including the subsurface-offset dimension in the model space (equation \ref{eqn:regugene}). The regularization in the subsurface offset complements the inversion while filling the illumination gaps in the reflection-angle dimension, as was shown with the image of the 32 scatterers (Figure \ref{fig:conv-Sis-inv-ang}).  

The filtered migrated image (Figure \ref{fig:migoff1-filt}) was chosen as the input to the inversion in order to reduce the detrimental multiples contamination in the shadow zones. The norm of the residuals show a convergence behavior similar to that in the image of the 32 scatterers (Figure \ref{fig:res-ir-off}). 

 Figure \ref{fig:invoff-zero} shows the inversion result in the subsurface-offset domain as the result of solving the linear system in equation \ref{eqn:regugene} after 20 iterations of the conjugate-gradients solver with no regularization. The front face corresponds to the image at zero subsurface offset, and the side face corresponds to the subsurface-offset gather at position $x=32300$ ft. Figure \ref{fig:invang-zero} shows the result of transforming the image in Figure \ref{fig:invoff-zero} to the reflection-angle domain. The front face corresponds to the image at a $25^\circ$ reflection angle, and the side face corresponds to the reflection-angle gather at a position below the tip of the salt ($x=32300$ ft). Comparing it to the migration (Figure \ref{fig:migang}), we see that some of the gaps in the angle coverage have been partially filled by the inversion.   

\plot{invoff-zero}{width=5.0in}{Sigsbee inversion without regularization in the subsurface-offset domain. The front face corresponds to the image at zero subsurface-offset, and the side face corresponds to the subsurface-offset gather at  $x=32300$ ft. }
\plot{invang-zero}{width=5.0in}{ Sigsbee inversion without regularization in the reflection angle. The front face corresponds to the image at a $25^\circ$ reflection angle, and the side face corresponds to the reflection-angle gather at  $x=32300$ ft. }

Adding the regularization in the subsurface offset makes the inversion result stable. Figure \ref{fig:invoff4}  shows the image in the subsurface-offset domain, obtained by solving the linear system in equation \ref{eqn:regugene} with differential semblance as the regularization operator and $\varepsilon=4.0$ (equation \ref{eqn:regusemb}). The front panel shows the zero subsurface offset, and the side panel shows the subsurface-offset gather at position $x=32300$ ft. Compare with the migration shown in Figure \ref{fig:migoff}; besides the usual benefit of balancing the amplitudes, the PSF has been reduced in all three dimensions. Also if we compare the inversion results with and without regularization (Figure  \ref{fig:invoff-zero} versus  \ref{fig:invoff4}), we see a reduction of the image in the subsurface-offset dimension. 

\plot{invoff4}{width=5.0in}{Sigsbee regularized inversion in the subsurface-offset domain  ($\varepsilon=4.0$). The front face corresponds to the image at zero subsurface-offset, and the side face corresponds to the subsurface-offset gather at  $x=32300$ ft. }

Figure \ref{fig:invang4} shows the inverted image in the reflection-angle domain. Comparing this result with the migration (Figure \ref{fig:migang}) and with the inversion with no regularization (Figure \ref{fig:invang-zero}), we see that the inversion is stable in the shadow zone. Also, the inversion-plus-regularization equalizes the amplitudes in the angle gathers, filling the illumination gaps. Figure \ref{fig:comp-ang} shows three angle gathers at three different locations that correspond to: the migration (top row), the inversion without regularization (middle row), and the regularized inversion (bottom row). Note that the illumination is compensated as we move from migration to regularized inversion. 

\plot{invang4}{width=5.0in}{Sigsbee regularized inversion in the reflection angle ($\varepsilon=4.0$). The front face corresponds to the image at $25^\circ$ reflection angle, and the side face corresponds to the reflection-angle gather at  $x=32300$ ft. }

\plot{comp-ang}{height=7in}{ Angle gathers corresponding to the migration (top row), the inversion without regularization (middle row), and the regularize inversion for $\varepsilon=4.0$ (bottom row) at three locations. a), b), and g) $x=32300$ ft; b), e), and h) $x=33725$ ft; and c), f), and i) $x=36500$ ft.}

\plot{inv-stack-filt0}{width=3.5in}{Angle stack of the inversion without regularization (the front panel of Figure \ref{fig:invang-zero} showed the result for just $25^\circ$ reflection angle).}
\plot{inv-stack-filt4}{width=3.5in}{Angle stack of the regularized inversion with $\varepsilon=4.0$ (the front panel of Figure \ref{fig:invang4} showed the result for just $25^\circ$ reflection angle).}

Finally, Figures \ref{fig:inv-stack-filt0} and  \ref{fig:inv-stack-filt4} show the angle stack without and with regularization. The regularized inversion give a less noisy result inside the shadow zone, since the regularization attenuates events that are not flat in the angle gathers (such as multiples). By comparing this images with the filtered migration (Figure \ref{fig:mig-stack-filt}), we see that the inversion produces a result with better resolution and more balanced amplitudes. As with inversion in the poststack image-domain, the reflectors can be better fallowed into the shadow zones with much-improved kinematics. Compared with the inversion in the poststack image-domain (Figure \ref{fig:invz-filt}), there is a better signal-to-noise ratio in the shadow zones. 

\section{Conclusions}

Regularization can be implemented in the zero subsurface offset, where the options for an operator are limited to a customary damping, or in the subsurface offset, where differential semblance can be applied. The differential semblance, although it is more expensive, can improve the continuity of the reflectors into the shadow zones with a higher signal-to-noise ratio than can the damping in the zero subsurface offset.

In different zones of the image, quality of the inversion result can vary depending on the regularization operators or parameters.
In my proposed methodology the cost of the computation of the Green functions and the Hessian matrix is higher than the inversion per se. This makes it suitable as an interpretation tool because, after computing the Hessian just once, the regularization parameters or the regularization operator can be changed to produce different inversion results with reasonable efficiency.


