MeteoInfoLab脚本示例：利用比湿、温度计算相对湿度(附中国领土完整的世界地图shp)

本文源码作者：王亚强 研究员 中国气象科学研究院

---

利用比湿和温度计算相对湿度的函数是qair2rh(qair, temp, press=1013.25)，三个参数分别是比湿、温度和气压，气压有一个缺省值1013.25，因此计算地面相对湿度的时候也可以不给气压参数。下面的例子打开CFS模式的地面结果，并从比湿和温度数组计算出相对湿度数组并绘图。

#Open and read data 
fn = 'F:/RMeteoInfo/data/data3/fnl_20170101_00_00.grib2' 
f = addfile(fn) 
temp = f['Temperature_surface'][0,::-1,:] 
sh = f['Specific_humidity_height_above_ground'][0,0,::-1,:] 
#Calculate relative humidity 
temp = temp - 273.15 
rh = meteo.qair2rh(sh, temp) 
figure(figsize=[700,400], newfig=False)
ax=axesm(tickfontsize=11,axis=True)
world = shaperead('F:/RMeteoInfo/data/map/country1.shp') 
lchina = shaperead('F:/RMeteoInfo/data/map/bou2_4p.shp')
lchina2 = shaperead('F:/RMeteoInfo/data/map/bou2_4l.shp')
geoshow(lchina,edgecolor='k',size=0.3)
bou1_layer=geoshow(lchina2,edgecolor='k',size=0.3)
china_layer = geoshow('china', visible=False)
levs = [40,50, 60, 70, 80,90,100,110]
layer = contourfm(rh, levs, cmap='grads_rainbow')
masklayer(lchina, [layer])
title(u'Relative humidity',bold=True,fontsize=11)
colorbar(layer,orientation='horizontal',ticklen=0,extendrect=False,shrink=1,pad=0.01,aspect=50,label='mm',bold=True,fontsize=11)
yaxis(tickvisible=True,location='left',tickwidth=2,linewidth=2,ticklength=3) #ticklength刻度线长度，tickwidth刻度线宽度，linewidth边框宽度
yaxis(tickvisible=False,location='right',tickwidth=2,linewidth=2,ticklength=4) #分别调试tick的宽度，边框线宽和tick的长度
xaxis(tickvisible=False,location='top',tickwidth=2,linewidth=2,ticklength=4)
xaxis(tickvisible=True,location='bottom',tickwidth=2,linewidth=2,ticklength=3)
xlim(70, 140)
ylim(15, 55)
xticks(arange(70, 140.1, 10),bold=True,fontsize=11)
yticks(arange(15, 55.1, 10),bold=True,fontsize=11)
xlabel('Longitude',bold=True,fontsize=11)
ylabel('Latitude',bold=True,fontsize=11)
ax.scale_bar(0.2,0.35,width=165,linewidth=1,bold=True,fontsize=10,bartype='alternating_bar')
ax.north_arrow(0.2,0.9,width=40,height=40,linewidth=1,bold=True,fontsize=10)
#Add south China Sea
axesm(position=[0.7305,0.238,0.15,0.2], axison=False, frameon=True)
geoshow(lchina,edgecolor='k',size=0.2)
bou1_layer=geoshow(lchina2,edgecolor='k',size=0.1)
levs = [40,50, 60, 70, 80,90,100,110]
layer = contourfm(rh, levs, cmap='grads_rainbow')
masklayer(lchina, [layer])
xlim(106, 123)
ylim(2, 23)
savefig('F:/RMeteoInfo/test/plot127.1.png',dpi=800)

---

众所周知，Matlab，R语言，Python，IDL等国外科研作图软件自带的中国地图对我们“十分的不友好”。本文分享的世界地图是以0°经线为中心，其中印边界藏南段、阿克赛钦段、台湾诸岛、克什米尔地区、科索沃地区等边界线与我国立场保持高度一致，使用时仍需要添加南海九段线及南海诸岛（2020年6月中国省、市、县行政区划shp）。

中国领土神圣不可侵犯，中国一点都不能少！

---

3 获取本文中国领土完整的世界地图shp数据的途径：

气象水文科研猫公众号后台回复：

“中国领土完整的世界地图shp”，

获取百度云免费下载链接。

---