I want to perform cubic interpolation of an array that contains some nan values using scipy.griddata. However, as soon as a single nan is present in the values parameter, the returned interpolation is filled with only nan. This is not the case when using the 'nearest' or the 'linear' interpolation methods.
What is the reason for this behavior and is there a simple way to ignore nans in the values input?
This is a minimal working example, adapted from griddata scipy interpolation not working (giving nan):
import numpy as np
def func(x, y):
return x*(1-x)*np.cos(4*np.pi*x) * np.sin(4*np.pi*y**2)**2
grid_x, grid_y = np.mgrid[0:1:10j, 0:1:10j]
points = np.random.rand(100, 2)
values = func(points[:,0], points[:,1])
values[0]=np.nan # now add a single nan value to the array
from scipy.interpolate import griddata
grid_z0 = griddata(points, values, (grid_x, grid_y), method='nearest') # no nans here
grid_z1 = griddata(points, values, (grid_x, grid_y), method='linear') # this has nans on the edges (as expected)
grid_z2 = griddata(points, values, (grid_x, grid_y), method='cubic') # this is filled only with nans.
A solution is to remove all nan from the points and values input arrays prior to interpolating the data. numpy can be used efficiently for doing so:
import numpy as np
def func(x, y):
return x*(1-x)*np.cos(4*np.pi*x) * np.sin(4*np.pi*y**2)**2
grid_x, grid_y = np.mgrid[0:1:10j, 0:1:10j]
points = np.random.rand(100, 2)
values = func(points[:,0], points[:,1])
values[0]=np.nan # now add a single nan value to the array
#Find all the indexes where there is no nan neither in values nor in points.
nonanindex=np.invert(np.isnan(points[:,0]))*np.invert(np.isnan(points[:,1]))*np.invert(np.isnan(values))
#Remove the nan using fancy indexing. griddata can now properly interpolate. The result will have nan only on the edges of the array
from scipy.interpolate import griddata
grid_z2 = riddata(np.stack((points[nonanindex,0],points[nonanindex,1]),axis=1), values[nonanindex], (grid_x, grid_y), method='cubic')
While this solves the problem, I did not answer yet why this issue of the griddata function only appears for a cubic interpolation.