2D Interpolation irregular grid fortran

Question

How can I implement 2-dimensional interpolation in FORTRAN where the data looks like as shown below. x and y are two co ordinates and z is a value dependent on them x is spaced evenly but y is not uniformly spaced and the maximum value of y corresponding to uniform values of x keep on increasing. Without losing much of accuracy-

• What is the simplest method to obtain a value of z based on a given x and y?
• What is the fastest method to obtain a value of z based on a given x and y?

Thanks SM

x    y    z
-----------
0   0     -
0   0.014 -
0   0.02  -

.....
....

0.1 0     -
0.1 0.02  -
0.1 0.03  - 

.......
.....

1.0  0     -
1.0  0.05  -
1.0  0.08  -

.......
.......

Answer 1

I am going to assume you have already read your data into an array N x 3, following the format you gave. I am assuming you don't know ahead of time what the X spacing is - you definitely don't know the Y spacing as it varies. Thus I would recommend the following strategy:

• Figure out the X spacing: start at the first row, and go through the X elements until you see a change in value. You now know XSTART and XSTEP - you will need those later.
• Do a binary search in your array for value X, until you find a value XFOUND such that XFOUND < X < XFOUND + XSTART
• Assuming you are pointing "somewhere in the list", you find the corresponding Y value - depending on whether it's bigger or smaller than the value you need, you step up or down the array until you find the first entry < Y. The corresponding values are X11, Y11, Z11. The next line in the array has X12 Y12 and Z12.
• You need two more points before you can do the interpolation - repeat this process, looking for the "next larger value of X". This will give you XYZ21 and XYZ22
• Now you can think about calculating the interpolated Z value. In general there are different techniques, with different accuracies:
• "Nearest neighbor": find the closest point, and use its value of Z (simplest, least accurate)
• "linear interpolation": find the three closest points, and do a linear interpolation of values based on the relative distances
• "higher order estimates": to do this, you typically need to create a full connectivity mapping of the grid points, so you can do spline interpolation and get a smooth interpolation that will typically be more accurate at points between the grid points (assuming that the function being described by the samples is in fact a smooth function!)

My FORTRAN is a bit rusty - hope this is some help.

PS - potentially a simpler approach is to use the fact that the X values are already evenly spaced. This allows you to make a better interpolation. See this picture:

Answer 2

Before finding the fastest way to solve this problem, I'd suggest finding a way to solve the problem. In brief, I'd suggest:

1) Find a Delaunay triangulation of the (x,y) points. Fortran code to do this is in GEOMPACK, for instance.

2) To interpolate, given the Delaunay triangulation, find the triangle that contains the point where the interpolation is to be done, then interpolate the z value, based on the location of the point relative to each of the triangle vertices.

(EDIT) I had forgotten the name of the method (if I ever knew it) but, thanks to @Floris, a good approach to interpolate in a triangle is called barycentric interpolation, which finds the interpolated value based on the ratios of the areas in the three smaller triangles into which a large triangle can be split by drawing lines from a point inside the large triangle to each of the corners of the triangle. The area of each small triangle can be found from the length of each side of the triangle using Heron's formula.

If speed improvements were needed, I think they could be obtained mainly from finding a quick way to find the triangle that contains the interpolation point but I'd want to profile some test runs before choosing which bits of the code to optimise.

Answer 3

I suggest using one of the TOMS interpolation algorithms: https://people.sc.fsu.edu/~jburkardt/f77_src/toms526/toms526.html https://people.sc.fsu.edu/~jburkardt/f77_src/toms660/toms660.html https://people.sc.fsu.edu/~jburkardt/f77_src/toms790/toms790.html