How to determine a point is inside or outside a cube?

Question

Given a cube with 8 vertex in 3D space. How could I determine the myPoint is inside or outside the cube?

cube[0] = (x0, y0, z0);
cube[1] = (x1, y1, z1);
cube[2] = (x2, y2, z2);
cube[3] = (x3, y3, z3);
cube[4] = (x4, y4, z4);
cube[5] = (x5, y5, z5);
cube[6] = (x6, y6, z6);
cube[7] = (x7, y7, z7);

myPoint = (x, y, z);

I am trying to implement this data filter technique in 3D

Answer 1

1) Determine plane equation (ax+by+cz+d=0) of each face (F1, F2, ...)

1.1 Calc the normal vector of F1, (point this vector to outside of the object!!!):

n1 = (n1x, n1y, n1z)

1.2 Set n1 on coefficients in plane equation of F1

n1xPx+ n1yPy+ n1z*Pz - d = val

in java it looks like this:

private Vector3f normal;
private float d;
public PlaneEquation(Vector3f normal, Vector v1){
    this.normal = normal;
    this.d = normal.x*v1x + normal.y*v1y+ normal.z*v1z;
}
public float relativePosition(Vector3f p){
    return  normal.x*p.x + normal.y*p.y+ normal.z*p.z - d;
}

This equation is a function F1(Px,Py,Pz), such that P is a point to test.

1.3 calc the value of d, by set a coplanar point of F1 on eq. above. I choose v1, thus:

d = n1xv1x+ n1yv1y+ n1z*v1z

2) Perform a loop over all faces of the object, using the plane function method relativePosition(Vector3f p), with the test point coordinates, F1(Px, Py, Pz), F2(Px, Py, Pz), etc...

2.2) if all results F1(P), F2(P), ... are negative, so the point is INSIDE of the obect, else is OUT or ON. Or if some result are positive, then the point is OUT.

2.3) if no result is positive, but at least one is zero, so the point is ON the object.

    boolean OUT = false;
    boolean IN = false;
    boolean ON = false;
    for(PlaneEquation plane : planeList){
       if(plane.relativePosition(point) > 0){
           OUT = true;
           ON = false;
           break;
       }
       if(plane.relativePosition(point) == 0){
          ON = true;
       }
    }
    if(OUT == false && ON == false){
        IN = true;
    }

METHOD 2 - BARICENTRIC COORDINATES

import org.joml.Vector3f;

public class HexaedronEquation {
    // a is origin point
    private Vector3f a,b,c,d;//vertices of parallelepiped
    private Vector3f u,v,w; //direction vectors
    private Vector3f vXu; //v cross u
    private Vector3f wXu; //w croos u
    
    public HexaedronEquation(Vector3f a,Vector3f b,Vector3f c,Vector3f d ) {
        this.a = a;
        this.b = b;
        this.c = c;
        this.d = d;
        this.u = this.b.sub(a,new Vector3f());  // b - a
        this.v = this.c.sub(a, new Vector3f()); // c - a
        this.w = this.d.sub(a, new Vector3f()); // d - a
        this.vXu = this.v.cross(u, new Vector3f());
        this.wXu = this.w.cross(u, new Vector3f());
    }
    
    public float[] getParams(Vector3f p) {
        /* baricentric coordinates
         * solve the system
         * px = ax +t1*ux +t2*vx + t3*wx
         * py = ay +t1*uy +t2*vy + t3*wy
         * pz = az +t1*uz +t2*vz + t3*wz
         */
        float t1 = 0; 
        float t2 = 0; 
        float t3 = 0; 
        Vector3f s = p.sub(a, new Vector3f());
        Vector3f sXu = s.cross(u, new Vector3f());
        
        if(vXu.length() != 0 && w.length()!=0) {
            t3 = s.dot(vXu)/w.dot(vXu);
        }
        
        if(vXu.x != 0) {
            t2 = (sXu.x - t3*wXu.x)/vXu.x;
        }
        if(vXu.y != 0) {
            t2 = (sXu.y - t3*wXu.y)/vXu.y;
        }
        if(vXu.z != 0) {
            t2 = (sXu.z - t3*wXu.z)/vXu.z;
        }
        
        if(u.x != 0) {
            t1 = (s.x - v.x*t2 -w.x*t3)/u.x;
        }
        if(u.y != 0) {
            t1 = (s.y - v.y*t2 -w.y*t3)/u.y;
        }
        if(u.z != 0) {
            t1 = (s.z - v.z*t2 -w.z*t3)/u.z;
        }
        
        return new float[] {t1, t2, t3};
    }
    
    public boolean isInside(Vector3f p) {
        float[] params = getParams(p);
        if(params[0]< 0 || params[0] > 1) return false;
        if(params[1]< 0 || params[1] > 1) return false;
        if(params[2]< 0 || params[2] > 1) return false;
        
        if(params[0]>0 || params[0] <1) return true;
        if(params[1]>0 || params[1] <1) return true;
        if(params[2]>0 || params[2] <1) return true;
        return false;
    }
    
    public boolean isOn(Vector3f p) {
        float[] params = getParams(p);
        if(params[0]< 0 || params[0] > 1) return false;
        if(params[1]< 0 || params[1] > 1) return false;
        if(params[2]< 0 || params[2] > 1) return false;
        
        if(params[0]==0 || params[0] ==1) return true;
        if(params[1]==0 || params[1] ==1) return true;
        if(params[2]==0 || params[2] ==1) return true;
        return false;
    }
    
    public boolean isOnInside(Vector3f p) {
        float[] params = getParams(p);
        System.out.println("params: "+ getParamsString(params));
        if(params[0]< 0 || params[0] > 1) return false;
        if(params[1]< 0 || params[1] > 1) return false;
        if(params[2]< 0 || params[2] > 1) return false;
        return true;
    }
    
    public boolean isOut(Vector3f p) {
        float[] params = getParams(p);
        if(params[0]<0 || params[0] >1) return true;
        if(params[1]<0 || params[1] >1) return true;
        if(params[2]<0 || params[2] >1) return true;
        return false;
    }

Answer 2

Special Case #1 (Axis Aligned Cube):

As maxim1000's answer shows, you can simply check if X, Y, Z coordinates of the point in consideration lie in the minimum and maximum of the X, Y, Z coordinates of the Cube.

X_min <= X <= X_max and Y_min <= Y <= Y_max  and Z_min <= Z <= Z_max

If the aforementioned condition is satisfied then, the point lies inside the cube, otherwise it does not.

General Case (Oriented Cube):

There are two approaches to solve this. First one brings the point in the local coordinate system of the cube and apply the aforementioned Special Case. Second case works on the concept of the projection of the vector. First case is little more complex than the second as you need to calculate the rotation matrix which transforms the point from the world coordinate system into the Cube's local coordinate system. Consider the cube as show in the following figure.

For both the approaches we would need to derive some basic information from the cube representation. Let us fix origin in the local coordinate system of the cube at the bottom left corner of the cube; in this case, it is Point D. Now calculate the unit direction vectors in three dimensions and extent of the cube in those directions. It can be done as follows.

The X_local, Y_local, and Z_local are illustrated in the figure with Blue, Red, Green colors. And X_length, Y_length, and Z_length are the extents along the axes.

Now lets get back to solving the problem.

Approach #1: Bring the point in consideration, in the local coordinate system of the Cube. To do that, we need to estimate the rotation matrix. Rotation matrix in this case is 3 x 3 matrix with X_local, Y_local, and Z_local as columns.

Using the rotation Matrix R, you can bring the point in the local coordinate system and then apply the special case of axis aligned cube.

Approach #2:

Construct the direction vector from the cube center to the point in consideration and project it onto each local axis and check if the projection spans beyond the extent of the cube along that axis. If the projection lies inside the extent along each axis, then point is inside, otherwise it is outside of the cube.

The center of the cube is I, as show in the figure. The direction vector from the center of the cube to the point P is V. The projection of the vector V on X_local, Y_local, and Z_local can be calculated as follows.

Now the point P is inside the cube only if all of the following conditions are satisfied.

Here's the quick implementation in python for the approach #2.

import numpy as np

def inside_test(points , cube3d):
    """
    cube3d  =  numpy array of the shape (8,3) with coordinates in the clockwise order. first the bottom plane is considered then the top one.
    points = array of points with shape (N, 3).

    Returns the indices of the points array which are outside the cube3d
    """
    b1,b2,b3,b4,t1,t2,t3,t4 = cube3d

    dir1 = (t1-b1)
    size1 = np.linalg.norm(dir1)
    dir1 = dir1 / size1

    dir2 = (b2-b1)
    size2 = np.linalg.norm(dir2)
    dir2 = dir2 / size2

    dir3 = (b4-b1)
    size3 = np.linalg.norm(dir3)
    dir3 = dir3 / size3

    cube3d_center = (b1 + t3)/2.0

    dir_vec = points - cube3d_center

    res1 = np.where( (np.absolute(np.dot(dir_vec, dir1)) * 2) > size1 )[0]
    res2 = np.where( (np.absolute(np.dot(dir_vec, dir2)) * 2) > size2 )[0]
    res3 = np.where( (np.absolute(np.dot(dir_vec, dir3)) * 2) > size3 )[0]

    return list( set().union(res1, res2, res3) )

Answer 3

Here is the algorithm that I used for solving a similar problem, where my cube wasn't aligned with the axis of the 3d space. It might not be the fastest option but it works.

1. Pick points a, b, c, d so that a will be the origin of your cube, and with b, c, d will form the the axis of the cube (a and b are connected through an edge, same for a and c, and a and d). It doesn't have to be aligned with the axis of your 3D space.

2. Calculate the edgeLength of the cube (ex. distance from a to b), using the following formula:

   distance = sqrt(sqr(a.X - b.X) + (a.Y - b.Y) + (a.Z - b.Z))
   

3. Determine if a point is inside the cube

   double latA = a DistanceTo(myPoint);
   double latB = b.DistanceTo(myPoint);
   if (Math.Abs(latA + latB - edgeLength) < 0.0001)
   {
       return true;
   }
   double angleA = CalculateAngleA(latA, latB, edgeLength);
   if (angleA > 90.0001) return false;
   double angleB = CalculateAngleA(latB, latA, edgeLength);
   if (angleB > 90.0001) return false;
   
   latB = e.DistanceTo(myPoint);
   if (Math.Abs(latA + latB - edgeLength) < 0.0001)
   {
       return true;
   }
   angleA = CalculateAngleA(latA, latB, edgeLength);
   if (angleA > 90.0001) return false;
   angleB = CalculateAngleA(latB, latA, edgeLength);
   if (angleB > 90.0001) return false;
   
   latB = d.DistanceTo(myPoint);
   if (Math.Abs(latA + latB - edgeLength) < 0.001)
   {
       return true;
   }
   angleA = CalculateAngleA(latA, latB, edgeLength);
   if (angleA > 90.0001) return false;
   angleB = CalculateAngleA(latB, latA, edgeLength);
   if (angleB > 90.0001) return false;
   
   //if all validations pass
   return true;
   

   And the formula for calculating the angless (in degrees):

   double CalculateAngleA(double latA, double latB, double latC)
   {
       return Math.Acos((latA*latA + latC*latC - latB*latB)/(2*latA*latC))*(180.0/Math.PI);
   }
   

Answer 4

1. Pick Any One Out of Four Diagonals Of Cube

2. Check Your Point is within the range of diagonal points.

example. you picked diagonal having start(0,9,1) and end(10,3,-1)

So for P1(5,6,0) 
P.x (5) > Start.x(0) & P.x(5) < End.x(10)  
P.y (6) < Start.x(9) & P.y(6) > End.x(3)  
P.z (0) < Start.x(1) & P.z(0) > End.x(-1)

P1 is inside

for P2(5,6,6)
P.x (5) > Start.x(0) & P.x(5) < End.x(10)  
P.y (6) < Start.x(9) & P.y(6) > End.x(3)  
P.z (6) > Start.x(1) & P.z(6) > End.x(-1) <-- Not Within Range

P2 is out side

Answer 5

If you want to implement idea from the linked post, it makes sense to consider axis aligned cubes (parallelepipeds, actually). In that case the check is xmin<=x<=xmax && ymin<=y<=ymax && zmin<=z<=zmax.

Answer 6

Probably the easiest way to do this is to calculate the plane equation for each of the 6 planes that bound the cube, plug the point into each one and make sure the resulting sign is positive (or negative, depending of whether you calculate your planes to face inwards or outwards). The plane equation is p * normal + k = 0, calculcate the normal by taking the cross product between two edges and then plug one of the points into the plane equation to get k.

A more advanced method would be to imagine the cube defining an X, Y and Z axis and an offset (defined by cube[0]) and plug these into a matrix to convert points between the two spaces. Transforming your point by the inverse of this matrix will put it in "cube space" where the cube is aligned to the X/Y/Z axis so you can then just do a magnitude comparison against the sides.