Fit plane to 3D data using least squares I have some samples of data of the form $x,y$ and $z=f(x,y)$. I wish to fit a plane (i.e. $z = Ax + By + C$) to the data with the smallest mean square errors. I have found an "answer" in section 3 of this document, and several other locations too but the answers always end with variations of "now solve these equations and you can find $A$, $B$ and $C$"
I just about have the ability to solve these equations, but the process gets so messy that the likelihood of me making a mistake is quite high (and I need a guaranteed-correct answer). Surely someone has written out the full solution longhand somewhere - i.e. in the form $A=\dots$, $B=\dots$ and $C=\dots$ anyone know where this has been done?
EDIT: I see two answers already which leave out the last step as being trivial... and indeed they don't require any advanced maths... but you do need quite a big whiteboard to work it all out and the scope for making a mistake along the way is quite large. Indeed I noticed that assorted tutorials going through worked examples always have rather neat whole numbers. I find it hard to believe that there is nowhere online where someone has worked out the general case. I.e. find $A$, $B$ and $C$ in the following set of equations...
$$Ad+Be+Cf=g$$
$$Ah+Bi+Cj=j$$
$$Al+Bm+Cn=p$$
...where $d$ to $p$ are all constants.
 A: If you don't feel confident with the resolution of a $3\times3$ system, work as follows:


*

*take the average of all equations,


$$\bar z=A\bar x+B\bar y+C$$


*

*subtract from all equations, giving


$$z_i-\bar z=A(x_i-\bar x)+B(y_i-\bar y)$$ or
$$\hat z_i=A\hat x_i+B\hat y_i.$$


*

*solve the least squares system


$$\sum \hat z_i\hat x_i=A\sum \hat x_i^2+B\sum \hat x_i\hat y_i\\
\sum \hat z_i\hat y_i=A\sum \hat x_i\hat y_i+B\sum \hat y_i^2$$
which is $2\times 2$. This gives you $A$ and $B$.
By Cramer,
$$A=\frac{\sum \hat z_i\hat x_i\sum \hat y_i^2-\sum \hat z_i\hat y_i\sum \hat x_i\hat y_i}{\sum \hat x_i^2\sum \hat y_i^2-\left(\sum \hat x_i\hat y_i\right)^2}\\
B=\frac{\sum \hat x_i^2\sum \hat z_i\hat x_i-\sum\hat z_i\hat x_i\sum \hat x_i\hat y_i}{\sum \hat x_i^2\sum \hat y_i^2-\left(\sum \hat x_i\hat y_i\right)^2}$$


*

*C is found from 


$$C=\bar z-A\bar x-B\bar y.$$

For convenient evaluation of the sums, notice that
$$\sum \hat u_i\hat v_i=\sum u_iv_i-N\bar u\bar v,$$ i.e.
$$\overline{\hat u\hat v}=\overline{uv}-\bar u\bar v.$$
Then
$$\begin{align}A&=\frac{(\overline{zx}-\bar z\bar x)(\overline{y^2}-\bar y^2)-(\overline{zy}-\bar z\bar y)(\overline{xy}-\bar x\bar y)}{(\overline{x^2}-\bar x^2)(\overline{y^2}-\bar y^2)-(\overline{xy}-\bar x\bar y)^2}\\
B&=\frac{(\overline{x^2}-\bar x^2)(\overline{zy}-\bar z\bar y)-(\overline{zx}-\bar z\bar x)(\overline{xy}-\bar x\bar y)}{(\overline{x^2}-\bar x^2)(\overline{y^2}-\bar y^2)-(\overline{xy}-\bar x\bar y)^2}\\
C&=\bar z-A\bar x-B\bar y.\end{align}$$
A: Say you have $n$ data sets:
$( x_0 , y_0 , z_0 ), ( x_1 , y_1 , z_1 ), ( x_2 , y_2 , z_2 ), ( x_3 , y_3 , z_3 ),...( x_n , y_n , z_n )$
Create a matrix $A\in R^{(n+1)x(n+1)}$ s.t.
$$
        \begin{matrix}
        1 & x_0 & y_0 \\
        1 & x_1 & y_1 \\
        1 & x_2 & y_2 \\
        1 & x_3 & y_3 \\
        . & . & . \\
        1 & x_n & y_n \\
        \end{matrix}
$$
Afterwards create a column (result) vector $b\in R^{(n+1)x1}$ s.t.
$$
        \begin{matrix}
        z_0 \\
        z_1 \\
        z_2 \\
        z_3 \\
        . \\
        z_n \\
        \end{matrix}
$$
Then we get the equation $A*u = b$ where $u$ is your unknown vector
$$
        \begin{matrix}
        C \\
        A \\
        B \\
        \end{matrix}
$$
Final step is to obtain the normal equation $A^T(A*u) = A^T(b)$ which is a 3x3 system of linear equations:
$(A^TA)*u = A^T(b)$
I think the rest is trivial.
A: To do least square fit, you simply follow these three steps.
1.) Define the least square error function, for example, you could do
$ e(A, B, C) = \sum(Ax + By + C - z)^2 $
2.) Assuming you have no special constraint you want on $ A $, $ B $ and $ C $, so it is simply an unconstrained optimization problem. You differentiate the error function and set zero.
$ \frac{\partial e}{\partial A} = \sum 2(Ax + By + C - z)x = 0 $
$ \frac{\partial e}{\partial B} = \sum 2(Ax + By + C - z)y = 0 $
$ \frac{\partial e}{\partial C} = \sum 2(Ax + By + C - z) = 0$
3.) At this point, I guess your reference just ask you to solve for $ A, B, C $ and you are unsure about how to do that. 
The key observation is that these are just linear equations!
Ley say, for example, that you have these 4 data points.
(12, 31, 27)
(22, 32, 37)
(13, 33, 17)
(0, 0, 0)
You put in to $ \frac{\partial e}{\partial A} $ and see?
That gives
$ 2(12A + 31B + C - 27) + 2(22A + 32B + C - 37) + 2(13A + 33B + C - 17) + 2(0A + 0B + C - 0) $
$ = 94A + 192B + 7C - 162 = 0 $
The other equations can be simplified similarly. Then the rest is just solving 3 linear equations with 3 unknowns, any matrix method would work.
