function omerf (x, y, sq)
c     ******************************************************************
c
c        this subroutine calculates chapman function as approximated by
c        john a. fitzmaurice, appl. opt. 3,640(1964).
c
c        method:
c                omerf(h,mu)= sqrt(h*pi/2) exp(h*(mu**2)/2)
c
c                where-  h=layer height from the center of earth
c                          (in units of scale height of attenuating
c                                          species)
c
c
c     ******************************************************************
      implicit integer*4(i-n),real*8 (a-h,o-z)
      real *8 a(5),dp,dy,dx,dt,dsum,dsq,domerf
      data dp/0.3275911d+00/
      data a    /
     1         1.061405429d+00,        -1.453152027d+00,
     2         1.421413741d+00,        -0.284496736d+00,
     3         0.254829592d+00     /
      dy = y
      dx = x
      dsq = sq
      dt = 1.0d+00 / (1.0d+00 + dp * dx)
      dsum = dt * a(1)
      if (y .ge. 0.2) go to 20
      do 10 i=2,5
   10 dsum = dt * (a(i) + dsum)
      domerf = dsum * dsq
      omerf = domerf
      return
c
c*****for all but large solar zenith angles use the series expansion
c*****       of the chapman function
c
   20 dt = 0.5d+00 / (dx**2)
      domerf=(1.0d+00-dt*(1.d+00-3.d+00*dt*(1.d+00-5.d+00*dt)))/dy
      omerf = domerf
      return
      end