I have a TriangleMesh in open3d and I would like to crop it using a bounding box.
Open3d has the crop function, which removes triangles if they are fully or partially outside the bounding box.
Is there a function that slices triangles instead if are partially outside the bounding box?
Here is a simple 2D example (see plot below). Given the bounding box and the input triangle, the open3d crop function would simply remove the triangle. I would like a function that takes this triangle that overlaps with the bounding box and slices it. Is there such a function?
I build a solution that can crop AxisAlignedBoundingBoxes (which is what I need):
Clean crop with the below code:
import open3d as o3d
import numpy as np
def sliceplane(mesh, axis, value, direction):
# axis can be 0,1,2 (which corresponds to x,y,z)
# value where the plane is on that axis
# direction can be True or False (True means remove everything that is
# greater, False means less
# than)
vertices = np.asarray(mesh.vertices)
triangles = np.asarray(mesh.triangles)
new_vertices = list(vertices)
new_triangles = []
# (a, b) -> c
# c refers to index of new vertex that sits at the intersection between a,b
# and the boundingbox edge
# a is always inside and b is always outside
intersection_edges = dict()
# find axes to compute
axes_compute = [0,1,2]
# remove axis that the plane is on
axes_compute.remove(axis)
def compute_intersection(vertex_in_index, vertex_out_index):
vertex_in = vertices[vertex_in_index]
vertex_out = vertices[vertex_out_index]
if (vertex_in_index, vertex_out_index) in intersection_edges:
intersection_index = intersection_edges[(vertex_in_index, vertex_out_index)]
intersection = new_vertices[intersection_index]
else:
intersection = [None, None, None]
intersection[axis] = value
const_1 = (value - vertex_in[axis])/(vertex_out[axis] - vertex_in[axis])
c = axes_compute[0]
intersection[c] = (const_1 * (vertex_out[c] - vertex_in[c])) + vertex_in[c]
c = axes_compute[1]
intersection[c] = (const_1 * (vertex_out[c] - vertex_in[c])) + vertex_in[c]
assert not (None in intersection)
# save new vertice and remember that this intersection already added an edge
new_vertices.append(intersection)
intersection_index = len(new_vertices) - 1
intersection_edges[(vertex_in_index, vertex_out_index)] = intersection_index
return intersection_index
for t in triangles:
v1, v2, v3 = t
if direction:
v1_out = vertices[v1][axis] > value
v2_out = vertices[v2][axis] > value
v3_out = vertices[v3][axis] > value
else:
v1_out = vertices[v1][axis] < value
v2_out = vertices[v2][axis] < value
v3_out = vertices[v3][axis] < value
bool_sum = sum([v1_out, v2_out, v3_out])
# print(f"{v1_out=}, {v2_out=}, {v3_out=}, {bool_sum=}")
if bool_sum == 0:
# triangle completely inside --> add and continue
new_triangles.append(t)
elif bool_sum == 3:
# triangle completely outside --> skip
continue
elif bool_sum == 2:
# two vertices outside
# add triangle using both intersections
vertex_in_index = v1 if (not v1_out) else (v2 if (not v2_out) else v3)
vertex_out_1_index = v1 if v1_out else (v2 if v2_out else v3)
vertex_out_2_index = v3 if v3_out else (v2 if v2_out else v1)
# print(f"{vertex_in_index=}, {vertex_out_1_index=}, {vertex_out_2_index=}")
# small sanity check if indices sum matches
assert sum([vertex_in_index, vertex_out_1_index, vertex_out_2_index]) == sum([v1,v2,v3])
# add new triangle
new_triangles.append([vertex_in_index, compute_intersection(vertex_in_index, vertex_out_1_index),
compute_intersection(vertex_in_index, vertex_out_2_index)])
elif bool_sum == 1:
# one vertice outside
# add three triangles
vertex_out_index = v1 if v1_out else (v2 if v2_out else v3)
vertex_in_1_index = v1 if (not v1_out) else (v2 if (not v2_out) else v3)
vertex_in_2_index = v3 if (not v3_out) else (v2 if (not v2_out) else v1)
# print(f"{vertex_out_index=}, {vertex_in_1_index=}, {vertex_in_2_index=}")
# small sanity check if outdices sum matches
assert sum([vertex_out_index, vertex_in_1_index, vertex_in_2_index]) == sum([v1,v2,v3])
new_triangles.append([vertex_in_1_index, compute_intersection(vertex_in_1_index, vertex_out_index), vertex_in_2_index])
new_triangles.append([compute_intersection(vertex_in_1_index, vertex_out_index),
compute_intersection(vertex_in_2_index, vertex_out_index), vertex_in_2_index])
else:
assert False
# TODO remap indices and remove unused
mesh = o3d.geometry.TriangleMesh()
mesh.vertices = o3d.utility.Vector3dVector(np.array(new_vertices))
mesh.triangles = o3d.utility.Vector3iVector(np.array(new_triangles))
return mesh
def clean_crop_xy(mesh, min_corner, max_corner):
min_x = min(min_corner[0], max_corner[0])
min_y = min(min_corner[1], max_corner[1])
max_x = max(min_corner[0], max_corner[0])
max_y = max(min_corner[1], max_corner[1])
# mesh = sliceplane(mesh, 0, min_x, False)
mesh_sliced = sliceplane(mesh, 0, max_x, True)
mesh_sliced = sliceplane(mesh_sliced, 0, min_x, False)
mesh_sliced = sliceplane(mesh_sliced, 1, max_y, True)
mesh_sliced = sliceplane(mesh_sliced, 1, min_y, False)
# mesh_sliced = mesh_sliced.paint_uniform_color([0,0,1])
return mesh_sliced
knot_mesh = o3d.data.KnotMesh()
mesh = o3d.io.read_triangle_mesh(knot_mesh.path)
v = np.asarray(mesh.vertices)
mesh_cropped = clean_crop_xy(mesh, (-200, 0), (0, -200))
bounding_box = o3d.geometry.AxisAlignedBoundingBox((-200, 0, -1000), (0,200,1000))
mesh_crop_bad = mesh.crop(bounding_box)
o3d.visualization.draw_geometries([mesh], mesh_show_back_face=True, mesh_show_wireframe=True)
o3d.visualization.draw_geometries([mesh_crop_bad], mesh_show_back_face=True, mesh_show_wireframe=True)
o3d.visualization.draw_geometries([mesh_cropped], mesh_show_back_face=True, mesh_show_wireframe=True)