source: trunk/SRC/ToBeReviewed/HOPE/computehopegrid.pro @ 226

Last change on this file since 226 was 226, checked in by pinsard, 17 years ago

corrections of some misspellings in some *.pro
  • Property svn:keywords set to Id
File size: 6.6 KB
Line 
1;------------------------------------------------------------
2;------------------------------------------------------------
3;------------------------------------------------------------
4;+
5;
6; @file_comments
7;
8; @categories
9; Grid
10;
11; @param XAXIS
12;
13; @param YAXIS
14;
15; @param ZAXIS
16;
17; @param LINETYPE
18;
19; @keyword FORTHEMASK
20;
21; @keyword WPOINT
22;
23; @keyword FIRSTS
24;
25; @keyword LASTS
26;
27; @keyword PTTYPE
28;
29; @uses
30; common.pro
31;
32; @history
33; Sebastien Masson (smasson\@lodyc.jussieu.fr)
34;                      2001-06
35;
36; @version
37; $Id$
38;
39; @todo seb
40;
41;-
42;------------------------------------------------------------
43;------------------------------------------------------------
44;------------------------------------------------------------
45PRO computehopegrid, xaxis, yaxis, zaxis, linetype, FORTHEMASK = forthemask, WPOINT = wpoint, FIRSTS = firsts, LASTS = lasts, PTTYPE = pttype
46;---------------------------------------------------------
47;
48  compile_opt idl2, strictarrsubs
49;
50@cm_4mesh
51@cm_4data
52  IF NOT keyword_set(key_forgetold) THEN BEGIN
53@updatenew
54  ENDIF
55;---------------------------------------------------------
56;------------------------------------------------------------
57   if not (keyword_set(scalar)*keyword_set(vector)) then scalar=1
58;------------------------------------------------------------
59   jpiglo = n_elements(xaxis)
60   jpjglo = n_elements(yaxis)
61   jpkglo = n_elements(zaxis)
62;
63   jpi = jpiglo
64   jpj = jpjglo
65   jpk = jpkglo
66;
67   if NOT keyword_set(firsts) then firsts = [0, 0, 0]
68   if NOT keyword_set(lasts) then lasts = [jpi-1, jpj-1, jpk-1]
69;
70; determination of the grid type and of the point type
71;
72   if keyword_set(pttype) then vargrid = pttype
73   if linetype EQ 'odd-even' then key_gridtype = 'e' ELSE key_gridtype = 'c'
74;
75; computation of the horizontal grid
76;
77   if key_gridtype EQ 'e' then begin
78      if vargrid EQ 'T' then begin
79         glamt=xaxis
80         glamt = glamt#replicate(1, jpj)
81         xstep=(glamt[1,0]-glamt[0,0])/2.
82         glamt[*,2*lindgen((jpj)/2)]=glamt[*,2*lindgen(jpj/2)]-xstep
83         glamu=glamt+xstep
84      ENDIF ELSE BEGIN
85         glamu=xaxis
86         glamu = glamu#replicate(1, jpj)
87         xstep=(glamu[1,0]-glamu[0,0])/2.
88         glamu[*,2*lindgen((jpj)/2)]=glamu[*,2*lindgen(jpj/2)]-xstep
89         glamt=glamu-xstep
90      ENDELSE
91;zoom
92      glamt = glamt[firsts[0]:lasts[0], firsts[1]:lasts[1]]
93      glamu = glamu[firsts[0]:lasts[0], firsts[1]:lasts[1]]
94      jpiglo = lasts[0]-firsts[0]+1
95      jpi = jpiglo
96      jpjglo = lasts[1]-firsts[1]+1
97      jpj = jpjglo
98      glamv = glamu
99      glamf = glamu
100      gphit = yaxis[firsts[1]:lasts[1]]
101      gphit = replicate(1, jpi)#gphit
102      gphif = gphit
103      gphiu = gphit
104      gphiv = gphif
105   ENDIF ELSE BEGIN
106      if vargrid eq 'T' then begin
107         glamt=xaxis
108         glamt = glamt#replicate(1, jpj)
109         glamu=glamt+(glamt[1,0]-glamt[0,0])/2.
110      ENDIF ELSE BEGIN
111         glamu=xaxis
112         glamu = glamu#replicate(1, jpj)
113         xstep=(glamu[1,0]-glamu[0,0])/2.
114         glamt=glamu-(glamu[1,0]-glamu[0,0])/2.
115      ENDELSE
116;zoom
117      glamt = glamt[firsts[0]:lasts[0], firsts[1]:lasts[1]]
118      glamu = glamu[firsts[0]:lasts[0], firsts[1]:lasts[1]]
119      jpiglo = lasts[0]-firsts[0]+1
120      jpi = jpiglo
121      jpjglo = lasts[1]-firsts[1]+1
122      jpj = jpjglo
123      glamv = glamt
124      glamf = glamu
125      gphit = yaxis[firsts[1]:lasts[1]]
126      gphit = replicate(1, jpi)#gphit
127      gphiu = gphit
128      if jpj GT 1 then begin
129         gphiv = gphit+(gphit[0, 1]-gphit[0, 0])/2.
130         gphif = gphit+(gphit[0, 1]-gphit[0, 0])/2.
131      ENDIF ELSE BEGIN
132         gphiv = gphit
133         gphif = gphit
134      ENDELSE
135   ENDELSE
136;
137; computation of the vertical grid
138;
139   if keyword_set(wpoint) then begin
140      gdepw = zaxis
141      if jpk ne 1 then begin
142         gdept=(shift(gdepw, -1)+gdepw)/2.
143         gdept[jpk-1]=gdepw[jpk-1]+.5*(gdepw[jpk-1]-gdepw[jpk-2])
144      endif else gdept=zaxis
145   endif else begin
146      gdept = zaxis
147      if jpk ne 1 then begin
148         gdepw=(shift(gdept, 1)+gdept)/2.
149         gdepw[0]=0
150      endif else gdepw = zaxis
151   endelse
152;
153; computation of the vertical scale factors
154;
155   if jpk ne 1 then begin
156      e3t = abs(shift(gdepw,-1)-gdepw)
157      e3t[jpk-1] = abs(gdept[jpk-1]-gdepw[jpk-1])
158      e3w = abs(gdept-shift(gdept,1))
159      e3w[0] = abs(gdept[0]-gdepw[0])
160   endif else begin
161      e3t=1
162      e3w=1
163   endelse
164; zoom
165   gdept = gdept[firsts[2]:lasts[2]]
166   gdepw = gdepw[firsts[2]:lasts[2]]
167   e3t = e3t[firsts[2]:lasts[2]]
168   e3w = e3w[firsts[2]:lasts[2]]
169   jpkglo = lasts[2]-firsts[2]+1
170   jpk = jpkglo
171;
172; computation of the horizontal scale factors
173;
174   e1t = replicate(1b,jpi,jpj)
175   e2t = replicate(1b,jpi,jpj)
176   e1u = e1t
177   e2u = e2t
178   e1v = e1t
179   e2v = e2t
180   e1f = e1t
181   e2f = e2t
182;
183; mask
184;
185   tmask = replicate(1b, jpi, jpj, jpk)
186   if keyword_set(forthemask) then BEGIN
187      land=where(forthemask ge valmask/10)
188      if land[0] ne -1 then tmask[land] = 0b
189   endif
190   umaskred = replicate(1, jpj, jpk)
191   vmaskred = replicate(1, jpi, jpk)
192   fmaskredy = replicate(1, jpj, jpk)
193   fmaskredx = replicate(1, jpi, jpk)
194   @updateold
195   domdef
196   if keyword_set(firsts) AND keyword_set(lasts) then BEGIN
197      if vargrid EQ 'T' then BEGIN
198         if jpj GT 1 then begin
199            lon1 = min(glamt[0, [0, 1]])
200            lon2 = max(glamt[jpi-1, [0, 1]])
201         endif ELSE BEGIN
202            lon1 = min(glamt[0, 0])
203            lon2 = max(glamt[jpi-1, [0]])
204         ENDELSE
205      ENDIF ELSE BEGIN
206         if jpj GT 1 then begin
207            lon1 = min(glamu[0, [0, 1]])
208            lon2 = max(glamu[jpi-1, [0, 1]])
209         endif ELSE BEGIN
210            lon1 = min(glamu[0, 0])
211            lon2 = max(glamu[jpi-1, 0])
212         ENDELSE
213      ENDELSE
214      lat1 = min([gphit[0, 0], gphit[0, jpj-1]])
215      lat2 = max([gphit[0, 0], gphit[0, jpj-1]])
216      domdef, lon1, lon2, lat1, lat2, gdepw[0], gdept[jpk-1], gridtype = vargrid
217   ENDIF
218;
219   ixminmesh  =0l
220   ixmaxmesh  =long(jpi-1)
221   iyminmesh  =0l
222   iymaxmesh  =long(jpj-1)
223   izminmesh  =0l
224   izmaxmesh  =long(jpk-1)
225;----------------------------------------------------------
226; for the triangulation
227;----------------------------------------------------------
228   key_periodic = glamt[0] EQ ((glamt[jpi-1]+glamt[+1]-glamt[0]) MOD 360)
229   if jpi gt 4 AND jpj GT 4 then begin
230      triangles_list = triangule(shifted = glamt[0, 0] LT glamt[0, 1])
231      twin_corners_up=-1
232      twin_corners_dn=-1
233   ENDIF ELSE BEGIN
234      triangles_list=-1
235      twin_corners_up=-1
236      twin_corners_dn=-1
237   ENDELSE
238;------------------------------------------------------------
239  IF NOT keyword_set(key_forgetold) THEN BEGIN
240   @updateold
241  ENDIF
242;------------------------------------------------------------
243   return
244end
245
246
Note: See TracBrowser for help on using the repository browser.