数字图像处理成长之路13：SIFT之极值点(子像元插值）

很多关于SIFT算法都会写到极值点的定位，这篇文章也 离不开这以问题，然而这篇文章的重点是发现在定位的过程中用到的知识，这也是我写SIFT的目的，我究竟学会了哪些知识，上一篇我初步领略了尺度空间这一全新概念，这一篇我将学到子像素差值或者叫做亚像素差值。
下面这个文章写的不错，
http://www.cnblogs.com/ronny/p/4028776.html
但是我们还是重新自己学习一下。

已知函数,求函数在[0,4]间的最大值。
如果函数是离散的，那么最大值是5，但是函数是连续的，所以我们只能估计当x=1附近时有最大值。

接下来在x=1附近做泰勒展开用一阶二阶差分，近似代替一阶二阶导数，
看这里的二阶导数差分近似，x点的二阶近似由x+1和x+2确定，偏离x点比较大，为了更贴近x点附近的值，我们对二阶导数稍微变换一下。由导数的定义，无论是向前差分还是向后差分，导数不变。

下面用程序求极值：

    double F [5] = {1,5,4,3,2};    double Fdy = F[2]-F[1];    double Fddy = F[2]+F[0]-2*F[1];    double maxY, maxX;    for (double x = 0; x < 4;x += 0.01)    {        double y ;        y = F[1] + Fdy*(x - 1) + Fddy* (x-1)*(x-1)/2;        if (y > maxY)        {            maxY = y;            maxX = x;            mydebug ("x: "<<maxX<<"  "<<"y: "<<maxY);        }    }

打印输出信息：

[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0    y:  3.5[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.01    y:  3.53975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.02    y:  3.579[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.03    y:  3.61775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.04    y:  3.656[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.05    y:  3.69375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.06    y:  3.731[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.07    y:  3.76775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.08    y:  3.804[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.09    y:  3.83975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.1    y:  3.875[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.11    y:  3.90975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.12    y:  3.944[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.13    y:  3.97775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.14    y:  4.011[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.15    y:  4.04375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.16    y:  4.076[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.17    y:  4.10775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.18    y:  4.139[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.19    y:  4.16975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.2    y:  4.2[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.21    y:  4.22975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.22    y:  4.259[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.23    y:  4.28775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.24    y:  4.316[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.25    y:  4.34375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.26    y:  4.371[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.27    y:  4.39775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.28    y:  4.424[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.29    y:  4.44975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.3    y:  4.475[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.31    y:  4.49975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.32    y:  4.524[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.33    y:  4.54775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.34    y:  4.571[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.35    y:  4.59375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.36    y:  4.616[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.37    y:  4.63775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.38    y:  4.659[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.39    y:  4.67975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.4    y:  4.7[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.41    y:  4.71975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.42    y:  4.739[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.43    y:  4.75775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.44    y:  4.776[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.45    y:  4.79375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.46    y:  4.811[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.47    y:  4.82775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.48    y:  4.844[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.49    y:  4.85975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.5    y:  4.875[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.51    y:  4.88975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.52    y:  4.904[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.53    y:  4.91775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.54    y:  4.931[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.55    y:  4.94375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.56    y:  4.956[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.57    y:  4.96775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.58    y:  4.979[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.59    y:  4.98975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.6    y:  5[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.61    y:  5.00975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.62    y:  5.019[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.63    y:  5.02775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.64    y:  5.036[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.65    y:  5.04375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.66    y:  5.051[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.67    y:  5.05775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.68    y:  5.064[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.69    y:  5.06975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.7    y:  5.075[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.71    y:  5.07975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.72    y:  5.084[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.73    y:  5.08775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.74    y:  5.091[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.75    y:  5.09375[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.76    y:  5.096[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.77    y:  5.09775[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.78    y:  5.099[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.79    y:  5.09975[debug- cqalgorithm.cpp - cqTest - 713 ]: x:  0.8    y:  5.1

求得当x= 0.8时，函数最大值为5.1

对于二元函数泰勒展开：

这是雅可比矩阵和黑塞矩阵的和。

写到这里，我大概理解SIFT算法为什么要弄尺度空间这样一个工具，并且强调其连续性，这一切准备都是为了要进行泰勒近似。让我们开看看SIFT的三元泰勒近似，即x,y,和sigma。

如何在编程中实现多元混合导数，即多元差分呢？有一篇文章写得挺好：
http://www.dataguru.cn/thread-744150-1-1.html
我们重新看一下差分，差分的三种方式为：
，向前差分
，向后差分
，中心差分
我们对中心差分进行泰勒展开，只展开到平方项：

当两式相减有：
把这个式子带入任一中心差分公式有：
这些还都比较简单，我们的关键是要找二阶混合偏导数的差分表示，先来看看二元函数的泰勒展开，只展开到平方项：
，仿照一阶差分的做法可得二阶差分中的混合偏导数：

现在回到我们的SIFT中，就可以用计算机来计算从离散的x，y，sigma空间到连续的空间的极值的逼近了。

写到这里大致把SIFT的中心思想总结完了，即类似子像元差值的思想，利用有限的离散数据来逼近或估计连续的数据中的极值，这是我学习SIFT算法最大的收获。

限于时间和精力，现在未能继续用c语言编写SIFT算法剩下的部分，待时间和精力充裕时或许再来完善它，相比于其他特征算法，SIFT还是相对复杂些，弄清楚了SIFT其他算法相比起来就容易一些了。