[AstroPy] LinearNDInterpolator

[Jonathan Slavin]

Erik,

Thanks for the response.  I didn't understand that the object returned
was a function which one then has to call with the desired coordinates.
I get it now.  

Unfortunately it seems to be not as fast as I'd hoped.  The simple
function that I wrote, which takes the "effective indices" for the
arrays and returns an array of interpolated values is much faster.
Mostly that's because it requires a 3-D grid for input so using the
effective indices is a simple vector operation.  I don't know how many
people would have use for my routine or if such a routine already exists
out there somewhere (it seems there must) but in any case I'd be happy
to make it available if anyone wants it.

Jon

> I'm not clear on what you're asking about LinearNDInterpolator in your
> second paragraph... what linear interpolator does is take in a list of
> coordinates and the corresponding values at those points (they can be
> gridded or not), and then can be used to give you the data value at
> any point inside the point values (although it does *not* extrapolate
> beyond the input points).  Is that what you want? If not, please
> clarify, or if so, see the example below for a sample usage - it makes
> a grid of points that represent a 3d gaussian and then uses the
> interpolator to do (I think) what you want.
> 
> 
> 
> 
> from numpy import mgrid,exp,array
> from scipy.interpolate import LinearNDInterpolator
> 
> #x,y,and z are 11x11x11 arrays with the x,y, and z values for a 3d box
> from -1 to 1
> x,y,z = mgrid[-1:1:11j,-1:1:11j,-1:1:11j]
> 
> #s is a gaussian evaluated at the grid points - replace this with your
> grid values
> s = exp(-(x**2+y**2+z**2)/2)
> 
> #this reshapes the points into the kind of input LinearNDInterpolator
> wants: Npts x Ndims
> pts = array((x.ravel(),y.ravel(),z.ravel())).T
> 
> #and the values should be a 1d array that matches the first dimension
> of pts
> sr = s.ravel()
> 
> i = LinearNDInterpolator(pts,sr)
> 
> 
> Now to get values at whatever points you want, you call the
> interpolator with an array of dimensions nyourpoints x ndims.  e.g.
> 
> >>> i([[0,0,0]])
> array([ 1.])
> >>> i([[0,0,0],[0,1,0],[1,1,0],[0,0,.5]])
> array([ 1.        ,  0.60653066,  0.36787944,  0.87919328])
> 
> Although note that if you go off the grid, you get NaN back:
> >>> i([[0,0,2]])
> array([ nan])
> 

> 
> On Fri, Apr 22, 2011 at 10:27 AM, Jonathan Slavin
> <jslavin at cfa.harvard.edu> wrote:
> > Hi all,
> >
> > I'm wondering if anyone has used the routine
> > scipy.interpolate.LinearNDInterpolator
> > >From what I can tell, it should do fast multidimensional linear
> > interpolations, which is something I could use.  My problem is that
> it
> > isn't well documented.
> >
> > What I have is a 3-D grid of data values and what I want is values
> at a
> > specific group of non-gridded points.  scipy.interpolate.griddata
> seems
> > to do the opposite to what I want -- given a set of (non-gridded)
> points
> > at which you know the data value, interpolate onto a regular grid.
>  It
> > seemed like LinearNDInterpolator was what I wanted, but it seems to
> want
> > the same input as griddata -- points at which the data value is
> known.
> > That is, it takes two inputs points and values where the number of
> > points must match the no. of values.  What I need is something more
> like
> > linearinterp(gridded_data, coords) -- that is, more like what
> > scipy.ndimage.map_coordinates does, but much faster and just linear
> > (rather than spline) interpolation.  I wrote a routine that does
> this,
> > and reasonably quickly, but it's pure python and not as fast as I'd
> > like.
> >
> > Jon
> >
> > _______________________________________________
> > AstroPy mailing list
> > AstroPy at scipy.org
> > http://mail.scipy.org/mailman/listinfo/astropy
> >
> 
> 
> 
> -- 
> Erik Tollerud
>