/[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.23 by jmc, Wed Nov 2 14:42:31 2005 UTC revision 1.57 by jmc, Mon Jun 27 22:27:23 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 :: temporary array; used to store a copy of diag. output field.
62    C     qtmp2 :: temporary array; used to store a copy of a 2nd diag. field.
63    C-  Note: local common block no longer needed.
64    c     COMMON /LOCAL_DIAGNOSTICS_OUT/ qtmp1
65          _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,NrMax,nSx,nSy)
66          _RL qtmp2(1-OLx:sNx+OLx,1-OLy:sNy+OLy,NrMax,nSx,nSy)
67    
68          INTEGER i, j, k, lm
69        INTEGER bi, bj        INTEGER bi, bj
70        INTEGER md, ndId, ip, im        INTEGER md, ndId, nn, ip, im
71        INTEGER mate, mVec        INTEGER mate, mDbl, mVec
72        CHARACTER*8 parms1        CHARACTER*10 gcode
73        CHARACTER*3 mate_index        _RL undefRL
74        _RL qtmp1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+Nrphys,nSx,nSy)        INTEGER nLevOutp, kLev
       _RL undef, getcon  
       EXTERNAL getcon  
       INTEGER ILNBLNK  
       EXTERNAL ILNBLNK  
       INTEGER ilen  
       INTEGER nlevsout  
75    
76          INTEGER iLen
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          _RL     tmpLoc
84    #ifdef ALLOW_MDSIO
85        LOGICAL glf        LOGICAL glf
86    #endif
87  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
       INTEGER ii  
88        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn
       INTEGER CW_DIMS, NLEN  
       PARAMETER ( CW_DIMS = 10 )  
       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  
89  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
90    
91  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
92    
93    C---  set file properties
94        ioUnit= standardMessageUnit        ioUnit= standardMessageUnit
95        undef = getcon('UNDEF')        undefRL = UNSET_RL
96        glf = globalFiles  #ifdef ALLOW_FIZHI
97          IF ( useFIZHI ) undefRL = getcon('UNDEF')
98    #endif
99          IF ( misvalFlt(listId).NE.UNSET_RL ) undefRL = misvalFlt(listId)
100    
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              tmpLoc = FLOAT(i)
130              IF ( timeRec(1).NE.tmpLoc ) 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            tmpLoc = FLOAT(i) - 0.5 _d 0
141            IF ( timeRec(1).EQ.tmpLoc ) 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  C--   Place the loop on lm (= averagePeriod) outside the loop on md (= field):
151        IF (useMNC .AND. diag_mnc) THEN        DO lm=1,averageCycle(listId)
         DO i = 1,MAX_LEN_FNAM  
           diag_mnc_bn(i:i) = ' '  
         ENDDO  
         DO i = 1,NLEN  
           dn_blnk(i:i) = ' '  
         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,1,1,'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       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  */  
152    
153        ENDIF  #ifdef ALLOW_MNC
154           IF (useMNC .AND. diag_mnc) THEN
155             CALL DIAGNOSTICS_MNC_SET(
156         I                    nLevOutp, listId, lm,
157         O                    diag_mnc_bn,
158         I                    undefRL, myTime, myIter, myThid )
159           ENDIF
160  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
161    
162        DO md = 1,nfields(listId)  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
163    
164           DO md = 1,nfields(listId)
165          ndId = jdiag(md,listId)          ndId = jdiag(md,listId)
166          parms1 = gdiag(ndId)(1:8)          gcode = gdiag(ndId)(1:10)
167          IF ( idiag(md,listId).NE.0 .AND. parms1(5:5).NE.'D' ) THEN          mate = 0
168            mVec = 0
169            mDbl = 0
170            IF ( gcode(5:5).EQ.'C' ) THEN
171    C-      Check for Mate of a Counter Diagnostic
172               mate = hdiag(ndId)
173            ELSEIF ( gcode(5:5).EQ.'P' ) THEN
174    C-      Also load the mate (if stored) for Post-Processing
175               nn = ndId
176               DO WHILE ( gdiag(nn)(5:5).EQ.'P' )
177                 nn = hdiag(nn)
178               ENDDO
179               IF ( mdiag(md,listId).NE.0 ) mDbl = hdiag(nn)
180            ELSEIF ( gcode(1:1).EQ.'U' .OR. gcode(1:1).EQ.'V' ) THEN
181    C-      Check for Mate of a Vector Diagnostic
182               mVec = hdiag(ndId)
183            ENDIF
184            IF ( idiag(md,listId).NE.0 .AND. gcode(5:5).NE.'D' ) THEN
185  C--     Start processing 1 Fld :  C--     Start processing 1 Fld :
186    
187            ip = ABS(idiag(md,listId))            ip = ABS(idiag(md,listId)) + kdiag(ndId)*(lm-1)
188            im = mdiag(md,listId)            im = mdiag(md,listId)
189              IF (mate.GT.0) im = im + kdiag(mate)*(lm-1)
190              IF (mDbl.GT.0) im = im + kdiag(mDbl)*(lm-1)
191              IF (mVec.GT.0) im = im + kdiag(mVec)*(lm-1)
192    
193            IF ( ndiag(ip,1,1).EQ.0 ) THEN            IF ( ndiag(ip,1,1).EQ.0 ) THEN
194  C-        Empty diagnostics case :  C-        Empty diagnostics case :
195    
# Line 204  C-        Empty diagnostics case : Line 198  C-        Empty diagnostics case :
198       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter       &        '- WARNING - from DIAGNOSTICS_OUT at iter=', myIter
199              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
200       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
201              WRITE(msgBuf,'(A,I4,3A,I3,2A)')              WRITE(msgBuf,'(A,I6,3A,I4,2A)')
202       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),       &       '- WARNING -   diag.#',ndId, ' : ',flds(md,listId),
203       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)       &       ' (#',md,' ) in outp.Stream: ',fnames(listId)
204              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
205       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
206              WRITE(msgBuf,'(A,I2,A)')              IF ( averageCycle(listId).GT.1 ) THEN
207       &       '- WARNING -   has not been filled (ndiag=',               WRITE(msgBuf,'(A,2(I3,A))')
208       &       ndiag(ip,1,1), ' )'       &        '- WARNING -   has not been filled (ndiag(lm=',lm,')=',
209         &                                            ndiag(ip,1,1), ' )'
210                ELSE
211                 WRITE(msgBuf,'(A,2(I3,A))')
212         &        '- WARNING -   has not been filled (ndiag=',
213         &                                            ndiag(ip,1,1), ' )'
214                ENDIF
215              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,              CALL PRINT_MESSAGE( msgBuf, errorMessageUnit,
216       &                          SQUEEZE_RIGHT, myThid)       &                          SQUEEZE_RIGHT, myThid)
217              WRITE(msgBuf,'(A)')              WRITE(msgBuf,'(A)')
# Line 221  C-        Empty diagnostics case : Line 221  C-        Empty diagnostics case :
221              _END_MASTER( myThid )              _END_MASTER( myThid )
222              DO bj = myByLo(myThid), myByHi(myThid)              DO bj = myByLo(myThid), myByHi(myThid)
223                DO bi = myBxLo(myThid), myBxHi(myThid)                DO bi = myBxLo(myThid), myBxHi(myThid)
224                  DO k = 1,nlevels(listId)                  DO k = 1,nLevOutp
225                    DO j = 1-OLy,sNy+OLy                    DO j = 1-OLy,sNy+OLy
226                      DO i = 1-OLx,sNx+OLx                      DO i = 1-OLx,sNx+OLx
227                        qtmp1(i,j,k,bi,bj) = 0. _d 0                        qtmp1(i,j,k,bi,bj) = 0. _d 0
# Line 234  C-        Empty diagnostics case : Line 234  C-        Empty diagnostics case :
234            ELSE            ELSE
235  C-        diagnostics is not empty :  C-        diagnostics is not empty :
236    
237              IF ( myThid.EQ.1 ) WRITE(ioUnit,'(A,I3,3A,I8,2A)')              IF ( debugLevel.GE.debLevB .AND. myThid.EQ.1 ) THEN
238                  IF ( gcode(5:5).EQ.'P' ) THEN
239                   WRITE(ioUnit,'(A,I6,7A,I8,2A)')
240         &         ' Post-Processing Diag # ', ndId, '  ', cdiag(ndId),
241         &         '   Parms: ',gdiag(ndId)
242                   IF ( mDbl.EQ.0 ) THEN
243                    WRITE(ioUnit,'(2(3A,I6,A,I8))') '   from diag: ',
244         &            cdiag(nn), ' (#', nn, ') Cnt=', ndiag(ip,1,1)
245                   ELSE
246                    WRITE(ioUnit,'(2(3A,I6,A,I8))') '   from diag: ',
247         &            cdiag(nn), ' (#', nn, ') Cnt=', ndiag(ip,1,1),
248         &          ' and diag: ',
249         &            cdiag(mDbl),' (#',mDbl,') Cnt=',ndiag(im,1,1)
250                   ENDIF
251                  ELSE
252                   WRITE(ioUnit,'(A,I6,3A,I8,2A)')
253       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),
254       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)
255                  ENDIF
256              IF ( parms1(5:5).EQ.'C' ) THEN                IF ( mate.GT.0 ) THEN
257  C             Check for Mate of a Counter Diagnostic                 WRITE(ioUnit,'(3A,I6,2A)')
 C             --------------------------------------  
               mate_index = parms1(6:8)  
               READ (mate_index,'(I3)') mate  
               IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,2A)')  
258       &         '       use Counter Mate for  ', cdiag(ndId),       &         '       use Counter Mate for  ', cdiag(ndId),
259       &         '     Diagnostic # ',mate, '  ', cdiag(mate)       &         '     Diagnostic # ',mate, '  ', cdiag(mate)
260                  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  
261                  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
262                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
263       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
264       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
265       &             ' exists '       &             ' exists '
266                  ELSE                  ELSE
267                   IF ( myThid.EQ.1 ) WRITE(ioUnit,'(3A,I3,3A)')                   WRITE(ioUnit,'(3A,I6,3A)')
268       &             '           Vector  Mate for  ', cdiag(ndId),       &             '           Vector  Mate for  ', cdiag(ndId),
269       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),       &             '     Diagnostic # ',mVec, '  ', cdiag(mVec),
270       &             ' not enabled'       &             ' not enabled'
# Line 269  C             -------------------------- Line 272  C             --------------------------
272                ENDIF                ENDIF
273              ENDIF              ENDIF
274    
275              DO bj = myByLo(myThid), myByHi(myThid)              IF ( fflags(listId)(2:2).EQ.' ' ) THEN
276               DO bi = myBxLo(myThid), myBxHi(myThid)  C-       get only selected levels:
277                DO k = 1,nlevels(listId)                DO bj = myByLo(myThid), myByHi(myThid)
278                  CALL GETDIAG(                 DO bi = myBxLo(myThid), myBxHi(myThid)
279       I                       levs(k,listId),undef,                  DO k = 1,nlevels(listId)
280       O                       qtmp1(1-OLx,1-OLy,k,bi,bj),                    kLev = NINT(levs(k,listId))
281       I                       ndId,mate,ip,im,bi,bj,myThid)                    CALL DIAGNOSTICS_GET_DIAG(
282         I                         kLev, undefRL,
283         O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
284         I                         ndId, mate, ip, im, bi, bj, myThid )
285                    ENDDO
286                   ENDDO
287                ENDDO                ENDDO
288               ENDDO                IF ( mDbl.GT.0 ) THEN
289              ENDDO                 DO bj = myByLo(myThid), myByHi(myThid)
290                    DO bi = myBxLo(myThid), myBxHi(myThid)
291  C-        end of empty diag / not empty block                   DO k = 1,nlevels(listId)
292            ENDIF                    kLev = NINT(levs(k,listId))
293                      CALL DIAGNOSTICS_GET_DIAG(
294            nlevsout = nlevels(listId)       I                         kLev, undefRL,
295         O                         qtmp2(1-OLx,1-OLy,k,bi,bj),
296  C-----------------------------------------------------------------------       I                         mDbl, 0, im, 0, bi, bj, myThid )
297  C         Check to see if we need to interpolate before output                   ENDDO
298  C-----------------------------------------------------------------------                  ENDDO
299           IF ( fflags(listId)(2:2).EQ.'P' ) THEN                 ENDDO
300  C-        Do vertical interpolation:                ENDIF
           CALL DIAGNOSTICS_INTERP_VERT(  
      I                     listId, md, ndId, ip, im,  
      U                     nlevsout,  
      U                     qtmp1,  
      I                     undef,  
      I                     myTime, myIter, myThid )  
          ENDIF  
   
 #ifdef ALLOW_MDSIO  
 C         Prepare for mdsio optionality  
           IF (diag_mdsio) THEN  
             IF (fflags(listId)(1:1) .EQ. 'R') THEN  
 C             Force it to be 32-bit precision  
               CALL MDSWRITEFIELD_NEW(fn,precFloat32,glf,.FALSE.,  
      &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)  
             ELSEIF (fflags(listId)(1:1) .EQ. 'D') THEN  
 C             Force it to be 64-bit precision  
               CALL MDSWRITEFIELD_NEW(fn,precFloat64,glf,.FALSE.,  
      &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)  
301              ELSE              ELSE
302  C             This is the old default behavior  C-       get all the levels (for vertical post-processing)
303                CALL MDSWRITEFIELD_NEW(fn,writeBinaryPrec,glf,.FALSE.,                DO bj = myByLo(myThid), myByHi(myThid)
304       &             'RL',Nr+Nrphys,nlevsout,qtmp1,md,myIter,myThid)                 DO bi = myBxLo(myThid), myBxHi(myThid)
305                      CALL DIAGNOSTICS_GET_DIAG(
306         I                         0, undefRL,
307         O                         qtmp1(1-OLx,1-OLy,1,bi,bj),
308         I                         ndId, mate, ip, im, bi, bj, myThid )
309                   ENDDO
310                  ENDDO
311                  IF ( mDbl.GT.0 ) THEN
312                   DO bj = myByLo(myThid), myByHi(myThid)
313                    DO bi = myBxLo(myThid), myBxHi(myThid)
314                     DO k = 1,nlevels(listId)
315                      CALL DIAGNOSTICS_GET_DIAG(
316         I                         0, undefRL,
317         O                         qtmp2(1-OLx,1-OLy,k,bi,bj),
318         I                         mDbl, 0, im, 0, bi, bj, myThid )
319                     ENDDO
320                    ENDDO
321                   ENDDO
322                  ENDIF
323              ENDIF              ENDIF
           ENDIF  
 #endif /*  ALLOW_MDSIO  */  
   
 #ifdef ALLOW_MNC  
           IF (useMNC .AND. diag_mnc) THEN  
   
             _BEGIN_MASTER( myThid )  
324    
325              DO ii = 1,CW_DIMS  C-----------------------------------------------------------------------
326                d_cw_name(1:NLEN) = dn_blnk(1:NLEN)  C--     Apply specific post-processing (e.g., interpolate) before output
327                dn(ii)(1:NLEN) = dn_blnk(1:NLEN)  C-----------------------------------------------------------------------
328              ENDDO              IF ( fflags(listId)(2:2).EQ.'P' ) THEN
329    C-          Do vertical interpolation:
330  C           Note that the "d_cw_name" variable is a hack that hides a               IF ( fluidIsAir ) THEN
331  C           subtlety within MNC.  Basically, each MNC-wrapped file is  C jmc: for now, this can only work in an atmospheric set-up (fluidIsAir);
332  C           caching its own concept of what each "grid name" (that is, a                CALL DIAGNOSTICS_INTERP_VERT(
333  C           dimension group name) means.  So one cannot re-use the same       I                         listId, md, ndId, ip, im, lm,
334  C           "grid" name for different collections of dimensions within a       U                         qtmp1, qtmp2,
335  C           given file.  By appending the "ndId" values to each name, we       I                         undefRL, myTime, myIter, myThid )
336  C           guarantee uniqueness within each MNC-produced file.               ELSE
337              WRITE(d_cw_name,'(a,i6.6)') 'd_cw_',ndId                 WRITE(msgBuf,'(2A)') 'DIAGNOSTICS_OUT: ',
338         &           'INTERP_VERT not allowed in this config'
339  C           XY dimensions                 CALL PRINT_ERROR( msgBuf , myThid )
340              dim(1)       = sNx + 2*OLx                 STOP 'ABNORMAL END: S/R DIAGNOSTICS_OUT'
341              dim(2)       = sNy + 2*OLy               ENDIF
             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', nlevels(listId)  
             IF ( (gdiag(ndId)(10:10) .EQ. 'R')  
      &           .AND. (gdiag(ndId)(9:9) .EQ. 'M') ) THEN  
               WRITE(dn(3),'(a,i6.6)') 'Zmd', nlevels(listId)  
             ENDIF  
             IF ( (gdiag(ndId)(10:10) .EQ. 'R')  
      &           .AND. (gdiag(ndId)(9:9) .EQ. 'L') ) THEN  
               WRITE(dn(3),'(a,i6.6)') 'Zld', nlevels(listId)  
342              ENDIF              ENDIF
343              IF ( (gdiag(ndId)(10:10) .EQ. 'R')              IF ( fflags(listId)(2:2).EQ.'I' ) THEN
344       &           .AND. (gdiag(ndId)(9:9) .EQ. 'U') ) THEN  C-          Integrate vertically: for now, output field has just 1 level:
345                WRITE(dn(3),'(a,i6.6)') 'Zud', nlevels(listId)                CALL DIAGNOSTICS_SUM_LEVELS(
346         I                         listId, md, ndId, ip, im, lm,
347         U                         qtmp1,
348         I                         undefRL, myTime, myIter, myThid )
349              ENDIF              ENDIF
350              dim(3) = Nr+Nrphys              IF ( gcode(5:5).EQ.'P' ) THEN
351              ib(3)  = 1  C-          Do Post-Processing:
352              ie(3)  = nlevels(listId)               IF ( flds(md,listId).EQ.'PhiVEL  '
353    c    &       .OR. flds(md,listId).EQ.'PsiVEL  '
354  C           Time dimension       &          ) THEN
355              dn(4)(1:1) = 'T'                CALL DIAGNOSTICS_CALC_PHIVEL(
356              dim(4) = -1       I                         listId, md, ndId, ip, im, lm,
357              ib(4)  = 1       U                         qtmp1, qtmp2,
358              ie(4)  = 1       I                         myTime, myIter, myThid )
359                 ELSE
360              CALL MNC_CW_ADD_GNAME(d_cw_name, 4,                 WRITE(msgBuf,'(2A)') 'DIAGNOSTICS_OUT: ',
361       &             dim, dn, ib, ie, myThid)       &           'unknown Processing method for diag="',cdiag(ndId),'"'
362              CALL MNC_CW_ADD_VNAME(cdiag(ndId), d_cw_name,                 CALL PRINT_ERROR( msgBuf , myThid )
363       &             4,5, myThid)                 STOP 'ABNORMAL END: S/R DIAGNOSTICS_OUT'
364              CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'description',               ENDIF
      &             tdiag(ndId),myThid)  
             CALL MNC_CW_ADD_VATTR_TEXT(cdiag(ndId),'units',  
      &             udiag(ndId),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)  
365              ENDIF              ENDIF
               
             CALL MNC_CW_DEL_VNAME(cdiag(ndId), myThid)  
             CALL MNC_CW_DEL_GNAME(d_cw_name, myThid)  
366    
367              _END_MASTER( myThid )  C--     End of empty diag / not-empty diag block
368              ENDIF
369    
370    C--     Ready to write field "md", element "lm" in averageCycle(listId)
371    
372    C-        write to binary file, using MDSIO pkg:
373              IF ( diag_mdsio ) THEN
374    c           nRec = lm + (md-1)*averageCycle(listId)
375                nRec = md + (lm-1)*nfields(listId)
376    C           default precision for output files
377                prec = writeBinaryPrec
378    C           fFlag(1)=R(or D): force it to be 32-bit(or 64) precision
379                IF ( fflags(listId)(1:1).EQ.'R' ) prec = precFloat32
380                IF ( fflags(listId)(1:1).EQ.'D' ) prec = precFloat64
381    C         a hack not to write meta files now: pass -nRec < 0 to MDS_WRITE S/R
382                CALL WRITE_REC_LEV_RL(
383         I                            fn, prec,
384         I                            NrMax, 1, nLevOutp,
385         I                            qtmp1, -nRec, myIter, myThid )
386              ENDIF
387    
388    #ifdef ALLOW_MNC
389              IF (useMNC .AND. diag_mnc) THEN
390                CALL DIAGNOSTICS_MNC_OUT(
391         I                       NrMax, nLevOutp, listId, ndId, mate,
392         I                       diag_mnc_bn, qtmp1,
393         I                       undefRL, myTime, myIter, myThid )
394            ENDIF            ENDIF
395  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
396    
397  C--     end of Processing Fld # md  C--     end of Processing Fld # md
398          ENDIF          ENDIF
399           ENDDO
400    
401    C--   end loop on lm counter (= averagePeriod)
402        ENDDO        ENDDO
403    
404    #ifdef ALLOW_MDSIO
405          IF (diag_mdsio) THEN
406    C-    Note: temporary: since it is a pain to add more arguments to
407    C     all MDSIO S/R, uses instead this specific S/R to write only
408    C     meta files but with more informations in it.
409                glf = globalFiles
410                nRec = averageCycle(listId)*nfields(listId)
411                CALL MDS_WR_METAFILES(fn, prec, glf, .FALSE.,
412         &              0, 0, nLevOutp, ' ',
413         &              nfields(listId), flds(1,listId), nTimRec, timeRec,
414         &              nRec, myIter, myThid)
415          ENDIF
416    #endif /*  ALLOW_MDSIO  */
417    
418        RETURN        RETURN
419        END        END
420    
Legend:
Removed from v.1.23  
changed lines
  Added in v.1.57
  ViewVC Help
Powered by ViewVC 1.1.22