
      SUBROUTINE hcorat(ear, ne, param, ifl, photar, photer)

      INTEGER ne, ifl
      REAL ear(0:ne), param(3), photar(ne), photer(ne)

c  based on hcoabspm 
c  coded by K. Hayashida 2008.5.30
c  in turn based on varabspm only HCO is considered
c
c  works out the redshifted transmission for material whose abundances 
c  are specified for the following params :
c    1   Carbon column  (1e18 cm2)
c    2   H/C ratio   (dimensionless)
c    3   O/C ratio   (dimensionless)

c Arguments :
c      ear     r        i: energy ranges
c      ne      i        i: number of energies
c      param   r        i: model parameters
c      ifl     i        i: file number
c      photar  r        o: transmitted fraction

      INTEGER NELTS
      PARAMETER (NELTS=18)
      REAL LYLIMIT
      PARAMETER (LYLIMIT=0.0135984)

      INTEGER atomic(NELTS)
      INTEGER ie, status, i

      REAL column(NELTS)
      REAL zfac, elow, ehi, xsect, pret

      CHARACTER*2 celts(NELTS)

c External references to cross section function

      REAL photo, gphoto, fgabnd
      CHARACTER fgxsct*4
      EXTERNAL photo, gphoto, fgabnd, fgxsct

      DATA atomic /1,    2,    6,    7,    8,    10,   11,   12,
     &             13,   14,   16,   17,   18,   20,   24,   26,
     &             27,   28/

      DATA celts /'H ', 'He', 'C ', 'N ', 'O ', 'Ne', 'Na', 'Mg',
     &            'Al', 'Si', 'S ', 'Cl', 'Ar', 'Ca', 'Cr', 'Fe',
     &            'Co', 'Ni'/


      status = 0

c Set the element columns

      DO i = 1, NELTS
c         column(i) = param(i) * 1.e22 * fgabnd(celts(i))
         column(i) = 0.0
      ENDDO
c     only C,H,O are considered...
      column(3) = param(1) * 1.e18
      column(1) = param(2) * 1.e18 * param(1)
      column(5) = param(3) * 1.e18 * param(1)

c      zfac = 1.0 + param(19)
      zfac = 1.0 
      elow = ear(0) * zfac

      DO ie = 1, ne

         ehi = ear(ie) * zfac
         xsect = 0.

c Set transmission to unity above 100 keV and below the Lyman limit

         IF ( elow .GT. 100. .OR. elow .LT. LYLIMIT ) THEN

            photar(ie) = 1.0

c else calculate the cross-section

         ELSE

            DO i = 1, NELTS
               IF ( fgxsct() .EQ. 'vern' ) THEN
                  pret = gphoto(elow, ehi, atomic(i), status) 
               ELSEIF ( fgxsct() .EQ. 'bcmc' ) THEN
                  pret = photo(elow, ehi, atomic(i), 3, status)
               ELSEIF ( fgxsct() .EQ. 'obcm' ) THEN
                  pret = photo(elow, ehi, atomic(i), 2, status)
               ENDIF
               IF ( column(i) .GT. 0. ) THEN
                  IF ( pret .LT. 1e30/column(i) ) THEN
                     xsect = xsect + pret * column(i)
                  ELSE
                     xsect = 1e30
                  ENDIF
               ELSE
c                negative column is allowed
                  xsect = xsect + pret * column(i)
               ENDIF

            ENDDO

            photar(ie) = EXP(-xsect)

         ENDIF

         elow = ehi

      ENDDO

      RETURN
      END



