[SciPy-User] [SciPy-user] 3D interpolation
denis
denis-bz-gg@t-online...
Mon Dec 28 12:26:31 CST 2009
Vincent, you're right, splines can overshoot, though not by a lot:
0 1 1 0 overshoots by 20 % in the code snippet below.
Burak, does order=1, just averaging 8 neighbors, work ?
There are many kinds of cubic splines, including
-- global cubic, which map_coordinates seems to use
-- local B-spline: does not interpolate, go through the input points,
but can't overshoot, stays in the convex hull of the input points
-- local Catmull-Rom: interpolates, but can overshoot (over shoots and
leaves :)
("Local" means that e.g. in 1d, out[ 3 to 4 ] is calculated from
the 4 nearest input points at 2 3 4 5,
or equivalently that in[3] affects only out[ 1 to 5 ]: a desirable
property.)
Different corners of Scipy have different spline routines, some with
doc.
To see exactly what they do, you can only plot their impulse response
like so.
cheers
-- denis
""" plot map_coordinates( 1d spike )
"""
from __future__ import division
import sys
import numpy as np
from scipy.ndimage import map_coordinates
import pylab as pl
N = 10
exec( "\n".join( sys.argv[1:] )) # N= ...
np.set_printoptions( 2, threshold=100, suppress=True ) # .2f
a = np.r_[ 5*[0], 1., 5*[0], 1,1, 5*[0], 5*[1] ] # 0 1 1 0 -> 1.2
na = len(a) - 1
x = np.linspace( 0, na, na*N + 1 )
print "a:", a.astype(int)
for order in (3,):
z = map_coordinates( a, [x], order=order )
print "z[%d] halfint:" % order, z[N//2::N]
pl.plot( x, z )
pl.show()
More information about the SciPy-User
mailing list