[cig-commits] r13344 - in seismo/3D/ADJOINT_TOMO/flexwin/latex: . figures/fig

alessia at geodynamics.org alessia at geodynamics.org
Wed Nov 19 07:22:24 PST 2008

Author: alessia
Date: 2008-11-19 07:22:24 -0800 (Wed, 19 Nov 2008)
New Revision: 13344

Modified:
Log:
Modifications to abstract, introduction, method.  Still need to look at results, discussion and appendix.

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/abstract.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/abstract.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -1,5 +1,5 @@
\begin{abstract}
%We present an algorithm for the automated selection of time windows on pairs of observed and synthetic seismograms.  {\bf The algorithm was designed specifically to automate window selection and measurement for adjoint tomography studies, but is sufficiently flexible to be adapted to many tomographic applications and seismological scenarios.}  Adjoint tomography utilizes 3D wavefield simulations that capture complex phases that do not necessarily exist in 1D simulations or traditional travel-time curves. {\bf It requires a data selection method that includes these new phases, maximizes the number of measurements made on each seismic record, while avoiding seismic noise. This selection method must also be automated in order to adapt to changes in the synthetic seismograms after each iteration of the tomographic inversion.  These considerations have led us to favor a signal processing approach to the time-window selection problem, and to the development of the FLEXWIN algorithm we present here.} We illustrate the algorithm using datasets from three distinct regions: the entire globe, the Japan subduction zone, and southern California.
-{\bf CHT modified version:}
-We present FLEXWIN, an algorithm for the automated selection of time windows on pairs of observed and synthetic seismograms.  The algorithm was designed specifically to accommodate synthetic seismograms produced from 3D wavefield simulations that capture complex phases that do not necessarily exist in 1D simulations or traditional travel-time curves. Relying on signal processing tools and several user-tuned parameters, the algorithm is able to include these new phases and to maximize the number of measurements made on each seismic record, while avoiding seismic noise.  Our motivation is to use the algorithm for an iterative tomographic inversion, whereby the synthetic seismograms will change from one iteration to the next. Hence, automation is needed to handle the sheer volume of measurements and to allow for the increasing number of windows for each model iteration. The algorithm is sufficiently flexible to be adapted to many tomographic applications and seismological scenarios, including those based on synthetics generated from 1D models.  We illustrate the algorithm using datasets from three distinct regions: the entire globe, the Japan subduction zone, and southern California.
+%{\bf CHT modified version:}
+We present FLEXWIN, an algorithm for the automated selection of time windows on pairs of observed and synthetic seismograms.  The algorithm was designed specifically to accommodate synthetic seismograms produced from 3D wavefield simulations, which capture complex phases that do not necessarily exist in 1D simulations or traditional travel-time curves. Relying on signal processing tools and several user-tuned parameters, the algorithm is able to include these new phases and to maximize the number of measurements made on each seismic record, while avoiding seismic noise.  Our motivation is to use the algorithm for an iterative tomographic inversion, in which the synthetic seismograms change from one iteration to the next. Hence, automation is needed to handle the volume of measurements and to allow for an increasing number of windows for each model iteration. The algorithm is sufficiently flexible to be adapted to many tomographic applications and seismological scenarios, including those based on synthetics generated from 1D models.  We illustrate the algorithm using datasets from three distinct regions: the entire globe, the Japan subduction zone, and southern California.
\end{abstract}

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/acknowledgements.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/acknowledgements.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -7,4 +7,4 @@
Additional global scale data were provided by the GEOSCOPE network.
We thank the Hi-net Data Center (NIED), especially Takuto Maeda and Kazushige Obara, for their help in providing the seismograms used in the Japan examples.
For the southern California examples, we used seismograms from the Southern California Seismic Network, operated by California Institute of Technology and U.S.G.S.
-{\bf The FLEXWIN code makes use of filtering and enveloping algorithms that are part of SAC (Seismic Analysis Code, Lawerence Livermore National Laboratory) provided for free to IRIS members.  We thank Brian Savage for adding interfaces to these algorithms in recent SAC distributions. }
+The FLEXWIN code makes use of filtering and enveloping algorithms that are part of SAC (Seismic Analysis Code, Lawerence Livermore National Laboratory) provided for free to IRIS members.  We thank Brian Savage for adding interfaces to these algorithms in recent SAC distributions.

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/appendix.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/appendix.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -1,5 +1,4 @@
\appendix
-{\bf
\section{Tuning considerations\label{ap:tuning}}
FLEXWIN is not a black-box application, and as such cannot be blindly applied
to any given dataset or tomographic scenario.  The data windowing required by
@@ -70,7 +69,9 @@

\subsection{Examples of user functions\label{ap:user_fn}}

-As concrete examples of how the time dependence of these tuning parameters can be used, we present here the functional forms of these time dependencies used for the three example tomographic scenarios described in the text (Windowing Examples, section~\ref{sec:results}).  CHT modified: In each example we use information (predicted arrival times) derived from 1D Earth models to help guide certain user functions in the windowing algorithm. Note, however, that the actual selection of individual windows is based primarily on the details of the waveforms, and not on information from 1D Earth models.
+As concrete examples of how the time dependence of these tuning parameters can be used, we present here the functional forms of these time dependencies used for the three example tomographic scenarios described in the text (Windowing Examples, section~\ref{sec:results}).
+%CHT modified:
+In each example we use information (predicted arrival times) derived from 1D Earth models to help guide certain user functions in the windowing algorithm. Note, however, that the actual selection of individual windows is based primarily on the details of the waveforms, and not on information from 1D Earth models.

\subsubsection{Global scenario\label{ap:user_global}}

@@ -157,7 +158,7 @@

In the following, $t_P$ and $t_S$ denote the start of the time windows for the crustal P wave and the crustal S wave, computed from a 1D layered model appropriate to Southern California \citep{Wald95}.  The start and end times for the surface-wave time window, $t_{R0}$ and $t_{R1}$, as well as the criteria for the time shifts $\Delta\tau_0(t)$, are derived from formulas in \cite{KomatitschEtal2004}. The source-receiver distance (in km) is denoted by $\Delta$.

-CHT modified
+%CHT modified

For the \trange{6}{40} and \trange{3}{40} data, we use constant values of $r_0(t)=r_0$, $\mathrm{CC}_0(t)=\mathrm{CC}_0$, $\Delta\tau_0(t)=\Delta\tau_0$, and $\Delta \ln A_0(t)=\Delta \ln A_0$. We exclude any arrivals before the $P$ wave and after the Rayleigh wave. This is achieved by the box-car function for $w_E(t)$:
%
@@ -201,5 +202,5 @@
%\Delta\tau_0(t) & = \Delta\tau_0.
\end{align}

-}
+
%-----------------------

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/discussion.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/discussion.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -1,4 +1,3 @@
-{\bf
\section{Using FLEXWIN for tomography}
\label{sec:discuss}

@@ -8,20 +7,20 @@
For this class of natural separation tomographic problems, the advantages of using FLEXWIN over manual or specifically designed automated windowing would be the encapsulation of the selection criteria entirely within the parameters of Table~\ref{tb:params} (and their time-dependent modulation), leading to greater clarity and portability between studies using different inversion methods.

%FLEXWIN is not indicated for tomographic problems in which the extraction and separation of information from overlapping portions of a single timeseries is required, for example studies of higher mode surface wave dispersion for which specific methods -- mode branch stripping \citep{vanHeijstWoodhouse1997}, separation of secondary observables \citep{CaraLeveque1987, Debayle1999}, partitioned waveform and automated multimode inversion \citep{Nolet1990, LebedevEtal2005}, non-linear direct search \citep{YoshizawaKennett2002b, VisserEtal2007}  -- have been developed.
-{\bf CHT modify}
+%{\bf CHT modify}
FLEXWIN is not intended for tomographic problems in which the extraction and separation of information from overlapping portions of a single timeseries is required, for example studies of higher mode surface wave dispersion for which specific methods have been developed, for example, mode branch stripping \citep{vanHeijstWoodhouse1997}, separation of secondary observables \citep{CaraLeveque1987, Debayle1999}, partitioned waveform and automated multimode inversion \citep{Nolet1990, LebedevEtal2005}, and non-linear direct search \citep{YoshizawaKennett2002b, VisserEtal2007}.

The full power of FLEXWIN can only be unleashed for problems -- such as adjoint tomography -- which do not require the separation (natural or otherwise) of seismic phases.  The specificity of adjoint tomography, among the 3D-3D tomographic methods, is to calculate the sensitivity kernels by interaction between the wavefield used to generate the synthetic seismograms and an adjoint wavefield whose source term is derived from measurements of misfit between the synthetic and observed seismograms \cite{TrompEtal2005, LiuTromp2006}.  The manner in which the adjoint sources are constructed is specific to each type of measurement (e.g. waveform difference, cross-correlation time-lag, multi-taper phase and amplitude anomaly), but once formulated can be applied indifferently to any part of the seismogram.  Adjoint methods have been used to calculate kernels of various body- and surface-wave phases with respect to isotropic elastic parameters and interface depths \citep{LiuTromp2006}, and with respect to anisotropic elastic parameters \citep{SieminskiEtal2007a,SieminskiEtal2007b}.  Adjoint methods allow us to calculate kernels for each and every wiggle on a given seismic record, thereby giving access to virtually all the information contained within.
-}

+
It is becoming clear, as more finite-frequency tomography models are published, that better kernels on their own are not the answer to the problem of improving the resolution of tomographic studies.  \cite{TrampertSpetzler2006} and \cite{BoschiEtal2007} investigate the factors limiting the quality of finite-frequency tomography images, and conclude that incomplete and inhomogeneous data coverage limit in practice the improvement in resolution that accurate finite-frequency kernels can provide.  The current frustration with the data-induced limitations to the improvements in wave-propagation theory is well summarized by \cite{Romanowicz2008}.  The ability of adjoint methods to deal with all parts of the seismogram indifferently means we can incorporate more information from each seismogram into a tomographic problem, thereby improving data coverage.

The computational cost of constructing an adjoint kernel is independent of the number of time windows on each seismogram we choose to measure, and also of the number of records of a given event we choose to work with.  It is therefore computationally advantageous to make measurements on as many records as possible for each event, while covering as much as possible of each record.  There are, however, certain limits we must be aware of.  As mentioned in the introduction, there is nothing in the adjoint method itself that prevents us from constructing a kernel from noise-dominated portions of the data.  As the purpose of 3D-3D tomography is to improve the fine details of Earth models, it would be counterproductive to pollute the inversion process with such kernels.
It is clear that the use of adjoint methods for tomography requires a strategy for selecting and windowing seismograms that avoids seismic noise while at the same time extracting as much information as possible from the signals.

-{\bf The adjoint kernels are only strictly valid for the 3D Earth model they were constructed in, and therefore need to be re-computed at each iteration of the tomographic inversion \citep{TapeEtal2007}.  At each iteration, the similarities between the synthetic and observed seismograms improve, such that for later iterations a greater proportion of the waveform is adequate for measurement.  In order to take advantage of this extra information, the windowing method used to isolate the portions of the waveform to be measured needs to be automated.}
+The adjoint kernels are only strictly valid for the 3D Earth model they were constructed in, and therefore need to be re-computed at each iteration of the tomographic inversion \citep{TapeEtal2007}.  At each iteration, the similarities between the synthetic and observed seismograms improve, such that for later iterations a greater proportion of the waveform is adequate for measurement.  In order to take advantage of this extra information, the windowing method used to isolate the portions of the waveform to be measured needs to be automated.
The method must also be adaptable to the features that exist in the seismograms themselves, because 3D wavefield simulations are able to synthesize phases that do not exist in 1D simulations or traditional travel-time curves.  All these considerations led us to favor a signal processing approach to the problem of data selection, an approach which in turn led to the development of the FLEXWIN algorithm we have presented here.

Finally, we note that the design of this algorithm is based on the desire {\em not} to use the entire time series of each event when making a measurement between data and synthetics. If one were to simply take the waveform difference between two time series, then there would be no need for selecting time windows of interest. However, this ideal approach \citep[e.g.,][]{GauthierEtal1986} may only work in real applications if the

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/figures/fig/window_rejection_global_data.fig	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/figures/fig/window_rejection_global_data.fig	2008-11-19 15:22:24 UTC (rev 13344)
@@ -78,5 +78,5 @@
2 5 0 1 0 -1 60 -1 20 0.000 0 0 -1 0 0 5
0 ../050295B.050-150/ABKT.II.LHZ.win_reject.eps
900 -190 10282 -190 10282 11700 900 11700 900 -190
-4 0 0 50 -1 16 20 0.0000 4 135 90 360 4950 B\001
-4 0 0 50 -1 16 20 0.0000 4 135 90 360 9225 C\001
+4 0 0 50 -1 16 20 0.0000 4 135 90 360 4950 C\001
+4 0 0 50 -1 16 20 0.0000 4 135 90 360 9225 D\001

===================================================================
(Binary files differ)

===================================================================
(Binary files differ)

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/introduction.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/introduction.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -10,7 +10,7 @@
in 1D media \citep[e.g.][]{LiTanimoto1993,DahlenBaig2002,DahlenZhou2006} and numeric sensitivity kernels in 3D media \citep[e.g.][]{Capdeville2005,TrompEtal2005,ZhaoEtal2005}.
The analytic kernels have been taken up rapidly by tomographers, and used to
produce new 3D Earth models \citep[e.g.][]{MontelliEtal2004a,ZhouEtal2006}.  The numeric kernels have
-opened up the possibility of 3D-3D' tomography, i.e.~seismic tomography based upon a 3D reference model, 3D numerical simulations of the seismic wavefield and finite-frequency sensitivity kernels \citep{TrompEtal2005,ChenEtal2007b}.
+opened up the possibility of 3D-3D' tomography, i.e.~seismic tomography based upon a 3D reference model, 3D numerical simulations of the seismic wavefield, and finite-frequency sensitivity kernels \citep{TrompEtal2005,ChenEtal2007b}.

It is common practice in tomography to work only with certain subsets of the
available seismic data. The choices made in selecting these subsets are
@@ -27,17 +27,17 @@
thereby reducing
the data restrictions required when using approximate forward modelling and simplified descriptions
of sensitivity.  These methods seem to be the best candidates for studying
-regions with complex 3D structure as they permit the use of a larger proportion of the information contained within each seismogram, including complex arrivals not predicted by 1D approximations of Earth structure.  In order to exploit the full power of 3D-3D tomographic methods, we require a new data selection strategy that does not exclude such complex arrivals.
+regions with complex 3D structure, as they permit the use of a larger proportion of the information contained within each seismogram, including complex arrivals not predicted by 1D approximations of Earth structure.  In order to exploit the full power of 3D-3D tomographic methods, we require a new data selection strategy that does not exclude such complex arrivals.

-{\bf As data selection strategies for tomography depend so closely on the tomographic technique, there are nearly as many such strategies as there are tomographic methods.  Furthermore, many of these strategies have been automated in some way, as larger and larger volumes of data have become available for tomographic studies.}
+As data selection strategies for tomography depend so closely on the tomographic technique, there are nearly as many such strategies as there are tomographic methods.  Furthermore, many of these strategies have been automated in some way, as larger and larger volumes of data have become available.
% CHT cut following sentence
%We shall not try to be exhaustive in the following short presentation of currently available data selection methods, but shall limit our discussion to a few representative examples.
Body wave studies that have moved away from using manual travel-time picks or catalog arrival times generally pick windows around specific arrivals defined by predicted arrival times, and include automated tests on arrival time separation and/or the fit of observed to synthetic pulses to reject inadequate data \cite[e.g.][]{RitsemaVanHeijst2002,LawrenceShearer2008}.  Partial automation of the \cite{vanDecarCrosson1990} multi-channel cross-correlation method has led to efficient methods for obtaining highly accurate travel-time \citep{SiglochNolet2006,HouserEtal2008} and even attenuation \citep{LawrenceEtal2006} measurements.  In the surface wave community, there has been much work done to automate methods for extracting dispersion characteristics of fundamental mode \citep{TrampertWoodhouse1995, LaskeMasters1996, EkstromEtal1997, LevshinRitzwoller2001} and higher mode \citep{vanHeijstWoodhouse1997, Debayle1999, YoshizawaKennett2002b, BeuclerEtal2003, LebedevEtal2005, VisserEtal2007} surface waves.  Recently, \cite{PanningRomanowicz2006} have described an algorithm to semi-automatically pick body and surface wavepackets based on the predicted traveltimes of several phases.

%{\bf The automated data selection method we present in this paper goes further than that of \cite{PanningRomanowicz2006}, as it is not tied to arrival time predictions of known phases.}
-{\bf CHT modified version:}
-Our algorithm is designed for tomographic applications with 3D earth reference models.
-Unlike the techniques discussed above, ours is not tied to arrival time predictions of known phases, and it therefore is able to accommodate complex phases due to 3D complexity in structure.
+%{\bf CHT modified version:}
+Our algorithm is designed for tomographic applications with 3D Earth reference models.
+Unlike the techniques discussed above, ours is not tied to arrival time predictions of known phases, and, therefore, is able to accommodate complex phases due to 3D structure.
One promising approach to 3D-3D tomography is based upon adjoint methods \citep{Tarantola1984,TrompEtal2005,LiuTromp2006,TapeEtal2007}.  In adjoint tomography'' the sensitivity kernels that tie variations
in Earth model parameters to variations in the misfit are obtained by
interaction between the wavefield used to generate the synthetic seismograms (the
@@ -47,13 +47,13 @@
The adjoint kernel calculation procedure allows us to measure and use for
tomographic inversion almost any part of the seismic signal.  We do not
need to identify specific seismic phases, as the kernel will take care of
-defining the relevant sensitivities.  However, there is nothing in the adjoint method itself that prevents us from constructing an adjoint kernel from noise-dominated data, and thereby polluting our inversion process.
-{\bf An appropriate data selection strategy for adjoint tomography should therefore define measurement time windows that cover as much of a given seismogram as possible, whilst avoiding portions of the waveform that are dominated by noise.}
+defining the relevant sensitivities.  However, there is nothing in the adjoint method itself that prevents us from constructing an adjoint kernel from noise-dominated data, and thereby polluting our inversion.
+An appropriate data selection strategy for adjoint tomography should therefore define measurement time windows that cover as much of a given seismogram as possible, whilst avoiding portions of the waveform that are dominated by noise.

From a signal processing point of view, the simplest way to avoid serious
contamination by noise is to select and measure strong signals, which in
seismology correspond to seismic arrivals.
-{\bf Our strategy is therefore to select time windows on the synthetic seismogram within which the waveform contains a distinct energy arrival, then require to an adequate correspondence between observed and synthetic waveforms within these windows.}
+Our strategy is therefore to select time windows on the synthetic seismogram within which the waveform contains a distinct energy arrival, then require an adequate correspondence between observed and synthetic waveforms within these windows.
This selection paradigm is general, and can be applied to synthetic seismograms regardless of how they have been obtained.
It is clear, however, that a synthetic seismogram obtained by 3D propagation through a good 3D Earth model will provide a better fit to the observed seismogram over a greater proportion of its length than will be the case for a more approximate synthetic seismogram.

@@ -61,12 +61,12 @@
associated with distinct energy arrivals, we need to analyse the character of the synthetic waveform itself.  This analysis is similar to that used on observed waveforms
in automated phase detection algorithms for the routine location of
earthquakes.
-{\bf In designing our time-window selection algorithm, we have taken a tool used in this detection process -- the long-term / short-term average ratio -- and applied it to the definition of time windows around distinct seismic phases.}
+In designing our time-window selection algorithm, we have taken a tool used in this detection process -- the long-term / short-term average ratio -- and applied it to the definition of time windows around distinct seismic phases.

The choices made in time-window selection for tomography are
interconnected with all aspects of the tomographic inversion process,
from the waveform simulation method (direct problem), through the choice of
measurement method, to the method
-used to obtain sensitivity kernels.  One of the major difficulties in defining a general data selection strategy is the great range of possible choices open to the tomographer.  We have  designed a configurable data selection process that can be adapted to different tomographic scenarios by tuning a handful of parameters (see Table~\ref{tb:params}).  Although we have designed our algorithm for use in adjoint tomography, its inherent flexibility should make it useful in many data-selection applications.
+used to obtain sensitivity kernels, and that used for the inversion itself.  One of the major difficulties in defining a general data selection strategy is the great range of possible choices open to the tomographer.  We have  designed a configurable data selection process that can be adapted to different tomographic scenarios by tuning a handful of parameters (see Table~\ref{tb:params}).  Although we have designed our algorithm for use in adjoint tomography, its inherent flexibility should make it useful in many data-selection applications.

We have successfully applied our windowing algorithm, the details of which are described in Section~\ref{sec:algorithm}, to diverse seismological scenarios: local and near regional tomography in Southern California, regional subduction-zone tomography in Japan, and global tomography.  We present examples from each of these scenarios in Section~\ref{sec:results}, and we discuss the use of the algorithm in the context of tomography in Section~\ref{sec:discuss}.

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/method.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/method.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -5,7 +5,7 @@
on the type of simulation used to generate the synthetics, though realistic
Earth models and more complete propagation theories yield waveforms that are more similar to the observed
seismograms, and thereby allow the definition of measurement windows
-covering more of the available data.  The input seismograms can be measures of displacement, velocity or acceleration, indifferently.  There is no requirement for horizontal signals to be rotated into radial and transverse directions.  All the synthetic seismograms presented in
+covering more of the available data.  The input seismograms can be measures of displacement, velocity, or acceleration, indifferently.  There is no requirement for horizontal signals to be rotated into radial and transverse directions.  All the synthetic seismograms presented in
this paper have been generated using the SPECFEM3D package \citep{KomatitschEtal2002}.

The window selection process has five stages, each of which is discussed individually
@@ -16,7 +16,7 @@
seismograms; {\em \stge:} resolution of preliminary window overlaps.  The parameters that permit tuning of the
window selection towards a specific tomographic scenario are all contained in a
simple parameter file (see Table~\ref{tb:params}).  More complexity and finer
-tuning can be obtained by rendering some of these parameters time dependent via user defined functions that can depend on the source parameters (e.g. event location or depth).
+tuning can be obtained by making some of these parameters time dependent via user defined functions that can depend on the source parameters (e.g. event location or depth).

%----------------------

@@ -28,7 +28,7 @@
%----------------------

We apply minimal and identical pre-processing to both observed and synthetic
-seismograms: {\bf removal of any linear trend, tapering, and} band-pass filtering with a non-causal Butterworth
+seismograms: removal of any linear trend, tapering, and band-pass filtering with a non-causal Butterworth
filter, whose
short and long period corners we denote by $T_0$ and $T_1$ respectively.
Values of these corner periods should reflect the information content of the data,
@@ -56,20 +56,20 @@

%----------------------

-{\bf Detection and identification of seismic phase arrivals is routinely performed by automated
+Detection and identification of seismic phase arrivals is routinely performed by automated
earthquake location algorithms \citep[e.g.][]{Allen1982, EarleShearer1994, AsterRowe2000, BaiKennett2000, SleemanVanEck2003}. We have taken a tool used in most implementations of the
automated detection process -- the short-term long-term average ratio \citep[STA:LTA, e.g.][]{WithersEtal1998,BaiKennett2001}
--- and adapted it to the task of defining time windows around seismic phases.  }
+-- and adapted it to the task of defining time windows around seismic phases.
Given a synthetic seismogram $s(t)$, we derive an
STA:LTA timeseries using an iterative algorithm applied to the envelope of the synthetic.
If we denote the Hilbert transform of the synthetic seismogram by
-$\mathcal{H}[s(t)]$, its envelope $e(t)$ is given by:
+$\mathcal{H}[s(t)]$, its envelope $e(t)$ is given by
$$e(t) = | s(t) + i \mathcal{H}[s(t)] |.$$
In order to create the STA:LTA waveform $E(t)$, we discretize the envelope time
series with timestep $\delta t$, calculate its short term average
-$S(t_i)$ and its long term average $L(t_i)$ {\bf recursively} as follows,
+$S(t_i)$ and its long term average $L(t_i)$ recursively
\begin{align}
S(t_i) & = C_S \; S(t_{i-1}) + e(t_i) , \label{eq:sta}\\
L(t_i) & = C_L \; L(t_{i-1}) + e(t_i) , \label{eq:lta}
@@ -91,15 +91,15 @@
where the use of $T_0$, the low-pass corner period of our band-pass filter,
substitutes that of the dominant period.

-{\bf
+
The STA:LTA of a constant signal converges to a constant value when the length
of the time-series is greater than the effective averaging length of the
long-term average.  An energy
-arrival in $e(t)$ causes the $E(t)$ to rise sharply, then drop to below the
+arrival in $e(t)$ causes $E(t)$ to rise sharply, then drop to below the
convergence value after the arrival, before stabilizing again.  The maximum
-height reached by the $E(t)$ for a given $T_0$ depends on the amplitude of the
+height reached by $E(t)$ for a given $T_0$ depends on the amplitude of the
arrival in $e(t)$ and on its duration: higher amplitudes and shorter durations
-cause higher $E(t)$ values.  The depth to which the $E(t)$ drops after the
+cause higher $E(t)$ values.  The depth to which $E(t)$ drops after the
end of the arrival depends on the same parameters: higher amplitudes and longer
durations cause deeper drops in $E(t)$ that take longer to return to the
convergence value for a constant signal.
@@ -113,8 +113,8 @@
synthetic, and that the local minima in $E(t)$ correspond to the
transitions between one phase and the next.  In the following sections we shall
explain how we use these correspondences to define time windows.
-}

+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\subsection{\stgb\ \label{sec:stageB}}
@@ -129,12 +129,13 @@
phases, we must determine rules for deciding when adjacent local maxima
should be part of a single window.  From an algorithmic point
of view, it is simpler to create all possible combinations of adjacent windows
-and subsequently reject the unacceptable ones, than to consider {\bf combining}
+and subsequently reject the unacceptable ones, than to consider combining
small, single-maximum windows into larger ones.

We start by defining a water level on $E(t)$ via the time dependent parameter
$w_E(t)$ in Table~\ref{tb:params}.
-{\bf All local maxima that lie above $w_E(t)$ are used for the creation of candidate time windows.}
+All local maxima that lie above $w_E(t)$ are considered acceptable,
+and are used for the creation of candidate time windows.
The water level shown in
Figure~\ref{fg:stalta} corresponds to $w_E=0.08$ for the duration of the main
seismic signal.  Once set for typical seismograms for a given
@@ -142,7 +143,7 @@
seismogram.  This is also true of all the other parameters in
Table~\ref{tb:params}: once the system has been tuned,
these parameters remain unchanged and are used for all seismic events in the same scenario.
-{\bf A summary of the main considerations the user should take into account to tune these parameters can be found in Appendix~\ref{ap:tuning}.}
+A summary of the main considerations the user should take into account to tune these parameters can be found in Appendix~\ref{ap:tuning}.
The functional forms of the time-dependent parameters are defined by the user, can depend on
remain unchanged once the system has been tuned.
@@ -167,7 +168,7 @@
\subsection{\stgc\ \label{sec:stageC}}
%{\em Parameters used: $T_0$, $w_E(t)$, $c_{0-4}$.}

-After having created a complete suite of candidate time windows in the manner
+After having created a complete set of candidate time windows in the manner
described above, we start the rejection process.  We reject windows based on
two sets of criteria concerning respectively the shape of the STA:LTA waveform $E(t)$,
and the similarity of the observed and synthetic waveforms
@@ -191,7 +192,8 @@
Thirdly, we reject windows whose seed maximum $E(t_M)$ rises by less than
$c_2 w_E(t_M)$ above either of its adjacent minima.  Subdued local maxima of
this kind represent minor changes in waveform character, and should not be used
-to anchor time windows.  They may, however, be considered as part of a time window with a more prominent maximum (see Figure~\ref{fg:win_composite}c).
+to anchor time windows.  They may, however, be included within a time window with a
+more prominent seed maximum (see Figure~\ref{fg:win_composite}c).
Lastly, we reject windows that contain at least
one strong phase arrival that is well separated in time from $t_M$.  The
rejection is performed using the following criterion:
@@ -213,13 +215,13 @@

If we take
as an example $c_{3a}=1$, this criterion leads to the automatic rejection of
-windows containing a local maximum that is higher than {\bf their} seed maximum; it also leads to the rejection of windows containing a local maximum that is
-lower than {\bf their} seed maximum if {\bf this local maximum} is sufficiently distant in time from
-$t_M$.  {\bf This criterion allows us to distinguish unseparable phase groups from distinct seismic phases.}
+windows containing a local maximum that is higher than their seed maximum; it also leads to the rejection of windows containing a local maximum that is
+lower than their seed maximum if this local maximum is sufficiently distant in time from
+$t_M$.  This criterion allows us to distinguish unseparable phase groups from distinct seismic phases.

The candidate windows that remain after application of these four rejection
-criteria are almost ready to be passed on to the next stage, in which they will
-be evaluated in terms of the similarity between observed and synthetic
+criteria are almost ready to be passed on to the next stage, in which we shall
+evaluate the similarity between observed and synthetic
waveforms within the window limits.  Special precautions may have to be taken,
however, in the case of windows that contain long coda waves: the
details of codas are often poorly matched by synthetic seismogram calculations,
@@ -256,12 +258,11 @@
A$. The signal-to-noise ratio for single windows is defined as an amplitude ratio,${\rm SNR}_W=A_{\rm window}/A_{\rm noise}$, where$A_{\rm window}$and$A_{\rm noise}$are the maximum absolute values of the observed seismogram$|d(t)|$in the window -and in the noise time-span respectively (see equation~\ref{eq:noise}). The cross-correlation value$\mathrm{CC}$is defined as the maximum value of the +and in the noise time-span respectively (the noise time-span is the same as that for equation~\ref{eq:noise}). The cross-correlation value$\mathrm{CC}$is defined as the maximum value of the cross-correlation function$\mathrm{CC}={\rm max}[\Gamma(t)]$, where $$-\Gamma(t) = \int s(t'-t)d(t')\,\rmd t', +\Gamma(t) = \frac {\int s(t'-t)d(t')\,\rmd t'}{ \left[\int s^2(t')\,\rmd t' \int d^2(t')\,\rmd t' \right]^{1/2}}$$ -and quantifies the similarity in shape between the$s(t)$and$d(t)$waveforms. The time lag$\Delta \tau$is defined as the value of$t$at which$\Gamma$is maximal, and quantifies the delay in time between a @@ -277,14 +278,14 @@ = 0.5 \ln \left[ \frac{\int d(t)^2 \,\rmd t}{\int s(t)^2 \,\rmd t} \right] \label{eq:dlnA_def} % -(Note that \citet[][Eq.~3]{DahlenBaig2002} is the first-order approximation of (\ref{eq:dlnA_def}).) +(note that \citet[][Eq.~3]{DahlenBaig2002} is the first-order approximation of equation~\ref{eq:dlnA_def}). The limits that trigger rejection of windows based on the values of these four -quantities are the time dependent parameters$r_0(t)$,$\mathrm{CC}_0(t)$,$\Delta
-\tau_0(t)$and$\Delta \ln A_0(t)$in Table~\ref{tb:params}. -As for the STA:LTA water level$w_E(t)$used in above, the functional form of -these parameters is defined by the user, and can depend on source and receiver +quantities are the parameters$r_0(t)$,$\mathrm{CC}_0(t)$,$\Delta \tau_{\rm ref}$,$\Delta
+\tau_0(t)$,$\Delta \ln A_{\rm ref}$and$\Delta \ln A_0(t)$in Table~\ref{tb:params}. +As for the STA:LTA water level$w_E(t)$used above, the functional forms of +the time dependent parameters are defined by the user, and can depend on source and receiver parameters such as epicentral distance and earthquake depth. -{\bf Examples of functional forms for these parameters can be found in Appendix~\ref{ap:user_fn}.} +Examples of functional forms for these parameters can be found in Appendix~\ref{ap:user_fn}. Figure~\ref{fg:criteria} shows the time dependence of$\mathrm{CC}_0$,$\Delta \tau_0$and$\Delta \ln A_0for the example seismogram of Figure~\ref{fg:win_rej_data}. @@ -296,7 +297,7 @@ %|\Delta\tau| & \leq \Delta\tau_0(t_M), \label{eq:tau} \\ %|\Delta\ln{A}| & \leq \Delta\ln{A}_0(t_M), \label{eq:dlnA} %\end{align} -{\bf CHT modify equation and following text} +%{\bf CHT modify equation and following text} \begin{align} {\rm SNR}_W & \geq r_0(t_M), \label{eq:snr_win} \\ {\rm CC} & \geq {\rm CC}_0(t_M), \label{eq:cc} \\ @@ -307,18 +308,19 @@ where % \begin{eqnarray} -\Delta\tau_{\rm min} &\equiv& \Delta\tau_{\rm ref} - \Delta\tau_0(t_M) \\ -\Delta\tau_{\rm max} &\equiv& \Delta\tau_{\rm ref} + \Delta\tau_0(t_M) \\ -\Delta\ln{A}_{\rm min} &\equiv& \Delta\ln{A}_{\rm ref} - \Delta\ln{A}_0(t_M) \\ -\Delta\ln{A}_{\rm max} &\equiv& \Delta\ln{A}_{\rm ref} + \Delta\ln{A}_0(t_M) +\Delta\tau_{\rm min} &\equiv& \Delta\tau_{\rm ref} - \Delta\tau_0(t_M), \\ +\Delta\tau_{\rm max} &\equiv& \Delta\tau_{\rm ref} + \Delta\tau_0(t_M), \\ +\Delta\ln{A}_{\rm min} &\equiv& \Delta\ln{A}_{\rm ref} - \Delta\ln{A}_0(t_M), \\ +\Delta\ln{A}_{\rm max} &\equiv& \Delta\ln{A}_{\rm ref} + \Delta\ln{A}_0(t_M), \end{eqnarray} % -wheret_M$is the time of the window's seed maximum. The reference measurement variables,$\Delta\tau_{\rm ref}$and$\Delta\ln{A}_{\rm ref}$, allow for systematic differences between data and synthetics. For example, if the 3D velocity model is on average too fast, then$\Delta\tau_{\rm ref}$should be set to a positive value. Or if the magnitudes of the synthetic sources lead to systematically non-zero$\Delta\ln{A}$values, then$\Delta\ln{A}_{\rm ref}$should be chosen accordingly. In essence, these reference values should designate the approximate center-value of the distribution of measurements (for example, Figure~\ref{fg:socal_rs_T06}c). - -In words, we only accept +and$t_M$is the time of the window's seed maximum. In words, we only accept windows in which the observed signal is sufficiently above the noise level, the observed and synthetic signals are reasonably similar in shape, their arrival time -differences are small, and their amplitudes are broadly compatible. When the synthetic and observed +differences are small, and their amplitudes are broadly compatible. +The parameters$\Delta\tau_{\rm ref}$and$\Delta\ln{A}_{\rm ref}$allow the algorithm to work efficiently when there are systematic differences between data and synthetics. For example, if the Earth model is on average too fast, then$\Delta\tau_{\rm ref}$should be set to a positive value. Or if the magnitudes of the synthetic sources lead to systematically non-zero$\Delta\ln{A}$values, then$\Delta\ln{A}_{\rm ref}$should be chosen accordingly. In practice, these reference values should designate the approximate center-value of the distribution of measurements (see, for example, Figure~\ref{fg:socal_rs_T06}c). + +When the synthetic and observed seismograms are similar, the fit-based criteria of equations~(\ref{eq:cc})-(\ref{eq:dlnA}) reject only a few of the candidate data windows (see lower portion of Figure~\ref{fg:win_rej_data}). They are @@ -353,8 +355,8 @@ The condition that optimal windows should have passed all previous tests removes the straightforward solution of merging overlapping windows. Indeed, given any two overlapping windows, we know that the window defined by their merger -existed in the complete list of candidate windows obtained at the end of -\stgb, and that its absence from the current list means it was rejected +existed in the complete set of candidate windows obtained at the end of +\stgb, and that its absence from the current set means it was rejected either because of the shape of its$E(t)$time-series (\stgc), or because of an inadequate similarity between observed and synthetic waveforms (\stgd). It would therefore be meaningless to re-instate such a window at this stage. @@ -367,9 +369,9 @@ seismogram covered by the windows, average cross-correlation value for the windows, and total number of windows. These criteria often work against each other. For example, a long window may have a lower$\mathrm{CC}$than two shorter ones, if the two -short ones have different time lags$\Delta\tau$. {\bf Weighting of the +short ones have different time lags$\Delta\tau$. Weighting of the three scores is necessary, and is controlled by the three parameters -$w_{\mathrm{CC}}$,$w_{\rm len}$and$w_{\rm nwin}$in Table~\ref{tb:params}. } +$w_{\mathrm{CC}}$,$w_{\rm len}$and$w_{\rm nwin}$in Table~\ref{tb:params}. As can be seen in Figure~\ref{fg:stageE}, the generation of subsets is facilitated by first grouping candidate windows such that no group overlaps @@ -391,9 +393,9 @@ S = \frac{w_{\mathrm{CC}}S_{\mathrm{CC}}+w_{\rm len}S_{\rm len}+w_{\rm nwin}S_{\rm nwin}}{w_{\mathrm{CC}}+w_{\rm len}+w_{\rm nwin}}. \label{eq:score} -{\bf The best subset of candidate windows within each group is the one with the +The best subset of candidate windows within each group is the one with the highest combined score$S\$.  The final set of windows is
given by concatenating the best subsets of candidate windows for each group.
Figure~\ref{fg:res_abkt} shows an example of final windows selected on real
-data.}
+data.

===================================================================
--- seismo/3D/ADJOINT_TOMO/flexwin/latex/results.tex	2008-11-19 06:49:04 UTC (rev 13343)
+++ seismo/3D/ADJOINT_TOMO/flexwin/latex/results.tex	2008-11-19 15:22:24 UTC (rev 13344)
@@ -26,7 +26,7 @@
scenario-dependent information is encapsulated in the tuning parameters of
Table~\ref{tb:params}.

-We tuned the windowing algorithm separately for each of the three scenarios we present here, and we present examples based on the events listed in Table~\ref{tb:events}.  Tuning parameter values for each scenario can be found in Table~\ref{tb:example_params}, while the functional forms of the time-dependent parameters can be found in Appendix~\ref{ap:user_fn}.  {\bf Once tuned for a given scenario, the algorithm is applied to all its events without further modification. }
+We tuned the windowing algorithm separately for each of the three scenarios we present here, and we present examples based on the events listed in Table~\ref{tb:events}.  Tuning parameter values for each scenario can be found in Table~\ref{tb:example_params}, while the functional forms of the time-dependent parameters can be found in Appendix~\ref{ap:user_fn}.  Once tuned for a given scenario, the algorithm is applied to all its events without further modification.

\subsection{Global tomography}