三维计算机视觉（五）--特征描述子

PPF机器人视觉中有一项重要人物就是从场景中提取物体的位置，姿态。图像处理算法借助Deep Learning 的东风已经在图像的物体标记领域耍的飞起了。而从三维场景中提取物体还有待研究。目前已有的思路是先提取关键点，再使用各种局部特征描述子对关键点进行描述，最后与待检测物体进行比对，得到点-点的匹配。个别文章在之后还采取了ICP对匹配结果进行优化。　对于缺乏表面纹理信息，或局部曲率变化很小

Eason.wxd

23662人浏览 · 2017-10-01 11:18:38

Eason.wxd · 2017-10-01 11:18:38 发布

PPF

机器人视觉中有一项重要人物就是从场景中提取物体的位置，姿态。图像处理算法借助Deep Learning 的东风已经在图像的物体标记领域耍的飞起了。而从三维场景中提取物体还有待研究。目前已有的思路是先提取关键点，再使用各种局部特征描述子对关键点进行描述，最后与待检测物体进行比对，得到点-点的匹配。个别文章在之后还采取了ICP对匹配结果进行优化。

　对于缺乏表面纹理信息，或局部曲率变化很小，或点云本身就非常稀疏的物体，采用局部特征描述子很难有效的提取到匹配对。所以就有了所谓基于Point Pair 的特征，该特征使用了一些全局的信息来进行匹配，更神奇的是，最终的位姿估计结果并不会陷入局部最小值。详细可参见论文：Model globally, match locally: Efficient and robust 3D object recognition. 与 Going further with point pair features。SLAM的重要研究方向object based Slam 也声称使用了Point Pair Feature进行匹配。

为了更好的理解这种方法，而在pcl中也没有找到现成的算法，所以我自己用matlab实现了一遍。

算法的思想很简单：

0、ppf 特征为[d,<d,n1>,<d,n2>,<n1,n2>].

1、针对目标模型，在两两点之间构造点对特征F，如果有N个点，那么就有N*N个特征（说明此算法是O(N2)的），N*N个特征形成特征集F_Set

2、在场景中任意取1定点a,再任意取1动点b，构造ppf特征，并从F_set中寻找对应的，那么理想情况下，如果找到了完全匹配的特征，则可获得点云匹配的结果。

3、此算法是一种投票算法，每次匹配都能得到一个旋转角度，如果m个b都投票给了某一旋转角度则可认为匹配成功

这个算法最大的问题就是不停的采样会导致极大的计算量。不过算法本身确实可以匹配物体和场景。

ppf 特征的构建

1 function obj = ppf(point1,point2)
2     d = point1.Location - point2.Location;
3     d_unit = d/norm(d);
4     apha1 = acos(point1.Normal*d_unit');
5     apha2 = acos(point2.Normal*d_unit');
6     apha3 = acos(point1.Normal*point2.Normal');
7     obj = [norm(d),apha1,apha2,apha3];
8 end

ppf 特征集的构建

 1 classdef modelFeatureSet < handle
 2     %MODELFEATURESET 此处显示有关此类的摘要
 3     %   此处显示详细说明
 4     
 5     properties
 6         FeatureTree
 7         ModelPointCloud
 8         Pairs
 9     end
10     
11     methods
12         function obj = modelFeatureSet(pt)
13             obj.ModelPointCloud = copy(pt.removeInvalidPoints());
14         end
15         function growTree(self)
16             self.ModelPointCloud = pcdownsample(self.ModelPointCloud,'GridAverage',.1);
17             pt_size = self.ModelPointCloud.Count;
18             idx = repmat(1:pt_size,pt_size,1);
19             tmp1 = reshape(idx,pt_size*pt_size,1);
20             tmp2 = reshape(idx',pt_size*pt_size,1);
21             pairs = [tmp1,tmp2];
22             rnd = randseed(1,1000,1,1,pt_size*pt_size);
23             pairs = pairs(rnd,:);
24             Features = zeros(size(pairs,1),4);
25             for i = 1:size(pairs,1)
26                 Features(i,:) = ppf(self.ModelPointCloud.select(pairs(i,1)),...
27                                     self.ModelPointCloud.select(pairs(i,2)));
28             end
29             self.FeatureTree = createns(Features);
30             self.Pairs = pairs;
31         end
32     end
33 end

PFH

正如点特征表示法所示，表面法线和曲率估计是某个点周围的几何特征基本表示法。虽然计算非常快速容易，但是无法获得太多信息，因为它们只使用很少的几个参数值来近似表示一个点的k邻域的几何特征。然而大部分场景中包含许多特征点，这些特征点有相同的或者非常相近的特征值，因此采用点特征表示法，其直接结果就减少了全局的特征信息。本小节介绍三维特征描述子中的一位成员：点特征直方图(Point Feature Histograms)，我们简称为PFH，本小节将介绍它的理论优势，从PCL实现的角度讨论其实施细节。PFH特征不仅与坐标轴三维数据有关，同时还与表面法线有关。

理论基础

PFH计算方式通过参数化查询点与邻域点之间的空间差异，并形成一个多维直方图对点的k邻域几何属性进行描述。直方图所在的高维超空间为特征表示提供了一个可度量的信息空间，对点云对应曲面的6维姿态来说它具有不变性，并且在不同的采样密度或邻域的噪音等级下具有鲁棒性。点特征直方图（PFH）表示法是基于点与其k邻域之间的关系以及它们的估计法线，简言之，它考虑估计法线方向之间所有的相互作用，试图捕获最好的样本表面变化情况，以描述样本的几何特征。因此，合成特征超空间取决于每个点的表面法线估计的质量。如图1所示，表示的是一个查询点(Pq)的PFH计算的影响区域，Pq 用红色标注并放在圆球的中间位置，半径为r，(Pq)的所有k邻元素（即与点Pq的距离小于半径r的所有点）全部互相连接在一个网络中。最终的PFH描述子通过计算邻域内所有两点之间关系而得到的直方图，因此存在一个O(k)的计算复杂性。

图1 查询点的PFH计算的影响区域

为了计算两点Pi和Pj及与它们对应的法线Ni和Nj之间的相对偏差，在其中的一个点上定义一个固定的局部坐标系，如图2所示。

图2 定义一个固定的局部坐标系

使用上图中uvw坐标系，法线和之间的偏差可以用一组角度来表示，如下所示：

d是两点Ps和Pt之间的欧氏距离,。计算k邻域内的每一对点的四组值，这样就把两点和它们法线相关的12个参数（xyz坐标值和法线信息）减少到4个。

为每一对点估计PFH四元组，可以使用：

computePairFeatures (const Eigen::Vector4f&p1,const Eigen::Vector4f&n1,

const Eigen::Vector4f&p2,const Eigen::Vector4f&n2,

float&f1,float&f2,float&f3,float&f4);

有关其他详细信息，请见API文件。为查询点创建最终的PFH表示，所有的四元组将会以某种统计的方式放进直方图中，这个过程首先把每个特征值范围划分为b个子区间，并统计落在每个子区间的点数目，因为四分之三的特征在上述中为法线之间的角度计量，在三角化圆上可以将它们的参数值非常容易地归一到相同的区间内。一个统计的例子是：把每个特征区间划分成等分的相同数目，为此在一个完全关联的空间内创建有个区间的直方图。在这个空间中，一个直方图中某一区间统计个数的增一对应一个点的四个特征值。如图3所示，就是点云中不同点的点特征直方图表示法的一个例子，在某些情况下，第四个特征量d在通常由机器人捕获的2.5维数据集中的并不重要，因为临近点间的距离从视点开始是递增的，而并非不变的，在扫描中局部点密度影响特征时，实践证明省略d是有益的。

图3 点云中不同点的点特征直方图表示法

注意：更多相关信息和数学推导，包括不同几何体表面点云的PFH特征分析，请见[RusuDissertation]。

估计PFH特征

点特征直方图（PFH）在PCL中的实现是pcl_features模块的一部分。默认PFH的实现使用5个区间分类（例如：四个特征值中的每个都使用5个区间来统计），其中不包括距离（在上文中已经解释过了——但是如果有需要的话，也可以通过用户调用computePairFeatures方法来获得距离值），这样就组成了一个125浮点数元素的特征向量（3⁵），其保存在一个pcl::PFHSignature125的点类型中。以下代码段将对输入数据集中的所有点估计其对应的PFH特征。

#include <pcl/point_types.h> //点类型头文件

#include <pcl/features/pfh.h> //pfh特征估计类头文件

...//其他相关操作

pcl::PointCloud<pcl::PointXYZ>::Ptrcloud(newpcl::PointCloud<pcl::PointXYZ>);

pcl::PointCloud<pcl::Normal>::Ptrnormals(newpcl::PointCloud<pcl::Normal>());

...//打开点云文件估计法线等

//创建PFH估计对象pfh，并将输入点云数据集cloud和法线normals传递给它

pcl::PFHEstimation<pcl::PointXYZ,pcl::Normal,pcl::PFHSignature125>pfh;

pfh.setInputCloud(cloud);

pfh.setInputNormals(normals);

//如果点云是类型为PointNormal,则执行pfh.setInputNormals (cloud);

//创建一个空的kd树表示法，并把它传递给PFH估计对象。

//基于已给的输入数据集，建立kdtree

pcl::KdTreeFLANN<pcl::PointXYZ>::Ptrtree(newpcl::KdTreeFLANN<pcl::PointXYZ>());

pfh.setSearchMethod(tree);

//输出数据集

pcl::PointCloud<pcl::PFHSignature125>::Ptrpfhs(newpcl::PointCloud<pcl::PFHSignature125>());

//使用半径在5厘米范围内的所有邻元素。

//注意：此处使用的半径必须要大于估计表面法线时使用的半径!!!

pfh.setRadiusSearch(0.05);

//计算pfh特征值

pfh.compute(*pfhs);

// pfhs->points.size ()应该与input cloud->points.size ()有相同的大小，即每个点都有一个pfh特征向量

PFHEstimation类的实际计算程序内部只执行以下：

对点云P中的每个点p

1.得到p点的最近邻元素

2.对于邻域内的每对点，计算其三个角度特征参数值

3.将所有结果统计到一个输出直方图中

使用下列代码，从一个k-邻域计算单一的PFH描述子：

computePointPFHSignature (const pcl::PointCloud<PointInT> &cloud,

const pcl::PointCloud<PointNT> &normals,

const std::vector<int> &indices,

int nr_split,

Eigen::VectorXf&pfh_histogram);

此处，cloud变量是包含点的输入点云，normals变量是包含对应cloud的法线的输入点云，indices代表输入点云（点与法线对应）中查询点的k-近邻元素集，nr_split是所分区间的数目，用于每个特征区间的统计过程，pfh_histogram是浮点数向量来存储输出的合成直方图。

FPFH

已知点云P中有n个点，那么它的点特征直方图（PFH）的理论计算复杂度是，其中k是点云P中每个点p计算特征向量时考虑的邻域数量。对于实时应用或接近实时应用中，密集点云的点特征直方图（PFH）的计算，是一个主要的性能瓶颈。本小节讲述PFH计算方式的简化形式，我们称为快速点特征直方图FPFH（Fast Point Feature Histograms）（更多详情，请见[RusuDissertation]），FPFH把算法的计算复杂度降低到了，但是任然保留了PFH大部分的识别特性。

理论基础

为了简化直方图的特征计算，我们执行以下过程：

第一步，对于每一个查询点，计算这个点和它的邻域点之间的一个元组（参考上一节PFH的介绍），第一步结果我们称之为简化的点特征直方图SPFH（Simple Point Feature Histograms）；

第二步，重新确定每个点的k邻域，使用邻近的SPFH值来计算的最终直方图（称为FPFH），如下所示：

上式中，权重在一些给定的度量空间中，表示查询点和其邻近点之间的距离，因此可用来评定一对点（ , ），但是如果需要的话，也可以把用另一种度量来表示。如图1所示可以帮助理解这个权重方式的重要性，它表示的是以点为中心的k邻域影响范围。

图1 以点为中心的k邻域影响范围图

因此，对于一个已知查询点，这个算法首先只利用和它邻域点之间对应对（上图中以红色线来说明），来估计它的SPFH值，很明显这样比PFH的标准计算少了邻域点之间的互联。点云数据集中的所有点都要执行这一计算获取SPFH，接下来使用它的邻近点的SPFH值和点的SPFH值重新权重计算，从而得到点的最终FPFH值。FPFH计算添加的计算连接对，在上图中以黑色线表示。如上图所示，一些重要对点（与直接相连的点）被重复计数两次（图中以粗线来表示），而其他间接相连的用细黑线表示。

PFH和FPFH的区别

PFH和FPFH计算方式之间的主要区别总结如下：

1.FPFH没有对全互连点的所有邻近点的计算参数进行统计，从图12-18中可以看到，因此可能漏掉了一些重要的点对，而这些漏掉的对点可能对捕获查询点周围的几何特征有贡献。

2.PFH特征模型是对查询点周围的一个精确的邻域半径内，而FPFH还包括半径r范围以外的额外点对（不过在2r内）；

3.因为重新权重计算的方式，所以FPFH结合SPFH值，重新捕获邻近重要点对的几何信息；

4.由于大大地降低了FPFH的整体复杂性，因此FPFH有可能使用在实时应用中；

5.通过分解三元组，简化了合成的直方图。也就是简单生成d分离特征直方图，对每个特征维度来单独绘制，并把它们连接在一起（见下2图）。

图2 PFH与FPFH示意图

估计FPFH特征

快速点特征直方图FPFH在点云库中的实现可作为pcl_features库的一部分。默认的FPFH实现使用11个统计子区间（例如：四个特征值中的每个都将它的参数区间分割为11个），特征直方图被分别计算然后合并得出了浮点值的一个33元素的特征向量，这些保存在一个pcl::FPFHSignature33点类型中。以下代码段将对输入数据集中的所有点估计一组FPFH特征值。

#include

#include //fpfh特征估计类头文件声明

...//其他相关操作

pcl::PointCloud<pcl::PointXYZ>::Ptrcloud(newpcl::PointCloud<pcl::PointXYZ>);

pcl::PointCloud<pcl::Normal>::Ptrnormals(newpcl::PointCloud<pcl::Normal>());

...//打开点云文件估计法线等

//创建FPFH估计对象fpfh，并把输入数据集cloud和法线normals传递给它。

pcl::FPFHEstimation<pcl::PointXYZ,pcl::Normal,pcl::FPFHSignature33>fpfh;

fpfh.setInputCloud(cloud);

fpfh.setInputNormals(normals);

//如果点云是类型为PointNormal，则执行fpfh.setInputNormals (cloud);

//创建一个空的kd树对象tree，并把它传递给FPFH估计对象。

//基于已知的输入数据集，建立kdtree

pcl::search::KdTree<PointXYZ>::Ptrtree(newpcl::search::KdTree<PointXYZ>);

fpfh.setSearchMethod(tree);

//输出数据集

pcl::PointCloud<pcl::FPFHSignature33>::Ptrfpfhs(newpcl::PointCloud<pcl::FPFHSignature33>());

//使用所有半径在5厘米范围内的邻元素

//注意：此处使用的半径必须要大于估计表面法线时使用的半径!!!

fpfh.setRadiusSearch(0.05);

//计算获取特征向量

fpfh.compute(*fpfhs);

//fpfhs->points.size ()应该和input cloud->points.size ()有相同的大小，即每个点有一个特征向量

FPFHEstimation类的实际计算内部只执行以下操作：

对点云P中的每个点p

第一步：

1.得到:math:`p`的邻域元素

2. 计算每一对:math:`p, p_k`的三个角度参数值（其中:math:`p_k`是:math:`p`的邻元素）

3.把所有结果统计输出到一个SPFH直方图

第二步:

1.得到:math:`p`的最近邻元素

2.使用:math:`p`的每一个SPFH和一个权重计算式，来计算最终:math:`p`的FPFH

利用OpenMP提高FPFH速度

对于计算速度要求苛刻的用户，PCL提供了一个FPFH估计的另一实现，它使用多核/多线程规范，利用OpenMP开发模式来提高计算速度。这个类的名称是pcl::FPFHEstimationOMP，并且它的应用程序接口（API）100%兼容单线程pcl::FPFHEstimation，这使它适合作为一个替换元件。在8核系统中，OpenMP的实现可以在6-8倍更快的计算时间内完全同样单核系统上的计算。