module;

#include <boost/geometry/geometries/geometries.hpp>
#include <cassert>

export module routemon:geo.utm;

import std;
import :geo;
import :geo.wgs84;

namespace routemon::geo::utm {

class zonable_wgs84_point
{
  wgs84::normalized_point p_;

  explicit inline zonable_wgs84_point(wgs84::normalized_point p) : p_{p} {}

public:
  inline zonable_wgs84_point() : p_{} {}

  inline static auto from(wgs84::normalized_point p)
      -> std::optional<zonable_wgs84_point>
  {
    if (p.lat() < -80 || p.lat() > 84)
      return std::nullopt;
    return zonable_wgs84_point{p};
  }

  inline auto lat() const -> double { return p_.lat(); }

  inline auto lon() const -> double { return p_.lon(); }

  inline auto lon(double lon) -> void { p_.lon(lon); }

  inline auto lat(double lat) -> void
  {
    if (lat < -80 || lat > 84)
      throw std::range_error{"latitude out of range for UTM"};
    p_.lat(lat);
  }
};

} // namespace routemon::geo::utm

namespace boost::geometry::traits {

namespace {

using zonable_wgs84_point = routemon::geo::utm::zonable_wgs84_point;

} // namespace

template <>
struct tag<zonable_wgs84_point>
{
  using type = point_tag;
};

template <>
struct dimension<zonable_wgs84_point> : boost::mpl::int_<2>
{
};

template <>
struct coordinate_type<zonable_wgs84_point>
{
  using type = double;
};

template <>
struct coordinate_system<zonable_wgs84_point>
{
  using type = routemon::geo::wgs84::cs;
};

template <std::size_t index>
struct access<zonable_wgs84_point, index>
{
  static inline auto get(zonable_wgs84_point const& p) -> double
  {
    if constexpr (index == 0)
      return p.lon();
    else if constexpr (index == 1)
      return p.lat();
    else
      static_assert(false, "Out of range");
  }

  static inline auto set(zonable_wgs84_point& p, double v) -> void
  {
    if constexpr (index == 0)
      p.lon(v);
    else if constexpr (index == 1)
      p.lat(v);
    else
      static_assert(false, "Out of range");
  }
};

} // namespace boost::geometry::traits

namespace routemon::geo::utm {

using zonable_wgs84_linestring =
    bgeo::model::linestring<zonable_wgs84_point, std::vector>;

// In the sense of the common WGS84 subdivisions by simple northing
// and easting (so no Norway/Svalbard exceptions).
class zone
{
  // Note: calculations heavily depend on the values of the variants.
  enum class hemisphere : std::uint8_t
  {
    northern = 0,
    southern = 1,
  };

  // Negative if in the southern hemisphere
  // Equal to (zone_no - 1) * 2 + hemisphere
  // In range [0, 119]
  std::uint8_t zone_;

  static constexpr auto zone_no_valid(std::uint8_t zone_no) -> bool
  {
    return 0 < zone_no && zone_no <= 60;
  }

  static constexpr auto from(std::uint8_t zone_no, hemisphere h) -> std::uint8_t
  {
    if (!zone_no_valid(zone_no))
      throw std::invalid_argument{"invalid UTM zone number"};
    return (zone_no - 1) * 2 + static_cast<std::uint8_t>(h);
  }

public:
  explicit constexpr zone(std::uint8_t zone_no, hemisphere h)
    : zone_{from(zone_no, h)}
  {
  }

  constexpr auto hemisphere() const -> enum hemisphere
  {
    return static_cast<enum hemisphere>(zone_ % 2);
  }

  // Return value in range [1, 60]
  constexpr auto zone_no() const -> std::uint8_t { return 1 + zone_ / 2; }

  // Return value in range [min(), max()]
  constexpr auto as_index() const -> std::uint8_t { return zone_; }

  static auto for_wgs84_point(zonable_wgs84_point p) noexcept -> zone
  {
    auto zone_no = static_cast<std::uint8_t>(1 + (p.lon() + 180.0) / 6.0);
    auto northern = p.lat() >= 0.0;
    return zone{
      zone_no, northern ? hemisphere::northern : hemisphere::southern
    };
  }

  constexpr auto wgs84_proj_epsg() const -> int
  {
    switch (hemisphere())
    {
    case hemisphere::northern:
      return 32600 + zone_no();
    case hemisphere::southern:
      return 32700 + zone_no();
    }
  }

  // Guarantee: min().as_index() == 0.
  static constexpr auto min() noexcept -> zone
  {
    return zone{1, hemisphere::northern};
  }

  static constexpr auto max() noexcept -> zone
  {
    return zone{60, hemisphere::southern};
  }

  constexpr auto next() const -> zone
  {
    assert(zone_ <= max().as_index());
    if (zone_ == max().as_index())
      throw std::range_error{"cannot take next of greatest UTM zone"};
    auto copy = zone{*this};
    copy.zone_++;
    return copy;
  }

  constexpr auto operator==(zone rhs) const -> bool
  {
    return zone_ == rhs.zone_;
  }
};
static_assert(zone::min().as_index() == 0);
static_assert(zone::max().as_index() == 119);

} // namespace routemon::geo::utm