QGpCoreWave/RefractionDippingModel.h

Go to the documentation of this file.

/***************************************************************************
**
**  This file is part of QGpCoreWave.
**
**  This library is free software; you can redistribute it and/or
**  modify it under the terms of the GNU Lesser General Public
**  License as published by the Free Software Foundation; either
**  version 2.1 of the License, or (at your option) any later version.
**
**  This file is distributed in the hope that it will be useful, but WITHOUT
**  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
**  FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
**  License for more details.
**
**  You should have received a copy of the GNU Lesser General Public
**  License along with this library; if not, write to the Free Software
**  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
**
**  See http://www.geopsy.org for more information.
**
**  Created : 2004-09-17
**  Authors :
**    Marc Wathelet
**    Marc Wathelet (ULg, Liège, Belgium)
**    Marc Wathelet (LGIT, Grenoble, France)
**
***************************************************************************/

#ifndef REFRACTIONDIPPINGMODEL_H
#define REFRACTIONDIPPINGMODEL_H

#include <math.h>
#include <QGpCoreTools.h>
#include "Seismic1DModel.h"
#include "QGpCoreWaveDLLExport.h"

namespace QGpCoreWave {

class TiltNode;
class TiltPath;

class QGPCOREWAVE_EXPORT RefractionDippingModel : public GeophysicalModel
{
  TRANSLATIONS("RefractionDippingModel");
public:
  RefractionDippingModel();
  RefractionDippingModel(int nLayers);
  RefractionDippingModel(const RefractionDippingModel& o);
  ~RefractionDippingModel();

  virtual GeophysicalModel * clone() const {return new RefractionDippingModel(*this);}

  int layerCount() const {return _layerCount;}
  virtual bool isEmpty() const {return layerCount()==0;}
  double depthLeft(int iLayer) const {return _depthL[iLayer];}
  double setDepthLeft(int iLayer, double d) {return _depthL[iLayer]=d;}
  double depthRight(int iLayer) const {return _depthR[iLayer];}
  double setDepthRight(int iLayer, double d) {return _depthR[iLayer]=d;}
  double slowness(int iLayer) const {return _slow[iLayer];}
  double setSlowness(int iLayer, double s) {return _slow[iLayer]=s;}

  void setXLeft(double x) {_xL=x;}
  double xLeft() const {return _xL;}

  void setXRight(double x) {_xR=x;}
  double xRight() const {return _xR;}

  void fromSeismic1DModel(const Seismic1DModel& ref, const QVector<double>& pitch,
                         Seismic1DModel::SlownessType slowType);

  enum Direction {LeftToRight, RightToLeft};

  void begin();
  inline bool isCriticalLayer(Direction dir, int iLayer) const;
  double propagate(Direction dir, int iLayer, int iDeepestLayer, double& x) const;
  void end();

  double travelTime(const TiltNode& src, const TiltNode& rec, int& minTimeLayer);
  Curve<Point2D> ray(const TiltNode& src, const TiltNode& rec);

  static QString formatHelp();
  virtual QString toString() const;
  virtual bool fromStream(QTextStream& s, QString * comments=0);
  virtual void toStream(QTextStream& s) const;
  void fromStream(QDataStream& s);
  void toStream(QDataStream& s) const;

  class RefractionContext : public GeophysicalContext
  {
  public:
    RefractionContext();
    RefractionDippingModel * model() const {return static_cast<RefractionDippingModel *>(_model);}
    virtual QString helpCode() const;
  };

  virtual GeophysicalContext * expressionContext() const {return new RefractionContext;}

  class RefractionStorage : public ExpressionStorage
  {
  public:
    RefractionStorage(RefractionContext * context) {_context=context;}
    virtual bool isReadOnly() const {return false;}
  protected:
    RefractionContext * _context;
  };

  class VariableN : public RefractionStorage
  {
  public:
    VariableN(RefractionContext * context) : RefractionStorage(context) {}
    virtual bool isReadOnly() const {return true;}
    virtual QVariant value(const QString&) const;
    virtual void setValue(const QString&, const QVariant& ) {}
  };

  class VariableDL : public RefractionStorage
  {
  public:
    VariableDL(RefractionContext * context) : RefractionStorage(context) {}
    virtual QVariant value(const QString& index) const;
    virtual void setValue(const QString&, const QVariant& );
  };

  class VariableDR : public RefractionStorage
  {
  public:
    VariableDR(RefractionContext * context) : RefractionStorage(context) {}
    virtual QVariant value(const QString& index) const;
    virtual void setValue(const QString&, const QVariant& );
  };

  class VariableV : public RefractionStorage
  {
  public:
    VariableV(RefractionContext * context) : RefractionStorage(context) {}
    virtual QVariant value(const QString& index) const;
    virtual void setValue(const QString&, const QVariant& );
  };
private:
  void initMembers();
  void allocateData();
  void deleteData();
  bool setValue(double& var, const StringSection& field, int iLayer, const QString& varName, bool optional);
  double travelTime(const TiltPath *& srcPaths, const TiltPath *& recPaths, int& minTimeLayer);
  bool node2paths(const TiltNode& src, const TiltNode& rec,
                   const TiltPath *& srcPaths, const TiltPath *& recPaths) const;

  // number of layers including the bottom half space
  int _layerCount;
  // pointer to depth vector (positive values) counted at profile Left
  double * _depthL;
  // pointer to depth vector (positive values) counted at profile Right
  double * _depthR;
  // pointer to slownesses vector
  double * _slow;
  double _xL, _xR;

  // Min and max x of rays along interfaces
  double * _minX;
  double * _maxX;

  // Cosine of the angle (to avoid multiple useless computations)
  double * _cosangle;
  // Sine of the angle (to avoid multiple useless computations)
  double * _sinangle;
  // ratio of slownesses at the interface
  double * _n;
  /* Abscissa to substract to abscissae measured from the intercept point
      to obtain the abcissa from the reference profile.
      The intercept is the projection of the intersection of the reference
      profile with the upper interface, made perpendiculary to the upper
      interface, and downto the current interface */
  double * _shift;
  /* Apparent thickness of the layer calculated prependiculary to upper interface
      and at reference profile */
  double * _h;
  /* For each layer a vector of relative angles may be allocated
    if the interface is a good candidate for refraction at critical angle */
  double ** _critl;
  double ** _critr;
};

inline bool RefractionDippingModel::isCriticalLayer(Direction dir, int iLayer) const
{
  TRACE;
  switch(dir) {
  case LeftToRight:
    return _critl[iLayer];
    break;
  default:
    return _critr[iLayer];
    break;
  }
}

} // namespace QGpCoreWave

#endif // REFRACTIONDIPPINGMODEL_H