【ORB-SLAM2源码梳理3】构造特征点图像金字塔，构造灰度质心圆（ORBextractor.cc）

---

特征提取器，ORBextractor构造函数

从Tracking.cc中Tracking构造函数中的构造ORB特征提取器mpORBextractorLeft = new ORBextractor()，进入ORBextractor.cc的ORBextractor构造函数。

mpORBextractorLeft = new ORBextractor(nFeatures,fScaleFactor,nLevels,fIniThFAST,fMinThFAST);

ORBextractor.cc：

ORBextractor::ORBextractor(int _nfeatures,		//指定要提取的特征点数目
						   float _scaleFactor,	//指定图像金字塔的缩放系数，即金字塔层之间的变化尺度
						   int _nlevels,		//指定图像金字塔的层数
						   int _iniThFAST,		//指定初始的FAST关键点提取时计算像素值差值的阈值
						   int _minThFAST):		//如果用上面那个阈值提取不到足够数量的fast点，则使用这个阈值
    nfeatures(_nfeatures), scaleFactor(_scaleFactor), nlevels(_nlevels),
    iniThFAST(_iniThFAST), minThFAST(_minThFAST)//初始化设置这些参数

随后设置金字塔的各项参数

重点1：逐层计算图像金字塔中图像相当于初始图像的缩放系数

   for(int i=1; i<nlevels; i++)  
   {
		// 累乘计算得出缩放因子
        mvScaleFactor[i]=mvScaleFactor[i-1]*scaleFactor;
		//sigma^2就是每层图像相对于初始图像缩放因子的平方
        mvLevelSigma2[i]=mvScaleFactor[i]*mvScaleFactor[i];
    }
    mvInvScaleFactor.resize(nlevels);		// 初始化缩放因子的倒数
    mvInvLevelSigma2.resize(nlevels);		// 初始化缩放因子平方的倒数
	// 获得实际上图像金字塔每一层的变化尺度
	for(int i=0; i<nlevels; i++)
    {
        mvInvScaleFactor[i]=1.0f/mvScaleFactor[i];
        mvInvLevelSigma2[i]=1.0f/mvLevelSigma2[i];
    }

重点2：分配每一层的特征点

来源：知乎@小葡萄

	// 第0层分配的特征点
    float nDesiredFeaturesPerScale = nfeatures*(1 - factor)/(1 - (float)pow((double)factor, (double)nlevels));

    int sumFeatures = 0;
    // 开始逐层计算要分配的特征点数目，除了最后一层
    for( int level = 0; level < nlevels-1; level++ )
    {
        mnFeaturesPerLevel[level] = cvRound(nDesiredFeaturesPerScale);		// mnFeaturesPerLevel[0]为第0层，即原图像
        sumFeatures += mnFeaturesPerLevel[level];							// 累计已经分配的特征点数目
        nDesiredFeaturesPerScale *= factor;									// 按比例取下一层的特征点数目
    }
    mnFeaturesPerLevel[nlevels-1] = std::max(nfeatures - sumFeatures, 0);	// 近似分配，最后把剩余的全给最顶层

重点3：描述子选点模式及灰度质心法计算FAST关键点方向

1. 选点模式（或者理解为：取点的方式）：

pattern是BRIEF描述子选取点对的模式，这个pattern一共有256个点对，即512个点，此处的pattern是由别人设计好的：

    const int npoints = 512;
    const Point* pattern0 = (const Point*)bit_pattern_31_;			// 获取用于计算BRIEF描述子的随机采样点点集头指针
    std::copy(pattern0, pattern0 + npoints, std::back_inserter(pattern));		// 使用std::back_inserter的目的是可以快覆盖掉这个容器pattern之前的数据

注：为了保持踩点的固定，工程上常采用特殊设计的固定pattern来计算BRIEF描述子。

2. 计算FAST关键点的方向

原始的fast关键点没有方向信息，即当图像发生旋转后，brief描述子也会发生变化，故特征点不具有旋转不变性。

解决方法：Orientated FAST，灰度质心法计算特征点方向。

ORB在计算BRIEF描述子时建立的坐标系：以关键点为圆心，以关键点和取点区域形心的连线作为X轴，建立二维坐标系。

一段不太流畅的解释：

传统的FAST关键点：
在当前关键点P周围以一定模式选取N个点对，组合这N个点对的T操作的结果就为最终的描述子。当我们选取点对的时候，是以当前关键点为原点，以水平方向为X轴，以垂直方向为Y轴建立坐标系（上图左）。当图片发生旋转时，坐标系不变，同样的取点模式取出来的点却不一样，计算得到的BRIEF描述子也不一样。

改进后的FAST关键点——加入方向：
我们从图像帧上取出含有关键点P的一小块区域（最外面的矩形），以关键点P为圆心，某一特定的长度为半径画圆，以此作为描述子的取点区域，用积分的方法（下面的公式）计算这不均匀（指的是灰度像素值不均匀）圆形区域的质心Q。以关键点为圆心，以关键点和取点区域形心的连线作为X轴，建立二维坐标系，这样有了参考方向固定的参考系，即使图像发生了旋转，在同一选点模式下也会得到一样BRIEF描述子。

计算方向的具体步骤如下（R为选点圆形区域的半径）：

亦或者：

来源：
CSDN@金木炎
CSDN@Mr.Silver

补充：ORB-SLAM2中关于BRIEF的选点

ORB-SLAM2中的BRIEF描述子是一种二进制描述子，总共具有256bit，即32个字节的长度，每位bit为0或1，即描述子是一个256维由0、1组成的向量。
根据一定的点对选取规则选择点对，该选取规则应当使点对与点对之间的相关性最低，换言之，点对与点对之间尽量垂直，并判断该各个点对两个像素点的灰度值大小（例如，像素p和像素q，p>q则取1,否则取0）。

重点4：为了计算关键点方向，做的准备，表述这些的代码有点复杂

    //This is for orientation
    // pre-compute the end of a row in a circular patch
    umax.resize(HALF_PATCH_SIZE + 1);						// +1：当v=0时，也有对应的u值存在。

	/**
	* 展开说说最大行数vmax：
	* 1. vmax先初始化为R*sin45°+1
	* 关于+1的意义：有人说是因为考虑了中间行，其实并没有什么关系。
	* 较为合适的解释是：为了在vmax和vmin在遍历的过程中产生交叉，防止漏掉。不然为什么不把vmax==vmin，还要分别上下取整。
	* 2. 最大行数：单纯是计算过程中的最大行数，而不是这个圆的最大行数。
	* 此时作者利用圆的对称性来快速计算v与u的关系，此时的vmax只是从0-45°这个过程的最大行数。
	*/
    int v, v0, vmax = cvFloor(HALF_PATCH_SIZE * sqrt(2.f) / 2 + 1);		// v,v0为辅助变量，vmax为最大行数（向下取整）
    
    int vmin = cvCeil(HALF_PATCH_SIZE * sqrt(2.f) / 2);		// vmin（向上取整）意义和vmax差不多，
    														// 此处指的是45-90°过程中的最小行数
    const double hp2 = HALF_PATCH_SIZE*HALF_PATCH_SIZE;		// R的平方
    for (v = 0; v <= vmax; ++v)								// 利用圆的方程计算每行像素值对应的u边界（umax）
        umax[v] = cvRound(sqrt(hp2 - v * v));				// v从0到11依次遍历，umax中存储12个数：15、15、15、15、14、14、
        													// 14、13、13、12、11、10

    // Make sure we are symmetric  确保对称性
    /**
    * 下面的是计算v = vmin至HALF_PATCH_SIZE时的umax[v]，那就是v从11到15之间的值，计算得umax[15]-umax[11]依次是3、6、8、9、10
    * 此时，算出了umax中所有16个数值，依次是：15 15 15 15 14 14 14 13 13 12 11 10 9 8 6 3
    * 这里其实是使用了对称的方式计算上四分之一的圆周上的umax，目的也是为了保持严格的对称
    * 如果按照常规的想法做，由于cvRound就会很容易出现不对称的情况
    */
    for (v = HALF_PATCH_SIZE, v0 = 0; v >= vmin; --v)
    {
        while (umax[v0] == umax[v0 + 1])
            ++v0;
        umax[v] = v0;
        ++v0;
    }
}

源码解释：
要以关键点keypoint像素坐标点为圆心、直径为PATCH_SIZE、半径为HALF_PATCH_SIZE的patch圆内计算关键点keypoint的方向。
那如何描述这个patch圆的范围呢？
这里选择的是：存储不同v所对应的的umax来描述这个patch圆的范围，即给定行数v，能得到对应的列数u。