#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Créé nov 2025
@auteur: David ALBERTO (www.astrolabe-science.fr) CC-BY-SA
Génère une carte du monde avec des bandes colorées pour identifier
les différents climats, tels que définis depuis l'Antiquité,
en fonction de la durée maximale du jour.
"""
import matplotlib.pyplot as plt
import numpy as np
import geopandas as gpd  # gestion de données géographiques
from shapely.ops import unary_union, Polygon # fonctions de traitement SIG
# import matplotlib as mpl # Normalisation du colormap
import matplotlib.patches as mpatches  # pour les légendes
from matplotlib.colors import ListedColormap # colormap personnalisé
from shapely.geometry import MultiPolygon # regrouper des polygones en multipolygones
# =============================================================================
# Personnalisation :
# =============================================================================

plt.rcParams["font.family"] = "Roboto Condensed"
plt.rcParams["font.size"] = 10
style_climats = dict(linewidth=0.1, edgecolor='k', alpha=0.5)
increment = 0.5  # différence de durée entre deux climats voisins (0.5 ou 1)
nmax = 12 / increment  # nombre de climats

#  ---------------------------------
eps = 23 + 26/60  # obliquité de l'écliptique
eps_rad = np.radians(eps)

# =============================================================================
# # import des contours des pays (chemin à adapter selon l'adresse des fichiers):
# =============================================================================
chemin ='world-administrative-boundaries/world-administrative-boundaries.shp'
df = gpd.read_file(chemin)
#  union de toutes les terres en une seule géométrie :
geom = df.geometry
# union des terres émergées :
landmass_gs = gpd.GeoSeries(unary_union(geom))  # géométrie de l'union
landmass = gpd.GeoDataFrame({
    'geometry': landmass_gs.geometry},
    crs = df.crs)

# =============================================================================
# Création d'un colormap qualitatif avec 24 couleurs,
# à partir des colormaps tab20 et tab10 :
# =============================================================================
colors1 = list(plt.get_cmap('tab20').colors)
colors2 = list(plt.get_cmap('tab10').colors[0:5])
colors = colors1 + colors2
custom_cmap = ListedColormap(colors)

# =============================================================================
# FONCTIONS
# =============================================================================

def lat_climat(duree):
    """
    duree : entier, durée du jour maximale, en heures.
    renvoie la latitude minimale pour la durée du jour donnée.
    """
    H0 = np.radians(duree / 2 * 15)
    lat = np.arctan2(- np.cos(H0), np.tan(eps_rad))
    lat = np.degrees(lat)
    return lat

def traceclimat(n):
    """
    n : entier (numéro du climat)
    trace la délimitation d'un climat.
    """
    if increment == 1:
        duree = 11 + n
    if increment == 0.5:
        duree = 11.5 + n/2
    latmin = lat_climat(duree)
    latmax = lat_climat(duree+increment)
    # polygones définis par les deux rubans du climat (nord et sud):
    polynord = Polygon([(-180, latmin), (180, latmin), (180, latmax), (-180,latmax)])
    polysud = Polygon([(-180, -latmax), (180, -latmax), (180, -latmin), (-180,-latmin)])
    # regroupement des deux polygones en une seule géométrie de type MultiPolygone :
    poly = MultiPolygon([polynord, polysud])
    # dataframe défini par le polygone :
    rect = gpd.GeoDataFrame({'geometry': [poly]}, crs=df.crs)
    intersec = landmass.intersection(rect)
    gdf = gpd.GeoDataFrame(geometry=intersec)
    label = f"climat {n:02d} - durée {duree} h".replace(".",",")
    gdf.plot(ax=ax, color=custom_cmap(n), **style_climats, aspect=1.5)
    legend_handle = mpatches.Patch(facecolor=custom_cmap(n),
                                   label=label, **style_climats)
    liste_leg.append(legend_handle)  # ajout des carac des légendes à list_leg

# =============================================================================
#
# définition de la figure et de l'axe :
# =============================================================================

fig, ax = plt.subplots(figsize=(9, 6), tight_layout=True)
ax.set_facecolor('lightsteelblue')
ax.set_title(
    'Sous quel climat vivez-vous ? Durée maximale de la journée, selon la latitude',
             fontsize=12)
ax.set_xticks([])
ax.set_yticks([])
ax.set_ylim(-60,86)

# =============================================================================
# # tracé des terres émergées :
# =============================================================================

df.plot(ax=ax, facecolor='white', edgecolor='k', lw=0.1)  # frontières

#  parallèles particuliers (tropiques, cercle polaire arctique):
for lat in [eps, -eps, 90-eps]:
    ax.axhline(lat, ls='--', c='gray', lw=0.2)

ax.axhline(0, ls='--', c='k', lw=0.5) # Équateur


liste_leg = []  # liste collectant les infos des légendes pour tous les climats
for nclim in range(1, int(nmax+1)):
    traceclimat(nclim)

ax.legend(handles=liste_leg, bbox_to_anchor=(1.02, 1.01), loc='upper left',
          fontsize=8)

# =============================================================================
# légendes des tropiques et cercle polaire
# =============================================================================
xmin, xmax = ax.get_xlim()

ax.text(0.98*xmin, eps, 'tropique du Cancer', c='dimgray', ha='left',
        fontsize=6)
ax.text(0.98*xmin, -eps, 'tropique du Capricorne', c='dimgray', ha='left',
        va='top', fontsize=6)
ax.text(0.98*xmin, 90-eps, 'cercle polaire\narctique', c='dimgray', ha='left',
        va='center', fontsize=6)

# =============================================================================
# # signature :
# =============================================================================
ax.text(1.01, 0.0, "D. Alberto (www.astrolabe-science.fr) CC-BY-SA",
        c='gray', ha='left', va='bottom', fontsize=6, rotation=90,
        transform=ax.transAxes)

# =============================================================================
# fichiers images :
# =============================================================================

fig.savefig('climats_carte.png', dpi=300)
fig.savefig('climats_carte.pdf',
            metadata={'Author': 'D. Alberto (www.astrolabe-science.fr)'}
            )

plt.show()
