% This texfile was generated by executing Latify on file paper_act.tex
\documentclass[12pt]{sep2e}
\begin{document}

\def\figdir{./Fig}
\chapter{Introduction}\label{chapter:Intro}
\def\w{\omega}
\def\xiv{\mbox{\boldmath$\xi$}}
\INT
Passive seismic imaging is a processing methodology which attempts
to extract subsurface information due to unconventional 
sources.  In this thesis, the goal is the production of subsurface
structural images by recording the ambient noisefield of the Earth
with surface arrays of seismometers or geophones.  The images produced
with this technique are directly analogous to those produced with
conventional reflection seismic data.  Manipulating wavefields
with such unconventional sources has been given the general descriptor
in the literature interferometry.  Passive imaging is not
new science.  This statement is true for two reasons.  
\par
First, the idea of imaging the subsurface by correlating the 
recordings from subsurface sources was first  
introduced by \longcite{GEO33-02-02640269}.  That work 
provided a one-dimensional proof that the auto-correlation of time
series collected by a receiver array can produce the
equivalent to a zero-offset time section.  Subsequently, 
\longcite{Zhang.sep.60.85.art}, used planewave decomposition 
to prove the result in 3D over a homogeneous medium.  \longcite{fink}
developed the Green's function of a heterogeneous medium
with acoustic waves via correlation and validated the observational 
theory with an ultrasonic experiment.  \longcite{kees04}, through 
one-way reciprocity, 
proved that cross-correlating traces of the transmission response of
an arbitrary medium synthesizes the complete reflection response,
in the form of shot-gathers, collected in a conventional active source
experiment.  In fact, in recent years, there has been a sudden
increase in interest around interferometry that has generated a wealth
of literature \cite{wapenaar:SI1}.  
\par
Second, interferometry is simply a generalization of the seismic 
problem formulated by earthquake and active-source 
exploration seismology.  On the of the continuum of purposefulness,
active experiments occupy an endpoint where humans 
initiate an experiment over a specific patch of Earth that is 
of interest.  In this case the source function is (hopefully) zero
phase, impulsive, and time zero is precisely known.  Toward the middle of the 
continuum, humans install receivers over an area of interest and 
capitalize on naturally occurring far-field (teleseismic) 
earthquake energy to probe the subsurface ADD BOSTOCK/SHRAGGE REFS.  Here, the source energy is 
areal, and time zero is defined by windowing data beginning at the 
onset of the earthquake arrival.  At the opposite end of 
this spectrum is this work, which asks if the complicated,
low amplitude ambient noisefield in the Earth, about which one has no 
knowledge, can provide sufficient energy to image the subsurface
\cite{artman:SI177}.  However, the physics, and mathematics used to 
manipulate these various data sets are the same.  The more purposeful 
the experiment, many of the details of seismic interferometry can be 
simplified or are implicitly included in the form of the 
acquisition.  The less one knows about the source energy probing the 
subsurface, the more precise the processing needs to be in terms of 
the rigorous formulation of seismic interferometry.
\par
Passive seismology is a very broad field that basically recapitulates 
the entire field of seismology with unconventional sources.  To date,
the literature concerning subsurface structural imaging at the reservoir 
scale has been largely relegated to theoretic and modeling pursuits due 
largely to a lack of suitable data acquisition.  This thesis adds to 
the further theoretic development of this topic by introducing an 
algorithm to directly migrate raw passive data without first correlating 
the traces to synthesize shot gathers.  Of equal importance is the 
presentation of results from the processing of a large volume of 
passive records with both the standard processing sequence and the 
direct migration methodology developed herein.

\section{Development of the correlation theory}
Figure~\ref{fig:sketchP}, simplified for clarity,
shows the basic kinematics exploited in processing transmission
data.  The figure includes two recording stations 
capturing an approximately planar wavefront emerging from a two-layer
subsurface.  Panel ({\tt a}) shows the raypaths associated with the direct
arrival and one reflected both at the free surface and a subsurface
interface.  The second travel path (labeled reflection ray) has the 
familiar kinematics of the reflection seismic experiment if a source
were excited at receiver one.  The transmission
wavefield is depicted in panel ({\tt b}).  Wavelet polarity is
appropriate for direct arrivals and reflection.  The three main
features of transmission wavefields can be appreciated here.  First, the
exact timing of the arrival is unknown.  Second, the phase 
and duration of transmitted energy are unknown and likely
complicated.  Third, if the incident wavetrain coda is long,
arrivals in the transmission record can interfere.
\par
%%\sideplot{figure1}{width=3.3in, height=2.4in}{ ({\tt a}) Approximately
\activesideplot{sketchP}{width=3.3in,height=2.4in}{NR}{({\tt a}) Approximately
  planar arrival with rays showing important propagation paths for
  passive imaging.  ({\tt b}) Idealized traces from a transmission
  wavefield. ({\tt c}) Shot-gather (reflection wavefield) synthesized
  using trace \textit{r1} as the source.  Many details, explained in
  the text, are neglected for simplicity.} 
Choosing trace \textit{r1} as the comparison (source) trace, panel ({\tt c})
depicts the correlation spikes associated with the arrivals in the
data panel ({\tt b}), where $\otimes$ is correlation.  A solid line
with linear move-out is super-imposed across the correlated traces
corresponding to the direct arrival recorded at each receiver
location.  The dashed line on panel ({\tt c}) has hyperbolic moveout.
However, no correlation peak exists on the \textit{r1}$\otimes$
\textit{r1} trace under the hyperbola.  Not drawn, the second arrival
on \textit{r2} will have a counterpart on \textit{r1} from a ray
reaching the free surface farther to the left of the model.  
\par
In fact, the correlations produced from a single plane wave will 
produce another plane wave.  However, each planar reflection is moved 
to the lag-time associated with
a two-way trip from the surface to the reflector.  Correlation 
removes the wait time for the initial arrival and maintains the 
time differences between the direct arrival and reflections.  Summing the 
correlations from a full suite of plane waves builds hyperbolic events 
through constructive and destructive interference.  Analyzing seismic 
data in terms of plane wave constituents is a commonly invoked tool in 
seismic processing.  Summing the correlations from incident plane waves 
is a plane wave superposition process that builds hyperbolas.
\par
%%Passive seismic imaging is predicated on all possible raypaths bouncing 
%%to and from every possible reflector or diffractor.  Cartoons depicting 
%%the experiment always leave something out that
%%causes an inconsistency that needs more raypaths and receivers to
%%explain.  Unfortunately the trend in complexity continues nearly 
%%forever.  Such a complication arises with the inclusion of a second 
%%reflector.  The two reflection rays correlate with each other with a
%%positive coefficient.  The two travel paths share the time through the
%%shallow layer, so they correlate at a lag (or time difference), equal
%%to the two-way travel time through the deep layer.  This correlation
%%would seem to cause an artifact in the synthesized shot gather.
%%\par
%%However, an upcoming interbed multiple event arrives after the direct arrival
%%delayed by the two-way travel-time through the second layer.  This event
%%has the opposite polarity compared to the direct arrival after once
%%changing its propagation direction from $\uparrow$ to $\downarrow$.  The 
%%delay of its arrival at receiver \textit{r2} compared to the direct
%%arrival at receiver \textit{r1} is the same lag as the correlation
%%between the reflectors and has opposite sign.  Therefore, the internal
%%multiple cancels potential artifacts of the correlation.  This shows
%%the importance of modeling transmission data with a two-way 
%%extrapolator.  Without all possible multiples, correlation artifacts
%%will quickly overwhelm the Earth structure within the correlated output. 
%%\par
Cross-correlation of each trace with every other trace handles the
three main difficulties of transmission data: timing, waveform, and
interference.  First, the output of the correlation is in lag units,
that when multiplied by the time sampling interval, provide the
time delays between like events on different traces.  The zero lag of the
correlation is the zero time for the synthesized 
shot-gathers.  Second, each trace records the character and duration
of the incident energy as it is reflected at the surface.  This
becomes the source wavelet analogous to a recorded vibrator 
sweep.  Third, overlapping wavelets are separated by correlation. 

\subsection{Natural phase encoding}
A simple mathematic explanation of how correlation of subsurface
sources returns the conventional active seismic data volume is in
terms of a phase encoding experiment.  Phase encoding shot records for
simultaneous processing was discussed in
\longcite{GEO65-02-04260436}.  The motivation for phase encoding
is to manipulate the independent shot gathers in a predictable manner
such that the entire data volume can be processed (in this case
migrated) simultaneously after a summation to reduce 
dimensionality.  The reduction in dimensionality similarly reduces the
cost of the processing.  To that end, each shot gather is 
convolved with a unique function that will minimize cross-talk
interference between the shots during a summation over the shot axes
in a 3D seismic survey.  Random time series are often used as effective
encoding functions.  Plane wave migration strategies can also be viewed
as a phase encoding migration operation \cite{CONOCO}.  Consider a set 
of arbitrary encoding
functions $a({\bf x};{\bf x}_s,\w)$ for each shot ${\bf x}_s$ 
over the map coordinates ${\bf x}$ with
which shot records are convolved.  Normally, source functions for
active seismic data are bandlimited (in time) impulses, $\delta({\bf x}-{\bf
  x}_s,\w)$, at the source location, and zero elsewhere.   After 
convolving the impulses with encoding functions the source 
wavefield, $\hat{S}$, required for simultaneous migration is
\begin{equation}\label{source}
\hat{S}_{z=0}({\bf x},\w)=\sum_{{\bf x}_s} 
a({\bf x};{\bf x}_s,\w)\delta({\bf x}-{\bf x}_s,\w)\;.
\end{equation}
The five dimensional active seismic data, $R({\bf x},{\bf x}_s,\w)$, 
are similarly convolved with encoding functions and summed to reduce 
its dimensionality
\begin{equation}\label{encode}
\hat{R}_{z=0}({\bf x},\w)=\sum_{{\bf x}_s} 
a({\bf x};{\bf x}_s,\w)R_{z=0}({\bf x};{\bf x}_s,\w)\;.
\end{equation}
If the encoding functions do not correlate at all, migrating the sum of 
the shot gathers, $\hat{R}$, would reproduce the subsurface image desired 
at a substantial cost savings.  Unfortunately, perfect encoding functions 
do not exist, and steps must be taken to minimize the cross-talk introduced
during summation.  This is normally done by summing the results of several 
encoded migrations with different encoding functions.  In this manner, 
the cross-talk incurred by each migration will be out of phase and 
stacking the images will reduce their effects.
\par
For passive seismic records, we can assume unknown and
arbitrarily complicated encoding functions arrive from the
subsurface.  Therefore a transmission record $T({\bf x},\w)$ recorded
by an array at the surface is $\hat{R}_{z=0}$.  To recover 
$R_{z=0}({\bf x};{\bf x}_s,\w)$, we need to invert 
equation~\ref{encode}.  Unfortunately, after the summation in the 
equation, this is not possible.  However, the rigorous development 
of the correlation theory developed in the literature shows that 
if a transmission wavefield due to a single subsurface source, 
$\xi$ (whose location need not be known), can be considered 
independently, it is possible to synthesize $R$ from $T$.  
\par
To calculate the Fourier transform of the reflection response, 
$R({\bf x}_r,{\bf x}_s, \omega)$, of a loss-less (to neglect amplitude
and attenuation issues) subsurface
\longcite{kees04}, proved
\begin{equation}\label{rigor}
2\Re [R({\bf x}_r,{\bf x}_s,\omega)]=
\delta({\bf x}_s-{\bf x}_r)-\int_{\partial\! D_m} 
T  ({\bf x}_r,\xi,\omega)
T^*({\bf x}_s,\xi,\omega)\;
\mbox{d}^2\xiv \;,
\end{equation}
where $^*$ represents conjugation and subscripts $r$ and
$s$ indicate different station locations from a transmission
wavefield.  After correlation $r$ and $s$ acquire the meaning of receiver
and source locations, respectively, associated with an active survey.  The 
RHS represents summing correlations of windows of passive data around
the occurrence of individual sources, $\xi$.  The
symbol $\partial \! D_m$ represents the domain boundary that surrounds
the subsurface region of interest on which the sources are 
located.  Alternatively, many impulses can be substituted for a 
full suite of plane waves emerging from all angles and azimuths 
as in the kinematic explanation above. 
\par
The two main assumptions within the equation are: 1) processing a 
complete set of transmission wavefields (integral), and 2) processing
wavefields due to individual sources, $\xi$.  If these requirements 
are met, the sum of the many correlation volumes does invert 
equation~\ref{encode} and provide a volume of
data that is at least a scaled version (acknowledging potential 
amplitude problems) of the active seismic reflection experiment.  If 
too few sources are captured (and the integral above not completed),
$\tilde{R}$ may not accurately reproduce the kinematics of 
$R({\bf x}_s)$.  For some applications, completing the integral of 
sources around the volume of interest may not be required.  If one is 
specifically interested in body waves probing subsurface structure, 
only point sources below the domain (or relatively flat plane waves)
are required to image the structure, but the synthesized gathers 
will contain some artifacts.  Alternatively, \longcite{louie}, 
described using surface waves generated by cultural noise for 
near-surface refraction velocity analysis.
\par
The phase encoding experiment requires convolving both $S$ and $R$
with encoding functions.  The source wavefield $\tilde{S}$ is the 
sum of all the encoding functions.  To numerically calculate 
$\tilde{S}$, it would be necessary to separate
the initial upcoming energy, $a({\bf x},\xi,\w)$, from the recorded data,
$T$, for correlation.  However, these
autocorrelations will produce the bandlimited impulse functions, 
$\delta({\bf x}-{\bf x}_s,\w)$ in equation~\ref{source}.  Therefore, this 
procedure need not be performed, but instead a zero-phase wavelet
can be modeled that shares the bandwidth of the recorded data.  Therefore,
if sufficient source energy is captured and the resulting transmission
wavefields are correlated and stacked, the appropriate source functions
for migration are impulses at every recording station.  This result 
recapitulates that produced by \longcite{gerry} where it is shown 
that the Kirchoff migration kernel for migrating correlated gathers 
is the same as that used to image prestack active data when one assumes 
that impulsive virtual sources are located at all receiver locations.
\par
Viewing passive seismic processing as an inverse phase-encoding 
problem reinforces the intuitive conception that a transmission 
wavefield is some convolution of the incoming energy from the 
subsurface and the resultant 
reflections from the subsurface after the energy has been mirrored 
at the free surface.  Correlation finds similar energy packets at
all times that correspond to the initial upcoming energy, and  
redefines time zero at each receiver location to correspond to the 
two-way traveltime to subsurface reflectors:  Initial 
waiting time for the arrival of energy at any location is removed
as the correlation calculates time differences between the initial
arrival and subsequent reflections.  Viewed in this light, passive
imaging is an identical operation, with a possibly more complicated
wavelet, to modeling primaries from multiples \cite{DELFT}.  The 
correlation removes the time in a record associated with waiting for 
the arrival of a multiple (or subsurface source), across the entire 
${\bf x}$-axis for every active shot gather (or subsurface source).

\subsection{Chapter description}
Chapter~\ref{chapter:TimeSynth} will analyze the effectiveness of
correlation processing using a two-layer, 2D synthetic model.  Issues
such as source character, distribution, and the inclusion of multiple
source functions within a transmission wavefield will be
illustrated with examples.  
\par
Chapter~\ref{chapter:TimeShot} introduces the Valhall field data and
describes and interprets the shot gathers generated by correlation and
spectral whitening.  Chapter~\ref{chapter:TimeCMP} extends the analysis
of the data by further processing the gathers with NMO and stacking.  
\par
Chapter~\ref{chapter:DepthSyn} Presents migrations of the correlated
gathers shown in Chapter~\ref{chapter:TimeSynth} to explore the
effects of the artifacts within synthesized shot gathers in the image
space.  I then introduce the concept of passive seismic processing as an
imaging problem within the course of shot profile migration: 
direct migration.  The 2D synthetic data is used to illustrate the 
procedure with the suite of modeled data used previously.
\par
Chapter~\ref{chapter:DepthVal} presents migration results from the
Valhall passive data set produced by both correlation/imaging and
direct migration strategies.  An image provided by BP, produced with
data from an active seismic survey, is used as ground truth to
interpret the subsurface content of the final results produced with
the passive data. 
\par
Chapter~\ref{chapter:DepthMult} extends the concept of data processing
at the imaging condition of shot profile migration presented in
Chapter~\ref{chapter:DepthSyn} where the transmission wavefield was
used as both the upcoming and downgoing wavefields.  In this case,
shot gathers are used for both cases, and the resulting image is a 
prediction of the surface-related multiples in the depth domain.  The
prediction volume is equivalent to migrating the multiple prediction
produced by surface-related multiple prediction (SRMP) in the data
domain. 

\clearpage

\end{document}
