restructure calib and plotMuso -hugh work

RBBGCMuso/R/calibMuso.R

@@ -16,13 +16,14 @@
 #' @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,parameters=NULL, timee="d", debugging=FALSE, logfilename=NULL, keepEpc=FALSE, export=FALSE, silent=FALSE, aggressive=FALSE, leapYear=FALSE,keepBinary=FALSE, binaryPlace="./", fileToChange="epc"){
+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){
 
 
 ##########################################################################
@@ -109,7 +110,7 @@ calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE,
 
     ##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(!skipSpinup) {
 
     ##Run the model for the spinup run.
 
@@ -158,7 +159,7 @@ calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE,
             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.
@@ -218,43 +219,60 @@ calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE,
     }
 
     
-    logfiles <- tryCatch(getLogs(outputLoc,outputNames,type="both"),
+    if(skipSpinup){
-                        error = function (e){
+       logfiles <- tryCatch(getLogs(outputLoc,outputNames,type="normal"),
-                            setwd(whereAmI)
+                                error = function (e){
-                            stop("Cannot find log files, something went wrong")})
+                                    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")
-    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
+        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(spincrash){
+            if(length(perror)==1){
                 errorsign <- 1
             } else {
-                errorsign <- 0
+                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)
 
         }
@@ -263,7 +281,8 @@ calibMuso <- function(settings=NULL,parameters=NULL, timee="d", debugging=FALSE,
 
 
     if(debugging){ #debugging is boolean
-                       logfiles <- file.path(outputLoc,logfiles)
+        logfiles <- file.path(outputLoc,logfiles)
+        
                        stampAndDir(stampDir=dirName, wrongDir=dirERROR, names=logfiles, type="general",errorsign=errorsign,logfiles=logfiles)}
   
     

RBBGCMuso/R/plotMuso.R

@@ -13,13 +13,15 @@
 #' @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 graphics
+#' @import ggplot2, dplyr
 #' @export
 
 plotMuso <- function(settings=NULL,
@@ -33,50 +35,165 @@ plotMuso <- function(settings=NULL,
                      logfilename=NULL,
                      aggressive=FALSE,
                      leapYear=FALSE,
-                     export=FALSE){
+                     plotName=NULL,
+                     plotType="cts",
+                     layerPlot=FALSE,
+                     colour="blue",
+                     skipSpinup=TRUE,
+                     fromData=FALSE,
+                     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()
     }
     
-    musoData <- rungetMuso(settings=settings,
+    numberOfYears <- settings$numYears
-                           silent=silent,
+    startYear <- settings$startYear
-                           timee=timee,
+    ## musoData <- rungetMuso(settings=settings,
-                           debugging=debugging,
+    ##                        silent=silent,
-                           keepEpc=keepEpc,
+    ##                        timee=timee,
-                           logfilename=logfilename,
+    ##                        debugging=debugging,
-                           export=export) 
+    ##                        keepEpc=keepEpc,
+    ##                        logfilename=logfilename,
+    ##                        export=export)
+    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 {
+            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)) %>%
+                tidyr::separate(date2,c("day","month","year"),sep="\\.")
+    }
 
-     xlab_muso<- switch(timee, "d"="days","y"="years","m"="months")
+    ## numVari <- ncol(musoData)
-     numVari <- ncol(musoData)
+     numVari <- ncol(musoData)-5
 
-     if(is.numeric(variable)){
+    pointOrLineOrPlot <- function(musoData, variableName, plotType="cts", expandPlot=FALSE, plotName=NULL){
-     if((variable>numVari)|(variable<1)){
+        if(!expandPlot){
-         return(print(paste("The variable parameter must be between 1 and ",numVari)))
+            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 <- baseData %>%
+                        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 <- baseData %>%
+                        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(colour=colour, aes_string("date",variableName))
+                    
+                } else {
+                    stop("you cannot add layers for multiple plots")
+                }
+            } else {
+                if(length(variableName)==1){
+                    geom_point(colour=colour, aes_string("date",variableName))
+                } else{
+                    stop("you cannot add layers for multiple plots")
+                }
+            }
+            
+        }
+        
+    }
     
 
-     plot(musoData[,variable],pch=20,col = "dark blue",xlab=xlab_muso,ylab=colnames(musoData)[variable])
+
+    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(variable=="all"){
          
-             musoData  <- rbind(1:numVari,musoData) #creating the column indexes
+         if(prod(sapply(variable,function(x){
-             par(mfrow = niceMatrixLayoutForPlots(numVari))
+             return(x >= 0 && x <= numVari)
-
+         }))){
-             apply(musoData, 2, function(x) plot(x[2:length(x)],pch=20,col="dark blue",xlab=xlab_muso,ylab = colnames(musoData)[x[1]]))
+             variableName <-  as.character(settings$outputVars[[1]])[variable]
-             par(mfrow=c(1,1))
-             return(print("Everything was Ok. ;)"))
          } else {
-             return(print("The variable option is the coloumn number of the output data-matrix, so it must be numeric, of if you want to plot the whole data matrix set it \"all\""))
+             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 
@@ -95,7 +212,7 @@ plotMuso <- function(settings=NULL,
 #' leapYear=FALSE, export=FALSE)
 #' @import ggplot2
 #' @export
-plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",", savePlot=NULL){
+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")
     }
@@ -121,8 +238,8 @@ plotMusoWithData <- function(csvFile, variable, NACHAR=NA, settings=NULL, sep=",
 
     p <- baseData %>%
         ggplot(aes_string("date",variable)) +
-        geom_line(colour = "blue") +
+        geom_line(colour=colour[1]) +
-        geom_point(aes(date,measuredData)) +
+        geom_point(colour=colour[2], aes(date,measuredData)) +
         labs(y = paste0(variable,"_measured"))+
         theme(axis.title.x = element_blank())
     if(!is.null(savePlot)){