第10章：特征检测实战：Harris 与 Canny

从"什么是特征"出发，拆开 Harris 角点检测与 Canny 边缘检测两大经典算法的原理，再用 cv::cornerHarris 与 cv::Canny 亲手跑通、调参、把特征画回原图

10.1 特征是什么：角点、边缘与"可重复性"

先问一个艺术问题：人眼扫过一张照片，第一眼抓住的是什么？不是每个像素的灰度，而是轮廓——楼房的边、天空和屋顶的分界线、窗户的四角。画家起稿时画的正是这些：几条线加上几个交点。计算机视觉复刻的就是这套直觉，它管这些"线条与交点"叫特征（feature）。官方文档把这一族算法收在 imgproc 模块的 Feature Detection（特征检测）分组里，其中最基础、最经典的两类就是本章主角：角点（corner，两条边交界、向任意方向移动窗口灰度都剧烈变化的点）与边缘（edge，灰度突变的边界线）。与之相对的是平坦区：窗口往哪挪，灰度都几乎不变。

为什么要在意这些点？因为工程上特征承担着"锚点"的角色：全景拼接要找出两幅图里同一批角点才能对齐，视频跟踪要在帧与帧之间找到同一个标记点，相机标定需要能稳定定位的棋盘格角点。所以特征有一个硬性要求——可重复性：同一场景换个角度、换个光照再拍，同一个角点还应该被检出来。这也是 Harris 与 Canny 从 1988 年和 1986 年用到现在依然是教学与工程地基的原因。动手前先做一件公共准备：第 9 章讲过 Mat 是"数字地图"，而 Harris 与 Canny 的输入都要求单通道灰度图，所以第一步永远是读图、转灰度、落盘——下面的 gray_prep.cpp 是本章所有程序的公共前奏。

// gray_prep.cpp —— 特征检测的公共前奏：读图、转灰度、落盘
#include <opencv2/opencv.hpp>

int main() {
    cv::Mat img = cv::imread("building.jpg", cv::IMREAD_COLOR);
    if (img.empty()) {
        return 1;                  // 读图失败，沿用 9.6 节的防护习惯
    }
    // 转灰度：Harris 与 Canny 都以单通道灰度图为输入
    cv::Mat gray;
    cv::cvtColor(img, gray, cv::COLOR_BGR2GRAY);
    cv::imwrite("gray.png", gray);
    return 0;
}

"平坦区、边缘、角点"的区别能不能量化？能——用一个小窗口里的灰度方差就能看出来：平坦区方差小、边缘中等、角点最大。下面这个程序用 cv::meanStdDev 对三个假设位置各取一个窗口，把三种区域的"数字手感"打出来——"窗口移动、灰度变化"的直觉，正是 10.2 节 Harris 公式的出发点。

// region_stats.cpp —— 三个窗口的灰度方差：平坦区小、边缘中等、角点大
#include <opencv2/opencv.hpp>
#include <iostream>

// 以 (x, y) 为中心、边长 2*half+1 的窗口，返回灰度方差
double region_var(const cv::Mat& gray, int x, int y, int half) {
    cv::Rect win(x - half, y - half, 2 * half + 1, 2 * half + 1);
    cv::Scalar m, sd;            // mean 与标准差一次算完
    cv::meanStdDev(gray(win), m, sd);
    return sd[0] * sd[0];     // 方差 = 标准差平方
}

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);
    std::cout << "平坦区: " << region_var(gray, 10, 10, 3) << std::endl;
    std::cout << "边缘上: " << region_var(gray, 100, 100, 3) << std::endl;
    std::cout << "角点处: " << region_var(gray, 200, 200, 3) << std::endl;
    return 0;
}

10.2 Harris 原理：窗口移动与角点响应

1988 年，Chris Harris 和 Mike Stephens 发表了论文 A Combined Corner and Edge Detector，把"角点"这个艺术直觉变成了一个可计算的公式。核心思想一句话：角点是"往哪个方向挪动窗口，灰度都剧烈变化"的点。设窗口位移为 (u, v)，灰度变化量写成 E(u, v)：E(u, v) = Σ w(x, y) [I(x+u, y+v) - I(x, y)]²，其中 w(x, y) 是窗口函数（矩形窗或高斯窗）。平坦区里 E 处处很小；在边缘上，沿边缘方向挪窗口 E 很小、垂直方向很大；而在角点处，所有方向的 E 都大。官方教程《Harris Corner Detection》的推导正是从最大化这个 E 出发的。

对 E 做泰勒展开，可以把它近似成一个二次型：E ≈ [u v] M [u v]ᵀ，其中 M 是 2×2 的结构张量（梯度协方差矩阵）：M = Σ w [Ix², IxIy; IxIy, Iy²]，Ix、Iy 是图像在 x、y 方向的梯度（用 Sobel 算子算）。M 的两个特征值 λ1、λ2 直接刻画局部形状：两个都大 → 角点；一大一小 → 边缘；两个都小 → 平坦区。但解特征值太贵，Harris 用一个响应函数代替：R = det(M) - k * trace(M)²，det 是行列式、trace 是迹，k 是自由参数（典型 0.04 到 0.06）。官方文档对 cornerHarris 的描述与此一致：每个像素在 blockSize×blockSize 邻域上算 M，再算 dst(x,y) = det M - k * (tr M)²，角点是响应图的局部极大值。

把公式落成代码最能消除"玄学感"：下面的 harris_math.cpp 不调用 cornerHarris，而是用 Sobel 求梯度、逐像素乘法和高斯模糊把 M 的四个分量算出来，再合成 R——每一步对应公式里的一项。跑完会看到 R 是一张和原图一样大的浮点图，值有正有负，正是 10.3 节要处理的"得分地图"。

// harris_math.cpp —— 把 Harris 公式亲手落成代码：Ix/Iy -> M -> R
#include <opencv2/opencv.hpp>
#include <iostream>

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);

    // 第一步：Sobel 求梯度，输出 CV_32F 防溢出
    cv::Mat Ix, Iy;
    cv::Sobel(gray, Ix, CV_32F, 1, 0, 3);   // x 方向一阶导
    cv::Sobel(gray, Iy, CV_32F, 0, 1, 3);   // y 方向一阶导

    // 第二步：结构张量 M 的四个分量，mul 是逐像素相乘
    cv::Mat Ix2 = Ix.mul(Ix);
    cv::Mat Iy2 = Iy.mul(Iy);
    cv::Mat Ixy = Ix.mul(Iy);

    // 窗口求和：高斯模糊扮演 w(x,y)
    cv::GaussianBlur(Ix2, Ix2, cv::Size(3, 3), 0);
    cv::GaussianBlur(Iy2, Iy2, cv::Size(3, 3), 0);
    cv::GaussianBlur(Ixy, Ixy, cv::Size(3, 3), 0);

    // 第三步：角点响应 R = det(M) - k * trace(M)^2
    double k = 0.04;
    cv::Mat R = Ix2.mul(Iy2) - Ixy.mul(Ixy)
               - k * (Ix2 + Iy2).mul(Ix2 + Iy2);

    double mn, mx;
    cv::minMaxLoc(R, &mn, &mx);
    std::cout << "R 的范围: " << mn << " ~ " << mx << std::endl;
    return 0;
}

R 是浮点图，直接存成 PNG 会丢失信息，所以看响应图之前要归一化。下面这段把 R 线性缩放到 0 到 255 再保存——你会看到一张"得分地图"：越亮的地方角点响应越强，窗户角、房檐交点都在图上发亮，而墙面、天空一片漆黑。这张图就是 Harris 的"素描稿"，10.3 节只需要在这张稿子上画阈值。

// harris_visual.cpp —— 响应图归一化：把浮点 R 变成能看的 8 位图
#include <opencv2/opencv.hpp>
#include <iostream>

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);

    // cornerHarris 一步出响应图（内部已含梯度计算与邻域求和）
    cv::Mat R;
    cv::cornerHarris(gray, R, 3, 3, 0.04);

    // 响应值有正有负，先线性归一化到 0~255 再保存
    cv::Mat Rn;
    cv::normalize(R, Rn, 0, 255, cv::NORM_MINMAX, CV_8UC1);
    cv::imwrite("harris_response.png", Rn);

    // 响应图里越亮 = 角点响应越强，它就是一张"得分地图"
    double mn, mx; cv::Point pmin, pmax;
    cv::minMaxLoc(R, &mn, &mx, &pmin, &pmax);
    std::cout << "R 范围: " << mn << " ~ " << mx << std::endl;
    return 0;
}

10.3 cv::cornerHarris 用法：blockSize、ksize 与 k

工程视角看 API。官方签名是 void cv::cornerHarris(InputArray src, OutputArray dst, int blockSize, int ksize, double k, int borderType = BORDER_DEFAULT)，参数含义按官方文档逐一核实过：src 是输入单通道 8 位或浮点图；dst 是输出响应图，类型恒为 CV_32FC1、尺寸与 src 相同；blockSize 是结构张量求和的邻域尺寸（越大越"宏观"，角点越少越稳）；ksize 是 Sobel 算子的孔径参数（梯度核大小，须为奇数）；k 是响应公式里的自由参数。borderType 是边界外推方式，绝大多数场景用默认值即可。使用流程永远三步：转灰度 → 算响应 → 阈值化。

阈值怎么定？响应值的绝对大小随图像灰度尺度和尺寸变化，不能写死一个数。官方 Python 教程的示例就是以最大响应的 1% 作为阈值（dst > 0.01 * dst.max()），这是最常用的经验起点：先用 cv::minMaxLoc 取到最大值，乘以 0.01 得到阈值，再把响应图里超过阈值的像素在原图副本上画红点。C++ 里没有 Python 那种向量化掩膜，直接用 9.4 节学的 ptr 逐行扫描即可。另外两个相关的进阶函数值得记住：想要亚像素级精度时用 cv::cornerSubPix 精化角点坐标（官方教程《Harris Corner Detection》专门有一节讲它）；想要"最多 N 个、彼此至少距离 d"的工程化角点列表，用 cv::goodFeaturesToTrack，它内部也可以选择 Harris 响应作为质量度量。

// cornerharris_demo.cpp —— 标准用法：响应图 + 阈值 + 画圈
#include <opencv2/opencv.hpp>
#include <iostream>

int main() {
    cv::Mat img = cv::imread("building.jpg", cv::IMREAD_COLOR);
    cv::Mat gray;
    cv::cvtColor(img, gray, cv::COLOR_BGR2GRAY);

    // 1) 算响应图：src 单通道 8U/浮点，dst 恒为 CV_32FC1
    cv::Mat R;
    cv::cornerHarris(gray, R, 3, 3, 0.04);

    // 2) 阈值：以最大响应的 1% 为界（经验起点，需按图微调）
    double mx;
    cv::minMaxLoc(R, nullptr, &mx);
    double t = 0.01 * mx;

    // 3) 在原图副本上，把响应超阈值的像素画成红点
    cv::Mat out = img.clone();
    for (int r = 0; r < R.rows; ++r) {
        float* p = R.ptr<float>(r);
        for (int c = 0; c < R.cols; ++c) {
            if (p[c] > t) {
                cv::circle(out, cv::Point(c, r), 3,
                             cv::Scalar(0, 0, 255), -1);   // 半径 3 的实心红点
            }
        }
    }
    cv::imwrite("harris_corners.png", out);
    std::cout << "阈值: " << t << std::endl;
    return 0;
}

参数	含义	典型值
blockSize	结构张量求和的邻域边长，越大越"宏观"	3 或 5
ksize	Sobel 算子的孔径（核）大小，须为奇数	3
k	响应公式里的自由参数，越大判定越严格	0.04 ~ 0.06

// harris_tune.cpp —— 参数扫描：blockSize 与 k 如何影响角点数量
#include <opencv2/opencv.hpp>
#include <iostream>

// 返回"响应超过最大响应 1%"的像素个数
int count_corners(const cv::Mat& gray, int block, double k) {
    cv::Mat R;
    cv::cornerHarris(gray, R, block, 3, k);
    double mx;
    cv::minMaxLoc(R, nullptr, &mx);
    return cv::countNonZero(R > 0.01 * mx);  // 比较表达式直接得二值掩膜
}

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);
    std::cout << "blockSize=3, k=0.04 -> "
              << count_corners(gray, 3, 0.04) << std::endl;
    std::cout << "blockSize=5, k=0.04 -> "
              << count_corners(gray, 5, 0.04) << std::endl;
    std::cout << "blockSize=3, k=0.10 -> "
              << count_corners(gray, 3, 0.10) << std::endl;
    // 规律：blockSize 越大越"宏观"（角点更少更稳），k 越大越严格
    return 0;
}

10.4 Canny 原理：高斯、梯度、非极大值抑制与双阈值

1986 年，John F. Canny 提出了以他命名的边缘检测算法，它的贡献是把"找边缘"从单一步骤变成一条多阶段流水线，每一阶段解决一个具体问题。官方教程《Canny Edge Detection》按顺序给出五步：① 噪声抑制——边缘检测对噪声极其敏感，先用 5×5 高斯滤波去噪；② 求梯度——用 Sobel 核分别算出水平方向 Gx 与垂直方向 Gy，再合成梯度幅度 G = sqrt(Gx² + Gy²) 和梯度方向 θ = atan2(Gy, Gx)；③ 非极大值抑制；④ 双阈值（hysteresis）；⑤ 连通性追踪。后三步是 Canny 区别于朴素"梯度阈值化"的关键。

先说非极大值抑制。梯度幅度大的像素连成一片，边缘看起来是"粗带"，而我们要的是单像素宽的细线。做法是：梯度方向始终垂直于边缘，把 θ 归到水平、垂直、两条对角线四个方向之一，然后看当前像素沿梯度方向的两个邻居——只有当前像素是局部最大才保留，否则置零。官方教程用 A、B、C 三点说明：A 在边缘上，B、C 是 A 沿梯度方向的邻居，A 必须大于 B、C 才留下，抑制后得到"细边缘"二值图。双阈值则回答"哪些细线是真的"：设低阈值 minVal 和高阈值 maxVal，高于 maxVal 的必是边缘，低于 minVal 的必不是；介于两者之间的，看它是否与"确定边缘"连通——连着的算边缘（保留完整曲线），孤立的丢弃，这一机制顺带清掉了小噪点。官方教程的图例里：A 是确定边缘，C 虽然低于 maxVal 但连着 A 所以保留，B 与任何确定边缘都不连通所以丢弃。

把前两阶段落成代码，中间产物一目了然：归一化后的幅度图就是"所有候选边缘的亮度图"。下面 canny_stages.cpp 只做高斯模糊和 Sobel 梯度，输出幅度图；接着用一段"流水线图"把五阶段串起来——这正是 10.5 节 cv::Canny 内部替你完成的事。

// canny_stages.cpp —— Canny 前两阶段手动实现：高斯降噪 + Sobel 梯度
#include <opencv2/opencv.hpp>
#include <iostream>

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);

    // 阶段 1：5x5 高斯模糊去噪（边缘检测对噪声极其敏感）
    cv::Mat smooth;
    cv::GaussianBlur(gray, smooth, cv::Size(5, 5), 0);

    // 阶段 2：Sobel 求 Gx、Gy，再合成幅度与方向
    cv::Mat Gx, Gy;
    cv::Sobel(smooth, Gx, CV_32F, 1, 0, 3);
    cv::Sobel(smooth, Gy, CV_32F, 0, 1, 3);
    cv::Mat mag, ang;
    cv::magnitude(Gx, Gy, mag);           // G = sqrt(Gx^2 + Gy^2)
    cv::phase(Gx, Gy, ang, true);        // theta = atan2(Gy, Gx)，角度制

    // 幅度图归一化保存：越亮 = 梯度越大
    cv::Mat mag8;
    cv::normalize(mag, mag8, 0, 255, cv::NORM_MINMAX, CV_8UC1);
    cv::imwrite("canny_magnitude.png", mag8);
    std::cout << "幅度图已保存（非极大值抑制之前的粗边缘）" << std::endl;
    return 0;
}

// canny 五阶段流水线（数字对应官方教程的顺序）
// 1) 5x5 高斯模糊      -> 去掉噪声，防止把噪点当边缘
// 2) Sobel 梯度 Gx/Gy  -> 算出幅度 G 与方向 theta
// 3) 非极大值抑制       -> 只保留梯度方向上的局部最大，边缘变细
// 4) 双阈值 hysteresis -> 高于 maxVal 必保，低于 minVal 丢弃
// 5) 连通性追踪         -> 介于两者之间：连着强边缘就保留
// 结果：一张二值边缘图（白 = 边缘，黑 = 背景）

10.5 cv::Canny 用法：threshold1 与 threshold2

工程视角看 API。官方签名是 void cv::Canny(InputArray image, OutputArray edges, double threshold1, double threshold2, int apertureSize = 3, bool L2gradient = false)。参数逐个说：image 是输入单通道图（8 位即可）；edges 是输出二值边缘图（CV_8UC1，白为边缘）；threshold1 是低阈值 minVal，threshold2 是高阈值 maxVal——这是最容易记反的一对参数；apertureSize 是求梯度用的 Sobel 核大小，默认 3；L2gradient 决定幅度公式：true 用精确的 sqrt(Gx² + Gy²)，false 用更快的 |Gx| + |Gy|，官方教程说 true 更准确、默认 false 更快。Canny 内部自动完成高斯模糊等前置步骤，输入一张灰度图、输出一张二值图，一步到位。

阈值怎么选？官方教程特别强调："必须把 minVal 和 maxVal 选对才能得到正确结果"。经验法则：两个阈值保持 2:1 到 3:1 的比例，比如 50/150、100/200——minVal 太低会把纹理和噪声都当成边缘，太高又会让边缘断成碎片。调参的工程化做法是阈值扫描：固定图片，跑几组阈值，对比边缘像素占比（cv::countNonZero 统计白点数）和目视效果。注意输出是二值图，直接 imwrite 保存即可；想叠加到原图上看，要等 10.6 节的颜色转换。

// canny_demo.cpp —— 标准用法：一步出二值边缘图
#include <opencv2/opencv.hpp>
#include <iostream>

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);

    // threshold1 = 低阈值 minVal，threshold2 = 高阈值 maxVal
    cv::Mat edges;
    cv::Canny(gray, edges, 50, 150);     // 经典 1:3 组合
    cv::imwrite("canny_50_150.png", edges);

    // 边缘像素占比：countNonZero 统计白色像素
    int n = cv::countNonZero(edges);
    double ratio = 100.0 * n / (edges.rows * edges.cols);
    std::cout << "边缘像素: " << n << " (" << ratio << "%)" << std::endl;
    return 0;
}

// canny_thresholds.cpp —— 阈值扫描：同样的图，三组阈值三张结果
#include <opencv2/opencv.hpp>
#include <iostream>

void run(const cv::Mat& gray, double t1, double t2, const char* name) {
    cv::Mat edges;
    cv::Canny(gray, edges, t1, t2);
    cv::imwrite(name, edges);
    std::cout << name << ": 边缘像素 " << cv::countNonZero(edges) << std::endl;
}

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);
    run(gray, 30, 90,   "canny_30_90.png");    // 低阈值：边缘多而杂
    run(gray, 50, 150,  "canny_50_150.png");   // 中等：常用起点
    run(gray, 100, 300, "canny_100_300.png");  // 高阈值：只留强边缘
    // L2gradient=true 用 sqrt(Gx^2+Gy^2)，更精确但稍慢
    cv::Mat edges2;
    cv::Canny(gray, edges2, 50, 150, 3, true);
    cv::imwrite("canny_l2.png", edges2);
    return 0;
}

10.6 可视化叠加：把特征画回原图

检测结果是一堆坐标和二值图，怎么"看"？答案是画回原图。角点已在 10.3 节用 cv::circle 画过红点；边缘是二值图，直接叠会盖住原图，标准做法是先用 cv::cvtColor（COLOR_GRAY2BGR）转成三通道，再用 cv::addWeighted 按权重混合——原图 0.9、边缘图 0.4，得到半透明绿色描边。把两段逻辑放进同一个程序，就得到一张"角点红圈 + 边缘绿线"的完整特征图——这正是笔记素材里的真实经历：第一次用 OpenCV，一小时就同时做出了 Harris 与 Canny 两个检测器。API 上手极快，功夫在原理与调参。

另一条更工程化的可视化路线是 features2d 模块的 cv::KeyPoint 与 cv::drawKeypoints：把超阈值角点包装成 KeyPoint 向量（参数是 x、y 和特征尺寸），一次调用画完所有点。KeyPoint 自带尺度信息，是 SIFT/ORB 等现代特征检测器的统一输出格式——你现在就把 Harris 的结果"翻译"成了那个格式。综合实战闭环：读图 → 灰度 → Harris 角点 + Canny 边缘 → 叠加可视化 → 保存。

// visualize_overlay.cpp —— 角点红圈 + 边缘绿色叠加，画回原图
#include <opencv2/opencv.hpp>

int main() {
    cv::Mat img = cv::imread("building.jpg", cv::IMREAD_COLOR);
    cv::Mat gray;
    cv::cvtColor(img, gray, cv::COLOR_BGR2GRAY);

    // 1) 角点：响应图 + 阈值 + 红圈
    cv::Mat R;
    cv::cornerHarris(gray, R, 3, 3, 0.04);
    double mx;
    cv::minMaxLoc(R, nullptr, &mx);
    for (int r = 0; r < R.rows; ++r) {
        float* p = R.ptr<float>(r);
        for (int c = 0; c < R.cols; ++c) {
            if (p[c] > 0.01 * mx) {
                cv::circle(img, cv::Point(c, r), 3,
                             cv::Scalar(0, 0, 255), -1);   // 红点
            }
        }
    }

    // 2) 边缘：Canny 出二值图，转 BGR 后按权重叠回（绿线半透明）
    cv::Mat edges, edges_bgr;
    cv::Canny(gray, edges, 50, 150);
    cv::cvtColor(edges, edges_bgr, cv::COLOR_GRAY2BGR);
    cv::addWeighted(img, 0.9, edges_bgr, 0.4, 0, img);
    cv::imwrite("overlay.png", img);
    return 0;
}

// draw_keypoints.cpp —— 可选的 KeyPoint 路线：drawKeypoints 一步画完
#include <opencv2/opencv.hpp>
#include <opencv2/features2d.hpp>
#include <vector>

int main() {
    cv::Mat gray = cv::imread("gray.png", cv::IMREAD_GRAYSCALE);
    cv::Mat R;
    cv::cornerHarris(gray, R, 3, 3, 0.04);
    double mx;
    cv::minMaxLoc(R, nullptr, &mx);

    // 超阈值像素装进 KeyPoint（x 列、y 行，尺寸 5）
    std::vector<cv::KeyPoint> kps;
    for (int r = 0; r < R.rows; ++r) {
        float* p = R.ptr<float>(r);
        for (int c = 0; c < R.cols; ++c) {
            if (p[c] > 0.01 * mx) {
                kps.push_back(cv::KeyPoint(c, r, 5));
            }
        }
    }
    // drawKeypoints：每个点画带半径的圆圈
    cv::Mat out;
    cv::drawKeypoints(gray, kps, out, cv::Scalar(0, 0, 255));
    cv::imwrite("keypoints.png", out);
    return 0;
}

这一章我们完成了从"看得见图"到"看得懂图"的第一步：10.1 厘清特征概念（角点与边缘、可重复性），10.2 拆解 Harris 的窗口移动思想与响应公式，10.3 上手 cv::cornerHarris 并学会调 blockSize、ksize、k，10.4 走完 Canny 五阶段流水线，10.5 用 cv::Canny 出二值边缘图并做阈值扫描，10.6 把两类特征叠加回原图。三个钩子留给后面：一是 cv::cornerSubPix 的亚像素精化与 cv::goodFeaturesToTrack 的工程化角点输出，是特征工程的进阶方向；二是 Harris 与 Canny 的"慢"——每个像素都要算梯度、每个窗口都要算结构张量、Canny 更是五遍全图扫描，全是互不依赖的像素级重复劳动，这种形状的计算天然适合并行；三是这个"像素级并行"的答案，第 11 章的 CUDA 会用 GPU 把"特征检测为什么慢"彻底解开——特征检测正是 GPU 加速的最佳练兵场。