From e146658a917e21b74ed660ebf9cc0662e55921a3 Mon Sep 17 00:00:00 2001 From: Jordan S Read Date: Mon, 14 Mar 2016 09:57:08 -0500 Subject: [PATCH] fixes #195 --- DESCRIPTION | 2 +- R/plot_temp.R | 5 +- R/plot_temp_compare.R | 14 ++-- R/plot_var.R | 5 +- R/plot_var_compare.R | 160 ++++++++++++++++++++------------------- R/timeseries_plots.R | 13 ++-- README.md | 6 -- man/plot_temp.Rd | 4 +- man/plot_temp_compare.Rd | 7 +- man/plot_var.Rd | 4 +- man/plot_var_compare.Rd | 4 +- 11 files changed, 118 insertions(+), 106 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 65ac051..1ca8499 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: glmtools Type: Package Title: glmtools -Version: 0.13.0 +Version: 0.13.1 Date: 2016-01-15 Authors@R: c( person("Jordan", "Read", role = c("aut","cre"), email = "jread@usgs.gov"), diff --git a/R/plot_temp.R b/R/plot_temp.R index d1e1753..c9fa675 100644 --- a/R/plot_temp.R +++ b/R/plot_temp.R @@ -2,6 +2,7 @@ #'@param file a string with the path to the netcdf output from GLM #'@param reference a string for 'surface' or 'bottom' #'@param fig_path F if plot to screen, string path if save plot as .png +#' @param col_lim range for heatmap (in units of the variable) #'@param ... additional arguments passed to \code{par()} #'@keywords methods #'@seealso \code{\link{get_temp}}, \code{\link{plot_var}} @@ -16,8 +17,8 @@ #'plot_temp(file = nc_file, fig_path = FALSE) #'plot_temp(file = nc_file, fig_path = 'test_figure.png', height = 3, reference = 'surface') #'@export -plot_temp <- function(file='output.nc', fig_path = FALSE, reference = 'surface', ...){ +plot_temp <- function(file='output.nc', fig_path = FALSE, reference = 'surface', col_lim, ...){ - plot_var(file, var_name = 'temp', fig_path, reference,...) + plot_var(file, var_name = 'temp', fig_path, reference, col_lim, ...) } diff --git a/R/plot_temp_compare.R b/R/plot_temp_compare.R index 2f57075..aa377cb 100644 --- a/R/plot_temp_compare.R +++ b/R/plot_temp_compare.R @@ -1,8 +1,10 @@ -#'@title Plot matching heatmaps for modeled and observed temp -#'@param nc_file Netcdf model output file -#'@param field_file CSV or TSV field data file (see \link{resample_to_field} for format) +#' @title Plot matching heatmaps for modeled and observed temp +#' @param nc_file Netcdf model output file +#' @param field_file CSV or TSV field data file (see \link{resample_to_field} for format) #' @param fig_path F if plot to screen, string path if save plot as .png -#'@param \dots additional arguments passed to \code{\link{resample_to_field}} +#' @param resample sample the model output to the same time points as the observations? +#' @param col_lim range for heatmap (in units of the variable) +#' @param \dots additional arguments passed to \code{\link{resample_to_field}} #' #'@seealso Internally uses \code{\link{plot_var_compare}}, \code{\link{get_temp}} and \code{\link{resample_to_field}} #' @@ -16,8 +18,8 @@ #' #'plot_temp_compare(nc_file, field_file) ##makes a plot! #'@export -plot_temp_compare = function(nc_file, field_file, fig_path=FALSE, ...){ +plot_temp_compare = function(nc_file, field_file, fig_path=FALSE, resample=TRUE, col_lim, ...){ - plot_var_compare(nc_file, field_file, var_name='temp', fig_path=fig_path, ...) + plot_var_compare(nc_file, field_file, var_name='temp', fig_path=fig_path, resample=resample, col_lim=col_lim, ...) } diff --git a/R/plot_var.R b/R/plot_var.R index a5847e5..d542a81 100644 --- a/R/plot_var.R +++ b/R/plot_var.R @@ -3,6 +3,7 @@ #'@param var_name a character vector of valid variable names (see \code{\link{sim_vars}}) #'@param fig_path F if plot to screen, string path if save plot as .png #'@param reference 'surface' or 'bottom'. Only used for heatmap plots. +#' @param col_lim range for heatmap (in units of the variable) #'@param ... additional arguments passed to \code{par()} #'@keywords methods #'@seealso \code{\link{get_temp}}, \code{\link{sim_var_longname}}, @@ -25,7 +26,7 @@ #'fig_path = 'aed_out.png') #'} #'@export -plot_var <- function(file='output.nc', var_name, fig_path = F, reference='surface', ...){ +plot_var <- function(file='output.nc', var_name, fig_path = F, reference='surface', col_lim, ...){ heatmaps <- .is_heatmap(file, var_name) num_divs <- length(var_name) @@ -42,7 +43,7 @@ plot_var <- function(file='output.nc', var_name, fig_path = F, reference='surfac # iterate through plots for (j in 1:num_divs){ if (heatmaps[j]){ - .plot_nc_heatmap(file, var_name[j], reference) + .plot_nc_heatmap(file, var_name[j], reference, col_lim=col_lim) } else { .plot_nc_timeseries(file, var_name[j]) if(is_heatmap) .plot_null() # to fill up the colormap div diff --git a/R/plot_var_compare.R b/R/plot_var_compare.R index 2d3ffa5..27c665f 100644 --- a/R/plot_var_compare.R +++ b/R/plot_var_compare.R @@ -1,79 +1,83 @@ -#' @title Plot matching heatmaps for modeled and observed temp -#' @param nc_file Netcdf model output file -#' @param field_file CSV or TSV field data file (see \link{resample_to_field} for format) -#' @param var_name a character vector of valid variable names (see \code{\link{sim_vars}}) -#' @param fig_path F if plot to screen, string path if save plot as .png -#' @param resample sample the model output to the same time points as the observations? -#' @param \dots additional arguments passed to \code{\link{resample_to_field}} -#' -#' @seealso Internally uses \link{get_var} and \link{resample_to_field} -#' -#' -#'@examples -#'sim_folder <- run_example_sim(verbose = FALSE) -#'nc_file <- file.path(sim_folder, 'output.nc') -#'nml_file <- file.path(sim_folder, 'glm2.nml') -#'field_file <- file.path(sim_folder, 'field_data.tsv') -#' -#'run_glm(sim_folder) -#' -#'plot_var_compare(nc_file, field_file, 'temp', resample=FALSE) ##makes a plot! -#' -#'@importFrom akima interp -#'@export -plot_var_compare = function(nc_file, field_file, var_name, fig_path = FALSE, resample = TRUE, ...){ - - heatmaps <- .is_heatmap(nc_file, var_name) - if (!heatmaps){ - warning('plot_var_compare not implemented for 1D variables') - return() - } - - start_par = par(no.readonly = TRUE) - #Create layout - - mod_temp = get_var(nc_file, var_name, reference='surface') - mod_depths = get.offsets(mod_temp) - - - data = resample_to_field(nc_file, field_file, var_name=var_name, ...) - if(resample){ - model_df <- resample_sim(mod_temp, t_out = unique(data$DateTime)) - }else{ - model_df = mod_temp - } - - #Pivot observed into table - x = as.numeric(as.POSIXct(data$DateTime)) - y = data$Depth - z = data[,paste0('Observed_', var_name)] - x_out = sort(unique(x)) - y_out = sort(unique(mod_depths)) - - #remove any NA values before the 2D interp - x = x[!is.na(z)] - y = y[!is.na(z)] - z = z[!is.na(z)] - - #Added a scaling factor to Y. Interp won't interpolate if X and Y are on vastly different scales. - # I don't use Y from here later, so it doesn't matter what the mangitude of the values is. - interped = interp(x, y*1e6, z, x_out, y_out*1e6) - - gen_default_fig(filename=fig_path, num_divs=2)#, omi = c(0.1, 0.5, 0, 0)) - .stacked_layout(heatmaps, num_divs=2) - obs_df <- data.frame(interped$z) - names(obs_df) <- paste('var_', y_out, sep='') - obs_df <- cbind(data.frame(DateTime=as.POSIXct(x_out, origin='1970-01-01')), obs_df) - - #Use model to define X-axis plotting extent for both graphs - xaxis <- get_xaxis(model_df[,1]) - - y.text = y_out[1]+diff(range(y_out))*0.05 # note, reference will ALWAYS be surface for compare to field data - .plot_df_heatmap(obs_df, bar_title = .unit_label(nc_file,var_name), overlays=c(points(x=x,y=y),text(x_out[1],y=y.text,'Observed', pos=4, offset = 1)), xaxis=xaxis) - - .plot_df_heatmap(model_df, bar_title = .unit_label(nc_file,var_name), overlays=text(x_out[1],y=y.text,'Modeled', pos=4, offset = 1), xaxis=xaxis) - - par(start_par)#set PAR back to what it started at - if(is.character(fig_path)) - dev.off() +#' @title Plot matching heatmaps for modeled and observed temp +#' @param nc_file Netcdf model output file +#' @param field_file CSV or TSV field data file (see \link{resample_to_field} for format) +#' @param var_name a character vector of valid variable names (see \code{\link{sim_vars}}) +#' @param fig_path F if plot to screen, string path if save plot as .png +#' @param resample sample the model output to the same time points as the observations? +#' @param col_lim range for heatmap (in units of the variable) +#' @param \dots additional arguments passed to \code{\link{resample_to_field}} +#' +#' @seealso Internally uses \link{get_var} and \link{resample_to_field} +#' +#' +#'@examples +#'sim_folder <- run_example_sim(verbose = FALSE) +#'nc_file <- file.path(sim_folder, 'output.nc') +#'nml_file <- file.path(sim_folder, 'glm2.nml') +#'field_file <- file.path(sim_folder, 'field_data.tsv') +#' +#'run_glm(sim_folder) +#' +#'plot_var_compare(nc_file, field_file, 'temp', resample=FALSE) ##makes a plot! +#' +#'@importFrom akima interp +#'@export +plot_var_compare = function(nc_file, field_file, var_name, fig_path = FALSE, resample = TRUE, col_lim, ...){ + + heatmaps <- .is_heatmap(nc_file, var_name) + if (!heatmaps){ + warning('plot_var_compare not implemented for 1D variables') + return() + } + + start_par = par(no.readonly = TRUE) + #Create layout + + mod_temp = get_var(nc_file, var_name, reference='surface') + mod_depths = get.offsets(mod_temp) + + + data = resample_to_field(nc_file, field_file, var_name=var_name, ...) + if(resample){ + model_df <- resample_sim(mod_temp, t_out = unique(data$DateTime)) + }else{ + model_df = mod_temp + } + + #Pivot observed into table + x = as.numeric(as.POSIXct(data$DateTime)) + y = data$Depth + z = data[,paste0('Observed_', var_name)] + x_out = sort(unique(x)) + y_out = sort(unique(mod_depths)) + + #remove any NA values before the 2D interp + x = x[!is.na(z)] + y = y[!is.na(z)] + z = z[!is.na(z)] + + #Added a scaling factor to Y. Interp won't interpolate if X and Y are on vastly different scales. + # I don't use Y from here later, so it doesn't matter what the mangitude of the values is. + interped = interp(x, y*1e6, z, x_out, y_out*1e6) + + gen_default_fig(filename=fig_path, num_divs=2)#, omi = c(0.1, 0.5, 0, 0)) + .stacked_layout(heatmaps, num_divs=2) + obs_df <- data.frame(interped$z) + names(obs_df) <- paste('var_', y_out, sep='') + obs_df <- cbind(data.frame(DateTime=as.POSIXct(x_out, origin='1970-01-01')), obs_df) + + #Use model to define X-axis plotting extent for both graphs + xaxis <- get_xaxis(model_df[,1]) + + y.text = y_out[1]+diff(range(y_out))*0.05 # note, reference will ALWAYS be surface for compare to field data + if (missing(col_lim)) + col_lim = range(data[, -1], na.rm = TRUE) + + .plot_df_heatmap(obs_df, bar_title = .unit_label(nc_file,var_name), overlays=c(points(x=x,y=y),text(x_out[1],y=y.text,'Observed', pos=4, offset = 1)), xaxis=xaxis, col_lim=col_lim) + + .plot_df_heatmap(model_df, bar_title = .unit_label(nc_file,var_name), overlays=text(x_out[1],y=y.text,'Modeled', pos=4, offset = 1), xaxis=xaxis, col_lim=col_lim) + + par(start_par)#set PAR back to what it started at + if(is.character(fig_path)) + dev.off() } \ No newline at end of file diff --git a/R/timeseries_plots.R b/R/timeseries_plots.R index c5c8ec5..46bca28 100644 --- a/R/timeseries_plots.R +++ b/R/timeseries_plots.R @@ -1,5 +1,5 @@ -.plot_nc_heatmap <- function(file, var_name, reference, num_cells=100, palette){ +.plot_nc_heatmap <- function(file, var_name, reference, num_cells=100, palette, ...){ surface <- get_surface_height(file) max_depth <- max(surface[, 2]) @@ -7,21 +7,20 @@ z_out <- seq(min_depth, max_depth,length.out = num_cells) data <- get_var(file, z_out = z_out, var_name = var_name, reference = reference) title = .unit_label(file, var_name) - .plot_df_heatmap(data, title, num_cells, palette) + .plot_df_heatmap(data, title, num_cells, palette, ...) } -.plot_df_heatmap <- function(data, bar_title, num_cells, palette, title_prefix=NULL, overlays=NULL, xaxis=NULL, ...){ +.plot_df_heatmap <- function(data, bar_title, num_cells, palette, title_prefix=NULL, overlays=NULL, xaxis=NULL, col_lim){ z_out <- rLakeAnalyzer::get.offsets(data) reference = ifelse(substr(names(data)[2],1,3) == 'elv', 'bottom', 'surface') + if (missing(col_lim)) + col_lim = range(data[, -1], na.rm = TRUE) if (missing(palette)) palette <- colorRampPalette(c("violet","blue","cyan", "green3", "yellow", "orange", "red"), bias = 1, space = "rgb") - col_lim <- range(data[, -1], na.rm = TRUE) - - col_subs <- head(pretty(col_lim, 6), -1) levels <- sort(unique(c(col_subs, pretty(col_lim, 15)))) colors <- palette(n = length(levels)-1) @@ -35,7 +34,7 @@ plot_layout(xaxis, yaxis, add=T) .filled.contour(x = dates, y = z_out, z =matrix_var, levels= levels, - col=colors, ...) + col=colors) overlays # will plot any overlay functions axis_layout(xaxis, yaxis) #doing this after heatmap so the axis are on top diff --git a/README.md b/README.md index 70c6602..95f5c5c 100644 --- a/README.md +++ b/README.md @@ -6,12 +6,6 @@ glmtools Tools for interacting with the [General Lake Model (GLM)](http://aed.see.uwa.edu.au/research/models/GLM/ "General Lake Model's website") in R. `glmtools` includes some basic functions for calculating physical derivatives and thermal properties of model output, and some plotting functionality (see example image below). -![alt tag](http://github.gleon.io/images/test_figure.png) - -`glmtools`, as of v0.2.5, can also call GLM using the `GLMr` package. Shown here running GLM from R with example driver data that is part of the package: - -![alt tag](http://github.gleon.io/images/glm-r.png) - `glmtools` Functions (as of v0.2.5.2) ===== | Function | Title | diff --git a/man/plot_temp.Rd b/man/plot_temp.Rd index fbc98a0..f723588 100644 --- a/man/plot_temp.Rd +++ b/man/plot_temp.Rd @@ -5,7 +5,7 @@ \title{plot water temperatures from a GLM simulation} \usage{ plot_temp(file = "output.nc", fig_path = FALSE, reference = "surface", - ...) + col_lim, ...) } \arguments{ \item{file}{a string with the path to the netcdf output from GLM} @@ -14,6 +14,8 @@ plot_temp(file = "output.nc", fig_path = FALSE, reference = "surface", \item{reference}{a string for 'surface' or 'bottom'} +\item{col_lim}{range for heatmap (in units of the variable)} + \item{...}{additional arguments passed to \code{par()}} } \note{ diff --git a/man/plot_temp_compare.Rd b/man/plot_temp_compare.Rd index 92fbc18..45e2a88 100644 --- a/man/plot_temp_compare.Rd +++ b/man/plot_temp_compare.Rd @@ -4,7 +4,8 @@ \alias{plot_temp_compare} \title{Plot matching heatmaps for modeled and observed temp} \usage{ -plot_temp_compare(nc_file, field_file, fig_path = FALSE, ...) +plot_temp_compare(nc_file, field_file, fig_path = FALSE, resample = TRUE, + col_lim, ...) } \arguments{ \item{nc_file}{Netcdf model output file} @@ -13,6 +14,10 @@ plot_temp_compare(nc_file, field_file, fig_path = FALSE, ...) \item{fig_path}{F if plot to screen, string path if save plot as .png} +\item{resample}{sample the model output to the same time points as the observations?} + +\item{col_lim}{range for heatmap (in units of the variable)} + \item{\dots}{additional arguments passed to \code{\link{resample_to_field}}} } \examples{ diff --git a/man/plot_var.Rd b/man/plot_var.Rd index 2509a27..9b5d438 100644 --- a/man/plot_var.Rd +++ b/man/plot_var.Rd @@ -5,7 +5,7 @@ \title{plot variables from a GLM simulation} \usage{ plot_var(file = "output.nc", var_name, fig_path = F, - reference = "surface", ...) + reference = "surface", col_lim, ...) } \arguments{ \item{file}{a string with the path to the netcdf output from GLM} @@ -16,6 +16,8 @@ plot_var(file = "output.nc", var_name, fig_path = F, \item{reference}{'surface' or 'bottom'. Only used for heatmap plots.} +\item{col_lim}{range for heatmap (in units of the variable)} + \item{...}{additional arguments passed to \code{par()}} } \note{ diff --git a/man/plot_var_compare.Rd b/man/plot_var_compare.Rd index 3e7f969..b87c7db 100644 --- a/man/plot_var_compare.Rd +++ b/man/plot_var_compare.Rd @@ -5,7 +5,7 @@ \title{Plot matching heatmaps for modeled and observed temp} \usage{ plot_var_compare(nc_file, field_file, var_name, fig_path = FALSE, - resample = TRUE, ...) + resample = TRUE, col_lim, ...) } \arguments{ \item{nc_file}{Netcdf model output file} @@ -18,6 +18,8 @@ plot_var_compare(nc_file, field_file, var_name, fig_path = FALSE, \item{resample}{sample the model output to the same time points as the observations?} +\item{col_lim}{range for heatmap (in units of the variable)} + \item{\dots}{additional arguments passed to \code{\link{resample_to_field}}} } \examples{