Inverse Square Weighted Interpolation

Suppose we wish to interpolate a surface through a set of non-uniformly 
distributed control points z = z(x, y).  Of course, there are 
infinitely many surfaces passing through any finite set of control 
points, so the result is not uniquely determined unless we specify 
some additional information, which as a criterion for "goodness of
fit", and that constrains the overall shape of the surface.

One common technique is multiple linear regression.  This requires you 
to specify the "form" of the surface.  For example, you could specify

        z(x,y) = A + Bx + Cy + Dxy + Ex^2 + Fy^2 + Gx^2y^2

and then use linear regression to determine the coefficients A,B,..,G
such that the sum of squares of the z errors at your control points is 
minimized.  To ensure that this surface actually passes through each of 
your N control points, you would need to choose a "form" with N
terms.  For example, the above form will fit any seven control points 
exactly.

Another interesting method, and one that automatically passes through 
each control point, is the "inverse square weighting" method.  Given 
the values of f(xi,yi) for N points (x1,y1), (x2,y2),...,(xN,yN), the
interpolated function is

               f(x1,y1)   f(x2,y2)         f(xN,yN)
               -------- + -------- + ... + --------
                 R1^2       R2^2             RN^2
    f(x,y) =  --------------------------------------
                    1         1               1
                   ----  +  ----  +  ...  +  ----
                   R1^2     R2^2             RN^2


where (Ri)^2 = (x-xi)^2 + (y-yi)^2.  This can obviously be extended
to any number of dimensions.  For points (x,y) far from any of the 
control points, the value of f approaches the average of the f values 
for the control points.

To illustrate, suppose we are given two control points f(0,0)=1 and 
f(1,0)=-1.  The surface generated by these two points is shown below.



For another example, suppose we are given the four points

              f(0,0) = -4          f(1,0) =  5
              f(0,1) =  5          f(1,1) = -1

The surface generated by these four points is shown below.