logo标签列表

连接两条线段

geometry
8

给定两条2D线段A和B，如何计算连接A和B的最短2D线段C的长度？

定义连接。您是指连接它们的端点，还是找到线段上任意一点之间的最短线段？ - strager
@strager，根据欧几里得几何学，它们要么是平行的，要么端点更近，因此您需要检查向量A1-B1，A1-B2，A2-B1和A2-B2。 - paxdiablo
@Pax：我没有看到任何相反的证据。 “无罪推定”去哪了？ - Konrad Rudolph
1
最接近的点必须是其中一条线段的端点之一。我是否可以通过使用从该端点到另一条线段的每个端点的距离来创建一个比值，以便沿着另一条线段进行插值？ - anon
可能是两条线段之间的最短距离的重复问题。 - Veedrac
显示剩余2条评论
7个回答
7
考虑将两条线段 A 和 B 表示为每个点：

线段 A 由 A1(x,y) 和 A2(x,y) 表示

线段 B 由 B1(x,y) 和 B2(x,y) 表示

首先使用此算法检查这两条线段是否相交。

如果它们相交，则两条线之间的距离为零，并且连接它们的线段是交点。

如果它们不相交，请使用此方法：http://paulbourke.net/geometry/pointlineplane/ 计算以下四个线段之间的最短距离：

  1. 点 A1 和线段 B 之间的距离
  2. 点 A2 和线段 B 之间的距离
  3. 点 B1 和线段 A 之间的距离
  4. 点 B2 和线段 A 之间的距离

其中最短的那个线段就是你要找的答案。

此外，首先检查 A 和 B 是否相交（A1 + vector (A1->A2) * a = B1 + vector (B1->B2) * b，其中 a 和 b 是非零实数）。其次，检查 A 和 B 是否平行（即矢量（A1->A2）* a = 矢量（B1->B2），其中 a 是非零实数）。 - Bodo Thiesen
我不相信这个答案是正确的。根据上面的内容，如果两条线相交，它们之间的最短距离为零。但是，任何一个端点到另一条线的最短距离都大于零。 - David Waters
你们两个都是对的：检查交集是必要的。然而，我认为检查线是否平行是不必要的。 - Alterlife
我认为需要进行交点检查的是线段是否相交，而不是延长线（通常会相交）。但是如何确定这一点并不清楚（算法没有给出？），即使在检查点到线（段）的确定性链接后也是如此。 - Rob Parker
@BodoThiesen 您的交点检查需要检查a和b是否在0到1之间，否则交点将超出线段的末端。 - Eyal
4
这个页面提供了您可能需要的信息，与IT技术有关。请点击此链接查看。
一个针对3D的答案同样适用于2D，2D只是3D的一种特殊情况，其中z始终等于0。因此，在下面的伪代码中，z（i）== z（j）== 0。 - David Waters
3

使用Afterlife和Rob Parker上面算法的一般思路，这里提供了一个C++版本的一组方法来获取两个任意2D线段之间的最小距离。它可以处理重叠线段、平行线段、相交和非相交线段。此外，它使用各种epsilon值来防止浮点精度误差。最后，除了返回最小距离外，该算法还将给出线段1中最靠近线段2的点（如果线段相交，则为交点）。如果需要，也很容易返回[p3，p4]上最接近[p1，p2]的点，但我会把它留给读者作为练习 :)

// minimum distance (squared) between vertices, i.e. minimum segment length (squared)
#define EPSILON_MIN_VERTEX_DISTANCE_SQUARED 0.00000001

// An arbitrary tiny epsilon.  If you use float instead of double, you'll probably want to change this to something like 1E-7f
#define EPSILON_TINY 1.0E-14

// Arbitrary general epsilon.  Useful for places where you need more "slop" than EPSILON_TINY (which is most places).
// If you use float instead of double, you'll likely want to change this to something like 1.192092896E-04
#define EPSILON_GENERAL 1.192092896E-07

bool AreValuesEqual(double val1, double val2, double tolerance)
{
    if (val1 >= (val2 - tolerance) && val1 <= (val2 + tolerance))
    {
        return true;
    }

    return false;
}


double PointToPointDistanceSquared(double p1x, double p1y, double p2x, double p2y)
{
    double dx = p2x - p1x;
    double dy = p2y - p1y;
    return (dx * dx) + (dy * dy);
}


double PointSegmentDistanceSquared( double px, double py,
                                    double p1x, double p1y,
                                    double p2x, double p2y,
                                    double& t,
                                    double& qx, double& qy)
{
    double dx = p2x - p1x;
    double dy = p2y - p1y;
    double dp1x = px - p1x;
    double dp1y = py - p1y;
    const double segLenSquared = (dx * dx) + (dy * dy);
    if (AreValuesEqual(segLenSquared, 0.0, EPSILON_MIN_VERTEX_DISTANCE_SQUARED))
    {
        // segment is a point.
        qx = p1x;
        qy = p1y;
        t = 0.0;
        return ((dp1x * dp1x) + (dp1y * dp1y));
    }
    else
    {
        t = ((dp1x * dx) + (dp1y * dy)) / segLenSquared;
        if (t <= EPSILON_TINY)
        {
            // intersects at or to the "left" of first segment vertex (p1x, p1y).  If t is approximately 0.0, then
            // intersection is at p1.  If t is less than that, then there is no intersection (i.e. p is not within
            // the 'bounds' of the segment)
            if (t >= -EPSILON_TINY)
            {
                // intersects at 1st segment vertex
                t = 0.0;
            }
            // set our 'intersection' point to p1.
            qx = p1x;
            qy = p1y;
            // Note: If you wanted the ACTUAL intersection point of where the projected lines would intersect if
            // we were doing PointLineDistanceSquared, then qx would be (p1x + (t * dx)) and qy would be (p1y + (t * dy)).
        }
        else if (t >= (1.0 - EPSILON_TINY))
        {
            // intersects at or to the "right" of second segment vertex (p2x, p2y).  If t is approximately 1.0, then
            // intersection is at p2.  If t is greater than that, then there is no intersection (i.e. p is not within
            // the 'bounds' of the segment)
            if (t <= (1.0 + EPSILON_TINY))
            {
                // intersects at 2nd segment vertex
                t = 1.0;
            }
            qx = p2x;
            qy = p2y;
            // Note: If you wanted the ACTUAL intersection point of where the projected lines would intersect if
            // we were doing PointLineDistanceSquared, then qx would be (p1x + (t * dx)) and qy would be (p1y + (t * dy)).
        }
        else
        {
            // The projection of the point to the point on the segment that is perpendicular succeeded and the point
            // is 'within' the bounds of the segment.  Set the intersection point as that projected point.
            qx = ((1.0 - t) * p1x) + (t * p2x);
            qy = ((1.0 - t) * p1y) + (t * p2y);
            // for debugging
            //ASSERT(AreValuesEqual(qx, p1x + (t * dx), EPSILON_TINY));
            //ASSERT(AreValuesEqual(qy, p1y + (t * dy), EPSILON_TINY));
        }
        // return the squared distance from p to the intersection point.
        double dpqx = px - qx;
        double dpqy = py - qy;
        return ((dpqx * dpqx) + (dpqy * dpqy));
    }
}


double SegmentSegmentDistanceSquared(   double p1x, double p1y,
                                        double p2x, double p2y,
                                        double p3x, double p3y,
                                        double p4x, double p4y,
                                        double& qx, double& qy)
{
    // check to make sure both segments are long enough (i.e. verts are farther apart than minimum allowed vert distance).
    // If 1 or both segments are shorter than this min length, treat them as a single point.
    double segLen12Squared = PointToPointDistanceSquared(p1x, p1y, p2x, p2y);
    double segLen34Squared = PointToPointDistanceSquared(p3x, p3y, p4x, p4y);
    double t = 0.0;
    double minDist2 = 1E+38;
    if (segLen12Squared <= EPSILON_MIN_VERTEX_DISTANCE_SQUARED)
    {
        qx = p1x;
        qy = p1y;
        if (segLen34Squared <= EPSILON_MIN_VERTEX_DISTANCE_SQUARED)
        {
            // point to point
            minDist2 = PointToPointDistanceSquared(p1x, p1y, p3x, p3y);
        }
        else
        {
            // point - seg
            minDist2 = PointSegmentDistanceSquared(p1x, p1y, p3x, p3y, p4x, p4y);
        }
        return minDist2;
    }
    else if (segLen34Squared <= EPSILON_MIN_VERTEX_DISTANCE_SQUARED)
    {
        // seg - point
        minDist2 = PointSegmentDistanceSquared(p3x, p3y, p1x, p1y, p2x, p2y, t, qx, qy);
        return minDist2;
    }

    // if you have a point class and/or methods to do cross products, you can use those here.
    // This is what we're actually doing:
    // Point2D delta43(p4x - p3x, p4y - p3y);    // dir of p3 -> p4
    // Point2D delta12(p1x - p2x, p1y - p2y);    // dir of p2 -> p1
    // double d = delta12.Cross2D(delta43);
    double d = ((p4y - p3y) * (p1x - p2x)) - ((p1y - p2y) * (p4x - p3x));
    bool bParallel = AreValuesEqual(d, 0.0, EPSILON_GENERAL);

    if (!bParallel)
    {
        // segments are not parallel.  Check for intersection.
        // Point2D delta42(p4x - p2x, p4y - p2y);    // dir of p2 -> p4
        // t = 1.0 - (delta42.Cross2D(delta43) / d);
        t = 1.0 - ((((p4y - p3y) * (p4x - p2x)) - ((p4y - p2y) * (p4x - p3x))) / d);
        double seg12TEps = sqrt(EPSILON_MIN_VERTEX_DISTANCE_SQUARED / segLen12Squared);
        if (t >= -seg12TEps && t <= (1.0 + seg12TEps))
        {
            // inside [p1,p2].   Segments may intersect.
            // double s = 1.0 - (delta12.Cross2D(delta42) / d);
            double s = 1.0 - ((((p4y - p2y) * (p1x - p2x)) - ((p1y - p2y) * (p4x - p2x))) / d);
            double seg34TEps = sqrt(EPSILON_MIN_VERTEX_DISTANCE_SQUARED / segLen34Squared);
            if (s >= -seg34TEps && s <= (1.0 + seg34TEps))
            {
                // segments intersect!
                minDist2 = 0.0;
                qx = ((1.0 - t) * p1x) + (t * p2x);
                qy = ((1.0 - t) * p1y) + (t * p2y);
                // for debugging
                //double qsx = ((1.0 - s) * p3x) + (s * p4x);
                //double qsy = ((1.0 - s) * p3y) + (s * p4y);
                //ASSERT(AreValuesEqual(qx, qsx, EPSILON_MIN_VERTEX_DISTANCE_SQUARED));
                //ASSERT(AreValuesEqual(qy, qsy, EPSILON_MIN_VERTEX_DISTANCE_SQUARED));
                return minDist2;
            }
        }
    }

    // Segments do not intersect.   Find closest point and return dist.   No other way at this
    // point except to just brute-force check each segment end-point vs opposite segment.  The
    // minimum distance of those 4 tests is the closest point.
    double tmpQx, tmpQy, tmpD2;
    minDist2 = PointSegmentDistanceSquared(p3x, p3y, p1x, p1y, p2x, p2y, t, qx, qy);
    tmpD2 = PointSegmentDistanceSquared(p4x, p4y, p1x, p1y, p2x, p2y, t, tmpQx, tmpQy);
    if (tmpD2 < minDist2)
    {
        qx = tmpQx;
        qy = tmpQy;
        minDist2 = tmpD2;
    }
    tmpD2 = PointSegmentDistanceSquared(p1x, p1y, p3x, p3y, p4x, p4y, t, tmpQx, tmpQy);
    if (tmpD2 < minDist2)
    {
        qx = p1x;
        qy = p1y;
        minDist2 = tmpD2;
    }
    tmpD2 = PointSegmentDistanceSquared(p2x, p2y, p3x, p3y, p4x, p4y, t, tmpQx, tmpQy);
    if (tmpD2 < minDist2)
    {
        qx = p2x;
        qy = p2y;
        minDist2 = tmpD2;
    }

    return minDist2;
}
这段代码无法编译，因为minDist2 = PointSegmentDistanceSquared(p1x, p1y, p3x, p3y, p4x, p4y);这一行需要额外的参数。 - Richard
啊，对不起 - 我显然忘记了包含不需要t、qx和qy的方法版本。由于它们只是返回参数，您可以在那里添加虚拟值并忽略返回的结果。 - devnullicus
2
2
Afterlife提到了“首先使用此算法检查两条线是否相交”，但他没有指出他所指的算法。很明显，重点是线段的交点而不是延长线；任何非平行线段（不包括不定义一条线的重合端点）都会相交，但线段之间的距离不一定为零。因此我认为他指的是“线段”而不是“线”。
Afterlife给出的链接是找到一条线（或线段、射线）上离另一个任意点最近的点的非常优雅的方法。这适用于找到每个端点到另一条线段的距离（将计算参数u约束为线段或射线的最小值为0，最大值为1），但它不能处理一个线段上的内部点比任何一个端点更接近的情况，因此需要额外的交点检查。
至于确定线段交点的算法，一种方法是找到延长线的交点（如果平行则结束），然后确定该点是否在两个线段内，例如通过取从交点T到两个端点的向量的点积：
((Tx - A1x) * (Tx - A2x)) + ((Ty - A1y) * (Ty - A2y))
如果这是负数（或“零”），那么T在A1和A2之间（或在一个端点）。类似地检查另一条线段。如果任何一个大于“零”，则线段不相交。当然，这取决于首先找到延长线的交点，这可能需要将每条线表示为一个方程，并通过高斯消元（等）解决系统。
但是可能有一种更直接的方法，而无需解出交点，即通过将向量（B1-A1）和（B2-A1）的叉积以及向量（B1-A2）和（B2-A2）的叉积进行比较。如果这些叉积方向相同，则A1和A2在线B的同侧；如果它们方向相反，则它们在线B的两侧（如果为0，则其中一个或两个在线B上）。类似地检查向量（A1-B1）和（A2-B1）以及（A1-B2）和（A2-B2）的叉积。如果这些叉积中的任何一个是“零”，或者如果两个线段的端点落在另一条线的两侧，则线段本身必须相交，否则它们不相交。
当然，您需要一个方便的公式来计算两个向量的叉积。或者如果您能确定角度（为正或负），您就不需要实际的叉积，因为我们实际关心的是向量之间的角度方向（或者说是角度的正弦值）。但我认为叉积的公式（在2-D中）很简单：
Cross(V1,V2) = (V1x * V2y) - (V2x * V1y)
这是3-D叉积向量的z轴分量（因为初始向量在平面z=0上，所以x和y分量必须为零），因此您可以简单地查看符号（或“零”）。
因此，在Afterlife描述的算法中（参考链接），您可以使用其中一种方法来检查线段交点。
2

小贴士：如果你想基于点来比较距离，那么不必进行平方根运算。

例如，要查看 P 到 Q 的距离是否比 Q 到 R 更短，只需检查以下伪代码：

square(P.x-Q.x) + square(P.y-Q.y) < square(Q.x-R.x) + square(Q.y-R.y)
0

这个页面有一个很好的简短描述，可以找到两条线之间的最短距离，尽管@strager的链接包括一些Fortran代码！

这是3D的，它指的是线而不是线段。在这里不相关。 - Jason S
原始问题没有指定2D。 - Ian Hopkinson
好的。抱歉。建议那些给这个帖子点踩的人撤销他们的踩。（我没有点踩。） - Jason S
没问题 - 作为一个副作用，我发现了如何查看帖子的编辑历史！ - Ian Hopkinson
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接