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

Last change on this file since 72 was 25, checked in by pinsard, 18 years ago

upgrade of CALCULS according to cerbere.lodyc.jussieu.fr: /usr/home/smasson/SAXO_RD/ : files

  • 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@cm_4mesh
73@cm_4data
74@cm_4cal
75  IF NOT keyword_set(key_forgetold) THEN BEGIN
76@updatenew
77@updatekwd
78  ENDIF
79;---------------------------------------------------------
80   tempsun = systime(1)         ; pour key_performance
81;
82   IF finite(glamu[0])*finite(gphiu[0])*finite(glamv[0])*finite(gphiv[0]) EQ 0 THEN $
83     return, report(['This version of norme is based on Arakawa C-grid.' $
84                     , 'U and V grids must therefore be defined'])
85;
86;------------------------------------------------------------
87  if keyword_set(boxzoom) then BEGIN
88    Case 1 Of
89      N_Elements(Boxzoom) Eq 1:bte = [lon1, lon2, lat1, lat2, 0., boxzoom[0]]
90      N_Elements(Boxzoom) Eq 2:bte = [lon1, lon2, lat1, lat2, boxzoom[0], boxzoom[1]]
91      N_Elements(Boxzoom) Eq 4:bte = [Boxzoom, vert1, vert2]
92      N_Elements(Boxzoom) Eq 5:bte = [Boxzoom[0:3], 0, Boxzoom[4]]
93      N_Elements(Boxzoom) Eq 6:bte = Boxzoom
94      Else: return, report('Mauvaise Definition de Boxzoom')
95    ENDCASE
96    domdef, boxzoom
97  ENDIF
98;------------------------------------------------------------
99   if NOT keyword_set(direc) then direc = 0
100; construction de u et v aux pts T
101   u = litchamp(composanteu)
102   v = litchamp(composantev)
103   date1 = time[0]
104   if n_elements(jpt) EQ 0 then date2 = date1 ELSE date2 = time[jpt-1]
105
106   if (size(u))[0] NE (size(v))[0] then return,  -1
107
108   vargrid='T'
109   varname = 'norme '
110   valmask = 1e20
111;
112   grilleu = litchamp(composanteu, /grid)
113   if grilleu EQ '' then grilleu = 'U'
114   grillev = litchamp(composantev, /grid)
115   if grillev EQ '' then grillev = 'V'
116   IF grilleu EQ 'V' AND grillev EQ 'U' THEN inverse = 1
117   IF grilleu EQ 'T' AND grillev EQ 'T' THEN BEGIN
118      interpolle  = 0
119      return, report('cas non code mais facile a faire!')
120   ENDIF ELSE interpolle = 1
121   if keyword_set(inverse) then begin
122      rien = u
123      u = v
124      v = rien
125   endif
126
127
128;------------------------------------------------------------
129; on trouve les points que u et v ont en communs
130;------------------------------------------------------------
131   indicexu = (lindgen(jpi))[firstxu:firstxu+nxu-1]
132   indicexv = (lindgen(jpi))[firstxv:firstxv+nxv-1]
133   indicex = inter(indicexu, indicexv)
134   indiceyu = (lindgen(jpj))[firstyu:firstyu+nyu-1]
135   indiceyv = (lindgen(jpj))[firstyv:firstyv+nyv-1]
136   indicey = inter(indiceyu, indiceyv)
137   nx = n_elements(indicex)
138   ny = n_elements(indicey)
139;----------------------------------------------------------------------------
140   case 1 of
141;----------------------------------------------------------------------------
142;----------------------------------------------------------------------------
143;xyz
144;----------------------------------------------------------------------------
145;----------------------------------------------------------------------------
146      (size(u))[0] EQ 3 AND date1 EQ date2 :BEGIN
147;----------------------------------------------------------------------------
148         indice2d = lindgen(jpi, jpj)
149         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
150         indice3d = lindgen(jpi, jpj, jpk)
151         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
152;------------------------------------------------------------
153; extraction de u et v sur le domaine qui convient
154;------------------------------------------------------------
155         case 1 of
156            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
157             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
158               case (size(u))[3] OF
159                  nzt:BEGIN
160                     if nxu NE nx then $
161                      if indicex[0] EQ firstxu then u = u[0:nx-1,*,*] ELSE u = u[1: nx, *,*]
162                     IF nxv NE nx THEN $
163                      if indicex[0] EQ firstxv then v = v[0:nx-1,*,*] ELSE v = v[1: nx, *,*]
164                     IF nyu NE ny THEN $
165                      if indicey[0] EQ firstyu then u = u[*,0:ny-1,*] ELSE u = u[*, 1: ny,*]
166                     IF nyv NE ny THEN $
167                      if indicey[0] EQ firstyv then v = v[*,0:ny-1,*] ELSE v = v[*, 1: ny,*]
168                  end
169                  jpk:BEGIN
170                     if nxu NE nx then $
171                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt] ELSE u = u[1: nx, *,firstzt:lastzt]
172                     IF nxv NE nx THEN $
173                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt] ELSE v = v[1: nx, *,firstzt:lastzt]
174                     IF nyu NE ny THEN $
175                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt] ELSE u = u[*, 1: ny,firstzt:lastzt]
176                     IF nyv NE ny THEN $
177                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt] ELSE v = v[*, 1: ny,firstzt:lastzt]
178                  end
179                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
180               endcase
181            END
182            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
183             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
184               u = u[indice3d]
185               v = v[indice3d]
186            END
187            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
188         endcase
189;------------------------------------------------------------------
190; on reform u et v pour s'assurer qu'aucune dimension n'a ete ecrasee
191;------------------------------------------------------------------
192         if nzt EQ 1 then begin
193            u = reform(u, nx, ny, nzt, /over)
194            v = reform(v, nx, ny, nzt, /over)
195         endif
196;------------------------------------------------------------------
197; construction de u et v aux pts T
198;-----------------------------------------------------------
199         a=u(0,*,*)
200         u=(u+shift(u,1,0,0))/2.
201         if NOT keyword_set(key_periodic) OR nx NE jpi then u(0,*,*)=a
202         a=v(*,0,*)
203         v=(v+shift(v,0,1,0))/2.
204         if NOT keyword_set(key_periodic) OR nx NE jpi then v(*,0,*)=a
205;----------------------------------------------------------------------------
206; attribution du mask et des tableau de longitude et latitude
207;----------------------------------------------------------------------------
208         mask = tmask[indice3d]
209         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
210;-----------------------------------------------------------
211         if n_elements(valmask) EQ 0 THEN valmask = 1e20
212         landu = where(u GE valmask/10)
213         if landu[0] NE -1 then u[landu] = 0
214         landv = where(v GE valmask/10)
215         if landv[0] NE -1 then v[landv] = 0
216         res=sqrt(u^2+v^2)
217         if NOT keyword_set(key_periodic) OR nx NE jpi then res(0,*, *)=!values.f_nan
218         res(*,0, *)=!values.f_nan
219         mask = where(mask eq 0)
220         IF mask[0] NE -1 THEN res(mask) = valmask
221; moyennes en tous genres
222         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
223         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
224;-----------------------------------------------------------
225      END
226;----------------------------------------------------------------------------
227;----------------------------------------------------------------------------
228;xyt
229;----------------------------------------------------------------------------
230;----------------------------------------------------------------------------
231      date1 NE date2 AND (size(u))[0] EQ 3 :BEGIN
232         indice2d = lindgen(jpi, jpj)
233         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
234;------------------------------------------------------------
235; extraction de u et v sur le domaine qui convient
236;------------------------------------------------------------
237         case 1 of
238            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
239             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
240               if nxu NE nx then $
241                if indicex[0] EQ firstxu then u = u[0:nx-1, *, *] ELSE u = u[1: nx, *, *]
242               IF nxv NE nx THEN $
243                if indicex[0] EQ firstxv then v = v[0:nx-1, *, *] ELSE v = v[1: nx, *, *]
244               IF nyu NE ny THEN $
245                if indicey[0] EQ firstyu then u = u[*, 0:ny-1, *] ELSE u = u[*, 1: ny, *]
246               IF nyv NE ny THEN $
247                if indicey[0] EQ firstyv then v = v[*, 0:ny-1, *] ELSE v = v[*, 1: ny, *]
248            END
249            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
250             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN
251               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
252               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1, *]
253            END
254            ELSE:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
255         endcase
256;------------------------------------------------------------------
257; construction de u et v aux pts T
258;-----------------------------------------------------------
259         a=u(0,*,*)
260         u=(u+shift(u,1,0,0))/2.
261         if NOT keyword_set(key_periodic) OR nx NE jpi then u(0,*,*)=a
262         a=v(*,0,*)
263         v=(v+shift(v,0,1,0))/2.
264         if NOT keyword_set(key_periodic) OR nx NE jpi then v(*,0,*)=a
265;----------------------------------------------------------------------------
266; attribution du mask et des tableau de longitude et latitude
267; on recupere la grille complette pour etablir un grand mask etendu ds les 4
268; directions pour couvrir les points pour lesquels un pt terre a ete pris en
269; compte (faire un petit dessin...)
270;----------------------------------------------------------------------------
271         mask = tmask[indice2d+jpi*jpj*firstzt]
272         if ny EQ 1 then mask = reform(mask, nx, ny, /over)
273;-----------------------------------------------------------
274; construction de terre qui contient tous les point a masquer
275;-----------------------------------------------------------
276         if n_elements(valmask) EQ 0 THEN valmask = 1e20
277         landu = where(u GE valmask/10)
278         if landu[0] NE -1 then u[landu] = 0
279         landv = where(v GE valmask/10)
280         if landv[0] NE -1 then v[landv] = 0
281         res=sqrt(u^2+v^2)
282         if NOT keyword_set(key_periodic) OR nx NE jpi then res(0,*, *)=!values.f_nan
283         res(*,0, *)=!values.f_nan
284         mask = where(mask eq 0)
285         IF mask[0] NE -1 THEN BEGIN
286            coeftps = lindgen(jpt)*nx*ny
287            coeftps = replicate(1, n_elements(mask))#coeftps
288            mask = (temporary(mask))[*]#replicate(1, jpt)
289            mask =temporary(mask[*]) + temporary(coeftps[*])
290            res(temporary(mask)) = valmask
291         ENDIF
292; moyennes en tous genres
293         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
294         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
295      END
296;----------------------------------------------------------------------------
297;----------------------------------------------------------------------------
298;xyzt
299;----------------------------------------------------------------------------
300;----------------------------------------------------------------------------
301      date1 NE date2 AND (size(u))[0] EQ 4:BEGIN
302         indice2d = lindgen(jpi, jpj)
303         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
304         indice3d = lindgen(jpi, jpj, jpk)
305         indice3d = indice3d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt]
306;------------------------------------------------------------
307; extraction de u et v sur le domaine qui convient
308;------------------------------------------------------------
309         case 1 of
310            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
311             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
312               case (size(u))[3] OF
313                  nzt:BEGIN
314                     if nxu NE nx then $
315                      if indicex[0] EQ firstxu then u = u[0:nx-1,*,*,*] ELSE u = u[1: nx, *,*,*]
316                     IF nxv NE nx THEN $
317                      if indicex[0] EQ firstxv then v = v[0:nx-1,*,*,*] ELSE v = v[1: nx, *,*,*]
318                     IF nyu NE ny THEN $
319                      if indicey[0] EQ firstyu then u = u[*,0:ny-1,*,*] ELSE u = u[*, 1: ny,*,*]
320                     IF nyv NE ny THEN $
321                      if indicey[0] EQ firstyv then v = v[*,0:ny-1,*,*] ELSE v = v[*, 1: ny,*,*]
322                  end
323                  jpk:BEGIN
324                     if nxu NE nx then $
325                      if indicex[0] EQ firstxu then u = u[0:nx-1, *,firstzt:lastzt,*] ELSE u = u[1: nx, *,firstzt:lastzt,*]
326                     IF nxv NE nx THEN $
327                      if indicex[0] EQ firstxv then v = v[0:nx-1, *,firstzt:lastzt,*] ELSE v = v[1: nx, *,firstzt:lastzt,*]
328                     IF nyu NE ny THEN $
329                      if indicey[0] EQ firstyu then u = u[*, 0:ny-1,firstzt:lastzt,*] ELSE u = u[*, 1: ny,firstzt:lastzt,*]
330                     IF nyv NE ny THEN $
331                      if indicey[0] EQ firstyv then v = v[*, 0:ny-1,firstzt:lastzt,*] ELSE v = v[*, 1: ny,firstzt:lastzt,*]
332                  end
333                  ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
334               endcase
335            END
336            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND (size(u))[3] EQ jpk AND $
337             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj AND (size(u))[3] EQ jpk :BEGIN
338               u = u[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
339               v = v[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1,firstzt:lastzt, *]
340            END
341            ELSE: return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
342         endcase
343;------------------------------------------------------------------
344; construction de u et v aux pts T
345;-----------------------------------------------------------
346         a=u(0,*,*,*)
347         u=(u+shift(u,1,0,0,0))/2.
348         if NOT keyword_set(key_periodic) OR nx NE jpi then u(0,*,*,*)=a
349         a=v(*,0,*,*)
350         v=(v+shift(v,0,1,0,0))/2.
351         if NOT keyword_set(key_periodic) OR nx NE jpi then v(*,0,*,*)=a
352;----------------------------------------------------------------------------
353; attribution du mask et des tableau de longitude et latitude
354;----------------------------------------------------------------------------
355         mask = tmask[indice3d]
356         if nzt EQ 1 then mask = reform(mask, nx, ny, nzt, /over)
357;-----------------------------------------------------------
358         if n_elements(valmask) EQ 0 THEN valmask = 1e20
359         landu = where(u GE valmask/10)
360         if landu[0] NE -1 then u[landu] = 0
361         landv = where(v GE valmask/10)
362         if landv[0] NE -1 then v[landv] = 0
363         res=sqrt(u^2+v^2)
364         if NOT keyword_set(key_periodic) OR nx NE jpi then res(0,*, *, *)=!values.f_nan
365         res(*,0, *, *)=!values.f_nan
366         mask = where(mask eq 0)
367         IF mask[0] NE -1 THEN BEGIN
368            coeftps = lindgen(jpt)*nx*ny*nzt
369            coeftps = replicate(1, n_elements(mask))#coeftps
370            mask = (temporary(mask))[*]#replicate(1, jpt)
371            mask =temporary(mask[*]) + temporary(coeftps[*])
372            res(temporary(mask)) = valmask
373         ENDIF
374; moyennes en tous genres
375         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
376         if keyword_set(direc) then res = grossemoyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
377      END
378;----------------------------------------------------------------------------
379;----------------------------------------------------------------------------
380;xy
381;----------------------------------------------------------------------------
382;----------------------------------------------------------------------------
383      ELSE:BEGIN                ;xy
384         indice2d = lindgen(jpi, jpj)
385         indice2d = indice2d[indicex[0]:indicex[0]+nx-1,indicey[0]:indicey[0]+ny-1]
386;------------------------------------------------------------
387; extraction de u et v sur le domaine qui convient
388;------------------------------------------------------------
389         case 1 of
390            (size(u))[1] EQ nxu AND (size(u))[2] EQ nyu AND $
391             (size(v))[1] EQ nxv AND (size(v))[2] EQ nyv:BEGIN
392               if nxu NE nx then $
393                if indicex[0] EQ firstxu then u = u[0:nx-1, *] ELSE u = u[1: nx, *]
394               IF nxv NE nx THEN $
395                if indicex[0] EQ firstxv then v = v[0:nx-1, *] ELSE v = v[1: nx, *]
396               IF nyu NE ny THEN $
397                if indicey[0] EQ firstyu then u = u[*, 0:ny-1] ELSE u = u[*, 1: ny]
398               IF nyv NE ny THEN $
399                if indicey[0] EQ firstyv then v = v[*, 0:ny-1] ELSE v = v[*, 1: ny]
400            END
401            (size(u))[1] EQ jpi AND (size(u))[2] EQ jpj AND $
402             (size(v))[1] EQ jpi AND (size(v))[2] EQ jpj:BEGIN
403               u = u[indice2d]
404               v = v[indice2d]
405            END
406            ELSE:return, report('problemes d''adequation entre la taille du domaine et la taille des matrices necessaires a tracer des vecteurs')
407         endcase
408;------------------------------------------------------------------
409; on reform u et v pour s'assurer qu'aucune dimension n'a ete ecrasee
410;------------------------------------------------------------------
411         if ny EQ 1 then begin
412            u = reform(u, nx, ny, /over)
413            v = reform(v, nx, ny, /over)
414         endif
415;------------------------------------------------------------------
416; construction de u et v aux pts T
417;-----------------------------------------------------------
418         a=u(0,*)
419         u=(u+shift(u,1,0))/2.
420         if NOT keyword_set(key_periodic) OR nx NE jpi then u(0,*)=a
421         a=v(*,0)
422         v=(v+shift(v,0,1))/2.
423         if NOT keyword_set(key_periodic) OR nx NE jpi then v(*,0)=a
424;----------------------------------------------------------------------------
425; attribution du mask et des tableau de longitude et latitude
426; on recupere la grille complette pour etablir un grand mask etendu ds les 4
427; directions pour couvrir les points pour lesquels un pt terre a ete pris en
428; compte (faire un petit dessin...)
429;----------------------------------------------------------------------------
430         mask = tmask[indice2d+jpi*jpj*firstzt]
431         if nyt EQ 1 THEN mask = reform(mask, nx, ny, /over)
432;-----------------------------------------------------------
433; construction de terre qui contient tous les point a masquer
434;-----------------------------------------------------------
435         if n_elements(valmask) EQ 0 THEN valmask = 1e20
436         landu = where(u GE valmask/10)
437         if landu[0] NE -1 then u[landu] = 0
438         landv = where(v GE valmask/10)
439         if landv[0] NE -1 then v[landv] = 0
440         res=sqrt(u^2+v^2)
441         if NOT keyword_set(key_periodic) OR nx NE jpi then res(0,*)=!values.f_nan
442         res(*,0)=!values.f_nan
443         mask = where(mask eq 0)
444         IF mask[0] NE -1 THEN res(mask) = valmask
445; moyennes en tous genres
446         domdef, (glamt[indice2d])[0, 0], (glamu[indice2d])[nx-1, 0],(gphit[indice2d])[0, 0], (gphiv[indice2d])[0, ny-1], vert1, vert2, /meme
447         if keyword_set(direc) then res = moyenne(res,direc,/nan, boxzoom = boxzoom, /nodomdef)
448      END
449;----------------------------------------------------------------------------
450   endcase
451;------------------------------------------------------------
452   if keyword_set(key_performance) THEN print, 'temps norme', systime(1)-tempsun
453   return, res
454end
Note: See TracBrowser for help on using the repository browser.