Python兰伯特投影中国区域等值线图(含南海小地图)

2021-01-04 14:52:36 浏览数 (1)


自定义兰伯特投影:

原作者:“坎坷”大佬


PlateCarree (无坐标转换)作图:

代码调试作者:气象水文科研猫

注:因小编时间有限,代码未进行精简。

代码语言:javascript复制
import numpy as np
import xarray as xr
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
import maskout2
import maskout
import lambert_ticks
from scipy.interpolate import Rbf
from cartopy.io.shapereader import Reader
from cartopy.mpl.gridliner import LATITUDE_FORMATTER, LONGITUDE_FORMATTER
from matplotlib import rcParams
config={"font.family":'Times New Roman',"font.size":16,"mathtext.fontset":'stix'}
rcParams.update(config)
plt.rcParams['figure.figsize']=(15,15)
nc_data = xr.open_dataset('F:/Rpython/lp28/data/r160.nc')
olon = nc_data.coords["station_lon"]
olat = nc_data.coords["station_lat"]
stid = nc_data.coords["station_id"]
stid_data = np.random.rand(len(olon)) * 500
nlon = np.linspace(70, 140, 140)
nlat = np.linspace(15, 55, 60)
nlon, nlat = np.meshgrid(nlon, nlat)
func = Rbf(olon, olat, stid_data, function='linear')
stid_pre = func(nlon, nlat)
# 地图
area_str = ["China"]
shp_add=r'F:/RMeteoInfo/data/map/bou2_4l.shp'
shp_mask=r'F:/RMeteoInfo/data/map/country1.shp'
proj = ccrs.LambertConformal(central_latitude=90, central_longitude=104)
fig, ax = plt.subplots(figsize=(15, 15), subplot_kw=dict(projection=proj))  # 建立页面
ax.set_extent([76, 132, 16, 54.5],ccrs.PlateCarree())  # 设置经纬度范围
ax.add_geometries(Reader(shp_add).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=0.7)
# 等值线图
con=ax.contourf(nlon,nlat,stid_pre,levels=np.arange(0,2000,200),cmap='gist_rainbow',extend='both',transform=ccrs.PlateCarree())
maskout.maskout_areas(originfig=con, ax=ax, shp_path=shp_mask, area_str=area_str, proj=proj)
b=plt.colorbar(con,shrink=0.6,orientation='vertical',extend='both',pad=0.035,aspect=20)
# 站点分布
point_locs_f = maskout.maskout_points(lon=olon, lat=olat, shp_path=shp_mask, area_str=area_str)
lons_f = np.array([olon[point_loc_f] for point_loc_f in point_locs_f])
lats_f = np.array([olat[point_loc_f] for point_loc_f in point_locs_f])
stids_f=np.array([stid[point_loc_f] for point_loc_f in point_locs_f])
ax.scatter(lons_f[::3], lats_f[::3], c='k', s=25, marker='o', transform=ccrs.PlateCarree())
# 南海小地图
ax_n = fig.add_axes([0.643,0.2645,0.15,0.09], projection=ccrs.PlateCarree())
ax_n.set_extent([105, 123, 2, 23], ccrs.PlateCarree())
ax_n.add_geometries(Reader(shp_add).geometries(),ccrs.PlateCarree(),facecolor='none',edgecolor='k',linewidth=0.7)
# 等值线图
con2=ax_n.contourf(nlon,nlat,stid_pre,levels=np.arange(0,2000,200),cmap='gist_rainbow',extend='both',transform=ccrs.PlateCarree())
clip2=maskout2.shp2clip(con2,ax_n,r'F:/Rpython/lp3/hls/china0') #白化1
# 坐标与经纬网格(兰伯特投影)
fig.canvas.draw()
xticks = list(range(-180, 180, 4))
yticks = list(range(-90, 90, 3))
ax.gridlines(xlocs=xticks, ylocs=yticks, linewidth=1.2, linestyle='--')
ax.xaxis.set_major_formatter(LONGITUDE_FORMATTER)
ax.yaxis.set_major_formatter(LATITUDE_FORMATTER)
lambert_ticks.lambert_xticks(ax, xticks)
lambert_ticks.lambert_yticks(ax, yticks)
plt.tick_params(labelsize=13)
font2={'family':'SimHei','size':24,'color':'k'}
ax.set_title("中国降水量分布",fontdict=font2)
plt.savefig('F:/Rpython/lp28/plot13.3.png',dpi=600)
plt.show()

0 人点赞