4

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. Bad Image

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

wbest
  • 611
  • 1
  • 6
  • 15

2 Answers2

7

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()

Interpolated missing values

Stefan van der Walt
  • 7,165
  • 1
  • 32
  • 41
  • This doesn't quite solve my problem (though it is in part my fault for not including an image), as I need to do more extrapolation than interpolation. But pre-calculating the grid coordinates is a good call. – wbest Feb 11 '14 at 16:03
  • What would a correct extrapolation look like? Why not simply mirror the image? Pretty much any extrapolation is "valid" outside the given image. – Stefan van der Walt Feb 11 '14 at 22:32
  • 1
    A nearest neighbor approach would be fine, except that that method results in striations that I don't want. I think I actually have a method working now that produces reasonable results wherein I fill the edges with the nanmean of all the square, and then use griddata to interpolate linearly. – wbest Feb 11 '14 at 23:09
4

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
wbest
  • 611
  • 1
  • 6
  • 15
  • 1
    This seems to be similar to one impainting method in OpenCV (INPAINT_TELEA): https://docs.opencv.org/3.4.0/df/d3d/tutorial_py_inpainting.html – moi Jan 26 '19 at 21:05