Index: /trunk/SRC/Computation/curl.pro
===================================================================
--- /trunk/SRC/Computation/curl.pro	(revision 313)
+++ /trunk/SRC/Computation/curl.pro	(revision 314)
@@ -2,5 +2,6 @@
 ;
 ; @file_comments
-; Calculate the vertical component of the curl of a field of horizontal vectors
+; Calculate the vertical component of the curl of a vectors field
+; located on Arakawa C-grid.
 ;
 ; @categories
@@ -8,5 +9,5 @@
 ;
 ; @param UU
-; Matrix representing the zonal coordinates (U point) of a field of vectors
+; Matrix representing the zonal coordinates (at U point) of a field of vectors
 ; A 2D (xy), 3D (xyz or yt) or a structure readable by <pro>litchamp</pro> and containing
 ; a 2D (xy), 3D (xyz or yt) array (4D case is not coded yet).
@@ -14,5 +15,5 @@
 ;
 ; @param VV
-; Matrix representing the meridional coordinates (V point) of a field of vectors
+; Matrix representing the meridional coordinates (at V point) of a field of vectors
 ; A 2D (xy), 3D (xyz or yt) or a structure readable by <pro>litchamp</pro> and containing
 ; a 2D (xy), 3D (xyz or yt) array (4D case is not coded yet).
@@ -24,11 +25,18 @@
 ; with a mean done in the DIREC direction
 ;
+; @keyword MILLION {default=0}{type=scalar: 0 or 1}
+; Activate to multiply returned array by 1.e6
+; 
+; @keyword _EXTRA
+; Used to declare that this routine accepts inherited keywords
+;
 ; @returns
-; the vertical component of the curl of the input data (with the same size)
+; the vertical component of the curl of the input data (with the same
+; size) located at F point (see restrictions)
 ;
 ; @uses
 ; cm_4cal
 ; cm_4data
-; cm_4mmesh
+; cm_4mesh
 ;
 ; @restrictions
@@ -48,5 +56,5 @@
 ;
 ; @examples
-; IDL> \@tst_initorca2
+; IDL> @tst_initorca2
 ; IDL> plt, curl(dist(jpi,jpj), dist(jpi,jpj))
 ;
@@ -65,5 +73,5 @@
 ;-
 ;
-FUNCTION curl, uu, vv, DIREC = DIREC
+FUNCTION curl, uu, vv, DIREC = direc, MILLION = million, _EXTRA = ex
 ;
   compile_opt idl2, strictarrsubs
@@ -73,5 +81,5 @@
 @cm_4mesh
 ;
-  tempsun = systime(1)          ; To key_performance
+  time1 = systime(1)          ; for key_performance
 ;
   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
@@ -79,4 +87,12 @@
                      , 'U and V grids must therefore be defined'])
 ;
+  gr = litchamp(uu, /grid)
+  IF gr NE '' THEN BEGIN 
+    IF gr NE 'U' THEN return, report('the first parameter is not located on U grid, but on '+ strtrim(gr, 2) +'grid')
+  ENDIF 
+  gr = litchamp(vv, /grid)
+  IF gr NE '' THEN BEGIN 
+    IF gr NE 'V' THEN return, report('the second parameter is not located on V grid, but on '+ strtrim(gr, 2) +'grid')
+  ENDIF 
   u = litchamp(uu)
   v = litchamp(vv)
@@ -103,8 +119,9 @@
   vargrid = 'F'
   varname = 'vorticity'
-  varunits = 's-1'
+  IF keyword_set(million) THEN varunits = '1.e6*'+varunit+'/m' ELSE varunits = varunit+'/m'
+  IF keyword_set(million) THEN scale = 1.e6 ELSE scale = 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
+  firstxf = indicex[0] & lastxf = indicex[0]+nx-1 & nxf = nx
+  firstyf = indicey[0] & lastyf = indicey[0]+ny-1 & nyf = ny
 ;----------------------------------------------------------------------------
 ;----------------------------------------------------------------------------
@@ -148,5 +165,5 @@
       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] $
+      tabf = scale * (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)
@@ -196,9 +213,9 @@
           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
+        ELSE:return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      ENDCASE
 ;----------------------------------------------------------------------------
 ; curl computation
@@ -214,5 +231,5 @@
       coefv = temporary(coefv[*])#replicate(1., jpt)
 ;
-      tabf = (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
+      tabf = scale * (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
       tabf = reform(temporary(tabf[*])#replicate(1., jpt), nx, ny, jpt, /overwrite)
       landf = where(tabf EQ 0)
@@ -279,5 +296,5 @@
       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])
+      tabf = scale * (fmask())[indice2d+jpi*jpj*firstzt]/(e1f[indice2d]*e2f[indice2d])
       landf =  where(tabf EQ 0)
 ;
@@ -305,5 +322,5 @@
   ENDCASE
 ;------------------------------------------------------------
-  if keyword_set(key_performance) THEN print, 'temps curl', systime(1)-tempsun
+  if keyword_set(key_performance) THEN print, 'time curl', systime(1)-time1
 
   return, psi
Index: /trunk/SRC/Computation/div.pro
===================================================================
--- /trunk/SRC/Computation/div.pro	(revision 313)
+++ /trunk/SRC/Computation/div.pro	(revision 314)
@@ -2,5 +2,5 @@
 ;
 ; @file_comments
-; compute the horizontal divergence of a vectors field
+; compute the horizontal divergence of a vectors field located on Arakawa C-grid.
 ;
 ; @categories
@@ -8,5 +8,5 @@
 ;
 ; @param UU
-; Matrix representing the zonal coordinates (U point) of a field of vectors
+; Matrix representing the zonal coordinates (at U point) of a field of vectors
 ; A 2D (xy), 3D (xyz or yt) or a structure readable by <pro>litchamp</pro> and containing
 ; a 2D (xy), 3D (xyz or yt) array (4D case is not coded yet).
@@ -14,5 +14,5 @@
 ;
 ; @param VV
-; Matrix representing the meridional coordinates (V point) of a field of vectors
+; Matrix representing the meridional coordinates (at V point) of a field of vectors
 ; A 2D (xy), 3D (xyz or yt) or a structure readable by <pro>litchamp</pro> and containing
 ; a 2D (xy), 3D (xyz or yt) array (4D case is not coded yet).
@@ -24,6 +24,13 @@
 ; (stupid ?) with a mean done in the DIREC direction
 ;
+; @keyword MILLION {default=0}{type=scalar: 0 or 1}
+; Activate to multiply returned array by 1.e6
+; 
+; @keyword _EXTRA
+; Used to declare that this routine accepts inherited keywords
+;
 ; @returns
-; the divergence of the input data (with the same size)
+; the divergence of the input data (with the same size) located at T
+;point (see restrictions).
 ;
 ; @uses
@@ -64,5 +71,5 @@
 ;-
 ;
-FUNCTION div, uu, vv, DIREC = direc
+FUNCTION div, uu, vv, DIREC = direc, MILLION = million, _EXTRA = ex
 ;
   compile_opt idl2, strictarrsubs
@@ -72,5 +79,5 @@
 @cm_4mesh
 ;
-  tempsun = systime(1)          ; For key_performance
+  time1 = systime(1)          ; for key_performance
 ;
   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
@@ -78,4 +85,12 @@
                      , 'U and V grids must therefore be defined'])
 ;
+  gr = litchamp(uu, /grid)
+  IF gr NE '' THEN BEGIN 
+    IF gr NE 'U' THEN return, report('the first parameter is not located on U grid, but on '+ strtrim(gr, 2) +'grid')
+  ENDIF 
+  gr = litchamp(vv, /grid)
+  IF gr NE '' THEN BEGIN 
+    IF gr NE 'V' THEN return, report('the second parameter is not located on V grid, but on '+ strtrim(gr, 2) +'grid')
+  ENDIF 
   u = litchamp(uu)
   v = litchamp(vv)
@@ -102,5 +117,6 @@
   vargrid = 'T'
   varname = 'div'
-  varunits = '1.e6*s-1'
+  IF keyword_set(million) THEN varunits = '1.e6*'+varunit+'/m' ELSE varunits = varunit+'/m'
+  IF keyword_set(million) THEN scale = 1.e6 ELSE scale = 1.
   if n_elements(valmask) EQ 0 THEN valmask = 1.e20
   firstxt = indicex[0] & lastxt = indicex[0]+nx-1 & nxt = nx
@@ -132,6 +148,9 @@
           v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, *]
         END
-        ELSE:return, -1
-      endcase
+        ELSE:return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      ENDCASE
 ;------------------------------------------------------------
 ; divergence computation
@@ -147,5 +166,5 @@
       zv = temporary(v) * temporary(zv)
 ;
-      zdiv = (e1t[indice2d]*e2t[indice2d])[*]#replicate(1.e6, nzt)
+      zdiv = (scale / (e1t[indice2d]*e2t[indice2d]))[*]#replicate(1., nzt)
       zdiv = ( zu - shift(zu, 1, 0, 0) + zv - shift(zv, 0, 1, 0) ) * temporary(zdiv)
 ;------------------------------------------------------------
@@ -205,5 +224,5 @@
       zv = temporary(v) * temporary(zv)
 ;
-      zdiv = (e1t[indice2d]*e2t[indice2d])[*]#replicate(1.e6, jpt)
+      zdiv = (scale / (e1t[indice2d]*e2t[indice2d]))[*]#replicate(1., jpt)
       zdiv = ( zu - shift(zu, 1, 0, 0) + zv - shift(zv, 0, 1, 0) ) * temporary(zdiv)
 ;------------------------------------------------------------
@@ -272,5 +291,5 @@
       zv = temporary(v) * temporary(zv)
 
-      zdiv = 1.e6 / (e1t[indice2d]*e2t[indice2d])
+      zdiv = scale / (e1t[indice2d]*e2t[indice2d])
       zdiv = ( zu - shift(zu, 1, 0) + zv - shift(zv, 0, 1) ) * temporary(zdiv)
 ;------------------------------------------------------------
@@ -293,5 +312,5 @@
   ENDCASE
 ;------------------------------------------------------------
-  if keyword_set(key_performance) THEN print, 'temps div', systime(1)-tempsun
+  if keyword_set(key_performance) THEN print, 'time div', systime(1)-time1
 
   return, zdiv
Index: /trunk/SRC/Computation/grad.pro
===================================================================
--- /trunk/SRC/Computation/grad.pro	(revision 313)
+++ /trunk/SRC/Computation/grad.pro	(revision 314)
@@ -2,5 +2,5 @@
 ;
 ; @file_comments
-; compute the gradient of a variable
+; compute the gradient of a variable located on Arakawa C-grid.
 ;
 ; @categories
@@ -10,6 +10,5 @@
 ; 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 
-; <pro>litchamp</pro>
-; and containing a 2D (xy), 3D (xyz or yt) or 4D (xyzt) array.
+; <pro>litchamp</pro> and containing a 2D (xy), 3D (xyz or yt) or 4D (xyzt) array.
 ; Note that the dimension of the array must suit the domain dimension.
 ;
@@ -17,11 +16,18 @@
 ; the gradient direction: 'x', 'y', 'z'
 ;
+; @keyword MILLION {default=0}{type=scalar: 0 or 1}
+; Activate to multiply returned array by 1.e6
+; 
+; @keyword _EXTRA
+; Used to declare that this routine accepts inherited keywords
+;
 ; @returns
-; the gradient of the input data with the same size 2D, 3D or 4D array
+; the gradient of the input data with the same size 2D, 3D or 4D
+; array, located on the appropriate grid (see restrictions)
 ;
 ; @uses
 ; cm_4cal
 ; cm_4data
-; cm_4mmesh
+; cm_4mesh
 ;
 ; @restrictions
@@ -47,5 +53,5 @@
 ;-
 ;
-FUNCTION grad, field, direc
+FUNCTION grad, field, direc, MILLION = million, _EXTRA = ex
 ;
   compile_opt idl2, strictarrsubs
@@ -54,4 +60,6 @@
 @cm_4data  ; for varname, vargrid, vardate, varunit, valmask
 @cm_4mesh
+;
+  time1 = systime(1)          ; for key_performance
 ;------------------------------------------------------------
 ;
@@ -67,5 +75,5 @@
   if n_elements(valmask) EQ 0 then valmask = 1.e20
   varname = 'grad of '+varname
-  varunit = varunit+'/m'
+    IF keyword_set(million) THEN varunit = '1.e6*'+varunit+'/m' ELSE varunit = varunit+'/m'
   case strupcase(vargrid) of
     'T':BEGIN
@@ -180,4 +188,5 @@
     ELSE:return, report('Bad definition of vargrid')
   ENDCASE
+  IF keyword_set(million) THEN divi = temporary(divi)*1.e-6
   res = fitintobox(temporary(res))
   IF n_elements(res) EQ 1 AND res[0] EQ -1 THEN return, res
@@ -311,4 +320,7 @@
     ELSE:return, report('input array must have 2, 3 or 4 dimensions')
   ENDCASE
+
+  if keyword_set(key_performance) THEN print, 'time curl', systime(1)-time1
+
   return, res
 END
Index: /trunk/SRC/Computation/norm.pro
===================================================================
--- /trunk/SRC/Computation/norm.pro	(revision 314)
+++ /trunk/SRC/Computation/norm.pro	(revision 314)
@@ -0,0 +1,453 @@
+;+
+;
+; @file_comments
+; calculate the norm of vectors field located on Arakawa C-grid
+;
+; @categories
+; Calculation
+;
+; @param UU {in}{required}
+; Matrix representing the zonal coordinates (at U/V point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) or a structure readable by
+; <pro>litchamp</pro> and containing a 2D (xy), 3D (xyz or yt), 4D (xyzt) array.
+; Note that the dimension of the array must suit the domain dimension.
+;
+; @param VV {in}{required}
+; Matrix representing the meridional coordinates (at V/U point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) or a structure readable by
+; <pro>litchamp</pro> and containing a 2D (xy), 3D (xyz or yt), 4D (xyzt) array.
+; Note that the dimension of the array must suit the domain dimension.
+;
+; @keyword DIREC
+; 't' 'x' 'y' 'z' 'xys' 'xz' 'yz' 'xyz' 'xt' 'yt' 'zt' 'xyt'
+;       'xzt' 'yzt' 'xyzt' Direction on which do averages
+;
+; @keyword _EXTRA
+; Used to declare that this routine accepts inherited keywords
+;
+; @returns
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) Array
+;
+; @uses
+; @cm_4mesh
+; @cm_4data
+; @cm_4cal
+;
+; @restrictions
+; The norm is calculated on points T. To do this calculation, we average
+; field U and V on points T before calculate the norm. At the edge of
+; coast and of domain, we can not calculate fields U and V at points T,
+; that is why these points are at value !values.f_nan.
+;
+; When we calculate on a reduce geographic domain, field U and V have not
+; necessarily the same number of point. In this case, we recut U and V to
+; keep only common points. We profit of this to redo a <pro>domdef</pro> which redefine
+; a geographic domain on which fields U and V are extracted on same points
+;
+; To know what type of array we work with, we  test its size and dates
+; gave by time[0] and time[jpt-1] to know if thee is a temporal dimension.
+; Before to start norm, make sure that time and jpt are defined how
+; they have to!
+;
+; @examples
+; To calculate the average of the norm of streams on all the domain
+; between 0 and 50:
+;      IDL> domdef, 0, 50
+;      IDL> res = norm(un, vn, dir = 'xyz')
+;
+; @history
+; Sebastien Masson (smasson\@lodyc.jussieu.fr)
+;                       9/6/1999
+; @version
+; $Id$
+;
+;-
+;
+FUNCTION norm, uu, vv, DIREC = direc, _EXTRA = ex
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4mesh
+@cm_4data
+@cm_4cal
+  IF NOT keyword_set(key_forgetold) THEN BEGIN
+@updatenew
+@updatekwd
+  ENDIF
+;---------------------------------------------------------
+  time1 = 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 norm is based on Arakawa C-grid.' $
+                     , 'U and V grids must therefore be defined'])
+;------------------------------------------------------------
+  if NOT keyword_set(direc) then direc = 0
+; construction of u and v at points T
+  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
+;
+  grilleu = litchamp(uu, /grid)
+  if grilleu EQ '' then grilleu = 'U'
+  grillev = litchamp(vv, /grid)
+  if grillev EQ '' then grillev = 'V'
+  IF grilleu EQ 'V' AND grillev EQ 'U' THEN inverse = 1
+  IF grilleu EQ 'T' AND grillev EQ 'T' THEN BEGIN
+    interpolle  = 0
+    return, report('Case not coded, but easy to do...!')
+  ENDIF ELSE interpolle = 1
+  if keyword_set(inverse) then begin
+    tmp = u
+    u = temporary(v)
+    v = temporary(tmp)
+  endif
+;------------------------------------------------------------
+; 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)
+;----------------------------------------------------------------------------
+  vargrid = 'T'
+  varname = 'norm'
+  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
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    (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]
+      indice3d = lindgen(jpi, jpj, jpk)
+      indice3d = indice3d[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
+;------------------------------------------------------------
+; 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
+          case (size(u))[3] OF
+            nzt: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
+            jpk:BEGIN
+              if nxu NE nx then $
+                 if indicex[0] EQ firstxu then u = u[0:nx-1, *, firstzt:lastzt] ELSE u = u[1: nx, *, firstzt:lastzt]
+              IF nxv NE nx THEN $
+                 if indicex[0] EQ firstxv then v = v[0:nx-1, *, firstzt:lastzt] ELSE v = v[1: nx, *, firstzt:lastzt]
+              IF nyu NE ny THEN $
+                 if indicey[0] EQ firstyu then u = u[*, 0:ny-1, firstzt:lastzt] ELSE u = u[*, 1: ny, firstzt:lastzt]
+              IF nyv NE ny THEN $
+                 if indicey[0] EQ firstyv then v = v[*, 0:ny-1, firstzt:lastzt] ELSE v = v[*, 1: ny, firstzt:lastzt]
+            end
+            ELSE: return, report(['the third dimension of u (' + strtrim((size(u))[3], 1) $
+                                  +') must be equal to nzt (' + strtrim(nzt, 1) + ') or jpk ('+strtrim(jpk, 1)+')'])
+          endcase
+        END
+        (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
+           (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
+          u = u[indice3d]
+          v = v[indice3d]
+        END
+        ELSE: return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      endcase
+;------------------------------------------------------------------
+; We reshape u and v to make sure that no dimension has been erased
+;------------------------------------------------------------------
+      if nzt EQ 1 then begin
+        u = reform(u, nx, ny, nzt, /over)
+        v = reform(v, nx, ny, nzt, /over)
+      endif
+;------------------------------------------------------------------
+; construction of u and v at points T
+;-----------------------------------------------------------
+      a = u[0, *, *]
+      u = (u+shift(u, 1, 0, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then u[0, *, *] = a
+      a = v[*, 0, *]
+      v = (v+shift(v, 0, 1, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then v[*, 0, *] = a
+;----------------------------------------------------------------------------
+; attribution of the mask and of longitude and latitude arrays
+;----------------------------------------------------------------------------
+      mask = tmask[indice3d]
+      if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
+;-----------------------------------------------------------
+      if n_elements(valmask) EQ 0 THEN valmask = 1e20
+      landu = where(u GE valmask/10)
+      if landu[0] NE -1 then u[landu] = 0
+      landv = where(v GE valmask/10)
+      if landv[0] NE -1 then v[landv] = 0
+      res = sqrt(u^2+v^2)
+      if NOT keyword_set(key_periodic) OR nx NE jpi then res[0, *, *] = !values.f_nan
+      res[*, 0, *] = !values.f_nan
+      mask = where(mask eq 0)
+      IF mask[0] NE -1 THEN res[mask] = valmask
+; All kind of average
+      domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0], (gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
+      if keyword_set(direc) then res = moyenne(res, direc, /nan, boxzoom = boxzoom, /nodomdef)
+;-----------------------------------------------------------
+    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:return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      endcase
+;------------------------------------------------------------------
+; construction of u and v at points T
+;-----------------------------------------------------------
+      a = u[0, *, *]
+      u = (u+shift(u, 1, 0, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then u[0, *, *] = a
+      a = v[*, 0, *]
+      v = (v+shift(v, 0, 1, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then v[*, 0, *] = a
+;----------------------------------------------------------------------------
+; attribution of the mask and of longitude and latitude arrays.
+; We recover the complete grid to establish a big mask extent in the four
+; direction to cover pointsfor which a land point has been
+; considerated (make a small drawing)
+;----------------------------------------------------------------------------
+      mask = tmask[indice2d+jpi*jpj*firstzt]
+      if ny EQ 1 then mask = reform(mask, nx, ny, /over)
+;-----------------------------------------------------------
+; construction of land containing all points to mask
+;-----------------------------------------------------------
+      if n_elements(valmask) EQ 0 THEN valmask = 1e20
+      landu = where(u GE valmask/10)
+      if landu[0] NE -1 then u[landu] = 0
+      landv = where(v GE valmask/10)
+      if landv[0] NE -1 then v[landv] = 0
+      res = sqrt(u^2+v^2)
+      if NOT keyword_set(key_periodic) OR nx NE jpi then res[0, *, *] = !values.f_nan
+      res[*, 0, *] = !values.f_nan
+      mask = where(mask eq 0)
+      IF mask[0] NE -1 THEN BEGIN
+        coeftps = lindgen(jpt)*nx*ny
+        coeftps = replicate(1, n_elements(mask))#coeftps
+        mask = (temporary(mask))[*]#replicate(1, jpt)
+        mask = temporary(mask[*]) + temporary(coeftps[*])
+        res[temporary(mask)] = valmask
+      ENDIF
+; moyennes en tous genres
+      domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0], (gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
+      if keyword_set(direc) then res = grossemoyenne(res, direc, /nan, boxzoom = boxzoom, /nodomdef)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xyzt
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    date1 NE date2 AND (size(u))[0] EQ 4:BEGIN
+      indice2d = lindgen(jpi, jpj)
+      indice2d = indice2d[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1]
+      indice3d = lindgen(jpi, jpj, jpk)
+      indice3d = indice3d[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt]
+;------------------------------------------------------------
+; 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
+          case (size(u))[3] OF
+            nzt: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
+            jpk:BEGIN
+              if nxu NE nx then $
+                 if indicex[0] EQ firstxu then u = u[0:nx-1, *, firstzt:lastzt, *] ELSE u = u[1: nx, *, firstzt:lastzt, *]
+              IF nxv NE nx THEN $
+                 if indicex[0] EQ firstxv then v = v[0:nx-1, *, firstzt:lastzt, *] ELSE v = v[1: nx, *, firstzt:lastzt, *]
+              IF nyu NE ny THEN $
+                 if indicey[0] EQ firstyu then u = u[*, 0:ny-1, firstzt:lastzt, *] ELSE u = u[*, 1: ny, firstzt:lastzt, *]
+              IF nyv NE ny THEN $
+                 if indicey[0] EQ firstyv then v = v[*, 0:ny-1, firstzt:lastzt, *] ELSE v = v[*, 1: ny, firstzt:lastzt, *]
+            end
+            ELSE: report, (['the third dimension of u (' + strtrim((size(u))[3], 1) $
+                            +') must be equal to nzt (' + strtrim(nzt, 1) + ') or jpk ('+strtrim(jpk, 1)+')'])
+          endcase
+        END
+        (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
+           (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
+          u = u[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt, *]
+          v = v[indicex[0]:indicex[0]+nx-1, indicey[0]:indicey[0]+ny-1, firstzt:lastzt, *]
+        END
+        ELSE: return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      endcase
+;------------------------------------------------------------------
+; construction of u and v at points T
+;-----------------------------------------------------------
+      a = u[0, *, *, *]
+      u = (u+shift(u, 1, 0, 0, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then u[0, *, *, *] = a
+      a = v[*, 0, *, *]
+      v = (v+shift(v, 0, 1, 0, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then v[*, 0, *, *] = a
+;----------------------------------------------------------------------------
+; attribution of the mask and of longitude and latitude arrays
+;----------------------------------------------------------------------------
+      mask = tmask[indice3d]
+      if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
+;-----------------------------------------------------------
+      if n_elements(valmask) EQ 0 THEN valmask = 1e20
+      landu = where(u GE valmask/10)
+      if landu[0] NE -1 then u[landu] = 0
+      landv = where(v GE valmask/10)
+      if landv[0] NE -1 then v[landv] = 0
+      res = sqrt(u^2+v^2)
+      if NOT keyword_set(key_periodic) OR nx NE jpi then res[0, *, *, *] = !values.f_nan
+      res[*, 0, *, *] = !values.f_nan
+      mask = where(mask eq 0)
+      IF mask[0] NE -1 THEN BEGIN
+        coeftps = lindgen(jpt)*nx*ny*nzt
+        coeftps = replicate(1, n_elements(mask))#coeftps
+        mask = (temporary(mask))[*]#replicate(1, jpt)
+        mask = temporary(mask[*]) + temporary(coeftps[*])
+        res[temporary(mask)] = valmask
+      ENDIF
+; All kind of average
+      domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0], (gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
+      if keyword_set(direc) then res = grossemoyenne(res, direc, /nan, boxzoom = boxzoom, /nodomdef)
+    END
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+;xy
+;----------------------------------------------------------------------------
+;----------------------------------------------------------------------------
+    ELSE:BEGIN                  ;xy
+      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[indice2d]
+          v = v[indice2d]
+        END
+        ELSE:return $
+           , report(['We cannot find re-cut on common indexes (see restriction un the header.).' $
+                     , 'To avoid this problem, read the full domain' $
+                     , 'or use the keyword /memeindice in domdef when defining the zoom area'])
+      endcase
+;------------------------------------------------------------------
+; We reshape u and v to make sure that no dimension has been erased
+;------------------------------------------------------------------
+      if ny EQ 1 then begin
+        u = reform(u, nx, ny, /over)
+        v = reform(v, nx, ny, /over)
+      endif
+;------------------------------------------------------------------
+; construction of u and v at points T
+;-----------------------------------------------------------
+      a = u[0, *]
+      u = (u+shift(u, 1, 0))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then u[0, *] = a
+      a = v[*, 0]
+      v = (v+shift(v, 0, 1))/2.
+      if NOT keyword_set(key_periodic) OR nx NE jpi then v[*, 0] = a
+;----------------------------------------------------------------------------
+; attribution of the mask and of longitude and latitude arrays.
+; We recover the complete grid to establish a big mask extent in the four
+; direction to cover pointsfor which a land point has been
+; considerated (make a small drawing)
+;----------------------------------------------------------------------------
+      mask = tmask[indice2d+jpi*jpj*firstzt]
+      if nyt EQ 1 THEN mask = reform(mask, nx, ny, /over)
+;-----------------------------------------------------------
+; construction of land containing all points to mask
+;-----------------------------------------------------------
+      if n_elements(valmask) EQ 0 THEN valmask = 1e20
+      landu = where(u GE valmask/10)
+      if landu[0] NE -1 then u[landu] = 0
+      landv = where(v GE valmask/10)
+      if landv[0] NE -1 then v[landv] = 0
+      res = sqrt(u^2+v^2)
+      if NOT keyword_set(key_periodic) OR nx NE jpi then res[0, *] = !values.f_nan
+      res[*, 0] = !values.f_nan
+      mask = where(mask eq 0)
+      IF mask[0] NE -1 THEN res[mask] = valmask
+; All kind of average
+      domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0], (gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
+      if keyword_set(direc) then res = moyenne(res, direc, /nan, boxzoom = boxzoom, /nodomdef)
+    END
+;----------------------------------------------------------------------------
+  endcase
+;------------------------------------------------------------
+  if keyword_set(key_performance) THEN print, 'time norm', systime(1)-time1
+  return, res
+end
Index: /trunk/SRC/Obsolete/norme.pro
===================================================================
--- /trunk/SRC/Obsolete/norme.pro	(revision 314)
+++ /trunk/SRC/Obsolete/norme.pro	(revision 314)
@@ -0,0 +1,90 @@
+;+
+;
+; @file_comments
+; calculate the norm of vectors field located on Arakawa C-grid
+;
+; @obsolete
+;
+; @param COMPOSANTEU {in}{required}
+; Matrix representing the zonal coordinates (at U/V point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) or a structure readable by
+; <pro>litchamp</pro> and containing a 2D (xy), 3D (xyz or yt), 4D (xyzt) array.
+; Note that the dimension of the array must suit the domain dimension.
+;
+; @param COMPOSANTEV {in}{required}
+; Matrix representing the meridional coordinates (at V/U point) of a field of vectors
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) or a structure readable by
+; <pro>litchamp</pro> and containing a 2D (xy), 3D (xyz or yt), 4D (xyzt) array.
+; Note that the dimension of the array must suit the domain dimension.
+;
+; @keyword BOXZOOM
+; boxzoom on which do the average (by default the domain selected
+; by the last <pro>domdef</pro> done)
+;
+; @keyword DIREC
+; 't' 'x' 'y' 'z' 'xys' 'xz' 'yz' 'xyz' 'xt' 'yt' 'zt' 'xyt'
+;       'xzt' 'yzt' 'xyzt' Direction on which do averages
+;
+; @returns
+; A 2D (xy), 3D (xyz or yt), 4D (xyzt) Array
+;
+; @uses
+; common.pro
+;
+; @restrictions
+; The norm is calculated on points T. To do this calculation, we average
+; field U and V on points T before calculate the norme. At the edge of
+; coast and of domain, we can not calculate fields U and V at points T,
+; that is why these points are at value !values.f_nan.
+;
+; When we calculate on a reduce geographic domain, field U and V have not
+; necessarily the same number of point. In this case, we recut U and V to
+; keep only common points. We profit of this to redo a <pro>domdef</pro> which redefine
+; a geographic domain on which fields U and V are extracted on same points
+;
+; To know what type of array we work with, we  test its size and dates
+; gave by time[0] and time[jpt-1] to know if thee is a temporal dimension.
+; Before to start norme, make sure that time and jpt are defined how
+; they have to!
+;
+; @examples
+; To calculate the average of the norme of streams on all the domain
+; between 0 and 50:
+;      IDL> res=norme(un,vn,boxzoom=[0,50],dir='xyz')
+;
+; @history
+; 9/6/1999 Sebastien Masson (smasson\@lodyc.jussieu.fr) 
+; December 2007: Sebastien Masson, move to obsolete
+;                       
+;
+; @version
+; $Id$
+;
+;-
+;
+FUNCTION norme, composanteu, composantev, BOXZOOM = boxzoom, DIREC = direc, _EXTRA = ex
+;
+  compile_opt idl2, strictarrsubs
+;
+@cm_4mesh
+@cm_4data
+@cm_4cal
+  IF NOT keyword_set(key_forgetold) THEN BEGIN
+@updatenew
+@updatekwd
+  ENDIF
+;------------------------------------------------------------
+  if keyword_set(boxzoom) then BEGIN
+    Case 1 Of
+      N_Elements(Boxzoom) Eq 1:bte = [lon1, lon2, lat1, lat2, 0., boxzoom[0]]
+      N_Elements(Boxzoom) Eq 2:bte = [lon1, lon2, lat1, lat2, boxzoom[0], boxzoom[1]]
+      N_Elements(Boxzoom) Eq 4:bte = [Boxzoom, vert1, vert2]
+      N_Elements(Boxzoom) Eq 5:bte = [Boxzoom[0:3], 0, Boxzoom[4]]
+      N_Elements(Boxzoom) Eq 6:bte = Boxzoom
+      Else: return, report('wrong number of element in Boxzoom')
+    ENDCASE
+    domdef, boxzoom
+  ENDIF
+;------------------------------------------------------------
+  return, norm(composanteu, composantev, DIREC = direc, _EXTRA = ex)
+end
Index: unk/SRC/ToBeReviewed/CALCULS/norme.pro
===================================================================
--- /trunk/SRC/ToBeReviewed/CALCULS/norme.pro	(revision 313)
+++ 	(revision )
@@ -1,451 +1,0 @@
-;+
-;
-; @file_comments
-; calculate the norm of a field of vectors, then make a possible average.
-;   Comment 1: The field of vector can be, 2d:xy, 3d: xyz or xyt,
-; 4d: xyzt
-;   Comment 2:
-; The calculation of the norm is made before the possible spatial or
-; temporal average because the average of the norm is not equal to the
-; norm of averages
-;
-; @categories
-; Calculation
-;
-; @param COMPOSANTEU {in}{required}
-; an 2d, 3d or 4d array
-;
-; @param COMPOSANTEV {in}{required}
-; an 2d, 3d or 4d array
-;
-; @keyword BOXZOOM
-; boxzoom on which do the average (by default the domain selected
-; by the last <pro>domdef</pro> done)
-;
-; @keyword DIREC
-; 't' 'x' 'y' 'z' 'xys' 'xz' 'yz' 'xyz' 'xt' 'yt' 'zt' 'xyt'
-;       'xzt' 'yzt' 'xyzt' Direction on which do averages
-;
-; @returns
-; Array to trace with <pro>plt</pro>, <pro>pltz</pro> or <pro>pltt</pro>.
-;
-; @uses
-; common.pro
-;
-; @restrictions
-; The norm is calculated on points T. To do this calculation, we average
-; field U and V on points T before calculate the norme. At the edge of
-; coast and of domain, we can not calculate fields U and V at points T,
-; that is why these points are at value !values.f_nan.
-;
-; When we calculate on a reduce geographic domain, field U and V have not
-; necessarily the same number of point. In this case, we recut U and V to
-; keep only common points. We profit of this to redo a <pro>domdef</pro> which redefine
-; a geographic domain on which fields U and V are extracted on same points
-;
-; To know what type of array we work with, we  test its size and dates
-; gave by time[0] and time[jpt-1] to know if thee is a temporal dimension.
-; Before to start norme, make sure that time and jpt are defined how
-; they have to!
-;
-; @examples
-; To calculate the average of the norme of streams on all the domain
-; between 0 and 50:
-;      IDL> res=norme(un,vn,boxzoom=[0,50],dir='xyz')
-;
-; @history
-; Sebastien Masson (smasson\@lodyc.jussieu.fr)
-;                       9/6/1999
-;
-; @version
-; $Id$
-;
-;-
-;
-FUNCTION norme, composanteu, composantev, BOXZOOM = boxzoom, DIREC = direc, _EXTRA = ex
-;
-  compile_opt idl2, strictarrsubs
-;
-@cm_4mesh
-@cm_4data
-@cm_4cal
-  IF NOT keyword_set(key_forgetold) THEN BEGIN
-@updatenew
-@updatekwd
-  ENDIF
-;---------------------------------------------------------
-   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 norme is based on Arakawa C-grid.' $
-                     , 'U and V grids must therefore be defined'])
-;
-;------------------------------------------------------------
-  if keyword_set(boxzoom) then BEGIN
-    Case 1 Of
-      N_Elements(Boxzoom) Eq 1:bte = [lon1, lon2, lat1, lat2, 0., boxzoom[0]]
-      N_Elements(Boxzoom) Eq 2:bte = [lon1, lon2, lat1, lat2, boxzoom[0], boxzoom[1]]
-      N_Elements(Boxzoom) Eq 4:bte = [Boxzoom, vert1, vert2]
-      N_Elements(Boxzoom) Eq 5:bte = [Boxzoom[0:3], 0, Boxzoom[4]]
-      N_Elements(Boxzoom) Eq 6:bte = Boxzoom
-      Else: return, report('Mauvaise Definition de Boxzoom')
-    ENDCASE
-    domdef, boxzoom
-  ENDIF
-;------------------------------------------------------------
-   if NOT keyword_set(direc) then direc = 0
-; construction of u and v at points T
-   u = litchamp(composanteu)
-   v = litchamp(composantev)
-   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
-
-   vargrid='T'
-   varname = 'norme '
-   valmask = 1e20
-;
-   grilleu = litchamp(composanteu, /grid)
-   if grilleu EQ '' then grilleu = 'U'
-   grillev = litchamp(composantev, /grid)
-   if grillev EQ '' then grillev = 'V'
-   IF grilleu EQ 'V' AND grillev EQ 'U' THEN inverse = 1
-   IF grilleu EQ 'T' AND grillev EQ 'T' THEN BEGIN
-      interpolle  = 0
-      return, report('cas non code mais facile a faire!')
-   ENDIF ELSE interpolle = 1
-   if keyword_set(inverse) then begin
-      rien = u
-      u = v
-      v = rien
-   endif
-
-
-;------------------------------------------------------------
-; 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]
-         indice3d = lindgen(jpi, jpj, jpk)
-         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
-;------------------------------------------------------------
-; 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
-               case (size(u))[3] OF
-                  nzt: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
-                  jpk:BEGIN
-                     if nxu NE nx then $
-                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt] ELSE u = u[1: nx, *,firstzt:lastzt]
-                     IF nxv NE nx THEN $
-                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt] ELSE v = v[1: nx, *,firstzt:lastzt]
-                     IF nyu NE ny THEN $
-                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt] ELSE u = u[*, 1: ny,firstzt:lastzt]
-                     IF nyv NE ny THEN $
-                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt] ELSE v = v[*, 1: ny,firstzt:lastzt]
-                  end
-                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-               endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
-               u = u[indice3d]
-               v = v[indice3d]
-            END
-            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-         endcase
-;------------------------------------------------------------------
-; We reshape u and v to make sure that no dimension has been erased
-;------------------------------------------------------------------
-         if nzt EQ 1 then begin
-            u = reform(u, nx, ny, nzt, /over)
-            v = reform(v, nx, ny, nzt, /over)
-         endif
-;------------------------------------------------------------------
-; construction of u and v at points T
-;-----------------------------------------------------------
-         a=u[0,*,*]
-         u=(u+shift(u,1,0,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*]=a
-         a=v[*,0,*]
-         v=(v+shift(v,0,1,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*]=a
-;----------------------------------------------------------------------------
-; attribution of the mask and of longitude and latitude arrays
-;----------------------------------------------------------------------------
-         mask = tmask[indice3d]
-         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
-;-----------------------------------------------------------
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         landu = where(u GE valmask/10)
-         if landu[0] NE -1 then u[landu] = 0
-         landv = where(v GE valmask/10)
-         if landv[0] NE -1 then v[landv] = 0
-         res=sqrt(u^2+v^2)
-         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *]=!values.f_nan
-         res[*,0, *]=!values.f_nan
-         mask = where(mask eq 0)
-         IF mask[0] NE -1 THEN res[mask] = valmask
-; All kind of average
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
-         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
-;-----------------------------------------------------------
-      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:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-         endcase
-;------------------------------------------------------------------
-; construction of u and v at points T
-;-----------------------------------------------------------
-         a=u[0,*,*]
-         u=(u+shift(u,1,0,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*]=a
-         a=v[*,0,*]
-         v=(v+shift(v,0,1,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*]=a
-;----------------------------------------------------------------------------
-; attribution of the mask and of longitude and latitude arrays.
-; We recover the complete grid to establish a big mask extent in the four
-; direction to cover pointsfor which a land point has been
-; considerated (make a small drawing)
-;----------------------------------------------------------------------------
-         mask = tmask[indice2d+jpi*jpj*firstzt]
-         if ny EQ 1 then mask = reform(mask, nx, ny, /over)
-;-----------------------------------------------------------
-; construction of land containing all points to mask
-;-----------------------------------------------------------
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         landu = where(u GE valmask/10)
-         if landu[0] NE -1 then u[landu] = 0
-         landv = where(v GE valmask/10)
-         if landv[0] NE -1 then v[landv] = 0
-         res=sqrt(u^2+v^2)
-         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *]=!values.f_nan
-         res[*,0, *]=!values.f_nan
-         mask = where(mask eq 0)
-         IF mask[0] NE -1 THEN BEGIN
-            coeftps = lindgen(jpt)*nx*ny
-            coeftps = replicate(1, n_elements(mask))#coeftps
-            mask = (temporary(mask))[*]#replicate(1, jpt)
-            mask =temporary(mask[*]) + temporary(coeftps[*])
-            res[temporary(mask)] = valmask
-         ENDIF
-; moyennes en tous genres
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
-         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xyzt
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      date1 NE date2 AND (size(u))[0] EQ 4:BEGIN
-         indice2d = lindgen(jpi, jpj)
-         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
-         indice3d = lindgen(jpi, jpj, jpk)
-         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
-;------------------------------------------------------------
-; 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
-               case (size(u))[3] OF
-                  nzt: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
-                  jpk:BEGIN
-                     if nxu NE nx then $
-                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt,*] ELSE u = u[1: nx, *,firstzt:lastzt,*]
-                     IF nxv NE nx THEN $
-                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt,*] ELSE v = v[1: nx, *,firstzt:lastzt,*]
-                     IF nyu NE ny THEN $
-                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt,*] ELSE u = u[*, 1: ny,firstzt:lastzt,*]
-                     IF nyv NE ny THEN $
-                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt,*] ELSE v = v[*, 1: ny,firstzt:lastzt,*]
-                  end
-                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-               endcase
-            END
-            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
-             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
-               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
-               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
-            END
-            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-         endcase
-;------------------------------------------------------------------
-; construction of u and v at points T
-;-----------------------------------------------------------
-         a=u[0,*,*,*]
-         u=(u+shift(u,1,0,0,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*,*]=a
-         a=v[*,0,*,*]
-         v=(v+shift(v,0,1,0,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*,*]=a
-;----------------------------------------------------------------------------
-; attribution of the mask and of longitude and latitude arrays
-;----------------------------------------------------------------------------
-         mask = tmask[indice3d]
-         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
-;-----------------------------------------------------------
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         landu = where(u GE valmask/10)
-         if landu[0] NE -1 then u[landu] = 0
-         landv = where(v GE valmask/10)
-         if landv[0] NE -1 then v[landv] = 0
-         res=sqrt(u^2+v^2)
-         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *, *]=!values.f_nan
-         res[*,0, *, *]=!values.f_nan
-         mask = where(mask eq 0)
-         IF mask[0] NE -1 THEN BEGIN
-            coeftps = lindgen(jpt)*nx*ny*nzt
-            coeftps = replicate(1, n_elements(mask))#coeftps
-            mask = (temporary(mask))[*]#replicate(1, jpt)
-            mask =temporary(mask[*]) + temporary(coeftps[*])
-            res[temporary(mask)] = valmask
-         ENDIF
-; All kind of average
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
-         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
-      END
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-;xy
-;----------------------------------------------------------------------------
-;----------------------------------------------------------------------------
-      ELSE:BEGIN                ;xy
-         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[indice2d]
-               v = v[indice2d]
-            END
-            ELSE:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
-         endcase
-;------------------------------------------------------------------
-; We reshape u and v to make sure that no dimension has been erased
-;------------------------------------------------------------------
-         if ny EQ 1 then begin
-            u = reform(u, nx, ny, /over)
-            v = reform(v, nx, ny, /over)
-         endif
-;------------------------------------------------------------------
-; construction of u and v at points T
-;-----------------------------------------------------------
-         a=u[0,*]
-         u=(u+shift(u,1,0))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*]=a
-         a=v[*,0]
-         v=(v+shift(v,0,1))/2.
-         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0]=a
-;----------------------------------------------------------------------------
-; attribution of the mask and of longitude and latitude arrays.
-; We recover the complete grid to establish a big mask extent in the four
-; direction to cover pointsfor which a land point has been
-; considerated (make a small drawing)
-;----------------------------------------------------------------------------
-         mask = tmask[indice2d+jpi*jpj*firstzt]
-         if nyt EQ 1 THEN mask = reform(mask, nx, ny, /over)
-;-----------------------------------------------------------
-; construction of land containing all points to mask
-;-----------------------------------------------------------
-         if n_elements(valmask) EQ 0 THEN valmask = 1e20
-         landu = where(u GE valmask/10)
-         if landu[0] NE -1 then u[landu] = 0
-         landv = where(v GE valmask/10)
-         if landv[0] NE -1 then v[landv] = 0
-         res=sqrt(u^2+v^2)
-         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*]=!values.f_nan
-         res[*,0]=!values.f_nan
-         mask = where(mask eq 0)
-         IF mask[0] NE -1 THEN res[mask] = valmask
-; All kind of average
-         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
-         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
-      END
-;----------------------------------------------------------------------------
-   endcase
-;------------------------------------------------------------
-   if keyword_set(key_performance) THEN print, 'temps norme', systime(1)-tempsun
-   return, res
-end
