Commit c28c3575 authored by Franssen, Wietse's avatar Franssen, Wietse
Browse files

Added selection of timesteps

parent 5112c518
......@@ -8,7 +8,13 @@
## "export": export the help and make it available
#' @export
plotje <-function(data, title = " ") {
image.plot(data$xyCoords$x,data$xyCoords$y,data$Data[,,1], asp = 1, main = title, xlab = '', ylab = '')
if (length(dim(data$Data)) > 2) {
dataTmp<-data$Data[,,1]
cat("Plotting first timestep\n")
} else {
dataTmp<-data$Data
}
image.plot(data$xyCoords$x,data$xyCoords$y,dataTmp, asp = 1, main = title, xlab = '', ylab = '')
world(add = TRUE)
}
......@@ -18,9 +24,8 @@ plotje <-function(data, title = " ") {
#' @author Wietse Franssen \email{wietse.franssen@@wur.nl}
#' @keywords internal
#' @export
ncPlot <-function(file, varName = NULL, lonlatbox = NULL) {
data<-ncLoad(file, varName, lonlatbox)
plotje(data, title = paste0(data$Variable$varName, " (" ,data$Dates$start[1], ")"))
ncPlot <-function(file, varName = NULL, lonlatbox = NULL, timesteps = NULL) {
data<-ncLoad(file, varName, lonlatbox, timesteps)
title=paste0(data$Variable$longName, " [",data$Variable$units,"]\n(" ,data$Dates$start[1], ")")
plotje(data, title = title)
}
......@@ -55,7 +60,7 @@ ncCheck <-function(ncFile, variable) {
result$dims[i]<-ncFile$dim[[i]]$name
result$dimids[i]<-ncFile$dim[[i]]$id
}
if (length(times) > 1) { result$tArray<-TRUE }
if ((lons[2]-lons[1]) < 0) { result$xFlip<-TRUE }
if ((lats[2]-lats[1]) < 0) { result$yFlip<-TRUE }
......@@ -76,7 +81,7 @@ ncCheck <-function(ncFile, variable) {
#' @author Wietse Franssen \email{wietse.franssen@@wur.nl}
#' @keywords internal
#' @export
ncLoad <-function(file, varName = NULL, lonlatbox = NULL) {
ncLoad <-function(file, varName = NULL, lonlatbox = NULL, timesteps = NULL) {
## first build an empty data structure
data<-rDataStructure()
......@@ -106,16 +111,21 @@ ncLoad <-function(file, varName = NULL, lonlatbox = NULL) {
LonIdx <- c(1: ncFile$dim$lon$len)
LatIdx <- c(1: ncFile$dim$lat$len)
} else {
# LonIdx <- which( ncFile$dim$lon$vals > lonlatbox[1] | ncFile$dim$lon$vals < lonlatbox[2])
LonIdx <- which( ncFile$dim$lon$vals > lonlatbox[1] & ncFile$dim$lon$vals < lonlatbox[2])
LatIdx <- which( ncFile$dim$lat$vals > lonlatbox[3] & ncFile$dim$lat$vals < lonlatbox[4])
# LonIdx <- which( ncFile$dim$lon$vals > lonlatbox[1] | ncFile$dim$lon$vals < lonlatbox[2])
LonIdx <- which( ncFile$dim$lon$vals >= lonlatbox[1] & ncFile$dim$lon$vals <= lonlatbox[2])
LatIdx <- which( ncFile$dim$lat$vals >= lonlatbox[3] & ncFile$dim$lat$vals <= lonlatbox[4])
}
data$xyCoords$x <- ncFile$dim$lon$vals[LonIdx]
data$xyCoords$y <- ncFile$dim$lat$vals[LatIdx]
#data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
## Fill time
NCtime <- ncvar_get( ncFile, "time" )
if (is.null(timesteps)) {
timeIdx<-c(1)
} else {
timeIdx<-timesteps
}
NCtime <- ncvar_get( ncFile, "time" ) [timeIdx]
#NCtime <- ncvar_get( ncFile, "time" )
NCtimeAtt <- ncatt_get( ncFile, "time", "units" )$value
firstTime<-unlist(strsplit(NCtimeAtt, split=' ', fixed=TRUE))[3]
firstTime<-strptime(firstTime, format = "%Y-%m-%d", tz = "GMT")
......@@ -124,7 +134,7 @@ ncLoad <-function(file, varName = NULL, lonlatbox = NULL) {
## Fill data
#data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
data$Data <- ncvar_get( ncFile, varName ) [ LonIdx, LatIdx, ]
data$Data <- ncvar_get( ncFile, varName ) [ LonIdx, LatIdx, timeIdx]
#data$Data <- ncvar_get( ncFile, varName )
attr(data$Data,"dimensions") <- ncCheckResult$dims
data$Variable$varName <- varName
......@@ -167,60 +177,60 @@ ncLoad <-function(file, varName = NULL, lonlatbox = NULL) {
# rm(list=ls())
# library(ncdf4)
# library(fields) # e.g: using the fields library
#
# plotje <-function(plottitle) {
# image.plot(NClon,NClat,data, asp = 1, main = plottitle, xlab = '', ylab = '')
# world(add = TRUE)
# }
#
# domainName<-c( "GHA", "EU")
# lonmin<- c( 27.75, -24.75)
# lonmax<- c( 49.25, 39.75)
# latmin<- c(-12.25, 33.25)
# latmax<- c( 18.25, 71.75)
#
# iDomain<-1
#
# ## READ WFD NETCDF
# ncFile <- nc_open( "~/Desktop/gg/wfd_pr_1974.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
# landmask<-data
# landmask[!is.na(landmask)] <- 1
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ WFDEI NETCDF
# ncFile <- nc_open( "~/Desktop/gg/wfd_pr_1979.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ SOIL NETCDF
# ncFile <- nc_open( "~/Desktop/gg/soil_GHA.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "stexture")[ LonIdx, LatIdx]
# landmask<-data
# landmask[!is.na(landmask)] <- 1
#
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ SOIL NETCDF
# nameFileNCin<-"~/Desktop/gg/soil_GHA.nc"
# ncid_in=nc_open(nameFileNCin)
# NCdata <- ncvar_get( ncid_in, "stexture")
# nc_close(ncid_in)
#
# ## Plot and add to LandMask
# data<-NCdata[,];plotje(plottitle = "WFDEI")
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
# library(ncdf4)
# library(fields) # e.g: using the fields library
#
# plotje <-function(plottitle) {
# image.plot(NClon,NClat,data, asp = 1, main = plottitle, xlab = '', ylab = '')
# world(add = TRUE)
# }
#
# domainName<-c( "GHA", "EU")
# lonmin<- c( 27.75, -24.75)
# lonmax<- c( 49.25, 39.75)
# latmin<- c(-12.25, 33.25)
# latmax<- c( 18.25, 71.75)
#
# iDomain<-1
#
# ## READ WFD NETCDF
# ncFile <- nc_open( "~/Desktop/gg/wfd_pr_1974.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
# landmask<-data
# landmask[!is.na(landmask)] <- 1
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ WFDEI NETCDF
# ncFile <- nc_open( "~/Desktop/gg/wfd_pr_1979.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "pr")[ LonIdx, LatIdx, 1]
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ SOIL NETCDF
# ncFile <- nc_open( "~/Desktop/gg/soil_GHA.nc" )
# LonIdx <- which( ncFile$dim$lon$vals > lonmin[iDomain] | ncFile$dim$lon$vals < lonmax[iDomain])
# LatIdx <- which( ncFile$dim$lat$vals > latmin[iDomain] & ncFile$dim$lat$vals < latmax[iDomain])
# data <- ncvar_get( ncFile, "stexture")[ LonIdx, LatIdx]
# landmask<-data
# landmask[!is.na(landmask)] <- 1
#
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
# nc_close(ncFile)
#
# ## READ SOIL NETCDF
# nameFileNCin<-"~/Desktop/gg/soil_GHA.nc"
# ncid_in=nc_open(nameFileNCin)
# NCdata <- ncvar_get( ncid_in, "stexture")
# nc_close(ncid_in)
#
# ## Plot and add to LandMask
# data<-NCdata[,];plotje(plottitle = "WFDEI")
# landmask[is.na(data)] <- NA
# sum( !is.na( landmask ) )
......@@ -4,7 +4,7 @@
\alias{ncLoad}
\title{ncLoad}
\usage{
ncLoad(file, varName = NULL, lonlatbox = NULL)
ncLoad(file, varName = NULL, lonlatbox = NULL, timesteps = NULL)
}
\arguments{
\item{file}{Name of the NetCDF file}
......
......@@ -4,7 +4,7 @@
\alias{ncPlot}
\title{Make a image of the NetCDF data}
\usage{
ncPlot(file, varName = NULL, lonlatbox = NULL)
ncPlot(file, varName = NULL, lonlatbox = NULL, timesteps = NULL)
}
\value{
Nothing
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment