source: trunk/SRC/Grid/checkperio.pro @ 213

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

add 2 routines + fix minor bugfix
  • Property svn:keywords set to Id
File size: 11.6 KB
Line 
1;+
2; @file_comments
3; check ORCA2 and ORCA05 periodicity (east-west and noth fold) for files related to the T GRID ONLY
4;
5; @categories
6; Grid
7;
8; @param file {in}{required}{type=scalar string}
9; A string giving the name of the NetCdf file
10;
11; @keyword WRITE {default=0}{type=scalar: 0 or 1}
12; Define to 1 to write the good periodicity in the file.
13; By default print only mesage of periodicity errors
14;
15; @restrictions
16; Works only for grid T points
17; if the WRITE keyword is activated, we need the writing access on the file!
18;
19; @examples
20;    IDL> checkperio, '/Users/smasson/OPA/modipsl/bin/runoff_1m_nomask.nc'
21;    checkperio, '/Users/smasson/OPA/modipsl/bin/runoff_1m_nomask.nc'
22;    nav_lon bad y periodicity (1)
23;    nav_lon bad y periodicity (2)
24;    nav_lon bad x periodicity (1)
25;    nav_lon bad x periodicity (2)
26;    nav_lat bad y periodicity (1)
27;    nav_lat bad y periodicity (2)
28;    nav_lat bad x periodicity (1)
29;    nav_lat bad x periodicity (2)
30;    socoefr bad y periodicity (1)
31;    socoefr bad x periodicity (1)
32;    socoefr bad x periodicity (2)
33;    sorunoff bad x periodicity (1)
34;    sorunoff bad x periodicity (2)
35;    IDL> checkperio, '/Users/smasson/OPA/modipsl/bin/runoff_1m_nomask.nc', /WRITE
36;    IDL> checkperio, '/Users/smasson/OPA/modipsl/bin/runoff_1m_nomask.nc'
37;
38; @history
39; October 2006: Sebastien Masson (smasson\@locean-ipsl.upmc.fr)
40;   
41; @version
42; $Id$
43;-
44PRO checkperio, file, WRITE = write
45
46  IF file_test(file) EQ 0 THEN BEGIN
47    print, 'file '+file+' not found'
48    RETURN
49  ENDIF
50
51  cdfid = ncdf_open(file, WRITE = write)
52  finq = ncdf_inquire(cdfid)
53;
54  dimsz = lonarr(finq.ndims)
55  FOR did = 0, finq.ndims-1 DO BEGIN
56    ncdf_diminq, cdfid, did, name, size
57    dimsz[did] = size
58  ENDFOR
59;
60  FOR vid = 0, finq.nvars-1 DO BEGIN
61    vinq =  ncdf_varinq(cdfid, vid)
62    IF vinq.ndims GE 2 THEN BEGIN
63      vdimsz = dimsz[vinq.dim]
64      jpi = vdimsz[0]
65      jpj = vdimsz[1]
66      CASE vinq.ndims OF
67;------------------------------
68; 2D VAR
69;------------------------------
70        2:BEGIN
71; y periodicity (1)
72          CASE 1 OF
73            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
74              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1], offset = [1, jpj-1]
75              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1], offset = [1, jpj-3]
76              part2 = reverse(shift(part2, -1))
77              IF array_equal(part1, part2) NE 1 THEN BEGIN
78                print, vinq.name+' bad y periodicity (1)'
79                IF keyword_set(write) THEN $
80                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1], offset = [1, jpj-1]
81              ENDIF
82            END
83            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
84              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1], offset = [1, jpj-1]
85              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1], offset = [1, jpj-2]
86              part2 = reverse(part2)
87              IF array_equal(part1, part2) NE 1 THEN BEGIN
88                print, vinq.name+' bad y periodicity (1)'
89                IF keyword_set(write) THEN $
90                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1], offset = [1, jpj-1]
91              ENDIF
92            END
93            ELSE:print, vinq.name+' nothing to check'
94          ENDCASE
95; y periodicity (2)
96          CASE 1 OF
97            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
98              ncdf_varget, cdfid, vid, part1, count = [(jpi-1)/2-2+1, 1], offset = [(jpi-1)/2+2, jpj-2]
99              ncdf_varget, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1], offset = [          2, jpj-2]
100              part2 = reverse(part2)
101              IF array_equal(part1, part2) NE 1 THEN BEGIN
102                print, vinq.name+' bad y periodicity (2)'
103                IF keyword_set(write) THEN $
104                   ncdf_varput, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1], offset = [(jpi-1)/2+2, jpj-2]
105              ENDIF
106            END
107            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
108            END
109            ELSE:print, vinq.name+' nothing to check'
110          ENDCASE
111; x periodicity (1)
112          ncdf_varget, cdfid, vid, part1, count = [1, jpj], offset = [0, 0]
113          ncdf_varget, cdfid, vid, part2, count = [1, jpj], offset = [jpi-2, 0]
114          IF array_equal(part1, part2) NE 1 THEN BEGIN
115            print, vinq.name+' bad x periodicity (1)'
116            IF keyword_set(write) THEN $
117               ncdf_varput, cdfid, vid, part2, count = [1, jpj], offset = [0, 0]
118          ENDIF
119; x periodicity (2)
120          ncdf_varget, cdfid, vid, part1, count = [1, jpj], offset = [jpi-1, 0]
121          ncdf_varget, cdfid, vid, part2, count = [1, jpj], offset = [1, 0]
122          IF array_equal(part1, part2) NE 1 THEN BEGIN
123            print, vinq.name+' bad x periodicity (2)'
124            IF keyword_set(write) THEN $
125               ncdf_varput, cdfid, vid, part2, count = [1, jpj], offset = [jpi-1, 0]
126          ENDIF
127        END
128;------------------------------
129; 3D VAR
130;------------------------------
131        3:BEGIN
132          jpk = vdimsz[2]
133; y periodicity (1)
134          CASE 1 OF
135            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
136              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1, jpk], offset = [1, jpj-1, 0]
137              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1, jpk], offset = [1, jpj-3, 0]
138              IF jpk EQ 1 THEN part2 = reform(part2, jpi-2, 1, jpk, /over)
139              part2 = reverse(shift(part2, -1, 0, 0), 1)
140              IF array_equal(part1, part2) NE 1 THEN BEGIN
141                print, vinq.name+' bad y periodicity (1)'
142                IF keyword_set(write) THEN $
143                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1, jpk], offset = [1, jpj-1, 0]
144              ENDIF
145            END
146            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
147              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1, jpk], offset = [1, jpj-1, 0]
148              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1, jpk], offset = [1, jpj-2, 0]
149              IF jpk EQ 1 THEN part2 = reform(part2, jpi-2, 1, jpk, /over)
150              part2 = reverse(part2, 1)
151              IF array_equal(part1, part2) NE 1 THEN BEGIN
152                print, vinq.name+' bad y periodicity (1)'
153                IF keyword_set(write) THEN $
154                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1, jpk], offset = [1, jpj-1, 0]
155              ENDIF
156            END
157            ELSE:print, vinq.name+' nothing to check'
158          ENDCASE
159; y periodicity (2)
160          CASE 1 OF
161            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
162              ncdf_varget, cdfid, vid, part1, count = [(jpi-1)/2-2+1, 1, jpk], offset = [(jpi-1)/2+2, jpj-2, 0]
163              ncdf_varget, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1, jpk], offset = [          2, jpj-2, 0]
164              part2 = reverse(part2, 1)
165              IF array_equal(part1, part2) NE 1 THEN BEGIN
166                print, vinq.name+' bad y periodicity (2)'
167                IF keyword_set(write) THEN $
168                   ncdf_varput, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1, jpk], offset = [(jpi-1)/2+2, jpj-2, 0]
169              ENDIF
170            END
171            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
172            END
173            ELSE:print, vinq.name+' nothing to check'
174          ENDCASE
175; x periodicity (1)
176          ncdf_varget, cdfid, vid, part1, count = [1, jpj, jpk], offset = [0, 0, 0]
177          ncdf_varget, cdfid, vid, part2, count = [1, jpj, jpk], offset = [jpi-2, 0, 0]
178          IF array_equal(part1, part2) NE 1 THEN BEGIN
179            print, vinq.name+' bad x periodicity (1)'
180            IF keyword_set(write) THEN $
181               ncdf_varput, cdfid, vid, part2, count = [1, jpj, jpk], offset = [0, 0, 0]
182          ENDIF
183; x periodicity (2)
184          ncdf_varget, cdfid, vid, part1, count = [1, jpj, jpk], offset = [jpi-1, 0, 0]
185          ncdf_varget, cdfid, vid, part2, count = [1, jpj, jpk], offset = [1, 0, 0]
186          IF array_equal(part1, part2) NE 1 THEN BEGIN
187            print, vinq.name+' bad x periodicity (2)'
188            IF keyword_set(write) THEN $
189               ncdf_varput, cdfid, vid, part2, count = [1, jpj, jpk], offset = [jpi-1, 0, 0]
190          ENDIF
191        END
192;------------------------------
193; 4D VAR
194;------------------------------
195        4:BEGIN
196          jpk = vdimsz[2]
197          jpt = vdimsz[3]
198; y periodicity (1)
199          CASE 1 OF
200            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
201              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-1, 0, 0]
202              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-3, 0, 0]
203              IF jpt EQ 1 THEN part2 = reform(part2, jpi-2, 1, jpk, jpt, /over)
204              part2 = reverse(shift(part2, -1, 0, 0, 0), 1)
205              IF array_equal(part1, part2) NE 1 THEN BEGIN
206                print, vinq.name+' bad y periodicity (1)'
207                IF keyword_set(write) THEN $
208                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-1, 0, 0]
209              ENDIF
210            END
211            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
212              ncdf_varget, cdfid, vid, part1, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-1, 0, 0]
213              ncdf_varget, cdfid, vid, part2, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-2, 0, 0]
214              IF jpt EQ 1 THEN part2 = reform(part2, jpi-2, 1, jpk, jpt, /over)
215              part2 = reverse(part2, 1)
216              IF array_equal(part1, part2) NE 1 THEN BEGIN
217                print, vinq.name+' bad y periodicity (1)'
218                IF keyword_set(write) THEN $
219                   ncdf_varput, cdfid, vid, part2, count = [jpi-2, 1, jpk, jpt], offset = [1, jpj-1, 0, 0]
220              ENDIF
221            END
222            ELSE:print, vinq.name+' nothing to check'
223          ENDCASE
224; y periodicity (2)
225          CASE 1 OF
226            array_equal([182, 149], vdimsz[0:1]):BEGIN ; ORCA2
227              ncdf_varget, cdfid, vid, part1, count = [(jpi-1)/2-2+1, 1, jpk, jpt], offset = [(jpi-1)/2+2, jpj-2, 0, 0]
228              ncdf_varget, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1, jpk, jpt], offset = [          2, jpj-2, 0, 0]
229              part2 = reverse(part2, 1)
230              IF array_equal(part1, part2) NE 1 THEN BEGIN
231                print, vinq.name+' bad y periodicity (2)'
232                IF keyword_set(write) THEN $
233                   ncdf_varput, cdfid, vid, part2, count = [(jpi-1)/2-2+1, 1, jpk, jpt], offset = [(jpi-1)/2+2, jpj-2, 0, 0]
234              ENDIF
235            END
236            array_equal([722, 511], vdimsz[0:1]):BEGIN ; ORCA05
237            END
238            ELSE:print, vinq.name+' nothing to check'
239          ENDCASE
240; x periodicity (1)
241          ncdf_varget, cdfid, vid, part1, count = [1, jpj, jpk, jpt], offset = [0, 0, 0, 0]
242          ncdf_varget, cdfid, vid, part2, count = [1, jpj, jpk, jpt], offset = [jpi-2, 0, 0, 0]
243          IF array_equal(part1, part2) NE 1 THEN BEGIN
244            print, vinq.name+' bad x periodicity (1)'
245            IF keyword_set(write) THEN $
246               ncdf_varput, cdfid, vid, part2, count = [1, jpj, jpk, jpt], offset = [0, 0, 0, 0]
247          ENDIF
248; x periodicity (2)
249          ncdf_varget, cdfid, vid, part1, count = [1, jpj, jpk, jpt], offset = [jpi-1, 0, 0, 0]
250          ncdf_varget, cdfid, vid, part2, count = [1, jpj, jpk, jpt], offset = [1, 0, 0, 0]
251          IF array_equal(part1, part2) NE 1 THEN BEGIN
252            print, vinq.name+' bad x periodicity (2)'
253            IF keyword_set(write) THEN $
254               ncdf_varput, cdfid, vid, part2, count = [1, jpj, jpk, jpt], offset = [jpi-1, 0, 0, 0]
255          ENDIF
256        END
257        ELSE:print, vinq.name+' nothing to check'
258      ENDCASE
259    ENDIF   
260  ENDFOR
261
262  ncdf_close, cdfid
263
264  RETURN
265END
Note: See TracBrowser for help on using the repository browser.