diff --git a/AGENTS.md b/AGENTS.md
index bb54896..50aa880 100644
--- a/AGENTS.md
+++ b/AGENTS.md
@@ -19,10 +19,10 @@
## Conventions
- **Un seul serveur web : `mapserve.py`** (image `lidar-maps`, port 8975 léger / 8973 worker). L'ancienne webapp (`webapp.py`, `export.py`, index.html/`_APP_JS`) a été supprimée : la génération de tuiles (portée de la webapp historique — `/api/preview`, `/api/generate`, `/api/status`, `/api/stop`, `/api/queue/clear`, `/api/cell`) vit dans `mapserve.py`, l'interface dans `mapui.py` (constantes `_MAP_HTML`/`_MAP_CSS`/`_MAP_JS`, écrites par `write_map_assets()` et bâchées dans les images). Sur l'image légère sans `LIDAR_GENERATION_URL`, `/api/status` répond `available: false` et l'interface masque les boutons.
-- **`index.py` = catalogue + registres partagés** (plus d'interface) : `VIZ_LABELS`/`VIZ_LEGENDS`, défauts d'affichage (`DEFAULT_LAYERS`/`DEFAULT_OPACITY`/`DEFAULT_BLEND`), `PANEL_VIZ`/`KEYWORD_TO_STEP`, `scan_tiles`/`cells_with_all_viz`, vignettes + sous-tuiles + inventaire `index_tiles.json` (`build_index`). L'inventaire est servi par `/api/tiles` de mapserve aux machines légères (`LIDAR_SOURCE_URL`).
+- **`index.py` = catalogue + registres partagés** (plus d'interface) : `VIZ_LABELS`/`VIZ_LEGENDS`, défauts d'affichage (`DEFAULT_VIZ`/`PRECISION_VIZ`/`VIEW_MODES`), `PANEL_VIZ`/`KEYWORD_TO_STEP`, `scan_tiles`/`cells_with_all_viz`, vignettes + sous-tuiles + inventaire `index_tiles.json` (`build_index`). L'inventaire est servi par `/api/tiles` de mapserve aux machines légères (`LIDAR_SOURCE_URL`).
- **Generation is 0.2 m only** (policy): `/api/generate` (`GENERATE_RESOLUTIONS` in `mapserve.py`), the compose `process` command and the CLI `-r` default all produce 0.2 m exclusively; 0.5 m stays available via explicit `-r 0.5`. Completeness detection (`complete_cells`) requires the viz at 0.2 m only.
- **Génération du nord au sud** : les tuiles sont traitées par ligne décroissante (row = nord en km), colonnes croissantes — `find_laz_files` (pipeline.py) pour les passes batch et `_resolve_request` (mapserve.py) pour les runs lancés depuis la carte. Les workers prennent les fichiers dans l'ordre de soumission : la carte se remplit de haut en bas pendant un run (`--file` explicite au CLI = ordre utilisateur préservé). Parallélisme de génération : `LIDAR_WORKERS` (10 dans les compose ; `auto` sinon).
-- **Sub-tuilage intégral** : `_CARTO_SUBTILED_VIZ` (vide dans `index.py`) découpe TOUTES les couches en quadrants 500 m à 0,2 m ; ortho/topo sont encodées en AVIF q75 (`_SUBTILE_DETAIL_VIZ`) contre q55 pour les rampes de couleur. Une couche qui échoue à la découpe retombe en dalle entière (`_fallback_full_dalle`) sans pénaliser les autres.
+- **Sub-tuilage intégral** : `_CARTO_SUBTILED_VIZ` (vide dans `index.py`) découpe TOUTES les couches en quadrants 500 m à 0,2 m ; toutes en AVIF 4:2:0 q75 (`_SUBTILE_AVIF_QUALITY`), encodées UNE fois depuis le raster d'origine : `tif_to_crop` appelle `index.write_subtiles` juste après la dalle (plus récentes qu'elle, `build_index` ne les réencode pas ; l'ancienne chaîne dalle q60 → sous-tuile q55 cumulait deux pertes : 18,1 dB/SSIM 0,84 contre 19,5 dB/0,93). Le 4:2:0 plafonne vers 19,5 dB sur le relief (teinte pixel par pixel moyennée par 2 × 2) : seul le 4:4:4 irait plus loin (q75 : 28 dB, ~1,6× plus lourd). Exception : les aplats de niveaux (`_SUBTILE_LOSSLESS_VIZ`, la précision) en WebP sans perte (`.webp`, accepté par `tiles.py`). **Pillow ignore `lossless=True` en AVIF** (q75 avec perte) ; en RGB aucun réglage AVIF n'est exact. Une couche qui échoue à la découpe retombe en dalle entière (`_fallback_full_dalle`) sans pénaliser les autres.
- **Tuiles XYZ réutilisables hors du projet** : `/tiles/{layer}/{z}/{x}/{y}.png` suit le schéma OpenStreetMap (256 px, EPSG:3857, PNG RGBA transparent hors emprise, CORS `*`) — JOSM, iD, QGIS, uMap et MapLibre le consomment tel quel, avec découverte via TileJSON / WMTS / `josm.imagery.xml`. Le 512 px (`@2x.webp`) est réservé à l'interface interne (moitié moins de requêtes en HTTP/1.1). Toute évolution du gabarit d'URL casse des configurations clientes : la changer demande une décision explicite.
- **Bilingual naming**: all code identifiers are English; every user-facing string, log message, argparse help, and comment is French.
- **Adding a visualization requires 4 edits**: (1) `generate_X()` in `visualizations.py`, (2) entry in `VIZ_STEPS` in `pipeline.py`, (3) entry in `COLORMAPS` in `rendering.py` (or `RGB_LEGENDS` for an RGB output), (4) entry in `VIZ_LEGENDS` in `index.py` (title/legend/description + sampled cmap gradient — single text source merged into `COLORMAPS` at import, also served in `/api/map/meta` and le TileJSON). Missing any one breaks the pipeline.
@@ -30,6 +30,7 @@
- **Return `None` on failure, never raise**: `dtm.py`, `visualizations.py`, and `ign.py` all return `None` to let the pipeline continue. Raising aborts the entire file.
- **Logger is always `logging.getLogger("lidar")`**, never `__name__`. All modules route through this single logger so worker processes can configure it.
- **Ajustement conjoint des lignes de balayage (3ᵉ passe du calage, `STRIP_ALIGN_VERSION` 3)** : chaque ligne de balayage (~150 lignes/s, découpée par `_scan_line_ids` : retour de la dent de scie de `scan_angle`) peut être décalée ET inclinée (roulis) de 1 à 15 cm — lignes en creux isolées, passe entière basculée visible en bande (mesuré sur LHD_FXX_0999_6882). `_joint_line_corrections` (`dtm.py`) recale chaque ligne (a + b·u) contre le consensus des AUTRES faisceaux de sa maille 1 m (plan local par la pente), moindres carrés tronqués à 5 MAD par `bincount`, pas amorti 0,5, arrêt sous 1 mm, recentrage (pas de dérive d'ensemble), estimation sur 1 point sur 3 ; CuPy si GPU actif, repli numpy. S'y ajoute un **profil d'étalonnage par faisceau et par degré d'angle** (`STRIP_ANGLE_BIN`, commun à toutes les lignes du faisceau, moyenne retirée, pente conservée car l'inclinaison par ligne est indéterminable quand la fauchée ne traverse la dalle qu'en partie ; classes pauvres = valeur voisine) : il retire les écarts non linéaires en travers de la fauchée (±1,5 cm mesurés au bord), sources de marches parallèles au vol au bord de fauchée. Lignes sans recouvrement : `_scan_line_corrections_beam` (contre la surface de leur propre faisceau, composante ligne à ligne seule, σ 3 lignes). Validé sur blocs jamais vus (damier 10 m). La gigue par fenêtres de temps (2ᵉ passe) n'est plus calculée quand `scan_angle` existe (redondante). Coût CPU ~26 s par dalle pour tout le calage (73 s avant). Sidecar : `lines`, `line_window`, `line_cell`, `line_model`.
+- **Comblement des vides entre points borné à l'enveloppe (`GAP_FILL_VERSION` 2, `_fill_small_gaps` dans `dtm.py`)** : à 0,2 m ~80 % des pixels n'ont aucun point. L'ancien `fillnodata` à 1 m comblait tout pixel à moins de 1 m d'un point — chaque point isolé devenait une pastille plate de 2 m (teinte `atan2(0,0)` = rose saturé dans le relief orienté) et chaque trou recevait une bande extrapolée de 1 m (liseré coloré). Désormais : fermeture morphologique des pixels mesurés (rien n'est étendu vers l'extérieur) dont le rayon suit l'espacement local des points (`GAP_RADIUS_K` 1,5 × espacement mesuré sur 5 m, paliers `GAP_RADII_M` 1/1,5/2/3 m, 1 m mini pour boucher les trous de voitures), puis îlots < `GAP_MIN_ISLAND_M2` (1 m²) retirés. Morphologie par tranches numpy (`_morph_step`, ~10× scipy) ; ~5–9 s par dalle 5000². Version dans le tag GeoTIFF `LIDAR_GAP_FILL` (3 : + fichier annexe de densité) : absent/différent ⇒ DTM régénéré (`_gap_fill_matches` dans `pipeline.py`). `rasterio.fill.fillnodata` écrit dans son entrée : toujours lui passer une copie.
- **Openness à échelle fixe** : `generate_openness` normalise par des références figées `OPENNESS_POS_REF` / `OPENNESS_NEG_REF` (degrés, médianes mesurées sur 15 dalles réparties sur le territoire) et plus par z-score de dalle — même ouverture = même couleur, mosaïque jointive.
- **Filename special-cases** in `_expected_output_path()`: `pos_open` → `positive_openness`, `neg_open` → `negative_openness`, `hillshade` → `hillshade_multi`.
- **Génération imposée** : la carte ne propose plus de réglage de classification ni de raccord — `_build_command` (`mapserve.py`) passe toujours `--ground-classification ign --ign-classes sol --edge-buffer 100` ; les champs `ground_class`/`ign_classes`/`reclassify`/`edge_buffer` éventuellement envoyés sont ignorés. Défauts CLI alignés (`ign`, `100`).
@@ -47,7 +48,8 @@
- **Default output is AVIF**, not WebP. Use `--format webp` for WebP. Quality default is 60 (visually lossless on smooth color ramps, ~÷3 vs q98).
- **Calage vertical des faisceaux de vol** : chaque tuile mélange plusieurs passes (1-2 `PointSourceId` par passe) parfois biaisées verticalement de quelques cm (±2,5 cm mesurés sur 1000_6882). `create_dtm_fast` mesure l'offset robuste de chaque faisceau (points sol, maille 1 m, surface médiane itérée 3×) et retranche les offsets ≥ 0,5 cm (`STRIP_ALIGN_THRESHOLD` dans `dtm.py`) avant rastérisation. Offsets calculés **par tuile** (ils dérivent le long d'une ligne de vol : jamais de table globale), mémoïsés par LAS sol, consignés dans `DTM/*_dtm*_stripalign.json` (sidecar de cache : absent, ou version/seuil/paramètres différents ⇒ régénération du DTM). Désactivable : `--no-strip-align`.
- **Gigue intra-faisceau (2ᵉ passe du calage)** : les lignes de balayage successives d'une MÊME passe peuvent être décalées verticalement de façon aléatoire (vibration capteur / bruit haute fréquence de trajectoire) — un offset constant par faisceau n'y suffit pas. `_strip_jitter_offsets` découpe chaque faisceau en fenêtres de temps GPS (`STRIP_JITTER_BIN` = 0,1 s, origine de temps propre à chaque faisceau), mesure l'offset robuste de chaque fenêtre contre la surface médiane des AUTRES faisceaux (maille 1 m partagée, ≥ `STRIP_JITTER_MIN_CELLS` = 40 cellules), lisse la série (médiane glissante `STRIP_JITTER_SMOOTH` = 5 fenêtres), la borne à ± `STRIP_JITTER_MAX` (10 cm) puis l'interpole au temps GPS de chaque point (`_apply_strip_jitter`) ; en recouvrement à deux faisceaux, chacun reçoit une série (chacun absorbe sa part). Requiert la dimension `gps_time` (silencieusement ignorée sinon). Sidecar version 2 (séries dans `jitter`), couverte par `--no-strip-align`.
-- **Couches produites et servies = `PANEL_VIZ` (`index.py`, aujourd'hui `('relief_oriente',)`)** : pipeline sans `--only`/`--skip` (`panel_steps()`), génération lancée depuis la carte (`_panel_viz_steps` dans `mapserve.py`) et couches servies (`tiles.available_layers` filtre : panneau, `/tiles/…`, TileJSON, WMTS, JOSM). Les autres visualisations restent calculables avec `--only`. Les tests de la carte lèvent la restriction via `_setup(..., panel=None)`.
+- **Couches produites et servies = `PANEL_VIZ` (`index.py`, aujourd'hui `('relief_oriente', 'densite_sol')`)** : pipeline sans `--only`/`--skip` (`panel_steps()`), génération lancée depuis la carte (`_panel_viz_steps` dans `mapserve.py`) et couches servies (`tiles.available_layers` filtre : panneau, `/tiles/…`, TileJSON, WMTS, JOSM). Les autres visualisations restent calculables avec `--only`. Les tests de la carte lèvent la restriction via `_setup(..., panel=None)`.
+- **Précision (`densite_sol`) et affichage de la carte** : `create_dtm_fast` écrit la densité des points sol retenus (bande de raccord comprise) dans `DTM/*_dtm*_density.tif` (mailles 1 m, moyenne 3 × 3, pts/m²) ; `generate_densite_sol` la quantifie en 16 niveaux (`density_levels` : échelle log fixe, niveau k dès 0,25 × 2^(k/2) pts/m²) et la garde à 1 m (dalle 1000², ~200 Ko en WebP sans perte, `rendering.LOSSLESS_GRAY_KEYWORDS`) ; tuiles au plus proche voisin (`tiles.NEAREST_LAYERS`) pour garder 16 gris nets aux zooms 18–19. Interface (`mapui.py`) : plus de pile de couches — une couche principale (`DEFAULT_VIZ`) et trois modes relief / précision / les deux (`VIEW_MODES`, précision en `multiply` à `DEFAULT_PRECISION_OPACITY` 0,6), touche P, lien `&M=…&P=mode:opacité`, défauts figés `{main, mode, precision_opacity, base}`.
- **Relief orienté (`relief_oriente`, seule couche de la carte)** : image RGB unique (GeoTIFF uint8 3 bandes, rendue telle quelle comme ortho/topo via `RGB_KEYWORDS` dans `rendering.py`) — clarté CIELAB = openness positive **locale** (MNT − gaussienne `RELIEF_DETREND_M` = 10 m, rayons `RELIEF_RADII_M` = 5/10/20 m, 16 directions) 65 % + ombrage 35 % ; teinte = aspect, chroma fixe (`RELIEF_CHROMA`). Échelle log fixe `RELIEF_OPEN_RANGE` (pas de statistique par dalle) et support total 40 m < bande de raccord 100 m : dalles jointives. Rapide : détendance + rayons sur grille décimée ~0,8 m (`RELIEF_GRID_M`), noyau dédié qui n'accumule que la moyenne des angles (`_mean_horizon_*` : CuPy `RawKernel` sur GPU, numba parallèle sur CPU, numpy en repli), colorisation fusionnée (numba) ou vectorisée sans trigonométrie (CuPy) via une table L* × teinte (`_relief_lut`). ~5 s de calcul hors préparation sur CPU 12 cœurs. Tout changement de constante change le rendu : régénérer les dalles (`--only relief_oriente --force`).
- **Openness sous-échantillonnée** : `generate_openness` calcule le lancé de rayons (l'étape la plus coûteuse : 532 s/tuile à 0,2 m sur CPU) sur une grille décimée par blocs (`OPENNESS_DOWNSAMPLE = 2` : max par bloc en positive, min en négative — préserve les reliefs qui bornent l'horizon) puis rééchantillonne en bilinéaire. Coût ÷ facteur³ : 532 s → 40 s (×13). Signal archéologique préservé (corr. 0,93 après lissage) ; la texture de bruit sub-métrique disparaît. `--openness-downsample 1` = pleine résolution. SVF et openness anisotrope ne sont PAS concernées.
- **Raccord des bords entre tuiles** : les rendus à grand noyau (openness/SVF : rayons 100 m ; LRM : 15 m) tronquent leur fenêtre au bord de dalle — bandes d'artefacts à chaque changement de tuile. `--edge-buffer N` (défaut 100 ; toujours appliqué par la génération depuis la carte, `EDGE_BUFFER_METERS` = 100 m dans `mapserve.py`, plus de case à cocher ; 0 = off en ligne de commande) fait rastériser le DTM sur la **dalle nominale 1 km alignée sur la grille** plus une bande de N m remplie avec les points sol des 8 LAZ voisines (`_neighbor_ground_points` dans `dtm.py` : PDAL en flux, découpe + filtre de classes IGN ; voisine absente = téléchargement automatique depuis le catalogue IGN avant le run, **isolée dans `input/edge_neighbors/`** pour ne pas gonfler le corpus des passes globales (dédupliqué sur tout le lot, `_fetch_edge_neighbors` dans `pipeline.py`) ; introuvable ou échec = bande vide). Les visualisations calculent sur l'emprise étendue puis `rendering.py` (`_core_tile_window`, via `tif_to_crop`/`tif_to_png`) recadre les sorties sur la dalle 1 km exacte lue dans le nom LHD — les AVIF restent des carrés 1 km alignés dans la mosaïque. Tampon consigné dans le tag GeoTIFF `LIDAR_EDGE_BUFFER` du DTM : changer `--edge-buffer` invalide le cache DTM automatiquement (tag absent = 0). Bandes voisines non calées par faisceaux (contexte seul, recadrée hors image finale). Coût : ~7 s de lecture par voisine + ~44 % de pixels en plus à 100 m/0,2 m. Nom hors pattern LHD : option ignorée (bornes d'en-tête, pas de recadrage).
diff --git a/docs/MAPS.md b/docs/MAPS.md
index 18f3278..b7772e9 100644
--- a/docs/MAPS.md
+++ b/docs/MAPS.md
@@ -120,34 +120,38 @@ Quand le relief orienté est affiché, une rose des vents donne la couleur de
chaque orientation de pente (même formule CIELAB que le rendu) ; la clarté
porte le relief local (clair = bosse, sombre = creux).
-## Pile de couches
+## Affichage : relief et précision
-La carte ne sert que les couches de `PANEL_VIZ` (`index.py`) : aujourd'hui le
-seul **relief orienté**, qui fusionne openness locale et orientation des pentes
-en une image. Les autres visualisations présentes sur disque ne sont ni listées
-ni servies en tuiles.
+La carte ne sert que les couches de `PANEL_VIZ` (`index.py`) : le **relief
+orienté** (couche d'affichage principal, qui fusionne openness locale et
+orientation des pentes) et la **précision** (`densite_sol` : densité des points
+sol retenus pour le MNT). Les autres visualisations présentes sur disque ne
+sont ni listées ni servies en tuiles.
-Chaque couche du panneau porte trois réglages, tous persistés et transportés
-par le lien de partage :
+Il n'y a plus de pile de couches. Le panneau propose trois modes, un clic
+chacun, ou la touche **P** pour passer au suivant :
-- **ordre** — poignée ⠿ en glisser-déposer, ou boutons ▲▼ (les seuls
- utilisables au doigt). La liste est encadrée par « Haut de pile — devant » et
- « Bas de pile — derrière ». Pendant un glisser, la destination est explicite :
- la ligne tirée s'estompe et une **barre d'insertion** marque le point de
- chute, au-dessus ou en dessous de la ligne survolée selon la moitié visée ;
- à l'arrivée, la couche déplacée clignote brièvement et est ramenée dans le
- champ de vision ;
-- **opacité** — curseur par couche ;
-- **mode de fusion** — `mix-blend-mode` CSS : Normal, Produit, Écran,
- Superposé, Doux, Lumière crue, Différence, Luminosité. La superposition ne se
- réduit donc pas à de la transparence : « Produit » pose un ombrage sur une
- rampe de couleur sans la délaver, « Superposé » creuse le contraste local,
- « Différence » fait ressortir les écarts entre deux couches.
+- **Relief** — le relief orienté seul ;
+- **Précision** — la densité seule, en 16 gris (échelle log fixe : niveau k à
+ partir de 0,25 × 2^(k/2) pts/m², noir ≤ 0,35 ou aucun point, blanc ≥ 45) ;
+ se lit comme une carte de fiabilité géométrique ;
+- **Les deux** — la précision en « produit » (`mix-blend-mode: multiply`) sur
+ le relief, opacité réglable (60 % par défaut) : les zones où le relief est
+ interpolé s'assombrissent sans masquer le relief.
+
+La légende de la précision (16 paliers, info-bulle en pts/m² sur chaque
+palier) s'affiche dès que la précision est visible. Les couches LiDAR vivent
+dans un conteneur **isolé** (`isolation: isolate`) : le produit n'agit que sur
+le relief, jamais sur le fond de carte.
+
+Le lien de partage transporte la couche principale, le mode et l'opacité :
+`#z/lat/lng&M=relief_oriente&P=both:60&B=1:85:1`. Les anciens liens de la pile
+(`&L=…`) s'ouvrent sans erreur sur l'état courant.
### Figer la configuration
-Le bouton **★ Définir par défaut** enregistre l'état courant — ordre, couches
-allumées, opacités, modes de fusion, fond de carte et fusion de la pile — dans
+Le bouton **★ Définir par défaut** enregistre l'affichage courant — couche
+principale, mode, opacité de la précision, fond de carte — dans
`output/.map-defaults.json`. Tout navigateur sans réglage local part alors de
cette configuration ; **↺ Réinitialiser** oublie l'état local et y revient.
@@ -157,14 +161,11 @@ curl -X DELETE http://localhost:8975/api/map/defaults # retour au registre
```
Sans fichier enregistré, les défauts viennent du registre du pipeline
-(`DEFAULT_LAYERS`, `DEFAULT_OPACITY`, `DEFAULT_BLEND` dans `index.py`). Les
-valeurs reçues sont filtrées : couches inconnues écartées, opacités bornées à
-0–1, modes de fusion validés.
-
-Les couches LiDAR vivent dans un conteneur **isolé** (`isolation: isolate`) :
-les fusions agissent entre elles, jamais sur le fond de carte à travers les
-zones sans donnée. La fusion de la **pile entière sur le fond** se règle à
-part, en bas du panneau.
+(`DEFAULT_VIZ`, `PRECISION_VIZ`, `DEFAULT_VIEW_MODE`,
+`DEFAULT_PRECISION_OPACITY` dans `index.py`). Les valeurs reçues sont
+filtrées : couche inconnue (ou la précision elle-même) refusée comme
+principale, mode validé, opacités bornées à 0–1. Un fichier de l'ancienne pile
+(`order`/`on`/`blend`) est ignoré, sauf le fond.
> **Licence** — LiDAR HD est diffusé sous **Licence Ouverte 2.0** :
> l'attribution IGN est obligatoire et doit rester visible chez le client.
diff --git a/lidar_pipeline/dtm.py b/lidar_pipeline/dtm.py
index 7407cb5..1e077e9 100644
--- a/lidar_pipeline/dtm.py
+++ b/lidar_pipeline/dtm.py
@@ -1446,6 +1446,200 @@ def _repair_laz_with_laspy(input_laz, output_las):
return False
+# ------------------------------------------------------------
+# Comblement des vides entre points (petits trous)
+# ------------------------------------------------------------
+# À 0,2 m, ~80 % des pixels ne reçoivent aucun point : le MNT est comblé entre
+# les points mesurés. Un comblement « à distance fixe » (tout pixel à moins de
+# 1 m d'un point) faisait grossir chaque point isolé en pastille plate de 2 m
+# et extrapolait une bande de 1 m dans chaque trou (plateaux recopiés, fausses
+# pentes : pastilles roses et liserés colorés dans le relief orienté). Le
+# comblement est désormais borné à l'ENVELOPPE des points : fermeture
+# morphologique dont le rayon suit l'espacement local des points (court en
+# zone dense, long sous couvert clairsemé), puis suppression des îlots isolés.
+GAP_FILL_TAG = "LIDAR_GAP_FILL" # tag GeoTIFF : version du comblement
+GAP_FILL_VERSION = 3 # 1 (absent) = distance fixe 1 m ; 3 = + densité
+GAP_RADIUS_K = 1.5 # rayon = K × espacement local des points
+GAP_RADII_M = (1.0, 1.5, 2.0, 3.0) # paliers de rayon (m) ; 1 m mini : trous < 2 m comblés
+GAP_DENSITY_WINDOW_M = 5.0 # fenêtre de mesure de la densité locale
+GAP_MIN_ISLAND_M2 = 1.0 # îlots de données plus petits : supprimés
+
+
+def _morph_step(mask, square, erode):
+ """Un pas 3×3 (carré ou croix) de dilatation/érosion binaire par tranches
+ numpy (~10× plus rapide que scipy.ndimage). Au bord, le voisin manquant
+ compte comme le pixel lui-même."""
+ op = np.logical_and if erode else np.logical_or
+ out = mask.copy()
+ op(out[1:], mask[:-1], out=out[1:])
+ op(out[:-1], mask[1:], out=out[:-1])
+ src = out.copy() if square else mask # carré : séparable (colonnes puis lignes)
+ op(out[:, 1:], src[:, :-1], out=out[:, 1:])
+ op(out[:, :-1], src[:, 1:], out=out[:, :-1])
+ return out
+
+
+def _morph_disk(mask, radius_px, erode=False, first_step=1):
+ """Dilatation/érosion par un disque approché : octogone, pas carrés
+ (impairs) et croix (pairs) alternés. `first_step` poursuit une chaîne."""
+ for step in range(first_step, first_step + radius_px):
+ mask = _morph_step(mask, step % 2 == 1, erode)
+ return mask
+
+
+def _gap_fill_envelope(valid, resolution):
+ """Masque des pixels à combler ou conserver : enveloppe des points mesurés.
+
+ Fermeture morphologique (dilatation puis érosion) de rayon variable : les
+ vides ENTRE points voisins sont comblés, rien n'est étendu vers l'extérieur
+ (un point isolé reste un pixel). Le rayon vaut GAP_RADIUS_K × l'espacement
+ local des points, mesuré sur GAP_DENSITY_WINDOW_M dans la zone couverte
+ (le bord d'un trou ne fait pas chuter la densité), arrondi au palier
+ GAP_RADII_M le plus proche. Les îlots de moins de GAP_MIN_ISLAND_M2 sont
+ retirés (retours épars sur l'eau, dans une cour…).
+ """
+ from scipy import ndimage as nd
+ radii_px = sorted({max(1, int(round(r / resolution))) for r in GAP_RADII_M})
+ r_max = radii_px[-1]
+ # Masque étendu en miroir : la bordure de la grille n'érode pas les
+ # données et ne comble pas les bandes vides. Une seule chaîne de
+ # dilatations sert tous les paliers (disque de rayon r = r premiers pas).
+ pad = 2 * r_max
+ dilated = {}
+ grown = np.pad(valid, pad, mode="symmetric")
+ for step in range(1, r_max + 1):
+ grown = _morph_disk(grown, 1, first_step=step)
+ if step in radii_px:
+ dilated[step] = grown
+ del grown
+
+ # Densité locale : part de pixels mesurés dans la zone couverte (points à
+ # moins du plus grand rayon), convertie en espacement moyen des points.
+ # Calcul sur une grille de ~1 m (fenêtre de 5 m : largement suffisant).
+ covered = dilated[r_max][pad:-pad, pad:-pad]
+ h, w = valid.shape
+ f = max(1, int(round(1.0 / resolution)))
+ hc, wc = -(-h // f), -(-w // f)
+
+ def _block_mean(a):
+ full = np.zeros((hc * f, wc * f), dtype=np.float32)
+ full[:h, :w] = a
+ return full.reshape(hc, f, wc, f).mean(axis=(1, 3))
+
+ win = max(3, int(round(GAP_DENSITY_WINDOW_M / (resolution * f))) | 1)
+ frac = nd.uniform_filter(_block_mean(valid), win, mode="nearest")
+ frac_cov = nd.uniform_filter(_block_mean(covered), win, mode="nearest")
+ with np.errstate(divide="ignore", invalid="ignore"):
+ # rayon voulu en pixels = K × espacement / résolution
+ wanted = GAP_RADIUS_K / np.sqrt(frac / np.maximum(frac_cov, 1e-6))
+ # Rayon voulu agrandi en bilinéaire (sinon changements de palier en
+ # créneaux de 1 m sur le bord des trous), puis palier le plus proche.
+ wanted = np.nan_to_num(np.clip(wanted, 0, 2 * r_max), nan=2 * r_max)
+ if f > 1:
+ wanted = nd.zoom(wanted.astype(np.float32), f, order=1, grid_mode=True,
+ mode="nearest")[:h, :w]
+ mids = [(a + b) / 2 for a, b in zip(radii_px, radii_px[1:])]
+ level = np.digitize(wanted, mids).astype(np.uint8)
+ del covered, frac, frac_cov, wanted
+
+ # Érosion par palier, restreinte à l'emprise de ses pixels (les grands
+ # rayons ne concernent que les zones clairsemées)
+ envelope = valid.copy()
+ for i, r in enumerate(radii_px):
+ sel = (level == i) & ~valid
+ rows = np.flatnonzero(sel.any(axis=1))
+ if rows.size == 0:
+ continue
+ cols = np.flatnonzero(sel.any(axis=0))
+ y0, y1 = rows[0], rows[-1] + 1
+ x0, x1 = cols[0], cols[-1] + 1
+ m = 2 * r # marge : l'érosion au bord du découpage reste hors emprise
+ crop = dilated[r][y0 + pad - m:y1 + pad + m, x0 + pad - m:x1 + pad + m]
+ closed = _morph_disk(crop, r, erode=True)[m:-m, m:-m]
+ envelope[y0:y1, x0:x1] |= sel[y0:y1, x0:x1] & closed
+ del dilated, level
+
+ # Îlots trop petits : retirés (8-connexité, surface en m²)
+ labels, n = nd.label(envelope, structure=np.ones((3, 3), dtype=bool))
+ if n:
+ area = np.bincount(labels.ravel()) * resolution * resolution
+ small = area < GAP_MIN_ISLAND_M2
+ small[0] = False
+ if small.any():
+ envelope &= ~small[labels]
+ return envelope
+
+
+def _fill_small_gaps(dtm, resolution):
+ """Comble les vides entre points dans l'enveloppe des données.
+
+ Returns:
+ Tuple (dtm, filled_count, removed_count) : pixels comblés et pixels
+ mesurés retirés avec les îlots isolés.
+ """
+ valid = ~np.isnan(dtm)
+ if valid.all() or not valid.any():
+ return dtm, 0, 0
+ envelope = _gap_fill_envelope(valid, resolution)
+ from rasterio.fill import fillnodata
+ # Tout pixel de l'enveloppe est à moins du plus grand rayon d'un point
+ max_px = max(1, int(round(max(GAP_RADII_M) / resolution)))
+ # fillnodata écrit dans le tableau fourni : copie
+ filled = fillnodata(dtm.copy(), mask=valid, max_search_distance=max_px)
+ to_fill = envelope & ~valid & ~np.isnan(filled)
+ removed = valid & ~envelope
+ out = dtm.copy()
+ out[to_fill] = filled[to_fill]
+ out[removed] = np.nan
+ return out, int(to_fill.sum()), int(removed.sum())
+
+
+# Densité de points sol (couche « densite_sol ») : comptage des points retenus
+# pour le MNT en mailles de DENSITY_CELL_M, moyenné sur DENSITY_SMOOTH ×
+# DENSITY_SMOOTH mailles (pts/m² sur 9 m² : moins de bruit de petits entiers).
+# Écrite à côté du DTM (*_dtm*_density.tif) : les visualisations ne reçoivent
+# que le MNT, plus les points.
+DENSITY_CELL_M = 1.0
+DENSITY_SMOOTH = 3
+
+
+def density_path(dtm_path):
+ """Fichier annexe de densité d'un DTM."""
+ dtm_path = Path(dtm_path)
+ return dtm_path.with_name(f"{dtm_path.stem}_density.tif")
+
+
+def _write_density(xs, ys, bounds, dtm_path):
+ """Écrit la densité de points sol (pts/m², GeoTIFF float32 à 1 m)."""
+ from scipy import ndimage as nd
+ min_x, min_y, max_x, max_y = bounds
+ cell = DENSITY_CELL_M
+ nx = max(1, int(np.ceil((max_x - min_x) / cell - 1e-6)))
+ ny = max(1, int(np.ceil((max_y - min_y) / cell - 1e-6)))
+ top = max_y
+ counts, _, _ = np.histogram2d(top - ys, xs - min_x, bins=(ny, nx),
+ range=[[0, ny * cell], [0, nx * cell]])
+ density = nd.uniform_filter(counts, DENSITY_SMOOTH, mode="nearest") / (cell * cell)
+ out = density_path(dtm_path)
+ with rasterio.open(
+ out, 'w', driver='GTiff', height=ny, width=nx, count=1,
+ dtype='float32', crs='EPSG:2154',
+ transform=from_bounds(min_x, top - ny * cell, min_x + nx * cell, top, nx, ny),
+ compress='deflate',
+ ) as dst:
+ dst.write(density.astype('float32'), 1)
+ return out
+
+
+def read_dtm_gap_fill(dtm_path):
+ """Version du comblement enregistrée dans un DTM (1 si absente/illisible)."""
+ try:
+ with rasterio.open(dtm_path) as src:
+ return int(src.tags().get(GAP_FILL_TAG, 1) or 1)
+ except Exception:
+ return 1
+
+
def _interpolate_holes(dtm, downsample=8):
"""Fill remaining NaN holes with a terrain-aware surface interpolation.
@@ -1834,6 +2028,13 @@ def create_dtm_fast(las_file, basename, dtm_dir, resolution, force=False,
ys = np.concatenate([ys, ny])
zs = np.concatenate([zs, nz])
+ # Densité des points sol retenus (couche densite_sol), best-effort
+ try:
+ _write_density(xs, ys, (min_x, min_y, max_x, max_y),
+ dtm_dir / f"{basename}_dtm{output_suffix}.tif")
+ except Exception as e:
+ logger.warning(f" Densité de points non écrite ({e})")
+
t_raster = time.perf_counter()
dtm = bin_mean_2d(xs, ys, zs, width, height,
(min_x, max_x), (min_y, max_y))
@@ -1861,24 +2062,18 @@ def create_dtm_fast(las_file, basename, dtm_dir, resolution, force=False,
# dense ce retour est la végétation, qui imprimerait les arbres dans
# le MNT.
- # Fill small gaps (< 1 m from data) precisely — comme avant
+ # Vides entre points comblés dans l'enveloppe des données (rayon
+ # adapté à la densité locale), îlots isolés retirés
nan_count = np.count_nonzero(np.isnan(dtm))
if nan_count > 0:
total = dtm.size
nan_pct = 100.0 * nan_count / total
logger.info(f" {nan_count:,} pixels sans données ({nan_pct:.1f}%)")
-
- max_gap_pixels = max(1, int(1.0 / resolution))
t_fill = time.perf_counter()
- from rasterio.fill import fillnodata
- valid_mask = ~np.isnan(dtm)
- dtm_filled = fillnodata(dtm, mask=valid_mask, max_search_distance=max_gap_pixels)
- small_gap_mask = np.isnan(dtm) & ~np.isnan(dtm_filled)
- filled_count = np.count_nonzero(small_gap_mask)
- if filled_count > 0:
- dtm = np.where(small_gap_mask, dtm_filled, dtm)
- logger.info(f" {filled_count:,} petits trous comblés "
- f"(< {max_gap_pixels}px, {time.perf_counter() - t_fill:.1f}s)")
+ dtm, filled_count, removed_count = _fill_small_gaps(dtm, resolution)
+ logger.info(f" {filled_count:,} pixels comblés entre points, "
+ f"{removed_count:,} pixels d'îlots isolés retirés "
+ f"({time.perf_counter() - t_fill:.1f}s)")
# Save as GeoTIFF
output_tif = dtm_dir / f"{basename}_dtm{output_suffix}.tif"
@@ -1894,6 +2089,9 @@ def create_dtm_fast(las_file, basename, dtm_dir, resolution, force=False,
compress='lzw'
) as dst:
dst.write(dtm.astype('float32'), 1)
+ # Version du comblement : un DTM d'une version antérieure est
+ # régénéré (pastilles et liserés de l'ancien comblement).
+ dst.update_tags(**{GAP_FILL_TAG: str(GAP_FILL_VERSION)})
if used_edge_buffer > 0:
# Tampon de raccord inscrit dans le fichier : changement de
# --edge-buffer ⇒ invalidation automatique du cache DTM.
diff --git a/lidar_pipeline/index.py b/lidar_pipeline/index.py
index d638637..69932fa 100644
--- a/lidar_pipeline/index.py
+++ b/lidar_pipeline/index.py
@@ -6,7 +6,7 @@ Ce module produit aujourd'hui les artefacts dont le pipeline et la carte ont
besoin :
- registres partagés : VIZ_LABELS/VIZ_LEGENDS (libellés, légendes), défauts
- d'affichage (DEFAULT_LAYERS/OPACITY/BLEND), correspondances mot-clé ↔
+ d'affichage (DEFAULT_VIZ, PRECISION_VIZ, VIEW_MODES), correspondances mot-clé ↔
étape du pipeline (KEYWORD_TO_STEP) ;
- vignettes (index_thumbs/) et sous-tuiles (index_subtiles/) servant de
paliers sources à la pyramide XYZ (tiles.py) ;
@@ -46,6 +46,7 @@ VIZ_LABELS = {
'solar': 'Éclairage solaire',
'anomaly': 'Carte d\'anomalies',
'relief_oriente': 'Relief orienté',
+ 'densite_sol': 'Précision (densité de points sol)',
'ortho': 'Orthophoto IGN',
'topo': 'Carte topographique IGN',
}
@@ -185,6 +186,15 @@ VIZ_LEGENDS = {
'gradient': None,
'ticks': None,
},
+ 'densite_sol': {
+ 'title': 'Précision géométrique (densité de points sol)',
+ 'legend': 'Points sol retenus pour le MNT, par m² (moyenne sur 3 × 3 m)\n16 gris, échelle log fixe : 2 niveaux = densité doublée\nNoir = ≤ 0,35 pt/m² ou aucun point (relief interpolé ou absent)\nBlanc = ≥ 45 pts/m²',
+ 'description': 'Où le relief est mesuré (clair) et où il est interpolé (sombre) — en « produit » sur le relief, assombrit les zones peu fiables',
+ 'cmap': 'gray',
+ 'gradient': ('#000000', '#202020', '#404040', '#606060', '#808080',
+ '#a0a0a0', '#c0c0c0', '#e0e0e0', '#ffffff'),
+ 'ticks': ('≤ 0,35 pt/m²', '≥ 45 pts/m²'),
+ },
'ortho': {
'title': 'Photographie Aérienne IGN',
'legend': 'Orthophotographie\nImage aérienne',
@@ -203,21 +213,14 @@ VIZ_LEGENDS = {
},
}
-# Couches activées par défaut à l'ouverture de la carte : seules celles-ci
-# sont allumées, toutes les autres sont désactivées (mais restent accessibles
-# via le panneau). Ordre = pile basse → haute : une couche opaque (opacité 1)
-# doit rester en bas pour ne pas masquer les autres. Le relief orienté
-# fusionne déjà openness et aspect en une seule image (plus d'empilement).
-DEFAULT_LAYERS = ('relief_oriente',)
-
-# Opacité par défaut des couches allumées (0–1) ; les absentes restent à 1.
-DEFAULT_OPACITY = {}
-
-# Mode de fusion (mix-blend-mode CSS) par couche : les couches d'ombrage en
-# niveaux de gris fondent en « multiply » — elles assombrissent le relief sans
-# délaver la rampe de couleurs en dessous, contrairement à la transparence.
-# Les couches absentes restent en « normal » (simple transparence).
-DEFAULT_BLEND = {'hillshade_multi': 'multiply', 'solar': 'multiply', 'svf': 'multiply'}
+# Affichage de la carte : UNE couche principale (DEFAULT_VIZ) et la couche
+# « précision » (densité de points sol), montrée seule ou superposée en
+# « produit » sur la principale (assombrit les zones interpolées sans masquer
+# le relief). Plus de pile de couches réordonnable.
+PRECISION_VIZ = 'densite_sol'
+VIEW_MODES = ('relief', 'precision', 'both') # principale / précision / les deux
+DEFAULT_VIEW_MODE = 'relief'
+DEFAULT_PRECISION_OPACITY = 0.6 # opacité de la précision en mode « les deux »
# Couche « principale » (mise en avant ; index_tiles.json + affichage en
# direct de la webapp) : le relief orienté, qui fusionne openness et aspect.
@@ -228,7 +231,7 @@ DEFAULT_VIZ = 'relief_oriente'
# carte (panneau, tuiles XYZ, TileJSON, WMTS). Les autres visualisations
# restent calculables avec --only mais ne sont plus proposées.
# None = toutes les visualisations présentes sur disque.
-PANEL_VIZ = ('relief_oriente',)
+PANEL_VIZ = ('relief_oriente', 'densite_sol')
# 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,
@@ -255,14 +258,13 @@ def step_to_keyword(step):
return STEP_TO_KEYWORD.get(step, step)
-def default_layers_present(all_viz_keys):
- """Couches activées par défaut (DEFAULT_LAYERS) présentes sur le disque.
-
- Ordre conservé (DEFAULT_LAYERS) ; sert de base à l'état initial du
- panneau (toutes les autres couches sont désactivées).
- """
- present = set(all_viz_keys)
- return [v for v in DEFAULT_LAYERS if v in present]
+def default_main_layer(all_viz_keys):
+ """Couche principale par défaut présente sur le disque (DEFAULT_VIZ, sinon
+ la première qui n'est pas la précision), ou None."""
+ keys = [k for k in all_viz_keys if k != PRECISION_VIZ]
+ if DEFAULT_VIZ in keys:
+ return DEFAULT_VIZ
+ return keys[0] if keys else None
def cells_with_all_viz(vis_dir, viz_keys, resolutions=(0.5,)):
@@ -289,11 +291,15 @@ def cells_with_all_viz(vis_dir, viz_keys, resolutions=(0.5,)):
# 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
+# Encodage AVIF des sous-tuiles : q75 en 4:2:0, encodé UNE fois depuis le
+# raster d'origine (write_subtiles appelé par tif_to_crop). Mesuré sur le
+# relief orienté (2 dalles réelles, galerie de compression) : l'ancienne
+# chaîne dalle q60 → sous-tuile q55 donnait 18,1 dB / SSIM 0,84 pour
+# 3,9 Mo par dalle ; q75 unique 19,5 dB / SSIM 0,93 pour 8,8 Mo. Au-delà,
+# le 4:2:0 plafonne (la teinte du relief, pixel par pixel, est moyennée par
+# 2 × 2) : seul le 4:4:4 irait plus loin (q75 : 28 dB, 14 Mo).
+# speed 9 : encodage rapide (cf. rendering.AVIF_SPEED).
+_SUBTILE_AVIF_QUALITY = 75
_SUBTILE_AVIF_SPEED = 9
# Vignette de sous-tuile (px) : à la vue d'ensemble chaque sous-tuile active
@@ -303,11 +309,15 @@ _SUBTILE_AVIF_SPEED = 9
_SUBTILE_THUMB_PX = 160
# 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.
+# qualité propre (égale aujourd'hui à celle des rampes de couleur).
_SUBTILE_AVIF_QUALITY_DETAIL = 75
_SUBTILE_DETAIL_VIZ = frozenset({'ortho', 'topo'})
+# Couches en aplats de niveaux codés (cf. rendering.LOSSLESS_GRAY_KEYWORDS) :
+# sous-tuiles en WebP sans perte (niveaux de gris ; 3× plus léger que l'AVIF
+# q100, seul AVIF exact) et vignette intermédiaire sans perte.
+_SUBTILE_LOSSLESS_VIZ = frozenset({'densite_sol'})
+
# 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.
@@ -759,6 +769,74 @@ def _fallback_full_dalle(entries, viz_key, info):
entry['viz'][viz_key] = dict(info)
+def _subtile_ext(viz_key):
+ """Extension des sous-tuiles pleine résolution d'une couche."""
+ return '.webp' if viz_key in _SUBTILE_LOSSLESS_VIZ else '.avif'
+
+
+def _save_subtile_thumb(quad, path):
+ """Écrit la vignette de sous-tuile (redimensionnée à _SUBTILE_THUMB_PX)."""
+ from PIL import Image as PILImage
+ scale = min(1.0, _SUBTILE_THUMB_PX / max(quad.size))
+ out = quad
+ if scale < 1.0:
+ out = quad.resize((max(1, int(quad.size[0] * scale)),
+ max(1, int(quad.size[1] * scale))), PILImage.LANCZOS)
+ out.save(str(path), format='WEBP', quality=80)
+
+
+def write_subtiles(output_dir, dir_name, viz_key, img, k, sub_dir_name='index_subtiles'):
+ """Découpe une image de dalle (PIL, nord en haut) en k × k sous-tuiles :
+ pleine résolution, vignette intermédiaire et vignette.
+
+ Appelée par build_index depuis l'image de dalle, et par le pipeline
+ (rendering.tif_to_crop) directement depuis le raster d'origine : la
+ sous-tuile ne subit alors qu'un seul encodage avec perte (réencoder la
+ dalle AVIF cumulait deux pertes). Encodage : AVIF 4:2:0 q75
+ (_SUBTILE_AVIF_QUALITY, q75 aussi pour ortho/topo), WebP sans perte pour
+ les aplats de niveaux (_SUBTILE_LOSSLESS_VIZ). Lève en cas d'échec.
+ """
+ from PIL import Image as PILImage
+ output_dir = Path(output_dir)
+ out_dir = output_dir / sub_dir_name
+ out_dir.mkdir(parents=True, exist_ok=True)
+ lossless = viz_key in _SUBTILE_LOSSLESS_VIZ
+ ext = _subtile_ext(viz_key)
+ other_ext = '.avif' if ext == '.webp' else '.webp'
+ quality = (_SUBTILE_AVIF_QUALITY_DETAIL if viz_key in _SUBTILE_DETAIL_VIZ
+ else _SUBTILE_AVIF_QUALITY)
+ if img.mode not in ('RGB', 'L'):
+ img = img.convert('RGB')
+ if lossless:
+ img = img.convert('L')
+ W, H = img.size
+ for j in range(k):
+ for i in range(k):
+ stem = f"{dir_name}_{viz_key}_{i}_{j}"
+ # 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))
+ if lossless:
+ quad.save(str(out_dir / (stem + ext)), format='WEBP', lossless=True)
+ else:
+ quad.save(str(out_dir / (stem + ext)), format='AVIF', quality=quality,
+ subsampling='4:2:0', speed=_SUBTILE_AVIF_SPEED)
+ # Anciens fichiers d'une génération précédente (autre format,
+ # vignette 256 px sans taille dans le nom)
+ (out_dir / (stem + other_ext)).unlink(missing_ok=True)
+ (out_dir / (stem + '_thumb.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))), PILImage.LANCZOS)
+ mid_img.save(str(out_dir / (stem + '_mid.webp')), format='WEBP',
+ **({'lossless': True} if lossless else {'quality': 82}))
+ _save_subtile_thumb(quad, out_dir / (stem + f"_thumb{_SUBTILE_THUMB_PX}.webp"))
+
+
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.
@@ -796,42 +874,32 @@ def _build_subtiles(tile, offered_viz_keys, output_dir, sub_dir_name):
}
try:
- resample = getattr(PILImage, 'LANCZOS', 1)
thumb_suffix = f"_thumb{_SUBTILE_THUMB_PX}.webp"
- def _save_thumb(quad, stem):
- """Écrit la vignette de sous-tuile (redimensionnée à _SUBTILE_THUMB_PX)."""
- scale = min(1.0, _SUBTILE_THUMB_PX / max(quad.size))
- out = quad
- if scale < 1.0:
- out = quad.resize((max(1, int(quad.size[0] * scale)),
- max(1, int(quad.size[1] * scale))), resample)
- out.save(str(output_dir / sub_dir_name / (stem + thumb_suffix)),
- format='WEBP', quality=80)
-
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)
+ ext = _subtile_ext(viz_key)
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.
+ # Des sous-tuiles écrites par le pipeline depuis le raster
+ # d'origine (tif_to_crop) sont plus récentes que la dalle : elles
+ # sont gardées telles quelles (un seul encodage avec perte).
src = output_dir / info['full'].split('?')[0]
src_mtime = _mtime(src)
def _fresh(name):
return _cached_file_fresh(output_dir / sub_dir_name / name, src_mtime)
- # AVIF + vignette intermédiaire d'un côté, vignettes de l'autre :
- # un changement de taille de vignette (dalle source inchangée)
- # re-découpe depuis les AVIF de sous-tuiles déjà présents, sans
- # ré-encoder ces derniers (des minutes de CPU par rebuild).
- heavy_fresh = all(_fresh(stem + '.avif') and _fresh(stem + '_mid.webp')
+ # Pleine résolution + vignette intermédiaire d'un côté, vignettes
+ # de l'autre : un changement de taille de vignette (dalle source
+ # inchangée) re-découpe depuis les sous-tuiles déjà présentes, sans
+ # ré-encoder ces dernières (des minutes de CPU par rebuild).
+ heavy_fresh = all(_fresh(stem + ext) and _fresh(stem + '_mid.webp')
for stem in stems.values())
thumbs_fresh = all(_fresh(stem + thumb_suffix) for stem in stems.values())
if heavy_fresh and not thumbs_fresh:
@@ -839,11 +907,11 @@ def _build_subtiles(tile, offered_viz_keys, output_dir, sub_dir_name):
f"{tile['dir_name']}/{viz_key} ({len(stems)} découpages)")
try:
for stem in stems.values():
- with PILImage.open(str(output_dir / sub_dir_name / (stem + '.avif'))) as q:
+ with PILImage.open(str(output_dir / sub_dir_name / (stem + ext))) as q:
q.load()
if q.mode not in ('RGB', 'L'):
q = q.convert('RGB')
- _save_thumb(q, stem)
+ _save_subtile_thumb(q, output_dir / sub_dir_name / (stem + thumb_suffix))
# Vignette 256 px d'une génération précédente
(output_dir / sub_dir_name / (stem + '_thumb.webp')).unlink(missing_ok=True)
except Exception as e:
@@ -871,38 +939,11 @@ def _build_subtiles(tile, offered_viz_keys, output_dir, sub_dir_name):
_fallback_full_dalle(entries, viz_key, info)
if img is None:
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 les anciens crops d'une génération précédente
- # (repli .webp, vignette 256 px sans taille dans le nom)
- (output_dir / sub_dir_name / (stem + '.webp')).unlink(missing_ok=True)
- (output_dir / sub_dir_name / (stem + '_thumb.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)
- _save_thumb(quad, stem)
- if cut_failed:
+ try:
+ write_subtiles(output_dir, tile['dir_name'], viz_key, img, k, sub_dir_name)
+ except Exception as e:
+ logger.warning(f"Encodage des sous-tuiles impossible ({tile['dir_name']}), "
+ f"repli dalle entière pour {viz_key} : {e}")
_fallback_full_dalle(entries, viz_key, info)
continue
# Chaque URL est versionnée par la mtime de SON fichier (cf.
@@ -916,7 +957,7 @@ def _build_subtiles(tile, offered_viz_keys, output_dir, sub_dir_name):
entries[(i, j)]['viz'][viz_key] = {
'thumb': f"{sub_dir_name}/{stem}{thumb_suffix}{_file_v(stem + thumb_suffix)}",
'mid': f"{sub_dir_name}/{stem}_mid.webp{_file_v(stem + '_mid.webp')}",
- 'full': f"{sub_dir_name}/{stem}.avif{_file_v(stem + '.avif')}",
+ 'full': f"{sub_dir_name}/{stem}{ext}{_file_v(stem + ext)}",
}
except Exception as e:
logger.warning(f"Sous-tuilage abandonné pour {tile['dir_name']}: {e}")
diff --git a/lidar_pipeline/mapserve.py b/lidar_pipeline/mapserve.py
index 3c095c9..5672d3e 100644
--- a/lidar_pipeline/mapserve.py
+++ b/lidar_pipeline/mapserve.py
@@ -878,11 +878,12 @@ def tilejson(layer: str, request: Request):
# ---------------------------------------------------------------------------
-# Configuration par défaut de la pile de couches
+# Configuration d'affichage par défaut
# ---------------------------------------------------------------------------
-# L'état réglé depuis la carte (ordre, visibilité, opacité, fusion) peut être
-# figé comme configuration servie à TOUT nouveau navigateur : le réglage utile
-# ne se perd pas dans un localStorage et se partage sans lien.
+# L'affichage réglé depuis la carte (couche principale, mode relief /
+# précision / les deux, opacité de la précision, fond) peut être figé comme
+# configuration servie à TOUT nouveau navigateur : le réglage utile ne se perd
+# pas dans un localStorage et se partage sans lien.
DEFAULTS_FILE = OUTPUT_DIR / ".map-defaults.json"
_defaults_lock = threading.Lock()
@@ -903,38 +904,26 @@ def _clamp01(value, fallback=1.0):
return fallback
-BLEND_MODES = ("normal", "multiply", "screen", "overlay", "soft-light",
- "hard-light", "difference", "luminosity")
-
-
def _sanitize_defaults(req):
- """Retient ce qui est connu et borné : couches réelles, opacités 0–1, fusions valides."""
- known = set(tiles_mod.available_layers(OUTPUT_DIR))
- order = [k for k in (req.order or []) if k in known]
- for key in sorted(known - set(order)):
- order.append(key)
- on = [k for k in (req.on or []) if k in known]
- opacity = {k: _clamp01(v) for k, v in (req.opacity or {}).items() if k in known}
- blend = {k: v for k, v in (req.blend or {}).items()
- if k in known and v in BLEND_MODES}
+ """Retient ce qui est connu et borné : couche réelle, mode valide, opacités 0–1."""
+ from .index import DEFAULT_PRECISION_OPACITY, PRECISION_VIZ, VIEW_MODES
+ known = set(tiles_mod.available_layers(OUTPUT_DIR)) - {PRECISION_VIZ}
base_in = req.base if isinstance(req.base, dict) else {}
base = {"on": bool(base_in.get("on", True)),
"opacity": _clamp01(base_in.get("opacity"), 0.85),
"dark": bool(base_in.get("dark", True))}
- stack = req.stack_blend if req.stack_blend in BLEND_MODES else "normal"
- return {"order": order, "on": on, "opacity": opacity, "blend": blend,
- "base": base, "stack_blend": stack, "saved_at": time.time()}
+ return {"main": req.main if req.main in known else None,
+ "mode": req.mode if req.mode in VIEW_MODES else "relief",
+ "precision_opacity": _clamp01(req.precision_opacity, DEFAULT_PRECISION_OPACITY),
+ "base": base, "saved_at": time.time()}
class DefaultsRequest(BaseModel):
- order: list = Field(default_factory=list,
- description="ordre de pile, du bas vers le haut")
- on: list = Field(default_factory=list, description="couches allumées")
- opacity: dict = Field(default_factory=dict, description="opacité par couche, 0–1")
- blend: dict = Field(default_factory=dict, description="mode de fusion par couche")
+ main: Optional[str] = Field(None, description="couche d'affichage principal")
+ mode: str = Field("relief", description="relief | precision | both")
+ precision_opacity: float = Field(0.6, description="opacité de la précision en mode both, 0–1")
base: dict = Field(default_factory=dict,
description="fond de carte : on, opacity, dark")
- stack_blend: str = Field("normal", description="fusion de la pile sur le fond")
@app.get("/api/map/defaults")
@@ -952,8 +941,8 @@ def set_defaults(req: DefaultsRequest):
DEFAULTS_FILE.parent.mkdir(parents=True, exist_ok=True)
tmp.write_text(json.dumps(data, ensure_ascii=False), encoding="utf-8")
os.replace(tmp, DEFAULTS_FILE)
- logger.info(f"Configuration de pile par défaut enregistrée : "
- f"{len(data['on'])} couche(s) allumée(s)")
+ logger.info(f"Configuration d'affichage par défaut enregistrée : "
+ f"{data['main'] or 'principale du registre'}, mode {data['mode']}")
return {"enregistré": True, "defaults": data}
@@ -973,28 +962,26 @@ def clear_defaults():
@app.get("/api/map/meta")
def map_meta():
"""Couches, zooms, emprise et version — tout l'état initial de l'interface."""
- from .index import DEFAULT_BLEND, DEFAULT_OPACITY, default_layers_present
+ from .index import (DEFAULT_PRECISION_OPACITY, DEFAULT_VIEW_MODE, PRECISION_VIZ,
+ VIEW_MODES, default_main_layer)
infos = _layer_infos()
keys = [i["key"] for i in infos]
- saved = _load_defaults()
+ saved = _load_defaults() or {}
# Configuration figée depuis la carte, sinon réglages du registre (index.py).
- # Couches figées toutes disparues (ex. passage au seul relief orienté) :
- # registre, sinon un nouveau navigateur arrive sur une carte sans LiDAR.
- saved_on = [k for k in saved["on"] if k in keys] if saved else []
- if saved and (saved_on or not saved["on"]):
- default_layers = saved_on
- else:
- default_layers = default_layers_present(keys)
- default_opacity = saved["opacity"] if saved else dict(DEFAULT_OPACITY)
- default_blend = saved["blend"] if saved else dict(DEFAULT_BLEND)
+ # Couche figée disparue, ou fichier de l'ancienne pile (clés order/on) :
+ # couche principale du registre.
+ main = saved.get("main")
+ if main not in keys or main == PRECISION_VIZ:
+ main = default_main_layer(keys)
+ mode = saved.get("mode") if saved.get("mode") in VIEW_MODES else DEFAULT_VIEW_MODE
return {
"layers": infos,
- "default_layers": default_layers,
- "default_opacity": default_opacity,
- "default_blend": default_blend,
- "default_order": saved["order"] if saved else None,
- "default_base": saved["base"] if saved else None,
- "default_stack_blend": saved["stack_blend"] if saved else "normal",
+ "default_main": main,
+ "precision_layer": PRECISION_VIZ if PRECISION_VIZ in keys else None,
+ "default_mode": mode,
+ "default_precision_opacity": _clamp01(saved.get("precision_opacity"),
+ DEFAULT_PRECISION_OPACITY),
+ "default_base": saved.get("base"),
"defaults_saved": bool(saved),
# L'interface consomme des tuiles 512 px (2× moins de requêtes qu'en
# 256 : décisif en HTTP/1.1) ; les clients OSM gardent le 256 canonique.
diff --git a/lidar_pipeline/mapui.py b/lidar_pipeline/mapui.py
index 6bec269..bf24031 100644
--- a/lidar_pipeline/mapui.py
+++ b/lidar_pipeline/mapui.py
@@ -67,19 +67,34 @@ _MAP_HTML = """
-
Haut de pile — devant
-
-
Bas de pile — derrière
-
- Fusion de la pile sur le fond
-
+
+ Affichage
+
+
+
+
+
+
+
+
+
+
+ Opacité de la précision
+
+
+
+
+
+
≤ 0,354≥ 45 pts/m²
+
Points sol par m² : sombre = relief interpolé, noir = aucun point
+
-
+
@@ -220,49 +235,36 @@ button:hover { background: rgba(255, 255, 255, 0.12); border-color: var(--accent
.icon-btn:hover { color: var(--text); background: none; }
.link-btn { flex: 1; font-size: 11.5px; }
-/* Une couche tient sur deux lignes : identité (poignée, œil, nom, ordre) puis
- réglages (fusion + opacité). Les boutons ▲▼ doublent le glisser-déposer —
- seuls utilisables au doigt. */
+/* Ligne du fond de carte : œil, nom, opacité. */
.layer-row { display: grid; grid-template-columns: 14px 24px minmax(0, 1fr) auto 34px;
align-items: center; gap: 6px; padding: 5px 4px; border-radius: 8px; user-select: none; }
.layer-ctl { grid-column: 1 / -1; display: grid; grid-template-columns: minmax(0, 1fr) 84px;
gap: 6px; align-items: center; padding: 2px 0 2px 38px; }
-.layer-move { display: inline-flex; gap: 2px; }
-.layer-move button { padding: 0 5px; font-size: 10px; line-height: 16px; color: var(--muted);
- background: rgba(255, 255, 255, 0.05); }
-.layer-move button:hover:not(:disabled) { color: var(--text); }
-.layer-move button:disabled { opacity: 0.25; cursor: default; }
-.blend { background: rgba(255, 255, 255, 0.06); color: var(--muted);
- border: 1px solid var(--border); border-radius: 6px; padding: 2px 4px;
- font: inherit; font-size: 11px; min-width: 0; width: 100%; height: 22px; }
-.blend:hover { color: var(--text); }
-.blend option { background-color: #10141b; color: var(--text); }
.layer-row:hover { background: rgba(255, 255, 255, 0.05); }
-/* Réorganisation : la destination doit être lisible AVANT de lâcher. Une
- barre d'insertion épaisse marque le point de chute, la ligne déplacée
- s'estompe, et à l'arrivée la couche clignote pour que l'œil la retrouve. */
-.layer-row.dragging { opacity: 0.35; outline: 1px dashed var(--accent); }
-.layer-row.drop-before { box-shadow: inset 0 3px 0 0 var(--accent); }
-.layer-row.drop-after { box-shadow: inset 0 -3px 0 0 var(--accent); }
-.layer-row.drop-before::before, .layer-row.drop-after::before {
- content: '▸'; position: absolute; margin-left: -12px; color: var(--accent);
- font-size: 11px; line-height: 1;
-}
-.layer-row.drop-before::before { margin-top: -14px; }
-.layer-row.drop-after::before { margin-top: 14px; }
-.layer-row { position: relative; }
-.layer-row.just-moved { animation: moved 1s ease-out; }
-@keyframes moved {
- 0% { background: rgba(90, 169, 230, 0.38); }
- 100% { background: transparent; }
-}
-body.reordering, body.reordering * { cursor: grabbing !important; }
-.stack-edge {
- display: flex; align-items: center; gap: 6px;
- font-size: 10px; letter-spacing: 0.04em; text-transform: uppercase;
- color: var(--muted); padding: 2px 6px;
-}
-.stack-edge::after { content: ''; flex: 1; height: 1px; background: var(--border); }
+/* Affichage : une couche principale, la précision seule ou superposée. */
+.view-row { display: grid; grid-template-columns: auto minmax(0, 1fr); align-items: center;
+ gap: 8px; padding: 8px 4px 4px; border-top: 1px solid var(--border); margin-top: 4px; }
+.view-label { font-size: 10.5px; font-weight: 700; letter-spacing: 0.06em; text-transform: uppercase; color: var(--muted); }
+.view-row .layer-name { color: var(--text); }
+.view-sel { background: rgba(255, 255, 255, 0.06); color: var(--text); border: 1px solid var(--border);
+ border-radius: 6px; padding: 2px 4px; font: inherit; font-size: 12px; min-width: 0; width: 100%; height: 24px; }
+.view-sel option { background-color: #10141b; color: var(--text); }
+.mode-seg { display: grid; grid-template-columns: repeat(3, 1fr); margin: 6px 4px 2px;
+ border: 1px solid var(--border); border-radius: 8px; overflow: hidden; }
+.mode-seg button { border: none; border-radius: 0; background: rgba(255, 255, 255, 0.04);
+ color: var(--muted); padding: 6px 4px; font-size: 12px; }
+.mode-seg button + button { border-left: 1px solid var(--border); }
+.mode-seg button:hover { color: var(--text); background: rgba(255, 255, 255, 0.1); }
+.mode-seg button[aria-pressed="true"] { background: var(--accent); color: #0b0e13; font-weight: 600; }
+.mode-seg button:focus-visible { outline: 2px solid var(--accent); outline-offset: -2px; }
+.prec-ctl { display: grid; grid-template-columns: auto minmax(0, 1fr) 34px; gap: 8px;
+ align-items: center; padding: 6px 4px 0; font-size: 11.5px; color: var(--muted); }
+.prec-legend { padding: 8px 4px 2px; }
+.prec-bar { display: grid; grid-template-columns: repeat(16, 1fr); height: 10px;
+ border: 1px solid var(--border); border-radius: 3px; overflow: hidden; }
+.prec-ticks { display: flex; justify-content: space-between; font-size: 10.5px; color: var(--muted);
+ font-variant-numeric: tabular-nums; margin-top: 2px; }
+.prec-cap { font-size: 10.5px; color: var(--muted); margin-top: 2px; }
.grip { color: var(--muted); cursor: grab; text-align: center; font-size: 12px; }
.grip::before { content: '⠿'; }
#baseRow .grip::before { content: ''; }
@@ -371,35 +373,29 @@ _MAP_JS = r"""'use strict';
d'image, aucune rotation CSS, aucune mosaïque maison. */
const UI_VERSION = "__UI_VERSION__";
-const LS_KEY = 'lidarMapLayers_v1';
+const LS_KEY = 'lidarMapView_v2';
const el = (id) => document.getElementById(id);
const SVG_EYE = '';
const SVG_EYE_OFF = '';
-// Modes de fusion proposés (CSS mix-blend-mode) : la superposition de couches
-// ne se réduit pas à de la transparence — « Produit » assombrit le relief sans
-// délaver la rampe du dessous, « Écran » éclaircit, « Superposé »/« Doux »
-// creusent le contraste local, « Différence » fait ressortir les écarts.
-const BLEND_MODES = [
- ['normal', 'Normal'],
- ['multiply', 'Produit'],
- ['screen', 'Écran'],
- ['overlay', 'Superposé'],
- ['soft-light', 'Doux'],
- ['hard-light', 'Lumière crue'],
- ['difference', 'Différence'],
- ['luminosity', 'Luminosité'],
-];
-const BLEND_KEYS = BLEND_MODES.map(m => m[0]);
-const validBlend = (v) => (BLEND_KEYS.includes(v) ? v : null);
-const blendOptions = (current) => BLEND_MODES.map(
- m => '').join('');
+// Affichage : une couche principale (relief orienté) et la couche
+// « précision » (densité de points sol), montrée seule ou superposée en
+// « produit » (mix-blend-mode multiply) sur la principale : les zones peu
+// fiables s'assombrissent sans masquer le relief. Touche P : mode suivant.
+const VIEW_MODES = [['relief', 'Relief'], ['precision', 'Précision'], ['both', 'Les deux']];
+const VIEW_KEYS = VIEW_MODES.map(m => m[0]);
+const validMode = (v) => (VIEW_KEYS.includes(v) ? v : null);
+const clamp01 = (v, d) => {
+ const n = Number(v);
+ return Number.isFinite(n) ? Math.min(1, Math.max(0, n)) : d;
+};
+// Paliers de la précision : niveau k à partir de 0,25 × 2^(k/2) pts/m²
+// (visualizations.density_levels), gris = colormap « gray » sur 0..15.
+const PREC_LEVELS = 16;
let META = null;
let STATE = null;
const LAYERS = new Map(); // clé → L.TileLayer
-let dragKey = null;
const map = L.map('map', {
zoomControl: false,
@@ -506,33 +502,28 @@ function toast(msg, ms) {
}
// --- état persisté --------------------------------------------------------
+function precisionKey() {
+ const p = META && META.precision_layer;
+ return p && META.layers.some(l => l.key === p) ? p : null;
+}
+
+function mainKeys(meta) {
+ return meta.layers.map(l => l.key).filter(k => k !== meta.precision_layer);
+}
+
+function validMain(meta, k) {
+ return mainKeys(meta).includes(k) ? k : null;
+}
+
function defaultState(meta) {
- // `default_order` n'existe que si une configuration a été figée depuis la
- // carte (POST /api/map/defaults) ; sinon on part du registre du pipeline.
- const order = (Array.isArray(meta.default_order) && meta.default_order.length)
- ? meta.default_order.filter(k => meta.layers.some(l => l.key === k))
- : meta.layers.map(l => l.key);
- for (const l of meta.layers) if (!order.includes(l.key)) order.push(l.key);
- const on = {}, opacity = {}, blend = {};
- const wanted = (meta.default_layers || []).filter(k => order.includes(k));
- // Ordre figé depuis la carte : respecté tel quel. Sinon, les couches
- // allumées par défaut passent en bas de pile (bas → haut).
- const stacked = (Array.isArray(meta.default_order) && meta.default_order.length)
- ? order
- : wanted.concat(order.filter(k => !wanted.includes(k)));
- for (const k of stacked) {
- on[k] = wanted.includes(k);
- const o = Number((meta.default_opacity || {})[k]);
- opacity[k] = Number.isFinite(o) ? Math.min(1, Math.max(0, o)) : 1;
- blend[k] = (meta.default_blend || {})[k] || 'normal';
- }
const b = meta.default_base || {};
return {
- order: stacked, on, opacity, blend,
+ main: validMain(meta, meta.default_main) || mainKeys(meta)[0] || null,
+ mode: validMode(meta.default_mode) || 'relief',
+ precOpacity: clamp01(meta.default_precision_opacity, 0.6),
base: { on: b.on !== undefined ? !!b.on : true,
opacity: b.opacity !== undefined ? Number(b.opacity) : 0.85,
dark: b.dark !== undefined ? !!b.dark : true },
- stackBlend: validBlend(meta.default_stack_blend) || 'normal',
};
}
@@ -543,33 +534,19 @@ function loadState(meta) {
const raw = localStorage.getItem(LS_KEY);
if (raw) {
const s = JSON.parse(raw);
- const order = (Array.isArray(s.order) ? s.order : []).filter(k => def.order.includes(k));
- for (const k of def.order) if (!order.includes(k)) order.push(k);
- const on = {}, opacity = {}, blend = {};
- for (const k of order) {
- on[k] = !!(s.on && s.on[k]);
- const o = Number(s.opacity ? s.opacity[k] : NaN);
- opacity[k] = Number.isFinite(o) ? Math.min(1, Math.max(0, o)) : def.opacity[k];
- blend[k] = (s.blend && s.blend[k]) || def.blend[k] || 'normal';
- }
- st = { order, on, opacity, blend,
- base: Object.assign({}, def.base, s.base || {}),
- stackBlend: validBlend(s.stackBlend) || 'normal' };
+ st = { main: validMain(meta, s.main) || def.main,
+ mode: validMode(s.mode) || def.mode,
+ precOpacity: clamp01(s.precOpacity, def.precOpacity),
+ base: Object.assign({}, def.base, s.base || {}) };
}
} catch (e) { st = def; }
+ // Lien partagé : prime sur l'état local (les anciens liens &L=… de la
+ // pile de couches sont ignorés, la vue s'ouvre alors sur l'état courant).
const hash = parseHash();
- if (hash && hash.layers) {
- const order = hash.layers.order.filter(k => def.order.includes(k));
- for (const k of def.order) if (!order.includes(k)) order.push(k);
- st = { order, on: hash.layers.on, opacity: hash.layers.opacity,
- blend: Object.assign({}, st.blend, hash.layers.blend || {}),
- base: hash.base || st.base,
- stackBlend: hash.stackBlend || st.stackBlend || 'normal' };
- for (const k of order) {
- if (st.on[k] === undefined) st.on[k] = false;
- if (st.opacity[k] === undefined) st.opacity[k] = 1;
- if (!st.blend[k]) st.blend[k] = def.blend[k] || 'normal';
- }
+ if (hash) {
+ if (hash.main && validMain(meta, hash.main)) st.main = hash.main;
+ if (hash.view) { st.mode = hash.view.mode; st.precOpacity = hash.view.opacity; }
+ if (hash.base) st.base = hash.base;
}
return st;
}
@@ -579,33 +556,23 @@ function saveState() {
}
// --- hash de vue ----------------------------------------------------------
-// #z/lat/lng&L=clé:opacité:visible,...&B=on:opacité:sombre
+// #z/lat/lng&M=principale&P=mode:opacité&B=on:opacité:sombre
function parseHash() {
const m = /^#([\d.]+)\/(-?[\d.]+)\/(-?[\d.]+)(.*)$/.exec(location.hash || '');
if (!m) return null;
const out = { zoom: parseFloat(m[1]), lat: parseFloat(m[2]), lng: parseFloat(m[3]),
- layers: null, base: null };
+ main: null, view: null, base: null };
const seg = {};
for (const part of m[4].split('&')) {
const i = part.indexOf('=');
if (i > 0) seg[part.slice(0, i)] = part.slice(i + 1);
}
- if (seg.L) {
- // clé:opacité:visible[:fusion] — le 4e champ est optionnel (liens anciens)
- const order = [], on = {}, opacity = {}, blend = {};
- for (const item of seg.L.split(',')) {
- const f = item.split(':');
- if (f.length < 3) continue;
- order.push(f[0]);
- on[f[0]] = f[2] === '1';
- const o = Number(f[1]);
- opacity[f[0]] = Number.isFinite(o) ? Math.min(100, Math.max(0, o)) / 100 : 1;
- const b = validBlend(f[3]);
- if (b) blend[f[0]] = b;
- }
- if (order.length) out.layers = { order, on, opacity, blend };
+ if (seg.M) out.main = seg.M;
+ if (seg.P) {
+ const f = seg.P.split(':');
+ const mode = validMode(f[0]);
+ if (mode) out.view = { mode, opacity: clamp01(Number(f[1]) / 100, 0.6) };
}
- if (seg.S) out.stackBlend = validBlend(seg.S) || 'normal';
if (seg.B) {
const f = seg.B.split(':');
if (f.length === 3) {
@@ -620,12 +587,10 @@ function viewHash(full) {
const c = map.getCenter();
let h = '#' + map.getZoom().toFixed(1) + '/' + c.lat.toFixed(6) + '/' + c.lng.toFixed(6);
if (full) {
- h += '&L=' + STATE.order.map(k =>
- k + ':' + Math.round(STATE.opacity[k] * 100) + ':' + (STATE.on[k] ? '1' : '0')
- + ':' + (STATE.blend[k] || 'normal')).join(',');
+ if (STATE.main) h += '&M=' + STATE.main;
+ h += '&P=' + STATE.mode + ':' + Math.round(STATE.precOpacity * 100);
h += '&B=' + (STATE.base.on ? '1' : '0') + ':' + Math.round(STATE.base.opacity * 100)
+ ':' + (STATE.base.dark ? '1' : '0');
- h += '&S=' + (STATE.stackBlend || 'normal');
}
return h;
}
@@ -662,11 +627,10 @@ const EvenLevelTileLayer = ZoomAnimTileLayer.extend({
});
function buildLayers() {
- // Conteneur isolé : `isolation: isolate` cantonne les mix-blend-mode à la
- // pile LiDAR. Sans lui, « Produit » assombrirait aussi le fond OSM à travers
- // les zones sans donnée. La fusion de la pile ENTIÈRE sur le fond reste
- // réglable à part (STATE.stackBlend). Les panneaux sont réutilisés si déjà
- // créés (rafraîchissement après un run de génération).
+ // Conteneur isolé : `isolation: isolate` cantonne le « produit » de la
+ // précision à la couche principale. Sans lui, il assombrirait aussi le fond
+ // OSM. Les panneaux sont réutilisés si déjà créés (rafraîchissement après
+ // un run de génération).
const stack = map.getPane('lidarStack') || map.createPane('lidarStack');
stack.style.zIndex = 400;
stack.style.isolation = 'isolate';
@@ -703,20 +667,30 @@ function buildLayers() {
applyLayers();
}
+// Mode réellement appliqué : sans couche précision servie, relief seul.
+function effectiveMode() {
+ return precisionKey() ? STATE.mode : 'relief';
+}
+
function applyLayers() {
- const stack = map.getPane('lidarStack');
- if (stack) stack.style.mixBlendMode = STATE.stackBlend || 'normal';
- let z = 400;
- for (const key of STATE.order) {
- const layer = LAYERS.get(key);
- if (!layer) continue;
+ const prec = precisionKey();
+ const mode = effectiveMode();
+ for (const [key, layer] of LAYERS) {
const pane = map.getPane('lidar-' + key);
- z += 1;
- pane.style.zIndex = String(z);
- pane.style.mixBlendMode = STATE.blend[key] || 'normal';
- if (STATE.on[key]) {
+ let show, opacity = 1, blend = 'normal';
+ if (key === prec) {
+ show = mode !== 'relief';
+ if (mode === 'both') { opacity = STATE.precOpacity; blend = 'multiply'; }
+ } else {
+ show = key === STATE.main && mode !== 'precision';
+ }
+ if (pane) {
+ pane.style.zIndex = key === prec ? '402' : '401';
+ pane.style.mixBlendMode = blend;
+ }
+ if (show) {
if (!map.hasLayer(layer)) layer.addTo(map);
- layer.setOpacity(STATE.opacity[key]);
+ layer.setOpacity(opacity);
} else if (map.hasLayer(layer)) {
map.removeLayer(layer);
}
@@ -729,156 +703,60 @@ function applyLayers() {
updateRose();
}
-function moveLayer(key, delta) {
- // STATE.order est bas → haut ; le panneau affiche haut en premier.
- const ord = STATE.order;
- const i = ord.indexOf(key);
- const j = i + delta;
- if (i < 0 || j < 0 || j >= ord.length) return;
- ord.splice(j, 0, ord.splice(i, 1)[0]);
- saveState(); renderPanel(key); applyLayers();
+function setMode(mode, announce) {
+ if (!validMode(mode)) return;
+ STATE.mode = mode;
+ saveState(); renderPanel(); applyLayers();
+ if (announce) toast('Affichage : ' + VIEW_MODES.find(m => m[0] === mode)[1], 1200);
}
-// Insère `src` avant ou après `target` dans l'ordre AFFICHÉ (haut → bas).
-function reorderRelative(src, target, before) {
- const top = STATE.order.slice().reverse();
- const from = top.indexOf(src);
- if (from < 0) return;
- top.splice(from, 1);
- let at = top.indexOf(target);
- if (at < 0) return;
- if (!before) at += 1;
- top.splice(at, 0, src);
- STATE.order = top.reverse();
- saveState(); renderPanel(src); applyLayers();
-}
-
-function clearDropMarks() {
- document.querySelectorAll('.layer-row.drop-before, .layer-row.drop-after')
- .forEach(r => r.classList.remove('drop-before', 'drop-after'));
-}
-
-function renderPanel(movedKey) {
- const host = el('layers');
- host.innerHTML = '';
- const labels = new Map(META.layers.map(l => [l.key, l]));
- // Haut de pile en premier (STATE.order est bas → haut).
- const topFirst = STATE.order.slice().reverse();
- topFirst.forEach((key, idx) => {
- const info = labels.get(key);
- if (!info) return;
- const on = !!STATE.on[key];
- const row = document.createElement('div');
- row.className = 'layer-row' + (on ? ' on' : '')
- + (key === movedKey ? ' just-moved' : '');
- row.dataset.layer = key;
- row.title = info.description || info.label;
- row.innerHTML =
- '' +
- '' +
- '' + info.label + '' +
- '' +
- '' +
- '' +
- '' +
- '' + Math.round(STATE.opacity[key] * 100) + '%' +
- (on ? '
' +
- '' +
- '' +
- '
' : '');
- row.draggable = false;
- row.addEventListener('mousedown', (e) => {
- row.draggable = !!(e.target.closest && e.target.closest('.grip'));
- });
- row.querySelector('.eye').addEventListener('click', () => {
- STATE.on[key] = !STATE.on[key];
- if (STATE.on[key] && STATE.opacity[key] === 0) STATE.opacity[key] = 1;
- saveState(); renderPanel(); applyLayers();
- });
- row.querySelector('.up').addEventListener('click', () => moveLayer(key, 1));
- row.querySelector('.down').addEventListener('click', () => moveLayer(key, -1));
- const blend = row.querySelector('.blend');
- if (blend) blend.addEventListener('change', () => {
- STATE.blend[key] = blend.value;
- saveState();
- const pane = map.getPane('lidar-' + key);
- if (pane) pane.style.mixBlendMode = blend.value;
- });
- const range = row.querySelector('input[type="range"]');
- const pct = row.querySelector('.pct');
- if (range) {
- range.addEventListener('input', () => {
- STATE.opacity[key] = Number(range.value) / 100;
- pct.textContent = range.value + '%';
- const layer = LAYERS.get(key);
- if (layer && map.hasLayer(layer)) layer.setOpacity(STATE.opacity[key]);
- });
- range.addEventListener('change', saveState);
- }
- row.addEventListener('dragstart', (e) => {
- dragKey = key;
- row.classList.add('dragging');
- document.body.classList.add('reordering');
- if (e.dataTransfer) {
- e.dataTransfer.effectAllowed = 'move';
- e.dataTransfer.setData('text/plain', key); // exigé par Firefox
- }
- });
- row.addEventListener('dragend', () => {
- dragKey = null;
- row.classList.remove('dragging');
- row.draggable = false;
- document.body.classList.remove('reordering');
- clearDropMarks();
- });
- row.addEventListener('dragover', (e) => {
- e.preventDefault();
- if (!dragKey || dragKey === key) return;
- // Moitié haute survolée → insertion au-dessus, sinon en dessous.
- const box = row.getBoundingClientRect();
- const before = (e.clientY - box.top) < box.height / 2;
- clearDropMarks();
- row.classList.add(before ? 'drop-before' : 'drop-after');
- });
- row.addEventListener('dragleave', () => {
- row.classList.remove('drop-before', 'drop-after');
- });
- row.addEventListener('drop', (e) => {
- e.preventDefault();
- if (!dragKey || dragKey === key) { clearDropMarks(); return; }
- const before = row.classList.contains('drop-before');
- const src = dragKey;
- clearDropMarks();
- document.body.classList.remove('reordering');
- reorderRelative(src, key, before);
- });
- host.appendChild(row);
- });
- // La couche qui vient d'être déplacée est ramenée dans le champ de vision
- // (le panneau défile dès qu'il y a beaucoup de couches).
- if (movedKey) {
- const moved = host.querySelector('.layer-row.just-moved');
- if (moved && moved.scrollIntoView) moved.scrollIntoView({ block: 'nearest' });
+function renderPrecLegend() {
+ const bar = el('precBar');
+ if (bar.childElementCount) return;
+ for (let k = 0; k < PREC_LEVELS; k++) {
+ const g = Math.round(255 * k / (PREC_LEVELS - 1));
+ const lo = 0.25 * Math.pow(2, k / 2);
+ const hi = 0.25 * Math.pow(2, (k + 1) / 2);
+ const cell = document.createElement('i');
+ cell.style.background = 'rgb(' + g + ',' + g + ',' + g + ')';
+ cell.title = k === 0 ? '≤ ' + hi.toFixed(2).replace('.', ',') + ' pt/m²'
+ : k === PREC_LEVELS - 1 ? '≥ ' + Math.round(lo) + ' pts/m²'
+ : lo.toFixed(lo < 10 ? 1 : 0).replace('.', ',') + '–'
+ + hi.toFixed(hi < 10 ? 1 : 0).replace('.', ',') + ' pts/m²';
+ bar.appendChild(cell);
}
}
+function renderPanel() {
+ const labels = new Map(META.layers.map(l => [l.key, l]));
+ const mains = mainKeys(META);
+ const info = labels.get(STATE.main);
+ const name = el('mainName'), sel = el('mainSel');
+ name.textContent = info ? info.label : 'Aucune couche';
+ name.title = info ? (info.description || info.label) : '';
+ // Plusieurs couches principales possibles : menu ; sinon simple libellé.
+ const many = mains.length > 1;
+ name.hidden = many;
+ sel.hidden = !many;
+ if (many) {
+ sel.innerHTML = mains.map(k => '').join('');
+ }
+ const prec = precisionKey();
+ el('precBlock').hidden = !prec;
+ if (!prec) return;
+ renderPrecLegend();
+ const mode = effectiveMode();
+ el('modeSeg').querySelectorAll('button').forEach(b =>
+ b.setAttribute('aria-pressed', String(b.dataset.mode === mode)));
+ el('precCtl').hidden = mode !== 'both';
+ el('precOpacity').value = Math.round(STATE.precOpacity * 100);
+ el('precPct').textContent = Math.round(STATE.precOpacity * 100) + '%';
+ el('precLegend').hidden = mode === 'relief';
+}
+
function renderBaseRow() {
const b = STATE.base;
- const sel = el('stackBlend');
- if (sel && !sel.options.length) {
- sel.innerHTML = blendOptions(STATE.stackBlend || 'normal');
- sel.addEventListener('change', () => {
- STATE.stackBlend = sel.value;
- saveState(); applyLayers();
- });
- } else if (sel) {
- sel.value = STATE.stackBlend || 'normal';
- }
const toggle = el('baseToggle');
toggle.classList.toggle('on', b.on);
toggle.innerHTML = b.on ? SVG_EYE : SVG_EYE_OFF;
@@ -1062,7 +940,8 @@ function drawRose() {
}
function updateRose() {
- const on = !!(STATE && STATE.on && STATE.on.relief_oriente && LAYERS.get('relief_oriente'));
+ const on = !!(STATE && STATE.main === 'relief_oriente' && effectiveMode() !== 'precision'
+ && LAYERS.get('relief_oriente'));
el('rose').hidden = !on;
if (on) drawRose();
}
@@ -1429,7 +1308,7 @@ window.addEventListener('hashchange', () => {
if (!h) return;
applyingHash = true;
map.setView([h.lat, h.lng], h.zoom);
- if (h.layers) {
+ if (h.main || h.view || h.base) {
STATE = loadState(META);
saveState(); renderPanel(); renderBaseRow(); applyLayers();
}
@@ -1463,21 +1342,37 @@ el('btnShare').addEventListener('click', () => {
});
el('btnDefault').addEventListener('click', () => {
const body = {
- order: STATE.order,
- on: STATE.order.filter(k => STATE.on[k]),
- opacity: Object.fromEntries(STATE.order.map(k => [k, STATE.opacity[k]])),
- blend: Object.fromEntries(STATE.order.map(k => [k, STATE.blend[k] || 'normal'])),
+ main: STATE.main,
+ mode: STATE.mode,
+ precision_opacity: STATE.precOpacity,
base: STATE.base,
- stack_blend: STATE.stackBlend || 'normal',
};
fetch('api/map/defaults', {
method: 'POST', headers: { 'Content-Type': 'application/json' },
body: JSON.stringify(body),
}).then(r => r.ok ? r.json() : Promise.reject('HTTP ' + r.status))
- .then(() => { META.default_order = STATE.order.slice(); toast('Configuration enregistrée comme défaut'); })
+ .then(() => {
+ META.default_main = STATE.main; META.default_mode = STATE.mode;
+ META.default_precision_opacity = STATE.precOpacity; META.default_base = STATE.base;
+ toast('Affichage enregistré comme défaut');
+ })
.catch(msg => toast('Échec de l\'enregistrement : ' + msg, 6000));
});
+el('modeSeg').querySelectorAll('button').forEach(b =>
+ b.addEventListener('click', () => setMode(b.dataset.mode)));
+el('mainSel').addEventListener('change', (e) => {
+ STATE.main = e.target.value;
+ saveState(); renderPanel(); applyLayers();
+});
+el('precOpacity').addEventListener('input', (e) => {
+ STATE.precOpacity = Number(e.target.value) / 100;
+ el('precPct').textContent = e.target.value + '%';
+ const layer = LAYERS.get(precisionKey());
+ if (layer && map.hasLayer(layer)) layer.setOpacity(STATE.precOpacity);
+});
+el('precOpacity').addEventListener('change', saveState);
+
el('btnReset').addEventListener('click', () => {
// Oublie l'état local et reprend la configuration servie par le serveur.
try { localStorage.removeItem(LS_KEY); } catch (e) { /* privé */ }
@@ -1529,6 +1424,14 @@ el('gps').addEventListener('click', () => {
}, { enableHighAccuracy: true, timeout: 15000, maximumAge: 60000 });
});
document.addEventListener('keydown', (e) => {
+ // P : relief → précision → les deux → relief (hors saisie de texte)
+ const typing = e.target && e.target.closest && e.target.closest('input, select, textarea');
+ if ((e.key === 'p' || e.key === 'P') && !typing && !e.ctrlKey && !e.metaKey && !e.altKey
+ && META && STATE && precisionKey()) {
+ const i = VIEW_KEYS.indexOf(STATE.mode);
+ setMode(VIEW_KEYS[(i + 1) % VIEW_KEYS.length], true);
+ return;
+ }
if (e.key === 'Escape') {
el('tilecard').hidden = true; el('usecard').hidden = true;
if (genDrawing) genSetDrawing(false);
diff --git a/lidar_pipeline/pipeline.py b/lidar_pipeline/pipeline.py
index eaf0326..cd50907 100644
--- a/lidar_pipeline/pipeline.py
+++ b/lidar_pipeline/pipeline.py
@@ -90,10 +90,11 @@ from .visualizations import (
generate_flow_accumulation,
generate_anomaly_mask,
generate_relief_oriente,
+ generate_densite_sol,
)
from .gpu import gpu_cleanup, num_gpus, available_gpu_ids, restrict_gpus, safe_gpu_call, gpu_worker_slots
from .ign import generate_ign_overlay
-from .rendering import tif_to_crop
+from .rendering import tif_to_crop, LOSSLESS_GRAY_KEYWORDS
# Ordered list of visualization steps.
@@ -114,6 +115,7 @@ VIZ_STEPS = [
('solar', generate_solar),
('anomaly', generate_anomaly_mask),
('relief_oriente', generate_relief_oriente),
+ ('densite_sol', generate_densite_sol),
('ortho', lambda d, b, v, r: generate_ign_overlay(
d, b, v, r,
layer='ORTHOIMAGERY.ORTHOPHOTOS',
@@ -262,6 +264,8 @@ class LidarArchaeoPipeline:
def _expected_output_path(name, basename, file_vis_dir, output_format='avif'):
"""Return the expected output filename for a visualization step."""
ext = 'avif' if output_format == 'avif' else 'webp'
+ if name in LOSSLESS_GRAY_KEYWORDS:
+ ext = 'webp' # aplats de niveaux : WebP sans perte (rendering.tif_to_crop)
if name == 'pos_open':
return file_vis_dir / f"{basename}_positive_openness.{ext}"
elif name == 'neg_open':
@@ -370,7 +374,7 @@ class LidarArchaeoPipeline:
logger.info(f" Conversion images {fmt_label}:")
for name, tif_file in vis_results.items():
if tif_file and isinstance(tif_file, Path) and tif_file.suffix == '.tif' and tif_file.exists():
- img_file = tif_to_crop(tif_file, file_vis_dir, resolution, keep_tif=self.keep_tif, quality=self.quality, output_format=self.output_format)
+ img_file = tif_to_crop(tif_file, file_vis_dir, resolution, keep_tif=self.keep_tif, quality=self.quality, output_format=self.output_format, subtiles_dir=self.output_dir)
if img_file:
logger.info(f" ✓ {img_file.name}")
@@ -468,6 +472,12 @@ class LidarArchaeoPipeline:
from .dtm import read_dtm_edge_buffer
return abs(read_dtm_edge_buffer(dtm_path) - self.edge_buffer) < 1e-6
+ def _gap_fill_matches(self, dtm_path):
+ """True si le DTM a été comblé par la version courante (tag
+ LIDAR_GAP_FILL, dtm.py) ; absent = ancien comblement à distance fixe."""
+ from .dtm import read_dtm_gap_fill, GAP_FILL_VERSION
+ return read_dtm_gap_fill(dtm_path) == GAP_FILL_VERSION
+
def _fetch_edge_neighbors(self, files):
"""Télécharge les dalles LAZ voisines manquantes (raccord des bords).
@@ -556,6 +566,10 @@ class LidarArchaeoPipeline:
if not self._strip_align_matches(basename, res_suffix):
logger.info(f" DTM{res_suffix} sans calage de faisceaux conforme — régénération (offsets verticaux mesurés et appliqués)")
dtm_path.unlink()
+ elif not self._gap_fill_matches(dtm_path):
+ logger.info(f" DTM{res_suffix} au comblement d'une version antérieure — régénération "
+ f"(vides comblés dans l'enveloppe des points, rayon selon la densité)")
+ dtm_path.unlink()
elif not self._edge_buffer_matches(dtm_path):
from .dtm import read_dtm_edge_buffer
recorded = read_dtm_edge_buffer(dtm_path)
diff --git a/lidar_pipeline/rendering.py b/lidar_pipeline/rendering.py
index 6beb02d..0ca6d54 100644
--- a/lidar_pipeline/rendering.py
+++ b/lidar_pipeline/rendering.py
@@ -160,8 +160,20 @@ COLORMAPS = {
'vmin_mode': 'fixed', 'vmin_val': 0,
'vmax_mode': 'fixed', 'vmax_val': 1,
},
+ # Densité de points sol : niveau 0..15 (visualizations.density_levels)
+ 'densite_sol': {
+ 'cmap': 'gray',
+ 'vmin_mode': 'fixed', 'vmin_val': 0,
+ 'vmax_mode': 'fixed', 'vmax_val': 15,
+ },
}
+# Couches en aplats de gris codés (niveaux discrets) : image niveaux de gris
+# (mode L) encodée SANS PERTE (AVIF q100 en L, seul réglage exact) — l'AVIF q55-60 déplace ~20 % des pixels jusqu'à 4
+# niveaux sur ces aplats. Elles restent à leur résolution propre (1 m pour
+# la densité : ~150–260 Ko par dalle).
+LOSSLESS_GRAY_KEYWORDS = ('densite_sol',)
+
# RGB entries (ortho/topo) are handled specially
RGB_LEGENDS = {
'ortho': {},
@@ -824,7 +836,24 @@ def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False, source_info=None,
return None
-def tif_to_crop(tif_file, vis_dir, resolution, keep_tif=False, quality=60, output_format='avif'):
+def _write_subtiles_from(img, tif_file, vis_dir, resolution, subtiles_dir):
+ """Sous-tuiles d'une dalle depuis son image d'origine (best-effort :
+ en cas d'échec, build_index les découpera depuis la dalle)."""
+ from .index import VIZ_LABELS, _subdivision_k, write_subtiles
+ try:
+ k = _subdivision_k(resolution)
+ stem = Path(tif_file).stem
+ viz_key = next((v for v in sorted(VIZ_LABELS, key=len, reverse=True)
+ if stem.endswith(f"_{v}")), None)
+ if k <= 1 or viz_key is None:
+ return
+ write_subtiles(subtiles_dir, Path(vis_dir).name, viz_key, img, k)
+ except Exception as e:
+ logger.warning(f" Sous-tuiles non écrites ({Path(tif_file).name}) : {e}")
+
+
+def tif_to_crop(tif_file, vis_dir, resolution, keep_tif=False, quality=60, output_format='avif',
+ subtiles_dir=None):
"""Convert GeoTIFF to a cropped visualization image (no legend, no overlay).
Applies colormap and saves the image as a pure 1×1 km square.
@@ -837,6 +866,9 @@ def tif_to_crop(tif_file, vis_dir, resolution, keep_tif=False, quality=60, outpu
keep_tif: If True, keep the source TIFF after conversion.
quality: Image quality (1-100). Use 100 for lossless.
output_format: Output format ('webp' or 'avif').
+ subtiles_dir: dossier de sortie du pipeline : les sous-tuiles de la
+ carte (index_subtiles) y sont écrites depuis l'image d'origine,
+ un seul encodage avec perte (build_index les trouve à jour).
Returns:
Path to output image file, or None on failure.
@@ -844,7 +876,8 @@ def tif_to_crop(tif_file, vis_dir, resolution, keep_tif=False, quality=60, outpu
if not tif_file or not tif_file.exists():
return None
- ext = 'avif' if output_format == 'avif' else 'webp'
+ lossless_gray = any(k in tif_file.stem for k in LOSSLESS_GRAY_KEYWORDS)
+ ext = 'webp' if lossless_gray else ('avif' if output_format == 'avif' else 'webp')
output_file = vis_dir / f"{tif_file.stem}.{ext}"
try:
@@ -892,12 +925,22 @@ def tif_to_crop(tif_file, vis_dir, resolution, keep_tif=False, quality=60, outpu
# Save as AVIF/WebP
img = PILImage.fromarray(rgb_data)
pil_format = 'AVIF' if output_format == 'avif' else 'WEBP'
- if quality >= 100:
+ if lossless_gray:
+ # WebP sans perte en niveaux de gris : exact et 3× plus léger que
+ # l'AVIF q100 (Pillow ignore `lossless=True` en AVIF : q75 avec perte)
+ img = PILImage.fromarray(rgb_data[:, :, 0])
+ img.save(str(output_file), format='WEBP', lossless=True)
+ elif quality >= 100:
img.save(str(output_file), format=pil_format, lossless=True)
else:
img.save(str(output_file), format=pil_format, quality=quality,
**({'speed': AVIF_SPEED} if pil_format == 'AVIF' else {}))
+ # Sous-tuiles de la carte depuis l'image d'origine (après la dalle :
+ # plus récentes qu'elle, build_index ne les ré-encode pas)
+ if subtiles_dir is not None:
+ _write_subtiles_from(img, tif_file, vis_dir, resolution, subtiles_dir)
+
# Delete source TIFF (unless --keep-tif)
if not keep_tif:
tif_file.unlink(missing_ok=True)
diff --git a/lidar_pipeline/tests/test_dtm.py b/lidar_pipeline/tests/test_dtm.py
index 8d45fd2..60583a4 100644
--- a/lidar_pipeline/tests/test_dtm.py
+++ b/lidar_pipeline/tests/test_dtm.py
@@ -140,6 +140,105 @@ class TestInterpolateHoles:
assert np.isnan(filled).all()
+class TestFillSmallGaps:
+ """Comblement borné à l'enveloppe des points (plus de pastilles ni de liseré)."""
+
+ @staticmethod
+ def _grid(n=200, step=2, value=10.0):
+ """Semis régulier de points (1 pixel sur `step`) sur n × n pixels."""
+ dtm = np.full((n, n), np.nan)
+ dtm[::step, ::step] = value
+ return dtm
+
+ def test_isolated_point_is_removed_not_grown(self):
+ from lidar_pipeline.dtm import _fill_small_gaps
+ dtm = np.full((100, 100), np.nan)
+ dtm[50, 50] = 5.0
+ out, filled, removed = _fill_small_gaps(dtm, 0.2)
+ assert np.isnan(out).all() # ni pastille, ni point seul
+ assert (filled, removed) == (0, 1)
+
+ def test_gaps_between_points_filled_without_edge_band(self):
+ from lidar_pipeline.dtm import _fill_small_gaps
+ dtm = self._grid()
+ dtm[:, 100:] = np.nan # grand trou à l'est
+ out, filled, _ = _fill_small_gaps(dtm, 0.2)
+ assert filled > 0
+ assert not np.isnan(out[10:190, 10:99]).any() # vides entre points comblés
+ assert np.isnan(out[:, 99:]).all() # rien d'extrapolé dans le trou
+
+ def test_radius_follows_local_density(self):
+ """Semis clairsemé (1 point / 1,8 m, vides de 2,5 m en diagonale)
+ comblé ; trou de 3 m en zone dense conservé ; trou de 1,6 m comblé."""
+ from lidar_pipeline.dtm import _fill_small_gaps
+ sparse = self._grid(step=9)
+ out, _, _ = _fill_small_gaps(sparse, 0.2)
+ assert not np.isnan(out[30:170, 30:170]).any()
+ dense = self._grid(step=1)
+ dense[90:105, 90:105] = np.nan # trou de 3 m dans un semis plein
+ dense[40:48, 40:48] = np.nan # trou de 1,6 m (voiture)
+ out, _, _ = _fill_small_gaps(dense, 0.2)
+ assert np.isnan(out[95:100, 95:100]).all()
+ assert not np.isnan(out[40:48, 40:48]).any()
+
+ def test_morph_disk_matches_scipy(self):
+ from scipy import ndimage as nd
+ from lidar_pipeline.dtm import _morph_disk
+ rng = np.random.default_rng(0)
+ mask = rng.random((60, 70)) < 0.05
+ square = np.ones((3, 3), dtype=bool)
+ cross = nd.generate_binary_structure(2, 1)
+ ref = mask
+ for i in range(5):
+ ref = nd.binary_dilation(ref, structure=square if i % 2 == 0 else cross)
+ assert (_morph_disk(mask, 5) == ref).all()
+ ref_e = ref
+ for i in range(5):
+ ref_e = nd.binary_erosion(ref_e, structure=square if i % 2 == 0 else cross,
+ border_value=1)
+ assert (_morph_disk(ref, 5, erode=True) == ref_e).all()
+
+ def test_dtm_records_gap_fill_version(self, tmp_output_dir):
+ import laspy
+ import rasterio
+ from lidar_pipeline.dtm import (create_dtm_fast, read_dtm_gap_fill,
+ GAP_FILL_VERSION)
+ hdr = laspy.LasHeader(version='1.2', point_format=0)
+ las = laspy.LasData(hdr)
+ las.x, las.y, las.z = [0.5, 1.5, 0.5, 1.5], [0.5, 0.5, 1.5, 1.5], [10.0] * 4
+ las.write(str(tmp_output_dir / "g.las"))
+ out = create_dtm_fast(tmp_output_dir / "g.las", "g", tmp_output_dir, 1.0,
+ force=True, strip_align=False)
+ assert read_dtm_gap_fill(out) == GAP_FILL_VERSION
+ with rasterio.open(str(out)) as src:
+ src.tags() # lisible
+ legacy = tmp_output_dir / "legacy.tif"
+ with rasterio.open(str(legacy), "w", driver="GTiff", width=2, height=2,
+ count=1, dtype="float32") as dst:
+ dst.write(np.zeros((1, 2, 2), dtype="float32"))
+ assert read_dtm_gap_fill(legacy) == 1
+
+
+class TestDensitySidecar:
+ def test_dtm_writes_ground_density(self, tmp_output_dir):
+ """4 points par m² sur 10 × 10 m : densité 4 au cœur, grille de 1 m."""
+ import laspy
+ import rasterio
+ from lidar_pipeline.dtm import create_dtm_fast, density_path
+ g = (np.arange(20) + 0.25) / 2.0 # pas de 0,5 m
+ xx, yy = np.meshgrid(g, g)
+ hdr = laspy.LasHeader(version='1.2', point_format=0)
+ las = laspy.LasData(hdr)
+ las.x, las.y, las.z = xx.ravel(), yy.ravel(), np.full(xx.size, 10.0)
+ las.write(str(tmp_output_dir / "d.las"))
+ out = create_dtm_fast(tmp_output_dir / "d.las", "d", tmp_output_dir, 0.5,
+ force=True, strip_align=False)
+ with rasterio.open(density_path(out)) as src:
+ dens = src.read(1)
+ assert abs(src.transform.a - 1.0) < 1e-9
+ assert np.allclose(dens[2:-2, 2:-2], 4.0)
+
+
class TestDetectGroundMethod:
def _make_mock_las(self, num_returns, z_values):
"""Create a mock laspy object with specified NumberOfReturns and z."""
diff --git a/lidar_pipeline/tests/test_index.py b/lidar_pipeline/tests/test_index.py
index 8151f1d..9e05d61 100644
--- a/lidar_pipeline/tests/test_index.py
+++ b/lidar_pipeline/tests/test_index.py
@@ -287,17 +287,15 @@ def test_subtiles_cover_all_layers(tmp_path):
def test_panel_restricted_to_kept_layers():
- """Le panneau est restreint aux couches conservées (PANEL_VIZ).
-
- Les couches allumées par défaut (DEFAULT_LAYERS) y figurent toutes — le
- panneau peut en proposer davantage (ex. pente éteinte par défaut), et
- chaque opacité par défaut vise une couche allumée.
- """
- from lidar_pipeline.index import (PANEL_VIZ, DEFAULT_LAYERS, DEFAULT_OPACITY,
- KEYWORD_TO_STEP)
+ """Le panneau est restreint aux couches conservées (PANEL_VIZ) : la couche
+ principale par défaut et la précision y figurent."""
+ from lidar_pipeline.index import (PANEL_VIZ, DEFAULT_VIZ, PRECISION_VIZ, VIEW_MODES,
+ DEFAULT_VIEW_MODE, KEYWORD_TO_STEP, default_main_layer)
from lidar_pipeline.pipeline import VIZ_STEPS
- assert set(DEFAULT_LAYERS) <= set(PANEL_VIZ)
- assert set(DEFAULT_OPACITY) <= set(DEFAULT_LAYERS)
+ assert DEFAULT_VIZ in PANEL_VIZ and PRECISION_VIZ in PANEL_VIZ
+ assert DEFAULT_VIEW_MODE in VIEW_MODES
+ assert default_main_layer(["densite_sol", "aspect"]) == "aspect"
+ assert default_main_layer(["densite_sol"]) is None
# Chaque couche conservée correspond à une étape --only valide
steps = {name for name, _ in VIZ_STEPS}
for key in PANEL_VIZ:
diff --git a/lidar_pipeline/tests/test_mapserve.py b/lidar_pipeline/tests/test_mapserve.py
index b074a6e..841a0e1 100644
--- a/lidar_pipeline/tests/test_mapserve.py
+++ b/lidar_pipeline/tests/test_mapserve.py
@@ -210,7 +210,9 @@ def test_map_meta(tmp_path, monkeypatch):
assert meta["tile_size"] == 512 and meta["zoom_offset"] == -1
assert meta["max_native_zoom"] == 19
assert meta["bounds"] and meta["stamp"] > 0
- assert meta["default_layers"] == ["relief_oriente"]
+ assert meta["default_main"] == "relief_oriente"
+ assert meta["precision_layer"] is None # densité absente de la fixture
+ assert meta["default_mode"] == "relief"
def test_only_panel_layers_are_served(tmp_path, monkeypatch):
@@ -256,31 +258,31 @@ def test_healthz_and_assets(tmp_path, monkeypatch):
mapserve.assets("introuvable.js")
-def test_ui_offers_reorder_and_blend_modes():
- """Pile de couches : réordonnancement au doigt ET modes de fusion.
+def test_ui_main_layer_and_precision_modes():
+ """Une couche d'affichage principal et la précision : relief seul,
+ précision seule, ou précision en « produit » sur le relief.
- La superposition ne se réduit pas à de la transparence : chaque couche
- porte un mix-blend-mode, cantonné à la pile LiDAR par un conteneur isolé
- (sinon « Produit » assombrirait aussi le fond de carte).
+ Plus de pile réordonnable ni de fusion par couche. Le produit est
+ cantonné à la pile LiDAR par un conteneur isolé (sinon il assombrirait
+ aussi le fond de carte). Bascule rapide au clavier (touche P).
"""
from lidar_pipeline.mapui import _MAP_CSS, _MAP_HTML, _MAP_JS
- # Modes proposés
- for mode in ("multiply", "screen", "overlay", "soft-light", "difference",
- "luminosity"):
- assert f"'{mode}'" in _MAP_JS, mode
- # Réordonnancement : glisser-déposer ET boutons (seuls utilisables au doigt)
- assert "function moveLayer(" in _MAP_JS
- assert "moveLayer(key, 1)" in _MAP_JS and "moveLayer(key, -1)" in _MAP_JS
- assert "dragstart" in _MAP_JS and "layer-move" in _MAP_CSS
- # Pile isolée + fusion par couche et fusion de la pile sur le fond
- assert "createPane('lidarStack')" in _MAP_JS
- assert "isolation = 'isolate'" in _MAP_JS
- assert "map.createPane(paneName, stack)" in _MAP_JS
- assert "pane.style.mixBlendMode" in _MAP_JS
- assert "STATE.stackBlend" in _MAP_JS and 'id="stackBlend"' in _MAP_HTML
- # Le lien partagé transporte la fusion (4e champ, rétrocompatible)
- assert "(STATE.blend[k] || 'normal')" in _MAP_JS
- assert "const b = validBlend(f[3]);" in _MAP_JS
+ for mode in ('data-mode="relief"', 'data-mode="precision"', 'data-mode="both"'):
+ assert mode in _MAP_HTML
+ assert 'id="precOpacity"' in _MAP_HTML and 'id="precLegend"' in _MAP_HTML
+ assert 'id="mainSel"' in _MAP_HTML
+ assert "blend = 'multiply'" in _MAP_JS
+ assert "createPane('lidarStack')" in _MAP_JS and "isolation = 'isolate'" in _MAP_JS
+ assert "e.key === 'p'" in _MAP_JS
+ # Légende en 16 paliers exacts (échelle log de la densité)
+ assert "const PREC_LEVELS = 16;" in _MAP_JS and "repeat(16, 1fr)" in _MAP_CSS
+ # Pile et fusions supprimées
+ for gone in ("function moveLayer(", "reorderRelative", "stackBlend", "dragstart", "blendOptions"):
+ assert gone not in _MAP_JS, gone
+ assert "Haut de pile" not in _MAP_HTML
+ # Lien partagé : principale, mode:opacité, fond
+ assert "'&M=' + STATE.main" in _MAP_JS
+ assert "'&P=' + STATE.mode" in _MAP_JS
def test_map_ui_tile_load_indicator_and_selection_pane():
@@ -301,85 +303,73 @@ def test_map_ui_tile_load_indicator_and_selection_pane():
def test_defaults_roundtrip(tmp_path, monkeypatch):
- """L'état réglé depuis la carte devient la configuration servie à tous."""
+ """L'affichage réglé depuis la carte devient la configuration servie à tous."""
import lidar_pipeline.mapserve as mapserve
- _setup(tmp_path, monkeypatch)
+ _setup(tmp_path, monkeypatch, layers=("relief_oriente", "aspect", "densite_sol"))
monkeypatch.setattr(mapserve, "DEFAULTS_FILE", tmp_path / ".map-defaults.json")
# Sans enregistrement : réglages du registre (index.py)
meta = mapserve.map_meta()
assert meta["defaults_saved"] is False
- assert meta["default_order"] is None
- assert meta["default_stack_blend"] == "normal"
+ assert meta["default_main"] == "relief_oriente"
+ assert meta["precision_layer"] == "densite_sol"
+ assert meta["default_mode"] == "relief"
+ assert meta["default_precision_opacity"] == 0.6
- req = mapserve.DefaultsRequest(
- order=["slope", "aspect"], on=["aspect"],
- opacity={"aspect": 0.6, "slope": 1.0},
- blend={"aspect": "multiply", "slope": "normal"},
- base={"on": True, "opacity": 0.4, "dark": False},
- stack_blend="soft-light")
+ req = mapserve.DefaultsRequest(main="aspect", mode="both", precision_opacity=0.4,
+ base={"on": True, "opacity": 0.4, "dark": False})
assert mapserve.set_defaults(req)["enregistré"] is True
meta = mapserve.map_meta()
assert meta["defaults_saved"] is True
- assert meta["default_order"] == ["slope", "aspect"]
- assert meta["default_layers"] == ["aspect"]
- assert meta["default_opacity"]["aspect"] == 0.6
- assert meta["default_blend"]["aspect"] == "multiply"
+ assert meta["default_main"] == "aspect"
+ assert meta["default_mode"] == "both"
+ assert meta["default_precision_opacity"] == 0.4
assert meta["default_base"] == {"on": True, "opacity": 0.4, "dark": False}
- assert meta["default_stack_blend"] == "soft-light"
# Persisté sur disque : survit au redémarrage du conteneur
assert (tmp_path / ".map-defaults.json").is_file()
- assert mapserve.get_defaults()["defaults"]["on"] == ["aspect"]
+ assert mapserve.get_defaults()["defaults"]["mode"] == "both"
# Retrait : retour aux réglages du registre
assert mapserve.clear_defaults()["supprimé"] is True
assert mapserve.map_meta()["defaults_saved"] is False
-def test_defaults_obsolete_layers_fall_back(tmp_path, monkeypatch):
- """Défauts figés sur des couches qui ne sont plus servies : registre.
-
- Sinon un nouveau navigateur arrive sur une carte sans aucune couche
- LiDAR allumée (cas de la prod après le passage au seul relief orienté).
- """
+def test_defaults_legacy_stack_file_falls_back(tmp_path, monkeypatch):
+ """Fichier de l'ancienne pile (order/on/blend) ou couche disparue :
+ couche principale du registre, mode relief ; le fond est conservé."""
import json
import lidar_pipeline.mapserve as mapserve
_setup(tmp_path, monkeypatch, layers=("relief_oriente",))
defaults = tmp_path / ".map-defaults.json"
monkeypatch.setattr(mapserve, "DEFAULTS_FILE", defaults)
defaults.write_text(json.dumps({
- "order": ["aspect", "positive_openness"],
- "on": ["aspect", "positive_openness"], "opacity": {}, "blend": {},
- "base": {"on": True, "opacity": 0.85, "dark": True},
- "stack_blend": "normal"}), encoding="utf-8")
+ "order": ["aspect"], "on": ["aspect"], "opacity": {}, "blend": {},
+ "base": {"on": True, "opacity": 0.5, "dark": True}, "stack_blend": "normal"}),
+ encoding="utf-8")
meta = mapserve.map_meta()
- assert meta["default_layers"] == ["relief_oriente"]
- # Tout éteint volontairement (couches connues) : respecté
- defaults.write_text(json.dumps({
- "order": ["relief_oriente"], "on": [], "opacity": {}, "blend": {},
- "base": {}, "stack_blend": "normal"}), encoding="utf-8")
- assert mapserve.map_meta()["default_layers"] == []
+ assert meta["default_main"] == "relief_oriente"
+ assert meta["default_mode"] == "relief"
+ assert meta["default_base"]["opacity"] == 0.5
+ defaults.write_text(json.dumps({"main": "aspect", "mode": "vaudou"}), encoding="utf-8")
+ meta = mapserve.map_meta()
+ assert meta["default_main"] == "relief_oriente" and meta["default_mode"] == "relief"
def test_defaults_sanitised(tmp_path, monkeypatch):
- """Couches inconnues écartées, opacités bornées, fusions validées."""
+ """Couche inconnue ou précision refusée comme principale, mode validé, opacités bornées."""
import lidar_pipeline.mapserve as mapserve
- _setup(tmp_path, monkeypatch)
+ _setup(tmp_path, monkeypatch, layers=("relief_oriente", "densite_sol"))
monkeypatch.setattr(mapserve, "DEFAULTS_FILE", tmp_path / ".map-defaults.json")
data = mapserve.set_defaults(mapserve.DefaultsRequest(
- order=["aspect", "inexistante"], on=["aspect", "inexistante"],
- opacity={"aspect": 5, "slope": "x", "inexistante": 0.5},
- blend={"aspect": "vaudou", "slope": "screen"},
- base={"opacity": -3}, stack_blend="vaudou"))["defaults"]
- assert "inexistante" not in data["order"] and "inexistante" not in data["on"]
- # Les couches présentes mais non citées complètent l'ordre
- assert set(data["order"]) == {"aspect", "slope"}
- assert data["opacity"]["aspect"] == 1.0 # borné à 1
- assert data["opacity"]["slope"] == 1.0 # valeur illisible → opaque
- assert data["blend"] == {"slope": "screen"} # mode inconnu écarté
+ main="inexistante", mode="vaudou", precision_opacity=5,
+ base={"opacity": -3}))["defaults"]
+ assert data["main"] is None
+ assert data["mode"] == "relief"
+ assert data["precision_opacity"] == 1.0 # borné à 1
assert data["base"]["opacity"] == 0.0 # borné à 0
- assert data["stack_blend"] == "normal" # repli
+ data = mapserve.set_defaults(mapserve.DefaultsRequest(main="densite_sol"))["defaults"]
+ assert data["main"] is None # la précision n'est pas une principale
def test_ui_applies_server_defaults():
@@ -387,38 +377,13 @@ def test_ui_applies_server_defaults():
from lidar_pipeline.mapui import _MAP_HTML, _MAP_JS
assert 'id="btnDefault"' in _MAP_HTML and 'id="btnReset"' in _MAP_HTML
assert "api/map/defaults" in _MAP_JS
- assert "meta.default_order" in _MAP_JS
+ assert "meta.default_main" in _MAP_JS and "meta.default_mode" in _MAP_JS
assert "meta.default_base" in _MAP_JS
- assert "validBlend(meta.default_stack_blend)" in _MAP_JS
+ assert "meta.default_precision_opacity" in _MAP_JS
# Réinitialiser oublie l'état local avant de reprendre celui du serveur
assert "localStorage.removeItem(LS_KEY)" in _MAP_JS
-def test_ui_reorder_shows_destination():
- """Réorganiser doit être lisible : repère de chute, fantôme, arrivée signalée.
-
- Sans repère, on lâche une couche sans savoir où elle atterrit — la pile
- n'est pas un détail, elle décide de ce qu'on voit.
- """
- from lidar_pipeline.mapui import _MAP_CSS, _MAP_HTML, _MAP_JS
- # Barre d'insertion selon la moitié survolée, effacée entre deux survols
- assert ".layer-row.drop-before" in _MAP_CSS and ".layer-row.drop-after" in _MAP_CSS
- assert "clearDropMarks" in _MAP_JS
- assert "box.height / 2" in _MAP_JS
- assert "row.classList.add(before ? 'drop-before' : 'drop-after')" in _MAP_JS
- # Ligne déplacée estompée pendant le glisser
- assert ".layer-row.dragging" in _MAP_CSS and "classList.add('dragging')" in _MAP_JS
- # Insertion relative (avant/après la ligne visée), pas un simple échange
- assert "function reorderRelative(" in _MAP_JS
- # Arrivée mise en évidence puis ramenée dans le champ de vision
- assert "just-moved" in _MAP_CSS and "@keyframes moved" in _MAP_CSS
- assert "scrollIntoView" in _MAP_JS
- assert "renderPanel(key)" in _MAP_JS and "renderPanel(src)" in _MAP_JS
- # Extrémités de la pile nommées
- assert 'id="edgeTop"' in _MAP_HTML and 'id="edgeBottom"' in _MAP_HTML
- assert "Haut de pile" in _MAP_HTML and "Bas de pile" in _MAP_HTML
-
-
def test_ui_uses_standard_tilelayer():
"""L'interface s'appuie sur le LOD natif de Leaflet, pas sur un palier maison."""
from lidar_pipeline.mapui import _MAP_JS, render_html, ui_version
diff --git a/lidar_pipeline/tests/test_pipeline.py b/lidar_pipeline/tests/test_pipeline.py
index d0ca77b..21e9d25 100644
--- a/lidar_pipeline/tests/test_pipeline.py
+++ b/lidar_pipeline/tests/test_pipeline.py
@@ -20,16 +20,16 @@ class TestVizSteps:
assert len(names) == len(set(names)), "VIZ_STEPS has duplicate names"
def test_expected_visualization_count(self):
- """Should have 16 visualizations (14 terrain + ortho + topo)."""
+ """17 visualisations : 14 terrain + densité de points + ortho + topo."""
from lidar_pipeline.pipeline import VIZ_STEPS
- assert len(VIZ_STEPS) == 16
+ assert len(VIZ_STEPS) == 17
- def test_default_run_produces_only_relief(self, tmp_path):
- """Sans --only : seule la couche affichée (relief orienté) est produite ;
- --only reste libre pour les autres visualisations."""
+ def test_default_run_produces_only_panel_layers(self, tmp_path):
+ """Sans --only : seules les couches affichées (relief orienté, densité
+ de points) sont produites ; --only reste libre pour les autres."""
from lidar_pipeline.pipeline import LidarArchaeoPipeline
p = LidarArchaeoPipeline(tmp_path, tmp_path / "out")
- assert [n for n, _ in p.viz_steps] == ["relief_oriente"]
+ assert [n for n, _ in p.viz_steps] == ["relief_oriente", "densite_sol"]
p = LidarArchaeoPipeline(tmp_path, tmp_path / "out2", only_viz=["slope"])
assert [n for n, _ in p.viz_steps] == ["slope"]
diff --git a/lidar_pipeline/tests/test_rendering.py b/lidar_pipeline/tests/test_rendering.py
index 0993fcc..41b0886 100644
--- a/lidar_pipeline/tests/test_rendering.py
+++ b/lidar_pipeline/tests/test_rendering.py
@@ -217,3 +217,52 @@ class TestCoreTileWindow:
(637900, 6626900, 639100, 6628100), 1200)
with rasterio.open(tif) as src:
assert _core_tile_window(tif, src) is None
+
+
+class TestDensiteSolCrop:
+ """Densité de points : WebP sans perte en niveaux de gris, sous-tuiles
+ écrites depuis l'image d'origine (un seul encodage)."""
+
+ def test_lossless_gray_levels_and_subtiles(self, tmp_path):
+ from PIL import Image as PILImage
+ from lidar_pipeline.rendering import tif_to_crop
+ base = "LHD_FXX_0660_6701_PTS_LAMB93_IGN69"
+ vis = tmp_path / "visualisations" / f"{base}_r0p2"
+ vis.mkdir(parents=True)
+ levels = np.tile(np.arange(16, dtype=np.float32).repeat(4), (64, 1)) # 64 × 64
+ tif = TestTifToCrop._write_named_tif(vis, f"{base}_densite_sol.tif", levels)
+ out = tif_to_crop(tif, vis, 0.2, output_format='avif', subtiles_dir=tmp_path)
+ assert out is not None and out.suffix == ".webp"
+ # WebP n'a pas de mode gris : relu en RGB à canaux égaux
+ rgb = np.asarray(PILImage.open(str(out)).convert("RGB")).astype(int)
+ assert (rgb[..., 0] == rgb[..., 1]).all() and (rgb[..., 1] == rgb[..., 2]).all()
+ img = PILImage.fromarray(rgb[..., 0].astype(np.uint8))
+ row = np.asarray(img)[0, ::4].astype(int)
+ assert len(set(row.tolist())) == 16 and np.all(np.diff(row) > 0)
+ # Sous-tuiles 2 × 2 (0,2 m/px) en WebP sans perte, identiques à la dalle
+ q = tmp_path / "index_subtiles" / f"{base}_r0p2_densite_sol_0_1.webp"
+ assert q.exists()
+ assert np.array_equal(np.asarray(PILImage.open(str(q)).convert("L")), np.asarray(img)[:32, :32])
+ assert (tmp_path / "index_subtiles" / f"{base}_r0p2_densite_sol_0_1_mid.webp").exists()
+
+ def test_relief_subtiles_encoded_from_source(self, tmp_path):
+ """Relief : sous-tuiles AVIF écrites par tif_to_crop, plus récentes que la dalle."""
+ import pytest
+ from lidar_pipeline.rendering import tif_to_crop
+ base = "LHD_FXX_0660_6701_PTS_LAMB93_IGN69"
+ vis = tmp_path / "visualisations" / f"{base}_r0p2"
+ vis.mkdir(parents=True)
+ rgb = np.random.default_rng(1).integers(0, 255, (3, 64, 64)).astype('uint8')
+ tif = vis / f"{base}_relief_oriente.tif"
+ with rasterio.open(tif, 'w', driver='GTiff', height=64, width=64, count=3,
+ dtype='uint8', crs='EPSG:2154',
+ transform=from_bounds(660000, 6700000, 661000, 6701000, 64, 64)) as dst:
+ dst.write(rgb)
+ try:
+ out = tif_to_crop(tif, vis, 0.2, output_format='avif', subtiles_dir=tmp_path)
+ except Exception:
+ pytest.skip("encodeur AVIF indisponible")
+ sub = tmp_path / "index_subtiles"
+ quads = sorted(sub.glob(f"{base}_r0p2_relief_oriente_?_?.avif"))
+ assert len(quads) == 4
+ assert all(q.stat().st_mtime_ns >= out.stat().st_mtime_ns for q in quads)
diff --git a/lidar_pipeline/tests/test_tiles.py b/lidar_pipeline/tests/test_tiles.py
index afaf62f..69675b0 100644
--- a/lidar_pipeline/tests/test_tiles.py
+++ b/lidar_pipeline/tests/test_tiles.py
@@ -687,3 +687,28 @@ def test_existing_cache_without_registry_is_refreshed_once(tmp_path):
tiles._seen_cache.clear()
tiles.source_index(tmp_path, force=True)
assert tiles.cached_tile(tmp_path, "slope", z, x, y)[1] == "fresh"
+
+
+def test_webp_subtiles_indexed_and_rendered_nearest(tmp_path):
+ """Couche densité : quadrants .webp sans perte reconnus comme palier fin,
+ et rendus au plus proche voisin (16 gris exacts même agrandis)."""
+ from PIL import Image
+ from lidar_pipeline import tiles
+ base = _make_dalle(tmp_path, 1054, 6882, ["densite_sol"])
+ sub = tmp_path / "index_subtiles"
+ sub.mkdir(parents=True, exist_ok=True)
+ for i in range(2):
+ for j in range(2):
+ im = Image.new("L", (8, 8), 0)
+ im.paste(255, (0, 0, 4, 8)) # moitié blanche, moitié noire
+ im.save(str(sub / f"{base}_densite_sol_{i}_{j}.webp"), format="WEBP", lossless=True)
+ tiers = tiles.source_index(tmp_path, force=True)["densite_sol"][(1054, 6882)]
+ quads = next(t for t in tiers if len(t) == 4)
+ assert all(q.path.suffix == ".webp" for q in quads)
+ assert "densite_sol" in tiles.NEAREST_LAYERS
+ # z17 : résolution visée (~0,8 m) atteinte par les quadrants (0,5 m ici)
+ x, y = _tile_of_cell(1054, 6882, 17)
+ img = tiles.render_tile(tmp_path, "densite_sol", 17, x, y)
+ assert img is not None
+ grays = {p[0] for p in img.getdata() if p[3] == 255}
+ assert grays <= {0, 255} # aucun gris intermédiaire inventé
diff --git a/lidar_pipeline/tests/test_visualizations.py b/lidar_pipeline/tests/test_visualizations.py
index 0cc4be4..e598d2b 100644
--- a/lidar_pipeline/tests/test_visualizations.py
+++ b/lidar_pipeline/tests/test_visualizations.py
@@ -553,3 +553,30 @@ class TestPriorityFlood:
result = _priority_flood(dem, nodata)
assert result[2, 2] == 999.0
+
+
+class TestDensiteSol:
+ def test_density_levels_log_scale(self):
+ from lidar_pipeline.visualizations import density_levels
+ d = np.array([0.0, np.nan, 0.25, 0.36, 0.5, 1.0, 45.0, 45.3, 1000.0])
+ assert density_levels(d).tolist() == [0, 0, 0, 1, 2, 4, 14, 15, 15]
+
+ def test_generate_reads_density_sidecar(self, tmp_path):
+ import rasterio
+ from rasterio.transform import from_bounds
+ from lidar_pipeline.dtm import density_path
+ from lidar_pipeline.visualizations import generate_densite_sol
+ dem = tmp_path / "T_dtm_r0p2.tif"
+ dem.touch()
+ with rasterio.open(density_path(dem), 'w', driver='GTiff', width=4, height=1,
+ count=1, dtype='float32', crs='EPSG:2154',
+ transform=from_bounds(0, 0, 4, 1, 4, 1)) as dst:
+ dst.write(np.array([[0.0, 1.0, 4.0, 64.0]], dtype='float32'), 1)
+ out = generate_densite_sol(dem, "T", tmp_path, 0.2)
+ with rasterio.open(out) as src:
+ assert src.read(1).tolist() == [[0, 4, 8, 15]]
+ assert src.width == 4 # grille de 1 m conservée
+
+ def test_missing_sidecar_returns_none(self, tmp_path):
+ from lidar_pipeline.visualizations import generate_densite_sol
+ assert generate_densite_sol(tmp_path / "X_dtm.tif", "X", tmp_path, 0.2) is None
diff --git a/lidar_pipeline/tiles.py b/lidar_pipeline/tiles.py
index 272daef..5c9c261 100644
--- a/lidar_pipeline/tiles.py
+++ b/lidar_pipeline/tiles.py
@@ -48,6 +48,9 @@ AVIF_SPEED = 9 # encodage rapide (cf. rendering.AVIF_SPEED : ×7, +3 %
# l'active quand la bande passante compte (consultation mobile).
PNG_PALETTE = os.environ.get("LIDAR_TILE_PNG_PALETTE", "") == "1"
+# Couches en aplats de niveaux codés, rééchantillonnées au plus proche voisin
+NEAREST_LAYERS = frozenset({'densite_sol'})
+
# Demi-circonférence équatoriale : emprise du Web Mercator (EPSG:3857).
ORIGIN = 20037508.342789244
@@ -452,8 +455,10 @@ def _build_index(output_dir):
quads = {}
for f in sub_dir.iterdir():
name = f.name
+ # « .webp » plein en dernier : les vignettes finissent aussi en .webp
for suffix, tier in ((".avif", "full"), ("_mid.webp", "mid"),
- (f"_thumb{_SUBTILE_THUMB_PX}.webp", "thumb")):
+ (f"_thumb{_SUBTILE_THUMB_PX}.webp", "thumb"),
+ (".webp", "full")):
if not name.endswith(suffix):
continue
m = _SUBTILE_RE.search(name[:-len(suffix)] + ".avif")
@@ -810,7 +815,9 @@ def render_tile(output_dir, layer, z, x, y, scale=1):
canvas = Image.new("RGBA", (size, size), (0, 0, 0, 0))
# Agrandissement (zoom natif) : bicubique ; réduction : bilinéaire suffit
# puisque le palier source est déjà calé sur la résolution de la tuile.
- resample = Image.BICUBIC
+ # Aplats de niveaux codés (densité, 1 m/px) : plus proche voisin, sinon
+ # l'agrandissement aux zooms 18–19 invente des gris intermédiaires.
+ resample = Image.NEAREST if layer in NEAREST_LAYERS else Image.BICUBIC
painted = False
for src in sources:
painted |= _paste_source(canvas, src, z, x, y, size, resample)
diff --git a/lidar_pipeline/visualizations.py b/lidar_pipeline/visualizations.py
index 5dfb823..d69161e 100644
--- a/lidar_pipeline/visualizations.py
+++ b/lidar_pipeline/visualizations.py
@@ -1129,6 +1129,48 @@ def generate_relief_oriente(dem_file, basename, vis_dir, resolution, shared=None
return None
+# Densité de points sol : 16 niveaux sur une échelle log fixe (même gris =
+# même densité partout, mosaïque jointive). Le niveau k commence à
+# DENSITY_BASE × 2^(k/2) pts/m² : 0 = ≤ 0,35 (noir, y compris sans aucun
+# point), 15 = ≥ 45 (blanc) ; deux niveaux = densité doublée.
+DENSITY_LEVELS = 16
+DENSITY_BASE = 0.25
+
+
+def density_levels(density):
+ """Niveau 0..15 d'une densité (pts/m²)."""
+ with np.errstate(divide="ignore", invalid="ignore"):
+ k = np.floor(2.0 * np.log2(np.asarray(density, dtype=np.float64) / DENSITY_BASE))
+ return np.clip(np.nan_to_num(k, nan=0.0, neginf=0.0), 0, DENSITY_LEVELS - 1).astype(np.uint8)
+
+
+def generate_densite_sol(dem_file, basename, vis_dir, resolution, shared=None):
+ """Densité des points sol retenus pour le MNT, en 16 niveaux de gris.
+
+ Lue dans le fichier annexe écrit à la rastérisation (dtm.density_path) :
+ grille de 1 m, gardée telle quelle (la densité n'a pas plus de détail —
+ image 25× plus légère qu'à 0,2 m). Sortie : niveau 0..15 (float32).
+ """
+ logger.info(" → Densité de points sol...")
+ t0 = time.time()
+ output = vis_dir / f"{basename}_densite_sol.tif"
+ try:
+ from .dtm import density_path
+ src_path = density_path(dem_file)
+ if not src_path.exists():
+ logger.error(f" ✗ Densité absente ({src_path.name}) : régénérer le DTM")
+ return None
+ with rasterio.open(src_path) as src:
+ density = src.read(1)
+ transform, crs = src.transform, src.crs
+ _save_tif(output, density_levels(density).astype(np.float32), transform, crs)
+ logger.info(f" ✓ Densité de points terminée ({time.time()-t0:.1f}s)")
+ return output
+ except Exception as e:
+ logger.error(f" ✗ Erreur densité de points: {e}", exc_info=True)
+ return None
+
+
def generate_mslrm(dem_file, basename, vis_dir, resolution, shared=None):
"""Multi-Scale Relief Model (MSRM) - LRM at adaptive scales combined (GPU if available).