利用NCEP数据做水汽通量散度图

MeteoInfo

这个例子展示了利用NCEP数据做水汽通量散度图的MeteoInfo脚本程序。需要更新MeteoInfo最新的文件（可以在置顶帖子中找到），这个没有经过验证，有兴趣的朋友可以用GrADS对同样的数据进行做图，对比一下。从图形来看2011年6月23日北京地区是水汽汇集中心，那天北京是大暴雨。

#---- Import class libraries
import clr
clr.AddReferenceByPartialName("System.Windows.Forms")
clr.AddReferenceByPartialName("System.Drawing")
from System.Windows.Forms import *
from System.Drawing import *
clr.AddReference("MeteoInfoC.dll")
from MeteoInfoC import *
from MeteoInfoC.Data import *
from MeteoInfoC.Data.MeteoData import *

#---- Define variables
BaseDir = "C:\\Program Files (x86)\\MeteoInfo\\"
MapDir = BaseDir + "Map\\"
DataDir = "E:\\Temp\\nc\\"

#---- Create MIApp object
myApp = MIApp()
#---- Open layers
myApp.OpenLayer(MapDir + "country1.shp")
myApp.SetLegendBreak("country1.shp",0,Color.Yellow,Color.Black,1,True,False,True)
#---- Open data
DataAir = MeteoDataInfo()
DataUwnd = MeteoDataInfo()
DataVwnd = MeteoDataInfo()
DataRhum = MeteoDataInfo()
DataAir.OpenNCData(DataDir + "air.2011.nc")
DataUwnd.OpenNCData(DataDir + "uwnd.2011.nc")
DataVwnd.OpenNCData(DataDir + "vwnd.2011.nc")
DataRhum.OpenNCData(DataDir + "rhum.2011.nc")
#---- Calculation
#---- Set date to June 23
tIdx = 171
DataAir.TimeIndex = tIdx;
DataUwnd.TimeIndex = tIdx;
DataVwnd.TimeIndex = tIdx;
DataRhum.TimeIndex = tIdx;
#---- Set level to 700hPa
lIdx = 3
DataAir.LevelIndex = lIdx;
DataUwnd.LevelIndex = lIdx;
DataVwnd.LevelIndex = lIdx;
DataRhum.LevelIndex = lIdx;
#---- Get grid data
air = DataAir.GetGridData("air")
uwnd = DataUwnd.GetGridData("uwnd")
vwnd = DataVwnd.GetGridData("vwnd")
rhum = DataRhum.GetGridData("rhum")
#---- Calculate
prs = 700
g = 9.8
es = 6.112*DataMath.Exp(17.67*(air-273.16)/(air-29.65))
qs = 0.62197*es/(prs-0.378*es)
q = qs*rhum/100
qhdivg = DataMath.Hdivg(q*uwnd/g,q*vwnd/g)
qv = rhum*es/100
uv = DataMath.Magnitude(uwnd, vwnd)
uvq = uv*qv/(9.8*1000)

#---- Plot data
myApp.SetDrawType("shaded")
myApp.Display(qhdivg)
myApp.MoveLayerToTop("country1.shp")
myApp.ZoomLonLatEx(0,360,-90.1,90.1)
aTime = DataAir.GetTime()
myApp.SetTitle("Water Vapor Flux Divergence - " + aTime.ToString("yyyy-MM-dd"))
myApp.MapLayout.PaintGraphics()
Application.Run(myApp)

mofangbao

顶起楼主，问个有点傻的问题啊，能不能用其他语言来写脚本？

MeteoInfo

也可以用IronRuby，重要的是脚本语言要能支持.Net。

MeteoInfo

用GrADS验证了一下，结果是一样的。



gs文件内容如下。GrADS中t, lev的设置是从1开始的，MeteoInfo脚本是从0开始，因此gs中t要大1。
'reinit'
'sdfopen e:/temp/nc/uwnd.2011.nc'
'sdfopen e:/temp/nc/vwnd.2011.nc'
'sdfopen e:/temp/nc/rhum.2011.nc'
'sdfopen e:/temp/nc/air.2011.nc'
t=172
'set lev 700'
'set t 't
'set gxout shaded'

'define prs=lev'
'define g=9.8'
'define es=(6.112*exp((17.67*(air.4-273.16))/(air.4-29.65)))'
'define qs=(0.62197*es/(prs-0.378*es))'
'define q=qs*rhum.3/100'
'define qhdivg=hdivg(q*uwnd.1/g,q*vwnd.2/g)'
'd qhdivg'

mofangbao

哦  明白了 哈哈  看来我用的话还是用c#吧

MeteoInfo

C#并非脚本语言，MeteoInfo目前的脚本语言用的是IronPython，MeteoInfo里面已经包含了IronPython的库，无须再单独安装。

mofangbao

我是打算在c#里用类库的，放到自己的软件里去，呵呵

MeteoInfo

好哇，有什么成果以后给大家秀秀

mofangbao

MeteoInfo 发表于 2011-7-21 10:17

                               
登录/注册后可看大图

好哇，有什么成果以后给大家秀秀

呵呵 那是当然的 还得多多向你请教呢

zjusutao

这段gs里公式的系数有错误，经验公式的那个17.2693822应该是17.67