/[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.51 by jmc, Sun Jun 12 13:58:33 2011 UTC
# Line 8  CBOP 0 Line 8  CBOP 0
8  C     !ROUTINE: DIAGNOSTICS_OUT  C     !ROUTINE: DIAGNOSTICS_OUT
9    
10  C     !INTERFACE:  C     !INTERFACE:
11        SUBROUTINE  DIAGNOSTICS_OUT(        SUBROUTINE DIAGNOSTICS_OUT(
12       I     listId,       I     listId,
      I     myIter,  
13       I     myTime,       I     myTime,
14         I     myIter,
15       I     myThid )       I     myThid )
16    
17  C     !DESCRIPTION:  C     !DESCRIPTION:
# 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     bi,bj :: tile indices
52    C     lm    :: loop index (averageCycle)
53  C     md    :: field number in the list "listId".  C     md    :: field number in the list "listId".
54  C     ndId  :: diagnostics  Id number (in available diagnostics list)  C     ndId  :: diagnostics  Id number (in available diagnostics list)
55  C     mate  :: counter mate Id number (in available diagnostics list)  C     mate  :: counter mate Id number (in available diagnostics list)
56  C     ip    :: diagnostics  pointer to storage array  C     ip    :: diagnostics  pointer to storage array
57  C     im    :: counter-mate pointer to storage array  C     im    :: counter-mate pointer to storage array
58        INTEGER i, j, k  C     nLevOutp :: number of levels to write in output file
59    C
60    C--   COMMON /LOCAL_DIAGNOSTICS_OUT/ local common block (for multi-threaded)
61    C     qtmp1 :: thread-shared temporary array (needs to be in common block):
62    C              to write a diagnostic field to file, copy it first from (big)
63    C              diagnostic storage qdiag into it.
64          COMMON /LOCAL_DIAGNOSTICS_OUT/ qtmp1
65          _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,NrMax,nSx,nSy)
66    
67          INTEGER i, j, k, lm
68        INTEGER bi, bj        INTEGER bi, bj
69        INTEGER md, ndId, ip, im        INTEGER md, ndId, ip, im
70        INTEGER mate, mVec        INTEGER mate, mVec
71        CHARACTER*8 parms1        CHARACTER*10 gcode
72        CHARACTER*3 mate_index        _RL undef
73        _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+Nrphys,nSx,nSy)        _RL tmpLev
74        _RL undef, getcon        INTEGER iLen
75        EXTERNAL getcon        INTEGER nLevOutp
       INTEGER ILNBLNK  
       EXTERNAL ILNBLNK  
       INTEGER ilen  
       INTEGER nlevsout  
76    
77        INTEGER ioUnit        INTEGER ioUnit
78        CHARACTER*(MAX_LEN_FNAM) fn        CHARACTER*(MAX_LEN_FNAM) fn
79        CHARACTER*(MAX_LEN_MBUF) suff        CHARACTER*(MAX_LEN_MBUF) suff
80        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
81          INTEGER prec, nRec, nTimRec
82          _RL     timeRec(2)
83    #ifdef ALLOW_MDSIO
84        LOGICAL glf        LOGICAL glf
85    #endif
86  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
87        INTEGER ii        INTEGER ll, llMx, jj, jjMx
88        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn
89        INTEGER CW_DIMS, NLEN        LOGICAL useMissingValue
90        PARAMETER ( CW_DIMS = 10 )        REAL*8 misValLoc
       PARAMETER ( NLEN    = 80 )  
       INTEGER dim(CW_DIMS), ib(CW_DIMS), ie(CW_DIMS)  
       CHARACTER*(NLEN) dn(CW_DIMS)  
       CHARACTER*(NLEN) d_cw_name  
       CHARACTER*(NLEN) dn_blnk  
 #ifdef DIAG_MNC_COORD_NEEDSWORK  
       CHARACTER*(5) ctmp  
       _RS ztmp(Nr+Nrphys)  
 #endif  
91  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
92    
93  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
94    
95    C---  set file properties
96        ioUnit= standardMessageUnit        ioUnit= standardMessageUnit
97        undef = getcon('UNDEF')        undef = UNSET_RL
98        glf = globalFiles  #ifdef ALLOW_FIZHI
99          IF ( useFIZHI ) undef = getcon('UNDEF')
100    #endif
101        WRITE(suff,'(I10.10)') myIter        WRITE(suff,'(I10.10)') myIter
102        ilen = ILNBLNK(fnames(listId))        iLen = ILNBLNK(fnames(listId))
103        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:ilen),'.',suff(1:10)        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:iLen),'.',suff(1:10)
104    C-    for now, if integrate vertically, output field has just 1 level:
105          nLevOutp = nlevels(listId)
106          IF ( fflags(listId)(2:2).EQ.'I' ) nLevOutp = 1
107    
108    C--   Set time information:
109          IF ( freq(listId).LT.0. ) THEN
110    C-    Snap-shot: store a unique time (which is consistent with State-Var timing)
111            nTimRec = 1
112            timeRec(1) = myTime
113          ELSE
114    C-    Time-average: store the 2 edges of the time-averaging interval.
115    C      this time is consitent with intermediate Var (i.e., non-state, e.g, flux,
116    C      tendencies) timing. For State-Var, this is shifted by + halt time-step.
117            nTimRec = 2
118    
119    C-    end of time-averaging interval:
120            timeRec(2) = myTime
121    
122    C-    begining of time-averaging interval:
123    c       timeRec(1) = myTime - freq(listId)
124    C     a) find the time of the previous multiple of output freq:
125            timeRec(1) = myTime-deltaTClock*0.5 _d 0
126            timeRec(1) = (timeRec(1)-phase(listId))/freq(listId)
127            i = INT( timeRec(1) )
128            IF ( timeRec(1).LT.0. ) THEN
129              tmpLev = FLOAT(i)
130              IF ( timeRec(1).NE.tmpLev ) i = i - 1
131            ENDIF
132            timeRec(1) = phase(listId) + freq(listId)*FLOAT(i)
133    c       if ( listId.eq.2 ) write(0,*) 'f',i,timeRec(1)/deltaTClock
134            timeRec(1) = MAX( timeRec(1), startTime )
135    
136    C     b) round off to nearest multiple of time-step:
137            timeRec(1) = (timeRec(1)-baseTime)/deltaTClock
138            i = NINT( timeRec(1) )
139    C     if just half way, NINT will return the next time-step: correct this
140            tmpLev = FLOAT(i) - 0.5 _d 0
141            IF ( timeRec(1).EQ.tmpLev ) i = i - 1
142            timeRec(1) = baseTime + deltaTClock*FLOAT(i)
143    c       if ( listId.eq.2 ) write(0,*) i,timeRec(1)/deltaTClock
144          ENDIF
145    C--   Convert time to iteration number (debug)
146    c     DO i=1,nTimRec
147    c       timeRec(i) = timeRec(i)/deltaTClock
148    c     ENDDO
149    
150  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
151    C-- this is a trick to reverse the order of the loops on md (= field)
152    C   and lm (= averagePeriod): binary output: lm loop inside md loop ;
153    C                                 mnc ouput: md loop inside lm loop.
154        IF (useMNC .AND. diag_mnc) THEN        IF (useMNC .AND. diag_mnc) THEN
155          DO i = 1,MAX_LEN_FNAM          jjMx = averageCycle(listId)
156            diag_mnc_bn(i:i) = ' '          llMx = 1
157          ENDDO        ELSE
158          DO i = 1,NLEN          jjMx = 1
159            dn_blnk(i:i) = ' '          llMx = averageCycle(listId)
         ENDDO  
         WRITE( diag_mnc_bn, '(A)' ) fnames(listId)(1:ilen)  
   
 C       Update the record dimension by writing the iteration number  
         CALL MNC_CW_SET_UDIM(diag_mnc_bn, -1, myThid)  
         CALL MNC_CW_RL_W_S('D',diag_mnc_bn,0,0,'T',myTime,myThid)  
         CALL MNC_CW_SET_UDIM(diag_mnc_bn, 0, myThid)  
         CALL MNC_CW_I_W_S('I',diag_mnc_bn,0,0,'iter',myIter,myThid)  
   
 C       NOTE: at some point it would be a good idea to add a time_bounds  
 C       variable that has dimension (2,T) and clearly denotes the  
 C       beginning and ending times for each diagnostics period  
   
         dn(1)(1:NLEN) = dn_blnk(1:NLEN)  
         WRITE(dn(1),'(a,i6.6)') 'Zmd', nlevels(listId)  
         dim(1) = nlevels(listId)  
         ib(1)  = 1  
         ie(1)  = nlevels(listId)  
   
         CALL MNC_CW_ADD_GNAME('diag_levels', 1,  
      &       dim, dn, ib, ie, myThid)  
         CALL MNC_CW_ADD_VNAME('diag_levels', 'diag_levels',  
      &       0,0, myThid)  
         CALL MNC_CW_ADD_VATTR_TEXT('diag_levels','description',  
      &       'Idicies of vertical levels within the source arrays',  
      &       myThid)  
           
         CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,  
      &       'diag_levels', levs(1,listId), myThid)  
   
         CALL MNC_CW_DEL_VNAME('diag_levels', myThid)  
         CALL MNC_CW_DEL_GNAME('diag_levels', myThid)  
   
 #ifdef DIAG_MNC_COORD_NEEDSWORK  
 C       This part has been placed in an #ifdef because, as its currently  
 C       written, it will only work with variables defined on a dynamics  
 C       grid.  As we start using diagnostics for physics grids, ice  
 C       levels, land levels, etc. the different vertical coordinate  
 C       dimensions will have to be taken into account.  
   
 C       20051021 JMC & EH3 : We need to extend this so that a few  
 C       variables each defined on different grids do not have the same  
 C       vertical dimension names so we should be using a pattern such  
 C       as: Z[uml]td000000 where the 't' is the type as specified by  
 C       gdiag(10)  
   
 C       Now define:  Zmdxxxxxx, Zudxxxxxx, Zldxxxxxx  
         ctmp(1:5) = 'mul  '  
         DO i = 1,3  
           dn(1)(1:NLEN) = dn_blnk(1:NLEN)  
           WRITE(dn(1),'(3a,i6.6)') 'Z',ctmp(i:i),'d',nlevels(listId)  
           CALL MNC_CW_ADD_GNAME(dn(1), 1, dim, dn, ib, ie, myThid)  
           CALL MNC_CW_ADD_VNAME(dn(1), dn(1), 0,0, myThid)  
   
 C         The following three ztmp() loops should eventually be modified  
 C         to reflect the fractional nature of levs(j,l) -- they should  
 C         do something like:  
 C            ztmp(j) = rC(INT(FLOOR(levs(j,l))))  
 C                      + ( rC(INT(FLOOR(levs(j,l))))  
 C                          + rC(INT(CEIL(levs(j,l)))) )  
 C                        / ( levs(j,l) - FLOOR(levs(j,l)) )  
 C         for averaged levels.  
           IF (i .EQ. 1) THEN  
             DO j = 1,nlevels(listId)  
               ztmp(j) = rC(NINT(levs(j,listId)))  
             ENDDO  
             CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',  
      &           'Dimensional coordinate value at the mid point',  
      &           myThid)  
           ELSEIF (i .EQ. 2) THEN  
             DO j = 1,nlevels(listId)  
               ztmp(j) = rF(NINT(levs(j,listId)) + 1)  
             ENDDO  
             CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',  
      &           'Dimensional coordinate value at the upper point',  
      &           myThid)  
           ELSEIF (i .EQ. 3) THEN  
             DO j = 1,nlevels(listId)  
               ztmp(j) = rF(NINT(levs(j,listId)))  
             ENDDO  
             CALL MNC_CW_ADD_VATTR_TEXT(dn(1),'description',  
      &           'Dimensional coordinate value at the lower point',  
      &           myThid)  
           ENDIF  
           CALL MNC_CW_RS_W('D',diag_mnc_bn,0,0, dn(1), ztmp, myThid)  
           CALL MNC_CW_DEL_VNAME(dn(1), myThid)  
           CALL MNC_CW_DEL_GNAME(dn(1), myThid)  
         ENDDO  
 #endif /*  DIAG_MNC_COORD_NEEDSWORK  */  
   
160        ENDIF        ENDIF
161          DO jj=1,jjMx
162    
163           IF (useMNC .AND. diag_mnc) THEN
164             misValLoc = undef
165             IF ( misvalFlt(listId).NE.UNSET_RL )
166         &        misValLoc = misvalFlt(listId)
167             CALL DIAGNOSTICS_MNC_SET(
168         I                    nLevOutp, listId, jj,
169         O                    diag_mnc_bn, useMissingValue,
170         I                    misValLoc, myTime, myIter, myThid )
171           ENDIF
172  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
173    
174        DO md = 1,nfields(listId)  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
175    
176           DO md = 1,nfields(listId)
177          ndId = jdiag(md,listId)          ndId = jdiag(md,listId)
178          parms1 = gdiag(ndId)(1:8)          gcode = gdiag(ndId)(1:10)
179          IF ( idiag(md,listId).NE.0 .AND. parms1(5:5).NE.'D' ) THEN          mate = 0
180            mVec = 0
181            IF ( gcode(5:5).EQ.'C' ) THEN
182    C-      Check for Mate of a Counter Diagnostic
183               mate = hdiag(ndId)
184            ELSEIF ( gcode(1:1).EQ.'U' .OR. gcode(1:1).EQ.'V' ) THEN
185    C-      Check for Mate of a Vector Diagnostic
186               mVec = hdiag(ndId)
187            ENDIF
188            IF ( idiag(md,listId).NE.0 .AND. gcode(5:5).NE.'D' ) THEN
189  C--     Start processing 1 Fld :  C--     Start processing 1 Fld :
190    #ifdef ALLOW_MNC
191             DO ll=1,llMx
192              lm = jj+ll-1
193    #else
194             DO lm=1,averageCycle(listId)
195    #endif
196    
197            ip = ABS(idiag(md,listId))            ip = ABS(idiag(md,listId)) + kdiag(ndId)*(lm-1)
198            im = mdiag(md,listId)            im = mdiag(md,listId)
199              IF (mate.GT.0) im = im + kdiag(mate)*(lm-1)
200              IF (mVec.GT.0) im = im + kdiag(mVec)*(lm-1)
201    
202            IF ( ndiag(ip,1,1).EQ.0 ) THEN            IF ( ndiag(ip,1,1).EQ.0 ) THEN
203  C-        Empty diagnostics case :  C-        Empty diagnostics case :
204    
# Line 210  C-        Empty diagnostics case : Line 207  C-        Empty diagnostics case :
207       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter
208              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
209       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
210              WRITE(msgBuf,'(A,I4,3A,I3,2A)')              WRITE(msgBuf,'(A,I6,3A,I4,2A)')
211       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),
212       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)
213              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
214       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
215              WRITE(msgBuf,'(A,I2,A)')              IF ( averageCycle(listId).GT.1 ) THEN
216       &       '- WARNING -   has not been filled (ndiag=',               WRITE(msgBuf,'(A,2(I3,A))')
217       &       ndiag(ip,1,1), ' )'       &        '- WARNING -   has not been filled (ndiag(lm=',lm,')=',
218         &                                            ndiag(ip,1,1), ' )'
219                ELSE
220                 WRITE(msgBuf,'(A,2(I3,A))')
221         &        '- WARNING -   has not been filled (ndiag=',
222         &                                            ndiag(ip,1,1), ' )'
223                ENDIF
224              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
225       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
226              WRITE(msgBuf,'(A)')              WRITE(msgBuf,'(A)')
# Line 227  C-        Empty diagnostics case : Line 230  C-        Empty diagnostics case :
230              _END_MASTER( myThid )              _END_MASTER( myThid )
231              DO bj = myByLo(myThid), myByHi(myThid)              DO bj = myByLo(myThid), myByHi(myThid)
232                DO bi = myBxLo(myThid), myBxHi(myThid)                DO bi = myBxLo(myThid), myBxHi(myThid)
233                  DO k = 1,nlevels(listId)                  DO k = 1,nLevOutp
234                    DO j = 1-OLy,sNy+OLy                    DO j = 1-OLy,sNy+OLy
235                      DO i = 1-OLx,sNx+OLx                      DO i = 1-OLx,sNx+OLx
236                        qtmp1(i,j,k,bi,bj) = 0. _d 0                        qtmp1(i,j,k,bi,bj) = 0. _d 0
# Line 240  C-        Empty diagnostics case : Line 243  C-        Empty diagnostics case :
243            ELSE            ELSE
244  C-        diagnostics is not empty :  C-        diagnostics is not empty :
245    
246              IF ( myThid.EQ.1 ) WRITE(ioUnit,'(A,I3,3A,I8,2A)')              IF ( debugLevel.GE.debLevB .AND. myThid.EQ.1 ) THEN
247                  WRITE(ioUnit,'(A,I6,3A,I8,2A)')
248       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),
249       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)
250                  IF ( mate.GT.0 ) THEN
251              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)')  
252       &         '       use Counter Mate for  ', cdiag(ndId),       &         '       use Counter Mate for  ', cdiag(ndId),
253       &         '     Diagnostic # ',mate, '  ', cdiag(mate)       &         '     Diagnostic # ',mate, '  ', cdiag(mate)
254                  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  
255                  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
256                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
257       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
258       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
259       &             ' exists '       &             ' exists '
260                  ELSE                  ELSE
261                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
262       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
263       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
264       &             ' not enabled'       &             ' not enabled'
# Line 275  C             -------------------------- Line 266  C             --------------------------
266                ENDIF                ENDIF
267              ENDIF              ENDIF
268    
269              DO bj = myByLo(myThid), myByHi(myThid)              IF ( fflags(listId)(2:2).NE.' ' ) THEN
270               DO bi = myBxLo(myThid), myBxHi(myThid)  C-       get all the levels (for vertical post-processing)
271                DO k = 1,nlevels(listId)                DO bj = myByLo(myThid), myByHi(myThid)
272                  CALL GETDIAG(                 DO bi = myBxLo(myThid), myBxHi(myThid)
273       I                       levs(k,listId),undef,                  DO k = 1,kdiag(ndId)
274       O                       qtmp1(1-OLx,1-OLy,k,bi,bj),                    tmpLev = k
275       I                       ndId,mate,ip,im,bi,bj,myThid)                    CALL GETDIAG(
276         I                         tmpLev,undef,
277         O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
278         I                         ndId,mate,ip,im,bi,bj,myThid)
279                    ENDDO
280                   ENDDO
281                ENDDO                ENDDO
282               ENDDO              ELSE
283              ENDDO  C-       get only selected levels:
284                  DO bj = myByLo(myThid), myByHi(myThid)
285  C-        end of empty diag / not empty block                 DO bi = myBxLo(myThid), myBxHi(myThid)
286            ENDIF                  DO k = 1,nlevels(listId)
287                      CALL GETDIAG(
288            nlevsout = nlevels(listId)       I                         levs(k,listId),undef,
289         O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
290         I                         ndId,mate,ip,im,bi,bj,myThid)
291                    ENDDO
292                   ENDDO
293                  ENDDO
294                ENDIF
295    
296  C-----------------------------------------------------------------------  C-----------------------------------------------------------------------
297  C         Check to see if we need to interpolate before output  C--     Apply specific post-processing (e.g., interpolate) before output
298  C-----------------------------------------------------------------------  C-----------------------------------------------------------------------
299           IF ( fflags(listId)(2:2).EQ.'P' ) THEN              IF ( fflags(listId)(2:2).EQ.'P' ) THEN
300  C-        Do vertical interpolation:  C-          Do vertical interpolation:
301            CALL DIAGNOSTICS_INTERP_VERT(               IF ( fluidIsAir ) THEN
302       I                     listId, md, ndId, ip, im,  C jmc: for now, this can only work in an atmospheric set-up (fluidIsAir);
303       U                     nlevsout,                CALL DIAGNOSTICS_INTERP_VERT(
304       U                     qtmp1,       I                         listId, md, ndId, ip, im, lm,
305       I                     undef,       U                         qtmp1,
306       I                     myTime, myIter, myThid )       I                         undef, myTime, myIter, myThid )
307           ENDIF               ELSE
308                   WRITE(msgBuf,'(2A)') 'DIAGNOSTICS_OUT: ',
309  #ifdef ALLOW_MDSIO       &           'INTERP_VERT not allowed in this config'
310  C         Prepare for mdsio optionality                 CALL PRINT_ERROR( msgBuf , myThid )
311            IF (diag_mdsio) THEN                 STOP 'ABNORMAL END: S/R DIAGNOSTICS_OUT'
312              IF (fflags(listId)(1:1) .EQ. 'R') THEN               ENDIF
313  C             Force it to be 32-bit precision              ENDIF
314                CALL MDSWRITEFIELD_NEW(fn,precFloat32,glf,.FALSE.,              IF ( fflags(listId)(2:2).EQ.'I' ) THEN
315       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)  C-          Integrate vertically: for now, output field has just 1 level:
316              ELSEIF (fflags(listId)(1:1) .EQ. 'D') THEN                CALL DIAGNOSTICS_SUM_LEVELS(
317  C             Force it to be 64-bit precision       I                         listId, md, ndId, ip, im, lm,
318                CALL MDSWRITEFIELD_NEW(fn,precFloat64,glf,.FALSE.,       U                         qtmp1,
319       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)       I                         undef, myTime, myIter, myThid )
             ELSE  
 C             This is the old default behavior  
               CALL MDSWRITEFIELD_NEW(fn,writeBinaryPrec,glf,.FALSE.,  
      &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)  
320              ENDIF              ENDIF
           ENDIF  
 #endif /*  ALLOW_MDSIO  */  
   
 #ifdef ALLOW_MNC  
           IF (useMNC .AND. diag_mnc) THEN  
   
             _BEGIN_MASTER( myThid )  
321    
322              DO ii = 1,CW_DIMS  C--     End of empty diag / not-empty diag block
323                d_cw_name(1:NLEN) = dn_blnk(1:NLEN)            ENDIF
               dn(ii)(1:NLEN) = dn_blnk(1:NLEN)  
             ENDDO  
324    
325  C           Note that the "d_cw_name" variable is a hack that hides a  C--     Ready to write field "md", element "lm" in averageCycle(listId)
 C           subtlety within MNC.  Basically, each MNC-wrapped file is  
 C           caching its own concept of what each "grid name" (that is, a  
 C           dimension group name) means.  So one cannot re-use the same  
 C           "grid" name for different collections of dimensions within a  
 C           given file.  By appending the "ndId" values to each name, we  
 C           guarantee uniqueness within each MNC-produced file.  
             WRITE(d_cw_name,'(a,i6.6)') 'd_cw_',ndId  
   
 C           XY dimensions  
             dim(1)       = sNx + 2*OLx  
             dim(2)       = sNy + 2*OLy  
             ib(1)        = OLx + 1  
             ib(2)        = OLy + 1  
             IF (gdiag(ndId)(2:2) .EQ. 'M') THEN  
               dn(1)(1:2) = 'X'  
               ie(1)      = OLx + sNx  
               dn(2)(1:2) = 'Y'  
               ie(2)      = OLy + sNy  
             ELSEIF (gdiag(ndId)(2:2) .EQ. 'U') THEN  
               dn(1)(1:3) = 'Xp1'  
               ie(1)      = OLx + sNx + 1  
               dn(2)(1:2) = 'Y'  
               ie(2)      = OLy + sNy  
             ELSEIF (gdiag(ndId)(2:2) .EQ. 'V') THEN  
               dn(1)(1:2) = 'X'  
               ie(1)      = OLx + sNx  
               dn(2)(1:3) = 'Yp1'  
               ie(2)      = OLy + sNy + 1  
             ELSEIF (gdiag(ndId)(2:2) .EQ. 'Z') THEN  
               dn(1)(1:3) = 'Xp1'  
               ie(1)      = OLx + sNx + 1  
               dn(2)(1:3) = 'Yp1'  
               ie(2)      = OLy + sNy + 1  
             ENDIF  
               
 C           Z is special since it varies  
             WRITE(dn(3),'(a,i6.6)') 'Zd', nlevsout  
             IF ( (gdiag(ndId)(10:10) .EQ. 'R')  
      &           .AND. (gdiag(ndId)(9:9) .EQ. 'M') ) THEN  
               WRITE(dn(3),'(a,i6.6)') 'Zmd', nlevsout  
             ENDIF  
             IF ( (gdiag(ndId)(10:10) .EQ. 'R')  
      &           .AND. (gdiag(ndId)(9:9) .EQ. 'L') ) THEN  
               WRITE(dn(3),'(a,i6.6)') 'Zld', nlevsout  
             ENDIF  
             IF ( (gdiag(ndId)(10:10) .EQ. 'R')  
      &           .AND. (gdiag(ndId)(9:9) .EQ. 'U') ) THEN  
               WRITE(dn(3),'(a,i6.6)') 'Zud', nlevsout  
             ENDIF  
             dim(3) = Nr+Nrphys  
             ib(3)  = 1  
             ie(3)  = nlevsout  
   
 C           Time dimension  
             dn(4)(1:1) = 'T'  
             dim(4) = -1  
             ib(4)  = 1  
             ie(4)  = 1  
   
             CALL MNC_CW_ADD_GNAME(d_cw_name, 4,  
      &             dim, dn, ib, ie, myThid)  
             CALL MNC_CW_ADD_VNAME(cdiag(ndId), d_cw_name,  
      &             4,5, myThid)  
             CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'description',  
      &             tdiag(ndId),myThid)  
             CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'units',  
      &             udiag(ndId),myThid)  
   
 C           Per the observations of Baylor, this has been commented out  
 C           until we have code that can write missing_value attributes  
 C           in a way thats compatible with most of the more popular  
 C           netCDF tools including ferret.  Using all-zeros completely  
 C           breaks ferret.  
   
 C           CALL MNC_CW_ADD_VATTR_DBL(cdiag(ndId),'missing_value',  
 C           &             0.0 _d 0,myThid)  
   
             IF ( ( (writeBinaryPrec .EQ. precFloat32)  
      &           .AND. (fflags(listId)(1:1) .NE. 'D')  
      &           .AND. (fflags(listId)(1:1) .NE. 'R') )  
      &           .OR. (fflags(listId)(1:1) .EQ. 'R')) THEN  
               CALL MNC_CW_RL_W('R',diag_mnc_bn,0,0,  
      &             cdiag(ndId), qtmp1, myThid)  
             ELSEIF ( (writeBinaryPrec .EQ. precFloat64)  
      &             .OR. (fflags(listId)(1:1) .EQ. 'D') ) THEN  
               CALL MNC_CW_RL_W('D',diag_mnc_bn,0,0,  
      &             cdiag(ndId), qtmp1, myThid)  
             ENDIF  
               
             CALL MNC_CW_DEL_VNAME(cdiag(ndId), myThid)  
             CALL MNC_CW_DEL_GNAME(d_cw_name, myThid)  
326    
327              _END_MASTER( myThid )  C-        write to binary file, using MDSIO pkg:
328              IF ( diag_mdsio ) THEN
329                nRec = lm + (md-1)*averageCycle(listId)
330    C           default precision for output files
331                prec = writeBinaryPrec
332    C           fFlag(1)=R(or D): force it to be 32-bit(or 64) precision
333                IF ( fflags(listId)(1:1).EQ.'R' ) prec = precFloat32
334                IF ( fflags(listId)(1:1).EQ.'D' ) prec = precFloat64
335    C         a hack not to write meta files now: pass -nRec < 0 to MDS_WRITE S/R
336                CALL WRITE_REC_LEV_RL(
337         I                            fn, prec,
338         I                            NrMax, 1, nLevOutp,
339         I                            qtmp1, -nRec, myIter, myThid )
340              ENDIF
341    
342    #ifdef ALLOW_MNC
343              IF (useMNC .AND. diag_mnc) THEN
344                CALL DIAGNOSTICS_MNC_OUT(
345         I                       NrMax, nLevOutp, listId, ndId,
346         I                       diag_mnc_bn,
347         I                       useMissingValue, misValLoc,
348         I                       qtmp1,
349         I                       myTime, myIter, myThid )
350            ENDIF            ENDIF
351  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
352    
353    C--      end loop on lm (or ll if ALLOW_MNC) counter
354             ENDDO
355  C--     end of Processing Fld # md  C--     end of Processing Fld # md
356          ENDIF          ENDIF
357           ENDDO
358    
359    #ifdef ALLOW_MNC
360    C--   end loop on jj counter
361        ENDDO        ENDDO
362    #endif
363    
364    #ifdef ALLOW_MDSIO
365          IF (diag_mdsio) THEN
366    C-    Note: temporary: since it is a pain to add more arguments to
367    C     all MDSIO S/R, uses instead this specific S/R to write only
368    C     meta files but with more informations in it.
369                glf = globalFiles
370                nRec = nfields(listId)*averageCycle(listId)
371                CALL MDS_WR_METAFILES(fn, prec, glf, .FALSE.,
372         &              0, 0, nLevOutp, ' ',
373         &              nfields(listId), flds(1,listId), nTimRec, timeRec,
374         &              nRec, myIter, myThid)
375          ENDIF
376    #endif /*  ALLOW_MDSIO  */
377    
378        RETURN        RETURN
379        END        END
Legend:
Removed from v.1.28  
changed lines
  Added in v.1.51
  ViewVC Help
Powered by ViewVC 1.1.22