Off your first post - you have a bunch of data points which are from some function f(x,y,z) right? You have asked for a program to recover the equation of f(x,y,z) from the set of ordered triples {(xi,yi,zi)}. The problem is that f(x,y,z) can be anything.
Look at the simpler example in 2D ... if I have ordered pairs {(1,1),(2,1),(3,1),...} we could say that the curve is y=1 ... but it could also be y=cos(kx): k=2pi, so which is it? There are an infinite number of possible cosine curves... then there are other periodic functions ... and, in this case, there are infinite data points. If there were finite data points, then I could also fit polynomials. Then there are all the possible irregular functions, peicewise functions and on and on and on. In 3D f(x,y,z) could even loop back on itself.
Programs like MATLAB have built-in assumptions that they use to do interpolation - they do not, in general, find the equation of the generating function.
gridfit, your example, for instance, assumes the data corresponds to a surface z=f(x,y)... and makes assumptions about the nature of the surface.
For the 2D case I could use the polyfit function - but I have to also input the order of the polynomial I want to fit.