在python中找到直线和圆的交点最有效的方法是什么？

Question

我有一个由许多点组成的多边形.我想找到多边形和圆的交点.提供[x0,y0]的圆心和r0的半径,我写了一个粗略的函数来简单地求解圆和线的二次方程.但是逐个找到多边形的每个线段的交点效率怎么样？有更有效的方法吗？

我知道同情已经提供了获得不同几何形状的交叉点的特征.但是,如果我想处理大量的多边形,那么如同通过我自己的函数计算外部库的效率相比如何？

def LineIntersectCircle(p,lsp,lep):
# p is the circle parameter, lsp and lep is the two end of the line
  x0,y0,r0 = p
  x1,y1 = lsp
  x2,y2 = esp
  if x1 == x2:
    if abs(r0) >= abs(x1 - x0):
        p1 = x1, y0 - sqrt(r0**2 - (x1-x0)**2)
        p2 = x1, y0 + sqrt(r0**2 - (x1-x0)**2)
        inp = [p1,p2]
        # select the points lie on the line segment
        inp = [p for p in inp if p[1]>=min(y1,y2) and p[1]<=max(y1,y2)]
    else:
        inp = []
  else:
    k = (y1 - y2)/(x1 - x2)
    b0 = y1 - k*x1
    a = k**2 + 1
    b = 2*k*(b0 - y0) - 2*x0
    c = (b0 - y0)**2 + x0**2 - r0**2
    delta = b**2 - 4*a*c
    if delta >= 0:
        p1x = (-b - sqrt(delta))/(2*a)
        p2x = (-b + sqrt(delta))/(2*a)
        p1y = k*x1 + b0
        p2y = k*x2 + b0
        inp = [[p1x,p1y],[p2x,p2y]]
        # select the points lie on the line segment
        inp = [p for p in inp if p[0]>=min(x1,x2) and p[0]<=max(x1,x2)]
    else:
        inp = []
  return inp

Answer 1

下面是计算圆与由两个 (x, y) 点定义的直线或线段的交点的解决方案：

def circle_line_segment_intersection(circle_center, circle_radius, pt1, pt2, full_line=True, tangent_tol=1e-9):
    """ Find the points at which a circle intersects a line-segment.  This can happen at 0, 1, or 2 points.

    :param circle_center: The (x, y) location of the circle center
    :param circle_radius: The radius of the circle
    :param pt1: The (x, y) location of the first point of the segment
    :param pt2: The (x, y) location of the second point of the segment
    :param full_line: True to find intersections along full line - not just in the segment.  False will just return intersections within the segment.
    :param tangent_tol: Numerical tolerance at which we decide the intersections are close enough to consider it a tangent
    :return Sequence[Tuple[float, float]]: A list of length 0, 1, or 2, where each element is a point at which the circle intercepts a line segment.

    Note: We follow: http://mathworld.wolfram.com/Circle-LineIntersection.html
    """

    (p1x, p1y), (p2x, p2y), (cx, cy) = pt1, pt2, circle_center
    (x1, y1), (x2, y2) = (p1x - cx, p1y - cy), (p2x - cx, p2y - cy)
    dx, dy = (x2 - x1), (y2 - y1)
    dr = (dx ** 2 + dy ** 2)**.5
    big_d = x1 * y2 - x2 * y1
    discriminant = circle_radius ** 2 * dr ** 2 - big_d ** 2

    if discriminant < 0:  # No intersection between circle and line
        return []
    else:  # There may be 0, 1, or 2 intersections with the segment
        intersections = [
            (cx + (big_d * dy + sign * (-1 if dy < 0 else 1) * dx * discriminant**.5) / dr ** 2,
             cy + (-big_d * dx + sign * abs(dy) * discriminant**.5) / dr ** 2)
            for sign in ((1, -1) if dy < 0 else (-1, 1))]  # This makes sure the order along the segment is correct
        if not full_line:  # If only considering the segment, filter out intersections that do not fall within the segment
            fraction_along_segment = [(xi - p1x) / dx if abs(dx) > abs(dy) else (yi - p1y) / dy for xi, yi in intersections]
            intersections = [pt for pt, frac in zip(intersections, fraction_along_segment) if 0 <= frac <= 1]
        if len(intersections) == 2 and abs(discriminant) <= tangent_tol:  # If line is tangent to circle, return just one point (as both intersections have same location)
            return [intersections[0]]
        else:
            return intersections

Answer 2

我想也许您的问题是理论上如何以最快的方式做到这一点。但是，如果您想快速执行此操作，则应真正使用以C / C ++编写的内容。

我已经习惯了Shapely，所以我将提供一个示例说明如何使用此库。有许多用于python的几何库。我将在此答案的结尾列出它们。

from shapely.geometry import LineString
from shapely.geometry import Point

p = Point(5,5)
c = p.buffer(3).boundary
l = LineString([(0,0), (10, 10)])
i = c.intersection(l)

print i.geoms[0].coords[0]
(2.8786796564403576, 2.8786796564403576)

print i.geoms[1].coords[0]
(7.121320343559642, 7.121320343559642)

在Shapely中有点反直觉的是，圆是具有缓冲区的点的边界。这就是为什么我这样做p.buffer(3).boundry

交集i也是一个几何形状的列表，在这种情况下，它们都是点，这就是为什么我必须从中获得它们的原因i.geoms[]

还有另一个Stackoverflow问题，针对那些感兴趣的人来介绍有关这些库的详细信息。

• 象征
• CGAL Python绑定
• 欧几里德
• PythonOCC
• 几何简单

编辑，因为评论：

Shapely基于GEOS（trac.osgeo.org/geos），它是用C ++构建的，并且比您用python本地编写的任何东西都要快得多。SymPy似乎基于mpmath（mpmath.org），后者也似乎是python，但是似乎其中集成了许多非常复杂的数学。实现自己可能需要很多工作，并且可能不会像GEOS C ++实现那样快。