Darwin PC - 1 year ago 102
R Question

# How to draw a shape (ellipse or oval) following some points and calculate its area?

I am trying to plot rings of trees and calculate their areas. However, I have noticed that in reality not all rings have symmetric radii like a circle. I have data measurements of 4 radii, and I would like to plot rings (or any similar shape) following each point of every radio like this example (this figure was done manually with vectors in PowerPoint):

the problem is that in R I found only the possibility to plot these rings with the

`circles`
option from the
`symbols()`
function, and I got this graph:

using this R script:

``````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)
``````

I have checked some packages with functions to draw ellipses (
`elliplot`
,
`ellipse`
,
`ellipseplot`
,
`car`
, etc), but I don't like their functions. I am not interested in use these packages, on the contrary I would like to write an own code.

My idea is to plot a shape which best meets the real figure of a ring with my data values of the four radii, it can be an ellipse, oval, etc.

With a circle I am using only data of one radio (in my example, the mean of all radii).
With a ellipse would be better, because I can use at least two values, the major-axis (A+B), and the minor-axis (C+D). But would be great to draw a shape that use the values of four radii (A, B, C, D) or even more radii.

Here a guy drew a very nice superellipse using a R script, and another one drew some ellipses likes rings also in R.

However, I don't know how to use their methods to my specific problem.

If somebody have idea how to start drawing at least an ellipse in R would be nice. But would be great to know how to draw a shape (oval, ellipse, etc.) using the values of four radii and finally calculate their area.

I would appreciate very much your help or any direction to do that.

If A & B are on y-axis and C & D are on x-axis, it isn't difficult to calculate the parameters of ellipses. I used `optim()` to get params (Note: this approach has tiny error, such as 2.439826e-12).

data manipulation
`````` # 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
``````
calculation of center coordinates, ox & oy
``````center <- df %>% group_by(Var2) %>% summarize(sum(x)/2, sum(y)/2) %>% as.data.frame()
``````
calculation of parameters of ellipse; semi-major and -minor axis, 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
``````
function to plot (Probably some packages have similar one)
``````radian <- function(degree) degree/180*pi
plot.ellipse <- function(ox, oy, ra, rb, phi=0, start=0, end=360, length=100, func=lines, ...) {
if (phi == 0) {
func(ra*cos(theta)+ox, rb*sin(theta)+oy, ...)
} else {
x <- ra*cos(theta)
y <- rb*sin(theta)
cosine <- cos(phi)
sine <- sin(phi)
func(cosine*x-sine*y+ox, sine*x+cosine*y+oy, ...)
}
}
``````
draw
``````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])
``````

Recommended from our users: Dynamic Network Monitoring from WhatsUp Gold from IPSwitch. Free Download