在計算機視覺中，特徵點的概念被大量用於解決物體識別、影象匹配、視覺跟蹤、三維重建等問題，比如影象中物體的角點，它們是在影象中可被輕易而精確地定位的二維特徵。顧名思義，特徵點檢測的思想是無需觀察整幅影象，而是通過選擇某些特殊點，然後對它們執行區域性分析。如果能檢測到足夠多的這種點，同時它們的區分度很高，並且可以精確定位穩定的特徵，那麼這個方法就很有效。這裡主要使用Harris特徵檢測器檢測影象角點。開發平臺為Qt5.3.2+OpenCV2.4.9。

在此之前，先給出OpenCV中cv::cornerHarris函式的呼叫方式：

    cv::cornerHarris(image, // 輸入影象
                     cornerStrength, //
 
 輸出為表示角點強度的32位浮點影象
                     3, // 導數平滑的相鄰畫素的尺寸
                     3, // 梯度計算的濾波器孔徑大小
                     0.01); // Harris引數

描述Harris運算元的經典論文可參考：

The article by C.Harris, M.J. Stephens, A Combined Corner and Edge Detector, Alvey Vision Conference, pp.147-152, 1988

The article by J. Shi and C. Tomasi, Good features to track, Int. Conference on Computer Vision and Pattern Recognition, pp. 593-600, 1994

The article by K. Mikolajczyk and C. Schmid, Scale and Affine invariant interest point detectors, International Journal of Computer Vision, vol 60, no 1, pp. 63-86, 2004

在第一篇論文中提到，Harris角點檢測最直觀的解釋是在任意兩個相互垂直的方向上，都有較大變化的點。為了定義一幅影象中的角點，Harris觀察一個假定的特徵點周圍小視窗內的方向性強度平均變化。

如圖所示，假設一個小視窗在影象上移動，在平滑區域如左圖所示，視窗在各個方向上均沒有變化。對於中間圖，小視窗在邊緣的方向上移動時也沒有變化。而知道小視窗移動到右圖的角點處，視窗在各個方向上均有明顯的變化。Harris角點檢測正是利用了這個直觀的物理現象，通過視窗在各個方向上的變化程度，決定是否存在著角點。這裡我們考慮偏移量(u,v)，則將影象視窗平移(u,v)產生的E(u,v)可表示為：

由以下公式可得到E(u,v)：

對於區域性微小的移動，E(u,v)可近似表達為：

其中M的詳細表達式為：

E(u,v)的橢圓形式如下：

E(u,v)是一個協方差矩陣，表現的是所有方向上強度的變化率。該定義涉及影象的一階導數，這通常是Sobel運算元的計算結果。而在OpenCV中cv::cornerHarris函式的第四個引數對應的正是用於計算Sobel濾波器的孔徑(aperture)。協方差的兩個特徵值給出了最大平均強度變化以及垂直方向上的平均強度變化，如果這兩個特徵值均較低，就認為當前是同質區域；如果其中一個特徵值較高，另外一個較低，則認為當前位於邊緣上；最後，若兩個特徵值均較高，則判定當前位於角點處。

因此，定義角點響應函式R，其中k為函式cv::cornerHarris中的最後一個引數；之後，對R進行閾值處理，設定若R大於閾值threshold，則提取出區域性極大值：

下面記錄一下harris角點檢測的幾種方案。

一、基本的Harris角點檢測實現

#include <QCoreApplication>
#include <opencv2/core/core.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/highgui/highgui.hpp>

int main(int argc, char *argv[])
{
    QCoreApplication a(argc, argv);
    // Harris法的簡單實現
    cv::Mat image = cv::imread("c:/031.jpg", 0);
    cv::Mat cornerStrength;
    cv::cornerHarris(image,
                     cornerStrength, // 輸出為表示角點強度的32位浮點影象
                     3, // 導數平滑的相鄰畫素的尺寸
                     3, // 梯度計算的濾波器孔徑大小
                     0.01); // Harris的相關引數

    // 要在視窗輸出，需要轉化為CV_8U格式，閾值化即可
    // 角點強度的閾值
    cv::Mat harrisCorner;
    double threshold = 0.0001;
    cv::threshold(cornerStrength,
                  harrisCorner,
                  threshold,
                  255,
                  cv::THRESH_BINARY_INV); // 輸出為翻轉的二值影象

    cv::namedWindow("Original Image");
    cv::imshow("Original Image", image);
    cv::namedWindow("Harris Corner");
    cv::imshow("Harris Corner", harrisCorner);

    return a.exec();
}

得到的結果為二值影象，可以看到影象中角點的位置包含許多圓圈，這與精確定位特徵點的目標相悖：

二、改進的Harris角點檢測實現(Shi-Tomasi演算法)

這裡通過封裝自定義類來改進角點檢測的效果。定義一個類HarrisDetector（其中已封裝了Harris引數和相關函式）：

#ifndef HARRISDETECTOR_H
#define HARRISDETECTOR_H
#include <opencv2/core/core.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/highgui/highgui.hpp>

class HarrisDetector
{
private:
    // 表示影象角點強度的32位浮點影象
    cv::Mat cornerStrength;
    // 將輸出影象閾值化後的32位浮點影象
    cv::Mat cornerThreshold;
    // 區域性極大值影象
    cv::Mat localMax;

    // 以下三個為cornerHarris函式的必要引數
    // 導數平滑的相鄰畫素的尺寸
    int neighbourhoodPixelSize;
    // 濾波器的孔徑大小
    int aperture;
    // Harris引數
    double k;

    // 閾值計算的最大強度
    double maxStrength;
    // 計算得到的閾值
    double threshold;
    // 非極大值抑制的相鄰畫素的尺寸
    int noneighbourhoodPixelSize;
    // 非極大值抑制的核
    cv::Mat kernel;

public:
    // 初始化引數
    HarrisDetector():neighbourhoodPixelSize(3),
                     aperture(3),
                     k(0.01),
                     maxStrength(0.0),
                     threshold(0.01),
                     noneighbourhoodPixelSize(3)
    {
        setLocalMaxWindowSize(noneighbourhoodPixelSize);
    }

    // 建立非極大值抑制的核
    void setLocalMaxWindowSize(int size);

    // 計算Harris角點
    void detect(const cv::Mat &image);

    // 由Harris的值獲取角點圖
    cv::Mat getCornerMap(double qualityLevel);

    // 由Harris的值獲取特徵點
    void getCorners(std::vector<cv::Point> &points, double qualityLevel);

    // 由角點圖獲取特徵點
    void getCorners(std::vector<cv::Point> &points, const cv::Mat& cornerMap);

    // 在特徵點的位置繪製圖
    void drawOnImage(cv::Mat &image,
                     const std::vector<cv::Point> &points,
                     cv::Scalar color = cv::Scalar(255, 255, 255),
                     int radius = 4, int thickness = 2);
};


#endif // HARRISDETECTOR_H

接著，在harrisdetector.cpp中定義各個函式和初始化的引數：

#include "harrisdetector.h"

// 建立非極大值抑制的核
void HarrisDetector::setLocalMaxWindowSize(int size)
{
    noneighbourhoodPixelSize = size;
    kernel.create(noneighbourhoodPixelSize, noneighbourhoodPixelSize, CV_8U);
}

// 計算Harris角點
void HarrisDetector::detect(const cv::Mat &image)
{
    // Harris計算
    cv::cornerHarris(image,cornerStrength,
                     neighbourhoodPixelSize,
                     aperture,
                     k);
    // 內部閾值計算
    double minStrength; // 未使用
    cv::minMaxLoc(cornerStrength,
                  &minStrength,
                  &maxStrength);
    // 區域性極大值檢測
    cv::Mat dilate; // 臨時影象
    cv::dilate(cornerStrength, dilate, cv::Mat());
    cv::compare(cornerStrength, dilate, localMax, cv::CMP_EQ);
}

// 由Harris的值獲取角點圖
cv::Mat HarrisDetector::getCornerMap(double qualityLevel)
{
    cv::Mat cornerMap;
    // 對角點影象進行閾值化
    threshold = qualityLevel * maxStrength;
    cv::threshold(cornerStrength, cornerThreshold,
                  threshold,255,cv::THRESH_BINARY);
    // 轉換為8點陣圖像
    cornerThreshold.convertTo(cornerMap, CV_8U);
    // 非極大值抑制
    cv::bitwise_and(cornerMap, localMax, cornerMap);
    return cornerMap;
}

// 由Harris的值獲取特徵點
void HarrisDetector::getCorners(std::vector<cv::Point> &points, double qualityLevel)
{
    // 得到角點圖
    cv::Mat cornerMap = getCornerMap(qualityLevel);
    getCorners(points, cornerMap);
}

// 由角點圖獲取特徵點
void HarrisDetector::getCorners(std::vector<cv::Point> &points, const cv::Mat& cornerMap)
{
    // 遍歷畫素得到所有特徵
    for( int y = 0; y < cornerMap.rows; y++ )
    {
        const uchar* rowPtr = cornerMap.ptr<uchar>(y);
        for( int x = 0; x < cornerMap.cols; x++ )
        {
            // 如果是特徵點
            if (rowPtr[x])
            {
                points.push_back(cv::Point(x,y));
            }
        }
    }
}

// 在特徵點的位置繪製圖
void HarrisDetector::drawOnImage(cv::Mat &image,
                 const std::vector<cv::Point> &points,
                 cv::Scalar color,
                 int radius, int thickness)
{
    std::vector<cv::Point>::const_iterator it = points.begin();
    // 對於所有角點，繪製白色圓圈
    while(it != points.end())
    {
        cv::circle(image, *it, radius, color, thickness);
        ++ it;
    }
}

最後，使用該類的步驟如下，直接修改main函式：

#include <QCoreApplication>
#include <opencv2/core/core.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <QDebug>
#include "harrisdetector.h"

int main(int argc, char *argv[])
{
    QCoreApplication a(argc, argv);

    cv::Mat image = cv::imread("c:/031.jpg", 0);
    // 建立Harris物件
    HarrisDetector harris;
    // 計算Harris的值
    harris.detect(image);
    // 檢測Harris的角點
    std::vector<cv::Point>pts;
    harris.getCorners(pts, 0.1);
    // 繪製角點圖
    harris.drawOnImage(image, pts);
    cv::namedWindow("Harris Corners Final");

    return a.exec();
}

生成的影象：

這裡為了改進特徵點檢測結果，添加了額外的非極大值抑制步驟，目的是移除彼此相鄰的Harris角點。這就要求Harris角點不只需要得分高於給定閾值，它還必須是區域性極大值。在檢測中使用了一個技巧，即將Harris得分的影象進行膨脹：

cv::dilate(cornerStrength, dilate, cv::Mat());

這是由於膨脹運算替換每個畫素值為相鄰範圍內的最大值，因此只有區域性極大值的點才會保留原樣，並通過以下函式進行測試：

cv::compare(cornerStrength, dilate, localMax, cv::CMP_EQ);

其中，localMax矩陣僅在區域性極大值的位置為真，因此又可以在getCornerMap函式中用它來抑制所有非極大值的特徵（基於cv::bitwise_and函式）。

三、引入適合跟蹤的優質特徵的Harris檢測實現

在浮點處理器的幫助下，為了避免特徵值分解而引入的數學上的簡化變得微不足道，因此Harris檢測可以基於計算而得的特徵值。原則上這個修改不會顯著影響檢測的結果，但是能夠避免使用任意的k引數。

以上第二中方法引入了局部極大值的條件，改善了部分效果。然而，特徵點傾向於影象中不均勻分佈、普遍集中在紋理豐富的部分。這裡記錄一種新的解決方案：

該方案利用兩個特徵點之間的最小距離，從Harris得分最高的點開始，僅接受離開有效特徵點距離大於特定值的那些點。在OpenCV中提供cv::goodFeaturesToTrack實現這一方法，它檢測到的特徵能用於視覺跟蹤應用中的優質特徵集合。其呼叫方式如下：

    // 計算適合跟蹤的優質特徵
    std::vector<cv::Point2f> corners;
    cv::goodFeaturesToTrack(image,corners,
        250,    // 返回的最大特徵點的數目
        0.01,   // 質量等級
        8); // 兩點之間的最小允許距離

除了質量等級閾值、特徵點之間的最小允許距離，該函式還需要指定返回的最大特徵點數目，這是因為特徵點是按照強度進行排序的。以下給出該方法的實現程式碼，直接在main函式中新增：

    // 輸入影象
    cv::Mat image= cv::imread("c:/031.jpg",0);

    // 計算適合跟蹤的優質特徵
    std::vector<cv::Point2f> corners;
    cv::goodFeaturesToTrack(image,corners,
        250,    // 返回的最大特徵點的數目
        0.01,   // 質量等級
        8); // 兩點之間的最小允許距離

    // 遍歷所有特徵點並畫圓圈
    std::vector<cv::Point2f>::const_iterator it= corners.begin();
    while (it!=corners.end())
    {
        cv::circle(image, *it, 3, cv::Scalar(255,255,255), 2);
        ++it;
    }

    // 顯示輸出結果
    cv::namedWindow("Good Features to Track");
    cv::imshow("Good Features to Track",image);

返回生成的結果：

可以看到，該方法顯著改進了特徵點的分佈情況，但是這樣也增加了檢測的複雜度，因為要求特徵點要安裝Harris的得分進行排序。該函式也可以指定一個可選的引數，使得按照經典的焦點分數定義進行計算。

其中，cv::goodFeaturesToTrack函式擁有一個封裝類cv::GoodFeatureToTrackDetector，它繼承自cv::FeaturesDetector類。其用法與以上的Harris類相類似：

    // 特徵點向量
    std::vector<cv::KeyPoint> keypoints;
    cv::goodFeaturesToTrackDetector gftt(
        250,    // 返回的最大特徵點的數目
        0.01,   // 質量等級
        8); // 兩點之間的最小允許距離
    // 使用FeatureDetector的函式進行檢測
    gftt.detect(image, keypoints);

結果與先前得到的結果是一樣的，因為它們呼叫的是同一個函式。

關於Harris的理論研究有待進一步研究……