检查地理点是在多边形内部还是外部

Question

我正在使用 python 并且我已经在地图上定义了多边形的纬度和经度（以度为单位）。我的目标是检查坐标 x,y 的通用点 P 是否落在这样的多边形内。因此，我希望有一个功能可以让我检查这种情况，如果该点在多边形内部或外部，return True 或 False。

在这个例子中，点在外面所以结果是 False

问题: 有library/package可以达到我的目标吗？如果是，你推荐哪一个？您能举个小例子说明如何使用它吗？

这是我到目前为止编写的代码：

import numpy as np

# Define vertices of polygon (lat/lon)
v0 = [7.5, -2.5] 
v1 = [2, 3.5]
v2 = [-2, 4]
v3 = [-5.5, -4]
v4 = [0, -10]
lats_vect = np.array([v0[0],v1[0],v2[0],v3[0],v4[0]])
lons_vect = np.array([v0[1],v1[1],v2[1],v3[1],v4[1]])

# Point of interest P
x, y = -6, 5 # x = Lat, y = Lon

## START MODIFYING FROM HERE; DO NOT MODIFY POLYGON VERTICES AND DATA TYPE
# Check if point of interest falls within polygon boundaries
# If yes, return True
# If no, return False

为了绘制多边形和兴趣点，我使用了 cartopy 并编写了以下代码行：

import cartopy.crs as ccrs
import matplotlib.pyplot as plt
ax = plt.axes(projection=ccrs.PlateCarree())
ax.stock_img() 

# Append first vertex to end of vector to close polygon when plotting
lats_vect = np.append(lats_vect, lats_vect[0])
lons_vect = np.append(lons_vect, lons_vect[0])
plt.plot([lons_vect[0:-1], lons_vect[1:]], [lats_vect[0:-1], lats_vect[1:]],
         color='black', linewidth=1, 
         transform=ccrs.Geodetic(),
         )   

plt.plot(y, x, 
        '*',          # marker shape
        color='blue',  # marker colour
        markersize=8  # marker size
        )  

plt.show()

注意：

点通过大圆相互连接！
我在实习期间进行了研究，最后发现了一些类似的问题 like this one 但我没有成功，因为他们都使用了我没有的 .shp 文件。

Answer 1

这是我的问题的可能解决方案。

必须正确存储地理坐标。示例 np.array([[Lon_A, Lat_A], [Lon_B, Lat_B], [Lon_C, Lat_C]])
创建多边形
创建要测试的点
使用 polygon.contains(point) 测试点是在多边形内部 (True) 还是外部 (False)。

这里是代码的缺失部分：

from shapely.geometry import Point
from shapely.geometry.polygon import Polygon

lons_lats_vect = np.column_stack((lons_vect, lats_vect)) # Reshape coordinates
polygon = Polygon(lons_lats_vect) # create polygon
point = Point(y,x) # create point
print(polygon.contains(point)) # check if polygon contains point
print(point.within(polygon)) # check if a point is in the polygon

注意：多边形不考虑大圆，因此需要将边分割成很多段，增加顶点数

特例：如果点位于多边形的边界上

例如print(Polygon([(0, 0), (1, 0), (1, 1)]).contains(Point(0, 0))) 会失败

所以可以用

print(polygon.touches(point)) # check if point lies on border of polygon

Answer 2

另一种方法是使用 link https://wrf.ecse.rpi.edu//Research/Short_Notes/pnpoly.html The python code is given in wikipedia https://en.wikipedia.org/wiki/Even–odd_rule

中解释的奇偶算法

伙计们，请记住使多边形重要的点的顺序！我的意思是，不同的顺序会产生不同的多边形。

Answer 3

还有一个新兴的 python 库 turfpy。用于地理空间分析。

PyPI

Github

示例：

from turfpy.measurement import boolean_point_in_polygon
from geojson import Point, Polygon, Feature

point = Feature(geometry=Point((-46.6318, -23.5523)))
polygon = Polygon(
    [
        [
            (-46.653, -23.543),
            (-46.634, -23.5346),
            (-46.613, -23.543),
            (-46.614, -23.559),
            (-46.631, -23.567),
            (-46.653, -23.560),
            (-46.653, -23.543),
        ]
    ]
)
boolean_point_in_polygon(point, polygon)

Answer 4

您可以使用 pygeodesy 包，它不依赖于系统级地理包，并且它使用 Kenneth Gade 的 n-vector 方法，

https://github.com/mrJean1/PyGeodesy

就pip install pygeodesy

示例代码

from pygeodesy.sphericalNvector import LatLon

p = LatLon(45.1, 1.1)
b = LatLon(45, 1), LatLon(45, 2), LatLon(46, 2), LatLon(46, 1)
print (p.isenclosedBy(b))

这应该输出 True

检查地理点是在多边形内部还是外部

Check if geo-point is inside or outside of polygon

python

polygon

coordinate-systems

point-in-polygon

shapely