diff --git a/DESCRIPTION b/DESCRIPTION index 92a566250..4d856dba7 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,5 +1,5 @@ Package: aqp -Version: 2.0.1 +Version: 2.0.2 Title: Algorithms for Quantitative Pedology Authors@R: c(person(given="Dylan", family="Beaudette", role = c("aut", "cre"), email = "dylan.beaudette@usda.gov"), person(given="Pierre", family="Roudier", email="roudierp@landcareresearch.co.nz", role = c("aut", "ctb")), person(given="Andrew", family="Brown", email="andrew.g.brown@usda.gov", role = c("aut", "ctb"))) Author: Dylan Beaudette [aut, cre], Pierre Roudier [aut, ctb], Andrew Brown [aut, ctb] diff --git a/NAMESPACE b/NAMESPACE index 1045acd4c..96351c839 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -23,6 +23,7 @@ export(brierScore) export(buntley.westin.index) export(checkHzDepthLogic) export(checkSPC) +export(col2Munsell) export(colorChart) export(colorContrast) export(colorContrastPlot) diff --git a/NEWS.md b/NEWS.md index f50a190cd..9964b1530 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,4 +1,8 @@ +# aqp 2.0.2 (2023-09-14) + * new function `col2Munsell()` generalizes and replaces `rgb2munsell()` (thanks Shawn Salley for the suggestion) + # aqp 2.0.1 (2023-09-03) + * CRAN release (CRAN check bugfix) * new function `flagOverlappingHz()` for identifying horizons with perfect overlap * `fillHzGaps()`, `dice()`, `slab()`, and several other functions now safely handle horizons with perfect overlap (#296) diff --git a/R/col2Munsell.R b/R/col2Munsell.R new file mode 100644 index 000000000..30e4c6124 --- /dev/null +++ b/R/col2Munsell.R @@ -0,0 +1,215 @@ + + +#' @title Convert colors into Munsell Notation +#' +#' @description Lookup the `n` closest Munsell chips from the `munsell` lookup table from various color notations. This function replaces `rgb2munsell()`. +#' +#' @param col character vector of colors, `data.frame` or `matrix` of color coordinates in sRGB or CIELAB color space +#' @param space character, one of `sRGB` or `CIELAB`, defines the input color system +#' @param nClosest integer, number of closest Munsell colors to return (valid range is 1-20) +#' + +#' +#' @note This function is fully vectorized and will pad output with NA-records when NA are present in \code{color}. +#' +#' @author D.E. Beaudette +#' +#' @references +#' \url{http://ncss-tech.github.io/AQP/} +#' \url{http://www.brucelindbloom.com/index.html?ColorCalcHelp.html} +#' \url{https://www.munsellcolourscienceforpainters.com/MunsellAndKubelkaMunkToolbox/MunsellAndKubelkaMunkToolbox.html} +#' http://www.cis.rit.edu/mcsl/online/munsell.php + + +#' @return an (NA-padded) `data.frame` containing `hue`, `value`, `chroma`, and CIE delta-E 2000 color contrast metric between source and nearest matching color(s). +#' +#' @export +#' +#' @examples +#' +#' # vector of named R colors +#' col2Munsell(c('red', 'green', 'blue')) +#' +#' # sRGB matrix in the range of 0-255 +#' col2Munsell(cbind(255, 0, 0)) +#' +#' # sRGB matrix in the range of 0-1 +#' col2Munsell(cbind(1, 0, 0)) +#' +#' # 10YR 5/6 in CIELAB +#' col2Munsell( +#' cbind(51.4337, 9.917916, 38.6889), +#' space = 'CIELAB' +#' ) +#' +#' # 2.5YR 6/8 in hex notation +#' col2Munsell("#D18158FF") +#' +#' # 7.5YR 8/1 in sRGB {0, 1} +#' col2Munsell( +#' cbind(0.8240707, 0.7856834, 0.7541048) +#' ) +#' +#' # 7.5YR 8/1 in sRGB {0, 255} +#' col2Munsell( +#' cbind(0.8240707, 0.7856834, 0.7541048) * 255 +#' ) +#' +#' # multple colors in CIELAB +#' col2Munsell( +#' parseMunsell(c('10BG 6/6', '2.5YR 4/6'), returnLAB = TRUE), +#' space = 'CIELAB' +#' ) +#' +#' # data.frame input +#' col2Munsell( +#' data.frame(r = 1, g = 0, b = 0), +#' space = 'sRGB' +#' ) +#' +#' # keep examples from using more than 2 cores +#' data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2)) +#' +#' # Munsell notation to sRGB triplets {0, 1} +#' color <- munsell2rgb( +#' the_hue = c('10YR', '2.5YR', '5YR'), +#' the_value = c(3, 5, 2.5), +#' the_chroma = c(5, 6, 2), +#' return_triplets = TRUE +#' ) +#' +#' # result is a data.frame of sRGB {0, 1} +#' color +#' +#' # back-transform sRGB -> closest Munsell color +#' # sigma is the dE00 color contrast metric +#' col2Munsell(color, space = 'sRGB') +#' +col2Munsell <- function(col, space = c('sRGB', 'CIELAB'), nClosest = 1) { + + # reasonable constraints on n-closest chips + if(nClosest < 1) { + message('setting `nClosest to 1`') + nClosest <- 1 + } + + if(nClosest > 20) { + message('setting `nClosest to 20`') + nClosest <- 20 + } + + # sacrafice to CRAN gods + munsell <- NULL + + # note: this is incompatible with LazyData: true + # load look-up table from our package + # This should be more foolproof than data(munsell) c/o PR + load(system.file("data/munsell.rda", package="aqp")[1]) + + # col can be either character vector or numeric matrix + + # convert numeric vector of colors -> sRGB -> CIELAB + if(inherits(col, 'character')) { + + # sRGB matrix [0,1] + col <- t(col2rgb(col) / 255) + + # sRGB [0,1] -> CIELAB + col <- convertColor(col, from = 'sRGB', to = 'Lab', from.ref.white = 'D65', to.ref.white = 'D65') + + } else { + + # convert to matrix as needed + if(inherits(col, 'data.frame')) { + col <- as.matrix(col) + } + + # interpret space argument + space <- match.arg(space) + + # sRGB input + # convert to CIELAB + if(space == 'sRGB') { + + # detect 0-1 vs. 0-255 sRGB range + if(range(col)[2] > 2) { + col <- col / 255 + } + + # sRGB [0,1] -> CIELAB + col <- convertColor(col, from = 'sRGB', to = 'Lab', from.ref.white = 'D65', to.ref.white = 'D65') + } + + # CIELAB input + # nothing left to do + # if(space == 'CIELAB') { + # message('no conversion required') + # } + } + + # vectorize via for-loop + n <- nrow(col) + res <- vector(length = n, mode = 'list') + + # accounting for the possibility of NA + # result should be an empty record + not.na.idx <- which(apply(col, 1, function(i) ! any(is.na(i)))) + + + # iterate over colors + for(i in not.na.idx) { + # convert current color to matrix, this will allow matrix and DF as input + this.color <- as.matrix(col[i, , drop = FALSE]) + dimnames(this.color)[[2]] <- c('L', 'A', 'B') + + # CIE dE00 + # fully-vectorized + sigma <- farver::compare_colour( + from = this.color, + to = munsell[, 7:9], + from_space = 'lab', + method = 'CIE2000', + white_from = 'D65' + ) + + # return the closest n-matches + idx <- order(sigma)[1:nClosest] + + # pack results, sorted by closest results + res[[i]] <- data.frame( + munsell[idx, 1:3], + sigma = sigma[idx] + ) + + } + + # pad records with NA in the sRGB input + # https://github.com/ncss-tech/aqp/issues/160 + na.idx <- which(sapply(res, is.null)) + + if(length(na.idx) > 0) { + for(i in na.idx){ + res[[i]] <- data.frame( + hue = NA, + value = NA, + chroma = NA, + sigma = NA, + stringsAsFactors = FALSE + ) + } + } + + + # convert to DF + res <- do.call('rbind', res) + row.names(res) <- NULL + + # save sigma units + attr(res, which = 'sigma') <- 'CIE delta-E 2000' + + return(res) + +} + + + diff --git a/R/munsell2rgb.R b/R/munsell2rgb.R index 679f776ac..d44fe4a76 100644 --- a/R/munsell2rgb.R +++ b/R/munsell2rgb.R @@ -1,7 +1,8 @@ #' @title sRGB to Munsell Color Conversion #' -#' @description Convert sRGB color coordinates to the closest `n` Munsell chips in the \code{munsell} lookup table. +#' @description Convert sRGB color coordinates to the closest `n` Munsell chips in the `munsell` lookup table. This function will be replaced by `col2Munsell()` in future versions of aqp. +#' #' #' @param color a `data.frame` or `matrix` object containing sRGB coordinates in the range of (0,1) #' @@ -19,7 +20,7 @@ #' \url{https://www.munsellcolourscienceforpainters.com/MunsellAndKubelkaMunkToolbox/MunsellAndKubelkaMunkToolbox.html} #' http://www.cis.rit.edu/mcsl/online/munsell.php #' -#' @return an (NA-padded) \code{data.frame} containing `hue`, `value`, `chroma`, and distance (dE00 when \code{colorSpace = 'CIE2000'}, Euclidean distance otherwise) to nearest matching color. +#' @return an (NA-padded) `data.frame` containing `hue`, `value`, `chroma`, and distance (dE00 when \code{colorSpace = 'CIE2000'}, Euclidean distance otherwise) to nearest matching color. #' #' @export #' @@ -44,10 +45,10 @@ #' rgb2munsell(color) #' rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nClosest = 1) { - + # argument check colorSpace <- match.arg(colorSpace) - + # reasonable constraints on n-closest chips if(nClosest < 1) { message('setting `nClosest to 1`') @@ -65,16 +66,16 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses # vectorize via for-loop n <- nrow(color) res <- vector(length = n, mode='list') - + # This is a hack to avoid munsell2rgb: "no visible binding for global variable munsell" at package R CMD check munsell <- NULL - + # note: this is incompatible with LazyData: true # load look-up table from our package # This should be more foolproof than data(munsell) c/o PR load(system.file("data/munsell.rda", package="aqp")[1]) - - + + # CIE2000 requires farver >= 2.0.3 if(colorSpace == 'CIE2000') { if( !requireNamespace('farver', quietly = TRUE) || packageVersion("farver") < '2.0.3' ) { @@ -82,17 +83,17 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses colorSpace <- 'LAB'; } } - + # accounting for the possibility of NA # result should be an empty record not.na.idx <- which(apply(color, 1, function(i) ! any(is.na(i)))) - - + + # iterate over colors for(i in not.na.idx) { # convert current color to matrix, this will allow matrix and DF as input this.color <- as.matrix(color[i, , drop=FALSE]) - + # TODO: there isn't any reason to use sRGB other than demonstration if(colorSpace == 'sRGB') { # euclidean distance (in sRGB space) is our metric for closest-color @@ -104,7 +105,7 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses # return the closest n-matches idx <- order(sigma)[1:nClosest] } - + if(colorSpace == 'LAB') { # euclidean distance (in CIELAB space) is our metric for closest-color # convert sRGB to LAB @@ -117,29 +118,42 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses # return the closest n-matches idx <- order(sigma)[1:nClosest] } - + # most accurate / efficient method as of farver >= 2.0.3 if(colorSpace == 'CIE2000') { # CIE dE00 # convert sRGB to LAB - this.color.lab <- convertColor(this.color, from='sRGB', to='Lab', from.ref.white='D65', to.ref.white = 'D65') + this.color.lab <- convertColor( + color = this.color, + from = 'sRGB', + to = 'Lab', + from.ref.white = 'D65', + to.ref.white = 'D65' + ) + dimnames(this.color.lab)[[2]] <- c('L', 'A', 'B') - + # fully-vectorized - sigma <- farver::compare_colour(this.color.lab, munsell[, 7:9], from_space='lab', method = 'CIE2000', white_from = 'D65') - + sigma <- farver::compare_colour( + from = this.color.lab, + to = munsell[, 7:9], + from_space = 'lab', + method = 'CIE2000', + white_from = 'D65' + ) + # return the closest n-matches idx <- order(sigma)[1:nClosest] } - + # pack results, sorted by closest results res[[i]] <- data.frame( munsell[idx, 1:3], sigma = sigma[idx] ) - + } - + # pad records with NA in the sRGB input # https://github.com/ncss-tech/aqp/issues/160 na.idx <- which(sapply(res, is.null)) @@ -148,15 +162,15 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses res[[i]] <- data.frame(hue=NA, value=NA, chroma=NA, sigma=NA, stringsAsFactors=FALSE) } } - - + + # convert to DF res <- do.call('rbind', res) row.names(res) <- as.character(1:nrow(res)) - + # save sigma units attr(res, which = 'sigma') <- switch(colorSpace, 'sRGB' = 'distance in sRGB', 'LAB' = 'distance in CIELAB', 'CIE2000' = 'dE00') - + return(res) } @@ -277,16 +291,16 @@ rgb2munsell <- function(color, colorSpace = c('CIE2000', 'LAB', 'sRGB'), nCloses munsell2rgb <- function(the_hue, the_value, the_chroma, alpha = 1, maxColorValue = 1, return_triplets = FALSE, returnLAB = FALSE) { - ## important: change the default behavior of data.frame and melt + ## important: change the default behavior of data.frame and melt opt.original <- options(stringsAsFactors = FALSE) - + # check for missing data - if(missing(the_hue) | missing(the_chroma) | missing(the_value)) - stop('Must supply a valid Munsell color.') - - # check to make sure that each vector is the same length - if(length(unique( c(length(the_hue), length(the_value), length(the_chroma)))) != 1) - stop('All inputs must be vectors of equal length.') + if(missing(the_hue) | missing(the_chroma) | missing(the_value)) + stop('Must supply a valid Munsell color.') + + # check to make sure that each vector is the same length + if(length(unique( c(length(the_hue), length(the_value), length(the_chroma)))) != 1) + stop('All inputs must be vectors of equal length.') # in case of factors, why would anyone do this? if(is.factor(the_hue)) { @@ -312,7 +326,7 @@ munsell2rgb <- function(the_hue, the_value, the_chroma, alpha = 1, maxColorValue # This should be moreover more foolproof than data(munsell) c/o PR munsell <- NULL load(system.file("data/munsell.rda", package="aqp")[1]) - + ## 2016-03-07: "fix" neutral hues ## they will typically be missing chroma or have some arbitrary number ## set it to 0 for correct matching @@ -366,51 +380,51 @@ munsell2rgb <- function(the_hue, the_value, the_chroma, alpha = 1, maxColorValue # reset options: options(opt.original) - + # syntax is now a little muddled, test for sRGB and LAB if(return_triplets & returnLAB) return(res[, c('r','g','b', 'L', 'A', 'B')]) - - # if the user wants the raw sRGB triplets, give those back - if(return_triplets) - return(res[, c('r','g','b')]) - + + # if the user wants the raw sRGB triplets, give those back + if(return_triplets) + return(res[, c('r','g','b')]) + if(returnLAB) return(res[, c('L','A','B')]) - - # keep track of NA values - rgb.na <- which(is.na(res$r)) - - # truncate alpha at maxColorValue, otherwise error - if(alpha > maxColorValue) { - alpha <- maxColorValue - } - - # convert to R color - res$soil_color <- NA # init an empy column - + + # keep track of NA values + rgb.na <- which(is.na(res$r)) + + # truncate alpha at maxColorValue, otherwise error + if(alpha > maxColorValue) { + alpha <- maxColorValue + } + + # convert to R color + res$soil_color <- NA # init an empy column + # account for missing values if present: we have to do this because rgb() doesn't like NA - if(length(rgb.na > 0)) { - res$soil_color[-rgb.na] <- with( - res[-rgb.na,], - rgb(red = r, green = g, blue = b, alpha = alpha, maxColorValue = maxColorValue) - ) - } else { - # no missing values - res$soil_color <- with( - res, - rgb(red = r, green = g, blue = b, alpha = alpha, maxColorValue = maxColorValue) - ) - } - - + if(length(rgb.na > 0)) { + res$soil_color[-rgb.na] <- with( + res[-rgb.na,], + rgb(red = r, green = g, blue = b, alpha = alpha, maxColorValue = maxColorValue) + ) + } else { + # no missing values + res$soil_color <- with( + res, + rgb(red = r, green = g, blue = b, alpha = alpha, maxColorValue = maxColorValue) + ) + } + + # default behavior, vector of colors is returned - return(res$soil_color) + return(res$soil_color) } # if (!isGeneric("munsell2spc")) - setGeneric("munsell2spc", function(object, ...) standardGeneric("munsell2spc")) +setGeneric("munsell2spc", function(object, ...) standardGeneric("munsell2spc")) #' @title Merge Munsell Hue, Value, Chroma converted to sRGB & CIELAB into a SoilProfileCollection #' @@ -455,70 +469,70 @@ setMethod("munsell2spc", signature(object = "SoilProfileCollection"), hue = "hue", value = "value", chroma = "chroma", .data = NULL, as.spc = TRUE) { - - # if .data vector/column/data.frame not specified - if (is.null(.data)) { - - # need hue, value, chroma as existing columns in horizon data - if (!all(c(hue, value, chroma) %in% horizonNames(object))) { - stop("arguments `hue` [character], `value` [numeric] and `chroma` [numeric] must specify column names in the horizon data", - call. = FALSE) - } else { - h <- horizons(object) - } - - } else { - # .data might be a data.frame, containing hue, value, chroma, or a character vector with munsell notation - # (e.g from parseMunsell(..., convertColors=FALSE)) - if (inherits(.data, 'data.frame')) { - - if (ncol(.data) == 1 && is.character(.data[[1]])) { - - h <- parseMunsell(.data[[1]], convertColors = FALSE) - - } else if (!all(c(hue, value, chroma) %in% colnames(.data))) { - stop("arguments `hue` [character], `value` [numeric] and `chroma` [numeric] must specify column names in `.data`", - call. = FALSE) - } else { - h <- .data - } - - } else { - - if (length(.data) == 1 && .data %in% horizonNames(object)) { - .data <- object[[.data]] - - # otherwise need munsell character columnname, or a vector with equal length to horizons - } else if (length(.data) != nrow(object)) { - stop("argument `.data` [character or data.frame], must specify either a character vector of equal length to number of horizons , a column name in the horizon data (containing Munsell notation) or a data.frame with three columns (names specified in `hue`, `value`, `chroma`)", - call. = FALSE) - } - - h <- parseMunsell(.data, convertColors = FALSE) - } - } - - # makes a data.frame - # return sRGB + CIELAB at the same time, note that hex notation of color is not returned - drgb <- munsell2rgb(h[[hue]], h[[value]], h[[chroma]], return_triplets = TRUE, returnLAB = TRUE) - colnames(drgb)[1:3] <- paste0("rgb_", c("R","G","B")) - colnames(drgb)[4:6] <- paste0("lab_", c("L","A","B")) - - # ID management - idn <- idname(object) - hidn <- hzidname(object) - h <- horizons(object) - - # munsell2rgb does not return ID names (not inherently aware of the SPC) - idcol <- data.frame(h[[idn]], h[[hidn]]) - colnames(idcol) <- c(idn, hidn) - - if (as.spc) { - # horizons<- will ensure merge.data.table triggers if @horizons is data.table - horizons(object) <- cbind(idcol, drgb) - - return(object) - } else { - return(.as.data.frame.aqp(cbind(idcol, drgb), aqp_df_class(object))) - } -}) + + # if .data vector/column/data.frame not specified + if (is.null(.data)) { + + # need hue, value, chroma as existing columns in horizon data + if (!all(c(hue, value, chroma) %in% horizonNames(object))) { + stop("arguments `hue` [character], `value` [numeric] and `chroma` [numeric] must specify column names in the horizon data", + call. = FALSE) + } else { + h <- horizons(object) + } + + } else { + # .data might be a data.frame, containing hue, value, chroma, or a character vector with munsell notation + # (e.g from parseMunsell(..., convertColors=FALSE)) + if (inherits(.data, 'data.frame')) { + + if (ncol(.data) == 1 && is.character(.data[[1]])) { + + h <- parseMunsell(.data[[1]], convertColors = FALSE) + + } else if (!all(c(hue, value, chroma) %in% colnames(.data))) { + stop("arguments `hue` [character], `value` [numeric] and `chroma` [numeric] must specify column names in `.data`", + call. = FALSE) + } else { + h <- .data + } + + } else { + + if (length(.data) == 1 && .data %in% horizonNames(object)) { + .data <- object[[.data]] + + # otherwise need munsell character columnname, or a vector with equal length to horizons + } else if (length(.data) != nrow(object)) { + stop("argument `.data` [character or data.frame], must specify either a character vector of equal length to number of horizons , a column name in the horizon data (containing Munsell notation) or a data.frame with three columns (names specified in `hue`, `value`, `chroma`)", + call. = FALSE) + } + + h <- parseMunsell(.data, convertColors = FALSE) + } + } + + # makes a data.frame + # return sRGB + CIELAB at the same time, note that hex notation of color is not returned + drgb <- munsell2rgb(h[[hue]], h[[value]], h[[chroma]], return_triplets = TRUE, returnLAB = TRUE) + colnames(drgb)[1:3] <- paste0("rgb_", c("R","G","B")) + colnames(drgb)[4:6] <- paste0("lab_", c("L","A","B")) + + # ID management + idn <- idname(object) + hidn <- hzidname(object) + h <- horizons(object) + + # munsell2rgb does not return ID names (not inherently aware of the SPC) + idcol <- data.frame(h[[idn]], h[[hidn]]) + colnames(idcol) <- c(idn, hidn) + + if (as.spc) { + # horizons<- will ensure merge.data.table triggers if @horizons is data.table + horizons(object) <- cbind(idcol, drgb) + + return(object) + } else { + return(.as.data.frame.aqp(cbind(idcol, drgb), aqp_df_class(object))) + } + }) diff --git a/man/col2Munsell.Rd b/man/col2Munsell.Rd new file mode 100644 index 000000000..182b6e5ba --- /dev/null +++ b/man/col2Munsell.Rd @@ -0,0 +1,94 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/col2Munsell.R +\name{col2Munsell} +\alias{col2Munsell} +\title{Convert colors into Munsell Notation} +\usage{ +col2Munsell(col, space = c("sRGB", "CIELAB"), nClosest = 1) +} +\arguments{ +\item{col}{character vector of colors, \code{data.frame} or \code{matrix} of color coordinates in sRGB or CIELAB color space} + +\item{space}{character, one of \code{sRGB} or \code{CIELAB}, defines the input color system} + +\item{nClosest}{integer, number of closest Munsell colors to return (valid range is 1-20)} +} +\value{ +an (NA-padded) \code{data.frame} containing \code{hue}, \code{value}, \code{chroma}, and CIE delta-E 2000 color contrast metric between source and nearest matching color(s). +} +\description{ +Lookup the \code{n} closest Munsell chips from the \code{munsell} lookup table from various color notations. This function replaces \code{rgb2munsell()}. +} +\note{ +This function is fully vectorized and will pad output with NA-records when NA are present in \code{color}. +} +\examples{ + +# vector of named R colors +col2Munsell(c('red', 'green', 'blue')) + +# sRGB matrix in the range of 0-255 +col2Munsell(cbind(255, 0, 0)) + +# sRGB matrix in the range of 0-1 +col2Munsell(cbind(1, 0, 0)) + +# 10YR 5/6 in CIELAB +col2Munsell( + cbind(51.4337, 9.917916, 38.6889), + space = 'CIELAB' +) + +# 2.5YR 6/8 in hex notation +col2Munsell("#D18158FF") + +# 7.5YR 8/1 in sRGB {0, 1} +col2Munsell( + cbind(0.8240707, 0.7856834, 0.7541048) +) + +# 7.5YR 8/1 in sRGB {0, 255} +col2Munsell( + cbind(0.8240707, 0.7856834, 0.7541048) * 255 +) + +# multple colors in CIELAB +col2Munsell( + parseMunsell(c('10BG 6/6', '2.5YR 4/6'), returnLAB = TRUE), + space = 'CIELAB' +) + +# data.frame input +col2Munsell( + data.frame(r = 1, g = 0, b = 0), + space = 'sRGB' +) + +# keep examples from using more than 2 cores +data.table::setDTthreads(Sys.getenv("OMP_THREAD_LIMIT", unset = 2)) + +# Munsell notation to sRGB triplets {0, 1} +color <- munsell2rgb( + the_hue = c('10YR', '2.5YR', '5YR'), + the_value = c(3, 5, 2.5), + the_chroma = c(5, 6, 2), + return_triplets = TRUE +) + +# result is a data.frame of sRGB {0, 1} +color + +# back-transform sRGB -> closest Munsell color +# sigma is the dE00 color contrast metric +col2Munsell(color, space = 'sRGB') + +} +\references{ +\url{http://ncss-tech.github.io/AQP/} +\url{http://www.brucelindbloom.com/index.html?ColorCalcHelp.html} +\url{https://www.munsellcolourscienceforpainters.com/MunsellAndKubelkaMunkToolbox/MunsellAndKubelkaMunkToolbox.html} +http://www.cis.rit.edu/mcsl/online/munsell.php +} +\author{ +D.E. Beaudette +} diff --git a/man/rgb2munsell.Rd b/man/rgb2munsell.Rd index 5c7cb5706..bef7a3292 100644 --- a/man/rgb2munsell.Rd +++ b/man/rgb2munsell.Rd @@ -17,7 +17,7 @@ rgb2munsell(color, colorSpace = c("CIE2000", "LAB", "sRGB"), nClosest = 1) an (NA-padded) \code{data.frame} containing \code{hue}, \code{value}, \code{chroma}, and distance (dE00 when \code{colorSpace = 'CIE2000'}, Euclidean distance otherwise) to nearest matching color. } \description{ -Convert sRGB color coordinates to the closest \code{n} Munsell chips in the \code{munsell} lookup table. +Convert sRGB color coordinates to the closest \code{n} Munsell chips in the \code{munsell} lookup table. This function will be replaced by \code{col2Munsell()} in future versions of aqp. } \note{ This function is fully vectorized and will pad output with NA-records when NA are present in \code{color}. diff --git a/tests/testthat/test-color-conversion.R b/tests/testthat/test-color-conversion.R index a564e7ec4..ed0338347 100644 --- a/tests/testthat/test-color-conversion.R +++ b/tests/testthat/test-color-conversion.R @@ -11,15 +11,34 @@ m.rgb <- munsell2rgb(x.p$hue, x.p$value, x.p$chroma, return_triplets = TRUE) # sRGB --> Munsell x.back <- rgb2munsell(color = m.rgb, colorSpace = 'LAB', nClosest = 1) # using truncated sRGB values -x.back.trunc <- rgb2munsell(data.frame(r=0.36, g=0.26, b=0.13)) +x.back.trunc <- rgb2munsell(data.frame(r = 0.36, g = 0.26, b = 0.13)) # neutral colors map to shades of gray -x.neutral <- parseMunsell('N 2/', return_triplets=TRUE) +x.neutral <- parseMunsell('N 2/', return_triplets = TRUE) ## tests -test_that("parseMunsell()", { +test_that("col2Munsell works as expected", { + + # 10YR 3/4 as hex notation sRGB + .res <- col2Munsell(m) + expect_equal(.res$hue, '10YR') + expect_equal(.res$value, 3L) + expect_equal(.res$chroma, 4L) + + + # N 2/ as sRGB [0,1] + .res <- col2Munsell(x.neutral, space = 'sRGB') + expect_equal(.res$hue, 'N') + expect_equal(.res$value, 2L) + expect_equal(.res$chroma, 0L) + +}) + + +test_that("parseMunsell()", { + # parsing bogus notation generates NA # will also generate a warning from munsell2rgb() expect_equal(suppressWarnings(parseMunsell('10YZ 4/5')), NA_character_) @@ -27,13 +46,13 @@ test_that("parseMunsell()", { expect_equal(suppressWarnings(parseMunsell('10YR ')), NA_character_) expect_equal(suppressWarnings(parseMunsell('10YR 4/')), NA_character_) expect_equal(suppressWarnings(parseMunsell('G1 6/N')), NA_character_) - + # parsing bogus notation without conversion bogus <- parseMunsell('G1 3/X', convertColors = FALSE) expect_equal(bogus$hue, NA_character_) expect_equal(bogus$value, NA_real_) expect_equal(bogus$chroma, NA_real_) - + # test NA some.NA <- parseMunsell(c(NA, '10YR 3/3')) expect_true(inherits(some.NA, 'character')) @@ -45,7 +64,7 @@ test_that("parseMunsell()", { # splitting of text into columns within data.frame expect_identical(x.p, data.frame(hue = "10YR", value = 3, chroma = 4, stringsAsFactors = FALSE)) - + # Test not using spaces expect_equal(suppressWarnings(parseMunsell('2.5YR 3/4')), suppressWarnings(parseMunsell('2.5YR3/4'))) @@ -57,14 +76,14 @@ test_that("parseMunsell()", { # addresses #66 (https://github.com/ncss-tech/aqp/issues/66) test_that("Munsell hue parsing", { - + # normal operation res <- .parseMunsellHue('10YR') expect_true(inherits(res, 'data.frame')) expect_equal(res$hue.numeric, 10L) expect_equal(res$hue.character, 'YR') expect_equal(nrow(res), 1) - + # white space is trimmed res <- .parseMunsellHue('10 YR') expect_true(inherits(res, 'data.frame')) @@ -89,7 +108,7 @@ test_that("Munsell hue parsing", { test_that("non-integer value and chroma are selectively rounded", { - + # rounding of value, throws warning expect_warning(res <- parseMunsell('10YR 3.3/4'), regexp = 'non-standard notation') # this will not throw a warning @@ -100,7 +119,7 @@ test_that("non-integer value and chroma are selectively rounded", { suppressWarnings(parseMunsell('10YR 3.3/4')), parseMunsell('10YR 3/4') ) - + # rounding of chroma, throws warning expect_warning(res <- parseMunsell('10YR 3/4.6'), regexp = 'non-standard notation') # this will not throw a warning @@ -111,7 +130,7 @@ test_that("non-integer value and chroma are selectively rounded", { suppressWarnings(parseMunsell('10YR 3/4.6')), parseMunsell('10YR 3/5') ) - + # no rounding of 2.5 values res <- parseMunsell('10YR 2.5/2') res.test <- rgb2munsell(t(col2rgb('#493A2BFF')/255)) @@ -121,27 +140,27 @@ test_that("non-integer value and chroma are selectively rounded", { test_that("Munsell <--> sRGB and back again", { - + # sRGB in hex notation expect_equal(m, '#5C4222FF') expect_equal(parseMunsell(x), m) - + # sRGB triplets expect_equal(m.rgb$r, 0.3618738, tolerance=0.0001) expect_equal(m.rgb$g, 0.2598939, tolerance=0.0001) expect_equal(m.rgb$b, 0.1337521, tolerance=0.0001) - + # neutral colors expect_equal(x.neutral$r, 0.03278, tolerance=0.001) expect_equal(x.neutral$g, 0.03305, tolerance=0.001) expect_equal(x.neutral$b, 0.03305, tolerance=0.001) - + # sRGB --> Munsell expect_equal(x.back$hue, '10YR') expect_equal(x.back$value, 3) expect_equal(x.back$chroma, 4) expect_equal(x.back$sigma, 0) - + expect_equal(x.back.trunc$hue, '10YR') expect_equal(x.back.trunc$value, 3) expect_equal(x.back.trunc$chroma, 4) @@ -149,7 +168,7 @@ test_that("Munsell <--> sRGB and back again", { test_that("missing data", { - + # data with missing sRGB coordinates color <- rbind( cbind(NA, NA, NA), @@ -157,16 +176,16 @@ test_that("missing data", { cbind(1, 1, 1), cbind(NA, NA, NA) ) - + # conversion should work without error res <- rgb2munsell(color) - + # same number of rows in / out expect_true(nrow(res) == nrow(color)) - + # row order preserved expect_true(is.na(res$hue[1]) & is.na(res$hue[4])) - + }) test_that("neutral hues", { @@ -201,7 +220,7 @@ test_that("neutral hues", { test_that("closest Munsell chip based on sRGB coordinates", { - + # closest chip in aqp LUT expect_equal(getClosestMunsellChip('10YR 3.3/5', convertColors = FALSE), '10YR 3/5') expect_equal(getClosestMunsellChip('9YR 3.8/3', convertColors = FALSE), '10YR 4/3') @@ -213,23 +232,23 @@ test_that("closest Munsell chip based on sRGB coordinates", { # https://github.com/ncss-tech/aqp/issues/69 test_that("Munsell --> LAB + sRGB coordinates", { - + # sRGB test.1 <- parseMunsell("10YR 3/5", return_triplets=TRUE) expect_equal(names(test.1), c('r', 'g', 'b')) - - + + # sRGB and LAB test.2 <- parseMunsell("10YR 3/5", return_triplets=TRUE, returnLAB=TRUE) expect_equal(names(test.2), c('r', 'g', 'b', 'L', 'A', 'B')) - + # LAB test.3 <- parseMunsell("10YR 3/5", return_triplets=FALSE, returnLAB=TRUE) expect_equal(names(test.3), c('L', 'A', 'B')) - + # test the LAB ---> sRGB is close test.4 <- grDevices::convertColor(test.3, from = 'Lab', to='sRGB') - + # sRGB (r) expect_equal(test.1[, 1], test.4[, 1], tolerance=0.1) # sRGB (g) @@ -240,10 +259,10 @@ test_that("Munsell --> LAB + sRGB coordinates", { test_that("similar colors result in same, closest chip", { - + cols <- t(col2rgb(c('#5F5335', '#5F5236'))) / 255 res <- rgb2munsell(cols) - + expect_equal(res$hue[1], res$hue[2]) expect_equal(res$value[1], res$value[2]) expect_equal(res$chroma[1], res$chroma[2]) @@ -251,36 +270,36 @@ test_that("similar colors result in same, closest chip", { test_that("munsell2spc wrapper method works as expected", { - + data(sp3) depths(sp3) <- id ~ top + bottom - + # inspect input data # horizons(sp3)[,c("hue","value","chroma")] - + # do color conversions to sRGB and LAB, join into horizon data expect_silent( {sp3 <- munsell2spc(sp3)}) expect_true(inherits(sp3, 'SoilProfileCollection')) - + # # plot rgb "R" coordinate by horizon # plot(sp3, color = "rgb_R") # # # plot lab "A" coordinate by horizon # plot(sp3, color = "lab_A") - + # test returning profile+horizon ID data.frame with results expect_silent( { dftest <- munsell2spc(sp3, as.spc = FALSE) } ) expect_true(inherits(dftest, 'data.frame')) - + # foo is not a column in horizons() expect_error( { err1 <- munsell2spc(sp3, hue = "foo") } ) - + # chip is not a column in horizons expect_error( { d1 <- horizons(munsell2spc(sp3, .data = "chip")) } ) - + # create chip as a combination of hue value/chroma sp3$chip <- with(horizons(sp3), sprintf("%s %s/%s", hue, value, chroma)) - + # calculate from: column name, vector, data.frame expect_silent( { d1 <- horizons(munsell2spc(sp3, .data = "chip")) } ) expect_silent( { d2 <- horizons(munsell2spc(sp3, .data = sp3$chip)) } )