Numpy inpaint nans interpolate AND extrapolate

Question

I'm working on a project using numpy and scipy and I need to fill in nanvalues. Currently I use scipy.interpolate.rbf, but it keeps causing python to crash so severely try/except won't even save it. However, after running it a few times, it seems as if it may keep failing in cases where there is data in the middle surrounded by all nans, like an island. Is there a better solution to this that won't keep crashing?

By the way, this is a LOT of data I need to extrapolate. Sometimes as much as half the image (70x70, greyscale), but it doesn't need to be perfect. It's part of an image stitching program, so as long as it's similar to the actual data, it'll work. I've tried nearest neighbor to fill in the nans, but the results are too different.

EDIT:

The image it always seems to fail on. Isolating this image allowed it to pass the image ONCE before crashing.

I'm using at least version NumPy 1.8.0 and SciPy 0.13.2.

Answer 1

Using SciPy's LinearNDInterpolator. If all images are of the same size, grid coordinates can be pre-computed and re-used.

import numpy as np
import matplotlib.pyplot as plt
from scipy import interpolate

x = np.linspace(0, 1, 500)
y = x[:, None]
image = x + y

# Destroy some values
mask = np.random.random(image.shape) > 0.7
image[mask] = np.nan

valid_mask = ~np.isnan(image)
coords = np.array(np.nonzero(valid_mask)).T
values = image[valid_mask]

it = interpolate.LinearNDInterpolator(coords, values, fill_value=0)

filled = it(list(np.ndindex(image.shape))).reshape(image.shape)

f, (ax0, ax1) = plt.subplots(1, 2)

ax0.imshow(image, cmap='gray', interpolation='nearest')
ax0.set_title('Input image')
ax1.imshow(filled, cmap='gray', interpolation='nearest')
ax1.set_title('Interpolated data')
plt.show()

Answer 2

This proved sufficient for my needs. It is actually quite fast and produces reasonable results:

ipn_kernel = np.array([[1,1,1],[1,0,1],[1,1,1]]) # kernel for inpaint_nans

def inpaint_nans(im):
    nans = np.isnan(im)
    while np.sum(nans)>0:
        im[nans] = 0
        vNeighbors = scipy.signal.convolve2d((nans==False),ipn_kernel,mode='same',boundary='symm')
        im2 = scipy.signal.convolve2d(im,ipn_kernel,mode='same',boundary='symm')
        im2[vNeighbors>0] = im2[vNeighbors>0]/vNeighbors[vNeighbors>0]
        im2[vNeighbors==0] = np.nan
        im2[(nans==False)] = im[(nans==False)]
        im = im2
        nans = np.isnan(im)
    return im