如何根据一些点绘制形状(椭圆或椭圆形)并计算其面积?
How to draw a shape (ellipse or oval) following some points and calculate its area?
我正在尝试绘制树木的年轮并计算它们的面积。但是,我注意到实际上并不是所有的环都像圆一样具有对称的半径。我有 4 个半径的数据测量值,我想像这个例子一样在每个收音机的每个点之后绘制环(或任何类似的形状)(这个图是用 PowerPoint 中的矢量手动完成的):
问题是,在 R 中,我发现只能使用 symbols()
函数中的 circles
选项来绘制这些环,我得到了这张图:
使用这个 R 脚本:
data <- data.frame(
a = c(1,4,5,8, 10),
b = c(1, 3,7,9, 10),
c = c(2, 6, 8, 9 ,10),
d = c(1, 3, 4, 7, 9) )
data$y <- (data$a - data$b)/2 # y position
data$x <- (data$d - data$c)/2 # x position
data$z <- rowMeans(data[,1:4]) # radio length
symbols(x = data$x, y = data$y, circles=data$z,
xlim = c(-10, 10)*1.5, ylim = c(-10, 10)*1.5, inches = F, fg = "orange", lwd = 2)
我检查了一些带有绘制椭圆函数的包(elliplot
、ellipse
、ellipseplot
、car
等),但我不喜欢它们职能。我对使用这些包不感兴趣,相反我想写一个自己的代码。
我的想法是用我的四个半径的数据值绘制一个最符合圆环真实图形的形状,可以是椭圆,椭圆等。
对于一个圆,我只使用一个无线电的数据(在我的例子中,所有半径的平均值)。
使用椭圆会更好,因为我可以使用至少两个值,major-axis (A+B) 和 minor-axis (C+D)。但是绘制一个使用四个半径(A、B、C、D)或更多半径的值的形状会很棒。
Here a guy drew a very nice superellipse using a R script, and another one drew some ellipses likes rings also in R。
但是,我不知道如何使用他们的方法来解决我的具体问题。
如果有人知道如何开始在 R 中至少绘制一个椭圆,那就太好了。但是如果知道如何使用四个半径的值绘制形状(椭圆形、椭圆形等)并最终计算它们的面积,那就太好了。
非常感谢你的帮助或任何指示。
更新:
感谢@cuttlefish44 的出色回答,这对向我的学生解释树木生长非常有用。然而,大多数热带树木的形状非常不规则,现在我想知道我是否可以使用额外的无线电 "E" 和不同位置的半径轴来绘制其他形状,就像这个方案:
任何方向对我都非常有用。
如果A&B在y轴上,C&D在x轴上,那么计算椭圆的参数就不难了。我使用 optim()
获取参数(注意:这种方法有微小的错误,例如 2.439826e-12)。
数据操作
# change all data into xy coordinates and make ring-factor
library(reshape2); library(dplyr)
data <- data.frame(
a = c(1, 4, 5, 8, 10),
b = c(1, 3, 7, 9, 10) * -1,
c = c(2, 6, 8, 9, 10) * -1,
d = c(1, 3, 4, 7, 9) )
data <- t(data)
colnames(data) <- LETTERS[1:ncol(data)] # ring-factor
df <- melt(data, value.name = "x") # change into long-form
df$y <- df$x # make xy coordinates
df[df$Var1=="a"|df$Var1=="b", "x"] <- 0
df[df$Var1=="c"|df$Var1=="d", "y"] <- 0
中心坐标的计算,ox & oy
center <- df %>% group_by(Var2) %>% summarize(sum(x)/2, sum(y)/2) %>% as.data.frame()
椭圆参数的计算;半长轴和半短轴,ra & rb
opt.f <- function(par, subset, center) { # target function
ox <- center[[1]] # par[1] and par[2] are ra and rb
oy <- center[[2]]
x <- subset$x
y <- subset$y
sum(abs((x - ox)^2/par[1]^2 + (y - oy)^2/par[2]^2 - 1)) # from ellipse equation
}
lev <- levels(df$Var2)
## search parameters
res <- sapply(1:length(lev), function(a)
optim(c(1,1), opt.f, subset = subset(df, Var2 == lev[a]),
center = center[a, 2:3], control = list(reltol = 1.0e-12)))
res # result. you can get detail by res[,1etc]. values are not 0 but much nearly 0
绘制函数(可能有些包有类似的)
radian <- function(degree) degree/180*pi
plot.ellipse <- function(ox, oy, ra, rb, phi=0, start=0, end=360, length=100, func=lines, ...) {
theta <- c(seq(radian(start), radian(end), length=length), radian(end))
if (phi == 0) {
func(ra*cos(theta)+ox, rb*sin(theta)+oy, ...)
} else {
x <- ra*cos(theta)
y <- rb*sin(theta)
phi <- radian(phi)
cosine <- cos(phi)
sine <- sin(phi)
func(cosine*x-sine*y+ox, sine*x+cosine*y+oy, ...)
}
}
画
plot(0, type="n", xlim=c(-10, 10), ylim =c(-10, 10), asp=1, xlab="x", ylab="y", axes = F)
axis(1, pos=0);axis(2, pos=0, las=2)
points(df$x, df$y)
for(a in 1:length(lev)) plot.ellipse(ox = center[a, 2], oy = center[a, 3],
ra = res[,a]$par[1], rb = res[,a]$par[2], length=300)
area <- sapply(res[1,], function(a) pi * a[1] * a[2])
我正在尝试绘制树木的年轮并计算它们的面积。但是,我注意到实际上并不是所有的环都像圆一样具有对称的半径。我有 4 个半径的数据测量值,我想像这个例子一样在每个收音机的每个点之后绘制环(或任何类似的形状)(这个图是用 PowerPoint 中的矢量手动完成的):
问题是,在 R 中,我发现只能使用 symbols()
函数中的 circles
选项来绘制这些环,我得到了这张图:
使用这个 R 脚本:
data <- data.frame(
a = c(1,4,5,8, 10),
b = c(1, 3,7,9, 10),
c = c(2, 6, 8, 9 ,10),
d = c(1, 3, 4, 7, 9) )
data$y <- (data$a - data$b)/2 # y position
data$x <- (data$d - data$c)/2 # x position
data$z <- rowMeans(data[,1:4]) # radio length
symbols(x = data$x, y = data$y, circles=data$z,
xlim = c(-10, 10)*1.5, ylim = c(-10, 10)*1.5, inches = F, fg = "orange", lwd = 2)
我检查了一些带有绘制椭圆函数的包(elliplot
、ellipse
、ellipseplot
、car
等),但我不喜欢它们职能。我对使用这些包不感兴趣,相反我想写一个自己的代码。
我的想法是用我的四个半径的数据值绘制一个最符合圆环真实图形的形状,可以是椭圆,椭圆等。
对于一个圆,我只使用一个无线电的数据(在我的例子中,所有半径的平均值)。 使用椭圆会更好,因为我可以使用至少两个值,major-axis (A+B) 和 minor-axis (C+D)。但是绘制一个使用四个半径(A、B、C、D)或更多半径的值的形状会很棒。
Here a guy drew a very nice superellipse using a R script, and another one drew some ellipses likes rings also in R。
但是,我不知道如何使用他们的方法来解决我的具体问题。
如果有人知道如何开始在 R 中至少绘制一个椭圆,那就太好了。但是如果知道如何使用四个半径的值绘制形状(椭圆形、椭圆形等)并最终计算它们的面积,那就太好了。
非常感谢你的帮助或任何指示。
更新:
感谢@cuttlefish44 的出色回答,这对向我的学生解释树木生长非常有用。然而,大多数热带树木的形状非常不规则,现在我想知道我是否可以使用额外的无线电 "E" 和不同位置的半径轴来绘制其他形状,就像这个方案:
任何方向对我都非常有用。
如果A&B在y轴上,C&D在x轴上,那么计算椭圆的参数就不难了。我使用 optim()
获取参数(注意:这种方法有微小的错误,例如 2.439826e-12)。
# change all data into xy coordinates and make ring-factor
library(reshape2); library(dplyr)
data <- data.frame(
a = c(1, 4, 5, 8, 10),
b = c(1, 3, 7, 9, 10) * -1,
c = c(2, 6, 8, 9, 10) * -1,
d = c(1, 3, 4, 7, 9) )
data <- t(data)
colnames(data) <- LETTERS[1:ncol(data)] # ring-factor
df <- melt(data, value.name = "x") # change into long-form
df$y <- df$x # make xy coordinates
df[df$Var1=="a"|df$Var1=="b", "x"] <- 0
df[df$Var1=="c"|df$Var1=="d", "y"] <- 0
中心坐标的计算,ox & oy
center <- df %>% group_by(Var2) %>% summarize(sum(x)/2, sum(y)/2) %>% as.data.frame()
椭圆参数的计算;半长轴和半短轴,ra & rb
opt.f <- function(par, subset, center) { # target function
ox <- center[[1]] # par[1] and par[2] are ra and rb
oy <- center[[2]]
x <- subset$x
y <- subset$y
sum(abs((x - ox)^2/par[1]^2 + (y - oy)^2/par[2]^2 - 1)) # from ellipse equation
}
lev <- levels(df$Var2)
## search parameters
res <- sapply(1:length(lev), function(a)
optim(c(1,1), opt.f, subset = subset(df, Var2 == lev[a]),
center = center[a, 2:3], control = list(reltol = 1.0e-12)))
res # result. you can get detail by res[,1etc]. values are not 0 but much nearly 0
绘制函数(可能有些包有类似的)
radian <- function(degree) degree/180*pi
plot.ellipse <- function(ox, oy, ra, rb, phi=0, start=0, end=360, length=100, func=lines, ...) {
theta <- c(seq(radian(start), radian(end), length=length), radian(end))
if (phi == 0) {
func(ra*cos(theta)+ox, rb*sin(theta)+oy, ...)
} else {
x <- ra*cos(theta)
y <- rb*sin(theta)
phi <- radian(phi)
cosine <- cos(phi)
sine <- sin(phi)
func(cosine*x-sine*y+ox, sine*x+cosine*y+oy, ...)
}
}
画
plot(0, type="n", xlim=c(-10, 10), ylim =c(-10, 10), asp=1, xlab="x", ylab="y", axes = F)
axis(1, pos=0);axis(2, pos=0, las=2)
points(df$x, df$y)
for(a in 1:length(lev)) plot.ellipse(ox = center[a, 2], oy = center[a, 3],
ra = res[,a]$par[1], rb = res[,a]$par[2], length=300)
area <- sapply(res[1,], function(a) pi * a[1] * a[2])