Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
937 views
in Technique[技术] by (71.8m points)

math - how to plot ellipse given a general equation in R?

Ellipse general equation:

a * x ^ 2 + b * y ^ 2 + c * x * y + d * x + e * y + f = 0

enter image description here

See Question&Answers more detail:os

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

We can start from the parametric equation of an ellipse (the following one is from wikipedia), we need 5 parameters: the center (xc, yc) or (h,k) in another notation, axis lengths a, b and the angle between x axis and the major axis phi or tau in another notation.

enter image description here

xc <- 1 # center x_c or h
yc <- 2 # y_c or k
a <- 5 # major axis length
b <- 2 # minor axis length
phi <- pi/3 # angle of major axis with x axis phi or tau

t <- seq(0, 2*pi, 0.01) 
x <- xc + a*cos(t)*cos(phi) - b*sin(t)*sin(phi)
y <- yc + a*cos(t)*cos(phi) + b*sin(t)*cos(phi)
plot(x,y,pch=19, col='blue')

enter image description here

Now if we want to start from the cartesian conic equation, it's a 2-step process.

  1. Convert the cartesian equation to the polar (parametric), form we can use the following equations to first obtain the 5 parameters using the 5 equations from the below figure (taken from http://www.cs.cornell.edu/cv/OtherPdf/Ellipse.pdf, detailed math can be found there).

  2. Plot the ellipse, by using the parameters obtained, as shown above.

enter image description here

For step (1) we can use the following code (when we have known A,B,C,D,E,F):

M0 <- matrix(c(F,D/2,E/2, D/2, A, B/2, E/2, B/2, C), nrow=3, byrow=TRUE)
M <- matrix(c(A,B/2,B/2,C), nrow=2)
lambda <- eigen(M)$values
abs(lambda - A)
abs(lambda - C) 

# assuming abs(lambda[1] - A) < abs(lambda[1] - C), if not, swap lambda[1] and lambda[2] in the following equations:

a <- sqrt(-det(M0)/(det(M)*lambda[1]))  
b <- sqrt(-det(M0)/(det(M)*lambda[2]))
xc <- (B*E-2*C*D)/(4*A*C-B^2)
yc <- (B*D-2*A*E)/(4*A*C-B^2)
phi <- pi/2 - atan((A-C)/B)*2

For step (2) use the following code:

t <- seq(0, 2*pi, 0.01) 
x <- xc + a*cos(t)*cos(phi) - b*sin(t)*sin(phi)
y <- yc + a*cos(t)*cos(phi) + b*sin(t)*cos(phi)
plot(x,y,pch=19, col='blue')

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...