Nothing
#' Aligns images using \code{\link[RNiftyReg]{niftyreg}} utilities for automated image registration
#' and extracts color pattern using watershed segmentation. This function works interactively by
#' allowing to pick a starting pixel within each pattern element from which the watershed will
#' extract the pattern. This function works best for patterns with sharp boundaries.
#'
#' @param sampleList List of RasterStack objects.
#' @param target Image imported as RasterStack used as target for registration.
#' @param resampleFactor Integer for downsampling image used by \code{\link{redRes}}.
#' @param crop Vector c(xmin, xmax, ymin, ymax) that specifies the pixel coordinates to crop the
#' original image.
#' @param useBlockPercentage Block percentage as used in \code{\link[RNiftyReg]{niftyreg}}
#' (default = 75).
#' @param removebgR Integer indicating the range RGB treshold to remove from image (e.g. 100 removes
#' pixels with average RGB > 100; default = NULL) for registration analysis. This works only to
#' remove a white background.
#' @param maskOutline When outline is specified, everything outside of the outline will be masked for
#' the color extraction (default = NULL).
#' @param cartoonID ID of the sample for which the cartoon was drawn and will be used for masking.
#' @param correct Correct image illumination using a linear model (default = FALSE).
#' @param blur Blur image for priority map extraction (default = TRUE).
#' @param sigma Size of sigma for Gaussian blurring (default = 5).
#' @param bucketfill Use a bucket fill on the background to fill holes (default = TRUE).
#' @param cleanP Integer to remove spurious areas with width smaller than cleanP (default = NULL).
#' @param splitC Integer to split selected patterns into connected components and remove ones with
#' areas smaller than splitC (default = NULL).
#' @param plotTransformed Plot transformed image (default = FALSE).
#' @param plotCorrect Plot corrected image, corrected for illumination using a linear model
#' (default = FALSE).
#' @param plotEdges Plot image gradient (default = FALSE).
#' @param plotPriority Plot priority map (default = FALSE).
#' @param plotWS Plot watershed result (default = FALSE).
#' @param plotBF Plot bucketfill (default = FALSE).
#' @param plotFinal Plot extracted patterns (default = FALSE).
#'
#' @return List of raster objects.
#'
#' @examples
#' \dontrun{
#' IDlist <- c('BC0077','BC0071','BC0050','BC0049','BC0004')
#' prepath <- system.file("extdata", package = 'patternize')
#' extension <- '.jpg'
#'
#' imageList <- makeList(IDlist, 'image', prepath, extension)
#'
#' target <- imageList[[1]]
#'
#' outline_BC0077 <- read.table(paste(system.file("extdata", package = 'patternize'),
#' '/BC0077_outline.txt', sep=''), header = FALSE)
#'
#' rasterList_regW <- patRegW(imageList, target, plotTransformed = FALSE, cartoonID = 'BC0077',
#' correct = TRUE, plotCorrect = FALSE, blur = FALSE, sigma = 2,
#' bucketfill = FALSE, cleanP = 0, splitC = 10, plotPriority = TRUE,
#' plotWS = FALSE, plotBF = FALSE, plotFinal = TRUE, removebgR = 100,
#' maskOutline = outline_BC0077)
#' }
#'
#' @export
#'
#' @import raster
#' @importFrom imager as.cimg imfill isoblur imgradient imsplit watershed bucketfill clean split_connected highlight enorm parany add
#' @importFrom purrr discard
#' @importFrom dplyr sample_n
#' @importFrom magrittr %>%
#' @importFrom graphics locator
#' @importFrom stats lm
patRegW <- function(sampleList,
target,
resampleFactor = NULL,
useBlockPercentage = 75,
crop = c(0,0,0,0),
removebgR = NULL,
maskOutline = NULL,
cartoonID = NULL,
correct = FALSE,
blur = TRUE,
sigma = 3,
bucketfill = TRUE,
cleanP = NULL,
splitC = NULL,
plotTransformed = FALSE,
plotCorrect = FALSE,
plotEdges = FALSE,
plotPriority = FALSE,
plotWS = FALSE,
plotBF = FALSE,
plotFinal = FALSE){
declare( . <- as_nonstandard())
rasterList <- list()
# Crop target image
if(!identical(crop, c(0,0,0,0))){
targetExtRaster <- crop
target <- raster::crop(target, targetExtRaster)
}
# Reduce resolution of target image
if(!is.null(resampleFactor)){
target <- redRes(target, resampleFactor)
}
# Average color channels to gray scale
targetA <- apply(raster::as.array(target), 1:2, mean)
# Remove bacground of target image
if(is.numeric(removebgR)){
targetA <- apply(targetA, 1:2, function(x) ifelse(x > removebgR, 0, x))
}
if(!is.null(maskOutline)){
if(is.null(cartoonID)){
print('Please specify cartoonID')
}
indx <- which(names(sampleList) == cartoonID)
maskOutlineNew <- maskOutline
extPicture <- raster::extent(sampleList[[indx]])
maskOutlineNew[,2] <- extPicture[4]-maskOutlineNew[,2]
}
n <- 1
todo <- 0
# Run the loop for each sample
while(n <= length(sampleList)){
if(todo != 'r' && todo != 'a'){
sStack <- sampleList[[n]]
extRaster <- raster::extent(sStack)
# Crop image
if(!identical(crop, c(0,0,0,0))){
extRaster <- crop
sStack <- crop(sStack, extRaster)
}
sourceRaster <- redRes(sStack, 1)
# Reduce resolution
if(!is.null(resampleFactor)){
sourceRaster <- redRes(sStack, resampleFactor)
}
sourceRasterK <- sourceRaster
# Average color channels to gray scale
sourceRaster <- apply(raster::as.array(sourceRaster), 1:2, mean)
# Remove bacground for registration
if(is.numeric(removebgR)){
sourceRaster <- apply(sourceRaster, 1:2, function(x) ifelse(x > removebgR, 0, x))
}
# Calculate transformation
result <- RNiftyReg::niftyreg(sourceRaster, targetA, useBlockPercentage=useBlockPercentage)
# Apply transformation
transformedMap <- RNiftyReg::applyTransform(RNiftyReg::forward(result), raster::as.array(sourceRasterK), interpolation=0)
r1 <- raster::raster(transformedMap[1:nrow(transformedMap),ncol(transformedMap):1,1])
# r1 <- focal(r1, w=matrix(1,nrow=3,ncol=3), fun=fill.na, pad = TRUE, na.rm = FALSE)
r2 <- raster::raster(transformedMap[1:nrow(transformedMap),ncol(transformedMap):1,2])
# r2 <- focal(r2, w=matrix(1,nrow=3,ncol=3), fun=fill.na, pad = TRUE, na.rm = FALSE)
r3 <- raster::raster(transformedMap[1:nrow(transformedMap),ncol(transformedMap):1,3])
# r3 <- focal(r3, w=matrix(1,nrow=3,ncol=3), fun=fill.na, pad = TRUE, na.rm = FALSE)
imageTr <-raster::stack(r1,r2,r3)
imageTr <- raster::flip(imageTr,'x')
raster::extent(imageTr) <- raster::extent(sourceRasterK)
# imageTr[imageTr == 0] <- NA
# Plot transformed raster
if(plotTransformed){
imageTr <- raster::flip(imageTr, 'y')
x <- as.array(imageTr)
cols <- rgb(x[,,1], x[,,2], x[,,3], maxColorValue=255)
uniqueCols <- unique(cols)
x2 <- match(cols, uniqueCols)
dim(x2) <- dim(x)[1:2]
raster::image(t(x2), col=uniqueCols, yaxt='n', xaxt='n', main = paste(names(sampleList)[n],'transformed', sep=' '))
imageTr <- raster::flip(imageTr, 'y')
}
# MaskOutline
if(!is.null(maskOutline)){
imageTr <- maskOutline(imageTr, maskOutline, refShape = 'target', flipOutline = 'y', crop = crop,
maskColor = 255, imageList = sampleList)
cropEx <- c(min(maskOutlineNew[,1]), max(maskOutlineNew[,1]), min(maskOutlineNew[,2]), max(maskOutlineNew[,2]))
imageTr <- raster::crop(imageTr, cropEx)
raster::extent(imageTr) <- cropEx
}
imageTr[is.na(imageTr)] <- 0
# Transform to imager format
imA <- raster::as.array(imageTr)
imR <- as.raster(imA, nrow = dim(imageTr)[1], ncol = dim(imageTr)[2], max = 255)
im <- imager::as.cimg(imR)
if(n == 1){
im1 <- imager::as.cimg(imR)
}
# Correct illumination using linear model
if(correct){
if(plotCorrect){
layout(t(1:2))
par(mar=c(2, 2, 2, 2))
plot(im, main = paste(names(sampleList)[n], 'original', sep=' '))
}
d <- as.data.frame(im)
# Subsamble, fit a linear model
m <- sample_n(d, 1e4) %>% lm(value ~ x*y, data=.)
# Correct by removing the trend
im <- im-predict(m,d)
if(plotCorrect) plot(im, main = paste(names(sampleList)[n], 'corrected', sep=' '))
}
# # Blur image for color extraction
# if(focal){
#
# gf <- focalWeight(image, sigma, "Gauss")
#
# rrr1 <- raster::focal(image[[1]], gf)
# rrr2 <- raster::focal(image[[2]], gf)
# rrr3 <- raster::focal(image[[3]], gf)
#
# image <- raster::stack(rrr1, rrr2, rrr3)
#
# image[is.na(image)] <- 255
# }
# Blur image for priority map extraction
if(blur){
edges <- isoblur(im, sigma) %>% imgradient("xy") %>% enorm %>% imsplit("c") %>% add %>% sqrt
}
else{
edges <- imgradient(im,"xy") %>% enorm
}
layout(t(1:1))
if(plotEdges) plot(edges)
priority <- 1/(255+edges)
if(plotPriority) plot(priority)
# Extract seed points from image
seed <- imager::imfill(dim=dim(im)) #Empty image
}
layout(t(1:1))
par(mar=c(2, 2, 2, 2))
plot(im)
print("Choose points to identify patterns. Click outside image area to stop.")
while(1){
xy <- locator(n=1)
if(as.numeric(xy)[1] > dim(im)[1] || as.numeric(xy)[1] < 0 || as.numeric(xy)[2] > dim(im)[2] || as.numeric(xy)[2] < 0){
break
}
print(paste('x: ', as.character(xy)[1], 'y: ', as.character(xy)[2]))
seed[as.numeric(xy[1]),as.numeric(xy[2]),1,1] <- 2
}
print("Choose at least one point to identify background. Click outside image area to stop.")
while(1){
xy <- locator(n=1)
if(any(c(as.numeric(xy)[1] > dim(im)[1], as.numeric(xy)[1] < 0, as.numeric(xy)[2] > dim(im)[2], as.numeric(xy)[2] < 0))){
break
}
print(paste('x: ', as.character(xy)[1], 'y: ', as.character(xy)[2]))
xyBF <- xy
seed[as.numeric(xy[1]),as.numeric(xy[2]),1,1] <- 1
}
# Run watershed transform
ws <- (imager::watershed(seed, priority) == 2) # multiply by 255 if RGB extracted
if(plotWS) plot(ws)
# Use a bucket fill on the background to fill holes (starts from last selected bg pixel)
if(bucketfill){
ws <- imager::bucketfill(ws, as.numeric(xyBF[1]), as.numeric(xyBF[2]), color=2) %>% {!( . == 2)}
}
if(plotBF) plot(ws)
# Remove spurious areas with width < cleanP
if(!is.null(cleanP)){
ws <- imager::clean(ws,cleanP)
}
# Split into connected components and remove ones with areas < splitC
if(!is.null(splitC)){
ws <- imager::split_connected(ws) %>% purrr::discard(~ sum(.) < splitC) %>% parany
}
# Transform to raster
wsR <- as.raster(ws)
wsM <- matrix(t(col2rgb(wsR))[,1]/255, nrow = dim(im)[1], ncol = dim(im)[2])
patternRaster <- raster::raster(t(wsM))
raster::extent(patternRaster) <- raster::extent(0, dim(im)[1], 0, dim(im)[2])
if(plotFinal){
layout(t(1:2))
plot(im, main="Watershed")
highlight(ws, col="red")
plot(patternRaster, legend = FALSE, bty = 'n', main = 'Extracted pattern')
}
todo <- readline(prompt="Press [enter] to continue, 'a' to add points 'r' to redo or 'q' to quit >>> ")
if(todo == 'q') break
if(todo == 'a') next
if(todo == 'r'){
seed <- imager::imfill(dim=dim(im))
next
}
else{
y <- raster::raster(ncol = dim(im1)[2], nrow = dim(im1)[1])
extent(y) <- raster::extent(0, dim(im)[1], 0, dim(im)[2])
patternRaster <- resample(patternRaster, y, method = 'ngb')
raster::extent(patternRaster) <- raster::extent(target)
rasterList[[names(sampleList)[n]]] <- patternRaster
print(paste('sample', names(sampleList)[n], 'done and added to rasterList', sep=' '))
n = n + 1
}
}
return(rasterList)
}
###
# The following code is added to suppress raising a R CMD check note resulting from magrittr '.'
###
# Placeholder function.
as_nonstandard <- function()
{
stop("This is a placeholder function, and should not be called directly.")
}
# Declare symbols for nonstandard evalution.
declare <- function(declaration)
{
d <- substitute(declaration)
p <- parent.frame()
if (!identical(d[[1]], quote(`<-`)))
stop("Declarations follow the <- syntax.", call. = FALSE)
if (!is.symbol(d[[2]]))
stop("LHS needs to be a symbol", call. = FALSE)
sym <- as.character(d[[2]])
msg <- sprintf("Symbol '%s' is declared for non-standard use.", sym)
eval(
substitute(
delayedAssign(x, stop(m, call. = FALSE), environment(), environment()),
list(m = msg, x = sym)), p, p)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.