1## Werner, Martin and Tim have added several useful
2## things. At the end of this e-mail there is our final result.
3##
4## As an example we reproduced a similar figure as Fig. 4.23 of Chambers et
5## al. (1983) "Graphical Methods For Data Analysis":
6##
7## ii_3:4
8## x <- matrix(aperm(iris[,ii,], perm =c(1,3,2)), ncol=2,
9##             dimnames=list(dimnames(iris)[[1]],dimnames(iris)[[2]][ii]))
10## xr <- round(2*x,1)/2
11## nam <- dimnames(xr)[[2]]
12## p.sunflowers(xr[,1],xr[,2], xlab=nam[1], ylab=nam[2], size= 1/16,
13##              main="Iris data")
14##
15##
16## Andreas Ruckstuhl <ruckstuhl@stat.math.ethz.ch>
17## Seminar fuer Statistik, SOL G5, ETH (Federal Institute of Technology)
18## 8092 Zurich	SWITZERLAND  	phone: x-41-1-256-5319  fax: x-41-1-252-3410
19##
20##
21##================================ S function ========================
22##
23
24if(!.R.) {
25  p.sunflowers <- function(x, y, number, size = 0.125, add = FALSE,
26                           pch = 16, ...)
27  {
28    ## Purpose: Produce a 'sunflower'-Plot
29    ## -------------------------------------------------------------------------
30    ## Arguments: x,y: coordinates;
31    ##    number[i] = number of times for (x[i],y[i])  [may be 0]
32    ##    size: in inches;  1 in := 2.54 cm
33    ##    add : (logical) Should I add to a previous plot ?
34    ##    further args: as for plot(..)
35    ## -------------------------------------------------------------------------
36    ## Authors: Andreas Ruckstuhl, Werner Stahel, Martin Maechler, Tim Hesterberg
37    ## Date   : Aug 89 / Jan 93,   March 92,      Jan 93,          Jan 93
38    ## Examples: p.sunflowers(x=sort(round(rnorm(100))), y= round(2*rnorm(100),0))
39    ## ~~~~~~~~  p.sunflowers(rnorm(100),rnorm(100), number=rpois(n=100,lambda=2),
40    ##                        main="Sunflower plot")
41
42    n <- length(x)
43    if(length(y) != n)
44      stop("x & y must have same length !")
45
46    if(missing(number)) {
47      orderxy <- order(x, y)
48      x <- x[orderxy]
49      y <- y[orderxy]
50      first <- c(TRUE, (x[-1] != x[ - n]) | (y[-1] != y[ - n]))
51      x <- x[first]
52      y <- y[first]
53      number <- diff(c((1:n)[first], n + 1))
54    } else {
55      if(length(number) != n)
56        stop("number must have same length as x & y !")
57
58      x <- x[number > 0]
59      y <- y[number > 0]
60      number <- number[number > 0]
61    }
62
63    n <- length(x)
64    if(!add) {
65      axislabels <- match(c("xlab", "ylab"), names(list(...)))
66      if(!is.na(axislabels[1]))
67        xlab <- list(...)[[axislabels[1]]]
68      else xlab <- deparse(substitute(x))
69
70      if(!is.na(axislabels[2]))
71        ylab <- list(...)[[axislabels[2]]]
72      else ylab <- deparse(substitute(y))
73
74      plot(x, y, xlab = xlab, ylab = ylab, type = "n", ...)
75    }
76
77    nequ1 <- number == 1
78    if(any(nequ1))
79      points(x[nequ1], y[nequ1], pch = pch, csi = size * 1.25)
80
81    if(any(!nequ1))
82      points(x[!nequ1], y[!nequ1], pch = pch, csi = size * 0.8)
83
84    i.multi <- (1:n)[number > 1]
85    if(length(i.multi)) {
86      ppin <- par()$pin
87      pusr <- par()$usr
88      xr <- (size * abs(pusr[2] - pusr[1]))/ppin[1]
89      yr <- (size * abs(pusr[4] - pusr[3]))/ppin[2]
90      i.rep <- rep(i.multi, number[number > 1])
91      z <- NULL
92      for(i in i.multi)
93        z <- c(z, 1:number[i])
94
95      deg <- (2 * pi * z)/number[i.rep]
96      segments(x[i.rep], y[i.rep], x[i.rep] + xr * sin(deg), y[i.rep] +
97               yr * cos(deg))
98    }
99
100    invisible()
101  }
102
103  NULL
104}
105