Improve visualizations: adaptive scales, revert z-score to std normalization

- MSRM/TPI/roughness/anomalies: revert z-score (x-mean)/std to std normalization x/std
  to preserve contrast and visibility of linear features (paths, ditches, trenches)
- MSRM: adaptive scales based on resolution, archaeological weight combination
- TPI: extend from 2 to 4 scales (3m/15m/50m/200m) with weighted combination
- Hillshade: 8 directions instead of 4, altitude 35° instead of 30°
- LRM: adaptive sigma based on resolution
- Openness: doubled radius (100m instead of 50m)
- Roughness: multi-scale (3m fine + 15m broad) instead of single 5x5 window
- Anomalies: uses MSRM multi-scale relief instead of single LRM 15m
- Wavelet: 8 adaptive scales, std normalization, archaeological weights
- Remove svf (Sky-View Factor) and local_dominance visualizations
- Add AVIF format support (default), quality 98
- Add multi-resolution support (-r 0.5,0.2)
- Improve Ctrl+C handling for immediate process termination
- Update rendering.py descriptions for all modified visualizations

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

.gitignore

@@ -45,3 +45,4 @@ htmlcov/
 
 # Éventuels fichiers de cache matplotlib
 matplotlibrc
+.playwright-mcp/

CLAUDE.md

@@ -4,7 +4,7 @@ This file provides guidance to Claude Code (claude.ai/code) when working with co
 
 ## Project Overview
 
-LiDAR archaeological processing pipeline that generates 17 terrain visualizations from LAZ/LAS point clouds. Runs in Docker with optional NVIDIA GPU acceleration (CuPy). Designed for French LiDAR HD data in Lambert 93 (EPSG:2154).
+LiDAR archaeological processing pipeline that generates 16 terrain visualizations from LAZ/LAS point clouds. Runs in Docker with optional NVIDIA GPU acceleration (CuPy). Designed for French LiDAR HD data in Lambert 93 (EPSG:2154).
 
 ## Commands
 
@@ -14,12 +14,15 @@ All commands run inside Docker. Use `./run.sh` as the primary interface.
 ./run.sh -g                                    # Standard run with GPU
 ./run.sh -g -w 4                               # GPU + 4 parallel workers
 ./run.sh -g -r 0.2                             # High resolution (0.2m/px)
+./run.sh -g -r 0.5,0.2                         # Multi-resolution (0.5m + 0.2m)
 ./run.sh --test                                 # Run unit tests
 ./run.sh -g --file LHD_FXX_1000_6882_PTS_LAMB93_IGN69.copc  # Single file
-./run.sh --ground-classification pmf            # Force PMF ground classification
-./run.sh -g --keep-tif                          # Keep intermediate TIFF files
-./run.sh -g --no-viewer                         # Skip web viewer generation
-./run.sh serve                                  # Start web map server
+./run.sh --ground-classification csf            # Force CSF ground classification (complex terrain)
+./run.sh -g --keep-tif                          # Keep TIFF files (allows WebP regeneration without recalculating DTM)
+./run.sh -g --only hillshade svf lrm            # Only generate specific visualizations
+./run.sh -g --skip ortho topo                   # Exclude specific visualizations
+./run.sh -g --quality 90                       # WebP quality 90 (default: 85)
+./run.sh -g --lossless                          # Lossless WebP compression
 ./run.sh                                        # Print help (no args)
 ```
 
@@ -34,34 +37,50 @@ docker run --rm --gpus all -v $(pwd)/input:/data/input:ro -v $(pwd)/output:/data
 ### Module responsibilities
 
 - **`cli.py`** — argparse + logging setup. Entry point via `python -m lidar_pipeline`.
-- **`pipeline.py`** — `LidarArchaeoPipeline` orchestrator. `VIZ_STEPS` registry maps names to generate functions. `FilePrefixFilter` for parallel logging.
-- **`dtm.py`** — PDAL ground classification (SMRF/PMF/CSF + auto-detection) and DTM generation via scipy `binned_statistic_2d`.
-- **`visualizations.py`** — 15 `generate_*` functions + 2 IGN overlay lambdas. All take `(dem_file, basename, vis_dir, resolution)` and return a TIF path or None.
+- **`pipeline.py`** — `LidarArchaeoPipeline` orchestrator. `VIZ_STEPS` registry maps names to generate functions. `FilePrefixFilter` for parallel logging. Creates `SharedDEM` once per file and passes it to all visualizations. Multi-resolution support: `self.resolutions` list, `_res_suffix()` for naming, `generate_all_visualizations()` accepts `vis_dir` override.
+- **`dtm.py`** — PDAL ground classification (SMRF/CSF + auto-detection) and DTM generation via scipy `binned_statistic_2d`. `create_dtm_fast()` accepts `output_suffix` for multi-resolution DTM naming.
+- **`visualizations.py`** — 13 `generate_*` functions + 2 IGN overlay lambdas. All take `(dem_file, basename, vis_dir, resolution, shared=None)` and return a TIF path or None. `SharedDEM` class pre-computes gradient, NaN mask, LRM to avoid redundant I/O and computation. Lazy evaluation: properties computed on first access.
 - **`gpu.py`** — CuPy/numpy abstraction: `HAS_GPU`, `to_gpu()`, `to_cpu()`, `xp_gaussian_filter()`, `xp_uniform_filter()`, `xp_minimum_filter()`, `gpu_cleanup()`. Falls back to CPU gracefully.
 - **`ign.py`** — IGN WMTS tile download + overlay generation for orthophoto and topographic maps.
-- **`rendering.py`** — `COLORMAPS` dict maps filename keywords to (cmap, title, legend, description). `tif_to_png()` converts TIF→WebP with legend/scale/north arrow. `convert_to_cog()` converts TIF→Cloud Optimized GeoTIFF. `generate_cog_metadata()` creates metadata JSON for web viewer. `generate_pdf_report()` creates A3 PDF.
-- **`viewer.py`** — Generates MapLibre GL JS HTML viewer with layer controls, opacity sliders, and IGN/OSM basemaps.
-- **`server.py`** — TiTiler-based Starlette server for serving COG tiles and viewer HTML. Entry point via `python -m lidar_pipeline.server`.
+- **`rendering.py`** — `COLORMAPS` dict maps filename keywords to (cmap, title, legend, description). `tif_to_png()` converts TIF→WebP with legend/scale/north arrow. Quality parameter controls WebP compression (default 85).
+
+### SharedDEM optimization
+
+`SharedDEM` pre-computes once per file:
+- DEM data (single I/O read)
+- NaN mask + filled DEM (single `_fill_nans` call, avoiding ~20 redundant calls)
+- Gradient components (dy, dx, slope, aspect) shared by hillshade, slope, aspect, curvature
+- LRM at 15m kernel (shared by lrm + anomalies)
+
+`_filter_nanaware_from_filled()` applies filters on the pre-filled DEM, skipping the expensive `_fill_nans` interpolation.
 
 ### Adding a visualization
 
 Three places must be updated:
-1. `visualizations.py` — add `generate_X(dem_file, basename, vis_dir, resolution)` function
+1. `visualizations.py` — add `generate_X(dem_file, basename, vis_dir, resolution, shared=None)` function
 2. `pipeline.py` `VIZ_STEPS` — add `('name', generate_X)` entry
 3. `rendering.py` `COLORMAPS` — add entry keyed by the output filename keyword
 
 ### Ground classification
 
 Auto-detection in `dtm.py` `detect_ground_method()`:
-- Single-return ratio > 0.6 → PMF (urban terrain)
+- Single-return ratio > 0.6 → CSF (urban terrain, cloth simulation)
 - Height std > 30m → CSF (complex/mountainous terrain)
 - Default → SMRF (natural terrain)
 
-Override with `--ground-classification {auto,smrf,pmf,csf}`.
+Override with `--ground-classification {auto,smrf,csf}`.
 
 ### NaN handling
 
-DTM small gaps (< 1m from existing data) are filled using `rasterio.fill.fillnodata`. Large gaps remain as NaN. Visualization functions use `_fill_nans()` and `_filter_nanaware()` to avoid NaN propagation through filters.
+DTM small gaps (< 1m from existing data) are filled using `rasterio.fill.fillnodata`. Large gaps remain as NaN. `SharedDEM` fills NaN once; `_filter_nanaware_from_filled()` applies filters on the pre-filled array and restores the NaN mask.
+
+### Flow accumulation
+
+Uses priority-flood algorithm (Wang & Liu 2006) for sink filling, which is O(n log n) instead of iterative minimum_filter. D8 accumulation uses numba JIT; falls back to pure Python if numba unavailable.
+
+### Multi-resolution
+
+`-r 0.5,0.2` processes each tile at both 0.5m and 0.2m. Ground classification is shared (done once per tile). Each resolution gets its own DTM (`_dtm.tif` / `_dtm_r0p2.tif`) and visualization subdirectory (`basename/` / `basename_r0p2/`).
 
 ### Parallel processing
 
@@ -72,7 +91,6 @@ Uses `ProcessPoolExecutor` with `'spawn'` start method (required for CUDA). Each
 - **Language**: UI messages and comments in French. Code identifiers in English.
 - **Logging**: Use `logger = logging.getLogger("lidar")`. Prefix per-file logs via `_file_filter.basename`.
 - **GPU pattern**: `arr_gpu = to_gpu(arr)` → compute → `result = to_cpu(arr_gpu)` → `gpu_cleanup()` between visualizations.
-- **Output format**: Visualizations saved as WebP (not PNG). TIFF intermediates deleted unless `--keep-tif` or viewer enabled. COGs generated for web viewer by default. PDF reports use `PILImage.open().convert('RGB')`.
-- **Web viewer**: MapLibre GL JS + TiTiler. COGs served as raster tiles. `./run.sh serve` starts server on port 8000.
-- **Flow accumulation**: Uses numba JIT for D8 accumulation loop. Falls back to pure Python if numba unavailable.
+- **Output format**: Visualizations saved as AVIF (quality 98 by default, best quality/size ratio). Use `--format webp` for WebP output. TIFF intermediates deleted by default. Use `--keep-tif` to keep DTM+TIF for regeneration with `--force`. No PDF reports, no COGs or viewer.
+- **Compression**: TIF intermediates use `deflate` compression (faster than LZW for float32 data).
 - **Tests**: Run only inside Docker via `./run.sh --test`. Synthetic DEM fixture in `tests/conftest.py`.

Dockerfile

@@ -3,17 +3,21 @@ FROM nvidia/cuda:12.4.0-devel-ubuntu22.04
 ENV DEBIAN_FRONTEND=noninteractive
 ENV TZ=Europe/Paris
 
-# Install PDAL and system packages
+# Install system packages + Miniforge for PDAL >= 2.5 (Ubuntu 22.04 ships PDAL 2.3 which can't read COPC v1.1)
 RUN apt-get update && apt-get install -y --no-install-recommends \
-    pdal \
-    liblaszip8 \
     gdal-bin \
     python3-gdal \
     python3-pip \
     python3-dev \
     build-essential \
     wget \
-    && rm -rf /var/lib/apt/lists/*
+    && rm -rf /var/lib/apt/lists/* \
+    && wget -q https://github.com/conda-forge/miniforge/releases/latest/download/Miniforge3-Linux-x86_64.sh -O /tmp/miniforge.sh \
+    && bash /tmp/miniforge.sh -b -p /opt/conda \
+    && rm /tmp/miniforge.sh \
+    && /opt/conda/bin/conda install -y -c conda-forge pdal \
+    && ln -sf /opt/conda/bin/pdal /usr/local/bin/pdal \
+    && /opt/conda/bin/conda clean -afy
 
 WORKDIR /app
 
@@ -37,7 +41,8 @@ RUN pip3 install --no-cache-dir \
     titiler.core \
     fastapi \
     uvicorn \
-    piexif
+    piexif \
+    pillow-avif-plugin
 
 # Install CuPy for GPU acceleration (optional - will fallback to numpy if not available)
 RUN pip3 install --no-cache-dir cupy-cuda12x || echo "CuPy not available - GPU acceleration disabled"

download_lidar.sh

@@ -0,0 +1,584 @@
+#!/bin/bash
+# download_lidar.sh — Téléchargement de fichiers LiDAR HD depuis l'IGN
+#
+# Usage:
+#   ./download_lidar.sh LHD_FXX_1049_6895_PTS_LAMB93_IGN69.copc.laz
+#   ./download_lidar.sh 1049_6895
+#   ./download_lidar.sh 1049_6895 1029_6884
+#   ./download_lidar.sh --list-zones
+#   ./download_lidar.sh --search 1049
+#   ./download_lidar.sh --zone RE --output input/
+#   ./download_lidar.sh --zone RE --dry-run
+#
+# L'API IGN expose les fichiers organisés par zones (RE, SE, AE, etc.)
+# Le script cherche automatiquement dans quelle zone se trouve chaque tuile.
+
+set -euo pipefail
+
+BASE_URL="https://data.geopf.fr/telechargement/download/LiDARHD-NUALID"
+API_URL="https://data.geopf.fr/telechargement/resource/LiDARHD-NUALID"
+OUTPUT_DIR="${OUTPUT_DIR:-.}"
+
+# Colors
+RED='\033[0;31m'
+GREEN='\033[0;32m'
+YELLOW='\033[1;33m'
+CYAN='\033[0;36m'
+NC='\033[0m'
+
+log_info()  { echo -e "${CYAN}[INFO]${NC} $*"; }
+log_ok()    { echo -e "${GREEN}[OK]${NC} $*"; }
+log_warn()  { echo -e "${YELLOW}[WARN]${NC} $*"; }
+log_err()   { echo -e "${RED}[ERR]${NC} $*"; }
+
+# Normalize a tile identifier to a full filename
+# Supports: "1049_6895", "LHD_FXX_1049_6895_PTS_LAMB93_IGN69", "LHD_FXX_1049_6895_PTS_LAMB93_IGN69.copc.laz"
+normalize_name() {
+    local input="$1"
+    # If it's already a full filename, return as-is
+    if [[ "$input" == LHD_FXX_*_PTS_LAMB93_IGN69.copc.laz ]]; then
+        echo "$input"
+        return
+    fi
+    if [[ "$input" == LHD_FXX_*_PTS_LAMB93_IGN69 ]]; then
+        echo "${input}.copc.laz"
+        return
+    fi
+    # Bare coordinates like "1049_6895"
+    echo "LHD_FXX_${input}_PTS_LAMB93_IGN69.copc.laz"
+}
+
+# List all available zone releases
+list_zones() {
+    log_info "Récupération de la liste des zones..."
+    local page=1
+    local total_pages=""
+
+    while true; do
+        local content
+        content=$(curl -s "${API_URL}?page=${page}&limit=50") || true
+
+        if [ -z "$content" ]; then
+            break
+        fi
+
+        # Parse zone codes and dates
+        echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    for entry in root.findall('atom:entry', ns):
+        title = entry.find('atom:title', ns)
+        if title is not None and title.text:
+            # Extract zone code and date
+            parts = title.text.split('_')
+            # NUALHD_1-0__LAZ_LAMB93_{ZONE}_{DATE}
+            zone = parts[-2] if len(parts) >= 2 else '?'
+            date = parts[-1] if len(parts) >= 1 else '?'
+            print(f'{zone}\t{date}\t{title.text}')
+except Exception:
+    pass
+" 2>/dev/null || true
+
+        # Check if more pages
+        if [ -z "$total_pages" ]; then
+            total_pages=$(echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom', 'gpf_dl': 'https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    pc = root.get('{https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd}pagecount', '1')
+    print(pc)
+except:
+    print('1')
+" 2>/dev/null)
+        fi
+
+        if [ "$page" -ge "${total_pages:-1}" ]; then
+            break
+        fi
+        page=$((page + 1))
+    done | sort
+}
+
+# Search for a tile across all zones
+search_tile() {
+    local tile_name="$1"
+    local filename
+    filename=$(normalize_name "$tile_name")
+
+    log_info "Recherche de ${filename} dans toutes les zones..."
+
+    local page=1
+    local found=0
+    local total_pages=""
+
+    while true; do
+        local content
+        content=$(curl -s "${API_URL}?page=${page}&limit=50") || true
+        [ -z "$content" ] && break
+
+        # Get sub-resources and search for the file
+        echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    for entry in root.findall('atom:entry', ns):
+        title = entry.find('atom:title', ns)
+        if title is not None and title.text:
+            print(title.text)
+except:
+    pass
+" 2>/dev/null | while read -r release; do
+            # Check if file exists in this release
+            local url="${BASE_URL}/${release}/${filename}"
+            local code
+            code=$(curl -sI -L -o /dev/null -w "%{http_code}" "$url" 2>/dev/null)
+            if [ "$code" = "200" ]; then
+                log_ok "Trouvé: ${filename} dans la zone ${release}"
+                echo "  URL: ${url}"
+                found=1
+            fi
+        done
+
+        if [ -z "$total_pages" ]; then
+            total_pages=$(echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom', 'gpf_dl': 'https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    pc = root.get('{https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd}pagecount', '1')
+    print(pc)
+except:
+    print('1')
+" 2>/dev/null)
+        fi
+
+        if [ "$page" -ge "${total_pages:-1}" ]; then
+            break
+        fi
+        page=$((page + 1))
+    done
+
+    if [ "$found" -eq 0 ]; then
+        log_err "${filename} non trouvé dans aucune zone"
+    fi
+}
+
+# Download a single file with resume support
+download_file() {
+    local url="$1"
+    local output_path="$2"
+    local filename
+    filename=$(basename "$output_path")
+
+    # Resume download with retries
+    local max_retries=15
+    local retry=0
+
+    while [ $retry -lt $max_retries ]; do
+        retry=$((retry + 1))
+        local result
+        result=$(curl -L --http1.1 -C - -o "$output_path" "$url" 2>&1)
+        local exit_code=$?
+
+        if [ $exit_code -eq 0 ]; then
+            # Verify it's a valid LAS/LAZ file
+            local magic
+            magic=$(head -c 4 "$output_path" 2>/dev/null | xxd -p 2>/dev/null || echo "")
+            if [ "$magic" = "4c415346" ]; then
+                local size
+                size=$(stat -c%s "$output_path" 2>/dev/null || echo "0")
+                log_ok "${filename} téléchargé (${size} octets)"
+                return 0
+            else
+                log_err "${filename} n'est pas un fichier LAZ valide (magic: ${magic})"
+                rm -f "$output_path"
+                return 1
+            fi
+        fi
+
+        # curl exit code 18 = partial download, resume
+        if [ $exit_code -eq 18 ]; then
+            local size
+            size=$(stat -c%s "$output_path" 2>/dev/null || echo "0")
+            log_warn "Téléchargement partiel (${size} octets), tentative ${retry}/${max_retries}..."
+            sleep 2
+            continue
+        fi
+
+        # HTTP 404
+        if echo "$result" | grep -q "404"; then
+            log_err "Fichier non trouvé (404): ${url}"
+            rm -f "$output_path"
+            return 1
+        fi
+
+        log_err "Erreur curl (${exit_code}): ${result}"
+        return 1
+    done
+
+    # If we exhausted retries, check if partial file is usable
+    local magic
+    magic=$(head -c 4 "$output_path" 2>/dev/null | xxd -p 2>/dev/null || echo "")
+    if [ "$magic" = "4c415346" ]; then
+        local size
+        size=$(stat -c%s "$output_path" 2>/dev/null || echo "0")
+        log_warn "Fichier partiel mais valide: ${filename} (${size} octets)"
+        return 0
+    fi
+
+    log_err "Échec après ${max_retries} tentatives"
+    rm -f "$output_path"
+    return 1
+}
+
+# Find which zone contains a given tile and download it
+find_and_download() {
+    local tile_name="$1"
+    local output_dir="$2"
+    local filename
+    filename=$(normalize_name "$tile_name")
+    local output_path="${output_dir}/${filename}"
+
+    # If file already exists and is valid, skip
+    if [ -f "$output_path" ]; then
+        local magic
+        magic=$(head -c 4 "$output_path" 2>/dev/null | xxd -p 2>/dev/null || echo "")
+        if [ "$magic" = "4c415346" ]; then
+            local size
+            size=$(stat -c%s "$output_path" 2>/dev/null || echo "0")
+            if [ "$size" -gt 1000000 ]; then
+                log_info "${filename} déjà présent (${size} octets) — ignoré"
+                return 0
+            fi
+        fi
+        # File exists but seems invalid/truncated, re-download
+        log_warn "${filename} existe mais semble incomplet, re-téléchargement..."
+    fi
+
+    log_info "Recherche de ${filename}..."
+
+    # Try all zone releases
+    local page=1
+    local total_pages=""
+
+    while true; do
+        local content
+        content=$(curl -s "${API_URL}?page=${page}&limit=50") || true
+        [ -z "$content" ] && break
+
+        echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    for entry in root.findall('atom:entry', ns):
+        title = entry.find('atom:title', ns)
+        if title is not None and title.text:
+            print(title.text)
+except:
+    pass
+" 2>/dev/null | while read -r release; do
+            local url="${BASE_URL}/${release}/${filename}"
+            local code
+            code=$(curl -sI -L -o /dev/null -w "%{http_code}" "$url" 2>/dev/null || echo "000")
+            if [ "$code" = "200" ]; then
+                log_ok "Trouvé dans la zone ${release}"
+                download_file "$url" "$output_path"
+                return 0
+            fi
+        done
+
+        if [ -z "$total_pages" ]; then
+            total_pages=$(echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom', 'gpf_dl': 'https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    pc = root.get('{https://data.geopf.fr/annexes/ressources/xsd/gpf_dl.xsd}pagecount', '1')
+    print(pc)
+except:
+    print('1')
+" 2>/dev/null)
+        fi
+
+        if [ "$page" -ge "${total_pages:-1}" ]; then
+            break
+        fi
+        page=$((page + 1))
+    done
+
+    log_err "${filename} non trouvé dans aucune zone disponible"
+    return 1
+}
+
+# Download all tiles from a specific zone
+download_zone() {
+    local zone="$1"
+    local output_dir="$2"
+    local dry_run="${3:-false}"
+
+    # Find the release name for this zone
+    local release=""
+    local page=1
+
+    while true; do
+        local content
+        content=$(curl -s "${API_URL}?page=${page}&limit=50") || true
+        [ -z "$content" ] && break
+
+        release=$(echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    for entry in root.findall('atom:entry', ns):
+        title = entry.find('atom:title', ns)
+        if title is not None and title.text:
+            # Extract zone code from title like NUALHD_1-0__LAZ_LAMB93_RE_2025-02-17
+            parts = title.text.split('_')
+            z = parts[-2] if len(parts) >= 2 else ''
+            if z == '${zone}':
+                print(title.text)
+                break
+except:
+    pass
+" 2>/dev/null)
+
+        if [ -n "$release" ]; then
+            break
+        fi
+        page=$((page + 1))
+        # Safety limit
+        if [ $page -gt 25 ]; then
+            break
+        fi
+    done
+
+    if [ -z "$release" ]; then
+        log_err "Zone '${zone}' non trouvée"
+        return 1
+    fi
+
+    log_info "Zone ${zone}: ${release}"
+
+    # Get list of files in this release
+    local sub_url="${API_URL}/${release}?limit=50"
+    local total_files=0
+    local page=1
+
+    while true; do
+        local content
+        content=$(curl -s "${API_URL}/${release}?page=${page}&limit=50") || true
+        [ -z "$content" ] && break
+
+        local files
+        files=$(echo "$content" | python3 -c "
+import sys, xml.etree.ElementTree as ET
+ns = {'atom': 'http://www.w3.org/2005/Atom'}
+try:
+    tree = ET.parse(sys.stdin)
+    root = tree.getroot()
+    for entry in root.findall('atom:entry', ns):
+        title = entry.find('atom:title', ns)
+        if title is not None and title.text:
+            print(title.text)
+except:
+    pass
+" 2>/dev/null)
+
+        if [ -z "$files" ]; then
+            break
+        fi
+
+        local count
+        count=$(echo "$files" | wc -l)
+        total_files=$((total_files + count))
+
+        if [ "$dry_run" = "true" ]; then
+            echo "$files" | head -10
+            echo "... (${count} fichiers sur cette page)"
+        else
+            echo "$files" | while read -r fname; do
+                local url="${BASE_URL}/${release}/${fname}"
+                download_file "$url" "${output_dir}/${fname}"
+            done
+        fi
+
+        page=$((page + 1))
+        # Safety limit
+        if [ $page -gt 500 ]; then
+            break
+        fi
+    done
+
+    log_info "Total: ${total_files} fichiers dans la zone ${zone}"
+}
+
+# Show usage
+usage() {
+    cat <<'EOF'
+download_lidar.sh — Téléchargement de fichiers LiDAR HD depuis l'IGN
+
+Usage:
+  ./download_lidar.sh <tuile> [tuile2...]    Télécharger une ou plusieurs tuiles
+  ./download_lidar.sh --search <motif>        Chercher une tuile dans toutes les zones
+  ./download_lidar.sh --list-zones            Lister les zones disponibles
+  ./download_lidar.sh --zone <ZONE> [--dry-run]  Télécharger une zone entière
+  ./download_lidar.sh --check [dossier/]      Vérifier l'intégrité des fichiers LAZ
+
+Formats de tuile acceptés:
+  1049_6895                                    Coordonnées seules
+  LHD_FXX_1049_6895_PTS_LAMB93_IGN69          Nom complet sans extension
+  LHD_FXX_1049_6895_PTS_LAMB93_IGN69.copc.laz Nom complet avec extension
+
+Options:
+  -o, --output DIR    Répertoire de destination (défaut: .)
+  -n, --dry-run       Afficher sans télécharger (avec --zone)
+  --check [DIR]       Vérifier l'intégrité des fichiers LAZ (défaut: ./)
+  --list-zones        Lister toutes les zones disponibles
+  --search MOTIF      Chercher un motif dans les noms de tuile
+  --zone ZONE         Télécharger tous les fichiers d'une zone
+  -h, --help          Afficher cette aide
+
+Exemples:
+  ./download_lidar.sh 1049_6895
+  ./download_lidar.sh LHD_FXX_0713_6347_PTS_LAMB93_IGN69.copc.laz -o input/
+  ./download_lidar.sh --search 1049
+  ./download_lidar.sh --list-zones
+  ./download_lidar.sh --zone RE -o input/
+  ./download_lidar.sh --check input/
+EOF
+}
+
+# Check integrity of LAZ files
+check_files() {
+    local dir="${1:-.}"
+    local ok=0
+    local fail=0
+
+    log_info "Vérification de l'intégrité des fichiers LAZ dans ${dir}/"
+
+    for f in "${dir}"/*.copc.laz "${dir}"/*.laz; do
+        [ -f "$f" ] || continue
+        local magic
+        magic=$(head -c 4 "$f" 2>/dev/null | xxd -p 2>/dev/null || echo "")
+        local size
+        size=$(stat -c%s "$f" 2>/dev/null || echo "0")
+
+        if [ "$magic" = "4c415346" ] && [ "$size" -gt 1000000 ]; then
+            log_ok "$(basename $f) (${size} octets)"
+            ok=$((ok + 1))
+        elif [ "$magic" = "7b227469" ]; then
+            log_err "$(basename $f) est une page HTML (404), pas un fichier LAZ"
+            fail=$((fail + 1))
+        elif [ "$size" -lt 1000 ]; then
+            log_err "$(basename $f) trop petit (${size} octets) — probablement corrompu"
+            fail=$((fail + 1))
+        else
+            log_warn "$(basename $f) magic=${magic}, taille=${size} — vérification nécessaire"
+            fail=$((fail + 1))
+        fi
+    done
+
+    log_info "Résultat: ${ok} OK, ${fail} problème(s)"
+    return 0
+}
+
+# Main
+CMD=""
+TILES=()
+ZONE=""
+DRY_RUN=false
+
+while [ $# -gt 0 ]; do
+    case "$1" in
+        -h|--help)
+            usage
+            exit 0
+            ;;
+        --list-zones)
+            CMD="list-zones"
+            shift
+            ;;
+        --search)
+            CMD="search"
+            shift
+            TILES+=("$1")
+            shift
+            ;;
+        --zone)
+            CMD="zone"
+            shift
+            ZONE="$1"
+            shift
+            ;;
+        --check)
+            CMD="check"
+            shift
+            OUTPUT_DIR="${1:-.}"
+            [ -d "$OUTPUT_DIR" ] && shift || true
+            ;;
+        -o|--output)
+            shift
+            OUTPUT_DIR="$1"
+            shift
+            ;;
+        -n|--dry-run)
+            DRY_RUN=true
+            shift
+            ;;
+        -*)
+            log_err "Option inconnue: $1"
+            usage
+            exit 1
+            ;;
+        *)
+            TILES+=("$1")
+            shift
+            ;;
+    esac
+done
+
+case "${CMD}" in
+    list-zones)
+        list_zones
+        ;;
+    search)
+        for tile in "${TILES[@]}"; do
+            search_tile "$tile"
+        done
+        ;;
+    zone)
+        if [ -z "$ZONE" ]; then
+            log_err "Zone non spécifiée. Utilisez --zone <CODE>"
+            exit 1
+        fi
+        mkdir -p "$OUTPUT_DIR"
+        download_zone "$ZONE" "$OUTPUT_DIR" "$DRY_RUN"
+        ;;
+    check)
+        check_files "$OUTPUT_DIR"
+        ;;
+    *)
+        # Default: download tiles
+        if [ ${#TILES[@]} -eq 0 ]; then
+            usage
+            exit 1
+        fi
+        mkdir -p "$OUTPUT_DIR"
+        for tile in "${TILES[@]}"; do
+            find_and_download "$tile" "$OUTPUT_DIR"
+        done
+        ;;
+esac

lidar_pipeline/cli.py

@@ -23,6 +23,12 @@ def setup_logging(verbose=False, debug=False):
         verbose: If True, include timestamps and level names.
         debug:    If True, set level to DEBUG and add file:line info.
     """
+    # Ensure UTF-8 output for French messages (é, è, ê, etc.)
+    if hasattr(sys.stdout, 'reconfigure'):
+        sys.stdout.reconfigure(encoding='utf-8', errors='replace')
+    if hasattr(sys.stderr, 'reconfigure'):
+        sys.stderr.reconfigure(encoding='utf-8', errors='replace')
+
     if debug:
         level = logging.DEBUG
         fmt = "%(asctime)s.%(msecs)03d %(levelname)-5s [%(filename)s:%(lineno)d] %(message)s"
@@ -91,9 +97,9 @@ Exemples:
     )
     parser.add_argument(
         "-r", "--resolution",
-        type=float,
-        default=0.5,
-        help="Résolution en mètres par pixel (défaut: 0.5)"
+        type=str,
+        default="0.5",
+        help="Résolution en m/px, ou multiples séparées par virgules (défaut: 0.5, ex: 0.5,0.2)"
     )
     parser.add_argument(
         "-w", "--workers",
@@ -114,18 +120,44 @@ Exemples:
     parser.add_argument(
         "--keep-tif",
         action="store_true",
-        help="Conserver les fichiers TIFF intermédiaires (sinon supprimés après conversion WebP)"
-    )
-    parser.add_argument(
-        "--no-viewer",
-        action="store_true",
-        help="Ne pas générer le viewer web (COGs + HTML MapLibre)"
+        help="Conserver les fichiers TIFF (DTM + visualisations) pour pouvoir régénérer les WebP sans recalculer"
     )
     parser.add_argument(
         "--ground-classification",
-        choices=["auto", "smrf", "pmf", "csf"],
+        choices=["auto", "smrf", "csf"],
         default="auto",
-        help="Méthode de classification du sol : auto (détection), smrf, pmf, csf (défaut: auto)"
+        help="Méthode de classification du sol : auto (détection), smrf, csf (défaut: auto)"
+    )
+    parser.add_argument(
+        "--quality",
+        type=int,
+        default=98,
+        help="Qualité image (1-100, défaut: 98). Utilisez 100 pour lossless."
+    )
+    parser.add_argument(
+        "--lossless",
+        action="store_true",
+        help="Forcer la compression lossless (équivalent à --quality 100)"
+    )
+    parser.add_argument(
+        "--format",
+        choices=["webp", "avif"],
+        default="avif",
+        help="Format de sortie : avif (défaut, meilleure qualité) ou webp"
+    )
+    parser.add_argument(
+        "--only",
+        nargs="+",
+        type=str,
+        default=None,
+        help="Générer uniquement ces visualisations (ex: --only hillshade svf lrm)"
+    )
+    parser.add_argument(
+        "--skip",
+        nargs="+",
+        type=str,
+        default=None,
+        help="Exclure ces visualisations (ex: --skip ortho topo)"
     )
     parser.add_argument(
         "--file",
@@ -167,6 +199,14 @@ Exemples:
     log_gpu_status()
 
     try:
+        quality = 100 if args.lossless else args.quality
+        # Parse --only and --skip: accept comma-separated values
+        only_viz = None
+        if args.only:
+            only_viz = [v.strip() for item in args.only for v in item.split(',')]
+        skip_viz = None
+        if args.skip:
+            skip_viz = [v.strip() for item in args.skip for v in item.split(',')]
         pipeline = LidarArchaeoPipeline(
             input_dir=args.input,
             output_dir=args.output,
@@ -176,7 +216,10 @@ Exemples:
             ground_method=args.ground_classification,
             force_classify=args.force_classification,
             keep_tif=args.keep_tif,
-            no_viewer=args.no_viewer
+            quality=quality,
+            only_viz=only_viz,
+            skip_viz=skip_viz,
+            output_format=args.format,
         )
 
         # If --file is specified, process only matching files
@@ -187,10 +230,16 @@ Exemples:
             # Also supports bare name without .copc (e.g. LHD_FXX_1000_6882_PTS_LAMB93_IGN69)
             selected_files = []
             for pattern in args.file:
-                matches = (list(input_dir.glob(f"{pattern}.laz"))
-                          + list(input_dir.glob(f"{pattern}.las"))
-                          + list(input_dir.glob(f"{pattern}.copc.laz"))
-                          + list(input_dir.glob(f"{pattern}.copc.las")))
+                # Try exact filename first (e.g. LHD_FXX_...copc.laz)
+                exact_match = input_dir / pattern
+                if exact_match.exists() and exact_match.is_file():
+                    matches = [exact_match]
+                else:
+                    # Try with added extensions (e.g. pattern=LHD_FXX_...IGN69)
+                    matches = (list(input_dir.glob(f"{pattern}.laz"))
+                              + list(input_dir.glob(f"{pattern}.las"))
+                              + list(input_dir.glob(f"{pattern}.copc.laz"))
+                              + list(input_dir.glob(f"{pattern}.copc.las")))
                 # Remove duplicates
                 matches = list(dict.fromkeys(matches))
                 if not matches:
@@ -235,9 +284,15 @@ def _kill_orphan_pdal(signum=None, frame=None):
     """Kill orphan PDAL processes on interrupt or exit."""
     import subprocess
     try:
-        subprocess.run(["pkill", "-f", "pdal"], capture_output=True, timeout=5)
+        subprocess.run(["pkill", "-9", "-f", "pdal"], capture_output=True, timeout=3)
     except Exception:
         pass
     if signum is not None:
-        logger.info("Interruption — nettoyage des processus PDAL")
+        logger.info("Interruption — nettoyage des processus")
+        # Force-kill all child processes immediately
+        try:
+            import os
+            os.killpg(os.getpgrp(), signal.SIGKILL)
+        except Exception:
+            pass
         sys.exit(130)

lidar_pipeline/dtm.py

@@ -1,6 +1,6 @@
 """DTM generation from classified LiDAR point clouds.
 
-Handles ground classification via PDAL/SMRF and DTM rasterisation
+Handles ground classification via PDAL (SMRF or CSF) and DTM rasterisation
 using scipy binned_statistic_2d. Zones without LiDAR data remain as NaN.
 """
 
@@ -27,13 +27,13 @@ def _create_ground_pipeline(input_laz, output_las, method):
     1. Reset Classification to 0
     2. ELM (Extended Local Minimum) — mark low outliers as noise (Classification=7)
     3. Statistical outlier removal
-    4. Ground classification (SMRF/PMF/CSF)
+    4. Ground classification (SMRF or CSF)
     5. Extract ground points (Classification=2)
 
     Args:
         input_laz: Path to input LAZ/LAS file.
         output_las: Path to output classified LAS file.
-        method: Ground classification method ('smrf', 'pmf', or 'csf').
+        method: Ground classification method ('smrf' or 'csf').
 
     Returns:
         JSON string of the PDAL pipeline.
@@ -84,15 +84,6 @@ def _create_ground_pipeline(input_laz, output_las, method):
             "threshold": 0.5,
             "scalar": 1.25
         }
-    elif method == 'pmf':
-        ground_step = {
-            "type": "filters.pmf",
-            "max_window": 33,
-            "slope": 0.15,
-            "max_distance": 2.5,
-            "initial_distance": 0.15,
-            "cell_size": 1.0
-        }
     elif method == 'csf':
         ground_step = {
             "type": "filters.csf",
@@ -128,22 +119,118 @@ def create_smrf_pipeline(input_laz, output_las):
     return _create_ground_pipeline(input_laz, output_las, 'smrf')
 
 
-def create_pmf_pipeline(input_laz, output_las):
-    """Create a PDAL pipeline JSON for PMF ground classification."""
-    return _create_ground_pipeline(input_laz, output_las, 'pmf')
-
-
 def create_csf_pipeline(input_laz, output_las):
     """Create a PDAL pipeline JSON for CSF ground classification."""
     return _create_ground_pipeline(input_laz, output_las, 'csf')
 
 
+def validate_laz(laz_file):
+    """Quick integrity check for a LAZ/LAS file.
+
+    Tries laspy first (fast header read), then PDAL as fallback for COPC files
+    that laspy cannot read. Also checks that the file contains points.
+
+    Returns:
+        True if file is readable and contains points, False otherwise.
+    """
+    # Try laspy first (fast)
+    import laspy
+    try:
+        with laspy.open(str(laz_file)) as f:
+            header = f.header
+            point_count = header.point_count
+            if point_count == 0:
+                logger.error(f"  ✗ Fichier vide (0 points): {laz_file.name}")
+                logger.error(f"    → Re-télécharger depuis https://ign.fr/lidar-hd")
+                return False
+        return True
+    except Exception:
+        pass
+
+    # Fallback: try PDAL (handles COPC v1.1 that laspy can't read)
+    import subprocess
+    try:
+        result = subprocess.run(
+            ["pdal", "info", str(laz_file), "--summary"],
+            capture_output=True, text=True, timeout=30
+        )
+        if result.returncode == 0:
+            # Check point count from PDAL info output
+            import json as _json
+            try:
+                info = _json.loads(result.stdout)
+                count = info.get('summary', {}).get('num_points', 0)
+                if count == 0:
+                    logger.error(f"  ✗ Fichier vide (0 points PDAL): {laz_file.name}")
+                    logger.error(f"    → Re-télécharger depuis https://ign.fr/lidar-hd")
+                    return False
+            except Exception:
+                pass  # Can't parse — assume valid
+            return True
+        logger.error(f"  ✗ Fichier illisible: {laz_file.name}")
+        logger.error(f"    PDAL: {result.stderr.strip()[:200]}")
+    except (subprocess.TimeoutExpired, FileNotFoundError):
+        logger.error(f"  ✗ Impossible de vérifier le fichier: {laz_file.name}")
+    logger.error(f"    → Re-télécharger depuis https://ign.fr/lidar-hd")
+    return False
+
+
+def _read_with_pdal(laz_file):
+    """Read a LAZ/LAS file via PDAL when laspy fails (e.g. COPC v1.1).
+
+    Returns a laspy.LasData object, or None on failure.
+    """
+    import subprocess
+    import tempfile
+    import os
+
+    try:
+        # Convert COPC to LAS via PDAL, then read with laspy
+        with tempfile.NamedTemporaryFile(suffix='.las', delete=False) as tmp:
+            tmp_path = tmp.name
+
+        pipeline = json.dumps({
+            "pipeline": [
+                str(laz_file),
+                {"type": "writers.las", "filename": tmp_path}
+            ]
+        })
+
+        result = subprocess.run(
+            ["pdal", "pipeline", "--stdin"],
+            input=pipeline, capture_output=True, text=True, timeout=300
+        )
+
+        if result.returncode != 0:
+            logger.warning(f"  PDAL conversion échouée: {result.stderr[:200]}")
+            try:
+                os.unlink(tmp_path)
+            except Exception:
+                pass
+            return None
+
+        import laspy
+        las = laspy.read(tmp_path)
+        try:
+            os.unlink(tmp_path)
+        except Exception:
+            pass
+        if len(las.points) == 0:
+            logger.warning(f"  PDAL: conversion réussie mais 0 points")
+            return None
+        return las
+
+    except Exception as e:
+        logger.warning(f"  PDAL fallback échoué: {e}")
+        return None
+
+
 def detect_ground_method(laz_file):
     """Detect the best ground classification method based on point cloud statistics.
 
-    Auto-selects between SMRF (natural terrain) and PMF (urban) only.
-    CSF is available only via --ground-classification csf (slow but robust
-    on complex terrain).
+    Auto-selects between SMRF and CSF:
+    - SMRF: fast, robust for most natural terrain (PDAL recommended default)
+    - CSF: cloth simulation, better for complex/urban terrain
 
     Falls back to SMRF if the file cannot be read or attributes are missing.
 
@@ -151,18 +238,28 @@ def detect_ground_method(laz_file):
         laz_file: Path to input LAZ/LAS file.
 
     Returns:
-        String: 'smrf', 'pmf', or 'csf'
+        String: 'smrf' or 'csf'
     """
     import laspy
 
+    # Try laspy first, then PDAL for COPC files
+    las = None
     try:
         las = laspy.read(str(laz_file))
     except Exception as e:
-        logger.warning(f"  Impossible de lire le nuage pour auto-détection: {e}")
+        logger.warning(f"  laspy: {e}")
+        logger.info(f"  → Lecture via PDAL pour auto-détection...")
+        las = _read_with_pdal(laz_file)
+
+    if las is None:
         logger.info(f"  → Méthode: SMRF (défaut — lecture impossible)")
         return 'smrf'
 
     total_points = len(las.points)
+    if total_points == 0:
+        logger.warning(f"  Nuage vide (0 points) — méthode par défaut: SMRF")
+        return 'smrf'
+
     z = np.array(las.z)
 
     # Height variance (always available)
@@ -182,13 +279,16 @@ def detect_ground_method(laz_file):
                 f"ratio_retours_uniques={single_return_ratio:.2f}, "
                 f"écart_Z={z_std:.1f}m, amplitude_Z={z_range:.1f}m")
 
-    # Decision logic (auto selects between SMRF and PMF only):
-    # - High single-return ratio (>0.6) → urban (buildings, roads) → PMF
+    # Decision logic:
+    # - High single-return ratio (>0.6) → urban (buildings, roads) → CSF (cloth simulation)
+    # - High elevation variance (>30m) → complex/mountainous terrain → CSF
     # - Default → SMRF (fast, robust for most natural terrain)
-    # Note: CSF is available only via --ground-classification csf (slow but robust on complex terrain)
     if single_return_ratio > 0.6:
-        method = 'pmf'
+        method = 'csf'
         reason = f"ratio retours uniques={single_return_ratio:.2f} > 0.6 → milieu urbain"
+    elif z_std > 30:
+        method = 'csf'
+        reason = f"écart_Z={z_std:.1f}m > 30m → terrain complexe"
     else:
         method = 'smrf'
         reason = f"terrain naturel standard"
@@ -203,7 +303,7 @@ def classify_ground(laz_file, temp_dir, method='auto', force=False):
     Args:
         laz_file: Path to input LAZ/LAS file.
         temp_dir: Directory for temporary files (pipeline.json, ground.las).
-        method: Ground classification method ('auto', 'smrf', 'pmf', or 'csf').
+        method: Ground classification method ('auto', 'smrf', or 'csf').
         force: If True, reclassify even if output file already exists.
 
     Returns:
@@ -243,14 +343,66 @@ def classify_ground(laz_file, temp_dir, method='auto', force=False):
             ["pdal", "pipeline", str(pipeline_file)],
             capture_output=True, check=True
         )
+        # Verify that ground file has points
+        if output_las.exists() and output_las.stat().st_size < 100:
+            logger.error(f"  ✗ Fichier ground vide (taille < 100 octets)")
+            output_las.unlink(missing_ok=True)
+            return None
         logger.info(f"  ✓ Classification sol {method.upper()} terminée")
         return output_las
     except subprocess.CalledProcessError as e:
-        logger.error(f"  ✗ Erreur classification PDAL ({method.upper()}): {e.stderr.decode()}")
+        error_msg = e.stderr.decode() if e.stderr else str(e)
+        logger.warning(f"  ✗ Erreur classification PDAL ({method.upper()}): {error_msg}")
+
+        # Try repairing file with laspy if PDAL fails on EVLR/VLR
+        if 'VLR' in error_msg or 'Invalid' in error_msg:
+            logger.info(f"  → Tentative de réparation du fichier avec laspy...")
+            repaired_las = temp_dir / f"{laz_base}_repaired.las"
+            if _repair_laz_with_laspy(laz_file, repaired_las):
+                # Retry PDAL pipeline with repaired file
+                pipeline_json = _create_ground_pipeline(repaired_las, output_las, method)
+                with open(pipeline_file, 'w') as f:
+                    f.write(pipeline_json)
+                try:
+                    subprocess.run(
+                        ["pdal", "pipeline", str(pipeline_file)],
+                        capture_output=True, check=True
+                    )
+                    logger.info(f"  ✓ Classification sol {method.upper()} terminée (fichier réparé)")
+                    return output_las
+                except subprocess.CalledProcessError as e2:
+                    error_msg2 = e2.stderr.decode() if e2.stderr else str(e2)
+                    logger.error(f"  ✗ Échec classification même après réparation: {error_msg2}")
+            else:
+                logger.error(f"  ✗ Impossible de réparer le fichier")
         return None
 
 
-def create_dtm_fast(las_file, basename, dtm_dir, resolution):
+def _repair_laz_with_laspy(input_laz, output_las):
+    """Try to repair a corrupt LAZ file by re-reading with laspy and saving as LAS.
+
+    Works around PDAL errors like 'Invalid Extended VLR size' by stripping
+    problematic VLR/EVLR metadata during re-save.
+
+    Args:
+        input_laz: Path to corrupt LAZ/LAS file.
+        output_las: Path for repaired LAS output.
+
+    Returns:
+        True if repair succeeded, False otherwise.
+    """
+    import laspy
+    try:
+        las = laspy.read(str(input_laz))
+        las.write(str(output_las))
+        logger.info(f"  ✓ Fichier réparé via laspy ({len(las.points):,} points)")
+        return True
+    except Exception as e:
+        logger.warning(f"  ✗ Réparation laspy échouée: {e}")
+        return False
+
+
+def create_dtm_fast(las_file, basename, dtm_dir, resolution, force=False, output_suffix=""):
     """Create DTM using fast binning method with gap filling.
 
     Args:
@@ -258,16 +410,38 @@ def create_dtm_fast(las_file, basename, dtm_dir, resolution):
         basename: Base name for output file.
         dtm_dir: Directory for output DTM GeoTIFF.
         resolution: Grid resolution in meters per pixel.
+        force: If True, regenerate even if DTM already exists.
+        output_suffix: Suffix for output filename (e.g. '_r0p2' for additional resolutions).
 
     Returns:
         Path to output DTM GeoTIFF, or None on failure.
     """
+    output_tif = dtm_dir / f"{basename}_dtm{output_suffix}.tif"
+
+    if output_tif.exists() and not force:
+        logger.info(f"  DTM déjà existant — fichier réutilisé: {output_tif.name}")
+        return output_tif
+
     import laspy
 
     logger.info("  → Génération DTM...")
 
     try:
         las = laspy.read(str(las_file))
+    except Exception as e:
+        # laspy can't read COPC v1.1 — try PDAL conversion
+        logger.warning(f"  laspy: {e}")
+        logger.info(f"  → Conversion via PDAL pour lecture COPC...")
+        las = _read_with_pdal(las_file)
+        if las is None:
+            logger.error(f"  ✗ Impossible de lire {las_file.name}")
+            return None
+
+    if len(las.points) == 0:
+        logger.error(f"  ✗ Fichier vide (0 points): {las_file.name}")
+        return None
+
+    try:
 
         min_x, max_x = float(las.header.min[0]), float(las.header.max[0])
         min_y, max_y = float(las.header.min[1]), float(las.header.max[1])
@@ -309,7 +483,7 @@ def create_dtm_fast(las_file, basename, dtm_dir, resolution):
             logger.info(f"    {remaining:,} pixels restent sans données (grands écarts)")
 
         # Save as GeoTIFF
-        output_tif = dtm_dir / f"{basename}_dtm.tif"
+        output_tif = dtm_dir / f"{basename}_dtm{output_suffix}.tif"
         transform = from_bounds(min_x, min_y, max_x, max_y, width, height)
 
         with rasterio.open(

lidar_pipeline/gpu.py

@@ -110,6 +110,16 @@ def xp_minimum_filter(arr, footprint=None, size=None):
     return ndimage.minimum_filter(arr, footprint=footprint, size=size)
 
 
+def xp_maximum_filter(arr, footprint=None, size=None):
+    """Maximum filter — uses GPU if array is on GPU, CPU otherwise."""
+    if _cp is not None and isinstance(arr, _cp.ndarray):
+        try:
+            return _cp_ndimage.maximum_filter(arr, footprint=footprint, size=size)
+        except Exception:
+            arr = to_cpu(arr)
+    return ndimage.maximum_filter(arr, footprint=footprint, size=size)
+
+
 def gpu_cleanup():
     """Free GPU memory. Call between visualizations to prevent OOM."""
     if _cp is not None:

lidar_pipeline/ign.py

@@ -76,6 +76,8 @@ def _lat_lon_to_px(lat, lon, zoom, tile_size=256):
 def download_ign_tiles(min_x, max_x, min_y, max_y, layer, zoom_level=15):
     """Download IGN WMTS tiles for the given bounds using Web Mercator (PM).
 
+    If the first tile returns 404, automatically retries at lower zoom levels.
+
     Args:
         min_x, max_x, min_y, max_y: Bounds in Lambert 93.
         layer: IGN WMTS layer name.
@@ -106,76 +108,100 @@ def download_ign_tiles(min_x, max_x, min_y, max_y, layer, zoom_level=15):
     tile_matrix_set = "PM"
     tile_size = 256
 
-    col_min, row_min = _lat_lon_to_tile(nw_lat, nw_lon, zoom_level)
-    col_max, row_max = _lat_lon_to_tile(se_lat, se_lon, zoom_level)
+    # Try downloading at the requested zoom level; fall back to lower zooms on 404
+    min_zoom = 15
+    for zoom in range(zoom_level, min_zoom - 1, -1):
+        col_min, row_min = _lat_lon_to_tile(nw_lat, nw_lon, zoom)
+        col_max, row_max = _lat_lon_to_tile(se_lat, se_lon, zoom)
 
-    nw_px_x, nw_px_y = _lat_lon_to_px(nw_lat, nw_lon, zoom_level)
-    se_px_x, se_px_y = _lat_lon_to_px(se_lat, se_lon, zoom_level)
+        nw_px_x, nw_px_y = _lat_lon_to_px(nw_lat, nw_lon, zoom)
+        se_px_x, se_px_y = _lat_lon_to_px(se_lat, se_lon, zoom)
 
-    out_width = int(se_px_x - nw_px_x)
-    out_height = int(se_px_y - nw_px_y)
+        out_width = int(se_px_x - nw_px_x)
+        out_height = int(se_px_y - nw_px_y)
 
-    if out_width <= 0 or out_height <= 0 or out_width > 10000 or out_height > 10000:
-        logger.warning(f"    Image IGN trop grande: {out_width}x{out_height}px — abandon")
-        return None
+        if out_width <= 0 or out_height <= 0 or out_width > 10000 or out_height > 10000:
+            logger.warning(f"    Image IGN trop grande: {out_width}x{out_height}px — zoom {zoom} abandon")
+            continue
 
-    total_tiles = (col_max - col_min + 1) * (row_max - row_min + 1)
-    logger.info(f"    Zoom {zoom_level}: {total_tiles} tuiles à télécharger ({out_width}x{out_height}px)")
+        total_tiles = (col_max - col_min + 1) * (row_max - row_min + 1)
+        if zoom < zoom_level:
+            logger.info(f"    ↓ Retry zoom {zoom_level}→{zoom}: {total_tiles} tuiles ({out_width}x{out_height}px)")
+        else:
+            logger.info(f"    Zoom {zoom}: {total_tiles} tuiles à télécharger ({out_width}x{out_height}px)")
 
-    composite = np.full((out_height, out_width, 3), 255, dtype=np.uint8)
+        composite = np.full((out_height, out_width, 3), 255, dtype=np.uint8)
 
-    tiles_downloaded = 0
-    fmt = "image/png" if 'PLAN' in layer else "image/jpeg"
+        tiles_downloaded = 0
+        tiles_404 = 0
+        fmt = "image/png" if 'PLAN' in layer else "image/jpeg"
 
-    for col in range(col_min, col_max + 1):
-        for row in range(row_min, row_max + 1):
-            url = (
-                f"{wmts_url}?SERVICE=WMTS&VERSION=1.0.0&REQUEST=GetTile"
-                f"&LAYER={layer}&STYLE=normal"
-                f"&TILEMATRIXSET={tile_matrix_set}"
-                f"&TILEMATRIX={zoom_level}&TILECOL={col}&TILEROW={row}"
-                f"&FORMAT={fmt}"
-            )
+        for col in range(col_min, col_max + 1):
+            for row in range(row_min, row_max + 1):
+                url = (
+                    f"{wmts_url}?SERVICE=WMTS&VERSION=1.0.0&REQUEST=GetTile"
+                    f"&LAYER={layer}&STYLE=normal"
+                    f"&TILEMATRIXSET={tile_matrix_set}"
+                    f"&TILEMATRIX={zoom}&TILECOL={col}&TILEROW={row}"
+                    f"&FORMAT={fmt}"
+                )
 
-            try:
-                req = urllib.request.Request(url, headers={'User-Agent': 'Mozilla/5.0'})
-                with urllib.request.urlopen(req, timeout=10) as response:
-                    tile_data = response.read()
-                    tile_img = PILImage.open(io.BytesIO(tile_data)).convert('RGB')
-                    tile_arr = np.array(tile_img)
+                try:
+                    req = urllib.request.Request(url, headers={'User-Agent': 'Mozilla/5.0'})
+                    with urllib.request.urlopen(req, timeout=10) as response:
+                        tile_data = response.read()
+                        tile_img = PILImage.open(io.BytesIO(tile_data)).convert('RGB')
+                        tile_arr = np.array(tile_img)
 
-                    tile_origin_x = col * tile_size
-                    tile_origin_y = row * tile_size
+                        tile_origin_x = col * tile_size
+                        tile_origin_y = row * tile_size
 
-                    px_x = int(tile_origin_x - nw_px_x)
-                    px_y = int(tile_origin_y - nw_px_y)
+                        px_x = int(tile_origin_x - nw_px_x)
+                        px_y = int(tile_origin_y - nw_px_y)
 
-                    x_off = max(0, -px_x)
-                    y_off = max(0, -px_y)
-                    dst_x_start = max(0, px_x)
-                    dst_y_start = max(0, px_y)
-                    dst_x_end = min(out_width, px_x + tile_size)
-                    dst_y_end = min(out_height, px_y + tile_size)
+                        x_off = max(0, -px_x)
+                        y_off = max(0, -px_y)
+                        dst_x_start = max(0, px_x)
+                        dst_y_start = max(0, px_y)
+                        dst_x_end = min(out_width, px_x + tile_size)
+                        dst_y_end = min(out_height, px_y + tile_size)
 
-                    src_x = x_off
-                    src_y = y_off
-                    src_w = dst_x_end - dst_x_start
-                    src_h = dst_y_end - dst_y_start
+                        src_x = x_off
+                        src_y = y_off
+                        src_w = dst_x_end - dst_x_start
+                        src_h = dst_y_end - dst_y_start
 
-                    if src_w > 0 and src_h > 0 and tile_arr.shape[0] >= src_y + src_h and tile_arr.shape[1] >= src_x + src_w:
-                        composite[dst_y_start:dst_y_end, dst_x_start:dst_x_end] = \
-                            tile_arr[src_y:src_y+src_h, src_x:src_x+src_w]
-                        tiles_downloaded += 1
+                        if src_w > 0 and src_h > 0 and tile_arr.shape[0] >= src_y + src_h and tile_arr.shape[1] >= src_x + src_w:
+                            composite[dst_y_start:dst_y_end, dst_x_start:dst_x_end] = \
+                                tile_arr[src_y:src_y+src_h, src_x:src_x+src_w]
+                            tiles_downloaded += 1
 
-            except Exception as e:
-                if tiles_downloaded == 0 and col == col_min and row == row_min:
-                    logger.error(f"    ✗ Erreur tuile IGN: {e}")
+                except urllib.error.HTTPError as e:
+                    if e.code == 404:
+                        tiles_404 += 1
+                        # If the very first tile is 404, this zoom level is unavailable
+                        if col == col_min and row == row_min:
+                            logger.info(f"    Zoom {zoom} non disponible (404) — essai zoom inférieur")
+                            break
+                    continue
+                except Exception:
+                    continue
+            else:
                 continue
+            # Only reach here if inner loop broke (first tile 404)
+            break
 
-    logger.info(f"    → {tiles_downloaded} tuiles IGN téléchargées ({layer})")
-    if tiles_downloaded == 0:
-        return None
-    return composite
+        if tiles_404 > 0 and tiles_downloaded == 0:
+            # No tiles at this zoom, try lower
+            continue
+
+        logger.info(f"    → {tiles_downloaded} tuiles IGN téléchargées ({layer})")
+        if tiles_downloaded == 0:
+            continue
+        return composite
+
+    logger.error("    ✗ Aucun zoom disponible pour cette zone")
+    return None
 
 
 def generate_ign_overlay(dem_file, basename, vis_dir, resolution, layer, title, legend_label, description, out_suffix):

lidar_pipeline/pipeline.py

@@ -12,6 +12,7 @@ import multiprocessing
 import shutil
 import time
 from concurrent.futures import ProcessPoolExecutor, as_completed
+from datetime import datetime
 from pathlib import Path
 import subprocess
 
@@ -55,16 +56,16 @@ _file_filter = FilePrefixFilter()
 
 from .dtm import classify_ground, create_dtm_fast
 from .visualizations import (
+    SharedDEM,
     generate_hillshade, generate_slope, generate_aspect, generate_curvature,
-    generate_lrm, generate_svf, generate_openness,
+    generate_lrm, generate_openness,
     generate_mslrm, generate_tpi, generate_sailore,
     generate_roughness, generate_anomalies, generate_wavelet,
     generate_flow,
 )
 from .gpu import gpu_cleanup
 from .ign import generate_ign_overlay
-from .rendering import tif_to_png, generate_pdf_report, convert_to_cog, generate_cog_metadata
-from .viewer import generate_viewer
+from .rendering import tif_to_png
 
 
 # Ordered list of visualization steps.
@@ -75,10 +76,9 @@ VIZ_STEPS = [
     ('slope',      generate_slope),
     ('aspect',     generate_aspect),
     ('curvature',  generate_curvature),
-    ('svf',        generate_svf),
     ('lrm',        generate_lrm),
-    ('pos_open',   lambda d, b, v, r: generate_openness(d, b, v, r, positive=True)),
-    ('neg_open',   lambda d, b, v, r: generate_openness(d, b, v, r, positive=False)),
+    ('pos_open',   lambda d, b, v, r, shared=None: generate_openness(d, b, v, r, positive=True, shared=shared)),
+    ('neg_open',   lambda d, b, v, r, shared=None: generate_openness(d, b, v, r, positive=False, shared=shared)),
     ('mslrm',      generate_mslrm),
     ('tpi',        generate_tpi),
     ('sailore',    generate_sailore),
@@ -106,16 +106,26 @@ VIZ_STEPS = [
 class LidarArchaeoPipeline:
     """Orchestrates the LiDAR archaeological analysis pipeline."""
 
-    def __init__(self, input_dir, output_dir, resolution=0.5, workers=1, force=False, ground_method='auto', force_classify=False, keep_tif=False, no_viewer=False):
+    def __init__(self, input_dir, output_dir, resolution=0.5, workers=1, force=False, ground_method='auto', force_classify=False, keep_tif=False, quality=98, only_viz=None, skip_viz=None, output_format='avif'):
         self.input_dir = Path(input_dir)
         self.output_dir = Path(output_dir)
-        self.resolution = resolution
+        # Accept single float or comma-separated string for multi-resolution
+        if isinstance(resolution, str):
+            self.resolutions = [float(r.strip()) for r in resolution.split(',')]
+        elif isinstance(resolution, (list, tuple)):
+            self.resolutions = [float(r) for r in resolution]
+        else:
+            self.resolutions = [float(resolution)]
+        self.resolution = self.resolutions[0]  # Primary resolution (backward compat)
         self.workers = workers
         self.force = force
         self.ground_method = ground_method
         self.force_classify = force_classify
         self.keep_tif = keep_tif
-        self.no_viewer = no_viewer
+        self.quality = quality
+        self.only_viz = only_viz
+        self.skip_viz = skip_viz
+        self.output_format = output_format
         self.temp_dir = self.output_dir / "temp"
 
         if not self.input_dir.exists():
@@ -126,21 +136,43 @@ class LidarArchaeoPipeline:
 
         self.dtm_dir = self.output_dir / "DTM"
         self.vis_dir = self.output_dir / "visualisations"
-        self.pdf_dir = self.output_dir / "rapports"
 
-        for d in [self.dtm_dir, self.vis_dir, self.pdf_dir]:
+        for d in [self.dtm_dir, self.vis_dir]:
             d.mkdir(exist_ok=True)
 
+        # Filter visualizations based on --only / --skip
+        all_viz_names = [name for name, _ in VIZ_STEPS]
+        if only_viz:
+            invalid = set(only_viz) - set(all_viz_names)
+            if invalid:
+                raise ValueError(f"Visualisations inconnues: {', '.join(invalid)}. Disponibles: {', '.join(all_viz_names)}")
+            self.viz_steps = [(n, f) for n, f in VIZ_STEPS if n in only_viz]
+        elif skip_viz:
+            invalid = set(skip_viz) - set(all_viz_names)
+            if invalid:
+                raise ValueError(f"Visualisations inconnues: {', '.join(invalid)}. Disponibles: {', '.join(all_viz_names)}")
+            self.viz_steps = [(n, f) for n, f in VIZ_STEPS if n not in skip_viz]
+        else:
+            self.viz_steps = VIZ_STEPS
+
         logger.info("Pipeline initialisé")
         logger.info(f"  Entrée      : {self.input_dir}")
         logger.info(f"  Sortie      : {self.output_dir}")
-        logger.info(f"  Résolution  : {resolution}m/px")
+        if len(self.resolutions) > 1:
+            logger.info(f"  Résolutions : {', '.join(f'{r}m/px' for r in self.resolutions)}")
+        else:
+            logger.info(f"  Résolution  : {self.resolution}m/px")
         logger.info(f"  Workers     : {workers}")
         logger.info(f"  Force       : {'OUI' if self.force else 'non (skip existing)'}")
         logger.info(f"  Classification sol : {self.ground_method}")
         logger.info(f"  Force classif.: {'OUI' if self.force_classify else 'non'}")
         logger.info(f"  Keep TIFF    : {'OUI' if self.keep_tif else 'non'}")
-        logger.info(f"  Viewer web   : {'non' if self.no_viewer else 'OUI'}")
+        logger.info(f"  Qualité {self.output_format.upper()}: {self.quality if self.quality < 100 else 'lossless'}")
+        if only_viz:
+            logger.info(f"  Visualisations: uniquement {', '.join(only_viz)}")
+        elif skip_viz:
+            logger.info(f"  Visualisations: tout sauf {', '.join(skip_viz)}")
+        logger.info(f"  Visualisations: {len(self.viz_steps)}/{len(VIZ_STEPS)}")
 
     def find_laz_files(self):
         """Find all LAZ/LAS files in input directory."""
@@ -162,27 +194,76 @@ class LidarArchaeoPipeline:
                 return False
         return True
 
-    def generate_all_visualizations(self, dtm_file, basename):
+    @staticmethod
+    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 == 'pos_open':
+            return file_vis_dir / f"{basename}_positive_openness.{ext}"
+        elif name == 'neg_open':
+            return file_vis_dir / f"{basename}_negative_openness.{ext}"
+        elif name == 'hillshade':
+            return file_vis_dir / f"{basename}_hillshade_multi.{ext}"
+        else:
+            return file_vis_dir / f"{basename}_{name}.{ext}"
+
+    def generate_all_visualizations(self, dtm_file, basename, resolution=None, vis_dir=None):
         """Generate all archaeological visualizations for one DTM file.
 
-        Returns a dict of {name: tif_path} for successful generations.
+        Optimisation: SharedDEM is only computed if at least one visualization
+        needs to be generated. When all WebP outputs exist, SharedDEM is
+        skipped entirely (saves ~2min per file on re-runs).
         """
+        if resolution is None:
+            resolution = self.resolution
         logger.info("  Génération visualisations:")
 
-        # Create per-file subdirectory
-        file_vis_dir = self.vis_dir / basename
+        # Use provided vis_dir (for multi-resolution subdirectories) or default
+        file_vis_dir = vis_dir if vis_dir else (self.vis_dir / basename)
         file_vis_dir.mkdir(exist_ok=True)
+        total = len(self.viz_steps)
 
+        # Phase 1: determine which visualizations need generation
+        needs_generation = {}  # name -> True/False
+        for name, func in self.viz_steps:
+            if self.force:
+                needs_generation[name] = True
+            else:
+                expected_webp = self._expected_output_path(name, basename, file_vis_dir, self.output_format)
+                needs_generation[name] = not expected_webp.exists()
+
+        to_generate = [n for n, needed in needs_generation.items() if needed]
+        ign_only = all(name in ('ortho', 'topo') for name in to_generate)
+        needs_shared = any(name not in ('ortho', 'topo') for name in to_generate)
+
+        if not to_generate:
+            logger.info("    Toutes les visualisations déjà existantes — ignorées")
+            # Still need to return results dict for PDF check
+            vis_results = {}
+            for name, func in self.viz_steps:
+                vis_results[name] = self._expected_output_path(name, basename, file_vis_dir, self.output_format)
+            return vis_results
+
+        # Phase 2: compute SharedDEM only if needed
+        shared = None
+        if needs_shared:
+            logger.info("    Pré-calcul données partagées (gradient, LRM)...")
+            t_shared = time.time()
+            shared = SharedDEM(dtm_file, resolution)
+            logger.info(f"    ✓ Données partagées prêtes ({time.time()-t_shared:.1f}s)")
+
+        # Phase 3: generate visualizations
         vis_results = {}
-        total = len(VIZ_STEPS)
-        elapsed_times = []
+        for idx, (name, func) in enumerate(self.viz_steps, 1):
+            if not needs_generation[name]:
+                logger.info(f"    [{idx}/{total}] {name}: déjà existant, ignoré")
+                vis_results[name] = self._expected_output_path(name, basename, file_vis_dir, self.output_format)
+                continue
 
-        for idx, (name, func) in enumerate(VIZ_STEPS, 1):
             # When --force, delete existing TIF to ensure clean regeneration
             if self.force:
                 for tif in file_vis_dir.glob(f"{basename}_{name}.tif"):
                     tif.unlink(missing_ok=True)
-                # Special cases for differently-named TIFs
                 if name == 'pos_open':
                     for tif in file_vis_dir.glob(f"{basename}_positive_openness.tif"):
                         tif.unlink(missing_ok=True)
@@ -193,40 +274,18 @@ class LidarArchaeoPipeline:
                     for tif in file_vis_dir.glob(f"{basename}_hillshade_multi.tif"):
                         tif.unlink(missing_ok=True)
 
-            # Check if output WebP already exists (skip unless --force)
-            if not self.force:
-                # Determine expected WebP filename from the viz name
-                # Special cases for openness and IGN overlays
-                if name == 'pos_open':
-                    expected_webp = file_vis_dir / f"{basename}_positive_openness.webp"
-                elif name == 'neg_open':
-                    expected_webp = file_vis_dir / f"{basename}_negative_openness.webp"
-                elif name == 'hillshade':
-                    expected_webp = file_vis_dir / f"{basename}_hillshade_multi.webp"
-                elif name in ('ortho', 'topo'):
-                    expected_webp = file_vis_dir / f"{basename}_{name}.webp"
-                else:
-                    expected_webp = file_vis_dir / f"{basename}_{name}.webp"
-
-                if expected_webp.exists():
-                    logger.info(f"    [{idx}/{total}] {name}: déjà existant, ignoré")
-                    vis_results[name] = expected_webp  # Track as existing file
-                    continue
-
             logger.info(f"    [{idx}/{total}] {name}...")
             t0 = time.time()
             try:
-                result = func(dtm_file, basename, file_vis_dir, self.resolution)
+                # IGN overlays don't use SharedDEM (they download external data)
+                if name in ('ortho', 'topo'):
+                    result = func(dtm_file, basename, file_vis_dir, resolution)
+                else:
+                    result = func(dtm_file, basename, file_vis_dir, resolution, shared=shared)
                 vis_results[name] = result
                 elapsed = time.time() - t0
-                elapsed_times.append(elapsed)
                 if result:
-                    eta = ""
-                    if len(elapsed_times) > 1:
-                        avg_time = sum(elapsed_times) / len(elapsed_times)
-                        remaining = (total - idx) * avg_time
-                        eta = f" — ETA: {remaining:.0f}s"
-                    logger.info(f"    [{idx}/{total}] ✓ {name} ({elapsed:.1f}s){eta}")
+                    logger.info(f"    [{idx}/{total}] ✓ {name} ({elapsed:.1f}s)")
                 else:
                     logger.warning(f"    [{idx}/{total}] ✗ {name} — no output ({elapsed:.1f}s)")
             except Exception as e:
@@ -236,25 +295,43 @@ class LidarArchaeoPipeline:
             # Free GPU memory between visualizations to prevent OOM
             gpu_cleanup()
 
-        # Convert to WebP (only newly generated TIFs, not skipped ones)
-        # Also generate COGs for web viewer if enabled
-        logger.info("  Conversion images WebP:")
-        cog_dir = file_vis_dir / "cog"
-        if not self.no_viewer:
-            cog_dir.mkdir(exist_ok=True)
+        # Convert to output format (only newly generated TIFs, not skipped ones)
+        fmt_label = self.output_format.upper()
+        logger.info(f"  Conversion images {fmt_label}:")
+        source_info = {
+            'method': self.ground_method,
+            'date': datetime.now().strftime('%Y-%m-%d'),
+            'basename': basename,
+        }
         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():
-                # Generate COG before WebP conversion (which may delete the TIF)
-                if not self.no_viewer:
-                    convert_to_cog(tif_file, cog_dir)
-                webp_file = tif_to_png(tif_file, file_vis_dir, self.resolution, keep_tif=self.keep_tif or not self.no_viewer)
-                if webp_file:
-                    logger.info(f"    ✓ {webp_file.name}")
+                img_file = tif_to_png(tif_file, file_vis_dir, resolution, keep_tif=self.keep_tif, source_info=source_info, quality=self.quality, output_format=self.output_format)
+                if img_file:
+                    logger.info(f"    ✓ {img_file.name}")
+
+        # Clean up remaining TIF files unless --keep-tif
+        if not self.keep_tif:
+            for tif in file_vis_dir.glob("*.tif"):
+                tif.unlink(missing_ok=True)
 
         return vis_results
 
+    @staticmethod
+    def _res_suffix(resolution):
+        """Return suffix for additional resolutions (empty string for primary)."""
+        if resolution == 0.5:
+            return ""  # Default resolution — no suffix
+        res_str = f"{resolution}".replace('.', 'p')
+        return f"_r{res_str}"
+
     def process_file(self, laz_file):
-        """Process a single LAZ file through the full pipeline."""
+        """Process a single LAZ file through the full pipeline.
+
+        If self.resolutions has multiple entries, processes each resolution:
+        - Primary resolution uses current naming (no suffix)
+        - Additional resolutions use _r0p2 suffix in directories/filenames
+        - Ground classification is done once and shared across resolutions
+        """
         basename = _file_basename(laz_file)
         _file_filter.basename = basename
         t_start = time.time()
@@ -263,64 +340,89 @@ class LidarArchaeoPipeline:
         logger.info(f"FICHIER : {basename}")
         logger.info("=" * 60)
 
-        # Step 1: Ground classification
-        logger.info("[1/6] Classification du sol...")
-        t1 = time.time()
-        las_file = classify_ground(laz_file, self.temp_dir, method=self.ground_method, force=self.force_classify)
-        t_classif = time.time() - t1
-        if not las_file:
-            logger.error(f"  ✗ Échec classification ({t_classif:.1f}s)")
+        # Validate file integrity before any processing
+        from .dtm import validate_laz
+        if not validate_laz(laz_file):
             return False
-        logger.info(f"  ✓ Classification terminée ({t_classif:.1f}s)")
 
-        # Step 2: Generate DTM
-        logger.info("[2/6] Génération DTM...")
-        t2 = time.time()
-        dtm_file = create_dtm_fast(las_file, basename, self.dtm_dir, self.resolution)
-        t_dtm = time.time() - t2
-        if not dtm_file:
-            logger.error(f"  ✗ Échec DTM ({t_dtm:.1f}s)")
-            return False
-        logger.info(f"  ✓ DTM terminé ({t_dtm:.1f}s)")
+        # Step 1: Ground classification (shared across all resolutions)
+        las_file = None
+        t_classif = 0
+        for i, res in enumerate(self.resolutions):
+            res_suffix = self._res_suffix(res)
+            dtm_path = self.dtm_dir / f"{basename}_dtm{res_suffix}.tif"
+            if dtm_path.exists():
+                import rasterio
+                try:
+                    with rasterio.open(dtm_path) as src:
+                        existing_res = abs(src.transform.a)
+                    if abs(existing_res - res) > 0.01:
+                        logger.info(f"  DTM{res_suffix} existant à {existing_res}m/px — résolution demandée {res}m/px → régénération")
+                        dtm_path.unlink()
+                    else:
+                        if i == 0:
+                            logger.info(f"[1/5] Classification du sol — sautée (DTM existant)")
+                            logger.info(f"[2/5] Génération DTM {res}m/px — sautée (DTM existant)")
+                        else:
+                            logger.info(f"  DTM {res}m/px déjà existant — ignoré")
+                        continue
+                except Exception:
+                    logger.warning(f"Impossible de lire le DTM existant — régénération")
+                    dtm_path.unlink()
 
-        # Step 3: Visualizations
-        logger.info("[3/6] Visualisations archéologiques...")
-        self.generate_all_visualizations(dtm_file, basename)
+            # Need to classify/generate DTM for this resolution
+            if las_file is None:
+                # First time: do ground classification
+                logger.info("[1/5] Classification du sol...")
+                t1 = time.time()
+                las_file = classify_ground(laz_file, self.temp_dir, method=self.ground_method, force=self.force_classify)
+                t_classif = time.time() - t1
+                if not las_file:
+                    logger.error(f"  ✗ Échec classification ({t_classif:.1f}s)")
+                    return False
+                logger.info(f"  ✓ Classification terminée ({t_classif:.1f}s)")
 
-        # Step 4: PDF report
-        file_vis_dir = self.vis_dir / basename
-        logger.info("[4/6] Rapport PDF A3...")
-        t4 = time.time()
-        generate_pdf_report(basename, file_vis_dir, self.pdf_dir, self.resolution)
-        t_pdf = time.time() - t4
-        logger.info(f"  ✓ Rapport PDF terminé ({t_pdf:.1f}s)")
+            # Generate DTM at this resolution
+            logger.info(f"{'[2/5]' if i == 0 else '  '} Génération DTM {res}m/px...")
+            t2 = time.time()
+            dtm_file = create_dtm_fast(las_file, basename, self.dtm_dir, res, force=self.force, output_suffix=res_suffix)
+            t_dtm = time.time() - t2
+            if not dtm_file:
+                logger.error(f"  ✗ Échec DTM {res}m/px ({t_dtm:.1f}s)")
+                if i == 0:
+                    return False  # Primary resolution failure is fatal
+                continue  # Additional resolution failure is non-fatal
+            logger.info(f"  ✓ DTM {res}m/px terminé ({t_dtm:.1f}s)")
 
-        # Step 5: COGs for web viewer
-        logger.info("[5/6] Génération métadonnées viewer web...")
-        t5 = time.time()
-        if not self.no_viewer:
-            # Convert DTM to COG as well
-            dtm_cog_dir = self.dtm_dir / "cog"
-            dtm_cog_dir.mkdir(exist_ok=True)
-            for dtm_file in sorted(self.dtm_dir.glob(f"{basename}_dtm.tif")):
-                convert_to_cog(dtm_file, dtm_cog_dir)
-            generate_cog_metadata(self.vis_dir, basename)
-        t_cog = time.time() - t5
-        logger.info(f"  ✓ Métadonnées viewer web terminées ({t_cog:.1f}s)")
+        # Process each resolution: visualizations + PDF
+        all_vis_results = {}
+        for res in self.resolutions:
+            res_suffix = self._res_suffix(res)
+            dtm_path = self.dtm_dir / f"{basename}_dtm{res_suffix}.tif"
 
-        # Step 6: Web viewer
-        if not self.no_viewer:
-            logger.info("[6/6] Génération viewer web...")
-            t6 = time.time()
-            generate_viewer(basename, file_vis_dir, self.vis_dir)
-            t_viewer = time.time() - t6
-            logger.info(f"  ✓ Viewer web terminé ({t_viewer:.1f}s)")
-        else:
-            logger.info("[6/6] Viewer web: ignoré (--no-viewer)")
+            if not dtm_path.exists():
+                logger.warning(f"  DTM {res}m/px manquant — visualisations ignorées")
+                continue
+
+            import rasterio
+            with rasterio.open(dtm_path) as src:
+                actual_res = abs(src.transform.a)
+
+            if len(self.resolutions) > 1:
+                logger.info(f"  --- Résolution {res}m/px ---")
+
+            # For additional resolutions, use suffixed subdirectory
+            if res_suffix:
+                vis_dir = self.vis_dir / f"{basename}{res_suffix}"
+            else:
+                vis_dir = self.vis_dir / basename
+
+            vis_dir.mkdir(exist_ok=True)
+
+            self.generate_all_visualizations(dtm_path, basename, actual_res, vis_dir=vis_dir)
 
         t_total = time.time() - t_start
         logger.info(f"✓ {basename} terminé en {t_total:.1f}s")
-        logger.debug(f"  Détails: classification={t_classif:.1f}s, DTM={t_dtm:.1f}s, PDF={t_pdf:.1f}s")
         _file_filter.basename = None
         return True
 
@@ -348,29 +450,42 @@ class LidarArchaeoPipeline:
             logger.info(f"Traitement parallèle avec {self.workers} workers...")
             logger.info(f"Fichiers: {len(files)}")
             with ProcessPoolExecutor(max_workers=self.workers) as executor:
+                # Pass resolutions as comma-separated string for multiprocessing serialization
+                resolutions_str = ','.join(str(r) for r in self.resolutions)
                 future_to_file = {
-                    executor.submit(_process_file_standalone, str(laz_file), str(self.input_dir), str(self.output_dir), self.resolution, self.force, self.ground_method, self.force_classify, self.keep_tif, self.no_viewer): laz_file
+                    executor.submit(_process_file_standalone, str(laz_file), str(self.input_dir), str(self.output_dir), resolutions_str, self.force, self.ground_method, self.force_classify, self.keep_tif, self.quality, self.only_viz, self.skip_viz, self.output_format): laz_file
                     for laz_file in files
                 }
                 done = 0
-                for future in as_completed(future_to_file):
-                    laz_file = future_to_file[future]
-                    done += 1
-                    try:
-                        success = future.result()
-                        results[laz_file.name] = success
-                        status = "✓" if success else "✗"
-                        logger.info(f"  [{done}/{len(files)}] {status} {laz_file.name}")
-                    except Exception as e:
-                        logger.error(f"  [{done}/{len(files)}] ✗ {laz_file.name}: {e}")
-                        logger.debug(f"  Traceback:", exc_info=True)
-                        results[laz_file.name] = False
+                try:
+                    for future in as_completed(future_to_file):
+                        laz_file = future_to_file[future]
+                        done += 1
+                        try:
+                            success = future.result()
+                            results[laz_file.name] = success
+                            status = "✓" if success else "✗"
+                            logger.info(f"  [{done}/{len(files)}] {status} {laz_file.name}")
+                        except Exception as e:
+                            logger.error(f"  [{done}/{len(files)}] ✗ {laz_file.name}: {e}")
+                            logger.debug(f"  Traceback:", exc_info=True)
+                            results[laz_file.name] = False
+                except KeyboardInterrupt:
+                    logger.info("Interruption — annulation des travaux en cours...")
+                    for f in future_to_file:
+                        f.cancel()
+                    executor.shutdown(wait=False, cancel_futures=True)
+                    logger.info("Travaux annulés.")
+                    return
         else:
             total = len(files)
             for idx, laz_file in enumerate(files, 1):
                 logger.info(f"--- Fichier {idx}/{total} ---")
                 try:
                     results[laz_file.name] = self.process_file(laz_file)
+                except KeyboardInterrupt:
+                    logger.info("Interruption — arrêt immédiat.")
+                    return
                 except Exception as e:
                     logger.error(f"✗ Erreur traitement {laz_file.name}: {e}")
                     logger.debug("Traceback:", exc_info=True)
@@ -384,18 +499,18 @@ class LidarArchaeoPipeline:
         logger.info("=" * 60)
         logger.info("RÉSUMÉ")
         logger.info("=" * 60)
+        for name, ok in results.items():
+            status = "✓" if ok else "✗"
+            logger.info(f"  {status} {name}")
+        logger.info("-" * 60)
         logger.info(f"  Succès : {success_count}/{len(results)}")
         if fail_count:
             logger.info(f"  Échecs : {fail_count}/{len(results)}")
-            for name, ok in results.items():
-                if not ok:
-                    logger.info(f"    ✗ {name}")
         logger.info(f"  Durée totale : {t_pipeline_total:.1f}s ({t_pipeline_total/60:.1f}min)")
 
         logger.info(f"\nRésultats dans: {self.output_dir}")
         logger.info(f"  • DTM           : {self.dtm_dir}")
         logger.info(f"  • Visualisations: {self.vis_dir}")
-        logger.info(f"  • Rapports PDF  : {self.pdf_dir}")
 
         # Clean up temporary files
         logger.info("Nettoyage des fichiers temporaires...")
@@ -411,7 +526,7 @@ class LidarArchaeoPipeline:
             logger.warning(f"  Note: Impossible de supprimer les fichiers temporaires: {e}")
 
 
-def _process_file_standalone(laz_file_str, input_dir, output_dir, resolution, force=False, ground_method='auto', force_classify=False, keep_tif=False, no_viewer=False):
+def _process_file_standalone(laz_file_str, input_dir, output_dir, resolution, force=False, ground_method='auto', force_classify=False, keep_tif=False, quality=98, only_viz=None, skip_viz=None, output_format='avif'):
     """Standalone function for multiprocessing — creates its own pipeline instance.
 
     Each worker gets its own temp directory to avoid file conflicts.
@@ -419,6 +534,11 @@ def _process_file_standalone(laz_file_str, input_dir, output_dir, resolution, fo
     # Configure logging in worker process (spawn doesn't inherit parent config)
     import logging
     import sys
+    # Ensure UTF-8 output — spawn workers may default to ASCII
+    if hasattr(sys.stdout, 'reconfigure'):
+        sys.stdout.reconfigure(encoding='utf-8', errors='replace')
+    if hasattr(sys.stderr, 'reconfigure'):
+        sys.stderr.reconfigure(encoding='utf-8', errors='replace')
     worker_logger = logging.getLogger("lidar")
     if not worker_logger.handlers:
         handler = logging.StreamHandler(sys.stdout)
@@ -427,7 +547,7 @@ def _process_file_standalone(laz_file_str, input_dir, output_dir, resolution, fo
         worker_logger.addHandler(handler)
     worker_logger.addFilter(_file_filter)
 
-    pipeline = LidarArchaeoPipeline(input_dir, output_dir, resolution=resolution, workers=1, force=force, ground_method=ground_method, force_classify=force_classify, keep_tif=keep_tif, no_viewer=no_viewer)
+    pipeline = LidarArchaeoPipeline(input_dir, output_dir, resolution=resolution, workers=1, force=force, ground_method=ground_method, force_classify=force_classify, keep_tif=keep_tif, quality=quality, only_viz=only_viz, skip_viz=skip_viz, output_format=output_format)
     basename = _file_basename(laz_file_str)
     pipeline.temp_dir = pipeline.output_dir / "temp" / basename
     pipeline.temp_dir.mkdir(exist_ok=True)

lidar_pipeline/rendering.py

@@ -1,13 +1,14 @@
-"""Rendering module: colormap registry, GeoTIFF-to-WebP conversion, and PDF report generation.
+"""Rendering module: colormap registry, GeoTIFF-to-image conversion, and PDF report generation.
 
 Contains:
 - COLORMAPS: registry mapping filename keywords to (cmap, title, legend, description)
-- tif_to_png(): convert a GeoTIFF to a WebP visualization with legend, scale bar, north arrow
+- tif_to_png(): convert a GeoTIFF to a WebP/AVIF visualization with legend, scale bar, north arrow
 - generate_pdf_report(): generate an A3 PDF report with all visualizations
 """
 
 import logging
 import time
+from datetime import datetime
 from pathlib import Path
 
 import numpy as np
@@ -73,18 +74,10 @@ COLORMAPS = {
         'description': 'Détecte les ruptures de pente — utile pour bords de terrasses et levées',
         'vmin_mode': 'symmetric', 'sym_pct': (5, 95),
     },
-    'svf': {
-        'cmap': 'viridis',
-        'title': 'Sky-View Factor (Ray-tracing 16 directions)',
-        'legend': 'Fraction de ciel visible depuis chaque point\nSombre = Encaissé (fossés, vallées, rues)\nClair = Dégagé (sommets, plateformes, plateaux)',
-        'description': 'Ray-tracing sur 16 azimuts — élimine l\'éclairage, détecte structures linéaires et enclos',
-        'vmin_mode': 'percentile', 'vmin_pct': 5,
-        'vmax_mode': 'percentile', 'vmax_pct': 95,
-    },
     'mslrm': {
         'cmap': 'RdBu_r',
-        'title': 'MSRM - Multi-Scale Relief Model (5 échelles)',
-        'legend': 'Relief combiné à 5 échelles (5m à 100m)\nRouge = Surélévation (mur, tumulus, levée)\nBleu = Dépression (fossé, douve, fossé)\n\nDifférence avec LRM:\nLRM = 1 échelle (15m)\nMSRM = 5 échelles combinées\nMSRM détecte du micro au macro',
+        'title': 'MSRM - Multi-Scale Relief Model (échelles adaptatives)',
+        'legend': 'Relief combiné multi-échelles (adapté à la résolution)\nRouge = Surélévation (mur, tumulus, levée)\nBleu = Dépression (fossé, douve)\n\nDifférence avec LRM:\nLRM = 1 échelle (15m)\nMSRM = échelles combinées pondérées\nMSRM détecte du micro au macro',
         'description': 'Combine LRM à 5 échelles — détecte structures de 5m à 100m simultanément',
         'vmin_mode': 'symmetric', 'sym_pct': (2, 98),
     },
@@ -113,8 +106,8 @@ COLORMAPS = {
     },
     'tpi': {
         'cmap': 'BrBG',
-        'title': 'TPI - Topographic Position Index (2 échelles)',
-        'legend': 'Position dans le paysage\nBrun/Sombre = Plus bas que le voisinage (fossé, vallée)\nVert/Clair = Plus haut que le voisinage (crête, plateau)\nCombine échelle fine 5m + large 100m',
+        'title': 'TPI - Topographic Position Index (4 échelles)',
+        'legend': 'Position dans le paysage\nBrun/Sombre = Plus bas que le voisinage (fossé, vallée)\nVert/Clair = Plus haut que le voisinage (crête, plateau)\nCombine 4 échelles : 3m, 15m, 50m, 200m',
         'description': 'Identifie la position topographique — utile pour repérer crêtes vs vallées à grande échelle',
         'vmin_mode': 'symmetric', 'sym_pct': (2, 98),
     },
@@ -127,23 +120,23 @@ COLORMAPS = {
     },
     'roughness': {
         'cmap': 'magma',
-        'title': 'Rugosité de Surface (Écart-type local 5m)',
-        'legend': 'Irrégularité du terrain dans un voisinage de 5m\nSombre = Surface lisse (route, mur, sol plat)\nClair = Surface rugueuse (végétation, ruines, pierres)\nMax: {vmax:.2f}m',
+        'title': 'Rugosité Multi-Échelle (3m + 15m)',
+        'legend': 'Irrégularité du terrain combinée fine + large\nSombre = Surface lisse (route, mur, sol plat)\nClair = Surface rugueuse (végétation, ruines, pierres)\nCombine rugosité fine 3m (70%) + large 15m (30%)',
         'description': 'Mesure la variabilité locale — surfaces anthropiques lisses vs naturelles rugueuses',
         'vmin_mode': 'fixed', 'vmin_val': 0,
         'vmax_mode': 'percentile', 'vmax_pct': 97,
     },
     'anomalies': {
         'cmap': 'coolwarm',
-        'title': 'Anomalies Statistiques (Z-score x Moran\'s I)',
-        'legend': 'Anomalies topographiques significatives\nRouge vif = Surélévation anormale (mur, tumulus)\nBleu vif = Dépression anormale (fossé, doline)\nBlanc/gris = Normal\n\nCombine z-score (intensité) et\nMoran\'s I (regroupement spatial)',
+        'title': 'Anomalies Statistiques (MSRM multi-échelle + Moran\'s I)',
+        'legend': 'Anomalies topographiques significatives\nRouge vif = Surélévation anormale (mur, tumulus)\nBleu vif = Dépression anormale (fossé, doline)\nBlanc/gris = Normal\n\nCombine MSRM normalisé (intensité) et\nMoran\'s I (regroupement spatial)',
         'description': 'Détecte uniquement les anomalies statistiquement significatives — filtre le bruit de fond',
         'vmin_mode': 'symmetric', 'sym_pct': (5, 95),
     },
     'wavelet': {
         'cmap': 'cividis',
         'title': 'Ondelette Mexican Hat (CWT multi-échelle)',
-        'legend': 'Réponse de la transformée en ondelette à 5 échelles\nÉchelles: 2m, 5m, 10m, 20m, 50m\n\nClair = Structure détectée à cette échelle\nSombre = Pas de structure\n\nOptimisé pour formes circulaires:\ntumulus, enclos, fossés annulaires',
+        'legend': 'Réponse de la transformée en ondelette\nÉchelles adaptées à la résolution\n\nClair = Structure détectée à cette échelle\nSombre = Pas de structure\n\nOptimisé pour formes circulaires:\ntumulus, enclos, fossés annulaires',
         'description': 'Transformée en ondelette 2D — excellente pour détecter structures circulaires',
         'vmin_mode': 'symmetric', 'sym_pct': (2, 98),
     },
@@ -238,26 +231,54 @@ def _apply_colormap(data, tif_file):
     return data, 'terrain', title, 'Altitude normalisée', '', False
 
 
-def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False):
-    """Convert GeoTIFF to visualization WebP with GPS coordinates, legend, and scale bar.
+def _nice_scale(extent_m):
+    """Choose a nice round scale distance that fits well in the image.
+
+    Returns (scale_m, label) where label is like '100 m' or '500 m' or '1 km'.
+    """
+    nice_scales = [50, 100, 200, 500, 1000, 2000, 5000, 10000]
+    # Pick the largest scale that is <= 15% of the extent
+    for s in nice_scales:
+        if s <= extent_m * 0.15:
+            continue
+        # s is the first one > 15% — take the one below
+        break
+    else:
+        s = nice_scales[-1]
+    # Actually pick the largest scale <= 20% of extent
+    chosen = nice_scales[0]
+    for s in nice_scales:
+        if s <= extent_m * 0.20:
+            chosen = s
+    if chosen >= 1000:
+        return chosen, f"{chosen // 1000} km"
+    return chosen, f"{chosen} m"
+
+
+def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False, source_info=None, quality=98, output_format='avif'):
+    """Convert GeoTIFF to visualization image (WebP or AVIF) with GPS coordinates, legend, and scale bar.
 
     Args:
         tif_file: Path to input GeoTIFF.
-        vis_dir: Output directory for the WebP file.
+        vis_dir: Output directory for the image file.
         resolution: Grid resolution in m/px.
         keep_tif: If True, keep the source TIFF after conversion.
+        source_info: Dict with method/date/basename for metadata.
+        quality: Image quality (1-100). Use 100 for lossless. Default 95.
+        output_format: Output format ('webp' or 'avif'). Default 'webp'.
 
     Returns:
-        Path to output WebP file, or None on failure.
+        Path to output image file, or None on failure.
     """
     if not tif_file or not tif_file.exists():
         return None
 
-    webp_file = vis_dir / f"{tif_file.stem}.webp"
+    ext = 'avif' if output_format == 'avif' else 'webp'
+    output_file = vis_dir / f"{tif_file.stem}.{ext}"
 
     try:
         with rasterio.open(tif_file) as src:
-            is_rgb = src.count >= 3 and ('ortho' in str(tif_file) or 'topo' in str(tif_file))
+            is_rgb = src.count >= 3 and any(k in str(tif_file) for k in ('ortho', 'topo'))
 
             if is_rgb:
                 data = src.read([1, 2, 3])
@@ -362,9 +383,9 @@ def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False):
             # Fixed data area position — identical for ALL visualization types
             # This ensures overlay/superposition works across all WebP images
             data_left = 0.08
-            data_bottom = 0.12
+            data_bottom = 0.19
             data_width_frac = 0.74
-            data_height_frac = 0.78
+            data_height_frac = 0.71
 
             ax = fig.add_axes([data_left, data_bottom, data_width_frac, data_height_frac])
             if is_rgba or is_rgb:
@@ -440,20 +461,32 @@ def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False):
                 spine.set_color('black')
                 spine.set_linewidth(0.8)
 
-            # North arrow
-            north_ax = inset_axes(ax, width="4%", height="7%", loc='upper right',
-                                 bbox_to_anchor=(-0.03, 0.12, 1, 1), bbox_transform=ax.transAxes)
-            north_ax.set_xlim(0, 1)
-            north_ax.set_ylim(0, 1)
+            # North arrow — compass rose style
+            north_ax = inset_axes(ax, width="5%", height="9%", loc='upper right',
+                                 bbox_to_anchor=(-0.03, 0.08, 1, 1), bbox_transform=ax.transAxes)
+            north_ax.set_xlim(-1.2, 1.2)
+            north_ax.set_ylim(-0.5, 1.5)
             north_ax.axis('off')
-            north_ax.plot([0.5, 0.5], [0.1, 0.65], color='black', linewidth=2.5, zorder=10)
-            north_ax.add_patch(MplPolygon([[0.5, 0.2], [0.35, 0.4], [0.5, 0.65], [0.65, 0.4]],
-                                    closed=True, facecolor='black', edgecolor='black', zorder=9))
-            north_ax.text(0.5, 0.95, 'N', ha='center', va='top',
-                         fontsize=13, fontweight='bold', color='black', zorder=11)
+            north_ax.set_aspect('equal')
+            # N arrow
+            north_ax.annotate('N', xy=(0, 1.3), fontsize=11, fontweight='bold',
+                             ha='center', va='bottom', color='#b22222')
+            north_ax.plot([0, 0], [0.0, 1.0], color='#b22222', linewidth=2.0, zorder=10)
+            north_ax.add_patch(MplPolygon([[0, 0.3], [-0.2, 0.7], [0, 1.0], [0.2, 0.7]],
+                                    closed=True, facecolor='#b22222', edgecolor='#b22222', zorder=9))
+            # Cardinal ticks
+            for angle, label in [(90, ''), (0, 'E'), (180, 'O'), (270, 'S')]:
+                rad = np.radians(angle)
+                r_text = 1.25
+                north_ax.plot([0.85*np.cos(rad), 1.05*np.cos(rad)],
+                             [0.85*np.sin(rad), 1.05*np.sin(rad)],
+                             color='#555555', linewidth=0.8, zorder=5)
+                if label:
+                    north_ax.text(r_text*np.cos(rad), r_text*np.sin(rad), label,
+                                 fontsize=7, ha='center', va='center', color='#555555')
 
-            # Bottom info bar — fixed position
-            info_ax = fig.add_axes([data_left, 0.02, data_width_frac + cbar_width + 0.02, 0.07])
+            # Bottom info bar — enriched with source, method, date
+            info_ax = fig.add_axes([data_left, 0.015, data_width_frac + cbar_width + 0.02, 0.09])
             info_ax.axis('off')
 
             extent_km_x = (max_x - min_x) / 1000
@@ -465,46 +498,77 @@ def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False):
                 alt_min = float(np.nanmin(valid_data)) if len(valid_data) > 0 else 0
                 alt_max = float(np.nanmax(valid_data)) if len(valid_data) > 0 else 0
 
+            # Build info lines
+            line1_parts = []
             if gps_coords:
                 nw_lat, nw_lon = gps_coords['NW']
                 se_lat, se_lon = gps_coords['SE']
-                info_text = (
-                    f"GPS: {nw_lat:.5f}N {nw_lon:.5f}E - {se_lat:.5f}N {se_lon:.5f}E  |  "
-                    f"EPSG:2154  |  Res: {resolution}m/px  |  "
-                    f"Emprise: {extent_km_x:.1f}x{extent_km_y:.1f}km"
-                )
+                line1_parts.append(f"GPS: {nw_lat:.5f}°N {nw_lon:.5f}°E — {se_lat:.5f}°N {se_lon:.5f}°E")
             else:
-                info_text = (
-                    f"EPSG:2154  |  X: {min_x:.0f}-{max_x:.0f}  Y: {min_y:.0f}-{max_y:.0f}  |  "
-                    f"Res: {resolution}m/px  |  Emprise: {extent_km_x:.1f}x{extent_km_y:.1f}km"
-                )
+                line1_parts.append(f"X: {min_x:.0f}–{max_x:.0f}  Y: {min_y:.0f}–{max_y:.0f}")
+            line1_parts.append(f"EPSG:2154")
+            line1_parts.append(f"Res: {resolution}m/px")
+            line1_parts.append(f"Emprise: {extent_km_x:.1f}×{extent_km_y:.1f}km")
+            if not is_rgb:
+                line1_parts.append(f"Alt: {alt_min:.1f}–{alt_max:.1f}m")
 
-            info_ax.text(0.01, 0.5, info_text,
-                        transform=info_ax.transAxes, fontsize=8.5,
+            line2_parts = []
+            line2_parts.append("Source: LiDAR HD IGN")
+            if source_info:
+                if source_info.get('method'):
+                    line2_parts.append(f"Classif.: {source_info['method'].upper()}")
+                if source_info.get('date'):
+                    line2_parts.append(f"Date: {source_info['date']}")
+            else:
+                line2_parts.append(datetime.now().strftime("Date: %Y-%m-%d"))
+
+            info_text_line1 = "  |  ".join(line1_parts)
+            info_text_line2 = "  |  ".join(line2_parts)
+
+            info_ax.text(0.01, 0.7, info_text_line1,
+                        transform=info_ax.transAxes, fontsize=8,
                         verticalalignment='center', family='monospace',
-                        bbox=dict(boxstyle='round,pad=0.3', facecolor='#f0f0f0',
+                        bbox=dict(boxstyle='round,pad=0.2', facecolor='#f0f0f0',
                                   edgecolor='#aaaaaa', alpha=0.95))
+            info_ax.text(0.01, 0.2, info_text_line2,
+                        transform=info_ax.transAxes, fontsize=7.5,
+                        verticalalignment='center', family='monospace',
+                        color='#444444',
+                        bbox=dict(boxstyle='round,pad=0.2', facecolor='#f8f8f8',
+                                  edgecolor='#cccccc', alpha=0.9))
 
-            # Scale bar
-            scale_m = 100
+            # Scale bar — adaptive with alternating black/white segments
+            extent_m_x = max_x - min_x
+            scale_m, scale_label = _nice_scale(extent_m_x)
             pixels_per_meter = 1.0 / pixel_size_x
             scale_px = int(scale_m * pixels_per_meter)
-            scale_bar_frac = scale_px / width
-            bar_left = 0.80
-            bar_bottom = 0.15
-            bar_width_frac = min(scale_bar_frac, 0.15)
-            bar_height = 0.35
+            n_segments = 5
+            segment_px = scale_px / n_segments
+            bar_bottom_y = 0.55
+            bar_top_y = 0.85
+            bar_height = bar_top_y - bar_bottom_y
+            scale_start_x = 0.80
 
-            info_ax.add_patch(RectPatch((bar_left, bar_bottom), bar_width_frac, bar_height,
-                                       facecolor='black', edgecolor='black', linewidth=0.5,
-                                       transform=info_ax.transAxes, clip_on=False))
-            info_ax.text(bar_left + bar_width_frac / 2, bar_bottom + bar_height + 0.12,
-                        f"{scale_m} m", ha='center', va='bottom', fontsize=9, fontweight='bold',
-                        transform=info_ax.transAxes)
+            for seg_i in range(n_segments):
+                color = 'black' if seg_i % 2 == 0 else 'white'
+                seg_left = scale_start_x + seg_i * segment_px / width
+                seg_width_frac = segment_px / width
+                info_ax.add_patch(RectPatch((seg_left, bar_bottom_y), seg_width_frac, bar_height,
+                                           facecolor=color, edgecolor='black', linewidth=0.5,
+                                           transform=info_ax.transAxes, clip_on=False))
+            info_ax.text(scale_start_x + scale_px / (2 * width), bar_top_y + 0.12,
+                         f"{scale_label}", ha='center', va='bottom', fontsize=8, fontweight='bold',
+                         transform=info_ax.transAxes)
+            # Scale end ticks
+            info_ax.plot([scale_start_x, scale_start_x], [bar_bottom_y - 0.05, bar_top_y + 0.05],
+                        color='black', linewidth=1, transform=info_ax.transAxes, clip_on=False)
+            info_ax.plot([scale_start_x + scale_px / width, scale_start_x + scale_px / width],
+                        [bar_bottom_y - 0.05, bar_top_y + 0.05],
+                        color='black', linewidth=1, transform=info_ax.transAxes, clip_on=False)
 
             fig.patch.set_facecolor('white')
 
-            # Save as PNG then convert to WebP — fixed layout, no bbox_inches='tight'
+            # Save as PNG then convert to final format — fixed layout, no bbox_inches='tight'
             save_dpi = 200 if width > 3000 else 150
             png_temp = vis_dir / f"{tif_file.stem}_temp.png"
             try:
@@ -513,159 +577,24 @@ def tif_to_png(tif_file, vis_dir, resolution, keep_tif=False):
                 plt.close()
 
             img = PILImage.open(str(png_temp))
-            img.save(str(webp_file), format='WEBP', lossless=True)
+            pil_format = 'AVIF' if output_format == 'avif' else 'WEBP'
+            if quality >= 100:
+                img.save(str(output_file), format=pil_format, lossless=True)
+            else:
+                img.save(str(output_file), format=pil_format, quality=quality)
             png_temp.unlink(missing_ok=True)
 
             # Delete source TIFF (unless --keep-tif)
             if not keep_tif:
                 tif_file.unlink(missing_ok=True)
 
-            return webp_file
+            return output_file
 
     except Exception as e:
-        logger.error(f"    Erreur conversion WebP: {e}", exc_info=True)
+        logger.error(f"    Erreur conversion {ext.upper()}: {e}", exc_info=True)
         return None
 
 
-def convert_to_cog(tif_file, cog_dir):
-    """Convert a GeoTIFF to Cloud Optimized GeoTIFF for web map serving.
-
-    Args:
-        tif_file: Path to input GeoTIFF.
-        cog_dir: Directory for output COG file.
-
-    Returns:
-        Path to output COG file, or None on failure.
-    """
-    cog_file = cog_dir / tif_file.name
-
-    if cog_file.exists():
-        logger.debug(f"    COG déjà existant: {cog_file.name}")
-        return cog_file
-
-    try:
-        from rio_cogeo.cogeo import cog_translate
-        from rio_cogeo.profiles import cog_profiles
-
-        with rasterio.open(tif_file) as src:
-            is_rgb = src.count >= 3
-
-        # Use deflate compression profile
-        profile = dict(cog_profiles["deflate"])  # Make a mutable copy
-
-        if not is_rgb:
-            # Single-band terrain data: keep float32 for continuous values
-            profile.update(dtype="float32")
-
-        cog_translate(str(tif_file), str(cog_file), profile, quiet=True)
-        logger.debug(f"    COG créé: {cog_file.name}")
-        return cog_file
-
-    except ImportError:
-        logger.warning("    rio-cogeo non disponible — COG non généré")
-        return None
-    except Exception as e:
-        logger.warning(f"    Erreur conversion COG: {e}")
-        return None
-
-
-def generate_cog_metadata(vis_dir, basename):
-    """Generate metadata JSON for all visualization layers.
-
-    Scans the COG directory for files and reads their geographic bounds.
-    Creates a metadata.json with WGS84 bounds and layer info for the web viewer.
-
-    Args:
-        vis_dir: Directory containing COG files (vis_dir/basename/cog/).
-        basename: Base name for the LiDAR tile.
-
-    Returns:
-        Path to metadata.json, or None on failure.
-    """
-    import json
-
-    cog_dir = vis_dir / basename / "cog"
-    if not cog_dir.exists():
-        return None
-
-    # Find the DTM to get geographic bounds
-    # The COG files inherit bounds from the DTM, so we can read any COG
-    cog_files = sorted(cog_dir.glob("*.tif"))
-    if not cog_files:
-        return None
-
-    # Read bounds from first COG file
-    try:
-        ref_cog = cog_files[0]
-        with rasterio.open(ref_cog) as src:
-            bounds = src.bounds
-            crs = src.crs
-            if crs and HAS_WARP:
-                l93_xs = [bounds.left, bounds.right]
-                l93_ys = [bounds.top, bounds.bottom]
-                lons, lats = warp_transform(crs, 'EPSG:4326', l93_xs, l93_ys)
-                bounds_wgs84 = {
-                    'west': float(min(lons)),
-                    'south': float(min(lats)),
-                    'east': float(max(lons)),
-                    'north': float(max(lats)),
-                }
-            else:
-                # Fallback: use Lambert 93 bounds directly
-                bounds_wgs84 = {
-                    'west': float(bounds.left),
-                    'south': float(bounds.bottom),
-                    'east': float(bounds.right),
-                    'north': float(bounds.top),
-                }
-
-            bounds_l93 = {
-                'min_x': float(bounds.left),
-                'min_y': float(bounds.bottom),
-                'max_x': float(bounds.right),
-                'max_y': float(bounds.top),
-            }
-            resolution = float(abs(src.transform.a))
-    except Exception as e:
-        logger.warning(f"    Erreur lecture bounds COG: {e}")
-        return None
-
-    # Build layer list with titles
-    layers = []
-    for cog_file in cog_files:
-        name = cog_file.stem.replace(basename + "_", "")
-        # Find matching title from COLORMAPS or RGB_LEGENDS
-        title = name.replace("_", " ").title()
-        for key in sorted(COLORMAPS.keys(), key=len, reverse=True):
-            if key in name:
-                title = COLORMAPS[key]['title']
-                break
-        for key in RGB_LEGENDS:
-            if key in name:
-                title = RGB_LEGENDS[key]['title']
-                break
-
-        layers.append({
-            'name': name,
-            'title': title,
-            'file': cog_file.name,
-        })
-
-    metadata = {
-        'basename': basename,
-        'bounds_l93': bounds_l93,
-        'bounds_wgs84': bounds_wgs84,
-        'resolution': resolution,
-        'layers': layers,
-    }
-
-    metadata_file = vis_dir / basename / "metadata.json"
-    with open(metadata_file, 'w') as f:
-        json.dump(metadata, f, indent=2, ensure_ascii=False)
-
-    logger.info(f"    Métadonnées: {len(layers)} couches, bounds={bounds_wgs84}")
-    return metadata_file
-
 
 def generate_pdf_report(basename, vis_dir, pdf_dir, resolution):
     """Generate A3 PDF report for a LiDAR file with all visualizations.

lidar_pipeline/server.py

@@ -1,117 +0,0 @@
-"""TiTiler-based web server for serving LiDAR COG visualizations.
-
-Provides:
-- COG tile serving via TiTiler API
-- Static file serving for viewer HTML
-- CORS headers for local development
-
-Usage:
-    python -m lidar_pipeline.server /path/to/output [--port 8000]
-"""
-
-import logging
-import sys
-from pathlib import Path
-
-logger = logging.getLogger("lidar")
-
-
-def create_app(output_dir):
-    """Create the FastAPI application with TiTiler and static file serving.
-
-    Args:
-        output_dir: Path to the output directory containing visualisations/ and viewer/.
-
-    Returns:
-        FastAPI application instance.
-    """
-    from fastapi import FastAPI
-    from fastapi.responses import HTMLResponse, FileResponse, JSONResponse
-    from fastapi.middleware.cors import CORSMiddleware
-    from titiler.core.factory import TilerFactory
-
-    output_dir = Path(output_dir).resolve()
-
-    # TiTiler COG endpoint
-    cog = TilerFactory(router_prefix="/cog")
-
-    app = FastAPI(title="LiDAR Archéologique", docs_url="/docs")
-
-    # CORS for local development
-    app.add_middleware(CORSMiddleware, allow_origins=["*"], allow_methods=["*"], allow_headers=["*"])
-
-    # Mount TiTiler routes
-    app.include_router(cog.router, prefix="/cog")
-
-    @app.get("/")
-    async def index():
-        """Landing page with links."""
-        return HTMLResponse(
-            '<html><head><title>LiDAR Server</title>'
-            '<style>body{font-family:sans-serif;max-width:600px;margin:40px auto;padding:0 20px}'
-            'h1{color:#1a1a2e}a{color:#4fc3f7}</style></head>'
-            '<body>'
-            '<h1>Serveur LiDAR Archéologique</h1>'
-            '<p><a href="/viewer">Visualisations</a></p>'
-            '<p>TiTiler API: <a href="/cog/">/cog/</a></p>'
-            '</body></html>'
-        )
-
-    @app.get("/viewer/{basename}")
-    @app.get("/viewer")
-    async def serve_viewer(basename: str = ""):
-        """Serve viewer HTML files."""
-        if not basename:
-            # List available viewers
-            viewer_dir = output_dir / 'visualisations' / 'viewer'
-            if viewer_dir.exists():
-                viewers = sorted(viewer_dir.glob('*.html'))
-                if viewers:
-                    links = ''.join(
-                        f'<li><a href="/viewer/{f.stem}">{f.stem}</a></li>'
-                        for f in viewers
-                    )
-                    return HTMLResponse(
-                        f'<html><head><title>LiDAR Viewer</title>'
-                        f'<style>body{{font-family:sans-serif;max-width:600px;margin:40px auto;padding:0 20px}}'
-                        f'h1{{color:#1a1a2e}}li{{margin:8px 0}}a{{color:#4fc3f7}}</style></head>'
-                        f'<body><h1>Visualisations LiDAR</h1><ul>{links}</ul></body></html>'
-                    )
-            return HTMLResponse('<html><body><p>Aucun viewer disponible</p></body></html>')
-
-        viewer_file = output_dir / 'visualisations' / 'viewer' / f"{basename}.html"
-        if viewer_file.exists():
-            return FileResponse(str(viewer_file), media_type='text/html')
-        return JSONResponse({'error': f'Viewer not found: {basename}'}, status_code=404)
-
-    return app
-
-
-def main():
-    """Entry point for the LiDAR web server."""
-    import argparse
-    import uvicorn
-
-    parser = argparse.ArgumentParser(description='Serveur cartographique LiDAR')
-    parser.add_argument('output_dir', help='Répertoire de sortie contenant les visualisations')
-    parser.add_argument('--host', default='0.0.0.0', help='Hôte (défaut: 0.0.0.0)')
-    parser.add_argument('--port', type=int, default=8000, help='Port (défaut: 8000)')
-    args = parser.parse_args()
-
-    output_dir = Path(args.output_dir)
-    if not output_dir.exists():
-        logger.error(f"Répertoire introuvable: {output_dir}")
-        sys.exit(1)
-
-    app = create_app(output_dir)
-
-    print(f"Serveur LiDAR Archéologique")
-    print(f"  Viewer:  http://{args.host}:{args.port}/viewer")
-    print(f"  TiTiler: http://{args.host}:{args.port}/cog/")
-    print(f"  Répertoire: {output_dir}")
-
-    uvicorn.run(app, host=args.host, port=args.port, log_level='info')
-
-
-if __name__ == '__main__':
-    main()

lidar_pipeline/tests/test_visualizations.py

@@ -198,4 +198,35 @@ class TestFlow:
             data = src.read(1)
             # log1p(x) >= 0 for x >= 0
             valid = data[~np.isnan(data)]
-            assert np.nanmin(valid) >= 0
+            assert np.nanmin(valid) >= 0
+
+
+class TestLocalDominance:
+    def test_generates_tif(self, synthetic_dem, tmp_output_dir):
+        from lidar_pipeline.visualizations import generate_local_dominance
+        result = generate_local_dominance(synthetic_dem, "test", tmp_output_dir, 5.0)
+        assert result is not None
+        assert result.exists()
+        assert result.suffix == ".tif"
+
+    def test_dominance_values_0_1(self, synthetic_dem, tmp_output_dir):
+        import rasterio
+        from lidar_pipeline.visualizations import generate_local_dominance
+        result = generate_local_dominance(synthetic_dem, "test", tmp_output_dir, 5.0)
+        with rasterio.open(result) as src:
+            data = src.read(1)
+            valid = data[~np.isnan(data)]
+            assert np.nanmin(valid) >= 0, "Local dominance should be >= 0"
+            assert np.nanmax(valid) <= 1, "Local dominance should be <= 1"
+
+    def test_dominance_nan_mask_preserved(self, synthetic_dem, tmp_output_dir):
+        """Check that NaN zones from original DEM are preserved."""
+        import rasterio
+        from lidar_pipeline.visualizations import generate_local_dominance
+        result = generate_local_dominance(synthetic_dem, "test", tmp_output_dir, 5.0)
+        # The synthetic DEM has no NaN, so this just verifies the output is valid
+        with rasterio.open(result) as src:
+            data = src.read(1)
+            # Shape should match input
+            assert data.shape[0] > 0
+            assert data.shape[1] > 0

lidar_pipeline/viewer.py

@@ -1,416 +0,0 @@
-"""Web map viewer generator for LiDAR visualizations.
-
-Generates a self-contained HTML file with MapLibre GL JS that displays
-all visualization layers with opacity controls and IGN/OSM basemaps.
-Served by TiTiler for COG tile access.
-"""
-
-import json
-import logging
-from pathlib import Path
-
-logger = logging.getLogger("lidar")
-
-# Layer ordering for the viewer panel (archaeological priority)
-LAYER_ORDER = [
-    'hillshade_multi',
-    'slope',
-    'aspect',
-    'curvature',
-    'svf',
-    'lrm',
-    'positive_openness',
-    'negative_openness',
-    'mslrm',
-    'tpi',
-    'sailore',
-    'roughness',
-    'anomalies',
-    'wavelet',
-    'flow',
-    'ortho',
-    'topo',
-]
-
-# Colormaps for TiTiler rendering of single-band COGs
-LAYER_COLORMAPS = {
-    'hillshade_multi': 'gray',
-    'slope': 'inferno',
-    'aspect': 'turbo',
-    'curvature': 'RdYlBu_r',
-    'svf': 'viridis',
-    'lrm': 'RdBu_r',
-    'positive_openness': 'YlOrBr',
-    'negative_openness': 'GnBu_r',
-    'mslrm': 'RdBu_r',
-    'tpi': 'BrBG',
-    'sailore': 'seismic',
-    'roughness': 'magma',
-    'anomalies': 'coolwarm',
-    'wavelet': 'cividis',
-    'flow': 'Blues',
-}
-
-
-def generate_viewer(basename, vis_dir, output_vis_dir, titiler_url='http://localhost:8000'):
-    """Generate a MapLibre GL JS viewer HTML file for the LiDAR tile.
-
-    Args:
-        basename: Base name for the LiDAR tile.
-        vis_dir: Per-file visualization directory (vis_dir/basename/).
-        output_vis_dir: Parent visualization directory for viewer output.
-        titiler_url: Base URL of the TiTiler server.
-
-    Returns:
-        Path to the generated HTML file, or None on failure.
-    """
-    # Read metadata.json
-    metadata_file = vis_dir / "metadata.json"
-    if not metadata_file.exists():
-        logger.error(f"    Métadonnées manquantes: {metadata_file}")
-        return None
-
-    with open(metadata_file) as f:
-        metadata = json.load(f)
-
-    bounds_wgs84 = metadata['bounds_wgs84']
-    resolution = metadata.get('resolution', 0.5)
-
-    # Determine which files are RGB (ortho/topo)
-    layers = []
-    for layer in metadata['layers']:
-        is_rgb = 'ortho' in layer['name'] or 'topo' in layer['name']
-        layers.append({
-            'name': layer['name'],
-            'title': layer['title'],
-            'file': layer['file'],
-            'is_rgb': is_rgb,
-        })
-
-    # Sort layers by archaeological priority
-    def layer_sort_key(l):
-        name = l['name']
-        for i, key in enumerate(LAYER_ORDER):
-            if key in name:
-                return i
-        return len(LAYER_ORDER)
-
-    layers.sort(key=layer_sort_key)
-
-    # COG path prefix for TiTiler (absolute path inside Docker container)
-    cog_path_prefix = f'/data/output/visualisations/{basename}/cog'
-
-    # Build layer controls HTML
-    layer_controls = []
-    for layer in layers:
-        checked = 'checked' if layer['name'] == 'hillshade_multi' else ''
-        initial_opacity = 100 if layer['name'] == 'hillshade_multi' else 0
-        layer_type = 'RGB' if layer['is_rgb'] else 'DEM'
-        layer_controls.append(
-            f'      <div class="layer-control">\n'
-            f'        <label>\n'
-            f'          <input type="checkbox" data-layer="{layer["name"]}" {checked}>\n'
-            f'          <span class="layer-title">{layer["title"]}</span>\n'
-            f'          <span class="layer-type">{layer_type}</span>\n'
-            f'        </label>\n'
-            f'        <input type="range" class="opacity-slider" data-layer="{layer["name"]}" '
-            f'min="0" max="100" value="{initial_opacity}">\n'
-            f'      </div>'
-        )
-    layer_controls_html = '\n'.join(layer_controls)
-
-    # Serialize data for JS
-    bounds_json = json.dumps(bounds_wgs84)
-    layers_json = json.dumps(layers, ensure_ascii=False)
-    colormaps_json = json.dumps(LAYER_COLORMAPS)
-    center_lng = (bounds_wgs84['west'] + bounds_wgs84['east']) / 2
-    center_lat = (bounds_wgs84['south'] + bounds_wgs84['north']) / 2
-
-    # Build the full HTML using string concatenation to avoid f-string issues with {z}/{x}/{y}
-    html_parts = []
-    html_parts.append(f'''<!DOCTYPE html>
-<html lang="fr">
-<head>
-<meta charset="utf-8">
-<meta name="viewport" content="width=device-width, initial-scale=1">
-<title>LiDAR Archéologique — {basename}</title>
-<link rel="stylesheet" href="https://unpkg.com/maplibre-gl@4.7.1/dist/maplibre-gl.css">
-<script src="https://unpkg.com/maplibre-gl@4.7.1/dist/maplibre-gl.js"></script>
-<style>
-  * {{ margin: 0; padding: 0; box-sizing: border-box; }}
-  body {{ font-family: -apple-system, BlinkMacSystemFont, 'Segoe UI', Roboto, sans-serif; }}
-  #map {{ position: absolute; top: 0; left: 280px; right: 0; bottom: 0; }}
-  #panel {{
-    position: absolute; top: 0; left: 0; width: 280px; bottom: 0;
-    background: #1a1a2e; color: #e0e0e0; overflow-y: auto;
-    z-index: 10; box-shadow: 2px 0 8px rgba(0,0,0,0.3);
-  }}
-  #panel h1 {{
-    font-size: 14px; padding: 12px 14px; background: #16213e;
-    border-bottom: 1px solid #0f3460; letter-spacing: 0.5px;
-  }}
-  #panel h2 {{
-    font-size: 11px; padding: 8px 14px; color: #a0a0a0;
-    text-transform: uppercase; letter-spacing: 1px;
-    border-bottom: 1px solid #2a2a4a;
-  }}
-  .layer-group {{ padding: 4px 0; border-bottom: 1px solid #2a2a4a; }}
-  .layer-control {{
-    padding: 6px 14px;
-    transition: background 0.15s;
-  }}
-  .layer-control:hover {{ background: #252550; }}
-  .layer-control label {{
-    display: flex; align-items: center; gap: 8px;
-    font-size: 13px; cursor: pointer; margin-bottom: 4px;
-  }}
-  .layer-control input[type="checkbox"] {{
-    width: 16px; height: 16px; accent-color: #4fc3f7;
-  }}
-  .layer-control .layer-title {{ flex: 1; }}
-  .layer-control .layer-type {{
-    font-size: 10px; color: #888; text-transform: uppercase;
-  }}
-  .opacity-slider {{
-    width: 100%; height: 4px; margin: 2px 0 0 24px;
-    -webkit-appearance: none; appearance: none;
-    background: #333; border-radius: 2px; outline: none;
-  }}
-  .opacity-slider::-webkit-slider-thumb {{
-    -webkit-appearance: none; appearance: none;
-    width: 14px; height: 14px; border-radius: 50%;
-    background: #4fc3f7; cursor: pointer;
-  }}
-  .basemap-section {{
-    padding: 10px 14px; border-bottom: 1px solid #2a2a4a;
-  }}
-  .basemap-section h3 {{
-    font-size: 11px; color: #a0a0a0; text-transform: uppercase;
-    letter-spacing: 1px; margin-bottom: 8px;
-  }}
-  .basemap-btn {{
-    display: inline-block; padding: 5px 10px; margin: 2px;
-    background: #2a2a4a; color: #ccc; border: 1px solid #444;
-    border-radius: 4px; cursor: pointer; font-size: 12px;
-    transition: all 0.15s;
-  }}
-  .basemap-btn:hover {{ background: #3a3a6a; }}
-  .basemap-btn.active {{ background: #4fc3f7; color: #000; border-color: #4fc3f7; }}
-  #coords {{
-    position: absolute; bottom: 8px; left: 288px;
-    background: rgba(26,26,46,0.9); color: #4fc3f7;
-    padding: 4px 10px; border-radius: 4px; font-size: 12px;
-    font-family: monospace; z-index: 5;
-  }}
-</style>
-</head>
-<body>
-
-<div id="panel">
-  <h1>LiDAR Archéologique</h1>
-  <div style="padding: 6px 14px; font-size: 12px; color: #888;">
-    {basename}<br>
-    <span style="font-size: 10px;">Res: {resolution}m/px | EPSG:2154</span>
-  </div>
-
-  <div class="basemap-section">
-    <h3>Fond de carte</h3>
-    <button class="basemap-btn active" onclick="setBasemap('osm')">OSM</button>
-    <button class="basemap-btn" onclick="setBasemap('ign-photo')">Photo IGN</button>
-    <button class="basemap-btn" onclick="setBasemap('ign-map')">Carte IGN</button>
-  </div>
-
-  <h2>Visualisations</h2>
-  <div class="layer-group" id="layers">
-{layer_controls_html}
-  </div>
-</div>
-
-<div id="map"></div>
-<div id="coords"></div>
-
-<script>
-const TITILER_URL = '{titiler_url}';
-const BASENAME = '{basename}';
-const BOUNDS = {bounds_json};
-const LAYERS = {layers_json};
-const COG_PATH_PREFIX = '{cog_path_prefix}';
-const LAYER_COLORMAPS = {colormaps_json};''')
-
-    # JavaScript code — use regular string to avoid f-string issues with {z}/{x}/{y}
-    html_parts.append(r'''
-// Map initialization
-const map = new maplibregl.Map({
-  container: 'map',
-  style: {
-    version: 8,
-    sources: {},
-    layers: []
-  },
-  center: [' + f'{center_lng}, {center_lat}' + r'],
-  zoom: 13,
-  maxZoom: 19,
-  minZoom: 8,
-});
-
-// Basemap layers
-const basemapStyles = {
-  'osm': {
-    type: 'raster',
-    tiles: ['https://tile.openstreetmap.org/{z}/{x}/{y}.png'],
-    tileSize: 256,
-    attribution: '&copy; OpenStreetMap',
-    maxzoom: 19
-  },
-  'ign-photo': {
-    type: 'raster',
-    tiles: ['https://data.geopf.fr/wmts?SERVICE=WMTS&VERSION=1.0.0&REQUEST=GetTile&LAYER=ORTHOIMAGERY.ORTHOPHOTOS&STYLE=normal&TILEMATRIXSET=PM&TILEMATRIX={z}&TILECOL={x}&TILEROW={y}&FORMAT=image/jpeg'],
-    tileSize: 256,
-    attribution: '&copy; IGN',
-    maxzoom: 20
-  },
-  'ign-map': {
-    type: 'raster',
-    tiles: ['https://data.geopf.fr/wmts?SERVICE=WMTS&VERSION=1.0.0&REQUEST=GetTile&LAYER=GEOGRAPHICALGRIDSYSTEMS.PLANIGNV2&STYLE=normal&TILEMATRIXSET=PM&TILEMATRIX={z}&TILECOL={x}&TILEROW={y}&FORMAT=image/png'],
-    tileSize: 256,
-    attribution: '&copy; IGN',
-    maxzoom: 19
-  }
-};
-
-let currentBasemap = 'osm';
-
-function setBasemap(name) {
-  if (map.getLayer('basemap')) map.removeLayer('basemap');
-  if (map.getSource('basemap')) map.removeSource('basemap');
-
-  currentBasemap = name;
-  const style = basemapStyles[name];
-  map.addSource('basemap', style);
-  map.addLayer({ id: 'basemap', type: 'raster', source: 'basemap' });
-
-  document.querySelectorAll('.basemap-btn').forEach(btn => {
-    const names = {'osm': 'OSM', 'ign-photo': 'Photo IGN', 'ign-map': 'Carte IGN'};
-    btn.classList.toggle('active', btn.textContent.trim() === names[name]);
-  });
-
-  // Re-add data layers on top of basemap
-  LAYERS.forEach(layer => {
-    const layerId = 'layer-' + layer.name;
-    if (map.getLayer(layerId)) {
-      map.moveLayer(layerId, 'basemap');
-    }
-  });
-}
-
-// Add visualization layers
-function addVisualizationLayer(layerInfo) {
-  const name = layerInfo.name;
-  const isRgb = layerInfo.is_rgb;
-  const cogUrl = COG_PATH_PREFIX + '/' + layerInfo.file;
-
-  let tileUrl;
-  if (isRgb) {
-    tileUrl = TITILER_URL + '/cog/tiles/{z}/{x}/{y}.png?url=' +
-      encodeURIComponent(cogUrl) + '&bidx=1&bidx=2&bidx=3';
-  } else {
-    const cmap = LAYER_COLORMAPS[name] || 'viridis';
-    tileUrl = TITILER_URL + '/cog/tiles/{z}/{x}/{y}.png?url=' +
-      encodeURIComponent(cogUrl) + '&colormap_name=' + cmap + '&rescale=-1,1';
-  }
-
-  const sourceId = 'source-' + name;
-  const layerId = 'layer-' + name;
-
-  map.addSource(sourceId, {
-    type: 'raster',
-    tiles: [tileUrl],
-    tileSize: 256,
-    minzoom: 8,
-    maxzoom: 19,
-  });
-
-  map.addLayer({
-    id: layerId,
-    type: 'raster',
-    source: sourceId,
-    paint: {
-      'raster-opacity': name === 'hillshade_multi' ? 1.0 : 0.0,
-      'raster-resampling': 'bilinear',
-    }
-  });
-}
-
-// Layer control event handlers
-document.querySelectorAll('.layer-control input[type="checkbox"]').forEach(cb => {
-  cb.addEventListener('change', () => {
-    const layerId = 'layer-' + cb.dataset.layer;
-    if (map.getLayer(layerId)) {
-      map.setLayoutProperty(layerId, 'visibility',
-        cb.checked ? 'visible' : 'none');
-    }
-    // When checking a layer, set its opacity to 70% if it was 0
-    const slider = document.querySelector('.opacity-slider[data-layer="' + cb.dataset.layer + '"]');
-    if (cb.checked && slider && parseInt(slider.value) === 0) {
-      slider.value = 70;
-      const layerId = 'layer-' + cb.dataset.layer;
-      if (map.getLayer(layerId)) {
-        map.setPaintProperty(layerId, 'raster-opacity', 0.7);
-      }}
-  });
-});
-
-document.querySelectorAll('.opacity-slider').forEach(slider => {
-  slider.addEventListener('input', () => {
-    const layerId = 'layer-' + slider.dataset.layer;
-    if (map.getLayer(layerId)) {
-      map.setPaintProperty(layerId, 'raster-opacity', slider.value / 100);
-    }
-  });
-});
-
-// Coordinate display
-map.on('mousemove', (e) => {
-  const { lng, lat } = e.lngLat;
-  document.getElementById('coords').textContent =
-    lat.toFixed(6) + 'N  ' + lng.toFixed(6) + 'E';
-});
-
-// Fit map to bounds on load
-map.on('load', () => {
-  setBasemap('osm');
-
-  // Add all visualization layers
-  LAYERS.forEach(layer => {
-    addVisualizationLayer(layer);
-  });
-
-  // Initially hide all except hillshade
-  LAYERS.forEach(layer => {
-    if (layer.name !== 'hillshade_multi') {
-      const layerId = 'layer-' + layer.name;
-      if (map.getLayer(layerId)) {
-        map.setLayoutProperty(layerId, 'visibility', 'none');
-      }
-    }
-  });
-
-  // Fit bounds
-  map.fitBounds([[BOUNDS.west, BOUNDS.south], [BOUNDS.east, BOUNDS.north]], { padding: 20 });
-});
-</script>
-</body>
-</html>''')
-
-    html = ''.join(html_parts)
-
-    # Write viewer HTML
-    viewer_dir = output_vis_dir / 'viewer'
-    viewer_dir.mkdir(exist_ok=True)
-    viewer_file = viewer_dir / f"{basename}.html"
-
-    with open(viewer_file, 'w', encoding='utf-8') as f:
-        f.write(html)
-
-    logger.info(f"    Viewer: {viewer_file}")
-    return viewer_file

run.sh

@@ -3,21 +3,22 @@
 # Utilisation: ./run.sh [options]
 #
 # Options:
-#   -r RESOLUTION   Résolution en m/px (défaut: 0.5)
+#   -r RESOLUTION   Résolution en m/px, ou multiples séparées par virgules (défaut: 0.5, ex: 0.5,0.2)
 #   -w WORKERS      Nombre de workers parallèles (défaut: 1)
 #   -g              Activer l'accélération GPU
-#   -v              Mode verbeux (timestamps + niveaux)
-#   --debug         Mode debug (détails internes fichier:ligne)
-#   -f / --force    Régénérer tous les fichiers même si existants
-#   --keep-tif       Conserver les fichiers TIFF intermédiaires
-#   --no-viewer     Ne pas générer le viewer web (COGs + HTML)
+#   -v              Mode verbeux
+#   --debug         Mode debug
+#   -f / --force    Régénérer tous les fichiers
+#   --keep-tif      Conserver les fichiers TIFF
 #   --force-classification
-#                   Reclassifier le sol même si le fichier .las existe déjà
-#   --ground-classification {auto,smrf,pmf,csf}
-#                   Méthode de classification du sol (défaut: auto)
-#   --file NOM...   Traiter un ou plusieurs fichiers LAZ spécifiques
+#   --ground-classification {auto,smrf,csf}
+#   --quality N     Qualité image 1-100 (défaut: 98)
+#   --lossless      Compression lossless
+#   --format FMT    Format de sortie : avif (défaut) ou webp
+#   --only VIZ...   Générer uniquement ces visualisations
+#   --skip VIZ...   Exclure ces visualisations
+#   --file NOM...   Traiter un ou plusieurs fichiers LAZ
 #   --test          Exécuter les tests unitaires
-#   serve           Démarrer le serveur cartographique
 #   -h              Afficher l'aide complète
 
 set -e
@@ -27,36 +28,14 @@ INPUT_DIR="${SCRIPT_DIR}/input"
 OUTPUT_DIR="${SCRIPT_DIR}/output"
 IMAGE_NAME="lidar-lidar"
 
-# Serve command — start the web map server
-if [ "$1" = "serve" ]; then
-    # Build l'image si elle n'existe pas
-    if ! docker image inspect "$IMAGE_NAME" >/dev/null 2>&1; then
-        echo "Build de l'image Docker..."
-        docker build -t "$IMAGE_NAME" "$SCRIPT_DIR"
-    fi
-    mkdir -p "$OUTPUT_DIR"
-    echo "============================================"
-    echo " Serveur cartographique LiDAR"
-    echo "============================================"
-    echo "  Viewer:  http://localhost:8000/viewer"
-    echo "  TiTiler: http://localhost:8000/cog/"
-    echo "============================================"
-    docker run --rm -p 8000:8000 \
-        -v "${OUTPUT_DIR}:/data/output" \
-        "$IMAGE_NAME" \
-        python3 -m lidar_pipeline.server /data/output
-    exit 0
-fi
-
 # Afficher l'aide si aucun argument
 if [ $# -eq 0 ]; then
     echo "Pipeline LiDAR Archéologique"
     echo ""
     echo "Usage: $0 [options]"
-    echo "  $0 serve                                  # Démarrer le serveur cartographique"
     echo ""
     echo "Options:"
-    echo "  -r RESOLUTION   Résolution en m/px (défaut: 0.5)"
+    echo "  -r RESOLUTION   Résolution en m/px, ou multiples (défaut: 0.5, ex: 0.5,0.2)"
     echo "  -w WORKERS      Nombre de workers CPU parallèles (défaut: 1)"
     echo "  -g              Activer l'accélération GPU NVIDIA"
     echo "  -v              Mode verbeux (timestamps + niveaux)"
@@ -64,13 +43,11 @@ if [ $# -eq 0 ]; then
     echo "  -f / --force    Régénérer tous les fichiers même si les WebP existent"
     echo "  --force-classification"
     echo "                  Reclassifier le sol même si le fichier .las existe"
-    echo "  --keep-tif      Conserver les fichiers TIFF intermédiaires"
-    echo "  --no-viewer     Ne pas générer le viewer web"
-    echo "  --ground-classification {auto,smrf,pmf,csf}"
+    echo "  --keep-tif      Conserver les TIFF pour régénérer les WebP"
+    echo "  --ground-classification {auto,smrf,csf}"
     echo "                  Méthode de classification du sol (défaut: auto)"
     echo "  --file NOM...   Traiter un ou plusieurs fichiers LAZ (nom complet sans .laz/.las)"
     echo "  --test          Exécuter les tests unitaires"
-    echo "  serve           Démarrer le serveur cartographique"
     echo "  -h              Afficher cette aide"
     echo ""
     echo "Exemples:"
@@ -78,11 +55,11 @@ if [ $# -eq 0 ]; then
     echo "  $0 -g -w 4                               # GPU + 4 workers"
     echo "  $0 -g -v                                 # GPU + verbeux"
     echo "  $0 -g -r 0.2                             # Haute résolution"
-    echo "  $0 -g --force                            # Tout régénérer (WebP + classification)"
+    echo "  $0 -g -r 0.5,0.2                          # Multi-résolution (0.5m + 0.2m)"
+    echo "  $0 -g --force                            # Régénérer WebP (DTM conservé si --keep-tif)"
     echo "  $0 -g --force-classification             # Reclassifier le sol seulement"
-    echo "  $0 -g --ground-classification pmf        # Forcer PMF"
+    echo "  $0 -g --ground-classification csf        # Forcer CSF (terrain complexe)"
     echo "  $0 -g --file LHD_...IGN69.copc           # Un fichier"
-    echo "  $0 serve                                  # Démarrer le serveur web"
     exit 0
 fi
 
@@ -95,7 +72,10 @@ FILE_ARGS=""
 GROUND_METHOD=""
 FORCE_CLASSIFY_FLAG=""
 KEEP_TIF_FLAG=""
-NO_VIEWER_FLAG=""
+QUALITY=""
+FORMAT_FLAG=""
+ONLY_FLAG=""
+SKIP_FLAG=""
 
 # Parse arguments manually (more robust than getopts for mixed short/long options)
 while [ $# -gt 0 ]; do
@@ -109,16 +89,19 @@ while [ $# -gt 0 ]; do
         --force) FORCE_FLAG="--force"; shift ;;
         --force-classification) FORCE_CLASSIFY_FLAG="--force-classification"; shift ;;
         --keep-tif) KEEP_TIF_FLAG="--keep-tif"; shift ;;
-        --no-viewer) NO_VIEWER_FLAG="--no-viewer"; shift ;;
         --ground-classification) GROUND_METHOD="$2"; shift 2 ;;
         --ground-classification=*) GROUND_METHOD="${1#--ground-classification=}"; shift ;;
+        --quality) QUALITY="--quality $2"; shift 2 ;;
+        --lossless) QUALITY="--lossless"; shift ;;
+        --format) FORMAT_FLAG="--format $2"; shift 2 ;;
+        --only) shift; ONLY_FLAG="--only"; while [ $# -gt 0 ] && [[ ! "$1" =~ ^- ]]; do ONLY_FLAG="$ONLY_FLAG $1"; shift; done ;;
+        --skip) shift; SKIP_FLAG="--skip"; while [ $# -gt 0 ] && [[ ! "$1" =~ ^- ]]; do SKIP_FLAG="$SKIP_FLAG $1"; shift; done ;;
         --file) shift; while [ $# -gt 0 ] && [[ ! "$1" =~ ^- ]]; do FILE_ARGS="$FILE_ARGS $1"; shift; done ;;
         --test) ;;  # Handled below
         -h|--help|-help)
             echo "Pipeline LiDAR Archéologique"
             echo ""
             echo "Usage: $0 [options]"
-            echo "  $0 serve                                  # Démarrer le serveur cartographique"
             echo ""
             echo "Options:"
             echo "  -r RESOLUTION   Résolution en m/px (défaut: 0.5)"
@@ -129,25 +112,35 @@ while [ $# -gt 0 ]; do
             echo "  -f / --force    Régénérer tous les fichiers même si les WebP existent"
             echo "  --force-classification"
             echo "                  Reclassifier le sol même si le fichier .las existe"
-            echo "  --keep-tif      Conserver les fichiers TIFF intermédiaires"
-            echo "  --no-viewer     Ne pas générer le viewer web"
-            echo "  --ground-classification {auto,smrf,pmf,csf}"
+            echo "  --keep-tif      Conserver les TIFF pour régénérer les WebP"
+            echo "  --ground-classification {auto,smrf,csf}"
             echo "                  Méthode de classification du sol (défaut: auto)"
-            echo "  --file NOM...   Traiter un ou plusieurs fichiers LAZ (nom complet sans .laz/.las)"
+            echo "  --quality N     Qualité image 1-100 (défaut: 98, 100=lossless)"
+            echo "  --lossless      Compression lossless (équivalent à --quality 100)"
+            echo "  --format FMT    Format de sortie : avif (défaut) ou webp"
+            echo "  --only VIZ...   Générer uniquement ces visualisations"
+            echo "  --skip VIZ...   Exclure ces visualisations"
+            echo "  --file NOM...   Traiter un ou plusieurs fichiers LAZ"
             echo "  --test          Exécuter les tests unitaires"
-            echo "  serve           Démarrer le serveur cartographique"
             echo "  -h              Afficher cette aide"
             echo ""
+            echo "Visualisations disponibles:"
+            echo "  hillshade slope aspect curvature lrm pos_open neg_open"
+            echo "  mslrm tpi sailore roughness anomalies wavelet flow ortho topo"
+            echo ""
             echo "Exemples:"
             echo "  $0 -g                                    # GPU, auto"
             echo "  $0 -g -w 4                               # GPU + 4 workers"
             echo "  $0 -g -v                                 # GPU + verbeux"
             echo "  $0 -g -r 0.2                             # Haute résolution"
-            echo "  $0 -g --force                            # Tout régénérer (WebP + classification)"
-            echo "  $0 -g --force-classification             # Reclassifier le sol seulement"
-            echo "  $0 -g --ground-classification pmf        # Forcer PMF"
+    echo "  $0 -g -r 0.5,0.2                          # Multi-résolution (0.5m + 0.2m)"
+            echo "  $0 -g --force                            # Régénérer WebP"
+            echo "  $0 -g --only hillshade svf lrm           # Seulement 3 visualisations"
+            echo "  $0 -g --skip ortho topo                  # Sans les overlays IGN"
+            echo "  $0 -g --lossless                         # WebP lossless"
+            echo "  $0 -g --quality 90                        # WebP qualité 90"
+            echo "  $0 -g --ground-classification csf        # Forcer CSF (terrain complexe)"
             echo "  $0 -g --file LHD_...IGN69.copc           # Un fichier"
-            echo "  $0 serve                                  # Démarrer le serveur web"
             exit 0
             ;;
         *) echo "Option invalide: $1" >&2; exit 1 ;;
@@ -167,7 +160,7 @@ if [[ " $* " == *" --test "* ]]; then
     echo "============================================"
     echo " Tests unitaires LiDAR Pipeline"
     echo "============================================"
-    docker run --rm $GPU_FLAG \
+    docker run --rm --init $GPU_FLAG \
         "$IMAGE_NAME" \
         python3 -m pytest -v --pyargs lidar_pipeline.tests
     exit $?
@@ -193,22 +186,34 @@ echo "  Verbeux    : $([ -n "$VERBOSE_FLAG" ] && echo 'OUI' || echo 'non')"
 echo "  Force      : $([ -n "$FORCE_FLAG" ] && echo 'OUI' || echo 'non')"
 echo "  Force classif.: $([ -n "$FORCE_CLASSIFY_FLAG" ] && echo 'OUI' || echo 'non')"
 echo "  Keep TIFF    : $([ -n "$KEEP_TIF_FLAG" ] && echo 'OUI' || echo 'non')"
-echo "  Viewer web   : $([ -n "$NO_VIEWER_FLAG" ] && echo 'non' || echo 'OUI')"
+echo "  Qualité image : $([ -n "$QUALITY" ] && echo "$QUALITY" || echo '98')"
+echo "  Format      : $([ -n "$FORMAT_FLAG" ] && echo "${FORMAT_FLAG#--format }" || echo 'avif')"
 echo "  Classification sol : $([ -n "$GROUND_METHOD" ] && echo "$GROUND_METHOD" || echo 'auto')"
+if [ -n "$ONLY_FLAG" ]; then
+    echo "  Visualisations: uniquement${ONLY_FLAG#--only}"
+elif [ -n "$SKIP_FLAG" ]; then
+    echo "  Visualisations: tout sauf${SKIP_FLAG#--skip}"
+fi
 if [ -n "$FILE_ARGS" ]; then
     echo "  Fichiers   :${FILE_ARGS}"
 fi
 echo "============================================"
 
-CMD_ARGS="-o /data/output -r $RESOLUTION -w $WORKERS $VERBOSE_FLAG $FORCE_FLAG $FORCE_CLASSIFY_FLAG $KEEP_TIF_FLAG $NO_VIEWER_FLAG"
+CMD_ARGS="-o /data/output -r $RESOLUTION -w $WORKERS $VERBOSE_FLAG $FORCE_FLAG $FORCE_CLASSIFY_FLAG $KEEP_TIF_FLAG $QUALITY $FORMAT_FLAG"
 if [ -n "$GROUND_METHOD" ]; then
     CMD_ARGS="$CMD_ARGS --ground-classification $GROUND_METHOD"
 fi
+if [ -n "$ONLY_FLAG" ]; then
+    CMD_ARGS="$CMD_ARGS $ONLY_FLAG"
+fi
+if [ -n "$SKIP_FLAG" ]; then
+    CMD_ARGS="$CMD_ARGS $SKIP_FLAG"
+fi
 if [ -n "$FILE_ARGS" ]; then
     CMD_ARGS="$CMD_ARGS --file $FILE_ARGS"
 fi
 
-docker run --rm $GPU_FLAG \
+docker run --rm --init $GPU_FLAG \
     --user 1000:1000 \
     -v "${INPUT_DIR}:/data/input:ro" \
     -v "${OUTPUT_DIR}:/data/output" \