minor updates

SCRIPTS/eenmaligeScripts/wfdei_GPCC.R

+rm(list=ls())
+library(fields) # e.g: using the fields library
+library(abind)
+library(ncdf4)
+source(file = "./functions/infoGeneral.R")
+source(file = "./functions/functionsGeneral.R")
+source(file = "./functions/functionR2Netcdf.R")
+source(file = "./functions/functionNetcdf2R.R")
+
+submitscript <- FALSE
+
+if (submitscript) {
+  targetMonths <- c(1:12)
+  targetYears  <- c(1981:2011)
+  locName    <- 'X'
+  resolution<-"0.50"
+  inPath<-sprintf("../WFDEI_nc/%sdeg/%s", resolution, locName)
+  outPath<-sprintf("../DATA/wfdei_rev4.0/%sdeg/%s", resolution, locName)
+} else {
+  targetMonths <-c(1:12)
+  targetYears  <-c(1981:2010)
+  # targetYears  <-c(1981:1982)
+  locName<-"GHA"
+  locName<-"EU"
+  resolution<-"0.50"
+  inPath<-sprintf("../WFDEI_nc/")
+  outPath<-sprintf("../DATA/wfdei_rev4.0/%sdeg/%s", resolution, locName)
+}
+
+variables<-names(variableInfo)
+variables<-c( "pr", "sfcWind", "rsds", "rlds", "huss", "tas", "tasmin", "tasmax")
+variables<-c( "pr")
+
+dir.create(outPath, recursive = TRUE, showWarnings = FALSE)
+variableName<-"pr"
+targetMonth<-1
+for (targetMonth in targetMonths) {
+  
+  print(sprintf("targetMonth: %s, years: %d-%d",month.name[targetMonth], targetYears[1], targetYears[length(targetYears)]))
+  iFile <- sprintf("%s/%s_daily_WFDEI_GPCC_grid_%4d-%4d_%s_%02d.nc4",
+                   inPath, "Precip", 
+                   targetYears[1], targetYears[length(targetYears)],
+                   locName, targetMonth)
+  print(iFile)
+  RData<-Netcdf2R(iFile,"Precip")
+  
+  oPrefix <- sprintf("%s/%s_wfdei_GPCC_%s_%4d-%4d_%02d",
+                     outPath, variableName, locName, 
+                     targetYears[1], targetYears[length(targetYears)],
+                     targetMonth)
+  print(oPrefix)
+  
+  if (resolution == "0.50") {
+    lonLim = c(locationInfo$res0.50[[locName]]$lonmin, locationInfo$res0.50[[locName]]$lonmax)
+    latLim = c(locationInfo$res0.50[[locName]]$latmin, locationInfo$res0.50[[locName]]$latmax) 
+  } else {
+    lonLim = c(locationInfo$res0.75[[locName]]$lonmin, locationInfo$res0.75[[locName]]$lonmax)
+    latLim = c(locationInfo$res0.75[[locName]]$latmin, locationInfo$res0.75[[locName]]$latmax) 
+  }
+  
+  RData$xyCoords$x[]<-round(RData$xyCoords$x[],2)
+  RData$xyCoords$y[]<-round(RData$xyCoords$y[],2)
+  
+  if (resolution == "0.75") {
+    newGrid<-getGrid(RData)
+    newGrid$x <- c(locationInfo$res0.75[[locName]]$lonmin, locationInfo$res0.75[[locName]]$lonmax)
+    newGrid$y <- c(locationInfo$res0.75[[locName]]$latmin, locationInfo$res0.75[[locName]]$latmax)
+    attr(newGrid,"resX")<-0.75
+    attr(newGrid,"resY")<-0.75
+    print("Regridding")
+    RData <- interpGridData(RData, new.grid = newGrid, method = "bilinear");
+  }
+  
+  if (resolution == "0.50") {
+    RData_new<-RData
+    iLonMin<-which(RData$xyCoords$x==locationInfo$res0.50[[locName]]$lonmin)
+    iLonMax<-which(RData$xyCoords$x==locationInfo$res0.50[[locName]]$lonmax)
+    iLatMin<-which(RData$xyCoords$y==locationInfo$res0.50[[locName]]$latmin)
+    iLatMax<-which(RData$xyCoords$y==locationInfo$res0.50[[locName]]$latmax)
+    RData_new$xyCoords$x<-RData$xyCoords$x[iLonMin:iLonMax]
+    RData_new$xyCoords$y<-RData$xyCoords$y[iLatMin:iLatMax]
+    RData_new$Data<-RData$Data[iLonMin:iLonMax,iLatMin:iLatMax,]
+    attr(RData_new$Data,"dimensions")<-c("lon", "lat", "time")
+    attr(RData_new$xyCoords,"projection")<-"+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0"
+    RData<-RData_new
+    rm(RData_new)
+  }
+    
+    # Check units
+    RData$Data[]<-RData$Data[]*86400
+    RData$Variable$varName<-variableName
+    attr(RData$Variable,"standard_name")<-variableInfo[[variableName]]$standardName
+    attr(RData$Variable,"long_name")<-variableInfo[[variableName]]$longName
+    attr(RData$Variable,"units")<-variableInfo[[variableName]]$unitsEcoms
+    
+    ## add some extra attributes
+    attr(RData,"contact") <- "Wietse Franssen (wietse.franssen@wur.nl)"
+    attr(RData,"based on") <- "GPCC"
+    
+    save(file = sprintf("%s.RData",oPrefix), RData)
+    # R2Netcdf(sprintf("%s.nc4",oPrefix), RData)
+}
+