module modis28

	use CKD_common

	implicit none
	private
	
	integer, parameter :: ni = 5
	integer, parameter :: nch4 = 5
	
	public :: ch28, init_ch28
	
contains
	subroutine init_ch28(modis_channel)
	
		type(ckd_type), intent(inout) :: modis_channel
		integer :: all_weights
		
		all_weights = ni*nch4
		
		modis_channel%num_weights = all_weights
		allocate(modis_channel%weights(all_weights))
		allocate(modis_channel%taus(all_weights, max_layers))
	
	end subroutine init_ch28


      subroutine ch28(u0, p0,t0, ux, fac,modis_channel)
!
!      This routine was developed for MODIS Channel 28: 1335--1395 cm^{-1}
!
!-----------------------------------------------------------------
!     Channel 28   7.175 - 7.475 microns
!-----------------------------------------------------------------
!  INPUTS:
!     u0 --> water vapor amount [g/cm2/km]
!     p0 --> pressure in atmospheres
!     t0 --> temp in K
!	  ux --> ozone (g/cm^2/km)
!     fac --> layer thickness in km
!-----------------------------------------------------------------

		real, dimension(:), intent(in) :: u0, p0, t0, ux, fac
		type(ckd_type), intent(inout) :: modis_channel
		
		integer, parameter :: nlev = max_levels
		integer, parameter :: nlay = max_layers
		integer, parameter :: mlv = nlay

		real :: p(nlay), t(nlay),u(nlay), xlog, &
				  uch4(nlay)
				
		real :: tau(ni,nlay),tauch4(nch4, nlay), &
				 f(ni),  fch4(nch4), twin
		
		integer :: i,m,ia, index
	  
      do  m=1,use_layers
        t(m)=(t0(m)+t0(m+1))/2.0
        p(m)=(p0(m)+p0(m+1))/2.0
!       u(m)=fac(m)*(u0(m)+u0(m+1))/2.0*6.023e+23/2.687e+19/18.01534
        u(m)=fac(m)*(u0(m)+u0(m+1))/2.0
        uch4(m)=p(m)*1.01325e+06/(1.3805e-16*t(m))*1.75e-06*fac(m)* &
     &       1.0e+05/6.023e+23*16.04303
      enddo
	  
      call ck_28(u,f,p,t,tau)
      call ckch4_28(uch4,fch4,p,t,tauch4)
	  	  
      do  m=1,use_layers
        call window_28(u(m),twin,p(m),t(m),fac(m))
          do ia=1,nch4 ! 5 k's for ch4
        do i=1,ni  ! 5 k's for h2o
			index = (ia-1)*ni + i		
            modis_channel%taus(index,m)=tau(i,m)+twin+tauch4(ia,m)
 			if (m == 1) modis_channel%weights(index) = f(i)*fch4(ia)
          end do
        end do
      enddo
	  
	
	end subroutine ch28

! *********************************************************************
!     h2o 1335--1395 cm^{-1}
      subroutine ck_28(u,f,p,t,tau)
		real, dimension(:), intent(in) :: p, t, u
		
		real, dimension(:), intent(inout) :: f
		real, dimension(:,:), intent(inout) :: tau

		real :: k(ni)
		real :: coefk(ni,3,num_pressures)
		integer :: i, jp, jt
      f(1)=0.526064
      f(2)=0.309565
      f(3)=0.117418
      f(4)=0.041242
      f(5)=0.005711
      k(1)=1.7591042
      do i=2,ni
        k(i)=8.0*k(i-1)
      end do
      data ( (coefk(1,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00062,0.00098,0.00154,0.00242,0.00381,0.00599, &
     & 0.00939,0.01471,0.02305,0.03606,0.05636,0.08828, &
     & 0.13929,0.22014,0.34626,0.53718,0.81621,1.22661, &
     & 1.80673, &
     & 6.500E-06,1.025E-05,1.613E-05,2.538E-05,3.963E-05,6.213E-05, &
     & 9.675E-05,1.514E-04,2.360E-04,3.688E-04,5.795E-04,9.074E-04, &
     & 1.431E-03,2.270E-03,3.581E-03,5.594E-03,8.534E-03,1.287E-02, &
     & 1.904E-02, &
     & 1.250E-08,1.875E-08,3.438E-08,5.938E-08,8.438E-08,1.281E-07, &
     & 2.000E-07,3.219E-07,4.875E-07,7.437E-07,1.212E-06,1.916E-06, &
     & 3.000E-06,4.819E-06,7.669E-06,1.271E-05,2.050E-05,3.285E-05, &
     & 4.776E-05/
      data ( (coefk(2,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00096,0.00147,0.00226,0.00347,0.00534,0.00825, &
     & 0.01276,0.01977,0.03076,0.04771,0.07370,0.11322, &
     & 0.17291,0.25854,0.37848,0.55473,0.82280,1.21548, &
     & 1.74685, &
     & 1.088E-05,1.637E-05,2.487E-05,3.762E-05,5.725E-05,8.800E-05, &
     & 1.355E-04,2.096E-04,3.231E-04,4.994E-04,7.725E-04,1.182E-03, &
     & 1.806E-03,2.692E-03,3.979E-03,5.992E-03,9.123E-03,1.369E-02, &
     & 1.967E-02, &
     & 3.438E-08,4.687E-08,5.937E-08,8.437E-08,1.312E-07,2.000E-07, &
     & 3.062E-07,4.781E-07,6.719E-07,9.906E-07,1.562E-06,2.372E-06, &
     & 3.675E-06,5.878E-06,9.709E-06,1.389E-05,1.807E-05,2.477E-05, &
     & 3.312E-05/
      data ( (coefk(3,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00196,0.00281,0.00407,0.00594,0.00870,0.01285, &
     & 0.01912,0.02843,0.04214,0.06186,0.08957,0.12728, &
     & 0.18341,0.26963,0.39774,0.58307,0.84273,1.16653, &
     & 1.47341, &
     & 2.925E-05,3.950E-05,5.437E-05,7.587E-05,1.071E-04,1.531E-04, &
     & 2.199E-04,3.173E-04,4.546E-04,6.528E-04,9.160E-04,1.298E-03, &
     & 1.931E-03,2.918E-03,4.393E-03,6.534E-03,9.599E-03,1.356E-02, &
     & 1.745E-02, &
     & 1.562E-07,1.812E-07,2.156E-07,2.531E-07,3.281E-07,4.469E-07, &
     & 5.531E-07,7.375E-07,9.406E-07,1.344E-06,1.763E-06,2.972E-06, &
     & 4.247E-06,5.406E-06,7.066E-06,1.065E-05,1.552E-05,2.320E-05, &
     & 3.507E-05/
      data ( (coefk(4,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.02474,0.02716,0.03081,0.03609,0.04349,0.05438, &
     & 0.06979,0.09164,0.12349,0.17027,0.23742,0.33297, &
     & 0.45910,0.60364,0.75125,0.88178,0.97072,1.02249, &
     & 1.00497, &
     & 4.100E-04,4.280E-04,4.606E-04,5.182E-04,6.022E-04,7.234E-04, &
     & 9.085E-04,1.172E-03,1.547E-03,2.079E-03,2.894E-03,4.029E-03, &
     & 5.545E-03,7.384E-03,9.206E-03,1.095E-02,1.171E-02,1.222E-02, &
     & 1.193E-02, &
     & 3.275E-06,3.306E-06,3.422E-06,3.612E-06,3.956E-06,4.078E-06, &
     & 4.463E-06,5.184E-06,5.691E-06,5.903E-06,7.572E-06,9.675E-06, &
     & 1.257E-05,1.655E-05,2.163E-05,2.821E-05,2.849E-05,2.838E-05, &
     & 2.990E-05/
      data ( (coefk(5,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 2.47130,2.46377,2.45640,2.44544,2.42681,2.40072, &
     & 2.36233,2.31177,2.24881,2.16870,2.07667,1.96404, &
     & 1.82806,1.66439,1.48253,1.29109,1.09475,0.90080, &
     & 0.71857, &
     & 2.732E-02,2.725E-02,2.721E-02,2.712E-02,2.697E-02,2.681E-02, &
     & 2.645E-02,2.605E-02,2.544E-02,2.467E-02,2.370E-02,2.244E-02, &
     & 2.090E-02,1.900E-02,1.687E-02,1.450E-02,1.260E-02,1.029E-02, &
     & 7.818E-03, &
     & 4.491E-05,4.539E-05,4.539E-05,4.364E-05,4.401E-05,4.425E-05, &
     & 4.409E-05,4.453E-05,4.300E-05,4.284E-05,4.264E-05,4.075E-05, &
     & 3.800E-05,3.592E-05,3.070E-05,2.342E-05,2.162E-05,1.951E-05, &
     & 1.032E-05/

		call apply_ckd(ni, k, u, f, p, t, tau, coefk, .false.)

      end subroutine ck_28
! *********************************************************************
!     ch4 1335--1395 cm^{-1}
      subroutine ckch4_28(u,f,p,t,tau)
		real, dimension(:), intent(in) :: p, t, u
		
		real, dimension(:), intent(inout) :: f
		real, dimension(:,:), intent(inout) :: tau

		real :: k(nch4)
		real :: coefk(nch4,3,num_pressures)
		integer :: i, jp, jt
      f(1)=0.426433
      f(2)=0.340613
      f(3)=0.166590
      f(4)=0.065309
      f(5)=0.001055
      k(1)=15.6509586
      do i=2,nch4
        k(i)=8.0*k(i-1)
      end do
      data ( (coefk(1,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00000,0.00000,0.00000,0.00000,0.00008,0.00029, &
     & 0.00082,0.00193,0.00425,0.00898,0.01875,0.04152, &
     & 0.08782,0.16846,0.31236,0.52235,0.80696,1.21861, &
     & 1.73129, &
     &  0.000E+00, 0.000E+00, 0.000E+00, 0.000E+00, 1.250E-06, 4.750E-06, &
     &  1.363E-05, 3.175E-05, 6.413E-05, 1.294E-04, 2.410E-04, 4.560E-04, &
     &  8.580E-04, 1.494E-03, 2.569E-03, 4.527E-03, 7.177E-03, 1.153E-02, &
     &  1.783E-02, &
     &  0.000E+00, 0.000E+00, 0.000E+00, 0.000E+00, 1.137E-15, 1.875E-08, &
     &  5.938E-08, 1.688E-07, 3.281E-07, 6.781E-07, 1.225E-06, 1.594E-06, &
     &  2.006E-06, 4.303E-06, 1.306E-06,-1.872E-06,-3.906E-06,-1.012E-05, &
     & -9.453E-06/
      data ( (coefk(2,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00027,0.00029,0.00033,0.00059,0.00108,0.00212, &
     & 0.00423,0.00838,0.01613,0.02932,0.05124,0.08652, &
     & 0.14456,0.23771,0.37767,0.57391,0.83532,1.19534, &
     & 1.66934, &
     &  1.062E-05, 1.138E-05, 1.200E-05, 1.562E-05, 2.362E-05, 3.738E-05, &
     &  6.275E-05, 1.058E-04, 1.706E-04, 2.501E-04, 3.873E-04, 5.951E-04, &
     &  9.106E-04, 1.358E-03, 1.935E-03, 2.680E-03, 3.623E-03, 4.858E-03, &
     &  6.160E-03, &
     &  1.281E-07, 1.406E-07, 1.500E-07, 1.469E-07, 1.969E-07, 2.594E-07, &
     &  3.813E-07, 5.125E-07, 6.094E-07, 8.094E-07, 9.125E-07, 1.116E-06, &
     &  1.478E-06, 2.119E-06, 2.588E-06, 1.216E-06, 5.719E-07,-1.653E-06, &
     & -6.644E-06/
      data ( (coefk(3,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.00245,0.00255,0.00279,0.00408,0.00655,0.01082, &
     & 0.01690,0.02530,0.03707,0.05464,0.08205,0.12517, &
     & 0.19330,0.29851,0.44980,0.64644,0.87803,1.11719, &
     & 1.33794, &
     &  6.575E-05, 6.825E-05, 6.750E-05, 7.600E-05, 9.987E-05, 1.271E-04, &
     &  1.564E-04, 1.872E-04, 2.242E-04, 2.740E-04, 3.414E-04, 4.508E-04, &
     &  6.183E-04, 8.915E-04, 1.301E-03, 1.913E-03, 2.686E-03, 3.347E-03, &
     &  3.775E-03, &
     &  5.750E-07, 6.125E-07, 6.500E-07, 5.375E-07, 6.094E-07, 6.094E-07, &
     &  6.469E-07, 7.313E-07, 8.750E-07, 9.188E-07, 7.344E-07, 3.000E-07, &
     & -6.187E-07,-1.894E-06,-4.491E-06,-6.503E-06,-8.113E-06,-7.216E-06, &
     & -9.734E-06/
      data ( (coefk(4,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 0.47004,0.47577,0.48301,0.49214,0.50136,0.51167, &
     & 0.52363,0.54155,0.57086,0.61648,0.68084,0.76979, &
     & 0.87775,0.98090,1.04569,1.05977,1.02796,0.96499, &
     & 0.87486, &
     &  3.241E-03, 3.261E-03, 3.252E-03, 3.223E-03, 3.161E-03, 3.051E-03, &
     &  2.861E-03, 2.616E-03, 2.324E-03, 2.063E-03, 1.916E-03, 1.867E-03, &
     &  1.926E-03, 2.047E-03, 2.058E-03, 1.941E-03, 1.730E-03, 1.466E-03, &
     &  1.213E-03, &
     &  3.438E-08,-1.562E-08, 9.313E-12,-5.969E-07,-5.719E-07,-6.656E-07, &
     & -8.813E-07,-6.219E-07,-6.031E-07,-1.522E-06,-2.325E-06,-3.216E-06, &
     & -4.716E-06,-6.494E-06,-6.875E-06,-6.991E-06,-7.169E-06,-8.038E-06, &
     & -6.550E-06/
      data ( (coefk(5,jt,jp), jp = 1, 19), jt = 1, 3)/ &
     & 8.90708,8.95982,9.00443,9.01234,8.93422,8.74311, &
     & 8.39020,7.87647,7.20785,6.40016,5.57208,4.67567, &
     & 3.78880,2.98249,2.26846,1.65120,1.18555,0.83498, &
     & 0.58795, &
     & -4.827E-03,-4.787E-03,-5.012E-03,-5.177E-03,-5.327E-03,-5.256E-03, &
     & -5.075E-03,-4.295E-03,-2.663E-03, 1.191E-04, 2.788E-03, 5.592E-03, &
     &  6.374E-03, 5.920E-03, 5.053E-03, 4.087E-03, 2.991E-03, 2.124E-03, &
     &  1.317E-03, &
     & -7.298E-05,-7.411E-05,-7.494E-05,-8.078E-05,-7.758E-05,-7.528E-05, &
     & -7.372E-05,-6.625E-05,-6.545E-05,-6.135E-05,-6.455E-05,-6.573E-05, &
     & -4.441E-05,-2.553E-05,-2.307E-05,-1.565E-05,-1.322E-05,-9.872E-06, &
     & -5.906E-06/
		call apply_ckd(nch4, k, u, f, p, t, tau, coefk, .false.)

      end subroutine ckch4_28
! *********************************************************************
      subroutine window_28 (amnt,twin,patm,temp,fac)
!     Parameterization of CKD_2.4 continuum over 1335 to 1395 cm-1 band
!     Fit inaccurate for amnt < 0.01 g/cm*2;
!          however optical depths < 0.001.
!     INPUT:
!     amnt = h2O layer amount (g/cm**2)
!     patm= pressure (atm)
!     temp = temperature (k)
!     OUTPUT:
!     twin = parameterized CKD_2.4 optical depth over band
!
	  real, intent(in) :: amnt, patm, temp, fac
	  real, intent(inout) :: twin

      integer, parameter :: ncoef=7

      real ::  aa(ncoef), ph2o, tau_log
!
!
!     data aa/0.8065,1.019,-0.003620,1.0180,15.940,-8.057e-03,-0.03242/
      data aa/0.0000,1.200,-0.003620,1.0180,15.940,-8.057e-03,-0.03242/
!
      ph2o = amnt*(8.314e+07*temp)/(fac*1.0e+05*18.01534*1.01325e+06)
!
          tau_log   = aa(1)              +  &
     &                 aa(2) * log(amnt)  + &
     &                 aa(3) * temp       + &
     &                 aa(4) * patm       + &
     &                 aa(5) * ph2o       + &
     &                 aa(6) * amnt       + &
     &                 aa(7) * log(ph2o)

        twin = exp ( tau_log )
!
        end subroutine window_28

end module modis28
