source: trunk/SRC/ToBeReviewed/CALCULS/norme.pro @ 114

Last change on this file since 114 was 114, checked in by smasson, 18 years ago

new compilation options (compile_opt idl2, strictarrsubs) in each routine
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 22.4 KB
Line 
1;------------------------------------------------------------
2;------------------------------------------------------------
3;------------------------------------------------------------
4;+
5; NAME:norme
6;
7; PURPOSE: calcule la norme d''un champ de vecteurs, puis fait une
8; moyenne eventuelle.
9;          Rq1: le champ de vecteur peut etre, 2d:xy, 3d: xyz ou xyt,
10;          4d: xyzt
11;          Rq2: le calcul de la norme est fait avant l''eventuelle
12;          moyenne spatiale ou temporelle car la moyenne de la norme
13;          n''est pas egale a la norme des moyennes.
14;
15; CATEGORY: calcul de post traitement
16;
17; CALLING SEQUENCE:res=norme(champ_de_vecteurs)
18;
19; INPUTS:un tableau 2d, 3d ou 4d
20;
21; KEYWORD PARAMETERS:
22;
23;       BOXZOOM: boxzoom sur laquelle moyenner (par defaut le domaine
24;       selectionner par le dernier domdef effectue)
25;
26;       DIREC:'t' 'x' 'y' 'z' 'xy' 'xz' 'yz' 'xyz' 'xt' 'yt' 'zt' 'xyt'
27;       'xzt' 'yzt' 'xyzt' directions selon lesquelles effectuer les
28;       moyennes
29;
30; OUTPUTS:tableau a tracer avec plt, pltz ou pltt.
31;
32; COMMON BLOCKS:
33;       common.pro
34;
35; SIDE EFFECTS:
36;
37;      La norme est calculee aux points T. Pour faire ce calcul, on
38;      moyenne les champs U et V aux points T avant de calculer la norme.
39;      Au bord des cotes et du domaine, on ne peut pas calculer les
40;      champs U et V aux points T, ces points sont donc a la valeur
41;      !values.f_nan.
42;
43;      lorsqu''on fait le calcul sur un domaine geographique reduit,
44;      les champs U et V ne comprennent pas forcement le meme nombre
45;      de points. Dans ce cas on redecoupe U et V pour ne garder que
46;      les points en commun. Au passage on refait un domdef qui
47;      redefinit un domaine geographique sur lequel les champs U et V
48;      sont extraits sur les meme points.
49;
50; RESTRICTIONS:
51;
52;      pour savoir a quel type de tableau on a a faire, on teste la
53;      taille de celui-ci et les dates donnees par time[0] et
54;      time[jpt-1] pour savoir si il y a une dimension
55;      temporelle. Avant de lancer norme s''assurer que time et jpt
56;      sont bien definis comme il faut!
57;
58; EXAMPLE:
59;
60;     pour calculer la moyenne de la norme des courants sur tout le
61;     dommaine entre 0 et 50:
62;      IDL> res=norme(un,vn,boxzoom=[0,50],dir='xyz')
63;
64; MODIFICATION HISTORY: Sebastien Masson (smasson@lodyc.jussieu.fr)
65;                       9/6/1999
66;-
67;------------------------------------------------------------
68;------------------------------------------------------------
69;------------------------------------------------------------
70FUNCTION norme, composanteu, composantev, BOXZOOM = boxzoom, DIREC = direc, _extra = ex
71;---------------------------------------------------------
72;
73  compile_opt idl2, strictarrsubs
74;
75@cm_4mesh
76@cm_4data
77@cm_4cal
78  IF NOT keyword_set(key_forgetold) THEN BEGIN
79@updatenew
80@updatekwd
81  ENDIF
82;---------------------------------------------------------
83   tempsun = systime(1)         ; pour key_performance
84;
85   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
86     return, report(['This version of norme is based on Arakawa C-grid.' $
87                     , 'U and V grids must therefore be defined'])
88;
89;------------------------------------------------------------
90  if keyword_set(boxzoom) then BEGIN
91    Case 1 Of
92      N_Elements(Boxzoom) Eq 1:bte = [lon1, lon2, lat1, lat2, 0., boxzoom[0]]
93      N_Elements(Boxzoom) Eq 2:bte = [lon1, lon2, lat1, lat2, boxzoom[0], boxzoom[1]]
94      N_Elements(Boxzoom) Eq 4:bte = [Boxzoom, vert1, vert2]
95      N_Elements(Boxzoom) Eq 5:bte = [Boxzoom[0:3], 0, Boxzoom[4]]
96      N_Elements(Boxzoom) Eq 6:bte = Boxzoom
97      Else: return, report('Mauvaise Definition de Boxzoom')
98    ENDCASE
99    domdef, boxzoom
100  ENDIF
101;------------------------------------------------------------
102   if NOT keyword_set(direc) then direc = 0
103; construction de u et v aux pts T
104   u = litchamp(composanteu)
105   v = litchamp(composantev)
106   date1 = time[0]
107   if n_elements(jpt) EQ 0 then date2 = date1 ELSE date2 = time[jpt-1]
108
109   if (size(u))[0] NE (size(v))[0] then return,  -1
110
111   vargrid='T'
112   varname = 'norme '
113   valmask = 1e20
114;
115   grilleu = litchamp(composanteu, /grid)
116   if grilleu EQ '' then grilleu = 'U'
117   grillev = litchamp(composantev, /grid)
118   if grillev EQ '' then grillev = 'V'
119   IF grilleu EQ 'V' AND grillev EQ 'U' THEN inverse = 1
120   IF grilleu EQ 'T' AND grillev EQ 'T' THEN BEGIN
121      interpolle  = 0
122      return, report('cas non code mais facile a faire!')
123   ENDIF ELSE interpolle = 1
124   if keyword_set(inverse) then begin
125      rien = u
126      u = v
127      v = rien
128   endif
129
130
131;------------------------------------------------------------
132; on trouve les points que u et v ont en communs
133;------------------------------------------------------------
134   indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
135   indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
136   indicex = inter(indicexu, indicexv)
137   indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
138   indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
139   indicey = inter(indiceyu, indiceyv)
140   nx = n_elements(indicex)
141   ny = n_elements(indicey)
142;----------------------------------------------------------------------------
143   case 1 of
144;----------------------------------------------------------------------------
145;----------------------------------------------------------------------------
146;xyz
147;----------------------------------------------------------------------------
148;----------------------------------------------------------------------------
149      (size(u))[0] EQ 3 AND date1 EQ date2 :BEGIN
150;----------------------------------------------------------------------------
151         indice2d = lindgen(jpi, jpj)
152         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
153         indice3d = lindgen(jpi, jpj, jpk)
154         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
155;------------------------------------------------------------
156; extraction de u et v sur le domaine qui convient
157;------------------------------------------------------------
158         case 1 of
159            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
160             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
161               case (size(u))[3] OF
162                  nzt:BEGIN
163                     if nxu NE nx then $
164                      if indicex[0] EQ firstxu then u = u[0:nx-1,*,*] ELSE u = u[1: nx, *,*]
165                     IF nxv NE nx THEN $
166                      if indicex[0] EQ firstxv then v = v[0:nx-1,*,*] ELSE v = v[1: nx, *,*]
167                     IF nyu NE ny THEN $
168                      if indicey[0] EQ firstyu then u = u[*,0:ny-1,*] ELSE u = u[*, 1: ny,*]
169                     IF nyv NE ny THEN $
170                      if indicey[0] EQ firstyv then v = v[*,0:ny-1,*] ELSE v = v[*, 1: ny,*]
171                  end
172                  jpk:BEGIN
173                     if nxu NE nx then $
174                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt] ELSE u = u[1: nx, *,firstzt:lastzt]
175                     IF nxv NE nx THEN $
176                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt] ELSE v = v[1: nx, *,firstzt:lastzt]
177                     IF nyu NE ny THEN $
178                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt] ELSE u = u[*, 1: ny,firstzt:lastzt]
179                     IF nyv NE ny THEN $
180                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt] ELSE v = v[*, 1: ny,firstzt:lastzt]
181                  end
182                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
183               endcase
184            END
185            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
186             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
187               u = u[indice3d]
188               v = v[indice3d]
189            END
190            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
191         endcase
192;------------------------------------------------------------------
193; on reform u et v pour s'assurer qu'aucune dimension n'a ete ecrasee
194;------------------------------------------------------------------
195         if nzt EQ 1 then begin
196            u = reform(u, nx, ny, nzt, /over)
197            v = reform(v, nx, ny, nzt, /over)
198         endif
199;------------------------------------------------------------------
200; construction de u et v aux pts T
201;-----------------------------------------------------------
202         a=u[0,*,*]
203         u=(u+shift(u,1,0,0))/2.
204         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*]=a
205         a=v[*,0,*]
206         v=(v+shift(v,0,1,0))/2.
207         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*]=a
208;----------------------------------------------------------------------------
209; attribution du mask et des tableau de longitude et latitude
210;----------------------------------------------------------------------------
211         mask = tmask[indice3d]
212         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
213;-----------------------------------------------------------
214         if n_elements(valmask) EQ 0 THEN valmask = 1e20
215         landu = where(u GE valmask/10)
216         if landu[0] NE -1 then u[landu] = 0
217         landv = where(v GE valmask/10)
218         if landv[0] NE -1 then v[landv] = 0
219         res=sqrt(u^2+v^2)
220         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *]=!values.f_nan
221         res[*,0, *]=!values.f_nan
222         mask = where(mask eq 0)
223         IF mask[0] NE -1 THEN res[mask] = valmask
224; moyennes en tous genres
225         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
226         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
227;-----------------------------------------------------------
228      END
229;----------------------------------------------------------------------------
230;----------------------------------------------------------------------------
231;xyt
232;----------------------------------------------------------------------------
233;----------------------------------------------------------------------------
234      date1 NE date2 AND (size(u))[0] EQ 3 :BEGIN
235         indice2d = lindgen(jpi, jpj)
236         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
237;------------------------------------------------------------
238; extraction de u et v sur le domaine qui convient
239;------------------------------------------------------------
240         case 1 of
241            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
242             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
243               if nxu NE nx then $
244                if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
245               IF nxv NE nx THEN $
246                if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
247               IF nyu NE ny THEN $
248                if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
249               IF nyv NE ny THEN $
250                if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
251            END
252            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
253             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN
254               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
255               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
256            END
257            ELSE:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
258         endcase
259;------------------------------------------------------------------
260; construction de u et v aux pts T
261;-----------------------------------------------------------
262         a=u[0,*,*]
263         u=(u+shift(u,1,0,0))/2.
264         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*]=a
265         a=v[*,0,*]
266         v=(v+shift(v,0,1,0))/2.
267         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*]=a
268;----------------------------------------------------------------------------
269; attribution du mask et des tableau de longitude et latitude
270; on recupere la grille complette pour etablir un grand mask etendu ds les 4
271; directions pour couvrir les points pour lesquels un pt terre a ete pris en
272; compte (faire un petit dessin...)
273;----------------------------------------------------------------------------
274         mask = tmask[indice2d+jpi*jpj*firstzt]
275         if ny EQ 1 then mask = reform(mask, nx, ny, /over)
276;-----------------------------------------------------------
277; construction de terre qui contient tous les point a masquer
278;-----------------------------------------------------------
279         if n_elements(valmask) EQ 0 THEN valmask = 1e20
280         landu = where(u GE valmask/10)
281         if landu[0] NE -1 then u[landu] = 0
282         landv = where(v GE valmask/10)
283         if landv[0] NE -1 then v[landv] = 0
284         res=sqrt(u^2+v^2)
285         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *]=!values.f_nan
286         res[*,0, *]=!values.f_nan
287         mask = where(mask eq 0)
288         IF mask[0] NE -1 THEN BEGIN
289            coeftps = lindgen(jpt)*nx*ny
290            coeftps = replicate(1, n_elements(mask))#coeftps
291            mask = (temporary(mask))[*]#replicate(1, jpt)
292            mask =temporary(mask[*]) + temporary(coeftps[*])
293            res[temporary(mask)] = valmask
294         ENDIF
295; moyennes en tous genres
296         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
297         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
298      END
299;----------------------------------------------------------------------------
300;----------------------------------------------------------------------------
301;xyzt
302;----------------------------------------------------------------------------
303;----------------------------------------------------------------------------
304      date1 NE date2 AND (size(u))[0] EQ 4:BEGIN
305         indice2d = lindgen(jpi, jpj)
306         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
307         indice3d = lindgen(jpi, jpj, jpk)
308         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
309;------------------------------------------------------------
310; extraction de u et v sur le domaine qui convient
311;------------------------------------------------------------
312         case 1 of
313            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
314             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
315               case (size(u))[3] OF
316                  nzt:BEGIN
317                     if nxu NE nx then $
318                      if indicex[0] EQ firstxu then u = u[0:nx-1,*,*,*] ELSE u = u[1: nx, *,*,*]
319                     IF nxv NE nx THEN $
320                      if indicex[0] EQ firstxv then v = v[0:nx-1,*,*,*] ELSE v = v[1: nx, *,*,*]
321                     IF nyu NE ny THEN $
322                      if indicey[0] EQ firstyu then u = u[*,0:ny-1,*,*] ELSE u = u[*, 1: ny,*,*]
323                     IF nyv NE ny THEN $
324                      if indicey[0] EQ firstyv then v = v[*,0:ny-1,*,*] ELSE v = v[*, 1: ny,*,*]
325                  end
326                  jpk:BEGIN
327                     if nxu NE nx then $
328                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt,*] ELSE u = u[1: nx, *,firstzt:lastzt,*]
329                     IF nxv NE nx THEN $
330                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt,*] ELSE v = v[1: nx, *,firstzt:lastzt,*]
331                     IF nyu NE ny THEN $
332                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt,*] ELSE u = u[*, 1: ny,firstzt:lastzt,*]
333                     IF nyv NE ny THEN $
334                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt,*] ELSE v = v[*, 1: ny,firstzt:lastzt,*]
335                  end
336                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
337               endcase
338            END
339            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
340             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
341               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
342               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
343            END
344            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
345         endcase
346;------------------------------------------------------------------
347; construction de u et v aux pts T
348;-----------------------------------------------------------
349         a=u[0,*,*,*]
350         u=(u+shift(u,1,0,0,0))/2.
351         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*,*,*]=a
352         a=v[*,0,*,*]
353         v=(v+shift(v,0,1,0,0))/2.
354         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0,*,*]=a
355;----------------------------------------------------------------------------
356; attribution du mask et des tableau de longitude et latitude
357;----------------------------------------------------------------------------
358         mask = tmask[indice3d]
359         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
360;-----------------------------------------------------------
361         if n_elements(valmask) EQ 0 THEN valmask = 1e20
362         landu = where(u GE valmask/10)
363         if landu[0] NE -1 then u[landu] = 0
364         landv = where(v GE valmask/10)
365         if landv[0] NE -1 then v[landv] = 0
366         res=sqrt(u^2+v^2)
367         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*, *, *]=!values.f_nan
368         res[*,0, *, *]=!values.f_nan
369         mask = where(mask eq 0)
370         IF mask[0] NE -1 THEN BEGIN
371            coeftps = lindgen(jpt)*nx*ny*nzt
372            coeftps = replicate(1, n_elements(mask))#coeftps
373            mask = (temporary(mask))[*]#replicate(1, jpt)
374            mask =temporary(mask[*]) + temporary(coeftps[*])
375            res[temporary(mask)] = valmask
376         ENDIF
377; moyennes en tous genres
378         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
379         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
380      END
381;----------------------------------------------------------------------------
382;----------------------------------------------------------------------------
383;xy
384;----------------------------------------------------------------------------
385;----------------------------------------------------------------------------
386      ELSE:BEGIN                ;xy
387         indice2d = lindgen(jpi, jpj)
388         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
389;------------------------------------------------------------
390; extraction de u et v sur le domaine qui convient
391;------------------------------------------------------------
392         case 1 of
393            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
394             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
395               if nxu NE nx then $
396                if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *]
397               IF nxv NE nx THEN $
398                if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
399               IF nyu NE ny THEN $
400                if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny]
401               IF nyv NE ny THEN $
402                if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
403            END
404            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
405             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN
406               u = u[indice2d]
407               v = v[indice2d]
408            END
409            ELSE:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
410         endcase
411;------------------------------------------------------------------
412; on reform u et v pour s'assurer qu'aucune dimension n'a ete ecrasee
413;------------------------------------------------------------------
414         if ny EQ 1 then begin
415            u = reform(u, nx, ny, /over)
416            v = reform(v, nx, ny, /over)
417         endif
418;------------------------------------------------------------------
419; construction de u et v aux pts T
420;-----------------------------------------------------------
421         a=u[0,*]
422         u=(u+shift(u,1,0))/2.
423         if NOT keyword_set(key_periodic) OR nx NE jpi then u[0,*]=a
424         a=v[*,0]
425         v=(v+shift(v,0,1))/2.
426         if NOT keyword_set(key_periodic) OR nx NE jpi then v[*,0]=a
427;----------------------------------------------------------------------------
428; attribution du mask et des tableau de longitude et latitude
429; on recupere la grille complette pour etablir un grand mask etendu ds les 4
430; directions pour couvrir les points pour lesquels un pt terre a ete pris en
431; compte (faire un petit dessin...)
432;----------------------------------------------------------------------------
433         mask = tmask[indice2d+jpi*jpj*firstzt]
434         if nyt EQ 1 THEN mask = reform(mask, nx, ny, /over)
435;-----------------------------------------------------------
436; construction de terre qui contient tous les point a masquer
437;-----------------------------------------------------------
438         if n_elements(valmask) EQ 0 THEN valmask = 1e20
439         landu = where(u GE valmask/10)
440         if landu[0] NE -1 then u[landu] = 0
441         landv = where(v GE valmask/10)
442         if landv[0] NE -1 then v[landv] = 0
443         res=sqrt(u^2+v^2)
444         if NOT keyword_set(key_periodic) OR nx NE jpi then res[0,*]=!values.f_nan
445         res[*,0]=!values.f_nan
446         mask = where(mask eq 0)
447         IF mask[0] NE -1 THEN res[mask] = valmask
448; moyennes en tous genres
449         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
450         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
451      END
452;----------------------------------------------------------------------------
453   endcase
454;------------------------------------------------------------
455   if keyword_set(key_performance) THEN print, 'temps norme', systime(1)-tempsun
456   return, res
457end
Note: See TracBrowser for help on using the repository browser.