SphericalHarmonic1.hpp

Go to the documentation of this file.

#if !defined(GEOGRAPHICLIB_SPHERICALHARMONIC1_HPP)
#define GEOGRAPHICLIB_SPHERICALHARMONIC1_HPP 1

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

namespace GeographicLib {

  class GEOGRAPHICLIB_EXPORT SphericalHarmonic1 {
  public:
    enum normalization {
      FULL = SphericalEngine::FULL,
      SCHMIDT = SphericalEngine::SCHMIDT,
    };

  private:
    typedef Math::real real;
    SphericalEngine::coeff _c[2];
    real _a;
    unsigned _norm;

  public:
    SphericalHarmonic1(const std::vector<real>& C,
                       const std::vector<real>& S,
                       int N,
                       const std::vector<real>& C1,
                       const std::vector<real>& S1,
                       int N1,
                       real a, unsigned norm = FULL)
      : _a(a)
      , _norm(norm) {
      if (!(N1 <= N))
        throw GeographicErr("N1 cannot be larger that N");
      _c[0] = SphericalEngine::coeff(C, S, N);
      _c[1] = SphericalEngine::coeff(C1, S1, N1);
    }

    SphericalHarmonic1(const std::vector<real>& C,
                       const std::vector<real>& S,
                       int N, int nmx, int mmx,
                       const std::vector<real>& C1,
                       const std::vector<real>& S1,
                       int N1, int nmx1, int mmx1,
                       real a, unsigned norm = FULL)
      : _a(a)
      , _norm(norm) {
      if (!(nmx1 <= nmx))
        throw GeographicErr("nmx1 cannot be larger that nmx");
      if (!(mmx1 <= mmx))
        throw GeographicErr("mmx1 cannot be larger that mmx");
      _c[0] = SphericalEngine::coeff(C, S, N, nmx, mmx);
      _c[1] = SphericalEngine::coeff(C1, S1, N1, nmx1, mmx1);
    }

    SphericalHarmonic1() {}

    Math::real operator()(real tau, real x, real y, real z) const {
      real f[] = {1, tau};
      real v = 0;
      real dummy;
      switch (_norm) {
      case FULL:
        v = SphericalEngine::Value<false, SphericalEngine::FULL, 2>
          (_c, f, x, y, z, _a, dummy, dummy, dummy);
        break;
      case SCHMIDT:
        v = SphericalEngine::Value<false, SphericalEngine::SCHMIDT, 2>
          (_c, f, x, y, z, _a, dummy, dummy, dummy);
        break;
      }
      return v;
    }

    Math::real operator()(real tau, real x, real y, real z,
                          real& gradx, real& grady, real& gradz) const {
      real f[] = {1, tau};
      real v = 0;
      switch (_norm) {
      case FULL:
        v = SphericalEngine::Value<true, SphericalEngine::FULL, 2>
          (_c, f, x, y, z, _a, gradx, grady, gradz);
        break;
      case SCHMIDT:
        v = SphericalEngine::Value<true, SphericalEngine::SCHMIDT, 2>
          (_c, f, x, y, z, _a, gradx, grady, gradz);
        break;
      }
      return v;
    }

    CircularEngine Circle(real tau, real p, real z, bool gradp) const {
      real f[] = {1, tau};
      switch (_norm) {
      case FULL:
        return gradp ?
          SphericalEngine::Circle<true, SphericalEngine::FULL, 2>
          (_c, f, p, z, _a) :
          SphericalEngine::Circle<false, SphericalEngine::FULL, 2>
          (_c, f, p, z, _a);
        break;
      case SCHMIDT:
      default:                  // To avoid compiler warnings
        return gradp ?
          SphericalEngine::Circle<true, SphericalEngine::SCHMIDT, 2>
          (_c, f, p, z, _a) :
          SphericalEngine::Circle<false, SphericalEngine::SCHMIDT, 2>
          (_c, f, p, z, _a);
        break;
      }
    }

    const SphericalEngine::coeff& Coefficients() const
    { return _c[0]; }
    const SphericalEngine::coeff& Coefficients1() const
    { return _c[1]; }
  };

} // namespace GeographicLib

#endif  // GEOGRAPHICLIB_SPHERICALHARMONIC1_HPP