subroutine jdayi(julday, iyear, imon, iday)
! This routine calculates the month (imon) and day (iday)
! from the Julian day (julday) and year (iyear).
! The appropriate correction is made for leap year.

! Modifications:
! 4/29/21 (J. Dostalek) Started this Fortran 90 version
!-------------------------------------------------------------------------------

  implicit none

  integer, intent(in) :: julday
  integer, intent(in) :: iyear

  integer, intent(out) :: imon
  integer, intent(out) :: iday

  integer :: i
  integer :: mxjul

  integer, dimension(12) :: ndmon
  integer, dimension(13) :: nsum

! Specify the number of days in each month
  ndmon(1)  = 31
  ndmon(2)  = 28
  ndmon(3)  = 31
  ndmon(4)  = 30
  ndmon(5)  = 31
  ndmon(6)  = 30
  ndmon(7)  = 31
  ndmon(8)  = 31
  ndmon(9)  = 30
  ndmon(10) = 31
  ndmon(11) = 30
  ndmon(12) = 31

! Correct for leap year
  if (mod(iyear,   4) == 0) ndmon(2) = 29
  if (mod(iyear, 100) == 0) ndmon(2) = 28
  if (mod(iyear, 400) == 0) ndmon(2) = 29

! Check for illegal input
  if (ndmon(2) == 29) then
    mxjul = 366
  else
    mxjul = 365
  end if

  if (julday < 1 .or. julday > mxjul) then
    imon = -1
    iday = -1
    return
  end if

! Calculate the month and day
  nsum(1) = 0
  do i = 1, 12
    nsum(i + 1) = nsum(i) + ndmon(i)
  end do

  do i = 2, 13
    if (julday <= nsum(i)) then
      imon = i - 1
      exit
    end if
  end do

  iday = julday - nsum(imon)

  return
end subroutine jdayi