在 Basemap 中使用 maskoceans 时,世界的一半被遮盖了
Half of the world masked when using maskoceans in Basemap
我想在绘制来自 netCDF 数据集的数据时遮盖海洋。我遵循了 给出的重要指示。它适用于半个世界,但不知何故,格林威治以西的一切都被掩盖了,包括海洋和陆地。
这是我的代码:
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()
我知道有人提出了一个有点类似的问题 here 但从未得到答复,最后看来,问题是将经纬度坐标与 x-y 坐标混淆了。我尝试切换到 x-y 坐标,但得到了相同的半图。知道这里会发生什么吗?
N.B. 当使用 datapc = map1.pcolormesh(lons, lats, data[season,:,:], vmin=vmin, vmax=vmax, cmap=cmap, latlon=True)
绘制未屏蔽数据时,绘制了整个世界(陆地 + 海洋)。
如您所见,未绘制经度为 -180 到 0 的点。假设它们在您的数据中,它们必须出于某种原因被屏蔽或丢弃。
我的直觉是数据集经度 运行 0-360 而不是 -180 到 180,即 in the .
快速解决这个问题的方法是添加
lons_1d[lons_1d>180]-=360
就在你从 nc
中取出 lons_1d
之后。这是有效的,因为 lons_1d
是一个 numpy 数组,它使用 numpy boolean array indexing(通常称为 "fancy" 索引)有条件地 select 大于 180 的经度值并从中减去 360。
正如您所注意到的,如果您省略掩码,pcolormesh
绘图会起作用,这看起来像是 maskoceans
函数中包含包装的错误,或者至少是意外行为。
供参考 - 我不认为你是第一个遇到类似 "wrapping" 类型问题的面具,我认为这个 issue on the matplotlib github 看起来很相似。
我想在绘制来自 netCDF 数据集的数据时遮盖海洋。我遵循了
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()
我知道有人提出了一个有点类似的问题 here 但从未得到答复,最后看来,问题是将经纬度坐标与 x-y 坐标混淆了。我尝试切换到 x-y 坐标,但得到了相同的半图。知道这里会发生什么吗?
N.B. 当使用 datapc = map1.pcolormesh(lons, lats, data[season,:,:], vmin=vmin, vmax=vmax, cmap=cmap, latlon=True)
绘制未屏蔽数据时,绘制了整个世界(陆地 + 海洋)。
如您所见,未绘制经度为 -180 到 0 的点。假设它们在您的数据中,它们必须出于某种原因被屏蔽或丢弃。
我的直觉是数据集经度 运行 0-360 而不是 -180 到 180,即
快速解决这个问题的方法是添加
lons_1d[lons_1d>180]-=360
就在你从 nc
中取出 lons_1d
之后。这是有效的,因为 lons_1d
是一个 numpy 数组,它使用 numpy boolean array indexing(通常称为 "fancy" 索引)有条件地 select 大于 180 的经度值并从中减去 360。
正如您所注意到的,如果您省略掩码,pcolormesh
绘图会起作用,这看起来像是 maskoceans
函数中包含包装的错误,或者至少是意外行为。
供参考 - 我不认为你是第一个遇到类似 "wrapping" 类型问题的面具,我认为这个 issue on the matplotlib github 看起来很相似。