如何实现rtree来计算相交面积?
How to implement rtree to calculate intersecting area?
这是一些示例数据:
import pandas as pd
import geopandas as gp
import shapely.geometry
from shapely.geometry import Polygon
from shapely.geometry import Point
import shapely.affinity
import matplotlib.pyplot as plt
df = gp.GeoDataFrame([['a', Polygon([(0,1), (1,1), (2,2), (1,2)])],
['b', Polygon([(1.5,0.75), (2, 1.25), (3,0.25)])]],
columns = ['name', 'geometry'])
df = gp.GeoDataFrame(df, geometry = 'geometry')
df['area'] = df.area
points = gp.GeoDataFrame([['box', Point(1.2, 1.115), 4],
['triangle', Point(2.5, 1.25), 8]],
columns = ['name', 'geometry', 'value'],
geometry = 'geometry')
buf = points.buffer(0.5, cap_style = 3)
points['buffer'] = buf
points = points.drop(['geometry'], axis = 1)
points = points.rename(columns = {'buffer': 'geometry'})
看起来像这样:
基本上我试图找到这些对象之间的交集区域。
到目前为止,我已经使用以下代码完成了此操作:
data = []
for index, geo in df.iterrows():
for index2, poin in points.iterrows():
if geo['geometry'].intersects(poin['geometry']):
data.append({'geometry':geo['geometry'].intersection(poin['geometry']), 'area': geo['geometry'].intersection(poin['geometry']).area})
df2 = gp.GeoDataFrame(data, columns = ['geometry', 'area'])
但是,我将使用它的真实数据有 100,000 个多边形,因此这段代码将非常耗时。我知道我可以通过使用 r-trees 来加快速度。但是,我似乎无法正确实施它。
我试过这样的事情:
spatial_index = df.sindex
results_list = []
for index, row in points.iterrows():
buffer = row['geometry']
possible_matches_index = list(spatial_index.intersection(buffer.bounds))
possible_matches = df.iloc[possible_matches_index]
results_list.append({'geometry':possible_matches['geometry'].intersection(row['geometry']), 'area': possible_matches['geometry'].intersection(row['geometry']).area})
df = gp.GeoDataFrame(results_list, columns = ['geometry', 'area'])
但这会将每个方块的所有交点都放在一条线上。
geometry area
0 name a POLYGON ((1.615 1.615, 1 1, 0.7 1, 0... name a 0.000037 b 0.000003 dtype: float64
1 name a ... name a 0.000000 b 0.000013 dtype: float64
我怎样才能生成一个数据框,每个交叉点都用一条线表示它的几何形状和面积?
要计算两个 GeoDataFrame 之间的交集,您可以使用 geopandas.overlay
函数:
geopandas.overlay(df, points, how='intersection')
这利用了引擎盖下的 rtree
空间索引(因此应该比蛮力加倍 for-loop 更有效),并将 return 所有组合的交集两个数据集的几何形状作为一个新的 GeoDataFrame(然后您可以计算面积)。
有关这方面的文档,请参阅 https://geopandas.readthedocs.io/en/latest/set_operations.html。
这是一些示例数据:
import pandas as pd
import geopandas as gp
import shapely.geometry
from shapely.geometry import Polygon
from shapely.geometry import Point
import shapely.affinity
import matplotlib.pyplot as plt
df = gp.GeoDataFrame([['a', Polygon([(0,1), (1,1), (2,2), (1,2)])],
['b', Polygon([(1.5,0.75), (2, 1.25), (3,0.25)])]],
columns = ['name', 'geometry'])
df = gp.GeoDataFrame(df, geometry = 'geometry')
df['area'] = df.area
points = gp.GeoDataFrame([['box', Point(1.2, 1.115), 4],
['triangle', Point(2.5, 1.25), 8]],
columns = ['name', 'geometry', 'value'],
geometry = 'geometry')
buf = points.buffer(0.5, cap_style = 3)
points['buffer'] = buf
points = points.drop(['geometry'], axis = 1)
points = points.rename(columns = {'buffer': 'geometry'})
看起来像这样:
基本上我试图找到这些对象之间的交集区域。
到目前为止,我已经使用以下代码完成了此操作:
data = []
for index, geo in df.iterrows():
for index2, poin in points.iterrows():
if geo['geometry'].intersects(poin['geometry']):
data.append({'geometry':geo['geometry'].intersection(poin['geometry']), 'area': geo['geometry'].intersection(poin['geometry']).area})
df2 = gp.GeoDataFrame(data, columns = ['geometry', 'area'])
但是,我将使用它的真实数据有 100,000 个多边形,因此这段代码将非常耗时。我知道我可以通过使用 r-trees 来加快速度。但是,我似乎无法正确实施它。
我试过这样的事情:
spatial_index = df.sindex
results_list = []
for index, row in points.iterrows():
buffer = row['geometry']
possible_matches_index = list(spatial_index.intersection(buffer.bounds))
possible_matches = df.iloc[possible_matches_index]
results_list.append({'geometry':possible_matches['geometry'].intersection(row['geometry']), 'area': possible_matches['geometry'].intersection(row['geometry']).area})
df = gp.GeoDataFrame(results_list, columns = ['geometry', 'area'])
但这会将每个方块的所有交点都放在一条线上。
geometry area
0 name a POLYGON ((1.615 1.615, 1 1, 0.7 1, 0... name a 0.000037 b 0.000003 dtype: float64
1 name a ... name a 0.000000 b 0.000013 dtype: float64
我怎样才能生成一个数据框,每个交叉点都用一条线表示它的几何形状和面积?
要计算两个 GeoDataFrame 之间的交集,您可以使用 geopandas.overlay
函数:
geopandas.overlay(df, points, how='intersection')
这利用了引擎盖下的 rtree
空间索引(因此应该比蛮力加倍 for-loop 更有效),并将 return 所有组合的交集两个数据集的几何形状作为一个新的 GeoDataFrame(然后您可以计算面积)。
有关这方面的文档,请参阅 https://geopandas.readthedocs.io/en/latest/set_operations.html。