Skip to content

Commit

Permalink
-Fixed the error of the AUC ploting function
Browse files Browse the repository at this point in the history
-Added a new plot function based on scans
-Added an intensity threshold based on the histogram of the whole dataset
-Minor fixes
  • Loading branch information
LlucSF committed Feb 14, 2019
1 parent 75b8278 commit 4f0104e
Show file tree
Hide file tree
Showing 4 changed files with 230 additions and 152 deletions.
21 changes: 11 additions & 10 deletions NAMESPACE
Original file line number Diff line number Diff line change
@@ -1,10 +1,11 @@
# Generated by roxygen2: do not edit by hand
import(raster)
import(wavelets)
import(ggplot2)
import(sp)
import(flux)
export(DataInput.txt)
export(PlotAUCImage)
export(SmoothSpectra)

# Generated by roxygen2: do not edit by hand
import(raster)
import(wavelets)
import(ggplot2)
import(sp)
import(flux)
export(DataInput.txt)
export(PlotAUCImage)
export(PlotScanImage)
export(SmoothSpectra)

270 changes: 139 additions & 131 deletions R/DataInput.R
Original file line number Diff line number Diff line change
@@ -1,131 +1,139 @@
#########################################################################
# rCRMI - R package for Confocal Raman Spectroscopy Imaging data processing and visualization
# Copyright (C) 2018 Lluc Sementé Fernàndez
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
############################################################################

#' DataInput.txt
#'
#' @param path full path to .txt file.
#' @param numBands Number of .txt files to process.
#'
#' @description Reads a .txt file containing Raman data from Reinshaw devices and convert it to an RamanR data object.
#'
#' @return a rCRSIObj data object.
#' @export
#'

DataInput.txt <- function(path = NULL, numBands = 1)
{
writeLines("Select files...")

#Data input read
PathToFile <- list()
for(i in 1:(numBands))
{
if(is.null(path))
{
PathToFile[[i]] <- file.choose()
}else
{
PathToFile[[i]] <- path[i]
}
}


writeLines("Starting the data structuring process...")

rawData <- list()
for(i in 1:(numBands))
{
tmp <- utils::read.table(file = PathToFile[[i]],header = FALSE)
if(i==1)
{
rawData <- tmp
} else
{
rawData <- cbind(rawData,tmp[,3:4])
}
}

#Raman Shift Axis, number of spectra & Raman band length
RamanShift <- array(dim = c(length(union(rawData[[3]],rawData[[3]])),numBands))
for(i in 1:numBands)
{
p1 <- 3+(i-1)*2
RamanShift[,i] <- union(rawData[[p1]],rawData[[p1]])
}
BandLength <- (length(RamanShift)/numBands)

#Coordenates
X <- union(rawData[[1]],rawData[[1]])
Y <- union(rawData[[2]],rawData[[2]])
numPixels <- length(X)*length(Y)
Coords <- list(X = X, Y = Y)

#Data Shaping
RamanData <- array(dim=c(numBands*numPixels,BandLength))
for(i in 1:numBands)
{
p <- 2*i+2
for(j in 1:numPixels)
{
RamanData[j+numPixels*(i-1),] <- rawData[[p]][((j-1)*BandLength+1):(j*BandLength)]
}
}

writeLines("Data structuring complete...")
return(rCRMIObj(RamanData, numPixels, numBands, RamanShift, BandLength, Coords))
}


rCRMIObj <- function(RamanData, numPixels, numBands = 1, RamanShift, BandsLenght, Coords)
{
rCRSIObj <- list()

rCRSIObj$numPixels <- numPixels
rCRSIObj$numBands <- numBands
rCRSIObj$numSpectr <- numPixels*numBands
rCRSIObj$BandsLength <- BandsLenght
rCRSIObj$AbsCoords <- Coords
rCRSIObj$RelCoords <- list(X = 1:length(Coords$X), Y = 1:length(Coords$Y))
rCRSIObj$PixelSize <- c(Coords$X[2]-Coords$X[1], Coords$Y[2]-Coords$Y[1])
rCRSIObj$Data <- RamanData

# ID <- gsub(x = as.character(Sys.time()),pattern = ("-"),replacement = "")
# ID <- gsub(x = ID,pattern = (":"),replacement = "")
# ID <- gsub(x = ID,pattern = (" "),replacement = "")
# rCRSIObj$ID <- ID
#
# if (numBands == 1 & BandsLenght > 391)
# {
# rCRSIObj$FullSpect <- TRUE
# }
# else rCRSIObj$FullSpect <- FALSE

rCRSIObj$RamanShiftAxis <- RamanShift
rCRSIObj$AvrgSpectr <- RamanShift

for (i in 1:numBands)
{
for (j in 1:numPixels)
{
rCRSIObj$AvrgSpectr[,i] <- 0
}
}
rCRSIObj$AvrgSpectr <- AverageSpectrum(numBands, numPixels, rCRSIObj)
return (rCRSIObj)
}

#########################################################################
# rCRMI - R package for Confocal Raman Spectroscopy Imaging data processing and visualization
# Copyright (C) 2018 Lluc Sementé Fernàndez
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
############################################################################

#' DataInput.txt
#'
#' @param path full path to .txt file.
#' @param numBands Number of .txt files to process.
#'
#' @description Reads a .txt file containing Raman data from Reinshaw devices and convert it to an RamanR data object.
#'
#' @return a rCRSIObj data object.
#' @export
#'

DataInput.txt <- function(path = NULL, numBands = 1, IntensityThreshold = T)
{
writeLines("Select files...")

#Data input read
PathToFile <- list()
for(i in 1:(numBands))
{
if(is.null(path))
{
PathToFile[[i]] <- file.choose()
}else
{
PathToFile[[i]] <- path[i]
}
}


writeLines("Starting the data structuring process...")

rawData <- list()
for(i in 1:(numBands))
{
tmp <- utils::read.table(file = PathToFile[[i]],header = FALSE)
if(i==1)
{
rawData <- tmp
} else
{
rawData <- cbind(rawData,tmp[,3:4])
}
}

#Raman Shift Axis, number of spectra & Raman band length
RamanShift <- array(dim = c(length(union(rawData[[3]],rawData[[3]])),numBands))
for(i in 1:numBands)
{
p1 <- 3+(i-1)*2
RamanShift[,i] <- union(rawData[[p1]],rawData[[p1]])
}
BandLength <- (length(RamanShift)/numBands)

#Coordenates
X <- union(rawData[[1]],rawData[[1]])
Y <- union(rawData[[2]],rawData[[2]])
numPixels <- length(X)*length(Y)
Coords <- list(X = X, Y = Y)

#Data Shaping
RamanData <- array(dim=c(numBands*numPixels,BandLength))
for(i in 1:numBands)
{
p <- 2*i+2
for(j in 1:numPixels)
{
RamanData[j+numPixels*(i-1),] <- rawData[[p]][((j-1)*BandLength+1):(j*BandLength)]
}
}

return(rCRMIObj(RamanData, numPixels, numBands, RamanShift, BandLength, Coords, IntensityThreshold))
}


rCRMIObj <- function(RamanData, numPixels, numBands = 1, RamanShift, BandsLenght, Coords, IntensityThreshold)
{
rCRSIObj <- list()

rCRSIObj$numPixels <- numPixels
rCRSIObj$numBands <- numBands
rCRSIObj$numSpectr <- numPixels*numBands
rCRSIObj$BandsLength <- BandsLenght
rCRSIObj$AbsCoords <- Coords
rCRSIObj$RelCoords <- list(X = 1:length(Coords$X), Y = 1:length(Coords$Y))
rCRSIObj$PixelSize <- c(Coords$X[2]-Coords$X[1], Coords$Y[2]-Coords$Y[1])
rCRSIObj$Data <- RamanData
rCRSIObj$RamanShiftAxis <- RamanShift
rCRSIObj$AvrgSpectr <- RamanShift

for (i in 1:numBands)
{
for (j in 1:numPixels)
{
rCRSIObj$AvrgSpectr[,i] <- 0
}
}


rCRSIObj$AvrgSpectr <- AverageSpectrum(numBands, numPixels, rCRSIObj)

if(IntensityThreshold)
{
dens <- density(rCRSIObj$Data)
aucvec <- vector(length = (length(dens$x)-1))
for(i in 2:length(dens$x))
{
aucvec[i-1] <- flux::auc(x = dens$x[1:i], y = dens$y[1:i])
}

intThr <- dens$x[which(aucvec>0.99*max(aucvec))[1]]

for(i in 1:rCRSIObj$BandsLength)
{
rCRSIObj$Data[which(rCRSIObj$Data[,i]>intThr),i] <- median(rCRSIObj$Data[,i])
}
}

writeLines("Data structuring complete...")
return (rCRSIObj)
}

24 changes: 13 additions & 11 deletions R/PlotAUCimg.R
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@
#' @export
#'

PlotAUCImage <- function(rCRSIObj, intervals, numImg = 1)
PlotAUCImage <- function(rCRSIObj, intervals, numImg = 1, interpolate = T)
{
for(i in 1:numImg)
{
Expand All @@ -48,18 +48,20 @@ PlotAUCImage <- function(rCRSIObj, intervals, numImg = 1)
from <- which.min(abs(rCRSIObj$RamanShiftAxis-intervals[2*i-1]))
to <- which.min(abs(rCRSIObj$RamanShiftAxis-intervals[2*i]))

AUCimg <- matrix(nrow = length(rCRSIObj$RelCoords$Y), ncol = length(rCRSIObj$RelCoords$X))
AUCvec <- rep(0, times = rCRSIObj$numPixels)

for(X in 1:ncol(AUCimg))
for(i in 1:rCRSIObj$numPixels)
{
for(Y in 1:nrow(AUCimg))
{
AUCimg[Y,X] <- flux::auc(x = rCRSIObj$RamanShiftAxis[from:to],
y = rCRSIObj$Data[(Y+(X-1)*length(rCRSIObj$RelCoords$Y)), from:to]
)
}
AUCvec[i] <- flux::auc(x = rCRSIObj$RamanShiftAxis[from:to],
y = rCRSIObj$Data[i,from:to]
)
}

AUCimg <- matrix(data = AUCvec,
nrow = length(rCRSIObj$RelCoords$X),
ncol = length(rCRSIObj$RelCoords$Y))


#Create the raster
img <- raster::raster(nrows = nrow(AUCimg),
ncols = ncol(AUCimg),
Expand All @@ -69,12 +71,12 @@ PlotAUCImage <- function(rCRSIObj, intervals, numImg = 1)
ymx = nrow(AUCimg)
)

raster::values(img) <- AUCimg
img <- raster::setValues(img,AUCimg)

raster::plot(x = img,
col = sp::bpy.colors(64),
useRaster = T,
interpolate = T,
interpolate = interpolate
)

graphics::title(main = paste("AUC: ",signif(rCRSIObj$RamanShiftAxis[from],5),"~",signif(rCRSIObj$RamanShiftAxis[to],5)),
Expand Down
Loading

0 comments on commit 4f0104e

Please sign in to comment.