Index: trunk/SRC/Computation/curl.pro
===================================================================
--- trunk/SRC/Computation/curl.pro	(revision 167)
+++ trunk/SRC/Computation/curl.pro	(revision 167)
@@ -0,0 +1,310 @@
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+;+
+;
+; @file_comments
+; Calculate the vertical component of the curl of a field of horizontal vectors
+;
+; @categories
+; Calculation
+; 
+; @param UU
+; Matrix representing the zonal coordinates (U point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt) or a structure readable by litchamp and containing
+; a 2D (xy), 3D (xyz or yt) array (4d case is not coded yet).
+; note that the dimension of the arry must suit the domain dimension.
+;
+; @param VV 
+; Matrix representing the meridional coordinates (V point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt) or a structure readable by litchamp and containing
+; a 2D (xy), 3D (xyz or yt) array (4d case is not coded yet).
+; note that the dimension of the arry must suit the domain dimension.
+;
+; @keyword DIREC {type=scalar string}
+; Use if you want to call moyenne or grossemoyenne after the div computation
+; with a mean done in the DIREC direction
+;
+; @returns RES
+; the vertical component of the curl of the input data (with the same size)
+;
+; @uses
+; cm_4cal, cm_4data, cm_4mmesh 
+;
+; @restrictions
+;
+; - Works only for Arakawa C-grid. 
+; - UU must be on U grid, VV must be on V grid
+; - 4d case is not coded yet
+; - the common variable jpt is used to differ xyz (jpt=1) and xyt (jpt\=1) cases.
+; - U and V arrays are cut in the same geographic domain. Because of the shift between 
+;   T, U, V and F grids, it is possible that these two arrays do not have the same 
+;   size and refer to different indexes. In this case, arrays are re-cut on 
+;   common indexes. To avoid these re-cuts, use the keyword /memeindice in domdef.pro
+; - When computing the divergence, we update, vargrid, varname, varunits and the
+;   grid position parameters (firstxf, lastxf, nxf, firstyf, lastyf, nyf).
+; - points that cannot be computed (domain bondaries, coastline) are set to NaN
+;
+; @examples
+; IDL> @tst_initorca2
+; IDL> plt, curl(dist(jpi,jpj), dist(jpi,jpj))
+;
+; @history 
+; Guillaume Roullet (grlod\@ipsl.jussieu.fr)
+; Sebastien Masson (smasson\@lodyc.jussieu.fr)
+; adaptation to work with a reduce domain
+; 21/5/1999: missing values at !values.f_nan
+; 
+; @version
+; $Id$
+;
+; @todo 
+; code the 4d case
+;-
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+FUNCTION curl, uu, vv, DIREC = DIREC
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4cal                        ; for jpt
+@cm_4data                       ; for varname, vargrid, vardate, varunit, valmask
+@cm_4mesh 
+;
+  tempsun = systime(1)          ; To key_performance
+;
+  IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
+     return, report(['This version of curl is based on Arakawa C-grid.' $
+                     , 'U and V grids must therefore be defined'])
+;
+  u = litchamp(uu)
+  v = litchamp(vv)
+
+  szu = size(u)
+  szv = size(v)
+  
+  if szu[0] NE szv[0] then return, report('U and V input data must have the same number of dimensions')
+
+;------------------------------------------------------------
+; We find common points between U and V
+;------------------------------------------------------------
+  indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
+  indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
+  indicex = inter(indicexu, indicexv)
+  indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
+  indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
+  indicey = inter(indiceyu, indiceyv)
+  nx = n_elements(indicex) 
+  ny = n_elements(indicey)
+  indice2d = lindgen(jpi, jpj)
+  indice2d = indice2d[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1]
+;----------------------------------------------------------------------------
+  vargrid = 'F'
+  varname = 'vorticity'
+  varunits = 's-1'
+  if n_elements(valmask) EQ 0 THEN valmask = 1e20
+  firstxt = indicex[0] & lastxt = indicex[0]+nx-1 & nxt = nx 
+  firstyt = indicey[0] & lastyt = indicey[0]+ny-1 & nyt = ny
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+  case 1 of
+;----------------------------------------------------------------------------
+;xyz
+;----------------------------------------------------------------------------
+    szu[0] EQ 3 AND jpt EQ 1:BEGIN 
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          case 1 of
+            nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
+            nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
+            nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
+            nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
+            ELSE :
+          endcase
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+          v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+        END
+        ELSE:return,  -1
+      endcase
+;------------------------------------------------------------
+; curl computation
+;------------------------------------------------------------
+      coefu = ((e1u[indice2d])[*]#replicate(1., nzt)) $
+              * (umask())[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
+      landu = where(coefu EQ 0)
+      if landu[0] NE -1 then coefu[temporary(landu)] = !values.f_nan
+      
+      coefv = ((e2v[indice2d])[*]#replicate(1., nzt)) $
+              *(vmask())[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
+      landv = where(coefv EQ 0)
+      if landv[0] NE -1 then coefv[temporary(landv)] = !values.f_nan
+
+      tabf = (fmask())[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt] $
+             / ((e1f[indice2d]*e2f[indice2d])[*]#replicate(1., nzt))
+      landf =  where(tabf EQ 0)
+;
+      zu = temporary(u) * temporary(coefu)
+      zv = temporary(v) * temporary(coefv)
+
+      psi = (shift(zv, -1, 0, 0)-zv) + (zu-shift(zu, 0, -1, 0))
+      psi = temporary(tabf) * temporary(psi)
+;------------------------------------------------------------
+; Edging put at !values.f_nan
+;------------------------------------------------------------
+      if NOT keyword_set(key_periodic)  OR nx NE jpi then begin
+        psi[0, *, *] = !values.f_nan
+        psi[nx-1, *, *] = !values.f_nan
+      endif
+      psi[*, 0, *] = !values.f_nan
+      psi[*, ny-1, *] = !values.f_nan
+;
+      if landf[0] NE -1 then psi[temporary(landf)] = valmask
+      if keyword_set(direc) then psi = moyenne(psi, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xyt
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 3 AND jpt GT 1:BEGIN 
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          if nxu NE nx then $
+             if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *] 
+          IF nxv NE nx THEN $
+             if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
+          IF nyu NE ny THEN $
+             if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *] 
+          IF nyv NE ny THEN $
+             if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+          v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+        END
+        ELSE:BEGIN 
+          print, 'problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs'
+          return, -1
+        end
+      endcase
+;----------------------------------------------------------------------------
+; curl computation
+;----------------------------------------------------------------------------
+      coefu = e1u[indice2d]*(umask())[indice2d+jpi*jpj*firstzt]
+      landu = where(coefu EQ 0)
+      if landu[0] NE -1 then coefu[temporary(landu)] = !values.f_nan
+      coefu = temporary(coefu[*])#replicate(1., jpt)
+;
+      coefv = e2v[indice2d]*(vmask())[indice2d+jpi*jpj*firstzt]
+      landv = where(coefv EQ 0)
+      if landv[0] NE -1 then coefv[temporary(landv)] = !values.f_nan
+      coefv = temporary(coefv[*])#replicate(1., jpt)
+;
+      tabf = (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
+      tabf = temporary(tabf[*])#replicate(1., jpt)
+      landf = where(tabf EQ 0)
+;
+      zu = temporary(u) * temporary(coefu)
+      zv = temporary(v) * temporary(coefv)
+;
+      psi = (shift(zv, -1, 0, 0)-zv) + (zu-shift(zu, 0, -1, 0))
+      psi = temporary(tabf) * temporary(psi)
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      if NOT keyword_set(key_periodic) OR nx NE jpi then begin
+        psi[0, *, *] = !values.f_nan
+        psi[nx-1, *, *] = !values.f_nan
+      endif
+      psi[*, 0, *] = !values.f_nan
+      psi[*, ny-1, *] = !values.f_nan
+      if landf[0] NE -1 then psi[temporary(landf)] = valmask
+      if keyword_set(direc) then psi = grossemoyenne(psi, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xyzt
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 4:BEGIN 
+      return, report('Case not coded contact saxo team or make a do loop!')
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xy
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 2:BEGIN 
+;------------------------------------------------------------
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          if nxu NE nx then $
+             if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *] 
+          IF nxv NE nx THEN $
+             if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
+          IF nyu NE ny THEN $
+             if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny] 
+          IF nyv NE ny THEN $
+             if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indice2d]
+          v = v[indice2d]
+        END
+        ELSE:return,  -1
+      endcase
+;------------------------------------------------------------
+; curl computation
+;------------------------------------------------------------
+      coefu = e1u[indice2d]*(umask())[indice2d+jpi*jpj*firstzt]
+      landu = where(coefu EQ 0)
+      if landu[0] NE -1 then coefu[temporary(landu)] = !values.f_nan
+      coefv = e2v[indice2d]*(vmask())[indice2d+jpi*jpj*firstzt]
+      landv = where(coefv EQ 0)
+      if landv[0] NE -1 then coefv[temporary(landv)] = !values.f_nan
+      tabf = (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
+      landf =  where(tabf EQ 0)
+;
+      zu = temporary(u) * temporary(coefu)
+      zv = temporary(v) * temporary(coefv)
+
+      psi = (shift(zv, -1, 0)-zv) + (zu-shift(zu, 0, -1))
+      psi = temporary(tabf) * temporary(psi)
+;------------------------------------------------------------
+; Edging put at !values.f_nan
+;------------------------------------------------------------
+      if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
+        psi[0, *] = !values.f_nan
+        psi[nx-1, *] = !values.f_nan
+      endif
+      psi[*, 0] = !values.f_nan
+      psi[*, ny-1] = !values.f_nan
+;
+      if landf[0] NE -1 then psi[temporary(landf)] = valmask
+      if keyword_set(direc) then psi = moyenne(psi, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    ELSE:return, report('U and V input arrays must have 2, 3 or 4 dimensions')
+  ENDCASE
+;------------------------------------------------------------
+  if keyword_set(key_performance) THEN print, 'temps curl', systime(1)-tempsun 
+
+  return, psi
+end
Index: trunk/SRC/Computation/div.pro
===================================================================
--- trunk/SRC/Computation/div.pro	(revision 167)
+++ trunk/SRC/Computation/div.pro	(revision 167)
@@ -0,0 +1,298 @@
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+;+
+;
+; @file_comments
+; compute the horizontal divergence of a vectors field
+;
+; @categories
+; Calculation
+;
+; @param UU 
+; Matrix representing the zonal coordinates (U point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt) or a structure readable by litchamp and containing
+; a 2D (xy), 3D (xyz or yt) array (4d case is not coded yet).
+; note that the dimension of the arry must suit the domain dimension.
+;
+; @param VV 
+; Matrix representing the meridional coordinates (V point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt) or a structure readable by litchamp and containing
+; a 2D (xy), 3D (xyz or yt) array (4d case is not coded yet).
+; note that the dimension of the arry must suit the domain dimension.
+;
+; @keyword DIREC {type=scalar string}
+; Use if you want to call moyenne or grossemoyenne after the div computation
+; (stupid ?) with a mean done in the DIREC direction
+;
+; @returns RES
+; the divergence of the input data (with the same size)
+;
+; @uses
+; cm_4cal, cm_4data, cm_4mmesh 
+;
+; @restrictions
+;
+; - Works only for Arakawa C-grid. 
+; - UU must be on U grid, VV must be on V grid
+; - 4d case is not coded yet
+; - the common variable jpt is used to differ xyz (jpt=1) and xyt (jpt\=1) cases.
+; - U and V arrays are cut in the same geographic domain. Because of the shift between 
+;   T, U, V and F grids, it is possible that these two arrays do not have the same 
+;   size and refer to different indexes. In this case, arrays are re-cut on 
+;   common indexes. To avoid these re-cuts, use the keyword /memeindice in domdef.pro
+; - When computing the divergence, we update, vargrid, varname, varunits and the
+;   grid position parameters (firstxt, lastxt, nxt, firstyt, lastyt, nyt).
+; - points that cannot be computed (domain bondaries, coastline) are set to NaN
+;
+; @examples
+; IDL> @tst_initorca2
+; IDL> plt, div(dist(jpi,jpj), dist(jpi,jpj))
+;
+; @history
+; Guillaume Roullet (grlod\@ipsl.jussieu.fr): creation; spring 1998
+; Sebastien Masson (smasson\@lodyc.jussieu.fr)
+; adaptation to work with a reduce domain; 12/1/2000
+;
+; @version
+; $Id$
+;
+; @todo 
+; code the 4d case
+;-
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+FUNCTION div, uu, vv, DIREC = DIREC
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4cal                        ; for jpt
+@cm_4data                       ; for varname, vargrid, vardate, varunit, valmask
+@cm_4mesh 
+;
+  tempsun = systime(1)          ; For key_performance
+;
+  IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
+     return, report(['This version of div is based on Arakawa C-grid.' $
+                     , 'U and V grids must therefore be defined'])
+;
+  u = litchamp(uu)
+  v = litchamp(vv)
+;
+  szu = size(u)
+  szv = size(v)
+  
+  if szu[0] NE szv[0] then return, report('U and V input data must have the same number of dimensions')
+
+;------------------------------------------------------------
+; We find common points between U and V
+;------------------------------------------------------------
+  indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
+  indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
+  indicex = inter(indicexu, indicexv)
+  indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
+  indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
+  indicey = inter(indiceyu, indiceyv)
+  nx = n_elements(indicex) 
+  ny = n_elements(indicey)
+  indice2d = lindgen(jpi, jpj)
+  indice2d = indice2d[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1]
+;----------------------------------------------------------------------------
+  vargrid = 'T'
+  varname = 'div'
+  varunits = '1.e6*s-1'
+  if n_elements(valmask) EQ 0 THEN valmask = 1.e20
+  firstxt = indicex[0] & lastxt = indicex[0]+nx-1 & nxt = nx 
+  firstyt = indicey[0] & lastyt = indicey[0]+ny-1 & nyt = ny
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+  case 1 of
+;----------------------------------------------------------------------------
+;xyz
+;----------------------------------------------------------------------------
+    szu[0] EQ 3 AND jpt EQ 1:BEGIN 
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          case 1 of
+            nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
+            nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
+            nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
+            nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
+            ELSE :
+          endcase
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+          v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+        END
+        ELSE:return, -1
+      endcase
+;------------------------------------------------------------
+; divergence computation
+;------------------------------------------------------------
+      zu = (e2u[indice2d])[*]#replicate(1., nzt)
+      landu = where((umask())[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt] EQ 0)
+      if landu[0] NE -1 then zu[temporary(landu)] = !values.f_nan
+      zu = temporary(u) * temporary(zu)
+;
+      zv = (e1v[indice2d])[*]#replicate(1., nzt)
+      landv = where((vmask())[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt] EQ 0)
+      if landv[0] NE -1 then zv[temporary(landv)] = !values.f_nan
+      zv = temporary(v) * temporary(zv)
+;
+      zdiv = (e1t[indice2d]*e2t[indice2d])[*]#replicate(1.e6, nzt)
+      zdiv = ( zu - shift(zu, 1, 0, 0) + zv - shift(zv, 0, 1, 0) ) * temporary(zdiv)
+;------------------------------------------------------------
+; Edging put at !values.f_nan
+;------------------------------------------------------------
+      if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
+        zdiv[0, *, *] = !values.f_nan
+        zdiv[nx-1, *, *] = !values.f_nan
+      endif
+      zdiv[*, 0, *] = !values.f_nan
+      zdiv[*, ny-1, *] = !values.f_nan
+;
+      land = where(tmask[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt] EQ 0)
+      if land[0] NE -1 then zdiv[temporary(land)] = valmask
+      if keyword_set(direc) then  zdiv = moyenne(zdiv, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xyt
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 3 AND jpt GT 1:BEGIN 
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          case 1 of
+            nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
+            nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
+            nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
+            nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
+            ELSE :
+          endcase
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+          v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
+        END
+        ELSE:return, -1
+      endcase
+;------------------------------------------------------------
+; divergence computation
+;------------------------------------------------------------
+      zu = e2u[indice2d]
+      landu = where((umask())[indice2d+jpi*jpj*firstzt] EQ 0)
+      if landu[0] NE -1 then zu[temporary(landu)] = !values.f_nan
+      zu = (temporary(zu))[*]#replicate(1., jpt)
+      zu = temporary(u) * temporary(zu)
+;
+      zv = e1v[indice2d]
+      landv = where((vmask())[indice2d+jpi*jpj*firstzt] EQ 0)
+      if landv[0] NE -1 then zv[temporary(landv)] = !values.f_nan
+      zv = (temporary(zv))[*]#replicate(1., jpt)
+      zv = temporary(v) * temporary(zv)
+;
+      zdiv = (e1t[indice2d]*e2t[indice2d])[*]#replicate(1.e6, jpt)
+      zdiv = ( zu - shift(zu, 1, 0, 0) + zv - shift(zv, 0, 1, 0) ) * temporary(zdiv)
+;------------------------------------------------------------
+; Edging put at !values.f_nan
+;------------------------------------------------------------
+      if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
+        zdiv[0, *, *] = !values.f_nan
+        zdiv[nx-1, *, *] = !values.f_nan
+      endif
+      zdiv[*, 0, *] = !values.f_nan
+      zdiv[*, ny-1, *] = !values.f_nan
+;
+      land = where(tmask[indice2d+jpi*jpj*firstzt] EQ 0, cnt)
+      if land[0] NE -1 then BEGIN 
+        land = (temporary(land))#replicate(1L, jpt) + replicate(1L, cnt)#(nx*ny*lindgen(jpt))
+        zdiv[temporary(land)] = valmask
+      ENDIF
+      if keyword_set(direc) then  zdiv = grossemoyenne(zdiv, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xyzt
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 4:BEGIN 
+      return, report('Case not coded contact saxo team or make a do loop!')
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xy
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    szu[0] EQ 2:BEGIN 
+;------------------------------------------------------------
+; extraction of U and V on the appropriated domain
+;------------------------------------------------------------
+      case 1 of
+        szu[1] EQ nxu AND szu[2] EQ nyu AND $
+           szv[1] EQ nxv AND szv[2] EQ nyv:BEGIN 
+          case 1 of
+            nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *]
+            nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
+            nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny]
+            nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
+            ELSE :
+          endcase
+        END
+        szu[1] EQ jpi AND szu[2] EQ jpj AND $
+           szv[1] EQ jpi AND szv[2] EQ jpj:BEGIN 
+          u = u[indice2d]
+          v = v[indice2d]
+        END
+        ELSE:return, -1
+      endcase
+;------------------------------------------------------------
+; divergence computation
+;------------------------------------------------------------
+      zu = e2u[indice2d]
+      landu = where((umask())[indice2d+jpi*jpj*firstzt] EQ 0)
+      if landu[0] NE -1 then zu[temporary(landu)] = !values.f_nan
+      zu = temporary(u) * temporary(zu)
+
+      zv = e1v[indice2d]
+      landv = where((vmask())[indice2d+jpi*jpj*firstzt] EQ 0)
+      if landv[0] NE -1 then zv[temporary(landv)] = !values.f_nan
+      zv = temporary(v) * temporary(zv)
+
+      zdiv = 1.e6 / (e1t[indice2d]*e2t[indice2d])
+      zdiv = ( zu - shift(zu, 1, 0) + zv - shift(zv, 0, 1) ) * temporary(zdiv)
+;------------------------------------------------------------
+; Edging put at !values.f_nan
+;------------------------------------------------------------
+      if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
+        zdiv[0, *] = !values.f_nan
+        zdiv[nx-1, *] = !values.f_nan
+      endif
+      zdiv[*, 0] = !values.f_nan
+      zdiv[*, ny-1] = !values.f_nan
+;
+      land =  where(tmask[indice2d+jpi*jpj*firstzt] EQ 0)
+      if land[0] NE -1 then zdiv[temporary(land)] = valmask
+      if keyword_set(direc) then zdiv = moyenne(zdiv, direc, /nan)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    ELSE:return, report('U and V input arrays must have 2, 3 or 4 dimensions')
+  ENDCASE   
+;------------------------------------------------------------
+  if keyword_set(key_performance) THEN print, 'temps div', systime(1)-tempsun 
+  
+  return, zdiv
+end
Index: trunk/SRC/Computation/grad.pro
===================================================================
--- trunk/SRC/Computation/grad.pro	(revision 167)
+++ trunk/SRC/Computation/grad.pro	(revision 167)
@@ -0,0 +1,327 @@
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+;+
+; @file_comments
+; compute the gradient of a variable
+;
+; @categories
+; Calculation
+;
+; @param FIELD
+; The field for which we want to compute the gradient.  A 2D (xy),
+; 3D (xyz or yt) or 4D (xyzt) array or a structure readable by litchamp
+; and containing a 2D (xy), 3D (xyz or yt) or 4D (xyzt) array.
+; note that the dimension of the arry must suit the domain dimension.
+;
+; @param DIREC {type=scalar string}
+; the gradient direction: 'x', 'y', 'z'
+;
+; @returns RES {type=2D, 3D or 4D array}
+; the gradient of the input data (with the same size)
+;
+; @uses
+; cm_4cal, cm_4data, cm_4mmesh 
+;
+; @restrictions
+; - Works only for Arakawa C-grid. 
+; - When computing the gradient, the result is not on the same grid point
+;   than the input data. In consequence, we update, vargrid and the grid position
+;   parameters (firstx[tuvf], lastx[tuvf], nx[tuvf], firsty[tuvf], lasty[tuvf],
+;   ny[tuvf], firstz[tw], lastz[tw], nz[tw]).
+; - points that cannot be computed (domain bondaries, coastline) are set to NaN
+; - the common variable jpt is used to differ xyz (jpt=1) and xyt (jpt\=1) cases.
+;
+; @examples
+; IDL> @tst_initorca2
+; IDL> plt, grad({arr:gphit,g:'T'}, 'x')
+; IDL> plt, grad({arr:gphit,g:'T'}, 'y')
+;
+; @history
+; Sebastien Masson (smasson\@lodyc.jussieu.fr)
+;
+; @version
+; $Id$
+;
+;-
+;------------------------------------------------------------
+;------------------------------------------------------------
+;------------------------------------------------------------
+FUNCTION grad, field, direc
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4cal   ; for jpt
+@cm_4data  ; for varname, vargrid, vardate, varunit, valmask
+@cm_4mesh 
+;------------------------------------------------------------
+;
+  IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
+     return, report(['This version of grad is based on Arakawa C-grid.' $
+                     , 'U and V grids must therefore be defined'])
+;
+  res = litchamp(field)
+  szres = size(res)
+  grille, mask, glam, gphi, gdep, nx, ny, nz $
+          , firstx, firsty, firstz, lastx, lasty, lastz     
+;
+  if n_elements(valmask) EQ 0 then valmask = 1.e20
+  varname = 'grad of '+varname
+  varunit = varunit+'/m'
+  case strupcase(vargrid) of
+    'T':BEGIN
+      case direc of
+        'x':BEGIN 
+          divi = e1u[firstx:lastx, firsty:lasty]
+          newmask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'U'
+          firstxu = firstxt & lastxu = lastxt & nxu = nxt 
+          firstyu = firstyt & lastyu = lastyt & nyu = nyt
+        END
+        'y':BEGIN
+          divi = e2v[firstx:lastx, firsty:lasty]
+          newmask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'V'
+          firstxv = firstxt & lastxv = lastxt & nxv = nxt 
+          firstyv = firstyt & lastyv = lastyt & nyv = nyt
+        END
+        'z':BEGIN
+          divi = e3w[firstz:lastz]
+          newmask = mask
+          vargrid = 'W'
+          firstzw = firstzt & lastzw = lastzt & nzw = nzt 
+        END
+        ELSE:return, report('Bad definition of direction argument')
+      ENDCASE
+    END
+    'W':BEGIN
+      case direc of
+        'x':BEGIN 
+          divi = e1u[firstx:lastx, firsty:lasty]
+          newmask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'U'
+          firstxu = firstxt & lastxu = lastxt & nxu = nxt 
+          firstyu = firstyt & lastyu = lastyt & nyu = nyt
+        END
+        'y':BEGIN
+          divi = e2v[firstx:lastx, firsty:lasty]
+          newmask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'V'
+          firstxv = firstxt & lastxv = lastxt & nxv = nxt 
+          firstyv = firstyt & lastyv = lastyt & nyv = nyt
+        END
+        'z':BEGIN
+          divi = e3t[firstz:lastz]
+          newmask = mask
+          vargrid = 'T'
+          firstzt = firstzw & lastzt = lastzw & nzt = nzw 
+        END
+        ELSE:return, report('Bad definition of direction argument')
+      endcase
+    END
+    'U':BEGIN
+      case direc of
+        'x':BEGIN
+          divi = (shift(e1t, -1, 0))[firstx:lastx, firsty:lasty]
+          newmask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'T'
+          firstxt = firstxu & lastxt = lastxu & nxt = nxu 
+          firstyt = firstyu & lastyt = lastyu & nyt = nyu
+        END
+        'y':BEGIN
+          divi = e2f[firstx:lastx, firsty:lasty]
+          newmask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'F'
+          firstxf = firstxu & lastxf = lastxu & nxf = nxu 
+          firstyf = firstyu & lastyf = lastyu & nyf = nyu
+        END
+        'z':BEGIN
+          divi = e3w[firstz:lastz]
+          newmask = mask
+          vargrid = 'W'
+          firstzw = firstzt & lastzw = lastzt & nzw = nzt 
+        END
+        ELSE:return, report('Bad definition of direction argument')
+      endcase
+    END
+    'V':BEGIN
+      case direc of
+        'x':BEGIN
+          divi = e1f[firstx:lastx, firsty:lasty]
+          newmask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'F'
+          firstxf = firstxv & lastxf = lastxv & nxf = nxv
+          firstyf = firstyv & lastyf = lastyv & nyf = nyv
+        END
+        'y':BEGIN
+          divi = (shift(e2t, 0, -1))[firstx:lastx, firsty:lasty]
+          newmask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
+          vargrid = 'T'
+          firstxt = firstxv & lastxt = lastxv & nxt = nxv
+          firstyt = firstyv & lastyt = lastyv & nyt = nyv
+        END
+        'z':BEGIN
+          divi = e3w[firstz:lastz]
+          newmask = mask
+          vargrid = 'W'
+          firstzw = firstzt & lastzw = lastzt & nzw = nzt 
+        END
+        ELSE:return, report('Bad definition of direction argument')
+      endcase
+    END
+    'F':BEGIN
+;          case direc of
+;             'x':divi = (shift(e1v, -1, 0))[firstx:lastx, firsty:lasty]
+;             'y':divi = (shift(e2u, 0, -1))[firstx:lastx, firsty:lasty]
+;             'z':divi = e3w[firstz:lastz]
+;             ELSE:return, report('Bad definition of direction argument')
+;          endcase
+      return, report('F grid: case not coded, please contact SAXO team...')
+    END
+    ELSE:return, report('Bad definition of vargrid')
+  ENDCASE
+  res = fitintobox(temporary(res))
+  IF n_elements(res) EQ 1 AND res[0] EQ -1 THEN return, res
+  case 1 of
+;----------------------------------------------------------------------------
+;xy
+;----------------------------------------------------------------------------
+    szres[0] EQ 2:BEGIN
+      land = where((temporary(mask))[*, *, firstz] EQ 0)
+      if land[0] NE -1 then res[temporary(land)] = !values.f_nan
+      case direc of
+        'x':BEGIN 
+          res = (shift(res, -1, 0)-res)/temporary(divi)
+          if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(temporary(res), 1, 0)
+        END
+        'y':BEGIN 
+          res = (shift(res, 0, -1)-res)/temporary(divi)
+          res[*, ny-1] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(temporary(res), 0, 1)           
+        END
+        ELSE:return,  report('Bad definition of direction argument for the type of array')
+      ENDCASE
+      land = where((temporary(newmask))[*, *, firstz] EQ 0)
+      if land[0] NE -1 then res[temporary(land)] = valmask
+    END
+;----------------------------------------------------------------------------
+;xyt
+;----------------------------------------------------------------------------
+    szres[0] EQ 3 AND jpt NE 1:BEGIN
+      land = where((temporary(mask))[*, *, firstz] EQ 0, cnt)
+      if land[0] NE -1 then BEGIN
+        land = (temporary(land))#replicate(1L, jpt) + replicate(1L, cnt)#(nx*ny*lindgen(jpt))
+        res[temporary(land)] = !values.f_nan
+      ENDIF
+      divi = (temporary(divi))[*]#replicate(1., jpt)
+      case direc of
+        'x':BEGIN 
+          res = (shift(res, -1, 0, 0)-res)/temporary(divi)
+          if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(temporary(res), 1, 0, 0)
+        END
+        'y':BEGIN 
+          res = (shift(res, 0, -1, 0)-res)/temporary(divi)
+          res[*, ny-1, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(temporary(res), 0, 1, 0)        
+        END
+        ELSE:return,  report('Bad definition of direction argument for the type of array')
+      ENDCASE
+      land = where((temporary(newmask))[*, *, firstz] EQ 0, cnt)
+      if land[0] NE -1 then BEGIN
+        land = (temporary(land))#replicate(1L, jpt) + replicate(1L, cnt)#(nx*ny*lindgen(jpt))
+        res[temporary(land)] = valmask
+      ENDIF
+    END
+;----------------------------------------------------------------------------
+;xyz
+;----------------------------------------------------------------------------
+    szres[0] EQ 3 AND jpt EQ 1:BEGIN
+      land = where(mask EQ 0)
+      if land[0] NE -1 then res[temporary(land)] = !values.f_nan
+      case direc OF
+        'x':BEGIN 
+          divi = (temporary(divi))[*]#replicate(1., nz)
+          res = (shift(res, -1, 0, 0)-res)/temporary(divi)
+          if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(temporary(res), 1, 0, 0)
+        END
+        'y':BEGIN 
+          divi = (temporary(divi))[*]#replicate(1., nz)
+          res = (shift(res, 0, -1, 0)-res)/temporary(divi)
+          res[*, ny-1, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(temporary(res), 0, 1, 0)        
+        END
+        'z':BEGIN
+          divi = replicate(1., nx*ny)#(temporary(divi))[*]
+          if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz, /overwrite)
+          if vargrid EQ 'W' THEN BEGIN 
+            res = (shift(res, 0, 0, 1)-res)/temporary(divi)
+            res[*, *, 0] = !values.f_nan
+          ENDIF ELSE BEGIN
+            res = (res-shift(res, 0, 0, -1))/temporary(divi)
+            res[*, *, nz-1] = !values.f_nan
+          ENDELSE
+        END
+      ENDCASE
+      land = where(temporary(newmask) EQ 0)
+      if land[0] NE -1 then res[temporary(land)] = valmask
+    END
+;----------------------------------------------------------------------------
+;xyzt
+;----------------------------------------------------------------------------
+    szres[0] EQ 4:BEGIN
+      land = where(temporary(mask) EQ 0, cnt)
+      if land[0] NE -1 then BEGIN
+        land = (temporary(land))#replicate(1L, jpt) + replicate(1L, cnt)#(nx*ny*nz*lindgen(jpt))
+        res[temporary(land)] = !values.f_nan
+      ENDIF
+      case direc OF
+        'x':BEGIN 
+          divi = (temporary(divi))[*]#replicate(1., nz*jpt)
+          res = (shift(res, -1, 0, 0, 0)-res)/temporary(divi)
+          if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(temporary(res), 1, 0, 0, 0)
+        END
+        'y':BEGIN 
+          divi = (temporary(divi))[*]#replicate(1., nz*jpt)
+          res = (shift(res, 0, -1, 0, 0)-res)/temporary(divi)
+          res[*, ny-1, *, *] = !values.f_nan
+          if vargrid EQ 'T' OR vargrid EQ 'U' then res = shift(temporary(res), 0, 1, 0, 0)        
+        END
+        'z':BEGIN
+          divi = replicate(1., nx*ny)#(temporary(divi))[*]
+          divi = (temporary(divi))[*]#replicate(1L, jpt)
+          if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz, jpt, /overwrite)
+          if vargrid EQ 'W' THEN BEGIN 
+            res = (shift(res, 0, 0, 1, 0)-res)/temporary(divi)
+            res[*, *, 0, *] = !values.f_nan
+          ENDIF ELSE BEGIN
+            res = (res-shift(res, 0, 0, -1, 0))/temporary(divi)
+            res[*, *, nz-1, *] = !values.f_nan
+          ENDELSE
+        END
+      ENDCASE
+      land = where(newmask EQ 0, cnt)
+      if land[0] NE -1 then BEGIN
+        land = (temporary(land))#replicate(1L, jpt) + replicate(1L, cnt)#(nx*ny*nz*lindgen(jpt))
+        res[temporary(land)] = valmask
+      ENDIF
+    END
+;------------------------------------------------------------
+;------------------------------------------------------------
+    ELSE:return, report('input array must have 2, 3 or 4 dimensions')
+  ENDCASE
+;------------------------------------------------------------
+
+
+;------------------------------------------------------------
+  return, res
+END
+
+
+
+
+
+
Index: trunk/SRC/Postscript/openps.pro
===================================================================
--- trunk/SRC/Postscript/openps.pro	(revision 166)
+++ trunk/SRC/Postscript/openps.pro	(revision 167)
@@ -128,5 +128,5 @@
            , LANDSCAPE = 1 - key_portrait, PORTRAIT = key_portrait $
            , xsize = xs, ysize = ys, xoffset = xoff, yoffset = yoff $
-           , bits_per_pixel = 8, _extra = ex
+           , bits_per_pixel = 8, language_level = 2, _extra = ex
 ; to make smaller postcripts
    IF NOT (keyword_set(keeppfont) OR keyword_set(keep_pfont)) $
Index: trunk/SRC/ToBeReviewed/CALCULS/curl.pro
===================================================================
--- trunk/SRC/ToBeReviewed/CALCULS/curl.pro	(revision 166)
+++ 	(revision )
@@ -1,318 +1,0 @@
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-;+
-;
-; @file_comments
-; Calculate the vertical component of the curl of a field of horizontal vectors
-;
-; @categories
-; Calculation
-; 
-; @param UU
-; Matrix representing coordinates of a field of vectors
-;
-; @param VV 
-; Matrix representing coordinates of a field of vectors
-;
-; @returns RES
-; A 2d matrix
-;
-; @uses
-; common.pro
-;
-; @restrictions
-; U and V matrices can be 2 or 4d.
-; Beware, to discern different configuration of U and V (xy, xyz, xyt, xyzt), 
-; we look at the variable of the common 
-;        -time which contain the calendar in IDL julian days to which U and 
-; V refered to, in the same way as the variable 
-;        -jpt which is the number of time's step to consider in time.
-; U and V arrays ae cut in the same geographic domain. Because of the gap of 
-; T, U, V and F grids, it is possible that these two arrays has not the same 
-; size and refered to different indexes. In this case, arrays are re-cut on 
-; common indexes and the domain is redefined to match with these common 
-; indexes. To avoid these re-cuts, use the keyword /memeindice in domdef.pro
-;
-;
-; Points on the drawing edge are at !values.f_nan 
-;
-; @history 
-; Guillaume Roullet (grlod\@ipsl.jussieu.fr)
-; 
-; Sebastien Masson (smasson\@lodyc.jussieu.fr)
-; adaptation to work with a reduce domain
-; 21/5/1999: missing values at !values.f_nan
-; 
-; @version
-; $Id$
-;
-;-
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-FUNCTION curl, uu, vv
-;
-  compile_opt idl2, strictarrsubs
-;
-@common
-   tempsun = systime(1)         ; To key_performance
-;
-   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
-     return, report(['This version of curl is based on Arakawa C-grid.' $
-                     , 'U and V grids must therefore be defined'])
-;
-   u = litchamp(uu)
-   v = litchamp(vv)
-
-   date1 = time[0]
-   if n_elements(jpt) EQ 0 then date2 = date1 ELSE date2 = time[jpt-1]
-   if (size(u))[0] NE (size(v))[0] then return,  -1
-
-;------------------------------------------------------------
-; We find common points between U and V
-;------------------------------------------------------------
-   indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
-   indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
-   indicex = inter(indicexu, indicexv)
-   indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
-   indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
-   indicey = inter(indiceyu, indiceyv)
-   nx = n_elements(indicex) 
-   ny = n_elements(indicey)
-   case 1 of
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyz
-;----------------------------------------------------------------------------
-      (size(u))[0] EQ 3 AND date1 EQ date2 :BEGIN 
-         indice2d = lindgen(jpi, jpj)
-         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         case 1 of
-            (size(u))[0] NE 3 OR (size(v))[0] NE 3: return,  -1
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-            case 1 of
-               nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
-               nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
-               nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
-               nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
-               ELSE :
-            endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-            END
-            ELSE:return,  -1
-         endcase
-;------------------------------------------------------------
-; calculation of the curl
-;------------------------------------------------------------
-         coefu = (e1u[indice2d])[*]#replicate(1, nzt)
-         coefu = reform(coefu, nx, ny, nzt, /over)
-         coefu = coefu*(umask())[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-         terreu = where(coefu EQ 0)
-         if terreu[0] NE -1 then coefu[temporary(terreu)] = !values.f_nan
-         
-         coefv = (e2v[indice2d])[*]#replicate(1, nzt)
-         coefv = reform(coefv, nx, ny, nzt, /over)
-         coefv = coefv*(vmask())[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-         terrev = where(coefv EQ 0)
-         if terrev[0] NE -1 then coefv[temporary(terrev)] = !values.f_nan
-
-         tabf = (fmask())[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-         div = (e1f[indice2d]*e2f[indice2d])[*]#replicate(1, nzt)
-         div = reform(div, nx, ny, nzt, /over)
-         tabf = tabf/div
-;
-         zu = u*temporary(coefu)
-         zv = v*temporary(coefv)
-
-         psi =  (shift(zv, -1, 0, 0)-zv) + (zu-shift(zu, 0, -1, 0))
-         psi = tabf*psi
-;------------------------------------------------------------
-; Edging put at !values.f_nan
-;------------------------------------------------------------
-         if NOT keyword_set(key_periodic)  OR nx NE jpi then begin
-            psi[0, *, *] = !values.f_nan
-            psi[nx-1, *, *] = !values.f_nan
-         endif
-         psi[*, 0, *] = !values.f_nan
-         psi[*, ny-1, *] = !values.f_nan
-;
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         terref =  where(tabf EQ 0)
-         if terref[0] NE -1 then psi[temporary(terref)] = valmask
-;------------------------------------------------------------
-; For the graphic drawing
-;------------------------------------------------------------
-         domdef, (glagmt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t', 'f']
-         if keyword_set(direc) then psi = moyenne(psi,direc,/nan)
-
-;----------------------------------------------------------------------------
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyt
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      date1 NE date2 AND (size(u))[0] EQ 3 :BEGIN 
-         indice2d = lindgen(jpi, jpj)
-         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         case 1 of
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-               if nxu NE nx then $
-                if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *] 
-               IF nxv NE nx THEN $
-                if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
-               IF nyu NE ny THEN $
-                if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *] 
-               IF nyv NE ny THEN $
-                if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-            END
-            ELSE:BEGIN 
-               print, 'problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs'
-               return, -1
-            end
-         endcase
-;----------------------------------------------------------------------------
-; Calculation of the curl
-;----------------------------------------------------------------------------
-         coefu = e1u[indice2d]*(umask())[indice2d+jpi*jpj*firstzt]
-         terreu = where(coefu EQ 0)
-         if terreu[0] NE -1 then coefu[temporary(terreu)] = !values.f_nan
-         coefu = temporary(coefu[*])#replicate(1, jpt)
-         coefu = reform(coefu, nx, ny, jpt, /over)
-;
-         coefv = e2v[indice2d]*(vmask())[indice2d+jpi*jpj*firstzt]
-         terrev = where(coefv EQ 0)
-         if terrev[0] NE -1 then coefv[temporary(terrev)] = !values.f_nan
-         coefv = temporary(coefv[*])#replicate(1, jpt)
-         coefv = reform(coefv, nx, ny, jpt, /over)
-;
-         tabf = (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
-         tabf = temporary(tabf[*])#replicate(1, jpt)
-         tabf = reform(tabf, nx, ny, jpt, /over)
-;------------------------------------------------------------
-; Calculation of the curl
-;------------------------------------------------------------
-         zu = u*temporary(coefu)
-         zv = v*temporary(coefv)
-;
-         psi =  (shift(zv, -1, 0, 0)-zv) + (zu-shift(zu, 0, -1, 0))
-         psi = tabf*psi
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         if NOT keyword_set(key_periodic) OR nx NE jpi then begin
-            psi[0, *, *] = !values.f_nan
-            psi[nx-1, *, *] = !values.f_nan
-         endif
-         psi[*, 0, *] = !values.f_nan
-         psi[*, ny-1, *] = !values.f_nan
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         terref =  where(tabf EQ 0)
-         if terref[0] NE -1 then psi[temporary(terref)] = valmask
-;----------------------------------------------------------------------------
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t', 'f']
-         if keyword_set(direc) then psi = grossemoyenne(psi,direc,/nan)
-;----------------------------------------------------------------------------
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyzt
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      date1 NE date2 AND (size(u))[0] EQ 4:BEGIN 
-         return, report('non code!')
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xy
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      ELSE:BEGIN                ;xy
-         indice2d = lindgen(jpi, jpj)
-         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
-;------------------------------------------------------------
-;------------------------------------------------------------
-         case 1 of
-            (size(u))[0] NE 2 OR (size(v))[0] NE 2: return,  -1
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-               if nxu NE nx then $
-                if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *] 
-               IF nxv NE nx THEN $
-                if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
-               IF nyu NE ny THEN $
-                if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny] 
-               IF nyv NE ny THEN $
-                if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indice2d]
-               v = v[indice2d]
-            END
-            ELSE:return,  -1
-         endcase
-;------------------------------------------------------------
-; Calculation of the curl
-;------------------------------------------------------------
-         coefu = e1u[indice2d]*(umask())[indice2d+jpi*jpj*firstzt]
-         terreu = where(coefu EQ 0)
-         if terreu[0] NE -1 then coefu[temporary(terreu)] = !values.f_nan
-         coefv = e2v[indice2d]*(vmask())[indice2d+jpi*jpj*firstzt]
-         terrev = where(coefv EQ 0)
-         if terrev[0] NE -1 then coefv[temporary(terrev)] = !values.f_nan
-         tabf = (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
-;
-         zu = u*temporary(coefu)
-         zv = v*temporary(coefv)
-
-         psi =  (shift(zv, -1, 0)-zv) + (zu-shift(zu, 0, -1))
-         psi = tabf*psi
-
-;------------------------------------------------------------
-; Edging put at !values.f_nan
-;------------------------------------------------------------
-         if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
-            psi[0, *] = !values.f_nan
-            psi[nx-1, *] = !values.f_nan
-         endif
-         psi[*, 0] = !values.f_nan
-         psi[*, ny-1] = !values.f_nan
-;
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         terref =  where(tabf EQ 0)
-         if terref[0] NE -1 then psi[temporary(terref)] = valmask
-;------------------------------------------------------------
-; for the graphic drawing
-;------------------------------------------------------------
-         domdef, (glamt[indice2d])[0, 0], (glamf[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphif[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t', 'f']
-         if keyword_set(direc) then psi = moyenne(psi,direc,/nan)
-
-      END
-;----------------------------------------------------------------------------
-   endcase
-;------------------------------------------------------------
-   if keyword_set(key_performance) THEN print, 'temps curl', systime(1)-tempsun 
-;------------------------------------------------------------
-   vargrid = 'F'
-   varname = 'vorticity'
-   return, psi
-end
Index: trunk/SRC/ToBeReviewed/CALCULS/div.pro
===================================================================
--- trunk/SRC/ToBeReviewed/CALCULS/div.pro	(revision 166)
+++ 	(revision )
@@ -1,319 +1,0 @@
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-;+
-;
-; @file_comments
-; calculation of the divergence of a 2d field
-;
-; @categories
-; Calculation
-;
-; @param UU 
-; Matrix representing coordinates of a field of vectors
-;
-; @param VV 
-; Matrix representing coordinates of a field of vectors
-;
-; @returns RES
-; A 2d matrix
-;
-; @uses
-; common.pro
-;
-; @restrictions
-; U and V matrices can be 2 or 4d.
-; Beware, to discern different configuration of U and V (xy, xyz, xyt, xyzt), 
-; we look at the variable of the common 
-;        -time which contain the calendar in IDL julian days to which U and 
-; V refered to, in the same way as the variable 
-;        -jpt which is the number of time's step to consider in time.
-; U and V arrays are cut in the same geographic domain. Because of the gap of 
-; T, U, V and F grids, it is possible that these two arrays has not the same 
-; size and refered to different indexes. In this case, arrays are re-cut on 
-; common indexes and the domain is redefined to match with these common 
-; indexes. To avoid these re-cuts, use the keyword /memeindice in domdef.pro
-;
-;
-; Points on the drawing edge are at !values.f_nan 
-;
-;
-; @history
-; Guillaume Roullet (grlod\@ipsl.jussieu.fr)
-;                      Creation : printemps 1998
-;                      Sebastien Masson (smasson\@lodyc.jussieu.fr)
-;                      adaptation to work with a reduce domain
-;                      12/1/2000;
-;
-; @version
-; $Id$
-;
-;-
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-FUNCTION div, uu, vv
-;
-  compile_opt idl2, strictarrsubs
-;
-   tempsun = systime(1)         ; For key_performance
-@common
-;
-   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
-     return, report(['This version of div is based on Arakawa C-grid.' $
-                     , 'U and V grids must therefore be defined'])
-;
-   u = litchamp(uu)
-   v = litchamp(vv)
-;
-   date1 = time[0]
-   if n_elements(jpt) EQ 0 then date2 = date1 ELSE date2 = time[jpt-1]
-   if (size(u))[0] NE (size(v))[0] then return,  -1
-
-;------------------------------------------------------------
-; We find common points between U and V
-;------------------------------------------------------------
-   indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
-   indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
-   indicex = inter(indicexu, indicexv)
-   indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
-   indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
-   indicey = inter(indiceyu, indiceyv)
-   nx = n_elements(indicex) 
-   ny = n_elements(indicey)
-   indice2d = lindgen(jpi, jpj)
-   indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
-   case 1 of
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyz
-;----------------------------------------------------------------------------
-      (size(u))[0] EQ 3 AND date1 EQ date2 :BEGIN 
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         case 1 of
-            (size(v))[0] NE 3: return,  -1
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-               case 1 of
-                  nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
-                  nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
-                  nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
-                  nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
-                  ELSE :
-               endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-            END
-            ELSE:BEGIN 
-               zdiv = -1
-               GOTO, sortie
-            end
-            
-         endcase
-;------------------------------------------------------------
-; calcul de la divergence
-;------------------------------------------------------------
-         zu = (e2u[indice2d])[*]#replicate(1, nzt)
-         zu = reform(zu, nx, ny, nzt, /over)
-         zu = temporary(u)*temporary(zu) $
-          *(umask())[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-         terreu = where(zu EQ 0)
-         if terreu[0] NE -1 then zu[temporary(terreu)] = !values.f_nan
-;
-         zv = (e1v[indice2d])[*]#replicate(1, nzt)
-         zv = reform(zv, nx, ny, nzt, /over)
-         zv = temporary(v)*temporary(zv) $
-          *(vmask())[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-         terrev = where(zv EQ 0)
-         if terrev[0] NE -1 then zv[temporary(terrev)] = !values.f_nan
-;
-         zdiv = 1e6/(e1t[indice2d]*e2t[indice2d])
-         zdiv = (zdiv)[*]#replicate(1, nzt)
-         zdiv = reform(zdiv, nx, ny, nzt, /over)
-         zdiv = temporary(zdiv)*(zu-shift(zu, 1, 0, 0)+zv-shift(zv, 0, 1, 0)) $
-          *tmask[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
-;------------------------------------------------------------
-; Edging put at !values.f_nan
-;------------------------------------------------------------
-         if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
-            zdiv[0, *, *] = !values.f_nan
-            zdiv[nx-1, *, *] = !values.f_nan
-         endif
-         zdiv[*, 0, *] = !values.f_nan
-         zdiv[*, ny-1, *] = !values.f_nan
-;
-         zdiv = temporary(zdiv)
-;
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         terre =  where(tmask[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt:lastzt] EQ 0)
-         if terre[0] NE -1 then zdiv[temporary(terre)] = valmask
-;------------------------------------------------------------
-; For the graphic drawing
-;------------------------------------------------------------
-         vargrid = 'T'
-         varname = 'div'
-         varunits = '1e6*s-1'
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t']
-         if keyword_set(direc) then  zdiv = moyenne(zdiv,direc,/nan)
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyt
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      date1 NE date2 AND (size(u))[0] EQ 3 :BEGIN 
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         case 1 of
-            (size(u))[0] NE 3 OR (size(v))[0] NE 3: return,  -1
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-               case 1 of
-                  nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
-                  nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
-                  nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
-                  nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
-                  ELSE :
-               endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
-            END
-            ELSE:return,  -1
-         endcase
-;------------------------------------------------------------
-; Calculation of the divergence
-;------------------------------------------------------------
-         zu = e2u[indice2d]*(umask())[indice2d+jpi*jpj*firstzt]
-         terreu = where(zu EQ 0)
-         if terreu[0] NE -1 then zu[temporary(terreu)] = !values.f_nan
-         zu = (zu)[*]#replicate(1, jpt)
-         zu = reform(zu, nx, ny, jpt, /over)
-         zu = temporary(u)*temporary(zu)
-;
-         zv = e1v[indice2d]*(vmask())[indice2d+jpi*jpj*firstzt]
-         terrev = where(zv EQ 0)
-         if terrev[0] NE -1 then zv[temporary(terrev)] = !values.f_nan
-         zv = (zv)[*]#replicate(1, jpt)
-         zv = reform(zv, nx, ny, jpt, /over)
-         zv = temporary(v)*temporary(zv)
-;
-         zdiv = 1e6*tmask[indice2d+jpi*jpj*firstzt]/(e1t[indice2d]*e2t[indice2d])
-         zdiv = (zdiv)[*]#replicate(1, jpt)
-         zdiv = reform(zdiv, nx, ny, jpt, /over)
-         terre =  where(zdiv EQ 0)
-         zdiv = temporary(zdiv)*(zu-shift(zu, 1, 0, 0)+zv-shift(zv, 0, 1, 0))
-;------------------------------------------------------------
-; Edging put at !values.f_nan
-;------------------------------------------------------------
-         if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
-            zdiv[0, *, *] = !values.f_nan
-            zdiv[nx-1, *, *] = !values.f_nan
-         endif
-         zdiv[*, 0, *] = !values.f_nan
-         zdiv[*, ny-1, *] = !values.f_nan
-;
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         if terre[0] NE -1 then zdiv[temporary(terre)] = valmask
-;------------------------------------------------------------
-; for the graphic drawing
-;------------------------------------------------------------
-         vargrid = 'T'
-         varname = 'div'
-         varunits = '1e6*s-1'
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t']
-         if keyword_set(direc) then  zdiv = grossemoyenne(zdiv,direc,/nan)
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyzt
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      date1 NE date2 AND (size(u))[0] EQ 4:BEGIN 
-         return, report('non code!')
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xy
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      ELSE:BEGIN                ;xy
-         indice3d = lindgen(jpi, jpj, jpk)
-         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, firstzt]
-;------------------------------------------------------------
-; extraction of U and V on the appropriated domain
-;------------------------------------------------------------
-         case 1 of
-            (size(u))[0] NE 2 OR (size(v))[0] NE 2: BEGIN 
-               zdiv = -1
-               GOTO, sortie
-            end
-            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
-             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN 
-               case 1 of
-                  nxu NE nx:if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *]
-                  nxv NE nx:if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
-                  nyu NE ny:if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny]
-                  nyv NE ny:if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
-                  ELSE :
-               endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN 
-               u = u[indice2d]
-               v = v[indice2d]
-            END
-            ELSE:return,  -1
-         endcase
-;------------------------------------------------------------
-; Calculation of the divergence
-;------------------------------------------------------------
-         zu = temporary(u)*e2u[indice2d]*(umask())[indice3d]
-         terreu = where(zu EQ 0)
-         if terreu[0] NE -1 then zu[temporary(terreu)] = !values.f_nan
-         zv = temporary(v)*e1v[indice2d]*(vmask())[indice3d]
-         terrev = where(zv EQ 0)
-         if terrev[0] NE -1 then zv[temporary(terrev)] = !values.f_nan
-         zdiv = zu-shift(zu, 1, 0)+zv-shift(zv, 0, 1)
-         zdiv = temporary(zdiv)*tmask[indice3d]/(e1t[indice2d]*e2t[indice2d])
-;------------------------------------------------------------
-; Edging put at !values.f_nan
-;------------------------------------------------------------
-         if  NOT keyword_set(key_periodic) OR nx NE jpi then begin
-            zdiv[0, *] = !values.f_nan
-            zdiv[nx-1, *] = !values.f_nan
-         endif
-         zdiv[*, 0] = !values.f_nan
-         zdiv[*, ny-1] = !values.f_nan
-;
-         zdiv = temporary(zdiv)*1e6
-;
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         terre =  where(tmask[indice3d] EQ 0)
-         if terre[0] NE -1 then zdiv[temporary(terre)] = valmask
-;------------------------------------------------------------
-; for the graphic drawing
-;------------------------------------------------------------
-         vargrid = 'T'
-         varname = 'div'
-         varunits = '1e6*s-1'
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, gridtype = ['t']
-         if keyword_set(direc) then  zdiv = moyenne(zdiv,direc,/nan)
-
-      END
-   endcase
-   
-;------------------------------------------------------------
-sortie:
-   if keyword_set(key_performance) THEN print, 'temps div', systime(1)-tempsun 
-   
-   return, zdiv
-end
Index: trunk/SRC/ToBeReviewed/CALCULS/grad.pro
===================================================================
--- trunk/SRC/ToBeReviewed/CALCULS/grad.pro	(revision 166)
+++ 	(revision )
@@ -1,302 +1,0 @@
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-;+
-; @file_comments
-;
-;
-; @categories
-;
-;
-; @param FIELD
-;
-;
-; @param DIREC
-;
-;
-; @returns
-;
-;
-; @uses
-;
-;
-; @restrictions
-;
-;
-; @examples
-;
-;
-; @history
-; Sebastien Masson (smasson\@lodyc.jussieu.fr)
-;
-; @version
-; $Id$
-;
-; @todo 
-; seb: remplir les truc!
-;-
-;------------------------------------------------------------
-;------------------------------------------------------------
-;------------------------------------------------------------
-FUNCTION grad, field, direc
-;
-  compile_opt idl2, strictarrsubs
-;
-@common
-;------------------------------------------------------------
-;
-   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
-     return, report(['This version of grad is based on Arakawa C-grid.' $
-                     , 'U and V grids must therefore be defined'])
-;
-   res = litchamp(field)
-   taille = size(res)
-   grille, mask, glam, gphi, gdep, nx, ny, nz $
-           , firstx, firsty, firstz, lastx, lasty, lastz     
-   if n_elements(valmask) EQ 0 then valmask = 1e20
-   case strupcase(vargrid) of
-      'T':BEGIN
-         case direc of
-            'x':BEGIN 
-               divi = e1u[firstx:lastx, firsty:lasty]
-               newmask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'U'
-               firstxu = firstxt & lastxu = lastxt & nxu = nxt 
-               firstyu = firstyt & lastyu = lastyt & nyu = nyt
-            END
-            'y':BEGIN
-               divi = e2v[firstx:lastx, firsty:lasty]
-               newmask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'V'
-               firstxv = firstxt & lastxv = lastxt & nxv = nxt 
-               firstyv = firstyt & lastyv = lastyt & nyv = nyt
-            END
-            'z':BEGIN
-               divi = e3w[firstz:lastz]
-               newmask = mask
-               vargrid = 'W'
-               firstzw = firstzt & lastzw = lastzt & nzw = nzt 
-            END
-            ELSE:return, report('Bad definition of direction argument')
-         ENDCASE
-      END
-      'W':BEGIN
-         case direc of
-            'x':BEGIN 
-               divi = e1u[firstx:lastx, firsty:lasty]
-               newmask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'U'
-               firstxu = firstxt & lastxu = lastxt & nxu = nxt 
-               firstyu = firstyt & lastyu = lastyt & nyu = nyt
-            END
-            'y':BEGIN
-               divi = e2v[firstx:lastx, firsty:lasty]
-               newmask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'V'
-               firstxv = firstxt & lastxv = lastxt & nxv = nxt 
-               firstyv = firstyt & lastyv = lastyt & nyv = nyt
-            END
-            'z':BEGIN
-               divi = e3t[firstz:lastz]
-               newmask = mask
-               vargrid = 'T'
-               firstzt = firstzw & lastzt = lastzw & nzt = nzw 
-            END
-            ELSE:return, report('Bad definition of direction argument')
-         endcase
-      END
-      'U':BEGIN
-         case direc of
-            'x':BEGIN
-               divi = (shift(e1t, -1, 0))[firstx:lastx, firsty:lasty]
-               newmask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'T'
-               firstxt = firstxu & lastxt = lastxu & nxt = nxu 
-               firstyt = firstyu & lastyt = lastyu & nyt = nyu
-            END
-            'y':BEGIN
-               divi = e2f[firstx:lastx, firsty:lasty]
-               newmask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'F'
-               firstxf = firstxu & lastxf = lastxu & nxf = nxu 
-               firstyf = firstyu & lastyf = lastyu & nyf = nyu
-            END
-            'z':BEGIN
-               divi = e3w[firstz:lastz]
-               newmask = mask
-               vargrid = 'W'
-               firstzw = firstzt & lastzw = lastzt & nzw = nzt 
-            END
-            ELSE:return, report('Bad definition of direction argument')
-         endcase
-      END
-      'V':BEGIN
-         case direc of
-            'x':BEGIN
-               divi = e1f[firstx:lastx, firsty:lasty]
-               newmask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'F'
-               firstxf = firstxv & lastxf = lastxv & nxf = nxv
-               firstyf = firstyv & lastyf = lastyv & nyf = nyv
-            END
-            'y':BEGIN
-               divi = (shift(e2t, 0, -1))[firstx:lastx, firsty:lasty]
-               newmask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
-               vargrid = 'T'
-               firstxt = firstxv & lastxt = lastxv & nxt = nxv
-               firstyt = firstyv & lastyt = lastyv & nyt = nyv
-            END
-            'z':BEGIN
-               divi = e3w[firstz:lastz]
-               newmask = mask
-               vargrid = 'W'
-               firstzw = firstzt & lastzw = lastzt & nzw = nzt 
-            END
-            ELSE:return, report('Bad definition of direction argument')
-         endcase
-      END
-;       'F':BEGIN
-;          case direc of
-;             'x':divi = (shift(e1v, -1, 0))[firstx:lastx, firsty:lasty]
-;             'y':divi = (shift(e2u, 0, -1))[firstx:lastx, firsty:lasty]
-;             'z':divi = e3w[firstz:lastz]
-;             ELSE:return, report('Bad definition of direction argument')
-;          endcase
-;       END
-      ELSE:return, report('Bad definition of vargrid')
-   ENDCASE
-   res = fitintobox(res)
-   IF n_elements(res) EQ 1 AND res[0] EQ -1 THEN return, res
-   case 1 of
-;----------------------------------------------------------------------------
-;xy
-;----------------------------------------------------------------------------
-      taille[0] EQ 2:BEGIN
-         earth = where(mask[*, *, firstz] EQ 0)
-         if earth[0] NE -1 then res[earth] = !values.f_nan
-         case direc of
-            'x':BEGIN 
-               res = (shift(res, -1, 0)-res)/divi
-               if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0)
-            END
-            'y':BEGIN 
-               res = (shift(res, 0, -1)-res)/divi
-               res[*, ny-1] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'U' then res =   shift(res, 0, 1)           
-            END
-            ELSE:return,  report('Bad definition of direction argument for the type of array')
-         ENDCASE
-         earth = where(newmask[*, *, firstz] EQ 0)
-         if earth[0] NE -1 then res[earth] = valmask
-      END
-;----------------------------------------------------------------------------
-;xyt
-;----------------------------------------------------------------------------
-      taille[0] EQ 3 AND jpt NE 1:BEGIN
-         earth = where(mask[*, *, firstz] EQ 0)
-         if earth[0] NE -1 then BEGIN
-            earth = earth#replicate(1, jpt)+replicate(1, n_elements(earth))#(nx*ny*lindgen(jpt))
-            res[earth[*]] = !values.f_nan
-         ENDIF
-         divi = divi[*]#replicate(1, jpt)
-         case direc of
-            'x':BEGIN 
-               res = (shift(res, -1, 0, 0)-res)/divi
-               if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'V'  then res = shift(res, 1, 0, 0)
-            END
-            'y':BEGIN 
-               res = (shift(res, 0, -1, 0)-res)/divi
-               res[*, ny-1, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'U'  then res = shift(res, 0, 1, 0)        
-            END
-            ELSE:return,  report('Bad definition of direction argument for the type of array')
-         ENDCASE
-         earth = where(newmask[*, *, firstz] EQ 0)
-         if earth[0] NE -1 THEN res[earth] = valmask
-      END
-;----------------------------------------------------------------------------
-;xyz
-;----------------------------------------------------------------------------
-      taille[0] EQ 3 AND jpt EQ 1:BEGIN
-         earth = where(mask EQ 0)
-         if earth[0] NE -1 then res[earth] = !values.f_nan
-         case direc OF
-            'x':BEGIN 
-               divi = divi[*]#replicate(1, nz)
-               res = (shift(res, -1, 0, 0)-res)/divi
-               if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0, 0)
-            END
-            'y':BEGIN 
-               divi = divi[*]#replicate(1, nz)
-               res = (shift(res, 0, -1, 0)-res)/divi
-               res[*, ny-1, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'U'  then res = shift(res, 0, 1, 0)        
-            END
-            'z':BEGIN
-               divi = reform(replicate(1, nx*ny)#divi, nx, ny, nz)
-               if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz)
-               if vargrid EQ 'W' THEN BEGIN 
-                  res = (shift(res, 0, 0, 1)-res)/divi
-                  res[*, *, 0] = !values.f_nan
-               ENDIF ELSE BEGIN
-                  res = (res-shift(res, 0, 0, -1))/divi
-                  res[*, *, nz-1] = !values.f_nan
-               ENDELSE
-               if earth[0] NE -1 then res[earth] = valmask
-            END
-         ENDCASE
-      END
-;------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyzt
-;----------------------------------------------------------------------------
-      taille[0] EQ 4:BEGIN
-         earth = where((mask[*]#replicate(1, jpt)) EQ 0)
-         if earth[0] NE -1 then res[earth] = !values.f_nan
-         case direc OF
-            'x':BEGIN 
-               divi = divi[*]#replicate(1, nz*jpt)
-               res = (shift(res, -1, 0, 0, 0)-res)/divi
-               if key_periodic EQ 0 OR nx NE jpi THEN res[nx-1, *, *, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'V' then res = shift(res, 1, 0, 0, 0)
-            END
-            'y':BEGIN 
-               divi = divi[*]#replicate(1, nz*jpt)
-               res = (shift(res, 0, -1, 0, 0)-res)/divi
-               res[*, ny-1, *, *] = !values.f_nan
-               if vargrid EQ 'T' OR vargrid EQ 'U'  then res = shift(res, 0, 1, 0, 0)        
-            END
-            'z':BEGIN
-               divi = replicate(1, nx*ny)#divi
-               divi = reform(divi[*]#replicate(1, jpt), nx, ny, nz, jpt, /over)
-               if nx EQ 1 OR ny EQ 1 then res = reform(res, nx, ny, nz, jpt)
-               if vargrid EQ 'W' THEN BEGIN 
-                  res = (shift(res, 0, 0, 1, 0)-res)/divi
-                  res[*, *, 0, *] = !values.f_nan
-               ENDIF ELSE BEGIN
-                  res = (res-shift(res, 0, 0, -1, 0))/divi
-                  res[*, *, nz-1, *] = !values.f_nan
-               ENDELSE
-            END
-         ENDCASE
-         if earth[0] NE -1 then res[earth] = valmask
-      END
-;------------------------------------------------------------
-;------------------------------------------------------------
-   endcase
-   varname = 'grad of '+varname
-   varunit = varunit+'/m'
-
-
-;------------------------------------------------------------
-   return, res
-end
-
-
-
-
-
-
Index: trunk/SRC/ToBeReviewed/INIT/initncdf.pro
===================================================================
--- trunk/SRC/ToBeReviewed/INIT/initncdf.pro	(revision 166)
+++ trunk/SRC/ToBeReviewed/INIT/initncdf.pro	(revision 167)
@@ -123,5 +123,6 @@
   if keyword_set(xaxisname) then xaxisname = strlowcase(xaxisname) ELSE xaxisname = 'x'
   xvarid = where(namevar EQ xaxisname OR namevar EQ 'longitude' $
-                 OR namevar EQ 'nav_lon' OR namevar EQ 'lon')
+                 OR namevar EQ 'nav_lon' OR namevar EQ 'lon' $
+                 OR namevar EQ 'NbLongitudes')
   xvarid = xvarid[0]
   if xvarid EQ -1 then begin
@@ -145,5 +146,7 @@
 ; find the yaxis
   if keyword_set(yaxisname) then yaxisname = strlowcase(yaxisname) ELSE yaxisname = 'y'
-  yvarid = where(namevar EQ yaxisname OR namevar EQ 'latitude' OR namevar EQ 'nav_lat' OR namevar EQ 'lat')
+  yvarid = where(namevar EQ yaxisname OR namevar EQ 'latitude' $
+                 OR namevar EQ 'nav_lat' OR namevar EQ 'lat' $
+                 OR namevar EQ 'NbLatitudes')
   yvarid = yvarid[0]
   if yvarid EQ -1 then begin
Index: trunk/SRC/ToBeReviewed/WIDGET/AUTOUR_de_XXX/scanfile.pro
===================================================================
--- trunk/SRC/ToBeReviewed/WIDGET/AUTOUR_de_XXX/scanfile.pro	(revision 166)
+++ trunk/SRC/ToBeReviewed/WIDGET/AUTOUR_de_XXX/scanfile.pro	(revision 167)
@@ -114,18 +114,50 @@
     namedim[dimiq] = strlowcase(tmpname)
   ENDFOR
-; we are looking for a x dimension...
-  dimidx = where(namedim EQ 'x' OR strmid(namedim, 0, 3) EQ 'lon' OR strmid(namedim, 0, 3) EQ 'xi_' OR namedim EQ 'xt_i7_156')
-  dimidx = dimidx[0]
-  if dimidx EQ -1 then begin
-    print, 'one of the dimensions must have the name: ''x'' or ''lon...'' or ''xi_...'' or ''xt_i7_156'''
-    stop
-  endif
-; we are looking for a y dimension...
-  dimidy = where(namedim EQ 'y' OR strmid(namedim, 0, 3) EQ 'lat' OR strmid(namedim, 4) EQ 'eta_' OR namedim EQ 'yt_j6_75')
-  dimidy = dimidy[0]
-  if dimidy EQ -1 then begin
-    print, 'one of the dimensions must have the name: ''y'' or ''lat...'' or ''eta_...'' or ''yt_j6_75'''
-    stop
-  endif
+; we are looking for a x dimension with a name matching one of the following regular expression:
+  testname = ['longitude', 'lon', 'x', 'longitude*', 'lon*', 'x*', '*longitude*', '*lon*', '*x*']
+  cnt = -1
+  ii = 0
+  WHILE cnt NE 1 AND ii LT n_elements(testname) DO BEGIN
+    dimidx = where(strmatch(namedim, testname[ii]) EQ 1, cnt)
+    ii = ii+1
+  ENDWHILE
+  CASE cnt OF
+    0:begin
+      dummy = report(['none of the dimensions name matches one of the following regular expression:' $
+                      , '''longitude'', ''lon'', ''x'', ''longitude*'', ''lon*'', ''x*'', ''*longitude*'', ''*lon*'', ''*x*''' $
+                      , ' => we cannot find the x dimension'])
+      stop
+    END
+    1:dimidx = dimidx[0]
+    ELSE:begin
+      dummy = report(['several (and not one unique) dimensions name matches the following regular expression:' $
+                      , '''longitude'', ''lon'', ''x'', ''longitude*'', ''lon*'', ''x*'', ''*longitude*'', ''*lon*'', ''*x*''' $
+                      , ' => we cannot find the x dimension'])
+      stop
+    END
+  ENDCASE
+; we are looking for a y dimension with a name matching one of the following regular expression:
+  testname = ['latitude', 'lat', 'y', 'latitude*', 'lat*', 'y*', 'eta_*', '*latitude*', '*lat*', '*y*']
+  cnt = -1
+  ii = 0
+  WHILE cnt NE 1 AND ii LT n_elements(testname) DO BEGIN
+    dimidy = where(strmatch(namedim, testname[ii]) EQ 1, cnt)
+    ii = ii+1
+  ENDWHILE
+  CASE cnt OF
+    0:begin
+      dummy = report(['none of the dimensions name matches one of the following regular expression:' $
+                      , '''latitude'', ''lat'', ''y'', ''latitude*'', ''lat*'', ''y*'', ''eta_*'', ''*latitude*'', ''*lat*'', ''*y*''' $
+                      , ' => we cannot find the y dimension'])
+      stop
+    END
+    1:dimidy = dimidy[0]
+    ELSE:begin
+      dummy = report(['several (and not one unique) dimensions name matches the following regular expression:' $
+                      , '''latitude'', ''lat'', ''y'', ''latitude*'', ''lat*'', ''y*'', ''eta_*'', ''*latitude*'', ''*lat*'', ''*y*''' $
+                      , ' => we cannot find the x dimension'])
+      stop
+    END
+  ENDCASE
 ;------------------------------------------------------------
 ; name of all variables
@@ -156,7 +188,9 @@
     varid = 0
     repeat BEGIN
-      invar = ncdf_varinq(cdfid, varid) 
-      varid = varid+1
-    endrep until n_elements(invar.dim) EQ 1 AND invar.dim[0] EQ infile.recdim
+      IF varid LT infile.nvars THEN BEGIN
+        invar = ncdf_varinq(cdfid, varid) 
+        varid = varid+1        
+      ENDIF ELSE varid = 0
+    endrep until (n_elements(invar.dim) EQ 1 AND invar.dim[0] EQ infile.recdim) OR (varid EQ 0)
     varid = varid-1
 ;
@@ -179,5 +213,5 @@
         for attiq = 0, invar.natts-1 do attnames[attiq] = ncdf_attname(cdfid, varid, attiq)
         if (where(attnames EQ 'units'))[0] EQ -1 then BEGIN
-          dummy = report('Attribut ''units'' not found for the variable '+varid.name+'!C we create a fake calendar ...')
+          dummy = report('Attribut ''units'' not found for the variable '+invar.name+'!C we create a fake calendar ...')
           fakecal = 1
           time = date0fk + lindgen(jpt)
@@ -218,8 +252,8 @@
           ENDELSE
 ;
-; BEWARE we have to recuperate the calendar attribute and ajust TIME by consequence...
-;
-;
-; We pass TIME in IDL julian days
+; BEWARE we have to get back the calendar attribute and ajust time by consequence...
+;
+;
+; We pass time in IDL julian days
 ;
           unite = strlowcase(unite)
Index: trunk/SRC/ToBeReviewed/WIDGET/COMPOUND_WIDGET/cw_domain.pro
===================================================================
--- trunk/SRC/ToBeReviewed/WIDGET/COMPOUND_WIDGET/cw_domain.pro	(revision 166)
+++ trunk/SRC/ToBeReviewed/WIDGET/COMPOUND_WIDGET/cw_domain.pro	(revision 167)
@@ -662,10 +662,7 @@
 ;
   min2 = gdep2[indice1]
-  if indice2 EQ jpk-1 then BEGIN 
-    max2 = max([gdept, gdepw])
-    max2 = strtrim(string(max2, format = '(e8.0)'), 1)
-    max2 = float('1'+strmid(max2, 1))+float(max2)
-  ENDIF ELSE max2 = gdep1[indice2+1]
-  if max2 EQ min2 then max2 = min2+1
+  if indice2 EQ jpk-1 then max2 = ceil(max([gdept, gdepw])) $
+  ELSE max2 = gdep1[indice2+1]
+  if floor(max2) LE floor(min2) then max2 = min2+1
   rien = cw_slider_pm(basez, value = (min2 > boxzoom[4]) > boxzoom[5] < max2 $
                       , uvalue = {name:'depth2'}, minimum = min2, maximum =  max2 $
