9

Representing plane intersections as a system of linear equations

 3 years ago
source link: https://gieseanw.wordpress.com/2017/04/27/representing-plane-intersections-as-a-system-of-linear-equations/
Go to the source link to view the article. You can view the picture content, updated content and better typesetting reading experience. If the link is broken, please click the button below to view the snapshot at that time.

Representing plane intersections as a system of linear equations – Andy G's BlogSkip to content Skip to menu

DoritoPlanet

When you are watching a digitally-rendered battle onscreen in the latest blockbuster movie, you don’t always think about the “camera” moving about that scene. In the real-world, cameras have a field of view that dictates how much of the world about them they can see. Virtual cameras have a similar concept (called the viewing frustum) whereby they can only show so much of the digital scene. Everything else gets chopped off, so to speak. Because rendering a digital scene is a laborious task, computer scientists are very interested in making it go faster. Understandably, they only want to spend time drawing pixels that you’ll see in the final picture and forget about everything else (outside the field of view, or occluded (hidden) behind something else).

Our friends in the digital graphics world make heavy use of planes everyday, and being able to test for plane-plane intersection is extremely important.

In this post, I’ll try to break down what a plane is in understandable terms, how we can create one given a triangle, and how we would go about testing for the intersection between two of them.

(This article was originally written in \LaTeX and transcribed here. To get a hard copy that you can also cite, grab the pdf here)

What is a plane, anyway?

Say you’re given a triangle in 3D space. It consists of three points, P_0, P_1, and P_2 that each have x, y, and z components:

 P_0   =  \left(x_0, y_0, z_0\right)
P_1   =  \left(x_1, y_1, z_1\right)
P_2   =  \left(x_2, y_2, z_2\right)

For simplicity’s sake, we’ll assume that, moving “counter-clockwise”, the triangle’s points go in the order P_0, P_1, and then P_2 (and then back to P_0).

A triangle is distinctly “flat”. It’s like a perfectly constructed Dorito. Also, no matter which way you turn it and rotate it, there’s always a definite front and back. That is, if you extended the corners of your Dorito off towards infinity, then you could completely bisect the room, your state, the planet, the universe!

DoritoPlanet
And feed a lot of hungry people, too

This ability to split a space in two is a salient feature of a plane, and we can construct one with the triangle we’ve just defined. In the next few sections we’ll go more into depth in exactly how to do that.

Normal vector

In the previous section I said that our gargantuan flat Dorito had two distinct sides, front and back. That means you could draw a straight line from the “back” that eventually comes out the “front”. Imagine that our Dorito is so large that it is its own planet orbiting the Sun. You could hop in a space ship, fly up to that Dorito, take one small step for man, and then plant your flag on top. In geometric terms, your flag pole would be the normal vector to the plane’s (Dorito’s) surface.

DoritoPlanet2
Our Dorito planet with a normal vector pointing “up” out of it

A normal vector \vec{n} to a plane is a line that is perfectly perpendicular to any other line drawn entirely within the plane. We call it a vector because it has a direction associated with it. That is, the flag you planted on the Dorito planet points “up”, but you could just have easily landed on the other side and planted the flag pointing “down”, and it would still be perfectly valid; all we need is just to choose one.

Our normal vector will also have x, y, and z components, but we’ll label them a, b, and c to distinguish them from the points in the triangle:

\vec{n} = \left(a, b, c \right)

Vector equation of a plane

At this point we have a triangle, but we only know about what a normal vector is. How do we compute the normal vector given a triangle’s points?

Our normal vector is perpendicular to any line we draw in the plane. In computer graphics terms, we say the normal vector is “orthogonal” to our plane. An nice example of orthogonality is that the Y-axis is perpendicular to the X-axis (and the Z-axis is perpendicular to both!).

In geometry, we know how to obtain a perpendicular vector (a normal vector) so long as we have 2 non-parallel vectors in the plane. But first we need to get those two vectors in the plane. So how do we do that?

Obtaining non-parallel vectors in a plane given a triangle

Remember from before that we defined our plane originally using three points, P_0, P_1, and P_2. We can use these points to compute our vectors.

Think about it this way. The center of the earth is position (0,0,0), and you’re standing at some random point (x_0, y_0, z_0). Your friend Andy is standing at position (x_1, y_1, z_1). Vectors imply a direction, so how do you get a direction from you to your friend? Subtraction! You can compute your friend’s offset from yourself by subtracting your own position from Andy’s: (x_1 - x_0, y_1 - y_0, z_1 - z_0). The resulting offset is a direction out from the center of the earth in the same direction as the direction from you to your friend.

If your other friend Bernhard is standing at position (x_2, y_2, z_2), you’d do the exact same thing to get a vector representing the direction from you to him. So long as all three of you are not standing in a line, you now have two non-parallel vectors. Fortunately, our original points P_0, P_1, and P_2 define a triangle, and none of the sides of that triangle are parallel, so we’re good.

In our example, you are standing at P_0, Andy is at P_1, and Bernhard is at P_2. We can compute two non-parallel vectors, \vec{V_0}, \vec{V_1} as

\vec{V_0}  =  \left(x_1 - x_0, y_1 - y_0, z_1 - z_0\right)
\vec{V_1}  =  \left(x_2 - x_0, y_2 - y_0, z_2 - z_0\right)

Computing a normal vector to a plane given two vectors in the plane

In geometry there is this fancy thing called a cross product that enables us to compute an orthogonal (perpendicular) vector so long as we are given two input vectors. What does this mean? Well, there are a couple different interpretations, and I’ll work up to them after first discussing something called a dot product.

The cross product is a special kind of multiplication. The multiplication you were taught in grade school resulted in numbers many times larger than they were before. For example, 9*8=72, which is quite bit larger than either 9 or 8. The “problem” here is that 72 doesn’t have any sense of direction to it!

Our grade school multiplication has a ready analog in the world of vectors (things that have direction) called the dot product: we start by  multiplying the components! Afterwards we add them all together to get a single value. The dot product can be thought of as “multiplication in the same direction” because we multiply the x components together, the y components together, as well as the z components.

\vec{V_0} \cdot \vec{V_1} = \left(x_0 * x_1\right) + \left(y_0 * y_1\right) + \left(z_0 * z_1\right)

An alternate interpretation of the dot product is that it’s a measure of just “how much” in the same direction your two vectors are:

\vec{V_0} \cdot \vec{V_1} = ||\vec{V_0}|| ||\vec{V_1}||\cos{\theta}

Where ||V_0|| and ||V_1|| are the magnitudes of the vectors:

||\vec{V_0}||= \sqrt(x^2 + y^2 + z^2)

and \theta is the angle between them.

VectorVector
Two vectors with angle Θ between them

When \theta is 0, the vectors are in the same direction (parallel), meaning the cosine is 1.0 (it’s largest value) and therefore the result is the largest possible dot product between the two vectors. (Data miners love to use this and call it “cosine similarity”). But if the vectors were orthogonal (perpendicular), the dot product would be zero.

The cross product is similar to the dot product except it’s more of a measure of how different two vectors are. The data miners might choose the following representation of a cross product:

\vec{V_0} \times \vec{V_1} = ||\vec{V_0}|| ||\vec{V_1}||\sin{\theta}

When the \theta term is 0, the vectors are in the same direction (parallel), which maximizes the dot product, but the cross product is zero! And when \theta is 90 degrees (\frac{\pi}{2} radians), the vectors are perpendicular and the cross product is maximized! (The dot product, sadly, becomes 0).

There exists another way to compute a cross product where the result is not a single number, but rather a vector that is orthogonal to both the input vectors:

\vec{V_0} \times \vec{V_1} = \left(y_0z_1 - z_0y_1, z_0x_1 - x_0z_1, x_0y_1 - y_0x_1\right) (1)

It’s not the most intuitive computation. Personally, I enjoyed the “xyzxyz” explanation at
Better Explained.

All we need to know is that the result of this equation is a shiny new orthogonal vector.

Putting it together

We now have two vectors in our plane, \vec{V_0} and \vec{V_1}, which we computed using our triangle’s points P_0, P_1, and P_2. We also know how to take two vectors and compute an orthogonal vector. Our normal vector is exactly that; an orthogonal vector to our plane, so when we apply the cross product to \vec{V_0} and \vec{V_1}, we obtain our normal vector \vec{n} = \left(a, b, c\right):

\left(a, b, c\right) = \vec{V_0} \times \vec{V_1} = \left(y_0z_1 - z_0y_1, z_0x_1 - x_0z_1, x_0y_1 - y_0x_1\right)

Which gives us

a  =  y_0z_1 - z_0y_1
b  =  z_0x_1 - x_0z_1
c  =  x_0y_1 - y_0x_1

Now that we have a normal vector, we can define our plane intuitively. We know that our normal vector \vec{n} = (a, b, c) is perpendicular to all vectors in the plane, and as we saw before, the dot product of any two vectors is zero if the vectors are orthogonal. We can therefore say that, given our point P_0 on the plane, and any other point P on the plane, the dot product between \vec{n} and \overrightarrow{\left(P - P_0\right)} (subtracting two points results in a vector) is 0:

\vec{n} \cdot \overrightarrow{\left(P - P_0\right)} = 0(2)

Which is our generic representation of a plane. This is a nice way to think about a plane, but alone it doesn’t help us find the intersection between two planes. Ultimately what we want is a system of linear equations like the title talked about.

Representing a plane as a linear equation

Given our plane equation \vec{n} \cdot \overrightarrow{\left(P - P_0\right)} = 0

Let’s expand the terms.

\left(a, b, c\right) \cdot \left(\left(x - x_0\right), \left(y - y_0\right), \left(z - z_0\right)\right) = 0

After performing the dot product:

a\left(x - x_0\right) + b\left(y - y_0\right) + c\left(z - z_0\right) = 0

Distributing a, b, and c:

ax - ax_0 + by - by_0 + cz - cz_0 = 0

We know what the values of a, b, c, x_0, y_0, and z_0 are (they’re constants) because they were given to us in the plane’s definition, so let’s move them to the other side of the = sign in order to have only variables on one side and only constants on the other.

ax + by + cz = ax_0 + by_0 + cz_0

and since ax_0 + by_0 + cz_0 is constant, it will be easier to relabel it as d instead of writing it out every time. That is,

d = ax_0 + by_0 + cz_0

Which results in a linear equation that looks like this:

ax + by + cz = d

If we have two planes, then we’ll distinguish between their definitions via subscripts on our d constant:

     a_1x  +  b_1y +  c_1z  = d_1(3)
a_2x  +  b_2y  +  c_2z  = d_2

Finally we have a system of linear equations. Great! We have ways to solve those! However, this particular one presents a little problem. (A tiny one that we’ll make go away)

A problem

In linear algebra, we are often provided with a number of equations and an equal number of unknowns (variables) we must solve for. However, in our system, we appear to have more variables than equations!

Let’s take another look.

a_1x  +  b_1y  +  c_1z  = d_1
a_2x  +  b_2y  +  c_2z  = d_2

We have three variables to solve for — x, y, z, and yet only two equations. This is called an Under-determined system meaning there are infinite solutions (assuming the planes intersect at all). Infinite? You see, when a line intersects a line, they intersect at a single point, but when a plane intersects another plane, that intersection is a line, which has infinitely many points. We would need three planes to intersect before we could find a single point of intersection.

Solving an under-determined system

There are more math-intensive approaches to solving under-determined systems that involve computing something called a pseudo-inverse, but for our purposes here, we can take advantage of our knowledge of the domain! (Also, a Moore-Penrose pseudo-inverse would only provide an approximate solution whereas we can compute an exact one here).

One solution

Let’s assume that the planes we’re intersecting are not parallel; they have a “nice” intersection. If this is true, then we’ll end up with a line. Because planes are infinitely large, the resulting line will be infinitely long. At some point, this line will cross at least one of the x, y or z axes.

How can I be so sure of that? Start by thinking of it this way: we can represent a line as a point plus a direction of travel. We can move forwards and backwards in that direction (this is called the vector equation of a line as shown below:

PointVector
A line represented as a point and a direction

If we don’t want our line to pass through x = 0, then we set x to some arbitrary constant like x = 1. In two dimensions, this would give us a vertical line. In three, it gives us a line with direction as a function of y and z. If we don’t want the line to cross the y axis, then we fix y to some arbitrary constant like y = 33. Now our line varies only in the z direction. However, if we didn’t want the line to cross the z axis, and tried to assign it to some arbitrary constant as well, say z=7, then we’d end up with a single point \left(1,33,7\right), not a line. If we allow the line’s direction to be a function of z, then z=0 is a feasible value, and if we allow it to be a function of y and z, then both y = 0 and z = 0 are feasible values. Thus we can always guarantee that our line passes through 0 for at least one of our axes. In summary, in a D dimensional system, if we fix D-1 dimensions to non-zero values and allow the last dimension to vary, 0 is included in the last dimensions’ possible values, otherwise you would not have a line. Graphically, we can show this like so:

(a) Fixing the x dimension (green plane)
(b) Fixing the z dimension (red plane)
(c) Fixing x and z dimensions
(d) Intersection of fixed x and z varies in y (white line)

When a line cross through 0 for one of our axes, this is great! It’s great because if a component at that point has a value of zero, then its associated term in equation 3 will drop out. This leaves us with a system that is *not* under-determined; we will have two variables and two linear equations. And as I discussed in my previous post, we know how to solve for that. At this point we’ve obtained a point on the line of intersection!

Two remaining problems

Okay, we know how to solve a 2D system of linear equations to get a single point, but what use is that? Indeed a point by itself isn’t a solution, we need a line. The point we obtained, however, is in fact the first part in determining the line of intersection between two planes.

We still have two problems to solve:

  • How do we discover which component will become zero?
  • Even if we obtained a single point on the line, how do we figure out the rest of the line?

Let’s solve each of these in turn.

Discovering which component will become zero

I said before that we could use our knowledge of the domain to solve this plane-plane intersection problem. That knowledge is that, if the planes intersect “nicely”, then eventually that line will pass through zero for at least one of \left(x, y, z\right). So how do we discover if it’s x, y, or z? The short answer is: why not try them all?

First we set x to 0 (causing its term to drop out) in equation 3, leaving us with two lines. We can test to see if the resulting two lines intersect at all by testing to see if they’re parallel via the cross-product. If that doesn’t work, we move onto y, and then z.

The result is that we’ll have a point on the line of intersection where one of the components is zero, like \left(x, y, 0\right). All that remains is to discover the direction of the line.

Finding the line’s direction

I actually hinted at this earlier, and you may have caught it if you were reading closely.

Up until now, we’ve assumed that the two planes intersect and are not co-planar. In other words, the planes are not parallel.

How can we make sure of that? Well, we already know our two planes have normal vectors \vec{n_1} and \vec{n_2}, respectively. Also, we know that the cross product (equation 1) between these normal vectors \vec{n_1} \times \vec{n_2} gives us a new vector that is orthogonal to both. What happens when you compute the cross product of a vector with another in the same direction? Well, since the angle between them is 0 (or 180) degrees, the resulting vector is the \vec{0} vector, with magnitude 0.

For our purposes, when we encounter co-planar planes, we will disregard this case (the intersection is a plane equal to one of the input planes).

After we’ve ensured that the cross product between the two planes’ normal vectors is non-zero, we are left with a vector that’s orthogonal to both. How is this useful?

Well, we know from our earlier discussion that a normal vector is orthogonal to any vector we could draw in a plane. So it follows that a vector that is orthogonal to two normal vectors must lie in both those planes! And since our planes intersect at exactly one line, our fancy new vector we got from the cross product of the normals must be in the same direction as the line of intersection between the planes!

Armed with this information, we now have a point and a direction to represent our line of intersection.

If you’ll allow me, let the point on the line be a, and the direction it travels in be \vec{d}. Now we can represent any point on our line as:

r = a + \lambda\vec{d}(4)

Where r is a resulting point, and \lambda is a value from -\infty to \infty. Let’s go one step further, set \lambda to 1 so that we can obtain a second point on the line, which we’ll call b. And now we have two points on the line. We can use those two points to transform our line representation into Two-Point Form (like I use in my previous post where subscript of 1 indicates a component of a (x, y, and z), and a subscript of 2 is a component of b:

y - y_1 = \frac{y_2 - y_1}{x_2 - x_1}\left(x - x_1\right)(5)

And thus we’ve computed the line of intersection between two planes.

Conclusions

Finding the line of intersection between two planes is generally done the same way as you’d intersect any two geometric objects — set their equations equal to each other and solve. We discovered here that the result was an under-determined system and we overcame that. This is exciting! From this understanding, we are ready to take the next steps and handle “real world” situations that involve line segments and triangles instead of mathematical lines and planes.


About Joyk


Aggregate valuable and interesting links.
Joyk means Joy of geeK