/[MITgcm]/MITgcm/pkg/diagnostics/diagnostics_out.F
ViewVC logotype

Diff of /MITgcm/pkg/diagnostics/diagnostics_out.F

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

revision 1.28 by edhill, Tue Feb 7 15:52:02 2006 UTC revision 1.41 by jmc, Sun Jan 25 20:22:57 2009 UTC
# Line 26  C     !USES: Line 26  C     !USES:
26  #include "DIAGNOSTICS_SIZE.h"  #include "DIAGNOSTICS_SIZE.h"
27  #include "DIAGNOSTICS.h"  #include "DIAGNOSTICS.h"
28    
29  #ifdef ALLOW_FIZHI        INTEGER NrMax
30  #include "fizhi_SIZE.h"        PARAMETER( NrMax = numLevels )
 #else  
       INTEGER Nrphys  
       PARAMETER (Nrphys=0)  
 #endif  
   
31    
32  C     !INPUT PARAMETERS:  C     !INPUT PARAMETERS:
33  C     listId  :: Diagnostics list number being written  C     listId  :: Diagnostics list number being written
# Line 43  C     myThid  :: my Thread Id number Line 38  C     myThid  :: my Thread Id number
38        INTEGER listId, myIter, myThid        INTEGER listId, myIter, myThid
39  CEOP  CEOP
40    
41    C     !FUNCTIONS:
42          INTEGER ILNBLNK
43          EXTERNAL ILNBLNK
44    #ifdef ALLOW_FIZHI
45          _RL   getcon
46          EXTERNAL getcon
47    #endif
48    
49  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
50  C     i,j,k :: loop indices  C     i,j,k :: loop indices
51    C     lm    :: loop index (averageCycle)
52  C     md    :: field number in the list "listId".  C     md    :: field number in the list "listId".
53  C     ndId  :: diagnostics  Id number (in available diagnostics list)  C     ndId  :: diagnostics  Id number (in available diagnostics list)
54  C     mate  :: counter mate Id number (in available diagnostics list)  C     mate  :: counter mate Id number (in available diagnostics list)
55  C     ip    :: diagnostics  pointer to storage array  C     ip    :: diagnostics  pointer to storage array
56  C     im    :: counter-mate pointer to storage array  C     im    :: counter-mate pointer to storage array
57        INTEGER i, j, k  C
58    C--   COMMON /LOCAL_DIAGNOSTICS_OUT/ local common block (for multi-threaded)
59    C     qtmp1 :: thread-shared temporary array (needs to be in common block):
60    C              to write a diagnostic field to file, copy it first from (big)
61    C              diagnostic storage qdiag into it.
62          COMMON /LOCAL_DIAGNOSTICS_OUT/ qtmp1
63          _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,NrMax,nSx,nSy)
64    
65          INTEGER i, j, k, lm, klev
66        INTEGER bi, bj        INTEGER bi, bj
67        INTEGER md, ndId, ip, im        INTEGER md, ndId, ip, im
68        INTEGER mate, mVec        INTEGER mate, mVec
69        CHARACTER*8 parms1        CHARACTER*10 gcode
70        CHARACTER*3 mate_index        _RL undef
71        _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+Nrphys,nSx,nSy)        _RL tmpLev
       _RL undef, getcon  
       EXTERNAL getcon  
       INTEGER ILNBLNK  
       EXTERNAL ILNBLNK  
72        INTEGER ilen        INTEGER ilen
       INTEGER nlevsout  
73    
74        INTEGER ioUnit        INTEGER ioUnit
75        CHARACTER*(MAX_LEN_FNAM) fn        CHARACTER*(MAX_LEN_FNAM) fn
76        CHARACTER*(MAX_LEN_MBUF) suff        CHARACTER*(MAX_LEN_MBUF) suff
77        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
78          INTEGER prec, nRec
79    #ifdef ALLOW_MDSIO
80        LOGICAL glf        LOGICAL glf
81    #endif
82  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
83          INTEGER ll, llMx, jj, jjMx
84        INTEGER ii        INTEGER ii
85        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn
86        INTEGER CW_DIMS, NLEN        INTEGER CW_DIMS, NLEN
# Line 81  C     im    :: counter-mate pointer to s Line 92  C     im    :: counter-mate pointer to s
92        CHARACTER*(NLEN) dn_blnk        CHARACTER*(NLEN) dn_blnk
93  #ifdef DIAG_MNC_COORD_NEEDSWORK  #ifdef DIAG_MNC_COORD_NEEDSWORK
94        CHARACTER*(5) ctmp        CHARACTER*(5) ctmp
95        _RS ztmp(Nr+Nrphys)        _RS ztmp(NrMax)
96  #endif  #endif
97          LOGICAL useMissingValue, useMisValForThisDiag
98          REAL*8 misvalLoc
99          REAL*8 misval_r8(2)
100          REAL*4 misval_r4(2)
101          INTEGER misvalIntLoc, misval_int(2)
102  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
103    
104  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
105    
106        ioUnit= standardMessageUnit        ioUnit= standardMessageUnit
107          undef = UNSET_RL
108    #ifdef ALLOW_FIZHI
109    c     IF ( useFIZHI ) undef = getcon('UNDEF')
110        undef = getcon('UNDEF')        undef = getcon('UNDEF')
111        glf = globalFiles  #endif
112        WRITE(suff,'(I10.10)') myIter        WRITE(suff,'(I10.10)') myIter
113        ilen = ILNBLNK(fnames(listId))        ilen = ILNBLNK(fnames(listId))
114        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:ilen),'.',suff(1:10)        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:ilen),'.',suff(1:10)
115    
116  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
117    C-- this is a trick to reverse the order of the loops on md (= field)
118    C   and lm (= averagePeriod): binary output: lm loop inside md loop ;
119    C                                 mnc ouput: md loop inside lm loop.
120        IF (useMNC .AND. diag_mnc) THEN        IF (useMNC .AND. diag_mnc) THEN
121          DO i = 1,MAX_LEN_FNAM          jjMx = averageCycle(listId)
122            diag_mnc_bn(i:i) = ' '          llMx = 1
123          ENDDO        ELSE
124          DO i = 1,NLEN          jjMx = 1
125            dn_blnk(i:i) = ' '          llMx = averageCycle(listId)
126          ENDDO        ENDIF
127          WRITE( diag_mnc_bn, '(A)' ) fnames(listId)(1:ilen)        DO jj=1,jjMx
128    
129           IF (useMNC .AND. diag_mnc) THEN
130    C     Handle missing value attribute (land points)
131             useMissingValue = .FALSE.
132    #ifdef DIAGNOSTICS_MISSING_VALUE
133             useMissingValue = .TRUE.
134    #endif /* DIAGNOSTICS_MISSING_VALUE */
135             IF ( misvalFlt(listId) .NE. UNSET_RL ) THEN
136              misvalLoc = misvalFlt(listId)
137             ELSE
138              misvalLoc = undef
139             ENDIF
140    C     Defaults to UNSET_I
141             misvalIntLoc = misvalInt(listId)
142             DO ii=1,2
143    C         misval_r4(ii)  = UNSET_FLOAT4
144    C         misval_r8(ii)  = UNSET_FLOAT8
145              misval_r4(ii)  = misvalLoc
146              misval_r8(ii)  = misvalLoc
147              misval_int(ii) = UNSET_I
148             ENDDO
149             DO i = 1,MAX_LEN_FNAM
150               diag_mnc_bn(i:i) = ' '
151             ENDDO
152             DO i = 1,NLEN
153               dn_blnk(i:i) = ' '
154             ENDDO
155             WRITE( diag_mnc_bn, '(A)' ) fnames(listId)(1:ilen)
156    
157  C       Update the record dimension by writing the iteration number  C       Update the record dimension by writing the iteration number
158          CALL MNC_CW_SET_UDIM(diag_mnc_bn, -1, myThid)           klev = myIter + jj - jjMx
159          CALL MNC_CW_RL_W_S('D',diag_mnc_bn,0,0,'T',myTime,myThid)           tmpLev = myTime + deltaTClock*(jj -jjMx)
160          CALL MNC_CW_SET_UDIM(diag_mnc_bn, 0, myThid)           CALL MNC_CW_SET_UDIM(diag_mnc_bn, -1, myThid)
161          CALL MNC_CW_I_W_S('I',diag_mnc_bn,0,0,'iter',myIter,myThid)           CALL MNC_CW_RL_W_S('D',diag_mnc_bn,0,0,'T',tmpLev,myThid)
162             CALL MNC_CW_SET_UDIM(diag_mnc_bn, 0, myThid)
163             CALL MNC_CW_I_W_S('I',diag_mnc_bn,0,0,'iter',klev,myThid)
164    
165  C       NOTE: at some point it would be a good idea to add a time_bounds  C       NOTE: at some point it would be a good idea to add a time_bounds
166  C       variable that has dimension (2,T) and clearly denotes the  C       variable that has dimension (2,T) and clearly denotes the
167  C       beginning and ending times for each diagnostics period  C       beginning and ending times for each diagnostics period
168    
169          dn(1)(1:NLEN) = dn_blnk(1:NLEN)           dn(1)(1:NLEN) = dn_blnk(1:NLEN)
170          WRITE(dn(1),'(a,i6.6)') 'Zmd', nlevels(listId)           WRITE(dn(1),'(a,i6.6)') 'Zmd', nlevels(listId)
171          dim(1) = nlevels(listId)           dim(1) = nlevels(listId)
172          ib(1)  = 1           ib(1)  = 1
173          ie(1)  = nlevels(listId)           ie(1)  = nlevels(listId)
174    
175          CALL MNC_CW_ADD_GNAME('diag_levels', 1,           CALL MNC_CW_ADD_GNAME('diag_levels', 1,
176       &       dim, dn, ib, ie, myThid)       &        dim, dn, ib, ie, myThid)
177          CALL MNC_CW_ADD_VNAME('diag_levels', 'diag_levels',           CALL MNC_CW_ADD_VNAME('diag_levels', 'diag_levels',
178       &       0,0, myThid)       &        0,0, myThid)
179          CALL MNC_CW_ADD_VATTR_TEXT('diag_levels','description',           CALL MNC_CW_ADD_VATTR_TEXT('diag_levels','description',
180       &       'Idicies of vertical levels within the source arrays',       &        'Idicies of vertical levels within the source arrays',
181       &       myThid)       &        myThid)
182            C     suppress the missing value attribute (iflag = 0)
183          CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,           IF (useMissingValue)
184       &       'diag_levels', levs(1,listId), myThid)       &       CALL MNC_CW_VATTR_MISSING('diag_levels', 0,
185         I       misval_r8, misval_r4, misval_int,
186         I       myThid )
187    
188             CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,
189         &        'diag_levels', levs(1,listId), myThid)
190    
191          CALL MNC_CW_DEL_VNAME('diag_levels', myThid)           CALL MNC_CW_DEL_VNAME('diag_levels', myThid)
192          CALL MNC_CW_DEL_GNAME('diag_levels', myThid)           CALL MNC_CW_DEL_GNAME('diag_levels', myThid)
193    
194  #ifdef DIAG_MNC_COORD_NEEDSWORK  #ifdef DIAG_MNC_COORD_NEEDSWORK
195  C       This part has been placed in an #ifdef because, as its currently  C       This part has been placed in an #ifdef because, as its currently
# Line 148  C       as: Z[uml]td000000 where the 't' Line 205  C       as: Z[uml]td000000 where the 't'
205  C       gdiag(10)  C       gdiag(10)
206    
207  C       Now define:  Zmdxxxxxx, Zudxxxxxx, Zldxxxxxx  C       Now define:  Zmdxxxxxx, Zudxxxxxx, Zldxxxxxx
208          ctmp(1:5) = 'mul  '           ctmp(1:5) = 'mul  '
209          DO i = 1,3           DO i = 1,3
210            dn(1)(1:NLEN) = dn_blnk(1:NLEN)             dn(1)(1:NLEN) = dn_blnk(1:NLEN)
211            WRITE(dn(1),'(3a,i6.6)') 'Z',ctmp(i:i),'d',nlevels(listId)             WRITE(dn(1),'(3a,i6.6)') 'Z',ctmp(i:i),'d',nlevels(listId)
212            CALL MNC_CW_ADD_GNAME(dn(1), 1, dim, dn, ib, ie, myThid)             CALL MNC_CW_ADD_GNAME(dn(1), 1, dim, dn, ib, ie, myThid)
213            CALL MNC_CW_ADD_VNAME(dn(1), dn(1), 0,0, myThid)             CALL MNC_CW_ADD_VNAME(dn(1), dn(1), 0,0, myThid)
214    
215  C         The following three ztmp() loops should eventually be modified  C         The following three ztmp() loops should eventually be modified
216  C         to reflect the fractional nature of levs(j,l) -- they should  C         to reflect the fractional nature of levs(j,l) -- they should
217  C         do something like:  C         do something like:
218  C            ztmp(j) = rC(INT(FLOOR(levs(j,l))))  C            ztmp(j) = rC(INT(FLOOR(levs(j,l))))
219  C                      + ( rC(INT(FLOOR(levs(j,l))))  C                      + ( rC(INT(FLOOR(levs(j,l))))
220  C                          + rC(INT(CEIL(levs(j,l)))) )  C                          + rC(INT(CEIL(levs(j,l)))) )
221  C                        / ( levs(j,l) - FLOOR(levs(j,l)) )  C                        / ( levs(j,l) - FLOOR(levs(j,l)) )
222  C         for averaged levels.  C         for averaged levels.
223            IF (i .EQ. 1) THEN             IF (i .EQ. 1) THEN
224              DO j = 1,nlevels(listId)               DO j = 1,nlevels(listId)
225                ztmp(j) = rC(NINT(levs(j,listId)))                 ztmp(j) = rC(NINT(levs(j,listId)))
226              ENDDO               ENDDO
227              CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',               CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
228       &           'Dimensional coordinate value at the mid point',       &            'Dimensional coordinate value at the mid point',
229       &           myThid)       &            myThid)
230            ELSEIF (i .EQ. 2) THEN             ELSEIF (i .EQ. 2) THEN
231              DO j = 1,nlevels(listId)               DO j = 1,nlevels(listId)
232                ztmp(j) = rF(NINT(levs(j,listId)) + 1)                 ztmp(j) = rF(NINT(levs(j,listId)) + 1)
233              ENDDO               ENDDO
234              CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',               CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
235       &           'Dimensional coordinate value at the upper point',       &            'Dimensional coordinate value at the upper point',
236       &           myThid)       &            myThid)
237            ELSEIF (i .EQ. 3) THEN             ELSEIF (i .EQ. 3) THEN
238              DO j = 1,nlevels(listId)               DO j = 1,nlevels(listId)
239                ztmp(j) = rF(NINT(levs(j,listId)))                 ztmp(j) = rF(NINT(levs(j,listId)))
240              ENDDO               ENDDO
241              CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',               CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',
242       &           'Dimensional coordinate value at the lower point',       &            'Dimensional coordinate value at the lower point',
243       &           myThid)       &            myThid)
244            ENDIF             ENDIF
245            CALL MNC_CW_RS_W('D',diag_mnc_bn,0,0, dn(1), ztmp, myThid)  C     suppress the missing value attribute (iflag = 0)
246            CALL MNC_CW_DEL_VNAME(dn(1), myThid)             IF (useMissingValue)
247            CALL MNC_CW_DEL_GNAME(dn(1), myThid)       &          CALL MNC_CW_VATTR_MISSING(dn(1), 0,
248          ENDDO       I          misval_r8, misval_r4, misval_int,
249         I          myThid )
250               CALL MNC_CW_RS_W('D',diag_mnc_bn,0,0, dn(1), ztmp, myThid)
251               CALL MNC_CW_DEL_VNAME(dn(1), myThid)
252               CALL MNC_CW_DEL_GNAME(dn(1), myThid)
253             ENDDO
254  #endif /*  DIAG_MNC_COORD_NEEDSWORK  */  #endif /*  DIAG_MNC_COORD_NEEDSWORK  */
255    
256        ENDIF         ENDIF
257  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
258    
259        DO md = 1,nfields(listId)  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
260    
261           DO md = 1,nfields(listId)
262          ndId = jdiag(md,listId)          ndId = jdiag(md,listId)
263          parms1 = gdiag(ndId)(1:8)          gcode = gdiag(ndId)(1:10)
264          IF ( idiag(md,listId).NE.0 .AND. parms1(5:5).NE.'D' ) THEN          mate = 0
265            mVec = 0
266            IF ( gcode(5:5).EQ.'C' ) THEN
267    C-      Check for Mate of a Counter Diagnostic
268               mate = hdiag(ndId)
269            ELSEIF ( gcode(1:1).EQ.'U' .OR. gcode(1:1).EQ.'V' ) THEN
270    C-      Check for Mate of a Vector Diagnostic
271               mVec = hdiag(ndId)
272            ENDIF
273            IF ( idiag(md,listId).NE.0 .AND. gcode(5:5).NE.'D' ) THEN
274  C--     Start processing 1 Fld :  C--     Start processing 1 Fld :
275    #ifdef ALLOW_MNC
276             DO ll=1,llMx
277              lm = jj+ll-1
278    #else
279             DO lm=1,averageCycle(listId)
280    #endif
281    
282            ip = ABS(idiag(md,listId))            ip = ABS(idiag(md,listId)) + kdiag(ndId)*(lm-1)
283            im = mdiag(md,listId)            im = mdiag(md,listId)
284              IF (mate.GT.0) im = im + kdiag(mate)*(lm-1)
285              IF (mVec.GT.0) im = im + kdiag(mVec)*(lm-1)
286    
287            IF ( ndiag(ip,1,1).EQ.0 ) THEN            IF ( ndiag(ip,1,1).EQ.0 ) THEN
288  C-        Empty diagnostics case :  C-        Empty diagnostics case :
289    
# Line 210  C-        Empty diagnostics case : Line 292  C-        Empty diagnostics case :
292       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter
293              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
294       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
295              WRITE(msgBuf,'(A,I4,3A,I3,2A)')              WRITE(msgBuf,'(A,I6,3A,I4,2A)')
296       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),
297       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)
298              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
299       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
300              WRITE(msgBuf,'(A,I2,A)')              IF ( averageCycle(listId).GT.1 ) THEN
301       &       '- WARNING -   has not been filled (ndiag=',               WRITE(msgBuf,'(A,2(I3,A))')
302       &       ndiag(ip,1,1), ' )'       &        '- WARNING -   has not been filled (ndiag(lm=',lm,')=',
303         &                                            ndiag(ip,1,1), ' )'
304                ELSE
305                 WRITE(msgBuf,'(A,2(I3,A))')
306         &        '- WARNING -   has not been filled (ndiag=',
307         &                                            ndiag(ip,1,1), ' )'
308                ENDIF
309              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
310       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
311              WRITE(msgBuf,'(A)')              WRITE(msgBuf,'(A)')
# Line 240  C-        Empty diagnostics case : Line 328  C-        Empty diagnostics case :
328            ELSE            ELSE
329  C-        diagnostics is not empty :  C-        diagnostics is not empty :
330    
331              IF ( myThid.EQ.1 ) WRITE(ioUnit,'(A,I3,3A,I8,2A)')              IF ( debugLevel.GE.debLevA .AND. myThid.EQ.1 ) THEN
332                  WRITE(ioUnit,'(A,I6,3A,I8,2A)')
333       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),
334       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)
335                  IF ( mate.GT.0 ) THEN
336              IF ( parms1(5:5).EQ.'C' ) THEN                 WRITE(ioUnit,'(3A,I6,2A)')
 C             Check for Mate of a Counter Diagnostic  
 C             --------------------------------------  
               mate_index = parms1(6:8)  
               READ (mate_index,'(I3)') mate  
               IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,2A)')  
337       &         '       use Counter Mate for  ', cdiag(ndId),       &         '       use Counter Mate for  ', cdiag(ndId),
338       &         '     Diagnostic # ',mate, '  ', cdiag(mate)       &         '     Diagnostic # ',mate, '  ', cdiag(mate)
339                  ELSEIF ( mVec.GT.0 ) THEN
             ELSE  
               mate = 0  
   
 C             Check for Mate of a Vector Diagnostic  
 C             -------------------------------------  
               IF ( parms1(1:1).EQ.'U' .OR. parms1(1:1).EQ.'V' ) THEN  
                 mate_index = parms1(6:8)  
                 READ (mate_index,'(I3)') mVec  
340                  IF ( im.GT.0 .AND. ndiag(MAX(1,im),1,1).GT.0 ) THEN                  IF ( im.GT.0 .AND. ndiag(MAX(1,im),1,1).GT.0 ) THEN
341                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
342       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
343       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
344       &             ' exists '       &             ' exists '
345                  ELSE                  ELSE
346                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
347       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
348       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
349       &             ' not enabled'       &             ' not enabled'
# Line 275  C             -------------------------- Line 351  C             --------------------------
351                ENDIF                ENDIF
352              ENDIF              ENDIF
353    
354              DO bj = myByLo(myThid), myByHi(myThid)              IF ( fflags(listId)(2:2).EQ.'P' ) THEN
355               DO bi = myBxLo(myThid), myBxHi(myThid)  C-       get all the levels (for vertical interpolation)
356                DO k = 1,nlevels(listId)                DO bj = myByLo(myThid), myByHi(myThid)
357                  CALL GETDIAG(                 DO bi = myBxLo(myThid), myBxHi(myThid)
358       I                       levs(k,listId),undef,                  DO k = 1,kdiag(ndId)
359       O                       qtmp1(1-OLx,1-OLy,k,bi,bj),                    tmpLev = k
360       I                       ndId,mate,ip,im,bi,bj,myThid)                    CALL GETDIAG(
361         I                         tmpLev,undef,
362         O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
363         I                         ndId,mate,ip,im,bi,bj,myThid)
364                    ENDDO
365                   ENDDO
366                ENDDO                ENDDO
367               ENDDO              ELSE
368              ENDDO  C-       get only selected levels:
369                  DO bj = myByLo(myThid), myByHi(myThid)
370                   DO bi = myBxLo(myThid), myBxHi(myThid)
371                    DO k = 1,nlevels(listId)
372                      CALL GETDIAG(
373         I                         levs(k,listId),undef,
374         O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
375         I                         ndId,mate,ip,im,bi,bj,myThid)
376                    ENDDO
377                   ENDDO
378                  ENDDO
379                ENDIF
380    
381  C-        end of empty diag / not empty block  C-        end of empty diag / not empty block
382            ENDIF            ENDIF
383    
           nlevsout = nlevels(listId)  
   
384  C-----------------------------------------------------------------------  C-----------------------------------------------------------------------
385  C         Check to see if we need to interpolate before output  C         Check to see if we need to interpolate before output
386  C-----------------------------------------------------------------------  C-----------------------------------------------------------------------
387           IF ( fflags(listId)(2:2).EQ.'P' ) THEN            IF ( fflags(listId)(2:2).EQ.'P' ) THEN
388  C-        Do vertical interpolation:  C-        Do vertical interpolation:
389            CALL DIAGNOSTICS_INTERP_VERT(             IF ( fluidIsAir ) THEN
390       I                     listId, md, ndId, ip, im,  C jmc: for now, this can only work in an atmospheric set-up (fluidIsAir);
391       U                     nlevsout,              CALL DIAGNOSTICS_INTERP_VERT(
392         I                     listId, md, ndId, ip, im, lm,
393       U                     qtmp1,       U                     qtmp1,
394       I                     undef,       I                     undef, myTime, myIter, myThid )
395       I                     myTime, myIter, myThid )             ELSE
396           ENDIF               WRITE(msgBuf,'(2A)') 'DIAGNOSTICS_OUT: ',
397         &         'INTERP_VERT not allowed in this config'
398                 CALL PRINT_ERROR( msgBuf , myThid )
399                 STOP 'ABNORMAL END: S/R DIAGNOSTICS_OUT'
400               ENDIF
401              ENDIF
402    
403  #ifdef ALLOW_MDSIO  C--    Ready to write field "md", element "lm" in averageCycle(listId)
404  C         Prepare for mdsio optionality  
405            IF (diag_mdsio) THEN  C-        write to binary file, using MDSIO pkg:
406              IF (fflags(listId)(1:1) .EQ. 'R') THEN            IF ( diag_mdsio ) THEN
407  C             Force it to be 32-bit precision              nRec = lm + (md-1)*averageCycle(listId)
408                CALL MDSWRITEFIELD_NEW(fn,precFloat32,glf,.FALSE.,  C           default precision for output files
409       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)              prec = writeBinaryPrec
410              ELSEIF (fflags(listId)(1:1) .EQ. 'D') THEN  C           fFlag(1)=R(or D): force it to be 32-bit(or 64) precision
411  C             Force it to be 64-bit precision              IF ( fflags(listId)(1:1).EQ.'R' ) prec = precFloat32
412                CALL MDSWRITEFIELD_NEW(fn,precFloat64,glf,.FALSE.,              IF ( fflags(listId)(1:1).EQ.'D' ) prec = precFloat64
413       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)  C         a hack not to write meta files now: pass -nRec < 0 to MDS_WRITE S/R
414              ELSE              CALL WRITE_REC_LEV_RL(
415  C             This is the old default behavior       I                            fn, prec,
416                CALL MDSWRITEFIELD_NEW(fn,writeBinaryPrec,glf,.FALSE.,       I                            NrMax, 1, nlevels(listId),
417       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)       I                            qtmp1, -nRec, myIter, myThid )
             ENDIF  
418            ENDIF            ENDIF
 #endif /*  ALLOW_MDSIO  */  
419    
420  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
421            IF (useMNC .AND. diag_mnc) THEN            IF (useMNC .AND. diag_mnc) THEN
# Line 347  C           XY dimensions Line 441  C           XY dimensions
441              dim(2)       = sNy + 2*OLy              dim(2)       = sNy + 2*OLy
442              ib(1)        = OLx + 1              ib(1)        = OLx + 1
443              ib(2)        = OLy + 1              ib(2)        = OLy + 1
444              IF (gdiag(ndId)(2:2) .EQ. 'M') THEN              IF (gdiag(ndId)(2:2) .EQ. 'M') THEN
445                dn(1)(1:2) = 'X'                dn(1)(1:2) = 'X'
446                ie(1)      = OLx + sNx                ie(1)      = OLx + sNx
447                dn(2)(1:2) = 'Y'                dn(2)(1:2) = 'Y'
# Line 368  C           XY dimensions Line 462  C           XY dimensions
462                dn(2)(1:3) = 'Yp1'                dn(2)(1:3) = 'Yp1'
463                ie(2)      = OLy + sNy + 1                ie(2)      = OLy + sNy + 1
464              ENDIF              ENDIF
465                
466  C           Z is special since it varies  C           Z is special since it varies
467              WRITE(dn(3),'(a,i6.6)') 'Zd', nlevsout              WRITE(dn(3),'(a,i6.6)') 'Zd', nlevels(listId)
468              IF ( (gdiag(ndId)(10:10) .EQ. 'R')              IF ( (gdiag(ndId)(10:10) .EQ. 'R')
469       &           .AND. (gdiag(ndId)(9:9) .EQ. 'M') ) THEN       &           .AND. (gdiag(ndId)(9:9) .EQ. 'M') ) THEN
470                WRITE(dn(3),'(a,i6.6)') 'Zmd', nlevsout                WRITE(dn(3),'(a,i6.6)') 'Zmd', nlevels(listId)
471              ENDIF              ENDIF
472              IF ( (gdiag(ndId)(10:10) .EQ. 'R')              IF ( (gdiag(ndId)(10:10) .EQ. 'R')
473       &           .AND. (gdiag(ndId)(9:9) .EQ. 'L') ) THEN       &           .AND. (gdiag(ndId)(9:9) .EQ. 'L') ) THEN
474                WRITE(dn(3),'(a,i6.6)') 'Zld', nlevsout                WRITE(dn(3),'(a,i6.6)') 'Zld', nlevels(listId)
475              ENDIF              ENDIF
476              IF ( (gdiag(ndId)(10:10) .EQ. 'R')              IF ( (gdiag(ndId)(10:10) .EQ. 'R')
477       &           .AND. (gdiag(ndId)(9:9) .EQ. 'U') ) THEN       &           .AND. (gdiag(ndId)(9:9) .EQ. 'U') ) THEN
478                WRITE(dn(3),'(a,i6.6)') 'Zud', nlevsout                WRITE(dn(3),'(a,i6.6)') 'Zud', nlevels(listId)
479              ENDIF              ENDIF
480              dim(3) = Nr+Nrphys              dim(3) = NrMax
481              ib(3)  = 1              ib(3)  = 1
482              ie(3)  = nlevsout              ie(3)  = nlevels(listId)
483    
484  C           Time dimension  C           Time dimension
485              dn(4)(1:1) = 'T'              dn(4)(1:1) = 'T'
# Line 393  C           Time dimension Line 487  C           Time dimension
487              ib(4)  = 1              ib(4)  = 1
488              ie(4)  = 1              ie(4)  = 1
489    
490              CALL MNC_CW_ADD_GNAME(d_cw_name, 4,              CALL MNC_CW_ADD_GNAME(d_cw_name, 4,
491       &             dim, dn, ib, ie, myThid)       &             dim, dn, ib, ie, myThid)
492              CALL MNC_CW_ADD_VNAME(cdiag(ndId), d_cw_name,              CALL MNC_CW_ADD_VNAME(cdiag(ndId), d_cw_name,
493       &             4,5, myThid)       &             4,5, myThid)
494              CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'description',              CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'description',
495       &             tdiag(ndId),myThid)       &             tdiag(ndId),myThid)
496              CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'units',              CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'units',
497       &             udiag(ndId),myThid)       &             udiag(ndId),myThid)
498    
499  C           Per the observations of Baylor, this has been commented out  C     Missing values only for scalar diagnostics at mass points (so far)
500  C           until we have code that can write missing_value attributes              useMisValForThisDiag = useMissingValue
501  C           in a way thats compatible with most of the more popular       &           .AND.gdiag(ndId)(1:2).EQ.'SM'
502  C           netCDF tools including ferret.  Using all-zeros completely              IF ( useMisValForThisDiag ) THEN
503  C           breaks ferret.  C     assign missing values and set flag for adding the netCDF atttibute
504                 CALL MNC_CW_VATTR_MISSING(cdiag(ndId), 2,
505  C           CALL MNC_CW_ADD_VATTR_DBL(cdiag(ndId),'missing_value',       I            misval_r8, misval_r4, misval_int,
506  C           &             0.0 _d 0,myThid)       I            myThid )
507    C     and now use the missing values for masking out the land points
508              IF ( ( (writeBinaryPrec .EQ. precFloat32)               DO bj = myByLo(myThid), myByHi(myThid)
509       &           .AND. (fflags(listId)(1:1) .NE. 'D')                DO bi = myBxLo(myThid), myBxHi(myThid)
510       &           .AND. (fflags(listId)(1:1) .NE. 'R') )                 DO k = 1,nlevels(listId)
511       &           .OR. (fflags(listId)(1:1) .EQ. 'R')) THEN                  klev = NINT(levs(k,listId))
512                    DO j = 1-OLy,sNy+OLy
513                     DO i = 1-OLx,sNx+OLx
514                      IF ( _hFacC(I,J,klev,bi,bj) .EQ. 0. )
515         &                 qtmp1(i,j,k,bi,bj) = misvalLoc
516                     ENDDO
517                    ENDDO
518                   ENDDO
519                  ENDDO
520                 ENDDO
521                ELSE
522    C     suppress the missing value attribute (iflag = 0)
523    C     Note: We have to call the following subroutine for each mnc that has
524    C     been created "on the fly" by mnc_cw_add_vname and will be deleted
525    C     by mnc_cw_del_vname, because all of these variables use the same
526    C     identifier so that mnc_cw_vfmv(indv) needs to be overwritten for
527    C     each of these variables
528                 CALL MNC_CW_VATTR_MISSING(cdiag(ndId), 0,
529         I            misval_r8, misval_r4, misval_int,
530         I            myThid )
531                ENDIF
532    
533                IF (  ((writeBinaryPrec .EQ. precFloat32)
534         &            .AND. (fflags(listId)(1:1) .NE. 'D'))
535         &             .OR. (fflags(listId)(1:1) .EQ. 'R') ) THEN
536                CALL MNC_CW_RL_W('R',diag_mnc_bn,0,0,                CALL MNC_CW_RL_W('R',diag_mnc_bn,0,0,
537       &             cdiag(ndId), qtmp1, myThid)       &             cdiag(ndId), qtmp1, myThid)
538              ELSEIF ( (writeBinaryPrec .EQ. precFloat64)              ELSEIF ( (writeBinaryPrec .EQ. precFloat64)
539       &             .OR. (fflags(listId)(1:1) .EQ. 'D') ) THEN       &             .OR. (fflags(listId)(1:1) .EQ. 'D') ) THEN
540                CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,                CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,
541       &             cdiag(ndId), qtmp1, myThid)       &             cdiag(ndId), qtmp1, myThid)
542              ENDIF              ENDIF
543                
544              CALL MNC_CW_DEL_VNAME(cdiag(ndId), myThid)              CALL MNC_CW_DEL_VNAME(cdiag(ndId), myThid)
545              CALL MNC_CW_DEL_GNAME(d_cw_name, myThid)              CALL MNC_CW_DEL_GNAME(d_cw_name, myThid)
546    
# Line 431  C           &             0.0 _d 0,myThi Line 549  C           &             0.0 _d 0,myThi
549            ENDIF            ENDIF
550  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
551    
552    C--      end loop on lm (or ll if ALLOW_MNC) counter
553             ENDDO
554  C--     end of Processing Fld # md  C--     end of Processing Fld # md
555          ENDIF          ENDIF
556           ENDDO
557    
558    #ifdef ALLOW_MNC
559    C--   end loop on jj counter
560        ENDDO        ENDDO
561    #endif
562    
563    #ifdef ALLOW_MDSIO
564          IF (diag_mdsio) THEN
565    C-    Note: temporary: since it's a pain to add more arguments to
566    C     all MDSIO S/R, uses instead this specific S/R to write only
567    C     meta files but with more informations in it.
568                glf = globalFiles
569                nRec = nfields(listId)*averageCycle(listId)
570                CALL MDS_WR_METAFILES(fn, prec, glf, .FALSE.,
571         &              0, 0, nlevels(listId), ' ',
572         &              nfields(listId), flds(1,listId), 1, myTime,
573         &              nRec, myIter, myThid)
574          ENDIF
575    #endif /*  ALLOW_MDSIO  */
576    
577        RETURN        RETURN
578        END        END
Legend:
Removed from v.1.28  
changed lines
  Added in v.1.41
  ViewVC Help
Powered by ViewVC 1.1.22