Reputation: 61
I am trying to approximate a given route (coordinates) with a three dimensional cubic spline. Example data:
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d
import numpy as np
%matplotlib inline
x = np.array([1, 2, 2.3, 3, 4, 5, 5.5, 8, 9, 9.5])
y = np.arange(0, 10)
z = np.sin(x) * np.cos(y^2) + x
fig = plt.figure(figsize=(10,6))
ax = axes3d.Axes3D(fig)
ax.stem(x, y, z)
I now approximate the data points with RBF Interpolator from scipy. Is this the right approach?
from scipy.interpolate import RBFInterpolator
coord_data = np.stack([x, y], -1)
spline = RBFInterpolator(coord_data, z, kernel = 'cubic')
How do I get the resulting spline now (the points it is following through)? And how do I access it's derivatives?
Upvotes: 1
Views: 1138
Reputation: 2066
You are trying to approximate a route, i.e., a curve in 3D, and not a surface. The approach you are trying results in a surface and is therefore not suitable for your case.
A suitable representation for a 3D curve is in parametric form as a tuple (x(u), y(u), z(u))
where u
is some parameter and each coordinate is a function of u
.
The curve fitting problem is reduced to three 2D problems of fitting (ui, xi)
, (ui, yi)
, and (ui, zi)
separately (see also my answer here).
So, in order to perform a curve fit you are required to provide the parameterization ui
for each input point.
A common parameterization in spline fitting is the chord-length parameterization. This parameterization is defined by the accumulated length of the distances between the ordered points (u0=0, u1=|p1-p0|, u2 = u1+|p2-p1|
... etc.).
The following code implements a spline interpolation of your data with the chord-length parameterization.
from scipy import interpolate
xyz = np.vstack([x, y, z]).T
u = np.cumsum(np.r_[[0], np.linalg.norm(np.diff(xyz, axis=0), axis=1)])
# u is the chord-legth parameterization for each xyz point
sx = interpolate.InterpolatedUnivariateSpline(u, x) # x(u) spline
sy = interpolate.InterpolatedUnivariateSpline(u, y) # y(u) spline
sz = interpolate.InterpolatedUnivariateSpline(u, z) # z(u) spline
The code below samples the resulting spline and plots the result (in blue) on your data (black polyline). The result looks like this:
uu = np.linspace(u[0], u[-1], 100)
xx = sx(uu)
yy = sy(uu)
zz = sz(uu)
plt.plot(xx, yy, zz, "b")
You can also use the RBF
functions you suggested, for univariate interpolation.
The following code is an example of how to do this:
from scipy.interpolate import Rbf
rbfi_x = Rbf(u, x, function='cubic')
rbfi_y = Rbf(u, y, function='cubic')
rbfi_z = Rbf(u, z, function='cubic')
Sampling the resulting functions in a similar way to the spline sample above, and plotting the results on the previous figure, we get the plot below. As can be seen the spline and RBF interpolations are similar but not the same (for example near the endpoints).
Upvotes: 3