2014-04-10 56 views
2

我收集了关于船舶对海鸟造成干扰的数据。我在船上配有测距望远镜和角板。对于我所调查的每只鸟,我都有一个相对于船舶航向的起始距离和方位。我也有鸟的反应距离和方位(或在某些情况下没有)。绘图距离和方位在R

我想做一个两个面板的情节显示一个起始距离和轴承位置,另一个终止距离和轴承。理想情况下,第二个小区将进行彩色编码(或pch编码)以显示不同的反应类型。

我的数据是这种格式

 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 

这里是一个格式的数据,你可以玩弄

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) 

正如你可以看到有一些ID的,即重复和一些有只有一行。我想绘制第一个dist,并在另一个情节上承担一个情节和最后一个dist并承担(每个id)。对于只有一次观察的鸟类,这些可能是相同的。根据“行为”列对第二个图中的点进行颜色编码会很好。也没有左或右的轴承名称,所以我可以把所有的点都放在中间线的另一侧,但是如果你知道将它们随机放置在中线的左侧还是右侧将会很酷。理想情况下,情节会看起来像这样。

Dist Bear Plot

UPDATE:从@jbaums使用他的代码从另一个问题的建议之后发现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) 

这给了我这个越来越接近我的目标的阴谋!谢谢@jbaums。

closer

我的问题是,我无法弄清楚如何做到绘制从0 90楔到90(因为这是我的数据收集。

我也还是需要一些指导只有当一个以上的观察已经收集,选择第一(后来最后)观察

+0

将两张图组合成一张,使用不同的颜色作为同一只鸟的起点和终点,并用直线连接起点和终点,会更好吗? – jinlong

+0

我可能会做的第一件事就是计算每次观察的X和Y距离。像'calcXY < - 函数(tdat)dist < - as.numeric(tdat [2]) bear < - as.numeric(tdat [3]) bear.rad < - bear * pi/180 Y < - cos(bear.rad)* dist X < - sin(bear.rad)* dist return(c(X,Y)) } XYPos < - apply(dat,1,calcXY)'Sorry about代码全部在一行上。应该清楚间距应该在哪里。我不认为这是合理的回答你的问题 – Adrian

+0

@Adrian:使用分号分隔注释中的代码行。 :) – jbaums

回答

1

如果你想更仔细地重新创建示例图,请尝试以下,使用您的datget.coords功能最初发布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') 

dist & bearing

请注意,我已经通过了由于您的示例图表明您正在从正Y轴逆时针计算轴承,因此您的指向get.coords的点的角度为-bear。函数get.coords函数期望从正x轴顺时针计算角度,负角度(如将出现在-bear中)将被解释为360减去角度。

+1

非常感谢。我用jinlong提供的代码修改了一些代码,并且非常接近你在这里完成的代码,但是肯定会包含一些你的代码。由于我有超过4000点的阴谋,我选择沿着X轴放置仪表标记,因为它们将被点覆盖。再次感谢。 – marcellt

1

不知道如果我理解您的所有要求,但下面是我对“出发点”的情节的解决方案:

#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) 
    } 
} 

如果这是你想要的,你可以适应它的“终点”情节。你可以添加一个col参数给point()函数,并使用“act”对点进行颜色编码。

+0

感谢@jinlong,该代码为我想要的完美工作。我很感谢帮助。 – marcellt

+1

@jbaums我会留下来自你以前的文章的代码,因为它让我获得了2/3的目标,我认为它可以帮助其他人。 – marcellt