Commit 3a2c2e8b authored by Daniel Armbruster's avatar Daniel Armbruster Committed by thomas.forbriger
Browse files

read and write any format

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: 4579
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent 5bd33a75
......@@ -41,6 +41,7 @@
* - provide exponential decay
* - 02/08/2007 V1.5 explain results line
* - 11/02/2011 V1.6 add definition of traces
* - 20/02/2012 V1.7 read and write any file format (damb)
*
* ============================================================================
*/
......@@ -55,9 +56,12 @@
#include <vector>
#include <cmath>
#include <tfxx/commandline.h>
#include <tfxx/xcmdline.h>
#include <tfxx/rangelist.h>
#include <tfxx/rangestring.h>
#include <tfxx/error.h>
#include <sffxx.h>
#include <sffostream.h>
#include <datrwxx/readany.h>
#include <datrwxx/writeany.h>
#include <tsxx/tsxx.h>
#include <tsxx/innerproduct.h>
#include <aff/array.h>
......@@ -65,10 +69,10 @@
#include <aff/seriesoperators.h>
#include <linearxx/lapackxx.h>
typedef std::list<std::string> Tnamelist;
// double precision time series bundle
//typedef std::list<std::string> Tnamelist;
// double precision time series bundle which includes a WID2 header
typedef ts::TDsfftimeseries Tbundle;
typedef Tbundle::Tseries Tseries;
typedef datrw::Tdseries Tseries;
typedef std::vector<Tbundle> Tbundlevec;
typedef aff::Array<Tbundle::Tseries::Tvalue> Tmatrix;
......@@ -88,6 +92,7 @@ struct Options {
bool equalsearch;
bool fittrend, fitoffset, doskip, fitexp;
double amptrend, ampoffset, skip, tcexp;
std::string itype, otype;
}; // struct Options
/*======================================================================*/
......@@ -142,7 +147,8 @@ int main(int iargc, char* argv[])
" [-Sramp[=v]] [-Sconst[=v]] [-Sexp[=v]]" "\n"
" [-residual f]" "\n"
" [-searchrange[=r]] [-equalsearch] [-skip n]" "\n"
" signal trial [trial ...]" "\n"
" [-itype F] [-otype F]" "\n"
" signal [t:T f:F] trial [t:T f:F] [trial [t:T f:F]...]" "\n"
" or: sigfit --help|-h" "\n"
};
......@@ -168,10 +174,18 @@ int main(int iargc, char* argv[])
" rms(signal)/rms(trial signals)" "\n"
"-equalsearch use same search range for each trial signal" "\n"
"-skip n skip n seconds at beginning of each trace" "\n"
"-itype F set input file formats to F (default: sff)" "\n"
"-otype F set output (residual) file format to F (default: sff)" "\n"
"\n"
"signal signal to by explained by a linear combination" "\n"
" of trial signals" "\n"
"trial any set of trial signals"
"trial any set of trial signals" "\n"
"\n"
"File specific options:" "\n"
"t:T select specific traces from input file" "\n"
" T can be a list of traces like '1,4,5' or" "\n"
" a range like '6-19' or mixed like '5,8,12-17,20'" "\n"
"f:F specific file format (overrides -itype setting)" "\n"
"\n"
"The last line of the output is a summary of all results." "\n"
"This line contains the following information separated by blanks:" "\n"
......@@ -226,9 +240,26 @@ int main(int iargc, char* argv[])
{"skip",arg_yes,"0"},
// 10: fit an offset
{"Sexp",arg_opt,"1."},
// 11: input file format
{"itype",arg_yes,"sff"},
// 12: output file format
{"otype",arg_yes,"sff"},
{NULL}
};
//! key to select traces
const char* const tracekey="t";
//! key to select file format
const char* const formatkey="f";
//! list of keys for filename specific parameters
static const char* keys[]={
//! select traces
tracekey,
//! set file format
formatkey,
0
}; // const char* keys[]
// no arguments? print usage...
if (iargc<2)
{
......@@ -242,8 +273,8 @@ int main(int iargc, char* argv[])
// help requested? print full help text...
if (cmdline.optset(0))
{
cerr << usage_text << endl;
cerr << help_text << endl;
cout << usage_text << endl;
cout << help_text << endl;
exit(0);
}
......@@ -265,78 +296,128 @@ int main(int iargc, char* argv[])
opt.doskip=cmdline.optset(9);
opt.fitexp=cmdline.optset(10);
opt.tcexp=cmdline.double_arg(10);
opt.itype=cmdline.string_arg(11);
opt.otype=cmdline.string_arg(12);
/*----------------------------------------------------------------------*/
// read file names
TFXX_assert(cmdline.extra(),"ERROR: missing signal file name!");
std::string signalname=cmdline.next();
TFXX_assert((cmdline.extra()||opt.fittrend||opt.fitoffset),
"ERROR: missing trial signal file name!");
Tnamelist namelist;
while (cmdline.extra()) { namelist.push_back(cmdline.next()); }
//std::string signalname=cmdline.next();
tfxx::cmdline::Tparsed infiles(tfxx::cmdline::parse_cmdline(cmdline, keys));
//TFXX_assert((cmdline.extra()||opt.fittrend||opt.fitoffset),
// "ERROR: missing trial signal file name!");
//Tnamelist namelist;
//while (cmdline.extra()) { namelist.push_back(cmdline.next()); }
// read signal file
tfxx::cmdline::Tparsed::const_iterator infile(infiles.begin());
TFXX_assert((infile!=infiles.end()||opt.fittrend||opt.fitoffset),
"ERROR (sigfit): missing trial signal file name");
std::string tracelist("1");
if (infile->haskey(tracekey)) { tracelist=infile->value(tracekey); }
tfxx::RangeListStepper<int> rls(tfxx::string::rangelist<int>(tracelist));
TFXX_assert(rls.valid(), "Illegal tracelist");
std::string signalname(infile->name);
int signaltrace = rls.current();
std::string signalformat(opt.itype);
if (infile->haskey(formatkey)) { signalformat=infile->value(formatkey); }
Tbundle signal;
sff::FREE signalfilefree, signaltracefree;
bool hasSignalFileFree;
{
if (opt.verbose)
{ cout << "read signal file \"" << signalname << "\"" << endl; }
std::ifstream is(signalname.c_str());
sff::FileHeader fileheader(is);
if (opt.verbose)
{ cout << " read trace" << endl; }
sff::InputWaveform<Tseries> inputwaveform;
try {
inputwaveform.read(is);
}
catch(...) {
cerr << "ERROR (sigfit): read failed!" << endl;
throw;
if (opt.verbose)
{
cout << "sigfit: open signal file " << signalname << " of format "
<< signalformat << endl;
}
sff::TraceHeader traceheader=inputwaveform.header();
signal.header=traceheader.wid2();
signal=inputwaveform.series();
bool last=traceheader.last();
std::string wid2line=signal.header.line();
if (opt.verbose)
std::ifstream ifs(signalname.c_str(),
datrw::ianystream::openmode(signalformat));
datrw::ianystream is(ifs, signalformat);
if (is.hasfree())
{
cout << " " << wid2line.substr(0,69) << endl;
if (!last)
{
cout << "there are more traces in that files" << endl
<< "ignoring them..." << endl;
}
is >> signalfilefree;
hasSignalFileFree = true;
}
int itrace=1;
while (itrace!=signaltrace) { is.skipseries(); ++itrace; }
TFXX_assert(is.good(), "Illegal trace for signal input");
if (opt.verbose) { cout << "read trace" << endl; }
is >> signal;
is >> signal.header;
if (is.hasfree()) { is >> signaltracefree; }
}
++infile;
// read trial files
if (opt.verbose)
{ cout << namelist.size() << " trial data files to read" << endl; }
Tbundlevec bundlevec;
for (Tnamelist::const_iterator i=namelist.begin(); i!=namelist.end(); i++)
while (infile != infiles.end())
{
if (opt.verbose)
{ cout << "read trial data file \"" << *i << "\"" << endl; }
std::ifstream is(i->c_str());
sff::FileHeader fileheader(is);
bool last=false;
while(!last)
std::string itype(opt.itype);
if (infile->haskey(formatkey)) { itype=infile->value(formatkey); }
if (opt.verbose)
{
cout << "sigfit: Open trial signal file " << infile->name << " of format "
<< itype << endl;
}
bool selectedTraces=infile->haskey(tracekey);
// no traces of file were selected -> read the entire file
if (! selectedTraces)
{
if (opt.verbose)
{ cout << " read trial signal #" << bundlevec.size()+1 << endl; }
sff::InputWaveform<Tseries> inputwaveform(is);
sff::TraceHeader traceheader=inputwaveform.header();
Tbundle bundle;
bundle.header=traceheader.wid2();
bundle=inputwaveform.series();
bundlevec.push_back(bundle);
last=traceheader.last();
std::string wid2line=bundle.header.line();
if (opt.verbose)
{ cout << " " << wid2line.substr(0,69) << endl; }
std::ifstream ifs(infile->name.c_str(),
datrw::ianystream::openmode(itype));
datrw::ianystream is(ifs, itype);
int traceNum = 1;
while (is.good())
{
if (opt.verbose)
{
cout << "sigfit: Reading trace #" << traceNum << " ..." << endl;
}
Tbundle bundle;
is >> bundle;
is >> bundle.header;
bundlevec.push_back(bundle);
if (opt.verbose)
{
cout << " " << bundle.header.line().substr(0,69) << endl;
}
++traceNum;
}
} else
// traces were selected by rangelist -> read only selected traces
{
tfxx::RangeListStepper<int> rls(tfxx::string::rangelist<int>(
infile->value(tracekey)));
while (rls.valid())
{
std::ifstream ifs(infile->name.c_str(),
datrw::ianystream::openmode(itype));
datrw::ianystream is(ifs, itype);
int traceNum = 1;
// spool to selected trace
while (traceNum!=rls.current()) { is.skipseries(); ++traceNum; }
TFXX_assert(is.good(), "Illegal trace selected.");
if (opt.verbose)
{
cout << "sigfit: Reading trace #" << traceNum << " ..." << endl;
}
Tbundle bundle;
is >> bundle;
is >> bundle.header;
bundlevec.push_back(bundle);
if (opt.verbose)
{
cout << " " << bundle.header.line().substr(0,69) << endl;
}
++rls;
}
}
++infile;
}
unsigned int ntracesread=bundlevec.size();
// add synthetic trial files
......@@ -575,17 +656,33 @@ int main(int iargc, char* argv[])
{
if (opt.verbose)
{
cout << "write residual waveform to "
<< opt.residualname << endl;
cout << "write residual waveform to " << opt.residualname
<< " in format " << opt.otype << endl;
}
std::ofstream ofs(opt.residualname.c_str(),
datrw::oanystream::openmode(opt.otype));
datrw::oanystream os(ofs, opt.otype);
if (os.handlesfilefree())
{
sff::FREE residualFileFree;
residualFileFree.append(SIGFIT_VERSION);
residualFileFree.append(SIGFIT_CVSID);
if (hasSignalFileFree)
{
residualFileFree.append("Signal file free block:");
residualFileFree.append(signalfilefree);
}
os << residualFileFree;
}
std::ofstream ofs(opt.residualname.c_str());
sff::SFFostream<Tseries> os(ofs);
os << signal;
os << synthetics;
os << residual;
os << signal.header;
os << signal.series();
os << synthetics.header;
os << synthetics.series();
os << residual.header;
os << residual.series();
for (Tbundlevec::const_iterator i=bundlevec.begin();
i!=bundlevec.end(); i++)
{ os << *i; }
{ os << i->header; os << *i; }
for (int i=1; i<=N; ++i)
{
Tbundle psyn=bundlevec[i-1];
......@@ -593,11 +690,15 @@ int main(int iargc, char* argv[])
char chan[6];
std::sprintf(chan,"s%1.1d", i);
psyn.header.channel=std::string(chan);
os << psyn.header;
os << psyn;
}
os << corrsignal;
os << corrsynthetics;
os << correction;
os << corrsignal.header;
os << corrsignal.series();
os << corrsynthetics.header;
os << corrsynthetics.series();
os << correction.header;
os << correction.series();
}
// calculate rms values
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment