Re: Draw coastal line in different sea depth levels.

[Mick Spillane <spillane@xxxxxxxxxxxxx>]

Here is a way to draw the boundary of a mask using straight line
segments for the top, bottom, left, and right edges of grid cells
lying along the "land-sea" border.

Example: The Gulf of Alaska

use etopo20
region/x=180:270/y=45:60
shade rose

let mask=if(rose[d=1] gt 0)then 1 else 0
go bounder

where bounder.jnl is the following script

 ! bounder.jnl : drawing the "land-sea" edge in a mask

def grid/like=mask gm
let xl=x[g=gm]-xbox[g=gm]/2 ; let xr=x[g=gm]+xbox[g=gm]/2
let yb=y[g=gm]-ybox[g=gm]/2 ; let yt=y[g=gm]+ybox[g=gm]/2
! the left edge
let ylb=if(mask[x=@ddb] ne 0)then yb else 999
let ylt=if(mask[x=@ddb] ne 0)then yt else 999
! the right edge
let yrb=if(mask[x=@ddf] ne 0)then yb else 999
let yrt=if(mask[x=@ddf] ne 0)then yt else 999
! the bottom edge
let xbl=if(mask[y=@ddb] ne 0)then xl else 999
let xbr=if(mask[y=@ddb] ne 0)then xr else 999
! the top edge
let xtl=if(mask[y=@ddf] ne 0)then xl else 999
let xtr=if(mask[y=@ddf] ne 0)then xr else 999

! list all edges to the file bounds.d ...
list/nohead/form=(4f8.2,' 999 999')/file=bounds.d/clob xl,ylb,xl,ylt
list/nohead/form=(4f8.2,' 999 999')/file=bounds.d/app  xr,yrb,xr,yrt
list/nohead/form=(4f8.2,' 999 999')/file=bounds.d/app  xbl,yb,xbr,yb
list/nohead/form=(4f8.2,' 999 999')/file=bounds.d/app  xtl,yt,xtr,yt
! ... then use grep to exclude the interior cells ...
sp rm bounds2.d                ! this may require the user to type y
sp grep -v 999.00 bounds.d > bounds2.d

file/form=free/var=xx,yy/col=6 bounds2.d
set var/bad=999 xx ; set var/bad=999 yy
can region
plot/o/nolab/vs/line=1 xx,yy

------
Note that if a very large region were set, or the mask had a very fine
resolution, the number of points in bounds2.d might exceed the 20000
or so default limit and the read statement would need to be replaced
with a gridded read.

Mick Spillane


|--****--****-*---*---***--***--|____spillane@pmel.noaa.gov____|
|-*__---*-----*--*-*--*--*-*--*-|_SCIENCE APPLICATIONS SUPPORT_|
|--***--*-----*-*---*-***--***--|____EPIC/Ferret/PlotPlus______|
|-----*-*-----*-*****-*----*----|__Room 2070 Bldg#3 NOAA/PMEL__|
|-****---****-*-*---*-*----*----|____Phone_:_(206)526-6780_____|