在之前的博客。python极地极地投影绘制中,我曾经提到过python cartopy在极地投影中难以添加坐标标签的问题,当时解决方法是使用文本的方式添加,但这种方式的麻烦且并不统一适用。而在之后的可视化过程,我发现python在极地投影时仍会出现一些奇怪的问题:比如,等值线扭曲,出现不规则多边形,风场分布不均匀等等。
然而,这些问题都在最新版的cartopy21中解决,因此被cartopy极地投影折磨了数月而不得不转用m_map的我得到了解脱,下面将以极地极地春季的10m风场为例,给大家看看效果。
安装
最新版的cartopy似乎无法通过conda upgrde与conda install的方式直接安装,但是没关系,我们直接去官网下载安装包:cartopy21.0,根据自己的python版本与系统选择即可,这里我选的是win64+py39。
下载了安装包后,进入anaconda Prompt:
conda uninstall cartopy --force#只删除cartipy,不删除cartopy依赖conda install cartopy-0.21.0-py39h4915f10_0.tar.bz2
即可。
使用
绘制风场的代码与其他并没有区别,只是由于bug的修复,原来用于普通投影的代码可以正常使用。
import osimport matplotlib.ticker as mtickerimport netCDF4 as ncimport matplotlib.path as mpathimport cmapsimport matplotlib.pyplot as plt###引入库包import numpy as npimport matplotlib as mplimport cartopy.crs as ccrsimport cartopy.feature as cfeaturefrom netCDF4 import Datasetfrom wrf import getvar, interplevel, vertcross,vinterp, ALL_TIMES, CoordPair, xy_to_ll, ll_to_xy, to_np, get_cartopy, latlon_coords, cartopy_xlim, cartopy_ylimfrom cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTERd=nc.Dataset('F:/ERA5/met_em/use_1.nc')u=np.array(d.variables['UU'])v=np.array(d.variables['VV'])os.chdir('F:/wrfout/wrfout/')filepath='F:/wrfout/wrfout/'file_list = os.listdir(filepath)#数据导入files = [os.path.join(filepath,x) for x in file_list]wrflist=[Dataset(d) for d in file_list]times = to_np(getvar(wrflist, "times", timeidx=ALL_TIMES)) #lat2D = to_np(getvar(wrflist[10], "lat" )) # units: decimal degreeslon2D = to_np(getvar(wrflist[10], "lon" )) windspeed=np.sqrt(u**2 + v**2)#绘图mpl.rcParams["font.family"] = 'Arial' #默认字体类型mpl.rcParams["mathtext.fontset"] = 'cm' #数学文字字体mpl.rcParams["font.size"] = 12 #字体大小mpl.rcParams["axes.linewidth"] = 1 proj =ccrs.NorthPolarStereo(central_longitude=0)#设置地图投影#在圆柱投影中proj = ccrs.PlateCarree(central_longitude=xx)leftlon, rightlon, lowerlat, upperlat = (-180,180,60,90)#经纬度范围img_extent = [leftlon, rightlon, lowerlat, upperlat]fig1 = plt.figure(figsize=(12,10))#设置画布大小f1_ax1 = fig1.add_axes([0.2, 0.3, 0.5, 0.5],projection = ccrs.NorthPolarStereo(central_longitude=0))#绘制地图位置#注意此处添加了projection = ccrs.NorthPolarStereo(),指明该axes为北半球极地投影f1_ax1.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,linewidth=1, color='grey',linestyle='--')f1_ax1.set_extent(img_extent, ccrs.PlateCarree())f1_ax1.add_feature(cfeature.COASTLINE.with_scale('110m'))theta = np.linspace(0, 2*np.pi, 100)center, radius = [0.5, 0.5], 0.5verts = np.vstack([np.sin(theta), np.cos(theta)]).Tcircle = mpath.Path(verts * radius + center)f1_ax1.set_boundary(circle, transform=f1_ax1.transAxes)quiver = f1_ax1.quiver(lon2D, lat2D, musp, mvsp,windspeed, pivot='tail',width=0.002, scale=50, cmap=plt.cm.jet, headwidth=4,regrid_shape=35,alpha=1,transform=ccrs.PlateCarree())f1_ax1.quiverkey(quiver, 0.91, 1.03, 3, "3m/s",labelpos='E', coordinates='axes', fontproperties={'size': 10,'family':'Times New Roman'})plt.show()
可以看出并无不同,唯一值得一提的是regrid_shape选项,它可以帮助我们绘制风场时不那么密集。
效果:
可以看出比之前的bug都有所修复。