[R] Draw ellipses in S-PLUS or R?

[Martin Maechler]

>>>>> "Paul" == Paul Y Peng <ypeng at math.mun.ca>
>>>>>     on Tue, 22 Oct 2002 16:59:52 -0230 writes:

    Paul> Dear S-PLUS/R users:
    Paul> Do you know any default function or a user contributed function that
    Paul> can draw an ellipse with given axes and origin? Thanks for any help.

Yes, here (written for R, should work in S-plus as well;
	   for the examples, and S-plus, you need to replace
	   col = "gray" by col = <number>}


### Copyright 2002 (C) Martin Maechler <maechler at stat.math.ethz.ch>

## Martin's Solution : radially equispaced ("radial gleichmässig") :

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)
## rotate by 30 degrees :
plot(ellipsePoints(2,5, alpha = 30), asp=1)
abline(h=0,v=0,col="gray")
abline(a=0,b= tan( 30 *pi/180), col=2, lty = 2)
abline(a=0,b= tan(120 *pi/180), col=3, lty = 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 {less nice to look at}
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)
}

###-------------------------------------------------



For R, there's also the  "ellipse" package on CRAN
which does other ellipses (confidence ...).
i.e.
	install.package("ellipse")
	library(ellipse)
if you're online.

Regards,
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			<><
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._