MeteoInfoLab脚本示例：垂直螺旋度

半颗麦芽糖

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:50
恩恩，整个meteoinfo安装包我都是重新传上去的，脚本也是按照您上午给的改了路径的，报错的line8也就是这 ...

跟那个没关系

MeteoInfo

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

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

有事秦天

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

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

半颗麦芽糖

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 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)
这两句应该不需要。

半颗麦芽糖

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-24 20:11
王老师，照您说的进行了修改，第十四行wd1 = vort[:,:,'92']报了错：
indices must be 2 dimensions!
T ...

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

1. fn = 'D:/Temp/grads/wrfout_d02.ctl'
   
2. f = addfile(fn)
   
3. print 'Calculate vertical vort...'
   
4. tidx = 0
   
5. level = '500:200'
   
6. lat = '39:49'
   
7. lon = '120:130'
   
8. u = f['UMET'][tidx,level,lat,lon]
   
9. v = f['VMET'][tidx,level,lat,lon]
   
10. toproj = projinfo()
    
11. u = u.project(toproj=toproj)
    
12. v = v.project(toproj=toproj)
    
13. vort = hcurl(u, v)
    
14. wd1 = vort[:,:,'125']
    
15. lev1 = wd1.dimvalue(0)
    
16. lev2 = meteo.p2h(lev1)
    
17. levels = []
    
18. for j in range(0, len(lev1)):
    
19.     levels.append('%i' % lev1[j])
    
20. wd1.setdimvalue(0, lev2)
    
21. print 'Plot...'
    
22. layer = contourf(wd1*100000,cmap='BlueRed')
    
23. title(u'Vertical vort (Longitude: 125°E)')
    
24. colorbar(layer,fontsize=15)
    
25. layer = contour(wd1*100000,colors='k')
    
26. clabel(layer)
    
27. #xticks(arange(29,40,2),fontsize=15)
    
28. yticks(lev2, levels,fontsize=15)
    
29. xlabel('Latitude')
    
30. ylabel('Pressure (hPa)')

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

半颗麦芽糖

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

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

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

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

MeteoInfo

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

其中 ...

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