C $Header: /u/gcmpack/MITgcm/pkg/thsice/thsice_output.F,v 1.4 2005/06/24 04:36:54 edhill Exp $
C $Name:  $

#include "THSICE_OPTIONS.h"

CBOP
C     !ROUTINE: THSICE_OUTPUT
C     !INTERFACE:
      SUBROUTINE THSICE_OUTPUT( myTime, myIter, myThid )

C     !DESCRIPTION: \bv
C     *==========================================================*
C     | S/R THSICE_OUTPUT
C     | o general routine for ThSIce output
C     *==========================================================*
C     | - write time-average & snap-shot output
C     | - call monitor to write global quantities
C     *==========================================================*
C     \ev

C     !USES:
      IMPLICIT NONE

C     === Global variables ===
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "THSICE_PARAMS.h"
#include "THSICE_VARS.h"
#include "THSICE_TAVE.h"

C     !INPUT/OUTPUT PARAMETERS:
C     == Routine arguments ==
C     myTime :: Current time of simulation ( s )
C     myIter :: Iteration number
C     myThid :: my Thread Id number
      _RL     myTime
      INTEGER myIter
      INTEGER myThid
CEOP

#ifdef ALLOW_THSICE

C     == Local variables ==
      LOGICAL  DIFFERENT_MULTIPLE
      EXTERNAL 
      INTEGER bi, bj, kl
      CHARACTER*(MAX_LEN_MBUF) suff, fn
      LOGICAL gf
#ifdef ALLOW_TIMEAVE
      INTEGER i
      CHARACTER*(MAX_LEN_MBUF) mncf
#endif

#ifdef ALLOW_TIMEAVE

c set arrays to zero if first timestep
      IF ( myIter .EQ. nIter0 ) THEN
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)
         CALL TIMEAVE_RESET(ice_fract_Ave,  1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_iceH_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_snowH_Ave,  1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_Tsrf_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_Tice1_Ave,  1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_Tice2_Ave,  1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_snowPr_Ave, 1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_flx2oc_Ave, 1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_frw2oc_Ave, 1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_salFx_Ave,  1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_flxAtm_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_frwAtm_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_albedo_Ave, 1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_tMxL_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_sMxL_Ave,   1, bi, bj, myThid)
         ice_timeAve(1,bi,bj)=0.
        ENDDO
       ENDDO

C     Dump files and restart average computation if needed
      ELSEIF ( 
     &  DIFFERENT_MULTIPLE(thSIce_taveFreq,myTime,deltaTClock)
     &     ) THEN

C--    Normalize by integrated time
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)
C-     area weighted average (with ice-fraction)
         CALL TIMEAVE_NORMA_2V( ice_iceH_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_snowH_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_Tsrf_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_Tice1_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_Tice2_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_snowPr_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )
         CALL TIMEAVE_NORMA_2V( ice_albedo_Ave,
     &                          ice_fract_Ave, 1, bi, bj, myThid )

C-     simple time average :
         CALL TIMEAVE_NORMALIZ(ice_fract_Ave,ice_timeAve,1 ,
     &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_iceH_Ave,   ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_snowH_Ave,  ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_Tsrf_Ave,   ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_Tice1_Ave,  ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_Tice2_Ave,  ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
c        CALL TIMEAVE_NORMALIZ(ice_snowPr_Ave, ice_timeAve, 1 ,
c    &                 bi,bj,myThid)
         CALL TIMEAVE_NORMALIZ(ice_flx2oc_Ave, ice_timeAve, 1 ,
     &                 bi,bj,myThid)
         CALL TIMEAVE_NORMALIZ(ice_frw2oc_Ave, ice_timeAve, 1 ,
     &                 bi,bj,myThid)
         CALL TIMEAVE_NORMALIZ(ice_salFx_Ave,  ice_timeAve, 1 ,
     &                 bi,bj,myThid)
         IF ( fluidIsWater ) THEN
          CALL TIMEAVE_NORMALIZ(ice_flxAtm_Ave,ice_timeAve, 1 ,
     &                 bi,bj,myThid)
          CALL TIMEAVE_NORMALIZ(ice_frwAtm_Ave,ice_timeAve, 1 ,
     &                 bi,bj,myThid)
         ENDIF
         IF ( stepFwd_oceMxL ) THEN
          CALL TIMEAVE_NORMALIZ(ice_tMxL_Ave,  ice_timeAve, 1 ,
     &                 bi,bj,myThid)
          CALL TIMEAVE_NORMALIZ(ice_sMxL_Ave,  ice_timeAve, 1 ,
     &                 bi,bj,myThid)
         ENDIF
        ENDDO
       ENDDO
         
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|

C--    Write 1 file per time-average variable:
c      WRITE(suff,'(I10.10)') myIter
c      CALL WRITE_FLD_XY_RL('ICE_fract-T.',  suff,
c    &                       ice_fract_Ave,  myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_iceH-T.',   suff,
c    &                       ice_iceH_Ave,   myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_snowH-T.',  suff,
c    &                       ice_snowH_Ave,  myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_Tsrf-T.',   suff,
c    &                       ice_Tsrf_Ave,   myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_Tice1-T.',  suff,
c    &                       ice_Tice1_Ave,  myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_Tice2-T.',  suff,
c    &                       ice_Tice2_Ave,  myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_snowPr-T.', suff,
c    &                       ice_snowPr_Ave, myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_albedo-T.', suff,
c    &                       ice_albedo_Ave, myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_flx2oc-T.', suff,
c    &                       ice_flx2oc_Ave, myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_frw2oc-T.', suff,
c    &                       ice_frw2oc_Ave, myIter,myThid)
c      CALL WRITE_FLD_XY_RL('ICE_salFx-T.',  suff,
c    &                       ice_salFx_Ave,  myIter,myThid)
c      IF ( fluidIsWater ) THEN
c       CALL WRITE_FLD_XY_RL('ICE_flxAtm-T.', suff,
c    &                       ice_flxAtm_Ave, myIter,myThid)
c       CALL WRITE_FLD_XY_RL('ICE_frwAtm-T.', suff,
c    &                       ice_frwAtm_Ave, myIter,myThid)
c      ENDIF
c      IF ( stepFwd_oceMxL ) THEN
c       CALL WRITE_FLD_XY_RL('ICE_tMxL-T.', suff,
c    &                       ice_tMxL_Ave, myIter,myThid)
c       CALL WRITE_FLD_XY_RL('ICE_sMxL-T.', suff,
c    &                       ice_sMxL_Ave, myIter,myThid)
c      ENDIF


C--    Write all time-average variables in 1 file :
       _BARRIER
       _BEGIN_MASTER( myThid )

       IF ( thSIce_tave_mdsio ) THEN

C        find wether we are writing globalFile or tile-files:
         CALL GET_WRITE_GLOBAL_FLD( gf )

         WRITE(fn,'(A,I10.10)') 'thSIce_tave.', myIter

         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_fract_Ave,  1, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_iceH_Ave,   2, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_snowH_Ave,  3, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_Tsrf_Ave,   4, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_Tice1_Ave,  5, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_Tice2_Ave,  6, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_snowPr_Ave, 7, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_albedo_Ave, 8, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_flx2oc_Ave, 9, myIter, myThid )
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_frw2oc_Ave,10, myIter, myThid )
         kl = 11
         CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &        ice_salFx_Ave, kl, myIter, myThid )
         IF ( fluidIsWater ) THEN
           kl = kl+1
           CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &          ice_flxAtm_Ave,kl, myIter, myThid )
           kl = kl+1
           CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &          ice_frwAtm_Ave,kl, myIter, myThid )
         ENDIF
         IF ( stepFwd_oceMxL ) THEN
           kl = kl+1
           CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &          ice_tMxL_Ave,  kl, myIter, myThid )
           kl = kl+1
           CALL MDSWRITEFIELD( fn, writeBinaryPrec, gf, 'RL', 1,
     &          ice_sMxL_Ave,  kl, myIter, myThid )
         ENDIF

       ENDIF

#ifdef ALLOW_MNC
       IF ( useMNC .AND. thSIce_tave_mnc ) THEN
         DO i = 1,MAX_LEN_MBUF
           mncf(i:i) = ' '
         ENDDO
C        C             123456789 1
         mncf(1:11) = 'thsice_tave'
         CALL MNC_CW_SET_UDIM(fn, -1, myThid)
         CALL MNC_CW_RL_W_S('D',fn,0,0,'T', myTime, myThid)
         CALL MNC_CW_SET_UDIM(fn, 0, myThid)
         CALL MNC_CW_I_W_S('I',fn,0,0,'iter', myIter, myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'fract_ave', ice_fract_Ave, myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'iceH_ave',  ice_iceH_Ave,  myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'snowH_ave', ice_snowH_Ave, myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'Tsrf_ave',  ice_Tsrf_Ave,  myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'Tice1_ave', ice_Tice1_Ave, myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'Tice2_ave', ice_Tice2_Ave, myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'snowPr_ave',ice_snowPr_Ave,myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'albedo_ave',ice_albedo_Ave,myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'flx2oc_ave',ice_flx2oc_Ave,myThid)
         CALL MNC_CW_RL_W('D',fn,0,0,'frw2oc_ave',ice_frw2oc_Ave,myThid)
         IF ( fluidIsWater ) THEN
           CALL MNC_CW_RL_W('D',fn,0,0,'flxAtm_ave',
     &          ice_flxAtm_Ave,myThid)
           CALL MNC_CW_RL_W('D',fn,0,0,'frwAtm_ave',
     &          ice_frwAtm_Ave,myThid)
         ENDIF
         IF ( stepFwd_oceMxL ) THEN
           CALL MNC_CW_RL_W('D',fn,0,0,'tMxL_ave',ice_tMxL_Ave,myThid)
           CALL MNC_CW_RL_W('D',fn,0,0,'sMxL_ave',ice_sMxL_Ave,myThid)
         ENDIF
       ENDIF
#endif

       _END_MASTER( myThid )
       _BARRIER
         
C--    Reset averages to zero
       DO bj = myByLo(myThid), myByHi(myThid)
        DO bi = myBxLo(myThid), myBxHi(myThid)
         CALL TIMEAVE_RESET(ice_fract_Ave,  1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_iceH_Ave,   1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_snowH_Ave,  1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_Tsrf_Ave,   1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_Tice1_Ave,  1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_Tice2_Ave,  1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_snowPr_Ave, 1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_flx2oc_Ave, 1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_frw2oc_Ave, 1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_salFx_Ave,  1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_flxAtm_Ave,   1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_frwAtm_Ave,   1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_albedo_Ave, 1, bi,bj, myThid)
         CALL TIMEAVE_RESET(ice_tMxL_Ave,   1, bi, bj, myThid)
         CALL TIMEAVE_RESET(ice_sMxL_Ave,   1, bi, bj, myThid)
         ice_timeAve(1,bi,bj) = 0.
        ENDDO
       ENDDO
         
      ENDIF

#endif /* ALLOW_TIMEAVE */

C     Dump a snap-shot of main state variables:
      IF ( 
     &  DIFFERENT_MULTIPLE(thSIce_diagFreq,myTime,deltaTClock)
     &     .OR. myTime.EQ.startTime
     &     .OR. myTime.EQ.endTime ) THEN

        IF ( thSIce_snapshot_mdsio ) THEN
          
          WRITE(suff,'(I10.10)') myIter
          
          CALL WRITE_FLD_XY_RL('ice_fract.',  suff,iceMask,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_iceH.',   suff,iceheight,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_snowH.',  suff,snowheight,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_Tsrf.',   suff,Tsrf,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_Tice1.',  suff,Tice1,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_Tice2.',  suff,Tice2,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_Qice1.',  suff,Qice1,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_Qice2.',  suff,Qice2,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_snowPrc.',suff,snowPrc,
     &         myIter,myThid)
          CALL WRITE_FLD_XY_RL('ice_snowAge.',suff,snowAge,
     &         myIter,myThid)
          IF ( stepFwd_oceMxL ) THEN
            CALL WRITE_FLD_XY_RL('ice_tOceMxL.', suff, tOceMxL,
     &           myIter,myThid)
            CALL WRITE_FLD_XY_RL('ice_sOceMxL.', suff, sOceMxL,
     &     myIter,myThid)
          ENDIF
          
        ENDIF

#ifdef ALLOW_MNC
        IF ( thSIce_snapshot_mnc ) THEN
          DO i = 1,MAX_LEN_MBUF
            mncf(i:i) = ' '
          ENDDO
C         C             123456789 12345
          mncf(1:15) = 'thsice_snapshot'
          CALL MNC_CW_SET_UDIM(mncf, -1, myThid)
          CALL MNC_CW_I_W_S('I',mncf,0,0,'iter', myIter, myThid)
          CALL MNC_CW_SET_UDIM(mncf, 0, myThid)
          CALL MNC_CW_RL_W_S('D',mncf,0,0,'T', myTime, myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'iceMask',   iceMask,   myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'iceHeight', iceHeight, myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'snowHeight',snowHeight,myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'Tsrf',      Tsrf,      myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'Tice1',     Tice1,     myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'Tice2',     Tice1,     myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'Qice1',     Qice1,     myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'Qice2',     Qice2,     myThid)
          CALL MNC_CW_RL_W('D',mncf,0,0,'snowAge',   snowAge,   myThid)
          IF ( stepFwd_oceMxL ) THEN
            CALL MNC_CW_RL_W('D',mncf,0,0,'tOceMxL',tOceMxL,myThid)
            CALL MNC_CW_RL_W('D',mncf,0,0,'sOceMxL',sOceMxL,myThid)
          ENDIF
        ENDIF
#endif  /* ALLOW_MNC */


      ENDIF
         
      IF ( thSIce_monFreq.GT. 0. _d 0 )
     &    CALL THSICE_MONITOR( myTime, myIter, myThid )

#endif /* ALLOW_THSICE */
      
      RETURN
      END