Project pipeline

For this project, we start with a .ply file containing vertices and triangular face indices. This makes it easy to render the meshes provided using MeshEdit. However, all of the .ply files are missing vertex normal data and face normal data. This data is necessary for the Ball-Pivoting Algorithm (BPA), so we had to figure out a way to extract vertex normals from a set of vertices and faces. Our approach was simple, construct the mesh with the provided faces, and use the MeshEdit code we used to compute vertex-based normals (for lighting) and instead use it for creating .txt files that contained the vertex-normal pairs we needed for BPA.

Mesh Reconstruction, Links and References

Project Proposal and Background Information
Presentation Slides
Source 1: Fausto Bernardini et al. The Ball-Pivoting Algorithm for Surface Reconstruction
Source 2: Julie Digne, An Analysis and Implementation of a Parallel Ball Pivoting Algorithm

We were unable to implement the Poisson Surface Reconstruction. However, we created a visual debugger for the Ball-Pivoting Algorithm.

Overview of the Ball-Pivoting Algorithm

The Ball-Pivoting Algorithm (BPA) is named for the simulated use of a virtual ball to help reconstruct a mesh from a point cloud.

To construct the mesh from the point cloud, one first assumes that the given point cloud consists of points sampled from the surface of an object. Points strictly represent a surface, therefore a mesh constructed from those points also represent a surface.

Using this assumption, image rolling a tiny ball across a "surface" of points. This tiny ball is dependent on the scale of the mesh, but it should be slightly larger than the average space between points. When you drop a ball onto the surface of points, the ball will get caught and settle upon three points. This is the Seed Triangle: the three points propping the ball up form a triangle, which becomes part of the mesh. From that location, the ball can roll in any direction. Roll the ball, pivoting along the triangle edge formed from two points. The ball then settles in a new location: a new triangle is formed from two of the previous vertices and one new point, the Expansion Triangle is added to the mesh. As we continue rolling and pivoting the ball, new triangles are formed and added to the mesh. The ball continues rolling and rolling until the mesh is fully formed.

The following is a 2D visualiziation of the BPA, wherein a ball rolls and settles between pairs of points. Whenever the ball settles, the two points are connected to fill in the mesh. In the 3D case, the ball settles between a triplet of points and the triplet is connected.

Nuances

The idea behind the Ball-Pivoting Algorithm is simple, but of course there are many caveats to the procedure:

• How is the ball radius chosen? The radius, $\rho$, is obtained empirically based on the size and scale of the input point cloud. In theory, the diameter of the ball should be slightly larger than the average distance between points. So, one might calculate $$\rho = 1.25\times(\dfrac{AverageDistance}{2})$$In the case brought up two bullet points below, for multiple radii, one could scale the original radius to obtain smaller or larger radii. For example, one might use a factor of $2$, $\{0.5\rho, 1.0\rho, 2.0\rho\}$ or a factor of $10$, $\{0.1\rho, 1.0\rho, 10\rho\}$.


• What if the sampled points are too far apart at some locations and the ball "falls through" the surface? When the ball pivots along an edge, it may miss the appropriate point on the surface and instead hit another point on the object or even exactly its three old points. In this case, we check that the normal of the new triangle Facet is consistently oriented with the point's Vertex normals. If it is not, then we reject that triangle.


• What if the surface has a crease or valley, such that the distance between the surface and itself is less than the size of the ball? In this case, the ball would just roll over the crease and ignore the points within the crease. But, this is not ideal behavior as the reconstructed mesh is not accurate to the object. To fix this, multiple ball radii must be used in BPA: the smaller balls are able to roll into the crease and construct that part of the surface; the larger balls would prevent the "fall throughs" noted in the above point and improve surface continuity.


• What if the surface is spaced into regions of points such that the ball cannot successfully roll between the regions? The virtual ball is dropped onto the surface multiple times at varying locations. (Actually, our implementation has the ball drop onto all unprocessed sets of three points.) This ensures that the ball captures the entire mesh, even when the points are inconsistently spaced out.

Additional peculiarities are discussed within the Algorithm sections below.

Data Structure: VoxelArray

The Ball-Pivoting Algorithm is run on a .txt file, which contains an unordered set of positions and normals for each respective point in the point cloud. A point's position and normal is stored in the Vertex structure; and, the set of Vertexs are stored in the "Cloud." In order for the BPA algorithm to work, we need a method to access all of the points within a given radius of any local point. To "roll" the ball we need quick access to all of the neighboring points to a local point. At most for the BPA, we need to obtain all points within $2\rho$ distance from the local point. To allow these queries for neighboring points, we created the VoxelArray data structure.

The VoxelArray is a 3D grid of cubic voxels, each voxel having side lengths of $2\rho$. With this voxel side length, we can successfully obtain all of the neighboring points: all points within $2\rho$ of the local point are contained within the local (center) voxel and it's 26 neighboring voxels. (Which constitutes a $3$x$3$x$3$ voxels, $6\rho$ length cube.) The points contained within the 27 local voxels can then be quickly filtered and sorted based on their distances from the local point.

The VoxelArray structure has both advantages and disadvantages.

Advantages

• Fast voxel access: The voxels are all indexed in the VoxelArray for constant time access. Each voxel is a vector container of contained points.
• Intuitive internal structure: Voxels' indices are a flattened value based upon the integer x-, y-, and z-coordinates of the voxel. Float coordinates can simply be truncated into integers for indexing. Neighboring voxels can be found by $\pm 1$ to the integer x-, y-, and z-coordinates.

Disadvantages

• Memory intensive: The entire space of the point cloud is contained within the VoxelArray cube, whose voxels are all indexed and reserved. The amount of memory reserved is directly proportional to the size of the mesh, and inversely prportional to the input radius $\rho$.
• Empty voxels: Voxels are still allocated space on the memory even when no points are contained within its volume.

Data Structures: Vertex, Edge, Facet

Our BPA implementation uses a mixture of Vertexs, Edges, and Facets in order to represent the reconstructed mesh. We do not explicitly need any "mesh" structure which encompasses Vertex, Edge, Facet, because of the nature of the .ply format. .ply files only need to store point positions and faces, which are triangles indexing the points, to represent a mesh. So our final output only requires that we keep track of a list of Vertexs, which contain point positions, and a list of Facets, which contain sets of 3 Vertexs.

Vertex

A Vertex contains a 3D vector containing a point's position and the point's normal. These two pieces of data are gathered from the input .txt file. All of the points are stored in the "Point Cloud." The Cloud is simply a std::vector of Vertexs. This also happens to be the list of points used in the final output.

A Vertex also has a flag indicating whether the it is an inner vertex, which signifies that it has already been included within the reconstructed parts of the mesh -- whether it has already been processed.

Edge

An Edge contains pointers to its two Vertexs and two Facets.

Additionally, it holds two flags of whether if (1) it is an inner edge, and if (2) it is on a boundary. An inner edge has already been processed and is part of the mesh. A boundary edge is on the border of the mesh: the surface is on one side of this edge and no further surface is on the other side of this edge. A boundary edge is essentially a "cliff" or "crack" in the mesh, when the physical object was improperly sampled or surface data is missing.

Facet

A Facet is a triangle which contains pointers to its three Vertexs.

Upon creation, the three Vertexs are oriented such that the Facet surface normal is in the same direction as each Vertex's normal. This ensures that the resulting mesh has consistently oriented triangles. Additionally upon creation, the circumcenter and circumradius of the Facet is calculated: these values are used in the virtual pivoting of the ball. The center of the $\rho$-circumsphere of the Facet on the positive normal side is also calculated.

Visual Debugger

In addition to implementing BPA, we created a visual debugger. It renders the creation of Seed Triangles and Expansion Triangles in real time. Additionally, the BPA function calls can be paused and stepped through triangle-by-triangle.

Video demo of BPA run on the Dragon point cloud. Multiple radii were used, this can be seen as the multiple passes of the BPA when holes in the mesh get filled in. During the demo, the camera zooms in at around 0:30: the green triangles are Seed Triangles, the yellow triangles are Expansion Triangles, the red triangles are Hole Filling Triangles, which are discussed in the Hole Filling Case below.

For the sake of the demo, smaller radii were used. A large input radius $\rho$ would cause the BPA to run extremely slowly under the visual debugger. This is because the VoxelArray creates Voxels of size $2\rho$, so larger $\rho$ values lead to a greater amount of neighboring points which must be processed while "rolling the ball."

In creating the render-debugger, we needed two threads to accomplish the asynchronous rendering and meshing. The main thread renders the model, and receives new triangles via three separate callbacks for Seed Triangles, Expansion Triangles, and Hole-filling Triangles. The meshing thread just runs the BPA meshing algorithm, and sends triangles to the main thread through callbacks given to the Mesher object from the main thread.

Vertex Normal Debugger

In addition to the BPA visual debugger, we created a vertex normal debugger, for the .ply.txt conversion part of the Project Pipeline.

We couldn't think of a way to test our code for creating vertex normals without viewing the normals, and our initial rendering code would just display the point cloud like this:

We decided to add some code to display the normals, but no matter what depth we viewed them at, it was difficult to tell whether or not they were oriented correctly.

Finally, we had the idea of adjusting a few parameters to make visual vertex normal debugging possible.

• Extend the length of the normal vector
• Decrease the far-clip of the camera
• Make the vertices red
• Make the normal vectors blue
• Make the end of the normal vectors green

Then, if we zoom in on the mesh, we should see green before we see red, which would tell us that the normals are indeed facing outwards.

Using this method, we were able to verify our code for calculating vertex normals and produce the .txt data we needed, which has this format:

        vx1 vy1 vz1 nx1 ny1 nz1
        vx2 vy2 vz2 nx2 ny2 nz2
        .   .   .   .   .   .
        .   .   .   .   .   .
        .   .   .   .   .   .
        vxm vym vzm nxm nym nzm
    

Here is an example:

        5.89458 11.7884 27.2832 -0.83141 -0.404449 0.381023
        -53.3251 67.1044 -57.4501 -0.878005 -0.293709 0.377946
        3.75049 16.5054 29.454 -0.736516 -0.579169 0.349439
        -40.0742 -33.2375 10.7422 -0.859218 -0.0742719 0.50619
        -52.8133 67.0694 -57.5725 -0.844274 -0.203288 0.495858
        -52.4255 66.8347 -57.3804 -0.820147 -0.109775 0.561524
    

The Ball-Pivoting Algorithm

The BPA takes in a oriented point Cloud of Vertexs. That is, each point has a position and a normal. Additionally, BPA requires a user-given radius or set of radii.

While the Find Seed Triangle algorithm completes successfully, we call Expand Triangulation on the returned Seed Triangle. Otherwise, we break and the BPA is completed for the given radius. If multiple radii are given, the Find Seed Triangle-Expand Triangulation loop is run again until all radii are used.

The conclusion to BPA is a set of Vertexs and a set of Facets, whose data will be outputted into a .ply file.

Algorithm: Find Seed Triangle

Find Seed Triangle iterates throughs the entire point cloud and attempts to find valid Seed Triangles. Iterating through the Cloud, this function searches for two other points in the neighborhood of the current point. The three points must be in Empty Ball Configuration in order to return a valid Facet for Expand Triangulation.

Empty Ball Configuration

Empty Ball Configuration signifies that a sphere of radius $\rho$ can be placed onto the three points with no other point contained within that sphere. The sphere is placed such that it is on the positive normal side of the Facet that would be formed by those three points. The intuitive reasoning behind the Empty Ball Configuration is that we want to drop the ball onto the surface of the object: if some other point is contained within the ball at any time, then the ball has essentially clipped through the physical surface of the object. We want the ball to be strictly on the surface, so that the created Facets do in fact represent the surface of the object.

Algorithm: Expand Triangulation

Expand Triangulation takes in a Facet returned by Seed Triangle and uses that Facet to start the Edge Front. The Edge Front is a queue of Edges that the ball can pivot over and create a new triangle Facet for the final mesh. Within a loop, an edge is popped off the Edge Front and a ball is pivoted over that edge: A valid Vertex may be returned by Find Candidate. If so, the new valid point and the two points of the pivot edge form a new Facet for the mesh. Then the two new Edges, formed by the first old point + the new point, and the second old point + the new point, are pushed on to the Edge Front queue (if they are still available). The loop continues until the Edge Front is empty: new Edges are not always added to the Edge Front.

Throughout the entire loop, we must ensure that the Edges being used are still valid and are available to be used. That is, we must make sure each Edge has at most two adjacent Facets: if so then the Edge is marked as an inner edge, its job is finished and cannot be used anymore. If Find Candidate does not return a valid Vertex upon which the ball pivots onto, the Edge is marked as a boundary edge -- it only has one adjacent Facet, but its job is finished nonetheless. The ball in Expand Triangulation never pivots over any edge already tagged as inner or boundary, these edges are removed from the Edge Front.

Case: Hole-Filling Triangle

Hole-Filling Triangles are a special variant of Expansion Triangles. These are triangles that essentially fill in the hole as two Edge Fronts collide with each other or an Edge Front with itself. A Hole-Filling Triangle results from a ball pivoting upon a popped pivot Edge, and the two resulting Edges adjacent to the new Facet are both on an Edge Front. This means that all of the three Edges adjacent to the Facet are marked as inner. The local region of the mesh is completely filled.

Hole-Filling Triangles also appear in a post-processing step. For a variety of reasons (points too far apart, bad normal processing, etc.) a mesh may contain unintentional holes in its surface. Thus we check that if any three boundary Edges form an oriented loop, they do in fact form a valid triangle. This triangle made from boundary Edges is also a Hole-Filling Triangle.

Algorithm: Find Candidate

Called by Expand Triangulation, given an Edge this function returns a valid Vertex upon which the ball can pivot onto. (Or returns null, indicating that there are no points available to pivot towards and the given Edge should be marked as boundary.) Essentially, this function simulates the rolling of the ball: the first point that the ball hits while pivoting is returned.

To find a valid candidate Vertex, let the given Edge be $e$, let $c = (center\ of\ \rho$-$circumsphere\ of\ existing\ Facet\ adjacent\ to\ e)$ and let $m = (midpoint\ of\ e)$. We calculate $\rho' = ||m - c|| + \rho$. This $\rho'$ is the distance from $m$ to the center of the ball as the ball pivots around $e$, at all orientations of the ball. Then, we query the VoxelArray for all points within $2\rho'$ distance to $m$, sorted in order of increasing distance. The returned points are all of the points the ball can possibly hit while pivoting upon $e$. We then find a Vertex $v$ such that $v$:

• $v$ is not an inner Vertex of the mesh, and is not either of $e$'s two Vertexs.
• $v$ is compatible with $e$. That is, the Facet $f_{new}$ formed by $v$ and $e$ has a normal oriented in the same direction as the normals of $v$ and $e$'s Vertexs. (All three dot products are positive.)
• the $f_{new}$ triangle has a circumcenter that exists: the square circumradius is less than $\rho$. Essentially, the ball fits onto the three new points and does not fall through.
• $v$ and $e$ are in Empty Ball Configuration.
• The angle between $e$'s old Facet and $f_{new}$ is minimized: the returned $v$ is the first point that the ball hits while pivoting.

Iterating through all of the neighboring points, the $v$ meeting the above conditions is returned. Otherwise, null is returned.

Multiple Radii and Post-Processing

The BPA is run multiple times with varying radii to improve the accuracy of the output mesh. We first run BPA with the smallest radii, which fills in the smallest creases of the mesh (as mentioned in section Nuances). Then BPA is run with increasing radii until all input radii have been processed. Each time, BPA is run on points not contained in an Edge AND the points that are contained within boundary edges. The latter is because the Edges marked as being on a boundary cannot find a valid Candidate Vertex due to $\rho$ being too small on the previous BPA runs.

After finishing up with running the list of radii, a Hole-Filling post-processing step is run. For a variety of reasons (points too far apart, bad normal processing, etc.) a mesh may contain unintentional holes in its surface. Thus we check that if any three boundary Edges form an oriented loop, they do in fact form a valid triangle. This triangle made from the boundary Edges is a Hole-Filling Triangle.

Results

Stanford Dragon

Couldn't get a watertight mesh for this one. It appears that the distance between vertices just changes too much.

Armadillo man

Here's an example of a really bad mesh reconstruction. We couldn't find a list of radii that would work with the armadillo man, so we just chose our own, and they didn't work very well.

Simple Sphere

This example worked really well. This was because it had a good radius chosen. Because the point cloud was small, it was really easy to try different radii quickly and get different results.

Conclusions

We had initially set out for implementing both BPA and the Poisson Reconstruction, though throughout the course of the workload we had shifted focus to implementing BPA and creating a way to visualize our BPA algorithm in real time. Our implementation is run faster on the command line with no visualizer.

As we increase the input radius $\rho$, the running of BPA slows down proportionally. This is because the size of the Voxels inside the VoxelArray is dependent on $\rho$, specifically their side length is $2\rho$. Larger $\rho$ values lead to a greater amount of neighboring points which must be processed per local point while rolling the ball.

So to create as complete of a mesh as possible, we must use a wide range of input radii. Processing time to create a complete mesh is rather significant, especially if arbitrary input radii are used.

Efficient data structures are key for faster BPA processing and for visual debug rendering. And visual demos on smaller data sets are great for illustrating what is occuring step-by-step within the algorithm, helping with debugging.

Contributions

We all worked on the project together, and had significant inputs during all parts of the project. For the project proposal, work was split pretty evenly.
For the coding, Michael did a majority of the actual typing and Allen typed a few of the data structures. While typing is an individual task, all team members were present and contributing intellectually to the code. All members worked together to understand the Ball-Pivoting algorithm and catch edge cases.
Brett created the presentation slides and designed the flow of the presentation. Michael and Allen helped fill in the presentation slides, adding especially to the parts they respectively presented.
Allen wrote the majority of the writeup's body. Brett and Michael contributed to the Debugger sections; and the Results and Conclusions sections especially.