Filtrować według obwiedni w geopandzie?

11

Mam geoprzestrzenną ramkę danych w EPSG: 4326 i zrobiłbym nową ramkę danych składającą się ze wszystkich wierszy, które mieszczą się w określonej ramce granicznej.

Najpierw otrzymuję obwiednię, na której mi zależy (która w rzeczywistości jest obwiednią innej ramki danych):

print df_sussex.total_bounds
[ -1.57239292  50.57467674   0.14528384  51.27465152]

Następnie tworzę ramkę danych składającą się tylko z tego obwiedni:

pts = gpd.GeoDataFrame(df_sussex.total_bounds)

I w końcu staram się uzyskać wszystkie funkcje, które przecinają się z obwiednią:

sac_sussex = gpd.overlay(pts, df_sac, how='intersection')

Ale to mi daje AttributeError: No geometry data set yet (expected in column 'geometry'.

Co ja robię źle?

Richard
źródło
Problem polega na tym, że używasz metody „total_bounds”. Wytwarza tylko krotkę z maksymalnymi i minimalnymi punktami obwiedni. Metodą, którą należy zastosować, jest „koperta”; poprzedni, aby zbudować odpowiednią GeoDataFrame .
xunilk

Odpowiedzi:

6

Problem polega na tym, że używasz metody „total_bounds”. Wytwarza tylko krotkę z maksymalnymi i minimalnymi punktami obwiedni. Metodą, którą należy zastosować, jest „koperta”; poprzedni, aby zbudować odpowiednią „GeoDataFrame”. Na przykład czytanie moich plików kształtów jako GeoDataFrame :

import geopandas as gpd
pol1 = gpd.GeoDataFrame.from_file("pyqgis_data/polygon1.shp")
pol8 = gpd.GeoDataFrame.from_file("pyqgis_data/polygon8.shp")

Budowanie obwiedni pol1 i tworzenie odpowiedniej GeoDataFrame :

bounding_box = pol1.envelope
df = gpd.GeoDataFrame(gpd.GeoSeries(bounding_box), columns=['geometry'])

Przecięcie obu GeoDataFrame :

intersections = gpd.overlay(df, pol8, how='intersection')

Wyniki wykreślania:

from matplotlib import pyplot as plt
plt.ion()
intersections.plot() 

wprowadź opis zdjęcia tutaj

Działało zgodnie z oczekiwaniami.

Uwaga do edycji:

Korzystając z metody „total_bounds” (ponieważ metoda „koperty” zwraca ramkę graniczną dla każdej cechy wielokątów), można zastosować następujące podejście:

from matplotlib import pyplot as plt
import geopandas as gpd
from shapely.geometry import Point, Polygon

pol1 = gpd.GeoDataFrame.from_file("pyqgis_data/polygon1.shp")
pol8 = gpd.GeoDataFrame.from_file("pyqgis_data/polygon8.shp")

bbox = pol1.total_bounds

p1 = Point(bbox[0], bbox[3])
p2 = Point(bbox[2], bbox[3])
p3 = Point(bbox[2], bbox[1])
p4 = Point(bbox[0], bbox[1])

np1 = (p1.coords.xy[0][0], p1.coords.xy[1][0])
np2 = (p2.coords.xy[0][0], p2.coords.xy[1][0])
np3 = (p3.coords.xy[0][0], p3.coords.xy[1][0])
np4 = (p4.coords.xy[0][0], p4.coords.xy[1][0])

bb_polygon = Polygon([np1, np2, np3, np4])

df2 = gpd.GeoDataFrame(gpd.GeoSeries(bb_polygon), columns=['geometry'])

intersections2 = gpd.overlay(df2, pol8, how='intersection')

plt.ion()
intersections2.plot()

a wynik jest identyczny.

Xunilk
źródło
21

Możesz użyć tej cxmetody na ramce geodata, aby wybrać wiersze w obwiedni. Dla przykładowych ramek:

xmin, ymin, xmax, ymax = df_sussex.total_bounds
sac_sussex = df_sac.cx[xmin:xmax, ymin:ymax]

Ze strony http://geopandas.org/indexing.html :

Oprócz standardowych metod pand GeoPandas zapewnia także indeksowanie oparte na współrzędnych za pomocą indeksu cx , który wycina za pomocą obwiedni. Zwrócone zostaną geometrie w GeoSeries lub GeoDataFrame, które przecinają obwiednię.

jdmcbr
źródło
To rozwiązanie działało dla mnie. Dziękuję Ci. Zastanawiałem się jednak, czy istnieje szybszy sposób wdrożenia. Filtrowanie zagospodarowania przestrzennego OSM i miejsc, które mieszczą się w polu granicznym prowincji.
EFL
Zauważ, że .cxrobi coś nieco innego niż gpd.overlayrozwiązanie: wybiera wiersze, które przecinają obwiednię, ale pozostawiają geometrie nienaruszone, podczas gdy gpd.overlayrozwiązanie zwróci tylko części geometrii w obwiedni. W zależności od sytuacji możesz chcieć jednego lub drugiego.
danvk