!+ Subroutine to find energy transfer coefficients for triple overlap.
!
! Method:
!	Energy transfer coefficients for upward and downward radiation
!       at the edges of the layers are calculated assuming maximal 
!       overlap of regions of the same nature and random overlap of
!       regions of a different nature.
!
! Current Owner of Code: J. M. Edwards
!
! History:
!	Version		Date			Comment
!	1.0		12-04-95		First Version under RCS
!						(J. M. Edwards)
!
! Description of Code:
!   FORTRAN 77  with extensions listed in documentation.
!
!- ---------------------------------------------------------------------
      SUBROUTINE OVERLAP_TRIPLE(N_PROFILE, N_LAYER, N_CLOUD_TOP
     &   , W_CLOUD, W_FREE, FRAC_REGION
     &   , CLOUD_OVERLAP
     &   )
!
!
!
      IMPLICIT NONE
!
!
!     INCLUDE HEADER FILES.
      INCLUDE 'incl/dimensions_field.ucf'
      INCLUDE 'incl/dimensions_fixed.pcf'
      INCLUDE 'incl/def_std_io.icf'
      INCLUDE 'incl/prec_machine.ucf'
      INCLUDE 'incl/precision.pcf'
      INCLUDE 'incl/error.pcf'
      INCLUDE 'incl/cloud_coeff_pointer.pcf'
      INCLUDE 'incl/cloud_type.pcf'
      INCLUDE 'incl/cloud_region.pcf'
!
!     DUMMY ARGUMENTS.
      INTEGER 	!, INTENT(IN)
     &     N_PROFILE
!             NUMBER OF PROFILES
     &   , N_LAYER
!             NUMBER OF LAYERS
     &   , N_CLOUD_TOP
!             TOPMOST CLOUDY LAYER
      REAL*8 	!, INTENT(IN)
     &     W_CLOUD(NPD_PROFILE, NPD_LAYER)
!             CLOUD AMOUNTS
     &   , FRAC_REGION(NPD_PROFILE, NPD_LAYER, NPD_REGION)
!             FRACTIONS OF TOTAL CLOUD AMOUNT OCCUPIED BY 
!             DIFFERENT REGIONS
!
      REAL*8 	!, INTENT(OUT)
     &     W_FREE(NPD_PROFILE, NPD_LAYER)
!             CLOUD-FREE AMOUNTS
     &   , CLOUD_OVERLAP(NPD_PROFILE, 0: NPD_LAYER, NPD_OVERLAP_COEFF)
!             COEFFICIENTS FOR TRANSFER OF ENERGY AT INTERFACE
!
!
!     LOCAL ARGUMENTS.
      INTEGER
     &     I
!             LOOP VARIABLE
     &   , J
!             LOOP VARIABLE
     &   , L
!             LOOP VARIABLE
     &   , K
!             LOOP VARIABLE
     &   , LL
!             LOOP VARIABLE
     &   , N_INDEX
!             NUMBER OF INDICES IN WHEN TEST
     &   , INDEX(NPD_PROFILE)
!             INDICES OF WHEN TEST
!
!
!     FIXED LOCAL VALUES:
      INTEGER
     &     N_REGION
!             NUMBER OF REGIONS
      REAL*8
     &     AREA_LOWER(NPD_PROFILE, NPD_REGION)
!             AREAS OF REGIONS IN LOWER LAYER
     &   , AREA_UPPER(NPD_PROFILE, NPD_REGION)
!             AREAS OF REGIONS IN LOWER LAYER
     &   , AREA_OVERLAP(NPD_PROFILE, NPD_REGION, NPD_REGION)
!             AREAS OF REGIONS IN LOWER LAYER
     &   , AREA_RANDOM(NPD_PROFILE)
!             AREAS OF RANDOM OVERLAP
!
!     SUBROUTINES CALLED:
      EXTERNAL
     &     WHENFGT
!
!     CRAY DIRECTIVES FOR THE WHOLE ROUTINE:
!     POINTS ARE NOT REPEATED IN THE INDEXING ARRAY, SO IT IS SAFE
!     TO VECTORIZE OVER INDIRECTLY ADDRESSED ARRAYS.
Cfpp$ NODEPCHK R
!
!
!
!     SET THE FREE FRACTIONS IN EACH LAYER.
      DO I=1, N_CLOUD_TOP-1
         DO L=1, N_PROFILE
            W_FREE(L, I)=1.0D+00
         ENDDO
      ENDDO
      DO I=N_CLOUD_TOP, N_LAYER
         DO L=1, N_PROFILE
            W_FREE(L, I)=1.0D+00-W_CLOUD(L, I)
         ENDDO
      ENDDO
!
!
!
!     SET NUMBER OF REGIONS FOR POSSIBLE FUTURE EXTENSION.
      N_REGION=3
!
!     WE CONSIDER EACH BOUNDARY IN TURN, COMPARING THE FRACTIONS
!     OF EACH REGION IN THE LAYERS ABOVE AND BELOW THE BOUNDARY.
!
!     INITIALIZE FOR THE LAYER ABOVE THE CLOUDS.
      DO L=1, N_PROFILE
         AREA_UPPER(L, IP_REGION_CLEAR)=1.0D+00
         AREA_UPPER(L, IP_REGION_STRAT)=0.0D+00
         AREA_UPPER(L, IP_REGION_CONV)=0.0D+00
      ENDDO
     
      DO I=N_CLOUD_TOP-1, N_LAYER
!
!        SET AREAS OF THE REGIONS IN THE LOWER LAYER.
         IF (I.LT.N_LAYER) THEN
            DO L=1, N_PROFILE
               AREA_LOWER(L, IP_REGION_CLEAR)=W_FREE(L, I+1)
               AREA_LOWER(L, IP_REGION_STRAT)=W_CLOUD(L, I+1)
     &            *FRAC_REGION(L, I+1, IP_REGION_STRAT)
               AREA_LOWER(L, IP_REGION_CONV)=W_CLOUD(L, I+1)
     &            *FRAC_REGION(L, I+1, IP_REGION_CONV)
            ENDDO
         ELSE
            DO L=1, N_PROFILE
               AREA_LOWER(L, IP_REGION_CLEAR)=1.0D+00
               AREA_LOWER(L, IP_REGION_STRAT)=0.0D+00
               AREA_LOWER(L, IP_REGION_CONV)=0.0D+00
            ENDDO
         ENDIF
!
!        SET THE AREAS OF OVERLAP BETWEEN LIKE REGIONS.
         DO K=1, N_REGION
            DO L=1, N_PROFILE
               AREA_OVERLAP(L, K, K)=MIN(AREA_LOWER(L, K)
     &            , AREA_UPPER(L, K))
            ENDDO
         ENDDO
!  
!        FIND THE AREAS OF OVERLAP BETWEEN UNLIKE REGIONS. THE OVERLAP 
!        BETWEEN UNLIKE REGIONS IS ASSUMED TO BE RANDOM. THESE AREAS
!        ARE SET EQUAL TO 0 FOR THE CASE WHERE THERE IS NO SUCH AREA
!        AND ARE RESET WHEN SUCH AN AREA IS PRESENT.
         DO K=1, N_REGION
            DO J=1, K-1
               DO L=1, N_PROFILE
                  AREA_OVERLAP(L, K, J)=0.0D+00
                  AREA_OVERLAP(L, J, K)=0.0D+00
               ENDDO
            ENDDO
         ENDDO
         DO L=1, N_PROFILE
            AREA_RANDOM(L)=1.0D+00
         ENDDO
         DO K=1, N_REGION
            DO L=1, N_PROFILE
               AREA_RANDOM(L)=AREA_RANDOM(L)-AREA_OVERLAP(L, K, K)
            ENDDO
         ENDDO
         CALL WHENFGT(N_PROFILE, AREA_RANDOM, 1 , TOL_DIV
     &      , INDEX, N_INDEX)
         DO K=1, N_REGION
            DO J=1, K-1
               DO L=1, N_INDEX
                  LL=INDEX(L)
                  AREA_OVERLAP(LL, K, J)
     &               =(AREA_UPPER(LL, K)-AREA_OVERLAP(LL, K, K))
     &               *(AREA_LOWER(LL, J)-AREA_OVERLAP(LL, J, J))
     &               /AREA_RANDOM(LL)
                  AREA_OVERLAP(LL, J, K)
     &               =(AREA_UPPER(LL, J)-AREA_OVERLAP(LL, J, J))
     &               *(AREA_LOWER(LL, K)-AREA_OVERLAP(LL, K, K))
     &               /AREA_RANDOM(LL)
               ENDDO
            ENDDO
         ENDDO
!
!        NOW PROCEED TO FIND THE ENERGY TRANSFER COEFFICIENTS
!        BETWEEN THE VARIOUS REGIONS.
!
!        COEFFICIENTS FOR THE DOWNWARD TRANSFER OF ENERGY:
!
!        TO AVOID DIVISION BY 0 WE INITIALIZE TO DEFAULT VALUES
!        AND RESET.
         DO L=1, N_PROFILE
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V11)=1.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V21)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V31)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V12)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V22)=1.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V32)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V13)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V23)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_V33)=1.0D+00
         ENDDO
!
!        TRANSFER FROM CLEAR-SKY REGION:
         CALL WHENFGT(N_PROFILE, AREA_UPPER(1, IP_REGION_CLEAR), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V11)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_CLEAR)
     &         /AREA_UPPER(LL, IP_REGION_CLEAR)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V21)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_STRAT)
     &         /AREA_UPPER(LL, IP_REGION_CLEAR)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V31)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_CONV)
     &         /AREA_UPPER(LL, IP_REGION_CLEAR)
         ENDDO
!
!        TRANSFER FROM STARTIFORM REGION:
         CALL WHENFGT(N_PROFILE, AREA_UPPER(1, IP_REGION_STRAT), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V12)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_CLEAR)
     &         /AREA_UPPER(LL, IP_REGION_STRAT)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V22)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_STRAT)
     &         /AREA_UPPER(LL, IP_REGION_STRAT)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V32)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_CONV)
     &         /AREA_UPPER(LL, IP_REGION_STRAT)
         ENDDO
!
!        TRANSFER FROM CONVECTIVE REGION:
         CALL WHENFGT(N_PROFILE, AREA_UPPER(1, IP_REGION_CONV), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V13)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_CLEAR)
     &         /AREA_UPPER(LL, IP_REGION_CONV)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V23)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_STRAT)
     &         /AREA_UPPER(LL, IP_REGION_CONV)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_V33)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_CONV)
     &         /AREA_UPPER(LL, IP_REGION_CONV)
         ENDDO
!
!
!        TRANSFER COEFFICIENTS FOR UPWARD FLOW OF ENERGY:
!
!        TO AVOID DIVISION BY 0 WE INITIALIZE TO DEFAULT VALUES
!        AND RESET.
         DO L=1, N_PROFILE
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U11)=1.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U21)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U31)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U12)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U22)=1.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U32)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U13)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U23)=0.0D+00
            CLOUD_OVERLAP(L, I, IP_CLOVLP_U33)=1.0D+00
         ENDDO
!
!        TRANSFER FROM CLEAR-SKY REGION:
         CALL WHENFGT(N_PROFILE, AREA_LOWER(1, IP_REGION_CLEAR), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U11)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_CLEAR)
     &         /AREA_LOWER(LL, IP_REGION_CLEAR)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U21)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_CLEAR)
     &         /AREA_LOWER(LL, IP_REGION_CLEAR)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U31)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_CLEAR)
     &         /AREA_LOWER(LL, IP_REGION_CLEAR)
         ENDDO
!
!        TRANSFER FROM STRATIFORM REGION:
         CALL WHENFGT(N_PROFILE, AREA_LOWER(1, IP_REGION_STRAT), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U12)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_STRAT)
     &         /AREA_LOWER(LL, IP_REGION_STRAT)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U22)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_STRAT)
     &         /AREA_LOWER(LL, IP_REGION_STRAT)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U32)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_STRAT)
     &         /AREA_LOWER(LL, IP_REGION_STRAT)
         ENDDO
!
!        TRANSFER FROM CONVECTIVE REGION:
         CALL WHENFGT(N_PROFILE, AREA_LOWER(1, IP_REGION_CONV), 1
     &      , TOL_DIV, INDEX, N_INDEX)
         DO L=1, N_INDEX
            LL=INDEX(L)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U13)
     &         =AREA_OVERLAP(LL, IP_REGION_CLEAR, IP_REGION_CONV)
     &         /AREA_LOWER(LL, IP_REGION_CONV)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U23)
     &         =AREA_OVERLAP(LL, IP_REGION_STRAT, IP_REGION_CONV)
     &         /AREA_LOWER(LL, IP_REGION_CONV)
            CLOUD_OVERLAP(LL, I, IP_CLOVLP_U33)
     &         =AREA_OVERLAP(LL, IP_REGION_CONV, IP_REGION_CONV)
     &         /AREA_LOWER(LL, IP_REGION_CONV)
         ENDDO
!
!        REASSIGN THE FRACTIONS IN THE UPPER LAYER TO STEP DOWN
!        THROUGH THE ATMOSPHERE.
         IF (I.LT.N_LAYER) THEN
            DO K=1, N_REGION
               DO L=1, N_PROFILE
                  AREA_UPPER(L, K)=AREA_LOWER(L, K)
               ENDDO
            ENDDO
         ENDIF
!
      ENDDO
!
!
!
      RETURN
      END
