👽 GeoParquet

GeoParquet¶

GeoParquet to format zoptymalizowany do przechowywania i przetwarzania danych geoprzestrzennych.

• Umożliwia szybkie odczytywanie i zapisywanie dużych zbiorów danych dzięki kolumnowej strukturze.
• Obsługuje równoległe przetwarzanie, co czyni go idealnym do analizy dużych zbiorów danych.
• Umożliwia efektywne filtrowanie i kompresję, zmniejszając ilość przesyłanych i zapisywanych danych.
• W przeciwieństwie do formatów takich jak Shapefile czy GeoJSON, pozwala na przechowywanie atrybutów o różnych typach danych bez ograniczeń.

In [ ]:

import geopandas as gpd
import folium
from shapely.geometry import box
from fsspec.implementations.http import HTTPFileSystem
import pyarrow.parquet as pq
from geopandas.io.arrow import _arrow_to_geopandas
from shapely.geometry import box
from folium.plugins import HeatMap

import geopandas as gpd import folium from shapely.geometry import box from fsspec.implementations.http import HTTPFileSystem import pyarrow.parquet as pq from geopandas.io.arrow import _arrow_to_geopandas from shapely.geometry import box from folium.plugins import HeatMap

In [ ]:

file_url = "https://data.source.coop/cholmes/eurocrops/unprojected/geoparquet/FR_2018_EC21.parquet"
filesystem = HTTPFileSystem()

file_url = "https://data.source.coop/cholmes/eurocrops/unprojected/geoparquet/FR_2018_EC21.parquet" filesystem = HTTPFileSystem()

Wczytanie pliku GeoParquet do GeoDataFrame i podgląd danych¶

GeoPandas to rozszerzenie dla Pandas, które umożliwia pracę z danymi przestrzennymi w Pythonie. Ułatwia analizę geometrii (punkty, linie, poligony), przekształcanie układów współrzędnych, operacje przestrzenne (przecięcia, łączenia, buforowanie) oraz wizualizację na mapach.

GeoPandas wykorzystuje Shapely do obsługi geometrii, Fiona do wczytywania danych wektorowych, a Matplotlib do prostych wizualizacji.

To tylko wczytuje metadane

In [ ]:

parquet_file = pq.ParquetFile(file_url, filesystem=filesystem)
print(f"Columns: {parquet_file.schema_arrow.names}")
print(f"Number of rows: {parquet_file.metadata.num_rows}")
print(f"Number of row groups: {parquet_file.num_row_groups}")

parquet_file = pq.ParquetFile(file_url, filesystem=filesystem) print(f"Columns: {parquet_file.schema_arrow.names}") print(f"Number of rows: {parquet_file.metadata.num_rows}") print(f"Number of row groups: {parquet_file.num_row_groups}")

Tutaj wczytujemy faktycznie dane, ale tylko dla jednej z grup wierszy

In [ ]:

pyarrow_table = parquet_file.read_row_group(0, columns=["ID_PARCEL", "SURF_PARC", "geometry"])
geopandas_gdf = _arrow_to_geopandas(pyarrow_table)

pyarrow_table = parquet_file.read_row_group(0, columns=["ID_PARCEL", "SURF_PARC", "geometry"]) geopandas_gdf = _arrow_to_geopandas(pyarrow_table)

In [ ]:

geopandas_gdf.head(10)

geopandas_gdf.head(10)

In [ ]:

# Przekształcenie do EPSG:4326 jeśli dane są w innym układzie
if geopandas_gdf.crs.to_epsg() != 4326:
    geopandas_gdf = geopandas_gdf.to_crs(epsg=4326)

# Przekształcenie do EPSG:4326 jeśli dane są w innym układzie if geopandas_gdf.crs.to_epsg() != 4326: geopandas_gdf = geopandas_gdf.to_crs(epsg=4326)

Naprawianie geometrii¶

In [ ]:

# Usunięcie pustych i nieprawidłowych geometrii
geopandas_gdf = geopandas_gdf[~geopandas_gdf.geometry.is_empty & geopandas_gdf.geometry.notnull()]

# Naprawa geometrii, aby uniknąć błędów topologicznych
geopandas_gdf["geometry"] = geopandas_gdf["geometry"].buffer(0)

# Usunięcie pustych i nieprawidłowych geometrii geopandas_gdf = geopandas_gdf[~geopandas_gdf.geometry.is_empty & geopandas_gdf.geometry.notnull()] # Naprawa geometrii, aby uniknąć błędów topologicznych geopandas_gdf["geometry"] = geopandas_gdf["geometry"].buffer(0)

Filtrowanie¶

In [ ]:

# Definicja granic Francji w EPSG:4326
france_bbox = box(-5.0, 41.0, 9.7, 51.1)

# Usunięcie geometrii znajdujących się poza granicami Francji
geopandas_gdf = geopandas_gdf[geopandas_gdf.geometry.within(france_bbox)]

# Definicja granic Francji w EPSG:4326 france_bbox = box(-5.0, 41.0, 9.7, 51.1) # Usunięcie geometrii znajdujących się poza granicami Francji geopandas_gdf = geopandas_gdf[geopandas_gdf.geometry.within(france_bbox)]

In [ ]:

# Konwersja kolumny SURF_PARC do float, jeśli zawiera decimal.Decimal
geopandas_gdf["SURF_PARC"] = geopandas_gdf["SURF_PARC"].astype(float)

# Obliczenie progu dla 10% największych wartości powierzchni
threshold = geopandas_gdf["SURF_PARC"].quantile(0.9)

# Konwersja kolumny SURF_PARC do float, jeśli zawiera decimal.Decimal geopandas_gdf["SURF_PARC"] = geopandas_gdf["SURF_PARC"].astype(float) # Obliczenie progu dla 10% największych wartości powierzchni threshold = geopandas_gdf["SURF_PARC"].quantile(0.9)

In [ ]:

# Filtrowanie 10% największych parceli
gdf_top_10 = geopandas_gdf[geopandas_gdf["SURF_PARC"] >= threshold]

# Filtrowanie 10% największych parceli gdf_top_10 = geopandas_gdf[geopandas_gdf["SURF_PARC"] >= threshold]

Mapa interaktywna¶

In [ ]:

# Tworzenie mapy interaktywnej
m = folium.Map(
    location=[gdf_top_10.geometry.centroid.y.mean(), gdf_top_10.geometry.centroid.x.mean()],
    zoom_start=10,
)


# Dodanie warstwy polygonów
def add_gdf_to_map(gdf, fmap):
    for _, row in gdf.iterrows():
        folium.GeoJson(row.geometry).add_to(fmap)


add_gdf_to_map(gdf_top_10, m)

# Wyświetlenie mapy
m

# Tworzenie mapy interaktywnej m = folium.Map( location=[gdf_top_10.geometry.centroid.y.mean(), gdf_top_10.geometry.centroid.x.mean()], zoom_start=10, ) # Dodanie warstwy polygonów def add_gdf_to_map(gdf, fmap): for _, row in gdf.iterrows(): folium.GeoJson(row.geometry).add_to(fmap) add_gdf_to_map(gdf_top_10, m) # Wyświetlenie mapy m

Heatmap¶

In [ ]:

# Tworzenie mapy heatmap
heatmap_data = [[point.y, point.x] for point in geopandas_gdf.geometry.centroid]

m_heatmap = folium.Map(
    location=[gdf_top_10.geometry.centroid.y.mean(), gdf_top_10.geometry.centroid.x.mean()],
    zoom_start=10,
)
HeatMap(heatmap_data).add_to(m_heatmap)

# Wyświetlenie mapy heatmap
m_heatmap

# Tworzenie mapy heatmap heatmap_data = [[point.y, point.x] for point in geopandas_gdf.geometry.centroid] m_heatmap = folium.Map( location=[gdf_top_10.geometry.centroid.y.mean(), gdf_top_10.geometry.centroid.x.mean()], zoom_start=10, ) HeatMap(heatmap_data).add_to(m_heatmap) # Wyświetlenie mapy heatmap m_heatmap