Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I am working with many 3D planes and looking for a Least square solution for below case.

IF I am having many number of 3D planes knowing only one point and the normal vector (for eg. O1 and N1), and all those planes intersect each other and make nearly very close 3d lines, then how to compute the least square adjusted one single 3d line to represent all those intersections.

To get a clear idea, I have inserted a figure.

  • Known: one point and the normal vector of each plane.
  • Find: Least Square fitted single line3d

enter image description here

By kowing the 3d lines like O1-B1 and O1-E1, I know this can be done taking pairwise planes and then by finding 3dline segments. i.e. end points and then by finding the best fitted 3dline. but I heard without doing this and considering many planes at once, can also do this. SO, looking for the solution for that with Least Square. thanks in advance

share|cite|improve this question
How do segments get into it? Everything in the question except that one point where you say "segments, i.e. endpoints" seems to refer to unbounded subspaces; how do endpoints arise? – joriki Mar 14 '13 at 14:22
@joriki: sorry, it is my mistake of simplifying the question. i will correct it. – gnp Mar 14 '13 at 14:40
up vote 2 down vote accepted

In order to get that "optimal" line, the key is to establish the objective functions. In order to determine a line in 3D, we need to determine a point $p\in\mathbb{R}^3$ on this line, and the vector $v\in\mathbb{R}^3$ which is parallel to the line. Next I will first determine the optimal $p$ and then the optimal $v$.

Step 1: determine $p$

It is natural to choose $p$ as the point that is closest to all of the planes. For plane $i$, denote $n_i\in\mathbb{R}^3$ and $o_i\in\mathbb{R}^3$ respectively as its normal vector and a point on the plane. Then the distance between a point $p$ and the plane can be expressed as

$$\left\|\left(\frac{n_i n_i^T}{n_i^Tn_i}\right)(p-o_i)\right\|$$

Hence to determine the point $p$, we can minimize the objective function $$J_p= \sum_{i=1}^n \left\|\left(\frac{n_i n_i^T}{n_i^Tn_i}\right)(p-o_i)\right\|^2 $$ Denote $P_i=\frac{n_i n_i^T}{n_i^Tn_i}$, then $J_p$ can be written as $$ J_p=\sum_{i=1}^n (p-o_i)^T P_i (p-o_i) $$ Note $P_i$ is an orthogonal projection matrix satisfying $P_i=P_i^T$ and $P_i^2=P_i$. A short calculation shows that $$ p^* = \left(\sum_{i=1}^n P_i\right)^{-1}\left(\sum_{i=1}^nP_io_i\right)$$ You can check the condition under which $\sum_{i=1}^n P_i$ is invertible.

Step 2: determine $v$

Intuitively speaking, the angle between the vector $v$ and a plane should be as small as possible. In other words, $v$ should be perpendicular to the normal vector $n_i$ as much as possible. Or we can say the orthogonal projection of $v$ on $n_i$ should as small as possible. Note we only need to determine the direction of $v$. Hence we can assume $v^Tv=1$. To determine $v$ we can minimize the objective function $$ J_v= \sum_{i=1}^n \left\|\left(\frac{n_i n_i^T}{n_i^Tn_i}\right) v \right\|^2 = v^T \left(\sum_{i=1}^n P_i\right) v$$ Denote $A=\sum_{i=1}^n P_i $. Clearly $A$ is positive definite (you need to check under what condition $A$ is invertible). Let the SVD of $A$ be $A=U\Lambda U^T$. Hence $$ J_v=v^T A v \ge \lambda_{\mathrm{min}}(A) v^T v = \lambda_{\mathrm{min}}(A)$$ where the equality holds when $v$ is the last column of $U$.

share|cite|improve this answer
thank you for the response. I am studying some basics and contact you again for the clarifications. – gnp Mar 17 '13 at 20:53
again, i am thanking for giving this valuable solution. To be honest, me i am poor in mathematics and lacking basics. if you can could you please emphasis how this describe with matrices, because my objective is to implement this in c++. please. – gnp Mar 19 '13 at 21:19
The solution given above are all described with matrices. For implementation, I think it would not be complicated. You only need to implement the $p^*$ equation and the last sentence. But for C++, I think you need to have some libraries regarding matrix operation. – Shiyu Mar 20 '13 at 1:59
to be honest, the implementation of second part is still confusing me. What should I actually tested? what do you mean by last sentence. please be kind enough to simplify this problem again. thank you – gnp Jul 12 '13 at 20:13

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.