python/todo/meshwav_2_subdivision_surfaces.ipynb

{"worksheets": [{"cells": [{"cell_type": "markdown", "metadata": {}, "source": ["Subdivision Surfaces", "====================", "", "*Important:* Please read the [installation page](http://gpeyre.github.io/numerical-tours/installation_python/) for details about how to install the toolboxes.", "$\\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}$"]}, {"cell_type": "markdown", "metadata": {}, "source": ["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."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["from __future__ import division", "import nt_toolbox as nt", "from nt_solutions import meshwav_2_subdivision_surfaces as solutions", "%matplotlib inline", "%load_ext autoreload", "%autoreload 2"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Subdivision of a Regular Polyedra", "---------------------------------", "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."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["[vertex1, face1] = compute_base_mesh('oct')", "[vertex0, face0] = compute_base_mesh('ico')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display it."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Initialize the subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["face = face0", "vertex = vertex0"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the set of edges."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["edge = compute_edges(face)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Number of vertex and edges."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["n = size(vertex, 2)", "ne = size(edge, 2)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the number of the three edges associated to each face."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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))"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the new faces, each old face generates 4 faces."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["face = [cat(1, face(1, : ), v12, v31), ...", "    cat(1, face(2, : ), v23, v12), ...", "    cat(1, face(3, : ), v31, v23), ...", "    cat(1, v12, v23, v31)]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Add new vertices at the edges center."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["vertex = [vertex, (vertex(: , edge(1, : )) + vertex(: , edge(2, : )))/ 2]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Project the points on the sphere."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["d = sqrt(sum(vertex.^2, 1))", "vertex = vertex ./ repmat(d, [size(vertex, 1) 1])"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display before/after subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["__Exercise 1__", "", "Perform the full subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo1()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}, {"cell_type": "markdown", "metadata": {}, "source": ["__Exercise 2__", "", "Try with other control meshes."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo2()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Triangulated Mesh Subdivision", "-----------------------------", "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."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["name = 'mannequin'", "[vertex0, face0] = read_mesh(name)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display it."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["options.name = name", "", "plot_mesh(vertex0, face0, options)", "shading('faceted'); lighting('flat'); axis('tight')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Initialize."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["face = face0", "vertex = vertex0"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Perform the subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["edge = compute_edges(face)", "n = size(vertex, 2)", "ne = size(edge, 2)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the number of the three edges associated to each face."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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))"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the new faces, each old face generates 4 faces."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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)]"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the vertex and face ring."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Compute the interpolated position using"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["for k in 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)"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display before/after subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["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')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["__Exercise 3__", "", "Perform several steps of subdivision."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo3()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display the new mesh."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["plot_mesh(vertex, face, options)", "shading('interp'); lighting('phong'); axis('tight')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display the new mesh faceted."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["plot_mesh(vertex, face, options)", "shading('faceted'); lighting('phong'); axis('tight')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["__Exercise 4__", "", "Try on different 3D models."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo4()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}, {"cell_type": "markdown", "metadata": {}, "source": ["__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."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo5()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display the new mesh."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["plot_mesh(vertex, face, options)", "shading('interp'); lighting('phong'); axis('tight')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["Display the new mesh faceted."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["plot_mesh(vertex, face, options)", "shading('faceted'); lighting('phong'); axis('tight')"]}, {"cell_type": "markdown", "metadata": {}, "source": ["__Exercise 6__", "", "Implement another subdivision scheme that does not perform a 1:4 split", "of each face, for instance the sqrt(3) scheme."]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["solutions.exo6()"]}, {"cell_type": "code", "language": "python", "metadata": {}, "outputs": [], "collapsed": false, "input": ["## Insert your code here."]}]}], "nbformat": 3, "metadata": {"name": ""}, "nbformat_minor": 0}