上文做了说明,我们用像素点的海塞矩阵来描述一副图像的关键特征点。下面就来讲解图像的海塞矩阵在图像处理中如何计算,海塞矩阵就是求X,Y方向的二阶偏导。
公式2.1
以上公式就是对每一个像素点求出一个海赛矩阵,对于X方向和Y方向的二阶偏导。特征点要具备尺度无关性(这样做的目的就是进一步筛选那些关键点,在不同尺度下都达到要求的点,才列为关键点),所以在计算海赛矩阵以前,要对图像进行高斯滤波,得到不同尺度的图像。并且计算不同尺度下的图像特征点。那是什么是图像的尺度?通俗的讲,一副图像在不同远近的情况下,清晰程度不同,越远越模糊,越近越清楚。这图像处理领域我们可以通过高斯滤波来模拟不同远近条件下图像的清晰程度不同,来模拟人眼在不同远近对物体的成像。
不同尺度的图像得到,通过高斯滤波,高斯滤波在二维空间的定义为:
公式2.2
图2.1 matlab中图像表示
在二维空间中,这个公式生成的曲面的等高线是从中心开始呈正态分布的同心圆。分布不为零的像素组成的卷积矩阵与原始图像做变换。每个像素的值都是周围相邻像素值的加权平均。原始像素的值有最大的高斯分布值,所以有最大的权重,相邻像素随着距离原始像素越来越远,其权重也越来越小。这样进行模糊处理比其它的均衡模糊滤波器更高地保留了边缘效果。其中,Sigma是正态分布的标准差,值越大,图像越模糊(平滑)。如二维模板大小为m*n,则模板上的元素(x,y)对应的高斯计算公式为:
公式2.3
通过以上公式我们可以计算出不同σ值下,不同大小的高斯模板。下图是5*5的高斯模板卷积计算示意图。高斯模板是中心对称的。
图2.2
图2.3
以上图片代表不同尺度因子σ去不同值的结果。由此得到不同尺度的图像,如图2.3
图像尺度不同实际就是图像的模糊程度不同,就像一副图像在人眼中表示远近不同,越近越清楚,越远越模糊。所以用高斯滤波可以很好的模拟不同尺度的同一幅图像。
图2.4
如图2.4所示,和SIFT算法不同的是SURF算法采用不同的高斯模板大小得到不同尺度的图像,一共分为四组,每组四个不同尺度的图像。就像最上面一层第一副图像的高斯模板为9X9,第二幅为15X15….以此类推。下面说明其中一副图像如何计算得到,其他图层就是采用不同大小的高斯模板,也可一样得到。
经过高斯滤波的海赛矩阵公式如下:
公式2.5
由于求海赛矩阵时先对图片进行高斯平滑,而后求二阶偏导,在图像处理中一般使用模板卷积形成。模板就是数学公式求值在图像上针对离散像素点的简化,模板就是在程序中就是一个矩阵。矩阵的成员数值通过数学计算得到。为了便于计算海赛矩阵行列式,还引入了积分图想的概念。对于一幅灰度的图像,积分图像中的任意一点(x,y)的值是指从图像的左上角到这个点的所构成的矩形区域内所有的点的灰度值之和。下面介绍积分图如何计算:
图2.4(灰度原图) 图2.5(积分图)
用公式表示就是:
图2.6
区域1 : = sum(A);
区域2 : = sum(A + B);
区域3 : = sum(A + C);
区域4 : = sum(A + B + C + D);
所以,如果需要计算D区域中的灰度和,则
sum(D) = 区域4 - 区域2 - 区域3 + 区域1 (都是灰度值)。
很明显,这里仅仅只需要通过查表得到 1、2、3、4点的积分图像的值即可得到。
下面给出积分图像的C语言程序片段:
void integral_matrix_Old(unsigned char* srcImage,unsigned int* sumMatrix,int width,int height,int line_length )
{int y,x,k;
int offset = line_length - width;
int srcstep = width;
int s = 0;//sum的第0列的值都是0
int sumstep = width+1;
int cn = 1;
memset( sumMatrix, 0,(width+1)*sizeof(sumMatrix[0]));/*将sum的第一行设为0*/
sumMatrix += sumstep + cn;
for( y = 0; y < height; y++, srcImage += srcstep - cn+offset, sumMatrix += sumstep - cn )
{
for( k = 0; k < cn; k++, srcImage++, sumMatrix++ )
{
s = sumMatrix[-cn] = 0;//sum的第0列的值都是0
for( x = 0; x <width; x += cn )
{
s += srcImage[x];
sumMatrix[x] = sumMatrix[x - sumstep] + s;//积分图像
}
}
}
}
在计算图像特征点的时候,surf算法是采用计算原图中每个像素的的海赛矩阵行列式的近似值构成,其行列式计算公式如下:
公式2.4
其中0.9是原算法作者给出的经验值,有一套理论计算,具体请看原文论文。算法是在不同尺度图像下进行求海赛矩阵的,不同尺度图像通过高斯平滑得到,然后求不同方向的二阶偏导。在离散图像像素中,用一个图像处理模板代替就可以了。比如X方向,Y方向,XY方向上的模板如图2.7,我们通过X方向的模板求出Dxx,通过Y方向模板求出Dyy,通过XY方向模板求出Dxy. 最后分别求出之后带入公式2.4就可以得到每个像素海塞矩阵行列式的值。
图2.7
图2.7分别表示在X方向,Y方向,XY方向的运算模板。
下面就Y方向的滤波器详细说明这些箱式滤波器的设计原理
图2.8
上图的左边即用高斯平滑然后在y方向上求二阶导数的模板,为了加快运算用了近似处理,其处理结果如右图所示,这样就简化了很多。请看简化的区域有两块2X9的灰色区域,两块3X5的白色区域,一块3X5的黑色区域组成。这个9X9的模板被划分成5块不同区域。对于高斯滤波我们知道,对于一个像素的高斯滤波计算,是对周围9X9区域的像素的加权平均得到。区域内不同位置的像素,权值不一样,权值的意义是当前点像素的值对最终滤波值的影响程度。且看简化图像五个不同区域的权值为白色区域为1,黑色区域为-2,灰色区域为0。我们规定了不同区域权值大小。那么这个区域的像素灰度值乘以权值,就是这个区域对最终高斯滤波值的贡献。那么这个区域的灰度值如何得到,前面计算了积分图像就很容易得到这个区域的灰度值的总大小。如此可以看出极大的简化了计算。其他方向的原理类似,这里不做解释。注意1,-2,1的设计,这个权值的设计让我们完成了二阶偏导和高斯滤波的计算。图像中离散数据的求导运算,简化为减法运算。
其滤波器初始数据如下:
const int dx_s[3][5] = { {0, 2, 3, 7, 1}, {3, 2, 6, 7, -2}, {6, 2, 9, 7, 1} };
const int dy_s[3][5] = { {2, 0, 7, 3, 1}, {2, 3, 7, 6, -2}, {2, 6, 7, 9, 1} };
const int dxy_s[4][5] = { {1, 1, 4, 4, 1}, {5, 1, 8, 4, -1}, {1, 5, 4, 8, -1}, {5, 5, 8, 8, 1} };
void ResizeHaarPattern( const int* src, SURF_HARRI_FILTER* dst, unsigned int n, unsigned int oldSize, unsigned int newSize, unsigned int widthStep )
{
float ratio = (float)newSize/oldSize;
unsigned int k=0;
for( k = 0; k < n; k++ )
{
int dx1 =(int) ( ratio*(src+k*5)[0]+0.5f);/*+0.5主要是为了四舍五入*/
int dy1 = (int)( ratio*(src+k*5)[1]+0.5f);
int dx2 =(int) ( ratio*(src+k*5)[2]+0.5f );
int dy2 = (int)( ratio*(src+k*5)[3]+0.5f);
dst[k].p0 = dy1*widthStep + dx1;//第一个矩形相对于图像原点的坐标
dst[k].p1 = dy2*widthStep + dx1;//第二个矩形相对于图像原点的坐标
dst[k].p2 = dy1*widthStep + dx2;//第三个矩形相对于图像原点的坐标
dst[k].p3 = dy2*widthStep + dx2;//第四个矩形相对于图像原点的坐标
dst[k].w = (src+k*5)[4]/((float)(dx2-dx1)*(dy2-dy1));//这个矩形内的数据所占的权重
}
}
上面的代码计算得到那些2X9,3X5的矩形区域相对当前像素点的坐标。知道这些子区域的位置坐标乘以他们所代表的权重,就是这块区域对当前像素点的贡献。下面举例计算在Y方向上的矩形区域,得到Dyy的过程:
分别设d1_y ,d2_y,d3_y为三个不同区域的像素灰度权值和,具体如图所示:
图2.9
d1_y=(s_ptr[pDy_0->p0]+s_ptr[pDy_0->p3] - s_ptr[pDy_0->p1]-s_ptr[pDy_0->p2])*pDy_0->w;
d2_y=(s_ptr[pDy_1->p0]+s_ptr[pDy_1->p3] - s_ptr[pDy_1->p1] - s_ptr[pDy_1->p2])*pDy_1->w;
d3_y=(s_ptr[pDy_2->p0]+s_ptr[pDy_2->p3] - s_ptr[pDy_2->p1] - s_ptr[pDy_2->p2])*pDy_2->w;
看上面一段程序s_ptr[pDy_0->p0]+s_ptr[pDy_0->p3] - s_ptr[pDy_0->p1]-s_ptr[pDy_0->p2]这就代d1_y这个矩形区域里面的灰度值的和是多少,就是根据前面计算的积分图得到的。那么DYY= d1_y + d2_y + d3_y,同理我们可以计算DXX,DXY,那么根据公式2.4,我们就可以计算出经过高斯滤波的每个像素点的海塞矩阵行列式的值大小。我们设定一个阈值来挑选大于某一行列式值的像素。把大于某一值的像素作为一个关键点的候选点。
图2.10
拿出这些关键点的候选点,在两个相邻尺度的的海赛矩阵行列式集合上的26个点做比较。注意因为我们设定一组由四副尺度不同的图像,那么注定第一幅和最后一幅图像向上或者向下缺少一个邻域,这里我们不去计算即可。如果这个候选关键点的海赛矩阵行列式的值确实大于其他26个点的海赛矩阵行列式的值,那么这个关键点就被列为关键点。从这里我们得出一个这个算法的原理,对于关键点就是一副图像上的角点,边缘点,暗区域的亮点,亮区域的暗点,这些点其实都对应一条数学原理,就是这些点的二阶编导较大。说白了就是图像在这个区域灰度值变化比较激烈。至于为何要进行高斯模糊,其实得到的是不同尺度下的图像,这样相邻的尺度比较得到的关键点,我想更为稳定一些,有效的去除那些干扰点。以上就是针对一副图像怎么计算关键点。两幅一样的图像或者一副图像包涵另一幅图像的情况,总会有位置一样的关键件。如果只在考虑图像的水平方向或者垂直方向,从这些关键点中匹配到目标图像的位置,也是可以了。但是我们要考虑,目标位置会发生旋转,角度偏移的情况。那么我们就给这些关键点,找到他们的主方向,让这些关键点有大小有方向,那么就可以应对目标位置发生角度偏移的情况了。
三 SURF算法关键点的描述符
3.1 什么是关键特征点的主方向
图3.1
原算法中规定,以前面得到关键点为中心,计算一半径6计算半径为6
(为特征点所在的尺度值)的邻域内的点在x,y方向的Haar小波(Haar小波边长取4)响应,并给这些响应值赋予高斯权重系数,使得靠近特征点的响应贡献大,而远离特征点的响应贡献小。在360度范围内,每次扫描一个60度的扇区,把这个扇区内的Haar小波响应相加形成新的矢量,选择最长的方向为该特征点的主方向(如图3.1)。
注:Harr小波特征:Haar特征分为三类:边缘特征、线性特征、中心特征和对角线特征,组合成特征模板。特征模板内有白色和黑色两种矩形,并定义该模板的特征值为白色矩形像素和减去黑色矩形像素和。Haar特征值反映了图像的灰度变化情况。
图3.2
所示为哈尔小波滤波器在x方向和y方向上描述,黑色部分的权重是-1,白色部分的权重是1。计算出图像在哈尔小波的x和y方向上的响应值之后,对两个值进行因子为的高斯加权,加权后的值分别表示在水平和垂直方向上的方向分量。Harr特征值反应了图像灰度变化的情况,那么这个主方向就是描述那些灰度变化特别剧烈的区域方向。就像图3.1第三幅图像的样子。
3.2 构造关键点的主方向的描述子
我们已经得到了关键点的主方向了,然后以这个关键点为中心,画一个以边长为20S(S为该特征点所在的尺度),以该关键点的主方向为转角的正方形。如图:
图 3.3
V =
图3.4
然后把这些方框分为16个子区域,每个区域统计25个像素的水平方向和垂直方向的Harr小波特征,当然这里的水平和垂直方向都是相对主方向而言。该haar小波特征为水平方向值之和,水平方向绝对值之和,垂直方向之和,垂直方向绝对值之和以及分别的负方向之和。所以每个特征点就是16*8=128维的向量.
3.3 关键点的匹配
当搜索区域和搜索目标的关键点的主方向和描述子都计算出来之后,下一步进行关键点匹配,就是说对于目标图片上的任意关键点A,在搜索区域中找到最相似的点,与之对应。每一个关键都有一个描述向量,我们采用关键点描述向量的欧式距离来作为两幅图像中关键点的相似性判定度量。取图像1中的某个关键点,并找出其与图像2中欧式距离最近的前两个关键点,在这两个关键点中,如果最近的距离除以次近的距离少于某个比例阈值,则接受这一对匹配点。调整这一比例阈值可以应对不同的匹配情况。这样一来我们可以得到,很多个关键点。
3.4区域匹配和旋转角度计算
上面的计算我们就得到了,各个关键点的特征向量。在工业检测里面,我们必须精确跟踪目标区域的位置,下面我讲解如何通过得到的关键点计算出目标区域的矩形坐标。
图3.5
如图3.5所示我们要在search的区域,跟踪object区域。方便理解我们构造了坐标系。
设object’为目标区域,search为搜索区域。经过前面的计算,我们得到了在object’里面的i,j,k和serach里面的is,js,ks这几个对应的关键点。那么我们要通过这些个对应的关键点,来在search计算出object的目标位置。通过这三个点,我们实际上构造出了,一对三角形。如果在object’里面的i,j,k关键点serach里面的is,js,ks对应的位置就是图像中对应的几何位置的话(就像i点对应object中M字母的一个角点,is点对应search中M字母的一个角点一样).那么着六个点所构造的两个三角形,就是相似三角形。我们可以通过几何位置关系求出他们各个对应边的比值。在意图图像中的关键点有很多这样的对应的三角形。我们设定一个阈值,如果是完全相似的话,各个边对应的比值是一个固定的值,如果三角形一模一样比值就是1.在图像处理中考虑计算误差等图像干扰。不可能比值如完全相等的。我们寻找多个这样的相对三角形,做计算,如果比例对应关系越好说明,我们寻找的这三个关键点的精确度越高。我们就把这三个点作为最佳匹配点。然后我们取这两个三角形的两条对应边,我们计算这两个边的角度,就可以得到目标区域的转角,就可以跟踪如果目标区域带角度偏移的情况。
以上就是整个算法的介绍,如有描述错误欢迎批评指正。