Rasterio: Numpy for Geospatial Data

Salut, c'est Alex de DEV STORIES DOT EU. Rasterio : Numpy pour les données géospatiales, épisode 1 sur 5. Tu veux extraire des valeurs de pixels d'une image satellite, mais les librairies d'images standards galèrent avec les formats géospatiaux, et les bindings officiels te donnent l'impression d'écrire du code C en Python. C'est exactement pour ça que Rasterio existe. Rasterio est une librairie qui lit et écrit des données géospatiales sous forme de grille, comme les fichiers GeoTIFF. Il est très probable que tu aies déjà entendu parler de GDAL, la librairie géospatiale de référence dans l'industrie. Une idée reçue courante est que Rasterio remplace GDAL. Ce n'est pas le cas. En fait, Rasterio utilise GDAL sous le capot. La différence, c'est l'interface. Si tu as déjà utilisé les bindings Python natifs de GDAL, tu sais qu'ils exposent une API C un peu lourde, pleine de gestion manuelle des ressources et avec une syntaxe qui semble totalement étrangère à Python. Rasterio laisse tout ça derrière lui. Il fournit une interface moderne, pensée d'abord pour Python. Il traite les datasets raster comme des fichiers normaux et gère leurs pixels comme des arrays Numpy standards. Supposons que tu aies un GeoTIFF d'images Landsat sur ton disque. Tu ouvres ce fichier en utilisant la fonction rasterio point open. Comme Rasterio adopte les patterns Python standards, tu le fais avec un statement with. Ça crée un context manager. Tout comme quand tu ouvres un simple fichier texte, le bloc with garantit que le dataset est proprement fermé et que la mémoire est libérée dès que ton code a fini de tourner. Tu n'as jamais besoin de détruire manuellement des objets ou de t'inquiéter des memory leaks causés par des fichiers non fermés. À l'intérieur de ce bloc with, la fonction open te donne un objet dataset reader. Avant même de lire les vrais pixels, tu peux inspecter le fichier. Tu peux demander au dataset son count, ce qui t'indique le nombre total de bandes, ou de couches, dans l'image. Pour une image Landsat, ce count pourrait être quelque chose comme sept ou onze, représentant différentes longueurs d'onde de lumière. Tu peux lire les propriétés width et height pour obtenir les dimensions spatiales en pixels. Tu peux aussi vérifier la propriété dtypes. Ça te donne le data type des pixels pour chaque bande, comme des entiers non signés 16 bits ou des floats 32 bits. Toutes ces metadata sont disponibles instantanément, sans avoir à charger les lourdes données de l'image dans la mémoire de ton système. Voici le point clé. Quand tu es prêt à regarder les vrais pixels, tu appelles la méthode read sur l'objet dataset. Tu peux passer un index spécifique pour charger juste une couche. Fais bien attention à ça : les bandes Rasterio sont indexées à partir de un, ce qui veut dire que tu demandes la bande un, et non la bande zéro. Si tu appelles read et que tu passes le chiffre un, la méthode te renvoie ces pixels sous la forme d'un array Numpy standard en deux dimensions. Si tu appelles read sans aucun argument, elle lit tout, et te renvoie un array Numpy en trois dimensions contenant toutes les bandes. Il n'y a pas d'objet pixel propriétaire spécial. Tu obtiens juste un array de nombres. Une fois que tu as cet array Numpy, tu es de retour en terrain connu. Tu peux le slicer, lancer des opérations statistiques, ou le passer directement à des modèles de machine learning en utilisant exactement le même code que tu utilises pour n'importe quelle autre matrice. Rasterio agit comme un pont. Il gère les formats de fichiers compliqués et les metadata sur le disque, et te donne un objet mathématique propre en mémoire. Si tu aimes l'émission et que tu veux soutenir ce qu'on fait, tu peux chercher DevStoriesEU sur Patreon — c'est une aide énorme. La vraie puissance de Rasterio n'est pas de réinventer le traitement géospatial, mais de le rendre banal ; en transformant des images satellites complexes en arrays Numpy standards, il te permet d'arrêter de te battre avec les formats GIS et de commencer à faire de la vraie data science. Merci d'avoir écouté, et happy coding à tous !

Bonjour, c'est Alex de DEV STORIES DOT EU. Rasterio : Numpy pour les données géospatiales, épisode 2 sur 5. Un array Numpy de base, c'est juste une grille plate de nombres. Si tu le visualises, tu obtiens simplement un rectangle coloré sur ton écran. Qu'est-ce qui transforme vraiment ces nombres abstraits en une carte précise du plateau du Colorado ? La réponse, c'est l'Affine Transform. On imagine souvent que les données géospatiales, c'est juste une image dont les coins sont vaguement épinglés à des coordonnées géographiques. Ce n'est pas le cas. Le transform est un lien mathématique strict qui calcule la position réelle précise de chaque pixel dans ton array. Avant de pouvoir mapper des pixels sur la Terre, tu dois savoir comment tu mesures la Terre. C'est le Coordinate Reference System. Tu lis ça depuis le dataset en accédant à dataset point crs. Cette propriété indique à Rasterio si tes coordonnées spatiales sont mesurées en degrés de longitude et de latitude, ou en mètres par rapport à un équateur et un méridien d'origine spécifiques. Une fois le système de référence défini, tu as besoin de la logique de mapping. Tu trouves ça en vérifiant dataset point transform. Ça te donne une matrice de transformation affine. Ça ressemble à de l'algèbre linéaire un peu lourde, mais c'est en fait un ensemble de six nombres qui décrivent l'empreinte physique de ta grille. La matrice définit les coordonnées spatiales x et y exactes du coin supérieur gauche du pixel en haut à gauche. Elle définit la largeur physique d'un pixel, par exemple, exactement trente mètres. Elle définit la hauteur physique d'un pixel, qui est généralement négative parce que les lignes de l'image comptent vers le bas alors que les coordonnées de la carte comptent vers le haut. La matrice contient aussi des valeurs de rotation au cas où le satellite était incliné quand l'image a été capturée. Voici le point clé. Tu n'as jamais besoin de multiplier ces valeurs de matrice à la main. Rasterio te donne des méthodes directes pour basculer entre l'espace des pixels et l'espace des coordonnées. Supposons que tu analyses une image et que tu trouves une anomalie à un endroit précis, disons la ligne cinq cents et la colonne mille. Pour savoir où ça se trouve sur la planète réelle, tu passes ces deux entiers à la méthode dataset point xy. Rasterio fait passer la ligne et la colonne à travers la matrice affine et retourne les coordonnées spatiales x et y exactes du point central de ce pixel spécifique. Tu peux aussi faire tourner cette machine à l'envers. Disons que tu as les coordonnées spatiales exactes d'un emplacement physique. C'est peut-être un canyon situé exactement à cent kilomètres à l'est et cinquante kilomètres au sud de l'origine de ton image. Tu dois extraire les données de pixels sous-jacentes pour ce canyon. Tu prends tes valeurs spatiales x et y et tu les passes à la méthode dataset point index. Rasterio inverse la matrice affine, traite tes coordonnées, et te renvoie les entiers exacts de la ligne et de la colonne. Tu utilises ensuite ces entiers pour slicer dans ton array Numpy et lire la mesure exacte enregistrée par le satellite. L'Affine Transform isole ton array de données brutes de ta géométrie spatiale, garantissant que peu importe comment tu recadres ou lis tes données, tu ne perds jamais ta position précise sur le globe. C'est tout pour aujourd'hui. Merci d'avoir écouté — va construire un truc cool.

Salut, c'est Alex de DEV STORIES DOT EU. Rasterio : Numpy pour les données géospatiales, épisode 3 sur 5. Si tu essaies de lire un GeoTIFF de dix gigaoctets directement dans la mémoire de ton ordinateur portable, ton script va presque certainement planter. Tu n'as peut-être besoin que d'une petite portion de pixels pour tester un algorithme, mais le système essaie d'engloutir tout le continent d'un coup. La solution à cette limite de mémoire, c'est le Windowed Processing. Le Windowed Processing te permet de lire ou d'écrire des sous-ensembles rectangulaires spécifiques d'un fichier raster, pendant que le reste du fichier reste en sécurité sur ton disque dur. C'est entièrement basé sur des index. On ne va pas manipuler de coordonnées géographiques aujourd'hui. On va gérer strictement des lignes et des colonnes de pixels. Pour lire un sous-ensemble, tu dois avoir un moyen d'indiquer à Rasterio quels pixels récupérer. Tu fais ça en important l'objet Window depuis le module rasterio windows. La façon standard de définir un Window, c'est de fournir quatre nombres : le column offset, le row offset, la width et la height. Fais bien attention à ça. L'ordre des arguments est important. Il faut d'abord le column offset, puis le row offset. Le column offset, c'est ta position horizontale de départ, en mesurant depuis le bord gauche. Le row offset, c'est ta position verticale de départ, en mesurant vers le bas depuis le bord supérieur. Après les offsets, tu spécifies la width et la height de ta sélection en pixels. Imagine un énorme fichier de test de la taille d'un gigaoctet. Tu veux extraire une petite tuile qui fait 256 pixels de large et 512 pixels de haut. Tu veux que cette tuile commence à mille pixels du bord gauche et à deux mille pixels du bord supérieur. Tu instancies ton Window en passant mille pour le column offset, deux mille pour le row offset, 256 pour la width, et 512 pour la height. Ensuite, tu ouvres ton dataset en utilisant la fonction open standard. Quand tu appelles la méthode read, tu ne laisses pas les arguments vides. Des arguments vides disent à Rasterio de charger tout le fichier en mémoire. À la place, tu assignes ton nouvel objet Window au keyword argument window dans la méthode read. Rasterio traduit ce window en byte offsets sur le disque, récupère uniquement ce chunk spécifique de 256 par 512, et le retourne sous forme de numpy array. Ton utilisation de la mémoire bouge à peine. Si tu bosses avec numpy toute la journée, passer des offsets et des dimensions peut te sembler un peu contre-nature. Tu as probablement l'habitude de définir des slices avec des index start et stop pour tes lignes et tes colonnes. Rasterio supporte ce workflow avec une méthode alternative appelée Window dot from slices. Cette méthode prend deux inputs. Le premier input est un tuple qui définit le row start index et le row stop index. Le deuxième input est un tuple qui définit le column start index et le column stop index. Remarque bien que l'ordre ici, c'est d'abord les lignes, puis les colonnes. Ça reflète exactement la façon dont tu slices un numpy array en deux dimensions. Sous le capot, cette méthode calcule les offsets et les dimensions pour toi, et retourne un objet Window standard que tu peux passer à la méthode read, exactement comme avant. En imposant des windowed reads dans tes data pipelines, tu découples de façon permanente la taille physique de tes fichiers d'input des limites matérielles de ta machine. C'est tout pour cet épisode. Merci d'avoir écouté, et continue de développer !

Bonjour, ici Alex de DEV STORIES DOT EU. Rasterio : Numpy pour les données géospatiales, épisode 4 sur 5. Tu viens de télécharger une immense image satellite couvrant un État entier, mais seuls les pixels situés à l'intérieur des limites d'une ferme spécifique t'intéressent. Charger toute la grille gaspille de la mémoire, et la lecture par décalage de pixels est inutile quand ta limite est un polygone de forme irrégulière. La solution, c'est le Masking avec des Vectors. Dans un épisode précédent, nous avons regardé le windowed reading. Cette approche nécessitait des indices de pixels exacts, comme demander la ligne cinquante et la colonne cent. Le Masking est totalement différent. Le Masking repose sur des coordonnées géographiques réelles. Au lieu de compter les pixels, tu fournis une géométrie vectorielle, comme un polygone défini par des longitudes et latitudes réelles. Rasterio détermine quels pixels tombent à l'intérieur de cette forme et gère la traduction complexe entre l'espace de coordonnées vectorielles continu et la grille raster discrète. L'outil principal pour ce workflow est la fonction mask, située dans le module rasterio dot mask. Tu lui passes deux arguments principaux. Premièrement, un dataset raster ouvert. Deuxièmement, un iterable de formes géométriques. Ces formes doivent être formatées comme des dictionnaires respectant la spécification GeoJSON. En général, tu lis un shapefile contenant le polygone de ta ferme, tu extrais le dictionnaire de géométrie pour cette feature spécifique, tu le places dans une list Python standard, et tu le passes à la fonction mask. Voici l'élément clé. La fonction mask possède deux comportements distincts contrôlés par un seul paramètre boolean appelé crop. Tu dois comprendre la différence entre simplement faire un masking sur un dataset et le cropper réellement. Si tu lances la fonction avec crop défini sur false, Rasterio regarde le polygone de ta ferme et définit chaque pixel en dehors de cette limite à une valeur nodata. Il lit la valeur nodata par défaut directement depuis les metadata raster d'origine. L'output array que tu reçois a exactement les mêmes dimensions que ton immense image d'origine couvrant tout l'État. Les pixels à l'intérieur de la ferme conservent leurs valeurs réelles, mais tout ce qui est à l'extérieur est effacé. L'étendue spatiale globale du fichier ne change pas. Si tu définis crop sur true, le comportement change. Rasterio définit toujours les pixels en dehors de ton polygone sur nodata. Cependant, il ignore ensuite les vastes zones d'espace vide. Il calcule la bounding box rectangulaire du polygone irrégulier de ta ferme et réduit l'array renvoyé pour qu'il s'insère exactement dans cette box spécifique. Ton array géant couvrant tout l'État est physiquement réduit à une petite grille gérable contenant juste la ferme et un peu de padding nodata autour de ses bords irréguliers. Parce que le cropping modifie physiquement les dimensions de la grille, il modifie aussi le point de départ géographique. Le coin supérieur gauche n'est plus le coin supérieur gauche de l'État. C'est maintenant le coin supérieur gauche de la bounding box de la ferme. Par conséquent, la fonction mask renvoie un tuple contenant deux éléments. Le premier élément est ton array Numpy fraîchement croppé. Le second élément est un tout nouvel objet affine transform. Ce transform mis à jour contient les nouvelles coordonnées d'origine. Si tu veux sauvegarder cet array fraîchement croppé dans un nouveau fichier, tu dois passer ce transform mis à jour dans ton write profile. Si tu réutilises par erreur le transform d'origine de tout l'État, l'image de ta petite ferme sera étirée à nouveau sur l'ensemble de l'État. La vraie valeur de la fonction mask est qu'elle gère tout ce recalcul spatial pour toi automatiquement. J'espère que ça t'a été utile. Merci d'avoir écouté, et profite bien du reste de ta journée.

Salut, c'est Alex de DEV STORIES DOT EU. Rasterio : Numpy pour les données géospatiales, épisode 5 sur 5. Parfois, tu as de magnifiques images satellites haute résolution, mais ton modèle climatique n'accepte que des inputs grossiers et basse résolution. Tu ne peux pas juste slicer l'array, parce que les pixels qui en résultent doivent toujours représenter exactement la même zone physique, juste sur une grille moins dense. Pour faire ça correctement, tu as besoin du resampling et du writing. Il ne s'agit pas de reprojection. Le coordinate reference system reste exactement le même. On change uniquement la résolution. Disons que tu dois downsampler un raster de moitié. Dans Rasterio, tu fais le resampling au moment précis où tu lis les données en mémoire. Quand tu appelles la méthode read sur ton dataset, tu passes un paramètre appelé out shape. Au lieu de lire tout l'array, Rasterio le lit et le scale immédiatement pour qu'il corresponde à la shape que tu as demandée. Pour un raster à une seule bande, ton out shape serait un tuple spécifiant une bande, puis la moitié de la hauteur d'origine, et la moitié de la largeur d'origine. Tu dois aussi indiquer à Rasterio comment calculer les valeurs de ces nouveaux pixels plus grands. Tu fais ça en passant un paramètre resampling. Rasterio fournit un enum Resampling pour ça. Si tu traites des données continues comme la température, tu peux passer Resampling dot bilinear, qui calcule une moyenne pondérée des pixels environnants. Si tu traites des données catégorielles comme l'occupation du sol, tu vas utiliser Resampling dot nearest pour conserver les valeurs de classe d'origine exactes. Voici le point clé. Tu as maintenant un nouvel array Numpy plus petit, mais il représente toujours exactement la même zone géographique sur Terre. Comme l'array contient moins de pixels, chaque pixel individuel couvre un espace physique plus grand. Ça veut dire que ton Affine transform d'origine est maintenant complètement fausse. Si tu écris le nouvel array sur le disque en utilisant l'ancienne transform, ton image va rétrécir sur la carte pour ne couvrir plus qu'un quart de sa surface d'origine. Tu dois scaler la transform. Tu fais ça en prenant la transform du dataset d'origine et en la multipliant par une matrice de scaling Affine. Tu trouves les ratios de scaling en divisant la largeur d'origine par ta nouvelle largeur, et la hauteur d'origine par ta nouvelle hauteur. Puisqu'on a divisé les dimensions par deux, ces deux ratios sont exactement de deux. Quand tu multiplies la transform d'origine par ces facteurs de scaling, tu obtiens une nouvelle transform où les dimensions en pixels sont doublées, mais le point d'origine en haut à gauche reste parfaitement ancré. Tu as maintenant ton array resamplé et ta transform scalée. L'étape finale, c'est d'écrire un GeoTIFF valide sur le disque. Pour faire ça, tu as besoin de metadata. L'approche la plus sûre est de copier la propriété profile de ton dataset source. Ça te donne un dictionnaire contenant le file driver, le coordinate reference system, le data type et la valeur nodata. Tu mets ensuite à jour ce dictionnaire profile avec ta nouvelle hauteur, ta nouvelle largeur et ta transform fraîchement calculée. Une fois le profile mis à jour, tu appelles rasterio dot open. Tu fournis le nouveau nom de fichier, tu mets le mode sur write, et tu unpack ton dictionnaire profile comme keyword arguments. Ça crée un dataset vide sur le disque, formaté selon tes spécifications exactes. Enfin, tu appelles la méthode write sur ce nouvel objet dataset et tu lui passes ton array Numpy downsamplé. L'array est écrit dans le fichier, et ton nouveau GeoTIFF basse résolution est prêt pour le modèle climatique. N'oublie pas, modifier des données raster est toujours un contrat en deux parties : si tu changes la shape de l'array Numpy, tu dois explicitement modifier l'Affine transform pour qu'elle corresponde, sinon ta spatial footprint est cassée. Ça nous amène à la fin de notre série. Je t'encourage à lire la documentation officielle de Rasterio pour explorer d'autres algorithmes de resampling et essayer de construire ces pipelines toi-même. Si tu veux suggérer des sujets pour une future série, visite devstories dot eu. Je voudrais prendre un moment pour te remercier de ton écoute, ça nous aide beaucoup. Passe une bonne journée !