<html><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; "><br><div><div>On Aug 19, 2009, at 2:27 PM, PiotrD wrote:</div><blockquote type="cite"><div><font class="Apple-style-span" color="#000000"><br></font>Does anybody have script that can make hurricane track from GFS ensembles?<br><font class="Apple-style-span" color="#000000"><font class="Apple-style-span" color="#144FAE"><br></font></font></div></blockquote><br></div><div>I wrote the script below to use TIGGE forecasts to draw forecasted hurricane tracks when Ike was entering the Gulf of Mexico in September 2008. An example of the final product is here:</div><div><br></div><div><a href="http://iges.org/grads/Ike_2008090912_12z09sep2008_ecmf.png">http://iges.org/grads/Ike_2008090912_12z09sep2008_ecmf.png</a></div><div><a href="http://iges.org/grads/Ike_2008090912_12z09sep2008_kwbc.png">http://iges.org/grads/Ike_2008090912_12z09sep2008_kwbc.png</a></div><div>&nbsp;</div><div>You'll have to modify the script a bit to allow for your data layout. Basically, the script take the minimum of SLP over a small area, using minloc to find the grid location, converting it to x,y coordinates, drawing a dot for each forecast interval, then connecting all the dots for each ensemble member. If you compare the two images, you'll see what a difference grid resolution can make. The ECMWF (ecmf) grid is ~0.45 degrees, NCEP (kwbc) is 1.0.&nbsp;</div><div><br></div><div><br></div><div><span class="Apple-style-span" style="font-family: Times; "><pre style="word-wrap: break-word; white-space: pre-wrap; ">function main(arg)
'reinit'
*initdate='2008090912'
initdate=subwrd(arg,1)
firsttime = '12z09sep2008'
lasttime  = '00z13sep2008'
path='/put_your_path_here/ctl/'initdate
models='ammc babj cwao ecmf egrr kwbc rksl sbsj'
colors='2    6    8    14   4    11   3    13'
nmodels=8
models='kwbc babj cwao ecmf sbsj'
colors='11   2    8    3    14  '
nmodels=4


minlat=20
maxlat=31
minlon=261
maxlon=281

'set rgb 98 66 110 44'
'set rgb 99 60 60 60'

verb=0
threshold=1008
marksize=0.04

'set warn off'
cnt=0
m=1
while(m&lt;=nmodels)
 model=subwrd(models,m)
 color=subwrd(colors,m)
 'open 'path'/'model'/sl.ctl'
 'set dfile 'm
 'set t 2'
 'set lon 'minlon' 'maxlon
 'set lat 'minlat' 'maxlat
 'q file 'm
 sizes = sublin(result,5)
 esize = subwrd(sizes,15)
 if (cnt=0)
  'clear'
  'set clab off'
  'set gxout contour'
  'set ccolor 0'
  'set grads off'
  'set grid on'
  'set map 1'
  'set xlint 3'
  'set ylint 2'
  'set mpdset hires'
  'd msl(e=1)'
  'basemap L 99 98 hires'
  'printim bg0.png x1100 y850'
  'q gxinfo'
  xlims=sublin(result,3)
  ylims=sublin(result,4)
  gxmin=subwrd(xlims,4)
 endif
 if (1)
  'clear'
  'set clab off'
  'set gxout contour'
  'set ccolor 0'
  'set grads off'
  'set grid off'
  'set annot 0'
  'set map 0'
  'd msl(e=1)'
 endif
 e=1
 while (e&lt;=esize)
  'set dfile 'm
  'set e 'e
  'set time 'firsttime' 'lasttime
  'q dims'
  times=sublin(result,5)
  t1=subwrd(times,11)
  t2=subwrd(times,13)
  t=t1
  slp='maskout(msl/100,'threshold'-msl/100)'
  cmd='draw line'
  while (t&lt;=t2)
    'set t 't
*   get the grid locations of the slp minimum
    'd minloc(min('slp',lon='minlon',lon='maxlon'),lat='minlat',lat='maxlat')'
    line1=sublin(result,1)
    line2=sublin(result,2)
    word1=subwrd(line1,1)
    if (word1="Notice:")
      yloc=subwrd(line2,4)
    else
      yloc=subwrd(line1,4)
    endif
    'd minloc(min('slp',lat='minlat',lat='maxlat'),lon='minlon',lon='maxlon')'
    line1=sublin(result,1)
    line2=sublin(result,2)
    word1=subwrd(line1,1)
    if (word1="Notice:")
      xloc=subwrd(line2,4)
    else
      xloc=subwrd(line1,4)
    endif
    if (xloc &lt; 9.999e+20 &amp; yloc &lt; 9.999e+20 &amp; xloc &gt; -9.999e+30 &amp; yloc &gt; -9.999e+30)
*     convert these to world coordinates
      'q gr2xy 'xloc' 'yloc
      line1=sublin(result,1)
      line2=sublin(result,2)
      word1=subwrd(line1,1)
      if (word1="Notice:")
        xpos = subwrd(line2,3)
        ypos = subwrd(line2,6)
      else
        xpos = subwrd(line1,3)
        ypos = subwrd(line1,6)
      endif
      if (xpos != gxmin) 
        'set line 'color
        'draw mark 3 'xpos' 'ypos' 'marksize
        cmd = cmd' 'xpos' 'ypos
        if (verb); say 'm='m' xloc='xloc' yloc='yloc'  xpos='xpos' ypos='ypos; endif
      endif
      cnt=cnt+1
    endif
    t=t+2
  endwhile
  cmd ;* draw the line
  e=e+1
 endwhile
 'printim Ike_'initdate'_'firsttime'_'model'.png -t 0 -b bg0.png x1100 y850' 
 m=m+1
endwhile
say 'cnt='cnt

</pre></span></div><div><br></div></body></html>