module ken_hart

    ! Given the colour, this interpolates
    ! the Kenyon and Hartmann relationship to get the bolometric correction
    ! and colour.

    ! Error flag.
    ! 00 = success.
    ! 12 = out of range for colour.

  use quad

  implicit none

  integer, parameter :: nteff=55
  real, private, save, dimension(nteff) :: teff=0.0
  real, private, save, dimension(nteff) :: bc_v=0.0, v_i=0.0
  logical, private, save :: first=.true.

  contains

  subroutine readfiles_kh()

    integer :: iteff
    real :: u_b, b_v, v_rc, v_rj, v_ij, v_j, v_h, v_k, v_l 
    character(len=5) :: sp 

    character(len=30) :: data_dir

    first=.false.
    ! Open up table 5.
    call get_environment_variable('CMDDATA', data_dir)
    open(unit=1, file=trim(data_dir)//'/KenyonHartmann/table5.dat', &
    action='read')
    read(1,'(/,/)')
    each_teff: do iteff=1, nteff
      read(1,*) sp, teff(iteff), bc_v(iteff), u_b, b_v, v_rc, v_rj, &
      v_i(iteff), v_ij, v_j, v_h, v_k, v_l 
    end do each_teff
    teff=log10(teff)
    close(1)

  end subroutine readfiles_kh

  real function ken_hart_col_vi(teff_in, cflag)

    real, intent(in) :: teff_in
    character(len=*), intent(out) :: cflag

    if (first) call readfiles_kh()
    ken_hart_col_vi=quadint(teff, v_i, teff_in, cflag)

  end function ken_hart_col_vi


  real function ken_hart_bc_v(teff_in, iflag)

    real, intent(in) :: teff_in
    integer, intent(out) :: iflag

    character(len=2) :: cflag

    if (first) call readfiles_kh()
    ken_hart_bc_v=quadint(teff, bc_v, teff_in, cflag)
    iflag=0
    if (cflag /= 'OK') iflag=1

  end function ken_hart_bc_v

end module ken_hart
