【PCL】—常用的点云滤波算法详解

参考：

1. https://www.freesion.com/article/6131168793/
2. https://www.cnblogs.com/ironstark/p/4991232.html

1. 点云滤波的概念

点云滤波是点云处理的基本步骤，也是进行 high level 三维图像处理之前必须要进行的预处理。其作用类似于信号处理中的滤波，但实现手段却和信号处理不一样。我认为原因有以下几个方面：

1. 点云不是函数，对于复杂三维外形其x,y,z之间并非以某种规律或某种数值关系定义，所以点云无法建立横纵坐标之间的联系。
2. 点云在空间中是离散的，和图像、信号不一样，并不定义在某个区域上，无法以某种模板的形式对其进行滤波。换言之，点云没有图像与信号那么明显的定义域。
3. 点云在空间中分布很广泛。遍历整个点云中的每个点，并建立点与点之间的相互位置关系成了最大难点。不像图像与信号，可以有迹可循。
4. 点云滤波依赖于几何信息，而不是数值关系。

综上所述，点云滤波只在抽象意义上与信号、图像滤波类似。因为滤波的功能都是突出需要的信息。

2. 点云滤波的方法

PCL对常规滤波手段均进行了很好的封装，对点云的滤波通过调用各个滤波器对象来完成。主要的滤波器有直通滤波器、体素滤波器、统计滤波器、半径滤波器等。不同特性的滤波器构成了较为完整的点云前处理族，并组合使用完成任务。实际上，滤波手段的选择和采集方式是密不可分的。

1. 如果使用线结构光扫描的方式采集点云，必然是物体沿z方向分布较广，但x、y方向的分布处于有限范围内。此时可使用直通滤波器，确定点云在x或y方向上的范围，可较快剪除离群点，达到第一步粗处理的目的。
2. 如果使用高分辨率相机等设备对点云进行采集，往往点云会较为密集。过多的点云数量会对后续分割工作带来困难。体素滤波器可以达到向下采样同时不破坏点云本身几何结构的功能。点云几何结构不仅是宏观的几何外形，也包括其微观的排列方式，比如横向相似的尺寸、纵向相同的距离。随机下采样虽然效率比体素滤波器高，但会破坏点云微观结构。
3. 统计滤波器用于去除明显离群点（离群点往往由测量噪声引入）。其特征是在空间中分布稀疏，可以理解为：每个点都表达一定信息量，某个区域点越密集则可能信息量越大。噪声信息属于无用信息，信息量较小。所以离群点表达的信息可以忽略不计。考虑到离群点的特征，则可以定义某处点云小于某个密度，即点云无效。计算每个点到其最近的k个点平均距离，则点云中所有点的距离应构成高斯分布。给定均值与方差，可剔除3∑之外的点。
4. 半径滤波器与统计滤波器相比更加简单粗暴。以某点为中心画一个圆计算落在该圆中点的数量，当数量大于给定值时，则保留该点，数量小于给定值则剔除该点。此算法运行速度快，依序迭代留下的点一定是最密集的，但是圆的半径和圆内点的数目都需要人工指定。

实际上点云滤波的手段和传统的信号滤波与图像滤波在自动化程度、滤波效果上还有很大的差距。学者大多关注图像识别与配准算法在点云处理方面的移植，而对滤波算法关注较少。其实点云前处理对测量精度与识别速度都有很大影响。

3. 点云库对滤波算法的实现

点云库中已经包含了上述所有的滤波算法。PCL滤波算法的实现是通过滤波器类来完成的，需要实现滤波功能时则新建一个滤波器对象并设置参数，从而保证可以针对不同的滤波任务，使用不同参数的滤波器对点云进行处理。

3.1 直通滤波器

// Create the filtering object
pcl::PassThrough<pcl::PointXYZ> pass;
pass.setInputCloud(cloud);
pass.setFilterFieldName("z");
pass.setFilterLimits(0.0, 1.0);
//pass.setFilterLimitsNegative(true);
pass.filter(*cloud_filtered);

3.2 体素滤波器

// 创建过滤对象 
  pcl::VoxelGrid<pcl::PCLPointCloud2> sor;
  sor.setInputCloud(cloud);
  sor.setLeafSize(0.01f，0.01f，0.01f);
  sor.filter(*cloud_filtered);

3.3 统计滤波器

  // 创建过滤对象 
  pcl::StatisticalOutlierRemoval<pcl::PointXYZ> sor;
  sor.setInputCloud(cloud);
  sor.setMeanK(50);
  sor.setStddevMulThresh(1.0);
  sor.filter(*cloud_filtered);

3.4 半径滤波器

// 构建过滤器 
pcl::RadiusOutlierRemoval<pcl::PointXYZ> outrem;
outrem.setInputCloud(cloud);
outrem.setRadiusSearch(0.8);
outrem.setMinNeighborsInRadius(2);
// apply filter 
outrem.filter(*cloud_filtered);

显然，不同的滤波器在滤波过程中，总是先创建一个对象，再设置对象参数，最后调用滤波函数对点云进行处理（点云为智能指针指向的一块区域）。

4. 基于点云频率的滤波方法

4.1 点云的频率

• 今天在阅读分割有关的文献时，惊喜的发现，点云和图像一样，有可能也存在频率的概念。但这个概念并未在文献中出现也未被使用，谨在本博文中滥用一下“高频”一词。点云表达的是三维空间中的一种信息，这种信息本身并没有一一对应的函数值。故点云本身并没有在讲诉一种变化的信号。但在抽象意义上，点云必然是在表达某种信号的，虽然没有明确的时间关系，但应该会存在某种空间关系（例如激光雷达点云）。我们可以人为的指定点云空间中的一个点（例如扫描的重心或LIDAR的“源”），基于此点来讨论点云在各个方向上所谓的频率。
• 在传统的信号处理中，高频信号一般指信号变化快，低频信号一般指信号变化缓慢。在图像处理中，高低频的概念被引申至不同方向上图像灰度的变化，傅里叶变换可以用于提取图像的周期成分滤除布纹噪声。在点云处理中定义：点云法线向量差为点云所表达的信号，换言之，如果某处点云曲率大，则点云表达的是一个变化的信号；如果点云曲率小，则其表达的是一个不变的信号。这和我们的直观感受也是相近的，地面曲率小，它表达的信息量也小；人的五官部分曲率大，其表达了整个Scan中最大的信息量。

4.2 基于点云频率的滤波原理

虽然点云频率之前并没有被讨论，但使用频率信息的思想已经被广泛的应用在了各个方面，最著名的莫过于DoN算法。

• DoN算法被作者归类于点云分割算法中，但我认为并不准确，本质上DoN只是一种前处理，应该算是一种比较先进的点云滤波算法。分割本质上还是由欧式分割算法完成的。
• DoN 是 Difference of Normal 的简写。算法的目的是在去除点云低频滤波，低频信息（例如建筑物墙面，地面）往往会对分割产生干扰，高频信息（例如建筑物窗框，路面障碍锥）往往尺度上很小，直接采用 基于临近信息 的滤波器会将此类信息合并至墙面或路面中，所以DoN算法利用了多尺度空间的思想，算法如下：

1. 在小尺度上计算点云法线1；
2. 在大尺度上计算点云法线2；
3. 法线1-法线2；
4. 滤去3中值较小的点；
5. 欧式分割。

显然，在小尺度上是可以对高频信息进行检测的，此算法可以很好的保留小尺度高频信息。其在大规模点云（如LiDar点云）中优势尤其明显。

4.3 PCL对该算法的实现

// Create a search tree, use KDTreee for non-organized data.
  pcl::search::Search<PointXYZRGB>::Ptr tree;
  if (cloud->isOrganized())
  {
    tree.reset(new pcl::search::OrganizedNeighbor<PointXYZRGB> ());
  }
  else
  {
    tree.reset(new pcl::search::KdTree<PointXYZRGB> (false));
  }
  // Set the input pointcloud for the search tree
  tree->setInputCloud(cloud);

  //生成法线估计器（OMP是并行计算，忽略）
  pcl::NormalEstimationOMP<PointXYZRGB, PointNormal> ne;
  ne.setInputCloud(cloud);
  ne.setSearchMethod(tree);
  //设定法线方向（要做差，同向很重要）
  ne.setViewPoint(std::numeric_limits<float>::max (), std::numeric_limits<float>::max (), std::numeric_limits<float>::max ());

  //计算小尺度法线
  pcl::PointCloud<PointNormal>::Ptr normals_large_scale (new pcl::PointCloud<PointNormal>);
  ne.setRadiusSearch(scale2);
  ne.compute(*normals_large_scale);
  //计算大尺度法线
  pcl::PointCloud<PointNormal>::Ptr normals_large_scale (new pcl::PointCloud<PointNormal>);
  ne.setRadiusSearch(scale2);
  ne.compute(*normals_large_scale);

  //生成DoN分割器
  pcl::DifferenceOfNormalsEstimation<PointXYZRGB, PointNormal, PointNormal> don;
  don.setInputCloud(cloud);
  don.setNormalScaleLarge(normals_large_scale);
  don.setNormalScaleSmall(normals_small_scale);

  //计算法线差
  PointCloud<PointNormal>::Ptr doncloud (new pcl::PointCloud<PointNormal>);
  copyPointCloud<PointXYZRGB, PointNormal>(*cloud, *doncloud);
  don.computeFeature(*doncloud);

  //生成滤波条件：把法线差和阈值比
  pcl::ConditionOr<PointNormal>::Ptr range_cond (new pcl::ConditionOr<PointNormal> ());
  range_cond->addComparison(pcl::FieldComparison<PointNormal>::ConstPtr (new pcl::FieldComparison<PointNormal> ("curvature", pcl::ComparisonOps::GT, threshold)));
  //生成条件滤波器,输入滤波条件和点云
  pcl::ConditionalRemoval<PointNormal> condrem (range_cond);
  condrem.setInputCloud(doncloud);

  //导出滤波结果
  pcl::PointCloud<PointNormal>::Ptr doncloud_filtered (new pcl::PointCloud<PointNormal>);
  condrem.filter(*doncloud_filtered);

  //欧式聚类~~~（略）

算法运行过程可用图表示为：