# 导入所需的库import xarray as xrimport numpy as npimport matplotlib.pyplot as pltimport cartopy.crs as ccrsimport cartopy.feature as cfeatureimport cartopy.mpl.ticker as ctickerimport cmapsfrom scipy.interpolate import make_interp_spline# 打开数据集dsu_file = r'E:\SXJ_Pre_data\ERA5_monthly\uwnd_12lev_mon_mean_0.5.nc'dsu = xr.open_dataset(dsu_file)# 筛选时间范围：1979-01-01 至 2020-12-31dsu = dsu.sel(time=slice('1979-01-01', '2020-12-31'))# 选择0-180E和0-90N范围的数据dsu = dsu.sel(longitude=slice(0, 180), latitude=slice(90, 0))# 提取200hPa层次的冬季（1月、2月、3月）风速数据（这里为了方便直接取1-3月，严格冬季为12-2月）u_200 = dsu['u'].loc[dsu.time.dt.month.isin([1, 2, 3])].loc['1979-01-01':'2020-03-31', 200]# 读取经纬度信息lonu = dsu['longitude']latu = dsu['latitude']# 对时间维度求平均u_200_mean = u_200.mean(axis=0)# 绘图部分fig = plt.figure(figsize=(10, 7))# 设置地图投影proj = ccrs.PlateCarree(central_longitude=0)leftlon, rightlon, lowerlat, upperlat = (0, 180, 0, 90)img_extent = [leftlon, rightlon, lowerlat, upperlat]# 创建地图ax1 = fig.add_subplot(111, projection=proj)ax1.set_extent(img_extent, crs=ccrs.PlateCarree())ax1.add_feature(cfeature.COASTLINE.with_scale('110m'), edgecolor='grey') # 添加海岸线# 设置经纬度刻度ax1.set_xticks(np.arange(leftlon, rightlon + 30, 30), crs=ccrs.PlateCarree())ax1.set_yticks(np.arange(lowerlat, upperlat + 10, 10), crs=ccrs.PlateCarree())# 设置经纬度格式lon_formatter = cticker.LongitudeFormatter()lat_formatter = cticker.LatitudeFormatter()ax1.xaxis.set_major_formatter(lon_formatter)ax1.yaxis.set_major_formatter(lat_formatter)plt.tick_params(labelsize=15)# 设置标题ax1.set_title('Clim. U200 ', loc='left', fontdict={'family': 'Times New Roman', 'size': 20})# 绘制风速场c1 = ax1.contourf(lonu, latu, u_200_mean, zorder=0, transform=ccrs.PlateCarree(), cmap=cmaps.WhiteYellowOrangeRed, extend='both', levels=np.arange(10, 90, 10))# 设置颜色条colorbar1 = fig.add_axes([ax1.get_position().x0, ax1.get_position().y0 - 0.075, ax1.get_position().width, 0.015])cb1 = fig.colorbar(c1, cax=colorbar1, orientation='horizontal', format='%d')cb1.ax.set_title('m/s', fontsize=13, pad=0.05)plt.tick_params(labelsize=15)# 绘制急流轴# 找到每个经度上风速最大的纬度max_speed_lat = u_200_mean.argmax(dim='latitude')x = lonu.valuesy = latu.values[max_speed_lat.values]# 生成平滑曲线model = make_interp_spline(x, y, k=3) # 三次样条插值x_new = np.linspace(x.min(), x.max(), 20)y_new = model(x_new)# 绘制平滑后的曲线ax1.plot(x_new, y_new, 'k-', lw=4)# 保存图像plt.savefig('drawmap.png', dpi=600, bbox_inches='tight')plt.show()
