10.9 Surface Normals
When a surface is not planar but is near planar, it can be difficult to perceive its curvature. In such cases, it is very useful to have a mechanism that allows us to better visualize this curvature without changing the shape of the surface. The superimposition of normal vectors on the surface is a simple process to improve this visualization, as can be seen in this figure. On the left, we see a surface where the curvature is not perfectly obvious. On the right, we see the same surface but with the normal vectors clearly showing the surface’s curvature.
A near planar surface (on the left) and the same surface with its normal vectors, allowing a better understanding of its curvature (on the right).
To superimpose the normal vectors, we have to calculate their position and direction. We can assume that the surface is discretized in a mesh composed by a succession of quadrilaterals that extend along two directions (as it is produced, for example, by the function map_division). To better perceive the surface’s curvature, we can position each vector at the center of each of these quadrilaterals, following the direction of the normal to each quadrilateral.
To calculate the center of a quadrilateral we can employ the approach shown in this figure. The center of the quadrilateral defined by the points \(P_0\), \(P_1\), \(P_2\) and \(P_3\) can be obtained by firstly determining the midpoints \(M_{02}\) and \(M_{13}\) of the quadrilateral’s diagonals and, secondly, the midpoint \(M\) of the segment defined by \(M_{02}\) and \(M_{13}\). Translating this into Julia we have:
midpoint(p0, p1) = xyz((p0.x+p1.x)/2, (p0.y+p1.y)/2, (p0.z+p1.z)/2)
quadrilateral_center(p0, p1, p2, p3) =
midpoint(midpoint(p0, p2), midpoint(p1, p3))
The center of a quadrilateral as the midpoint of the segment defined by the diagonals’ midpoints.
The calculation of the normal to a quadrilateral is straightforward when the quadrilateral is planar, as we only need to produce the cross product of two adjacent edges (i.e., edges that share a vertex), as shown on the left of this figure. However, in the case of a non-planar quadrilateral (on the right of this figure), that logic can no longer applicable, because the normals at each vertex are different from each other.
Normal at the center of a planar quadrilateral (on the left) and a non-planar quadrilateral (on the right).
However, it is possible to compute an approximate normal using Newell’s method [Sutherland et al 1974], which is based on the fact that the areas of a polygon’s projections onto the Cartesian planes are proportional to the components of the normal vector to that polygon. Newell’s method calculates these areas and, from them, derives the normal vector \(\vec{n}\). Mathematically speaking, given a polygon with \(n\) vertices \(V_0,V_1,\ldots,V_n\), the method consists in calculating the components of the normal vector \(\vec{n}\) through the summation of areas:
\[n_x=\sum_{i=0}^{n}(V_{i_y}-V_{{i+1}_y})\cdot(V_{i_z}+V_{{i+1}_z})\] \[n_y=\sum_{i=0}^{n}(V_{i_z}-V_{{i+1}_z})\cdot(V_{i_x}+V_{{i+1}_x})\] \[n_z=\sum_{i=0}^{n}(V_{i_x}-V_{{i+1}_x})\cdot(V_{i_y}+V_{{i+1}_y})\]
To transform the resulting vector into a unit vector (also called normalized vector, we simply have to normalize it, i.e., divide each component by its magnitude. Translated into Julia, we have:
cross_product(v0, v1) =
vxyz((v0.y-v1.y)*(v0.z+v1.z),
(v0.z-v1.z)*(v0.x+v1.x),
(v0.x-v1.x)*(v0.y+v1.y))
cross_products(vs) =
length(vs) == 1 ?
vxyz(0, 0, 0) :
cross_product(vs[1], vs[2])+cross_products(vs[2:end])
normalized_vector(v) =
let m = sqrt(v.x^2+v.y^2+v.z^2)
vxyz(v.x/m, v.y/m, v.z/m)
end
polygon_normal(vs) = normalized_vector(cross_products([vs..., [vs[1]]...]))
In the case of a quadrilateral, we can simply define:
quadrilateral_normal(p0, p1, p2, p3) =
polygon_normal([p1-p0, p2-p1, p3-p2, p0-p3])
Note that some of the operations previously outlined are predefined in Khepri, particularly, the cross_product (called cross), the normalized_vector (called unitized), the polygon_normal, the quadrilateral_center (called quad_center) and the quadrilateral_normal (called quad_normal).
To visually represent the normal to a quadrilateral, we are going to use a very thin cylinder whose size will be proportional to the size of the quadrilateral, so that it automatically adapts to different dimensions:
normal(p0, p1, p2, p3) =
let c = quadrilateral_center(p0, p1, p2, p3),
n = quadrilateral_normal(p0, p1, p2, p3),
d = distance(p0, p2)
cylinder(c, d/40, c+n*d*1.5)
end
The next step is to go through all the mesh quadrilaterals and construct the normal to each one:
normals(ptss) =
map((pts0, pts1) -> map((p0, p1, p2, p3) -> normal(p0, p1, p2, p3), pts0, pts1, pts1[2:end], pts0[2:end]),
ptss,
ptss[2:end])
Even though the visualization of the normals is complementary to the visualization of the surface, it does not make much sense to view the first without also viewing the second. Therefore, we will define a function that combines the two visualizations:
surface_and_normals(ptss) = begin
normals(ptss)
surface_grid(ptss)
end
To test this function, we can try to apply it to the generation of normals to the Möbius strip (defined in section The Möbius Strip):
surface_and_normals(mobius_strip(0, 4*pi, 160, 0, 0.3, 5))
This figure shows the result of evaluating the previous expression.
The Möbius strip surface with the superimposed normals.
