氣象編程 | Python基礎地圖構建（網轉）

Python-basemap-中國南海小地圖：

import numpy as np 
import matplotlib.pyplot as plt 
from mpl_toolkits.basemap import Basemap
import cmaps
import shapefile
from matplotlib.path import Path
from matplotlib.patches import PathPatch
import os
import maskout2
from netCDF4 import Dataset #netCDFを扱うため
nc = Dataset('F:/Rpython/lp14/MERRA2_400.tavgM_2d_slv_Nx.201912.nc4.nc4', mode='r')
lon = nc.variables["lon"][:].data
lat = nc.variables["lat"][:].data
t2m0= nc.variables["T2M"][:,:,:].data
t2m = t2m0.mean(axis=0)
grid_z = t2m -273.15
fig=plt.figure(figsize=(16,9))
ax=fig.add_subplot(111)
mp=Basemap(llcrnrlon=70,llcrnrlat=15,urcrnrlon=140,urcrnrlat=55,projection='cyl')
mp.drawparallels(np.arange(15, 55 + 1, 10), labels = [1, 0, 0, 0],fontsize=14,linewidth='0.2',color='black')
mp.drawmeridians(np.arange(70, 140 + 1, 10), labels = [0, 0, 0, 1],fontsize=14,linewidth='0.2',color='black')
CHN='F:/Rpython/lp27/data/china-shapefiles-master/'
mp.readshapefile(CHN+'china_nine_dotted_line','china_nine_dotted_line',drawbounds=True,linewidth=1.5)  
mp.readshapefile(CHN+'china','china',drawbounds=True)
clevs = np.arange(-50,50,10)
cf=plt.contourf(lon,lat,grid_z,clevs,cmap='Spectral_r',extend='both')
cbar=mp.colorbar(cf,location='right',size=0.3)
clip=maskout2.shp2clip(cf,ax,r'F:/Rpython/lp27/data/china-shapefiles-master/china_country')
#添加南海，實際上就是新建一個子圖覆蓋在之前子圖的右下角
sub_ax = fig.add_axes([0.6888,0.1183,0.25,0.25])    #[*left*, *bottom*, *width*,*height*]
mf=Basemap(llcrnrlon=107,llcrnrlat=2,urcrnrlon=122,urcrnrlat=22,projection='cyl')
clevs2 = np.arange(-50,50,10)
cf2=sub_ax.contourf(lon, lat, grid_z,clevs2,cmap='Spectral_r')
CHN='F:/Rpython/lp27/data/china-shapefiles-master/'
mf.readshapefile(CHN+'china','china',drawbounds=True)
CHN='F:/Rpython/lp27/data/china-shapefiles-master/'
mf.readshapefile(CHN+'china_nine_dotted_line','china_nine_dotted_line',drawbounds=True,linewidth=1.5)
clip=maskout2.shp2clip(cf2,sub_ax,r'F:/Rpython/lp27/data/china-shapefiles-master/china_country')
plt.savefig(r"F:/Rpython/lp27/plot42.3.png",dpi=600)
plt.show()

 

Python-cartopy-中國南海小地圖：

import maskout2
import os
import xarray as xr
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeat
import matplotlib.colors as mcolors
from cartopy.io.shapereader import Reader
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import matplotlib.pyplot as plt
import matplotlib.ticker as mticker
import cartopy.feature as cfeature
import cartopy.io.shapereader as shpreader
from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter
from matplotlib import rcParams
config = {"font.family":'Times New Roman', "font.size": 16, "mathtext.fontset":'stix'}
rcParams.update(config)
SHP = 'F:/Rpython/lp12/es'
# 數據讀取
ds = xr.open_dataset(r'F:/Rpython/lp27/data/pr_Amon_FGOALS-f3-L_ssp370_r1i1p1f1_gr_201501-210012.nc')
pr = ds['pr']
#取一段時間內的降水作為繪圖的數據yuan
pr_0=pr[793:1032]
pr_ave=pr_0.mean(dim='time')*86400*30*12
lon=ds.lon
lat=ds.lat
lon_range=lon[(lon>70)&(lon<140)]
lat_range=lat[(lat>0)&(lat<60)]
pr_region=pr_ave.sel(lon=lon_range,lat=lat_range)
# 創建畫圖空間
proj = ccrs.PlateCarree()  #創建投影
fig = plt.figure(figsize=(16,12))  #創建頁面
ax = fig.subplots(1, 1, subplot_kw={'projection': proj})  #主圖
# 設置網格點屬性
region=[70, 140, 0, 60]
ax.set_extent(region,crs=proj)
ax.set_xticks(np.arange(region[0],region[1]+1,10),crs=proj)
ax.set_yticks(np.arange(region[-2],region[-1]+1,10),crs=proj)
#ax.grid(linestyle = '--')
ax.grid(linewidth=1.2,color='r',alpha=0.2,linestyle='--')
ax.xaxis.set_major_formatter(LongitudeFormatter(zero_direction_label=False))
ax.yaxis.set_major_formatter(LatitudeFormatter())
# 畫圖
levels = np.arange(0,3000,500)
font2={'family':'SimHei','size':20,'color':'k'}
ax.set_title('2080-2100年均降水總量unit:mm',loc='left',fontdict=font2)
cs=ax.contourf(lon_range,lat_range,pr_region,levels=levels,cmap='Spectral_r',transform=ccrs.PlateCarree(),extend='both')
b=plt.colorbar(cs,shrink=0.88,orientation='vertical',extend='both',pad=0.035,aspect=20)
#白化
ax.add_geometries(Reader(os.path.join(SHP,'bou2_4l.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1)
clip=maskout2.shp2clip(cs,ax,r'F:/Rpython/lp3/hls/china0')
#添加南海，實際上就是新建一個子圖覆蓋在之前子圖的右下角
f2_ax2=fig.add_axes([0.6566,0.133,0.1,0.17],projection=proj)
f2_ax2.set_extent([105,125,0,25],crs=ccrs.PlateCarree())
f2_ax2.add_feature(cfeature.COASTLINE.with_scale('50m'))
china=shpreader.Reader('F:/Rpython/lp12/es/bou2_4l.dbf').geometries()
f2_ax2.add_geometries(china,ccrs.PlateCarree(),facecolor='none',edgecolor='black',zorder=1)
cs3=f2_ax2.contourf(lon_range,lat_range,pr_region,range(0,3000,500),cmap='Spectral_r')
#白化
f2_ax2.add_geometries(Reader(os.path.join(SHP,'bou2_4l.shp')).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=1)
clip=maskout2.shp2clip(cs3,f2_ax2,r'F:/Rpython/lp3/hls/china0')
plt.savefig('F:/Rpython/lp27/plot8.7.png',dpi=300)
plt.show()#顯示圖片