Latitude.cpp

 
/* SPDX-License-Identifier: CC0-1.0 */
#include "Latitude.h"
 
#include <cmath>
 
#include "Constants.h"
#include "Distance.h"
#include "IException.h"
#include "IString.h"
#include "PvlGroup.h"
#include "QString"
#include "SpecialPixel.h"
#include "Target.h"
 
namespace Isis {
  
  Latitude::Latitude() : Angle() {
    
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_errors = AllowPastPole;
  }
 
 
  Latitude::Latitude(double latitude, Angle::Units latitudeUnits,
                     ErrorChecking errors) : Angle() {
 
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_errors = errors;
 
    setPlanetocentric(latitude, latitudeUnits);
  }
 
 
  Latitude::Latitude(Angle latitude, ErrorChecking errors) : Angle() {
    
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_errors = errors;
 
    setPlanetocentric(latitude.radians(), Radians);
  }
 
 
  Latitude::Latitude(Angle latitude, PvlGroup mapping,
                     ErrorChecking errors) : Angle(latitude) {
    
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    if (mapping.hasKeyword("EquatorialRadius") && mapping.hasKeyword("PolarRadius")) {
      m_equatorialRadius = new Distance(toDouble(mapping["EquatorialRadius"][0]),
                                        Distance::Meters);
      m_polarRadius = new Distance(toDouble(mapping["PolarRadius"][0]),
                                   Distance::Meters);
    }
    else {
      try {
        PvlGroup radiiGrp = Target::radiiGroup(mapping["TargetName"][0]);
        m_equatorialRadius = new Distance(toDouble(radiiGrp["EquatorialRadius"][0]),
                                          Distance::Meters);
        m_polarRadius = new Distance(toDouble(radiiGrp["PolarRadius"][0]),
                                     Distance::Meters);
      }
      catch (IException &e) {
        QString msg = "Unable to create Latitude object from given mapping group.";
        throw IException(e, IException::Unknown, msg, _FILEINFO_);
      }
    }
 
    m_errors = errors;
 
    if (mapping["LatitudeType"][0] == "Planetographic") {
      setPlanetographic(latitude.radians(), Radians);
    }
    else if (mapping["LatitudeType"][0] == "Planetocentric") {
      setPlanetocentric(latitude.radians(), Radians);
    }
    else {
      QString msg = "Latitude type [" + mapping["LatitudeType"][0] +
        "] is not recognized";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
  }
 
 
  Latitude::Latitude(double latitude,
                     PvlGroup mapping,
                     Angle::Units latitudeUnits,
                     ErrorChecking errors) : Angle(latitude, latitudeUnits) {
 
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    if (mapping.hasKeyword("EquatorialRadius") && mapping.hasKeyword("PolarRadius")) {
      m_equatorialRadius = new Distance(toDouble(mapping["EquatorialRadius"][0]),
                                        Distance::Meters);
      m_polarRadius = new Distance(toDouble(mapping["PolarRadius"][0]),
                                   Distance::Meters);
    }
    else {
      try {
        PvlGroup radiiGrp = Target::radiiGroup(mapping["TargetName"][0]);
        m_equatorialRadius = new Distance(toDouble(radiiGrp["EquatorialRadius"][0]),
                                          Distance::Meters);
        m_polarRadius = new Distance(toDouble(radiiGrp["PolarRadius"][0]),
                                     Distance::Meters);
      }
      catch (IException &e) {
        QString msg = "Unable to create Latitude object from given mapping group.";
        throw IException(e, IException::Unknown, msg, _FILEINFO_);
      }
    }
 
    m_errors = errors;
 
    if (mapping["LatitudeType"][0] == "Planetographic") {
      setPlanetographic(latitude, latitudeUnits);
    }
    else if (mapping["LatitudeType"][0] == "Planetocentric") {
      setPlanetocentric(latitude, latitudeUnits);
    }
    else {
      QString msg = "Latitude type [" + mapping["LatitudeType"][0] +
        "] is not recognized";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
  }
 
 
  Latitude::Latitude(double latitude,
                     Distance equatorialRadius, Distance polarRadius,
                     CoordinateType latType,
                     Angle::Units latitudeUnits,
                     ErrorChecking errors) : Angle(latitude, latitudeUnits) {
    
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_equatorialRadius = new Distance(equatorialRadius);
    m_polarRadius = new Distance(polarRadius);
 
    m_errors = errors;
 
    if (latType == Planetocentric) {
      setPlanetocentric(latitude, latitudeUnits);
    }
    else if (latType == Planetographic) {
      setPlanetographic(latitude, latitudeUnits);
    }
    else {
      QString msg = "Enumeration value [" + toString(latType) + "] is not a valid CoordinateType";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
  }
 
 
  Latitude::Latitude(const Latitude &latitudeToCopy) : Angle(latitudeToCopy) {
    
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_errors = latitudeToCopy.m_errors;
 
    if (latitudeToCopy.m_equatorialRadius) {
      m_equatorialRadius = new Distance(*latitudeToCopy.m_equatorialRadius);
    }
 
    if (latitudeToCopy.m_polarRadius) {
      m_polarRadius = new Distance(*latitudeToCopy.m_polarRadius);
    }
  }
 
 
  Latitude::~Latitude() {
    if (m_equatorialRadius) {
      delete m_equatorialRadius;
      m_equatorialRadius = NULL;
    }
 
    if (m_polarRadius) {
      delete m_polarRadius;
      m_polarRadius = NULL;
    }
  }
 
 
  double Latitude::planetocentric(Angle::Units units) const {
    return angle(units);
  }
 
 
  void Latitude::setPlanetocentric(double latitude, Angle::Units units) {
    setAngle(latitude, units);
  }
 
 
  double Latitude::planetographic(Angle::Units units) const {
    
    if (m_equatorialRadius == NULL || m_polarRadius == NULL) {
      QString msg = "Latitude [" + toString(true) + "] cannot "
          "be converted to Planetographic without the planetary radii, please "
          "use the other Latitude constructor.";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    if (*this > Angle(90.0, Angle::Degrees) ||
        *this < Angle(-90.0, Angle::Degrees)) {
      QString msg = "Latitudes outside of the -90/90 range cannot be converted "
          "between Planetographic and Planetocentric";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    if (!isValid()) {
      QString msg = "Invalid planetographic latitudes are not currently "
          "supported";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    double ographicLatitude = atan(tan(radians()) *
        (*m_equatorialRadius / *m_polarRadius) *
        (*m_equatorialRadius / *m_polarRadius));
 
    // This theoretically should just be an angle, but make it a Latitude so
    //   we can access angle
    return Latitude(ographicLatitude, Angle::Radians).angle(units);
  }
 
 
  void Latitude::setPlanetographic(double latitude, Angle::Units units) {
    
    if (m_equatorialRadius == NULL || m_polarRadius == NULL) {
      QString msg = "Latitude [" + Isis::toString(latitude) + " degrees] cannot be "
          "converted to Planetocentic without the planetary radii, please use "
          "the other Latitude constructor.";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    Angle inputAngle(latitude, units);
 
    if (inputAngle > Angle(90.0, Angle::Degrees) ||
        inputAngle < Angle(-90.0, Angle::Degrees)) {
      QString msg = "Latitudes outside of the -90/90 range cannot be converted "
          "between Planetographic and Planetocentric";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    // Since the Angle constructor handles special pixels, we will never get to this line of code
    // when passing in a special pixel.
    // Left this here just in case the functionality in Angle changes.
    if (IsSpecial(latitude)) {  
      QString msg = "Invalid planetographic latitudes are not currently "
          "supported";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    double ocentricLatitude = atan(tan(inputAngle.radians()) *
        (*m_polarRadius / *m_equatorialRadius) *
        (*m_polarRadius / *m_equatorialRadius));
 
    // Sometimes the trig functions return the negative of the expected value at the pole.
    if ((ocentricLatitude > 0) != (inputAngle.radians() > 0)) {
      ocentricLatitude *= -1;
    }
 
    setAngle(ocentricLatitude, Angle::Radians);
  }
 
 
  Latitude::ErrorChecking Latitude::errorChecking() const {
    return m_errors;
  }
 
 
  void Latitude::setErrorChecking(ErrorChecking errors) {
    m_errors = errors;
  }
 
 
  bool Latitude::inRange(Latitude min, Latitude max) const {
    
    // Validity check on the range
    if (min > max) {
      QString msg = "Minimum latitude [" + min.toString(true) + 
                    "] is greater than maximum latitude [" + 
                    max.toString(true) + "]";
      throw IException(IException::User, msg, _FILEINFO_);
    }
 
    // Provide a little wriggle room for precision problems
    Angle epsilon(DBL_EPSILON, Angle::Degrees);
    Latitude adjustedMin = min - epsilon;
    Latitude adjustedMax = max + epsilon;
 
    // Is this latitude between the min and the max
    return *this >= adjustedMin && *this <= adjustedMax;
  }
 
 
  Latitude& Latitude::operator=(const Latitude & latitudeToCopy) {
    
    if (this == &latitudeToCopy) return *this;
 
    m_equatorialRadius = NULL;
    m_polarRadius = NULL;
 
    m_errors = latitudeToCopy.m_errors;
 
    if (latitudeToCopy.m_equatorialRadius) {
      m_equatorialRadius = new Distance(*latitudeToCopy.m_equatorialRadius);
    }
 
    if (latitudeToCopy.m_polarRadius) {
      m_polarRadius = new Distance(*latitudeToCopy.m_polarRadius);
    }
 
    setPlanetocentric(latitudeToCopy.planetocentric());
 
    return *this;
  }
 
 
  Latitude Latitude::add(Angle angleToAdd, PvlGroup mapping) {
 
    CoordinateType latType;
 
    Distance equatorialRadius;
    Distance polarRadius;
    if (mapping.hasKeyword("EquatorialRadius") && mapping.hasKeyword("PolarRadius")) {
      equatorialRadius = Distance(toDouble(mapping["EquatorialRadius"][0]),
                                  Distance::Meters);
      polarRadius = Distance(toDouble(mapping["PolarRadius"][0]),
                             Distance::Meters);
    }
    else {
      try {
        PvlGroup radiiGrp = Target::radiiGroup(mapping["TargetName"][0]);
        equatorialRadius = Distance(toDouble(radiiGrp["EquatorialRadius"][0]),
                                    Distance::Meters);
        polarRadius = Distance(toDouble(radiiGrp["PolarRadius"][0]),
                               Distance::Meters);
      }
      catch (IException &e) {
        QString msg = "Unable to add angle to Latitude object from given mapping group.";
        throw IException(e, IException::Unknown, msg, _FILEINFO_);
      }
    }
 
    if (mapping["LatitudeType"][0] == "Planetocentric")
      latType = Planetocentric;
    else if (mapping["LatitudeType"][0] == "Planetographic")
      latType = Planetographic;
    else {
      QString msg = "Latitude type [" + mapping["LatitudeType"][0] + "] is not recognized";
      throw IException(IException::Programmer, msg, _FILEINFO_);
    }
 
    return add(angleToAdd, equatorialRadius, polarRadius, latType);
  }
 
 
  Latitude Latitude::add(Angle angleToAdd, 
                         Distance equatorialRadius, 
                         Distance polarRadius,
                         CoordinateType latType) {
    Latitude result;
 
    switch (latType) {
      case Planetocentric:
        result = Latitude(planetocentric() + angleToAdd.radians(), equatorialRadius, polarRadius,
                          latType, Radians, m_errors);
        break;
 
      case Planetographic:
        result = Latitude(planetographic() + angleToAdd.radians(), equatorialRadius, polarRadius,
                          latType, Radians, m_errors);
        break;
    }
 
    return result;
  }
 
 
  void Latitude::setAngle(const double &angle, 
                          const Angle::Units &units) {
    
    // Check for passing 90 degrees if that error checking is on
    if (!IsSpecial(angle) && (m_errors & AllowPastPole) != AllowPastPole) {
      Angle tmpAngle(angle, units);
      if (tmpAngle > Angle(90, Angle::Degrees) ||
          tmpAngle < Angle(-90, Angle::Degrees)) {
        QString msg = "Latitudes past 90 degrees are not valid. The latitude [" 
                      + Isis::toString(tmpAngle.degrees(), 8) + "] is not allowed";
        throw IException(IException::Programmer, msg, _FILEINFO_);
      }
    }
    Angle::setAngle(angle, units);
  }
}