R/wasa_modify_pars.R
In tpilz/WasaEchseTools: Title: Handy tools for preprocessing and analysis of the hydrological models WASA-SED and ECHSE (limited to the WASA engine). Currently only meant for internal use! If you want to try it, handle with care!
#' Modify input parameters of the WASA-SED model
#'
#' Function replaces specified parameter values within prepared WASA-SED model input
#' files by given substitutes.
#'
#' @param pars A named vector of type numeric with the values of selected parameters.
#'
#' @param dir_run Character string of the path to the WASA-SED input directory.
#'
#' @details IMPORTANT: Note that some parameters are factors or summands and the
#' parameters are updated by multiplication or summation rather than mere replacement.
#' Read notes below!
#'
#' Parameter names below can be modified by adding prefix 'log_*'. In this case,
#' the value will be replaced by exp(value) internally for better handling of highly
#' uncertain parameters. For instance, 'kf_bedrock_f' might vary between 0.01 and 100
#' which can cause some numerical issues to the calibration algorithm (i.e. underrepresenting
#' small numbers < 1). In such a case it is better to specify, e.g. log(100) or log(0.01)
#' (i.e. 4.6 or -4.6, respectively) as limits for the calibration algorithm. Note that
#' only the natural logarithm (R function log()) is supported herein!
#'
#'
#' Currently supported are the following parameters (use these as names for the input
#' vector \code{pars}):
#'
#' File Hillslope/vegetation.dat (equally applied to all vegetation types and node points):
#'
#' - stomr_f: Factor, somatal resistance\cr
#' - rootd_f: Factor, rootdepth\cr
#' - albedo_f: Factor, albedo
#'
#' File Hillslope/soil.dat (equally applied to all types):
#'
#' - n_f: Factor, porosity\cr
#' - cfr: Summand, added to coarse fragments\cr
#'
#' File Hillslope/frac_direct_gw.dat:
#'
#' - f_direct_gw: partitioning of groundwater into river and alluvia
#'
#' File River/response.dat:
#'
#' - uhg_f: Factor, applied to both lag and response times
#'
#' File Hillslope/soter.dat:
#'
#' - gw_delay_f: Factor, groundwater delay\cr
#' - soildepth_f: Factor, soil depth\cr
#' - kf_bedrock_f: Factor, kf of bedrock\cr
#' - riverdepth_f: Factor, depth of riverbed
#'
#' File Reservoir/lake_maxvol.dat:
#'
#' - f_lakemaxvol: Factor, water storage capacity for the reservoir size classes
#'
#' File Others/calib_wind.dat:
#'
#' - f_wind: Factor, modification of wind speed (hard-coded in WASA-SED; calibration of evapotranspiration)
#'
#' File do.dat:
#'
#' - kfcorr: Replacement of kfcorr value\cr
#' - kfcorr0: Replacement of kfcorr0 value\cr
#' - kfcorr_a: Replacement of kfcorr_a value\cr
#' - intcf: Replacement of intcf value\cr
#' - dosediment: Replacement of dosediment flag
#'
#' File calibration.dat:
#'
#' - ksat_factor: Replacement of ksat factor (factor internally used in WASA-SED)
#'
#'
#' @return Function returns nothing.
#'
#' @author Tobias Pilz \email{tpilz@@uni-potsdam.de}
#'
#' @export
wasa_modify_pars <- function(
pars = NULL,
dir_run = NULL
) {
# replace log values by normal values for use in WASA (log is just a matter of sampling)
if(any(grepl(names(pars), pattern = "^log_"))) {
log_trans = grepl(names(pars), pattern = "^log_") #find log-transformed parameters
pars[log_trans]=exp(pars[log_trans]) #transform back to non-log scale
names(pars) = sub(names(pars), pattern = "^log_", rep="") #remove "log_" from name
}
pars_t <- pars
# parameters in vegetation.dat
if (any(names(pars_t) %in% c("stomr_f", "rootd_f", "albedo_f"))) {
# file that hold the parameters to be changed
target_file <- paste(dir_run,"Hillslope/vegetation.dat",sep="/")
# read data
file_content <- read.table(target_file, skip=2, header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply somatal resistance in vegetation.dat by factor
nn= which(names(pars_t)=="stomr_f")
if (length(nn)>0) {
file_content[,2]=file_content[,2]*pars_t[nn]
pars_t <- pars_t[-nn]
}
# multiply rootdepth in vegetation.dat by factor (all 4 values)
nn= which(names(pars_t)=="rootd_f")
if (length(nn)>0) {
file_content[,c(9:12)]=file_content[,c(9:12)]*pars_t[nn]
pars_t <- pars_t[-nn]
}
# multiply albedo in vegetation.dat by factor (all 4 values)
nn= which(names(pars_t)=="albedo_f")
if (length(nn)>0) {
file_content[,c(17:20)]=pmin(1,file_content[,c(17:20)]*pars_t[nn])
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of vegetation parameters\nVeg-ID Stomata_Resistance[s/m] minsuction[hPa] maxsuction[hPa] height1[m] height2[m] height3[m] height4[m] rootdepth1[m] rootdepth2[m] rootdepth3[m] rootdepth4[m] LAI1[-] LAI2[-] LAI3[-] LAI4[-] albedo1[-] albedo2[-] albedo3[-] albedo4[-]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content, 2), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# parameters in soil.dat
if (any(names(pars_t) %in% c("n_f", "cfr"))) {
# file that hold the parameters to be changed
target_file=paste(dir_run,"Hillslope/soil.dat",sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply porosity in soil.dat by factor (equally for every horizon and soil)
nn= which(names(pars_t)=="n_f")
if (length(nn)>0) {
# iterate over soils (lines in data)
for (s in 1:nrow(file_content)) {
# get number of horizons
nhor = file_content[s,2]
# update parameter
file_content[s,12+(1:nhor-1)*13] = pmin(1,file_content[s,12+(1:nhor-1)*13]*pars_t[nn])
}
pars_t <- pars_t[-nn]
}
# add 'cfr' to coarse fragments in soil.dat (additive!, equally for every horizon and soil)
nn= which(names(pars_t)=="cfr")
if (length(nn)>0) {
# iterate over soils (lines in data)
for (s in 1:nrow(file_content)) {
# get number of horizons
nhor = file_content[s,2]
# update parameter
file_content[s,14+(1:nhor-1)*13] = pmin(1,pmax(0,file_content[s,14+(1:nhor-1)*13]+pars_t[nn]))
}
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of soil parameters\nSoil-ID[-] number(horizons)[-] (n_res[Vol-] n_PWP[-] n_FK2.6[-] n_FK1.8[-] n_nFK[-] n_saturated[-] n_thickness[mm] n_ks[mm/d] n_suction[mm] n_pore-size-index[-] n_bubblepressure[cm] n_coarse_frag[-]*n n_shrinks[0/1]) bedrock[0/1] alluvial[0/1]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content,3), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# params in frac_direct_gw.dat
if (any(names(pars_t) %in% c("f_direct_gw"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Hillslope/frac_direct_gw.dat",sep="/")
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t")
nn= which(names(pars_t)=="f_direct_gw")
if (length(nn)>0) {
file_content=pars_t[nn]
pars_t <- pars_t[-nn]
}
write.table(round(file_content,3), file = target_file, quote = F,row.names=F,col.names=F,sep="\t")
}
# parameters in response.dat
if (any(names(pars_t) %in% c("uhg_f"))){
# file that hold the parameters to be changed
target_file=paste(dir_run,"River/response.dat",sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply lag and response times in response.dat by factor
nn= which(names(pars_t)=="uhg_f")
if (length(nn)>0) {
file_content[,2:3] <- file_content[,2:3]*pars_t[nn]
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of routing parameter\nSubbasin-ID lag time [d] retention [d]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content,2), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# parameters in soter.dat
if (any(names(pars_t) %in% c("gw_delay_f","soildepth_f","kf_bedrock_f","riverdepth_f"))) {
# file that hold the parameters to be changed
target_file=paste(dir_run,"Hillslope/soter.dat", sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# consider shorter lines (LUs with less TCs) and bring fields to consistent position in matrix
max_n_tcs=max(file_content[,2])
shorter_lines=which(file_content[,2] != max_n_tcs)
n_fields=ncol(file_content)-max_n_tcs -2 #number of fields after specification of TC-IDs
if(length(shorter_lines) > 0) {
for (ll in shorter_lines) {
n_tcs = file_content[ll,2]
file_content[ll, 2+max_n_tcs+(1:n_fields)] = file_content[ll, 2+n_tcs +(1:n_fields)]
}
}
nn= which(names(pars_t)=="gw_delay_f")
if (length(nn)>0) {
file_content[,ncol(file_content)]=file_content[,ncol(file_content)]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="soildepth_f")
if (length(nn)>0) {
not_minus1=file_content[,ncol(file_content)-4]!=-1 #only modify entries not having flag=-1
file_content[not_minus1,ncol(file_content)-4]=file_content[not_minus1,ncol(file_content)-4]*pars_t[nn]
not_minus1=file_content[,ncol(file_content)-5]!=-1
file_content[not_minus1,ncol(file_content)-5]=file_content[not_minus1,ncol(file_content)-5]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="kf_bedrock_f")
if (length(nn)>0) {
file_content[,ncol(file_content)-7]=file_content[,ncol(file_content)-7]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="riverdepth_f")
if (length(nn)>0) {
file_content[,ncol(file_content)-3]=file_content[,ncol(file_content)-3]*pars_t[nn]
pars_t <- pars_t[-nn]
}
#consider shorter lines (LUs with less TCs) and bring fields to "spars_te" representation (unequal number of fields)
if(length(shorter_lines) > 0) {
for (ll in shorter_lines) {
file_content[ll,]
n_tcs = file_content[ll,2]
file_content[ll, 2+n_tcs +(1:n_fields)] =
file_content[ll, 2+max_n_tcs+(1:n_fields)]
file_content[ll, ncol(file_content)+1-1:(max_n_tcs-n_tcs)] = NA #mask obsolete fields with NA
}
}
#re-write file
content=paste("Specification of landscape units\nLU-ID[id] No._of_TC[-] TC1[id] TC2[id] TC3[id] kfsu[mm/d] length[m] meandep[mm] maxdep[mm] riverbed[mm] gwflag[0/1] gw_dist[mm] frgw_delay[day]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content, 1), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# volumes in lake.dat / lake_maxvol.dat
if (any(names(pars_t) %in% c("f_lakemaxvol"))) {
nn= which(names(pars_t)=="f_lakemaxvol")
# file that holds the parameters to be changed
target_file=paste(dir_run,"Reservoir/lake_maxvol.dat",sep="/")
if(file.exists(target_file)) {
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t", skip=2)
# update parameters
file_content[,2:ncol(file_content)] <- file_content[,2:ncol(file_content)] * pars_t[nn]
# write updated parameters
write(c("Specification of water storage capacity for the reservoir size classes",
"Sub-basin-ID, maxlake[m**3] (five reservoir size classes)"),
file = target_file, sep="\n")
write.table(round(file_content, 2), target_file, append = T, row.names=F, col.names=F, quote=F, sep="\t")
}
# file that holds the parameters to be changed
target_file=paste(dir_run,"Reservoir/lake.dat",sep="/")
if(file.exists(target_file)) {
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t", skip=2)
# update parameters
file_content[,2] <- file_content[,2] * pars_t[nn]
# write updated parameters
write(c("Specification of parameters for the reservoir size classes",
"Reservoir_class-ID, maxlake0[m**3], lake_vol0_factor[-], lake_change[-], alpha_Molle[-], damk_Molle[-], damc_hrr[-], damd_hrr[-]"),
file = target_file, sep="\n")
write.table(round(file_content, 2), target_file, append = T, row.names=F, col.names=F, quote=F, sep="\t")
}
pars_t <- pars_t[-nn]
}
# param in calib_wind.dat
if (any(names(pars_t) %in% c("f_wind"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Others/calib_wind.dat",sep="/")
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t")
nn= which(names(pars_t)=="f_wind")
if (length(nn)>0) {
file_content=pars_t[nn]
pars_t <- pars_t[-nn]
}
write.table(round(file_content,3), file = target_file, quote = F,row.names=F,col.names=F,sep="\t")
}
# params in do.dat
if (any(names(pars_t) %in% c("kfcorr","kfcorr0","kfcorr_a", "intcf", "dosediment"))) {
# adjust kfcorr in do.dat
target_file=paste(dir_run,"do.dat",sep="/") #file that hold the parameters to be changed
# read data
file_content = read.table(target_file, skip=0,header = FALSE, sep = "$",stringsAsFactors=FALSE)
nn= which(names(pars_t)=="kfcorr")
if (length(nn)>0) {
file_content[24,1]=paste(round(pars_t[nn],2)," //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr)",sep="")
pars_t <- pars_t[-nn]
}
if(any(names(pars_t)=="kfcorr0") || any(names(pars_t)=="kfcorr_a")) {
nn= which(names(pars_t)=="kfcorr0")
if (length(nn)>0) {
file_content[24,1]=paste(round(pars_t[nn],2)," ",sep="")
pars_t <- pars_t[-nn]
} else { # default value
file_content[24,1]="1 "
}
nn= which(names(pars_t)=="kfcorr_a")
if (length(nn)>0) {
file_content[24,1]=paste(file_content[24,1],round(pars_t[nn],2)," 0 //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr0) [optional: a <tab> b for kfcorr=kfcorr0*(a*1/daily_precip+b+1) ",sep="")
pars_t <- pars_t[-nn]
} else { # default value
file_content[24,1]=paste(file_content[24,1],"0 0 //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr0) [optional: a <tab> b for kfcorr=kfcorr0*(a*1/daily_precip+b+1) ",sep="")
}
}
nn= which(names(pars_t)=="intcf")
if (length(nn)>0) {
file_content[25,1]=paste(round(pars_t[nn], 2)," //intcf: interception capacity per unit LAI (mm)",sep="")
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="dosediment")
if (length(nn)>0) {
file_content[31,1]=paste(pars_t[nn]," //dosediment")
pars_t <- pars_t[-nn]
}
#rewrite file
write.table(file_content, file = target_file, append = F, quote = F,row.names=F,col.names=F,sep="\t")
}
# params in calibration.dat
if (any(names(pars_t) %in% c("ksat_factor"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Others/calibration.dat",sep="/")
# read data
file_content = read.table(target_file, skip=1,row.names=NULL, header = FALSE, sep = "\t", dec = ".")
nn= which(names(pars_t)=="ksat_factor")
if (length(nn)>0) {
file_content[,2]=pars_t[nn]
pars_t <- pars_t[-nn]
}
content=paste("Soil-ID","Ksat-calib-factor",sep="\t")
write(content, file = target_file) #write header
write.table(round(file_content,3), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t")
}
# all parameters used?
if(length(pars_t) > 0)
stop(paste("The following parameter(s) was/were not used (no instruction implemented):", paste(names(pars_t), collapse = ", ")))
} # EOF
tpilz/WasaEchseTools documentation built on May 5, 2019, 12:33 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
#' Modify input parameters of the WASA-SED model
#'
#' Function replaces specified parameter values within prepared WASA-SED model input
#' files by given substitutes.
#'
#' @param pars A named vector of type numeric with the values of selected parameters.
#'
#' @param dir_run Character string of the path to the WASA-SED input directory.
#'
#' @details IMPORTANT: Note that some parameters are factors or summands and the
#' parameters are updated by multiplication or summation rather than mere replacement.
#' Read notes below!
#'
#' Parameter names below can be modified by adding prefix 'log_*'. In this case,
#' the value will be replaced by exp(value) internally for better handling of highly
#' uncertain parameters. For instance, 'kf_bedrock_f' might vary between 0.01 and 100
#' which can cause some numerical issues to the calibration algorithm (i.e. underrepresenting
#' small numbers < 1). In such a case it is better to specify, e.g. log(100) or log(0.01)
#' (i.e. 4.6 or -4.6, respectively) as limits for the calibration algorithm. Note that
#' only the natural logarithm (R function log()) is supported herein!
#'
#'
#' Currently supported are the following parameters (use these as names for the input
#' vector \code{pars}):
#'
#' File Hillslope/vegetation.dat (equally applied to all vegetation types and node points):
#'
#' - stomr_f: Factor, somatal resistance\cr
#' - rootd_f: Factor, rootdepth\cr
#' - albedo_f: Factor, albedo
#'
#' File Hillslope/soil.dat (equally applied to all types):
#'
#' - n_f: Factor, porosity\cr
#' - cfr: Summand, added to coarse fragments\cr
#'
#' File Hillslope/frac_direct_gw.dat:
#'
#' - f_direct_gw: partitioning of groundwater into river and alluvia
#'
#' File River/response.dat:
#'
#' - uhg_f: Factor, applied to both lag and response times
#'
#' File Hillslope/soter.dat:
#'
#' - gw_delay_f: Factor, groundwater delay\cr
#' - soildepth_f: Factor, soil depth\cr
#' - kf_bedrock_f: Factor, kf of bedrock\cr
#' - riverdepth_f: Factor, depth of riverbed
#'
#' File Reservoir/lake_maxvol.dat:
#'
#' - f_lakemaxvol: Factor, water storage capacity for the reservoir size classes
#'
#' File Others/calib_wind.dat:
#'
#' - f_wind: Factor, modification of wind speed (hard-coded in WASA-SED; calibration of evapotranspiration)
#'
#' File do.dat:
#'
#' - kfcorr: Replacement of kfcorr value\cr
#' - kfcorr0: Replacement of kfcorr0 value\cr
#' - kfcorr_a: Replacement of kfcorr_a value\cr
#' - intcf: Replacement of intcf value\cr
#' - dosediment: Replacement of dosediment flag
#'
#' File calibration.dat:
#'
#' - ksat_factor: Replacement of ksat factor (factor internally used in WASA-SED)
#'
#'
#' @return Function returns nothing.
#'
#' @author Tobias Pilz \email{tpilz@@uni-potsdam.de}
#'
#' @export
wasa_modify_pars <- function(
pars = NULL,
dir_run = NULL
) {
# replace log values by normal values for use in WASA (log is just a matter of sampling)
if(any(grepl(names(pars), pattern = "^log_"))) {
log_trans = grepl(names(pars), pattern = "^log_") #find log-transformed parameters
pars[log_trans]=exp(pars[log_trans]) #transform back to non-log scale
names(pars) = sub(names(pars), pattern = "^log_", rep="") #remove "log_" from name
}
pars_t <- pars
# parameters in vegetation.dat
if (any(names(pars_t) %in% c("stomr_f", "rootd_f", "albedo_f"))) {
# file that hold the parameters to be changed
target_file <- paste(dir_run,"Hillslope/vegetation.dat",sep="/")
# read data
file_content <- read.table(target_file, skip=2, header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply somatal resistance in vegetation.dat by factor
nn= which(names(pars_t)=="stomr_f")
if (length(nn)>0) {
file_content[,2]=file_content[,2]*pars_t[nn]
pars_t <- pars_t[-nn]
}
# multiply rootdepth in vegetation.dat by factor (all 4 values)
nn= which(names(pars_t)=="rootd_f")
if (length(nn)>0) {
file_content[,c(9:12)]=file_content[,c(9:12)]*pars_t[nn]
pars_t <- pars_t[-nn]
}
# multiply albedo in vegetation.dat by factor (all 4 values)
nn= which(names(pars_t)=="albedo_f")
if (length(nn)>0) {
file_content[,c(17:20)]=pmin(1,file_content[,c(17:20)]*pars_t[nn])
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of vegetation parameters\nVeg-ID Stomata_Resistance[s/m] minsuction[hPa] maxsuction[hPa] height1[m] height2[m] height3[m] height4[m] rootdepth1[m] rootdepth2[m] rootdepth3[m] rootdepth4[m] LAI1[-] LAI2[-] LAI3[-] LAI4[-] albedo1[-] albedo2[-] albedo3[-] albedo4[-]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content, 2), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# parameters in soil.dat
if (any(names(pars_t) %in% c("n_f", "cfr"))) {
# file that hold the parameters to be changed
target_file=paste(dir_run,"Hillslope/soil.dat",sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply porosity in soil.dat by factor (equally for every horizon and soil)
nn= which(names(pars_t)=="n_f")
if (length(nn)>0) {
# iterate over soils (lines in data)
for (s in 1:nrow(file_content)) {
# get number of horizons
nhor = file_content[s,2]
# update parameter
file_content[s,12+(1:nhor-1)*13] = pmin(1,file_content[s,12+(1:nhor-1)*13]*pars_t[nn])
}
pars_t <- pars_t[-nn]
}
# add 'cfr' to coarse fragments in soil.dat (additive!, equally for every horizon and soil)
nn= which(names(pars_t)=="cfr")
if (length(nn)>0) {
# iterate over soils (lines in data)
for (s in 1:nrow(file_content)) {
# get number of horizons
nhor = file_content[s,2]
# update parameter
file_content[s,14+(1:nhor-1)*13] = pmin(1,pmax(0,file_content[s,14+(1:nhor-1)*13]+pars_t[nn]))
}
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of soil parameters\nSoil-ID[-] number(horizons)[-] (n_res[Vol-] n_PWP[-] n_FK2.6[-] n_FK1.8[-] n_nFK[-] n_saturated[-] n_thickness[mm] n_ks[mm/d] n_suction[mm] n_pore-size-index[-] n_bubblepressure[cm] n_coarse_frag[-]*n n_shrinks[0/1]) bedrock[0/1] alluvial[0/1]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content,3), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# params in frac_direct_gw.dat
if (any(names(pars_t) %in% c("f_direct_gw"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Hillslope/frac_direct_gw.dat",sep="/")
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t")
nn= which(names(pars_t)=="f_direct_gw")
if (length(nn)>0) {
file_content=pars_t[nn]
pars_t <- pars_t[-nn]
}
write.table(round(file_content,3), file = target_file, quote = F,row.names=F,col.names=F,sep="\t")
}
# parameters in response.dat
if (any(names(pars_t) %in% c("uhg_f"))){
# file that hold the parameters to be changed
target_file=paste(dir_run,"River/response.dat",sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# multiply lag and response times in response.dat by factor
nn= which(names(pars_t)=="uhg_f")
if (length(nn)>0) {
file_content[,2:3] <- file_content[,2:3]*pars_t[nn]
pars_t <- pars_t[-nn]
}
#re-write file
content=paste("Specification of routing parameter\nSubbasin-ID lag time [d] retention [d]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content,2), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# parameters in soter.dat
if (any(names(pars_t) %in% c("gw_delay_f","soildepth_f","kf_bedrock_f","riverdepth_f"))) {
# file that hold the parameters to be changed
target_file=paste(dir_run,"Hillslope/soter.dat", sep="/")
# read data
file_content = read.table(target_file, skip=2,header = FALSE, sep = "\t", dec = ".", fill = TRUE)
if (all(!is.finite(file_content[,ncol(file_content)]))) file_content[,ncol(file_content)]=NULL #discard last column if empty
# consider shorter lines (LUs with less TCs) and bring fields to consistent position in matrix
max_n_tcs=max(file_content[,2])
shorter_lines=which(file_content[,2] != max_n_tcs)
n_fields=ncol(file_content)-max_n_tcs -2 #number of fields after specification of TC-IDs
if(length(shorter_lines) > 0) {
for (ll in shorter_lines) {
n_tcs = file_content[ll,2]
file_content[ll, 2+max_n_tcs+(1:n_fields)] = file_content[ll, 2+n_tcs +(1:n_fields)]
}
}
nn= which(names(pars_t)=="gw_delay_f")
if (length(nn)>0) {
file_content[,ncol(file_content)]=file_content[,ncol(file_content)]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="soildepth_f")
if (length(nn)>0) {
not_minus1=file_content[,ncol(file_content)-4]!=-1 #only modify entries not having flag=-1
file_content[not_minus1,ncol(file_content)-4]=file_content[not_minus1,ncol(file_content)-4]*pars_t[nn]
not_minus1=file_content[,ncol(file_content)-5]!=-1
file_content[not_minus1,ncol(file_content)-5]=file_content[not_minus1,ncol(file_content)-5]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="kf_bedrock_f")
if (length(nn)>0) {
file_content[,ncol(file_content)-7]=file_content[,ncol(file_content)-7]*pars_t[nn]
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="riverdepth_f")
if (length(nn)>0) {
file_content[,ncol(file_content)-3]=file_content[,ncol(file_content)-3]*pars_t[nn]
pars_t <- pars_t[-nn]
}
#consider shorter lines (LUs with less TCs) and bring fields to "spars_te" representation (unequal number of fields)
if(length(shorter_lines) > 0) {
for (ll in shorter_lines) {
file_content[ll,]
n_tcs = file_content[ll,2]
file_content[ll, 2+n_tcs +(1:n_fields)] =
file_content[ll, 2+max_n_tcs+(1:n_fields)]
file_content[ll, ncol(file_content)+1-1:(max_n_tcs-n_tcs)] = NA #mask obsolete fields with NA
}
}
#re-write file
content=paste("Specification of landscape units\nLU-ID[id] No._of_TC[-] TC1[id] TC2[id] TC3[id] kfsu[mm/d] length[m] meandep[mm] maxdep[mm] riverbed[mm] gwflag[0/1] gw_dist[mm] frgw_delay[day]",sep = "")
write(content, file = target_file) #write header
write.table(round(file_content, 1), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t", na = "")
}
# volumes in lake.dat / lake_maxvol.dat
if (any(names(pars_t) %in% c("f_lakemaxvol"))) {
nn= which(names(pars_t)=="f_lakemaxvol")
# file that holds the parameters to be changed
target_file=paste(dir_run,"Reservoir/lake_maxvol.dat",sep="/")
if(file.exists(target_file)) {
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t", skip=2)
# update parameters
file_content[,2:ncol(file_content)] <- file_content[,2:ncol(file_content)] * pars_t[nn]
# write updated parameters
write(c("Specification of water storage capacity for the reservoir size classes",
"Sub-basin-ID, maxlake[m**3] (five reservoir size classes)"),
file = target_file, sep="\n")
write.table(round(file_content, 2), target_file, append = T, row.names=F, col.names=F, quote=F, sep="\t")
}
# file that holds the parameters to be changed
target_file=paste(dir_run,"Reservoir/lake.dat",sep="/")
if(file.exists(target_file)) {
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t", skip=2)
# update parameters
file_content[,2] <- file_content[,2] * pars_t[nn]
# write updated parameters
write(c("Specification of parameters for the reservoir size classes",
"Reservoir_class-ID, maxlake0[m**3], lake_vol0_factor[-], lake_change[-], alpha_Molle[-], damk_Molle[-], damc_hrr[-], damd_hrr[-]"),
file = target_file, sep="\n")
write.table(round(file_content, 2), target_file, append = T, row.names=F, col.names=F, quote=F, sep="\t")
}
pars_t <- pars_t[-nn]
}
# param in calib_wind.dat
if (any(names(pars_t) %in% c("f_wind"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Others/calib_wind.dat",sep="/")
# read data
file_content = read.table(target_file, header = FALSE, sep = "\t")
nn= which(names(pars_t)=="f_wind")
if (length(nn)>0) {
file_content=pars_t[nn]
pars_t <- pars_t[-nn]
}
write.table(round(file_content,3), file = target_file, quote = F,row.names=F,col.names=F,sep="\t")
}
# params in do.dat
if (any(names(pars_t) %in% c("kfcorr","kfcorr0","kfcorr_a", "intcf", "dosediment"))) {
# adjust kfcorr in do.dat
target_file=paste(dir_run,"do.dat",sep="/") #file that hold the parameters to be changed
# read data
file_content = read.table(target_file, skip=0,header = FALSE, sep = "$",stringsAsFactors=FALSE)
nn= which(names(pars_t)=="kfcorr")
if (length(nn)>0) {
file_content[24,1]=paste(round(pars_t[nn],2)," //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr)",sep="")
pars_t <- pars_t[-nn]
}
if(any(names(pars_t)=="kfcorr0") || any(names(pars_t)=="kfcorr_a")) {
nn= which(names(pars_t)=="kfcorr0")
if (length(nn)>0) {
file_content[24,1]=paste(round(pars_t[nn],2)," ",sep="")
pars_t <- pars_t[-nn]
} else { # default value
file_content[24,1]="1 "
}
nn= which(names(pars_t)=="kfcorr_a")
if (length(nn)>0) {
file_content[24,1]=paste(file_content[24,1],round(pars_t[nn],2)," 0 //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr0) [optional: a <tab> b for kfcorr=kfcorr0*(a*1/daily_precip+b+1) ",sep="")
pars_t <- pars_t[-nn]
} else { # default value
file_content[24,1]=paste(file_content[24,1],"0 0 //kfcorr: hydraulic conductivity factor (for daily model version) (kfcorr0) [optional: a <tab> b for kfcorr=kfcorr0*(a*1/daily_precip+b+1) ",sep="")
}
}
nn= which(names(pars_t)=="intcf")
if (length(nn)>0) {
file_content[25,1]=paste(round(pars_t[nn], 2)," //intcf: interception capacity per unit LAI (mm)",sep="")
pars_t <- pars_t[-nn]
}
nn= which(names(pars_t)=="dosediment")
if (length(nn)>0) {
file_content[31,1]=paste(pars_t[nn]," //dosediment")
pars_t <- pars_t[-nn]
}
#rewrite file
write.table(file_content, file = target_file, append = F, quote = F,row.names=F,col.names=F,sep="\t")
}
# params in calibration.dat
if (any(names(pars_t) %in% c("ksat_factor"))) {
# file that holds the parameters to be changed
target_file=paste(dir_run,"Others/calibration.dat",sep="/")
# read data
file_content = read.table(target_file, skip=1,row.names=NULL, header = FALSE, sep = "\t", dec = ".")
nn= which(names(pars_t)=="ksat_factor")
if (length(nn)>0) {
file_content[,2]=pars_t[nn]
pars_t <- pars_t[-nn]
}
content=paste("Soil-ID","Ksat-calib-factor",sep="\t")
write(content, file = target_file) #write header
write.table(round(file_content,3), file = target_file, append = TRUE, quote = F,row.names=F,col.names=F,sep="\t")
}
# all parameters used?
if(length(pars_t) > 0)
stop(paste("The following parameter(s) was/were not used (no instruction implemented):", paste(names(pars_t), collapse = ", ")))
} # EOF
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.