## Friday, April 27, 2012

### Create polygons from a matrix

The following function matrix.poly allows for the addition of polygons to a plot based on a matrix and defined matrix positions. I have used this function on occasion to highlight specific matrix locations (e.g. in the above figure). You can do the same by overlaying another image (left in above plot) but with this function you will have all other polygon plotting possibilities (e.g. borders etc.).
Another useful application is in the creation of polygons from a matrix representing a data field on a regular grid. The generated polygons can then be easily drawn onto maps of various projection.

matrix.poly function...
```matrix.poly <- function(x, y, z=mat, n=NULL){
if(missing(z)) stop("Must define matrix 'z'")
if(missing(n)) stop("Must define at least 1 grid location 'n'")
if(missing(x)) x <- seq(0,1,,dim(z)[1])
if(missing(y)) y <- seq(0,1,,dim(z)[2])
poly <- vector(mode="list", length(n))
for(i in seq(n)){
ROW <- ((n[i]-1) %% dim(z)[1]) +1
COL <- ((n[i]-1) %/% dim(z)[1]) +1

dist.left <- (x[ROW]-x[ROW-1])/2
dist.right <- (x[ROW+1]-x[ROW])/2
if(ROW==1) dist.left <- dist.right
if(ROW==dim(z)[1]) dist.right <- dist.left

dist.down <- (y[COL]-y[COL-1])/2
dist.up <- (y[COL+1]-y[COL])/2
if(COL==1) dist.down <- dist.up
if(COL==dim(z)[2]) dist.up <- dist.down

xs <- c(x[ROW]-dist.left, x[ROW]-dist.left, x[ROW]+dist.right, x[ROW]+dist.right)
ys <- c(y[COL]-dist.down, y[COL]+dist.up, y[COL]+dist.up, y[COL]-dist.down)
poly[[i]] <- data.frame(x=xs, y=ys)
}
return(poly)
}```
Created by Pretty R at inside-R.org

to reproduce the above figures...
```###
#required packages
library(maps)
library(mapproj)

source("matrix.poly.R")
source("val2col.R")

png("matrix.poly.ex.png", width=9, height=3, units="in", res=200)
par(mfcol=c(1,3), mar=c(5,5,0.5,0.5))
set.seed(1)

#Ex 1 - add another image layer
m=8
n=10
x <- seq(m)
y <- seq(n)
z <- matrix(runif(m*n), nrow=m, ncol=n)
image(x,y,z, col=grey.colors(20))
N <- sample(1:(m*n),20)
z2 <- NaN*z
z2[N] <- 1
box()

m=8
n=10
x <- seq(m)
y <- seq(n)
z <- matrix(runif(m*n), nrow=m, ncol=n)
image(x,y,z, col=grey.colors(20))
N <- sample(1:(m*n),20)
poly <- matrix.poly(x,y,z=z,n=N)
sapply(poly, function(X){polygon(X, col=rgb(1,1,0,0.4), border=1)})
box()

#Ex 2 - add polygons to unequal grid
x <- cumsum(round(runif(m, min=1, max=10)))
y <- cumsum(round(runif(n, min=1, max=10)))
z <- matrix(runif(m*n), nrow=m, ncol=n)
image(x,y,z, col=grey.colors(20))
N <- sample(1:(m*n),20)
poly <- matrix.poly(x,y,z=z,n=N)
sapply(poly, function(X){polygon(X, col=rgb(1,0,0,0.4), border=1)})
box()

dev.off()

####
png("matrix.poly.ex2.png", width=6, height=3, units="in", res=400)
par(mfcol=c(1,2), mar=c(5,5,0.5,0.5))
set.seed(1)

#Ex 1 - simple image
m=10
n=10
x <- seq(from=-60, to=60,, m)
y <- seq(from=30, to=85,, n)
z <- matrix(runif(m*n), nrow=m, ncol=n)
image(x,y,z, col=rainbow(20))

#Ex.2 project polygons
par(mar=c(0,0,0,0))
map("world", proj="orthographic", orient=c(60,0,20), par=NULL, t="n")
poly <- matrix.poly(x,y,z=z,n=seq(length(z)))
COL <- val2col(z,col=rainbow(20))
for(i in seq(poly)){
polygon(mapproject(poly[[i]]), col=COL[i], border=COL[i], lwd=0.3)
}