使用 Numpy 找到两个线性方程的交点

Question

目标是找到两个线性方程的交点。这两个线性方程是使用 Numpy

polyfit

函数导出的。

给定两个时间序列 (

xLeft

,

yLeft

) 和 (

xRight

,

yRight

)，使用

polyfit

计算适合每个时间序列的线性最小suqares，如下所示：

xLeft = [
    6168, 6169, 6170, 6171, 6172, 6173, 6174, 6175, 6176, 6177,
    6178, 6179, 6180, 6181, 6182, 6183, 6184, 6185, 6186, 6187
]
yLeft = [
    0.98288751, 1.3639959, 1.7550986, 2.1539073, 2.5580614,
    2.9651523, 3.3727503, 3.7784295, 4.1797948, 4.5745049,
    4.9602985, 5.3350167, 5.6966233, 6.0432272, 6.3730989,
    6.6846867, 6.9766307, 7.2477727, 7.4971657, 7.7240791
]

xRight = [
    6210, 6211, 6212, 6213, 6214, 6215, 6216, 6217, 6218, 6219,
    6220, 6221, 6222, 6223, 6224, 6225, 6226, 6227, 6228, 6229,
    6230, 6231, 6232, 6233, 6234, 6235, 6236, 6237, 6238, 6239,
    6240, 6241, 6242, 6243, 6244, 6245, 6246, 6247, 6248, 6249,
    6250, 6251, 6252, 6253, 6254, 6255, 6256, 6257, 6258, 6259,
    6260, 6261, 6262, 6263, 6264, 6265, 6266, 6267, 6268, 6269,
    6270, 6271, 6272, 6273, 6274, 6275, 6276, 6277, 6278, 6279,
    6280, 6281, 6282, 6283, 6284, 6285, 6286, 6287, 6288]
yRight=[
    7.8625913, 7.7713094, 7.6833806, 7.5997391, 7.5211883,
    7.4483986, 7.3819046, 7.3221073, 7.2692747, 7.223547,
    7.1849418, 7.1533613, 7.1286001, 7.1103559,  7.0982385,
    7.0917811, 7.0904517, 7.0936642, 7.100791, 7.1111741,
    7.124136, 7.1389918, 7.1550579, 7.1716633, 7.1881566,
    7.2039142, 7.218349, 7.2309117, 7.2410989, 7.248455,
    7.2525721, 7.2530937, 7.249711, 7.2421637, 7.2302341,
    7.213747, 7.1925621, 7.1665707, 7.1356878, 7.0998487,
    7.0590014, 7.0131001, 6.9621005, 6.9059525, 6.8445964,
    6.7779589, 6.7059474, 6.6284504, 6.5453324, 6.4564347,
    6.3615761, 6.2605534, 6.1531439, 6.0391097, 5.9182019,
    5.7901659, 5.6547484, 5.5117044, 5.360805, 5.2018456,
    5.034656, 4.8591075, 4.6751242, 4.4826899, 4.281858,
    4.0727611, 3.8556159, 3.6307325, 3.3985188, 3.1594861,
    2.9142516, 2.6635408, 2.4081881, 2.1491354, 1.8874279,
    1.6242117,1.3607255,1.0982931,0.83831298
]

left_line = np.polyfit(xleft, yleft, 1)
right_line = np.polyfit(xRight, yRight, 1)

在这种情况下，

polyfit

分别输出

m

的系数

b

和

y = mx + b

。

两个线性方程的交集可以计算如下：

x0 = -(left_line[1] - right_line[1]) / (left_line[0] - right_line[0])
y0 = x0 * left_line[0] + left_line[1]

但是，我想知道是否存在 Numpy 内置方法来计算最后两步？

Answer 1

不完全是内置方法，但您可以简化问题。假设我有给定

y = m1 * x + b1

和

y = m2 * x + b2

的台词。您可以轻松地找到差值的方程，它也是一条线：

y = (m1 - m2) * x + (b1 - b2)

请注意，如果两条原始线相交，这条线将在两条原始线的交点处有一个根。您可以使用

numpy.polynomial.Polynomial

类来执行以下操作：

>>> (np.polynomial.Polynomial(left_line[::-1]) - np.polynomial.Polynomial(right_line[::-1])).roots()
array([6192.0710885])

请注意，我必须交换系数的顺序，因为

Polynomial

期望从最小到最大，而

np.polyfit

返回相反的值。其实

np.polyfit

不推荐。相反，您可以使用 Polynomial

 类方法直接获取 

np.polynomial.Polynomial.fit 对象。您的代码将如下所示：

left_line = np.polynomial.Polynomial.fit(xLeft, yLeft, 1, domain=[-1, 1])
right_line = np.polynomial.Polynomial.fit(xRight, yRight, 1, domain=[-1, 1])
x0 = (left_line - right_line).roots()
y0 = left_line(x0)

domain

映射到窗口

[-1, 1]

。如果您不指定

domain

，则将使用 x 值的峰峰值。您不希望这样做，因为它将导致输入值的映射。相反，我们明确指定域

[-1, 1]

映射到同一窗口。另一种方法是使用默认域并设置例如

window=[xLeft.min(), xLeft.max()]

。这种方法的问题在于，它会为多项式创建不同的域，从而阻止运算

left_line - right_line

。

请参阅 https://numpy.org/doc/stable/reference/routines.polynomials.classes.html 了解更多信息。

Answer 2

您可以将其建模为线性系统并使用简单的线性代数：

def get_intersection(m1,b1,m2,b2):
    A = np.array([[-m1, 1], [-m2, 1]])
    b = np.array([[b1], [b2]])
    # you have to solve linear System AX = b where X = [x y]'
    X = np.linalg.pinv(A) @ b
    x, y = np.round(np.squeeze(X), 4)
    return x, y # returns point of intersection (x,y) with 4 decimal precision

m1,b1,m2,b2 = left_line(0), left_line(1), right_line(0), right_line(1)
print(get_intersection(m1,b1,m2,b2))

例如，对于线 y - x = 1 和 y + x = 1，我们期望交集为 (0,1)：

m1,b1,m2,b2 = 1, 1, -1, 1
print(get_intersection(m1,b1,m2,b2))

输出：

(0.0, 1.0)

符合预期。