Simple way to interpolate between points in 3D space to form a smooth surface

Question

I'm trying to come up with a simple and efficient way to create a smooth surface which intersects a number of given "sample" points.

For any X,Y point on the surface, I identify up to 4 sample points in each of the 4 directions (the next higher and lower points on the X, and then the Y axes). Given this point, I want a way to compute a Z value that interpolates between the 4 sample points.

Of course the algorithm, given the X, Y position of any of the 4 sample points, should output the Z value for that point. Note also that there may be less than 4 sample points.

I guess some function of the Z values for the 4 sample points, somehow inversely biased by the distance to the sample point, but I can't figure out how to do this.

Anyone got any ideas about a simple way to do this?

Answer 1

Use catmull-rom patches

Answer 2

You can utilize bilinear/bicubic interpolation, but in three directions (trilinear/tricubic, respectively). It is quite trivial if you understand how these forms of interpolation work. See Tricubic Interpolation on Wikipedia for more information.

Answer 3

Try NURBS.

Answer 4

If you want a simple linear interpolation of that point, then the Z value of the center point is just the mean of the 4 neighboring Z-values, given the distances are symmetric in both Y and X.

If the distances are not symmetric, but the center point is always on the same X and Y lines, you can calculate both Y and X interpolations and the final value is the mean of those.

So Zc would be: Zc=(Zx1+x*(Zx2-Zx1)/(x2-x1)+Zy1+y*(Zy2-Zy1)/(y2-y1))/2, where x and y are distances from x1 and y1.

Answer 5

Are you looking for a surface interpolation or would a grid be enough?

For a surface interpolation I see that other have suggested using triangulations (for example use this: http://en.wikipedia.org/wiki/Delaunay_triangulation)

For creating a grid: one of my colleagues used the heat equation (http://en.wikipedia.org/wiki/Heat_equation) to compute the values for pixels outside the given sample points. This produced extremely realistic looking terrain surfaces, and it was trivial to parallelize.

Answer 6

You can do this by constructing patches from Catmull-Rom splines. These splines will hit each of the control points and they are continuous in the first derivative (though not the second). I also find them to be extremely easy to work with. The math is straightforward and they behave intuitively with slight changes in the control points.

At the highest level, you'll need 16 points per patch (at the edge of your dataset, you can use the corner and edge points twice in the same spline).

First, you'll need to interpolate across the points p[i][j] in each row in your 4x4 matrix to create a set of four intermediate control points q[i]. Here's a rough ASCII sketch of what I mean.

p00 p01 q0 p02 p03
p10 p11 q1 p12 p13
p20 p21 q2 p22 p23
p30 p31 q3 p32 p33

Now you can interpolate between each of those four intermediate control points to find a final splined point on your surface.

Here is a construction of the Catmull-Rom spline across four points. In this example, you are interpolating between points p[i-1] and p[i], using control points on either side p[i-2] and p[i+1]. u is the interpolation factor ranging from zero to one. τ is defined as the tension on the spline and will affect how tightly your splined surface conforms to your control points.

                 | 0   1   0    0 | | p[i−2] |
                 |−τ   0   τ    0 | | p[i−1] |
p(u) = 1 u u2 u3 | 2τ τ−3 3−2τ −τ | | p[i]   |
                 |−τ  2−τ τ−2   τ | | p[i+1] |

NOTE: it's not immediately obvious how to lay this out in Stackoverflow's gui but u2 and u3 are supposed to represent u squared and u cubed, respectively.

Answer 7

A problem when interpolating using the scheme suggested in the question, wherein some subset of nearest neighbors are chosen from a scattered set, is that the result need not be continuous.

Think about it. Suppose I were to move along some smooth, continuous path through the (x,y) plane. As long as the 4 nearest neighbors don't change, the interpolant will be smooth, defined however you chose to do so. However, at some point, that subset of nearest neighbors will suddenly change. At that point, you must have the interpolant consistent across the boundary. Best is to use a triangulation of the independent variable set. This assures a continuous (linear) interpolant within the convex hull of the data. With more work, higher order interpolations can also be achieved with a triangulation.

Radial basis functions are also commonly also be used for interpolation, or Kriging, if you happen to be into geostatistics. Since you were looking into distance based methods, consider radial basis functions. For example, search for "inverse multiquadric interpolation".