气象编程 | Python基础地图构建（网转）
Python-basemap-中国南海小地图：
1 import numpy as np 2 import matplotlib.pyplot as plt 3 from mpl_toolkits.basemap import Basemap 4 import cmaps 5 import shapefile 6 from matplotlib.path import Path 7 from matplotlib.patches import PathPatch 8 import os 9 import maskout2 10 from netCDF4 import Dataset #netCDFを扱うため 11 nc = Dataset('F:/Rpython/lp14/MERRA2_400.tavgM_2d_slv_Nx.201912.nc4.nc4', mode='r') 12 lon = nc.variables["lon"][:].data 13 lat = nc.variables["lat"][:].data 14 t2m0= nc.variables["T2M"][:,:,:].data 15 t2m = t2m0.mean(axis=0) 16 grid_z = t2m -273.15 17 fig=plt.figure(figsize=(16,9)) 18 ax=fig.add_subplot(111) 19 mp=Basemap(llcrnrlon=70,llcrnrlat=15,urcrnrlon=140,urcrnrlat=55,projection='cyl') 20 mp.drawparallels(np.arange(15, 55 + 1, 10), labels = [1, 0, 0, 0],fontsize=14,linewidth='0.2',color='black') 21 mp.drawmeridians(np.arange(70, 140 + 1, 10), labels = [0, 0, 0, 1],fontsize=14,linewidth='0.2',color='black') 22 CHN='F:/Rpython/lp27/data/china-shapefiles-master/' 23 mp.readshapefile(CHN+'china_nine_dotted_line','china_nine_dotted_line',drawbounds=True,linewidth=1.5) 24 mp.readshapefile(CHN+'china','china',drawbounds=True) 25 clevs = np.arange(-50,50,10) 26 cf=plt.contourf(lon,lat,grid_z,clevs,cmap='Spectral_r',extend='both') 27 cbar=mp.colorbar(cf,location='right',size=0.3) 28 clip=maskout2.shp2clip(cf,ax,r'F:/Rpython/lp27/data/china-shapefiles-master/china_country') 29 #添加南海，实际上就是新建一个子图覆盖在之前子图的右下角 30 sub_ax = fig.add_axes([0.6888,0.1183,0.25,0.25]) #[*left*, *bottom*, *width*,*height*] 31 mf=Basemap(llcrnrlon=107,llcrnrlat=2,urcrnrlon=122,urcrnrlat=22,projection='cyl') 32 clevs2 = np.arange(-50,50,10) 33 cf2=sub_ax.contourf(lon, lat, grid_z,clevs2,cmap='Spectral_r') 34 CHN='F:/Rpython/lp27/data/china-shapefiles-master/' 35 mf.readshapefile(CHN+'china','china',drawbounds=True) 36 CHN='F:/Rpython/lp27/data/china-shapefiles-master/' 37 mf.readshapefile(CHN+'china_nine_dotted_line','china_nine_dotted_line',drawbounds=True,linewidth=1.5) 38 clip=maskout2.shp2clip(cf2,sub_ax,r'F:/Rpython/lp27/data/china-shapefiles-master/china_country') 39 plt.savefig(r"F:/Rpython/lp27/plot42.3.png",dpi=600) 40 plt.show()
Python-cartopy-中国南海小地图：
1 import maskout2 2 import os 3 import xarray as xr 4 import numpy as np 5 import cartopy.crs as ccrs 6 import cartopy.feature as cfeat 7 import matplotlib.colors as mcolors 8 from cartopy.io.shapereader import Reader 9 from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER 10 import matplotlib.pyplot as plt 11 import matplotlib.ticker as mticker 12 import cartopy.feature as cfeature 13 import cartopy.io.shapereader as shpreader 14 from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter 15 from matplotlib import rcParams 16 config = {"font.family":'Times New Roman', "font.size": 16, "mathtext.fontset":'stix'} 17 rcParams.update(config) 18 SHP = 'F:/Rpython/lp12/es' 19 # 数据读取 20 ds = xr.open_dataset(r'F:/Rpython/lp27/data/pr_Amon_FGOALS-f3-L_ssp370_r1i1p1f1_gr_201501-210012.nc') 21 pr = ds['pr'] 22 #取一段时间内的降水作为绘图的数据yuan 23 pr_0=pr[793:1032] 24 pr_ave=pr_0.mean(dim='time')*86400*30*12 25 lon=ds.lon 26 lat=ds.lat 27 lon_range=lon[(lon>70)&(lon<140)] 28 lat_range=lat[(lat>0)&(lat<60)] 29 pr_region=pr_ave.sel(lon=lon_range,lat=lat_range) 30 # 创建画图空间 31 proj = ccrs.PlateCarree() #创建投影 32 fig = plt.figure(figsize=(16,12)) #创建页面 33 ax = fig.subplots(1, 1, subplot_kw={'projection': proj}) #主图 34 # 设置网格点属性 35 region=[70, 140, 0, 60] 36 ax.set_extent(region,crs=proj) 37 ax.set_xticks(np.arange(region[0],region[1]+1,10),crs=proj) 38 ax.set_yticks(np.arange(region[-2],region[-1]+1,10),crs=proj) 39 #ax.grid(linestyle = '--') 40 ax.grid(linewidth=1.2,color='r',alpha=0.2,linestyle='--') 41 ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=False)) 42 ax.yaxis.set_major_formatter(LatitudeFormatter()) 43 # 画图 44 levels = np.arange(0,3000,500) 45 font2={'family':'SimHei','size':20,'color':'k'} 46 ax.set_title('2080-2100年均降水总量unit:mm',loc='left',fontdict=font2) 47 cs=ax.contourf(lon_range,lat_range,pr_region,levels=levels,cmap='Spectral_r',transform=ccrs.PlateCarree(),extend='both') 48 b=plt.colorbar(cs,shrink=0.88,orientation='vertical',extend='both',pad=0.035,aspect=20) 49 #白化 50 ax.add_geometries(Reader(os.path.join(SHP,'bou2_4l.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1) 51 clip=maskout2.shp2clip(cs,ax,r'F:/Rpython/lp3/hls/china0') 52 #添加南海，实际上就是新建一个子图覆盖在之前子图的右下角 53 f2_ax2=fig.add_axes([0.6566,0.133,0.1,0.17],projection=proj) 54 f2_ax2.set_extent([105,125,0,25],crs=ccrs.PlateCarree()) 55 f2_ax2.add_feature(cfeature.COASTLINE.with_scale('50m')) 56 china=shpreader.Reader('F:/Rpython/lp12/es/bou2_4l.dbf').geometries() 57 f2_ax2.add_geometries(china,ccrs.PlateCarree(),facecolor='none',edgecolor='black',zorder=1) 58 cs3=f2_ax2.contourf(lon_range,lat_range,pr_region,range(0,3000,500),cmap='Spectral_r') 59 #白化 60 f2_ax2.add_geometries(Reader(os.path.join(SHP,'bou2_4l.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1) 61 clip=maskout2.shp2clip(cs3,f2_ax2,r'F:/Rpython/lp3/hls/china0') 62 plt.savefig('F:/Rpython/lp27/plot8.7.png',dpi=300) 63 plt.show()#显示图片
