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

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

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

revision 1.6 by edhill, Fri Mar 10 16:09:31 2006 UTC revision 1.14 by jmc, Tue May 12 19:56:36 2009 UTC
# Line 21  C     Write a CF-convention coordinate v Line 21  C     Write a CF-convention coordinate v
21    
22  C     !USES:  C     !USES:
23        implicit none        implicit none
24  #include "netcdf.inc"  #include "MNC_COMMON.h"
 #include "mnc_common.h"  
25  #include "SIZE.h"  #include "SIZE.h"
26  #include "EEPARAMS.h"  #include "EEPARAMS.h"
27  #include "EESUPPORT.h"  #include "EESUPPORT.h"
28  #include "PARAMS.h"  #include "PARAMS.h"
29  #include "GRID.h"  #include "GRID.h"
30  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
31    #include "W2_EXCH2_SIZE.h"
32  #include "W2_EXCH2_TOPOLOGY.h"  #include "W2_EXCH2_TOPOLOGY.h"
 #include "W2_EXCH2_PARAMS.h"  
33  #endif  #endif
34    #include "netcdf.inc"
35    
36  C     Functions  C     Functions
37        integer IFNBLNK, ILNBLNK        integer IFNBLNK, ILNBLNK
# Line 44  C     !INPUT PARAMETERS: Line 44  C     !INPUT PARAMETERS:
44  CEOP  CEOP
45    
46  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
47        integer i,j, vid, nnf,nnl, doit, err        integer i, vid, nnf, nnl, doit, err
48        integer nids, cv_did(1), xtmin,ytmin        integer nids, cv_did(1), xtmin,ytmin
49        character*(MAX_LEN_MBUF) msgbuf        character*(MAX_LEN_MBUF) msgbuf
50        integer cv_start(1), cv_count(1)        integer cv_start(1), cv_count(1)
51        _RS rtmp(sNx + 2*OLx + sNy + 2*OLy + Nr)        _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)        nnf = IFNBLNK(cvname)
64        nnl = ILNBLNK(cvname)        nnl = ILNBLNK(cvname)
# Line 58  C     !LOCAL VARIABLES: Line 68  C     !LOCAL VARIABLES:
68  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
69        xtmin = exch2_tbasex(W2_myTileList(bi))        xtmin = exch2_tbasex(W2_myTileList(bi))
70        ytmin = exch2_tbasey(W2_myTileList(bi))        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  #endif
81        doit = 1        doit = 1
82        nids = 1        nids = 1
# Line 68  C     Check all the coordinate variables Line 87  C     Check all the coordinate variables
87    
88          cv_start(1) = 1          cv_start(1) = 1
89          cv_count(1) = sNx          cv_count(1) = sNx
         DO i = cv_start(1),cv_count(1)  
90  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
91            rtmp(i) = xtmin + i          DO i = cv_start(1),cv_count(1)
92             rtmp(i) = xtmin + i
93            ENDDO
94  #else  #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)            rtmp(i) = xC(i,1,bi,bj)
102             ENDDO
103            ENDIF
104  #endif  #endif
105          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Xp1') THEN
123    
124          cv_start(1) = 1          cv_start(1) = 1
125          cv_count(1) = sNx + 1          cv_count(1) = sNx + 1
         DO i = cv_start(1),cv_count(1)  
126  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
127            rtmp(i) = xtmin + i          DO i = cv_start(1),cv_count(1)
128             rtmp(i) = xtmin + i
129            ENDDO
130  #else  #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)            rtmp(i) = xG(i,1,bi,bj)
138             ENDDO
139            ENDIF
140  #endif  #endif
141          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Xwh') THEN
159    
160          cv_start(1) = 1          cv_start(1) = 1
161          cv_count(1) = sNx + 2*OLx          cv_count(1) = sNx + 2*OLx
         DO i = cv_start(1),cv_count(1)  
162  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
163            rtmp(i) = xtmin - OLx + i          DO i = cv_start(1),cv_count(1)
164             rtmp(i) = xtmin + i
165            ENDDO
166  #else  #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)            rtmp(i) = xC(i,1,bi,bj)
174    CML????          rtmp(i) = xC(i-Olx,1,bi,bj)
175             ENDDO
176            ENDIF
177  #endif  #endif
178          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Y') THEN
197    
198          cv_start(1) = 1          cv_start(1) = 1
199          cv_count(1) = sNy          cv_count(1) = sNy
         DO i = cv_start(1),cv_count(1)  
200  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
201            rtmp(i) = ytmin + i          DO i = cv_start(1),cv_count(1)
202             rtmp(i) = ytmin + i
203            ENDDO
204  #else  #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)            rtmp(i) = yC(1,i,bi,bj)
212             ENDDO
213            ENDIF
214  #endif  #endif
215          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Yp1') THEN
233    
234          cv_start(1) = 1          cv_start(1) = 1
235          cv_count(1) = sNy + 1          cv_count(1) = sNy + 1
         DO i = cv_start(1),cv_count(1)  
236  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
237            rtmp(i) = ytmin + i          DO i = cv_start(1),cv_count(1)
238             rtmp(i) = ytmin + i
239            ENDDO
240  #else  #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)            rtmp(i) = yG(1,i,bi,bj)
248             ENDDO
249            ENDIF
250  #endif  #endif
251          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Ywh') THEN
269    
270          cv_start(1) = 1          cv_start(1) = 1
271          cv_count(1) = sNy + 2*OLy          cv_count(1) = sNy + 2*OLy
         DO i = cv_start(1),cv_count(1)  
272  #ifdef ALLOW_EXCH2  #ifdef ALLOW_EXCH2
273            rtmp(i) = ytmin - OLy + i          DO i = cv_start(1),cv_count(1)
274             rtmp(i) = ytmin + i
275            ENDDO
276  #else  #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)            rtmp(i) = yC(1,i-OLy,bi,bj)
284             ENDDO
285            ENDIF
286  #endif  #endif
287          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Z') THEN
306    
# Line 143  C     Check all the coordinate variables Line 309  C     Check all the coordinate variables
309          DO i = cv_start(1),cv_count(1)          DO i = cv_start(1),cv_count(1)
310            rtmp(i) = rC(i)            rtmp(i) = rC(i)
311          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Zp1') THEN
325    
# Line 151  C     Check all the coordinate variables Line 328  C     Check all the coordinate variables
328          DO i = cv_start(1),cv_count(1)          DO i = cv_start(1),cv_count(1)
329            rtmp(i) = rF(i)            rtmp(i) = rF(i)
330          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Zu') THEN
344    
# Line 159  C     Check all the coordinate variables Line 347  C     Check all the coordinate variables
347          DO i = cv_start(1),cv_count(1)          DO i = cv_start(1),cv_count(1)
348            rtmp(i) = rF(i + 1)            rtmp(i) = rF(i + 1)
349          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Zl') THEN
364    
# Line 167  C     Check all the coordinate variables Line 367  C     Check all the coordinate variables
367          DO i = cv_start(1),cv_count(1)          DO i = cv_start(1),cv_count(1)
368            rtmp(i) = rF(i)            rtmp(i) = rF(i)
369          ENDDO          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        ELSEIF (cvname(nnf:nnl) .EQ. 'Zm1') THEN
384    
# Line 175  C     Check all the coordinate variables Line 387  C     Check all the coordinate variables
387          DO i = cv_start(1),cv_count(1)          DO i = cv_start(1),cv_count(1)
388            rtmp(i) = rF(i + 1)            rtmp(i) = rF(i + 1)
389          ENDDO          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        ELSE
404    
# Line 185  C     Check all the coordinate variables Line 409  C     Check all the coordinate variables
409        IF ( doit .EQ. 1 ) THEN        IF ( doit .EQ. 1 ) THEN
410    
411          CALL MNC_FILE_REDEF(fname, myThid)          CALL MNC_FILE_REDEF(fname, myThid)
412    #ifdef REAL4_IS_SLOW
413          err = NF_DEF_VAR(fid, cvname, NF_DOUBLE,          err = NF_DEF_VAR(fid, cvname, NF_DOUBLE,
414       &       nids, cv_did, vid)       &       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 )          i = ILNBLNK( fname )
420          write(msgbuf,'(5a)') 'defining coordinate variable ''',          write(msgbuf,'(5a)') 'defining coordinate variable ''',
421       &       cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''       &       cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
422          CALL MNC_HANDLE_ERR(err, msgbuf, myThid)          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)          CALL MNC_FILE_ENDDEF(fname, myThid)
450    #ifdef REAL4_IS_SLOW
451          err = NF_PUT_VARA_DOUBLE(fid, vid,          err = NF_PUT_VARA_DOUBLE(fid, vid,
452       &       cv_start, cv_count, rtmp)       &       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 ''',          write(msgbuf,'(5a)') 'writing coordinate variable ''',
458       &       cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''       &       cvname(nnf:nnl), ''' in file ''', fname(1:i), ''''
459          CALL MNC_HANDLE_ERR(err, msgbuf, myThid)          CALL MNC_HANDLE_ERR(err, msgbuf, myThid)

Legend:
Removed from v.1.6  
changed lines
  Added in v.1.14

  ViewVC Help
Powered by ViewVC 1.1.22