1.原理
对于点云中的每个点云点P，使用一大一小两个半径计算该点的两个法向量，然后用两个单位法向量的差异计算Don特征。
如下图：
图1为大半径计算得到的法向量；图2为小半径计算得到的法向量；图3为一大一小两个法向量计算得出的Don特征。
步骤：
（1）建KD树，计算点云平均密度（即点云点之间的平均邻近距离）。
（2）用户依据平均密度设置两个尺度的Don半径，然后分别计算两个Don半径下的法向量。
（3）通过公式：
计算两个尺度下Don法向量的差异。
（4）得到对应的法线微分差异，得出Don结果。
2.使用场景
（1）寻找特征点（例如边缘点，角点等；因为法线微分分割之后，变化差异大的位置即边缘位置）
（2）分割点云，得到点云的各个部分【往往用于分割具有多个平面的点云，如建筑物，道路面等】（在剔除Don计算出的边缘点之后，对其他点进行聚类，然后得到点云的不同部分）
3.注意事项
代码中的pcl::PointNormal不要写为pcl::Normal（因为我就写了...后面滤波报错了）
4.关键函数
（1）生成点云法线（在代码最初，生成点云法线时，最好使用K邻近搜索）
setKSearch()（2）计算平均半径，用于辅助分割（以下代码可以搜索与当前点云点最邻近的点；注意，参数设置为2，因为索引=0的点，是当前点云点本身）
std::vector<int> pointIdxNKNSearch(2); // 保存最邻近点的索引（该方法最低设置为2，其中第一个值为当前点）
std::vector<float> pointNKNSquaredDistance(2); // 保存最邻近点的距离的平方（同上）
tree->nearestKSearch(cloud->points[i], 2, pointIdxNKNSearch, pointNKNSquaredDistance);
double dis = sqrt(pointNKNSquaredDistance[1]);（3）大小尺度半径（大小尺度半径通过这个函数设置，以基于这两个函数，设置两个不同的值）
setRadiusSearch()（4）Don法线微分分割（设置大小半径的参数）
pcl::DifferenceOfNormalsEstimation<pcl::PointXYZRGB, pcl::PointNormal, pcl::PointNormal> feature;
feature.setNormalScaleLarge(ne_sm); // 设置小尺寸半径
feature.setNormalScaleSmall(ne_lg); // 设置大尺寸半径（5）条件滤波相关（用于过滤边缘点【保留边缘点同理】）条件滤波
pcl::ConditionOr<pcl::PointNormal>::Ptr filtersOr_keypoint(new pcl::ConditionOr<pcl::PointNormal>());
filtersOr_keypoint->addComparison(pcl::FieldComparison<pcl::PointNormal>::ConstPtr(new pcl::FieldComparison<pcl::PointNormal>("curvature", pcl::ComparisonOps::GT, 0.3))); // 条件选择器，曲率值大于0.3
pcl::ConditionalRemoval<pcl::PointNormal> filters_keypoint;
filters_keypoint.setInputCloud(donCloud);
filters_keypoint.setCondition(filtersOr_keypoint); // 设置条件选择器（或条件）
pcl::PointCloud<pcl::PointNormal>::Ptr cloud_filtered_keypoint(new pcl::PointCloud<pcl::PointNormal>);
filters_keypoint.filter(*cloud_filtered_keypoint);（6）聚类相关（Don之后聚类可将去掉边缘的点云分成几部分；聚类得到的结果为聚类前点云点的索引）
// 对去掉角点的内部点按照欧式聚类进行聚类，实现对物体不同面的分割功能
pcl::search::KdTree<pcl::PointNormal>::Ptr segtree(new pcl::search::KdTree<pcl::PointNormal>);
segtree->setInputCloud(cloud_filtered_inside);
std::vector<pcl::PointIndices> cluster_indices; // 聚类得到的索引集群（每个vector是一个集群）
pcl::EuclideanClusterExtraction<pcl::PointNormal> seg; // 欧式聚类
seg.setClusterTolerance(meanDis*5); // 设置欧式聚类的阈值
seg.setMinClusterSize(50); // 最小聚类点个数
seg.setMaxClusterSize(100000); // 最大聚类点个数
seg.setSearchMethod(segtree);
seg.setInputCloud(cloud_filtered_inside);
seg.extract(cluster_indices); // 得到聚类结果5.代码
#include <pcl/point_types.h>
#include <pcl/io/pcd_io.h>
#include <pcl/search/kdtree.h>
#include <pcl/features/normal_3d_omp.h>
#include <pcl/filters/conditional_removal.h> //条件滤波
#include <pcl/segmentation/extract_clusters.h>
#include <pcl/segmentation/impl/extract_clusters.hpp>
#include <pcl/features/don.h>
#include <pcl/visualization/pcl_visualizer.h>
#include <boost/thread/thread.hpp>
#include <random>
int main()
{
/****************Don法线微分分割********************/
// 原始点云
pcl::PointCloud<pcl::PointXYZRGB>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZRGB>);
pcl::io::loadPCDFile("D:/code/csdn/data/two_tree.pcd", *cloud); // 加载原始点云数据
// 生成点云法线
pcl::NormalEstimationOMP<pcl::PointXYZRGB, pcl::PointNormal> normEst; // 基于邻域的法线估计器(OMP型，可多线程处理)
pcl::PointCloud<pcl::PointNormal>::Ptr normals(new pcl::PointCloud<pcl::PointNormal>); // 法向量保存
pcl::search::KdTree<pcl::PointXYZRGB>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZRGB>()); // KD树
normEst.setNumberOfThreads(4);
normEst.setInputCloud(cloud);
normEst.setSearchMethod(tree);
normEst.setKSearch(50); // 使用K领域搜索
normEst.setRadiusSearch(0.3); // 进行法向估计的领域半径（该半径范围内的点云拟合一个平面，生成点云法向量）
normEst.compute(*normals);
// 计算平均半径，用于辅助分割
double totalDis = 0;
for (int i = 0; i < cloud->size(); i++) // 搜索最邻近的两个点
{
std::vector<int> pointIdxNKNSearch(2); // 保存最邻近点的索引（该方法最低设置为2，其中第一个值为当前点）
std::vector<float> pointNKNSquaredDistance(2); // 保存最邻近点的距离的平方（同上）
tree->nearestKSearch(cloud->points[i], 2, pointIdxNKNSearch, pointNKNSquaredDistance);
totalDis += sqrt(pointNKNSquaredDistance[1]);
}
double meanDis = totalDis / cloud->size(); // 平均距离
double scale_sm = meanDis * 2.0; // 小尺度辅助半径
double scale_lg = meanDis * 8.0; // 大尺度辅助半径
// 小尺度don半径
pcl::NormalEstimationOMP<pcl::PointXYZRGB, pcl::PointNormal> normEst_sm; // 基于邻域的法线估计器(OMP型，可多线程处理)
pcl::search::KdTree<pcl::PointXYZRGB>::Ptr tree_sm(new pcl::search::KdTree<pcl::PointXYZRGB>());
pcl::PointCloud<pcl::PointNormal>::Ptr ne_sm(new pcl::PointCloud<pcl::PointNormal>);
normEst_sm.setNumberOfThreads(4);
normEst_sm.setInputCloud(cloud);
normEst_sm.setSearchMethod(tree_sm);
normEst_sm.setRadiusSearch(scale_sm);
normEst_sm.compute(*ne_sm);
// 大尺度don半径
pcl::NormalEstimationOMP<pcl::PointXYZRGB, pcl::PointNormal> normEst_lg; // 基于邻域的法线估计器(OMP型，可多线程处理)
pcl::search::KdTree<pcl::PointXYZRGB>::Ptr tree_lg(new pcl::search::KdTree<pcl::PointXYZRGB>());
pcl::PointCloud<pcl::PointNormal>::Ptr ne_lg(new pcl::PointCloud<pcl::PointNormal>);
normEst_lg.setNumberOfThreads(4);
normEst_lg.setInputCloud(cloud);
normEst_lg.setSearchMethod(tree_lg);
normEst_lg.setRadiusSearch(scale_lg);
normEst_lg.compute(*ne_lg);
// Don法线微分分割
pcl::DifferenceOfNormalsEstimation<pcl::PointXYZRGB, pcl::PointNormal, pcl::PointNormal> feature;
feature.setInputCloud(cloud);
feature.setNormalScaleLarge(ne_sm); // 设置小尺寸半径
feature.setNormalScaleSmall(ne_lg); // 设置大尺寸半径
pcl::PointCloud<pcl::PointNormal>::Ptr donCloud(new pcl::PointCloud<pcl::PointNormal>);
pcl::copyPointCloud<pcl::PointXYZRGB, pcl::PointNormal>(*cloud, *donCloud); // 准备结果点云
feature.computeFeature(*donCloud);
// 应用一（寻找特征点）
// 按照曲率进行条件滤波，删除小于给定曲率阈值的点云，得到特征点(曲率变化大的点，即角点)
pcl::ConditionOr<pcl::PointNormal>::Ptr filtersOr_keypoint(new pcl::ConditionOr<pcl::PointNormal>());
filtersOr_keypoint->addComparison(pcl::FieldComparison<pcl::PointNormal>::ConstPtr(new pcl::FieldComparison<pcl::PointNormal>("curvature", pcl::ComparisonOps::GT, 0.3))); // 条件选择器，曲率值大于0.3
pcl::ConditionalRemoval<pcl::PointNormal> filters_keypoint;
filters_keypoint.setInputCloud(donCloud);
filters_keypoint.setCondition(filtersOr_keypoint); // 设置条件选择器（或条件）
pcl::PointCloud<pcl::PointNormal>::Ptr cloud_filtered_keypoint(new pcl::PointCloud<pcl::PointNormal>);
filters_keypoint.filter(*cloud_filtered_keypoint);
// 应用二（分割点云）
// 按照曲率进行条件滤波，删除大于给定曲率阈值的点云，得到内部点
pcl::ConditionOr<pcl::PointNormal>::Ptr filtersOr_inside(new pcl::ConditionOr<pcl::PointNormal>());
filtersOr_inside->addComparison(pcl::FieldComparison<pcl::PointNormal>::ConstPtr(new pcl::FieldComparison<pcl::PointNormal>("curvature", pcl::ComparisonOps::LT, 0.3))); // 条件选择器，曲率值小于0.3
pcl::ConditionalRemoval<pcl::PointNormal> filters_inside;
filters_inside.setInputCloud(donCloud);
filters_inside.setCondition(filtersOr_inside); // 设置条件选择器（或条件）
pcl::PointCloud<pcl::PointNormal>::Ptr cloud_filtered_inside(new pcl::PointCloud<pcl::PointNormal>);
filters_inside.filter(*cloud_filtered_inside);
// 对去掉角点的内部点按照欧式聚类进行聚类，实现对物体不同面的分割功能
pcl::search::KdTree<pcl::PointNormal>::Ptr segtree(new pcl::search::KdTree<pcl::PointNormal>);
segtree->setInputCloud(cloud_filtered_inside);
std::vector<pcl::PointIndices> cluster_indices; // 聚类得到的索引集群（每个vector是一个集群）
pcl::EuclideanClusterExtraction<pcl::PointNormal> seg; // 欧式聚类
seg.setClusterTolerance(meanDis*5); // 设置欧式聚类的阈值
seg.setMinClusterSize(50); // 最小聚类点个数
seg.setMaxClusterSize(100000); // 最大聚类点个数
seg.setSearchMethod(segtree);
seg.setInputCloud(cloud_filtered_inside);
seg.extract(cluster_indices); // 得到聚类结果
/****************展示********************/
// 整体曲率展示
boost::shared_ptr<pcl::visualization::PCLVisualizer> pointView(new pcl::visualization::PCLVisualizer("Curvature"));
pcl::visualization::PointCloudColorHandlerGenericField<pcl::PointNormal> handler(donCloud, "curvature");
pointView->setBackgroundColor(1, 1, 1);
pointView->addPointCloud(donCloud, handler);
pointView->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 3);
pointView->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_OPACITY, 0.5);
/*
while (!pointView->wasStopped())
{
pointView->spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(100000));
}
*/
// 内部点展示
boost::shared_ptr<pcl::visualization::PCLVisualizer> insideView(new pcl::visualization::PCLVisualizer("KeyPoint"));
pcl::visualization::PointCloudColorHandlerGenericField<pcl::PointNormal> handler_inside(cloud_filtered_inside, "curvature");
insideView->setBackgroundColor(1, 1, 1);
insideView->addPointCloud(cloud_filtered_inside, handler_inside);
insideView->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 3);
insideView->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_OPACITY, 0.5);
/*
while (!insideView->wasStopped())
{
insideView->spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(100000));
}
*/
// 聚类展示
pcl::PointCloud<pcl::PointXYZRGB>::Ptr clusters_cloud(new pcl::PointCloud<pcl::PointXYZRGB>);
clusters_cloud->resize(cloud_filtered_inside->points.size());
for (size_t i = 0; i < cloud_filtered_inside->points.size(); i++) // 设置初始值，默认为白色
{
clusters_cloud->points[i].x = cloud_filtered_inside->points[i].x;
clusters_cloud->points[i].y = cloud_filtered_inside->points[i].y;
clusters_cloud->points[i].z = cloud_filtered_inside->points[i].z;
clusters_cloud->points[i].r = 255;
clusters_cloud->points[i].g = 255;
clusters_cloud->points[i].b = 255;
}
for (size_t i = 0; i < cluster_indices.size(); i++)
{
std::random_device rand; // 生成随机颜色，每个聚类颜色一致
std::mt19937 gen(rand());
std::uniform_int_distribution<>dis(1, 254);
int r = dis(gen);
int g = dis(gen);
int b = dis(gen);
for (size_t j = 0; j < cluster_indices[i].indices.size(); j++)
{
clusters_cloud->points[cluster_indices[i].indices[j]].x = cloud_filtered_inside->points[cluster_indices[i].indices[j]].x;
clusters_cloud->points[cluster_indices[i].indices[j]].y = cloud_filtered_inside->points[cluster_indices[i].indices[j]].y;
clusters_cloud->points[cluster_indices[i].indices[j]].z = cloud_filtered_inside->points[cluster_indices[i].indices[j]].z;
clusters_cloud->points[cluster_indices[i].indices[j]].r = r;
clusters_cloud->points[cluster_indices[i].indices[j]].g = g;
clusters_cloud->points[cluster_indices[i].indices[j]].b = b;
}
}
boost::shared_ptr<pcl::visualization::PCLVisualizer> view_raw(new pcl::visualization::PCLVisualizer("cluster"));
view_raw->addPointCloud<pcl::PointXYZRGB>(clusters_cloud, "cluster cloud");
view_raw->setPointCloudRenderingProperties(pcl::visualization::PCL_VISUALIZER_POINT_SIZE, 1, "cluster cloud");
while (!view_raw->wasStopped())
{
view_raw->spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(100000));
}
return 0;
}
6.结果展示
（1）整体曲率展示
（2）去除边缘点之后的内部点展示
（3）聚类之后的内部点展示（将内部点分为几部分，用于分割平面，每个颜色是一个平面）
