T639细网格湿螺旋度绘图gs

njlo235

学了一个月的湿螺旋度绘制，终于勉强完成，总感觉有点啰嗦，请高手指点一二
'reinit'
'open D:\grads\T639\T639.ctl'
'set parea 1 9.6 0 8.3'
'set csmooth on'
start_time=2014090220
rct.1='000'
rct.2='003'
rct.3='006'
rct.4='009'
rct.5='012'
rct.6='015'
rct.7='018'
rct.8='021'
rct.9='024'
rct.10='027'
rct.11='030'
rct.12='033'
rct.13='036'
rct.14='039'
rct.15='042'
rct.16='045'
rct.17='048'
rct.18='051'
rct.19='054'
rct.20='057'
rct.21='060'
rct.22='063'
rct.23='066'
rct.24='069'
rct.25='072'
rct.26='075'
rct.27='078'
rct.28='081'
rct.29='084'
rct.30='087'
rct.31='090'
rct.32='093'
rct.33='096'
rct.34='099'
rct.35='102'
rct.36='105'
rct.37='108'
rct.38='111'
rct.39='114'
rct.40='117'
rct.41='120'

rtt.1=1
rtt.2=2
rtt.3=3
rtt.4=4
rtt.5=5
rtt.6=6
rtt.7=7
rtt.8=8
rtt.9=9
rtt.10=10
rtt.11=11
rtt.12=12
rtt.13=13
rtt.14=14
rtt.15=15
rtt.16=16
rtt.17=17
rtt.18=18
rtt.19=19
rtt.20=20
rtt.21=21
rtt.22=22
rtt.23=23
rtt.24=24
rtt.25=25
rtt.26=26
rtt.27=27
rtt.28=28
rtt.29=29
rtt.30=30
rtt.31=31
rtt.32=32
rtt.33=33
rtt.34=34
rtt.35=35
rtt.36=36
rtt.37=37
rtt.38=38
rtt.39=39
rtt.40=40
rtt.41=41

it=1
while(it<=41)
'set t 'rtt.it
'set lat 26.125 36.25'
'set lon 97 110.5'
'set lev 500'
'set xlopts 1 4 0.12'
'set ylopts 1 4 0.12'
'SET XLevs  97 100  103 106  109 '
'SET YLEVS 27 29 31 33 35 37'
'set grads off'
'set grid  off'
'set mpdset cnworld'
'set map 15 1 1'
'draw map'
'set mpdset sichuan sichuan_q '
'draw map'

'set gxout contour'
'set csmooth on'
'set cthick 10'
'set cint 10'
'set clab on'
'set clevs -40 -60  -120   -180  -240 -300 -360 -420  10 20 30 40 60 80  120 160  200 240  320 400 480 580 620  '
'set clab forced'

'define prs=500'
'define es=(6.112*exp((17.27*(tmpprs-273.16))/(tmpprs-35.86)))'
'define qs=0.62197*es/(prs-0.378*es)'
'define q=1000*rhprs*qs/100'
'define fuq=1/9.80665*ugrdprs*q'
'define fvq=1/9.80665*vgrdprs*q'
'define divf=hdivg(fuq,fvq)'
'define hlx=VVELprs*RELVprs*divf'  
'd hlx*100000000000'

'set gxout barb'
'set lev 500'
'set cthick 6'
'set digsiz 0.05'
'set ccolor 1'
'd 2.5*skip(UGRDprs,3,3);2.5*VGRDprs'
'q time'
rst=substr(result,8,12)
'set strsiz 0.16 0.16'
'draw string 0.7 8.2 T639028 500hPa RHLX[unit:10e-11pa*s-3]&wind on 'start_time'_'rct.it''
'printim  d:\grads\micaps\T639_'start_time'\rhlxgif\500\'start_time'_'rct.it'.gif     white'

'clear'
it=it+1
endwhile

river

挺好的，感谢分享

天边的彩虹

感谢分享。

sunflower

感谢楼主分享，就是前面赋了那么多值，有木有更简洁些的办法·

吃柠檬的猫

可以问一下楼主这个gs计算出来的湿螺旋度单位是什么吗？很多文献里单位是pa/s^3,对这个单位有点疑惑啊，楼主算出来的单位也是这个吗？

xgm324

看不太懂，感谢分享。