[SciPy-User] interpolation with inverse-distance weighting + KDTree

[denis]

Folks,
  here's a tiny class Invdisttree for interpolation with inverse-
distance weighting + KDTree.
It's solid, pretty fast, local, works for scattered data in any number
of dimensions,
leverages the excellent KDTree module.
Comments would be welcome; real test cases, 3+d, most welcome.

(For interpolating 2d data to a fine uniform grid,
matplotlib._delaunay.nn_interpolate_grid is ~ 10 times faster (on my
old mac ppc
One reason is that the dot( 1/dist, z[ix] ) takes over half the time;
another may be that nn_grid caches the current triangle ? )

cheers
  -- denis


    import numpy as np
    from scipy.spatial import cKDTree as KDTree

    class Invdisttree:
        """ inverse-distance-weighted interpolation using KDTree:
        invdisttree = Invdisttree( X, z )  -- points, values
        interpol = invdisttree( q, k=6, eps=0 )
            -- interpolate z from the 6 points nearest each q;
               q may be one point, or a batch of points
        """
        def __init__( self, X, z, leafsize=10 ):
            self.tree = KDTree( X, leafsize=leafsize )  # build the
tree
            self.z = z

        def __call__( self, q, k=6, eps=0 ):
                # k nearest neighbours of each query point --
            self.distances, self.ix = self.tree.query( q, k=k,
eps=eps )
            interpol = []  # np.zeros( (len(self.distances),) +
np.shape(z[0]) )
            for dist, ix in zip( self.distances, self.ix ):
                if dist[0] > 1e-10:
                    w = 1 / dist
                    wz = np.dot( w, self.z[ix] ) / np.sum(w)  # weight
z s by 1/dist
                else:
                    wz = self.z[ix[0]]
                interpol.append( wz )
            return interpol