Rob Hess sift源码详解（二）

最新推荐文章于 2021-06-21 02:23:56 发布
原创 最新推荐文章于 2021-06-21 02:23:56 发布 · 1.4k 阅读 · 1 ·
CC 4.0 BY-SA版权
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
文章标签：

#OpenCV #源码 #sift

本文详细介绍了SIFT算法中关键点检测的过程，包括关键点的粗定位与精确定位，并解释了如何通过插值和Hessian矩阵计算进行关键点位置与强度的优化。

摘要生成于 C知道 ，由 DeepSeek-R1 满血版支持， 前往体验 >

接上回继续，这部分为关键点定位。

static CvSeq* scale_space_extrema( IplImage*** dog_pyr, int octvs, int intvls,
				   double contr_thr, int curv_thr,
				   CvMemStorage* storage )//关键点定位
{
  CvSeq* features;
  double prelim_contr_thr = 0.5 * contr_thr / intvls;
  struct feature* feat;
  struct detection_data* ddata;
  int o, i, r, c;

  features = cvCreateSeq( 0, sizeof(CvSeq), sizeof(struct feature), storage );
  for( o = 0; o < octvs; o++ )
    for( i = 1; i <= intvls; i++ )
      for(r = SIFT_IMG_BORDER; r < dog_pyr[o][0]->height-SIFT_IMG_BORDER; r++)
	for(c = SIFT_IMG_BORDER; c < dog_pyr[o][0]->width-SIFT_IMG_BORDER; c++)//加入一个长度为SIFT_IMG_BORDER的边界，目的是排除边界的影响
	  /* perform preliminary check on contrast */
	  if( ABS( pixval32f( dog_pyr[o][i], r, c ) ) > prelim_contr_thr )//不理解为什么要大于某个阈值，希望知道者赐教
	    if( is_extremum( dog_pyr, o, i, r, c ) )//关键点初始定位
	      {
		feat = interp_extremum(dog_pyr, o, i, r, c, intvls, contr_thr);//关键点精确定位
		if( feat )
		  {
		    ddata = feat_detection_data( feat );
                    //去除边缘响应
		    if( ! is_too_edge_like( dog_pyr[ddata->octv][ddata->intvl],
					    ddata->r, ddata->c, curv_thr ) )//阈值判断
		      {
			cvSeqPush( features, feat );
		      }
		    else
		      free( ddata );
		    free( feat );
		  }
	      }
  
  return features;
}
static void interp_step( IplImage*** dog_pyr, int octv, int intvl, int r, int c,double* xi, double* xr, double* xc )
{
  CvMat* dD, * H, * H_inv, X;
  double x[3] = { 0 };
  
  dD = deriv_3D( dog_pyr, octv, intvl, r, c );//求偏导
  H = hessian_3D( dog_pyr, octv, intvl, r, c );//计算Hessian矩阵
  H_inv = cvCreateMat( 3, 3, CV_64FC1 );
  cvInvert( H, H_inv, CV_SVD );//矩阵求逆
  cvInitMatHeader( &X, 3, 1, CV_64FC1, x, CV_AUTOSTEP );
  cvGEMM( H_inv, dD, -1, NULL, 0, &X, 0 );
  
  cvReleaseMat( &dD );
  cvReleaseMat( &H );
  cvReleaseMat( &H_inv );

  *xi = x[2];
  *xr = x[1];
  *xc = x[0];
}
static int is_extremum( IplImage*** dog_pyr, int octv, int intvl, int r, int c )//关键点初始定位
{
  double val = pixval32f( dog_pyr[octv][intvl], r, c );
  int i, j, k;
//为寻找DoG函数的极值点，每一个像素点都要和它所有的相邻点相比较，看其是否比它的图像域和尺度域的相邻点大或者小。包括和它同尺度的8个相邻点和上下相邻尺度对应的9*2个点共计26个点比较
  /* check for maximum */  
if( val > 0 )    
{      
   for( i = -1; i <= 1; i++ )
      for( j = -1; j <= 1; j++ )  
         for( k = -1; k <= 1; k++ )    
             if( val < pixval32f( dog_pyr[octv][intvl+i], r + j, c + k ) )      
                return 0;    
}  
/* check for minimum */  
else    
{      
   for( i = -1; i <= 1; i++ )
      for( j = -1; j <= 1; j++ )  
          for( k = -1; k <= 1; k++ )    
              if( val > pixval32f( dog_pyr[octv][intvl+i], r + j, c + k ) )      
                 return 0;   
} 
 return 1;
}

//关键点精确定位
static struct feature* interp_extremum( IplImage*** dog_pyr, int octv,
					int intvl, int r, int c, int intvls,
					double contr_thr )
{
  struct feature* feat;
  struct detection_data* ddata;
  double xi, xr, xc, contr;
  int i = 0;
  
  while( i < SIFT_MAX_INTERP_STEPS )//lowe选择迭代SIFT_MAX_INTERP_STEPS次

    {
      interp_step( dog_pyr, octv, intvl, r, c, &xi, &xr, &xc );//通过三维二次函数来精确定位
      if( ABS( xi ) < 0.5  &&  ABS( xr ) < 0.5  &&  ABS( xc ) < 0.5 )
      //当相对插值中心的偏移量在任一维度大于0.5时代表插值中心已经偏移到其邻近点上了，故在新位置上反复插值直到收敛
	break;
      
      //四舍五入叠加到原位置上
      c += cvRound( xc );//cols      
      r += cvRound( xr );//rows      
      intvl += cvRound( xi );//intvls            
      if( intvl < 1  ||  
          intvl > intvls  ||  
          c < SIFT_IMG_BORDER  ||  
          r < SIFT_IMG_BORDER  ||  
          c >= dog_pyr[octv][0]->width - SIFT_IMG_BORDER  ||  
          r >= dog_pyr[octv][0]->height - SIFT_IMG_BORDER )
          {  
             return NULL;
          }           
       i++;    
    }    
    /* ensure convergence of interpolation */  
    if( i >= SIFT_MAX_INTERP_STEPS )    
       return NULL;    
    contr = interp_contr( dog_pyr, octv, intvl, r, c, xi, xr, xc );  
    if( ABS( contr ) < contr_thr / intvls )//过小的点易受噪声的干扰而变得不稳定，所以将其小于contr_thr / intvls的极值点删除。 
       return NULL;  
    feat = new_feature();  
    ddata = feat_detection_data( feat );  
    //获取特征点的精确位置（原位置上加上拟合的偏移量）以及尺度
    feat->img_pt.x = feat->x = ( c + xc ) * pow( 2.0, octv );  
    feat->img_pt.y = feat->y = ( r + xr ) * pow( 2.0, octv );  
    ddata->r = r;  
    ddata->c = c;  
    ddata->octv = octv;  
    ddata->intvl = intvl;  
    ddata->subintvl = xi;  
    return feat;
}
static CvMat* deriv_3D( IplImage*** dog_pyr, int octv, int intvl, int r, int c )
{
  CvMat* dI;
  double dx, dy, ds;

  //通过有限差分求导代替微分方程中的微分，将连续变量离散化
  dx = ( pixval32f( dog_pyr[octv][intvl], r, c+1 ) -
	 pixval32f( dog_pyr[octv][intvl], r, c-1 ) ) / 2.0;
  dy = ( pixval32f( dog_pyr[octv][intvl], r+1, c ) -
	 pixval32f( dog_pyr[octv][intvl], r-1, c ) ) / 2.0;
  ds = ( pixval32f( dog_pyr[octv][intvl+1], r, c ) -
	 pixval32f( dog_pyr[octv][intvl-1], r, c ) ) / 2.0;
  
  dI = cvCreateMat( 3, 1, CV_64FC1 );
  cvmSet( dI, 0, 0, dx );
  cvmSet( dI, 1, 0, dy );
  cvmSet( dI, 2, 0, ds );


  return dI;
}
static CvMat* hessian_3D( IplImage*** dog_pyr, int octv, int intvl, int r,int c )
{
  CvMat* H;
  double v, dxx, dyy, dss, dxy, dxs, dys;
  
  v = pixval32f( dog_pyr[octv][intvl], r, c );
  dxx = ( pixval32f( dog_pyr[octv][intvl], r, c+1 ) + 
	  pixval32f( dog_pyr[octv][intvl], r, c-1 ) - 2 * v );
  dyy = ( pixval32f( dog_pyr[octv][intvl], r+1, c ) +
	  pixval32f( dog_pyr[octv][intvl], r-1, c ) - 2 * v );
  dss = ( pixval32f( dog_pyr[octv][intvl+1], r, c ) +
	  pixval32f( dog_pyr[octv][intvl-1], r, c ) - 2 * v );
  dxy = ( pixval32f( dog_pyr[octv][intvl], r+1, c+1 ) -
	  pixval32f( dog_pyr[octv][intvl], r+1, c-1 ) -
	  pixval32f( dog_pyr[octv][intvl], r-1, c+1 ) +
	  pixval32f( dog_pyr[octv][intvl], r-1, c-1 ) ) / 4.0;
  dxs = ( pixval32f( dog_pyr[octv][intvl+1], r, c+1 ) -
	  pixval32f( dog_pyr[octv][intvl+1], r, c-1 ) -
	  pixval32f( dog_pyr[octv][intvl-1], r, c+1 ) +
	  pixval32f( dog_pyr[octv][intvl-1], r, c-1 ) ) / 4.0;
  dys = ( pixval32f( dog_pyr[octv][intvl+1], r+1, c ) -
	  pixval32f( dog_pyr[octv][intvl+1], r-1, c ) -
	  pixval32f( dog_pyr[octv][intvl-1], r+1, c ) +
	  pixval32f( dog_pyr[octv][intvl-1], r-1, c ) ) / 4.0;
  
  H = cvCreateMat( 3, 3, CV_64FC1 );
  cvmSet( H, 0, 0, dxx );
  cvmSet( H, 0, 1, dxy );
  cvmSet( H, 0, 2, dxs );
  cvmSet( H, 1, 0, dxy );
  cvmSet( H, 1, 1, dyy );
  cvmSet( H, 1, 2, dys );
  cvmSet( H, 2, 0, dxs );
  cvmSet( H, 2, 1, dys );
  cvmSet( H, 2, 2, dss );
  /// Ixx  Ixy  Ixs \ 
  //| Ixy  Iyy  Iys | 
  //\ Ixs  Iys  Iss /
  return H;
}
static double interp_contr( IplImage*** dog_pyr, int octv, int intvl, int r,
			    int c, double xi, double xr, double xc )
{
  CvMat* dD, X, T;
  double t[1], x[3] = { xc, xr, xi };

  cvInitMatHeader( &X, 3, 1, CV_64FC1, x, CV_AUTOSTEP );
  cvInitMatHeader( &T, 1, 1, CV_64FC1, t, CV_AUTOSTEP );
  dD = deriv_3D( dog_pyr, octv, intvl, r, c );
  cvGEMM( dD, &X, 1, NULL, 0, &T,  CV_GEMM_A_T );
  cvReleaseMat( &dD );

  return pixval32f( dog_pyr[octv][intvl], r, c ) + t[0] * 0.5;
}
static struct feature* new_feature( void )
{
  struct feature* feat;
  struct detection_data* ddata;

  feat = malloc( sizeof( struct feature ) );
  memset( feat, 0, sizeof( struct feature ) );
  ddata = malloc( sizeof( struct detection_data ) );
  memset( ddata, 0, sizeof( struct detection_data ) );
  feat->feature_data = ddata;
  feat->type = FEATURE_LOWE;

  return feat;
}
static int is_too_edge_like( IplImage* dog_img, int r, int c, int curv_thr )//消除边缘响应
{
  double d, dxx, dyy, dxy, tr, det;

  /* principal curvatures are computed using the trace and det of Hessian */
  d = pixval32f(dog_img, r, c);
  dxx = pixval32f( dog_img, r, c+1 ) + pixval32f( dog_img, r, c-1 ) - 2 * d;
  dyy = pixval32f( dog_img, r+1, c ) + pixval32f( dog_img, r-1, c ) - 2 * d;
  dxy = ( pixval32f(dog_img, r+1, c+1) - pixval32f(dog_img, r+1, c-1) -
	  pixval32f(dog_img, r-1, c+1) + pixval32f(dog_img, r-1, c-1) ) / 4.0;
  tr = dxx + dyy;//矩阵H的对角线元素之和
  det = dxx * dyy - dxy * dxy;//矩阵H的行列式

  /* negative determinant -> curvatures have different signs; reject feature */
  if( det <= 0 )
    return 1;

  if( tr * tr / det < ( curv_thr + 1.0 )*( curv_thr + 1.0 ) / curv_thr )//检测主曲率是否小于阈值curv_thr

    return 0;
  return 1;
}
未完待续

Rob Hess sift源码详解（一）
低调做人高调做事
10-23 2727
网上已经有许多关于sift算法详细解释。本文适合于已经对sift算法有一定了解的读者，写本文的目的在于对Rob Hess的sift库中的c++源代码进行注释，从代码的角度对提升对sift算法的理解。闲话少说，进入正文。 首先从网址http://robwhess.github.io/opensift/上下载sift源代码。解压之后不能直接运行，因为作者是在linux下编译的，所以需要自己在wind
OpenCV2马拉松第27圈——SIFT论文，原理及源码解读
热门推荐
abcd1992719g的专栏
06-07 1万+
计算机视觉讨论群162501053 转载请注明：http://blog.youkuaiyun.com/abcd1992719g/article/details/28424579
RobHess的SIFT代码解析步骤二
weixin_30670925的博客
06-03 300
平台：win10 x64 +VS 2015专业版 +opencv-2.4.11 + gtk_-bundle_2.24.10_win32 主要参考：1.代码：RobHess的SIFT源码 2.书：王永明 王贵锦 《图像局部不变性特征与描述》 SIFT四步骤和特征匹配及筛选： 步骤一：建立尺度空间，即建立高斯差分(DoG)金字塔dog_pyr 步骤二：在尺度空间中检测极值点，并进行精确...
OpenCV-跟我一起学数字图像处理之拉普拉斯算子
weixin_30721077的博客
09-26 1470
Laplace算子和Sobel算子一样，属于空间锐化滤波操作。起本质与前面的Spatial Filter操作大同小异，下面就通过Laplace算子来介绍一下空间锐化滤波，并对OpenCV中提供的Laplacian函数进行一些说明。 数学原理 离散函数导数 离散函数的导数退化成了差分，一维一阶差分公式和二阶差分公式分别为， Laplace算子的差分形式 分别...
图像处理中的一阶偏导和二阶偏导
xiaoxin_ling的专栏
02-11 7772
1. 一阶差分： 2. 二阶偏导数的推导和近似： 3. 上式以点（i+1，j）为中心，用i代换i+1可得以（i,j）为中心的二阶偏导数则有： 4. 同理： 5. 进而可推导： 6. 这样我们就可以很好的运用其他的一阶偏导的定义，如SIFT特征OpenCV实现版本中的一阶以及二阶偏导： [cpp] view plaincopypri
图像处理中的一阶偏导数和二阶偏导数
槑宝包
04-15 4337
1. 一阶差分：   2. 二阶偏导数的推导和近似：   3. 上式以点（i+1，j）为中心，用i代换i+1可得以（i,j）为中心的二阶偏导数则有：   4. 同理：   5. 进而可推导：   6. 这样我们就可以很好的运用其他的一阶偏导的定义，如SIFT特征OpenCV实现版本中的一阶以及二阶偏
sift算法c语言实现
weixin_34112208的博客
03-01 218
前段时间在做三维測量方面的研究。须要得到物体表面三维数据。sift算法是立体匹配中的经典算法。以下是对RobHess的SIFT源码的凝视。部分内容參考网上，在这里向各位大神表示感谢。 http://blog.youkuaiyun.com/lsh_2013/article/details/46826141 头文件及函数声明 #include "sift.h" #include "i...
Rob Hess sift源码详解（三）
低调做人高调做事
10-23 1082
这部分为关键点分配方向和计算尺度 static void calc_feature_scales( CvSeq* features, double sigma, int intvls )//计算特征点尺度 { struct feature* feat; struct detection_data* ddata; double intvl; int i, n; n = f
SIFT算法实现理解及注释详解(基于RobHess源码).pdf
06-19
在Rob Hess的SIFT算法实现中，主要包括以下几个步骤： 1. **图像初始化**： 首先，将输入的图像转换为32位浮点灰度图像，并通过三次插值放大图像，以便在不同的尺度上检测特征。接着，应用高斯滤波器对图像进行...
SIFT算法详解与RobHess源码分析
"SIFT算法实现理解及注释详解(基于RobHess源码).pdf" SIFT（尺度不变特征变换）算法是一种经典的图像处理技术，由David G. Lowe在1999年提出，主要用于图像匹配和物体识别。由于其对尺度、旋转和亮度变化的不变性，...
Rob Hess的SIFT算法
08-26
**Rob Hess的SIFT算法详解** SIFT（Scale-Invariant Feature Transform，尺度不变特征变换）是一种在计算机视觉领域中广泛使用的图像特征检测方法，由David G. Lowe在1999年提出。Rob Hess对SIFT算法进行了进一步的...
如何用计算机编程实现泰勒展开,如何在opencv中对图像进行泰勒级数展开？
weixin_39918128的博客
06-21 315
static void interp_step( IplImage*** dog_pyr, int octv, int intvl, int r, int c,double* xi, double* xr, double* xc ){CvMat* dD, * H, * H_inv, X;double x[3] = { 0 };dD = deriv_3D( dog_pyr, octv, intvl,...
尺度的概念和sift算法的解释
lihuajie1003的专栏
11-21 3356
首先了解一下多尺度的概念。 尺度空间也就是分辨率，就好比我们看一个物体，在不同的距离看到的大小不一样，也就是从不同的尺度空间观察物体。       如果想要描述现实世界的结构，或者将三维物体映射到二维的图像上去，多尺度表示将会至关重要。多尺度表示的概念很容易理解，举例说明，绘制地图时会有比例尺的概念。世界地图中就只能够显示大洲大洋，以及较大的地域和国家；而一个城市地图，甚至可以详细的显
Rob Hess的SIFT程序详细解释
zhuoyue08的专栏
11-12 4018
最近也一直在研究SIFT（尺度不变特征匹配算法），先理解算法的数学原理，然后再一步步分析Rob Hess的程序代码： /** default number of sampled intervals per octave sub-level 每个尺度空间的采样间隔数*/ //初始值为3 #define SIFT_INTVLS 3 /** default sigma for initial g
SIFT算法原理详解
monk1992的博客
05-07 3286
通过《图像局部不变性特征与描述》学习SIFT，遇到各种Issue，总结了这篇博客和另外九篇博客。感谢关注，希望可以互相学习，不断提升。转载请注明链接：https://www.cnblogs....
Opencv--图像处理之一阶和二阶偏导数
Jensen Lee的博客
04-02 4369
 1. 一阶差分：   2. 二阶偏导数的推导和近似：   3. 上式以点（i+1，j）为中心，用i代换i+1可得以（i,j）为中心的二阶偏导数则有：   4. 同理：   5. 进而可推导：   6. 这样我们就可以很好的运用其他的一阶偏导的定义，如SIFT特征OpenCV实现版本中的一阶以及二阶偏导： ``` /* Comp...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值