Mesh deformation

If we want to deform a triangle mesh according to a small number of constraints, we can use mesh deformation algorithms. Open3D implements the as-rigid-as-possible method by [SorkineAndAlexa2007] that optimizes the following energy function

$$\begin{equation} \sum_i \sum_{j \in \mathcal{N}(i)} w_{ij} || (\mathbf{p}'_i - \mathbf{p}'_j) - \mathbf{R}_i (\mathbf{p}_i - \mathbf{p}_j)||^2 \,, \end{equation}$$

where $\mathbf{R}_i$ are the rotation matrices that we want to optimize for, and $\mathbf{p}_i$ and $\mathbf{p}'_i$ are the vertex positions before and after the optimization, respectively. $\mathcal{N}(i)$ is the set of neighbors of vertex $i$. The weights $w_{ij}$ are cotangent weights.

Open3D implements this method in deform_as_rigid_as_possible. The first argument to this method is a set of constraint_ids that refer to the vertices in the triangle mesh. The second argument constrint_pos defines at which position those vertices should be after the optimization. The optimization process is an iterative scheme. Hence, we also can define the number of iterations via max_iter.

mesh = o3dtut.get_armadillo_mesh()

vertices = np.asarray(mesh.vertices)
static_ids = [idx for idx in np.where(vertices[:, 1] < -30)[0]]
static_pos = []
for id in static_ids:
    static_pos.append(vertices[id])
handle_ids = [2490]
handle_pos = [vertices[2490] + np.array((-40, -40, -40))]
constraint_ids = o3d.utility.IntVector(static_ids + handle_ids)
constraint_pos = o3d.utility.Vector3dVector(static_pos + handle_pos)

with o3d.utility.VerbosityContextManager(
        o3d.utility.VerbosityLevel.Debug) as cm:
    mesh_prime = mesh.deform_as_rigid_as_possible(constraint_ids,
                                                  constraint_pos,
                                                  max_iter=50)

[Open3D DEBUG] [DeformAsRigidAsPossible] setting up S'
[Open3D DEBUG] [DeformAsRigidAsPossible] done setting up S'
[Open3D DEBUG] [DeformAsRigidAsPossible] setting up system matrix L
[Open3D DEBUG] [DeformAsRigidAsPossible] done setting up system matrix L
[Open3D DEBUG] [DeformAsRigidAsPossible] setting up sparse solver
[Open3D DEBUG] [DeformAsRigidAsPossible] done setting up sparse solver
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=0, energy=4.199704e+03
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=1, energy=1.902613e+03
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=2, energy=1.315890e+03
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=3, energy=1.004800e+03
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=4, energy=8.093276e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=5, energy=6.758837e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=6, energy=5.799615e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=7, energy=5.084362e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=8, energy=4.535641e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=9, energy=4.104631e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=10, energy=3.759097e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=11, energy=3.477058e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=12, energy=3.243155e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=13, energy=3.046429e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=14, energy=2.878899e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=15, energy=2.734642e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=16, energy=2.609177e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=17, energy=2.499063e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=18, energy=2.401625e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=19, energy=2.314760e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=20, energy=2.236803e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=21, energy=2.166421e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=22, energy=2.102539e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=23, energy=2.044279e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=24, energy=1.990918e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=25, energy=1.941850e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=26, energy=1.896564e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=27, energy=1.854625e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=28, energy=1.815662e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=29, energy=1.779354e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=30, energy=1.745426e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=31, energy=1.713637e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=32, energy=1.683781e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=33, energy=1.655677e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=34, energy=1.629167e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=35, energy=1.604115e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=36, energy=1.580397e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=37, energy=1.557908e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=38, energy=1.536553e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=39, energy=1.516247e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=40, energy=1.496916e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=41, energy=1.478493e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=42, energy=1.460920e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=43, energy=1.444141e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=44, energy=1.428109e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=45, energy=1.412779e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=46, energy=1.398114e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=47, energy=1.384075e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=48, energy=1.370629e+02
[Open3D DEBUG] [DeformAsRigidAsPossible] iter=49, energy=1.357746e+02

print('Original Mesh')
R = mesh.get_rotation_matrix_from_xyz((0, np.pi, 0))
o3d.visualization.draw_geometries([mesh.rotate(R, center=mesh.get_center())])
print('Deformed Mesh')
mesh_prime.compute_vertex_normals()
o3d.visualization.draw_geometries(
    [mesh_prime.rotate(R, center=mesh_prime.get_center())])

Original Mesh

Deformed Mesh

Smoothed ARAP

Open3D also implements a smoothed version of the ARAP objective defined as

$$\begin{equation} \sum_i \sum_{j \in \mathcal{N}(i)} w_{ij} || (\mathbf{p}'_i - \mathbf{p}'_j) - \mathbf{R}_i (\mathbf{p}_i - \mathbf{p}_j)||^2 + \alpha A ||\mathbf{R}_i - \mathbf{R}_j||^2\,, \end{equation}$$

that penalizes a deviation of neighboring rotation matrices. $\alpha$ is a trade-off parameter for the regularization term and $A$ is the surface area.

The smoothed objective can be used in deform_as_rigid_as_possible by using the argument energy with the parameter Smoothed.