in=read(desc.dat) val=sublin(in,2) starttime=getgdtg(val) in=read(desc.dat) val=sublin(in,2) endtime=getgdtg(val) in=read(desc.dat) name=sublin(in,2) in=read(desc.dat) ctl=sublin(in,2) year=substr(starttime,9,4) "open tmp.ctl" "open "ctl trackfile="track.tmp" "set rgb 30 200 200 200" smooth=2 "set dfile 1" "set t 1" "q dims" rec=sublin(result,5) minit=subwrd(rec,6) ndate=getndate(minit) "set dfile 1" "q file" rec=sublin(result,5) tmax=subwrd(rec,12) "set t 1 "tmax "define var1=ave(vtu,t-"smooth",t+"smooth")" "define var2=ave(vtl,t-"smooth",t+"smooth")" "define var3=ave(size,t-"smooth",t+"smooth")" "set t 1" "set display color white" "set parea 1.5 9.5 1.25 7.5" "clear" xmin=-600 xmax=300 ymin=-600 ymax=300 "set mproj scaled" "set x "xmin " " xmax "set y "ymin " " ymax "set xaxis "xmin " " xmax " 100 " "set yaxis "ymin " " ymax " 100 " "set xlopts 1 3 0.15" "set ylopts 1 3 0.15" "set clevs 5000" "set mpdraw off" "set grads off" "set grid off" "d lat" "set strsiz 0.15" "set string 2 c 5 0" "draw string 4.0 0.85 Cold Core" "draw string 8.5 0.85 Warm Core" "set strsiz 0.2" "set string 4 c 10 0" "draw string 5.5 0.5 -V`bT`n`aL`n [900hPa-600hPa Thermal Wind]" "set string 4 c 10 90" "draw string 0.25 4.25 -V`bT`n`aU`n [600hPa-300hPa Thermal Wind]" "set strsiz 0.15" "set string 2 c 5 90" "draw string 0.75 3.0 Cold Core" "draw string 0.75 6.5 Warm Core" "set string 1 c 5 0" "q gr2xy "xmin " -10 " xloc1=subwrd(result,3) yloc1=subwrd(result,6) "q gr2xy "xmax " 10 " xloc2=subwrd(result,3) yloc2=subwrd(result,6) "set line 30 1 10" "draw recf "xloc1 " " yloc1" " xloc2 " " yloc2 "q gr2xy -10 " ymin xloc3=subwrd(result,3) yloc3=subwrd(result,6) "q gr2xy 10 " ymax xloc4=subwrd(result,3) yloc4=subwrd(result,6) "set line 30 1 10" "draw recf "xloc3 " " yloc3 " " xloc4 " " yloc4 "set line 1 1 10" *"draw line "xloc1 " " yloc1 " " xloc2 " " yloc2 "set line 1 1 10" *"draw line "xloc3 " " yloc3 " " xloc4 " " yloc4 "set strsiz 0.1" "set string 1 c 5 0" "q gr2xy -300 125" xloc=subwrd(result,3) yloc=subwrd(result,6) *"draw string "xloc " " yloc " SHALLOW COLD CORE" *"draw string "xloc " " yloc-0.2 " (e.g. EXTRATROPICAL CYCLONE)" "q gr2xy 150 150" xloc=subwrd(result,3) yloc=subwrd(result,6) "draw string "xloc " " yloc " DEEP WARM CORE" *"draw string "xloc " " yloc-0.2 " (e.g. STRONG TC)" "q gr2xy 150 20" xloc=subwrd(result,3) yloc=subwrd(result,6) "draw string "xloc " " yloc " MODERATE WARM CORE" *"draw string "xloc " " yloc-0.2 " (e.g. WEAK TC," *"draw string "xloc " " yloc-0.4 " HYBRID CYCLONE," *"draw string "xloc " " yloc-0.6 " OCEANIC MIDLAT CYCLONE," *"draw string "xloc " " yloc-0.8 " WARM SECLUSION)" "q gr2xy -300 -300" xloc=subwrd(result,3) yloc=subwrd(result,6) "draw string "xloc " " yloc " DEEP COLD CORE" *"draw string "xloc " " yloc-0.2 " (e.g. EXTRATROPICAL CYCLONE)" "q gr2xy 150 -100" xloc=subwrd(result,3) yloc=subwrd(result,6) "draw string "xloc " " yloc " SHALLOW WARM-CORE" tt=1 "set gxout stat" minlat=9999 maxlat=-9999 minlon=9999 maxlon=-9999 while (tt <= tmax) "set t "tt "q dims" rec=sublin(result,5) time.tt=subwrd(rec,6) time2.tt=getndate(time.tt) "d "var1 rec=sublin(result,8) thick.tt=subwrd(rec,4) "d "var2 rec=sublin(result,8) slope.tt=subwrd(rec,4) "d "var3 rec=sublin(result,8) size.tt=subwrd(rec,4) data=getstorm(time.tt,trackfile) ilat.tt=subwrd(data,1) ilon.tt=subwrd(data,2) if (ilat.tt < minlat & ilat.tt > -900) minlat=ilat.tt endif if (ilat.tt > maxlat & ilat.tt > -900) maxlat=ilat.tt endif if (ilon.tt < minlon & ilon.tt > -900) minlon=ilon.tt endif if (ilon.tt > maxlon & ilon.tt > -900) maxlon=ilon.tt endif tt=tt+1 endwhile "set string 3 r 10 0" "set strsiz 0.175" timea=getndate(time.1) timeb=getndate(minit) daya=getday(time.1) dayb=getday(minit) dayc=getday(time.tmax) "draw string 11 8.25 "name"("year") [0.5`a0`n GFS Analysis]" "set string 1 r 5 0" "set strsiz 0.15" "q gr2xy 300 375" xloc=subwrd(result,3) yloc=subwrd(result,6) "draw string "xloc " " yloc-0.2 " Start (A): "time.1 " ("daya")" anltime=time.1 "draw string "xloc " " yloc-0.4 " End (Z): "time.tmax" ("dayc")" tt=1 xlocold=-999999 ylocold=-999999 while (tt <= tmax) if (ilon.tt > -100) data=getstorm(time.tt,trackfile) pres=subwrd(data,3) col=getcol(pres) hour=substr(time.tt,1,2) day=substr(time.tt,4,2) "q gr2xy "slope.tt " " thick.tt xloc=subwrd(result,3) yloc=subwrd(result,6) if (time2.tt > ndate) ltype=1 else ltype=3 endif if (tt = 1) xloca=xloc yloca=yloc endif if (tt = tmax) xlocz=xloc ylocz=yloc endif if (time.tt = minit) xlocc=xloc ylocc=yloc endif "set line "col " 1 10" if (tt > 1 & math_abs(xloc) < 15 & math_abs(xlocold) <15 & math_abs(yloc) < 15 & math_abs(ylocold) < 15) "set line "colold " "ltype" 10" "draw line "xloc " " yloc " " xlocold " " ylocold endif xlocold=xloc ylocold=yloc colold=col endif tt=tt+1 endwhile tt=1 xlocold=-999999 ylocold=-999999 while (tt <= tmax) if (ilon.tt > -100) data=getstorm(time.tt,trackfile) pres=subwrd(data,3) col=getcol(pres) hour=substr(time.tt,1,2) day=substr(time.tt,4,2) "q gr2xy "slope.tt " " thick.tt xloc=subwrd(result,3) yloc=subwrd(result,6) if (time2.tt > ndate) ltype=1 else ltype=3 endif if (tt = 1) xloca=xloc yloca=yloc endif if (tt = tmax) xlocz=xloc ylocz=yloc endif if (time.tt = minit) xlocc=xloc ylocc=yloc endif "set line "col " 1 10" msize=getsize(size.tt) "draw mark 3 " xloc " " yloc " "msize if (hour = 00) "set string 1 c 8 0" "set strsiz 0.15" "draw string "xloc " " yloc+0.1 " "day endif xlocold=xloc ylocold=yloc colold=col endif tt=tt+1 endwhile "set string 1 c 10 0" "set strsiz 0.35" "draw string "xloca " " yloca " A" "draw string "xlocz " " ylocz " Z" * Produce cheap legend for pressure labeling "set strsiz 0.11" "set string 1 c 3 0" "set line 1 1 10" "draw string 10.25 6.50 Intensity (hPa):" "draw line 9.60 6.40 10.90 6.40" pps = "1015 1010 1005 1000 995 990 985" dd=1 while (subwrd(pps,dd) != "") pp=subwrd(pps,dd) col=getcol(pp) xloc=10.00 yloc=6.25-(dd-1)*0.25 "set string "col " c 3 0" "draw string "xloc " " yloc " "pp dd=dd+1 endwhile pps="980 975 970 965 960 955 950" dd=1 while (subwrd(pps,dd) != "") pp=subwrd(pps,dd) col=getcol(pp) xloc=10.50 yloc=6.25-(dd-1)*0.25 "set string "col " c 3 0" "draw string "xloc " " yloc " "pp dd=dd+1 endwhile "set strsiz 0.11" "set string 1 c 3 0" "draw string 10.25 4.50 Mean radius of" "draw string 10.25 4.25 925hPa gale" "draw string 10.25 4.00 force wind (km):" "draw line 9.535 3.90 11.00 3.90" rrs="100 200 300 500 750" dd=1 while (subwrd(rrs,dd) != "") rr=subwrd(rrs,dd) msize=getsize(rr) xloc=10.00 yloc=3.75-(dd-1)*0.35 "set line 1 1 10" "draw mark 3 "xloc " " yloc " " msize "set string "col " c 3 0" if (rr > 100) "draw string "xloc+0.75 " " yloc " "rr else "draw string "xloc+0.75 " " yloc " <100" endif dd=dd+1 endwhile * Produce inset frame of storm track and Reynolds weekly SST field meanlat=(maxlat+minlat)/2 meanlon=(maxlon+minlon)/2 dlat=maxlat-minlat+10 dlon=maxlon-minlon+10 if (dlat > dlon) dist=dlat else dist=dlon endif lat1=meanlat-dist/2 lat2=meanlat+dist/2 lon1=meanlon-1.5*(dist/2) lon2=meanlon+1.5*(dist/2) yinc=10 if (dlon > 90) xinc=20 else xinc=10 endif "set vpage 1.0 5.5 5.50 8.50" "set parea 1.2 9 0 6.25" if (dist < 30) "set mpdset hires" endif "set gxout contour" "set grid off" "set grads off" "set clevs 500" "set mpdraw on" "set mproj latlon" "set map 3 1 10" "set xlab off" "set ylab off" "set xlint 5" "set ylint 5" "set lat "lat1 " " lat2 "set lon "lon1 " " lon2 "set ylpos 0 r" "d lat" "q w2xy "lon1 " " lat1 xloc1=subwrd(result,3) yloc1=subwrd(result,6) "q w2xy "lon2 " " lat2 xloc2=subwrd(result,3) yloc2=subwrd(result,6) xloc1=xloc1+0.05 xloc2=xloc2-0.05 yloc1=yloc1+0.05 yloc2=yloc2-0.05 "set line 0" "draw recf "xloc1 " " yloc1 " " xloc2 " " yloc2 "set gxout contour" "set grid off" "set grads off" "set clevs 500" "set mpdraw on" "set map 3 1 10" "set xlab on" "set ylab on" "set xlpos 0 t" "set ylpos 0 r" "set xlopts 1 3 0.25" "set ylopts 1 3 0.25" "set xlint "xinc "set ylint "yinc "set lat "lat1 " " lat2 "set lon "lon1 " " lon2 "d lat" "set dfile 2" "set time "starttime "q dims" rec=sublin(result,5) sstdate=subwrd(rec,6) "define sst=tmpsfc.2-273.15" "set gxout stat" "d sst" rec=sublin(result,7) valid=subwrd(rec,8) if (valid > 0) "set mpdraw on" "set map 15 1 10" "set gxout shaded" "set rgb 48 140 140 140" "set rgb 49 150 150 150" "set rgb 50 160 160 160" "set rgb 51 170 170 170" "set rgb 52 180 180 180" "set rgb 53 190 190 190" "set rgb 54 200 200 200" "set rgb 55 210 210 210" "set rgb 56 220 220 220" "set rgb 57 230 230 230" "set rgb 58 240 240 240" "set rgb 59 250 250 250" "set clevs 2 5 8 11 14 17 20 23 26 28 30" "set ccols 59 58 57 56 55 54 53 52 51 50 49 48" "set ylpos 0 r" "d sst" * "run cbarnlg.gs 1.2 0 5.5 7.10" "set gxout contour" "set ccolor 0" "set clab on" "set clopts 1 3 0.2" "set clevs 2 5 8 11 14 17 20 23 26 28 30" "set ylpos 0 r" "d sst" "set rgb 77 210 180 140" "basemap.gs L 77 1" "set clevs 500" "set mpdraw off" "d lat" endif "set line 0" "draw recf "xloc1 " " yloc2-0.5 " " xloc2 " "yloc2 "set line 1" "draw rec "xloc1 " " yloc2-0.5 " " xloc2 " "yloc2 "set strsiz 0.25" "set string 1 l 3 0" "draw string "xloc1+0.1 " " yloc2-0.25 " "substr(sstdate,1,12) " GFS SST (shaded)" rc=draw(minit,trackfile) "printim phase2-smooth.png x1280 y1024" "quit" function getstorm(reqtime,trackfile) rc=1 data=" " while (rc != 2 & tstamp!=reqtime) out=read(trackfile) rc=sublin(out,1) if (rc != 2) data=sublin(out,2) date=subwrd(data,1) year=substr(date,1,4) month=substr(date,5,2) mo2=getmonth(month) day=substr(date,7,2) hour=substr(date,9,2) tstamp=hour"Z" day mo2 year lat=subwrd(data,3) lon=subwrd(data,4) pres=subwrd(data,5) endif endwhile rc=close(trackfile) if (tstamp != reqtime) lat=-9999 lon=-9999 pres=-9999 endif return(lat " " lon " "pres) function getmonth(month) monstr = "JAN FEB MAR APR MAY JUN JUL AUG SEP OCT NOV DEC" if (month > 0 & month <= 12) mo2=subwrd(monstr,month) else mo2 = "???" endif return(mo2) function getnum(string) monstr = "JAN FEB MAR APR MAY JUN JUL AUG SEP OCT NOV DEC" ii=1 while (ii <= 12) if (subwrd(monstr,ii) = string) num=ii endif ii=ii+1 endwhile if (num < 10) num="0"num endif return(num) function getcol(pres) pres2=5*math_int(0.5+pres/5) if (pres2 < 850 | pres2 > 1020) col=1 else rbcols = "9 14 4 11 5 13 3 10 7 12 8 2 6" val=math_int((1015-pres2)/5) if (val > 13) val=13 endif if (val < 1) col=1 else col=subwrd(rbcols,val) endif endif return(col) function draw(init,trackfile) ndate=getndate(init) rc=1 data=" " ii=0 while (rc != 2) ii=ii+1 out=read(trackfile) rc=sublin(out,1) if (rc != 2) data=sublin(out,2) date=subwrd(data,1) lat=subwrd(data,3) lon=subwrd(data,4) day=substr(date,7,2) hr=substr(date,9,2) press=subwrd(data,5) col=getcol(press) "q w2xy "lon " " lat xloc=subwrd(result,3) yloc=subwrd(result,6) if (ndate < date) ltype=1 else ltype=6 endif "set line "colold" "ltype" 10" "draw line "xloc " " yloc " " xlocold " " ylocold xlocold=xloc ylocold=yloc colold=col endif endwhile rc=close(trackfile) rc=1 data=" " ii=0 while (rc != 2) ii=ii+1 out=read(trackfile) rc=sublin(out,1) if (rc != 2) data=sublin(out,2) date=subwrd(data,1) lat=subwrd(data,3) lon=subwrd(data,4) day=substr(date,7,2) hr=substr(date,9,2) press=subwrd(data,5) col=getcol(press) "q w2xy "lon " " lat xloc=subwrd(result,3) yloc=subwrd(result,6) if (ndate < date) ltype=1 else ltype=6 endif if (hr = 00) "set line "col" "ltype" 10" "draw mark 5 "xloc " " yloc " 0.6" "set string 0 c 3 0" "set strsiz 0.3" if (day < 10) day=substr(day,2,1) endif "draw string "xloc " " yloc " "day "set line 1 1 3" "draw mark 4 "xloc " " yloc " 0.6" endif if (ii = 1) xloca=xloc yloca=yloc endif xlocold=xloc ylocold=yloc colold=col endif endwhile rc=close(trackfile) "set string 1 c 10 0" "set strsiz 1.00" "draw string "xloca " " yloca " A" "draw string "xloc " " yloc " Z" return function getndate(string) hour=substr(string,1,2) day=substr(string,4,2) mon=substr(string,6,3) year=substr(string,9,4) mons = "JAN FEB MAR APR MAY JUN JUL AUG SEP OCT NOV DEC" ii=1 while (subwrd(mons,ii) != mon) ii=ii+1 endwhile if (ii < 10) fii="0"ii else fii=ii endif string2=year fii day hour return(string2) function getsize(radius) size=0.25*radius/500 if (size < 0.075) size=0.075 endif return(size) function getday(time) "set time "time "q time" day=subwrd(result,6) return(day) function getgdtg(dtg) year=substr(dtg,1,4) month=substr(dtg,5,2) day=substr(dtg,7,2) hr=substr(dtg,9,2) month2=getmonth(month) string=hr"Z"day""month2""year return(string)