Skip to content

Commit

Permalink
carb_chem and post_plot fn updates
Browse files Browse the repository at this point in the history
  • Loading branch information
dustintharper committed Jan 9, 2024
1 parent 1449309 commit c349ccd
Show file tree
Hide file tree
Showing 5 changed files with 74 additions and 43 deletions.
58 changes: 37 additions & 21 deletions R/carbchem.R
Original file line number Diff line number Diff line change
Expand Up @@ -9,49 +9,65 @@
#' @param ccparm2 User indicates which 2nd carbonate chemistry variable to use: 'dic', 'alk', 'pco2', 'co3', 'hco3', or 'omegac'.
#' Defaults to 'dic'.
#'
#' @param ccout Specify which parameter to compute: 'dic', 'alk', 'pco2', 'co3', 'hco3', or 'omegac'. Defaults to 'co2a'.
#' @param ccout Specify which parameter to compute: 'dic', 'alk', 'pco2', 'co3', 'hco3', or 'omegac'. Defaults to 'pco2'.
#'
#' @param inv_out Use inversion output, returned object from 'run_inversion' function, to compute carbonate chemistry. Otherwise, specify the
#' vectors with individual variable arguments. Defaults to NULL.
#'
#' @param pH Provide a vector of pH data in total scale. Defaults to 'pH'.
#' @param pH Provide a vector of pH data in total scale if 'inv_out' is not used. Defaults to NULL.
#'
#' @param ccparm2.vec Provide a vector of data of the type specified by 'ccparm2'. Defaults to 'ccparm2.vec'.
#' @param ccparm2.vec Provide a vector of data of the type specified by 'ccparm2' if 'inv_out' is not used. Defaults to NULL.
#'
#' @param temp Provide a vector of temperature data in degrees C. Defaults to 'temp'.
#' @param temp Provide a vector of temperature data in degrees C if 'inv_out' is not used. Defaults to NULL.
#'
#' @param sal Provide a vector of salinity data. Defaults to 'sal'.
#' @param sal Provide a vector of salinity data if 'inv_out' is not used. Defaults to NULL.
#'
#' @param press Provide a vector of pressure data in bar. Defaults to 'press'.
#' @param press Provide a vector of pressure data in bar if 'inv_out' is not used. Defaults to NULL.
#'
#' @param xca Provide a vector of Ca ion concentration data in mmol/kg. Defaults to 'xca'.
#' @param xca Provide a vector of Ca ion concentration data in mmol/kg if 'inv_out' is not used. Defaults to NULL.
#'
#' @param xmg Provide a vector of Mg ion concentration data in mmol/kg. Defaults to 'xmg'.
#' @param xmg Provide a vector of Mg ion concentration data in mmol/kg if 'inv_out' is not used. Defaults to NULL.
#'
#' @param xso4 Provide a vector of SO4 ion concentration data in mmol/kg. Defaults to 'xso4'.
#' @param xso4 Provide a vector of SO4 ion concentration data in mmol/kg if 'inv_out' is not used. Defaults to NULL.
#'
#' @returns Returns vector of values for parameter of user interest specified by 'ccout' argument.
#'
#' @examples
#' carbchem(ccparm2 = "dic", ccout = "pco2", pH = pH, ccparm2.vec = ccparm2.vec, temp = temp, sal = sal, press = press,
#' xca = xca, xmg = xmg, xso4 = xso4)
#' carbchem(ccparm2 = "dic", ccout = "alk", inv_out = inv_out)
#'
#' @export
carbchem <- function(ccparm2 = "dic",
ccout = "pco2",
pH = pH,
ccparm2.vec = ccparm2.vec,
temp = temp,
sal = sal,
press = press,
xca = xca,
xmg = xmg,
xso4 = xso4){
ccout = "alk",
inv_out,
pH,
ccparm2.vec,
temp,
sal,
press,
xca,
xmg,
xso4){

if(!(is.null(inv_out))){
io_med <- inv_out$jout$BUGSoutput$median
pH <- as.numeric(io_med$pH)
ccparm2.vec <- as.numeric(io_med[[ccparm2]])
temp <- as.numeric(io_med$tempC)
sal <- as.numeric(io_med$sal)
press <- rep_len(as.numeric(io_med$press), length.out = as.integer(length(sal)))
xca <- as.numeric(io_med$xca)
xmg <- as.numeric(io_med$xmg)
xso4 <- as.numeric(io_med$xso4)
}


if(length(pH) != length(ccparm2.vec) | length(pH) != length(temp) | length(pH) != length(sal) | length(pH) != length(press) |
length(pH) != length(xca) | length(pH) != length(xmg) | length(pH) != length(xso4)){
stop("Supplied vectors must be equal length")
}

ccparm1 <- "pH"
cc_in <- data.frame(c(pH, ccparm2.vec, temp, sal, press, xca, xmg, xso4))
cc_in <- data.frame(pH, ccparm2.vec, temp, sal, press, xca, xmg, xso4)
colnames(cc_in) <- c("pH","ccparm2.vec", "temp", "sal", "press", "xca", "xmg", "xso4")


Expand Down
10 changes: 9 additions & 1 deletion R/post_plot.R
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,9 @@
#' @param type Indicate which plot style desired. Either 'CI' for 95 percent credible interval, or 'draws' for 500 draws of time
#' series posteriors. Defaults to 'CI'.
#'
#' @param n.draws Indicate the number of posterior draws if something other than 500 is wanted and 'draws' was selected for
#' 'type' argument. Defaults to NULL.
#'
#' @param show.legend Logical. Specify TRUE if you would like a legend in the plot. Defaults to TRUE.
#'
#' @param leg.pos Option to include a character string to adjust the legend position. Position options are: 'bottomright',
Expand All @@ -25,6 +28,7 @@
post_plot <- function(inv_out = inv_out,
parm = "dic",
type = "CI",
n.draws,
show.legend = TRUE,
leg.pos = "topleft"){

Expand Down Expand Up @@ -71,10 +75,14 @@ post_plot <- function(inv_out = inv_out,
stop("Must input a time series parameter (i.e., 'parm' argument) that is in the 'save.parms' list used as argument for 'inv_out' function")
}

if(type =='draws' & is.null(n.draws)){
n.draws = 500
}

# generate time series plot (either random draws or confidence intervals) of parameter of interest
if(type == "draws"){
post_plot <- plot(ages_prox, parm_out[as.integer(runif(1,1,nrow(parm_out))),], type="l", xlab = "Age (kyr)", ylab = paste(parm, units), xlim = rev(range(ages_prox)), ylim = range(parm_out), col=rgb(red=0, green=0, blue=0, alpha=0.05), lwd=0.3)
for (i in as.integer(runif(500,1,nrow(parm_out)))){
for (i in as.integer(runif(n.draws,1,nrow(parm_out)))){
lines(ages_prox, parm_out[i,], col=rgb(red=0, green=0, blue=0, alpha=0.05), lwd=0.3)
}
lines(ages_prox, parm_med, col=rgb(red=1, green=0, blue=0), lwd=1.5)
Expand Down
2 changes: 1 addition & 1 deletion R/post_plot_ind.R
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,7 @@ post_plot_ind <- function(inv_out = inv_out,

if(ncol(parm_out) > 1 & !(tstep_age %in% ages_prox)){
tstep_age <- ages_prox[(which.min(abs(ages_prox-tstep_age)))]
warning("'tstep_age' is not a value in 'ages_prox'. The nearest time step to 'tstep_age'in 'ages_prox' vector has been used.")
warning("'tstep_age' is not a value in 'ages_prox'. The closest time step to 'tstep_age' in 'ages_prox' vector has been used.")
}

# set units for parameters
Expand Down
43 changes: 23 additions & 20 deletions man/carbchem.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 4 additions & 0 deletions man/post_plot.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

0 comments on commit c349ccd

Please sign in to comment.