Skip to content
Anima Mundi Anima Mundi

距离查询

Open3D 中的 RaycastingScene 类提供了一组距离查询功能，可用于将三角网格转换为隐式函数、查询到表面的距离，或判断某个点是否位于网格内部。

本教程展示如何生成这些查询，并从网格生成几何处理与三维机器学习中常用的各种隐式表示。

将网格转换为隐式表示

Section titled “将网格转换为隐式表示”

distance_field_animation.gif tutorial_geometry_distance_queries_18_1.png

初始化

Section titled “初始化”

第一步，用一个（闭合的）三角网格初始化 RaycastingScene：

# Load mesh and convert to open3d.t.geometry.TriangleMesh
armadillo_data = o3d.data.ArmadilloMesh()
mesh = o3d.io.read_triangle_mesh(armadillo_data.path)
mesh = o3d.t.geometry.TriangleMesh.from_legacy(mesh)


# Create a scene and add the triangle mesh
scene = o3d.t.geometry.RaycastingScene()
_ = scene.add_triangles(mesh)  # we do not need the geometry ID for mesh

计算单个点的距离与占据状态

Section titled “计算单个点的距离与占据状态”

RaycastingScene 直接提供了若干函数，用于计算从一个点到网格表面的无符号距离和有符号距离，还提供了一个用于计算查询点占据状态（occupancy）的函数。

query_point = o3d.core.Tensor([[10, 10, 10]], dtype=o3d.core.Dtype.Float32)


# Compute distance of the query point from the surface
unsigned_distance = scene.compute_distance(query_point)
signed_distance = scene.compute_signed_distance(query_point)
occupancy = scene.compute_occupancy(query_point)

无符号距离始终可以计算，而有符号距离和占据状态只有在网格是水密（watertight）的、且内部与外部有明确定义时才有效。如果查询点位于网格内部，则有符号距离为负值。占据状态在点位于网格外部时为 0，位于内部时为 1。

在本例中，网格是水密的，从有符号距离和占据状态的值可以看出查询点位于网格内部。

print("unsigned distance", unsigned_distance.numpy())
print("signed_distance", signed_distance.numpy())
print("occupancy", occupancy.numpy())

输出：

unsigned distance [7.8694816]
signed_distance [-7.8694816]
occupancy [1.]

计算多个点和网格的距离

Section titled “计算多个点和网格的距离”

RaycastingScene 支持一次性进行多次查询。例如，我们可以传入一个 [N,3] 的 Tensor，包含 N 个查询点，这可用于在机器学习中随机采样一个体积，以训练隐式神经表示。

min_bound = mesh.vertex.positions.min(0).numpy()
max_bound = mesh.vertex.positions.max(0).numpy()


N = 256
query_points = np.random.uniform(low=min_bound, high=max_bound,
                                 size=[N, 3]).astype(np.float32)


# Compute the signed distance for N random points
signed_distance = scene.compute_signed_distance(query_points)

RaycastingScene 允许以任意数量的前导维度来组织查询点。要查询一个网格的有符号距离，可以这样做：

xyz_range = np.linspace(min_bound, max_bound, num=32)


# query_points is a [32,32,32,3] array ..
query_points = np.stack(np.meshgrid(*xyz_range.T), axis=-1).astype(np.float32)


# signed distance is a [32,32,32] array
signed_distance = scene.compute_signed_distance(query_points)


# We can visualize a slice of the distance field directly with matplotlib
plt.imshow(signed_distance.numpy()[:, :, 15])

同样的流程也适用于 RaycastingScene.compute_distance 和 RaycastingScene.compute_occupancy，可用于生成无符号距离场和占据场。

基于最近点查询计算距离

Section titled “基于最近点查询计算距离”

距离函数是建立在 compute_closest_points() 函数之上的。本部分将重新实现有符号距离的计算，并展示如何利用 compute_closest_points() 函数返回的额外信息。

初始化

Section titled “初始化”

首先用两个三角网格初始化一个 RaycastingScene。两个网格都是水密的，并且我们将它们放置为彼此不相交。

cube = o3d.t.geometry.TriangleMesh.from_legacy(
    o3d.geometry.TriangleMesh.create_box().translate([-1.2, -1.2, 0]))
sphere = o3d.t.geometry.TriangleMesh.from_legacy(
    o3d.geometry.TriangleMesh.create_sphere(0.5).translate([0.7, 0.8, 0]))


scene = o3d.t.geometry.RaycastingScene()
# Add triangle meshes and remember ids
mesh_ids = {}
mesh_ids[scene.add_triangles(cube)] = 'cube'
mesh_ids[scene.add_triangles(sphere)] = 'sphere'

计算表面上的最近点

Section titled “计算表面上的最近点”

RaycastingScene.compute_closest_points() 可以计算相对于查询点的表面上最近点。

query_point = o3d.core.Tensor([[0, 0, 0]], dtype=o3d.core.Dtype.Float32)


# We compute the closest point on the surface for the point at position [0,0,0].
ans = scene.compute_closest_points(query_point)


# Compute_closest_points provides the point on the surface, the geometry id,
# and the primitive id.
# The dictionary keys are
#.    points
#.    geometry_ids
#.    primitive_ids
print('The closest point on the surface is', ans['points'].numpy())
print('The closest point is on the surface of the',
      mesh_ids[ans['geometry_ids'][0].item()])
print('The closest point belongs to triangle', ans['primitive_ids'][0].item())

输出：

The closest point on the surface is [[-0.2 -0.2  0. ]]
The closest point is on the surface of the cube
The closest point belongs to triangle 0

要判断点位于内部还是外部，可以从查询点出发投射一条射线，并统计交点数量：

rays = np.concatenate(
    [query_point.numpy(),
     np.ones(query_point.shape, dtype=np.float32)],
    axis=-1)
intersection_counts = scene.count_intersections(rays).numpy()
# A point is inside if the number of intersections with the scene is even
# This assumes that inside and outside are well-defined for the scene.
is_inside = intersection_counts % 2 == 1

我们可以把这些步骤组合成一个函数，从而创建一个返回额外信息的特殊有符号距离函数：

def compute_signed_distance_and_closest_goemetry(query_points: np.ndarray):
    closest_points = scene.compute_closest_points(query_points)
    distance = np.linalg.norm(query_points - closest_points['points'].numpy(),
                              axis=-1)
    rays = np.concatenate([query_points, np.ones_like(query_points)], axis=-1)
    intersection_counts = scene.count_intersections(rays).numpy()
    is_inside = intersection_counts % 2 == 1
    distance[is_inside] *= -1
    return distance, closest_points['geometry_ids'].numpy()

我们可以用这个函数生成一个包含距离和几何体 id 信息的网格：

# compute range
xyz_range = np.linspace([-2, -2, -2], [2, 2, 2], num=32)
# query_points is a [32,32,32,3] array ..
query_points = np.stack(np.meshgrid(*xyz_range.T), axis=-1).astype(np.float32)


sdf, closest_geom = compute_signed_distance_and_closest_goemetry(query_points)


# We can visualize a slice of the grids directly with matplotlib
fig, axes = plt.subplots(1, 2)
axes[0].imshow(sdf[:, :, 16])
axes[1].imshow(closest_geom[:, :, 16])

tutorial_geometry_distance_queries_33_1.png