# Find extra arbitrary two points for a plane, given the normal and a point that lies on the plane

For a plane, I have the normal $n$, and also a point $P$ that lies on the plane.

Now, how am I going to find extra arbitrary two points ($P_1$ and $P_2$) for the plane so that these three points $P$, $P_1$ and $P_2$ completely define the plane?

The solution here suggests that one assumes a certain $x$ and $y$ to substitute into the plane equation and find the remaining $z$. But this method is only suitable for hand calculation; it breaks down for plane $z=0$. As such, it is not suitable for computer implementation.

I would need an algorithm that is robust and can handle all the cases, any idea how to construct it?

There is a similar question here, and the answer suggests me to use Gram-Schmidt, but I don't see how it can be applied in my case here.

-
Your question is flawed in at least two ways: 1. You say "the basis" for a plane, which suggests that the plane only has one basis. Perhaps this is just a language problem. 2. As Shiyu pointed out, if the plane doesn't include the origin, then it is not a vector space, so it doesn't have a basis. So $P$ plays no role. –  TonyK Sep 14 '11 at 13:11
@TonyK, if the plane doesn't pass through the origin, then it still has two basis-- in the sense that cross product these two basis will give you the normal. –  Graviton Sep 14 '11 at 14:35
@TonyK, I've updated the question so that it is clearer (hopefully!) –  Graviton Sep 14 '11 at 14:43
Forget about $P$. When you have two basis vectors $\{e_1,e_2\}$ for the plane through the origin, you can just add $P$ to both of them to get your points $P_1=P+e_1,P_2=P+e_2$. –  TonyK Sep 14 '11 at 16:57
@TonyK, that's the problem, how to get that two basis vectors ${e_1, e_2}$? –  Graviton Sep 15 '11 at 1:38

The fundamental problem here is to find two vectors $u, v$ that are orthogonal to the vector $n = (n_x, n_y, n_z)^T$. Morevover, the procedure should be numerically stable. I assume that $|n| = 1$.

Approach 1. Find the component of $n$ that has the smallest absolute value, say $n_z = \epsilon_z$. The vector $n$ is hence "almost" in the x/y-plane (z = 0). Take the vector $w = (0, 0, 1)^T$, with all components zero except that one on the place of the smallest component. Compute $u = w \times n$. $u$ is orthogonal to $n$ and is computed numerically stable because it is "almost" orthogonal to $n$. Finalize the procedure by computing $v = n \times u$.

Approach 2. Often Householder orthogonalization is recommended for ill-conditioned orthogonalization. It could work better than the Approach 1.

1. Build that "mirror" vector $w = (n_x + 1, n_y, n_z)$ or $w = (n_x - 1, n_y, n_z)$ whose first component has bigger absolute value.
2. Compute the 3x3 Householder matrix $H = I - 2 w w^T / (w^T w)$
3. The first row of H will be a unit vector parallel to $n$, and the other two rows will be unit vectors orthogonal to $n$ and to each other, i.e. these are the vectors $u^T, v^T$.

After $u, v$ are computed, we set $P_1 = P + u, P_2 = P + v$.

-

Case 1: If the plane does not pass the origin, then there is no basis of this plane because the plane does not form a subspace.

Case 2: If the plane passes the origin, then only one normal vector $n$ is already sufficient to determine this plane, and there is no need to give the point $P$.

In the case 2, let $x$ be an point on the plane, then $$n^Tx=0$$ Using Singular Value Decomposition (SVD), you can find an orthogonal basis of the plane. The procedure is as below. Suppose one SVD of $n\in\mathbb{R}^d$ is $$n=U\Sigma V^T$$ where $U\in\mathbb{R}^{d\times d}$ and $V\in\mathbb{R}$. Denote $U=(u_1,u_2,\cdots,u_d)$, then $$u_i^Tn=n^Tu_i=0$$ for all $i=2,\cdots,d$. So one orthogonal basis of the plane is $\{u_2,\cdots,u_d\}$.

-
I think you err in Case 1: if the plane doesn't pass the origin, there is still basis of this plane, think $z=6$ –  Graviton Sep 14 '11 at 12:18
I don't get it. Bases only exist in subspaces. If a plane doesn't pass the origin, it is not a subspace. If the plane doesn't pass the origin, then the point $P$ on the plane is required to determine the plane. In this case, the equation of the plane is $n^T(x-P)=0$. –  Shiyu Sep 14 '11 at 13:30
In the case 1, any point $x$ in the plane satisfies $n^T(x-p)=0$. So one "basis" is $\{u_2+p, \cdots, u_d+p\}$, if you would like to call it as a basis. –  Shiyu Sep 14 '11 at 13:49
just want to clarify, the plane $z=6$, what is the basis for this plane? and what is the algorithm to compute it? Note that this plane doesn't pass through the origin, but there are still basis in the sense that crossing these two basis will give you the normal vector. –  Graviton Sep 14 '11 at 14:34
I've updated the question so that it is clearer (hopefully!) –  Graviton Sep 14 '11 at 14:43

Let $n=(a,b,c)$. If $a=b=0$, take $(1,0,0)$ and $(0,1,0)$ as the basis. Otherwise, take $u=(b,-a,0)$ and $v=n \times u$.

-
I think this is not a good solution. if $a$ and $b$ is small, then the $u$ and $v$ will be inaccurate due to rounding off error. Also, in the second case, you assume that there must be one basis vector lying at $z=0$, but this may not be true for all cases –  Graviton Sep 14 '11 at 12:49

Here is a naive "hack" that will likely work, and will avoid complex computations. My apologies to the purists.

First, that the plane may not pass through the origin is irrelevant, as others have mentioned. Just imagine moving it so it does pass through the origin (e.g., by subtracting $p$ from every point on the plane), and after computing the $p_1$ and $p_2$ you want for the shifted plane, add $p$ to both points to shift the plane back.

Second, you only need find $p_1$ robustly, because once you have $p_1$ not too close to the origin (the former $p$), you can just rotate $p_1$ about $n$ by $90^\circ$.

So it all comes down to finding a point $p_1$ on a plane through the origin that is not itself too close to the origin, given the normal $n$. Let $n=(a,b,c)$ as per lhf. Any point $(x,y,z)$ on the plane satisfies $a x + b y + c z = 0$.

Here's the hack. Choose a random point on the plane, say by selecting two coordinates, generating random numbers, and solving for the third coordinate. Of course you have to be a bit intelligent about the cases where one or two coordinates are zero, but that is not too difficult. If you end up too close to the origin for your robustness requirements, try again.

Eventually you find $p_1$. Rotate about $n$ to get $p_2$. Finally shift back to $\{ p, p_1 + p, p_2 + p \}$.

-

If you any component of your normal is non-zero (and we know at least ONE component won't be), any value of the other two components will fall SOMEWHERE on the plane. Say you have a normal with a positive or negative y component, you know any pair of x, z values will fall somewhere on its plane. Plug in any two values for it and solve for the last using the plane equation.

vec3 get_point_on_plane(vec3 const& point, vec3 const& norm) {

//get random x, y, z points
float x = point.x + 1;
float y = point.y + 1;
float z = point.z + 1;

/*
Plane equation:
Ax + By + Cz = D  ==
norm.x(x - point.x) + norm.y(y - point.y) + norm.z(z - point.z) = 0
*/
if (! nearly_equal(norm.x, 0)) //any y, z pair will fall SOMEWHERE on the plane
{
//solve for x with rand_y and rand_z
x = point.x - ( norm.y * (y - point.y) + norm.z * (z - point.z) ) / norm.x;

}
else if (! nearly_equal(norm.y, 0)) //any x, z pair will fall SOMEWHERE on the plane
{
//solve for y with rand_x and rand_z
y = point.y - ( norm.x * (x - point.x) + norm.z * (z - point.z) ) / norm.y;
}
else //any x, y pair will fall SOMEWHERE on the plane
{
//solve for z with rand_x and rand_y
z = point.z - ( norm.x * (x - point.x) + norm.y * (y - point.y) ) / norm.z;
}
return vec3(x, y, z);


}

-