CGAL-5.4.1三角剖分和点云分割简单案例

1、二维Delaunay三角剖分

#include <CGAL/CORE_Expr.h>               // 包含CORE库，用于精确的实数计算
#include <CGAL/Simple_cartesian.h>        // 包含CGAL的简单笛卡尔坐标系统
#include <CGAL/Delaunay_triangulation_2.h>// 包含CGAL的二维Delaunay三角剖分库

// 使用CORE库中的Expr，这是一个用于表达精确实数的类
typedef CORE::Expr Real;
// 使用CGAL的Simple_cartesian模板，设定坐标类型为Real
typedef CGAL::Simple_cartesian<Real> K;
// 定义Delaunay三角剖分的类型，使用上面定义的坐标系统
typedef CGAL::Delaunay_triangulation_2<K> DT;

// 定义二维点的类型
typedef K::Point_2 Point_2;

int main() {
    // 创建一个Delaunay三角剖分的实例
    DT dt;
    // 定义一个双精度变量，并赋值为2，用于后续计算
    double two = 2;
    // 创建三个点，分别位于(0,0)，(sqrt(2),1)，(0,1)
    Point_2 p(0, 0), q(std::sqrt(two), 1), r(0, 1);

    // 将这三个点插入到Delaunay三角剖分中
    dt.insert(p);
    dt.insert(q);
    dt.insert(r);

    // 输出当前Delaunay三角剖分的状态
    std::cout << dt << std::endl;

    return 0;
}

2、点云分割

左边是分割前，右边是分割后

  

 

example_classification.cpp 

这段代码实现了点云数据的分类。它首先读取点云数据，然后计算了一系列的特征，并对数据应用了分类。最后，它使用标签为每个点着色，并将结果保存为PLY格式的文件。 

// 包含必要的CGAL库头文件
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Classification.h>
#include <CGAL/bounding_box.h>
#include <CGAL/tags.h>
#include <CGAL/IO/read_points.h>
#include <CGAL/IO/write_ply_points.h>
#include <CGAL/Real_timer.h>

#include <cstdlib>
#include <fstream>
#include <iostream>
#include <string>

// 定义并行或可用的并发标签
typedef CGAL::Parallel_if_available_tag Concurrency_tag;

// 设置CGAL的内核和相关类型
typedef CGAL::Simple_cartesian<double> Kernel;
typedef Kernel::Point_3 Point;
typedef Kernel::Iso_cuboid_3 Iso_cuboid_3;
typedef std::vector<Point> Point_range;
typedef CGAL::Identity_property_map<Point> Pmap;

// 设置CGAL的分类模块的命名空间
namespace Classification = CGAL::Classification;

// 定义分类器和相关特征类型
typedef Classification::Sum_of_weighted_features_classifier Classifier;
typedef Classification::Planimetric_grid<Kernel, Point_range, Pmap> Planimetric_grid;
typedef Classification::Point_set_neighborhood<Kernel, Point_range, Pmap> Neighborhood;
typedef Classification::Local_eigen_analysis Local_eigen_analysis;

// 定义标签和特征处理类型
typedef Classification::Label_handle Label_handle;
typedef Classification::Feature_handle Feature_handle;
typedef Classification::Label_set Label_set;
typedef Classification::Feature_set Feature_set;

// 定义具体特征
typedef Classification::Feature::Distance_to_plane<Point_range, Pmap> Distance_to_plane;
typedef Classification::Feature::Elevation<Kernel, Point_range, Pmap> Elevation;
typedef Classification::Feature::Vertical_dispersion<Kernel, Point_range, Pmap> Dispersion;

// 主函数开始
int main(int argc, char** argv)
{
    // 读取输入文件路径，如果没有提供则使用CGAL的数据文件路径
    const std::string filename = (argc > 1) ? argv[1] : CGAL::data_file_path("meshes/b9.ply");

    // 读取点云数据
    std::vector<Point> pts;
    if (!(CGAL::IO::read_points(filename, std::back_inserter(pts),
        CGAL::parameters::use_binary_mode(false))))
    {
        std::cerr << "Error: cannot read " << filename << std::endl;
        return EXIT_FAILURE;
    }

    // 设置格网分辨率和邻域数量
    float grid_resolution = 0.34f;
    unsigned int number_of_neighbors = 6;

    // 计算点云的边界盒
    Iso_cuboid_3 bbox = CGAL::bounding_box(pts.begin(), pts.end());

    // 初始化网格和邻域分析工具
    Planimetric_grid grid(pts, Pmap(), bbox, grid_resolution);  // 创建一个平面网格，用于辅助分类特征的计算
    Neighborhood neighborhood(pts, Pmap());  // 创建一个邻域对象，用于查询点云中点的邻近点
    Local_eigen_analysis eigen = Local_eigen_analysis::create_from_point_set(pts, Pmap(), neighborhood.k_neighbor_query(number_of_neighbors));  // 进行局部特征分析，计算每个点的局部特征

    // 计算特征
    float radius_neighbors = 1.7f;  // 设置邻域查询的半径
    float radius_dtm = 15.0f;  // 设置地形模型(DTM)计算的半径
    Feature_set features;  // 创建一个特征集合
    features.begin_parallel_additions();  // 开始并行添加特征（在并行模式下有效）
    Feature_handle distance_to_plane = features.add<Distance_to_plane>(pts, Pmap(), eigen);  // 添加“距离到平面”的特征
    Feature_handle dispersion = features.add<Dispersion>(pts, Pmap(), grid, radius_neighbors);  // 添加“垂直分散”的特征
    Feature_handle elevation = features.add<Elevation>(pts, Pmap(), grid, radius_dtm);  // 添加“高程”的特征
    features.end_parallel_additions();  // 结束特征的并行添加

    // 初始化标签
    Label_set labels;  // 创建一个标签集
    Label_handle ground = labels.add("ground");  // 添加地面标签
    Label_handle vegetation = labels.add("vegetation", CGAL::IO::Color(0, 255, 0));  // 添加植被标签，绿色表示
    Label_handle roof = labels.add("roof", CGAL::IO::Color(255, 0, 0), 6);  // 添加屋顶标签，红色表示

    // 设置分类器的权重和效果
    Classifier classifier(labels, features);  // 创建分类器，指定标签和特征
    classifier.set_weight(distance_to_plane, 6.75e-2f);  // 为“距离到平面”特征设置权重
    classifier.set_weight(dispersion, 5.45e-1f);  // 为“垂直分散”特征设置权重
    classifier.set_weight(elevation, 1.47e1f);  // 为“高程”特征设置权重

    // 设置每个标签在不同特征下的影响
    classifier.set_effect(ground, distance_to_plane, Classifier::NEUTRAL);  // 对于地面标签，“距离到平面”特征的影响是中性的
    classifier.set_effect(ground, dispersion, Classifier::NEUTRAL);  // 对于地面标签，“垂直分散”特征的影响是中性的
    classifier.set_effect(ground, elevation, Classifier::PENALIZING);  // 对于地面标签，“高程”特征的影响是惩罚性的

    classifier.set_effect(vegetation, distance_to_plane, Classifier::FAVORING);  // 对于植被标签，“距离到平面”特征的影响是有利的
    classifier.set_effect(vegetation, dispersion, Classifier::FAVORING);  // 对于植被标签，“垂直分散”特征的影响是有利的
    classifier.set_effect(vegetation, elevation, Classifier::NEUTRAL);  // 对于植被标签，“高程”特征的影响是中性的

    classifier.set_effect(roof, distance_to_plane, Classifier::NEUTRAL);  // 对于屋顶标签，“距离到平面”特征的影响是中性的
    classifier.set_effect(roof, dispersion, Classifier::NEUTRAL);  // 对于屋顶标签，“垂直分散”特征的影响是中性的
    classifier.set_effect(roof, elevation, Classifier::FAVORING);  // 对于屋顶标签，“高程”特征的影响是有利的


    // 执行分类
    std::vector<int> label_indices(pts.size(), -1);  // 初始化一个标签索引数组，用于存储每个点的标签索引
    CGAL::Real_timer t;  // 实例化一个计时器
    t.start();  // 开始计时
    // 使用分类器对点云进行分类，结果存储在label_indices中
    Classification::classify<CGAL::Parallel_if_available_tag>(pts, labels, classifier, label_indices);
    t.stop();  // 停止计时
    t.reset();  // 重置计时器

    // 应用局部平滑
    t.start();  // 开始计时
    // 使用局部平滑算法对分类结果进行平滑处理，以提高分类的准确性和连贯性
    Classification::classify_with_local_smoothing<CGAL::Parallel_if_available_tag>(
        pts, Pmap(), labels, classifier,
        neighborhood.sphere_neighbor_query(radius_neighbors),
        label_indices);
    t.stop();  // 停止计时
    t.reset();  // 重置计时器

    // 应用图割算法进行分类
    t.start();  // 开始计时
    // 使用图割算法进一步优化分类结果，这通常可以提供更精确和更一致的分类
    Classification::classify_with_graphcut<CGAL::Parallel_if_available_tag>(
        pts, Pmap(), labels, classifier,
        neighborhood.k_neighbor_query(12),
        0.2f, 4, label_indices);
    t.stop();  // 停止计时

    // 保存输出结果到PLY文件
    std::vector<unsigned char> red, green, blue;  // 初始化颜色数组
    // 为每个点根据其分类标签设置颜色
    for (std::size_t i = 0; i < pts.size(); ++i)
    {
        Label_handle label = labels[std::size_t(label_indices[i])];  // 获取每个点的标签
        unsigned r = 0, g = 0, b = 0;  // 初始化颜色值
        // 根据不同的标签设置不同的颜色
        if (label == ground) { r = 245; g = 180; b = 0; }
        else if (label == vegetation) { r = 0; g = 255; b = 27; }
        else if (label == roof) { r = 255; g = 0; b = 170; }
        red.push_back(r);
        green.push_back(g);
        blue.push_back(b);
    }

    std::ofstream f("classification.ply");  // 打开文件流准备写入
    // 将分类结果写入PLY文件，包括点的位置和对应的颜色
    CGAL::IO::write_PLY_with_properties(
        f, CGAL::make_range(boost::counting_iterator<std::size_t>(0), boost::counting_iterator<std::size_t>(pts.size())),
        CGAL::make_ply_point_writer(CGAL::make_property_map(pts)),
        std::make_pair(CGAL::make_property_map(red), CGAL::PLY_property<unsigned char>("red")),
        std::make_pair(CGAL::make_property_map(green), CGAL::PLY_property<unsigned char>("green")),
        std::make_pair(CGAL::make_property_map(blue), CGAL::PLY_property<unsigned char>("blue")));

    return EXIT_SUCCESS;
}