<think>我们
构建2.5D
Delaunay三角网，通常是指对三维空间中的点集在xy平面上的投影进行二维
Delaunay三角剖分，但每个顶点保留其原始z坐标（高程），从而形成具有地形起伏特征的
三角网。这种
三角网常用于地形建模。在
CGAL中，我们可以使用`
Delaunay_triangulation_2`类来实现。以下是详细步骤和代码示例：###步骤1.**定义内核和
Delaunay三角剖分类**：选择合适的内核（如精确谓词但不精确构造的内核）和三角剖分数据结构。2.**插入点集**：将三维点投影到二维（忽略z坐标），插入到三角剖分中。3.**遍历
三角网**：遍历所有面，获取每个面的顶点信息（包括原始z坐标）。4.**可视化**：可以使用
Open3D或其他工具可视化结果。###示例代码```cpp#include<
CGAL/Exact_predicates_inexact_constructions_kernel.h>#include<
CGAL/
Delaunay_triangulation_2.h>#include<vector>typedef
CGAL::Exact_predicates_inexact_constructions_kernelKernel;typedefKernel
::Point_2Point_2;//二维点类型（用于三角剖分）typedefKernel
::Point_3Point_3;//原始三维点类型（用于输入）typedef
CGAL::Delaunay_triangulation_2<Kernel>
Delaunay;intmain(){//示例点集（2.5D：x,y为平面坐标，z为高程）std
::vector<Point_3>points={Point_3(0.0,0.0,1.0),Point_3(1.0,0.0,2.0),Point_3(0.0,1.0,3.0),Point_3(1.0,1.0,4.0),Point_3(0.5,0.5,5.0)};//创建
Delaunay三角剖分对象
Delaunaydt;//将三维点投影到二维（忽略z坐标），并插入三角剖分for(constauto&p
:points){Point_2p2(p.x(),p.y());//投影到xy平面dt.insert(p2);}//输出三角剖分信息（顶点和面）std
::cout<<"顶点数量
:"<<dt.number_of_vertices()<<std
::endl;std
::cout<<"面数量
:"<<dt.number_of_faces()<<std
::endl;//遍历所有面，输出每个面的顶点（包括原始z坐标）for(
Delaunay::Finite_faces_iteratorfit=dt.finite_faces_begin();fit!=dt.finite_faces_end();++fit){
Delaunay::Face_handleface=fit;//获取面的三个顶点（在三角剖分中的二维点）Point_2p0=face->vertex(0)->point();Point_2p1=face->vertex(1)->point();Point_2p2=face->vertex(2)->point();//查找原始三维点（通过xy坐标匹配）//注意：实际应用中，如果点集很大，应建立映射关系以提高效率Point_3v0,v1,v2;for(constauto&p
:points){if(p.x()==p0.x()&&p.y()==p0.y())v0=p;if(p.x()==p1.x()&&p.y()==p1.y())v1=p;if(p.x()==p2.x()&&p.y()==p2.y())v2=p;}//输出三角形（三维点）std
::cout<<"三角形
:"<<"("<<v0<<")-"<<"("<<v1<<")-"<<"("<<v2<<")"<<std
::endl;}return0;}```###优化建议1.**点映射优化**：上述代码中通过循环查找三维点效率较低。实际应用中，可以在插入点之前建立从`Point_2`到`Point_3`的映射（例如`std
::map`）：```cppstd
::map<Point_2,Point_3>point_map;for(constauto&p
:points){Point_2p2(p.x(),p.y());point_map[p2]=p;//注意：如果有重复点（xy相同），需要处理}```然后在遍历面时，直接通过`Point_2`查找对应的`Point_3`。2.**可视化**：
CGAL本身不提供可视化功能，可以将结果导出为OFF格式，用MeshLab查看，或使用其他库（如
Open3D）进行可视化。以下是导出为OFF文件的示例：```cppstd
::ofstreamout("terrain.off");out<<"OFF\n";out<<dt.number_of_vertices()<<""<<dt.number_of_faces()<<"0\n";//首先输出所有顶点（三维点）for(
Delaunay::Vertex_iteratorvit=dt.vertices_begin();vit!=dt.vertices_end();++vit){Point_2p2=vit->point();Point_3p3=point_map[p2];//从映射中获取三维点out<<p3.x()<<""<<p3.y()<<""<<p3.z()<<"\n";}//输出所有面（三角形）for(
Delaunay::Finite_faces_iteratorfit=dt.finite_faces_begin();fit!=dt.finite_faces_end();++fit){//注意：OFF文件中面的索引从0开始out<<"3";for(inti=0;i<3;++i){//计算顶点索引（需要建立顶点到索引的映射）//这里省略了索引映射的建立，实际需要遍历顶点并分配索引}out<<"\n";}out.close();```###注意事项-**重复点处理**：在插入点之前，应确保没有重复的二维点（xy坐标相同），否则会导致三角剖分失败。可以使用`std
::set`或
CGAL的`Unique_hash_map`进行去重。-**性能**：插入$n$个点的时间复杂度为$O(n\logn)$，适合大规模点集。-**约束条件**：如果需要考虑地形特征线（如河流、山脊），需使用约束
Delaunay三角剖分（Constrained
DelaunayTriangulation,CDT）。###数学原理
Delaunay三角剖分满足最小角最大化原则，即任意三角形的外接圆内不包含其他点（空圆准则）。对于点集$P$，其
Delaunay三角剖分$DT(P)$满足：$$\forallt\inDT(P),\text{圆}(t)\capP=\emptyset$$其中$\text{圆}(t)$表示三角形$t$的外接圆[^1]。###相关扩展1.**带约束的三角剖分**：使用`
CGAL::Constrained_
Delaunay_triangulation_2`类。2.**地形插值**：在
三角网基础上，可以通过重心坐标等方法进行高程插值。§§相关问题§§1.如何处理带断裂线的地形三角剖分？2.
CGAL中如何进行三维
点云的
Delaunay四面体剖分？3.如何从
Delaunay三角网生成等高线？</think>### 使用
CGAL 构建 2.5D
Delaunay 三角网
2.5D
Delaunay 三角网（也称为 TIN - 不规则
三角网）是地形建模中的关键技术，它在 xy 平面进行
Delaunay 三角剖分，同时保留 z 坐标（高程）信息。以下是完整实现方法：
#### 核心原理
1. **投影处理**：在 xy 平面进行 2D
Delaunay 三角剖分
2. **高程保留**：使用原始
3D 点的高程数据
3. **
Delaunay 准则**：所有三角形的外接圆内不包含其他点[^1]
#### 代码实现
```cpp
#include <
CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <
CGAL/
Delaunay_triangulation_2.h>
#include <
CGAL/Projection_traits_xy_3.h>
#include <fstream>
#include <vector>
typedef
CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel
::Point_3 Point_3;
typedef
CGAL::Projection_traits_xy_3<Kernel> Projection_traits;
typedef
CGAL::Delaunay_triangulation_2<Projection_traits>
Delaunay;
typedef
Delaunay::Vertex_handle Vertex_handle;
typedef
Delaunay::Face_handle Face_handle;
int main() {
// 1. 读取
点云数据 (x,y,z)
std
::vector<Point_3> points;
points.push_back(Point_3(0, 0, 100)); // 点1
points.push_back(Point_3(10, 0, 105)); // 点2
points.push_back(Point_3(0, 10, 110)); // 点3
points.push_back(Point_3(10, 10, 115)); // 点4
points.push_back(Point_3(5, 5, 120)); // 点5
// 2. 创建
三角网对象
Delaunay dt;
// 3. 插入点集（自动进行
Delaunay三角剖分）
dt.insert(points.begin(), points.end());
// 4. 验证三角剖分质量
std
::cout << "顶点数
: " << dt.number_of_vertices() << std
::endl;
std
::cout << "三角形数
: " << dt.number_of_faces() << std
::endl;
// 5. 遍历并输出所有三角形
std
::cout << "\n
三角网结构
:" << std
::endl;
for (Face_handle fh
: dt.finite_face_handles()) {
Point_3 p0 = fh->vertex(0)->point();
Point_3 p1 = fh->vertex(1)->point();
Point_3 p2 = fh->vertex(2)->point();
std
::cout << "三角形
: ("
<< p0 << "), ("
<< p1 << "), ("
<< p2 << ")" << std
::endl;
}
// 6. 导出为OFF格式（可用MeshLab查看）
std
::ofstream offFile("terrain_tin.off");
offFile << "OFF\n";
offFile << dt.number_of_vertices() << " "
<< dt.number_of_faces() << " 0\n";
// 输出顶点
for (Vertex_handle vh
: dt.finite_vertex_handles()) {
offFile << vh->point().x() << " "
<< vh->point().y() << " "
<< vh->point().z() << "\n";
}
// 输出三角形面
for (Face_handle fh
: dt.finite_face_handles()) {
offFile << "3 ";
for (int i = 0; i < 3; ++i) {
offFile << dt.triangle(fh).vertex(i).index() << " ";
}
offFile << "\n";
}
offFile.close();
return 0;
}
```
#### 关键组件说明
1. **投影特征**：`Projection_traits_xy_3`
仅使用 xy 坐标进行三角剖分，但保留 z 坐标信息
2. **
Delaunay 准则验证**
自动确保每个三角形满足：
$$ \forall p \in \text{点集}, p \notin \text{circumcircle}(ABC) $$
其中 ABC 是任意三角形
3. **数据结构**
- `Vertex_handle`
: 访问顶点属性
- `Face_handle`
: 访问三角形面
#### 处理真实数据
```cpp
// 从文件读取
点云（XYZ格式）
std
::ifstream input("terrain.xyz");
std
::vector<Point_3> points;
double x, y, z;
while (input >> x >> y >> z) {
points.push_back(Point_3(x, y, z));
}
input.close();
// 添加边界约束点（可选）
points.push_back(Point_3(min_x, min_y, 0)); // 边界点1
points.push_back(Point_3(max_x, max_y, 0)); // 边界点2
```
#### 可视化建议
1. **MeshLab**：打开导出的 OFF 文件
2. **ParaView**：支持大规模地形渲染
3. **
CGAL 内置查看器**（需安装 Qt）：
```cpp
#include <
CGAL/draw_triangulation_2.h>
CGAL::draw(dt); // 交互式查看
```
#### 性能优化
1. **空间索引**：使用 `Spatial_sort` 加速插入
```cpp
#include <
CGAL/spatial_sort.h>
CGAL::spatial_sort(points.begin(), points.end());
dt.insert(points.begin(), points.end());
```
2. **并行处理**（
CGAL ≥ 5.0）：
```cpp
dt.insert(points.begin(), points.end(), true); // 启用并行
```
### 应用场景
1. 数字高程模型（DEM）生成
2. 地形等高线提取
3. 洪水淹没分析
4. 地质结
构建模
本文介绍了如何结合Open3D和CGAL库构建点云的Delaunay三角网，包括导入点云数据、使用CGAL的Delaunay三角化算法以及通过Open3D进行可视化。通过示例代码，读者可以学习到点云处理和Delaunay三角网构建的方法。
订阅专栏 解锁全文
扫一扫
8100
到【灌水乐园】发言
