[SciPy-User] Quick interpolation question

[David Baddeley]

Hi Jose,

am I correct in assuming that you are trying to recover missing data points on 
an existing regular grid? Depending on what fraction of your data is missing, 
you might be able to do something based on convolution - for example (2D, but 
easily generalisable):

#set up a convolution kernel
kernel = np.array([[1,1,1],[1,0,1],[1,1,1]], 'f')  # could also use a less 
connected version, eg [[0,1,0],[1,0,1], [0,1,0]]

#normalise so integral is 1
kernel = kernel/kernel.sum()

#start with 0 for all missing values
f = data*flags

#convolve with the kernel - will replace each point with the mean of it's 
neighbours
f_smooth = scipy.ndimage.convolve(f, kernel)

#obtain a new estimate for the data by replacing the missing values with those 
from the smoothed version
f = data*flags + f_smooth*(1-flags)

this will effectively fix any isolated missing values. If you have larger 
missing regions, you should be able to reach a solution by iterating these last 
two steps - stop when the solution doesn't change. As convolution with such a 
small kernel is fast, the method should be fast for small missing regions, but 
will tend to blow up if the regions get too large & the iteration number goes 
up. 

cheers,
David




________________________________
From: Jose Gomez-Dans <jgomezdans at gmail.com>
To: SciPy Users List <scipy-user at scipy.org>
Sent: Tue, 5 April, 2011 6:59:22 AM
Subject: [SciPy-User] Quick interpolation question

Hi,

I just want to interpolate, in the simplest possible terms, a 3D dataset. Linear 
interpolation, nearest neighbour, all that would suffice (this is to start off 
some algorithm, so no accurate estimate is required). In new scipy versions, 
things like griddata would be useful, but currently I only have scipy 0.8. So I 
have a "cube" (data[:,:,:], (NixNjxNk)) array, and an array of flags 
(flags[:,:,:,], True or False) of the same size. I want to interpolate my data 
for the elements of data where the corresponding element of  flag is False, 
using eg the nearest valid datapoint in data, or some linear combination of 
"close by" points.

I sort of think that map_coordinates would work, but I would need to filter 
data, and recalculate all the locations, so not efficient. Other ideas are to 
use np.interp1d and do it dimension by dimension. Would work, but clumsy.

In essence, what I need is some "gapfilling" simple multi-dimensional 
interpolation that works fast (Ni, Nj and Nk are pretty large, >500 each). Any 
suggestions?
Thanks!

J



      
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.scipy.org/pipermail/scipy-user/attachments/20110404/903a99e6/attachment.html>