湿位涡与雷达反射率叠加

nikolan

function b( args  )
dat = subwrd( args,1)
'open f:\ziliao\'dat'.ctl'
'set lon 90 140 1'
'set lat 10 60 1'
'set lev 1000 100'
'define prs=lev'
'set ccolor 1'
'set grid off'
'set cthick 1'
'define  l1=TMPprs'
if(l1>263.0);
'define a1=0.622*6.11*EXP(17.26*(TMPprs-273.16)/(TMPprs-35.86))'
'define b1=lev-0.278*EXP(17.26*(TMPprs-273.16)/(TMPprs-35.86))'
'define  qs1=a1/b1'
endif
if(l1<=263.0);
'define a1=0.622*6.11*EXP(21.87*(TMPprs-273.16)/(TMPprs-7.66))'
'define b1=lev-0.278*EXP(21.87*(TMPprs-273.16)/(TMPprs-7.66))'
'define  qs1=a1/b1'
endif
'define  q1=qs1*RHprs'
'define e1=lev*q1/100./(0.62197+q1/100.0)'
'define tk1=55.0+2840.0/(3.5*log(TMPprs)-log(e1)-4.805)'
'define pot10=TMPprs*pow((1000/lev),(0.2854*(1.0-0.28*q1/100.0)))'
'define ept1=pot10*exp(((3376./tk1)-2.54)*q1/100.0*(1.0+0.81*q1/100.0))'
'define  l1=TMPprs(z-1)'
if(l1>263.0);
'define a1=0.622*6.11*EXP(17.26*(TMPprs(z-1)-273.16)/(TMPprs(z-1)-35.86))'
'define b1=lev(z-1)-0.278*EXP(17.26*(TMPprs(z-1)-273.16)/(TMPprs(z-1)-35.86))'
'define  qs1=a1/b1'
endif
if(l1<=263.0);
'define a1=0.622*6.11*EXP(21.87*(TMPprs(z-1)-273.16)/(TMPprs(z-1)-7.66))'
'define b1=lev(z-1)-0.278*EXP(21.87*(TMPprs(z-1)-273.16)/(TMPprs(z-1)-7.66))'
'define  qs1=a1/b1'
endif
'define  q1=qs1*RHprs(z-1)'
'define e1=lev(z-1)*q1/100./(0.62197+q1/100.0)'
'define tk1=55.0+2840.0/(3.5*log(TMPprs(z-1))-log(e1)-4.805)'
'define pot11=TMPprs(z-1)*pow((1000/lev(z-1)),(0.2854*(1.0-0.28*q1/100.0)))'
'define ept2=pot11*exp(((3376./tk1)-2.54)*q1/100.0*(1.0+0.81*q1/100.0))'
'define  l1=TMPprs(z+1)'
if(l1>263.0);
'define a1=0.622*6.11*EXP(17.26*(TMPprs(z+1)-273.16)/(TMPprs(z+1)-35.86))'
'define b1=lev(z+1)-0.278*EXP(17.26*(TMPprs(z+1)-273.16)/(TMPprs(z+1)-35.86))'
'define  qs1=a1/b1'
endif
if(l1<=263.0);
'define a1=0.622*6.11*EXP(21.87*(TMPprs(z+1)-273.16)/(TMPprs(z+1)-7.66))'
'define b1=lev(z+1)-0.278*EXP(21.87*(TMPprs(z+1)-273.16)/(TMPprs(z+1)-7.66))'
'define  qs1=a1/b1'
endif
'define  q1=qs1*RHprs(z+1)'
'define e1=lev(z+1)*q1/100./(0.62197+q1/100.0)'
'define tk1=55.0+2840.0/(3.5*log(TMPprs(z+1))-log(e1)-4.805)'
'define pot12=TMPprs(z+1)*pow((1000/lev(z+1)),(0.2854*(1.0-0.28*q1/100.0)))'
'define ept3=pot12*exp(((3376./tk1)-2.54)*q1/100.0*(1.0+0.81*q1/100.0))'
'define coriol = 2*7.29e-5 * sin(lat*3.1415/180)'
'define dudy = cdiff(UGRDprs,y)/(111177*cdiff(lat,y))'
'define dvdx = cdiff(VGRDprs,x)/(111177*cdiff(lon,x)*cos(lat*3.1415/180))'
'define t=ept1'
'define dt = ept2-ept3'
'define dp = 100*(lev(z-1)-lev(z+1))'
'define dtdp = dt/dp'
'define du = UGRDprs(z-1)*pow(1000/lev(z-1),0.286)-UGRDprs(z+1)*pow(1000/lev(z+1),0.286)'
'define dv = VGRDprs(z-1)*pow(1000/lev(z-1),0.286)-VGRDprs(z+1)*pow(1000/lev(z+1),0.286)'
'define dudp=du/dp'
'define dvdp=dv/dp'
'define dtdy = cdiff(t,y)/(111177*cdiff(lat,y))'
'define dtdx = cdiff(t,x)/(111177*cdiff(lon,x)*cos(lat*3.1415/180))'
'define mpv2=9.8*(dvdp*dtdx-dudp*dtdy)'
'define mpv1= -9.8*(coriol+dvdx-dudy)*dtdp'
'define mpv=mpv1+mpv2'
'set lon 105 125 1'
'set lat 30 50 1'
'set lev 400'
'set mpdset cnworld'
'set grid off'
'set grads off'
'set lin 0 1 6'
'set yaxis 30  50 1'
'set xaxis 105 125 1'
'set xlopts 1 4 0.14'
'set ylopts 1 4 0.14'
'set cthick 1'
'set cint 0.05'
'd smth9(mpv2*1000000)'
'printim d:\fan\wz\mpv2\400\'dat'.gif gif white x800 y600'
'reinit'

良辰

请问雷达反射率的gs部分那？

ayzqs

求雷达反射率读写程序！

cuhum

求雷达反射率读写程序！

天边的那边

好长的程序{:5_213:}

彩虹依旧

求雷达反射率读写程序！

刘大壮壮

kzpyezi

多谢楼主分享!

rossbyman

这个程序不错，拿走好好研究下

qxtlyf

这个程序不错，拿走好好研究下