      SUBROUTINE DISTRIB (C, H, C BAR, H BAR, ILG, IL1, IL2,
     1                    I TEST, C FCT, METHOD, FRACTN)
  
      IMPLICIT REAL (A-H,O-Z),
     +INTEGER (I-N)
      DIMENSION H(ILG), C(ILG), C FCT(ILG)
C 
C     * M.LAZARE. - MAY 30/88. - SPEED-UP METHODE 2.
C     * J.P. BLANCHET. - NOV 13/87. - MODIFIED FOR TOPO FIX.
C 
C     * NOTE THAT THE WORK ARRAY "C FCT" IS USED TWICE IN THIS ROUTINE, 
C     * ONCE IN THE SECOND PASS AND ONCE IN METHODE 2.
C---------------------------------------------------------------------- 
      NLEN=IL2-IL1+1
  
C........ LIMITING CASE OF C BAR = 1
  
      IF (C BAR .GT. 0.99E0*HBAR) THEN
            DO 100 I = IL1, IL2 
                  C(I) = MIN (C BAR,0.999E0)
C                 WRITE (*,600) I, H(I), C(I) 
  100       CONTINUE
            RETURN
      END IF
  
C     IF (C BAR .LT. 0.01) METHOD = 2 
  
C........ RESCALE H(I) IN RANGE (0,1) (OPTIMIZED BUNCHING) (OPTIONAL).
  
C     H MX = 0. 
C     H MN = 1. 
C     DO 160 I = IL1, IL2 
C         H MX = MAX (H(I), H MX) 
C         H MN = MIN (H(I), H MN) 
C 160 CONTINUE
  
C     D H = H MX - H MN 
C     DO 170 I = IL1, IL2 
C         H(I) = (H(I) - H MN) / D H
C 170 CONTINUE
  
C........ FUNCTION CONTROLING THE DEGRE OF BUNCHING (OPTIONAL). 
  
C     DO 200 I = IL1, IL2 
C         H(I) = H(I) ** 3
C         H(I) = H(I) * H(I)
C         H(I) = H(I) 
C         H(I) = H(I) ** 0.5
C 200 CONTINUE
  
C........ FIND MAX RELATIVE HUMIDITIES ALONG THE LONGITUDE CIRCLE.
  
      H MIN = 1.E+10
      H MAX =-1.E+10
      DO 210 I = IL1, IL2 
            H MAX = MAX (H MAX, H(I)) 
            H MIN = MIN (H MIN, H(I)) 
  210 CONTINUE
  
C........ COMPUTE PARAMETERS (CASE A=0) 
  
      SUM 1 = FLOAT(IL2 - IL1 + 1)
      SUM 1 I = 1.E0 / SUM 1
      IF (METHOD .GE. 3) GO TO 280
  215 SUM 2 = 0.E0
      DO 220 I = IL1, IL2 
            SUM 2 = SUM 2 + H(I)
  220 CONTINUE
  
C........ CASE C(H MAX) = 1.
  
C     HO = (SUM 2 - C BAR * SUM 1 * H MAX) / (SUM 1 * (1. - C BAR)) 
C     B = 1. / (H MAX - HO) 
  
C........ CASE C(H MAX) = H MAX  (UPPER PHYSICAL LIMIT: H BAR(CLEAR)=0.)
  
      HO = H MAX * (SUM 2 - C BAR * SUM 1) / (SUM 1 * (H MAX - C BAR))
      B  = H MAX / (H MAX - HO) 
  
C........ DISTRIBUTE CLOUDS ALONG THE LONGITUDE CIRCLE. FIRST GUESS.
  
      DO 230 I = IL1, IL2 
            C(I) =  B * (H(I) - HO) 
  230 CONTINUE
  
C........ 2ND PASS: SUM POSITIVE AND NEGATIVE CLOUDS AND RESET LATER TO ZERO. 
C         NOTE: THE FUNCTION "C FCT" ALSO CONTROLES THE DEGRE OF BUNCHING.
  
      SUM CP = 0.E0 
      SUM CN = 0.E0 
      DO 240 I = IL1, IL2 
          C FCT(I) = C(I) 
C         C FCT(I) = C(I) * (1. - C(I)) 
          SUM CP = MERGE(SUM CP + C FCT(I), SUM CP, C(I) .GE. 0.E0)
          SUM CN = MERGE(SUM CN - C(I), SUM CN, C(I) .LT. 0.E0)
          C(I)   = MAX (C(I), 0.E0) 
  240 CONTINUE
  
C........ REDISTRIBUTE NEGATIVE PORTION. (TAKE 1 OF 2 OPTIONS)
C . . . . METHODE 1 : DEPLETE MAXIMA (SPREADS THE CLOUDS).
  
      IF (METHOD .EQ. 1) FRACTN = 1.E0
      IF (METHOD .EQ. 2) FRACTN = 0.E0
  
      IF (METHOD .NE. 2) THEN 
      RATIO C = FRACTN * SUM CN / SUM CP
      DO 250 I = IL1, IL2 
            C(I) = C(I) * (1.E0 - RATIO C)
  250 CONTINUE
      END IF
  
C . . . . METHODE 2 : DEPLETE MINIMA (CONFINEMENTS OF CLOUDS).
  
      IF (METHOD .NE. 1) THEN 
          RESIDU = (1.E0 - FRACTN) * SUM CN 
  260     CONTINUE
          DO 270 I = IL1, IL2 
              C FCT(I)=MERGE(C(I),1.E+10,C(I).GT.0.001E0) 
  270     CONTINUE
          I MIN = ISMIN(NLEN,C FCT,1) 
          C MIN = C(I MIN)
          C(I MIN) = MAX (C MIN - RESIDU, 0.E0) 
          RESIDU = RESIDU - C MIN 
          IF (RESIDU .GT. 0.E0) GO TO 260 
      END IF
      GO TO 299 
  
C . . . . METHOD 3: MORE ACCURATE ITERATION METHOD. 
  
  280 CONTINUE
      ACCURCY = 0.001E0 
      ITER = 0
      H1 = 0.E0 
      H2 = MIN (H MAX, 0.999999E0)
  282 HP = (H1 + H2) * 0.5E0
      ITER = ITER + 1 
      IF (ITER .GT. 20) THEN
          METHOD = 2
          GO TO 215 
      END IF
      HP M 1 I = SUM 1 I / (1.E0 - HP)
      C PRIME = 0.E0
      DO 290 I = IL1, IL2 
          C(I) = MAX ((H(I) - HP) * HP M 1 I, 0.E0) 
          C PRIME = C PRIME + C(I)
  290 CONTINUE
      IF (ABS (C PRIME - C BAR) .LT. ACCURCY) GO TO 299 
      IF (C PRIME .GT. C BAR) THEN
          H1 = HP 
      ELSE
          H2 = HP 
      END IF
      GO TO 282 
  299 CONTINUE
  
C........ FOR TEST ONLY.
      IF (I TEST .EQ. 0) RETURN 
  
      SUM CLD = 0.E0
      DO 300 I = IL1, IL2 
C           WRITE (*,600) I, H(I), C(I) 
            SUM CLD = SUM CLD + C(I)
            IF(C(I).GT.1.001E0) CALL ABORT
  300 CONTINUE
      C BAR CM = SUM CLD / SUM 1
      IF (ABS(C BAR CM - C BAR) .GT. 0.01E0) THEN 
          WRITE (*,601) HO, C BAR, C BAR CM 
C         DO 400 I = IL1, IL2 
C             WRITE (*,600) I, H(I), C(I) 
C 400     CONTINUE
      END IF
  
      IF (H MAX .GT. 1.E0) WRITE (*,602) H MAX
      IF (H MIN .LT. 0.E0) WRITE (*,603) H MIN
  
      RETURN
C---------------------------------------------------------------------- 
  
  600 FORMAT (5X,'I, H(I), C(I) = ',I5,2F10.4)
  601 FORMAT (/5X,'<<< BUNCH >>> HUMIDITY HO = ',F10.2,5X,
     *    'C BAR PRESCRIBED = ',F10.4,5X,'C BAR COMPUTED =   ',F10.4) 
  602 FORMAT (5X,'<<< BUNCH >>> WARNING! THE MAXIMUM RELATIVE', 
     1 ' HUMIDITY IS ',F10.4) 
  603 FORMAT (5X,'<<< BUNCH >>> WARNING! THE MINIMUM RELATIVE', 
     2 ' HUMIDITY IS ',F10.4) 
  
      END 
