vectorrenderingraytracingnormalsheightmap

Generating normals from a voxel heightmap within a regular grid


I am making a voxel game, where each level is represented by a heightmap.

I am implementing precomputed raytraced ambient occlusion for it, by computing a 3D regular grid of occlusion terms (one per point in the world, ranging from 0 to the max X and Z values, and from 0 to the max voxel height for Y). I'm not baking the occlusion terms into the vertices so that all other objects in the scene can read from this 3D texture.

Here's what an example scene looks like: enter image description here

In order to compute am ambient occlusion term for each vertex, I need a surface normal at each point, to cast rays from the hemisphere defined by that normal. In some cases, where a point on the regular grid is either below or above the heightmap, there will be no normal; but that's just an edge case.

At the moment, I'm struggling with creating a function to compute this normal. Here it is at the moment:

// The inputs are x, y (height), and z.

const byte left_x = (x == 0) ? 0 : (x - 1), top_z = (z == 0) ? 0 : (z - 1);
    
#define SIGN_DIFF(a, b) sign_of((GLfloat) (a) - (GLfloat) (b))

/*
| a | b |
| c | d |
*/

const byte
    a = sample_heightmap(heightmap, left_x, top_z),
    b = sample_heightmap(heightmap, x, top_z),
    c = sample_heightmap(heightmap, left_x, z),
    d = sample_heightmap(heightmap, x, z);

vec3 normal = {
    SIGN_DIFF(c, d),
    y == a || y == b || y == c || y == d,
    SIGN_DIFF(b, d)
};

normalize(normal);

#undef SIGN_DIFF

Here's how it works: First, I compute the sign difference between the current y and adjacent points b and c, and use those gradients as the initial x and z components of the normal. Then, if the y height equals any of the 4 sampled heights, the y-component is set to point upwards (i.e. it's set to 1); otherwise it points straight ahead (i.e. it's set to 0).

When visualizing the normals, you can see that most of them are correct (ignore the ones over the heightmap; I'm not as concerned about the incorrectness of those yet).

enter image description here

Other normals are not correct, though:

enter image description here

The bottom incorrect normal here is pointing in -x, +y, and 0 for z, for reference. The top normal is pointing in -x, and 0 for y and z.

For someone who has worked with raytracing in a voxel environment like this before, how did you solve this problem of finding the right normal on a heightmap? And with that, do you see what is wrong with my normal-calculating algorithm?


Solution

  • I managed to solve this about a week ago like this:

    These are some wrapper types and helper functions used:

    typedef unsigned char byte;
    typedef signed char signed_byte;
    typedef float vec3[3];
    
    signed_byte sign_between_bytes(const byte a, const byte b) {
        if (a > b) return 1;
        else if (a < b) return -1;
        else return 0;
    }
    
    signed_byte clamp_signed_byte_to_directional_range(const signed_byte x) {
        if (x > 1) return 1;
        else if (x < -1) return -1;
        else return x;
    }
    
    byte sample_map_point(const byte* const map, const byte x, const byte z, const byte map_width) {
        return map[z * map_width + x];
    }
    

    Then, I defined this enum, to represent possible normal states:

    typedef enum {OnMap, AboveMap, BelowMap} NormalLocationStatus;
    

    The normal can only be valid when it's on the map. Otherwise, there is no normal, because that point in space isn't touching the map as a surface.

    Next, I defined the normal getter like this:

    NormalLocationStatus get_normal_at(
        const byte* const heightmap, const byte map_width,
        const byte x, const byte y, const byte z, vec3 normal) {
    
        const byte left_x = (x == 0) ? 0 : (x - 1), top_z = (z == 0) ? 0 : (z - 1);
    
        // 't' = top, 'b' = bottom, 'l' = left, and 'r' = right
        const struct {const signed_byte tl, tr, bl, br;} diffs = {
            sign_between_bytes(sample_map_point(heightmap, left_x, top_z, map_width), y),
            sign_between_bytes(sample_map_point(heightmap, x, top_z, map_width), y),
            sign_between_bytes(sample_map_point(heightmap, left_x, z, map_width), y),
            sign_between_bytes(sample_map_point(heightmap, x, z, map_width), y)
        };
    
        //////////
    
        normal[0] = clamp_signed_byte_to_directional_range((diffs.bl - diffs.br) + (diffs.tl - diffs.tr));
        normal[1] = !diffs.tl || !diffs.tr || !diffs.bl || !diffs.br;
        normal[2] = clamp_signed_byte_to_directional_range((diffs.tl - diffs.bl) + (diffs.tr - diffs.br));
    
        if (normal[0] == 0.0f && normal[1] == 0.0f && normal[2] == 0.0f)
            return (diffs.tl == 1) ? BelowMap : AboveMap;
    
        normalize(normal);
        return OnMap;
    }
    

    I am still not completely sure how this works. I managed to make this function by writing out a lookup table mapping possible flow states to unnormalized normals (this was around 81 entries, I think, since 3 states per flow tile, and 3^4 = 81), and then writing another helper function that told me what percentage of locations on the heightmap yielded the correct normal with my non-lookup function. Iterating more and more until I reached a very high correct percentage lead me to where I am now.

    I hope this helped you if you had the same issue as me!