mirror of https://github.com/CGAL/cgal
394 lines
22 KiB
Plaintext
394 lines
22 KiB
Plaintext
namespace CGAL {
|
|
/*!
|
|
|
|
\mainpage User Manual
|
|
\anchor Chapter_SurfaceModeling
|
|
|
|
\cgalAutoToc
|
|
\authors Olga Sorkine, Yin Xu and Ilker %O. Yaz
|
|
|
|
\section Surface_Modeling_Intro Introduction
|
|
Surface mesh deformation is the process of transforming the shape of the surface mesh under control constraints without requiring any additional structure other than the surface itself.
|
|
|
|
In this package, we implement the algorithm described in \cite Sorkine2007AsRigidAs,
|
|
where a nonlinear deformation energy is minimized under constraints to preserve rigidity as much as possible.
|
|
|
|
\cgalFigureBegin{Main_image_suggestion, main_image_suggestion_4_resized.png}
|
|
Deformation results for dinosaur model.
|
|
\cgalFigureEnd
|
|
|
|
\section Surface_Modeling_Overview Overview
|
|
Since there is a close relationship between differential representations and as-rigid-as possible deformation, we first
|
|
give a brief overview of related topics.
|
|
|
|
\subsection Surface_Modeling_Overview_Laplacian Laplacian Representation
|
|
Laplacian coordinates are a way of representing the mesh differentially, using the local coordinate of a vertex
|
|
with respect to its one-ring neighbors (i.e. cell) rather than its global position. The Laplacian coordinate of a vertex \f$ \mathbf{v}_i \f$ is:
|
|
|
|
\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$ is the one-ring neighbor vertices of \f$\mathbf{v}_i\f$;
|
|
-\f$w_{ij}\f$ is the weight between \f$\mathbf{v}_i\f$ and the neighbor \f$\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 vertex
|
|
\f$\mathbf{v}_j\f$ in the one-ring neighborhood. For the computation of the Laplacian,
|
|
this actually corresponds to subtracting the average of neighbor vertices from the current vertex (\cgalFigureRef{Simple_laplacian}).
|
|
|
|
We use cotangent weights \cite Pinkall1993Cotangent to eliminate the effect of non-uniform discretization of the mesh model.
|
|
|
|
\cgalFigureBegin{Simple_laplacian,simple_mesh_with_laplacian.png}
|
|
Laplacian with uniform weights: the red square vertex is the average position of the neighbor vertices of the yellow circle vertex. Blue vector represents the Laplacian coordinate of the yellow vertex under uniform weights.
|
|
\cgalFigureEnd
|
|
|
|
Considering the mesh as a whole, it is possible to generate a linear system which results in Laplacian coordinates of all vertices:
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\mathbf{L}\mathbf{V} = \delta,
|
|
\label{eq:lap_system}
|
|
\end{equation}
|
|
\f]
|
|
|
|
where
|
|
- \f$n\f$ is the number of vertices in the mesh,
|
|
- \f$\mathbf{L}\f$ is an \f$n \times n\f$ matrix, called the <i>Laplacian matrix</i>,
|
|
in which each row is filled according to Eq.\f$\eqref{eq:lap_open}\f$, that is the element \f$ m_{ij} \f$ of \f$ L \f$ is equal to \f$ -w_{ij} \f$
|
|
if \f$ \mathbf{v}_j \in N(\mathbf{v}) \f$, \f$ m_{ii} \f$ is equal to \f$ \sum_{\mathbf{v}_j \in N(\mathbf{v}_i)}w_{ij} \f$ and \f$ 0 \f$ elsewhere,
|
|
- \f$\mathbf{V}\f$ is an \f$n \times 1\f$ vector consisting of
|
|
global positions of vertices; and \f$\delta\f$ represents Laplacian coordinates of the mesh.
|
|
|
|
|
|
\subsection Surface_Modeling_Overview_Laplacian_Deformation Naive Laplacian Deformation
|
|
A simple mesh deformation system consists of a mesh, a few selected control vertices (i.e. <i>handles</i>)
|
|
and target positions of the handles as deformation constraints.
|
|
|
|
The main idea behind Laplacian deformation is to preserve the Laplacian coordinates under deformation constraints.
|
|
In other words, Laplacian coordinates are treated as a representative form for the shape, and the deformation
|
|
process must follow the deformation constraints while preserving the Laplacian coordinates as much as possible.
|
|
|
|
There are different ways to incorporate deformation constraints into the system \cite Botsch2008OnLinearVariational.
|
|
We provide support for hard constraints, which means constrained positions of handles are preserved as is after deformation.
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\left[
|
|
\begin{array}{ccc}
|
|
\mathbf{L}_f \\
|
|
\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{L}_f\f$ is Laplacian matrix of free (i.e. not constrained) vertices, \f$\mathbf{I}_c\f$ contains a row for each constrained vertex
|
|
with one in corresponding position of the vertex and zeros elsewhere, \f${\delta}_f\f$ is Laplacian coordinates of free vertices,
|
|
and \f$\mathbf{V}_c\f$ holds the target positions of constrained vertices. Solving the system under deformation constraints
|
|
preserves Laplacians of free vertices (upper side of the system) while satisfying the constraints (lower side of the system).
|
|
|
|
The constructed matrix in Eq.\f$\eqref{eq:lap_energy_system}\f$ is square nonsymmetric sparse matrix and an appropriate solver should be used.
|
|
|
|
\subsection Surface_Modeling_Overview_ARAP As-Rigid-As Possible Deformation
|
|
|
|
The energy function proposed to be minimized in \cite Sorkine2007AsRigidAs under modeling constraints for a mesh structure \f$M\f$ is:
|
|
|
|
\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 rotation matrix for \f$\mathbf{v}_i\f$, and \f$\mathbf{v}'_i\f$ is deformed version of \f$\mathbf{v}_i\f$.
|
|
|
|
The intuitive idea behind this energy function is that allowing cells to
|
|
have individual rotations, and at the same time using overlapping nature of the
|
|
cells to prevent sharing \cgalFigureRef{Overlapping_cells}.
|
|
|
|
\cgalFigureBegin{Overlapping_cells,overlapping_cells.png}
|
|
Illustration of overlapping cells
|
|
\cgalFigureEnd
|
|
|
|
Minimization of the energy function in Eq.\f$\eqref{eq:arap_energy}\f$ for every vertex under deformation constraints provides a
|
|
as-rigid-as possible deformation. In the equation, there are two unknowns: the new positions of vertices (\f$\mathbf{v}'_i\f$)
|
|
and the rotation matrices (\f$\mathbf{R}_i\f$). The problem is solved 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 \f$\mathbf{R}_i\f$ is the only unknown.
|
|
Given the covariance matrix \f$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]
|
|
|
|
\f$\mathbf{R}_i\f$ matrices minimizing Eq.\f$\eqref{eq:arap_energy}\f$ can be found using a singular value decomposition (SVD ) of \f$\mathbf{S}_i\f$
|
|
and by also guaranteeing \f$det(\mathbf{R}_i)\f$ to be positive by changing the sign of the eigen vector in \f$\mathbf{U}_i\f$ which corresponds to smallest eigen value \cite Sorkine2009LeastSquaresRigid.
|
|
|
|
\f[
|
|
\begin{equation}
|
|
\mathbf{R}_i = \mathbf{V}_i \mathbf{U}_i^T \:where\: \mathbf{S}_i = \mathbf{U}_i\Sigma_i \mathbf{V}_i^T.
|
|
\label{eq:optimal_rot}
|
|
\end{equation}
|
|
\f]
|
|
|
|
In the second step, the computed rotation matrices are replaced into the partial derivative of Eq.\f$\eqref{eq:arap_energy}\f$
|
|
with respect to \f$\mathbf{v}'_i\f$. In order to find \f$\mathbf{v}'_i\f$
|
|
which minimizes the deformation energy, 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 Eq.\f$\eqref{eq:lap_open}\f$,
|
|
so that we can fill \f$\delta\f$ using the right-hand side of the equation which consists only of known values in our system.
|
|
Then, the linear system in Eq. \f$\eqref{eq:lap_energy_system}\f$ can be solved in order to find \f$\mathbf{V}'\f$.
|
|
This two-step optimization can be applied iteratively for several times to obtain a better result.
|
|
|
|
An important observation is the linear system in Eq. \f$\eqref{eq:lap_energy_system}\f$ only depends on the initial mesh structure and which vertices are handle vertices.
|
|
So that, once handle vertices are determined, we can use a direct solver to factorize sparse matrix in Eq. \f$\eqref{eq:lap_energy_system}\f$, and use this factorization
|
|
in each iteration of the optimization procedure with the cost of back-substitution.
|
|
|
|
\note The original algorithm in \cite Sorkine2007AsRigidAs assumes that weights are symmetric (i.e. \f$w_{ij} = w_{ji}\f$) while deriving the Eq. \f$\eqref{eq:lap_ber}\f$.
|
|
In our implementation, in order to support asymmetric weights, we slightly changed the 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]
|
|
|
|
A small problem with this approach is that when meshing is bad (in a way of leading many negative edge weights), the energy function might become negative and thus optimization fails.
|
|
As a workaround, we eliminate negative values by setting them to zero in our default weights.
|
|
|
|
\subsection Surface_Modeling_Overview_ARAP_Rims Spokes and Rims Version
|
|
|
|
\cgalFigureBegin{Spoke_and_rim_edges, spoke_and_rim_edges_2.png}
|
|
Considered edges \f$E(\mathbf{v}_i)\f$ in energy function for a vertex \f$\mathbf{v}_i\f$.
|
|
\cgalFigureEnd
|
|
|
|
The modified energy function proposed by \cite Chao2010SimpleGeomModel take into account all edges around the facets incident to a vertex. The energy function 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$ is set of edges around facets incident to \f$\mathbf{v}_i\f$ (see \cgalFigureRef{Spoke_and_rim_edges}).
|
|
While constructing the solution, we can follow the same approach (i.e. two-step optimization) explained in \ref Surface_Modeling_Overview_ARAP.
|
|
For the first step, corresponding Eq. \f$\eqref{eq:cov_matrix}\f$ is simply modified to include all edges in \f$E(\mathbf{v}_i)\f$.
|
|
For the second step, setting partial derivative of Eq. \f$\eqref{eq:arap_energy_rims}\f$ zero with respect to \f$\mathbf{v}_i\f$ results in:
|
|
|
|
\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 rotation matrices of \f$\mathbf{v}_m\f$,
|
|
\f$\mathbf{v}_n\f$ which are opposite vertices of edge \f$\mathbf{v}_i, \mathbf{v}_j\f$.
|
|
|
|
The main difference compared to \ref Surface_Modeling_Overview_ARAP is that the energy function is guaranteed to stay nonnegative.
|
|
We are using cotangent weights as it is proposed in \cite Chao2010SimpleGeomModel.
|
|
However, produced results can be asymmetric as demonstrated in \cgalFigureRef{Arap_spokes_comparison}.
|
|
|
|
\cgalFigureBegin{Arap_spokes_comparison, arap_spokes_comparison.png}
|
|
Comparison between \ref Surface_Modeling_Overview_ARAP and \ref Surface_Modeling_Overview_ARAP_Rims. \ref Surface_Modeling_Overview_ARAP_Rims is affected from meshing bias and produces asymmetric results.
|
|
\cgalFigureEnd
|
|
|
|
\section Surface_Modeling_API API
|
|
|
|
\cgalFigureBegin{Flow_of_API, usage_flow_of_API.png}
|
|
The deformation work flow.
|
|
\cgalFigureEnd
|
|
|
|
\todo <b>SL</b>: what is the meaning of the different box style? what software did you use? are the source inside the git repo? Not sure that it is obvious that <i>clear</i> is optional.
|
|
|
|
\subsection Preprocess_section Preprocess Section
|
|
|
|
We present the API by considering the user point of view during a deformation session (see \cgalFigureRef{Flow_of_API} for an illustration).
|
|
|
|
The deformation is handled thanks to a function object that must be initialized:
|
|
- `Deform_mesh deform = Deform_mesh(Polyhedron& polyhedron, VertexIndexMap vertex_index_map, EdgeIndexMap edge_index_map, unsigned int iterations = 5, double tolerance = 1e-4, WeightCalculator weight_calculator = WeightCalculator())`
|
|
- `vertex_index_map` and `edge_index_map` are property maps used for assigning indexes to vertices and edges.
|
|
- For details on `iterations` and `tolerance` parameters see the \ref deform_method "deform method".
|
|
- For details on `weight_calculator` see \ref Example_4
|
|
- `deformation_type` can be used as nontype template parameter for selecting deformation algorithm. The default algorithm is \ref Surface_Modeling_Overview_ARAP_Rims.
|
|
|
|
The region the algorithm operates on (called the region-of-interest and abbreviated ROI) must be defined.
|
|
The API is convenient enough for either using a subset or the whole mesh as region of interest.
|
|
- `Deform_mesh::insert_roi(vertex_descriptor vd)` insert vertices one by one
|
|
- `Deform_mesh::insert_roi(InputIterator begin, InputIterator end)` insert vertices in the range
|
|
|
|
The ROI should also include the handle vertices \cgalFigureRef{Simple_roi}.
|
|
Note that the algorithm automatically finds the exterior-boundary of the ROI and treat them as <i>fixed handles</i>.
|
|
|
|
\cgalFigureBegin{Simple_roi, roi_example.png}
|
|
Red vertices represent handles, and yellow vertices together with red vertices represent ROI.
|
|
\cgalFigureEnd
|
|
|
|
Similarly, the handle vertices which will constrain the deformation should be added.
|
|
At this point vertices should be grouped that are expected to translate or rotate together.
|
|
- `Deform_mesh::create_handle_group()` create an empty group for vertices and return a representative `Deform_mesh::Handle_group` of the group
|
|
- `Deform_mesh::insert_handle(Handle_group handle_group, vertex_descriptor vd)` insert a vertex in a group of handles
|
|
- `Deform_mesh::insert_handle(Handle_group handle_group, InputIterator begin, InputIterator end)` insert a set of vertices to a group of handles
|
|
|
|
Vertices in the same group can be connected or disconnected. The reason why vertices are gathered in groups is to be able to apply the same translation and rotation to the vertices of a group.
|
|
|
|
The status of a vertex (roi or handle) can be reinitialized :
|
|
- `Deform_mesh::erase_roi(vertex_descriptor vd)` erase the vertex from ROI
|
|
- `Deform_mesh::erase_handle(Handle_group handle_group)` erase the whole group
|
|
- `Deform_mesh::erase_handle(Handle_group handle_group, vertex_descriptor vd)` erase the vertex which is inside given group
|
|
- `Deform_mesh::erase_handle(vertex_descriptor vd)` erase the vertex by searching through all handle groups
|
|
|
|
Iteration over ROI and handle vertices are possible using:
|
|
- Deform_mesh::roi_vertices() iterate over ROI vertices
|
|
- Deform_mesh::handle_groups() iterate over the groups of handles
|
|
- Deform_mesh::handles(Handle_group handle_group) iterate over the vertices inside a given group of handles
|
|
|
|
Also constant time querying the state of a vertex can be done using:
|
|
- Deform_mesh::is_roi(vertex_descriptor vd) const check whether the vertex is inside ROI
|
|
- Deform_mesh::is_handle(vertex_descriptor vd) const check whether a vertex is inside any group of handles
|
|
|
|
After having set the ROI and the handle groups, a preprocess function needs to be called.
|
|
- `Deform_mesh::preprocess()` preprocess ROI and handles, make internal stuff ready for deformation
|
|
|
|
Preprocess function returns true if factorization of the linear system in Eq. \f$\eqref{eq:lap_energy_system}\f$ is successful, false in case of a failure.
|
|
|
|
Rank deficiency is the main reason for the failure, which happens if there is no path between a free vertex and a handle vertex (i.e. both fixed and user-inserted).
|
|
Some cases for such failure:
|
|
- Inserting whole mesh as ROI and inserting no handle vertices,
|
|
- Using custom weights which includes too many zero weights, so that the connectivity is broken.
|
|
|
|
Calling preprocess function is optional, which means it is called by functions in \ref Deform_section, if there is any need to.
|
|
Of course, it might create a delay in corresponding function caused by preprocess computations.
|
|
|
|
\note The ROI does not have to be a connected component, however considering performance issues it is best to use an individual deformation object per each connected component ROI.
|
|
The main performance gain in such case is related to factorization process (i.e. `Deform_mesh::preprocess`).
|
|
In deformation, if we assume that fill-ins produced by factorization is negligible, and cost of back-substitution is linearly proportional to non-zero elements in the system,
|
|
then the performance hit should not be remarkable.
|
|
|
|
|
|
\subsection Deform_section Deform Section
|
|
|
|
In order to deform the mesh, the final position of handle vertices must be provided. The first method provided is
|
|
- `Deform_mesh::assign(vertex_descriptor vd, const Point& target_position)` that directly assigns the target position for the handle vertex
|
|
|
|
If the target positions for the handles are provided directly, then handles can be simply gathered into a single group
|
|
and `Deform_mesh::assign` can be used for direct assignment of target positions.
|
|
|
|
\anchor rot
|
|
The second possibility offered by the API is to set translation and rotation for vertex handle groups.
|
|
- `Deform_mesh::translate(Const_handle_group handle_group, const Vector& translation)` translate the vertices in `handle_group` by `translation`
|
|
- `Deform_mesh::rotate(Const_handle_group handle_group, const Point& rotation_center, const Quaternion& quat, const Vect& translation)`
|
|
rotate the `handle_group` around `rotation_center` using the quaternion `quat` and then translate it by `translation`
|
|
|
|
`rotate` and `translate` operate on <i>original positions</i> (i.e. positions of vertices at the time of construction) of handle vertices.
|
|
So that, they are not accumulative (i.e. any call to `rotate`, `translate`, or `assign` overrides previous calls).
|
|
|
|
Note that changing target positions of handle vertices has no immediate effect on actual positions.
|
|
Also unconstrained handles will be constrained to their last target positions.
|
|
|
|
\anchor deform_method
|
|
The final action is invoking the deformation algorithm:
|
|
- `Deform_mesh::deform(unsigned int iterations, double tolerance)` deform mesh with one-time parameters.
|
|
- Each iteration consists two steps of optimization procedure explained in \ref Surface_Modeling_Overview_ARAP
|
|
- In case of tolerance is smaller than or equal to zero, no energy based termination occurs.
|
|
Otherwise it is compared against \f$|energy(m_i) - energy(m_{i-1})| / energy(m_i)\f$ where \f$m_i\f$ is the resulting surface mesh at iteration \f$i\f$.
|
|
- `Deform_mesh::deform()` deform the mesh using global iterations and tolerance provided in constructor or using `Deform_mesh::set_iterations`, and `Deform_mesh::set_tolerance`
|
|
|
|
After deformation, handle vertices will be placed to their target positions and free vertices will be deformed accordingly.
|
|
The steps in \ref Deform_section can be repeated for further deformations.
|
|
|
|
\warning The transformation provided to the functions `Deform_mesh::rotate` and `Deform_mesh::translate` always applies on the original coordinates of the input graph.
|
|
That is in an interactive deformation session, the user is responsible for composing the transformations given to `Deform_mesh::rotate` and `Deform_mesh::translate`.
|
|
In particular, any vertex which is no longer inside the ROI will be assigned to its original position in `Deform_mesh::preprocess`.
|
|
The original positions can be updated by calling `Deform_mesh::overwrite_original_positions`. This behavior is illustrated in \ref Surface_Modeling_Demo_Tutorial.
|
|
|
|
\subsection Surface_modeling_examples Examples
|
|
|
|
\subsubsection Example_1 Using the Whole Mesh as Region-of-interest
|
|
In this example, the whole mesh is used as ROI and a few vertices are added as handles. `Deform_mesh::assign` is used for constraining handle vertices.
|
|
|
|
\cgalExample{Surface_modeling/all_roi_assign_example.cpp}
|
|
|
|
\cgalFigureBegin{Example_1_results, example_1_results.png}
|
|
Deformation results from example, deform_1.off and deform_2.off.
|
|
\cgalFigureEnd
|
|
|
|
\subsubsection Example_2 Using Different Groups of Handles
|
|
In this example, we use handle grouping feature of the API together with `Deform_mesh::translate` and `Deform_mesh::rotate`.
|
|
|
|
\cgalExample{Surface_modeling/k_ring_roi_translate_rotate_example.cpp}
|
|
|
|
\cgalFigureBegin{Example_2_results, example_2_results.png}
|
|
Deformation results from example, deform_1.off and deform_2.off.
|
|
\cgalFigureEnd
|
|
|
|
\subsubsection Example_3 Using an Enriched Polyhedron with Ids
|
|
In this example, we use an enriched polyhedron whose edges and vertices support an id field to store indexes.
|
|
The main advantage is to decrease the complexity of accessing associated indexes with edges and vertices from logarithmic to constant.
|
|
|
|
\cgalExample{Surface_modeling/enriched_polyhedron_with_custom_pmap_example.cpp}
|
|
|
|
\subsubsection Example_4 Using a Custom Edge Weighting Scheme
|
|
Using custom weights for edges is also possible via providing a model of SurfaceModelingWeightCalculator.
|
|
In this example, we assume having precomputed weights for each edge and use an internal map for store them in custom weight calculator.
|
|
|
|
You can also check the code example in SurfaceModelingWeightCalculator for direct computation.
|
|
|
|
\cgalExample{Surface_modeling/custom_weight_for_edges_example.cpp}
|
|
|
|
\section Surface_Modeling_Demo How to use the demo
|
|
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/288N_SKrubY" 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 Handle / ROI and Overwrite
|
|
|
|
@htmlonly
|
|
<center>
|
|
<iframe width="420" height="315" src="http://www.youtube.com/embed/PWXmViSzD8M" frameborder="0" allowfullscreen></iframe>
|
|
</center>
|
|
@endhtmlonly
|
|
|
|
\section Surface_Modeling_History Implementation History
|
|
The 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.
|
|
*/
|
|
|
|
} /* namespace CGAL */
|
|
|