Commit 138f1cd4 authored by thomas.forbriger's avatar thomas.forbriger

libstfinv [WP][FIX]: use correct (padded) number of samples

parent 4fca620a
......@@ -576,6 +576,8 @@ namespace stfinv {
*/
void STFFourierDomainEngine::taperstf()
{
STFINV_debug(Mdebug&16, "STFFourierDomainEngine::taperstf()",
STFINV_value(Mapplystftaper));
if (this->Mapplystftaper)
{
// transform correction filter to time domain
......@@ -592,8 +594,10 @@ namespace stfinv {
* t2. Behind this is the concept that time series are implicitely
* periodic in discrete Fourier transform.
*/
// samples in series
int nsamples=thestfseries.size(0);
// duration of time series
double T=this->dt()*this->nsamples();
double T=this->dt()*nsamples;
// make sure taper definition is not longer than time series.
STFINV_assert((this->Mtt4-this->Mtt1)<T,
"Taper width is larger than length of time series");
......@@ -603,15 +607,26 @@ namespace stfinv {
double tt1=this->Mtt1+T;
double tt2=this->Mtt2+T;
// sample to start at
int l3=int(std::floor(tt3/T));
int l2=int(std::ceil(tt2/T));
int l3=int(std::floor(tt3/this->dt()));
int l2=int(std::ceil(tt2/this->dt()));
STFINV_debug(Mdebug&16, "STFFourierDomainEngine::taperstf()",
STFINV_value(T) << "\n " <<
STFINV_value(this->dt()) << "\n " <<
STFINV_value(nsamples) << "\n " <<
STFINV_value(tt3) << "\n " <<
STFINV_value(tt4) << "\n " <<
STFINV_value(tt1) << "\n " <<
STFINV_value(tt2) << "\n " <<
STFINV_value(l3) << "\n " <<
STFINV_value(l2));
for (int l=l3; l<=l2; ++l)
{
// time of sample
double t=l*this->dt();
// index to series
int rl = l%this->nsamples();
if (rl < 0) { rl += this->nsamples(); }
int rl = l%nsamples;
int dl = rl;
if (rl < 0) { rl += nsamples; }
STFINV_assert(((rl >= thestfseries.f(0)) &&
(rl <= thestfseries.l(0))),
"Index out of range. Internal programming error! "
......@@ -632,6 +647,12 @@ namespace stfinv {
}
// apply to series
thestfseries(rl) = thestfseries(rl)*factor;
STFINV_debug(Mdebug&16, "STFFourierDomainEngine::taperstf()",
STFINV_value(t) << " " <<
STFINV_value(l) << " " <<
STFINV_value(dl) << " " <<
STFINV_value(rl) << " " <<
STFINV_value(factor));
}
// transform correction filter to Fourier domain
Mfftenginestf.r2c();
......
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