--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/nmr_bucketing/NmrBucketing_wrapper.R	Thu Apr 20 08:53:46 2017 -0400
@@ -0,0 +1,295 @@
+#!/usr/local/public/bin/Rscript --vanilla --slave --no-site-file
+
+## 070115_NmrBucketing2galaxy_v1.R
+## Marie Tremblay-Franco
+## MetaboHUB: The French Infrastructure for Metabolomics and Fluxomics
+## www.metabohub.fr/en
+## marie.tremblay-franco@toulouse.inra.fr
+
+runExampleL <- FALSE
+
+if(runExampleL) {
+##------------------------------
+## Example of arguments
+##------------------------------
+argLs <- list(StudyDir = "Tlse_BPASourisCerveau",
+              upper = "10.0",
+              lower = "0.50",
+              bucket.width = "0.01",
+              exclusion = "TRUE",
+              exclusion.zone = list(c(6.5,4.5)),
+              graph="Overlay")
+
+argLs <- c(argLs,
+           list(dataMatrixOut = paste(directory,"_NmrBucketing_dataMatrix.tsv",sep=""),
+                sampleMetadataOut = paste(directory,"_NmrBucketing_sampleMetadata.tsv",sep=""),
+                variableMetadataOut = paste(directory,"_NmrBucketing_variableMetadata.tsv",sep=""),
+                graphOut = paste(directory,"_NmrBucketing_graph.pdf",sep=""),
+                logOut = paste(directory,"_NmrBucketing_log.txt",sep="")))
+}
+
+##------------------------------
+## Options
+##------------------------------
+strAsFacL <- options()$stringsAsFactors
+options(stringsAsFactors = FALSE)
+
+
+##------------------------------
+## Libraries laoding
+##------------------------------
+# For parseCommandArgs function
+library(batch)
+# For cumtrapz function
+library(pracma)
+
+# R script call
+source_local <- function(fname)
+{
+	argv <- commandArgs(trailingOnly = FALSE)
+	base_dir <- dirname(substring(argv[grep("--file=", argv)], 8))
+	source(paste(base_dir, fname, sep="/"))
+}
+#Import the different functions
+source_local("NmrBucketing_script.R")
+source_local("DrawSpec.R")
+
+##------------------------------
+## Errors ?????????????????????
+##------------------------------
+
+
+##------------------------------
+## Constants
+##------------------------------
+topEnvC <- environment()
+flagC <- "\n"
+
+
+##------------------------------
+## Script
+##------------------------------
+if(!runExampleL)
+    argLs <- parseCommandArgs(evaluate=FALSE)
+
+
+## Parameters Loading
+##-------------------
+  # Inputs
+if (!is.null(argLs[["zipfile"]])){
+	fileType="zip"
+	zipfile= argLs[["zipfile"]]
+	directory=unzip(zipfile, list=F)
+	directory=paste(getwd(),strsplit(directory[1],"/")[[1]][2],sep="/")
+} else if (!is.null(argLs[["tsvfile"]])){
+	fileType="tsv"
+	directory <- read.table(argLs[["tsvfile"]],check.names=FALSE,header=TRUE,sep="\t")
+}
+
+leftBorder <- argLs[["left_border"]]
+rightBorder <- argLs[["right_border"]]
+bucketSize <- argLs[["bucket_width"]]
+exclusionZones <- argLs[["zone_exclusion_choices.choice"]]
+
+exclusionZonesBorders <- NULL
+if (!is.null(argLs$zone_exclusion_left))
+{
+   for(i in which(names(argLs)=="zone_exclusion_left"))
+   {
+     exclusionZonesBorders <- c(exclusionZonesBorders,list(c(argLs[[i]],argLs[[i+1]])))
+   }
+}
+
+graphique <- argLs[["graphType"]]
+
+  # Outputs
+nomGraphe <- argLs[["graphOut"]]
+dataMatrixOut <- argLs[["dataMatrixOut"]]
+logFile <- argLs[["logOut"]]
+if (fileType=="zip")
+{
+  sampleMetadataOut <- argLs[["sampleOut"]]
+  variableMetadataOut <- argLs[["variableOut"]]
+}
+
+## Checking arguments
+##-------------------
+error.stock <- "\n"
+
+if(length(error.stock) > 1)
+  stop(error.stock)
+
+
+## Computation
+##------------
+outputs <- NmrBucketing(fileType=fileType, fileName=directory, leftBorder=leftBorder, rightBorder=rightBorder, bucketSize=bucketSize,
+						exclusionZones=exclusionZones, exclusionZonesBorders=exclusionZonesBorders, graph=graphique, nomFichier=nomGraphe,
+						savLog.txtC=logFile)
+data_bucket <- outputs[[1]]
+data_sample <- outputs[[2]]
+data_variable <- outputs[[3]]
+ppm <- outputs[[4]]
+ppm <- round(ppm,2)
+
+## Graphical outputs
+##------------------
+if (graphique != "None")
+{
+  excludedZone <- NULL
+  for (c in 1:length(exclusionZonesBorders))
+  {
+    excludedZone <- c(excludedZone,exclusionZonesBorders[[c]])
+    excludedZone <- sort(excludedZone)
+  }
+  nbZones <- length(excludedZone)/2
+  n <- length(excludedZone)
+  
+  # Graphic Device opening
+  pdf(nomGraphe,onefile=TRUE)
+  
+  if (graphique == "Overlay")
+  {
+    # Global spectral window
+    spectra <- data.frame(t(data_bucket))
+    drawSpec(spectra,xlab="", ylab="Intensity", main="")
+    
+    ## Zoomed spectral window depending on exclusion zone(s)
+    if (nbZones != 0)
+    {
+      BInf <- excludedZone[n]
+      if (round(BInf,1) == BInf)
+      {
+        BInf <- BInf+0.01
+      }
+      spectra <- data.frame(t(data_bucket[1:(which(ppm == BInf)[[1]]),]))
+      drawSpec(spectra,xlab="", ylab="Intensity", main="")			
+      n <- n - 1
+      
+      while (n >= nbZones & nbZones > 1)
+      {
+        BInf <- excludedZone[n-1]
+        if (round(BInf,1) > BInf)
+        {
+          BInf <- BInf+0.01
+        }
+        spectra <- data.frame(t(data_bucket[(which(ppm == excludedZone[n])[[1]]):(which(ppm == BInf)[[1]]),]))
+        drawSpec(spectra,xlab="", ylab="Intensity", main="")
+        n <- n - 2
+      }
+      
+      BInf <- excludedZone[1]
+      if (round(BInf,1) <= BInf)
+      {
+        BInf <- BInf+0.01
+      }
+      spectra <- data.frame(t(data_bucket[(which(ppm == BInf)[[1]]):nrow(data_bucket),]))
+      drawSpec(spectra,xlab="", ylab="Intensity", main="")
+    }
+  }
+  else
+  {
+    for (i in 1:ncol(data_bucket))
+    {
+      par(mfrow=c((nbZones+2),1))
+      n <- length(excludedZone)
+      spectra <- t(data_bucket[,i])
+	  names(spectra) <- rownames(data_bucket)
+      plot(1:length(spectra), spectra, type='l', xlab="", ylab="Intensity", main=colnames(data_bucket)[i], xaxt = "n")
+	  xPos <- 1
+	  nAxisPos <- 4
+	  startP <- length(nAxisPos) 
+	  endP <- nrow(data_bucket)
+	  GraphRange <- c(startP:endP)
+	  tempVal = trunc(length(GraphRange)/nAxisPos)
+	  xPos = c(0:nAxisPos) * tempVal
+	  axis(1, at = xPos, labels = rownames(data_bucket)[xPos + startP])
+     
+      ## Zoomed spectral window depending on exclusion zone(s)
+      if (nbZones != 0)
+      {
+        BInf <- excludedZone[n]
+        if (round(BInf,1) == BInf)
+        {
+          BInf <- BInf+0.01
+        }
+        spectra <- t(data_bucket[1:(which(ppm == BInf)[[1]]),i])
+		names(spectra) <- rownames(data_bucket)[1:(which(ppm == BInf)[[1]])]
+		plot(1:length(spectra), spectra, type='l',xlab="", ylab="Intensity", main="", xaxt = "n")			
+		xPos <- 1
+		nAxisPos <- 4
+		startP <- length(nAxisPos) 
+		endP <- length(spectra)
+		GraphRange <- c(startP:endP)
+		tempVal = trunc(length(GraphRange)/nAxisPos)
+		xPos = c(0:nAxisPos) * tempVal
+		axis(1, at = xPos, labels = rownames(data_bucket)[xPos + startP])
+        n <- n - 1
+        
+        while (n >= nbZones & nbZones > 1)
+        {
+          BInf <- excludedZone[n-1]
+          if (round(BInf,1) > BInf)
+          {
+            BInf <- BInf+0.01
+          }
+          spectra <- t(data_bucket[(which(ppm == excludedZone[n])[[1]]):(which(ppm == BInf)[[1]]),i])
+		  names(spectra) <- rownames(data_bucket)[(which(ppm == excludedZone[n])[[1]]):(which(ppm == BInf)[[1]])]
+          plot(1:length(spectra), spectra, type='l',xlab="", ylab="Intensity", main="", xaxt = "n")
+		  xPos <- 1
+		  nAxisPos <- 4
+		  startP <- length(nAxisPos) 
+		  endP <- length(spectra)
+		  GraphRange <- c(startP:endP)
+		  tempVal = trunc(length(GraphRange)/nAxisPos)
+		  xPos = c(0:nAxisPos) * tempVal
+		  axis(1, at = xPos, labels = rownames(data_bucket)[xPos + startP])
+          n <- n - 2
+        }
+        
+        BInf <- excludedZone[1]
+        if (round(BInf,1) <= BInf)
+        {
+          BInf <- BInf+0.01
+        }
+        spectra <- t(data_bucket[(which(ppm == BInf)[[1]]):nrow(data_bucket),i])
+		names(spectra) <- rownames(data_bucket)[(which(ppm == BInf)[[1]]):nrow(data_bucket)]
+        plot(1:length(spectra), spectra, type='l',xlab="", ylab="Intensity", main="", xaxt = "n")
+		xPos <- 1
+		nAxisPos <- 4
+		startP <- length(nAxisPos) 
+		endP <- length(spectra)
+		GraphRange <- c(startP:endP)
+		tempVal = trunc(length(GraphRange)/nAxisPos)
+		xPos = c(0:nAxisPos) * tempVal
+		axis(1, at = xPos, labels = rownames(data_bucket)[xPos + startP])
+      }
+    }
+  }
+  dev.off()
+}
+## Saving
+##-------
+  # Data
+data_bucket <- cbind(rownames(data_bucket),data_bucket)
+colnames(data_bucket) <- c("Bucket",colnames(data_bucket)[-1])
+write.table(data_bucket,file=argLs$dataMatrixOut,quote=FALSE,row.names=FALSE,sep="\t")
+  # Sample
+data_sample <- cbind(rownames(data_sample),data_sample)
+colnames(data_sample) <- c("Sample",colnames(data_sample)[-1])
+write.table(data_sample,file=argLs$sampleOut,quote=FALSE,row.names=FALSE,sep="\t")
+  # Variable
+data_variable <- cbind(rownames(data_variable),data_variable)
+colnames(data_variable) <- c("Bucket",colnames(data_variable)[-1])
+write.table(data_variable,file=argLs$variableOut,quote=FALSE,row.names=FALSE,sep="\t")
+
+
+## Ending
+##---------------------
+
+cat("\nEnd of 'NMR bucketing' Galaxy module call: ", as.character(Sys.time()), sep = "")
+
+## sink(NULL)
+
+options(stringsAsFactors = strAsFacL)
+
+rm(list = ls())