I would like to mask the oceans while building data from the netCDF dataset. I followed the excellent instructions provided in the answer to this question . It works great for half the world, but for some reason everything west of Greenwich disguises both the ocean and the earth. Here is my code:
import netCDF4
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import matplotlib.cm as cm
import mpl_toolkits
from mpl_toolkits import basemap
from mpl_toolkits.basemap import Basemap, maskoceans
filename = 'myfile.nc'
vmin = 0.
vmax = 1
nc = netCDF4.Dataset(filename, 'r')
data = nc.variables['sum'][:]
lats_1d = nc.variables['lat'][:]
lons_1d = nc.variables['lon'][:]
lons, lats = np.meshgrid(lons_1d, lats_1d)
labels = ['DJF', 'MAM', 'JJA', 'SON']
cmap = cm.RdYlBu
cmap.set_over('#00FF00')
my_dpi = 96
fig = plt.figure(figsize=(1200/my_dpi, 800./my_dpi))
for season in range(4):
ax = fig.add_subplot(2, 2, season+1)
map1 = basemap.Basemap(resolution='c', projection='kav7', lon_0=0)
map1.drawcoastlines()
map1.drawcountries()
nc_new = maskoceans(lons,lats,data[season,:,:],resolution='c', grid = 1.25)
datapc = map1.pcolormesh(lons, lats, nc_new, vmin=vmin, vmax=vmax, cmap=cmap, latlon=True)
plt.title(labels[season])
fig.tight_layout(pad=1, w_pad=1, h_pad=4)
ax = fig.add_axes([0.05, 0.52, 0.9, 0.025])
cb = plt.colorbar(cax=ax, orientation='horizontal', cmap=cmap,
extend='max', format="%.2f",
ticks=[0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1])
plt.show()
I know that a somewhat similar problem was raised here , but did not receive an answer, and it seems that in the end the problem was mixing lat-long coordinates with xy. I tried switching to xy coordinates, but got the same half of the map. Any idea what could happen here?
N.B. datapc = map1.pcolormesh(lons, lats, data[season,:,:], vmin=vmin, vmax=vmax, cmap=cmap, latlon=True) ( + ).