我正在尝试使用matplotlib绘制CMC grib2压力预测文件来绘制压力等高线。可以在这里找到grib2网格的描述:https://weather.gc.ca/grib/grib2_reg_10km_e.html。 grib2文件位于以下目录中:http://dd.weather.gc.ca/model_gem_regional/10km/grib2/00/000/,以CMC_reg_PRMSL_MSL_0_ps10km开头,后跟日期。它是一个包含平均海平面压力的grib文件。
我的问题是,我最终得到的一些直线轮廓遵循实际压力轮廓顶部的纬度线。我想这可能是因为我在PlateCarree中绘制而不是大地测量,但等高线图不允许使用大地测量。我的情节结果是:
代码如下:
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import datetime as dt
import cartopy
import cartopy.crs as ccrs
import Nio
gr = Nio.open_file('./data/CMC_reg_PRMSL_MSL_0_ps10km_2018111800_P000.grib2', 'r')
print(gr)
names = gr.variables.keys()
print("Variable Names:", names)
dims = gr.dimensions
print("Dimensions: ", dims)
attr = gr.attributes.keys()
print("Attributes: ", attr)
obs = gr.variables['PRMSL_P0_L101_GST0'][:]
lats = gr.variables["gridlat_0"][:]
lons = gr.variables["gridlon_0"][:]
fig = plt.figure(figsize=(15, 2))
intervals = range(95000, 105000, 400)
ax=plt.axes([0.,0.,1.,1.],projection=ccrs.PlateCarree())
obsobj = plt.contour(lons, lats, obs, intervals, cmap='jet',transform=ccrs.PlateCarree())
states_provinces = cartopy.feature.NaturalEarthFeature(
category='cultural',
name='admin_1_states_provinces_lines',
scale='50m',
facecolor='none')
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
ax.add_feature(states_provinces,edgecolor='gray')
obsobj.clabel()
colbar =plt.colorbar(obsobj)
任何建议,将不胜感激。
UPDATE
对于没有PyNIO的任何人,可以使用以下内容来使用评论部分中的转储文件进行重现。
只需删除所有对NIO的引用,并用以下内容替换lats,lons,obs赋值。
lats = np.load('lats.dump')
lons = np.load('lons.dump')
obs = np.load('obs.dump')
问题是网格缠绕在地球上。因此,在-180°的网格上将存在点,其近邻位于+ 180°,即网格环绕反射影。下面绘制了沿两个方向的网格索引。可以看到第一个网格行(黑色)出现在图的两侧。
因此,太平洋西部之后的等高线需要直接穿过地块,继续朝向地块另一侧的日本。这将导致不希望的线条
解决方案是屏蔽PlateCarree的外部点。那些发生在网格的中间。在经度大于179°或小于-179°的坐标处切割网格,以及将北极移出看起来像
其中蓝色表示切出点。
将其应用于等高线图给出:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
如果您将经度总和+180以避免负坐标,那么您的代码应该正在运行。从我的观点来看,坐标转换应该是合法的。