rfvel.cpp

Go to the documentation of this file.

#include "RFVelCalc.h"
#include "CalcRecConf.h"
#include "SeismicDataComp.h"
#include "FatalException.h"
#include <iostream>
#include <fstream>
#include <string>
#include "Util.h"

using namespace std;
using namespace gplib;

string version = "$Id: rfvel.cpp 1839 2010-03-05 17:04:34Z mmoorkamp $";

int main()
  {
    cout
        << " This is rfvel: Calculate absolute velocties from receiver functions"
        << endl; // write some info
    cout << " Reads in radial and vertical components in SAC format" << endl;
    cout
        << " Writes a file with ending .vel containing apparent velocities and corresponding periods"
        << endl;
    cout
        << " The settings for the calculation method, omega, sigma and shift are taken from calcrec.conf"
        << endl;
    cout << " This is Version: " << version << endl << endl;
    string modelfilename;

    CalcRecConf Config;
    SeismicDataComp Radial, Vertical;
    string outfilename;
    double slowness;
    try
      {
        string radfilename, verfilename;
        //read the configuration file
        Config.GetData("calcrec.conf");
        //the default calculation method is spectral division
        RecCalc::trfmethod rfmethod = RecCalc::specdiv;
        if (Config.recmethod == "iterdecon")
          {
            rfmethod = RecCalc::iterdecon;
          }
        //get the names for the radial and vertical component seismograms
        radfilename = AskFilename("Radial Component: ");
        verfilename = AskFilename("Vertical Component: ");
        outfilename = radfilename + ".vel";
        //we have to ask for the slowness of the incoming P-wave
        cout << "Slowness in s/km: ";
        cin >> slowness;
        //read in the seismograms
        Radial.ReadData(radfilename);
        Vertical.ReadData(verfilename);
        ofstream outfile(outfilename.c_str());
        RFVelCalc RFVel(Config.sigma, Config.cc, rfmethod);

        ttsdata AppVel;
        //do the calculation and write to a file
        RFVel.CalcRFVel(slowness, Radial, Vertical, AppVel);
        RFVel.WriteVelocities(outfilename);
      } catch (FatalException &e)
      {
        cerr << e.what() << endl; // if something fails print error
        return -1; // and stop execution
      }
  }