\title{Theory of differential offset continuation}\author{Sergey Fomel}\maketitle\begin{abstract} I introduce a partial differential equation to describe the  process of prestack reflection data transformation in the offset,  midpoint, and time coordinates. The equation is proved theoretically  to provide correct kinematics and amplitudes on the transformed  constant-offset sections. Solving an initial-value problem with the  proposed equation leads to integral and frequency-domain offset  continuation operators, which reduce to the known forms of dip  moveout operators in the case of continuation to zero offset.  \end{abstract}\section{Introduction}  The Earth subsurface is three-dimensional, while seismic reflectiondata from a multi-coverage acquisition belong to a five-dimensionalspace (time, 2-D offset, and 2-D midpoint coordinates). This factalone indicates the additional connection that existsin the data space. I show in this paper that it is possible, undercertain assumptions, to express this connection in a concisemathematical form of a partial differential equation. The theoreticalanalysis of this equation allows us to explain and predict the datatransformation between different offsets.  The partial differential equation, introduced in thispaper\footnote{To my knowledge, the first derivation of the revisedoffset continuation equation was accomplished by Joseph Higginbothamof Texaco in 1989.  Unfortunately, Higginbotham's derivation neverappeared in the open literature.}, describes the process of\emph{offset continuation}, which is a transformation of common-offsetseismic gathers from one constant offset to another\cite[]{GPR30-06-08130828}. \cite{GEO61-06-18461858} identified offsetcontinuation (OC) with a whole family of prestack continuationoperators, such as shot continuation \cite[]{SEG-1993-0673}, dipmoveout as a continuation to zero offset \cite[]{DMObook}, andthree-dimensional azimuth moveout \cite[]{GEO63-02-05740588}. Anintuitive introduction to the concept of offset continuation ispresented by \cite{TLE20-01-02100213}. A general data mappingprospective is developed by \cite{GPR48-01-01350162}.   As early as in 1982, Bolondi et al.  came up with the idea ofdescribing offset continuation and dip moveout (DMO) as a continuousprocess by means of a partial differential equation\cite[]{GPR30-06-08130828}.  However, their approximate differentialoperator, built on the results of Deregowski and Rocca's classic paper\cite[]{GPR29-03-03740406}, failed in the cases of steep reflector dipsor large offsets.  \cite{Hale.sepphd.36} writes:\begin{quote}  The differences between this algorithm [DMO by Fourier transform]  and previously published finite-difference DMO algorithms are  analogous to the differen\-ces be\-t\-ween  fre\-qu\-en\-cy-wave\-num\-ber  \cite[]{GEO43-01-00230048,GEO43-07-13421351} and  fi\-nite-dif\-fe\-rence \cite[]{Claerbout.fgdp.76} algorithms for  migration. For example, just as finite-difference migration  algorithms require approximations that break down at steep dips,  finite-difference DMO algorithms are inaccurate for large offsets  and steep dips, even for constant velocity.\end{quote}Continuing this analogy, we can observe that both finite-differenceand frequency-domain migration algorithms share a common origin: thewave equation. The new OC equation, presented in this paper and validfor all offsets and dips, plays a role analogous to that of the waveequation for offset continuation and dip moveout algorithms. Amultitude of seismic migration algorithms emerged from the fundamentalwave-propagation theory that is embedded in the wave equation.Likewise, the fundamentals of DMO algorithms can be traced to the OCdifferential equation.In the first part of the paper, I prove that the revised equationis, under certain assumptions, kinematically valid. This means thatwavefronts of the offset continuation process correspond to thereflection wave traveltimes and correctly transform between differentoffsets.  Moreover, the wave amplitudes are also propagated correctlyaccording to the \emph{true-amplitude} criterion \cite[]{GEO58-01-00470066}. In the second part of the paper, I relate the offset continuationequation to different methods of dip moveout. Considering DMO as acontinuation to zero offset, I show that DMO operators can be obtainedby solving a special initial value problem for the OCequation.  Different known forms of DMO \cite[]{DMObook} appear as specialcases of more general offset continuation operators.The companion paper \cite[]{GEO68-02-07330744} demonstrates apractical application of differential offset continuation to seismicdata interpolation.\section{Introducing the offset continuation equation}%%%%%%%%%%%%%%%%%%%%%%%%%%%% Most of the contents of this paper refer to the following linearpartial differential equation:\begin{equation}h \, \left( {\partial^2 P \over \partial y^2} - {\partial^2 P \over \partialh^2} \right) \, = \, t_n \, {\partial^2 P \over {\partial t_n \,\partial h}} \,\,\, . \label{eqn:OCequation} \end{equation}Equation~(\ref{eqn:OCequation}) describes an {\em artificial}(non-physical) process of transforming reflection seismic data$P(y,h,t_n)$ in the offset-midpoint-time domain. Inequation~(\ref{eqn:OCequation}), $h$ stands for the half-offset($h=(r-s)/2$, where $s$ and $r$ are the source and the receiversurface coordinates), $y$ is the midpoint ($y=(r+s)/2$), and $t_n$ is the timecoordinate after normal moveout correction is applied:\begin{equation}\label{eqn:tnmo}t_n=\sqrt{t^2-{4 \, h^2 \over v^2}}\;.\end{equation}The velocity $v$ is assumed to be known a priori.Equation~(\ref{eqn:OCequation}) belongs to the class of linearhyperbolic equations, with the offset $h$ acting as a time-likevariable. It describes a wave-like propagation in the offsetdirection.\subsection{Proof of validity}A simplified version of the ray method technique \cite[]{cerveny,babich}can allow us to prove the theoretical validity ofequation~(\ref{eqn:OCequation}) for all offsets and reflector dips byderiving two equations that describe separately wavefront (traveltime)and amplitude transformation.  According to the formal ray theory,the leading term of the high-frequency asymptotics for a reflectedwave recorded on a seismogram takes the form\begin{equation}   P\left(y,h,t_n\right) \approxA_n(y,h)\,R_n\left(t_n-\tau_n(y,h)\right) \;,\label{eqn:raymethod} \end{equation}  where $A_n$ stands for the amplitude, $R_n$ is the wavelet shape ofthe leading high-frequency term, and $\tau_n$ is the traveltime curveafter normal moveout. Inserting~(\ref{eqn:raymethod}) as a trialsolution for~(\ref{eqn:OCequation}), collecting terms that have thesame asymptotic order (correspond to the same-order derivatives of thewavelet $R_n$), and neglecting low-order terms, we arrive at the set oftwo first-order partial differential equations:\begin{equation}h \, \left[     {\left( \partial \tau_n \over \partial y \right)}^2 -                 {\left( \partial \tau_n \over \partial h \right)}^2     \right] = \, - \, \tau_n \, {\partial \tau_n \over \partial h} \,\,\,,  \label{eqn:eikonal} \end{equation}   \begin{equation}\left( \tau_n - 2h \, {\partial \tau_n \over {\partial h}} \right)\, {\partial A_n \over \partial h} + 2h {\partial \tau_n \over \partialy}   {\partial A_n \over \partial y} + h A_n \left( {\partial^2 \tau_n\over {\partial y^2}} - {\partial^2 \tau_n \over {\partial h^2}} \right) \, = \, 0 \,\,\,.\label{eqn:transport} \end{equation}Equation~(\ref{eqn:eikonal}) describes the transformation oftraveltime curve geometry in the OC process analogously to how theeikonal equation describes the front propagation in the classic wavetheory.  What appear to be wavefronts of the wave motion described byequation~(\ref{eqn:OCequation}) are traveltime curves of reflectedwaves recorded on seismic sections.  The law of amplitudetransformation for high-frequency wave components related to thosewavefronts is given by equation~(\ref{eqn:transport}).  In terms ofthe theory of partial differential equations,equation~(\ref{eqn:eikonal}) is the characteristic equationfor~(\ref{eqn:OCequation}).\subsection{Proof of kinematic equivalence}%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%In order to prove the validity of equation~(\ref{eqn:eikonal}), it isconvenient to transform it to the coordinates of the initial shotgathers: $s=y-h$, $r=y+h$, and $\tau = \sqrt{\tau_n^2+{{4h^2} \over    {v^2}}}$. The transformed equation takes the form\begin{equation}\left( \tau^2 + {{(r-s)^2} \over {v^2}} \right) \left( {\partial \tau\over \partial r} -   {\partial \tau \over \partial s} \right) = 2 \, (r-s) \,\tau \left( {1 \over {v^2}} - {\partial \tau \over \partial r}{\partial \tau \over \partial s} \right) \,\,\,.\label{eqn:SCeikonal} \end{equation}Now the goal is to prove that any reflection traveltime function$\tau(r,s)$ in a constant velocity medium satisfies equation~(\ref{eqn:SCeikonal}). Let $S$ and $R$ be the source and the receiver locations, and $O$ be areflection point for that pair.  Note that the incident ray $SO$ andthe reflected ray $OR$ form a triangle with the basis on the offset$SR$ ($l=|SR|=|r-s|$).  Let $\alpha_1$ be the angle of $SO$ from thevertical axis, and $\alpha_2$ be the analogous angle of $RO$ (Figure\ref{fig:ocoray}). The law of sines gives us the following explicitrelationships between the sides and the angles of the triangle $SOR$:\begin{eqnarray}|SO|\,=\,|SR|\, {\cos{\alpha_2} \over\sin{\left(\alpha_2-\alpha_1\right)}} \,\,\,, \label{eqn:triangle1} \\ |RO|\,=\,|SR|\, {\cos{\alpha_1} \over\sin{\left(\alpha_2-\alpha_1\right)}} \,\,\,.\label{eqn:triangle2} \end{eqnarray} Hence, the total length of the reflected ray satisfies\begin{equation}v \tau = |SO|+|RO|=|SR|\,  {{\cos{\alpha_1}+ \cos{\alpha_2}} \over\sin{\left(\alpha_2-\alpha_1\right)}} = |r-s|\,{\cos{\alpha} \over\sin{\gamma}} \,\,\,.\label{eqn:length} \end{equation}Here $\gamma$ is the reflection angle ($\gamma = (\alpha_2 -\alpha_1)/2$), and $\alpha$ is the central ray angle ($\alpha =(\alpha_2 + \alpha_1)/2$), which coincides with the local dip angle ofthe reflector at the reflection point.  Recalling the well-knownrelationships between the ray angles and the first-order traveltimederivatives\begin{eqnarray}{{\partial \tau} \over {\partial s}} \,=\,{ {\sin{\alpha_1}} \over{v}} \,\,\,,\label{eqn:snell1}\\{{\partial \tau} \over {\partial r}} \,=\, {{\sin{\alpha_2}} \over {v}}  \label{eqn:snell2}\,\,\,,\end{eqnarray}  we can substitute~(\ref{eqn:length}),~(\ref{eqn:snell1}), and~(\ref{eqn:snell2}) into (\ref{eqn:SCeikonal}), which  leads to the simple trigonometric equality\begin{equation}\cos^2{\left( {\alpha_1 + \alpha_2} \over 2 \right)} +  \sin^2{\left( {\alpha_1 - \alpha_2} \over 2 \right)}\, = \, 1 -\sin{\alpha_1} \sin{\alpha_2} \,\,\,. \label{eqn:equality} \end{equation}It is now easy to show that equality~(\ref{eqn:equality}) is true for any$\alpha_1$ and $\alpha_2$, since\[\sin^2{a} - \sin^2 = \sin{(a+b)}\,\sin{(a-b)}\;.\]\inputdir{XFig}\sideplot{ocoray}{height=2.5in}{Reflection rays in a constantvelocity medium (a scheme).}Thus we have proved that equation(\ref{eqn:SCeikonal}), equivalent to~(\ref{eqn:eikonal}), is valid in constantvelocity media independently of the reflector geometry and the offset.This means that high-frequency asymptotic components of the waves,described by the OC equation,are located on the true reflection traveltime curves.  The theory of characteristics can provide other ways to prove thekinematic validity of equation~(\ref{eqn:eikonal}), as described by\cite{me} and \cite{plag}.%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\subsection{Comparison with Bolondi's OC equation}Equation~(\ref{eqn:OCequation}) and the previously published OCequation \cite[]{GPR30-06-08130828} differ only with respect to thesingle term $\partial^2 P \over {\partial h^2}$. However, thisdifference is substantial.  From the offset continuation characteristic equation(\ref{eqn:eikonal}), we can conclude that the first-order traveltimederivative with respect to offset decreases with decreasingoffset. The derivative equals zero at the zero offset, as predicted by theprinciple of reciprocity (the reflection traveltime has to be an {\em  even} function of offset). Neglecting $\left({\partial \tau_n} \over  {\partial h}\right)^2$ in (\ref{eqn:eikonal}) leads to thecharacteristic equation\begin{equation} h \,   {\left( \partial \tau_n \over \partial y \right)}^2   = \, - \, \tau_n \, {\partial \tau_n \over \partial h}\;, \label{eqn:ITeikonal} \end{equation}which corresponds to the approximate OC equation of\cite{GPR30-06-08130828}. The approximate equation has the form\begin{equation}h \, {\partial^2 P \over \partial y^2} \, = \, t_n \, {\partial^2 P\over {\partial t_n \, \partial h}}\;.\label{eqn:bolondi} \end{equation}Comparing equations~(\ref{eqn:ITeikonal}) and (\ref{eqn:eikonal}), wecan note that approximation (\ref{eqn:ITeikonal}) is valid only if\begin{equation}{\left( \partial \tau_n \over \partial h \right)}^2 \, \ll\, {\left(\partial \tau_n \over \partial y \right)}^2 \,\,\,. \label{eqn:condition} \end{equation}To find the geometric constraints implied by inequality(\ref{eqn:condition}), we can express the traveltime derivatives ingeometric terms. As follows from expressions (\ref{eqn:snell1}) and(\ref{eqn:snell2}),\begin{eqnarray}\label{eqn:snells1}{{\partial \tau} \over {\partial y}} & = & {{\partial \tau} \over{\partial r}} + {{\partial \tau} \over {\partial s}} \,=\, { {2\sin{\alpha} \cos{\gamma}} \over {v}}\;, \\{{\partial \tau} \over{\partial h}} & = & {{\partial \tau} \over {\partial r}} - {{\partial\tau} \over {\partial s}} \,=\, { {2 \cos{\alpha} \sin{\gamma}} \over{v}}\;.\label{eqn:snells2}\end{eqnarray}Expression (\ref{eqn:length}) allows transforming