关于图像特征提取的几点小结：
基于颜色的特征提取：RGB颜色空间——HSV颜色空间
基于纹理的特征提取：  LBP（Local Binary
 Pattern，局部二值模式）局部纹理算子；gabor滤波提取全局纹理，忽略表情局部变化
基于形状的特征提取：伦敦特征Hough (http://www.open-open.com/lib/view/open1440832074794.html) 从2D的图像特征提取应用到3D图像的形状特征提取、傅里叶形状描述法
基于空间关系的特征提取：指图像中分割出来的多个目标之间的相互的空间位置或相对方向关系，这些关系也可以分为连接/邻接关系、交叠/重叠关系和包含/包容关系等。

网上虽然出现了很多讲解ORB特征提取和描述的方法，但都不够详尽。为了搞明白到底是怎么回事，只能结合别人的博客和原著对ORB的详细原理做一个研究和学习。哪里有不对的地方，请多多指教

 

1、算法介绍

ORB（Oriented FAST and Rotated BRIEF）是一种快速特征点提取和描述的算法。这个算法是由Ethan Rublee, Vincent Rabaud, Kurt Konolige以及Gary R.Bradski在2011年一篇名为“ORB：An Efficient Alternative to SIFTor SURF”的文章中提出。ORB算法分为两部分，分别是特征点提取和特征点描述。特征提取是由FAST（Features from  Accelerated Segment Test）算法发展来的，特征点描述是根据BRIEF（Binary Robust IndependentElementary Features）特征描述算法改进的。ORB特征是将FAST特征点的检测方法与BRIEF特征描述子结合起来，并在它们原来的基础上做了改进与优化。据说，ORB算法的速度是sift的100倍，是surf的10倍。

1.1 Fast特征提取

ORB算法的特征提取是由FAST算法改进的，这里成为oFAST（FASTKeypoint Orientation）。也就是说，在使用FAST提取出特征点之后，给其定义一个特征点方向，以此来实现特征点的旋转不变形。FAST算法是公认的最快的特征点提取方法。FAST算法提取的特征点非常接近角点类型。oFAST算法如下：



 

图1  FAST特征点判断示意图

步骤一：粗提取。该步能够提取大量的特征点，但是有很大一部分的特征点的质量不高。下面介绍提取方法。从图像中选取一点P，如上图1。我们判断该点是不是特征点的方法是，以P为圆心画一个半径为3pixel的圆。圆周上如果有连续n个像素点的灰度值比P点的灰度值大或者小，则认为P为特征点。一般n设置为12。为了加快特征点的提取，快速排出非特征点，首先检测1、9、5、13位置上的灰度值，如果P是特征点，那么这四个位置上有3个或3个以上的的像素值都大于或者小于P点的灰度值。如果不满足，则直接排出此点。

步骤二：机器学习的方法筛选最优特征点。简单来说就是使用ID3算法训练一个决策树，将特征点圆周上的16个像素输入决策树中，以此来筛选出最优的FAST特征点。

步骤三：非极大值抑制去除局部较密集特征点。使用非极大值抑制算法去除临近位置多个特征点的问题。为每一个特征点计算出其响应大小。计算方式是特征点P和其周围16个特征点偏差的绝对值和。在比较临近的特征点中，保留响应值较大的特征点，删除其余的特征点。

步骤四：特征点的尺度不变形。建立金字塔，来实现特征点的多尺度不变性。设置一个比例因子scaleFactor（opencv默认为1.2）和金字塔的层数nlevels（pencv默认为8）。将原图像按比例因子缩小成nlevels幅图像。缩放后的图像为：I’= I/scaleFactork(k=1,2,…, nlevels)。nlevels幅不同比例的图像提取特征点总和作为这幅图像的oFAST特征点。

步骤五：特征点的旋转不变性。ORB算法提出使用矩（moment）法来确定FAST特征点的方向。也就是说通过矩来计算特征点以r为半径范围内的质心，特征点坐标到质心形成一个向量作为该特征点的方向。矩定义如下：



 

其中，I(x,y)为图像灰度表达式。该矩的质心为：



 

假设角点坐标为O，则向量的角度即为该特征点的方向。计算公式如下：



 

1.2 rBRIEF特征描述

rBRIEF特征描述是在BRIEF特征描述的基础上加入旋转因子改进的。下面先介绍BRIEF特征提取方法，然后说一说是怎么在此基础上修改的。

BRIEF算法描述

BRIEF算法计算出来的是一个二进制串的特征描述符。它是在一个特征点的邻域内，选择n对像素点pi、qi（i=1,2,…,n）。然后比较每个点对的灰度值的大小。如果I(pi)> I(qi)，则生成二进制串中的1，否则为0。所有的点对都进行比较，则生成长度为n的二进制串。一般n取128、256或512，opencv默认为256。另外，值得注意的是为了增加特征描述符的抗噪性，算法首先需要对图像进行高斯平滑处理。在ORB算法中，在这个地方进行了改进，在使用高斯函数进行平滑后，又用了其他操作，使其更加的具有抗噪性。具体方法下面将会描述。

关于在特征点SxS的区域内选取点对的方法，BRIEF论文(附件2)中测试了5种方法：

1）在图像块内平均采样；

2）p和q都符合(0,S2/25)的高斯分布；

3）p符合(0,S2/25)的高斯分布，而q符合(0,S2/100)的高斯分布；

4）在空间量化极坐标下的离散位置随机采样；

5）把p固定为(0,0)，q在周围平均采样。

五种采样方法的示意图如下：



论文指出，第二种方法可以取得较好的匹配结果。在旋转不是非常厉害的图像里，用BRIEF生成的描述子的匹配质量非常高，作者测试的大多数情况中都超越了SURF。但在旋转大于30°后，BRIEF的匹配率快速降到0左右。BRIEF的耗时非常短，在相同情形下计算512个特征点的描述子时，SURF耗时335ms,BRIEF仅8.18ms；匹配SURF描述子需28.3ms，BRIEF仅需2.19ms。在要求不太高的情形下，BRIEF描述子更容易做到实时。

改进BRIEF算法—rBRIEF（Rotation-AwareBrief）

（1）steered BRIEF（旋转不变性改进）

在使用oFast算法计算出的特征点中包括了特征点的方向角度。假设原始的BRIEF算法在特征点SxS（一般S取31）邻域内选取n对点集。



 

经过旋转角度θ旋转，得到新的点对



 

在新的点集位置上比较点对的大小形成二进制串的描述符。这里需要注意的是，在使用oFast算法是在不同的尺度上提取的特征点。因此，在使用BRIEF特征描述时，要将图像转换到相应的尺度图像上，然后在尺度图像上的特征点处取SxS邻域，然后选择点对并旋转，得到二进制串描述符。

（2）rBRIEF-改进特征点描述子的相关性

使用steeredBRIEF方法得到的特征描述子具有旋转不变性，但是却在另外一个性质上不如原始的BRIEF算法。是什么性质呢，是描述符的可区分性，或者说是相关性。这个性质对特征匹配的好坏影响非常大。描述子是特征点性质的描述。描述子表达了特征点不同于其他特征点的区别。我们计算的描述子要尽量的表达特征点的独特性。如果不同特征点的描述子的可区分性比较差，匹配时不容易找到对应的匹配点，引起误匹配。ORB论文中，作者用不同的方法对100k个特征点计算二进制描述符，对这些描述符进行统计，如下表所示：



图2 特征描述子的均值分布.X轴代表距离均值0.5的距离，y轴是相应均值下的特征点数量统计

我们先不看rBRIEF的分布。对BRIEF和steeredBRIEF两种算法的比较可知，BRIEF算法落在0上的特征点数较多，因此BRIEF算法计算的描述符的均值在0.5左右，每个描述符的方差较大，可区分性较强。而steeredBRIEF失去了这个特性。至于为什么均值在0.5左右，方差较大，可区分性较强的原因，这里大概分析一下。这里的描述子是二进制串，里面的数值不是0就是1，如果二进制串的均值在0.5左右的话，那么这个串有大约相同数目的0和1，那么方差就较大了。用统计的观点来分析二进制串的区分性，如果两个二进制串的均值都比0.5大很多，那么说明这两个二进制串中都有较多的1时，在这两个串的相同位置同时出现1的概率就会很高。那么这两个特征点的描述子就有很大的相似性。这就增大了描述符之间的相关性，减小之案件的可区分性。

下面我们介绍解决上面这个问题的方法：rBRIEF。

原始的BRIEF算法有5中去点对的方法，原文作者使用了方法2。为了解决描述子的可区分性和相关性的问题，ORB论文中没有使用5种方法中的任意一种，而是使用统计学习的方法来重新选择点对集合。

首先建立300k个特征点测试集。对于测试集中的每个点，考虑其31x31邻域。这里不同于原始BRIEF算法的地方是，这里在对图像进行高斯平滑之后，使用邻域中的某个点的5x5邻域灰度平均值来代替某个点对的值，进而比较点对的大小。这样特征值更加具备抗噪性。另外可以使用积分图像加快求取5x5邻域灰度平均值的速度。

从上面可知，在31x31的邻域内共有(31-5+1)x(31-5+1)=729个这样的子窗口，那么取点对的方法共有M=265356种，我们就要在这Ｍ种方法中选取256种取法，选择的原则是这256种取法之间的相关性最小。怎么选取呢？

1）在300k特征点的每个31x31邻域内按M种方法取点对，比较点对大小，形成一个300kxM的二进制矩阵Q。矩阵的每一列代表300k个点按某种取法得到的二进制数。

2）对Q矩阵的每一列求取平均值，按照平均值到0.5的距离大小重新对Q矩阵的列向量排序，形成矩阵T。

3）将T的第一列向量放到R中。

4）取T的下一列向量和R中的所有列向量计算相关性，如果相关系数小于设定的阈值，则将T中的该列向量移至R中。

5）按照4）的方式不断进行操作，直到R中的向量数量为256。

通过这种方法就选取了这256种取点对的方法。这就是rBRIEF算法。

 

2、ORB特征提取实验

实验代码用opencv中的

2.1ORB特征提取和匹配实验



 

（1）



 

（2）



（3）



 

（4）

2-1 ORB特征匹配

从上图（1）（2）（3）可以看出，ORB算法的特征匹配效果比较理想，并且具有较稳定的旋转不变性。但是通过（4）看出，ORB算法在尺度方面效果较差，在增加算法的尺度变换的情况下仍然没有取得较好的结果。

ORB是一种快速的特征提取和匹配的算法。它的速度非常快，但是相应的算法的质量较差。和sift相比，ORB使用二进制串作为特征描述，这就造成了高的误匹配率。

 

SIFT（Scale-invariant feature transform）是一种检测局部特征的算法，该算法通过求一幅图中的特征点（interest points,or corner points）及其有关scale 和 orientation 的描述子得到特征并进行图像特征点匹配，获得了良好效果，详细解析如下：
算法描述
SIFT特征不只具有尺度不变性，即使改变旋转角度，图像亮度或拍摄视角，仍然能够得到好的检测效果。整个算法分为以下几个部分：
1. 构建尺度空间
这是一个初始化操作，尺度空间理论目的是模拟图像数据的多尺度特征。
高斯卷积核是实现尺度变换的唯一线性核，于是一副二维图像的尺度空间定义为：

其中 G(x,y,σ) 是尺度可变高斯函数 
（x，y）是空间坐标，是尺度坐标。σ大小决定图像的平滑程度，大尺度对应图像的概貌特征，小尺度对应图像的细节特征。大的σ值对应粗糙尺度(低分辨率)，反之，对应精细尺度(高分辨率)。为了有效的在尺度空间检测到稳定的关键点，提出了高斯差分尺度空间（DOG scale-space）。利用不同尺度的高斯差分核与图像卷积生成。
下图所示不同σ下图像尺度空间：


关于尺度空间的理解说明：2kσ中的2是必须的，尺度空间是连续的。在  Lowe的论文中 ，将第0层的初始尺度定为1.6（最模糊），图片的初始尺度定为0.5（最清晰）. 在检测极值点前对原始图像的高斯平滑以致图像丢失高频信息，所以 Lowe 建议在建立尺度空间前首先对原始图像长宽扩展一倍，以保留原始图像信息，增加特征点数量。尺度越大图像越模糊。 

图像金字塔的建立：对于一幅图像I,建立其在不同尺度(scale)的图像，也成为子八度（octave），这是为了scale-invariant，也就是在任何尺度都能够有对应的特征点，第一个子八度的scale为原图大小，后面每个octave为上一个octave降采样的结果，即原图的1/4（长宽分别减半），构成下一个子八度（高一层金字塔）。


尺度空间的所有取值，i为octave的塔数（第几个塔），s为每塔层数
由图片size决定建几个塔，每塔几层图像(S一般为3-5层)。0塔的第0层是原始图像(或你double后的图像)，往上每一层是对其下一层进行Laplacian变换（高斯卷积，其中σ值渐大，例如可以是σ, k*σ, k*k*σ…），直观上看来越往上图片越模糊。塔间的图片是降采样关系，例如1塔的第0层可以由0塔的第3层down sample得到，然后进行与0塔类似的高斯卷积操作。

2. LoG近似DoG找到关键点<检测DOG尺度空间极值点>
为了寻找尺度空间的极值点，每一个采样点要和它所有的相邻点比较，看其是否比它的图像域和尺度域的相邻点大或者小。如图所示，中间的检测点和它同尺度的8个相邻点和上下相邻尺度对应的9×2个点共26个点比较，以确保在尺度空间和二维图像空间都检测到极值点。 一个点如果在DOG尺度空间本层以及上下两层的26个领域中是最大或最小值时，就认为该点是图像在该尺度下的一个特征点,如图所示。

同一组中的相邻尺度（由于k的取值关系，肯定是上下层）之间进行寻找

s=3的情况
 在极值比较的过程中，每一组图像的首末两层是无法进行极值比较的，为了满足尺度变化的连续性（下面有详解）
，我们在每一组图像的顶层继续用高斯模糊生成了 3 幅图像，高斯金字塔有每组S+3层图像。DOG金字塔每组有S+2层图像.

==========================================
这里有的童鞋不理解什么叫“为了满足尺度变化的连续性”，现在做仔细阐述：
假设s=3，也就是每个塔里有3层，则k=21/s=21/3，那么按照上图可得Gauss Space和DoG space 分别有3个（s个）和2个（s-1个）分量，在DoG space中，1st-octave两项分别是σ,kσ; 2nd-octave两项分别是2σ,2kσ;由于无法比较极值，我们必须在高斯空间继续添加高斯模糊项，使得形成σ,kσ,k2σ,k3σ,k4σ这样就可以选择DoG space中的中间三项kσ,k2σ,k3σ（只有左右都有才能有极值），那么下一octave中（由上一层降采样获得）所得三项即为2kσ,2k2σ,2k3σ，其首项2kσ=24/3。刚好与上一octave末项k3σ=23/3尺度变化连续起来，所以每次要在Gaussian space添加3项，每组（塔）共S+3层图像，相应的DoG金字塔有S+2层图像。
==========================================

使用Laplacian of Gaussian能够很好地找到找到图像中的兴趣点，但是需要大量的计算量，所以使用Difference of Gaussian图像的极大极小值近似寻找特征点.DOG算子计算简单，是尺度归一化的LoG算子的近似,有关DOG寻找特征点的介绍及方法详见http://blog.csdn.net/abcjennifer/article/details/7639488，极值点检测用的Non-Maximal Suppression。

3. 除去不好的特征点
这一步本质上要去掉DoG局部曲率非常不对称的像素。

通过拟和三维二次函数以精确确定关键点的位置和尺度（达到亚像素精度），同时去除低对比度的关键点和不稳定的边缘响应点(因为DoG算子会产生较强的边缘响应)，以增强匹配稳定性、提高抗噪声能力，在这里使用近似Harris Corner检测器。
①空间尺度函数泰勒展开式如下：，对上式求导,并令其为0,得到精确的位置, 得
②在已经检测到的特征点中,要去掉低对比度的特征点和不稳定的边缘响应点。去除低对比度的点：把公式(2)代入公式(1)，即在DoG Space的极值点处D(x)取值，只取前两项可得：

若   ，该特征点就保留下来，否则丢弃。
③边缘响应的去除
一个定义不好的高斯差分算子的极值在横跨边缘的地方有较大的主曲率，而在垂直边缘的方向有较小的主曲率。主曲率通过一个2×2 的Hessian矩阵H求出:
导数由采样点相邻差估计得到。
D的主曲率和H的特征值成正比，令α为较大特征值，β为较小的特征值，则
令α=γβ，则
 (r + 1)2/r的值在两个特征值相等的时候最小，随着r的增大而增大，因此，为了检测主曲率是否在某域值r下，只需检测
if (α+β)/ αβ> (r+1)2/r, throw it out.   在Lowe的文章中，取r＝10。


4. 给特征点赋值一个128维方向参数
上一步中确定了每幅图中的特征点，为每个特征点计算一个方向，依照这个方向做进一步的计算， 利用关键点邻域像素的梯度方向分布特性为每个关键点指定方向参数，使算子具备旋转不变性。

为(x,y)处梯度的模值和方向公式。其中L所用的尺度为每个关键点各自所在的尺度。至此，图像的关键点已经检测完毕，每个关键点有三个信息：位置，所处尺度、方向，由此可以确定一个SIFT特征区域。

梯度直方图的范围是0～360度，其中每10度一个柱，总共36个柱。随着距
      中心点越远的领域其对直方图的贡献也响应减小.Lowe论文中还提到要使用高斯函数对直方图进行平滑，减少突变的影响。

在实际计算时，我们在以关键点为中心的邻域窗口内采样，并用直方图统计邻域像素的梯度方向。梯度直方图的范围是0～360度，其中每45度一个柱，总共8个柱, 或者每10度一个柱，总共36个柱。Lowe论文中还提到要使用高斯函数对直方图进行平滑，减少突变的影响。直方图的峰值则代表了该关键点处邻域梯度的主方向，即作为该关键点的方向。

直方图中的峰值就是主方向，其他的达到最大值80%的方向可作为辅助方向
由梯度方向直方图确定主梯度方向
该步中将建立所有scale中特征点的描述子（128维）
Identify peak and assign orientation and sum of magnitude to key point.
  The user may choose a threshold to exclude key points based on their assigned sum of magnitudes.


关键点描述子的生成步骤

 通过对关键点周围图像区域分块，计算块内梯度直方图，生成具有独特性的向量，这个向量是该区域图像信息的一种抽象，具有唯一性。

5. 关键点描述子的生成
首先将坐标轴旋转为关键点的方向，以确保旋转不变性。以关键点为中心取8×8的窗口。
Figure.16*16的图中其中1/4的特征点梯度方向及scale，右图为其加权到8个主方向后的效果。
图左部分的中央为当前关键点的位置，每个小格代表关键点邻域所在尺度空间的一个像素，利用公式求得每个像素的梯度幅值与梯度方向，箭头方向代表该像素的梯度方向，箭头长度代表梯度模值，然后用高斯窗口对其进行加权运算。

图中蓝色的圈代表高斯加权的范围（越靠近关键点的像素梯度方向信息贡献越大）。然后在每4×4的小块上计算8个方向的梯度方向直方图，绘制每个梯度方向的累加值，即可形成一个种子点，如图右部分示。此图中一个关键点由2×2共4个种子点组成，每个种子点有8个方向向量信息。这种邻域方向性信息联合的思想增强了算法抗噪声的能力，同时对于含有定位误差的特征匹配也提供了较好的容错性。

计算keypoint周围的16*16的window中每一个像素的梯度，而且使用高斯下降函数降低远离中心的权重。
在每个4*4的1/16象限中，通过加权梯度值加到直方图8个方向区间中的一个，计算出一个梯度方向直方图。
这样就可以对每个feature形成一个4*4*8=128维的描述子，每一维都可以表示4*4个格子中一个的scale/orientation. 将这个向量归一化之后，就进一步去除了光照的影响。
5. 根据SIFT进行Match
生成了A、B两幅图的描述子，（分别是k1*128维和k2*128维），就将两图中各个scale（所有scale）的描述子进行匹配，匹配上128维即可表示两个特征点match上了。

实际计算过程中，为了增强匹配的稳健性，Lowe建议对每个关键点使用4×4共16个种子点来描述，这样对于一个关键点就可以产生128个数据，即最终形成128维的SIFT特征向量。此时SIFT特征向量已经去除了尺度变化、旋转等几何变形因素的影响，再继续将特征向量的长度归一化，则可以进一步去除光照变化的影响。 当两幅图像的SIFT特征向量生成后，下一步我们采用关键点特征向量的欧式距离来作为两幅图像中关键点的相似性判定度量。取图像1中的某个关键点，并找出其与图像2中欧式距离最近的前两个关键点，在这两个关键点中，如果最近的距离除以次近的距离少于某个比例阈值，则接受这一对匹配点。降低这个比例阈值，SIFT匹配点数目会减少，但更加稳定。为了排除因为图像遮挡和背景混乱而产生的无匹配关系的关键点,Lowe提出了比较最近邻距离与次近邻距离的方法,距离比率ratio小于某个阈值的认为是正确匹配。因为对于错误匹配,由于特征空间的高维性,相似的距离可能有大量其他的错误匹配,从而它的ratio值比较高。Lowe推荐ratio的阈值为0.8。但作者对大量任意存在尺度、旋转和亮度变化的两幅图片进行匹配，结果表明ratio取值在0. 4~0. 6之间最佳，小于0. 4的很少有匹配点，大于0. 6的则存在大量错误匹配点。(如果这个地方你要改进，最好给出一个匹配率和ration之间的关系图，这样才有说服力)作者建议ratio的取值原则如下:
ratio=0. 4　对于准确度要求高的匹配；
ratio=0. 6　对于匹配点数目要求比较多的匹配； 
ratio=0. 5　一般情况下。
也可按如下原则:当最近邻距离<200时ratio=0. 6，反之ratio=0. 4。ratio的取值策略能排分错误匹配点。


当两幅图像的SIFT特征向量生成后，下一步我们采用关键点特征向量的欧式距离来作为两幅图像中关键点的相似性判定度量。取图像1中的某个关键点，并找出其与图像2中欧式距离最近的前两个关键点，在这两个关键点中，如果最近的距离除以次近的距离少于某个比例阈值，则接受这一对匹配点。降低这个比例阈值，SIFT匹配点数目会减少，但更加稳定。
 
实验结果：






Python+opencv实现：

import cv2
import numpy as np
#import pdb
#pdb.set_trace()#turn on the pdb prompt

#read image
img = cv2.imread('D:\privacy\picture\little girl.jpg',cv2.IMREAD_COLOR)
gray = cv2.cvtColor(img,cv2.COLOR_BGR2GRAY)
cv2.imshow('origin',img);

#SIFT
detector = cv2.SIFT()
keypoints = detector.detect(gray,None)
img = cv2.drawKeypoints(gray,keypoints)
#img = cv2.drawKeypoints(gray,keypoints,flags = cv2.DRAW_MATCHES_FLAGS_DRAW_RICH_KEYPOINTS)
cv2.imshow('test',img);
cv2.waitKey(0)
cv2.destroyAllWindows()



C实现：

// FeatureDetector.cpp : Defines the entry point for the console application.
//  
//  Created by Rachel on 14-1-12.  
//  Copyright (c) 2013年 ZJU. All rights reserved.  
//  

#include "stdafx.h"
#include "highgui.h"
#include "cv.h"
#include "vector"
#include "opencv\cxcore.hpp"
#include "iostream"
#include "opencv.hpp"
#include "nonfree.hpp"
#include "showhelper.h"

using namespace cv;
using namespace std;

int _tmain(int argc, _TCHAR* argv[])
{
	//Load Image 
	Mat c_src1 =  imread( "..\\Images\\3.jpg");
	Mat c_src2 = imread("..\\Images\\4.jpg");
	Mat src1 = imread( "..\\Images\\3.jpg", CV_LOAD_IMAGE_GRAYSCALE);
	Mat src2 = imread( "..\\Images\\4.jpg", CV_LOAD_IMAGE_GRAYSCALE);
	if( !src1.data || !src2.data )
	{ std::cout<< " --(!) Error reading images " << std::endl; return -1; }

	//sift feature detect
	SiftFeatureDetector detector;
	std::vector<KeyPoint> kp1, kp2;

	detector.detect( src1, kp1 );
	detector.detect( src2, kp2 );
	SiftDescriptorExtractor extractor;
	Mat des1,des2;//descriptor
	extractor.compute(src1,kp1,des1);
	extractor.compute(src2,kp2,des2);	
	Mat res1,res2; 
	int drawmode = DrawMatchesFlags::DRAW_RICH_KEYPOINTS;
	drawKeypoints(c_src1,kp1,res1,Scalar::all(-1),drawmode);//在内存中画出特征点
	drawKeypoints(c_src2,kp2,res2,Scalar::all(-1),drawmode);
	cout<<"size of description of Img1: "<<kp1.size()<<endl;
	cout<<"size of description of Img2: "<<kp2.size()<<endl;

	//write the size of features on picture
	CvFont font;    
	double hScale=1;   
	double vScale=1;    
	int lineWidth=2;// 相当于写字的线条    
	cvInitFont(&font,CV_FONT_HERSHEY_SIMPLEX|CV_FONT_ITALIC, hScale,vScale,0,lineWidth);//初始化字体，准备写到图片上的   
	// cvPoint 为起笔的x，y坐标   
	IplImage* transimg1 = cvCloneImage(&(IplImage) res1);
	IplImage* transimg2 = cvCloneImage(&(IplImage) res2);
		
	char str1[20],str2[20];
	sprintf(str1,"%d",kp1.size());
	sprintf(str2,"%d",kp2.size());


	const char* str = str1;
	cvPutText(transimg1,str1,cvPoint(280,230),&font,CV_RGB(255,0,0));//在图片中输出字符 

	str = str2;
	cvPutText(transimg2,str2,cvPoint(280,230),&font,CV_RGB(255,0,0));//在图片中输出字符 

	//imshow("Description 1",res1);
	cvShowImage("descriptor1",transimg1);
	cvShowImage("descriptor2",transimg2);

	BFMatcher matcher(NORM_L2);
	vector<DMatch> matches;
	matcher.match(des1,des2,matches);
	Mat img_match;
	drawMatches(src1,kp1,src2,kp2,matches,img_match);//,Scalar::all(-1),Scalar::all(-1),vector<char>(),drawmode);
	cout<<"number of matched points: "<<matches.size()<<endl;
	imshow("matches",img_match);
	cvWaitKey();
	cvDestroyAllWindows();

	return 0;
}





===============================
基本概念及一些补充
什么是局部特征？
　　•局部特征从总体上说是图像或在视觉领域中一些有别于其周围的地方
　　•局部特征通常是描述一块区域，使其能具有高可区分度
　　•局部特征的好坏直接会决定着后面分类、识别是否会得到一个好的结果
 局部特征需具备的特性
　　•重复性
　　•可区分性
　　•准确性
　　•数量以及效率
　　•不变性
 局部特征提取算法-sift
　　•SIFT算法由D.G.Lowe 1999年提出，2004年完善总结。后来Y.Ke将其描述子部分用PCA代替直方图的方式，对其进行改进。
 　　•SIFT算法是一种提取局部特征的算法，在尺度空间寻找极值点，提取位置，尺度，旋转不变量
　　•SIFT特征是图像的局部特征，其对旋转、尺度缩放、亮度变化保持不变性，对视角变化、仿射变换、噪声也保持一定程度的稳定性。
　　•独特性好，信息量丰富，适用于在海量特征数据库中进行快速、准确的匹配。
　　•多量性，即使少数的几个物体也可以产生大量SIFT特征向量。
　　•可扩展性，可以很方便的与其他形式的特征向量进行联合。
尺度空间理论
　　•尺度空间理论目的是模拟图像数据的多尺度特征
　　•其基本思想是在视觉信息图像信息处理模型中引入一个被视为尺度的参数, 通过连续变化尺度参数获得不同尺度下的视觉处理信息, 然后综合这些信息以深入地挖掘图像的本质特征。
描述子生成的细节
　　•以极值点为中心点，并且以此点所处于的高斯尺度sigma值作为半径因子。对于远离中心点的梯度值降低对其所处区域的直方图的贡献，防止一些突变的影响。
　　•每个极值点对其进行三线性插值，这样可以把此极值点的贡献均衡的分到直方图中相邻的柱子上
归一化处理
　　•在求出4*4*8的128维特征向量后，此时SIFT特征向量已经去除了尺度变化、旋转等几何变形因素的影响。而图像的对比度变化相当于每个像素点乘上一个因子，光照变化是每个像素点加上一个值，但这些对图像归一化的梯度没有影响。因此将特征向量的长度归一化，则可以进一步去除光照变化的影响。
　　•对于一些非线性的光照变化，SIFT并不具备不变性，但由于这类变化影响的主要是梯度的幅值变化，对梯度的方向影响较小，因此作者通过限制梯度幅值的值来减少这类变化造成的影响。
PCA-SIFT算法
　　•PCA-SIFT与标准SIFT有相同的亚像素位置，尺度和主方向。但在第4步计算描述子的设计，采用的主成分分析的技术。
　　•下面介绍一下其特征描述子计算的部分：
　　　　•用特征点周围的41×41的像斑计算它的主元，并用PCA-SIFT将原来的2×39×39维的向量降成20维，以达到更精确的表示方式。
　　　　•它的主要步骤为，对每一个关键点：在关键点周围提取一个41×41的像斑于给定的尺度，旋转到它的主方向 ；计算39×39水平和垂直的梯度，形成一个大小为3042的矢量；用预先计算好的投影矩阵n×3042与此矢量相乘；这样生成一个大小为n的PCA-SIFT描述子。
===============================
辅助资料：




===============================

Reference:
Lowe SIFT 原文：http://www.cs.ubc.ca/~lowe/papers/ijcv04.pdf
SIFT 的C实现：https://github.com/robwhess/opensift/blob/master/src
MATLAB 应用Sift算子的模式识别方法：http://blog.csdn.net/abcjennifer/article/details/7372880
http://blog.csdn.net/abcjennifer/article/details/7365882
http://en.wikipedia.org/wiki/Scale-invariant_feature_transform#David_Lowe.27s_method
http://blog.sciencenet.cn/blog-613779-475881.html
http://www.cnblogs.com/linyunzju/archive/2011/06/14/2080950.html
http://www.cnblogs.com/linyunzju/archive/2011/06/14/2080951.html
http://blog.csdn.net/ijuliet/article/details/4640624
http://www.cnblogs.com/cfantaisie/archive/2011/06/14/2080917.html  (部分图片有误，以本文中的图片为准)




关于computer vision的更多讨论与交流，敬请关注本博客和新浪微博Rachel____Zhang。

颜色特征提取
颜色直方图
一般颜色直方图    cv2.calcHist
全局累加直方图      np.cumsum()
主色调直方图
总结
颜色矩
颜色聚合向量
颜色相关图


图像特征是指可以对图像的特点或内容进行表征的一系列属性的集合，主要包括图像自然特征（如亮度、色彩、纹理等）和图像人为特征（如图像频谱、图像直方图等）。


图像特征提取根据其相对尺度可分为全局特征提取和局部特征提取两类。全局特征提取关注图像的整体表征。常见的全局特征包括颜色特征、纹理特征、形状特征、空间位置关系特征等。局部特征提取关注图像的某个局部区域的特殊性质。一幅图像中往往包含若干兴趣区域，从这些区域中可以提取数量不等的若干个局部特征。


特征特点


代表性（可区分性）
稳定性
独立性


颜色特征是比较简单但是应用较广泛的一种视觉特征。颜色特征往往和图像中包含的对象或场景相关。与其他图像特征相比，颜色特征对图像的尺寸、方向、视角变化的依赖性较小，即相对于图像的尺寸、方向、视角变化具有较好的健壮性。


目前使用的颜色特征主要包括颜色直方图、颜色矩、颜色集、颜色聚合向量以及颜色相关图。

颜色直方图

颜色直方图用于描述图像中像素颜色的数值分布情况，可以反映图像颜色的统计分布和图像的基本色调。颜色直方图仅可表征图像中某一颜色值出现的频数，无法描述图像像素分布的空间位置信息。


任意一幅图像都能唯一给出一幅与它对应的颜色直方图，但不同的图像可能有相同的颜色直方图，因此直方图与图像存在一对多的关系。

一般颜色直方图    cv2.calcHist

颜色直方图是最基本的颜色特征，它反映的是图像中像素颜色值的组成分布，即出现了哪些颜色以及各种颜色出现的概率。假设s(xi)为图像F中某一特定颜色xi的像素个数，则xi像素出现的频率为：



整个图像F的一般颜色直方图可以表示为：H(F)=[h(x1),h(x}xb2),…,h(xn)]


方法说明


cv2.calcHist(image,channels,mask,histSize,ranges)

参数说明：

1、image：原始图像，需要[]括起来

2、channel：通道编号，需要[]括起来

3、mask：掩膜图像，当统计全局图像，值设为None，统计局部图像的直方图，需要使用到它

4、histSize：表示直方图的灰度级数，例如[0,255]一共256级

5、range：表示像素值范围[a,b]，表示a到b

6、accumulate：表示累计标识，默认False，如果设为True，则直方图在开始计算时不会被清零，计算的是多个直方图的累积结果，用于对一组图像计算直方图。参数可选

import matplotlib.pyplot as plt
import cv2 as cv
from skimage import data
img = data.coffee()
hist0 = cv.calcHist([img],[0],None,[256],[0,255])
hist1 = cv.calcHist([img],[1],None,[256],[0,255])
hist2 = cv.calcHist([img],[2],None,[256],[0,255])


全局累加直方图      np.cumsum()

当图像中的颜色值不能取遍所有可能的颜色值时，一般颜色直方图中就会出现一些零值。这些零值的出现会影响相似性的度量，进而会使计算出的相似度不能准确反映图像之间的颜色分布差异。

为了弥合一般颜色直方图法的上述缺陷，在一般直方图方法的基础之上，通过对直方图元素进行累加，消除零值影响，形成全局累加直方图。



img = data.coffee()
hist0 = np.cumsum(cv.calcHist([img],[0],None,[256],[0,255]))
hist1 = np.cumsum(cv.calcHist([img],[1],None,[256],[0,255]))
hist2 = np.cumsum(cv.calcHist([img],[2],None,[256],[0,255]))
plt.figure(figsize=(6,6))
plt.plot(range(256),hist0,label = 'R')
plt.plot(range(256),hist1,label = 'G')
plt.plot(range(256),hist2,label = 'B')
plt.legend()                 
plt.title("累积直方图")


主色调直方图

在一幅图像中，不同颜色值出现的概率不尽相同，且通常情况下少数几种颜色就能涵盖整幅图像的主色调。基于该思想，主色调直方图法会计算出图像中每种颜色出现的频率，选择出现频率最高的几种颜色并将其作为主色调。

总结

优点


计算简单，对图像的平移和旋转变换不敏感，能简单描述图像中颜色的全局分布情况。


缺点


无法捕捉（即会丢失）颜色组成之间的空间位置关系。

颜色矩

矩是非常重要的统计量，用于表征数据分布的特点。在统计中，一阶矩表示数据分布的均值，二阶矩表示数据分布的方差，三阶矩表示数据分布的偏移度。


图像的颜色矩用于对图像内的颜色分布进行表征，是比较重要的一种全局图像特征表示。图像的颜色矩特征提取时主要瞄准图像颜色矩中的一阶矩、二阶矩和三阶矩，对于图像而言，这3种统计特征已经足以表达数字图像的颜色分布。相对于颜色直方图特征提取，颜色矩特征提取的优点是无须对颜色特征进行提前量化。


原理


对于数字图像P，其一阶颜色矩的定义为：



其中Pij表示数字图像P的第i个图像通道的第j个像素的像素值，N表示图像中像素的个数。


二阶颜色矩的定义为：




三阶颜色矩的定义为：




其中一阶矩可以表征该颜色通道的平均响应强度，二阶矩可以表示该颜色通道的响应方差，三阶矩可以表征该颜色通道数据分布的偏移度。针对彩色图像，图像的颜色矩一共有9个分量，每个颜色通道均有3个低阶矩。

def color_comment(img):
    r,g,b = cv.split(img)
    color_featrue = []
    # 一阶矩
    r_mean = np.mean(r)
    g_mean = np.mean(g)
    b_mean = np.mean(b)
    # 二阶矩
    r_std = np.std(r)
    g_std = np.std(g)
    b_std = np.std(b)
    #三阶矩
    r_offset = (np.mean(np.abs((r - r_mean)**3)))**(1./3)
    g_offset = (np.mean(np.abs((g - g_mean)**3)))**(1./3)
    b_offset = (np.mean(np.abs((b - b_mean)**3)))**(1./3)
    color_featrue.extend([r_mean,g_mean,b_mean,r_std,g_std,b_std,r_offset,g_offset,b_offset])
    return color_featrue
print(color_comment(img))

[158.5690875, 85.794025, 51.48475, 62.9728671221504, 60.958103707650785, 52.93569362069574, 73.2849737717244, 69.97427088321201, 67.60422322535385]

颜色矩的特点：图像的颜色矩有9个分量（3个颜色通道，每个通道上3个低阶矩）；与其他颜色特征相比非常简洁；分辨力较弱；颜色矩一般和其他特征结合使用，可以起到缩小范围的作用。

颜色聚合向量
参考

针对颜色直方图和颜色矩无法表达图像色彩的空间位置的缺点，有人提出了图像的颜色聚合向量。


其核心思想是将属于颜色直方图的每个颜色量化区间的像素分为两部分，如果该颜色量化区间中的某些像素占据的连续区域的面积大于指定阈值，则将该区域内的像素作为聚合像素，否则作为非聚合像素。


颜色聚合向量可表示为<(α1,β1),(αn,βn)>，其中αi与βi分别代表颜色直方图的第i个颜色量化区间中的聚合像素和非聚合像素的数量。颜色聚合向量除了包含颜色频率信息外，也包含颜色的部分空间分布信息，因此其可以获得比颜色直方图更好的表示效果。


步骤

1、量化

颜色聚合向量算法的第一步与求普通的颜色直方图类似，即对图像进行量化处理。一般采用均匀量化处理方法，量化的目标是使图像中只保留有限个颜色区间。

2、连通区域划分

针对重新量化后的像素值矩阵，根据像素间的连通性把图像划分成若干个连通区域。

3、判断聚合性

统计每个连通区域中的像素数目，根据设定的阈值判断该区域中的像素是聚合的，还是非聚合，得出每个颜色区间中聚合像素和非聚合像素的数量αi和βi。

4、 聚合向量形成

图像的聚合向量可以表示为〈(α1,β1),(αn,βn)〉

颜色相关图

颜色相关图是图像颜色分布的另外一种表达方式。颜色相关图不仅可以显示像素在图像中的占比，也可以反映不同颜色对间的空间位置的相关性。颜色相关图利用颜色对间的相对距离分布来描述空间位置信息。


原理


颜色相关图是一张用颜色对<i,j>索引的表，其中<i,j>的第k个分量表示颜色为c(i)的像素和颜色为c(j)的像素之间的距离等于k的概率。设I表示整张图像的全部像素，Ic(i)表示颜色为c(i)的所有像素，则图像的颜色相关图可以表达为：



|p1-p2|表示像素p1和p2之间的距离

def is_vaild(X,Y,point): #判断像素分布点是否超出图像范围，超出返回False
    if point[0] < 0 or point[0] >= X:
        return False
    if point[1] < 0 or point[1] >= Y:
        return False
    return True
    
def getNeighbors(X,Y,x,y,dist):  # 输入图片的一个像素点的位置，返回它的8邻域
    cn1 = (x+dist,y+dist)
    cn2 = (x+dist,y)
    cn3 = (x+dist,y-dist)
    cn4 = (x,y-dist)
    cn5 = (x-dist,y-dist)
    cn6 = (x-dist,y)
    cn7 = (x-dist,y+dist)
    cn8 = (x,y+dist)
    point = (cn1,cn2,cn3,cn4,cn5,cn6,cn7,cn8)
    Cn = []
    for i in point:
        if is_vaild(X,Y,i):
            Cn.append(i)
    return Cn
        
def corrlogram(img,dist):
    xx,yy,tt = img.shape
    cgram = np.zeros((256,256),np.uint8)
    for x in range(xx):
        for y in range(yy):
            for t in range(tt):
                color_i = img[x,y,t]   # X的某一个通道的像素值
                neighbors_i = getNeighbors(xx,yy,x,y,dist)
            for j in neighbors_i:
                j0 = j[0]
                j1 = j[1]
                color_j = img[j0,j1,t]   #X的某一个邻域像素点的某一个通道的像素值
                cgram[color_i,color_j] +=  1  #统计像素值i核像素值j的个数
    return cgram
crgam = corrlogram(img,4)  # 4 
plt.imshow(crgam)
plt.show()