Python-气象-大气科学-可视化绘图系列(二)——利用basemap叠加地图,并添加白化效果(代码+示例)
本文为原创链接: https:////www.cnblogs.com/zhanling/p/12193031.html
1 白化单图代码: 2 import numpy as np 3 import xarray as xr 4 from mpl_toolkits.basemap import Basemap 5 import matplotlib.pyplot as plt 6 from matplotlib.patches import Polygon 7 import matplotlib.patches as mpatches 8 9 ds = xr.open_dataset('2019072300.006.nc') 10 t = ds['value'] 11 lons = ds.lon.data 12 lats = ds.lat.data 13 temp = xr.DataArray(t.data.T, coords=[lats,lons], dims=['latitude','longitude']) 14 # 创建画图空间 15 fig, ax = plt.subplots() 16 m = Basemap(projection='cyl',resolution='i',llcrnrlon=lons.min(),llcrnrlat=lats.min(), 17 urcrnrlon=lons.max(),urcrnrlat=lats.max(),lon_0=120.,lat_0=90) 18 Lon,Lat = np.meshgrid(lons[:],lats[:]) 19 X,Y = m(Lon,Lat) 20 21 shp_info3 = m.readshapefile("CHN_adm_shp\\CHN_adm3",'states',drawbounds=False,linewidth = 0.4,zorder=10) 22 for info, shp in zip(m.states_info, m.states): 23 proid = info['NAME_1'] # 可以用notepad打开CHN_adm1.csv文件,可以知道'NAME_1'代表各省的名称 24 if proid == 'Hubei': 25 poly = Polygon(shp,facecolor='None',edgecolor='b', lw=0.8) 26 ax.add_patch(poly) 27 else: 28 poly = Polygon(shp,facecolor='w',edgecolor='w', lw=0.8) 29 ax.add_patch(poly) 30 31 #市一级底图 32 shp_info2 = m.readshapefile("CHN_adm_shp\\CHN_adm2",'states',drawbounds=False,linewidth = 0.4,zorder=10) 33 for info, shp in zip(m.states_info, m.states): 34 proid = info['NAME_1'] 35 if proid == 'Hubei': 36 poly = Polygon(shp,facecolor='None',edgecolor='black', lw=0.8) 37 ax.add_patch(poly) 38 cs=m.contourf(Lon,Lat,t.data.T,np.arange(0,1.1,0.1),cmap='rainbow') 39 legend_elements = [mpatches.Patch(facecolor='#a50026',label='0-0.1'), 40 mpatches.Patch(facecolor='#da362a',label='0.1-0.2'), 41 mpatches.Patch(facecolor='#f67a49',label='0.2-0.3'), 42 mpatches.Patch(facecolor='#fdbf6f',label='0.3-0.4'), 43 mpatches.Patch(facecolor='#feeda1',label='0.4-0.5'), 44 mpatches.Patch(facecolor='#ebf7a3',label='0.5-0.6'), 45 mpatches.Patch(facecolor='#b7e075',label='0.6-0.7'), 46 mpatches.Patch(facecolor='#75c465',label='0.7-0.8'), 47 mpatches.Patch(facecolor='#249d53',label='0.8-0.9'), 48 mpatches.Patch(facecolor='#006837',label='0.9-1')] 49 plt.legend(handles=legend_elements, loc='lower right',title='NDVI') 50 # cbar = plt.colorbar(cs,orientation='horizontal',label='Potential')
个例效果:
1 多图代码: 2 import numpy as np 3 import xarray as xr 4 from mpl_toolkits.basemap import Basemap 5 import matplotlib.pyplot as plt 6 from matplotlib.patches import Polygon 7 import matplotlib.patches as mpatches 8 9 ds = xr.open_dataset('2019072300.006.nc') 10 t = ds['value'] 11 lons = ds.lon.data 12 lats = ds.lat.data 13 temp = xr.DataArray(t.data.T, coords=[lats,lons], dims=['latitude','longitude']) 14 # 创建画图空间 15 fig, ax = plt.subplots(2,2,figsize=(15,15)) 16 for ii in np.arange(4): 17 print('ii:::::',ii) 18 plt.sca(ax[np.unravel_index(ii,(2,2))]) # 获取一维索引在二维矩阵中的索引 19 m = Basemap(projection='cyl',resolution='i',llcrnrlon=lons.min(),llcrnrlat=lats.min(), 20 urcrnrlon=lons.max(),urcrnrlat=lats.max(),lon_0=120.,lat_0=90) 21 Lon,Lat = np.meshgrid(lons[:],lats[:]) 22 X,Y = m(Lon,Lat) 23 24 shp_info3 = m.readshapefile("CHN_adm_shp\\CHN_adm3",'states',drawbounds=False,linewidth = 0.4,zorder=10) 25 for info, shp in zip(m.states_info, m.states): 26 proid = info['NAME_1'] # 可以用notepad打开CHN_adm1.csv文件,可以知道'NAME_1'代表各省的名称 27 if proid == 'Hubei': 28 poly = Polygon(shp,facecolor='None',edgecolor='b', lw=0.8) 29 ax[np.unravel_index(ii,(2,2))].add_patch(poly) 30 else: 31 poly = Polygon(shp,facecolor='w',edgecolor='w', lw=0.8) 32 ax[np.unravel_index(ii,(2,2))].add_patch(poly) 33 34 #市一级底图 35 shp_info2 = m.readshapefile("CHN_adm_shp\\CHN_adm2",'states',drawbounds=False,linewidth = 0.4,zorder=10) 36 for info, shp in zip(m.states_info, m.states): 37 proid = info['NAME_1'] 38 if proid == 'Hubei': 39 poly = Polygon(shp,facecolor='None',edgecolor='black', lw=0.8) 40 ax[np.unravel_index(ii,(2,2))].add_patch(poly) 41 42 cs=m.contourf(Lon,Lat,t.data.T,np.arange(0,1.1,0.1),cmap='rainbow') 43 44 cbar = plt.colorbar(cs, ax=ax.ravel().tolist(), orientation='horizontal') 45 plt.savefig('multi_masked.jpg',bbox_inches = 'tight')
效果如下: