7 o. ]. K8 u& I
摘自《MeteoInfo气象GIS、科学计算与可视化平台》182页。 , o, b0 H" T+ m
水汽通量散度是分析水气循环和降水常用的指标之一,需要用到风场、温度和湿度进行计算。下面这个例子从多个NetCDF中分别读取700hPa温度、风场U/V风量、相关湿度数据。 用meteolib包的saturation_mixing_ratio函数来计算饱和比湿,再通过计算比湿,再用divergence函数计算散度。 2 c/ r3 a* l8 F, b2 T) B' E% n
; V+ F# Q7 E4 |2 b6 g/ u
3 t: @7 }2 z6 |8 C[Python] 纯文本查看 复制代码 print 'Open data files...'
f_air = addfile('/Users//Documents/MeteoInfo/milab_sample/temp_adv//air.2011.nc')
f_uwnd = addfile('/Users//Documents/MeteoInfo/milab_sample/temp_adv//uwnd.2011.nc')
f_vwnd = addfile('/Users//Documents/MeteoInfo/milab_sample/temp_adv//vwnd.2011.nc')
f_rhum = addfile('/Users//Documents/MeteoInfo/milab_sample/temp_adv//rhum.2011.nc')
print 'Read data array...'
tidx = 173 # Jun 23, 2011
t = f_air.gettime(tidx)
lidx = 3 # 700 hPa
air = f_air['air'][tidx,lidx,:,:]
uwnd = f_uwnd['uwnd'][tidx,lidx,:,:]
vwnd = f_vwnd['vwnd'][tidx,lidx,:,:]
rhum = f_rhum['rhum'][tidx,lidx,:,:]
# Calculate
print 'Calculate...'
prs = 700
g = 9.8
qs = meteolib.saturation_mixing_ratio(prs, air-273.15)
q = qs * rhum / 100
qhdivg = meteolib.divergence(q*uwnd/g, q*vwnd/g) * 1000
#Plot
print 'Plot...'
axesm()
geoshow('coastline', color='k')
layer = contourf(qhdivg, 20)
title(u'Water Vapor Flux Divergency (气象可视化)')
colorbar(layer)
xlim(0, 360)
ylim(-90, 90)
savefig('/Users//Desktop/MeteoInfo/script/wxj/water_vapor_flux_divg.png', 600, 300, dpi=300) , ^% e% M/ ?% r' l
3 Q) W+ q( q* M! f o& a+ @
欢迎大家加入meteoinfo官方qq交流群,一起交流。王老师和众多气科院的老师们都在的。只要申请,我都会通过的。
7 Q+ {1 H' [0 n0 a' p* y5 {/ I |