Tutoriel 4 : CostGrow et comparaison de modèles#
Ce tutoriel prolonge le tutoriel 2 en exécutant les mêmes données d’exemple avec le modèle intégré CostGrow_Terrain et le modèle téléchargé ResUNet_16x_DEM.
À la fin, vous aurez :
un raster de sortie CostGrow
un raster de sortie ResUNet
un raster de comparaison contenant
CostGrow - ResUNetun graphique côte à côte des deux sorties de modèle
Avant de commencer#
Ce tutoriel nécessite l’installation étendue, car CostGrow_Terrain dépend de PCRaster.
CLI (
pipx) : utilisez l’environnement étendu de Ligne de commande (CLI) - étendue.carnet local (Jupyter) : utilisez l’environnement étendu de Carnet local (Jupyter) - étendu.
carnet hébergé (Colab) : non recommandé pour ce tutoriel, car la pile étendue est plus lourde et moins fiable.
Installer des paquets supplémentaires#
# %pip install -q matplotlib rasterio
Vérifier l’installation#
Vérifions maintenant les versions :
import matplotlib
import rasterio
print(f"matplotlib=={matplotlib.__version__}")
print(f"rasterio=={rasterio.__version__}")
matplotlib==3.10.8
rasterio==1.5.0
!floodsr --version
Importations#
from pathlib import Path
from urllib.request import urlretrieve
import matplotlib.pyplot as plt
import numpy as np
import rasterio
Télécharger les mêmes données de test que dans le tutoriel 2#
Nous réutilisons les petits rasters d’exemple des tutoriels précédents afin que la comparaison reste rapide à exécuter et facile à examiner.
urlretrieve(
"https://github.com/cefect/floodsr/releases/download/v0.0.3/hires002_dem.tif",
"hires002_dem.tif",
)
urlretrieve(
"https://github.com/cefect/floodsr/releases/download/v0.0.3/lowres032.tif",
"lowres032.tif",
)
('lowres032.tif', <http.client.HTTPMessage at 0x7374613018e0>)
Résoudre les chemins d’entrée#
lowres_fp = Path("lowres032.tif").resolve()
dem_fp = Path("hires002_dem.tif").resolve()
assert lowres_fp.is_file(), f"missing low-res raster\n {lowres_fp}"
assert dem_fp.is_file(), f"missing DEM raster\n {dem_fp}"
costgrow_fp = Path("lowres032_costgrow_sr.tif").resolve()
resunet_fp = Path("lowres032_resunet_sr.tif").resolve()
comparison_fp = Path("lowres032_costgrow_minus_resunet.tif").resolve()
Récupérer les poids ResUNet#
CostGrow_Terrain est intégré, il n’a donc pas besoin de téléchargement de poids. ResUNet_16x_DEM en a besoin, donc nous les récupérons une fois ici avant de lancer la comparaison.
!floodsr models fetch --no-progress ResUNet_16x_DEM
version=ResUNet_16x_DEM stored=/tmp/.cache/floodsr/ResUNet_16x_DEM/model_infer.onnx retrieved_from=cache
Tracer les entrées du tutoriel#
def read_raster_band_1(fp):
"""Read band 1 as float and promote masked nodata pixels to NaN."""
assert Path(fp).is_file(), f"missing raster\n {fp}"
with rasterio.open(fp) as src:
return src.read(1, masked=True).astype(float).filled(np.nan)
def write_float_raster_like(reference_fp, out_fp, arr):
"""Write one float32 raster using the profile from a reference raster."""
with rasterio.open(reference_fp) as src:
profile = src.profile.copy()
profile.update(dtype="float32", count=1, nodata=np.nan)
with rasterio.open(out_fp, "w", **profile) as dst:
dst.write(arr.astype("float32"), 1)
lowres_arr = read_raster_band_1(lowres_fp)
dem_arr = read_raster_band_1(dem_fp)
fig, axes = plt.subplots(1, 2, figsize=(10, 4))
im0 = axes[0].imshow(np.where(lowres_arr > 0, lowres_arr, np.nan), cmap="Blues")
axes[0].set_title("Low-res flood depth")
fig.colorbar(im0, ax=axes[0])
axes[0].set_axis_off()
im1 = axes[1].imshow(dem_arr, cmap="terrain")
axes[1].set_title("High-res DEM")
fig.colorbar(im1, ax=axes[1])
axes[1].set_axis_off()
fig.tight_layout()
Exécuter CostGrow#
Exécutez maintenant tohr avec le modèle intégré CostGrow_Terrain. Nous définissons explicitement --window-method hard afin que le tutoriel reflète le chemin CostGrow actuel adapté aux grands rasters.
!floodsr -q tohr --in lowres032.tif --dem hires002_dem.tif --model-version CostGrow_Terrain --window-method hard --tile-overlap 0 --out lowres032_costgrow_sr.tif
Exécuter la référence ResUNet#
Exécutez ensuite le modèle appris ResUNet_16x_DEM sur les mêmes entrées afin de comparer directement les deux sorties.
!floodsr -q tohr --in lowres032.tif --dem hires002_dem.tif --model-version ResUNet_16x_DEM --out lowres032_resunet_sr.tif
Écrire un raster de comparaison#
Pour une comparaison numérique rapide, nous enregistrons un raster de CostGrow - ResUNet sur la même grille que les sorties de modèle. Les valeurs positives indiquent des prédictions plus profondes avec CostGrow, et les valeurs négatives indiquent des prédictions plus profondes avec ResUNet.
costgrow_arr = read_raster_band_1(costgrow_fp)
resunet_arr = read_raster_band_1(resunet_fp)
comparison_arr = costgrow_arr - resunet_arr
write_float_raster_like(costgrow_fp, comparison_fp, comparison_arr)
comparison_arr
Tracer les résultats côte à côte#
Les trois premiers panneaux comparent l’entrée et les deux sorties haute résolution. Le quatrième panneau montre le raster de différence signée afin que vous puissiez repérer rapidement les zones où les sorties des modèles divergent.
fig, axes = plt.subplots(1, 4, figsize=(16, 4))
plot_l = [
(axes[0], np.where(lowres_arr > 0, lowres_arr, np.nan), "Input low-res", "Blues"),
(axes[1], np.where(costgrow_arr > 0, costgrow_arr, np.nan), "CostGrow", "Blues"),
(axes[2], np.where(resunet_arr > 0, resunet_arr, np.nan), "ResUNet", "Blues"),
]
for ax, arr, title, cmap in plot_l:
wet_pct = 100.0 * float(np.count_nonzero(np.nan_to_num(arr, nan=0.0) > 0.0)) / float(arr.size)
im = ax.imshow(arr, cmap=cmap)
ax.set_title(title)
ax.axis("off")
ax.text(0.03, 0.03, f"wet={wet_pct:.1f}%", transform=ax.transAxes, color="black", fontsize=10, bbox={"facecolor": "white", "alpha": 0.8, "edgecolor": "none"})
fig.colorbar(im, ax=ax, fraction=0.046, pad=0.04)
delta_vmax = float(np.nanmax(np.abs(comparison_arr)))
if not np.isfinite(delta_vmax) or delta_vmax == 0.0:
delta_vmax = 1.0
im_delta = axes[3].imshow(comparison_arr, cmap="PiYG", vmin=-delta_vmax, vmax=delta_vmax)
axes[3].set_title("CostGrow - ResUNet")
axes[3].axis("off")
fig.colorbar(im_delta, ax=axes[3], fraction=0.046, pad=0.04)
fig.tight_layout()
