      SUBROUTINE PTETAE(PT,KBIMX,LMIMX,NBCEX
     1,                 IMI,JMI,LMI0
     2,                 DETAO,IUDETI
     3,                 PDBO,TBO,QBO,UBO,VBO,Q2BO,CWMBO
     4,                 HTMBO,USTARBO
     5,                 PDBI,TBI,QBI,UBI,VBI,Q2BI,CWMBI)
C
C$$$  SUBPROGRAM DOCUMENTATION BLOCK
C                .      .    .
C   SUBROUTINE:  PTETAE      ETA-TO-ETA INTERPOLATION PROGRAM
C   PRGRMMR: BLACK           ORG: W/NP22     DATE: 99-05-25
C
C  ABSTRACT:  VERTICALLY INTERPOLATES ETA DATA AT ONE VERTICAL
C             RESOLUTION TO ETA LAYERS AT ANOTHER RESOLUTION.
C             PROGRAM ASSUMES HORIZONTAL INTERPOLATION HAS BEEN
C             PERFORMED EARLIER AND THAT THE TOP PRESSURE OF THE
C             INITIAL AND FINAL VERTICAL GRIDS IS THE SAME
C
C PROGRAM HISTORY LOG:
C    95-07-14  T BLACK - WROTE CODE FOR ETA RESTRT-->ETA RESTRT ABILITY
C    96-10-29  T BLACK - MODIFIED FOR GENERAL NESTING CAPABILITY
C    99-05-25  T BLACK - MODIFIED FOR BACKGROUND PRE-POST JOB
C
C USAGE:  CALL PTETAE FROM SUBROUTINE BCEX
C
C   INPUT ARGUMENT LIST:
C         PT - THE TOP PRESSURE OF THE DOMAIN
C      KBIMX - THE MAXIMUM LENGTH OF ANY NEST BOUNDARY'S HORIZONTAL
C              EXTENT
C      LMIMX - THE MAXIMUM NUMBER OF MODEL LEVELS FOR ANY NEST
C      NBCEX - THE NUMBER OF NESTS
C        IMI - THE ARRAY OF IM's FOR ALL NESTS
C        JMI - THE ARRAY OF JM's FOR ALL NESTS
C       LMI0 - THE ARRAY OF LM's FOR ALL NESTS
C     IUDETI - THE BASE UNIT NUMBER OF THE NEST DETA's
C       PDBO - THE PD ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C        TBO - THE TEMPERATURE ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C        QBO - THE SPEC HUMIDITY ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C        UBO - THE U COMPONENT ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C        VBO - THE V COMPONENT ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C       Q2BO - THE TKE ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C      CWMBO - THE CLOUD WATER ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C      HTMBO - THE HEIGHT MASK ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C    USTARBO - THE USTAR ARRAY FOR ALL NEST BOUNDARIES AFTER
C              HORIZONTAL INTERPOLATION FROM THE OPERATIONAL GRID
C
C   OUTPUT ARGUMENT LIST:
C       PDBI - THE FINAL PD ARRAY FOR ALL NEST BOUNDARIES
C        TBI - THE FINAL TEMP ARRAY FOR ALL NEST BOUNDARIES
C        QBI - THE FINAL SPEC HUMIDITY ARRAY FOR ALL NEST BOUNDARIES
C        UBI - THE FINAL U COMPONENT ARRAY FOR ALL NEST BOUNDARIES 
C        VBI - THE FINAL V COMPONENT ARRAY FOR ALL NEST BOUNDARIES
C       Q2BI - THE FINAL TKE ARRAY FOR ALL NEST BOUNDARIES 
C      CWMBI - THE FINAL CLOUD WATER ARRAY FOR ALL NEST BOUNDARIES
C
C   SUBPROGRAMS CALLED:
C     UNIQUE:  NONE
C
C----------------------------------------------------------------------
      INCLUDE "parmeta"
C----------------------------------------------------------------------
                             P A R A M E T E R
     1 (LMOP=LM+1,LMOM=LM-1)
C
                             P A R A M E T E R
     1 (CM1=2937.4,CM2=4.9283,CM3=23.5518,EPS=0.622
     2, RD=287.04,G=9.80,RG=1./G,EPSQ=2.E-12,EPSQ2=0.2
     3, B1=11.87799326)
C--------------------------------------------------------------------
                             P A R A M E T E R
     & (K15=SELECTED_REAL_KIND(15))
C
                             R E A L
     & (KIND=K15) B,C,ALPET,ALPETK
     &,           ALP(KBIMX,LM+1),ALPETA(KBIMX),ALPSQ(KBIMX,LM+1)
C--------------------------------------------------------------------
                               R E A L
     1 PDBI(KBIMX,NBCEX)
     2,TBI(KBIMX,LMIMX,NBCEX),QBI(KBIMX,LMIMX,NBCEX)
     3,UBI(KBIMX,LMIMX,NBCEX),VBI(KBIMX,LMIMX,NBCEX)
     4,Q2BI(KBIMX,LMIMX,NBCEX),CWMBI(KBIMX,LMIMX,NBCEX)
C
                               R E A L
     1 PDBO(KBIMX,NBCEX)
     2,TBO(KBIMX,LM,NBCEX),QBO(KBIMX,LM,NBCEX)
     3,UBO(KBIMX,LM,NBCEX),VBO(KBIMX,LM,NBCEX)
     4,Q2BO(KBIMX,LM,NBCEX),CWMBO(KBIMX,LM,NBCEX)
     5,HTMBO(KBIMX,LM,NBCEX),USTARBO(KBIMX,NBCEX)
C----------------------------------------------------------------------
                               R E A L
     1  ETAO(LM+1),ETAI(LMIMX+1)
     2, DETAO(LM),AETAO(LM),DETAI(LMIMX),AETAI(LMIMX),PINTO(KBIMX,LM+1)
     3, ALPR(KBIMX,LM)
     4, PDSLO(KBIMX),ZETAO(KBIMX,LM+1),ZETAI(KBIMX,LMIMX)
     5, DLT(KBIMX,LMIMX),PMIDO(KBIMX,LM)
     6, PINTI(LMIMX),ZSFCI(KBIMX),REFI(KBIMX),DFLO(LM+1)
     7, RESO(KBIMX),PSFCO(KBIMX),TX(KBIMX)
     8, PDIJ(KBIMX),PETA(KBIMX)
     9, PDVP(KBIMX),PMVP(KBIMX,2),PMIDI(KBIMX)
     O, UTOP(KBIMX),VTOP(KBIMX),Q2TOP(KBIMX),ALPTU(KBIMX),PHUB(KBIMX)
     1, ALPSFC(KBIMX),Q2SFC(KBIMX)
C-----------------------------------------------------------------------
                               I N T E G E R
     1  LDT1(KBIMX,LMIMX),LMHO(KBIMX)
     2, IDX(KBIMX),IDXX(KBIMX),KNTN(KBIMX),K1(KBIMX),K2(KBIMX),LL(KBIMX)
     3, KSTRTX(KBIMX)
C
                               I N T E G E R
     1 IMI(9),JMI(9),LMI0(9)
C
      CHARACTER DETAFIL*50
C-----------------------------------------------------------------------
C&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
C-----------------------------------------------------------------------
C***
C***  DEFINITION OF CONSTANTS NEEDED FOR VERTICAL INTERPOLATION
C***
      CP=1004.6
      P00=1.E5
      GAMMA=0.0065
      PRF0=101325.
      T0=288.
      GORG=G/(RD*GAMMA)
      RGOG=1./GORG
      ROG=RD/G
      CPOG=CP/G
      CAPA=RD/CP
      ALPT=ALOG(PT)
C-------------------------------------------------------------
      IUDETI0=IUDETI
C-------------------------------------------------------------
C-------------------------------------------------------------
C***
C***  THE GRAND LOOP OVER EACH OF THE NESTS
C***
C-------------------------------------------------------------
C-------------------------------------------------------------
      DO 1000 NB=1,NBCEX
C-------------------------------------------------------------
      KBI=2*IMI(NB)+JMI(NB)-3
      LMI=LMI0(NB)
C-------------------------------------------------------------
ccccc IUDETI0=IUDETI0+1
      CLOSE(IUDETI0)
      WRITE(DETAFIL,10)NB
   10 FORMAT('deta.',I2.2)
      OPEN(UNIT=IUDETI0,FILE=DETAFIL,FORM='UNFORMATTED',IOSTAT=IER)
      REWIND IUDETI0
      READ(IUDETI0)(DETAI(L),L=1,LMI)
C
      DO K=1,KBI
        ZSFCI(K)=0.
        REFI(K)=1.
      ENDDO
C***
C***  SINCE WE ARE ON THE NEST'S BOUNDARY, LMH ALWAYS EQUALS LMI
C***
      LMHI=LMI
C
C----------------------------------------------------------------------
C***  COMPUTE ETA AT THE INTERFACES OF THE MODEL LAYERS
C----------------------------------------------------------------------
C***
C***  FIRST ON THE OUTER GRID
C***
      ETAO(1)=0.
C
      DO L=2,LM
        ETAO(L)=ETAO(L-1)+DETAO(L-1)
      ENDDO
C
      ETAO(LM+1)=1.
      DETAO(LM)=1.-ETAO(LM)
C
C***
C***  THEN ON THE INNER GRID
C***
C
      ETAI(1)=0.
C
      DO L=2,LMI
        ETAI(L)=ETAI(L-1)+DETAI(L-1)
      ENDDO
C
      ETAI(LMI+1)=1.
      DETAI(LMI)=1.-ETAI(LMI)
C***
C***  COMPUTE AETA IN THE MIDDLES OF THE LAYERS
C***
      DO L=1,LMI
         AETAI(L)=0.5*(ETAI(L)+ETAI(L+1))
      ENDDO
      DO L=1,LM
         AETAO(L)=0.5*(ETAO(L)+ETAO(L+1))
      ENDDO
C----------------------------------------------------------------
C***
C***  COMPUTE THE HEIGHTS OF THE STANDARD INTERFACES
C***  FOR THE OUTER GRID
C***
      DO L=1,LM+1
        DFLO(L)=G*T0*(1.-((PT+ETAO(L)*(PRF0-PT))/PRF0)**RGOG)/GAMMA
      ENDDO
C***
C***  GENERATE THE LMH AND THE SURFACE ETA ARRAYS
C***
      DO K=1,KBI
        LMHO(K)=0
      ENDDO
C
      DO L=2,LM
        DO K=1,KBI
          IF(HTMBO(K,L,NB).LT.0.1.AND.LMHO(K).EQ.0)THEN
            LMHO(K)=L-1
            RESO(K)=1./ETAO(L)
          ENDIF
        ENDDO
      ENDDO
C
      DO K=1,KBI
        IF(LMHO(K).EQ.0)THEN
          LMHO(K)=LM
          RESO(K)=1.
        ENDIF
      ENDDO
C----------------------------------------------------------------
C***
C***  COMPUTE HEIGHTS OF THE OUTER GRID ETA INTERFACES
C***
      DO K=1,KBI
        PDSLO(K)=RESO(K)*PDBO(K,NB)
        ZETAO(K,LM+1)=0.0
C       ZETAO(K,LM+1)=FISO(K)*RG
      ENDDO
C       
      DO L=LM,1,-1
        DO K=1,KBI
          DPDE=DETAO(L)*PDSLO(K)
          APEL=PT+AETAO(L)*PDSLO(K)
          RTOP=RD*TBO(K,L,NB)*(1.+0.608*QBO(K,L,NB))/APEL
          FI=G*ZETAO(K,L+1)+RTOP*DPDE
          ZETAO(K,L)=((FI-DFLO(L))*HTMBO(K,L,NB)+DFLO(L))*RG
        ENDDO
      ENDDO
C***
C***  FIND THE PRESSURE ON THE OUTER GRID INTERFACES
C***
      DO K=1,KBI
        PINTO(K,1)=PT
      ENDDO
C
      DO L=2,LM+1
        DO K=1,KBI
          PINTO(K,L)=PDSLO(K)*ETAO(L)+PT
        ENDDO
      ENDDO
C***
C***  FIND THE PRESSURE AT THE OUTER GRID MIDLAYER LOCATIONS
C***
      DO L=1,LM
        DO K=1,KBI
          PMIDO(K,L)=PINTO(K,L)+0.5*DETAO(L)*PDSLO(K)
        ENDDO
      ENDDO
C***
C***  COMPUTE THE LOGS OF PRESSURE FOR UPCOMING INTERPOLATIONS
C***
      DO LD=1,LM
        DO K=1,KBI
          ALPR(K,LD)=ALOG(PINTO(K,LD+1)/PINTO(K,LD))
        ENDDO
      ENDDO
C
      DO LD=1,LM+1
        DO K=1,KBI
          ALP(K,LD)=ALOG(PINTO(K,LD))
          ALPSQ(K,LD)=ALP(K,LD)**2
        ENDDO
      ENDDO
C-------------------------------------------------------------
C***
C***  VERTICAL INTERPOLATION: QUADRATIC INTERPOLATION OF
C***  HEIGHTS (EQUIVALENT TO TEMPERATURE BEING LINEAR IN LN P)
C***  AND LINEAR INTERPOLATION OF WINDS
C***
C-------------------------------------------------------------
C***
C***  ASSIGN TO PD ON THE NEST BOUNDARY THE VALUE OF PDSL AT THE
C***  NEAREST H POINT IN THE PARENT ETA RUN
C***
C----------------------------------------------------------------------
      DO K=1,KBI
        PDIJ(K)=PDSLO(K)
C
C**********************************************************************
C***  THE FOLLOWING LINE IS ONLY TRUE IF PT IS THE SAME
C***  FOR BOTH THE INNER AND OUTER DOMAINS
C**********************************************************************
C
        ZETAI(K,1)=ZETAO(K,1)
C
      ENDDO
C----------------------------------------------------------------------
C***
C***  COMPUTATION OF HEIGHTS AT UPPER BOUNDARIES OF ETA LAYERS 
C***  AND TEMPERATURE AND SPEC HUMIDITY AT MIDDLE OF ETA LAYERS
C***
C***  LDT1 HOLDS THE INDEX OF THE OUTER GRID ETA INTERFACE THAT IS
C***  DIRECTLY ABOVE THE INNER GRID ETA INTERFACE IN QUESTION.
C***
      DO K=1,KBI
        PETA(K)=PT
        LDT1(K,1)=1
      ENDDO
C
C----------------------------------------------------------------------
      MNKNT=1
      DO 100 L=2,LMI
      DETA1=DETAI(L-1)
C
      NUMK=KBI
      KNTX=MNKNT
      KSTRT=KNTX+1
C
      DO K=1,KBI
        PETA(K)=PETA(K)+DETA1*PDIJ(K)
        IDXX(K)=K
        KSTRTX(K)=1
      ENDDO
C
   50 NK=NUMK
      DO NN=1,NK
        IDX(NN)=IDXX(NN)
      ENDDO
      KNTIN=KNTX
      NUMK=0
C
cdir$ ivdep
      DO N=1,NK
        K=IDX(N)
        PETAK=PETA(K)
        PINTOA=PINTO(K,KNTIN)
        PINTOB=PINTO(K,KNTIN+1)
        IF(PETAK.GT.PINTOA.AND.PETAK.LE.PINTOB.OR.
     1     PETAK.LT.PINTO(K,2))THEN
          LDT1(K,L)=MIN(KNTIN,LM-1)
          IF(PETAK+DETAI(L)*PDIJ(K).LE.PINTOB)THEN
            KSTRTX(K)=MIN(KSTRT,KNTIN)
          ELSE
            KSTRTX(K)=MIN(KSTRT,KNTIN+1)
          ENDIF
        ELSE
          NUMK=NUMK+1
          IDXX(NUMK)=K
          KNTX=KNTIN+1
        ENDIF
      ENDDO
C
      IF(NUMK.GT.0)GO TO 50
      IF(NUMK.EQ.0)THEN
        MNKNT=MNKNT+1
        DO K=1,KBI
          MNKNT=MIN(MNKNT,KSTRTX(K))
        ENDDO
      ENDIF
  100 CONTINUE
C
      DO 110 L=1,LMI
      DO K=1,KBI
        L1=LDT1(K,L)
        DLT(K,L)=ALPR(K,L1+1)*(-ALPR(K,L1))*ALOG(PINTO(K,L1+2)
     1                                           /PINTO(K,L1))
      ENDDO    
  110 CONTINUE
C----------------------------------------------------------------------
C***
C***  NOW COMPUTE THE HEIGHTS OF THE INNER GRID INTERFACES
C***
C----------------------------------------------------------------------
                               DO 200 L=2,LMI
C----------------------------------------------------------------------
      DO 190 K=1,KBI
C
      L1=LDT1(K,L)
      L2=L1+1
      L3=L1+2
      H1=ZETAO(K,L1)
      H2=ZETAO(K,L2)
      H3=ZETAO(K,L3)
      H3M2=H3-H2
      H3M1=H3-H1
C***
C***  THE COEFFICIENTS FOR QUADRATIC INTERPOLATION OF HEIGHTS.
C***
      B=(H3M2*(ALPSQ(K,L3)-ALPSQ(K,L1))-H3M1*(ALPSQ(K,L3)-ALPSQ(K,L2)))
     1        /DLT(K,L)
      C=(H3M1*(ALP(K,L3)-ALP(K,L2))-H3M2*(ALP(K,L3)-ALP(K,L1)))
     1        /DLT(K,L)
C
      ARGLOG=PT+PDIJ(K)*ETAI(L)
      ALPET=ALOG(ARGLOG)
C***
C***  THE HEIGHTS OF THE INNER GRID ETA INTERFACES.
C***
      ZETAI(K,L)=H2+B*(ALPET-ALP(K,L2))+C*(ALPET**2-ALPSQ(K,L2))
C
  190 CONTINUE
C
  200 CONTINUE
C
      DO K=1,KBI
        PDBI(K,NB)=PDIJ(K)
      ENDDO
C--------------------------------------------------------------
C***
C***  COMPUTATION OF WINDS WITHIN ETA LAYERS
C***
C--------------------------------------------------------------
      IMT=2*IMI(NB)-1
      JMT=JMI(NB)/2+1
C***
C***  ARRANGEMENT OF LOCATIONS IN THE V BOUNDARY FILE
C***
C                       2ND
C           IBBP--------->----------ITB
C           ILB                     KBI
C            |                       |
C       3RD /|\                     /|\ 4TH
C            |                       |
C           IMT                     ILBP
C            1----------->----------IBB
C                       1ST
      IBB=IMT/2
      IBBP=IBB+1
      ITB=IMT-1
      ILB=IMT+JMT-2
      ILBP=ILB+1
      ILBP2=ILB+2
      ILBM=ILB-1
      IMTP=IMT+1
      KBM=KBI-1
C***
C***  FIND AVERAGE PD OVER VELOCITY BOUNDARY POINTS
C***
      DO K=1,IBB
        K1(K)=K
        K2(K)=K+1
      ENDDO
      DO K=IBBP,ITB
        K1(K)=K+1
        K2(K)=K+2
      ENDDO
      DO K=IMTP,ILBM
        K1(K)=K+1
        K2(K)=K+2
      ENDDO
      DO K=ILBP2,KBM
        K1(K)=K
        K2(K)=K+1
      ENDDO
      K1(IMT)=1
      K2(IMT)=IMT+2
      K1(ILB)=ILB+1
      K2(ILB)=IBB+2
      K1(ILBP)=IBB+1
      K2(ILBP)=ILBP2
      K1(KBI)=KBI
      K2(KBI)=ITB+2
C***
C***  AVERAGE THE PRESSURE AT VELOCITY POINTS.
C***
      DO 210 K=1,KBI
      K1K=K1(K)
      K2K=K2(K)
      PDVP(K)=0.5*(PDBI(K1K,NB)+PDBI(K2K,NB))
      PMVP(K,1)=0.5*(PMIDO(K1K,1)+PMIDO(K2K,1))
      ALP(K,1)=ALOG(PMVP(K,1))
  210 CONTINUE
C
      DO 225 L=2,LM
      DO K=1,KBI
        K1K=K1(K)
        K2K=K2(K)
        PMVP(K,2)=0.5*(PMIDO(K1K,L)+PMIDO(K2K,L))
        ALPR(K,L-1)=ALOG(PMVP(K,2)/PMVP(K,1))
        ALP(K,L)=ALOG(PMVP(K,2))
        PMVP(K,1)=PMVP(K,2)
      ENDDO   
  225 CONTINUE
C-------------------------------------------------------------------
C***
C***  INTERPOLATE THE WINDS
C***
C-------------------------------------------------------------------
C***  FOR EACH MODEL LAYER IN THE NEST DOMAIN, SEARCH THROUGH
C***  LEVELS IN THE PARENT DOMAIN.  WHEN TWO PARENT LEVELS
C***  SANDWICH POINTS ON A NEST LEVEL, DO THE VERTICAL
C***  INTERPOLATION TO THE NEST LEVEL THEN REMOVE THOSE
C***  NEST LEVEL POINTS FROM CONSIDERATION AND CONTINUE SEARCHING
C***  PARENT LEVELS FOR THE REMAINING POINTS ON THE GIVEN
C***  NEST LEVEL.
C***
      MNKNT=1
      DO 270 L=1,LMI
C
C***  ALPETA IS THE LN OF PRESSURE ON THE 
C***  INNER GRID MID-LEVEL
C
      DO K=1,KBI
        ALPETA(K)=ALOG(PT+PDVP(K)*AETAI(L))
        IDXX(K)=K
      ENDDO
C
      NUMK=KBI
      KNTX=MNKNT
C
  230 NK=NUMK
      DO NN=1,NK
        IDX(NN)=IDXX(NN)
      ENDDO
      KNTIN=KNTX
      KP1=MIN(LM,KNTIN+1)
      NUMK=0
C
cdir$ ivdep
      DO N=1,NK
        K=IDX(N)
        ALPETK=ALPETA(K)
        IF(ALPETK.GT.ALP(K,KNTIN).AND.ALPETK.LE.ALP(K,KP1)
     1    .AND.ABS(UBO(K,KP1,NB)).GT.0..OR.ALPETK.LT.ALP(K,1))THEN
          ULD=UBO(K,KP1,NB)
          VLD=VBO(K,KP1,NB)
          CF=(ALP(K,KP1)-ALPETK)/ALPR(K,KNTIN)
          UKB=ULD+(UBO(K,KNTIN,NB)-ULD)*CF
          VKB=VLD+(VBO(K,KNTIN,NB)-VLD)*CF
          UTOP(K)=UKB
          VTOP(K)=VKB
          ALPTU(K)=ALPETK
          LL(K)=L
          UBI(K,L,NB)=UKB
          VBI(K,L,NB)=VKB
C***
C***  IF THE NEW ETA GROUND SURFACE IS BELOW THE OLD GROUND, THEN
C***  ASSUME THE WIND COMPONENTS DECREASE TO ZERO AT THE SURFACE.
C***  THE SAME IS TRUE FOR ANY NEW ETA LEVEL BELOW THE LOWEST
C***  OLD ETA SURFACE.
C***
        ELSEIF(ALPETK.GT.ALP(K,LM).OR.UBO(K,KP1,NB).EQ.0.)THEN
          UDIF=UTOP(K)
          VDIF=VTOP(K)
          ALPDIF=ALPTU(K)-ALOG(PT+PDVP(K))
          LL(K)=LL(K)+1
          LLK=LL(K)
          UBI(K,LLK,NB)=UTOP(K)-UDIF*(ALPTU(K)-ALPETK)/ALPDIF
          VBI(K,LLK,NB)=VTOP(K)-VDIF*(ALPTU(K)-ALPETK)/ALPDIF
        ELSEIF(ALPETK.GT.ALP(K,KP1))THEN
          NUMK=NUMK+1
          IDXX(NUMK)=K
          KNTX=KNTIN+1
        ENDIF
      ENDDO
      IF(NUMK.GT.0)THEN
        IF(NUMK.EQ.KBI)MNKNT=MNKNT+1
        GO TO 230
      ENDIF
c     IF(NUMK.EQ.0)MNKNT=KNTX
C
  270 CONTINUE
C
C--------------------------------------------------------------------
C***
C***  INTERPOLATE Q AND CLD WATER LINEARLY IN LN(P)
C***
cfpp$ noconcur l
      DO 300 L=LMI,1,-1
      DO K=1,KBI
        IF(L.EQ.LMI)PHUB(K)=0.
        PHLB=PHUB(K)
        PHUB(K)=G*ZETAI(K,L)
        TBI(K,L,NB)=-(PHLB-PHUB(K))*(PT+AETAI(L)*PDBI(K,NB))
     1             /(RD*DETAI(L)*PDBI(K,NB))
        IF(TBI(K,L,NB).LT.150.)TBI(K,L,NB)=150.
      ENDDO
  300 CONTINUE
C
      DO K=1,KBI
        ALP(K,1)=ALOG(PMIDO(K,1))
        PSFCO(K)=PDBO(K,NB)+PT
      ENDDO
C
      DO LD=2,LM
      DO K=1,KBI
        ALPR(K,LD-1)=ALOG(PMIDO(K,LD)/PMIDO(K,LD-1))
        ALP(K,LD)=ALOG(PMIDO(K,LD))
      ENDDO
      ENDDO
C
      MNKNT=1
      DO 340 L=1,LMHI
C
      DO K=1,KBI
        PMIDI(K)=AETAI(L)*PDBI(K,NB)+PT
        ALPETA(K)=ALOG(PMIDI(K))
        IDXX(K)=K
      ENDDO
C
      NUMK=KBI
      KNTX=MNKNT
C
  310 NK=NUMK
      DO NN=1,NK
        IDX(NN)=IDXX(NN)
      ENDDO
      KNTIN=KNTX
      KP1=KNTIN+1
      NUM0=0
      NUMK=0
C
cdir$ ivdep
      DO N=1,NK
        K=IDX(N)
        LKO=LMHO(K)
        ALPETK=ALPETA(K)
        IF(ALPETK.GT.ALP(K,KNTIN).AND.ALPETK.LE.ALP(K,KP1)
     1    .AND.KP1.LE.LKO.OR.ALPETK.LT.ALP(K,1)
     2    .OR.ALPETK.GT.ALP(K,KP1).AND.PMIDI(K).LE.PSFCO(K)
     3    .AND.KP1.EQ.LKO)THEN
          QLD=QBO(K,KP1,NB)
          CWMLD=CWMBO(K,KP1,NB)
          CF=(ALP(K,KP1)-ALPETA(K))/ALPR(K,KNTIN)
          CWMBI(K,L,NB)=CWMLD+(CWMBO(K,KNTIN,NB)-CWMLD)*CF
          QETAIJ=QLD+(QBO(K,KNTIN,NB)-QLD)*CF
          CLOGES=-CM1/TBI(K,L,NB)-CM2*ALOG10(TBI(K,L,NB))+CM3
          ESE=10.**(CLOGES+2.)
          QSX=EPS*ESE/(PMIDI(K)-ESE*(1.-EPS))
          QSMX=0.98*QSX
          QBI(K,L,NB)=AMIN1(QETAIJ,QSMX,20.E-3)
C***
C***  IF THE NEW ETA GROUND SURFACE IS BELOW THE OLD GROUND THEN
C***  ASSUME THE RH (GIVEN T) IS EQUAL TO THAT ON THE LOWEST
C***  OLD SURFACE AND THAT THE CLOUD WATER IS ZERO.
C***  NOTE:  TBI INSIDE THE NEXT BLOCK IS ACTUALLY THE VIRTUAL
C***         TEMPERATURE YET IT IS BEING USED AS THE TRUE T
C***
        ELSEIF(PMIDI(K).GT.PSFCO(K))THEN
          CLOGES=-CM1/TBO(K,LKO,NB)-CM2*ALOG10(TBO(K,LKO,NB))+CM3
          ESS=10.**(CLOGES+2.)
          QSLDM=EPS*ESS/(PMIDO(K,LKO)-ESS*(1.-EPS))
          RHLDM=QBO(K,LKO,NB)/QSLDM
          RHLDM=AMIN1(RHLDM,0.98)
          CLOGES=-CM1/TBI(K,L,NB)-CM2*ALOG10(TBI(K,L,NB))+CM3
          ESL=10.**(CLOGES+2.)
          QSL=EPS*ESL/(AETAI(L)*PDBI(K,NB)+PT-ESL*(1.-EPS))
          QBI(K,L,NB)=RHLDM*QSL
          QBI(K,L,NB)=AMIN1(QBI(K,L,NB),20.E-3)
          CWMBI(K,L,NB)=0.
        ELSEIF(ALPETK.GT.ALP(K,KP1))THEN
          NUMK=NUMK+1
          IDXX(NUMK)=K
          KNTX=KNTIN+1
        ENDIF
      ENDDO
      IF(NUMK.GT.0)THEN
        IF(NUMK.EQ.KBI)MNKNT=MNKNT+1
        GO TO 310
      ENDIF
c     IF(NUMK.EQ.0)MNKNT=KNTX
  340 CONTINUE
C
      DO L=1,LMI
      DO K=1,KBI
        QBI(K,L,NB)=AMAX1(QBI(K,L,NB),EPSQ)
        CWMBI(K,L,NB)=AMAX1(CWMBI(K,L,NB),0.)
      ENDDO   
      ENDDO   
C
C***
C***  CONVERT FROM VIRTUAL TO TRUE TEMPERATURE.
C***
      DO 375 L=1,LMI
C
      DO K=1,KBI
        TBI(K,L,NB)=AMIN1(TBI(K,L,NB),350.)
        TBI(K,L,NB)=AMAX1(TBI(K,L,NB),150.)
        TX(K)=TBI(K,L,NB)
        IDXX(K)=K
      ENDDO
      NUMK=KBI
C
      ITER=0
  360 ITER=ITER+1
      NK=NUMK
      NUMK=0
      DO NN=1,NK
        IDX(NN)=IDXX(NN)
      ENDDO
C
cdir$ ivdep
      DO N=1,NK
        K=IDX(N)
        TBI(K,L,NB)=TX(K)/(1.+0.608*QBI(K,L,NB))
        CLOGES=-CM1/TBI(K,L,NB)-CM2*ALOG10(TBI(K,L,NB))+CM3
        ESE=10.**(CLOGES+2.)
        DENOM2=AETAI(L)*PDBI(K,NB)+PT-ESE*(1.-EPS)
        QSX=EPS*ESE/DENOM2
        QSMX=0.98*QSX
        IF(QBI(K,L,NB).GT.QSMX)THEN
          QBI(K,L,NB)=QSMX
          NUMK=NUMK+1
          IDXX(NUMK)=K
        ENDIF
      ENDDO
      IF(NUMK.GT.0.AND.ITER.LE.2)GO TO 360
C
  375 CONTINUE
C--------------------------------------------------------------------
C***
C***  INTERPOLATE Q2 LINEARLY IN LN(P) (BUT ON INTERFACES)
C***
      DO K=1,KBI
        Q2SFC(K)=AMAX1(B1**(2./3.)*USTARBO(K,NB)**2,EPSQ2)
        ALPSFC(K)=ALOG(PDBI(K,NB)*REFI(K)+PT)
        ALP(K,LMOP)=ALOG(PINTO(K,LMOP))
      ENDDO
C
      DO 400 LD=2,LM
      DO K=1,KBI
        ALPR(K,LD)=ALOG(PINTO(K,LD+1)/PINTO(K,LD))
        ALP(K,LD)=ALOG(PINTO(K,LD))
      ENDDO
  400 CONTINUE
C
      MNKNT=2
      LKI=LMHI-1
      DO 450 L=1,LKI
C
      DO K=1,KBI
        ALPETA(K)=ALOG(ETAI(L+1)*PDBI(K,NB)+PT)
        IDXX(K)=K
      ENDDO
C
      NUMK=KBI
      KNTX=MNKNT
C
  410 NK=NUMK
      DO NN=1,NK
        IDX(NN)=IDXX(NN)
      ENDDO
C
      KNTIN=KNTX
      KP1=KNTIN+1
      NUMK=0
C
cdir$ ivdep
      DO N=1,NK
        K=IDX(N)
        LKO=LMHO(K)
        ALPETK=ALPETA(K)
        IF(ALPETK.GT.ALP(K,KNTIN).AND.ALPETK.LE.ALP(K,KP1)
     1     .AND.KP1.LE.LKO.OR.ALPETK.LT.ALP(K,2))THEN
          Q2LD=Q2BO(K,KP1,NB)
          CF=(ALP(K,KP1)-ALPETK)/ALPR(K,KNTIN)
          Q2IJ=Q2LD+(Q2BO(K,KNTIN,NB)-Q2LD)*CF
          Q2TOP(K)=AMAX1(Q2IJ,EPSQ2)
          ALPTU(K)=ALPETK
          LL(K)=L
          Q2BI(K,L,NB)=Q2TOP(K)
C***
C***  IF THE NEW Q2 INTERFACE IS BELOW THE TOP OF THE LOWEST 
C***  OUTER GRID LAYER, ASSUME Q2 DECREASES LINEARLY TO Q2 AT
C***  THE NEW GROUND.
C***
        ELSEIF(ALPETK.GT.ALP(K,LKO))THEN
          Q2DIF=Q2TOP(K)-Q2SFC(K)
          ALPDIF=ALPTU(K)-ALPSFC(K)
          LL(K)=LL(K)+1
          LLK=LL(K)
          Q2BI(K,LLK,NB)=Q2TOP(K)-Q2DIF*(ALPTU(K)-ALPETK)/ALPDIF
        ELSEIF(ALPETK.GT.ALP(K,KP1))THEN
          NUMK=NUMK+1
          IDXX(NUMK)=K
          KNTX=KNTIN+1
        ENDIF
      ENDDO
      IF(NUMK.GT.0)THEN
        IF(NUMK.EQ.KBI)MNKNT=MNKNT+1
        GO TO 410
      ENDIF
c     IF(NUMK.EQ.0)MNKNT=KNTX
C
  450 CONTINUE
C-----------------------------------------------------------
C***
C***  NOW REDEFINE PD TO HAVE IT EQUAL TO PSFC-PT.
C***
      DO K=1,KBI
        PDBI(K,NB)=PDBI(K,NB)*REFI(K)
        Q2BI(K,LMHI,NB)=Q2SFC(K)
      ENDDO
C------------------------------------------------------------
 1000 CONTINUE
C------------------------------------------------------------
      RETURN
      END
