点到多边形的切线

最新推荐文章于 2020-03-08 16:32:26 发布

SimpleFelix

最新推荐文章于 2020-03-08 16:32:26 发布

阅读量2k

点赞数 1

分类专栏：计算几何算法与实现文章标签：切线

本文链接：https://blog.csdn.net/u013279723/article/details/104239723

版权

计算几何算法与实现专栏收录该内容

10 篇文章 16 订阅

订阅专栏

前言

切线的构造属于计算几何中的基础问题。存在高效的算法可以计算凸集的切线。对于凹集，其切线与对应凸集的切线相同，所以也能高效计算。因此，我们这里仅仅讨论凸多边形的切线。

寻找切线的算法与寻找极值点的算法（凸多边形的极值点）类似。通过运用相似的二分搜索方法，我们能构造O(logn)复杂度的切线寻找算法。

定义

直线L与多边形S相切定义为：L与S接触，但不穿过S的边界。也就是L仅仅与S的一个点相交或者与一条边重合。

根据上面的定义，S中的所有点要么在L上，要么都在L的同一侧。从外部任意一点，存在两条切线与S相切。

不妨设 $\mathbf{L}=\{V_0, V_1, ..., V_{n-1}, V_n}\}$ ,以逆时针方向给出。 $\mathbf{e_i}$ 为第i条边， $\mathbf{ev}_i=V_{i+1}-V_i$ 为边向量，由Vi指向Vi+1。想象一下，站在Vi点面朝 $\mathbf{e_i}$ 的方向，如果P点在左手边，那么称P点在 $\mathbf{e_i}$ 左侧。

如上图所示。记P点为多边形S外部任意一点。对于任意一个顶点Vi，可以用下面的方法很方便地检查直线PVi是否是切线：

考虑Vi前后的两条边 $\mathbf{e_{i-1}}$ 和 $\mathbf{e_i}$ 。如果点P在这两条边的同一侧（同在左侧或者同在右侧），那么Vi-1和Vi+1一定在直线PVi的不同侧，PVi就一定不是切线。反之，若点P在 $\mathbf{e_{i-1}}$ 的左侧，而在 $\mathbf{e_i}$ 的右侧，或者相反，那么点P就是一条切线。如上图，点P在 $\mathbf{e_{i-1}}$ 右侧而在 $\mathbf{e_i}$ 左侧，那么点P是最右切线；同样地，点P在 $\mathbf{e_{k-1}}$ 左侧而在 $\mathbf{e_k}$ 右侧，那么点P是最左切线。

对于凸多边形，可以用上述方法求出切线。但对于凹多边形，需要先求凸包，然后再用上述方法求其凸包的切线，其凸包的切线必然是原多边形的切线。相对于原始多边形，凸包往往点数大大减少，因此求解起来很快。多边形凸包算法可以参考之前的博文Melkman Algorithm。

对于凹多边形还有一种方法是：使用暴力搜索，先求出所有的局部切线，然后从这些切线中找到最右侧的 $T_R$ 和最左侧的 $T_L$ ，就是我们要求的两条切线了。如下图：

暴力算法

寻找任意多边形切线的暴力算法很容易实现，通常当多边形点数n较小时，采用这一算法比较高效。

其C++代码如下：

// Assume that classes are already given for the objects:
//    Point with coordinates {float x, y;}
//===================================================================
 

// isLeft(): test if a point is Left|On|Right of an infinite line.
//    Input:  three points P0, P1, and P2
//    Return: >0 for P2 left of the line through P0 and P1
//            =0 for P2 on the line
//            <0 for P2 right of the line
//    See: Algorithm 1 on Area of Triangles
inline float
isLeft( Point P0, Point P1, Point P2 )
{
    return (P1.x - P0.x)*(P2.y - P0.y) - (P2.x - P0.x)*(P1.y - P0.y);
}


// tests for polygon vertex ordering relative to a fixed point P
#define above(P,Vi,Vj)  (isLeft(P,Vi,Vj) > 0)   // true if Vi is above Vj
#define below(P,Vi,Vj)  (isLeft(P,Vi,Vj) < 0)   // true if Vi is below Vj
//===================================================================
 

// tangent_PointPoly(): find any polygon's exterior tangents
//    Input:  P = a 2D point (exterior to the polygon)
//            n = number of polygon vertices
//            V = array of vertices for any 2D polygon with V[n]=V[0]
//    Output: *rtan = index of rightmost tangent point V[*rtan]
//            *ltan = index of leftmost tangent point V[*ltan]
void
tangent_PointPoly( Point P, int n, Point* V, int* rtan, int* ltan )
{
    float  eprev, enext;        // V[i] previous and next edge turn direction

    *rtan = *ltan = 0;          // initially assume V[0] = both tangents
    eprev = isLeft(V[0], V[1], P);
    for (int i=1; i<n; i++) {
        enext = isLeft(V[i], V[i+1], P);
        if ((eprev <= 0) && (enext > 0)) {
            if (!below(P, V[i], V[*rtan]))
                 *rtan = i;
        }
        else if ((eprev > 0) && (enext <= 0)) {
            if (!above(P, V[i], V[*ltan]))
                 *ltan = i;
        }
        eprev = enext;
    }
    return;
}
//===================================================================

二分搜索

对于凸多边形，存在速度更快的二分搜索算法。只需要在多边形的点之间定义一个次序，就可以使用类似于凸多边形的极值点中的相同的二分搜索算法。例如，如下图，我们定义 $\mathbf{e_i}$ 是增大的, 当Vi+1相对于点P在Vi之上。或者说点P在 $\mathbf{e_i}$ 右侧。

可以看到对于凸多边形而言，从最下面的切点开始，先单调增大，到最上面大的切点之后，再单调减小。因此与凸多边形的极值点算法相比，唯一的差别就是顺序关系判断函数发生了变化，其他部分基本完全一样。下面直接给出一个C++实现：

// tangent_PointPolyC(): fast binary search for tangents to a convex polygon
//    Input:  P = a 2D point (exterior to the polygon)
//            n = number of polygon vertices
//            V = array of vertices for a 2D convex polygon with V[n] = V[0]
//    Output: *rtan = index of rightmost tangent point V[*rtan]
//            *ltan = index of leftmost tangent point V[*ltan]
void
tangent_PointPolyC( Point P, int n, Point* V, int* rtan, int* ltan )
{
    *rtan = Rtangent_PointPolyC(P, n, V);
    *ltan = Ltangent_PointPolyC(P, n, V);
}


// Rtangent_PointPolyC(): binary search for convex polygon right tangent
//    Input:  P = a 2D point (exterior to the polygon)
//            n = number of polygon vertices
//            V = array of vertices for a 2D convex polygon with V[n] = V[0]
//    Return: index "i" of rightmost tangent point V[i]
int
Rtangent_PointPolyC( Point P, int n, Point* V )
{
    // use binary search for large convex polygons
    int     a, b, c;            // indices for edge chain endpoints
    int     upA, dnC;           // test for up direction of edges a and c

    // rightmost tangent = maximum for the isLeft() ordering
    // test if V[0] is a local maximum
    if (below(P, V[1], V[0]) && !above(P, V[n-1], V[0]))
        return 0;               // V[0] is the maximum tangent point

    for (a=0, b=n;;) {          // start chain = [0,n] with V[n]=V[0]
        c = (a + b) / 2;        // midpoint of [a,b], and 0<c<n
        dnC = below(P, V[c+1], V[c]);
        if (dnC && !above(P, V[c-1], V[c]))
            return c;          // V[c] is the maximum tangent point

        // no max yet, so continue with the binary search
        // pick one of the two subchains [a,c] or [c,b]
        upA = above(P, V[a+1], V[a]);
        if (upA) {                       // edge a points up
            if (dnC)                         // edge c points down
                 b = c;                           // select [a,c]
            else {                           // edge c points up
                 if (above(P, V[a], V[c]))     // V[a] above V[c]
                     b = c;                       // select [a,c]
                 else                          // V[a] below V[c]
                     a = c;                       // select [c,b]
            }
        }
        else {                           // edge a points down
            if (!dnC)                        // edge c points up
                 a = c;                           // select [c,b]
            else {                           // edge c points down
                 if (below(P, V[a], V[c]))     // V[a] below V[c]
                     b = c;                       // select [a,c]
                 else                          // V[a] above V[c]
                     a = c;                       // select [c,b]
            }
        }
    }
}
 


// Ltangent_PointPolyC(): binary search for convex polygon left tangent
//    Input:  P = a 2D point (exterior to the polygon)
//            n = number of polygon vertices
//            V = array of vertices for a 2D convex polygon with V[n]=V[0]
//    Return: index "i" of leftmost tangent point V[i]
int
Ltangent_PointPolyC( Point P, int n, Point* V )
{
    // use binary search for large convex polygons
    int     a, b, c;            // indices for edge chain endpoints
    int     dnA, dnC;           // test for down direction of edges a and c

    // leftmost tangent = minimum for the isLeft() ordering
    // test if V[0] is a local minimum
    if (above(P, V[n-1], V[0]) && !below(P, V[1], V[0]))
        return 0;               // V[0] is the minimum tangent point

    for (a=0, b=n;;) {          // start chain = [0,n] with V[n] = V[0]
        c = (a + b) / 2;        // midpoint of [a,b], and 0<c<n
        dnC = below(P, V[c+1], V[c]);
        if (above(P, V[c-1], V[c]) && !dnC)
            return c;          // V[c] is the minimum tangent point

        // no min yet, so continue with the binary search
        // pick one of the two subchains [a,c] or [c,b]
        dnA = below(P, V[a+1], V[a]);
        if (dnA) {                       // edge a points down
            if (!dnC)                        // edge c points up
                 b = c;                           // select [a,c]
            else {                           // edge c points down
                 if (below(P, V[a], V[c]))     // V[a] below V[c]
                     b = c;                       // select [a,c]
                 else                          // V[a] above V[c]
                     a = c;                       // select [c,b]
            }
        }
        else {                           // edge a points up
            if (dnC)                         // edge c points down
                 a = c;                           // select [c,b]
            else {                           // edge c points up
                 if (above(P, V[a], V[c]))     // V[a] above V[c]
                     b = c;                       // select [a,c]
                 else                          // V[a] below V[c]
                     a = c;                       // select [c,b]
            }
        }
    }
}
//===================================================================