initialize()

void stfinv::STFFourierDomainEngine::initialize ( )

private

initialize work space

Create FFT engines.

Number of samples

For the FFT we usually use a number of samples larger than that of the underlying time series. This way we at least partially avoid the undesired effect of the cyclic discrete Fourier transform, which can lead to wrap-around. Currently there are three ways to modify the number of samples actually used:

1. Option fpad specifies an padding factor. The number of samples simply is larger than the number of samples of the underlying time series by this factor.
2. Option fpow2: If set, the next power of two larger than the fpad*nsamples is used.
3. Option fdiv specifies a divisor. If set and if fpow2 is not set the number of samples will be the next integer multiple of the divisor larger than fpad*nsamples.

Workspace

Two FFT engines will be created. They are only used one-way, since input data is not altered during processing input time series only have to be transformed once to Fourier domain. Results of processing have to be transformed to time domain.

1. One engine (STFFourierDomainEngine::Mfftengineinput) being shared by recorded data and synthetic data, because both have to be transformed to Fourier domain at once. Transformation to Fourier domain takes place in STFFourierDomainEngine::fftinput.
2. One engine (STFFourierDomainEngine::Mfftengineoutput) being shared by the stf and the convolved synthetics, because both have to be transformed to the time domain at once. Transformation to time domain takes place in STFFourierDomainEngine::fftoutput.

Definition at line 191 of file stfinvfourier.cc.

References Mapplyshift, Mapplystftaper, stfinv::STFBaseEngine::Mdebug, Mfftengineinput, Mfftengineoutput, Mfftenginestf, Mtshift, Mtt1, Mtt2, Mtt3, Mtt4, stfinv::STFBaseEngine::npairs(), stfinv::STFBaseEngine::nreceivers(), stfinv::STFBaseEngine::nsamples(), stfinv::STFBaseEngine::parameter(), stfinv::STFBaseEngine::parameterisset(), stfinv::tools::secomtospace(), STFINV_assert, STFINV_debug, STFINV_value, stfseries(), and stfspec().

Referenced by STFFourierDomainEngine().

  {
    // extract parameter values
    // ------------------------
    // padding factor
    double padfactor;
    {
      std::istringstream is(this->parameter("fpad","1.5"));
      is >> padfactor;
    }
    STFINV_assert(padfactor >= 1.,
                  "ERROR: parameter for option \"fpad\" not larger or equal 1");
    // flag: use power of two
    bool poweroftwo=(this->parameter("fpow2","false")=="true");
    // flag: use integer multiples of divisor
    bool divisorset=this->parameterisset("fdiv");
    // number of samples shall be integer power of divisor
    unsigned int divisor=1;
    if (divisorset)
    {
      std::istringstream is (this->parameter("fdiv","100")); 
      is >> divisor;
      STFINV_assert(divisor > 0,
                  "ERROR: parameter for option \"fdiv\" not larger than 0");
    }

    // define number of samples to be used by Fourier engine
    // -----------------------------------------------------
    // use at least the number of samples sepcified by the padfactor
    unsigned int nsamples=padfactor*this->nsamples();
    if (nsamples<this->nsamples()) { nsamples=this->nsamples(); }
    // power of two has precendence
    if (poweroftwo)
    {
      unsigned int n=2;
      while (n<nsamples) { n *= 2; }
      nsamples=n;
    }
    else if (divisorset)
    {
      unsigned int rest=nsamples % divisor;
      nsamples -= rest;
      nsamples += divisor;
    } // if (poweroftwo) or (divisorset)

    // allocate workspace by creating engines
    // --------------------------------------
    Mfftengineinput=fourier::fft::DRFFTWAFFArrayEngine(
                                     2*this->nreceivers()+this->npairs(),
                                                       nsamples);
    Mfftengineoutput=fourier::fft::DRFFTWAFFArrayEngine(
                                     1+this->nreceivers()+this->npairs(),
                                                        nsamples);
    Mfftenginestf=fourier::fft::DRFFTWAFFArrayEngine(this->stfseries(),
                                                     this->stfspec());

    // set time shift for STF if requested
    // -----------------------------------
    {
      std::istringstream is(this->parameter("tshift","0."));
      is >> Mtshift;
    }
    Mapplyshift=this->parameterisset("tshift");
      
    // taper definition
    // ----------------
    {
      std::istringstream is(stfinv::tools::secomtospace(
                            this->parameter("irtap","0.,1.,2.,3.")));
      is >> Mtt1 >> Mtt2 >> Mtt3 >> Mtt4;
    }
    Mapplystftaper=this->parameterisset("irtap");
    STFINV_assert((Mtt1<=Mtt2) && (Mtt2<=Mtt3) && (Mtt3<=Mtt4),
                  "ERROR: taper definition times not in increasing order");

    STFINV_debug(Mdebug&1, "STFFourierDomainEngine::initialize()",
                 STFINV_value(this->nsamples()) << "\n  " <<
                 STFINV_value(padfactor) << "\n  " <<
                 STFINV_value(divisor) << "\n  " <<
                 STFINV_value(nsamples) << "\n  " <<
                 STFINV_value(Mapplyshift) << "\n  " <<
                 STFINV_value(Mtshift) << "\n  " <<
                 STFINV_value(Mapplystftaper) << "\n  " <<
                 STFINV_value(Mtt1) << "\n  " <<
                 STFINV_value(Mtt2) << "\n  " <<
                 STFINV_value(Mtt3) << "\n  " <<
                 STFINV_value(Mtt4));
  } // void STFFourierDomainEngine::initialize() 

Here is the call graph for this function:

Here is the caller graph for this function: