      subroutine shap_indx(npt,nxp,nyp,mask,isx,nfxk,nfpxk,nfyk,ifx,ifpx,ify)
c-------------------------------------------------------------------------

c     nfx     : number of contiguous latitudinal segments 
c     ifx(1,i): starting point (x-compressed index) of i'th segment  
c     ifx(2,i): length of i-th segment
c     ifx(3,i): highest order filter used in i'th segment

c     nfpx    : number of segments with periodic overlap
c     ifx(1,i): starting point western portion of segment 
c     ifx(2,i): length of western portion of segment
c     ifx(3,i): starting point eastern portion of segment 
c     ifx(4,i): length of eastern portion of segment
c     ifx(5,i): highest order filter used in segment

c     nfy     : number of contiguous longitudinal segments 
c     ify(1,j): starting point (y-compressed index) of j'th segment
c     ify(2,j): length of j'th segment
c     ify(3,j): highest order filter used in j'th segment

      include "comm_new.h"
      include 'comm_para.h'
      dimension mask(nxp,1),isx(1),ifx(3,1),ifpx(5,1),ify(3,1)
      logical prev, curr, peri, ZER
      common /shap_c25/ cs25(10), ZER, s_coef
      common /new_filt/  MAXFO, nbx, nby, nfx, nfpx, nfy

      do i = 1, MAXFO
         cs25(i) = 1./(4.**i)
      enddo
 
      nfxk  = 0
      nfpxk = 0
      do irow = 1, nyp

         ista = 1
         prev = (mask(1, irow) .eq. 0)
         peri = (.not. prev) .and. (iglob .ne. 0)  

         do icol = 2, nxp
            
            curr = (mask(icol,irow) .eq. 0)
            if (curr .ne. prev) then
                if (peri .and. ista .eq. 1) then
                  nfpxk = nfpxk + 1
                  ifpx(1,nfpxk) = mask(1,irow)
                  ifpx(2,nfpxk) = icol - 1
               elseif (curr) then
                  nfxk = nfxk + 1
                  ifx(1,nfxk) = mask(ista,irow)
                  ifx(2,nfxk) = icol-ista
                  ifx(3,nfxk) = min((icol-ista-1)/2,MAXFO)
               endif
               ista = icol
               prev = curr
            endif
         enddo

         if (.not. curr) then
            if (peri) then
               if (ista .eq. 1) then
                  nfpxk = nfpxk + 1
                  ifpx(1,nfpxk) = mask(1,irow)
                  ifpx(2,nfpxk) = nxp
                  ifpx(3,nfpxk) = 0
                  ifpx(4,nfpxk) = 0
                  ifpx(5,nfpxk) = MAXFO
               else
                  ifpx(3,nfpxk) = mask(ista,irow)
                  ifpx(4,nfpxk) = mask(nxp,irow) - mask(ista,irow) + 1
                  ifpx(5,nfpxk) = min((ifpx(2,nfpxk)+ifpx(4,nfpxk)-1)/2,MAXFO)
               endif
            else
               nfxk = nfxk + 1
               ifx(1,nfxk) = mask(ista,irow)
               ifx(2,nfxk) = nxp-ista+1
               ifx(3,nfxk) = min((nxp-ista)/2,MAXFO)
            endif
         endif
      enddo

      nfyk  = 0
      do icol = 1, nxp

         jsta = 1
         prev = (mask(icol, 1) .eq. 0)
         do irow = 2, nyp
            
            curr = (mask(icol,irow) .eq. 0)
            if (curr .ne. prev) then
               if (curr) then
                  nfyk = nfyk + 1
                  ify(1,nfyk) = isx(mask(icol,jsta))
                  ify(2,nfyk) = irow-jsta
                  ify(3,nfyk) = min((irow-jsta-1)/2,MAXFO)
               endif
               jsta = irow
               prev = curr
            endif
         enddo
         if (.not. curr) then
            nfyk = nfyk + 1
            ify(1,nfyk) = isx(mask(icol,jsta))
            ify(2,nfyk) = nyp-jsta+1
            ify(3,nfyk) = min((nyp-jsta)/2,MAXFO)
         endif
      enddo

      return
      end

      subroutine out_mean
c--------------------------
      include 'comm_data.h'
      include 'comm_new.h'
      include 'comm_bio.h'
      include 'comm_para.h'
      include 'comm_main.h'

      common /mean_comm/ nmcount, nmcount_bio, nmcount_mosf, nmcount_conv

      coef_dens = 1. / (1000. + SITUD_BOT)

      coef = 1./real(max(1,nmcount))
      nmcount = 0

      do i = 1, nz*nxyc
         um(i) = coef*um(i) 
         vm(i) = coef*vm(i) 
         wm(i) = coef*wm(i) 
         hm(i) = coef*hm(i) 
         tm(i) = coef*tm(i) 
         salm(i)  = coef*salm(i) 
         densm(i) = coef*densm(i) 
      enddo

      do i = 1, nxyc
         hflx_mean(i) = coef*hflx_mean(i) 
         sflx_mean(i) = coef*sflx_mean(i) 
         plat_mean(i) = coef*plat_mean(i) 
         defi_mean(i) = coef*defi_mean(i) 
         plon_mean(i) = coef*plon_mean(i) 
         sens_mean(i) = coef*sens_mean(i) 
         psim(i) = coef*psim(i) 
         hc_top(i) = hm(i)*tm(i)
         sc_top(i) = hm(i)*salm(i)
         uc_top(i) = hm(i)*um(i)
         vc_top(i) = hm(i)*vm(i)
         dh_top(i) = hm(i)*densm(i)* coef_dens
         do k = 2, nzi(i)
            ik = (k-1)*nxyc + i
            ikm = ik - nxyc
            hc_top(ik) = hc_top(ikm)  + hm(ik) * tm(ik)
            sc_top(ik) = sc_top(ikm)  + hm(ik) * salm(ik)
            uc_top(ik) = uc_top(ikm)  + hm(ik) * um(ik)
            vc_top(ik) = vc_top(ikm)  + hm(ik) * vm(ik)
            dh_top(ik) = dh_top(ikm)  + hm(ik) * densm(ik)* coef_dens
         enddo
         do k = nzi(i)+ 1, nz
            ik = (k-1)*nxyc + i
            hc_top(ik) = rnan
            sc_top(ik) = rnan
            uc_top(ik) = rnan
            vc_top(ik) = rnan
            dh_top(ik) = rnan
         enddo
      enddo

      if (initq .eq. 8) then
         do i = 1, nxp*nyp
            amhum_mean(i)= coef*amhum_mean(i) 
            amth_mean(i) = coef*amth_mean(i) 
         enddo
      endif
      if (use_ice) then
         do i = 1, nxyc*ntrac_sur
            trsm(i) = coef*trsm(i) 
         enddo
      endif

      do in = 1, nz*nxyc
         do k = 1, ntrac
            i = (k-1)*nz*nxyc + in
            trm(i) = coef*trm(i) 
         enddo
      enddo

c  works to here

      call zonal_int(nxp,nyp,1,nxyc,mask,hflx_mean,hflx_y,sxm)
      call zonal_int(nxp,nyp,1,nxyc,mask,sflx_mean,sflx_y,sxm)

c  if means are being saved, then only make restart when means have been
c     written 

      nskip = 0
      call keep_rstrt(1,nstep, nstep)
      call keep_rstrt(0,nstep, nstep)

      return
      end

      subroutine keep_rstrt (itype, nstep, nskip)
c--------------------------------------------------------------------
      include 'comm_new.h'
      include 'comm_bio.h'
      include 'comm_data.h'
      include 'diffiso.h'
      common /new_save/ iword, nruns, nscpu, nswll, inf(100), rnf(100)
      common /mean_comm/ nmcount, nmcount_bio, nmcount_mosf, nmcount_conv
      equivalence (inf(20),kd_sta),
     *            (inf(21),kd_xy),  (inf(22),kd_z),  (inf(23),kd_time), 
     *            (inf(24),kd_seg), (inf(25),kd_uv), (inf(26),kd_temp), 
     *            (inf(27),kd_salt),(inf(28),kd_tr), (inf(29),kd_phi), 
     *            (inf(30),kd_mask),(inf(32),kd_psi), (inf(33),kd_ice), 
     *            (inf(100),kd_end)

      if (nskip.eq.0) return
      if (mod(nstep, nskip) .ne. 0) goto 101

      nptz  = inf(4)*inf(6)
      npt   = inf(6)
      nptzp  = inf(4)*inf(6) +inf(6)

      kd = kd_uv
      call wrda (itype, iword, kd, nptz, u) 
      kd = iwrda (iword, nptz, v) 
      kd = iwrda (iword, nptz, uc) 
      kd = iwrda (iword, nptz, vc) 
      kd = iwrda (iword, nptz, h) 

      kd_temp = kd
      kd = iwrda (iword, nptz, t)
      kd_salt = kd
      kd   = iwrda (iword, nptz, sal) 
      kd   = iwrda (iword, nptz, dens) 

      if ( ibaro.ne.0 ) then
         kd_psi = kd
         kd = iwrda (iword, npt, psi) 
      endif
      if ( use_ice ) then
         kd_ice = kd
         kd = iwrda (iword, npt, trs(i_hice)) 
         kd = iwrda (iword, npt, trs(i_cice)) 
         kd = iwrda (iword, npt, trs(i_thice)) 
      endif

      if ( use_trac.or.use_bio ) then
         kd_tr = kd
         kd    = iwrda (iword, nptz*ntrac_tot, tr) 
      endif

      kd_end = kd

      inf(11) = nstep
      inf(19) = nmcount

      rnf(5) = f0

      rnf(10) = ekf1 
      rnf(11) = epf1
      rnf(12) = hcf1
      rnf(13) = wcf1
      rnf(14) = vlf1
      rnf(15) = budget_sp
      rnf(16) = budget_ri
      rnf(17) = budget_eb
      rnf(18) = budget_pr
      rnf(19) = coef_bud
      rnf(20) = budget_re 

      do k = 1, inf(4)+1
         rnf(20+k) = hsave(k) 
      enddo

  101 continue

      inf(18) = nstep
      inf(13) = nscpu + ipast_scpu()
      inf(14) = nswll + ipast_swll()

      kd = 0
      call wrda (itype, iword, kd, 100, inf) 
      kd = iwrda (iword, 100, rnf) 

      call flda (itype)
      return
      end

