/[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.52 by jmc, Sun Jun 12 19:16:09 2011 UTC revision 1.61 by jmc, Wed Feb 6 21:25:26 2013 UTC
# Line 41  CEOP Line 41  CEOP
41  C     !FUNCTIONS:  C     !FUNCTIONS:
42        INTEGER ILNBLNK        INTEGER ILNBLNK
43        EXTERNAL ILNBLNK        EXTERNAL ILNBLNK
 #ifdef ALLOW_FIZHI  
       _RL   getcon  
       EXTERNAL getcon  
 #endif  
44    
45  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
46  C     i,j,k :: loop indices  C     i,j,k :: loop indices
# Line 52  C     bi,bj :: tile indices Line 48  C     bi,bj :: tile indices
48  C     lm    :: loop index (averageCycle)  C     lm    :: loop index (averageCycle)
49  C     md    :: field number in the list "listId".  C     md    :: field number in the list "listId".
50  C     ndId  :: diagnostics  Id number (in available diagnostics list)  C     ndId  :: diagnostics  Id number (in available diagnostics list)
 C     mate  :: counter mate Id number (in available diagnostics list)  
51  C     ip    :: diagnostics  pointer to storage array  C     ip    :: diagnostics  pointer to storage array
52  C     im    :: counter-mate pointer to storage array  C     im    :: counter-mate pointer to storage array
53    C     mate  :: counter mate Id number (in available diagnostics list)
54    C     mDbl  :: processing mate Id number (in case processing requires 2 diags)
55    C     mVec  :: vector mate Id number
56    C     ppFld :: post-processed diag or not (=0): =1 stored in qtmp1 ; =2 in qtmp2
57    C   isComputed :: previous post-processed diag (still available in qtmp)
58  C     nLevOutp :: number of levels to write in output file  C     nLevOutp :: number of levels to write in output file
59  C  C
60  C--   COMMON /LOCAL_DIAGNOSTICS_OUT/ local common block (for multi-threaded)  C--   COMMON /LOCAL_DIAGNOSTICS_OUT/ local common block (for multi-threaded)
# Line 67  c     COMMON /LOCAL_DIAGNOSTICS_OUT/ qtm Line 67  c     COMMON /LOCAL_DIAGNOSTICS_OUT/ qtm
67    
68        INTEGER i, j, k, lm        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          INTEGER ppFld, isComputed
73        CHARACTER*10 gcode        CHARACTER*10 gcode
74        _RL undefRL        _RL undefRL
75        INTEGER nLevOutp, kLev        INTEGER nLevOutp, kLev
# Line 85  c     COMMON /LOCAL_DIAGNOSTICS_OUT/ qtm Line 86  c     COMMON /LOCAL_DIAGNOSTICS_OUT/ qtm
86        LOGICAL glf        LOGICAL glf
87  #endif  #endif
88  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
       INTEGER ll, llMx, jj, jjMx  
89        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn        CHARACTER*(MAX_LEN_FNAM) diag_mnc_bn
       LOGICAL useMissingValue  
       REAL*8 misValLoc  
90  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
91    
92  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
93    
94  C---  set file properties  C---  set file properties
95        ioUnit= standardMessageUnit        ioUnit= standardMessageUnit
96        undefRL = UNSET_RL        undefRL = misValFlt(listId)
97  #ifdef ALLOW_FIZHI  
       IF ( useFIZHI ) undefRL = getcon('UNDEF')  
 #endif  
98        WRITE(suff,'(I10.10)') myIter        WRITE(suff,'(I10.10)') myIter
99        iLen = ILNBLNK(fnames(listId))        iLen = ILNBLNK(fnames(listId))
100        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:iLen),'.',suff(1:10)        WRITE( fn, '(A,A,A)' ) fnames(listId)(1:iLen),'.',suff(1:10)
# Line 148  c     DO i=1,nTimRec Line 144  c     DO i=1,nTimRec
144  c       timeRec(i) = timeRec(i)/deltaTClock  c       timeRec(i) = timeRec(i)/deltaTClock
145  c     ENDDO  c     ENDDO
146    
147  #ifdef ALLOW_MNC  C--   Place the loop on lm (= averagePeriod) outside the loop on md (= field):
148  C-- this is a trick to reverse the order of the loops on md (= field)        DO lm=1,averageCycle(listId)
 C   and lm (= averagePeriod): binary output: lm loop inside md loop ;  
 C                                 mnc ouput: md loop inside lm loop.  
       IF (useMNC .AND. diag_mnc) THEN  
         jjMx = averageCycle(listId)  
         llMx = 1  
       ELSE  
         jjMx = 1  
         llMx = averageCycle(listId)  
       ENDIF  
       DO jj=1,jjMx  
149    
150    #ifdef ALLOW_MNC
151         IF (useMNC .AND. diag_mnc) THEN         IF (useMNC .AND. diag_mnc) THEN
          misValLoc = undefRL  
          IF ( misvalFlt(listId).NE.UNSET_RL )  
      &        misValLoc = misvalFlt(listId)  
152           CALL DIAGNOSTICS_MNC_SET(           CALL DIAGNOSTICS_MNC_SET(
153       I                    nLevOutp, listId, jj,       I                    nLevOutp, listId, lm,
154       O                    diag_mnc_bn, useMissingValue,       O                    diag_mnc_bn,
155       I                    misValLoc, myTime, myIter, myThid )       I                    undefRL, myTime, myIter, myThid )
156         ENDIF         ENDIF
157  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
158    
159  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|  C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
160    
161           isComputed = 0
162         DO md = 1,nfields(listId)         DO md = 1,nfields(listId)
163          ndId = jdiag(md,listId)          ndId = jdiag(md,listId)
164          gcode = gdiag(ndId)(1:10)          gcode = gdiag(ndId)(1:10)
165          mate = 0          mate = 0
166          mVec = 0          mVec = 0
167            mDbl = 0
168            ppFld = 0
169          IF ( gcode(5:5).EQ.'C' ) THEN          IF ( gcode(5:5).EQ.'C' ) THEN
170  C-      Check for Mate of a Counter Diagnostic  C-      Check for Mate of a Counter Diagnostic
171             mate = hdiag(ndId)             mate = hdiag(ndId)
172            ELSEIF ( gcode(5:5).EQ.'P' ) THEN
173               ppFld = 1
174               IF ( gdiag(hdiag(ndId))(5:5).EQ.'P' ) ppFld = 2
175    C-      Also load the mate (if stored) for Post-Processing
176               nn = ndId
177               DO WHILE ( gdiag(nn)(5:5).EQ.'P' )
178                 nn = hdiag(nn)
179               ENDDO
180               IF ( mdiag(md,listId).NE.0 ) mDbl = hdiag(nn)
181    c          write(0,*) ppFld,' ndId=', ndId, nn, mDbl, isComputed
182          ELSEIF ( gcode(1:1).EQ.'U' .OR. gcode(1:1).EQ.'V' ) THEN          ELSEIF ( gcode(1:1).EQ.'U' .OR. gcode(1:1).EQ.'V' ) THEN
183  C-      Check for Mate of a Vector Diagnostic  C-      Check for Mate of a Vector Diagnostic
184             mVec = hdiag(ndId)             mVec = hdiag(ndId)
185          ENDIF          ENDIF
186          IF ( idiag(md,listId).NE.0 .AND. gcode(5:5).NE.'D' ) THEN          IF ( idiag(md,listId).NE.0 .AND. gcode(5:5).NE.'D' ) THEN
187  C--     Start processing 1 Fld :  C--     Start processing 1 Fld :
 #ifdef ALLOW_MNC  
          DO ll=1,llMx  
           lm = jj+ll-1  
 #else  
          DO lm=1,averageCycle(listId)  
 #endif  
188    
189            ip = ABS(idiag(md,listId)) + kdiag(ndId)*(lm-1)            ip = ABS(idiag(md,listId)) + kdiag(ndId)*(lm-1)
190            im = mdiag(md,listId)            im = mdiag(md,listId)
191            IF (mate.GT.0) im = im + kdiag(mate)*(lm-1)            IF (mate.GT.0) im = im + kdiag(mate)*(lm-1)
192              IF (mDbl.GT.0) im = im + kdiag(mDbl)*(lm-1)
193            IF (mVec.GT.0) im = im + kdiag(mVec)*(lm-1)            IF (mVec.GT.0) im = im + kdiag(mVec)*(lm-1)
194    
195            IF ( ndiag(ip,1,1).EQ.0 ) THEN            IF ( ppFld.EQ.2 .AND. isComputed.EQ.hdiag(ndId) ) THEN
196    C-        Post-Processed diag from an other Post-Processed diag -and-
197    C         both of them have just been calculated and are still stored in qtmp:
198    C         => skip computation and just write qtmp2
199                IF ( debugLevel.GE.debLevB .AND. myThid.EQ.1 ) THEN
200                   WRITE(ioUnit,'(A,I6,3A,I6)')
201         &         '  get Post-Proc. Diag # ', ndId, '  ', cdiag(ndId),
202         &         ' from previous computation of Diag # ', isComputed
203                ENDIF
204                isComputed = 0
205              ELSEIF ( ndiag(ip,1,1).EQ.0 ) THEN
206  C-        Empty diagnostics case :  C-        Empty diagnostics case :
207                isComputed = 0
208    
209              _BEGIN_MASTER( myThid )              _BEGIN_MASTER( myThid )
210              WRITE(msgBuf,'(A,I10)')              WRITE(msgBuf,'(A,I10)')
# Line 243  C-        Empty diagnostics case : Line 246  C-        Empty diagnostics case :
246    
247            ELSE            ELSE
248  C-        diagnostics is not empty :  C-        diagnostics is not empty :
249                isComputed = 0
250    
251              IF ( debugLevel.GE.debLevB .AND. myThid.EQ.1 ) THEN              IF ( debugLevel.GE.debLevB .AND. myThid.EQ.1 ) THEN
252                WRITE(ioUnit,'(A,I6,3A,I8,2A)')                IF ( ppFld.GE.1 ) THEN
253                   WRITE(ioUnit,'(A,I6,7A,I8,2A)')
254         &         ' Post-Processing Diag # ', ndId, '  ', cdiag(ndId),
255         &         '   Parms: ',gdiag(ndId)
256                   IF ( mDbl.EQ.0 ) THEN
257                    WRITE(ioUnit,'(2(3A,I6,A,I8))') '   from diag: ',
258         &            cdiag(nn), ' (#', nn, ') Cnt=', ndiag(ip,1,1)
259                   ELSE
260                    WRITE(ioUnit,'(2(3A,I6,A,I8))') '   from diag: ',
261         &            cdiag(nn), ' (#', nn, ') Cnt=', ndiag(ip,1,1),
262         &          ' and diag: ',
263         &            cdiag(mDbl),' (#',mDbl,') Cnt=',ndiag(im,1,1)
264                   ENDIF
265                  ELSE
266                   WRITE(ioUnit,'(A,I6,3A,I8,2A)')
267       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),       &         ' Computing Diagnostic # ', ndId, '  ', cdiag(ndId),
268       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)       &         '     Counter:',ndiag(ip,1,1),'   Parms: ',gdiag(ndId)
269                  ENDIF
270                IF ( mate.GT.0 ) THEN                IF ( mate.GT.0 ) THEN
271                 WRITE(ioUnit,'(3A,I6,2A)')                 WRITE(ioUnit,'(3A,I6,2A)')
272       &         '       use Counter Mate for  ', cdiag(ndId),       &         '       use Counter Mate for  ', cdiag(ndId),
# Line 276  C-       get only selected levels: Line 295  C-       get only selected levels:
295                    CALL DIAGNOSTICS_GET_DIAG(                    CALL DIAGNOSTICS_GET_DIAG(
296       I                         kLev, undefRL,       I                         kLev, undefRL,
297       O                         qtmp1(1-OLx,1-OLy,k,bi,bj),       O                         qtmp1(1-OLx,1-OLy,k,bi,bj),
298       I                         ndId,mate,ip,im,bi,bj,myThid)       I                         ndId, mate, ip, im, bi, bj, myThid )
299                  ENDDO                  ENDDO
300                 ENDDO                 ENDDO
301                ENDDO                ENDDO
302                  IF ( mDbl.GT.0 ) THEN
303                   DO bj = myByLo(myThid), myByHi(myThid)
304                    DO bi = myBxLo(myThid), myBxHi(myThid)
305                     DO k = 1,nlevels(listId)
306                      kLev = NINT(levs(k,listId))
307                      CALL DIAGNOSTICS_GET_DIAG(
308         I                         kLev, undefRL,
309         O                         qtmp2(1-OLx,1-OLy,k,bi,bj),
310         I                         mDbl, 0, im, 0, bi, bj, myThid )
311                     ENDDO
312                    ENDDO
313                   ENDDO
314                  ENDIF
315              ELSE              ELSE
316  C-       get all the levels (for vertical post-processing)  C-       get all the levels (for vertical post-processing)
317                DO bj = myByLo(myThid), myByHi(myThid)                DO bj = myByLo(myThid), myByHi(myThid)
# Line 287  C-       get all the levels (for vertica Line 319  C-       get all the levels (for vertica
319                    CALL DIAGNOSTICS_GET_DIAG(                    CALL DIAGNOSTICS_GET_DIAG(
320       I                         0, undefRL,       I                         0, undefRL,
321       O                         qtmp1(1-OLx,1-OLy,1,bi,bj),       O                         qtmp1(1-OLx,1-OLy,1,bi,bj),
322       I                         ndId,mate,ip,im,bi,bj,myThid)       I                         ndId, mate, ip, im, bi, bj, myThid )
323                 ENDDO                 ENDDO
324                ENDDO                ENDDO
325                  IF ( mDbl.GT.0 ) THEN
326                   DO bj = myByLo(myThid), myByHi(myThid)
327                    DO bi = myBxLo(myThid), myBxHi(myThid)
328                      CALL DIAGNOSTICS_GET_DIAG(
329         I                         0, undefRL,
330         O                         qtmp2(1-OLx,1-OLy,1,bi,bj),
331         I                         mDbl, 0, im, 0, bi, bj, myThid )
332                    ENDDO
333                   ENDDO
334                  ENDIF
335              ENDIF              ENDIF
336    
337  C-----------------------------------------------------------------------  C-----------------------------------------------------------------------
# Line 317  C-          Integrate vertically: for no Line 359  C-          Integrate vertically: for no
359       U                         qtmp1,       U                         qtmp1,
360       I                         undefRL, myTime, myIter, myThid )       I                         undefRL, myTime, myIter, myThid )
361              ENDIF              ENDIF
362                IF ( ppFld.GE.1 ) THEN
363    C-          Do Post-Processing:
364                 IF ( flds(md,listId).EQ.'PhiVEL  '
365         &       .OR. flds(md,listId).EQ.'PsiVEL  '
366         &          ) THEN
367                  CALL DIAGNOSTICS_CALC_PHIVEL(
368         I                         listId, md, ndId, ip, im, lm,
369         I                         NrMax,
370         U                         qtmp1, qtmp2,
371         I                         myTime, myIter, myThid )
372                  isComputed = ndId
373                 ELSE
374                   WRITE(msgBuf,'(2A)') 'DIAGNOSTICS_OUT: ',
375         &           'unknown Processing method for diag="',cdiag(ndId),'"'
376                   CALL PRINT_ERROR( msgBuf , myThid )
377                   STOP 'ABNORMAL END: S/R DIAGNOSTICS_OUT'
378                 ENDIF
379                ENDIF
380    
381  C--     End of empty diag / not-empty diag block  C--     End of empty diag / not-empty diag block
382            ENDIF            ENDIF
# Line 325  C--     Ready to write field "md", eleme Line 385  C--     Ready to write field "md", eleme
385    
386  C-        write to binary file, using MDSIO pkg:  C-        write to binary file, using MDSIO pkg:
387            IF ( diag_mdsio ) THEN            IF ( diag_mdsio ) THEN
388              nRec = lm + (md-1)*averageCycle(listId)  c          nRec = lm + (md-1)*averageCycle(listId)
389  C           default precision for output files             nRec = md + (lm-1)*nfields(listId)
390              prec = writeBinaryPrec  C         default precision for output files
391  C           fFlag(1)=R(or D): force it to be 32-bit(or 64) precision             prec = writeBinaryPrec
392              IF ( fflags(listId)(1:1).EQ.'R' ) prec = precFloat32  C         fFlag(1)=R(or D): force it to be 32-bit(or 64) precision
393              IF ( fflags(listId)(1:1).EQ.'D' ) prec = precFloat64             IF ( fflags(listId)(1:1).EQ.'R' ) prec = precFloat32
394               IF ( fflags(listId)(1:1).EQ.'D' ) prec = precFloat64
395  C         a hack not to write meta files now: pass -nRec < 0 to MDS_WRITE S/R  C         a hack not to write meta files now: pass -nRec < 0 to MDS_WRITE S/R
396               IF ( ppFld.LE.1 ) THEN
397              CALL WRITE_REC_LEV_RL(              CALL WRITE_REC_LEV_RL(
398       I                            fn, prec,       I                            fn, prec,
399       I                            NrMax, 1, nLevOutp,       I                            NrMax, 1, nLevOutp,
400       I                            qtmp1, -nRec, myIter, myThid )       I                            qtmp1, -nRec, myIter, myThid )
401               ELSE
402                CALL WRITE_REC_LEV_RL(
403         I                            fn, prec,
404         I                            NrMax, 1, nLevOutp,
405         I                            qtmp2, -nRec, myIter, myThid )
406               ENDIF
407            ENDIF            ENDIF
408    
409  #ifdef ALLOW_MNC  #ifdef ALLOW_MNC
410            IF (useMNC .AND. diag_mnc) THEN            IF (useMNC .AND. diag_mnc) THEN
411               IF ( ppFld.LE.1 ) THEN
412                CALL DIAGNOSTICS_MNC_OUT(
413         I                       NrMax, nLevOutp, listId, ndId, mate,
414         I                       diag_mnc_bn, qtmp1,
415         I                       undefRL, myTime, myIter, myThid )
416               ELSE
417              CALL DIAGNOSTICS_MNC_OUT(              CALL DIAGNOSTICS_MNC_OUT(
418       I                       NrMax, nLevOutp, listId, ndId,       I                       NrMax, nLevOutp, listId, ndId, mate,
419       I                       diag_mnc_bn,       I                       diag_mnc_bn, qtmp2,
420       I                       useMissingValue, misValLoc,       I                       undefRL, myTime, myIter, myThid )
421       I                       qtmp1,             ENDIF
      I                       myTime, myIter, myThid )  
422            ENDIF            ENDIF
423  #endif /*  ALLOW_MNC  */  #endif /*  ALLOW_MNC  */
424    
 C--      end loop on lm (or ll if ALLOW_MNC) counter  
          ENDDO  
425  C--     end of Processing Fld # md  C--     end of Processing Fld # md
426          ENDIF          ENDIF
427         ENDDO         ENDDO
428    
429  #ifdef ALLOW_MNC  C--   end loop on lm counter (= averagePeriod)
 C--   end loop on jj counter  
430        ENDDO        ENDDO
 #endif  
431    
432  #ifdef ALLOW_MDSIO  #ifdef ALLOW_MDSIO
433        IF (diag_mdsio) THEN        IF (diag_mdsio) THEN
# Line 366  C-    Note: temporary: since it is a pai Line 435  C-    Note: temporary: since it is a pai
435  C     all MDSIO S/R, uses instead this specific S/R to write only  C     all MDSIO S/R, uses instead this specific S/R to write only
436  C     meta files but with more informations in it.  C     meta files but with more informations in it.
437              glf = globalFiles              glf = globalFiles
438              nRec = nfields(listId)*averageCycle(listId)              nRec = averageCycle(listId)*nfields(listId)
439              CALL MDS_WR_METAFILES(fn, prec, glf, .FALSE.,              CALL MDS_WR_METAFILES(fn, prec, glf, .FALSE.,
440       &              0, 0, nLevOutp, ' ',       &              0, 0, nLevOutp, ' ',
441       &              nfields(listId), flds(1,listId), nTimRec, timeRec,       &              nfields(listId), flds(1,listId),
442         &              nTimRec, timeRec, undefRL,
443       &              nRec, myIter, myThid)       &              nRec, myIter, myThid)
444        ENDIF        ENDIF
445  #endif /*  ALLOW_MDSIO  */  #endif /*  ALLOW_MDSIO  */
Legend:
Removed from v.1.52  
changed lines
  Added in v.1.61
  ViewVC Help
Powered by ViewVC 1.1.22