plot distance and bearing in R

Question

I have data collected on seabird disturbance from ships. I was on board ships with range finder binoculars and an angle board. For each bird I surveyed I have a starting distance and bearing relative to the ships course. I also have the distance and bearing at which the bird reacted (or didn't in some cases).

I would like to make a two panel plot showing on one the starting distance and bearing positions and on the other the terminating distance and bearings. Ideally the second plot will be color coded (or pch coded) to show the different reaction type.

My data is in this format

      date_id dist bear act
550 40711_027  200   30   f
551 40711_028  500   45   n
552 40711_028  450   60   n
553 40711_028  400   75   n
554 40711_028  371   80   f
555 40711_029  200    5   f
556 40711_030  200   10   d
557 40711_031  400   30   n
558 40711_031  350   30   d

Here is the data in a format you can play around with

id <- c(1,2,2,2,2,3,4,5,5)
dist <- c(200,500,450,400,371,200,200,400,350)
bear <- c(30,45,60,75,80,5,10,30,30)
act <- c("f","n","n","n","f","f","d","n","d")

dat <- data.frame(id, dist, bear, act)

As you can see there are some id's that repeat and some that have only one row. I would like to plot the first dist and bear on one plot and the last dist and bear (per id) on another plot. These may be the same for birds with only one observation. It would be nice to color code the points in the second plot based on the 'act' column. Also there is no left or right designation for bearing so I am okay with all the points being on one side of the middle line or the other but if you know how it would be cool to randomly place them left or right of the center line. Ideally the plots will look something like this.

UPDATE: Following suggestions from @jbaums using his code from another question found here.

get.coords <- function(a, d, x0, y0) {
  a <- ifelse(a <= 90, 90 - a, 450 - a)
  data.frame(x = x0 + d * cos(a / 180 * pi), 
             y = y0+ d * sin(a / 180 * pi))
}

rotatedAxis <- function(x0, y0, a, d, symmetrical=FALSE, tickdist, ticklen, ...) {
  if(isTRUE(symmetrical)) {
    axends <- get.coords(c(a, a + 180), d, x0, y0)    
    tick.d <- c(seq(0, d, tickdist), seq(-tickdist, -d, -tickdist))      
  } else {
    axends <- rbind(get.coords(a, d, x0, y0), c(x0, y0))
    tick.d <- seq(0, d, tickdist)
  }
  invisible(lapply(apply(get.coords(a, d=tick.d, x0, y0), 1, function(x) {
    get.coords(a + 90, c(-ticklen, ticklen), x[1], x[2])
  }), function(x) lines(x$x, x$y, ...)))
  lines(axends$x, axends$y, ...)
}

plot.new()
plot.window(xlim=c(-1000,1000),ylim=c(-1000, 1000), asp=1) 
polygon(get.coords(seq(0,180, length.out=1000),1000,0,0),lwd=2)
polygon(get.coords(seq(0,180, length.out=750),750,0,0),lwd=2)
polygon(get.coords(seq(0,180, length.out=500),500,0,0),lwd=2)
polygon(get.coords(seq(0,180, length.out=250),250,0,0),lwd=2)

rotatedAxis(0, 0, a=90, d=1000, tickdist=100, ticklen=1)
rotatedAxis(0, 0, a=45, d=1000, tickdist=100, ticklen=1)
rotatedAxis(0, 0, a=135, d=1000, tickdist=100, ticklen=1)

obs <- with(dat, get.coords(bear, dist, 0, 0))
points(obs)

This gives me this plotted figure which is getting closer to my goal! Thanks @jbaums.

My issue is that I cannot figure out how to just plot the 90 wedge from 0 to 90 (as this is where my data was collected in.

I also still need some guidance on only selecting the first (and later the last) observation when more than one observations have been collected.

Answer 1

If you want to recreate your example plot more closely, try the following, which uses your dat, and the get.coords function originally posted here:

# Define function to calculate coordinates given distance and bearing
get.coords <- function(a, d, x0, y0) {
  a <- ifelse(a <= 90, 90 - a, 450 - a)
  data.frame(x = x0 + d * cos(a / 180 * pi), 
             y = y0+ d * sin(a / 180 * pi))
}

# Set up plotting device
plot.new()
par(mar=c(2, 0, 0, 0), oma=rep(0, 4))
plot.window(xlim=c(-1100, 1100), ylim=c(-100, 1100), asp=1)

# Semicircles with radii = 100 through 1000
sapply(seq(100, 1000, 100), function(x) {
  lines(get.coords(seq(270, 450, length.out=1000), x, 0, 0))
})

# Horizontal line
segments(-1000, 0, 1000, 0)

# 45-degree lines
apply(get.coords(c(360-45, 45), 1000, 0, 0), 1, 
      function(x) lines(rbind(x, c(0, 0)), lwd=2))

# Plot white curves over black curves and add text
sapply(seq(100, 1000, 100), function(x) {
  txt <- paste0(x, 'm')
  w <- strwidth(txt, cex=0.9)/2
  a <- atan(w/x)/pi*180
  lines(get.coords(seq(-a, a, length=100), x, 0, 0), 
        lwd=2.5, col='white')
  text(0, x, txt, cex=0.8)
})

# Add points
points(with(dat, get.coords(-bear, dist, 0, 0)), pch=20)

# Add triangle
polygon(c(0, -30, 30), c(-5, -55, -55), col='black')

Note that I've passed the angles of your points to get.coords as -bear, since your example figure suggests you are calculating bearings counter-clockwise from the positive y-axis. The get.coords function expects angles to be calculated clockwise from the positive x-axis, and negative angles (as will arise with -bear) will be interpreted as 360 minus the angles.

Answer 2

Not sure if I understand all your requirements but below is my solution for the "starting points" plot:

#install.packages("plotrix")
library("plotrix")

id <- c(1,2,2,2,2,3,4,5,5)
dist <- c(200,500,450,400,371,200,200,400,350)
bear <- c(30,45,60,75,80,5,10,30,30)
act <- c("f","n","n","n","f","f","d","n","d")

dat <- data.frame(id, dist, bear, act)

##Define a function that converts degrees to radians
#NOTE: Authored by Fabio Marroni
#URL: http://fabiomarroni.wordpress.com/2010/12/23/r-function-to-convert-degrees-to-radians/
degrees.to.radians<-function(degrees=45,minutes=30)
{
  if(!is.numeric(minutes)) stop("Please enter a numeric value for minutes!\n")
  if(!is.numeric(degrees)) stop("Please enter a numeric value for degrees!\n")
  decimal<-minutes/60
  c.num<-degrees+decimal
  radians<-c.num*pi/180
  return(radians)
}

#Plot the canvas
plot(0, 0, type = "n", xaxt = "n", yaxt = "n", asp=1,
     xlim = c(0, max(dat$dist)), ylim =  c(0, max(dist)), 
     bty="n", xlab = "", ylab = "", 
     main = "Whatever observations (starting points only)")

#Plot x/y axes
segments(0, 0, max(dat$dist), 0)
segments(0, 0, 0, max(dat$dist))

#Plot axes labels
axis(1, at = seq(0, max(dat$dist), 100), labels = seq(0, max(dat$dist), 100))

#Plot the equal-distance arcs
dist = 100
while(dist < max(dat$dist)){
  draw.arc(0, 0, radius = dist, deg1 = 0, deg2 = 90, n = 100, col = "blue")
  dist <- dist + 100
}

#Plot the 1st point (cause it's always an starting point)
x <- dat[1, ]$dist*sin(degrees.to.radians(dat[1, ]$bear))
y <- dat[1, ]$dist*cos(degrees.to.radians(dat[1, ]$bear))
points(x, y, pch = 21)

for(i in 2:nrow(dat)){
  #Only plot starting points
  if(dat[i, ]$id != dat[i-1, ]$id){
    #Determin the x and y for each point
    x <- dat[i, ]$dist*sin(degrees.to.radians(dat[i, ]$bear))
    y <- dat[i, ]$dist*cos(degrees.to.radians(dat[i, ]$bear))

    #Adding starting points
    points(x, y, pch = 21)
  }
}

If this is what you want, you can adapt it for the "ending points" plot. And you can add a col parameter to the point() function and use "act" to color code the points.