Mesh Geometry
Meshes are collections of simple polygons, typically triangles or quads, used for representing surfaces and solids. They are are motivated by their capability for accelerated processing via GPUs and for overcoming the limitations of mathematical surfaces which require well-behaved functions in the sense of continuity and derivatives.
Representations
Triangle Lists
The simplest way to define a mesh is via a collection of 3D points, every three thereof describing a triangle. A table ( n, 3 ) of points M = [[a1, b1, c1], ... [an, bn, cn]] is often used, however, it is also common to use a ( n, 3, 3 ) tensors in numpy and pytorch or flat lists of ( n * 9 ) coordinates when uploading the data to the GPU.
class TriangleList:
def __init__( self, triangles ):
self.Triangles = triangles
def ToFloats( self ):
coordinates = []
for triangle in self.Triangles:
for point in triangle:
coordinates.append( point.X )
coordinates.append( point.Y )
coordinates.append( point.Z )
return coordinates
def ToNumPy( self ):
array = np.zeros( ( len( self.Triangles ), 3, 3 ) )
for i, triangle in enumerate( self.Triangles ):
for i, point in enumerate( triangle ):
array[i, j, 0] = point.X
array[i, j, 1] = point.Y
array[i, j, 2] = point.Z
return array
mesh = TriangleList( [
[Point3d( 0.0, 0.0, 0.0 ),
Point3d( 1.0, 0.0, 0.0 ),
Point3d( 1.0, 1.0, 0.0 )],
[Point3d( 1.0, 1.0, 0.0 ),
Point3d( 0.0, 1.0, 0.0 ),
Point3d( 0.0, 0.0, 0.0 )]] )
This representation is the most unstructured way to describe meshes also known as a soup of triangles. It is used when each point is associated with a unique color [r, g, b, a] or a texture coordinate [u, v]. Otherwise, it is considered wasteful because often triangles share points which are duplicated here. For instance a cube has 8 vertices and 6 * 2 = 12 triangles, but with this representation we need 12 * 3 = 36 vertices! Note that this definition may be trivially extended to include quads or polygons instead of only triangles.
Indexed Meshes
The common way to represent meshes is using a list of points P = [P0, P1, ..., Pn-1], also known as vertices, and a list of indices F = [F0, F1, ..., Fn-1], also known as faces. A face F[i] is a just list of integers [a, b, c], where a, b and c describe which vertices to select from the vertex-list, by index, to form a triangle, namely P[a], P[b] and P[c]. For example, the face [4, 1, 6] represents a triangle with points P[4], P[1], and P[6]. The motivation for this representation is primarily for saving memory. A cube only requires 8 vertices as opposed to 36 with the equivalent triangle soup representation.
class IndexedMesh:
def __init__( self, vertices, faces ):
self.Vertices = vertices
self.Faces = indices
''' Face Index -> Points
'''
def TriangleAt( self, face ):
ia, ib, ic = self.Faces[face]
pa = self.Vertices[ia]
pb = self.Vertices[ib]
pc = self.Vertices[ic]
return [pa, pb, pc]
mesh = IndexedMesh(
[Point3d( 0.0, 0.0, 0.0 ),
Point3d( 1.0, 0.0, 0.0 ),
Point3d( 1.0, 1.0, 0.0 ),
Point3d( 0.0, 1.0, 0.0 )],
[[0, 1, 2], [2, 3, 0]] )
Rhino's meshes are based on the latter representation. The collection of points is accessed via the Vertices property and the triangulation indices via the Faces property. Note that Rhino supports mixed meshes containing not only triangles but also quads simultaneously. Newer versions also allow meshes with polygonal faces.
""" Mesh Construction
"""
mesh = Mesh( )
a = mesh.Vertices.Add( 0.0, 0.0, 0.0 )
b = mesh.Vertices.Add( 1.0, 0.0, 0.0 )
c = mesh.Vertices.Add( 1.0, 1.0, 0.0 )
d = mesh.Vertices.Add( 0.0, 1.0, 0.0 )
mesh.Faces.AddFace( a, b, c )
mesh.Faces.AddFace( c, d, a )
''' add quadrilateral '''
mesh.Faces.Add( a, b, c, d )
Note that mesh points in Rhino are represented by Point3f using single precision floating point numbers instead of Point3d which use doubles. This is because graphics cards work faster with lower precision values. Nevertheless, Point3f does not implement the range of methods Point3d does, thus for performing computations conversion is required.
""" Retrieve Mesh Vertex
"""
vertex = mesh.Vertices[0]
""" Convert to Point3d
"""
point = Point3d( vertex )
""" Perform Computation
"""
point += Vector3d.ZAxis
""" Convert to Point3f
"""
mesh.Vertices[0] = Point3f( point.X, point.Y, point.Z )
''' equivalent '''
mesh.Vertices.SetVertex( 0, point )
Conversions
Converting between indexed meshes and triangle lists is straight forward application of the indexing logic to select the appropriate vertices and package them into a list.
class IndexedMesh:
def __init__( self, vertices, faces ):
self.Vertices = vertices
self.Faces = indices
def ToTriangleList( self ):
triangles = [self.TriangleAt( index )
for index in range( len( self.Faces ) )]
return TriangleList( triangles )
The opposite conversion, from a triangle list to an indexed mesh, requires flattening all triangles into a list of vertices and serially numbering the faces [0, 1, 2], [3, 4, 5] .... Note that this process will not merge duplicate vertices and provide a compact the triangle soup.
class TriangleList:
def __init__( self, triangles ):
self.Triangles = triangles
def ToIndexedMesh( self ):
vertices, faces = [], []
for triangle in self.Triangles:
index = len( vertices )
faces.append( [index + 0,
index + 1,
index + 2] )
for point in triangle:
vertices.append( point )
return IndexedMesh( vertices, faces )
Mesh Normals
The concept of mesh normals is associated with ideas presented in the analytical curves and surfaces section. Mesh normals are defined per mesh point and express the locally tangent plane. Therefore, the number of normals is the same with the number of points. Mesh normals can be assigned arbitrarily or computed in the approximate sense, since there are no surface derivatives for meshes as they are fundamentally non-differentiable.
This approximation is performed by averaging the plane normals from the triangles connected with each point. There are two common methods for computing mesh normals, namely weighted and unweighted. The difference is in accumulating unnormalized or normalized normals for each triangle. In other words, using the triangle's area as a weight encoded in the cross product of its sides or not. Weighted normals is motivated by the potential presence of tiny or sliver triangles which may affect the normal direction adversely.
class IndexedMesh:
def __init__( self, vertices, faces, normals ):
self.Vertices = vertices
self.Faces = faces
self.Normals = normals
if( normals is None ):
self.ComputeNormals( )
def ComputeNormals( self, weighted = False ):
self.Normals = [Vector3d.Zero
for _ in self.Vertices]
for a, b, c in self.Faces:
normal = Vector3d.CrossProduct(
self.Vertices[b] - self.Vertices[a],
self.Vertices[c] - self.Vertices[b] )
if( not weighted ): normal.Unitize( )
self.Normals[a] += normal
self.Normals[b] += normal
self.Normals[c] += normal
for normal in self.Normals:
normal.Unitize( )
return self
In Rhino mesh normals have the same semantics. They can be assigned manually or computed automatically. Note that mesh normals are Vector3f instead of Vector3d, so conversions may be required.
""" Mesh Construction
"""
mesh = Mesh( )
a = mesh.Vertices.Add( 0.0, 0.0, 0.0 )
b = mesh.Vertices.Add( 1.0, 0.0, 0.0 )
c = mesh.Vertices.Add( 1.0, 1.0, 0.0 )
d = mesh.Vertices.Add( 0.0, 1.0, 0.0 )
mesh.Faces.AddFace( a, b, c )
mesh.Faces.AddFace( c, d, a )
''' manually '''
for _ in range( mesh.Vertices.Count ):
mesh.Normals.Add( Vector3d.ZAxis )
''' computed '''
mesh.Normals.ComputeNormals( )
Face Winding
The order by which triangle points are used to define their normals is important because the cross product is not commutative. This implies that a triangle [a, b, c] will have its normal flipped if we use [a, c, b]. Meshes with inconsistent triangle winding can cause problems from visualization to 3d printing because mesh normals are used for lighting, hidden surface removal and solid closure determination.
""" Inconsistent Winding
"""
mesh = IndexedMesh(
[Point3d( 0.0, 0.0, 0.0 ),
Point3d( 1.0, 0.0, 0.0 ),
Point3d( 1.0, 1.0, 0.0 ),
Point3d( 0.0, 1.0, 0.0 )],
[[0, 1, 2], [0, 2, 3]] ) #-- vs. [2, 3, 0]
Flat Shading
Mesh normals are used for smoothly shading triangles against light sources. However, a box expressed as an indexed mesh with 8 vertices and 12 triangles arguably looks weird. This is because the normals are averaged across faces with sharp angles.
To avoid this artifact and/or in general to force a mesh to appear with each triangle flat shaded the following approach can be used. The method below copies every triangle but the points are duplicated or unshared such that the normals will not wrap around corners.
class IndexedMesh:
def __init__( self, vertices, faces, normals ):
self.Vertices = vertices
self.Faces = faces
self.Normals = normals
if( normals is None ):
self.ComputeNormals( )
def ToFlatMesh( self ):
points = []
planes = []
for index, ( a, b, c ) in enumerate( self.Faces ):
points.append( self.Vertices[a] )
points.append( self.Vertices[b] )
points.append( self.Vertices[c] )
planes.append( [index * 3 + 0,
index * 3 + 1,
index * 3 + 2] )
return IndexedMesh( points, planes )
Mesh Colors
Assigning colors to mesh vertices is common from old-school video games to 3d printing and scientific visualization. Supporting per vertex colors is usually implemented with a list containing as many colors as the number of points. Colors may be encoded as either RGB [r, g, b] or RGBA [r, g, b, a] values, where the color components r, g, b and a, are either expressed as bytes [0, 255] or floats in [0.0, 1.0].
class IndexedMesh:
def __init__( self, vertices, faces, normals ):
self.Vertices = vertices
self.Faces = faces
self.Normals = normals
self.Colors = colors
def PaintNormals( self ):
self.Colors = []
for normal in self.Normals:
''' Float RGB values [0, 1]
'''
r = ( 1.0 + normal.X ) / 2.0
g = ( 1.0 + normal.Y ) / 2.0
b = ( 1.0 + normal.Z ) / 2.0
''' Byte RGB values [0, 255]
'''
r = int( 255 * r ) & 0xFF
g = int( 255 * g ) & 0xFF
b = int( 255 * b ) & 0xFF
self.Colors.append( [r, g, b] )
return self
def PaintZ( self ):
''' Min and Max Z-coordinates
'''
zmin, zmax = 1e10, -1e10
for point in self.Vertices:
zmin = min( zmin, point.Z )
zmax = max( zmax, point.Z )
dz = zmax - zmin
dz = 0.0 if( dz == 0.0 ) else 1.0 / dz
self.Colors = []
for point in self.Vertices:
''' Float Gray Tone
'''
gray = ( point.Z - zmin ) * dz
''' Byte Gray Tone
'''
gray = int( 255 * gray ) & 0xFF
self.Colors.append( [gray, gray, gray] )
return self
Rhino uses the color format representation of .NET / Windows where alpha, red, green and blue are packed into a single 32-bit integer.
from System.Drawing import Color
mesh = Mesh( )
a = mesh.Vertices.Add( 0.0, 0.0, 0.0 )
b = mesh.Vertices.Add( 1.0, 0.0, 0.0 )
c = mesh.Vertices.Add( 1.0, 1.0, 0.0 )
d = mesh.Vertices.Add( 0.0, 1.0, 0.0 )
mesh.VertexColors.Add( Color.FromArbg( 0xFF, 0xFF, 0xFF ) ) #-- R, G, B
mesh.VertexColors.Add( Color.FromArbg( 0xFF0000 ) ) #-- RGB
mesh.VertexColors.Add( Color.FromArbg( 0x80, 0x00, 0xFF, 0x00 ) ) #-- A, R, G, B
mesh.VertexColors.Add( Color.FromArbg( 0x800000FF ) ) #-- ARGB
mesh.Faces.AddFace( a, b, c )
mesh.Faces.AddFace( c, d, a )
Mesh Textures
Associating colors per vertex and letting the GPU smoothly bend them in the interior of triangles can only go as far as painting meshes with gradients. Texture mapping enables applying raster images with much more complex patterns. This is performed by associating each vertex with a texture coordinate, that is the pixel location in the image space.
Texture coordinates u and v for the horizontal and vertical image direction, respectively, are not expressed in pixels but in the unit range [0, 1]. When the letters s and t are used instead of u and v, the vertical direction of the image is flipped, that is u and s are always left-to-right but u is bottom-to-top while t is top-to-bottom, that is v = 1 - t.
class IndexedMesh:
def __init__( self, vertices, faces, textures ):
self.Vertices = vertices
self.Faces = faces
self.Textures = textures
mesh = IndexedMesh(
[Point3d( 0.0, 0.0, 0.0 ),
Point3d( 1.0, 0.0, 0.0 ),
Point3d( 1.0, 1.0, 0.0 ),
Point3d( 0.0, 1.0, 0.0 )],
[[0, 1, 2], [2, 3, 0]],
[[0.0, 0.0],
[1.0, 0.0],
[1.0, 1.0],
[0.0, 1.0]] )
Rhino uses the same semantics via the TextureCoordinates property. Note that texture coordinates use the s and t convention and they are single precision 2D points Point2f.
""" Texture Coordinates
"""
mesh = Mesh( )
a = mesh.Vertices.Add( 0.0, 0.0, 0.0 )
b = mesh.Vertices.Add( 1.0, 0.0, 0.0 )
c = mesh.Vertices.Add( 1.0, 1.0, 0.0 )
d = mesh.Vertices.Add( 0.0, 1.0, 0.0 )
mesh.TextureCoordinates.Add( 0.0, 0.0 )
mesh.TextureCoordinates.Add( 1.0, 0.0 )
mesh.TextureCoordinates.Add( 1.0, 1.0 )
mesh.TextureCoordinates.Add( 0.0, 1.0 )
mesh.Faces.AddFace( a, b, c )
mesh.Faces.AddFace( c, d, a )
Barycentrics
Vertex normals and colors are interpolated by the GPU for shading the interior pixels of triangles in visually smooth way. This process is known as barycentric interpolation and aims to smoothly blend between three quantities, which may be vertices, normals, colors etc.
A point p within a triangle [a, b, c], where a, b and c are its vertices, may be expressed as the ratio of areas of the triangles formed by a, b, c and p. The areas are computed using the convention A = |p, b, c|, B = |p, c, a| and C = |p, a, b|, where the area for vertex a involves point p and the opposite side of the triangle bc. The sum of these areas A + B + C equals the area of the triangle T = |a, b, c|, that is A + B + C = T. By dividing with the area of the triangle produces A/T + B/T + C/T = 1.
The terms wa = A/T, wb = B/T and wc = C/T are the barycentric coordinates or the weights of the interpolation. This because we can address every point p in the interior of the triangle by the weighted sum of the vertices p = a * wa + b * wb + c * wc. The values of the weights are in the unit range [0, 1]. When the weights are [1, 0, 0], [0, 1, 0]and [0, 0, 1], then the interpolated points are on the vertices a, b and c, respectively. Moreover, when the weights are [1/3, 1/3, 1/3] the interpolated point is at the centroid of the triangle ( a + b + c ) / 3, hence the name barycentric.
Barycentric interpolation enables blending normals and colors by replacing the triangles' vertices with the associated vectors and RGB values. This technique produces smooth gradients and light shading in 3D graphics.
""" Barycentric Interpolation
"""
def ToBarycentric( triangle, point ):
a, b, c = triangle
area = Vector3d.CrossProduct( b - a, c - b ).Length
return Point3d(
Vector3d.CrossProduct( b - point, c - b ).Length / area,
Vector3d.CrossProduct( c - point, a - c ).Length / area,
Vector3d.CrossProduct( a - point, b - a ).Length / area )
def FromBarycentric( triangle, point ):
a, b, c = triangle
return ( a * point.X +
b * point.Y +
c * point.Z )
Operations
Surface Area
The surface area of a mesh is computed by the sum of the areas of its triangles. The area of triangles can be computed using the cross product.
class IndexedMesh:
def __init__( self, vertices, faces ):
self.Vertices = vertices
self.Faces = faces
@property
def Triangles( self ):
for ia, ib, ic in self.Faces:
pa = self.Vertices[ia]
pb = self.Vertices[ib]
pc = self.Vertices[ic]
yield [pa, pb, pc]
@property
def Area( self ):
area = 0.0
for pa, pb, pc in self.Triangles:
area += Vector3d.CrossProduct(
pb - pa, pc - pb ).Length
return area / 2
In Rhino the AreaMassProperties class may be used for computing mesh volumes.
""" Computing Mesh Areas
"""
props = AreaMassProperties.Create( mesh )
print( mesh.Area )
Mesh Volume
Mesh volume is computed as the sum of signed volumes of tetrahedra formed by each triangle and the world origin, see detailed derivation. The volume associated with a triangle [a, b, c] and the world origin o = [0, 0, 0] is surprisingly easy to compute using a × b · c. Note that only properly closed meshes have meaningful volume.
class IndexedMesh:
def __init__( self, vertices, faces ):
self.Vertices = vertices
self.Faces = faces
@property
def Volume( self ):
volume = 0.0
for pa, pb, pc in self.Triangles:
volume += Vector3d.CrossProduct( pa, pb ) * pc
return abs( volume / 6 )
In Rhino the VolumeMassProperties class may be used for computing mesh volumes.
""" Computing Mesh Volumes
"""
props = VolumeMassProperties.Create( mesh )
print( mesh.Volume )
Deduplication
Removing duplicate mesh vertices or converting a triangle soup to an indexed mesh with shared vertices across faces is fairly complicated. This is because it requires finding all the unique or in other words shared vertices.
A simple approach is presented below, where an initially empty list of unique vertices and a list of faces is created. For each vertex from each triangle, first we look up whether the vertex exists or not in the unique vertex list. If the vertex exists, then we know its index, otherwise we append it at the end of the list and record the index as the length of unique vertices before insertion.
This algorithm is slow so it should not be used for large meshes.
''' Checks if a point existing in a list of points
Returns the index if present, or -1 otherwise.
'''
def FindPoint( target, points, epsilon = 1e-5 ):
for index, source in enumerate( points ):
if( source.DistanceTo( target ) <= epsilon ):
return index
return -1
class IndexedMesh:
def __init__( self, vertices, faces ):
self.Vertices = vertices
self.Faces = faces
''' Removes duplicate vertices and returns a
new compacted mesh.
'''
def RemoveDuplicateVertices( self, epsilon = 1e-5 ):
vertices = [] #-- Unique vertices list
faces = [] #-- Reindexed faces list
for triangle in self.Triangles:
face = []
for vertex in triangle:
index = FindPoint( vertex, vertices )
if( index < 0 ):
index = len( vertices )
vertices.append( vertex )
face.append( index )
faces.append( face )
return IndexedMesh( vertices, faces )