/[MITgcm]/MITgcm/pkg/mnc/mnc_cw_cvars.F
ViewVC logotype

Contents of /MITgcm/pkg/mnc/mnc_cw_cvars.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.13 - (show annotations) (download)
Tue Mar 3 10:01:24 2009 UTC (15 years, 3 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint61n, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k
Changes since 1.12: +10 -1 lines
fix the coordinate variables for multi-tile curvi-linear grids that
are not cubed-sphere grids

1 C $Header: /u/gcmpack/MITgcm/pkg/mnc/mnc_cw_cvars.F,v 1.12 2008/10/25 20:36:34 mlosch Exp $
2 C $Name: $
3
4 #include "MNC_OPTIONS.h"
5
6 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7 CBOP 1
8 C !ROUTINE: MNC_CW_WRITE_CVAR
9
10 C !INTERFACE:
11 SUBROUTINE MNC_CW_WRITE_CVAR(
12 I fname,
13 I cvname,
14 I fid,
15 I did,
16 I bi, bj,
17 I myThid )
18
19 C !DESCRIPTION:
20 C Write a CF-convention coordinate variable (a vector).
21
22 C !USES:
23 implicit none
24 #include "MNC_COMMON.h"
25 #include "SIZE.h"
26 #include "EEPARAMS.h"
27 #include "EESUPPORT.h"
28 #include "PARAMS.h"
29 #include "GRID.h"
30 #ifdef ALLOW_EXCH2
31 #include "W2_EXCH2_TOPOLOGY.h"
32 #include "W2_EXCH2_PARAMS.h"
33 #endif
34 #include "netcdf.inc"
35
36 C Functions
37 integer IFNBLNK, ILNBLNK
38
39 C !INPUT PARAMETERS:
40 character*(*) fname
41 character*(*) cvname
42 integer fid, did, bi,bj
43 integer myThid
44 CEOP
45
46 C !LOCAL VARIABLES:
47 integer i, vid, nnf, nnl, doit, err
48 integer nids, cv_did(1), xtmin,ytmin
49 character*(MAX_LEN_MBUF) msgbuf
50 integer cv_start(1), cv_count(1)
51 _RS rtmp(sNx + 2*OLx + sNy + 2*OLy + Nr)
52 C variables for text attributes
53 integer MAX_LEN_NAME, ia
54 PARAMETER ( MAX_LEN_NAME = 128 )
55 character*(MAX_LEN_NAME) units, long_name, positive
56
57 DO i=1,MAX_LEN_NAME
58 units(i:i) = ' '
59 long_name(i:i) = ' '
60 positive(i:i) = ' '
61 ENDDO
62
63 nnf = IFNBLNK(cvname)
64 nnl = ILNBLNK(cvname)
65
66 xtmin = 0
67 ytmin = 0
68 #ifdef ALLOW_EXCH2
69 xtmin = exch2_tbasex(W2_myTileList(bi))
70 ytmin = exch2_tbasey(W2_myTileList(bi))
71 #else
72 IF ( .NOT. useCubedSphereExchange ) THEN
73 C make sure for a non-cubed-sphere curvi-linear grid,
74 C that the X/Y coordinate variables are monotonous
75 C bi+(myXGlobalLo-1)/sNx is the i-tile number
76 C bj+(myYGlobalLo-1)/sNy is the j-tile number
77 xtmin = sNx * ( bi+(myXGlobalLo-1)/sNx - 1 )
78 ytmin = sNy * ( bj+(myYGlobalLo-1)/sNy - 1 )
79 ENDIF
80 #endif
81 doit = 1
82 nids = 1
83 cv_did(1)= did
84
85 C Check all the coordinate variables that we know about
86 IF (cvname(nnf:nnl) .EQ. 'X') THEN
87
88 cv_start(1) = 1
89 cv_count(1) = sNx
90 #ifdef ALLOW_EXCH2
91 DO i = cv_start(1),cv_count(1)
92 rtmp(i) = xtmin + i
93 ENDDO
94 #else
95 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
96 DO i = cv_start(1),cv_count(1)
97 rtmp(i) = xtmin + i
98 ENDDO
99 ELSE
100 DO i = cv_start(1),cv_count(1)
101 rtmp(i) = xC(i,1,bi,bj)
102 ENDDO
103 ENDIF
104 #endif
105 IF ( usingCartesianGrid ) THEN
106 long_name = 'X-coordinate of cell center'
107 units = 'meters'
108 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
109 long_name = 'i-index of cell center'
110 units = 'none'
111 ELSEIF ( usingSphericalPolarGrid ) THEN
112 long_name = 'longitude of cell center'
113 units = 'degrees_east'
114 ELSEIF ( usingCylindricalGrid ) THEN
115 long_name = 'polar angle coordinate of cell center'
116 units = 'degrees'
117 ELSE
118 C unknown grid type
119 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
120 ENDIF
121
122 ELSEIF (cvname(nnf:nnl) .EQ. 'Xp1') THEN
123
124 cv_start(1) = 1
125 cv_count(1) = sNx + 1
126 #ifdef ALLOW_EXCH2
127 DO i = cv_start(1),cv_count(1)
128 rtmp(i) = xtmin + i
129 ENDDO
130 #else
131 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
132 DO i = cv_start(1),cv_count(1)
133 rtmp(i) = xtmin + i
134 ENDDO
135 ELSE
136 DO i = cv_start(1),cv_count(1)
137 rtmp(i) = xG(i,1,bi,bj)
138 ENDDO
139 ENDIF
140 #endif
141 IF ( usingCartesianGrid ) THEN
142 long_name = 'X-Coordinate of cell corner'
143 units = 'meters'
144 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
145 long_name = 'i-index of cell corner'
146 units = 'none'
147 ELSEIF ( usingSphericalPolarGrid ) THEN
148 long_name = 'longitude of cell corner'
149 units = 'degrees_east'
150 ELSEIF ( usingCylindricalGrid ) THEN
151 long_name = 'polar angle of cell corner'
152 units = 'degrees'
153 ELSE
154 C unknown grid type
155 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
156 ENDIF
157
158 ELSEIF (cvname(nnf:nnl) .EQ. 'Xwh') THEN
159
160 cv_start(1) = 1
161 cv_count(1) = sNx + 2*OLx
162 #ifdef ALLOW_EXCH2
163 DO i = cv_start(1),cv_count(1)
164 rtmp(i) = xtmin + i
165 ENDDO
166 #else
167 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
168 DO i = cv_start(1),cv_count(1)
169 rtmp(i) = xtmin - OLx + i
170 ENDDO
171 ELSE
172 DO i = cv_start(1),cv_count(1)
173 rtmp(i) = xC(i,1,bi,bj)
174 CML???? rtmp(i) = xC(i-Olx,1,bi,bj)
175 ENDDO
176 ENDIF
177 #endif
178 IF ( usingCartesianGrid ) THEN
179 long_name = 'X-Coordinate of cell center including overlaps'
180 units = 'meters'
181 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
182 long_name = 'i-index of cell center including overlaps'
183 units = 'none'
184 ELSEIF ( usingSphericalPolarGrid ) THEN
185 long_name = 'longitude of cell center including overlaps'
186 units = 'degrees_east'
187 ELSEIF ( usingCylindricalGrid ) THEN
188 long_name =
189 & 'polar angle coordinate of cell center including overlaps'
190 units = 'degrees'
191 ELSE
192 C unknown grid type
193 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
194 ENDIF
195
196 ELSEIF (cvname(nnf:nnl) .EQ. 'Y') THEN
197
198 cv_start(1) = 1
199 cv_count(1) = sNy
200 #ifdef ALLOW_EXCH2
201 DO i = cv_start(1),cv_count(1)
202 rtmp(i) = ytmin + i
203 ENDDO
204 #else
205 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
206 DO i = cv_start(1),cv_count(1)
207 rtmp(i) = ytmin + i
208 ENDDO
209 ELSE
210 DO i = cv_start(1),cv_count(1)
211 rtmp(i) = yC(1,i,bi,bj)
212 ENDDO
213 ENDIF
214 #endif
215 IF ( usingCartesianGrid ) THEN
216 long_name = 'Y-Coordinate of cell center'
217 units = 'meters'
218 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
219 long_name = 'j-index of cell center'
220 units = 'none'
221 ELSEIF ( usingSphericalPolarGrid ) THEN
222 long_name = 'latitude of cell center'
223 units = 'degrees_north'
224 ELSEIF ( usingCylindricalGrid ) THEN
225 long_name = 'radial coordinate of cell center'
226 units = 'meters'
227 ELSE
228 C unknown grid type
229 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
230 ENDIF
231
232 ELSEIF (cvname(nnf:nnl) .EQ. 'Yp1') THEN
233
234 cv_start(1) = 1
235 cv_count(1) = sNy + 1
236 #ifdef ALLOW_EXCH2
237 DO i = cv_start(1),cv_count(1)
238 rtmp(i) = ytmin + i
239 ENDDO
240 #else
241 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
242 DO i = cv_start(1),cv_count(1)
243 rtmp(i) = ytmin + i
244 ENDDO
245 ELSE
246 DO i = cv_start(1),cv_count(1)
247 rtmp(i) = yG(1,i,bi,bj)
248 ENDDO
249 ENDIF
250 #endif
251 IF ( usingCartesianGrid ) THEN
252 long_name = 'Y-Coordinate of cell corner'
253 units = 'meters'
254 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
255 long_name = 'j-index of cell corner'
256 units = 'none'
257 ELSEIF ( usingSphericalPolarGrid ) THEN
258 long_name = 'latitude of cell corner'
259 units = 'degrees_north'
260 ELSEIF ( usingCylindricalGrid ) THEN
261 long_name = 'radial coordinate of cell corner'
262 units = 'meters'
263 ELSE
264 C unknown grid type
265 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
266 ENDIF
267
268 ELSEIF (cvname(nnf:nnl) .EQ. 'Ywh') THEN
269
270 cv_start(1) = 1
271 cv_count(1) = sNy + 2*OLy
272 #ifdef ALLOW_EXCH2
273 DO i = cv_start(1),cv_count(1)
274 rtmp(i) = ytmin + i
275 ENDDO
276 #else
277 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
278 DO i = cv_start(1),cv_count(1)
279 rtmp(i) = ytmin - OLy + i
280 ENDDO
281 ELSE
282 DO i = cv_start(1),cv_count(1)
283 rtmp(i) = yC(1,i-OLy,bi,bj)
284 ENDDO
285 ENDIF
286 #endif
287 IF ( usingCartesianGrid ) THEN
288 long_name = 'Y-Coordinate of cell center including overlaps'
289 units = 'meters'
290 ELSEIF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
291 long_name = 'j-index of cell center including overlaps'
292 units = 'none'
293 ELSEIF ( usingSphericalPolarGrid ) THEN
294 long_name = 'latitude of cell center including overlaps'
295 units = 'degrees_north'
296 ELSEIF ( usingCylindricalGrid ) THEN
297 long_name =
298 & 'radial coordinate of cell center including overlaps'
299 units = 'meters'
300 ELSE
301 C unknown grid type
302 print *, 'S/R MNC_CW_CVARS: Ooops, unknown horizontal grid!'
303 ENDIF
304
305 ELSEIF (cvname(nnf:nnl) .EQ. 'Z') THEN
306
307 cv_start(1) = 1
308 cv_count(1) = Nr
309 DO i = cv_start(1),cv_count(1)
310 rtmp(i) = rC(i)
311 ENDDO
312 C
313 long_name = 'vertical coordinate of cell center'
314 IF ( usingZCoords ) THEN
315 units = 'meters'
316 positive = 'up'
317 ELSEIF ( usingPCoords ) THEN
318 units = 'pascal'
319 ELSE
320 C unknown grid type
321 print *, 'S/R MNC_CW_CVARS: Ooops, unknown vertical grid!'
322 ENDIF
323
324 ELSEIF (cvname(nnf:nnl) .EQ. 'Zp1') THEN
325
326 cv_start(1) = 1
327 cv_count(1) = Nr + 1
328 DO i = cv_start(1),cv_count(1)
329 rtmp(i) = rF(i)
330 ENDDO
331 C
332 long_name = 'vertical coordinate of cell interface'
333 IF ( usingZCoords ) THEN
334 units = 'meters'
335 positive = 'up'
336 ELSEIF ( usingPCoords ) THEN
337 units = 'pascal'
338 ELSE
339 C unknown grid type
340 print *, 'S/R MNC_CW_CVARS: Ooops, unknown vertical grid!'
341 ENDIF
342
343 ELSEIF (cvname(nnf:nnl) .EQ. 'Zu') THEN
344
345 cv_start(1) = 1
346 cv_count(1) = Nr
347 DO i = cv_start(1),cv_count(1)
348 rtmp(i) = rF(i + 1)
349 ENDDO
350 C
351 IF ( usingZCoords ) THEN
352 long_name = 'vertical coordinate of lower cell interface'
353 units = 'meters'
354 positive = 'up'
355 ELSEIF ( usingPCoords ) THEN
356 long_name = 'vertical coordinate of upper cell interface'
357 units = 'pascal'
358 ELSE
359 C unknown grid type
360 print *, 'S/R MNC_CW_CVARS: Ooops, unknown vertical grid!'
361 ENDIF
362
363 ELSEIF (cvname(nnf:nnl) .EQ. 'Zl') THEN
364
365 cv_start(1) = 1
366 cv_count(1) = Nr
367 DO i = cv_start(1),cv_count(1)
368 rtmp(i) = rF(i)
369 ENDDO
370 C
371 IF ( usingZCoords ) THEN
372 long_name = 'vertical coordinate of upper cell interface'
373 units = 'meters'
374 positive = 'up'
375 ELSEIF ( usingPCoords ) THEN
376 long_name = 'vertical coordinate of lower cell interface'
377 units = 'pascal'
378 ELSE
379 C unknown grid type
380 print *, 'S/R MNC_CW_CVARS: Ooops, unknown vertical grid!'
381 ENDIF
382
383 ELSEIF (cvname(nnf:nnl) .EQ. 'Zm1') THEN
384
385 cv_start(1) = 1
386 cv_count(1) = Nr - 1
387 DO i = cv_start(1),cv_count(1)
388 rtmp(i) = rF(i + 1)
389 ENDDO
390 C
391 IF ( usingZCoords ) THEN
392 long_name = 'vertical coordinate of lower cell interface'
393 units = 'meters'
394 positive = 'up'
395 ELSEIF ( usingPCoords ) THEN
396 long_name = 'vertical coordinate of upper cell interface'
397 units = 'pascal'
398 ELSE
399 C unknown grid type
400 print *, 'S/R MNC_CW_CVARS: Ooops, unknown vertical grid!'
401 ENDIF
402
403 ELSE
404
405 doit = 0
406
407 ENDIF
408
409 IF ( doit .EQ. 1 ) THEN
410
411 CALL MNC_FILE_REDEF(fname, myThid)
412 #ifdef REAL4_IS_SLOW
413 err = NF_DEF_VAR(fid, cvname, NF_DOUBLE,
414 & nids, cv_did, vid)
415 #else
416 err = NF_DEF_VAR(fid, cvname, NF_FLOAT,
417 & nids, cv_did, vid)
418 #endif /* REAL4_IS_SLOW */
419 i = ILNBLNK( fname )
420 write(msgbuf,'(5a)') 'defining coordinate variable ''',
421 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
422 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
423 C add attributes if set
424 ia = ILNBLNK(long_name)
425 IF ( ia .GT. 0 ) THEN
426 err = NF_PUT_ATT_TEXT(fid, vid, 'long_name', ia, long_name)
427 write(msgbuf,'(5a)')
428 & 'adding attribute ''long_name'' to coordinate variable ''',
429 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
430 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
431 ENDIF
432 ia = ILNBLNK(units)
433 IF ( ia .GT. 0 ) THEN
434 err = NF_PUT_ATT_TEXT(fid, vid, 'units', ia, units)
435 write(msgbuf,'(5a)')
436 & 'adding attribute ''units'' to coordinate variable ''',
437 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
438 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
439 ENDIF
440 ia = ILNBLNK(positive)
441 IF ( ia .GT. 0 ) THEN
442 err = NF_PUT_ATT_TEXT(fid, vid, 'positive', ia, positive)
443 write(msgbuf,'(5a)')
444 & 'adding attribute ''positive'' to coordinate variable ''',
445 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
446 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
447 ENDIF
448 C
449 CALL MNC_FILE_ENDDEF(fname, myThid)
450 #ifdef REAL4_IS_SLOW
451 err = NF_PUT_VARA_DOUBLE(fid, vid,
452 & cv_start, cv_count, rtmp)
453 #else
454 err = NF_PUT_VARA_REAL(fid, vid,
455 & cv_start, cv_count, rtmp)
456 #endif /* REAL4_IS_SLOW */
457 write(msgbuf,'(5a)') 'writing coordinate variable ''',
458 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
459 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
460
461 ENDIF
462
463 RETURN
464 END
465
466 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
467

  ViewVC Help
Powered by ViewVC 1.1.22