SXXXXXXX_FlightMonitor/flightmonitor/map/map_utils.py
2025-06-03 15:31:23 +02:00

508 lines
24 KiB
Python

# FlightMonitor/map/map_utils.py
import logging
import math
from typing import Tuple, Optional, List, Set, Dict, Any
try:
import pyproj
PYPROJ_MODULE_LOCALLY_AVAILABLE = True
except ImportError:
pyproj = None # type: ignore
PYPROJ_MODULE_LOCALLY_AVAILABLE = False
logging.warning("MapUtils: 'pyproj' not found. Geographic calculations will be impaired.")
try:
import mercantile
MERCANTILE_MODULE_LOCALLY_AVAILABLE = True
except ImportError:
mercantile = None # type: ignore
MERCANTILE_MODULE_LOCALLY_AVAILABLE = False
logging.warning("MapUtils: 'mercantile' not found. Tile operations will fail.")
try:
from ..utils.logger import get_logger
logger = get_logger(__name__)
except ImportError:
logger = logging.getLogger(__name__)
if not logger.hasHandlers():
logging.basicConfig(level=logging.DEBUG)
logger.warning("MapUtils using fallback standard Python logger.")
try:
from . import map_constants
except ImportError:
logger.error("MapUtils: map_constants not found. Using mock constants.")
class MockMapConstants:
COORDINATE_DECIMAL_PLACES = 4
DMS_DEGREE_SYMBOL = " deg"
DMS_MINUTE_SYMBOL = "' "
DMS_SECOND_SYMBOL = "''"
MIN_ZOOM_LEVEL = 0
DEFAULT_MAX_ZOOM_FALLBACK = 20
map_constants = MockMapConstants() # type: ignore
class MapCalculationError(Exception):
"""Custom exception for map calculations."""
pass
def _is_valid_bbox_dict(bbox_dict: Dict[str, Any]) -> bool:
if not isinstance(bbox_dict, dict):
# logger.debug("BBox validation failed: not dict.") # Troppo verboso per un check comune
return False
required_keys = ["lat_min", "lon_min", "lat_max", "lon_max"]
if not all(key in bbox_dict for key in required_keys):
# logger.debug(f"BBox validation failed: missing keys {required_keys}.")
return False
try:
lat_min = float(bbox_dict["lat_min"])
lon_min = float(bbox_dict["lon_min"])
lat_max = float(bbox_dict["lat_max"])
lon_max = float(bbox_dict["lon_max"])
if not (-90.0 <= lat_min <= 90.0 and -90.0 <= lat_max <= 90.0 and
-180.0 <= lon_min <= 180.0 and -180.0 <= lon_max <= 180.0):
# logger.debug(f"BBox validation failed: coordinates out of range {bbox_dict}.")
return False
if lat_min >= lat_max :
# logger.debug(f"BBox validation failed: lat_min ({lat_min}) >= lat_max ({lat_max}).")
return False
if lon_min >= lon_max: # Stretta per OpenSky, ma potrebbe essere valida per mappe che wrappano
# logger.debug(f"BBox validation failed: lon_min ({lon_min}) >= lon_max ({lon_max}).")
return False
if not all(math.isfinite(coord) for coord in [lat_min, lon_min, lat_max, lon_max]):
# logger.debug(f"BBox validation failed: non-finite values {bbox_dict}.")
return False
return True
except (ValueError, TypeError):
# logger.debug(f"BBox validation failed: invalid types/values {bbox_dict}. Error: {e}")
return False
except Exception as e_val: # Dovrebbe essere raro
logger.error(f"Unexpected error validating BBox: {bbox_dict}. Error: {e_val}", exc_info=True)
return False
def get_bounding_box_from_center_size(
center_latitude_deg: float, center_longitude_deg: float, area_size_km: float
) -> Optional[Tuple[float, float, float, float]]:
if not PYPROJ_MODULE_LOCALLY_AVAILABLE or pyproj is None:
logger.error("'pyproj' required for BBox from center/size. Not available.")
return None
if not (isinstance(area_size_km, (int, float)) and area_size_km > 0):
logger.error(f"Invalid area_size_km: {area_size_km}.")
return None
if not (isinstance(center_latitude_deg, (int, float)) and -90.0 <= center_latitude_deg <= 90.0):
logger.error(f"Invalid center_latitude_deg: {center_latitude_deg}.")
return None
if not (isinstance(center_longitude_deg, (int, float)) and -180.0 <= center_longitude_deg <= 180.0):
logger.error(f"Invalid center_longitude_deg: {center_longitude_deg}.")
return None
try:
geodetic_calculator = pyproj.Geod(ellps="WGS84")
half_side_m = (area_size_km / 2.0) * 1000.0
_, north_lat, _ = geodetic_calculator.fwd(center_longitude_deg, center_latitude_deg, 0.0, half_side_m)
_, south_lat, _ = geodetic_calculator.fwd(center_longitude_deg, center_latitude_deg, 180.0, half_side_m)
east_lon, _, _ = geodetic_calculator.fwd(center_longitude_deg, center_latitude_deg, 90.0, half_side_m)
west_lon, _, _ = geodetic_calculator.fwd(center_longitude_deg, center_latitude_deg, 270.0, half_side_m)
north_lat = min(90.0, max(-90.0, north_lat))
south_lat = min(90.0, max(-90.0, south_lat))
east_lon = (east_lon + 180) % 360 - 180
west_lon = (west_lon + 180) % 360 - 180
if east_lon == -180.0 and center_longitude_deg > 0 : east_lon = 180.0
if west_lon == 180.0 and center_longitude_deg < 0 : west_lon = -180.0
if not all(math.isfinite(coord) for coord in [west_lon, south_lat, east_lon, north_lat]):
logger.error("Calculated BBox from center/size contains non-finite values.")
return None
if west_lon > east_lon and area_size_km < 20000:
logger.warning(f"Calculated BBox has west_lon ({west_lon}) > east_lon ({east_lon}) for non-global size. Swapping.")
west_lon, east_lon = east_lon, west_lon
return (west_lon, south_lat, east_lon, north_lat)
except Exception as e_bbox_calc:
logger.exception(f"Error calculating BBox from center/size: {e_bbox_calc}")
return None
def get_tile_ranges_for_bbox(
bounding_box_deg: Tuple[float, float, float, float],
zoom_level: int,
) -> Optional[Tuple[Tuple[int, int], Tuple[int, int]]]:
if not MERCANTILE_MODULE_LOCALLY_AVAILABLE or mercantile is None:
logger.error("'mercantile' required for tile ranges, but not available.")
return None
west_lon, south_lat, east_lon, north_lat = bounding_box_deg
clamped_south_lat = max(-85.05112877, min(85.05112877, south_lat))
clamped_north_lat = max(-85.05112877, min(85.05112877, north_lat))
if clamped_south_lat >= clamped_north_lat:
logger.warning(f"Clamped latitude range invalid for tile calculation: S={clamped_south_lat:.4f}, N={clamped_north_lat:.4f}.")
center_lat = (south_lat + north_lat) / 2.0
center_lon = (west_lon + east_lon) / 2.0
if west_lon > east_lon: center_lon = (west_lon + (east_lon + 360.0)) / 2.0;_ = 0;_ = (0 if center_lon <= 180 else (center_lon := center_lon - 360.0))
try:
tile = mercantile.tile(center_lon, center_lat, zoom_level)
return ((tile.x, tile.x), (tile.y, tile.y))
except Exception as e_single_tile:
logger.error(f"Could not get single tile for invalid lat range: {e_single_tile}")
return None
if not (isinstance(zoom_level, int) and 0 <= zoom_level <= 25):
logger.error(f"Invalid zoom level for tile calculation: {zoom_level}.")
return None
try:
tiles_in_bbox_generator = mercantile.tiles(west_lon, clamped_south_lat, east_lon, clamped_north_lat, zooms=[zoom_level])
list_of_tiles = list(tiles_in_bbox_generator)
if not list_of_tiles:
logger.warning(f"No tiles found by mercantile for BBox {bounding_box_deg} at zoom {zoom_level}. Attempting center tile.")
center_lon = (west_lon + east_lon) / 2.0
if west_lon > east_lon: center_lon = (west_lon + (east_lon + 360.0)) / 2.0; _ = (0 if center_lon <= 180 else (center_lon := center_lon - 360.0))
center_lat = (clamped_south_lat + clamped_north_lat) / 2.0
try:
center_tile = mercantile.tile(center_lon, center_lat, zoom_level)
min_x, max_x = center_tile.x, center_tile.x
min_y, max_y = center_tile.y, center_tile.y
except Exception as e_center_tile:
logger.error(f"Fallback failed: Error getting tile for BBox center: {e_center_tile}", exc_info=True)
return None
else:
x_coords = [tile.x for tile in list_of_tiles]
y_coords = [tile.y for tile in list_of_tiles]
min_x, max_x = min(x_coords), max(x_coords)
min_y, max_y = min(y_coords), max(y_coords)
return ((min_x, max_x), (min_y, max_y))
except Exception as e_tile_range_calc:
logger.exception(f"Error calculating tile ranges for BBox {bounding_box_deg} at zoom {zoom_level}: {e_tile_range_calc}")
return None
def calculate_meters_per_pixel(
latitude_degrees: float, zoom_level: int, tile_pixel_size: int = 256
) -> Optional[float]:
try:
if not (isinstance(latitude_degrees, (int, float)) and -90.0 <= latitude_degrees <= 90.0): return None
if not (isinstance(zoom_level, int) and 0 <= zoom_level <= 25): return None
if not (isinstance(tile_pixel_size, int) and tile_pixel_size > 0): return None
EARTH_CIRCUMFERENCE_METERS = 40075016.686
clamped_lat = max(-85.05112878, min(85.05112878, latitude_degrees))
lat_radians = math.radians(clamped_lat)
denominator = tile_pixel_size * (2**zoom_level)
if abs(denominator) < 1e-9: return None
resolution_m_px = (EARTH_CIRCUMFERENCE_METERS * math.cos(lat_radians)) / denominator
if not math.isfinite(resolution_m_px) or resolution_m_px <= 0: return None
return resolution_m_px
except Exception as e_mpp_calc:
logger.exception(f"Error calculating meters per pixel: {e_mpp_calc}")
return None
def calculate_geographic_bbox_size_km(
bounding_box_deg: Tuple[float, float, float, float],
) -> Optional[Tuple[float, float]]:
if not PYPROJ_MODULE_LOCALLY_AVAILABLE or pyproj is None:
logger.error("'pyproj' required for BBox size calculation, but not available.")
return None
west_lon, south_lat, east_lon, north_lat = bounding_box_deg
if not (isinstance(west_lon, (int, float)) and isinstance(south_lat, (int, float)) and
isinstance(east_lon, (int, float)) and isinstance(north_lat, (int, float))):
logger.warning(f"Invalid coordinate types for BBox size calculation: {bounding_box_deg}.")
return None
if not (-90.0 <= south_lat <= north_lat <= 90.0):
logger.warning(f"Invalid latitude range for BBox size calculation: S={south_lat}, N={north_lat}.")
south_lat = max(-90.0, south_lat); north_lat = min(90.0, north_lat)
if south_lat >= north_lat : return None
try:
geod = pyproj.Geod(ellps="WGS84")
avg_lon_for_height = (west_lon + east_lon) / 2.0
if west_lon > east_lon:
avg_lon_for_height = (west_lon + (east_lon + 360.0)) / 2.0
if avg_lon_for_height > 180.0: avg_lon_for_height -= 360.0
_, _, height_m = geod.inv(avg_lon_for_height, south_lat, avg_lon_for_height, north_lat)
height_m = abs(height_m)
avg_lat_for_width = (south_lat + north_lat) / 2.0
_, _, width_m = geod.inv(west_lon, avg_lat_for_width, east_lon, avg_lat_for_width)
width_m = abs(width_m)
approx_width_km = width_m / 1000.0
approx_height_km = height_m / 1000.0
if approx_width_km < 1e-6: approx_width_km = 0.001
if approx_height_km < 1e-6: approx_height_km = 0.001
if not all(math.isfinite(val) for val in [approx_width_km, approx_height_km]):
logger.warning(f"Calculated BBox size has non-finite values: W={approx_width_km} km, H={approx_height_km} km.")
return None
return (approx_width_km, approx_height_km)
except Exception as e_size_calc:
logger.exception(f"Error calculating BBox size for {bounding_box_deg}: {e_size_calc}")
return None
def get_hgt_tile_geographic_bounds(
lat_coord: int, lon_coord: int
) -> Optional[Tuple[float, float, float, float]]:
if not (isinstance(lat_coord, int) and isinstance(lon_coord, int)): return None
south_lat = float(lat_coord); west_lon = float(lon_coord)
north_lat = south_lat + 1.0; east_lon = west_lon + 1.0
south_lat = max(-90.0, min(south_lat, 89.0)); north_lat = max(-89.0, min(north_lat, 90.0))
west_lon = max(-180.0, min(west_lon, 179.0)); east_lon = max(-179.0, min(east_lon, 180.0))
if south_lat >= north_lat or west_lon >= east_lon: return None
return (west_lon, south_lat, east_lon, north_lat)
def calculate_zoom_level_for_geographic_size(
latitude_degrees: float,
geographic_span_meters: float,
target_pixel_span: int,
tile_pixel_size: int = 256
) -> Optional[int]:
if not (isinstance(latitude_degrees, (int, float)) and -90.0 <= latitude_degrees <= 90.0): return None
if not (isinstance(geographic_span_meters, (int, float)) and geographic_span_meters > 0): return None
if not (isinstance(target_pixel_span, int) and target_pixel_span > 0): return None
if not (isinstance(tile_pixel_size, int) and tile_pixel_size > 0): return None
try:
required_resolution_m_px = geographic_span_meters / target_pixel_span
if required_resolution_m_px <= 1e-9 or not math.isfinite(required_resolution_m_px):
logger.warning(f"Calculated non-positive/non-finite required resolution ({required_resolution_m_px}). Cannot determine zoom.")
return None
EARTH_CIRCUMFERENCE_METERS = 40075016.686
clamped_lat = max(-85.05112878, min(85.05112878, latitude_degrees))
lat_radians = math.radians(clamped_lat)
cos_lat = math.cos(lat_radians)
if abs(cos_lat) < 1e-9:
logger.warning(f"Latitude {latitude_degrees} too close to pole for reliable zoom calculation based on span.")
return map_constants.MIN_ZOOM_LEVEL
term_for_log = (EARTH_CIRCUMFERENCE_METERS * cos_lat) / (tile_pixel_size * required_resolution_m_px)
if term_for_log <= 1e-9:
logger.warning(f"Calculated non-positive term for log2 ({term_for_log:.2e}) in zoom calculation.")
return map_constants.MIN_ZOOM_LEVEL
precise_zoom = math.log2(term_for_log)
# MODIFIED: Usare round() per uno zoom più aderente. Il margine sarà gestito altrove.
integer_zoom = int(round(precise_zoom))
min_zoom = getattr(map_constants, "MIN_ZOOM_LEVEL", 0)
# Usa un max_zoom più realistico per i servizi di tile comuni come OSM
max_zoom = getattr(map_constants, "DEFAULT_MAX_ZOOM_FALLBACK", 19) # OSM va fino a 19-20
clamped_zoom = max(min_zoom, min(integer_zoom, max_zoom))
logger.debug(f"CalcZoomForGeoSize: precise={precise_zoom:.2f}, rounded={integer_zoom}, clamped={clamped_zoom} for span {geographic_span_meters}m in {target_pixel_span}px at lat {latitude_degrees:.2f}")
return clamped_zoom
except Exception as e_zoom_calc:
logger.exception(f"Error calculating zoom level for geographic size: {e_zoom_calc}")
return None
def deg_to_dms_string(degree_value: float, coord_type: str) -> str:
if not isinstance(degree_value, (int, float)) or not math.isfinite(degree_value):
return "N/A"
if coord_type.lower() not in ("lat", "lon"):
return "N/A (Invalid Type)"
val_for_suffix = degree_value
if coord_type.lower() == "lat":
if not (-90.0 <= degree_value <= 90.0): val_for_suffix = max(-90.0, min(90.0, degree_value))
elif coord_type.lower() == "lon":
val_for_suffix = (degree_value + 180) % 360 - 180
if val_for_suffix == -180.0 and degree_value != -180.0: val_for_suffix = 180.0
abs_deg_val = abs(degree_value)
degrees = int(abs_deg_val)
minutes_decimal = (abs_deg_val - degrees) * 60.0
minutes = int(minutes_decimal)
seconds = (minutes_decimal - minutes) * 60.0
seconds = round(seconds, 2)
if seconds >= 60.0: seconds = 0.0; minutes += 1
if minutes >= 60: minutes = 0; degrees +=1
if coord_type.lower() == "lat" and degrees > 90: degrees = 90
elif coord_type.lower() == "lon" and degrees > 180: degrees = 180
suffix = ""
if coord_type.lower() == "lat": suffix = "N" if val_for_suffix >= 0 else "S"
elif coord_type.lower() == "lon":
if abs(val_for_suffix) < 1e-9 : suffix = ""
else: suffix = "E" if val_for_suffix >= 0 else "W"
try:
dms_str = (f"{degrees}{map_constants.DMS_DEGREE_SYMBOL} "
f"{minutes:02d}{map_constants.DMS_MINUTE_SYMBOL} "
f"{seconds:05.2f}{map_constants.DMS_SECOND_SYMBOL}")
except AttributeError:
dms_str = f"{degrees}deg {minutes:02d}m {seconds:05.2f}s"
return f"{dms_str} {suffix}".strip()
def get_combined_geographic_bounds_from_tile_info_list(
tile_info_list: List[Dict],
) -> Optional[Tuple[float, float, float, float]]:
if not tile_info_list: return None
min_lon_overall, min_lat_overall, max_lon_overall, max_lat_overall = 180.0, 90.0, -180.0, -90.0
initialized = False
for tile_info in tile_info_list:
lat_c = tile_info.get("latitude_coord")
lon_c = tile_info.get("longitude_coord")
if lat_c is None or lon_c is None: continue
try:
tile_bounds_wesn = get_hgt_tile_geographic_bounds(int(lat_c), int(lon_c))
if tile_bounds_wesn:
w, s, e, n = tile_bounds_wesn
if not initialized:
min_lon_overall, min_lat_overall, max_lon_overall, max_lat_overall = w, s, e, n
initialized = True
else:
min_lon_overall = min(min_lon_overall, w)
min_lat_overall = min(min_lat_overall, s)
max_lon_overall = max(max_lon_overall, e)
max_lat_overall = max(max_lat_overall, n)
except (ValueError, TypeError): continue
except Exception: continue # Catch any other error during tile bound processing
if not initialized: return None
# Basic check for validity, though individual HGT bounds should be valid
if min_lat_overall >= max_lat_overall or min_lon_overall >= max_lon_overall : return None
return (min_lon_overall, min_lat_overall, max_lon_overall, max_lat_overall)
def calculate_geographic_bbox_from_pixel_size_and_zoom(
center_latitude_deg: float,
center_longitude_deg: float,
target_pixel_width: int,
target_pixel_height: int,
zoom_level: int,
tile_pixel_size: int = 256
) -> Optional[Tuple[float, float, float, float]]:
if not MERCANTILE_MODULE_LOCALLY_AVAILABLE or mercantile is None:
logger.error("Mercantile missing, cannot calculate BBox from pixel size.")
return None
if not (isinstance(center_latitude_deg, (int, float)) and -90.0 <= center_latitude_deg <= 90.0 and
isinstance(center_longitude_deg, (int, float)) and -180.0 <= center_longitude_deg <= 180.0 and
isinstance(target_pixel_width, int) and target_pixel_width > 0 and
isinstance(target_pixel_height, int) and target_pixel_height > 0 and
isinstance(zoom_level, int) and 0 <= zoom_level <= 25 and
isinstance(tile_pixel_size, int) and tile_pixel_size > 0):
logger.error(f"Invalid parameters for calculate_geographic_bbox_from_pixel_size_and_zoom: center=({center_latitude_deg},{center_longitude_deg}), target_px=({target_pixel_width}x{target_pixel_height}), zoom={zoom_level}")
return None
try:
center_merc_x, center_merc_y = mercantile.xy(center_longitude_deg, center_latitude_deg, truncate=False)
resolution_m_px = calculate_meters_per_pixel(center_latitude_deg, zoom_level, tile_pixel_size)
if resolution_m_px is None: return None
half_width_meters = (target_pixel_width / 2.0) * resolution_m_px
half_height_meters = (target_pixel_height / 2.0) * resolution_m_px
ul_merc_x = center_merc_x - half_width_meters
ul_merc_y = center_merc_y + half_height_meters
lr_merc_x = center_merc_x + half_width_meters
lr_merc_y = center_merc_y - half_height_meters
west_lon, north_lat = mercantile.lnglat(ul_merc_x, ul_merc_y, truncate=False)
east_lon, south_lat = mercantile.lnglat(lr_merc_x, lr_merc_y, truncate=False)
MAX_MERC_LAT = 85.05112878
north_lat = max(-MAX_MERC_LAT, min(MAX_MERC_LAT, north_lat))
south_lat = max(-MAX_MERC_LAT, min(MAX_MERC_LAT, south_lat))
west_lon = (west_lon + 180) % 360 - 180
east_lon = (east_lon + 180) % 360 - 180
if west_lon == 180 and center_longitude_deg < 0: west_lon = -180.0
if east_lon == -180 and center_longitude_deg > 0: east_lon = 180.0
if abs(west_lon - east_lon) < 1e-9 or abs(south_lat - north_lat) < 1e-9 : return None
if not all(math.isfinite(coord) for coord in [west_lon, south_lat, east_lon, north_lat]): return None
return (west_lon, south_lat, east_lon, north_lat)
except Exception as e_bbox_calc:
logger.exception(f"Error calculating geo BBox from pixel size/zoom: {e_bbox_calc}")
return None
def _pixel_to_geo(
canvas_x: int,
canvas_y: int,
current_map_geo_bounds: Optional[Tuple[float, float, float, float]],
current_stitched_map_pixel_shape: Optional[Tuple[int, int]],
) -> Tuple[Optional[float], Optional[float]]:
if not MERCANTILE_MODULE_LOCALLY_AVAILABLE or mercantile is None:
logger.error("_pixel_to_geo: Mercantile missing or not available.")
return None, None
if current_map_geo_bounds is None or len(current_map_geo_bounds) != 4:
logger.warning("_pixel_to_geo: Map geo bounds missing or invalid.")
return None, None
if current_stitched_map_pixel_shape is None or len(current_stitched_map_pixel_shape) != 2:
logger.warning("_pixel_to_geo: Map pixel shape missing or invalid.")
return None, None
img_w, img_h = current_stitched_map_pixel_shape
if img_w <= 0 or img_h <= 0:
logger.warning(f"_pixel_to_geo: Map pixel dims invalid ({img_w}x{img_h}).")
return None, None
map_west_lon, map_south_lat, map_east_lon, map_north_lat = current_map_geo_bounds
try:
map_ul_merc_x, map_ul_merc_y = mercantile.xy(map_west_lon, map_north_lat, truncate=False)
map_lr_merc_x, map_lr_merc_y = mercantile.xy(map_east_lon, map_south_lat, truncate=False)
total_map_width_merc = map_lr_merc_x - map_ul_merc_x
total_map_height_merc = map_ul_merc_y - map_lr_merc_y
if map_east_lon < map_west_lon:
map_east_lon_equiv = map_east_lon + 360.0
temp_lr_merc_x, _ = mercantile.xy(map_east_lon_equiv, map_south_lat, truncate=False)
total_map_width_merc = temp_lr_merc_x - map_ul_merc_x
if abs(total_map_width_merc) < 1e-9 or abs(total_map_height_merc) < 1e-9:
logger.warning(f"_pixel_to_geo: Map mercator dims zero/near-zero ({total_map_width_merc:.2e}x{total_map_height_merc:.2e}).")
return None, None
relative_x_on_map_pixels = canvas_x / img_w
relative_y_on_map_pixels = canvas_y / img_h
target_merc_x = map_ul_merc_x + (relative_x_on_map_pixels * total_map_width_merc)
target_merc_y = map_ul_merc_y - (relative_y_on_map_pixels * total_map_height_merc)
clicked_lon, clicked_lat = mercantile.lnglat(target_merc_x, target_merc_y, truncate=False)
MAX_MERCATOR_LATITUDE = 85.05112878
clicked_lat = max(-MAX_MERCATOR_LATITUDE, min(MAX_MERCATOR_LATITUDE, clicked_lat))
clicked_lon = (clicked_lon + 180) % 360 - 180
if clicked_lon == -180.0 and target_merc_x > 0: clicked_lon = 180.0
return clicked_lon, clicked_lat
except Exception as e_pixel_to_geo:
logger.error(f"Error converting pixel ({canvas_x}, {canvas_y}) to geo. Error: {e_pixel_to_geo}", exc_info=True)
return None, None