mode_read_date_netcdf.F90

Go to the documentation of this file.

MODULE mode_trip_date_netcdf
!
! Module to read correctly time in a netcdf file
! Author : Matthieu Lafaysse
! Creation : 2012-11-12
! Modifications
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
INTEGER,PARAMETER::iyear_gregoire=1582 !Année où le pape Grégoire XIII décida de passer du calendrier julien au calendrier grégorien
!Avant cette date tous les années multiples de 4 sont des années bissextiles
!Au-delà de cette date les années multiples de 100 ne sont pas bissextiles sauf les années multiples de 400 qui le sont
!De plus les journées du 5 au 14 octobre 1582 ont été supprimées pour compenser le retard accumulé
!
 CONTAINS
!
!----------------------------------------------------------------------------------------------------------------
!
SUBROUTINE read_date(PTIME,HUNITS,KYEAR,KMONTH,KDAY,PTIME_DATE)
!
!*      0.1   Declarations of arguments
!
REAL,DIMENSION(:),INTENT(IN)  :: PTIME
 CHARACTER(*),     INTENT(IN)  :: HUNITS
!
INTEGER,          INTENT(OUT) :: KYEAR
INTEGER,          INTENT(OUT) :: KMONTH
INTEGER,          INTENT(OUT) :: KDAY
REAL,             INTENT(OUT) :: PTIME_DATE
!
!*      0.2   Declarations of local variables
!
INTEGER,DIMENSION(SIZE(PTIME)) :: ITIMEHOURS
REAL,   DIMENSION(SIZE(PTIME)) :: ZREST
!
INTEGER :: IYEAR
INTEGER :: IMONTH
INTEGER :: IDAY
INTEGER :: IHOUR
REAL    :: ZTIME
!
 CHARACTER(LEN=14) :: YHEADER
 CHARACTER(LEN=100) :: YFMT 
 CHARACTER :: YC1,YC2,YC3
!
LOGICAL :: GHEADER
!
INTEGER,DIMENSION(4),PARAMETER :: ITV_DATEREF = (/1900,1,1,0/)
INTEGER,DIMENSION(4),PARAMETER :: ITV_DATEREFBIS = (/1850,1,1,0/)
!
INTEGER,DIMENSION(4) :: ITV_DATETEMPO
INTEGER :: ICARACUNITS,ICARACHOUR,ICARACDAY,ICARACMONTH,ICARACYEAR
INTEGER :: IERROR
INTEGER :: IV_DIFFREF
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:READ_DATE',0,zhook_handle)
!
SELECT CASE (hunits(1:4))
  CASE ("days")
    icaracunits = 11
    itimehours(:) = int(ptime(:)*24.)
    zrest(:) = ptime(:)*24.-itimehours(:)
  CASE ("hour")
    icaracunits = 12
    itimehours(:) = int(ptime(:))
    zrest(:) = ptime(:)-itimehours(:)
  CASE ("minu")
    icaracunits = 14
    itimehours(:) = int(ptime(:)/60.)
    zrest(:) = (ptime(:)/60.)-itimehours(:)
  CASE ("seco")
    icaracunits = 14
    itimehours(:) = int(ptime(:)/3600.)
    zrest(:) = (ptime(:)/3600.)-itimehours(:)
  CASE DEFAULT
    print*, "ERROR date_netcdf.F90 : Can't read time units :"
    print*,trim(hunits)
    stop "Error units time"
END SELECT
!
! * Read reference date
!
SELECT CASE (hunits(1:4))
  CASE ("days","hour","minu","seco")
    boucles_hours:DO icarachour=1,2 !Hour written by 1 or 2 char
      DO icaracday=1,2 !Day written by 1 or 2 char
        DO icaracmonth=1,2 !Month written by 1 or 2 char
          DO icaracyear=1,4 !Year written by 1 to 4 char
            WRITE(yfmt,fmt='("(A",I2,",",2("I",I1,",A1,"),"I",I1,",X,I",I1,",A1)")')&
              icaracunits,icaracyear,icaracmonth,icaracday,icarachour
            READ(hunits,fmt=yfmt,iostat=ierror)yheader,iyear,&
              yc1,imonth,yc2,iday,ihour,yc3
            IF (ierror==0) THEN
              gheader=(trim(yheader)=='hours since') .OR. (trim(yheader)=='minutes since') .OR. &
                      (trim(yheader)=='seconds since') .OR. (trim(yheader)=='days since')
              IF (gheader.AND.(yc1=='-').AND.(yc2=='-').AND.(yc3==':')) THEN
                 EXIT boucles_hours
              END IF
            END IF
          END DO
        END DO
      END DO
    END DO boucles_hours
  CASE DEFAULT
    stop "Error units time"
END SELECT
!
! * Check successful reading
!
IF (ierror>0) THEN
  print*, "ERROR date_netcdf.F90 : Can't read time units :"
  print*,trim(hunits)
  stop "Error units time"
END IF
!
!Initialiaze the date to the reference date
!
ztime = REAL(ihour)+ZREST(1)
!
!Add the number of hours to a date
!
 CALL addtime(itimehours(1),iyear,imonth,iday,ztime)
!
! Current date
!
kyear      = iyear
kmonth     = imonth
kday       = iday
ptime_date = ztime
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:READ_DATE',1,zhook_handle)

END SUBROUTINE read_date
!
!----------------------------------------------------------------------------------------------------------------
!
SUBROUTINE addtime (KNHOURS,KYEAR,KMONTH,KDAY,PTIME)
!
!*      0.1   Declarations of arguments
!
INTEGER, INTENT(IN)  :: KNHOURS ! number of hours
!
INTEGER, INTENT(OUT) :: KYEAR
INTEGER, INTENT(OUT) :: KMONTH
INTEGER, INTENT(OUT) :: KDAY
REAL,    INTENT(OUT) :: PTIME
!
!*      0.2   Declarations of local variables
!
INTEGER,DIMENSION(12) :: INBDM          !Number of days per months
INTEGER               :: IREMAININGDAYS !Number of remaining days to add
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!-------------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:ADDTIME',0,zhook_handle)
!
!special case of year 1582
!
IF (kyear==iyear_gregoire) THEN
  inbdm=(/31,28,31,30,31,30,31,31,30,21,30,31/)
ELSE
  IF (leapyear(kyear)) THEN
    inbdm=(/31,29,31,30,31,30,31,31,30,31,30,31/)
  ELSE
    inbdm=(/31,28,31,30,31,30,31,31,30,31,30,31/)
  ENDIF
ENDIF
!
IF (knhours>=0) THEN
  !Increase hour
  IF ((int(ptime)+knhours)<=23) THEN
    ptime=ptime+knhours
  ELSE
    iremainingdays=(int(ptime)+knhours)/24 !Number of remaining days
    ptime=mod(int(ptime)+knhours,24)+ptime-int(ptime) ! new hour (+minutes)
    !Increase days
    DO
      IF ((kday+iremainingdays)<=inbdm(kmonth)) THEN
        kday=kday+iremainingdays
        EXIT
      ELSE
        iremainingdays=iremainingdays-(inbdm(kmonth)-kday+1)
        IF (kmonth<12) THEN
          kmonth=kmonth+1 !Month change
          kday=1
        ELSE
          !Year change
          kyear=kyear+1
          kmonth=1
          kday=1
          !Update february number of days
          IF (leapyear(kyear)) inbdm(2)=29
          IF (leapyear(kyear-1)) inbdm(2)=28
          !Update october month for 1582 and 1583
          IF (kyear==iyear_gregoire) inbdm(10)=21
          IF (kyear==iyear_gregoire+1) inbdm(10)=31
        ENDIF
      ENDIF
    ENDDO
  ENDIF
ELSE
  IF ((int(ptime)+knhours)>=0) THEN
    ptime=ptime+knhours
  ELSE
    iremainingdays=(int(ptime)-knhours)/24+1
    ptime=mod(int(ptime)-knhours,24)+ptime-int(ptime) ! new hour
    ! decrease days
    DO
      IF  ((kday-iremainingdays)>=1) THEN
        kday=kday-iremainingdays
        EXIT
      ELSE
        iremainingdays=iremainingdays-kday
        IF (kmonth>=1) THEN
          kmonth=kmonth-1 !Month change
          kday=inbdm(kmonth)
        ELSE
          !Year change
          kyear=kyear-1
          kmonth=12
          kday=inbdm(kmonth)
          !Update february number of days
          IF (leapyear(kyear)) inbdm(2)=29
          IF (leapyear(kyear+1)) inbdm(2)=28
          !Update october month for 1582 and 1583
          IF (kyear==iyear_gregoire) inbdm(10)=21
          IF (kyear==iyear_gregoire-1) inbdm(10)=31
        END IF
      END IF
    END DO
  END IF
END IF
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:ADDTIME',1,zhook_handle)
!
END SUBROUTINE addtime
!
!--------------------------------------------------------------------------------------------------------------
!
LOGICAL FUNCTION leapyear (PYEAR)
INTEGER, INTENT(IN) :: PYEAR ! Is the year a leap year ?
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:LEAPYEAR',0,zhook_handle)
!
IF (pyear>iyear_gregoire) THEN
  leapyear = (((mod(pyear,4)==0).AND.(mod(pyear,100)/=0)).OR.(mod(pyear,400)==0))
ELSE
  leapyear = (mod(pyear,4)==0)
ENDIF
!
IF (lhook) CALL dr_hook('MODE_TRIP_DATE_NETCDF:LEAPYEAR',1,zhook_handle)
!
END FUNCTION leapyear
!
!--------------------------------------------------------------------------------------------------------------
!
END MODULE mode_trip_date_netcdf