Python用Basemap画位势高度距平图和5880线

import xarray as xr import matplotlib.pyplot as plt import numpy as np from mpl_toolkits.basemap import Basemap plt.rcParams['font.sans-serif'] = ['KaiTi','SimHei'] plt.rcParams['font.size']= 12 plt.rcParams['axes.unicode_minus']=False hgt=r"D:\论文\论文数据\nc文件1998hgt.mon.mean.nc"  #阅读1998年目标文件 def DrawSeaLPJuly():     Hgt = xr.open_dataset(hgt)["hgt"]     time_label = (Hgt["time"].dt.month == [6,7,8])    选择6、7、8三个月     Hgt_main=Hgt[time_label].loc[dict(level=500)].mean("time").loc[60:0,70:180]  # 选定气压为500hpa，经纬度为0~60N，70~180E     lat = np.array(Hgt_main['lat'])  # 设置经纬度     lon=np.array(Hgt_main['lon'])  # 设置经纬度    #  算气候状态(1981~2010）30年：     i=1981     hgt1=r"D:\论文\论文数据\nc文件1998hgt.mon.mean.nc"     Hgt1= xr.open_dataset(hgt1)["hgt"]     time_label = (Hgt1["time"].dt.month == [6, 7, 8])     Hgt_total = Hgt1[time_label].loc[dict(level=500)].mean("time").loc[60:0,70:180]     i=i 1     while i<=2010:         hgt1=r"D:\论文\论文数据\nc文件" str(i) "hgt.mon.mean.nc"         Hgt1=xr.open_dataset(hgt1)["hgt"]         time_label = (Hgt1["time"].dt.month == [6,7,8])         Hgt_total = Hgt1[time_label].loc[dict(level=500)].mean("time").loc[60:0,70:180] Hgt_total         i=i 1     Hgt_total=Hgt_total/30  #求平均     con=Hgt_main-Hgt_total     X, Y = np.meshgrid(lon, lat)     fig = plt.figure(figsize=[12, 18])     ax = fig.add_subplot(111)     m = Basemap(lat_0=40, lon_0=130, projection="cyl",                 llcrnrlon=70,  # 地图左侧的经度                 llcrnrlat=0.0,  # 地图下方的纬度                 urcrnrlon=180.0,  # 地图上方的经度                 urcrnrlat=60.0,  # 地图上方的纬度                 resolution=None,  # 分辨率，这里设置为无                 )     m.readshapefile(r"D:\论文\论文数据\Shp\世界国家SHP\世界国家", "World Map", color="k", linewidth=1.2)  # 读取shp文件，想放多少就写多少这个代码     m.readshapefile(r"D:\论文\论文数据\Shp\国家水系矢量数据\hyd1_4l", "river Map", color="lightskyblue", linewidth=1)  #river Map是名字，可以随意更改，但必须有     levels = np.arange(-30, 70, 12)     c = m.contourf(X, Y, con,cmap="RdBu_r",alpha=0.7)  # 绘制等值线     a1=m.contour(X, Y, con, colors="k", linewidths=0.5)     # ax.clabel(a2, fmt="%2.0f", fontsize=10, colors="k")  #色标     a2=m.contour(X,Y, Hgt_main, range(5880, 5890, 10),linewidths=0.5,colors='r')   #画出特定的等值线，range()函数画出实际和关键线，间隔多少？     ax.clabel(a2, fmt="%2.0f", fontsize=10, colors="r")     ax.clabel(a1, fmt="%2.0f", fontsize=10, colors="g")  #标出值     # fig.colorbar(c, ax=ax, orientation='horizontal')     plt.yticks([0, 10, 20, 30, 40, 50, 60], ["EQ", "10°N", "20°N", "30°N", "40°N", "50°N", "60°N"], fontsize=16)  #设置x和y坐标轴     plt.xticks([70, 80, 90, 100, 110, 120, 130, 140, 150, 160, 170, 180],                ["70°E", "80°E", "90°E", "100°E", "110°E", "120°E", "130°E", "140°E", "150°E", "160°E", "170°E", "180°E"],                fontsize=16)     plt.ylabel("纬度", fontsize=16)     plt.xlabel("经度", fontsize=16)     m.drawparallels(np.arange(0, 60, 10))     m.drawmeridians(np.arange(70, 181, 10))  # 画出格线     plt.title("500hpa位置高度距平场", fontsize=18)  # 写上title      plt.show()   if __name__ == '__main__':     DrawSeaLPJuly()  

本文数据为NCEP/NCAR 再分析数据，感觉这真的是比较ERA5易于使用，即分辨率与内容类型的比较ERA5差。链接如下：https://psl.noaa.gov/data/gridded/data.ncep.reanalysis.html。本文比较初级，就是画了等值线还有一个画出目标区域的线，shp网上有很多文件。如果你想在评论区发邮箱，我可以给你我收集的一切。在文章中，我把一个大文件剪成30个小文件，没有必要，这个可以再优化一次。希望能给大家提供帮助。