PROGRAM VIIRS_COMPUTE_COEFF

c ... Given a VIIRS spectral response data file in Hal Woolf format,
c ... compute central wavenumber, temperature correction slope,
c ... and temperature correction intercept

      implicit none

      integer lun      
      character*80 filename, header
      integer band, npts
      real w1, w2
      real srf(10000), wav(10000)
      integer index
      real center
      real t1, t2
      real tcs, tci
      integer nbands
      integer band_store(100)
      real center_store(100), tcs_store(100), tci_store(100)
      real tlower(16), tupper(16)
      
      real cenwav
      external cenwav

c ... Lower and upper temperature limits for VIIRS
      data tlower /16 * 180.0/
      data tupper /16 * 340.0/
      
c ... Open the spectral response data file and read header
      write(*, '(''Enter name of VIIRS spectral response file: ''$)')
      read(*, '(a80)') filename
      lun = 50
      open(lun, file=filename, status='old')
      read(lun, '(a80)') header
      write(*, '(a80,/)') header

c ... Start loop over VIIRS spectral bands
      nbands = 0
20    continue

c ... Read the data for this spectral band
      read(lun, '(i3, i5, 2f8.2)', end=40) band, npts, w1, w2
      read(lun, '(8f9.6)') (srf(index), index = 1, npts)
      nbands = nbands + 1
      band_store(nbands) = band
      
c ... Compute the effective central wavenumber
      center = cenwav(npts, w1, w2, srf)
      center_store(nbands) = center

c ... Construct wavenumber array
      do index = 1, npts
        wav(index) = float(index - 1) * 0.1 + w1
      end do

c ... Compute the temperature correction slope and intercept
      t1 = tlower(band)
      t2 = tupper(band)
      call bndfit(npts, wav, srf, center, t1, t2, tcs, tci)
      tcs_store(nbands) = tcs
      tci_store(nbands) = tci

      write(*,'(a,i2,a,f7.2,a,f7.2,a)')
     &  'BAND ', band, ' TEMPERATURE RANGE WAS ', t1, ' K TO ', t2, ' K'

      goto 20
         
40    continue

c ... Close spectral response file
      close(lun)

      write(*,'(/,a)') 'BANDS'
      write(*,'(4(i4,'',''))')
     &  (band_store(index), index= 1, nbands)
      
      write(*,'(/,a)') 'EFFECTIVE CENTRAL WAVENUMBERS'
      write(*,'(4(1pe13.6,'',''))')
     &  (center_store(index), index= 1, nbands)

      write(*,'(/,a)') 'TEMPERATURE CORRECTION SLOPES'
      write(*,'(4(1pe13.6,'',''))')
     &  (tcs_store(index), index= 1, nbands)
      
      write(*,'(/,a)') 'TEMPERATURE CORRECTION INTERCEPTS'
      write(*,'(4(1pe13.6,'',''))')
     &  (tci_store(index), index= 1, nbands)
      
      END