SIFT特征--构造DOG尺度空间

    最近打算做全景图拼接,尝试过hugin和opencv stitch,一直没有很满意的效果。打算深入研究,恰好,在github看到清华一小哥的项目,故逐步分析他的代码。

    sift特征是全景拼接的第一步,相当基础,但sift特征提取与匹配从0开始实现起来又挺折腾的,故花点时间整理下。

    清华小哥的博客:http://ppwwyyxx.com/2013/SIFT-and-Image-Stiching/ 

    我的github连接:https://github.com/haopo2005/panorama.git

    在ubuntu平台下跑,eigen安装后,记得添加makefile中eigen头文件所在目录,其他没啥了。

    

一、特征提取与sift特征简介

1.      整体特征与局部特征

在做目标识别与分类时,传统的方法大都将目标作为一个整体,从大量包含目标的图片集中学习并提取整体特征,如面积、周长、不变矩和傅立叶描绘子等,并进行统计分类。这种识别方法有很大缺陷:对结构复杂的图像,识别效果受到图像分割精度的制约,或者需要大量的时间学习大量的数据,又比如目标的形状发生较大变化时(目标遮挡),目标的整体特征会突然变化,这些对于目标识别不利。

  而人类的视觉系统可以将物体分解为多个许多有意义的小块,并通过这些局部信息进行目标的辨识。利用局部特征对图像进行描述,在图像内容复杂、噪声干扰大、存在局部遮挡、目标姿态发生较大变化等情况下,利用局部信息进行目标识别是非常有效的。

 

2.      局部特征区域检测方法分类

a.      密集选取方法:

研究者认为,在模式识别的低层处理中,所有图像区域都有一定的作用,丢失任何细节都可能对最终效果产生很大的影响。常见做法:

(1)    将图像密集地分为互不重叠的特征窗口,每个特征窗口都当作一个局部特征区域

(2)    以整幅图像的每一个像素点为中心,选取周围的区域作为局部特征区域

(3)    在检测窗口的每个像素位置、不同尺度下提取大量的特征区域,供进一步使用

密集选取方法在滑动窗口模型中应用很多,其优点是基本没有丢失图像的细节,可以得到非常丰富局部特征。但是其中很大一部分特征区域信息量过小,对后期的识别没有作用甚至起到干扰作用,加重了下一步特征优化工作。

 

b.      稀疏选取方法:

通过特征检测,选取具有显著特点的图像区域作为局部特征。检测算子一般可以分为基于形状的和基于外观的。

(1)    基于形状的检测算子是根据图像的边界、直线、弧线等来确定特定区域的位置。主要用于外形区分度明显的目标识别,比如canny特征

(2)    基于外观的检测算子是在图像的灰度模式下,搜索具有某种稳定性和不变性的特征点或关键区域。比如Hessian矩阵、Harris角点、高斯差分DOG、SIFT特征

稀疏选取法检测的特征区域数量一般在200~3000,其主要优点是简洁、紧致,图像的关键点远少于图像的像素,使得后面的识别过程能大大加速。但很多特征区域检测算法往往和图像的特性有关,应用到通用的目标识别时,可能会有一定的局限。

 

c.      其他方法

有人在研究向量空间模型的取样策略时发现,当训练集的样本数足够多时,随机取样法能达到和某些稀疏取样相近甚至更好的结果。比如使用显著性映射在分类过程中动态选取图像块的方法。

 

3.      sift特征简介

SIFT算子(Scale-Invariant FeatureTransform),是一种图像的局部描述子,具有尺度、旋转、平移的不变性,而且对光照变化、仿射变换和3维投影变换具有一定的鲁棒性。另外,它还具有以下有点:

a.很好的独特性,适于在海量特征数据库中进行快速、准确的匹配;

b.算法产生的特征点在图像中的密度很大,速度可以达到实时要求;

c.由于 SIFT 特征描述子是向量的形式,它可以与其他形式的特征向量进行联合。

  Sift在一定程度上可以解决:

n  目标的旋转、缩放、平移(RST)

n  图像仿射/投影变换(视点viewpoint)

n  光照影响(illumination)

n  目标遮挡(occlusion)

n  杂物场景(clutter)

n  噪声

 

二、sift特征提取步骤

SIFT算法的实质可以归为在不同尺度空间上查找特征点(关键点)的问题。SIFT算法实现物体识别主要有三大工序:

1.      构造DOG尺度空间

2.      关键点搜索与精确定位

3.      方向赋值

4.      关键点描述

5.      特征匹配

 

三、DOG尺度空间构造

在讲解DOG尺度空间构造之前,我们先弄明白以下几个问题:

1.      什么是尺度空间

从认知学的角度讲,在一幅图像中,即使对一个事物没有概念,或者并不熟悉它,人仍然能够感知此物体的结构。对于计算机而言,想要得知图像中哪些是有意义的,必须先要明确这样一个问题:在一幅图像中,只有在一定的尺度范围内,一个物体才有意义。举一个例子,树枝这个概念,只有在几厘米到几米的距离去观察它,才能感知到它的确是树枝;如果在微米级或者千米级去观察,就不能感知到树枝这个概念了,这样的话可以感知到的是细胞或者是森林的概念。

因而,尺度空间中各尺度图像的模糊程度逐渐变大,能够模拟人在距离目标由近到远时目标在视网膜上的形成过程。并且尺度越大,图像越模糊。

注意,高斯核是唯一可以产生多尺度空间的核函数

 

2.      哪些点是SIFT中要查找的关键点(特征点)

所谓关键点,就是在不同尺度空间的图像下检测出的具有方向信息的局部极值点。这些点不会因光照条件的改变而消失,比如角点、边缘点、暗区域的亮点以及亮区域的暗点,既然两幅图像中有相同的景物,那么使用某种方法分别提取各自的稳定点,这些点之间会有相互对应的匹配点。

 

3. 多尺度和多分辨率

尺度空间表达和金字塔多分辨率表达之间最大的不同是:

l  尺度空间表达是由不同高斯核平滑卷积得到,在所有尺度上有相同的分辨率;

l  而金字塔多分辨率表达每层分辨率减少固定比率。

所以,金字塔多分辨率生成较快,且占用存储空间少;而多尺度表达随着尺度参数的增加冗余信息也变多。多尺度表达的优点在于图像的局部特征可以用简单的形式在不同尺度上描述;而金字塔表达没有理论基础,难以分析图像局部特征。

通过上面的图可以形象的感觉到多尺度表示和多分辨率表示的区别。多尺度在所有尺度上像素数是相同的,细节通过平滑逐步丢失。多分辨率(金字塔)是对图像进行降采样(通常是水平,竖直方向1/2),从而得到一系列尺寸缩小的图像(上图是对缩小图像插值放大的结果)

 

DOG高斯差分金字塔的构建分为两部:高斯金字塔构建和差分图像计算。

(1)高斯金子塔构建

多次进行如下步骤:

a.      对图像做高斯滤波平滑

b.      对图像做将采样

这样,一幅图像可以产生O组(octave)图像,一组图像包括S层(interval)图像。

参数设置啥的,可以看github上的config.cfg文件。具体实现:

a.给定一张图片,按照组数Oi,定义每组图片的长和宽,δ为缩放因子

b.将原始图像用双线性插值的方式缩小到指定大小

c.将缩放后的图像转换为灰度图

d.然后根据组内层数Si,指定高斯模糊核函数

为了加速计算,采用一维高斯核函数分别对图像的行和列进行卷积,达到高斯模糊效果,具体实现见代码。

缩放的代码:

ScaleSpace::ScaleSpace(const Mat32f& mat, int num_octave, int num_scale):
	noctave(num_octave), nscale(num_scale),
	origw(mat.width()), origh(mat.height())
{
	// #pragma omp parallel for schedule(dynamic)
	REP(i, noctave) {
		if (!i)
			pyramids.emplace_back(mat, nscale);//拷入某组图像O,及该组图像的层数S
		else {
			float factor = pow(SCALE_FACTOR, -i);
			int neww = ceil(origw * factor),
					newh = ceil(origh * factor);
			Mat32f resized(newh, neww, 3);
			resize(mat, resized);				//双线性插值,两组图像之间
			pyramids.emplace_back(resized, nscale);
		}
	}
}

resize双线性插值代码,自己github上看吧。结果:

原图:

第一组:

第二组:

第三组:

一维卷积高斯模糊代码:

/*该构造函数在ScaleSpace函数内,由emplace_back阶段发起调用*/
GaussianPyramid::GaussianPyramid(const Mat32f& m, int num_scale):
	nscale(num_scale),
	data(num_scale), mag(num_scale), ort(num_scale),	/*为每个vector指定成员数量?*/
	w(m.width()), h(m.height())
{
	TotalTimer tm("build pyramid");
	if (m.channels() == 3)
		data[0] = rgb2grey(m);
	else
		data[0] = m.clone();

	Filter blurer(nscale, GAUSS_SIGMA, SCALE_FACTOR);	//组内层间进行不同种类的高斯模糊,不再做图像缩放,生成s个高斯核函数
	for (int i = 1; i < nscale; i ++) {
		data[i] = blurer.GaussianBlur(data[0], i);	// sigma needs a better one
	}
}

GaussCache::GaussCache(float sigma) {
	// TODO decide window size ?
	/*
	 *const int kw = round(GAUSS_WINDOW_FACTOR * sigma) + 1;
	 */
	kw = ceil(0.3 * (sigma / 2 - 1) + 0.8) * GAUSS_WINDOW_FACTOR;		//ceil向上取整,高斯核大小(6*sigma +1)*(6*sigma +1),核中央为1,往外依次减小
	//cout << kw << " " << sigma << endl;
	if (kw % 2 == 0) kw ++;
	kernel_buf = create_auto_buf<float>(kw);
	const int center = kw / 2;
	kernel = kernel_buf.get() + center;

	kernel[0] = 1;

	float exp_coeff = -1.0 / (sigma * sigma * 2),
				 wsum = 1;
	for (int i = 1; i <= center; i ++)									//生成一维的高斯核函数,wsum为了归一化
		wsum += (kernel[i] = exp(i * i * exp_coeff)) * 2;

	float fac = 1.0 / wsum;
	kernel[0] = fac;
	for (int i = 1; i <= center; i ++)
		kernel[-i] = (kernel[i] *= fac);
}

Filter::Filter(int nscale,
		float gauss_sigma,
		float scale_factor) {
	REP(k, nscale - 1) {						//每组有s种高斯核函数
		gcache.emplace_back(gauss_sigma);
		gauss_sigma *= scale_factor;
	}
}

// TODO faster convolution,先调用filter.hh中的Mat32f GaussianBlur(const Mat32f& img, int n)
Mat32f Filter::GaussianBlur(
		const Mat32f& img, const GaussCache& gauss) const {
	m_assert(img.channels() == 1);
	TotalTimer tm("gaussianblur");
	const int w = img.width(), h = img.height();
	Mat32f ret(h, w, 1);

	const int kw = gauss.kw;
	const int center = kw / 2;
	float * kernel = gauss.kernel;

	auto cur_line_mem = create_auto_buf<float>(
			center * 2 + max(w, h), true);
	float *cur_line = cur_line_mem.get() + center;

	// apply to columns
	REP(j, w){
		const float* src = img.ptr(0, j);		//某列第一个元素地址
		// copy a column of src
		REP(i, h) {
			cur_line[i] = *src;					//每列的一个元素对应高斯核的中央位置
			src += w;
		}

		// pad the border with border value,构造边界外的图像值,上边界外的全用第一个列元素表示,下边界外的全用最后一个列元素表示
		float v0 = cur_line[0];
		for (int i = 1; i <= center; i ++)
			cur_line[-i] = v0;
		v0 = cur_line[h - 1];
		for (int i = 0; i < center; i ++)
			cur_line[h + i] = v0;

		float *dest = ret.ptr(0, j);
		REP(i, h) {
			float tmp = 0;
			for (int k = -center; k <= center; k ++)			//对列方向做图像和高斯核的卷积
				tmp += kernel[k] * cur_line[i + k];
			*dest = tmp;
			dest += w;
		}
	}

	// apply to rows,原理同上,只是输入图像从dest继续,行方向做卷积
	REP(i, h) {
		float *dest = ret.ptr(i);
		memcpy(cur_line, dest, sizeof(cur_line[0]) * w);
		{	// pad the border
			float v0 = cur_line[0];
			for (int j = 1; j <= center; j ++)
				cur_line[-j] = v0;
			v0 = cur_line[w - 1];
			for (int j = 0; j < center; j ++)
				cur_line[center + j] = v0;
		}
		REP(j, w) {
			float tmp = 0;
			for (int k = -center; k <= center; k ++)
				tmp += kernel[k] * cur_line[j + k];
			*(dest ++) = tmp;
		}
	}
	return ret;
}

刚才有4组分辨率,现在每个分辨率下有5种尺度,供20张图,为了示意,我只贴几组了(按尺度层数递增的顺序)

第二组:

第五组:

(2)高斯差分图像计算

  计算同一组内,相邻两层之间平滑图像之差,并把查的绝对值作为差分图像。

/*计算同一组内相邻两个尺度(高斯模糊)图片的差*/
Mat32f DOGSpace::diff(const Mat32f& img1, const Mat32f& img2) const {
	int w = img1.width(), h = img1.height();
	m_assert(w == img2.width() && h == img2.height());
	Mat32f ret(h, w, 1);
	REP(i, h) {
		// speed up
		const float *p1 = img1.ptr(i),
					      *p2 = img2.ptr(i);
		float* p = ret.ptr(i);
		REP(j, w)
			p[j] = fabs(p1[j] - p2[j]);
	}
	return ret;
}

DOGSpace::DOGSpace(ScaleSpace& ss):
	noctave(ss.noctave), nscale(ss.nscale),
	origw(ss.origw), origh(ss.origh),
	dogs(noctave)
{
#pragma omp parallel for schedule(dynamic)
	REP(i, noctave) {
		auto& o = ss.pyramids[i];
		int ns = o.get_len();	//获取某一组图像内层数S
		REP(j, ns - 1)
			dogs[i].emplace_back(diff(o.get(j), o.get(j+1)));
	}
}

结果同样会很多,我只贴部分结果:

第一组:

第四组:

我们可以通过高斯差分图像看出图像上的像素值变化情况。(如果没有变化,也就没有特征。特征必须是变化尽可能多的点。)DOG图像描绘的是目标的轮廓。

时间: 2015-11-25

SIFT特征--构造DOG尺度空间的相关文章

sift特征--关键点搜索与定位

     前面介绍了如何生成高斯图像金字塔,并计算了每组图像的高斯差分图像.现在介绍如何进行关键点搜索与定位(都在灰度图上搞的). 一.极值点计算    关键点是由DOG空间的局部极值点组成的.为了寻找DoG函数的极值点,每一个像素点要和它所有的相邻点比较,看其是否比它的图像域和尺度域的相邻点大或者小. 每一幅高斯差分图像中的一个像素点,要和它同尺度的8个相邻点和上下相邻尺度对应的9×2个点共26个点比较,以确保在尺度空间和二维图像空间都检测到极值点. 具体算法: (1)   遍历不同分辨率和尺

SIFT特征--方向赋值与关键点描述

  一个SIFT关键点拥有三个信息:位置,尺度和方向.前面已经介绍了如何精确定位关键点的位置,通过尺度不变性求极值点,可以使其具有缩放不变的性质.现在来谈谈为特征点指定方向参数,使提取的特征对图像旋转具有不变性,从而实现匹配时图像的旋转无关性.最后,再介绍该用什么样的描述符来表达sift特征.   一.关键点方向分配   SIFT特征的一个关键的特性是旋转不变性,实现旋转不变的基本思想是采用"相对"的概念.利用关键点邻域像素的梯度方向分布特性,我们可以为每个关键点指定方向参数方向.而后

sift-已经用opencv对2副图片进行了SIFT特征点匹配,这一对一对的特征点存储在哪里呢,网上找的代码

问题描述 已经用opencv对2副图片进行了SIFT特征点匹配,这一对一对的特征点存储在哪里呢,网上找的代码 已经用opencv对2副图片进行了SIFT特征点匹配,那么这一对一对的特征点存储在哪里呢,怎么提取其中一对特征点的坐标? #include "stdafx.h" #include #include #include "opencv2/core/core.hpp"//因为在属性中已经配置了opencv等目录,所以把其当成了本地目录一样 #include &qu

矩阵分解-局部特征(比如gabor或sift特征)的低秩矩阵恢复可行吗?有没有实际意义?

问题描述 局部特征(比如gabor或sift特征)的低秩矩阵恢复可行吗?有没有实际意义? 我想在原始图像的基础上提取局部的特征(入gabor),并降维,然后以这个作为初始字典进行低秩矩阵恢复或低秩矩阵分解,不知道这样会不会比用全局特征效果好?可行吗?

OpenCV教程(47) sift特征和surf特征

     在前面三篇教程中的几种角检测方法,比如harris角检测,都是旋转无关的,即使我们转动图像,依然能检测出角的位置,但是图像缩放后,harris角检测可能会失效,比如下面的图像,图像放大之前可以检测出为harris角,但是图像放大后,则变成了边,不能检测出角了.所以,harris角是缩放相关的.      在paper Distinctive Image Features from Scale-Invariant Keypoints中,D.Lowe提出了SIFT算法,该算法是缩 放无关的

sift

简介 SIFT算法是用来提取图像局部特征的经典算法,SIFT算法的实质是在不同的尺度空间上查找关键点(特征点),并计算出关键点的方向.SIFT所查找到的关键点是一些十分突出,不会因光照,仿射变换和噪音等因素而变化的点,如角点.边缘点.暗区的亮点及亮区的暗点等.  主要步骤 1.构建DOG尺度空间 2.关键点搜索和定位 3.方向赋值 4.关键点描述子的生成 算法详细流程 ----------------------------------------------------------------

OpenCV特征点检测------Surf(特征点篇)

Surf(Speed Up Robust Feature) Surf算法的原理                                                                           1.构建Hessian矩阵构造高斯金字塔尺度空间 其实surf构造的金字塔图像与sift有很大不同,就是因为这些不同才加快了其检测的速度.Sift采用的是DOG图像,而surf采用的是Hessian矩阵行列式近似值图像.Hessian矩阵是Surf算法的核心,为了方便运算

理解SIFT

理解SIFT 理解SIFT MarkdownHTML ## SIFT步骤### 找尺度空间极值 可以用DoG找 需要讲理论依据和DoG结合进行推导 s值和$\sigma$值的选取问题### 关键点定位 和泰勒展开相关### 方向赋值### 局部图像描述符## SIFT用于目标识别的步骤 目标图像提取SIFT特征:搭建SIFT特征数据库 用快速kNN将目标特征在数据库中查找 用霍夫变换寻找属于同一单目标的簇(cluster) 用least-square(最小二乘法)验证## SIFT算法步

SIFT 特征提取算法总结

原文链接:http://www.cnblogs.com/cfantaisie/archive/2011/06/14/2080917.html   主要步骤 1).尺度空间的生成: 2).检测尺度空间极值点: 3).精确定位极值点: 4).为每个关键点指定方向参数: 5).关键点描述子的生成. L(x,y,σ), σ= 1.6 a good tradeoff     D(x,y,σ), σ= 1.6 a good tradeoff 关于尺度空间的理解说明:图中的2是必须的,尺度空间是连续的.在