一文详解点云分割算法-轻识

点击下方卡片，关注“新机器视觉”公众号

重磅干货，第一时间送达

作者丨书生意封侯@知乎

来源丨https://zhuanlan.zhihu.com/p/470782623

编辑丨3D视觉工坊

从某种意义上说，地面点剔除（分割）也属于点云分割的一种，但两者技术路线有所不同，因此本节主要是针对地面点剔除后的点云分割。

传统点云分割的方式，通过查阅文档后发现网上文档大都雷同，但由于时间关系且并不是本次学习的重点，因此本次记录暂不详细探究，但对相关论文进行下载。

1.1基于边缘的分割方法

边缘是描述点云物体形状的基本特征，这种方法检测点云一些区域的边界来获取分割区域，这些方法的原理是定位出边缘点的强度变化。论文[1]提出了一种边缘检测技术，通过计算梯度，检测表面上单位法向量方向的变化来拟合线段。论文[2]是基于扫描线的分组进行快速分割。

基于边缘的方法虽然分割速度比较快但是准确度不能保证，因为边缘对于噪声和不均匀的或稀疏的点云非常敏感。

[1]Range data processing:Representation of surfaces by edges. In proc.int. Pattern recognition conf, 1995——范围数据处理：用边界表示表面

[2] Fast range image segmentation using high-level segmentation primitives, In 3rd IEEE Workshop on Applications of Computer Vision, USA, 1996——基础扫描线分组的快速分割

1.2基于区域的分割方法

基于区域的方法使用邻域信息来将具有相似属性的附近点归类，以获得到分割区域，并区分出不同区域之间的差异性。基于区域的方法比基于边缘的方法更准确。但是他们在分割过度或不足以及在如何准确确定区域边界方面存在问题。研究者们将基于区域的方法分为两类：种子区域（自下而上）方法和非种子区域（自上而下）方法。

1.2.1种子区域方法

基于种子的区域分割通过选择多个种子点来开始做分割，从这些种子点为起始点，通过添加种子的邻域点的方式逐渐形成点云区域，算法主要包含了两个步骤（论文[3]）：

（1）基于每个点的曲率识别种子点，

（2）根据预定标准，该标准可以是点的相似度和点云的表面的相似度来生长种子点。

这种方法对噪声点也非常敏感，并且耗时。但后续有很多基于这种方法的改进，比如对于激光雷达数据的区域增长的方法，提出了基于种子点的法向量和与生长平面的距离来生长种子点。种子区域方法高度依赖于选定的种子点。不准确选择种子点会影响分割过程，并可能导致分割不足或过度。选择种子点以及控制生长过程是耗时的。分割结果可能对所选的兼容性阈值敏感。另一个困难是决定是否在给定区域中添加点，因为这种方法对点云的噪声也很敏感。该方法详情由5.2.4介绍。

pcl中提供了基于欧式聚类的区域增长方式，以第一个点为标准作为种子点，候选其周边的点作为它的对比或者比较的对象，如果满足条件就加入到聚类的对象中。其主要的缺点：该算法没有初始化种子系统，没有过度分割或者分割不足的控制，还有就是从主循环运算中调用条件函数时，效率比较低。该方法详情由5.3.2介绍。

1.2.2非种子区域方法

这种方法时基于自上而下的方法。首先，所有点都分为一个区域。然后细分过程开始将其划分为更小的区域。使用这种方法（论文[4]）指导聚类平面区域的过程，以重建建筑物的完整几何形状。该工作引入了基于局部区域的置信率为平面的分割方法。这种方法的局限性在于它也会可能过度分割，并且在分割其他对象（例如树）时它不能很好地执行。非种子区域方法的主要困难是决定细分的位置和方式。这些方法的另一个限制是它们需要大量的先验知识（例如，对象模型，区域数量等），然后这些未知的先验知识在复杂场景中通常是未知的。

[3]P.J. Besl, R.C. Jain, Segmentation through variable order surface fitting, IEEE Transaction on Pattern Analysis and Machine Intelligence 10, 1988.——变阶曲面拟合分割

[4]Architectural Modeling from Sparsely Scanned Range Data——稀疏扫描范围内的建筑建模

1.2.3基于图像区域增长方法

该方法主要是通过将点云转换成二值图像，再通过图像方法中的区域增长进行聚类，再转换成点云。其优点为速度快，而缺点是存在过度分割以及分割不足。

1.将激光点云数据放入栅格中，并生成对应的深度图像，每格图像中的像素值用转换后的高度表示，分别用像素值中的r、g、b中的0-b表示高度最大值，1-g-表示高度最小值，2-r表示有无聚类的标志（0-未聚类，1-255-聚类ID）

/* @brief:将点云转换为vector,俯视图栅格化

 * @param [in]: in_cloud,输入点云
 * @param [out]: out_cloud,转换的点云vector
 * @return NONE
   */
   void FrontLidarAlg::convetCloudToVector(pcl::PointCloud::Ptr in_cloud, \
                                       cv::Mat& out_img, \
                                       std::vector>* out_cloud)
   {
     out_img = cv::Mat::zeros(cv::Size(img_width_,img_height_), CV_8UC3);

  volatile int row;
  volatile int col;

  for(int i=0;isize();i++)
  {
    //将横向在检测范围之外的去掉
    if(in_cloud->points[i].x <= cloud_x_min_ || in_cloud->points[i].x >= cloud_x_max_)
    {
      continue;
    }
    //将纵向在检测范围之外的去掉
    if(in_cloud->points[i].y >= cloud_y_max_ || in_cloud->points[i].y <= 0)
    {
      continue;
    }

    //限制最高和最低点大小
    if(in_cloud->points[i].z < min_dect_height_)//高度低于最低点,则赋值为最低点
    {
      in_cloud->points[i].z = min_dect_height_;
    }
    
    if(in_cloud->points[i].z > max_dect_height_)//高度高于最高点,则赋值最高点
    {
      in_cloud->points[i].z = max_dect_height_;
    }
    
    //计算点云所在图像的行数和列数,四舍五入
    col = int((in_cloud->points[i].x - cloud_x_min_)/img_res_);
    row = int(in_cloud->points[i].y/img_res_);
    out_cloud->at(col+row*img_width_).points.push_back(in_cloud->points[i]);
    
    //输入点云高度值转换到图像坐标系下的数值
    int value = (int)((in_cloud->points[i].z - min_dect_height_)/dect_height_res_);
    
    if(out_img.at(row,col)[0] < value)//b-最大高度
    {
      out_img.at(row,col)[0] = value;
      if(0 == out_img.at(row,col)[1])
      {
        out_img.at(row,col)[1] = value;
      }
    }
    else if(out_img.at(row,col)[1] > value)//b-最小高度
    {
      out_img.at(row,col)[1] = value;
    }

  }//for in_cloud->size
}

2、判断图像的边界是否越界；

3、对生成的二值图像做形态学滤波操作，进行闭运算，将临近的小目标进行合并，一定程度上优化过分割；

4、采用递归的方式对图像中的任意一个像素值附近的8邻域进行标记操作，即对所有的像素对应的r值进行赋值（ID）；

5、将相同ID（r值）对应的像素和栅格内的点云取出，即可得到不同类簇的点云。

6、滤波操作，滤除点云较少的类簇。

#include "object_segment.h"

namespace front_lidar {

/* @brief:点云分割构造函数,初始化参数

 * @param [in]: NONE
 * @param [out]: NONE
 * @return NONE
   */
   ObjectSegment::ObjectSegment(string config_file_str)
   {
     config_parser_ptr_.reset(new ConfigParser(config_file_str));

  use_seg_raw_cloud_ = config_parser_ptr_->getInt("use_seg_raw_cloud");//聚类后是否输出原始数据
  img_res_ = config_parser_ptr_->getDouble("img_res");//图像分辨率
  cloud_x_min_ = config_parser_ptr_->getDouble("cloud_x_min");//x轴最小值
  cloud_x_max_ = config_parser_ptr_->getDouble("cloud_x_max");//x轴最大值
  cloud_y_min_ = config_parser_ptr_->getDouble("cloud_y_min");//y轴最小值
  cloud_y_max_ = config_parser_ptr_->getDouble("cloud_y_max");//y轴最大值
  min_cluster_size_ = config_parser_ptr_->getDouble("object_detect_min_cluster_size");//最小聚类个数
  max_cluster_size_ = config_parser_ptr_->getDouble("object_detect_max_cluster_size");//最大聚类个数
  cell_size_ = config_parser_ptr_->getDouble("cell_size");//聚类临域尺寸

  min_dect_height_ = -10;//最小检测高度
  max_dect_height_ = 10;//最大检测高度
  dect_height_res_ = (max_dect_height_ - min_dect_height_)/255;//检测高度分辨率

  region_growing_clusters_.clear();
  region_growing_img_.resize(0);//区域增长后的图像

}

/* @brief:点云分割处理,将点云分割为不同障碍物点云

  * @param [in]: in_img,输入图像栅格, in_vector-输入栅格化的点云
  * @param [out]: cloud_clusters_ptr,输出聚类结果
  * @return NONE
    */
    void ObjectSegment::process(const std::vector>* in_vector, \
                          const cv::Mat& in_img, \
                          std::vector::Ptr>* cloud_clusters_ptr)
    {
      cloud_clusters_ptr->clear();

  cluster_num_ = 0;

  //区域增长,图像聚类分割
  pcl::StopWatch regin_grow_timer;
  reginGrowing(in_img, region_growing_img_);
  LOG_Debug()<<"regin_grow time:"<
  //将图像转换为点云
  pcl::StopWatch vector_to_cloud_timer;
  if(1 == use_seg_raw_cloud_)
  {
    convetVectorToCloud(in_vector, region_growing_img_, ®ion_growing_clusters_);
  }
  else
  {
    convetImageToCloud(region_growing_img_, ®ion_growing_clusters_);
  }
  LOG_Debug()<<"vector_to_cloud time:"<
  //滤除聚类点数少于设定阈值的类别
  pcl::StopWatch filter_timer;
  filterObjectsCloud(region_growing_clusters_, cloud_clusters_ptr);
  LOG_Debug()<<"filter time:"<
}


/* @brief:区域增长,8邻域聚类

 * @param [in]: in_img-输入图像
 * @param [out]: out_img,转换的图像,0-b-最大值,1-g-最小值,2-r-有无聚类标志(0-未聚类,1-255-已聚类ID)
 * @return NONE
   */
   void ObjectSegment::reginGrowing(const cv::Mat &in_img, cv::Mat &out_img)
   {
     unsigned short class_id = 0;
     unsigned short element_class_id = 0;
     out_img = in_img.clone();

  for(int row = 1; row < out_img.rows -1; row++)
  {
    for(int col = 1; col < out_img.cols - 1; col++)
    {
      //像素的高度最大值为0,则该像素为空,无效像素
      if(0 == out_img.at(row,col)[0])
      {
        continue;
      }
      //像素的类别标记为空,即为未标记类别的像素,则分配类别ID
      if(0 == out_img.at(row,col)[2])
      {
        if(class_id > 250)//超出限定类别总数,则返回
        {
          return;
        }
        out_img.at(row,col)[2] = ++class_id;//给该像素赋值类别ID
      }
      element_class_id = out_img.at(row,col)[2];
      //根据输入种子像素,递归区域增长
      elementReginGrowing(out_img, row, col, element_class_id);
    }//col
  }//row
  cluster_num_ = class_id;
}

/* @brief:单个元素区域增长,cell_size_邻域聚类

 * @param [in]: 输入图像-in_img,当前元素的r和c,输出图像-in_img
 * @param [out]: out_img,转换的图像,0-b-最大值,1-g-最小值,2-r-有无聚类标志(0-未聚类,1-255-聚类ID)
 * @return NONE
   */
   void ObjectSegment::elementReginGrowing(cv::Mat &in_img, int row, int col, unsigned short class_id)
   {
     int start_row;//起始行
     int end_row;//截止行

  int start_col;//起始列
  int end_col;//截止列

  //判断起始行 截止行 起始列 截止列,防止越界
  start_row = row - cell_size_;
  end_row = row + cell_size_;
  start_col = col - cell_size_;
  end_col = col + cell_size_;

  if(start_row < 0) 
    start_row = 0;
  if(end_row > in_img.rows - 1)
    end_row = in_img.rows - 1;
  if(start_col < 0)
    start_col = 0;
  if(end_col > in_img.cols - 1)
    end_col = in_img.cols - 1; 

  for(int m = start_row; m <= end_row; m++)
  {
    for(int n = start_col; n <= end_col; n++)
    {
      //该像素对应最大高度为0,则为无效像素
      if(0 == in_img.at(m,n)[0])
      {
         continue;
      }
      //该元素为初始输入,已经标记过
      if(m==row && n==col)
      {
        continue;
      }
      //未标记过的元素,标记该元素,以该元素为种子,进行区域增长
      if(0 == in_img.at(m,n)[2])
      {
         in_img.at(m,n)[2] = class_id;//标记元素
         elementReginGrowing(in_img, m, n, class_id);//以该元素为初始值,进行区域增长
      }
    }//for col
  }//for row
}

/* @brief:将vector转化为聚类结果,2D到3D

 * @param [in]: in_img,输入的栅格图像
 * @param [out]: out_cloud_ptr,转换的聚类结果
 * @return NONE
   */
   void ObjectSegment::convetImageToCloud(const cv::Mat& in_img, \
                                     std::vector::Ptr>* out_cloud_ptr)
   {
     out_cloud_ptr->clear();
     pcl::PointXYZI point;
     int class_id = 0;//类ID,1开始递增

  //根据类别分配向量空间
  out_cloud_ptr->resize(cluster_num_);//扩充输出向量大小
  //分配点云空间
  for(int i=0;i  {
    out_cloud_ptr->at(i).reset(new pcl::PointCloud);
  }
  //将图像转换为点云
  for(int row=0;row  {
    for(int col=0;col    {
      //如果分类标记为0,则为无效数据
      if(0 == in_img.at(row,col)[2])
      {
        continue;
      }

      class_id = in_img.at(row,col)[2];//获取该区域的类别ID
      
      //图像坐标转换为点云坐标
      point.x = (col+0.5)*img_res_ + cloud_x_min_;
      point.y = (row+0.5)*img_res_;
    
      //点云坐标高度最大值
      point.z = in_img.at(row,col)[0]*dect_height_res_ + min_dect_height_;
      out_cloud_ptr->at(class_id-1)->points.push_back(point);
    
      //点云坐标高度最小值
      point.z = in_img.at(row,col)[1]*dect_height_res_ + min_dect_height_;
      out_cloud_ptr->at(class_id-1)->points.push_back(point);
    }

  }
}

/* @brief:将vector转化为聚类结果,2D到3D

 * @param [in]: in_vector,输入的栅格化后点云，in_img,输入的聚类标志位的图像
 * @param [out]: out_cloud_ptr,转换的聚类结果
 * @return NONE
   */
   void ObjectSegment::convetVectorToCloud(const std::vector>* in_vector, \
                                       const cv::Mat& in_img, \
                                       std::vector::Ptr>* out_cloud_ptr)
   {
     out_cloud_ptr->clear();

  int class_id = 0;//类ID,1开始递增

  //根据类别分配向量空间
  out_cloud_ptr->resize(cluster_num_);//扩充输出向量大小
  //分配点云空间
  for(int i=0;i  {
    out_cloud_ptr->at(i).reset(new pcl::PointCloud);
  }
  //将图像转换为点云
  for(int row=0;row  {
    for(int col=0;col    {
      //如果分类标记为0,则为无效数据
      if(0 == in_img.at(row,col)[2])
      {
        continue;
      }

      class_id = in_img.at(row,col)[2];//获取该区域的类别ID
      
      int count = in_vector->at(col+row*in_img.cols).points.size();
      if(count > 25)
        count = 25;
    
      for(int j = 0; j < count; j++)
      {
        out_cloud_ptr->at(class_id-1)->points.push_back(in_vector->at(col+row*in_img.cols).points[j]);
      }
    }

  }
}

/* @brief:滤除聚类后点云数量小于设定值的类别

 * @param [in]: in_objects_cloud-聚类后的目标点云
 * @param [out]: out_objects_cloud_ptr-滤除类别中点数量少于设定阈值的类别
 * @return NONE
   */
   void ObjectSegment::filterObjectsCloud(const std::vector::Ptr> in_objects_cloud, \
                                         std::vector::Ptr> *out_objects_cloud_ptr)
   {
     out_objects_cloud_ptr->clear();

  for(int i=0;i  {
    if(in_objects_cloud[i]->size() < min_cluster_size_)
    {
      continue;
    }

    out_objects_cloud_ptr->push_back(in_objects_cloud[i]);

  }
}

}//namespace front_lidar

1.2.4基于点云区域增长方法

算法思想：

首先依据点的曲率值对点进行排序，之所以排序是因为，区域生长算法是从曲率最小的点开始生长的，这个点就是初始种子点，初始种子点所在的区域即为最平滑的区域，从最平滑的区域开始生长可减少分割片段的总数，提高效率，设置一空的种子点序列和空的聚类区域，选好初始种子后，将其加入到种子点序列中，并搜索邻域点，对每一个邻域点，比较邻域点的法线与当前种子点的法线之间的夹角，小于平滑阀值的将当前点加入到当前区域，然后检测每一个邻域点的曲率值，小于曲率阀值的加入到种子点序列中，删除当前的种子点，循环执行以上步骤，直到种子序列为空。

算法步骤：

1.种子周围的临近点和种子点云相比较

2.法线的方向是否足够相近

3.曲率是否足够小

3.如果满足1，2则该点可用做种子点

4.如果只满足1，则归类而不做种

5.从某个种子出发，其“子种子”不再出现则一类聚集完成

6.类的规模既不能太大也不能太小

上述算法是针对小曲率变化面设计的。尤其适合对连续阶梯平面进行分割：比如SLAM算法所获得的建筑走廊。

算法源码：

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

int
main (int argc, char** argv)
{ 
  //点云的类型
  pcl::PointCloud::Ptr cloud (new pcl::PointCloud);
  //打开点云
  if ( pcl::io::loadPCDFile  ("region_growing_tutorial.pcd", *cloud) == -1)
  {
    std::cout << "Cloud reading failed." << std::endl;
    return (-1);
  }
 //设置搜索的方式或者说是结构
  pcl::search::Search::Ptr tree = boost::shared_ptr > (new pcl::search::KdTree);
   //求法线
  pcl::PointCloud ::Ptr normals (new pcl::PointCloud );
  pcl::NormalEstimation normal_estimator;
  normal_estimator.setSearchMethod (tree);
  normal_estimator.setInputCloud (cloud);
  normal_estimator.setKSearch (50);
  normal_estimator.compute (*normals);
   //直通滤波在Z轴的0到1米之间
  pcl::IndicesPtr indices (new std::vector );
  pcl::PassThrough pass;
  pass.setInputCloud (cloud);
  pass.setFilterFieldName ("z");
  pass.setFilterLimits (0.0, 1.0);
  pass.filter (*indices);
  //聚类对象<点，法线>
  pcl::RegionGrowing reg;
  reg.setMinClusterSize (50);  //最小的聚类的点数
  reg.setMaxClusterSize (1000000);  //最大的
  reg.setSearchMethod (tree);    //搜索方式
  reg.setNumberOfNeighbours (30);    //设置搜索的邻域点的个数
  reg.setInputCloud (cloud);         //输入点
  //reg.setIndices (indices);
  reg.setInputNormals (normals);     //输入的法线
  reg.setSmoothnessThreshold (3.0 / 180.0 * M_PI);  //设置平滑度
  reg.setCurvatureThreshold (1.0);     //设置曲率的阀值

  std::vector  clusters;
  reg.extract (clusters);

  std::cout << "Number of clusters is equal to " << clusters.size () << std::endl;
  std::cout << "First cluster has " << clusters[0].indices.size () << " points." << endl;
  std::cout << "These are the indices of the points of the initial" <<
    std::endl << "cloud that belong to the first cluster:" << std::endl;
 
 int counter = 0;
  while (counter < clusters[0].indices.size ())
  {
    std::cout << clusters[0].indices[counter] << ", ";
    counter++;
    if (counter % 10 == 0)
      std::cout << std::endl;
  }
  std::cout << std::endl;
  
  //可视化聚类的结果
  pcl::PointCloud ::Ptr colored_cloud = reg.getColoredCloud ();
  pcl::visualization::CloudViewer viewer ("Cluster viewer");
  viewer.showCloud(colored_cloud);
  while (!viewer.wasStopped ())
  {
  }

  return (0);
}

1.2.5基于颜色区域增长方法

基于颜色的区域生长分割原理上和基于曲率，法线的分割方法是一致的。只不过比较目标换成了颜色，去掉了点云规模上限的限制。可以认为，同一个颜色且挨得近，是一类的可能性很大，不需要上限来限制。所以这种方式比较适合用于室内场景分割。

算法步骤：

算法分为两步：

（1）分割，当前种子点和领域点之间色差小于色差阀值的视为一个聚类

（2）合并，聚类之间的色差小于色差阀值和并为一个聚类，且当前聚类中点的数量小于聚类点数量的与最近的聚类合并在一起。

pcl中的函数调用：

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

int
main (int argc, char** argv)
{
  pcl::search::Search ::Ptr tree = boost::shared_ptr > (new pcl::search::KdTree);

  pcl::PointCloud ::Ptr cloud (new pcl::PointCloud );
  if ( pcl::io::loadPCDFile  ("region_growing_rgb_tutorial.pcd", *cloud) == -1 )
  {
    std::cout << "Cloud reading failed." << std::endl;
    return (-1);
  }
  //存储点云的容器
  pcl::IndicesPtr indices (new std::vector );
  //滤波
  pcl::PassThrough pass;
  pass.setInputCloud (cloud);
  pass.setFilterFieldName ("z");
  pass.setFilterLimits (0.0, 1.0);
  pass.filter (*indices);
  
 //基于颜色的区域生成的对象
  pcl::RegionGrowingRGB reg;
  reg.setInputCloud (cloud);
  reg.setIndices (indices);   //点云的索引
  reg.setSearchMethod (tree);
  reg.setDistanceThreshold (10);  //距离的阀值
  reg.setPointColorThreshold (6);  //点与点之间颜色容差
  reg.setRegionColorThreshold (5);  //区域之间容差
  reg.setMinClusterSize (600);       //设置聚类的大小

  std::vector  clusters;
  reg.extract (clusters);

  pcl::PointCloud ::Ptr colored_cloud = reg.getColoredCloud ();
  pcl::visualization::CloudViewer viewer ("Cluster viewer");
  viewer.showCloud (colored_cloud);
  while (!viewer.wasStopped ())
  {
    boost::this_thread::sleep (boost::posix_time::microseconds (100));
  }

  return (0);
}

1.3基于属性的分割方法

该方法是基于点云数据的属性的一种鲁棒性较好的分割方法，这种方法一般包括了两个单独的步骤：

第一步，基于属性的计算。

第二步，将根据计算点的属性进行聚类，这种聚类方法一般能适应空间关系和点云的各种属性，最终将不同的属性的点云分割出来，但是这种方法局限性在于他们高度依赖派生属性的质量所以要求第一步能够精确的计算点云数据的属性，这样才会在第二步中根据属性的类别分割出最佳的效果。

论文[5]则是这种方法实现的，提出了一种基于特征空间聚类分析方法，在该方法中，使用一种自适应斜率的邻域系统导出法向量，使用点云数据的属性，例如距离，点密度，点在水平或者垂直方向的分布，来定义测量点之间的领域，然后将每个方向上的法向量的斜率和点邻域的数据之差作为聚类的属性，这种方法可以消除异常值和噪声的影响，基于属性的方法是将点云分割相同属性区域的高效方法，并且分割的结果灵活而准确。然而，这些方法依赖于点之间邻域的定义和点云数据的点密度。当处理大量输入点的多维属性时，这种方法的另一个限制是比较耗时。

论文[5]Segmentation of airborne data using a slope adaptive filter未能找到原文。

目前使用的关于属性的分割方法有欧式聚类和密度聚类：

1.3.1欧式聚类

https://www.cnblogs.com/li-yao7758258/p/6602342.html

对于欧式聚类来说，距离判断准则为欧氏距离。对于空间某点P，通过KD-Tree近邻搜索算法找到k个离p点最近的点，这些点中距离小于设定阈值的便聚类到集合Q中。如果Q中元素的数目不在增加，整个聚类过程便结束；否则须在集合Q中选取p点以外的点，重复上述过程，直到Q中元素的数目不再增加为止。

在pcl中，欧式聚类的代码如下：

实现原理：

pcl::extractEuclideanClusters (const PointCloud &cloud,
                               const typename search::Search::Ptr &tree,
                               float tolerance, std::vector &clusters,
                               unsigned int min_pts_per_cluster,
                               unsigned int max_pts_per_cluster)
{
  if (tree->getInputCloud ()->points.size () != cloud.points.size ())   // 点数量检查
  {
    PCL_ERROR ("[pcl::extractEuclideanClusters] Tree built for a different point cloud dataset (%lu) than the input cloud (%lu)!\n", tree->getInputCloud ()->points.size (), cloud.points.size ());
    return;
  }
  // Check if the tree is sorted -- if it is we don't need to check the first element
  int nn_start_idx = tree->getSortedResults () ? 1 : 0;
  // Create a bool vector of processed point indices, and initialize it to false
  std::vector processed (cloud.points.size (), false);
  std::vector nn_indices;
  std::vector nn_distances;   // 定义需要的变量


  // Process all points in the indices vector
  for (int i = 0; i < static_cast (cloud.points.size ()); ++i)   //遍历点云中的每一个点
  {
    if (processed[i])   //如果该点已经处理则跳过
      continue;

    std::vector seed_queue;   //定义一个种子队列

    int sq_idx = 0;
    seed_queue.push_back (i);  //加入一个种子


    processed[i] = true;

    while (sq_idx < static_cast (seed_queue.size ()))    //遍历每一个种子
    {
      // Search for sq_idx  kdtree 树的近邻搜索
      if (!tree->radiusSearch (seed_queue[sq_idx], tolerance, nn_indices, nn_distances))  
      {
        sq_idx++;
        continue;   //没找到近邻点就继续
      }

      for (size_t j = nn_start_idx; j < nn_indices.size (); ++j)             // can't assume sorted (default isn't!)
      {
        if (nn_indices[j] == -1 || processed[nn_indices[j]])        // Has this point been processed before ?
          continue;   // 种子点的近邻点中如果已经处理就跳出此次循环继续

        // Perform a simple Euclidean clustering
        seed_queue.push_back (nn_indices[j]);   //将此种子点的临近点作为新的种子点。入队操作
        processed[nn_indices[j]] = true;  // 该点已经处理，打标签
      }

      sq_idx++;
    }

    // If this queue is satisfactory, add to the clusters    最大点数和最小点数的类过滤
    if (seed_queue.size () >= min_pts_per_cluster && seed_queue.size () <= max_pts_per_cluster)
    {
      pcl::PointIndices r;
      r.indices.resize (seed_queue.size ());
      for (size_t j = 0; j < seed_queue.size (); ++j)
        r.indices[j] = seed_queue[j];

      // These two lines should not be needed: (can anyone confirm?) -FF
      std::sort (r.indices.begin (), r.indices.end ());
      r.indices.erase (std::unique (r.indices.begin (), r.indices.end ()), r.indices.end ());

      r.header = cloud.header;
      clusters.push_back (r);   // We could avoid a copy by working directly in the vector
    }
  }
}

调用方式及可视化：

// Creating the KdTree object for the search method of the extraction
  pcl::search::KdTree::Ptr tree (new pcl::search::KdTree);
  tree->setInputCloud (cloud_filtered);
 
  std::vector cluster_indices;
  pcl::EuclideanClusterExtraction ec;
  ec.setClusterTolerance (0.02); //设置近邻搜索的搜索半径为2cm
  ec.setMinClusterSize (100);    //设置一个聚类需要的最少点数目为100
  ec.setMaxClusterSize (25000);  //设置一个聚类需要的最大点数目为25000
  ec.setSearchMethod (tree);     //设置点云的搜索机制
  ec.setInputCloud (cloud_filtered); //设置原始点云 
  ec.extract (cluster_indices);  //从点云中提取聚类
 
                                 // 可视化部分
  pcl::visualization::PCLVisualizer viewer("segmention");
  // 我们将要使用的颜色
  float bckgr_gray_level = 0.0;  // 黑色
  float txt_gray_lvl = 1.0 - bckgr_gray_level;
  int num = cluster_indices.size();
 
  int j = 0;
  for (std::vector::const_iterator it = cluster_indices.begin (); it != cluster_indices.end (); ++it)
  {
    pcl::PointCloud::Ptr cloud_cluster (new pcl::PointCloud);
    for (std::vector::const_iterator pit = it->indices.begin (); pit != it->indices.end (); pit++)
      cloud_cluster->points.push_back (cloud_filtered->points[*pit]); //*
    cloud_cluster->width = cloud_cluster->points.size ();
    cloud_cluster->height = 1;
    cloud_cluster->is_dense = true;
 
    std::cout << "PointCloud representing the Cluster: " << cloud_cluster->points.size () << " data points." << std::endl;
    std::stringstream ss;
    ss << "cloud_cluster_" << j << ".pcd";
    writer.write (ss.str (), *cloud_cluster, false); //*
 
    pcl::visualization::PointCloudColorHandlerCustom cloud_in_color_h(cloud,
        color_bar[j][0],
        color_bar[j][1],
        color_bar[j][2]);//赋予显示点云的颜色
    viewer.addPointCloud(cloud_cluster, cloud_in_color_h, std::to_string(j));
 
    j++;
  }

1.3.2条件欧式聚类

条件欧几里德分割的工作方式与（1）所示的标准的欧几里德分割方法基本一样，条件分割除了要距离检查，点还需要满足一个特殊的，可以自定义的要求的限制，这样它们被添加到一个集群。此条件是用户指定的。它归结为如下：对于每一对点（第一个点，作为种子点，第二个点，作为候选点，是一个临近点的选择与计算比较等操作）将会有自定义函数。此函数具有一定的特性：这两个点具有副本，以便我们可以执行我们自己的选择计算的平方距离函数，并返回布尔值。如果值为true，则可以将候选添加到群集。如果假，它不会被添加，即使它通过距离检查。

其主要的缺点：该算法没有初始化种子系统，没有过度分割或者分割不足的控制，还有就是从主循环运算中调用条件函数时，效率比较低。

pcl中提供的函数调用方式：

#include 
#include 

#include 

// 如果这个函数返回的是真，这这个候选点将会被加入聚类中
bool
customCondition(const pcl::PointXYZ& seedPoint, const pcl::PointXYZ& candidatePoint, float squaredDistance)
{
    // Do whatever you want here.做你想做的条件的筛选
    if (candidatePoint.y < seedPoint.y)  //如果候选点的Y的值小于种子点的Y值（就是之前被选择为聚类的点），则不满足条件，返回假
        return false;

    return true;
}

int
main(int argc, char** argv)
{
    pcl::PointCloud::Ptr cloud(new pcl::PointCloud);

    if (pcl::io::loadPCDFile(argv[1], *cloud) != 0)
    {
        return -1;
    }

    pcl::ConditionalEuclideanClustering clustering;
    clustering.setClusterTolerance(0.02);
    clustering.setMinClusterSize(100);
    clustering.setMaxClusterSize(25000);
    clustering.setInputCloud(cloud);
        //设置每次检测一对点云时的函数
    clustering.setConditionFunction(&customCondition);
    std::vector clusters;
    clustering.segment(clusters);

    int currentClusterNum = 1;
    for (std::vector::const_iterator i = clusters.begin(); i != clusters.end(); ++i)
    {
        pcl::PointCloud::Ptr cluster(new pcl::PointCloud);
        for (std::vector::const_iterator point = i->indices.begin(); point != i->indices.end(); point++)
            cluster->points.push_back(cloud->points[*point]);
        cluster->width = cluster->points.size();
        cluster->height = 1;
        cluster->is_dense = true;

        if (cluster->points.size() <= 0)
            break;
        std::cout << "Cluster " << currentClusterNum << " has " << cluster->points.size() << " points." << std::endl;
        std::string fileName = "cluster" + boost::to_string(currentClusterNum) + ".pcd";
        pcl::io::savePCDFileASCII(fileName, *cluster);

        currentClusterNum++;
    }
}

1.3.3密度聚类

https://blog.csdn.net/jjmjim/article/details/107064271

1.4基于模型的分割方法

该方法时基于几何的形状比如球形，圆锥，平面和圆柱形来对点云进行分组，那么根据这些几个形状，具有相同的数学表示的点将会被分割为同一组点，论文[6]中引入了一种众所周知的算法RANSAC（RANdom SAmple Consensus），RANSAC是强大的模型，用于检测直线，圆等数学特征，这种应用极为广泛且可以认为是模型拟合的最先进技术，在3D点云的分割中需要改进的方法都是继承了这种方法。基于模型的方法具有纯粹的数学原理，快速且强大，具有异值性，这种方法的主要局限性在于处理不同点云是的不准确性。这种方法在点云库中已经实现了基于线，平面，圆等各种模型。

论文[6]在ACM官网上，名为Random sample consensus: a paradigm for model fitting with applications to image analysis and automated cartography Communications of the ACM

pcl中提供的关于ransac进行分割的方法可以参考4.2。

1.5基于图割的分割方法

图优化的方法在机器人的应用中十分流行，众所周知的方法是FH算法[7]，该方法简单且高效，并且像Kruskal算法一样用于在图中查找最小生成树。许多基于图的方法的工作被投入到概率推理模型中，例如条件随机场（CRF），使用CRF标记具有不同几何表面基元的点的方法。基于图优化的方法在复杂的城市环境中成功地分割点云，具有接近实时的性能。为了与其他方法进行比较，基于图形的方法可以对点云数据中的复杂场景进但是，这些方法通常无法实时运行。其中一些可能需要离线训练等步骤。

[7]Efficient Graph-Based Image Segmentation——一种快速的基于图的图像分割方法

1.6关于pcl中提供的关于点云分割的模块

1.6.1最小分割算法

https://blog.csdn.net/lemonxiaoxiao/article/details/106211294

原理：

最小割算法是图论中的一个概念，其作用是以某种方式，将两个点分开，当然这两个点中间可能是通过无数的点再相连的。如下图所示：

如果要分开最左边的点和最右边的点，红绿两种割法都是可行的，但是红线跨过了三条线，绿线只跨过了两条。单从跨线数量上来论可以得出绿线这种切割方法更优的结论。但假设线上有不同的权值，那么最优切割则和权值有关了。它到底是怎么找到那条绿线的暂且不论。总而言之，就是有那么一个算法，当你给出了点之间的 “图” （广义的），以及连线的权值时，最小割算法就能按照你的要求把图分开。

显而易见，切割有两个非常重要的因素，第一个是获得点与点之间的拓扑关系，也就是生成一张“图”。第二个是给图中的连线赋予合适的权值。只要这两个要素合适，最小割算法就会办好剩下的事情。

算法流程：

连接算法如下：1.找到每个点最近的n个点；

1.将这n个点和父点连接；

2.找到距离最小的两个块（A块中某点与B块中某点距离最小），并连接；

3.重复3，直至只剩一个块。

可以用点与点之间的欧式距离来构造权值。所有线的权值可映射为线长的函数：

目标需要人为指定(center)，尺寸需要提前给出(radius)。

让除此对象之外的物体被保护起来，不受到打击。保护的方法就是人为加重目标范围之外的权值（罚函数）：

pcl中对最小分割算法的实现：

//生成分割器
pcl::MinCutSegmentation seg;
//分割输入分割目标
seg.setInputCloud(cloud);
//指定打击目标（目标点）
pcl::PointCloud::Ptr foreground_points(new pcl::PointCloud ());
pcl::PointXYZ point;
point.x = 68.97;
point.y = -18.55;
point.z = 0.57;
foreground_points->points.push_back(point);
seg.setForegroundPoints(foreground_points);
//指定权函数sigma
seg.setSigma(0.25);
//物体大概范围
seg.setRadius(3.0433856);
//用多少生成图
seg.setNumberOfNeighbours(14);
//和目标点相连点的权值(至少有14个）
seg.setSourceWeight(0.8);
//分割结果
std::vector  clusters;
seg.extract(clusters);

1.6.2超体聚类

https://www.cnblogs.com/ironstark/p/5013968.html

原理：

超体聚类实际上是一种特殊的区域生长算法，和无限制的生长不同，超体聚类首先需要规律的布置区域生长“晶核”。晶核在空间中实际上是均匀分布的,并指定晶核距离（Rseed)。再指定粒子距离(Rvoxel)。再指定最小晶粒(MOV)，过小的晶粒需要融入最近的大晶粒。超体聚类之前，必须以八叉树对点云进行划分，获得不同点团之间的邻接关系。

pcl中的超体聚类调用方式及可视化

#include 
#include 
#include 
#include 
#include 
#include 

//VTK include needed for drawing graph lines
#include 

// 数据类型
typedef pcl::PointXYZRGBA PointT;
typedef pcl::PointCloud PointCloudT;
typedef pcl::PointNormal PointNT;
typedef pcl::PointCloud PointNCloudT;
typedef pcl::PointXYZL PointLT;
typedef pcl::PointCloud PointLCloudT;

//可视化
void addSupervoxelConnectionsToViewer (PointT &supervoxel_center,
                                       PointCloudT &adjacent_supervoxel_centers,
                                       std::string supervoxel_name,
                                       boost::shared_ptr & viewer);


int
main (int argc, char ** argv)
{
//解析命令行
  if (argc < 2)
  {
    pcl::console::print_error ("Syntax is: %s  \n "
                                "--NT Dsables the single cloud transform \n"
                                "-v \n-s \n"
                                "-c  \n-z  \n"
                                "-n \n", argv[0]);
    return (1);
  }

  //打开点云
  PointCloudT::Ptr cloud = boost::shared_ptr  (new PointCloudT ());
  pcl::console::print_highlight ("Loading point cloud...\n");
  if (pcl::io::loadPCDFile (argv[1], *cloud))
  {
    pcl::console::print_error ("Error loading cloud file!\n");
    return (1);
  }


  bool disable_transform = pcl::console::find_switch (argc, argv, "--NT");

  float voxel_resolution = 0.008f;  //分辨率
  bool voxel_res_specified = pcl::console::find_switch (argc, argv, "-v");
  if (voxel_res_specified)
    pcl::console::parse (argc, argv, "-v", voxel_resolution);

  float seed_resolution = 0.1f;
  bool seed_res_specified = pcl::console::find_switch (argc, argv, "-s");
  if (seed_res_specified)
    pcl::console::parse (argc, argv, "-s", seed_resolution);

  float color_importance = 0.2f;
  if (pcl::console::find_switch (argc, argv, "-c"))
    pcl::console::parse (argc, argv, "-c", color_importance);

  float spatial_importance = 0.4f;
  if (pcl::console::find_switch (argc, argv, "-z"))
    pcl::console::parse (argc, argv, "-z", spatial_importance);

  float normal_importance = 1.0f;
  if (pcl::console::find_switch (argc, argv, "-n"))
    pcl::console::parse (argc, argv, "-n", normal_importance);

//如何使用SupervoxelClustering函数
  pcl::SupervoxelClustering super (voxel_resolution, seed_resolution);
  if (disable_transform)//如果设置的是参数--NT  就用默认的参数
  super.setUseSingleCameraTransform (false);
  super.setInputCloud (cloud);
  super.setColorImportance (color_importance); //0.2f
  super.setSpatialImportance (spatial_importance); //0.4f
  super.setNormalImportance (normal_importance); //1.0f

  std::map ::Ptr > supervoxel_clusters;

  pcl::console::print_highlight ("Extracting supervoxels!\n");
  super.extract (supervoxel_clusters);
  pcl::console::print_info ("Found %d supervoxels\n", supervoxel_clusters.size ());

  boost::shared_ptr viewer (new pcl::visualization::PCLVisualizer ("3D Viewer"));
  viewer->setBackgroundColor (0, 0, 0);

  PointCloudT::Ptr voxel_centroid_cloud = super.getVoxelCentroidCloud ();//获得体素中心的点云
  viewer->addPointCloud (voxel_centroid_cloud, "voxel centroids");
  viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_POINT_SIZE,2.0, "voxel centroids");     //渲染点云
  viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_OPACITY,0.95, "voxel centroids");

  PointLCloudT::Ptr labeled_voxel_cloud = super.getLabeledVoxelCloud ();
  viewer->addPointCloud (labeled_voxel_cloud, "labeled voxels");
  viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_OPACITY,0.8, "labeled voxels");

  PointNCloudT::Ptr sv_normal_cloud = super.makeSupervoxelNormalCloud (supervoxel_clusters);

  //We have this disabled so graph is easy to see, uncomment to see supervoxel normals
  //viewer->addPointCloudNormals (sv_normal_cloud,1,0.05f, "supervoxel_normals");

  pcl::console::print_highlight ("Getting supervoxel adjacency\n");

  std::multimap supervoxel_adjacency;
  super.getSupervoxelAdjacency (supervoxel_adjacency);
  //To make a graph of the supervoxel adjacency, we need to iterate through the supervoxel adjacency multimap
  //为了使整个超体形成衣服图，我们需要遍历超体的每个临近的个体
  std::multimap::iterator label_itr = supervoxel_adjacency.begin ();
  for ( ; label_itr != supervoxel_adjacency.end (); )
  {
    //First get the label
    uint32_t supervoxel_label = label_itr->first;
    //Now get the supervoxel corresponding to the label
    pcl::Supervoxel::Ptr supervoxel = supervoxel_clusters.at (supervoxel_label);

    //Now we need to iterate through the adjacent supervoxels and make a point cloud of them
    PointCloudT adjacent_supervoxel_centers;
    std::multimap::iterator adjacent_itr = supervoxel_adjacency.equal_range (supervoxel_label).first;
    for ( ; adjacent_itr!=supervoxel_adjacency.equal_range (supervoxel_label).second; ++adjacent_itr)
    {
      pcl::Supervoxel::Ptr neighbor_supervoxel = supervoxel_clusters.at (adjacent_itr->second);
      adjacent_supervoxel_centers.push_back (neighbor_supervoxel->centroid_);
    }
    //Now we make a name for this polygon
    std::stringstream ss;
    ss << "supervoxel_" << supervoxel_label;
    //This function is shown below, but is beyond the scope of this tutorial - basically it just generates a "star" polygon mesh from the points given
//从给定的点云中生成一个星型的多边形，
    addSupervoxelConnectionsToViewer (supervoxel->centroid_, adjacent_supervoxel_centers, ss.str (), viewer);
    //Move iterator forward to next label
    label_itr = supervoxel_adjacency.upper_bound (supervoxel_label);
  }

  while (!viewer->wasStopped ())
  {
    viewer->spinOnce (100);
  }
  return (0);
}

//VTK可视化构成的聚类图
void
addSupervoxelConnectionsToViewer (PointT &supervoxel_center,
                                  PointCloudT &adjacent_supervoxel_centers,
                                  std::string supervoxel_name,
                                  boost::shared_ptr & viewer)
{
  vtkSmartPointer points = vtkSmartPointer::New ();
  vtkSmartPointer cells = vtkSmartPointer::New ();
  vtkSmartPointer polyLine = vtkSmartPointer::New ();

  //Iterate through all adjacent points, and add a center point to adjacent point pair
  PointCloudT::iterator adjacent_itr = adjacent_supervoxel_centers.begin ();
  for ( ; adjacent_itr != adjacent_supervoxel_centers.end (); ++adjacent_itr)
  {
    points->InsertNextPoint (supervoxel_center.data);
    points->InsertNextPoint (adjacent_itr->data);
  }
  // Create a polydata to store everything in
  vtkSmartPointer polyData = vtkSmartPointer::New ();
  // Add the points to the dataset
  polyData->SetPoints (points);
  polyLine->GetPointIds  ()->SetNumberOfIds(points->GetNumberOfPoints ());
  for(unsigned int i = 0; i < points->GetNumberOfPoints (); i++)
    polyLine->GetPointIds ()->SetId (i,i);
  cells->InsertNextCell (polyLine);
  // Add the lines to the dataset
  polyData->SetLines (cells);
  viewer->addModelFromPolyData (polyData,supervoxel_name);
}

1.6.3基于凹凸性的点云分割

https://www.cnblogs.com/ironstark/p/5027269.html

1.6.3.1lccp方法

算法原理：

超体聚类存在明显的过分割现象。点云完成超体聚类之后，对于过分割的点云需要计算不同的块之间凹凸关系。凹凸关系通过 CC（Extended Convexity Criterion）和 SC （Sanity criterion）判据来进行判断。其中 CC 利用相邻两片中心连线向量与法向量夹角来判断两片是凹是凸。显然，如果图中a1>a2则为凹，反之则为凸。

考虑到测量噪声等因素，需要在实际使用过程中引入门限值（a1需要比a2大出一定量）来滤出较小的凹凸误判。此外，为去除一些小噪声引起的误判，还需要引入“第三方验证”，如果某块和相邻两块都相交，则其凹凸关系必相同。CC 判据最终如CCe：

如果相邻两面中，有一个面是单独的，cc判据是无法将其分开的。举个简单的例子，两本厚度不同的书并排放置，视觉算法应该将两本书分割开。如果是台阶，则视觉算法应该将台阶作为一个整体。本质上就是因为厚度不同的书存在surface-singularities。为此需要引入SC判据，来对此进行区分。

如图所示，相邻两面是否真正联通，是否存在单独面，与θ角有关，θ角越大，则两面真的形成凸关系的可能性就越大。据此，可以设计SC判据：

其中S(向量)为两平面法向量的叉积。

　　最终，两相邻面之间凸边判据为：

在标记完各个小区域的凹凸关系后，则采用区域增长算法将小区域聚类成较大的物体。此区域增长算法受到小区域凹凸性限制，即只允许区域跨越凸边增长。

//生成LCCP分割器
pcl::LCCPSegmentation::LCCPSegmentation LCCPseg;
//输入超体聚类结果
seg.setInputSupervoxels(supervoxel_clusters,supervoxel_adjacency);
//CC效验beta值
seg.setConcavityToleranceThreshold (concavity_tolerance_threshold);
//CC效验的k邻点
seg.setKFactor (k_factor_arg)
//
seg.setSmoothnessCheck (bool_use_smoothness_check_arg,voxel_res_arg,seed_res_arg,smoothness_threshold_arg = 0.1);
//SC效验
seg.setSanityCheck (bool_use_sanity_criterion_arg);
//最小分割尺寸
seg.setMinSegmentSize (min_segment_size_arg)

seg.segment();
seg.relabelCloud (pcl::PointCloud &labeled_cloud_arg);

1.6.3.2cpc方法

出自论文：*Constrained Planar Cuts - Object Partitioning for Point Clouds* 。

算法原理：

和其他基于凹凸性的方法相同，本方法也需要先进行超体聚类。在完成超体聚类之后，采用和LCCP相同的凹凸性判据获得各个块之间的凹凸关系。在获得凹凸性之后，CPC方法所采取的措施是不同的。其操作称为半全局分割

　　在分割之前，首先需要生成 EEC(Euclidean edge cloud)， EEC的想法比较神奇，因为凹凸性定义在相邻两个”片“上，换言之，定义在连接相邻两“片”的edge上。将每个edge抽象成一个点云，则得到了附带凹凸信息的点云。如图所示，左图是普通点云，但附带了邻接和凹凸信息。右边是EEC，对凹边赋权值1，其他为0。

　显而易见，某处如果蓝色的点多，那么就越凹，就越应该切开（所谓切开实际上是用平面划分）。问题就转化为利用蓝点求平面了。利用点云求一个最可能的平面当然需要请出我们的老朋友 RanSaC . 但此处引入一个评价函数，用于评价此次分割的优良程度Sm,**Pm 是EEC中的点.

单纯的weighted RanSac算法并不够。其会导致对某些图形的错误分割，所以作者对此做了第一次“修补".错误的分割如下图所示

方法的原理很简单，垂直于凹边表面的点具有更高的权重,显然，对于EEC中的凹点，只要取其少量邻点即可估计垂直方向。

　　这种修补后还有一个问题，如果这个分割面过长的情况下，有可能会误伤。如图所示：

这种修补方法的原理就更加简单粗暴了，对凹点先进行欧式分割（限制增长上限），之后再分割所得的子域里进行分割。

　　在修修补补之后，CPC算法终于可以投入使用了，从测试集的结果来看，效果还是很好的。

pcl中cpc算法的调用：

//生成CPC分割器
pcl::CPCSegmentation::CPCSegmentation seg;
//输入超体聚类结果
seg.setInputSupervoxels(supervoxel_clusters,supervoxel_adjacency);
//设置分割参数
setCutting (max_cuts = 20,
                 cutting_min_segments = 0,
                 cutting_min_score = 0.16,
                 locally_constrained = true,
                 directed_cutting = true,
                 clean_cutting = false)；
seg.setRANSACIterations (ransac_iterations);
seg.segment();
seg.relabelCloud (pcl::PointCloud &labeled_cloud_arg);

1.6.4基于形态学的点云分割

本方法主要用于航空测量中。

　出自论文：　A Progressive Morphological Filter for Removing Nonground Measurements From Airborne LIDAR Data

算法原理：

对于图像而言，形态学运算一般是针对二值图像而言的。当然也有针对灰度的形态学运算，其原理应该和针对点云的形态学运算类似（我猜的）。形态学算子的设计实际上非常简单，只要能设计出基础的膨胀和腐蚀算子就可以组合得到一系列的处理。

　　其中，d表示膨胀算子，e表示腐蚀算子。算子的原理有些像中值滤波，通过选取一个窗w中最高点或最低点来完成图像的膨胀和腐蚀，其效果如下图所示：

　　在航拍图的横截面上可以很清楚的看出膨胀与腐蚀的效果。对于房子和树可以用不同尺度窗（从小到大）先腐蚀至地面。但是这会导致一个巨大的问题。。。如果地面上有个土包（比如秦始皇陵），那么这个土包也会在一次次的腐蚀中被消耗。那岂不是秦始皇陵就发现不了？所以还有一个补偿算法用于解决这个问题，称为线性补偿算法。

　　建筑物和土包有一个巨大的区别，建筑物往往相对比较陡峭，而土包却是变化比较平缓的。这个可以作为一个判据，用于判断物体是否需要被腐蚀，也作为窗收敛的判据。

　　式中k称为斜率，代表下一个窗的大小是上一个窗的2^k倍

　　s是一个因子

　　dh是切深判据，每一次腐蚀大于切深判据才认为是有效的，小于切深判据则是土包。

pcl中本算法的调用：

//生成形态滤波器
  pcl::ProgressiveMorphologicalFilter pmf;
  pmf.setInputCloud (cloud);
  //设置窗的大小以及切深，斜率信息
  pmf.setMaxWindowSize (20);
  pmf.setSlope (1.0f);
  pmf.setInitialDistance (0.5f);
  pmf.setMaxDistance (3.0f);
  //提取地面
  pmf.extract (ground->indices);

  // 从标号到点云
  pcl::ExtractIndices extract;
  extract.setInputCloud (cloud);
  extract.setIndices (ground);
  extract.filter (*cloud_filtered);

1.7总结

一般来说，点云的分割有两种基本方法。

第一种方法使用纯数学模型和几何推理技术，如区域增长或模型拟合，将线性和非线性模型拟合到点云数据。这种方法允许快速运行时间能实现良好的结果。这种方法的局限性在于在拟合物体时难以选择模型的大小，对噪声敏感并且在复杂场景中不能很好地工作。

第二种方法使用特征描述子的方法从点云数据中提取3D特征，并使用机器学习技术来学习不同类别的对象类型，然后使用结果模型对所获取的数据进行分类。在复杂场景中，机器学习技术将优于纯粹基于几何推理的技术。原因是由于噪声，密度不均匀，点云数据中的遮挡，很难找到并将复杂的几何图元拟合到物体上。虽然机器学习技术可以提供更好的结果，但它们通常很慢并且依赖于特征提取过程的结果。

本文仅做学术分享，如有侵权，请联系删文。

—THE END—