import numpy as np
from scipy.interpolate import interp1d, interp2d, UnivariateSpline, griddata
import matplotlib.pyplot as plt
np.random.seed(114514)
x0 = np.random.uniform(-3, 3, 50)
y0 = np.random.uniform(-4, 4, 50)
f = lambda x, y: (x**2 - 2*x)*np.exp(-x**2 - y**2 - x*y)
z0 = f(x0, y0)
xy0 = np.vstack([x0, y0]).T
x = np.linspace(x0.min(), x0.max(), 300)
y = np.linspace(y0.min(), y0.max(), 300)
X, Y = np.meshgrid(x, y)
z = griddata(xy0, z0, (X,Y), 'cubic')
zl = griddata(xy0, z0, (X,Y), 'linear')
z[np.isnan(z)] = zl[np.isnan(z)]
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(X, Y, z, cmap='hot')
fig.show()