計算幾何常用演算法及numpy模擬

在很久之前的一篇的文章點乘和叉乘及其物理意義（C++STL實現），我們用C++（STL）實現了對向量內積和叉積的定義與簡單計算，最後演示瞭如何用幾何的方法計算點到直線的距離，計算任意三角形的面積等問題。把這些放在更大的範圍內，其實就是今天的主角，計算幾何（Computational geometry）。

本文所有的模擬程式均在numpy（python第三方的科學計算庫）下進行：

import numpy as np

1. 向量減法

設二維向量 P=(x 1 ,y 1 ) ，Q=(x 2 ,y 2 )
則向量減法定義（對應位相減）為： P−Q=(x 1 −x 2 ,y 1

−y 2 )
顯然有性質 P−Q=−(Q−P)
如不加說明，下面所有的點都看作向量，兩點的減法就是向量相減；

def sub(p, q):
    return p-q
p, q = np.array([1, 2]), np.array([3, 5])
print(sub(p, q))
print(-sub(q, p))

2. 向量叉積

設向量P=(x 1 ,y 1 ),Q=(x 2 ,y 2 ) （本文以二維為例，可自然推廣到三維）
則向量叉積定義為： P×Q=x 1 ∗y 2 −x 2 ∗y 1 得到的是一個標量
顯然有性質 P×Q=−(Q×P),P×(−

Q)=−(P×Q)
如不加說明，下面所有的點都看作向量，點的乘法看作向量叉積；

叉乘的重要性質：

若 P × Q > 0 , 則P 在Q的順時針方向
若 P × Q < 0 , 則P 在Q的逆時針方向（與性質1等價）
若 P × Q = 0 , 則P 與Q共線，但可能同向也可能反向

def cross_prod(p, q):
    assert len(p)==2 and len(q)==2
    return p[0]*q[1]-p[1]*q[0]
p, q = np.array([1, 0]), np.array([0, 1])
                        # p在q的順時針方向 

print(cross_prod(p, q))
                        # 1>0
p, q = np.array([1, 2], [2, 4])         
                        # p, q 同向
print(cross_prod(p, q))
                        # 0

3. 判斷點線上段上

設點為Q ，線段為P 1 P 2   ，判斷點Q 在該線段上的依據是：
(Q−P 1 )×(P 2 −P 1 )=0 且 Q 在以 P 1 ,P 2   為對角頂點的矩形內；
這是叉積性質的直接應用，也即(Q−P 1 )×(P 2 −P 1 )=0 ，則 Q−P 1  與P 2 −P 1  共線；

def is_point_on_line(q, p1, p2):
    return True if cross_prod(q-p1, p2-p1)==0 else False
p1, p2 = np.array([1, 0]), np.array([0, 1])
q = np.array([1/2, 1/2])
print(is_point_on_line(q, p1, p2))
                    # True
q = np.array([1/2, 1/3])
print(is_point_on_line(q, p1, p2))
                    # False

4. 判斷兩線段是否相交

我們分兩步確定兩條線段是否相交：

快速排斥試驗

設以線段 P 1 P 2 為對角線的矩形為R ，設以線段 Q 1 Q 2 為對角線的矩形為T ，如果R 和T 不相交，顯然兩線段不會相交；

def rect_overlap(r1, r2):
    return not (((r1[1][0] < r2[0][0])|(r1[0][1] > r2[1][1]))
               |((r2[1][0] < r1[0][0])|(r2[0][1] > r1[1][1]))

跨立試驗

如果兩線段相交，則兩線段必然相互跨立對方，如圖1所示。在圖1中，P 1 P 2 跨立Q 1 Q 2 ，則向量 (P1−Q1) 和(P 2 −Q 1 ) 位於向量(Q 2 −Q 1 ) 的兩側，即(P1−Q1)×(Q2−Q1)∗(P2−Q1)×(Q2−Q1)<0 上式可改寫成 ( P1 - Q1 ) × ( Q2 - Q1 ) * ( Q2 - Q1 ) × ( P2 - Q1 ) > 0，
當 ( P1 - Q1 ) × ( Q2 - Q1 ) = 0 時，說明 ( P1 - Q1 ) 和 ( Q2 - Q1 )共線，但是因為已經通過快速排斥試驗，所以 P1 一定線上段 Q1Q2上；同理，( Q2 - Q1 ) ×(P2 - Q1 ) = 0 說明 P2 一定線上段 Q1Q2上。所以判斷P1P2跨立Q1Q2的依據是： ( P1 - Q1 ) × ( Q2 - Q1 ) * ( Q2 - Q1 ) × ( P2 - Q1 ) ≥ 0 同理判斷Q1Q2跨立P1P2的依據是： ( Q1 - P1 ) × ( P2 - P1 ) * ( P2 - P1 ) × ( Q2 - P1 ) ≥ 0至此已經完全解決判斷線段是否相交的問題。

def line_intersection(p1, p2, q1, q2):
    return False if not rect_overlap([p1, p2], [q1, q2]) else \             # 不相交直接返回
            (True if cross_prod(p1-q1, q2-q1)*cross_prod(q2-q1, p2-q1) >= 0 else False)

5. 判斷線段和直線是否相交

從標題即可看出，這裡的直線和線段是不相同的；
如果線段 P 1 P 2 和直線 Q 1 Q 2 相交，則 P 1 P 2 跨立 Q 1 Q 2 ，根據4的結論：

計算幾何常用演算法及numpy模擬

1. 向量減法

2. 向量叉積

3. 判斷點線上段上

4. 判斷兩線段是否相交

5. 判斷線段和直線是否相交

計算幾何常用演算法及numpy模擬

BZOJ 1100 POI2007 對稱軸osi 計算幾何+KMP演算法

計算幾何工具演算法-求任意多邊形的面積

關於計算幾何一些演算法

ACM常用演算法及練習

計算幾何-常用幾何函式（模板）

深入淺出PID控制演算法（一）————連續控制系統的PID演算法及MATLAB模擬

計算機幾何常用演算法

Codeforces Round #558 (Div. 2)C(計算幾何，排列組合，模擬)

計算幾何基礎演算法幾何C++實現

2017級演算法模擬上機準備篇（計算幾何初階_1)

計算幾何與圖形學有關的幾種常用演算法（二）

算法系列之九：計算幾何與圖形學有關的幾種常用演算法（二）

計算幾何與圖形學有關的幾種常用演算法（一）

計算幾何與圖形學有關的幾種常用演算法

計算幾何演算法

POJ-2420 A Star not a Tree? 計算幾何模擬退火

HDU-3644 A Chocolate Manufacturer's Problem 計算幾何模擬退火

POJ-1379 Run Away 計算幾何模擬退火

計算幾何平面最近點對 nlogn分治演算法求平面中距離最近的兩點

計算幾何常用演算法及numpy模擬

1. 向量減法

2. 向量叉積

3. 判斷點線上段上

4. 判斷兩線段是否相交

5. 判斷線段和直線是否相交

相關推薦