number_to_date Subroutine

    subroutine number_to_date(date_num, date_array, ref_year)
        !! Returns array with date and time from a date number
        !!
        !! 'date-num' can be any non-negative number.
        real(dp), intent(in) :: date_num !! Date number in seconds since ref_year
        integer, intent(out) :: date_array(6) !! Datetime [y,m,d,h,m,s]
        integer, intent(in) :: ref_year !! Reference year

        ! Local variables
        integer :: y, m, d
        real(dp) :: day_fraction
        integer :: day_int
        integer :: day_count, days_in_year
        integer :: rest_seconds
        integer :: daysinmonth(12) = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]

        ! Check that date number if positive
        if (date_num < 0) then
            print "(a)", "ERROR! In number_to_date, 'date_num' must be >= 0.0"
            stop 1
        end if

        ! Set day fraction to the nearest second. Avoiding round off errors
        day_int = idint(date_num)
        day_fraction = date_num - day_int

        ! Determine hours, minutes and seconds
        date_array = 0
        rest_seconds = int(day_fraction*24.0*3600.0 + 0.5) ! Rounded off
        date_array(4) = int(rest_seconds/3600.0)
        date_array(5) = int((rest_seconds/60.0 - date_array(4)*60.0))
        date_array(6) = int((rest_seconds - date_array(4)*3600.0 - date_array(5)*60.0))

        ! Count up days keeping track of the year, month and day of month
        
        ! Determine year
        y = ref_year - 1
        day_count = 0
        do while (day_count .le. day_int)
            y = y + 1
            days_in_year = 365
            if (((mod(y,4) .eq. 0) .and. (mod(y,100) .ne. 0)) .or. (mod(y,400) .eq. 0)) days_in_year = 366 ! Leap year
            day_count = day_count + days_in_year
        end do
        date_array(1) = y
        day_count = day_count - days_in_year

        ! Determine month given the year
        daysinmonth(2) = 28
        if (((mod(date_array(1),4) .eq. 0) .and. (mod(date_array(1),100) .ne. 0)) .or. (mod(date_array(1),400) .eq. 0)) daysinmonth(2) = 29 ! Leap month
        m = 0

        do while (day_count .le. day_int)
            m = m + 1
            day_count = day_count + daysinmonth(m)
        end do
        date_array(2) = m
        day_count = day_count - daysinmonth(m)

        ! Determine day
        d = 0
        do while (day_count .le. day_int)
            d = d + 1
            day_count = day_count + 1
        end do
        date_array(3) = d
    end subroutine number_to_date