特征算法之Moravec算法

---

特征算法之Moravec算法

前言

最近一直在学习图像特征提取方面的算法，打算做一个博客系列作为总结。在特征算法的学习中，我认为主要关注以下几个方面。

• 特征的定义：要搞懂特征的定义，首先明白特征用做什么？我理解特征就是图像内容的一个关键信息，在图像检索、检测、匹配等后续操作中，不可能使用全部图像信息，所以必须找到图像中的一些关键信息。比如边缘、纹理、角点。大多数情况下，角点和特征点、特征可以作为等价概念。举个例子，比如要表达一个三角形，我们只需要找到这个三角形的三个顶点就能得到这个三角形的很多信息。（由于大量特征算法关注的都是角点特征，后续再提到特征的时候我默认为特征点、角点。不再提及边缘和纹理特征，当需要用到边缘和纹理特征我会指明边缘特征或纹理特征。）

  那么边缘、纹理、角点都可以作为特征，显然对特征的定义是不同的。对特征的定义直接决定了特征算法的原理和实现。

• 特征的检测：提出特征的检测（detection）是为了区分特征的描述（descriptor），在开始学习的时候，容易混淆这两个关键的概念。学习一个算法，首先要明白它是特征的检测算法，还是特征的描述算法，还是同时实现特征检测和描述。它们有什么区别呢？特征检测是提取图像中的特征点，即定位图像中角点的位置，理解为仅找到那些关键点在哪儿。

  故评价一个特征检测算法，需要从以下几个方面：速度、准确率（检测的角点中有没有伪角点、有没有是角点没有被检测到）、鲁棒性（关于鲁棒性，下文再介绍）

• 特征描述：特征描述算法很多时候能看到同意义的其他表达，如特征描述符、描述子。它说的是如何描述这个特征点，即可以认为特征检测算法找到特征点（的位置），而特征描述算法实现对这个特征点的内容的描述。

  那么评价一个特征描述算法，也需要从几个方面：速度、特征描述的区分度（我理解为可靠性，比如让不同角点的描述分布呈现高方差，可以先不深究这个东西，具体算法提到它的时候再细说）、鲁棒性。

• 特征匹配：说到特征匹配其实更多的是整个匹配的过程。以两张图片为例，图A、B内容可能近似，可能不同。首先在AB上进行特征检测算法，分别找出两个图像的特征点，然后用特征描述算法，分别对AB图像的特征点进行定量的描述（比如将每个特征点转换为一个数值的特征向量）。之后，将两个图像的所有特征点进行比对（匹配），然后实现目标的跟踪、图像的匹配等。

• 鲁棒性：特征算法的鲁棒性主要包括

  • 尺度不变性：不同图像尺度下，算法的结果。比如一幅图像，特征检测算法能检测到的角点，在图像放大10倍后，任然可以检测到那些角点吗。比如特征描述算法，在图像A中某角点，描述算法得到该角点的描述向量为(x1,x2,x3)，那么图像A放大10倍这个角点任然被检测到，但是它的描述还是(x1,x3,x3)吗。
  • 旋转不变性：即图像在经过旋转之后，角点还是能被检测到吗？角点的描述依然不变吗？
  • 抗噪：即噪声对算法的影响。
  • 仿射不变性：仿射变换后算法的结果是否受影响。至于仿射变换可以借用一句简单的定义：仿射变换=线性变换+平移。如何通俗地讲解「仿射变换」这个概念？

我认为特征算法的学习，主要是从其任务出发，着眼于如何解决几个鲁棒性问题。后续的博客计划探讨以下几个算法：

• Moravec(本文的后续内容将学习该算法)
• Harris
• FAST
• SIFT
• HOG
• Haar
• SUSAN
• ORB

---

Moravec算法原理

在介绍Moravec算法之前，先介绍一下角点的概念。

什么是角点？

角点，corner，目前没有明确的数学定义。广泛的理解为，角点就是极值点，在某方面（这个某方面的规则是自定义的）特别突出的点，根据这一特征可以自己定义角点，在这种定义下，“角点”、“兴趣点”、“特征点”通常指同一内容。狭隘理解定义到某一具体的研究对象，可以认为角点就是图像中的结构（轮廓）的拐点，如直线交点、线段终点，曲线局部曲率最大点。

角点的作用

角点是图像的很重要的局部特征，它决定图像中目标的形状。常用于三维场景重建，运动估计，目标跟踪，目标识别，图像配准等。

---

Moravec原理

Moravec角点检测算法是最早的角点检测算法之一。该算法将角点定义为具有低“自相关性”的点。算法会检测图像的每一个像素，将像素周边的一个邻域作为一个patch，并检测这个patch和周围其他patch的相关性。这种相关性通过两个patch间的平方差之和(SSD)来衡量，SSD值越小则相似性越高。

• 如果像素位于平滑图像区域内，周围的patch都会非常相似。SSD会很小。
• 如果像素在边缘上，则周围的patch在与边缘正交的方向上会有很大差异，在与边缘平行的方向上则较为相似。某一个方向SSD很小，其他方向存在很大差异。
• 如果像素是各个方向上都有变化的特征点，则周围所有的patch都不会很相似。所有方向SSD都会很大。

总结一下Moravec对于角点的定义

角点定义为具有低“自相关性”的点。

相关性：两个图块之间的平方差之和。

那么寻找角点的过程就是计算各个像素的各个方向上的SSD。

步骤

• 依次以每个像素点为锚点计算其5x5的窗口内4个方向的SSD，记为S1,S2,S3,S4
• 设定该点角点响应为：$CFR_{x,y} = min\{S1,S2,S3,S4\}$
• 设定一个阈值R，则角点响应值大于R的认为是候选角点。
• 局部非最大值抑制。（防止局部窗口内出现过各角点）

缺陷

Moravec会计算每个像素patch和周围patch的SSD最小值作为强度值，取局部强度最大的点作为特征点。
Moravec角点检测算法有几个很明显的缺陷：

• 强度值的计算并不是各向同性的，只有离散的8个45度角方向被考虑。因为patch的评议比较最多只有8个方向；
• 由于窗口是方形并且二元的，因此相应函数会有噪声；
• 对边缘的相应太简单，因为强度值尽取SSD的最小值；

---

Moravec算法实现

代码源地址

Moravec算子学习(http://blog.chinaunix.net/uid-13321460-id-3088421.html)

#include "opencv2/highgui/highgui.hpp"
#include "opencv2/imgproc/imgproc_c.h"
#include "opencv2/objdetect/objdetect.hpp"

#include <stdio.h>

#define BORDER_VALUE 2
#define THRESHOLD_VALUE 1000
#define NMS_OFFSET_VALUE 3

#define NMS     1
#define THRESHOLD   1
//#define SHOW_CORNERNESS_PIC   1

unsigned char * gray_image_mat;
int * cornerness_map_mat;
int cornerness_value_window[8];
int img_width;
int img_height;

struct offset_map {
    int x;
    int y;
};

struct offset_map cornerness_window_offset[9] = {
    {-1 , -1},
    { 0 , -1},
    { 1 , -1},
    {-1 ,  0},
    { 0 ,  0},
    { 1 ,  0},
    {-1 ,  1},
    { 0 ,  1},
    { 1 ,  1},
};

struct offset_map shift_offset[8] = {
    {-1 , -1},
    { 0 , -1},
    { 1 , -1},
    {-1 ,  0},
    { 1 ,  0},
    {-1 ,  1},
    { 0 ,  1},
    { 1 ,  1},
};

static CvScalar colors[] =
{
    {{0,0,255}},
    {{0,128,255}},
    {{0,255,255}},
    {{0,255,0}},
    {{255,128,0}},
    {{255,255,0}},
    {{255,0,0}},
    {{255,0,255}},
    {{255,255,255}}
};

int get_cornerness_value(int x , int y)
{
    int u,v;
    int a,b;
    int shift_count;
    int offset_count;
    int origin_value;
    int shift_value;
    int cornerness_value;
    int tmp;

    //printf("\n");

    for(shift_count = 0 ; shift_count < 8 ; shift_count++){
        u = shift_offset[shift_count].x;
        v = shift_offset[shift_count].y;

        cornerness_value = 0;

        for(offset_count = 0 ; offset_count < 9 ; offset_count++){
            int tmp = 0;

            a = cornerness_window_offset[offset_count].x;
            b = cornerness_window_offset[offset_count].y;

            origin_value = gray_image_mat[(y+b)*img_width + (x+a)];
            shift_value  = gray_image_mat[(y+v+b)*img_width + (x+u+a)];

            tmp = shift_value - origin_value;
            tmp = tmp*tmp;

            cornerness_value += tmp;
            //printf("(x->%d u->%d a->%d)(y->%d v->%d b->%d), origin->%d shift->%d cornerness_value->%d \n" , x , u , a , y , v , b , origin_value , shift_value , cornerness_value);
        }
        //printf("cornerness_value_window[%d][%d]->%d \n" , u , v , cornerness_value);
        cornerness_value_window[shift_count] = cornerness_value;
    }

    tmp = cornerness_value_window[0];

    for(shift_count = 0 ; shift_count < 8 ; shift_count++){
        if(cornerness_value_window[shift_count] < tmp){
        //if(cornerness_value_window[shift_count] > tmp){
            //printf("tmp->%d cornerness_value->%d \n" , tmp , cornerness_value_window[count]);
            tmp = cornerness_value_window[shift_count];
        }
    }

#ifdef THRESHOLD
    if((tmp != 0) && (tmp < THRESHOLD_VALUE)){
        //printf("tmp->%d \n" , tmp);
        tmp = 0;
    }
#endif
    return tmp;
}

void calc_cornerness_map()
{
    int x,y;
    int x_offset , y_offset;
    int count;
    int tmp , x_tmp , y_tmp;

    for(y = BORDER_VALUE ; y < (img_height - BORDER_VALUE) ; y++){
        for(x = BORDER_VALUE ; x < (img_width - BORDER_VALUE) ; x++){
            cornerness_map_mat[y*img_width + x] = get_cornerness_value(x , y);
            //printf("cornerness_map_mat[%d][%d]->%d \n" , x , y , cornerness_map_mat[y*img_width + x]);
            //printf("[%3d][%3d]->%3d " , x , y , cornerness_map_mat[y*img_width + x]);
        }
        //printf("\n");
    }

#ifdef NMS
    for(y = BORDER_VALUE ; y < (img_height - BORDER_VALUE) ; y++){
        for(x = BORDER_VALUE ; x < (img_width - BORDER_VALUE) ; x++){

            tmp = cornerness_map_mat[y*img_width + x];
            x_tmp = x;
            y_tmp = y;

            for(x_offset = -NMS_OFFSET_VALUE ; x_offset <= NMS_OFFSET_VALUE ; x_offset++){
                if((x_offset == 0) || ((x + x_offset) < 0) || (x + x_offset) > img_width){
                    //printf("x0 [%d][%d] pass\n" , (x + x_offset) , y);
                    continue;
                }

                count = y*img_width + x + x_offset; 

                if(cornerness_map_mat[count] <= tmp){
                    //printf("x1 [%d][%d]->%d \n" , (x + x_offset) , y , cornerness_map_mat[count]);
                    cornerness_map_mat[count] = 0;
                }else{
                    //printf("x2 [%d][%d]->%d \n" , (x + x_offset) , y , cornerness_map_mat[count]);
                    tmp = cornerness_map_mat[count];
                    cornerness_map_mat[y_tmp * img_width + x_tmp] = 0;
                    x_tmp = x + x_offset;
                    y_tmp = y;
                }
            }

            for(y_offset = -NMS_OFFSET_VALUE ; y_offset <= NMS_OFFSET_VALUE ; y_offset++){
                if((y_offset == 0) || ((y + y_offset) < 0) || (y + y_offset) > img_height){
                    //printf("y0 [%d][%d] pass\n" , (x + x_offset) , y);
                    continue;
                }

                count = (y + y_offset)*img_width + x; 

                if(cornerness_map_mat[count] <= tmp){
                    //printf("y1 [%d][%d]->%d \n" , (x + x_offset) , y , cornerness_map_mat[count]);
                    cornerness_map_mat[count] = 0;
                }else{
                    //printf("y2 [%d][%d]->%d \n" , (x + x_offset) , y , cornerness_map_mat[count]);
                    tmp = cornerness_map_mat[count];
                    cornerness_map_mat[y_tmp * img_width + x_tmp] = 0;
                    x_tmp = x;
                    y_tmp = y + y_offset;
                }
            }
        }
    }
#endif
}

int main(int argc , char** argv)
{
    int count;
    char * p;
    const char* gray_image_name = argc == 2 ? argv[1] : "test.png";

    IplImage* gray_image = cvLoadImage( gray_image_name , CV_LOAD_IMAGE_GRAYSCALE );
    IplImage* cornerness_image = cvCreateImage(cvGetSize(gray_image), 8, 3);
    cvCvtColor( gray_image , cornerness_image , CV_GRAY2BGR );

    if(!gray_image){
        printf("load gray image failed \n");
        return -1;
    }

    cvNamedWindow("gray_image", 1);
    cvNamedWindow("cornerness_image", 1);

    printf("-----image-----\n");
    printf("align->%d \n" , gray_image->align);
    printf("depth->%d \n" , gray_image->depth);
    printf("width->%d \n" , gray_image->width);
    printf("height->%d \n" , gray_image->height);
    printf("imageSize->%d \n" , gray_image->imageSize);
    printf("widthStep->%d \n" , gray_image->widthStep);

    if((gray_image->width < 5) || (gray_image->height < 5)){
        printf("picture's width and height must greater than 5 \n");
        return -1;
    }

    gray_image_mat = (unsigned char *)calloc(gray_image->width * gray_image->height , 1);
    cornerness_map_mat = (int *)calloc(gray_image->width * gray_image->height , 4);

    p = (char *)gray_image->imageData;

    for(img_height = 0 ; img_height < gray_image->height ; img_height++){
        for(img_width = 0 ; img_width < gray_image->width ; img_width++){
            count = img_height*gray_image->width + img_width;
            gray_image_mat[count] = p[count];
            cornerness_map_mat[count] = 0;
            //printf("[%3d][%3d]->%3d " , img_width , img_height , p[count] & 0xff);
        }
        //printf("\n");
    }

    img_height = gray_image->height;
    img_width = gray_image->width;

    printf("-----start-----\n");

    calc_cornerness_map();

    p = (char *)gray_image->imageData;
    for(img_height = 0 ; img_height < gray_image->height ; img_height++){
        for(img_width = 0 ; img_width < gray_image->width ; img_width++){
            count = img_height*gray_image->width + img_width;
            if(cornerness_map_mat[count] == 0){
#ifdef SHOW_CORNERNESS_PIC
                p[count] = 0;
#endif
            }else{
                printf("[%d][%d]->%d \n" , img_width , img_height , cornerness_map_mat[count]);
#ifdef SHOW_CORNERNESS_PIC  
                p[count] = 255;
#endif
                CvPoint center;
                int radius;
                center.x = cvRound(img_width);
                center.y = cvRound(img_height);
                radius = cvRound(9*1.2/9.*2);
                cvCircle( cornerness_image , center, radius, colors[0], 1, 8, 0 );
            }
            //p[count] = cornerness_map_mat[count];
        }
    }
    cvSaveImage("output_gray_image.png" , gray_image);
    cvSaveImage("output_cornerness_image.png" , cornerness_image);

    cvShowImage("gray_image", gray_image);
    cvShowImage("cornerness_image", cornerness_image);

    cvWaitKey(0);

    cvDestroyWindow("gray_image");

    return 0;
}

---

参考

如何通俗地讲解「仿射变换」这个概念

Moravec算子学习(http://blog.chinaunix.net/uid-13321460-id-3088421.html)