MeteoInfoLab脚本示例：垂直螺旋度

半颗麦芽糖 · 发表于 2017-4-24 14:50:48

MeteoInfo 发表于 2017-4-24 14:42
很有可能是版本问题

恩恩，整个meteoinfo安装包我都是重新传上去的，脚本也是按照您上午给的改了路径的，报错的line8也就是这句：u = f['u-component_of_wind_isobaric'][tidx,level,lat,lon][::-1,:,:]，所以考虑是不是milab.sh的版本问题。

MeteoInfo · 发表于 2017-4-24 14:52:04

半颗麦芽糖发表于 2017-4-24 14:50
恩恩，整个meteoinfo安装包我都是重新传上去的，脚本也是按照您上午给的改了路径的，报错的line8也就是这 ...

跟那个没关系

MeteoInfo · 发表于 2017-4-24 14:59:26

半颗麦芽糖发表于 2017-4-24 14:50
恩恩，整个meteoinfo安装包我都是重新传上去的，脚本也是按照您上午给的改了路径的，报错的line8也就是这 ...

刚在大型机上试了一下，应该没问题。或者你下载最新的1.4.2R6版本试试（http://www.meteothinker.com/downloads/index.html）

有事秦天 · 发表于 2017-4-24 16:38:56

MeteoInfo 发表于 2017-4-24 10:01
参考下面的脚本：

感谢王老师，已经成功出图。

半颗麦芽糖 · 发表于 2017-4-24 17:16:34

本帖最后由半颗麦芽糖于 2017-4-24 17:18 编辑

MeteoInfo 发表于 2017-4-24 14:59
刚在大型机上试了一下，应该没问题。或者你下载最新的1.4.2R6版本试试（http://www.meteothinker.com/dow ...

王老师，下载了最新版本的MeteoInfo成功出图，现在考虑对WRF跑出来的数据进行同样的作图，考虑到投影问题，加入了转化投影的语句，但是最后报错了，请教王老师。

脚本如下：
fn = 'E:/gongzuo/3.8.1/test38101.ctl'
f = addfile(fn)
print 'Calculate vertical helicity...'
tidx = 12
level = '500:200'
lat = '29:39'
lon = '75:105'
u = f['U'][tidx,level,lat,lon][:,:,:]
v = f['V'][tidx,level,lat,lon][:,:,:]
toproj = projinfo()
u = u.project(toproj=toproj)
u.dims[0].setDimType(DimensionType.Y)
u.dims[1].setDimType(DimensionType.X)
v = v.project(toproj=toproj)
vort = hcurl(u, v)
wd1 = vort[:,:,'92']
lev1 = wd1.dimvalue(0)
lev1 = lev1[::-1]
lev2 = meteo.p2h(lev1)
levels = []
for j in range(0, len(lev1)):
levels.append('%i' % lev1[j])
wd1.setdimvalue(0, lev2)
print 'Plot...'
layer = contourf(wd1*100000,cmap='BlueRed')
#title(u'Vertical helicity (Longitude: 90°E)')
colorbar(layer,fontsize=15)
layer = contour(wd1*100000,colors='k')
clabel(layer)
xticks(arange(29,40,2),fontsize=15)
yticks(lev2, levels,fontsize=15)
xlabel('Latitude')
ylabel('Pressure (hPa)')
savefig('F:/tuer/fnl20140816/test1.png', width=600, height=600)

报错如下：
Traceback (most recent call last):
File "<iostream>", line 12, in <module>
NameError: name 'DimensionType' is not defined

MeteoInfo · 发表于 2017-4-24 19:23:27

半颗麦芽糖发表于 2017-4-24 17:16
王老师，下载了最新版本的MeteoInfo成功出图，现在考虑对WRF跑出来的数据进行同样的作图，考虑到投影问 ...

u = f['U'][tidx,level,lat,lon][:,:,:]
v = f['V'][tidx,level,lat,lon][:,:,:]
后面的[:,:,:]是多余的。

u.dims[0].setDimType(DimensionType.Y)
u.dims[1].setDimType(DimensionType.X)
这两句应该不需要。

半颗麦芽糖 · 发表于 2017-4-24 20:11:54

MeteoInfo 发表于 2017-4-24 19:23
u = f['U'][tidx,level,lat,lon][:,:,:]
v = f['V'][tidx,level,lat,lon][:,:,:]
后面的[:,:,:]是多余 ...

王老师，照您说的进行了修改，第十四行wd1 = vort[:,:,'92']报了错：
indices must be 2 dimensions!
Traceback (most recent call last):
  File "<iostream>", line 14, in <module>
  File "F:\11\MeteoInfo\pylib\mipylib\numeric\dimarray.py", line 73, in __getitem__
raise IndexError()
IndexError
将其改成二维的，另一个IndexError仍然存在。请王老师指点。
ctl前面部分如下：
dset E:\gongzuo\3.8.1\test38101.dat
options  byteswapped
undef 1.e30
title  OUTPUT FROM WRF V3.8.1 MODEL
pdef  195 132 lcc  35.000 94.000 98.000 66.500  60.00000  30.00000 94.00000  30000.000  30000.000
xdef  625 linear 51.75831 0.13513513
ydef  310 linear 12.61555 0.13513513
zdef 19 levels
1000.00000
950.00000
900.00000
850.00000
800.00000
750.00000
700.00000
650.00000
600.00000
550.00000
500.00000
450.00000
400.00000
350.00000
300.00000
250.00000
200.00000
150.00000
100.00000
tdef 41 linear 00Z25OCT2014    180MN

MeteoInfo · 发表于 2017-4-25 13:30:05

半颗麦芽糖发表于 2017-4-24 20:11
王老师，照您说的进行了修改，第十四行wd1 = vort[:,:,'92']报了错：
indices must be 2 dimensions!
T ...

参考这个脚本，MeteoInfo做了更新，需要下载最新版本：

fn = 'D:/Temp/grads/wrfout_d02.ctl'
f = addfile(fn)
print 'Calculate vertical vort...'
tidx = 0
level = '500:200'
lat = '39:49'
lon = '120:130'
u = f['UMET'][tidx,level,lat,lon]
v = f['VMET'][tidx,level,lat,lon]
toproj = projinfo()
u = u.project(toproj=toproj)
v = v.project(toproj=toproj)
vort = hcurl(u, v)
wd1 = vort[:,:,'125']
lev1 = wd1.dimvalue(0)
lev2 = meteo.p2h(lev1)
levels = []
for j in range(0, len(lev1)):
levels.append('%i' % lev1[j])
wd1.setdimvalue(0, lev2)
print 'Plot...'
layer = contourf(wd1*100000,cmap='BlueRed')
title(u'Vertical vort (Longitude: 125°E)')
colorbar(layer,fontsize=15)
layer = contour(wd1*100000,colors='k')
clabel(layer)
#xticks(arange(29,40,2),fontsize=15)
yticks(lev2, levels,fontsize=15)
xlabel('Latitude')
ylabel('Pressure (hPa)')

半颗麦芽糖 · 发表于 2017-4-25 18:59:08

MeteoInfo 发表于 2017-4-25 13:30
参考这个脚本，MeteoInfo做了更新，需要下载最新版本：

王老师，您好，我在超算上使用milab.sh参照您上面的脚本画图报了错，报错如下，还弹出了小窗口：

其中line25是layer = contour(wd1*100000,colors='k')

当我把等值线的指令注释掉以后，有报了错，显示的是保存图片错误，但是我用相同的语句画fnl资料的图却是对的，请老师指点，谢谢！

MeteoInfo · 发表于 2017-4-25 21:02:27

半颗麦芽糖发表于 2017-4-25 18:59
王老师，您好，我在超算上使用milab.sh参照您上面的脚本画图报了错，报错如下，还弹出了小窗口：

其中 ...

1、需要下载1.4.2R7版本
2、我是根据自己的数据写的经纬度范围，你需要根据自己数据来确定。

		自动登录	找回密码
密码			立即注册