单应性 Homography

相机把三维世界压成二维图像，深度丢了，一般无法反过来。但有一种特殊情形，反过来完全没有歧义：被看的点都落在同一个平面上。地面、桌面、墙、棋盘格、一页文档、AprilTag 的黑白方块，都是平面。同一个平面被两台相机拍到，两幅图像里对应点之间的关系就是一个 3×3 矩阵 \(H\)：\(\lambda\tilde{\mathbf x}'=H\tilde{\mathbf x}\)，这个矩阵叫单应性（homography）。

另一种同样"免费"的情形是相机只转动、不平移（原地环顾、云台转动、手机拍全景）。这时场景有多深都无所谓，两幅图像也只差一个单应性。这个 3×3 矩阵撑起了机器人视觉里一长串应用：

• 鸟瞰图（逆透视映射）：把前视相机的地面像素映成俯视的度量地图，用于车道线、AMR 导航和可通行区域（FIG 5）；
• 障碍物检测：地面上的像素服从 \(H\)，不服从的就是"高出地面"的东西（FIG 2）；
• 全景拼接：旋转相机拍的多张图合成一张（FIG 3）；
• 平面标记的位姿：AprilTag / ArUco 的四个角点 → \(H\) → 标记相对相机的位姿；文档与白板的矫正也是同一件事；
• 平面标定：10.4 节中每张棋盘格照片给出一个 \(H\)，Zhang 方法从多个 \(H\) 里解出内参 \(K\)；
• SLAM 初始化：场景接近平面或相机几乎纯旋转时，对极几何（11.1 节）退化，要改用 \(H\)。

本节先用四个手柄建立直觉（FIG 1），看 DLT 如何从四个点对解出矩阵；再从几何推出 \(H=K'(R+\mathbf t\mathbf n\T/d)K^{-1}\) 并让平面之外的点"露馅"（FIG 2）；然后看纯旋转（FIG 3）、噪声与归一化（FIG 4）；最后在一个合成的公路场景上做鸟瞰图（FIG 5）。

想象把一张印好的图平放在桌上，从斜上方拍照：矩形变成了梯形，近大远小，直线仍是直线，平行线却可能汇聚。左图是原图，右图是你"拍到"的样子。拖动右图的四个角，每个角对应原图的一个角。四个点对 \((\mathbf x_i\mapsto\mathbf x_i')\) 恰好够解一个 \(H\)：每个点对给两个方程，\(H\) 有 9 个元素，但整体乘一个常数不改变映射，只剩 8 个自由度，\(4\times2=8\)。页面上每一帧都在运行真正的 DLT 算法（见下面的数学部分），再用 \(H^{-1}\) 反向映射每个输出像素到原图并做双线性采样，得到右图。

注意 \(H\) 的列有几何意义：第 1、2 列是原图 \(x\) 轴、\(y\) 轴方向无穷远点的像，也就是消失点（10.2 节）；第 3 列是原点的像。当四边形变成平行四边形，\(H\) 的最后一行变成 \((0,0,1)\)，消失点跑到无穷远，映射退化为仿射。

FIG 1 里的 \(H\) 是"拟合"出来的。真实的相机对里，它由几何完全决定。设相机 1 坐标系中有平面 \(\mathbf n\T\mathbf X_1=d\)（\(\mathbf n\) 单位法向，\(d\) 是相机 1 光心到平面的距离），相机 2 与相机 1 的关系为 \(\mathbf X_2=R\mathbf X_1+\mathbf t\)。则（见下面的数学部分） $$H=K'\Big(R+\frac{\hl{t}{\mathbf t}\,\hl{n}{\mathbf n}\T}{\hl{d}{d}}\Big)K^{-1}.$$ 下面的实验里，相机 1、2 看着一块平面（可以倾斜），平面中心立着一根杆。你可以移动相机 2、倾斜平面、改变杆高，页面做三件事：用上面的公式由几何算出 \(H\)；用平面上 25 个点的投影，用 DLT 独立解出 \(H\)；然后用 \(H\) 把相机 1 的点搬到相机 2。平面上的点严丝合缝（误差 \(10^{-12}\) 像素量级）；杆顶不在平面上，被 \(H\) 搬到的位置是错的，这个误差就是视差。

几何解释在 3D 视图里：\(H\) 把相机 1 的射线"压"到平面上，再从相机 2 看这个平面点。杆顶的射线先落在平面的 \(\mathbf Q\)（紫色），\(H\) 预测相机 2 在 \(\mathbf Q\) 的位置看到杆顶，而真正的杆顶在别处。这就是"地面单应性残差"检测障碍物的原理：地面像素残差为零，高出地面的物体残差正比于高度和基线。

让 \(\mathbf t=\mathbf 0\)，公式里的 \(\mathbf t\mathbf n\T/d\) 整项消失：\(H=K'RK^{-1}\)，\(\mathbf n\) 和 \(d\) 都不出现，场景可以是任意的三维结构。直观地说，绕光心旋转时每条射线只是转了个方向，物体多远都沿同一条射线成像。下面的场景里有 60 个深度从 3 m 到 25 m 的点。滑块让相机绕竖直轴转 \(\theta\)，并（或不）向右平移 \(t_x\)。紫色圆圈是用 DLT 在所有点上拟合的最优单应性的预测；纯旋转时它与真实点重合，哪怕深度差了 8 倍；一旦平移，近处的点移动得更多（\(\approx f\,t_x/Z\)），任何一个 \(H\) 都解释不了。

反过来，由 \(H\) 还能读出旋转：\(R=K^{-1}HK\)（先把 \(H\) 乘上 \(\det^{-1/3}\) 使 \(\det=1\)），再用 \(\theta=\arccos\frac{\operatorname{tr}R-1}2\) 或 5.6 节的 \(\Log\)。全景拼接就是这样估计相邻图像的相对旋转的（Brown–Lowe，2007）。

真实数据有噪声，点数多于 4，还有错配。下面把图像 2 中的点加上高斯噪声 \(\sigma\)，分别用Hartley 归一化的 DLT（青色）和直接用像素坐标的 DLT（红色）估计 \(H\)，用整幅 640×480 图像上一个 9×7 测试网格上的均方根误差（灰色是真值）评价。有三个旋钮：点数、点占图像的比例（"分布"）、噪声。要看到的现象是：点铺满图像时两者差不多；点挤在一个小区域且离坐标原点很远时，未归一化的 DLT 条件数 \(\kappa(A)\) 暴涨，误差大一个数量级；而不论哪个，外推到没有点的区域都是危险的。"加外点"按钮把个别点偏移 50 像素量级，最小二乘立刻被拖偏；"RANSAC"用 4 点采样、重复拟合找出内点（11.5 节详细讲）。

用 10.2 节的齐次坐标 \(\tilde{\mathbf x}=(x,y,1)\T\)，单应性是一个可逆的 3×3 矩阵作用在 \(\mathbb P^2\) 上：

\(cH\)（\(c\neq0\)）给出同一个映射，所以 9 个元素只有 8 个自由度，每个点对提供 2 个方程，至少 4 个点对。分母 \(w=h_{31}x+h_{32}y+h_{33}\) 就是 FIG 1 里的"\(w\) 比值"，除以 \(w\) 这一步让映射成为非线性的。\(H\) 的最后一行为 \((0,0,1)\) 时是仿射（6 个自由度，3 个点对，平行线保持平行）；再要求左上 2×2 是旋转乘缩放时是相似变换（4 个自由度，2 个点对）。复合就是矩阵乘法：\(H_{13}=H_{23}H_{12}\)，反向是 \(H^{-1}\)。

\(\lambda\) 未知，所以不写 \(\tilde{\mathbf x}'=H\tilde{\mathbf x}\)，而写"两个向量平行"：\(\tilde{\mathbf x}'\times H\tilde{\mathbf x}=\mathbf 0\)（叉积 1.8 节）。记 \(\tilde{\mathbf x}'=(u',v',1)\T\)，\(H\) 的三行为 \(\mathbf h_1\T,\mathbf h_2\T,\mathbf h_3\T\)，展开：

第三行等于前两行的线性组合（\(u'\)×第一行 + \(v'\)×第二行 \(=-\)第三行），所以每个点对只有 2 个独立方程，这就是"8 = 4×2"。把 \(n\) 个点对的 \(A_i\) 叠成 \(A\in\R^{2n\times9}\)：\(A\mathbf h=\mathbf 0\)。恰好 4 个点（任意三点不共线）时 \(A\) 秩为 8，零空间（1.6 节）是一维的，\(\mathbf h\) 唯一（相差一个常数）；\(n>4\) 时方程无精确解，取

（1.12 节的 SVD。）若三点在两幅图里都共线，\(A\) 的秩降到 7，\(\sigma_8\approx0\)，零空间多出一维，\(H\) 不唯一；只在一幅图里共线时 \(A\) 仍满秩，但解出的 \(H\) 是奇异的（FIG 1 读数里的 \(\sigma_3(H)/\sigma_1(H)\)）。\(\|A\mathbf h\|\) 叫代数误差，它没有几何含义，只是方便；本节末尾会回到几何误差。

像素坐标在百的量级，\(A\) 的一行里同时出现 \(1\)、\(x\sim10^2\)、\(u'x\sim10^4\) 三种尺度，\(\kappa(A)\) 高达 \(10^6\) 以上，SVD 的最小奇异向量对噪声极其敏感（FIG 4 读数）。办法：对两幅图各做一次相似变换 \(T,T'\)——平移使质心在原点，缩放使到原点的平均距离是 \(\sqrt2\)——在归一化坐标里解出 \(\hat H\)，再换回来：

它不增加点数，也不改变理想数据的解（4 个无噪声点时两种 DLT 一样准），只改善条件数，代价是几行代码。Hartley（1997）为八点算法提出，同一思想适用于所有 DLT。

相机 2 的坐标 \(\mathbf X_2=R\mathbf X_1+\mathbf t\)（与 5.3 节一致：\(R,\mathbf t\) 是相机 1 → 相机 2 的相对位姿），平面 \(\mathbf n\T\mathbf X_1=d\)，\(\mathbf n\) 是相机 1 坐标系里的单位法向，\(d>0\)。平面上的点满足 \(\mathbf n\T\mathbf X_1/d=1\)，所以可以把 \(\mathbf t\) 乘上这个"1"：

用的是 \(\mathbf X_1=\lambda K^{-1}\tilde{\mathbf x}_1\)，\(\lambda\) 被"\(\sim\)"吞掉了。几点说明。

• 符号：这里平面写成 \(\mathbf n\T\mathbf X=d\)，所以是加号。Hartley–Zisserman 把平面写成 \(\mathbf n\T\mathbf X+d=0\)，得到 \(R-\mathbf t\mathbf n\T/d\)。FIG 2 对加号做了数值验证；换号会差几百像素。
• 只出现 \(\mathbf t/d\)：把整个场景和基线同时放大 2 倍，\(H\) 不变，这就是单目的尺度不确定。
• 欧氏单应性：\(K'^{-1}HK=R+\mathbf t\mathbf n\T/d\) 只含运动与平面。分解它可得 \((R,\mathbf t/d,\mathbf n)\)，一般有至多 4 组数学解，要求点在两个相机前方后通常剩 2 组，第三个视角或已知法向再消歧（Malis–Vargas 2007；OpenCV 的 decomposeHomographyMat）。
• 平面外的点：若 \(\mathbf n\T\mathbf X_1=d'\neq d\)，则 \(\mathbf X_2=(R+\mathbf t\mathbf n\T/d)\mathbf X_1+\mathbf t\,(1-d'/d)\)。多出来的 \(\mathbf t(1-d'/d)\) 就是视差：正比于基线 \(\mathbf t\) 和"离平面有多远"，方向沿着极线（指向极点，11.1 节）。\(\mathbf t=\mathbf 0\) 时为零。

纯旋转：\(d\to\infty\)（平面在无穷远）与 \(\mathbf t=\mathbf 0\) 效果相同，\(H_\infty=K'RK^{-1}\)。相机绕光心转 \(\theta\)，\(K\) 相同时，图像中心移动约 \(f\tan\theta\) 像素（小角时 \(\approx f\theta\)），与深度无关。地面：世界平面 \(Z=0\) 上的点 \(\mathbf X=(X,Y,0)\)，\(\mathbf x_c=R\mathbf X+\mathbf t=\mathbf r_1X+\mathbf r_2Y+\mathbf t\)，\(R\) 的第三列 \(\mathbf r_3\) 乘 0 消失，所以 \(H_g=K[\mathbf r_1\ \mathbf r_2\ \mathbf t]\)，这里 \((R,\mathbf t)\) 是世界 → 相机的外参 \({}^cT_w\)（10.3 节）。对高 \(h\)、俯仰角 \(\alpha\) 向前看的相机（世界 \(X\) 向右、\(Y\) 向前、\(Z\) 向上）：

\(Y\to\infty\) 得到地平线 \(v=c_y-f\tan\alpha\)。反过来，从 \(H\) 求相机位姿（平面目标的位姿，AprilTag 与 10.4 节标定的第一步）：\(K^{-1}H_g\propto[\mathbf r_1\ \mathbf r_2\ \mathbf t]\)，取 \(\lambda=1/\|K^{-1}\mathbf h_1\|\) 使 \(\mathbf r_1\) 为单位向量，\(\mathbf r_3=\mathbf r_1\times\mathbf r_2\)，再用 SVD 把 \([\mathbf r_1\ \mathbf r_2\ \mathbf r_3]\) 投影到最近的旋转矩阵（\(R=UV\T\)，1.12 节），选符号使 \(t_z>0\)（目标在相机前方）。FIG 5 的读数用 4 个点的 DLT 恢复出相机的高度和俯仰角，就是这一段。有噪声时更好的做法是 PnP（11.4 节）或迭代优化。

DLT 最小化的 \(\|A\mathbf h\|\) 是代数量。真正想要的是像素误差。转移误差 \(\sum_i\|\mathbf x_i'-\pi(H\tilde{\mathbf x}_i)\|^2\)（\(\pi\) 是除以 \(w\)）只算图像 2 的误差；对称转移误差再加上 \(\|\mathbf x_i-\pi(H^{-1}\tilde{\mathbf x}_i')\|^2\)；若两幅图的点都有噪声，最大似然是重投影误差：同时估计 \(H\) 和精确点 \(\hat{\mathbf x}_i,\hat{\mathbf x}_i'\)（Sampson 距离是它的一阶近似）。后两种是非线性最小二乘，用 6.6、6.9 节的高斯–牛顿 / LM，初值取 DLT 结果。OpenCV 的 findHomography 在 RANSAC 之后也会在内点上做 LM 精化。

前视相机在路上看到的车道线汇聚到消失点，这对"在地图上走多远、车道多宽"毫无帮助。如果假设地面是平面，则 \(H_g\) 把地面坐标 \((X,Y)\)（米）映到像素，它的逆 \(H_g^{-1}\) 把每个像素映回地面米制坐标。对鸟瞰图的每个输出像素 \((X,Y)\)，用 \(H_g\) 找到它在相机图像中的位置采样，就得到一张俯视的度量图：车道线平行、宽度恒为 3.5 m、距离刻度均匀。这就是反向映射（输出 → 输入），不会留下空洞。下面是一个合成的两车道公路（所有像素由射线与地面求交渲染），相机高 \(h\)、俯仰 \(\alpha\) 可调。\(H_g\) 有两种来源：公式 \(K[\mathbf r_1\ \mathbf r_2\ \mathbf t]\)（已知安装高度和俯仰），或者四个手柄：在图中点出地面已知矩形（车道线与两条横向标线的交点，3.5 m × 10 m）的四个角，用 DLT 求 \(H\)，这是没有标定数据时的做法。

但这只对地面上的东西成立。高出地面的物体（琥珀色箱子）违反了平面假设：射线穿过它的顶部，落在地面上更远的位置，于是箱子被沿着视线方向拖长。侧视图里的相似三角形给出顶端落点 \(d'=d\,\dfrac{h}{h-z_b}\)（\(d\) 是水平距离，\(z_b\) 是高度）。