;+
;
; @file_comments
; Draw the grid
;
; @categories
; Grid
;
; @param GLAMIN {in}{optional}{type=1d or 2d array}{default=glam specified by vargrid, on the domain defined by domdef}
; Longitude of points of the grid.
;
; @param GPHIIN {in}{optional}{type=1d or 2d array}{default=gphi specified by vargrid, on the domain defined by domdef}
; Latitude of points of the grid.
;
; @keyword XSTRIDE {type=integer}{default=1}
; It specify that we want to trace only one line of
; constant i every xstride points
;
; @keyword YSTRIDE {type=integer}{default=1}
; It specify that we want to trace only one line of
; constant j every ystride points
;
; @keyword OCEAN
; To trace the grid only on ocean points.
;
; @keyword EARTH
; To trace the grid only on land points.
;
; @keyword RMOUT
; Select to remove all cell having one corner out of the
; plot boundaries (!x.range, !y.range)
;
; @keyword _EXTRA
; Used to pass keywords
;
; @uses
; common
;
; @examples
;
; IDL> plt,indgen(jpi,jpj),/nocontour,/nofill
; IDL> vargrid='T'
; IDL> tracegrille,/ocean,color=20
; IDL> tracegrille,/earth,color=80
;
; @history
; Sebastien Masson (smasson\@lodyc.jussieu.fr)
;
; @version
; $Id$
;
;-
PRO tracegrille, glamin, gphiin, OCEAN=ocean, EARTH=earth $
, XSTRIDE=xstride, YSTRIDE=ystride, RMOUT=rmout $
, _EXTRA=extra
;
compile_opt idl2, strictarrsubs
;
@cm_4mesh
@cm_4data
IF NOT keyword_set(key_forgetold) THEN BEGIN
@updatenew
ENDIF
;---------------------------------------------------------
tempsun = systime(1) ; For key_performance
; to avoid warning message
oldexcept = !except
!except = 0
if n_elements(key_gridtype) EQ 0 then key_gridtype = 'c'
;
if n_elements(glamin) * n_elements(gphiin) EQ 0 then BEGIN
grille, mask, glam, gphi, gdep, nx, ny, nz, firstx, firsty, firstz, lastx, lasty, lastz
IF keyword_set(ocean) AND strmid(key_gridtype, 0, 1) EQ 'c' THEN BEGIN
; we reduce the mask to take into account the point located ON the coastline.
CASE vargrid OF
'U':BEGIN
mask = tmask[firstx:lastx, firsty:lasty]
IF NOT keyword_set(key_periodic) OR nx NE jpi $
THEN tmpx = mask[nx-1, *]
mask = (mask+shift(mask, -1, 0)) < 1
IF NOT keyword_set(key_periodic) OR nx NE jpi $
THEN mask[nx-1, *] = temporary(tmpx)
END
'V':BEGIN
mask = tmask[firstx:lastx, firsty:lasty]
tmpy = mask[*, ny-1]
mask = (mask+shift(mask, 0, -1)) < 1
mask[*, ny-1] = temporary(tmpy)
END
'F':BEGIN
mask = tmask[firstx:lastx, firsty:lasty]
IF NOT keyword_set(key_periodic) OR nx NE jpi $
THEN tmpx = mask[nx-1, *]
tmpy = mask[*, ny-1]
mask = (mask+shift(mask, -1, 0)+shift(mask, 0, -1)+shift(mask, -1, -1)) < 1
mask[*, ny-1] = temporary(tmpy)
IF NOT keyword_set(key_periodic) OR nx NE jpi $
THEN mask[nx-1, *] = temporary(tmpx)
END
ELSE:
ENDCASE
ENDIF
ENDIF ELSE BEGIN
glam = glamin
gphi = gphiin
IF (size(glam))[0] EQ 1 AND (size(gphi))[0] EQ 1 THEN BEGIN
nx = n_elements(glam)
ny = n_elements(gphi)
glam = glam#replicate(1, ny)
gphi = replicate(1, nx)#gphi
ENDIF ELSE BEGIN
nx = (size(glam))[1]
ny = (size(glam))[2]
ENDELSE
ENDELSE
if n_elements(mask) EQ 0 then mask = replicate(1b, nx, ny)
if (size(mask))[0] EQ 3 then mask = mask[*, *, 0]
;
IF keyword_set(RMOUT) THEN BEGIN
out = where(glam GT max(!x.range) OR glam LT min(!x.range) $
OR gphi GT max(!y.range) OR gphi LT min(!y.range))
IF out[0] NE -1 THEN BEGIN
glam[out] = !values.f_nan
gphi[out] = !values.f_nan
ENDIF
ENDIF
;
IF keyword_set(ocean) then BEGIN
earth = where(mask EQ 0)
if earth[0] NE -1 then begin
glam[earth] = !values.f_nan
gphi[earth] = !values.f_nan
ENDIF
earth = 0
ENDIF
;
IF keyword_set(earth) THEN BEGIN
ocean = where(mask EQ 1)
if ocean[0] NE -1 then begin
glam[ocean] = !values.f_nan
gphi[ocean] = !values.f_nan
ENDIF
ocean = 0
ENDIF
;
if NOT keyword_set(xstride) then xstride = 1
if NOT keyword_set(ystride) then ystride = 1
case strmid(key_gridtype, 0, 1) of
'c':BEGIN
for i = 0, ny-1, ystride do begin
plots, glam[*, i], gphi[*, i], _extra = extra
endfor
for i = 0, nx-1, xstride do begin
plots, glam[i, *], gphi[i, *], _extra = extra
endfor
END
'e':BEGIN
shifted = glam[0, 0] LT glam[0, 1]
glam2 = glam+(glam[1]-glam[0])/2.
if shifted then begin
for i = 0, ny-2 do BEGIN
xx = (transpose([[glam[*, i]], [glam2[*, i]]]))[*]
yy = (transpose([[gphi[*, i]], [gphi[*, i+1]]]))[*]
plots, xx[0:2*nx-2], yy[0:2*nx-2], _extra = extra
ENDFOR
ENDIF ELSE BEGIN
for i = 1, ny-1 do BEGIN
xx = (transpose([[glam[*, i]], [glam2[*, i]]]))[*]
yy = (transpose([[gphi[*, i]], [gphi[*, i-1]]]))[*]
plots, xx[0:2*nx-2], yy[0:2*nx-2], _extra = extra
ENDFOR
ENDELSE
for i = 1, (ny-1)/2 do $
plots, [glam[0, 2*i-1], glam[0, 2*i]] $
, [gphi[0, 2*i-1], gphi[0, 2*i]], _extra = extra
for i = 0, (ny-2)/2 do $
plots, [glam[nx-1, 2*i], glam[nx-1, 2*i+1]] $
, [gphi[nx-1, 2*i], gphi[nx-1, 2*i+1]], _extra = extra
END
endcase
if keyword_set(key_performance) THEN print, 'temps trace grille', systime(1)-tempsun
!except = oldexcept
return
end