I'm working on visualizing some climate model output. The computation is done on a projected latitude/longitude grid. Since the model is simulating sea ice, all land grid cells are masked. The standard tools for plotting geographical information in Python are Basemap and Cartopy, both of which use matplotlib routines. In particular, pcolormesh is the obvious choice for plotting. If there was no land mask, it would be simple:
X = longitude
Y = latitude
C = variable
fig, ax = plt.subplots()
plt.pcolormesh(X,Y,C)
Although C is allowed to be a masked array, pcolormesh can't handle masked arrays on X and Y. So how can I get around this?
For a simple example, take:
n = 100
X,Y = np.meshgrid(np.linspace(1,5,n),np.linspace(1,5,n))
C = np.sin(X*Y)
fig, ax = plt.subplots()
plt.pcolormesh(X,Y,C)
Now imagine we have a mask:
X[50:60,:] = np.nan
X[:,50:60] = np.nan
Y[50:60,:] = np.nan
Y[:,50:60] = np.nan
C[50:60,:] = np.nan
C[:,50:60] = np.nan
The first idea I had to get around this was to just select the valid entries and reshape X, Y, and C:
M = np.isnan(X)
X_valid = X[~M]
Y_valid = Y[~M]
C_valid = C[~M]
X_valid.shape = (81,100)
Y_valid.shape = (81,100)
C_valid.shape = (81,100)
plt.pcolormesh(X_valid, Y_valid, C_valid)
Like many naive approaches, this doesn't work.
Ideally the resulting plot would just be blank where the mask was. How can this be done?