Flat, 3D triangle, made out of voxels

Question

I have a problem that I can't seem to get a working algorithm for, I've been trying to days and get so close but yet so far.

I want to draw a triangle defined by 3 points (p0, p1, p2). This triangle can be any shape, size, and orientation. The triangle must also be filled inside.

Here's a few things I've tried and why they've failed:

1

• Drawing lines along the triangle from side to side
  • Failed because the triangle would have holes and would not be flat due to the awkwardness of drawing lines across the angled surface with changing locations

2

• Iterate for an area and test if the point falls past the plane parallel to the triangle and 3 other planes projected onto the XY, ZY, and XZ plane that cover the area of the triangle
  • Failed because for certain triangles (that have very close sides) there would be unpredictable results, e.g. voxels floating around not connected to anything

3

• Iterate for an area along the sides of the triangle (line algorithm) and test to see if a point goes past a parallel plane
  • Failed because drawing a line from p0 to p1 is not the same as a line from p1 to p0 and any attempt to rearrange either doesn't help, or causes more problems. Asymmetry is the problem with this one.

This is all with the intent of making polygons and flat surfaces. 3 has given me the most success and makes accurate triangles, but when I try to connect these together everything falls apart and I get issues with things not connecting, asymmetry, etc. I believe 3 will work with some tweaking but I'm just worn out from trying to make this work for so long and need help.

There's a lot of small details in my algorithms that aren't really relevant so I left them out. For number 3 it might be a problem with my implementation and not the algorithm itself. If you want code I'll try and clean it up enough to be understandable, it will take me a few minutes though. But I'm looking for algorithms that are known to work. I can't seem to find any voxel shape making algorithms anywhere, I've been doing everything from scratch.

EDIT:

Here's the third attempt. It's a mess, but I tried to clean it up.

// Point3i is a class I made, however the Vector3fs you'll see are from lwjgl

public void drawTriangle (Point3i r0, Point3i r1, Point3i r2)
{
    // Util is a class I made with some useful stuff inside

    // Starting values for iteration

    int sx = (int) Util.min(r0.x, r1.x, r2.x);
    int sy = (int) Util.min(r0.y, r1.y, r2.y);
    int sz = (int) Util.min(r0.z, r1.z, r2.z);

    // Ending values for iteration

    int ex = (int) Util.max(r0.x, r1.x, r2.x);
    int ey = (int) Util.max(r0.y, r1.y, r2.y);
    int ez = (int) Util.max(r0.z, r1.z, r2.z);

    // Side lengths

    float l0 = Util.distance(r0.x, r1.x, r0.y, r1.y, r0.z, r1.z);
    float l1 = Util.distance(r2.x, r1.x, r2.y, r1.y, r2.z, r1.z);
    float l2 = Util.distance(r0.x, r2.x, r0.y, r2.y, r0.z, r2.z);

    // Calculate the normal vector

    Vector3f nn = new Vector3f(r1.x - r0.x, r1.y - r0.y, r1.z - r0.z);
    Vector3f n   = new Vector3f(r2.x - r0.x, r2.y - r0.y, r2.z - r0.z);
    Vector3f.cross(nn, n, n);

    // Determines which direction we increment for

    int iz = n.z >= 0 ? 1 : -1;
    int iy = n.y >= 0 ? 1 : -1;
    int ix = n.x >= 0 ? 1 : -1;

    // Reorganize for the direction of iteration

    if (iz < 0) {
        int tmp = sz;
        sz = ez;
        ez = tmp;
    }
    if (iy < 0) {
        int tmp = sy;
        sy = ey;
        ey = tmp;
    }
    if (ix < 0) {
        int tmp = sx;
        sx = ex;
        ex = tmp;
    }

    // We're we want to iterate over the end vars so we change the value
    // by their incrementors/decrementors

    ex += ix;
    ey += iy;
    ez += iz;

    // Maximum length

    float lmax = Util.max(l0, l1, l2);

    // This is a class I made which manually iterates over a line, I already
    // know that this class is working

    GeneratorLine3d g0, g1, g2;

    // This is a vector for the longest side

    Vector3f v = new Vector3f();

    // make the generators

    if (lmax == l0) {
        v.x = r1.x - r0.x;
        v.y = r1.y - r0.y;
        v.z = r1.z - r0.z;

        g0 = new GeneratorLine3d(r0, r1);
        g1 = new GeneratorLine3d(r0, r2);
        g2 = new GeneratorLine3d(r2, r1);
    }
    else if (lmax == l1) {
        v.x = r1.x - r2.x;
        v.y = r1.y - r2.y;
        v.z = r1.z - r2.z;

        g0 = new GeneratorLine3d(r2, r1);
        g1 = new GeneratorLine3d(r2, r0);
        g2 = new GeneratorLine3d(r0, r1);
    }
    else {
        v.x = r2.x - r0.x;
        v.y = r2.y - r0.y;
        v.z = r2.z - r0.z;

        g0 = new GeneratorLine3d(r0, r2);
        g1 = new GeneratorLine3d(r0, r1);
        g2 = new GeneratorLine3d(r1, r2);
    }

    // Absolute values for the normal

    float anx = Math.abs(n.x);
    float any = Math.abs(n.y);
    float anz = Math.abs(n.z);

    int i, o;
    int si, so;
    int ii, io;
    int ei, eo;

    boolean maxx, maxy, maxz,
    midy, midz, midx,
    minx, miny, minz;

    maxx = maxy = maxz =
    midy = midz = midx =
    minx = miny = minz = false;

    // Absolute values for the longest side vector      

    float rnx = Math.abs(v.x);
    float rny = Math.abs(v.y);
    float rnz = Math.abs(v.z);

    int rmid = Util.max(rnx, rny, rnz);

    if (rmid == rnz) midz = true;
    else if (rmid == rny) midy = true;

    midx = !midz && !midy;

    // Determine the inner and outer loop directions

    if (midz) {
        if (any > anx) 
        {
            maxy = true;
            si = sy;
            ii = iy;
            ei = ey;
        }
        else {
            maxx = true;
            si = sx;
            ii = ix;
            ei = ex;
        }
    }
    else {
        if (anz > anx) {
            maxz = true;
            si = sz;
            ii = iz;
            ei = ez;
        }
        else {
            maxx = true;
            si = sx;
            ii = ix;
            ei = ex;
        }
    }

    if (!midz && !maxz) {
        minz = true;
        so = sz;
        eo = ez;
    }
    else if (!midy && !maxy) {
        miny = true;
        so = sy;
        eo = ey;
    }
    else {
        minx = true;
        so = sx;
        eo = ex;
    }

    // GeneratorLine3d is iterable
    Point3i p1;
    for (Point3i p0 : g0) {
        // Make sure the two 'mid' coordinate correspond for the area inside the triangle

        if (midz)
            do p1 = g1.hasNext() ? g1.next() : g2.next();
            while (p1.z != p0.z);
        else if (midy)
            do p1 = g1.hasNext() ? g1.next() : g2.next();
            while (p1.y != p0.y);
        else
            do p1 = g1.hasNext() ? g1.next() : g2.next();
            while (p1.x != p0.x);

        eo = (minx ? p0.x : miny ? p0.y : p0.z);
        so = (minx ? p1.x : miny ? p1.y : p1.z);
        io = eo - so >= 0 ? 1 : -1;

        for (o = so; o != eo; o += io) {
            for (i = si; i != ei; i += ii) {
                int x = maxx ? i : midx ? p0.x : o;
                int y = maxy ? i : midy ? p0.y : o;
                int z = maxz ? i : midz ? p0.z : o;

                // isPassing tests to see if a point goes past a plane
                // I know it's working, so no code

                // voxels is a member that is an arraylist of Point3i

                if (isPassing(x, y, z, r0, n.x, n.y, n.z)) {
                    voxels.add(new Point3i(x, y, z));
                    break;
                }
            }
        }
    }
}

Answer 1

You could use something like Besenham's line algorithm, but extended into three dimensions. The two main ideas we want to take from it are:

• rotate the initial line so its slope isn't too steep.
• for any given x value, find an integer value that is closest to the ideal y value.

Just as Bresenham's algorithm prevents gaps by performing an initial rotation, we'll avoid holes by performing two initial rotations.

1. Get the normal vector and point that represent the plane your triangle lies on. Hint: use the cross product of (line from p0 to p1) and (line from p0 to p2) for the vector, and use any of your corner points for the point.

2. You want the plane to be sufficiently not-steep, to avoid holes. You must satisfy these conditions:

   • -1 >= norm.x / norm.y >= 1
   • -1 >= norm.z / norm.y >= 1

   Rotate your normal vector and initial points 90 degrees about the x axis and 90 degrees about the z axis until these conditions are satisfied. I'm not sure how to do this in the fewest number of rotations, but I'm fairly sure you can satisfy these conditions for any plane.

3. Create a function f(x,z) which represents the plane your rotated triangle now lies on. It should return the Y value of any pair of X and Z values.
4. Project your triangle onto the XZ plane (i.e., set all the y values to 0), and use your favorite 2d triangle drawing algorithm to get a collection of x-and-z coordinates.
5. For each pixel value from step 4, pass the x and z values into your function f(x,z) from step 3. Round the result to the nearest integer, and store the x, y, and z values as a voxel somewhere.
6. If you performed any rotations in step 2, perform the opposite of those rotations in reverse order on your voxel collection.

Answer 2

Start with a function that checks for triangle/voxel intersection. Now you can scan a volume and find the voxels that intersect the triangle - these are the ones you're interested in. This is a lousy algorithm but is also a regression test for anything else you try. This test is easy to implement using SAT (separating axis theorem) and considering the triangle a degenerate volume (1 face, 3 edges) and considering the voxels symmetry (only 3 face normals).

I use octtrees, so my preferred method is to test a triangle against a large voxel and figure out which of the 8 child octants it intersects. Then use recursion on the intersected children until the desired level of subdivision is attained. Hint: at most 6 of the children can be intersected by the triangle and often fewer than that. This is tricky but will produce the same results as the first method but much quicker.

Rasterization in 3d is probably fastest, but IMHO is even harder to guarantee no holes in all cases. Again, use the first method for comparison.