Praca z danymi rastrowymi

Praca z danymi geoprzestrzennymi w Pythonie

Joris Van den Bossche

Open source software developer and teacher, GeoPandas maintainer

Ilustracja siatki wartości

Źródło obrazu: dokumentacja QGIS

Dane rastrowe

Mapa szansy na opady

Dane rastrowe z wieloma kanałami

Pakiet rasterio

import rasterio

Interfejs „pythoniczny" dla GDAL
Odczyt i zapis plików rastrowych
Narzędzia przetwarzania (maskowanie, reprojekcja, resampling, ..)

https://rasterio.readthedocs.io/en/latest/

Otwieranie pliku rastrowego

import rasterio

src = rasterio.open("DEM_world.tif")

Metadane:

src.count

src.width, src.height

(4320, 2160)

Dane rastrowe = tablica numpy

array = src.read()

Standardowa tablica numpy:

array

array([[[-4290, -4290, -4290, ..., -4290, -4290, -4290],
        [-4278, -4278, -4278, ..., -4278, -4278, -4278],
        [-4269, -4269, -4269, ..., -4269, -4269, -4269],
        ...,
        [ 2804,  2804,  2804, ...,  2804,  2804,  2804],
        [ 2804,  2804,  2804, ...,  2804,  2804,  2804],
        [ 2804,  2804,  2804, ...,  2804,  2804,  2804]]], dtype=int16)

Wizualizacja danych rastrowych

Metoda rasterio.plot.show():

import rasterio.plot

rasterio.plot.show(src, cmap='terrain')

Ekstrakcja informacji na podstawie danych wektorowych

rasterstats: Statystyki podsumowujące rastrowych danych geoprzestrzennych na podstawie geometrii wektorowych (https://github.com/perrygeo/python-rasterstats)

Ekstrakcja wartości rastrowych z rasterstats

Dla geometrii punktowych:

rasterstats.point_query(geometries, "path/to/raster", 
                        interpolation='nearest'|'bilinear')

Dla geometrii poligonowych:

rasterstats.zonal_stats(geometries, "path/to/raster",
                        stats=['min', 'mean', 'max'])

Ekstrakcja wartości rastrowych z rasterstats

result = rasterstats.zonal_stats(countries.geometry, "DEM_gworld.tif", 
                                 stats=['mean'])

countries['mean_elevation'] = pd.DataFrame(result)

countries.sort_values('mean_elevation', ascending=False).head()

            name    continent                     geometry  mean_elevation
157   Tajikistan         Asia  POLYGON ((74.98 37.41, ...      3103.231105
85    Kyrgyzstan         Asia  POLYGON ((80.25 42.34, ...      2867.717142
24        Bhutan         Asia  POLYGON ((91.69 27.77, ...      2573.559846
119        Nepal         Asia  POLYGON ((81.11 30.18, ...      2408.907816
6     Antarctica   Antarctica  (POLYGON ((-59.57 -80.04...     2374.075028
..           ...          ...                          ...             ...

Czas na ćwiczenia!

Praca z danymi geoprzestrzennymi w Pythonie