Most Efficient Way to Calculate an SDF to a Triangle Mesh

Question

Hi.

I've been gathering info from various sources consistently over the past month, yet without any luck of finding an idea that would suit my particular issue. So here's the formulation of the problem:

Given a mesh in a form of buffer geometry (32-bit arrays of vertex coords and vertex indices + additional arrays, such as vertex normals, uvs, or tangents) calculate a signed distance function (SDF) for a uniform grid of points surrounding the geometry.

To be more specific, I intend to create something resembling a MetaBall object in 3D engines such as Maxon's Cinema4D or Blender. I've successfully implemented distance functions for all geometric primitives, but an arbitrary mesh SDF required me to implement a brute force approach - testing the distance for each triangle of the mesh for each grid point - which, of course, gets really slow for complex meshes.

Now, I recalled that most of these problems require one to build a tree-like structure, such as an Octree, KD-tree, BSP-tree, or an AABB-tree. Then I found a few articles on the so-called Fast-Sweeping algorithm (for solving an Eikonal equation) which requires one first to fill the grid points laying on the boundary (in my case the mesh, or closest to the mesh) with 0 and the remaining with large values (Infinity) and then solve the non-linear hyperbolic boundary value problem iteratively (Gauss-Seidel). I also found an open-source implementation of a mesh SDF method in the CGAL library. Alternatively, I also thought about using some shader library (like GLSL) and perhaps try building the trees using GPU, but I have never used shaders in a JS or TS project.

The step I keep getting stuck on isn't just picking the best option, but actually practically using, at least one of these methods efficiently. For example:

• If I wanted to implement the Fast-Marching Method, I'd have to iterate through all triangles and then for each one, iterate through all grid points Gijk, and using something similar to a Marching Cubes lookup table for grid cell intersections (but with even more options), I'd interpolate values close to 0 for intersected cell vertices. I have a feeling that this would take unnecessarily long, and prove to be unsuitable for real-time updating.

• I've managed to find some examples of Ray Marching SDF computation in Unity. Also, since I've never tried computing anything directly on GPU, I have no idea what the limits for example to parallel computing on it actually are, nor do I understand how such computations get carried out. Can I parallelize the computation of distance to each triangle and then quick-sort all the distances for each grid point Gijk ? If so, how would I include it into a TypeScript project?

• Say I build an AABB-tree around all the triangles in the mesh (which should be O(n * log(n)) ), then given a grid point Gijk (assuming the tree root is a box containing both Gijk and all triangles), I'd search for the closest leaf and compute the Euclidean distance to the triangle contained within (or if there's more than one, pick the). (Please, correct me if I got this wrong) The search should be O(log(n)), and if the user moves the mesh, an update would take O(n), where n is the number of triangles. So if m is the grid size, then an entire SDF computation would be O(m * n * log(n)) ? That doesn't look like an improvement against an O(n * m) brute force approach, but perhaps I've got the complexity wrong.

I thought about combining multiple approaches, but it all seems very time consuming. I also thought of using the CGAL library, or perhaps refactoring it for the purposes of my project, but I find it quite difficult to understand the C++ code because of all the dependencies inside the library. Has anyone of you done something similar to this? What would you recommend?

Thank you for all the insights.

Answer 1

I solved this using an approach outlined by three (four) main steps:

1. Generate an AABB tree of the mesh triangle soup

2. Use the AABB tree to subdivide regular bounding cubes whenever a cube intersects the mesh (fast intersection using AABB tree) forming an Octree. When subdivision reaches leaves, exact squared distances (cube centroid to triangle) are computed and a square root is taken of the smallest, writing it into a regular grid based on cube min-max coordinates.

3. Using exact distance values on mesh-intersecting voxels and a large value elsewhere, run the Fast Sweeping algorithm 8 times and fill the scalar grid with distance values.

4. (optional) flood fill the outside voxels of a negated distance grid to make sure voxels inside an outline of marked initial voxels (with exact values) have negative sign (this is, unfortunately, the slowest part for grids with resolution > 100^3 )

For a more detailed look on my implementation with computation times, feel free to read my blog posts: https://mshgrid.com/blog/

Thanks for the upvotes and comments :) .