%\documentclass{segabs}
%\documentclass[paper,twocolumn,twoside]{geophysics}
%\documentclass[manuscript,revised]{geophysics}
%\newcommand{\rs}[1]{\mathstrut\mbox{\scriptsize\rm #1}}
%\newcommand{\rr}[1]{\mbox{\rm #1}}
%\usepackage{amsmath}
%\usepackage{multicol}
%\usepackage{placeins}
% % figure directory
%\begin{document}
%\renewcommand{\thefootnote}{\fnsymbol{footnote}} 
 
\lefthead{}
\righthead{}
\footer{July--2009}
\title{Migration as an adjoint of up-ward continuation modeling in shot-profile domain}
\author{}
\maketitle

\section{upward-continuation modeling in the shot-profile domain}
Starting from the cross-correlation imaging condition in the shot-profile migration
\begin{eqnarray}
r({\bf x},{\bf h}) = \sum_{{\bf x}_s}\sum_{\omega} p^*({\bf x}-{\bf h},{\bf x}_s,\omega) u({\bf x}+{\bf h},{\bf x}_s,\omega), \label{eq:imgcond}
\end{eqnarray}
where $p({\bf x},{\bf x}_s,\omega)$ is the source wavefield for a single frequency $\omega$ at image point ${\bf x}=(x,y,z)$ with the source located at ${\bf x}_s=(x_s,y_s,0)$; $u({\bf x},{\bf x}_s,\omega)$ is the receiver wavefield and ${\bf h}=(h_x,h_y,h_z)$ is the subsurface half-offset, and `$^*$' stands for the complex-conjugate. We can compute the adjoint of equation \ref{eq:imgcond} as:

\begin{subequations} 
\label{subeqn:list} 
\begin{align} 
%\begin{eqnarray}
p({\bf x},{\bf x}_s,\omega) = \sum_{\bf h} u({\bf x}+{\bf h},{\bf x}_s,\omega) r({\bf x},{\bf h}) \label{subeqn:adjimgcond1} \\
u({\bf x},{\bf x}_s,\omega) = \sum_{\bf h} p({\bf x}-{\bf h},{\bf x}_s,\omega) r({\bf x},{\bf h}) . \label{subeqn:adjimgcond2}
%\end{eqnarray}
\end{align} 
\end{subequations} 
In a more compact notation, not explicitly writing the dependencies on $x$,  $y$ and $\bf h$, and for a given depth level, $z$, equation \ref{subeqn:adjimgcond2} can be written as:
\begin{eqnarray}
u^{z}(\omega) = {\bf{P}}^z(\omega) r^z, 
\label{eq:imgconcompact}
\end{eqnarray}
where ${\bf{P}}^z$ is a convolutional matrix. If $r=r({\bf x},h_x)$, for the case of 5 $x$ positions and $3$ subsurface offsets and subsurface-offset interval equals to the CMP interval, equation \ref{eq:imgconcompact} can be explicitly written as:
\begin{equation} 
\left[\begin{array}{l}
u_1\\ 
u_2\\ 
u_3\\ 
u_4\\ 
u_5\\ 
u_6\\ 
u_7\\ 
\end{array}\right]=
\left[\begin{array}{lllllllllllllll}
0 & p_1& 0  & p_3& 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0 \\ 
0 & 0  & 0  & 0  & p_2& 0  & p_4& 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0 \\ 
0 & 0  & 0  & 0  & 0  & p_1& 0  & p_3& 0  & p_5& 0  & 0  & 0  & 0  & 0 \\ 
0 & 0  & 0  & 0  & 0  & 0  & 0  & 0  & p_2& 0  & p_4& 0  & 0  & 0  & 0 \\ 
0 & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & 0  & p_3& 0  & p_5& 0 \\ 
\end{array}\right]
\left[\begin{array}{l}
r_{1,1}\\
r_{1,2}\\
r_{1,3}\\
r_{2,1}\\
r_{2,2}\\
r_{2,3}\\
r_{3,1}\\
r_{3,2}\\
r_{3,3}\\
r_{4,1}\\
r_{4,2}\\
r_{4,3}\\
r_{5,1}\\
r_{5,2}\\
r_{5,3}\\
\end{array}\right]
\end{equation}
where the subscripts of $u$ and $p$ represent $x$ positions. The first and second subscripts of $r$ represent $x$ positions and subsurface offset, respectively.
\par
The recursive upward propagation of $u^z(\omega)$ and recording at the surface are given by the following equations
\begin{eqnarray}
\left\{ \begin{array}{l}
u^{z-1}(\omega)  = T^z_{_\uparrow}(\omega,s)P^z(\omega) r^z\\ 
d(\omega) = u^1(\omega)
\label{eq:upwmod}
\end{array}
\right.
\end{eqnarray}
where $T^z_{_\uparrow}$ is the upward continuation operator, which is function of the slowness, $s$.
\par
bla, bla ...
The matrix representation of equation \ref{eq:upwmod} is
\begin{eqnarray}
\begin{array}{l}
{\underline{\bf u}} = {\bf T}_\uparrow {\bf P} {\bf r}\\ 
{\bf d} = {\bf E}{\underline{\bf u}}
\label{eq:upwmodmatrix}
\end{array}
\end{eqnarray}

\section{Weighted least-squares data fitting}

The optimum reflctivity, $\bf r$, minimizes the objective function
\begin{eqnarray}
J({\bf r,\underline{u}})=h({\bf r,\underline{u}})=\frac{1}{2} \left | \left | {\bf W} \left ( {\bf d - d_{obs}}) \right ) \right | \right | ^2. \\ 
\label{eq:objfun}
\end{eqnarray}

The forward equation reads:
\begin{eqnarray}
F({\bf r,\underline{u}})={\bf {\underline u}} - {\bf T}_\uparrow {\bf P} {\bf r} = 0.
\label{eq:forwmod}
\end{eqnarray}

\begin{itemize}
\item Build the augmented functional
\begin{eqnarray}
\begin{array}{l}
\mathcal{L}(\underline{\lambda},{\bf \underline{u}},{\bf r}) = h({\bf r,\underline{u}}) - \langle \underline{\lambda}, F({\bf r,\underline{u}}) \rangle \\
\mathcal{L}(\underline{\lambda},{\bf \underline{u}},{\bf r}) = \frac{1}{2} \left | \left | {\bf W} \left ( {\bf d - d_{obs}}) \right ) \right | \right | ^2 - \langle \underline{\lambda},{\bf {\underline u}} - {\bf T}_\uparrow {\bf P} {\bf r} \rangle.
\end{array}
\end{eqnarray}
\item Define the adjoint state equations
\begin{eqnarray}
\frac{\partial \mathcal{L}(\underline{\lambda},{\bf \underline{u}},{\bf r})}{\partial {\bf \underline{u}}} = {\bf E}^T {\bf W}^2({\bf d - d_{obs}}) - \lambda = 0, 
\end{eqnarray}
which gives
\begin{eqnarray}
\lambda = {\bf E}^T {\bf W}^2({\bf d - d_{obs}}).
\end{eqnarray}
\item Compute the gradient with respect to $\bf r$
\begin{eqnarray}
\frac{\partial \mathcal{L}(\underline{\lambda},{\bf \underline{u}},{\bf r})}{\partial {\bf r}} = \frac{\partial J}{\partial {\bf r}} = \lambda ^* {\bf T}_\uparrow {\bf P}
\end{eqnarray}
and
\begin{eqnarray}
\nabla J = {\bf P}^* {\bf T}_\downarrow \lambda,
\end{eqnarray}
\end{itemize}
that corresponds to the downward propagation of the residuals and cross-correlation imaging condition of the downward propagated  residuals with the source wavefield.

