Open3D 最小二乘法拟合二次曲面

目录

一、概述

1.1原理

1.2实现步骤

二、代码实现

2.1关键代码

2.2完整代码

三、实现效果

3.1原始点云

3.2拟合后点云


前期试读,后续会将博客加入该专栏,欢迎订阅

Open3D点云算法与点云深度学习案例汇总(长期更新)-CSDN博客

一、概述

1.1原理

        最小二乘法是一种统计方法,通过最小化观测值与模型之间的平方误差来估计模型参数。在拟合二次曲面时,我们假设曲面方程为:

1.2实现步骤

1.生成或读取点云数据:
        - 生成带有噪声的二次曲面点云数据。
        - 或者从文件中读取点云数据。
2.构建设计矩阵𝐴和观测向量 𝐵:
        - 设计矩阵 A 包含点云数据的各项多项式项。
        - 观测向量 𝐵 包含点云数据的 𝑧坐标。
3.解线性方程 𝐴𝑋=𝐵:
        - 使用最小二乘法求解该方程,获得拟合的二次曲面参数。
4.创建二次曲面的 Mesh:
        - 使用拟合的参数,计算指定范围和分辨率下的二次曲面点。
        - 生成三角形网格,用于可视化二次曲面。
5.可视化点云和拟合的二次曲面:
        - 使用 Open3D 可视化原始点云数据和拟合的二次曲面。

二、代码实现

2.1关键代码

def fit_quadratic_surface(points):
    """
    使用最小二乘法拟合二次曲面。

    参数:
    points (numpy.ndarray): 点云数据,形状为 (N, 3)。

    返回:
    numpy.ndarray: 拟合的二次曲面的系数。
    """
    A = np.c_[points[:, 0]**2, points[:, 1]**2, points[:, 0] * points[:, 1],
              points[:, 0], points[:, 1], np.ones(points.shape[0])]
    B = points[:, 2]

    coeff, _, _, _ = np.linalg.lstsq(A, B, rcond=None)
    return coeff

2.2完整代码

import open3d as o3d
import numpy as np

def generate_noisy_paraboloid(a, b, c, d, e, f, g, num_points=1000, noise_level=0.1):
    """
    生成带有噪声的二次曲面点云数据。

    参数:
    a, b, c, d, e, f, g: 二次曲面的系数。
    num_points (int): 点的数量。
    noise_level (float): 噪声水平。

    返回:
    numpy.ndarray: 生成的点云数据。
    """
    x = np.random.uniform(-1, 1, num_points)
    y = np.random.uniform(-1, 1, num_points)
    z = a * x**2 + b * y**2 + c * x * y + d * x + e * y + f + np.random.normal(0, noise_level, num_points)
    points = np.vstack((x, y, z)).T

    return points

def fit_quadratic_surface(points):
    """
    使用最小二乘法拟合二次曲面。

    参数:
    points (numpy.ndarray): 点云数据,形状为 (N, 3)。

    返回:
    numpy.ndarray: 拟合的二次曲面的系数。
    """
    A = np.c_[points[:, 0]**2, points[:, 1]**2, points[:, 0] * points[:, 1],
              points[:, 0], points[:, 1], np.ones(points.shape[0])]
    B = points[:, 2]

    coeff, _, _, _ = np.linalg.lstsq(A, B, rcond=None)
    return coeff

def create_quadratic_surface_mesh(coeff, x_range, y_range, resolution=100):
    """
    创建二次曲面的 Mesh,用于可视化。

    参数:
    coeff (numpy.ndarray): 二次曲面的系数。
    x_range (tuple): x 轴的范围。
    y_range (tuple): y 轴的范围。
    resolution (int): 网格分辨率。

    返回:
    open3d.geometry.TriangleMesh: 二次曲面 Mesh 对象。
    """
    a, b, c, d, e, f = coeff

    x = np.linspace(x_range[0], x_range[1], resolution)
    y = np.linspace(y_range[0], y_range[1], resolution)
    X, Y = np.meshgrid(x, y)
    Z = a * X**2 + b * Y**2 + c * X * Y + d * X + e * Y + f

    points = np.vstack((X.ravel(), Y.ravel(), Z.ravel())).T
    faces = []

    for i in range(resolution - 1):
        for j in range(resolution - 1):
            idx = i * resolution + j
            faces.append([idx, idx + 1, idx + resolution])
            faces.append([idx + 1, idx + resolution + 1, idx + resolution])

    mesh = o3d.geometry.TriangleMesh()
    mesh.vertices = o3d.utility.Vector3dVector(points)
    mesh.triangles = o3d.utility.Vector3iVector(faces)
    mesh.compute_vertex_normals()

    return mesh

# 生成带有噪声的二次曲面点云数据
a, b, c, d, e, f = 1, 1, 0, 0, 0, 0
num_points = 1000
noise_level = 0.1
points = generate_noisy_paraboloid(a, b, c, d, e, f, 0, num_points, noise_level)

# 使用最小二乘法拟合二次曲面
coeff = fit_quadratic_surface(points)

# 创建点云对象
pcd = o3d.geometry.PointCloud()
pcd.points = o3d.utility.Vector3dVector(points)
o3d.visualization.draw_geometries([pcd], window_name="pcd",
                                  width=800, height=600, left=50, top=50)
# 创建拟合的二次曲面 Mesh 对象
x_range = (-1, 1)
y_range = (-1, 1)
surface_mesh = create_quadratic_surface_mesh(coeff, x_range, y_range)

# 可视化点云和拟合的二次曲面
o3d.visualization.draw_geometries([pcd, surface_mesh], window_name="Least Squares Quadratic Surface Fitting",
                                  width=800, height=600, left=50, top=50)

三、实现效果

3.1原始点云

3.2拟合后点云

相关推荐

  1. 点云乘法直线 Matlab

    2024-07-18 07:26:03       59 阅读
  2. 乘法-平面方程

    2024-07-18 07:26:03       25 阅读

最近更新

  1. docker php8.1+nginx base 镜像 dockerfile 配置

    2024-07-18 07:26:03       101 阅读
  2. Could not load dynamic library ‘cudart64_100.dll‘

    2024-07-18 07:26:03       109 阅读
  3. 在Django里面运行非项目文件

    2024-07-18 07:26:03       87 阅读
  4. Python语言-面向对象

    2024-07-18 07:26:03       96 阅读

热门阅读

  1. uniapp小程序项目解决键盘问题

    2024-07-18 07:26:03       24 阅读
  2. C# 类型的默认值

    2024-07-18 07:26:03       23 阅读
  3. [PostgreSql]获取表结构数据

    2024-07-18 07:26:03       22 阅读
  4. 设计模式-工厂设计

    2024-07-18 07:26:03       23 阅读
  5. 构建完成,通知我:在Gradle中配置构建通知

    2024-07-18 07:26:03       21 阅读