      SUBROUTINE CONVEC6 (TH,QH,PCP,TCV,PRESSG, 
     1                   SFCT,QSCRN,DEL,SHJ,SHTJ,SHXKJ,QADD,
     2                   ILEV,ILEVM,LEV,ILG,IL1,IL2,MSG,NUPS,NSUPS, 
     3                   Q1,RGOCP,
     4                   T,Q,P,ST,TSCRN,
     5                   PCPDH,SPHI,SPHJ,GAC,GAM, 
     6                   SUMQ,DQR,TDEEP,DTDEEP,SUMDEEP, 
     7                   DTF,DTH,DQF,DQH,DTFMAX,
     8                   QSI,DQI,QSJ,DQJ,HT,
     9                   H,HS,EPH,Y,HZ, 
     A                   EI,EJ,DQEX,GAMSAT,ICON,
     B                   ITER,IADJ,IWET)
  
C     * JAN 30/89 - M.LAZARE. ICE-PHASE PARAMETERS PASSED FROM SPWCONH
C     *                       COMMON BLOCK EPSICE. IN-LINE FUNCTIONS
C     *                       HTVOCP AND ESTEFF USED TO EVALUATE L/CP 
C     *                       AND ES(T) FOR ICE/WATER PHASES, WITH ICE
C     *                       PHASE INITITATING AT 0 DEG LINEARLY TO
C     *                       COMPLETE ICE-PHASE AT -15 DEG (LIKE OLD)
C     *                       THEREFORE ICE-PHASE SECTION AT END
C     *                       REMOVED.
C     *                       NEW CALCULATION OF MOISTURE REQUIRED TO 
C     *                       FILL HOLES OF NEGATIVE Q. 
C     * DEC 3/88 - M.LAZARE. DERIVED FROM CONVEC5 BUT VECTORIZED. 
C 
C     * NOV 9/88 - M.LAZARE,N.MCFARLANE. VERSION CONVEC5. LIKE CONVC3V
C     *            EXCEPT:  
C     *                    A.) RH ARRAY REMOVED.
C     *                    B.) SUPERSATURATION W.R.T. ICE (BELOW -40C)
C     *                        REMOVED AT END.
C     *                    C.) LTOP DEFINED BASED ON POTENTIAL FOR
C     *                        CONVECTION (USING TSCRN,QSCRN).
C 
C     *****   NOV 1975  -  CLAUDE GIRARD, ROGER DALEY   ***** 
  
C     *    PERFORMS MOIST AND DRY CONVECTIVE ADJUSTMENTS AND
C     *    COMPUTES BOTH LARGE AND SMALL SCALE PRECIPITATION
  
C     **************** INDEX OF VARIABLES *********************** 
C     *  I    => INPUT ARRAYS.
C     * I/O   => INPUT/OUTPUT ARRAYS. 
C     *  W    => WORK ARRAYS. 
C     * IC    => INPUT DATA CONSTANTS.
C     *  C    => DATA CONSTANTS PERTAINING TO SUBROUTINE ITSELF.
  
C  I  * DEL      LOCAL SIGMA HALF-LEVEL THICKNESS (I.E. DSHJ).
C  W  * DQEX     MOISTURE REMOVED BY CONDENSATION.
C  W  * DQF      MOISTURE CHANGE CAUSED BY A CONVECTIVE ADJUSTMENT. 
C  W  * DQH      MOISTURE CHANGE CAUSED BY A RELEASE OF LATENT HEAT.
C  W  * DQI,DQJ  SATURATION DEFICITS IN "I" AND "J" ADJACENT LEVELS.
C  W  * DTDEEP   DIFFERENCE BETWEEN AMBIENT TEMPERATURE AND TDEEP.
C  C  * DTDSCR   ASYMPTOTIC VALUE FOR MULTIPLE ADJUSTMENT CRITERIA. 
C  W  * DTF      TEMPERATURE CHANGE CAUSED BY A CONVECTIVE ADJUSTMENT.
C  W  * DTFMAX   MAXIMUM MASS-WEIGHTED DTF IN A COLUMN. 
C  W  * DTH      TEMPERATURE CHANGE CAUSED BY A RELEASE OF LATENT HEAT. 
C  W  * DQR      MOISTURE EQUIVALENT TO PRECIPITATION AMOUNT. 
C  W  * EI,EJ    VAPOUR PRESSURE FOR EACH ADJACENT LAYER PAIR.
C  W  * EPH      FACTOR TO CALCULATE MOISTURE RELEASE UPON CONDENSATION.
C  W  * GAM,GAC  STABILITY FACTOR AND ITS CRITICAL VALUE. 
C  W  * GAMSAT   SATURATED ADIABATIC LAPSE RATE.
C  W  * H        RELATIVE HUMIDITY. 
C  W  * HS       CRITICAL SATURATION RELATIVE HUMIDITY. 
C  W  * HT       RATIO OF LATENT HEAT OF VAPOUR TO SPECIFIC HEAT. 
C  W  * HZ       FACTOR TO CALCULATE MOISTURE RELEASE UPON CONDENSATION.
C  W  * IADJ     INDICATES WHETHER FURTHER ADJUSTMENT REQUIRED(=1) OR NOT.
C  W  * ICON     HOLDS POSITION OF GATHERED POINTS VS LONGITUDE INDEX.
C  W  * IWET     INDICATES WHETHER POINT WAS SATURATED UNSTABLE.
C IC  * ILEV     NUMBER OF MODEL LEVELS.
C IC  * ILEVM    ILEV-1.
C IC  * ILG      LON+2 = SIZE OF GRID SLICE.
C IC  * IL1,IL2  START AND END LONGITUDE INDICES FOR LATITUDE CIRCLE. 
C  W  * ITER     NUMBER OF ITERATIONS AT EACH LONGITUDE POINT.
C IC  * LEV      ILEV+1.
C  C  * MAD      LOGICAL SWITCH TO INDICATE WHETHER TO MULTIPLE ADJUST. 
C IC  * MSG      NUMBER OF MISSING MOISTURE LEVELS AT THE TOP OF MODEL. 
C  C  * NSUPS    NUMBER OF SATURATED UNSTABLE POINTS. 
C  C  * NUPS     NUMBER OF DRY UNSTABLE POINTS. 
C  W  * P        GRID SLICE OF AMBIENT PRESSURE IN MBS. 
C I/O * PCP      ROW OF PRECIPITATION RATE. 
C  W  * PCPDH    SCALED SURFACE PRESSURE. 
C  I  * PRESSG   ROW OF SURFACE PRESSURE IN PA. 
C  W  * Q        GATHERED SLICE OF SPECIFIC HUMIDITY. 
C  O  * QADD     "FICTITIOUS" SOURCE OF SPECIFIC HUMIDITY.
C I/O * QH       GRID SLICE OF SPECIFIC HUMIDITY. 
C  I  * QSCRN    ROW OF SCREEN-LEVEL SPECIFIC HUMIDITY. 
C  W  * QSI,QSJ  CRITICAL SPECIFIC HUMIDITY FOR EACH ADJACENT LAYER PAIR
C IC  * Q1       SPECIFIED CONSTANT SPECIFIC HUMIDITY IN "MSG" LEVELS.
C IC  * RGOCP    RGAS/CP = "KAPA".
C  W  * SFCT     ROW OF SCREEN-LEVEL TEMPERATURE. 
C  I  * SHJ      GRID SLICE OF LOCAL HALF-LEVEL SIGMA VALUES. 
C  I  * SHTJ     GRID SLICE OF LOCAL HALF-LEVEL SIGMA INTERFACE VALUES. 
C  I  * SHXKJ    GRID SLICE OF LOCAL SH**KAPA.
C  W  * SPHI,J   SPECIFIC HUMIDITY FOR EACH ADJACENT LAYER PAIR.
C  W  * ST       STABILITY MATRICES.
C  W  * SUMDEEP  VERTICAL SUM OF POSITIVE AREA BETWEEN "T" AND "TDEEP". 
C  W  * SUMQ     MOISTURE EQUIVALENT TO PRECIP AMOUNT FOR ICE PHASE.
C  W  * T        GATHERED GRID SLICE OF TEMPERATURE.
C I/O * TCV      ROW OF INDICES OF TOP MOST CONVECTIVE LAYER. 
C  W  * TDEEP    TEMPERATURES ON "PARCEL" CURVE FROM SCREEN-LEVEL.
C I/O * TH       GRID SLICE OF TEMPERATURE. 
C  W  * TSCRN    SCREEN-LEVEL TEMPERATURE ADIABATIC EXTRAPOLATE TO SFC. 
C  W  * Y        FACTOR TO CALCULATE MOISTURE RELEASE UPON CONDENSATION.
  
C     *********************************************************** 
  
C     * MULTI-LEVEL I/O FIELDS: 
  
      IMPLICIT REAL (A-H,O-Z),
     +INTEGER (I-N)
      REAL TH(ILG,ILEV), QH(ILG,ILEV), DEL(ILG,ILEV),   SHJ(ILG,ILEV),
     1   SHTJ(ILG, LEV), SHXKJ(ILG,ILEV), QADD(ILG) 
  
C     * SINGLE-LEVEL I/O FIELDS:  
  
      REAL     PCP(ILG),     TCV(ILG),  PRESSG(ILG),    SFCT(ILG) 
      REAL   QSCRN(ILG) 
  
C     * WORK FIELDS:  
  
      REAL  T(ILG,ILEV),  Q(ILG,ILEV),  P(ILG,ILEV),    ST(ILG,ILEVM,6) 
  
      REAL   TSCRN(ILG),   PCPDH(ILG),    SPHI(ILG),    SPHJ(ILG) 
      REAL     GAC(ILG),     GAM(ILG),    SUMQ(ILG),     DQR(ILG) 
      REAL   TDEEP(ILG),  DTDEEP(ILG), SUMDEEP(ILG),     DTF(ILG) 
      REAL     DTH(ILG),     DQF(ILG),     DQH(ILG),  DTFMAX(ILG) 
      REAL     QSI(ILG),     DQI(ILG),     QSJ(ILG),     DQJ(ILG) 
      REAL      HT(ILG),       H(ILG),      HS(ILG),     EPH(ILG) 
      REAL       Y(ILG),      HZ(ILG),      EI(ILG),      EJ(ILG) 
      REAL    DQEX(ILG),  GAMSAT(ILG) 
  
      INTEGER ICON(ILG),    ITER(ILG),    IADJ(ILG),    IWET(ILG) 
  
      LOGICAL MAD 
  
      COMMON/ADJPCP/ HC,HF,HM,AA,DEPTH,LHEAT,MOIADJ,MOIFLX
      COMMON/EPS   / A, B, EPS1, EPS2 
      COMMON/EPSICE/ AICE,BICE,TICE,QMIN
      COMMON/GAMS  / EPSS,CAPA
      COMMON/HTCP  / T1S,T2S,AI,BI,AW,BW,SLP
  
C---------------------------------------------------------------------- 
C     * STATEMENT FUNCTIONS DEFENITIONS.
  
C     * A) LINEAR INTERPOLATION BETWEEN HC AND 1.0 .
  
      CHIC(HHH,HC) = MERGE( 1.0E0, (HHH-HC)/(1.E0-HC), HHH.GE.1.E0)
  
C     * B) CRITICAL SATURATION RELATIVE HUMIDITY. 
  
      CR(HM,HHH,AA)     =     MERGE( AA*(HHH-HM)**3,
     &  AA*(2.E0-HM-HHH)**3+HHH-1.E0,
     &  HHH.LE.1.E0) 
      CRIRLH(HM,HHH,AA) = HHH -
     &  MERGE( CR(HM,HHH,AA), HHH-1.E0, HM+HHH.LE.2)
  
C     * C) COMPUTES THE RATIO OF LATENT HEAT OF VAPORIZATION OF 
C     *    WATER OR ICE TO THE SPECIFIC HEAT OF AIR AT CONSTANT 
C     *    PRESSURE CP. 
  
      TW(TTT)     = AW-BW*TTT 
      TI(TTT)     = AI-BI*TTT 
      HTVOCP(TTT) = MERGE( TW(TTT), 
     1              MERGE( TI(TTT), 
     2                     SLP*((TTT-T2S)*TW(TTT)+(T1S-TTT)*TI(TTT)), 
     3                     TTT.LE.T2S), 
     4                     TTT.GE.T1S)
  
C     * D) COMPUTES THE SATURATION VAPOUR PRESSURE OVER WATER OR ICE. 
  
      ESW(TTT)    = EXP(A-B/TTT)
      ESI(TTT)    = EXP(AICE-BICE/TTT)
      ESTEFF(TTT) = MERGE( ESW(TTT),
     1              MERGE( ESI(TTT),
     2                     SLP*((TTT-T2S)*ESW(TTT)+(T1S-TTT)*ESI(TTT)), 
     3                     TTT.LE.T2S), 
     4                     TTT.GE.T1S)
C---------------------------------------------------------------------- 
      DATA MAD   /.TRUE./ 
      DATA DTDSCR/0.002E0 / 
C-----------------------------------------------------------------------
C     * CALCULATE ROW OF STABILITY MATRICES.
C 
      CALL STAWCL2 (ST,SHTJ(1,2),SHTJ,SHJ,
     1              ILEVM,ILEV,ILEV+1,IL1,IL2,ILG,RGOCP)
C 
C     * INITIALIZE NECESSARY ARRAYS.
C 
      DO 50 IL=IL1,IL2
         ICON(IL)  = IL 
         IWET(IL)  = 0
         IADJ(IL)  = 0
         ITER(IL)  = 0
         SUMQ(IL)  = 0.E0 
         PCP(IL)   = 0.E0 
         TCV(IL)   = FLOAT(ILEV)
         TSCRN(IL) = SFCT(IL)*SHXKJ(IL,ILEV)
         PCPDH(IL) = PRESSG(IL)*DEPTH 
   50 CONTINUE
C 
      DO 60 L=1,ILEV
      DO 60 IL=IL1,IL2
         P(IL,L)   = SHJ(IL,L) * PRESSG(IL) * 0.01E0
   60 CONTINUE
C 
      LEN     = IL2 - IL1 + 1 
      LENGATH = LEN 
C 
C     ****************************************************************
C     * BEGIN ITERATION LOOP *****************************************
C     ****************************************************************
  
  100 CONTINUE
C 
C     * CALCULATE GATHERED TEMPERATURE AND SPECIFIC HUMIDITY GRID 
C     * SLICES. DO-LOOP INDEX "K" REFERS TO GATHERED LONGITUDE POINTS 
C     * FROM NOW ON, WHILE INDEX "IL" REFERS TO SCATTERED OUT VALUES. 
C     * THE ARRAY "ICON" CONTAINS THE GATHERED POINTS' LOCATIONS. 
C 
      DO 110 L=1,ILEV 
      DO 110 K=1,LENGATH
          Q(K,L) = QH(ICON(K),L)
          T(K,L) = TH(ICON(K),L)
  110 CONTINUE
C 
C     * UPDATE ITERATION COUNTER ARRAY. 
C 
      DO 115 K=1,LENGATH
          ITER(ICON(K)) = ITER(ICON(K)) + 1 
  115 CONTINUE
C 
C     * CALCULATE INDICATORS OF HOW HIGH CONVECTION MAY GO. 
C 
      DO 120 K=1,LENGATH
          DQR(K)        = 0.0E0 
          DTFMAX(K)     = 0.0E0 
          EILEV         = ESTEFF(T(K,ILEV)) 
          SPHI(K)       = EPS1 * EILEV /(P(ICON(K),ILEV) - EPS2 * EILEV)
          HT(K)         = HTVOCP(T(K,ILEV)) 
          TBASE         = MAX(T(K,ILEV),TSCRN(ICON(K))) 
          QBASE         = MAX(Q(K,ILEV),QSCRN(ICON(K))) 
          EBASE         = ESTEFF(TBASE) 
          SPHBASE       = EPS1*EBASE/(P(ICON(K),ILEV)-EPS2*EBASE) 
          TDEEP(K)      = TBASE-HT(K)*(SPHBASE-QBASE) 
          SUMDEEP(K)    = DEL(ICON(K),ILEV)*(TDEEP(K)-T(K,ILEV))
  120 CONTINUE
C 
C     * CONVECTIVE ADJUSTMENT "BOOTSTRAPPING" LOOP. 
C     * INDEX "J" REFERS TO LOWER LEVEL AND INDEX "I" TO UPPER LEVEL. 
C 
      DO 200 J=ILEV,2,-1
          I    = J-1
          DO 130 K=1,LENGATH
              GAC(K)     = 0.0E0
              DTF(K)     = 0.0E0
              DTH(K)     = 0.0E0
              DQF(K)     = 0.0E0
              DQH(K)     = 0.0E0
              EJ(K)      = ESTEFF(T(K,J)) 
              SPHJ(K)    = EPS1*EJ(K)/(P(ICON(K),J)-EPS2*EJ(K)) 
              EI(K)      = ESTEFF(T(K,I)) 
              SPHI(K)    = EPS1*EI(K)/(P(ICON(K),I)-EPS2*EI(K)) 
C 
C             ********************************************************
C             * CONVECTIVE HEAT FLUX *********************************
C             ********************************************************
C 
C             * COMPUTE  GAM, H, HS.
C 
              H(K)       = MAX(Q(K,J),QMIN)/SPHJ(K) 
              IF(H(K).GT.HM.AND.ITER(ICON(K)).EQ.1)              THEN 
                 HS(K)   = CRIRLH(HM,H(K),AA) 
              ELSE
                 HS(K)   = MIN(H(K),1.E0) 
              ENDIF 
              TT         = ST(ICON(K),I,1)*T(K,I)+ST(ICON(K),I,2)*T(K,J)
              GAM(K)     = TT+ST(ICON(K),I,3)*(T(K,I)-T(K,J)) 
              EST        = ESTEFF(TT) 
              QST        = EPS1*EST/(P(ICON(K),J)*ST(ICON(K),I,6) 
     1                     -EPS2*EST) 
              HT(K)      = HTVOCP(TT) 
              GA         = HT(K)*QST/(CAPA*TT)
              GB         = GA*EPSS*HT(K)/TT 
              GAMSAT(K)  = TT*(GB-GA)/(1.E0+GB) 
              RATIO      = (1.E0+GA)/(1.E0+GB)
              R1         = ST(ICON(K),I,1)*RATIO + ST(ICON(K),I,3)
              R2         = ST(ICON(K),I,3)-RATIO * ST(ICON(K),I,2)
              TDEEP(K)   = R2 * TDEEP(K) / R1 
              DTDEEP(K)  = TDEEP(K)-T(K,I)
              IF(DTDEEP(K).LE.0.E0.AND.SUMDEEP(K).GT.0.E0.AND.I.GT.MSG) 
     1                     TCV(ICON(K)) = I 
              SUMDEEP(K) = SUMDEEP(K) + DEL(ICON(K),I) * DTDEEP(K)
C 
C             * COMPUTE GAC AND UPDATE WET UNSTABLE POINT INDICATOR 
C             * ARRAY IF APPROPRIATE. 
C 
              IF (H(K)*MOIADJ.GE.HC)                             THEN 
                  GAC(K) = CHIC(HS(K),HC)*GAMSAT(K) 
              ENDIF 
              IF (GAC(K).GT.0.0E0.AND.GAC(K).GE.GAM(K))
     &          IWET(ICON(K)) = 1 
C 
C             * NOW DRY ADJUST. 
C 
              IF (GAM(K).LT.GAC(K))                              THEN 
                  DTF(K)    = ST(ICON(K),I,4)*(GAC(K)-GAM(K)) 
                  DTFMAX(K) = MAX( DTFMAX(K), 
     1                               ABS(DTF(K)) * DEL(ICON(K),J))
                  T(K,I)    = T(K,I) + DTF(K) * ST(ICON(K),I,5) 
                  T(K,J)    = T(K,J) + DTF(K) 
                  EJ(K)     = ESTEFF(T(K,J))
                  SPHJ(K)   = EPS1*EJ(K)/(P(ICON(K),J)-EPS2*EJ(K))
                  EI(K)     = ESTEFF(T(K,I))
                  SPHI(K)   = EPS1*EI(K)/(P(ICON(K),I)-EPS2*EI(K))
              ENDIF 
              QSJ(K)     = HS(K) * SPHJ(K)
              DQJ(K)     = QSJ(K) - Q(K,J)
  130     CONTINUE
C 
C         ********************************************************
C         * MOISTURE FLUX ****************************************
C         ********************************************************
C 
          IF(I.GT.MSG)                                           THEN 
              DO 140 K=1,LENGATH
                  IF (DQJ(K).LT.0.0E0.AND.DTF(K)*MOIFLX.NE.0.0E0) THEN
                      QSI(K) = MAX(HS(K),HF) * SPHI(K)
                      DQI(K) = MAX( QSI(K)-Q(K,I), 0.E0 ) 
                      DQF(K) = MAX( DQI(K)/ST(ICON(K),I,5), DQJ(K) )
                      Q(K,I) = Q(K,I)+DQF(K)*ST(ICON(K),I,5)
                      Q(K,J) = Q(K,J)+DQF(K)
                      DQJ(K) = QSJ(K)-Q(K,J)
                  ENDIF 
  140         CONTINUE
          ENDIF 
C 
C         ************************************************
C         * CONVECTIVE OR STABLE HEATING BY CONDENSATION *
C         ************************************************
C 
          DO 150 K=1,LENGATH
              IF (J.GT.MSG .AND. DQJ(K).LT.0.0E0)                  THEN 
                  HT(K)     = HTVOCP(T(K,J))
                  EPH(K)    = QSJ(K)*EPS2/(HS(K)*EPS1)
                  Y(K)      = B/(T(K,J)*T(K,J))*QSJ(K)*(1.E0+EPH(K))
                  HZ(K)     = HT(K)*Y(K)/(T(K,J)*(1.E0+HT(K)*Y(K))) 
                  DQH(K)    = DQJ(K)/(1.E0+HT(K)*Y(K))
     1                        * (1.E0+DQJ(K)/Y(K)*HZ(K)*HZ(K) 
     2                        * (B*(EPH(K)+0.5E0)-T(K,J)))
                  DTH(K)    = -HT(K)*DQH(K) 
                  T(K,J)    = T(K,J)+DTH(K) 
                  Q(K,J)    = Q(K,J)+DQH(K) 
C 
C                 * REMOVES LEFT OVER SUPERSATURATION TO INSURE 
C                 * MAX(R.H.)<=1. 
C 
                  EJ(K)     = ESTEFF(T(K,J))
                  SPHJ(K)   = EPS1*EJ(K)/(P(ICON(K),J)-EPS2*EJ(K))
                  DQEX(K)   = MAX(Q(K,J) - SPHJ(K), 0.E0) 
                  DQH(K)    = DQH(K) - DQEX(K)
                  Q(K,J)    = Q(K,J) - DQEX(K)
                  T(K,J)    = T(K,J) + HT(K) * DQEX(K)
C 
C                 * CONVECTIVE OR STABLE PRECIPITATION. 
C 
                  DQR(K) = DQR(K) + DEL(ICON(K),J) * DQH(K) 
              ENDIF 
  150     CONTINUE
C 
C         * UPDATE WHETHER GATHERED LONGITUDE POINT NEEDS FURTHER 
C         * ADJUSTING.
C         * THIS LOOP WILL NOT VECTORIZE. 
C 
          DO 160 K=1,LENGATH
              IF (J.GT.MSG .AND. DQJ(K).LT.0.0E0)                  THEN 
                  DTFTEMP   = MERGE(0.E0, ABS(DTH(K))*DEL(ICON(K),J), 
     1                              DTF(K).EQ.0.0E0)
                  DTFMAX(K) = MAX( DTFMAX(K), DTFTEMP)
              ENDIF 
              IF (DTFMAX(K).GT.DTDSCR)    IADJ(ICON(K)) = 1 
  160     CONTINUE
  
  200 CONTINUE
C 
C     ***************************************************************** 
C     * CONDENSATION IN TOP LAYER (IF MOISTURE VARIABLE CARRIED THERE). 
C     ***************************************************************** 
C 
      IF(MSG.LT.1)                                               THEN 
          DO 210 K=1,LENGATH
              H(K)     = Q(K,1)/SPHI(K) 
              HS(K)    = CRIRLH(HM,H(K),AA) 
              QSI(K)   = HS(K)*SPHI(K)
              DQI(K)   = QSI(K)-Q(K,1)
              IF (DQI(K).LT.0.0E0)                                 THEN 
                  HT(K)     = HTVOCP(T(K,1))
                  EPH(K)    = QSI(K)*EPS2/(HS(K)*EPS1)
                  Y(K)      = B/(T(K,1)*T(K,1))*QSI(K)*(1.E0+EPH(K))
                  HZ(K)     = HT(K)*Y(K)/(T(K,1)*(1.E0+HT(K)*Y(K))) 
                  DQH(K)    = DQI(K)/(1.E0+HT(K)*Y(K))
     1                        * (1.E0+DQI(K)/Y(K)*HZ(K)*HZ(K) 
     2                        * (B*(EPH(K)+0.5E0)-T(K,1)))
                  DTH(K)    = -HT(K)*DQH(K) 
                  T(K,1)    = T(K,1)+DTH(K) 
                  Q(K,1)    = Q(K,1)+DQH(K) 
C 
C                 * REMOVES LEFT OVER SUPERSATURATION TO INSURE 
C                 * MAX(R.H.)<=1. 
C 
                  EI(K)     = ESTEFF(T(K,1))
                  SPHI(K)   = EPS1*EI(K)/(P(ICON(K),1)-EPS2*EI(K))
                  DQEX(K)   = MAX(Q(K,1) - SPHI(K), 0.E0) 
                  DQH(K)    = DQH(K) - DQEX(K)
                  Q(K,1)    = Q(K,1) - DQEX(K)
                  T(K,1)    = T(K,1) + HT(K) * DQEX(K)
C 
C                 * CONVECTIVE OR STABLE PRECIPITATION. 
C 
                  DQR(K) = DQR(K) + DEL(ICON(K),1) * DQH(K) 
              ENDIF 
  210     CONTINUE
      ENDIF 
C 
C     * UPDATE PRECIPITATION. 
C 
      DO 220 K=1,LENGATH
          PCP(ICON(K)) = PCP(ICON(K)) - PCPDH(ICON(K))*DQR(K) 
  220 CONTINUE
C 
C     ****************************************************************
C     * END OF ITERATION LOOP ****************************************
C     ****************************************************************
C 
C     * GATHER BACK TEMPERATURE AND SPECIFIC HUMIDITY.
C 
      DO 230 L=1,ILEV 
      DO 230 K=1,LENGATH
          QH(ICON(K),L) = Q(K,L)
          TH(ICON(K),L) = T(K,L)
  230 CONTINUE
C 
      IF(MAD)                                                    THEN 
C 
C         * RECALCULATE GATHERED INDICES. 
C 
          JYES = 0
          JNO  = LEN + 1
          DO 240 IL=IL1,IL2 
              IF(IADJ(IL).NE.0)                                  THEN 
                  JYES       = JYES + 1 
                  ICON(JYES) = IL 
              ELSE
                  JNO        = JNO  - 1 
                  ICON(JNO)  = IL 
              ENDIF 
  240     CONTINUE
          LENGATH = JYES
          IF(LENGATH.GT.0)                                       THEN 
C 
C             * RE-INITIALIZE ADJUSTMENT INDICATOR ARRAY AND START OVER.
C 
              DO 245 IL=IL1,IL2 
                  IADJ(IL) = 0.E0 
  245         CONTINUE
              GO TO 100 
          ENDIF 
      ENDIF 
C 
C     * CALCULATE NUMBER OF WET/DRY UNSTABLE POINTS IN THIS LATITUDE. 
C 
      DO 250 IL=IL1,IL2 
          IF(ITER(IL).NE.1)                     NUPS  = NUPS + 1
          IF(ITER(IL).NE.1.AND.IWET(IL).NE.0)   NSUPS = NSUPS + 1 
  250 CONTINUE
C 
C     * FILL IN "HOLES" IN MOISTURE FIELD (I.E. VALUES IN CARRIED 
C     * MOISTURE LEVELS HAVING QH<QMIN).
C 
      DO 390 I=IL1,IL2
          SUMQEX=0.0E0
          SUMQDEF=0.0E0 
          DO 360 L=MSG+1,ILEV 
              QEXC=MAX(QH(I,L)-QMIN,0.E0) 
              QDEF=MAX(QMIN-QH(I,L),0.E0) 
              SUMQEX=SUMQEX+DEL(I,L)*QEXC 
              SUMQDEF=SUMQDEF+DEL(I,L)*QDEF 
              QH(I,L)=MAX(QH(I,L),QMIN) 
  360     CONTINUE
          IF(SUMQEX.GT.0.E0) THEN 
              RATIO=MAX(SUMQEX-SUMQDEF,0.E0)/SUMQEX 
              DO 370 L=MSG+1,ILEV 
                  QH(I,L)=QMIN+RATIO*(QH(I,L)-QMIN) 
  370         CONTINUE
          ENDIF 
          QADD(I)=MAX(SUMQDEF-SUMQEX,0.E0)
  390 CONTINUE
  
C 
C     * OBTAIN FINAL PRECIPITATION RATE.
C 
      DO 400 IL=IL1,IL2 
          PCP(IL) = PCP(IL) + PCPDH(IL) * MAX(SUMQ(IL),0.E0)
  400 CONTINUE
  
      RETURN
      END 
