Python气象绘图教程—（十九）剖面图

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.io.shapereader import Reader as shpreader
import xarray as xr
plt.rcParams['font.sans-serif']=['SimHei']#显示中文
filename=r'E:aaaaworld_geo.nc'#地形文件储存地址
f=xr.open_dataset(filename)#读取文件
print(f)#打印其属性

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.io.shapereader import Reader as shpreader
import xarray as xr
plt.rcParams['font.sans-serif']=['SimHei']#显示中文
#######################################
filename=r'E:aaaaworld_geo.nc'#地形文件地址
proj=ccrs.PlateCarree()#缩写投影
extent=[70,140,5,75]#绘图范围
def truncate_colormap(cmap, minval=0.0, maxval=1.0, n=128):
    new_cmap = mpl.colors.LinearSegmentedColormap.from_list('trunc({n},{a:.2f},{b:.2f})'.format(n=cmap.name,a=minval, b=maxval),cmap(np.linspace(minval, maxval, n)))
    return new_cmap#新的等高图色条，比较符合地理地形图的样子，来源气象家园巨佬
############################################
f=xr.open_dataset(filename)#读取文件
lon=f['x'][:]#将文件中的x变量赋值为经度
lat=f['y'][:]#赋值为纬度
height=f['z'][:]#将z变量赋值为高度
fig=plt.figure(figsize=(10,9),dpi=500)
ax=fig.add_subplot(projection=proj)
cmap_new = truncate_colormap(plt.cm.terrain, 0.23, 1.0) #截取colormap，要绿色以上的（>=0.23）
cmap_new.set_under([198/255,234/255,250/255]) #低于0的填色为海蓝
lev=np.arange(0,4000,200)
norm3 = mpl.colors.BoundaryNorm(lev, cmap_new.N) #标准化level，映射色标
cf=ax.contourf(lon[7500:9600],lat[2850:4950],height[2850:4950,7500:9600],levels=lev,norm=norm3,cmap=cmap_new,extend='both')
ax.set_extent(extent)
ax.set_title('一个框框地形图',fontsize=20)
plt.show()

ax.contourf(lon[7500:9600],lat[2850:4960],height[2850:4960,7500:9600],levels=lev,norm=norm3,cmap=cmap_new,extend='both')
ax.set_extent(extent)
############################################
ax.contourf(lon,lat,height,levels=lev,norm=norm3,cmap=cmap_new,extend='both')
ax.set_extent(extent)

import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.io.shapereader import Reader as shpreader
import matplotlib.ticker as mticker
import xarray as xr
plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus']=False
#######################################
filename=r'E:aaaaworld_geo.nc'
f=xr.open_dataset(filename)
lat=f['y'][3591:3621]
height=f['z'][3591:3621,8669]
fig=plt.figure(figsize=(4,1.5),dpi=700)
ax=fig.add_axes([0,0,1,1])
ax.plot(lat,height,c='k',lw=1)
ax.fill_between(lat,height,facecolor='white',hatch='///')#填充阴影
ax.set_xlim(29.7,30.6)
ax.set_xlabel('北纬(N)',fontsize=7)
ax.set_ylim(700,1650)
ax.set_ylabel('海拔高度(m)',fontsize=7)
ax.tick_params(which='both',labelsize=5)
plt.show()

lat=f['y'][3591:3621]
height=f['z'][3591:3621,8669]
print(lat)
print(height)
print(len(lat))
print(len(height))

import xarray as xr
ds = xr.open_dataset(r'E:aaaadatanchefnl_20200715_12_00.nc')
print(ds)

import xarray as xr
ds = xr.open_dataset(r'E:aaaadatanchefnl_20200715_12_00.nc')
single_temp=ds['TMP_P0_L1_GLL0'][:][:]#这是二维的
many_temp=ds['TMP_P0_L100_GLL0'][:][:][:]#这是三维的
single_many_temp=many_temp[0][:][:]#这就变为二维的，只取了一层次

import numpy as np
import matplotlib.pyplot as plt
import xarray as xr
import cartopy.crs as ccrs
import cartopy.feature as cf
import cartopy.io.shapereader as shpreader
from cartopy.mpl.ticker import LongitudeFormatter,LatitudeFormatter
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import matplotlib.ticker as mticker
plt.rcParams['font.sans-serif']=['SimHei']
filename=r'E:aaaadatanchefnl_20200715_18_00.nc'
f=xr.open_dataset(filename)
lon=f['lon_0'][:]#读取经度，一维的
lat=f['lat_0'][:]#读取纬度，一维的
RH_P0_L100_GLL0=f['RH_P0_L100_GLL0'][:][:][:]#读取相对湿度，三维的
UGRD_P0_L100_GLL0=f['UGRD_P0_L100_GLL0'][:][:][:]#读取纬向风，三维的
VGRD_P0_L100_GLL0=f['VGRD_P0_L100_GLL0'][:][:][:]#读取经向风，三维的
pres= f['lv_ISBL5'][:]*0.01#读取气压层数，一维的
fig=plt.figure(figsize=(5,3),dpi=700)#添加画布
ax=fig.add_axes([0,0,1,1])#添加子图
ax.invert_yaxis()#反转纵轴，使1000hPa作为起点
ax.set_yticks([1000,925,850,700,500,300])#突出我们感兴趣的气压层
ax.set_xticks(np.arange(28,36,1))
ax.set_xticklabels([r'28$^degree$N',r'29$^degree$N',r'30$^degree$N',r'31$^degree$N',r'32$^degree$N',
                  r'33$^degree$N',r'34$^degree$N',r'35$^degree$N'])#转换为纬度格式
ax.set(ylim=(1000,300))#设置气压层绘图范围，此处我们只显示到300hPa
ax.set_xlabel('纬度',fontsize=7)
ax.set_ylabel('层次（hPa）',fontsize=7)
ax.tick_params(axis='both',which='both',labelsize=7)
##################################################################################
ac=ax.contourf(lat[55:63],pres[:],RH_P0_L100_GLL0[:,55:63,109],
               cmap='Greens',levels=np.arange(60,101,10),extend='both',alpha=0.75)
ax.barbs(lat[55:63],pres[:],UGRD_P0_L100_GLL0[:,55:63,109],VGRD_P0_L100_GLL0[:,55:63,109],
         barb_increments={'half':2,'full':4,'flag':20},length=5)
cb=fig.colorbar(ac,extend='both',shrink=1,label='相对湿度（%）',pad=0.01)
cb.ax.tick_params(axis='both',which='both',length=1,labelsize=6)
ax.set_title('2020年7月15日20时相对湿度与风场剖面图',fontsize=15)

ax.contourf(lat[55:63],pres[:],RH_P0_L100_GLL0[:,55:63,109],
               cmap='Greens',levels=np.arange(60,101,10),extend='both',alpha=0.75)

ax.barbs(lat[55:63],pres[:],UGRD_P0_L100_GLL0[:,55:63,109],VGRD_P0_L100_GLL0[:,55:63,109],
         barb_increments={'half':2,'full':4,'flag':20},length=5)

import xarray as xr
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['font.sans-serif']=['SimHei']
filename=r'E:aaaaz1.nc'
f=xr.open_dataset(filename)
lat=f['lat'][:]
level=f['level'][:]
z=f['z'][:][:][:][:][:]
fig=plt.figure(figsize=(5,4),dpi=500)
ax=fig.add_axes([0,0,1,1])
ax.invert_yaxis()
ca=ax.contourf(lat[90:181],level[:],z[1,1,:,90:181,100],levels=np.arange(0,320000,10000))
fig.colorbar(ca)
ax.set(xlabel='北纬',ylabel='气压层（百帕）',title='多维数据的剖面图')
plt.show()