气象编程 | 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()#显示图片