13.2 泊松曲面重建

泊松重建是曲面重建系列方法中最常用的算法之一，它通过求解泊松方程从有向点云重建水密表面。

0. 算法动机：为什么需要泊松重建？

现实问题引入

想象你是一位文物保护专家，收到了一尊古希腊雕像的激光扫描数据：

数百万个点，带有位置信息和法向方向
但这些都是离散的点，没有表面
你需要生成一个封闭的水密网格用于3D打印

或者你是一位医学影像工程师：

从CT扫描获得了骨骼的点云数据
点云有噪声，有缺失区域
需要重建出光滑的骨骼表面用于手术规划

传统方法的局限

基于Delaunay的方法：

对噪声敏感，会产生不规则的三角形
难以处理非均匀采样
可能产生自相交

基于局部拟合的方法：

在复杂拓扑处容易失败
难以保证水密性
对噪声鲁棒性差

推进前沿法：

需要高质量的点云
难以处理大孔洞
可能产生自相交

泊松重建的优势

全局优化：不是局部决策，而是全局求解 噪声鲁棒：通过求解泊松方程自然平滑噪声 水密保证：输出总是封闭的流形 理论保证：有坚实的数学基础

如果不使用泊松重建

表面不封闭：出现孔洞，无法3D打印
噪声放大：高频噪声产生表面褶皱
拓扑错误：把手连接到主体，或产生不连通组件
细节丢失：过度平滑导致重要特征消失

1. 直观理解：3D扫描的轮廓

1.1 “指示函数就像剪影”

想象你正在拍摄一个人的剪影：

2D剪影：

人物内部：黑色（1）
背景：白色（0）
轮廓：黑白交界处

3D指示函数：

物体内部：1
外部：0
表面：0.5等值面（轮廓）

2D剪影：                      3D指示函数截面：
  000000000000                 000000000000
  000111110000                 000111110000
  001111111000                 001111111000
  001111111000                 001111111000
  000111110000                 000111110000
  000000000000                 000000000000
  
  轮廓：黑/白边界              表面：函数值=0.5的位置

1.2 “梯度就像坡度指示”

指示函数的梯度指向变化最剧烈的方向：

在表面处，梯度垂直于表面
梯度大小反映了从内部到外部的变化率
理想情况下，梯度等于表面法向

类比：想象一座山

指示函数：海拔高度（内部=山下，外部=山上）
梯度：坡度方向（垂直于等高线）
表面：海平面（高度=0）

1.3 “泊松方程就像平滑填充”

我们有：

输入：稀疏的法向采样（只在点云位置知道法向）
目标：找到一个平滑的指示函数，其梯度匹配输入法向

类比：拼图游戏

法向是边缘碎片，告诉你边界在哪里
泊松方程帮你填充内部，形成完整的图像
结果是平滑过渡，不是突兀变化

2. 数学基础

2.1 指示函数与梯度

泊松重建的核心思想是构造一个指示函数（indicator function） $χ : R^{3} \to {0, 1}$ ，该函数在物体内部取值为 1，外部为 0：

$χ (x) = {10 if x \in Ω otherwise$

其中 $Ω \subset R^{3}$ 是待重建的实体区域。

指示函数的梯度在表面边界处形成一个向量场： $\nabla χ (x) = {n (x) \cdot δ (x) 0 on \partial Ω otherwise$

其中 $n (x)$ 是表面法向量， $δ (x)$ 是 Dirac delta 函数。

2.2 泊松方程

由于直接计算指示函数较困难，算法转而求解一个平滑的近似函数 $f$ ，使其梯度最佳匹配输入法向量场 $V$ ：

$Δ f = \nabla \cdot V$

这是泊松方程的标准形式，其中：

$Δ f = \nabla^{2} f$ 是 Laplacian 算子
$\nabla \cdot V$ 是向量场的散度

2.3 变分形式

泊松方程等价于最小化以下能量泛函：

$E (f) = \int_{R^{3}} ∥\nabla f - V ∥^{2} d V$

通过变分法，最小化该能量得到欧拉-拉格朗日方程正是泊松方程。

3. 算法流程详解

3.1 八叉树结构

泊松重建使用八叉树（Octree）来组织空间：

八叉树层次结构：

Level 0:          Level 1:          Level 2:
+--------+        +----+----+       +--+--+--+--+ 
|        |        |    |    |       |##|  |  |  |
|        |   →    +----+----+   →   +--+--+--+--+ 
|        |        |    |    |       |  |##|##|  |
+--------+        +----+----+       +--+--+--+--+ 
                                     |  |##|  |  |
                                     +--+--+--+--+ 
                                     |  |  |  |##|
                                     +--+--+--+--+ 

## 表示包含点云的区域（需要细分）

八叉树的作用：

自适应细分：在点密集处细分更深
空间索引：快速定位相邻点
多分辨率：支持不同精度的重建

3.2 从点云到表面的重建过程

步骤可视化：

1. 输入点云：
   •     •     •
      •     •
   •     •     •
   （带有法向的离散点）

2. 构建八叉树：
   +---+---+---+---+
   |   |   |   |   |
   +---+---+---+---+
   |   |•••|•••|   |
   +---+---+---+---+
   |   |•••|•••|   |
   +---+---+---+---+
   |   |   |   |   |
   +---+---+---+---+
   （自适应细分）

3. 定义向量场：
   ↖  ↑  ↗
   ←  •  →    （每个节点有插值法向）
   ↙  ↓  ↘

4. 计算散度：
   [ 0.2] [ 0.5] [-0.1]
   [ 0.3] [ 1.0] [-0.2]   （标量场）
   [-0.1] [-0.3] [ 0.0]

5. 求解泊松方程：
   [0.1] [0.3] [0.5]
   [0.2] [0.6] [0.8]    （指示函数）
   [0.4] [0.7] [0.9]

6. 提取等值面：
       __________
      /          \
     /            \
    |              |   （三角网格）
     \\            /
      \\__________/

3.3 向量场构造

从输入点云 ${(p_{i}, n_{i})}$ 构造一个平滑的向量场：

$V (x) = \sum_{i} w_{i} (x) \cdot n_{i}$

其中 $w_{i} (x)$ 是基于八叉树的插值权重。

3.4 散度计算

计算向量场的散度：

$\nabla \cdot V (x) = \sum_{i} \nabla w_{i} (x) \cdot n_{i}$

3.5 泊松方程求解

在离散化后的空间中求解线性系统。CGAL 使用 Delaunay 三角剖分代替原始论文中的八叉树：

$L \cdot f = b$

其中：

$L$ 是离散 Laplacian 矩阵
$f$ 是待求的隐函数值
$b$ 是散度向量

3.6 等值面提取

使用 Marching Cubes 或 Delaunay 细化提取 $f (x) = iso$ 的等值面。等值通常选择为输入点处函数值的中位数。

4. CGAL 架构分析

4.1 核心类层次

Poisson_reconstruction_function<Gt>
    |
    +-- Reconstruction_triangulation_3<BaseGt, Gt, Tds>
            |
            +-- Delaunay_triangulation_3<Gt, Tds>
            |
            +-- Reconstruction_vertex_base_3<Gt, Vb>
            |       +-- Triangulation_vertex_base_3<Gt>
            |
            +-- Delaunay_triangulation_cell_base_3<Gt, Cb>
                    +-- Triangulation_cell_base_with_info_3<int, Gt>

4.2 关键文件位置

文件	路径	功能
`Poisson_reconstruction_function.h`	`Poisson_surface_reconstruction_3/include/CGAL/`	核心重建类
`Reconstruction_triangulation_3.h`	`Poisson_surface_reconstruction_3/include/CGAL/`	重建专用三角剖分
`Poisson_mesh_domain_3.h`	`Poisson_surface_reconstruction_3/include/CGAL/`	Mesh_3 域适配器
`poisson_surface_reconstruction.h`	`Poisson_surface_reconstruction_3/include/CGAL/`	便捷函数接口

4.3 Poisson_reconstruction_function 类

template <class Gt>
class Poisson_reconstruction_function
{
public:
    // 几何类型定义
    typedef Gt Geom_traits;
    typedef typename Geom_traits::FT FT;
    typedef typename Geom_traits::Point_3 Point;
    typedef typename Geom_traits::Vector_3 Vector;
    
    // 三角剖分类型
    typedef Reconstruction_triangulation_3<...> Triangulation;
    typedef typename Triangulation::Cell_handle Cell_handle;
    typedef typename Triangulation::Vertex_handle Vertex_handle;
    
    // 构造函数：从点云创建
    template <typename InputIterator, typename PointPMap, typename NormalPMap>
    Poisson_reconstruction_function(
        InputIterator first, InputIterator beyond,
        PointPMap point_pmap, NormalPMap normal_pmap);
    
    // 核心操作：计算隐函数
    template <class SparseLinearAlgebraTraits_d>
    bool compute_implicit_function(SparseLinearAlgebraTraits_d solver, 
                                   bool smoother_hole_filling = false);
    
    // 隐函数求值（实现 ImplicitFunction 概念）
    FT operator()(const Point& p) const;
    
    // 获取内部点和包围球
    Point get_inner_point() const;
    Sphere bounding_sphere() const;
};

5. 实现细节

5.1 Delaunay 细化

CGAL 的实现使用自适应 Delaunay 细化来离散化空间，而非固定八叉树：

// Delaunay 细化参数
const FT radius_edge_ratio_bound = 2.5;  // 半径-边比上界
const unsigned int max_vertices = 10000000;  // 最大顶点数
const FT enlarge_ratio = 1.5;  // 包围盒扩展比例
 
// 执行细化
unsigned int delaunay_refinement(FT radius_edge_ratio_bound,
                                 Sizing_field sizing_field,
                                 unsigned int max_vertices,
                                 FT enlarge_ratio);

细化过程迭代地插入 Steiner 点，直到满足以下条件：

所有四面体的半径-边比小于给定阈值
四面体尺寸符合尺寸场要求
顶点数不超过上限

5.2 线性系统组装

泊松方程离散化为稀疏线性系统：

template <class SparseLinearAlgebraTraits_d>
bool solve_poisson(SparseLinearAlgebraTraits_d solver, double lambda)
{
    // 初始化单元索引
    initialize_cell_indices();
    initialize_barycenters();
    
    // 约束一个顶点以固定常数偏移
    constrain_one_vertex_on_convex_hull();
    m_tr->index_unconstrained_vertices();
    
    unsigned int nb_variables = m_tr->number_of_vertices() - 1;
    
    // 组装矩阵 A 和右端项 B
    typename SparseLinearAlgebraTraits_d::Matrix A(nb_variables);
    typename SparseLinearAlgebraTraits_d::Vector X(nb_variables), B(nb_variables);
    
    // 对每个顶点组装行
    for (auto v = m_tr->finite_vertices_begin(); v \!= m_tr->finite_vertices_end(); ++v)
    {
        if (\!m_tr->is_constrained(v)) {
            B[v->index()] = div(v);  // 计算散度
            assemble_poisson_row<SparseLinearAlgebraTraits_d>(A, v, B, lambda);
        }
    }
    
    // 求解系统
    double D;
    if (\!solver.linear_solver(A, B, X, D))
        return false;
    
    // 将解复制回顶点
    // ...
    return true;
}

5.3 散度计算

散度计算有两种模式：

非归一化模式（默认）：

FT div(Vertex_handle v)
{
    FT div = 0.0;
    for (auto cell : incident_cells(v))
    {
        if (m_tr->is_infinite(cell)) continue;
        
        const int index = cell->index(v);
        const Point& a = cell->vertex(m_tr->vertex_triple_index(index, 0))->point();
        const Point& b = cell->vertex(m_tr->vertex_triple_index(index, 1))->point();
        const Point& c = cell->vertex(m_tr->vertex_triple_index(index, 2))->point();
        const Vector nn = CGAL::cross_product(b-a, c-a);
        
        // 累加相邻顶点的法向
        div += nn * (cell->vertex((index+1)%4)->normal() +
                     cell->vertex((index+2)%4)->normal() +
                     cell->vertex((index+3)%4)->normal());
    }
    return div;
}

归一化模式：考虑 Voronoi 面积权重，更精确但计算量更大。

5.4 Laplacian 矩阵组装

使用余切公式（cotangent formula）组装离散 Laplacian：

template <class SparseLinearAlgebraTraits_d>
void assemble_poisson_row(typename SparseLinearAlgebraTraits_d::Matrix& A,
                          Vertex_handle vi,
                          typename SparseLinearAlgebraTraits_d::Vector& B,
                          double lambda)
{
    double diagonal = 0.0;
    
    // 遍历所有邻边
    for (auto edge : incident_edges(vi))
    {
        Vertex_handle vj = edge.first->vertex(edge.third);
        if (vj == vi) vj = edge.first->vertex(edge.second);
        if (m_tr->is_infinite(vj)) continue;
        
        // 计算余切权重
        double cij = cotan_geometric(edge);
        
        if (m_tr->is_constrained(vj)) {
            B[vi->index()] -= cij * vj->f();  // 移到右端项
        } else {
            A.set_coef(vi->index(), vj->index(), -cij, true);
        }
        diagonal += cij;
    }
    
    // 设置对角元
    if (vi->type() == Triangulation::INPUT)
        A.set_coef(vi->index(), vi->index(), diagonal + lambda, true);
    else
        A.set_coef(vi->index(), vi->index(), diagonal, true);
}

5.5 等值面提取

CGAL 使用 Mesh_3 包进行等值面提取：

// 创建网格域
Mesh_domain domain = Mesh_domain::create_Poisson_mesh_domain(
    function, 
    Sphere(inner_point, sm_sphere_radius),
    CGAL::parameters::relative_error_bound(sm_dichotomy_error / sm_sphere_radius)
);
 
// 定义网格标准
Mesh_criteria criteria(
    CGAL::parameters::facet_angle = sm_angle,
    CGAL::parameters::facet_size = sm_radius * spacing,
    CGAL::parameters::facet_distance = sm_distance * spacing
);
 
// 生成网格
C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(
    domain, criteria,
    CGAL::parameters::surface_only()
                     .manifold_with_boundary()
);

6. 使用示例

6.1 完整示例代码

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Polyhedron_3.h>
#include <CGAL/Poisson_reconstruction_function.h>
#include <CGAL/Mesh_triangulation_3.h>
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
#include <CGAL/Mesh_criteria_3.h>
#include <CGAL/Poisson_mesh_domain_3.h>
#include <CGAL/make_mesh_3.h>
#include <CGAL/facets_in_complex_3_to_triangle_mesh.h>
#include <CGAL/property_map.h>
#include <CGAL/IO/read_points.h>
#include <CGAL/compute_average_spacing.h>
 
#include <vector>
#include <fstream>
 
// 类型定义
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel::FT FT;
typedef Kernel::Point_3 Point;
typedef Kernel::Vector_3 Vector;
typedef std::pair<Point, Vector> Point_with_normal;
typedef CGAL::First_of_pair_property_map<Point_with_normal> Point_map;
typedef CGAL::Second_of_pair_property_map<Point_with_normal> Normal_map;
typedef CGAL::Polyhedron_3<Kernel> Polyhedron;
 
// 泊松重建相关类型
typedef CGAL::Poisson_reconstruction_function<Kernel> Poisson_function;
typedef CGAL::Poisson_mesh_domain_3<Kernel> Mesh_domain;
typedef CGAL::Mesh_triangulation_3<Mesh_domain>::type Tr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3t3;
typedef CGAL::Mesh_criteria_3<Tr> Mesh_criteria;
 
int main(int argc, char* argv[])
{
    // 1. 读取点云数据
    std::vector<Point_with_normal> points;
    const char* input_file = (argc > 1) ? argv[1] : "input.xyz";
    
    if (\!CGAL::IO::read_points(input_file, std::back_inserter(points),
                               CGAL::parameters::point_map(Point_map())
                                                .normal_map(Normal_map())))
    {
        std::cerr << "Error: cannot read file " << input_file << std::endl;
        return EXIT_FAILURE;
    }
    
    std::cout << "Loaded " << points.size() << " points" << std::endl;
    
    // 2. 创建泊松重建函数
    Poisson_function function(points.begin(), points.end(),
                              Point_map(), Normal_map());
    
    // 3. 计算隐函数（使用默认的 Eigen 求解器）
    if (\!function.compute_implicit_function())
    {
        std::cerr << "Error: cannot solve Poisson equation" << std::endl;
        return EXIT_FAILURE;
    }
    
    // 4. 计算平均点间距
    FT average_spacing = CGAL::compute_average_spacing<CGAL::Sequential_tag>(
        points, 6, CGAL::parameters::point_map(Point_map()));
    
    // 5. 设置网格参数
    FT sm_angle = 20.0;      // 最小三角形角度（度）
    FT sm_radius = 30.0;     // 最大三角形大小（相对于平均间距）
    FT sm_distance = 0.375;  // 表面逼近误差（相对于平均间距）
    
    Sphere bsphere = function.bounding_sphere();
    FT radius = CGAL::approximate_sqrt(bsphere.squared_radius());
    
    // 6. 创建网格域并生成网格
    Mesh_domain domain = Mesh_domain::create_Poisson_mesh_domain(
        function, bsphere,
        CGAL::parameters::relative_error_bound(sm_distance * average_spacing / (1000.0 * 5.0 * radius))
    );
    
    Mesh_criteria criteria(
        CGAL::parameters::facet_angle = sm_angle,
        CGAL::parameters::facet_size = sm_radius * average_spacing,
        CGAL::parameters::facet_distance = sm_distance * average_spacing
    );
    
    C3t3 c3t3 = CGAL::make_mesh_3<C3t3>(
        domain, criteria,
        CGAL::parameters::surface_only()
                         .manifold_with_boundary()
    );
    
    // 7. 提取并保存结果
    Polyhedron output_mesh;
    CGAL::facets_in_complex_3_to_triangle_mesh(c3t3, output_mesh);
    
    std::ofstream out("output.off");
    out << output_mesh;
    out.close();
    
    std::cout << "Reconstruction completed. Output saved to output.off" << std::endl;
    
    return EXIT_SUCCESS;
}

6.2 便捷函数接口

CGAL 还提供了更简单的便捷函数：

#include <CGAL/poisson_surface_reconstruction.h>
 
// 使用便捷函数进行重建
bool success = CGAL::poisson_surface_reconstruction_delaunay(
    points.begin(), points.end(),           // 点云范围
    Point_map(), Normal_map(),              // 属性映射
    output_mesh,                            // 输出网格
    spacing,                                // 间距参数
    20.0,                                   // 角度参数
    30.0,                                   // 半径参数
    0.375,                                  // 距离参数
    CGAL::Manifold_with_boundary_tag()      // 流形标签
);

6.3 实际应用：文物数字化

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Surface_mesh.h>
#include <CGAL/Poisson_reconstruction_function.h>
#include <CGAL/IO/read_points.h>
#include <CGAL/IO/OBJ.h>
#include <CGAL/compute_average_spacing.h>
#include <CGAL/edge_aware_upsampling_point_set.h>
#include <CGAL/jet_estimate_normals.h>
#include <CGAL/mst_orient_normals.h>
#include <iostream>
#include <vector>
 
typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel;
typedef Kernel::Point_3 Point;
typedef Kernel::Vector_3 Vector;
typedef std::pair<Point, Vector> Point_with_normal;
typedef CGAL::Surface_mesh<Point> Mesh;
 
// 文物数字化完整流程
void digitize_cultural_heritage(const std::string& input_file, 
                                const std::string& output_file)
{
    std::vector<Point_with_normal> points;
    
    // 1. 读取扫描数据
    if (\!CGAL::IO::read_points(input_file, std::back_inserter(points),
                               CGAL::parameters::point_map(CGAL::First_of_pair_property_map<Point_with_normal>())
                                                .normal_map(CGAL::Second_of_pair_property_map<Point_with_normal>())))
    {
        std::cerr << "Failed to read input" << std::endl;
        return;
    }
    
    std::cout << "Loaded " << points.size() << " points" << std::endl;
    
    // 2. 预处理：上采样（如果点云太稀疏）
    if (points.size() < 10000) {
        std::vector<Point_with_normal> upsampled;
        CGAL::edge_aware_upsampling_point_set(
            points,
            std::back_inserter(upsampled),
            CGAL::parameters::point_map(CGAL::First_of_pair_property_map<Point_with_normal>())
                             .normal_map(CGAL::Second_of_pair_property_map<Point_with_normal>())
                             .sharpness_angle(25.0)
                             .edge_sensitivity(0.5)
                             .neighbor_radius(2.0)
                             .number_of_output_points(points.size() * 4)
        );
        points.swap(upsampled);
        std::cout << "Upsampled to " << points.size() << " points" << std::endl;
    }
    
    // 3. 预处理：法向估计（如果没有）
    // 法向定向
    CGAL::mst_orient_normals(
        points,
        12,  // K-近邻数
        CGAL::parameters::point_map(CGAL::First_of_pair_property_map<Point_with_normal>())
                         .normal_map(CGAL::Second_of_pair_property_map<Point_with_normal>())
    );
    
    // 4. 泊松重建
    CGAL::Poisson_reconstruction_function<Kernel> function(
        points.begin(), points.end(),
        CGAL::First_of_pair_property_map<Point_with_normal>(),
        CGAL::Second_of_pair_property_map<Point_with_normal>()
    );
    
    if (\!function.compute_implicit_function()) {
        std::cerr << "Poisson reconstruction failed" << std::endl;
        return;
    }
    
    // 5. 参数设置（根据文物细节调整）
    FT average_spacing = CGAL::compute_average_spacing<CGAL::Sequential_tag>(
        points, 6, CGAL::parameters::point_map(CGAL::First_of_pair_property_map<Point_with_normal>()));
    
    // 对于精细文物，使用更小的间距
    FT spacing = average_spacing * 0.5;  // 更精细的网格
    FT angle = 20.0;                      // 标准角度
    FT radius = 20.0;                     // 更小的三角形
    FT distance = 0.2;                    // 更精确的表面逼近
    
    // 6. 网格生成
    Mesh mesh;
    bool success = CGAL::poisson_surface_reconstruction_delaunay(
        points.begin(), points.end(),
        CGAL::First_of_pair_property_map<Point_with_normal>(),
        CGAL::Second_of_pair_property_map<Point_with_normal>(),
        mesh,
        spacing,
        angle,
        radius,
        distance,
        CGAL::Manifold_with_boundary_tag()
    );
    
    if (success) {
        CGAL::IO::write_OBJ(output_file, mesh);
        std::cout << "Saved to " << output_file << std::endl;
        std::cout << "Vertices: " <> mesh.number_of_vertices() << std::endl;
        std::cout << "Faces: " << mesh.number_of_faces() << std::endl;
    }
}
 
int main(int argc, char* argv[])
{
    std::string input = (argc > 1) ? argv[1] : "artifact.xyz";
    std::string output = (argc > 2) ? argv[2] : "artifact.obj";
    
    digitize_cultural_heritage(input, output);
    
    return 0;
}

7. 复杂度分析

7.1 时间复杂度

步骤	复杂度	说明
三角剖分构建	$O (n lo g n)$	$n$ 为输入点数
Delaunay 细化	$O (m lo g m)$	$m$ 为细化后顶点数
线性系统组装	$O (m)$	稀疏矩阵操作
线性系统求解	$O (m \cdot k)$	$k$ 为迭代次数（共轭梯度）
等值面提取	$O (t)$	$t$ 为输出三角形数

总体复杂度约为 $O (n lo g n + m \cdot k)$ 。

7.2 空间复杂度

三角剖分： $O (m)$
稀疏矩阵： $O (m)$ （平均每个顶点约 15 个非零元）
临时存储： $O (m)$

总体空间复杂度为 $O (m)$ ，其中 $m$ 通常与 $n$ 同阶或略大。

7.3 参数影响

spacing：控制输出网格分辨率。较小的值产生更精细的网格，但增加计算时间和内存使用。
sm_angle：控制三角形最小角度。较大的值产生更规则的三角形，但可能增加顶点数。
sm_distance：控制表面逼近精度。较小的值产生更精确的表面，但需要更多三角形。

8. 应用建议

8.1 适用场景

泊松重建最适合以下场景：

封闭物体的扫描数据
具有法向信息的点云
需要水密输出的应用
允许输出网格不经过输入点

8.2 预处理建议

法向估计：如果输入没有法向，使用 PCA 或 jet 拟合估计
法向定向：确保法向一致朝外（可使用 MST 传播）
去噪：移除离群点，减少噪声
采样：对于非常密集的点云，可适当下采样

8.3 常见问题

问题 1：重建结果有孔洞

原因：输入点云不完整或法向不一致
解决：启用 smoother_hole_filling 选项，或修复输入数据

问题 2：表面过于平滑，丢失细节

原因：spacing 参数过大
解决：减小 spacing 值，或使用更小的 sm_distance

问题 3：重建失败（线性系统求解失败）

原因：输入点过少或法向全为零
解决：检查输入数据质量

9. 参数调优指南

9.1 根据数据特征选择参数

数据特征	推荐参数	说明
高精度扫描	spacing=0.3, distance=0.2	保留细节
噪声数据	spacing=0.8, distance=0.5	平滑噪声
大规模数据	spacing=1.0, max_vertices=1M	控制内存
精细特征	spacing=0.2, angle=15	保留尖锐特征

9.2 多分辨率重建策略

// 粗粒度快速预览
void coarse_reconstruction(const Point_set& points) {
    FT spacing = CGAL::compute_average_spacing(points, 6) * 2.0;
    poisson_reconstruction(points, spacing, 20.0, 30.0, 0.5);
}
 
// 细粒度最终输出
void fine_reconstruction(const Point_set& points) {
    FT spacing = CGAL::compute_average_spacing(points, 6) * 0.5;
    poisson_reconstruction(points, spacing, 20.0, 20.0, 0.2);
}

参考文献

Kazhdan, M., Bolitho, M., & Hoppe, H. (2006). Poisson surface reconstruction. In Proceedings of the fourth Eurographics symposium on Geometry processing (pp. 61-70).
Kazhdan, M., & Hoppe, H. (2013). Screened Poisson surface reconstruction. ACM Transactions on Graphics, 32(3), 1-13.
CGAL Documentation: Poisson Surface Reconstruction. https://doc.cgal.org/latest/Poisson_surface_reconstruction_3/

Chunibyo

Explorer

13.2 泊松曲面重建

13.2 泊松曲面重建 §

0. 算法动机：为什么需要泊松重建？ §

现实问题引入 §

传统方法的局限 §

泊松重建的优势 §

如果不使用泊松重建 §

1. 直观理解：3D扫描的轮廓 §

1.1 “指示函数就像剪影” §

1.2 “梯度就像坡度指示” §

1.3 “泊松方程就像平滑填充” §

2. 数学基础 §

2.1 指示函数与梯度 §

2.2 泊松方程 §

2.3 变分形式 §

3. 算法流程详解 §

3.1 八叉树结构 §

3.2 从点云到表面的重建过程 §

3.3 向量场构造 §

3.4 散度计算 §

3.5 泊松方程求解 §

3.6 等值面提取 §

4. CGAL 架构分析 §

4.1 核心类层次 §

4.2 关键文件位置 §

4.3 Poisson_reconstruction_function 类 §

5. 实现细节 §

5.1 Delaunay 细化 §

5.2 线性系统组装 §

5.3 散度计算 §

5.4 Laplacian 矩阵组装 §

5.5 等值面提取 §

6. 使用示例 §

6.1 完整示例代码 §

6.2 便捷函数接口 §

6.3 实际应用：文物数字化 §

7. 复杂度分析 §

7.1 时间复杂度 §

7.2 空间复杂度 §

7.3 参数影响 §

8. 应用建议 §

8.1 适用场景 §

8.2 预处理建议 §

8.3 常见问题 §

9. 参数调优指南 §

9.1 根据数据特征选择参数 §

9.2 多分辨率重建策略 §

参考文献 §

Graph View

Table of Contents

Backlinks

13.2 泊松曲面重建

0. 算法动机：为什么需要泊松重建？

现实问题引入

传统方法的局限

泊松重建的优势

如果不使用泊松重建

1. 直观理解：3D扫描的轮廓

1.1 “指示函数就像剪影”

1.2 “梯度就像坡度指示”

1.3 “泊松方程就像平滑填充”

2. 数学基础

2.1 指示函数与梯度

2.2 泊松方程

2.3 变分形式

3. 算法流程详解

3.1 八叉树结构

3.2 从点云到表面的重建过程

3.3 向量场构造

3.4 散度计算

3.5 泊松方程求解

3.6 等值面提取

4. CGAL 架构分析

4.1 核心类层次

4.2 关键文件位置

4.3 Poisson_reconstruction_function 类

5. 实现细节

5.1 Delaunay 细化

5.2 线性系统组装

5.3 散度计算

5.4 Laplacian 矩阵组装

5.5 等值面提取

6. 使用示例

6.1 完整示例代码

6.2 便捷函数接口

6.3 实际应用：文物数字化

7. 复杂度分析

7.1 时间复杂度

7.2 空间复杂度

7.3 参数影响

8. 应用建议

8.1 适用场景

8.2 预处理建议

8.3 常见问题

9. 参数调优指南

9.1 根据数据特征选择参数

9.2 多分辨率重建策略

参考文献