[Bio] / FigKernelPackages / gjonewicklib.pm Repository:
ViewVC logotype

Diff of /FigKernelPackages/gjonewicklib.pm

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.8, Sun Feb 11 18:25:48 2007 UTC revision 1.14, Tue Jul 7 22:14:22 2009 UTC
# Line 77  Line 77 
77  #  Tree format interconversion:  #  Tree format interconversion:
78  #===============================================================================  #===============================================================================
79  #  #
80    #  $bool      = is_overbeek_tree( $tree )
81    #  $bool      = is_gjonewick_tree( $tree )
82    #
83  #  $gjonewick = overbeek_to_gjonewick( $overbeek )  #  $gjonewick = overbeek_to_gjonewick( $overbeek )
84  #  $overbeek  = gjonewick_to_overbeek( $gjonewick )  #  $overbeek  = gjonewick_to_overbeek( $gjonewick )
85  #  #
# Line 107  Line 110 
110  #  set_newick_c4( $noderef, $listref )  #  set_newick_c4( $noderef, $listref )
111  #  set_newick_c5( $noderef, $listref )  #  set_newick_c5( $noderef, $listref )
112  #  set_newick_desc_list( $noderef, @desclist )  #  set_newick_desc_list( $noderef, @desclist )
113  #  set_newick_desc_i( $noderef1, $i, $noderef2 )  #  set_newick_desc_i( $noderef1, $i, $noderef2 )  # 1-based numbering
114  #  #
115  #  $bool    = newick_is_valid( $noderef )       # verify that tree is valid  #  $bool    = newick_is_valid( $noderef )       # verify that tree is valid
116  #  #
# Line 118  Line 121 
121  #  #
122  #  $n       = newick_tip_count( $noderef )  #  $n       = newick_tip_count( $noderef )
123  #  @tiprefs = newick_tip_ref_list( $noderef )  #  @tiprefs = newick_tip_ref_list( $noderef )
124    # \@tiprefs = newick_tip_ref_list( $noderef )
125  #  @tips    = newick_tip_list( $noderef )  #  @tips    = newick_tip_list( $noderef )
126    # \@tips    = newick_tip_list( $noderef )
127    #
128  #  $tipref  = newick_first_tip_ref( $noderef )  #  $tipref  = newick_first_tip_ref( $noderef )
129  #  $tip     = newick_first_tip( $noderef )  #  $tip     = newick_first_tip( $noderef )
130    #
131  #  @tips    = newick_duplicated_tips( $noderef )  #  @tips    = newick_duplicated_tips( $noderef )
132    # \@tips    = newick_duplicated_tips( $noderef )
133    #
134  #  $bool    = newick_tip_in_tree( $noderef, $tipname )  #  $bool    = newick_tip_in_tree( $noderef, $tipname )
135    #
136  #  @tips    = newick_shared_tips( $tree1, $tree2 )  #  @tips    = newick_shared_tips( $tree1, $tree2 )
137    # \@tips    = newick_shared_tips( $tree1, $tree2 )
138  #  #
139  #  $length  = newick_tree_length( $noderef )  #  $length  = newick_tree_length( $noderef )
140    #
141    #  %tip_distances = newick_tip_distances( $noderef )
142    # \%tip_distances = newick_tip_distances( $noderef )
143    #
144  #  $xmax    = newick_max_X( $noderef )  #  $xmax    = newick_max_X( $noderef )
145  #  ( $tipref,  $xmax ) = newick_most_distant_tip_ref( $noderef )  #  ( $tipref,  $xmax ) = newick_most_distant_tip_ref( $noderef )
146  #  ( $tipname, $xmax ) = newick_most_distant_tip( $noderef )  #  ( $tipname, $xmax ) = newick_most_distant_tip_name( $noderef )
147  #  #
148  #  Tree tip insertion point (tip is on branch of length x that  #  Tree tip insertion point (tip is on branch of length x that
149  #  is inserted into branch connecting node1 and node2, a distance  #  is inserted into branch connecting node1 and node2, a distance
150  #  x1 from node1 and x2 from node2):  #  x1 from node1 and x2 from node2):
151  #  #
152  #  [ $node1, $x1, $node2, $x2, $x ]  #  [ $node1, $x1, $node2, $x2, $x ] = newick_tip_insertion_point( $tree, $tip )
 #           = newick_tip_insertion_point( $tree, $tip )  
153  #  #
154  #  Standardized label for a node in terms of intersection of 3 lowest sorting  #  Standardized label for a node in terms of intersection of 3 lowest sorting
155  #  tips (sort is lower case):  #  tips (sort is lower case):
156  #  #
157  #  @TipOrTips = std_node_name( $Tree, $Node )  #  @TipOrTips = std_node_name( $tree, $node )
158  #  #
159  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
160  #  Paths from root of tree:  #  Paths from root of tree:
# Line 206  Line 220 
220  #  $newtree = reroot_newick_next_to_tip( $tree, $tip )  #  $newtree = reroot_newick_next_to_tip( $tree, $tip )
221  #  $newtree = reroot_newick_to_node( $tree, @node )  #  $newtree = reroot_newick_to_node( $tree, @node )
222  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )
223  #  $newtree = reroot_newick_to_node_ref( $tree, $noderef )  #  $newtree = reroot_newick_between_nodes( $tree, $node1, $node2, $fraction )
224    #  $newtree = reroot_newick_to_midpoint( $tree )           # unweighted
225    #  $newtree = reroot_newick_to_midpoint_w( $tree )         # weight by tips
226  #  $newtree = reroot_newick_to_approx_midpoint( $tree )    # unweighted  #  $newtree = reroot_newick_to_approx_midpoint( $tree )    # unweighted
227  #  $newtree = reroot_newick_to_approx_midpoint_w( $tree )  # weight by tips  #  $newtree = reroot_newick_to_approx_midpoint_w( $tree )  # weight by tips
228  #  $newtree = uproot_tip_rooted_newick( $tree )  #  $newtree = uproot_tip_rooted_newick( $tree )
229  #  $newtree = uproot_newick( $tree )  #  $newtree = uproot_newick( $tree )
230  #  #
231  #  $newtree = prune_from_newick( $tree, $tip )  #  $newtree = prune_from_newick( $tree, $tip )
232    #  $newtree = rooted_newick_subtree( $tree,  @tips )
233    #  $newtree = rooted_newick_subtree( $tree, \@tips )
234  #  $newtree = newick_subtree( $tree,  @tips )  #  $newtree = newick_subtree( $tree,  @tips )
235  #  $newtree = newick_subtree( $tree, \@tips )  #  $newtree = newick_subtree( $tree, \@tips )
236    #  $newtree = newick_covering_subtree( $tree,  @tips )
237    #  $newtree = newick_covering_subtree( $tree, \@tips )
238  #  #
239  #  $newtree = collapse_zero_length_branches( $tree )  #  $newtree = collapse_zero_length_branches( $tree )
240  #  #
# Line 222  Line 242 
242  #  $tree = newick_insert_between_nodes( $tree, $subtree, $node1, $node2, $fraction )  #  $tree = newick_insert_between_nodes( $tree, $subtree, $node1, $node2, $fraction )
243  #  #
244  #===============================================================================  #===============================================================================
245    #  Tree neighborhood: subtree of n tips to represent a larger tree.
246    #===============================================================================
247    #
248    #  Focus around root:
249    #
250    #  $subtree = root_neighborhood_representative_tree( $tree, $n, \%tip_priority )
251    #  $subtree = root_neighborhood_representative_tree( $tree, $n )
252    #  @tips    = root_neighborhood_representative_tips( $tree, $n, \%tip_priority )
253    #  @tips    = root_neighborhood_representative_tips( $tree, $n )
254    # \@tips    = root_neighborhood_representative_tips( $tree, $n, \%tip_priority )
255    # \@tips    = root_neighborhood_representative_tips( $tree, $n )
256    #
257    #  Focus around a tip insertion point (the tip is not in the subtree):
258    #
259    #  $subtree = tip_neighborhood_representative_tree( $tree, $tip, $n, \%tip_priority )
260    #  $subtree = tip_neighborhood_representative_tree( $tree, $tip, $n )
261    #  @tips    = tip_neighborhood_representative_tips( $tree, $tip, $n, \%tip_priority )
262    #  @tips    = tip_neighborhood_representative_tips( $tree, $tip, $n )
263    # \@tips    = tip_neighborhood_representative_tips( $tree, $tip, $n, \%tip_priority )
264    # \@tips    = tip_neighborhood_representative_tips( $tree, $tip, $n )
265    #
266    #===============================================================================
267  #  Tree reading and writing:  #  Tree reading and writing:
268  #===============================================================================  #===============================================================================
269    #  Write machine-readable trees:
270  #  #
271  #   writeNewickTree( $tree )  #   writeNewickTree( $tree )
272  #   writeNewickTree( $tree, $file )  #   writeNewickTree( $tree, $file )
# Line 231  Line 274 
274  #  fwriteNewickTree( $file, $tree )  # Matches the C arg list for f... I/O  #  fwriteNewickTree( $file, $tree )  # Matches the C arg list for f... I/O
275  #  $treestring = swriteNewickTree( $tree )  #  $treestring = swriteNewickTree( $tree )
276  #  $treestring = formatNewickTree( $tree )  #  $treestring = formatNewickTree( $tree )
277    #
278    #  Write human-readable trees:
279    #
280  #  @textlines  = text_plot_newick( $node, $width, $min_dx, $dy )  #  @textlines  = text_plot_newick( $node, $width, $min_dx, $dy )
281  #   printer_plot_newick( $node, $file, $width, $min_dx, $dy )  #   printer_plot_newick( $node, $file, $width, $min_dx, $dy )
282  #  #
283    #  Read trees:
284    #
285  #  $tree  = read_newick_tree( $file )  # reads to a semicolon  #  $tree  = read_newick_tree( $file )  # reads to a semicolon
286  #  @trees = read_newick_trees( $file ) # reads to end of file  #  @trees = read_newick_trees( $file ) # reads to end of file
287  #  $tree  = parse_newick_tree_str( $string )  #  $tree  = parse_newick_tree_str( $string )
# Line 243  Line 291 
291    
292  use Carp;  use Carp;
293  use Data::Dumper;  use Data::Dumper;
294    use strict;
295    
296  require Exporter;  require Exporter;
297    
298  our @ISA = qw(Exporter);  our @ISA = qw(Exporter);
299  our @EXPORT = qw(  our @EXPORT = qw(
300            is_overbeek_tree
301            is_gjonewick_tree
302          overbeek_to_gjonewick          overbeek_to_gjonewick
303          gjonewick_to_overbeek          gjonewick_to_overbeek
   
304          newick_is_valid          newick_is_valid
305          newick_is_rooted          newick_is_rooted
306          newick_is_unrooted          newick_is_unrooted
307          tree_rooted_on_tip          tree_rooted_on_tip
308          newick_is_bifurcating          newick_is_bifurcating
309          newick_tip_count          newick_tip_count
310            newick_tip_ref_list
311          newick_tip_list          newick_tip_list
312    
313          newick_first_tip          newick_first_tip
314          newick_duplicated_tips          newick_duplicated_tips
315          newick_tip_in_tree          newick_tip_in_tree
316          newick_shared_tips          newick_shared_tips
317    
318          newick_tree_length          newick_tree_length
319            newick_tip_distances
320          newick_max_X          newick_max_X
321          newick_most_distant_tip_ref          newick_most_distant_tip_ref
322          newick_most_distant_tip_name          newick_most_distant_tip_name
323    
324            newick_tip_insertion_point
325    
326          std_newick_name          std_newick_name
327    
328          path_to_tip          path_to_tip
# Line 305  Line 360 
360          reroot_newick_next_to_tip          reroot_newick_next_to_tip
361          reroot_newick_to_node          reroot_newick_to_node
362          reroot_newick_to_node_ref          reroot_newick_to_node_ref
363            reroot_newick_between_nodes
364            reroot_newick_to_midpoint
365            reroot_newick_to_midpoint_w
366          reroot_newick_to_approx_midpoint          reroot_newick_to_approx_midpoint
367          reroot_newick_to_approx_midpoint_w          reroot_newick_to_approx_midpoint_w
368          uproot_tip_rooted_newick          uproot_tip_rooted_newick
369          uproot_newick          uproot_newick
370    
371          prune_from_newick          prune_from_newick
372            rooted_newick_subtree
373          newick_subtree          newick_subtree
374            newick_covering_subtree
375          collapse_zero_length_branches          collapse_zero_length_branches
376    
377          newick_insert_at_node          newick_insert_at_node
378          newick_insert_between_nodes          newick_insert_between_nodes
379    
380            root_neighborhood_representative_tree
381            root_neighborhood_representative_tips
382            tip_neighborhood_representative_tree
383            tip_neighborhood_representative_tips
384    
385          writeNewickTree          writeNewickTree
386          fwriteNewickTree          fwriteNewickTree
387          strNewickTree          strNewickTree
# Line 362  Line 427 
427          );          );
428    
429    
 use gjolists qw(  
         common_prefix  
         unique_suffixes  
   
         duplicates  
         random_order  
   
         intersection  
         set_difference  
         );  
   
   
 use strict;  
   
   
430  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
431  #  Internally used definitions  #  Internally used definitions
432  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
# Line 386  Line 436 
436    
437    
438  #===============================================================================  #===============================================================================
439  #  Interconvert Overbeek and gjonewick trees:  #  Interconvert overbeek and gjonewick trees:
440  #===============================================================================  #===============================================================================
441    
442    sub is_overbeek_tree { array_ref( $_[0] ) && array_ref( $_[0]->[2] ) }
443    
444    sub is_gjonewick_tree { array_ref( $_[0] ) && array_ref( $_[0]->[0] ) }
445    
446  sub overbeek_to_gjonewick  sub overbeek_to_gjonewick
447  {  {
448      return () unless ref( $_[0] ) eq 'ARRAY';      return () unless ref( $_[0] ) eq 'ARRAY';
# Line 408  Line 462 
462      return $node;      return $node;
463  }  }
464    
465    
466  #===============================================================================  #===============================================================================
467  #  Extract tree structure values:  #  Extract tree structure values:
468  #===============================================================================  #===============================================================================
# Line 428  Line 483 
483  #  #
484  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
485    
486  sub newick_desc_ref { $_[0]->[0] }  # = ${$_[0]}[0]  sub newick_desc_ref { ref($_[0]) ? $_[0]->[0] : Carp::confess() }  # = ${$_[0]}[0]
487  sub newick_lbl      { ref($_[0]) ? $_[0]->[1] : Carp::confess() }  sub newick_lbl      { ref($_[0]) ? $_[0]->[1] : Carp::confess() }
488  sub newick_x        { $_[0]->[2] }  sub newick_x        { ref($_[0]) ? $_[0]->[2] : Carp::confess() }
489  sub newick_c1       { $_[0]->[3] }  sub newick_c1       { ref($_[0]) ? $_[0]->[3] : Carp::confess() }
490  sub newick_c2       { $_[0]->[4] }  sub newick_c2       { ref($_[0]) ? $_[0]->[4] : Carp::confess() }
491  sub newick_c3       { $_[0]->[5] }  sub newick_c3       { ref($_[0]) ? $_[0]->[5] : Carp::confess() }
492  sub newick_c4       { $_[0]->[6] }  sub newick_c4       { ref($_[0]) ? $_[0]->[6] : Carp::confess() }
493  sub newick_c5       { $_[0]->[7] }  sub newick_c5       { ref($_[0]) ? $_[0]->[7] : Carp::confess() }
494    
495  sub newick_desc_list {  sub newick_desc_list {
496      my $node = $_[0];      my $node = $_[0];
# Line 641  Line 696 
696  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
697  #  List of tip nodes:  #  List of tip nodes:
698  #  #
699  #  @tips = newick_tip_ref_list( $node )  #  @tips = newick_tip_ref_list( $noderef )
700    # \@tips = newick_tip_ref_list( $noderef )
701  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
702  sub newick_tip_ref_list {  sub newick_tip_ref_list {
703      my ( $node, $not_root ) = @_;      my ( $node, $not_root ) = @_;
# Line 658  Line 714 
714          push @list, newick_tip_ref_list( $_, 1 );          push @list, newick_tip_ref_list( $_, 1 );
715      }      }
716    
717      @list;      wantarray ? @list : \@list;
718  }  }
719    
720    
# Line 666  Line 722 
722  #  List of tips:  #  List of tips:
723  #  #
724  #  @tips = newick_tip_list( $node )  #  @tips = newick_tip_list( $node )
725    # \@tips = newick_tip_list( $node )
726  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
727  sub newick_tip_list {  sub newick_tip_list {
728      map { newick_lbl( $_ ) } newick_tip_ref_list( $_[0] );      my @tips = map { newick_lbl( $_ ) } newick_tip_ref_list( $_[0] );
729        wantarray ? @tips : \@tips;
730  }  }
731    
732    
# Line 707  Line 765 
765  #  List of duplicated tip labels.  #  List of duplicated tip labels.
766  #  #
767  #  @tips = newick_duplicated_tips( $node )  #  @tips = newick_duplicated_tips( $node )
768    # \@tips = newick_duplicated_tips( $node )
769  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
770  sub newick_duplicated_tips {  sub newick_duplicated_tips {
771      gjolists::duplicates( newick_tip_list( $_[0] ) );      my @tips = &duplicates( newick_tip_list( $_[0] ) );
772        wantarray ? @tips : \@tips;
773  }  }
774    
775    
# Line 740  Line 800 
800  #  Tips shared between 2 trees.  #  Tips shared between 2 trees.
801  #  #
802  #  @tips = newick_shared_tips( $tree1, $tree2 )  #  @tips = newick_shared_tips( $tree1, $tree2 )
803    # \@tips = newick_shared_tips( $tree1, $tree2 )
804  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
805  sub newick_shared_tips {  sub newick_shared_tips {
806      my ( $Tree1, $Tree2 ) = @_;      my ( $tree1, $tree2 ) = @_;
807      my ( @Tips1 ) = newick_tip_list( $Tree1 );      my $tips1 = newick_tip_list( $tree1 );
808      my ( @Tips2 ) = newick_tip_list( $Tree2 );      my $tips2 = newick_tip_list( $tree2 );
809      gjolists::intersection( \@Tips1, \@Tips2 );      my @tips = &intersection( $tips1, $tips2 );
810        wantarray ? @tips : \@tips;
811  }  }
812    
813    
# Line 767  Line 829 
829    
830    
831  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
832    #  Hash of tip nodes and corresponding distances from root:
833    #
834    #   %tip_distances = newick_tip_distances( $node )
835    #  \%tip_distances = newick_tip_distances( $node )
836    #-------------------------------------------------------------------------------
837    sub newick_tip_distances
838    {
839        my ( $node, $x, $hash ) = @_;
840        my $root = ! $hash;
841        ref( $hash ) eq 'HASH' or $hash = {};
842    
843        $x ||= 0;
844        $x  += newick_x( $node ) || 0;
845    
846        #  Is it a tip?
847    
848        my $n_desc = newick_n_desc( $node );
849        if ( ! $n_desc )
850        {
851            $hash->{ newick_lbl( $node ) } = $x;
852            return $hash;
853        }
854    
855        #  Tree rooted on tip?
856    
857        if ( ( $n_desc == 1 ) && $root && ( newick_lbl( $node ) ) )
858        {
859            $hash->{ newick_lbl( $node ) } = 0;  # Distance to root is zero
860        }
861    
862        foreach ( newick_desc_list( $node ) )
863        {
864            newick_tip_distances( $_, $x, $hash );
865        }
866    
867        wantarray ? %$hash : $hash;
868    }
869    
870    
871    #-------------------------------------------------------------------------------
872  #  Tree max X.  #  Tree max X.
873  #  #
874  #  $xmax = newick_max_X( $node )  #  $xmax = newick_max_X( $node )
# Line 910  Line 1012 
1012    
1013      else      else
1014      {      {
1015          my ( $n1, $x1 ) = describe_desc( $dl->[0] );          my ( $n1, $x1 ) = describe_descendant( $dl->[0] );
1016          my ( $n2, $x2 ) = describe_desc( $dl->[1] );          my ( $n2, $x2 ) = describe_descendant( $dl->[1] );
1017    
1018          if ( @$n1 == 2 ) { push @$n1, $n2->[0] }          if ( @$n1 == 2 ) { push @$n1, $n2->[0] }
1019          if ( @$n2 == 2 )          if ( @$n2 == 2 )
# Line 926  Line 1028 
1028  }  }
1029    
1030    
1031  sub describe_desc  sub describe_descendant
1032  {  {
1033      my $node = shift;      my $node = shift;
1034    
# Line 951  Line 1053 
1053      #  Return the two lowest of those (the third will come from the      #  Return the two lowest of those (the third will come from the
1054      #  other side of the original node).      #  other side of the original node).
1055    
     else  
     {  
1056          my @rep_tips = sort { lc $a cmp lc $b }          my @rep_tips = sort { lc $a cmp lc $b }
1057                         map  { ( sort { lc $a cmp lc $b } newick_tip_list( $_ ) )[0] }                         map  { ( sort { lc $a cmp lc $b } newick_tip_list( $_ ) )[0] }
1058                         @$dl;                         @$dl;
1059          return ( [ @rep_tips[0,1] ], $x );          return ( [ @rep_tips[0,1] ], $x );
1060      }      }
 }  
1061    
1062    
1063  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
# Line 967  Line 1066 
1066  #     Three sorted tip labels intersecting at node, each being smallest  #     Three sorted tip labels intersecting at node, each being smallest
1067  #           of all the tips of their subtrees  #           of all the tips of their subtrees
1068  #  #
1069  #  @TipOrTips = std_node_name( $Tree, $Node )  #  @TipOrTips = std_node_name( $tree, $node )
1070  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1071  sub std_node_name {  sub std_node_name {
1072      my $tree = $_[0];      my $tree = $_[0];
# Line 985  Line 1084 
1084      #  @rest, and keeping the best tip for each subtree.      #  @rest, and keeping the best tip for each subtree.
1085    
1086      my @rest = newick_tip_list( $tree );      my @rest = newick_tip_list( $tree );
1087      my @best = map {      my @best = map
1088              {
1089              my @tips = sort { lc $a cmp lc $b } newick_tip_list( $_ );              my @tips = sort { lc $a cmp lc $b } newick_tip_list( $_ );
1090              @rest = gjolists::set_difference( \@rest, \@tips );              @rest = &set_difference( \@rest, \@tips );
1091              $tips[0];              $tips[0];
1092          } newick_desc_list( $noderef );          } newick_desc_list( $noderef );
1093    
# Line 1075  Line 1175 
1175      my $imax = newick_n_desc( $node );      my $imax = newick_n_desc( $node );
1176      for ( my $i = 1; $i <= $imax; $i++ ) {      for ( my $i = 1; $i <= $imax; $i++ ) {
1177         @path = path_to_node_ref( newick_desc_i( $node, $i ), $noderef, ( @path0, $i ) );         @path = path_to_node_ref( newick_desc_i( $node, $i ), $noderef, ( @path0, $i ) );
1178         if ( @path ) { return @path }          return @path if @path;
1179      }      }
1180    
1181      ();  #  Not found      ();  #  Not found
# Line 1106  Line 1206 
1206      @p2 && @p3 || return ();                             #  Were they found?      @p2 && @p3 || return ();                             #  Were they found?
1207    
1208      # Find the common prefix for each pair of paths      # Find the common prefix for each pair of paths
1209      my @p12 = gjolists::common_prefix( \@p1, \@p2 );      my @p12 = &common_prefix( \@p1, \@p2 );
1210      my @p13 = gjolists::common_prefix( \@p1, \@p3 );      my @p13 = &common_prefix( \@p1, \@p3 );
1211      my @p23 = gjolists::common_prefix( \@p2, \@p3 );      my @p23 = &common_prefix( \@p2, \@p3 );
1212    
1213      # Return the longest common prefix of any two paths      # Return the longest common prefix of any two paths
1214      ( @p12 >= @p13 && @p12 >= @p23 ) ? @p12 :      ( @p12 >= @p13 && @p12 >= @p23 ) ? @p12 :
# Line 1159  Line 1259 
1259      @p1 && @p2 || return undef;                          # Were they found?      @p1 && @p2 || return undef;                          # Were they found?
1260    
1261      # Find the unique suffixes of the two paths      # Find the unique suffixes of the two paths
1262      my ( $suf1, $suf2 ) = gjolists::unique_suffixes( \@p1, \@p2 ); # Common node is lost      my ( $suf1, $suf2 ) = &unique_suffixes( \@p1, \@p2 ); # Common node is lost
1263      my $d1 = @$suf1 ? distance_along_path_2( @$suf1 ) : 0;      my $d1 = @$suf1 ? distance_along_path_2( @$suf1 ) : 0;
1264      my $d2 = @$suf2 ? distance_along_path_2( @$suf2 ) : 0;      my $d2 = @$suf2 ? distance_along_path_2( @$suf2 ) : 0;
1265    
# Line 1184  Line 1284 
1284      my @p2 = path_to_node( $node, $node2 ) or return undef;      my @p2 = path_to_node( $node, $node2 ) or return undef;
1285    
1286      # Find the unique suffixes of the two paths      # Find the unique suffixes of the two paths
1287      my ( $suf1, $suf2 ) = gjolists::unique_suffixes( \@p1, \@p2 ); # Common node is lost      my ( $suf1, $suf2 ) = &unique_suffixes( \@p1, \@p2 ); # Common node is lost
1288      my $d1 = @$suf1 ? distance_along_path_2( @$suf1 ) : 0;      my $d1 = @$suf1 ? distance_along_path_2( @$suf1 ) : 0;
1289      my $d2 = @$suf2 ? distance_along_path_2( @$suf2 ) : 0;      my $d2 = @$suf2 ? distance_along_path_2( @$suf2 ) : 0;
1290    
# Line 1621  Line 1721 
1721  #  #
1722  #  $tree = unaesthetic_newick_tree( $treeref, $dir )  #  $tree = unaesthetic_newick_tree( $treeref, $dir )
1723  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1724  sub unaesthetic_newick_tree {  sub unaesthetic_newick_tree
1725    {
1726      my ( $tree, $dir ) = @_;      my ( $tree, $dir ) = @_;
1727      my %cnt;      my %cnt;
1728    
# Line 1641  Line 1742 
1742  #           = 0 for no change, and  #           = 0 for no change, and
1743  #           > 0 for downward branch (small group first).  #           > 0 for downward branch (small group first).
1744  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1745  sub reorder_against_tip_count {  sub reorder_against_tip_count
1746    {
1747      my ( $node, $cntref, $dir ) = @_;      my ( $node, $cntref, $dir ) = @_;
1748    
1749      my $nd = newick_n_desc( $node );      my $nd = newick_n_desc( $node );
# Line 1680  Line 1782 
1782  #  #
1783  #  $tree = random_order_newick_tree( $tree )  #  $tree = random_order_newick_tree( $tree )
1784  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1785  sub random_order_newick_tree {  sub random_order_newick_tree
1786    {
1787      my ( $node ) = @_;      my ( $node ) = @_;
1788    
1789      my $nd = newick_n_desc( $node );      my $nd = newick_n_desc( $node );
# Line 1689  Line 1792 
1792      #  Reorder this subtree:      #  Reorder this subtree:
1793    
1794      my $dl_ref = newick_desc_ref( $node );      my $dl_ref = newick_desc_ref( $node );
1795      @$dl_ref = gjolists::random_order( @$dl_ref );      @$dl_ref = &random_order( @$dl_ref );
1796    
1797      #  Reorder descendants:      #  Reorder descendants:
1798    
# Line 1704  Line 1807 
1807  #  #
1808  #  $newtree = reroot_newick_by_path( @path )  #  $newtree = reroot_newick_by_path( @path )
1809  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1810  sub reroot_newick_by_path {  sub reroot_newick_by_path
1811    {
1812      my ( $node1, $path1, @rest ) = @_;      my ( $node1, $path1, @rest ) = @_;
1813      array_ref( $node1 ) || return undef;      #  Always expect a node      array_ref( $node1 ) || return undef;      #  Always expect a node
1814    
# Line 1812  Line 1916 
1916    
1917    
1918  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
1919    #  Reroot a newick tree along the path between 2 nodes:
1920    #
1921    #  $tree = reroot_newick_between_nodes( $tree, $node1, $node2, $fraction )
1922    #-------------------------------------------------------------------------------
1923    sub reroot_newick_between_nodes
1924    {
1925        my ( $tree, $node1, $node2, $fraction ) = @_;
1926        array_ref( $tree ) or return undef;
1927        $fraction >= 0 && $fraction <= 1 or return undef;
1928    
1929        #  Find the paths to the nodes:
1930    
1931        my @path1 = path_to_node( $tree, $node1 ) or return undef;
1932        my @path2 = path_to_node( $tree, $node2 ) or return undef;
1933    
1934        reroot_newick_between_nodes_by_path( $tree, \@path1, \@path2, $fraction )
1935    }
1936    
1937    
1938    #-------------------------------------------------------------------------------
1939    #  Reroot a newick tree along the path between 2 nodes:
1940    #
1941    #  $tree = reroot_newick_between_node_refs( $tree, $node1, $node2, $fraction )
1942    #-------------------------------------------------------------------------------
1943    sub reroot_newick_between_node_refs
1944    {
1945        my ( $tree, $node1, $node2, $fraction ) = @_;
1946        array_ref( $tree ) or return undef;
1947    
1948        #  Find the paths to the nodes:
1949    
1950        my @path1 = path_to_node_ref( $tree, $node1 ) or return undef;
1951        my @path2 = path_to_node_ref( $tree, $node2 ) or return undef;;
1952    
1953        reroot_newick_between_nodes_by_path( $tree, \@path1, \@path2, $fraction )
1954    }
1955    
1956    
1957    #-------------------------------------------------------------------------------
1958    #  Reroot a newick tree along the path between 2 nodes defined by paths:
1959    #
1960    #  $tree = reroot_newick_between_nodes_by_path( $tree, $path1, $path2, $fraction )
1961    #-------------------------------------------------------------------------------
1962    sub reroot_newick_between_nodes_by_path
1963    {
1964        my ( $tree, $path1, $path2, $fraction ) = @_;
1965        array_ref( $tree ) and array_ref( $path1 ) and  array_ref( $path2 )
1966           or return undef;
1967        $fraction >= 0 && $fraction <= 1 or return undef;
1968    
1969        my @path1 = @$path1;
1970        my @path2 = @$path2;
1971    
1972        #  Trim the common prefix, saving it:
1973    
1974        my @prefix = ();
1975        while ( $path1[1] == $path2[1] )
1976        {
1977            push @prefix, splice( @path1, 0, 2 );
1978            splice( @path2, 0, 2 );
1979        }
1980    
1981        my ( @path, $dist );
1982        if    ( @path1 < 3 )
1983        {
1984            @path2 >= 3 or return undef;              # node1 = node2
1985            $dist = $fraction * newick_path_length( @path2 );
1986            @path = @path2;
1987        }
1988        elsif ( @path2 < 3 )
1989        {
1990            $dist = ( 1 - $fraction ) * newick_path_length( @path1 );
1991            @path = @path1;
1992        }
1993        else
1994        {
1995            my $dist1 = newick_path_length( @path1 );
1996            my $dist2 = newick_path_length( @path2 );
1997            $dist = $fraction * ( $dist1 + $dist2 ) - $dist1;
1998            @path = ( $dist <= 0 ) ? @path1 : @path2;
1999            $dist = abs( $dist );
2000        }
2001    
2002        #  Descend tree until we reach the insertion branch:
2003    
2004        my $x;
2005        while ( ( $dist > ( $x = newick_x( $path[2] ) ) ) && ( @path > 3 ) )
2006        {
2007            $dist -= $x;
2008            push @prefix, splice( @path, 0, 2 );
2009        }
2010    
2011        #  Insert the new node:
2012    
2013        my $newnode = [ [ $path[2] ], undef, $dist ];
2014        set_newick_desc_i( $path[0], $path[1], $newnode );
2015        set_newick_x( $path[2], ( ( $x > $dist ) ? ( $x - $dist ) : 0 ) );
2016    
2017        #  We can now build the path from root to the new node
2018    
2019        reroot_newick_by_path( @prefix, @path[0,1], $newnode );
2020    }
2021    
2022    
2023    #-------------------------------------------------------------------------------
2024  #  Move root of tree to an approximate midpoint.  #  Move root of tree to an approximate midpoint.
2025  #  #
2026  #  $newtree = reroot_newick_to_approx_midpoint( $tree )  #  $newtree = reroot_newick_to_approx_midpoint( $tree )
# Line 1823  Line 2032 
2032    
2033      my $dists1 = average_to_tips_1( $tree );      my $dists1 = average_to_tips_1( $tree );
2034    
2035      #  Compile average tip to node distances descending, returning midpoint node      #  Compile average tip to node distances descending, returning midpoint
2036        #  cadidates as a list of [ $node1, $node2, $fraction ]
2037    
2038        my @mids = average_to_tips_2( $dists1, undef, undef );
2039    
2040        #  Reroot to first midpoint candidate
2041    
2042        return $tree if ! @mids;
2043        my ( $node1, $node2, $fraction ) = @{ $mids[0] };
2044        reroot_newick_to_node_ref( $tree, $fraction >= 0.5 ? $node2 : $node1 );
2045    }
2046    
2047    
2048    #-------------------------------------------------------------------------------
2049    #  Move root of tree to a midpoint.
2050    #
2051    #  $newtree = reroot_newick_to_midpoint( $tree )
2052    #-------------------------------------------------------------------------------
2053    sub reroot_newick_to_midpoint {
2054        my ( $tree ) = @_;
2055    
2056        #  Compile average tip to node distances assending
2057    
2058        my $dists1 = average_to_tips_1( $tree );
2059    
2060      my $node = average_to_tips_2( $dists1, undef, undef );      #  Compile average tip to node distances descending, returning midpoint
2061        #  [ $node1, $node2, $fraction ]
2062    
2063      #  Reroot      my @mids = average_to_tips_2( $dists1, undef, undef );
2064    
2065      $node ? reroot_newick_to_node_ref( $tree, $node ) : $tree      @mids ? reroot_newick_between_node_refs( $tree, @{ $mids[0] } ) : $tree;
2066  }  }
2067    
2068    
2069    #-------------------------------------------------------------------------------
2070    #  Compile average tip to node distances assending
2071    #-------------------------------------------------------------------------------
2072  sub average_to_tips_1 {  sub average_to_tips_1 {
2073      my ( $node ) = @_;      my ( $node ) = @_;
2074    
# Line 1843  Line 2079 
2079          foreach ( @desc_dists ) { $x_below += $_->[0] }          foreach ( @desc_dists ) { $x_below += $_->[0] }
2080          $x_below /= @desc_dists;          $x_below /= @desc_dists;
2081      }      }
2082    
2083      my $x = newick_x( $node ) || 0;      my $x = newick_x( $node ) || 0;
2084      my $x_net = $x_below + $x;      my $x_net = $x_below + $x;
2085    
# Line 1850  Line 2087 
2087  }  }
2088    
2089    
2090    #-------------------------------------------------------------------------------
2091    #  Compile average tip to node distances descending, returning midpoint as
2092    #  [ $node1, $node2, $fraction_of_dist_between ]
2093    #-------------------------------------------------------------------------------
2094  sub average_to_tips_2 {  sub average_to_tips_2 {
2095      my ( $dists1, $x_above, $anc_node ) = @_;      my ( $dists1, $x_above, $anc_node ) = @_;
2096      my ( undef, $x, $x_below, $desc_list, $node ) = @$dists1;      my ( undef, $x, $x_below, $desc_list, $node ) = @$dists1;
2097    
2098      #  Are we done?  Root is in this node's branch, or "above"?      #  Are we done?  Root is in this node's branch, or "above"?
2099    
2100      # defined( $x_above ) and print STDERR "x_above = $x_above\n";      my @mids = ();
     # print STDERR "x       = $x\n";  
     # print STDERR "x_below = $x_below\n";  
     # print STDERR "n_desc  = ", scalar @$desc_list, "\n\n";  
   
2101      if ( defined( $x_above ) && ( ( $x_above + $x ) >= $x_below ) )      if ( defined( $x_above ) && ( ( $x_above + $x ) >= $x_below ) )
2102      {      {
2103          #  At this point the root can only be in this node's branch,          #  At this point the root can only be in this node's branch,
# Line 1872  Line 2109 
2109    
2110          if ( ( $x_below + $x ) >= $x_above )          if ( ( $x_below + $x ) >= $x_above )
2111          {          {
2112              return ( $x_above >= $x_below ) ? $anc_node : $node;              #  We will need to make a new node for the root, $fract of
2113          }              #  the way from $node to $anc_node:
2114          else              my $fract = ( $x > 0 ) ? 0.5 * ( ( $x_above - $x_below ) / $x + 1 )
2115          {                                     : 0.5;
2116              return undef;              push @mids, [ $node, $anc_node, $fract ];
2117          }          }
2118      }      }
2119    
2120      #  The root must be somewhere below this node:      #  The root might be somewhere below this node:
2121    
2122      my $n_1      =   @$desc_list - ( $anc_node ? 0 : 1 );      my $n_1      =   @$desc_list - ( $anc_node ? 0 : 1 );
2123      my $ttl_dist = ( @$desc_list * $x_below ) + ( defined( $x_above ) ? ( $x_above + $x ) : 0 );      my $ttl_dist = ( @$desc_list * $x_below ) + ( defined( $x_above ) ? ( $x_above + $x ) : 0 );
# Line 1890  Line 2127 
2127          #  If input tree is tip_rooted, $n-1 can be 0, so:          #  If input tree is tip_rooted, $n-1 can be 0, so:
2128    
2129          my $above2 = $n_1 ? ( ( $ttl_dist - $_->[0] ) / $n_1 ) : 0;          my $above2 = $n_1 ? ( ( $ttl_dist - $_->[0] ) / $n_1 ) : 0;
2130          my $root = average_to_tips_2( $_, $above2, $node );          push @mids, average_to_tips_2( $_, $above2, $node );
         if ( $root ) { return $root }  
2131      }      }
2132    
2133      #  Was not anywhere below this node (oh-oh):      return @mids;
   
     return undef;  
2134  }  }
2135    
2136    
# Line 1908  Line 2142 
2142  sub reroot_newick_to_approx_midpoint_w {  sub reroot_newick_to_approx_midpoint_w {
2143      my ( $tree ) = @_;      my ( $tree ) = @_;
2144    
2145        #  Compile average tip to node distances assending from tips
2146    
2147        my $dists1 = average_to_tips_1_w( $tree );
2148    
2149        #  Compile average tip to node distances descending, returning midpoints
2150    
2151        my @mids = average_to_tips_2_w( $dists1, undef, undef, undef );
2152    
2153        #  Reroot to first midpoint candidate
2154    
2155        return $tree if ! @mids;
2156        my ( $node1, $node2, $fraction ) = @{ $mids[0] };
2157        reroot_newick_to_node_ref( $tree, $fraction >= 0.5 ? $node2 : $node1 );
2158    }
2159    
2160    
2161    #-------------------------------------------------------------------------------
2162    #  Move root of tree to an approximate midpoint.  Weight by tips.
2163    #
2164    #  $newtree = reroot_newick_to_midpoint_w( $tree )
2165    #-------------------------------------------------------------------------------
2166    sub reroot_newick_to_midpoint_w {
2167        my ( $tree ) = @_;
2168    
2169      #  Compile average tip to node distances assending      #  Compile average tip to node distances assending
2170    
2171      my $dists1 = average_to_tips_1_w( $tree );      my $dists1 = average_to_tips_1_w( $tree );
2172    
2173      #  Compile average tip to node distances descending, returning midpoint node      #  Compile average tip to node distances descending, returning midpoint node
2174    
2175      my $node = average_to_tips_2_w( $dists1, undef, undef, undef );      my @mids = average_to_tips_2_w( $dists1, undef, undef, undef );
2176    
2177      #  Reroot      #  Reroot at first candidate midpoint
2178    
2179      $node ? reroot_newick_to_node_ref( $tree, $node ) : $tree      @mids ? reroot_newick_between_node_refs( $tree, @{ $mids[0] } ) : $tree;
2180  }  }
2181    
2182    
# Line 1939  Line 2197 
2197          }          }
2198          $x_below /= $n_below;          $x_below /= $n_below;
2199      }      }
2200    
2201      my $x = newick_x( $node ) || 0;      my $x = newick_x( $node ) || 0;
2202      my $x_net = $x_below + $x;      my $x_net = $x_below + $x;
2203    
# Line 1952  Line 2211 
2211    
2212      #  Are we done?  Root is in this node's branch, or "above"?      #  Are we done?  Root is in this node's branch, or "above"?
2213    
2214      # defined( $x_above ) and print STDERR "x_above = $x_above\n";      my @mids = ();
     # print STDERR "x       = $x\n";  
     # print STDERR "x_below = $x_below\n";  
     # print STDERR "n_desc  = ", scalar @$desc_list, "\n\n";  
   
2215      if ( defined( $x_above ) && ( ( $x_above + $x ) >= $x_below ) )      if ( defined( $x_above ) && ( ( $x_above + $x ) >= $x_below ) )
2216      {      {
2217          #  At this point the root can only be in this node's branch,          #  At this point the root can only be in this node's branch,
# Line 1964  Line 2219 
2219          #  would mean that the midpoint is actually down a different          #  would mean that the midpoint is actually down a different
2220          #  path from the root of the current tree).          #  path from the root of the current tree).
2221          #          #
2222          #  Is the root in the current branch?          #  Is their a root in the current branch?
2223    
2224          if ( ( $x_below + $x ) >= $x_above )          if ( ( $x_below + $x ) >= $x_above )
2225          {          {
2226              return ( $x_above >= $x_below ) ? $anc_node : $node;              #  We will need to make a new node for the root, $fract of
2227          }              #  the way from $node to $anc_node:
2228          else              my $fract = ( $x > 0 ) ? 0.5 * ( ( $x_above - $x_below ) / $x + 1 )
2229          {                                     : 0.5;
2230              return undef;              push @mids, [ $node, $anc_node, $fract ];
2231          }          }
2232      }      }
2233    
# Line 1992  Line 2247 
2247    
2248          my $x_above2 = $n_above2 ? ( ( $ttl_w_dist - $n_2 * $_->[0] ) / $n_above2 )          my $x_above2 = $n_above2 ? ( ( $ttl_w_dist - $n_2 * $_->[0] ) / $n_above2 )
2249                                   : 0;                                   : 0;
2250          my $root = average_to_tips_2_w( $_, $x_above2, $n_above2 || 1, $node );          push @mids, average_to_tips_2_w( $_, $x_above2, $n_above2 || 1, $node );
         if ( $root ) { return $root }  
2251      }      }
2252    
2253      #  Was not anywhere below this node (oh-oh):      return @mids;
   
     return undef;  
2254  }  }
2255    
2256    
# Line 2313  Line 2565 
2565    
2566    
2567  #-------------------------------------------------------------------------------  #-------------------------------------------------------------------------------
2568    #  Produce a potentially rooted subtree with the desired tips:
2569    #
2570    #     Except for (some) tip nodes, the tree produced is a copy.
2571    #     There is no check that requested tips exist.
2572    #
2573    #  $newtree = rooted_newick_subtree( $tree,  @tips )
2574    #  $newtree = rooted_newick_subtree( $tree, \@tips )
2575    #-------------------------------------------------------------------------------
2576    sub rooted_newick_subtree {
2577        my ( $tr, @tips ) = @_;
2578        if ( @tips == 1 && ref( $tips[0] ) eq "ARRAY" ) { @tips = @{ $tips[0] } }
2579    
2580        if ( @tips < 2 ) { return undef }
2581        my $keephash = { map { ( $_, 1 ) } @tips };
2582        my $tr2 = subtree1( $tr, $keephash );
2583        $tr2->[2] = undef if $tr2;                   # undef root branch length
2584        $tr2;
2585    }
2586    
2587    
2588    #-------------------------------------------------------------------------------
2589  #  Produce a subtree with the desired tips:  #  Produce a subtree with the desired tips:
2590  #  #
2591  #     Except for (some) tip nodes, the tree produced is a copy.  #     Except for (some) tip nodes, the tree produced is a copy.
# Line 2386  Line 2659 
2659  }  }
2660    
2661    
2662    #-------------------------------------------------------------------------------
2663    #  The smallest subtree of rooted tree that includes @tips:
2664    #
2665    #    $node = newick_covering_subtree( $tree,  @tips )
2666    #    $node = newick_covering_subtree( $tree, \@tips )
2667    #-------------------------------------------------------------------------------
2668    
2669    sub newick_covering_subtree {
2670        my $tree = shift;
2671        my %tips = map { $_ => 1 } ( ( ref( $_[0] ) eq 'ARRAY' ) ? @{ $_[0] } : @_ );
2672    
2673        #  Return smallest covering node, if any:
2674    
2675        ( newick_covering_subtree( $tree, \%tips ) )[ 0 ];
2676    }
2677    
2678    
2679    sub newick_covering_subtree_1 {
2680        my ( $node, $tips ) = @_;
2681        my $n_cover = 0;
2682        my @desc = newick_desc_list( $node );
2683        if ( @desc )
2684        {
2685            foreach ( @desc )
2686            {
2687                my ( $subtree, $n ) = newick_covering_subtree_1( $_, $tips );
2688                return ( $subtree, $n ) if $subtree;
2689                $n_cover += $n;
2690            }
2691        }
2692        elsif ( $tips->{ newick_lbl( $node ) } )
2693        {
2694            $n_cover++;
2695        }
2696    
2697        #  If all tips are covered, return node
2698    
2699        ( $n_cover == keys %$tips ) ? ( $node, $n_cover ) : ( undef, $n_cover );
2700    }
2701    
2702    
2703    #===============================================================================
2704    #
2705    #  Representative subtrees
2706    #
2707    #===============================================================================
2708    #  Find subtree of size n representating vicinity of the root:
2709    #
2710    #   $subtree = root_neighborhood_representative_tree( $tree, $n, \%tip_priority )
2711    #   $subtree = root_neighborhood_representative_tree( $tree, $n )
2712    #
2713    #  Note that if $tree is rooted, then the subtree will also be.  This can have
2714    #  consequences on downstream programs.
2715    #-------------------------------------------------------------------------------
2716    sub root_neighborhood_representative_tree
2717    {
2718        my ( $tree, $n, $tip_priority ) = @_;
2719        array_ref( $tree ) && ( $n >= 2 ) or return undef;
2720        if ( newick_tip_count( $tree ) <= $n ) { return $tree }
2721    
2722        $tip_priority ||= default_tip_priority( $tree );
2723        my @tips = map { representative_tip_of_newick_node( $_, $tip_priority ) }
2724                   root_proximal_newick_subtrees( $tree, $n );
2725    
2726        newick_subtree( copy_newick_tree( $tree ), \@tips );
2727    }
2728    
2729    
2730    #-------------------------------------------------------------------------------
2731    #  Find n tips to represent tree lineages in vicinity of another tip.
2732    #  Default tip priority is short total branch length.
2733    #
2734    #  \@tips = root_neighborhood_representative_tips( $tree, $n, \%tip_priority )
2735    #   @tips = root_neighborhood_representative_tips( $tree, $n, \%tip_priority )
2736    #  \@tips = root_neighborhood_representative_tips( $tree, $n )
2737    #   @tips = root_neighborhood_representative_tips( $tree, $n )
2738    #-------------------------------------------------------------------------------
2739    sub root_neighborhood_representative_tips
2740    {
2741        my ( $tree, $n, $tip_priority ) = @_;
2742        array_ref( $tree ) && ( $n >= 2 ) or return undef;
2743    
2744        my @tips;
2745        if ( newick_tip_count( $tree ) <= $n )
2746        {
2747            @tips = newick_tip_list( $tree );
2748        }
2749        else
2750        {
2751            $tip_priority ||= default_tip_priority( $tree );
2752            @tips = map { representative_tip_of_newick_node( $_, $tip_priority ) }
2753                    root_proximal_newick_subtrees( $tree, $n );
2754        }
2755    
2756        wantarray ? @tips : \@tips;
2757    }
2758    
2759    
2760    #-------------------------------------------------------------------------------
2761    #  Find subtree of size n representating vicinity of a tip:
2762    #
2763    #   $subtree = tip_neighborhood_representative_tree( $tree, $tip, $n, \%tip_priority )
2764    #   $subtree = tip_neighborhood_representative_tree( $tree, $tip, $n )
2765    #-------------------------------------------------------------------------------
2766    sub tip_neighborhood_representative_tree
2767    {
2768        my ( $tree, $tip, $n, $tip_priority ) = @_;
2769        array_ref( $tree ) && $tip && ( $n >= 2 ) or return undef;
2770        newick_tip_in_tree( $tree, $tip ) or return undef;
2771    
2772        my $tree1 = copy_newick_tree( $tree );
2773        if ( newick_tip_count( $tree1 ) - 1 <= $n )
2774        {
2775            return prune_from_newick( $tree1, $tip )
2776        }
2777    
2778        $tree1 = reroot_newick_to_tip( $tree1, $tip );
2779        $tree1 = newick_desc_i( $tree1, 1 );        # Node immediately below tip
2780        my @tips = root_neighborhood_representative_tips( $tree1, $n, $tip_priority );
2781        newick_subtree( copy_newick_tree( $tree ), \@tips );
2782    }
2783    
2784    
2785    #-------------------------------------------------------------------------------
2786    #  Find n tips to represent tree lineages in vicinity of another tip.
2787    #  Default tip priority is short total branch length.
2788    #
2789    #  \@tips = tip_neighborhood_representative_tips( $tree, $tip, $n, \%tip_priority )
2790    #   @tips = tip_neighborhood_representative_tips( $tree, $tip, $n, \%tip_priority )
2791    #  \@tips = tip_neighborhood_representative_tips( $tree, $tip, $n )
2792    #   @tips = tip_neighborhood_representative_tips( $tree, $tip, $n )
2793    #-------------------------------------------------------------------------------
2794    sub tip_neighborhood_representative_tips
2795    {
2796        my ( $tree, $tip, $n, $tip_priority ) = @_;
2797        array_ref( $tree ) && $tip && ( $n >= 2 ) or return undef;
2798        newick_tip_in_tree( $tree, $tip ) or return undef;
2799    
2800        my @tips = newick_tip_list( $tree );
2801        if ( newick_tip_count( $tree ) - 1 <= $n )
2802        {
2803            @tips = grep { $_ ne $tip } @tips;
2804        }
2805        else
2806        {
2807            my $tree1 = copy_newick_tree( $tree );
2808            $tree1 = reroot_newick_to_tip( $tree1, $tip );
2809            $tree1 = newick_desc_i( $tree1, 1 );        # Node immediately below tip
2810            @tips = root_neighborhood_representative_tips( $tree1, $n, $tip_priority );
2811        }
2812    
2813        wantarray ? @tips : \@tips;
2814    }
2815    
2816    
2817    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2818    #  Anonymous hash of the negative distance from root to each tip:
2819    #
2820    #   \%tip_priority = default_tip_priority( $tree )
2821    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2822    sub default_tip_priority
2823    {
2824        my ( $tree ) = @_;
2825        my $tip_distances = newick_tip_distances( $tree ) || {};
2826        return { map { $_ => -$tip_distances->{$_} } keys %$tip_distances };
2827    }
2828    
2829    
2830    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2831    #  Select a tip from a subtree base on a priority value:
2832    #
2833    #    $tip = representative_tip_of_newick_node( $node, \%tip_priority )
2834    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2835    sub representative_tip_of_newick_node
2836    {
2837        my ( $node, $tip_priority ) = @_;
2838        my ( $tip ) = sort { $b->[1] <=> $a->[1] }   # The best
2839                      map  { [ $_, $tip_priority->{ $_ } ] }
2840                      newick_tip_list( $node );
2841        $tip->[0];                                   # Label from label-priority pair
2842    }
2843    
2844    
2845    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2846    #  Find n subtrees focused around the root of a tree.  Typically each will
2847    #  then be reduced to a single tip to make a representative tree:
2848    #
2849    #   @subtrees = root_proximal_newick_subtrees( $tree, $n )
2850    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2851    sub root_proximal_newick_subtrees
2852    {
2853        my ( $tree, $n ) = @_;
2854        my $node_start_end = newick_branch_intervals( $tree );
2855        n_representative_branches( $n, $node_start_end );
2856    }
2857    
2858    
2859    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2860    #   @node_start_end = newick_branch_intervals( $tree )
2861    #  \@node_start_end = newick_branch_intervals( $tree )
2862    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2863    sub newick_branch_intervals
2864    {
2865        my ( $node, $parent_x ) = @_;
2866        $parent_x ||= 0;
2867        my ( $desc, undef, $dx ) = @$node;
2868        my $x = $parent_x + $dx;
2869        my $interval = [ $node, $parent_x, $desc && @$desc ? $x : 1e100 ];
2870        my @intervals = ( $interval,
2871                          map { &newick_branch_intervals( $_, $x ) } @$desc
2872                        );
2873        return wantarray ? @intervals : \@intervals;
2874    }
2875    
2876    
2877    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2878    #   @ids = n_representative_branches( $n,  @id_start_end )
2879    #   @ids = n_representative_branches( $n, \@id_start_end )
2880    #  \@ids = n_representative_branches( $n,  @id_start_end )
2881    #  \@ids = n_representative_branches( $n, \@id_start_end )
2882    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2883    sub n_representative_branches
2884    {
2885        my $n = shift;
2886        #  Sort intervals by start point:
2887        my @unprocessed = sort { $a->[1] <=> $b->[1] }
2888                          ( @_ == 1 ) ? @{ $_[0] } : @_;
2889        my @active = ();
2890        my ( $interval, $current_point );
2891        foreach $interval ( @unprocessed )
2892        {
2893            $current_point = $interval->[1];
2894            #  Filter out intervals that have ended.  This is N**2 in the number
2895            #  of representatives.  Fixing this would require maintaining a sorted
2896            #  active list.
2897            @active = grep { $_->[2] > $current_point } @active;
2898            push @active, $interval;
2899            last if ( @active >= $n );
2900        }
2901    
2902        my @ids = map { $_->[0] } @active;
2903        return wantarray() ? @ids : \@ids;
2904    }
2905    
2906    
2907  #===============================================================================  #===============================================================================
2908  #  #
2909  #  Tree writing and reading  #  Tree writing and reading
# Line 2553  Line 3071 
3071      my ( $fh, $close ) = open_input( $file );      my ( $fh, $close ) = open_input( $file );
3072      my $tree;      my $tree;
3073      my @lines = ();      my @lines = ();
3074      while ( defined( $_ = <$fh> ) )      foreach ( <$fh> )
3075      {      {
3076          chomp;          chomp;
3077          push @lines, $_;          push @lines, $_;
# Line 2575  Line 3093 
3093      my ( $fh, $close ) = open_input( $file );      my ( $fh, $close ) = open_input( $file );
3094      my @trees = ();      my @trees = ();
3095      my @lines = ();      my @lines = ();
3096      while ( defined( $_ = <$fh> ) )      foreach ( <$fh> )
3097      {      {
3098          chomp;          chomp;
3099          push @lines, $_;          push @lines, $_;
# Line 2762  Line 3280 
3280      #  Loop while it is a comment:      #  Loop while it is a comment:
3281      while ( substr( $s, $ind, 1 ) eq "[" ) {      while ( substr( $s, $ind, 1 ) eq "[" ) {
3282          $ind++;          $ind++;
3283            my $depth = 1;
3284            my $ind2  = $ind;
3285    
3286          #  Find end          #  Find end
3287          if ( substr( $s, $ind ) !~ /^([^]]*)\]/ ) {          while ( $depth > 0 )
3288            {
3289                if ( substr( $s, $ind2 ) =~ /^([^][]*\[)/ )     # nested [ ... ]
3290                {
3291                    $ind2 += length( $1 );  #  Points at char just past [
3292                    $depth++;               #  If nested comments are allowed
3293                }
3294                elsif ( substr( $s, $ind2 ) =~ /^([^][]*\])/ )  # close bracket
3295                {
3296                    $ind2 += length( $1 );  #  Points at char just past ]
3297                    $depth--;
3298                }
3299                else
3300                {
3301              treeParseError( "comment missing closing bracket '["              treeParseError( "comment missing closing bracket '["
3302                             . substr( $s, $ind ) . "'" )                             . substr( $s, $ind ) . "'" )
3303          }          }
3304          my $comment = $1;          }
3305    
3306          #  Save if it includes any "text"          my $comment = substr( $s, $ind, $ind2-$ind-1 );
3307          if ( $comment =~ m/\S/ ) { push @clist, $comment }          if ( $comment =~ m/\S/ ) { push @clist, $comment }
3308    
3309          $ind += length( $comment ) + 1;     #  Comment plus closing bracket          $ind = $ind2;
3310    
3311          #  Skip white space          #  Skip white space
3312          if ( substr( $s, $ind ) =~ /^(\s+)/ ) { $ind += length( $1 ) }          if ( substr( $s, $ind ) =~ /^(\s+)/ ) { $ind += length( $1 ) }
# Line 2792  Line 3325 
3325  #===============================================================================  #===============================================================================
3326  #  Make a printer plot of a tree:  #  Make a printer plot of a tree:
3327  #  #
3328  #     $node   newick tree root node  #  printer_plot_newick( $node, $file, $width, $min_dx, $dy )
3329  #     $file   undef (= \*STDOUT), \*STDOUT, \*STDERR, or a file name.  #  printer_plot_newick( $node, $file, \%options )
3330  #     $width  the approximate characters for the tree without labels  #
3331  #     $min_dx the minimum horizontal branch length  #     $node   # newick tree root node
3332  #     $dy     the vertical space per taxon  #     $file   # undef = \*STDOUT, \*FH, or a file name.
3333    #     $width  # the approximate characters for the tree without labels (D = 68)
3334    #     $min_dx # the minimum horizontal branch length (D = 2)
3335    #     $dy     # the vertical space per taxon (D = 1, most compressed)
3336    #
3337    #  Options:
3338    #
3339    #    dy     => nat_number    # the vertical space per taxon
3340    #    chars  => key           # line drawing character set:
3341    #                            #     html_unicode
3342    #                            #     text (default)
3343    #    min_dx => whole_number  # the minimum horizontal branch length
3344    #    width  => whole_number  # approximate tree width without labels
3345  #  #
 #  printer_plot_newick( $node, $file (D=\*STDOUT), $width (D=68), $min_dx (D=2), $dy (D=1) )  
3346  #===============================================================================  #===============================================================================
3347  sub printer_plot_newick {  sub printer_plot_newick
3348      my ( $node, $file, $width, $min_dx, $dy ) = @_;  {
3349        my ( $node, $file, @opts ) = @_;
3350    
3351      my ( $fh, $close ) = open_output( $file );      my ( $fh, $close ) = open_output( $file );
3352      $fh or return;      $fh or return;
3353    
3354      print $fh join( "\n", text_plot_newick( $node, $width, $min_dx, $dy ) ), "\n";      my $html = $opts[0] && ref($opts[0]) eq 'HASH'
3355                            && $opts[0]->{ chars }
3356                            && $opts[0]->{ chars } =~ /html/;
3357        print $fh '<PRE>' if $html;
3358        print $fh join( "\n", text_plot_newick( $node, @opts ) ), "\n";
3359        print $fh "</PRE>\n" if $html;
3360    
3361      if ( $close ) { close $fh }      if ( $close ) { close $fh }
3362  }  }
3363    
3364    
3365  #===============================================================================  #===============================================================================
3366    #  Character sets for printer plot trees:
3367    #-------------------------------------------------------------------------------
3368    
3369    my %char_set =
3370      ( text1     => { space  => ' ',
3371                       horiz  => '-',
3372                       vert   => '|',
3373                       el_d_r => '/',
3374                       el_u_r => '\\',
3375                       el_d_l => '\\',
3376                       el_u_l => '/',
3377                       tee_l  => '+',
3378                       tee_r  => '+',
3379                       tee_u  => '+',
3380                       tee_d  => '+',
3381                       half_l => '-',
3382                       half_r => '-',
3383                       half_u => '|',
3384                       half_d => '|',
3385                       cross  => '+',
3386                     },
3387        text2     => { space  => ' ',
3388                       horiz  => '-',
3389                       vert   => '|',
3390                       el_d_r => '+',
3391                       el_u_r => '+',
3392                       el_d_l => '+',
3393                       el_u_l => '+',
3394                       tee_l  => '+',
3395                       tee_r  => '+',
3396                       tee_u  => '+',
3397                       tee_d  => '+',
3398                       half_l => '-',
3399                       half_r => '-',
3400                       half_u => '|',
3401                       half_d => '|',
3402                       cross  => '+',
3403                     },
3404        html_box  => { space  => '&nbsp;',
3405                       horiz  => '&#9472;',
3406                       vert   => '&#9474;',
3407                       el_d_r => '&#9484;',
3408                       el_u_r => '&#9492;',
3409                       el_d_l => '&#9488;',
3410                       el_u_l => '&#9496;',
3411                       tee_l  => '&#9508;',
3412                       tee_r  => '&#9500;',
3413                       tee_u  => '&#9524;',
3414                       tee_d  => '&#9516;',
3415                       half_l => '&#9588;',
3416                       half_r => '&#9590;',
3417                       half_u => '&#9589;',
3418                       half_d => '&#9591;',
3419                       cross  => '&#9532;',
3420                     },
3421        utf8_box  => { space  => ' ',
3422                       horiz  => chr(226) . chr(148) . chr(128),
3423                       vert   => chr(226) . chr(148) . chr(130),
3424                       el_d_r => chr(226) . chr(148) . chr(140),
3425                       el_u_r => chr(226) . chr(148) . chr(148),
3426                       el_d_l => chr(226) . chr(148) . chr(144),
3427                       el_u_l => chr(226) . chr(148) . chr(152),
3428                       tee_l  => chr(226) . chr(148) . chr(164),
3429                       tee_r  => chr(226) . chr(148) . chr(156),
3430                       tee_u  => chr(226) . chr(148) . chr(180),
3431                       tee_d  => chr(226) . chr(148) . chr(172),
3432                       half_l => chr(226) . chr(149) . chr(180),
3433                       half_r => chr(226) . chr(149) . chr(182),
3434                       half_u => chr(226) . chr(149) . chr(181),
3435                       half_d => chr(226) . chr(149) . chr(183),
3436                       cross  => chr(226) . chr(148) . chr(188),
3437                     },
3438      );
3439    
3440    %{ $char_set{ html1 } } = %{ $char_set{ text1 } };
3441    $char_set{ html1 }->{ space } = '&nbsp;';
3442    
3443    %{ $char_set{ html2 } } = %{ $char_set{ text2 } };
3444    $char_set{ html2 }->{ space } = '&nbsp;';
3445    
3446    #  Define some synonyms
3447    
3448    $char_set{ html } = $char_set{ html_box };
3449    $char_set{ line } = $char_set{ utf8_box };
3450    $char_set{ symb } = $char_set{ utf8_box };
3451    $char_set{ text } = $char_set{ text1 };
3452    $char_set{ utf8 } = $char_set{ utf8_box };
3453    
3454    #  Define tree formats and synonyms
3455    
3456    my %tree_format =
3457        ( text         => 'text',
3458          tree_tab_lbl => 'tree_tab_lbl',
3459          tree_lbl     => 'tree_lbl',
3460          chrlist_lbl  => 'chrlist_lbl',
3461          raw          => 'chrlist_lbl',
3462        );
3463    
3464    #===============================================================================
3465  #  Make a text plot of a tree:  #  Make a text plot of a tree:
3466  #  #
3467  #     $node   newick tree root node  #  @lines = text_plot_newick( $node, $width, $min_dx, $dy )
3468  #     $width  the approximate characters for the tree without labels  #  @lines = text_plot_newick( $node, \%options )
3469  #     $min_dx the minimum horizontal branch length  #
3470  #     $dy     the vertical space per taxon  #     $node   # newick tree root node
3471    #     $width  # the approximate characters for the tree without labels (D = 68)
3472    #     $min_dx # the minimum horizontal branch length (D = 2)
3473    #     $dy     # the vertical space per taxon (D = 1, most compressed)
3474    #
3475    #  Options:
3476    #
3477    #    chars  => keyword       # the output character set for the tree
3478    #    dy     => nat_number    # the vertical space per taxon
3479    #    format => keyword       # output format of each line
3480    #    min_dx => whole_number  # the minimum horizontal branch length
3481    #    width  => whole_number  # approximate tree width without labels
3482    #
3483    #  Character sets:
3484    #
3485    #    html       #  synonym of html1
3486    #    html_box   #  html encoding of unicode box drawing characters
3487    #    html1      #  text1 with nonbreaking spaces
3488    #    html2      #  text2 with nonbreaking spaces
3489    #    line       #  synonym of utf8_box
3490    #    raw        #  pass out the internal representation
3491    #    symb       #  synonym of utf8_box
3492    #    text       #  synonym of text1 (Default)
3493    #    text1      #  ascii characters: - + | / \ and space
3494    #    text2      #  ascii characters: - + | + + and space
3495    #    utf8       #  synonym of utf8_box
3496    #    utf8_box   #  utf8 encoding of unicode box drawing characters
3497    #
3498    #  Formats for row lines:
3499    #
3500    #    text           #    $textstring              # Default
3501    #    tree_tab_lbl   #    $treestr \t $labelstr
3502    #    tree_lbl       # [  $treestr,  $labelstr ]
3503    #    chrlist_lbl    # [ \@treechar, $labelstr ]   # Forced with raw chars
3504    #    raw            #  synonym of chrlist_lbl
3505  #  #
 #  @textlines = text_plot_newick( $node, $width (D=68), $min_dx (D=2), $dy (D=1) )  
3506  #===============================================================================  #===============================================================================
3507  sub text_plot_newick {  sub text_plot_newick
3508      my ( $node, $width, $min_dx, $dy ) = @_;  {
3509        my $node = shift @_;
3510      array_ref( $node ) || die "Bad node passed to text_plot_newick\n";      array_ref( $node ) || die "Bad node passed to text_plot_newick\n";
3511      defined( $min_dx ) and ( $min_dx >=  0 ) or $min_dx =  2;  
3512      defined(     $dy ) and (     $dy >=  1 ) or     $dy =  1;      my ( $opts, $width, $min_dx, $dy, $chars, $fmt );
3513      defined( $width  )                       or  $width = 68;      if ( $_[0] && ref $_[0] eq 'HASH' )
3514        {
3515            $opts = shift;
3516        }
3517        else
3518        {
3519            ( $width, $min_dx, $dy ) = @_;
3520            $opts = {};
3521        }
3522    
3523        $chars = $opts->{ chars } || '';
3524        my $charH;
3525        $charH = $char_set{ $chars } || $char_set{ 'text1' } if ( $chars ne 'raw' );
3526        my $is_box = $charH eq $char_set{ html_box }
3527                  || $charH eq $char_set{ utf8_box }
3528                  || $chars eq 'raw';
3529    
3530        $fmt = ( $chars eq 'raw' ) ? 'chrlist_lbl' : $opts->{ format };
3531        $fmt = $tree_format{ $fmt || '' } || 'text';
3532    
3533        $dy    ||= $opts->{ dy     } ||  1;
3534        $width ||= $opts->{ width  } || 68;
3535        $min_dx  = $opts->{ min_dx } if ( ! defined $min_dx || $min_dx < 0 );
3536        $min_dx  = $is_box ? 1 : 2   if ( ! defined $min_dx || $min_dx < 0 );
3537    
3538        #  Layout the tree:
3539    
3540      $min_dx = int( $min_dx );      $min_dx = int( $min_dx );
3541      $dy     = int( $dy );      $dy     = int( $dy );
# Line 2835  Line 3544 
3544      my $hash = {};      my $hash = {};
3545      layout_printer_plot( $node, $hash, 0, -0.5 * $dy, $x_scale, $min_dx, $dy );      layout_printer_plot( $node, $hash, 0, -0.5 * $dy, $x_scale, $min_dx, $dy );
3546    
3547      # dump_tree_hash( $node, $hash ); exit;      #  Generate the lines of the tree-one by-one:
   
     #  Generate the lines of the tree one by one:  
3548    
3549      my ( $y1, $y2 ) = @{ $hash->{ $node } };      my ( $y1, $y2 ) = @{ $hash->{ $node } };
3550      map { text_tree_row( $node, $hash, $_, "", "+" ) } ( $y1 .. $y2 );      my @lines;
3551        foreach ( ( $y1 .. $y2 ) )
3552        {
3553            my $line = text_tree_row( $node, $hash, $_, [], 'tee_l' );
3554            my $lbl  = '';
3555            if ( @$line )
3556            {
3557                if ( $line->[-1] eq '' ) { pop @$line; $lbl = pop @$line }
3558                #  Translate tree characters
3559                @$line = map { $charH->{ $_ } } @$line if $chars ne 'raw';
3560            }
3561    
3562            # Convert to requested output format:
3563    
3564            push @lines, $fmt eq 'text'         ? join( '', @$line, ( $lbl ? " $lbl" : () ) )
3565                       : $fmt eq 'text_tab_lbl' ? join( '', @$line, "\t", $lbl )
3566                       : $fmt eq 'tree_lbl'     ? [ join( '', @$line ), $lbl ]
3567                       : $fmt eq 'chrlist_lbl'  ? [ $line, $lbl ]
3568                       :                          ();
3569  }  }
3570    
3571        # if ( $cells )
3572        # {
3573        #     my $nmax = 0;
3574        #     foreach ( @lines ) { $nmax = @$_ if @$_ > $nmax }
3575        #     foreach ( @lines )
3576        #     {
3577        #         @$_ = map { "<TD>$_</TD>" } @$_;
3578        #         my $span = $nmax - @$_ + 1;
3579        #         $_->[-1] =~ s/^<TD>/<TD NoWrap ColSpan=$span>/;
3580        #     }
3581        # }
3582        # elsif ( $tables )
3583        # {
3584        #     my $nmax = 0;
3585        #     foreach ( @lines ) { $nmax = @$_ if @$_ > $nmax }
3586        #     foreach ( @lines )
3587        #     {
3588        #         @$_ = map { "<TD>$_</TD>" } @$_;
3589        #         my $span = $nmax - @$_ + 1;
3590        #         $_->[-1] =~ s/^<TD>/<TD NoWrap ColSpan=$span>/;
3591        #     }
3592        # }
3593    
3594        wantarray ? @lines : \@lines;
3595    }
3596    
3597    
3598  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3599  #  ( $xmax, $ymax, $root_y ) = layout_printer_plot( $node, $hash, $x0, $y0, $x_scale, $min_dx, $dy )  #  ( $xmax, $ymax, $root_y ) = layout_printer_plot( $node, $hash, $x0, $y0, $x_scale, $min_dx, $dy, $yrnd )
3600  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3601  sub layout_printer_plot {  sub layout_printer_plot
3602      my ( $node, $hash, $x0, $y0, $x_scale, $min_dx, $dy ) = @_;  {
3603        my ( $node, $hash, $x0, $y0, $x_scale, $min_dx, $dy, $yrnd ) = @_;
3604      array_ref( $node ) || die "Bad node ref passed to layout_printer_plot\n";      array_ref( $node ) || die "Bad node ref passed to layout_printer_plot\n";
3605      hash_ref(  $hash ) || die "Bad hash ref passed to layout_printer_plot\n";      hash_ref(  $hash ) || die "Bad hash ref passed to layout_printer_plot\n";
3606    
3607      my $dx = newick_x( $node );      my $dx = newick_x( $node );
3608      if ( defined( $dx ) ) {      if ( defined( $dx ) ) {
3609          $dx *= $x_scale;          $dx *= $x_scale;
3610          $dx >= $min_dx or $dx = $min_dx;          $dx = $min_dx if $dx < $min_dx;
3611      }      }
3612      else {      else {
3613          $dx = ( $x0 > 0 ) ? $min_dx : 0;          $dx = ( $x0 > 0 ) ? $min_dx : 0;
# Line 2881  Line 3634 
3634          $ymax = $y0;          $ymax = $y0;
3635    
3636          foreach ( @dl ) {          foreach ( @dl ) {
3637              ( $xmaxi, $ymax, $yi ) = layout_printer_plot( $_, $hash, $x, $ymax, $x_scale, $min_dx, $dy );              ( $xmaxi, $ymax, $yi ) = layout_printer_plot( $_, $hash, $x, $ymax, $x_scale, $min_dx, $dy,
3638                                                              ( 2*@ylist < @dl ? 0.5001 : 0.4999 )
3639                                                            );
3640              push @ylist, $yi;              push @ylist, $yi;
3641              if ( $xmaxi > $xmax ) { $xmax = $xmaxi }              if ( $xmaxi > $xmax ) { $xmax = $xmaxi }
3642          }          }
# Line 2891  Line 3646 
3646    
3647          $yn1 = $ylist[ 0];          $yn1 = $ylist[ 0];
3648          $yn2 = $ylist[-1];          $yn2 = $ylist[-1];
3649          $y = int( 0.5 * ( $yn1 + $yn2 ) + 0.4999 );          $y = int( 0.5 * ( $yn1 + $yn2 ) + ( $yrnd || 0.4999 ) );
3650      }      }
3651    
3652      $y2 = int( $ymax - 0.5 * $dy + 0.4999 );      $y2 = int( $ymax - 0.5 * $dy + 0.4999 );
# Line 2901  Line 3656 
3656  }  }
3657    
3658    
3659  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #  What symbol do we get if we add a leftward line to some other symbol?
 #  Debug routine  
 #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  
 sub dump_tree {  
     my ( $node, $prefix ) = @_;  
     defined( $prefix ) or $prefix = "";  
     print STDERR $prefix, join(", ", @$node), "\n";  
     my @dl = $node->[0] ? @{$node->[0]} : ();  
     foreach ( @dl ) { dump_tree( $_, $prefix . "  " ) }  
     $prefix or print STDERR "\n";  
 }  
   
   
 #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  
 #  Debug routine  
 #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  
 sub dump_tree_hash {  
     my ( $node, $hash, $prefix ) = @_;  
     defined( $prefix ) or print STDERR "node; [ y1, y2, x0, x, y, yn1, yn2 ]\n" and $prefix = "";  
     print STDERR $prefix, join(", ", @$node), "; ", join(", ", @{ $hash->{ $node } } ), "\n";  
     my @dl = $node->[0] ? @{$node->[0]} : ();  
     foreach (@dl) { dump_tree_hash( $_, $hash, $prefix . "  " ) }  
 }  
3660    
3661    my %with_left_line = ( space  => 'half_l',
3662                           horiz  => 'horiz',
3663                           vert   => 'tee_l',
3664                           el_d_r => 'tee_d',
3665                           el_u_r => 'tee_u',
3666                           el_d_l => 'el_d_l',
3667                           el_u_l => 'el_u_l',
3668                           tee_l  => 'tee_l',
3669                           tee_r  => 'cross',
3670                           tee_u  => 'tee_u',
3671                           tee_d  => 'tee_d',
3672                           half_l => 'half_l',
3673                           half_r => 'horiz',
3674                           half_u => 'el_u_l',
3675                           half_d => 'el_d_l',
3676                           cross  => 'cross',
3677                         );
3678    
3679  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3680  #  $line = text_tree_row( $node, $hash, $row, $line, $symb )  #  Produce a description of one line of a printer plot tree.
3681    #
3682    #  \@line = text_tree_row( $node, $hash, $row, \@line, $symb )
3683    #
3684    #     \@line is the character descriptions accumulated so far, one per array
3685    #          element, except for a label, which can be any number of characters.
3686    #          Labels are followed by an empty string, so if $line->[-1] eq '',
3687    #          then $line->[-2] is a label. The calling program translates the
3688    #          symbol names to output characters.
3689    #
3690    #     \@node is a newick tree node
3691    #     \%hash contains tree layout information
3692    #      $row  is the row number (y value) that we are building
3693    #      $symb is the plot symbol proposed for the current x and y position
3694    #
3695  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3696  sub text_tree_row {  sub text_tree_row
3697    {
3698      my ( $node, $hash, $row, $line, $symb ) = @_;      my ( $node, $hash, $row, $line, $symb ) = @_;
3699    
3700      my ( $y1, $y2, $x0, $x, $y, $yn1, $yn2 ) = @{ $hash->{ $node } };      my ( $y1, $y2, $x0, $x, $y, $yn1, $yn2 ) = @{ $hash->{ $node } };
3701      if ( $row < $y1 || $row > $y2 ) { return $line }      if ( $row < $y1 || $row > $y2 ) { return $line }
3702    
3703      if ( length( $line ) < $x0 ) { $line .= " " x ( $x0 - length( $line ) ) }      if ( @$line < $x0 ) { push @$line, ('space') x ( $x0 - @$line ) }
3704    
3705      if ( $row == $y ) {      if ( $row == $y ) {
3706          $line = substr( $line, 0, $x0 ) . $symb . (( $x > $x0 ) ? "-" x ($x - $x0) : "");          while ( @$line > $x0 ) { pop @$line }  # Actually 0-1 times
3707            push @$line, $symb,
3708                         ( ( $x > $x0 ) ? ('horiz') x ($x - $x0) : () );
3709      }      }
3710    
3711      elsif ( $row > $yn1 && $row < $yn2 ) {      elsif ( $row > $yn1 && $row < $yn2 ) {
3712          if ( length( $line ) < $x ) { $line .= " " x ( $x - length( $line ) ) . "|" }          if ( @$line < $x ) { push @$line, ('space') x ( $x - @$line ), 'vert' }
3713          else { substr( $line, $x ) = "|" }          else               { $line->[$x] = 'vert' }
3714      }      }
3715    
3716      my @dl = newick_desc_list( $node );      my @dl = newick_desc_list( $node );
3717    
3718      if ( @dl < 1 ) {      if ( @dl < 1 ) { push @$line, ( newick_lbl( $node ) || '' ), '' }
         $line .= " " . $node->[1];  
     }  
3719    
3720      else {      else {
3721          my @list = map { [ $_, "+" ] } @dl;  #  Print symbol for line          my @list = map { [ $_, 'tee_r' ] } @dl;  # Line to the right
3722          $list[ 0]->[1] = "/";          if ( @list > 1 ) { #  Fix top and bottom sympbols
3723          $list[-1]->[1] = "\\";              $list[ 0]->[1] = 'el_d_r';
3724                $list[-1]->[1] = 'el_u_r';
3725            }
3726            elsif ( @list ) {  # Only one descendent
3727                $list[ 0]->[1] = 'half_r';
3728            }
3729          foreach ( @list ) {          foreach ( @list ) {
3730              my ( $n, $s ) = @$_;              my ( $n, $s ) = @$_;
3731              if ( $row >= $hash->{ $n }->[0] && $row <= $hash->{ $n }->[1] ) {              if ( $row >= $hash->{ $n }->[0] && $row <= $hash->{ $n }->[1] ) {
# Line 2964  Line 3733 
3733              }              }
3734           }           }
3735    
3736          if ( $row == $y ) { substr( $line, $x, 1 ) = "+" }          if ( $row == $y ) {
3737                $line->[$x] = ( $line->[$x] eq 'horiz' ) ? 'tee_l'
3738                                                         : $with_left_line{ $line->[$x] };
3739            }
3740      }      }
3741    
3742      return $line;      return $line;
3743  }  }
3744    
3745    
3746    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3747    #  Debug routine
3748    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3749    sub dump_tree {
3750        my ( $node, $prefix ) = @_;
3751        defined( $prefix ) or $prefix = "";
3752        print STDERR $prefix, join(", ", @$node), "\n";
3753        my @dl = $node->[0] ? @{$node->[0]} : ();
3754        foreach ( @dl ) { dump_tree( $_, $prefix . "  " ) }
3755        $prefix or print STDERR "\n";
3756    }
3757    
3758    
3759    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3760    #  Debug routine
3761    #- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3762    sub dump_tree_hash {
3763        my ( $node, $hash, $prefix ) = @_;
3764        defined( $prefix ) or print STDERR "node; [ y1, y2, x0, x, y, yn1, yn2 ]\n" and $prefix = "";
3765        print STDERR $prefix, join(", ", @$node), "; ", join(", ", @{ $hash->{ $node } } ), "\n";
3766        my @dl = $node->[0] ? @{$node->[0]} : ();
3767        foreach (@dl) { dump_tree_hash( $_, $hash, $prefix . "  " ) }
3768    }
3769    
3770    
3771  #===============================================================================  #===============================================================================
3772  #  Open an input file stream:  #  Open an input file stream:
3773  #  #
# Line 3012  Line 3809 
3809      return undef;      return undef;
3810  }  }
3811    
3812    
3813    #===============================================================================
3814    #  Some subroutines copied from gjolists
3815    #===============================================================================
3816    #  Return the common prefix of two lists:
3817    #
3818    #  @common = common_prefix( \@list1, \@list2 )
3819    #-----------------------------------------------------------------------------
3820    sub common_prefix
3821    {
3822        my ($l1, $l2) = @_;
3823        ref($l1) eq "ARRAY" || die "common_prefix: arg 1 is not an array ref\n";
3824        ref($l2) eq "ARRAY" || die "common_prefix: arg 2 is not an array ref\n";
3825    
3826        my $i = 0;
3827        my $l1_i;
3828        while ( defined( $l1_i = $l1->[$i] ) && $l1_i eq $l2->[$i] ) { $i++ }
3829    
3830        return @$l1[ 0 .. ($i-1) ];  # perl handles negative range
3831    }
3832    
3833    
3834    #-----------------------------------------------------------------------------
3835    #  Return the unique suffixes of each of two lists:
3836    #
3837    #  ( \@suffix1, \@suffix2 ) = unique_suffixes( \@list1, \@list2 )
3838    #-----------------------------------------------------------------------------
3839    sub unique_suffixes
3840    {
3841        my ($l1, $l2) = @_;
3842        ref($l1) eq "ARRAY" || die "common_prefix: arg 1 is not an array ref\n";
3843        ref($l2) eq "ARRAY" || die "common_prefix: arg 2 is not an array ref\n";
3844    
3845        my $i = 0;
3846        my @l1 = @$l1;
3847        my @l2 = @$l2;
3848        my $l1_i;
3849        while ( defined( $l1_i = $l1[$i] ) && $l1_i eq $l2[$i] ) { $i++ }
3850    
3851        splice @l1, 0, $i;
3852        splice @l2, 0, $i;
3853        return ( \@l1, \@l2 );
3854    }
3855    
3856    
3857    #-------------------------------------------------------------------------------
3858    #  List of values duplicated in a list (stable in order by second occurance):
3859    #
3860    #  @dups = duplicates( @list )
3861    #-------------------------------------------------------------------------------
3862    sub duplicates
3863    {
3864        my %cnt = ();
3865        grep { ++$cnt{$_} == 2 } @_;
3866    }
3867    
3868    
3869    #-------------------------------------------------------------------------------
3870    #  Randomize the order of a list:
3871    #
3872    #  @random = random_order( @list )
3873    #-------------------------------------------------------------------------------
3874    sub random_order
3875    {
3876        my ( $i, $j );
3877        for ( $i = @_ - 1; $i > 0; $i-- )
3878        {
3879            $j = int( ($i+1) * rand() );
3880            ( $_[$i], $_[$j] ) = ( $_[$j], $_[$i] ); # Interchange i and j
3881        }
3882    
3883       @_;
3884    }
3885    
3886    
3887    #-----------------------------------------------------------------------------
3888    #  Intersection of two or more sets:
3889    #
3890    #  @intersection = intersection( \@set1, \@set2, ... )
3891    #-----------------------------------------------------------------------------
3892    sub intersection
3893    {
3894        my $set = shift;
3895        my @intersection = @$set;
3896    
3897        foreach $set ( @_ )
3898        {
3899            my %set = map { $_ => 1 } @$set;
3900            @intersection = grep { exists $set{ $_ } } @intersection;
3901        }
3902    
3903        @intersection;
3904    }
3905    
3906    
3907    #-----------------------------------------------------------------------------
3908    #  Elements in set 1, but not set 2:
3909    #
3910    #  @difference = set_difference( \@set1, \@set2 )
3911    #-----------------------------------------------------------------------------
3912    sub set_difference
3913    {
3914        my ($set1, $set2) = @_;
3915        my %set2 = map { $_ => 1 } @$set2;
3916        grep { ! ( exists $set2{$_} ) } @$set1;
3917    }
3918    
3919    
3920  1;  1;

Legend:
Removed from v.1.8  
changed lines
  Added in v.1.14

MCS Webmaster
ViewVC Help
Powered by ViewVC 1.0.3