R/scales.R

   1 ## scales.R (2012-12-19)
   2
   3 ##   Add a Scale Bar or Axis to a Phylogeny Plot
   4
   5 ## add.scale.bar: add a scale bar to a phylogeny plot
   6 ## axisPhylo: add a scale axis on the side of a phylogeny plot
   7
   8 ## Copyright 2002-2011 Emmanuel Paradis
   9
  10 ## This file is part of the R-package `ape'.
  11 ## See the file ../COPYING for licensing issues.
  12
  13 add.scale.bar <- function(x, y, length = NULL, ask = FALSE,
  14                           lwd = 1, lcol = "black", ...)
  15 {
  16     lastPP <- get("last_plot.phylo", envir = .PlotPhyloEnv)
  17     direc <- lastPP$direction
  18     if (is.null(length)) {
  19         nb.digit <-
  20           if (direc %in% c("rightwards", "leftwards")) diff(range(lastPP$xx))
  21           else diff(range(lastPP$yy))
  22         nb.digit <- ceiling(log10(nb.digit)) - 2
  23         length <- eval(parse(text = paste("1e", nb.digit, sep = "")))
  24     }
  25
  26     if (ask) {
  27         cat("\nClick where you want to draw the bar\n")
  28         x <- unlist(locator(1))
  29         y <- x[2]
  30         x <- x[1]
  31     } else if (missing(x) || missing(y)) {
  32         if (lastPP$type %in% c("phylogram", "cladogram")) {
  33             switch(direc,
  34                    "rightwards" = {
  35                        x <- 0
  36                        y <- 1
  37                    },
  38                    "leftwards" = {
  39                        x <- max(lastPP$xx)
  40                        y <- 1
  41                    },
  42                    "upwards" = {
  43                        x <- max(lastPP$xx)
  44                        y <- 0
  45                    },
  46                    "downwards" = {
  47                        x <- 1
  48                        y <- max(lastPP$yy)
  49                    })
  50         } else {
  51             direc <- "rightwards" # just to be sure for below
  52             x <- lastPP$x.lim[1]
  53             y <- lastPP$y.lim[1]
  54         }
  55     }
  56
  57     switch(direc,
  58            "rightwards" = {
  59                segments(x, y, x + length, y, col = lcol, lwd = lwd)
  60                text(x + length * 1.1, y, as.character(length), adj = c(0, 0.5), ...)
  61            },
  62            "leftwards" = {
  63                segments(x - length, y, x, y, col = lcol, lwd = lwd)
  64                text(x - length * 1.1, y, as.character(length), adj = c(1, 0.5), ...)
  65            },
  66            "upwards" = {
  67                segments(x, y, x, y + length, col = lcol, lwd = lwd)
  68                text(x, y + length * 1.1, as.character(length), adj = c(0, 0.5), srt = 90, ...)
  69            },
  70            "downwards" = {
  71                segments(x, y - length, x, y, col = lcol, lwd = lwd)
  72                text(x, y - length * 1.1, as.character(length), adj = c(0, 0.5), srt = 270, ...)
  73            })
  74 }
  75
  76 axisPhylo <- function(side = 1, ...)
  77 {
  78     lastPP <- get("last_plot.phylo", envir = .PlotPhyloEnv)
  79     type <- lastPP$type
  80
  81     if (type == "unrooted")
  82         stop("axisPhylo() not available for unrooted plots; try add.scale.bar()")
  83     if (type == "radial")
  84         stop("axisPhylo() not meaningful for this type of plot")
  85
  86     if (type %in% c("phylogram", "cladogram")) {
  87         if (lastPP$direction %in% c("rightwards", "leftwards")) {
  88             x <- pretty(lastPP$xx)
  89             if (lastPP$direction == "rightwards") maxi <- max(lastPP$xx)
  90             else {
  91                 maxi <- min(lastPP$xx)
  92                 x <- -x
  93             }
  94         } else {
  95             x <- pretty(lastPP$yy)
  96             if (lastPP$direction == "upwards") maxi <- max(lastPP$yy)
  97             else {
  98                 maxi <- min(lastPP$yy)
  99                 x <- -x
 100             }
 101         }
 102         axis(side = side, at = c(maxi - x), labels = abs(x), ...)
 103     } else { # type == "fan"
 104         n <- lastPP$Ntip
 105         xx <- lastPP$xx[1:n]; yy <- lastPP$yy[1:n]
 106         r0 <- max(sqrt(xx^2 + yy^2))
 107         firstandlast <- c(1, n)
 108         theta0 <- mean(atan2(yy[firstandlast], xx[firstandlast]))
 109         x0 <- r0 * cos(theta0); y0 <- r0 * sin(theta0)
 110         inc <- diff(pretty(c(0, r0))[1:2])
 111         srt <- 360*theta0/(2*pi)
 112         coef <- -1
 113         if (abs(srt) > 90) {
 114             srt <- srt + 180
 115             coef <- 1
 116         }
 117         len <- 0.025 * r0 # the length of tick marks
 118         r <- r0
 119         while (r > 1e-8) {
 120             x <- r * cos(theta0); y <- r * sin(theta0)
 121             if (len/r < 1) {
 122                 ra <- sqrt(len^2 + r^2); thetaa <- theta0 + coef * asin(len/r)
 123                 xa <- ra * cos(thetaa); ya <- ra * sin(thetaa)
 124                 segments(xa, ya, x, y)
 125                 text(xa, ya, r0 - r, srt = srt, adj = c(0.5, 1.1), ...)
 126             }
 127             r <- r - inc
 128         }
 129         segments(x, y, x0, y0)
 130     }
 131 }