Upcoming Events
Unite 2010
11/10 - 11/12 @ Montrιal, Canada

GDC China
12/5 - 12/7 @ Shanghai, China

Asia Game Show 2010
12/24 - 12/27  

GDC 2011
2/28 - 3/4 @ San Francisco, CA

More events...
Quick Stats
115 people currently visiting GDNet.
2406 articles in the reference section.

Help us fight cancer!
Join SETI Team GDNet!
Link to us Events 4 Gamers
Intel sponsors gamedev.net search:

Normal Computations for Heightfield Lighting


Normal Calculations for Heightfields

Method 1: Mean Weighted Equally (MWE)

Mean Weighted Equally is a two stage algorithm for computing the normal at a vertex in a heightfield (vertex normal). The first stage in the process is to compute the surface normals of each of the triangles in the heightfield using cross-product calculations. The computed surface normals can either be normalized or not, depending on the quality of vertex normals you are looking for. The second stage of the process is to add up all of the surface normals attached to a given vertex and then either average the resultant vertex normal or normalize it. The following two equations are given as the mathematical definitions for the surface and vertex normals.

Surface Normal: Ns ═ A Χ B
Vertex Normal: Nnwe = Σ (Nsi) / 6.0f | i = 0…5

The equation for the surface normal above can be read as "The Surface Normals is equal to the vector A crossed with the vector B" where A and B are two vectors that lie in the plane of the surface.

The equation for the vertex normal above can be read as "The vertex normal is the sum of all surface normals which this vertex is a part of, divided by 6, where the surface normal index is in the range of 0 to 5, inclusive. I should take a moment to note that in Gouraud's algorithm he does not strictly state that there are up to 6 surfaces connected to a single vertex, however this is a mathematical certainty for a heightfield where all of the triangles are oriented the same direction. As well, based on the algorithm we use to implement the above equation each vertex will have exactly 6 surfaces connected to it. I do this as a simplification to help in the computation of vertex normals and to speed up the processing. Now that we've got the equations defined, let's take a closer look at each one, see how they apply to our heightfield, and see if we can't make the equations more efficient.

Computation of surface normals
The equation above for surface normals requires first and foremost two vectors A and B which lie in the plane of the surface. For the purpose of this article we can assume we're using a heightfield with vertices laid out as in the diagram to the right.

In the diagram to the right each quad is made up of two triangles, one in the lower left of the quad and one in the top right. For the sake of this article you can assume that we're using DirectX's default left-handed coordinate system and that triangles are oriented positively in a clock-wise ordering. You can also assume that our heightfield lies in the xz-plane, and that y is the coordinate of height. From the diagram you can see that quad 0 is made up of triangles 0-3-1 and 1-3-4. Mathematically speaking a plane or surface can be computed from three points using the cross-product of the vectors formed by two of the edges of the triangle.

Equation for a Normal Vector from a Left-Handed Triangle:

A = P1 – P0; B = P2 – P0
N = A Χ B

As an example, given triangles 031, 134, and height values as a function of (x,z) pairs you can compute the normals to the triangles with the following equations:

// Triangle 0 Normal
A = 3 – 0; A = ( 0, h3, 1 ) – ( 0, h0, 0 )
B = 1 – 0; B = ( 1, h1, 0 ) – ( 0, h0, 0 )
N0 = [ ( 0, h3, 1 ) – ( 0, h0, 0 ) ] Χ [( 1, h1, 0 ) – ( 0, h0, 0 ) ]
Optionally normalize N0
// Triangle 1 Normal
A = 3 – 1; A = ( 0, h3, 1 ) – ( 1, h1, 0 )
B = 4 – 1; B = ( 1, h4, 1 ) – ( 1, h1, 0 )
N1 = [ ( 0, h3, 1 ) – ( 1, h1, 0 ) ] Χ [ ( 1, h4, 1 ) – ( 1, h1, 0 ) ]
Optionally normalize N1

The above equations illustrate how surface normals can be computed for a single quad of the heightfield, however, this is a special case and we need a general form of the equation if we're going to apply it to the terrain as a whole. In general, we'd like an equation that uses i, and j notation, instead of specific x,z pairs to allow us to iterate over the entire terrain and maintain the same equation. So given that i is an index in the range of 0 to numVertsWide – 1 and j is an index in the range of 0 to numVertsDeep – 1, let's see if we can't establish a general equation that works in code such as:

for( j = 0; j < numVertsDeep; j++ )
  for( i = 0; i < numVertsWide; i++ )
  {
    ...
  }

For the triangle in the bottom left, hereafter referred to as Triangle 0, let P0 be the point in the bottom left, P1 be the point in the top left and P2 be the point in the bottom right. If P0, P1, and P2 were a function of the current [i, j] pair, the points could be represented as follows:

P0 = [ i,     h0,  j ]
P1 = [ i,     h1,  j+1 ]
P2 = [ i+1, h2,  j ]

// Triangle 0
A = P1 – P0
B = P2 – P0

Now that we can represent any point on the heightfield as a function of ‘i' and ‘j', let's re-examine the equation for the normal. Even though we are rephrasing the equations with (i, j) coordinates, the pairs will drop out leaving us with simpler calculations at the end.

// Compute the two vectors used in the cross product
A = P1 – P0; A = (i, h1, j+1) – (i, h0, j); A = ( 0, h1 - h0, 1 )
B = P2 – P0; B = (i+1, h2, j) – (i, h0, j); B = ( 1, h2 - h0, 0 )

// Compute the normal to Triangle 0, using the generic form i,j and the full definition of cross-product
N0 = A Χ B
N0 = [ AyBz − AzBy, AzBx – AxBz, AxBy – AyBx ]

N0 = [ ( 0, h1 - h0, 1 ) ] Χ ( 1, h2 - h0, 0 ) ]
N0 = [ 0( h1 - h0 ) – 1( h2 - h0 ), 1(1) – 0(0), 0( h1 - h0 ) – 1( h1 - h0 ) ]
N0 = [ h0 - h2, 1.0f, h0 - h1 ]

In the equation above the cross product of two orthogonal vectors can be simplified dramatically. The calculation of a cross product normally requires 6 multiplications and 3 subtractions, however once we've computed it to the end we see that many of the terms drop out and that only 2 subtractions are required. Before you reach for the D3DXVec3Cross() function remember that DirectX knows nothing about the nature of your vectors. In this case, based on the orthogonal nature of our vectors it is far faster to compute the cross product ourselves than to let DirectX do it for us. Now that we've computed the cross product of Triangle 0 in general form, let's quickly compute the cross product of Triangle 1 using generic form. Once we're done we'll have the two equations necessary to iterate over the entire heightfield and compute surface normals.

// The 3 points used to compute the top right triangle
P0 = [ i+1, j ]
P1 = [ i, j+1 ]
P2 = [ i+1, j+1 ]

// Compute the two vectors used in the cross product
A = P1 – P0; A = (     i, h1, j+1 ) – ( i+1, h0, j ); A = ( -1, h1 - h0, 1 )
B = P2 – P0; B = ( i+1, h2, j+1 ) – ( i+1, h0, j ); B = ( 0, h2 - h0, 1 )

N0 = [ AyBz − AzBy, AzBx – AxBz, AxBy – AyBx ]
N0 = [( -1, h1 - h0, 1 ) ] Χ ( 0, h2 - h0, 1 ) ]
N0 = [ 1( h1 - h0 ) – 1( h2 - h0 ), 1(0) – (-1)(1), (-1)( h2 - h0 ) – (0)( h1 - h0) ]
N0 = [ h1 - h2, 1.0f, h0 - h2 ]

As with the equation for Triangle 0, Triangle 1 also simplifies down to 2 subtractions as one might expect. You are now armed with everything you need to compute the surface normals for every triangle on the heightfield. Simply iterate over the quads, and for each quad compute the two surface normals using the equations above. Before we look at the code for computing the surface normals, lets quickly evaluate the range over which we will iterate.

A heightfield will have numVertsWide vertices in the x-direction and numVertsDeep in the z-direction. You could iterate from 0 to numVertsWide – 1 and from 0 to numVertsDeep -1, and that would give you every triangle on the heightfield. Although this would initially seem like a good idea, there is a better way. Instead of computing those ranges, imagine that there is a heightfield that is one size larger than your own in all directions. That is, if you were to place your heightfield inside of this larger heightfield, as indicated in the diagram below, you would see that there is one extra row and column of quads on each edge. The purpose of those extra quads will become clear when we get into vertex normal calculations. For now, just make sure that you allocate enough space to store the extra surface normals and initialize them to the up ( 0.0, 1.0, 0.0 ) vector. Listed below is the code from my demo for computing the surface normals. You will see that although I iterate over a larger range of surfaces I do not actually calculate the normals for those extra quads. I want those quad's surfaces to remain with their normals facing the Up direction. I still iterate through the range, however, because I want my indices to update so that I store my surface normals in the correct place in the array. Follow along:

voidHeightfield::ComputeSurfaceNormals( void )
{
  INT normalIndex = 0;

  // Iterate through the z components
  for( int z = -1; z <= m_NumVertsDeep - 1; z++ )
  {
    // Iterate through the x components
    for( int x = -1; x <= m_NumVertsWide - 1; x++ )
    {
      // Does the current x,z lie on the terrain, or just outside?
      if( !IsValidCoord( x, z )   || !IsValidCoord( x+1, z ) ||
          !IsValidCoord( x, z+1 ) || !IsValidCoord( x+1, z+1 ) )
      {
        normalIndex += 2;
        continue;
      }

      // Set/Compute the normal for triangle 0
      float height0 = GetHeight( x,   z );
      float height1 = GetHeight( x, z+1 );
      float height2 = GetHeight( x+1, z );
      D3DXVECTOR3 normal1( height0 - height2, 1.0f, height0 - height1 );

      // Set/Compute the normal for triangle 1
      height0 = height2;
      height2 = GetHeight( x+1, z+1 );
      D3DXVECTOR3 normal2( height1 - height2, 1.0f, height0 - height2 );

      // Set the surface normals in the array
      m_pSurfaceNormals[normalIndex++] = normal1;
      m_pSurfaceNormals[normalIndex++] = normal2;
    }
  }
}

Computation of vertex normals
The second stage of Gouraud's normal algorithm is calculation of vertex normals. The second stage is far easier than the first stage and is less based on mathematics than on logic. Computing the vertex normals is as simple as taking the adjacent surfaces to the current vertex, adding them up, and dividing by 6. The hard part here is not the division, but instead determining which of the surface normals we computed above are adjacent to the current triangle. The equation for computing vertex normals is:

Vertex Normal: Nnwe = Σ (Nsi) / 6.0f | i = 0…5

In order to determine which of the surfaces are adjacent to the current vertex let's once again take a look at the image to the right. In the drawing to the right each of the vertices has a different number of surrounding surfaces. Vertex 0 for example has only one surface adjacent to it. Vertex 1 has three surfaces adjacent, two in the top left, and one in the top right. Vertex 2 has two adjacent surfaces to its upper left. And vertex 4 has six adjacent surfaces – three above it, and three below it. In order to create an algorithm in this example that accurately computes the vertex normal of every vertex we would have to take into account each of the special cases. Is the vertex on an edge – top, left, bottom, right? Is it on a corner – top left, top right, bottom left, or bottom right? Such an algorithm would first have to compute conditional statements to determine which of the above statements is true. Once the condition of the vertex was known it would likely progress through a series of if-statements where certain adjacent surfaces would be summed. There's got to be a cleaner way!

The primary reason we have this problem is because the heightfield has edges. If the heightfield had no edges then all vertices would have 6 surface normals such as the middle vertex to the image to the right. In which case, it might be possible to create an algorithm that reliably determined the six surrounding surfaces for any given vertex. So let's remove the edges from the heightfield! If you recall from the section above on computing the surface normals I illustrated how you should add an extra row and column of quads to each edge. The quads would not be updated but would merely have the UP vector as its normal. The reason we did that is so we could change the previous image to the one just below.

In the image to the right we no longer have edges to our heightfield because of the additional sets of quads. Every vertex on our heightfield now has 6 surrounding surfaces. Based on the direction of our triangles, each vertex has two surfaces to the top left, one to the top right, one to the bottom left, and two to the bottom right. Now imagine there was a small square around any given vertex. The square would be a part of 6 surfaces. If you slide the square to the left, right, up, or down, it is always located in six triangles. Therefore, we will call this algorithm The Sliding Six Algorithm.

The idea behind the algorithm is to put the rectangle in the bottom left corner, determine which surfaces it's part of, and then add those surfaces to a list to be summed for the current vertex index. Let's try it with a few vertices to help you get the hang of it.

Vertex 0
The triangle to the bottom left is triangle 1
The triangle just below it is triangle 2
The triangle to the bottom right is triangle 3
The triangle to the top right is triangle 8
The triangle just above it is triangle 9
The triangle to the top right is triangle 10.

Vertex 1
The triangles on the bottom from left to right are 3, 4, and 5
The triangles on the top from left to right are 10, 11, 12

Vertex 2
The triangles on the bottom from left to right are 5, 6, 7
The triangles on the top from left to right are 12, 13, 14

I wanted to get three vertices above because I wanted you to see the pattern evolve. Each row of triangles is sequential from one vertex to the next. Vertex 0 had surfaces 1, 2, and 3, where vertex 1 had surfaces 3, 4, 5 and so on. This would have continued all the way across the heightfield. The same was true for the top row. So if you could get a starting index for each vertex, we'll call it the triangleIndex, then we'd just need to start with that index and move along two to the right, and the corresponding three to the top and we'd be done. The only question would be how to get the surface which is immediately above the triangle index. That can be obtained from the following equation:

rowOffset = 2( numVertsWide + 1 ) - 1;

Here's the algorithm as implemented in the attached demo:

void Heightfield::ComputeVertexNormals( void )
{
  D3DXVECTOR3 vec[6];
  D3DXVECTOR3 vertexNormal;
  INT triIndex    = 1;
  INT indexPlusOne = 2;
  INT indexPlusTwo = 3;
  INT vertIndex   = 0;

  // Computes the offset based on the number of triangles across
  INT rowOffset = ( ( m_NumVertsWide + 1 ) * 2 ) - 1;
  for( INT z = 0; z < m_NumVertsDeep; z++ )
  {
    for( INT x = 0; x < m_NumVertsWide; x++ )
    {
      indexPlusOne = triIndex + 1;
      indexPlusTwo = triIndex + 2;

      // Get the three triangles below the vertex
      vec[0] = m_pSurfaceNormals[ triIndex ];
      vec[1] = m_pSurfaceNormals[ indexPlusOne ];
      vec[2] = m_pSurfaceNormals[ indexPlusTwo ];

      // Get the three triangles above the vertex
      vec[3] = m_pSurfaceNormals[ rowOffset + triIndex ];
      vec[4] = m_pSurfaceNormals[ rowOffset + indexPlusOne ];
      vec[5] = m_pSurfaceNormals[ rowOffset + indexPlusTwo ];

      // Sum the vectors and then divide by 6 to average them
      vertexNormal = vec[0] + vec[1] + vec[2] + vec[3] + vec[4] + vec[5];
      vertexNormal /= 6.0f;

      m_pVertexNormals[vertIndex] = vertexNormal;

      triIndex += 2;
      vertIndex++;
    }

    triIndex += 2;
  }
}




Method 2


Contents
  Introduction
  Quick History of Normal Calculations of Polyhedra
  Normal Calculations for Heightfields
  Method 2
  Analysis of the performance of each approach
  Future research

  Printable version
  Discuss this article