Reputation: 401
Apparently scipy
doesn't provide a method similar to MATLAB's spapi
that allows for setting the derivatives only at specific points, ensuring smooth transitions all along the interpolated spline. I'd like to find a way to do that in Python.
My goal is to make spline interpolation using Python 3, but setting the desired derivatives only at specific points. For example, say the array X defines position of an object, point to point, I'd want a spline to represent a feasible trajectory, but also ensuring that the speed (first derivative) is zero at both starting and ending points, other points having no constraints on the derivatives. In this example, I could also desire to have the acceleration (second derivative) equal to zero at the same points.
What I want is, in other words, an implementation of spline interpolation similar to the MATLAB's spapi
function, as referenced here:
spline = spapi(knots,x,y)
returns the spline f (if any) of order
k = length(knots) - length(x)
with knot sequence knots for which
(*) f(x(j)) = y(:,j), all j.
If some of the entries of
x
are the same, then this is taken in the osculatory sense, i.e., in the sense thatDm(j)f(x(j)) = y(:, j)
, withm(j) : = #{ i < j : x(i) = x(j) }
, andDmf
themth
derivative of f. Thus r-fold repetition of a sitez
inx
corresponds to the prescribing of value and the firstr – 1
derivatives of f atz
.
I know of the method from_derivatives
of the BPoly
class of scipy.interpolate
, but that poses a critical problem, being that when derivatives are not specified at any point, the algorithm does not guarantee smooth transition, as noted here. I also tried what was proposed here, but as expected, the same problem arises.
So, next I present simple reproductions of what I'm trying to achieve, successful or not.
spapi
without derivativesHere, example using spapi
method, in MATLAB, without setting any derivatives. Not as what I want since the derivatives are high at start and end points:
xi = linspace(0, 10, 6);
xnew = linspace(0, 10, 100);
yi = [0 2 1 4 2 0];
knots = optknt(xi, order);
ref_spline = spapi(knots, xi, yi);
spline = fnval(xnew, ref_spline);
der_spline = gradient(spline);
And the corresponding plot of reference points along with the interpolated spline:
And here the 1st order derivative of that spline:
spapi
with derivativesHere, example using spapi
method, in MATLAB, setting derivatives as 0
at both start and end points. Exactly the expected result, with smooth transitions all along the spline, and derivatives equal to 0
at start and end points:
xi = linspace(0, 10, 6);
xnew = linspace(0, 10, 100);
xder = [0 xi 10];
yi = [0 2 1 4 2 0];
ynew = [0 yi 0];
knots = optknt(xder, order);
ref_spline = spapi(knots, xder, ynew);
spline = fnval(xnew, ref_spline);
der_spline = gradient(spline);
And the plot of reference points along with the spline:
And here the 1st order derivative of that spline:
BPoly.from_derivatives
with derivativesHere, example using BPoly.from_derivatives
setting derivatives as 0
at both start and end points. Unsuccessful, even though the derivatives are 0
at start and end, smooth transition is not guaranteed along the spline:
ref_points = [0, 2, 1, 4, 2, 0]
time_vector = np.linspace(0, 10, 100)
time_points = np.linspace(0, 10, 6)
ref_complete = [[ref_points[j] if (i == 0) else 0 for i in range(2)] if (
(j == 0) or (j == len(ref_points) - 1)) else [ref_points[j]] for j in range(len(ref_points))]
ref_spline = BPoly.from_derivatives(time_points, ref_complete)
spline = ref_spline(time_vector)
der_spline = ref_spline.derivative(1)
der_y = der_spline(time_vector)
To clarify, the line that defines ref_complete
simply replaces the first and last elements of ref_points
by an array containing the original value at index 0
, and a 0
at index 1
:
>>> ref_points
[0, 2, 1, 4, 2, 0]
>>> ref_complete
[[0, 0], [2], [1], [4], [2], [0, 0]]
And the plot of reference points along with the spline:
And here the 1st order derivative of that spline:
Upvotes: 3
Views: 1547
Reputation: 19
If you just want to do cubic spline interpolation with zero derivatives at both ends this is enough (docs here)
from scipy.interpolate import CubicSpline
CubicSpline(time_points, ref_points, bc_type="clamped")
Upvotes: 1