Reputation: 9762
My goal is to calculate a smooth trajectory passing through a set of points as shown below. I have had a look at the available methods of scipy.interpolate
, and also the scipy user guide. However, the choice of the right method is not quite clear to me.
What is the difference between
BSpline
,
splprep
,
splrep
,
UnivariateSpline
,
interp1d
,
make_interp_spline
and
CubicSpline
?
According to the documentation, all functions compute some polynomial spline function through a sequence of input points. Which function should I select? A) What is the difference between a cubic spline and a B-spline of order 3? B) What is the difference between splprep()
and splrep()
? C) Why is interp1d()
going to be deprecated?
I know I'm asking different questions here, but I don't see the point in splitting the questions up as I assume the answers will be related.
All in all, I find that the scipy.interpolate module is organized a little confusingly. I thought maybe I'm not the only one who has this impression, which is why I'm reaching out to SO.
Here's how far I've come. Below is some code that runs the different spline functions for some test data. It creates the figure below.
CublicSpline
and any of the B-spline methods?BSpline
object from the output of splprep
and splrep
, such that the results of splev()
and BSpline()
are equivalent. That way, we can convert the output of splrep
and splprep
into the object-oriented interface of scipy.interpolate()
.splrep
, UnivariateSpline
and make_interp_spline
lead to the same result. In my 2D data, I need to apply the interpolation independently per data dimension that it works. The convenience function interp1d
yields the same result, too. Related SO-question: Linksplprep
and splrep
seem unrelated. Even if I compute splprep
twice for every data axis independently (see p0_new), the result looks differently. I see in the docs that splprep computes the B-spline representation of an n-D curve. But should splrep
and splprep
not be related?splprep
, splrep
and UnivariateSpline
have a smoothing parameter, while other interpolators have no such parameter.splrep
pairs with UnivariateSpline
. However, I couldn't find a matching object-oriented counterpart for splprep
. Is there one?import numpy as np
from scipy.interpolate import *
import matplotlib.pyplot as plt
points = [[0, 0], [4, 4], [-1, 9], [-4, -1], [-1, -9], [4, -4], [0, 0]]
points = np.asarray(points)
n = 50
ts = np.linspace(0, 1, len(points))
ts_new = np.linspace(0, 1, n)
(t0_0,c0_0,k0_0), u = splprep(points[:,[0]].T, s=0, k=3)
(t0_1,c0_1,k0_1), u = splprep(points[:,[1]].T, s=0, k=3)
p0_new = np.r_[np.asarray(splev(ts_new, (t0_0,c0_0,k0_0))),
np.asarray(splev(ts_new, (t0_1,c0_1,k0_1))),
].T
# splprep/splev
(t1,c1,k1), u = splprep(points.T, s=0, k=3)
p1_new = splev(ts_new, (t1,c1,k1))
# BSpline from splprep
p2_new = BSpline(t1, np.asarray(c1).T, k=k1)(ts_new)
# splrep/splev (per dimension)
(t3_0,c3_0,k3_0) = splrep(ts, points[:,0].T, s=0, k=3)
(t3_1,c3_1,k3_1) = splrep(ts, points[:,1].T, s=0, k=3)
p3_new = np.c_[splev(ts_new, (t3_0,c3_0,k3_0)),
splev(ts_new, (t3_1,c3_1,k3_1)),
]
# Bspline from splrep
p4_new = np.c_[BSpline(t3_0, np.asarray(c3_0), k=k3_0)(ts_new),
BSpline(t3_1, np.asarray(c3_1), k=k3_1)(ts_new),
]
# UnivariateSpline
p5_new = np.c_[UnivariateSpline(ts, points[:,0], s=0, k=3)(ts_new),
UnivariateSpline(ts, points[:,1], s=0, k=3)(ts_new),]
# make_interp_spline
p6_new = make_interp_spline(ts, points, k=3)(ts_new)
# CubicSpline
p7_new = CubicSpline(ts, points, bc_type="clamped")(ts_new)
# interp1d
p8_new = interp1d(ts, points.T, kind="cubic")(ts_new).T
fig, ax = plt.subplots()
ax.plot(*points.T, "o-", label="Original points")
ax.plot(*p1_new, "o-", label="1: splprep/splev")
ax.plot(*p2_new.T, "x-", label="1: BSpline from splprep")
ax.plot(*p3_new.T, "o-", label="2: splrep/splev")
ax.plot(*p4_new.T, "x-", label="2: BSpline from splrep")
ax.plot(*p5_new.T, "*-", label="2: UnivariateSpline")
ax.plot(*p6_new.T, "+-", label="2: make_interp_spline")
ax.plot(*p7_new.T, "x-", label="3: CubicSpline")
#ax.plot(*p8_new.T, "k+-", label="3: interp1d")
#ax.plot(*p0_new.T, "k+-", label="3: CubicSpline")
ax.set_aspect("equal")
ax.grid("on")
ax.legend(loc='center left', bbox_to_anchor=(1, 0.5))
plt.show()
Upvotes: 3
Views: 1491
Reputation: 26030
What is the difference between BSpline, splprep, splrep, UnivariateSpline, interp1d, make_interp_spline and CubicSpline?
a BSpline
object represents a spline function in terms of knots t
, coefficients c
and degree k
. It does not know anything about the data x, y
. It is a low-level implementation object, on par with PPoly
--- think of PPoly vs BSpline as a change of basis.
make_interp_spline
constructs a spline (BSpline
) which passes through the data x
and y
-- it's an interpolating spline, so that spl(x) == y
. It can handle batches of data: if y is two-dimensional and the second dimension has length n
, make_interp_spline(x, y)
represents a stack of functions y_1(x), y_2(x), ..., y_n(x)
.
CubicSpline
is similar but more limited: it only allows cubic splines, k=3
. If you ever need to switch from cubics to e.g. linear interpolation, with make_interp_spline
you just change k=3
to k=1
at the call site. OTOH, CubicSpline
has some features BSpline does not: e.g. root-finding. It is a PPoly instance, not BSpline.
interp1d is not deprecated, it is legacy. A useful part of interp1d is nearest/previous/next modes. The rest just delegates to make_interp_spline, so better use that directly.
splrep
constructs a smoothing spline function given data. The amount of smoothing is controlled by the s
parameter, with s=0
being interpolation. It returns not a BSpline but a tck tuple (knots, coefficients and degree).
splprep
constructs a smoothing spline curve in a parametric form. That is (x(u), y(u))
not y(x)
. Also returns tck-tuples.
UnivariateSpline
is equivalent to splrep
.
Note that scipy.interpolate
has grown organically. With at least four generations of developers over almost quarter century. And the FITPACK library which splrep
, splprep
and UnivariateSpline
delegate to is from 1980s. CubicSpline
, BSpline
and make_interp_spline
are not using FITPACK.
So in new code I personally would recommend make_interp_spline
+ BSpline
; if you only need cubics, CubicSpline
works, too.
But why do I end up in different results if I use CublicSpline and any of the B-spline methods?
Boundary conditions. Both make_interp_spline
and CubicSpline
default to not-a-knot
, and you can change it.
splrep
et al only use not-a-knot
.
Upvotes: 3