;+
;
; @file_comments
; Choose the grid which must be used to do the graph in function of
; vargrid and send back corresponding parameters calculated in
; domdef and reduced at the domain defined by
; domdef (contrarily to
; grandegrille)
; BEWARE!! The choice of the grid is made from the value of the
; global variable vargrid, which can be equal to 'T', 'U', 'V', 'W' ou 'F'.
;
; @categories
; Grid
;
; @keyword TRI
; This keyword serve to obtain the triangulation which
; refer to the grid but only on the part of the zoom. This array of triangulation
; is passed in the variable we have equate at TRI.
; For example: grille,...,tri=triangulation_reduite.
; This keyword is used in plt.
;
; @keyword WDEPTH
; To specify that the field is at W depth instead of T
; depth (automatically activated if vargrid eq 'W')
;
; @keyword FORPLT
; In plt, we want that land points,
; glam and gphi, be equal to glamt and
; gphit regardless of the grid.
;
; @keyword NOTRI
; Useful only when TRI is activated. In this case, grille send back -1 in the
; variable tri even if the variable of the common triangles_list is defined
; and different of -1
;
; @keyword _EXTRA
; Used to pass keywords
;
; @keyword TOUT
;
; @param MASK {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param GLAM {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param GPHI {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param GDEP {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param NX {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param NY {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param NZ {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param FIRSTX {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param FIRSTY {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param FIRSTZ {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param LASTX {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param LASTY {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param LASTZ {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param E1 {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param E2 {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @param E3 {out}{optional}
; For the definition, see domdef and the management of subdomains on the web.
;
; @uses
; cm_4mesh
; cm_4data
;
; @restrictions
; Use the variable vargrid
;
; @restrictions
; Vargrid must be 'T', 'W', 'U', 'V' or 'F'
;
; @history
; Sebastien Masson (smasson\@lodyc.jussieu.fr)
; 12/2/1999
; 10/11/1999 /forplt
;
; @version
; $Id$
;
; @todo
; Comment ecrire la remarque sur les inputs?
;
;-
PRO grille, mask, glam, gphi, gdep, nx, ny, nz $
, firstx, firsty, firstz, lastx, lasty, lastz $
, e1, e2, e3 $
, TRI=tri, NOTRI=notri, TOUT=tout, FORPLT=forplt $
, TYPE=type, WDEPTH=wdepth, _EXTRA=ex
;
compile_opt idl2, strictarrsubs
;
@cm_4mesh
@cm_4data
@cm_4cal ; for jpt
IF NOT keyword_set(key_forgetold) THEN BEGIN
@updatenew
ENDIF
;---------------------
tempsun = systime(1) ; For key_performance
;------------------------------------------------------------
vargrid = strupcase(strmid(vargrid,0,/reverse_offset))
;
if vargrid eq 'W' then wdepth = 1
if keyword_set(tout) then begin
savedbox = 1b
saveboxparam, 'boxparam4grille.dat'
domdef, gridtype = vargrid, _EXTRA = ex
endif
tempdeux = systime(1) ; For key_performance =2
;------------------------------------------------------------
;------------------------------------------------------------
IF keyword_set(wdepth) THEN BEGIN
firstz = firstzw
lastz = lastzw
nz = nzw
ENDIF ELSE BEGIN
firstz = firstzt
lastz = lastzt
nz = nzt
ENDELSE
;------------------------------------------------------------
;------------------------------------------------------------
CASE 1 OF
;------------------------------------------------------------
; grid T and W
;------------------------------------------------------------
vargrid eq 'T' OR vargrid eq 'W' : begin
;scalars
nx = nxt
ny = nyt
firstx = firstxt
firsty = firstyt
lastx = lastxt
lasty = lastyt
;2d vectors
IF arg_present(glam) THEN glam = glamt[firstx:lastx, firsty:lasty]
IF arg_present(gphi) THEN gphi = gphit[firstx:lastx, firsty:lasty]
IF arg_present(e1) THEN e1 = e1t[firstx:lastx, firsty:lasty]
IF arg_present(e2) THEN e2 = e2t[firstx:lastx, firsty:lasty]
;3d vectors
IF keyword_set(forplt) THEN mask = tmask[firstx:lastx, firsty:lasty, firstz] $
ELSE IF arg_present(mask) THEN mask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
end
;------------------------------------------------------------
; grid U
;------------------------------------------------------------
vargrid eq 'U': begin
;scalars
nx = nxu
ny = nyu
firstx = firstxu
firsty = firstyu
lastx = lastxu
lasty = lastyu
;2d vectors
IF arg_present(glam) THEN glam = glamu[firstx:lastx, firsty:lasty]
IF arg_present(gphi) THEN gphi = gphiu[firstx:lastx, firsty:lasty]
if keyword_set(forplt) then BEGIN
mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
eastboarder = mask-shift(mask, 1, 0)*mask
westboarder = mask-shift(mask, -1, 0)*mask
if key_periodic NE 1 OR nx NE jpi then westboarder[nx-1, *] = 0b
tmp1 = shift(eastboarder, 0, 1)
tmp1[*, 0] = 0b
tmp2 = shift(eastboarder, 0, -1)
tmp2[*, ny-1] = 0b
add = (temporary(tmp1)+temporary(tmp2))*(1b-eastboarder)*(1b-temporary(westboarder))
eastboarder = temporary(eastboarder)+temporary(add)
tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
tmp1[*, ny-1] = 1b
tmp1[*, 0] = 1b
tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
if key_periodic NE 1 OR nx NE jpi then begin
tmp2[nx-1, *] = 1b
tmp2[0, *] = 0b
endif
no1 = temporary(tmp1)*temporary(tmp2)
tmp = temporary(eastboarder)*temporary(no1)*mask
mask[0:nx-2, *] = 0b
tmp = temporary(tmp)+temporary(mask)
tmp = where(tmp GE 1)
if tmp[0] NE -1 then begin
glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
endif
ENDIF
IF arg_present(e1) THEN e1 = e1u[firstx:lastx, firsty:lasty]
IF arg_present(e2) THEN e2 = e2u[firstx:lastx, firsty:lasty]
;3d vectors
IF keyword_set(forplt) THEN mask = (umask())[firstx:lastx, firsty:lasty, firstz] $
ELSE IF arg_present(mask) THEN mask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
end
;------------------------------------------------------------
; grid V
;------------------------------------------------------------
vargrid eq 'OPAPTDHV' or vargrid eq 'OPAPT3DV' $
or vargrid eq 'V': begin
;scalars
nx = nxv
ny = nyv
firstx = firstxv
firsty = firstyv
lastx = lastxv
lasty = lastyv
;2d vectors
IF arg_present(glam) THEN glam = glamv[firstx:lastx, firsty:lasty]
IF arg_present(gphi) THEN gphi = gphiv[firstx:lastx, firsty:lasty]
if keyword_set(forplt) then BEGIN
mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
northboarder = mask-shift(mask, 0, 1)*mask
southboarder = mask-shift(mask, 0, -1)*mask
southboarder[*, ny-1] = 0b
tmp1 = shift(northboarder, -1, 0)
if key_periodic NE 1 OR nx NE jpi then tmp1[nx-1, *] = 0b
tmp2 = shift(northboarder, 1, 0)
if key_periodic NE 1 OR nx NE jpi then tmp2[0, *] = 0b
add = (temporary(tmp1)+temporary(tmp2))*(1b-northboarder)*(1b-southboarder)
northboarder = temporary(northboarder)+temporary(add)
tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
tmp1[*, ny-1] = 1b
tmp1[*, 0] = 0b
tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
if key_periodic NE 1 OR nx NE jpi then begin
tmp2[nx-1, *] = 1b
tmp2[0, *] = 1b
endif
no1 = temporary(tmp1)*temporary(tmp2)
tmp = temporary(northboarder)*mask*temporary(no1)
mask[*, 0:ny-2] = 0b
tmp = temporary(tmp)+temporary(mask)
tmp = where(tmp GE 1)
if tmp[0] NE -1 then begin
glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
endif
ENDIF
IF arg_present(e1) THEN e1 = e1v[firstx:lastx, firsty:lasty]
IF arg_present(e2) THEN e2 = e2v[firstx:lastx, firsty:lasty]
;3d vecteurs
IF keyword_set(forplt) THEN mask = (vmask())[firstx:lastx, firsty:lasty, firstz] $
ELSE IF arg_present(mask) THEN mask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
end
;------------------------------------------------------------
; grid F
;------------------------------------------------------------
vargrid eq 'OPAPTDHF' or vargrid eq 'OPAPT3DF' $
or vargrid eq 'F': begin
;scalars
nx = nxf
ny = nyf
firstx = firstxf
firsty = firstyf
lastx = lastxf
lasty = lastyf
;2d vectors
IF arg_present(glam) THEN glam = glamf[firstx:lastx, firsty:lasty]
IF arg_present(gphi) THEN gphi = gphif[firstx:lastx, firsty:lasty]
if keyword_set(forplt) then BEGIN
mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
eastboarder = mask-shift(mask, 1, 0)*mask
westboarder = mask-shift(mask, -1, 0)*mask
westboarder[nx-1, *] = 0b
northboarder = mask-shift(mask, 0, 1)*mask
southboarder = mask-shift(mask, 0, -1)*mask
southboarder[*, ny-1] = 0b
tmp1 = shift(northboarder, -1, 0)
if key_periodic NE 1 OR nx NE jpi then tmp1[nx-1, *] = 0b
tmp2 = shift(northboarder, 1, 0)
if key_periodic NE 1 OR nx NE jpi then tmp2[0, *] = 0b
add = (temporary(tmp1)+temporary(tmp2))*(1b-northboarder)*(1b-southboarder)
northboarder = temporary(northboarder)+temporary(add)
tmp1 = shift(eastboarder, 0, 1)
tmp1[*, 0] = 0b
tmp2 = shift(eastboarder, 0, -1)
tmp2[*, ny-1] = 0b
add = (temporary(tmp1)+temporary(tmp2))*(1b-eastboarder)*(1b-temporary(westboarder))
eastboarder = temporary(eastboarder)+temporary(add)
tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
tmp1[*, ny-1] = 1b
tmp1[*, 0] = 1b
tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
if key_periodic NE 1 OR nx NE jpi then begin
tmp2[nx-1, *] = 1b
tmp2[0, *] = 1b
endif
no1 = temporary(tmp1)*temporary(tmp2)
tmp = (temporary(northboarder)+temporary(eastboarder))*mask*temporary(no1)
mask[0:nx-2, *] = 0b
mask[*, 0:ny-2] = 0b
tmp = temporary(tmp)+temporary(mask)
tmp = where(tmp GE 1)
if tmp[0] NE -1 then begin
glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
endif
ENDIF
IF arg_present(e1) THEN e1 = e1f[firstx:lastx, firsty:lasty]
IF arg_present(e2) THEN e2 = e2f[firstx:lastx, firsty:lasty]
;3d vectors
IF keyword_set(forplt) THEN mask = (fmask())[firstx:lastx, firsty:lasty, firstz] $
ELSE IF arg_present(mask) THEN mask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
END
;------------------------------------------------------------
ELSE:BEGIN
ras = report('Wrong definition of vargrid = '+vargrid+'. Only T, U, V, W or F are acceptable')
stop
END
ENDCASE
IF testvar(var = key_performance) EQ 2 THEN $
print, 'temps grille: attribution des scalaires, vecteurs et tableaux ', systime(1)-tempdeux
;
;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
; Variables refering to the vertical dimension
;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
;
;
tempdeux = systime(1) ; For key_performance =2
if keyword_set(wdepth) then begin
gdep = gdepw[firstz:lastz]
e3 = e3w[firstz:lastz]
endif else begin
gdep = gdept[firstz:lastz]
e3 = e3t[firstz:lastz]
ENDELSE
; for the vertical sections with partial steps
IF keyword_set(type) AND keyword_set(key_partialstep) THEN BEGIN
CASE 1 OF
type EQ 'xz' AND ny EQ 1:BEGIN
bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
good = where(bottom NE 0 AND bottom NE nz+keyword_set(wdepth))
bottom = lindgen(nx)+(bottom-1l+keyword_set(wdepth))*nx
IF good[0] NE -1 THEN BEGIN
bottom = bottom[good]
IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
gdep = replicate(1, nx)#gdep
if keyword_set(wdepth) THEN $
truegdep = hdepw[firstx:lastx, firsty:lasty] $
ELSE truegdep = hdept[firstx:lastx, firsty:lasty]
gdep[bottom] = truegdep[good]
ENDIF
END
type EQ 'yz' AND nx EQ 1:BEGIN
bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
good = where(bottom NE 0 AND bottom NE nz+keyword_set(wdepth))
bottom = lindgen(ny)+(bottom-1l+keyword_set(wdepth))*ny
IF good[0] NE -1 THEN BEGIN
bottom = bottom[good]
IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
gdep = replicate(1, ny)#gdep
if keyword_set(wdepth) THEN $
truegdep = hdepw[firstx:lastx, firsty:lasty] $
ELSE truegdep = hdept[firstx:lastx, firsty:lasty]
gdep[bottom] = truegdep[good]
ENDIF
END
ELSE:
ENDCASE
ENDIF
; for the vertical sections with roms
IF keyword_set(type) AND n_elements(romszinfos) NE 0 THEN BEGIN
romsdp = romsdepth()
IF romsdp[0] NE -1 THEN BEGIN
IF jpt EQ 1 THEN BEGIN
CASE type OF
'z':gdep = moyenne(temporary(romsdp), 'xy')
'xz':gdep = moyenne(temporary(romsdp), 'y')
'yz':gdep = moyenne(temporary(romsdp), 'x')
'zt':gdep = moyenne(temporary(romsdp), 'xy')
ELSE:
ENDCASE
ENDIF ELSE BEGIN
CASE type OF
'z':gdep = moyenne(temporary(romsdp), 'xyt')
'xz':gdep = grossemoyenne(temporary(romsdp), 'yt')
'yz':gdep = grossemoyenne(temporary(romsdp), 'xt')
'zt':gdep = grossemoyenne(temporary(romsdp), 'xy')
ELSE:
ENDCASE
ENDELSE
IF size(gdep, /n_dimensions) EQ 2 THEN $
gdep = remplit(gdep, niter = 32000, mask = gdep LT 1.E19, /basique, /fillxdir)
ENDIF
ENDIF
;
IF testvar(var = key_performance) EQ 2 THEN $
print, 'temps grille: Variables se rapportant a la dimension verticale ', systime(1)-tempdeux
;------------------------------------------------------------
; Triangulation vector when TRI is activated.
;------------------------------------------------------------
if arg_present(TRI) then $
if triangles_list[0] EQ -1 OR keyword_set(notri) then tri = -1 ELSE BEGIN
tempdeux = systime(1) ; pour key_performance =2
msk = bytarr(jpi, jpj)
msk[firstx:lastx, firsty:lasty] = 1
ind = where( msk[triangles_list[0, *]]*msk[triangles_list[1, *]]*msk[triangles_list[2, *]] EQ 1 )
tri = triangles_list[*, ind]-(firstx+firsty*jpi)
y = tri/jpi
x = tri-y*jpi
tri = x+y*nx
IF testvar(var = key_performance) EQ 2 THEN $
print, 'temps grille: decoupage de la triangulation ', systime(1)-tempdeux
ENDELSE
;------------------------------------------------------------------
; To make sure there is not any degenerated dimension (=1)
;------------------------------------------------------------------
; mask=reform(mask, /over)
; glam=reform(glam, /over)
; gphi=reform(gphi, /over)
; gdep=reform(gdep, /over)
; e1=reform(e1, /over)
; e2=reform(e2, /over)
; e3=reform(e3, /over)
if keyword_set(savedbox) THEN restoreboxparam, 'boxparam4grille.dat'
if keyword_set(key_performance) THEN print, 'temps grille', systime(1)-tempsun
;------------------------------------------------------------
IF NOT keyword_set(key_forgetold) THEN BEGIN
@updateold
ENDIF
;---------------------
return
end