Combien d’écoles françaises en zone inondable ? (ou comment manipuler des fichiers .tif avec Python et Xarray)

Début février, le sud-ouest de la France a connu des crues qui ont nécessité la fermeture de plusieurs écoles. A cette occasion la climatologue Valérie Masson-Delmotte s’est interrogé sur Twitter :

Dans ce tutoriel, nous allons voir comment on pourrait répondre à cette question en utilisant Python et Xarray.

Sources de données

Pour déterminer combien d’établissements scolaires français sont situés en zone inondable, il nous faut :

1. La localisation des établissements,
2. Les limites des zones inondables,

Il suffira ensuite de prendre les établissements un par un et de vérifier s’ils sont ou non dans une zone menacée d’inondation.

Le premier jeu de données est facile à trouver : le Ministère de l’éducation publie une liste de toutes écoles, collèges et lycées français avec leurs géolocalisation. Cette liste est accessible ici au format .csv.

C’est un peu plus compliqué pour le deuxième jeu de données. On pourrait essayer de rechercher les “zones bleues” des plans de prévention du risque inondation (PPRI) mais ces plans sont réalisés à l’échelle de la commune et ne sont pas consolidés nationalement, d’ailleurs de nombreuses communes n’en ont toujours pas. Le même problème se pose avec les atlas des zones inondables, établis par région et bassin versant… bref il n’existe pas vraiment de données officielles françaises se prêtant à une analyse rapide.

Comme pour beaucoup de projets de data science, l’accès aux données et leur compréhension est l’étape la plus délicate.

Nous allons plutôt utiliser une cartographie mondiale des zones inondables réalisée par la Commission Européenne. Il existe plusieurs versions en fonction du temps de retour et de la zone couverte : de 10 ans à 500 ans. Les cartes mondiales ont une résolution d’environ 1km contre 100 mètres pour les versions européennes.

Pour simplifier et alléger les calculs, dans ce tutoriel, nous n’utiliserons que la carte mondiale de crue centennalle.

Manipuler un GeoTIFF avec Xarray

Ces données sont disponible sous format GeoTIFF. Un fichier TIFF n’est en fait rien d’autre qu’une image en 2D avec des métadonnées, le GeoTIFF ajoutant simplement des métadonnées de géoréférencement.

Nous avons déjà parlé dans des tutos précédents de Xarray, la librarie de référence en Python pour traiter des données climatiques au format NetCDF. Ce serait bien pratique si on pouvait l’utiliser aussi pour l’exploitation de notre .tif, non ?

Et bien il se trouve que c’est possible. Pour cela, il faudra commencer par installer Rasterio et si ce n’est pas déjà fait Xarray, par exemple :

conda install -c conda-forge Xarray
conda install -c conda-forge Rasterio

Une fois que c’est fait vous pourrez ouvrir le fichier .tif sous forme de datarray très simplement :

Nous sommes maintenant en terrain connu avec un datarray qui contient pour chaque point la hauteur atteinte lors de l’inondation centenale. Pour les zones non indondable cette valeur n’est pas définie et le dataarray contient -3.402823e+38 qui remplace ici un NaN.

Commençons par réduire le fichier à la zone qui nous intéresse, c’est à dire la France métropolitaine.

Les coordonnées du fichier sont appelées x et y. x correspond à la longitude, on va donc sélectionner les x compris entre -5 (un peu à l’ouest de Brest) et 10 (un peu à l’est de Bastia). Même logique avec les y/latitudes que l’on sélectionne entre 41 et 52 :

A ce stade, on peut représenter nos données pour s’assurer qu’on ne s’est pas trompé. Pour obtenir une carte sommaire, il suffit de faire :

da.plot()

Déterminer si un établissement est en zone inondable

Ouvrons maintenant le fichier contenant la liste des établissements scolaires :

Comme c’est souvent le cas pour les csv qui sont passés par une version française d’Excel, le séparateur n’est pas une virgule comme le voudrait le nom (Comma-Separated Values) mais un point virgule. Il suffit de le préciser avec l’argument sep = ‘;’ :

import pandas as pd
data = pd.read_csv("etablissements.csv", sep = ';')

Commençons par récupérer le nom et la commune du premier établissement :

etablissement = data.iloc[0]
print(etablissement.appellation_officielle)
print(etablissement.libelle_commune)

Il s’agit du “Lycée général et technologique privé Jules Froment” à Aubenas. Est-il en zone inondable ? Pour le savoir, nous allons simplement récupérer ses coordonnées et chercher la hauteur d’eau correspondant à l’inondation centenale pour ce point dans le datarray que nous avons déjà ouvert :

lat = etablissement.latitude
lon = etablissement.longitude
float(da.sel(x = lat, y = lon, method = 'nearest'))

On obtient -3.402823e+38, ce qui signifie que ce point n’est pas inondable. Au contraire une valeur supérieure à 0 indiquerait que l’établissement se trouve en zone inondable.

Synthèse

Il nous suffit de reproduire ce calcul pour l’ensemble des établissements. Mais auparavant nous allons supprimer de la liste ceux qui ne sont pas en métropole et ceux pour lesquels les coordonnées sont manquantes :

On voit que sur les 65694 établissements contenus dans le fichier, 63217 sont en métropole et parmi ceux-ci 63162 sont géolocalisés.

On va parcourir cette liste et chercher pour chaque établissement la hauteur d’eau correspondant à l’inondation centennale puis ajouter le résultat dans une nouvelle colonne :

Il ne reste plus qu’à compter les lignes pour lesquelles la hauteur d’eau lors de l’inondation centennale est supérieure à 0 :

Et voilà le résultat :

3926 établissements scolaires français, soit 6.2% du total, se trouvent dans une zone inondable avec un temps de retour de 100 ans.

A propos : Callendar est une start-up spécialisée dans le développement de solutions innovantes pour l’évaluation des risques climatiques. Conscients du défi que représente l’adaptation au changement climatique, nous nous efforçons de partager notre expertise au travers d’outils gratuits ou de tutoriels comme celui-ci.