一文详解点云分割算法

作者丨书生意封侯@知乎

来源丨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)


    
    
  1. /* @brief:将点云转换为vector,俯视图栅格化
  2. * @param [in]: in_cloud,输入点云
  3. * @param [out]: out_cloud,转换的点云vector
  4. * @return NONE
  5. */
  6. void FrontLidarAlg::convetCloudToVector(pcl::PointCloud<pcl::PointXYZI>::Ptr in_cloud, \
  7. cv::Mat& out_img, \
  8. std::vector<pcl::PointCloud<pcl::PointXYZI>>* out_cloud)
  9. {
  10. out_img = cv::Mat::zeros(cv::Size(img_width_,img_height_), CV_8UC3);
  11. volatile int row;
  12. volatile int col;
  13. for( int i= 0;i<in_cloud->size();i++)
  14. {
  15. //将横向在检测范围之外的去掉
  16. if(in_cloud->points[i].x <= cloud_x_min_ || in_cloud->points[i].x >= cloud_x_max_)
  17. {
  18. continue;
  19. }
  20. //将纵向在检测范围之外的去掉
  21. if(in_cloud->points[i].y >= cloud_y_max_ || in_cloud->points[i].y <= 0)
  22. {
  23. continue;
  24. }
  25. //限制最高和最低点大小
  26. if(in_cloud->points[i].z < min_dect_height_) //高度低于最低点,则赋值为最低点
  27. {
  28. in_cloud->points[i].z = min_dect_height_;
  29. }
  30. if(in_cloud->points[i].z > max_dect_height_) //高度高于最高点,则赋值最高点
  31. {
  32. in_cloud->points[i].z = max_dect_height_;
  33. }
  34. //计算点云所在图像的行数和列数,四舍五入
  35. col = int((in_cloud->points[i].x - cloud_x_min_)/img_res_);
  36. row = int(in_cloud->points[i].y/img_res_);
  37. out_cloud->at(col+row*img_width_).points.push_back(in_cloud->points[i]);
  38. //输入点云高度值转换到图像坐标系下的数值
  39. int value = ( int)((in_cloud->points[i].z - min_dect_height_)/dect_height_res_);
  40. if(out_img.at<cv::Vec3b>(row,col)[ 0] < value) //b-最大高度
  41. {
  42. out_img.at<cv::Vec3b>(row,col)[ 0] = value;
  43. if( 0 == out_img.at<cv::Vec3b>(row,col)[ 1])
  44. {
  45. out_img.at<cv::Vec3b>(row,col)[ 1] = value;
  46. }
  47. }
  48. else if(out_img.at<cv::Vec3b>(row,col)[ 1] > value) //b-最小高度
  49. {
  50. out_img.at<cv::Vec3b>(row,col)[ 1] = value;
  51. }
  52. } //for in_cloud->size
  53. }

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

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

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

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

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


    
    
  1. #include "object_segment.h"
  2. namespace front_lidar {
  3. /* @brief:点云分割构造函数,初始化参数
  4. * @param [in]: NONE
  5. * @param [out]: NONE
  6. * @return NONE
  7. */
  8. ObjectSegment::ObjectSegment( string config_file_str)
  9. {
  10. config_parser_ptr_.reset( new ConfigParser(config_file_str));
  11. use_seg_raw_cloud_ = config_parser_ptr_->getInt( "use_seg_raw_cloud"); //聚类后是否输出原始数据
  12. img_res_ = config_parser_ptr_->getDouble( "img_res"); //图像分辨率
  13. cloud_x_min_ = config_parser_ptr_->getDouble( "cloud_x_min"); //x轴最小值
  14. cloud_x_max_ = config_parser_ptr_->getDouble( "cloud_x_max"); //x轴最大值
  15. cloud_y_min_ = config_parser_ptr_->getDouble( "cloud_y_min"); //y轴最小值
  16. cloud_y_max_ = config_parser_ptr_->getDouble( "cloud_y_max"); //y轴最大值
  17. min_cluster_size_ = config_parser_ptr_->getDouble( "object_detect_min_cluster_size"); //最小聚类个数
  18. max_cluster_size_ = config_parser_ptr_->getDouble( "object_detect_max_cluster_size"); //最大聚类个数
  19. cell_size_ = config_parser_ptr_->getDouble( "cell_size"); //聚类临域尺寸
  20. min_dect_height_ = -10; //最小检测高度
  21. max_dect_height_ = 10; //最大检测高度
  22. dect_height_res_ = (max_dect_height_ - min_dect_height_)/ 255; //检测高度分辨率
  23. region_growing_clusters_.clear();
  24. region_growing_img_.resize( 0); //区域增长后的图像
  25. }
  26. /* @brief:点云分割处理,将点云分割为不同障碍物点云
  27. * @param [in]: in_img,输入图像栅格, in_vector-输入栅格化的点云
  28. * @param [out]: cloud_clusters_ptr,输出聚类结果
  29. * @return NONE
  30. */
  31. void ObjectSegment::process( const std::vector<pcl::PointCloud<pcl::PointXYZI>>* in_vector, \
  32. const cv::Mat& in_img, \
  33. std::vector<pcl::PointCloud<pcl::PointXYZI>::Ptr>* cloud_clusters_ptr)
  34. {
  35. cloud_clusters_ptr->clear();
  36. cluster_num_ = 0;
  37. //区域增长,图像聚类分割
  38. pcl::StopWatch regin_grow_timer;
  39. reginGrowing(in_img, region_growing_img_);
  40. LOG_Debug()<< "regin_grow time:"<<regin_grow_timer.getTime()<< "ms";
  41. //将图像转换为点云
  42. pcl::StopWatch vector_to_cloud_timer;
  43. if( 1 == use_seg_raw_cloud_)
  44. {
  45. convetVectorToCloud(in_vector, region_growing_img_, &region_growing_clusters_);
  46. }
  47. else
  48. {
  49. convetImageToCloud(region_growing_img_, &region_growing_clusters_);
  50. }
  51. LOG_Debug()<< "vector_to_cloud time:"<<vector_to_cloud_timer.getTime()<< "ms";
  52. //滤除聚类点数少于设定阈值的类别
  53. pcl::StopWatch filter_timer;
  54. filterObjectsCloud(region_growing_clusters_, cloud_clusters_ptr);
  55. LOG_Debug()<< "filter time:"<<filter_timer.getTime()<< "ms";
  56. }
  57. /* @brief:区域增长,8邻域聚类
  58. * @param [in]: in_img-输入图像
  59. * @param [out]: out_img,转换的图像,0-b-最大值,1-g-最小值,2-r-有无聚类标志(0-未聚类,1-255-已聚类ID)
  60. * @return NONE
  61. */
  62. void ObjectSegment::reginGrowing( const cv::Mat &in_img, cv::Mat &out_img)
  63. {
  64. unsigned short class_id = 0;
  65. unsigned short element_class_id = 0;
  66. out_img = in_img.clone();
  67. for( int row = 1; row < out_img.rows -1; row++)
  68. {
  69. for( int col = 1; col < out_img.cols - 1; col++)
  70. {
  71. //像素的高度最大值为0,则该像素为空,无效像素
  72. if( 0 == out_img.at<cv::Vec3b>(row,col)[ 0])
  73. {
  74. continue;
  75. }
  76. //像素的类别标记为空,即为未标记类别的像素,则分配类别ID
  77. if( 0 == out_img.at<cv::Vec3b>(row,col)[ 2])
  78. {
  79. if(class_id > 250) //超出限定类别总数,则返回
  80. {
  81. return;
  82. }
  83. out_img.at<cv::Vec3b>(row,col)[ 2] = ++class_id; //给该像素赋值类别ID
  84. }
  85. element_class_id = out_img.at<cv::Vec3b>(row,col)[ 2];
  86. //根据输入种子像素,递归区域增长
  87. elementReginGrowing(out_img, row, col, element_class_id);
  88. } //col
  89. } //row
  90. cluster_num_ = class_id;
  91. }
  92. /* @brief:单个元素区域增长,cell_size_邻域聚类
  93. * @param [in]: 输入图像-in_img,当前元素的r和c,输出图像-in_img
  94. * @param [out]: out_img,转换的图像,0-b-最大值,1-g-最小值,2-r-有无聚类标志(0-未聚类,1-255-聚类ID)
  95. * @return NONE
  96. */
  97. void ObjectSegment::elementReginGrowing(cv::Mat &in_img, int row, int col, unsigned short class_id)
  98. {
  99. int start_row; //起始行
  100. int end_row; //截止行
  101. int start_col; //起始列
  102. int end_col; //截止列
  103. //判断起始行 截止行 起始列 截止列,防止越界
  104. start_row = row - cell_size_;
  105. end_row = row + cell_size_;
  106. start_col = col - cell_size_;
  107. end_col = col + cell_size_;
  108. if(start_row < 0)
  109. start_row = 0;
  110. if(end_row > in_img.rows - 1)
  111. end_row = in_img.rows - 1;
  112. if(start_col < 0)
  113. start_col = 0;
  114. if(end_col > in_img.cols - 1)
  115. end_col = in_img.cols - 1;
  116. for( int m = start_row; m <= end_row; m++)
  117. {
  118. for( int n = start_col; n <= end_col; n++)
  119. {
  120. //该像素对应最大高度为0,则为无效像素
  121. if( 0 == in_img.at<cv::Vec3b>(m,n)[ 0])
  122. {
  123. continue;
  124. }
  125. //该元素为初始输入,已经标记过
  126. if(m==row && n==col)
  127. {
  128. continue;
  129. }
  130. //未标记过的元素,标记该元素,以该元素为种子,进行区域增长
  131. if( 0 == in_img.at<cv::Vec3b>(m,n)[ 2])
  132. {
  133. in_img.at<cv::Vec3b>(m,n)[ 2] = class_id; //标记元素
  134. elementReginGrowing(in_img, m, n, class_id); //以该元素为初始值,进行区域增长
  135. }
  136. } //for col
  137. } //for row
  138. }
  139. /* @brief:将vector转化为聚类结果,2D到3D
  140. * @param [in]: in_img,输入的栅格图像
  141. * @param [out]: out_cloud_ptr,转换的聚类结果
  142. * @return NONE
  143. */
  144. void ObjectSegment::convetImageToCloud( const cv::Mat& in_img, \
  145. std::vector<pcl::PointCloud<pcl::PointXYZI>::Ptr>* out_cloud_ptr)
  146. {
  147. out_cloud_ptr->clear();
  148. pcl::PointXYZI point;
  149. int class_id = 0; //类ID,1开始递增
  150. //根据类别分配向量空间
  151. out_cloud_ptr->resize(cluster_num_); //扩充输出向量大小
  152. //分配点云空间
  153. for( int i= 0;i<cluster_num_;i++)
  154. {
  155. out_cloud_ptr->at(i).reset( new pcl::PointCloud<pcl::PointXYZI>);
  156. }
  157. //将图像转换为点云
  158. for( int row= 0;row<in_img.rows;row++)
  159. {
  160. for( int col= 0;col<in_img.cols;col++)
  161. {
  162. //如果分类标记为0,则为无效数据
  163. if( 0 == in_img.at<cv::Vec3b>(row,col)[ 2])
  164. {
  165. continue;
  166. }
  167. class_id = in_img.at<cv::Vec3b>(row,col)[ 2]; //获取该区域的类别ID
  168. //图像坐标转换为点云坐标
  169. point.x = (col+ 0.5)*img_res_ + cloud_x_min_;
  170. point.y = (row+ 0.5)*img_res_;
  171. //点云坐标高度最大值
  172. point.z = in_img.at<cv::Vec3b>(row,col)[ 0]*dect_height_res_ + min_dect_height_;
  173. out_cloud_ptr->at(class_id -1)->points.push_back(point);
  174. //点云坐标高度最小值
  175. point.z = in_img.at<cv::Vec3b>(row,col)[ 1]*dect_height_res_ + min_dect_height_;
  176. out_cloud_ptr->at(class_id -1)->points.push_back(point);
  177. }
  178. }
  179. }
  180. /* @brief:将vector转化为聚类结果,2D到3D
  181. * @param [in]: in_vector,输入的栅格化后点云,in_img,输入的聚类标志位的图像
  182. * @param [out]: out_cloud_ptr,转换的聚类结果
  183. * @return NONE
  184. */
  185. void ObjectSegment::convetVectorToCloud( const std::vector<pcl::PointCloud<pcl::PointXYZI>>* in_vector, \
  186. const cv::Mat& in_img, \
  187. std::vector<pcl::PointCloud<pcl::PointXYZI>::Ptr>* out_cloud_ptr)
  188. {
  189. out_cloud_ptr->clear();
  190. int class_id = 0; //类ID,1开始递增
  191. //根据类别分配向量空间
  192. out_cloud_ptr->resize(cluster_num_); //扩充输出向量大小
  193. //分配点云空间
  194. for( int i= 0;i<cluster_num_;i++)
  195. {
  196. out_cloud_ptr->at(i).reset( new pcl::PointCloud<pcl::PointXYZI>);
  197. }
  198. //将图像转换为点云
  199. for( int row= 0;row<in_img.rows;row++)
  200. {
  201. for( int col= 0;col<in_img.cols;col++)
  202. {
  203. //如果分类标记为0,则为无效数据
  204. if( 0 == in_img.at<cv::Vec3b>(row,col)[ 2])
  205. {
  206. continue;
  207. }
  208. class_id = in_img.at<cv::Vec3b>(row,col)[ 2]; //获取该区域的类别ID
  209. int count = in_vector->at(col+row*in_img.cols).points.size();
  210. if(count > 25)
  211. count = 25;
  212. for( int j = 0; j < count; j++)
  213. {
  214. out_cloud_ptr->at(class_id -1)->points.push_back(in_vector->at(col+row*in_img.cols).points[j]);
  215. }
  216. }
  217. }
  218. }
  219. /* @brief:滤除聚类后点云数量小于设定值的类别
  220. * @param [in]: in_objects_cloud-聚类后的目标点云
  221. * @param [out]: out_objects_cloud_ptr-滤除类别中点数量少于设定阈值的类别
  222. * @return NONE
  223. */
  224. void ObjectSegment::filterObjectsCloud( const std::vector<pcl::PointCloud<pcl::PointXYZI>::Ptr> in_objects_cloud, \
  225. std::vector<pcl::PointCloud<pcl::PointXYZI>::Ptr> *out_objects_cloud_ptr)
  226. {
  227. out_objects_cloud_ptr->clear();
  228. for( int i= 0;i<in_objects_cloud.size();i++)
  229. {
  230. if(in_objects_cloud[i]->size() < min_cluster_size_)
  231. {
  232. continue;
  233. }
  234. out_objects_cloud_ptr->push_back(in_objects_cloud[i]);
  235. }
  236. }
  237. } //namespace front_lidar

1.2.4基于点云区域增长方法

算法思想:

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

算法步骤:

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

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

3.曲率是否足够小

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

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

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

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

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

算法源码:


    
    
  1. #include <iostream>
  2. #include <vector>
  3. #include <pcl/point_types.h>
  4. #include <pcl/io/pcd_io.h>
  5. #include <pcl/search/search.h>
  6. #include <pcl/search/kdtree.h>
  7. #include <pcl/features/normal_3d.h>
  8. #include <pcl/visualization/cloud_viewer.h>
  9. #include <pcl/filters/passthrough.h>
  10. #include <pcl/segmentation/region_growing.h>
  11. int
  12. main ( int argc, char** argv)
  13. {
  14. //点云的类型
  15. pcl::PointCloud<pcl::PointXYZ>::Ptr cloud ( new pcl::PointCloud<pcl::PointXYZ>);
  16. //打开点云
  17. if ( pcl::io::loadPCDFile <pcl::PointXYZ> ( "region_growing_tutorial.pcd", *cloud) == -1)
  18. {
  19. std::cout << "Cloud reading failed." << std::endl;
  20. return ( -1);
  21. }
  22. //设置搜索的方式或者说是结构
  23. pcl::search::Search<pcl::PointXYZ>::Ptr tree = boost::shared_ptr<pcl::search::Search<pcl::PointXYZ> > ( new pcl::search::KdTree<pcl::PointXYZ>);
  24. //求法线
  25. pcl::PointCloud <pcl::Normal>::Ptr normals ( new pcl::PointCloud <pcl::Normal>);
  26. pcl::NormalEstimation<pcl::PointXYZ, pcl::Normal> normal_estimator;
  27. normal_estimator.setSearchMethod (tree);
  28. normal_estimator.setInputCloud (cloud);
  29. normal_estimator.setKSearch ( 50);
  30. normal_estimator.compute (*normals);
  31. //直通滤波在Z轴的0到1米之间
  32. pcl::IndicesPtr indices ( new std::vector < int>);
  33. pcl::PassThrough<pcl::PointXYZ> pass;
  34. pass.setInputCloud (cloud);
  35. pass.setFilterFieldName ( "z");
  36. pass.setFilterLimits ( 0.0, 1.0);
  37. pass.filter (*indices);
  38. //聚类对象<点,法线>
  39. pcl::RegionGrowing<pcl::PointXYZ, pcl::Normal> reg;
  40. reg.setMinClusterSize ( 50); //最小的聚类的点数
  41. reg.setMaxClusterSize ( 1000000); //最大的
  42. reg.setSearchMethod (tree); //搜索方式
  43. reg.setNumberOfNeighbours ( 30); //设置搜索的邻域点的个数
  44. reg.setInputCloud (cloud); //输入点
  45. //reg.setIndices (indices);
  46. reg.setInputNormals (normals); //输入的法线
  47. reg.setSmoothnessThreshold ( 3.0 / 180.0 * M_PI); //设置平滑度
  48. reg.setCurvatureThreshold ( 1.0); //设置曲率的阀值
  49. std::vector <pcl::PointIndices> clusters;
  50. reg.extract (clusters);
  51. std::cout << "Number of clusters is equal to " << clusters.size () << std::endl;
  52. std::cout << "First cluster has " << clusters[ 0].indices.size () << " points." << endl;
  53. std::cout << "These are the indices of the points of the initial" <<
  54. std::endl << "cloud that belong to the first cluster:" << std::endl;
  55. int counter = 0;
  56. while (counter < clusters[ 0].indices.size ())
  57. {
  58. std::cout << clusters[ 0].indices[counter] << ", ";
  59. counter++;
  60. if (counter % 10 == 0)
  61. std::cout << std::endl;
  62. }
  63. std::cout << std::endl;
  64. //可视化聚类的结果
  65. pcl::PointCloud <pcl::PointXYZRGB>::Ptr colored_cloud = reg.getColoredCloud ();
  66. pcl::visualization::CloudViewer viewer ( "Cluster viewer");
  67. viewer.showCloud(colored_cloud);
  68. while (!viewer.wasStopped ())
  69. {
  70. }
  71. return ( 0);
  72. }

1.2.5基于颜色区域增长方法

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

算法步骤:

算法分为两步:

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

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

pcl中的函数调用:


    
    
  1. #include <iostream>
  2. #include <vector>
  3. #include <pcl/point_types.h>
  4. #include <pcl/io/pcd_io.h>
  5. #include <pcl/search/search.h>
  6. #include <pcl/search/kdtree.h>
  7. #include <pcl/visualization/cloud_viewer.h>
  8. #include <pcl/filters/passthrough.h>
  9. #include <pcl/segmentation/region_growing_rgb.h>
  10. int
  11. main ( int argc, char** argv)
  12. {
  13. pcl::search::Search <pcl::PointXYZRGB>::Ptr tree = boost::shared_ptr<pcl::search::Search<pcl::PointXYZRGB> > ( new pcl::search::KdTree<pcl::PointXYZRGB>);
  14. pcl::PointCloud <pcl::PointXYZRGB>::Ptr cloud ( new pcl::PointCloud <pcl::PointXYZRGB>);
  15. if ( pcl::io::loadPCDFile <pcl::PointXYZRGB> ( "region_growing_rgb_tutorial.pcd", *cloud) == -1 )
  16. {
  17. std::cout << "Cloud reading failed." << std::endl;
  18. return ( -1);
  19. }
  20. //存储点云的容器
  21. pcl::IndicesPtr indices ( new std::vector < int>);
  22. //滤波
  23. pcl::PassThrough<pcl::PointXYZRGB> pass;
  24. pass.setInputCloud (cloud);
  25. pass.setFilterFieldName ( "z");
  26. pass.setFilterLimits ( 0.0, 1.0);
  27. pass.filter (*indices);
  28. //基于颜色的区域生成的对象
  29. pcl::RegionGrowingRGB<pcl::PointXYZRGB> reg;
  30. reg.setInputCloud (cloud);
  31. reg.setIndices (indices); //点云的索引
  32. reg.setSearchMethod (tree);
  33. reg.setDistanceThreshold ( 10); //距离的阀值
  34. reg.setPointColorThreshold ( 6); //点与点之间颜色容差
  35. reg.setRegionColorThreshold ( 5); //区域之间容差
  36. reg.setMinClusterSize ( 600); //设置聚类的大小
  37. std::vector <pcl::PointIndices> clusters;
  38. reg.extract (clusters);
  39. pcl::PointCloud <pcl::PointXYZRGB>::Ptr colored_cloud = reg.getColoredCloud ();
  40. pcl::visualization::CloudViewer viewer ( "Cluster viewer");
  41. viewer.showCloud (colored_cloud);
  42. while (!viewer.wasStopped ())
  43. {
  44. boost::this_thread::sleep (boost::posix_time::microseconds ( 100));
  45. }
  46. return ( 0);
  47. }

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

401a3385c05d72b2687cdc0b201a71f2.png

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

在pcl中,欧式聚类的代码如下:

实现原理:


    
    
  1. pcl::extractEuclideanClusters ( const PointCloud<PointT> &cloud,
  2. const typename search::Search<PointT>::Ptr &tree,
  3. float tolerance, std::vector<PointIndices> &clusters,
  4. unsigned int min_pts_per_cluster,
  5. unsigned int max_pts_per_cluster)
  6. {
  7. if (tree->getInputCloud ()->points.size () != cloud.points.size ()) // 点数量检查
  8. {
  9. 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 ());
  10. return;
  11. }
  12. // Check if the tree is sorted -- if it is we don't need to check the first element
  13. int nn_start_idx = tree->getSortedResults () ? 1 : 0;
  14. // Create a bool vector of processed point indices, and initialize it to false
  15. std::vector< bool> processed (cloud.points.size (), false);
  16. std::vector< int> nn_indices;
  17. std::vector<float> nn_distances; // 定义需要的变量
  18. // Process all points in the indices vector
  19. for ( int i = 0; i < static_cast< int> (cloud.points.size ()); ++i) //遍历点云中的每一个点
  20. {
  21. if (processed[i]) //如果该点已经处理则跳过
  22. continue;
  23. std::vector< int> seed_queue; //定义一个种子队列
  24. int sq_idx = 0;
  25. seed_queue.push_back (i); //加入一个种子
  26. processed[i] = true;
  27. while (sq_idx < static_cast< int> (seed_queue.size ())) //遍历每一个种子
  28. {
  29. // Search for sq_idx kdtree 树的近邻搜索
  30. if (!tree->radiusSearch (seed_queue[sq_idx], tolerance, nn_indices, nn_distances))
  31. {
  32. sq_idx++;
  33. continue; //没找到近邻点就继续
  34. }
  35. for (size_t j = nn_start_idx; j < nn_indices.size (); ++j) // can't assume sorted (default isn't!)
  36. {
  37. if (nn_indices[j] == -1 || processed[nn_indices[j]]) // Has this point been processed before ?
  38. continue; // 种子点的近邻点中如果已经处理就跳出此次循环继续
  39. // Perform a simple Euclidean clustering
  40. seed_queue.push_back (nn_indices[j]); //将此种子点的临近点作为新的种子点。入队操作
  41. processed[nn_indices[j]] = true; // 该点已经处理,打标签
  42. }
  43. sq_idx++;
  44. }
  45. // If this queue is satisfactory, add to the clusters 最大点数和最小点数的类过滤
  46. if (seed_queue.size () >= min_pts_per_cluster && seed_queue.size () <= max_pts_per_cluster)
  47. {
  48. pcl::PointIndices r;
  49. r.indices.resize (seed_queue.size ());
  50. for (size_t j = 0; j < seed_queue.size (); ++j)
  51. r.indices[j] = seed_queue[j];
  52. // These two lines should not be needed: (can anyone confirm?) -FF
  53. std::sort (r.indices.begin (), r.indices.end ());
  54. r.indices.erase (std::unique (r.indices.begin (), r.indices.end ()), r.indices.end ());
  55. r.header = cloud.header;
  56. clusters.push_back (r); // We could avoid a copy by working directly in the vector
  57. }
  58. }
  59. }

调用方式及可视化:


    
    
  1. // Creating the KdTree object for the search method of the extraction
  2. pcl::search::KdTree<pcl::PointXYZ>::Ptr tree ( new pcl::search::KdTree<pcl::PointXYZ>);
  3. tree->setInputCloud (cloud_filtered);
  4. std::vector<pcl::PointIndices> cluster_indices;
  5. pcl::EuclideanClusterExtraction<pcl::PointXYZ> ec;
  6. ec.setClusterTolerance ( 0.02); //设置近邻搜索的搜索半径为2cm
  7. ec.setMinClusterSize ( 100); //设置一个聚类需要的最少点数目为100
  8. ec.setMaxClusterSize ( 25000); //设置一个聚类需要的最大点数目为25000
  9. ec.setSearchMethod (tree); //设置点云的搜索机制
  10. ec.setInputCloud (cloud_filtered); //设置原始点云
  11. ec.extract (cluster_indices); //从点云中提取聚类
  12. // 可视化部分
  13. pcl::visualization::PCLVisualizer viewer( "segmention");
  14. // 我们将要使用的颜色
  15. float bckgr_gray_level = 0.0; // 黑色
  16. float txt_gray_lvl = 1.0 - bckgr_gray_level;
  17. int num = cluster_indices.size();
  18. int j = 0;
  19. for (std::vector<pcl::PointIndices>::const_iterator it = cluster_indices.begin (); it != cluster_indices.end (); ++it)
  20. {
  21. pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_cluster ( new pcl::PointCloud<pcl::PointXYZ>);
  22. for (std::vector< int>::const_iterator pit = it->indices.begin (); pit != it->indices.end (); pit++)
  23. cloud_cluster->points.push_back (cloud_filtered->points[*pit]); //*
  24. cloud_cluster->width = cloud_cluster->points.size ();
  25. cloud_cluster->height = 1;
  26. cloud_cluster->is_dense = true;
  27. std::cout << "PointCloud representing the Cluster: " << cloud_cluster->points.size () << " data points." << std::endl;
  28. std::stringstream ss;
  29. ss << "cloud_cluster_" << j << ".pcd";
  30. writer.write<pcl::PointXYZ> (ss.str (), *cloud_cluster, false); //*
  31. pcl::visualization::PointCloudColorHandlerCustom<pcl::PointXYZ> cloud_in_color_h(cloud,
  32. color_bar[j][ 0],
  33. color_bar[j][ 1],
  34. color_bar[j][ 2]); //赋予显示点云的颜色
  35. viewer.addPointCloud(cloud_cluster, cloud_in_color_h, std::to_string(j));
  36. j++;
  37. }

1.3.2条件欧式聚类

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

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

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


    
    
  1. #include <pcl/io/pcd_io.h>
  2. #include <pcl/segmentation/conditional_euclidean_clustering.h>
  3. #include <iostream>
  4. // 如果这个函数返回的是真,这这个候选点将会被加入聚类中
  5. bool
  6. customCondition( const pcl::PointXYZ& seedPoint, const pcl::PointXYZ& candidatePoint, float squaredDistance)
  7. {
  8. // Do whatever you want here.做你想做的条件的筛选
  9. if (candidatePoint.y < seedPoint.y) //如果候选点的Y的值小于种子点的Y值(就是之前被选择为聚类的点),则不满足条件,返回假
  10. return false;
  11. return true;
  12. }
  13. int
  14. main( int argc, char** argv)
  15. {
  16. pcl::PointCloud<pcl::PointXYZ>::Ptr cloud( new pcl::PointCloud<pcl::PointXYZ>);
  17. if (pcl::io::loadPCDFile<pcl::PointXYZ>(argv[ 1], *cloud) != 0)
  18. {
  19. return -1;
  20. }
  21. pcl::ConditionalEuclideanClustering<pcl::PointXYZ> clustering;
  22. clustering.setClusterTolerance( 0.02);
  23. clustering.setMinClusterSize( 100);
  24. clustering.setMaxClusterSize( 25000);
  25. clustering.setInputCloud(cloud);
  26. //设置每次检测一对点云时的函数
  27. clustering.setConditionFunction(&customCondition);
  28. std::vector<pcl::PointIndices> clusters;
  29. clustering.segment(clusters);
  30. int currentClusterNum = 1;
  31. for (std::vector<pcl::PointIndices>::const_iterator i = clusters.begin(); i != clusters.end(); ++i)
  32. {
  33. pcl::PointCloud<pcl::PointXYZ>::Ptr cluster( new pcl::PointCloud<pcl::PointXYZ>);
  34. for (std::vector< int>::const_iterator point = i->indices.begin(); point != i->indices.end(); point++)
  35. cluster->points.push_back(cloud->points[*point]);
  36. cluster->width = cluster->points.size();
  37. cluster->height = 1;
  38. cluster->is_dense = true;
  39. if (cluster->points.size() <= 0)
  40. break;
  41. std::cout << "Cluster " << currentClusterNum << " has " << cluster->points.size() << " points." << std::endl;
  42. std:: string fileName = "cluster" + boost::to_string(currentClusterNum) + ".pcd";
  43. pcl::io::savePCDFileASCII(fileName, *cluster);
  44. currentClusterNum++;
  45. }
  46. }

1.3.3密度聚类

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

相关概念:

DBSCAN中为了增加抗噪声的能力,引入了核心对象等概念。

ε:参数,邻域距离。

minPts:参数,核心点领域内最少点数。

核心点:在 ϵ \epsilonϵ 邻域内有至少 m i n P t s minPtsminPt**s 个邻域点的点为核心点。

直接密度可达:对于样本集合 D DD,如果样本点 q qq 在 p pp 的 ϵ \epsilonϵ 邻域内,并且 p pp 为核心对象,那么对象 q qq 从对象 p pp 直接密度可达。

密度可达:对于样本集合 D DD,给定一串样本点p 1 p_1p1, p 2 p_2p2, …, p n p_np**n,p = p 1 p = p_1p=p1, q = p n q = p_nq=p**n, 假如对象 p i p_ip**i 从 p i − 1 p_{i-1}p**i−1 直接密度可达,那么对象 q qq 从对象 p pp 密度可达。

密度相连:存在样本集合 D DD 中的一点 o oo ,如果对象 o oo 到对象 p pp 和对象 q qq 都是密度可达的,那么 p pp 和 q qq 密度相联。

DBSCAN 算法核心是找到密度相连对象的最大集合。

b99a063a494ba9f4036df252d52faf4a.png

算法:

逐点遍历,如果该点非核心点,则认为是噪声点并忽视(噪声点可能在后续被核心点归入聚类中),若为核心点则新建聚类,并将所有邻域点加入聚类。对于邻域点中的核心点,还要递归地把其邻域点加入聚类。依此类推直到无点可加入该聚类,并开始考虑新的点,建立新的聚类。利用kdtree查询邻域点。

实现代码:


    
    
  1. #ifndef DBSCAN_H
  2. #define DBSCAN_H
  3. #include <pcl/point_types.h>
  4. #define UN_PROCESSED 0
  5. #define PROCESSING 1
  6. #define PROCESSED 2
  7. inline bool comparePointClusters ( const pcl::PointIndices &a, const pcl::PointIndices &b) {
  8. return (a.indices.size () < b.indices.size ());
  9. }
  10. template <typename PointT>
  11. class DBSCANSimpleCluster {
  12. public:
  13. typedef typename pcl::PointCloud<PointT>::Ptr PointCloudPtr;
  14. typedef typename pcl::search::KdTree<PointT>::Ptr KdTreePtr;
  15. virtual void setInputCloud(PointCloudPtr cloud) {
  16. input_cloud_ = cloud;
  17. }
  18. void setSearchMethod(KdTreePtr tree) {
  19. search_method_ = tree;
  20. }
  21. void extract(std::vector<pcl::PointIndices>& cluster_indices) {
  22. std::vector< int> nn_indices;
  23. std::vector<float> nn_distances;
  24. std::vector< bool> is_noise(input_cloud_->points.size(), false);
  25. std::vector< int> types(input_cloud_->points.size(), UN_PROCESSED);
  26. for ( int i = 0; i < input_cloud_->points.size(); i++) {
  27. if (types[i] == PROCESSED) {
  28. continue;
  29. }
  30. int nn_size = radiusSearch(i, eps_, nn_indices, nn_distances);
  31. if (nn_size < minPts_) {
  32. is_noise[i] = true;
  33. continue;
  34. }
  35. std::vector< int> seed_queue;
  36. seed_queue.push_back(i);
  37. types[i] = PROCESSED;
  38. for ( int j = 0; j < nn_size; j++) {
  39. if (nn_indices[j] != i) {
  40. seed_queue.push_back(nn_indices[j]);
  41. types[nn_indices[j]] = PROCESSING;
  42. }
  43. } // for every point near the chosen core point.
  44. int sq_idx = 1;
  45. while (sq_idx < seed_queue.size()) {
  46. int cloud_index = seed_queue[sq_idx];
  47. if (is_noise[cloud_index] || types[cloud_index] == PROCESSED) {
  48. // seed_queue.push_back(cloud_index);
  49. types[cloud_index] = PROCESSED;
  50. sq_idx++;
  51. continue; // no need to check neighbors.
  52. }
  53. nn_size = radiusSearch(cloud_index, eps_, nn_indices, nn_distances);
  54. if (nn_size >= minPts_) {
  55. for ( int j = 0; j < nn_size; j++) {
  56. if (types[nn_indices[j]] == UN_PROCESSED) {
  57. seed_queue.push_back(nn_indices[j]);
  58. types[nn_indices[j]] = PROCESSING;
  59. }
  60. }
  61. }
  62. types[cloud_index] = PROCESSED;
  63. sq_idx++;
  64. }
  65. if (seed_queue.size() >= min_pts_per_cluster_ && seed_queue.size () <= max_pts_per_cluster_) {
  66. pcl::PointIndices r;
  67. r.indices.resize(seed_queue.size());
  68. for ( int j = 0; j < seed_queue.size(); ++j) {
  69. r.indices[j] = seed_queue[j];
  70. }
  71. // These two lines should not be needed: (can anyone confirm?) -FF
  72. std::sort (r.indices.begin (), r.indices.end ());
  73. r.indices.erase (std::unique (r.indices.begin (), r.indices.end ()), r.indices.end ());
  74. r.header = input_cloud_->header;
  75. cluster_indices.push_back (r); // We could avoid a copy by working directly in the vector
  76. }
  77. } // for every point in input cloud
  78. std::sort (cluster_indices.rbegin (), cluster_indices.rend (), comparePointClusters);
  79. }
  80. void setClusterTolerance(double tolerance) {
  81. eps_ = tolerance;
  82. }
  83. void setMinClusterSize ( int min_cluster_size) {
  84. min_pts_per_cluster_ = min_cluster_size;
  85. }
  86. void setMaxClusterSize ( int max_cluster_size) {
  87. max_pts_per_cluster_ = max_cluster_size;
  88. }
  89. void setCorePointMinPts( int core_point_min_pts) {
  90. minPts_ = core_point_min_pts;
  91. }
  92. protected:
  93. PointCloudPtr input_cloud_;
  94. double eps_ { 0.0};
  95. int minPts_ { 1}; // not including the point itself.
  96. int min_pts_per_cluster_ { 1};
  97. int max_pts_per_cluster_ {std::numeric_limits< int>::max()};
  98. KdTreePtr search_method_;
  99. virtual int radiusSearch(
  100. int index, double radius, std::vector< int> &k_indices,
  101. std::vector<float> &k_sqr_distances) const
  102. {
  103. return this->search_method_->radiusSearch(index, radius, k_indices, k_sqr_distances);
  104. }
  105. }; // class DBSCANCluster
  106. #endif // DBSCAN_H

使用示例:


    
    
  1. pcl::search::KdTree<pcl::PointXYZ>::Ptr tree( new pcl::search::KdTree<pcl::PointXYZ>);
  2. tree->setInputCloud(point_cloud_input);
  3. std::vector<pcl::PointIndices> cluster_indices;
  4. DBSCANKdtreeCluster<pcl::PointXYZ> ec;
  5. ec.setCorePointMinPts( 20);
  6. ec.setClusterTolerance( 0.05);
  7. ec.setMinClusterSize( 100);
  8. ec.setMaxClusterSize( 25000);
  9. ec.setSearchMethod(tree);
  10. ec.setInputCloud(point_cloud_input);
  11. ec.extract(cluster_indices);

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

原理:

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

65e17ba3dd6c92ca39c8ab31a0e1c570.png

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

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

算法流程:

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

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

2.找到距离最小的两个块(A块中某点与B块中某点距离最小),并连接;

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

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

c99fe77139dd83e81433c8661059df49.png

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

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

a1a353ca75e03d38879f159c07636145.png

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


    
    
  1. //生成分割器
  2. pcl::MinCutSegmentation<pcl::PointXYZ> seg;
  3. //分割输入分割目标
  4. seg.setInputCloud(cloud);
  5. //指定打击目标(目标点)
  6. pcl::PointCloud<pcl::PointXYZ>::Ptr foreground_points( new pcl::PointCloud<pcl::PointXYZ> ());
  7. pcl::PointXYZ point;
  8. point.x = 68.97;
  9. point.y = -18.55;
  10. point.z = 0.57;
  11. foreground_points->points.push_back(point);
  12. seg.setForegroundPoints(foreground_points);
  13. //指定权函数sigma
  14. seg.setSigma( 0.25);
  15. //物体大概范围
  16. seg.setRadius( 3.0433856);
  17. //用多少生成图
  18. seg.setNumberOfNeighbours( 14);
  19. //和目标点相连点的权值(至少有14个)
  20. seg.setSourceWeight( 0.8);
  21. //分割结果
  22. std::vector <pcl::PointIndices> clusters;
  23. seg.extract(clusters);

1.6.2超体聚类

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

原理:

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

6dd4c2f4119a2629ae63f6832f8e828d.png

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


    
    
  1. #include <pcl/console/parse.h>
  2. #include <pcl/point_cloud.h>
  3. #include <pcl/point_types.h>
  4. #include <pcl/io/pcd_io.h>
  5. #include <pcl/visualization/pcl_visualizer.h>
  6. #include <pcl/segmentation/supervoxel_clustering.h>
  7. //VTK include needed for drawing graph lines
  8. #include <vtkPolyLine.h>
  9. // 数据类型
  10. typedef pcl::PointXYZRGBA PointT;
  11. typedef pcl::PointCloud<PointT> PointCloudT;
  12. typedef pcl::PointNormal PointNT;
  13. typedef pcl::PointCloud<PointNT> PointNCloudT;
  14. typedef pcl::PointXYZL PointLT;
  15. typedef pcl::PointCloud<PointLT> PointLCloudT;
  16. //可视化
  17. void addSupervoxelConnectionsToViewer (PointT &supervoxel_center,
  18. PointCloudT &adjacent_supervoxel_centers,
  19. std:: string supervoxel_name,
  20. boost::shared_ptr<pcl::visualization::PCLVisualizer> & viewer);
  21. int
  22. main ( int argc, char ** argv)
  23. {
  24. //解析命令行
  25. if (argc < 2)
  26. {
  27. pcl::console::print_error ( "Syntax is: %s <pcd-file> \n "
  28. "--NT Dsables the single cloud transform \n"
  29. "-v <voxel resolution>\n-s <seed resolution>\n"
  30. "-c <color weight> \n-z <spatial weight> \n"
  31. "-n <normal_weight>\n", argv[ 0]);
  32. return ( 1);
  33. }
  34. //打开点云
  35. PointCloudT::Ptr cloud = boost::shared_ptr <PointCloudT> ( new PointCloudT ());
  36. pcl::console::print_highlight ( "Loading point cloud...\n");
  37. if (pcl::io::loadPCDFile<PointT> (argv[ 1], *cloud))
  38. {
  39. pcl::console::print_error ( "Error loading cloud file!\n");
  40. return ( 1);
  41. }
  42. bool disable_transform = pcl::console::find_switch (argc, argv, "--NT");
  43. float voxel_resolution = 0.008f; //分辨率
  44. bool voxel_res_specified = pcl::console::find_switch (argc, argv, "-v");
  45. if (voxel_res_specified)
  46. pcl::console::parse (argc, argv, "-v", voxel_resolution);
  47. float seed_resolution = 0.1f;
  48. bool seed_res_specified = pcl::console::find_switch (argc, argv, "-s");
  49. if (seed_res_specified)
  50. pcl::console::parse (argc, argv, "-s", seed_resolution);
  51. float color_importance = 0.2f;
  52. if (pcl::console::find_switch (argc, argv, "-c"))
  53. pcl::console::parse (argc, argv, "-c", color_importance);
  54. float spatial_importance = 0.4f;
  55. if (pcl::console::find_switch (argc, argv, "-z"))
  56. pcl::console::parse (argc, argv, "-z", spatial_importance);
  57. float normal_importance = 1.0f;
  58. if (pcl::console::find_switch (argc, argv, "-n"))
  59. pcl::console::parse (argc, argv, "-n", normal_importance);
  60. //如何使用SupervoxelClustering函数
  61. pcl::SupervoxelClustering<PointT> super (voxel_resolution, seed_resolution);
  62. if (disable_transform) //如果设置的是参数--NT 就用默认的参数
  63. super.setUseSingleCameraTransform ( false);
  64. super.setInputCloud (cloud);
  65. super.setColorImportance (color_importance); //0.2f
  66. super.setSpatialImportance (spatial_importance); //0.4f
  67. super.setNormalImportance (normal_importance); //1.0f
  68. std:: map <uint32_t, pcl::Supervoxel<PointT>::Ptr > supervoxel_clusters;
  69. pcl::console::print_highlight ( "Extracting supervoxels!\n");
  70. super.extract (supervoxel_clusters);
  71. pcl::console::print_info ( "Found %d supervoxels\n", supervoxel_clusters.size ());
  72. boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer ( new pcl::visualization::PCLVisualizer ( "3D Viewer"));
  73. viewer->setBackgroundColor ( 0, 0, 0);
  74. PointCloudT::Ptr voxel_centroid_cloud = super.getVoxelCentroidCloud (); //获得体素中心的点云
  75. viewer->addPointCloud (voxel_centroid_cloud, "voxel centroids");
  76. viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 2.0, "voxel centroids"); //渲染点云
  77. viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_OPACITY, 0.95, "voxel centroids");
  78. PointLCloudT::Ptr labeled_voxel_cloud = super.getLabeledVoxelCloud ();
  79. viewer->addPointCloud (labeled_voxel_cloud, "labeled voxels");
  80. viewer->setPointCloudRenderingProperties (pcl::visualization::PCL_VISUALIZER_OPACITY, 0.8, "labeled voxels");
  81. PointNCloudT::Ptr sv_normal_cloud = super.makeSupervoxelNormalCloud (supervoxel_clusters);
  82. //We have this disabled so graph is easy to see, uncomment to see supervoxel normals
  83. //viewer->addPointCloudNormals<PointNormal> (sv_normal_cloud,1,0.05f, "supervoxel_normals");
  84. pcl::console::print_highlight ( "Getting supervoxel adjacency\n");
  85. std::multimap<uint32_t, uint32_t> supervoxel_adjacency;
  86. super.getSupervoxelAdjacency (supervoxel_adjacency);
  87. //To make a graph of the supervoxel adjacency, we need to iterate through the supervoxel adjacency multimap
  88. //为了使整个超体形成衣服图,我们需要遍历超体的每个临近的个体
  89. std::multimap<uint32_t,uint32_t>::iterator label_itr = supervoxel_adjacency.begin ();
  90. for ( ; label_itr != supervoxel_adjacency.end (); )
  91. {
  92. //First get the label
  93. uint32_t supervoxel_label = label_itr->first;
  94. //Now get the supervoxel corresponding to the label
  95. pcl::Supervoxel<PointT>::Ptr supervoxel = supervoxel_clusters.at (supervoxel_label);
  96. //Now we need to iterate through the adjacent supervoxels and make a point cloud of them
  97. PointCloudT adjacent_supervoxel_centers;
  98. std::multimap<uint32_t,uint32_t>::iterator adjacent_itr = supervoxel_adjacency.equal_range (supervoxel_label).first;
  99. for ( ; adjacent_itr!=supervoxel_adjacency.equal_range (supervoxel_label).second; ++adjacent_itr)
  100. {
  101. pcl::Supervoxel<PointT>::Ptr neighbor_supervoxel = supervoxel_clusters.at (adjacent_itr->second);
  102. adjacent_supervoxel_centers.push_back (neighbor_supervoxel->centroid_);
  103. }
  104. //Now we make a name for this polygon
  105. std::stringstream ss;
  106. ss << "supervoxel_" << supervoxel_label;
  107. //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
  108. //从给定的点云中生成一个星型的多边形,
  109. addSupervoxelConnectionsToViewer (supervoxel->centroid_, adjacent_supervoxel_centers, ss.str (), viewer);
  110. //Move iterator forward to next label
  111. label_itr = supervoxel_adjacency.upper_bound (supervoxel_label);
  112. }
  113. while (!viewer->wasStopped ())
  114. {
  115. viewer->spinOnce ( 100);
  116. }
  117. return ( 0);
  118. }
  119. //VTK可视化构成的聚类图
  120. void
  121. addSupervoxelConnectionsToViewer (PointT &supervoxel_center,
  122. PointCloudT &adjacent_supervoxel_centers,
  123. std:: string supervoxel_name,
  124. boost::shared_ptr<pcl::visualization::PCLVisualizer> & viewer)
  125. {
  126. vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New ();
  127. vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New ();
  128. vtkSmartPointer<vtkPolyLine> polyLine = vtkSmartPointer<vtkPolyLine>::New ();
  129. //Iterate through all adjacent points, and add a center point to adjacent point pair
  130. PointCloudT::iterator adjacent_itr = adjacent_supervoxel_centers.begin ();
  131. for ( ; adjacent_itr != adjacent_supervoxel_centers.end (); ++adjacent_itr)
  132. {
  133. points->InsertNextPoint (supervoxel_center.data);
  134. points->InsertNextPoint (adjacent_itr->data);
  135. }
  136. // Create a polydata to store everything in
  137. vtkSmartPointer<vtkPolyData> polyData = vtkSmartPointer<vtkPolyData>::New ();
  138. // Add the points to the dataset
  139. polyData->SetPoints (points);
  140. polyLine->GetPointIds ()->SetNumberOfIds(points->GetNumberOfPoints ());
  141. for(unsigned int i = 0; i < points->GetNumberOfPoints (); i++)
  142. polyLine->GetPointIds ()->SetId (i,i);
  143. cells->InsertNextCell (polyLine);
  144. // Add the lines to the dataset
  145. polyData->SetLines (cells);
  146. viewer->addModelFromPolyData (polyData,supervoxel_name);
  147. }
  • 1
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值