code and commentary cleanup. caught some interpolated spectra < 0

ncss-tech · Oct 3, 2024 · 32f7bef · 32f7bef
1 parent 2e75115
commit 32f7bef
Show file tree

Hide file tree

Showing 8 changed files with 199 additions and 153 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -15,6 +15,6 @@ URL: https://github.com/ncss-tech/aqp, https://ncss-tech.github.io/AQP/
 BugReports: https://github.com/ncss-tech/aqp/issues
 Language: en-US
 Encoding: UTF-8
-RoxygenNote: 7.3.1
+RoxygenNote: 7.3.2
 Roxygen: list(markdown = TRUE)
 VignetteBuilder: knitr
diff --git a/data/munsell.spectra.rda b/data/munsell.spectra.rda
diff --git a/data/munsell.spectra.wide.rda b/data/munsell.spectra.wide.rda
diff --git a/misc/utils/Munsell/interpolate-spectra.R b/misc/utils/Munsell/interpolate-spectra.R
@@ -14,10 +14,18 @@ library(tactile)
 library(purrr)
 library(reshape2)
 
+source('local-functions.R')
+
 # load simplified spectra
 m.rel <- readRDS('simplified-Munsell-spectra.rds')
 
-## investigate slices -- can interpolate reflectance vs chroma (by wavelengths) for odd chroma
+# review original range
+range(m.rel$reflectance)
+min.reflectance <- min(m.rel$reflectance)
+
+hist(m.rel$reflectance, breaks = 50)
+
+## interpolate spectra for odd Munsell chroma
 
 idx <- which(m.rel$hue %in% c('7.5YR') & m.rel$value == 3)
 s <- m.rel[idx, ]
@@ -32,75 +40,6 @@ xyplot(reflectance ~ chroma | factor(wavelength), data=s,
 # split by hue/value/wavelength
 m <- split(m.rel, list(m.rel$hue, m.rel$value, m.rel$wavelength))
 
-## TODO: clamp to original range of Munsell chroma
-
-# interpolation of odd chroma
-interpolateOddChromaSpectra <- function(i) {
-
-  # 0-row input
-  if(nrow(i) < 1)
-    return(NULL)
-
-  # chroma stats
-  u.chroma <- unique(i$chroma)
-  r.chroma <- range(u.chroma)
-  n.chroma <- length(u.chroma)
-
-  # reflectance stats
-  r.reflectance <- range(i$reflectance)
-
-  # sequence of candidate chroma
-  s <- seq(from = r.chroma[1], to = r.chroma[2], by = 1)
-  s.chroma <- setdiff(s, u.chroma)
-
-  # short circuit: single chroma, interpolation impossible
-  if(n.chroma < 2)
-    return(NULL)
-
-  # short circuit: 0 candidates for interpolation
-  if(length(s.chroma) < 1)
-    return(NULL)
-
-
-  # setup interpolation function: natural splines
-  # fit is exact at training points
-  af <- splinefun(i$chroma, i$reflectance, method = 'natural')
-
-  # check: fit should be exact at points
-  if(sum(af(i$chroma) - i$reflectance) > 0.001){
-    message('spline not fitting at training data!')
-  }
-
-  # interpolate candidate chroma
-  s.reflectance <- af(s.chroma)
-
-  # # check for over / undershoots
-  # if( ! (s.reflectance < r.reflectance[1] | s.reflectance > r.reflectance[2])){
-  #   message('spline exceeds original range!')
-  # }
-
-  # re-assemble into original format
-  res <- data.frame(
-    munsell = sprintf("%s %s/%s", i$hue[1], i$value[1], s.chroma),
-    hue = i$hue[1],
-    value = i$value[1],
-    chroma = s.chroma,
-    wavelength = i$wavelength[1],
-    reflectance = s.reflectance,
-    stringsAsFactors = FALSE
-  )
-
-
-  # debugging: graphical check
-  # OK
-  # plot(reflectance ~ chroma, data = i )
-  # lines(seq(r.chroma[1], r.chroma[2], by = 0.1), af(seq(r.chroma[1], r.chroma[2], by = 0.1)), col = 'red')
-  # points(s.chroma, s.reflectance, pch = 15)
-
-  return(res)
-
-}
-
 # do interpolation
 mm <- map(m, .f = interpolateOddChromaSpectra, .progress = TRUE)
 
@@ -113,7 +52,7 @@ m.final <- m.final[order(m.final$hue, m.final$value, m.final$chroma), ]
 
 
 # graphical check
-idx <- which(m.final$hue %in% c('7.5YR') & m.final$value == 5)
+idx <- which(m.final$hue %in% c('7.5YR') & m.final$value == 3)
 s <- m.final[idx, ]
 
 xyplot(reflectance ~ chroma | factor(wavelength), data=s, 
@@ -134,7 +73,8 @@ xyplot(reflectance ~ chroma | factor(wavelength), data=s,
 
 
 # check for reflectance <= 0
-m.final[m.final$reflectance <= 0, ]
+# 6 rows
+nrow(m.final[m.final$reflectance <= 0, ])
 
 # hmm
 idx <- which(m.final$hue %in% c('2.5R') & m.final$value == 2)
@@ -151,8 +91,18 @@ xyplot(reflectance ~ chroma | factor(wavelength), groups = reflectance <= 0, dat
 idx <- which(m.final$reflectance <= 0)
 m.final[idx, ]
 
-# replace with minimum reflectance, ignoring these values
-m.final$reflectance[idx] <- min(m.final$reflectance[-idx])
+# replace with original minimum reflectance, ignoring these values
+m.final$reflectance[idx] <- min.reflectance
+
+xyplot(reflectance ~ chroma | factor(wavelength), groups = reflectance <= 0, 
+       data = m.final,
+       subset = hue == '2.5R' & value == 2,
+       type='b', as.table=TRUE,
+       scales = list(y = list(tick.number = 10)),
+       auto.key=list(lines=FALSE, points=TRUE, cex=1, space='top'),
+       par.settings = tactile.theme()
+)
+
 
 
 ## check: OK
@@ -177,58 +127,12 @@ xyplot(reflectance ~ wavelength, data = s,
 
 
 
-## interpolate half-chip Munsell values
-
-# # just 2/3 value
-# m.sub <- subset(m.final, subset = value %in% 2:3)
-# 
-# head(m.sub)
-# 
-# # check
-# s <- subset(m.sub, subset = hue == '2.5YR' & chroma == 3)
-# 
-# xyplot(reflectance ~ wavelength, data = s, 
-#        groups = munsell, type='b',
-#        scales = list(y = list(tick.number = 10)),
-#        auto.key=list(lines=TRUE, points=FALSE, cex=1, space='top', columns = 2),
-#        par.settings = tactile.theme()
-# )
-
+## interpolate spectra for select half-chip Munsell values
 
 # split by hue/chroma/wavelength
 m <- split(m.final, list(m.final$hue, m.final$chroma, m.final$wavelength))
 
 
-## TODO: clamp to original range of Munsell value
-
-interpolateValueSpectra <- function(i) {
-
-  # 0 or 1 row input: no interpolation possible
-  if(nrow(i) < 2)
-    return(NULL)
-
-  # setup interpolation function: natural splines
-  # fit is exact at training points
-  a.fun <- splinefun(i$value, i$reflectance, method = 'natural')
-
-  # new Munsell values
-  v.target <- c(2.5, 8.5, 9.5)
-
-  # re-assemble into original format
-  res <- data.frame(
-    munsell = sprintf("%s %s/%s", i$hue[1], v.target, i$chroma[1]),
-    hue = i$hue[1],
-    value = v.target,
-    chroma = i$chroma[1],
-    wavelength = i$wavelength[1],
-    reflectance = a.fun(v.target),
-    stringsAsFactors = FALSE
-  )
-
-  return(res)
-}
-
-
 # do interpolation
 mm <- map(m, .f = interpolateValueSpectra, .progress = TRUE)
 
@@ -240,10 +144,35 @@ mm <- do.call('rbind', mm)
 m.final <- rbind(m.final, mm)
 m.final <- m.final[order(m.final$hue, m.final$value, m.final$chroma), ]
 
-
 # check: OK
 str(m.final)
 
+
+# check for reflectance <= 0
+# 3403 rows, all very close to 0
+# most of these are very low value + low chroma | low value + high chroma
+nrow(m.final[m.final$reflectance <= 0, ])
+
+# hmm
+idx <- which(m.final$munsell == '7.5YR 2.5/14')
+s <- m.final[idx, ]
+
+xyplot(reflectance ~ chroma | factor(wavelength), groups = reflectance <= 0, data=s, 
+       type='b', as.table=TRUE,
+       scales = list(y = list(tick.number = 10)),
+       auto.key=list(lines=FALSE, points=TRUE, cex=1, space='top'),
+       par.settings = tactile.theme()
+)
+
+# probably spline undershoots
+idx <- which(m.final$reflectance <= 0)
+
+# replace with minimum reflectance, ignoring these values
+m.final$reflectance[idx] <- min.reflectance
+
+
+
+
 s <- subset(m.final, subset = hue == '10YR' & chroma == 4 & value %in% c(2, 2.5, 3, 4))
 
 xyplot(reflectance ~ wavelength, data = s, 

diff --git a/misc/utils/Munsell/local-functions.R b/misc/utils/Munsell/local-functions.R
@@ -1,4 +1,103 @@
 
+
+
+## TODO: clamp to original range of Munsell chroma
+
+# interpolation of odd chroma
+interpolateOddChromaSpectra <- function(i) {
+
+  # 0-row input
+  if(nrow(i) < 1) {
+    return(NULL)
+  }
+
+  # chroma stats
+  u.chroma <- unique(i$chroma)
+  r.chroma <- range(u.chroma)
+  n.chroma <- length(u.chroma)
+
+  # reflectance stats
+  r.reflectance <- range(i$reflectance)
+
+  # sequence of candidate chroma
+  s <- seq(from = r.chroma[1], to = r.chroma[2], by = 1)
+  s.chroma <- setdiff(s, u.chroma)
+
+  # short circuit: single chroma, interpolation impossible
+  if(n.chroma < 2)
+    return(NULL)
+
+  # short circuit: 0 candidates for interpolation
+  if(length(s.chroma) < 1)
+    return(NULL)
+
+
+  # setup interpolation function: natural splines
+  # fit is exact at training points
+  sf <- splinefun(i$chroma, i$reflectance, method = 'natural')
+
+  # check: fit should be exact at points
+  if(sum(sf(i$chroma) - i$reflectance) > 0.001){
+    message('spline not fitting at training data!')
+  }
+
+  # interpolate candidate chroma
+  s.reflectance <- sf(s.chroma)
+
+  # re-assemble into original format
+  res <- data.frame(
+    munsell = sprintf("%s %s/%s", i$hue[1], i$value[1], s.chroma),
+    hue = i$hue[1],
+    value = i$value[1],
+    chroma = s.chroma,
+    wavelength = i$wavelength[1],
+    reflectance = s.reflectance,
+    stringsAsFactors = FALSE
+  )
+
+
+  # debugging: graphical check
+  # OK
+  # plot(reflectance ~ chroma, data = i )
+  # lines(seq(r.chroma[1], r.chroma[2], by = 0.1), af(seq(r.chroma[1], r.chroma[2], by = 0.1)), col = 'red')
+  # points(s.chroma, s.reflectance, pch = 15)
+
+  return(res)
+
+}
+
+
+interpolateValueSpectra <- function(i) {
+
+  # 0 or 1 row input: no interpolation possible
+  if(nrow(i) < 2)
+    return(NULL)
+
+  # setup interpolation function: natural splines
+  # fit is exact at training points
+  a.fun <- splinefun(i$value, i$reflectance, method = 'natural')
+
+  # new Munsell values
+  v.target <- c(2.5, 8.5, 9.5)
+
+  # re-assemble into original format
+  res <- data.frame(
+    munsell = sprintf("%s %s/%s", i$hue[1], v.target, i$chroma[1]),
+    hue = i$hue[1],
+    value = v.target,
+    chroma = i$chroma[1],
+    wavelength = i$wavelength[1],
+    reflectance = a.fun(v.target),
+    stringsAsFactors = FALSE
+  )
+
+  return(res)
+}
+
+
+
+
+
 # 2022-03-29
 # interpolate odd chroma from Munsell renotation data
 # m.i: subset renotation data.frame, for a single hue/value

diff --git a/misc/utils/Munsell/main.R b/misc/utils/Munsell/main.R
@@ -1,5 +1,5 @@
 ## Code / Data related to preparation of Munsell color interpretation in aqp
-## 2024-09-26
+## 2024-10-03
 ## D.E. Beaudette, A.G. Brown
 
 # make Munsell and related LUT