如何轉換大熊貓經緯度點並查看它們是否落在某些邊界多邊形中？

我有一個熊貓數據幀df這樣的：如何轉換大熊貓經緯度點並查看它們是否落在某些邊界多邊形中？

id lat lon 
jhg 2.7 3.5 
ytr 3.1 3.5 
...

我也有一個Geopandas數據幀poly一些多邊形。現在，我想只繪製點df中的內部的一些多邊形。所以我應該可以做點像poly.intersects(p)，其中p是Shapely Point。但我做錯了什麼;

from shapely.geometry import Point 
for index, row in df.iterrows(): 
    t = poly.intersects(Point(row.lon, row.lat))

什麼是通過與經/緯度點的數據框並繪製它們疊加到poly最好的方法是什麼？請注意，我可以定義一個最小/最大緯度/經度範圍，但也可以打印poly以外的點，但位於（較大）邊界框內。

來源

2016-09-14 Dervin Thunk

你的出發點：

import pandas as pd 
from shapely.geometry import box 
import matplotlib.pyplot as plt 

from matplotlib.collections import PatchCollection 
from matplotlib.patches import Polygon 
from shapely.geometry import Point 
import seaborn as sns 
import numpy as np 

# some pretend data 
data = {'lat':[2.7,3.5,1.4,2.3,.9,1.9], 'lon':[1.2,.9,1.9,2.2,3,1.1]} 
df = pd.DataFrame(data) 

# the 'bounding' polygon 
poly = box(1,1,2,2) 
patches = PatchCollection([Polygon(poly.exterior)], facecolor='red', linewidth=.5, alpha=.5) 


# plot the bounding box 
fig, ax = sns.plt.subplots(1, figsize=(4,4)) 
ax.add_collection(patches, autolim=True) 

# plot the lat/lon points 
df.plot(x='lat',y='lon', kind='scatter',ax=ax) 
plt.show()

的數字看起來是這樣的：

一種方式來擺脫不必要的點是使用一個布爾面膜：

#probably more efficient ways to do this, but this works 
mask = [poly.intersects(Point(lat,lon)) for lat,lon in zip(df.lat,df.lon)] 
df = df[mask] 

# make new plot (must make a new 'patch' object) 
patches1 = PatchCollection([Polygon(poly.exterior)], facecolor='red', linewidth=.5, alpha=.5) 
fig1, ax1 = sns.plt.subplots(1, figsize=(4,4)) 
ax1.add_collection(patches1, autolim=True) 

# make the axis bounds the same 
ax1.set_xlim(ax.get_xlim()) 
ax1.set_ylim(ax.get_ylim()) 

# plot the lat/lon points 
df.plot(x='lat',y='lon', kind='scatter',ax=ax1) 
plt.show()

將此圖像給我。

請注意，你可以做其他，更快的方式一個布爾面具等是否LAT是上述多邊形中的最高點。這些可能並不完美，但可以減少問題，因此您沒有多次致電intersects()。

[編輯：如果你的多邊形是一個矩形，]另一種方式（就像你在你的問題中所建議的那樣）將只是圍繞邊界多邊形「裁剪」圖像。這是一個更快的解決方案，因爲您不必一次又一次調用該函數。要修剪基於邊界多邊形的圖像，你可以plt.plot()之前插入這一權利：

ax.set_xlim((np.min(poly.exterior.xy[0]),np.max(poly.exterior.xy[0]))) 
ax.set_ylim((np.min(poly.exterior.xy[1]),np.max(poly.exterior.xy[1])))

提供了以下：

來源

2016-09-15 01:16:55 benten

這個解決方案不是假設矩形多邊形，這與OP的意圖不符，因爲我讀它... – beroe

@beroe：你是對的。問題是我從一些地理信息中得到不規則的多邊形，所以最後的解決方案（構建一個邊界框）在這一點上對我沒有用處。我目前還沒有能夠完成第一個工作，但這要比我的自己的限制做得更多，而不是OP的答案:) –

謝謝@beroe，相應地改進了我的答案。德文，讓我知道你是否真的陷入了某個地方。 – benten

This tutorial似乎做你想要什麼。它還利用geopandas內置的Rtree空間索引進行快速交叉查詢。

spatial_index = gdf.sindex 
possible_matches_index = list(spatial_index.intersection(polygon.bounds)) 
possible_matches = gdf.iloc[possible_matches_index] 
precise_matches = possible_matches[possible_matches.intersects(polygon)]

然後它繪製多邊形和它內部的點和外部它以不同的顏色，就像你想要做的一樣。

來源

2016-12-12 20:30:36

如何轉換大熊貓經緯度點並查看它們是否落在某些邊界多邊形中？

回答

相關問題