计算几何问题汇总--点与线的位置关系

点与点之间， 线与线之间，点与线之间的位置关系是一类非常重要的问题。它不仅是平面几何学的基石，也常常应用于LBS(Location Based Service)，社交网络，以及数据库查询等领域。

本文中，我将给出判断这些关系的相关算法，作为参考。需要说明的是，我给出的这些问题的解法，都是建立在二维平面空间之上。有关多维空间的位置关系，大家可以仿照二维空间中问题的思路，做相应的拓展。

语言上，我用的当然还是Python.

点与点之间的距离

先从最简单的点与点的位置关系说起。一般情况下，我们只关心点与点之间的距离。

1. 点类的定义

为使算法思路更加清晰，先定义点类 Point，既然是在二维空间上，那么每个点都应该有两个属性：x, y分别代表点的横纵坐标。

class Point(object):
    """Point are two-dimension"""

    def __init__(self, x, y):
        self.x = x
        self.y = y

接下来就看看如何计算两点之间距离：当然可以用初中学的欧氏距离最基本的计算方法。但是考虑到代码编写的效率，以及方便以后向高维空间拓展。我在本文中将尽量使用向量计算。

而为了简化代码，我们使用对于向量运算已经相当成熟的库numpy

2. 两点之间距离的计算

显然，两点可以构成向量，而向量的长度则是其内积的开方。空间中，点$A$与点$B$的距离可以用向量$\overrightarrow{AB}$的模$|\overrightarrow{AB}|$表示。所以，现在需要做的，就是写一个函数，以两点为参数，计算由这两点构成的向量的模。

为了和本文之后的问题保持编码风格上一致，同时简化代码编写。我使用对向量运算已经极为成熟的库numpy帮助计算。并且定义了一个新的类 Vector ，类 Vector 以向量的起点和终点作为输入，生成一个只拥有属性x和y的向量对象。

最后，和前面定义的类放在一起，代码如下：

import numpy as np
# numpy help us do some vector calculation


class Point(object):
    """Point are two-dimension"""

    def __init__(self, x, y):
        self.x = x
        self.y = y

class Vector(object):
    """start and end are two points"""

    def __init__(self, start, end):
        self.x = end.x - start.x
        self.y = end.y - start.y

def pointDistance(p1, p2):
    """calculate the distance between point p1 and p2"""

    # v: a Vector object
    v = Vector(p1, p2)

    # translate v to a ndarray object
    t = np.array([v.x, v.y])

    # calculate the inner product of ndarray t
    return float(np.sqrt(t @ t))

说明一下，在Python3.5以后的版本中，使用numpy库时，ndarray对象之间的乘法可以用 @ ，代替之前的 v1.dot(v2) 这样的形式。

点与线之间的位置关系

1. 线的分类

点与线之间的位置关系就要稍微复杂一些了，复杂之处在于线分线段和直线两种情况。但是，在定义类的时候我都用两点来代表线段（直线）的两个属性。于是，至少代码看上去是没什么分别的。

不同之处在于，线段的两个点事两个端点，而直线的两个点是直线上任意两点。

class Segment(object):
    """the 2 points p1 and p2 are unordered"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2


class Line(object):
    """p1 and p2 are 2 points in straight line"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2

需要注意的是，这里并没有说线段的两个点是什么顺序（不一定说左边的点就是p1，右边就是p2）

2. 点与线的位置关系

(1) 计算点到直线的距离

如Fig.1(a)所示，现要求点$C$到直到直线$AB$的距离。还是向量法，据向量知识可知：

$$\begin{equation} cos \angle CAB = \frac{\overrightarrow{AC} \cdot \overrightarrow{AB}}{|\overrightarrow{AC}| \cdot |\overrightarrow{AB}|} \end{equation}$$

再由三角形知识可知，线段$AD$的长度为：

$$\begin{equation} |\overrightarrow{AC}| \cdot cos \angle CAB \end{equation}$$

所以，$\overrightarrow{AD}$可以这样计算：

$$\begin{equation} \overrightarrow{AD} = \frac{\overrightarrow{AB}}{|\overrightarrow{AB}|} \cdot |\overrightarrow{AD}| = \frac{\overrightarrow{AB}}{|\overrightarrow{AB}|} \cdot |\overrightarrow{AC}| \cdot cos \angle CAB = \frac{\overrightarrow{AB} \cdot \overrightarrow{AC}}{|\overrightarrow{AB}|^2} \overrightarrow{AB} \end{equation}$$

当$\overrightarrow{AD}$计算完成之后，可以根据$\overrightarrow{AD}$相应的坐标值得到点$D$的坐标，再由上面点和点之间的距离，即可得到线段$CD$的长度。

给出完整的代码如下：

import numpy as np
# numpy help us do some vector calculation


class Point(object):
    """Point are two-dimension"""

    def __init__(self, x, y):
        self.x = x
        self.y = y


class Segment(object):
    """the 2 points p1 and p2 are unordered"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2


class Line(object):
    """p1 and p2 are 2 points in straight line"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2


class Vector(object):
    """start and end are two points"""

    def __init__(self, start, end):
        self.x = end.x - start.x
        self.y = end.y - start.y


def pointDistance(p1, p2):
    """calculate the distance between point p1 and p2"""

    # v: a Vector object
    v = Vector(p1, p2)

    # translate v to a ndarray object
    t = np.array([v.x, v.y])

    # calculate the inner product of ndarray t
    return float(np.sqrt(t @ t))


def pointToLine(C, AB):
    """calculate the shortest distance between point C and straight line AB, return: a float value"""

    # two Vector object
    vector_AB = Vector(AB.p1, AB.p2)
    vector_AC = Vector(AB.p1, C)

    # two ndarray object
    tAB = np.array([vector_AB.x, vector_AB.y])
    tAC = np.array([vector_AC.x, vector_AC.y])

    # vector AD, type: ndarray
    tAD = ((tAB @ tAC) / (tAB @ tAB)) * tAB

    # get point D
    Dx, Dy = tAD[0] + AB.p1.x, tAD[1] + AB.p1.y
    D = Point(Dx, Dy)

    return pointDistance(D, C)

(2) 判断点是否在直线上
既然已经能够计算点到直线的距离了，那么，只需要看点到直线的距离是否为0即可知道这个点在不在直线上。

接着上面的代码，可以写出如下函数：

def pointInLine(C, AB):
    """determine whether a point is in a straight line"""

    return pointToLine(C, AB) < 1e-9

(3) 判断点是否在线段上

处理完了点和直线的位置关系，我们接着来看点与线段的位置关系。其实，最常用的就是一点：判断点是否在线段上。这和判断点是否在直线上最大的区别在于线段有起点、终点。

如Fig.1(b)所示，判断点$C$在不在线段$AB$上，可以这样解决：

1. 计算点$C$到线段$AB$所在直线的距离
2. 若这个距离为0，继续第3步；否则，返回False
3. 若点$C$的横坐标在点$A$与点$B$的横坐标之间，则返回True，否则返回False

函数如下：

def pointInSegment(C, AB):
    """determine whether a point is in a segment"""

    # if C in segment AB, it first in straight line AB
    if pointInLine(C, Line(AB.p1, AB.p2)):
        return min(AB.p1.x, AB.p2.x) <= C.x <= max(AB.p1.x, AB.p2.x)
    return False

还是需要结合上面的代码的，这里省略，只是写出函数

线与线之间的位置关系

1. 直线与直线

主要就是判断两条之间是否相交。很容易解决：先把两条直线做成向量，再判断两个向量是否平行即可。

如Fig.1(c)所示，判断直线$AB$与$CD$是否平行，可以通过向量平行的判别公式：

$$\begin{equation} \frac{B.y - A.y}{B.x - A.x} = \frac{D.y - C.y}{D.x - C.x} \end{equation}$$

也就是判断两个向量的斜率是否相等。若相等，则平行；若不等，则不平行。

def linesAreParallel(l1, l2):
    """determine whether 2 straight lines l1, l2 are parallel"""

    v1 = Vector(l1.p1, l1.p2)
    v2 = Vector(l2.p1, l2.p2)

    return abs((v1.y / v1.x) - (v2.y / v2.x)) < 1e-9

说明两点：

1. 这个函数既能判断直线是否平行，也能判断线段是否平行
2. 考虑到斜率在大多情况下是浮点数，所以不能用操作符 == 判断两个浮点数是否相等，而是通过 abs(k1 - k2) < 1e-9 判断。

2. 线段与线段

线段与线段是否平行通过上面的函数可以实现，那现在主要要解决的问题是如何判断两直线是否相交。这也是一个非常经典的算法。

一种比较标准的做法是通过“跨立实验”。基本逻辑是这样，倘若两条线段满足以下两个条件之一，则这两条直线一定相交：

1. 相互跨越对方线段所在的直线
2. 一条线段的端点在另一条直线上

上面的这两个条件包含了所有可能的线段相交的情况，如Fig.2所示。

基于此，可以先想想如何判断一条线段是“跨越”对方所在直线的。不过在此之前，先补充一个知识：向量叉积。

向量叉积是用来判断两个向量之间方向的

他的计算方法如下：

假设两个向量分别为$\vec a = (x_1, y_1), \vec b = (x_2, y_2)$，则$\vec a$与$\vec b$的叉积如下定义：

$$\begin{equation} \vec a \times \vec b = x_1 \cdot y_2 - x_2 \cdot y_1 \end{equation}$$

其中，符号$\times$用来表示向量叉积运算。而$\vec a$和$\vec b$都是以坐标原点为起点，以上面给出的坐标：$(x_1, y_1), (x_2, y_2)$为终点的两个向量。

根据这个叉积的算法，我们不难得到以下规律：

1. 若$\vec a \times \vec b > 0$，则$\vec a$在$\vec b$的顺时针方向
2. 若$\vec a \times \vec b < 0$，则$\vec a$在$\vec b$的逆时针方向
3. 若$\vec a \times \vec b = 0$，则$\vec a$与$\vec b$共线

其实向量叉积表示的是由这两个向量组成的平行四边形的有向面积。具体的推导我这里省略了。感兴趣的话，可以去学习一下行列式的几何意义。以后如果有机会，我可能也会写成博客。

写出计算叉积的函数：当然要结合上面定义的函数和类，这里只是给出函数

def crossProduct(v1, v2):
    """calculate the cross product of 2 vectors"""

    # v1, v2 are two Vector object
    return v1.x * v2.y - v1.y * v2.x

了解了叉积的性质以及实现，我们回过头再看线段的相互跨越问题。如图Fig.2，若线段$Q_1Q_2$跨越线段$P_1P_2$所在直线，则无论是Fig.2(a)或Fig.2(b)中的哪种情况，一定有下面的结论成立：

$$\begin{equation} \overrightarrow{Q_1P_1} \times \overrightarrow{Q_1Q_2} * \overrightarrow{Q_1P_2} \times \overrightarrow{Q_1Q_2} < 0 \end{equation}$$

更进一步，如果两条线段相互跨越，如Fig.3(a)所示的那样（Fig.3(b)的情况是线段$P_1P_2$跨越了线段$Q_1Q_2$所在的直线，但是反过来，线段$Q_1Q_2$却没有跨越线段$P_1P_2$所在的直线）。则相互跨越的两条线段一定同时满足下面的两个条件：

$$\begin{equation} \overrightarrow{Q_1P_1} \times \overrightarrow{Q_1Q_2} * \overrightarrow{Q_1P_2} \times \overrightarrow{Q_1Q_2} < 0 \end{equation}$$

$$\begin{equation} \overrightarrow{P_1Q_1} \times \overrightarrow{P_1P_2} * \overrightarrow{P_1Q_2} \times \overrightarrow{P_1P_2} < 0 \end{equation}$$

此外，别忘了上面红字标出的两线段相交的另一个条件，就是一条线段的一个端点在另一条线段上的情况。比如Fig.3(c)，此时，

$$\begin{equation} \overrightarrow{P_2Q_1} \times \overrightarrow{P_2P_1} = 0 \end{equation}$$

所以，还需要检测点$Q_1$是否在线段$P_1P_2$上。因为只要$Q_1$在$P_1P_2$所在的直线上，上面的式子就会成立。

综上所述，判断两条线段是否相交的算法思路就很清晰了：

先计算相应向量的叉积，一共是4个结果，根据上面的讲解，分为一下几种情况：

1. 4个结果分为两组；若两组组的结果都是异号的，则一定相交；
2. 4个结果中有0存在，则检查相应的点是否在另一条线段上，在，则相交；不在，则不相交
3. 以上两个条件都不成立，则肯定不相交

最后，我将上面讲的：计算两点的距离；计算点到直线的距离；判断点是否在直线上；判断点是否在线段上；判断两直线（线段）是否平行；判断两条线段是否相交等等计算几何中关于点与线的问题的完整代码写在下面，供大家参考：

点与线之间位置关系的完整代码

import numpy as np
# numpy help us do some vector calculation


class Point(object):
    """Point are two-dimension"""

    def __init__(self, x, y):
        self.x = x
        self.y = y


class Segment(object):
    """the 2 points p1 and p2 are unordered"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2


class Line(object):
    """p1 and p2 are 2 points in straight line"""

    def __init__(self, p1, p2):
        self.p1 = p1
        self.p2 = p2


class Vector(object):
    """start and end are two points"""

    def __init__(self, start, end):
        self.x = end.x - start.x
        self.y = end.y - start.y


def pointDistance(p1, p2):
    """calculate the distance between point p1 and p2"""

    # v: a Vector object
    v = Vector(p1, p2)

    # translate v to a ndarray object
    t = np.array([v.x, v.y])

    # calculate the inner product of ndarray t
    return float(np.sqrt(t @ t))


def pointToLine(C, AB):
    """calculate the shortest distance between point C and straight line AB, return: a float value"""

    # two Vector object
    vector_AB = Vector(AB.p1, AB.p2)
    vector_AC = Vector(AB.p1, C)

    # two ndarray object
    tAB = np.array([vector_AB.x, vector_AB.y])
    tAC = np.array([vector_AC.x, vector_AC.y])

    # vector AD, type: ndarray
    tAD = ((tAB @ tAC) / (tAB @ tAB)) * tAB

    # get point D
    Dx, Dy = tAD[0] + AB.p1.x, tAD[1] + AB.p1.y
    D = Point(Dx, Dy)

    return pointDistance(D, C)


def pointInLine(C, AB):
    """determine whether a point is in a straight line"""

    return pointToLine(C, AB) < 1e-9


def pointInSegment(C, AB):
    """determine whether a point is in a segment"""

    # if C in segment AB, it first in straight line AB
    if pointInLine(C, Line(AB.p1, AB.p2)):
        return min(AB.p1.x, AB.p2.x) <= C.x <= max(AB.p1.x, AB.p2.x)
    return False


def linesAreParallel(l1, l2):
    """determine whether 2 straight lines l1, l2 are parallel"""

    v1 = Vector(l1.p1, l1.p2)
    v2 = Vector(l2.p1, l2.p2)

    return abs((v1.y / v1.x) - (v2.y / v2.x)) < 1e-9


def crossProduct(v1, v2):
    """calculate the cross product of 2 vectors"""

    # v1, v2 are two Vector object
    return v1.x * v2.y - v1.y * v2.x


def segmentsIntersect(s1, s2):
    """determine whether 2 segments s1, s2 intersect with each other"""

    v1 = Vector(s1.p1, s1.p2)
    v2 = Vector(s2.p1, s2.p2)

    t1 = Vector(s1.p1, s2.p1)
    t2 = Vector(s1.p1, s2.p2)

    d1 = crossProduct(t1, v1)
    d2 = crossProduct(t2, v1)

    t3 = Vector(s2.p1, s1.p1)
    t4 = Vector(s2.p1, s1.p2)

    d3 = crossProduct(t3, v2)
    d4 = crossProduct(t4, v2)

    if d1 * d2 < 0 and d3 * d4 < 0:
        return True

    if d1 == 0:
        return pointInSegment(s2.p1, s1)
    elif d2 == 0:
        return pointInSegment(s2.p2, s1)
    elif d3 == 0:
        return pointInSegment(s1.p1, s2)
    elif d4 == 0:
        return pointInSegment(s1.p2, s2)

    return False 

在上面完成类和函数定义的基础上，给出一个测试脚本，方便检验：

if __name__ == "__main__":
    p1 = Point(0, 0)
    p2 = Point(2, 2)

    # 计算点p1, p2之间的距离
    print(pointDistance(p1, p2))  # >>> 2

    # 通过p1, p2分别建立一个线段和一个直线
    l1 = Line(p1, p2)
    s1 = Segment(p1, p2)

    # 设点p3，显然p3在l1上，却不在l2上
    p3 = Point(3, 3)

    print(pointInLine(p3, l1))  # >>> True
    print(pointInSegment(p3, s1))  # >>> False

    # 设点p4, p5得到一条与l1平行的直线l2
    p4 = Point(0, 1)
    p5 = Point(2, 3)

    l2 = Line(p4, p5)

    print(linesAreParallel(l1, l2))  # >>> True

    # 计算p4到l1的距离
    print(pointToLine(p4, l1))  # >>> 0.7071067...

    # 设两条线段s2, s3
    s2 = Segment(Point(0, 2), Point(5, -1))
    s3 = Segment(Point(1, 0.7), Point(5, -1))

    # s2与s1相交；s3与s1不相交
    print(segmentsIntersect(s2, s1))  # >>> True
    print(segmentsIntersect(s3, s1))  # >>> False