1

I've been tasked to develop an algorithm that, given a set of sparse points representing measurements of an existing surface, would allow us to compute the z coordinate of any point on the surface. The challenge is to find a suitable interpolation method that can recreate the 3D surface given only a few points and extrapolate values also outside of the range containing the initial measurements (a notorious problem for many interpolation methods).

After trying to fit many analytic curves to the points I've decided to use RBF interpolation as I thought this will better reproduce the surface given that the points should all lie on it (I'm assuming the measurements have a negligible error).

The first results are quite impressive considering the few points that I'm using.

Interpolation results

In the picture that I'm showing the blue points are the ones used for the RBF interpolation which produces the shape represented in gray scale. The red points are instead additional measurements of the same shape that I'm trying to reproduce with my interpolation algorithm.

Unfortunately there are some outliers, especially when I'm trying to extrapolate points outside of the area where the initial measurements were taken (you can see this in the upper right and lower center insets in the picture). This is to be expected, especially in RBF methods, as I'm trying to extract information from an area that initially does not have any.

Apparently the RBF interpolation is trying to flatten out the surface while I would just need to continue with the curvature of the shape. Of course the method does not know anything about that given how it is defined. However this causes a large discrepancy from the measurements that I'm trying to fit.

That's why I'm asking if there is any way to constrain the interpolation method to keep the curvature or use a different radial basis function that doesn't smooth out so quickly only on the border of the interpolation range. I've tried different combination of the epsilon parameters and distance functions without luck. This is what I'm using right now:

from scipy import interpolate
import numpy as np

spline = interpolate.Rbf(df.X.values, df.Y.values, df.Z.values,
                            function='thin_plate')
X,Y = np.meshgrid(np.linspace(xmin.round(), xmax.round(), precision),
                      np.linspace(ymin.round(), ymax.round(), precision))
Z = spline(X, Y)

I was also thinking of creating some additional dummy points outside of the interpolation range to constrain the model even more, but that would be quite complicated.

I'm also attaching an animation to give a better idea of the surface.

Animation

Amessihel
  • 5,891
  • 3
  • 16
  • 40
Droid
  • 441
  • 1
  • 3
  • 18
  • Did you try also bivariate spline fitting (for instance [LSQBivariateSpline](https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.LSQBivariateSpline.html#scipy.interpolate.LSQBivariateSpline)? I think it will do a better job for describing the surface – xdze2 Oct 09 '19 at 14:02
  • First of all, there is a bug in scipy's RBF implementation when using polyharmonic splines. See https://stackoverflow.com/questions/56820251/unexpected-results-from-scipy-interpolate-rbf/57094397#57094397 Either try Gaussian RBF or a better RBF library altogether, like one mentioned in that answer. – rych Oct 09 '19 at 15:02
  • I've tried with a Gaussian RBF but results were even worse. I will try with the https://docs.scipy.org/doc/scipy/reference/generated/scipy.interpolate.LSQBivariateSpline.html#scipy.interpolate.LSQBivariateSpline and let everyone know the result in this thread. – Droid Oct 09 '19 at 15:23

1 Answers1

1

Just wanted to post my solution in case someone has the same problem. The issue was indeed with scipy implementation of the RBF interpolation. I tried instead to adopt a more flexible library, https://rbf.readthedocs.io/en/latest/index.html#. The results are pretty cool! Using the following options 

from rbf.interpolate import RBFInterpolant
spline = RBFInterpolant(X_obs, U_obs, phi='phs5', order=1, sigma=0.0, eps=1.)

I was able to get the right shape even at the edge.

Surface interpolation

I've played around with the different phi functions and here is the boxplot of the spread between the interpolated surface and the points that I'm testing the interpolation against (the red points in the picture).

Boxplot

With phs5 I get the best result with an average spread of about 0.5 mm on the upper surface and 0.8 on the lower surface. Before I was getting a similar average but with many outliers > 15 mm. Definitely a success :)

Droid
  • 441
  • 1
  • 3
  • 18
  • Thank you for trying out my recommendations and posting your results that validate my conviction that SciPy's RBF is simply wrong. Not sure why it had to be the `phs5` that agreed best with your test points, but it's great the new library even gives us all these options to investigate, even higher order polyharmonic splines! RBF interpolation remains an esoteric knowledge, but at least we now have a good python library for it. Maybe `phs5` somehow best approaches the smoothness class of the original surface whence the test points are taken from. – rych Oct 19 '19 at 13:32