如何检查一个点是否在三个点的外接圆内？

Question

有没有简单的解决方案？或者有没有任何一个实现的例子？

谢谢，乔纳斯。

- J. Meier

---

为什么不直接计算外接圆，然后检查是否包含？ - Nico Schertler

Answer 1

我们称

• a、b、c为三个点，
• C是（a、b、c）的外接圆，
• d是另一个点。

判断d是否在C内的快速方法是计算此行列式：

      | ax-dx, ay-dy, (ax-dx)² + (ay-dy)² |
det = | bx-dx, by-dy, (bx-dx)² + (by-dy)² |
      | cx-dx, cy-dy, (cx-dx)² + (cy-dy)² |

如果a、b、c按逆时针顺序排列，则：

• 如果行列式等于0，则d在C上。
• 如果行列式大于0，则d在C内部。
• 如果行列式小于0，则d在C外部。

这是一个执行上述操作的JavaScript函数：

function inCircle (ax, ay, bx, by, cx, cy, dx, dy) {
    let ax_ = ax-dx;
    let ay_ = ay-dy;
    let bx_ = bx-dx;
    let by_ = by-dy;
    let cx_ = cx-dx;
    let cy_ = cy-dy;
    return (
        (ax_*ax_ + ay_*ay_) * (bx_*cy_-cx_*by_) -
        (bx_*bx_ + by_*by_) * (ax_*cy_-cx_*ay_) +
        (cx_*cx_ + cy_*cy_) * (ax_*by_-bx_*ay_)
    ) > 0;
}

你可能还需要检查一下你的点是否按逆时针顺序排列：

function ccw (ax, ay, bx, by, cx, cy) {
    return (bx - ax)*(cy - ay)-(cx - ax)*(by - ay) > 0;
}

我没有在inCircle函数中放置ccw检查，因为你不应该每次都进行检查。

这个过程不需要任何除法或平方根运算。

你可以在这里看到代码的实际使用：https://titouant.github.io/testTriangle/

以及源代码在这里：https://github.com/TitouanT/testTriangle

- TitouanT

Answer 2

（如果您对非明显/“疯狂”的解决方案感兴趣。）

Delaunay三角剖分的一个等价属性是：如果您在三角剖分中建立任何三角形的外接圆，它保证不包含任何其他三角剖分的顶点。

Delaunay三角剖分的另一个等价属性是：它最大化了最小三角形角度（即在相同点集上的所有三角剖分中将其最大化）。

这提示了您测试的算法：

1.考虑在原始3个点上构建的三角形ABC。

2.如果测试点P位于三角形内部，则它肯定在圆内。

3.如果测试点P属于其中一个“角落”区域（请参见下图中的阴影区域），则它肯定在圆外。

否则（假设P在区域1中），考虑四边形ABCP的两个三角剖分：原始三角剖分包含原始三角形（红色对角线），而“翻转”对角线的替代三角剖分为蓝色对角线。

1. 通过测试“翻转”条件，比较α = min(∠1,∠4,∠5,∠8)和β = min(∠2,∠3,∠6,∠7)来确定这些三角剖分中哪一个是Delaunay三角剖分。
2. 如果原始三角剖分是Delaunay三角剖分（α > β），则P在圆外。如果备选三角剖分是Delaunay三角剖分（α < β），则P在圆内。

完成。

---

（一段时间后重新访问这个答案。）

这种解决方案可能并不像乍一看那么“疯狂”。请注意，在步骤5和6中比较角度时，没有必要计算实际的角度值。只需知道它们的余弦值即可（即无需涉及三角函数）。

代码的C++版本

#include <cmath>
#include <array>
#include <algorithm>

struct pnt_t
{
  int x, y;

  pnt_t ccw90() const
    { return { -y, x }; }

  double length() const
    { return std::hypot(x, y); }

  pnt_t &operator -=(const pnt_t &rhs)
  {
    x -= rhs.x;
    y -= rhs.y;
    return *this;
  }

  friend pnt_t operator -(const pnt_t &lhs, const pnt_t &rhs)
    { return pnt_t(lhs) -= rhs; }

  friend int operator *(const pnt_t &lhs, const pnt_t &rhs)
    { return lhs.x * rhs.x + lhs.y * rhs.y; }
};

int side(const pnt_t &a, const pnt_t &b, const pnt_t &p)
{
  int cp = (b - a).ccw90() * (p - a);
  return (cp > 0) - (cp < 0);
}

void make_ccw(std::array<pnt_t, 3> &t)
{
  if (side(t[0], t[1], t[2]) < 0)
    std::swap(t[0], t[1]);
}

double ncos(pnt_t a, const pnt_t &o, pnt_t b)
{
  a -= o;
  b -= o;
  return -(a * b) / (a.length() * b.length());
}

bool inside_circle(std::array<pnt_t, 3> t, const pnt_t &p)
{
  make_ccw(t);

  std::array<int, 3> s = 
    { side(t[0], t[1], p), side(t[1], t[2], p), side(t[2], t[0], p) };

  unsigned outside = std::count(std::begin(s), std::end(s), -1);
  if (outside != 1)
    return outside == 0;

  while (s[0] >= 0)
  {
    std::rotate(std::begin(t), std::begin(t) + 1, std::end(t));
    std::rotate(std::begin(s), std::begin(s) + 1, std::end(s));
  }

  double 
    min_org = std::min({
      ncos(t[0], t[1], t[2]), ncos(t[2], t[0], t[1]), 
      ncos(t[1], t[0], p), ncos(p, t[1], t[0]) }),
    min_alt = std::min({
      ncos(t[1], t[2], p), ncos(p, t[2], t[0]), 
      ncos(t[0], p, t[2]), ncos(t[2], p, t[1]) });

  return min_org <= min_alt;
}

以任意选择的三角形和大量随机点进行一些测试。

---

当然，整个问题可以轻松地重新表述，甚至不必提及 Delaunay 三角剖分。从第四步开始，此解决方案基于循环四边形的对顶角性质，其总和必须为180°。

- AnT stands with Russia

---

我正在尝试理解这个问题，因为它回答了我一个问题，但我遇到了一个问题。您提到您的方法基于翻转条件，其中角度之和<= 180度。但在您的示例代码中，您没有使用角度，而是计算余弦值，而不是将它们相加，而是只取所有角度的最小值。这是如何工作的？某些余弦角始终为负数，这不会影响计算吗？取最小余弦如何完成与求和角度相同的任务，以确定哪个是Delaunay？ - Nairou

为了比较角度，我使用 ncos 函数，它代表“负余弦”，即实际余弦的否定。 它在整个 [0, Pi) 范围内均匀地行为：角度越大，则其 ncos 值越大。（虽然在这个问题中，如果我没有漏掉什么，我们只处理锐角。） - AnT stands with Russia

关于翻转条件，它有多个等价的表述方式。维基百科链接中的那个是基于角度之和的。我使用了另一种（但同样等价的）翻转条件形式，我在我的回答的第5点中进行了描述：翻转以最大化最小角度。这个版本的翻转条件不需要对角度求和，甚至不需要知道角度值。它只需要能够比较角度。这就是我选择这个翻转条件形式的全部原因。 - AnT stands with Russia

如果你更深入地研究 Delaunay 三角剖分，你会发现它有很多有趣的特性，这些特性实际上可以相互推导。这些特性可以用来建立不同（但等效）的反转条件。在这里我们看到了最为人熟知的两个条件。 - AnT stands with Russia

Answer 3

在我发布的这篇数学SE帖子中，我包括了一个检查四个点是否共圆的方程式，通过计算4x4行列式来实现。将该方程式转化为不等式可以进行内部检查。

如果您想知道不等式的方向，考虑到远离点的情况。在这种情况下，x²+y²项将支配所有其他项。因此，您可以假设对于所询问的点，此项为一，而其他三项为零。然后选择不等式的符号，以使该值不满足它，因为该点肯定在外面，但您想表述其内部。

如果数字精度是一个问题，Shewchuk教授的这个页面描述了如何获得使用常规双精度浮点数表示的点的一致谓词。

- MvG

Answer 4

假设有三个点(x1,y1)、(x2,y2)、(x3,y3)，并且你想要检查的点(x,y)是否在由上述三个点定义的圆内，你可以这样做：

/**
 *
 * @param x coordinate of point want to check if inside
 * @param y coordinate of point want to check if inside
 * @param cx center x
 * @param cy center y
 * @param r radius of circle
 * @return whether (x,y) is inside circle
 */
static boolean g(double x,double y,double cx,double cy,double r){
    return Math.sqrt((x-cx)*(x-cx)+(y-cy)*(y-cy))<r;
}

// check if (x,y) is inside circle defined by (x1,y1),(x2,y2),(x3,y3)
static boolean isInside(double x,double y,double x1,double y1,double x2,double y2,double x3,double y3){
    double m1 = (x1-x2)/(y2-y1);
    double m2 = (x1-x3)/(y3-y1);
    double b1 = ((y1+y2)/2) - m1*(x1+x2)/2;
    double b2 = ((y1+y3)/2) - m2*(x1+x3)/2;
    double xx = (b2-b1)/(m1-m2);
    double yy = m1*xx + b1;
    return g(x,y,xx,yy,Math.sqrt((xx-x1)*(xx-x1)+(yy-y1)*(yy-y1)));
}

public static void main(String[] args) {
    // if (0,1) is inside the circle defined by (0,0),(0,2),(1,1)
    System.out.println(isInside(0,1,0,0,0,2,1,1));
}

从3个点获取圆心表达式的方法是，找到2条线段垂线中线的交点，我选择了上面的(x1,y1)-(x2,y2)和(x1,y1)-(x3,y3)。既然你知道每个垂线中线上的一点，即(x1+x2)/2和(x1+x3)/2，并且你也知道每个垂线中线的斜率，即(x1-x2)/(y2-y1)和(x1-x3)/(y3-y1)，分别从上述2条线段得出，那么你就可以解出它们相交的(x,y)。

- Bobas_Pett