[SciPy-Dev] Splines

Pablo Winant pablo.winant@gmail....
Thu Feb 21 04:33:33 CST 2013


Hi,

It's good to see people working in that direction. I will be more than 
happy to participate if I can.

As a user : I really wish I had a cubic spline interpolation with 
natural bounding conditions (which implies linear extrapolation). This 
is what is implemented in Matlab's griddedInterpolator and is missing 
from the existing options in scipy (right?). Being able to evaluate the 
derivatives is also a big advantage.

As a developper : there is the Einspline library that has a 
straightforward implementation in C. In particular, the representation 
of spline objects (C structs) and the low-level API are very 
straightforward and may be a good of inspiration. The library is 
currently GPL, but the author told me it could be made BSD if some m4 
macros are removed.
For what it's worth, I have made a simple Cython wrapper around some of 
its functions 
(https://github.com/albop/dolo/tree/master/dolo/numeric/interpolation) 
and was considering repackaging it. I had some other plans (like writing 
code for higher dimensions using some kind of templating, and updating 
SIMD evaluation). If there is a more elegant alternative, I would be 
happy to jump on the bandwagon.
I agree completely about the remarks on the spline format: it should be 
left as simple as possible. Having it isomorphic to C structs would be a 
good thing too as it permits easy experiments with extensions (for 
instance PyCuda kernels)

Best,

Pablo

On 21/02/2013 10:41, Pauli Virtanen wrote:
> Hi,
>
> Charles R Harris <charlesr.harris <at> gmail.com> writes:
>> There have been several threads on the list about splines
>> and consolidation of splines. For instance, there are several
>> uniform spline implementations for images and signal processing,
>> various low level functions in fitlib that are unexposed, and
>> perhaps useful altenatives to b-splines for some applications
>> like straight interpolation. For myself, I've started implementing
>> several functions in pure Python with an eye to converting them
>> to Cython once the interface and documentation is in place,
>> mostly for doing things that fitpack doesn't do because it is
>> very integrated, as opposed to supplying a basic toolset.
>> As part of this project, I'd like to get some feedback
>> on which functions people use most and what features they
>> would like to see. I'm not interested in the high level classes
>> at this point, either the current classes or combo functions
>> like fpcurf, but rather a collection of good lower level function
>> that could be used to implement the higher level functions
>> in a more flexible way. Thoughts?
> Great! I was going to start on this for 0.13.0, but this should
> speed things up considerably :)
>
> Overall, I think what we need is are (i) a well-specified spline
> format, and (ii) solid basic functions for evaluating and
> manipulating them, (iii) ditto for tensor product splines.
>
> How to construct the splines (interpolation, smoothing, etc.) should
> be considered as a separate problem. We can turn to FITPACK for
> smoothing splines, but it should not be used for interpolating
> splines.
>
>      ***
>
> Some misc thoughts on this:
>
> * The spline data format should be documented and set in
> stone as a first step. Users (and future developers) will
> want to toy around with it.
>
> Also, the data format for tensor product N-dim splines needs
> to be set. They are what we are missing, and what people are
> constantly asking for. We don't want them turn to
> `ndimage.map_coordinates`, which is clunky to use.
>
> The Fitpack tck format looks like this:
>
> https://github.com/pv/scipy-work/commits/spline-unify
>
> Currently, there's also a second B-spline data format used in
> scipy.interpolate with different padding, but we may want to
> stick with the FITPACK one, it's probably as good as any.
>
> * Functions for splines with varying x-coordinates are needed.
> Uniform grid splines would be a nice bonus as a speed-gain
> optimization.
>
> * The 1-D spline routines should be able to work over an
> arbitrary axis of multidimensional data. Even better if this
> can be done without reshaping and copying the input data
> (e.g. with Numpy iterators).
>
> This sounds like providing strided 1-D loops for heavy lifting,
> and bolting array iterators on top.
>
> * Functions for integration + differentiation of splines
> as as abstract objects would be useful. For efficiency,
> evaluation of derivatives & integrals probably might need
> to be provided separately.
>
> * For tensor product splines, evaluation on a scattered point
> set + on a grid would be useful.
>
> * It will probably be easiest to start from a clean slate,
> rather than trying to reuse scipy.interpolate.
>
> * FITPACK should not be used as a basis for this work, no
> need to use ancient FORTRAN 77 code for simple stuff. We can
> use its routines for generating smoothing splines, though.
>
> * Routines for constructing interpolating splines --- I think most
> of the time people will use these for simple gridded data interpolation
> rather than smoothing. FITPACK's knot selection is nice when it works,
> but often it doesn't (or requires careful fiddling), so we should have
> something simple and robust as a default algorithm.
>
> * I'm not sure what to do with the various boundary conditions
> and out-of-bounds value handling. It's probably best to leave room
> for various ways to do this...
>
> * Scipy also has two implementations of piecewise polynomials
> --- these should be consolidated into one, too.
>
> Cheers,
> Pauli
>
>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev@scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev



More information about the SciPy-Dev mailing list