3D点云——NARF特征

NARF: 3D Range Image Features for Object Recognition

NARF 全称 normal aligned radial feature(法线对齐的径向特征) 。是一种3D特征点检测和描述算法。下面将分别详述 NARF 特征点的提取和描述。

NARF 特征点提取

NARF 特征点提取算法的两个目标是：1). 提取出的特征点应位于表面稳定的区域(where the surface is stable)以确保对法线的估计具有较强的鲁邦性，同时在其邻域点的法线应具有明显的变化。2). 提取出的点应尽可能的位于物体的边缘处，从而使得在不同的视点下都能提取到大部分相同的特征点。

特性

1. 能提取到 borders points 和其对应的 shadow points
2. 提取到的特征点都是位于物体的外轮廓(边缘)上，更加具有稳定性和可重复性
3. 对特征点Normal的估计只使用共面点，因而对Normal的估计更加稳定
4. 所有操作都是基于2D的range image，计算量相比于直接操作点云要小

步骤

1. 边缘提取

论文将点分为三种类型：object borders , shadow borders , veil points .

论文提出“提取边缘的最鲁棒的指导原则是：相邻点之间的距离变化(is a change in the distance between neighboring points)” 。因此边缘提取的算法就是基于这个指导原则进行的。

对于每一个 range image 中的点 ${\mathbf{p}}_{\mathbf{i}}$ ：

1. 使用启发式方法寻找不跨越边缘的相邻点的典型3D距离
2. 使用 (1) 得到的距离信息，计算该点是边缘点的可信度(score)
3. 如果是边缘点，将该点分为上述3类点中的一类
4. 非极大抑制以获得精确的边缘点位置

实现细节 ：

1. 使用启发式方法寻找不跨越边缘的相邻点的典型3D距离

   1. 对每一个 range image 中的点 $p_i$ ，选择其2D平面内的一个方形窗内的邻域点，记为 N p i = { n 1 , n ) 2 , . . . , n s 2 } N_{p_i} = \{n_1, n)2, ..., n_{s^2}\} Npi​​={n1​,n)2,...,ns2​} ，其中 $s$ 为窗宽

   2. 计算 $p_i$ 与 $N_{p_i}$ 内每个点的3D距离(欧氏距离)，得到 D p i = { d 0 , d 1 , . . . , d s 2 } D_{p_i} = \{d_0, d_1, ..., d_{s^2}\} Dpi​​={d0​,d1​,...,ds2​} ，对 $D_{p_i}$ 按升序排列，得到 D p i ′ = { d 0 ′ , d 1 ′ , . . . , d s 2 ′ } D_{p_i}^{'} = \{d_0^{'}, d_1^{'}, ..., d_{s^2}^{'}\} Dpi​′​={d0′​,d1′​,...,ds2′​}

   3. 假设 $N_{p_i}$ 中至少有 $M$ 个点与 $p_i$ 位于同一表面，我们选取 $\delta =d_M^{{}^{\prime }}$ 作为 $p_i$ 的 $p_i-s$ 邻域内的、不考虑邻域内非共面点的典型3D距离(typical 3D distance)

      Note ：关于 $M$ 和 $s$ 的选取，论文选择 $s=5$， $M=(\frac{s+1}{2}{)}^2$ 。为了自适应，可以选择距离 $d_j$ 跳变的点作为 $M$

2. 使用 (1) 得到的距离信息，计算该点是边缘点的可信度(score)

   对点 $p_i$ ，计算四个分数值，分别表示点的上、左、右、下方存在边缘的可能(下面仅以 right 方向为例，其他方向计算方法相似)：
   notation : $p_{x,y}$ 表示图片中 $(x,y)$ 位置处的点

   1. 首先计算右边点的平均距离
      $$p_{right}=\frac{1}{m_p}\sum _{i=1}^{m_p}{p}_{x+i,y}\text{\hspace{1em}(1)}$$
      其中 $m_p$ 为考虑的右边点的个数，论文中采用 $m_p=3$

   2. 计算3D距离 $d_{right}=||p_{x,y}-p_{right}||$

   3. 计算 $s_{right}\in [0,1)$
      $$s_{right}=\max (0,1-\frac{\delta }{d_{right}})\text{\hspace{1em}(2)}$$
      $s_{right}$ 值越大，表明点 $p_i$ 右边距离变化越明显，则存在边缘的概率越大。

   4. 对得到的 $S_{right}$ 图进行平滑操作，以得到连续的边缘同时提高对噪声的鲁棒性

3. 如果是边缘点，将该点分为上述3类点中的一类

   1. 确定 $p_{x,y}$ 所属的点的类型：根据 $p_{x,y}$ 与 $p_{right}$ 的3D距离的比较，$p_{x,y}$ 距离更小则表明可能是 obstacle border 否则可能是 shadow border 。
   2. 对所有可能是 obstacle border 的点，在2D平面内搜索其右侧对应的 shadow border 。选择2D平面上右边 $m_p$ 范围内3D距离最大的点作为对应的 shadow border 并计算 $s_{shadow}$ 。对 $s_{right}$ 做如下调整：
      $$s_{right}^{{}^{\prime }}=\max (0.9,1-(1-s_{shadow}{)}^3)\times s_{right}\text{\hspace{1em}(3)}$$

4. 非极大抑制以获得精确的边缘点位置

   如果 $s_{right}^{{}^{\prime }}$ 大于阈值(论文中取 0.8) ，并且是 $(p_{x-1,y},p_{x+1,y})$ 邻域内的极大值，则认为 $p_{x,y}$ 为 obstacle border ，(4.6) 中找到的对应的 shadow point 即其对应的 shadow border ，位于两者之间的点就是 veil point

2. 特征点提取

特征点提取的指导原则如下：

1. 必须考虑边缘和表面结构信息
2. 必须是能在不同视角下都能被检测到的点
3. 必须位于稳定的区域，从而可以得到的稳定的法向量

步骤如下：

1. 对 range image 中的每一个点，计算一个表示其邻域表面变化的分数和这种变化的一个主方向(dominant direction)
2. 对每个点，根据邻域点的dominant direction 计算一个 interest value，该值表达了两个性质 1). 邻域内点的dominant directions 之间的差异；2). 该点所处的表面的局部变化量(表面的稳定程度) 。
3. 对 interest value 图进行平滑
4. 非极大抑制，找到最终的 interest points

实现细节 ：

由 边缘提取 步骤，我们将 range image 中的点分成了两类：borders point 和 others point 。

notation ： $p_{x,y}$ 表示 range image 中 $(x,y)$ 处的像素点

1. 估计法向量 ：对每个点 $p_i\in RangeImage$ ，利用其2D邻域内欧氏距离在 $2\delta$ 内的点的3D位置和PCA 算法，估计其法向量 。

2. 对 range image 中的每一个点，计算一个表示其邻域表面变化的分数和这种变化的一个主方向(dominant direction)

   1. 计算主方向 $v$ 和幅值 $\lambda$ ：对于 borders point ，$v$ 为边缘方向；对于 others point ，$v$ 为主曲率方向(main direction of the curvature)，$\lambda$ 为方向的幅值。

   2. 对每个 $v$ ，都赋一个权值 $w\in [0,1)$ ：对 borders point ， $w=1$ ；对 others point ，$w=1-(1-\lambda {)}^3$ 。

   3. 指定一个参数 $\sigma$ 表示 support size 。对图像中的每个点 $p$ ，找到其所有欧氏距离在 $\frac{\sigma }{2}$ 内且中间不存在边缘(两点位于同一表面上)的点，记为 N p = { p 0 , p 1 , . . . , p N } N_p = \{p_0, p_1, ..., p_N\} Np​={p0​,p1​,...,pN​} 。

   4. 对 $\forall p_i\in N_p$ ，将 $v_{p_i}$ 投影到一个平面 $P$ 上 ($P$ 为垂直于 $p$ 到原点连线的一个平面) ，得到一个角度 ${\alpha }_{p_i}$ 。对 ${\alpha }_{p_i}$ 作如下调整

      $${\alpha }^{{}^{\prime }}=\{\begin{array}{l}2\times (\alpha -180)for\alpha >90\\ 2\times (\alpha +180)for\alpha \le -90\end{array}\text{\hspace{1em}(4)}$$

      使得 ${\alpha }^{{}^{\prime }}\in (-180,180]$ 。

   5. 使用 bounded Gaussian kernel 对 $w$ 和 $\alpha$ 进行平滑

   6. 计算 interest value：

      $$I_1(p)=\underset{i}{\min }(1-w_{p_i}\times \max (1-\frac{10\times ||p-p_i||}{\sigma },0))\text{\hspace{1em}(5)}$$

      $$f(p_i)=\sqrt{w_{p_i}\times (1-|\frac{2\times ||p-p_i||}{\sigma }-\frac{1}{2}|)}\text{\hspace{1em}(6)}$$

      KaTeX parse error: Got function '\max' as argument to '\underset' at position 24: … \underset{i,j}\̲m̲a̲x̲(f(p_i)f(p_j)(1…

      $$I(p)=I_1(p)I_2(p)\text{\hspace{1em}(8)}$$

3. 非极大抑制，找到最终的 interest points

NARF 描述子计算

特性

1. 具有 6DOF 的不变性
2. 使用距离函数即可进行特征之间的比较，匹配简单
3. 可以根据需要灵活设置描述子的大小，且描述子的大小具有明确的物理意义

步骤

计算 NARF 描述子的指导原则如下：

1. Captures the existence of occupied and free space, so that parts on the surface and also the outer shape of an object can be described.
2. 对噪声具有鲁棒性
3. 能在特征点周围提取一个唯一的局部坐标系

步骤如下：

1. 计算 interest point $p$ 附近的法向量对齐的 range value patch：即视角为沿着 $p$ 的法向量的local range image
2. 将一个星形模板(star pattern)覆盖在这个patch上，模板的每一个beam(梁、束)都是最终NARF描述子中的一个值，这个值表示了local range image 中沿着该beam方向的像素变化了多少
3. 从(2)中的描述子中提取一个方向(类似于 SIFT 中的主方向)，记为 d
4. 将描述子旋转(3) d 方向，使得最终的描述子具有旋转不变性

实现细节 ：

1. 计算Local range image

   1. 建立 Local coordinate system：原点在特征点 $p$ 上，z 轴指向Normal的方向，y轴指向世界坐标系的 upright(垂直的，直立的) 方向，x轴由右手准则定义。
   2. 将点云中和 $p$ 点距离在 $\frac{\sigma }{2}$ 内的点都转换到 local 坐标系下表达。转换后的 (x, y)表示该点落在描述子的哪一个cell中，而所有转换后落在该cell中的点的的z值的最小值则是该cell的值。没有点落入的cell则被置为最大值 $\frac{\sigma }{2}$ 。
   3. 论文中选取的 Local range value patch 大小为 10x10
   4. 为了避免在patch中点云分辨率的分布不一致问题，使用插值或 ray tracing 技术来处理cell之间的区域
   5. 使用 Gaussian 模糊对 patch 进行平滑处理

2. 计算NARF-descriptor描述子

   1. 将一个具有 $n$ 个beams 的星形模板覆盖到patch上，$n$ 就是 NARF-descriptor 的大小：
      
   2. 对每一个 beam $b_i$ ，选择patch上所有落在 $b_i$ 上的cell 的值 C = { c 0 , c 1 , . . . , c m } C = \{c_0, c_1, ..., c_m\} C={c0​,c1​,...,cm​} 。其中 $c_0$ 为patch的中心点的值，其余的都按到 $c_0$ 的距离进行升序排列 。
   3. NARF-descriptor 的第 $i$ 个描述子 $D_I$ 的值由如下公式计算：
      
      其中 $w(c_j)$ 是一个基于距离的加权系数。 随着 $c_j$ 到 $c_0$ 的距离的增加，$w(c_j)$ 由2逐渐降低到1。公式(10)中 $m$ 为落在 $b_i$ 上的cell的个数。公式(11)是为了将所有的 $D_i^{{}^{\prime }}$ 归一化到 $[-0.5,0.5]$ 内。
   4. 由图(b)可以看到，在最终的 NARF-descriptor 中，平坦区域的 D 值很小，而跨越边缘区域的 D 值很大

3. 添加对于Normal方向的旋转不变性

   1. 建立一个旋转方向(360°)的离散直方图

   2. 直方图中对应角度 $\beta$ 的bin中存放的值为：

      $$h(\beta )=\frac{1}{2}+\frac{1}{n}\sum _{i=1}^n{D}_i\times (1-\frac{|\beta -{\gamma }_i|}{180}{)}^2\text{\hspace{1em}(12)}$$
      其中 ${\gamma }_i$ 为NARF-descriptor 中第 $i$ 个cell(即 $b_i$ ) 所对应的角度。

   3. 选择直返图中值最大的bin所对应的角度作为 patch 的dominant orientation。如果存在值为最大值 $80\%$ 的bin，也保留下来，并为其创建另外一个特征描述子。

   4. 将NARF-descriptor shift(移动？) dominant orientation ，使之具有对 Normal 方向的旋转不变性。

关于特征比较(匹配) ：

1. 使用距离函数就可以对特征进行比较，论文选择 Manhattan 距离作为距离函数，并将距离除以描述符中cell的个数，使得度量出来的距离位于 $[0,1]$ 内。

个人总结

NARF描述子，整体看下来给我的感觉很像SIFT描述子。尤其是建立方向直方图那一步。通过建立 Normal aligned 的range image patch， 使得描述子具有姿态的不变性，而通过主方向的选取，又使得描述子具有沿着Normal方向旋转的不变性。因此，从理论上来看，NARF描述子具有很好的6DOF的不变性。但是这都是基于一个前提——Normal的估计具有稳定性和准确定。因此将NARF描述子和NARF特征点提取结合使用，应该具有不错的效果。