Numpy 实现透视变换

其实就是用numpy.linalg.solve解方程啦

透视变换就是三维空间上的一个点乘上变换矩阵转移到别的坐标（因为图片本身是二维的，所以新旧坐标的 z 值都为 1），即：

\begin{bmatrix} XW \\ YW \\ W \end{bmatrix} = \begin{bmatrix} a & b & c \\ d & e & f \\ g & h & 1 \end{bmatrix} \begin{bmatrix} x \\ y \\ 1 \end{bmatrix}

解得新坐标：

X = \frac {ax+by+c} W = \frac {ax+by+c} {gx+hy+1}
Y = \frac {dx+ey+f} W = \frac {dx+ey+f} {gx+hy+1}

但其实 a-h 是未知的，所以得先通过几对给定的新旧点来解出它们，才能得到新旧图片上所有点的映射。

解 a-h 的步骤如下：首先把分母与左边的 X/Y 相乘，再整理等式把单独的 X/Y 拿到等号右边可以得到两个方程：

ax+by+\phantom 1 c-0d+0e+0f-Xxg-Xyh=X
0a+0b+0c+xd+ye+\phantom 1 f-Yxg-Yyh=Y

转换成矩阵形式：

\begin{bmatrix} x & y & 1 & 0 & 0 & -1 & -Xx & -Xy \\ 0 & 0 & 0 & x & y & 1 & -Yx & -Yy \\ & & & & & \vdots \end{bmatrix} \begin{bmatrix} a \\ b \\ c \\ d \\ e \\ f \\ g \\ h \end{bmatrix} = \begin{bmatrix} X \\ Y \\ \vdots \end{bmatrix}

8 个未知数需要 8 个方程来解，每对新旧点提供两个方程，所以需要 4 对点来填满矩阵，求解 a-h 相当于解形如 Ax=b 的方程。 numpy.linalg.solve就是用来解 Ax=b 的函数，用 A 与 b 当参数调用就可以得到 x（a-h），这样我们就得到了变换矩阵。

构造 A、b 与解变换矩阵的代码如下：

A_list = []
b_list = []
for xy, XY in zip(old_points, new_points):
    A_list.append([xy[0], xy[1], 1, 0, 0, 0, -XY[0] * xy[0], -XY[0] * xy[1]])
    A_list.append([0, 0, 0, xy[0], xy[1], 1, -XY[1] * xy[0], -XY[1] * xy[1]])
    b_list.extend(XY)

a, b, c, d, e, f, g, h = np.linalg.solve(A_list, b_list)
H = np.array([[a, b, c], [d, e, f], [g, h, 1]])

拿到变换矩阵遍历新图的每一处 (X, Y) 坐标，再通过numpy.linalg.solve就可以算出变换前的 (x, y) 坐标：

\begin{bmatrix} a & b & c \\ d & e & f \\ g & h & 1 \end{bmatrix} \begin{bmatrix} x/W \\ y/W \\ 1/W \end{bmatrix} = \begin{bmatrix} X \\ Y \\ 1 \end{bmatrix}
xy1 = np.linalg.solve(H, [X, Y, 1])
x, y, _ = (1 / xy1[2] * xy1).astype(int)

然后把原图 (x, y) 处的灰度值拿来放到新图的 (X, Y) 处即可（是的，这里没考虑插值啥的）：

if 0 <= x < im.shape[0] and 0 <= y < im.shape[1]:
    newim[X, Y] = im[x, y]

效果图：

perspective-transform

完整代码在这里

参考了这篇文章