点云-L1骨架算法

最新推荐文章于 2025-03-02 11:19:11 发布

我是瓦力

最新推荐文章于 2025-03-02 11:19:11 发布

阅读量2.1k

点赞数 29

分类专栏： 2024最新点云处理文章标签：算法 numpy python 计算机视觉 3d

本文链接：https://blog.csdn.net/u012901740/article/details/144107523

版权

2024最新点云处理专栏收录该内容

18 篇文章

订阅专栏

L1骨架算法通过最小化点到骨架的距离，把点云的“外壳”剥开，留下最核心的“线条”或“骨架”，像用最少的笔画画出物体的大致形状。

点云中的L1骨架算法可以简单地理解为一种提取点云数据主要骨架结构的方法，就像从一团毛线中理出它的中心线。下面用通俗的方式解释它的核心思路和步骤：

L1骨架的优点
- 抗噪声能力强：不需要高质量的点云，也能得到合理的骨架。
- 简洁高效：直接提炼出结构化的骨架模型，方便后续的分析和处理。
应用场景
- 三维建模：将复杂点云简化为骨架用于建模。
- 形状分析：提取物体的几何特征。
- 机器人导航：骨架可以用作导航路径。

什么是L1骨架算法？

L1骨架算法的核心目标是找到一个点云的“骨架”——一种简单而紧凑的几何结构，能够描述点云的整体形状。它强调通过最小化点到骨架的距离（L1范数），逐步提炼出骨架线条。

结果是什么样的？

如果点云是树干形状，L1骨架会提取出一条接近树干中心的线。
如果点云是网状或者复杂结构，L1骨架会形成类似树状或网状的骨架图。

通俗步骤讲解

点云就像撒在地上的豆子

想象一团点云就像散在桌子上的一堆豆子。这些豆子是某个物体表面上的采样点。

提取中心点（骨架）

我们的目标是找到这堆豆子中最能代表它整体形状的“主线”（骨架），比如一条线或树状结构。

怎么找到这条骨架？

用吸引点的方法：
假设骨架是由一组“种子点”构成的，这些种子点会“吸引”点云中的豆子。吸引的规则是通过最小化距离来决定。
L1范数的作用：
在这个过程中，我们使用L1范数（点到种子点的直线距离）来优化骨架的位置和形状。相对于其他距离范数，L1范数对点云中的噪声和稀疏分布更鲁棒。

逐步优化

算法从初始的随机骨架点开始，通过不断调整骨架的位置，让它更贴合点云数据，最终提炼出最优的骨架结构。

具体实现步骤

L1骨架算法的实现可以分为以下几个主要步骤，每一步都旨在提取点云的核心骨架。我们从直观的思路逐步深入技术细节：

1. 点云数据初始化

输入：三维空间中的点云数据（一个包含很多点的位置坐标的集合）。
目的：为骨架的提取提供原始数据。

2. 种子点初始化（骨架候选点）

思路：在点云范围内，初始化一组种子点（候选骨架点）。
- 这些种子点可以通过均匀采样点云空间生成，也可以基于密度较高的区域优先选点。
目标：这些种子点将逐步优化，最终收敛为骨架点。

3. 分配吸引关系（Voronoi分区）

概念：将点云中的每个点分配给其最近的种子点（骨架候选点）。
- 这类似于构建一个基于L1距离的Voronoi图，即每个种子点吸引一定范围内的点云点。
实现：
- 计算点云中每个点到所有种子点的L1距离。
- 每个点归属于最近的种子点，形成分区。

4. 重心计算（骨架点更新）

思路：根据分配关系，将每个种子点更新为其吸引的点云点的“重心”。
- 重心的计算可以通过直接取分配点的坐标平均值完成。
公式：
$x_i^{(t+1)} = \frac{1}{|S_i|} \sum_{p \in S_i} p$
其中， $x_i^{(t+1)}$ 是第 $i$ 个种子点的新位置， $S_i$ 是该种子点吸引的点集合。

5. 骨架点的稀疏化

问题：种子点可能过于密集，不利于提取简洁的骨架。
解决方案：
- 使用一种稀疏化策略移除过于靠近的骨架点。
- 常用方法是设置一个最小距离阈值，将距离小于阈值的骨架点合并。

6. 迭代优化

目的：逐步优化种子点的位置，直到收敛为稳定的骨架。
流程：
- 每次迭代完成后，检查种子点位置的变化是否小于某个设定的阈值。
- 如果变化很小（如小于0.01），或者达到最大迭代次数，则停止迭代。

7. 后处理（生成连接结构）

生成骨架连线：
- 为了使骨架更具结构性，将骨架点用线段连接。
- 可以通过邻近点连接或基于最小生成树的方法生成连线。
去噪：移除骨架中孤立的点或线段。

8. 输出结果

最终输出：点云的骨架模型，包括骨架点和它们的连接关系。
结果可以用来可视化、简化建模、形状分析等。

算法伪代码

L1骨架算法能够提取点云的骨架结构，并适用于多种复杂场景。

分配吸引关系：点云点被分配到最近的骨架点。
重心更新：骨架点根据点云点的位置不断调整。
稀疏化和收敛：控制骨架点的数量并优化形状。
后处理：将骨架点连线以生成结构化结果。

# 输入：点云数据点集 P
# 输出：骨架点集 S

初始化种子点集合 S
repeat:
    # 1. 分配吸引关系
    for 每个点 p ∈ P:     
        找到离 p 最近的种子点 s ∈ S     
        将 p 分配给 s 的吸引域     

    # 2. 更新种子点位置      
    for 每个种子点 s ∈ S:      
        计算其吸引域中所有点的重心      
        更新 s 为重心位置      

    # 3. 种子点稀疏化      
    移除距离小于阈值的多余种子点      

until 骨架点的位置收敛或达到最大迭代次数      
return 骨架点集合 S

代码实现

import numpy as np 
from scipy.spatial import distance 
from scipy.sparse.csgraph import minimum_spanning_tree

import matplotlib.pyplot as plt 
from mpl_toolkits.mplot3d import Axes3D 


def extract_l1_skeleton(point_cloud, num_seeds=10, epsilon=0.01, max_iterations=100):

    """
    提取点云的L1骨架。 
    
    参数:
        point_cloud (ndarray): 输入点云数据 (n, 3) 
        num_seeds (int): 初始骨架点数量
        epsilon (float): 收敛阈值
        max_iterations (int): 最大迭代次数
    
    返回:
        skeleton_points (ndarray): 提取的骨架点 (m, 3)
    """
    # 初始化骨架点 (随机选择点云中的一些点)
    np.random.seed(42)
    indices = np.random.choice(point_cloud.shape[0], num_seeds, replace=False)
    skeleton_points = point_cloud[indices]
    
    for iteration in range(max_iterations):
        # 分配吸引关系：计算每个点到骨架点的最近距离
        distances = distance.cdist(point_cloud, skeleton_points, metric='cityblock')  # L1距离
        nearest_indices = np.argmin(distances, axis=1)  # 每个点的最近骨架点索引
        
        # 更新骨架点：计算每个骨架点的几何中位数
        new_skeleton_points = []
        for j in range(num_seeds):
            assigned_points = point_cloud[nearest_indices == j]
            if len(assigned_points) > 0:
                median = np.median(assigned_points, axis=0)
                new_skeleton_points.append(median)
        new_skeleton_points = np.array(new_skeleton_points)
        
        # 检查收敛条件：骨架点的移动是否小于阈值
        movement = np.linalg.norm(new_skeleton_points - skeleton_points, ord=1, axis=1).max()
        skeleton_points = new_skeleton_points
        if movement < epsilon:
            break

    return skeleton_points


def sparsify_skeleton(skeleton_points, min_distance=0.05):
    """
    稀疏化骨架点，移除距离小于指定阈值的点。
    
    参数:
        skeleton_points (ndarray): 骨架点集合 (m, 3)
        min_distance (float): 最小距离阈值
    
    返回:
        sparse_skeleton (ndarray): 稀疏化后的骨架点
    """
    sparse_skeleton = []
    for point in skeleton_points:
        if all(distance.cityblock(point, p) > min_distance for p in sparse_skeleton):
            sparse_skeleton.append(point)
    return np.array(sparse_skeleton)


def connect_skeleton_points_with_mst(skeleton_points):
    """
    使用最小生成树 (MST) 生成骨架点的连接线。
    
    参数:
        skeleton_points (ndarray): 骨架点集合 (m, 3)
    
    返回:
        lines (list of tuples): 每一对连接点的索引
    """
    # 计算骨架点之间的欧几里得距离矩阵
    dist_matrix = distance.cdist(skeleton_points, skeleton_points, metric='euclidean')
    
    # 构建最小生成树
    mst = minimum_spanning_tree(dist_matrix)
    
    # 提取MST的连接关系
    connections = np.array(mst.nonzero()).T  # 非零元素为连接关系
    return [(i, j) for i, j in connections]


def visualize_point_cloud_and_skeleton(point_cloud, skeleton_points, connections=None):
    """
    可视化点云、骨架点和骨架连线。
    """
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.scatter(point_cloud[:, 0], point_cloud[:, 1], point_cloud[:, 2], c='gray', s=1, label="Point Cloud")
    ax.scatter(skeleton_points[:, 0], skeleton_points[:, 1], skeleton_points[:, 2], c='red', s=20, label="Skeleton")

    # 绘制骨架连线
    if connections:
        for (i, j) in connections:
            p1 = skeleton_points[i]
            p2 = skeleton_points[j]
            ax.plot([p1[0], p2[0]], [p1[1], p2[1]], [p1[2], p2[2]], c='blue', linewidth=1)

    ax.legend()
    plt.show()


# 主程序
if __name__ == "__main__":
    # 生成接近直线的点云数据
    np.random.seed(42)
    num_points = 500

    # 基础直线：从 (0, 0, 0) 到 (1, 1, 1)
    t = np.linspace(0, 1, num_points)
    line = np.array([t, t, t]).T  # 直线数据

    # 添加随机噪声
    noise = np.random.normal(scale=0.05, size=(num_points, 3))  # 控制噪声强度
    point_cloud = line + noise

    # 提取骨架
    skeleton_points = extract_l1_skeleton(point_cloud, num_seeds=20, epsilon=0.01, max_iterations=100)

    # 稀疏化骨架
    sparse_skeleton = sparsify_skeleton(skeleton_points, min_distance=0.05)

    # 生成骨架点连接线（使用MST）
    connections = connect_skeleton_points_with_mst(sparse_skeleton)

    # 可视化结果
    visualize_point_cloud_and_skeleton(point_cloud, sparse_skeleton, connections)