Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
import numpy as np
class TriangulatedSurface(object):
"""
center: center of mass [m^3]
m: total mass [kg]
im: 3, 3 | inertia matrix relative to center of mass [kg m^2]
"""
def __init__(self, rho, vertices=None, triangles=None):
"""
:param rho: constant density of the solid object [kg/m^3]
:param vertices: n, 3 | vertices of 3d object [m]
:param triangles: m, 3 | triangles composing the object's surface
Assume that these are vertices are consistently ordered
For a closed body, all surface normals pointing outward is sufficient
"""
self.rho = rho
self.vertices = vertices
self.triangles = triangles
self.m = None
self.center = None
self.im = None
def calculate_mass_properties(vertices, triangles, rho):
"""
Calculate mass properties of solid object via surface integrals
following supplemental material for:
M. Bacher, et al., ”Spin-It: Optimizing Moment of
Inertia for Spinnable Objects", ACM Trans. Graphics, 2014.
:param vertices: n, 3 | vertices of 3d object [m]
:param triangles: m, 3 | triangles composing the object's surface
Assume that these are vertices are consistently ordered
:param rho: constant density of the solid object [kg/m^3]
:return:
"""
s = np.zeros((10,))
for ijk in triangles:
a, b, c = vertices[ijk]
bars = vertices[ijk, :].copy()
bars = bars[:, [1, 2, 0]]
u = b - a
v = c - a
n = np.cross(u, v)
h1 = a + b + c
h2 = a**2 + b * (a + b)
h3 = h2 + c * h1
h4 = a**3 + b * h2 + c * h3
h5 = h3 + a * (h1 + a)
h6 = h3 + b * (h1 + b)
h7 = h3 + c * (h1 + c)
h8 = bars[0] * h5 + bars[1] * h6 + bars[2] * h7
s[0] += (n * h1)[0]
s[1:4] += n * h3
s[4:7] += n * h8
s[7:] += n * h4
s[0] /= 6
s[1:4] /= 24
s[4:7] /= 120
s[7:] /= 60
s *= rho
m = s[0]
center = s[1:4] / m
im = np.array([
[s[8] + s[9], -s[4], -s[6]],
[-s[4], s[7] + s[9], -s[5]],
[-s[6], -s[5], s[7] + s[8]],
])
return m, center, im
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
def remove_unused_vertices(vertices, triangles):
"""
Remove vertices that do not form a triangle, and relabel triangle vertices
to their new indices.
:param vertices: n_points, 3 | includes all points forming triangle faces
:param triangles: n_faces, 3 | indices of vertices that form each triangle
:return: vertices, triangles
"""
used_inds = np.unique(triangles.ravel())
vertices_new = vertices[used_inds]
triangles_new = triangles.copy()
for i in range(used_inds.size):
triangles_new[triangles_new == used_inds[i]] = i
return vertices_new, triangles_new
def make_convex_shape_triangles_consistent(vertices, triangles):
"""
For a convex shape described by a triangulated surface, ensure
each triangle's normal faces outward.
:param vertices: n_points, 3 | includes all points forming triangle faces
:param triangles: n_faces, 3 | indices of vertices that form each triangle
:return: triangles: n_faces, 3 | triangles made consistent
"""
used_inds = np.unique(triangles.ravel())
center = vertices[used_inds].mean(axis=0)
triangle_centers = np.array([
vertices[triangles[i]].mean(axis=0) for i in range(len(triangles))
])
outward_vec = triangle_centers - center
triangle_normals = np.cross(
vertices[triangles[:, 1]] - vertices[triangles[:, 0]],
vertices[triangles[:, 2]] - vertices[triangles[:, 0]],
)
is_mismatch = np.einsum('ij,ij->i', outward_vec, triangle_normals) < 0
triangles[is_mismatch, :] = triangles[is_mismatch, :][:, [0, 2, 1]]
return triangles