program COMPUTE_CHI2

c  Author:  Glen Cowan
c  Date:    7-OCT-1997

c  Program reads in a data histogram and a theoretical prediction
c  and computes a chi2 statistic.

c  The program also demonstrates some handy code for reading
c  data from ASCII files.  If the first character in a line is !, then
c  the line is skipped.  This can be used for adding comments in
c  data files or for commenting out data lines.

      integer         max_points
      parameter       (max_points = 200)      ! increase if necessary

      character*80    infile
      character*80    line

      integer         i
      integer         ios
      integer         lun
      integer         n (max_points)
      integer         num_points

      real            chi2
      real            nu (max_points)
      real            x_min (max_points)
      real            x_max (max_points)

c  Read data histogram from file.  Lines starting with ! are skipped.

      write (*, *) 'enter file name for data histogram'
      read (*, '(a80)') infile
      lun = 30
      open (unit=lun, file = infile,
     1  form = 'formatted', status = 'old')

      i = 0
      ios = 0
      do while ( ios .eq. 0 )
        read (lun, fmt = '(a80)', iostat = ios) line
        if (ios .eq. 0 .and. .not. line(1:1) .eq. '!') then
          i = i + 1
          read (line, *) x_min(i), x_max(i), n(i)
        endif
      end do
      close (lun)
      num_points = i

c  Read theoretical prediction from file.

      write (*, *) 'enter file name for theoretical prediction'
      read (*, '(a80)') infile
      lun = 30
      open (unit=lun, file = infile,
     1  form = 'formatted', status = 'old')

      i = 0
      ios = 0
      do while ( ios .eq. 0 )
        read (lun, fmt = '(a80)', iostat = ios) line
        if (ios .eq. 0 .and. .not. line(1:1) .eq. '!') then
          i = i + 1
          read (line, *) x_min(i), x_max(i), nu(i)
        endif
      end do
      close (lun)

      if ( i .ne. num_points ) then
        write (*, *) 'Number of points in data and theory not equal'
        write (*, *) 'Data file has   ', num_points, ' lines'
        write (*, *) 'Theory file has ', num_points, ' lines'
      else
        write (*, *) num_points, ' lines read in for data and theory'
      endif

c  Compute chi2.  (Note:  all nu(i) should be nonzero.)

      chi2 = 0.
      do i = 1, num_points
        chi2 = chi2 + (FLOAT(n(i)) - nu(i))**2 / nu(i)
      end do

      write (*, *) 'chi2 = ', chi2
      write (*, *) 'number of data points = ', num_points

      stop
      END