!
! 19981120  add sections in the main program to read ECMWF verion 0 GRIB
!
! New main program to solve 'GRIB' in BDS problem   19971105
!
! This GRIB decoder has been tested on SGI/iris4, HP/ux, SUN/sunOS4.1.3 and
!*********************************************************************
      FUNCTION julday(id,mm,iyyy)
! See numerical recipes 2nd ed. The order of month and day have been swapped!
!*********************************************************************
      IMPLICIT NONE
	   INTEGER id,mm,iyyy
      INTEGER IGREG 
      INTEGER*4 ja,jm,jy
      INTEGER*4 julday
      IGREG=15+31*(10+12*1582)
      jy=iyyy
      IF (jy.EQ.0) WRITE(6,*) 'There is no zero year !!'
      IF (jy.LT.0) jy=jy+1
      IF (mm.GT.2) THEN
         jm=mm+1
      ELSE
         jy=jy-1
         jm=mm+13
      ENDIF
      julday=INT(365.25*jy)+int(30.6001*jm)+id+1720995
      IF (id+31*(mm+12*iyyy).GE.IGREG) THEN
         ja=INT(0.01*jy)
         julday=julday+2-ja+INT(0.25*ja)
      ENDIF
      RETURN
      END 

!*********************************************************************
      SUBROUTINE caldat(julian,id,mm,iyyy)
! See numerical recipes 2nd ed. The order of month and day have been swapped!
!*********************************************************************
      IMPLICIT NONE
	   INTEGER*4 julian
      INTEGER id,mm,iyyy
      INTEGER IGREG
      INTEGER*4 ja,jalpha,jb,jc,jd,je

      IGREG=2299161
      if (julian.GE.IGREG) THEN
         jalpha=INT(((julian-1867216)-0.25)/36524.25)
         ja=julian+1+jalpha-INT(0.25*jalpha)
      ELSE
         ja=julian
      ENDIF
      jb=ja+1524
      jc=INT(6680.+((jb-2439870)-122.1)/365.25)
      jd=365*jc+INT(0.25*jc)
      je=INT((jb-jd)/30.6001)
      id=jb-jd-INT(30.6001*je)
      mm=je-1
      IF (mm.GT.12) mm=mm-12
      iyyy=jc-4715
      IF (mm.GT.2) iyyy=iyyy-1
         IF (iyyy.LE.0) iyyy=iyyy-1
      RETURN
      END 

PROGRAM jm1
INTEGER iy,im,id,iy2,im1,im2,id1,id2
INTEGER(KIND=4) jday,jday2
 
  OPEN(2,file='date1')
  READ(2,'(I4,I2,I2)') iy,im,id
  CLOSE(2)
  jday=julday(id,im,iy)
  jday2=jday-1
  call caldat(jday2,id2,im2,iy2)
  OPEN(2,file='date2')
         im21=INT(im2/10)
         im22=MOD(im2,10)
         id21=INT(id2/10)
         id22=MOD(id2,10)
  WRITE(2,'(I4,4I1)') iy2,im21,im22,id21,id22
  CLOSE(2)
  jday2=jday+1
  call caldat(jday2,id2,im2,iy2)
  OPEN(2,file='date3')
         im21=INT(im2/10)
         im22=MOD(im2,10)
         id21=INT(id2/10)
         id22=MOD(id2,10)
  WRITE(2,'(I4,4I1)') iy2,im21,im22,id21,id22
  CLOSE(2)
END
