import logging
import re
from math import floor

import geopandas
from pyproj import Transformer
from shapely.geometry import Point, Polygon

TRANSFORMER_OS_GRID_TO_WGS84 = Transformer.from_crs("EPSG:27700", "EPSG:4326")
TRANSFORMER_IRISH_GRID_TO_WGS84 = Transformer.from_crs("EPSG:29903", "EPSG:4326")
TRANSFORMER_CI_UTM_GRID_TO_WGS84 = Transformer.from_crs("+proj=utm +zone=30 +ellps=WGS84", "EPSG:4326")

cq_zone_data = geopandas.GeoDataFrame.from_features(geopandas.read_file("datafiles/cqzones.geojson"))
itu_zone_data = geopandas.GeoDataFrame.from_features(geopandas.read_file("datafiles/ituzones.geojson"))


# Finds out which CQ zone a lat/lon point is in.
def lat_lon_to_cq_zone(lat, lon):
    for index, row in cq_zone_data.iterrows():
        polygon = Polygon(row["geometry"])
        test_point = Point(lon, lat)
        if polygon.contains(test_point):
            return int(row["name"])
    return None


# Finds out which ITU zone a lat/lon point is in.
def lat_lon_to_itu_zone(lat, lon):
    for index, row in itu_zone_data.iterrows():
        polygon = Polygon(row["geometry"])
        test_point = Point(lon, lat)
        if polygon.contains(test_point):
            return int(row["name"])
    return None


# Convert a Maidenhead grid reference of arbitrary precision to the lat/long of the centre point of the square.
# Returns None if the grid format is invalid.
def lat_lon_for_grid_centre(grid):
    lat, lon, lat_cell_size, lon_cell_size = lat_lon_for_grid_sw_corner_plus_size(grid)
    if lat is not None and lon is not None and lat_cell_size is not None and lon_cell_size is not None:
        return [lat + lat_cell_size / 2.0, lon + lon_cell_size / 2.0]
    else:
        return None


# Convert a Maidenhead grid reference of arbitrary precision to the lat/long of the southwest corner of the square.
# Returns None if the grid format is invalid.
def lat_lon_for_grid_sw_corner(grid):
    lat, lon, lat_cell_size, lon_cell_size = lat_lon_for_grid_sw_corner_plus_size(grid)
    if lat is not None and lon is not None:
        return [lat, lon]
    else:
        return None

# Convert a Maidenhead grid reference of arbitrary precision to the lat/long of the northeast corner of the square.
# Returns None if the grid format is invalid.
def lat_lon_for_grid_ne_corner(grid):
    lat, lon, lat_cell_size, lon_cell_size = lat_lon_for_grid_sw_corner_plus_size(grid)
    if lat is not None and lon is not None and lat_cell_size is not None and lon_cell_size is not None:
        return [lat + lat_cell_size, lon + lon_cell_size]
    else:
        return None


# Convert a Maidenhead grid reference of arbitrary precision to lat/long, including in the result the size of the
# lowest grid square. This is a utility method used by the main methods that return the centre, southwest, and
# northeast coordinates of a grid square.
# The return type is always a tuple of size 4. The elements in it are None if the grid format is invalid.
def lat_lon_for_grid_sw_corner_plus_size(grid):
    # Make sure we are in upper case so our maths works. Case is arbitrary for Maidenhead references
    grid = grid.upper()

    # Return None if our Maidenhead string is invalid or too short
    length = len(grid)
    if length <= 0 or (length % 2) != 0:
        return (None, None, None, None)

    lat = 0.0  # aggregated latitude
    lon = 0.0  # aggregated longitude
    lat_cell_size = 10.0  # Size in degrees latitude of the current cell. Starts at 10 and gets smaller as the calculation progresses
    lon_cell_size = 20.0  # Size in degrees longitude of the current cell. Starts at 20 and gets smaller as the calculation progresses

    # Iterate through blocks (two-character sections)
    block = 0
    while block * 2 < length:
        if block % 2 == 0:
            # Letters in this block
            lon_cell_no = ord(grid[block * 2]) - ord('A')
            lat_cell_no = ord(grid[block * 2 + 1]) - ord('A')
            # Bail if the values aren't in range. Allowed values are A-R (0-17) for the first letter block, or
            # A-X (0-23) thereafter.
            max_cell_no = 17 if block == 0 else 23
            if lat_cell_no < 0 or lat_cell_no > max_cell_no or lon_cell_no < 0 or lon_cell_no > max_cell_no:
                return (None, None, None, None)
        else:
            # Numbers in this block
            try:
                lon_cell_no = int(grid[block * 2])
                lat_cell_no = int(grid[block * 2 + 1])
            except ValueError:
                return (None, None, None, None)
            # Bail if the values aren't in range 0-9
            if lat_cell_no < 0 or lat_cell_no > 9 or lon_cell_no < 0 or lon_cell_no > 9:
                return (None, None, None, None)

        # Aggregate the angles
        lat += lat_cell_no * lat_cell_size
        lon += lon_cell_no * lon_cell_size

        # Reduce the cell size for the next block, unless we are on the last cell.
        if block * 2 < length - 2:
            # Still have more work to do, so reduce the cell size
            if block % 2 == 0:
                # Just dealt with letters, next block will be numbers so cells will be 1/10 the current size
                lat_cell_size = lat_cell_size / 10.0
                lon_cell_size = lon_cell_size / 10.0
            else:
                # Just dealt with numbers, next block will be letters so cells will be 1/24 the current size
                lat_cell_size = lat_cell_size / 24.0
                lon_cell_size = lon_cell_size / 24.0

        block += 1

    # Offset back to (-180, -90) where the grid starts
    lon -= 180.0
    lat -= 90.0

    # Return None values on maths errors
    if any(x != x for x in [lat, lon, lat_cell_size, lon_cell_size]):  # NaN check
        return None, None, None, None

    return lat, lon, lat_cell_size, lon_cell_size


# Convert a Worked All Britain or Worked All Ireland reference to a lat/lon point.
def wab_wai_square_to_lat_lon(ref):
    # First check we have a valid grid square, and based on what it looks like, use either the Ordnance Survey, Irish,
    # or UTM grid systems to perform the conversion.
    if re.match(r"^[HNOST][ABCDEFGHJKLMNOPQRSTUVWXYZ][0-9]{2}$", ref):
        return os_grid_square_to_lat_lon(ref)
    elif re.match(r"^[ABCDEFGHJKLMNOPQRSTUVWXYZ][0-9]{2}$", ref):
        return irish_grid_square_to_lat_lon(ref)
    elif re.match(r"^W[AV][0-9]{2}$", ref):
        return utm_grid_square_to_lat_lon(ref)
    else:
        logging.warning("Invalid WAB/WAI square: " + ref)
        return None


# Get a lat/lon point for the centre of an Ordnance Survey grid square
def os_grid_square_to_lat_lon(ref):
    # Convert the letters into multipliers for the 500km squares and 100km squares
    offset_500km_multiplier = ord(ref[0]) - 65
    offset_100km_multiplier = ord(ref[1]) - 65

    # The letter "I" is not used in the grid, so any offset of 8 or more needs to be reduced by 1.
    if offset_500km_multiplier >= 8:
        offset_500km_multiplier = offset_500km_multiplier - 1
    if offset_100km_multiplier >= 8:
        offset_100km_multiplier = offset_100km_multiplier - 1

    # Convert the offsets into increments of 100km from the false origin (grid square SV):
    easting_100km = ((offset_500km_multiplier - 2) % 5) * 5 + (offset_100km_multiplier % 5)
    northing_100km = (19 - floor(offset_500km_multiplier / 5) * 5) - floor(offset_100km_multiplier / 5)

    # Take the numeric parts of the grid square and multiply by 10000 to get metres, then combine with the 100km
    # box offsets
    easting = int(ref[2]) * 10000 + easting_100km * 100000
    northing = int(ref[3]) * 10000 + northing_100km * 100000

    # Add 5000m to each value to get the middle of the box rather than the south-west corner
    easting = easting + 5000
    northing = northing + 5000

    # Reproject to WGS84 lat/lon
    lat, lon = TRANSFORMER_OS_GRID_TO_WGS84.transform(easting, northing)
    return lat, lon


# Get a lat/lon point for the centre of an Irish Grid square.
def irish_grid_square_to_lat_lon(ref):
    # Convert the letters into multipliers for the 100km squares
    offset_100km_multiplier = ord(ref[0]) - 65

    # The letter "I" is not used in the grid, so any offset of 8 or more needs to be reduced by 1.
    if offset_100km_multiplier >= 8:
        offset_100km_multiplier = offset_100km_multiplier - 1

    # Convert the offsets into increments of 100km from the false origin:
    easting_100km = offset_100km_multiplier % 5
    northing_100km = 4 - floor(offset_100km_multiplier / 5)

    # Take the numeric parts of the grid square and multiply by 10000 to get metres, then combine with the 100km
    # box offsets
    easting = int(ref[1]) * 10000 + easting_100km * 100000
    northing = int(ref[2]) * 10000 + northing_100km * 100000

    # Add 5000m to each value to get the middle of the box rather than the south-west corner
    easting = easting + 5000
    northing = northing + 5000

    # Reproject to WGS84 lat/lon
    lat, lon = TRANSFORMER_IRISH_GRID_TO_WGS84.transform(easting, northing)
    return lat, lon


# Get a lat/lon point for the centre of a UTM grid square (supports only squares WA & WV for the Channel Islands, nothing else implemented)
def utm_grid_square_to_lat_lon(ref):
    # Take the numeric parts of the grid square and multiply by 10000 to get metres from the corner of the letter-based grid square
    easting = int(ref[2]) * 10000
    northing = int(ref[3]) * 10000

    # Apply the appropriate offset based on whether the square is WA or WV
    easting = easting + 500000
    if ref[1] == "A":
        northing = northing + 5500000
    else:
        northing = northing + 5400000

    # Add 5000m to each value to get the middle of the box rather than the south-west corner
    easting = easting + 5000
    northing = northing + 5000

    # Reproject to WGS84 lat/lon
    lat, lon = TRANSFORMER_CI_UTM_GRID_TO_WGS84.transform(easting, northing)
    return lat, lon