Help on view culling algo

Hi guys, I have a little maths problem and I’d like some help.

I decided to have a go at LWJGL, and after some basic setting up stuff, I’m now trying to implement a view culling algorithm. I have some OpenGL tutorials in my hands, but no code, so this is what I tried to do:


private final static int RIGHT_PLANE = 0;
private final static int LEFT_PLANE = 4;
private final static int TOP_PLANE = 8;
private final static int BOTTOM_PLANE = 12;
private final static int NEAR_PLANE = 16;
private final static int FAR_PLANE = 20;
      
private double[] frustrum;

// Computes the frustrum plane equations      
public void computeFrustrum() {
      gl.loadIdentity();
            
      // Camera Transformation
      gl.translatef(0.0f, 0.0f, camZ);
      gl.rotatef(camRotX, 1.0f, 0.0f, 0.0f);
            
      // Interest Transformation
      gl.rotatef(intRotY, 0.0f, 1.0f, 0.0f);
      gl.translatef(intX, intY, intZ);
            
      // Get ModelView Matrix
      BufferCache.loadMatrix(GL.MODELVIEW_MATRIX);
      matrix.set(BufferCache.matrix);
            
      // Multiply by Projection Matrix
      matrix.mul(EngineGraph.PROJECT_MATRIX);
            
      // Right
      frustrum[0] = matrix.m03-matrix.m00;
      frustrum[1] = matrix.m13-matrix.m10;
      frustrum[2] = matrix.m23-matrix.m20;
      frustrum[3] = matrix.m33-matrix.m30;
            
      // Left
      frustrum[4] = matrix.m03+matrix.m00;
      frustrum[5] = matrix.m13+matrix.m10;
      frustrum[6] = matrix.m23+matrix.m20;
      frustrum[7] = matrix.m33+matrix.m30;
            
      // Top
      frustrum[8] = matrix.m03-matrix.m01;
      frustrum[9] = matrix.m13-matrix.m11;
      frustrum[10] = matrix.m23-matrix.m21;
      frustrum[11] = matrix.m33-matrix.m31;
            
      // Bottom
      frustrum[12] = matrix.m01+matrix.m03;
      frustrum[13] = matrix.m11+matrix.m13;
      frustrum[14] = matrix.m21+matrix.m23;
      frustrum[15] = matrix.m31+matrix.m33;
            
      // Near
      frustrum[16] = matrix.m02-matrix.m03;
      frustrum[17] = matrix.m12-matrix.m13;
      frustrum[18] = matrix.m22-matrix.m23;
      frustrum[19] = matrix.m32-matrix.m33;
      
      // Far
      frustrum[20] = matrix.m03-matrix.m02;
      frustrum[21] = matrix.m13-matrix.m12;
      frustrum[22] = matrix.m23-matrix.m22;
      frustrum[23] = matrix.m33-matrix.m32;
}

// Returns true when point is within the view      
public boolean isVisible(float x, float y, float z) {
      if ( visibleInPlane(RIGHT_PLANE, x, y, z)
      && visibleInPlane(LEFT_PLANE, x, y, z)
      && visibleInPlane(TOP_PLANE, x, y, z)
      && visibleInPlane(BOTTOM_PLANE, x, y, z)
      && !visibleInPlane(NEAR_PLANE, x, y, z)
      && visibleInPlane(FAR_PLANE, x, y, z) )
            return true;
      else
            return false;
}

// Returns true if the point is "inside" the plane      
private boolean visibleInPlane(int plane, float x, float y, float z) {
      double xc = frustrum[plane++] * x;
      double yc = frustrum[plane++] * y;
      double zc = frustrum[plane++] * z;
      double wc = frustrum[plane];
      
      if ( -wc < xc && xc < wc && -wc < yc && yc < wc && -wc < zc && zc < wc )
            return true;
      else
            return false;
}

It works pretty nice, but I don’t know why I have to put that ! in front of the near plane comparison (it works that way, I just can’t understand why it doesn’t without it ???). And there is a problem with the bottom plane. When the view looks straight down the y axis or straight forward, down the z axis, it works correctly, but at other angles it fails. It returns false even when the point is just below the middle of the screen, or lower depending on the view angle.

I’m not very good at math, but I think I implemented what I read in the tutorials properly. I just can’t see where the problem is. Any hints people?

  • Any comments on the algo itself will be appreciated (I’m not even sure if it is a “decent” way to do it). And of course suggestions for making it better.

Spasi

I think most likely you just have your normal for the near plane reversed. It is probably in the same direction as your far plane, thus needing the ! for your inside/outside test. Am I correct in assuming each plane is defined in the form Ax + By + Cz + D = 0?
plane[i] = A
plane[i+1] = B… etc

  • Tristan

heh, view frustum clipping is somewhere on my todo list for my next project, but I most of your code looks good to me.

However, I’m wondering what happens in your code when an object is sufficiantly large enough to cover the entire screen? Individual points may be outside the view, but the object is still visible. I was thinking of doing a 3d version of the Cohen-Sutherland algorithm, which should be nice and fast…

[rant]
Oh, how I hated this kind of problems when I did D3D earlier… Holes in the terrain, object popping in and out … why why why… ?

I’m so glad there Java3D around now…

[/rant]

I solve my culling problems last week and I thought I should share my solution. It’s quite different from what I started with and it works wonderfully! I’ll just paste the whole code and hope it might help someone here…

I’ll post it in two messages, because it’s too long for YABB…


public final class Frustum {

    private final double[] frustum;

    private final Point3d corner;
    private final Point3d lower;
    private final Point3d upper;

    public Frustum() {
        frustum = new double[24];

        corner = new Point3d();
        lower = new Point3d();
        upper = new Point3d();
    }

    public final void calculate(final Matrix4d matrix) {
        // Left
        frustum[0] = matrix.m03 + matrix.m00;
        frustum[1] = matrix.m13 + matrix.m10;
        frustum[2] = matrix.m23 + matrix.m20;
        frustum[3] = matrix.m33 + matrix.m30;

        double norm = Math.sqrt(frustum[0] * frustum[0] + frustum[1] * frustum[1] + frustum[2] * frustum[2]);
        frustum[0] /= norm;
        frustum[1] /= norm;
        frustum[2] /= norm;
        frustum[3] /= norm;

        // Right
        frustum[4] = matrix.m03 - matrix.m00;
        frustum[5] = matrix.m13 - matrix.m10;
        frustum[6] = matrix.m23 - matrix.m20;
        frustum[7] = matrix.m33 - matrix.m30;

        norm = Math.sqrt(frustum[4] * frustum[4] + frustum[5] * frustum[5] + frustum[6] * frustum[6]);
        frustum[4] /= norm;
        frustum[5] /= norm;
        frustum[6] /= norm;
        frustum[7] /= norm;

        // Top
        frustum[8] = matrix.m03 - matrix.m01;
        frustum[9] = matrix.m13 - matrix.m11;
        frustum[10] = matrix.m23 - matrix.m21;
        frustum[11] = matrix.m33 - matrix.m31;

        norm = Math.sqrt(frustum[8] * frustum[8] + frustum[9] * frustum[9] + frustum[10] * frustum[10]);
        frustum[8] /= norm;
        frustum[9] /= norm;
        frustum[10] /= norm;
        frustum[11] /= norm;

        // Bottom
        frustum[12] = matrix.m01 + matrix.m03;
        frustum[13] = matrix.m11 + matrix.m13;
        frustum[14] = matrix.m21 + matrix.m23;
        frustum[15] = matrix.m31 + matrix.m33;

        norm = Math.sqrt(frustum[12] * frustum[12] + frustum[13] * frustum[13] + frustum[14] * frustum[14]);
        frustum[12] /= norm;
        frustum[13] /= norm;
        frustum[14] /= norm;
        frustum[15] /= norm;

        // Near
        frustum[16] = matrix.m02 + matrix.m03;
        frustum[17] = matrix.m12 + matrix.m13;
        frustum[18] = matrix.m22 + matrix.m23;
        frustum[19] = matrix.m32 + matrix.m33;

        norm = Math.sqrt(frustum[16] * frustum[16] + frustum[17] * frustum[17] + frustum[18] * frustum[18]);
        frustum[16] /= norm;
        frustum[17] /= norm;
        frustum[18] /= norm;
        frustum[19] /= norm;

        // Far
        frustum[20] = matrix.m03 - matrix.m02;
        frustum[21] = matrix.m13 - matrix.m12;
        frustum[22] = matrix.m23 - matrix.m22;
        frustum[23] = matrix.m33 - matrix.m32;

        norm = Math.sqrt(frustum[20] * frustum[20] + frustum[21] * frustum[21] + frustum[22] * frustum[22]);
        frustum[20] /= norm;
        frustum[21] /= norm;
        frustum[22] /= norm;
        frustum[23] /= norm;

        lower.x = lower.y = lower.z = Double.MAX_VALUE;
        upper.x = upper.y = upper.z = Double.MIN_VALUE;

        calcPoint(0, 8, 16);
        calcPoint(0, 12, 16);
        calcPoint(4, 12, 16);
        calcPoint(4, 8, 16);
        calcPoint(4, 8, 20);
        calcPoint(0, 8, 20);
        calcPoint(0, 12, 20);
        calcPoint(4, 12, 20);
    }

    private void calcPoint(final int a, final int b, final int c) {
        MarathonMath.plane3Cut(frustum, a, b, c, corner);

        if ( corner.x < lower.x )
            lower.x = corner.x;
        if ( corner.x > upper.x )
            upper.x = corner.x;

        if ( corner.y < lower.y )
            lower.y = corner.y;
        if ( corner.y > upper.y )
            upper.y = corner.y;

        if ( corner.z < lower.z )
            lower.z = corner.z;
        if ( corner.z > upper.z )
            upper.z = corner.z;
    }


    public final boolean isVisible(final float x, final float y, final float z) {
        int plane = 0;

        // LEFT PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane++] < 0 )
            return false;

        // RIGHT PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane++] < 0 )
            return false;

        // TOP PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane++] < 0 )
            return false;

        // BOTTOM PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane++] < 0 )
            return false;

        // NEAR PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane++] < 0 )
            return false;

        // FAR PLANE
        if ( frustum[plane++] * x + frustum[plane++] * y + frustum[plane++] * z + frustum[plane] < 0 )
            return false;

        return true;
    }

    public final boolean isVisible(final BoundingSphere sphere) {
        int plane = 0;
        final double invRadius = -sphere.radius;

        // LEFT PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane++] < invRadius ) {
            return false;
        }

        // RIGHT PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane++] < invRadius ) {
            return false;
        }

        // TOP PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane++] < invRadius ) {
            return false;
        }

        // BOTTOM PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane++] < invRadius ) {
            return false;
        }

        // NEAR PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane++] < invRadius ) {
            return false;
        }

        // FAR PLANE
        if ( frustum[plane++] * sphere.x + frustum[plane++] * sphere.y + frustum[plane++] * sphere.z + frustum[plane] < invRadius ) {
            return false;
        }

        return true;
    }

    public final int isVisible(final BoundingBox box) {
        // First check intersection with the frustum's bounding box
        if ( box.upper.x > lower.x && box.lower.x < upper.x && box.upper.y > lower.y && box.lower.y < upper.y && box.upper.z > lower.z && box.lower.z < upper.z ) {
            double x, y, z, w;
            boolean inside = true;
            // For each of the six frustum planes
            for ( int index = 0; index <= 20; index += 4 ) {
                if ( inside ) { // There's still a chance the box is entirely inside the frustum
                    boolean pointIn = false;

                    // Check if any of the boxes' points is visible
                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.lower.y;
                    z = frustum[index + 2] * box.lower.z;
                    w = frustum[index + 3];
                    if ( x + y + z + w > 0 )
                        pointIn = true;

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 ) { // The point is visible
                        if ( !pointIn ) { // If previous points were not visible
                            inside = false; // No chance for the box to be entirely inside frustum, so continue
                            continue;
                        }
                    } else { // The point is not visible
                        inside = false;
                        if ( pointIn ) // If there was a previous visible point
                            continue;
                    }

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.upper.y;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.lower.y;
                    z = frustum[index + 2] * box.upper.z;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.upper.y;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 ) {
                        if ( !pointIn ) {
                            inside = false;
                            continue;
                        }
                    } else {
                        inside = false;
                        if ( pointIn )
                            continue;
                    }

                    if ( !pointIn )
                        return 0;
                } else { // If we get here, the box is not entirely inside the frustum
                    // Check if any of the boxes' points is visible
                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.lower.y;
                    z = frustum[index + 2] * box.lower.z;
                    w = frustum[index + 3];
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.upper.y;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.lower.y;
                    z = frustum[index + 2] * box.upper.z;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.lower.x;
                    y = frustum[index + 1] * box.upper.y;
                    if ( x + y + z + w > 0 )
                        continue;

                    x = frustum[index] * box.upper.x;
                    if ( x + y + z + w > 0 )
                        continue;

                    return 0;
                }
            }
            return inside ? 2 : 1;
        } else // Does not intersect the frustum's bounding box
            return 0;
    }
}

The method MarathonMath.plane3Cut(…) is this:


public static void plane3Cut(final double[] planes, final int a, final int b, final int c, final Point3d p) {
        System.arraycopy(planes, a, system[0], 0, 3);
        System.arraycopy(planes, b, system[1], 0, 3);
        System.arraycopy(planes, c, system[2], 0, 3);

        system[0][3] = -planes[a + 3];
        system[1][3] = -planes[b + 3];
        system[2][3] = -planes[c + 3];

        double l;
        for ( int i = 0; i < 3; i++ ) {
            for ( int j = 0; j < 3; j++ ) {
                if ( i == j )
                    continue;

                l = system[j][i] / system[i][i];
                for ( int k = i + 1; k < 4; k++ ) {
                    system[j][k] = system[j][k] - l * system[i][k];
                }
            }
        }

        p.x = system[0][3] / system[0][0];
        p.y = system[1][3] / system[1][1];
        p.z = system[2][3] / system[2][2];
}

I’ll try to explain some of it:

  • I call calculate(…) whenever the view changes. I pass MODELVIEW_MATRIX * PROJECT_MATRIX. From this matrix I get the plane equations and normalize the coefficients (it’s in the form Ax + By + Cz + D = 0). Normalization is necessary for the BoundingSphere culling to work.

  • From there I find the intersections of three planes for each point of the frustum (calling calcPoint with the starting indices of the proper planes as parameters). calcPoint(…) calls plane3cut(…) which returns to the intersection of the three planes. plane3cut(…) is an implementation of the Gauss-Jordan method, which solves multiple equations with multiple variables. In this case we have three plane equations and want to find the values of three variables (x,y,z). From the points that I find, I just hold the maximum and minimum values, that define a “BoundingBox” for the frustum.

  • I believe point and sphere culling are simple. I just find the distances from each plane.

  • The BoundingBox culling is a little mess. It’s the method I’m using most in my engine, that’s why it’s more improved than the others. It doesn’t return boolean, but rather 0 for not visible at all, 1 for visible and 2 for completely inside the view. I first try an intersection with the frustum bounding box (which saves a lot of calculations) and then I go plane by plane. The mess begins here because I have split the method in two ( if ( inside ) … else … ). In the first part there’s still a chance for the BoundingBox to be completely inside the frustum and in the second there’s not. I did it this way because I believe it saves me a lot of calculation.

  • I could make the sphere culling to work as the box (0,1,2), I just don’t need it yet. It should be really easy.

That’s it! If anyone takes the time to study it, I would really appreciate any comments.

Spasi

I just implemented frustum clipping in my software renderer and I had the same questions myself when I started. It actually turned out to be fairly easy using a 3D version of the Sutherland-Hodgman re-entrant polygon clipper. When clipping against a canonical view volume you can really optimize it. Culling is next on my list. I think gamasutra has a nice article on bounding box and sphere frustum culling in their archives.