c     CERES Window 2   980--1100 cm^{-1}
      real wi(8),xi(8),p0(35),p(34),t0(35),t(34),bnu(34),u0(35),u(34),
     *uco2(34),ux(35),uo3(34),uf11(34),uf12(34),uf22(34),tau(4,34),
     *tauco2(1,34),tauo3(5,34),taucfc(1,34),tauf(4,1,5,34),
     *sx(4,1,5,34),sy(4,1,5,34),x(34,8),y(34,8),f(4),fco2(1),fo3(5),
     *fcfc(2),fac(34),tsx(35),tsy(35)
      data xi /.01986,.10167,.23723,.40828,
     *         .59172,.76277,.89833,.98014/
      data wi /.05062,.11119,.15685,.18134,
     *         .18134,.15685,.11119,.05062/
      open(unit=9,file='/usr1/aquila/kratz/ozone/mls70.dat',
     *status='old')
      read(9,*)(u0(m),ux(m),p0(m),t0(m),m=35,1,-1)
      omega=1040.0
      delw=120.0
      mlv=34
      do 5 m=1,34
        t(m)=(t0(m)+t0(m+1))/2.0
        xlog=(log10(p0(m))+log10(p0(m+1)))/2.0
        p(m)=10.0**xlog
        call flux(bnu(m),omega,t(m))
        if(m.ge.26)then
         fac(m)=5.0
        else
         fac(m)=1.0
        end if
        xlog=(log10(u0(m))+log10(u0(m+1)))/2.0
        u(m)=fac(m)*10.0**xlog
        xlog=(log10(ux(m))+log10(ux(m+1)))/2.0
        uo3(m)=fac(m)*10.0**xlog
        uco2(m)=p(m)*1.01325e+06/(1.3805e-16*t(m))*3.50e-04*fac(m)*
     *          1.0e+05/6.023e+23*44.00995
        uf11(m)=p(m)*1.01325e+06/(1.3805e-16*t(m))*0.268e-09*fac(m)*
     *          1.0e+05/6.023e+23*137.3685
        if(m.ge.26)uf11(m)=uf11(m)*exp((-7.5-5.0*(m-25))/3.0)
        if(m.le.25.and.m.ge.16)uf11(m)=uf11(m)*exp((15.5-m)/3.0)
        uf12(m)=p(m)*1.01325e+06/(1.3805e-16*t(m))*0.503e-09*fac(m)*
     *          1.0e+05/6.023e+23*120.9139
        if(m.ge.26)uf12(m)=uf12(m)*exp((-7.5-5.0*(m-25))/5.0)
        if(m.le.25.and.m.ge.16)uf12(m)=uf12(m)*exp((15.5-m)/5.0)
        uf22(m)=p(m)*1.01325e+06/(1.3805e-16*t(m))*0.105e-09*fac(m)*
     *          1.0e+05/6.023e+23*86.4689
        if(m.ge.26)uf22(m)=uf22(m)*exp((-14.5-5.0*(m-25))/18.0)
        if(m.le.25.and.m.ge.9)uf22(m)=uf22(m)*exp((8.5-m)/18.0)
 5    continue
      do 7 m=1,35
        tsx(m)=0.0
        tsy(m)=0.0
 7    continue
      call flux(back,omega,t0(1))
      call ck(ni,mlv,u,delw,f,p,t,tau)
      call ckco2(nco2,mlv,uco2,delw,fco2,p,t,tauco2)
      call cko3(no3,mlv,uo3,delw,fo3,p,t,tauo3)
      call ckcfc(ncfc,mlv,uf11,uf12,uf22,delw,fc,p,t,taucfc)
      do 10 m=1,34
        call window(u(m),twin,p(m),t(m),fac(m))
        do i=1,ni
         do ia=1,nco2
          do ib=1,no3
            tauf(i,ia,ib,m)=tau(i,m)+twin+tauco2(ia,m)+tauo3(ib,m)
     *                      +taucfc(1,m)
          end do
         end do
        end do
 10   continue
      do 90 i=1,ni
       do 89 ia=1,nco2
       do 85 ib=1,no3
        do 30 m=1,34
          sx(i,ia,ib,m)=0.0
          sy(i,ia,ib,m)=0.0
          do 20 j=1,8
           if(m.eq.1)then
            x(m,j)=back*exp(-tauf(i,ia,ib,m)/xi(j))+
     *      bnu(m)*(1.0-exp(-tauf(i,ia,ib,m)/xi(j)))
           else
            x(m,j)=x(m-1,j)*exp(-tauf(i,ia,ib,m)/xi(j))+
     *      bnu(m)*(1.0-exp(-tauf(i,ia,ib,m)/xi(j)))
           end if
 20       continue
 30     continue
        do 50 m=34,1,-1
          do 40 j=1,8
           if(m.eq.34)then
            y(m,j)=bnu(m)*(1.0-exp(-tauf(i,ia,ib,m)/xi(j)))
           else
            y(m,j)=y(m+1,j)*exp(-tauf(i,ia,ib,m)/xi(j))+
     *      bnu(m)*(1.0-exp(-tauf(i,ia,ib,m)/xi(j)))
           end if
 40       continue
 50     continue
        do 70 m=1,34
          do 60 j=1,8
            sx(i,ia,ib,m)=sx(i,ia,ib,m)+2.0*x(m,j)*wi(j)*xi(j)
            sy(i,ia,ib,m)=sy(i,ia,ib,m)+2.0*y(m,j)*wi(j)*xi(j)
 60       continue
 70     continue
        do 80 m=1,34
           tsx(m+1)=tsx(m+1)+sx(i,ia,ib,m)*delw*f(i)*fco2(ia)
     *     *fo3(ib)
           tsy(m)=tsy(m)+sy(i,ia,ib,m)*delw*f(i)*fco2(ia)
     *     *fo3(ib)
 80     continue
        tb=tb+back*delw*f(i)*fco2(ia)*fo3(ib)
        tsx(1)=tb
 85    continue
 89    continue
 90   continue
      filter=0.7595234
      do 95 m=1,35
        write(6,100)m-1,filter*tsx(m),filter*(tb-tsx(m)),filter*tsy(m)
 95   continue
 100  format(i3,3f15.7)
      stop
      end
c *********************************************************************
      subroutine flux(bnu,omega,t)
      bnu=0.0
      z=3.14159
      bnu=z/1000.0*((2.*6.625e-27)*(2.997925e+10)**2
     *    *omega**3)/(exp(1.4388*omega/t)-1.0)
      return
      end
c *********************************************************************
c     h2o 980--1100 cm^{-1}
      subroutine ck(ni,mlv,u,delw,f,p,t,tau)
      real f(4),k(34),tau(4,34),delw,p(34),t(34),u(34),coefk(4,3,19),
     *     stp(19),fkg(4,34),pmb(34)
      ni=4
      f(1)=0.914825
      f(2)=0.045491
      f(3)=0.036386
      f(4)=0.003298
      k(1)=0.00238376099
      do i=2,ni
        k(i)=10.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.00000,0.00000,
     *0.00009,0.00078,0.00452,0.01483,0.03426,0.06639,
     *0.11794,0.19817,0.32318,0.51457,0.80714,1.23416,
     *1.80137,
     *0.000E+00,0.000E+00,0.000E+00,0.000E+00,0.000E+00,0.000E+00,
     *3.087E-06,2.678E-05,1.455E-04,4.626E-04,1.070E-03,2.083E-03,
     *3.708E-03,6.252E-03,1.015E-02,1.610E-02,2.492E-02,3.792E-02,
     *5.561E-02,
     *0.000E+00,0.000E+00,0.000E+00,0.000E+00,0.000E+00,0.000E+00,
     *3.219E-08,2.975E-07,1.519E-06,4.734E-06,1.090E-05,2.135E-05,
     *3.804E-05,6.464E-05,1.046E-04,1.658E-04,2.524E-04,3.817E-04,
     *5.685E-04/
      data ( (coefk(2,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.00000,0.00000,0.00000,0.00000,0.00005,0.00057,
     *0.00374,0.01246,0.02747,0.04804,0.07634,0.11725,
     *0.17851,0.26999,0.40956,0.60833,0.85850,1.14913,
     *1.39638,
     *0.000E+00,0.000E+00,0.000E+00,0.000E+00,1.700E-06,1.995E-05,
     *1.177E-04,3.789E-04,8.418E-04,1.517E-03,2.429E-03,3.740E-03,
     *5.666E-03,8.523E-03,1.261E-02,1.852E-02,2.637E-02,3.494E-02,
     *4.316E-02,
     *0.000E+00,0.000E+00,0.000E+00,0.000E+00,1.688E-08,2.206E-07,
     *1.194E-06,3.778E-06,8.425E-06,1.556E-05,2.510E-05,3.895E-05,
     *5.897E-05,8.877E-05,1.274E-04,1.836E-04,2.658E-04,3.456E-04,
     *4.354E-04/
      data ( (coefk(3,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.02572,0.02678,0.02958,0.03599,0.04565,0.06011,
     *0.08145,0.11230,0.15518,0.21712,0.30136,0.40518,
     *0.52890,0.66382,0.79336,0.90048,0.97465,1.01230,
     *0.98289,
     *1.362E-03,1.393E-03,1.450E-03,1.638E-03,1.924E-03,2.335E-03,
     *2.941E-03,3.848E-03,5.123E-03,6.871E-03,9.322E-03,1.241E-02,
     *1.614E-02,2.032E-02,2.429E-02,2.793E-02,3.071E-02,3.224E-02,
     *3.241E-02,
     *2.040E-05,2.072E-05,2.110E-05,2.297E-05,2.579E-05,2.946E-05,
     *3.465E-05,4.291E-05,5.505E-05,7.050E-05,9.342E-05,1.237E-04,
     *1.601E-04,2.030E-04,2.433E-04,2.842E-04,3.160E-04,3.324E-04,
     *3.490E-04/
      data ( (coefk(4,jt,jp), jp = 1, 19), jt = 1, 3)/
     *2.71173,2.70906,2.69995,2.68591,2.66765,2.64146,
     *2.59903,2.53778,2.45614,2.34626,2.21089,2.05945,
     *1.88700,1.69715,1.50251,1.30286,1.10293,0.90151,
     *0.71955,
     *8.568E-02,8.552E-02,8.527E-02,8.498E-02,8.437E-02,8.368E-02,
     *8.244E-02,8.084E-02,7.844E-02,7.541E-02,7.146E-02,6.686E-02,
     *6.154E-02,5.576E-02,4.970E-02,4.319E-02,3.651E-02,3.006E-02,
     *2.363E-02,
     *8.699E-04,8.665E-04,8.643E-04,8.640E-04,8.569E-04,8.513E-04,
     *8.400E-04,8.293E-04,8.071E-04,7.826E-04,7.463E-04,7.006E-04,
     *6.477E-04,5.909E-04,5.275E-04,4.576E-04,3.863E-04,3.228E-04,
     *2.508E-04/
	data stp / 0.251, 0.398, 0.631, 1.000, 1.58, 2.51, 
     *	             3.98, 6.31, 10.0, 15.8, 25.1, 39.8, 63.1,
     *	             100.0, 158.0, 251.0, 398.0, 631.0, 1000.0/
      do i=1,ni
      do m=1,mlv
        ml=1
        pmb(m)=p(m)*1013.25
        if(pmb(m).lt.stp(1))then
         x1=coefk(i,1,1)+coefk(i,2,1)*(t(m)-250.0)
     *   +coefk(i,3,1)*(t(m)-250.0)**2
c     Line shapes assumed to remain dominated by pressure broadening 
c     even at low pressures (default).
c        fkg(i,m)=x1*pmb(m)/stp(1)
c     Line shapes assumed to be dominated by Doppler broadening at 
c     low pressures.
         fkg(i,m)=x1
        else if (pmb(m).gt.stp(19)) then
         x1=coefk(i,1,18)+coefk(i,2,18)*(t(m)-250.0)
     *   +coefk(i,3,18)*(t(m)-250.0)**2
         x2=coefk(i,1,19)+coefk(i,2,19)*(t(m)-250.0)
     *   +coefk(i,3,19)*(t(m)-250.0)**2
         fkg(i,m)=x1+(x2-x1)/(stp(19)-stp(18))
     *  *(pmb(m)-stp(18))
        else
         do while(pmb(m).ge.stp(ml))
           ml=ml+1
         end do
         x1=coefk(i,1,ml-1)+coefk(i,2,ml-1)*(t(m)-250.0)
     *   + coefk(i,3,ml-1)*(t(m)-250.0)**2
         x2=coefk(i,1,ml)+coefk(i,2,ml)*(t(m)-250.0)
     *   + coefk(i,3,ml)*(t(m)-250.0)**2
         fkg(i,m)=x1+(x2-x1)/(stp(ml)-stp(ml-1))
     *   *(pmb(m)-stp(ml-1))
        end if
      end do
      end do
      do i=1,ni
        do m=1,mlv
          tau(i,m)=k(i)*u(m)*fkg(i,m)
c         write(6,*)m,i,u(m),k(i),fkg(i,m),tau(i,m)
        end do
      end do
      return
      end
c *********************************************************************
      subroutine window (amnt,twin,patm,temp,fac)
c     Parameterization of CKD_2.1 continuum over 980 to 1100 cm-1 band
c     Fit inaccurate for amnt < 0.01 g/cm*2;
c          however optical depths < 0.001.
c     INPUT:
c     amnt = h2O layer amount (g/cm**2)
c     patm= pressure (atm)
c     temp = temperature (k)
c     OUTPUT:
c     twin = parameterized CKD_2.1 optical depth over band
c
      parameter (ncoef=7)
      real aa(ncoef)
c     data aa/8.408,1.000,-0.02319,0.02111,-0.8137,6.080e-06,0.9933/
      data aa/8.436,1.000,-0.02324,0.02087,-0.8128,1.466e-05,0.9934/
c
c
      ph2o = amnt*(8.314e+07*temp)/(fac*1.0e+05*18.01534*1.01325e+06)
c
          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 )
c
	return
	end
c *********************************************************************
c     co2 980--1100 cm^{-1}
      subroutine ckco2(ni,mlv,u,delw,f,p,t,tau)
      real f(1),k(34),tau(1,34),delw,p(34),t(34),u(34)
      ni=1
      f(1)=1.000000
      k(1)=0.03479137
      do m=1,mlv
        tau(1,m)=0.837764*k(1)*(1.0000+3.4268e-02*(t(m)-250.0)+
     *           3.7401e-04*(t(m)-250.0)**2)*u(m)
      end do
      return
      end
c *********************************************************************
c     o3 980--1100 cm^{-1}
      subroutine cko3(ni,mlv,u,delw,f,p,t,tau)
      real f(5),k(34),tau(5,34),delw,p(34),t(34),u(34),coefk(5,3,19),
     *     stp(19),fkg(5,34),pmb(34)
      ni=5
      f(1)=0.248933
      f(2)=0.457202
      f(3)=0.216992
      f(4)=0.062813
      f(5)=0.014060
      k(1)=6.71032
      do i=2,ni
        k(i)=10.0*k(i-1)
      end do
      data ( (coefk(1,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.00320,0.01284,0.03445,0.07584,0.14230,0.24621,
     *0.40450,0.64840,1.00000,1.47430,2.11962,2.97938,
     *4.03962,5.22412,6.41274,7.49490,8.35317,9.01560,
     * 9.57730,
     * 6.713E-05, 2.180E-04, 5.124E-04, 1.001E-03, 1.788E-03, 3.011E-03,
     * 4.707E-03, 7.103E-03, 1.051E-02, 1.496E-02, 2.094E-02, 2.883E-02,
     * 3.855E-02, 4.917E-02, 5.841E-02, 6.480E-02, 6.702E-02, 6.669E-02,
     * 6.519E-02,
     * 4.906E-07, 1.025E-06, 2.322E-06, 3.903E-06, 7.034E-06, 1.048E-05,
     * 1.567E-05, 2.058E-05, 2.782E-05, 4.112E-05, 5.237E-05, 6.130E-05,
     * 7.315E-05, 8.429E-05, 8.187E-05, 5.862E-05, 3.568E-05, 2.766E-06,
     * 6.475E-06/
      data ( (coefk(2,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.04389,0.06296,0.09425,0.14304,0.21227,0.30795,
     *0.45567,0.67711,1.00000,1.42635,1.95939,2.62854,
     *3.45876,4.48749,5.73926,7.15974,8.60381,9.93792,
     *10.99805,
     * 9.019E-04, 1.101E-03, 1.390E-03, 1.740E-03, 2.215E-03, 2.860E-03,
     * 3.749E-03, 4.932E-03, 6.600E-03, 8.888E-03, 1.187E-02, 1.549E-02,
     * 1.921E-02, 2.273E-02, 2.638E-02, 3.084E-02, 3.554E-02, 4.023E-02,
     * 4.487E-02,
     * 6.734E-06, 6.888E-06, 7.319E-06, 7.753E-06, 8.350E-06, 9.003E-06,
     * 9.425E-06, 9.863E-06, 1.010E-05, 1.298E-05, 1.863E-05, 2.550E-05,
     * 3.188E-05, 3.390E-05, 2.814E-05, 1.474E-05, 5.072E-06,-3.215E-06,
     *-3.543E-05/
      data ( (coefk(3,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.20985,0.23325,0.26466,0.30880,0.36614,0.44217,
     *0.56267,0.74252,1.00000,1.34751,1.80220,2.38159,
     *3.03788,3.69088,4.25143,4.66915,4.93614,5.06337,
     * 5.12396,
     * 2.802E-03, 2.790E-03, 2.812E-03, 2.893E-03, 3.082E-03, 3.363E-03,
     * 3.718E-03, 4.052E-03, 4.358E-03, 4.363E-03, 4.042E-03, 3.636E-03,
     * 3.105E-03, 2.555E-03, 1.738E-03, 6.342E-04,-5.584E-04,-1.655E-03,
     *-2.769E-03,
     * 1.499E-05, 1.422E-05, 1.358E-05, 1.303E-05, 1.228E-05, 1.200E-05,
     * 1.229E-05, 1.230E-05, 1.307E-05, 1.499E-05, 1.212E-05, 7.131E-06,
     * 3.791E-06, 6.062E-07, 2.347E-06, 3.931E-06, 5.284E-06, 6.159E-06,
     * 1.432E-05/
      data ( (coefk(4,jt,jp), jp = 1, 19), jt = 1, 3)/
     *0.49106,0.49922,0.51093,0.53265,0.56184,0.60461,
     *0.69186,0.82488,1.00000,1.18415,1.32702,1.40221,
     *1.39564,1.31679,1.19139,1.05403,0.93148,0.84829,
     * 0.78764,
     * 4.379E-03, 4.323E-03, 4.217E-03, 4.049E-03, 3.851E-03, 3.597E-03,
     * 3.167E-03, 2.562E-03, 1.776E-03, 1.053E-03, 6.011E-04, 2.164E-04,
     *-1.765E-04,-4.361E-04,-5.357E-04,-5.249E-04,-5.462E-04,-5.599E-04,
     *-5.494E-04,
     * 4.275E-06, 4.425E-06, 4.631E-06, 5.059E-06, 5.156E-06, 5.506E-06,
     * 4.609E-06, 3.256E-06, 2.106E-06, 2.719E-07,-3.094E-07, 1.147E-06,
     * 8.563E-07, 8.969E-07, 1.925E-06, 1.522E-06, 2.088E-06, 1.428E-06,
     * 6.344E-07/
      data ( (coefk(5,jt,jp), jp = 1, 19), jt = 1, 3)/
     *1.42231,1.41001,1.38601,1.35629,1.30044,1.22096,
     *1.15796,1.08959,1.00000,0.87232,0.72065,0.57595,
     *0.44530,0.33968,0.26222,0.20970,0.17386,0.14736,
     * 0.12964,
     *-2.145E-03,-2.036E-03,-1.836E-03,-1.651E-03,-1.195E-03,-6.566E-04,
     *-5.384E-04,-6.499E-04,-6.870E-04,-5.831E-04,-4.473E-04,-3.281E-04,
     *-2.534E-04,-1.961E-04,-1.825E-04,-1.681E-04,-1.485E-04,-1.250E-04,
     *-1.093E-04,
     * 2.531E-07,-4.312E-07,-1.747E-06,-2.100E-06,-4.587E-06,-6.747E-06,
     *-5.572E-06,-3.597E-06,-2.581E-06,-8.531E-07, 2.250E-07,-4.907E-07,
     * 5.531E-07, 5.031E-07, 4.500E-07, 4.406E-07, 4.125E-07, 4.875E-07,
     * 3.562E-07/
	data stp / 0.251, 0.398, 0.631, 1.000, 1.58, 2.51, 
     *	             3.98, 6.31, 10.0, 15.8, 25.1, 39.8, 63.1,
     *	             100.0, 158.0, 251.0, 398.0, 631.0, 1000.0/
      do i=1,ni
      do m=1,mlv
        ml=1
        pmb(m)=p(m)*1013.25
        if(pmb(m).lt.stp(1))then
         x1=coefk(i,1,1)+coefk(i,2,1)*(t(m)-250.0)
     *   +coefk(i,3,1)*(t(m)-250.0)**2
c     Line shapes assumed to remain dominated by pressure broadening 
c     even at low pressures (default).
c        fkg(i,m)=x1*pmb(m)/stp(1)
c     Line shapes assumed to be dominated by Doppler broadening at 
c     low pressures.
         fkg(i,m)=x1
        else if (pmb(m).gt.stp(19)) then
         x1=coefk(i,1,18)+coefk(i,2,18)*(t(m)-250.0)
     *   +coefk(i,3,18)*(t(m)-250.0)**2
         x2=coefk(i,1,19)+coefk(i,2,19)*(t(m)-250.0)
     *   +coefk(i,3,19)*(t(m)-250.0)**2
         fkg(i,m)=x1+(x2-x1)/(stp(19)-stp(18))
     *  *(pmb(m)-stp(18))
        else
         do while(pmb(m).ge.stp(ml))
           ml=ml+1
         end do
         x1=coefk(i,1,ml-1)+coefk(i,2,ml-1)*(t(m)-250.0)
     *   + coefk(i,3,ml-1)*(t(m)-250.0)**2
         x2=coefk(i,1,ml)+coefk(i,2,ml)*(t(m)-250.0)
     *   + coefk(i,3,ml)*(t(m)-250.0)**2
         fkg(i,m)=x1+(x2-x1)/(stp(ml)-stp(ml-1))
     *   *(pmb(m)-stp(ml-1))
        end if
      end do
      end do
      do i=1,ni
        do m=1,mlv
          tau(i,m)=k(i)*u(m)*fkg(i,m)
        end do
      end do
      return
      end
c *********************************************************************
      subroutine ckcfc(ni,mlv,uf11,uf12,uf22,delw,f,p,t,tau)
      real f(2),k(34),tau(1,34),delw,p(34),t(34),uf11(34),uf12(34),
     *uf22(34)
      ni=1
c     f(1)=ratio of spectral interval used vs. cfc spectral interval
c     correction for spectral interval offset
      cor11=0.891109
      cor12=0.869315
      cor22=0.902223
      f(1)=1.0
      do m=1,mlv
        fk11=cor11*891.794*(1.0000+1.9494e-04*(t(m)-250.0)+
     *       7.0829e-06*(t(m)-250.0)**2)
        fk12=cor12*800.771*(1.0000+1.3996e-03*(t(m)-250.0)+
     *       1.7760e-06*(t(m)-250.0)**2)
        fk22=cor22*417.009*(1.0000+2.8908e-03*(t(m)-250.0)+
     *       7.1388e-06*(t(m)-250.0)**2)
        tau(1,m)=fk11*uf11(m)+fk12*uf12(m)+fk22*uf22(m)
      end do
      return
      end
