R/display_2D_histo.R

In espadon: Easy Study of Patient DICOM Data in Oncology

#' Display of a 2D histogram
#' @description The \code{display.2D.histo} function displays the density map of
#' a "histo2D" class object.
#' @param histo.2D "histo2D" class object.
#' @param add Boolean indicating whether to display the background image.
#' @param main Title of the background image. If \code{main = NULL},
#' the title indicates "2D histogram".
#' @param x.lab Label of the x-axis of the background image. If \code{x.lab = NULL},
#' this label is \code{histo.2D$x.file.src}
#' @param y.lab Label of the y-axis of the background image. If \code{y.lab = NULL},
#' this label is \code{histo.2D$y.file.src}.
#' @param x.lim Vector, representing the range of the x-axis.
#' @param y.lim Vector, representing the range of the y-axis.
#' @param bg Background color of the image. By default, this color is black.
#' @param i.rng Vector of 2 elements giving the minimum and the maximum intensity
#' of the image. If \code{i.rng = NULL}, then the minimum is 0 and the maximum the
#' maximum density.
#' @param display.mode function display mode. See Details.
#' @param col Color of the color gradient in \code{display.mode = "mono.color"}, 
#' or of the level lines in \code{display.mode = "line"}. By default, this color 
#' is white.
#' @param alpha A number from 0 to 1, indicating the opacity of the image in
#' \code{"rainbow.color"} mode. Default alpha = 1.
#' @param line.pc.levels Vector of level lines in percent of maximum density
#' in \code{display.mode = "line"}. By default lines 1% and 100% are displayed.
#' @param line.lwd Line thickness of the level lines in \code{display.mode = "line"}.
#' @param line.lty Type of lines for level lines in \code{display.mode = "line"}.
#' @details The \code{display.mode} argument can be set to three values: 
#' \code{"mono.color"}, \code{"rainbow.color"}, or \code{"line"}. The 2D 
#' histogram graph is displayed by default in \code{"mono.color"} mode.
#' \itemize{
#' \item The \code{"mono.color"} mode displays a gradient of the color defined
#' by the col argument, depending on the intensity of \code{$density.map} 2-dimensional
#' array.
#' \item The \code{"rainbow.color"} mode makes a display according to the
#' \code{"rainbow"} palette, while managing the opacity of the colors.
#' \item The \code{"line"} mode draws level lines defined in percent by the
#' \code{line.pc.levels} argument.
#' }
#' @return Returns a display of the density map of \code{histo.2D}. This one must 
#' be an object of class "histo2D". See \link[espadon]{espadon.class} for 
#' class definitions.
#' @seealso \link[espadon]{histo.2D}.
#' @examples
#' # loading of toy-patient objects (decrease dxyz for better result)
#' step <- 4
#' patient <- toy.load.patient (modality = c("ct", "mr", "rtstruct"), 
#'                              roi.name =  "brain", 
#'                              dxyz = rep (step, 3))
#' CT <- patient$ct[[1]]
#' MR <- patient$mr[[1]]
#' S <- patient$rtstruct[[1]]
#' T.MAT <- patient$T.MAT
#'
#' # restriction of the volume around the RoI
#' CT.on.roi <- nesting.roi (CT, S, roi.name = "brain", vol.restrict = TRUE,
#'                           xyz.margin = c (1, 1, 1), alias = CT$description)
#' MR.on.CT <- vol.regrid (vol = MR, back.vol = CT.on.roi, interpolate = TRUE,
#'                         T.MAT = T.MAT, alias = CT$description,
#'                         description = NULL)

#' # selection of voxels included in the RoI.
#' roi.bin <- bin.from.roi (vol = CT.on.roi, struct = S, roi.sname = "brain",
#'                          verbose = FALSE)
#' MR.select <- vol.from.bin (MR.on.CT, roi.bin, alias = MR$description)
#' CT.select <- vol.from.bin (CT.on.roi, roi.bin, alias = CT$description)

#' # 2D histogram
#' H2D <- histo.2D (MR.select, CT.select, x.breaks = seq (50, 400, 10),
#' 			  y.breaks = seq (50, 400, 10), alias = "H2D MR1 MR2")
#' display.2D.histo (H2D, display.mode = "mono.color", col = "#ffff00", 
#'                   main ="mono color mode")
#' display.2D.histo (H2D, display.mode = "rainbow.color", main ="rainbow mode")
#' display.2D.histo (H2D, display.mode = "line", main ="level lines mode",
#'                   line.pc.levels = c (0, 25, 50, 75, 100), col = "#ff0000")

#' @importFrom stats ecdf
#' @importFrom grDevices contourLines rainbow
#' @importFrom colorspace hex2RGB
#' @importFrom methods is
#' @export
display.2D.histo  <- function (histo.2D, add=TRUE, main= NULL, 
                               x.lab = NULL, y.lab = NULL, 
                               x.lim=NULL, y.lim = NULL,
                               bg="#000000", i.rng = NULL, 
                               display.mode=c("mono.color","rainbow.color","line"),
                               col= "#ffffff", alpha = 1,
                               line.pc.levels = c (1, 100),
                               line.lwd = 2, line.lty = 1){
  
  
  .apply.color.map <- function (M, alpha) {
    p <- rev(rainbow(256, s = seq(256, 1, by=-1)/256, v = 1, start = 0, end = 4/6, alpha = alpha * seq (256, 1, by=-1)/256))
    p[M * 255 + 1]
  }
  
  if (!is (histo.2D, "histo2D")) stop ("histo2D must be a histo2D class object.")

  

  scale1 <- histo.2D$x.mids[c(1,length(histo.2D$x.mids))]
  
  scale2 <- histo.2D$y.mids[c(1,length(histo.2D$y.mids))]
  
  if(is.null(main)) main <-"2D histogram"
  if(is.null(x.lab)) x.lab <- histo.2D$x.file.src
  if(is.null(y.lab)) y.lab <- histo.2D$y.file.src
  abs.left <- histo.2D$x.breaks[1]
  ord.bottom <- histo.2D$y.breaks[1]
  abs.right <- histo.2D$x.breaks[length(histo.2D$x.breaks)]
  ord.top <- histo.2D$y.breaks[length(histo.2D$y.breaks)]
 
  xlim <- suppressWarnings(range (x.lim, na.rm=T))
  ylim <- suppressWarnings(range (y.lim, na.rm=T))
  if (any(is.infinite(xlim)))  xlim <- c(abs.left,abs.right)
  if (any(is.infinite(ylim)))  ylim <- c(ord.bottom,ord.top)
  if (add){
    plot (scale1, scale2, type="n", xlim = xlim, ylim = ylim, #asp = 1,
          main=main ,xlab=x.lab, ylab = y.lab, xaxs="i", yaxs="i")
    rect(par("usr")[1], par("usr")[3], par("usr")[2], par("usr")[4], col = bg)
  }
  
  
  if (!is.null(i.rng)) {
    imin <- i.rng[1]
    imax <- i.rng[2]
  } else {
    imin <- 0
    imax <- max(histo.2D$density.map, na.rm=TRUE)
  }
  
  
  map.layer <- (t(histo.2D$density.map) -imin)/(imax - imin)
  map.layer[map.layer>1] <- 1
  map.layer[map.layer<0] <- 0
  
  display.mode  <- display.mode[1]
  if (display.mode[1]=="rainbow.color"){
    
    map.layer <- matrix (.apply.color.map(map.layer, alpha), nrow=nrow(map.layer))
    rasterImage (map.layer[nrow(map.layer):1,], xleft = abs.left,
                 ybottom = ord.bottom, xright = abs.right, ytop = ord.top,
                 interpolate=FALSE)
  } else if (display.mode[1]=="line") {
    
    # z <- histo.2D
    imax <- max(histo.2D$density.map, na.rm=TRUE)
    z <-(histo.2D$density.map/imax) *100
    z[is.na(z)]<- 0
    ctr <- contourLines(x = histo.2D$x.mids, y= histo.2D$y.mids ,
                        z=z, levels = line.pc.levels)
    if (length(ctr)>0) for (i in 1:length(ctr)) {
      lines(c (ctr[[i]]$x,ctr[[i]]$x[1]), 
            c (ctr[[i]]$y,ctr[[i]]$y[1]), col=col, lwd=line.lwd, lty=line.lty)
    }
  } else {
    
    
    ramp<-seq(0,1,1/256)
    full.scale <- hex2RGB(col[1])@coords - hex2RGB(bg[1])@coords
    col.pal <- sapply(ramp, function(k) round((full.scale*k + hex2RGB(bg[1])@coords)*255))
    
    pal <-sapply(1:ncol(col.pal), function (c.idx) {
      vect <- as.character (as.hexmode(col.pal[,c.idx]))
      paste ("#",paste(sapply(vect, function(e)
        ifelse(nchar(e)==1,paste("0",e,sep=""),e)),collapse=""),sep="")
    })
    map.layer <- matrix (pal[round(map.layer*255 +1)], nrow=nrow(map.layer))
    rasterImage (map.layer[nrow(map.layer):1,], xleft = abs.left,
                 ybottom = ord.bottom, xright = abs.right, ytop = ord.top,
                 interpolate=FALSE)
    
  }
  
  
  
}