[Phylobase-commits] r615 - in pkg: R src

[noreply at r-forge.r-project.org]

Author: regetz
Date: 2009-09-03 01:47:29 +0200 (Thu, 03 Sep 2009)
New Revision: 615

Added:
   pkg/src/descendants.c
Modified:
   pkg/R/treewalk.R
Log:
Reimplemented descendants() in C, leaving original R code commented out.


Modified: pkg/R/treewalk.R
===================================================================
--- pkg/R/treewalk.R	2009-09-02 22:52:59 UTC (rev 614)
+++ pkg/R/treewalk.R	2009-09-02 23:47:29 UTC (rev 615)
@@ -66,21 +66,54 @@
 {
     ## FIXME: allow vector of nodes? (or just let people lapply?)
     type <- match.arg(type)
-    if (type=="children") return(children(phy,node))
-    node <- getNode(phy,node)
-    if (is.na(node)) stop("node ",node," not found in tree")
+    if (type=="children") {
+        return(children(phy, node))
+    }
+    node <- getNode(phy, node)
+    if (is.na(node)) stop("node ", node, " not found in tree")
     n <- nTips(phy)
-    if (node <= n) return(node)
-    l <- numeric()
-    d <- children(phy, node)
-    for (j in d) {
-        if (j <= n)
-          l <- c(l,j)
-        else if (type=="all") l <- c(l,j,
-                   descendants(phy,j,type="all"))
-        else l <- c(l, descendants(phy,j,type=type))
+    if (node <= n) {
+        return(node)
     }
-    return(getNode(phy,l))
+
+    ## edge matrix must be in preorder for the C function!
+    if (phy at order=="preorder") {
+        edge <- phy at edge
+    } else {
+        edge <- reorder(phy, order="preorder")@edge
+    }
+    ## deal with NA root node
+    edge[is.na(edge)] <- 0
+    ## extract edge columns
+    ancestor <- as.integer(edge[, 1])
+    descendant <- as.integer(edge[, 2])
+    ## create indicator vector and seed it based on children
+    isChild <- rep(0L, length=nrow(edge))
+    isChild[ancestor==node] <- 1L
+    numEdges <- as.integer(length(isChild))
+
+    ## returned value of isChild will indicate *all* descendants 
+    isDescendant <- .C("descendants", isChild, ancestor, descendant,
+        numEdges)[[1]]
+
+    l <- descendant[as.logical(isDescendant)]
+    if (type=="tips") {
+        l <- l[l<=n]
+    }
+
+    ## Original pure R implementation of the above
+    ## (note that it does not require preorder ordering)
+    ##l <- numeric()
+    ##d <- children(phy, node)
+    ##for (j in d) {
+    ##    if (j <= n)
+    ##      l <- c(l,j)
+    ##    else if (type=="all") l <- c(l,j,
+    ##               descendants(phy,j,type="all"))
+    ##    else l <- c(l, descendants(phy,j,type=type))
+    ##}
+
+    return(getNode(phy, l))
 }
 
 siblings <- function(phy, node, include.self=FALSE) {

Added: pkg/src/descendants.c
===================================================================
--- pkg/src/descendants.c	                        (rev 0)
+++ pkg/src/descendants.c	2009-09-02 23:47:29 UTC (rev 615)
@@ -0,0 +1,24 @@
+/*
+  descendants.c:
+    Identify all descendants of a given node. Function inputs are
+  derived from a phylo4 edge matrix, which *must* be in preorder order.
+  The isDescendant input vector should contain 1 for the immediate
+  children of the node, and 0 otherwise. The function returns this
+  vector updated to include all further descendants.
+*/
+
+#include <R.h>
+
+void descendants(int *isDescendant, int *ancestor, int *descendant, 
+    int *numEdges) {
+
+    int child=0;
+    for (int i=0; i<*numEdges; i++) {
+        if (isDescendant[i]==1) {
+            child = descendant[i];
+            for (int j=i+1; j<*numEdges; j++) {
+                if (ancestor[j]==child) isDescendant[j]=1; 
+            }
+        }
+    }
+}


Property changes on: pkg/src/descendants.c
___________________________________________________________________
Name: svn:eol-style
   + native