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