]> git.donarmstrong.com Git - ape.git/blobdiff - R/plot.phylo.R
some updates for ape 3.0-7
[ape.git] / R / plot.phylo.R
index da5ebec9d1fb58e74b314da9f24dc6723dd46c28..b9309a004882249ce7444517e1eb1bdb1157e63f 100644 (file)
@@ -1,8 +1,8 @@
-## plot.phylo.R (2011-06-14)
+## plot.phylo.R (2013-01-11)
 
 ##   Plot Phylogenies
 
-## Copyright 2002-2011 Emmanuel Paradis
+## Copyright 2002-2013 Emmanuel Paradis
 
 ## This file is part of the R-package `ape'.
 ## See the file ../COPYING for licensing issues.
@@ -15,11 +15,12 @@ plot.phylo <-
              adj = NULL, srt = 0, no.margin = FALSE, root.edge = FALSE,
              label.offset = 0, underscore = FALSE, x.lim = NULL,
              y.lim = NULL, direction = "rightwards", lab4ut = "horizontal",
-             tip.color = "black", plot = TRUE, ...)
+             tip.color = "black", plot = TRUE, rotate.tree = 0,
+             open.angle = 0, ...)
 {
     Ntip <- length(x$tip.label)
-    if (Ntip == 1) {
-        warning("found only one tip in the tree")
+    if (Ntip < 2) {
+        warning("found less than 2 tips in the tree")
         return(NULL)
     }
     if (any(tabulate(x$edge[, 1]) == 1))
@@ -96,17 +97,7 @@ plot.phylo <-
     }
     ## 'z' is the tree in pruningwise order used in calls to .C
     z <- reorder(x, order = "pruningwise")
-###    edge.color <- rep(edge.color, length.out = Nedge)
-###    edge.width <- rep(edge.width, length.out = Nedge)
-###    edge.lty <- rep(edge.lty, length.out = Nedge)
-###    ## fix from Li-San Wang (2007-01-23):
-###    xe <- x$edge
-###    x <- reorder(x, order = "pruningwise")
-###    ereorder <- match(x$edge[, 2], xe[, 2])
-###    edge.color <- edge.color[ereorder]
-###    edge.width <- edge.width[ereorder]
-###    edge.lty <- edge.lty[ereorder]
-###    ## end of fix
+
     if (phyloORclado) {
         if (is.null(node.pos)) {
             node.pos <- 1
@@ -131,12 +122,16 @@ plot.phylo <-
         } else  {
             xx <- .nodeDepthEdgelength(Ntip, Nnode, z$edge, Nedge, z$edge.length)
         }
-    } else switch(type, "fan" = {
+    } else {
+    twopi <- 2 * pi
+    rotate.tree <- twopi * rotate.tree/360
+    switch(type, "fan" = {
         ## if the tips are not in the same order in tip.label
         ## and in edge[, 2], we must reorder the angles: we
         ## use `xx' to store temporarily the angles
         TIPS <- x$edge[which(x$edge[, 2] <= Ntip), 2]
-        xx <- seq(0, 2*pi*(1 - 1/Ntip), 2*pi/Ntip)
+        xx <- seq(0, twopi * (1 - 1/Ntip) - twopi * open.angle/360,
+                  length.out = Ntip)
         theta <- double(Ntip)
         theta[TIPS] <- xx
         theta <- c(theta, numeric(Nnode))
@@ -147,14 +142,15 @@ plot.phylo <-
             r <- .nodeDepth(Ntip, Nnode, z$edge, Nedge)
             r <- 1/r
         }
-        xx <- r*cos(theta)
-        yy <- r*sin(theta)
+        theta <- theta + rotate.tree
+        xx <- r * cos(theta)
+        yy <- r * sin(theta)
     }, "unrooted" = {
         nb.sp <- .nodeDepth(Ntip, Nnode, z$edge, Nedge)
         XY <- if (use.edge.length)
-            unrooted.xy(Ntip, Nnode, z$edge, z$edge.length, nb.sp)
+            unrooted.xy(Ntip, Nnode, z$edge, z$edge.length, nb.sp, rotate.tree)
         else
-            unrooted.xy(Ntip, Nnode, z$edge, rep(1, Nedge), nb.sp)
+            unrooted.xy(Ntip, Nnode, z$edge, rep(1, Nedge), nb.sp, rotate.tree)
         ## rescale so that we have only positive values
         xx <- XY$M[, 1] - min(XY$M[, 1])
         yy <- XY$M[, 2] - min(XY$M[, 2])
@@ -164,11 +160,11 @@ plot.phylo <-
         ## radius:
         X <- 1 - X/Ntip
         ## angle (1st compute the angles for the tips):
-        yy <- c((1:Ntip)*2*pi/Ntip, rep(0, Nnode))
+        yy <- c((1:Ntip)*twopi/Ntip, rep(0, Nnode))
         Y <- .nodeHeight(Ntip, Nnode, z$edge, Nedge, yy)
-        xx <- X * cos(Y)
-        yy <- X * sin(Y)
-    })
+        xx <- X * cos(Y + rotate.tree)
+        yy <- X * sin(Y + rotate.tree)
+    })}
     if (phyloORclado) {
         if (!horizontal) {
             tmp <- yy
@@ -267,12 +263,12 @@ plot.phylo <-
         y.lim <- c(0, y.lim)
         if (phyloORclado && horizontal) y.lim[1] <- 1
         if (type %in% c("fan", "unrooted") && show.tip.label)
-          y.lim[1] <- -max(nchar(x$tip.label) * 0.018 * max(yy) * cex)
+            y.lim[1] <- -max(nchar(x$tip.label) * 0.018 * max(yy) * cex)
         if (type == "radial")
-          y.lim[1] <- if (show.tip.label) -1 - max(nchar(x$tip.label) * 0.018 * max(yy) * cex) else -1
+            y.lim[1] <- if (show.tip.label) -1 - max(nchar(x$tip.label) * 0.018 * max(yy) * cex) else -1
     }
     ## mirror the yy:
-    if (phyloORclado && direction == "downwards") yy <- y.lim[2] - yy
+    if (phyloORclado && direction == "downwards") yy <- max(yy) - yy
     if (phyloORclado && root.edge) {
         if (direction == "leftwards") x.lim[2] <- x.lim[2] + x$root.edge
         if (direction == "downwards") y.lim[2] <- y.lim[2] + x$root.edge
@@ -382,7 +378,6 @@ if (plot) {
         if (type %in% c("fan", "radial")) {
             xx.tips <- xx[1:Ntip]
             yy.tips <- yy[1:Ntip]
-            ## using atan2 considerably facilitates things compared to acos...
             angle <- atan2(yy.tips, xx.tips) # in radians
             if (label.offset) {
                 xx.tips <- xx.tips + label.offset * cos(angle)
@@ -428,7 +423,7 @@ phylogram.plot <- function(edge, Ntip, Nnode, xx, yy, horizontal,
         yy <- xx
         xx <- tmp
     }
-    ## un trait vertical à chaque noeud...
+    ## un trait vertical a chaque noeud...
     x0v <- xx[nodes]
     y0v <- y1v <- numeric(Nnode)
     ## store the index of each node in the 1st column of edge:
@@ -558,7 +553,7 @@ circular.plot <- function(edge, Ntip, Nnode, xx, yy, theta,
     }
 }
 
-unrooted.xy <- function(Ntip, Nnode, edge, edge.length, nb.sp)
+unrooted.xy <- function(Ntip, Nnode, edge, edge.length, nb.sp, rotate.tree)
 {
     foo <- function(node, ANGLE, AXIS) {
         ind <- which(edge[, 1] == node)
@@ -573,7 +568,7 @@ unrooted.xy <- function(Ntip, Nnode, edge, edge.length, nb.sp)
             yy[sons[i]] <<- h*sin(beta) + yy[node]
         }
         for (i in sons)
-          if (i > Ntip) foo(i, angle[i], axis[i])
+            if (i > Ntip) foo(i, angle[i], axis[i])
     }
     Nedge <- dim(edge)[1]
     yy <- xx <- numeric(Ntip + Nnode)
@@ -581,7 +576,7 @@ unrooted.xy <- function(Ntip, Nnode, edge, edge.length, nb.sp)
     ## `axis': the axis of each branch
     axis <- angle <- numeric(Ntip + Nnode)
     ## start with the root...
-    foo(Ntip + 1L, 2*pi, 0)
+    foo(Ntip + 1L, 2*pi, 0 + rotate.tree)
 
     M <- cbind(xx, yy)
     axe <- axis[1:Ntip] # the axis of the terminal branches (for export)
@@ -602,14 +597,63 @@ node.depth <- function(phy)
        as.integer(N), double(n + m), DUP = FALSE, PACKAGE = "ape")[[6]]
 }
 
+node.depth.edgelength <- function(phy)
+{
+    n <- length(phy$tip.label)
+    m <- phy$Nnode
+    N <- dim(phy$edge)[1]
+    phy <- reorder(phy, order = "pruningwise")
+    .C("node_depth_edgelength", as.integer(n), as.integer(n),
+       as.integer(phy$edge[, 1]), as.integer(phy$edge[, 2]),
+       as.integer(N), as.double(phy$edge.length), double(n + m),
+       DUP = FALSE, PACKAGE = "ape")[[7]]
+}
+
+node.height <- function(phy)
+{
+    n <- length(phy$tip.label)
+    m <- phy$Nnode
+    N <- dim(phy$edge)[1]
+    phy <- reorder(phy, order = "pruningwise")
+
+    e1 <- phy$edge[, 1]
+    e2 <- phy$edge[, 2]
+
+    yy <- numeric(n + m)
+    TIPS <- e2[e2 <= n]
+    yy[TIPS] <- 1:n
+
+    .C("node_height", as.integer(n), as.integer(m),
+       as.integer(e1), as.integer(e2), as.integer(N),
+       as.double(yy), DUP = FALSE, PACKAGE = "ape")[[6]]
+}
+
+node.height.clado <- function(phy)
+{
+    n <- length(phy$tip.label)
+    m <- phy$Nnode
+    N <- dim(phy$edge)[1]
+    phy <- reorder(phy, order = "pruningwise")
+
+    e1 <- phy$edge[, 1]
+    e2 <- phy$edge[, 2]
+
+    yy <- numeric(n + m)
+    TIPS <- e2[e2 <= n]
+    yy[TIPS] <- 1:n
+
+    .C("node_height_clado", as.integer(n), as.integer(m),
+       as.integer(e1), as.integer(e2), as.integer(N),
+       double(n + m), as.double(yy), DUP = FALSE,
+       PACKAGE = "ape")[[7]]
+}
+
 plot.multiPhylo <- function(x, layout = 1, ...)
 {
-    if (layout > 1)
-      layout(matrix(1:layout, ceiling(sqrt(layout)), byrow = TRUE))
-    else layout(matrix(1))
-    if (!par("ask")) {
-        par(ask = TRUE)
-        on.exit(par(ask = FALSE))
+    layout(matrix(1:layout, ceiling(sqrt(layout)), byrow = TRUE))
+    if (!devAskNewPage() && interactive()) {
+        devAskNewPage(TRUE)
+        on.exit(devAskNewPage(FALSE))
     }
     for (i in 1:length(x)) plot(x[[i]], ...)
 }
@@ -675,9 +719,9 @@ kronoviz <- function(x, layout = length(x), horiz = TRUE, ...)
         h <- 1
     }
     layout(matrix(1:layout, nrow), widths = w, heights = h)
-    if (layout > Ntree && !par("ask")) {
-        par(ask = TRUE)
-        on.exit(par(ask = FALSE))
+    if (layout < Ntree && !devAskNewPage() && interactive()) {
+        devAskNewPage(TRUE)
+        on.exit(devAskNewPage(FALSE))
     }
     if (horiz) {
         for (i in 1:Ntree)