R/NumRate7BenpTable.R

In RafaelSdeSouza/nuclear: Functions to estimate S-factor for nuclear reactions

# Rafael de Souza, UNC
#This program is free software: you can redistribute it and/or modify
#it under the terms of the GNU General Public License version 3 as published by
#the Free Software Foundation.

#This program is distributed in the hope that it will be useful,
#but WITHOUT ANY WARRANTY; without even the implied warranty of
#MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#GNU General Public License for more details.

#  A copy of the GNU General Public License is available at
#  http://www.r-project.org/Licenses/
#
#' @title  Estimate reaction rates over a temperature grid
#' @description Provides a reaction rate table
#' @aliases NumRateHe3dpTable
#' @usage NumRateHe3dpTable (mat = mat, vars=vars, N = N, T9  =T9)
#' @format \describe{
#' \item{x}{
#' The function has 4 arguments: mat, vars, N, T9}
#' }
#' @param mat mat
#' @param vars  vars
#' @param N    N
#' @param T9  T9
#' @return NumRate7BenpTable
#' @import gsl
#' @author Rafael de Souza, UNC,  and Christian Illiadis, UNC
#'
#' @keywords NumRate7BenpTable
#' @export
#'
#'
#'
#'


NumRate7BenpTable <- function(mat, N = 1000,T9){

  index <- sample(1:nrow(mat),size=N,replace=FALSE)
  mcdat_I  <- mat[index,]

  gdat <-  sapply(Tgrid,function(Tgrid){apply(mcdat_I,1,NumRate7Benp,T9=Tgrid)})

  gg <-  as.data.frame(gdat)

  gg2 <- apply(gg, 2, quantile, probs=c(0.16, 0.5, 0.84), na.rm=TRUE)

  fu <- function(x){exp(sqrt(log(1+var(x)/mean(x)^2)))}

  fu_I<-apply(gg, 2, fu)

  gg2data <- data.frame(T9 =T9, lower = gg2["16%",], mean = gg2["50%",], upper = gg2["84%",] )
  gg2data$fu <- fu_I
  rownames(gg2data) <- c()
  return(gg2data)
}