Issue with behavior of ray-plane intersection algorithm

[STENDEC]

I'm trying to determine the point, in 3D space, where an arbitrary line/ray intersects with an infinite plane.

Using an article on Wikipedia, I tried to reproduce the presented formulas in code. This seems to work fine so long the ray is emitted from the origin of the coordinate system (0, 0, 0). As soon as I move the origin, the intersection point quickly starts to diverge from the plane in a curved path. My understanding of math is not good enough to figure out whether this is to be expected (i.e. a limitation of the formula), or an error in my code.

I know this is not a programming forum, but maybe you can spot the mathematical error in my code, if there is one. The input variables p0, p1, l0, l1 have been verified to carry the correct values.

/* plane normal */
vector p0; // origin position or normal (or point on plane, if you will)
vector p1; // tip position of normal

/* line/ray */
vector l0;
vector l1;

/* compute ray-plane intersection */
vector n = p1 - p0; // isolated plane normal
float t = dot(p0 - l0, n) / dot(l1, n);
vector hit = l0 + (l1 * t); // projects ray (e.g. rotations) correctly so long l0 is (0, 0, 0)

Thanks.

 

[jedishrfu]

All point vectors here emanate from the origin.

If P1 and P0 are two distinct points both in the plane then the resultant vector (P1 - P0) resides in the plane and that means the resultant vector would be normal to N hence (P1 - P0) . N = 0

In your example you're saying that (P1 - P0) is normal to the plane and that may be the problem.

If you're trying to compute the plane's normal then you need a third distinct point P2 and use it to compute two vectors V1 = (P1 - P0) and V2 = (P2 - P0) and then the normal would be N = V1 x V2

 

[Daz]

You've forgotten to subtract l0 when calculating the dot product of the ray's "velocity" with the plane normal.

Try this:

/* plane normal */
vector p0; // origin position or normal (or point on plane, if you will)
vector p1; // tip position of normal

/* line/ray */
vector l0;
vector l1;

/* compute ray-plane intersection */
vector n = p1 - p0; // isolated plane normal
float t = dot(p0 - l0, n) / dot((l1-l0), n);
vector hit = l0 + ((l1-l0) * t); // projects ray (e.g. rotations) correctly so long l0 is (0, 0, 0)

 

[STENDEC]

All point vectors here emanate from the origin.

If P1 and P0 are two distinct points both in the plane then the resultant vector (P1 - P0) resides in the plane and that means the resultant vector would be normal to N hence (P1 - P0) . N = 0

In your example you're saying that (P1 - P0) is normal to the plane and that may be the problem.

If you're trying to compute the plane's normal then you need a third distinct point P2 and use it to compute two vectors V1 = (P1 - P0) and V2 = (P2 - P0) and then the normal would be N = V1 x V2

My setup might have sent a confusing message, sorry! P0 and P1 are the two vertices of a line representing the plane normal, where P0 sits on the plane and P1 offset from P0 by the plane's "normal". I did this because I wanted to work with geometry (a line primitive needing two end points) in a 3D viewport to gauge whether I'm on the right track.

You've forgotten to subtract l0 when calculating the dot product of the ray's "velocity" with the plane normal.

Try this:

Unfortunately this only amplified the error in projection, considerably so.

Could you guys confirm whether the formula from the article should work with a ray emanating from a point not at the origin? The use of l and l0 in the article gave me the impression that an origin, combined with an offset, are a feature; maybe I'm misinterpreting the meaning.

As an end result, I want to be able to "slice" randomly positioned lines using planes, hence the desire to find a way to determine collision points.

 

[Daz]

I assume that your ray is represented in the same way as the normal. i.e. The ray is represented as a vector from l0 to l1 in the same way as the normal is represented as a vector from p0 to p1? If so, I can’t see any problem.

p0 is a point in the plane and p1 is such that p1-p0=n the normal to the plane?

l0 is the starting point of the ray and l1 is such that l1-l0 is the direction of the ray? (Usually called the direction cosines vector.)

On the other hand, if l0 is the starting point of the ray and l1 the direction cosines, then your original code should have worked.

Can you give a concrete examples listing the x,y and z components of each vector?

 

[STENDEC]

I assume that your ray is represented in the same way as the normal. i.e. The ray is represented as a vector from l0 to l1 in the same way as the normal is represented as a vector from p0 to p1?

Correct, same approach for both.

p0 is a point in the plane and p1 is such that p1-p0=n the normal to the plane?

Yes, I visualized n to confirm that it's the normal of the plane offset from the origin of the coordinate system, just as one would expect.

l0 is the starting point of the ray and l1 is such that l1-l0 is the direction of the ray? (Usually called the direction cosines vector.)

On the other hand, if l0 is the starting point of the ray and l1 the direction cosines, then your original code should have worked.

Can you give a concrete examples listing the x,y and z components of each vector?

I've created two images that hopefully help illustrate the behavior: Ray at origin, ray offset from origin (along the y axis).

While preparing the images, I noticed that the intersection point (8) no longer moves in curved paths (as described in my original post). It now stays on the plane at all times, just increasingly far away from the ray as it is moved. I don't know what I changed to affect this (code looks the same), but I imagine it's closer to the solution.

The corresponding vector values for the two images are:
(points 0-3 are inaccessible to the code, hence cannot be a source of error)

Ray at origin:
Point 4 = p0 is: {0.000000, 0.000000, 1.000000}
Point 5 = p1 is: {-0.263435, -0.674302, 1.689868}
Point 6 = l0 is: {0.000000, 0.000000, 0.000000}
Point 7 = l1 is: {0.000000, -0.199539, 0.979890}
Point 8 (hit) is: {0.000000, -0.169831, 0.834001}

Ray offset from origin:
Point 4 = p0 is: {0.000000, 0.000000, 1.000000}
Point 5 = p1 is: {-0.263435, -0.674302, 1.689868}
Point 6 = l0 is: {0.000000, -0.175000, 0.000000}
Point 7 = l1 is: {0.000000, -0.374539, 0.979890}
Point 8 (hit) is: {0.000000, -0.405668, 0.603486}

 

[Daz]

Well, this is odd. Just plugging your numbers into a pocket calculator, and going through the lines of code one by one, I don't get the same values as you. For the first case, I get the intersection at (0, -0.1699742, 0.8338610) and for the second case it's at (0, -0.3158999, 0.6912280). There's a small discreoancy, even in the first case. What's more you can easily check that my values are exactly right by taking the dot product of (hit-p0) with (p1-p0) which comes out to be zero showing that the hit lies in the plane. Then take the dot product of (hit-l0) with (l1-l0) and divide by the norm of both vectors and that comes out to be exactly 1, showing that the hit lies on the ray.

I don't see how you obtained the values you did because if your code really was as shown above, then it would have to give the values I calculate. The only thing I can think of is that you've mis-transcribed one of the numbers (or I've mis-read them) and you have omitted the second correction that I made to your code. Did you notice I made two changes? (The last two lines.)

 

[STENDEC]

The only thing I can think of is that you've mis-transcribed one of the numbers (or I've mis-read them) and you have omitted the second correction that I made to your code. Did you notice I made two changes? (The last two lines.)

Ah I'm dumb, I had only seen the change in the last line and missed the first l1 - l0 subtraction in the assignment to t. That also explains why your code seemed more "off" to me than what I had before.

Now that both lines have this operation, the projection works correctly no matter how I offset the ray's origin. Thanks a lot for you help!