圆与圆的交点

Question

如何计算两个圆的交点。我希望在所有情况下都会有两个、一个或没有交叉点。

我有中心点的 x 和 y 坐标以及每个圆的半径。

python 中的答案是首选，但任何工作算法都是可以接受的。

Answer 1

两个圆的交点

作者：Paul Bourke

下面的注释描述了如何找到交点
在平面上的两个圆之间，使用以下符号。这
目的是找到两个点 P₃ = (x₃,
y₃)（如果存在）。

首先计算中心之间的距离d
的圈子。 d = ||P₁ - P₀||。

• 如果d> r₀ + r₁ 则无解，
  圆圈是分开的。

• 如果 d < |r<子>0 -
  r<子>1|那么就没有解，因为一个圆是
  包含在另一个中。

• 如果 d = 0 且 r₀ =
  r₁ 则两个圆重合，并且有
  无限多个解决方案。

考虑两个三角形 P₀P₂P₃
和 P₁P₂P₃ 我们可以写

a² + h² = r₀² 并且
b² + h² = r₁²

使用 d = a + b 我们可以求解 a，

a =
(r₀² - r₁² +
d² ) / (2 d)

可以很容易地证明这可以简化为
r₀ 当两个圆相交于一点时，即： d =
r₀ + r₁

通过将 a 代入第一个来求解 h
方程，h² = r₀² - a²

所以

P₂ = P₀ + a ( P₁ -
P₀ ) / d

最后，P₃ =
(x₃,y₃) 表示为 P₀ =
(x₀,y₀), P₁ =
(x₁,y₁) 和 P₂ =
(x₂,y₂), 是

x₃ =
x₂ +- h ( y₁ - y₀ ) / d

y₃ = y₂ -+ h ( x₁ - x₀ ) /
d

来源：http://paulbourke.net/geometry/circlesphere/

Answer 2

为什么不直接使用您最喜欢的程序语言（或可编程计算器！）的 7 行代码，如下所示。

假设给定 P0 坐标 (x0,y0)、P1 坐标 (x1,y1)、r0 和 r1，并且您想要找到 P3 坐标 (x3,y3)：

d=sqr((x1-x0)^2 + (y1-y0)^2)
a=(r0^2-r1^2+d^2)/(2*d)
h=sqr(r0^2-a^2)
x2=x0+a*(x1-x0)/d   
y2=y0+a*(y1-y0)/d   
x3=x2+h*(y1-y0)/d       // also x3=x2-h*(y1-y0)/d
y3=y2-h*(x1-x0)/d       // also y3=y2+h*(x1-x0)/d

Answer 3

这是我基于 Paul Bourke 的文章 的 C++ 实现。它仅在有两个交集时才有效，否则它可能返回 NaN NAN NAN NAN。

class Point{
    public:
        float x, y;
        Point(float px, float py) {
            x = px;
            y = py;
        }
        Point sub(Point p2) {
            return Point(x - p2.x, y - p2.y);
        }
        Point add(Point p2) {
            return Point(x + p2.x, y + p2.y);
        }
        float distance(Point p2) {
            return sqrt((x - p2.x)*(x - p2.x) + (y - p2.y)*(y - p2.y));
        }
        Point normal() {
            float length = sqrt(x*x + y*y);
            return Point(x/length, y/length);
        }
        Point scale(float s) {
            return Point(x*s, y*s);
        }
};

class Circle {
    public:
        float x, y, r, left;
        Circle(float cx, float cy, float cr) {
            x = cx;
            y = cy;
            r = cr;
            left = x - r;
        }
        pair<Point, Point> intersections(Circle c) {
            Point P0(x, y);
            Point P1(c.x, c.y);
            float d, a, h;
            d = P0.distance(P1);
            a = (r*r - c.r*c.r + d*d)/(2*d);
            h = sqrt(r*r - a*a);
            Point P2 = P1.sub(P0).scale(a/d).add(P0);
            float x3, y3, x4, y4;
            x3 = P2.x + h*(P1.y - P0.y)/d;
            y3 = P2.y - h*(P1.x - P0.x)/d;
            x4 = P2.x - h*(P1.y - P0.y)/d;
            y4 = P2.y + h*(P1.x - P0.x)/d;

            return pair<Point, Point>(Point(x3, y3), Point(x4, y4));
        }

};

Answer 4

这是使用向量的 Javascript 实现。该代码有详细记录，您应该能够遵循它。这是原始来源

查看现场演示此处：

// Let EPS (epsilon) be a small value
var EPS = 0.0000001;

// Let a point be a pair: (x, y)
function Point(x, y) {
  this.x = x;
  this.y = y;
}

// Define a circle centered at (x,y) with radius r
function Circle(x,y,r) {
  this.x = x;
  this.y = y;
  this.r = r;
}

// Due to double rounding precision the value passed into the Math.acos
// function may be outside its domain of [-1, +1] which would return
// the value NaN which we do not want.
function acossafe(x) {
  if (x >= +1.0) return 0;
  if (x <= -1.0) return Math.PI;
  return Math.acos(x);
}

// Rotates a point about a fixed point at some angle 'a'
function rotatePoint(fp, pt, a) {
  var x = pt.x - fp.x;
  var y = pt.y - fp.y;
  var xRot = x * Math.cos(a) + y * Math.sin(a);
  var yRot = y * Math.cos(a) - x * Math.sin(a);
  return new Point(fp.x+xRot,fp.y+yRot);
}

// Given two circles this method finds the intersection
// point(s) of the two circles (if any exists)
function circleCircleIntersectionPoints(c1, c2) {

  var r, R, d, dx, dy, cx, cy, Cx, Cy;

  if (c1.r < c2.r) {
    r  = c1.r;  R = c2.r;
    cx = c1.x; cy = c1.y;
    Cx = c2.x; Cy = c2.y;
  } else {
    r  = c2.r; R  = c1.r;
    Cx = c1.x; Cy = c1.y;
    cx = c2.x; cy = c2.y;
  }

  // Compute the vector <dx, dy>
  dx = cx - Cx;
  dy = cy - Cy;

  // Find the distance between two points.
  d = Math.sqrt( dx*dx + dy*dy );

  // There are an infinite number of solutions
  // Seems appropriate to also return null
  if (d < EPS && Math.abs(R-r) < EPS) return [];

  // No intersection (circles centered at the 
  // same place with different size)
  else if (d < EPS) return [];

  var x = (dx / d) * R + Cx;
  var y = (dy / d) * R + Cy;
  var P = new Point(x, y);

  // Single intersection (kissing circles)
  if (Math.abs((R+r)-d) < EPS || Math.abs(R-(r+d)) < EPS) return [P];

  // No intersection. Either the small circle contained within 
  // big circle or circles are simply disjoint.
  if ( (d+r) < R || (R+r < d) ) return [];

  var C = new Point(Cx, Cy);
  var angle = acossafe((r*r-d*d-R*R)/(-2.0*d*R));
  var pt1 = rotatePoint(C, P, +angle);
  var pt2 = rotatePoint(C, P, -angle);
  return [pt1, pt2];

}

Answer 5

试试这个；

    def ri(cr1,cr2,cp1,cp2):
        int1=[]
        int2=[]
        ori=0
        if cp1[0]<cp2[0] and cp1[1]!=cp2[1]:
            p1=cp1
            p2=cp2
            r1=cr1
            r2=cr2
            if cp1[1]<cp2[1]:
                ori+=1
            elif cp1[1]>cp2[1]:
                ori+=2        
        elif cp1[0]>cp2[0] and cp1[1]!=cp2[1]:
            p1=cp2
            p2=cp1
            r1=cr2
            r2=cr1
            if p1[1]<p2[1]:
                ori+=1
            elif p1[1]>p2[1]:
                ori+=2
        elif cp1[0]==cp2[0]:
            ori+=4
            if cp1[1]>cp2[1]:
                p1=cp1
                p2=cp2
                r1=cr1
                r2=cr2
            elif cp1[1]<cp2[1]:
                p1=cp2
                p2=cp1
                r1=cr2
                r2=cr1
        elif cp1[1]==cp2[1]:
            ori+=3
            if cp1[0]>cp2[0]:
                p1=cp2
                p2=cp1
                r1=cr2
                r2=cr1
            elif cp1[0]<cp2[0]:
                p1=cp1
                p2=cp2
                r1=cr1
                r2=cr2
        if ori==1:#+
            D=calc_dist(p1,p2)
            tr=r1+r2
            el=tr-D
            a=r1-el
            b=r2-el
            A=a+(el/2)
            B=b+(el/2)
            thta=math.degrees(math.acos(A/r1))
            rs=p2[1]-p1[1]
            rn=p2[0]-p1[0]
            gd=rs/rn
            yint=p1[1]-((gd)*p1[0])
            dty=calc_dist(p1,[0,yint])

            aa=p1[1]-yint
            bb=math.degrees(math.asin(aa/dty))
            d=90-bb
            e=180-d-thta
            g=(dty/math.sin(math.radians(e)))*math.sin(math.radians(thta))
            f=(g/math.sin(math.radians(thta)))*math.sin(math.radians(d))
            oty=yint+g
            h=f+r1
            i=90-e
            j=180-90-i
            l=math.sin(math.radians(i))*h
            k=math.cos(math.radians(i))*h
            iy2=oty-l
            ix2=k
            int2.append(ix2)
            int2.append(iy2)

            m=90+bb
            n=180-m-thta
            p=(dty/math.sin(math.radians(n)))*math.sin(math.radians(m))
            o=(p/math.sin(math.radians(m)))*math.sin(math.radians(thta))
            q=p+r1
            r=90-n
            s=math.sin(math.radians(r))*q
            t=math.cos(math.radians(r))*q
            otty=yint-o
            iy1=otty+s
            ix1=t
            int1.append(ix1)
            int1.append(iy1)
        elif ori==2:#-
            D=calc_dist(p1,p2)
            tr=r1+r2
            el=tr-D
            a=r1-el
            b=r2-el
            A=a+(el/2)
            B=b+(el/2)
            thta=math.degrees(math.acos(A/r1))
            rs=p2[1]-p1[1]
            rn=p2[0]-p1[0]
            gd=rs/rn
            yint=p1[1]-((gd)*p1[0])
            dty=calc_dist(p1,[0,yint])

            aa=yint-p1[1]
            bb=math.degrees(math.asin(aa/dty))
            c=180-90-bb
            d=180-c-thta
            e=180-90-d
            f=math.tan(math.radians(e))*p1[0]
            g=math.sqrt(p1[0]**2+f**2)
            h=g+r1
            i=180-90-e
            j=math.sin(math.radians(e))*h
            jj=math.cos(math.radians(i))*h
            k=math.cos(math.radians(e))*h
            kk=math.sin(math.radians(i))*h
            l=90-bb
            m=90-e
            tt=l+m+thta
            n=(dty/math.sin(math.radians(m)))*math.sin(math.radians(thta))
            nn=(g/math.sin(math.radians(l)))*math.sin(math.radians(thta))
            oty=yint-n
            iy1=oty+j
            ix1=k
            int1.append(ix1)
            int1.append(iy1)

            o=bb+90
            p=180-o-thta
            q=90-p
            r=180-90-q
            s=(dty/math.sin(math.radians(p)))*math.sin(math.radians(o))
            t=(s/math.sin(math.radians(o)))*math.sin(math.radians(thta))
            u=s+r1
            v=math.sin(math.radians(r))*u
            vv=math.cos(math.radians(q))*u
            w=math.cos(math.radians(r))*u
            ww=math.sin(math.radians(q))*u
            ix2=v
            otty=yint+t
            iy2=otty-w
            int2.append(ix2)
            int2.append(iy2)

        elif ori==3:#y
            D=calc_dist(p1,p2)
            tr=r1+r2
            el=tr-D
            a=r1-el
            b=r2-el
            A=a+(el/2)
            B=b+(el/2)
            b=math.sqrt(r1**2-A**2)
            int1.append(p1[0]+A)
            int1.append(p1[1]+b)
            int2.append(p1[0]+A)
            int2.append(p1[1]-b)
        elif ori==4:#x
            D=calc_dist(p1,p2)
            tr=r1+r2
            el=tr-D
            a=r1-el
            b=r2-el
            A=a+(el/2)
            B=b+(el/2)
            b=math.sqrt(r1**2-A**2)
            int1.append(p1[0]+b)
            int1.append(p1[1]-A)
            int2.append(p1[0]-b)
            int2.append(p1[1]-A)
        return [int1,int2]
    def calc_dist(p1,p2):
        return math.sqrt((p2[0] - p1[0]) ** 2 +
                         (p2[1] - p1[1]) ** 2)

Answer 6

确定交点 P1 和 P2 的快速方法是获取圆心 A 和 B 之间的向量 w。

现在想象一个由点 A 和 P1 之间的两个向量 a 和 b 组成的矩形，我们可以说

P1 = A + a + b
P2 = A + a - b

唯一剩下的问题是，向量 a 和 b 有多长：

我们知道 |a|^2 + |b|^2 = r_A^2 和 (|w| - |a|)^2 + |b|^2 = r_B^2，将它们设置为相等并求解 < code>|a| 产生

|a| = (r_A^2 - r_b^2 + |w|^2) / (2|w|)
|b| = |b| = +- sqrt(r_A^2 - |a|^2)

现在向量的长度可用于通过使用归一化向量 w 来构造实际向量 a 和 b。

当实现这个想法时，解决方案非常简单

w = {
    x: B.x - A.x,
    y: B.y - A.y
}
    
d = hypot(w.x, w.y)

if (d <= A.r + B.r && abs(B.r - A.r) <= d) {

    w.x/= d;
    w.y/= d;

    a = (A.r * A.r - B.r * B.r + d * d) / (2 * d);
    b = Math.sqrt(A.r * A.r - a * a);

    P1 = {
        x: A.x + a * w.x - b * w.y,
        y: A.y + a * w.y + b * w.x
    }

    P2 = {
        x: A.x + a * w.x + b * w.y,
        y: A.y + a * w.y - b * w.x
    }

} else {
    P1 = P2 = null
}