使用Shapely查找多边形的最大内接矩形

Question

使用Shapely查找多边形的最大内接矩形

12

我正在尝试在六边形中定位数百万个点。这是我的代码：

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"

我已经阅读了其他关于使用shapely改善点多边形查询性能的问题，他们建议使用Rtree。然而，这似乎只对有很多多边形的情况有用（一个问题中有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)

是否也检查一下最大内接矩形（或者可能是三角形）这样的形状的简化版本会有用呢？

查看了Shapely文档，似乎没有针对此操作的函数。也许可以通过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%的点查询。

EDIT 2：GIS StackExchange上另一个相关的问题：Python效率 —— 需要关于如何更有效地使用OGR和Shapely的建议。这涉及在约3,000个多边形中处理1.5百万个点。

- Martin Burch

1

这篇文章介绍了如何将大型多边形分割成小块以加快计算交点的速度。 - Gary Liao

感谢@GaryLiao，Joshua Arnott的“katana”函数（沿最短维度递归地将每个多边形分成两半）是一个非常令人愉悦的想法。如果您有一个实现，在其中展示如何在交集查询（特别是点与多边形）中使用它，请分享，因为该博客并没有清楚地说明如何使用结果拆分的多边形。 - Martin Burch

通过katana，生成了两种矩形：（A）在多边形“内部”的矩形和（B）与多边形边缘“相交”的矩形。通过索引纬度和经度，很容易找到矩形内的点。在几秒钟内可以找到A和B中都包含的点。我们确信A中的点也在多边形内，所以我们要做的最后一件事就是顺序检查B中的点。 - Gary Liao

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Markus Schumann · Accepted Answer

我会使用R-Tree。但我会将所有的点（而不是多边形的边界框）插入到R-Tree中。

例如使用R-Tree索引：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树上得到半打查询，而不是在小R树上进行数百万次查询。