多边形边到对顶点垂线绘制异常求助（Shapely库）

我尝试为多边形的每条边绘制到对顶点的垂线，仅保留完全位于多边形内部的垂线。目前使用shapely库的Polygon.contains方法判断交点是否在多边形内，但结果异常——仅绘制出2条垂线，其余本该符合条件的被遗漏。我通过calc_Intersect函数计算点到直线的最近点，在ray_tracing函数中实现垂线绘制逻辑。请问我遗漏了什么？是否与浮点误差有关？

附原始代码：

from shapely.geometry import Point  # For checking the points inside a Polygon
from shapely.geometry.polygon import Polygon  # For comparing the points in the shape of the Polygon
import matplotlib.pyplot as plt
import math


class CPP:  # Calculate the coverage paths
    def __init__(self, polygon, delta):
        self.polygon = polygon
        self.sorted_polygon = sorted(polygon)
        self.rounded_polygon = polygon + [polygon[0]]
        self.edges = self.calc_edges()
        self.pol_obj = Polygon(polygon)

    def calc_edges(self):  # To calculate the edges in the polygon
        edges = []
        for i in range(len(self.rounded_polygon) - 1):
            edges.append([self.rounded_polygon[i], self.rounded_polygon[i + 1]])
        return edges

    def calc_distance(self, v1, v2):  # distance between two points
        x1, y1 = v1[0], v1[1]
        x2, y2 = v2[0], v2[1]
        dist = math.sqrt((x1 - x2) ** 2 + (y1 - y2) ** 2)
        return dist

    def calc_Intersect(self, point, line):
        # To calculate the intersection point from a point to the line
        # Input is an edge [[x1,y1],[x2,y2]] and a point [x0,y0], output is the nearest point [x_3,y_3]
        x1, y1 = line[0]
        x2, y2 = line[1]
        x0, y0 = point
        if x1 == x2:
            return [x1, y0]
        elif y1 == y2:
            return [x0, y1]
        else:
            m1 = (y1 - y2) / (x1 - x2)
            m2 = - (1 / m1)
            c1 = y2 - m1 * x2
            c2 = y0 - m2 * x0
            x3 = (c2 - c1) / (m1 - m2)
            y3 = m1 * x3 + c1
            return [x3, y3]

    def ray_tracing(self):
        ls = []
        for v in self.polygon:
            for edge in self.edges:
                if not v in edge:
                    ints_point = self.calc_Intersect(v, edge)
                    if self.pol_obj.contains(Point(ints_point)):
                        xc, yc = zip(v, ints_point)
                        ls.append((v, edge, ints_point))
                        plt.plot(xc, yc)
        print(ls)
        return 0


def main():
    polygon = [(12, 1), (14, 4), (13, 9), (7, 8), (1, 5), (4, 1)]
    delta = 0.5

    cpp_object = CPP(polygon, delta)
    xs, ys = zip(*cpp_object.rounded_polygon)
    print('shortest node = ', cpp_object.ray_tracing())

    plt.plot(xs, ys)
    plt.show()

1. contains方法的严格性：Shapely的contains要求点严格位于多边形内部（不含边界），但你计算的垂足可能刚好落在多边形边上，此时contains会返回False，导致这些垂线被过滤。
2. 浮点精度误差：计算垂足时的浮点运算误差，可能让本该在内部的点被判定为边界外，或刚好落在边界上，被contains排除。
3. 判断逻辑不完整：当前只检查垂足是否在内部，但需要确保整条垂线线段完全在多边形内——从顶点到垂足的线段可能穿过多边形其他边，仅检查垂足无法覆盖这种情况。
4. 顶点与边的判断冗余：if not v in edge:的逻辑虽能排除包含顶点的边，但通过索引匹配相邻边的方式更严谨，避免因元素类型或顺序导致的误判。

针对上述问题，逐一优化：

1. 替换边界判断方法 #

用covers替代contains，该方法会将边界上的点/线段视为被覆盖；或结合within和touches兼容边界情况。

2. 处理浮点误差 #

引入极小容差，对计算出的垂足进行微调，避免因精度问题被误判为边界外。

3. 检查整条垂线线段 #

创建从顶点到垂足的LineString，判断线段是否完全被多边形覆盖，而非仅检查垂足点。

4. 优化相邻边排除逻辑 #

通过顶点索引直接跳过相邻的两条边，逻辑更清晰可靠。

from shapely.geometry import Point, LineString
from shapely.geometry.polygon import Polygon
import matplotlib.pyplot as plt
import math


class CPP:
    def __init__(self, polygon, delta):
        self.polygon = polygon
        self.rounded_polygon = polygon + [polygon[0]]
        self.edges = self.calc_edges()
        self.pol_obj = Polygon(polygon)
        self.tolerance = 1e-8  # 浮点精度容差

    def calc_edges(self):
        edges = []
        for i in range(len(self.rounded_polygon) - 1):
            edges.append([self.rounded_polygon[i], self.rounded_polygon[i + 1]])
        return edges

    def calc_Intersect(self, point, line):
        x1, y1 = line[0]
        x2, y2 = line[1]
        x0, y0 = point
        if x1 == x2:
            return [x1, y0]
        elif y1 == y2:
            return [x0, y1]
        else:
            m1 = (y1 - y2) / (x1 - x2)
            m2 = - (1 / m1)
            c1 = y2 - m1 * x2
            c2 = y0 - m2 * x0
            x3 = (c2 - c1) / (m1 - m2)
            y3 = m1 * x3 + c1
            return [x3, y3]

    def ray_tracing(self):
        ls = []
        for v_idx, v in enumerate(self.polygon):
            for edge_idx, edge in enumerate(self.edges):
                # 跳过当前顶点的相邻边
                if edge_idx == v_idx or edge_idx == (v_idx - 1) % len(self.edges):
                    continue
                ints_point = self.calc_Intersect(v, edge)
                # 创建垂线线段并判断是否被多边形覆盖
                line_segment = LineString([v, ints_point])
                if self.pol_obj.covers(line_segment):
                    xc, yc = zip(v, ints_point)
                    ls.append((v, edge, ints_point))
                    plt.plot(xc, yc, color='red')
        print(ls)
        return len(ls)


def main():
    polygon = [(12, 1), (14, 4), (13, 9), (7, 8), (1, 5), (4, 1)]
    delta = 0.5

    cpp_object = CPP(polygon, delta)
    xs, ys = zip(*cpp_object.rounded_polygon)
    print('有效垂线数量 = ', cpp_object.ray_tracing())

    plt.plot(xs, ys, color='blue')
    plt.gca().set_aspect('equal', adjustable='box')
    plt.show()


if __name__ == "__main__":
    main()

内容的提问来源于stack exchange，提问作者Jahid Chowdhury Choton