!----------------------------------------------------------------------------
SUBROUTINE RRTM_TAUMOL16 (KLEV,P_TAU,&
 & P_TAUAERL,P_FAC00,P_FAC01,P_FAC10,P_FAC11,K_JP,K_JT,K_JT1,P_ONEMINUS,&
 & P_COLH2O,P_COLCH4,K_LAYTROP,P_SELFFAC,P_SELFFRAC,K_INDSELF,PFRAC)  

!     BAND 16:  2600-3000 cm-1 (low - H2O,CH4; high - nothing)

! Modifications
!        M.Hamrud      01-Oct-2003 CY28 Cleaning

!     D Salmond 1999-07-14 speed-up

USE PARKIND1  ,ONLY : JPIM     ,JPRB
USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK

USE PARRRTM  , ONLY : JPLAY  ,JPBAND ,JPGPT  ,NGS15
USE YOERRTWN , ONLY :      NSPA   
USE YOERRTA16, ONLY : ABSA   ,FRACREFA,SELFREF,STRRAT

IMPLICIT NONE

!  Output
INTEGER(KIND=JPIM),INTENT(IN)    :: KLEV 
REAL(KIND=JPRB)   ,INTENT(OUT)   :: P_TAU(JPGPT,JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_TAUAERL(JPLAY,JPBAND) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC00(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC01(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC10(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_FAC11(JPLAY) 
INTEGER(KIND=JPIM),INTENT(IN)    :: K_JP(JPLAY) 
INTEGER(KIND=JPIM),INTENT(IN)    :: K_JT(JPLAY) 
INTEGER(KIND=JPIM),INTENT(IN)    :: K_JT1(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_ONEMINUS 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_COLH2O(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_COLCH4(JPLAY) 
INTEGER(KIND=JPIM),INTENT(IN)    :: K_LAYTROP 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_SELFFAC(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: P_SELFFRAC(JPLAY) 
INTEGER(KIND=JPIM),INTENT(IN)    :: K_INDSELF(JPLAY) 
REAL(KIND=JPRB)   ,INTENT(OUT)   :: PFRAC(JPGPT,JPLAY) 
!- from AER
!- from INTFAC      
!- from INTIND
!- from PRECISE             
!- from PROFDATA             
!- from SELF             
!- from SP             
INTEGER(KIND=JPIM) :: IG, IND0, IND1, INDS, JS, I_LAY

REAL(KIND=JPRB) :: Z_FAC000, Z_FAC001, Z_FAC010, Z_FAC011, Z_FAC100, Z_FAC101,&
 & Z_FAC110, Z_FAC111, Z_FS, Z_SPECCOMB, Z_SPECMULT, Z_SPECPARM  
REAL(KIND=JPRB) :: ZHOOK_HANDLE

!  Input
!#include "yoeratm.h"

!      REAL TAUAER(JPLAY)
!      EQUIVALENCE (TAUAERL(1,16),TAUAER)

!     Compute the optical depth by interpolating in ln(pressure), 
!     temperature, and appropriate species.  Below LAYTROP, the water
!     vapor self-continuum is interpolated (in temperature) separately. 
 
IF (LHOOK) CALL DR_HOOK('RRTM_TAUMOL16',0,ZHOOK_HANDLE)
DO I_LAY = 1, K_LAYTROP
  Z_SPECCOMB = P_COLH2O(I_LAY) + STRRAT*P_COLCH4(I_LAY)
  Z_SPECPARM = P_COLH2O(I_LAY)/Z_SPECCOMB
  Z_SPECPARM = MIN(Z_SPECPARM,P_ONEMINUS)
  Z_SPECMULT = 8._JPRB*(Z_SPECPARM)
  JS = 1 + INT(Z_SPECMULT)
  Z_FS = MOD(Z_SPECMULT,1.0_JPRB)
!----jjm         
  Z_FAC000 = (1.0_JPRB - Z_FS) * P_FAC00(I_LAY)
  Z_FAC010 = (1.0_JPRB - Z_FS) * P_FAC10(I_LAY)
  Z_FAC100 = Z_FS * P_FAC00(I_LAY)
  Z_FAC110 = Z_FS * P_FAC10(I_LAY)
  Z_FAC001 = (1.0_JPRB - Z_FS) * P_FAC01(I_LAY)
  Z_FAC011 = (1.0_JPRB - Z_FS) * P_FAC11(I_LAY)
  Z_FAC101 = Z_FS * P_FAC01(I_LAY)
  Z_FAC111 = Z_FS * P_FAC11(I_LAY)
!-----         
  IND0 = ((K_JP(I_LAY)-1)*5+(K_JT(I_LAY)-1))*NSPA(16) + JS
  IND1 = (K_JP(I_LAY)*5+(K_JT1(I_LAY)-1))*NSPA(16) + JS
  INDS = K_INDSELF(I_LAY)
!         DO IG = 1, NG16
!-- DS_990714
  IG=1
  P_TAU (NGS15+IG,I_LAY) = Z_SPECCOMB *&
   !     &((1. - Z_FS)*(P_FAC00(I_LAY) * ABSA(IND0   ,IG) +
   !     &            P_FAC10(I_LAY) * ABSA(IND0+ 9,IG) +
   !     &            P_FAC01(I_LAY) * ABSA(IND1   ,IG) + 
   !     &            P_FAC11(I_LAY) * ABSA(IND1+ 9,IG))+
   !     &    Z_FS * (  P_FAC01(I_LAY) * ABSA(IND1+ 1,IG) +
   !     &            P_FAC10(I_LAY) * ABSA(IND0+10,IG) +
   !     &            P_FAC00(I_LAY) * ABSA(IND0+ 1,IG) +
   !     &            P_FAC11(I_LAY) * ABSA(IND1+10,IG))) +
   & (Z_FAC000 * ABSA(IND0   ,IG) +&
   & Z_FAC100 * ABSA(IND0+ 1,IG) +&
   & Z_FAC010 * ABSA(IND0+ 9,IG) +&
   & Z_FAC110 * ABSA(IND0+10,IG) +&
   & Z_FAC001 * ABSA(IND1   ,IG) +&
   & Z_FAC101 * ABSA(IND1+ 1,IG) +&
   & Z_FAC011 * ABSA(IND1+ 9,IG) +&
   & Z_FAC111 * ABSA(IND1+10,IG))+&
   & P_COLH2O(I_LAY) * &
   & P_SELFFAC(I_LAY) * (SELFREF(INDS,IG) + &
   & P_SELFFRAC(I_LAY) *&
   & (SELFREF(INDS+1,IG) - SELFREF(INDS,IG)))&
   & + P_TAUAERL(I_LAY,16)  
  PFRAC(NGS15+IG,I_LAY) = FRACREFA(IG,JS) + Z_FS *&
   & (FRACREFA(IG,JS+1) - FRACREFA(IG,JS))  
  IG=2
  P_TAU (NGS15+IG,I_LAY) = Z_SPECCOMB *&
   !     &((1. - Z_FS)*(P_FAC00(I_LAY) * ABSA(IND0   ,IG) +
   !     &            P_FAC10(I_LAY) * ABSA(IND0+ 9,IG) +
   !     &            P_FAC01(I_LAY) * ABSA(IND1   ,IG) +
   !     &            P_FAC11(I_LAY) * ABSA(IND1+ 9,IG))+
   !     &    Z_FS * (  P_FAC01(I_LAY) * ABSA(IND1+ 1,IG) +
   !     &            P_FAC10(I_LAY) * ABSA(IND0+10,IG) +
   !     &            P_FAC00(I_LAY) * ABSA(IND0+ 1,IG) +
   !     &            P_FAC11(I_LAY) * ABSA(IND1+10,IG))) +
   & (Z_FAC000 * ABSA(IND0   ,IG) +&
   & Z_FAC100 * ABSA(IND0+ 1,IG) +&
   & Z_FAC010 * ABSA(IND0+ 9,IG) +&
   & Z_FAC110 * ABSA(IND0+10,IG) +&
   & Z_FAC001 * ABSA(IND1   ,IG) +&
   & Z_FAC101 * ABSA(IND1+ 1,IG) +&
   & Z_FAC011 * ABSA(IND1+ 9,IG) +&
   & Z_FAC111 * ABSA(IND1+10,IG))+&
   & P_COLH2O(I_LAY) *&
   & P_SELFFAC(I_LAY) * (SELFREF(INDS,IG) +&
   & P_SELFFRAC(I_LAY) *&
   & (SELFREF(INDS+1,IG) - SELFREF(INDS,IG)))&
   & + P_TAUAERL(I_LAY,16)  
  PFRAC(NGS15+IG,I_LAY) = FRACREFA(IG,JS) + Z_FS *&
   & (FRACREFA(IG,JS+1) - FRACREFA(IG,JS))  

!         END DO
!-- DS_990714
ENDDO

DO I_LAY = K_LAYTROP+1, KLEV
!         DO IG = 1, NG16
!-- DS_990714
  IG=1
  P_TAU (NGS15+IG,I_LAY) = P_TAUAERL(I_LAY,16)
  PFRAC(NGS15+IG,I_LAY) = 0.0_JPRB
  IG=2
  P_TAU (NGS15+IG,I_LAY) = P_TAUAERL(I_LAY,16)
  PFRAC(NGS15+IG,I_LAY) = 0.0_JPRB
!-- DS_990714
!         END DO
ENDDO

IF (LHOOK) CALL DR_HOOK('RRTM_TAUMOL16',1,ZHOOK_HANDLE)
END SUBROUTINE RRTM_TAUMOL16