subroutine TemperaturesFromEnergy(theta, Tgap, Ebf, Ebb, nlayer, nperr, ErrorMessage)
!***********************************************************************
! this subroutine computes the new temperature distribution
!***********************************************************************
integer, intent(in) :: nlayer
real(r64), dimension(maxlay), intent(in) :: Ebf, Ebb
real(r64), dimension(maxlay2), intent(inout) :: theta
real(r64), dimension(maxlay1), intent(inout) :: Tgap
!real(r64), intent(out) :: dtmax
!integer, intent(out) :: MaxIndex
integer, intent(inout) :: nperr
character(len=*), intent(inout) :: ErrorMessage
real(r64), dimension(maxlay2) :: told
integer :: i, j
!dtmax = 0.0
!MaxIndex = 0
! first check for energy values. They cannot be negative because power to 0.25
! will crash application
do i = 1, nlayer
if ((Ebf(i).lt.0).and.(Ebb(i).lt.0)) then
nperr = 2 !this is flag for convergence error
ErrorMessage = 'Tarcog failed to converge.'
return ! stop execution
end if
end do
do i=1, nlayer
j = 2*i
told(j) = theta(j)
told(j-1) = theta(j-1)
theta(j-1) = (Ebf(i)/StefanBoltzmann)**0.25d0
theta(j) = (Ebb(i)/StefanBoltzmann)**0.25d0
if (i.ne.1) then
Tgap(i) = (theta(j-1) + theta(j-2)) / 2
end if
end do
end subroutine TemperaturesFromEnergy