MeteoInfoLab脚本示例：水汽通量散度计算

qxsw2016

瑛雩 发表于 2019-6-8 17:51
老师您好，我有一个问题想请教一下。
我在计算水汽通量散度时，前面的计算都有结果，但是在调用函数hdivg ...

你把代码贴出来

woaichimian666

http://bbs.06climate.com/forum.p ... p;extra=&page=1
想问这个hdivg里面的函数，也是用m计算的吗？

woaichimian666

就是函数hdivg的地球半径是m还是cm为单位？

920daibingbing

赞赞！已经成功出图

cc123

王老师，我计算出来的结果与预期的量级上有出入，参考其他帖子，简单计算了下单位的量级，觉得将计算公式改成如下的似乎更合理一点：
g = 9.8
es = 6.112*exp(17.67*(air-273.16)/(air-29.65))
qs = 621.97*es/(prs-0.378*es)                    # 克/千克
q = qs*rhum/100
qhdivg = hdivg(q*uwnd/g,q*vwnd/g) / 100        #考虑计算散度时，地球半径用的单位是m，变为cm，所以要除以100

上述公式计算出来的单位是，g*cm^-2*s^-1*hPa^-1

参考的帖子：
http://bbs.06climate.com/forum.p ... 24090&extra=&page=1
https://bbs.06climate.com/forum.php?mod=viewthread&tid=95245

Linda.W

王老师，望指教在计算水汽通量散度的基础上如何进一步用Meteoinfo画出某点水汽通量散度的时间-高度剖面图

MeteoInfo

Linda.W 发表于 2024-12-16 11:48
王老师，望指教在计算水汽通量散度的基础上如何进一步用Meteoinfo画出某点水汽通量散度的时间-高度剖面图{: ...

参考这个：

1. print('Open data files...')
   
2. f_air = addfile('D:/Temp/nc/air.2011.nc')
   
3. f_uwnd = addfile('D:/Temp/nc/uwnd.2011.nc')
   
4. f_vwnd = addfile('D:/Temp/nc/vwnd.2011.nc')
   
5. f_rhum = addfile('D:/Temp/nc/rhum.2011.nc')
   
6. 
   
7. print('Read data array...')
   
8. stidx = 0
   
9. etidx = 20
   
10. levels = f_rhum['level'][:]    # rhum has less levels
    
11. nlev = len(levels)
    
12. lonlim = '60:140'
    
13. latlim = '10:60'
    
14. air = f_air['air'][stidx:etidx,0:nlev,latlim,lonlim]
    
15. uwnd = f_uwnd['uwnd'][stidx:etidx,0:nlev,latlim,lonlim]
    
16. vwnd = f_vwnd['vwnd'][stidx:etidx,0:nlev,latlim,lonlim]
    
17. rhum = f_rhum['rhum'][stidx:etidx,0:nlev,latlim,lonlim]
    
18. 
    
19. # Calculate
    
20. print('Calculate...')
    
21. prs = 700
    
22. g = 9.8
    
23. es = 6.112*exp(17.67*(air-273.16)/(air-29.65))
    
24. qs = 0.62197*es/(prs-0.378*es)
    
25. q = qs*rhum/100
    
26. qhdivg = meteolib.divergence(q*uwnd/g, q*vwnd/g)
    
27. 
    
28. # Get data at a location
    
29. lon = 120
    
30. lat = 40
    
31. data = qhdivg[:,:,str(lat),str(lon)]
    
32. data = np.rot90(data)
    
33. tt = data.dimvalue(1)
    
34. height = meteolib.pressure_to_height_std(levels)
    
35. 
    
36. # Plot
    
37. print('Plot...')
    
38. contourf(tt, height, data, 20, extend='both')
    
39. colorbar()
    
40. xaxis(axistype='time')
    
41. yticks(height, levels)
    
42. title('Water Vapor Flux Divergency (Lon: {}, Lat: {})'.format(lon, lat))
    

print('Open data files...')<br /> f_air = addfile('D:/Temp/nc/air.2011.nc')<br /> f_uwnd = addfile('D:/Temp/nc/uwnd.2011.nc')<br /> f_vwnd = addfile('D:/Temp/nc/vwnd.2011.nc')<br /> f_rhum = addfile('D:/Temp/nc/rhum.2011.nc')<br /> <br /> print('Read data array...')<br /> stidx = 0<br /> etidx = 20<br /> levels = f_rhum['level'][:]    # rhum has less levels<br /> nlev = len(levels)<br /> lonlim = '60:140'<br /> latlim = '10:60'<br /> air = f_air['air'][stidx:etidx,0:nlev,latlim,lonlim]<br /> uwnd = f_uwnd['uwnd'][stidx:etidx,0:nlev,latlim,lonlim]<br /> vwnd = f_vwnd['vwnd'][stidx:etidx,0:nlev,latlim,lonlim]<br /> rhum = f_rhum['rhum'][stidx:etidx,0:nlev,latlim,lonlim]<br /> <br /> # Calculate<br /> print('Calculate...')<br /> prs = 700<br /> g = 9.8<br /> es = 6.112*exp(17.67*(air-273.16)/(air-29.65))<br /> qs = 0.62197*es/(prs-0.378*es)<br /> q = qs*rhum/100<br /> qhdivg = meteolib.divergence(q*uwnd/g, q*vwnd/g)<br /> <br /> # Get data at a location<br /> lon = 120<br /> lat = 40<br /> data = qhdivg[:,:,str(lat),str(lon)]<br /> data = np.rot90(data)<br /> tt = data.dimvalue(1)<br /> height = meteolib.pressure_to_height_std(levels)<br /> <br /> # Plot<br /> print('Plot...')<br /> contourf(tt, height, data, 20, extend='both')<br /> colorbar()<br /> xaxis(axistype='time')<br /> yticks(height, levels)<br /> title('Water Vapor Flux Divergency (Lon: {}, Lat: {})'.format(lon, lat))<br />