Coverage for trimesh/smoothing.py: 93%
112 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-24 04:40 +0000
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-24 04:40 +0000
1import numpy as np
3from .typed import ArrayLike
5try:
6 from scipy.sparse import coo_matrix, eye
7 from scipy.sparse.linalg import spsolve
8except ImportError as E:
9 from .exceptions import ExceptionWrapper
11 wrapper = ExceptionWrapper(E)
12 eye, spsolve, coo_matrix = wrapper, wrapper, wrapper
14from . import graph, triangles
15from .base import Trimesh
16from .geometry import index_sparse
17from .triangles import mass_properties
18from .util import unitize
21def filter_laplacian(
22 mesh,
23 lamb=0.5,
24 iterations=10,
25 implicit_time_integration=False,
26 volume_constraint=True,
27 laplacian_operator=None,
28):
29 """
30 Smooth a mesh in-place using laplacian smoothing.
31 Articles
32 1 - "Improved Laplacian Smoothing of Noisy Surface Meshes"
33 J. Vollmer, R. Mencl, and H. Muller
34 2 - "Implicit Fairing of Irregular Meshes using Diffusion
35 and Curvature Flow". M. Desbrun, M. Meyer,
36 P. Schroder, A.H.B. Caltech
37 Parameters
38 ------------
39 mesh : trimesh.Trimesh
40 Mesh to be smoothed in place
41 lamb : float
42 Diffusion speed constant
43 If 0.0, no diffusion
44 If > 0.0, diffusion occurs
45 implicit_time_integration: boolean
46 if False: explicit time integration
47 -lamb <= 1.0 - Stability Limit (Article 1)
48 if True: implicit time integration
49 -lamb no limit (Article 2)
50 iterations : int
51 Number of passes to run filter
52 laplacian_operator : None or scipy.sparse.coo.coo_matrix
53 Sparse matrix laplacian operator
54 Will be autogenerated if None
55 """
57 # if the laplacian operator was not passed create it here
58 if laplacian_operator is None:
59 laplacian_operator = laplacian_calculation(mesh)
61 # save initial volume
62 if volume_constraint:
63 vol_ini = mesh.volume
65 # get mesh vertices and faces as vanilla numpy array
66 vertices = mesh.vertices.copy().view(np.ndarray)
67 faces = mesh.faces.copy().view(np.ndarray)
69 # Set matrix for linear system of equations
70 if implicit_time_integration:
71 dlap = laplacian_operator.shape[0]
72 AA = eye(dlap) + lamb * (eye(dlap) - laplacian_operator)
74 # Number of passes
75 for _index in range(iterations):
76 # Classic Explicit Time Integration - Article 1
77 if not implicit_time_integration:
78 dot = laplacian_operator.dot(vertices) - vertices
79 vertices += lamb * dot
81 # Implicit Time Integration - Article 2
82 else:
83 vertices = spsolve(AA, vertices)
85 # volume constraint
86 if volume_constraint:
87 # find the volume with new vertex positions
88 vol_new = triangles.mass_properties(vertices[faces], skip_inertia=True)[
89 "volume"
90 ]
91 # scale by volume ratio
92 vertices *= (vol_ini / vol_new) ** (1.0 / 3.0)
94 # assign modified vertices back to mesh
95 mesh.vertices = vertices
96 return mesh
99def filter_humphrey(mesh, alpha=0.1, beta=0.5, iterations=10, laplacian_operator=None):
100 """
101 Smooth a mesh in-place using laplacian smoothing
102 and Humphrey filtering.
103 Articles
104 "Improved Laplacian Smoothing of Noisy Surface Meshes"
105 J. Vollmer, R. Mencl, and H. Muller
106 Parameters
107 ------------
108 mesh : trimesh.Trimesh
109 Mesh to be smoothed in place
110 alpha : float
111 Controls shrinkage, range is 0.0 - 1.0
112 If 0.0, not considered
113 If 1.0, no smoothing
114 beta : float
115 Controls how aggressive smoothing is
116 If 0.0, no smoothing
117 If 1.0, full aggressiveness
118 iterations : int
119 Number of passes to run filter
120 laplacian_operator : None or scipy.sparse.coo.coo_matrix
121 Sparse matrix laplacian operator
122 Will be autogenerated if None
123 """
124 # if the laplacian operator was not passed create it here
125 if laplacian_operator is None:
126 laplacian_operator = laplacian_calculation(mesh)
128 # get mesh vertices as vanilla numpy array
129 vertices = mesh.vertices.copy().view(np.ndarray)
130 # save original unmodified vertices
131 original = vertices.copy()
133 # run through iterations of filter
134 for _index in range(iterations):
135 vert_q = vertices.copy()
136 vertices = laplacian_operator.dot(vertices)
137 vert_b = vertices - (alpha * original + (1.0 - alpha) * vert_q)
138 vertices -= beta * vert_b + (1.0 - beta) * laplacian_operator.dot(vert_b)
140 # assign modified vertices back to mesh
141 mesh.vertices = vertices
142 return mesh
145def filter_taubin(mesh, lamb=0.5, nu=0.5, iterations=10, laplacian_operator=None):
146 """
147 Smooth a mesh in-place using laplacian smoothing
148 and taubin filtering.
149 Articles
150 "Improved Laplacian Smoothing of Noisy Surface Meshes"
151 J. Vollmer, R. Mencl, and H. Muller
152 Parameters
153 ------------
154 mesh : trimesh.Trimesh
155 Mesh to be smoothed in place.
156 lamb : float
157 Controls shrinkage, range is 0.0 - 1.0
158 nu : float
159 Controls dilation, range is 0.0 - 1.0
160 Nu shall be between 0.0 < 1.0/lambda - 1.0/nu < 0.1
161 iterations : int
162 Number of passes to run the filter
163 laplacian_operator : None or scipy.sparse.coo.coo_matrix
164 Sparse matrix laplacian operator
165 Will be autogenerated if None
166 """
167 # if the laplacian operator was not passed create it here
168 if laplacian_operator is None:
169 laplacian_operator = laplacian_calculation(mesh)
171 # get mesh vertices as vanilla numpy array
172 vertices = mesh.vertices.copy().view(np.ndarray)
174 # run through multiple passes of the filter
175 for index in range(iterations):
176 # do a sparse dot product on the vertices
177 dot = laplacian_operator.dot(vertices) - vertices
178 # alternate shrinkage and dilation
179 if index % 2 == 0:
180 vertices += lamb * dot
181 else:
182 vertices -= nu * dot
184 # assign updated vertices back to mesh
185 mesh.vertices = vertices
186 return mesh
189def filter_mut_dif_laplacian(
190 mesh, lamb=0.5, iterations=10, volume_constraint=True, laplacian_operator=None
191):
192 """
193 Smooth a mesh in-place using laplacian smoothing using a
194 mutable diffusion laplacian.
196 Articles
197 Barroqueiro, B., Andrade-Campos, A., Dias-de-Oliveira,
198 J., and Valente, R. (January 21, 2021).
199 "Bridging between topology optimization and additive
200 manufacturing via Laplacian smoothing." ASME. J. Mech. Des.
203 Parameters
204 ------------
205 mesh : trimesh.Trimesh
206 Mesh to be smoothed in place
207 lamb : float
208 Diffusion speed constant
209 If 0.0, no diffusion
210 If > 0.0, diffusion occurs
211 iterations : int
212 Number of passes to run filter
213 laplacian_operator : None or scipy.sparse.coo.coo_matrix
214 Sparse matrix laplacian operator
215 Will be autogenerated if None
216 """
218 # if the laplacian operator was not passed create it here
219 if laplacian_operator is None:
220 laplacian_operator = laplacian_calculation(mesh)
222 # Set volume constraint
223 if volume_constraint:
224 v_ini = mesh.volume
226 # get mesh vertices as vanilla numpy array
227 vertices = mesh.vertices.copy().view(np.ndarray)
228 faces = mesh.faces.copy().view(np.ndarray)
229 eps = 0.01 * (np.max(mesh.area_faces) ** 0.5)
231 # Number of passes
232 for _index in range(iterations):
233 # Mutable diffusion
234 normals = get_vertices_normals(mesh)
235 qi = laplacian_operator.dot(vertices)
236 pi_qi = vertices - qi
238 adil = np.abs((normals * pi_qi).dot(np.ones((3, 1))))
239 adil = 1.0 / np.maximum(1e-12, adil)
240 lamber = np.maximum(0.2 * lamb, np.minimum(1.0, lamb * adil / np.mean(adil)))
242 # Filter
243 dot = laplacian_operator.dot(vertices)
244 vertices += lamber * (dot - vertices)
246 # Volume constraint
247 if volume_constraint:
248 vol = mass_properties(vertices[faces], skip_inertia=True)["volume"]
249 if _index == 0:
250 slope = dilate_slope(vertices, faces, normals, vol, eps)
251 vertices += normals * slope * (v_ini - vol)
253 # assign modified vertices back to mesh
254 mesh.vertices = vertices
256 return mesh
259def laplacian_calculation(
260 mesh: Trimesh,
261 equal_weight: bool = True,
262 pinned_vertices: ArrayLike | None = None,
263):
264 """
265 Calculate a sparse matrix for laplacian operations.
267 Note that setting equal_weight to False significantly hampers performance.
269 Parameters
270 -------------
271 mesh : trimesh.Trimesh
272 Input geometry
273 equal_weight : bool
274 If True, all neighbors will be considered equally
275 If False, all neighbors will be weighted by inverse distance
276 pinned_vertices : None or list of ints
277 If None, no vertices are pinned
278 If list, vertices will be pinned, such that they will not be moved
279 Returns
280 ----------
281 laplacian : scipy.sparse.coo.coo_matrix
282 Laplacian operator
283 """
284 if equal_weight:
285 laplacian = graph.edges_to_coo(mesh.edges)
287 if pinned_vertices is not None:
288 # Set pinned vertices to only have themselves as neighbours
289 laplacian.data[np.isin(laplacian.row, pinned_vertices)] = False
290 laplacian.row = np.concatenate((laplacian.row, pinned_vertices))
291 laplacian.col = np.concatenate((laplacian.col, pinned_vertices))
292 laplacian.data = np.concatenate(
293 (laplacian.data, np.ones(len(pinned_vertices), dtype=bool))
294 )
296 laplacian = laplacian / laplacian.sum(axis=1)
298 else:
299 # get the vertex neighbors from the cache
300 neighbors = mesh.vertex_neighbors
302 # if a node is pinned, it will average his coordinates by himself
303 # in practice it will not move
304 if pinned_vertices is not None:
305 for i in pinned_vertices:
306 neighbors[i] = [i]
308 # avoid hitting crc checks in loops
309 vertices = mesh.vertices.view(np.ndarray)
311 # stack neighbors to 1D arrays
312 col = np.concatenate(neighbors)
313 row = np.concatenate([[i] * len(n) for i, n in enumerate(neighbors)])
315 # umbrella weights, distance-weighted
316 # use dot product of ones to replace array.sum(axis=1)
317 ones = np.ones(3)
318 # the distance from verticesex to neighbors
319 norms = [
320 1.0
321 / np.maximum(1e-6, np.sqrt(np.dot((vertices[i] - vertices[n]) ** 2, ones)))
322 for i, n in enumerate(neighbors)
323 ]
324 # normalize group and stack into single array
325 data = np.concatenate([i / i.sum() for i in norms])
327 # create the sparse matrix
328 laplacian = coo_matrix((data, (row, col)), shape=[len(vertices)] * 2)
330 return laplacian
333def get_vertices_normals(mesh):
334 """
335 Compute Vertex normals using equal weighting of neighbors faces.
336 Parameters
337 -------------
338 mesh : trimesh.Trimesh
339 Input geometry
340 Returns
341 ----------
342 vertices_normals: array
343 Vertices normals
344 """
346 # get mesh vertices and faces
347 vertices = mesh.vertices
348 faces = mesh.faces
350 # get face normals
351 face_normals = mesh.face_normals
353 # Compute Vert normals
354 vert_normals = index_sparse(len(vertices), faces).dot(face_normals)
356 return unitize(vert_normals)
359def dilate_slope(vertices, faces, normals, v, eps):
360 """
361 Get the derivate of dilation scalar by the volume variation by finite differences
362 Thus, Vertices += vertex_normals*dilate_slope*(Initial_Volume - Srinked_Volume)
363 Parameters
364 -------------
365 mesh : trimesh.Trimesh
366 Input geometry
367 vertices: mesh.vertices
368 faces: mesh.faces
369 normals: array
370 vertices normals
371 Returns
372 ----------
373 dilate_slope: float
374 derivative
375 """
377 # finite difference derivative
378 vertices2 = vertices + normals * eps
379 v2 = mass_properties(vertices2[faces], skip_inertia=True)["volume"]
381 return (eps) / (v2 - v)