$\newcommand{\dotp}[2]{\langle #1, #2 \rangle}$ $\newcommand{\enscond}[2]{\lbrace #1, #2 \rbrace}$ $\newcommand{\pd}[2]{ \frac{ \partial #1}{\partial #2} }$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\umax}[1]{\underset{#1}{\max}\;}$ $\newcommand{\umin}[1]{\underset{#1}{\min}\;}$ $\newcommand{\uargmin}[1]{\underset{#1}{argmin}\;}$ $\newcommand{\norm}[1]{\|#1\|}$ $\newcommand{\abs}[1]{\left|#1\right|}$ $\newcommand{\choice}[1]{ \left\{ \begin{array}{l} #1 \end{array} \right. }$ $\newcommand{\pa}[1]{\left(#1\right)}$ $\newcommand{\diag}[1]{{diag}\left( #1 \right)}$ $\newcommand{\qandq}{\quad\text{and}\quad}$ $\newcommand{\qwhereq}{\quad\text{where}\quad}$ $\newcommand{\qifq}{ \quad \text{if} \quad }$ $\newcommand{\qarrq}{ \quad \Longrightarrow \quad }$ $\newcommand{\ZZ}{\mathbb{Z}}$ $\newcommand{\CC}{\mathbb{C}}$ $\newcommand{\RR}{\mathbb{R}}$ $\newcommand{\EE}{\mathbb{E}}$ $\newcommand{\Zz}{\mathcal{Z}}$ $\newcommand{\Ww}{\mathcal{W}}$ $\newcommand{\Vv}{\mathcal{V}}$ $\newcommand{\Nn}{\mathcal{N}}$ $\newcommand{\NN}{\mathcal{N}}$ $\newcommand{\Hh}{\mathcal{H}}$ $\newcommand{\Bb}{\mathcal{B}}$ $\newcommand{\Ee}{\mathcal{E}}$ $\newcommand{\Cc}{\mathcal{C}}$ $\newcommand{\Gg}{\mathcal{G}}$ $\newcommand{\Ss}{\mathcal{S}}$ $\newcommand{\Pp}{\mathcal{P}}$ $\newcommand{\Ff}{\mathcal{F}}$ $\newcommand{\Xx}{\mathcal{X}}$ $\newcommand{\Mm}{\mathcal{M}}$ $\newcommand{\Ii}{\mathcal{I}}$ $\newcommand{\Dd}{\mathcal{D}}$ $\newcommand{\Ll}{\mathcal{L}}$ $\newcommand{\Tt}{\mathcal{T}}$ $\newcommand{\si}{\sigma}$ $\newcommand{\al}{\alpha}$ $\newcommand{\la}{\lambda}$ $\newcommand{\ga}{\gamma}$ $\newcommand{\Ga}{\Gamma}$ $\newcommand{\La}{\Lambda}$ $\newcommand{\si}{\sigma}$ $\newcommand{\Si}{\Sigma}$ $\newcommand{\be}{\beta}$ $\newcommand{\de}{\delta}$ $\newcommand{\De}{\Delta}$ $\newcommand{\phi}{\varphi}$ $\newcommand{\th}{\theta}$ $\newcommand{\om}{\omega}$ $\newcommand{\Om}{\Omega}$
Subdvision methods progressively refine a discrete mesh and converge to a smooth surface. This allows to perform an interpolation or approximation of a given coarse dataset.
addpath('toolbox_signal')
addpath('toolbox_general')
addpath('toolbox_graph')
addpath('toolbox_wavelet_meshes')
addpath('solutions/meshwav_2_subdivision_surfaces')
Starting from a control mesh which is a regular polyhedra, one can construct a sequence of mesh that converge to a sphere by subdividing each edge into two edges, and each triangle into four smaller triangles. The position of the mid points are projected onto the sphere.
Compute two examples of initial base mesh.
[vertex1,face1] = compute_base_mesh('oct');
[vertex0,face0] = compute_base_mesh('ico');
Display it.
clf;
subplot(1,2,1);
plot_mesh(vertex1,face1);
shading('faceted'); lighting('flat'); view(3); axis('tight');
subplot(1,2,2);
plot_mesh(vertex0,face0);
shading('faceted'); lighting('flat'); view(3); axis('tight');
Initialize the subdivision.
face = face0;
vertex = vertex0;
Compute the set of edges.
edge = compute_edges(face);
Number of vertex and edges.
n = size(vertex,2);
ne = size(edge,2);
Compute the number of the three edges associated to each face.
A = sparse([edge(1,:);edge(2,:)],[edge(2,:);edge(1,:)],[n+(1:ne);n+(1:ne)],n,n);
v12 = full( A( face(1,:) + (face(2,:)-1)*n ) );
v23 = full( A( face(2,:) + (face(3,:)-1)*n ) );
v31 = full( A( face(3,:) + (face(1,:)-1)*n ) );
Compute the new faces, each old face generates 4 faces.
face = [ cat(1,face(1,:),v12,v31),...
cat(1,face(2,:),v23,v12),...
cat(1,face(3,:),v31,v23),...
cat(1,v12,v23,v31) ];
Add new vertices at the edges center.
vertex = [vertex, (vertex(:,edge(1,:))+vertex(:,edge(2,:)))/2 ];
Project the points on the sphere.
d = sqrt( sum(vertex.^2,1) );
vertex = vertex ./ repmat( d, [size(vertex,1) 1]);
Display before/after subdivision.
clf;
subplot(1,2,1);
plot_mesh(vertex0,face0);
shading('faceted'); lighting('flat'); view(3); axis('tight');
subplot(1,2,2);
plot_mesh(vertex,face);
shading('faceted'); lighting('flat'); view(3); axis('tight');
Exercise 1
Perform the full subdivision.
exo1()
%% Insert your code here.
Exercise 2
Try with other control meshes.
exo2()
%% Insert your code here.
The same method can be applied to an arbitrary control mesh, but without the projection on the sphere. More clever interpolations should be used to avoid having a simple piecewise linear surface.
Load the base control mesh.
name = 'mannequin';
[vertex0,face0] = read_mesh(name);
Display it.
options.name = name;
clf;
plot_mesh(vertex0,face0,options);
shading('faceted'); lighting('flat'); axis('tight');
Initialize.
face = face0;
vertex = vertex0;
Perform the subdivision.
edge = compute_edges(face);
n = size(vertex,2);
ne = size(edge,2);
Compute the number of the three edges associated to each face.
A = sparse([edge(1,:);edge(2,:)],[edge(2,:);edge(1,:)],[n+(1:ne);n+(1:ne)],n,n);
v12 = full( A( face(1,:) + (face(2,:)-1)*n ) );
v23 = full( A( face(2,:) + (face(3,:)-1)*n ) );
v31 = full( A( face(3,:) + (face(1,:)-1)*n ) );
Compute the new faces, each old face generates 4 faces.
face_old = face;
face = [ cat(1,face(1,:),v12,v31),...
cat(1,face(2,:),v23,v12),...
cat(1,face(3,:),v31,v23),...
cat(1,v12,v23,v31) ];
Compute the vertex and face ring.
global vring e2f fring facej;
vring = compute_vertex_ring(face);
e2f = compute_edge_face_ring(face_old);
fring = compute_face_ring(face_old);
facej = face_old;
Compute the interpolated position using
for k=n+1:n+ne
[e,v,g] = compute_butterfly_neighbors(k, n);
vertex(:,k) = 1/2*sum(vertex(:,e),2) + 1/8*sum(vertex(:,v),2) - 1/16*sum(vertex(:,g),2);
end
Display before/after subdivision.
clf;
subplot(1,2,1);
plot_mesh(vertex0,face0,options);
shading('faceted'); lighting('flat'); axis('tight');
subplot(1,2,2);
plot_mesh(vertex,face,options);
shading('faceted'); lighting('flat'); axis('tight');
Exercise 3
Perform several steps of subdivision.
exo3()
%% Insert your code here.
Display the new mesh.
clf;
plot_mesh(vertex,face,options);
shading('interp'); lighting('phong'); axis('tight');
Display the new mesh faceted.
clf;
plot_mesh(vertex,face,options);
shading('faceted'); lighting('phong'); axis('tight');
Exercise 4
Try on different 3D models.
exo4()
%% Insert your code here.
Exercise 5
Implement another subdivision scheme that is not interpolating, for instance the loop scheme. Be careful about the handling of points that does not have valence 6.
exo5()
%% Insert your code here.
Display the new mesh.
clf;
plot_mesh(vertex,face,options);
shading('interp'); lighting('phong'); axis('tight');
Display the new mesh faceted.
clf;
plot_mesh(vertex,face,options);
shading('faceted'); lighting('phong'); axis('tight');
Exercise 6
Implement another subdivision scheme that does not perform a 1:4 split of each face, for instance the sqrt(3) scheme.
exo6()
%% Insert your code here.