我怎样才能让我的 point_in_polygon 函数更加 Pythonic？

问题描述

我从 2 个月开始学习 Python，由于我来自更结构化的编程语言背景（C/Java、JavaScript），我倾向于以结构化的方式而不是 Pythonic 的方式编写所有代码。

我使用 here 中的“光线投射算法”复制了 point_in_polygon 函数，该算法是用 JavaScript 编写的，并遵循完美的“结构化编程风格”。

不过，我把它改成了 Python：

def point_in_polygon(lat,lon,polygon):
    x = lon
    y = lat

    inside = False
    for edge in polygon:
        # i vertex
        xi = edge[1]
        yi = edge[0]

        # j vertex
        xj = edge[3]
        yj = edge[2]

        intersect = ((yi > y) != (yj > y)) and (
            x < (xj - xi) * (y - yi) / (yj - yi) + xi
        )

        # if intersections are an odd number,then the point is inside the polygon
        if intersect is True:
            inside = not inside

    return inside

但它并不能满足我对 Python 的需求。有什么办法可以让这个函数看起来更pythonic？甚至更好地使用 numpy 包？

编辑*

多边形输入示例：

polygon = [
    [39.393037,-0.391826,39.403516,-0.335521],[39.403516,-0.335521,39.424338,-0.298442],[39.424338,-0.298442,39.456420,-0.279903],[39.456420,-0.279903,39.467950,-0.321788],[39.467950,-0.321788,39.458806,-0.339469],[39.458806,-0.339469,39.393037,-0.391826],]

列表内的每个子列表代表多边形的一条边，子列表内有 [lat1,long1,lat2,long2]。

解决方法

您可以使用 numpy 来“Pythonize”您的函数：

def point_in_polygon(lon,lat,polygon):
    p = np.array(polygon)
    inside = ((p[:,1] > lat) != (p[:,3] > lat)) &\
        (lon < ((p[:,2] - p[:,0]) * (lat - p[:,1]) /
         (p[:,3] - p[:,1]) + p[:,0]))

    return (np.sum(inside) % 2 == 1)

请注意，我反转了您的 lat 和 lon 参数，因为给出 (x,y) 而不是 (y,x) 更自然（至少对我而言）。此外，我认为您在函数中颠倒了 x 和 y（或者我可能对经度和纬度感到困惑）。

这是一个 MWE :

import numpy as np
import matplotlib.pyplot as plt


def point_in_polygon(lon,0]))

    return (np.sum(inside) % 2 == 1)


polygon = [
    [39.393037,-0.391826,39.403516,-0.335521],[39.403516,-0.335521,39.424338,-0.298442],[39.424338,-0.298442,39.456420,-0.279903],[39.456420,-0.279903,39.467950,-0.321788],[39.467950,-0.321788,39.458806,-0.339469],[39.458806,-0.339469,39.393037,-0.391826],]

tests = [(round(x,4),round(y,4)) 
         for x in np.linspace(39.39,39.47,25)
         for y in np.linspace(-0.28,-0.40,25)]

polygons = polygon[:]
polygons.append(polygons[0])  # To have a closed polygon when plotting
p = np.array(polygons)  # To make plotting easier
plt.plot(p[:,0],p[:,1])  # plot polygon
for tt in tests:
    inside = point_in_polygon(*tt,polygon)
    # plot test points,with color whether they are within the polygon or not:
    plt.plot(*tt,ls="None",marker="o",color="green" if inside else "red")
plt.show()

coding-style python