GRADS对NCEP grib2资料求风暴相对水平螺旋度

戎马秀才

用于求NCEP FNL GRIB2资料的水平螺旋度，希望能方便大家。
'reinit'
'open D:\2015\2015_zw\grib21\fnl_20140423_06_00.ctl'
'enable print D:\2015\2015_zw\gmf\Hrs12.gmf'
t1=2;
*while(t1<=2)
'reset'
'set mpdset cnhimap'
'set map 1 1 7'
'set t 't1
'set lon 80 110'
'set lat 30 50'
***********************************************************
*hh=6
*while(hh<=11)
*'set lev 850 600'
***********************************************************
***********************************************************
*'define usum=usum+u'
*'define vsum=vsum+u'
*hh=hh+1
*endwhile
*'define uc=usum/6.'
*'define vc=vsum/6.'
*'undefine usum'
*'undefine vsum'
***********************************************************
*计算850hpa到600hpa的平均风
'uc=mean(UGRDprs,lev=850,lev=600)'
'vc=mean(VGRDprs,lev=850,lev=600)'
'speed=mag(uc,vc)'
'direc=atan2(vc,uc)'
'cdirec=direc-40./180.*3.14159'
if(cdirec<-3.14159)
cdirec=cdirec+3.14159*2
endif
'cuc=speed*cos(cdirec)'
'cvc=speed*sin(cdirec)'
***********************************************************
*此处验证一下速度矢量是否有错误
*'d cuc;cvc;mag(cuc,cvc)'
***********************************************************
'Hrs=0.'
'zHrs=0.'
hh=6;
while(hh<=10)
*'zHrs=cuc*v(z='hh+1')'
*上式证明下式中第一项与第四项相当，其余项为零
*'zHrs=-u(z='hh')*v(z='hh+1')+cuc*v(z='hh+1')-cuc*v(z='hh')+v(z='hh')*u(z='hh+1')-cvc*u(z='hh+1')+cvc*u(z='hh')'
'zHrs=(UGRDprs(z='hh+1')-cuc)*(VGRDprs(z='hh')-cvc)-(UGRDprs(z='hh')-cuc)*(VGRDprs(z='hh+1')-cvc)'
say hh+1
'Hrs=Hrs+zHrs'
hh=hh+1
endwhile
***********************************************************
*'define Hrs=UGRDprs(lev=850)-UGRDprs(lev=600)'
'set xlopts 1 6 0.17'
'set ylopts 1 6 0.17'
'set grid off'
'set timelab off'
'set grads off'
'set gxout shaded'
*'set cint 3'
'set csmooth on'
'd Hrs'
'set gxout contour'
*'set cint 3'
*'set cstyle 1'
'set csmooth on'
'set clopts 0 -1 0.08'
'd Hrs'
'print'
* t1=t1+1;
*endwhile
'disable print'
'c'

river

配个图就更好了

yf431

非常感谢，配两张图就更好了

青衣

多谢楼主的分享。

hanxinlei12

怎么打开的，处理过程能给以下吗啊啊

戎马秀才

hanxinlei12 发表于 2015-8-6 13:27
怎么打开的，处理过程能给以下吗啊啊

你好！
1.你说的“打开的”是指打开什么文件？
2.处理过程就是这个gs文件，你意思是指？

ENSO522

好东西

青衣

感谢楼主，很不错

sunflower

楼主给力，真是帮了我的大忙呀，非常感谢

momoxiangge

谢谢楼主，一看就是个编程高手