我收集了有关船舶海鸟干扰的数据。我在船上带着测距双筒望远镜和角度板。对于我调查的每只鸟,我都有一个相对于船舶航向的起始距离和方位。我还有鸟做出反应(或在某些情况下没有反应)的距离和方位。
我想制作一个两面板图,其中一个显示起始距离和轴承位置,另一个显示终止距离和轴承。理想情况下,第二个图将进行颜色编码(或 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
Run Code Online (Sandbox Code Playgroud)
这是您可以使用的格式的数据
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)
Run Code Online (Sandbox Code Playgroud)
正如您所看到的,有些 id 是重复的,有些则只有一行。我想在一个地块上绘制第一个距离和熊,在另一个地块上绘制最后一个距离和熊(每个 ID)。对于只有一次观察的鸟类来说,这些可能是相同的。最好根据“act”列对第二个图中的点进行颜色编码。此外,方位没有左侧或右侧的指定,所以我可以接受所有点都位于中线的一侧或另一侧,但如果您知道将它们随机放置在中心线的左侧或右侧会很酷。理想情况下,绘图看起来像这样。

更新:遵循@jbaums 的建议,使用他在此处找到的另一个问题中的代码。
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)
Run Code Online (Sandbox Code Playgroud)
这给了我这个绘制的数字,它越来越接近我的目标!谢谢@jbaums。

我的问题是,我无法弄清楚如何绘制从 0 到 90 的 90 楔形(因为这是我收集数据的地方。
当收集了多个观察值时,我还需要一些指导来仅选择第一个(以及后来的最后一个)观察值。
如果您想更仔细地重新创建示例图,请尝试以下操作,它使用您的dat, 以及get.coords最初发布于此处的函数:
# 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')
Run Code Online (Sandbox Code Playgroud)

请注意,我已将您的点的角度传递给get.coordsas -bear,因为您的示例图表明您正在从正 y 轴逆时针计算方位角。该get.coords函数期望从正 x 轴顺时针计算角度,负角度(如 中出现的-bear)将被解释为 360 减去角度。
| 归档时间: |
|
| 查看次数: |
2682 次 |
| 最近记录: |