R/classifyQA.R

In RStoolbox: Remote Sensing Data Analysis

#' Classify Landsat QA bands
#' 
#' extracts five classes from QA band: background, cloud, cirrus, snow and water.
#' 
#' @param img SpatRaster. Landsat 8 OLI QA band.
#' @param type Character. Classes which should be returned. One or more of c("background", "cloud", "cirrus","snow", "water").
#' @param confLayers Logical. Return one layer per class classified by confidence levels, i.e. cloud:low, cloud:med, cloud:high.
#' @param ... further arguments passed to \link[terra]{writeRaster}
#' @inheritParams encodeQA
#' 
#' @seealso \link{encodeQA} \link{decodeQA}
#' @details 
#' By default each class is queried for *high* confidence. See \link{encodeQA} for details. To return the different confidence levels per condition use \code{confLayers=TRUE}.
#' This approach corresponds to the way LandsatLook Quality Images are produced by the USGS.
#' @return Returns a SpatRaster with maximal five classes:
#' \tabular{rr}{
#' class \tab value \cr
#' background \tab 1L \cr 
#' cloud  \tab 2L \cr
#' cirrus \tab 3L \cr
#' snow   \tab 4L \cr
#' water  \tab 5L \cr 
#' }
#' Values outside of these classes are returned as NA.
#' If \code{confLayers = TRUE} then a RasterStack with one layer per condition (except 'background') is returned, whereby each layer contains the confidence level of the condition.
#' \tabular{rr}{
#' Confidence \tab value \cr
#' low    \tab 1L \cr
#' med    \tab 2L \cr
#' high   \tab 3L \cr
#' }
#' @export 
#' @examples
#' library(terra)
#' qa <- rast(ncol = 100, nrow=100, val = sample(1:2^14,  10000))
#' 
#' ## QA classes
#' qacs <- classifyQA(img = qa)
#' ## Confidence levels
#' qacs_conf <- classifyQA(img = qa, confLayers = TRUE)
classifyQA <- function(img, type = c("background", "cloud", "cirrus","snow", "water"), confLayers = FALSE, sensor = "OLI", legacy = "collection1", ...){
  img <- .toTerra(img)

  ## Input checks
  if(legacy == "pre_collection" & !any(sensor %in% c("OLI", "TIRS"))) stop("For argument legacy = 'pre_collection', argument sensor can only be 'OLI' or 'TIRS'.", call.=FALSE)
  
  if(any(!type %in% c("background", "cloud", "cirrus","snow", "water")) | !length(type)) stop("type must be element of c('background', 'cloud', 'cirrus','snow', 'water')")
  if(nlyr(img) != 1) stop("img should be a single layered SpatRaster")

  if(!confLayers){
    rclx <- rbind(
      if("background" %in% type) cbind(is = encodeQA(fill = "yes", sensor = sensor, legacy = legacy), becomes = 1L),
      if("cloud" %in% type)  cbind(is = encodeQA(cloud = "high", sensor = sensor, legacy = legacy),   becomes = 2L),
      if("cirrus" %in% type) cbind(is = encodeQA(cirrus = "high", sensor = sensor, legacy = legacy),  becomes = 3L),
      if("snow" %in% type)   cbind(is = encodeQA(snow = "high", sensor = sensor, legacy = legacy),    becomes = 4L),
      if("water" %in% type)  cbind(is = encodeQA(water = "high", sensor = sensor, legacy = legacy),   becomes = 5L)
    )

    out <- .paraRasterFun(
      img,
      rasterFun = app,
      args = list(
        fun = function(xi, na.rm = FALSE) classQA(x = xi, rcl = rclx)
      ),
      wrArgs = list(...)
    )

    names(out) <- "QAclass"
    return(out)
  } else {
    
    rclxList <- list(
      cloud = rbind(            
        cbind(is = encodeQA(cloud = "low", sensor = sensor, legacy = legacy),  becomes = 1L),
        cbind(is = encodeQA(cloud = "med", sensor = sensor, legacy = legacy),  becomes = 2L),
        cbind(is = encodeQA(cloud = "high", sensor = sensor, legacy = legacy),  becomes = 3L)),
      
      cirrus = rbind(
        cbind(is = encodeQA(cirrus = "low", sensor = sensor, legacy = legacy), becomes = 1L),
        cbind(is = encodeQA(cirrus = "med", sensor = sensor, legacy = legacy), becomes = 2L),
        cbind(is = encodeQA(cirrus = "high", sensor = sensor, legacy = legacy), becomes = 3L)),
      
      snow =  rbind(        
        cbind(is = encodeQA(snow = "low", sensor = sensor, legacy = legacy), becomes = 1L),
        cbind(is = encodeQA(snow = "med", sensor = sensor, legacy = legacy), becomes = 2L),        
        cbind(is = encodeQA(snow = "high", sensor = sensor, legacy = legacy), becomes = 3L)),
      
      water = rbind(
        cbind(is = encodeQA(water = "low", sensor = sensor, legacy = legacy),  becomes = 1L),
        cbind(is = encodeQA(water = "med", sensor = sensor, legacy = legacy),  becomes = 2L),
        cbind(is = encodeQA(water = "high", sensor = sensor, legacy = legacy),  becomes = 3L)))

    out <- lapply(type[type != "background"], function(i){
      .paraRasterFun(
        img,
        rasterFun = app,
        args = list(
          fun = function(xi, na.rm = FALSE) classQA(x = xi, rcl = rclxList[[i]])
        )
      )
    })

    names(out) <- type[type != "background"]
    return(out)
  }
    
}