"""Carte continue interactive des tuiles LiDAR traitées.
Génère l'interface web de consultation (servie par webapp.py via --serve) :
- output/index.html : coquille HTML minimaliste embarquant les données
- output/assets/app.css : styles (panneaux flottants type SIG)
- output/assets/app.js : logique carte (Leaflet), couches, infos tuiles
Fonctionnalités :
- Carte continue zoom/pan (Leaflet), tuiles rotées selon la projection
Lambert 93, vignettes ↔ images pleine résolution selon le zoom.
- Panneau de couches superposables avec opacité individuelle et
réordonnancement par glisser-déposer (persisté en localStorage).
- Panneau d'infos par tuile : résolution, méthode de classification du sol
(lue depuis output/DTM/*_dtm_method.txt), dates et tailles des fichiers.
- Génération de nouvelles zones via l'API web (./run.sh --serve) :
dessin d'un rectangle → téléchargement IGN + traitement.
- Onglet Export : sélection de dalles sur la carte → mosaïque jointive
(image ou PDF multi-couches) téléchargeable pour consultation sur
téléphone (/api/export de la webapp, module export.py).
Intégration:
- Appelé automatiquement à la fin de process_all() dans pipeline.py
- Régénération solo via --rebuild-index dans cli.py
"""
import json
import logging
import re
import time
from datetime import datetime
from pathlib import Path
logger = logging.getLogger("lidar")
# Noms d'affichage (français) pour le panneau de couches.
# Clé = mot-clé dans le nom de fichier de sortie (post-basename).
VIZ_LABELS = {
'hillshade_multi': 'Hillshade multidirectionnel',
'slope': 'Pente',
'aspect': 'Aspect',
'mslrm': 'MSRM (relief multi-échelle)',
'sailore': 'SAILORE (LRM adaptatif)',
'positive_openness': 'Openness positive',
'negative_openness': 'Openness négative',
'svf': 'Sky-View Factor',
'roughness': 'Rugosité',
'wavelet': 'Ondelette',
'flow_acc': 'Accumulation d\'écoulement',
'solar': 'Éclairage solaire',
'anomaly': 'Carte d\'anomalies',
'ortho': 'Orthophoto IGN',
'topo': 'Carte topographique IGN',
}
# Légendes des couches : titre, lecture du rendu (couleurs), méthode de calcul
# et dégradé de colormap. Source unique partagée entre rendering.py (images
# unitaires) et export.py (bandeau légende des mosaïques multi-dalles) — ce
# module reste volontairement sans dépendance lourde (webapp légère : ni
# matplotlib, ni GDAL).
# 'gradient' : 9 arrêts échantillonnés sur la colormap matplotlib (plt.get_cmap
# à i/8) pour dessiner une barre de dégradé sans matplotlib ; la cohérence
# avec rendering.COLORMAPS est vérifiée par les tests (test_export.py).
# 'ticks' : libellés des extrémités du dégradé (None = pas de barre).
VIZ_LEGENDS = {
'hillshade_multi': {
'title': 'Hillshade Multidirectionnel',
'legend': 'Illumination combinée de 8 directions (échelle fixe 0–1)\nBlanc = Face éclairée | Noir = Zone d\'ombre\nCouleurs homogènes entre tuiles',
'description': 'Ombres portées révélant micro-relief (murs, fossés, terrasses)',
'cmap': 'gray',
'gradient': ('#000000', '#202020', '#404040', '#606060', '#808080',
'#a0a0a0', '#c0c0c0', '#e0e0e0', '#ffffff'),
'ticks': ('0', '1'),
},
'slope': {
'title': 'Pente (Inclinaison du terrain)',
'legend': 'Inclinaison en degrés\nÉchelle fixe 0–30° — couleurs homogènes entre tuiles\nJaune = Forte pente | Violet foncé = Terrain plat',
'description': 'Murs, talus et bords ressortent en jaune — terrain plat en sombre',
'cmap': 'inferno',
'gradient': ('#000004', '#210c4a', '#57106e', '#8a226a', '#bc3754',
'#e45a31', '#f98e09', '#f9cb35', '#fcffa4'),
'ticks': ('0°', '30°'),
},
'aspect': {
'title': 'Aspect (Direction des pentes)',
'legend': 'Direction vers laquelle le terrain descend\nCycle continu : Nord→Est→Sud→Ouest→Nord\nCouleurs perceptuellement uniformes (pas de saut de teinte)',
'description': 'Orientation des pentes — utile pour distinguer structures des formes naturelles',
'cmap': 'twilight',
'gradient': ('#e2d9e2', '#95b5c7', '#6276ba', '#592a8f', '#2f1436',
'#741e4f', '#b25652', '#cca389', '#e2d9e2'),
'ticks': ('Nord 0°', 'Nord 360°'),
},
'mslrm': {
'title': 'MSRM - Multi-Scale Relief Model (échelles adaptatives)',
'legend': 'Relief combiné multi-échelles (σ locales, échelle fixe ±3σ)\nRouge = Surélévation (mur, tumulus, levée)\nBleu = Dépression (fossé, douve)\n\nÉchelles : 2 à 200m, pondérées 5–20m\nCouleurs homogènes entre tuiles\nDétecte du micro au macro',
'description': 'Combine LRM à 5–7 échelles — détecte structures de 5m à 100m simultanément',
'cmap': 'seismic',
'gradient': ('#00004c', '#0000a6', '#0101ff', '#8181ff', '#fffdfd',
'#ff7d7d', '#fe0000', '#be0000', '#800000'),
'ticks': ('-3σ', '+3σ'),
},
'sailore': {
'title': 'SAILORE - LRM Auto-Adaptatif',
'legend': 'Relief local adaptatif (σ locales, échelle fixe ±3σ)\nRouge = Surélévation | Bleu = Dépression\nCouleurs homogènes entre tuiles\n\nNoyau adapté à la pente locale\nPlat = grand noyau (25m) | Pente = petit noyau (2m)',
'description': 'Noyau qui s\'adapte à la pente locale — terrain plat=grand noyau, pente=petit noyau',
'cmap': 'seismic',
'gradient': ('#00004c', '#0000a6', '#0101ff', '#8181ff', '#fffdfd',
'#ff7d7d', '#fe0000', '#be0000', '#800000'),
'ticks': ('-3σ', '+3σ'),
},
'positive_openness': {
'title': 'Openness Positive (ouverture vers le haut)',
'legend': 'Angle d\'ouverture vers le ciel (z-score local)\nClair = Vue dégagée vers le ciel (sommets, plateaux)\nSombre = Vue bloquée (vallées encaissées)\nÉchelle fixe ±3σ — couleurs homogènes entre tuiles',
'description': 'Ray-tracing 8 directions, multi-rayon — détecte crêtes et sommets',
'cmap': 'YlOrBr',
'gradient': ('#ffffe5', '#fff7bc', '#fee390', '#fec34f', '#fe9829',
'#eb6f14', '#cb4b02', '#983404', '#662506'),
'ticks': ('-3σ', '+3σ'),
},
'negative_openness': {
'title': 'Openness Negative (ouverture vers le bas)',
'legend': 'Angle d\'ouverture vers le bas (z-score local)\nClair = Surplomb (bords de fossé, grottes)\nSombre = Terrain plat (fonds de vallée)\nÉchelle fixe ±3σ — couleurs homogènes entre tuiles\nMeilleur détecteur de cavités et dolines',
'description': 'Ray-tracing 8 directions, multi-rayon — détecte fossés, dolines, souterrains',
'cmap': 'PuBu',
'gradient': ('#fff7fb', '#ece7f2', '#d0d1e6', '#a5bddb', '#73a9cf',
'#358fc0', '#056faf', '#04598c', '#023858'),
'ticks': ('-3σ', '+3σ'),
},
'svf': {
'title': 'Sky-View Factor (fraction de ciel visible)',
'legend': 'Proportion de ciel visible (échelle physique fixe 0–1)\nBlanc/jaune = Ciel masqué (vallée, fossé, tranchée)\nNoir = Ciel dégagé (sommet, plateau)\nCouleurs homogènes entre tuiles\nLes fossés ressortent en vif — excellent pour structures linéaires',
'description': 'Détection de micro-relief — fossés en jaune/blanc, levées en sombre',
'cmap': 'hot_r',
'gradient': ('#ffffff', '#ffff81', '#ffff03', '#ffad00', '#ff5900',
'#ff0500', '#b00000', '#5c0000', '#0b0000'),
'ticks': ('0', '1'),
},
'roughness': {
'title': 'Rugosité Multi-Échelle (3m + 15m)',
'legend': 'Irrégularité du terrain combinée fine + large\nViolet foncé = Surface lisse (route, mur, sol plat)\nJaune vif = Surface rugueuse (végétation, ruines, pierres)\nCombine rugosité fine 3m (70%) + large 15m (30%)',
'description': 'Mesure la variabilité locale — surfaces anthropiques lisses vs naturelles rugueuses',
'cmap': 'plasma',
'gradient': ('#0d0887', '#4c02a1', '#7e03a8', '#aa2395', '#cc4778',
'#e66c5c', '#f89540', '#fdc527', '#f0f921'),
'ticks': ('lisse', 'rugueux'),
},
'wavelet': {
'title': 'Ondelette Mexican Hat (CWT multi-échelle)',
'legend': 'Indice RMS multi-échelles, centré sur la médiane de la\ntuile (1 = niveau moyen, plus haut = structure)\n\nGrands volumes retirés (moyenne locale 35 m) :\nun fossé en sommet ou en pente ne ressort\npas plus qu\'un fossé à plat\n\nÉtirement quantile global figé (calibré sur un\néchantillon de tuiles) : même valeur = même couleur\nsur toutes les tuiles et résolutions\n\nOptimisé petites structures :\nchemins, fossés, ramparts',
'description': 'Transformée en ondelette 2D : détection des petites structures (chemins, fossés, ramparts)',
'cmap': 'inferno',
'gradient': ('#000004', '#210c4a', '#57106e', '#8a226a', '#bc3754',
'#e45a31', '#f98e09', '#f9cb35', '#fcffa4'),
'ticks': ('bruit', 'structure'),
},
'flow_acc': {
'title': 'Accumulation d\'Écoulement (Flow Accumulation)',
'legend': 'Log10 du nombre de cellules amont\nVert foncé = Forte accumulation (fossé, chenal, drainage)\nJaune = Faible accumulation (terrain plat)\n\nDétection fossés et linéaires hydrologiques',
'description': 'Priority-flood + D8 — détecte fossés archéologiques et drainages',
'cmap': 'YlGn',
'gradient': ('#ffffe5', '#f7fcb9', '#d9f0a3', '#acdd8e', '#77c679',
'#40aa5c', '#228343', '#006737', '#004529'),
'ticks': ('faible', 'forte'),
},
'solar': {
'title': 'Éclairage Solaire',
'legend': "Illumination solaire (azimut 90°, altitude 30°)\nClair = Face éclairée | Sombre = Zone d'ombre",
'description': "Simulation de l'éclairage solaire matinal",
'cmap': 'gray',
'gradient': ('#000000', '#202020', '#404040', '#606060', '#808080',
'#a0a0a0', '#c0c0c0', '#e0e0e0', '#ffffff'),
'ticks': ('0', '1'),
},
'anomaly': {
'title': 'Carte d\'Anomalies (détection automatique)',
'legend': 'Score d\'anomalie composite (0–1)\nRouge = Haute anomalie (structures suspectes)\nJaune = Anomalie modérée\nBlanc = Aucun signal (terrain naturel)\n\nSeuil auto: pixels > 2σ de la moyenne locale\nCombiné: MSRM, SVF, Ondelette, Openness, Rugosité',
'description': 'Détection automatique — cible à vérifier sur terrain',
'cmap': 'YlOrRd',
'gradient': ('#ffffcc', '#ffeda0', '#fed976', '#feb24c', '#fd8c3c',
'#fc4d2a', '#e2191c', '#bb0026', '#800026'),
'ticks': ('0', '1'),
},
'ortho': {
'title': 'Photographie Aérienne IGN',
'legend': 'Orthophotographie\nImage aérienne',
'description': 'Photographie aérienne IGN (Orthophoto)',
'cmap': None,
'gradient': None,
'ticks': None,
},
'topo': {
'title': 'Carte Topographique IGN',
'legend': 'Carte IGN\nPlan topographique',
'description': 'Carte topographique IGN (Plan IGN)',
'cmap': None,
'gradient': None,
'ticks': None,
},
}
# Couche activée par défaut à l'ouverture de la carte.
DEFAULT_VIZ = 'hillshade_multi'
# Restriction optionnelle du panneau de couches : None = proposer toutes
# les visualisations présentes sur disque (aucune restriction). Le sélecteur
# de génération/régénération est piloté par le même registre (voir _render_html).
PANEL_VIZ = None
# Correspondance mot-clé de fichier de sortie → nom d'étape --only du pipeline
# (les trois visualisations dont le nom de sortie diffère du nom d'étape,
# cf. _expected_output_path dans pipeline.py).
KEYWORD_TO_STEP = {
'hillshade_multi': 'hillshade',
'positive_openness': 'pos_open',
'negative_openness': 'neg_open',
}
# Correspondance inverse : nom d'étape --only → mot-clé de fichier de sortie.
STEP_TO_KEYWORD = {step: kw for kw, step in KEYWORD_TO_STEP.items()}
def step_to_keyword(step):
"""Nom d'étape du pipeline (ex: 'pos_open') → mot-clé de fichier ('positive_openness')."""
return STEP_TO_KEYWORD.get(step, step)
def cells_with_all_viz(vis_dir, viz_keys, resolutions=(0.5,)):
"""Cellules (col, row) disposant de TOUTES les visualisations demandées.
Une cellule est complète si, pour chaque résolution de `resolutions`, un
dossier de visualisations lui correspond et contient chaque mot-clé de
`viz_keys` (ex: 'aspect', 'hillshade_multi'). Sert à distinguer les tuiles
réellement terminées de celles à compléter : une tuile existante mais
incomplète (visualisation ou résolution manquante) reste à traiter.
Returns:
Ensemble des (col, row) complets.
"""
by_cell = {}
for t in scan_tiles(vis_dir):
per_res = by_cell.setdefault((t['col'], t['row']), {})
per_res.setdefault(t['resolution'], set()).update(t['viz'].keys())
return {cell for cell, per_res in by_cell.items()
if all(kw in per_res.get(res, ()) for res in resolutions for kw in viz_keys)}
# Toutes les visualisations sont découpées en sous-tuiles (quadrants 500 m)
# pour alléger la carte. Renseigner un tuple pour restreindre la découpe —
# les viz exclues retombent en repli dalle entière.
_CARTO_SUBTILED_VIZ = ()
# Encodage AVIF des sous-tuiles (benchmark sur dalles d'aspect réelles :
# q55 ≈ −45 % vs WebP q82, visuellement propre sur rampes de couleur).
# speed 0 (lent/meilleur) → 10 (rapide) ; 6+8 quasi identique en taille,
# on garde 8 (2× plus rapide) pour 2500×2500 px.
_SUBTILE_AVIF_QUALITY = 55
_SUBTILE_AVIF_SPEED = 8
# Couches à contenu photographique ou traits fins (orthophoto, carte topo) :
# q55, calé sur des rampes de couleur, bave textes et linéaires — qualité
# relevée pour ces couches uniquement.
_SUBTILE_AVIF_QUALITY_DETAIL = 75
_SUBTILE_DETAIL_VIZ = frozenset({'ortho', 'topo'})
# Vignette intermédiaire (px) : pallier entre la vignette 256 px et l'image
# pleine résolution, pour éviter d'étirer la vignette ou de décoder l'AVIF
# complet dès qu'une tuile dépasse ~300 px à l'écran.
_MID_THUMB_SIZE = 640
# Ordre préféré des couches (bas → haut de pile) et choix de la vignette de repli.
_VIZ_FALLBACK_ORDER = [
'hillshade_multi', 'svf', 'slope', 'mslrm', 'positive_openness',
'negative_openness', 'aspect', 'sailore', 'roughness',
'wavelet', 'flow_acc', 'solar', 'anomaly', 'ortho', 'topo',
]
# Regex pour parser les coordonnées tuile dans le basename LHD.
# LHD_FXX_{COL}_{ROW}_PTS_LAMB93_IGN69 (COL/ROW en km, Lambert 93)
_RE_LHD_COORDS = re.compile(r'^LHD_FXX_(\d+)_(\d+)_PTS_LAMB93')
def parse_basename_coords(name):
"""Extrait les coordonnées tuile (col, row en km) depuis un basename.
Args:
name: basename potentiel (ex: 'LHD_FXX_1000_6881_PTS_LAMB93_IGN69')
ou nom de dossier avec suffixe résolution ('..._r0p2').
Returns:
(col_km, row_km) ou None si le nom ne correspond pas au pattern LHD.
"""
m = _RE_LHD_COORDS.match(name)
if not m:
return None
return int(m.group(1)), int(m.group(2))
def _strip_res_suffix(dirname):
"""Sépare le basename de base et la résolution d'un nom de dossier de viz.
'LHD_FXX_1000_6881_PTS_LAMB93_IGN69' → (basename, 0.5)
'LHD_FXX_1000_6881_PTS_LAMB93_IGN69_r0p2' → (basename, 0.2)
Returns:
(basename_without_suffix, resolution_float) ou (dirname, 0.5) si pas de suffixe.
"""
m = re.match(r'^(.+?)_r(\d+p\d+)$', dirname)
if m:
res_str = m.group(2).replace('p', '.')
try:
return m.group(1), float(res_str)
except ValueError:
pass
return dirname, 0.5
def _res_suffix_str(resolution):
"""Suffixe de nommage d'une résolution (délègue à pipeline._res_suffix)."""
from .pipeline import LidarArchaeoPipeline
return LidarArchaeoPipeline._res_suffix(resolution)
def scan_tiles(vis_dir):
"""Scanne le dossier des visualisations pour inventorier les tuiles traitées.
Args:
vis_dir: Path vers output/visualisations/
Returns:
Liste de dictionnaires:
{basename, col, row, resolution, dir_path,
viz: {viz_key: {filename, ext}}, dir_name}
Triée par (resolution, row, col).
"""
vis_dir = Path(vis_dir)
if not vis_dir.is_dir():
return []
tiles = []
for entry in sorted(vis_dir.iterdir()):
if not entry.is_dir():
continue
coords = parse_basename_coords(entry.name)
if coords is None:
continue
col, row = coords
basename, resolution = _strip_res_suffix(entry.name)
# Liste les fichiers image de visualisation dans le dossier.
viz = {}
for f in sorted(entry.iterdir()):
if not f.is_file():
continue
# Détection extension AVIF/WebP
ext = None
low = f.name.lower()
for e in ('.avif', '.webp'):
if low.endswith(e):
ext = e.lstrip('.')
break
if ext is None:
continue
# viz_key = nom sans le préfixe basename_ ni l'extension
stem = f.name[:-len('.' + ext)]
prefix = basename + '_'
if not stem.startswith(prefix):
continue
viz_key = stem[len(prefix):]
viz[viz_key] = {'filename': f.name, 'ext': ext}
if not viz:
# Dossier vide ou sans image valide → ignoré
continue
tiles.append({
'basename': basename,
'col': col,
'row': row,
'resolution': resolution,
'dir_name': entry.name,
'dir_path': str(entry),
'viz': viz,
})
tiles.sort(key=lambda t: (t['resolution'], -t['row'], t['col']))
return tiles
def compute_bbox(tiles):
"""Calcule la bounding box (en km) couverte par les tuiles.
Returns:
Dict {min_col, max_col, min_row, max_row} ou None si aucune tuile.
"""
if not tiles:
return None
cols = [t['col'] for t in tiles]
rows = [t['row'] for t in tiles]
return {
'min_col': min(cols),
'max_col': max(cols),
'min_row': min(rows),
'max_row': max(rows),
}
def compute_zones(tiles, proximity_threshold=15):
"""Regroupe les tuiles en zones géographiques par clustering de proximité.
Les tuiles à moins de proximity_threshold km les unes des autres
sont regroupées dans la même zone. Les zones sont triées par taille
décroissante (plus grande zone en premier).
Args:
tiles: Liste de dictionnaires de tuiles.
proximity_threshold: Distance maximale en km pour regrouper deux tuiles.
Returns:
Liste de dictionnaires de zones:
{label, tiles, bbox}
"""
if not tiles:
return []
# Union-Find pour le clustering
parent = list(range(len(tiles)))
def find(x):
while parent[x] != x:
parent[x] = parent[parent[x]]
x = parent[x]
return x
def union(x, y):
px, py = find(x), find(y)
if px != py:
parent[px] = py
# Regrouper les tuiles proches
for i in range(len(tiles)):
for j in range(i + 1, len(tiles)):
dc = abs(tiles[i]['col'] - tiles[j]['col'])
dr = abs(tiles[i]['row'] - tiles[j]['row'])
if max(dc, dr) <= proximity_threshold:
union(i, j)
# Construire les zones
zone_members = {}
for i in range(len(tiles)):
root = find(i)
if root not in zone_members:
zone_members[root] = []
zone_members[root].append(tiles[i])
zones = []
for i, zone_tiles in enumerate(zone_members.values(), 1):
zone_bbox = compute_bbox(zone_tiles)
zones.append({
'label': f'Zone {i} ({len(zone_tiles)} tuiles)',
'tiles': zone_tiles,
'bbox': zone_bbox,
})
# Trier par taille décroissante
zones.sort(key=lambda z: len(z['tiles']), reverse=True)
return zones
def _approx_l93_to_wgs84(x_m, y_m):
"""Approximation affine Lambert 93 → WGS84 (fallback sans rasterio).
Origine exacte: (700000, 6600000) L93 ↔ (3.0°E, 46.5°N).
Précision de l'ordre du km — utilisée seulement si rasterio/PROJ
est indisponible (jamais le cas dans le Docker).
"""
import math
lat = 46.5 + (y_m - 6600000.0) / 111320.0
lon = 3.0 + (x_m - 700000.0) / (111320.0 * math.cos(math.radians(47.0)))
return lon, lat
def _approx_wgs84_to_l93(lon, lat):
"""Approximation affine WGS84 → Lambert 93 (inverse exacte du précédent).
Précision de l'ordre du km — repli sans rasterio ni pyproj (webapp
légère, cf. bbox_to_cells dans webapp.py).
"""
import math
y = (lat - 46.5) * 111320.0 + 6600000.0
x = (lon - 3.0) * (111320.0 * math.cos(math.radians(47.0))) + 700000.0
return x, y
def attach_gps_bounds(tiles):
"""Attache à chaque tuile ses coins GPS pour l'affichage Leaflet.
Chaque tuile 1×1 km est définie par son coin nord-ouest en km L93
(col, row) → X ∈ [col, col+1] km, Y ∈ [row-1, row] km.
(Vérifié sur les bounds des DTM : X_min = col×1000, Y_max = row×1000.)
Utilise rasterio.warp (conversion PROJ exacte) si disponible, sinon
pyproj, sinon l'approximation affine _approx_l93_to_wgs84 (webapp
légère sans GDAL — précision ~km).
Ajoute à chaque tuile:
corners: [[lat, lon] × 4] dans l'ordre SW, SE, NE, NW
bounds : [[lat_sud, lon_ouest], [lat_nord, lon_est]]
"""
# Coins SW, SE, NE, NW — Y du bord sud = (row-1)×1000, bord nord = row×1000
xs = []
ys = []
for t in tiles:
xs.extend([t['col'] * 1000, (t['col'] + 1) * 1000,
(t['col'] + 1) * 1000, t['col'] * 1000])
ys.extend([(t['row'] - 1) * 1000, (t['row'] - 1) * 1000,
t['row'] * 1000, t['row'] * 1000])
try:
from rasterio.warp import transform as warp_transform
lons, lats = warp_transform('EPSG:2154', 'EPSG:4326', xs, ys)
ok = True
except ImportError:
try:
from pyproj import Transformer
transformer = Transformer.from_crs('EPSG:2154', 'EPSG:4326',
always_xy=True)
lons, lats = transformer.transform(xs, ys)
ok = True
except ImportError:
logger.debug("Coins GPS approximatifs (rasterio et pyproj indisponibles)")
lons = None
lats = None
ok = False
except Exception as e:
logger.debug(f"Coins GPS approximatifs (rasterio indisponible: {e})")
lons = None
lats = None
ok = False
for i, t in enumerate(tiles):
if ok:
corners = [[lats[4 * i], lons[4 * i]],
[lats[4 * i + 1], lons[4 * i + 1]],
[lats[4 * i + 2], lons[4 * i + 2]],
[lats[4 * i + 3], lons[4 * i + 3]]]
else:
corners = []
for cx in (t['col'] * 1000, (t['col'] + 1) * 1000):
for cy in ((t['row'] - 1) * 1000, t['row'] * 1000):
lon, lat = _approx_l93_to_wgs84(cx, cy)
corners.append([lat, lon])
# Reordonner SW, SE, NE, NW (la boucle donne SW, NW, SE, NE)
corners = [corners[0], corners[2], corners[3], corners[1]]
t['corners'] = corners
t['bounds'] = [[min(c[0] for c in corners), min(c[1] for c in corners)],
[max(c[0] for c in corners), max(c[1] for c in corners)]]
return ok
def _mtime(path):
"""Mtime d'un fichier, ou None si inaccessible."""
try:
return Path(path).stat().st_mtime
except OSError:
return None
def _cached_file_fresh(path, src_mtime):
"""True si un fichier cache existe et est plus récent que sa source.
Sert à invalider vignettes et sous-tuiles quand une tuile est recalculée :
l'image source (AVIF/WebP) étant réécrite, sa mtime devient plus récente
que celle du cache, qui doit alors être régénéré.
"""
cached_mtime = _mtime(path)
if cached_mtime is None:
return False
return src_mtime is None or cached_mtime >= src_mtime
def _url_version(mtime):
"""Suffixe d'invalidation de cache pour une URL image, ou '' si inconnu.
Les images de la carte sont servies par des montages statiques sans
Cache-Control : après une régénération, le cache heuristique du navigateur
peut continuer d'afficher l'ancien rendu (même URL = image considérée
fraîche). En suffixant chaque URL de la mtime de la source (?v=ms), un
recalcul change l'URL et force le rechargement — y compris en direct
pendant un run : pollLiveTiles fusionne les nouvelles données et
updateImages() voit un src différent.
"""
return f"?v={int(mtime * 1000)}" if mtime is not None else ""
def generate_thumbnail(src_path, thumb_path, max_size=256, mid_path=None,
mid_size=640):
"""Génère une vignette JPEG depuis une image AVIF/WebP existante.
Args:
src_path: chemin de l'image source (AVIF/WebP).
thumb_path: chemin de sortie JPEG.
max_size: taille maximale (côté le plus grand) en pixels.
mid_path: vignette intermédiaire optionnelle (JPEG, taille mid_size).
mid_size: taille maximale de la vignette intermédiaire.
Returns:
True si la vignette principale est OK, False si échec.
"""
try:
from PIL import Image as PILImage
except ImportError:
logger.warning("PIL indisponible — impossible de générer les vignettes")
return False
try:
try:
resample = PILImage.Resampling.LANCZOS
except AttributeError:
resample = getattr(PILImage, 'LANCZOS', 1)
def resized(source, target):
scale = min(1.0, target / max(source.size))
if scale >= 1.0:
return source
return source.resize((max(1, int(source.size[0] * scale)),
max(1, int(source.size[1] * scale))), resample)
img = PILImage.open(str(src_path))
img = img.convert('RGB')
Path(thumb_path).parent.mkdir(parents=True, exist_ok=True)
if mid_path is not None:
try:
resized(img, mid_size).save(str(mid_path), format='JPEG', quality=82)
except Exception as e:
logger.debug(f"Vignette intermédiaire ignorée {src_path}: {e}")
resized(img, max_size).save(str(thumb_path), format='JPEG', quality=80)
return True
except Exception as e:
logger.debug(f"Vignette ignorée {src_path}: {e}")
return False
def _pick_display_viz(viz_keys):
"""Choisit la visualisation par défaut pour une tuile.
Privilégie hillshade_multi, sinon la première disponible selon l'ordre de repli.
"""
for v in _VIZ_FALLBACK_ORDER:
if v in viz_keys:
return v
return sorted(viz_keys)[0]
def _subdivision_k(resolution, tile_m=1000, target_px=2500):
"""Facteur de découpage k (grille k×k) pour alléger le rendu carte.
Au 0,2 m/px une dalle de 1 km fait 5000×5000 px (~100 Mo décodés dans le
navigateur) : on la découpe en quadrants de 500 m (k=2, 2500×2500 px).
À 0,5 m/px (2000 px) la dalle reste entière (k=1).
"""
px = max(1, int(round(tile_m / resolution)))
return max(1, int(round(px / target_px)))
def _subtile_corners(corners, i, j, k):
"""Coins WGS84 [SW, SE, NE, NW] de la sous-tuile (i, j) d'un découpage k×k.
i : indice vers l'est (0..k-1), j : indice vers le nord (0..k-1).
Interpolation bilinéaire des coins de la dalle — le quadrilatère projeté
est quasi un parallélogramme à cette échelle (erreur écran < 1 px).
"""
sw, se, ne, nw = corners
def lerp(p, q, u):
return [p[0] + (q[0] - p[0]) * u, p[1] + (q[1] - p[1]) * u]
def at(u, v):
return lerp(lerp(sw, se, u), lerp(nw, ne, u), v)
u0, u1 = i / k, (i + 1) / k
v0, v1 = j / k, (j + 1) / k
return [at(u0, v0), at(u1, v0), at(u1, v1), at(u0, v1)]
def _fallback_full_dalle(entries, viz_key, info):
"""Repli dalle entière pour une couche non découpable : elle reste
superposable en couche (images plus lourdes, mais fonctionnelles)."""
for entry in entries.values():
entry['viz'][viz_key] = dict(info)
def _build_subtiles(tile, offered_viz_keys, output_dir, sub_dir_name):
"""Découpe une dalle en sous-tuiles (crops AVIF) pour la carte interactive.
Ne découpe que les visualisations proposées dans le panneau de couches. Retourne
la liste des entrées display (une par sous-tuile), ou None si le découpage
n'est pas nécessaire/possible (la dalle entière sera alors affichée).
"""
k = _subdivision_k(tile['resolution'])
if k <= 1:
return None
try:
from PIL import Image as PILImage
except ImportError:
return None
sub_dir = output_dir / sub_dir_name
sub_dir.mkdir(parents=True, exist_ok=True)
entries = {}
for j in range(k):
for i in range(k):
corners = _subtile_corners(tile['corners'], i, j, k)
entries[(i, j)] = {
'col': tile['col'], 'row': tile['row'],
'name': tile['name'], 'dir_name': tile['dir_name'],
'resolution': tile['resolution'],
'bounds': [[min(c[0] for c in corners), min(c[1] for c in corners)],
[max(c[0] for c in corners), max(c[1] for c in corners)]],
'corners': corners,
'display_viz': tile['display_viz'],
'viz': {},
'meta': tile.get('meta'),
'sub_i': i, 'sub_j': j, 'sub_k': k,
'size_km': round(1.0 / k, 3),
}
try:
resample = getattr(PILImage, 'LANCZOS', 1)
for viz_key in offered_viz_keys:
info = tile['viz'].get(viz_key)
if not info:
continue
avif_quality = (_SUBTILE_AVIF_QUALITY_DETAIL
if viz_key in _SUBTILE_DETAIL_VIZ
else _SUBTILE_AVIF_QUALITY)
stems = {key: f"{tile['dir_name']}_{viz_key}_{key[0]}_{key[1]}"
for key in entries}
# Régénère si au moins un fichier manque ou est périmé (dalle
# source recalculée depuis — comme les vignettes). L'URL pleine
# porte un suffixe ?v= d'invalidation : on l'ôte pour le chemin.
src = output_dir / info['full'].split('?')[0]
src_mtime = _mtime(src)
all_exist = all(
_cached_file_fresh(output_dir / sub_dir_name / (stem + '.avif'), src_mtime)
and _cached_file_fresh(output_dir / sub_dir_name / (stem + '_thumb.webp'), src_mtime)
and _cached_file_fresh(output_dir / sub_dir_name / (stem + '_mid.webp'), src_mtime)
for stem in stems.values())
if not all_exist:
logger.info(f" Sous-tuiles recalculées : {tile['dir_name']}/{viz_key} "
f"({len(stems)} découpages)")
try:
img = PILImage.open(str(src))
img.load()
except Exception as e:
logger.debug(f"Sous-tuilage {viz_key} impossible ({src.name}): {e}")
_fallback_full_dalle(entries, viz_key, info)
continue
if img.mode not in ('RGB', 'L'):
img = img.convert('RGB')
W, H = img.size
cut_failed = False
for (i, j), stem in stems.items():
# Image : ligne 0 = nord → la sous-tuile j (nord) est en haut
left, right = round(W * i / k), round(W * (i + 1) / k)
top = round(H * (1 - (j + 1) / k))
bottom = round(H * (1 - j / k))
quad = img.crop((left, top, right, bottom))
try:
quad.save(str(output_dir / sub_dir_name / (stem + '.avif')),
format='AVIF', quality=avif_quality,
speed=_SUBTILE_AVIF_SPEED)
except Exception as e:
logger.warning(f"Encodage AVIF impossible ({stem}), "
f"repli dalle entière pour {viz_key} : {e}")
cut_failed = True
break
# Supprime l'ancien crop .webp d'une génération précédente
(output_dir / sub_dir_name / (stem + '.webp')).unlink(missing_ok=True)
mid_scale = min(1.0, _MID_THUMB_SIZE / max(quad.size))
mid_img = quad
if mid_scale < 1.0:
mid_img = quad.resize((max(1, int(quad.size[0] * mid_scale)),
max(1, int(quad.size[1] * mid_scale))), resample)
mid_img.save(str(output_dir / sub_dir_name / (stem + '_mid.webp')),
format='WEBP', quality=82)
scale = min(1.0, 256 / max(quad.size))
if scale < 1.0:
quad = quad.resize((max(1, int(quad.size[0] * scale)),
max(1, int(quad.size[1] * scale))), resample)
quad.save(str(output_dir / sub_dir_name / (stem + '_thumb.webp')),
format='WEBP', quality=80)
if cut_failed:
_fallback_full_dalle(entries, viz_key, info)
continue
v = _url_version(src_mtime)
for (i, j), stem in stems.items():
entries[(i, j)]['viz'][viz_key] = {
'thumb': f"{sub_dir_name}/{stem}_thumb.webp{v}",
'mid': f"{sub_dir_name}/{stem}_mid.webp{v}",
'full': f"{sub_dir_name}/{stem}.avif{v}",
}
except Exception as e:
logger.warning(f"Sous-tuilage abandonné pour {tile['dir_name']}: {e}")
return None
usable = [e for e in entries.values() if e['viz']]
if not usable:
return None
# Repli dalle entière pour les visualisations hors sélection : elles
# restent superposables en couche (images plus lourdes, mais fonctionnelles).
sub_keys = set(offered_viz_keys)
for viz_key, info in tile['viz'].items():
if viz_key not in sub_keys:
_fallback_full_dalle(entries, viz_key, info)
return usable
def _viz_src_dir(tile, info):
"""Dossier réel du fichier d'une visualisation.
Les couches fusionnées d'une autre résolution (run interrompu entre les
deux passes) vivent dans leur dossier d'origine — info['dir_name'] — et
non dans dir_path de la dalle affichée.
"""
d = info.get('dir_name') or tile.get('dir_name')
if d and d != Path(tile['dir_path']).name:
return Path(tile['dir_path']).parent / d
return Path(tile['dir_path'])
def _collect_tile_metadata(tile, dtm_dir):
"""Rassemble les métadonnées de génération d'une tuile.
Lit la méthode de classification du sol depuis le sidecar DTM
(output/DTM/{basename}_dtm{suffix}_method.txt, écrit par pipeline.py),
et les dates/tailles des fichiers de visualisation.
Returns:
{method: str|None, generated: str|None,
viz: {viz_key: {date: str, size: int}}}
"""
meta = {'method': None, 'generated': None, 'viz': {}}
suffix = _res_suffix_str(tile['resolution'])
method_file = Path(dtm_dir) / f"{tile['basename']}_dtm{suffix}_method.txt"
if not method_file.exists() and suffix:
# La classification du sol est partagée entre résolutions : repli sur
# le sidecar de la résolution primaire si le spécifique manque.
method_file = Path(dtm_dir) / f"{tile['basename']}_dtm_method.txt"
try:
if method_file.exists():
method = method_file.read_text(encoding='utf-8').strip()
if method:
meta['method'] = method
# Le sidecar est écrit juste après la création du DTM :
# sa date ≈ date de génération de la dalle.
meta['generated'] = datetime.fromtimestamp(
method_file.stat().st_mtime).strftime('%Y-%m-%d %H:%M')
except OSError as e:
logger.debug(f"Métadonnées illisibles {method_file.name}: {e}")
for viz_key, info in tile['viz'].items():
try:
viz_dir = _viz_src_dir(tile, info)
st = (viz_dir / info['filename']).stat()
meta['viz'][viz_key] = {
'date': datetime.fromtimestamp(st.st_mtime).strftime('%Y-%m-%d %H:%M'),
'size': st.st_size,
}
except OSError:
continue
if meta['generated'] is None and meta['viz']:
dates = [v['date'] for v in meta['viz'].values()]
meta['generated'] = min(dates)
return meta
def build_index(output_dir, output_format='avif'):
"""Génère la carte interactive HTML des tuiles traitées, organisée par zones.
Scanne output_dir/visualisations/, collecte les métadonnées de génération,
génère les vignettes JPEG, écrit les assets (CSS/JS) puis
output_dir/index.html.
Args:
output_dir: dossier de sortie racine (contient visualisations/).
output_format: format des images ('avif' ou 'webp') — pour info.
Returns:
Path vers index.html si succès, None si échec ou aucune tuile.
"""
output_dir = Path(output_dir)
vis_dir = output_dir / 'visualisations'
dtm_dir = output_dir / 'DTM'
t_start = time.time()
tiles = scan_tiles(vis_dir)
if not tiles:
logger.info("Aucune tuile traitée trouvée — index global non généré")
return None
# Bounds GPS par tuile (géoréférencement exact pour la carte Leaflet)
attach_gps_bounds(tiles)
# Une seule tuile par position (col, row) : on garde la résolution la plus
# fine disponible. Sinon les versions 0,5 m et 0,2 m d'une même dalle se
# superposent exactement sur la carte et celle ajoutée en dernier dans le
# DOM (la moins résolue, tri croissant) masque l'autre.
best_by_pos = {}
tiles_by_pos = {}
for t in tiles:
key = (t['col'], t['row'])
tiles_by_pos.setdefault(key, []).append(t)
if key not in best_by_pos or t['resolution'] < best_by_pos[key]['resolution']:
best_by_pos[key] = t
# Complète chaque dalle affichée (résolution la plus fine) avec les
# visualisations produites uniquement à l'autre résolution : sans cela,
# une couche en cours de génération (passe 0,5 m faite, 0,2 m pas encore)
# reste invisible sur la carte et absente du panneau de couches.
# dir_name mémorise le dossier d'origine : vignettes, URLs et métadonnées
# doivent lire le fichier là où il existe réellement.
for key, best in best_by_pos.items():
for other in tiles_by_pos[key]:
if other is best:
continue
for viz_key, viz_info in other['viz'].items():
if viz_key not in best['viz']:
best['viz'][viz_key] = dict(viz_info, dir_name=other['dir_name'])
tiles = sorted(best_by_pos.values(),
key=lambda t: (t['resolution'], -t['row'], t['col']))
# Détecter les zones géographiques (pour info uniquement)
zones = compute_zones(tiles)
logger.info(f" {len(zones)} zone(s) détectée(s)")
thumb_dir = output_dir / 'index_thumbs'
thumb_dir.mkdir(parents=True, exist_ok=True)
# Collecte toutes les visualisations disponibles (pour le panneau de couches).
all_viz_keys = set()
for t in tiles:
all_viz_keys.update(t['viz'].keys())
# Restreint le découpage en sous-tuiles aux visualisations choisies
if _CARTO_SUBTILED_VIZ:
sub_viz = [v for v in _CARTO_SUBTILED_VIZ if v in all_viz_keys]
if sub_viz:
logger.info(f" Sous-tuilage limité à : {', '.join(sub_viz)}")
else:
sub_viz = list(all_viz_keys)
# Génère les vignettes et construit les données pour le HTML.
zone_records = []
thumbs_generated = 0
thumbs_failed = 0
tile_idx = 0
n_tiles = len(tiles)
logger.info(f" Vignettes : {n_tiles} tuile(s) × {len(all_viz_keys)} visualisation(s)")
for zone in zones:
zone_tile_records = []
for t in zone['tiles']:
tile_idx += 1
viz_thumbs = {}
regen = 0
for viz_key, info in t['viz'].items():
# Les couches fusionnées d'une autre résolution vivent dans
# leur dossier d'origine (info['dir_name']), pas dir_path.
src = _viz_src_dir(t, info) / info['filename']
src_mtime = _mtime(src)
thumb_name = f"{t['dir_name']}_{viz_key}.jpg"
thumb_path = thumb_dir / thumb_name
mid_name = f"{t['dir_name']}_{viz_key}_mid.jpg"
mid_path = thumb_dir / mid_name
# Régénère si manquante ou périmée (tuile recalculée depuis)
if (not _cached_file_fresh(thumb_path, src_mtime)
or not _cached_file_fresh(mid_path, src_mtime)):
if generate_thumbnail(src, thumb_path, mid_path=mid_path,
mid_size=_MID_THUMB_SIZE):
thumbs_generated += 1
regen += 1
else:
thumbs_failed += 1
continue
else:
thumbs_generated += 1
viz_dir_name = _viz_src_dir(t, info).name
v = _url_version(src_mtime)
viz_thumbs[viz_key] = {
'thumb': f"index_thumbs/{thumb_name}{v}",
'full': f"visualisations/{viz_dir_name}/{info['filename']}{v}",
}
if mid_path.is_file():
viz_thumbs[viz_key]['mid'] = f"index_thumbs/{mid_name}{v}"
if regen:
logger.info(f" [{tile_idx}/{n_tiles}] {t['dir_name']} — "
f"{regen} vignette(s) régénérée(s)")
if not viz_thumbs:
continue
display_viz = _pick_display_viz(viz_thumbs.keys())
tile_meta = _collect_tile_metadata(t, dtm_dir)
zone_tile_records.append({
'col': t['col'],
'row': t['row'],
'name': t['basename'],
'dir_name': t['dir_name'],
'resolution': t['resolution'],
'bounds': t.get('bounds'),
'corners': t.get('corners'),
'display_viz': display_viz,
'viz': viz_thumbs,
'meta': tile_meta,
})
if zone_tile_records:
zone_records.append({
'label': zone['label'],
'tiles': zone_tile_records,
'bbox': zone['bbox'],
})
if not zone_records:
logger.warning("Aucune vignette générée — index global abandonné")
return None
# Calculer le bbox global (pour la stat bar)
global_bbox = compute_bbox(tiles)
# HTML avec données intégrées : liste plate des quads affichables.
# Les dalles 0,2 m (5000×5000 px) sont découpées en sous-tuiles 500 m
# (quadrants 2500×2500 px) pour alléger mémoire et chargements navigateur.
sub_dir_name = 'index_subtiles'
display_tiles = []
n_dalles = 0
n_sous = 0
for zr in zone_records:
for t in zr['tiles']:
n_dalles += 1
subs = _build_subtiles(t, sub_viz, output_dir, sub_dir_name)
if subs:
display_tiles.extend(subs)
n_sous += len(subs)
else:
display_tiles.append(t)
if n_sous:
logger.info(f" {n_sous} sous-tuiles générée(s) pour {n_dalles} dalle(s)")
# Export des données de carte (index_tiles.json) : la webapp les sert via
# /api/tiles pour afficher les tuiles terminées en direct, sans recharger
# la page — le pipeline le réécrit après chaque tuile en mode incrémental.
ordered_viz, default_viz = _panel_selection(all_viz_keys)
(output_dir / 'index_tiles.json').write_text(json.dumps({
'tiles': display_tiles,
'viz_meta': {k: {'label': VIZ_LABELS.get(k, k)} for k in ordered_viz},
'stats': {'default_viz': default_viz, 'n_tiles': len(display_tiles)},
}, ensure_ascii=False), encoding='utf-8')
_write_assets(output_dir)
html = _render_html(display_tiles, global_bbox, all_viz_keys, output_format)
html_path = output_dir / 'index.html'
html_path.write_text(html, encoding='utf-8')
logger.info(f"Index global généré : {html_path} ({time.time() - t_start:.1f}s)")
logger.info(f" {len(tiles)} tuile(s) • {thumbs_generated} vignette(s) générée(s)"
+ (f" • {thumbs_failed} échec(s)" if thumbs_failed else ""))
logger.info(f" Grille : {global_bbox['min_col']}-{global_bbox['max_col']} km E × "
f"{global_bbox['min_row']}-{global_bbox['max_row']} km N")
return html_path
def _write_assets(output_dir):
"""Écrit les fichiers statiques de l'interface (CSS + JS) dans output/assets/.
Inclut Leaflet vendorisé (assets/vendor/leaflet) : la carte reste
fonctionnelle hors ligne, sans dépendance au CDN unpkg.
"""
import shutil
assets_dir = Path(output_dir) / 'assets'
assets_dir.mkdir(parents=True, exist_ok=True)
(assets_dir / 'app.css').write_text(_APP_CSS, encoding='utf-8')
(assets_dir / 'app.js').write_text(_APP_JS, encoding='utf-8')
src_vendor = Path(__file__).parent / 'assets' / 'vendor'
dst_vendor = assets_dir / 'vendor'
if src_vendor.is_dir():
shutil.rmtree(dst_vendor, ignore_errors=True)
shutil.copytree(src_vendor, dst_vendor)
def _panel_selection(all_viz_keys):
"""Couches proposées dans le panneau + couche activée par défaut.
Partagé par la coquille HTML et l'export index_tiles.json (la webapp doit
décrire le panneau exactement comme la page au moment du rebuild).
"""
# Panneau complet : toutes les couches du registre (VIZ_LABELS), même
# absentes du disque — activables dès leur apparition — plus les couches
# historiques éventuellement trouvées sur disque (libellé = clé brute).
# PANEL_VIZ permet toujours de restreindre volontairement la sélection.
panel_keys = set(VIZ_LABELS) | set(all_viz_keys)
if PANEL_VIZ is not None:
panel_keys &= set(PANEL_VIZ)
# Ordre des couches dans le panneau (selon ordre préféré puis alpha).
ordered_viz = [v for v in _VIZ_FALLBACK_ORDER if v in panel_keys]
for v in sorted(panel_keys):
if v not in ordered_viz:
ordered_viz.append(v)
# Couche activée par défaut : préférence globale puis aspect, à condition
# d'exister sur le disque ; sinon la 1re proposée présente, sinon la 1re.
default_viz = next(
(v for v in (DEFAULT_VIZ, 'aspect') if v in all_viz_keys),
next((v for v in ordered_viz if v in all_viz_keys),
ordered_viz[0] if ordered_viz else DEFAULT_VIZ))
return ordered_viz, default_viz
def _render_html(tiles, global_bbox, all_viz_keys, output_format):
"""Construit la coquille HTML (données intégrées, styles/JS dans assets/)."""
ordered_viz, default_viz = _panel_selection(all_viz_keys)
tiles_json = json.dumps(tiles, ensure_ascii=False)
# Options du sélecteur de génération/régénération : TOUTES les couches du
# registre (mêmes libellés que le panneau ; valeurs = noms d'étapes --only
# du pipeline), y compris celles absentes du disque — pour pouvoir les
# générer. Toutes pré-sélectionnées : une régénération refait la dalle
# complète, l'utilisateur pouvant désélectionner pour cibler.
gen_keys = [v for v in _VIZ_FALLBACK_ORDER if v in VIZ_LABELS]
gen_keys += sorted(k for k in VIZ_LABELS if k not in gen_keys)
gen_viz_options = "\n".join(
f' '
for v in gen_keys)
viz_meta_json = json.dumps({
k: {'label': VIZ_LABELS.get(k, k)}
for k in ordered_viz
}, ensure_ascii=False)
stats_json = json.dumps({'default_viz': default_viz, 'n_tiles': len(tiles)},
ensure_ascii=False)
n_tiles = len(tiles)
grid_w = global_bbox['max_col'] - global_bbox['min_col'] + 1
grid_h = global_bbox['max_row'] - global_bbox['min_row'] + 1
return _HTML_TEMPLATE.format(
tiles_json=tiles_json,
viz_meta_json=viz_meta_json,
stats_json=stats_json,
gen_viz_options=gen_viz_options,
n_tiles=n_tiles,
grid_w=grid_w,
grid_h=grid_h,
output_format=output_format.upper(),
)
_HTML_TEMPLATE = """
Carte LiDAR
Zoom
Molette : zoom · Clic-glisser : déplacer
Clic sur tuile : infos de génération