曲面重建技术在逆向工程、数据可视化、机器视觉、虚拟现实、医疗技术等领域中得到了广泛的应用 。 例如，在汽车、航空等工业领域中，复杂外形产品的设计仍需要根据手工模型，采用逆向工程的手段建立产品的数字化模型，根据测量数据建立人体以及骨骼和器官的计算机模型，在医学、定制生产等方面都有重要意义 。
除了上述传统的行业，随着新兴的廉价 RGBD 获取设备在数字娱乐行业的病毒式扩展，使得更多人开始使用点云来处理对象并进行工程应用 。 根据重建曲面和数据点云之间的关系，可将曲面重建分为两大类：插值法和逼近法。前者得到的重建曲面完全通过原始数据点，而后者则是用分片线性曲面或其他形式的曲面来逼近原始数据点，从而使得得到的重建曲面是原始点集的一个逼近曲面。
关联算法：Search、KdTree、Octree
基于多项式重构的平滑和法线估计
在测量较小的数据时会产生一些误差，这些误差所造成的不规则数据如果直接拿来曲面重建的话，会使得重建的曲面不光滑或者有漏洞，可以采用对数据重采样来解决这样问题，通过对周围的数据点进行高阶多项式插值来重建表面缺少的部分。
上图左侧，在包含两个配准点云的数据集中看到了配准后的效果及表面法线估计。由于对齐错误，所产生的法线有噪声。在右侧，使用移动最小二乘法对表面法线估计进行平滑处理后，在同一数据集中看到了该效果。绘制每个点的曲率，作为重新采样前后特征值关系的度量，得到：
代码实现
#include <pcl/point_types.h>
#include <pcl/io/pcd_io.h>
#include <pcl/search/kdtree.h> //kd-tree搜索对象的类定义的头文件
#include <pcl/surface/mls.h> //最小二乘法平滑处理类定义头文件
#include <pcl/visualization/cloud_viewer.h>
#include <boost/thread/thread.hpp>
int main(int argc, char** argv)
{
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>());
// 从文件读取点云图
pcl::PCDReader reader;
reader.read("G:/vsdata/PCLlearn/PCDdata/bun0.pcd", *cloud);
// 创建 KD-Tree
pcl::search::KdTree<pcl::PointXYZ>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZ>);
// 输出具有PointNormal类型，以便存储MLS计算的法线
pcl::PointCloud<pcl::PointNormal>::Ptr mls_points(new pcl::PointCloud<pcl::PointNormal>);
// Init对象（第二种点类型用于法线，即使未使用）
pcl::MovingLeastSquares<pcl::PointXYZ, pcl::PointNormal> mls;
mls.setComputeNormals(true);
// 设置参数
mls.setInputCloud(cloud);
mls.setPolynomialOrder(2);
mls.setSearchMethod(tree);
mls.setSearchRadius(0.1);
// 重建
mls.process(*mls_points);
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer1(new pcl::visualization::PCLVisualizer("before"));
viewer1->addPointCloud<pcl::PointXYZ>(cloud, "before");
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer(new pcl::visualization::PCLVisualizer("smoothed"));
viewer->addPointCloud<pcl::PointNormal>(mls_points, "smoothed");
while (!viewer->wasStopped()) {
viewer->spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(1000));
}
// 保存
pcl::io::savePCDFile("bun0-mls.pcd", *mls_points);
}实现效果
在平面模型上提取凸（凹）多边形
该实例先从点云中提取平面模型，再通过该估计的平面模型系数从滤波后的点云投影一组点集形成点云，最后为投影后的点云计算其对应的二维凸多边形。
代码实现
#include <pcl/ModelCoefficients.h> //采样一致性模型相关类头文件
#include <pcl/io/pcd_io.h>
#include <pcl/point_types.h>
#include <pcl/sample_consensus/method_types.h>
#include <pcl/sample_consensus/model_types.h>
#include <pcl/filters/passthrough.h>
#include <pcl/filters/project_inliers.h> //滤波相关类头文件
#include <pcl/segmentation/sac_segmentation.h> //基于采样一致性分割类定义的头文件
#include <pcl/surface/concave_hull.h> //创建凹多边形类定义头文件
#include <pcl/visualization/pcl_visualizer.h>
#include <boost/thread/thread.hpp>
int main(int argc, char** argv)
{
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>), cloud_filtered(new pcl::PointCloud<pcl::PointXYZ>), cloud_projected(new pcl::PointCloud<pcl::PointXYZ>);
pcl::PCDReader reader;
reader.read("G:/vsdata/PCLlearn/PCDdata/table_scene_mug_stereo_textured.pcd", *cloud);
// 创建过滤器消除杂散的NaN
pcl::PassThrough<pcl::PointXYZ> pass;
pass.setInputCloud(cloud); //设置输入点云
pass.setFilterFieldName("z"); //设置分割字段为z坐标
pass.setFilterLimits(0, 1.1); //设置分割阀值为(0, 1.1)
pass.filter(*cloud_filtered);
std::cerr << "PointCloud after filtering has: "
<< cloud_filtered->points.size() << " data points." << std::endl;
pcl::ModelCoefficients::Ptr coefficients(new pcl::ModelCoefficients);
pcl::PointIndices::Ptr inliers(new pcl::PointIndices); //inliers存储分割后的点云
// 创建分割对象
pcl::SACSegmentation<pcl::PointXYZ> seg;
// 设置优化系数，该参数为可选参数
seg.setOptimizeCoefficients(true);
// 强制性的
seg.setModelType(pcl::SACMODEL_PLANE);
seg.setMethodType(pcl::SAC_RANSAC);
seg.setDistanceThreshold(0.01);
seg.setInputCloud(cloud_filtered);
seg.segment(*inliers, *coefficients);
std::cerr << "PointCloud after segmentation has: "
<< inliers->indices.size() << " inliers." << std::endl;
// 投影模型inliers
pcl::ProjectInliers<pcl::PointXYZ> proj;//点云投影滤波模型
proj.setModelType(pcl::SACMODEL_PLANE); //设置投影模型
proj.setIndices(inliers);
proj.setInputCloud(cloud_filtered);
proj.setModelCoefficients(coefficients); //将估计得到的平面coefficients参数设置为投影平面模型系数
proj.filter(*cloud_projected); //得到投影后的点云
std::cerr << "PointCloud after projection has: "
<< cloud_projected->points.size() << " data points." << std::endl;
// 存储提取多边形上的点云
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud_hull(new pcl::PointCloud<pcl::PointXYZ>);
pcl::ConcaveHull<pcl::PointXYZ> chull; //创建多边形提取对象
chull.setInputCloud(cloud_projected); //设置输入点云为提取后点云
chull.setAlpha(0.1);
chull.reconstruct(*cloud_hull); //创建提取创建凹多边形
std::cerr << "Concave hull has: " << cloud_hull->points.size()
<< " data points." << std::endl;
pcl::PCDWriter writer;
writer.write("table_scene_mug_stereo_textured_hull.pcd", *cloud_hull, false);
// 可视化
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer(new pcl::visualization::PCLVisualizer("before"));
viewer->addPointCloud<pcl::PointXYZ>(cloud, "before");
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer1(new pcl::visualization::PCLVisualizer("filtered"));
viewer1->addPointCloud<pcl::PointXYZ>(cloud_filtered, "filtered");
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer2(new pcl::visualization::PCLVisualizer("projected"));
viewer2->addPointCloud<pcl::PointXYZ>(cloud_projected, "projected");
boost::shared_ptr<pcl::visualization::PCLVisualizer> viewer3(new pcl::visualization::PCLVisualizer("hull"));
viewer3->addPointCloud<pcl::PointXYZ>(cloud_hull, "hull");
while (!viewer->wasStopped())
{
viewer->spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(1000));
}
return (0);
}结果：
原始点云
滤波后点云（直通滤波）
分割平面
多边形
无序点云的快速三角化
使用贪婪投影三角化算法对有向点云进行三角化，具体方法是：
（1）先将有向点云投影到某一局部二维坐标平面内
（2）在坐标平面内进行平面内的三角化
（3）根据平面内三位点的拓扑连接关系获得一个三角网格曲面模型
贪婪投影三角化算法原理：是处理一系列可以使网格“生长扩大”的点（边缘点）延伸这些点直到所有符合几何正确性和拓扑正确性的点都被连上，该算法可以用来处理来自一个或者多个扫描仪扫描到得到并且有多个连接处的散乱点云但是算法也是有很大的局限性，它更适用于采样点云来自表面连续光滑的曲面且点云的密度变化比较均匀的情况。
代码实现
#include <pcl/point_types.h>
#include <pcl/io/pcd_io.h>
#include <pcl/kdtree/kdtree_flann.h>
#include <pcl/features/normal_3d.h>
#include <pcl/surface/gp3.h> //贪婪投影三角化算法
#include <pcl/visualization/pcl_visualizer.h>
#include <boost/thread/thread.hpp>
int main(int argc, char** argv)
{
// 将一个XYZ点类型的PCD文件打开并存储到对象中
pcl::PointCloud<pcl::PointXYZ>::Ptr cloud(new pcl::PointCloud<pcl::PointXYZ>);
pcl::PCLPointCloud2 cloud_blob;
pcl::io::loadPCDFile("G:/vsdata/PCLlearn/PCDdata/bun0.pcd", cloud_blob);
pcl::fromPCLPointCloud2(cloud_blob, *cloud);
// 正态估计
pcl::NormalEstimation<pcl::PointXYZ, pcl::Normal> n; //法线估计对象
pcl::PointCloud<pcl::Normal>::Ptr normals(new pcl::PointCloud<pcl::Normal>); //存储估计的法线
pcl::search::KdTree<pcl::PointXYZ>::Ptr tree(new pcl::search::KdTree<pcl::PointXYZ>); //定义kd树指针
tree->setInputCloud(cloud); //用cloud构建tree对象
n.setInputCloud(cloud);
n.setSearchMethod(tree);
n.setKSearch(20);
n.compute(*normals); //估计法线存储到其中
//* 法线不应包含点法线+曲面曲率
// 连接XYZ字段和法线字段*
pcl::PointCloud<pcl::PointNormal>::Ptr cloud_with_normals(new pcl::PointCloud<pcl::PointNormal>);
pcl::concatenateFields(*cloud, *normals, *cloud_with_normals); //连接字段
//* cloud_with_normals = cloud + normals
// 定义搜索树对象
pcl::search::KdTree<pcl::PointNormal>::Ptr tree2(new pcl::search::KdTree<pcl::PointNormal>);
tree2->setInputCloud(cloud_with_normals); //点云构建搜索树
// 初始化对象
pcl::GreedyProjectionTriangulation<pcl::PointNormal> gp3; //定义三角化对象
pcl::PolygonMesh triangles; //存储最终三角化的网络模型
// 设置连接点之间的最大距离（最大边长）
gp3.setSearchRadius(0.025); //设置连接点之间的最大距离，（即是三角形最大边长）
// 设置各参数值
gp3.setMu(2.5); //设置被样本点搜索其近邻点的最远距离为2.5，为了使用点云密度的变化
gp3.setMaximumNearestNeighbors(100); //设置样本点可搜索的邻域个数
gp3.setMaximumSurfaceAngle(M_PI / 4); // 设置某点法线方向偏离样本点法线的最大角度45
gp3.setMinimumAngle(M_PI / 18); // 设置三角化后得到的三角形内角的最小的角度为10
gp3.setMaximumAngle(2 * M_PI / 3); // 设置三角化后得到的三角形内角的最大角度为120
gp3.setNormalConsistency(false); //设置该参数保证法线朝向一致
// 获取结果
gp3.setInputCloud(cloud_with_normals); //设置输入点云为有向点云
gp3.setSearchMethod(tree2); //设置搜索方式
gp3.reconstruct(triangles); //重建提取三角化
// 附加顶点信息
std::vector<int> parts = gp3.getPartIDs();
std::vector<int> states = gp3.getPointStates();
// 可视化
// 创建可视化对象
pcl::visualization::PCLVisualizer viewer("PolygonMesh Viewer");
// 添加网格到可视化窗口
viewer.addPolygonMesh(triangles, "mesh");
// 设置观察点和方向
viewer.setCameraPosition(0, 0, -2, 0, -1, 0);
while (!viewer.wasStopped())
{
viewer.spinOnce(100);
boost::this_thread::sleep(boost::posix_time::microseconds(1000));
}
return (0);
}结果
