!+ Subroutine to calculate transmission and reflection coefficients.
!
! Method:
!	 Straightforward.
!
! 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 TRANS_SOURCE_COEFF(N_PROFILE
     &   , I_LAYER_FIRST, I_LAYER_LAST
     &   , ISOLIR, L_IR_SOURCE_QUAD
     &   , TAU, SUM, DIFF, LAMBDA, SEC_0
     &   , GAMMA_UP, GAMMA_DOWN
     &   , TRANS, REFLECT, TRANS_0, SOURCE_COEFF
     &   )
!
!
      IMPLICIT NONE
!
!
!     INCLUDE HEADER FILES.
      INCLUDE 'incl/dimensions_field.ucf'
      INCLUDE 'incl/dimensions_fixed.pcf'
      INCLUDE 'incl/spectral_region.pcf'
      INCLUDE 'incl/prec_machine.ucf'
      INCLUDE 'incl/source_coeff_pointer.pcf'
!
!     DUMMY VARIABLES.
      INTEGER	!, INTENT(IN)
     &     N_PROFILE
!             NUMBER OF PROFILES
     &   , I_LAYER_FIRST
!             FIRST LAYER TO CONSIDER
     &   , I_LAYER_LAST
!             LAST LAYER TO CONSIDER
!
!     ALGORITHMIC CONTROL
      LOGICAL	!, INTENT(IN)
     &     L_IR_SOURCE_QUAD
!             QUADRATIC SOURCE IN INFRA-RED
      INTEGER	!, INTENT(IN)
     &     ISOLIR
!             SPECTRAL REGION
!
!     OPTICAL PROPERTIES OF THE LAYER
      REAL*8	!, INTENT(IN)
     &     TAU(NPD_PROFILE, NPD_LAYER)
!             OPTICAL DEPTHS OF LAYERS
     &   , SUM(NPD_PROFILE, NPD_LAYER)
!             SUM OF ALPHA_1 AND ALPHA_2
     &   , DIFF(NPD_PROFILE, NPD_LAYER)
!             DIFFERENCE OF ALPHA_1 AND ALPHA_2
     &   , LAMBDA(NPD_PROFILE, NPD_LAYER)
!             LAMBDA
     &   , SEC_0(NPD_PROFILE)
!             SECANT OF SOLAR ZENITH ANGLE
     &   , GAMMA_UP(NPD_PROFILE, NPD_LAYER)
!             BASIC SOLAR COEFFICIENT FOR UPWARD RADIATION
     &   , GAMMA_DOWN(NPD_PROFILE, NPD_LAYER)
!             BASIC SOLAR COEFFICIENT FOR DOWNWARD RADIATION
!
!     TRANSMISSION AND REFLECTION COEFFICIENTS AND COEFFICIENTS FOR
!     SOURCE TERMS.
      REAL*8	!, INTENT(OUT)
     &     TRANS(NPD_PROFILE, NPD_LAYER)
!             DIFFUSE TRANSMISSION COEFFICIENT
     &   , REFLECT(NPD_PROFILE, NPD_LAYER)
!             DIFFUSE REFLECTION COEFFICIENT
     &   , TRANS_0(NPD_PROFILE, NPD_LAYER)
!             DIRECT TRANSMISSION COEFFICIENT
     &   , SOURCE_COEFF(NPD_PROFILE, NPD_LAYER, NPD_SOURCE_COEFF)
!             SOURCE COEFFICIENTS
!
!
!     LOCAL VARIABLES
      INTEGER
     &     I
!             LOOP VARIABLE
     &   , L
!             LOOP VARIABLE
      REAL*8
     &     GAMMA
!             GAMMA
     &   , EXPONENTIAL
!             EXPONENTIAL OF SCALED OPTICAL DEPTH
!
!
!
!
!     DETERMINE THE DIFFUSE TRANSMISSION AND REFLECTION COEFFICIENTS.
!
      DO I=I_LAYER_FIRST, I_LAYER_LAST
         DO L=1, N_PROFILE
c            print*,i,l,tau(l,i)
            EXPONENTIAL=EXP(-LAMBDA(L, I)*TAU(L, I))
            GAMMA=(SUM(L, I)-LAMBDA(L, I))
     &         /(SUM(L, I)+LAMBDA(L, I))
            TRANS(L, I)=(EXPONENTIAL*(1.0D+00-GAMMA**2)
     &         /(1.0D+00-(EXPONENTIAL*GAMMA)**2))
            REFLECT(L, I)=GAMMA*(1.0D+00-EXPONENTIAL**2)
     &         /(1.0D+00-(EXPONENTIAL*GAMMA)**2)
         ENDDO
      ENDDO
!
!
!
      IF (ISOLIR.EQ.IP_SOLAR) THEN
!
!        CALCULATE THE DIRECT TRANSMISSION AND THE SOURCE COEFFICIENTS 
!        FOR THE SOLAR BEAM: IN THE SOLAR CASE THESE ARE 
!        THE COEFFICIENTS WHICH WILL MULTIPLY THE DIRECT FLUX AT THE
!        TOP OF THE LAYER TO GIVE THE SOURCE TERMS FOR THE UPWARD 
!        DIFFUSE FLUX AND THE TOTAL DOWNWARD FLUX.
!
         DO I=I_LAYER_FIRST, I_LAYER_LAST
            DO L=1, N_PROFILE
            
               TRANS_0(L, I)=EXP(-TAU(L, I)*SEC_0(L))
               SOURCE_COEFF(L, I, IP_SCF_SOLAR_UP)
     &            =(GAMMA_UP(L, I)-REFLECT(L, I)
     &            *(1.0D+00+GAMMA_DOWN(L, I)))
     &            -GAMMA_UP(L, I)*TRANS(L, I)*TRANS_0(L, I)
               SOURCE_COEFF(L, I, IP_SCF_SOLAR_DOWN)=TRANS_0(L, I)
     &            *(1.0D+00+GAMMA_DOWN(L, I)
     &            -GAMMA_UP(L, I)*REFLECT(L, I))
     &            -(1.0D+00+GAMMA_DOWN(L, I))*TRANS(L, I)
            ENDDO
         ENDDO
!
!
      ELSE IF (ISOLIR.EQ.IP_INFRA_RED) THEN
!
!        IN THE CASE OF INFRA-RED RADIATION, THE FIRST SOURCE
!        COEFFICIENT HOLDS THE MULTIPLIER FOR THE FIRST DIFFERENCE
!        OF THE PLANCK FUNCTION ACROSS THE LAYER, AND THE SECOND
!        THAT FOR THE SECOND DIFFERENCE.
!
         DO I=I_LAYER_FIRST, I_LAYER_LAST
            DO L=1, N_PROFILE
!
!              A TOLERANCE IS ADDED TO THE NUMERATOR AND THE DENOMIATOR 
!              TO AVOID ILL-CONDITIONING AT SMALL OPTICAL DEPTHS.
!
               SOURCE_COEFF(L, I, IP_SCF_IR_1D)
     &            =((1.0D+00-TRANS(L, I)+REFLECT(L, I))
     &            +(TOL_MACHINE/(SQRT_TOL_MACHINE+TAU(L, I))))
     &            /((TOL_MACHINE/(SQRT_TOL_MACHINE+TAU(L, I)))
     &            +TAU(L, I)*SUM(L, I))
!
            ENDDO
         ENDDO
!
!
         IF (L_IR_SOURCE_QUAD) THEN
!
!           QUADRATIC CORRECTION TO SOURCE FUNCTION.
!           THIS CORRECTION IS VERY ILL-CONDITIONED FOR
!           SMALL OPTICAL DEPTHS SO THE ASYMPTOTIC FORM IS THEN USED.
!
            DO I=I_LAYER_FIRST, I_LAYER_LAST
               DO L=1, N_PROFILE
                  IF (TAU(L, I).GT.EXP(3.3D-01*LOG(TOL_MACHINE))) THEN
                     SOURCE_COEFF(L, I, IP_SCF_IR_2D)
     &                  =-2.0D+00*(1.0D+00-TRANS(L, I)-REFLECT(L, I)
     &                  +SQRT_TOL_MACHINE)
     &                  /(DIFF(L, I)*TAU(L, I)+SQRT_TOL_MACHINE)
                  ELSE
                     SOURCE_COEFF(L, I, IP_SCF_IR_2D)
     &                  =-2.0D+00+DIFF(L, I)*TAU(L, I)
                  ENDIF
                  SOURCE_COEFF(L, I, IP_SCF_IR_2D)
     &               =-(1.0D+00+REFLECT(L, I)+TRANS(L, I)
     &               +SOURCE_COEFF(L, I, IP_SCF_IR_2D))
     &               /(SUM(L, I)*TAU(L, I)+SQRT_TOL_MACHINE)
               ENDDO
            ENDDO
!
         ENDIF
!
      ENDIF
!
!
!
      RETURN
      END
