Index: /trunk/SRC/Computation/bsf.pro
===================================================================
--- /trunk/SRC/Computation/bsf.pro	(revision 388)
+++ /trunk/SRC/Computation/bsf.pro	(revision 388)
@@ -0,0 +1,131 @@
+;+
+;
+; @file_comments
+; compute the barotropic (along the model grid) stream function in Sv using
+; the following equation: total(total(e2u*e3u*un),3), 2, /cumulative)/1.e6
+;
+; @categories
+; diagnostics
+;
+; @param z3d {in} {type=3D xyz array}
+; zonal velocity
+;
+; @keyword NOSTRUCTURE {default=0}{type=scalar: 0 or 1}
+; activate if you do not want that msf returns a structure
+; but only the array referring to the field.
+;
+; @keyword REFPOINT {type=2 element vector}
+; a 2 element vector [x,y] giving the (approximative) position of the
+; point that should be taken as a reference. This position should be
+; given with the same coordinates system as the one used in
+; glam/gphi. see example
+;
+; @keyword REFVALUE {type=scalar} {default=0.}
+; the bsf value that we want to speficy at the position defined by refpoint
+;
+; @returns {type=2D xy array}
+; barotropic stream function in Sv
+; 
+; @uses
+; <pro>cm_4mesh</pro>
+; <pro>cm_4data</pro>
+; <pro>grille</pro>
+; <pro>litchamp</pro>
+; <pro>fitintobox</pro>
+;
+; @restrictions
+; - T and U boxzoom parameters must be the same
+; - to be valid, mask must be equal to 0 at the bottom and on each
+;   side of the domain along x direction
+; - define the common variables (of cm_4data) 
+;       varname = 'BSF'
+;       varunit = 'Sv'
+;       vargrid = 'F' 
+;
+; @examples
+;
+;   IDL> initorca05
+;   IDL> file = '/Volumes/pim/F31/oce/1y/F31_1y_0040_0040_grid_U.nc'
+;   IDL> un = read_ncdf('vozocrtx', 0, 0, /timestep, file = file)
+;   IDL> bar = bsf(un, refvalue = 0., refpoint = [0, 20])
+;   IDL> plt, bar, -150, 150, int = 10, style = 'so0so', format = '(i4)', lct = 64
+;
+; @history
+; Gurvan Madec, Christian Ethe, Claude Talandier and Sebastien Masson, Dec 2008
+;
+; @version
+; $Id$
+;
+;-
+;
+FUNCTION bsf, z3d, NOSTRUCTURE = nostructure, REFPOINT = refpoint, REFVALUE = refvalue
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4mesh
+@cm_4data
+;
+  CASE 1 OF
+    firstxt NE firstxu OR lastxt NE lastxu: return, report('T and U box must be defined over the same x indexes, use /memeindices when calling domdef')
+    firstyt NE firstyu OR lastyt NE lastyu: return, report('T and U box must be defined over the same y indexes, use /memeindices when calling domdef')
+    ELSE:
+  ENDCASE
+;
+  grille, -1, -1, -1, -1, nx, ny, nz, firstx, firsty, firstz, lastx, lasty, lastz, grid = 'T'
+  umsk = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
+  IF total(umsk[*, *, nz-1]) NE 0 THEN return, report('The bottom of the box (defined by lastzt) must be set to land value (=0)')
+  IF total(umsk[*, 0, *]) NE 0 THEN return, report('southern boudary must be land point if you want to compute the bsf with this function')
+;
+  un = litchamp(z3d)
+  IF size(un, /n_dimensions) NE 3 THEN return, report('input data must be a 3D array')
+  un = fitintobox(temporary(un))
+;
+  IF keyword_set(key_partialstep) THEN BEGIN
+; 2D e3u at the bottom of the ocean
+; see zgr_zps: e3u (ji,jj,jk) = MIN( e3t(ji,jj,jk), e3t(ji+1,jj,jk))
+    e3u_ps = [[[e3t_ps[firstx:lastx, firsty:lasty]]], [[shift(e3t_ps[firstx:lastx, firsty:lasty], -1, 0)]]]
+    e3u_ps = min(temporary(e3u_ps), dimension = 3)
+    flagdata = NOT (keyword_set(key_periodic) AND nx EQ jpi) AND total(umsk[nx-1, *, *]) NE 0
+; level of the bottom of the ocean
+    bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
+    bottom = 0L > ( long(temporary(bottom)) - 1L )
+; the bottom of the ocean in 3D index is:
+    bottom = lindgen(nx*ny) + temporary(bottom)*nx*ny
+; 3D e3t array
+    e33 = replicate(1., nx*ny) # e3t[firstz:lastz]
+; apply e3u_ps to e33 at the bottom of the ocean
+    e33[temporary(bottom)] = e3u_ps
+; 3D e2u array
+    e23d = (e2u[firstx:lastx, firsty:lasty])[*] # replicate(1., nzt)
+; e2u*e3u
+    e23 = temporary(e23d) * temporary(e33) 
+  ENDIF ELSE BEGIN 
+    e23 = (e2u[firstx:lastx, firsty:lasty])[*] # e3t[firstz:lastz]
+  ENDELSE
+;
+; mask the array
+  un = temporary(umsk) * temporary(un)
+; compute the bsf
+  bsf = 1.e-6 * total(total(temporary(un) * temporary(e23), 3), 2, /cumulative)
+; set bsf to 0 in the largest continent... no done...
+  IF keyword_set(refpoint) THEN BEGIN
+    refind = neighbor(refpoint[0], refpoint[1], glamf[firstx:lastx, firsty:lasty], gphif[firstx:lastx, firsty:lasty], SPHERE = key_onearth)
+    IF n_elements(refvalue) EQ 0 THEN refval = - bsf[refind[0]] ELSE refval = refvalue - bsf[refind[0]] 
+    bsf = temporary(bsf) + refval
+  ENDIF
+  IF keyword_set(flagdata) THEN bsf[nx-1, *] = !values.f_nan
+;
+; update data informations
+  varname = 'BSF'
+  varunit = 'Sv'
+  vargrid = 'F'
+;
+  IF keyword_set(nostructure) THEN return, bsf
+
+  IF keyword_set(key_forgetold) THEN BEGIN
+    return, {arr:temporary(bsf), grid:vargrid, unit:varunit, experiment:varexp, name:varname}
+  ENDIF ELSE BEGIN
+    return, {tab:temporary(bsf), grille:vargrid, unite:varunit, experience:varexp, nom:varname}
+  ENDELSE
+
+END
Index: /trunk/SRC/Computation/msf.pro
===================================================================
--- /trunk/SRC/Computation/msf.pro	(revision 388)
+++ /trunk/SRC/Computation/msf.pro	(revision 388)
@@ -0,0 +1,173 @@
+;+
+;
+; @file_comments
+; compute the meridional (along the model grid) stream function in Sv using
+; the following equation: total(total(e1v*e3v*vn),1), 2, /cumulative)/1.e6
+;
+; @categories
+; diagnostics
+;
+; @param z3d {in} {type=3D xyz array}
+; meridional velocity
+;
+; @param mask2d {in} {type=2D xy array} {optional} {default=1 everywhere}
+; land/sea (0/1) mask to be applied (in addtion to the 3D mask
+; of the model) to mask some parts of the domain (and keep for example
+; the northend atlantic).
+;
+; @keyword NOSTRUCTURE {default=0}{type=scalar: 0 or 1}
+; activate if you do not want that msf returns a structure
+; but only the array referring to the field.
+;
+; @keyword INDEXBOXZOOM
+; Set this keyword to a named variable in which msf will return the x
+; and y indexes definig the zoom-box that should be used by pltz to do
+; the plot (see example)
+;
+; @keyword MASKOUT
+; Set this keyword to a named variable in which msf will return the 2d
+; array (yz) should be used by pltz to do plot the land/sea mask(see example)
+;
+; @returns {type=2D yz array}
+; meridional stream function in Sv
+; 
+; @uses
+; <pro>cm_4mesh</pro>
+; <pro>cm_4data</pro>
+; <pro>grille</pro>
+; <pro>litchamp</pro>
+; <pro>fitintobox</pro>
+;
+; @restrictions
+; - T and V boxzoom parameters must be the same
+; - to be valid, mask must be equal to 0 at the bottom and on each
+;   side of the domain along x direction
+; - define the common variables (of cm_4data) 
+;       varname = 'MSF'
+;       varunit = 'Sv'
+;       vargrid = 'W' ; -> should introduce WF point...
+;
+; @examples
+;
+;   IDL> initorca05
+;   IDL> file = '/Volumes/pim/F31/oce/1y/F31_1y_0040_0040_grid_V.nc'
+;   IDL> vn = read_ncdf('vomecrty', 0, 0, /timestep, file = file)
+;   IDL> mm = msf(vn, mskatl, indexboxzoom = ind, maskout = ma)
+;   IDL> pltz, mm, -20, 20, int = 1, boxzoom = [ind, 0, 5200], /xindex  $
+;            , /yindex, zoom = 1000, maskdta = ma, /no_partial, style = 'so0so'
+;
+; @history
+; Gurvan Madec, Christian Ethe, Claude Talandier and Sebastien Masson, Dec 2008
+;
+; @version
+; $Id$
+;
+;-
+;
+FUNCTION msf, z3d, mask2d, INDEXBOXZOOM = indexboxzoom, MASKOUT = maskout, NOSTRUCTURE = nostructure
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4mesh
+@cm_4data
+;
+  CASE 1 OF
+    firstxt NE firstxv OR lastxt NE lastxv: return, report('T and V box must be defined over the same x indexes, use /memeindices when calling domdef')
+    firstyt NE firstyv OR lastyt NE lastyv: return, report('T and V box must be defined over the same y indexes, use /memeindices when calling domdef')
+    ELSE:
+  ENDCASE
+;
+  grille, -1, -1, -1, -1, nx, ny, nz, firstx, firsty, firstz, lastx, lasty, lastz, grid = 'T'
+;
+  vn = litchamp(z3d)
+  IF size(vn, /n_dimensions) NE 3 THEN return, report('input data must be a 3D array')
+  vn = fitintobox(temporary(vn))
+;
+  IF n_elements(mask2d) EQ 0 THEN mask2d = replicate(1b, nx, ny)
+  IF size(mask2d, /n_dimensions) NE 2 THEN return, report('input mask2d must be a 2D array')
+  msk = fitintobox(mask2d)
+  msk = (temporary(msk))[*]#replicate(1b, nz)
+  msk = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz] * temporary(msk)
+;
+  IF total(msk[*, *, nz-1]) NE 0 THEN return, report('The bottom of the box (defined by lastzt) must be set to land value (=0)')
+  IF (total(msk[0, *, *]) NE 0 OR total(msk[nx-1, *, *]) NE 0) AND NOT (keyword_set(key_periodic) AND nx EQ jpi) THEN $
+     return, report('for non-periodic domains, eastern and western bondaries must be land point')
+  IF keyword_set(key_partialstep) AND total(msk[*, ny-1, *]) NE 0 THEN flagdata = 1
+;
+;
+  IF keyword_set(key_partialstep) THEN BEGIN
+
+; Rebuild the T-point 3D partial steps array from 2D bottom e3t_ps
+;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
+; level of the bottom of the ocean
+     bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
+     bottom = 0L > ( long(temporary(bottom)) - 1L )
+; the bottom of the ocean in 3D index is:
+     bottom = lindgen(nx*ny) + temporary(bottom)*nx*ny
+; 3D e3t array
+     e3t_3D = replicate(1, nx*ny) # e3t[firstz:lastz]
+; apply e3t_ps to e3t_3D at the bottom of the ocean
+     e3t_3D[temporary(bottom)] = e3t_ps[firstx:lastx, firsty:lasty]
+; 3D e3t_3D array
+     e3t_3D = reform(e3t_3D, nx, ny, nz, /overwrite)
+     
+; Rebuild the V-point 3D partial steps array from T-point 3D e3t_3D array 
+     tmp = shift(e3t_3D, 0, -1, 0)
+     e3v_3D = [ [ (temporary(e3t_3D))[*] ], [ (temporary(tmp))[*] ] ]
+     e3v_3D = min(temporary(e3v_3D), dimension = 2)
+
+; 3D e1v array
+     e1_3D = (e1v[firstx:lastx, firsty:lasty])[*] # replicate(1., nz)
+; e1v*e3v
+     e13 = temporary(e1_3D) * temporary(e3v_3D) 
+  ENDIF ELSE BEGIN 
+     e13 = (e1v[firstx:lastx, firsty:lasty])[*] # e3t[firstz:lastz]
+  ENDELSE
+;
+; mask the array
+  vn = temporary(vn) * msk
+  IF arg_present(maskout) THEN BEGIN
+    msk = total(temporary(msk), 1)
+    maskout = temporary(msk) NE 0
+    maskout = shift(maskout, 0, 1) ; bottom value of W should be taken into account
+    maskout[*, 0] = maskout[*, 1]
+    maskout = reverse(temporary(maskout), 2)
+  ENDIF
+;  
+; defaut computation:
+;  msf = 1.e-6 * total(total(temporary(e13) * temporary(vn), 1), 2, /cumulative)
+; but we force the computation from bottom to up.
+  msf = total(temporary(vn) * temporary(e13), 1) ; -> yz array
+  msf = reverse(temporary(msf), 2)
+  msf = 1.e-6 * total(temporary(msf), 2, /cumulative)
+  msf = -1.*reverse(temporary(msf), 2)
+  IF keyword_set(flagdata) THEN msf[ny-1, *] = !values.f_nan
+
+; set msf to 0 in the largest continent... no done...
+
+; update data informations
+  varname = 'MSF'
+  varunit = 'Sv'
+  vargrid = 'W' ; -> should introduce WF point...
+;
+; we must define a domain related to the section
+;
+; are we using an ORCA grid???
+
+  IF (jpiglo EQ 182 AND jpjglo EQ 149) OR (jpiglo EQ 722 AND jpjglo EQ 511) $
+     OR (jpiglo EQ 1442 AND jpjglo EQ 1021) OR (jpiglo EQ 4322 AND jpjglo EQ 3059) THEN BEGIN 
+; find the i index with the northest latitude
+    ind = (where(gphit[firstx:lastx, firsty:lasty] EQ max(gphit[firstx:lastx, firsty:lasty])))[0] mod nx
+  ENDIF ELSE ind = 0
+;
+  indexboxzoom = [ind, ind, firsty, lasty]
+;
+  IF keyword_set(nostructure) THEN return, msf
+
+  IF keyword_set(key_forgetold) THEN BEGIN
+    return, {arr:temporary(msf), grid:vargrid, unit:varunit, experiment:varexp, name:varname}
+  ENDIF ELSE BEGIN
+    return, {tab:temporary(msf), grille:vargrid, unite:varunit, experience:varexp, nom:varname}
+  ENDELSE
+
+END
Index: /trunk/SRC/Grid/restoreboxparam.pro
===================================================================
--- /trunk/SRC/Grid/restoreboxparam.pro	(revision 387)
+++ /trunk/SRC/Grid/restoreboxparam.pro	(revision 388)
@@ -73,5 +73,12 @@
      nzt = boxzoomparam.indexes[28]
      nzw = boxzoomparam.indexes[29]
-     key_irregular = boxzoomparam.key
+     key_shift    = boxzoomparam.key[0]
+     key_periodic = boxzoomparam.key[1]
+     key_stride   = boxzoomparam.key[2]
+     key_gridtype = boxzoomparam.key[3]
+     key_yreverse = boxzoomparam.key[4]
+     key_zreverse = boxzoomparam.key[5]
+     key_partialstep = boxzoomparam.key[6]
+     key_onearth  = boxzoomparam.key[7]
 ;
       boxzoomparam = -1
Index: /trunk/SRC/Grid/saveboxparam.pro
===================================================================
--- /trunk/SRC/Grid/saveboxparam.pro	(revision 387)
+++ /trunk/SRC/Grid/saveboxparam.pro	(revision 388)
@@ -42,10 +42,11 @@
       boxzoomparam = {bound:[lon1, lon2, lat1, lat2, vert1, vert2] $
        , indexes:[firstxt, lastxt, firstyt, lastyt $
-       , firstxu, lastxu, firstyu, lastyu $
-       , firstxv, lastxv, firstyv, lastyv $
-       , firstxf, lastxf, firstyf, lastyf $
-       , firstzt, lastzt, firstzw, lastzw $
-       , nxt, nyt, nxu, nyu, nxv, nyv, nxf, nyf, nzt, nzw] $
-       , key:key_irregular}
+                , firstxu, lastxu, firstyu, lastyu $
+                , firstxv, lastxv, firstyv, lastyv $
+                , firstxf, lastxf, firstyf, lastyf $
+                , firstzt, lastzt, firstzw, lastzw $
+                , nxt, nyt, nxu, nyu, nxv, nyv, nxf, nyf, nzt, nzw] $
+       , key:[key_shift, key_periodic, key_stride, key_gridtype $
+            , key_yreverse, key_zreverse, key_partialstep, key_onearth]}
 ;
    ENDIF ELSE BEGIN
@@ -58,5 +59,7 @@
        , firstzt, lastzt, firstzw, lastzw $
        , nxt, nyt, nxu, nyu, nxv, nyv, nxf, nyf, nzt, nzw $
-       , key_irregular, filename = myuniquetmpdir + filename
+       , key_shift, key_periodic, key_stride, key_gridtype $
+       , key_yreverse, key_zreverse, key_partialstep, key_onearth $
+       , filename = myuniquetmpdir + filename
 ;
    ENDELSE
Index: /trunk/SRC/ToBeReviewed/GRILLE/decoupeterre.pro
===================================================================
--- /trunk/SRC/ToBeReviewed/GRILLE/decoupeterre.pro	(revision 387)
+++ /trunk/SRC/ToBeReviewed/GRILLE/decoupeterre.pro	(revision 388)
@@ -34,4 +34,8 @@
 ; @keyword USETRI
 ;
+; @keyword MASKDTA
+; 2d array defining the mask that must be applied to the data instead
+; of the default mask which is defined according to the grid (t/u/v/f mask)
+;
 ; @keyword _EXTRA
 ; Used to pass keywords
@@ -64,5 +68,5 @@
 PRO decoupeterre, mask, glam, gphi, gdep $
                 , TYPE=type, TRI=tri, INDICEZOOM=indicezoom $
-                , COINMONTE=coinmonte, COINDESCEND=coindescend $
+                , COINMONTE=coinmonte, COINDESCEND=coindescend, MASKDTA = maskdta $
                 , WDEPTH=wdepth, REALSECTION=realsection, USETRI=usetri $
                 , _EXTRA=ex
@@ -193,14 +197,25 @@
 ;------------------------------------------------------------
   if keyword_set(wdepth)  then begin
-    firstz = 0 > (min([firstzt, firstzw])-1)
-    lastz = (max([lastzt, lastzw])+1) < (jpk-1)
+    firstzw = 0 > (min([firstzt, firstzw])-1)
+    lastzw = (max([lastzt, lastzw])+1) < (jpk-1)
+    nzw = lastzw - firstzw + 1
+    firstz = firstzw
+    lastz = lastzw
+    nz = nzw
   ENDIF ELSE BEGIN
     firstz = firstzt
     lastz = lastzt
+    nz = nzt
   ENDELSE
-  nz = lastz-firstz+1
 ;------------------------------------------------------------
 ; mask
 ;------------------------------------------------------------
+  IF keyword_set(maskdta) THEN BEGIN
+    mask = byte(fitintobox(maskdta))
+    IF keyword_set(wdepth) THEN BEGIN
+      mask = shift(temporary(mask), 0, 1)
+      mask[*, 0] = 0b
+    ENDIF
+  ENDIF ELSE BEGIN 
   case type of
     'xy':BEGIN
@@ -221,5 +236,6 @@
     END
     ELSE:mask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
-  endcase
+  ENDCASE
+ENDELSE
 ;------------------------------------------------------------
 ; vertical axis
@@ -240,11 +256,12 @@
 ; in that case, the z axis is a 2d array and we modify the depth of
 ; the last level ocean with hdepw that is the real depth of the bottom.
+; in domzgr:  hdepw(ji,jj) = gdepw(ji,jj,ik+1)
     CASE 1 OF
-      keyword_set(key_partialstep) and type EQ 'xz' $
+      keyword_set(key_partialstep) AND type EQ 'xz' $
         AND ny EQ 1 AND keyword_set(realsection):BEGIN
-        bottom = total(mask, 3)
+        bottom = total(mask, 3 - keyword_set(maskdta))
         good = where(bottom NE 0 AND bottom NE nz)
-        bottom = lindgen(nx)+(bottom)*nx
         IF good[0] NE -1 THEN BEGIN
+          bottom = lindgen(nx)+(bottom)*nx
           bottom = bottom[good]
           gdep = replicate(1, nx)#gdep
@@ -253,10 +270,10 @@
         ENDIF
       END
-      keyword_set(key_partialstep) and type EQ 'yz' $
+      keyword_set(key_partialstep) AND type EQ 'yz' $
         AND nx EQ 1 AND keyword_set(realsection):BEGIN
-        bottom = total(mask, 3)
+        bottom = total(mask, 3 - keyword_set(maskdta))
         good = where(bottom NE 0 AND bottom NE nz)
-        bottom = lindgen(ny)+(bottom)*ny
         IF good[0] NE -1 THEN BEGIN
+          bottom = lindgen(ny)+(bottom)*ny
           bottom = bottom[good]
           gdep = replicate(1, ny)#gdep
Index: /trunk/SRC/ToBeReviewed/GRILLE/grille.pro
===================================================================
--- /trunk/SRC/ToBeReviewed/GRILLE/grille.pro	(revision 387)
+++ /trunk/SRC/ToBeReviewed/GRILLE/grille.pro	(revision 388)
@@ -21,5 +21,8 @@
 ; This keyword is used in <pro>plt</pro>.
 ;
-; @keyword WDEPTH
+; @keyword GRID {default=vargrid defined in <pro>cm_4data</pro>} 
+; To specify on which grid are located the data
+;
+; @keyword WDEPTH {default=0} 
 ; To specify that the field is at W depth instead of T
 ; depth (automatically activated if vargrid eq 'W')
@@ -114,5 +117,5 @@
           , e1, e2, e3 $
           , TRI=tri, NOTRI=notri, TOUT=tout, FORPLT=forplt $
-          , TYPE=type, WDEPTH=wdepth, _EXTRA=ex
+          , TYPE=type, WDEPTH=wdepth, GRID = grid, _EXTRA=ex
 ;
   compile_opt idl2, strictarrsubs
@@ -129,4 +132,5 @@
   vargrid = strupcase(strmid(vargrid,0,/reverse_offset))
 ;
+  IF size(grid, /type) EQ 7 THEN vargrid = strupcase(grid)
   if vargrid eq 'W' then wdepth = 1
   if keyword_set(tout) then begin
@@ -356,11 +360,13 @@
   ENDELSE
 ; for the vertical sections with partial steps
+; from domzgh: hdept(ji,jj) = gdept(ji,jj,ik  )
+;              hdepw(ji, jj) = gdepw(ji, jj, ik+1)
   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
+        good = where(bottom NE 0 AND bottom NE nz-1+keyword_set(wdepth))
         IF good[0] NE -1 THEN BEGIN
+          bottom = lindgen(nx)+(bottom-(1l+keyword_set(wdepth)))*nx
           bottom = bottom[good]
           IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
@@ -374,7 +380,7 @@
       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
+        good = where(bottom NE 0 AND bottom NE nz-1+keyword_set(wdepth))
         IF good[0] NE -1 THEN BEGIN
+          bottom = lindgen(ny)+(bottom-(1l+keyword_set(wdepth)))*ny
           bottom = bottom[good]
           IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
Index: /trunk/SRC/ToBeReviewed/PLOTS/DESSINE/pltz.pro
===================================================================
--- /trunk/SRC/ToBeReviewed/PLOTS/DESSINE/pltz.pro	(revision 387)
+++ /trunk/SRC/ToBeReviewed/PLOTS/DESSINE/pltz.pro	(revision 388)
@@ -93,4 +93,7 @@
 ; For these ones traced by a trait, use CONTLABEL.
 ;
+; @keyword MASKFILL
+; set to 0 to avoid to fill the masked values
+;
 ; @keyword MAXIN
 ; to specify the max value we want to plot with a keyword instead of with the
@@ -109,4 +112,10 @@
 ; @keyword NOFILL
 ; To make just isolines with no filling
+;
+; @keyword NOMASK
+; activate to supress the use of the land-sea mask when plotting the field
+; 
+; @keyword NO_PARTIALSTEP
+; activate to supress the use of partial steps (force key_partialstep = 0 within pltz)
 ;
 ; @keyword NOTRI
@@ -116,4 +125,8 @@
 ; the keyword CELL_FILL=2.
 ; Comment: if the field contain points !values.f_nan, then we even do a triangulation.
+;
+; @keyword MASKDTA
+; 2d array defining the mask that must be applied to the data instead
+; of the default mask which is defined according to the grid (t/u/v/f mask)
 ;
 ; @keyword OVERPLOT
@@ -172,8 +185,9 @@
         , STYLE=style, CONTMAX=contmax, SIN=sin, TYPEIN=typein $
         , CONTMIN=contmin, NLEVEL=nlevel, CONTNLEVEL=contnlevel $
-        , NOFILL=nofill, NOTRI=notri, USETRI=usetri, FILLXDIR=fillxdir $
+        , NOFILL=nofill, NOMASK = nomask, NO_PARTIALSTEP = no_partialstep, NOTRI=notri $
+        , USETRI = usetri, FILLXDIR = fillxdir $
         , ZOOM=zoom, XZ=xz, YZ=yz, MININ=minin, MAXIN=maxin $
-        , STRICTFILL=strictfill, OVERPLOT=overplot, MASKFILL=maskfill $
-        , WDEPTH=wdepth, REALSECTION=realsection $
+        , STRICTFILL=strictfill, OVERPLOT=overplot $
+        , WDEPTH=wdepth, REALSECTION=realsection, MASKFILL=maskfill, MASKDTA = maskdta $
         , _EXTRA=ex
 ;
@@ -197,4 +211,5 @@
 ;
   if n_elements(contour) ne 4 THEN saveboxparam, 'boxparam4pltz.dat'
+  key_partialstep = keyword_set(key_partialstep) * (1b - keyword_set(no_partialstep)) 
 ;--------------------------------------------------------------
 ;  Reading of the field.
@@ -262,4 +277,7 @@
 ; to draw from bottom to top (avoid using cell_fill)
   z2d = reverse(z2d, 2)
+  szmsk = size(mask, /dimensions)
+  IF keyword_set(nomask) THEN mask = replicate(1b, szmsk[0], szmsk[1]) 
+  IF keyword_set(maskdta) THEN mask = byte(fitintobox(maskdta))
 ;-----------------------------------------------------------------------------
 ; Determination of the mi:min and of the ma:max of z2d in the same way
@@ -355,6 +373,8 @@
     if keyword_set(nan) then trinan = triangule(masknan, /basic, coinmonte = coinmontenan, coindescend = coindescendnan)
     maskorg = mask
-    decoupeterre, mask, glammsk, gphimsk, gdepmsk, type = type, WDEPTH = wdepth, REALSECTION = realsection
+    decoupeterre, mask, glammsk, gphimsk, gdepmsk, type = type, WDEPTH = wdepth, REALSECTION = realsection, MASKDTA = maskdta
     axis4pltz, type, mask, glammsk, gphimsk, gdepmsk, XXAXIS = xmask, ZZAXIS = zmask, SIN = sin, ZRATIO = zratio, ZOOM = zoom, PROFMAX = profmax, PROFMIN = profmin, _extra = ex
+    szmsk = size(mask, /dimensions)
+    IF keyword_set(nomask) THEN mask = replicate(1b, szmsk[0], szmsk[1]) 
   ENDIF ELSE BEGIN
     xmask = xxaxis
Index: /trunk/SRC/Utilities/fitintobox.pro
===================================================================
--- /trunk/SRC/Utilities/fitintobox.pro	(revision 387)
+++ /trunk/SRC/Utilities/fitintobox.pro	(revision 388)
@@ -89,6 +89,9 @@
 ; @param lastz {in}{optional}{default=define by <pro>grille</pro>}
 ;
-; @keyword WDEPTH
-; To specify that we are at W level
+; @keyword GRID {default=vargrid defined in <pro>cm_4data</pro>} 
+; To specify on which grid are located the data, see <pro>grille</pro>
+;
+; @keyword WDEPTH {default=0}
+; To specify that we are at W level, see <pro>grille</pro>
 ;
 ; @returns
@@ -120,5 +123,5 @@
 ;-
 FUNCTION fitintobox, field, nx, ny, nz, firstx, firsty $
-                   , firstz, lastx, lasty, lastz, WDEPTH=wdepth
+                     , firstz, lastx, lasty, lastz, GRID = GRID, WDEPTH = wdepth
 ;
   compile_opt idl2, strictarrsubs
@@ -132,5 +135,5 @@
   arr = litchamp(field)
   IF n_params() EQ 1 THEN grille, -1, -1, -1, -1, nx, ny, nz $
-    , firstx, firsty, firstz, lastx, lasty, lastz, WDEPTH = wdepth
+    , firstx, firsty, firstz, lastx, lasty, lastz, GRID = GRID, WDEPTH = wdepth
 ;--------------------------------------------------------------
 ;--------------------------------------------------------------
