mirror of https://github.com/CGAL/cgal
491 lines
26 KiB
Plaintext
491 lines
26 KiB
Plaintext
namespace CGAL {
|
|
/*!
|
|
|
|
\mainpage User Manual
|
|
\anchor Chapter_SurfaceModeling
|
|
|
|
\cgalAutoToc
|
|
\authors Sébastien Loriot, Olga Sorkine-Hornung, Yin Xu and Ilker %O. Yaz
|
|
|
|
\image latex main_image_suggestion_4_resized.png
|
|
\image html main_image_suggestion_4_resized.png
|
|
|
|
\section Surface_Modeling_Intro Introduction
|
|
This package offers surface mesh deformation algorithms which compute new vertex positions of a surface mesh
|
|
under positional constraints of some of its vertices, without requiring any additional structure other than the surface mesh itself
|
|
|
|
This package implements the algorithm described in \cgalCite{Sorkine2007AsRigidAs} together with an
|
|
alternative energy function \cgalCite{Chao2010SimpleGeomModel}.
|
|
The algorithm minimizes a nonlinear deformation energy under positional constraints to preserve rigidity as much as possible.
|
|
The minimization of the energy relies on solving sparse linear systems and finding closest rotation matrices.
|
|
|
|
\section Surface_Modeling_Def Definitions
|
|
|
|
A <em>surface mesh deformation system</em> consists of:
|
|
- a triangulated surface mesh (<em>surface mesh</em> in the following),
|
|
- a set of vertices defining the region to deform (referred to as the <em>region-of-interest</em> and abbreviated <em>ROI</em>),
|
|
- a subset of vertices from the ROI that the user wants to move (referred to as the <em>control vertices</em>),
|
|
- a target position for each control vertex (defining the <em>deformation constraints</em>).
|
|
|
|
A vertex from the ROI that is not a control vertex is called an <em>unconstrained vertex</em>. These definitions are depicted in \cgalFigureRef{Simple_roi}.
|
|
|
|
\cgalFigureBegin{Simple_roi, roi_example.png}
|
|
The ROI features the green vertices (the unconstrained vertices) and the red vertices (the control vertices).
|
|
(Left) The initial surface mesh; (Right) A target position is defined for each control vertex,
|
|
the coordinates of the unconstrained vertices being updated by the deformation algorithm.
|
|
\cgalFigureEnd
|
|
|
|
In this package, two algorithms are implemented:
|
|
- The <em>As-Rigid-As-Possible</em> (<em>ARAP</em>) method described in \cgalCite{Sorkine2007AsRigidAs};
|
|
- The <em>Spokes and Rims</em> method \cgalCite{Chao2010SimpleGeomModel}.
|
|
|
|
Given an edge weighting scheme, both methods iteratively minimize an energy function and produce a different surface mesh at each step until convergence is reached.
|
|
|
|
<em>Spokes and Rims</em> is the default method proposed by the package. It provides an unconditional
|
|
convergence while the <em>ARAP</em> method requires the edge weights to be positive. However,
|
|
the results obtained using the <em>Spokes and Rims</em> method are more dependent on the discretization
|
|
of the deformed surface (See \cgalFigureRef{Arap_spokes_comparison}).
|
|
|
|
More details on these algorithms are provided in section \ref Surface_Modeling_Overview.
|
|
|
|
\cgalFigureBegin{Arap_spokes_comparison, arap_spokes_comparison.png}
|
|
Comparison between the As-Rigid-As-Possible and the Spokes and Rims deformation methods.
|
|
On the surface mesh of a square with spikes depicted on the left, the ROI consists of the green vertices. The control vertices are the red ones.
|
|
We translate the control vertices along the normal to the plane and observe the result produced by
|
|
the As-Rigid-As-Possible (center) and the Spokes and Rims (right) methods from the same view point. The
|
|
latter method provides unconditional convergence does not produce a symmetric result.
|
|
\cgalFigureEnd
|
|
|
|
|
|
\section Surface_Modeling_API User Interface Description
|
|
|
|
The deformation methods implemented rely on solving a sparse linear system.
|
|
The sparse matrix definition depends on the weighting scheme
|
|
and on the unconstrained and control vertices.
|
|
The right term depends only on the target positions of the control vertices.
|
|
|
|
The deformation process is handled by the class `Deform_mesh` and the
|
|
surface mesh is represented as a halfedge data structure that must
|
|
be a model of the `HalfedgeGraph` concept.
|
|
|
|
The class `Deform_mesh` provides two groups of functions for
|
|
the preprocessing (sparse matrix definition) and the deformation (right hand term definition).
|
|
|
|
\subsection Preprocess_section Preprocessing
|
|
|
|
The preprocessing consists in computing a factorization of the aforementioned sparse
|
|
matrix to speed up the linear system resolution.
|
|
It requires the ROI to be defined.
|
|
The following conventions are used for the definition of the ROI:
|
|
- A vertex inserted in the set of control vertices is inserted in the ROI;
|
|
- A control vertex erased from the ROI is no longer considered as a control vertex;
|
|
- A control vertex that is erased is not erased from the ROI.
|
|
|
|
Each time the ROI is modified, the preprocessing function `preprocess()`
|
|
must be called.
|
|
Note that if it is not done, the first deformation step calls this function automatically
|
|
and has a longer runtime compared to subsequent deformation steps.
|
|
|
|
<!--- This is a pure implementation details as it does not shows up in the ROI
|
|
This function also extends the ROI by adding as <em>fixed control vertex</em> any vertex not in the ROI that is adjacent
|
|
to a vertex of the ROI. The target position of a fixed control vertex is its current position and cannot be changed.
|
|
-->
|
|
|
|
The function `Deform_mesh::preprocess()` returns `true` if the factorization is successful, and `false` otherwise.
|
|
|
|
Rank deficiency is the main reason for failure. Typical failure cases are:
|
|
- All the vertices are in the ROI and no control vertices are set
|
|
- The weighting scheme used to fill the sparse matrix (model of `SurfaceModelingWeights`) features too many zeros and breaks the connectivity information
|
|
|
|
\cgalAdvancedBegin
|
|
The choice of the weighting scheme provides a mean to adjust the way the control vertices
|
|
influences the unconstrained vertices. The defaults provides satisfactory results in general
|
|
but other weighting schemes may be selected or designed to experiment or improve the results
|
|
in specific cases.
|
|
\cgalAdvancedEnd
|
|
|
|
\note The ROI does not have to be a connected component of the graph of the surface mesh. However, for better performances
|
|
it is better to use an individual instance of the deformation object for each connected component.
|
|
|
|
\subsection Deform_section Deformation
|
|
|
|
The deformation of the surface mesh is triggered by the displacement of the control vertices.
|
|
This is achieved through setting the target positions of the control vertices (directly or by using an affine transformation
|
|
to be applied to a control vertex or a range of control vertices).
|
|
|
|
Note that a rotation or a translation of a control vertex is always applied on its last target position set:
|
|
they are cumulative.
|
|
|
|
The deformation of the surface mesh happens when calling the function `Deform_mesh::deform()`. The number of optimization iterations
|
|
varies depending on whether the user chooses a fixed number of iterations or a stopping criterion based on the energy variation.
|
|
|
|
After the call to the deformation function, the input surface mesh is updated and the control vertices are at
|
|
their target positions and the unconstrained vertices are moved accordingly.
|
|
The function `Deform_mesh::deform()` can be called several times consecutively, in particular
|
|
if the convergence has not been reached yet (otherwise it has no effect).
|
|
|
|
\warning Vertices can be inserted into or erased from the ROI and the set of control vertices at any time.
|
|
In particular, any vertex that is no longer inside the ROI will be assigned to its original position when
|
|
`Deform_mesh::preprocess()` is first called.
|
|
The original positions can be updated by calling `Deform_mesh::overwrite_initial_geometry()` (
|
|
which will also require a new preprocessing step).
|
|
This behavior is illustrated in the video \ref Surface_Modeling_Demo_Tutorial.
|
|
|
|
\subsection Surface_modeling_arap_or_spokes_and_rims As-Rigid-As-Possible and Spokes-and-Rims Deformation Techniques
|
|
Two deformation techniques are provided by this package. This section summarizes from the user point of view what is
|
|
explained in details in the section \ref Surface_Modeling_Overview.
|
|
|
|
The As-Rigid-As-Possible deformation technique requires the use of a positive weighting scheme to garantee
|
|
the correct minimization of the energy. When using the default cotangent weighting scheme, this means that
|
|
the input surface mesh must be <em>clean</em>. That is, that for all edges in the surface mesh
|
|
the sum of the angles opposite to the edge in the incident triangles is less that \f$ \pi \f$.
|
|
If this is not the case and the targeted application allows the modification of the surface mesh connectivity,
|
|
a solution (amongst other) is to bissect (possibly recursively) the problematic edges.
|
|
See \cgalFigureRef{Cotangent_bisect_fig}.
|
|
|
|
\cgalFigureBegin{Cotangent_bisect_fig, cotangent_bisect.png}
|
|
Cotangent weight of edges. (Left) The sum of the opposite angles to the edge \f$ e \f$, \f$ \alpha + \beta \f$
|
|
is greater than \f$ \pi \f$. The cotangent weight of \f$ e \f$ is negative.
|
|
(Right) The edge \f$ e \f$ was split so that the sum of the angles opposite
|
|
to each sub-edge \f$ e_1 \f$ and \f$ e_2 \f$ is less than \f$ \pi \f$. The corresponding
|
|
cotangent weights are positive.
|
|
\cgalFigureEnd
|
|
|
|
If the mesh connectivity must be preserved, the Spokes and Rims deformation technique
|
|
is guaranteed to always correctly minimize the energy even if the weights are negative.
|
|
However, this technique is more dependent on the discretization of the deformed surface
|
|
(See \cgalFigureRef{Arap_spokes_comparison}).
|
|
|
|
\subsection Surface_modeling_examples Examples
|
|
|
|
\subsubsection SModelingExample_1 Using the Whole Surface Mesh as Region-of-Interest
|
|
In this example, the whole surface mesh is used as ROI and a few vertices are added as control vertices.
|
|
`Deform_mesh::set_target_position()` is used for setting the target positions of the control vertices.
|
|
|
|
\cgalExample{Surface_modeling/all_roi_assign_example.cpp}
|
|
|
|
\cgalFigureBegin{SModelingExample_1_results, example_1_results.png}
|
|
Deformation results when running example \ref SModelingExample_1 : `deform_1.off` and `deform_2.off`.
|
|
\cgalFigureEnd
|
|
|
|
\subsubsection SModelingExample_2 Using an Affine Transformation on a Range of Vertices
|
|
In this example, we use the functions `translate()` and `rotate()` on a range of control vertices.
|
|
Note that the translations and the rotations are defined using a 3D vector type and a quaternion type from the \ref thirdpartyEigen "Eigen library".
|
|
|
|
\cgalExample{Surface_modeling/k_ring_roi_translate_rotate_example.cpp}
|
|
|
|
\cgalFigureBegin{SModelingExample_2_results, example_2_results.png}
|
|
Deformation results when running example \ref SModelingExample_2 : `deform_1.off` and `deform_2.off`.
|
|
\cgalFigureEnd
|
|
|
|
\subsubsection SModelingExample_3 Using an Enriched Polyhedron with Ids
|
|
In the previous examples, we used an <i>enriched</i> polyhedron storing an ID
|
|
in its halfedges and vertices together with the default property maps
|
|
in the deformation object to access them.
|
|
In the following example, we show how we can use alternative property maps.
|
|
|
|
For practical performance however we recommend relying upon
|
|
the former examples instead, as using a `std::map`
|
|
to access indices increases the complexity from constant to logarithmic.
|
|
\cgalExample{Surface_modeling/deform_polyhedron_with_custom_pmap_example.cpp}
|
|
|
|
\subsubsection SModelingExample_4 Using a Custom Edge Weighting Scheme
|
|
Using a custom weighting scheme for edges is also possible if one provides a model of `SurfaceModelingWeights`.
|
|
In this example, the weight of each edge is pre-computed and an internal map is used for storing and accessing them.
|
|
|
|
|
|
Another example is given in the manual page of the concept `::SurfaceModelingWeights`.
|
|
|
|
\cgalExample{Surface_modeling/custom_weight_for_edges_example.cpp}
|
|
|
|
\section Surface_Modeling_Demo How to Use the Demo
|
|
\todo Update these videos when the demo is stable
|
|
|
|
A plugin for the polyhedron demo is available to test the algorithm. The following video tutorials explain how to use it.
|
|
|
|
\subsection Surface_Modeling_Demo_Showcase An Introduction and Real World Example
|
|
@htmlonly
|
|
<center>
|
|
<iframe width="420" height="315" src="http://www.youtube.com/embed/WvC47fVSPGU" frameborder="0" allowfullscreen></iframe>
|
|
|
|
<iframe width="420" height="315" src="http://www.youtube.com/embed/0oTi5aKdcmg" frameborder="0" allowfullscreen></iframe>
|
|
</center>
|
|
@endhtmlonly
|
|
|
|
\subsection Surface_Modeling_Demo_Tutorial Demonstration for Deleting Control Vertices / ROI and Overwrite
|
|
|
|
@htmlonly
|
|
<center>
|
|
<iframe width="420" height="315" src="http://www.youtube.com/embed/xPAG_AOLvU0" frameborder="0" allowfullscreen></iframe>
|
|
</center>
|
|
@endhtmlonly
|
|
|
|
|
|
\section Surface_Modeling_Overview Deformations Techniques, Energies and Weighting Scheme
|
|
This section gives the theorical background to make the user manual self-contained
|
|
and at the same time explains where the weights comes in. This allows advanced users
|
|
of this package to tune the weighting scheme by developing a model of the concept
|
|
`SurfaceModelingWeights` used in the class `Deform_mesh`.
|
|
|
|
\subsection Surface_Modeling_Overview_Preliminaries Preliminaries
|
|
\subsubsection Surface_Modeling_Overview_Laplacian Laplacian Representation
|
|
The <em>Laplacian representation</em> (referred to as <em>Laplace coordinates</em> in \cgalCite{botsch2010polygon})
|
|
of a vertex in a surface mesh is one way to <em>encode</em> the local neighborhood of a vertex in the surface mesh.
|
|
In this representation, a vertex \f$ \mathbf{v}_i \f$ is associated a 3D vector defined as:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
L(\mathbf{v}_i) = \sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} w_{ij}(\mathbf{v}_i - \mathbf{v}_j),
|
|
\label{eq:lap_open}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where:
|
|
- \f$N(\mathbf{v}_i)\f$ denotes the set of vertices adjacent to \f$\mathbf{v}_i\f$;
|
|
-\f$w_{ij}\f$ denotes a weight for the directed edge \f$\mathbf{v}_i\ \mathbf{v}_j\f$.
|
|
|
|
The simplest choice for the weights is the uniform scheme where \f$ w_{ij}=1/|N(\mathbf{v}_i)| \f$ for each adjacent vertex
|
|
\f$\mathbf{v}_j\f$. In this case, the Laplacian representation of a vertex is the vector between this vertex
|
|
and the centroid of its adjacent vertices (\cgalFigureRef{Simple_laplacian}).
|
|
|
|
In the surface mesh deformation context, a popular choice is the cotangent weight scheme
|
|
that derives from the discretization of the Laplace operator \cgalCite{Pinkall1993Cotangent} :
|
|
Given an edge of the surface mesh, its corresponding cotangent weight is the mean of the
|
|
cotangents of the angles opposite to the edge. It was shown to produce results that are not biased by the surface mesh of the approximated surface
|
|
\cgalCite{Sorkine2007AsRigidAs}.
|
|
|
|
\cgalFigureBegin{Simple_laplacian,simple_mesh_with_laplacian.png}
|
|
Laplacian representation of \f$ v_i \f$ with uniform weights: the red square vertex is the centroid of the vertices adjacent to \f$ v_i \f$. The Laplacian representation \f$ L(v_i) \f$ is represented as a blue vector.
|
|
\cgalFigureEnd
|
|
|
|
Considering a surface mesh with \f$n\f$ vertices, it is possible to define its <i>Laplacian representation</i> \f$\Delta\f$ as a \f$n \times 3\f$ matrix:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\mathbf{L}\mathbf{V} = \Delta,
|
|
\label{eq:lap_system}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where:
|
|
- \f$\mathbf{L}\f$ is a \f$n \times n\f$ sparse matrix, referred to as the <em>Laplacian matrix</em>. Its elements \f$ m_{ij} \f$, \f$i,j \in \{1 \dots n\} \f$ are defined as follows:
|
|
- \f$ m_{ii} \f$ = \f$ \sum_{\mathbf{v}_j \in N(\mathbf{v}_i)}w_{ij} \f$,
|
|
- \f$ m_{ij} = -w_{ij} \f$ if \f$ \mathbf{v}_j \in N(\mathbf{v_i}) \f$,
|
|
- \f$ 0 \f$ otherwise.
|
|
- \f$\mathbf{V}\f$ is a \f$n \times 3\f$ matrix made of the %Cartesian coordinates of the vertices.
|
|
|
|
|
|
\subsubsection Surface_Modeling_Overview_Laplacian_Deformation Laplacian Deformation
|
|
|
|
This section is an introduction to provide the background for the next two sub-sections describing the algorithms implemented
|
|
in this package. A system relying only on the approach described below results in non-smooth transitions in the neighborhood of
|
|
the control vertices. For a survey on different Laplacian-based editing techniques we refer to
|
|
\cgalCite{Botsch2008OnLinearVariational}.
|
|
|
|
The main idea behind Laplacian-based deformation techniques is to preserve the Laplacian representation
|
|
under deformation constraints. The Laplacian representation of a surface mesh is treated as a representative form of the discretized surface,
|
|
and the deformation process must follow the deformation constraints while preserving the Laplacian representation as much as possible.
|
|
|
|
There are different ways to incorporate deformation constraints into the deformation system \cgalCite{Botsch2008OnLinearVariational}.
|
|
This package supports hard constraints, that is, target positions of control vertices are preserved after the deformation.
|
|
|
|
Given a surface mesh deformation system with a ROI made of \f$ n \f$ vertices and \f$ k \f$ control vertices, we consider the following linear system:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\left[
|
|
\begin{array}{ccc}
|
|
\mathbf{L}_f\\
|
|
0 \; \mathbf{I}_c
|
|
\end{array}
|
|
\right]
|
|
\mathbf{V} =
|
|
\left[
|
|
\begin{array}{ccc}
|
|
{\Delta}_f \\
|
|
\mathbf{V}_c
|
|
\end{array}
|
|
\right],
|
|
\label{eq:lap_energy_system}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where:
|
|
- \f$\mathbf{V}\f$ is a \f$n \times 3\f$ matrix denoting the unknowns of the system that represent the vertex coordinates after deformation. The system is built so that the \f$ k \f$ last rows correspond to the control vertices.
|
|
- \f$\mathbf{L}_f\f$ denotes the Laplacian matrix of the unconstrained vertices. It is a \f$ (n-k) \times n \f$ matrix as defined in Eq. \f$\eqref{eq:lap_system}\f$ but removing the rows corresponding to the control vertices.
|
|
- \f$\mathbf{I}_c\f$ is the \f$k \times k\f$ identity matrix.
|
|
- \f${\Delta}_f\f$ denotes the Laplacian representation of the unconstrained vertices as defined in Eq. \f$\eqref{eq:lap_system}\f$ but removing the rows corresponding to the control vertices.
|
|
- \f$\mathbf{V}_c\f$ is a \f$k \times 3\f$ matrix containing the %Cartesian coordinates of the target positions of the control vertices.
|
|
|
|
The left-hand side matrix of the system of Eq.\f$\eqref{eq:lap_energy_system}\f$ is a square non-symmetric sparse matrix. To
|
|
solve the aforementioned system, an appropriate solver (e.g. LU solver) needs to be used. Note that solving this system
|
|
preserves the Laplacian representation of the surface mesh restricted to the unconstrained vertices while satisfying the deformation constraints.
|
|
|
|
\subsection Surface_Modeling_Overview_ARAP As-Rigid-As Possible Deformation
|
|
|
|
Given a surface mesh \f$M\f$ with \f$ n \f$ vertices \f$ \{\mathbf{v}_i\} i \in \{1 \dots n \} \f$ and some deformation
|
|
constraints, we consider the following energy function:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\sum_{\mathbf{v}_i \in M}
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} w_{ij}
|
|
\left\| (\mathbf{v}'_i - \mathbf{v}'_j) - \mathbf{R}_i(\mathbf{v}_i - \mathbf{v}_j) \right\|^2,
|
|
\label{eq:arap_energy}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where:
|
|
- \f$\mathbf{R}_i\f$ is a \f$ 3 \times 3 \f$ rotation matrix
|
|
- \f$w_{ij}\f$ denotes a weight
|
|
- \f$N(\mathbf{v}_i)\f$ denotes the set of vertices adjacent to \f$\mathbf{v}_i\f$ in \f$M\f$
|
|
- \f$N(\mathbf{v}'_i)\f$ denotes a new position of the vertex \f$N(\mathbf{v}_i)\f$ after a given deformation
|
|
|
|
An as-rigid-as possible surface mesh deformation \cgalCite{Sorkine2007AsRigidAs} is defined by minimizing
|
|
this energy function under the deformation constraints, i.e. the assigned position
|
|
\f$ {v}'_i\f$ for each vertex \f$ \mathbf{v}_i\f$ in the set of control vertices.
|
|
Defining the <i>one-ring neighborhood</i> of a vertex as its set of adjacent vertices,
|
|
the intuitive idea behind this energy function is to allow each one-ring neighborhood of
|
|
vertices to have an individual rotation, and at the same time to prevent shearing by taking advantage of the
|
|
overlapping of one-ring neighborhoods of adjacent vertices (see \cgalFigureRef{Overlapping_cells}).
|
|
|
|
|
|
\cgalFigureBegin{Overlapping_cells,overlapping_cells.png}
|
|
Overlaps of one-ring neighborhoods of vertices.
|
|
The one-ring neighborhoods of four vertices
|
|
are drawn with different colors, the corresponding vertex is colored accordingly.
|
|
\cgalFigureEnd
|
|
|
|
There are two unknowns per vertex in Eq. \f$\eqref{eq:arap_energy}\f$: the new positions (\f$\mathbf{v}'_k\f$) of the unconstrained
|
|
vertices and the rotation matrices (\f$\mathbf{R}_i\f$).
|
|
If the energy contribution of each vertex is positive, this boils down to minimizing the energy contribution of each vertex \f$\mathbf{v}_i\f$.
|
|
|
|
Each such term of the energy is minimized by using a two-step optimization approach (also called local-global approach).
|
|
|
|
In the first step, the positions of the vertices are considered as fixed so that the rotation matrices are the only unknowns.
|
|
|
|
For the vertex \f$\mathbf{v}_i\f$, we consider the covariance matrix \f$\mathbf{S}_i\f$:
|
|
\f[
|
|
\begin{equation}
|
|
\mathbf{S}_i = \sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} w_{ij} (\mathbf{v}_i - \mathbf{v}_j)(\mathbf{v}'_i - \mathbf{v}'_j)^T,
|
|
\label{eq:cov_matrix}
|
|
\end{equation}
|
|
\f]
|
|
|
|
It was shown \cgalCite{Sorkine2009LeastSquaresRigid} that minimizing the energy contribution of
|
|
\f$\mathbf{v}_i\f$ in Eq. \f$\eqref{eq:arap_energy}\f$ is equivalent to maximizing the trace of the matrix
|
|
\f$\mathbf{R}_i \mathbf{S}_i\f$. \f$\mathbf{R}_i \f$ is the transpose of the unitary matrix in the polar decomposition
|
|
of \f$\mathbf{S}_i\f$.
|
|
|
|
|
|
In the second step, the rotation matrices are substituted into the partial derivative of Eq.\f$\eqref{eq:arap_energy}\f$
|
|
with respect to \f$\mathbf{v}'_i\f$. Assuming the weights are symmetric, setting the derivative to zero results in the
|
|
following equation:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} w_{ij}(\mathbf{v}'_i - \mathbf{v}'_j) =
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} w_{ij} \frac{(\mathbf{R}_i + \mathbf{R}_j)}{2} (\mathbf{v}_i - \mathbf{v}_j).
|
|
\label{eq:lap_ber}
|
|
\end{equation}
|
|
\f]
|
|
|
|
The left-hand side of this equation corresponds to the one of Eq.\f$\eqref{eq:lap_open}\f$,
|
|
and we can set \f$\Delta\f$ to be the right-hand side.
|
|
Solving the linear system in Eq. \f$\eqref{eq:lap_energy_system}\f$ gives the new positions of the unconstrained vertices.
|
|
|
|
This two-step optimization can be applied several times iteratively to obtain a better result.
|
|
|
|
\note The matrix built with the Laplacian matrix of the unconstrained vertices in the left-hand side of
|
|
Eq. \f$\eqref{eq:lap_energy_system}\f$ depends only on the initial surface mesh structure and on which vertices
|
|
are control vertices. Once the control vertices are set, we can use a direct solver to factorize
|
|
the sparse matrix in Eq. \f$\eqref{eq:lap_energy_system}\f$, and reuse this factorization during
|
|
each iteration of the optimization procedure.
|
|
|
|
The original algorithm \cgalCite{Sorkine2007AsRigidAs} we described assumes that:
|
|
|
|
- the weight between two vertices is symmetric. In order to support asymmetric weights in our implementation,
|
|
we slightly change Eq. \f$\eqref{eq:lap_ber}\f$ to:
|
|
\f[
|
|
\begin{equation}
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} (w_{ij} + w_{ji})(\mathbf{v}'_i - \mathbf{v}'_j) =
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} (w_{ij}\mathbf{R}_i + w_{ji}\mathbf{R}_j)(\mathbf{v}_i - \mathbf{v}_j).
|
|
\label{eq:lap_ber_asym}
|
|
\end{equation}
|
|
\f]
|
|
|
|
- The energy contribution of each vertex is positive. If the weight between two vertices is always positive, this is always the case.
|
|
However, when using the cotangent weighting scheme (the default in our implementation), if the sum of the angles opposite to an edge is greater than \f$ \pi \f$,
|
|
its cotangent weight is negative. As a workaround for bad quality meshes, we eliminate those negative weights
|
|
by setting them to zero.
|
|
|
|
A method minimizing another energy function is described next to avoid the latter issue.
|
|
|
|
\subsection Surface_Modeling_Overview_ARAP_Rims Spokes and Rims Version
|
|
|
|
The elastic energy function proposed by \cgalCite{Chao2010SimpleGeomModel} additionally takes into account
|
|
all the opposite edges in the facets incident to a vertex. The energy function to minimize becomes:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\sum_{\mathbf{v}_i \in M}
|
|
\sum_{(\mathbf{v}_j, \mathbf{v}_k) \in E(\mathbf{v}_i)} w_{jk}
|
|
\left\| (\mathbf{v}'_j - \mathbf{v}'_k) - \mathbf{R}_i(\mathbf{v}_j - \mathbf{v}_k) \right\|^2,
|
|
\label{eq:arap_energy_rims}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where \f$E(\mathbf{v}_i)\f$ consists of the set of edges incident to \f$\mathbf{v}_i\f$ (the <em>spokes</em>) and
|
|
the set of edges in the link (the <em>rims</em>) of \f$\mathbf{v}_i\f$ in the surface mesh \f$M\f$
|
|
(see \cgalFigureRef{Spoke_and_rim_edges}).
|
|
|
|
\cgalFigureBegin{Spoke_and_rim_edges, spoke_and_rim_edges_2.png}
|
|
The vertices \f$ \mathbf{v}_n\f$ and \f$ \mathbf{v}_m\f$ are the opposite vertices to the edge
|
|
\f$ \mathbf{v}_i \mathbf{v}_j\f$.
|
|
\cgalFigureEnd
|
|
|
|
The method to get the new positions of the unconstrained vertices is similar to the two-step optimization
|
|
method explained in \ref Surface_Modeling_Overview_ARAP.
|
|
For the first step, the Eq. \f$\eqref{eq:cov_matrix}\f$ is modified to take into account the edges in \f$E(\mathbf{v}_i)\f$:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\mathbf{S}_i = \sum_{(\mathbf{v}_j, \mathbf{v}_k) \in E(\mathbf{v}_i)} w_{jk} (\mathbf{v}_j - \mathbf{v}_k)(\mathbf{v}'_j - \mathbf{v}'_k)^T,
|
|
\label{eq:cov_matrix_sr}
|
|
\end{equation}
|
|
\f]
|
|
|
|
For the second step, setting partial derivative of Eq. \f$\eqref{eq:arap_energy_rims}\f$ to zero with
|
|
respect to \f$\mathbf{v}_i\f$ gives the following equation:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} (w_{ij} + w_{ji})(\mathbf{v}'_i - \mathbf{v}'_j) =
|
|
\sum_{\mathbf{v}_j \in N(\mathbf{v}_i)} \frac{w_{ij}(\mathbf{R}_i + \mathbf{R}_j + \mathbf{R}_m) + w_{ji}(\mathbf{R}_i + \mathbf{R}_j + \mathbf{R}_n)}{3} (\mathbf{v}_i - \mathbf{v}_j).
|
|
\label{eq:lap_ber_rims}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where \f$\mathbf{R}_m\f$ and \f$\mathbf{R}_n\f$ are the rotation matrices of the vertices \f$\mathbf{v}_m\f$,
|
|
\f$\mathbf{v}_n\f$ which are the opposite vertices of the edge \f$\mathbf{v}_i \mathbf{v}_j\f$
|
|
(see \cgalFigureRef{Spoke_and_rim_edges}). Note that if the edge \f$ \mathbf{v}_i \mathbf{v}_j \f$ is on
|
|
the boundary of the surface mesh, then \f$ w_{ij} \f$ must be 0 and \f$ \mathbf{v}_m
|
|
\f$ does not exist.
|
|
|
|
An important property of this approach compared to \ref Surface_Modeling_Overview_ARAP is that the contribution to the global energy
|
|
of each vertex is guaranteed to be non-negative when using the cotangent weights \cgalCite{Chao2010SimpleGeomModel}.
|
|
Thus even with negative weights, the minimization of the energy with the iterative method presented is always guaranteed.
|
|
However, this method is more dependent on the discretization of the deformed surface (See \cgalFigureRef{Arap_spokes_comparison}).
|
|
|
|
The implementation in this package uses the cotangent weights by default (negative values included) as proposed in \cgalCite{Chao2010SimpleGeomModel}.
|
|
|
|
\section Surface_Modeling_History Design and Implementation History
|
|
An initial version of this package has been implemented during the 2011 Google Summer of Code by Yin Xu under the guidance of Olga Sorkine and Andreas Fabri.
|
|
Ilker O. Yaz took over the finalization of the package with the help of Sébastien Loriot for the documentation and the API.
|
|
The authors are grateful to Gaël Guennebaud for his great help on using the Eigen library and for providing the code to compute
|
|
the closest rotation.
|
|
*/
|
|
|
|
} /* namespace CGAL */
|
|
|