SeismicModel.h

Go to the documentation of this file.

#ifndef CSEISMICMODEL_H
#define CSEISMICMODEL_H

#include "types.h"

namespace gplib
  {
    /** \addtogroup seistools Seismic data analysis and modeling */
    /* @{ */

    //! The class SeismicModel is the base class for some of the model format for seismic codes
    /*!
     * We store the basic seismic parameters here and implement some of the common functionailty.
     *
     */
    class SeismicModel
      {
    public:
      //! Do we want to calculate the arrival of a direct S-Wave or a P-wave
      enum tArrivalType
        {
        DirectS, DirectP
        };
    private:
      trealdata SVelErrors;
      trealdata ThickErrors;
      bool
      FuzzComp(const double elem1, const double elem2, const double tolerance);
      trealdata PVelocity;
      trealdata SVelocity;
      trealdata Density;
      trealdata Thickness;
      trealdata Qp;
      trealdata Qs;
      double SourceDepth;
      double dt;
      unsigned int npts;
    public:
      //! Get the number of points for synthetic seismogram calculation
      unsigned int GetNpts() const
        {
          return npts;
        }
      //! Set the number of points for synthetic seismogram calculation
      void SetNpts(const unsigned int s)
        {
          npts = s;
        }
      //! Get the depth to the seismic source that generates the wavefield
      double GetSourceDepth() const
        {
          return SourceDepth;
        }
      //! Set the depth to the seismic source
      void SetSourceDepth(const double s)
        {
          SourceDepth = s;
        }
      //! Get the time between two samples in s, this is for synthetic forward calculations
      double GetDt() const
        {
          return dt;
        }
      //! Set the time between two samples in s, this is for synthetic forward calculations
      void SetDt(const double s)
        {
          dt = s;
        }
      //! Read-only access to the vector of P-velocities in km/s in each layer
      const trealdata &GetPVelocity() const
        {
          return PVelocity;
        }
      //! Read-write access to the vector of P-velocities in km/s in each layer
      trealdata &SetPVelocity()
        {
          return PVelocity;
        }
      //! Read-only access to the vector of S-velocities in km/s in each layer
      const trealdata &GetSVelocity() const
        {
          return SVelocity;
        }
      //! Read-write access to the vector of S-velocities in km/s in each layer
      trealdata &SetSVelocity()
        {
          return SVelocity;
        }
      //! Read-only access to the vector of densities in g/cm^3 in each layer
      const trealdata &GetDensity() const
        {
          return Density;
        }
      //! Read-write access to the vector of densities in g/cm^3 in each layer
      trealdata &SetDensity()
        {
          return Density;
        }
      //! Read-only access to the vector of thicknesses in km in each layer
      const trealdata &GetThickness() const
        {
          return Thickness;
        }
      //! Read-write access to the vector of thicknesses in km in each layer
      trealdata &SetThickness()
        {
          return Thickness;
        }
      //! Read-only access to the vector of P quality factors for each layer
      const trealdata &GetQp() const
        {
          return Qp;
        }
      //! Read-write access to the vector of P quality factors for each layer
      trealdata &SetQp()
        {
          return Qp;
        }
      //! Read-only access to the vector of S quality factors for each layer
      const trealdata &GetQs() const
        {
          return Qs;
        }
      //! Read-write access to the vector of S quality factors for each layer
      trealdata &SetQs()
        {
          return Qs;
        }
      //! Set error bars on S-velocities for plotting
      void SetSVelErrors(const trealdata &sve)
        {
          SVelErrors = sve;
        }
      //! Set error bars on Thicknesses for plotting
      void SetThickErrors(const trealdata &te)
        {
          ThickErrors = te;
        }
      //! Returns the index of the layer that correponds to depth in km
      int FindLayer(const double depth);
      //! Given a slowness in s/km and a wave type calculate the distance that matches this slowness
      double MatchSlowness(const double slowness, const tArrivalType mode);
      double CalcArrival(const tArrivalType mode, const double recdist);
      double CalcTravelTime(const tArrivalType mode, const double sdepth,
          const double rdepth, const double p);
      //! Read the model in its native format from a file
      virtual void ReadModel(const std::string filename)=0;
      //! Write the model in its native format to a file
      virtual void WriteModel(const std::string filename)=0;
      //! Write out a script that performs a forward calculation for the current model
      virtual void WriteRunFile(const std::string &filename)=0;
      //! Write out an ascii file for plotting with xmgrace etc.
      void WritePlot(const std::string &filename);
      //! Write out an ascii file with error bars for plotting with xmgrace etc.
      void PlotVelWithErrors(const std::string &filename);
      //! Init provides a convenient way to allocate memory in all structures for a given number of layers
      void Init(const int nlayers);
      SeismicModel(const int nlayers = 0);
      SeismicModel(const SeismicModel& source);
      SeismicModel& operator=(const SeismicModel& source);
      virtual ~SeismicModel();
      };
  /* @} */
  }
#endif // CSEISMICMODEL_H