Home

Awesome

Dierckx.jl

Julia library for 1-d and 2-d splines

Build Status Coverage Status

This is a Julia wrapper for the dierckx Fortran library, the same library underlying the spline classes in scipy.interpolate. Some of the functionality here overlaps with Interpolations.jl, a pure-Julia interpolation package. Take a look at it if you have a use case not covered here.

All new development on Dierckx.jl will be for Julia v1.3 and above. The master branch is therefore incompatible with earlier versions of Julia.

Features

Install (Julia 1.3 and later)

(v1.3) pkg> add Dierckx

(Type ] to enter package mode.)

Example Usage

using Dierckx

Fit a 1-d spline to some input data (points can be unevenly spaced):

x = [0., 1., 2., 3., 4.]
y = [-1., 0., 7., 26., 63.]  # x.^3 - 1.
spl = Spline1D(x, y)

Evaluate the spline at some new points:

spl([1.5, 2.5])  # result = [2.375, 14.625]
spl(1.5)  # result = 2.375

Equivalent to the above:

evaluate(spl, [1.5, 2.5])
evaluate(spl, 1.5)

Evaluate derivative, integral, or roots:

derivative(spl, 1.5)  # derivate at x=1.5; result is 5.75
integrate(spl, 0., 4.)  # integrate from x=0 to x=4; result is 60.0
roots(spl)  # result is [1.0]

Note that roots() only works for cubic splines (k=3).

Fit a 2-d spline to data on a (possibly irregular) grid:

x = [0.5, 2., 3., 4., 5.5, 8.]
y = [0.5, 2., 3., 4.]
z = [1. 2. 1. 2.;  # size is (length(x), length(y))
     1. 2. 1. 2.;
     1. 2. 3. 2.;
     1. 2. 2. 2.;
     1. 2. 1. 2.;
     1. 2. 3. 1.]

spline = Spline2D(x, y, z)

Note that if you consider z as a matrix, x refers to row coordinates and y refers to column coordinates.

Evaluate at element-wise points:

xi = [1., 1.5, 2.3, 4.5, 3.3, 3.2, 3.]
yi = [1., 2.3, 5.3, 0.5, 3.3, 1.2, 3.]
zi = spline(xi, yi)  # 1-d array of length 7
zi = evaluate(spline, xi, yi)  # equivalent to previous line

Evaluate at grid spanned by input arrays:

xi = [1., 1.5, 2.3, 4.5]
yi = [1., 2.3, 5.3]
zi = evalgrid(spline, xi, yi)  # 2-d array of size (4, 3)

Reference

1-d Splines

Spline1D(x, y; w=ones(length(x)), k=3, bc="nearest", s=0.0)
Spline1D(x, y, xknots; w=ones(length(x)), k=3, bc="nearest")
evaluate(spl, x)
derivative(spl, x; nu=1)
integrate(spl, a, b)
roots(spl; maxn=8)

Parametric Splines

These are the B-spline representation of a curve through N-dimensional space.

ParametricSpline(X; s=0.0, ...)
ParametricSpline(u, X; s=0.0, ...)
ParametricSpline(X, knots, ...)
ParametricSpline(u, X, knots, ...)

Keyword arguments common to all constructor methods:

2-d Splines

Spline2D(x, y, z; w=ones(length(x)), kx=3, ky=3, s=0.0)
Spline2D(x, y, z; kx=3, ky=3, s=0.0)
evaluate(spl, x, y)
evalgrid(spl, x, y)
derivative(spl, x, y; nux = 1, nuy = 1)
integrate(spl, x0, x1, y0, y1)

Translation from scipy.interpolate

The Spline classes in scipy.interpolate are also thin wrappers for the Dierckx Fortran library. The performance of Dierckx.jl should be similar or better than the scipy.interpolate classes. (Better for small arrays where Python overhead is more significant.) The equivalent of a specific classes in scipy.interpolate:

scipy.interpolate classDierckx.jl constructor method
UnivariateSplineSpline1D(x, y; s=length(x))
InterpolatedUnivariateSplineSpline1D(x, y; s=0.0)
LSQUnivariateSplineSpline1D(x, y, xknots)
SmoothBivariateSplineSpline2D(x, y, z; s=length(x))
LSQBivariateSpline
RectBivariateSplineSpline2D(x, y, z; s=0.0) (z = 2-d array)
SmoothSphereBivariateSpline
LSQSphereBivariateSpline
RectSphereBivariateSpline

Parametric splines:

scipy.interpolate functionDierckx.jl constructor method
splprep(X)ParametricSpline(X)
splprep(X, u=...)ParametricSpline(u, X)
splprep(X, t=...)ParametricSpline(X, t) (t = knots)
splprep(X, u=..., t=...)ParametricSpline(u, X, t)

License

Dierckx.jl is distributed under a 3-clause BSD license. See LICENSE.md for details. The real*8 version of the Dierckx Fortran library as well as some test cases and error messages are copied from the scipy package, which is distributed under this license.

Citation

If you use this package in a pulication and wish to cite it, you may want to cite the original Fortran dierckx library:

Paul Dierckx, Curve and Surface Fitting with Splines, Oxford University Press, 1993

If convenient, you can also include a link to the github repository for this package: https://github.com/kbarbary/Dierckx.jl.