CS 184: Computer Graphics and Imaging, Spring 2018
Project 2: Mesh Editor

Part 1: Bézier curves with 1D de Casteljau subdivision

De Casteljau's algorithm is one method to evaluate a point on a Bézier curve with recursion. It starts with an initial series of $n$ control points $\{p_0, p_1, ... , p_{n-1}\}$, and recursively linearly interpolates each segment $\{\{p_0, p_1\}, \{p_1, p_2\}, ... , \{p_{n-2}, p_{n-1}\}\}$ by a factor of $t$ (like convolution in one dimension), reducing the number of control points by one, resulting in $\{(1-t) * p_0 + t * p_1, (1-t) * p_1 + t * p_2\, ... , (1-t) * p_{n-2} + t * p_{n-1}\}$, until there be one point left in the set. Following this scheme, we can naively interpolate a Bézier curve by sampling multiple points on the curve at different $t$'s. Below shows a progression from the original series of six control points to the resulting curve.

By moving those control points around, we can create a loopy shape with Bézier curve too.

Part 2: Bézier surfaces with separable 1D de Casteljau subdivision

Beyond Bézier curves, we can interpolate Bézier patches with de Casteljau's algorithm as well. To give an example of a surface defined by $16$ control points, $\{p_{0,0}, p_{0,1}, p_{0,2}, p_{0,3}, p_{1,0}, p_{1,1}, p_{1,2}, p_{1,3}, p_{2,0}, p_{2,1}, p_{2,2}, p_{2,3}, p_{3,0}, p_{3,1}, p_{3,2}, p_{3,3}\}$. We can first interpolate parallel lines on the surface defined by $\{p_{0,0}, p_{0,1}, p_{0,2}, p_{0,3}\}$, $\{p_{1,0}, p_{1,1}, p_{1,2}, p_{1,3}\}$, $\{p_{2,0}, p_{2,1}, p_{2,2}, p_{2,3}\}$ and $\{p_{3,0}, p_{3,1}, p_{3,2}, p_{3,3}\}$; at factor $u$, we get $4$ intermediate control points $\{p_{0,0}^{0,3}, p_{1,0}^{1,3}, p_{2,0}^{2,3}, p_{3,0}^{3,3}\}$ (like on a cross section of the surface). With those intermediate control points, we can further interpolate the curve at that cross section by varying another factor $v$ to get a point on the Bézier surface $p_{0,0}^{3,3}$. The $v$ here does a similar job to $t$ in part one.

Also, it does not matter in which order the two-part-interpolation is carried out. We can choose to first interpolate with $\{p_{0,0}, p_{1,0}, p_{2,0}, p_{3,0}\}$, $\{p_{0,1}, p_{1,1}, p_{2,1}, p_{3,1}\}$, $\{p_{0,1}, p_{1,2}, p_{2,2}, p_{3,2}\}$ and $\{p_{0,3}, p_{1,3}, p_{2,3}, p_{3,3}\}$ at $v$ to get intermediate control points $\{p_{0,0}^{3,0}, p_{0,1}^{3,1}, p_{0,2}^{3,2}, p_{0,3}^{3,3}\}$. Afterwards, with those intermediate points we can further interpolate a point on the surface $p_{0,0}^{3,3}$ at a factor $u$.

Below is an example with Bézier surfaces interpolated with the above algorithm.

Part 3: Average normals for half-edge meshes

In this part, I implemented a vertex normal approximation based off the vertex's neighboring triangular faces that would allow smooth surface shading. In finer detail, it takes a weighted average of the neighboring face normals with respect to their areas. At a vertex $\vec{v}$, we can discover one of its outgoing half edge $h_0$; following this half edge we can find the vertices $\{\vec{v}_0, \vec{v}_1, \vec{v}_2\}$ that define the triangle, to which the half edge $h_0$ belongs, and can calculate the face normal by taking a cross product of $\vec{v}_0 - \vec{v}_1$ and $\vec{v}_2 - \vec{v}_1$. The area of the triangle is half of the magnitude of the cross product. After iterating through the triangles around the vertex $v$ we can find the sum of face normals and that of face areas too. For completeness of explanation, we can move from one outgoing half edge to a next by finding its twin's next half edge.

To optimize the computation slightly, we can postpone the intermediate vector normalization and weighing for each face until after processing all the faces:

The produced result is shown below:

Part 4: Half-edge flip

The edge flip function takes a non-boundary edge and reconnects the edge to the opposing vertices against the edge provided instead. Below is an illustration of the change in relationships between half-edges $\{h_0, ... h_5\}$, vertices and faces $\{f_0, f_1\}$. To give an example of each type of change, for the half-edge $h_0$, its new next half-edge becomes $h_4$ and its vertex becomes the same as that of $h_5$, its twin, edge and face staying the name. The vertex for $h_0$ before flipping needs $h_3$ as its new half-edge. For face $f_0$, its halfedge needs to be one of $\{h_0, h_4, h_3\}$ after edge flipping; before, it may be one of $\{h_0, h_2, h_4\}$.

For the implementation, I first exhaustively set pointers to all of the half-edges $\{h_0, ... , h_5\}$. Then, I just reassigned each of the properties for the half-edge to its new state, also taking care of updating the vertices that $h_0$ and $h_1$ connects to before flipping and updating the half-edges for the two faces.

A trick I used while implementing the function is to utilize Halfedge::setNeighbors(HalfedgeIter, HalfedgeIter, VertexIter, EdgeIter, FaceIter) to set every possible properties of the half-edges; it is better to cover all the property updates than to prematurely optimize the code. To debug the integrity of the mesh after edge flipping, I would manually counter-check the element properties (pointer values) and see, for example, whether the half-edge the face points to is actually on the face itself. The rendering result is shown below.

I didn't embark onto an eventful debugging journey while initially implementing this feature; however, later when working on part six, I found myself changing two lines of the code in this part that I thought I had made a typo. Afterwards I returned to working on part six. When the result from part six began to look weird, like faces crossing over each other, it took me about an hour to realize that I did not verify the changes made in this function which caused the issue. What I learnt: The vertex of a half-edge is its originating vertex, not the one it points to.

Part 5: Half-edge split

Below is a diagram I used while implementing the half-edge split. I initially set pointers to the new half-edges $\{h_6, ... , h_{11}\}$, to the new vertex inserted around the midpoint of the edge being split, to each of the three new edges and to the existing and the new faces $\{f_0, ... , f_3\}$. Then I went through each of these elements and update its relations to the neighboring ones.

For instance, the half-edge $h_8$'s vertex will point to the newly inserted vertex at the center, its next half-edge will be $h_2$, its twin half-edge will be $h_1$, its edge will point to a newly inserted edge and its face will point to $f_2$. For the newly inserted vertex, its halfedge could be any of $\{h_8, h_9, h_{10}, h_{11}\}$. Again, an implementation trick I used is to update as comprehensively as possible all attributes before optimizing. Otherwise, it is likely that something ends up not updated. With that, I didn't go on an epic debugging quest for this part.

Below is a pair of screenshots of the result after some edge splits.

Edge splits are compatible with edge flips. Below is a pair of screenshots of the result after a combination of both edge splits and flips.

I also implemented edge splits along boundary edges. Below is a pair of screenshots of this feature in action.

Part 6: Loop subdivision for mesh upsampling

My implementation of Loop subdivision can be broken down into four distinct steps. Each step is demonstrated in the slides below.

1. Pre-compute positions for the existing vertices and for the new ones;
2. Split all existing edges;
3. Flip new edges that touches an old vertex;
4. Assign all vertices with their new positions.

For an existing or new non-boundary vertex, its new position will be a weighted average of the positions of the surrounding vertices. An existing vertex of degree $n$ will a weight of $u = \frac{3}{16}$ if $n = 3$ or $u = \frac{3}{8n}$ otherwise for its neighboring vertex positions, and $1 - n*u$ for its current position. A newly inserted vertex along a shared edge of two triangles will have a weight of $\frac{3}{8}$ for positions of vertices along the edge, and $\frac{1}{8}$ for the furthest vertices' positions.

An implementation trick I used is to step through each of these steps and make sure the result is consistent with what's expected. It is a little easier to debug this part as its implementation details are heavily dependent on the previous two and it's a little easier to see the bigger picture.

On the result it produces, quite noticeably, after Loop subdivision, the sharp corners and edges disappear while the mesh shrinks slightly in size.

To preserve some of the sharp corners and edges, I found it quite effective to pre-subdivide areas that need sharper details. The bottom right result presents a flat surface from the pre-subdivided triangular surface.

An issue with Loop subdivision comes in when the topological model of the mesh doesn't well represent the geometry. For the symmetric mesh in dae/cube.dae, due to the asymmetric topology, the repeatedly subdivided mesh can result in a loss of symmetry. Below is the mesh before and after 4 times Loop subdivision.

To alleviate the asymmetry in the result geometry, we can preprocess the topology to be more symmetric while retaining the same geometry. The pair of screenshots below shows the mesh before and after 4 times Loop subdivision. On the left is the preprocessed mesh (by splitting the edge on each face once), with each face having 4 rotationally symmetric triangles. Now, the subdivided mesh is in better symmetry.

Special case at mesh boundary (extra credit)

For the completeness of Loop subdivision, the implementation is also expected work with meshes with boundaries (creases), which add a special case for updating vertex positions.

Part 7: Draw something interesting!

I thought making a pumpkin would be quite fun, so I ended up designing a single pumpkin, topologically equivalent to a UV-sphere. Then I used particle system to generate a bunch of those pumpkins.

I customized two additional shaders following modeling the scene. (However, they are not really meant for generic use.) And the results are shown below: