用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