source: trunk/SRC/ToBeReviewed/GRILLE/grille.pro @ 163

Last change on this file since 163 was 163, checked in by navarro, 18 years ago

header improvements : type of parameters and keywords, default values, spell checking + idldoc assistant (IDL online_help)
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 16.3 KB
Line 
1;------------------------------------------------------------
2;------------------------------------------------------------
3;------------------------------------------------------------
4;+
5;
6; @file_comments
7; Choose the grid which must be used to do the graph in function of
8; vargrid and send back corresponding parameters calculated in
9; domdef.pro and reduced at the domain defined by domdef (contrarily
10; to grandegrille.pro)
11; BEWARE!! The choice of the grid is made from the value of the
12; global variable vargrid, which can be equal to 'T', 'U', 'V', 'W' ou 'F'.
13;
14; @categories
15;
16; @keyword TRI
17; This keyword serve to obtain, thanks to grille, the triangulation which
18; refer to the grid but only on the part of the zoom. This of triangulation
19; is passed in the variable we have equate at TRI.
20; For example: grille,...,tri=triangulation_reduite.
21; This keyword is used in plt.pro
22;
23; @keyword WDEPTH
24; To specify that the field is at W depth instead of T
25; depth (automatically activated if vargrid eq 'W')
26;
27; @keyword FORPLT
28; In plt, we want that land points, glam and gphi, be equal to glamt and
29; gphit regardless of the grid.
30;
31; @keyword NOTRI
32; Useful only when TRI is activated. In this case, grill send back -1 in the
33; variable tri even if the variable of the common triangles_list is defined
34; and different of-1
35;
36; @keyword _EXTRA
37; Used to pass your keywords
38;
39; @keyword TOUT
40;
41; @param MASK {out}{optional}
42; For the definition, see domdef and the management of subdomains on the web.
43;
44; @param  GLAM {out}{optional}
45; For the definition, see domdef and the management of subdomains on the web.
46;
47; @param  GPHI {out}{optional}
48; For the definition, see domdef and the management of subdomains on the web.
49;
50; @param  GDEP {out}{optional}
51; For the definition, see domdef and the management of subdomains on the web.
52;
53; @param  NX {out}{optional}
54; For the definition, see domdef and the management of subdomains on the web.
55;
56; @param  NY {out}{optional}
57; For the definition, see domdef and the management of subdomains on the web.
58;
59; @param  NZ {out}{optional}
60; For the definition, see domdef and the management of subdomains on the web.
61;
62; @param  FIRSTX {out}{optional}
63; For the definition, see domdef and the management of subdomains on the web.
64;
65; @param  FIRSTY {out}{optional}
66; For the definition, see domdef and the management of subdomains on the web.
67;
68; @param  FIRSTZ {out}{optional}
69; For the definition, see domdef and the management of subdomains on the web.
70;         
71; @param LASTX {out}{optional}
72; For the definition, see domdef and the management of subdomains on the web.
73;
74; @param  LASTY {out}{optional}
75; For the definition, see domdef and the management of subdomains on the web.
76;
77; @param  LASTZ {out}{optional}
78; For the definition, see domdef and the management of subdomains on the web.
79;
80; @param  E1 {out}{optional}
81; For the definition, see domdef and the management of subdomains on the web.
82;
83; @param  E2 {out}{optional}
84; For the definition, see domdef and the management of subdomains on the web.
85;
86; @param  E3 {out}{optional}
87; For the definition, see domdef and the management of subdomains on the web.
88;
89; @uses
90; cm_4mesh
91; cm_4data
92;
93; @restrictions
94; Use the variable vargrid
95;
96; @restrictions
97; Vargrid must be 'T', 'W', 'U', 'V' or 'F'
98;
99; @history
100;  Sebastien Masson (smasson\@lodyc.jussieu.fr)
101;                       12/2/1999
102;                       10/11/1999 /forplt
103;
104; @version
105; $Id$
106;
107; @todo Comment ecrire la remarque sur les inputs?
108;
109;-
110;------------------------------------------------------------
111;------------------------------------------------------------
112;------------------------------------------------------------
113pro grille, mask, glam, gphi, gdep, nx, ny, nz, firstx, firsty, firstz, lastx, lasty, lastz, e1, e2, e3, TRI = tri, NOTRI = notri, TOUT = tout, FORPLT = forplt, IFPLTZ = ifpltz, WDEPTH = wdepth, _EXTRA = ex
114;------------------------------------------------------------
115; include commons
116;
117  compile_opt idl2, strictarrsubs
118;
119@cm_4mesh
120@cm_4data
121  IF NOT keyword_set(key_forgetold) THEN BEGIN
122@updatenew
123  ENDIF
124;---------------------
125  tempsun = systime(1)          ; For key_performance
126;------------------------------------------------------------
127  vargrid = strupcase(strmid(vargrid,0,/reverse_offset))
128;
129  if vargrid eq 'W' then wdepth = 1
130  if keyword_set(tout) then begin
131    savedbox = 1b
132    saveboxparam, 'boxparam4grille.dat'
133    domdef, gridtype = vargrid, _EXTRA = ex
134  endif
135  tempdeux = systime(1)         ; For key_performance =2
136;------------------------------------------------------------
137;------------------------------------------------------------
138  IF keyword_set(wdepth) THEN BEGIN
139    firstz = firstzw
140    lastz = lastzw
141    nz = nzw
142  ENDIF ELSE BEGIN
143    firstz = firstzt
144    lastz = lastzt
145    nz = nzt
146  ENDELSE
147;------------------------------------------------------------
148;------------------------------------------------------------
149  CASE 1 OF
150;------------------------------------------------------------
151; grid T and W
152;------------------------------------------------------------
153    vargrid eq 'T' OR vargrid eq 'W' : begin
154;scalars
155      nx = nxt
156      ny = nyt
157      firstx = firstxt
158      firsty = firstyt
159      lastx = lastxt
160      lasty = lastyt
161;2d vectors
162      IF arg_present(glam) THEN glam = glamt[firstx:lastx, firsty:lasty]
163      IF arg_present(gphi) THEN gphi = gphit[firstx:lastx, firsty:lasty]
164      IF arg_present(e1) THEN e1 = e1t[firstx:lastx, firsty:lasty]
165      IF arg_present(e2) THEN e2 = e2t[firstx:lastx, firsty:lasty]
166;3d vectors
167      IF keyword_set(forplt) THEN mask = tmask[firstx:lastx, firsty:lasty, firstz] $
168      ELSE IF arg_present(mask) THEN mask = tmask[firstx:lastx, firsty:lasty, firstz:lastz]
169    end
170;------------------------------------------------------------
171; grid U
172;------------------------------------------------------------
173    vargrid eq 'U': begin
174;scalars
175      nx = nxu
176      ny = nyu
177      firstx = firstxu
178      firsty = firstyu
179      lastx = lastxu
180      lasty = lastyu
181;2d vectors
182      IF arg_present(glam) THEN glam = glamu[firstx:lastx, firsty:lasty]
183      IF arg_present(gphi) THEN gphi = gphiu[firstx:lastx, firsty:lasty]
184      if keyword_set(forplt) then BEGIN
185        mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
186        eastboarder = mask-shift(mask, 1, 0)*mask
187        westboarder = mask-shift(mask, -1, 0)*mask
188        if key_periodic NE 1 OR nx NE jpi then westboarder[nx-1, *] = 0b
189        tmp1 = shift(eastboarder, 0, 1)
190        tmp1[*, 0] = 0b
191        tmp2 = shift(eastboarder, 0, -1)
192        tmp2[*, ny-1] = 0b
193        add = (temporary(tmp1)+temporary(tmp2))*(1b-eastboarder)*(1b-temporary(westboarder))
194        eastboarder = temporary(eastboarder)+temporary(add)
195        tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
196        tmp1[*, ny-1] = 1b
197        tmp1[*, 0] = 1b
198        tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
199        if key_periodic NE 1 OR nx NE jpi then begin
200          tmp2[nx-1, *] = 1b
201          tmp2[0, *] = 0b
202        endif
203        no1 = temporary(tmp1)*temporary(tmp2)
204        tmp = temporary(eastboarder)*temporary(no1)*mask
205        mask[0:nx-2, *] = 0b
206        tmp = temporary(tmp)+temporary(mask)
207        tmp = where(tmp GE 1)
208        if tmp[0] NE -1 then begin
209          glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
210          gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
211        endif
212      ENDIF
213      IF arg_present(e1) THEN e1  = e1u[firstx:lastx, firsty:lasty]
214      IF arg_present(e2) THEN e2  = e2u[firstx:lastx, firsty:lasty]
215;3d vectors
216      IF keyword_set(forplt) THEN mask = (umask())[firstx:lastx, firsty:lasty, firstz] $
217      ELSE IF arg_present(mask) THEN mask = (umask())[firstx:lastx, firsty:lasty, firstz:lastz]
218    end
219;------------------------------------------------------------
220; grid V
221;------------------------------------------------------------
222    vargrid eq 'OPAPTDHV' or vargrid eq 'OPAPT3DV' $
223      or vargrid eq 'V': begin
224;scalars
225      nx = nxv
226      ny = nyv
227      firstx = firstxv
228      firsty = firstyv
229      lastx = lastxv
230      lasty = lastyv
231;2d vectors
232      IF arg_present(glam) THEN glam = glamv[firstx:lastx, firsty:lasty]
233      IF arg_present(gphi) THEN gphi = gphiv[firstx:lastx, firsty:lasty]
234      if keyword_set(forplt) then BEGIN
235        mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
236        northboarder = mask-shift(mask, 0, 1)*mask
237        southboarder = mask-shift(mask, 0, -1)*mask
238        southboarder[*, ny-1] = 0b
239        tmp1 = shift(northboarder, -1, 0)
240        if key_periodic NE 1 OR nx NE jpi then tmp1[nx-1, *] = 0b
241        tmp2 = shift(northboarder, 1, 0)
242        if key_periodic NE 1 OR nx NE jpi then tmp2[0, *] = 0b
243        add = (temporary(tmp1)+temporary(tmp2))*(1b-northboarder)*(1b-southboarder)
244        northboarder = temporary(northboarder)+temporary(add)
245        tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
246        tmp1[*, ny-1] = 1b
247        tmp1[*, 0] = 0b
248        tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
249        if key_periodic NE 1 OR nx NE jpi then begin
250          tmp2[nx-1, *] = 1b
251          tmp2[0, *] = 1b
252        endif
253        no1 = temporary(tmp1)*temporary(tmp2)
254        tmp = temporary(northboarder)*mask*temporary(no1)
255        mask[*, 0:ny-2] = 0b
256        tmp = temporary(tmp)+temporary(mask)
257        tmp = where(tmp GE 1)
258        if tmp[0] NE -1 then begin
259          glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
260          gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
261        endif
262      ENDIF
263      IF arg_present(e1) THEN e1  = e1v[firstx:lastx, firsty:lasty]
264      IF arg_present(e2) THEN e2  = e2v[firstx:lastx, firsty:lasty]
265;3d vecteurs
266      IF keyword_set(forplt) THEN mask = (vmask())[firstx:lastx, firsty:lasty, firstz] $
267      ELSE IF arg_present(mask) THEN mask = (vmask())[firstx:lastx, firsty:lasty, firstz:lastz]
268    end
269;------------------------------------------------------------
270; grid F
271;------------------------------------------------------------
272    vargrid eq 'OPAPTDHF' or vargrid eq 'OPAPT3DF' $
273      or vargrid eq 'F': begin
274;scalars
275      nx = nxf
276      ny = nyf
277      firstx = firstxf
278      firsty = firstyf
279      lastx = lastxf
280      lasty = lastyf
281;2d vectors
282      IF arg_present(glam) THEN glam = glamf[firstx:lastx, firsty:lasty]
283      IF arg_present(gphi) THEN gphi = gphif[firstx:lastx, firsty:lasty]
284      if keyword_set(forplt) then BEGIN
285        mask = 1b-tmask[firstx:lastx, firsty:lasty, firstz]
286        eastboarder = mask-shift(mask, 1, 0)*mask
287        westboarder = mask-shift(mask, -1, 0)*mask
288        westboarder[nx-1, *] = 0b
289        northboarder = mask-shift(mask, 0, 1)*mask
290        southboarder = mask-shift(mask, 0, -1)*mask
291        southboarder[*, ny-1] = 0b
292        tmp1 = shift(northboarder, -1, 0)
293        if key_periodic NE 1 OR nx NE jpi then tmp1[nx-1, *] = 0b
294        tmp2 = shift(northboarder, 1, 0)
295        if key_periodic NE 1 OR nx NE jpi then tmp2[0, *] = 0b
296        add = (temporary(tmp1)+temporary(tmp2))*(1b-northboarder)*(1b-southboarder)
297        northboarder = temporary(northboarder)+temporary(add)
298        tmp1 = shift(eastboarder, 0, 1)
299        tmp1[*, 0] = 0b
300        tmp2 = shift(eastboarder, 0, -1)
301        tmp2[*, ny-1] = 0b
302        add = (temporary(tmp1)+temporary(tmp2))*(1b-eastboarder)*(1b-temporary(westboarder))
303        eastboarder = temporary(eastboarder)+temporary(add)
304        tmp1 = (mask+shift(mask, 0, -1)+shift(mask, 0, 1)) NE 1b
305        tmp1[*, ny-1] = 1b
306        tmp1[*, 0] = 1b
307        tmp2 = (mask+shift(mask, -1, 0)+shift(mask, 1, 0)) NE 1b
308        if key_periodic NE 1 OR nx NE jpi then begin
309          tmp2[nx-1, *] = 1b
310          tmp2[0, *] = 1b
311        endif
312        no1 = temporary(tmp1)*temporary(tmp2)
313        tmp = (temporary(northboarder)+temporary(eastboarder))*mask*temporary(no1)
314        mask[0:nx-2, *] = 0b
315        mask[*, 0:ny-2] = 0b
316        tmp = temporary(tmp)+temporary(mask)
317        tmp = where(tmp GE 1)
318        if tmp[0] NE -1 then begin
319          glam[tmp] = (glamt[firstx:lastx, firsty:lasty])[tmp]
320          gphi[tmp] = (gphit[firstx:lastx, firsty:lasty])[tmp]
321        endif
322      ENDIF
323      IF arg_present(e1) THEN e1  = e1f[firstx:lastx, firsty:lasty]
324      IF arg_present(e2) THEN e2  = e2f[firstx:lastx, firsty:lasty]
325;3d vectors
326      IF keyword_set(forplt) THEN mask = (fmask())[firstx:lastx, firsty:lasty, firstz] $
327      ELSE IF arg_present(mask) THEN mask = (fmask())[firstx:lastx, firsty:lasty, firstz:lastz]
328    END
329;------------------------------------------------------------
330    ELSE:BEGIN
331      ras = report('Wrong definition of vargrid = '+vargrid+'. Only T, U, V, W or F are acceptable')
332      stop
333    END
334  ENDCASE
335  IF testvar(var = key_performance) EQ 2 THEN $
336    print, 'temps grille: attribution des scalaires, vecteurs et tableaux ', systime(1)-tempdeux
337;
338;------------------------------------------------------------
339;------------------------------------------------------------
340;------------------------------------------------------------
341; Variables refering to the vertical dimension
342;------------------------------------------------------------
343;------------------------------------------------------------
344;------------------------------------------------------------
345;
346;
347  tempdeux = systime(1)         ; For key_performance =2
348  if keyword_set(wdepth) then begin
349    gdep = gdepw[firstz:lastz]
350    e3 = e3w[firstz:lastz]
351  endif else begin
352    gdep = gdept[firstz:lastz]
353    e3 = e3t[firstz:lastz]
354  ENDELSE
355; for the vertical sections with partial steps
356  IF keyword_set(ifpltz) AND keyword_set(key_partialstep) THEN BEGIN
357    CASE 1 OF
358      ifpltz EQ 'xz' AND ny EQ 1:BEGIN
359        bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
360        good = where(bottom NE 0 AND bottom NE nz+keyword_set(wdepth))
361        bottom = lindgen(nx)+(bottom-1l+keyword_set(wdepth))*nx
362        IF good[0] NE -1 THEN BEGIN
363          bottom = bottom[good]
364          IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
365          gdep = replicate(1, nx)#gdep
366          if keyword_set(wdepth) THEN $
367            truegdep = hdepw[firstx:lastx, firsty:lasty] $
368          ELSE truegdep = hdept[firstx:lastx, firsty:lasty]
369          gdep[bottom] = truegdep[good]
370        ENDIF
371      END
372      ifpltz EQ 'yz' AND nx EQ 1:BEGIN
373        bottom = total(tmask[firstx:lastx, firsty:lasty, firstz:lastz], 3)
374        good = where(bottom NE 0 AND bottom NE nz+keyword_set(wdepth))
375        bottom = lindgen(ny)+(bottom-1l+keyword_set(wdepth))*ny
376        IF good[0] NE -1 THEN BEGIN
377          bottom = bottom[good]
378          IF lastz EQ jpk-1 THEN gdep[nz-1] = max(hdepw)
379          gdep = replicate(1, ny)#gdep
380          if keyword_set(wdepth) THEN $
381            truegdep = hdepw[firstx:lastx, firsty:lasty] $
382          ELSE truegdep = hdept[firstx:lastx, firsty:lasty]
383          gdep[bottom] = truegdep[good]
384        ENDIF
385      END
386      ELSE:
387    ENDCASE
388  ENDIF
389  IF testvar(var = key_performance) EQ 2 THEN $
390    print, 'temps grille: Variables se rapportant a la dimension verticale ', systime(1)-tempdeux
391;------------------------------------------------------------
392; Triangulation vector when TRI is activated.
393;------------------------------------------------------------
394  if arg_present(TRI) then $
395    if triangles_list[0] EQ -1 OR keyword_set(notri) then tri = -1 ELSE BEGIN
396    tempdeux = systime(1)       ; pour key_performance =2
397    msk = bytarr(jpi, jpj)
398    msk[firstx:lastx, firsty:lasty] = 1
399    ind = where( msk[triangles_list[0, *]]*msk[triangles_list[1, *]]*msk[triangles_list[2, *]] EQ 1 )
400    tri = triangles_list[*, ind]-(firstx+firsty*jpi)
401    y = tri/jpi
402    x = tri-y*jpi
403    tri = x+y*nx
404    IF testvar(var = key_performance) EQ 2 THEN $
405      print, 'temps grille: decoupage de la triangulation ', systime(1)-tempdeux
406  ENDELSE
407;------------------------------------------------------------------
408; To make sure there is not any degenerated dimension (=1)
409;------------------------------------------------------------------
410;    mask=reform(mask, /over)
411;    glam=reform(glam, /over)
412;    gphi=reform(gphi, /over)
413;    gdep=reform(gdep, /over)
414;    e1=reform(e1, /over)
415;    e2=reform(e2, /over)
416;    e3=reform(e3, /over)
417
418  if keyword_set(savedbox) THEN restoreboxparam, 'boxparam4grille.dat'
419  if keyword_set(key_performance) THEN print, 'temps grille', systime(1)-tempsun
420
421;------------------------------------------------------------
422  IF NOT keyword_set(key_forgetold) THEN BEGIN
423@updateold
424  ENDIF
425;---------------------
426  return
427
428end
429
430
431
432
433
434
435
Note: See TracBrowser for help on using the repository browser.