判断地理坐标点是否在多边形内的Python实现

在地学问题中经常会遇到判断点是否在多边形内的问题，例如成像问题中射线是否穿过某个网格即为此类问题。在Python中我们对处理这样的问题需要满足以下条件：

输入多边形经纬度坐标和待判断点的经纬度坐标，输出True/False
能与numpy完美衔接
输入为多个点也可以处理

用 Shapely 定义多边形#

shapely 是一个处理几何对象的Python库，可以用来定义多边形并判断点是否在多边形内。我们可以使用shapely.geometry模块中的Point和Polygon类来实现。

1
from shapely.geometry import Point
2
from shapely.geometry.polygon import Polygon
3
from io import StringIO
4
import numpy as np
5

6
poly = '''
7
98.2 27
8
101 27
9
102.8 29.2
10
98.2 29.2
11
'''
12

13
def points(lats, lons):
14
    pg = Polygon(np.genfromtxt(StringIO(poly))) # 生成Polygon实例
15
    return np.array([pg.contains(Point(lo, la)) for lo, la in zip(lons, lats)]) # 判断所有输入的点是否在poly内，返回np.array

测试代码#

下面进行一个测试，我们定义两个点分别是(28, 100)和(44, 120)

1
p = np.array([[28, 100], [44, 120]])
2
print(points(p[:, 0], p[:, 1]))

结果为：

1
[ True False]

这样结合np.where函数可以获取原先列表里在或不再多边形内的点。