      SUBROUTINE GETBC2
Status: R

C
#include "slgcm.slgprm"
#include "slgcm.slgwrk"
#include "slgcm.slgcon"
C
      DIMENSION PHIX(IMP3,JM)
      EQUIVALENCE (PHIX(1,1),PHIS(1))
C
      REWIND (KBC)
      READ (KBC) ((PHIX(I,J),I=1,IMP1), J=1,JM)
C
      DO 22 I=1,IJM0M0
C     PHIS(I) = 0.0
      IF (PHIS(I) .LT. 0.0) PHIS(I) = 0.0
   22 CONTINUE
      PHIX(1,:) = PHIX(IMP1,:)
      CALL WRAP(1,IJM0M0,PHIS)
C
      CALL MAP(PHIS,IMP3,JM,'PHIS',0.0)
C     CALL SMTH2D(PHIS,WK1,IMP3,JM)
C     CALL MAP(PHIS,IMP3,JM,'PHIS',0.0)
C
      RETURN
      END
      SUBROUTINE CLKOFF(TIME)
C
C*********************************************************************
C************************* SLSI GCM **********************************
C********************* SUBROUTINE CLKOFF *****************************
C********************** 11 AUGUST 1989  ******************************
C*********************************************************************
C
C
      REAL CLKSTK(100)
      DATA JSTK/0/
C     SAVE
C
      IF (JSTK .LE. 0) RETURN
      TIME = TIME + (SECOND() - CLKSTK(JSTK))
      JSTK = JSTK - 1
      RETURN
C
      ENTRY CLKON
      JSTK = JSTK + 1
      CLKSTK(JSTK) = SECOND()
      RETURN
C
      END
      SUBROUTINE WRITE1(KTO, LEN, FLD)
C
C
C*********************************************************************
C************************* SLSI GCM **********************************
C********************* SUBROUTINE WRITE1 *****************************
C********************** 22 AUGUST  1989 ******************************
C*********************************************************************
C
C
      DIMENSION FLD(LEN)
C
      WRITE(KTO) FLD
C
      RETURN
      END
      SUBROUTINE TRNSDA(IM, IMP3, JM, FIN, FOU)
C
C
C*********************************************************************
C************************* SLSI GCM **********************************
C********************* SUBROUTINE TRNSDA *****************************
C********************** 22 AUGUST  1989 ******************************
C*********************************************************************
C
C
      DIMENSION FIN(IMP3,JM), FOU(IM,JM)
C     SAVE
C
      DO 10 J=1,JM
      DO 10 I=1,IM
      FOU(I,J) = FIN(I,J)
   10 CONTINUE
C
      RETURN
      END
      SUBROUTINE CONSL
C
#include "slgcm.slgprm"
#include "slgcm.slgclm"
#include "slgcm.slgwrk"
#include "slgcm.slgcon"
C
C
      DIMENSION SGS(17)
      DIMENSION TX1(LM*LMP2),  TX2(LMP1), REI(LM*2), REC(LM*LM*2)
     *,         ETEM(LM,LM)
      COMPLEX   ECP(LM,LM), EIG(LM)
      EQUIVALENCE (REI(1),EIG(1)), (REC(1),ECP(1,1))
C
      DATA SGS/1.000,  20.0,  40.0,  70.0, 105.0, 145.0,  200.0
     *,        270.0, 350.0, 435.0, 525.0, 615.0, 705.0
     *,        790.0, 875.0, 950.0, 1000.0/
c     DATA SGS/1.000,  20.0,  40.0,  70.0, 110.0, 160.0,  220.0
c    *,        290.0, 370.0, 460.0, 555.0, 650.0, 740.0
c    *,        820.0, 890.0, 950.0, 1000.0/
C
C     SAVE
C
      DO 220 J=1,JM
      PH(J,2) = -PIO2  + (J-1)*DPH
      PH(J,1) = PH(J,2) + 0.5*DPH
  220 CONTINUE
C
      DO 225 K=1,2
      DO 225 J=1,JM
      SPH(J,K) = SIN(PH(J,K))
      CPH(J,K) = COS(PH(J,K))
      DX(J,K) = AE*DLM*CPH(J,K)
  225 CONTINUE
      CPH(1,2)  = 0.0
      CPH(JM,2) = 0.0
      DX(1,2)   = 0.0
      DX(JM,2)  = 0.0
C
      DO 226 J=1,JMM1
      RH1JK(J,1) = ONE / (AE*CPH(J,1))
  226 CONTINUE
      DO 227 J=2,JMM1
      RH1JK(J,2) = ONE / (AE*CPH(J,2))
  227 CONTINUE
      RH1JK(1,2)  = 0.0
      RH1JK(JM,1) = 0.0
      RH1JK(JM,2) = 0.0
C
      DO 5010 J=1,JM
      II = (J-1)*IMP3
      DO 5009 I=1,IMP3
      PHV(II+I) =  PH(J,1)
      PHU(II+I) =  PH(J,2)
      SPHV(II+I) = SPH(J,1)
      SPHU(II+I) = SPH(J,2)
      CPHV(II+I) = CPH(J,1)
      CPHU(II+I) = CPH(J,2)
      RH1V(II+I) = RH1JK(J,1)
      RH1U(II+I) = RH1JK(J,2)
      WK1(II+I) =  DX(J,2)
 5009 CONTINUE
 5010 CONTINUE
C
      DO 215 I=1,IJM2M0
      TEM       = 1.0 / ((DY*DY)  * CPHU(I1+I))
      CD1(I1+I) = 1.0 / (WK1(I1+I)*WK1(I1+I))
      CD2(I1+I) = TEM * CPHV(I1+I)
      CD3(I1+I) = TEM * CPHV(I)
  215 CONTINUE
C
C  K=1 GIVES THE ANGLE LAMDA AT THE U-POINTS
C  K=2 GIVES THE ANGLE LAMDA AT THE H-POINTS AND THE V-POINTS
C
      DO 221 I=1,IM
      BLMIK(I,2) = (I-1)*DLM
      BLMIK(I,1) = BLMIK(I,2) +  0.5*DLM
  221 CONTINUE
      DO 222 I=1,3
      BLMIK(IM+I,1) = BLMIK(I,1)
      BLMIK(IM+I,2) = BLMIK(I,2)
  222 CONTINUE
      DO 223 K=1,2
      DO 223 I=1,IM
      CLM(I,K) = COS(BLMIK(I,K))
      SLM(I,K) = SIN(BLMIK(I,K))

      CLM(I,K+2) = COS(2.0*BLMIK(I,K))
      SLM(I,K+2) = SIN(2.0*BLMIK(I,K))

      CLM(I,K+4) = COS(3.0*BLMIK(I,K))
      SLM(I,K+4) = SIN(3.0*BLMIK(I,K))
  223 CONTINUE
C
C
      DO 230 J=2,JMM1
      CDIVU(J)  = ONE  / DX(J,2)
      CDIVV(J)  = ONE  / (DY*CPH(J,2))
  230 CONTINUE
C
C
      AREA   = 0.0
      DO 2005 J=2,JMM1
      DO 2005 I=1,IM
      AREA   = AREA + (DX(J,2)*DY)
 2005 CONTINUE
      ARE2   = AREA + IM * DY * DX(1,1) * 0.5
      AREA   = AREA + 2.0 * PI * DY*DY*0.25
c     PRINT *,'ARE2=',ARE2
C
      DO 2007 J=2,JMM1
      II = (J-1)*IMP3
      DO 2007 I=1,IMP3
      CDU(II+I) = CDIVU(J)
      CDV(II+I) = CDIVV(J)
 2007 CONTINUE
C
      DO 8888 J=1,JM
      II = (J-1)*IMP3
      DO 8888 I=1,IMP3
C8888 INDI(II+I) = I
      xcr(II+I,1) = (I - 1.0) * dlm
      xcr(II+I,2) = (I - 0.5) * dlm
 8888 continue
C
C      VERTICAL LEVELS
C
      DTBDSG(:) = 80.0
c     DTBDSG(:) =  0.0
      TTOP   = 210.0
c     TTOP   = 300.0
C     DTBDSG(:) = 70.0
C     TTOP   = 240.0
C     TTOP   = 200.0
C     DTBDSG(:) = 100.0
C     TTOP   = 220.0
C     DTBDSG(:) = 0.0
C     TTOP   = 250.0
C     SIG(1) = 0.000001
C
      IF (EQSIG) THEN
         SIG(1) = 0.01
         DO 490 L=1,LM

C        TE1 = FLOAT(L) / FLOAT(LM)
C        SIG(L+1) = 0.75 * TE1 + 1.75 * TE2 * TE1 - 1.5 * TE2 * TE2

         SIG(L+1) = SIG(L) + (1.0-SIG(1)) / FLOAT(LM)
  490    CONTINUE
      ELSE
         DO 491 L=1,LMP1
         SIG(L) = SGS(L) * 0.001
  491    CONTINUE 
      ENDIF
C
C     SIG(2)  = 0.125
C     SIG(LMP1) = 1.0
C     DO 490 L=2,LMM1
C     SIG(L+1) = SIG(L) + (0.875-SIG(2)) / FLOAT(LMM2)
C 490 CONTINUE
C
      DO 494 L=1,LMP1
      PRJ(L) = (SIG(L)/1000.0) ** RKAP
  494 CONTINUE
C
      DO 500 L=1,LM
      DSG(L) = SIG(L+1) - SIG(L)
      SGB(L) = (SIG(L+1) + SIG(L)) * HALF
      TBR(L) = TTOP  + DTBDSG(L) * (SGB(L) - SIG(1))
      GMA(L) = RKAP / SGB(L) - DTBDSG(L) / TBR(L)
      DSL(L) = ALOG(SIG(L+1)/SIG(L))
C     VAD(L) = (1.0-VADVF) / (DSG(L)*TBR(L))
  500 CONTINUE
C
      DO 495 L=2,LM
      SGX(L-1) = 1.0 / (SGB(L) - SGB(L-1))
      VDF(L) = VDIFC * 2.0 * (SIG(L) * GRAV / RGAS) ** 2
     *       /             (DSG(L)+DSG(L-1))
      TEM     = VDIFB
      IF (SIG(L) .GT. 0.75)
     *   TEM = VDIFB + (VDIFQ-VDIFB)*(SIG(L)-0.75)/(SIG(LM)-0.75)

      VDQ(L) = TEM * 8.0 * (SIG(L) * GRAV / (RGAS*PRJ(L))) ** 2
     *       /             (DSG(L)+DSG(L-1))
  495 CONTINUE
C
      SGX(LM) = 1.0
      VDF(1)    = 0.0
      VDF(LMP1) = 0.0
      VDQ(1)    = 0.0
      VDQ(LMP1) = 0.0
C
      DO 502 L=1,LM
c     PRH(L) = (SGB(L)/1000.0) ** RKAP
      PRH(L  ) = (SIG(L+1)*PRJ(L+1)-SIG(L)*PRJ(L)) /
     *         ( (SIG(L+1)-SIG(L))*(RKAP+ONE) )
      FPK(L  ) = CP * (PRJ(L+1)-PRJ(L))
      HPK(L  ) = CP * (PRJ(L+1)-PRH(L))
C
      HFBDSG(L)= HALF / DSG(L)
      CNV(L) = (CMB2PA*DSG(L)/GRAV)
      ODP(L) = 1./(DSG(L)*PRH(L))
      ODS(L) = 1./(DSG(L)       )

  502 CONTINUE
C
      DO 504 L=1,LMM1
      APM(L) = PRJ(L  )*DSG(L  )/(PRJ(L)*DSG(L)+PRJ(L+1)*DSG(L+1))
      AMP(L) = PRJ(L+1)*DSG(L+1)/(PRJ(L)*DSG(L)+PRJ(L+1)*DSG(L+1))
      QPM(L) = DSG(L  )/(DSG(L)+DSG(L+1))
      QMP(L) = DSG(L+1)/(DSG(L)+DSG(L+1))
      CPM(L) = 1.0 - CVADJ*(PRJ(L)-PRJ(L+1))*QMP(L) / PRJ(L)
      CMP(L) = 1.0 + CVADJ*(PRJ(L)-PRJ(L+1))*QPM(L) / PRJ(L+1)
      AKL(L)  = (PRJ(L+1)-PRH(L  ))/(PRH(L+1)-PRH(L))
C
      EFQ(L) = (ODS(L) + ODS(L+1))
      EFP(L) = (ODP(L) + ODP(L+1))
C
  504 CONTINUE
C
      PRINT *,' AMP=',AMP
      PRINT *,' APM=',APM
      PRINT *,' CMP=',CMP
      PRINT *,' CPM=',CPM
C
      TEM = 0.0
      DO 503 L=LM,1,-1
      TEM1    = DSL(L)*TBR(L) * RGAS
      GBRM(L) = TEM + TEM1*0.5
      TEM     = TEM + TEM1
  503 CONTINUE
C
      PRINT *,'SIG:',SIG
      PRINT *,'SGB:',SGB
      PRINT *,'DSG:',DSG
      PRINT *,'DSL:',DSL
      PRINT *,'TBR:',TBR
      PRINT *,'GMA:',GMA
      PRINT *,'PRJ:',PRJ
      PRINT *,'PRH:',PRH
      PRINT *,'GBRM:',GBRM
      PRINT *,' IM=',IM,' JM=',JM,' LM=',LM, ' PHMAXD=',PHMAXD
C
C COMPUTATION OF S(K)
C
      TE1 = ZERO
      DO 510 L=LM,1,-1
      TE2    = (RKAP - GMA(L)*(SGB(L)-SIG(1))) * DSL(L) * TBR(L)
      HFC(L) = TE1 + TE2 * HALF
      TE1    = TE1 + TE2
  510 CONTINUE
      DO 515 L=1,LM
      HFC(L) = (TBR(L) + HFC(L)) / (SIG(LMP1)-SIG(1))
      TX1(L) = GMA(L) * TBR(L) * DSL(L)
  515 CONTINUE
      TX2(LM) = TX1(LM)
      DO 520 L=LMM1,1,-1
      TX2(L) = TX2(L+1) + TX1(L)
  520 CONTINUE
      TX2(LMP1) = 0.0
C
C     COMPUTATION OF MATRIX M(K,KP)
C
      DO 550 K=1,LM
      KM1 = K - 1
      KP1 = K + 1
      IF (KM1 .GT. 0) THEN
         DO 530 L=1,KM1
         GDM(K,L) = (HFC(K) + TX1(K)*HALF + TX2(K+1)) * DSG(L)
  530    CONTINUE
      ENDIF
C
      GDM(K,K) = (HFC(K) + TX1(K)*0.25 + TX2(KP1)) * DSG(K)
      IF (KP1 .LE. LM) THEN
         DO 540 L=KP1,LM
         GDM(K,L) = (HFC(K) + TX1(L)*HALF + TX2(L+1)) * DSG(L)
  540    CONTINUE
      ENDIF
C
C     GDM(K,LM) = (HFC(LM) + TX1(LM)*HALF) * DSG(LM)
  550 CONTINUE

      DO 551 L=1,LM
      DO 551 K=1,LM
      GDM(K,L) = GDM(K,L) * RGAS
  551 CONTINUE

C
      DO 552 L=1,LM
      DO 552 K=1,LM
      ETEM(K,L) = GDM(K,L)
  552 CONTINUE
C
C IMSL ROUTINE FOR EIGENVALUES AND EIGENVECTORS
C
C     CALL EIGRF(ETEM, LM, LM, 2, REI, REC, LM, TX1, IER)
      CALL EVCRG(LM, ETEM, LM,  REI, REC, LM)
C
      DO 555 L=1,LM
      ELM(L) = REAL(EIG(L))
      DO 555 K=1,LM
      E(K,L)    = REAL(ECP(K,L))
      ETEM(K,L) = E(K,L)
  555 CONTINUE
C
      DO 570 KL=1,LM
      DO 560 L=1,LM
      TX1(L) = 0.0
      DO 556 K=1,LM
      TX1(L) = TX1(L) + GDM(L,K) * E(K,KL)
  556 CONTINUE
  560 CONTINUE
c     PRINT *,' AT KL=',KL,' EIGEN VAL=',ELM(KL)
c     PRINT *,' ME=',(TX1(K),K=1,LM)
c     PRINT *,' E=',(E(K,KL),K=1,LM)
  570 CONTINUE
C
C     CALCULATE INVERSE OF ETEM; OUTPUT IN EI
C
C     CALL LINV1F(ETEM, LM, LM, EI, 0, TX1, IER1)
      CALL LINRG(LM, ETEM, LM, EI, LM)
C
C
      IF (LM .GE. 4) THEN
         DO 2010 L=2,LMM2
         S1A(L,1) = 1.0 / ((SGB(L-1)-SGB(L))  *(SGB(L-1)-SGB(L+1))
     *                                        *(SGB(L-1)-SGB(L+2)))
         S1A(L,2) = 1.0 / ((SGB(L)-SGB(L-1))  *(SGB(L)-SGB(L+1))
     *                                        *(SGB(L)-SGB(L+2)))
         S1A(L,3) = 1.0 / ((SGB(L+1)-SGB(L-1))*(SGB(L+1)-SGB(L))
     *                                        *(SGB(L+1)-SGB(L+2)))
         S1A(L,4) = 1.0 / ((SGB(L+2)-SGB(L-1))*(SGB(L+2)-SGB(L))
     *                                        *(SGB(L+2)-SGB(L+1)))
 2010    CONTINUE
      ENDIF
C
      RETURN
      END
      SUBROUTINE GLBMN(IMP3, JM, U, V, T, SGD, CPH)
C
      DIMENSION U(IMP3,JM), V(IMP3,JM), T(IMP3,JM), SGD(IMP3,JM)
     *,         CPH(JM,2), GLBM(4)
C
      IM   = IMP3 - 3
      JMM1 = JM - 1
      DO 10 I=1,4
      GLBM(I) = 0.0
   10 CONTINUE
      DXPU = 0.0
      DXPV = 0.0
C
      DO 30 J=2,JMM1
      DO 20 I=1,IM
      GLBM(1) = GLBM(1) + U(I,J)  *CPH(J,2)
      GLBM(3) = GLBM(3) + T(I,J)  *CPH(J,2)
      GLBM(4) = GLBM(4) + SGD(I,J)*CPH(J,2)
      DXPU    = DXPU    +          CPH(J,2)
   20 CONTINUE
   30 CONTINUE
      DO 50 J=1,JMM1
      DO 40 I=1,IM
      GLBM(2) = GLBM(2) + V(I,J)*CPH(J,1)
      DXPV    = DXPV    +        CPH(J,1)
   40 CONTINUE
   50 CONTINUE
C
      DXPUI = 1.0 / DXPU
      GLBM(1)  = GLBM(1) * DXPUI
      GLBM(2)  = GLBM(2) / DXPV
      GLBM(3)  = GLBM(3) * DXPUI
      GLBM(4)  = GLBM(4) * DXPUI
C
      PRINT *,' GLBU=',GLBM(1),' GLBV=',GLBM(2),' GLBT=',GLBM(3)
     *,       ' GLBS=',GLBM(4)
C
      RETURN
      END
      SUBROUTINE GLBMNP(IMP3, JM, SRP, DX, DY, AREA, PI, GLBM)
C
      DIMENSION SRP(IMP3,JM), DX(JM)
C
      IM   = IMP3 - 3
      JMM1 = JM - 1
      GLBM   = 0.0
      DO 3005 J=2,JMM1
      DO 3005 I=1,IM
      GLBM = GLBM + SRP(I,J) * (DX(J)*DY)
 3005 CONTINUE
      GLBM = GLBM + (SRP(1,1) + SRP(1,JM)) * PI * DY*DY*0.25
C
      GLBM = GLBM / AREA
C
      RETURN
      END
      SUBROUTINE GLMNP2(IMP3, JM, SRP, DX, DY, ARE2, DXV, GLBM)
C
      DIMENSION SRP(IMP3,JM), DX(JM)
C
      IM   = IMP3 - 3
      JMM1 = JM - 1
      DYBA = DY / ARE2
C
      GLBM = (SRP(1,1) + SRP(1,JM)) * (DYBA*DXV*0.25*IM)
      DO 3005 J=2,JMM1
      DO 3005 I=1,IM
      GLBM = GLBM + SRP(I,J) * (DX(J)*DYBA)
 3005 CONTINUE
C
      RETURN
      END
      SUBROUTINE CHKENG(IMP3, JM, LM, U, V, T, PS, CPH, DSG)
      PARAMETER (CP=1003.5)
C
      DIMENSION U(IMP3,JM,LM), V(IMP3,JM,LM), T(IMP3,JM,LM)
     *,         PS(IMP3,JM),   CPH(JM,2)
     *,         WK1(IMP3,JM), WK2(IMP3,JM), WK3(IMP3,JM), WK4(IMP3,JM)
     *,         TKE(50), TPE(50), DSG(LM)
      LOGICAL FIRST
      DATA FIRST/.TRUE./, TKE0/0.0/, TPE0/0.0/, TOE0/0.0/
C
      IM     = IMP3 - 3
      JMM1   = JM - 1
      IJM0M0 = IMP3*JM
      IJM1M0 = IMP3*JMM1
      IJM0M1 = IJM0M0 - 1
C
      TKE(:) = 0.0
      TPE(:) = 0.0
      TKEM   = 0.0
      TPEM   = 0.0
C
      WK1(:,:) = 0.0
      WK2(:,:) = 0.0
      DO 5 I=1,IJM0M1
      WK1(I,1) = PS(I,1) + PS(I+1,1)
    5 CONTINUE
      DO 6 I=1,IJM1M0
      WK2(I,1) = PS(I,1) + PS(I,2)
    6 CONTINUE
C
      DO 100 L=1,LM
C
      DO 10 I=1,IJM0M0
      WK3(I,1) = U(I,1,L) * U(I,1,L) * WK1(I,1)
      WK4(I,1) = V(I,1,L) * V(I,1,L) * WK2(I,1)
   10 CONTINUE
C
      DO 30 J=2,JMM1
      DO 20 I=1,IM
      TKE(L) = TKE(L) + WK3(I,J) * CPH(J,2)
      TPE(L) = TPE(L) + T(I,J,L) * CPH(J,2) * PS(I,J)
   20 CONTINUE
   30 CONTINUE
      DO 50 J=1,JMM1
      DO 40 I=1,IM
      TKE(L) = TKE(L) + WK4(I,J) * CPH(J,1)
   40 CONTINUE
   50 CONTINUE
C
      TKEM = TKEM + 0.5 * TKE(L) * DSG(L)
      TPEM = TPEM + CP  * TPE(L) * DSG(L)
  100 CONTINUE
      TOTE = TKEM + TPEM
C
      IF (.NOT. FIRST) THEN
          TKEM = TKEM / TKE0
          TPEM = TPEM / TPE0
          TOTE = TOTE / TOE0
          PRINT *,' TKEM=',TKEM,' TPEM=',TPEM,' TOTE=',TOTE
      ELSE
          TKE0  = TKEM
          TPE0  = TPEM
          TOE0  = TOTE
          FIRST = .FALSE.
          PRINT *,' TKE0=',TKE0,' TPE0=',TPE0,' TOE0=',TOE0
      ENDIF
C
      RETURN
      END
      SUBROUTINE PTOSIG(US,  VS,  GS,  RS, PS,  PSTR, LR
     *,                 UOV, VOV, TOV,  QOV, PIV, TSR
     *,                 TX1, TX2, TX3, TX4, TX5)
C
#include "slgcm.slgprm"
#include "slgcm.slgcon"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE PTOSIG *****************************
C********************** 13 AUGUST 1990 *******************************
C*********************************************************************
C
C
      DIMENSION TX1(IMP3,JM), TX2(IMP3,JM), TX3(IMP3,JM)
     *,         TX4(IMP3,JM), TX5(IMP3,JM)
C
      REAL  US(IMP3,JM,LR), VS(IMP3,JM,LR)
     *,     GS(IMP3,JM,LR), PS(IMP3,JM), RS(IMP3,JM,LR)
C
      REAL UOV(IMP3,JM,LM), VOV(IMP3,JM,LM), TOV(IMP3,JM,LM)
     *,    QOV(IMP3,JM,LM), PIV(IMP3,JM   ), PHIX(IMP3,JM)
     *,    TSR(IMP3,JM)
C
      REAL  PSTR(LR), PKSTR(100)
C     EQUIVALENCE (PHIS(1),PHIX(1,1))
C
C     CALL GETBC2
C     CALL GETBCS(SEASNL)
C
      DO 5 L=1,LR
      PKSTR(L) = (PSTR(L)/1000.0)**RKAP
5     CONTINUE
C
      DO 15 L=1,LM
      DO 10 I=1,IJM0M0
      UOV(I,1,L) = ZERO
      VOV(I,1,L) = ZERO
      TOV(I,1,L) = ZERO
      QOV(I,1,L) = ZERO
   10 CONTINUE
   15 CONTINUE
C
      do 17 l=1,lr
      do 17 i=1,ijm0m0
      if (us(i,1,l) .gt. 999.0) us(i,1,l) = 0.0
      if (vs(i,1,l) .gt. 999.0) vs(i,1,l) = 0.0
      if (rs(i,1,l) .gt. 200.0) rs(i,1,l) = 0.8
   17 continue
C
      DO 20 I=1,IJM0M0
      PIV(I,1)  = ZERO
      PHIX(I,1) = PHIS(I)
      TSR(I,1)  = ZERO
   20 CONTINUE
C
      CALL SMTHX(GS(1,JMM3,LR),TX1,IMP3,3)
      CALL SMTHX(GS(1,JMM3,LR-1),TX1,IMP3,3)

      CALL SMTHX(PS(1,JMM3   ),TX1,IMP3,3)
      CALL SMTHX(PS(1,JMM3   ),TX1,IMP3,3)

      CALL SMTHX(GS(1,2,LR),TX1,IMP3,3)
      CALL SMTHX(GS(1,2,LR-1),TX1,IMP3,3)

      CALL SMTHX(PS(1,2   ),TX1,IMP3,3)
      CALL SMTHX(PS(1,2   ),TX1,IMP3,3)
      CALL SMTHX(PS(1,2   ),TX1,IMP3,3)

C
      DO 100 L=1,LR-1
      DO 30 I=1,IJM0M0
      IF (PHIX(I,1).LT.GS(I,1,L) .AND. PHIX(I,1).GE.GS(I,1,L+1)) THEN
         TX1(I,1) = PHIX(I,1) - GS(I,1,L+1)
         TX2(I,1) = GS(I,1,L) - GS(I,1,L+1)
         TX1(I,1) = TX1(I,1)  / TX2(I,1)
         TX2(I,1) = PSTR(L)   * TX1(I,1)
         TX1(I,1) = ONE       - TX1(I,1)
         PIV(I,1) = PSTR(L+1) * TX1(I,1) + TX2(I,1)
      ENDIF
   30 CONTINUE
  100 CONTINUE
C
      DO 110 I=1,IJM0M0
      IF (PHIX(I,1).LT.GS(I,1,LR) .AND. PHIX(I,1).GT.ZERO) THEN
         TX1(I,1) = GS(I,1,LR) - PHIX(I,1)
         TX1(I,1) = TX1(I,1)   / GS(I,1,LR)
         TX2(I,1) = TX1(I,1)   * PS(I,1)
         TX1(I,1) = ONE        - TX1(I,1)
         PIV(I,1) = TX2(I,1)   + PSTR(LR)*TX1(I,1)
      ENDIF
  110 CONTINUE
C
      DO 120 I=1,IJM0M0
      IF (PHIX(I,1) .LE. ZERO) PIV(I,1) = PS(I,1)
  120 CONTINUE
C
      CALL WRAP(1,IJM0M0,PIV(1,1))
      CALL SMTHX(PIV(1,JMM3),TX1,IMP3,3)
      CALL SMTHX(PIV(1,2   ),TX1,IMP3,3)
      CALL WRAP(1,IJM0M0,PIV(1,1))
C
      TEM1 = 0.0
      TEM2 = 0.0
      DO 122 I=1,IM
      TEM1 = TEM1 + PIV(I,1)
      TEM2 = TEM2 + PIV(I,JM)
  122 CONTINUE
C
      TEM1 = TEM1 / IM
      TEM2 = TEM2 / IM
C
      DO 123 I=1,IMP3
      PIV(I,1)  = TEM1
      PIV(I,JM) = TEM2
  123 CONTINUE
C
      CALL GLMNP2(IMP3, JM, PIV, DX(1,2), DY, ARE2, DX(1,1), GLBP)
      PRINT *,' GLOBAL MEAN PRESSURE =',GLBP
C
      DO 130 I=1,IJM0M0
      PIV(I,1) = PIV(I,1) -GLBP + SRFPRS
      TX5(I,1) = PIV(I,1) ** RKAP
  130 CONTINUE
C
      DO 200 L=1,LM
      DO 140 I=1,IJM0M0
      TX1(I,1) = TX5(I,1)*PRH(L)
      uov(i,1,l) = us(i,1,1)
      vov(i,1,l) = vs(i,1,1)
      qov(i,1,l) = 0.0
  140 CONTINUE

      DO 300 K=1,LR-1
      DO 150 I=1,IJM0M0
      IF (TX1(I,1).GT.PKSTR(K) .AND. TX1(I,1).LE.PKSTR(K+1)) THEN
         TX2(I,1)  = PKSTR(K+1) - TX1(I,1)
         TX2(I,1)  = TX2(I,1) / (PKSTR(K+1)-PKSTR(K))
         TX3(I,1)  = ONE - TX2(I,1)
C
         TX4(I,1)   = US(I,1,K)   * TX2(I,1)
         UOV(I,1,L) = US(I,1,K+1) * TX3(I,1)
         UOV(I,1,L) = UOV(I,1,L)  + TX4(I,1)
C
         TX4(I,1)   = VS(I,1,K)   * TX2(I,1)
         VOV(I,1,L) = VS(I,1,K+1) * TX3(I,1)
         VOV(I,1,L) = VOV(I,1,L)  + TX4(I,1)
C
         TX4(I,1)   = RS(I,1,K)   * TX2(I,1)
         QOV(I,1,L) = RS(I,1,K+1) * TX3(I,1)
         QOV(I,1,L) = QOV(I,1,L)  + TX4(I,1)

      ENDIF
  150 CONTINUE
  300 CONTINUE
      DO 160 I=1,IJM0M0
      IF (TX1(I,1) .GT. PKSTR(LR)) THEN
         UOV(I,1,L) = US(I,1,LR)
         VOV(I,1,L) = VS(I,1,LR)
         QOV(I,1,L) = RS(I,1,LR)
      ENDIF
  160 CONTINUE
C
      DO 170 I=1,IJM0M1
      TX1(I,1)   = UOV(I,1,L) + UOV(I+1,1,L)
      UOV(I,1,L) = TX1(I,1)*HALF
  170 CONTINUE
      DO 180 I=1,IJM2M0
      TX1(I,1)   = VOV(I,1,L) + VOV(I,2,L)
      VOV(I,1,L) = TX1(I,1)*HALF
  180 CONTINUE
C
      CALL WRAP(1,IJM0M0,UOV(1,1,L))
      CALL WRAP(1,IJM0M0,VOV(1,1,L))
      CALL FILPU(UOV(1,1,L))
      CALL FILPV(VOV(1,1,L))
  200 CONTINUE
C
      DO 400 L=1,LM
      IF (L .GT. 1) THEN
         DO 310 I=1,IJM0M0
         TX1(I,1) = TX5(I,1) * PRJ(L)
  310    CONTINUE
      ELSE
         DO 320 I=1,IJM0M0
         TX1(I,1) = TX5(I,1) * PRH(L)
  320    CONTINUE
      ENDIF
C
      DO 500 K=1,LR-1
      DO 410 I=1,IJM0M0
      IF (TX1(I,1).GT.PKSTR(K) .AND. TX1(I,1).LE.PKSTR(K+1)) THEN
         TX2(I,1)   = PKSTR(K+1)  - TX1(I,1)
         TX2(I,1)   = TX2(I,1)    / (PKSTR(K+1)-PKSTR(K))
         TX3(I,1)   = ONE         - TX2(I,1)
C
         TX4(I,1)   = GS(I,1,K)   * TX2(I,1)
         TOV(I,1,L) = GS(I,1,K+1) * TX3(I,1)
         TOV(I,1,L) = TOV(I,1,L)  + TX4(I,1)
      ENDIF
  410 CONTINUE
C
  500 CONTINUE
C
      DO 415 I=1,IJM0M0
      IF (TX1(I,1).GT.PKSTR(LR-1) .AND. GS(I,1,LR).LT.PHIX(I,1)) THEN
         TX2(I,1)   = TX5(I,1) * PRJ(LMP1)
         TX3(I,1)   = (TX2(I,1) - TX1(I,1)) / (TX2(I,1) - PKSTR(LR-1))
         TOV(I,1,L) = GS(I,1,LR-1) * TX3(I,1) + PHIX(I,1)*(1.0-TX3(I,1))
      ENDIF
  415 CONTINUE
      DO 420 I=1,IJM0M0
      IF (TX1(I,1) .GT. PKSTR(LR)) THEN
         TX2(I,1)   = PRJ(LMP1) * TX5(I,1) - TX1(I,1)
         TX4(I,1)   = PRJ(LMP1) * TX5(I,1) - PKSTR(LR)
         TX2(I,1)   = TX2(I,1) / TX4(I,1)
         TX3(I,1)   = ONE - TX2(I,1)
C
         TX4(I,1)   = GS(I,1,LR) * TX2(I,1)
         TOV(I,1,L) = PHIX(I,1)  * TX3(I,1)
         TOV(I,1,L) = TOV(I,1,L) + TX4(I,1)
      ENDIF
C
  420 CONTINUE

c     if (l .eq. lm) then
c        irt = i3+143
c        print *,' tov=',tov(irt,lm), ' tx1=',tx1(irt,1),' pkstr=',pkstr
c    *,   ' gs=',(gs(irt,1,ll),ll=1,lr)
c     endif

  400 CONTINUE
C
      DO 430 I=1,IJM0M0
      TOV(I,1,1) = TOV(I,1,1) - TOV(I,1,2)
      TOV(I,1,1) = TOV(I,1,1) / (RGAS*ALOG(SIG(2)/SGB(1)))
  430 CONTINUE
      IF (LM .GT. 2) THEN
         DO 600 L=2,LMM1
         DO 510 I=1,IJM0M0
         TOV(I,1,L) = TOV(I,1,L) - TOV(I,1,L+1)
         TOV(I,1,L) = TOV(I,1,L) / (RGAS * DSL(L))
  510    CONTINUE

  600    CONTINUE
      ENDIF
C
      DO 610 I=1,IJM0M0
      TOV(I,1,LM) = TOV(I,1,LM) - PHIX(I,1)
      TOV(I,1,LM) = TOV(I,1,LM) / (RGAS * DSL(LM))
  610 CONTINUE
c        print *,' tov=',tov(irt,lm), ' phix=',phix(irt,1),' rgas='
c    *,   rgas,' dsl=',dsl(LM)
C*****
C*****
      IF (LM .GT. 3) THEN
         DO 650 L=LMM2,2,-1
         DO 620 I=1,IJM0M0
         TX1(I,1) = TX5(I,1) * PRJ(L)
         TX2(I,1) = TX5(I,1) * PRH(L)
  620    CONTINUE
         TEM = ALOG(SGB(L+1)/SGB(L)) / ALOG(SGB(L+2)/SGB(L+1))
         DO 630 I=1,IJM0M0
         IF (TX1(I,1) .LT. PKSTR(1)) THEN
c           TEM1       = TOV(I,1,L+1) * PRH(L+1)
c           TEM2       = TEM1 - TOV(I,1,L+2) * PRH(L+2)
C
c           TOV(I,1,L) = (TEM1 + TEM2 * TEM) * (1.0/PRH(L))
            TX3(I,1)   = TOV(I,1,L+1) - TOV(I,1,L+2)
            TOV(I,1,L) = TOV(I,1,L+1) + TX3(I,1) * TEM
         ENDIF
  630    CONTINUE
         DO 640 I=1,IJM0M0
         IF (TX2(I,1) .LT. PKSTR(1)) THEN
            UOV(I,1,L) = UOV(I,1,L+1)
            VOV(I,1,L) = VOV(I,1,L+1)
         ENDIF
  640    CONTINUE
  650    CONTINUE
      ENDIF
      DO 660 I=1,IJM0M0
      TX2(I,1) = TX5(I,1) * PRH(1)
  660 CONTINUE
      TEM = ALOG(SGB(2)/SGB(1)) / ALOG(SGB(3)/SGB(2))
      DO 670 I=1,IJM0M0
      IF (TX2(I,1) .LT. PKSTR(1)) THEN
c        TEM1       = TOV(I,1,2) * PRH(2)
c        TEM2       = TEM1 - TOV(I,1,3) * PRH(3)
C
c        TOV(I,1,1) = (TEM1 + TEM2 * TEM) * (1.0/PRH(1))
c     IF (TX2(I,1) .LT. PKSTR(1)) THEN
         TX3(I,1)   = TOV(I,1,2) - TOV(I,1,3)
         TOV(I,1,1) = TOV(I,1,2) + TX3(I,1) * TEM
         UOV(I,1,1) = UOV(I,1,2)
         VOV(I,1,1) = VOV(I,1,2)
      ENDIF
  670 CONTINUE
C*****
C*****
   
      CALL SMTHX(TOV(1,JMM4,LM),TX1,IMP3,4)
      CALL SMTHX(TOV(1,JMM4,LM),TX1,IMP3,4)
      CALL SMTHX(TOV(1,JMM4,LM),TX1,IMP3,4)
      CALL SMTHX(TOV(1,   2,LM),TX1,IMP3,4)
      CALL SMTHX(TOV(1,   2,LM),TX1,IMP3,4)
      CALL SMTHX(TOV(1,   2,LM),TX1,IMP3,4)
      DO 700 L=1,LM
      DO 710 I=1,IJM0M0
      TX1(I,1) = TOV(I,1,L)
      TX2(I,1) = PIV(I,1) * SGB(L)
  710 CONTINUE
      CALL VQSAT(TX1,TX2,TX3,TX4,IJM0M0,.FALSE.)
c     print *,' qov=',(qov(I,2,L), I=1,IM,5)
      DO 720 I=1,IJM0M0
      QOV(I,1,L) = QOV(I,1,L) * TX3(I,1)
  720 CONTINUE
c     print *,' tx3=',(tx3(I,2), I=1,IM,5)
c     print *,' qov=',(qov(I,2,L), I=1,IM,5)
      call wrap(1,ijm0m0,qov(1,1,l))
C
      tem1 = 0.0
      tem2 = 0.0
      tem3 = 0.0
      tem4 = 0.0
      do 721 I=1,IM
      tem1 = tem1 + qov(i,2,l)
      tem3 = tem3 + tov(i,2,l)
      tem2 = tem2 + qov(i,jmm1,l)
      tem4 = tem4 + tov(i,jmm1,l)
  721 continue
      tem1 = tem1 / im
      tem2 = tem2 / im
      tem3 = tem3 / im
      tem4 = tem4 / im
      do 723 i=1,imp3
      qov(i,1,l)  = tem1
      tov(i,1,l)  = tem3
      qov(i,jm,l) = tem2
      tov(i,jm,l) = tem4
  723 continue
  700 CONTINUE
C
      DO 730 I=1,IJM0M0
      TSR(I,1) = TX5(I,1) * PRJ(LMP1)
      TSR(I,1) = TSR(I,1) * TOV(I,1,LM)
  730 CONTINUE
      tem1 = 0.0
      tem2 = 0.0
      do 731 I=1,IM
      tem1 = tem1 + tsr(i,2)
      tem2 = tem2 + tsr(i,jmm1)
  731 continue
      tem1 = tem1 / im
      tem2 = tem2 / im
      do 733 i=1,imp3
      tsr(i,1) = tem1
      tsr(i,jm) = tem2
  733 continue
C
c     print *,' tsrs=',(tsr(i,1),i=1,10)
c     print *,' tovs=',(tov(i,1),i=1,10)
c     print *,' tx5s=',(tx5(i,1),i=1,10)
c     print *,' tsrn=',(tsr(i,jm),i=1,10)
c     print *,' tovn=',(tov(i,jm),i=1,10)
c     print *,' tx5n=',(tx5(i,jm),i=1,10)
C
C
      DO 750 L=1,LM
      qmin = 1000.0
      qmax = -1000.0
      DO 745 I=1,IJM0M0
      TOV(I,1,L) = TOV(I,1,L) / (TX5(I,1)*PRH(L))
C     TOV(I,1,L) = TOV(I,1,L) / TBR(L) - 1.0
      qmin = amin1(qov(i,1,l),qmin)
      qmax = amax1(qov(i,1,l),qmax)
  745 CONTINUE
c     print *,' at level = ',l,' qmin=',qmin,' qmax=',qmax
  750 CONTINUE
C
      RETURN
      END
      SUBROUTINE INTER(IEM, JEM, IMP3, JM, FEC, FSL, WRK)
C
      DIMENSION FEC(IEM,JEM), FSL(IMP3,JM)
      DIMENSION WRK(IMP3,JEM)
C
      IM = IMP3 - 3
      JMM1 = JM - 1
C
      if (iem-1 .eq. im .and. jem .eq. jmm1) then
         do 5 j=2,jmm1
         do 5 i=1,im
         fsl(i,j) = fec(i,j-1)
    5    continue
         do 6 i=1,im
         wrk(i,1)     = fec(i,1)
         wrk(i,jem-1) = fec(i,jem-1)
    6    continue
      else
C
         DLON = 360.0 / FLOAT(IM)
         DLAT = 180.0 / FLOAT(JMM1)
         DLOE = FLOAT(IEM-1) / FLOAT(IM)
         DLAE = FLOAT(JEM)   / FLOAT(JMM1)
         FEC(1,:) = FEC(IEM,:)
C
C     FEC(IEM,:) = FEC(1,:)
         DO 20 J=1,JEM
         DO 10 I=1,IM
         XLON = (I-1) * DLOE + 1
         ILON = XLON
         XLON = XLON - ILON
         WRK(I,J) = (1.0-XLON)*FEC(ILON,J) + XLON*FEC(ILON+1,J)
   10    CONTINUE
   20    CONTINUE
C
         DO 50 J=2,JMM1
         YLAT =  (J-1) * DLAE
         JLAT = YLAT
         YLAT = YLAT - JLAT
         IF (JLAT .GT. 0 .AND. JLAT .LT. JEM-1) THEN
            DO 40 I=1,IM
            FSL(I,J) = (1.0-YLAT)*WRK(I,JLAT) + YLAT*WRK(I,JLAT+1)
   40       CONTINUE
         ELSEIF (JLAT .LE. 0) THEN
            DO 42 I=1,IM
            FSL(I,J) = WRK(I,1)
   42       CONTINUE
         ELSE
            DO 43 I=1,IM
            FSL(I,J) = WRK(I,JEM-1)
   43       CONTINUE
         ENDIF
C
   50    CONTINUE
      endif
C
      TEMS = 0.0
      TEMN = 0.0
      DO 60 I=1,IM
      TEMS = TEMS + WRK(I,1)
      TEMN = TEMN + WRK(I,JEM-1)
   60 CONTINUE
      TEMS = TEMS / FLOAT(IM)
      TEMN = TEMN / FLOAT(IM)
C
      DO 65 I=1,IM
      FSL(I,1)  = TEMS
      FSL(I,JM) = TEMN
   65 CONTINUE
C
      CALL WRAP (1,IMP3*JM,FSL(1,1))
C
      RETURN
      END
      SUBROUTINE INTRJ(JEM, JM, LM, FEC, FSL)
C
      DIMENSION FEC(JEM,LM), FSL(JM,LM)
C
      JMM1 = JM - 1
C
      if (jem .eq. jmm1) then
         do 5 j=2,jmm1
         do 5 l=1,lm
         fsl(j,l) = fec(j-1,l)
    5    continue
      else
C
         DLAT = 180.0 / FLOAT(JMM1)
         DLAE = FLOAT(JEM)   / FLOAT(JMM1)
C
         DO 50 J=2,JMM1
         YLAT =  (J-1) * DLAE
         JLAT = YLAT
         YLAT = YLAT - JLAT
         IF (JLAT .GT. 0 .AND. JLAT .LT. JEM-1) THEN
            DO 40 L=1,LM
            FSL(J,L) = (1.0-YLAT)*FEC(JLAT,L) + YLAT*FEC(JLAT+1,L)
   40       CONTINUE
         ELSEIF (JLAT .LE. 0) THEN
            DO 42 L=1,LM
            FSL(J,L) = FEC(1,L)
   42       CONTINUE
         ELSE
            DO 43 L=1,LM
            FSL(J,L) = FSL(JEM-1,L)
   43       CONTINUE
         ENDIF
C
   50    CONTINUE
      endif
C
      DO 65 L=1,LM
      FSL(1,L)  = FEC(1,L)
      FSL(JM,L) = FEC(JEM-1,L)
   65 CONTINUE
C
      RETURN
      END
      SUBROUTINE DCADJ(LEN, LENC, LM, CMP, CPM, AMP, APM
     *,                     POI, TX1, TX2, DT)
C
C
C*********************************************************************
C*********************** slgcm   MODEL *******************************
C********************* SUBROUTINE  DCADJ *****************************
C********************** 16 MARCH   1988 ******************************
C*********************************************************************
C
      PARAMETER (ZERO=0.0, ITRM=2)
C
      DIMENSION  POI(LEN,LM), TX1(LENC),   TX2(LENC)
      DIMENSION  CMP(LM),  CPM(LM), AMP(LM), APM(LM)
c     INTEGER DT
C
C     DRY CONVECTIVE ADJUSTMENT
C
      LMM1 = LM - 1
      DO 500 IRT=1,ITRM
      DO 30 L=1,LMM1
      DO 10 I=1,LENC
      TX1(I) = POI(I,L) * CMP(L) - POI(I,L+1)*CPM(L)
   10 CONTINUE
C
      DO 20 I=1,LENC
      IF (TX1(I) .LT. ZERO) THEN
         POI(I,L)   = POI(I,L)   - TX1(I)*AMP(L)
         POI(I,L+1) = POI(I,L+1) + TX1(I)*APM(L)
      ENDIF
   20 CONTINUE
C
   30 CONTINUE
  500 CONTINUE

      RETURN
      END
      SUBROUTINE DCADJ2(LEN, LENC, LM, CMP, CPM, AMP, APM
     *,                      QMP, QPM, POI, QOI, TX1, TX2, DT)
C
C
C*********************************************************************
C*********************** slgcm   MODEL *******************************
C********************* SUBROUTINE  DCADJ *****************************
C********************** 16 MARCH   1988 ******************************
C*********************************************************************
C
      PARAMETER (ZERO=0.0, ITRM=2)
C
      DIMENSION  POI(LEN,LM), QOI(LEN,LM), TX1(LENC),   TX2(LENC)
      DIMENSION  CMP(LM),  CPM(LM), AMP(LM), APM(LM)
     *,          QMP(LM),  QPM(LM)
c     INTEGER DT
C
C     DRY CONVECTIVE ADJUSTMENT
C
      LMM1 = LM - 1
      DO 500 IRT=1,ITRM
      DO 30 L=1,LMM1
      DO 10 I=1,LENC
      TX1(I) = POI(I,L) * CMP(L) - POI(I,L+1)*CPM(L)
      TX2(I) = QOI(I,L) - QOI(I,L+1)
   10 CONTINUE
C
      DO 20 I=1,LENC
      IF (TX1(I) .LT. ZERO) THEN
         POI(I,L)   = POI(I,L)   - TX1(I)*AMP(L)
         POI(I,L+1) = POI(I,L+1) + TX1(I)*APM(L)
         QOI(I,L)   = QOI(I,L)   - TX2(I)*QMP(L)
         QOI(I,L+1) = QOI(I,L+1) + TX2(I)*QPM(L)
      ENDIF
   20 CONTINUE
C
   30 CONTINUE
  500 CONTINUE

      RETURN
      END
      SUBROUTINE GETBCS(SEASNL)
C
#include "slgcm.slgprm"
#include "slgcm.slgcon"
#include "slgcm.slgprg"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE GETBCS *****************************
C********************** 29 January 1991 ******************************
C*********************************************************************
C
C
C    ARGUMENTS
C
      PARAMETER (JMLM13=JM*LM*13)
      PARAMETER (CAPFAC = (4.184E4)**2 * 86400.*.002/(2.*3.1416))
      PARAMETER (C1=.386, C2=.15)
      PARAMETER (ALBMIN=0.4, ALBMAX=1.0, ALBSNO=0.4)
      PARAMETER (SSTMIN=273.16, SSTMAX=310.0)

      LOGICAL SEASNL
C
C    LOCALS
C
      INTEGER     DAYS(13)
      DATA        DAYS/31,28,31,30,31,30,31,31,30,31,30,31,30/
C
      REAL  SSTM(IMP3,JM), SSTP(IMP3,JM), PHIM(IMP3,JM)
      REAL  WETM(IMP3,JM), WETP(IMP3,JM)
      REAL  ALBM(IMP3,JM), ALBP(IMP3,JM), GWDM(IMP3,JM)
      REAL  CAPM(IMP3,JM), PSGM(IMP3,JM)
      REAL  THKM(IMP3,JM), THKP(IMP3,JM)
      REAL  CDRM(IMP3,JM), CDRP(IMP3,JM)
      REAL  CDSM(IMP3,JM), CDSP(IMP3,JM), WK1(IMP3,JM)
      INTEGER LWI(IMP3,JM)
      character*3 mont
      character*4 yr
      character*4 nam1, nam2, nam3, nam4, nam5
      EQUIVALENCE (PHIS(1),PHIM(1,1))
C
      SAVE
      LOGICAL     FIRST, CHANGE
      DATA        FIRST/.TRUE./
      DATA K1/0/, K2/0/, MIDM/0/, MIDP/0/
C     DATA MONTH/01/, DAY/15/, HOUR/00/
C
      WRITE (6,988)
      WRITE (6,989) MONTH(), DAY(), HOUR()
C
  988 FORMAT(' MADE IT TO GETBCS ')
  989 FORMAT(' UPDATING BOUNDARY CONDITIONS FOR ', I3,I3,I3)
C
      IDAY   = DAY()
C     IDAY   = DAY
      IMON   = MONTH()
C     IMON   = MONTH
      MIDMON = DAYS(IMON)/2 + 1
      CHANGE = FIRST .OR.
c    *        ( (IDAY .EQ. MIDMON) .AND. (HOUR  .EQ.0) )
     *        ( (IDAY .EQ. MIDMON) .AND. (HOUR().EQ.0) )
C
      IF (CHANGE) THEN
         IF (SEASNL) THEN
            IF (IDAY .LT. MIDMON) THEN
               K1   = MOD(IMON+10,12) + 1
               MIDM = DAYS(K1)/2 + 1
               K2   = IMON
               MIDP = DAYS(K1) + MIDMON
               K1G  = K1 + 1
               K2G  = K2 + 1
            ELSE
               K1   = IMON
               MIDM = MIDMON
               K2   = MOD(IMON,12) + 1
               MIDP = DAYS(K2)/2 + 1 + DAYS(K1)
               K1G  = K1 + 1
               K2G  = K2 + 1
            ENDIF
         ELSE
            K1   = 13
            K2   = 13
            K1G  = K1
            K2G  = K2
            MIDM = 15
            MIDP = 45
         ENDIF
         REWIND KBG
         IF (FIRST) THEN
C
c           READ(KBG) ZIMJM
            READ(KBG) IMJMZ, yr
c           IMJMZ = ZIMJM + 0.00001
            IF (IMJMZ .NE. IM*JM) THEN
               PRINT *,' BOUNDARY CON. NOT COMPATIBLE'
               STOP
             ENDIF

             READ (KBG) nam1, ((PHIM(I,J),I=1,IM),J=1,JM)
             CALL WRAP(1,IJM0M0,PHIM)
C            PHIS(:) = 0.0
C            CALL MAP(PHIM,IMP3,JM,'PHIM',0.0)
C            CALL SMTH2D(PHIS,WK1,IMP3,JM)
C            CALL MAP(PHIS,IMP3,JM,'PHIS',0.0)
C
C            READ (KBG) CDR1
             READ (KBG) nam2, ((GWDM(I,J),I=1,IM),J=1,JM)
             READ (KBG) nam3, ((CAPM(I,J),I=1,IM),J=1,JM)
C
         ELSE
            READ(KBG) IMJMZ, yr
            READ (KBG) nam1, ((PHIM(I,J),I=1,IM),J=1,JM)
            READ (KBG) nam2, ((GWDM(I,J),I=1,IM),J=1,JM)
            READ (KBG) nam3, ((CAPM(I,J),I=1,IM),J=1,JM)
         ENDIF
         FIRST  = .FALSE.
         DO 100 K=1,K1G
         PRINT *,' IN 100 LOOP :, K=',K
         READ (KBG) mont, nam1, ((SSTM(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam2, ((ALBM(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam3, ((WETM(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam4, ((THKM(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam5, ((CDRM(I,J),I=1,IM),J=1,JM)
 100     CONTINUE
c        print *,' mont1=',mont
         REWIND KBG
         CALL WRAP(1,IJM0M0,SSTM)
         CALL WRAP(1,IJM0M0,ALBM)
         CALL WRAP(1,IJM0M0,WETM)
         CALL WRAP(1,IJM0M0,THKM)
         CALL WRAP(1,IJM0M0,CDRM)
C
c        print *,' in getbcs bef 20: cdrm=', (cdrm(i,1),i=1,10) 
         DO 20 I=1,IJM0M0
         CDSM(I,1) = PBLZET / CDRM(I,1) + ONE
         CDRM(I,1) = VKRM   / ALOG(CDSM(I,1))
         CDRM(I,1) = CDRM(I,1) * CDRM(I,1)
         CDSM(I,1) = SQRT(CDSM(I,1)) * CDRM(I,1) * BSF15
   20    CONTINUE
C
         READ(KBG) IMJMZ, yr
         READ (KBG) nam1, ((PHIM(I,J),I=1,IM),J=1,JM)
         READ (KBG) nam3, ((GWDM(I,J),I=1,IM),J=1,JM)
         READ (KBG) nam3, ((CAPM(I,J),I=1,IM),J=1,JM)
         CALL WRAP(1,IJM0M0,PHIM)
         CALL WRAP(1,IJM0M0,GWDM)
         CALL WRAP(1,IJM0M0,CAPM)
         DO 200 K=1,K2G
c        PRINT *,' IN 200 LOOP :, K=',K
         READ (KBG) mont, nam1, ((SSTP(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam2, ((ALBP(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam3, ((WETP(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam4, ((THKP(I,J),I=1,IM),J=1,JM)
         READ (KBG) mont, nam5, ((CDRP(I,J),I=1,IM),J=1,JM)
 200     CONTINUE
c        print *,' mont2=',mont
         REWIND KBG
         CALL WRAP(1,IJM0M0,SSTP)
         CALL WRAP(1,IJM0M0,ALBP)
         CALL WRAP(1,IJM0M0,WETP)
         CALL WRAP(1,IJM0M0,THKP)
         CALL WRAP(1,IJM0M0,CDRP)
C
C
         DO 30 I=1,IJM0M0
         CDSP(I,1) = PBLZET / CDRP(I,1) + 1.0
         CDRP(I,1) = VKRM   / ALOG(CDSP(I,1))
         CDRP(I,1) = CDRP(I,1) * CDRP(I,1)
         CDSP(I,1) = SQRT(CDSP(I,1)) * CDRP(I,1) * BSF15
   30    CONTINUE
C
      ENDIF
C
      IF (IDAY .LT. MIDMON) THEN
         ID = IDAY + DAYS(K1)
      ELSE
         ID = IDAY
      ENDIF
C
C     FAC = (REAL(ID  -MIDM)*DAYLEN + SECOFDY()) /
      FAC = (REAL(ID  -MIDM)*DAYLEN ) /
     *      (REAL(MIDP-MIDM)*DAYLEN            )
C
      DO 220 I=1,IJM0M0
      ALBIJ(I) = ALBP(I,1)*FAC + ALBM(I,1)*(ONE-FAC)
      WETIJ(I) = WETP(I,1)*FAC + WETM(I,1)*(ONE-FAC)
      THKIJ(I) = THKP(I,1)*FAC + THKM(I,1)*(ONE-FAC)
      SSTIJ(I) = SSTP(I,1)*FAC + SSTM(I,1)*(ONE-FAC)
      CDRIJ(I) = CDRP(I,1)*FAC + CDRM(I,1)*(1.0-FAC)
      CDSIJ(I) = CDSP(I,1)*FAC + CDSM(I,1)*(ONE-FAC)
C
      LWI(I,1) = CAPM(I,1) + EPSLN
      IF (LWI(I,1) .EQ. 1) THEN
         CAPIJ(I) = ZERO
         WETIJ(I) = 1.0
      ELSE
         CAPIJ(I) = ONE / SQRT((C1 + C2*WETIJ(I))*(ONE+WETIJ(I))*CAPFAC)
      ENDIF
  220 CONTINUE
C
      DO 260 I=1,IJM0M0
      IF (THKIJ(I).LT.THKMIN) THEN
         WK1(I,1) = THKMIN
      ELSE
         WK1(I,1) = THKIJ(I)
      ENDIF
  260 CONTINUE
      DO 280 I=1,IJM0M0
      IF(THKIJ(I).NE.ZERO) THEN
         CAPIJ(I) = CAPI
         THKIJ(I) = CICE / WK1(I,1)
      ENDIF
C     IF (CAPIJ(I,1) .EQ. 0.0 .AND. SSTIJ(I,1) .LT. 0.0)
C    *   PRINT *,' I=',I,' SST=',SSTIJ(I,1),' CAPIJ=',CAPIJ(I,1)
  280 CONTINUE
c
      DO 285 I=1,IJM0M0
      IF (LWI(I,1) .EQ. 3 .OR. THKIJ(I) .GT. 0.0) THEN
          IF (ALBIJ(I) .GT. ALBMAX) ALBIJ(I) = ALBMAX
          IF (ALBIJ(I) .LT. ALBMIN) ALBIJ(I) = ALBMIN
      ENDIF
c
      IF (LWI(I,1) .EQ. 1 .AND. THKIJ(I) .EQ. 0.0 ) THEN
         IF (SSTIJ(I) .LT. SSTMIN) SSTIJ(I) = SSTMIN
         IF (SSTIJ(I) .GT. SSTMAX) SSTIJ(I) = SSTMAX
      ENDIF
  285 CONTINUE
C
c     IF (.NOT. SURFR) THEN
         DO 222 I=1,IJM0M0
         CDRG(I) = 0.0
  222    CONTINUE
c     ENDIF
C
C     WRITE(6,221) (CDRG(I),I=1,IJM0M0)
C 221 FORMAT(' CDRG=',10E12.3)
C
      if (notopg) then
         PHIS(:) = 0.0
      endif
      CALL SMTHX(PHIS(I1+1),WK1,IMP3,3)
c     CALL SMTHX(PHIS(I1+1),WK1,IMP3,3)
C
      RETURN
      END
      SUBROUTINE WRTPBP(LEN, LGP, DIV, ISEL, NPT, PBR, NSTEP, PSAV,DSAV)
C
      REAL LGP(LEN), DIV(LEN), PSAV(NPT), DSAV(NPT)
      DIMENSION ISEL(NPT)
C
      DO 10 N=1,NPT
      PSAV(N)  = EXP(LGP(ISEL(N)) + PBR) * 0.01
      DSAV(N) = DIV(ISEL(N))
   10 CONTINUE
C
      WRITE(6,20) (PSAV(N),N=9,NPT), (DSAV(N), N=9,NPT), NSTEP
   20 FORMAT(' PS=',5F8.2,' DIV=',5E12.4,' NST=',I5)
C
      RETURN
      END
      SUBROUTINE DELSQR (CPHU, CPHV, POI, DSQ)
C    *,                WK1, WK2, WK3, WK4, WK5, WK6, WK7, WK8, WK9)
C
#include "slgcm.slgprm"
C
      DIMENSION POI(IMP3,JM), DSQ(IMP3,JM)
      DIMENSION CPHU(IMP3,JM),   CPHV(IMP3,JM)
C
      DIMENSION WK1(IMP3,JM), WK2(IMP3,JM), WK3(IMP1,JM)
C
      SAVE

      LOGICAL FIRST
      DATA FIRST /.TRUE./
C
      IF (FIRST) THEN
         DIFAC  = 5.0E3
         DIF2   = DIFAC * DY  / DXP
         DIF3   = DIFAC * DXP / DY
         FIRST = .FALSE.
C
         PRINT *,' IN DELSQR: DIFAC=',DIFAC,' DIF2=',DIF2,' DIF3=',DIF3
C
      ENDIF
C
C     DIFFUSION OF POTENTIAL TEMPERATURE AT TOP & BOTTOM LEVELS
C
C
      DO 130 I=1,IJM0M1
      WK1(I,1) =  POI(I+1,1) - POI(I,1)
  130 CONTINUE
      CALL WRAP(1,IJM0M0,WK1)
C
      DO 140 I=I1P1,IJM1M0
      WK1(I,1) = DIF2 * WK1(I,1) / CPHU(I,1)
  140 CONTINUE

      DO 150 I=I1P1,IJM1M1
      WK2(I+1,1) = WK1(I+1,1) - WK1(I,1)
  150 CONTINUE
      CALL WRAP(I1P1,IJM1M0,WK2)
C
      DO 160 I=1,IJM0M0
      WK1(I,1) = 0.0
      DSQ(I,1) = 0.0
  160 CONTINUE
      DO 170 I=1,IJM1M0
      WK3(I,1) =  POI(I,2) - POI(I,1)
C
      WK1(I,1) = DIF3 *  WK3(I,1) * CPHV(I,1)
  170 CONTINUE
      DO 180 I=1,IJM2M0
      WK2(I,2) = WK2(I,2) + WK1(I,2) - WK1(I,1)
  180 CONTINUE
C
      DO 260 I=I1P1,IJM1M0
      DSQ(I,1)  = WK2(I,1) / ((DXP*DY) * CPHU(I,1))
  260 CONTINUE
C
      RETURN
      END
      SUBROUTINE DELSQV (CPHU, CPHV, VOI, DSQ)
C    *,                WK1, WK2, WK3, WK4, WK5, WK6, WK7, WK8, WK9)
C
#include "slgcm.slgprm"
C
      DIMENSION VOI(IMP3,JM), DSQ(IMP3,JM)
      DIMENSION CPHU(IMP3,JM),   CPHV(IMP3,JM)
C
      DIMENSION WK1(IMP3,JM), WK2(IMP3,JM), WK3(IMP1,JM)
C
      SAVE

      LOGICAL FIRST
      DATA FIRST /.TRUE./
C
      IF (FIRST) THEN
         DIFAC  = 5.0E3
         DIF2   = DIFAC * DY  / DXP
         DIF3   = DIFAC * DXP / DY
         FIRST = .FALSE.
C
         PRINT *,' IN DELSQR: DIFAC=',DIFAC,' DIF2=',DIF2,' DIF3=',DIF3
C
      ENDIF
C
C     DIFFUSION OF V
C
C
      DO 130 I=1,IJM1M1
      WK1(I,1) =  VOI(I+1,1) - VOI(I,1)
  130 CONTINUE
      CALL WRAP(1,IJM1M0,WK1)
C
      DO 140 I=1,IJM1M0
      WK1(I,1) = DIF2 * WK1(I,1) / CPHV(I,1)
  140 CONTINUE

      DO 150 I=1,IJM1M1
      WK2(I+1,1) = WK1(I+1,1) - WK1(I,1)
  150 CONTINUE
      CALL WRAP(1,IJM1M0,WK2)
C
      DO 160 I=1,IJM0M0
      WK1(I,1) = 0.0
      DSQ(I,1) = 0.0
  160 CONTINUE
      DO 170 I=1,IJM2M0
      WK3(I,2) =  VOI(I,2) - VOI(I,1)
C
      WK1(I,2) = DIF3 *  WK3(I,2) * CPHU(I,2)
  170 CONTINUE
      DO 180 I=1,IJM1M0
      WK2(I,1) = WK2(I,1) + WK1(I,2) - WK1(I,1)
  180 CONTINUE
C
      DO 260 I=1,IJM1M0
      DSQ(I,1)  = WK2(I,1) / ((DXP*DY) * CPHV(I,1))
  260 CONTINUE
C
      RETURN
      END
      SUBROUTINE MAP(A,IM,JM,NAME,UND)
      DIMENSION A(IM,JM)
      CHARACTER*1 AI(IM,JM), CS(12)
      CHARACTER*(*) NAME
      DATA CS/'-','2','3','4','5','6','7','8','9','+','+','U'/

      AMAX = -1.0E20
      AMIN = -AMAX
      DO 100 I=1,IM*JM
      IF (A(I,1) .NE. UND) THEN
         IF (A(I,1) .GT. AMAX) AMAX = A(I,1)
         IF (A(I,1) .LT. AMIN) AMIN = A(I,1)
      ENDIF
100   CONTINUE
C
      DO 200 I=1,IM*JM
      IF (A(I,1) .NE. UND) THEN
         IX = 1 + 10.*(A(I,1)-AMIN)/(AMAX-AMIN)
      ELSE
         IX = 12
      ENDIF
      AI(I,1) = CS(IX)
200   CONTINUE

      PRINT 40
      PRINT 40
      PRINT *, NAME, '  AMAX = ', AMAX, '  AMIN = ', AMIN
      PRINT 40
      DO 300 J=JM,1,-1
      PRINT 30, (AI(I,J),I=1,IM)
300   CONTINUE

30    FORMAT(2X,130A1)
40    FORMAT(//)


      RETURN
      END

      SUBROUTINE SMTH2D(A,TEMP,IMP1,JMP1)
      REAL A(IMP1,JMP1), TEMP(IMP1,JMP1)
C
      IM   = IMP1 - 1
      JM   = JMP1 - 1
      JMM1 = JM   - 1
C
      DO 40 J=2,JM
         DO 30 I=2,IM
            TEMP(I,J) = .25*(A(I-1,J)+A(I+1,J)+2.*A(I,J))
   30    CONTINUE
         TEMP(   1,J) = A(   1,J)
         TEMP(IMP1,J) = A(IMP1,J)
   40 CONTINUE
      DO 60 I=1,IMP1
         DO 50 J=3,JMM1
            A(I,J) = .25*(TEMP(I,J-1)+TEMP(I,J+1)+2.*TEMP(I,J))
   50    CONTINUE
         A(I, 2) = TEMP(I, 2)
         A(I,JM) = TEMP(I,JM)
   60 CONTINUE
C
      RETURN
      END
      SUBROUTINE WRAP(IB,IE,XX)
#include "slgcm.slgprm"
      DIMENSION XX(IJM0M0)
C
      DO 10 I=IB,IE,IMP3
      XX(IM+I)   = XX(I)
      XX(IMP1+I) = XX(I+1)
      XX(IMP2+I) = XX(I+2)
   10 CONTINUE
C
      RETURN
      END
      SUBROUTINE CONDT(DT)
C
#include "slgcm.slgprm"
#include "slgcm.slgclm"
#include "slgcm.slgwrk"
#include "slgcm.slgcon"
C
      INTEGER DT
C
C     SAVE
C
      CDA = CDF * CDAMP
      CDB = (2.0 - CDF) * CDAMP
C
      PRINT *,' CDAMP=',CDAMP, ' CDA=',CDA,' CDB=',CDB
C
      DTOA   = DT/AE
      DTHF   = DT*0.5
      DTHFSQ = DTHF*DTHF
      DL1DT  = DL1 * DT
      DL2DT  = DL2 * DT
      HFCDDT = CDA*DL1DT
C
C PRECOMPUTE ARRAYS AT V POINTS
C
        DO 1226 J=1,JM
C
        FCOR(J,1)  = TWOM * SPH(J,1)
        FDTHF(J,1) = FCOR(J,1) * DL1DT
        FVSQ(J)    = FDTHF(J,1) * FDTHF(J,1)
        ARR1(J,1)  = ONE / (ONE + FVSQ(J))
        ARR2(J,1)  = FDTHF(J,1) / (ONE + FVSQ(J))
        ARR5(J)    = DL1DT * ARR2(J,1) * RH1JK(J,1)
        ARR6(J)    = DL1DT * ARR1(J,1) * RA
        F1M(J,1)   =(ONE - FVSQ(J))   * ARR1(J,1)
C
C PRECOMPUTE ARRAYS AT U POINTS
C
        FCOR(J,2)  = TWOM*SPH(J,2)
        FDTHF(J,2) = FCOR(J,2)  * DL1DT
        FUSQ(J)    = FDTHF(J,2) * FDTHF(J,2)
        ARR1(J,2)  = ONE / (ONE + FUSQ(J))
        ARR2(J,2)  = FDTHF(J,2) / (ONE + FUSQ(J))
        ARR3(J)    = DL1DT * ARR1(J,2) * RH1JK(J,2)
        ARR4(J)    = DL1DT * ARR2(J,2) * RA
        F1M(J,2)   = (ONE - FUSQ(J))  * ARR1(J,2)
C
 1226 CONTINUE
C
C
      TEM  = DPH * DPH * DL1DT * TWOM * 0.5 / (DLM*DY*DY)
      AESQ = AE * AE
      DO 240 J=2,JMM1
      TEM1 = 1.0 / (DX(J,2) * DX(J,2))
      TEM2 = TEM * (1.0 - FUSQ(J)) * ARR1(J,2)
      TEM5 = ONE / (ARR1(J,2) * DY * DY *CPH(J,2))
      TEM3 = ARR1(J,1)   * CPH(J,1)   * TEM5
      TEM4 = ARR1(J-1,1) * CPH(J-1,1) * TEM5
      II = (J-1) * IMP3
      DO 235 I=1,IMP3
      CDW(II+I) = TEM1 + TEM2
      CDE(II+I) = TEM1 - TEM2
      CDN(II+I) = TEM3
      CDS(II+I) = TEM4
C
  235 CONTINUE
C
      CDN2(J) = ARR1(J,1) * CPH(J,1) / (CPH(J,2) *  DY*DY)
      CDS2(J) = ARR1(J-1,1) * CPH(J-1,1) / (CPH(J,2) * DY*DY)
      CDC1(J) = ARR1(J,2) / (CPH(J,2) * CPH(J,2) *AESQ)
      CDC2(J) = ARR1(J,2)*ARR1(J,2) *(1.0-FUSQ(J))*(DL1DT*TWOM) / AESQ
C
  240 CONTINUE
      CDC1(1)  = 4.0 * ARR1(1,1) / (DY*DY)
      CDC1(JM) = 4.0 * ARR1(JMM1,1) / (DY*DY)
C
c     PRINT *,' ARR1=',ARR1
C     PRINT *,' CDW=',CDW
C     PRINT *,' CDE=',CDE
C     PRINT *,' CDN=',CDN
C     PRINT *,' CDS=',CDS
c     PRINT *,' CDN2=',CDN2
c     PRINT *,' CDS2=',CDS2
c     PRINT *,' CDC1=',CDC1
c     PRINT *,' CDC2=',CDC2
C
C
      TEM = 0.0
      DO 503 L=LM,1,-1
      TEM1    = DSL(L)*TBR(L) * RGAS
      GBRM(L) = TEM + TEM1*0.5
      TEM     = TEM + TEM1
  503 CONTINUE
C
C
      DO 600 L=1,LM
      TEM = DT*0.5
C     HFDTLK(L) = TEM * ELM(L)
C     SQMU(L)   = 1.0 / (TEM*(HFDTLK(L) + CDA))
      HFDTLK(L) = DL1DT * ELM(L)
      SQMU(L)   = 1.0 / (DL1DT*(HFDTLK(L) + CDA))
  600 CONTINUE
C
      PRINT *,' SQMU=',SQMU
C
      RETURN
      END
      SUBROUTINE SMTHX(A,TX1,IMP3,LM)
      REAL A(IMP3,LM), TX1(IMP3)
C
      IMP2 = IMP3 - 1
C
      DO 40 L=1,LM
      call wrap(1,imp3,A(1,L))
      DO 30 I=2,IMP2
      TX1(I) = .25*(A(I-1,L)+A(I+1,L)+2.*A(I,L))
   30 CONTINUE
      DO 35 I=2,IMP2
      A(I,L) = TX1(I)
   35 CONTINUE
      A(1,L)    = A(IMP2-1,L)
      A(IMP3,L) = A(3,L)
   40 CONTINUE
C
      RETURN
      END
      SUBROUTINE AVRX(U,V)
#include "slgcm.slgprm"
#include "slgcm.slgwrk"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE AVRX *******************************
C********************** 08 OCTO    1992 ******************************
C*********************************************************************
C
C
      PARAMETER (NLATS=JM/8)
C
      INTEGER   TOTLATS, TEMPSIZE, FORWARD, BACKWARD
C
      PARAMETER (TOTLATS=4*NLATS*LM)
      PARAMETER (TEMPSIZE=TOTLATS*IMP2, ITRG=3*IMP2+1)
      PARAMETER (FORWARD=-1,BACKWARD=1)
C
C  ARGUMENTS
C
      DIMENSION    U(IMP3,JM,LM), V(IMP3,JM,LM)
     *,            S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
C   LOCALS
C
C     COMMON /AVRX1/ S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
      LOGICAL        FIRST
C
C   SCRATCH SPACE ALLOCATION
C
      DIMENSION      A(IMP2,TOTLATS), B(IM*2,TOTLATS)
      DIMENSION      F(IMP2,JMM1),    Z(IMP3,TOTLATS)
C
      DATA           FIRST/.TRUE./
      DATA S/TEMPSIZE*0.0/, TR/ITRG*0.0/, IX/IM*0/
      SAVE
C
C   INITIALIZATION ON FIRST CALL
C
      IF (FIRST) THEN
         FIRST = .FALSE.
         print *,' NLATS = ', NLATS
C
         DO 10 I=1,IMP2*JMM1
         F(I,1) = 1.
   10    CONTINUE
         CONS   = ONE / COS(-PI/TWO+(NLATS+1)*(PI/JMM1))
C
         DO 100 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 100 J=2,JMM1
         TEM      = CONS * COS(-PI/TWO + J*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
100      CONTINUE
C
         DO 200 I=1,IMP2
         DO 200 J=1,NLATS
         S(I,J      ) = F(I,         1 +J)
         S(I,J+NLATS) = F(I,JM-NLATS-1 +J)
200      CONTINUE
C
         DO 300 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 300 J=1,JMM1
         TEM    = CONS * COS(-PI/TWO+(J-HALF)*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
300      CONTINUE
C
         DO 400 I=1,IMP2
         DO 400 J=1,NLATS
         S(I,J+2*NLATS) = F(I,            J)
         S(I,J+3*NLATS) = F(I,JM-NLATS-1 +J)
400      CONTINUE
         DO 405 L=1,LMM1
         DO 407 IQ=1,4*NLATS
         DO 409 I=1,IMP2
         S(I,4*NLATS*L+IQ) = S(I,4*NLATS*(L-1)+IQ)
409      CONTINUE
407      CONTINUE
405      CONTINUE
C
         CALL FFTFAX(IM, IX, TR)
C
      ENDIF
C
C       END INITIALIZATION
C**********************************************************************
C       BEGIN STANDARD CASE
C
C   COPY LATITUDES TO BE SMOOTHED INTO CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1005 L=1,LM
      DO 1003 I=1,IMP3*NLATS
      Z(I,IOFF)         = U(I,2,L)
      Z(I,NLATS+IOFF)   = U(I,JM-NLATS,L)
      Z(I,2*NLATS+IOFF) = V(I,1,L)
      Z(I,3*NLATS+IOFF) = V(I,JM-NLATS,L)
1003  CONTINUE
      IOFF = IOFF + 4*NLATS
1005  CONTINUE
C
C   TRANSFORM AND FILTER
C
      DO 1100 J=1,TOTLATS
      DO 1100 I=1,IMP2
      A(I,J) = Z(I,J)
1100  CONTINUE
C
c     CALL FFT991(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
C
      DO 1110 I=1,IMP2*TOTLATS
      A(I,1) = A(I,1)*S(I,1)
 1110 CONTINUE
C
c     CALL FFT991(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
C
      DO 1500 J=1,TOTLATS
      DO 1450 I=1,IM
      Z(I,J) = A(I,J)
1450  CONTINUE
      Z(IMP1,J) = Z(1,J)
      Z(IMP2,J) = Z(2,J)
      Z(IMP3,J) = Z(3,J)
1500  CONTINUE
C
C   COPY BACK FROM CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1650 L=1,LM
      DO 1645 I=1,IMP3*NLATS
      U(I,2,L)        = Z(I,IOFF)
      U(I,JM-NLATS,L) = Z(I,NLATS+IOFF)
      V(I,1,L)        = Z(I,2*NLATS+IOFF)
      V(I,JM-NLATS,L) = Z(I,3*NLATS+IOFF)
1645  CONTINUE
      IOFF = IOFF + 4*NLATS
1650  CONTINUE
C
      RETURN
      END
      SUBROUTINE AVRX1(U,V)
#include "slgcm.slgprm"
#include "slgcm.slgwrk"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE AVRX *******************************
C********************** 08 OCTO    1992 ******************************
C*********************************************************************
C
C
      PARAMETER (NLATS=JM/8)
C
      INTEGER   TOTLATS, TEMPSIZE, FORWARD, BACKWARD
C
      PARAMETER (TOTLATS=4*NLATS)
      PARAMETER (TEMPSIZE=TOTLATS*IMP2, ITRG=3*IMP2+1)
      PARAMETER (FORWARD=-1,BACKWARD=1)
C
C  ARGUMENTS
C
      DIMENSION    U(IMP3,JM), V(IMP3,JM)
     *,            S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
C   LOCALS
C
C     COMMON /AVRX1/ S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
      LOGICAL        FIRST
C
C   SCRATCH SPACE ALLOCATION
C
      DIMENSION      A(IMP2,TOTLATS), B(IM*2,TOTLATS)
      DIMENSION      F(IMP2,JMM1),    Z(IMP3,TOTLATS)
C
      DATA           FIRST/.TRUE./
      DATA S/TEMPSIZE*0.0/, TR/ITRG*0.0/, IX/IM*0/
      SAVE
C
C   INITIALIZATION ON FIRST CALL
C
      IF (FIRST) THEN
         FIRST = .FALSE.
C
         DO 10 I=1,IMP2*JMM1
         F(I,1) = 1.
   10    CONTINUE
         CONS   = ONE / COS(-PI/TWO+(NLATS+1)*(PI/JMM1))
C
         DO 100 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 100 J=2,JMM1
         TEM      = CONS * COS(-PI/TWO + J*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
100      CONTINUE
C
         DO 200 I=1,IMP2
         DO 200 J=1,NLATS
         S(I,J      ) = F(I,         1 +J)
         S(I,J+NLATS) = F(I,JM-NLATS-1 +J)
200      CONTINUE
C
         DO 300 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 300 J=1,JMM1
         TEM    = CONS * COS(-PI/TWO+(J-HALF)*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
300      CONTINUE
C
         DO 400 I=1,IMP2
         DO 400 J=1,NLATS
         S(I,J+2*NLATS) = F(I,            J)
         S(I,J+3*NLATS) = F(I,JM-NLATS-1 +J)
400      CONTINUE
C
         CALL FFTFAX(IM, IX, TR)
C
      ENDIF
C
C       END INITIALIZATION
C**********************************************************************
C       BEGIN STANDARD CASE
C
C   COPY LATITUDES TO BE SMOOTHED INTO CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1003 I=1,IMP3*NLATS
      Z(I,IOFF)         = U(I,2)
      Z(I,NLATS+IOFF)   = U(I,JM-NLATS)
      Z(I,2*NLATS+IOFF) = V(I,1)
      Z(I,3*NLATS+IOFF) = V(I,JM-NLATS)
1003  CONTINUE
C
C   TRANSFORM AND FILTER
C
      DO 1100 J=1,TOTLATS
      DO 1100 I=1,IMP2
      A(I,J) = Z(I,J)
1100  CONTINUE
C
c     CALL FFT991(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
C
      DO 1110 I=1,IMP2*TOTLATS
      A(I,1) = A(I,1)*S(I,1)
 1110 CONTINUE
C
c     CALL FFT991(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
C
      DO 1500 J=1,TOTLATS
      DO 1450 I=1,IM
      Z(I,J) = A(I,J)
1450  CONTINUE
      Z(IMP1,J) = Z(1,J)
      Z(IMP2,J) = Z(2,J)
      Z(IMP3,J) = Z(3,J)
1500  CONTINUE
C
C   COPY BACK FROM CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1645 I=1,IMP3*NLATS
      U(I,2)        = Z(I,IOFF)
      U(I,JM-NLATS) = Z(I,NLATS+IOFF)
      V(I,1)        = Z(I,2*NLATS+IOFF)
      V(I,JM-NLATS) = Z(I,3*NLATS+IOFF)
1645  CONTINUE
C
      RETURN
      END
      SUBROUTINE AVRX2(U)
#include "slgcm.slgprm"
#include "slgcm.slgwrk"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE AVRX *******************************
C********************** 08 OCTO    1992 ******************************
C*********************************************************************
C
C
      PARAMETER (NLATS=JM/8)
C
      INTEGER   TOTLATS, TEMPSIZE, FORWARD, BACKWARD
C
      PARAMETER (TOTLATS=2*NLATS)
      PARAMETER (TEMPSIZE=TOTLATS*IMP2, ITRG=3*IMP2+1)
      PARAMETER (FORWARD=-1,BACKWARD=1)
C
C  ARGUMENTS
C
      DIMENSION    U(IMP3,JM)
     *,            S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
C   LOCALS
C
C     COMMON /AVRX1/ S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
      LOGICAL        FIRST
C
C   SCRATCH SPACE ALLOCATION
C
      DIMENSION      A(IMP2,TOTLATS), B(IM*2,TOTLATS)
      DIMENSION      F(IMP2,JMM1),    Z(IMP3,TOTLATS)
C
      DATA           FIRST/.TRUE./
      DATA S/TEMPSIZE*0.0/, TR/ITRG*0.0/, IX/IM*0/
      SAVE
C
C   INITIALIZATION ON FIRST CALL
C
      IF (FIRST) THEN
         FIRST = .FALSE.
C
         DO 10 I=1,IMP2*JMM1
         F(I,1) = 1.
   10    CONTINUE
         CONS   = ONE / COS(-PI/TWO+(NLATS+1)*(PI/JMM1))
C
         DO 100 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 100 J=2,JMM1
         TEM      = CONS * COS(-PI/TWO + J*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
100      CONTINUE
C
         DO 200 I=1,IMP2
         DO 200 J=1,NLATS
         S(I,J      ) = F(I,         1 +J)
         S(I,J+NLATS) = F(I,JM-NLATS-1 +J)
200      CONTINUE
C
         CALL FFTFAX(IM, IX, TR)
C
      ENDIF
C
C       END INITIALIZATION
C**********************************************************************
C       BEGIN STANDARD CASE
C
C   COPY LATITUDES TO BE SMOOTHED INTO CONTIGUOUS BUFFER
C
      DO 1003 I=1,IMP3*NLATS
      Z(I,1)         = U(I,2)
      Z(I,NLATS+1)   = U(I,JM-NLATS)
1003  CONTINUE
C
C   TRANSFORM AND FILTER
C
      DO 1100 J=1,TOTLATS
      DO 1100 I=1,IMP2
      A(I,J) = Z(I,J)
1100  CONTINUE
C
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
C
      DO 1110 I=1,IMP2*TOTLATS
      A(I,1) = A(I,1)*S(I,1)
 1110 CONTINUE
C
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
C
      DO 1500 J=1,TOTLATS
      DO 1450 I=1,IM
      Z(I,J) = A(I,J)
1450  CONTINUE
      Z(IMP1,J) = Z(1,J)
      Z(IMP2,J) = Z(2,J)
      Z(IMP3,J) = Z(3,J)
1500  CONTINUE
C
C   COPY BACK FROM CONTIGUOUS BUFFER
C
      DO 1645 I=1,IMP3*NLATS
      U(I,2)        = Z(I,1)
      U(I,JM-NLATS) = Z(I,NLATS+1)
1645  CONTINUE
C
      RETURN
      END
      SUBROUTINE AVRX3(U,V)
#include "slgcm.slgprm"
#include "slgcm.slgwrk"
C
C*********************************************************************
C*********************** SLGCM   MODEL *******************************
C********************* SUBROUTINE AVRX *******************************
C********************** 08 OCTO    1992 ******************************
C*********************************************************************
C
C
      PARAMETER (NLATS=JM/8)
C
      INTEGER   TOTLATS, TEMPSIZE, FORWARD, BACKWARD
C
      PARAMETER (TOTLATS=4*NLATS)
      PARAMETER (TEMPSIZE=TOTLATS*IMP2, ITRG=3*IMP2+1)
      PARAMETER (FORWARD=-1,BACKWARD=1)
C
C  ARGUMENTS
C
      DIMENSION    U(IMP3,JM), V(IMP3,JM)
     *,            S(IMP2,TOTLATS), TR(ITRG), IX(IM)
C
C   LOCALS
C
      LOGICAL        FIRST
C
C   SCRATCH SPACE ALLOCATION
C
      DIMENSION      A(IMP2,TOTLATS), B(IM*2,TOTLATS)
      DIMENSION      F(IMP2,JMM1),    Z(IMP3,TOTLATS)
C
      DATA           FIRST/.TRUE./
      DATA S/TEMPSIZE*0.0/, TR/ITRG*0.0/, IX/IM*0/
      SAVE
C
C   INITIALIZATION ON FIRST CALL
C
      IF (FIRST) THEN
         FIRST = .FALSE.
C
         DO 10 I=1,IMP2*JMM1
         F(I,1) = 1.
   10    CONTINUE
         CONS   = ONE / COS(-PI/TWO+(NLATS+1)*(PI/JMM1))
C
         DO 100 I=4,IM
         ANGI = SIN((I/2)*PI/IM)
         DO 100 J=2,JMM1
         TEM      = CONS * COS(-PI/TWO + (J-1)*(PI/JMM1))
         F(I+1,J) = MIN(TEM*TEM/ANGI,ONE)
100      CONTINUE
C
         DO 200 I=1,IMP2
         DO 200 J=1,NLATS
         S(I,J      ) = F(I,         1 +J)
         S(I,J+NLATS) = F(I,JM-NLATS-1 +J)
         S(I,J+2*NLATS) = S(I,J)
         S(I,J+3*NLATS) = S(I,J+NLATS)
200      CONTINUE
C
         CALL FFTFAX(IM, IX, TR)
C
      ENDIF
C
C       END INITIALIZATION
C**********************************************************************
C       BEGIN STANDARD CASE
C
C   COPY LATITUDES TO BE SMOOTHED INTO CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1003 I=1,IMP3*NLATS
      Z(I,IOFF)         = U(I,2)
      Z(I,NLATS+IOFF)   = U(I,JM-NLATS)
      Z(I,2*NLATS+IOFF) = V(I,2)
      Z(I,3*NLATS+IOFF) = V(I,JM-NLATS)
1003  CONTINUE
C
C   TRANSFORM AND FILTER
C
      DO 1100 J=1,TOTLATS
      DO 1100 I=1,IMP2
      A(I,J) = Z(I,J)
1100  CONTINUE
C
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,FORWARD)
C
      DO 1110 I=1,IMP2*TOTLATS
      A(I,1) = A(I,1)*S(I,1)
 1110 CONTINUE
C
      CALL RFFTMLT(A,B,TR,IX,1,IMP2,IM,TOTLATS,BACKWARD)
C
      DO 1500 J=1,TOTLATS
      DO 1450 I=1,IM
      Z(I,J) = A(I,J)
1450  CONTINUE
      Z(IMP1,J) = Z(1,J)
      Z(IMP2,J) = Z(2,J)
      Z(IMP3,J) = Z(3,J)
1500  CONTINUE
C
C   COPY BACK FROM CONTIGUOUS BUFFER
C
      IOFF = 1
      DO 1645 I=1,IMP3*NLATS
      U(I,2)        = Z(I,IOFF)
      U(I,JM-NLATS) = Z(I,NLATS+IOFF)
      V(I,2)        = Z(I,2*NLATS+IOFF)
      V(I,JM-NLATS) = Z(I,3*NLATS+IOFF)
1645  CONTINUE
C
      RETURN
      END

