read_nuclear_rates Subroutine

subroutine read_nuclear_rates(fid, grp, rates)
    !+ Reads in a nuclear reaction rates table from file
    !+ and puts it into a [[NuclearRates]] type
    integer(HID_T), intent(in)              :: fid
        !+ HDF5 file ID
    character(len=*), intent(in)            :: grp
        !+ HDF5 group to read from
    type(NuclearRates), intent(inout)       :: rates
        !+ Atomic reaction rates

    integer(HSIZE_T), dimension(1) :: dim1
    integer(HSIZE_T), dimension(3) :: dim3
    logical :: path_valid, err
    integer :: i, j, error
    real(Float64) :: emin, emax, tmin, tmax, rmin

    err = .False.
    call h5ltpath_valid_f(fid, grp, .True., path_valid, error)
    if(.not.path_valid) then
        if(inputs%verbose.ge.0) then
            write(*,'(a,a)') 'READ_NUCLEAR_RATES: Unknown nuclear interaction: ', trim(grp)
            write(*,'(a)') 'Continuing without neutron calculation'
        endif
        inputs%calc_neutron=0
        inputs%calc_cfpd=0
        return
    endif

    dim1 = [2]
    call h5ltread_dataset_double_f(fid, grp//"/bt_amu", rates%bt_amu, dim1, error)

    if(abs(beam_mass-rates%bt_amu(1)).gt.0.2) then
        if(inputs%verbose.ge.0) then
            write(*,'(a,f6.3,a,f6.3,a)') 'READ_NUCLEAR_RATES: Unexpected beam species mass. Expected ',&
                rates%bt_amu(1),' amu got ', beam_mass, ' amu'
        endif
        err = .True.
    endif

    if(abs(thermal_mass(1)-rates%bt_amu(2)).gt.0.2) then
        if(inputs%verbose.ge.0) then
            write(*,'(a,f6.3,a,f6.3,a)') 'READ_NUCLEAR_RATES: Unexpected thermal species mass. Expected ',&
                 rates%bt_amu(2),' amu got ', thermal_mass(1), ' amu'
        endif
     !!!err = .True.
    endif

    if(err) then
        if(inputs%verbose.ge.0) then
            write(*,'(a)') 'Continuing without neutron calculation'
        endif
        inputs%calc_neutron=0
        inputs%calc_cfpd=0
        return
    endif

    call h5ltread_dataset_int_scalar_f(fid, grp//"/nbranch", rates%nbranch, error)
    call h5ltread_dataset_int_scalar_f(fid, grp//"/nenergy", rates%nenergy, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/emin", emin, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/emax", emax, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/dlogE", rates%dlogE, error)
    call h5ltread_dataset_int_scalar_f(fid, grp//"/ntemp", rates%ntemp, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/tmin", tmin, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/tmax", tmax, error)
    call h5ltread_dataset_double_scalar_f(fid, grp//"/dlogT", rates%dlogT, error)
    rates%logemin = log10(emin)
    rates%logemax = log10(emax)
    rates%logtmin = log10(tmin)
    rates%logtmax = log10(tmax)

    allocate(rates%log_rate(rates%nenergy, &
                            rates%ntemp,   &
                            rates%nbranch))

    dim3 = [rates%nenergy, rates%ntemp, rates%nbranch]
    call h5ltread_dataset_double_f(fid, grp//"/fusion", rates%log_rate, dim3, error)

    rmin = minval(rates%log_rate, rates%log_rate.gt.0.d0)
    where (rates%log_rate.le.0.d0)
        rates%log_rate = 0.9*rmin
    end where
    rates%minlog_rate = log10(rmin)
    rates%log_rate = log10(rates%log_rate)

end subroutine read_nuclear_rates