From be8eaa3bf8c1ec4e08511d548f3e3ef6cec01631 Mon Sep 17 00:00:00 2001 From: Johannes Oberpriller Date: Thu, 7 Jan 2021 12:50:20 +0100 Subject: [PATCH 1/3] major update to be able to run all model version --- rLPJGUESS/R/AdjustTemplates.R | 101 ++++++ rLPJGUESS/R/GetRunAbleParameters.R | 48 +++ rLPJGUESS/R/InferParameterAndDesignList.R | 399 ++++++++++++++++++++++ rLPJGUESS/R/callLPJ.R | 6 +- rLPJGUESS/R/checkDesign.R | 2 +- rLPJGUESS/R/checkParameters.R | 20 +- rLPJGUESS/R/createRunParameters.R | 25 +- rLPJGUESS/R/createSingleObject.R | 66 ++-- rLPJGUESS/R/getData.R | 129 ++++--- rLPJGUESS/R/getDesign.R | 8 +- rLPJGUESS/R/getParameterList.R | 5 +- rLPJGUESS/R/plotLPJData.R | 125 ++++--- rLPJGUESS/R/readTableHeaderLPJ.R | 2 +- rLPJGUESS/R/runLPJ.R | 22 +- rLPJGUESS/R/runLPJWrapper.R | 60 ++-- rLPJGUESS/R/sysdata.rda | Bin 15417 -> 7364 bytes rLPJGUESS/R/writeTemplate.R | 49 +-- 17 files changed, 818 insertions(+), 249 deletions(-) create mode 100644 rLPJGUESS/R/AdjustTemplates.R create mode 100644 rLPJGUESS/R/GetRunAbleParameters.R create mode 100644 rLPJGUESS/R/InferParameterAndDesignList.R diff --git a/rLPJGUESS/R/AdjustTemplates.R b/rLPJGUESS/R/AdjustTemplates.R new file mode 100644 index 0000000..8f9e0ed --- /dev/null +++ b/rLPJGUESS/R/AdjustTemplates.R @@ -0,0 +1,101 @@ +#' @title This function writes the two main templates in the correct order +#' @description This function takes the defaultparameters and a defaultlist of +#' of files and writes the main and pft instruction file into the defined destinations +#' which are MainTemplateDestination and PftTemplateDestination +#' @param defaultparameters a matrix with the defaultparameters which should be +#' used in rlpjguess runs +#' @param defaultlist a matrix with the design parameters which should be used +#' in rlpjguess runs +#' @param MainTemplateDestination a character string indicating the space where +#' one should write the main file +#' @param PftTemplateDestination a character string indicating the space where +#' one should write the pft file +#' @param NameMainFile a character string with the name of the main file +#' which was produced with the default parameters before +#' @param NamePftFile a character string with the name of the pft file +#' which was produced with the default parameters before +#' @return writes the two files which serve as templates +#' @author Johannes Oberpriller +#' @export + + +AdjustTemplates <- function(defaultparameters,defaultlist, + MainTemplateDestination = NULL, + PftTemplateDestination = NULL, + NameMainFile = NULL, + NamePftFile = NULL){ + + # check if there is a destination to save the main template + + if(is.null(MainTemplateDestination)){ + stop("Please provide a destination for the Main template") + } + + #check if there is a destination to save the pft/species template + + else if(is.null(PftTemplateDestination)){ + stop("Please provide a destination for the PFT template") + } + + else if(is.null(NameMainFile)){ + stop("Please provide a valid main template") + } + + else if(is.null(NamePftFile)){ + stop("Please provide a valid pft template") + } + + ## TO-DO: Dokumentation + ## TO-DO: Maybe give another return than the function does now + + + + tx <- readLines(NamePftFile) + for(i in 1:nrow(defaultparameters)){ + if(defaultparameters[i,"group"] == "run"){ + tx <- gsub(paste0(defaultparameters[i,"name"]," "), + paste0(defaultparameters[i,"name"]," ", + defaultparameters[i,"rlpjname"]," !"), x = tx) + } + else{ + linenumber_group = grep(paste0("group \"",defaultparameters[i,"group"],"\""), + x = tx, fixed = T) + linenumber_pft = grep(paste0("pft \"",defaultparameters[i,"group"],"\""), + x = tx, fixed = T) + linenumber_st = grep(paste0("st \"",defaultparameters[i,"group"],"\""), + x = tx, fixed = T) + linenumber = c(linenumber_group, linenumber_pft, linenumber_st) + for(j in 1:length(linenumber)){ + endofgroup = grep(")", x = substring(tx[(linenumber[j]+1):length(tx)],1,1))[1] + exchange = grep(paste0(as.character(defaultparameters[i,"name"])), + x = tx[(linenumber[j]+1):(linenumber[j]+1+endofgroup+1)], fixed = T) + tx[(linenumber[j])+exchange] <- gsub(paste0(defaultparameters[i,"name"]," "), + paste0(defaultparameters[i,"name"]," ", + defaultparameters[i,"rlpjname"]," ","!"), + x = tx[(linenumber[j])+exchange]) + } + } + } + + + + writeLines(tx, con = PftTemplateDestination) + + + files <- defaultlist + #change the actual file names with the rlpjnames + + maintemplate <- readLines(NameMainFile) + + for(i in 1:nrow(files)){ + maintemplate <- gsub(paste0("\"",files[i,"name"],"\""), + paste0("\"",files[i,"name"],"\"", + " ","(str ", "\"_",files[i,"name"],"_\")","!"), + x = maintemplate) + } + + #writes the Template to the given MainTemplateDestination + writeLines(maintemplate, con = MainTemplateDestination) + + +} diff --git a/rLPJGUESS/R/GetRunAbleParameters.R b/rLPJGUESS/R/GetRunAbleParameters.R new file mode 100644 index 0000000..cfce3bc --- /dev/null +++ b/rLPJGUESS/R/GetRunAbleParameters.R @@ -0,0 +1,48 @@ +#' @title A function to separate all default values into parameters, design parameters and input files +#' @description This function organizes the parameters extracted from the instruction files in order to run the model +#' @param defaultparameters the default parameters extracted from the instruction files +#' @param PFTs the species one wants to run the model with +#' @keywords rLPJGUESS +#' @author Johannes Oberpriller +#' @return list of runable parameters, design parameters and the files containing the default inputs +#' @export + +GetRunAbleParameters <- function(defaultparameters,PFTs){ + + + LPJparameters_PFT <- matrix(defaultparameters$defaultparameters[,4]) + rownames(LPJparameters_PFT) = defaultparameters$defaultparameters[,3] + LPJrunParameters = matrix(LPJparameters_PFT[-which(substr(rownames(LPJparameters_PFT),1,3)== "run"),]) + rownames(LPJrunParameters) = rownames(LPJparameters_PFT)[-which(substr(rownames(LPJparameters_PFT),1,3)== "run")] + + designLPJ = LPJparameters_PFT[which(substr(rownames(LPJparameters_PFT),1,3)== "run"),] + designLPJ = designLPJ[-which(names(designLPJ) == "run_outputdirectory")] + + lpjvalues = as.matrix(defaultparameters$defaultlist[,c(2,3)]) + + filelist = split(t(lpjvalues),rep(1:nrow(lpjvalues),each = ncol(lpjvalues))) + + names(filelist) = gsub("_", ".",defaultparameters$defaultlist[,1]) + + spp <- as.matrix(LPJrunParameters[grep("_include", rownames(LPJrunParameters )),]) + PFTsRows <- c() + for(i in PFTs) PFTsRows <- c(PFTsRows, grep(i, rownames(spp))) + spp[-PFTsRows, 1] <- 0 + spp[PFTsRows,1] <- 1 + print(spp) + parameterStandard_PFT <- as.matrix(LPJrunParameters[-grep("_include", + rownames(LPJrunParameters)),]) + + parameterStandard_PFT <- rbind(parameterStandard_PFT, spp) + + + tmp <- as.list(parameterStandard_PFT) + namestmp <- rownames(parameterStandard_PFT) + names(tmp) <- namestmp + + parameterStandard_PFT <- tmp + + + return(list(runParameters = parameterStandard_PFT, design = designLPJ, defaultfiles = filelist)) +} + diff --git a/rLPJGUESS/R/InferParameterAndDesignList.R b/rLPJGUESS/R/InferParameterAndDesignList.R new file mode 100644 index 0000000..2dc3d62 --- /dev/null +++ b/rLPJGUESS/R/InferParameterAndDesignList.R @@ -0,0 +1,399 @@ +#' @title A function to infer the parameters and design parameters from the +#' instruction file +#' @description This function returns a list of parameters and design parameters +#' used for the later model runs. It also allows the user to see a summary +#' of all his parameters +#' @param MainInputFile a list with one element named "main" which is the main +#' instruction file given to lpj-guess when running it in the command line +#' @param NameMainFile a character string with the name of the main file. In +#' this file this function will write all input parameters. When no value is +#' provided the default main instruction file while be main.ins +#' @param NamePftFile a character string with the name of the pft file. In this +#' file the function will write all parameters specific for pfts and design. +#' When no value is provided the default main instruction file while be pft.ins +#' @param vectorvaluedparams a vector of character strings providing the vector +#' valued parameters. rlpjguess can in the moment not handel vectorvalued +#' parameters, which can not be computed from the first element. Typically +#' for this is the rootdist, which depends on the different layers in the soil. +#' @return a list with a matrix for the defaultparameter values and a matrix +#' with defaultdesign values +#' @author Johannes Oberpriller +#' @export + + +InferParameterAndDesignList <- function(MainInputFile,NameMainFile = "main.ins", + NamePftFile = "pft.ins", + vectorvaluedparams){ + + ## get the main and pft template from the origin main instruction files + + runlist <- write_main_templates(MainInputFile, NameMainFile = NameMainFile, + NamePftFile = NamePftFile) + + # set up matrix to save all defaultparameters and give them names + + defaultparameters <- data.frame(matrix(nrow = 1, ncol = 5)) + colnames(defaultparameters) = c("type","name", "rlpjname","value", "group") + + for(i in 1:length(runlist$design)){ + parametersnew = vector(length = 4) + parametersnew[1] = "design" + parametersnew[2] = strsplit(runlist$design[[i]]," ")[[1]][1] + parametersnew[3] = paste0("run_",parametersnew[2]) + parametersnew[4] = gsub("[[:space:]]", "", strsplit(runlist$design[[i]]," ")[[1]][2]) + parametersnew[4] = gsub("!","",parametersnew[4]) + parametersnew[5] = "run" + defaultparameters = rbind(defaultparameters, parametersnew) + } + defaultparameters = defaultparameters[-1,] + + runlist_pfts = unlist(runlist$pfts) + # handel all possible ways the instruciton files are + # ordered or the code is saved + for(i in 1:length(runlist$pfts)){ + linenumber = c() + linenumber_group = c() + linenumber_pft =c() + linenumber_st = c() + # get if line is beginning of group,pft or stand + linenumber_group = grep("group ", x = substr(runlist_pfts[i], + start = 1, stop = 6)) + linenumber_pft = grep("pft ", x = substr(runlist_pfts[i], + start = 1, stop = 4)) + linenumber_st = grep("st ", x = substr(runlist_pfts[i], + start = 1, stop = 3)) + linenumber = c(linenumber_group, linenumber_pft, linenumber_st) + linename = strsplit(runlist_pfts[i]," ")[[1]][2] + if(length(linenumber)!=0){ + endoufgroup = c() + # get the end of the group + endofgroup = grep(")", x = substr(runlist_pfts[(i+1):length(runlist_pfts)],1,2), fixed = T)[1] + # put all lines between beginning and end into the parametermatrix + for(j in 1:(endofgroup-1)){ + parametersnew = vector(length = 4) + parametersnew[1] = "pft" + parametersnew[2] = gsub("[[:space:]]","",strsplit(runlist_pfts[i+j]," ")[[1]][1]) + parametersnew[3] = paste0(gsub("\"","",linename),"_", parametersnew[2]) + parameterposition = 2 + parametersnew[4] = "" + while(grepl("^\\s*$", parametersnew[4])){ + parametersnew[4] = gsub("[[:space:]]"," ", + strsplit(runlist_pfts[i+j]," ", fixed = T)[[1]][parameterposition]) + parameterposition = parameterposition + 1 + if(grepl("^!", parametersnew[4])){ + parametersnew[4] = "" + break + } + } + parametersnew[4] = gsub("^!","",parametersnew[4]) + parametersnew[5] = gsub("\"","",linename) + defaultparameters = rbind(defaultparameters, parametersnew) + } + } + } + #print(defaultparameters) + # handle possible errors + defaultparameters = defaultparameters[which(defaultparameters[,2] != ""),] + defaultparameters = defaultparameters[which(defaultparameters[,2] != "!"),] + defaultparameters = defaultparameters[which(defaultparameters[,4] != "NA"),] + # get ride of spaces and comments for possible multiplications later + for(i in 1:nrow(defaultparameters)){ + startofcomments = gregexpr("!",defaultparameters[i,4])[[1]][1] + if(startofcomments != "-1"){ + defaultparameters[i,4] = sub("[[:space:]]+","",substr(defaultparameters[i,4], start = 1, stop = startofcomments -1)) + } + } + # exclude all vectorvalued parameters from the list + for(i in 1:length(vectorvaluedparams)){ + defaultparameters = defaultparameters[which(defaultparameters[,2] != vectorvaluedparams[i]),] + } + + + + ### Infer the Design list + + # unpack all variables and climate input files + + runlist_variables = unlist(runlist$variables,use.names = F) + runlist_inputfiles = unlist(runlist$files, use.names = F) + + # define matrix to store them + + defaultlist = data.frame(matrix(nrow = 1, ncol = 3)) + colnames(defaultlist) = c("name","rlpjnames","value") + + # handle all possible ways stored in file + for(i in 1:length(runlist_variables)){ + new_variable = vector(length = 3) + new_variable[1] = gsub("\"","",strsplit(runlist_variables[i]," ")[[1]][2]) + beginn = gregexpr("\\(",runlist_variables[i]) + end = gregexpr("\\)", runlist_variables[i]) + parameterpart = gsub("\"", "",strsplit(substr(runlist_variables[i], start = beginn[[1]], stop = end[[1]]-1), " ")[[1]][2]) + new_variable[2] = paste0("_",new_variable[1],"_") + new_variable[3] = parameterpart + defaultlist = rbind(defaultlist, new_variable) + } + + defaultlist = defaultlist[-1,] + # handle all possible ways stored in file + for(i in 1:length(runlist_inputfiles)){ + new_variable = vector(length = 3) + new_variable[1] = gsub("\"","",strsplit(runlist_inputfiles[i]," ")[[1]][2]) + beginn = gregexpr("\\(",runlist_inputfiles[i]) + end = gregexpr("\\)", runlist_inputfiles[i]) + parameterpart = gsub("\"", "",strsplit(substr(runlist_inputfiles[i], start = beginn[[1]], stop = end[[1]]-1), " ")[[1]][2]) + new_variable[2] = paste0("_",new_variable[1],"_") + new_variable[3] = parameterpart + defaultlist = rbind(defaultlist, new_variable) + } + + # get ride of empty lines + #print(defaultparameters) + defaultlist = defaultlist[which(defaultlist[,1] != ""),] + + return(list(defaultparameters = defaultparameters, defaultlist = defaultlist)) +} + + +# @description This function takes the main input input files and desired file +# names and returns a list with files, variables, parameters and design parameters +# @param MainInputFile is the main input file see InferParameterAndDesignList +# @param NameMainFile is the main file after reordering see InferParameterAndDesignList +# @param NamePftFile is the pft file after reordering see InferParameterAndDesignList +# @return returns a list with the entries parameters,files, variables and design +# which are all themselves matrices +# @author Johannes Oberpriller + +write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ + + # Infer the correct order of instruction files + + insfilelist = infer_instruction_order(MainInputFile = MainInputFile) + + # Create the temporary files to get all lines + + file.create("./main_tryout.ins") + outcon <- file("./main_tryout.ins", "w") + for(i in 1:length(insfilelist)){ + lines <- readLines(insfilelist[[i]]) + writeLines(lines,outcon) + } + close(outcon) + + # Create the Main and Ins file + + file.create(paste0("./",NameMainFile)) + file.create(paste0("./",NamePftFile)) + lines <- readLines("./main_tryout.ins") + + # delete all import files + + importlines = grepl("import",substr(lines, start = 1, stop = 7)) + lines = lines[-which(importlines)] + + # Grep all "param" lines and the ndep_timeseries, delete duplicates + # and write them to the parameterfiles + + paramlines = grepl("param", substr(lines, start = 1, stop = 7)) + parameters = lines[which(paramlines)] + + parameternames = sapply(strsplit(parameters," "), function (x) x[2]) + duplicates_parameters = which(duplicated(parameternames) ==T) + + while(length(duplicates_parameters) != 0 ){ + first_apperance = c() + for(i in duplicates_parameters){ + first_apperance <- c(first_apperance, which(parameternames == parameternames[i])[1]) + } + parameters = parameters[-first_apperance] + parameternames = sapply(strsplit(parameters," "), function (x) x[2]) + duplicates_parameters = which(duplicated(parameternames) ==T) + } + + ndep_time = grepl("!ndep_timeseries", substr(lines, start = 1, stop = 16)) + if(length(which(ndep_time ==T)) ==0){ + ndep_time = grepl("ndep_timeseries", substr(lines, start =1, stop =15)) + } + parameters = c(parameters, lines[which(ndep_time)]) + parameternames = sapply(strsplit(parameters," "), function (x) x[2]) + parameterfiles = grepl("file", parameternames) + parametervariables = parameters[-which(parameterfiles)] + parameterfiles = parameters[parameterfiles] + + # Write the main Files and delete all lines which are written there + + FileCon = file(NameMainFile) + mainimport = paste("import ",paste0("\"", "path_to_globalTemplate", "\"", "\n"), sep =" ") + writeLines(c(mainimport,parameters), FileCon) + close(FileCon) + lines = lines[-c(which(paramlines),which(ndep_time))] + headingfileline <- paste("\n", "!Which files to include in the output","\n") + + ## Get the Title, the Outputfiles and the outcommented lines + + titleline <- grepl("title",substr(lines, start = 1, stop = 10)) + filelines <- grepl("file", substr(lines, start = 1, stop = 5)) + outcommented <- grepl("!", substr(lines, start = 1, stop = 1)) + + ## get the lines where pfts and/or species parameters are specified + + pftlines = c() + for(i in 1:length(lines)){ + linenumber = c() + linenumber_group = c() + linenumber_pft =c() + linenumber_st = c() + linenumber_group = grep("group ", x = substr(lines[i], start = 1, stop = 6)) + linenumber_pft = grep("pft ", x = substr(lines[i], start = 1, stop = 4)) + linenumber_st = grep("st ", x = substr(lines[i], start = 1, stop = 3)) + linenumber = c(linenumber_group, linenumber_pft, linenumber_st) + if(length(linenumber)!=0){ + endoufgroup = c() + endofgroup = grep(")", x = substr(lines[(i+1):length(lines)],1,2), fixed = T)[1] + pftlines = c(pftlines, lines[i:(i+endofgroup)],"\n") + } + } + + titel = lines[titleline] + files = lines[filelines] + files = gsub("^file","!file", files) + + # delete lines which are titlelines, filelines ( which we alter to be + # all outcommented), and outcommented lines and pftlines + + lines = lines[-c(which(titleline),which(filelines), + which(outcommented),which(lines %in% pftlines))] + + # grep empty and tabed lines and delete them + + tablines = grepl("\t", x = substr(lines, start =1, stop = 2)) + emptylines = grepl(" ", substr(lines, start = 1, stop =2)) + lines = lines[-c(which(tablines), which(emptylines))] + + # get the unique designnames from the parameters + + designnames = sapply(strsplit(lines," "), function (x) x[1]) + duplicates_design = which(duplicated(designnames) ==T) + + first_apperance = c() + for(i in duplicates_design){ + first_apperance <- c(first_apperance, which(designnames == designnames[i])[1]) + } + design = lines[-first_apperance] + + # write the PFT or Species file which includes the title, the files, + # the designvariables and all pft/speciesspecific lines into the + # main Pft file + + Pftcon <- file(NamePftFile) + writeLines(c(titel,headingfileline,files,"\n",design, "\n",pftlines),Pftcon) + close(Pftcon) + + # remove the main file written in the beginning because no more need + + file.remove("./main_tryout.ins") + + return(list(files = parameterfiles, variables = parametervariables, + pfts = pftlines, design = design )) + +} + +# @description This function infers the correct order of the instruction files +# and returns it +# @param MainInputFile is the main input file see InferParameterAndDesignList +# @return returns a vector with the names of correct ordered instruction files +# @author Johannes Oberpriller + +infer_instruction_order <- function(MainInputFile){ + + # get the main directory where the instruction file is stored + + dashes = gregexpr("/",MainInputFile[["main"]]) + + lastdash = dashes[[1]][length(dashes[[1]])] + + maindirectory = substr(MainInputFile[["main"]],start = 1, + stop = lastdash) + + # get all imports from the file + + imports = get_imports(MainInputFile[["main"]],maindirectory = maindirectory) + + # order the imports like they are in the main instruction file + + ordered_imports = imports + i = 1 + while(i <= length(ordered_imports)){ + new_imports = get_imports(ordered_imports[[i]], maindirectory = maindirectory) + if(any((new_imports %in% ordered_imports))){ + i = i+1 + } + else{ + if(length(new_imports) == 0){ + i = i+1 + } + else{ + number_new = length(new_imports) + number_old = length(ordered_imports) + if((i+1)>number_old){ + ordered_imports[number_new + number_old] = ordered_imports[number_old] + ordered_imports[i:(i+number_new-1)] = new_imports + } + else{ + ordered_imports[(i+number_new):(number_new+number_old)] = ordered_imports[i:number_old] + ordered_imports[(i):(i+number_new-1)] = new_imports + } + } + } + } + # return the correct order + return(c(MainInputFile, ordered_imports)) +} + +# @description this function takes the main file and the main directroy +# and returns the imported files +# @param a file with all instruction files below each other +# @param maindirectory the main directory for the analysis +# @return returns the ordered imports +# @author Johannes Oberpriller + + +get_imports <- function(file,maindirectory){ + + # read main file + + linesfile = readLines(file) + imports = list() + import_number = c() + + # loop over files and get instruction order + for(i in 1:length(linesfile)){ + is_in = grepl("import", linesfile[i]) + out_commented = grepl("!import",substr(linesfile[i], start = 1, stop = 8)) + if(is_in == T){ + if(out_commented == T){ + next + } + else{ + import_number =c(import_number, i) + } + } + else + { + next + } + } + # check that the imports have all the same apperance in code + ins_order = 1 + for(i in import_number){ + imports_file = substr(x=linesfile[i], + start = gregexpr("\"",linesfile[i])[[1]][1]+1, + stop = gregexpr("\"",linesfile[i])[[1]][2]-1) + if(substr(imports_file,1,2) =="./"){ + imports_file = substr(imports_file,3,stop =100) + } + imports[ins_order] = paste0(maindirectory,imports_file) + ins_order = ins_order +1 + } + return(imports) +} diff --git a/rLPJGUESS/R/callLPJ.R b/rLPJGUESS/R/callLPJ.R index 78ad299..4c6f20e 100644 --- a/rLPJGUESS/R/callLPJ.R +++ b/rLPJGUESS/R/callLPJ.R @@ -28,9 +28,9 @@ callLPJ <- function(mainDir, runDir , template2, mode){ if (is.null(template2) || !file.exists(file.path(runDir, template2))){ stop ("Please provide a valid template name") } - if (is.null(mode) || mode != "cf" & mode != "cru"){ - stop("Please provide a valid cluster type: cf or cru") - } + # if (is.null(mode) || mode != "cf" & mode != "cru"){ + # stop("Please provide a valid cluster type: cf or cru") + # } #----------------------------------------------------------------------------# # CALL MODELL diff --git a/rLPJGUESS/R/checkDesign.R b/rLPJGUESS/R/checkDesign.R index 3b13ef9..01bf5fe 100644 --- a/rLPJGUESS/R/checkDesign.R +++ b/rLPJGUESS/R/checkDesign.R @@ -1,7 +1,7 @@ # @title A check design function # @description This function checks the provided design against the -# default values, and returns a complete design. If wrong desing is passed, +# default values, and returns a complete design. If wrong design is passed, # the fucntion will raise an error. # @param scale a character string indicating whether the model runs global or for europe. # @param design a named list or matrix holding the design diff --git a/rLPJGUESS/R/checkParameters.R b/rLPJGUESS/R/checkParameters.R index 48afaa8..80971d2 100644 --- a/rLPJGUESS/R/checkParameters.R +++ b/rLPJGUESS/R/checkParameters.R @@ -17,9 +17,9 @@ # } checkParameters <- function(scale, parameterList = NULL, type = "serial"){ # include check - if ( scale != "global" & scale != "europe"){ - stop("Invalid scale: neither global nor europe") - } + # if ( scale != "global" & scale != "europe"){ + # stop("Invalid scale: neither global nor europe") + # } if (type != "parallel" & type != "serial" ){ stop("Invalid check type") } @@ -36,6 +36,7 @@ checkParameters <- function(scale, parameterList = NULL, type = "serial"){ stop("Please provide a valid parameterList") } # Throw an error if wrong parameter names + dummyCheck <- !names(parameterList) %in% names(default) if(any(dummyCheck)){ warning(paste("Wrong parameterList in ", paste(names(parameterList)[dummyCheck], collapse = ", " ))) @@ -64,15 +65,16 @@ checkParameters <- function(scale, parameterList = NULL, type = "serial"){ checkParameters.names <- function(scale, parameterNames){ # include check - if ( scale != "global" & scale != "europe"){ - stop("Invalid scale: neither global nor europe") - } + # if ( scale != "global" & scale != "europe"){ + # stop("Invalid scale: neither global nor europe") + # } # call the default template + default <- getParameterList(scale, list = T) - if(is.null(parameterNames)){ - parameterList <- default - } + # if(is.null(parameterNames)){ + # parameterList <- default + # } # Throw an error if wrong parameter names dummyCheck <- !parameterNames %in% names(default) if(any(dummyCheck)){ diff --git a/rLPJGUESS/R/createRunParameters.R b/rLPJGUESS/R/createRunParameters.R index cacc039..1e3f1d7 100644 --- a/rLPJGUESS/R/createRunParameters.R +++ b/rLPJGUESS/R/createRunParameters.R @@ -12,7 +12,8 @@ createRunParameters <- function(x, singleRun, parameterList){ # Check the parameters: length, type and names - parameterList <- try(checkParameters.matrix(singleRun$scale, parameterList), FALSE) + #parameterList <- try(checkParameters.matrix(singleRun$scale, parameterList), FALSE) + parameterList <- checkParameters.matrix(singleRun$scale, parameterList) if ('try-error' %in% class(parameterList)){ stop("Invalid parameterList provided") } if (class(parameterList[[1]])== "list"){ @@ -32,18 +33,20 @@ createRunParameters <- function(x, singleRun, parameterList){ # Based on parameter Names write the general template # and if no parameter is present raise and error? - parameterCommon <- checkParameters.names(scale= singleRun$scale, parameterNames) + #parameterCommon <- checkParameters.names(scale = singleRun$scale, parameterNames) + #print(parameterList) + # parameterCommon <- parameterList[[1]] + # + # if(length(parameterCommon) > 0){ + # #parameterCommon <- checkParameters.rootDist(parameterCommon) + # #write common template + # parameterCommonNames <- names(parameterCommon) + # for(i in 1:length(parameterCommon)){ + # singleRun$template1Mem <- sub(paste0(" ",parameterCommonNames[i]," "), paste0(" ",parameterCommon[[i]]," "), singleRun$template1Mem) + # } + # } - if(length(parameterCommon) > 0){ - parameterCommon <- checkParameters.rootDist(parameterCommon) - # write common template - parameterCommonNames <- names(parameterCommon) - for(i in 1:length(parameterCommon)) { - singleRun$template1Mem <- sub(parameterCommonNames[i], parameterCommon[[i]], singleRun$template1Mem) - } - - } # Check the grids gridListCell <- readLines(file.path(singleRun$mainDir,singleRun$gridList)) diff --git a/rLPJGUESS/R/createSingleObject.R b/rLPJGUESS/R/createSingleObject.R index 5b80e93..57ddf3a 100644 --- a/rLPJGUESS/R/createSingleObject.R +++ b/rLPJGUESS/R/createSingleObject.R @@ -6,7 +6,7 @@ # Default value is all outputs. # @param settings additional parameters # @seealso \code{\link{runLPJ}} -# @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig +# @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig, Johannes Oberpriller # @return TODO createSingleObject <- function(mainDir, typeList, settings){ @@ -14,24 +14,18 @@ createSingleObject <- function(mainDir, typeList, settings){ file.co2 = NULL, file.cru = NULL, file.cru.misc = NULL, file.ndep= NULL, file.temp = NULL, file.prec = NULL, file.insol = NULL, file.wetdays = NULL, file.minTemp = NULL, - file.maxTemp = NULL, variable.temp = NULL, + file.soildata = NULL, + file.maxTemp = NULL, variable.temp = NULL, variable.ndep = NULL, variable.prec = NULL, variable.insol = NULL, variable.wetdays = NULL, variable.minTemp = NULL, variable.maxTemp = NULL, template1 = NULL, template2=NULL, plot.data = FALSE, save.plots = FALSE, processing = FALSE, delete = TRUE, save= TRUE, runID = "", parallel = "auto", - checkParameters = "serial", design = NULL) + checkParameters = "serial", design = NULL, defaultlist = NULL) #, fun = NULL) # This would be to allow havin own functions in parallel. settings <- c(settings[names(settings) %in% names(defaultSettings)], defaultSettings[ !names(defaultSettings) %in% names(settings)]) - # mode - if (is.null(settings[["mode"]]) || settings[["mode"]] != "cf" & settings[["mode"]] != "cru"){ - stop("Please provide a valid cluster type: cf or cru") - } - if ( is.null(settings[["scale"]]) || settings[["scale"]] != "global" & settings[["scale"]] != "europe"){ # this is relevant if getting template - stop("Please provide a valid scale: global or europe") - } if (is.null(typeList) || !class(typeList) == "character"){ settings$typeList <- typelist.default message("Output typeList has not been provided") @@ -87,50 +81,54 @@ createSingleObject <- function(mainDir, typeList, settings){ # to replace in the template # Go throught the files and check whether they provided, if so add them to # default list, otherwise stop the function + + # Johannes: Depending on the version of LPJ different parameters are + # required to run the model. + # I think, every LPJ user has working instruction files, + # which allow to infer all required files and it is enough to provide + # the functionality to substitute them in the main file if they want to change + # this from the wrapper + files <- settings[grepl("file", names(settings))] - files.default <- files.parameters[[settings[["mode"]]]] - files.names <- names(files.default) + variables <- settings[grepl("variable", names(settings))] + filesandvariables = c(files,variables) + + files.names <- names(settings$defaultlist) + for (i in 1:length(files.names)){ - if (is.null(files[[files.names[i]]])){ - warning(paste("The", files.names[i], "has not been provided", sep = " ")) - }else if(!file.exists(files[[files.names[i]]])){ + if (is.null(filesandvariables[[files.names[i]]])){ + next + }else if(!file.exists(filesandvariables[[files.names[i]]])){ warning(paste("The", files.names[i], "does not exist", sep = " ")) }else{ - files.default[[files.names[i]]][2] <- files[[files.names[i]]] - } - } - #files.default <- files.default[keep] + #print(i) + #print(filesandvariables[[files.names[i]]]) + settings$defaultlist[[files.names[i]]][2] <- filesandvariables[[files.names[i]]] - variables <- settings[grepl("variable", names(settings))] - variables.default <- variables.parameters[[settings[["mode"]]]] - variables.names <- names(variables.default) - for (i in 1:length(variables.names)){ - if (is.null(variables[[variables.names[i]]])){ - warning(paste("The", variables.names[i], "has not been provided", sep = " ")) - }else{ - variables.default[[variables.names[i]]][2] <- variables[[variables.names[i]]] } } - #variables.default <- variables.default[keep] + singleObject <- settings[!grepl("file", names(settings))] singleObject <- singleObject[!grepl("variable", names(singleObject))] - singleObject$filesNames <- files.default - singleObject$variablesNames <- variables.default + singleObject$filesNames <- settings$defaultlist singleObject$mainDir <- mainDir singleObject$runInfoDir <- file.path(singleObject$mainDir, paste("runInfo", format(Sys.time(), "%Y_%m_%d_%H%M%S"), sep = "_")) - # Read template one and replace desing + # Read template one and replace design singleObject$template1Mem <- readLines(file.path(singleObject$mainDir, singleObject$template1)) # Check the design - settings$design <- checkDesign(settings$scale , settings$design) + #settings$design <- checkDesign(settings$scale , settings$design) designNames <- names(settings$design) - for(i in 1:length(settings$design)) { - singleObject$template1Mem <- sub(designNames[i], settings$design[[i]], singleObject$template1Mem) + + for(i in 1:length(settings$design)){ + singleObject$template1Mem <- gsub(paste0(" ",designNames[[i]]," "), + paste0(" ",settings$design[[i]]," "), + singleObject$template1Mem) } if(settings$design[["run_ifcalcsla"]]==as.character(0)){ singleObject$template1Mem <- sub("!sla", "sla", singleObject$template1Mem) diff --git a/rLPJGUESS/R/getData.R b/rLPJGUESS/R/getData.R index 7a909c4..972453d 100644 --- a/rLPJGUESS/R/getData.R +++ b/rLPJGUESS/R/getData.R @@ -23,18 +23,18 @@ #' @export #' @example /inst/examples/getDataHelp.R getLPJData <- function(x, typeList = NULL, runInfo=NULL, processing = FALSE){ - #, fun = NULL){ + #, fun = NULL){ # other options which could be included: - #lon.extent=c(-180, 180), lat.extent=c(-90, 90), - #area.weighted=FALSE, year.offset=0 ) { - # @param lon.extent a numeric vector containing the min and max values used - # for west-east extent(default -180 to 180) - # @param lat.extent a numeric vector containing the max and min values used - # for north-south extent (default 90 to -90). - # @param area.weighted a boolean indicating whether the gridcells should be - # weighted by their size (regular grids only, default FALSE). - # @param year.offset a integer indicating the value to be added to the 'Year' - # column in the LPJ-GUESS output. + #lon.extent=c(-180, 180), lat.extent=c(-90, 90), + #area.weighted=FALSE, year.offset=0 ) { + # @param lon.extent a numeric vector containing the min and max values used + # for west-east extent(default -180 to 180) + # @param lat.extent a numeric vector containing the max and min values used + # for north-south extent (default 90 to -90). + # @param area.weighted a boolean indicating whether the gridcells should be + # weighted by their size (regular grids only, default FALSE). + # @param year.offset a integer indicating the value to be added to the 'Year' + # column in the LPJ-GUESS output. # checking provided parameters @@ -63,91 +63,78 @@ getLPJData <- function(x, typeList = NULL, runInfo=NULL, processing = FALSE){ keep <- rep(FALSE, length(typeList)) for (i in 1:length(typeList)){ if (file.exists(file.path(x, paste(typeList[[i]], ".out", sep="")))){ - if ( file.info( file.path(x, paste(typeList[[i]], ".out", sep="")) )[['size']] == 0){ + if ( file.info( file.path(x, paste(typeList[[i]], ".out", sep="")) )[['size']] == 0){ warning( paste("The ", typeList[[i]], ".out is empty!", sep = "") ) }else{ keep[i] <- TRUE } }else{ warning( paste("There is no ", typeList[[i]], ".out", sep = "") ) - } } - if (any(keep)){ + } + if(any(keep)){ typeList.valid <- typeList[keep] # Creating a list to hold data listData <- vector(mode="list", length=length(typeList.valid)) names(listData) <- typeList.valid }else{ - stop("There are not model outputs. Please check the guess.log files") + stop("There are no model outputs. Please check the guess.log files") } - # storing run info + # storing run info #----------------------------------------------------------------------------# # OBTAIN OUTPUTS: #----------------------------------------------------------------------------# - if (processing == FALSE){ - # Adding Data to listData - # looping over data types, reading files, no processing data and adding it to the Data Class - # list append data, probably will have to use the name() function to give it the right name - # in the end, make the list of class LPJData - # Add data to LPJout Data class - for (j in 1:length(typeList.valid)) { - data <- try(read.table(file.path(x, paste( typeList.valid[j], ".out", sep="")),header=T), TRUE) - if ('try-error' %in% class(data)){ - data <- try(data <-readTableHeaderLPJ(file.path(x, paste( typeList.valid[j], ".out", sep=""))), TRUE) - if ('try-error' %in% class(data)){ - stop("Model outputs are not readable") - } + if (processing == FALSE){ + # Adding Data to listData + # looping over data types, reading files, no processing data and adding it to the Data Class + # list append data, probably will have to use the name() function to give it the right name + # in the end, make the list of class LPJData + # Add data to LPJout Data class + for(j in 1:length(typeList.valid)){ + data <- try(read.table(file.path(x, paste(typeList.valid[[j]], ".out", sep="")),header=T), TRUE) + if('try-error' %in% class(data)){ + data <- try(data <-readTableHeaderLPJ(file.path(x, paste(typeList.valid[[j]], ".out", sep=""))), TRUE) + if('try-error' %in% class(data)){ + stop("Model outputs are not readable") } - listData[[typeList.valid[[j]]]] <- data } - }else{ - for (j in 1:length(typeList.valid)) { - # Chekc how many grids in output - data <- try(read.table(file.path(x, paste( typeList.valid[j], ".out", sep="")),header=T), TRUE) + listData[[typeList.valid[[j]]]] <- data + } + }else{ + for(j in 1:length(typeList.valid)) { + data <- try(read.table(file.path(x, paste( typeList.valid[j], ".out", sep="")),header=T), TRUE) + if ('try-error' %in% class(data)){ + data <- try(data <-readTableHeaderLPJ(file.path(x, paste( typeList.valid[j], ".out", sep=""))), TRUE) if ('try-error' %in% class(data)){ - data <- try(data <-readTableHeaderLPJ(file.path(x, paste( typeList.valid[j], ".out", sep=""))), TRUE) - if ('try-error' %in% class(data)){ - stop("Model outputs are not readable") - } + stop("Model outputs are not readable") } - coordinates <- unique(paste(data$Lat, data$Lon, sep = "_")) + } + listData[[typeList.valid[[j]]]] <- data + } + coordinates <- unique(paste(data$Lat, data$Lon, sep = "_")) + if(length(coordinates) > 1 ){ + coordinates <- lapply(coordinates, function(x){as.numeric(unlist(strsplit(x, "_")))}) + listData <- rep(list(listData), length(coordinates)) + names(listData) <- paste0("gridcell",coordinates) # will need better name + #Adding Data to listData + #looping over data types, reading files, processing data and adding it to the Data Class + #list append data, probably will have to use the name() function to give it the right name + #in the end, make the list of class LPJData - if(length(coordinates) > 1 ){ - # Now we stop but eventually we want to do something out of this - # Below there is code that would do it - # This would also affect plotLPJData - stop("Processing is not supported for more than one grid") - # if only one grid, simplify the list - # if (length(listData) == 1){ - # listData <- listData[[1]] - # Adding Data to listData - # looping over data types, reading files, processing data and adding it to the Data Class - # list append data, probably will have to use the name() function to give it the right name - # in the end, make the list of class LPJData - # Add data to LPJout Data class - #coordinates <- lapply(coordinates, function(x){as.numeric(unlist(strsplit(x, "_")))}) - #listData <- rep(list(listData), length(coordinates)) - #names(listData) <- paste("grid", coordinates, sep = "_") # will need better named - # Adding Data to listData - # looping over data types, reading files, processing data and adding it to the Data Class - # list append data, probably will have to use the name() function to give it the right name - # in the end, make the list of class LPJData - #keep <- rep(TRUE, length(coord)) - #sub_data <- coord[coord>=min(lon.extent) & Lon<=max(lon.extent) & Lat <=max(lat.extent) & Lat>=min(lat.extent)) - #for (k in 1:length(listData)){ - # data <- data[data$Lat==coordinates[[k]][1] & data$Lon==coordinates[[k]][2],] + for(k in 1:length(listData)){ - # for (j in 1:length(typeList.valid)) { - # # reading output - # data <- read.table(file.path(x, paste(typeList.valid[[j]], ".out", sep="")),header=T) - # data.ts <- convertTS(data) - # listData[[k]][[typeList.valid[[j]]]] <- data.ts - # } - }else{ + for(j in 1:length(typeList.valid)){ + # reading output + data <- read.table(file.path(x, paste(typeList.valid[[j]], ".out", sep="")),header=T) + data <- data[data$Lat==coordinates[[k]][1] & data$Lon==coordinates[[k]][2],] data.ts <- convertTS(data) - listData[[typeList.valid[[j]]]] <- data.ts + listData[[k]][[typeList.valid[[j]]]] <- data.ts + } } + }else{ + data.ts <- convertTS(data) + listData[[typeList.valid[[j]]]] <- data.ts } } # add it to the data class diff --git a/rLPJGUESS/R/getDesign.R b/rLPJGUESS/R/getDesign.R index b495b68..a284576 100644 --- a/rLPJGUESS/R/getDesign.R +++ b/rLPJGUESS/R/getDesign.R @@ -8,10 +8,10 @@ #' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara #' @example /inst/examples/getDesignHelp.R getDesign <- function(scale, list = F){ - if ( is.null(scale) || scale != "global" & scale != "europe"){ - stop("Please provide a valid scale: global or europe") - } - tmp <- design.default[[scale]] + # if ( is.null(scale) || scale != "global" & scale != "europe"){ + # stop("Please provide a valid scale: global or europe") + # } + tmp <- design.default$scale if (list){ tmp.names <- rownames(tmp) tmp <- as.vector(tmp, mode = "list") diff --git a/rLPJGUESS/R/getParameterList.R b/rLPJGUESS/R/getParameterList.R index 820b2b2..b91e0e3 100644 --- a/rLPJGUESS/R/getParameterList.R +++ b/rLPJGUESS/R/getParameterList.R @@ -11,11 +11,12 @@ #' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara #' @example /inst/examples/getParameterListHelp.R getParameterList <- function(scale, list = TRUE){ - if ( is.null(scale) || scale != "global" & scale != "europe"){ + if (is.null(scale)){ stop("Please provide a valid scale: global or europe") } - tmp <- parameterList.default[[scale]] + #tmp <- parameterList.default[1]$normal + tmp <- parameterList.default if (list){ tmp.names <- rownames(tmp) diff --git a/rLPJGUESS/R/plotLPJData.R b/rLPJGUESS/R/plotLPJData.R index a907857..45d2d5e 100644 --- a/rLPJGUESS/R/plotLPJData.R +++ b/rLPJGUESS/R/plotLPJData.R @@ -17,89 +17,102 @@ #' \url{https://cran.r-project.org/web/packages/zoo/zoo.pdf} #' @export #' @keywords rLPJGUESS -#' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig +#' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig, Johannes Oberpriller #' @example /inst/examples/plotLPJDataHelp.R -plotLPJData <- function(x, typeList = NULL, outDir= NULL, save.plots = FALSE, prefix = ""){ +plotLPJData <- function(x, typeList = NULL, gridlist = NULL, outDir= NULL, + save.plots = FALSE, prefix = "", plot_avg = F){ # checking input parameters - if (is.null(x)){ + if(is.null(x)){ stop("No data has been provided") } - if (!class(x)=="LPJData"){ + if(!class(x)=="LPJData"){ stop("Invalid data has been provided") } if (save.plots){ - if ( is.null(outDir) || !file.exists(outDir)){ + if( is.null(outDir) || !file.exists(outDir)){ stop("No outDir has been provided") } } - if (!requireNamespace("zoo", quietly = TRUE)){ + if(!requireNamespace("zoo", quietly = TRUE)){ stop("Can't load required library 'zoo'") } data <- x@dataTypes # Plot from - typeList.available <- sort(names(data)) - - if (is.null(typeList) || !class(typeList) == "character"){ + gridlist.available <- sort(names(data)) + typeList.available <- names(data[[1]][which(names(data[[1]]) %in% + gridlist.available == F)]) + if(is.null(typeList) || !class(typeList) == "character"){ message("No typeList has been provided. Plotting all data") typeList.valid <- typeList.available }else{ keep <- rep(FALSE, length(typeList)) - for (i in 1:length(typeList)){ - if (typeList[i] %in% typeList.available){ + for(i in 1:length(typeList)){ + if(typeList[i] %in% typeList.available){ keep[i] <- TRUE } } - if (any(keep)){ + if(any(keep)){ typeList.valid <-typeList[keep] }else{ stop("None of the requested output types exists") } } - ## if plotting true, start the fucntion - # Check if data is a df - # Checking existentce of data types # in theory do not need if your plotting from class object - for (i in 1:length(typeList.valid)){ - df <- data[[typeList.valid[[i]] ]] + + if(is.null(gridlist)){ + message("No gridlist has been provided. Plotting all gridcells") + gridlist.valid <- gridlist.available + }else{ + gridlist.check = apply(X = gridlist,MARGIN = 1, FUN = function (x) { + y = paste0("gridcellc(",x[1], ", ", x[2],")") + return(y) + }) + keep <- rep(FALSE, length(gridlist.available)) + for(i in 1:length(gridlist.check)){ + if(gridlist.check[i] %in% gridlist.available){ + keep[i] <- TRUE + } + } + if(any(keep)){ + gridlist.valid <-gridlist.check[keep] + }else{ + stop("None of the requested gridcells exists") + } + } + average = matrix(data = 0,ncol = length(gridlist.valid), + nrow = nrow(data[[gridlist.valid[[1]]]][[typeList.valid[[1]]]])) + colnames(average) = typeList.valid + for (i in 1:length(gridlist.valid)){ + df <- data[[gridlist.valid[[i]]]] + coordinates = get_gridlist_values(gridlist.valid[[i]]) if(zoo::is.zoo(df) == FALSE){ - # check how many coordinates - coordinates <- unique(paste(df$Lat, df$Lon, sep = "_")) - #if(length(coordinates) > 1 ){ - # if only one grid, simplify the list - # if (length(listData) == 1){ - # listData <- listData[[1]] - # Adding Data to listData - # looping over data types, reading files, processing data and adding it to the Data Class - # list append data, probably will have to use the name() function to give it the right name - # in the end, make the list of class LPJData - # Add data to LPJout Data class - coordinates <- lapply(coordinates, function(x){as.numeric(unlist(strsplit(x, "_")))}) - #keep <- rep(TRUE, length(coord)) - #sub_data <- coord[coord>=min(lon.extent) & Lon<=max(lon.extent) & Lat <=max(lat.extent) & Lat>=min(lat.extent)) - for (k in 1:length(coordinates)){ - values <- df[df$Lat==coordinates[[k]][1] & df$Lon==coordinates[[k]][2],] - values <- convertTS(values) - if (save.plots){ - pdf(file.path(outDir, paste(prefix, typeList.valid[[i]], ".pdf", sep="")))#width=1000,height=750 - if(length(colnames(values))==1){ - plot(values, main =paste("Grid", coordinates[[k]][1], coordinates[[k]][2], - "Variable:", typeList.valid[[i]]),xlab="Years") - }else{ - plot(values, main =paste("Grid", coordinates[[k]][1], coordinates[[k]][2], - "Variable:", typeList.valid[[i]]),xlab="Years") - } - dev.off() + + for (k in 1:length(typeList.valid)){ + values <- df[[typeList.valid[[k]]]] + average[,k] = average[,k] + values[,ncol(values)] + if(save.plots){ + pdf(file.path(outDir, paste(prefix, gridlist.valid[[i]],typeList.valid[[k]], ".pdf", sep="")))#width=1000,height=750 + if(length(colnames(values))==1){ + plot(values, main =paste("Gridcell Lon:", coordinates[1],"Lat", coordinates[2], + "Variable:", typeList.valid[[k]]),xlab="Years") }else{ - plot(values, main =paste("Grid", coordinates[[k]][1], coordinates[[k]][2], - "Variable:", typeList.valid[[i]]),xlab="Years") + plot(values, main =paste("Gridcell Lon:", coordinates[1],"Lat", coordinates[2], + "Variable:", typeList.valid[[k]]),xlab="Years") } + dev.off() + }else{ + plot(values, main =paste("Gridcell Lon:", coordinates[1],"Lat", coordinates[2], + "Variable:", typeList.valid[[k]]),xlab="Years") } - }else{ + + } + } + else{ values <- df - # something like is zoo + # something like is zoo if (save.plots){ pdf(file.path(outDir, paste(prefix, typeList.valid[[i]], ".pdf", sep="")))#width=1000,height=750 plot(values, main =paste("Variable:", typeList.valid[[i]]),xlab="Years") @@ -108,11 +121,27 @@ plotLPJData <- function(x, typeList = NULL, outDir= NULL, save.plots = FALSE, pr plot(values, main =paste("Variable:", typeList.valid[[i]]),xlab="Years") } } + } + if(plot_avg == T){ + for(k in 1:length(typeList.valid)){ + plot(y = average[,typeList.valid[[k]]]/length(gridlist.valid), x = as.numeric(rownames(values)), + main = paste("Average of", "Variable:",typeList.valid[[k]]), + xlab = "Years",ylab = "total", type = "l") + } + + } } +get_gridlist_values <- function(gridcell){ + coordinates = strsplit(gridcell, "c")[[1]][3] + coordinates = unlist(strsplit(coordinates, ",")) + coordinates = sub("\\(","", coordinates) + coordinates = sub("\\)","", coordinates) + return(coordinates) +} diff --git a/rLPJGUESS/R/readTableHeaderLPJ.R b/rLPJGUESS/R/readTableHeaderLPJ.R index e5de9b5..788b485 100644 --- a/rLPJGUESS/R/readTableHeaderLPJ.R +++ b/rLPJGUESS/R/readTableHeaderLPJ.R @@ -14,7 +14,7 @@ readTableHeaderLPJ <- function(x){ columnNames <- unlist(strsplit(splitted_clean, " ")) df <- read.table(x, header =F, skip = 1) - if (ncol(df) != length(columnNames)){ + if(ncol(df) != length(columnNames)){ stop("Model output is not readable") }else{ names(df) <- columnNames diff --git a/rLPJGUESS/R/runLPJ.R b/rLPJGUESS/R/runLPJ.R index dc823e7..3bc5536 100644 --- a/rLPJGUESS/R/runLPJ.R +++ b/rLPJGUESS/R/runLPJ.R @@ -105,7 +105,7 @@ #' \linkS4class{LPJSetup}, \code{\link{getParameterList}}, \code{\link{getDesign}} #' @export #' @keywords rLPJGUESS -#' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig +#' @author Ramiro Silveyra Gonzalez, Maurizio Bagnara, Florian Hartig, Johannes Oberpriller #' @example /inst/examples/runLPJHelp.R runLPJ <- function(x, settings, typeList=NULL, parameterList=NULL){ @@ -150,7 +150,7 @@ runLPJ <- function(x, settings, typeList=NULL, parameterList=NULL){ #singleRun$template1Mem <- readLines(file.path(singleRun$mainDir, singleRun$template1)) # template 2: the cru or cf template #singleRun$template2Mem <- readLines(file.path(singleRun$mainDir,singleRun$template2)) - result <- try(runLPJWrapper(singleRun), FALSE) + result <- runLPJWrapper(singleRun) if ('try-error' %in% class(result)){ stop("Error when running the model") } @@ -192,19 +192,18 @@ runLPJ <- function(x, settings, typeList=NULL, parameterList=NULL){ #----------------------------------------------------------------------------# message("\n\nReading the parallel object structure") # do the settings check - runParameters <- try(createRunParameters(x, singleRun, parameterList), FALSE) + #runParameters <- try(createRunParameters(x, singleRun, parameterList), FALSE) + runParameters <- createRunParameters(x, singleRun, parameterList) if ('try-error' %in% class(runParameters)){ stop("Invalid settings provided") } # SOCK CLUSTER #----------------------------------------------------------------------------# - # Initialisation of snowfall. - #message("\n");message("\n");str(runParameters[[1]]) # Create cluster if (x@clusterType =="SOCK"){ message( paste ("Creating a", x@clusterType, "cluster with", x@numCores, " cores", sep = " " )) - cl <- snow::makeSOCKcluster(x@numCores) + cl <- snow::makeSOCKcluster(x@numCores,useXDR = F) # Exporting needed data and loading required # packages on workers. --> If daa is loaded firs it can be exporte to all workers snow::clusterEvalQ(cl, library(rLPJGUESS)) @@ -213,15 +212,12 @@ runLPJ <- function(x, settings, typeList=NULL, parameterList=NULL){ message ("Sending tasks to the cores") # Try catch prevent the package for crashing # the implemented try catch in snow is not satisfactory - #result <- try(snow::clusterMap(cl, runLPJWrapper, runParameters ), FALSE) - #if ('try-error' %in% class(result)){ - # stop("Error when running the model") - #} - result <- snow::clusterMap(cl, runLPJWrapper, runParameters ) - #result <- snow::clusterApply(cl, runParameters, runLPJWrapper )resul + result <- try(snow::clusterMap(cl, runLPJWrapper, runParameters), FALSE) + if ('try-error' %in% class(result)){ + stop("Error when running the model") + } # Destroy cluster snow::stopCluster(cl) - # deliver data to clusters # Snow's close command, shuts down and quits from script # MPI CLUSTER diff --git a/rLPJGUESS/R/runLPJWrapper.R b/rLPJGUESS/R/runLPJWrapper.R index 15cf256..12ed14a 100644 --- a/rLPJGUESS/R/runLPJWrapper.R +++ b/rLPJGUESS/R/runLPJWrapper.R @@ -27,24 +27,24 @@ runLPJWrapper <- function(runObject){ stop("Please provide a valid output directory") } if (is.null(runObject[["template1"]])){ - stop("Please provide a valid template1 name") + stop("Please provide a valid template1 name") } if (is.null(runObject[["template2"]])){ - stop("Please provide a valid template2 name") + stop("Please provide a valid template2 name") } #if (is.null(runObject[["parameterList"]])){ No because it could run with default values # stop("Please provide a valid parameter list.") #} - if (is.null(runObject[["plot.data"]])){ + if(is.null(runObject[["plot.data"]])){ warning("The plot.data boolean has not been provided. It will be set to FALSE") runObject$plot.data <- FALSE runObject$save.plots <- FALSE } - if ( is.null(runObject[["save.plots"]])){ + if(is.null(runObject[["save.plots"]])){ warning("The save.plots boolean has not been provided. It will be set to FALSE") runObject$save.plots <- FALSE } - if (is.null(runObject[["typeList"]])){ + if(is.null(runObject[["typeList"]])){ runObject$typeList <- typelist.default warning("The output type list has not been provided") warning("Setting type list to default values") @@ -60,24 +60,28 @@ runLPJWrapper <- function(runObject){ previouswd <- getwd() setwd(runObject$runDir) # write out files - runObject$template1Mem <- sub("path_to_output/", - paste(runObject$outDir, "/", sep =""), runObject$template1Mem) - for ( j in 1:length(runObject$typeList)) { - runObject$template1Mem <- sub(paste("! file", runObject$typeList[j], sep="_"), - paste("file", runObject$typeList[j], sep="_") , runObject$template1Mem) + runObject$template1Mem <- sub("run_outputdirectory", + paste("\"",runObject$outDir, "/","\"", sep =""), runObject$template1Mem) + + ## we have to do an exact search, this here has the problem of + ## cmass and cmass_st would both be activated + for (j in 1:length(runObject$typeList)) { + runObject$template1Mem <- sub(paste("!file", paste0(runObject$typeList[j]," ") ,sep="_"), + paste("file", paste0(runObject$typeList[j], " "),sep="_") , runObject$template1Mem) } writeLines(runObject$template1Mem,file.path(runObject$runDir,runObject$template1)) #runObject$template2Mem <- readLines(file.path(runObject$mainDir,runObject$template2)) runObject$template2Mem <- gsub("path_to_globalTemplate", - paste(runObject$runDir, "/", runObject$template1, sep=""), - runObject$template2Mem ) - runObject$template2Mem <- gsub("_file_gridlist_", - paste(runObject$runDir,"/", runObject$gridList, sep=""), + paste(runObject$runDir, "/", runObject$template1, sep=""), + runObject$template2Mem) + runObject$template2Mem <- sub("_file_gridlist_", + paste(runObject$runDir,"/", runObject$gridList,sep=""), runObject$template2Mem ) - for ( j in 1:length(runObject$filesNames)){ - if (is.na(runObject$filesNames[[j]][2])){ + + for (j in 1:length(runObject$filesNames)){ + if(is.na(runObject$filesNames[[j]][2])){ runObject$template2Mem <- gsub(paste("\\<",runObject$filesNames[[j]][1], "\\>", sep=""), "",runObject$template2Mem) }else{ @@ -86,17 +90,17 @@ runLPJWrapper <- function(runObject){ runObject$template2Mem) } } - for ( j in 1:length(runObject$variablesNames)){ - if (!is.na(runObject$variablesNames[[j]][2])){ - runObject$template2Mem[grep(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), - runObject$template2Mem, value=F)] <- gsub("! param", "param", - grep(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), - runObject$template2Mem, value=T)) - runObject$template2Mem <- gsub(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), - runObject$variablesNames[[j]][2], - runObject$template2Mem) - } - } + # for ( j in 1:length(runObject$variablesNames)){ + # if (!is.na(runObject$variablesNames[[j]][2])){ + # runObject$template2Mem[grep(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), + # runObject$template2Mem, value=F)] <- gsub("! param", "param", + # grep(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), + # runObject$template2Mem, value=T)) + # runObject$template2Mem <- gsub(paste("\\<",runObject$variablesNames[[j]][1], "\\>", sep=""), + # runObject$variablesNames[[j]][2], + # runObject$template2Mem) + # } + # } writeLines(runObject$template2Mem,file.path(runObject$runDir,runObject$template2)) writeLines(runObject$gridListCell, file.path(runObject$runDir, runObject$gridList)) @@ -145,7 +149,7 @@ runLPJWrapper <- function(runObject){ outDir = runObject$outDir, save.plots = runObject$save.plots, prefix = paste("run",runObject$runID, "_", sep="")) } - message(paste("Finished run ", runObject$runID, sep = "")) + message(paste("Finished run ", runObject$runID, sep = "")) #----------------------------------------------------------------------------# # END #----------------------------------------------------------------------------# diff --git a/rLPJGUESS/R/sysdata.rda b/rLPJGUESS/R/sysdata.rda index 137759a958104a3a27e20757001f3faef84e1649..73adf1098ebdb80da9447d7bfb7a0493349ac561 100644 GIT binary patch literal 7364 zcmXwd2{_c<7w~V)VC*xDH8jT9mn>tA-kGt@U}TxG6C%50OHtXfg%~@LeHr^&gb-y9 zS+XV}SrR4G+xve0@ArM@Ip;j*+;g6L?tRXE&b{a2Z9Oy;aW>)>Hq=h+KL9%H!N30x z{=e|&ONPe3GuSD#FvwW+n9Tcp$J!9EL6Jxicrva^z_Un+4Fg~>7{)5b zvUj37#<^gmwtg~)>7>A2>yB&Im$)K{1 z8nL+Ap+ox})F}g1i0CgVP6HR4pJI?D44h^z=6qom`<0T`6g*s^0jgW9+X_qYWO5%D=S7|2~+wahL5Ow@_5sjoz5$ z`u8^(e*2%C z!clbabyTqD`=f>nctVc0z1fbd|26 z_(aox8_<$edT=~-p(IUdB%e98Q5s_Yb}TWnC}-5pjBb%~G9@z(E z&%5!AA@C3e+!C|B^{-16$6b{YwI!76={PvM=!NmYz#If(g&r-J!F(j=jb~Q-@NJF- zgn;m>$5)$i{YBxp`M6^gU&Kup-_NN?y?Q&@2p_Ffb34PYGB)&T+DnOn{`qK=2g`v62pSVgseBfc=jN&G5|A!DpMvG8SX`3ywA*>hw>{)_A^Ru9KY8RR4H* z{4-@`yFs0;Ky0p{|1s!I6o7FFn0i-KFOT;h_e>$7fQTLF&UBmc0X zeY(-qbXe}pHJ7Fo`%ko? zehp@zN!qf4IE0jKJ|1+70>QAC-xShhWg2L9VGmOJXq`Kie!9tPoP_tf2F6qM_4OtW zQ2z}x@_KxHlWiWSYgJj52O%~lg_}Y#w_eCaWa3MxRb+F@IWoso(6|wqY)m|R{JBvj zR9XqOcs*lw#bLbHf-<*or5pN;b2}@!eZ2kBeV1pv3m?~_mlK5o&^mYmIf{%ylgx;C zniMH|KtP!)mqsK{rp=y5A>NrqWQriDY)l2Y)_MGhObVjoj68P6U!0#?r0mITTmzY? zxN$@G>x9AUC7F4Ly0gC8(J`r~#^6K^IV=C=6dOoqTj$wMY(mLw87E{Hzz7)vzP8b`+H9t8ki8N1%2Rogb10Nc*1QLmX^x52T&Fiu?o{BJBRkc(+ zJmXQ4apiYsH$LzBx2wVxOi4v|C3?#iv70)>>ZWINj0H%cv1sMPe{ot%?EwA?mGdH# z<|U*1*U+VW+21Q(U9W;IZ@o@eR({SAwz`kfx-TEQvVcGsnrDdc)z z(Nm_Gy7HIZdwp!)(dopUC=Fo8Ch)0)s|@0GJ!v;t@%Mk6n8wf@bu~174Q>h#Q>jtg z`>2Z~RM+pCR!I#z+JvoQnX%G?iiaoT2Dk(WENbcY(M%oEm!i1Yg5IhUW?fAhKy~qK zNmEa1m&eTRN$H=+)-z1y)WQcL@BRt2ag5UY8ptfLQNhHqz1eN0VZOa+R%3=guK6N} zari~4(Tgd9`f#?NKExyKXrO-f$dR&ClgKoRaam+;M42v{5TK&3CdsL!1W2d$nZL0Y z$~Awr#r`3P(LU7fZI;WE0Wsa;1}8dM2AxWB$raw-^q zX~Y`CjJz}Pc}AVCoXC?IP78rY-0Sye$EQYnAM`Fy@IAZa_4AT=k#CN54HJ!%oL06J z3JLP}FRFHmi!8UX)ogXW$QJj&cb$Wy-{$B?_^+e+Mo^1YP=EcngVWwWrBBH}Gd{Ur zfNNyP6qcxTCm+o|s~+xU>JVndN!%IQd#-XBNDUtnA%&)Vxyt)S^}fJK1Xv_y@I-{m zzaZU|oKTM0 z`UDeF$w&P$i(+e<(}jjb8xD#YN*9ui-(p}DEvzllSM_3vW{^xoLfM5YQ zmJR^;OdKme6Db-Ud{$Ozh*9Gr7kq4JHrsZV zxXlCIagp}6ajPR2MtTL~xK(PK=PD}J#7wGCwX<1wXJ%F-kvmcUd5}nhoRws_rG9Pa z9PNDJnCf`+P6VRPu3oTQ%u?t@ zFF$0H&^X2yU0rXccPpJXicNJE;wy;isuPOotZijt5PU!9lDO_*Sw$~fWpgYfJgLMU zS#A790rQqtw_ai=vzq$J-Sp+NO=zb@-RX9DId7RkG>iy3qi2v0Cppeo2^Lt|LZaDy zycPAP>7DcHP(`nr+JfDZgsatGfoh&9Ffto)&t303#hp4fC8WMW@5J5~ zeG+reJ8xvP;#S?fgUaM?tEPOVb9?H=-P zMOhbDg%^4-7CG_vvu+Z0`@TZ0lN)5uay(Z~Es^Ah2M~ja$&?u&B?UQozPO;8=Q9`bAuv}u>ADtMVN^#|)tVHJ&|b7D??ot@l3@N`o z&9PrRlk?8Mv0r$u@^9rPEuN|Tr+V|{GMg1MPRa#ll8ddWMDf9#RxDGfD*sN3{*3spqM|3Mke- z3;Jl`jZ3PIp+57z@9m#~MVQ!J_uj6~;PxKsvz1eAF2Cb=SorrH^sAdqbB3i4_Y0Ft zyiTiy8hBzc<`k4N1~4nhfocwMc#rQ@haqG!_|iUq>*cX=uD|Yu8BGFgJ<^1z*SK zAv3dc$+;-NoJT;!9N9oRU4}HiMu$7=GVZmpnTfj3d@$BYpvU1l(%`Tym`1632hJ_v z85RZ(xP?bT<$(>m!v-;9`(o!3Zeiw9=cIOcve#0e{Od(c_qVb9gF#0`T%WqP7WXB3 zCwg`E3p ze0Lz%^~Og)=^Pyd&)V14k}r20IKD3Q&6`-vz!5$co~>~oGTz5R zv}IVz_ZDb~)0}5YQdJ10R=vKM&i|CYyZ_tAUM<%^Uq9RFs9a|7`_@oI(drFgD6z$yP^n_84cm!oC} zFE^S`hi#_-zG1pctY)W15B*FGf>|03Aw0j$)xEmjqG!U_pCihHSP<1meP%-WBit)CbXYuy5KM+^ZKPcP2TYb}? z^t_FZ&HhpMc6Gn}M~8E$F*+h?<`MqisOo1r-9KP-LX^j0_h%a;{6( z3)Hz%iO5Z^lIEY5qoTyiyFSfB(DH;%qDWD>tF>%Gk{yF(4>4?UpWCa3LQZaDim7c7 zx3PrU&W0@wMhD|1o`TC_%zEmIQ}e#F^WQlAZN(ZAAimb|1_Fi6HJH7H5@M+I5h~z~ zy}D+(I)2O5>aY@;U}iIda?J)e zl_Z1G*rOIY95z)ChW7U)P25zZsw@N*puQFn<+{3GLRic;2-B$@f^ed((7Akjt=VhU5jE z>feFnwJ6x8FJb89i_LMW9@L`cmw#L;E|W9W&N=J0)_>IYRjQBa`auBx4XI%;c~P zXA52Xh?4BY1gBa|wQ)U}E@H{ID99M5nA{Ghn{!zCydNldt9)>{n_`hwH6?+fDVS85 za$FMrUByG*wkpj|I80M94H;AcAz|nwQBm|b!#7;d*J-J)O)F*_MnCoQZjb9p1QT1@ zS~ZYJ1^W9jq9l0p1vJ{3F9t;wn6E=(QwIyu7vg37I1t*mdwzZT-kz84tl}qsB;^El z`lv=>)(XoAee2)orPSXn!mUtJho>f@D3734lysCOj&(RHR*CK=fy&fXot z(!pJ(MJ1fD5?Fi$EMLtPnUN~~qmKy=fPe@XaC<~x^iB#7yNbkLx{?{@FI0DMi;n`- zWQs)~dXhykBIew-8`i~4uvBy>I5olB`l5WIPOVgZ!kf`YFJK{F^iYl((}}eM(;>(P zlM?IS{SNcBQO&^3V4mYj?2gGuE|?HI3l+Y#d$#u|UuB>C5!S&8yZnh6vtTkyvaPr! zFso4AkUbU0GiFX~ts%N%ku24+$PODsF743J%K$KLw>&^6kCxOI9*TI4-=c{>H0N1Ud5{?nU`rq zCZSzC=Hw<22+g|}^>-{l&vdM*iby?lTz-F~E1y(_Va}n)R=kAv@{iT^UsAOjkjre} z+y^@sGuT1|qlbpnq_0B>%5-AzfP+M~Um7}GwZ4orG%dy4YVY-J_SXQsgOg$nz^FJX zg^Z>KK(%dXx6gRlMxNw2p62<1A?VARF(kjnK(OXs2P)RL=Zu_0xD};#=jg<3HO0tp6)b~V&GF)1G;4gdfZk{dWgUa(4w=b> zBFuDSV&JJUv7eX2KB0zoMe-4!KiNt=sNlNBnR(;Y2!8}aIg-8Phv)dq$0-|Qpij#a z5|OoPOSf{K2_2G7&C1myBkEEyhFG<*g;OUp77Tyvy03_F3_9rua-q-s7s~nk#}b-} zM{rA_IT4})u`k1PSG&b#Wp3qfA=oMswO(;=0!>6#e#U0Qo=7>DxHew`r`Q{<^iY?+ zNU(D)xU_5`^)e_((nxY});C1FS1wqte<$nt!>J)sQCyips{#qG1kL5t<)FJIWXOON z)S8s_RPa+xg?o{3vIL4%1Oqf)(>C^vLH&wa#QXghC*Z_CH`X-_c&pY{dhWLQf0(BI z5=|EiXg+!0YUn@}=HyaMbi9;(Pt{%~KWgXF&E=KCmE3C+q3V05Q0;XwIKDC(JOEz% zcKa&dx8ASV4!rrpuF%z0ZZnsY+Nz$hx&r>Z?=uJj9~DnT4Y;}z#Ka+kkKDE3x4|oM zF;`3yZS=nYY>OG4`G_^2-enY9VK!gLZ|R?#8)-B%NN`ZlaBW3>X~}} z_*G&7t7aqJ>{V4A!M^?OO?#@c-g<)uAe`#_9?Xe7K=~OF&~&;kN}z>^1QBC=4!g=g zmC2Zt3SGycwUG6iAGz8qR#!>MjQ#k?PIx%qVaF%W5&#s)RB_IIl!Jr02CpviQK-O| za;QPR^lc3p{caUvY*j?(i`RS13jJl4e9q9emE~-74V#zJK8^MaQLC`KPKVD>tIkUi zps)@qliUX60-r?Y$b23NyDvKH6f|k!+2Vm1D4Ed}4B`%N7Yq;o$t(Nu$&;Kp??9v+ z9*-xolQDkPBiolA_@t8&Uw__B`6ZEiDE~ee9%78_t!^@m6F zpMoQ`Xi;{~NBgnc2F?k(nTm>@ZWP6w+3xEH->4J8|0>&1V?A%2nYZ~m~PpNRo7!usp=n#h{&Rw8wfE>fVxF|7Vt@`bq77%(g!iiV5Gp$o@I-Mc5U~-edF0 z=yoRLYIU|vaZWKr7MsFyBRz-31N{zh_CfZT5)|F8}u>8dFt&xo%J zoh_sB)6>&)&U8;tP4#rkn>j%Q(W6u3_g_v18SqUOzon=gul|h$9Xy0x$|RuRqBA{1p@=$)p{i z*CD;dDVcgQB$M9vNYaWy&(_HHMy~1`>UJ+vmL!xo6|aQ$1c-rG9!cheTb(+$>qjE{ z5Q^c1L0dwLK^P%VtCxUPWQl>H99ha;L5iP%g~1xawot^~w)98owOVO)g28;bR)y2j zgEVbkQKmDAZb-^fNL-Oxl`=C1#vH~IU>nQ|6(r*|frxkMzI$(hCYgmMK} z6jj#{Nz2~q=#`c)-A7ir2s@WJ%{gIavXl!)55G2CkJJOn$}-1Amh!9Gx*H4?)&4As zQ4t`P4EE~xJxT7ahy@B+^0v#yWC^*IYN^zN*vjSaDElSK)3s7vsX@ej9m9KdH!KO%OB|v-M*7&IJU>?2jujF)# z7H(&gE1_2=HZ10NMJ~jEi+@KhlNo}eLZMEoflFGrF?o@Tof%1d<9_~LMPAN^rRN#8Q*lk{8OW7*&&KX1n7ih93io2aO+OS08tc- zsMRgt#}6l3BK2~mo1y~e_muz$eg!H}4V2A^7rOlbLa05Jl5D8C2DIbCXSYDpXNk>9 z1DQs_aHtkmRQ*aOs>qnRNxQ62$|Fcg$_|SM=lmr#*qKa?8yrN?l0`x6iHw4hiS*h2 z&VK`|`f#4~x$T_?T+>OYp_d@#hn*i==^?kgehk{XR0U_JP!4i?vnJzF6-X6U58MQq zz&A)zic4y_H#d^``jt%#n_fRNQ zPmdA^VKjLs_rgIPG>PlrNob2jMS)2%6ciPQ*<%A6B1>+Jb}WRxvQJw1F&mECBUAi9p34>vd~Zt?sl5qMG^v z4*vf6U-EccsYJc9D=hLg5LY+~>+flg4l7+54C+xEo7{Rbo2EJ%E zB$HsJ)wbzA?&l>1>uP&PVUCcB=PO;L8I#r`x;5>po+#Av05g;vr4~?0q;$KsEe4*M z@}9Rp8MLB2`7p5r2Vs$QwXImX9owqSI3}psmXN6#D@*Nkkj4)R#(a5`CnPfC)dB>l z03?wH_)*Vnv#YwXwd$>zRcToZKYhfKdt=H1XsBWlXbgQIXzx9|Hk^|Tp!}_*sKGk#2h^eG1ixUOrc}E(BQJn|oysDReP>gmPCs{iyX8|pkWQd(bQ+@~;C7@Ks2w5lC zCWTw@k!k99tI*hiY1qy2p|2Ba1(U(`<|0UT!(k3X!7UhVi=tr_mAQSTnY|gz)}d{= zp;=84N){I4)ACvk7JpVtm1v(yoi%1Bh(ws0a_HMD<*C~X>YL=4IKUjB0B{R1iH8cD zLY)AQU;}Vf0F(e+AX`iTx;bIonYA9@o@|v;w`HtZ%>~MC*w*okl5806%myf@B2Z1l zOc+aH zTIMT>IJL}l7YaF^8~NchFt3ZYu@vmE3| zr#`XJWJtC0N00JdRMIw{wcWv$3$xV?p-|^z)?{!fN$*1^e+KAws^o6|Rp^GpqC!CA znh>f@UN+Nhf@t09TgUiVL&?VbTCFeIg^)9wdzUe2B`zR zYJ2@`N~l@zSlL7rrJ=y~)|$Gf8-7${})}i z+!FQUTi9Jzx7e!O`rjP?C;2)gcx}8L;THM0zpBNS^DjkNP6D{4~Noo6ocX zc*U{}WZ6RQ8c|f9T^m=UlBMnH{fHsBa{FyuXW@;UL3)@zWvX)#cNuk=g2 zTPVLx^HJI^a{y7&U|iyHjXzB94T&QAohDY(JgI3%6s&douYl^;X@DI@+y`;Hcc{o9tQu#Gn^zJP?eYq(B0?zu z>CUUe2-xX<~b4NbPn*CoxJKL?*w=kr4}9i>{aGJ^|}4> zxl@F(dem;+VBs+Fa5_Dp!1Onr)b{bXZEaZ4B=5lSkJVA)ctU}bU4gT~?^~SX$455f zV~jP9XcI?g3%%3DfsZ|6A+(?WOhITgnZ1u(fU|d%|Hi~DeB|ISQT^vnEI`c{eQpF8nfC1O1_td+~V>JuMo|!#Rez~9zMBbAK z@=X&SGhNXFQWH6fkk#b+h^PgwnEQT-A9alxLL`w>NABI3&4+mG;{5|CPneZh>ya2F z&9i4-WkPkBjY%ldw|l>bi&GwlPv9suD2V7XOo$Rw4K_d$M>zJ9ax4O!qA5iEA&Yg7 zJ~Df=^q@xVG5=q27_C?`o)i$fh65+;h3|%pXNxsQDHnT-nU+O7N23mMKCbt775m5% zSqku&ca8x-|I?F%2b#-z!+)B$eSdwA%|2Btq0gx++H6i&YXv`L;p-%{WFWD6|qMYU~z zCbt{QZIC642?PKrNF{Jc^uC6-PL8G!2KU>2N<28~zCPnbwDJvg2L^ilQJ+G`&2++lg3 z9TjP`piSdiYW_lBi}p1ZUVG+otCp>9#zSr@mps6$`R`lxTWR0=1UKFn*xyLsJY&@lvXNe_sZ>7Lr%H=<4M})zt4(FL zw;34dFL#GJI`{N?u_L_m+Q^0(SuM*|UY#J4&=g2PxsySmu>H(bPQK!q{yfw-VZIWB z{5XTRf-Z_R!_zvJq7U$~9;om=ncS6?if%0f@mfpK))c!an4J<&KPbr$CmE_319-73 z9iOT@!u%^UT=xs<;|OcId+;a8s6*PX;=!NKM1S*MiIFD{P73UY`g{sDjdSo(=BGXt zqPhLHGe=XY7f16PjC%WhqwrNlpj5h$^dI{9Z+42FSOW!8ix5ACJQAASkS!wG>D9-s z-uMq^+iLTwt8lBTYy0%Y9K*KYO@?&zR&~cW^RDuLR3&;lgb#k_yw2Jgx?+a$v#QS6 z!g%m4cs0fn`jzZ$#@t8`mKG0y%PrVP19q$KSs7 zSukp=G>lKDX+eu3<$NjM@PFW`h=9F(1yY1F^V1MHSp76754=3F(o7wM371bN3gQf7 zK3`&pBNR?^KNS5$JAPa-ESvjY%sY{sV>$Q!Ijw)b_*xaRAPnQgjt^TkEsjk!Mu3L> z(Ol7)3eK<%n)IkzAvW}_SdKPyNpmmg^M_n>tZ*nJi3|%ztN$54KMJq9zS!1ov;>{o zIwV;U&T9&*8sgW@DauMouidL-WlfvT)NqZAr;9MUhy~C(()}q97T4>gD_E>s7%VRx zzR!32tNszuNdGV5$}LpK*50UI>sQ5fOTZe3URI#nxM+(~=`AnPAe0j?N3EUCltV8^ zEi2j9+sErL#nn<3Fl8_8J<0w)X@vL5anT%vi=Ga(C~V=L#AuQS!KQn#Dh5GQpaBkM zJ=pZP0HvraIv5REHG7%v7&I={EY!jK!;(R?kbY9iU1hxDTkOS_vcvk&#E=gqE&aD@ z%_T=30m+fXpykc5c%Qd?d|}%YkTtp>(Qee1eo^(Jn06zokUr=0huu%?Zna;|gCO)@ z4028f3S@k0d3nAajOA;y&Ls}%)2PS4B+w^4#LQQCC8-TQgZhaVWAz{u0**UVrL?8! z{Do*;Ky>*0d*D}e?RsRo>Me8}RMd-E(^_vi)=C_>lsIJj51(tOKbOuXS?yuL+Vz|9 z_-oJMO=pioWvzhHrb*uB)r`<=eKsc`fES?p+Db=r_fnH`-K>tJ17LzR!8Z4}P!ZJmZ$p+>wROs^z6TB?iQpe>oUiO6Bz)tv`q z0~S}MTi*?(eKvk~>u=04<(Ka)?U+pdo*jpz9*>aKg};QCiQHu7j!y_sQc=(0_oZ&B z@sMJI6`0PE@wCH3=ua8?{W(jA~&lDx2b)O$a}2|snX z8CGu=hZG1hoC#M{b8wO&2WZdnJ|Abc_mqooC%>1FPNTHYma%E+J;u{oj5Z5Z9?7(e z_#~Q*-z=JhPSm}6Yx(;i`FTGbbBYLqz`K_}H-zuSIg_)l{p{@hWqkZwoWF*xU6jaU z#j2Y*fBS$0R~izB3G~-%LkRcq2xTDP9*5n`ei^Mdgzq*1cxBAFy$o^RkBfs2|HXmt z_}Se^aa_F_h3Q;=o)233OD4%wQ#CD}d!y;rkd-eMb{B}77hXUc7J!;Pc>VGphU!Nh z-oTY(x5q=wP2x33Ib-x+uc#3-a&0DW5;%75sWR3Nol^z?GO`ZDlphR>nsvbSxf=yd ziA~YHe%%y27UJmNAEvMpXz}tp!$MpnsYT5gJc3+Z*8o0xF}=2q6b$<<_nk@SYtwDx z*aA3oMrPr}eKIa(zoBu2V>xg>e|Yq`_gJgoZ}0Ff@YM-ZgR$!l-k10kMUMSg^LU3- z(9)hqk~EfpdKhB<@Xtk10qU#8k;i_(g~g zbNz6f({JCbngtlHbJjHtav3MZ0k@5Co%HYV4lUpvD;8bF<>riapm&lHd!I{8+0S4z z2B%ZFsSTsKx;maz&md+i|Aw6Dr%kyA{#CIqsD@ih4(*6phmOV%JC85cR}bBJ&fQIV5&Pq3W4TI(V<*jBe>YU)SjkUWz&Y_*Yceu9D7PxxC2seXhh-l)8J~K^0b* z%O!vJWa_lf7k$>(4Rvk$>$VE9V1#dL(Qi>V#B#06{xeFpbSfyGcQ;Y$Jn2{X zhG>%@`MpNPjeOWv_<1taCe~IiB#VRfNXUb!{2eF%w>B!nmTnqKokX=nt>vT@89Fy( zTu-JqM$hYTP z-Rh+``(wuO3-#VyP`|ZDv~{xYJUJ^&q)ykw!G=TLhz5qF%EvWn&S_c6st=D^Hlxcf zC#|9+Cv`!Z(Qqr82RQ47`FArRk!eLN#R537q5xq+rz zRrZs!nnEO8I5TT>S_Mm3)2mut&fW|0Z|{?F2K5`ypgjJDlY_i^PEC#jb7_|f+z5PX zFxd&74}g>CBFE;uGj?K2Q*@~1evJP*t{hOjH1aKc<;!N83eI_zy3*3<`fkHjbn1^c zk_R4H-0ZQO)KRYCHkzY9!@3P&|^9v%7WDCDVHL<)-#=vn=fz^ic-ixkj zw%gRIc4P5YqYdBF_D!S7x1d|;Yc`Zx@7J5*CK*9D`>|0WO*kiDyYrAuW}9DaH{iMRM1*h&(Wmt6`sy4n7EF3|0T{eqmj+(VZIuV!1ZWDsr0 zv6ljlQ5?>nSEJX4zEfm0CHDI{YXYYw?4s6tS@ldHf+=k=Oeo($%P2F9T%(U6tzODg zvzV&+rfUaAz_#?Bu;|G#vV%kQ%AMiyFj<+(;QQZ93u^e_K#O=SOh_Qg%K;q-wwP#J`)QeX4$ zCm&oP{JalJZtqe`$KryzvepeZwk7zN_7nL+Pa5C!FZ($vH8>MV2%W`*i#>4ISNysA zfupD(t0_WR)-jl)C7XOq*Ot{G+{!c<J)>+Xjq^PNQR-?!bm z`38%gFRBnXbou!^H&0h{7MTWp#lS&e;)4)$Q$ zCZ?B^6n9>wmnsNwAMW?S^6)7-5Hv8|Jbh1x69g{`}uFtm_H#S1A!v=NU(n}w0|9X;Ny(i#v#_nKi zYg}RE;CQFQX1r0c#H_Mo)$iPczdY+rb*+_r;uski7@!yojIHUAhA?~kFwbF>q=tMc}<3GHPmG^oybi@kcimZueF4ILf<1+_6UrPtjZ%X!bId5VPw%YYP zeX+;m&%>&(*VWZ-wU^h$#(yO3t@Zlu^zir0i+XL6u8o!NNzTe?t`@<>vP2jL1@vprd%!s(6P4|i z+YQ&StzOSiPn+Jw(AWqJ7JzVNTPkZsYD-5>MxuLVb4y3ITSKKekF#4w%W6l)_FKva z=$UaT)ap~O`cL#0)M}x_Xkra&&Rh8T7Hn5JyCY?_cqm%*iVb)oof8!agd8dXt_{KYzdD^73+z)V7n8(@B*8)RviVaLYH| z`?oKaHTLV@P|g^qTm2plEIqpKaz4c)G4M2hQlhT!G;zfW5wEI(9fjaG0L0z^6R$_& z?q{}S=1h@wA<=_!Pv7db+0+C7+@4$9kT@QO7VRr8+W80eF+yS!>`~Hiniv>+bO6_l zy6Oc-XCppN4GkI3u{QVq_k!@hWFsT=%rujIPKITAP(2WWz9jP9M8J=u?f~K=-0&v( zEF)jkAZt}+xgI)L?^`qu!F_reKvwm(p2Pa6Yrlvgx;wH_Rn;k_B>{eNm*4wuZ9E3WA*8+14;m_tb+Vg?+E`3xATq- z67v6FOJ1u$CRfc+J$qt`ZCNG%RNjN7x67vFIOF+gwhmHd+SI5I`qE4)#^j2TzUGF# z`kz4Uzkk13Us*el?;wK9KR=U@7BbbIk|^HW#sEx23#Fe2zNM zE@LQG9W+y)nLrco$YCThS@B1RRiE+t;LFAU*i6YOVj6j4Uxnqa#R(>Ib7$LeG9YNI z=azcpm21I+%Ta)b)&Bk;%oIKn#KEnzbm4^z(FwGt%kHQH;iS!BjZ1~Zt+j`~HiX)L zaw23gFn~|B7d^+owc*9(_9&H~k5Gj? z{?sNTP{|;T$s>sKEq#Cg<)xjG%PWc9Rmr>;BXihVc{9#FSBJ9_J_o03tH+=rR@Y?P z&2H0Hjl3yRH7+3U#LQLU(=r9hIUNb*S3hY6`Vwo&r+eDqT_8WWhUeUD&RuqzH6ck) zZf8*0oZJZ>YIRmq+mJLd`sF(iste*_7FE@Vq5 z`xXZsd_FUXuO6ob?R`3%-+_h`<8iWV3Opvr<&E=(3POHBrMJmrwJie_-cV zuO%XBAabTnF6_oJvumZMeoJs9C$C3fXYJOn`3j@fp3X!w@_CeghLQ(FKK8@5XOQN~ z!^G76XL4o@{NR-avyQTQ!?GL=jTB)Wi=?!Ri5?@}|5W@R@D6KmIPl}`Rd?n2-Ek)|xhV`o{>TyT`JA~H$wjk%p8Q0e4r>W7L&8^lSL(n0B zbpdupFpo1ECT_&j&1O(sTr)-%ASs>zqGaweP)Wn8oWgZDmJHic!;-~+I*9&fZIO=Z z%RYyG+lyLa4I&fmNLcq6l9fh)Yxgo&t)1SO1eb#WMZ>RoJcP?dXftGkiKB_1!6FEY zN;#V+203LsOdD%!c#Lx);E|1WWa}$&?WL+Mv2?b#+b23~za*WrK?tV^iMHJ6Sw~>LxdD6DM#*WYJy#NXglRLyUrs%b)C5X#JP2 zgNQ)Badmyg2fcT{UPP(R;1H!@;$TS)oQ7Bl$WWU7T<*Vq3$Tg#YVJW@{3_S8^6+(+ zO&YZy&{q70&su(s#m|{{sXt?EQJqyti;c&qipV~v-l&qIgE!w&OUO77UcE%wnoKLa z5%}JjAFldkrq+nv2=Dkd2m;K|cHb^j39%ei#ETk?t|%tdyo|AB+qi>5J>cZN7YOlj zmgO&-*Ma4GUTS@pTQ@Sd)$z|zn$isDCaP(q33GOtN3hb11qX*{?UEWU)(W`_@~m3! zi^Md4->EJtM=Bm8z5qC>4Vg{-GzM*}Y56L87#X~6=`%brE?|&vc+Tpk`c-=c=ooj6 zJE!7x5%Ae*B$@Wsla-DExMpx`xky%d?JYq5F;}Rxd2913j1MH{C#WN47t`1`MbLgT z!;rRAM_z}w1v1fs+N$K7GXE11UR#O{@PeAiex!miQC+7jsR`Ut-uF4koT1Sl zf#xGWj)b+;ep&V{^fh$mqY+7Y)izBcPE3Mfql6oH`Dgns?$bV_;@@^6;*daq>=#l;UH2auqzb zow(q%rQRAA)rkX1A923Vhi3nZ?>J`K+V7p0;7rMkE*x5|)W*lDUHuw!X7V+k__fVX zi^8;K)_~J$TUN%-UQsoUIGaYPvbl=$B?as8xWx0cAdohY9yIGDdr2YYygGlExf7(A z9!$O@O&cI}p{+q)rWEHbKLr>dTr?D4)OJ%bteW$-5mcpHDA&_pn6TP*o|l+4Z|>P& zjBFWCfhE=Q5Oi4BSg>na^6_Ri!yI6F*^bE^?8!Z>Zmx`qs%0+hcyk4uVqKNQ4Xj}} zs{1RPoCrd*7P1_+XonotV2+s&jLys31Vw+X*J2rP04=Z~k7}z0!l^tOV-o*mBJLJjTBqry}bnUW_hARnwfk9W`ajz4r|c%>j8(58J4* zp+9`+auc1;aodjc18+qr7Lpbq01cFOgx2KwvWk41m)0ya5(Gm3JRFL`QblgmCbw8les*_9{hBpKdo zt$Q2eky|#^;7dO`9gB{4hLIddXVx>tS_PPnG(7G ziy+%&?6Xk*yFFN)+)^6Q5#*qZE#x%BCWAfYQ%?+Rs_F<3~P> zdSXioQ4kY$UzmVmGUzoA%lC;ynThgZx8i@c)6K>@R@*wgrSEGV$>R~)*Pat9FRFc4 zy~{$a&dJfCr&Cwb?%S`z+&e&28zAn`-NjbQ`9w%ffpOQE?XsE?6s_KPzjY;+Mym3q zg)=4|uM6cS@mFv4R(t7ux(+HX0Ed9X4M=Y@f|*;M!>Ht742jzAY3X8nKuR&>O}1UUgaqT(5c ze^XLpj^cMQdNx%z=YhtRNohSD`}Lj>2udJY)4e%myl1gWTwwD{4$5_}m435O)Z14J zQsA$_IP#TPEM@IDk=8Jw&QP|}pjJRyK!2z3bo=aNVoF;)UN};9DgNqLVvOG{RRaB* z?tICrnIos#ctUZRDq6KxjGpxOtV$$|o2~tD`9~_`Bj{^Tc%bT2b}b^kiS4W7chR7P zqK{7SozPjkoe|q*sCvhwQfM8yA(br|Qc2LncBm;g=dxqTb~wo$+A^IlN()WOs5MeA zHLB&Uv+W~J=EUI`)*QKjH;|JFyHZb>9>-O*px0X_sq|CabWUggWod<#RZV!S{v=rl z`A7PMZ}&vXfZA5SijcpKRORUrz!{o)s)6*HSGje=tIp;n2T>)Db(w*vtOd5I$VPp1 zHwpTyHfbRXFDG|w!hqpfu`7e6O5y5S3T${FBXP)AD_!GF0BCMy~qzWDt(dBqfH zG)Y$fGU58f8{{7!!z-4I&1^!0w3|h>XOI%g^o@5GWe7y2un}g-~KYEh` z!JNYT%M8WOphgEA-Qq>QpPTxc4D`QO?s&ZU2b5D0{<3NdHxy#k_7VXlZwp}8O3yMP z2(kRw{QXM8EELb-AD5t^B&RrB(bN%pnhjMUhdl{i2AfT59pQ@YkF*8g3Gld9af+I7 zafg_*E#oJ})rT7yMaZ7V^# zayN0#tq_IqfXu874<&pBr7&dt7F1c#dGo43*qXBlf5dn-@nKOpz2@!e=$yQx4JYpG z6wZu8UJQ`|f(S>+)*HK;Op$m{cOV%on#@wpe`0dK$I10-Kdf7;Jj5zSy_ zeAcbZOr0Qx!}!s8a0|s2^UxJ09(Xb6Cx$~`NVKB*gGuk(8#0D$<}!4xl1U9(7d#9F zCsWkQwRxb z3)d$>Q^{#3gPLl7R+o7P^|9Th+bmxd{ABb!6dpRelw@j=0@i{#qu3R1Q{a0?S{p+j znok)%R5$g}XIrDuu{CJXU}@(~snaVY8E3a-eO94Y%<6+0Tsf`xb&ziD&|f!baCAV0^IGNe9#wfv>@RI%XK ztHcCG-^j1uJTJzt#$k8$Hns~BW^h~PiYPc4fhsqUxguf6mcU>=@+8M*pVJ?HFvP1li2gpfYzyH}xjWRLR#7bWqOvTCn>}uh*p-)tjaTJ7{AA z61KFMQ!t)Gg%c7I_#+$u?D=(OWK#n`YFSk>oC(_6G98!hR>Co47DmLAGl6T$J1(^L zJ;&6qul4AI`j1H!--X8|ba39__Li830ORSvYHU(If-aCbFlPQ&w#-8Dw1Il~sm|Ar z-Fx{`Q=F#J1KLOYErUuf4H$XPkd(1-2YmtPU2xKZ;yXCFYzm$%s0=ik*d+}lN_S~Q z=4}5Sv1m71@0zNgn#^&^)|b);E78uoEXfPQ;RxM+F)GRl$!D%YaxX}~2&l|lK_tn) zz^`H2RMioSdGqUm^PRNl#`{CcSAbK=IZ1CL2pYybG@?q0zxAe;ReJEOKHFg1TN?+J zQA~yhBd&%KMi#&l2SlpRgE{VOJ0(f)$14sg`E&qqQMg3uGQ(Z}Tb-a!ZN)Dr;!X`T z_mcFIp)B9ebZCHQ3XfaYu?Kri!UvgyeRU23pH;Y|CP8&TaP_C3X&?^xokros+U43= zCHaWUQpSitT+g_-2>DC^_rT3KCnrTk_?NozqY=@fsb4QwzX@BuGVOq2Puol}zXpn~9oh*d2` z!RKP_1{A!ovWrXGk?gss^=Ml`w+fW#TB|_XXBTw6^VV(#Hw?`xkN`~-zN@$mjn!`p zUW{#aOE6_s(bi5?G5i@S>LBpzgPrxCw0sAF?bNJ2b1E7)lcQEpiv}gD8sU_U09TbE zVyPX;(O||}#%jyK%F6X*RtX8FLIkn;Ey}jz*+Tg09K8DAr-UaB0fbBOB-eUVVjoqN zY?)m_tN>7mLj{2ETagu6o3CTw%rP`?xKETwjWs)%;+t^6TE#%ZTj zA)vQfOJcZp$XZRc@+z63`I84xWsl(xo)NX6TJ9*@aGu|WIYyMbBh0&^5RHLj*=DBJ zwN@1HSR8-?eotlGfBU-b3~JmtKR;KP@BDp1`))PV=>jgP65@$}rkaWTYrK51*nAbN zjy@87jXki@{!*FU>&i~Kfcz0q$ndqkIFCw*BWV!r-EX$L4$mvtOI3193-#B zM07mBmT79!qB6(3Ck1>#8-j{B0eiFBvPqKqMEdvuX2XYyf0qXzaZj(9i(S|*IBD>s zt$$9BvNNc6)|~2kr!ZSMUC0VHV_oy&#{GuYB64_tEA8VltCC|VhpES!Sxw-y|9(UXH_h#Ey~$71p-Nxx=R$b?DX_-VE+!(i|Os=YwYA);q!Yn znf-^vnQGA+Rv$H+Z?;NN0Yux}vFxnqQCC9m1G2IJ7Tso+U#&Bjl6H6Ypo7eC(!^;C zqD;~|$x&w${n`pxi-N-T8C=R%`oVTJvBPE_7Br;74bzdby9C9Zew8UF)GmT)Cq7$b z6_50gvkv|p|C7uq+Q{CuR=%y$MFHmHbSKc3ZwBuP+i+`3rf68p8?s4Bd`j1aZTO6h zi1Vn4n7N%f4@1rA$f-(+E1XsvESZ5vA+uc6TgO;+`69 z^ZVcVMYRWmh;NL)P02s}nq-UEt6tgEx>voI=IV>#Jop#;&h^E$d7QaTt*hzhO$Aiq z7gkfUG1RZoU{3GXz{K$l&W`&~M28(?z|sCc-S58|fopX_!fhO}yrgw-Y%AWYf3gMg zO0np^t1&U;Ekd!uMh31u zumjXqa9A~ zymSI-mfWB_qB3wzLVX}mU5$55B#Nd7c7r4oGus6jY4JV>FZOOgdNEIJRur{A1$KHG zN6EwZO)Ps5PJe|UjAR_AwLtMtI*RUOQsm0B8{a{52b;U&_u(o+ww3$LSBja7AnjA& z7;4*7_R;7?ummO4hFV^AIE6-_sJ^VxB4S0tw$PgTpwO`$teiQ*DmBjmGEfl=`T=f$ zn&fDZ*lbQLpz~^DDS0Snn12G&xQA}TA2$Cy3uE_qdp+xXTX(-)>HJ2^1l_dD;|~qb z@2w2852$OljIZp``ONb$fnN9+)Tq7}ez>}Tt{X6~^SzRyb^6QyVdhRLfBFG+yVVodw6KU3BL6!p4u9dwv1D+P&GW&tEYH%}-lq zI!rVOr0E=J9~IKtb*5aJcDK@bS1<2V4bcsr2*3Fv$qxT+K&%cG%XoYHoOy4(!^?ln zxy9Id`GeE0MOml@WV!vNfZriq%J#~|k7duUrZAJIXp8YkQ=~Z;EkJDHEowF{V_i=$4lnl-tex z2s5?T(m2Xrj^y&iy@r3E$ArQXdrGHgpy>G>Z8DZLeG(Tr7(G)2{8S7=uYZ-9uQA~f zH1xCkq3Nt?-ciLRVlJU}?Y)BP&s1y!oJ3%=TuE)dYQ-X(T$#-q_Y_hz+%45;7T4Tm ze56QoZDe{uh3>hnVjYY%lDZZPtjye;B#%#BMb078;qyi+I@+JRKV_{Fd~3Lr8?&e! znlo{-Cg}AR<##C2UX6#XR~~3IGSR?8*O~7ZLWz>JeDf0mit{=^>42B(Drc6rkA1Tz z4xx<^WTcT0EexDExkWxBquZ%Y5tI6H(>gUG(*_JR$dQR1_>zHWPsibHqTYC!$_r}r zUOPyF}g~Z`mtMH*4pF{@Wi6&eB2P>uXh5i1_impXL{M*h=H|EGzp7nA2;gNS za$|AGNIPI2Ej2JILmz7>nsqq2&>(C>cKcam@8md8`0h{i3dyF-M3qW6F9kST9)T!> z))D9@8rEEr>2#^F{YH_rOWgqPH`D&9drJUSF%s8GSIjSpIf;|TDQ=P+i@z>qSCp<+ zCWmJf1#nv-vkSm^pJX6x_*#1dQ#FWCg*(_%l(V=wL_V;HSX~`Kfi}L$rmo#dhe^W7 zp(N~_+~ksW`{8fKuVf%^7xzv31oFm>pxKV2^GjzTdy~9TP9DErM>C8P1jdNBfcAJ0 zK{L2VpE!?h1Eif30wgL>VC|~e!1Mh>CinTTQm z35xlu(V06e{Rv}KHb^}P0cK-!QoELZWwVC;&9G}xl}Y$kTLGB523kDy@_H8vJ08;ee}bAM%BFv;X$( zLWTM#^qTDEnGI&Q=*B0~b>#u?CW}jN5=zbQsmC$FOX5XiQkX9CARWQTVa}f;n6iHX zNleEB(s)u$%0e=2v}E$6+<5_;Af+uL1!b>HalN>X7xtm|ZzH0(H<3qbDe)Bmh7u@A zf4KgDp2F97snd8~(v}wtH1afuv_a^9UQ~nB*%pX`tde)}Ho>kIGd~N-D;`TXm{hxt zb_S!XeX56&G){w~A|o=~)j5qk5+aX_OziQ5TZ+7}za;Mab4C))A=$?GkmOa$k~- Date: Tue, 24 Aug 2021 15:04:07 +0200 Subject: [PATCH 2/3] fixes bugs of whitespaces in the ins file and removes unncessary outcommented parameters --- rLPJGUESS/R/AdjustTemplates.R | 19 ++-- rLPJGUESS/R/GetRunAbleParameters.R | 21 ++-- rLPJGUESS/R/InferParameterAndDesignList.R | 111 ++++++++++++++-------- rLPJGUESS/R/runLPJ.R | 4 +- 4 files changed, 97 insertions(+), 58 deletions(-) diff --git a/rLPJGUESS/R/AdjustTemplates.R b/rLPJGUESS/R/AdjustTemplates.R index 8f9e0ed..d4e9918 100644 --- a/rLPJGUESS/R/AdjustTemplates.R +++ b/rLPJGUESS/R/AdjustTemplates.R @@ -2,13 +2,12 @@ #' @description This function takes the defaultparameters and a defaultlist of #' of files and writes the main and pft instruction file into the defined destinations #' which are MainTemplateDestination and PftTemplateDestination -#' @param defaultparameters a matrix with the defaultparameters which should be -#' used in rlpjguess runs +#' @param defaultobjects R object produced by \code{\link{InferParameterAndDesignList}}. Its first list contains the data frame with default parameters extracted from the instruction files #' @param defaultlist a matrix with the design parameters which should be used #' in rlpjguess runs -#' @param MainTemplateDestination a character string indicating the space where +#' @param MainTemplateDestination a character string indicating the location to which #' one should write the main file -#' @param PftTemplateDestination a character string indicating the space where +#' @param PftTemplateDestination a character string indicating the location to which #' one should write the pft file #' @param NameMainFile a character string with the name of the main file #' which was produced with the default parameters before @@ -19,12 +18,14 @@ #' @export -AdjustTemplates <- function(defaultparameters,defaultlist, +AdjustTemplates <- function(defaultobjects, MainTemplateDestination = NULL, PftTemplateDestination = NULL, NameMainFile = NULL, NamePftFile = NULL){ + + # check if there is a destination to save the main template if(is.null(MainTemplateDestination)){ @@ -45,9 +46,13 @@ AdjustTemplates <- function(defaultparameters,defaultlist, stop("Please provide a valid pft template") } - ## TO-DO: Dokumentation - ## TO-DO: Maybe give another return than the function does now + ## define the thinks used here + + defaultparameters = defaultobjects$defaultparameters + defaultlist = defaultobjects$defaultlist + + ## read in text and adjust the templates tx <- readLines(NamePftFile) diff --git a/rLPJGUESS/R/GetRunAbleParameters.R b/rLPJGUESS/R/GetRunAbleParameters.R index cfce3bc..3f1b32b 100644 --- a/rLPJGUESS/R/GetRunAbleParameters.R +++ b/rLPJGUESS/R/GetRunAbleParameters.R @@ -1,28 +1,28 @@ -#' @title A function to separate all default values into parameters, design parameters and input files -#' @description This function organizes the parameters extracted from the instruction files in order to run the model -#' @param defaultparameters the default parameters extracted from the instruction files -#' @param PFTs the species one wants to run the model with +#' @title A function to separate all default values into parameters, design parameters and input files. +#' @description This function organizes the parameters extracted from the instruction files in order to run the model. +#' @param defaultobjects R object produced by \code{\link{InferParameterAndDesignList}}. Its first list contains the data frame with default parameters extracted from the instruction files +#' @param PFTs the PFTs or species, which should be included into the model run #' @keywords rLPJGUESS #' @author Johannes Oberpriller -#' @return list of runable parameters, design parameters and the files containing the default inputs +#' @return list consisting of a data.frame for runable parameters, a data.frame for design parameters and a data.frame for the default input files. #' @export -GetRunAbleParameters <- function(defaultparameters,PFTs){ +GetRunAbleParameters <- function(defaultobjects,PFTs){ - LPJparameters_PFT <- matrix(defaultparameters$defaultparameters[,4]) - rownames(LPJparameters_PFT) = defaultparameters$defaultparameters[,3] + LPJparameters_PFT <- matrix(defaultobjects$defaultparameters[,4]) + rownames(LPJparameters_PFT) = defaultobjects$defaultparameters[,3] LPJrunParameters = matrix(LPJparameters_PFT[-which(substr(rownames(LPJparameters_PFT),1,3)== "run"),]) rownames(LPJrunParameters) = rownames(LPJparameters_PFT)[-which(substr(rownames(LPJparameters_PFT),1,3)== "run")] designLPJ = LPJparameters_PFT[which(substr(rownames(LPJparameters_PFT),1,3)== "run"),] designLPJ = designLPJ[-which(names(designLPJ) == "run_outputdirectory")] - lpjvalues = as.matrix(defaultparameters$defaultlist[,c(2,3)]) + lpjvalues = as.matrix(defaultobjects$defaultlist[,c(2,3)]) filelist = split(t(lpjvalues),rep(1:nrow(lpjvalues),each = ncol(lpjvalues))) - names(filelist) = gsub("_", ".",defaultparameters$defaultlist[,1]) + names(filelist) = gsub("_", ".",defaultobjects$defaultlist[,1]) spp <- as.matrix(LPJrunParameters[grep("_include", rownames(LPJrunParameters )),]) PFTsRows <- c() @@ -30,6 +30,7 @@ GetRunAbleParameters <- function(defaultparameters,PFTs){ spp[-PFTsRows, 1] <- 0 spp[PFTsRows,1] <- 1 print(spp) + parameterStandard_PFT <- as.matrix(LPJrunParameters[-grep("_include", rownames(LPJrunParameters)),]) diff --git a/rLPJGUESS/R/InferParameterAndDesignList.R b/rLPJGUESS/R/InferParameterAndDesignList.R index 2dc3d62..7e68b40 100644 --- a/rLPJGUESS/R/InferParameterAndDesignList.R +++ b/rLPJGUESS/R/InferParameterAndDesignList.R @@ -3,20 +3,21 @@ #' @description This function returns a list of parameters and design parameters #' used for the later model runs. It also allows the user to see a summary #' of all his parameters -#' @param MainInputFile a list with one element named "main" which is the main -#' instruction file given to lpj-guess when running it in the command line -#' @param NameMainFile a character string with the name of the main file. In +#' @param MainInputFile a list with one element named "main" which contains the name of the main +#' instruction file given to lpj-guess when running it on the command line +#' @param NameMainFile a character string with the name of the new main file. Into #' this file this function will write all input parameters. When no value is -#' provided the default main instruction file while be main.ins -#' @param NamePftFile a character string with the name of the pft file. In this +#' provided the default main instruction file name while be main.ins +#' @param NamePftFile a character string with the name of the pft file. Into this #' file the function will write all parameters specific for pfts and design. -#' When no value is provided the default main instruction file while be pft.ins +#' When no file is provided the default main instruction file name will be pft.ins #' @param vectorvaluedparams a vector of character strings providing the vector -#' valued parameters. rlpjguess can in the moment not handel vectorvalued -#' parameters, which can not be computed from the first element. Typically -#' for this is the rootdist, which depends on the different layers in the soil. -#' @return a list with a matrix for the defaultparameter values and a matrix -#' with defaultdesign values +#' valued parameters. rlpjguess can in the moment not handle vectorvalued +#' parameters. An example for this type of parameters is the "rootdist"-parameters, +#' which has for each soil layer an individual parameter. +#' @return a list a containing matrix with the default parameter values (called defaultparameters) +#' read from the instruction file and a matrix +#' with default design values (called defaultdesign) read from the instruction file #' @author Johannes Oberpriller #' @export @@ -72,6 +73,11 @@ InferParameterAndDesignList <- function(MainInputFile,NameMainFile = "main.ins", for(j in 1:(endofgroup-1)){ parametersnew = vector(length = 4) parametersnew[1] = "pft" + if(nchar(trimws(runlist_pfts[i+j])) != 0){ + while(nchar(trimws(substring(runlist_pfts[i+j],1,1)))==0){ + runlist_pfts[i+j] = substring(runlist_pfts[i+j],2,nchar(runlist_pfts[i+j])) + } + } parametersnew[2] = gsub("[[:space:]]","",strsplit(runlist_pfts[i+j]," ")[[1]][1]) parametersnew[3] = paste0(gsub("\"","",linename),"_", parametersnew[2]) parameterposition = 2 @@ -87,7 +93,9 @@ InferParameterAndDesignList <- function(MainInputFile,NameMainFile = "main.ins", } parametersnew[4] = gsub("^!","",parametersnew[4]) parametersnew[5] = gsub("\"","",linename) - defaultparameters = rbind(defaultparameters, parametersnew) + if(substring(runlist_pfts[i+j],1,1) != "!"){ + defaultparameters = rbind(defaultparameters, parametersnew) + } } } } @@ -157,7 +165,7 @@ InferParameterAndDesignList <- function(MainInputFile,NameMainFile = "main.ins", # @description This function takes the main input input files and desired file # names and returns a list with files, variables, parameters and design parameters -# @param MainInputFile is the main input file see InferParameterAndDesignList +# @param MainInputFile is a list with the first element being the name of the main input file see InferParameterAndDesignList # @param NameMainFile is the main file after reordering see InferParameterAndDesignList # @param NamePftFile is the pft file after reordering see InferParameterAndDesignList # @return returns a list with the entries parameters,files, variables and design @@ -170,26 +178,39 @@ write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ insfilelist = infer_instruction_order(MainInputFile = MainInputFile) - # Create the temporary files to get all lines + # If the main input file imports information from other input files, create + # one large temporary .ins file - file.create("./main_tryout.ins") - outcon <- file("./main_tryout.ins", "w") - for(i in 1:length(insfilelist)){ - lines <- readLines(insfilelist[[i]]) - writeLines(lines,outcon) - } - close(outcon) + if(length(insfilelist) > 1){ + + # Create the temporary files to get all lines + + file.create("./main_tryout.ins") + outcon <- file("./main_tryout.ins", "w") + for(i in 1:length(insfilelist)){ + lines <- readLines(insfilelist[[i]]) + writeLines(lines,outcon) + } + close(outcon) + + # Create the Main and Ins file + + file.create(paste0("./",NameMainFile)) + file.create(paste0("./",NamePftFile)) + lines <- readLines("./main_tryout.ins") + + # delete all import files - # Create the Main and Ins file + importlines = grepl("import",substr(lines, start = 1, stop = 7)) # AHES I would increase stop, you never know whether the 'import' statement is somehwere further down. As we have learned today, the .ins files can be quite different. + lines = lines[-which(importlines)] - file.create(paste0("./",NameMainFile)) - file.create(paste0("./",NamePftFile)) - lines <- readLines("./main_tryout.ins") + }else{ - # delete all import files + # Else, if there is only one .ins file, read in all its lines. - importlines = grepl("import",substr(lines, start = 1, stop = 7)) - lines = lines[-which(importlines)] + lines <- readLines(MainInputFile[['main']]) + + } # Grep all "param" lines and the ndep_timeseries, delete duplicates # and write them to the parameterfiles @@ -200,21 +221,24 @@ write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ parameternames = sapply(strsplit(parameters," "), function (x) x[2]) duplicates_parameters = which(duplicated(parameternames) ==T) + while(length(duplicates_parameters) != 0 ){ - first_apperance = c() + first_appearance = c() for(i in duplicates_parameters){ - first_apperance <- c(first_apperance, which(parameternames == parameternames[i])[1]) + first_appearance <- c(first_appearance, which(parameternames == parameternames[i])[1]) } - parameters = parameters[-first_apperance] + parameters = parameters[-first_appearance] parameternames = sapply(strsplit(parameters," "), function (x) x[2]) duplicates_parameters = which(duplicated(parameternames) ==T) } - ndep_time = grepl("!ndep_timeseries", substr(lines, start = 1, stop = 16)) + + ndep_time = grepl("!ndep_timeseries", substr(lines, start = 1, stop = 15)) if(length(which(ndep_time ==T)) ==0){ ndep_time = grepl("ndep_timeseries", substr(lines, start =1, stop =15)) } parameters = c(parameters, lines[which(ndep_time)]) + parameternames = sapply(strsplit(parameters," "), function (x) x[2]) parameterfiles = grepl("file", parameternames) parametervariables = parameters[-which(parameterfiles)] @@ -275,11 +299,17 @@ write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ designnames = sapply(strsplit(lines," "), function (x) x[1]) duplicates_design = which(duplicated(designnames) ==T) - first_apperance = c() - for(i in duplicates_design){ - first_apperance <- c(first_apperance, which(designnames == designnames[i])[1]) + + if(length(duplicates_design) != 0){ + first_apperance = c() + for(i in duplicates_design){ + first_apperance <- c(first_apperance, which(designnames == designnames[i])[1]) + } + design = lines[-first_apperance] + }else{ + design = lines } - design = lines[-first_apperance] + # write the PFT or Species file which includes the title, the files, # the designvariables and all pft/speciesspecific lines into the @@ -289,9 +319,11 @@ write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ writeLines(c(titel,headingfileline,files,"\n",design, "\n",pftlines),Pftcon) close(Pftcon) + # If there were multiple .insfiles, which first had to be merged, a temporary merged file + # had been created. Delete this file now. # remove the main file written in the beginning because no more need - file.remove("./main_tryout.ins") + if(length(insfilelist) > 1) file.remove("./main_tryout.ins") return(list(files = parameterfiles, variables = parametervariables, pfts = pftlines, design = design )) @@ -301,7 +333,7 @@ write_main_templates <- function(MainInputFile,NameMainFile, NamePftFile){ # @description This function infers the correct order of the instruction files # and returns it # @param MainInputFile is the main input file see InferParameterAndDesignList -# @return returns a vector with the names of correct ordered instruction files +# @return returns a vector with the names of correctly ordered instruction files # @author Johannes Oberpriller infer_instruction_order <- function(MainInputFile){ @@ -317,7 +349,7 @@ infer_instruction_order <- function(MainInputFile){ # get all imports from the file - imports = get_imports(MainInputFile[["main"]],maindirectory = maindirectory) + imports = get_imports(MainInputFile[["main"]], maindirectory = maindirectory) # order the imports like they are in the main instruction file @@ -397,3 +429,4 @@ get_imports <- function(file,maindirectory){ } return(imports) } + diff --git a/rLPJGUESS/R/runLPJ.R b/rLPJGUESS/R/runLPJ.R index 3bc5536..75eb51e 100644 --- a/rLPJGUESS/R/runLPJ.R +++ b/rLPJGUESS/R/runLPJ.R @@ -6,8 +6,8 @@ #' character string indicating the path to the directory where #' the model link and template are located, and in which the function will create #' the directory structure for the outputs. -#' @param parameterList either a named list containing the parameters to be calibrated -#' or a matrix. If running in parallel, parameter list should be either a list of of list or +#' @param parameterList a named list containing the parameters (when not run in parallel) to be calibrated +#' or a matrix (when run in parallel). If running in parallel, parameter list should be either a list of of list or #' a matrix where each row is a parameter combination and the column names should be named #' after the parameters. See fucntion \code{\link{getParameterList}}) for default values. #' @param typeList a character vector with the outputs to be analyzed. From 07ecf104e04989193adcf6b30b6487cd8b413bad Mon Sep 17 00:00:00 2001 From: JohannesOberpriller <46445144+JohannesOberpriller@users.noreply.github.com> Date: Tue, 1 Feb 2022 09:25:30 +0100 Subject: [PATCH 3/3] Update README.md Update with actual workflow and placement of the files required to run the model --- README.md | 47 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 47 insertions(+) diff --git a/README.md b/README.md index ef6fc54..f104e09 100644 --- a/README.md +++ b/README.md @@ -54,3 +54,50 @@ In Rstudio, the vignette may not be built per default. You will turn this on in ```{r} devtools::build_vignettes("rLPJGUESS") ``` +### Workflow when using the package: + +#### 1. Define a folder to run the simulations, e.g. +```{r} +mainDir <- file.path(getwd(), "LPJrunTest") +``` +and place the binary (exectuable) of LPJ-GUESS into this folder (the binary does not come with this package). + +#### 2. Provide a valid instruction file that you usually use to define the set-up of LPJ-GUESS to the InferParameterAndDesignList function. +```{r} +defaultparameters <- InferParameterAndDesignList(list(main = paste0(mainDir,"LPJ_instruction_file.ins")), + NameMainFile = "main.ins", NamePftFile = "pft.ins", + vectorvaluedparams = c("rootdist","eps_mon", + "storfrac_mon","photo", + "fertdates","fertrate")) +``` + +This function extracts the parameters, design and driver-files used to run LPJ-GUESS and splits the instruction into a main file containing inputs and a PFT file, which has all the settings. Because we cannot handle vectorvalued parameters at the moment, they have to be provided to this function as well. + +#### 3. Call the adjust template function, which rewrites the templates to new destinations (IMPORTANT: the files need to be in the same directory as the binary, here: LPJrunTest) +```{r} +AdjustTemplates(defaultparameters = defaultparameters$defaultparameters, + defaultlist = defaultparameters$defaultlist, + MainTemplateDestination = "./LPJrunTest/main_new.ins", + PftTemplateDestination = "./LPJrunTest/pft_new.ins", + NameMainFile = "main.ins", NamePftFile = "pft.ins") +``` + +#### 4. Generate a set of parameters for the species (GetRunAbleParameters), we want to simulate, here Fagus sylvatica (Fag_syl): +```{r} +parameters <- GetRunAbleParameters(defaultparameters = defaultparameters, PFTs = c("Fag_syl")) +``` +#### 5. Change some of these parameters and provide it as matrix with a set of parameters per row and parameter names as colnames (here: new parameter matrix is called parameters_new) + +#### 6. Define the settings and input files, we want to use as input for the simulations in a list (the required files depend on the LPJ-Version, here: the list is called LPJsettings) + +#### 7. Define the setup (i.e. sequential, parallel, which kind of parallelisation, here a SOCK cluster with numCores cores) +```{r} +LPJsetup <- setupLPJParallel(numCores = numCores, clusterType = "SOCK", + mainDir = mainDir) +``` + +#### 8. Run the simulations with the model +```{r} +results <- runLPJ(x = LPJsetup, parameterList = parameters_new, + typeList = typeList, settings = LPJsettings) +```