subroutine rdate_to_jdate &
         (indate,dateout,outyear,outday,outhour)
      
! Change day in date from 'normal' (month/day) to julian day

      implicit none

      character*14 indate
      character*13 dateout
      integer :: year,month,date,hour
      integer :: outyear,outday,outhour
      integer :: iy,jday,julianday
      
      call rdate_unmake_big(year,month,date,hour,indate)
      
      iy = year - 1900
      jday=julianday(month,date,iy)

      outyear=year
      outday=jday
      outhour=hour

      write(dateout(1:4),10) year
      write(dateout(5:7),11) jday
      write(dateout(8:13),12) hour
   10 format (i4.4)      
   11 format (i3.3)      
   12 format (i6.6)      

      end subroutine rdate_to_jdate
!========================================================================

                              subroutine jdate_to_rdate &
         (datein,outdate)
      
! Change day in date from julian day to 'normal' (month/day)
!  Input:  datein  - julian date
! Output: outdate  - regular date

      implicit none

      character*13 datein
      character*14 outdate
      integer :: jday
      integer :: year,month,day,hour
      integer :: iy
      
      call jdate_unmake_big(year,jday,hour,datein)
      
      iy = year - 1900
      call unmake_julianday(jday,iy,month,day)

      write(outdate(1:4),10) year
      write(outdate(5:6),11) month
      write(outdate(7:8),12) day
      write(outdate(9:14),13) hour
   10 format (i4.4)      
   11 format (i2.2)      
   12 format (i2.2)      
   13 format (i6.6)      

      end subroutine jdate_to_rdate
!========================================================================

                              subroutine add_to_jdate &
         (year1,jday1,hour1,diff,year2,jday2,hour2)
      
! Add (subtract) diff (sec) from jdate (year1/jday1/hour1) 
! and get new year2-jday2-hour2
!  for diff>0 => add
!  for diff<0 => subtract

      implicit none

      integer :: year1,jday1,hour1
      integer :: year2,jday2,hour2
      integer :: yday
      integer :: diff
      integer :: diff_d,diff_h,diff_m,diff_s
      integer :: xhourin,xminin,xsecin
      integer :: newsec,newmin,newhrs,newday,newyrs
      integer :: addsec,addmin,addhrs,addday,addyrs
      
      xhourin=hour1/10000
      xminin=mod(hour1,10000)/100
      xsecin=mod(hour1,100)

      diff_d=float(diff)/(24.*3600.)
      diff_h=(diff-diff_d*24*3600)/3600.
      diff_m=(diff-diff_d*24*3600-diff_h*3600)/60.
      diff_s= diff-diff_d*24*3600-diff_h*3600-diff_m*60

!-- correct seconds -----
      newsec=xsecin+diff_s
      addmin=0

      if(newsec.ge.60) then
       newsec=newsec-60
       addmin=1
      elseif(newsec.lt.0) then
       newsec=60+newsec
       addmin=-1
      end if

!-- correct minutes -----
      newmin=xminin+diff_m+addmin
      addhrs=0

      if(newmin.ge.60) then
       newmin=newmin-60
       addhrs=1
      elseif(newmin.lt.0) then
       newmin=60+newmin
       addhrs=-1
      end if

!-- correct hours -------
      newhrs=xhourin+diff_h+addhrs
      addday=0

      if(newhrs.ge.24) then
       newhrs=newhrs-24
       addday=1
      elseif(newhrs.lt.0) then
       newhrs=24+newhrs
       addday=-1
      end if

!-- correct julian days -------
      newday=jday1+diff_d+addday
      addyrs=0

      yday=365+(4-mod(year1,4))/4
      
      if(newday.gt.yday) then
       newday=newday-yday
       addyrs=1
      elseif(newday.le.0) then
       newday=365+(4-mod(year1-1,4))/4 + newday
       addyrs=-1
      end if

!--- output ----------------------

      year2=year1+addyrs
      jday2=newday
      hour2=newhrs*10000 + newmin*100 + newsec

!---------------------------------

      end subroutine add_to_jdate
!========================================================================

                              subroutine diff_date1_date2 & 
         (year1,jday1,hour1,year2,jday2,hour2,diff)

! Substract dates: date1-date2 = diff (in seconds)

      implicit none

      integer :: year1,jday1,hour1
      integer :: year2,jday2,hour2
      integer :: yday1,yday2
      integer :: diff
      integer :: diff_d,diff_h,diff_m,diff_s
      integer :: xhourin,xminin,xsecin
      integer :: yhourin,yminin,ysecin
      integer :: newsec,newmin,newhrs,newday,newyrs
      integer :: addmin,addhrs,addday,addyrs
      
      xhourin=hour1/10000
      xminin=mod(hour1,10000)/100
      xsecin=mod(hour1,100)

      yhourin=hour2/10000
      yminin=mod(hour2,10000)/100
      ysecin=mod(hour2,100)

      diff_d=float(diff)/(24.*3600.)
      diff_h=(diff-diff_d*24*3600)/3600.
      diff_m=(diff-diff_d*24*3600-diff_h*3600)/60.
      diff_s= diff-diff_d*24*3600-diff_h*3600-diff_m*60

      addmin=0
      addhrs=0
      addday=0
      addyrs=0

!-- correct seconds -------------------------

        newsec=xsecin-ysecin
      if(newsec.lt.0) then
        newsec=60+newsec
        addmin=-1
      end if

!-- correct minutes -------------------------

      newmin=xminin-yminin+addmin
      if(newmin.lt.0) then
       newmin=(60+newmin)*60
       addhrs=-1
      else
       newmin=newmin*60
      end if

!-- correct hours -------------------------

      newhrs=xhourin-yhourin+addhrs
      if(newhrs.lt.0) then
       newhrs=(24+newhrs)*3600
       addday=-1
      else
       newhrs=newhrs*3600
      end if

!-- corect days and years -------------------------

      yday1=365+(4-mod(year1,4))/4
      yday2=365+(4-mod(year2,4))/4
       if(year1.gt.year2) then
      newday=(jday1+yday1-jday2+addday)*24*3600
       elseif(year1.eq.year2) then
      newday=(jday1-jday2+addday)*24*3600
       else
      write(*,*) "DATE ERROR: year1,year2=",year1,year2
      RETURN
       end if

!---------------------------------

      diff=newsec+newmin+newhrs+newday

!---------------------------------

      return
      end subroutine diff_date1_date2
!========================================================================

                              subroutine jdate_make_big &
         (inyear,inday,inhour,jdate)

      implicit none

      integer :: inyear,inday,inhour
      character*13 jdate

      write(jdate(1:4),10) inyear
      write(jdate(5:7),11) inday
      write(jdate(8:13),12) inhour
   10 format (i4.4)      
   11 format (i3.3)      
   12 format (i6.6)      

      end subroutine jdate_make_big
!========================================================================

                              subroutine rdate_make_big & 
         (inyear,inmonth,indate,inhour,outdate)

      implicit none

      integer :: inyear,inmonth,indate,inhour
      character*14 outdate

      write(outdate(1:4),10) inyear
      write(outdate(5:6),11) inmonth
      write(outdate(7:8),11) indate
      write(outdate(9:14),12) inhour
   10 format (i4.4)      
   11 format (i2.2)      
   12 format (i6.6)      

      end subroutine rdate_make_big
!========================================================================

                              subroutine jdate_unmake_big &
         (inyear,inday,inhour,dateout)

      implicit none

!  input: outdate - julian date (character)
! output: inyear  - year (integer)
! output: inday   - jukian day (integer)
! output: inhour  - hour,min,sec (integer_6)

      integer :: inyear,inmonth,indate,inhour
      integer :: inday
      character*13 dateout

      read(dateout(1:4),10) inyear
      read(dateout(5:7),11) inday
      read(dateout(8:13),12) inhour
   10 format (i4)      
   11 format (i3)      
   12 format (i6)      
   
      end subroutine jdate_unmake_big
!========================================================================

                              subroutine rdate_unmake_big &
         (inyear,inmonth,indate,inhour,outdate)

      implicit none

      integer :: inyear,inmonth,indate,inhour
      character*14 outdate

      read(outdate(1:4),10) inyear
      read(outdate(5:6),11) inmonth
      read(outdate(7:8),11) indate
      read(outdate(9:14),12) inhour
   10 format (i4)      
   11 format (i2)      
   12 format (i6)      
   
      end subroutine rdate_unmake_big
!========================================================================

                              subroutine unmake_julianday &
         (jday,iyear,imonth,iday)

      implicit none

! Get the month and day from a julian day

      integer :: jday,iyear,imonth,iday
      integer :: idiff,jdiff,n,idiff_old,molength
      integer :: imody(12),icorrect(12)
      data imody /31,28,31,30,31,30,31,31,30,31,30,31/

      icorrect(:)=0
      icorrect(2)=1-min(1,mod(iyear,4))

      idiff=0
      jdiff=999
      n=1
         do while (n.le.12.and.jdiff.gt.0)
      molength=imody(n)+icorrect(n)
      idiff_old=idiff
      idiff=idiff+molength
      jdiff=jday-idiff
      imonth=n
      iday=jday-idiff_old
      n=n+1
        end do

      imonth=n-1
      iday=jday-idiff_old

      end subroutine unmake_julianday
!========================================================================

      function julianday(imonth,iday,iyear)

!     compute the julian day from a normal date

      julianday= iday                                          &
          +min(1,max(0,imonth-1))*31                           &
          +min(1,max(0,imonth-2))*(28+(1-min(1,mod(iyear,4)))) &
          +min(1,max(0,imonth-3))*31                           &
          +min(1,max(0,imonth-4))*30                           &
          +min(1,max(0,imonth-5))*31                           &
          +min(1,max(0,imonth-6))*30                           &
          +min(1,max(0,imonth-7))*31                           &
          +min(1,max(0,imonth-8))*31                           &
          +min(1,max(0,imonth-9))*30                           &
          +min(1,max(0,imonth-10))*31                          &
          +min(1,max(0,imonth-11))*30                          &
          +min(1,max(0,imonth-12))*31
     
      end