Commit f806ca9f authored by thomas.forbriger's avatar thomas.forbriger Committed by thomas.forbriger
Browse files

reintegrate merge of branch lisousi

This is a legacy commit from before 2015-03-01.
It may be incomplete as well as inconsistent.
See COPYING.legacy and README.history for details.


SVN Path:     http://gpitrsvn.gpi.uni-karlsruhe.de/repos/TFSoftware/trunk
SVN Revision: 5277
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent 333dec5a
......@@ -5,6 +5,22 @@
# Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
#
# lisousi: Line Source Simulation
#
# ----
# lisousi is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
# ----
#
# REVISIONS and CHANGES
# 17/04/2013 V1.0 Thomas Forbriger
......@@ -110,8 +126,7 @@ lisousi.cc: description_text.h
lisousi.cc: help_text.h
lisousi: \
%: description_text.o help_text.o %.o \
filterresponse.o \
applytaper.o
$(patsubst %.cc,%.o,$(wildcard *.cc))
$(CXX) -o $@ $^ -I$(LOCINCLUDEDIR) -lfourierxx -lfftw3 -lm \
-lsffxx -ldatrwxx -llinearxx -lgsl -lgslcblas \
-ltsxx -ltfxx -lsffxx -lgsexx -ltime++ -laff \
......
......@@ -5,14 +5,129 @@ lisousi: Line Source Simulation
$Id$
============================================================================
Purpose of lisousi
------------------
Line-source simulation: Transform field data to apparent line-source generated
waveforms in preparation for Cartesian 2D full-waveform inversion.
Shallow seismic field data is excited by point sources (e.g. hammer blows).
Full waveform inversion (FWI) approaches which make use of Cartesian 2D
forward modeling implicitly use a line source to fit the observed data.
Therefore recorded waveforms must be transformed to simulate equivalent
line-source generated data prior to application of 2D FWI.
lisousi offers several single-trace approaches for vertical component and
radial component data excited by vertical hammer blows or explosions.
Single-trace approaches can be applied to each seismic trace individually in
contrast to integral transform approaches, which require data from a complete
profile. Single-trace approaches have the benefit, that they are independent
of geophone layout and that they perform reasonably well for data recorded on
laterally heterogeneous structures. On the downside they must estimate the
wave number from sample time and receiver offset. They inherently are
approximations, most of them being derived from the acoustic wave Green's
function. Nevertheless, they perform surprisingly well on data from
visco-elastic wave propagation in heterogeneous structures.
For shallow seismic data we recommend a simple but effective procedure:
1. scale waveform by r*sqrt(2) (offset times square root of 2)
2. convolve with 1/sqrt(t) (fractional half integration)
3. taper samples with 1/sqrt(t)
Where r is source-to-receiver offset and t is sample time.
Features
--------
- Fourier domain or time domain application of convolution filter (selectable)
- Selectable scaling optimized for either direct waves, non-dispersive waves,
or reflected waves
- Selectable delay of time-domain taper to match source wavelet centroid
delay
- Hybrid transformation for cases where near-offset traces suffer from
delayed taper
- Ability to read a variety of input formats including SeismicUn*x and raw
ASCII (through libdatrwxx )
References
----------
Forbriger, T., Groos, L. and Schäfer, M., 2013. Appropriate line source
simulation procedure for shallow seismic field data. 73rd Annual Meeting of
the German Geophysical Society (DGG), Leipzig.
(http://www.opentoast.de/Data_analysis_code_246.php)
Forbriger, T., Groos, L. and Schäfer, M., 2013. Line source-simulation for
shallow-seismic data. Part 1: Theoretical background. Geophys. J. Int.,
submitted
Schäfer, M., Groos, L., Forbriger, T. & Bohlen, T., 2013. Line-source
simulation for shallow-seismic data. Part 2: Waveform inversion – a 2D case
study, Geophys. J. Int., submitted.
Program documentation
---------------------
See files
usage_text.txt
help_text.txt
description_text.txt
and
for usage information or execute
description_text.txt
lisousi -help
============================================================================
Installation
------------
For compilation instructions see README.1st in the root directory of the
tar-ball or
http://gpitrsvn.gpi.uni-karlsruhe.de:8000/TFSoftware/wiki/docs/installation
The Makefile supports program compilation and linking. Some binary libraries
are required.
The command
make all
will compile and install the binary executables as well as the doxygen
documentation.
Environment variables control where the results are stored and where
libraries and library header files are expected:
LOCLIBDIR defines location of binary libraries
LOCINCLUDEDIR defines location of C/C++ header files (prototypes)
LOCBINDIR defines location of binary executables
TF_WWWBASEDIR defines location of doxygen output
Dependencies:
Compilers required to build the programs:
C++ compiler
C/C++ preprocessor
doxygen (required to process source code documentation)
TFSoftware libraries required to compile the code:
libaff
libtime
libseife
libfourier
liblinearxx
libgsexx
libdatrwxx
libsffxx
libtsxx
libtfxx
============================================================================
The home of this software suite is
http://gpitrsvn.gpi.uni-karlsruhe.de:8000/TFSoftware/wiki/trunk/src/ts/lisousi
for online documentation.
Please send bug reports and suggestions to
Thomas.Forbriger@kit.edu
----- END OF README -----
......@@ -12,12 +12,12 @@
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* This program is free software; you can redistribute it and/or modify
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
......@@ -39,6 +39,7 @@
"$Id$"
#include "lisousi.h"
#include "functions.h"
Tseries applytaper(const Tseries::Tcoc& input, const Tseries::Tcoc& taper,
const double& factor, const double& dt,
......
......@@ -37,7 +37,12 @@ the 1/sqrt(t) filter and the way this filter is applied. Application
of the filter (convolution with 1/sqrt(t)) and the application of
the scaling (tapering with 1/sqrt(t)) are not commutative in a
mathematical sense and the results may change for a different
order of the operations.
order of the operations, although no residuals where observed in
practical application yet. Some of the options are not really of
practical relevance due to more efficient and equally accurate
alternatives. For instance a time domain convolution is not
recommended. These options are primarily provided for experimental
purposes.
Available approaches:
---------------------
......@@ -105,7 +110,7 @@ Default:
-tlim f
-tfac f
They have a lesser (or even unnoticeable) impact when compared
to their impact whne constructing the filter resonse. This is
to their impact when constructing the filter resonse. This is
because the time series to be tapered usually are small or
vanish close to t=0.
......
......@@ -34,14 +34,23 @@
\author Thomas Forbriger
\since July 2012
\date April 2013
\date October 2013
\version V1.0
$Id$
<HR>
\section main_sec_toc Contents of this page
- \ref main_sec_readme
- \ref main_sec_invocation
- \ref main_subsec_help
- \ref main_subsec_description
<HR>
\section main_sec_readme Purpose
\include README
\section main_sec_invocation Invocation
\subsection main_subsec_help Program Options and Parameters
......
/*! \file fcsingleforce.cc
* \brief Fourier coefficients for a single force full space wavefield (implementation)
*
* ----------------------------------------------------------------------------
*
* $Id$
* \author Thomas Forbriger
* \date 17/04/2013
*
* Fourier coefficients for a single force full space wavefield (implementation)
*
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
* ----
*
*
* REVISIONS and CHANGES
* - 17/04/2013 V1.0 Thomas Forbriger
*
* ============================================================================
*/
#define TF_FCSINGLEFORCE_CC_VERSION \
"TF_FCSINGLEFORCE_CC V1.0 "
#define TF_FCSINGLEFORCE_CC_CVSID \
"$Id$"
#include "lisousi.h"
#include "functions.h"
/*! Fourier coefficients of the wave field of a vertical single point source
* in homogeneous full space
*/
TFourier::Tspectrum zfpointfc(const int& n, const double& dt,
const double& offset,
const double& vs,
const double& vp,
const bool& debug)
{
TFourier::Tseries gPN(0,n-1);
const double tp=offset/vp;
const double ts=offset/vs;
TFXX_debug(debug, "zfpointfc", "go..."
<< " " << TFXX_value(n)
<< " " << TFXX_value(vp)
<< " " << TFXX_value(vs)
<< " " << TFXX_value(offset)
<< " " << TFXX_value(dt)
<< " " << TFXX_value(tp)
<< " " << TFXX_value(ts)
);
// prepare near-field response
for (int i=0; i<n; ++i)
{
double t=i*dt;
if ((tp <= t) && (t >= ts))
{
gPN(i)=-t/(offset*offset*offset);
}
else
{
gPN(i)=0.;
}
}
TFXX_debug(debug, "zfpointfc", "Fourier Transformation");
TFourier::Tspectrum FCgPN=Fourier(gPN, dt);
TFourier::Tspectrum retval(FCgPN.shape());
TFXX_debug(debug, "zfpointfc", "total field "
<< TFXX_value(retval.l()) << " "
<< TFXX_value(FCgPN.l()));
const double df=1./(n*dt);
// prepare total response
for (int i=FCgPN.f(); i<=FCgPN.l(); ++i)
{
double f=(i-FCgPN.f())*df;
// far-field
retval(i) = exp(-IME*2.*M_PI*f*ts)/(offset*vs*vs);
TFXX_debug(debug, "zfpointfc",
TFXX_value(i) << " " <<
TFXX_value(f) << " " <<
TFXX_value(abs(retval(i))/abs(FCgPN(i))));
// near-field
retval(i) += FCgPN(i);
}
TFXX_debug(debug, "zfpointfc", "finished");
return(retval);
}
/*----------------------------------------------------------------------*/
/*! Helper function to save sign problems in sqrt
* function must not be called where sqrt argument would be negative
*/
inline
double sqrtfct(const double& psample, const double& v)
{
double retval;
double arg=(v*v*psample*psample)-1.;
TFXX_assert(arg>=0., "negative argument of sqrt: "
"programming error - check source code!");
retval=sqrt(arg)/v;
return (retval);
}
/*----------------------------------------------------------------------*/
/*! Fourier coefficients of the wave field of a vertical single line source
* in homogeneous full space
*/
TFourier::Tspectrum zflinefc(const int& n, const double& dt,
const double& offset,
const double& vs,
const double& vp,
const bool& debug)
{
TFourier::Tseries gLN(0,n-1);
const double tp=offset/vp;
const double ts=offset/vs;
// prepare near-field response
for (int i=0; i<n; ++i)
{
double t=i*dt;
if (t >= tp)
{
gLN(i) = sqrtfct(t/offset,vp);
if (t>=ts)
{
gLN(i) -= sqrtfct(t/offset,vs);
}
gLN(i) *= (-2./offset);
}
else
{
gLN(i)=0.;
}
}
TFourier::Tspectrum FCgLN=Fourier(gLN, dt);
TFourier::Tspectrum retval(FCgLN.shape());
const double df=1./(n*dt);
// prepare total response
for (int i=FCgLN.f(); i<=FCgLN.l(); ++i)
{
double f=(i-FCgLN.f())*df;
// far-field
retval(i) = -IME*M_PI*hankel(ts*2.*M_PI*f)/(vs*vs);
TFXX_debug(debug, "zflinefc",
TFXX_value(i) << " " <<
TFXX_value(f) << " " <<
TFXX_value(abs(retval(i))/abs(FCgLN(i))));
// near-field
retval(i) += FCgLN(i);
}
return(retval);
}
/* ----- END OF fcsingleforce.cc ----- */
......@@ -12,12 +12,12 @@
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* This program is free software; you can redistribute it and/or modify
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
......@@ -39,6 +39,7 @@
"$Id$"
#include "lisousi.h"
#include "functions.h"
/* create a 1/sqrt(t) filter response */
Tseries filterresponse(const int& size, const double& dt, const Options& opt)
......
/*! \file functions.h
* \brief lisousi functions (prototypes)
*
* ----------------------------------------------------------------------------
*
* $Id$
* \author Thomas Forbriger
* \date 19/04/2013
*
* lisousi functions (prototypes)
*
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
* ----
*
*
* REVISIONS and CHANGES
* - 19/04/2013 V1.0 Thomas Forbriger
*
* ============================================================================
*/
// include guard
#ifndef TF_FUNCTIONS_H_VERSION
#define TF_FUNCTIONS_H_VERSION \
"TF_FUNCTIONS_H V1.0 "
#define TF_FUNCTIONS_H_CVSID \
"$Id$"
#include "lisousi.h"
#include "wnintegration.h"
/*======================================================================*/
/* functions */
Tseries applytaper(const Tseries::Tcoc& input, const Tseries::Tcoc& taper,
const double& factor, const double& dt,
const double& offset, const Options& opt);
Tseries filterresponse(const int& size, const double& dt, const Options& opt);
TFourier::Tseries padseries(const TFourier::Tseries& series,
const Options& opt);
TFourier::Tseries singlevelocitytransformation(const TFourier::Tseries& series,
const Parameters& par,
const Exco& ec,
const Options& opt);
TFourier::Tseries tdtapertransformation(TFourier::Tseries series,
Parameters& par,
const Options& opt);
TFourier::Tseries transitionmixer(TFourier::Tseries singlevelocityseries,
TFourier::Tseries directwaveseries,
const Parameters& par,
const Options& opt);
TFourier::Tspectrum zfpointfc(const int& n, const double& dt,
const double& offset,
const double& vs,
const double& vp,
const bool& debug=false);
TFourier::Tspectrum zflinefc(const int& n, const double& dt,
const double& offset,
const double& vs,
const double& vp,
const bool& debug=false);
TFourier::Tcoeff hankel(const double& arg);
#endif // TF_FUNCTIONS_H_VERSION (includeguard)
/* ----- END OF functions.h ----- */
/*! \file globaldata.cc
* \brief a programming unit to create instances of global data (implementation)
*
* ----------------------------------------------------------------------------
*
* $Id$
* \author Thomas Forbriger
* \date 17/04/2013
*
* a programming unit to create instances of global data (implementation)
*
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
* ----
*
*
* REVISIONS and CHANGES
* - 17/04/2013 V1.0 Thomas Forbriger
*
* ============================================================================
*/
#define TF_GLOBALDATA_CC_VERSION \
"TF_GLOBALDATA_CC V1.0 "
#define TF_GLOBALDATA_CC_CVSID \
"$Id$"
#include "lisousi.h"
// global data
// other programming units will make use of this processor
// my Fourier transformation processor
TFourier Fourier;
// provide series for tapering and filtering
Tseries filter, taper, workseries;
/* ----- END OF globaldata.cc ----- */
/*! \file hankel.cc
* \brief Hankel function (implementation)
*
* ----------------------------------------------------------------------------
*
* $Id$
* \author Thomas Forbriger
* \date 17/04/2013
*
* Hankel function (implementation)
*
* Copyright (c) 2013 by Thomas Forbriger (BFO Schiltach)
*
* ----
* lisousi is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* lisousi is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
* ----
*
*
* REVISIONS and CHANGES
* - 17/04/2013 V1.0 Thomas Forbriger
*
* ============================================================================
*/
#define TF_HANKEL_CC_VERSION \
"TF_HANKEL_CC V1.0 "
#define TF_HANKEL_CC_CVSID \
"$Id$"
#include "lisousi.h"
#include "functions.h"
TFourier::Tcoeff hankel(const double& arg)
{