nccf_times.R 2.68 KB
Newer Older
Burkhardt Rockel's avatar
Burkhardt Rockel committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71
`nccf_times` <-
function(nc) {
    # this function converts netcdf times in a netCDF file 
    # with CF conventions to the 
    # R date-time format
        
    timei   <- which(names(nc$dim) %in% c("time", "TIME", "tim", "TIM"))
    units   <- nc$dim[[timei]]$units
    vals    <- nc$dim[[timei]]$vals
    tmp     <- ncatt_get(nc, names(nc$dim)[timei], "calendar")
    if (tmp$hasatt) {
        calendar <- trimws(tmp$value)
    } else {
        calendar <- "standard"
        print(paste("Warning: Calendar is missing in", nc$filename))
        print(paste('         "standard" calendar is used'))
    }    
#    if (calendar == "proleptic_gregorian") calendar <- "standard"
#... valid calendar names
    cals <- c("standard", "gregorian", "proleptic_gregorian", "360_day", "360", "365_day", "365", "noleap")
    itype.cals <- c(0, 0, 0, 1, 1, 2, 2, 2)
    if (identical(calendar %in% cals, FALSE)){
        print(paste("Warning: ",calendar," is not a valid calendar name" ))
        print(paste("valid names are: ", cals))
        print(paste('         "standard" calendar is used'))     
        calendar <- "standard"
    }

    units.dt <- c("seconds since", "minutes since", "hours since","days since","months since","years since")
    mul.dt <- c(1, 60, 3600, 86400, 1, 1)
#... find the index of matching units.td in units
    index.dt <- match(1,sapply (units.dt,function(x) charmatch(x, units, nomatch=0)),nomatch=0)
    if (index.dt == 0) {
         print(paste("ERROR: ",units,"is no valid unit"))
         print(paste("accepted units are:"))
         print(units.dt)
         stop
#         quit()
    }

    ref.time     <- trimws(gsub(units.dt[index.dt],"",units))
    ref.time.pcict <- as.PCICt(ref.time, calendar)

    if (index.dt == 5) {
      ymd<-format(as.Date(ref.time), "%Y-%m-01")
      ref.time.pcict.month<-as.PCICt(ymd,calendar)
      y<-sapply (vals,function(x) seq(ref.time.pcict.month, by = "month", length = x+1))
      z<-lapply(seq(1,length(y),1), function(x) tail(y[[x]],-(length(y[[x]])-1)))
      zz<-do.call(c,z)
      zzz<-as.PCICt(zz, calendar) - ref.time.pcict.month
      vals<-do.call(c,lapply(seq(1,length(zzz),1),function(x) zzz[[x]]))
    }

    if (index.dt == 6) {
      ymd<-format(as.Date(ref.time), "%Y-01-01")
      ref.time.pcict.year<-as.PCICt(ymd,calendar)
      y<-sapply (vals,function(x) seq(ref.time.pcict.year, by = "year", length = x+1))
      z<-lapply(seq(1,length(y),1), function(x) tail(y[[x]],-(length(y[[x]])-1)))
      zz<-do.call(c,z)
      zzz<-as.PCICt(zz, calendar) - ref.time.pcict.year
      vals<-do.call(c,lapply(seq(1,length(zzz),1),function(x) zzz[[x]]))
    }

    mul<-mul.dt[index.dt]
    times       <- vals * mul
    time <- ref.time.pcict + times

    time

}