Drawing a rectangle with lon/lat coordinates

GrADS has intrinsic functions to draw filled rectangles and unfilled rectangles. You probably know these functions:

draw rec xlo ylo xhi yhi
draw recf xlo ylo xhi yhi

The problem with these functions is that they use XY coordinates, and most of our plots are in world coordinates. Of course it is possible to convert between these two, but doing that every time you need to draw a rectangle is really annoying.

So I took sometime to write a script that does that:

draw rec lon_min lon_max lat_min lat_max <fill> <warn on|off>

It needs at least the four coordinates. There are optional parameters to draw a filled or unfilled rectangle, and to turn on or off warning messages.

Here's an example:

ga-> open vegetation.ctl
ga-> d mask
ga-> drawrec 285 310 -12 5
ga-> drawrec 180 240 -5 5 1 

And here's the code:

*
* drawrec.gs
* 
* This function draws a rectangle over a map using World Coordinates,
* i.e., longitude and latitude.
*
* Usage: drawrec lon_min lon_max lat_min lat_max <fill> <warn on|off>
*        <fill> can be 0 or 1
*        <warn on|off> can be on or off
*
* Written by Henrique Barbosa July 2004
* Last update September 05
*
function drawrec(args)

* Get arguments
    Rxa = subwrd(args,1)
    Rxb = subwrd(args,2)
    Rya = subwrd(args,3)
    Ryb = subwrd(args,4)
    fil = subwrd(arqs,5)
    opt = subwrd(args,6)

* Test for parameters
  if (Rxa='' | Rxb='' | Rya='' | Ryb='') 
    say 'usage: drawrec lon_min lon_max lat_min lat_max  
    return 
  endif

* Find allowed dimensions
  'q dims'
  lonlin =sublin(result,2)
  Bxa = subwrd(lonlin,6)
  Bxb = subwrd(lonlin,8)
  latlin =sublin(result,3)
  Bya = subwrd(latlin,6)
  Byb = subwrd(latlin,8)

* Test if user input within allowed dimensions
  if (opt!='off')
    if (RyaByb); say 'warn: lat_min='Rya' is outside range ['Bya','Byb']'; endif
    if (RybByb); say 'warn: lat_max='Ryb' is outside range ['Bya','Byb']'; endif

    if (RxaBxb); say 'warn: lon_min='Rxa' is outside range ['Bxa','Bxb']'; endif
    if (RxbBxb); say 'warn: lon_max='Rxb' is outside range ['Bxa','Bxb']'; endif
  endif 

* Transform world coordinates to xy
  'q w2xy 'Rxa' 'Rya
  x1=subwrd(result,3)
  y1=subwrd(result,6)

  'q w2xy 'Rxb' 'Ryb
  x2=subwrd(result,3)
  y2=subwrd(result,6)

* Draw a Filled/Empty rectangle
  if (fil=1)
    'draw recf 'x1' 'y1' 'x2' 'y2
  else
    'draw rec 'x1' 'y1' 'x2' 'y2''
  endif
return

GrADS pages you should read:

set strsiz
set string
draw string