How to interpolate scattered data to a regular grid in Python?
It is straightforward to do so with numpy, scipy.interpolate.griddata, and matplotlib. Here is an example:

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

# data coordinates and values
x = np.random.random(100)
y = np.random.random(100)
z = np.random.random(100)

# target grid to interpolate to
xi = yi = np.arange(0,1.01,0.01)
xi,yi = np.meshgrid(xi,yi)

# set mask
mask = (xi &gt; 0.5) &amp; (xi &lt; 0.6) &amp; (yi &gt; 0.5) &amp; (yi &lt; 0.6)

# interpolate
zi = griddata((x,y),z,(xi,yi),method='linear')

# mask out the field
zi[mask] = np.nan

# plot
fig = plt.figure()
ax = fig.add_subplot(111)
plt.contourf(xi,yi,zi,np.arange(0,1.01,0.01))
plt.plot(x,y,'k.')
plt.xlabel('xi',fontsize=16)
plt.ylabel('yi',fontsize=16)
plt.savefig('interpolated.png',dpi=100)
plt.close(fig)


Result:



How to use this:


x and y are locations of points - these correspond to lon and lat values of your stations;
z are the values of points - this corresponds to your temperature observations from stations;
xi and yi are target grid axes - these will be your target longitude and latitude coordinates, which must match your landmask field;
zi is the result;
This example includes a simple way to mask the field. You should replace this mask with the landmask on your grid.


Notice also the method argument to griddata. Besides linear, this can also be cubic or nearest. I suggest you play with each to see what yields the best result for your dataset.