文章目录

透视变换简介
预备知识
透视变换公式推导
- 投影
- 二维坐标向齐次坐标的变换
- 公式求解
例子：A4纸视角校正
- 原始图片和坐标变换模板
- 计算程序
- 结果
透视变换的限制
前向映射与后向映射矩阵不互逆

透视变换简介

真实的世界是三维的，而图像是二维的（至少目前是），如果要以二维图像描述三维世界，并且看起来足够真实，那么三维世界向二维图像转化的过程就需要满足一定的几何投影关系，即透视关系，用非常简单的话来讲就是近处的物体在图像中显得比较大，远处的物体在图像中显得比较小。比如下图中的铁轨，越往远处看，图像中铁轨的间距越小，但实际上铁轨的间距是不变的。

透视变换就是对图像中的物体进行空间坐标变换时，使变换结果满足一定的透视关系。透视变换包含以下三个过程：

二维坐标向齐次坐标的变换
齐次坐标投影
图像插值

说到这里可能仍然非常抽象，继续向下看应该能使各位对透视变换的概念逐渐变得清晰一些。

预备知识

在继续向下了解透视变换的公式推导之前，建议先详细了解一下仿射变换和图像插值的相关知识，这部分知识在其他地方已经写过了，所以这里不会再写，不过请务必详细了解，否则再往下理解起来可能会有点费劲。
图像坐标空间变换：仿射变换（Affine Transformation）
图像插值：最邻近（nearest）与双线性（bilinear）

其中特别注意以下内容：

符号约定：(u,v)用来表示原始图像中的坐标，(x,y)用来表示变换后图像的坐标
齐次坐标的初步概念
仿射变换通式
前向映射与后向映射，并理解为了图像插值方便，空间坐标变换通常使用后向映射方式实现

透视变换公式推导

下面我们先讲投影的概念，然后再讲二维坐标向齐次坐标的变换（请注意，这个讲述顺序与计算顺序相反）

投影

在仿射变换（Affine Transformation）中曾经简单讲过，(x,y,1)(x,y,1)(x,y,1)就是(x,y)(x,y)(x,y)对应的齐次坐标。现在我们要扩展一下齐次坐标的概念。当我们使用相机进行拍摄时，可以认为三维空间中的物体都投影在了一个二维图像上面，如下图所示，红色的平面就可以认为是投影平面，我们将其定义为z=1z=1z=1的平面，那么(x,y,1)(x,y,1)(x,y,1)就是二维图像上的点(x,y)(x,y)(x,y)在三维空间中的坐标。

接下来，对(x,y,1)(x,y,1)(x,y,1)乘以一个系数，即α(x,y,1)\alpha(x,y,1)α(x,y,1)，改变α\alphaα就相当于在红色虚线上移动(x,y,1)(x,y,1)(x,y,1)点，这一点可以通过相似三角形的原理进行证明。红色虚线通过连接观察者（或相机）与点(x,y,1)(x,y,1)(x,y,1)得到。

我们将α(x,y,1)\alpha(x,y,1)α(x,y,1)记为(X,Y,Z)(X,Y,Z)(X,Y,Z)，即(X=αx,Y=αy,Z=α)(X=\alpha x, Y=\alpha y, Z=\alpha)(X=αx,Y=αy,Z=α)。投影就是将(X,Y,Z)(X,Y,Z)(X,Y,Z)投影至(x,y,1)(x,y,1)(x,y,1)的过程。若已知(X,Y,Z)(X,Y,Z)(X,Y,Z)，可以通过投影计算出(x,y)(x,y)(x,y)，其计算过程比较简单，根据α(x,y,1)=(X,Y,Z)\alpha(x,y,1)=(X,Y,Z)α(x,y,1)=(X,Y,Z)的定义，(X,Y)(X,Y)(X,Y)同除以ZZZ即可：
{x=X/Zy=Y/Z\begin{cases} x = X/Z \\[2ex] y = Y/Z \end{cases} ⎩⎨⎧x=X/Zy=Y/Z
上述定义和计算过程表明，所有由α(x,y,1)\alpha(x,y,1)α(x,y,1)定义的点有着一个共同的投影坐标，即(x,y,1)(x,y,1)(x,y,1)。

二维坐标向齐次坐标的变换

此处我们直接推导后向映射过程，后向映射是编程时真正使用的计算方法。
记原始图像中的坐标为(u,v)(u,v)(u,v)，对应的齐次坐标为(U,V,W)(U,V,W)(U,V,W)，变换后图像的坐标为(x,y)(x,y)(x,y)。
那么有：
[UVW]=[xy1][t11t12t13t21t22t23t31t32t33]\begin{bmatrix} U & V & W \end{bmatrix} = \begin{bmatrix} x & y & 1 \end{bmatrix} \begin{bmatrix} t_{11} & t_{12} & t_{13} \\ t_{21} & t_{22} & t_{23} \\ t_{31} & t_{32} & t_{33} \end{bmatrix} [UVW]=[xy1]⎣⎡t11t21t31t12t22t32t13t23t33⎦⎤
请注意，这个式子不好求解，因为如果将此式化为方程组形式，你会发现没有常系数项。接下来，我们要运用上面投影部分讲过的概念变换此式，使其容易求解。

根据投影的概念，当我们求出(U,V,W)(U,V,W)(U,V,W)后，我们接着通过(u=U/W,v=V/W)(u = U/W,v=V/W)(u=U/W,v=V/W)来求解(u,v)(u,v)(u,v)，请特别注意，在这个过程中，我们可以对(U,V,W)(U,V,W)(U,V,W)乘以任意一个常系数而不会影响计算结果，也就是说，我们可以对上面的矩阵方程两边同乘以一个系数。

那么我们两边同除以t33t_{33}t33，得到的结果除了将t33t_{33}t33改写为1以外，其他符号保持不变，所以可得：
[UVW]=[xy1][t11t12t13t21t22t23t31t321]\begin{bmatrix} U & V & W \end{bmatrix} = \begin{bmatrix} x & y & 1 \end{bmatrix} \begin{bmatrix} t_{11} & t_{12} & t_{13} \\ t_{21} & t_{22} & t_{23} \\ t_{31} & t_{32} & 1 \end{bmatrix} [UVW]=[xy1]⎣⎡t11t21t31t12t22t32t13t231⎦⎤
这个式子就容易求解了。
现在可以去跟仿射变换（Affine Transformation）对比一下，然后会发现，这里的变换矩阵多了两个参数：t13t_{13}t13和t23t_{23}t23。

公式求解

将上述矩阵方程式展开得：
{U=t11∗x+t21∗y+t31V=t12∗x+t22∗y+t32W=t13∗x+t23∗y+1\begin{cases} U = t_{11}*x + t_{21}*y + t_{31} \\ V = t_{12}*x + t_{22}*y + t_{32} \\ W = t_{13}*x + t_{23}*y + 1 \end{cases} ⎩⎪⎨⎪⎧U=t11∗x+t21∗y+t31V=t12∗x+t22∗y+t32W=t13∗x+t23∗y+1
由投影过程可得：
{u=UW=t11∗x+t21∗y+t31t13∗x+t23∗y+1v=VW=t12∗x+t22∗y+t32t13∗x+t23∗y+1\begin{cases} u = \frac {U}{W} = \frac {t_{11}*x + t_{21}*y + t_{31}} {t_{13}*x + t_{23}*y + 1} \\[2ex] v = \frac {V}{W} = \frac {t_{12}*x + t_{22}*y + t_{32}} {t_{13}*x + t_{23}*y + 1} \end{cases} ⎩⎨⎧u=WU=t13∗x+t23∗y+1t11∗x+t21∗y+t31v=WV=t13∗x+t23∗y+1t12∗x+t22∗y+t32
进一步可得：
{u(t13∗x+t23∗y+1)=t11∗x+t21∗y+t31v(t13∗x+t23∗y+1)=t12∗x+t22∗y+t32\begin{cases} u(t_{13}*x + t_{23}*y + 1) = t_{11}*x + t_{21}*y + t_{31} \\[2ex] v(t_{13}*x + t_{23}*y + 1) = t_{12}*x + t_{22}*y + t_{32} \end{cases} ⎩⎨⎧u(t13∗x+t23∗y+1)=t11∗x+t21∗y+t31v(t13∗x+t23∗y+1)=t12∗x+t22∗y+t32
整理得：
{x∗t11+y∗t21+1∗t31+0∗t12+0∗t22+0∗t32−ux∗t13−uy∗t23=u0∗t11+0∗t21+0∗t31+x∗t12+y∗t22+1∗t32−vx∗t13−vy∗t23=v\begin{cases} x*t_{11} + y*t_{21} + 1*t_{31} + 0*t_{12} + 0*t_{22} + 0*t_{32} - ux*t_{13} - uy*t_{23} = u \\[2ex] 0*t_{11} + 0*t_{21} + 0*t_{31} + x*t_{12} + y*t_{22} + 1*t_{32} - vx*t_{13} - vy*t_{23} = v \end{cases} ⎩⎨⎧x∗t11+y∗t21+1∗t31+0∗t12+0∗t22+0∗t32−ux∗t13−uy∗t23=u0∗t11+0∗t21+0∗t31+x∗t12+y∗t22+1∗t32−vx∗t13−vy∗t23=v
假设我们现在有n组对应的点：(u1,v1)(u_1,v_1)(u1,v1)对应(x1,y1)(x_1,y_1)(x1,y1)，(u2,v2)(u_2,v_2)(u2,v2)对应(x2,y2)(x_2,y_2)(x2,y2), …, (un,vn)(u_n,v_n)(un,vn)对应(xn,yn)(x_n,y_n)(xn,yn)，那么可以得到矩阵形式的线性方程组：
[x1y11000−u1x1−u1y1000x1y11−v1x1−v1y1x2y21000−u2x2−u2y2000x2y21−v2x2−v2y2...xnyn1000−unxn−unyn000xnyn1−vnxn−vnyn][t11t21t31t12t22t32t13t23]=[u1v1u2v2...unvn]\begin{bmatrix} x_1 & y_1 & 1 & 0 & 0 & 0 & -u_1x_1 & -u_1y_1 \\ 0 & 0 & 0 & x_1 & y_1 & 1 & -v_1x_1 & -v_1y_1 \\ x_2 & y_2 & 1 & 0 & 0 & 0 & -u_2x_2 & -u_2y_2 \\ 0 & 0 & 0 & x_2 & y_2 & 1 & -v_2x_2 & -v_2y_2 \\ ... \\ x_n & y_n & 1 & 0 & 0 & 0 & -u_nx_n & -u_ny_n \\ 0 & 0 & 0 & x_n & y_n & 1 & -v_nx_n & -v_ny_n \\ \end{bmatrix} \begin{bmatrix} t_{11}\\ t_{21} \\ t_{31} \\ t_{12} \\ t_{22} \\ t_{32} \\ t_{13} \\ t_{23} \end{bmatrix} = \begin{bmatrix} u_1 \\ v_1 \\ u_2 \\ v_2 \\ ... \\ u_n \\ v_n \end{bmatrix} ⎣⎢⎢⎢⎢⎢⎢⎢⎢⎡x10x20...xn0y10y20yn01010100x10x20xn0y10y20yn010101−u1x1−v1x1−u2x2−v2x2−unxn−vnxn−u1y1−v1y1−u2y2−v2y2−unyn−vnyn⎦⎥⎥⎥⎥⎥⎥⎥⎥⎤⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡t11t21t31t12t22t32t13t23⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎡u1v1u2v2...unvn⎦⎥⎥⎥⎥⎥⎥⎥⎥⎤
记：
A=[x1y11000−u1x1−u1y1000x1y11−v1x1−v1y1x2y21000−u2x2−u2y2000x2y21−v2x2−v2y2...xnyn1000−unxn−unyn000xnyn1−vnxn−vnyn]A = \begin{bmatrix} x_1 & y_1 & 1 & 0 & 0 & 0 & -u_1x_1 & -u_1y_1 \\ 0 & 0 & 0 & x_1 & y_1 & 1 & -v_1x_1 & -v_1y_1 \\ x_2 & y_2 & 1 & 0 & 0 & 0 & -u_2x_2 & -u_2y_2 \\ 0 & 0 & 0 & x_2 & y_2 & 1 & -v_2x_2 & -v_2y_2 \\ ... \\ x_n & y_n & 1 & 0 & 0 & 0 & -u_nx_n & -u_ny_n \\ 0 & 0 & 0 & x_n & y_n & 1 & -v_nx_n & -v_ny_n \\ \end{bmatrix} A=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎡x10x20...xn0y10y20yn01010100x10x20xn0y10y20yn010101−u1x1−v1x1−u2x2−v2x2−unxn−vnxn−u1y1−v1y1−u2y2−v2y2−unyn−vnyn⎦⎥⎥⎥⎥⎥⎥⎥⎥⎤

T=[t11t21t31t12t22t32t13t23]T = \begin{bmatrix} t_{11}\\ t_{21} \\ t_{31} \\ t_{12} \\ t_{22} \\ t_{32} \\ t_{13} \\ t_{23} \end{bmatrix} T=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎡t11t21t31t12t22t32t13t23⎦⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎤

B=[u1v1u2v2...unvn]B = \begin{bmatrix} u_1 \\ v_1 \\ u_2 \\ v_2 \\ ... \\ u_n \\ v_n \end{bmatrix} B=⎣⎢⎢⎢⎢⎢⎢⎢⎢⎡u1v1u2v2...unvn⎦⎥⎥⎥⎥⎥⎥⎥⎥⎤
那么可得最小二乘解为：
T=(ATA)−1ATBT = (A^TA)^{-1}A^TB T=(ATA)−1ATB
最小二乘解的推导过程是：
AT=BATAT=ATB(ATA)−1ATAT=(ATA)−1ATBT=(ATA)−1ATBAT = B \\[2ex] A^TAT = A^TB \\[2ex] (A^TA)^{-1}A^TAT = (A^TA)^{-1}A^TB \\[2ex] T = (A^TA)^{-1}A^TB AT=BATAT=ATB(ATA)−1ATAT=(ATA)−1ATBT=(ATA)−1ATB
注意：由于上述方线性程组中有8个未知数，所以最少需要4组对应的点才能得到方程组的解，如果刚好有4组点，那么得到的解就是精确解，如果多于4组点，那么得到的解就称为最小二乘解。

例子：A4纸视角校正

原始图片和坐标变换模板

有的时候我们拍摄角度不太精准的情况下，A4纸在拍摄的图像中不是呈现出标准的长方形，而是不规则的四边形，这时我们可以通过透视变换进行视角校正。

比如下图中有一张A4纸，但是很明显它不是长方形，现在我们就用透视变换来校正它，步骤如下：

拾取A4纸四个顶点的坐标，从左上角开始，顺时针依次是(1208, 1192), (2814, 1116), (3557, 1866), (1160, 2089)
定义校正后的模板，这个模板通过A4纸的标准尺寸得来，其标准尺寸是210mm x 297mm。我这里以标准尺寸为基础定义了两个模板：模板1仅包含A4纸本身；模板2除了A4纸外，还包含周围的一些物体。
- 模板1：此处我将210和297两个数字乘以5，得到像素坐标模板，从左上角开始，顺时针依次为(0, 0), (1485, 0), (1485, 1050), (0, 1050)。
- 模板2：在模板1的基础上，左上角偏移(1000, 1000)，然后右下角边界扩展(500, 500)，得到的坐标为：(1000, 1000), (2485, 1000), (2485, 2050), (1000, 2050)
然后用后面的程序来做计算

下图在程序中的名字是side.jpg

计算程序

# -*- coding: utf-8 -*-import cv2
import numpy as npdef nearest_sampler(image, coords):"""Nearest Neighbour samplerParameters----------image: ndarraysource image. shape is [height, width, channels]coords: ndarraycoordinates to be interpolated, the length of last axis should be 2,meaning 2D coordinateReturns-------output: ndarraythe interpolated image, same shape as coords except the last axis"""height, width, channels = image.shape[0:3]output_shape = list(coords.shape)coords = np.reshape(coords, (-1, output_shape[-1]))output_shape[-1] = channelscoords = np.round(coords).astype(np.int32)idx = (coords[:, 0] >= 0) & (coords[:, 0] < width) & \(coords[:, 1] >= 0) & (coords[:, 1] < height)output = np.zeros((coords.shape[0], channels), dtype=np.uint8)output[idx] = image[coords[idx, 1], coords[idx, 0], :]output = np.reshape(output, output_shape)return outputdef perspective_transform_matrix(uv, xy):"""Compute perspective transform matrixParameters----------uv: ndarraycoordinates of feature points in original image. shape is [n, 2], n isthe number of pointsxy: ndarraycoordinates of feature points in perspective transformed image. shapeis same as uvReturns-------transform_matrix: ndarraytransform matrix. shape is [3, 3]"""A = np.zeros((2 * xy.shape[0], 8))B = np.zeros((2 * xy.shape[0], 1))for i in range(xy.shape[0]):A[2 * i, 0] = xy[i, 0]A[2 * i, 1] = xy[i, 1]A[2 * i, 2] = 1.0A[2 * i, 3] = 0.0A[2 * i, 4] = 0.0A[2 * i, 5] = 0.0A[2 * i, 6] = -uv[i, 0] * xy[i, 0]A[2 * i, 7] = -uv[i, 0] * xy[i, 1]A[2 * i + 1, 0] = 0.0A[2 * i + 1, 1] = 0.0A[2 * i + 1, 2] = 0.0A[2 * i + 1, 3] = xy[i, 0]A[2 * i + 1, 4] = xy[i, 1]A[2 * i + 1, 5] = 1.0A[2 * i + 1, 6] = -uv[i, 1] * xy[i, 0]A[2 * i + 1, 7] = -uv[i, 1] * xy[i, 1]B[2 * i] = uv[i, 0]B[2 * i + 1] = uv[i, 1]T = np.linalg.inv(A.T @ A) @ A.T @ Btransform_matrix = np.append(T, [1.0]).reshape([3, 3]).Treturn transform_matrixdef perspective_coordinates(transform_matrix, height, width):"""Compute perspective coordinates acoording to the transform matrix and thetransformed image shapeParameters----------transform_matrix: ndarrayperspective transform matrix. shape is [3, 3]height: intheight of transformed imageheight: intwidth of transformed imageReturns-------coords: ndarrayperspective coordinates. shape is [height, width, 2]"""coords = np.meshgrid(np.arange(0, width), np.arange(0, height))# shape = [height, width, 2] after transposecoords = np.array(coords).transpose([1, 2, 0])ones = np.ones([height, width, 1])# homogeneous  coordinatescoords = np.concatenate((coords, ones), axis=2)# transformed coordinatescoords = coords @ transform_matrix# projection of coordinatescoords[:, :, 0:2] /= coords[:, :, 2:]return coords[:, :, 0:2]if __name__ == '__main__':# feature pointsuv = np.array([[1208, 1192],[2814, 1116],[3557, 1866],[1160, 2089]])# model pointsxy1 = np.array([[0, 0],[1485, 0],[1485, 1050],[0, 1050]])xy2 = np.array([[1000, 1000],[2485, 1000],[2485, 2050],[1000, 2050]])transform_matrix1 = perspective_transform_matrix(uv, xy1)coords1 = perspective_coordinates(transform_matrix1, 1050, 1485)transform_matrix2 = perspective_transform_matrix(uv, xy2)coords2 = perspective_coordinates(transform_matrix2, 2550, 2985)# samplerimage = cv2.imread('side.jpg')img1 = nearest_sampler(image, coords1)img2 = nearest_sampler(image, coords2)cv2.namedWindow('projection1', 0)cv2.namedWindow('projection2', 0)cv2.imshow('projection1', img1)cv2.imshow('projection2', img2)cv2.waitKey(0)cv2.destroyAllWindows()

结果

运行上述程得到的结果是：

projection1：

projection2

透视变换的限制

从上述例子的projection2可以看出，除了A4纸本身以外，其他地方产生了比较奇怪的形变，这种形变不符合真实视角下的情况。下图是我从比较正的角度拍摄的图，各位可以对比看看。（只对比几何形变，不要在意光照变化）

这也正是透视变换的限制。这种只由一次矩阵变换+投影所做的空间坐标变换，只能对单一平面有准确的物理意义，而其他的平面则会产生非真实的形变。在实际应用时要特别注意这一特性，不是任何情况都可以随意使用透视变换，不然可能带来难以预料的结果。

从这个角度上讲，我们处理视频时通过透视变换所做的多帧对齐，严格来讲大多存在形变误差（除非镜头中只有单一平面，比如一堵墙）。因为我们通常会在图像中广泛地检测特征点，这些特征点往往分布在不同的平面中，所以这种形变误差可以认为是一种平均效果。当帧间变化较小时，可以称之为误差；当帧间变化较大时，就不能称之为误差而是错误了，此时大概率不能由透视变换来做好对齐工作。

前向映射与后向映射矩阵不互逆

参考仿射变换（Affine Transformation）的最后一部分，透视变换同样有此性质。

因为本文档中的例子刚好是4组特征点，得到的变换矩阵是精确解，所以不适合在此处作为例子呈现，因此此处不再举额外的例子进行说明。

图像坐标空间变换：透视变换（Perspective Transformation），或称为单应性（Homography）变换相关推荐

透视变换单应性矩阵怎么求 matlab,单应性(homography)变换的推导
矩阵的一个重要作用是将空间中的点变换到另一个空间中.这个作用在国内的<线性代数>教学中基本没有介绍.要能形像地理解这一作用,比较直观的方法就是图像变换,图像变换的方法很多,单应性变换是其中 ...
Opencv学习笔记透视变换/单应性变换
1.透视变换概述透视变换(Perspective Transformation),又称Homography Transformation. 在计算机视觉领域,空间中同一平面的任意两幅图像通过单应性关 ...
OpenCV中的「透视变换 / 投影变换 / 单应性」—cv.warpPerspective、cv.findHomography
文章目录引言透视变换(projective transform) 单应性(Homography) opencv代码仿射变换相关函数投影变换相关的函数鸟瞰图代码示例小结引言图像的几何变换 ...
单应性（Homography）变换
我们已经得到了像素坐标系和世界坐标系下的坐标映射关系: 其中,u.v表示像素坐标系中的坐标,s表示尺度因子,fx.fy.u0.v0.γ(由于制造误差产生的两个坐标轴偏斜参数,通常很小)表示5个相机内参 ...
透视变换（perspective transformation）和射影（投影）变换（projective transformation）
透视变换(perspective transformation)和射影(投影)变换(projective transformation) 0几点说明投影和变换的区别. 投影分为:平行投影或者正射投影 ...
【OpenCV】透视变换 Perspective Transformation（续）
透视变换的原理和矩阵求解请参见前一篇<透视变换 Perspective Transformation>.在OpenCV中也实现了透视变换的公式求解和变换函数. 求解变换公式的函数: [cp ...
【OpenCV】透视变换 Perspective Transformation
透视变换的原理和矩阵求解请参见前一篇<透视变换 Perspective Transformation>.在OpenCV中也实现了透视变换的公式求解和变换函数. 求解变换公式的函数: Mat ...
【图像处理】透视变换 Perspective Transformation
透视变换(Perspective Transformation)是将图片投影到一个新的视平面(Viewing Plane),也称作投影映射(Projective Mapping).通用的变换公式为: ...
【CG】透视变换(Perspective Transformation)
起源透视与消失点单点透视 - 1 个消失点两点透视 - 2个消失点三点透视 - 3个消失点透视投影变换透视变换矩阵透视变换矩阵与 Homography.相机内参的关系透视与车载环境的联 ...

图像坐标空间变换：透视变换（Perspective Transformation），或称为单应性（Homography）变换