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

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

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


Revision 1.13 - (hide 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 mlosch 1.13 C $Header: /u/gcmpack/MITgcm/pkg/mnc/mnc_cw_cvars.F,v 1.12 2008/10/25 20:36:34 mlosch Exp $
2 edhill 1.1 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 mlosch 1.9 #include "MNC_COMMON.h"
25 edhill 1.1 #include "SIZE.h"
26     #include "EEPARAMS.h"
27     #include "EESUPPORT.h"
28     #include "PARAMS.h"
29     #include "GRID.h"
30 edhill 1.4 #ifdef ALLOW_EXCH2
31     #include "W2_EXCH2_TOPOLOGY.h"
32     #include "W2_EXCH2_PARAMS.h"
33     #endif
34 utke 1.11 #include "netcdf.inc"
35 edhill 1.4
36     C Functions
37     integer IFNBLNK, ILNBLNK
38 edhill 1.1
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 mlosch 1.7 integer i, vid, nnf, nnl, doit, err
48 edhill 1.4 integer nids, cv_did(1), xtmin,ytmin
49 edhill 1.1 character*(MAX_LEN_MBUF) msgbuf
50     integer cv_start(1), cv_count(1)
51 edhill 1.3 _RS rtmp(sNx + 2*OLx + sNy + 2*OLy + Nr)
52 mlosch 1.10 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 edhill 1.1
63     nnf = IFNBLNK(cvname)
64     nnl = ILNBLNK(cvname)
65    
66 edhill 1.4 xtmin = 0
67     ytmin = 0
68     #ifdef ALLOW_EXCH2
69     xtmin = exch2_tbasex(W2_myTileList(bi))
70     ytmin = exch2_tbasey(W2_myTileList(bi))
71 mlosch 1.13 #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 edhill 1.4 #endif
81     doit = 1
82 edhill 1.1 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 mlosch 1.7 #ifdef ALLOW_EXCH2
91 edhill 1.1 DO i = cv_start(1),cv_count(1)
92 mlosch 1.7 rtmp(i) = xtmin + i
93     ENDDO
94     #else
95 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
96 mlosch 1.7 DO i = cv_start(1),cv_count(1)
97 edhill 1.4 rtmp(i) = xtmin + i
98 mlosch 1.7 ENDDO
99     ELSE
100     DO i = cv_start(1),cv_count(1)
101 edhill 1.1 rtmp(i) = xC(i,1,bi,bj)
102 mlosch 1.7 ENDDO
103     ENDIF
104 edhill 1.4 #endif
105 mlosch 1.10 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 edhill 1.1
122     ELSEIF (cvname(nnf:nnl) .EQ. 'Xp1') THEN
123    
124     cv_start(1) = 1
125     cv_count(1) = sNx + 1
126 mlosch 1.7 #ifdef ALLOW_EXCH2
127 edhill 1.1 DO i = cv_start(1),cv_count(1)
128 mlosch 1.7 rtmp(i) = xtmin + i
129     ENDDO
130     #else
131 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
132 mlosch 1.7 DO i = cv_start(1),cv_count(1)
133 edhill 1.4 rtmp(i) = xtmin + i
134 mlosch 1.7 ENDDO
135     ELSE
136     DO i = cv_start(1),cv_count(1)
137 edhill 1.1 rtmp(i) = xG(i,1,bi,bj)
138 mlosch 1.7 ENDDO
139     ENDIF
140 edhill 1.4 #endif
141 mlosch 1.10 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 edhill 1.1
158 edhill 1.3 ELSEIF (cvname(nnf:nnl) .EQ. 'Xwh') THEN
159    
160     cv_start(1) = 1
161     cv_count(1) = sNx + 2*OLx
162 mlosch 1.7 #ifdef ALLOW_EXCH2
163 edhill 1.3 DO i = cv_start(1),cv_count(1)
164 mlosch 1.7 rtmp(i) = xtmin + i
165     ENDDO
166     #else
167 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
168 mlosch 1.7 DO i = cv_start(1),cv_count(1)
169 edhill 1.4 rtmp(i) = xtmin - OLx + i
170 mlosch 1.7 ENDDO
171     ELSE
172     DO i = cv_start(1),cv_count(1)
173 edhill 1.3 rtmp(i) = xC(i,1,bi,bj)
174 mlosch 1.10 CML???? rtmp(i) = xC(i-Olx,1,bi,bj)
175 mlosch 1.7 ENDDO
176     ENDIF
177 edhill 1.4 #endif
178 mlosch 1.10 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 edhill 1.4
196 edhill 1.1 ELSEIF (cvname(nnf:nnl) .EQ. 'Y') THEN
197    
198     cv_start(1) = 1
199     cv_count(1) = sNy
200 mlosch 1.7 #ifdef ALLOW_EXCH2
201 edhill 1.1 DO i = cv_start(1),cv_count(1)
202 mlosch 1.7 rtmp(i) = ytmin + i
203     ENDDO
204     #else
205 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
206 mlosch 1.7 DO i = cv_start(1),cv_count(1)
207 edhill 1.4 rtmp(i) = ytmin + i
208 mlosch 1.7 ENDDO
209     ELSE
210     DO i = cv_start(1),cv_count(1)
211 edhill 1.1 rtmp(i) = yC(1,i,bi,bj)
212 mlosch 1.7 ENDDO
213     ENDIF
214 edhill 1.4 #endif
215 mlosch 1.10 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 edhill 1.1
232     ELSEIF (cvname(nnf:nnl) .EQ. 'Yp1') THEN
233    
234     cv_start(1) = 1
235     cv_count(1) = sNy + 1
236 mlosch 1.7 #ifdef ALLOW_EXCH2
237 edhill 1.1 DO i = cv_start(1),cv_count(1)
238 mlosch 1.7 rtmp(i) = ytmin + i
239     ENDDO
240     #else
241 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
242 mlosch 1.7 DO i = cv_start(1),cv_count(1)
243 edhill 1.4 rtmp(i) = ytmin + i
244 mlosch 1.7 ENDDO
245     ELSE
246     DO i = cv_start(1),cv_count(1)
247 edhill 1.1 rtmp(i) = yG(1,i,bi,bj)
248 mlosch 1.7 ENDDO
249     ENDIF
250 edhill 1.4 #endif
251 mlosch 1.10 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 edhill 1.1
268 edhill 1.3 ELSEIF (cvname(nnf:nnl) .EQ. 'Ywh') THEN
269    
270     cv_start(1) = 1
271     cv_count(1) = sNy + 2*OLy
272 mlosch 1.7 #ifdef ALLOW_EXCH2
273 edhill 1.3 DO i = cv_start(1),cv_count(1)
274 mlosch 1.7 rtmp(i) = ytmin + i
275     ENDDO
276     #else
277 mlosch 1.8 IF ( usingCurvilinearGrid .OR. rotateGrid ) THEN
278 mlosch 1.7 DO i = cv_start(1),cv_count(1)
279 edhill 1.4 rtmp(i) = ytmin - OLy + i
280 mlosch 1.7 ENDDO
281     ELSE
282     DO i = cv_start(1),cv_count(1)
283 edhill 1.3 rtmp(i) = yC(1,i-OLy,bi,bj)
284 mlosch 1.7 ENDDO
285     ENDIF
286 edhill 1.4 #endif
287 mlosch 1.10 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 edhill 1.3
305 edhill 1.1 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 mlosch 1.10 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 edhill 1.1
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 mlosch 1.10 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 edhill 1.4
343 edhill 1.5 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 mlosch 1.10 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 edhill 1.5
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 mlosch 1.10 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 edhill 1.5
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 mlosch 1.10 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 edhill 1.5
403 edhill 1.4 ELSE
404    
405     doit = 0
406 edhill 1.1
407     ENDIF
408    
409     IF ( doit .EQ. 1 ) THEN
410    
411     CALL MNC_FILE_REDEF(fname, myThid)
412 mlosch 1.12 #ifdef REAL4_IS_SLOW
413 edhill 1.1 err = NF_DEF_VAR(fid, cvname, NF_DOUBLE,
414     & nids, cv_did, vid)
415 mlosch 1.12 #else
416     err = NF_DEF_VAR(fid, cvname, NF_FLOAT,
417     & nids, cv_did, vid)
418     #endif /* REAL4_IS_SLOW */
419 edhill 1.6 i = ILNBLNK( fname )
420 edhill 1.1 write(msgbuf,'(5a)') 'defining coordinate variable ''',
421 edhill 1.6 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
422 edhill 1.1 CALL MNC_HANDLE_ERR(err, msgbuf, myThid)
423 mlosch 1.10 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 edhill 1.1 CALL MNC_FILE_ENDDEF(fname, myThid)
450 mlosch 1.12 #ifdef REAL4_IS_SLOW
451 edhill 1.1 err = NF_PUT_VARA_DOUBLE(fid, vid,
452     & cv_start, cv_count, rtmp)
453 mlosch 1.12 #else
454     err = NF_PUT_VARA_REAL(fid, vid,
455     & cv_start, cv_count, rtmp)
456     #endif /* REAL4_IS_SLOW */
457 edhill 1.1 write(msgbuf,'(5a)') 'writing coordinate variable ''',
458 edhill 1.6 & cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
459 edhill 1.1 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