Ray and Ellipsoid intersect

[PrisonLove]

How can I find the points of intersection between a ray and an ellispoid centered at an arbitrary point (i.e. not the origin)?

I simply just cannot do the math or form it into Java code.

So far I know that there is some 4x4 Matrix Q that represents the quadric surface for an Ellispoid, and there is some ray x = P + tu, where P is the ray start point, u is the direction vector (normalized) and t is some variable for a point along the ray. I also know that to get the intersection I need the equation xTQx = 0, where xT is the transpose of x. Expanding this equation out we get:

(P+tu)TQ(P+tu)

My problem is, I can't figure out what Q should be, and even if I did, I simply can't do the math to get the closest intersection point to the ray origin. Any help?

 

[saw792]

Fortunately, there is a simpler way to go about doing this. An ellipsoid can be represented as a sphere that has had some transformations applied to it:

1. Scaled in the x, y and z directions (for the different radii)
2. Rotated about the origin in x, y and z directions (requires two angles)
3. Translated to a new, non-origin location

These transformations can be represented by a transformation matrix, let's call it K. Since we are interested in affine transformations as well (i.e. translations) this will be a 4x4 matrix. The translation data is held in the 4th column.

Thus the transformation of each point (vector) x on the original sphere would be Kx.

The approach you are currently taking is to find when your ray, P + tu intersects this ellipsoid Kx. The other approach is to see when the ray K(-1)P + tK(-1)u intersects the sphere x. K(-1) means the inverse of matrix K (no superscript ).

So the first step in finding intersections this way is to construct the transformation matrix K. The order of transformations is important here. When changing a sphere to an ellipsoid the order is Scale, Rotate, Translate, so K can be constructed by multiplying these 3 matrices together: K = TRS (note the order of multiplication - it's the opposite of what you might think). The scaling matrix is quite simple, it's just taken from the radius along each major axis, eg:

[ rx  0  0  0 ]
[ 0  ry  0  0 ]
[ 0   0 rz  0 ]
[ 0   0  0  1 ]

The rotation matrix is more complicated. It's actually composed of the multiplication of 3 rotation matrices that deal with only one axis at a time. Where a, b and c are the rotation angles for the x, y and z axes respectively the combined rotation matrix R is as follows:

[ cos(b)cos(c)    -cos(a)sin(c) + sin(a)sin(b)cos(c)    sin(a)sin(c) + cos(a)sin(b)cos(c)        0 ]
[ cos(b)sin(c)     cos(a)cos(c) + sin(a)sin(b)sin(c)    -sin(a)cos(b) + cos(a)sin(b)sin(c)       0 ]
[   -sin(b)                   sin(a)cos(b)                          cos(a)cos(b)                 0 ]
[      0                           0                                     0                       1 ]

Translation matrix is also pretty simple:

[ 1  0  0 Tx ]
[ 0  1  0 Ty ]
[ 0  0  1 Tz ]
[ 0  0  0  1 ]

So you need to multiply these 3 matrices together in order K = T x R x S, then calculate the inverse and apply to the both the position vector and the direction vector of the ray. Then you can simply use your collision code for spheres with the new ray. If you need to recover the original ray you can do so by multiplying by K once again. Let me know if you need any more help.

 

[PrisonLove]

Wow, that really helps a lot, thank you so much. I do have a couple of questions though.

1. Given only these two fields for the Ellipsoidal class:

  public double x,y,z;	     // Center of the ellpsoid.
  public double rx,ry,rz;   // Radii of the ellipsoid.

How would I find the Ellipsoids rotation?

2. Where finding the Inverse matrix is concerned, is there a function to do this in java, or do I just have to use the reciprocal of all the values so that the matrix becomes (I'm purposely leaving out the rotation here)
M^-1

[ 1/rx  0   0    1/x ]
[  0  1/ry  0    1/y ]
[  0    0  1/rz  1/z ]
[  0    0   0     1  ]

 

[saw792]

Your description of an ellipsoid does not take into account rotation at all, which means that you can skip the rotation matrix completely or (perhaps more sensibly) allow for rotations in the future by leaving a placeholder there (which would be the identity matrix). If you want to add rotations in the future you can just add fields for x, y, and z rotation angle and use the formulae I gave above.

The inverse of a matrix is not simply the inverse of all it's values. It is defined as the matrix M^-1 such that M^-1M = I, the identity matrix (which can be thought of as a matrix equivalent of 1). Similarly the multiplication of two matrices is not simply the product of the elements in the same position.
See http://en.wikipedia.org/wiki/Matrix_(mathematics) for a better description of these operations.

There is a java package that can handle matrix operations called JAMA (google it ) which should do the job properly for you. It's not included in the standard library so you will have to import it. An example of intersect code in Java:

public boolean intersect(Matrix raypos, Matrix raydir) {
		Matrix S = new Matrix({{rx,  0, 0, 0}, 
								{0, ry, 0, 0}, 
								{0, 0, rz, 0}, 
								{0, 0, 0, 1}});
		Matrix R = new Matrix({{1, 0, 0, 0}, 
							   {0, 1, 0, 0}, 
							   {0, 0, 1, 0}, 
							   {0, 0, 0, 1}});
		Matrix T = new Matrix({1, 0, 0, x}, 
							  {0, 1, 0, y}, 
							  {0, 0, 1, z}, 
							  {0, 0, 0, 1}});
		Matrix K = T.times(R).times(S);
		Matrix InvK = K.inverse();
		
		Matrix kraypos = InvK.times(raypos);
		Matrix kraydir = InvK.times(raydir);
		
		Sphere es = new Sphere(...);
		
		return es.intersect(kraypos, kraydir);
	}

Good luck