Surface Geometry
The objective of this section is to introduce concepts associated with flat and curved surfaces. There are three very important operations related to surfaces:
- Evaluation: Forward mapping from R2 → R3, that is
[u, v] → [x, y, z]. - Derivatives: Geometric properties associated with differential concepts.
- Projection: Inverse mapping from R3 → R2, that is
[x, y, z] → [u, v].
Representation
Surfaces can be expressed in one of several ways, namely explicit, implicit and parametric, in the same sense there are multiple representations for curves. Those offer different kinds of intuition and capabilities. Nevertheless, the parametric representation of surfaces is the one most commonly used.
Point On Surface
Parametric surfaces can be considered as a direct extension of curves in two, namely u and v instead of one parameter, that is t. As such as surface S( u, v ) = [x, y, z] renders points in space p = [x, y, z].
In Rhino the Surface interface, which all surfaces comply with, captures the notion of a point-on-surface via the PointAt( ) method, see documentation.
""" Various Surfaces
"""
surface = Plane.WorldXY
surface = NurbsSurface.CreateFromPoints(
[Point3d( 0.0, 0.0, 0.0 ),
Point3d( 1.0, 0.0, 0.0 ),
Point3d( 0.0, 1.0, 1.0 ),
Point3d( 1.0, 1.0, 0.0 )], 2, 2, 1, 1 )
""" The PointAt Method
"""
point = surface.PointAt( 0.5, 0.5 )
Computing points on a surface along the parametric directions is commonly expressed as a grid or flat list of points as seen below.
""" Points on Normalized Surface
"""
nu, nv = 2, 3
flat = [surface.PointAt(
iu / ( nu - 1 ), iv / ( nv - 1 ) )
for iu in range( nu )
for iv in range( nv )]
grid = [[surface.PointAt(
iu / ( nu - 1 ), iv / ( nv - 1 ) )
for iu in range( nu )]
for iv in range( nv )]
Parameterization
The range of values associated with the parameters u and v follow the same conventions as with parametric curves. However, the domain of a surface requires two intervals, one for each parameter u: [umin, umax] and v: [vmin, vmax]. The minimum and maximum values are typically either associated with the unit range [0, 1] or arc-length [0, length].
In Rhino, the domain of a surface can be obtained using the Domain( ) method, where the parameter specifies which interval to select.
""" Surface Domain Intervals
"""
ud = surface.Domain( 0 )
vd = surface.Domain( 1 )
print( 'u: [{:.3f}, {:.3f}]'.format(
ud.Min, ud.Max ) )
print( 'v: [{:.3f}, {:.3f}]'.format(
vd.Min, vd.Max ) )
""" Points on Surface
"""
nu, nv = 2, 3
flat = [surface.PointAt(
ud.ParameterAt( iu / ( nu - 1 ) ),
vd.ParameterAt( iv / ( nv - 1 ) ) )
for iu in range( nu )
for iv in range( nv )]
grid = [[surface.PointAt(
ud.ParameterAt( iu / ( nu - 1 ) ),
vd.ParameterAt( iv / ( nv - 1 ) ) )
for iu in range( nu )]
for iv in range( nv )]
Derivatives
The concepts associated with parametric curves extend to surfaces but the terminology and mechanics change slightly because surfaces involve two parameters.
Tangent Plane
The tangent of a surface a point with parameters u and v is a plane instead of vector. It is defined via the process of computing the first derivative at S( u, v ). However, the result is two vectors or partial derivatives Su = dS / du and Sv = dS / dv. These define the basis vectors of the tangent plane, which are in general not orthogonal with one another.
""" Finite deltas
"""
du, dv = 1e-5, 1e-5
""" Point at u, v
"""
o = surface.PointAt( u, v )
""" Forward points
"""
pu = surface.PointAt( u + du, v )
pv = surface.PointAt( u, v + dv )
""" Finite partial derivatives
"""
Su = ( p - o ) / du
Sv = ( q - o ) / du
""" Finite tangent plane
"""
tangent = Plane( o, Su, Sv )
In Rhino the Evaluate( ) method can be used to compute partial derivative vectors of surfaces. Note that higher order partial derivatives come in expanding combinations, see Hessian matrix.
""" Surface Derivatives
"""
_, o, dS = surface.Evaluate( u, v, 2 )
""" First Order
"""
Su = dS[0]
Sv = dS[1]
""" Second Order
"""
Suu = dS[2]
Suv = dS[3]
Svv = dS[4]
Normal Vector
The normal of a surface at S( u, v ) is a vector perpendicular to the surface. It is defined as the normal of the tangent plane and computed via the cross product Su × Sv.
In Rhino the surface normal can be evaluated using the NormalAt( ) method.
""" Principle
"""
_, o, dS = surface.Evaluate( u, v, 1 )
normal = Vector3d.CrossProduct( dS[0], dS[1] )
normal.Unitize( )
""" Shortcut
"""
normal = surface.NormalAt( u, v )
Surface Frame
The surface frame at S( u, v ) is just the tangent plane as presented earlier. In Rhino, it is computed using the FrameAt( ) method. Note that planes in Rhino are orthonormal so the partial derivatives' directions and lengths are not preserved. Nevertheless, the X-axis is indeed in the direction of Su.
""" Principle
"""
_, o, dS = surface.Evaluate( u, v, 1 )
plane = Plane( o, dS[0], dS[1] )
""" Shortcut
"""
plane = surface.FrameAt( u, v )
Fundamentals
Measuring the length of curves as well as the area of regions within surfaces is associated with a concept from calculous known as the First Fundamental Form, see additional details. It is defined by three terms, E, F and G forming a 2x2 matrix computed from the first order partial derivatives' dot products as Su · Su, Su · Sv, and Sv · Sv, respectively.
Intuitively, the dot products capture the rate by which length is changing along each direction. However, because the partial derivative vectors Su and Sv are not necessarily orthogonal with one another, their combination Su · Sv needs to be considered, otherwise it would have been 0.
""" Surface Derivatives
"""
_, o, dS = surface.Evaluate( u, v, 1 )
Su = dS[0]
Sv = dS[1]
""" First Fundamental Form
"""
E = Su * Su
F = Su * Sv
G = Sv * Sv
I = [[ E, F ],
[ F, G ]]
While this concept is of theoretical significance, in practice lengths and areas are computed by approximation as seen in the case of parametric curves. Curves are discretized and the sum of the linear segment lengths is computed, while areas are measured using triangulations thereof.
The Second Fundamental Form is associated with the second order derivatives of surfaces and expresses the notion surface curvature, see additional details. Again the dot products between the normal n and the second order partial derivatives Suu, Svv and Suv, express projection lengths onto the normal, which measure the rate by which the surface is bending in each direction separately and in combinations.
""" Surface Derivatives
"""
_, o, dS = surface.Evaluate( u, v, 2 )
Su = dS[0]
Sv = dS[1]
n = Vector3d.CrossProduct( Su, Sv )
n.Unitize( )
""" Second Fundamental Form
"""
Suu = dS[2]
Suv = dS[3]
Svv = dS[4]
L = Suu * n
M = Suv * n
N = Svv * n
II = [[ L, M ],
[ M, N ]]
Normal Curvature
The complication with computing the curvature of surfaces is that unlike curves it is not a single value. For every point on a surface S( u, v ) = p there is an infinite number of curves, within the surface, passing through the point, and they all typically have different curvature values.
The way to investigate curvature therefore first requires selecting a particular direction w in the tangent plane. Then using p as the origin and w and n as basis vectors, we define a plane that intersects the surface producing a curve for which we can compute the curvature. This concept is also known as the normal curvature.
In practice, we do not need to perform a geometric intersection but only express the direction w in terms of the tangent plane's basis vectors Su and Sv. This requires solving a 2x2 linear system for s and t because the vectors are almost never orthogonal with one another. Finally, we can use the first and second fundamental forms to compute the normal curvature.
""" Direction in Tangent Plane
"""
w = Vector2d( wx, wy )
""" Express w as t * Su + s * Sv
"""
t = ( ( Sv.X * w.Y - Sv.Y * w.X ) /
( Sv.X * Su.Y - Sv.Y * Su.X ) )
s = ( ( Su.X * w.Y - Su.Y * w.X ) /
( Su.X * Sv.Y - Su.Y * Sv.X ) )
tt = t * t
ts = t * s
ss = s * s
""" Compute Normal Curvature
"""
k = ( ( L * tt + 2 * M * ts + N * ss ) /
( E * tt + 2 * F * ts + G * ss ) )
Principal Curvatures
Normal curvature has a minimum and a maximum value, namely k_min and k_max about the w_min and w_max directions in the tangent plane, respectively. Those are also known as the principal directions and they are orthogonal with one another and the normal, see derivation. Additionally, their mean H = ( k_min + k_max ) / 2 and product K = k_min * k_max are known as the Mean and Gaussian curvature of the surface at S( u, v ).
""" Mean and Gaussian Curvature
"""
K = ( ( E * N + G * L - 2 * F * M ) /
( 2 * ( E * G - F * F ) ) )
H = ( ( L * N - M * M ) /
( E * G - F * F ) )
""" Minimum and Maximum Curvature
"""
k_min = H - math.sqrt( H * H - K )
k_max = H + math.sqrt( H * H - K )
""" Principal Directions
"""
w_min = Sv * ( M - k_min * F ) / ( N - k_min * G ) - Su
w_max = Sv * ( M - k_max * F ) / ( N - k_max * G ) - Su
w_min.Unitize( )
w_max.Unitize( )
Rhino provides the CurvatureAt( ) method which returns a SurfaceCurvature object that contains the above computed surface properties.
""" Surface Curvature
"""
curvature = surface.CurvatureAt( u, v )
H = curvature.Mean
K = curvature.Gaussian
k_min = curvature.Kappa( 0 )
k_max = curvature.Kappa( 1 )
w_min = curvature.Direction( 0 )
w_max = curvature.Direction( 1 )
Projection
Projection, also known as finding the closest point, from a point p onto a surface S requires finding a point o = S( u, v ) such that p = o + n * s, where n is the surface normal at parameters [u, v] and s is a scalar. In other words, the point of intersection q between the surface and a line from p with direction perpendicular to the surface. Projection of points on simple surfaces, such as planes and spheres, can be computed analytically. For more complex curves it requires numerical computation.
In Rhino, surfaces implement the ClosestPoint( ) method, which accepts the point to project and return the parameter of the surface.
""" Surface Closest Point
"""
p = Point3d( px, py, pz )
_, u, v = surface.ClosestPoint( p )
q = surface.PointAt( u, v )