OpenCV 单应性矩阵（Homography）实战教程：原理精讲与五种演示代码逐行解析

单应性矩阵（Homography）是 OpenCV 计算机视觉中最核心的几何工具之一，它刻画了两个平面之间的投影变换关系，广泛用于相机位姿估计、透视校正、图像拼接等场景。本教程以 OpenCV 官方文档 doc/tutorials/features/homography/homography.markdown 为骨架，结合仓库内 samples/cpp/tutorial_code/features/Homography/ 下五个完整演示程序展开，帮助读者理解单应性的数学定义、findHomography / warpPerspective / decomposeHomographyMat 等关键 API 的调用方式，并掌握从平面点集恢复相机位姿、透视拉正、由相机位移构造单应性、单应性分解以及纯旋转全景拼接的完整实现。

适用版本：本仓库为 OpenCV 5.x 开发分支（文档声明兼容 OpenCV >= 3.0），演示代码可在 C++、Python、Java 三种语言中运行；文中所引用文件与命令行均以当前仓库实际内容为准。

---

1. 单应性矩阵的基本理论

1.1 什么是单应性矩阵？

简单地说，**平面单应性（planar homography）**刻画了两个平面之间、相差一个尺度因子（scale factor）的变换关系：

s · [x'  y'  1]ᵀ = H · [x  y  1]ᵀ

H = | h11  h12  h13 |
    | h21  h22  h23 |
    | h31  h32  h33 |

单应性矩阵是一个 3x3 矩阵，但它只有 8 个自由度（DoF），因为 H 是“尺度无关”的（乘以任意非零常数后，投影结果不变）。因此求解后通常要做归一化：约定 h33 = 1，或约定各元素平方和等于 1（h11² + h12² + ... + h33² = 1）。

1.2 单应性对应的三种典型“平面到平面”变换

并非只有真实平面结构才存在单应性。以下三种情况都可以用单应性描述：

1. 三维平面目标与其在图像平面上的投影：例如地面上一个标定板图案被相机拍摄，属于“平面表面 → 图像平面”的映射；
2. 同一个平面被两台不同位置的相机观察：属于“平面 → 平面”映射（经典双目/多视几何情形）；
3. 相机绕光心（投影中心）纯旋转拍摄任意场景：此时可以把场景中的所有点视为处于“无穷远平面”上，两幅图之间仍然存在单应性。

1.3 单应性有什么用？

文档给出三张典型应用场景图（存放于 doc/tutorials/features/homography/images/，即 homography_transformation_example1/2/3.jpg 与 homography_pose_estimation.jpg、homography_perspective_correction.jpg、homography_panorama_stitching.jpg）：

• 基于共面点的相机位姿估计（Camera pose estimation）：例如 AR 应用中通过平面 Marker 反算相机相对位姿；
• 透视消除/校正（Perspective removal / correction）：把倾斜拍摄的文档、棋盘等“拉正”成正视视角；
• 全景拼接（Panorama stitching）：把旋转拍摄的重叠图像对齐融合成宽幅图像。

---

2. 演示代码总览与运行方式

官方文档对应的演示程序位于仓库中的三个目录：

语言	位置
C++	samples/cpp/tutorial_code/features/Homography/（5 个 .cpp）
Python	samples/python/tutorial_code/features/Homography/
Java	samples/java/tutorial_code/features/Homography/

共五个演示，对应文件（C++ 版）：

Demo	主题	C++ 源文件
Demo 1	共面点相机位姿估计	pose_from_homography.cpp
Demo 2	透视校正	perspective_correction.cpp
Demo 3	由相机位移计算单应性	homography_from_camera_displacement.cpp
Demo 4	分解单应性矩阵	decompose_homography.cpp
Demo 5	纯旋转全景拼接	panorama_stitching_rotating_camera.cpp

演示使用的棋盘图像（left*.jpg）与相机内参文件（left_intrinsics.yml）都存放在 samples/data/ 目录；官方还提供了 C++ 之外的全套可编译示例工程框架，见 samples/CMakeLists.txt。运行任一 C++ 示例时均可用命令行传参控制输入，例如：

./pose_from_homography --image=left04.jpg --intrinsics=left_intrinsics.yml --width=9 --height=6 --square_size=0.025

下面逐一深入每个 Demo 的实现细节。

---

3. Demo 1：从共面点估计相机位姿

提示：本 Demo 旨在演示“如何从单应性快速恢复位姿”的原理，属于教学示例。若要在实际工程中估计平面或任意物体的相机位姿，应当使用 OpenCV 的 solvePnP。

3.1 整体思路

物体是平面时，“物体坐标系下的三维点 → 归一化相机坐标系下的图像点”之间的变换就是一个单应性（因为目标共面，Z=0 恒成立）。因此只要相机内参已知，就能从单应性矩阵中反推出相机位姿（旋转 + 平移）。本 Demo 用棋盘作为平面目标，借助 findChessboardCorners() 检测角点。

第一步：检测棋盘角点。棋盘尺寸 patternSize 在本例中为 9x6：

vector<Point2f> corners;
bool found = findChessboardCorners(img, patternSize, corners);

检测失败时程序直接返回；成功则用 drawChessboardCorners() 可视化标注，原文档配套结果图为 homography_pose_chessboard_corners.jpg。

3.2 构造平面上的三维物点

已知棋盘每格边长 squareSize（默认 0.025 m），物体坐标系下的三维角点可以轻易生成——注意只保留 XY、Z 恒为 0（源码见 pose_from_homography.cpp 的 calcChessboardCorners 与 [compute-chessboard-object-points] 段）：

for( int i = 0; i < boardSize.height; i++ )
    for( int j = 0; j < boardSize.width; j++ )
        corners.push_back(Point3f(float(j*squareSize),
                                  float(i*squareSize), 0));

由于平面目标 Z=0，估算单应性时要把 Z 坐标去掉、只留二维（见 [compute-object-points] 段）：把 Point3f 转成 Point2f(x, y)。

3.3 由像素角点得到归一化相机坐标下的像点

读取相机内参与畸变系数（用 FileStorage 从 left_intrinsics.yml 装载键 camera_matrix 与 distortion_coefficients），再用 undistortPoints() 反透视、去畸变，得到归一化相机坐标下的像点（源码见 [load-intrinsics] 与 [compute-image-points] 段）：

FileStorage fs( samples::findFile( intrinsicsPath ), FileStorage::READ);
Mat cameraMatrix, distCoeffs;
fs["camera_matrix"] >> cameraMatrix;
fs["distortion_coefficients"] >> distCoeffs;

vector<Point2f> imagePoints;
undistortPoints(corners, imagePoints, cameraMatrix, distCoeffs);

3.4 用 DLT 估算单应性

findHomography() 内部采用 Direct Linear Transform（DLT）等算法求解 H（源码见 [estimate-homography] 段）：

Mat H = findHomography(objectPointsPlanar, imagePoints);

3.5 从单应性快速恢复位姿

文档给出了从 H 恢复位姿的数学推导：设物点 X = (X, Y, 0, 1)，则

x = P·X = K·[r1  r2  r3  t]·(X,Y,0,1)ᵀ
        = K·[r1  r2  t]·(X,Y,1)ᵀ
        = H·(X,Y,1)ᵀ

从而 H = λ·K·[r1 r2 t]，即 K⁻¹·H = λ·[r1 r2 t]，完整位姿矩阵构造为：

P = K·[r1  r2  (r1 × r2)  t]

代码实现（源码 [pose-from-homography] 段）做了三步：先对第一列取范数归一化保证 ‖c1‖=1；取 c1、c2 为前两列，c3 = c1 × c2（即 r3）；第三列即为平移向量 tvec：

double norm = sqrt(H.at<double>(0,0)*H.at<double>(0,0) +
                   H.at<double>(1,0)*H.at<double>(1,0) +
                   H.at<double>(2,0)*H.at<double>(2,0));
H /= norm;
Mat c1  = H.col(0);
Mat c2  = H.col(1);
Mat c3 = c1.cross(c2);
Mat tvec = H.col(2);
Mat R(3, 3, CV_64F);
// 将 c1、c2、c3 填入 R 的列

文档特别指出：这是快速但不严谨的方案——它不能保证 R 是正交矩阵（带有反射/缩放误差），尺度也只是粗略地把第一列归一为 1。

3.6 用极分解获得“真正的”旋转矩阵

为了让 R 满足旋转矩阵的正交性，文档推荐对旋转部分做极分解 / 正交化，最实用的实现是利用 SVD（源码 [polar-decomposition-of-the-rotation-matrix] 段）：

SVDecomp(R, W, U, Vt);
R = U*Vt;
double det = determinant(R);
if (det < 0)   // 保证 det=+1（右手系旋转）
{
    Vt.at<double>(2,0) *= -1;
    Vt.at<double>(2,1) *= -1;
    Vt.at<double>(2,2) *= -1;
    R = U*Vt;
}

程序中会分别打印极分解前后的 R 与 det(R)（正交化后行列式应恰为 1）。最后用 Rodrigues(R, rvec) 转成旋转向量，用 drawFrameAxes() 把三维坐标系轴画回图像上验证结果（源码 [display-pose] 段），效果如 homography_pose.jpg 所示——物体坐标系的三色坐标轴被投影到图像中，且与棋盘贴合，说明位姿估计正确。

---

4. Demo 2：透视校正（Perspective Correction）

这个 Demo 演示如何计算把源图像映射到期望视角的单应性并做透视变换。左侧是源图像（斜视角拍摄的棋盘 left02.jpg），右侧是期望得到的正视角棋盘视图（left01.jpg），两者都是平面目标，因此存在单应关系（原图见 homography_source_desired_images.jpg）。

4.1 检测两组角点

对源图和目标图各做一次角点检测（C++ 源码 perspective_correction.cpp 的 [find-corners] 段）：

vector<Point2f> corners1, corners2;
bool found1 = findChessboardCorners(img1, patternSize, corners1);
bool found2 = findChessboardCorners(img2, patternSize, corners2);

4.2 估计单应性

直接把两组像素角点喂给 findHomography()，即得到 “源平面 → 目标平面” 的单应矩阵：

Mat H = findHomography(corners1, corners2);

4.3 用 warpPerspective 做透视拉正

单应矩阵 H 确定后，用 warpPerspective() 把源图像逐像素投影到目标视角：

Mat img1_warp;
warpPerspective(img1, img1_warp, H, img1.size());

在目标窗口里把“期望棋盘视图”与“拉正后的源棋盘视图”并排（hconcat）显示，肉眼应几乎看不出差异（结果见 homography_perspective_correction_chessboard_warp.jpg）。

4.4 验证：把源角点投影到目标图

为验证 H 的正确性，把源图每个角点用齐次坐标经 H 变换后反投影回目标图坐标，再与目标图对应角点连线（源码 [compute-transformed-corners] 段）。关键点：齐次坐标要除以第三个分量归一化：

for (size_t i = 0; i < corners1.size(); i++)
{
    Mat pt1 = Mat_<double>({3,1}, {corners1[i].x, corners1[i].y, 1});
    Mat pt2 = H * pt1;
    pt2 /= pt2.at<double>(2);   // 齐次坐标归一化
    // 在并排图上从源角点连到投影后的目标位置画线
    Point end( (int) (img1.cols + pt2.at<double>(0)), (int) pt2.at<double>(1) );
    line(img_draw_matches, corners1[i], end, randomColor(rng), 2);
}

若所有连线都精确地落在目标棋盘的同名角点上，说明单应性估算准确（匹配连线图见 homography_perspective_correction_chessboard_matches.jpg）。

---

5. Demo 3：由相机位移计算单应性

单应性既然描述两个平面视图之间的关系，那么只要知道两个相机位姿及目标平面的几何信息，就能直接算出两张视图之间的单应性（无需特征点匹配）。本节用公式推导 + 源码逐行讲解（详见 homography_from_camera_displacement.cpp，理论基础见 Malis & Vargas 关于单应性分解的论文）。

5.1 预备知识：相机位姿与齐次变换

solvePnP() 可以依据“三维物点 ↔ 二维像点”对应关系、内参与畸变系数解出相机位姿（即标定棋盘相对相机的位姿）：

s·(u, v, 1)ᵀ = K · [r11 r12 r13 t_x; r21 r22 r23 t_y; r31 r32 r33 t_z] · (Xo,Yo,Zo,1)ᵀ
             = K · ᶜMₒ · (Xo,Yo,Zo,1)ᵀ

其中 K 为内参矩阵，ᶜMₒ 表示“物体坐标系到相机坐标系”的刚体变换（4x4 齐次形式）；solvePnP 返回的 rvec（Rodrigues 旋转向量）与 tvec 即对应此位姿。把点从坐标系 A 变换到坐标系 B 只需做矩阵乘法：

ᶜ²M_c1 = ᶜ²Mₒ · (ᶜ¹Mₒ)⁻¹

利用 (ᶜMₒ)⁻¹ = [Rᵀ -Rᵀ·t; 0 1] 的性质即可展开计算。源码中 computeC2MC1() 用绝对位姿求两相机之间的相对旋转与平移（[compute-c2Mc1] 段）：

R_1to2 = R2 * R1.t();
tvec_1to2 = R2 * (-R1.t()*tvec1) + tvec2;

5.2 求解步骤：先算两个相机位姿

对两张棋盘图分别执行“检测角点 → 生成物点 → solvePnP”得到位姿 1 与位姿 2（源码 [compute-poses] 段，源文件 homography_from_camera_displacement.cpp）。drawFrameAxes 会将两个位姿的三维坐标系画在各自图像上（见 homography_camera_displacement_poses.jpg）。

5.3 平面法向量 n 与距离 d

如图（homography_camera_displacement.png）所示，n 是平面单位法向量、d 是相机坐标系原点沿法向到平面的距离。本案例棋盘平面即物体 XY 平面，其法向量在物体系为 (0,0,1)，转到相机 1 系只需乘以旋转矩阵 R1（源码 [compute-plane-normal-at-camera-pose-1] 段）：

Mat normal = Mat_<double>({3,1}, {0, 0, 1});
Mat normal1 = R1*normal;

平面到相机 1 坐标系原点的距离，可由法向量与平面上任意一点（这里取物体系原点经 R1、tvec1 变换到相机 1 系的点 origin1）做点积得到，代码中直接求其倒数 d_inv1（源码 [compute-plane-distance-to-the-camera-frame-1] 段）：

Mat origin(3, 1, CV_64F, Scalar(0));
Mat origin1 = R1*origin + tvec1;
double d_inv1 = 1.0 / normal1.dot(origin1);

5.4 欧氏单应性公式与“正负号”陷阱

文档给出由相机位移构造单应性的标准公式（平面到平面方程）：

²H₁ = ²R₁ − (²t₁ · ¹nᵀ) / ¹d

其中 ²R₁ 是两个相机坐标系之间的旋转、²t₁ 是平移。但是，本案例中棋盘坐标系的 Z 轴指向物体内部，而推导所用的平面图约定 Z 轴指向外部，两者方向相反，所以符号取反，实际实现为加号（源码 computeHomography 与 [compute-homography-from-camera-displacement] 段）：

Mat computeHomography(const Mat &R_1to2, const Mat &tvec_1to2, const double d_inv, const Mat &normal)
{
    Mat homography = R_1to2 + d_inv * tvec_1to2*normal.t();
    return homography;
}

Mat homography_euclidean = computeHomography(R_1to2, t_1to2, d_inv1, normal1);
Mat homography = cameraMatrix * homography_euclidean * cameraMatrix.inv();
homography /= homography.at<double>(2,2);

最后一行的 G = γ·K·H·K⁻¹ 是把欧氏单应性 H 转成投影单应性 G（需要内参矩阵 K，此处假设两幅视图由同一相机拍摄），并以 h33=1 归一化。

5.5 与 findHomography 估算结果对比

把“由相机位移推导出的单应性”与直接用角点对 findHomography() 估算的单应性打印对比，文档给出的参考输出高度接近：

findHomography H:
[0.32903393332201, -1.244138808862929, 536.4769088231476;
 0.6969763913334046, -0.08935909072571542, -80.34068504082403;
 0.00040511729592961, -0.001079740100565013, 0.9999999999999999]

homography from camera displacement:
[0.4160569997384721, -1.306889006892538, 553.7055461075881;
 0.7917584252773352, -0.06341244158456338, -108.2770029401219;
 0.0005926357240956578, -0.001020651672127799, 1]

两者数值相近（存在合理误差）。分别用两个 H 对图 1 做 warpPerspective 并排对比（见 homography_camera_displacement_compare.jpg），肉眼几乎无法区分两种方式得到的拉正图像——这从侧面印证了推导的正确性。

5.6 练习：位姿插值

文档为加深理解布置了一个练习：把一次大位移拆成 N 次中间单应性，逐次做 N 次 warping，观察中间视角的渐进变换过程，并对比最终帧与直接目标帧的“误差图”。参考输出如 homography_camera_poses_interpolation.jpg 所示：前三个子图是源图在三个内插相机视角下的 warp 结果，第 4 幅是最终视角 warp 结果与期望图像的误差图。

---

6. Demo 4：分解单应性矩阵

自 OpenCV 3 起，decomposeHomographyMat() 可以把单应性矩阵分解为若干组候选的旋转、平移与平面法向量（每组即一种几何解释）。官方文档说明其返回的可能解是：

• 平移向量 只能恢复到尺度（该尺度因子恰为平面距离 d，因为法向量是单位长度）；
• 若提供点对应关系，至少有两个解会因“正深度约束”（所有物点必须在相机前方）被剔除。

6.1 对“相机位移推导的单应性”做分解

先用上一节方法由两个 solvePnP 位姿构造投影单应性，再调用 decomposeHomographyMat（源码 decompose_homography.cpp 的 [decompose-homography-from-camera-displacement] 段）：

vector<Mat> Rs_decomp, ts_decomp, normals_decomp;
int solutions = decomposeHomographyMat(homography, cameraMatrix, Rs_decomp, ts_decomp, normals_decomp);

对每一组解，把旋转矩阵经 Rodrigues 转回旋转向量与“由位移算出的 rvec/tvec”逐项对照。文档打印的参考输出（节选 Solution 0）显示：四组解中恰有一组与真实相机位移几乎完全吻合：

Solution 0:
rvec from homography decomposition: [-0.0919829920641369, -0.5372581036567992, 1.310868863540717]
rvec from camera displacement:      [-0.09198299206413783, -0.5372581036567995, 1.310868863540717]
tvec from homography decomposition: [-0.7747961019053186, ...] and scaled by d: [-0.1578091561210742, ...]
tvec from camera displacement:      [0.1578091561210745, 0.005603443652993617, 0.1383378976078466]
plane normal from homography decomposition: [-0.1973513139420648, 0.6283451996579074, -0.7524857267431757]
plane normal at camera 1 pose:      [0.1973513139420654, -0.6283451996579068, 0.752485726743176]

注意：解出的 tvec 与真实位移方向相反、幅值差一个 1/d 因子，这正是“符号歧义 + 尺度歧义”的体现；而解出的法向量与相机 1 位姿下的真实法向量也互为相反方向。这正是正深度约束要剔除的解，而 Solution 1 则与真实值符号一致。对 findHomography() 直接估算的 H 做同样分解，同样存在一组吻合解（源码 [decompose-homography-estimated-by-findHomography] 段）。

6.2 如何选出唯一正确的解？

文档给出两个可操作的筛选思路：

1. 已知初始位姿 ᶜ¹Mₒ 时，用候选位移算 ᶜ²Mₒ = ᶜ²M_c1 · ᶜ¹Mₒ，把属于平面的三维物点投影回图像，验证是否都位于相机前方（正深度约束）；
2. 若知道平面在相机 1 位姿下的法向量，直接保留法向量最接近的那组解。

---

7. Demo 5：纯旋转相机的简易全景拼接

该示例仅为阐释“纯旋转（rotating camera）运动下任意场景可视为平面到平面映射”这一概念，不应替代正式的全景拼接方案——OpenCV 提供了完整的 stitching 模块（stitcher）用于生产级拼接。

7.1 为什么纯旋转下总有单应性？

单应性原则上只对平面结构成立。但当相机只绕光心旋转、没有平移时，可以把任意三维场景视为位于无穷远平面上的点，因此两幅图之间依然存在单应性。此时单应性可仅由旋转 R 与内参 K 直接构造：

s·(x', y', 1)ᵀ = K·R·K⁻¹·(x, y, 1)ᵀ

7.2 Blender 建模数据与代码

文档用免费开源三维软件 Blender 生成 Suzanne 模型的两次纯旋转相机视角（两视角原图对比见 homography_stitch_compare.jpg），并把相机位姿 c1Mo、c2Mo 与内参直接硬编码进示例（C++ 源码 panorama_stitching_rotating_camera.cpp，Python 版见 panorama_stitching_rotating_camera.py）。示例内参取 f = 700、主点 (320, 240)。

实现分四步：

① 提取旋转分量（[extract-rotation] 段）：

Mat R1 = c1Mo(Range(0,3), Range(0,3));
Mat R2 = c2Mo(Range(0,3), Range(0,3));

② 求相对旋转位移（[compute-rotation-displacement] 段，第二幅相对第一幅）：

// c1Mo * oMc2
Mat R_2to1 = R1*R2.t();

③ 构造投影单应性（[compute-homography] 段）：

Mat H = cameraMatrix * R_2to1 * cameraMatrix.inv();
H /= H.at<double>(2,2);

④ 拼接（[stitch] 段）：把第二幅图 warp 到第一幅的坐标系，画布加宽一倍，再把第一幅拷贝覆盖到左半部分即可：

Mat img_stitch;
warpPerspective(img2, img_stitch, H, Size(img2.cols*2, img2.rows));
Mat half = img_stitch(Rect(0, 0, img1.cols, img1.rows));
img1.copyTo(half);

拼接结果如 homography_stitch_Suzanne.jpg 所示，两视图在重叠区自然对齐，说明 H = K·R·K⁻¹ 的正确性。

---

8. 小结与延伸阅读

通过五个自底向上的 Demo，本教程完整覆盖了单应性的核心闭环：

1. 理论层面：3x3、8 自由度、尺度归一化、三种平面变换情形与三大应用方向；
2. 求解层面：findHomography()（DLT）估算、由位姿/平面几何推导、decomposeHomographyMat() 分解并用正深度约束选出唯一解；
3. 应用层面：位姿估计（solvePnP 是正式替代方案）、透视校正与 warpPerspective()、纯旋转全景拼接；
4. 源码细节：齐次坐标归一化、SVD 极分解恢复正交旋转、平面法向量/距离的点积与叉积求法、棋盘物点构造等工程技巧。

希望进一步深入可继续阅读本仓库的相关资源：

• 关于相机标定（内参、畸变、棋盘）的基础，参考 modules/calib 与棋盘数据 samples/data/left_intrinsics.yml；
• 生产级图像拼接的完整管线位于 modules/stitching；
• 位姿求解的正式 API solvePnP 可在 C++ 头文件 modules/calib/include/opencv2/calib.hpp 中查看；
• 全部演示源码与配套图像、参数，见 samples/cpp/tutorial_code/features/Homography/ 与 samples/data。

如果希望系统学习单应性的数学基础，官方文档建议参阅 Hartley & Zisserman 的《Multiple View Geometry in Computer Vision》、Szeliski 的《Computer Vision: Algorithms and Applications》等计算机视觉经典教材（本仓库根目录 README.md 亦列出了项目整体结构与构建方式）。