edges.R

## Jon-Michael Deldin
## Pattern Recognition
## Edge Detection
## Spring 2013
##
##
## USAGE
## ------
##  Rscript edges.R
##
## DEPENDENCIES
## ------------
## You'll need to install the libtiff-dev package on Linux and then install
## the biOps package in R:
##
##  install.packages('biOps')
##

## clear existing vars
rm(list = ls(all=TRUE))

source('utilities.R')
source('functions.R')

library(biOps)

mkdirp('data')
mkdirp('fig')

############
## VIOLET ##
############
violet <- readJpeg(system.file("samples", "violet.jpg", package="biOps"))

plotSobel(violet, 'violet-sobel-white.jpg')
plotCanny(violet, 'violet-canny.jpg')

####################
## ORIGINAL LENNA ##
####################

maybeDownload('data/lenna.jpg', 'http://www.cs.umt.edu/~dougr/PattRecFiles/projects/edgeDetection/Lenna.jpg')
rgb <- readJpeg('data/lenna.jpg')
greyImg <- imgRGB2Grey(rgb)

plotImage(greyImg, 'lenna.jpg')

####################
## LAPLACIAN MASK ##
####################
lpMask <- getLaplacianMask(5)
lpImg <- applyMask(lpMask, greyImg)
plotImage(lpImg, 'lenna-laplacian.jpg')

###################
## GAUSSIAN BLUR ##
###################

## quantiles for a 5x5 matrix
qVals <- distanceFromCenter(2)

## gaussian weights with a stddev of 1.4
blurMask <- round(15 * gaussWeights(qVals, 1.4))

## normalize the weights by the sum of the values
blurMask <- blurMask/sum(blurMask)

smoothedImg <- applyMask(blurMask, greyImg)
plotImage(smoothedImg, 'lenna-gaussian.jpg')

smoothedLpImg <- applyMask(lpMask, smoothedImg)
plotImage(smoothedLpImg, 'lenna-gaussian-laplacian.jpg')

higherThresholdImg <- applyMask(lpMask, smoothedImg, lower=100)
plotImage(higherThresholdImg, 'lenna-gaussian-100.jpg')

###########
## CANNY ##
###########
plotImage(imgCanny(greyImg, 1.4), 'lenna-canny.jpg')

functions.R

## Jon-Michael Deldin
## Pattern Recognition
## Edge Detection: Functions
## Spring 2013
##

library(biOps)
source('utilities.R')

## plot sobel detection with a white background
plotSobel <- function(image, filename) {
  plotImage(imgNegative(imgSobel(image)), filename)
}

## plot canny detection with a stddev of 1.4
plotCanny <- function(image, filename, stddev=1.4) {
  plotImage(imgCanny(image, stddev, 1, 255), filename)
}

## Returns the midpoint
getCenter <- function(size) {
  ceiling(size/2)
}

## Returns a Laplacian convolution mask of a given size.
getLaplacianMask <- function(size) {
  mask <- matrix(-1, size, size)
  center <- getCenter(size)
  mask[center, center] <- size * size - 1
  return(mask)
}

## first two rows, last two rows, and anything that will put us out of bounds
isBoundary <- function(coord, maxBorder, interval) {
  bad <- coord == 1 || coord == 2 || coord == maxBorder - 2 || coord == maxBorder
  out <- (coord + interval - 1) > maxBorder

  ## printf("...coord = %d, maxb = %d, interval=%d", coord, maxBorder, interval)

  return(bad || out)
}

## Normalize the color values
normalize <- function(x, upper, lower) {
  if (x > upper) { 255 }
  else if (x < lower) { 0 }
  else { x }
}

## Apply a convolution mask to an image
applyMask <- function(mask, imageData, upper=255, lower=0) {
  print("Applying mask...")
  interval <- ncol(mask)

  numRows <- nrow(imageData)
  numCols <- ncol(imageData)
  copy <- imageData
  for (row in 1:numRows) {
    for (col in 1:numCols) {
      ## if we're looking at one of the ignoreable zones, just zero out the pixel
      if (isBoundary(row, numRows, interval) || isBoundary(col, numCols, interval)) {
        copy[row, col] <- 0
      } else {
        ## otherwise, do convolution
        slice <- imageData[row:(row+interval-1), col:(col+interval-1)]
        ## highlight low-intensity edges near white
        pixel <- upper - normalize(sum(mask * slice), upper, lower)
        copy[row, col] <- pixel
      }
    }
  }
  return(imagedata(copy))
}

## Distance between two points
distance <- function(x1, y1, x2, y2) {
  sqrt((x2 - x1)^2 + (y2 - y1)^2)
}

## return a square distance matrix from a center point
distanceFromCenter <- function(distFromCenter) {
  n <- distFromCenter * 2 + 1
  mat <- matrix(ncol=n, nrow=n)
  center <- getCenter(n)

  for (i in 1:n) {
    for (j in 1:n) {
      mat[i, j] <- distance(i, j, center, center)
    }
  }
  return(mat)
}

## Returns the approximate Gaussian weights for a given radius and stddev
gaussWeights <- function(r, sigma) {
  exp(-(r**2) / (2 * sigma**2))
}

## z <- applyMask(GWs/sum(GWs), greyImg)
## zi <- z
## plot(applyMask(lpmask, zi, lower=100))

utilities.R

##
## Handy utilities to share across assignments.
##

library(MASS)

## create the directory to download the data (like `mkdir -p` on *nix)
mkdirp <- function(directory) {
  dir.create(directory, recursive=TRUE, showWarnings=FALSE)
}

## download a file if it doesn't exist
maybeDownload <- function(filename, url) {
  if (!file.exists(filename)) {
    download.file(url, filename)
  }
}

## drop a specific element from a list
dropElement <- function(elt, lst) {
  lst[-which(lst == elt)]
}

## runs lda() and returns an LDA object
getLda <- function(data, train, classes, tolerance) {
  return(lda(data[train,], classes[train], tol=tolerance))
}

## returns both projections
projectLda <- function(data, ldaObj) {
  mat <- as.matrix(data)
  first <- mat %*% ldaObj$scaling[,1]
  second <- mat %*% ldaObj$scaling[,2]
  return(list(first=first, second=second))
}

## Prints a formatted string
printf <- function(str, ...) {
  cat(sprintf(str, ...), "\n")
}

## Returns the Euclidean distance between two columns
eucDist <- function(x1, x2) {
  sqrt(sum((x1 - x2) ^ 2))
}

## sets the initial width, height, and margins
initPlot <- function(bottom=2.5, left=2.5, top=0.5, right=0.1) {
  X11(width=4, height=4)
  par(mar=c(bottom, left, top, right), mgp=c(1.5, .5, 0))
}

## save the plot to PDF
savePlot <- function(filename, type=pdf) {
  dev.copy(type, paste('fig/', filename, sep=''))
  dev.off()
}

## set sane defaults for plotting an image
initImagePlot <- function() {
  initPlot()
  par(mar=c(0,0,0,0))
}

## save a plot as a JPEG
saveJpeg <- function(filename) {
  savePlot(filename, jpeg)
}

## plot an imagedata() instance
plotImage <- function(image, filename) {
  initImagePlot()
  plot(image)
  saveJpeg(filename)
}