simulate_LM_OSL.R at master ⋅ R-Lum/RLumModel

#' sequence step LM-OSL-simulation
#'
#' This function simulates the LM-OSL measurement of quartz in the energy-band-model.
#'
#' @param temp \code{\link{numeric}} (\bold{required}): set temperature [deg. C] of the LM-OSL simulation
#'
#' @param duration \code{\link{numeric}} (\bold{required}): duration of the LM-OSL simulation
#'
#' @param start_power \code{\link{numeric}}: % of the power density at the start of the measurement.
#'
#' @param end_power \code{\link{numeric}}: % of the power density at the end of the measurement.
#' 100 % equates to 20 mW/cm^2
#'
#' @param RLumModel_ID \code{\link{numeric}} (optional): A ID-number for the LM-OSL-step. This ID
#' is pass down to \link{calc_concentrations} so all concentrations had the same ID as the
#' sequence step they were calculated from. This ID is identic to the sequence step in "sequence".
#'
#' @param n \code{\link{numeric}} or \code{\linkS4class{RLum.Results}} (\bold{required}):
#' concentration of electron-/holetraps, valence- and conduction band
#' from step before. This is necessary to get the boundary condition for the ODEs.
#'
#' @param parms \code{\linkS4class{RLum.Results}} (\bold{required}): The specific model parameters are used to simulate
#' numerical quartz luminescence results.
#'
#' @return This function returns an RLum.Results object from the LM-OSL simulation.
#'
#' @section Function version: 0.1.1
#'
#' @author Johannes Friedrich, University of Bayreuth (Germany),
#'
#' @references
#'
#' @seealso \code{\link{plot}}
#'
#' @examples
#'
#' #so far no example available
#'
#' @noRd
.simulate_LM_OSL <- function(
  temp,
  duration,
  start_power = 0,
  end_power = 100,
  RLumModel_ID = NULL,
  n,
  parms
){

# check input arguments ---------------------------------------------------

  ##check if temperature is > 0 K (-273 degree celsius)
  if(temp < -273){
    stop("\n [.simulate_LM_OSL()] Argument 'temp' has to be > 0 K!")
  }

  ##check if duration is > 0s
  if(duration < 0){
    stop("\n [.simulate_LM_OSL()] Argument 'duration' has to be a positive number!")
  }

  ##check if start_power is > 0
  if(start_power < 0){
    stop("\n [.simulate_LM_OSL()] Argument 'start_power' has to be a positive number!")
  }

  ##check if end_power > start_power
  if(start_power > end_power){
    stop("\n [.simulate_LM_OSL()] Argument 'start_power' has to be smaller than 'end_power'!")
  }

  ##check if object is of class RLum.Data.Curve
  if(class(n) != "RLum.Results"){
    n <- n
  } else {
    n <- n$n
  }

# Set parameters for ODE ---------------------------------------------------


  ##============================================================================##
  # SETTING PARAMETERS FOR ILLUMINATION
  #
  # R: electron-hole-production-rate = 0
  # P: Photonflux (in Bailey 2004: wavelength [nm]) = 1
  # b: heating rate [deg. C/s] = 0
  # a: rate of stimulationintensity, P*20, because in Bailey2001 P = 1 equates to 20 mW cm^(-2)
  ##============================================================================##

  if(parms$model == "Bailey2004" || parms$model == "Bailey2002"){
    P <- 0.02/(1.6*10^(-19)*(1240/470))
    a <- 20*((end_power - start_power)/100)/duration
  }
  else{
    P <- 2
    a <- P*20*((end_power - start_power)/100)/duration
  }

  b <- 0
  R <- 0

  ##============================================================================##
  # SETTING PARAMETERS FOR ODE
  ##============================================================================##

  times <- seq(from = 0, to = duration, by = 0.1)
  parameters.step <- .extract_pars(parameters.step = list(
    a = a,
    R = R,
    P = P,
    temp = temp,
    b = b,
    times = times,
    parms = parms))
  
  ##============================================================================##
  # SOLVING ODE (deSolve requiered)
  ##============================================================================##
  out <- deSolve::ode(y = n, times = times, parms = parameters.step, func = .set_ODE_Rcpp_LM_OSL, rtol=1e-3, atol=1e-3, maxsteps=1e5, method = "bdf")
  ##============================================================================##

  ##============================================================================##
  # CALCULATING RESULTS FROM ODE SOLVING
  ##============================================================================##

  signal <- .calc_signal(object = out, parameters = parameters.step)

  ##============================================================================##
  # CALCULATING CONCENTRATIONS FROM ODE SOLVING
  ##============================================================================##

  name <- c("LM-OSL")
  concentrations <- .calc_concentrations(
    data = out,
    times = times,
    name = name,
    RLumModel_ID = RLumModel_ID)

  ##============================================================================##
  # TAKING THE LAST LINE OF "OUT" TO COMMIT IT TO THE NEXT STEP
  ##============================================================================##

  return(Luminescence::set_RLum(
                  class = "RLum.Results",
                  data = list(
                    n = out[length(times),-1],
                    LM_OSL.data = Luminescence::set_RLum(
                      class = "RLum.Data.Curve",
                      data = matrix(data = c(times[2:length(times)], signal[2:length(signal)]),ncol = 2),
                      recordType = "LM-OSL",
                      curveType = "simulated",
                      info = list(
                        curveDescripter = NA_character_),
                      .pid = as.character(RLumModel_ID)
                      ),
                    temp = temp,
                    concentrations = concentrations)
                  )
         )


}#end function

Read our documentation on viewing source code .

Navigation	Overlay
`t` Navigate files	`h` Toggle hits
`y` Change url to tip of branch	`m` Toggle misses
`b / v` Jump to prev/next hit line	`p` Toggle partial
`z / x` Jump to prev/next missed or partial line	`1..9` Toggle flags
`shift + o` Open current page in GitHub	`a` Toggle all on
`/ or ?` Show keyboard shortcuts dialog	`c` Toggle context lines or commits

1		#' sequence step LM-OSL-simulation
2		#'
3		#' This function simulates the LM-OSL measurement of quartz in the energy-band-model.
4		#'
5		#' @param temp \code{\link{numeric}} (\bold{required}): set temperature [deg. C] of the LM-OSL simulation
6		#'
7		#' @param duration \code{\link{numeric}} (\bold{required}): duration of the LM-OSL simulation
8		#'
9		#' @param start_power \code{\link{numeric}}: % of the power density at the start of the measurement.
10		#'
11		#' @param end_power \code{\link{numeric}}: % of the power density at the end of the measurement.
12		#' 100 % equates to 20 mW/cm^2
13		#'
14		#' @param RLumModel_ID \code{\link{numeric}} (optional): A ID-number for the LM-OSL-step. This ID
15		#' is pass down to \link{calc_concentrations} so all concentrations had the same ID as the
16		#' sequence step they were calculated from. This ID is identic to the sequence step in "sequence".
17		#'
18		#' @param n \code{\link{numeric}} or \code{\linkS4class{RLum.Results}} (\bold{required}):
19		#' concentration of electron-/holetraps, valence- and conduction band
20		#' from step before. This is necessary to get the boundary condition for the ODEs.
21		#'
22		#' @param parms \code{\linkS4class{RLum.Results}} (\bold{required}): The specific model parameters are used to simulate
23		#' numerical quartz luminescence results.
24		#'
25		#' @return This function returns an RLum.Results object from the LM-OSL simulation.
26		#'
27		#' @section Function version: 0.1.1
28		#'
29		#' @author Johannes Friedrich, University of Bayreuth (Germany),
30		#'
31		#' @references
32		#'
33		#' @seealso \code{\link{plot}}
34		#'
35		#' @examples
36		#'
37		#' #so far no example available
38		#'
39		#' @noRd
40		.simulate_LM_OSL <- function(
41		temp,
42		duration,
43		start_power = 0,
44		end_power = 100,
45		RLumModel_ID = NULL,
46		n,
47		parms
48		){
49
50		# check input arguments ---------------------------------------------------
51
52		##check if temperature is > 0 K (-273 degree celsius)
53	5	if(temp < -273){
54	5	stop("\n [.simulate_LM_OSL()] Argument 'temp' has to be > 0 K!")
55		}
56
57		##check if duration is > 0s
58	5	if(duration < 0){
59	5	stop("\n [.simulate_LM_OSL()] Argument 'duration' has to be a positive number!")
60		}
61
62		##check if start_power is > 0
63	5	if(start_power < 0){
64	5	stop("\n [.simulate_LM_OSL()] Argument 'start_power' has to be a positive number!")
65		}
66
67		##check if end_power > start_power
68	5	if(start_power > end_power){
69	5	stop("\n [.simulate_LM_OSL()] Argument 'start_power' has to be smaller than 'end_power'!")
70		}
71
72		##check if object is of class RLum.Data.Curve
73	5	if(class(n) != "RLum.Results"){
74	5	n <- n
75		} else {
76	5	n <- n$n
77		}
78
79		# Set parameters for ODE ---------------------------------------------------
80
81
82		##============================================================================##
83		# SETTING PARAMETERS FOR ILLUMINATION
84		#
85		# R: electron-hole-production-rate = 0
86		# P: Photonflux (in Bailey 2004: wavelength [nm]) = 1
87		# b: heating rate [deg. C/s] = 0
88		# a: rate of stimulationintensity, P*20, because in Bailey2001 P = 1 equates to 20 mW cm^(-2)
89		##============================================================================##
90
91	5	if(parms$model == "Bailey2004" \|\| parms$model == "Bailey2002"){
92	5	P <- 0.02/(1.610^(-19)(1240/470))
93	5	a <- 20*((end_power - start_power)/100)/duration
94		}
95		else{
96	5	P <- 2
97	5	a <- P20((end_power - start_power)/100)/duration
98		}
99
100	5	b <- 0
101	5	R <- 0
102
103		##============================================================================##
104		# SETTING PARAMETERS FOR ODE
105		##============================================================================##
106
107	5	times <- seq(from = 0, to = duration, by = 0.1)
108	5	parameters.step <- .extract_pars(parameters.step = list(
109	5	a = a,
110	5	R = R,
111	5	P = P,
112	5	temp = temp,
113	5	b = b,
114	5	times = times,
115	5	parms = parms))
116
117		##============================================================================##
118		# SOLVING ODE (deSolve requiered)
119		##============================================================================##
120	5	out <- deSolve::ode(y = n, times = times, parms = parameters.step, func = .set_ODE_Rcpp_LM_OSL, rtol=1e-3, atol=1e-3, maxsteps=1e5, method = "bdf")
121		##============================================================================##
122
123		##============================================================================##
124		# CALCULATING RESULTS FROM ODE SOLVING
125		##============================================================================##
126
127	5	signal <- .calc_signal(object = out, parameters = parameters.step)
128
129		##============================================================================##
130		# CALCULATING CONCENTRATIONS FROM ODE SOLVING
131		##============================================================================##
132
133	5	name <- c("LM-OSL")
134	5	concentrations <- .calc_concentrations(
135	5	data = out,
136	5	times = times,
137	5	name = name,
138	5	RLumModel_ID = RLumModel_ID)
139
140		##============================================================================##
141		# TAKING THE LAST LINE OF "OUT" TO COMMIT IT TO THE NEXT STEP
142		##============================================================================##
143
144	5	return(Luminescence::set_RLum(
145	5	class = "RLum.Results",
146	5	data = list(
147	5	n = out[length(times),-1],
148	5	LM_OSL.data = Luminescence::set_RLum(
149	5	class = "RLum.Data.Curve",
150	5	data = matrix(data = c(times[2:length(times)], signal[2:length(signal)]),ncol = 2),
151	5	recordType = "LM-OSL",
152	5	curveType = "simulated",
153	5	info = list(
154	5	curveDescripter = NA_character_),
155	5	.pid = as.character(RLumModel_ID)
156		),
157	5	temp = temp,
158	5	concentrations = concentrations)
159		)
160		)
161
162
163		}#end function