# [R] drawing ellipses

Martin Maechler maechler at stat.math.ethz.ch
Wed Jun 19 09:45:01 CEST 2002

```>>>>> "Yves" == Yves Brostaux <brostaux.y at fsagx.ac.be> writes:

Yves> Hello again,

Yves> First I want to thank all the people who answered my
Yves> question about line width in graphs. I promise I will
Yves> learn the 'par' help page by heart for the end of the
Yves> month !

Yves> I now want to trace some ellipses to emphasize groups
Yves> of data. I found how to trace circles with
Yves> 'symbols()', but no ellipse. I'm planning on writing
Yves> my own function based on 'polygon()' and the ellipse
Yves> equation. Does anybody already have done similar work
Yves> or have another solution ?

yes. several.

1) The ellipse  package  does some.
2) The cluster  package (which is recommended, hence should be
available everywhere)
has an ellipsoidhull() function and a plot method for its results
(in the 2D case, i.e when the ellipsoid is an ellipse)

3) A `client' (researcher in another department here) needed
something  based on geometric parameters rather than
covariance matrices.
So I ended up doing the following, I hope it helps.

Martin Maechler <maechler at stat.math.ethz.ch>	http://stat.ethz.ch/~maechler/
Seminar fuer Statistik, ETH-Zentrum  LEO C16	Leonhardstr. 27
ETH (Federal Inst. Technology)	8092 Zurich	SWITZERLAND
phone: x-41-1-632-3408		fax: ...-1228			<><

-------------------------------------------------------------------

ellipsePoints <- function(a,b, alpha = 0, loc = c(0,0), n = 201)
{
## Purpose: ellipse points,radially equispaced, given geometric par.s
## -------------------------------------------------------------------------
## Arguments: a, b : length of half axes in (x,y) direction
##            alpha: angle (in degrees) for rotation
##            loc  : center of ellipse
##            n    : number of points
## -------------------------------------------------------------------------
## Author: Martin Maechler, Date: 19 Mar 2002, 16:26

B <- min(a,b)
A <- max(a,b)
## B <= A
d2 <- (A-B)*(A+B)                   #= A^2 - B^2
phi <- 2*pi*seq(0,1, len = n)
sp <- sin(phi)
cp <- cos(phi)
r <- a*b / sqrt(B^2 + d2 * sp^2)
xy <- r * cbind(cp, sp)
## xy are the ellipse points for alpha = 0 and loc = (0,0)
al <- alpha * pi/180
ca <- cos(al)
sa <- sin(al)
xy %*% rbind(c(ca, sa), c(-sa, ca)) + cbind(rep(loc[1],n),
rep(loc[2],n))
}

## Example:

ep <- ellipsePoints(2,5)
str(ep)
plot(ep, type="n",asp=1) ; polygon(ep, col = 2)

## Movie : rotating ellipse  :
nTurns <- 4 # #{full 360 deg turns}
for(al in 1:(nTurns*360)) {
ep <- ellipsePoints(3,6, alpha=al, loc = c(5,2))
plot(ep,type="l",xlim=c(-1,11),ylim=c(-4,8),
asp=1, axes = FALSE, xlab="", ylab="")
}

## Movie : rotating _filled_ ellipse  :
for(al in 1:180) {
ep <- ellipsePoints(3,6, alpha=al, loc = c(5,2))
plot(ep,type="n",xlim=c(-1,11),ylim=c(-4,8),
asp=1, axes = FALSE, xlab="", ylab="")
polygon(ep,col=2,border=3,lwd=2.5)
}
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._

```