使用Open3D空间雕刻算法进行实体体素化结果异常的原因分析

问题描述 #

使用Open3D 0.18.0实现空间雕刻算法，对斯坦福数据库的Dragon三角网格模型进行实体体素化，运行代码后出现多条窗口尺寸不匹配警告，最终生成的VoxelGrid仅包含2184个体素，结果严重不符合预期。

体素化结果图：

实现代码如下：

# ----------------------------------------------------------------------------
# -                        Open3D: www.open3d.org                            -
# ----------------------------------------------------------------------------
# The MIT License (MIT)
#
# Copyright (c) 2018-2021 www.open3d.org
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in
# all copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
# IN THE SOFTWARE.
# ----------------------------------------------------------------------------

import open3d as o3d
import numpy as np


def xyz_spherical(xyz):
    x = xyz[0]
    y = xyz[1]
    z = xyz[2]
    # 计算得到球体的半径
    r = np.sqrt(x * x + y * y + z * z)
    # 半径与x轴的夹角
    r_x = np.arccos(y / r)
    # 半径在y轴的夹角
    r_y = np.arctan2(z, x)
    return [r, r_x, r_y]


def get_rotation_matrix(r_x, r_y):
    rot_x = np.asarray([[1, 0, 0], [0, np.cos(r_x), -np.sin(r_x)],
                        [0, np.sin(r_x), np.cos(r_x)]])
    rot_y = np.asarray([[np.cos(r_y), 0, np.sin(r_y)], [0, 1, 0],
                        [-np.sin(r_y), 0, np.cos(r_y)]])
    return rot_y.dot(rot_x)


def get_extrinsic(xyz):
    rvec = xyz_spherical(xyz)
    # 计算该相机位姿下的旋转矩阵和平移向量并拼接成T矩阵
    r = get_rotation_matrix(rvec[1], rvec[2])
    t = np.asarray([0, 0, 2]).transpose()
    trans = np.eye(4)
    trans[:3, :3] = r
    trans[:3, 3] = t
    return trans


def preprocess(model):
    min_bound = model.get_min_bound()
    max_bound = model.get_max_bound()
    center = min_bound + (max_bound - min_bound) / 2.0
    scale = np.linalg.norm(max_bound - min_bound) / 2.0
    vertices = np.asarray(model.vertices)
    vertices -= center
    model.vertices = o3d.utility.Vector3dVector(vertices / scale)
    return model


def voxel_carving(mesh, cubic_size, voxel_resolution, w=300, h=300):
    # 计算mesh的顶点法向量
    mesh.compute_vertex_normals()
    # 创建球体
    camera_sphere = o3d.geometry.TriangleMesh.create_sphere(radius=1.0,
                                                            resolution=10)

    # o3d.visualization.draw_geometries([camera_sphere], mesh_show_back_face=True)

    # Setup dense voxel grid.
    voxel_carving = o3d.geometry.VoxelGrid.create_dense(
        width=cubic_size,
        height=cubic_size,
        depth=cubic_size,
        voxel_size=cubic_size / voxel_resolution,
        origin=[-cubic_size / 2.0, -cubic_size / 2.0, -cubic_size / 2.0],
        color=[1.0, 0.7, 0.0])

    # Rescale geometry.
    camera_sphere = preprocess(camera_sphere)
    mesh = preprocess(mesh)

    # Setup visualizer to render depthmaps.
    vis = o3d.visualization.Visualizer()
    vis.create_window(width=w, height=h, visible=False)
    vis.add_geometry(mesh)
    vis.get_render_option().mesh_show_back_face = True
    ctr = vis.get_view_control()
    param = ctr.convert_to_pinhole_camera_parameters()

    # Carve voxel grid.
    centers_pts = np.zeros((len(camera_sphere.vertices), 3))
    for cid, xyz in enumerate(camera_sphere.vertices):
        # Get new camera pose.
        trans = get_extrinsic(xyz)
        param.extrinsic = trans
        c = np.linalg.inv(trans).dot(np.asarray([0, 0, 0, 1]).transpose())
        centers_pts[cid, :] = c[:3]
        # 转换相机的参数到成open3d中的相机内外参
        ctr.convert_from_pinhole_camera_parameters(param)

        # Capture depth image and make a point cloud.
        vis.poll_events()
        vis.update_renderer()
        # 根据当前的位姿来进行渲染拍摄得到深度图
        depth = vis.capture_depth_float_buffer(False)

        # Depth map carving method.
        voxel_carving.carve_depth_map(o3d.geometry.Image(depth), param)
        print("Carve view %03d/%03d" % (cid + 1, len(camera_sphere.vertices)))
    vis.destroy_window()

    return voxel_carving

"""
流程：
1 先创建一个固定大小的稠密（dense）voxleGrid对象
2 创建一个球形用于虚拟相机的位姿来拍摄拍摄深度图
3 根据拍摄得到的深度图与相机位姿使用carve_depth_map融合到dense voxelGrid中
"""
if __name__ == "__main__":
    mesh = o3d.io.read_triangle_mesh("DragonTriangleMesh.ply")
    cubic_size = 2.0
    voxel_resolution = 128.0

    carved_voxels = voxel_carving(mesh, cubic_size, voxel_resolution)
    print("Carved voxels ...")
    print(carved_voxels)
    o3d.visualization.draw_geometries([carved_voxels])

原因分析 #

1. 虚拟相机视角数量不足 #

代码中创建的camera_sphere使用resolution=10，生成的顶点数量过少（仅约320个顶点），导致虚拟相机的拍摄视角覆盖不全面，大量体素未被多视角深度图验证，最终被错误剔除。

2. 相机内参未匹配窗口尺寸 #

创建Visualizer时指定了w=300, h=300，但未显式设置相机内参（焦距、主点），默认内参可能与窗口尺寸不匹配，引发窗口尺寸警告，同时导致深度图的像素坐标与相机投影关系错位，雕刻逻辑失效。

3. 深度图获取未对齐 #

调用capture_depth_float_buffer(False)时未启用对齐模式，获取的深度图可能与当前窗口尺寸不完全匹配，进一步加剧内参和深度图的不匹配问题。

4. 相机与模型的空间位置匹配问题 #

预处理后模型被缩放到直径为2的立方体，而相机平移向量t = [0,0,2]，相机与模型距离过近，可能导致模型部分区域超出相机视野，深度图无法完整捕捉模型信息。

解决建议 #

1. 增加相机视角数量 #

将camera_sphere的resolution调整为20或更高，例如：

camera_sphere = o3d.geometry.TriangleMesh.create_sphere(radius=1.0, resolution=20)

更多的顶点意味着更多拍摄视角，能更全面地覆盖模型表面。

2. 显式设置相机内参 #

在创建Visualizer后，手动设置与窗口尺寸匹配的内参，确保投影关系正确：

param = ctr.convert_to_pinhole_camera_parameters()
# 设置焦距和主点，这里焦距可根据需求调整，示例用500
param.intrinsic.set_intrinsics(w, h, 500, 500, w/2, h/2)
ctr.convert_from_pinhole_camera_parameters(param)

3. 启用深度图对齐 #

修改深度图获取代码，启用对齐模式：

depth = vis.capture_depth_float_buffer(True)

确保深度图与窗口尺寸完全对齐，解决尺寸不匹配警告。

4. 调整相机与模型的相对位置 #

可以适当增大相机与模型的距离，例如将t向量的z值调整为3：

t = np.asarray([0, 0, 3]).transpose()

确保模型完整处于相机视野内，提升深度图的完整性。

5. 验证模型加载正确性 #

在预处理前添加模型可视化代码，确认模型加载完整：

o3d.visualization.draw_geometries([mesh])

排除模型文件损坏或加载错误的可能。

内容的提问来源于stack exchange，提问作者user23620257