用Python Shapely查找多边形中最大的内接矩形
我试图在六个多边形内找到数百万个点。 这是我的代码:
def find_shape(longitude,latitude):
if longitude != 0 and latitude != 0:
point = shapely.geometry.Point(longitude,latitude)
else:
return "Unknown"
for current_shape in all_shapes:
if current_shape['bounding_box'].contains(point):
if current_shape['shape'].contains(point):
return current_shape['properties']['ShapeName']
break
return "Unknown"
我已经阅读了其他一些涉及改善点状多边形查询性能的问题。 他们建议Rtrees。 但是,对于有多个多边形的情况(一个问题中有36,000个,另一个问题中有100,000个)以及不希望遍历所有多边形的情况,这似乎很有用。
正如你所看到的,我已经设置了一个边界框。 这是我的形状设置代码:
with fiona.open(SHAPEFILE) as f_col:
all_shapes = []
for shapefile_record in f_col:
current_shape = {}
current_shape['shape'] = shapely.geometry.asShape(shapefile_record['geometry'])
minx, miny, maxx, maxy = current_shape['shape'].bounds
current_shape['bounding_box'] = shapely.geometry.box(minx, miny, maxx, maxy)
current_shape['properties'] = shapefile_record['properties']
all_shapes.append(current_shape)
检查另一个非常简化的形状 ,也就是由最大的内接矩形(或者可能是三角形) 构成的形状是否有用 ?
检查造型文档,它似乎没有这个功能。 也许一些simplify()
设置? 当然,我总是想确保新的简化形状不会超出原始形状的范围,所以我不必在实际形状上调用contains()
。 我也认为我想尽可能简化新的简化形状,以提高速度。
任何其他建议也赞赏。 谢谢!
编辑 :在等待回复时,我创建了一个符合我的要求的简化形状:
current_shape['simple_shape'] = current_shape['shape'].simplify(.5)
current_shape['simple_shape'] = current_shape['simple_shape'].intersection(current_shape['shape'])
以下是测试每个点时的使用方法:
if current_shape['simple_shape'].contains(point):
return current_shape['properties']['ShapeName']
elif current_shape['shape'].contains(point):
return current_shape['properties']['ShapeName']
这并不完美,因为形状并不像完成必要的intersection()
后那么简单。 尽管如此,这种方法已经使处理时间减少了60%。 在我的测试中,85%的点查询使用简单的多边形。
编辑2 :关于GIS StackExchange上的另一个相关问题:Python效率 - 需要关于如何以更有效的方式使用OGR和Shapely的建议。 这处理约3,000个多边形中的150万个点。
我会使用R-Tree。 但是我会将所有点(而不是多边形的边界框)插入到R-Tree中。
例如使用r树索引:http://toblerity.org/rtree/
from rtree import index
from rtree.index import Rtree
idx = index.Index();
//插入一个点,即left == right && top == bottom,将实质上插入一个单一的条目到索引中
for current_point in all_points:
idx.insert(current_point.id, (current_point.x, current_point.y, current_point.x, current_point.y))
//现在循环你的多边形
for current_shape in all_shapes:
for n in idx.intersect( current_shape['bounding_box'] )
if current_shape['shape'].contains(n):
# now we know that n is inside the current shape
所以你最终会在你的大型R-Tree上找到6个查询,而不是在一个小型R-Tree上进行数百万次查询。
链接地址: http://www.djcxy.com/p/84101.html上一篇: Finding largest inscribed rectangle of a polygon with Python Shapely
下一篇: Git "missing" commit