Geocentric.hpp

Go to the documentation of this file.

#if !defined(GEOGRAPHICLIB_GEOCENTRIC_HPP)
#define GEOGRAPHICLIB_GEOCENTRIC_HPP 1

#include <vector>
#include <GeographicLib/Constants.hpp>

namespace GeographicLib {

  class GEOGRAPHICLIB_EXPORT Geocentric {
  private:
    typedef Math::real real;
    friend class LocalCartesian;
    friend class MagneticCircle; // MagneticCircle uses Rotation
    friend class MagneticModel;  // MagneticModel uses IntForward
    friend class GravityCircle;  // GravityCircle uses Rotation
    friend class GravityModel;   // GravityModel uses IntForward
    friend class NormalGravity;  // NormalGravity uses IntForward
    static const size_t dim_ = 3;
    static const size_t dim2_ = dim_ * dim_;
    real _a, _f, _e2, _e2m, _e2a, _e4a, _maxrad;
    static void Rotation(real sphi, real cphi, real slam, real clam,
                         real M[dim2_]);
    static void Rotate(real M[dim2_], real x, real y, real z,
                       real& X, real& Y, real& Z) {
      // Perform [X,Y,Z]^t = M.[x,y,z]^t
      // (typically local cartesian to geocentric)
      X = M[0] * x + M[1] * y + M[2] * z;
      Y = M[3] * x + M[4] * y + M[5] * z;
      Z = M[6] * x + M[7] * y + M[8] * z;
    }
    static void Unrotate(real M[dim2_], real X, real Y, real Z,
                         real& x, real& y, real& z)  {
      // Perform [x,y,z]^t = M^t.[X,Y,Z]^t
      // (typically geocentric to local cartesian)
      x = M[0] * X + M[3] * Y + M[6] * Z;
      y = M[1] * X + M[4] * Y + M[7] * Z;
      z = M[2] * X + M[5] * Y + M[8] * Z;
    }
    void IntForward(real lat, real lon, real h, real& X, real& Y, real& Z,
                    real M[dim2_]) const;
    void IntReverse(real X, real Y, real Z, real& lat, real& lon, real& h,
                    real M[dim2_]) const;

  public:

    Geocentric(real a, real f);

    Geocentric() : _a(-1) {}

    void Forward(real lat, real lon, real h, real& X, real& Y, real& Z)
      const {
      if (Init())
        IntForward(lat, lon, h, X, Y, Z, NULL);
    }

    void Forward(real lat, real lon, real h, real& X, real& Y, real& Z,
                 std::vector<real>& M)
      const {
      if (!Init())
        return;
      if (M.end() == M.begin() + dim2_) {
        real t[dim2_];
        IntForward(lat, lon, h, X, Y, Z, t);
        std::copy(t, t + dim2_, M.begin());
      } else
        IntForward(lat, lon, h, X, Y, Z, NULL);
    }

    void Reverse(real X, real Y, real Z, real& lat, real& lon, real& h)
      const {
      if (Init())
        IntReverse(X, Y, Z, lat, lon, h, NULL);
    }

    void Reverse(real X, real Y, real Z, real& lat, real& lon, real& h,
                 std::vector<real>& M)
      const {
      if (!Init())
        return;
      if (M.end() == M.begin() + dim2_) {
        real t[dim2_];
        IntReverse(X, Y, Z, lat, lon, h, t);
        std::copy(t, t + dim2_, M.begin());
      } else
        IntReverse(X, Y, Z, lat, lon, h, NULL);
    }


    bool Init() const { return _a > 0; }
    Math::real MajorRadius() const
    { return Init() ? _a : Math::NaN(); }

    Math::real Flattening() const
    { return Init() ? _f : Math::NaN(); }

    static const Geocentric& WGS84();
  };

} // namespace GeographicLib

#endif  // GEOGRAPHICLIB_GEOCENTRIC_HPP