big changes

RBBGCMuso/DESCRIPTION

@@ -22,6 +22,7 @@ Imports:
 	tibble,
 	tidyr
 LinkingTo: Rcpp
+SystemRequirements: C++11
 Maintainer: Roland Hollo's <hollorol@gmail.com>
 RoxygenNote: 6.1.0
 Suggests: knitr,

RBBGCMuso/DESCRIPTION~

@@ -0,0 +1,29 @@
+Package: RBBGCMuso
+Title: An R package for BiomeBGC-MuSo ecosystem modelling
+Version: 0.5.0.0-0
+Authors@R: person("Roland", "Hollo's", , "hollorol@gmail.com", role = c("aut", "cre"))
+Description: What the package does (one paragraph).
+Depends: R (>= 2.10)
+License: GPL-2
+LazyData: true
+NeedsCompilation: no
+Packaged: 2017-07-19 14:00:04 UTC; hollorol
+Author: Roland Hollo's [aut, cre]
+Imports:
+	grDevices,
+	stats,
+	utils,
+	graphics,
+	Rcpp,
+	magrittr,
+	dplyr,
+	ggplot2,
+	rmarkdown,
+	tibble,
+	tidyr
+LinkingTo: Rcpp
+Maintainer: Roland Hollo's <hollorol@gmail.com>
+RoxygenNote: 6.1.0
+Suggests: knitr,
+    rmarkdown,
+VignetteBuilder: knitr

RBBGCMuso/NAMESPACE

@@ -2,15 +2,18 @@
 
 export(OtableMaker)
 export(calibMuso)
+export(changeMusoC)
 export(cleanupMuso)
-export(copyMusoExamleTo)
+export(compareMuso)
 export(corrigMuso)
 export(getyearlycum)
 export(getyearlymax)
 export(musoDate)
 export(musoMapping)
 export(musoMappingFind)
+export(musoMont)
 export(musoMonte)
+export(musoRand)
 export(musoRandomizer)
 export(musoSensi)
 export(normalMuso)

RBBGCMuso/R/#musoFilter.R#

@@ -0,0 +1,4 @@
+musoFilter <- function(dB=.,text){
+    eval(parse(text = paste0("filter(",dB,",.,",text,")"))) %>%
+        tbl_df
+}

RBBGCMuso/R/#musoRand.R#

@@ -0,0 +1,172 @@
+#' musoRand
+#'
+#' This funtion samples uniformly from the choosen parameters of the BiomeBGC-Muso model, which parameters are constrained by the model logic.
+#' @author Roland Hollos
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the input-file, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param constrains This is a matrics wich specify the constrain rules for the sampling. Further informations coming son.
+#' @param iteration The number of sample-s. It is adviced to use at least 3000 iteration, because it is generally fast and it can be subsampled later at any time. 
+#' @export
+
+musoRand <- function(parameters, constrains = NULL, iterations=3000){
+
+    if(!is.null(constrains)){
+        constMatrix <- constrains
+    }
+    
+    parameters <- parameters[,-1]
+    constMatrix <- constMatrix[,-1]
+    
+    depTableMaker <- function(constMatrix,parameters){
+	parameters <- parameters[order(parameters[,1]),]
+	constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],c(5,6)]<-parameters[,c(2,3)]
+	logiConstrain <- (constMatrix[,"GROUP"] %in% constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],"GROUP"] &
+			  (constMatrix[,"GROUP"]!=0)) | ((constMatrix[,"INDEX"] %in% parameters[,1]) & (constMatrix[,"GROUP"] == 0))
+	constMatrix<-constMatrix[logiConstrain,]
+	constMatrix <- constMatrix[order(apply(constMatrix[,7:8],1,function(x){x[1]/10+abs(x[2])})),]
+	constMatrix
+    }
+
+    genMat0 <- function(dep){
+	numberOfVariable <- nrow(dep)
+	G <- rbind(diag(numberOfVariable), -1*diag(numberOfVariable))
+	h <- c(dependences[,5], -1*dependences[,6])
+	return(list(G=G,h=h))
+    }
+
+    genMat1 <- function(dep, N){
+
+	## Range <- sapply(list(min,max),function(x){
+  
+	##   x(as.numeric(rownames(dep)))
+	## }) It is more elegant, more general, but slower
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	numberOfVariables <- nrow(dep)
+	G<- -1*diag(numberOfVariables)
+
+	for(i in 1:numberOfVariables){
+	    if(dep[i,4]!=0){
+		G[i,dep[i,4]] <- 1
+	    }
+
+	}
+
+	G<-G[dep[,4]!=0,]
+
+	if(Range[1]==1){
+	    G<-cbind(G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	} else{
+	    if(Range[2]==N){
+		G<-cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G)
+	    } else {
+		G <- cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	    }
+	}
+	return(list(G=G,h=rep(0,nrow(G))))
+    }
+
+    genMat2 <- function(dep, N){
+	G <- rep(1,nrow(dep))
+
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	if(Range[1]==1){
+	    G<-c(G, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		G<-c(numeric(Range[1]-1), G)
+	    } else {
+		G <- c(numeric(Range[1]-1), G, numeric(N-Range[2]))
+	    }
+	}
+
+	G <- t(matrix(sign(dep[2,4])*G))
+	h <- abs(dep[1,4])
+
+	return(list(G=G,h=h))
+    }
+
+    genMat3 <- function(dep, N){
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	E <- rep(1,nrow(dep))
+
+	if(Range[1]==1){
+	    E<-c(E, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		E<-c(numeric(Range[1]-1), E)
+	    } else {
+		E <- c(numeric(Range[1]-1), E, numeric(N-Range[2]))
+	    }
+	}
+
+
+	E <- t(matrix(E))
+	f <- dep[1,4]
+	return(list(E=E,f=f))
+    }
+
+
+    applyRandTypeG <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(type,
+	       invisible(Gh <- genMat1(dep, N)),
+	       invisible(Gh <- genMat2(dep, N)))
+	return(Gh)
+    }
+
+    applyRandTypeE <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(-type,
+	       stop("Not implemented yet"),
+	       stop("Not implemented yet"),
+	       invisible(Ef <- genMat3(dep, N)))
+	return(Ef)
+    }
+
+    dependences <- depTableMaker(constMatrix, parameters)
+    dependences <- cbind(dependences,1:nrow(dependences))
+    colnames(dependences)[ncol(dependences)] <- "rowIndex"
+    numberOfVariable <- nrow(dependences)
+    nonZeroDeps<-dependences[dependences[,"TYPE"]!=0,]
+    if(nrow(nonZeroDeps)!=0){
+        splitedDeps<- split(nonZeroDeps,nonZeroDeps[,"GROUP"])
+        Gh <- list()
+        Ef <- list()
+
+        for(i in 1:length(splitedDeps)){
+            print(splitedDeps[[i]][1,"TYPE"])
+            if(splitedDeps[[i]][1,"TYPE"]>0){
+                Gh[[i]]<-applyRandTypeG(splitedDeps[[i]],nrow(dependences))
+            } else {
+                Ef[[i]] <- applyRandTypeE(splitedDeps[[i]],nrow(dependences))
+            }
+        }
+        
+        Gh0<- genMat0(dependences)
+        G <- do.call(rbind,lapply(Gh,function(x){x$G}))
+        G<- rbind(Gh0$G,G)
+        h <- do.call(c,lapply(Gh,function(x){x$h}))
+        h <- c(Gh0$h,h)
+        E <- do.call(rbind,lapply(Ef,function(x){x$E}))
+        f <- do.call(c,lapply(Ef,function(x){x$f}))
+        randVal <- suppressWarnings(limSolve::xsample(G=G,H=h,E=E,F=f,iter = iterations))$X
+    } else{
+        Gh0<-genMat0(dependences)
+        randVal <- suppressWarnings(limSolve::xsample(G=Gh0$G,H=Gh0$h, iter = iterations))$X
+    }
+    
+    results <- list(INDEX =dependences$INDEX, randVal=randVal)
+    return(results)
+}

RBBGCMuso/R/.#musoFilter.R

@@ -0,0 +1 @@
+hollorol@meteor22.13528:1542837570

RBBGCMuso/R/.#musoRand.R

@@ -0,0 +1 @@
+hollorol@meteor22.13528:1542837570

RBBGCMuso/R/RcppExports.R

@@ -1,10 +1,30 @@
 # Generated by using Rcpp::compileAttributes() -> do not edit by hand
 # Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
 
+getWritePositions <- function(a) {
+    .Call('_RBBGCMuso_getWritePositions', PACKAGE = 'RBBGCMuso', a)
+}
+
+#' changeMusoC
+#' 
+#' This function is fastly randomize values based on min and max values
+#' @importFrom Rcpp evalCpp
+#' @useDynLib RBBGCMuso
+#' @param inFile is the big matrix
+#' @param outFile is the small matrix
+#' @export
+changeMusoC <- function(inFile, outFile, inMat) {
+    invisible(.Call('_RBBGCMuso_changeMusoC', PACKAGE = 'RBBGCMuso', inFile, outFile, inMat))
+}
+
 randTypeOne <- function(m) {
     .Call('_RBBGCMuso_randTypeOne', PACKAGE = 'RBBGCMuso', m)
 }
 
+randTypeTwo <- function(m) {
+    .Call('_RBBGCMuso_randTypeTwo', PACKAGE = 'RBBGCMuso', m)
+}
+
 #' musoRandomizer
 #' 
 #' This function is fastly randomize values based on min and max values,

RBBGCMuso/R/calibMuso.R

@@ -23,9 +23,14 @@
 #' @import utils
 #' @export
 
-calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE, logfilename=NULL, keepEpc=FALSE, export=FALSE, silent=FALSE, aggressive=FALSE, leapYear=FALSE,keepBinary=FALSE, binaryPlace="./", fileToChange="epc", skipSpinup = TRUE){
+calibMuso <- function(settings=NULL, calibrationPar=NULL,
-
+                      parameters=NULL, outVars = NULL, timee="d",
-
+                      debugging=FALSE, logfilename=NULL,
+                      keepEpc=FALSE, export=FALSE,
+                      silent=FALSE, aggressive=FALSE,
+                      leapYear=FALSE,keepBinary=FALSE,
+                      binaryPlace="./", fileToChange="epc",
+                      skipSpinup = TRUE, modifyOriginal =FALSE){
 ##########################################################################
 ###########################Set local variables and places########################
 ########################################################################
@@ -90,25 +95,46 @@ calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE,
          cleanupMuso(location=outputLoc,deep = TRUE)
      }
 
+    toModif<-c(epc[2],iniInput[2])
 
+    if(!modifyOriginal & (!is.null(parameters) | !is.null(outVars)))
+    {
+        
+        toModif <- sapply(toModif, function (x){
+            paste0(tools::file_path_sans_ext(basename(x)),"-tmp.",tools::file_ext(x))
+        })
+        
+    }
     
     ##change the epc file if and only if there are given parameters
     if(!is.null(parameters)){
-
+        changemulline(filePaths=c(epc[2],iniInput[2]), calibrationPar = calibrationPar,
-        switch(fileToChange,
+                                 contents = parameters, fileOut=toModif, fileToChange=fileToChange, modifyOriginal=modifyOriginal)
-               "epc"=tryCatch(changemulline(filename=epc[2],calibrationPar,parameters),
-                              error= function (e){
-                                  setwd(whereAmI)
-                                  stop("Cannot change the epc file")}),
-               "ini"=tryCatch(changemulline(filename=iniInput[2],calibrationPar,parameters),
-                              error= function (e){
-                                  setwd((whereAmI))
-                                  stop("Cannot change the ini file")}),
-               "both"=(stop("This option is not implemented yet, please choose epc or ini"))
-               )
     }
     
+
     ##We change the working directory becase of the model, but we want to avoid sideeffects, so we save the current location and after that we will change everything to it.
+    if(!modifyOriginal & (!is.null(parameters) | !is.null(outVars))){
+    epc[2]<-file.path(dirname(epc[2]),toModif[1]) # Writing back the lost path
+    toModif[2]<-file.path(dirname(iniInput[2]),toModif[2]) #for the Initmp, also
+    if((!is.null(outVars) | !file.exists(toModif[2])) & !modifyOriginal){
+        file.copy(iniInput[2],toModif[2],overwrite = TRUE)
+    }
+
+    iniInput[2] <- toModif[2]}
+
+    if(!is.null(parameters) & ((fileToChange == "epc") | (fileToChange == "both")) & !modifyOriginal){
+        tmp<-readLines(iniInput[2])
+        tmpInd<-grep("EPC_FILE",tmp)+1
+        tmp[tmpInd]<-file.path(dirname(tmp[tmpInd]),basename(epc[2]))
+        writeLines(tmp,iniInput[2])
+        rm(list=c("tmp","tmpInd"))
+    }
+    if(!is.null(outVars)){
+        outputVarChanges <- putOutVars(iniInput[2], outputVars = outVars, modifyOriginal = !modifyOriginal)
+        settings$outputVars[[1]]<-outputVarChanges[[1]]
+        settings$numData <- round(settings$numData*outputVarChanges[[2]])
+    }
     
    if(!skipSpinup) {
 

RBBGCMuso/R/calibMuso.R~

@@ -0,0 +1,348 @@
+#' calibMuso 
+#'
+#' This function changes the epc file and after that  runs the BBGC-MuSo model and reads in its outputfile in a well-structured way.
+#' 
+#' @author Roland Hollos
+#' @param settings You have to run the setupMuso function before calibMuso. It is its output which contains all of the necessary system variables. It sets the whole running environment
+#' @param timee The required timesteps in the modell output. It can be "d", if it is daily, "m", if it's monthly, "y", it it is yearly. I recommend to use daily data, the yearly and monthly data is not well-tested yet. 
+#' @param debugging If it is TRUE, it copies the log file to a Log directory to store it, if it is stamplog it contatenate a number before the logfile, which is one more than the maximum of the represented ones in the LOG directory. If it is true or stamplog it collects the "wrong" logfiles
+#' @param keepEpc If TRUE, it keeps the epc file and stamp it, after these copies it to the EPCS directory. If debugging True or false, it copies the wrong epc files to the wrong epc directory.
+#' @param export if it is yes or you give a filename here, it converts the ouxtput to the specific extension. For example, if you set export to "example.csv", it converts the output to "csv", if you set it to "example.xls" it converts to example.xls with the xlsx package. If it is not installed it gives back a warning message and converts it to csv.
+#' @param silent If you set it TRUE all off the modells output to the screen will be suppressed. It can be usefull, because it increases the model-speed.
+#' @param aggressive It deletes every possible modell-outputs from the previous modell runs.
+#' @param parameters In the settings variable you have set the row indexes of the variables, you wish to change. In this parameter you can give an exact value for them in a vector like: c(1,2,3,4)
+#' @param logfilename If you want to set a specific name for your logfiles you can set this via logfile parameter
+#' @param leapYear  Should the function do a leapyear correction on the outputdata? If TRUE, then the 31.12 day will be doubled.
+#' @param keepBinary In default RBBGCMuso to keep  working area as clean as possible, deletes all the regular output files. The results are directly printed to the standard output, but you can redirect it, and save it to a variable, or you can export your results to the desired destination in a desired format. Whith this variable you can enable to keep the binary output files. If you want to set the location of the binary output, please take a look at the binaryPlace argument.
+#' @param binaryPlace The place of the binary output files.
+#' @param fileToChange You can change any line of the epc or the ini file, you just have to specify with this variable which file you van a change. Two options possible: "epc", "ini"
+#' @param skipSpinup If TRUE, calibMuso wont do spinup simulation
+#' @return No return, outputs are written to file 
+#' @usage calibMuso(settings,parameters=NULL, timee="d", debugging=FALSE, logfilename=NULL,
+#' keepEpc=FALSE, export=FALSE, silent=FALSE, aggressive=FALSE, leapYear=FALSE)
+#' @import utils
+#' @export
+
+calibMuso <- function(settings=NULL, calibrationPar=NULL,
+                      parameters=NULL, timee="d",
+                      debugging=FALSE, logfilename=NULL,
+                      keepEpc=FALSE, export=FALSE,
+                      silent=FALSE, aggressive=FALSE,
+                      leapYear=FALSE,keepBinary=FALSE,
+                      binaryPlace="./", fileToChange="epc",
+                      skipSpinup = TRUE, modifyOriginal =FALSE){
+##########################################################################
+###########################Set local variables and places########################
+########################################################################
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    Linuxp <-(Sys.info()[1]=="Linux")
+    ##Copy the variables from settings
+    inputLoc <- settings$inputLoc
+    outputLoc <- settings$outputLoc
+    outputNames <- settings$outputNames
+    executable <- settings$executable
+    iniInput <- settings$iniInput
+    epc <- settings$epcInput
+    calibrationPar <- settings$calibrationPar
+    binaryPlace <- normalizePath(binaryPlace)
+    whereAmI<-getwd()
+
+
+    ## Set the working directory to the inputLoc temporarly.
+    setwd(inputLoc)
+
+      
+    if(debugging){#If debugging option turned on
+        #If log or ERROR directory does not exists create it!
+        dirName<-file.path(inputLoc,"LOG")
+        dirERROR<-file.path(inputLoc,"ERROR")
+        
+        if(!dir.exists(dirName)){
+            dir.create(dirName)
+        }
+
+        if(!dir.exists(dirERROR)){
+            dir.create(dirERROR)
+        }
+    }
+    
+    if(keepEpc) {
+        epcdir <- dirname(epc[1])
+        print(epcdir)
+        WRONGEPC<-file.path(inputLoc,"WRONGEPC")
+        EPCS<-file.path(inputLoc,"EPCS")
+        
+        if(!dir.exists(WRONGEPC)){
+            dir.create(WRONGEPC)
+        }
+        
+        if(!dir.exists(EPCS)){
+            dir.create(EPCS)
+        }
+    }   
+
+#############################################################
+############################spinup run############################
+   ########################################################## 
+
+    
+    
+
+     if(aggressive==TRUE){
+         cleanupMuso(location=outputLoc,deep = TRUE)
+     }
+
+    toModif<-c(epc[2],iniInput[2])
+
+    if(!modifyOriginal & !is.null(parameters))
+    {
+        
+        toModif <- sapply(toModif, function (x){
+            paste0(tools::file_path_sans_ext(basename(x)),"-tmp.",tools::file_ext(x))
+        })
+        
+    }
+    
+    ##change the epc file if and only if there are given parameters
+    if(!is.null(parameters)){
+        changemulline(filePaths=c(epc[2],iniInput[2]), calibrationPar = calibrationPar,
+                                 contents = parameters, fileOut=toModif, fileToChange=fileToChange, modifyOriginal=modifyOriginal)
+    }
+    
+
+    ##We change the working directory becase of the model, but we want to avoid sideeffects, so we save the current location and after that we will change everything to it.
+    if(!modifyOriginal & !is.null(parameters)){
+    epc[2]<-file.path(dirname(epc[2]),toModif[1]) # Writing back the lost path
+    toModif[2]<-file.path(dirname(iniInput[2]),toModif[2]) #for the Initmp, also
+    if(!file.exists(toModif[2]) & !modifyOriginal){
+        file.copy(iniInput[2],toModif[2],overwrite = TRUE)
+    }
+
+    iniInput[2] <- toModif[2]}
+
+    if(!is.null(parameters) & ((fileToChange == "epc") | (fileToChange == "both")) & !modifyOriginal){
+        tmp<-readLines(iniInput[2])
+        tmpInd<-grep("EPC_FILE",tmp)+1
+        tmp[tmpInd]<-file.path(dirname(tmp[tmpInd]),basename(epc[2]))
+        writeLines(tmp,iniInput[2])
+        rm(list=c("tmp","tmpInd"))
+    }
+    
+   if(!skipSpinup) {
+
+    ##Run the model for the spinup run.
+
+    if(silent){#silenc mode
+        if(Linuxp){
+            #In this case, in linux machines
+            tryCatch(system(paste(executable,iniInput[1],"> /dev/null",sep=" ")),
+                     error= function (e){
+                         setwd((whereAmI))
+                         stop("Cannot run the modell-check the executable!")})
+        } else {
+            #In windows machines there is a show.output.on.console option
+            tryCatch(system(paste(executable,iniInput[1],sep=" "),show.output.on.console = FALSE),
+                     error= function (e){
+                         setwd((whereAmI))
+                         stop("Cannot run the modell-check the executable!")})
+        }
+        
+    } else {
+        system(paste(executable,iniInput[1],sep=" "))
+    }
+
+    
+    logspinup <- getLogs(outputLoc,outputNames,type="spinup")
+    ## logspinup <- grep(paste0(outputNames[1],".log"), list.files(outputLoc),value = TRUE)
+    ## logspinup <- list.files(outputLoc)[grep("log$",list.files(outputLoc))]#load the logfiles
+    
+    if(length(logspinup)==0){
+        if(keepEpc){
+            stampnum<-stamp(EPCS)
+            lapply(epc,function (x) file.copy(from = x ,to=paste(EPCS,"/",(stampnum+1),"-", basename(x),sep="")))
+            lapply(epc, function (x) file.copy(from = paste(EPCS,"/",(stampnum+1),"-",basename(x),sep=""), to=WRONGEPC))
+            setwd(whereAmI)
+            stop("Modell Failure")
+        }
+        setwd(whereAmI)
+        stop("Modell Failure") #in that case the modell did not create even a logfile
+    }
+
+    if(length(logspinup)>1){
+        spincrash<-TRUE
+    } else {
+        if(identical(tail(readLines(paste(outputLoc,logspinup,sep="/"),-1),1),character(0))){
+            spincrash<-TRUE
+        } else {
+            spincrash <- (tail(readLines(paste(outputLoc,logspinup,sep="/"),-1),1)!=1)
+        }
+    }
+    } else {spincrash <- FALSE}
+                                        #If the last line in the logfile is 0 There are mistakes so the spinup crashes
+    
+    if(!spincrash){##If spinup did not crashed, run the normal run.
+        
+       ##################################################################### 
+       ###########################normal run#########################
+      #################################################################
+
+        ##for the sake of safe we set the location again
+        setwd(inputLoc)
+
+        if(silent){
+            if(Linuxp){
+                tryCatch(system(paste(executable,iniInput[2],"> /dev/null",sep=" ")),
+                         error =function (e){
+                             setwd((whereAmI))
+                             stop("Cannot run the modell-check the executable!")})
+            } else {
+                tryCatch(system(paste(executable,iniInput[2],sep=" "),show.output.on.console = FALSE),
+                         error =function (e){
+                             setwd((whereAmI))
+                             stop("Cannot run the modell-check the executable!")} )
+            }
+            
+        } else {
+            tryCatch(system(paste(executable,iniInput[2],sep=" ")),
+                     error =function (e){
+                         setwd((whereAmI))
+                         stop("Cannot run the modell-check the executable!")})
+        }
+
+
+        ##read the output
+         
+        switch(timee,
+               "d"=(Reva <- tryCatch(getdailyout(settings), #(:INSIDE: getOutput.R )
+                                    error = function (e){
+                                        setwd((whereAmI))
+                                        stop("Cannot read binary output, please check if the output type is set 2 in the ini files!")})),
+               "m"=(Reva <- tryCatch(getmonthlyout(settings), #(:INSIDE: getOutput.R )
+                                    error = function (e){
+                                        setwd((whereAmI))
+                                        stop("Cannot read binary output, please check if the output type is set 2 in the ini files!")})),
+               "y"=(Reva <- tryCatch(getyearlyout(settings), #(:INSIDE: getOutput.R )
+                                    error = function (e){
+                                        setwd((whereAmI))
+                                        stop("Cannot read binary output, please check if the output type is set 2 in the ini files!")}))
+               )
+        
+        if(keepBinary){
+            possibleNames <- tryCatch(getOutFiles(outputLoc = outputLoc,outputNames = outputNames),
+                                     error=function (e){
+                                         setwd((whereAmI))
+                                         stop("Cannot find output files")})
+            stampAndDir(outputLoc = outputLoc,names = possibleNames,stampDir=binaryPlace,type="output")
+        }
+    }
+
+    
+    if(skipSpinup){
+       logfiles <- tryCatch(getLogs(outputLoc,outputNames,type="normal"),
+                                error = function (e){
+                                    setwd(whereAmI)
+                                    stop("Cannot find log files, something went wrong")})
+    } else {
+        logfiles <- tryCatch(getLogs(outputLoc,outputNames,type="both"),
+                                error = function (e){
+                                    setwd(whereAmI)
+                                    stop("Cannot find log files, something went wrong")})
+    }
+    ## list.files(outputLoc)[grep("log$",list.files(outputLoc))]#creating a vector for logfilenames
+
+###############################################    
+#############LOG SECTION#######################
+###############################################
+
+    if(skipSpinup){
+        errorsign <- readErrors(outputLoc=outputLoc,logfiles=logfiles,type="normal")
+    } else {
+
+        perror <- readErrors(outputLoc=outputLoc,logfiles=logfiles)                                    #vector of spinup and normalrun error
+        
+        
+        ##if errorsign is 1 there is error, if it is 0 everything ok
+        perror[is.na(perror)]<-0
+        if(length(perror)>sum(perror)){
+        errorsign <- 1
+        } else {
+            if(length(perror)==1){
+                errorsign <- 1
+            } else {
+                if(spincrash){
+                    errorsign <- 1
+                } else {
+                    errorsign <- 0
+                } }
+        }
+        
+        
+        
+    }
+        
+    
+    
+
+    if(keepEpc){#if keepepc option turned on
+
+        if(length(unique(dirname(epc)))>1){
+            stop("Why are you playing with my nervs? Seriously? You hold your epc-s in different folders?")
+        } else {
+            if(skipSpinup){
+                stampAndDir(stampDir=EPCS, wrongDir=WRONGEPC, names=epc[2], type="general", errorsign=errorsign, logfiles=logfiles)
+            }
+            stampAndDir(stampDir=EPCS, wrongDir=WRONGEPC, names=epc, type="general", errorsign=errorsign, logfiles=logfiles)
+
+        }
+    }
+    
+
+
+    if(debugging){ #debugging is boolean
+        logfiles <- file.path(outputLoc,logfiles)
+        
+                       stampAndDir(stampDir=dirName, wrongDir=dirERROR, names=logfiles, type="general",errorsign=errorsign,logfiles=logfiles)}
+  
+    
+    #cleanupMuso(location=outputLoc,deep = FALSE)
+    if(errorsign==1){
+        return("Modell Failure")
+    }
+
+    if(timee=="d"){
+        colnames(Reva) <- unlist(settings$outputVars[[1]])
+    } else {
+        if(timee=="y")
+            colnames(Reva) <- unlist(settings$outputVars[[2]])
+    }
+
+
+
+    if(leapYear){
+        Reva <- corrigMuso(settings,Reva)
+        rownames(Reva) <- musoDate(settings$startYear,settings$numYears)
+    } else { 
+        rownames(Reva) <- musoDate(settings$startYear, settings$numYears, corrigated=FALSE)
+    }
+
+
+    
+    if(export!=FALSE){
+        setwd(whereAmI)
+
+        ## switch(fextension(export),
+        ##        "csv"=(write.csv(Reva,export)),
+        ##        "xlsx"=(),
+        ##        "odt"=
+                  
+            
+        ## )
+        write.csv(Reva,export)
+        
+    } else{
+        setwd(whereAmI)
+        return(Reva)}
+}

RBBGCMuso/R/changeMuso.R

@@ -5,23 +5,45 @@
 #' 
 #' @author Roland Hollos
 #' @keywords internal
+#' 
 
-changemulline <- function(filename,calibrationPar,contents){
+changemulline <- function(filePaths, calibrationPar, contents, fileOut, fileToChange, modifyOriginal=FALSE){
-    ##This is the function which is capable change multiple specific lines to other using their row numbers.
+    selectFileToWrite <- function(filePaths, fileTochange){
-    ##The function uses the previous changspecline function to operate.
+        if(fileToChange == "epc"){
-    ##From now changespecline is in the forarcheologist file, because its no longer needed
+            return(1)
-    varnum <- length(calibrationPar)
+        } else{
-    contents <- as.list(contents)
+            return(2)            
-    if(length(contents)!=varnum)
-    {
-        stop(" number of the values is not the same as the number of the changed parameters")
         }
 
-    readedFile = readLines(filename,-1)
-
-    for(i in 1:varnum){
-        readedFile[calibrationPar[i]] <- paste(contents[[i]],collapse = " ")
     }
     
-    writeLines(unlist(readedFile),filename)
+    if(xor(is.list(calibrationPar), is.list(contents)) ){
+        stop("If you change epc and ini files also, you have to use list for calibrationPar, and paramateters.")
     }
+     
+    if(!is.element(fileToChange,c("ini","epc","both"))){
+        stop("RBBGCMuso can only change ini or epc file, so fileToChange can be 'epc/ini/both'")
+    }
+
+    if(fileToChange == "epc" | fileToChange == "ini"){
+        parMat<-cbind(calibrationPar, contents)
+        parMat[order(parMat[,1]),]
+        changeMusoC(inFile = filePaths[selectFileToWrite(filePaths, fileToChange)],
+                    outFile = fileOut[selectFileToWrite(filePaths, fileToChange)],
+                    parMat)
+    }
+    
+    if(fileToChange == "both"){
+        parMat<-list()
+        parMat[[1]]<-cbind(calibrationPar[[1]], contents[[1]])
+        parMat[[1]][order(parMat[[1]][,1]),]
+        parMat[[2]]<-cbind(calibrationPar[[2]], contents[[2]])
+        parMat[[2]][order(parMat[[2]][,1]),]
+        
+        changeMusoC(filePaths[1],fileOut[1],parMat[[1]] )
+        changeMusoC(filePaths[2],fileOut[2],parMat[[2]] )
+    }
+
+}
+
+

RBBGCMuso/R/changeMuso.R~

@@ -0,0 +1,50 @@
+#' This is the function which is capable change multiple specific lines to other using their row numbers.
+#'
+#' he function uses the previous changspecline function to operate.
+         ##From now changespecline is in the forarcheologist file, because its no longer needed
+#' 
+#' @author Roland Hollos
+#' @keywords internal
+#' 
+
+changemulline <- function(filePaths, calibrationPar, contents, fileOut=NULL, fileToChange, modifyOriginal=FALSE){
+    selectFileToWrite <- function(filePaths, fileTochange){
+        if(fileToChange == "epc"){
+            return(filePaths[1])
+        } else{
+            return(filePaths[2])            
+        }
+
+    }
+
+    if(modifyOriginal){
+        fileOut <- filePaths
+        mutated <- TRUE
+    } else {
+        fileOut <- sapply(filePaths, function (x){
+            paste0(tools::file_path_sans_ext(basename(x)),"-tmp.",tools::file_ext(x))
+        })
+        mutated <- FALSE
+    }
+    
+    if(xor(is.list(calibrationPar), is.list(contents)) ){
+        stop("If you change epc and ini files also, you have to use list for calibrationPar, and paramateters.")
+    }
+    
+    if(!is.element(fileToChange,c("ini","epc","both"))){
+        stop("RBBGCMuso can only change ini or epc file, so fileToChange can be 'epc/ini/both'")
+    }
+
+    if(fileToChange == "epc" | fileToChange == "ini"){
+        changeMusoC(selectFileToWrite(filePaths, fileToChange),fileOut[1], cbind(calibrationPar, contents))
+    }
+    
+    if(fileToChange == "both"){
+        changeMusoC(filePaths[1],fileOut[1], cbind(calibrationPar[[1]],contents[[1]]))
+        changeMusoC(filePaths[2],fileOut[2], cbind(calibrationPar[[1]],contents[[2]]))
+    }
+
+    return(mutated)
+}
+
+

RBBGCMuso/R/musoExamples.R

@@ -1,56 +0,0 @@
-#'copyMusoExampleTo
-#'
-#'With this function you can copy RBBGCMuso example library whereever you want
-#'
-#'@param example The name of the example file, if it is NULL tcl/tk menu will pop up to select.
-#'@param destination The destination where the example files will be copied. 
-#'@export
-
-copyMusoExamleTo <- function(example = NULL, destination = NULL){
-    WindowsP <- Sys.info()[1] == "Windows"
-
-    chooseExample <- function(){
-            choiceWin <- tcltk::tktoplevel()
-            tcltk::tclRequire("BWidget")
-            tcltk::tktitle(choiceWin) <- "Choose an example!"
-            tcltk::tcl("wm","geometry",choiceWin,"200x50")
-            tcltk::tcl("wm", "attributes", choiceWin, topmost=TRUE)
-            choiceValues <-  basename(list.dirs(system.file("examples","",package = "RBBGCMuso"),recursive = FALSE))
-            choices <- tcltk::tkwidget(choiceWin,"ComboBox",
-                                       editable = FALSE, values = choiceValues,
-                                       textvariable = tcltk::tclVar(choiceValues[1]))
-            tcltk::tkpack(choices)
-            choiceValue <- NA
-            closeSelection <- tcltk::tkwidget(choiceWin,"button",text ="Select", command =function (){
-                choiceValue <<- tcltk::tclvalue(tcltk::tcl(choices,"get"))
-                tcltk::tkdestroy(choiceWin)
-            })
-            
-            tcltk::tkpack(closeSelection)
-            while(as.numeric(tcltk::tclvalue(tcltk::tcl("winfo","exists",choiceWin)))){
-                
-            }
-            return(choiceValue)
-    }
-
-
-    
-    if(is.null(example)){
-        cExample<-paste0(system.file("examples","",package = "RBBGCMuso"),"/",chooseExample())       
-    }
-    
-    if(is.null(destination)){
-        destination<-tcltk::tk_choose.dir(getwd(), "Choose folder to copy the examples!")
-    }
-    
-    currDir <- getwd()
-    setwd(cExample)
-    if(!WindowsP){
-        file.copy("./bin/muso", destination)
-    } else {
-        file.copy("./bin/muso.exe", destination)
-        file.copy("./bin/cygwin1.dll", destination)
-    }
-        file.copy(grep("bin", list.files(), value = TRUE, invert = TRUE),destination)
-     setwd(currDir)
-}

RBBGCMuso/R/musoFilter.R

@@ -0,0 +1,4 @@
+musoFilter <- function(text){
+    eval(parse(paste0("filter(.,",text,")"))) %>%
+        tbl_df
+}

RBBGCMuso/R/musoMont.R

@@ -0,0 +1,210 @@
+#' musoMont
+#'
+#' This function does monteCarlo on BiomeBGC-MuSo. It samples specified modell variables in given rangge from conditional multivariate uniform distribution, and runs the modell for each run.
+#' @author Roland Hollos
+#' @param settings A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically. 
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param inputDir The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file. 
+#' @param iterations number of the monteCarlo run.
+#' @param preTag It will be the name of the output files. For example preTag-1.csv, pretag-2csv...
+#' @param outputType  This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet. 
+#' @param fun If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.
+#' @param varIndex This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.
+#' @param debugging If you set this parameter, you can save every logfile, and RBBGCMuso will select those which contains errors.
+#' @param keepEpc if you set keepEpc also true, it will save every selected epc file, and put the wrong ones in the WRONGEPC directory.
+#' @export
+
+musoMont <- function(settings=NULL,
+                     parameters=NULL,
+                     inputDir = "./",
+                     outLoc = "./calib",
+                     iterations = 10,
+                     preTag = "mont-",
+                     outputType = "moreCsv",
+                     fun=mean,
+                     varIndex = 1,
+                     outVars = NULL,
+                     silent = TRUE,
+                     skipSpinup = TRUE,
+                     debugging = FALSE,
+                     keepEpc = FALSE,
+                     constrains = NULL,
+                     ...){
+
+
+    readValuesFromEpc  <- function(epc, linums){
+        epcFile <- readLines(epc)
+        rows <- numeric(2)
+        values <- sapply(linums, function(x){
+            rows[1] <- as.integer(x)
+            rows[2] <- as.integer(round(100*x)) %% 10 + 1 
+            epcFile <- readLines(epc)
+            selRow <- unlist(strsplit(epcFile[rows[1]], split= "[\t ]"))
+            selRow <- selRow[selRow!=""]
+            return(as.numeric(selRow[rows[2]]))
+            
+        })
+        
+        return(values)
+    } 
+    
+    if(is.null(parameters)){
+        parameters <- tryCatch(read.csv("parameters.csv", stringsAsFactor=FALSE), error = function (e) {
+            stop("You need to specify a path for the parameters.csv, or a matrix.")
+        })
+    } else {
+        if((!is.list(parameters)) & (!is.matrix(parameters))){
+             parameters <- tryCatch(read.csv(parameters, stringsAsFactor=FALSE), error = function (e){
+                                         stop("Cannot find neither parameters file neither the parameters matrix")
+                                     })
+        }}
+    
+    outLocPlain <- basename(outLoc)
+    currDir <- getwd()
+
+    if(!dir.exists(outLoc)){
+        dir.create(outLoc)
+        warning(paste(outLoc," is not exists, so it was created"))
+    }
+
+    outLoc <- normalizePath(outLoc)
+
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    
+
+    if(is.null(outVars)){
+        numVars <- length(settings$outputVars[[1]])
+        outVarNames <- settings$outputVars[[1]]
+    } else {
+        numVars <- length(outVars)
+        outVarNames <- sapply(outVars, musoMapping)
+    }
+    
+    parameterNames <- parameters[,1]
+    # settings$calibrationPar <- A[,1] #:LATER:
+    pretag <- file.path(outLoc,preTag)
+    npar <- length(settings$calibrationPar)
+    
+    ##reading the original epc file at the specified
+    ## row numbers
+    if(iterations < 3000){
+        randVals <- musoRand(parameters = parameters,constrains = constrains, iterations = 3000)
+        randVals[[2]]<- randVals[[2]][sample(1:3000,iterations),]
+    } else {
+        randVals <- musoRand(parameters = parameters,constrains = constrains, iterations = iterations)
+    }
+    
+    origEpc <- readValuesFromEpc(settings$epc[2],parameters[,2])
+    
+    ## Prepare the preservedEpc matrix for the faster
+    ##  run.
+
+    pretag <- file.path(outLoc,preTag)
+    
+    ## Creating function for generating separate
+    ## csv files for each run
+
+    progBar <- txtProgressBar(1,iterations,style=3)
+
+    
+    moreCsv <- function(){
+        randValues <- randVals[[2]]
+        settings$calibrationPar <- randVals[[1]]
+        ## randValues <- randValues[,randVals[[1]] %in% parameters[,2]][,rank(parameters[,2])]
+        modellOut <- matrix(ncol = numVars, nrow = iterations + 1)
+
+        origModellOut <- calibMuso(silent=TRUE)
+        write.csv(x=origModellOut, file=paste0(pretag,1,".csv"))
+
+        if(!is.list(fun)){
+            funct <- rep(list(fun), numVars)
+        }
+            
+        tmp2 <- numeric(numVars)
+
+        for(j in 1:numVars){
+            tmp2[j]<-funct[[j]](origModellOut[,j])
+        }
+        modellOut[1,]<- tmp2
+        
+        for(i in 2:(iterations+1)){
+            tmp <- calibMuso(settings = settings,
+                             parameters = randValues[(i-1),],
+                             silent= TRUE,
+                             skipSpinup = skipSpinup,
+                             keepEpc = keepEpc,
+                             debugging = debugging,
+                             outVars = outVars)
+            
+
+            for(j in 1:numVars){
+                tmp2[j]<-funct[[j]](tmp[,j])
+            }
+            
+            modellOut[i,]<- tmp2
+            write.csv(x=tmp, file=paste0(pretag,(i+1),".csv"))
+            setTxtProgressBar(progBar,i)
+        }
+
+        paramLines <- parameters[,2]
+        paramLines <- order(paramLines)
+        randInd <- randVals[[1]][(randVals[[1]] %in% parameters[,2])]
+        randInd <- order(randInd)
+        
+       
+        epcStrip <- rbind(origEpc[order(parameters[,2])],
+                          randValues[,randVals[[1]] %in% parameters[,2]][,randInd])
+        
+        
+            preservedEpc <- cbind(epcStrip,
+                                  modellOut)
+        colnames(preservedEpc) <- c(parameterNames[paramLines], sapply(outVarNames, function (x) paste0("mod.", x)))
+        return(preservedEpc)
+    }
+    
+    ## Creating function for generating one
+    ## csv files for each run
+    
+    oneCsv <- function () {
+        stop("This function is not implemented yet")
+        ## numDays <- settings$numdata[1]
+        ## if(!onDisk){
+        ##     for(i in 1:iterations){
+                
+        ##         parVar <- apply(parameters,1,function (x) {
+        ##             runif(1, as.numeric(x[3]), as.numeric(x[4]))})
+                
+        ##         preservedEpc[(i+1),] <- parVar
+        ##         exportName <- paste0(preTag,".csv")
+        ##         write.csv(parvar,"preservedEpc.csv",append=TRUE)
+        ##         calibMuso(settings,debugging = "stamplog",
+        ##                   parameters = parVar,keepEpc = TRUE) %>%
+        ##             {mutate(.,iD = i)} %>%
+        ##             {write.csv(.,file=exportName,append=TRUE)}
+        ##     }
+            
+        ##     return(preservedEpc)
+        ## } else {
+            
+        ## }
+    }
+    
+    netCDF <- function () {
+        stop("This function is not inplemented yet")
+    }
+    
+    ## Call one function according to the outputType
+    switch(outputType,
+           "oneCsv" = (a <- oneCsv()),
+           "moreCsv" = (a <- moreCsv()),
+           "netCDF" = (a <- netCDF()))
+    write.csv(a,"preservedEpc.csv")
+   
+    setwd(currDir)
+    return(a)
+}
+

RBBGCMuso/R/musoMont.R

@@ -0,0 +1,205 @@
+#' musoMont
+#'
+#' This function does monteCarlo on BiomeBGC-MuSo. It samples specified modell variables in given rangge from conditional multivariate uniform distribution, and runs the modell for each run.
+#' @author Roland Hollos
+#' @param settings A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically. 
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param inputDir The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file. 
+#' @param iterations number of the monteCarlo run.
+#' @param preTag It will be the name of the output files. For example preTag-1.csv, pretag-2csv...
+#' @param outputType  This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet. 
+#' @param fun If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.
+#' @param varIndex This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.
+#' @param debugging If you set this parameter, you can save every logfile, and RBBGCMuso will select those which contains errors.
+#' @param keepEpc if you set keepEpc also true, it will save every selected epc file, and put the wrong ones in the WRONGEPC directory.
+#' @export
+
+musoMonte <- function(settings=NULL,
+                     parameters=NULL,
+                     inputDir = "./",
+                     outLoc = "./calib",
+                     iterations = 10,
+                     preTag = "mont-",
+                     outputType = "moreCsv",
+                     fun=mean,
+                     varIndex = 1,
+                     silent = TRUE,
+                     skipSpinup = FALSE,
+                     debugging = FALSE,
+                     keepEpc = FALSE,
+                     ...){
+
+    getEpcValue  <- function(epc, linum){
+        numcord <- numeric(3)
+        numcord[1] <- as.integer(linNum)
+        linNum <- as.integer(round(linNum * 100))
+        numcord[3] <-linNum %% 10 +1 
+        numcord[2] <- (linNum %/% 10) %% 10 + 1
+        numcord 
+    } 
+    
+    
+    if(is.null(parameters)){
+        parameters <- tryCatch(read.csv("parameters.csv"), error = function (e) {
+            stop("You need to specify a path for the parameters.csv, or a matrix.")
+        })
+    } else {
+        if((!is.list(parameters)) & (!is.matrix(parameters))){
+             parameters <- tryCatch(read.csv(parameters), error = function (e){
+                                         stop("Cannot find neither parameters file neither the parameters matrix")
+                                     })
+        }}
+    
+    outLocPlain <- basename(outLoc) #Where to put the csv outputs 
+    currDir <- getwd() # just to go back, It is unlikely to be used
+    inputDir <- normalizePath(inputDir) # Where are the model files.
+
+    if(!dir.exists(outLoc)){
+        dir.create(outLoc)
+        warning(paste(outLoc," is not exists, so it was created"))
+    }
+
+    outLoc <- normalizePath(outLoc)
+
+    
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    parameterNames <- parameters[,1]
+    parReal <- parameters[,-1]
+    Otable <- OtableMaker(parReal)
+    A <- as.matrix(Otable[[1]][,c(2,4,5,6)])
+    B <- as.matrix(Otable[[2]])
+    settings$calibrationPar <- A[,1] 
+    pretag <- file.path(outLoc,preTag)
+    npar <- length(settings$calibrationPar)
+    
+    ##reading the original epc file at the specified
+    ## row numbers
+    
+    origEpcFile <- readLines(settings$epcInput[2])
+    
+    origEpc <- unlist(lapply(settings$calibrationPar, function (x) {
+        as.numeric(unlist(strsplit(origEpcFile[x],split="[\t ]"))[1])
+    }))
+    
+    ## Prepare the preservedEpc matrix for the faster
+    ##  run.
+    preservedEpc <- matrix(nrow = (iterations +1 ), ncol = npar)
+    preservedEpc[1,] <- origEpc
+    Otable[[1]][,1]  <- (as.character(Otable[[1]][,1]))
+    for(i in parameters[,2]){
+        Otable[[1]][Otable[[1]][,2]==i,1] <- as.character(parameters[parameters[,2]==i,1])
+    }
+    
+    colnames(preservedEpc) <- Otable[[1]][,1]
+    preservedEpc <- cbind(preservedEpc,rep(NA,(iterations+1)))
+    colnames(preservedEpc)[(npar+1)] <- "y"
+    ## Save the backupEpc, while change the settings
+    ## variable and set the output.
+    file.copy(settings$epc[2],"savedEpc",overwrite = TRUE) # do I need this? 
+    pretag <- file.path(outLoc,preTag)
+    
+    ## Creating function for generating separate
+    ## csv files for each run
+
+    progBar <- txtProgressBar(1,iterations,style=3)
+
+    modelRun <- function(settings, debugging, parameters, keepEpc, silent, skipSpinup){
+        if(!skipSpinup){
+            calibMuso(settings, debugging = debugging, parameters = parameters, keepEpc = keepEpc, silent = silent)        
+        } else {
+            normalMuso(settings, debugging = debugging, parameters = parameters, keepEpc = keepEpc, silent = silent)
+        } 
+        
+    }
+    
+    moreCsv <- function(){
+
+        if(skipSpinup){#skipSpinup is boolean
+            spinupMuso(settings = settings , silent = silent)
+        }
+        a <- numeric(iterations+1)
+        tempData <- modelRun(settings=settings,
+                             debugging = debugging,
+                             parameters = origEpc,
+                             keepEpc = keepEpc,
+                             silent = silent,
+                             skipSpinup = skipSpinup)
+        ## tempData <- calibMuso(settings, debugging = "stamplog", parameters = origEpc,keepEpc = TRUE,silent = silent)
+        a[1] <- tryCatch(fun(tempData[,varIndex]),error=function(e){return(NA)})
+        preservedEpc[1,(npar+1)] <- a[1]
+        write.table(t(preservedEpc[1,]),row.names = FALSE,"preservedEpc.csv",sep=",")
+        write.csv(x=tempData, file=paste0(preTag,1,".csv"))
+        for(i in 1:iterations){
+            parVar <- musoRandomizer(A,B)[,2]
+            preservedEpc[(i+1),] <- c(parVar,NA)
+            exportName <- paste0(preTag,(i+1),".csv")
+            tempData <- modelRun(settings = settings,
+                                 debugging = debugging,
+                                 parameters = parVar,
+                                 keepEpc = keepEpc,
+                                 silent=silent,
+                                 skipSpinup =skipSpinup)
+            write.csv(x=tempData,file=exportName)
+            
+            preservedEpc[(i+1),(npar+1)] <- a[i+1]<- tryCatch(fun(tempData[,varIndex]),error=function(e){return(NA)})
+            write.table(t(preservedEpc[(i+1),]),file="preservedEpc.csv",row.names=FALSE,col.names=FALSE, append=TRUE,sep=",")
+            setTxtProgressBar(progBar,i)
+        }
+        cat("\n")
+        return(preservedEpc)
+    }
+    
+    ## Creating function for generating one
+    ## csv files for each run
+    
+    oneCsv <- function () {
+        numDays <- settings$numdata[1]
+        if(!onDisk){
+            for(i in 1:iterations){
+                
+                parVar <- apply(parameters,1,function (x) {
+                    runif(1, as.numeric(x[3]), as.numeric(x[4]))})
+                
+                preservedEpc[(i+1),] <- parVar
+                exportName <- paste0(preTag,".csv")
+                write.csv(parvar,"preservedEpc.csv",append=TRUE)
+                calibMuso(settings,debugging = "stamplog",
+                          parameters = parVar,keepEpc = TRUE) %>%
+                    {mutate(.,iD = i)} %>%
+                    {write.csv(.,file=exportName,append=TRUE)}
+            }
+            
+            return(preservedEpc)
+        } else {
+            
+        }
+    }
+    
+    netCDF <- function () {
+        stop("This function is not inplemented yet")
+    }
+    
+    ## Call one function according to the outputType
+    switch(outputType,
+           "oneCsv" = (a <- oneCsv()),
+           "moreCsv" = (a <- moreCsv()),
+           "netCDF" = (a <- netCDF()))
+    
+    ## Change back the epc file to the original
+    for(i in file.path("./",grep(outLocPlain, list.files(inputDir), invert = TRUE, value = TRUE))){
+        file.remove(i,recursive=TRUE)
+    }
+    for(i in list.files()){
+        file.copy(i,outLoc,recursive=TRUE,overwrite = TRUE)
+    }
+    
+    unlink(tmp,recursive = TRUE)
+    setwd(currDir)
+    file.copy("savedEpc",settings$epc[2],overwrite = TRUE)
+    return(a)
+}
+

RBBGCMuso/R/musoMont.R ~

@@ -0,0 +1,208 @@
+#' musoMont
+#'
+#' This function does monteCarlo on BiomeBGC-MuSo. It samples specified modell variables in given rangge from conditional multivariate uniform distribution, and runs the modell for each run.
+#' @author Roland Hollos
+#' @param settings A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically. 
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param inputDir The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file. 
+#' @param iterations number of the monteCarlo run.
+#' @param preTag It will be the name of the output files. For example preTag-1.csv, pretag-2csv...
+#' @param outputType  This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet. 
+#' @param fun If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.
+#' @param varIndex This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.
+#' @param debugging If you set this parameter, you can save every logfile, and RBBGCMuso will select those which contains errors.
+#' @param keepEpc if you set keepEpc also true, it will save every selected epc file, and put the wrong ones in the WRONGEPC directory.
+#' @export
+
+musoMonte <- function(settings=NULL,
+                     parameters=NULL,
+                     inputDir = "./",
+                     outLoc = "./calib",
+                     iterations = 10,
+                     preTag = "mont-",
+                     outputType = "moreCsv",
+                     fun=mean,
+                     varIndex = 1,
+                     silent = TRUE,
+                     skipSpinup = FALSE,
+                     debugging = FALSE,
+                     keepEpc = FALSE,
+                     ...){
+
+    if(is.null(parameters)){
+        parameters <- tryCatch(read.csv("parameters.csv"), error = function (e) {
+            stop("You need to specify a path for the parameters.csv, or a matrix.")
+        })
+    } else {
+        if((!is.list(parameters)) & (!is.matrix(parameters))){
+             parameters <- tryCatch(read.csv(parameters), error = function (e){
+                                         stop("Cannot find neither parameters file neither the parameters matrix")
+                                     })
+        }}
+    
+    outLocPlain <- basename(outLoc) #Where to put the csv outputs 
+    currDir <- getwd() # just to go back, It is likely to not to be used
+    inputDir <- normalizePath(inputDir) # Where are the model files.
+    tmp <- file.path(outLoc,"tmp/") 
+
+    if(!dir.exists(outLoc)){
+        dir.create(outLoc)
+        warning(paste(outLoc," is not exists, so it was created"))
+    }
+
+    if(dir.exists(tmp)){
+        stop("There is a tmp directory inside the output location, please replace it. tmp is an important temporary directory for the function")
+    }
+    dir.create(tmp)
+    outLoc <- normalizePath(outLoc)
+    tmp <- normalizePath(tmp)
+
+    inputFiles <- file.path(inputDir,grep(basename(outLoc),list.files(inputDir),invert = TRUE,value = TRUE))
+
+ 
+    for(i in inputFiles){
+        file.copy(i,tmp)
+    }
+    setwd(tmp)
+    
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    parameterNames <- parameters[,1]
+    parReal <- parameters[,-1]
+    Otable <- OtableMaker(parReal)
+    A <- as.matrix(Otable[[1]][,c(2,4,5,6)])
+    B <- as.matrix(Otable[[2]])
+    settings$calibrationPar <- A[,1] 
+    pretag <- file.path(outLoc,preTag)
+    npar <- length(settings$calibrationPar)
+    
+    ##reading the original epc file at the specified
+    ## row numbers
+    
+    origEpcFile <- readLines(settings$epcInput[2])
+    
+    origEpc <- unlist(lapply(settings$calibrationPar, function (x) {
+        as.numeric(unlist(strsplit(origEpcFile[x],split="[\t ]"))[1])
+    }))
+    
+    ## Prepare the preservedEpc matrix for the faster
+    ##  run.
+    preservedEpc <- matrix(nrow = (iterations +1 ), ncol = npar)
+    preservedEpc[1,] <- origEpc
+    Otable[[1]][,1]  <- (as.character(Otable[[1]][,1]))
+    for(i in parameters[,2]){
+        Otable[[1]][Otable[[1]][,2]==i,1] <- as.character(parameters[parameters[,2]==i,1])
+    }
+    
+    colnames(preservedEpc) <- Otable[[1]][,1]
+    preservedEpc <- cbind(preservedEpc,rep(NA,(iterations+1)))
+    colnames(preservedEpc)[(npar+1)] <- "y"
+    ## Save the backupEpc, while change the settings
+    ## variable and set the output.
+    file.copy(settings$epc[2],"savedEpc",overwrite = TRUE) # do I need this? 
+    pretag <- file.path(outLoc,preTag)
+    
+    ## Creating function for generating separate
+    ## csv files for each run
+
+    progBar <- txtProgressBar(1,iterations,style=3)
+
+    modelRun <- function(settings, debugging, parameters, keepEpc, silent, skipSpinup){
+        if(!skipSpinup){
+            calibMuso(settings, debugging = debugging, parameters = parameters, keepEpc = keepEpc, silent = silent)        
+        } else {
+            normalMuso(settings, debugging = debugging, parameters = parameters, keepEpc = keepEpc, silent = silent)
+        } 
+        
+    }
+    
+    moreCsv <- function(){
+
+        if(skipSpinup){#skipSpinup is boolean
+            spinupMuso(settings = settings , silent = silent)
+        }
+        a <- numeric(iterations+1)
+        tempData <- modelRun(settings=settings,
+                             debugging = debugging,
+                             parameters = origEpc,
+                             keepEpc = keepEpc,
+                             silent = silent,
+                             skipSpinup = skipSpinup)
+        ## tempData <- calibMuso(settings, debugging = "stamplog", parameters = origEpc,keepEpc = TRUE,silent = silent)
+        a[1] <- tryCatch(fun(tempData[,varIndex]),error=function(e){return(NA)})
+        preservedEpc[1,(npar+1)] <- a[1]
+        write.table(t(preservedEpc[1,]),row.names = FALSE,"preservedEpc.csv",sep=",")
+        write.csv(x=tempData, file=paste0(preTag,1,".csv"))
+        for(i in 1:iterations){
+            parVar <- musoRandomizer(A,B)[,2]
+            preservedEpc[(i+1),] <- c(parVar,NA)
+            exportName <- paste0(preTag,(i+1),".csv")
+            tempData <- modelRun(settings = settings,
+                                 debugging = debugging,
+                                 parameters = parVar,
+                                 keepEpc = keepEpc,
+                                 silent=silent,
+                                 skipSpinup =skipSpinup)
+            write.csv(x=tempData,file=exportName)
+            
+            preservedEpc[(i+1),(npar+1)] <- a[i+1]<- tryCatch(fun(tempData[,varIndex]),error=function(e){return(NA)})
+            write.table(t(preservedEpc[(i+1),]),file="preservedEpc.csv",row.names=FALSE,col.names=FALSE, append=TRUE,sep=",")
+            setTxtProgressBar(progBar,i)
+        }
+        cat("\n")
+        return(preservedEpc)
+    }
+    
+    ## Creating function for generating one
+    ## csv files for each run
+    
+    oneCsv <- function () {
+        numDays <- settings$numdata[1]
+        if(!onDisk){
+            for(i in 1:iterations){
+                
+                parVar <- apply(parameters,1,function (x) {
+                    runif(1, as.numeric(x[3]), as.numeric(x[4]))})
+                
+                preservedEpc[(i+1),] <- parVar
+                exportName <- paste0(preTag,".csv")
+                write.csv(parvar,"preservedEpc.csv",append=TRUE)
+                calibMuso(settings,debugging = "stamplog",
+                          parameters = parVar,keepEpc = TRUE) %>%
+                    {mutate(.,iD = i)} %>%
+                    {write.csv(.,file=exportName,append=TRUE)}
+            }
+            
+            return(preservedEpc)
+        } else {
+            
+        }
+    }
+    
+    netCDF <- function () {
+        stop("This function is not inplemented yet")
+    }
+    
+    ## Call one function according to the outputType
+    switch(outputType,
+           "oneCsv" = (a <- oneCsv()),
+           "moreCsv" = (a <- moreCsv()),
+           "netCDF" = (a <- netCDF()))
+    
+    ## Change back the epc file to the original
+    for(i in file.path("./",grep(outLocPlain, list.files(inputDir), invert = TRUE, value = TRUE))){
+        file.remove(i,recursive=TRUE)
+    }
+    for(i in list.files()){
+        file.copy(i,outLoc,recursive=TRUE,overwrite = TRUE)
+    }
+    
+    unlink(tmp,recursive = TRUE)
+    setwd(currDir)
+    file.copy("savedEpc",settings$epc[2],overwrite = TRUE)
+    return(a)
+}
+

RBBGCMuso/R/musoMont.R~

@@ -0,0 +1,200 @@
+#' musoMont
+#'
+#' This function does monteCarlo on BiomeBGC-MuSo. It samples specified modell variables in given rangge from conditional multivariate uniform distribution, and runs the modell for each run.
+#' @author Roland Hollos
+#' @param settings A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically. 
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param inputDir The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file. 
+#' @param iterations number of the monteCarlo run.
+#' @param preTag It will be the name of the output files. For example preTag-1.csv, pretag-2csv...
+#' @param outputType  This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet. 
+#' @param fun If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.
+#' @param varIndex This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.
+#' @param debugging If you set this parameter, you can save every logfile, and RBBGCMuso will select those which contains errors.
+#' @param keepEpc if you set keepEpc also true, it will save every selected epc file, and put the wrong ones in the WRONGEPC directory.
+#' @export
+
+musoMont <- function(settings=NULL,
+                     parameters=NULL,
+                     inputDir = "./",
+                     outLoc = "./calib",
+                     iterations = 10,
+                     preTag = "mont-",
+                     outputType = "moreCsv",
+                     fun=mean,
+                     varIndex = 1,
+                     outVars = NULL,
+                     silent = TRUE,
+                     skipSpinup = TRUE,
+                     debugging = FALSE,
+                     keepEpc = FALSE,
+                     constrains = NULL,
+                     ...){
+
+
+    readValuesFromEpc  <- function(epc, linums){
+        rows <- numeric(2)
+        values <- sapply(linums, function(x){
+            rows[1] <- as.integer(x)
+            rows[2] <- as.integer(round(100*x)) %% 10 + 1 
+            epcFile <- readLines(epc)
+            
+            return(as.numeric(unlist(strsplit(epcFile[rows[1]], split= "[\t ]"))[rows[2]]))
+            
+        })
+        
+        return(values)
+    } 
+    
+    if(is.null(parameters)){
+        parameters <- tryCatch(read.csv("parameters.csv", stringsAsFactor=FALSE), error = function (e) {
+            stop("You need to specify a path for the parameters.csv, or a matrix.")
+        })
+    } else {
+        if((!is.list(parameters)) & (!is.matrix(parameters))){
+             parameters <- tryCatch(read.csv(parameters, stringsAsFactor=FALSE), error = function (e){
+                                         stop("Cannot find neither parameters file neither the parameters matrix")
+                                     })
+        }}
+    
+    outLocPlain <- basename(outLoc)
+    currDir <- getwd()
+
+    if(!dir.exists(outLoc)){
+        dir.create(outLoc)
+        warning(paste(outLoc," is not exists, so it was created"))
+    }
+
+    outLoc <- normalizePath(outLoc)
+
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    
+
+    if(is.null(outVars)){
+        numVars <- length(settings$outputVars[[1]])
+        outVarNames <- settings$outputVars[[1]]
+    } else {
+        numVars <- length(outVars)
+        outVarNames <- sapply(outVars, musoMapping)
+    }
+    
+    parameterNames <- parameters[,1]
+    # settings$calibrationPar <- A[,1] #:LATER:
+    pretag <- file.path(outLoc,preTag)
+    npar <- length(settings$calibrationPar)
+    
+    ##reading the original epc file at the specified
+    ## row numbers
+    if(iterations < 3000){
+        randVals <- musoRand(parameters = parameters,constrains = constrains, iterations = 3000)
+        randVals[[2]]<- randVals[[2]][sample(1:3000,iterations),]
+    } else {
+        randVals <- musoRand(parameters = parameters,constrains = constrains, iterations = iterations)
+    }
+    
+    origEpc <- readValuesFromEpc(settings$epc[2],parameters[,2])
+    
+    ## Prepare the preservedEpc matrix for the faster
+    ##  run.
+
+    pretag <- file.path(outLoc,preTag)
+    
+    ## Creating function for generating separate
+    ## csv files for each run
+
+    progBar <- txtProgressBar(1,iterations,style=3)
+
+    
+    moreCsv <- function(){
+        randValues <- randVals[[2]]
+        settings$calibrationPar <- randVals[[1]]
+        ## randValues <- randValues[,randVals[[1]] %in% parameters[,2]][,rank(parameters[,2])]
+        modellOut <- matrix(ncol = numVars, nrow = iterations + 1)
+
+        origModellOut <- calibMuso(silent=TRUE)
+        write.csv(x=origModellOut, file=paste0(pretag,1,".csv"))
+
+        if(!is.list(fun)){
+            funct <- rep(list(fun), numVars)
+        }
+            
+        tmp2 <- numeric(numVars)
+
+        for(j in 1:numVars){
+            tmp2[j]<-funct[[j]](origModellOut[,j])
+        }
+        modellOut[1,]<- tmp2
+        
+        for(i in 2:(iterations+1)){
+            tmp <- calibMuso(settings = settings,
+                             parameters = randValues[(i-1),],
+                             silent= TRUE,
+                             skipSpinup = skipSpinup,
+                             keepEpc = keepEpc,
+                             debugging = debugging,
+                             outVars = outVars)
+            
+
+            for(j in 1:numVars){
+                tmp2[j]<-funct[[j]](tmp[,j])
+            }
+            
+            modellOut[i,]<- tmp2
+            write.csv(x=tmp, file=paste0(pretag,(i+1),".csv"))
+            setTxtProgressBar(progBar,i)
+        }
+
+
+        
+            preservedEpc <- cbind(rbind(parameters[,2], randValues[,randVals[[1]] %in% parameters[,2]][,rank(parameters[,2])]),
+                                  modellOut)
+        colnames(preservedEpc) <- c(parameterNames, sapply(outVarNames, function (x) paste0("mod.", x)))
+        return(preservedEpc)
+    }
+    
+    ## Creating function for generating one
+    ## csv files for each run
+    
+    oneCsv <- function () {
+        stop("This function is not implemented yet")
+        ## numDays <- settings$numdata[1]
+        ## if(!onDisk){
+        ##     for(i in 1:iterations){
+                
+        ##         parVar <- apply(parameters,1,function (x) {
+        ##             runif(1, as.numeric(x[3]), as.numeric(x[4]))})
+                
+        ##         preservedEpc[(i+1),] <- parVar
+        ##         exportName <- paste0(preTag,".csv")
+        ##         write.csv(parvar,"preservedEpc.csv",append=TRUE)
+        ##         calibMuso(settings,debugging = "stamplog",
+        ##                   parameters = parVar,keepEpc = TRUE) %>%
+        ##             {mutate(.,iD = i)} %>%
+        ##             {write.csv(.,file=exportName,append=TRUE)}
+        ##     }
+            
+        ##     return(preservedEpc)
+        ## } else {
+            
+        ## }
+    }
+    
+    netCDF <- function () {
+        stop("This function is not inplemented yet")
+    }
+    
+    ## Call one function according to the outputType
+    switch(outputType,
+           "oneCsv" = (a <- oneCsv()),
+           "moreCsv" = (a <- moreCsv()),
+           "netCDF" = (a <- netCDF()))
+    write.csv(a,"preservedEpc.csv")
+   
+    setwd(currDir)
+    return(a)
+}
+

RBBGCMuso/R/musoRand.R

@@ -0,0 +1,171 @@
+#' musoRand
+#'
+#' This funtion samples uniformly from the choosen parameters of the BiomeBGC-Muso model, which parameters are constrained by the model logic.
+#' @author Roland Hollos
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the input-file, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param constrains This is a matrics wich specify the constrain rules for the sampling. Further informations coming son.
+#' @param iteration The number of sample-s. It is adviced to use at least 3000 iteration, because it is generally fast and it can be subsampled later at any time. 
+#' @export
+
+musoRand <- function(parameters, constrains = NULL, iterations=3000){
+
+    if(!is.null(constrains)){
+        constMatrix <- constrains
+    }
+    
+    parameters <- parameters[,-1]
+    constMatrix <- constMatrix[,-1]
+    
+    depTableMaker <- function(constMatrix,parameters){
+	parameters <- parameters[order(parameters[,1]),]
+	constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],c(5,6)]<-parameters[,c(2,3)]
+	logiConstrain <- (constMatrix[,"GROUP"] %in% constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],"GROUP"] &
+			  (constMatrix[,"GROUP"]!=0)) | ((constMatrix[,"INDEX"] %in% parameters[,1]) & (constMatrix[,"GROUP"] == 0))
+	constMatrix<-constMatrix[logiConstrain,]
+	constMatrix <- constMatrix[order(apply(constMatrix[,7:8],1,function(x){x[1]/10+abs(x[2])})),]
+	constMatrix
+    }
+
+    genMat0 <- function(dep){
+	numberOfVariable <- nrow(dep)
+	G <- rbind(diag(numberOfVariable), -1*diag(numberOfVariable))
+	h <- c(dependences[,5], -1*dependences[,6])
+	return(list(G=G,h=h))
+    }
+
+    genMat1 <- function(dep, N){
+
+	## Range <- sapply(list(min,max),function(x){
+	##   x(as.numeric(rownames(dep)))
+	## }) It is more elegant, more general, but slower
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	numberOfVariables <- nrow(dep)
+	G<- -1*diag(numberOfVariables)
+
+	for(i in 1:numberOfVariables){
+	    if(dep[i,4]!=0){
+		G[i,dep[i,4]] <- 1
+	    }
+
+	}
+
+	G<-G[dep[,4]!=0,]
+
+	if(Range[1]==1){
+	    G<-cbind(G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	} else{
+	    if(Range[2]==N){
+		G<-cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G)
+	    } else {
+		G <- cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	    }
+	}
+	return(list(G=G,h=rep(0,nrow(G))))
+    }
+
+    genMat2 <- function(dep, N){
+	G <- rep(1,nrow(dep))
+
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	if(Range[1]==1){
+	    G<-c(G, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		G<-c(numeric(Range[1]-1), G)
+	    } else {
+		G <- c(numeric(Range[1]-1), G, numeric(N-Range[2]))
+	    }
+	}
+
+	G <- t(matrix(sign(dep[2,4])*G))
+	h <- abs(dep[1,4])
+
+	return(list(G=G,h=h))
+    }
+
+    genMat3 <- function(dep, N){
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	E <- rep(1,nrow(dep))
+
+	if(Range[1]==1){
+	    E<-c(E, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		E<-c(numeric(Range[1]-1), E)
+	    } else {
+		E <- c(numeric(Range[1]-1), E, numeric(N-Range[2]))
+	    }
+	}
+
+
+	E <- t(matrix(E))
+	f <- dep[1,4]
+	return(list(E=E,f=f))
+    }
+
+
+    applyRandTypeG <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(type,
+	       invisible(Gh <- genMat1(dep, N)),
+	       invisible(Gh <- genMat2(dep, N)))
+	return(Gh)
+    }
+
+    applyRandTypeE <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(-type,
+	       stop("Not implemented yet"),
+	       stop("Not implemented yet"),
+	       invisible(Ef <- genMat3(dep, N)))
+	return(Ef)
+    }
+
+    dependences <- depTableMaker(constMatrix, parameters)
+    dependences <- cbind(dependences,1:nrow(dependences))
+    colnames(dependences)[ncol(dependences)] <- "rowIndex"
+    numberOfVariable <- nrow(dependences)
+    nonZeroDeps<-dependences[dependences[,"TYPE"]!=0,]
+    if(nrow(nonZeroDeps)!=0){
+        splitedDeps<- split(nonZeroDeps,nonZeroDeps[,"GROUP"])
+        Gh <- list()
+        Ef <- list()
+
+        for(i in 1:length(splitedDeps)){
+            print(splitedDeps[[i]][1,"TYPE"])
+            if(splitedDeps[[i]][1,"TYPE"]>0){
+                Gh[[i]]<-applyRandTypeG(splitedDeps[[i]],nrow(dependences))
+            } else {
+                Ef[[i]] <- applyRandTypeE(splitedDeps[[i]],nrow(dependences))
+            }
+        }
+        
+        Gh0<- genMat0(dependences)
+        G <- do.call(rbind,lapply(Gh,function(x){x$G}))
+        G<- rbind(Gh0$G,G)
+        h <- do.call(c,lapply(Gh,function(x){x$h}))
+        h <- c(Gh0$h,h)
+        E <- do.call(rbind,lapply(Ef,function(x){x$E}))
+        f <- do.call(c,lapply(Ef,function(x){x$f}))
+        randVal <- suppressWarnings(limSolve::xsample(G=G,H=h,E=E,F=f,iter = iterations))$X
+    } else{
+        Gh0<-genMat0(dependences)
+        randVal <- suppressWarnings(limSolve::xsample(G=Gh0$G,H=Gh0$h, iter = iterations))$X
+    }
+    
+    results <- list(INDEX =dependences$INDEX, randVal=randVal)
+    return(results)
+}

RBBGCMuso/R/musoRand.R~

@@ -0,0 +1,165 @@
+#' musoRand
+#'
+#' This funtion samples uniformly from the choosen parameters of the BiomeBGC-Muso model, which parameters are constrained by the model logic.
+#' @author Roland Hollos
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the input-file, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param constrains This is a matrics wich specify the constrain rules for the sampling. Further informations coming son.
+#' @param iteration The number of sample-s. It is adviced to use at least 3000 iteration, because it is generally fast and it can be subsampled later at any time. 
+#' @export
+
+musoRand <- function(parameters, constrains = NULL, iterations=3000){
+
+    if(!is.null(constrains)){
+        constMatrix <- constrains
+    }
+    
+    parameters <- parameters[,-1]
+    constMatrix <- constMatrix[,-1]
+    
+    depTableMaker <- function(constMatrix,parameters){
+	parameters <- parameters[order(parameters[,1]),]
+	constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],c(5,6)]<-parameters[,c(2,3)]
+	logiConstrain <- (constMatrix[,"GROUP"] %in% constMatrix[constMatrix[,"INDEX"] %in% parameters[,1],"GROUP"] &
+			  (constMatrix[,"GROUP"]!=0)) | ((constMatrix[,"INDEX"] %in% parameters[,1]) & (constMatrix[,"GROUP"] == 0))
+	constMatrix<-constMatrix[logiConstrain,]
+	constMatrix <- constMatrix[order(apply(constMatrix[,7:8],1,function(x){x[1]/10+abs(x[2])})),]
+	constMatrix
+    }
+
+    genMat0 <- function(dep){
+	numberOfVariable <- nrow(dep)
+	G <- rbind(diag(numberOfVariable), -1*diag(numberOfVariable))
+	h <- c(dependences[,5], -1*dependences[,6])
+	return(list(G=G,h=h))
+    }
+
+    genMat1 <- function(dep, N){
+
+	## Range <- sapply(list(min,max),function(x){
+	##   x(as.numeric(rownames(dep)))
+	## }) It is more elegant, more general, but slower
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	numberOfVariables <- nrow(dep)
+	G<- -1*diag(numberOfVariables)
+
+	for(i in 1:numberOfVariables){
+	    if(dep[i,4]!=0){
+		G[i,dep[i,4]] <- 1
+	    }
+
+	}
+
+	G<-G[dep[,4]!=0,]
+
+	if(Range[1]==1){
+	    G<-cbind(G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	} else{
+	    if(Range[2]==N){
+		G<-cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G)
+	    } else {
+		G <- cbind(matrix(ncol=(Range[1]-1),nrow=nrow(G),data=0),G,matrix(ncol=(N-Range[2]),nrow=nrow(G),data=0))
+	    }
+	}
+	return(list(G=G,h=rep(0,nrow(G))))
+    }
+
+    genMat2 <- function(dep, N){
+	G <- rep(1,nrow(dep))
+
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	if(Range[1]==1){
+	    G<-c(G, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		G<-c(numeric(Range[1]-1), G)
+	    } else {
+		G <- c(numeric(Range[1]-1), G, numeric(N-Range[2]))
+	    }
+	}
+
+	G <- t(matrix(sign(dep[2,4])*G))
+	h <- abs(dep[1,4])
+
+	return(list(G=G,h=h))
+    }
+
+    genMat3 <- function(dep, N){
+	Range <- (function(x){
+	    c(min(x), max(x))
+	})(as.numeric(dep[,"rowIndex"]))
+
+	E <- rep(1,nrow(dep))
+
+	if(Range[1]==1){
+	    E<-c(E, numeric(N-Range[2]))
+	} else{
+	    if(Range[2]==N){
+		E<-c(numeric(Range[1]-1), E)
+	    } else {
+		E <- c(numeric(Range[1]-1), E, numeric(N-Range[2]))
+	    }
+	}
+
+
+	E <- t(matrix(E))
+	f <- dep[1,4]
+	return(list(E=E,f=f))
+    }
+
+
+    applyRandTypeG <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(type,
+	       invisible(Gh <- genMat1(dep, N)),
+	       invisible(Gh <- genMat2(dep, N)))
+	return(Gh)
+    }
+
+    applyRandTypeE <- function(dep,N){
+	type <- unique(dep[,"TYPE"])
+	minR <- min(dep[,"rowIndex"])
+	maxR <- max(dep[,"rowIndex"])
+	switch(-type,
+	       stop("Not implemented yet"),
+	       stop("Not implemented yet"),
+	       invisible(Ef <- genMat3(dep, N)))
+	return(Ef)
+    }
+
+    dependences <- depTableMaker(constMatrix, parameters)
+    dependences <- cbind(dependences,1:nrow(dependences))
+    colnames(dependences)[ncol(dependences)] <- "rowIndex"
+    numberOfVariable <- nrow(dependences)
+    nonZeroDeps<-dependences[dependences[,"TYPE"]!=0,]
+    splitedDeps<- split(nonZeroDeps,nonZeroDeps[,"GROUP"])
+    Gh <- list()
+    Ef <- list()
+
+    for(i in 1:length(splitedDeps)){
+	print(splitedDeps[[i]][1,"TYPE"])
+	if(splitedDeps[[i]][1,"TYPE"]>0){
+	    Gh[[i]]<-applyRandTypeG(splitedDeps[[i]],nrow(dependences))
+	} else {
+	    Ef[[i]] <- applyRandTypeE(splitedDeps[[i]],nrow(dependences))
+	}
+    }
+
+    Gh0<- genMat0(dependences)
+    G <- do.call(rbind,lapply(Gh,function(x){x$G}))
+    G<- rbind(Gh0$G,G)
+    h <- do.call(c,lapply(Gh,function(x){x$h}))
+    h <- c(Gh0$h,h)
+    E <- do.call(rbind,lapply(Ef,function(x){x$E}))
+    f <- do.call(c,lapply(Ef,function(x){x$f}))
+    randVal <- suppressWarnings(limSolve::xsample(G=G,H=h,E=E,F=f,iter = iteration))$X
+    results <- list(INDEX =dependences$INDEX, randVal=randVal)
+    return(results)
+}

RBBGCMuso/R/musoSensi.R

@@ -26,6 +26,7 @@ musoSensi <- function(monteCarloFile = NULL,
                      settings = NULL,
                      inputDir = "./",
                      outLoc = "./calib",
+                     outVars = NULL,
                      iterations = 30,
                      preTag = "mont-",
                      outputType = "moreCsv",
@@ -34,23 +35,23 @@ musoSensi <- function(monteCarloFile = NULL,
                      outputFile = "sensitivity.csv",
                      plotName = "sensitivity.png",
                      plotTitle = "Sensitivity",
-                     skipSpinup = FALSE,
+                     skipSpinup = TRUE,
                      dpi=300){
 
     if(is.null(parameters)){
-        parameters <- tryCatch(read.csv("parameters.csv"), error = function (e) {
+        parameters <- tryCatch(read.csv("parameters.csv",stringsAsFactor=FALSE), error = function (e) {
             stop("You need to specify a path for the parameters.csv, or a matrix.")
         })
     } else {
         if((!is.list(parameters)) & (!is.matrix(parameters))){
-             parameters <- tryCatch(read.csv(parameters), error = function (e){
+             parameters <- tryCatch(read.csv(parameters,stringsAsFactor=FALSE), error = function (e){
                                          stop("Cannot find neither parameters file neither the parameters matrix")
                                      })
         }}
 
     doSensi <- function(M){
         npar <- ncol(M)-1
-        M <- M[which(!is.na(M[,"y"])),]
+        M <- M[which(!is.na(M[,ncol(M)])),]
         y <- M[,(npar+1)]
         M <- M[,colnames(M) %in% parameters[,1]]
         npar <- ncol(M)
@@ -84,21 +85,29 @@ musoSensi <- function(monteCarloFile = NULL,
     
 
     if(is.null(monteCarloFile)){
-        M <- musoMonte(parameters = parameters,
+        M <- musoMont(parameters = parameters,
                       settings = settings,
                       inputDir = inputDir,
                       outLoc = outLoc,
                       iterations = iterations,
                       preTag = preTag,
                       outputType = outputType,
+                      outVars = outVars,
                       fun = fun,
                       varIndex = varIndex,
                       skipSpinup = skipSpinup
                       )
+        yInd <-  grep("mod.", colnames(M))[varIndex]
+        parNames <- grep("mod.",colnames(M), invert=TRUE, value = TRUE)
+        M <- M[,c(grep("mod.", colnames(M),invert=TRUE),yInd)]
+        
         return(doSensi(M))
         
     } else {
         M <- read.csv(monteCarloFile)
+         yInd <-  grep("mod.", colnames(M))[varIndex]
+        parNames <- grep("mod.",colnames(M), invert=TRUE, value = TRUE)
+        M <- M[,c(grep("mod.", colnames(M),invert=TRUE),yInd)]
         return(doSensi(M))        
     }
 }

RBBGCMuso/R/musoSensi.R~

@@ -0,0 +1,108 @@
+#' musoSensi
+#'
+#' This function does regression based sensitivity analysis based on the output of musoMonte. 
+#' @author Roland Hollos
+#' @param monteCarloFile If you run musoMonte function previously, you did not have to rerun the monteCarlo, just provide the preservedEpc.csv file with its path. If you do not set this parameter, musoSensi will fun the musoMonte function to get all of the information.
+#' @param outputFile The filename in which the output of musoSensi function will be saved. It's default value is: "sensitivity.csv"
+#' @param plotName The name of the output barplot. It's default value is: "sensitivity.jpg"
+#' @param settings A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically. 
+#' @param parameters This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param inputDir The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file. 
+#' @param iterations number of the monteCarlo run.
+#' @param preTag It will be the name of the output files. For example preTag-1.csv, pretag-2csv...
+#' @param outputType This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet. 
+#' @param fun If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.
+#' @param varIndex This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.
+#' @param skipSpinup With this parameter, you can turn of the spinup phase after the first spinup. I will decrease the time significantly.
+#' @import dplyr
+#' @import graphics
+#' @import grDevices
+#' @import ggplot2
+#' @export
+
+musoSensi <- function(monteCarloFile = NULL,
+                     parameters = NULL,
+                     settings = NULL,
+                     inputDir = "./",
+                     outLoc = "./calib",
+                     iterations = 30,
+                     preTag = "mont-",
+                     outputType = "moreCsv",
+                     fun = mean,
+                     varIndex = 1,
+                     outputFile = "sensitivity.csv",
+                     plotName = "sensitivity.png",
+                     plotTitle = "Sensitivity",
+                     skipSpinup = TRUE,
+                     dpi=300){
+
+    if(is.null(parameters)){
+        parameters <- tryCatch(read.csv("parameters.csv",stringsAsFactor=FALSE), error = function (e) {
+            stop("You need to specify a path for the parameters.csv, or a matrix.")
+        })
+    } else {
+        if((!is.list(parameters)) & (!is.matrix(parameters))){
+             parameters <- tryCatch(read.csv(parameters,stringsAsFactor=FALSE), error = function (e){
+                                         stop("Cannot find neither parameters file neither the parameters matrix")
+                                     })
+        }}
+
+    doSensi <- function(M){
+        npar <- ncol(M)-1
+        M <- M[which(!is.na(M[,ncol(M)])),]
+        y <- M[,(npar+1)]
+        M <- M[,colnames(M) %in% parameters[,1]]
+        npar <- ncol(M)
+        M <- apply(M[,1:npar],2,function(x){x-mean(x)})
+        varNames<- colnames(M)[1:npar]
+        w <- lm(y~M)$coefficients[-1]
+        Sv <- apply(M,2,var)
+        overalVar <- sum(Sv*w^2)
+        S=numeric(npar)
+        
+        for(i in 1:npar){
+            S[i] <- ((w[i]^2*Sv[i])/(overalVar))*100
+        }
+        
+        S <- round(S,digits=2)
+        names(S)<-varNames
+        write.csv(file = outputFile, x = S)
+
+        sensiPlot <- ggplot(data.frame(name=varNames,y=S/100),aes(x=name,y=y))+
+            geom_bar(stat = 'identity')+
+            theme(axis.text.x = element_text(angle = 45, hjust = 1))+
+            xlab(NULL)+
+            ylab(NULL)+
+            ggtitle("Sensitivity")+
+            scale_y_continuous(labels = scales::percent,limits=c(0,1))
+        print(sensiPlot)
+        ggsave(plotName,dpi=dpi)
+        return(S)
+    }
+
+    
+
+    if(is.null(monteCarloFile)){
+        M <- musoMont(parameters = parameters,
+                      settings = settings,
+                      inputDir = inputDir,
+                      outLoc = outLoc,
+                      iterations = iterations,
+                      preTag = preTag,
+                      outputType = outputType,
+                      fun = fun,
+                      varIndex = varIndex,
+                      skipSpinup = skipSpinup
+                      )
+        yInd <-  grep("mod.", colnames(M))[varIndex]
+        parNames <- grep("mod.",colnames(M), invert=TRUE, value = TRUE)
+        M <- M[,c(grep("mod.", colnames(M),invert=TRUE),yInd)]
+        
+        return(doSensi(M))
+        
+    } else {
+        M <- read.csv(monteCarloFile)
+        return(doSensi(M))        
+    }
+}

RBBGCMuso/R/plotMuso.R

@@ -32,8 +32,11 @@ plotMuso <- function(settings=NULL,
                      ##plotname,
                      timee="d",
                      silent=TRUE,
+                     calibrationPar=NULL,
+                     parameters=NULL,
                      debugging=FALSE,
                      keepEpc=FALSE,
+                     fileToChange="epc",
                      logfilename=NULL,
                      aggressive=FALSE,
                      leapYear=FALSE,
@@ -68,7 +71,7 @@ plotMuso <- function(settings=NULL,
     
     groupByTimeFrame <- function(data, timeFrame, groupFun){
         datas <- data %>%
-            group_by(year) %>%
+            group_by(eval(parse(text=timeFrame))) %>%
             summarize(variable=groupFun(eval(parse(text=variable))))
         datas[,1]<-as.numeric(unlist(datas[,1]))
         colnames(datas) <- c("date",variable)
@@ -97,9 +100,9 @@ plotMuso <- function(settings=NULL,
                 tidyr::separate(date2,c("day","month","year"),sep="\\.")
             if(timeFrame!="day"){
                 musoData <- tryCatch(groupByTimeFrame(data=musoData, timeFrame = timeFrame, groupFun = groupFun),
-                                     error=function(e){stop("The timeframe or the gropFun is not found")})
+                                     error=function(e){stop("The timeFrame or the gropFun is not found")})
             }} else {
-                 musoData <- calibMuso(settings,silent = TRUE,skipSpinup=skipSpinup) %>%
+                 musoData <- calibMuso(settings,silent = TRUE,skipSpinup=skipSpinup,parameters = parameters, calibrationPar = calibrationPar,fileToChange = fileToChange) %>%
                      as.data.frame() %>%
                      tibble::rownames_to_column("date") %>%
                      mutate(date2=date,date=as.Date(date,"%d.%m.%Y"),
@@ -236,6 +239,8 @@ plotMuso <- function(settings=NULL,
 #' @param variable The name of the output variable to plot
 #' @param NACHAR This is not implemented yet
 #' @param csvFile The file of the measurement. It must contain a header.
+#' @param calibrationPar documentation in setupMuso()
+#' @param parameters documentation in calibMuso() 
 #' @usage plotMuso(settings, variable,
 #' timee="d", silent=TRUE,
 #' debugging=FALSE, keepEpc=FALSE,
@@ -243,7 +248,7 @@ plotMuso <- function(settings=NULL,
 #' leapYear=FALSE, export=FALSE)
 #' @import ggplot2
 #' @export
-plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",", savePlot=NULL,colour=c("black","blue")){
+plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",", savePlot=NULL,colour=c("black","blue"), calibrationPar=NULL, parameters=NULL){
     if(!is.na(NACHAR)){
         warning("NACHAR is not implemented yet")
     }
@@ -282,3 +287,29 @@ plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",
     
 }
 
+#' compareMuso 
+#'
+#' This function runs the modell, change one of it's input, and plot both in one plot. 
+#' 
+#' @author Roland Hollos
+#' @param settings You have to run the setupMuso function before rungetMuso. It is its output which contains all of the necessary system variables. It sets the whole environment
+#' @param parameters In the settings variable you have set the row indexes of the variables, you wish to change. In this parameter you can give an exact value for them in a vector like: c(1,2,3,4)
+#' @param variable The name of the output variable to plot
+#' @param calibrationPar  in the help of setupMuso function.
+#' @param fileToChange You can change any line of the epc or the ini file, you just have to specify with this variable which file you van a change. Two options possible: "epc", "ini", "both"
+#' @import ggplot2
+#' @export
+compareMuso <- function(settings=NULL,parameters, variable=1, calibrationPar=NULL, fileToChange="epc", skipSpinup=TRUE, timeFrame="day"){
+
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+
+    p1 <- plotMuso(settings = settings,variable = variable,timeFrame = timeFrame)
+    p2 <- p1+plotMuso(settings = settings,variable = variable, timeFrame = timeFrame,fileToChange=fileToChange,layerPlot=TRUE)
+    p2
+    
+}
+
+

RBBGCMuso/R/plotMuso.R~

@@ -0,0 +1,284 @@
+#'plot the BBGCMuso output 
+#'
+#' This function runs the BBGC-MuSo model and reads in its outputfile in a very structured way, and after that plot the results automaticly 
+#' 
+#' @author Roland Hollos, Dora Hidy
+#' @param settings You have to run the setupMuso function before rungetMuso. It is its output which contains all of the necessary system variables. It sets the whole environment
+#' @param timee The required timesteps in the modell output. It can be "d", if it is daily, "m", if it's monthly, "y", it it is yearly
+#' @param debugging If it is TRUE, it copies the log file to a Log directory to store it, if it is stamplog it contatenate a number before the logfile, which is one more than the maximum of the represented ones in the LOG directory. If it is true or stamplog it collects the "wrong" logfiles
+#' @param keepEpc If TRUE, it keeps the epc file and stamp it, after these copies it to the EPCS directory. If debugging True or false, it copies the wrong epc files to the wrong epc directory.
+#' @param export if it is yes or you give a filename here, it converts the output to the specific extension. For example, if you set export to "example.csv", it converts the output to "csv", if you set it to "example.xls" it converts to example.xls with the xlsx package. If it is not installed it gives back a warning message and converts it to csv.
+#' @param silent If you set it TRUE all off the modells output to the screen will be suppressed. It can be usefull, because it increases the model-speed.
+#' @param aggressive It deletes every possible modell-outputs from the previous modell runs.
+#' @param variable column number of the variable which should be plottedor "all" if you have less than 10 variables. In this case it will plot everything in a matrix layout 
+#' @param leapYear Should the function do a leapyear correction on the outputdata? If TRUE, then the 31.12 day will be doubled.
+#' @param logfilename If you want to set a specific name for your logfiles you can set this via logfile parameter
+#' @param plotType There are two options now: continious time series("cts")  or disctrete time series("dts")
+#' @param skipSpinup If TRUE, calibMuso wont do spinup simulation
+#' @return It depends on the export parameter. The function returns with a matrix with the modell output, or writes this in a file, which is given previously
+#' @usage plotMuso(settings, variable,
+#' timee="d", silent=TRUE,
+#' debugging=FALSE, keepEpc=FALSE,
+#' logfilename=NULL, aggressive=FALSE,
+#' leapYear=FALSE, export=FALSE)
+#' @import ggplot2
+#' @import dplyr
+#' @import tidyr
+#' @export
+
+plotMuso <- function(settings=NULL,
+                     variable=1,
+                     ##compare,
+                     ##plotname,
+                     timee="d",
+                     silent=TRUE,
+                     debugging=FALSE,
+                     keepEpc=FALSE,
+                     logfilename=NULL,
+                     aggressive=FALSE,
+                     leapYear=FALSE,
+                     plotName=NULL,
+                     plotType="cts",
+                     layerPlot=FALSE,
+                     colour="blue",
+                     skipSpinup=TRUE,
+                     fromData=FALSE,
+                     timeFrame="day",
+                     groupFun=mean,
+                     dpi=300){
+
+    if( plotType!="cts" && plotType != "dts"){
+        warning(paste0("The plotType ", plotType," is not implemented, plotType is set to cts"))
+        plotType <- "cts"
+    }
+    
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    numberOfYears <- settings$numYears
+    startYear <- settings$startYear
+    ## musoData <- rungetMuso(settings=settings,
+    ##                        silent=silent,
+    ##                        timee=timee,
+    ##                        debugging=debugging,
+    ##                        keepEpc=keepEpc,
+    ##                        logfilename=logfilename,
+    ##                        export=export)
+    
+    groupByTimeFrame <- function(data, timeFrame, groupFun){
+        datas <- data %>%
+            group_by(year) %>%
+            summarize(variable=groupFun(eval(parse(text=variable))))
+        datas[,1]<-as.numeric(unlist(datas[,1]))
+        colnames(datas) <- c("date",variable)
+        datas
+    }
+    
+    if(fromData){
+        Reva <- tryCatch(getdailyout(settings), #(:INSIDE: getOutput.R )
+                                    error = function (e){
+                                        setwd((whereAmI))
+                                        stop("Cannot read binary output, please check if the output type is set 2 in the ini files!")})
+        colnames(Reva) <- unlist(settings$outputVars[[1]])
+        rownames(Reva) <- NULL
+        musoData <- cbind(musoDate(startYear = startYear,numYears = numberOfYears,combined = TRUE,corrigated=FALSE),
+              rep(1:365,numberOfYears),
+              musoDate(startYear = startYear,numYears = numberOfYears,combined = FALSE,corrigated=FALSE),as.data.frame(Reva))
+        colnames(musoData)[1:5]<-c("date","yearDay","year","day","month")
+        musoData <-musoData %>%
+            mutate(date=as.Date(as.character(date),"%d.%m.%Y"))
+    } else {
+        if(!is.element("cum_yieldC_HRV",unlist(settings$outputVars[[1]]))){
+            musoData <- calibMuso(settings,silent = TRUE,skipSpinup=skipSpinup) %>%
+                as.data.frame() %>%
+                tibble::rownames_to_column("date") %>%
+                mutate(date2=date,date=as.Date(date,"%d.%m.%Y")) %>%
+                tidyr::separate(date2,c("day","month","year"),sep="\\.")
+            if(timeFrame!="day"){
+                musoData <- tryCatch(groupByTimeFrame(data=musoData, timeFrame = timeFrame, groupFun = groupFun),
+                                     error=function(e){stop("The timeFrame or the gropFun is not found")})
+            }} else {
+                 musoData <- calibMuso(settings,silent = TRUE,skipSpinup=skipSpinup) %>%
+                     as.data.frame() %>%
+                     tibble::rownames_to_column("date") %>%
+                     mutate(date2=date,date=as.Date(date,"%d.%m.%Y"),
+                            yearDay=rep(1:365,numberOfYears), cum_yieldC_HRV=cum_yieldC_HRV*22.22) %>%
+                     tidyr::separate(date2,c("day","month","year"),sep="\\.")
+                 if(timeFrame!="day"){
+                     musoData <- tryCatch(groupByTimeFrame(data=musoData, timeFrame = timeFrame, groupFun = groupFun),
+                                          error=function(e){stop("The timeframe or the gropFun is not found")})
+                 }
+                 
+             }
+    }
+
+    ## numVari <- ncol(musoData)
+     numVari <- ncol(musoData)-5
+
+    pointOrLineOrPlot <- function(musoData, variableName, plotType="cts", expandPlot=FALSE, plotName=NULL){
+        if(!expandPlot){
+            if(plotType=="cts"){
+                if(length(variableName)==1){
+                   p <- ggplot(musoData,aes_string("date",variableName))+geom_line(colour=colour)+theme(axis.title.x=element_blank())
+                   if(!is.null(plotName)){
+                       ggsave(as.character(plotName), plot = p)
+                    p
+                   }
+                   p
+                } else {
+                    p <- musoData %>%
+                        select(c("date", variableName))%>%
+                        gather(., key= outputs, value = bla, variableName) %>%
+                                        # head  %>%
+                        ggplot(aes(x=date,y=bla))+
+                        facet_wrap(~ outputs, scales = "free_y",ncol=1) +
+                        geom_line(colour=colour)+
+                        theme(
+                            axis.title.y = element_blank()
+                        )
+                    if(!is.null(plotName)){
+                        ggsave(as.character(plotName), plot = p)
+                    }
+                    p
+                }
+            } else {
+                if(length(variableName)==1){
+                    p <- ggplot(musoData,aes_string("date",variableName))+geom_point(colour=colour)+theme(axis.title.x=element_blank())
+                    if(!is.null(plotName)){
+                        ggsave(as.character(plotName),p)
+                    }
+                    p
+                } else{
+                    p <- musoData %>%
+                        select(c("date",variableName))%>%
+                        gather(., key= outputs, value = bla,variableName) %>%
+                                        # head  %>%
+                        ggplot(aes(x=date,y=bla))+
+                        facet_wrap(~ outputs, scales = "free_y",ncol=1) +
+                        geom_line(colour=colour)+
+                        theme(
+                            axis.title.y = element_blank()
+                        )
+                    if(!is.null(plotName)){
+                        ggsave(as.character(plotName),p)
+                    }
+                    p
+                }
+            }
+        } else {
+            if(!is.null(plotName)){
+                stop("Cannot save a single plot layer to a graphics device")
+            }
+            
+            if(plotType=="cts"){
+                if(length(variableName)==1){
+                    geom_line(data=musoData, colour=colour, aes_string("date",variableName))
+                    
+                } else {
+                    stop("you cannot add layers for multiple plots")
+                }
+            } else {
+                if(length(variableName)==1){
+                    geom_point(data=musoData, colour=colour, aes_string("date",variableName))
+                } else{
+                    stop("you cannot add layers for multiple plots")
+                }
+            }
+            
+        }
+        
+    }
+    
+
+
+    variableName <-  as.character(settings$outputVars[[1]])[variable]
+    if(is.character(variable)){
+        if(identical(variable,"all")){
+            variableName <- as.character(settings$outputVars[[1]])
+        } else {
+            if(identical(character(0),setdiff(variable,as.character(settings$outputVars[[1]])))){
+                variableName <- variable
+            } else {
+                stop("The symmetric difference of the set of the output variables specified in the ini files and the set specified with your variable parameter is not the empty set.")
+            }   
+        }
+        
+        if(length(variableName)>8){
+            warning("Too many variables to plot, the output quality can be poor")
+        }
+        
+     } else {
+         
+         if(prod(sapply(variable,function(x){
+             return(x >= 0 && x <= numVari)
+         }))){
+             variableName <-  as.character(settings$outputVars[[1]])[variable]
+         } else {
+             stop("Not all members of the variable parameter are among the output variables")
+         }}
+    
+    pointOrLineOrPlot(musoData = musoData,
+                      variableName = variableName,
+                      plotType = plotType,
+                      expandPlot = layerPlot,
+                      plotName = plotName)
+     }
+
+#'plot the BBGCMuso output with data 
+#'
+#' This function runs the BBGC-MuSo model and reads in its outputfile in a very structured way, and after that plot the results automaticly along with a given measurement 
+#' 
+#' @author Roland Hollos, Dora Hidy
+#' @param settings You have to run the setupMuso function before rungetMuso. It is its output which contains all of the necessary system variables. It sets the whole environment
+#' @param sep This is the separator used in the measurement file
+#' @param savePlot It it is specified, the plot will be saved in a format specified with the immanent extension
+#' @param variable The name of the output variable to plot
+#' @param NACHAR This is not implemented yet
+#' @param csvFile The file of the measurement. It must contain a header.
+#' @usage plotMuso(settings, variable,
+#' timee="d", silent=TRUE,
+#' debugging=FALSE, keepEpc=FALSE,
+#' logfilename=NULL, aggressive=FALSE,
+#' leapYear=FALSE, export=FALSE)
+#' @import ggplot2
+#' @export
+plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",", savePlot=NULL,colour=c("black","blue")){
+    if(!is.na(NACHAR)){
+        warning("NACHAR is not implemented yet")
+    }
+    if(is.null(settings)){
+        settings <- setupMuso()
+    }
+    
+    numberOfYears <- settings$numYears
+    startYear <- settings$startYear
+    yearVec <- seq(from = startYear, length=numberOfYears,by=1)
+
+    
+    data <- read.table(csvFile,header = TRUE, sep = ",") %>%
+        select(variable)
+    
+    baseData <- calibMuso(settings,silent=TRUE) %>%
+        as.data.frame() %>%
+        tibble::rownames_to_column("date") %>%
+        dplyr::mutate(date2=date,date=as.Date(date,"%d.%m.%Y"),yearDay=rep(1:365,numberOfYears)) %>%
+        tidyr::separate(date2,c("day","month","year"),sep="\\.")
+    baseData <- cbind(baseData,data)
+    colnames(baseData)[ncol(baseData)] <- "measuredData"
+
+    p <- baseData %>%
+        ggplot(aes_string("date",variable)) +
+        geom_line(colour=colour[1]) +
+        geom_point(colour=colour[2], aes(date,measuredData)) +
+        labs(y = paste0(variable,"_measured"))+
+        theme(axis.title.x = element_blank())
+    if(!is.null(savePlot)){
+        ggsave(savePlot,p)
+        return(p)
+    } else {
+        return(p)
+    }
+    
+}
+

RBBGCMuso/R/putOutVars.R

@@ -0,0 +1,32 @@
+#' putOutVars
+#'
+#'This function is for adding variables in the inifiles.
+#'
+#' @author Roland Hollos
+#' @param IniFile The name of the normal ini file.
+#' @param outputVars List of the output codes
+#' @keywords internal 
+putOutVars <- function(iniFile,outputVars,modifyOriginal = FALSE){
+    ini <- readLines(iniFile)
+    numVarsOriginal <- as.numeric(ini[grep("DAILY_OUTPUT",ini)+1])
+    if(!modifyOriginal){
+        iniOut <- paste0(tools::file_path_sans_ext(basename(iniFile)),"-tmp.",tools::file_ext(iniFile))
+    } else {
+        iniOut <- iniFile
+    }
+    
+    outNames <- sapply(outputVars,musoMapping)
+    partOne <- ini[1:grep("DAILY_OUTPUT",ini)]
+    partTwo <- ini[grep("ANNUAL_OUTPUT",ini):(length(ini))]
+    numVars <- length(outputVars)
+    fileContent <- c(partOne,
+                     as.character(numVars),
+                     sapply(outputVars,function (x) {
+                         paste(as.character(x),musoMapping(x),sep = " ")
+                     }),
+                     "",
+                     partTwo)
+    writeLines(fileContent,iniOut)
+    return(list(names=outNames,ratio=numVars/numVarsOriginal))
+}
+

RBBGCMuso/R/putOutVars.R~

@@ -0,0 +1,13 @@
+#' putOutVars
+#'
+#'This function is for adding variables in the inifiles.
+#'
+#' @author Roland Hollos
+#' @param IniFile The name of the normal ini file.
+#' @param outputVars List of the output codes
+#' @keywords internal 
+putOutVars <- function(IniFile,outputVars){
+    ini <- readLines(IniFile)
+    sapply(outputVars,musoMapping)
+}
+

RBBGCMuso/R/setupMuso.R

@@ -21,6 +21,7 @@
 #' @param nitInput Via this parameter, you can tell the modell where are the NO2 data files. As default it reads this from the iniFiles.
 #' @param iniInput Via this parameter, you can tell the modell where are the ini files. As default it reads this from the iniFiles.
 #' @param epcInput Via this parameter, you can tell the modell where are the epc data files. As default it reads this from the iniFiles.
+#' @param modelOutputs This parameter is the list of the codes of outpu 
 #' @usage setupMuso(executable=NULL, parallel = F, calibrationPar =c(1),
 #' outputLoc=NULL, inputLoc=NULL,
 #' metInput=NULL, CO2Input=NULL,
@@ -30,12 +31,13 @@
 #' fertInput=NULL, irrInput=NULL,
 #' nitInput=NULL, iniInput=NULL, epcInput=NULL)
 #' @return The output is a the modell setting list wich contains the following elements:
-#' executable, calibrationPar, outputLoc, outputName, inputLoc, iniInput, metInput, epcInput,thinInput,CO2Input, mowInput, grazInput, harvInput, plougInput, fertInput, irrInput, nitInput, inputFiles, numData, startyear, numYears, outputVars
+#' executable, calibrationPar, outputLoc, outputName, inputLoc, iniInput, metInput, epcInput,thinInput,CO2Input, mowInput, grazInput, harvInput, plougInput, fertInput,rrInput, nitInput, inputFiles, numData, startyear, numYears, outputVars
 #' @export
 setupMuso <- function(executable=NULL,
                       parallel = F,
                       calibrationPar =c(1),
                       outputLoc=NULL,
+                      modelOutputs=NULL,
                       inputLoc=NULL,
                       metInput=NULL,
                       CO2Input=NULL,
@@ -195,7 +197,7 @@ setupMuso <- function(executable=NULL,
     ##THIS IS AN UGLY SOLUTION, WHICH NEEDS AN UPGRADE!!! FiXED (2017.09.11)
     ## outputName <- unlist(strsplit(grep("prefix for output files",iniFiles[[2]],value=TRUE),"[\ \t]"))[1]
     if(is.null(outputName)){
-        stop("I cannot find outputName in your default ini file \n Please make sure that the line wich contains the name also contains the prefix and the outmut keywords!")
+        stop("I cannot find outputName in your default ini file \n Please make sure that the line wich contains the name also contains the prefix and the output keywords!")
         
     }
     ##    outputName<-unlist(read.table(iniInput[2],skip=93,nrows = 1))[1]
@@ -216,12 +218,12 @@ setupMuso <- function(executable=NULL,
     
     inputFiles<-c(iniInput,epcInput,metInput)
     numData<-rep(NA,3)
-    numYears <-  as.numeric(unlist(strsplit(grep("simulation years",iniFiles[[2]],value=TRUE),"[\ \t]"))[1])
+    numYears <-  as.numeric(unlist(strsplit(iniFiles[[2]][grep("TIME_DEFINE",iniFiles[[2]])+2],"[\ \t]"))[1])
     ##    numYears<-unlist(read.table(iniInput[2],skip = 14,nrows = 1)[1])
     numValues <-  as.numeric(unlist(strsplit(iniFiles[[2]][grep("DAILY_OUTPUT",iniFiles[[2]])+1],"[\ \t]"))[1])
     ## numValues will be replaced to numVar
     ## numValues<-unlist(read.table(iniInput[2],skip=102,nrows = 1)[1])
-    startYear <- as.numeric(unlist(strsplit(grep("first simulation year",iniFiles[[2]],value=TRUE),"[\ \t]"))[1])
+    startYear <- as.numeric(unlist(strsplit(iniFiles[[2]][grep("TIME_DEFINE",iniFiles[[2]])+3],"[\ \t]"))[1])
     numData[1] <- numValues * sumDaysOfPeriod(startYear,numYears,corrigated=leapYear)
 
     numData[2] <- numYears * numValues*12
@@ -232,6 +234,13 @@ setupMuso <- function(executable=NULL,
     writeLines(iniFiles[[1]],iniInput[1])
     writeLines(iniFiles[[2]],iniInput[2])
 
+    if(!is.null(modelOutputs)){
+        outVarChanges <- putOutVars(iniFile = iniInput[2],outputVars = modelOutputs, modifyOriginal = TRUE)
+        numData <- round(numDate*outVarChanges[[2]])
+        outputVars[[1]] <-outVarChanges[[1]] 
+    }
+    
+
     suppressWarnings(file.remove(paste0(file.path(outputLoc,outputName[1]),".log")))
     ## I use file.path additionally because We do not know if outputLoc ends or not to "/"
     suppressWarnings(file.remove(paste0(file.path(outputLoc,outputName[2]),".log")))

RBBGCMuso/R/setupMuso.R~

@@ -0,0 +1,274 @@
+#' setupMuso
+#'
+#' This funcion is fundamental for the BiomBGC-MuSo modell related functions like spinupMuso, normalMuso, rungetMuso, because it sets the modells environment.
+#'
+#' @author Roland Hollos
+#' @param parallel Do you want to run multiple modell paralelly, if yes, set this variable to TRUE
+#' @param executable This parameter stores the place of the modell-executable file. In normal usage, you don't have to be set this, because a RBBgcmuso package contains allways the latest modell executable. In spite of this, if you would like to use this package for modell development or just want to use different models (for example for comparison), you will find it useful
+#' @param calibrationPar You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)
+#' @param outputLoc Where should the modell puts its outputs. You should give a location for it via this variable, for example: outputLoc="/place/of/the/outputs/"
+#' @param inputLoc Usually it is the root directory, where you put the iniFiles for the modell
+#' @param metInput Via metInput parameter, you can tell the modell where are the meteorological files. As default it reads this from the iniFiles.
+#' @param CO2Input Via CO2 parameter, you can tell the modell where are the CO2 data files. As default it reads this from the iniFiles.
+#' @param plantInput Via planting parameter, you can tell the modell where are the data files, which contains the planting informations. As default it reads this from the iniFiles.
+#' @param thinInput Via thining parameter, you can tell the modell where are the data files, which contains the thining informations. As default it reads this from the iniFiles.
+#' @param mowInput Via mowing parameter, you can tell the modell where are the data files, which contains the mowing informations. As default it reads this from the iniFiles.
+#' @param grazInput Via grazing parameter, you can tell the modell where are the data files, which contains the grazing informations. As default it reads this from the iniFiles.
+#' @param harvInput Via harvesting parameter, you can tell the modell where are the data files, which contains the harvesting informations. As default it reads this from the iniFiles.
+#' @param plougInput Via ploughing parameter, you can tell the modell where are the data files, which contains the ploughing informations. As default it reads this from the iniFiles.
+#' @param fertInput Via fertilizing parameter, you can tell the modell where are the fertilizing data files, which contains the fertilizing informations. As default it reads this from the iniFiles.
+#' @param irrInput Via irrigation parameter, you can tell the modell where are the data files, which contains the irrigation informations. As default it reads this from the iniFiles.
+#' @param nitInput Via this parameter, you can tell the modell where are the NO2 data files. As default it reads this from the iniFiles.
+#' @param iniInput Via this parameter, you can tell the modell where are the ini files. As default it reads this from the iniFiles.
+#' @param epcInput Via this parameter, you can tell the modell where are the epc data files. As default it reads this from the iniFiles.
+#' @param modelOutputs This parameter is the list of the codes of outpu 
+#' @usage setupMuso(executable=NULL, parallel = F, calibrationPar =c(1),
+#' outputLoc=NULL, inputLoc=NULL,
+#' metInput=NULL, CO2Input=NULL,
+#' plantInput=NULL, thinInput=NULL,
+#' mowInput=NULL, grazInput=NULL,
+#' harvInput=NULL, plougInput=NULL,
+#' fertInput=NULL, irrInput=NULL,
+#' nitInput=NULL, iniInput=NULL, epcInput=NULL)
+#' @return The output is a the modell setting list wich contains the following elements:
+#' executable, calibrationPar, outputLoc, outputName, inputLoc, iniInput, metInput, epcInput,thinInput,CO2Input, mowInput, grazInput, harvInput, plougInput, fertInput,rrInput, nitInput, inputFiles, numData, startyear, numYears, outputVars
+#' @export
+setupMuso <- function(executable=NULL,
+                      parallel = F,
+                      calibrationPar =c(1),
+                      outputLoc=NULL,
+                      modelOutputs=NULL,
+                      inputLoc=NULL,
+                      metInput=NULL,
+                      CO2Input=NULL,
+                      plantInput=NULL,
+                      thinInput=NULL,
+                      mowInput=NULL,
+                      grazInput=NULL,
+                      harvInput=NULL,
+                      plougInput=NULL,
+                      fertInput=NULL,
+                      irrInput=NULL,
+                      nitInput=NULL,
+                      iniInput=NULL,
+                      epcInput=NULL,
+                      mapData=NULL,
+                      leapYear=FALSE,
+                      version=5
+                      ){
+
+    Linuxp <-(Sys.info()[1]=="Linux")
+    writep <- 0
+
+    if(is.null(mapData)&version==4){
+        mData <- mMapping4
+    }
+    
+    inputParser <- function(string,fileName,counter,value=TRUE){
+        unlist(strsplit(grep(string,fileName,value=TRUE),"[\ \t]"))[counter]
+    }
+
+    outMaker <- function(inputVar,grepString,filep){
+        tempVar <- eval(parse(text=inputVar))
+        if(is.null(tempVar)){
+            writep <<- writep+1
+            if(filep)
+            {
+                tempVar["spinup"] <- file.path(inputLoc,inputParser(string=grepString,fileName=iniFiles$spinup,counter=1,value=TRUE))
+                tempVar["normal"] <- file.path(inputLoc,inputParser(string=grepString,fileName=iniFiles$normal,counter=1,value=TRUE))
+            } else {
+                tempVar["spinup"] <- inputParser(string=grepString,fileName=iniFiles$spinup,counter=1,value=TRUE)
+                tempVar["normal"] <- inputParser(string=grepString,fileName=iniFiles$normal,counter=1,value=TRUE)
+                
+            }
+
+        } else {
+            iniFiles$spinup[grep(grepString,iniFiles$spinup)] <<- paste0(tempVar[1],"\t ",grepString)
+            
+            if(length(tempVar)==2){
+                iniFiles$normal[grep(" grepString",iniFiles$normal)] <<- paste0(tempVar[2],"\t ",grepString)
+            }
+        }
+        return(tempVar)
+    }
+    
+    if(is.null(inputLoc)){
+        inputLoc<- normalizePath("./")
+    } else{
+        inputLoc <- normalizePath(inputLoc)
+    }
+    
+    #iniChangedp <- FALSE
+
+    if(is.null(iniInput)){
+        spinups<-grep("s.ini$",list.files(inputLoc),value=TRUE)
+        normals<-grep("n.ini$",list.files(inputLoc),value=TRUE)
+
+        if(length(spinups)==1){
+            iniInput[1]<-file.path(inputLoc,spinups)
+        } else {stop("There are multiple or no spinup ini files, please choose")}
+        
+
+        if(length(normals)==1){
+            iniInput[2]<-file.path(inputLoc,normals)
+        } else {stop("There are multiple or no normal ini files, please choose")}
+
+    }
+
+    ##read the ini files for the further changes
+
+    iniFiles<-lapply(iniInput, function (x) readLines(x,-1))
+    iniFiles[[1]] <- gsub("\\","/", iniFiles[[1]],fixed=TRUE) #replacing \ to /
+    iniFiles[[2]] <- gsub("\\","/", iniFiles[[2]],fixed=TRUE) #replacing \ to /
+    names(iniFiles) <- c("spinup","normal")
+
+
+    
+    inputs <- lapply(1:nrow(grepHelper), function (x) {
+
+        outMaker(grepHelper[x,1],grepHelper[x,2],grepHelper[x,3])
+        
+    })   
+    names(inputs) <- grepHelper$inputVar   
+    ## grepHelper is in sysdata.rda it is a table like this:
+    ##
+    ## inputVar         string      isFile  
+    ## 1    epcInput  EPC file name TRUE 
+    ## 2    metInput  met file name TRUE
+    ## 3    CO2Input       CO2 file TRUE
+    ## 4    nitInput     N-dep file TRUE
+    ## 5   thinInput    do THINNING FALSE
+    ## 6  plantInput    do PLANTING FALSE
+    ## 7    mowInput      do MOWING FALSE
+    ## 8   grazInput     do GRAZING FALSE
+    ## 9   harvInput  do HARVESTING FALSE
+    ## 10 plougInput   do PLOUGHING FALSE
+    ## 11  fertInput do FERTILIZING FALSE
+    ## 12   irrInput  do IRRIGATION FALSE
+   # return(inputs) debug element 
+
+    
+    if(is.null(mapData)){
+        
+        outIndex<-grep("DAILY_OUTPUT",iniFiles[[2]])+1
+        numVar<-as.numeric(unlist(strsplit(iniFiles[[2]][outIndex],"[\ \t]"))[1])
+        dailyVarCodes<-tryCatch(iniFiles[[2]][(outIndex+1):(outIndex+numVar)],
+                                error = function(e){
+        stop("Cannot read indexes of output variables from the normal ini file, please make sure you have not skiped a line after the flag: \"DAILY_OUTPUT\"")    
+        })
+        dailyVarnames<-lapply(dailyVarCodes, function(x) musoMapping(unlist(strsplit(x,"[\ \t]"))[1]))
+
+        outIndex<-grep("ANNUAL_OUTPUT",iniFiles[[2]])+1
+        numVar<-as.numeric(unlist(strsplit(iniFiles[[2]][outIndex],"[\ \t]"))[1])
+        annualVarCodes<-iniFiles[[2]][(outIndex+1):(outIndex+numVar)]
+        annualVarnames<-lapply(annualVarCodes, function(x) musoMapping(unlist(strsplit(x,"[\ \t]"))[1]))
+        outputVars<-list(dailyVarnames,annualVarnames)} else {
+
+                                                          c<-grep("DAILY_OUTPUT",iniFiles[[2]])+1
+                                                          numVar<-as.numeric(unlist(strsplit(iniFiles[[2]][c],"[\ \t]"))[1])
+                                                          dailyVarCodes<-iniFiles[[2]][(c+1):(c+numVar)]
+                                                          dailyVarnames<-lapply(dailyVarCodes, function(x) musoMapping(unlist(strsplit(x,"[\ \t]"))[1],mapData))
+
+                                                          c<-grep("ANNUAL_OUTPUT",iniFiles[[2]])+1
+                                                          numVar<-as.numeric(unlist(strsplit(iniFiles[[2]][c],"[\ \t]"))[1])
+                                                          annualVarCodes<-iniFiles[[2]][(c+1):(c+numVar)]
+                                                          annualVarnames<-lapply(annualVarCodes, function(x) musoMapping(unlist(strsplit(x,"[\ \t]"))[1],mapData))
+                                                          outputVars<-list(dailyVarnames,annualVarnames)
+
+
+                                                          
+                                                      }
+    
+    
+    
+    if(is.null(executable)){
+        if(Linuxp){
+            executable<-file.path(inputLoc,"muso")
+        } else {
+            executable<-file.path(inputLoc,"muso.exe")
+        }
+    } else {        
+        file.copy(executable,inputLoc)
+    }
+    outputName <- character(2)
+    outputName[1] <- basename(unlist(strsplit(iniFiles[[1]][grep("OUTPUT_CONTROL",iniFiles[[1]])+1],"[\ \t]"))[1])
+    outputName[2] <- basename(unlist(strsplit(iniFiles[[2]][grep("OUTPUT_CONTROL",iniFiles[[2]])+1],"[\ \t]"))[1])
+    ##  outputName <- unlist(strsplit(grep("output",grep("prefix",iniFiles[[2]],value=TRUE),value=TRUE),"[\ \t]"))[1]
+    ##THIS IS AN UGLY SOLUTION, WHICH NEEDS AN UPGRADE!!! FiXED (2017.09.11)
+    ## outputName <- unlist(strsplit(grep("prefix for output files",iniFiles[[2]],value=TRUE),"[\ \t]"))[1]
+    if(is.null(outputName)){
+        stop("I cannot find outputName in your default ini file \n Please make sure that the line wich contains the name also contains the prefix and the output keywords!")
+        
+    }
+    ##    outputName<-unlist(read.table(iniInput[2],skip=93,nrows = 1))[1]
+
+
+    if(is.null(outputLoc)){
+        ##  outputLoc<-paste((rev(rev(unlist(strsplit(outputName,"/")))[-1])),collapse="/")
+        outputLoc <- dirname(unlist(strsplit(iniFiles[[2]][grep("OUTPUT_CONTROL",iniFiles[[2]])+1],"[\ \t]"))[1])
+        if(substr(outputLoc,start = 1,stop = 1)!="/"){
+            ##if the outputName is not absolute path make it absolute
+            outputLoc <- file.path(inputLoc,outputLoc)
+        } 
+    } else {
+        outputLoc <- normalizePath(outputLoc)
+    }
+
+    
+    
+    inputFiles<-c(iniInput,epcInput,metInput)
+    numData<-rep(NA,3)
+    numYears <-  as.numeric(unlist(strsplit(iniFiles[[2]][grep("TIME_DEFINE",iniFiles[[2]])+2],"[\ \t]"))[1])
+    ##    numYears<-unlist(read.table(iniInput[2],skip = 14,nrows = 1)[1])
+    numValues <-  as.numeric(unlist(strsplit(iniFiles[[2]][grep("DAILY_OUTPUT",iniFiles[[2]])+1],"[\ \t]"))[1])
+    ## numValues will be replaced to numVar
+    ## numValues<-unlist(read.table(iniInput[2],skip=102,nrows = 1)[1])
+    startYear <- as.numeric(unlist(strsplit(iniFiles[[2]][grep("TIME_DEFINE",iniFiles[[2]])+3],"[\ \t]"))[1])
+    numData[1] <- numValues * sumDaysOfPeriod(startYear,numYears,corrigated=leapYear)
+
+    numData[2] <- numYears * numValues*12
+    numData[3] <- numYears * numValues
+
+    ##Writing out changed ini-file
+
+    writeLines(iniFiles[[1]],iniInput[1])
+    writeLines(iniFiles[[2]],iniInput[2])
+
+    suppressWarnings(file.remove(paste0(file.path(outputLoc,outputName[1]),".log")))
+    ## I use file.path additionally because We do not know if outputLoc ends or not to "/"
+    suppressWarnings(file.remove(paste0(file.path(outputLoc,outputName[2]),".log")))
+    
+    settings = list(executable = executable,
+                    calibrationPar = calibrationPar,
+                    outputLoc=outputLoc,
+                    outputNames=outputName,
+                    inputLoc=inputLoc,
+                    iniInput=iniInput,
+                    metInput=inputs$metInput,
+                    epcInput=inputs$epcInput,
+                    thinInput=inputs$thinInput,
+                    CO2Input=inputs$CO2Input,
+                    mowInput=inputs$mowInput,
+                    grazInput=inputs$grazInput,
+                    harvInput=inputs$harvInput,
+                    plougInput=inputs$plougInput,
+                    fertInput=inputs$fertInput,
+                    irrInput=inputs$irrInput,
+                    nitInput=inputs$nitInput,
+                    inputFiles=inputFiles,
+                    numData=numData,
+                    startYear=startYear,
+                    numYears=numYears,
+                    outputVars=outputVars
+                    )
+    
+    if(writep!=nrow(grepHelper)){
+        writeLines(iniFiles[[1]],iniInput[[1]])
+        if(inputs$epcInput[1]!=inputs$epc$Input[2]){ #Change need here
+            writeLines(iniFiles[[2]],iniInput[[2]])      
+        }
+    }
+    return(settings)
+
+}
+

RBBGCMuso/inst/examples/hegyhatsal/CO2_from1961.txt

@@ -1,55 +0,0 @@
- 1961  317.33
- 1962  318.08
- 1963  318.70
- 1964  319.36
- 1965  320.02
- 1966  321.09
- 1967  321.99
- 1968  322.93
- 1969  324.21
- 1970  325.24
- 1971  326.06
- 1972  327.18
- 1973  328.84
- 1974  329.73
- 1975  330.73
- 1976  331.83
- 1977  333.25
- 1978  334.60
- 1979  336.85
- 1980  338.69
- 1981  339.93
- 1982  341.13
- 1983  342.78
- 1984  344.42
- 1985  345.90
- 1986  347.15
- 1987  348.93
- 1988  351.48
- 1989  352.91
- 1990  354.19
- 1991  355.59
- 1992  356.37
- 1993  357.04
- 1994  358.88
- 1995  360.88
- 1996  362.64
- 1997  363.76
- 1998  366.63
- 1999  368.31
- 2000  369.48
- 2001  372.59
- 2002  374.37
- 2003  378.04
- 2004  380.88
- 2005  383.88
- 2006  385.64
- 2007  385.76
- 2008  386.13
- 2009  387.37
- 2010  389.85
- 2011  391.62
- 2012  393.82
- 2013  396.48
- 2014  398.61
- 2015  400.00

RBBGCMuso/inst/examples/hegyhatsal/CO2_from2007.txt

@@ -1,9 +0,0 @@
- 2007  385.76
- 2008  386.13
- 2009  387.37
- 2010  389.85
- 2011  391.62
- 2012  393.82
- 2013  396.48
- 2014  398.61
- 2015  400.00

RBBGCMuso/inst/examples/hegyhatsal/Ndep_from1961.txt

@@ -1,55 +0,0 @@
-1961    0.0008407
-1962    0.0008601
-1963    0.0008795
-1964    0.0008989
-1965    0.0009183
-1966    0.0009378
-1967    0.0009572
-1968    0.0009766
-1969    0.0009960
-1970    0.0010348
-1971    0.0010591
-1972    0.0010465
-1973    0.0010524
-1974    0.0010582
-1975    0.0010641
-1976    0.0010699
-1977    0.0010758
-1978    0.0010816
-1979    0.0010875
-1980    0.0010933
-1981    0.0010992
-1982    0.0011050
-1983    0.0011109
-1984    0.0011167
-1985    0.0011226
-1986    0.0011284
-1987    0.0011343
-1988    0.0011401
-1989    0.0011460
-1990    0.0011519
-1991    0.0011577
-1992    0.0011636
-1993    0.0011694
-1994    0.0011753
-1995    0.0011811
-1996    0.0011870
-1997    0.0011928
-1998    0.0011987
-1999    0.0012045
-2000    0.0012104
-2001	0.0012239
-2002	0.0012347
-2003	0.0012654
-2004	0.0012762
-2005	0.0012870
-2006	0.0012977
-2007	0.0013085
-2008	0.0013192
-2009	0.0013200
-2010	0.0013407
-2011	0.0013515
-2012	0.0013722
-2013	0.0013830
-2014	0.0013938
-2015	0.0014000

RBBGCMuso/inst/examples/hegyhatsal/Ndep_from2007.txt

@@ -1,9 +0,0 @@
-2007	0.0013085
-2008	0.0013192
-2009	0.0013200
-2010	0.0013407
-2011	0.0013515
-2012	0.0013722
-2013	0.0013830
-2014	0.0013938
-2014	0.0014000

RBBGCMuso/inst/examples/hegyhatsal/c3grass_apriori_MuSo5.epc

@@ -1,180 +0,0 @@
-ECOPHYS FILE -    C3 grass muso5 
-----------------------------------------------------------------------------------------
-FLAGS
-0	(flag)		biome type flag (1 = WOODY 0 = NON-WOODY)
-0	(flag)		woody type flag (1 = EVERGREEN 0 = DECIDUOUS)
-1	(flag)		photosyn. type flag (1 = C3 PSN 0 = C4 PSN)
-1	(flag)		phenology flag (1 = MODEL PHENOLOGY 0 = USER-SPECIFIED PHENOLOGY)
-0	(flag)		transferGDD flag (1= transfer calc. from GDD 0 = transfer calc. from EPC)
-1       (flag)  	q10 flag (1 = temperature dependent q10 value; 0= constans q10 value)
-0       (flag)  	acclimation flag (1 = acclimation 0 = no acclimation)
-1	(flag)  	CO2 conductance reduction flag (0: no effect, 1: multiplier)
-1       (flag)  	use GSI index to calculate growing season (0: no GSI, 1: GSI)
-0       (flag)  	soil temperature calculation method (0: Zheng, 1: DSSAT)
-1       (flag)  	soil hydrological calculation method (0: Richards, 1: tipping DSSAT)
-0       (int)   	discretization level of SWC calculation (0: low, 1: medium, 2: high)
-0	(flag)          photosynthesis calculation method (0: Farquhar, 1: DSSAT)
-0	(flag)          evapotranspiration calculation method (0: Penman-Montieth, 1: Priestly-Taylor)
-0	(flag)          radiation calculation method (0: SWabs, 1: Rn)
-----------------------------------------------------------------------------------------
-PLANT FUNCTIONING PARAMETERS
-0	(yday)		yearday to start new growth (when phenology flag = 0)
-364	(yday)		yearday to end litterfall (when phenology flag = 0)
-0.5	(prop.)		transfer growth period as fraction of growing season (when transferGDD_flag = 0)
-0.5	(prop.)		litterfall as fraction of growing season (when transferGDD_flag = 0)
-5	(Celsius) 	base temperature
-0	(Celsius) 	minimum  temperature for growth displayed on current day (-9999: no T-dependence of allocation)
-10	(Celsius) 	optimal1 temperature for growth displayed on current day (-9999: no T-dependence of allocation) 
-30	(Celsius) 	optimal2 temperature for growth displayed on current day (-9999: no T-dependence of allocation) 
-40	(Celsius) 	maxmimum temperature for growth displayed on current day (-9999: no T-dependence of allocation)
-0	(Celsius) 	minimum  temperature for carbon assimilation displayed on current day (-9999: no limitation)
-10	(Celsius) 	optimal1 temperature for carbon assimilation displayed on current day (-9999: no limitation) 
-30	(Celsius) 	optimal2 temperature for carbon assimilation displayed on current day (-9999: no limitation) 
-50	(Celsius) 	maxmimum temperature for carbon assimilation displayed on current day (-9999: no limitation)
-1.0	(1/yr)		annual leaf and fine root turnover fraction
-0.00	(1/yr)		annual live wood turnover fraction
-0.05	(1/yr)		annual whole-plant mortality fraction
-0.0	(1/yr)		annual fire mortality fraction
-25.0	(kgC/kgN)	C:N of leaves
-45.0	(kgC/kgN)	C:N of leaf litter, after retranslocation
-50.0	(kgC/kgN)	C:N of fine roots
-25.0   *(kgC/kgN)       C:N of fruit  
-25.0    (kgC/kgN)  	C:N of soft stem 
-0.0    *(kgC/kgN)	C:N of live wood
-0.0    *(kgC/kgN)	C:N of dead wood
-0.68	(DIM)		leaf litter labile proportion
-0.23	(DIM)		leaf litter cellulose proportion
-0.34	(DIM)		fine root labile proportion
-0.44	(DIM)		fine root cellulose proportion
-0.68   *(DIM)		fruit litter labile proportion
-0.23   *(DIM)		fruit litter cellulose proportion
-0.68    (DIM)		soft stem litter labile proportion
-0.23    (DIM)		soft stem litter cellulose proportion
-0.00   *(DIM)		dead wood cellulose proportion
-0.01	(1/LAI/d)	canopy water interception coefficient
-0.50	(DIM)		canopy light extinction coefficient
-2.0     (g/MJ)          potential radiation use efficiency
-0.781   (DIM)           radiation parameter1 (Jiang et al.2015)
--13.596 (DIM)           radiation parameter2 (Jiang et al.2015)
-2.0	(DIM)		all-sided to projected leaf area ratio
-2.0	(DIM)		ratio of shaded SLA:sunlit SLA
-0.2	(DIM)		fraction of leaf N in Rubisco
-0.03	(DIM)		fraction of leaf N in PEP Carboxylase 
-0.004	(m/s)		maximum stomatal conductance (projected area basis)
-0.00006	(m/s)		cuticular conductance (projected area basis)
-0.04	(m/s)		boundary layer conductance (projected area basis)
-1.00	(prop)		relative SWC (prop. to FC)  to calc. soil moisture limit 1 (-9999: not used)
-0.99	(prop)		relative SWC (prop. to SAT) to calc. soil moisture limit 2 (-9999: not used)
--9999	(prop)		relative PSI (prop. to FC)  to calc. soil moisture limit 1 (-9999: not used)
--9999	(prop)		relative PSI (prop. to SAT) to calc. soil moisture limit 2 (-9999: not used)
-1000	(Pa)		vapor pressure deficit: start of conductance reduction
-5000	(Pa)		vapor pressure deficit: complete conductance reduction
-1.5     (m) 		maximum height of plant
-250	*(kgC/m2) 	stem weight at which maximum height attended
-1.0     (m) 		maximum depth of rooting zone
-3.67	(DIM) 		root distribution parameter
-0.4	(kgC/m2)	rootlenght parameter 1 (estimated max root weight)
-0.5	(prop.)		rootlenght parameter 2 (slope)
-0.3     (prop.) 	growth resp per unit of C grown
-0.218   (kgC/kgN/d) 	maintenance respiration in kgC/day per kg of tissue N 
-0.1     (DIM)           theoretical maximum prop. of non-structural and structural carbohydrates
-0.3	(DIM)		prop. of non-structural carbohydrates available for maintanance respiration
-0.0005  (kgN/m2/yr) 	symbiotic+asymbiotic fixation of N
-----------------------------------------------------------------------------------------
-CROP SPECIFIC PARAMETERS
-0   	(DIM)           number of phenophase of germination (from 1 to 7; 0: NO specific)
-0   	(DIM)           number of phenophase of emergence (from 1 to 7; 0: NO specific)
-0.5     (prop.)         critical relative SWC (prop. to FC) in germination
--9999   (DIM)           number of phenophase of photoperiodic slowing effect (from 1 to 7; -9999: NO effect)
-20      (hour) 	        critical photoslow daylength
-0.005   (DIM) 	        slope of relative photoslow development rate 
--9999   (DIM)           number of phenophase of vernalization (from 1 to 7; -9999: NO)
-0       (Celsius) 	critical vernalization temperature 1
-5       (Celsius) 	critical vernalization temperature 2
-8       (Celsius) 	critical vernalization temperature 3
-15      (Celsius) 	critical vernalization temperature 4
-0.04    (DIM) 	        slope of relative vernalization development rate 
-50      (n)             required vernalization days (in vernalization development rate)
--9999	(DIM)           number of flowering phenophase (from 1 to 7; -9999: NO effect)
-35	(Celsius)       critical flowering heat stress temperature 1
-40	(Celsius)       critical flowering heat stress temperature 2
-0.2	(prop.)         mortality parameter of flowering thermal stress
-----------------------------------------------------------------------------------------
-SENESCENCE AND SOIL PARAMETERS
-0.05    (prop.) 	maximum senescence mortality coefficient of aboveground plant material
-0.05    (prop.) 	maximum senescence mortality coefficient of belowground plant material
-0.0     (prop.) 	maximum senescence mortality coefficient of non-structured plant material
-2	(prop) 		effect of extreme high temperature on senesncene mortality
-30	(Celsius) 	lower limit extreme high temperature effect on senescence mortality
-40	(Celsius) 	upper limit extreme high temperature effect on senescence mortality
--9999   (Celsius) 	maximal lifetime of plant tissue  (-9999: no genetically programmed senescence)
-0.01    (prop.) 	turnover rate of wilted standing biomass to litter
-0.05    (prop.) 	turnover rate of non-woody cut-down biomass to litter
-0.01    (prop.) 	turnover rate of woody cut-down biomass to litter
-30	(prop.)		drought tolerance parameter (critical value of DSWS)
-0.02    (prop.) 	denitrification rate per g of CO2 respiration of SOM
-0.2     (prop.) 	nitrification coefficient 1 
-0.1	(prop.) 	nitrification coefficient 2
-0.02    (prop.)		coefficient of N2O emission of nitrification
-0.5     (prop.) 	proprortion of NH4 flux of N-deposition
-0.1     (prop.) 	NH4 mobilen proportion
-1.0     (prop.) 	NO3 mobilen proportion
-10      (m) 		e-folding depth of decomposition rate's depth scalar
-0.001	(prop.) 	fraction of dissolved part of SOIL1 organic matter
-0.001	(prop.) 	fraction of dissolved part of SOIL2 organic matter
-0.001	(prop.) 	fraction of dissolved part of SOIL3 organic matter
-0.001	(prop.) 	fraction of dissolved part of SOIL4 organic matter
-0.6	(DIM)           ratio of bare soil evaporation and pot.evaporation
-----------------------------------------------------------------------------------------
-RATE SCALARS
-0.39 	(DIM) 		respiration fractions for fluxes between compartments (l1s1)
-0.55 	(DIM) 		respiration fractions for fluxes between compartments (l2s2)
-0.29 	(DIM) 		respiration fractions for fluxes between compartments (l4s3)
-0.28 	(DIM) 		respiration fractions for fluxes between compartments (s1s2)
-0.46 	(DIM) 		respiration fractions for fluxes between compartments (s2s3)
-0.55 	(DIM) 		respiration fractions for fluxes between compartments (s3s4)
-0.7     (DIM) 		rate constant scalar of  labile litter pool
-0.07    (DIM) 		rate constant scalar of cellulose litter pool
-0.014   (DIM) 		rate constant scalar of lignin litter pool
-0.07    (DIM) 		rate constant scalar of  fast microbial recycling pool
-0.014   (DIM) 		rate constant scalar of  medium microbial recycling pool 
-0.0014  (DIM) 		rate constant scalar of  slow microbial recycling pool
-0.0001  (DIM) 		rate constant scalar of  recalcitrant SOM (humus) pool
-0.001   (DIM) 		rate constant scalar of  physical fragmentation of coarse woody debris
-----------------------------------------------------------------------------------------
-GROWING SEASON PARAMETERS
-5       (kg/m2) 	crit. amount of snow limiting photosyn.
-20      (Celsius) 	limit1 (under:full constrained) of HEATSUM index
-60      (Celsius) 	limit2 (above:unconstrained) of HEATSUM index
-0       (Celsius) 	limit1 (under:full constrained) of TMIN index
-5       (Celsius) 	limit2 (above:unconstrained) of TMIN index
-4000    (Pa) 		limit1 (above:full constrained) of VPD index
-1000    (Pa) 		limit2 (under:unconstrained) of VPD index
-0       (s) 		limit1 (under:full constrained) of DAYLENGTH index
-0       (s) 		limit2 (above:unconstrained) of DAYLENGTH index
-10      (day) 		moving average (to avoid the effects of extreme events)
-0.10    (dimless) 	GSI limit1 (greater that limit -> start of vegper)
-0.01    (dimless) 	GSI limit2 (less that limit -> end of vegper)
-----------------------------------------------------------------------------------------
-CH4 PARAMETERS
-212.5	(DIM) 		param1 for CH4 calculations (empirical function of BD)
-1.81	(DIM) 		param2 for CH4 calculations (empirical function of BD)
--1.353	(DIM) 		param1 for CH4 calculations (empirical function of VWC)
-0.2	(DIM) 		param2 for CH4 calculations (empirical function of VWC)
-1.781	(DIM) 		param3 for CH4 calculations (empirical function of VWC)
-6.786	(DIM) 		param4 for CH4 calculations (empirical function of VWC)
-0.010	(DIM) 		param1 for CH4 calculations (empirical function of Tsoil)
-----------------------------------------------------------------------------------------
-PHENOLOGICAL (ALLOCATION) PARAMETERS (7 phenological phases)
-phase1  phase2  phase3  phase4 	phase5 	phase6 	phase7  (text) name of the phenophase
-500	200	500	200	400	200	100	(Celsius) length of phenophase (GDD)
-0.4	0.4	0.4	0.4	0.4	0.4  	0.4	(ratio) leaf ALLOCATION 	                                                                                                                                                                                                                                                           
-0.2	0.2	0.2	0.2	0.2	0.2  	0.2	(ratio) fine root ALLOCATION
-0.2	0.2	0.2	0.2	0.2	0.2  	0.2	(ratio) fruit ALLOCATION 	                                                                                                                                                                                                                                             
-0.2	0.2	0.2	0.2	0.2	0.2  	0.2	(ratio) soft stem ALLOCATION
-0	0	0	0	0	0	0	(ratio) live woody stem ALLOCATION 
-0	0	0	0	0	0	0	(ratio) dead woody stem ALLOCATION 
-0	0	0	0	0	0	0	(ratio) live coarse root ALLOCATION 
-0	0	0	0	0	0	0	(ratio) dead coarse root ALLOCATION 
-49	49	49	49	49	49	49	(m2/kgC) canopy average specific leaf area (projected area basis)
-0.5	0.5	0.5	0.5	0.5	0.5     0.5	(prop.)	current growth proportion

RBBGCMuso/inst/examples/hegyhatsal/mowing_hhs_2007-2015.txt

@@ -1,30 +0,0 @@
-(yday) MOWING day
-168     -9999    -9999   -9999   -9999   -9999   -9999 
-151     231      -9999   -9999   -9999   -9999   -9999   
-159     219      -9999   -9999   -9999   -9999   -9999  
-163     269      -9999   -9999   -9999   -9999   -9999     
-152     233      -9999   -9999   -9999   -9999   -9999     
-145     230      -9999   -9999   -9999   -9999   -9999     
-167     272      -9999   -9999   -9999   -9999   -9999     
-160     272      -9999   -9999   -9999   -9999   -9999     
-164     273      -9999   -9999   -9999   -9999   -9999        
-(m2/m2) value of the LAI after MOWING (fixday method)                                     
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999    
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999   
-0.52    0.52    -9999   -9999   -9999   -9999   -9999    
-(%) transported part of plant material
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999    
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999  
-90.0    90.0    -9999   -9999   -9999   -9999   -9999        

RBBGCMuso/inst/examples/hegyhatsal/n.ini

@@ -1,204 +0,0 @@
-BBGC_MuSo simulation (missing data: -9999)
-
-MET_INPUT
-hhs_2007-2015.mtc43                                                                            (filename) met file name
-4                                                                                      (int) number of header lines in met file
-365                                                             (int) number of simdays in last simyear (truncated year: < 365)
-
-RESTART
-1                                                                                (flag) 1 = read restart; 0 = dont read restart
-0                                                                              (flag) 1 = write restart; 0 = dont write restart
-0                                                                             (flag) 1 = use restart metyear; 0 = reset metyear
-hhs_apriori_MuSo5.endpoint                                                    (filename) name of the input restart file
-hhs_apriori_MuSo5.endpoint                                                   (filename) name of the output restart file
-
-TIME_DEFINE
-9                                                                                     (int) number of meteorological data years
-9                                                                                             (int) number of simulation years 
-2007                                                                                                (int) first simulation year
-0                                                                                         (flag) 1 = spinup run; 0 = normal run
-6000                                                                                       (int) maximum number of spinup years
-
-CLIM_CHANGE
-0.0                                                                                                    (degC) - offset for Tmax
-0.0                                                                                                    (degC) - offset for Tmin
-1.0                                                                                                (degC) - multiplier for PRCP
-1.0                                                                                                 (degC) - multiplier for VPD
-1.0                                                                                                 (degC) - multiplier for RAD
-
-CO2_CONTROL
-1                                                                                           (flag) 0=constant; 1=vary with file
-395.0                                                                              (ppm) constant atmospheric CO2 concentration
-CO2_from2007.txt                                                                            (filename) name of the CO2 file
-
-NDEP_CONTROL
-1                                                                                           (flag) 0=constant; 1=vary with file
-0.001400                                                                        (kgN/m2/yr) wet+dry atmospheric deposition of N
-Ndep_from2007.txt                                                                   (filename)  name of the N-dep file
-
-SITE
-30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0                               (%) sand percentage by volume in rock-free soil
-50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0                               (%) silt percentage by volume in rock-free soil
-7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0                                                                      (dimless) soil pH 
-248.0                                                                                                       (m)  site elevation
-46.95                                                                                    (degrees) site latitude (- for S.Hem.)
-0.20                                                                                                (DIM) site shortwave albedo
-9.00                                                                                      (Celsius) mean annual air temperature
-10.15                                                                               (Celsius) mean annual air temperature range
-20.00                                                                                         (mm) maximum height of pond water
--9999                                                                                    (dimless) measured runoff curve number
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (g/cm3) bulk density
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                           (m3/m3) SWC at saturation
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                       (m3/m3) SWC at field capacity
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (m3/m3) SWC at wilting point
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                            (m3/m3) SWC at hygroscopic water content
-
-EPC_FILE
-c3grass_apriori_MuSo5.epc                                                                          (filename) EPC file name
-
-W_STATE
-0.0                                                                                            (kg/m2) water stored in snowpack
-1.0                                                                  (DIM) initial soil water as a proportion of field capacity
-
-CN_STATE
-0.001                                                                                  (kgC/m2)  first-year maximum leaf carbon
-0.001                                                                             (kgC/m2)  first-year maximum fine root carbon
-0.001                                                                                 (kgC/m2)  first-year maximum fruit carbon
-0.001                                                                              (kgC/m2)  first-year maximum softstem carbon
-0.001                                                                       (kgC/m2)  first-year maximum live woody stem carbon
-0.001                                                                      (kgC/m2)  first-year maximum live coarse root carbon
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  coarse woody debris carbon
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  litter carbon, labile pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                    (kgC/m2)  litter carbon, unshielded cellulose pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                      (kgC/m2)  litter carbon, shielded cellulose pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  litter carbon, lignin pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                  (kgC/m2)  soil carbon, fast microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                (kgC/m2)  soil carbon, medium microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                  (kgC/m2)  soil carbon, slow microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                     (kgC/m2)  soil carbon, recalcitrant SOM (slowest)
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                (kgN/m2)  litter nitrogen, labile pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                         (kgN/m2)  soil mineralized nitrogen, NH4 pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                         (kgN/m2)  soil mineralized nitrogen, NO3 pool
-
-OUTPUT_CONTROL
-hhs_apriori_MuSo5                                                                              (filename) output prefix
-1                                                                  (flag)  writing daily output (0 = no; 1 = binary; 2 = ascii)
-0                                               (flag)  writing monthly average of daily output (0 = no; 1 = binary; 2 = ascii)
-0                                                (flag)  writing annual average of daily output (0 = no; 1 = binary; 2 = ascii)
-2                                                                 (flag)  writing annual output (0 = no; 1 = binary; 2 = ascii)
-1                                                                                      (flag)  for on-screen progress indicator
-
-DAILY_OUTPUT
-14                                                                                             number of daily output variables
-49                                                                                                                    tsoil_2cm
-52                                                                                                                   tsoil_20cm
-2599                                                                                                                    vwc_2cm
-2602                                                                                                                   vwc_20cm
-2751                                                                                                               m_soilstress
-671                                                                                                                vegc_to_SNSC
-171                                                                                                                 evapotransp
-2520                                                                                                                   proj_LAI
-407                                                                                                                 STDBc_above
-3004                                                                                                                  daily_n2o
-3008                                                                                                                  daily_nbp
-3009                                                                                                                  daily_gpp
-3014                                                                                                                   daily_tr
-2502                                                                                                                  n_actphen
-
-ANNUAL_OUTPUT
-12                                                                                            number of annual output variables
-0                                                                                                                         onday
-1                                                                                                                        offday
-3000                                                                                                                    annprcp
-3001                                                                                                                    anntavg
-3002                                                                                                                  annrunoff
-2765                                                                                                                ytd_maxplai
-3024                                                                                                                    cum_nee
-3033                                                                                                            cum_Closs_THN_w
-3035                                                                                                              cum_Closs_MOW
-3037                                                                                                             cum_yieldC_HRV
-3059                                                                                                                  NH4_top10
-3060                                                                                                                  NO3_top10
-
--------------------
-MANAGEMENT_SECTION
--------------------
-PLANTING
-0                                                                   (flag) do PLANTING? 0=no; 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday)  PLANTING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (m) germination depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                          (n/m2) number of seedlings
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                       (g/1000n) weight of 1000-seed
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                              (%)  C content of seed
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                  (%) emergence rate
-
-THINNING
-0                                                                   (flag) do THINNING? 0=no; 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday)  THINNING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                        (prop) thinning rate of woody plant material
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                    (prop) thinning rate of non-woody plant material
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (%) transported part of stem
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (%) transported part of leaf
-
-MOWING
-mowing_hhs_2007-2015.txt                                      (flag) do MOWING? 0=no; 1=yes; filepath=reading from file
-0                                                                  (flag) mowing method? 0 - fixday method, 1 - fixvalue method
--9999                                                            (m2/m2) fixed value of the LAI before MOWING (fixvalue method)
--9999                                                             (m2/m2) fixed value of the LAI after MOWING (fixvalue method)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                   (yday) MOWING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999               (m2/m2) value of the LAI after MOWING (fixday method)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (%) transported part of plant material
-
-GRAZING
-0                                                                    (flag) do GRAZING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                         (yday) first day of GRAZING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                          (yday) last day of GRAZING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (kg/LSU) weight equivalent of one unit
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999          (LSU/ha) animal stocking rate: Livestock Units per hectare
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kg dry matter/LSU) daily ingested dry matter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999             (prop) trampling effect (standing dead biome to litter)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                             (%) ratio of DM intake formed excrement
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                          (%) ratio of excrement returning to litter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                    (%) carbon content of dry matter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                             (%) N content of manure
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                             (%) C content of manure
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999        (kgN2O-N:kgN) manure emission factor for direct N2O (T10.21)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999     (kgN/1000 kg animal mass/day) default N excretion rate (T10.19)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999             (kgCH4/LSU/yr*) manure emission factor for CH4 (T10.14)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999       (kgCH4/LSU/yr*) fermentation emission factor for CH4 (T10.11)
-
-HARVESTING
-0                                                                 (flag) do HARVESTING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (yday) HARVESTING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kgC/m2) soft stem C content after HARVESTING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (%) transported part of plant material
-
-PLOUGHING
-0                                                                  (flag) do PLOUGHING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday) PLOUGHING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                 (m) PLOUGHING depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999           (prop) dissolving coefficient of ploughed biome to litter
-
-FERTILIZING
-0                                                                (flag) do FERTILIZING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                              (yday) FERTILIZING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (m) fertilizing depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                        (kgN fertilizer/ha/day) amount of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                (%) dry matter content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                   (kgN/kg fertilizer) nitrate content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                  (kgN/kg fertilizer) ammonium content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999          (kgN/kg fertilizer) organic nitrogen content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999            (kgC/kg fertilizer) organic carbon content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                   (%) labile fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                     (%) unshielded cellulose fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (%) shielded cellulose fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                  (%) lignin  fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kgN2O-N:kgN) emission factor for N-additions
-
-IRRIGATION
-0                                                                 (flag) do IRRIGATION? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (yday) IRRIGATION day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                      (kgH2O/m2/day) amount of water
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                (%) utilization coefficient of water
-
-END_INIT

RBBGCMuso/inst/examples/hegyhatsal/s.ini

@@ -1,204 +0,0 @@
-BBGC_MuSo simulation (missing data: -9999)
-
-MET_INPUT
-hhs_1961-2014.mtc43                                                                            (filename) met file name
-4                                                                                      (int) number of header lines in met file
-365                                                             (int) number of simdays in last simyear (truncated year: < 365)
-
-RESTART
-0                                                                                (flag) 1 = read restart; 0 = dont read restart
-1                                                                              (flag) 1 = write restart; 0 = dont write restart
-0                                                                             (flag) 1 = use restart metyear; 0 = reset metyear
-hhs_apriori_MuSo5.endpoint                                                    (filename) name of the input restart file
-hhs_apriori_MuSo5.endpoint                                                   (filename) name of the output restart file
-
-TIME_DEFINE
-54                                                                                    (int) number of meteorological data years
-54                                                                                            (int) number of simulation years 
-1961                                                                                                (int) first simulation year
-1                                                                                         (flag) 1 = spinup run; 0 = normal run
-6000                                                                                       (int) maximum number of spinup years
-
-CLIM_CHANGE
-0.0                                                                                                    (degC) - offset for Tmax
-0.0                                                                                                    (degC) - offset for Tmin
-1.0                                                                                                (degC) - multiplier for PRCP
-1.0                                                                                                 (degC) - multiplier for VPD
-1.0                                                                                                 (degC) - multiplier for RAD
-
-CO2_CONTROL
-1                                                                                           (flag) 0=constant; 1=vary with file
-290.0                                                                              (ppm) constant atmospheric CO2 concentration
-CO2_from1961.txt                                                                            (filename) name of the CO2 file
-
-NDEP_CONTROL
-1                                                                                           (flag) 0=constant; 1=vary with file
-0.000200                                                                        (kgN/m2/yr) wet+dry atmospheric deposition of N
-Ndep_from1961.txt                                                                   (filename)  name of the N-dep file
-
-SITE
-30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0 30.0                               (%) sand percentage by volume in rock-free soil
-50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0 50.0                               (%) silt percentage by volume in rock-free soil
-7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0 7.0                                                                      (dimless) soil pH 
-248.0                                                                                                       (m)  site elevation
-46.95                                                                                    (degrees) site latitude (- for S.Hem.)
-0.20                                                                                                (DIM) site shortwave albedo
-9.00                                                                                      (Celsius) mean annual air temperature
-10.15                                                                               (Celsius) mean annual air temperature range
-20.00                                                                                         (mm) maximum height of pond water
--9999                                                                                    (dimless) measured runoff curve number
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (g/cm3) bulk density
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                           (m3/m3) SWC at saturation
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                       (m3/m3) SWC at field capacity
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (m3/m3) SWC at wilting point
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                            (m3/m3) SWC at hygroscopic water content
-
-EPC_FILE
-c3grass_apriori_MuSo5.epc                                                                          (filename) EPC file name
-
-W_STATE
-0.0                                                                                            (kg/m2) water stored in snowpack
-1.0                                                                  (DIM) initial soil water as a proportion of field capacity
-
-CN_STATE
-0.001                                                                                  (kgC/m2)  first-year maximum leaf carbon
-0.001                                                                             (kgC/m2)  first-year maximum fine root carbon
-0.001                                                                                 (kgC/m2)  first-year maximum fruit carbon
-0.001                                                                              (kgC/m2)  first-year maximum softstem carbon
-0.001                                                                       (kgC/m2)  first-year maximum live woody stem carbon
-0.001                                                                      (kgC/m2)  first-year maximum live coarse root carbon
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  coarse woody debris carbon
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  litter carbon, labile pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                    (kgC/m2)  litter carbon, unshielded cellulose pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                      (kgC/m2)  litter carbon, shielded cellulose pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                  (kgC/m2)  litter carbon, lignin pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                  (kgC/m2)  soil carbon, fast microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                (kgC/m2)  soil carbon, medium microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                  (kgC/m2)  soil carbon, slow microbial recycling pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                     (kgC/m2)  soil carbon, recalcitrant SOM (slowest)
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                                (kgN/m2)  litter nitrogen, labile pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                         (kgN/m2)  soil mineralized nitrogen, NH4 pool
-0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0   0.0                         (kgN/m2)  soil mineralized nitrogen, NO3 pool
-
-OUTPUT_CONTROL
-hhs_apriori_MuSo5_S                                                                            (filename) output prefix
-0                                                                  (flag)  writing daily output (0 = no; 1 = binary; 2 = ascii)
-0                                               (flag)  writing monthly average of daily output (0 = no; 1 = binary; 2 = ascii)
-0                                                (flag)  writing annual average of daily output (0 = no; 1 = binary; 2 = ascii)
-0                                                                 (flag)  writing annual output (0 = no; 1 = binary; 2 = ascii)
-1                                                                                      (flag)  for on-screen progress indicator
-
-DAILY_OUTPUT
-14                                                                                             number of daily output variables
-49                                                                                                                    tsoil_2cm
-52                                                                                                                   tsoil_20cm
-2599                                                                                                                    vwc_2cm
-2602                                                                                                                   vwc_20cm
-2751                                                                                                               m_soilstress
-671                                                                                                                vegc_to_SNSC
-171                                                                                                                 evapotransp
-2520                                                                                                                   proj_LAI
-407                                                                                                                 STDBc_above
-3004                                                                                                                  daily_n2o
-3008                                                                                                                  daily_nbp
-3009                                                                                                                  daily_gpp
-3014                                                                                                                   daily_tr
-2502                                                                                                                  n_actphen
-
-ANNUAL_OUTPUT
-12                                                                                            number of annual output variables
-0                                                                                                                         onday
-1                                                                                                                        offday
-3000                                                                                                                    annprcp
-3001                                                                                                                    anntavg
-3002                                                                                                                  annrunoff
-2765                                                                                                                ytd_maxplai
-3024                                                                                                                    cum_nee
-3033                                                                                                            cum_Closs_THN_w
-3035                                                                                                              cum_Closs_MOW
-3037                                                                                                             cum_yieldC_HRV
-3059                                                                                                                  NH4_top10
-3060                                                                                                                  NO3_top10
-
--------------------
-MANAGEMENT_SECTION
--------------------
-PLANTING
-0                                                                   (flag) do PLANTING? 0=no; 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday)  PLANTING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (m) germination depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                          (n/m2) number of seedlings
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                       (g/1000n) weight of 1000-seed
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                              (%)  C content of seed
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                  (%) emergence rate
-
-THINNING
-0                                                                   (flag) do THINNING? 0=no; 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday)  THINNING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                        (prop) thinning rate of woody plant material
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                    (prop) thinning rate of non-woody plant material
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (%) transported part of stem
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                        (%) transported part of leaf
-
-MOWING
-0                                                                     (flag) do MOWING? 0=no; 1=yes; filepath=reading from file
-0                                                                  (flag) mowing method? 0 - fixday method, 1 - fixvalue method
--9999                                                            (m2/m2) fixed value of the LAI before MOWING (fixvalue method)
--9999                                                             (m2/m2) fixed value of the LAI after MOWING (fixvalue method)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                   (yday) MOWING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999               (m2/m2) value of the LAI after MOWING (fixday method)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (%) transported part of plant material
-
-GRAZING
-0                                                                    (flag) do GRAZING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                         (yday) first day of GRAZING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                          (yday) last day of GRAZING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (kg/LSU) weight equivalent of one unit
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999          (LSU/ha) animal stocking rate: Livestock Units per hectare
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kg dry matter/LSU) daily ingested dry matter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999             (prop) trampling effect (standing dead biome to litter)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                             (%) ratio of DM intake formed excrement
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                          (%) ratio of excrement returning to litter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                    (%) carbon content of dry matter
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                             (%) N content of manure
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                             (%) C content of manure
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999        (kgN2O-N:kgN) manure emission factor for direct N2O (T10.21)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999     (kgN/1000 kg animal mass/day) default N excretion rate (T10.19)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999             (kgCH4/LSU/yr*) manure emission factor for CH4 (T10.14)
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999       (kgCH4/LSU/yr*) fermentation emission factor for CH4 (T10.11)
-
-HARVESTING
-0                                                                 (flag) do HARVESTING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (yday) HARVESTING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kgC/m2) soft stem C content after HARVESTING
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                              (%) transported part of plant material
-
-PLOUGHING
-0                                                                  (flag) do PLOUGHING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                (yday) PLOUGHING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                                 (m) PLOUGHING depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999           (prop) dissolving coefficient of ploughed biome to litter
-
-FERTILIZING
-0                                                                (flag) do FERTILIZING? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                              (yday) FERTILIZING day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (m) fertilizing depth
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                        (kgN fertilizer/ha/day) amount of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                (%) dry matter content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                   (kgN/kg fertilizer) nitrate content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                  (kgN/kg fertilizer) ammonium content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999          (kgN/kg fertilizer) organic nitrogen content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999            (kgC/kg fertilizer) organic carbon content of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                   (%) labile fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                     (%) unshielded cellulose fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (%) shielded cellulose fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                  (%) lignin  fraction of fertilizer
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                       (kgN2O-N:kgN) emission factor for N-additions
-
-IRRIGATION
-0                                                                 (flag) do IRRIGATION? 0=no, 1=yes; filepath=reading from file
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                               (yday) IRRIGATION day
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                      (kgH2O/m2/day) amount of water
--9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999 -9999                                (%) utilization coefficient of water
-
-END_INIT

RBBGCMuso/man/changeMusoC.Rd

@@ -0,0 +1,16 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/RcppExports.R
+\name{changeMusoC}
+\alias{changeMusoC}
+\title{changeMusoC}
+\usage{
+changeMusoC(inFile, outFile, inMat)
+}
+\arguments{
+\item{inFile}{is the big matrix}
+
+\item{outFile}{is the small matrix}
+}
+\description{
+This function is fastly randomize values based on min and max values
+}

RBBGCMuso/man/changemulline.Rd

@@ -4,7 +4,8 @@
 \alias{changemulline}
 \title{This is the function which is capable change multiple specific lines to other using their row numbers.}
 \usage{
-changemulline(filename, calibrationPar, contents)
+changemulline(filePaths, calibrationPar, contents, fileOut, fileToChange,
+  modifyOriginal = FALSE)
 }
 \description{
 he function uses the previous changspecline function to operate.

RBBGCMuso/man/compareMuso.Rd

@@ -0,0 +1,26 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/plotMuso.R
+\name{compareMuso}
+\alias{compareMuso}
+\title{compareMuso}
+\usage{
+compareMuso(settings, parameters, variable, calibrationPar = NULL,
+  fileToChange = NULL, skipSpinup = TRUE, timeFrame = "day")
+}
+\arguments{
+\item{settings}{You have to run the setupMuso function before rungetMuso. It is its output which contains all of the necessary system variables. It sets the whole environment}
+
+\item{parameters}{In the settings variable you have set the row indexes of the variables, you wish to change. In this parameter you can give an exact value for them in a vector like: c(1,2,3,4)}
+
+\item{variable}{The name of the output variable to plot}
+
+\item{calibrationPar}{in the help of setupMuso function.}
+
+\item{fileToChange}{You can change any line of the epc or the ini file, you just have to specify with this variable which file you van a change. Two options possible: "epc", "ini", "both"}
+}
+\description{
+This function runs the modell, change one of it's input, and plot both in one plot.
+}
+\author{
+Roland Hollos
+}

RBBGCMuso/man/copyMusoExamleTo.Rd

@@ -1,16 +0,0 @@
-% Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/musoExamples.R
-\name{copyMusoExamleTo}
-\alias{copyMusoExamleTo}
-\title{copyMusoExampleTo}
-\usage{
-copyMusoExamleTo(example = NULL, destination = NULL)
-}
-\arguments{
-\item{example}{The name of the example file, if it is NULL tcl/tk menu will pop up to select.}
-
-\item{destination}{The destination where the example files will be copied.}
-}
-\description{
-With this function you can copy RBBGCMuso example library whereever you want
-}

RBBGCMuso/man/musoMont.Rd

@@ -0,0 +1,41 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/musoMont.R
+\name{musoMont}
+\alias{musoMont}
+\title{musoMont}
+\usage{
+musoMont(settings = NULL, parameters = NULL, inputDir = "./",
+  outLoc = "./calib", iterations = 10, preTag = "mont-",
+  outputType = "moreCsv", fun = mean, varIndex = 1, outVars = NULL,
+  silent = TRUE, skipSpinup = TRUE, debugging = FALSE,
+  keepEpc = FALSE, constrains = NULL, ...)
+}
+\arguments{
+\item{settings}{A list of montecarlos environmental variables. It is generated by the setupMuso() function. In default the settings parameter is generated automatically.}
+
+\item{parameters}{This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the epc-fie, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.}
+
+\item{inputDir}{The location of the input directory, this directory must content a viable pack of all inputfiles and the executable file.}
+
+\item{iterations}{number of the monteCarlo run.}
+
+\item{preTag}{It will be the name of the output files. For example preTag-1.csv, pretag-2csv...}
+
+\item{outputType}{This parameter can be "oneCsv", "moreCsv", and "netCDF". If "oneCsv" is choosen the function create 1 big csv file for all of the runs, if "moreCsv" is choosen, every modell output goes to separate files, if netCDF is selected the outputs will be put in a netCDF file. The default value of the outputTypes is "moreCsv". netCDF is not implemented yet.}
+
+\item{fun}{If you select a variable from the possible outputs (with specify the varIndex parameter), you have to provide a function which maps to a subset of real numbers. The most frequent possibilities are: mean, min, max, var, but you can define any function for your need.}
+
+\item{varIndex}{This parameter specify which parameter of the output will be used. You can extract this information from the ini-files. At the output parameter specifications, the parameters order will determine this number. For example, if you have set these output parameters: 412, 874, 926, 888, and you want to use 926, you should address varIndex with 3.}
+
+\item{debugging}{If you set this parameter, you can save every logfile, and RBBGCMuso will select those which contains errors.}
+
+\item{keepEpc}{if you set keepEpc also true, it will save every selected epc file, and put the wrong ones in the WRONGEPC directory.}
+
+\item{calibrationPar}{You may want to change some parameters in your epc file, before you run the modell. You have to select the appropirate modell parameters. You can refence to these with the number of the line in the epc file where the variables are. It indexes from one. You should use a vector for this, like: c(1,5,8)}
+}
+\description{
+This function does monteCarlo on BiomeBGC-MuSo. It samples specified modell variables in given rangge from conditional multivariate uniform distribution, and runs the modell for each run.
+}
+\author{
+Roland Hollos
+}

RBBGCMuso/man/musoRand.Rd

@@ -0,0 +1,21 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/musoRand.R
+\name{musoRand}
+\alias{musoRand}
+\title{musoRand}
+\usage{
+musoRand(parameters, constrains = NULL, iterations = 3000)
+}
+\arguments{
+\item{parameters}{This is a dataframe (heterogen data-matrix), which first column is the name of the parameters, the second is a numeric vector of the rownumbers of the given variable in the input-file, the last two column consist the endpont of the parameter-ranges, where the parameters will be randomized.}
+
+\item{constrains}{This is a matrics wich specify the constrain rules for the sampling. Further informations coming son.}
+
+\item{iteration}{The number of sample-s. It is adviced to use at least 3000 iteration, because it is generally fast and it can be subsampled later at any time.}
+}
+\description{
+This funtion samples uniformly from the choosen parameters of the BiomeBGC-Muso model, which parameters are constrained by the model logic.
+}
+\author{
+Roland Hollos
+}

RBBGCMuso/man/musoSensi.Rd

@@ -5,11 +5,11 @@
 \title{musoSensi}
 \usage{
 musoSensi(monteCarloFile = NULL, parameters = NULL, settings = NULL,
-  inputDir = "./", outLoc = "./calib", iterations = 30,
+  inputDir = "./", outLoc = "./calib", outVars = NULL,
-  preTag = "mont-", outputType = "moreCsv", fun = mean,
+  iterations = 30, preTag = "mont-", outputType = "moreCsv",
-  varIndex = 1, outputFile = "sensitivity.csv",
+  fun = mean, varIndex = 1, outputFile = "sensitivity.csv",
   plotName = "sensitivity.png", plotTitle = "Sensitivity",
-  skipSpinup = FALSE, dpi = 300)
+  skipSpinup = TRUE, dpi = 300)
 }
 \arguments{
 \item{monteCarloFile}{If you run musoMonte function previously, you did not have to rerun the monteCarlo, just provide the preservedEpc.csv file with its path. If you do not set this parameter, musoSensi will fun the musoMonte function to get all of the information.}

RBBGCMuso/man/plotMusoWithData.Rd

@@ -22,6 +22,10 @@ leapYear=FALSE, export=FALSE)
 \item{sep}{This is the separator used in the measurement file}
 
 \item{savePlot}{It it is specified, the plot will be saved in a format specified with the immanent extension}
+
+\item{calibrationPar}{documentation in setupMuso()}
+
+\item{parameters}{documentation in calibMuso()}
 }
 \description{
 This function runs the BBGC-MuSo model and reads in its outputfile in a very structured way, and after that plot the results automaticly along with a given measurement

RBBGCMuso/man/putOutVars.Rd

@@ -0,0 +1,20 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/putOutVars.R
+\name{putOutVars}
+\alias{putOutVars}
+\title{putOutVars}
+\usage{
+putOutVars(iniFile, outputVars, modifyOriginal = FALSE)
+}
+\arguments{
+\item{outputVars}{List of the output codes}
+
+\item{IniFile}{The name of the normal ini file.}
+}
+\description{
+This function is for adding variables in the inifiles.
+}
+\author{
+Roland Hollos
+}
+\keyword{internal}

RBBGCMuso/man/setupMuso.Rd

@@ -22,6 +22,8 @@ nitInput=NULL, iniInput=NULL, epcInput=NULL)
 
 \item{outputLoc}{Where should the modell puts its outputs. You should give a location for it via this variable, for example: outputLoc="/place/of/the/outputs/"}
 
+\item{modelOutputs}{This parameter is the list of the codes of outpu}
+
 \item{inputLoc}{Usually it is the root directory, where you put the iniFiles for the modell}
 
 \item{metInput}{Via metInput parameter, you can tell the modell where are the meteorological files. As default it reads this from the iniFiles.}
@@ -52,7 +54,7 @@ nitInput=NULL, iniInput=NULL, epcInput=NULL)
 }
 \value{
 The output is a the modell setting list wich contains the following elements:
-executable, calibrationPar, outputLoc, outputName, inputLoc, iniInput, metInput, epcInput,thinInput,CO2Input, mowInput, grazInput, harvInput, plougInput, fertInput, irrInput, nitInput, inputFiles, numData, startyear, numYears, outputVars
+executable, calibrationPar, outputLoc, outputName, inputLoc, iniInput, metInput, epcInput,thinInput,CO2Input, mowInput, grazInput, harvInput, plougInput, fertInput,rrInput, nitInput, inputFiles, numData, startyear, numYears, outputVars
 }
 \description{
 This funcion is fundamental for the BiomBGC-MuSo modell related functions like spinupMuso, normalMuso, rungetMuso, because it sets the modells environment.

RBBGCMuso/src/RcppExports.cpp

@@ -5,6 +5,29 @@
 
 using namespace Rcpp;
 
+// getWritePositions
+IntegerVector getWritePositions(double a);
+RcppExport SEXP _RBBGCMuso_getWritePositions(SEXP aSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< double >::type a(aSEXP);
+    rcpp_result_gen = Rcpp::wrap(getWritePositions(a));
+    return rcpp_result_gen;
+END_RCPP
+}
+// changeMusoC
+void changeMusoC(std::string inFile, std::string outFile, NumericMatrix inMat);
+RcppExport SEXP _RBBGCMuso_changeMusoC(SEXP inFileSEXP, SEXP outFileSEXP, SEXP inMatSEXP) {
+BEGIN_RCPP
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< std::string >::type inFile(inFileSEXP);
+    Rcpp::traits::input_parameter< std::string >::type outFile(outFileSEXP);
+    Rcpp::traits::input_parameter< NumericMatrix >::type inMat(inMatSEXP);
+    changeMusoC(inFile, outFile, inMat);
+    return R_NilValue;
+END_RCPP
+}
 // randTypeOne
 NumericMatrix randTypeOne(NumericMatrix m);
 RcppExport SEXP _RBBGCMuso_randTypeOne(SEXP mSEXP) {
@@ -16,6 +39,17 @@ BEGIN_RCPP
     return rcpp_result_gen;
 END_RCPP
 }
+// randTypeTwo
+NumericMatrix randTypeTwo(NumericMatrix m);
+RcppExport SEXP _RBBGCMuso_randTypeTwo(SEXP mSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< NumericMatrix >::type m(mSEXP);
+    rcpp_result_gen = Rcpp::wrap(randTypeTwo(m));
+    return rcpp_result_gen;
+END_RCPP
+}
 // musoRandomizer
 NumericMatrix musoRandomizer(NumericMatrix A, NumericMatrix B);
 RcppExport SEXP _RBBGCMuso_musoRandomizer(SEXP ASEXP, SEXP BSEXP) {
@@ -30,7 +64,10 @@ END_RCPP
 }
 
 static const R_CallMethodDef CallEntries[] = {
+    {"_RBBGCMuso_getWritePositions", (DL_FUNC) &_RBBGCMuso_getWritePositions, 1},
+    {"_RBBGCMuso_changeMusoC", (DL_FUNC) &_RBBGCMuso_changeMusoC, 3},
     {"_RBBGCMuso_randTypeOne", (DL_FUNC) &_RBBGCMuso_randTypeOne, 1},
+    {"_RBBGCMuso_randTypeTwo", (DL_FUNC) &_RBBGCMuso_randTypeTwo, 1},
     {"_RBBGCMuso_musoRandomizer", (DL_FUNC) &_RBBGCMuso_musoRandomizer, 2},
     {NULL, NULL, 0}
 };

RBBGCMuso/src/musoIO.cpp

@@ -0,0 +1,138 @@
+#include <Rcpp.h>
+#include <numeric>
+#include <iostream>
+#include <fstream>
+#include <string>
+#include <new>
+#include <algorithm>
+#include <numeric>
+#include <ctime>
+#include <math.h>
+
+
+using namespace Rcpp;
+using namespace std;
+// [[Rcpp::plugins(cpp11)]]
+
+// [[Rcpp::export]]
+IntegerVector getWritePositions(double a){
+  //getWritePositions returns abstract rownumbers to rownumbers and other indexek
+  // getWritePositions(173.62) = c(173,7,3) // it supports up to 10 subvalues
+  IntegerVector outVec(3); // outVec is vector of rowIndex, colNumber, choosen Index
+  a = a * 100;
+  a = round(a); //without this line 155.92 ~= 155.9199999999, (int) 155.9199999999*100 = 15591
+  outVec[0] = (int)a / 100;
+  outVec[1] = ((int)a /10) % 10 + 1;
+  outVec[2] = (int)a % 10;
+  
+  
+  return outVec;
+}
+
+IntegerMatrix getPositions(NumericVector v){
+
+  int numVari = v.size();
+  IntegerMatrix indexek(numVari,3);
+  IntegerVector positions(3);
+  for(int i = 0; i < numVari; ++i){
+    positions = getWritePositions(v[i]);
+    indexek(i,_) = positions;
+  }
+  return indexek;
+}
+
+
+void goNextLine(std::ifstream& fin){
+  char c='a';
+  while((c!='\n')&&(fin.get(c))){}
+}
+
+#define NEXT goNextLine(fin)
+int fileChanger(std::string inFile, std::string outFile, IntegerVector linum, NumericVector num, IntegerVector colnum, IntegerVector colindex){
+  std::ifstream fin(inFile);
+  if (!fin.is_open()) {
+    stop("Cannot open " + inFile + " for read");
+  }
+  std::ofstream fot(outFile);
+
+  if (!fot.is_open()) {
+    stop("Cannot open " + outFile + " for read");    
+  }
+  string tempString;
+
+  int counter = 1;
+  int counterV = 0;
+  while (!fin.eof()) {
+
+
+    if(counter == linum[counterV]){
+      if(colnum[counterV]==1){
+	fot << num[counterV] << "\n";
+	NEXT;
+      } else {
+	double * elements;
+	elements = new double [colnum[counterV]];
+	if(linum[counterV]!=linum[counterV+1]){
+	  for(int i=0;i<colnum[counterV];++i){
+
+	    fin >> elements[i];
+	    
+	    if(i==colindex[counterV]){
+	      elements[i]=num[counterV];
+	    }
+	    // std::cout << colnum[counterV] << " " << colindex[counterV] << " " << elements[i] << " "<< i <<"\n";
+	    // std::cout << colindex[counterV] << getWritePositions(155.92) <<"\n";
+	    // std::cout << "======================== \n";
+	    fot << elements[i] << '\t';
+	  }
+	} else {
+	  int k=0;
+	  for(int i=0;i<colnum[counterV];++i){
+
+	    fin >> elements[i];
+	    if(i==colindex[counterV + k]){
+	      elements[i]=num[counterV + k];
+	      if(linum[counterV +k]==linum[counterV + k + 1]){
+		++k; 
+	      }
+	    }
+	    fot << elements[i] << '\t';
+
+	  }
+	  counterV = counterV + k;
+	}
+
+	fot << "\n";
+	delete [] elements;
+	NEXT;
+      }
+      ++counterV;
+    } else {
+      getline(fin,tempString);
+      fot << tempString << "\n";
+    }
+    ++counter;
+    if(counter > 1000){
+      stop("You modified a line which has not as many columns as you specified.");
+    }
+  }
+  fin.close();
+  fot.close();
+  return 0;
+}
+
+//' changeMusoC
+//' 
+//' This function is fastly randomize values based on min and max values
+//' @importFrom Rcpp evalCpp
+//' @useDynLib RBBGCMuso
+//' @param inFile is the big matrix
+//' @param outFile is the small matrix
+//' @export
+// [[Rcpp::export]]
+void changeMusoC(std::string inFile, std::string outFile, NumericMatrix inMat){
+  int  numChanges = inMat.nrow();
+  IntegerMatrix indexes(numChanges,3);
+  indexes = getPositions(inMat(_,0));
+  fileChanger(inFile,outFile,indexes(_,0),inMat(_,1),indexes(_,1),indexes(_,2));
+}

RBBGCMuso/src/musoIO.cpp~

@@ -0,0 +1,128 @@
+#include <Rcpp.h>
+#include <numeric>
+#include <iostream>
+#include <fstream>
+#include <string>
+#include <new>
+#include <algorithm>
+#include <numeric>
+#include <ctime>
+
+
+using namespace Rcpp;
+using namespace std;
+// [[Rcpp::plugins(cpp11)]]
+
+// [[Rcpp::export]]
+IntegerVector getWritePositions(double a){
+  //getWritePositions returns abstract rownumbers to rownumbers and other indexek
+  // getWritePositions(173.72) = c(173,7,3) // it supports up to 10 subvalues
+  IntegerVector outVec(3); // outVec is vector of rowIndex, colNumber, choosen Index
+  a = a * 100;
+  outVec[0] = (int)a / 100;
+  outVec[1] = ((int)a /10) % 10 + 1;
+  outVec[2] = (int)a % 10;
+  
+  return outVec;
+}
+
+IntegerMatrix getPositions(NumericVector v){
+
+  int numVari = v.size();
+  IntegerMatrix indexek(numVari,3);
+  IntegerVector positions(3);
+  for(int i = 0; i < numVari; ++i){
+    positions = getWritePositions(v[i]);
+    indexek(i,_) = positions;
+  }
+  return indexek;
+}
+
+
+void goNextLine(std::ifstream& fin){
+  char c='a';
+  while((c!='\n')&&(fin.get(c))){}
+}
+
+#define NEXT goNextLine(fin)
+int fileChanger(std::string inFile, std::string outFile, IntegerVector linum, NumericVector num, IntegerVector colnum, IntegerVector colindex){
+  std::ifstream fin(inFile);
+  if (!fin.is_open()) {
+    stop("Cannot open " + inFile + " for read");
+  }
+  std::ofstream fot(outFile);
+
+  if (!fot.is_open()) {
+    stop("Cannot open " + outFile + " for read");    
+  }
+  string tempString;
+
+  int counter = 1;
+  int counterV = 0;
+  while (!fin.eof()) {
+
+
+    if(counter == linum[counterV]){
+      if(colnum[counterV]==0){
+	fot << num[counterV] << "\n";
+	NEXT;
+      } else {
+	double * elements;
+	elements = new double [colnum[counterV]];
+	if(linum[counterV]!=linum[counterV+1]){
+	  for(int i=0;i<colnum[counterV];++i){
+
+	    fin >> elements[i];
+	    if(i==colindex[counterV]){
+	      elements[i]=num[counterV];
+	    }
+	    fot << elements[i] << '\t';
+	  }
+	} else {
+	  int k=0;
+	  for(int i=0;i<colnum[counterV];++i){
+
+	    fin >> elements[i];
+	    if(i==colindex[counterV + k]){
+	      elements[i]=num[counterV + k];
+	      if(linum[counterV +k]==linum[counterV + k + 1]){
+		++k;
+	      }
+	    }
+	    fot << elements[i] << '\t';
+
+	  }
+	  counterV = counterV + k;
+	}
+
+	fot << "\n";
+	delete [] elements;
+	NEXT;
+      }
+      ++counterV;
+    } else {
+      getline(fin,tempString);
+      fot << tempString << "\n";
+    }
+    ++counter;
+  }
+  fin.close();
+  fot.close();
+  return 0;
+}
+
+//' changeMusoC
+//' 
+//' This function is fastly randomize values based on min and max values
+//' @importFrom Rcpp evalCpp
+//' @useDynLib RBBGCMuso
+//' @param inFile is the big matrix
+//' @param outFile is the small matrix
+//' @export
+// [[Rcpp::export]]
+void changeMusoC(std::string inFile, std::string outFile, NumericMatrix inMat){
+  int  numChanges = inMat.nrow();
+  IntegerMatrix indexes(numChanges,3);
+  indexes = getPositions(inMat(_,0));
+  fileChanger(inFile,outFile,indexes(_,0),inMat(_,1),indexes(_,1),indexes(_,2));
+}

RBBGCMuso/src/musoRandomizer.cpp

@@ -13,6 +13,22 @@ using namespace std;
 
 
 NumericMatrix randTypeZero(NumericMatrix m){
+  /*
+    A typical matrix like m:
+
+    | INDEX | DEPENDENCE | MIN  | MAX  |
+    |-------+------------+------+------|
+    |    21 |          0 | 0    | 364  |
+    |    57 |          0 | bkla | asdf |
+   
+   This randomization type is the easiest,
+   the function produces a matrix which first
+   column contains the indexes, and the second contains
+   random numbers, which was drawn from uniform distribution
+   with the corresponding min and max parameters, specified
+   by the matrix m.
+
+   */
   int n=m.nrow()-1;
   NumericMatrix M(n+1,2);
   M(_,0)=m(_,0);
@@ -26,6 +42,15 @@ NumericMatrix randTypeZero(NumericMatrix m){
 
 // [[Rcpp::export]]
 NumericMatrix randTypeOne(NumericMatrix m){
+  /*
+    The stucture of matrix m here is the same
+    as in randTypeZero.
+
+
+
+   */
+
+  
   NumericVector dependence=m(_,2);
   int n=m.nrow()-1;
   NumericMatrix M(n+1,2);
@@ -43,11 +68,88 @@ NumericMatrix randTypeOne(NumericMatrix m){
 
 
 IntegerVector orderDec(NumericVector v){
+  //This function is order a vector decreasingly
   Function f("order");
   return f(v,_["decreasing"]=1);
 }
 
+// NumericMatrix randTypeTwo(NumericMatrix m){
+//   int n=m.nrow()-1;
+//   int N=n-1;
+//   NumericMatrix mv=m(Range(0,(n-1)),_);
+//   NumericVector dependence=m(_,2);
+//   NumericMatrix M(n+1,2);
+//   M(_,0)=m(_,0);
+//   IntegerVector indexes=orderDec(mv(_,2));
+//   NumericVector sorban=mv(_,2);
+//   sorban.sort(true);
+//   NumericVector sor=cumsum(sorban);
+//   sor.sort(true);
+//   for(int i=0;i<=N;++i){
+//     if(i!=N){
+//       mv((indexes[i]-1),3)-= sor[i+1];
+//     }
+//   }
+  
+//   double rollingNumber=0;
+  
+//   for(int i=0;i<=N;++i){
+//     double minimum=mv((indexes[i]-1),2);
+//     double maximum=mv((indexes[i]-1),3)-rollingNumber;
+//     M(i,1)=runif(1,minimum,maximum)[0];
+//     rollingNumber+=M(i,1);
+//     // cout << "minimum:\t" << minimum << endl;
+//     // cout << "maximum:\t" << maximum << endl;
+//     // cout << "indexes:\t" << indexes[i] << endl;
+//     // cout << "rollingNumber:\t" << rollingNumber << endl;
+//     // cout << "choosen:\t" << M(i,1) <<endl;
+//     // cout << "sor:\t" << sor <<endl;
+//     // cout << "\n\n\n" << mv;
+//   }
+//   M(n,1)=1-rollingNumber;
+//   return M;
+// }
+
+
+//[[Rcpp::export]]
 NumericMatrix randTypeTwo(NumericMatrix m){
+  int n = m.nrow()-1; //Just for optimalization (indexing from 0)
+  int N = n-1; //Just for optimalization (indexing from 0)
+  NumericMatrix mv = m(Range(0,(n-1)),_); // We take off the first n-1 element. What about the last element?????????
+  NumericVector dependence = m(_,2); // Th dependence vector represent the matrix and other information for function based on the elements.
+  NumericMatrix M(n+1,2); // M is the randomizated matrix. It will have n+1 rows (identical to m).
+  M(_,0) = m(_,0); //Insert the indexes
+  IntegerVector indexes = orderDec(mv(_,2)); // Why not use the dependence variable?
+  NumericVector sorban = mv(_,2); // Original ordering of the dependence variable? 
+  sorban.sort(true); // With the previous line, it sorts the dependence variable. 
+  NumericVector sor=cumsum(sorban); // ?????
+  sor.sort(true); // ??????
+  for(int i=0;i<=N;++i){
+    if(i!=N){
+      mv((indexes[i]-1),3)-= sor[i+1];
+    }
+  }
+  
+  double rollingNumber=0;
+  
+  for(int i=0;i<=N;++i){
+    double minimum=mv((indexes[i]-1),2);
+    double maximum=mv((indexes[i]-1),3)-rollingNumber;
+    M(i,1)=runif(1,minimum,maximum)[0];
+    rollingNumber+=M(i,1);
+    // cout << "minimum:\t" << minimum << endl;
+    // cout << "maximum:\t" << maximum << endl;
+    // cout << "indexes:\t" << indexes[i] << endl;
+    // cout << "rollingNumber:\t" << rollingNumber << endl;
+    // cout << "choosen:\t" << M(i,1) <<endl;
+    // cout << "sor:\t" << sor <<endl;
+    // cout << "\n\n\n" << mv;
+  }
+  M(n,1)=1-rollingNumber;
+  return M;
+}
+
+NumericMatrix randTypeThree(NumericMatrix m){
   int n=m.nrow()-1;
   int N=n-1;
   NumericMatrix mv=m(Range(0,(n-1)),_);
@@ -86,7 +188,6 @@ NumericMatrix randTypeTwo(NumericMatrix m){
 
 
 
-
 NumericMatrix copyMatSpec(NumericMatrix A, NumericMatrix B, int u, int v){
   int k=0;
   for(int i=u;i<=v;++i){
@@ -96,6 +197,8 @@ NumericMatrix copyMatSpec(NumericMatrix A, NumericMatrix B, int u, int v){
   return A;
 }
 
+
+
 //' musoRandomizer
 //' 
 //' This function is fastly randomize values based on min and max values,

RBBGCMuso/src/musoRandomizer.cpp~

@@ -0,0 +1,167 @@
+#include <Rcpp.h>
+#include <numeric>
+#include <iostream>
+#include <algorithm>
+#include <numeric>
+#include <ctime>
+
+
+using namespace Rcpp;
+using namespace std;
+// [[Rcpp::plugins(cpp11)]]
+
+
+
+NumericMatrix randTypeZero(NumericMatrix m){
+  /*
+    A typical matrix like m:
+
+    | INDEX | DEPENDENCE | MIN  | MAX  |
+    |-------+------------+------+------|
+    |    21 |          0 | 0    | 364  |
+    |    57 |          0 | bkla | asdf |
+   
+   This randomization type is the easiest,
+   the function produces a matrix which first
+   column contains the indexes, and the second contains
+   random numbers, which was drawn from uniform distribution
+   with the corresponding min and max parameters, specified
+   by the matrix m.
+
+   */
+  int n=m.nrow()-1;
+  NumericMatrix M(n+1,2);
+  M(_,0)=m(_,0);
+  for(int i=0;i<=n;++i){
+    double min=m(i,2);
+    double max=m(i,3);
+    M(i,1)=runif(1,min,max)[0];
+  }
+  return M;
+}
+
+// [[Rcpp::export]]
+NumericMatrix randTypeOne(NumericMatrix m){
+  /*
+    The stucture of matrix m here is the same
+    as in randTypeZero.
+
+
+
+   */
+
+  
+  NumericVector dependence=m(_,2);
+  int n=m.nrow()-1;
+  NumericMatrix M(n+1,2);
+  M(_,0)=m(_,0);
+  M(0,1)=runif(1,m(0,2),m(0,3))[0];
+  for(int i=1;i<=n;++i){
+    int dep=m(i,1)-1;
+    double min=max(M(dep,1),m(i,2));
+    double max=m(i,3);
+    M(i,1)=runif(1,min,max)[0];
+  }
+  return M;
+}
+
+
+
+IntegerVector orderDec(NumericVector v){
+  Function f("order");
+  return f(v,_["decreasing"]=1);
+}
+
+NumericMatrix randTypeTwo(NumericMatrix m){
+  int n=m.nrow()-1;
+  int N=n-1;
+  NumericMatrix mv=m(Range(0,(n-1)),_);
+  NumericVector dependence=m(_,2);
+  NumericMatrix M(n+1,2);
+  M(_,0)=m(_,0);
+  IntegerVector indexes=orderDec(mv(_,2));
+  NumericVector sorban=mv(_,2);
+  sorban.sort(true);
+  NumericVector sor=cumsum(sorban);
+  sor.sort(true);
+  for(int i=0;i<=N;++i){
+    if(i!=N){
+      mv((indexes[i]-1),3)-= sor[i+1];
+    }
+  }
+  
+  double rollingNumber=0;
+  
+  for(int i=0;i<=N;++i){
+    double minimum=mv((indexes[i]-1),2);
+    double maximum=mv((indexes[i]-1),3)-rollingNumber;
+    M(i,1)=runif(1,minimum,maximum)[0];
+    rollingNumber+=M(i,1);
+    // cout << "minimum:\t" << minimum << endl;
+    // cout << "maximum:\t" << maximum << endl;
+    // cout << "indexes:\t" << indexes[i] << endl;
+    // cout << "rollingNumber:\t" << rollingNumber << endl;
+    // cout << "choosen:\t" << M(i,1) <<endl;
+    // cout << "sor:\t" << sor <<endl;
+    // cout << "\n\n\n" << mv;
+  }
+  M(n,1)=1-rollingNumber;
+  return M;
+}
+
+
+
+
+NumericMatrix copyMatSpec(NumericMatrix A, NumericMatrix B, int u, int v){
+  int k=0;
+  for(int i=u;i<=v;++i){
+    A(i,_)=B(k,_);
+    k+=1;
+  }
+  return A;
+}
+
+//' musoRandomizer
+//' 
+//' This function is fastly randomize values based on min and max values,
+//' and row indexes.
+//' @importFrom Rcpp evalCpp
+//' @useDynLib RBBGCMuso
+//' @param A is the big matrix
+//' @param B is the small matrix
+//' @export
+// [[Rcpp::export]]
+NumericMatrix musoRandomizer(NumericMatrix A, NumericMatrix B){
+  NumericMatrix M(A.nrow(),2);
+  int nGroup = B.nrow()-1;
+  int k=0; 
+  for(int i=0;i<=nGroup;++i)
+  {
+    int b=B(i,0)-1;
+    int till=b+k;
+    int t=B(i,1);
+    // cout << b << "\t" << t <<endl;
+    switch(t){
+    case 0:
+      M=copyMatSpec(M,randTypeZero(A(Range(k,till),_)),k,till);
+      // cout << M << endl;
+      break;
+    case 1:
+      M=copyMatSpec(M,randTypeOne(A(Range(k,till),_)),k,till);
+      // cout << M << endl;
+      break;
+    case 2:
+      M=copyMatSpec(M,randTypeTwo(A(Range(k,till),_)),k,till);
+      // cout << M << endl;
+      break;
+    }
+    k=till+1;
+    // cout << k << endl;
+  }
+  return M;
+}
+
+std::string concatenate(std::string A, std::string B){
+  std::string C = A + B;
+  return C;
+}