白化_python_气象_

# coding=utf-8 import matplotlib.pyplot as plt import cartopy.crs as ccrs import cartopy.feature as cfeat # 世界地图信息 from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER import cartopy.mpl.ticker as cticker from matplotlib.font_manager import FontProperties import numpy import netCDF4 as nc import shapefile from matplotlib.path import Path from matplotlib.patches import PathPatch def shp2clip(originfig, ax, shpfile, region): sf = shapefile.Reader(shpfile) for shape_rec in sf.shapeRecords(): if shape_rec.record[3] == region: ####这里需要找到和region匹配的唯一标识符，record[]中必有一项是对应的。 vertices = [] codes = [] pts = shape_rec.shape.points prt = list(shape_rec.shape.parts) + [len(pts)] for i in range(len(prt) - 1): for j in range(prt[i], prt[i + 1]): vertices.append((pts[j][0], pts[j][1])) codes += [Path.MOVETO] codes += [Path.LINETO] * (prt[i + 1] - prt[i] - 2) codes += [Path.CLOSEPOLY] clip = Path(vertices, codes) clip = PathPatch(clip, transform=ax.transData) for contour in originfig.collections: contour.set_clip_path(clip) return clip font = FontProperties(fname=r"C:\\Windows\\Fonts\\impact.ttf") font1 = FontProperties(fname=u"C:\\Windows\\Fonts\\simkai.ttf") fig = plt.figure(figsize=(16, 9.6)) ax = fig.add_subplot(111) # ncdata=nc.Dataset(r'2004_500_grtuvw.nc') ncdata = nc.Dataset('pres.mon.ltm.nc') # data=ncdata.variables['t'][0,:,:] data = ncdata.variables['pres'][0, :, :] # exit() # lat=ncdata.variables['latitude'][:] # lon=ncdata.variables['longitude'][:] lat = ncdata.variables['lat'][:] lon = ncdata.variables['lon'][:] lon1, lon2 = lon[0], lon[-1] lat1, lat2 = lat[-1], lat[0] nx = data.shape[1] ny = data.shape[0] ax = fig.add_axes([0.1, 0.1, 0.8, 0.8], projection=ccrs.PlateCarree()) box = [lon1, lon2, lat1, lat2] # 画图区域 ax.set_extent(box, crs=ccrs.PlateCarree()) ax.add_geometries(Reader(os.path.join(SHP,'cnhimap.shp')).geometries(), ccrs.PlateCarree(),facecolor='none',edgecolor='k', linewidth=1) xx, yy = numpy.meshgrid(nx, ny) yy = yy[::-1] # 原来是这样的吗？放在同一文件夹下 ax.readshapefile('country1', 'whatevername', color='gray') minval, maxval = int(numpy.amin(data)), int(numpy.amax(data)) + 1 cs = ax.contourf(xx, yy, data, range(minval, maxval), cmap=plt.cm.get_cmap('jet')) bar = plt.colorbar(cs) bar.set_ticks(range(minval - 1, maxval, 40)) bar.set_ticklabels(range(minval - 1, maxval, 40)) clip = shp2clip(cs, ax, 'country1', 'China') plt.title(u'Python Super Mask', fontproperties=font, fontsize=40) lon1, lon2 = ax.set_xlim(70, 140) lat1, lat2 = ax.set_ylim(15, 55) signature = u"by 平流层的萝卜" plt.text(lon2 - (lon2 - lon1) * 3.0 / 10, lat2 + (lat2 - lat1) * 0.1 / 10, signature, fontproperties=font1, fontsize=25) plt.show() plt.savefig('111.png')