source: trunk/SRC/Obsolete/meshlec.pro @ 97

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

start to modify headers of Obsolete *.pro files for better idldoc output

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 14.8 KB
Line 
1;------------------------------------------------------------
2;------------------------------------------------------------
3;------------------------------------------------------------
4;+
5;
6; @file_comments lecture du mask de des sorties d'OPA. les sources se trouvent ds les
7; repertoires sur maia du type:
8;  /nom_exp/RESTARTS
9;
10; @obsolete
11; @examples
12; IDL> meshmask[,' nomfich']
13;
14; @param  nomfich {in}{required} string, c''est le nom du fichier a lire. Par defaut, c''est meshmask.
15;
16;
17;  @keyword GLAMBOUNDARY {in} un vecteur de 2 elements specifaint le min et le
18;    max qui doivent etre imposes en longitude (obligatoire si le
19;    tableau depasse 360 degres)
20;
21; @keyword /pasblabla {in} pour suprimer les blablas
22;
23; @keyword /DOUBLE {in} pour forcer a lire les tableaux en double
24;        precision. ce Mot clef est maintenant active
25;        automatiquement.
26;
27; @uses common.pro
28;
29; @restrictions
30;  La definition de ixminmesh,ixmaxmesh,iyminmesh,iymaxmesh
31;  ,izminmesh,izmaxmesh doit etre faite avant l''entree dans cette
32;  routine. pour attribuer automatiquement ces valeurs au maximum
33;  possible les mettre toutes a -1 et meshlec les calculera.
34;
35; @history Sebastien Masson (smasson\@lodyc.jussieu.fr)
36;     Marina Levy : lecture en double precision (cas calcul sur shine)
37;-
38;------------------------------------------------------------
39;------------------------------------------------------------
40;------------------------------------------------------------
41pro meshlec, nomfich, PASBLABLA = pasblabla,DOUBLE = double, GLAMBOUNDARY = glamboundary, GETDIMENSIONS = GETDIMENSIONS
42@common
43   tempsun = systime(1)         ; pour key_performance
44;---------------------------------------------------------
45;---------------------------------------------------------
46   jpiglo      = 0L
47   jpjglo      = 0L
48   jpkglo      = 0L
49   tab    = 'aaaaa'
50;---------------------------------------------------------
51;  definition du domaine de la grille surlequel sont
52;  effectuees les sorties.les indices des tableaux commencant
53;  a 1: cf le fichier wrivr2.F ds WKOPA sur le cray
54;----------------------------------------------------------
55; LECTURE DU MASK trouve ds les fichiers restart
56;----------------------------------------------------------
57;----------------------------------------------------------
58;
59;----------------------------------------------------------
60;       constitution de l'adresse s_fichier et
61;     ouverture du fichier a l'adresse s_fichier
62;-------------------------------------------------------
63
64   IF n_params() EQ 0 then nomfich = 'meshmask'
65   s_fichier = isafile(file = nomfich, iodir = iodir)
66   if not keyword_set(pasblabla) then print, ' '
67   if not keyword_set(pasblabla) then print,'adresse du fichier: ',s_fichier
68
69   openr, numlec, s_fichier, /get_lun, /f77_unformatted, /swap_if_little_endian
70   filepamameters = fstat(numlec)
71;---------------------------------------------------------
72;     lecture
73;---------------------------------------------------------
74   readu, numlec, jpiglo, jpjglo, jpkglo
75   if not keyword_set(pasblabla) then print, 'taille de la grille d''origine: ',jpiglo,', ',jpjglo,', ',jpkglo
76;
77   if keyword_set(getdimensions) then begin
78      free_lun,numlec
79      close, numlec
80      return
81   endif
82
83;---------------------------------------------------------
84; on determine si le fichier a ete ecrit en double precision on non...
85   sizenumber = 8l
86   sizefile8 = 4l+3l*4l+4l $
87    +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
88    +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
89    +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
90    +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
91    +4l*(4l+jpiglo*jpjglo*jpkglo*sizenumber+4l) $
92    +1l*(4l+jpiglo*jpjglo*sizenumber+4l) $
93    +4l*(4l+jpkglo*sizenumber+4l)
94   if filepamameters.size GE sizefile8 THEN double = 1
95;    sizenumber = 4l
96;    sizefile4 = 4l+3l*4l+4l $
97;     +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
98;     +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
99;     +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
100;     +4l*(4l+jpiglo*jpjglo*sizenumber+4l) $
101;     +4l*(4l+jpiglo*jpjglo*jpkglo*sizenumber+4l) $
102;     +1l*(4l+jpiglo*jpjglo*sizenumber+4l) $
103;     +4l*(4l+jpkglo*sizenumber+4l)
104; print, filepamameters.size,  sizefile4,  sizefile8
105;    case filepamameters.size of
106;       sizefile8:double = 1
107;       sizefile4:double = 0
108;       ELSE:BEGIN
109;          nothing = report('The OPA Mesh file as not the good size!')
110;          free_lun,numlec
111;          close, numlec
112;          return
113;       END
114;    endcase
115;
116;
117   if n_elements(ixminmesh) EQ 0 THEN ixminmesh = 0
118   if n_elements(ixmaxmesh) EQ 0 then ixmaxmesh = jpiglo-1
119   if ixminmesh EQ -1 THEN ixminmesh = 0
120   IF ixmaxmesh EQ -1 then ixmaxmesh = jpiglo-1
121   if n_elements(iyminmesh) EQ 0 THEN iyminmesh = 0
122   IF n_elements(iymaxmesh) EQ 0 then iymaxmesh = jpjglo-1
123   if iyminmesh EQ -1 THEN iyminmesh = 0
124   IF iymaxmesh EQ -1 then iymaxmesh = jpjglo-1
125   if n_elements(izminmesh) EQ 0 THEN izminmesh = 0
126   IF n_elements(izmaxmesh) EQ 0 then izmaxmesh = jpkglo-1
127   if izminmesh EQ -1 THEN izminmesh = 0
128   IF izmaxmesh EQ -1 then izmaxmesh = jpkglo-1
129;
130;
131
132   jpi    = long(ixmaxmesh-ixminmesh+1)
133   jpj    = long(iymaxmesh-iyminmesh+1)
134   jpk    = long(izmaxmesh-izminmesh+1)
135;-------------------------------------------------------
136; doit-on reellement lire la grille?
137;-------------------------------------------------------
138   meshparameters = {jpiglo:jpiglo, jpjglo:jpjglo, jpkglo:jpkglo $
139                     , jpi:jpi, jpj:jpj, jpk:jpk $
140                     , ixminmesh:ixminmesh, ixmaxmesh:ixmaxmesh $
141                     , iyminmesh:iyminmesh, iymaxmesh:iymaxmesh $
142                     , izminmesh:izminmesh, izmaxmesh:izmaxmesh $
143                     , key_shift:key_shift}
144;
145;-------------------------------------------------------
146   noticebase=xnotice('Lecture du fichier !C '+s_fichier+'!C ...')
147;-------------------------------------------------------
148
149   IF NOT keyword_set(double) THEN BEGIN
150      z3d    = fltarr(jpiglo, jpjglo, jpkglo)
151      z2d    = fltarr(jpiglo, jpjglo)
152      z1d    = fltarr(jpkglo)
153   ENDIF ELSE BEGIN
154      z3d    = dblarr(jpiglo, jpjglo, jpkglo)
155      z2d    = dblarr(jpiglo, jpjglo)
156      z1d    = dblarr(jpkglo)
157   ENDELSE
158
159   if not keyword_set(pasblabla) then print, ' '
160   readu, numlec, tab,z2d
161   GLAMT=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
162   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GLAMT[25,31]: ',GLAMT[25,31]
163   readu, numlec, tab,z2d
164   GLAMU=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
165   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GLAMU[25,31]: ',GLAMU[25,31]
166   readu, numlec, tab,z2d
167   GLAMV=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
168   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GLAMV[25,31]: ',GLAMV[25,31]
169   readu, numlec, tab,z2d
170   GLAMF=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
171   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GLAMF[25,31]: ',z2d[25,31]
172
173   if not keyword_set(pasblabla) then print, ' '
174   readu, numlec, tab, z2d 
175   GPHIT=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
176   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GPHIT[25,31]: ',GPHIT[25,31]
177   readu, numlec, tab, z2d 
178   GPHIU=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
179   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GPHIU[25,31]: ',GPHIU[25,31]
180   readu, numlec, tab, z2d 
181   GPHIV=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
182   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GPHIV[25,31]: ',GPHIV[25,31]
183   readu, numlec, tab, z2d 
184   GPHIF=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
185   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GPHIF[25,31]: ',z2d[25,31]
186
187   if not keyword_set(pasblabla) then print, ' '
188   readu, numlec, tab, z2d
189   E1T=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
190   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E1T[25,5]: ',z2d[25,5]
191   readu, numlec, tab, z2d
192   E1U=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
193   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E1U[25,5]: ',z2d[25,5]
194   readu, numlec, tab, z2d
195   E1V=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
196   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E1V[25,5]: ',z2d[25,5]
197   readu, numlec, tab, z2d
198   E1F=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
199   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E1F[25,5]: ',z2d[25,5]
200
201   if not keyword_set(pasblabla) then print, ' '
202   readu, numlec, tab, z2d
203   E2T=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
204   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E2T[25,5]: ',z2d[25,5]
205   readu, numlec, tab, z2d
206   E2U=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
207   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E2U[25,5]: ',z2d[25,5]
208   readu, numlec, tab, z2d
209   E2V=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
210   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E2V[25,5]: ',z2d[25,5]
211   readu, numlec, tab, z2d
212   E2F=float(z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh])
213   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E2F[25,5]: ',z2d[25,5]
214
215   if not keyword_set(pasblabla) then print, ' '
216   readu, numlec, tab, z3d     
217   TMASK=byte(z3d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh,izminmesh:izmaxmesh])
218   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur TMASK[25,5,0]: ',TMASK[25,5,0]
219   readu, numlec, tab, z3d   
220   UMASKred=byte(z3d[ixmaxmesh,iyminmesh:iymaxmesh,izminmesh:izmaxmesh])
221   umaskred = reform(umaskred)
222   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur UMASK[25,5,0]: ',z3d[25,5,0]
223   readu, numlec, tab, z3d   
224   VMASKred=byte(z3d[ixminmesh:ixmaxmesh,iymaxmesh,izminmesh:izmaxmesh])
225   vmaskred = reform(vmaskred)
226   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur VMASK[25,5,0]: ',z3d[25,5,0]
227   readu, numlec, tab, z3d   
228   fmaskredy=byte(z3d[ixmaxmesh,iyminmesh:iymaxmesh,izminmesh:izmaxmesh])
229   coast = where(fmaskredy NE 0 and fmaskredy NE 1)
230   IF coast[0] NE -1 THEN fmaskredy[coast] = 0b
231   fmaskredx=byte(z3d[ixminmesh:ixmaxmesh,iymaxmesh,izminmesh:izmaxmesh])
232   coast = where(fmaskredx NE 0 and fmaskredx NE 1)
233   IF coast[0] NE -1 THEN fmaskredx[coast] = 0b
234   fmaskredx = reform(fmaskredx)
235   fmaskredy = reform(fmaskredy)
236   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur FMASK[25,5,0]: ',z3d[25,5,0]
237;
238   if not keyword_set(pasblabla) then print, ' '
239   readu, numlec, tab, z2d   
240;FF=z2d[ixminmesh:ixmaxmesh,iyminmesh:iymaxmesh]
241   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur FF[25,5]: ',z2d[25,5]
242   readu, numlec, tab, z1d
243   GDEPT=float(z1d[izminmesh:izmaxmesh])
244   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GDEPT[1]: ',GDEPT[1]
245   readu, numlec, tab, z1d   
246   GDEPW=float(z1d[izminmesh:izmaxmesh])
247   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur GDEPW[1]: ',GDEPW[1]
248   readu, numlec, tab, z1d   
249   E3T=float(z1d[izminmesh:izmaxmesh])
250   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E3T[3]: ',E3T[3]
251   readu, numlec, tab, z1d   
252   E3W=float(z1d[izminmesh:izmaxmesh])
253   if not keyword_set(pasblabla) then print, 'tableau: ',tab,' exemple de valeur E3W[3]: ',E3W[3]
254;-------------------------------------------------------
255   free_lun,numlec
256   close, numlec
257;-------------------------------------------------------
258;-------------------------------------------------------
259; bornes de glam qui ne doivent pas depasser 360 degres...
260;-------------------------------------------------------
261;    minglam = min(glamt, max = maxglam)
262;    if maxglam-minglam GE 360 AND NOT keyword_set(glamboundary) then $
263;     nothing = execute('glamboundary = '+xquestion('What are the longitudes boundary?', '[-180,180]',/chkwidget))
264;-------------------------------------------------------
265   if keyword_set(glamboundary) then begin
266      if glamboundary[0] NE glamboundary[1] then begin
267         glamt = glamt MOD 360
268         smaller = where(glamt LT glamboundary[0])
269         if smaller[0] NE -1 then glamt[smaller] = glamt[smaller]+360
270         bigger = where(glamt GE glamboundary[1])
271         if bigger[0] NE -1 then glamt[bigger] = glamt[bigger]-360
272         glamu = glamu MOD 360
273         smaller = where(glamu LT glamboundary[0])
274         if smaller[0] NE -1 then glamu[smaller] = glamu[smaller]+360
275         bigger = where(glamu GE glamboundary[1])
276         if bigger[0] NE -1 then glamu[bigger] = glamu[bigger]-360
277         glamv = glamv MOD 360
278         smaller = where(glamv LT glamboundary[0])
279         if smaller[0] NE -1 then glamv[smaller] = glamv[smaller]+360
280         bigger = where(glamv GE glamboundary[1])
281         if bigger[0] NE -1 then glamv[bigger] = glamv[bigger]-360
282         glamf = glamf MOD 360
283         smaller = where(glamf LT glamboundary[0])
284         if smaller[0] NE -1 then glamf[smaller] = glamf[smaller]+360
285         bigger = where(glamf GE glamboundary[1])
286         if bigger[0] NE -1 then glamf[bigger] = glamf[bigger]-360
287      endif
288   endif
289;-------------------------------------------------------
290; shift en x
291;-------------------------------------------------------
292   if keyword_set(key_shift) AND jpi NE 1 then begin
293      glamt = shift(glamt,key_shift , 0)
294      gphit = shift(gphit,key_shift , 0)
295      e1t = shift(e1t, key_shift, 0)
296      e2t = shift(e2t,key_shift , 0)
297      glamu = shift(glamu, key_shift, 0)
298      gphiu = shift(gphiu, key_shift, 0)
299      e1u = shift(e1u,key_shift , 0)
300      e2u = shift(e2u, key_shift, 0)
301      glamv = shift(glamv, key_shift, 0)
302      gphiv = shift(gphiv, key_shift, 0)
303      e1v = shift(e1v,key_shift , 0)
304      e2v = shift(e2v, key_shift, 0)
305      glamf = shift(glamf, key_shift, 0)
306      gphif = shift(gphif, key_shift, 0)
307      e1f = shift(e1f, key_shift , 0)
308      e2f = shift(e2f, key_shift, 0)
309      if jpk EQ 1 then begin
310         tmask = shift(tmask, key_shift, 0)
311         vmaskred = shift(vmaskred, key_shift)
312         fmaskredx = shift(fmaskredx, key_shift)
313      ENDIF ELSE BEGIN
314         tmask = shift(tmask, key_shift, 0,0)
315         vmaskred = shift(vmaskred, key_shift, 0)
316         fmaskredx = shift(fmaskredx, key_shift, 0)
317      ENDELSE
318   endif
319;-------------------------------------------------------
320   key_yreverse = 0
321   key_zreverse = 0
322   key_partialstep = 0
323   key_stride = [1, 1, 1]
324   key_gridtype = 'c'
325;
326   if not keyword_set(pasblabla) then print,'lecture '+nomfich+' finie'
327   widget_control, noticebase, bad_id = toto, /destroy
328   if keyword_set(key_performance) THEN print, 'temps meshlec', systime(1)-tempsun
329
330   return
331end
332
333
334
335
Note: See TracBrowser for help on using the repository browser.