diff --git a/.Rbuildignore b/.Rbuildignore
index e5e9728e..c6eacc40 100644
--- a/.Rbuildignore
+++ b/.Rbuildignore
@@ -11,4 +11,9 @@ mkdocs
 docs
 logo.png.github
 .github
+<<<<<<< HEAD
+^.*\.Rproj$
+^\.Rproj\.user$
+=======
 ggtree_sticker.R
+>>>>>>> 6c3d0a213ebba123a58af768ef6d7b4ec93b85d1
diff --git a/.gitignore b/.gitignore
index 2a4b674c..6a7c648e 100644
--- a/.gitignore
+++ b/.gitignore
@@ -8,3 +8,6 @@ __init__.py
 __pycache__
 __init__.pyc
 
+.Rproj.user
+ggtree.Rproj
+.RData
diff --git a/R/ggtree.R b/R/ggtree.R
index 5844bb2c..6c9f855b 100644
--- a/R/ggtree.R
+++ b/R/ggtree.R
@@ -5,6 +5,7 @@
 ##' @param tr phylo object
 ##' @param mapping aes mapping
 ##' @param layout one of 'rectangular', 'slanted', 'fan', 'circular', 'radial' or 'unrooted'
+##' @param layout.method of 'equal_angle', 'daylight'.
 ##' @param open.angle open angle, only for 'fan' layout
 ##' @param mrsd most recent sampling date
 ##' @param as.Date logical whether using Date class in time tree
@@ -33,6 +34,7 @@
 ggtree <- function(tr,
                    mapping        = NULL,
                    layout         = "rectangular",
+                   layout.method  = "equal_angle",
                    open.angle     = 0,
                    mrsd           = NULL,
                    as.Date        = FALSE,
@@ -44,8 +46,10 @@ ggtree <- function(tr,
                    ndigits        = NULL,
                    ...) {
 
+    # Check if layout string is valid.
     layout %<>% match.arg(c("rectangular", "slanted", "fan", "circular", "radial", "unrooted"))
-
+    layout.method %<>% match.arg(c("equal_angle", "daylight"))
+  
     if (is(tr, "r8s") && branch.length == "branch.length") {
         branch.length = "TREE"
     }
@@ -60,8 +64,11 @@ ggtree <- function(tr,
     } else {
         mapping <- modifyList(aes_(~x, ~y), mapping)
     }
-    p <- ggplot(tr, mapping=mapping,
+  
+    p <- ggplot(tr, 
+                mapping       = mapping,
                 layout        = layout,
+                layout.method = layout.method,
                 mrsd          = mrsd,
                 as.Date       = as.Date,
                 yscale        = yscale,
diff --git a/R/method-fortify.R b/R/method-fortify.R
index 3230e79e..96b31cdf 100644
--- a/R/method-fortify.R
+++ b/R/method-fortify.R
@@ -110,6 +110,7 @@ rm.singleton.newick <- function(nwk, outfile = NULL) {
 ##' @export
 fortify.beast <- function(model, data,
                           layout        = "rectangular",
+                          layout.method = "equal_angle",
                           yscale        = "none",
                           ladderize     = TRUE,
                           right         = FALSE,
@@ -119,7 +120,8 @@ fortify.beast <- function(model, data,
 
     model <- set_branch_length(model, branch.length)
     phylo <- model@phylo
-    df    <- fortify(phylo, layout=layout, branch.length=branch.length,
+    df    <- fortify(phylo, layout=layout, layout.method=layout.method,
+                     branch.length=branch.length,
                      ladderize=ladderize, right=right, mrsd = mrsd, ...)
 
     stats <- model@stats
@@ -220,6 +222,7 @@ fortify.beast <- function(model, data,
 ##' @export
 fortify.codeml <- function(model, data,
                            layout        = "rectangular",
+                           layout.method = "equal_angle",
                            yscale        = "none",
                            ladderize     = TRUE,
                            right         = FALSE,
@@ -265,6 +268,7 @@ fortify.codeml <- function(model, data,
 ##' @export
 fortify.codeml_mlc <- function(model, data,
                                layout        = "rectangular",
+                               layout.method = "equal_angle",
                                yscale        = "none",
                                ladderize     = TRUE,
                                right         = FALSE,
@@ -306,6 +310,7 @@ merge_phylo_anno.codeml_mlc <- function(df, dNdS, ndigits = NULL) {
 
 fortify.codeml_mlc_ <- function(model, data,
                                 layout        = "rectangular",
+                                layout.method = "equal_angle",
                                 ladderize     = TRUE,
                                 right         = FALSE,
                                 branch.length = "branch.length",
@@ -317,8 +322,12 @@ fortify.codeml_mlc_ <- function(model, data,
 
 ##' @method fortify paml_rst
 ##' @export
-fortify.paml_rst <- function(model, data, layout = "rectangular", yscale="none",
-                             ladderize=TRUE, right=FALSE, mrsd=NULL, ...) {
+fortify.paml_rst <- function(model, data, 
+                             layout = "rectangular", 
+                             yscale="none",
+                             ladderize=TRUE, 
+                             right=FALSE, 
+                             mrsd=NULL, ...) {
     df <- fortify.phylo(model@phylo, data, layout, ladderize, right, mrsd=mrsd, ...)
     df <- merge_phylo_anno.paml_rst(df, model)
     df <- scaleY(model@phylo, df, yscale, layout, ...)
@@ -353,8 +362,11 @@ fortify.hyphy <- fortify.paml_rst
 ##' @importFrom treeio get.placements
 ##' @export
 fortify.jplace <- function(model, data,
-                           layout="rectangular", yscale="none",
-                           ladderize=TRUE, right=FALSE, mrsd=NULL, ...) {
+                           layout="rectangular",
+                           yscale="none",
+                           ladderize=TRUE, 
+                           right=FALSE, 
+                           mrsd=NULL, ...) {
     df <- extract.treeinfo.jplace(model, layout, ladderize, right, mrsd=mrsd, ...)
     place <- get.placements(model, by="best")
 
@@ -403,8 +415,13 @@ fortify.phylo4 <- function(model, data, layout="rectangular", yscale="none",
 
 ##' @method fortify phylo4d
 ##' @export
-fortify.phylo4d <- function(model, data, layout="rectangular", yscale="none",
-                            ladderize=TRUE, right=FALSE, branch.length="branch.length",
+fortify.phylo4d <- function(model, data, 
+                            layout="rectangular",
+                            layout.method = "equal_angle",
+                            yscale="none",
+                            ladderize=TRUE, 
+                            right=FALSE, 
+                            branch.length="branch.length",
                             mrsd=NULL, ...) {
     ## model <- set_branch_length(model, branch.length)
     ## phylo <- as.phylo.phylo4(model)
@@ -413,7 +430,7 @@ fortify.phylo4d <- function(model, data, layout="rectangular", yscale="none",
     ## tdata <- model@data[match(res$node, rownames(model@data)), , drop=FALSE]
     ## df <- cbind(res, tdata)
     ## scaleY(as.phylo.phylo4(model), df, yscale, layout, ...)
-    fortify(as.treedata(model), data, layout, yscale, ladderize, right, branch.length, mrsd, ...)
+    fortify(as.treedata(model), data, layout, layout.method, yscale, ladderize, right, branch.length, mrsd, ...)
 }
 
 
@@ -438,8 +455,12 @@ fortify.phylo4d <- function(model, data, layout="rectangular", yscale="none",
 ##' @method fortify phylo
 ##' @export
 ##' @author Yu Guangchuang
-fortify.phylo <- function(model, data, layout="rectangular",
-                          ladderize=TRUE, right=FALSE, mrsd=NULL, as.Date=FALSE, ...) {
+fortify.phylo <- function(model, data, 
+                          layout="rectangular",
+                          ladderize=TRUE, 
+                          right=FALSE, 
+                          mrsd=NULL, 
+                          as.Date=FALSE, ...) {
     tree <- reorder.phylo(model, 'postorder')
 
     if (ladderize == TRUE) {
@@ -605,18 +626,26 @@ fortify.multiPhylo <-  function(model, data, layout="rectangular",
 
 ##' @method fortify phylip
 ##' @export
-fortify.phylip <- function(model, data, layout="rectangular",
-                        ladderize=TRUE, right=FALSE,
-                        branch.length = "TREE", mrsd=NULL, ...) {
+fortify.phylip <- function(model, data, 
+                           layout="rectangular",
+                           layout.method = "equal_angle",
+                           ladderize=TRUE, 
+                           right=FALSE,
+                           branch.length = "TREE", 
+                           mrsd=NULL, ...) {
     trees <- get.tree(model)
     fortify(trees, layout=layout, ladderize = ladderize, right=right, mrsd=mrsd, ...)
 }
 
 ##' @method fortify r8s
 ##' @export
-fortify.r8s <- function(model, data, layout="rectangular",
-                        ladderize=TRUE, right=FALSE,
-                        branch.length = "TREE", mrsd=NULL, ...) {
+fortify.r8s <- function(model, data, 
+                        layout="rectangular",
+                        layout.method = "equal_angle",
+                        ladderize=TRUE, 
+                        right=FALSE,
+                        branch.length = "TREE", 
+                        mrsd=NULL, ...) {
     trees <- get.tree(model)
     branch.length %<>% match.arg(names(trees))
     phylo <- trees[[branch.length]]
@@ -625,8 +654,11 @@ fortify.r8s <- function(model, data, layout="rectangular",
 
 ##' @method fortify obkData
 ##' @export
-fortify.obkData <- function(model, data, layout="rectangular",
-                            ladderize=TRUE, right=FALSE, mrsd = NULL, ...) {
+fortify.obkData <- function(model, data, 
+                            layout="rectangular",
+                            ladderize=TRUE, 
+                            right=FALSE, 
+                            mrsd = NULL, ...) {
 
     df <- fortify(model@trees[[1]], layout=layout, ladderize=ladderize, right=right, mrsd=mrsd, ...)
 
@@ -644,8 +676,11 @@ fortify.obkData <- function(model, data, layout="rectangular",
 
 ##' @method fortify phyloseq
 ##' @export
-fortify.phyloseq <- function(model, data, layout="rectangular",
-                             ladderize=TRUE, right=FALSE, mrsd=NULL, ...) {
+fortify.phyloseq <- function(model, data, 
+                             layout="rectangular",
+                             ladderize=TRUE, 
+                             right=FALSE, 
+                             mrsd=NULL, ...) {
 
     df <- fortify(model@phy_tree, layout=layout, ladderize=ladderize, right=right, mrsd=mrsd, ...)
     phyloseq <- "phyloseq"
diff --git a/R/tidytree.R b/R/tidytree.R
index 7255b2c2..f8d9f977 100644
--- a/R/tidytree.R
+++ b/R/tidytree.R
@@ -159,79 +159,625 @@ reroot_node_mapping <- function(tree, tree2) {
 
 
 ##' @importFrom ape reorder.phylo
-layout.unrooted <- function(tree, branch.length="branch.length", ...) {
-    N <- getNodeNum(tree)
-    root <- getRoot(tree)
+layout.unrooted <- function(tree, branch.length="branch.length", layout.method="equal.angle", ...) {
 
-    df <- as.data.frame.phylo_(tree)
-    df$x <- NA
-    df$y <- NA
-    df$start <- NA
-    df$end   <- NA
-    df$angle <- NA
-    df[root, "x"] <- 0
-    df[root, "y"] <- 0
-    df[root, "start"] <- 0
-    df[root, "end"]   <- 2
-    df[root, "angle"] <- 0
+    switch(layout.method,
+           equal_angle = { df <- layoutEqualAngle(tree, branch.length) },
+           daylight = { df <- layoutDaylight(tree, branch.length) }
+    )
 
-    nb.sp <- sapply(1:N, function(i) length(get.offspring.tip(tree, i)))
+    return(df)
+}
 
-    nodes <- getNodes_by_postorder(tree)
 
-    for(curNode in nodes) {
-        curNtip <- nb.sp[curNode]
-        children <- getChild(tree, curNode)
+##' 'Equal-angle layout algorithm for unrooted trees'
+##'
+##' @references 
+##' "Inferring Phylogenies" by Joseph Felsenstein.
+##' 
+##' @title layoutEqualAngle
+##' @param tree phylo object
+##' @param branch.length set to 'none' for edge length of 1. Otherwise the phylogenetic tree edge length is used.
+##' @return tree as data.frame with equal angle layout.
+layoutEqualAngle <- function(tree, branch.length ){
+  root <- getRoot(tree)
+  # Convert Phylo tree to data.frame.
+  df <- as.data.frame.phylo_(tree)
+  
+  # NOTE: Angles (start, end, angle) are in half-rotation units (radians/pi or degrees/180)
+  
+  # create and assign NA to the following fields.
+  df$x <- NA
+  df$y <- NA
+  df$start <- NA # Start angle of segment of subtree.
+  df$end   <- NA # End angle of segment of subtree
+  df$angle <- NA # Orthogonal angle to beta ... for labels??
+  # Initialize root node position and angles.
+  df[root, "x"] <- 0
+  df[root, "y"] <- 0
+  df[root, "start"] <- 0 # 0-degrees
+  df[root, "end"]   <- 2 # 360-degrees
+  df[root, "angle"] <- 0 # Angle label.
+  
+  N <- getNodeNum(tree)
+
+  # Get number of tips for each node in tree.
+  nb.sp <- sapply(1:N, function(i) length(get.offspring.tip(tree, i)))
+  # Get list of node id's.
+  nodes <- getNodes_by_postorder(tree)
+  
+  for(curNode in nodes) { 
+    # Get number of tips for current node.
+    curNtip <- nb.sp[curNode]
+    # Get array of child node indexes of current node.
+    children <- getChild(tree, curNode)
+    
+    # Get "start" and "end" angles of a segment for current node in the data.frame.
+    start <- df[curNode, "start"]
+    end <- df[curNode, "end"]
+    
+    if (length(children) == 0) {
+      ## is a tip
+      next
+    }
+    
+    for (i in seq_along(children)) {
+      child <- children[i]
+      # Get the number of tips for child node.
+      ntip.child <- nb.sp[child] 
+      
+      # Calculated in half radians.
+      # alpha: angle of segment for i-th child with ntips_ij tips.
+      # alpha = (left_angle - right_angle) * (ntips_ij)/(ntips_current)
+      alpha <- (end - start) * ntip.child / curNtip
+      # beta = angle of line from parent node to i-th child.
+      beta <- start + alpha / 2
+      
+      if (branch.length == "none") {
+        length.child <- 1
+      } else {
+        length.child <- df[child, "length"]
+      }
+      
+      # update geometry of data.frame.
+      # Calculate (x,y) position of the i-th child node from current node.
+      df[child, "x"] <- df[curNode, "x"] + cospi(beta) * length.child
+      df[child, "y"] <- df[curNode, "y"] + sinpi(beta) * length.child
+      # Calculate orthogonal angle to beta. 
+      df[child, "angle"] <- -90 - 180 * beta * sign(beta - 1)
+      # Update the start and end angles of the childs segment.
+      df[child, "start"] <- start
+      df[child, "end"] <- start + alpha
+      start <- start + alpha
+    } 
+
+  }
+  
+  return(df)
+  
+}
 
-        start <- df[curNode, "start"]
-        end <- df[curNode, "end"]
+##' Equal daylight layout method for unrooted trees.
+##' 
+##' #' @title
+##' @param tree phylo object
+##' @param branch.length set to 'none' for edge length of 1. Otherwise the phylogenetic tree edge length is used.
+##' @return tree as data.frame with equal angle layout.
+##' @references 
+##' The following aglorithm aims to implement the vague description of the "Equal-daylight Algorithm"
+##' in "Inferring Phylogenies" pp 582-584 by Joseph Felsenstein.
+##' 
+##' ```
+##' Leafs are subtrees with no children
+##' Initialise tree using equal angle algorithm
+##' tree_df = equal_angle(tree)
+##'
+##' nodes = get list of nodes in tree_df breadth-first
+##' nodes = remove tip nodes.
+##' 
+##' ```
+layoutDaylight <- function( tree, branch.length ){
+
+  # How to set optimal
+  MAX_COUNT <- 100
+  MINIMUM_AVERAGE_ANGLE_CHANGE <- 0.01
+  
+
+  # Initialize tree.
+  tree_df <- layoutEqualAngle(tree, branch.length)
+
+  # nodes = get list of nodes in tree_df
+  # Get list of node id's.
+  #nodes <- getNodes_by_postorder(tree)
+  # nodes <- GetSubtree.df(tree_df, root)
+  
+  # Get list of internal nodes
+  #nodes <- tree_df[tree_df$IsTip != TRUE]$nodes
+  
+  nodes <- getNodesBreadthFirst.df(tree_df)
+  # select only internal nodes
+  internal_nodes <- tree_df[!tree_df$isTip,]$node
+  # Remove tips from nodes list, but keeping order.
+  nodes <- intersect(nodes, internal_nodes)
+  
+  i <- 1
+  ave_change <- 1.0
+  while( i <= MAX_COUNT & ave_change > MINIMUM_AVERAGE_ANGLE_CHANGE ){
+    cat('Iteration: ', i, '\n')
+
+    # Reset max_change after iterating over tree.
+    total_max <- 0.0
+    
+    # for node in nodes {
+    for( j in seq_along(nodes)){
+      currentNode_id <- nodes[j]
+      
+      result <- applyLayoutDaylight(tree_df, currentNode_id)
+      tree_df <- result$tree
+      total_max <- total_max + result$max_change
+      
+    }
+    
+    ave_change <- total_max / length(nodes)
+    
+    cat('Average angle change', ave_change,'\n')
+    
+    i <- i + 1
+  }
+  
+  return(tree_df)
+  
+}
 
-        if (length(children) == 0) {
-            ## is a tip
-            next
+##' Apply the daylight alorithm to adjust the spacing between the subtrees and tips of the
+##' specified node.
+##' 
+##' @title applyLayoutDaylight
+##' @param df tree data.frame
+##' @param node_id is id of the node from which daylight is measured to the other subtrees.
+##' @return list with tree data.frame with updated layout using daylight algorithm and max_change angle .
+##' 
+##'
+##' ```
+##' for node in nodes {
+##'   if node is a leaf {
+##'     next
+##'   }
+##'   
+##'   subtrees = get subtrees of node
+##'   
+##'   for i-th subtree in subtrees {
+##'     [end, start] = get left and right angles of tree from node id.
+##'     angle_list[i, 'left'] = end
+##'     angle_list[i, 'beta'] = start - end  # subtree arc angle
+##'     angle_list[i, 'index'] = i-th # index of subtree/leaf 
+##'   }
+##'
+##'   sort angle_list by 'left' column in ascending order.
+##'
+##'   D = 360 - sum( angle_list['beta'] ) # total daylight angle
+##'   d = D / |subtrees| # equal daylight angle.
+##' 
+##'   new_L = left angle of first subtree.
+##'
+##'   for n-th row in angle_list{
+##'     # Calculate angle to rotate subtree/leaf to create correct daylight angle.
+##'     new_left_angle = new_left_angle + d + angle_list[n, 'beta']
+##'     Calculate the difference between the old and new left angles.
+##'     adjust_angle = new_left_angle - angle_list[n, 'left'] 
+##'     
+##'     index = angle_list['index']
+##'     rotate subtree[index] wrt n-th node by adjust_angle
+##'     }
+##'   }
+##' }
+##' ```
+applyLayoutDaylight <- function(df, node_id ){
+  
+  max_change <- 0.0
+  
+  # Get lists of node ids for each subtree, including  rest of unrooted tree.
+  subtrees <- getSubtreeUnrooted.df(df, node_id)
+  angle_list <- data.frame(left=numeric(0), beta=numeric(0), subtree_id=integer(0) )
+
+  # Return tree if only 2 or less subtrees to adjust.
+  if(length(subtrees) <= 2){
+    return( list(tree = df, max_change = max_change) )
+  }
+  
+  # Find start and end angles for each subtree.
+  #   subtrees = get subtrees of node
+  #   for i-th subtree in subtrees {
+  for (i in seq_along(subtrees) ) {
+    subtree <- subtrees[[i]]
+    # [end, start] = get start and end angles of tree.
+    
+    angles <- getTreeArcAngles(df, node_id, subtree)
+    angle_list[ i, 'subtree_id'] <- i
+    angle_list[ i, 'left'] <- angles['left']
+    angle_list[ i, 'beta'] <- angles['left'] - angles['right'] # subtree arc angle
+    # If subtree arc angle is -ve, then + 2 (360).
+    if(angle_list[ i, 'beta'] < 0 ){
+      angle_list[ i, 'beta'] <- angle_list[ i, 'beta'] + 2
+    } 
+  }
+  #   sort angle_list by 'left angle' column in ascending order.
+  angle_list <- angle_list[with(angle_list, order(left)), ]
+  #   D = 360 - sum( angle_list['beta'] ) # total day
+  #   d = D / |subtrees| # equal daylight angle.
+  total_daylight <- 2 - colSums(angle_list['beta']) 
+  d <- total_daylight / length(subtrees)
+   
+  # Initialise new left-angle as first subtree left-angle.
+  new_left_angle <- angle_list[1, 'left']
+
+  # Adjust angles of subtrees and tips connected to current node.
+  # for n-th row in angle_list{
+  # Skip the first subtree as it is not adjusted.
+  for (i in 2:nrow(angle_list) ) {
+    # Calculate angle to rotate subtree/leaf to create correct daylight angle.
+    new_left_angle <- new_left_angle + d + angle_list[i, 'beta']
+    # Calculate the difference between the old and new left angles.
+    adjust_angle <- new_left_angle - angle_list[i, 'left']
+    
+    max_change <- max(max_change, abs(adjust_angle))
+    #cat('Adjust angle:', abs(adjust_angle), ' Max change:', max_change ,'\n')
+    
+    # rotate subtree[index] wrt current node
+    subtree_id <- angle_list[i, 'subtree_id']
+    subtree_nodes <- subtrees[[subtree_id]]$subtree
+    # update tree_df for all subtrees with rotated points.
+    df <- rotateTreePoints.df(df, node_id, subtree_nodes, adjust_angle)
+  }
+  
+  return( list(tree = df, max_change = max_change) ) 
+  
+}
+
+
+##' Find the right (clockwise rotation, angle from +ve x-axis to furthest subtree nodes) and 
+##' left (anti-clockwise angle from +ve x-axis to subtree)
+##' 
+##' @title getTreeArcAngles
+##' @param df tree data.frame
+##' @param origin_id node id from which to calculate left and right hand angles of subtree.
+##' @param subtree named list of root id of subtree and list of node ids for given subtree.
+##' @return named list with right and left angles in range [0,2] i.e 1 = 180 degrees, 1.5 = 270 degrees.
+getTreeArcAngles <- function(df, origin_id, subtree) {
+  # Initialise variables
+  theta_child <- 0.0
+  subtree_root_id <- subtree$node
+  subtree_node_ids <- subtree$subtree
+  
+  # Initialise angle from origin node to parent node.
+  # If subtree_root_id is child of origin_id
+  if( any(subtree_root_id == getChild.df(df, origin_id)) ){
+    theta_parent <- getNodeAngle.df(df, origin_id, subtree_root_id)
+  }else{
+    # get the real root of df tree to initialise left and right angles.
+    theta_parent <- getNodeAngle.df(df, origin_id, getRoot.df(df))
+  }
+    
+  # create vector with named columns
+  # left-hand and right-hand angles between origin node and the extremities of the tree nodes.
+  arc <- c('left' = theta_parent, 'right' = theta_parent)
+    
+  # Subtree has to have 1 or more nodes to compare.
+  if (length(subtree_node_ids) == 0 ){
+    return(0)
+  } 
+
+  
+  # Remove tips from nodes list, but keeping order.
+  # internal_nodes <- df[!df$isTip,]$node
+  # subtree_node_ids <- intersect(subtree_node_ids, internal_nodes)
+  
+  
+  # Calculate the angle from the origin node to each child node.
+  # Moving from parent to children in depth-first traversal.
+  for( i in seq_along(subtree_node_ids) ){
+    parent_id <- subtree_node_ids[i]
+    # Get angle from origin node to parent node.
+    # Skip if parent_id is a tip.
+    if(isTip.df(df, parent_id) ){ next }
+
+    theta_parent <- getNodeAngle.df(df, origin_id, parent_id)
+  
+    children_ids <- getChild.df(df, parent_id)
+    
+    for( j in seq_along(children_ids)){
+      #delta_x <- df[subtree_node_id, 'x'] - df[origin_id, 'x']
+      #delta_y <- df[subtree_node_id, 'y'] - df[origin_id, 'y']
+      #angles[i] <- atan2(delta_y, delta_x) / pi
+      child_id <- children_ids[j]
+      # Skip if child is parent node of subtree.
+      if( child_id == origin_id ){
+        next
+      }
+      
+      theta_child <- getNodeAngle.df(df, origin_id, child_id)
+
+      # Skip if child node is already inside arc.
+      # if left < right angle (arc crosses 180/-180 quadrant) and child node is not inside arc of tree.
+      # OR if left > right angle (arc crosses 0/360 quadrant) and child node is inside gap
+      if ( (arc['left'] < arc['right'] & !( theta_child > arc['left'] & theta_child < arc['right'])) |
+        (arc['left'] > arc['right'] & ( theta_child < arc['left'] & theta_child > arc['right'])) ){
+        # child node inside arc.
+        next
+      }
+      
+      
+      delta <- theta_child - theta_parent
+      delta_adj <- delta
+      # Correct the delta if parent and child angles cross the 180/-180 half of circle.
+      # If delta > 180
+      if( delta > 1){ # Edge between parent and child cross upper and lower quadrants of cirlce on 180/-180 side.
+        delta_adj <- delta - 2 # delta' = delta - 360
+      # If delta < -180
+      }else if( delta < -1){ # Edge between parent and child cross upper and lower quadrants of cirlce
+        delta_adj <- delta + 2 # delta' = delta - 360
+      }
+      
+      theta_child_adj <- theta_child
+
+      # If angle change from parent to node is positive (anti-clockwise), check left angle
+      if(delta_adj > 0){
+        # If child/parent edges cross the -180/180 quadrant (angle between them is > 180), 
+        # check if right angle and child angle are different signs and adjust if needed.
+        if( abs(delta) > 1){
+          if( arc['left'] > 0 & theta_child < 0){
+            theta_child_adj <- theta_child + 2
+          }else if (arc['left'] < 0 & theta_child > 0){
+            theta_child_adj <- theta_child - 2
+          }
+        }
+        
+          # check if left angle of arc is less than angle of child. Update if true.
+        if( arc['left'] < theta_child_adj ){
+          arc['left'] <- theta_child
         }
+      # If angle change from parent to node is negative (clockwise), check right angle  
+      }else if(delta_adj < 0){
+        # If child/parent edges cross the -180/180 quadrant (angle between them is > 180), 
+        # check if right angle and child angle are different signs and adjust if needed.
+        if( abs(delta) > 1){
+          # Else change in angle from parent to child is negative, then adjust child angle if right angle is a different sign.
+          if( arc['right'] > 0 & theta_child < 0){
+            theta_child_adj <- theta_child + 2
+          }else if (arc['right'] < 0 & theta_child > 0){
+            theta_child_adj <- theta_child - 2
+          }
+        }
+        # check if right angle of arc is greater than angle of child. Update if true.
+        if( arc['right'] > theta_child_adj  ){
+          arc['right'] <- theta_child
+        }  
+        
+      }
+    }
 
-        for (i in seq_along(children)) {
-            child <- children[i]
-            ntip.child <- nb.sp[child]
+  } 
+  # Convert arc angles of [1, -1] to [2,0] domain.
+  arc[arc<0] <- arc[arc<0] + 2
+  return(arc)
+  
+}
 
-            alpha <- (end - start) * ntip.child/curNtip
-            beta <- start + alpha / 2
+##' Rotate the points in a tree data.frame around a pivot node by the angle specified.
+##' 
+##' @title rotateTreePoints.data.fram
+##' @param df tree data.frame
+##' @param pivot_node is the id of the pivot node.
+##' @param nodes list of node numbers that are to be rotated by angle around the pivot_node
+##' @param angle in range [0,2], ie degrees/180, radians/pi
+##' @return updated tree data.frame with points rotated by angle
+rotateTreePoints.df <- function(df, pivot_node, nodes, angle){
+  # Rotate nodes around pivot_node.
+  # x' = cos(angle)*delta_x - sin(angle)*delta_y + delta_x
+  # y' = sin(angle)*delta_x + cos(angle)*delta_y + delta_y
+  
+  cospitheta <- cospi(angle)
+  sinpitheta <- sinpi(angle)
+  for(node in nodes){
+    # Update (x,y) of node
+    delta_x <- df[node, 'x'] - df[pivot_node, 'x']
+    delta_y <- df[node, 'y'] - df[pivot_node, 'y']
+    df[node, 'x'] <- cospitheta * delta_x - sinpitheta * delta_y + df[pivot_node, 'x']
+    df[node, 'y'] <- sinpitheta * delta_x + cospitheta * delta_y + df[pivot_node, 'y']
+    
+    # Update label angle if not root node.
+    # get parent
+    parent_id <- getParent.df(df, node)
+    # if 'node' is not root, then update label angle.
+    if( parent_id != 0){
+      theta_parent_child <- getNodeAngle.df(df, parent_id, node)
+      if(!is.na(theta_parent_child)){
+        # Update label angle
+        df[node, 'angle'] <- -90 - 180 * theta_parent_child * sign(theta_parent_child - 1)
+      }
+    }
+    
+  }
+  return(df)
+}
 
-            if (branch.length == "none") {
-                length.child <- 1
-            } else {
-                length.child <- df[child, "length"]
-            }
+##' Get the angle between the two nodes specified.
+##' 
+##' @title getNodeAngle.df
+##' @param df tree data.frame
+##' @param origin_node_id origin node id number
+##' @param node_id end node id number
+##' @return angle in range [-1, 1], i.e. degrees/180, radians/pi
+getNodeAngle.df <- function(df, origin_node_id, node_id){
+  if(origin_node_id != node_id){
+    delta_x <- df[node_id, 'x'] - df[origin_node_id, 'x']
+    delta_y <- df[node_id, 'y'] - df[origin_node_id, 'y']
+    angle <- atan2(delta_y, delta_x) / pi
+    return( angle )
+  }else{
+    return(NA)
+  }
+}
 
-            df[child, "x"] <- df[curNode, "x"] + cospi(beta) * length.child
-            df[child, "y"] <- df[curNode, "y"] + sinpi(beta) * length.child
-            df[child, "angle"] <- -90 -180 * beta * sign(beta - 1)
-            df[child, "start"] <- start
-            df[child, "end"] <- start + alpha
-            start <- start + alpha
-        }
 
-    }
 
-    return(df)
+##' Get all children of node from tree, including start_node.
+##' 
+##' @title getSubtree
+##' @param tree ape phylo tree object
+##' @param node is the tree node id from which the tree is derived.
+##' @return list of all child node id's from starting node.
+getSubtree <- function(tree, node){
+  
+  subtree <- c(node)
+  i <- 1
+  while( i <= length(subtree)){
+    subtree <- c(subtree, getChild(tree, subtree[i]))
+    # remove any '0' root nodes
+    subtree <- subtree[subtree != 0]
+    i <- i + 1
+  }
+  return(subtree)
 }
 
+# Get all children of node from df tree using breath-first.
+#
+##' @title GetSubtree.df
+##' @param df tree data.frame
+##' @param node id of starting node.
+##' @return list of all child node id's from starting node.
+##' @examples
+GetSubtree.df <- function(df, node){
+  subtree <- c(node)
+  i <- 1
+  while( i <= length(subtree)){
+    subtree <- c(subtree, getChild.df(df, subtree[i]))
+    # remove any '0' root nodes
+    subtree <- subtree[subtree != 0]
+    i <- i + 1
+  }
+  return(subtree)
+}
+
+##' Get all subtrees of specified node. This includes all ancestors and relatives of node and 
+##' return named list of subtrees.
+##' 
+##' @title GetSubtreeUnrooted
+##' @param tree ape phylo tree object
+##' @param node is the tree node id from which the subtrees are derived.
+##' @return named list of subtrees with the root id of subtree and list of node id's making up subtree.
+GetSubtreeUnrooted <- function(tree, node){
+  # if node leaf, return nothing.
+  if( isTip(tree, node) ){
+    # return NA
+    return(NA)
+  }
+  
+  subtrees <- list()
+  
+  # get subtree for each child node.
+  children_ids <- getChild(tree, node)
+  
+  remaining_nodes <- getNodes_by_postorder(tree)
+  # Remove current node from remaining_nodes list.
+  remaining_nodes <- setdiff(remaining_nodes, node)
+  
+  
+  for( child in children_ids ){
+    # Append subtree nodes to list if not 0 (root).
+    subtree <- getSubtree(tree, child)
+    subtrees[[length(subtrees)+1]] <- list( node = child, subtree = subtree)
+    # remove subtree nodes from remaining nodes.
+    remaining_nodes <- setdiff(remaining_nodes, as.integer(unlist(subtrees[[length(subtrees)]]['subtree']) ))
+  }
+  
+  # The remaining nodes that are not found in the child subtrees are the remaining subtree nodes.
+  # ie, parent node and all other nodes. We don't care how they are connect, just their ids.
+  parent_id <- getParent(tree, node)
+  # If node is not root, add remainder of tree nodes as subtree.
+  if( parent_id != 0 & length(remaining_nodes) >= 1){
+    subtrees[[length(subtrees)+1]] <- list( node = parent_id, subtree = remaining_nodes)
+  }
+  
+  return(subtrees)
+}
+
+
+##' Get all subtrees of node, as well as remaining branches of parent (ie, rest of tree structure as subtree)
+##' return named list of subtrees with list name as starting node id.
+##' @title GetSubtreeUnrooted
+##' @param df tree data.frame
+##' @param node is the tree node id from which the subtrees are derived.
+##' @return named list of subtrees with the root id of subtree and list of node id's making up subtree.
+getSubtreeUnrooted.df <- function(df, node){
+  # if node leaf, return nothing.
+  if( isTip.df(df, node) ){
+    return(NA)
+  }
+  
+  subtrees <- list()
+  
+  # get subtree for each child node.
+  children_ids <- getChild.df(df, node)
+  
+  # remaining_nodes <- getNodes_by_postorder(tree)
+  remaining_nodes <- df$node
+  
+  # Remove current node from remaining_nodes list.
+  remaining_nodes <- setdiff(remaining_nodes, node)
+  
+  for( child in children_ids ){
+    subtree <- GetSubtree.df(df, child)
+    # Append subtree nodes to list if more than 1 node in subtree (i.e. not a tip)
+    #if(length(subtree) >= 2){
+      subtrees[[length(subtrees)+1]] <- list( node = child, subtree = subtree)
+      # remove subtree nodes from remaining nodes.
+      remaining_nodes <- setdiff(remaining_nodes, as.integer(unlist(subtrees[[length(subtrees)]]['subtree']) ))
+    #}else{
+      # remove remaining nodes
+    #  remaining_nodes <- setdiff(remaining_nodes, subtree)
+    #}
+  }
+  
+  # The remaining nodes that are not found in the child subtrees are the remaining subtree nodes.
+  # ie, parent node and all other nodes. We don't care how they are connected, just their id.
+  parent_id <- getParent.df(df, node)
+  # If node is not root.
+  if( parent_id != 0 & length(remaining_nodes) >= 1){
+    subtrees[[length(subtrees)+1]] <- list( node = parent_id, subtree = remaining_nodes)
+  }
+  
+  return(subtrees)
+}
+
+
+getRoot.df <- function(df, node){
+  root <- which(is.na(df$parent))
+  return(root)
+}
+
+##' Get parent node id of child node.
+##'
+##' @title getParent.df
+##' @param df tree data.frame
+##' @param node is the node id of child in tree.
+##' @return integer node id of parent
 getParent.df <- function(df, node) {
     i <- which(df$node == node)
-    res <- df$parent[i]
-    if (res == node) {
+    parent_id <- df$parent[i]
+    if (parent_id == node | is.na(parent_id)) {
         ## root node
         return(0)
     }
-    return(res)
+    return(parent_id)
 }
 
 getAncestor.df <- function(df, node) {
     anc <- getParent.df(df, node)
     anc <- anc[anc != 0]
     if (length(anc) == 0) {
-        stop("selected node is root...")
+        # stop("selected node is root...")
+      return(0)
     }
     i <- 1
     while(i<= length(anc)) {
@@ -243,6 +789,12 @@ getAncestor.df <- function(df, node) {
 }
 
 
+##' Get list of child node id numbers of parent node
+##'
+##' @title getChild.df
+##' @param df tree data.frame
+##' @param node is the node id of child in tree.
+##' @return list of child node ids of parent
 getChild.df <- function(df, node) {
     i <- which(df$parent == node)
     if (length(i) == 0) {
@@ -257,7 +809,8 @@ get.offspring.df <- function(df, node) {
     sp <- getChild.df(df, node)
     sp <- sp[sp != 0]
     if (length(sp) == 0) {
-        stop("input node is a tip...")
+        #stop("input node is a tip...")
+      return(0)
     }
 
     i <- 1
@@ -305,7 +858,6 @@ get.offspring.tip <- function(tr, node) {
 ##     N <- Ntip + Nnode
 ##     return(N)
 ## }
-
 getParent <- function(tr, node) {
     if ( node == getRoot(tr) )
         return(0)
@@ -323,7 +875,9 @@ getParent <- function(tr, node) {
 }
 
 getChild <- function(tr, node) {
+    # Get edge matrix from phylo object.
     edge <- tr[["edge"]]
+    # Select all rows that match "node".
     res <- edge[edge[,1] == node, 2]
     ## if (length(res) == 0) {
     ##     ## is a tip
@@ -363,6 +917,16 @@ isRoot <- function(tr, node) {
     getRoot(tr) == node
 }
 
+isTip <- function(tr, node) {
+  children_ids <- getChild(tr, node)
+  length(children_ids) == 0
+}
+
+isTip.df <- function(df, node) {
+  return(df[node, 'isTip'])
+}
+
+
 getNodeName <- function(tr) {
     if (is.null(tr$node.label)) {
         n <- length(tr$tip.label)
@@ -400,7 +964,6 @@ getNodeName <- function(tr) {
 ##     }
 ##     return(root)
 ## }
-
 get.trunk <- function(tr) {
     root <- getRoot(tr)
     path_length <- sapply(1:(root-1), function(x) get.path_length(tr, root, x))
@@ -459,10 +1022,54 @@ get.path_length <- function(phylo, from, to, weight=NULL) {
 }
 
 getNodes_by_postorder <- function(tree) {
-    tree <- reorder.phylo(tree, "postorder")
+  tree <- reorder.phylo(tree, "postorder")
     unique(rev(as.vector(t(tree$edge[,c(2,1)]))))
 }
 
+#getNodes_by_postorder.df <- function(df) {
+    #tree <- reorder.phylo(tree, "postorder")
+    #unique(rev(as.vector(t(tree$edge[,c(2,1)]))))
+#}
+
+##' Get the nodes of tree from root in breadth-first order.
+##' 
+##' @title getNodesBreadthFirst.df
+##' @param df tree data.frame
+##' @return list of node id's in breadth-first order.
+getNodesBreadthFirst.df <- function(df){
+
+  root <- getRoot.df(df)
+  if(isTip.df(df, root)){
+    return(root)
+  }
+
+  tree_size <- nrow(df)
+  # initialise list of nodes
+  res <- root
+  
+  i <- 1
+  while(length(res) < tree_size){
+    parent <- res[i]
+    i <- i + 1
+    
+    # Skip if parent is a tip.
+    if(isTip.df(df, parent)){
+      next
+    }
+    
+    # get children of current parent.
+    children <- getChild.df(df,parent)
+
+    # add children to result
+    res <- c(res, children)
+    
+  }
+  
+  return(res)
+  
+}
+
+
 getXcoord2 <- function(x, root, parent, child, len, start=0, rev=FALSE) {
     x[root] <- start
     x[-root] <- NA  ## only root is set to start, by default 0
@@ -819,6 +1426,7 @@ add_angle_slanted <- function(res) {
     theta <- atan(dy/dx)
     theta[is.na(theta)] <- 0 ## root node
     res$angle <- theta/pi * 180
+    
     branch.y <- (res[res$parent, "y"] + res[, "y"])/2
     idx <- is.na(branch.y)
     branch.y[idx] <- res[idx, "y"]
@@ -929,8 +1537,6 @@ set_branch_length <- function(tree_object, branch.length) {
 
 ##     return(phylo)
 ## }
-
-
 re_assign_ycoord_df <- function(df, currentNode) {
     while(anyNA(df$y)) {
         pNode <- with(df, parent[match(currentNode, node)]) %>% unique
@@ -959,7 +1565,6 @@ re_assign_ycoord_df <- function(df, currentNode) {
 ## is.ggtree <- function(x) inherits(x, 'ggtree')
 
 
-
 calculate_angle <- function(data) {
     data$angle <- 360/(diff(range(data$y)) + 1) * data$y
     return(data)